Fixes for Intel and container build fixes
authorPaul Bauer <paul.bauer.q@gmail.com>
Mon, 23 Nov 2020 20:51:28 +0000 (20:51 +0000)
committerArtem Zhmurov <zhmurov@gmail.com>
Mon, 23 Nov 2020 20:51:28 +0000 (20:51 +0000)
Generalizes the TSAN build stage to a compiler build stage, and uses
it for containing the installation process for the free Intel tools

Fixed or worked around several bugs

Updated install guide section on tested platforms to describe
what we actually do

Refs #3459, #3620

Change-Id: I28e4b816923395edead9d425cec5d85581e22b5e

admin/containers/scripted_gmx_docker_builds.py
admin/containers/utility.py
docs/install-guide/index.rst
src/gromacs/gmxana/nsfactor.cpp
src/gromacs/hardware/detecthardware.cpp
src/gromacs/mdrun/md.cpp
src/gromacs/mdrun/mimic.cpp
src/gromacs/mdrun/rerun.cpp
src/gromacs/selection/CMakeLists.txt
src/gromacs/utility/basedefinitions.h

index 1efd05b5dca874499c26c136c0d9924f308d77a9..26a20c720a16eb6eb74a0f9aac51945e98ed9672 100644 (file)
@@ -48,6 +48,7 @@ Authors:
     * Paul Bauer <paul.bauer.q@gmail.com>
     * Eric Irrgang <ericirrgang@gmail.com>
     * Joe Jordan <e.jjordan12@gmail.com>
+    * Mark Abraham <mark.j.abraham@gmail.com>
 
 Usage::
 
@@ -179,7 +180,7 @@ def get_llvm_packages(args) -> typing.Iterable[str]:
         return []
 
 
-def get_compiler(args, tsan_stage: hpccm.Stage = None) -> bb_base:
+def get_compiler(args, compiler_build_stage: hpccm.Stage = None) -> bb_base:
     # Compiler
     if args.icc is not None:
         raise RuntimeError('Intel compiler toolchain recipe not implemented yet')
@@ -187,15 +188,15 @@ def get_compiler(args, tsan_stage: hpccm.Stage = None) -> bb_base:
     if args.llvm is not None:
         # Build our own version instead to get TSAN + OMP
         if args.tsan is not None:
-            if tsan_stage is not None:
-                compiler = tsan_stage.runtime(_from='tsan')
+            if compiler_build_stage is not None:
+                compiler = compiler_build_stage.runtime(_from='tsan')
             else:
-                raise RuntimeError('No TSAN stage!')
+                raise RuntimeError('No TSAN compiler build stage!')
         # Build the default compiler if we don't need special support
         else:
             compiler = hpccm.building_blocks.llvm(extra_repository=True, version=args.llvm)
 
-    elif (args.gcc is not None):
+    elif args.gcc is not None:
         compiler = hpccm.building_blocks.gnu(extra_repository=True,
                                              version=args.gcc,
                                              fortran=False)
@@ -218,6 +219,8 @@ def get_mpi(args, compiler):
                 raise RuntimeError('compiler is not an HPCCM compiler building block!')
 
         elif args.mpi == 'impi':
+            # TODO also consider hpccm's intel_mpi package if that doesn't need
+            # a license to run.
             raise RuntimeError('Intel MPI recipe not implemented yet.')
         else:
             raise RuntimeError('Requested unknown MPI implementation.')
@@ -260,7 +263,7 @@ def get_clfft(args):
         return None
 
 
