2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2015,2016,2017, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
36 * \brief Defines utility functionality for dividing resources and
37 * checking for consistency and usefulness.
39 * \author Mark Abraham <mark.j.abraham@gmail.com>
40 * \ingroup module_taskassignment
45 #include "resourcedivision.h"
54 #include "gromacs/ewald/pme.h"
55 #include "gromacs/hardware/cpuinfo.h"
56 #include "gromacs/hardware/detecthardware.h"
57 #include "gromacs/hardware/hardwaretopology.h"
58 #include "gromacs/hardware/hw_info.h"
59 #include "gromacs/mdlib/gmx_omp_nthreads.h"
60 #include "gromacs/mdtypes/commrec.h"
61 #include "gromacs/mdtypes/inputrec.h"
62 #include "gromacs/mdtypes/md_enums.h"
63 #include "gromacs/topology/mtop_util.h"
64 #include "gromacs/topology/topology.h"
65 #include "gromacs/utility/baseversion.h"
66 #include "gromacs/utility/fatalerror.h"
67 #include "gromacs/utility/gmxassert.h"
68 #include "gromacs/utility/logger.h"
69 #include "gromacs/utility/stringutil.h"
72 /* DISCLAIMER: All the atom count and thread numbers below are heuristic.
73 * The real switching points will depend on the system simulation,
74 * the algorithms used and the hardware it's running on, as well as if there
75 * are other jobs running on the same machine. We try to take into account
76 * factors that have a large influence, such as recent Intel CPUs being
77 * much better at wide multi-threading. The remaining factors should
78 * (hopefully) have a small influence, such that the performance just before
79 * and after a switch point doesn't change too much.
82 //! Constant used to help minimize preprocessed code
83 static const bool bHasOmpSupport
= GMX_OPENMP
;
85 /*! \brief The minimum number of atoms per thread-MPI thread when GPUs
86 * are present. With fewer atoms than this, the number of thread-MPI
87 * ranks will get lowered.
89 static const int min_atoms_per_mpi_thread
= 90;
90 /*! \brief The minimum number of atoms per GPU with thread-MPI
91 * active. With fewer atoms than this, the number of thread-MPI ranks
94 static const int min_atoms_per_gpu
= 900;
97 /*! \brief Constants for implementing default divisions of threads */
99 /* TODO choose nthreads_omp based on hardware topology
100 when we have a hardware topology detection library */
101 /* First we consider the case of no MPI (1 MPI rank).
102 * In general, when running up to 8 threads, OpenMP should be faster.
103 * Note: on AMD Bulldozer we should avoid running OpenMP over two dies.
104 * On Intel>=Nehalem running OpenMP on a single CPU is always faster,
105 * even on two CPUs it's usually faster (but with many OpenMP threads
106 * it could be faster not to use HT, currently we always use HT).
107 * On Nehalem/Westmere we want to avoid running 16 threads over
108 * two CPUs with HT, so we need a limit<16; thus we use 12.
109 * A reasonable limit for Intel Sandy and Ivy bridge,
110 * not knowing the topology, is 16 threads.
111 * Below we check for Intel and AVX, which for now includes
112 * Sandy/Ivy Bridge, Has/Broadwell. By checking for AVX instead of
113 * model numbers we ensure also future Intel CPUs are covered.
115 const int nthreads_omp_faster_default
= 8;
116 const int nthreads_omp_faster_Nehalem
= 12;
117 const int nthreads_omp_faster_Intel_AVX
= 16;
118 const int nthreads_omp_faster_AMD_Ryzen
= 16;
119 /* For CPU only runs the fastest options are usually MPI or OpenMP only.
120 * With one GPU, using MPI only is almost never optimal, so we need to
121 * compare running pure OpenMP with combined MPI+OpenMP. This means higher
122 * OpenMP threads counts can still be ok. Multiplying the numbers above
123 * by a factor of 2 seems to be a good estimate.
125 const int nthreads_omp_faster_gpu_fac
= 2;
127 /* This is the case with MPI (2 or more MPI PP ranks).
128 * By default we will terminate with a fatal error when more than 8
129 * OpenMP thread are (indirectly) requested, since using less threads
130 * nearly always results in better performance.
131 * With thread-mpi and multiple GPUs or one GPU and too many threads
132 * we first try 6 OpenMP threads and then less until the number of MPI ranks
133 * is divisible by the number of GPUs.
135 #if GMX_OPENMP && GMX_MPI
136 const int nthreads_omp_mpi_ok_max
= 8;
137 const int nthreads_omp_mpi_ok_min_cpu
= 1;
139 const int nthreads_omp_mpi_ok_min_gpu
= 2;
140 const int nthreads_omp_mpi_target_max
= 6;
144 /*! \brief Returns the maximum OpenMP thread count for which using a single MPI rank
145 * should be faster than using multiple ranks with the same total thread count.
147 static int nthreads_omp_faster(const gmx::CpuInfo
&cpuInfo
, gmx_bool bUseGPU
)
151 if (cpuInfo
.vendor() == gmx::CpuInfo::Vendor::Intel
&&
152 cpuInfo
.feature(gmx::CpuInfo::Feature::X86_Avx
))
154 nth
= nthreads_omp_faster_Intel_AVX
;
156 else if (gmx::cpuIsX86Nehalem(cpuInfo
))
159 nth
= nthreads_omp_faster_Nehalem
;
161 else if (cpuInfo
.vendor() == gmx::CpuInfo::Vendor::Amd
&& cpuInfo
.family() >= 23)
164 nth
= nthreads_omp_faster_AMD_Ryzen
;
168 nth
= nthreads_omp_faster_default
;
173 nth
*= nthreads_omp_faster_gpu_fac
;
176 nth
= std::min(nth
, GMX_OPENMP_MAX_THREADS
);
181 /*! \brief Returns that maximum OpenMP thread count that passes the efficiency check */
182 gmx_unused
static int nthreads_omp_efficient_max(int gmx_unused nrank
,
183 const gmx::CpuInfo
&cpuInfo
,
186 #if GMX_OPENMP && GMX_MPI
189 return nthreads_omp_mpi_ok_max
;
194 return nthreads_omp_faster(cpuInfo
, bUseGPU
);
198 /*! \brief Return the number of thread-MPI ranks to use.
199 * This is chosen such that we can always obey our own efficiency checks.
201 gmx_unused
static int get_tmpi_omp_thread_division(const gmx_hw_info_t
*hwinfo
,
202 const gmx_hw_opt_t
&hw_opt
,
207 const gmx::CpuInfo
&cpuInfo
= *hwinfo
->cpuInfo
;
209 GMX_RELEASE_ASSERT(nthreads_tot
> 0, "There must be at least one thread per rank");
211 /* There are no separate PME nodes here, as we ensured in
212 * check_and_update_hw_opt that nthreads_tmpi>0 with PME nodes
213 * and a conditional ensures we would not have ended up here.
214 * Note that separate PME nodes might be switched on later.
220 /* When the user sets nthreads_omp, we can end up oversubscribing CPU cores
221 * if we simply start as many ranks as GPUs. To avoid this, we start as few
222 * tMPI ranks as necessary to avoid oversubscription and instead leave GPUs idle.
223 * If the user does not set the number of OpenMP threads, nthreads_omp==0 and
224 * this code has no effect.
226 GMX_RELEASE_ASSERT(hw_opt
.nthreads_omp
>= 0, "nthreads_omp is negative, but previous checks should have prevented this");
227 while (nrank
*hw_opt
.nthreads_omp
> hwinfo
->nthreads_hw_avail
&& nrank
> 1)
232 if (nthreads_tot
< nrank
)
234 /* #thread < #gpu is very unlikely, but if so: waste gpu(s) */
235 nrank
= nthreads_tot
;
237 else if (nthreads_tot
> nthreads_omp_faster(cpuInfo
, ngpu
> 0) ||
238 (ngpu
> 1 && nthreads_tot
/ngpu
> nthreads_omp_mpi_target_max
))
240 /* The high OpenMP thread count will likely result in sub-optimal
241 * performance. Increase the rank count to reduce the thread count
242 * per rank. This will lead to GPU sharing by MPI ranks/threads.
246 /* Increase the rank count as long as have we more than 6 OpenMP
247 * threads per rank or the number of hardware threads is not
248 * divisible by the rank count. Don't go below 2 OpenMP threads.
256 while (nthreads_tot
/nrank
> nthreads_omp_mpi_target_max
||
257 (nthreads_tot
/(ngpu
*(nshare
+ 1)) >= nthreads_omp_mpi_ok_min_gpu
&& nthreads_tot
% nrank
!= 0));
260 else if (hw_opt
.nthreads_omp
> 0)
262 /* Here we could oversubscribe, when we do, we issue a warning later */
263 nrank
= std::max(1, nthreads_tot
/hw_opt
.nthreads_omp
);
267 if (nthreads_tot
<= nthreads_omp_faster(cpuInfo
, ngpu
> 0))
269 /* Use pure OpenMP parallelization */
274 /* Don't use OpenMP parallelization */
275 nrank
= nthreads_tot
;
282 //! Return whether hyper threading is enabled.
284 gmxSmtIsEnabled(const gmx::HardwareTopology
&hwTop
)
286 return (hwTop
.supportLevel() >= gmx::HardwareTopology::SupportLevel::Basic
&& hwTop
.machine().sockets
[0].cores
[0].hwThreads
.size() > 1);
292 //! Handles checks for algorithms that must use a single rank.
293 class SingleRankChecker
297 SingleRankChecker() : value_(false), reasons_() {}
298 /*! \brief Call this function for each possible condition
299 under which a single rank is required, along with a string
300 describing the constraint when it is applied. */
301 void applyConstraint(bool condition
, const char *description
)
306 reasons_
.push_back(gmx::formatString("%s only supports a single rank.", description
));
309 //! After applying any conditions, is a single rank required?
310 bool mustUseOneRank() const
314 /*! \brief Return a formatted string to use when writing a
315 message when a single rank is required, (or empty if no
316 constraint exists.) */
317 std::string
getMessage() const
319 return formatAndJoin(reasons_
, "\n", gmx::IdentityFormatter());
323 std::vector
<std::string
> reasons_
;
328 /* Get the number of MPI ranks to use for thread-MPI based on how many
329 * were requested, which algorithms we're using,
330 * and how many particles there are.
331 * At the point we have already called check_and_update_hw_opt.
332 * Thus all options should be internally consistent and consistent
333 * with the hardware, except that ntmpi could be larger than #GPU.
335 int get_nthreads_mpi(const gmx_hw_info_t
*hwinfo
,
336 gmx_hw_opt_t
*hw_opt
,
337 const std::vector
<int> &gpuIdsToUse
,
340 const t_inputrec
*inputrec
,
341 const gmx_mtop_t
*mtop
,
342 const gmx::MDLogger
&mdlog
,
345 int nthreads_hw
, nthreads_tot_max
, nrank
, ngpu
;
346 int min_atoms_per_mpi_rank
;
348 const gmx::CpuInfo
&cpuInfo
= *hwinfo
->cpuInfo
;
349 const gmx::HardwareTopology
&hwTop
= *hwinfo
->hardwareTopology
;
353 GMX_RELEASE_ASSERT((EEL_PME(inputrec
->coulombtype
) || EVDW_PME(inputrec
->vdwtype
)) &&
354 pme_gpu_supports_input(inputrec
, nullptr),
355 "PME can't be on GPUs unless we are using PME");
357 // A single rank is all that is supported with PME on GPUs
358 if (hw_opt
->nthreads_tmpi
< 1)
362 if (hw_opt
->nthreads_tmpi
> 1)
364 gmx_fatal(FARGS
, "PME on GPUs is only supported with a single rank");
369 /* Check if an algorithm does not support parallel simulation. */
370 // TODO This might work better if e.g. implemented algorithms
371 // had to define a function that returns such requirements,
372 // and a description string.
373 SingleRankChecker checker
;
374 checker
.applyConstraint(inputrec
->eI
== eiLBFGS
, "L-BFGS minimization");
375 checker
.applyConstraint(inputrec
->coulombtype
== eelEWALD
, "Plain Ewald electrostatics");
376 checker
.applyConstraint(doMembed
, "Membrane embedding");
377 bool useOrientationRestraints
= (gmx_mtop_ftype_count(mtop
, F_ORIRES
) > 0);
378 checker
.applyConstraint(useOrientationRestraints
, "Orientation restraints");
379 if (checker
.mustUseOneRank())
381 std::string message
= checker
.getMessage();
382 if (hw_opt
->nthreads_tmpi
> 1)
384 gmx_fatal(FARGS
, "%s However, you asked for more than 1 thread-MPI rank, so mdrun cannot continue. Choose a single rank, or a different algorithm.", message
.c_str());
386 GMX_LOG(mdlog
.warning
).asParagraph().appendTextFormatted("%s Choosing to use only a single thread-MPI rank.", message
.c_str());
391 if (hw_opt
->nthreads_tmpi
> 0)
393 /* Trivial, return the user's choice right away */
394 return hw_opt
->nthreads_tmpi
;
397 // Now implement automatic selection of number of thread-MPI ranks
398 nthreads_hw
= hwinfo
->nthreads_hw_avail
;
400 if (nthreads_hw
<= 0)
402 /* This should normally not happen, but if it does, we handle it */
403 gmx_fatal(FARGS
, "The number of available hardware threads can not be detected, please specify the number of MPI ranks and the number of OpenMP threads (if supported) manually with options -ntmpi and -ntomp, respectively");
406 /* How many total (#tMPI*#OpenMP) threads can we start? */
407 if (hw_opt
->nthreads_tot
> 0)
409 nthreads_tot_max
= hw_opt
->nthreads_tot
;
413 nthreads_tot_max
= nthreads_hw
;
416 /* nonbondedOnGpu might be false e.g. because this simulation uses
417 * the group scheme, or is a rerun with energy groups. */
418 ngpu
= (nonbondedOnGpu
? static_cast<int>(gpuIdsToUse
.size()) : 0);
420 if (inputrec
->cutoff_scheme
== ecutsGROUP
)
422 /* We checked this before, but it doesn't hurt to do it once more */
423 GMX_RELEASE_ASSERT(hw_opt
->nthreads_omp
== 1, "The group scheme only supports one OpenMP thread per rank");
427 get_tmpi_omp_thread_division(hwinfo
, *hw_opt
, nthreads_tot_max
, ngpu
);
429 if (inputrec
->eI
== eiNM
|| EI_TPI(inputrec
->eI
))
431 /* Dims/steps are divided over the nodes iso splitting the atoms.
432 * With NM we can't have more ranks than #atoms*#dim. With TPI it's
433 * unlikely we have fewer atoms than ranks, and if so, communication
434 * would become a bottleneck, so we set the limit to 1 atom/rank.
436 min_atoms_per_mpi_rank
= 1;
442 min_atoms_per_mpi_rank
= min_atoms_per_gpu
;
446 min_atoms_per_mpi_rank
= min_atoms_per_mpi_thread
;
450 if (mtop
->natoms
/nrank
< min_atoms_per_mpi_rank
)
454 /* the rank number was chosen automatically, but there are too few
455 atoms per rank, so we need to reduce the rank count */
456 nrank_new
= std::max(1, mtop
->natoms
/min_atoms_per_mpi_rank
);
458 /* Avoid partial use of Hyper-Threading */
459 if (gmxSmtIsEnabled(hwTop
) &&
460 nrank_new
> nthreads_hw
/2 && nrank_new
< nthreads_hw
)
462 nrank_new
= nthreads_hw
/2;
465 /* If the user specified the total thread count, ensure this is
466 * divisible by the number of ranks.
467 * It is quite likely that we have too many total threads compared
468 * to the size of the system, but if the user asked for this many
469 * threads we should respect that.
471 while (hw_opt
->nthreads_tot
> 0 &&
472 hw_opt
->nthreads_tot
% nrank_new
!= 0)
477 /* Avoid large prime numbers in the rank count */
480 /* Use only 6,8,10 with additional factors of 2 */
484 while (3*fac
*2 <= nrank_new
)
489 nrank_new
= (nrank_new
/fac
)*fac
;
493 /* Avoid 5, since small system won't fit 5 domains along
494 * a dimension. This might lead to waisting some cores, but this
495 * will have a small impact in this regime of very small systems.
505 /* We reduced the number of tMPI ranks, which means we might violate
506 * our own efficiency checks if we simply use all hardware threads.
508 if (bHasOmpSupport
&& hw_opt
->nthreads_omp
<= 0 && hw_opt
->nthreads_tot
<= 0)
510 /* The user set neither the total nor the OpenMP thread count,
511 * we should use all hardware threads, unless we will violate
512 * our own efficiency limitation on the thread count.
516 nt_omp_max
= nthreads_omp_efficient_max(nrank
, cpuInfo
, ngpu
>= 1);
518 if (nrank
*nt_omp_max
< hwinfo
->nthreads_hw_avail
)
520 /* Limit the number of OpenMP threads to start */
521 hw_opt
->nthreads_omp
= nt_omp_max
;
525 fprintf(stderr
, "\n");
526 fprintf(stderr
, "NOTE: Parallelization is limited by the small number of atoms,\n");
527 fprintf(stderr
, " only starting %d thread-MPI ranks.\n", nrank
);
528 fprintf(stderr
, " You can use the -nt and/or -ntmpi option to optimize the number of threads.\n\n");
535 void check_resource_division_efficiency(const gmx_hw_info_t
*hwinfo
,
537 bool willUsePhysicalGpu
,
538 gmx_bool bNtOmpOptionSet
,
540 const gmx::MDLogger
&mdlog
)
542 #if GMX_OPENMP && GMX_MPI
543 int nth_omp_min
, nth_omp_max
;
546 const char *mpi_option
= " (option -ntmpi)";
548 const char *mpi_option
= "";
551 /* This function should be called after thread-MPI (when configured) and
552 * OpenMP have been initialized. Check that here.
555 GMX_RELEASE_ASSERT(nthreads_omp_faster_default
>= nthreads_omp_mpi_ok_max
, "Inconsistent OpenMP thread count default values");
557 GMX_RELEASE_ASSERT(gmx_omp_nthreads_get(emntDefault
) >= 1, "Must have at least one OpenMP thread");
559 nth_omp_min
= gmx_omp_nthreads_get(emntDefault
);
560 nth_omp_max
= gmx_omp_nthreads_get(emntDefault
);
562 bool anyRankIsUsingGpus
= willUsePhysicalGpu
;
563 /* Thread-MPI seems to have a bug with reduce on 1 node, so use a cond. */
564 if (cr
->nnodes
+ cr
->npmenodes
> 1)
566 int count
[3], count_max
[3];
568 count
[0] = -nth_omp_min
;
569 count
[1] = nth_omp_max
;
570 count
[2] = willUsePhysicalGpu
;
572 MPI_Allreduce(count
, count_max
, 3, MPI_INT
, MPI_MAX
, cr
->mpi_comm_mysim
);
574 /* In case of an inhomogeneous run setup we use the maximum counts */
575 nth_omp_min
= -count_max
[0];
576 nth_omp_max
= count_max
[1];
577 anyRankIsUsingGpus
= count_max
[2] > 0;
580 int nthreads_omp_mpi_ok_min
;
582 if (!anyRankIsUsingGpus
)
584 nthreads_omp_mpi_ok_min
= nthreads_omp_mpi_ok_min_cpu
;
588 /* With GPUs we set the minimum number of OpenMP threads to 2 to catch
589 * cases where the user specifies #ranks == #cores.
591 nthreads_omp_mpi_ok_min
= nthreads_omp_mpi_ok_min_gpu
;
594 if (DOMAINDECOMP(cr
) && cr
->nnodes
> 1)
596 if (nth_omp_max
< nthreads_omp_mpi_ok_min
||
597 nth_omp_max
> nthreads_omp_mpi_ok_max
)
599 /* Note that we print target_max here, not ok_max */
600 sprintf(buf
, "Your choice of number of MPI ranks and amount of resources results in using %d OpenMP threads per rank, which is most likely inefficient. The optimum is usually between %d and %d threads per rank.",
602 nthreads_omp_mpi_ok_min
,
603 nthreads_omp_mpi_target_max
);
607 GMX_LOG(mdlog
.warning
).asParagraph().appendTextFormatted("NOTE: %s", buf
);
611 /* This fatal error, and the one below, is nasty, but it's
612 * probably the only way to ensure that all users don't waste
613 * a lot of resources, since many users don't read logs/stderr.
615 gmx_fatal(FARGS
, "%s If you want to run with this setup, specify the -ntomp option. But we suggest to change the number of MPI ranks%s.", buf
, mpi_option
);
621 const gmx::CpuInfo
&cpuInfo
= *hwinfo
->cpuInfo
;
623 /* No domain decomposition (or only one domain) */
624 if (nth_omp_max
> nthreads_omp_faster(cpuInfo
, anyRankIsUsingGpus
))
626 /* To arrive here, the user/system set #ranks and/or #OMPthreads */
630 bEnvSet
= (getenv("OMP_NUM_THREADS") != nullptr);
632 if (bNtOmpOptionSet
|| bEnvSet
)
634 sprintf(buf2
, "You requested %d OpenMP threads", nth_omp_max
);
638 sprintf(buf2
, "Your choice of %d MPI rank%s and the use of %d total threads %sleads to the use of %d OpenMP threads",
639 cr
->nnodes
+ cr
->npmenodes
,
640 cr
->nnodes
+ cr
->npmenodes
== 1 ? "" : "s",
641 numTotalThreads
> 0 ? numTotalThreads
: hwinfo
->nthreads_hw_avail
,
642 hwinfo
->nphysicalnode
> 1 ? "on a node " : "",
645 sprintf(buf
, "%s, whereas we expect the optimum to be with more MPI ranks with %d to %d OpenMP threads.",
646 buf2
, nthreads_omp_mpi_ok_min
, nthreads_omp_mpi_target_max
);
648 /* We can not quit with a fatal error when OMP_NUM_THREADS is set
649 * with different values per rank or node, since in that case
650 * the user can not set -ntomp to override the error.
652 if (bNtOmpOptionSet
|| (bEnvSet
&& nth_omp_min
!= nth_omp_max
))
654 GMX_LOG(mdlog
.warning
).asParagraph().appendTextFormatted("NOTE: %s", buf
);
658 gmx_fatal(FARGS
, "%s If you want to run with this many OpenMP threads, specify the -ntomp option. But we suggest to increase the number of MPI ranks%s.", buf
, mpi_option
);
662 #else /* GMX_OPENMP && GMX_MPI */
663 /* No OpenMP and/or MPI: it doesn't make much sense to check */
664 GMX_UNUSED_VALUE(bNtOmpOptionSet
);
665 GMX_UNUSED_VALUE(numTotalThreads
);
666 GMX_UNUSED_VALUE(willUsePhysicalGpu
);
667 GMX_UNUSED_VALUE(cr
);
668 /* Check if we have more than 1 physical core, if detected,
669 * or more than 1 hardware thread if physical cores were not detected.
671 if (!GMX_OPENMP
&& !GMX_MPI
&& hwinfo
->hardwareTopology
->numberOfCores() > 1)
673 GMX_LOG(mdlog
.warning
).asParagraph().appendText("NOTE: GROMACS was compiled without OpenMP and (thread-)MPI support, can only use a single CPU core");
675 #endif /* GMX_OPENMP && GMX_MPI */
679 //! Dump a \c hw_opt to \c fp.
680 static void print_hw_opt(FILE *fp
, const gmx_hw_opt_t
*hw_opt
)
682 fprintf(fp
, "hw_opt: nt %d ntmpi %d ntomp %d ntomp_pme %d gpu_id '%s' gputasks '%s'\n",
683 hw_opt
->nthreads_tot
,
684 hw_opt
->nthreads_tmpi
,
685 hw_opt
->nthreads_omp
,
686 hw_opt
->nthreads_omp_pme
,
687 hw_opt
->gpuIdsAvailable
.c_str(),
688 hw_opt
->userGpuTaskAssignment
.c_str());
691 void check_and_update_hw_opt_1(gmx_hw_opt_t
*hw_opt
,
695 /* Currently hw_opt only contains default settings or settings supplied
696 * by the user on the command line.
698 if (hw_opt
->nthreads_omp
< 0)
700 gmx_fatal(FARGS
, "The number of OpenMP threads supplied on the command line is %d, which is negative and not allowed", hw_opt
->nthreads_omp
);
703 /* Check for OpenMP settings stored in environment variables, which can
704 * potentially be different on different MPI ranks.
706 gmx_omp_nthreads_read_env(&hw_opt
->nthreads_omp
, SIMMASTER(cr
));
708 /* Check restrictions on the user supplied options before modifying them.
709 * TODO: Put the user values in a const struct and preserve them.
712 if (hw_opt
->nthreads_tot
> 0)
714 gmx_fatal(FARGS
, "Setting the total number of threads is only supported with thread-MPI and GROMACS was compiled without thread-MPI");
716 if (hw_opt
->nthreads_tmpi
> 0)
718 gmx_fatal(FARGS
, "Setting the number of thread-MPI ranks is only supported with thread-MPI and GROMACS was compiled without thread-MPI");
724 /* Check restrictions on PME thread related options set by the user */
726 if (hw_opt
->nthreads_omp_pme
> 0 && hw_opt
->nthreads_omp
<= 0)
728 gmx_fatal(FARGS
, "You need to specify -ntomp in addition to -ntomp_pme");
731 if (hw_opt
->nthreads_omp_pme
>= 1 &&
732 hw_opt
->nthreads_omp_pme
!= hw_opt
->nthreads_omp
&&
735 /* This can result in a fatal error on many MPI ranks,
736 * but since the thread count can differ per rank,
737 * we can't easily avoid this.
739 gmx_fatal(FARGS
, "You need to explicitly specify the number of PME ranks (-npme) when using different number of OpenMP threads for PP and PME ranks");
744 /* GROMACS was configured without OpenMP support */
746 if (hw_opt
->nthreads_omp
> 1 || hw_opt
->nthreads_omp_pme
> 1)
748 gmx_fatal(FARGS
, "More than 1 OpenMP thread requested, but GROMACS was compiled without OpenMP support");
750 hw_opt
->nthreads_omp
= 1;
751 hw_opt
->nthreads_omp_pme
= 1;
754 if (hw_opt
->nthreads_tot
> 0 && hw_opt
->nthreads_omp_pme
<= 0)
756 /* We have the same number of OpenMP threads for PP and PME ranks,
757 * thus we can perform several consistency checks.
759 if (hw_opt
->nthreads_tmpi
> 0 &&
760 hw_opt
->nthreads_omp
> 0 &&
761 hw_opt
->nthreads_tot
!= hw_opt
->nthreads_tmpi
*hw_opt
->nthreads_omp
)
763 gmx_fatal(FARGS
, "The total number of threads requested (%d) does not match the thread-MPI ranks (%d) times the OpenMP threads (%d) requested",
764 hw_opt
->nthreads_tot
, hw_opt
->nthreads_tmpi
, hw_opt
->nthreads_omp
);
767 if (hw_opt
->nthreads_tmpi
> 0 &&
768 hw_opt
->nthreads_tot
% hw_opt
->nthreads_tmpi
!= 0)
770 gmx_fatal(FARGS
, "The total number of threads requested (%d) is not divisible by the number of thread-MPI ranks requested (%d)",
771 hw_opt
->nthreads_tot
, hw_opt
->nthreads_tmpi
);
774 if (hw_opt
->nthreads_omp
> 0 &&
775 hw_opt
->nthreads_tot
% hw_opt
->nthreads_omp
!= 0)
777 gmx_fatal(FARGS
, "The total number of threads requested (%d) is not divisible by the number of OpenMP threads requested (%d)",
778 hw_opt
->nthreads_tot
, hw_opt
->nthreads_omp
);
782 if (hw_opt
->nthreads_tot
> 0)
784 if (hw_opt
->nthreads_omp
> hw_opt
->nthreads_tot
)
786 gmx_fatal(FARGS
, "You requested %d OpenMP threads with %d total threads. Choose a total number of threads that is a multiple of the number of OpenMP threads.",
787 hw_opt
->nthreads_omp
, hw_opt
->nthreads_tot
);
790 if (hw_opt
->nthreads_tmpi
> hw_opt
->nthreads_tot
)
792 gmx_fatal(FARGS
, "You requested %d thread-MPI ranks with %d total threads. Choose a total number of threads that is a multiple of the number of thread-MPI ranks.",
793 hw_opt
->nthreads_tmpi
, hw_opt
->nthreads_tot
);
799 print_hw_opt(debug
, hw_opt
);
802 /* Asserting this simplifies the hardware resource division later
804 GMX_RELEASE_ASSERT(!(hw_opt
->nthreads_omp_pme
>= 1 && hw_opt
->nthreads_omp
<= 0),
805 "PME thread count should only be set when the normal thread count is also set");
808 void check_and_update_hw_opt_2(gmx_hw_opt_t
*hw_opt
,
811 if (cutoff_scheme
== ecutsGROUP
)
813 /* We only have OpenMP support for PME only nodes */
814 if (hw_opt
->nthreads_omp
> 1)
816 gmx_fatal(FARGS
, "OpenMP threads have been requested with cut-off scheme %s, but these are only supported with cut-off scheme %s",
817 ecutscheme_names
[cutoff_scheme
],
818 ecutscheme_names
[ecutsVERLET
]);
820 hw_opt
->nthreads_omp
= 1;
824 void checkAndUpdateRequestedNumOpenmpThreads(gmx_hw_opt_t
*hw_opt
,
825 const gmx_hw_info_t
&hwinfo
,
827 PmeRunMode pmeRunMode
,
828 const gmx_mtop_t
&mtop
)
831 GMX_RELEASE_ASSERT(hw_opt
->nthreads_tmpi
>= 1, "Must have at least one thread-MPI rank");
833 /* If the user set the total number of threads on the command line
834 * and did not specify the number of OpenMP threads, set the latter here.
836 if (hw_opt
->nthreads_tot
> 0 && hw_opt
->nthreads_omp
<= 0)
838 hw_opt
->nthreads_omp
= hw_opt
->nthreads_tot
/hw_opt
->nthreads_tmpi
;
840 if (!bHasOmpSupport
&& hw_opt
->nthreads_omp
> 1)
842 gmx_fatal(FARGS
, "You (indirectly) asked for OpenMP threads by setting -nt > -ntmpi, but GROMACS was compiled without OpenMP support");
847 /* With both non-bonded and PME on GPU, the work left on the CPU is often
848 * (much) slower with SMT than without SMT. This is mostly the case with
849 * few atoms per core. Thus, if the number of threads is set to auto,
850 * we turn off SMT in that case. Note that PME on GPU implies that also
851 * the non-bonded are computed on the GPU.
852 * We only need to do this when the number of hardware theads is larger
853 * than the number of cores. Note that a queuing system could limit
854 * the number of hardware threads available, but we are not trying to be
855 * too smart here in that case.
857 /* The number 10000 atoms per core is a reasonable threshold
858 * for Intel x86 and AMD Threadripper. Note that around this regime
859 * the performance is not very sensitive to the SMT choice.
861 constexpr int c_numAtomsPerCoreSmtThreshold
= 10000;
863 if (bHasOmpSupport
&&
864 hw_opt
->nthreads_omp
<= 0 &&
865 hwinfo
.nthreads_hw_avail
> hwinfo
.hardwareTopology
->numberOfCores() &&
866 pmeRunMode
== PmeRunMode::GPU
)
868 /* Note that the queing system might have limited us from using
869 * all detected ncore_tot physical cores. We are currently not
870 * checking for that here.
872 int numAtomsPerPhysicalCore
= mtop
.natoms
/hwinfo
.ncore_tot
;
873 if (numAtomsPerPhysicalCore
< c_numAtomsPerCoreSmtThreshold
)
875 int numRanksInThisNode
= (cr
? cr
->nrank_intranode
: 1);
876 /* Choose one OpenMP thread per physical core */
877 hw_opt
->nthreads_omp
= std::max(1, hwinfo
.hardwareTopology
->numberOfCores()/numRanksInThisNode
);
881 GMX_RELEASE_ASSERT(bHasOmpSupport
|| hw_opt
->nthreads_omp
== 1, "Without OpenMP support, only one thread per rank can be used");
883 /* We are done with updating nthreads_omp, we can set nthreads_omp_pme */
884 if (hw_opt
->nthreads_omp_pme
<= 0 && hw_opt
->nthreads_omp
> 0)
886 hw_opt
->nthreads_omp_pme
= hw_opt
->nthreads_omp
;
891 print_hw_opt(debug
, hw_opt
);