-def add_tsan_stage(input_args, output_stages: typing.Mapping[str, hpccm.Stage]):
+def add_tsan_compiler_build_stage(input_args, output_stages: typing.Mapping[str, hpccm.Stage]):
     """Isolate the expensive TSAN preparation stage.
 
     This is a very expensive stage, but has few and disjoint dependencies, and
@@ -295,8 +298,7 @@ def add_tsan_stage(input_args, output_stages: typing.Mapping[str, hpccm.Stage]):
                      'ln -s /usr/local/bin/clang-format /usr/local/bin/clang-format-' + str(input_args.llvm),
                      'ln -s /usr/local/bin/clang-tidy /usr/local/bin/clang-tidy-' + str(input_args.llvm),
                      'ln -s /usr/local/libexec/c++-analyzer /usr/local/bin/c++-analyzer-' + str(input_args.llvm)])
-    output_stages['tsan'] = tsan_stage
-
+    output_stages['compiler_build'] = tsan_stage
 
 def prepare_venv(version: StrictVersion) -> typing.Sequence[str]:
     """Get shell commands to set up the venv for the requested Python version."""
@@ -415,16 +417,18 @@ def build_stages(args) -> typing.Iterable[hpccm.Stage]:
     # object early in this function.
     stages = collections.OrderedDict()
 
-    # If we need the TSAN compilers, the early build is more involved.
+    # If we need TSAN or oneAPI support the early build is more complex,
+    # so that our compiler images don't have all the cruft needed to get those things
+    # installed.
     if args.llvm is not None and args.tsan is not None:
-        add_tsan_stage(input_args=args, output_stages=stages)
+        add_tsan_compiler_build_stage(input_args=args, output_stages=stages)
 
     # Building blocks are chunks of container-builder instructions that can be
     # copied to any build stage with the addition operator.
     building_blocks = collections.OrderedDict()
 
     # These are the most expensive and most reusable layers, so we put them first.
-    building_blocks['compiler'] = get_compiler(args, tsan_stage=stages.get('tsan'))
+    building_blocks['compiler'] = get_compiler(args, compiler_build_stage=stages.get('compiler_build'))
     building_blocks['mpi'] = get_mpi(args, building_blocks['compiler'])
 
     # Install additional packages early in the build to optimize Docker build layer cache.
index e413a19be0c7f1fb661ddf9cd54f33dfd8b62e7b..c79de2a6085508d41d08b218cf057dce232af513 100644 (file)
@@ -40,6 +40,7 @@ Authors:
     * Paul Bauer <paul.bauer.q@gmail.com>
     * Eric Irrgang <ericirrgang@gmail.com>
     * Joe Jordan <e.jjordan12@gmail.com>
+    * Mark Abraham <mark.j.abraham@gmail.com>
 
 """
 
index c0e7238c99b5e059bb5068b2c4a1f6871cf810d9..6f332dff78f33cc2670539335a7a2f9ed7c18d2b 100644 (file)
@@ -1289,17 +1289,18 @@ Tested platforms
 
 While it is our best belief that |Gromacs| will build and run pretty
 much everywhere, it is important that we tell you where we really know
-it works because we have tested it. We do test on Linux, Windows, and
-Mac with a range of compilers and libraries for a range of our
-configuration options. Every commit in our git source code repository
-is currently tested on x86 with a number of gcc versions ranging from 5.1
-through 9.1, version 19 of the Intel compiler, and Clang
-versions 3.6 through 8. For this, we use a variety of GNU/Linux
-flavors and versions as well as Windows (where we test only MSVC 2017).
+it works because we have tested it.
+Every commit in our git source code repository
+is currently tested with a range of configuration options on x86 with
+gcc versions 6 and 7,
+clang versions 3.6 and 8,
+and
+For this testing, we use Ubuntu 16.04 or 18.04 operating system.
 Other compiler, library, and OS versions are tested less frequently.
 For details, you can
 have a look at the `continuous integration server used by GROMACS`_,
-which runs Jenkins_.
+which uses GitLab runner on a local k8s x86 cluster with NVIDIA and
+AMD GPU support.
 
 We test irregularly on ARM v7, ARM v8, Cray, Fujitsu
 PRIMEHPC, Power8, Power9,
index 3c53901f132c91a27c4a9246220eeaf73c3abc11..e56b568d756502b9ce707536d2b2b89d6a50ff6a 100644 (file)
@@ -1,7 +1,8 @@
 /*
  * This file is part of the GROMACS molecular simulation package.
  *
- * Copyright (c) 2012,2013,2014,2015,2016,2017,2018,2019, by the GROMACS development team, led by
+ * Copyright (c) 2012,2013,2014,2015,2016, The GROMACS development team.
+ * Copyright (c) 2017,2018,2019,2020, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -45,6 +46,7 @@
 #include "gromacs/random/threefry.h"
 #include "gromacs/random/uniformintdistribution.h"
 #include "gromacs/topology/topology.h"
+#include "gromacs/utility/basedefinitions.h"
 #include "gromacs/utility/cstringutil.h"
 #include "gromacs/utility/exceptions.h"
 #include "gromacs/utility/fatalerror.h"
@@ -208,7 +210,7 @@ gmx_radial_distribution_histogram_t* calc_radial_distribution_histogram(gmx_sans
     int                       nthreads;
     gmx::DefaultRandomEngine* trng = nullptr;
 #endif
-    int64_t                  mc = 0, mc_max;
+    int64_t                  mc_max;
     gmx::DefaultRandomEngine rng(seed);
 
     /* allocate memory for pr */
@@ -249,13 +251,14 @@ gmx_radial_distribution_histogram_t* calc_radial_distribution_histogram(gmx_sans
             snew(tgr[i], pr->grn);
             trng[i].seed(rng());
         }
-#    pragma omp parallel shared(tgr, trng, mc) private(tid, i, j)
+#    pragma omp parallel shared(tgr, trng) private(tid, i, j)
         {
             gmx::UniformIntDistribution<int> tdist(0, isize - 1);
             tid = gmx_omp_get_thread_num();
-/* now starting parallel threads */
+            /* now starting parallel threads */
+            INTEL_DIAGNOSTIC_IGNORE(593)
 #    pragma omp for
-            for (mc = 0; mc < mc_max; mc++)
+            for (int64_t mc = 0; mc < mc_max; mc++)
             {
                 try
                 {
@@ -270,6 +273,7 @@ gmx_radial_distribution_histogram_t* calc_radial_distribution_histogram(gmx_sans
                 GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
             }
         }
+        INTEL_DIAGNOSTIC_RESET
         /* collecting data from threads */
         for (i = 0; i < pr->grn; i++)
         {
@@ -287,7 +291,7 @@ gmx_radial_distribution_histogram_t* calc_radial_distribution_histogram(gmx_sans
         delete[] trng;
 #else
         gmx::UniformIntDistribution<int> dist(0, isize - 1);
-        for (mc = 0; mc < mc_max; mc++)
+        for (int64_t mc = 0; mc < mc_max; mc++)
         {
             i = dist(rng); // [0,isize-1]
             j = dist(rng); // [0,isize-1]
index c24a8b2a4cec6c4d0f9c761d97d368c952981e01..e5c754e78f20b7d69bd7aa58e861d69246661fd0 100644 (file)
@@ -1,7 +1,8 @@
 /*
  * This file is part of the GROMACS molecular simulation package.
  *
- * Copyright (c) 2012,2013,2014,2015,2016,2017,2018,2019, by the GROMACS development team, led by
+ * Copyright (c) 2012,2013,2014,2015,2016, The GROMACS development team.
+ * Copyright (c) 2017,2018,2019,2020, by the GROMACS development team, led by
  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
  * and including many others, as listed in the AUTHORS file in the
  * top-level source directory and at http://www.gromacs.org.
@@ -300,9 +301,17 @@ static void gmx_collect_hardware_mpi(const gmx::CpuInfo&              cpuInfo,
  *  2 seconds, or until all cores have come online. This can be used prior to
  *  hardware detection for platforms that take unused processors offline.
  *
- *  This routine will not throw exceptions.
+ *  This routine will not throw exceptions. In principle it should be
+ *  declared noexcept, but at least icc 19.1 and 21-beta08 with the
+ *  libstdc++-7.5 has difficulty implementing a std::vector of
+ *  std::thread started with this function when declared noexcept. It
+ *  is a known compiler bug that should be fixed after 19.1.
+ *  Fortunately, this function is not performance sensitive,
+ *  and only runs on platforms other than x86 and POWER (ie ARM),
+ *  so the possible overhead introduced by omitting noexcept is not
+ *  important.
  */
-static void spinUpCore() noexcept
+static void spinUpCore()
 {
 #if defined(HAVE_SYSCONF) && defined(_SC_NPROCESSORS_CONF) && defined(_SC_NPROCESSORS_ONLN)
     float dummy           = 0.1;
index 5ca5064229efa6182a30e4e2782847d41c16de22..9e519d44daffc17591d92f82307a15ab93159128 100644 (file)
@@ -159,7 +159,7 @@ void gmx::LegacySimulator::do_md()
     int64_t      step, step_rel;
     double       t, t0 = ir->init_t, lam0[efptNR];
     gmx_bool     bGStatEveryStep, bGStat, bCalcVir, bCalcEnerStep, bCalcEner;
-    gmx_bool     bNS, bNStList, bStopCM, bFirstStep, bInitStep, bLastStep = FALSE;
+    gmx_bool     bNS = FALSE, bNStList, bStopCM, bFirstStep, bInitStep, bLastStep = FALSE;
     gmx_bool     bDoDHDL = FALSE, bDoFEP = FALSE, bDoExpanded = FALSE;
     gmx_bool     do_ene, do_log, do_verbose;
     gmx_bool     bMasterState;
@@ -677,6 +677,9 @@ void gmx::LegacySimulator::do_md()
     bExchanged       = FALSE;
     bNeedRepartition = FALSE;
 
+    step     = ir->init_step;
+    step_rel = 0;
+
     auto stopHandler = stopHandlerBuilder->getStopHandlerMD(
             compat::not_null<SimulationSignal*>(&signals[eglsSTOPCOND]), simulationsShareState,
             MASTER(cr), ir->nstlist, mdrunOptions.reproducible, nstSignalComm,
@@ -695,9 +698,6 @@ void gmx::LegacySimulator::do_md()
 
     const DDBalanceRegionHandler ddBalanceRegionHandler(cr);
 
-    step     = ir->init_step;
-    step_rel = 0;
-
     // TODO extract this to new multi-simulation module
     if (MASTER(cr) && isMultiSim(ms) && !useReplicaExchange)
     {
index 0894cfb76168667f9ce36e676aed372a3b207cd2..ebdbd3031c2dc85b4fe1d93380033ebad62c2ef7 100644 (file)
@@ -347,6 +347,9 @@ void gmx::LegacySimulator::do_mimic()
                     "MiMiC does not report kinetic energy, total energy, temperature, virial and "
                     "pressure.");
 
+    step     = ir->init_step;
+    step_rel = 0;
+
     auto stopHandler = stopHandlerBuilder->getStopHandlerMD(
             compat::not_null<SimulationSignal*>(&signals[eglsSTOPCOND]), false, MASTER(cr),
             ir->nstlist, mdrunOptions.reproducible, nstglobalcomm, mdrunOptions.maximumHoursToRun,
@@ -357,9 +360,6 @@ void gmx::LegacySimulator::do_mimic()
 
     const DDBalanceRegionHandler ddBalanceRegionHandler(cr);
 
-    step     = ir->init_step;
-    step_rel = 0;
-
     /* and stop now if we should */
     isLastStep = (isLastStep || (ir->nsteps >= 0 && step_rel > ir->nsteps));
     while (!isLastStep)
index 956bd1ab1cbdd33e3ca7e8f673e5ebe4586607b9..310c8516b176a4eac89dd8e383486482cf26d3a3 100644 (file)
@@ -467,6 +467,9 @@ void gmx::LegacySimulator::do_rerun()
         calc_shifts(rerun_fr.box, fr->shift_vec);
     }
 
+    step     = ir->init_step;
+    step_rel = 0;
+
     auto stopHandler = stopHandlerBuilder->getStopHandlerMD(
             compat::not_null<SimulationSignal*>(&signals[eglsSTOPCOND]), false, MASTER(cr),
             ir->nstlist, mdrunOptions.reproducible, nstglobalcomm, mdrunOptions.maximumHoursToRun,
@@ -477,9 +480,6 @@ void gmx::LegacySimulator::do_rerun()
 
     const DDBalanceRegionHandler ddBalanceRegionHandler(cr);
 
-    step     = ir->init_step;
-    step_rel = 0;
-
     /* and stop now if we should */
     isLastStep = (isLastStep || (ir->nsteps >= 0 && step_rel > ir->nsteps));
     while (!isLastStep)
index 9c83ba1dfaa8f717b4a8dbcc7745361f89229771..43f8530ceccfe31acbc877a8057bc5ebd98543e7 100644 (file)
@@ -55,7 +55,6 @@ else()
     gmx_target_warning_suppression(scanner -Wno-unused-parameter HAS_NO_UNUSED_PARAMETER)
     gmx_target_warning_suppression(scanner -Wno-missing-declarations HAS_NO_MISSING_DECLARATIONS)
     gmx_target_warning_suppression(scanner -Wno-null-conversion HAS_NO_NULL_CONVERSIONS)
-    gmx_target_warning_suppression(scanner -wd1419 HAS_DECL_IN_SOURCE)
 endif()
 list(APPEND libgromacs_object_library_dependencies scanner)
 set(libgromacs_object_library_dependencies ${libgromacs_object_library_dependencies} PARENT_SCOPE)
index d1ae0c6ce71f26f60b9398f866daf8008be71ff3..939a7a27f89c104d338029c576e853f1a4ddae32 100644 (file)
@@ -169,6 +169,16 @@ index ssize(const T& t)
 #    define MSVC_DIAGNOSTIC_RESET
 #endif
 
+#ifdef __INTEL_COMPILER
+// Ignore unused loop variable warning - it was used until the compiler removes the use!
+#    define DO_PRAGMA(x) _Pragma(#    x)
+#    define INTEL_DIAGNOSTIC_IGNORE(id) DO_PRAGMA(warning push) DO_PRAGMA(warning(disable : id))
+#    define INTEL_DIAGNOSTIC_RESET DO_PRAGMA(warning pop)
+#else
+#    define INTEL_DIAGNOSTIC_IGNORE(id)
+#    define INTEL_DIAGNOSTIC_RESET
+#endif
+
 namespace gmx
 {
 namespace internal