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.
38 #include "resource-division.h"
47 #include "gromacs/gmxlib/md_logging.h"
48 #include "gromacs/hardware/cpuinfo.h"
49 #include "gromacs/hardware/detecthardware.h"
50 #include "gromacs/hardware/gpu_hw_info.h"
51 #include "gromacs/hardware/hardwaretopology.h"
52 #include "gromacs/hardware/hw_info.h"
53 #include "gromacs/mdlib/gmx_omp_nthreads.h"
54 #include "gromacs/mdtypes/commrec.h"
55 #include "gromacs/mdtypes/inputrec.h"
56 #include "gromacs/mdtypes/md_enums.h"
57 #include "gromacs/topology/mtop_util.h"
58 #include "gromacs/topology/topology.h"
59 #include "gromacs/utility/fatalerror.h"
60 #include "gromacs/utility/gmxassert.h"
61 #include "gromacs/utility/stringutil.h"
64 /* DISCLAIMER: All the atom count and thread numbers below are heuristic.
65 * The real switching points will depend on the system simulation,
66 * the algorithms used and the hardware it's running on, as well as if there
67 * are other jobs running on the same machine. We try to take into account
68 * factors that have a large influence, such as recent Intel CPUs being
69 * much better at wide multi-threading. The remaining factors should
70 * (hopefully) have a small influence, such that the performance just before
71 * and after a switch point doesn't change too much.
74 static const bool bHasOmpSupport = GMX_OPENMP;
77 /* The minimum number of atoms per tMPI thread. With fewer atoms than this,
78 * the number of threads will get lowered.
80 static const int min_atoms_per_mpi_thread = 90;
81 static const int min_atoms_per_gpu = 900;
82 #endif /* GMX_THREAD_MPI */
84 /* TODO choose nthreads_omp based on hardware topology
85 when we have a hardware topology detection library */
86 /* First we consider the case of no MPI (1 MPI rank).
87 * In general, when running up to 8 threads, OpenMP should be faster.
88 * Note: on AMD Bulldozer we should avoid running OpenMP over two dies.
89 * On Intel>=Nehalem running OpenMP on a single CPU is always faster,
90 * even on two CPUs it's usually faster (but with many OpenMP threads
91 * it could be faster not to use HT, currently we always use HT).
92 * On Nehalem/Westmere we want to avoid running 16 threads over
93 * two CPUs with HT, so we need a limit<16; thus we use 12.
94 * A reasonable limit for Intel Sandy and Ivy bridge,
95 * not knowing the topology, is 16 threads.
96 * Below we check for Intel and AVX, which for now includes
97 * Sandy/Ivy Bridge, Has/Broadwell. By checking for AVX instead of
98 * model numbers we ensure also future Intel CPUs are covered.
100 const int nthreads_omp_faster_default = 8;
101 const int nthreads_omp_faster_Nehalem = 12;
102 const int nthreads_omp_faster_Intel_AVX = 16;
103 /* For CPU only runs the fastest options are usually MPI or OpenMP only.
104 * With one GPU, using MPI only is almost never optimal, so we need to
105 * compare running pure OpenMP with combined MPI+OpenMP. This means higher
106 * OpenMP threads counts can still be ok. Multiplying the numbers above
107 * by a factor of 2 seems to be a good estimate.
109 const int nthreads_omp_faster_gpu_fac = 2;
111 /* This is the case with MPI (2 or more MPI PP ranks).
112 * By default we will terminate with a fatal error when more than 8
113 * OpenMP thread are (indirectly) requested, since using less threads
114 * nearly always results in better performance.
115 * With thread-mpi and multiple GPUs or one GPU and too many threads
116 * we first try 6 OpenMP threads and then less until the number of MPI ranks
117 * is divisible by the number of GPUs.
119 #if GMX_OPENMP && GMX_MPI
120 const int nthreads_omp_mpi_ok_max = 8;
121 const int nthreads_omp_mpi_ok_min_cpu = 1;
123 const int nthreads_omp_mpi_ok_min_gpu = 2;
124 const int nthreads_omp_mpi_target_max = 6;
127 /* Returns the maximum OpenMP thread count for which using a single MPI rank
128 * should be faster than using multiple ranks with the same total thread count.
130 static int nthreads_omp_faster(const gmx::CpuInfo &cpuInfo, gmx_bool bUseGPU)
134 if (cpuInfo.vendor() == gmx::CpuInfo::Vendor::Intel &&
135 cpuInfo.feature(gmx::CpuInfo::Feature::X86_Avx))
137 nth = nthreads_omp_faster_Intel_AVX;
139 else if (gmx::cpuIsX86Nehalem(cpuInfo))
142 nth = nthreads_omp_faster_Nehalem;
146 nth = nthreads_omp_faster_default;
151 nth *= nthreads_omp_faster_gpu_fac;
154 nth = std::min(nth, GMX_OPENMP_MAX_THREADS);
159 /* Returns that maximum OpenMP thread count that passes the efficiency check */
160 static int nthreads_omp_efficient_max(int gmx_unused nrank,
161 const gmx::CpuInfo &cpuInfo,
164 #if GMX_OPENMP && GMX_MPI
167 return nthreads_omp_mpi_ok_max;
172 return nthreads_omp_faster(cpuInfo, bUseGPU);
176 /* Return the number of thread-MPI ranks to use.
177 * This is chosen such that we can always obey our own efficiency checks.
179 static int get_tmpi_omp_thread_division(const gmx_hw_info_t *hwinfo,
180 const gmx_hw_opt_t *hw_opt,
185 const gmx::CpuInfo &cpuInfo = *hwinfo->cpuInfo;
187 GMX_RELEASE_ASSERT(nthreads_tot > 0, "There must be at least one thread per rank");
189 /* There are no separate PME nodes here, as we ensured in
190 * check_and_update_hw_opt that nthreads_tmpi>0 with PME nodes
191 * and a conditional ensures we would not have ended up here.
192 * Note that separate PME nodes might be switched on later.
198 /* When the user sets nthreads_omp, we can end up oversubscribing CPU cores
199 * if we simply start as many ranks as GPUs. To avoid this, we start as few
200 * tMPI ranks as necessary to avoid oversubscription and instead leave GPUs idle.
201 * If the user does not set the number of OpenMP threads, nthreads_omp==0 and
202 * this code has no effect.
204 GMX_RELEASE_ASSERT(hw_opt->nthreads_omp >= 0, "nthreads_omp is negative, but previous checks should have prevented this");
205 while (nrank*hw_opt->nthreads_omp > hwinfo->nthreads_hw_avail && nrank > 1)
210 if (nthreads_tot < nrank)
212 /* #thread < #gpu is very unlikely, but if so: waste gpu(s) */
213 nrank = nthreads_tot;
215 else if (gmx_gpu_sharing_supported() &&
216 (nthreads_tot > nthreads_omp_faster(cpuInfo, ngpu > 0) ||
217 (ngpu > 1 && nthreads_tot/ngpu > nthreads_omp_mpi_target_max)))
219 /* The high OpenMP thread count will likely result in sub-optimal
220 * performance. Increase the rank count to reduce the thread count
221 * per rank. This will lead to GPU sharing by MPI ranks/threads.
225 /* Increase the rank count as long as have we more than 6 OpenMP
226 * threads per rank or the number of hardware threads is not
227 * divisible by the rank count. Don't go below 2 OpenMP threads.
235 while (nthreads_tot/nrank > nthreads_omp_mpi_target_max ||
236 (nthreads_tot/(ngpu*(nshare + 1)) >= nthreads_omp_mpi_ok_min_gpu && nthreads_tot % nrank != 0));
239 else if (hw_opt->nthreads_omp > 0)
241 /* Here we could oversubscribe, when we do, we issue a warning later */
242 nrank = std::max(1, nthreads_tot/hw_opt->nthreads_omp);
246 if (nthreads_tot <= nthreads_omp_faster(cpuInfo, ngpu > 0))
248 /* Use pure OpenMP parallelization */
253 /* Don't use OpenMP parallelization */
254 nrank = nthreads_tot;
262 static int getMaxGpuUsable(FILE *fplog, const t_commrec *cr, const gmx_hw_info_t *hwinfo,
263 int cutoff_scheme, gmx_bool bUseGpu)
265 /* This code relies on the fact that GPU are not detected when GPU
266 * acceleration was disabled at run time by the user.
268 if (cutoff_scheme == ecutsVERLET &&
270 hwinfo->gpu_info.n_dev_compatible > 0)
272 if (gmx_multiple_gpu_per_node_supported())
274 return hwinfo->gpu_info.n_dev_compatible;
278 if (hwinfo->gpu_info.n_dev_compatible > 1)
280 md_print_warn(cr, fplog, "More than one compatible GPU is available, but GROMACS can only use one of them. Using a single thread-MPI rank.\n");
295 gmxSmtIsEnabled(const gmx::HardwareTopology &hwTop)
297 return (hwTop.supportLevel() >= gmx::HardwareTopology::SupportLevel::Basic && hwTop.machine().sockets[0].cores[0].hwThreads.size() > 1);
303 class SingleRankChecker
307 SingleRankChecker() : value_(false), reasons_() {}
308 /*! \brief Call this function for each possible condition
309 under which a single rank is required, along with a string
310 describing the constraint when it is applied. */
311 void applyConstraint(bool condition, const char *description)
316 reasons_.push_back(gmx::formatString("%s only supports a single rank.", description));
319 //! After applying any conditions, is a single rank required?
320 bool mustUseOneRank() const
324 /*! \brief Return a formatted string to use when writing a
325 message when a single rank is required, (or empty if no
326 constraint exists.) */
327 std::string getMessage() const
329 return formatAndJoin(reasons_, "\n", gmx::IdentityFormatter());
333 std::vector<std::string> reasons_;
338 /* Get the number of MPI ranks to use for thread-MPI based on how many
339 * were requested, which algorithms we're using,
340 * and how many particles there are.
341 * At the point we have already called check_and_update_hw_opt.
342 * Thus all options should be internally consistent and consistent
343 * with the hardware, except that ntmpi could be larger than #GPU.
345 int get_nthreads_mpi(const gmx_hw_info_t *hwinfo,
346 gmx_hw_opt_t *hw_opt,
347 const t_inputrec *inputrec,
348 const gmx_mtop_t *mtop,
354 int nthreads_hw, nthreads_tot_max, nrank, ngpu;
355 int min_atoms_per_mpi_rank;
357 const gmx::CpuInfo &cpuInfo = *hwinfo->cpuInfo;
358 const gmx::HardwareTopology &hwTop = *hwinfo->hardwareTopology;
361 /* Check if an algorithm does not support parallel simulation. */
362 // TODO This might work better if e.g. implemented algorithms
363 // had to define a function that returns such requirements,
364 // and a description string.
365 SingleRankChecker checker;
366 checker.applyConstraint(inputrec->eI == eiLBFGS, "L-BFGS minimization");
367 checker.applyConstraint(inputrec->coulombtype == eelEWALD, "Plain Ewald electrostatics");
368 checker.applyConstraint(doMembed, "Membrane embedding");
369 bool useOrientationRestraints = (gmx_mtop_ftype_count(mtop, F_ORIRES) > 0);
370 checker.applyConstraint(useOrientationRestraints, "Orientation restraints");
371 if (checker.mustUseOneRank())
373 std::string message = checker.getMessage();
374 if (hw_opt->nthreads_tmpi > 1)
376 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());
378 md_print_warn(cr, fplog, "%s Choosing to use only a single thread-MPI rank.", message.c_str());
383 if (hw_opt->nthreads_tmpi > 0)
385 /* Trivial, return right away */
386 return hw_opt->nthreads_tmpi;
389 nthreads_hw = hwinfo->nthreads_hw_avail;
391 if (nthreads_hw <= 0)
393 /* This should normally not happen, but if it does, we handle it */
394 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");
397 /* How many total (#tMPI*#OpenMP) threads can we start? */
398 if (hw_opt->nthreads_tot > 0)
400 nthreads_tot_max = hw_opt->nthreads_tot;
404 nthreads_tot_max = nthreads_hw;
407 ngpu = getMaxGpuUsable(fplog, cr, hwinfo, inputrec->cutoff_scheme, bUseGpu);
409 if (inputrec->cutoff_scheme == ecutsGROUP)
411 /* We checked this before, but it doesn't hurt to do it once more */
412 GMX_RELEASE_ASSERT(hw_opt->nthreads_omp == 1, "The group scheme only supports one OpenMP thread per rank");
416 get_tmpi_omp_thread_division(hwinfo, hw_opt, nthreads_tot_max, ngpu);
418 if (inputrec->eI == eiNM || EI_TPI(inputrec->eI))
420 /* Dims/steps are divided over the nodes iso splitting the atoms.
421 * With NM we can't have more ranks than #atoms*#dim. With TPI it's
422 * unlikely we have fewer atoms than ranks, and if so, communication
423 * would become a bottleneck, so we set the limit to 1 atom/rank.
425 min_atoms_per_mpi_rank = 1;
431 min_atoms_per_mpi_rank = min_atoms_per_gpu;
435 min_atoms_per_mpi_rank = min_atoms_per_mpi_thread;
439 if (mtop->natoms/nrank < min_atoms_per_mpi_rank)
443 /* the rank number was chosen automatically, but there are too few
444 atoms per rank, so we need to reduce the rank count */
445 nrank_new = std::max(1, mtop->natoms/min_atoms_per_mpi_rank);
447 /* Avoid partial use of Hyper-Threading */
448 if (gmxSmtIsEnabled(hwTop) &&
449 nrank_new > nthreads_hw/2 && nrank_new < nthreads_hw)
451 nrank_new = nthreads_hw/2;
454 /* If the user specified the total thread count, ensure this is
455 * divisible by the number of ranks.
456 * It is quite likely that we have too many total threads compared
457 * to the size of the system, but if the user asked for this many
458 * threads we should respect that.
460 while (hw_opt->nthreads_tot > 0 &&
461 hw_opt->nthreads_tot % nrank_new != 0)
466 /* Avoid large prime numbers in the rank count */
469 /* Use only 6,8,10 with additional factors of 2 */
473 while (3*fac*2 <= nrank_new)
478 nrank_new = (nrank_new/fac)*fac;
482 /* Avoid 5, since small system won't fit 5 domains along
483 * a dimension. This might lead to waisting some cores, but this
484 * will have a small impact in this regime of very small systems.
494 /* We reduced the number of tMPI ranks, which means we might violate
495 * our own efficiency checks if we simply use all hardware threads.
497 if (bHasOmpSupport && hw_opt->nthreads_omp <= 0 && hw_opt->nthreads_tot <= 0)
499 /* The user set neither the total nor the OpenMP thread count,
500 * we should use all hardware threads, unless we will violate
501 * our own efficiency limitation on the thread count.
505 nt_omp_max = nthreads_omp_efficient_max(nrank, cpuInfo, ngpu >= 1);
507 if (nrank*nt_omp_max < hwinfo->nthreads_hw_avail)
509 /* Limit the number of OpenMP threads to start */
510 hw_opt->nthreads_omp = nt_omp_max;
514 fprintf(stderr, "\n");
515 fprintf(stderr, "NOTE: Parallelization is limited by the small number of atoms,\n");
516 fprintf(stderr, " only starting %d thread-MPI ranks.\n", nrank);
517 fprintf(stderr, " You can use the -nt and/or -ntmpi option to optimize the number of threads.\n\n");
522 #endif /* GMX_THREAD_MPI */
525 void check_resource_division_efficiency(const gmx_hw_info_t *hwinfo,
526 const gmx_hw_opt_t *hw_opt,
527 gmx_bool bNtOmpOptionSet,
531 #if GMX_OPENMP && GMX_MPI
532 int nth_omp_min, nth_omp_max, ngpu;
535 const char *mpi_option = " (option -ntmpi)";
537 const char *mpi_option = "";
540 /* This function should be called after thread-MPI (when configured) and
541 * OpenMP have been initialized. Check that here.
544 GMX_RELEASE_ASSERT(nthreads_omp_faster_default >= nthreads_omp_mpi_ok_max, "Inconsistent OpenMP thread count default values");
545 GMX_RELEASE_ASSERT(hw_opt->nthreads_tmpi >= 1, "Must have at least one thread-MPI rank");
547 GMX_RELEASE_ASSERT(gmx_omp_nthreads_get(emntDefault) >= 1, "Must have at least one OpenMP thread");
549 nth_omp_min = gmx_omp_nthreads_get(emntDefault);
550 nth_omp_max = gmx_omp_nthreads_get(emntDefault);
551 ngpu = hw_opt->gpu_opt.n_dev_use;
553 /* Thread-MPI seems to have a bug with reduce on 1 node, so use a cond. */
554 if (cr->nnodes + cr->npmenodes > 1)
556 int count[3], count_max[3];
558 count[0] = -nth_omp_min;
559 count[1] = nth_omp_max;
562 MPI_Allreduce(count, count_max, 3, MPI_INT, MPI_MAX, cr->mpi_comm_mysim);
564 /* In case of an inhomogeneous run setup we use the maximum counts */
565 nth_omp_min = -count_max[0];
566 nth_omp_max = count_max[1];
570 int nthreads_omp_mpi_ok_min;
574 nthreads_omp_mpi_ok_min = nthreads_omp_mpi_ok_min_cpu;
578 /* With GPUs we set the minimum number of OpenMP threads to 2 to catch
579 * cases where the user specifies #ranks == #cores.
581 nthreads_omp_mpi_ok_min = nthreads_omp_mpi_ok_min_gpu;
584 if (DOMAINDECOMP(cr) && cr->nnodes > 1)
586 if (nth_omp_max < nthreads_omp_mpi_ok_min ||
587 (!(ngpu > 0 && !gmx_gpu_sharing_supported()) &&
588 nth_omp_max > nthreads_omp_mpi_ok_max))
590 /* Note that we print target_max here, not ok_max */
591 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.",
593 nthreads_omp_mpi_ok_min,
594 nthreads_omp_mpi_target_max);
598 md_print_warn(cr, fplog, "NOTE: %s\n", buf);
602 /* This fatal error, and the one below, is nasty, but it's
603 * probably the only way to ensure that all users don't waste
604 * a lot of resources, since many users don't read logs/stderr.
606 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);
612 const gmx::CpuInfo &cpuInfo = *hwinfo->cpuInfo;
614 /* No domain decomposition (or only one domain) */
615 if (!(ngpu > 0 && !gmx_gpu_sharing_supported()) &&
616 nth_omp_max > nthreads_omp_faster(cpuInfo, ngpu > 0))
618 /* To arrive here, the user/system set #ranks and/or #OMPthreads */
622 bEnvSet = (getenv("OMP_NUM_THREADS") != NULL);
624 if (bNtOmpOptionSet || bEnvSet)
626 sprintf(buf2, "You requested %d OpenMP threads", nth_omp_max);
630 sprintf(buf2, "Your choice of %d MPI rank%s and the use of %d total threads %sleads to the use of %d OpenMP threads",
631 cr->nnodes + cr->npmenodes,
632 cr->nnodes + cr->npmenodes == 1 ? "" : "s",
633 hw_opt->nthreads_tot > 0 ? hw_opt->nthreads_tot : hwinfo->nthreads_hw_avail,
634 hwinfo->nphysicalnode > 1 ? "on a node " : "",
637 sprintf(buf, "%s, whereas we expect the optimum to be with more MPI ranks with %d to %d OpenMP threads.",
638 buf2, nthreads_omp_mpi_ok_min, nthreads_omp_mpi_target_max);
640 /* We can not quit with a fatal error when OMP_NUM_THREADS is set
641 * with different values per rank or node, since in that case
642 * the user can not set -ntomp to override the error.
644 if (bNtOmpOptionSet || (bEnvSet && nth_omp_min != nth_omp_max))
646 md_print_warn(cr, fplog, "NOTE: %s\n", buf);
650 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);
654 #else /* GMX_OPENMP && GMX_MPI */
655 /* No OpenMP and/or MPI: it doesn't make much sense to check */
656 GMX_UNUSED_VALUE(hw_opt);
657 GMX_UNUSED_VALUE(bNtOmpOptionSet);
658 /* Check if we have more than 1 physical core, if detected,
659 * or more than 1 hardware thread if physical cores were not detected.
661 #if !GMX_OPENMP && !GMX_MPI
662 if (hwinfo->hardwareTopology->numberOfCores() > 1)
664 md_print_warn(cr, fplog, "NOTE: GROMACS was compiled without OpenMP and (thread-)MPI support, can only use a single CPU core\n");
667 GMX_UNUSED_VALUE(hwinfo);
668 GMX_UNUSED_VALUE(cr);
669 GMX_UNUSED_VALUE(fplog);
672 #endif /* GMX_OPENMP && GMX_MPI */
676 static void print_hw_opt(FILE *fp, const gmx_hw_opt_t *hw_opt)
678 fprintf(fp, "hw_opt: nt %d ntmpi %d ntomp %d ntomp_pme %d gpu_id '%s'\n",
679 hw_opt->nthreads_tot,
680 hw_opt->nthreads_tmpi,
681 hw_opt->nthreads_omp,
682 hw_opt->nthreads_omp_pme,
683 hw_opt->gpu_opt.gpu_id != NULL ? hw_opt->gpu_opt.gpu_id : "");
686 /* Checks we can do when we don't (yet) know the cut-off scheme */
687 void check_and_update_hw_opt_1(gmx_hw_opt_t *hw_opt,
691 /* Currently hw_opt only contains default settings or settings supplied
692 * by the user on the command line.
694 if (hw_opt->nthreads_omp < 0)
696 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);
699 /* Check for OpenMP settings stored in environment variables, which can
700 * potentially be different on different MPI ranks.
702 gmx_omp_nthreads_read_env(&hw_opt->nthreads_omp, SIMMASTER(cr));
704 /* Check restrictions on the user supplied options before modifying them.
705 * TODO: Put the user values in a const struct and preserve them.
708 if (hw_opt->nthreads_tot > 0)
710 gmx_fatal(FARGS, "Setting the total number of threads is only supported with thread-MPI and GROMACS was compiled without thread-MPI");
712 if (hw_opt->nthreads_tmpi > 0)
714 gmx_fatal(FARGS, "Setting the number of thread-MPI ranks is only supported with thread-MPI and GROMACS was compiled without thread-MPI");
720 /* Check restrictions on PME thread related options set by the user */
722 if (hw_opt->nthreads_omp_pme > 0 && hw_opt->nthreads_omp <= 0)
724 gmx_fatal(FARGS, "You need to specify -ntomp in addition to -ntomp_pme");
727 if (hw_opt->nthreads_omp_pme >= 1 &&
728 hw_opt->nthreads_omp_pme != hw_opt->nthreads_omp &&
731 /* This can result in a fatal error on many MPI ranks,
732 * but since the thread count can differ per rank,
733 * we can't easily avoid this.
735 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");
740 /* GROMACS was configured without OpenMP support */
742 if (hw_opt->nthreads_omp > 1 || hw_opt->nthreads_omp_pme > 1)
744 gmx_fatal(FARGS, "More than 1 OpenMP thread requested, but GROMACS was compiled without OpenMP support");
746 hw_opt->nthreads_omp = 1;
747 hw_opt->nthreads_omp_pme = 1;
750 if (hw_opt->nthreads_tot > 0 && hw_opt->nthreads_omp_pme <= 0)
752 /* We have the same number of OpenMP threads for PP and PME ranks,
753 * thus we can perform several consistency checks.
755 if (hw_opt->nthreads_tmpi > 0 &&
756 hw_opt->nthreads_omp > 0 &&
757 hw_opt->nthreads_tot != hw_opt->nthreads_tmpi*hw_opt->nthreads_omp)
759 gmx_fatal(FARGS, "The total number of threads requested (%d) does not match the thread-MPI ranks (%d) times the OpenMP threads (%d) requested",
760 hw_opt->nthreads_tot, hw_opt->nthreads_tmpi, hw_opt->nthreads_omp);
763 if (hw_opt->nthreads_tmpi > 0 &&
764 hw_opt->nthreads_tot % hw_opt->nthreads_tmpi != 0)
766 gmx_fatal(FARGS, "The total number of threads requested (%d) is not divisible by the number of thread-MPI ranks requested (%d)",
767 hw_opt->nthreads_tot, hw_opt->nthreads_tmpi);
770 if (hw_opt->nthreads_omp > 0 &&
771 hw_opt->nthreads_tot % hw_opt->nthreads_omp != 0)
773 gmx_fatal(FARGS, "The total number of threads requested (%d) is not divisible by the number of OpenMP threads requested (%d)",
774 hw_opt->nthreads_tot, hw_opt->nthreads_omp);
778 if (hw_opt->nthreads_tot == 1)
780 hw_opt->nthreads_tmpi = 1;
782 if (hw_opt->nthreads_omp > 1)
784 gmx_fatal(FARGS, "You requested %d OpenMP threads with %d total threads",
785 hw_opt->nthreads_tmpi, hw_opt->nthreads_tot);
787 hw_opt->nthreads_omp = 1;
788 hw_opt->nthreads_omp_pme = 1;
791 /* Parse GPU IDs, if provided.
792 * We check consistency with the tMPI thread count later.
794 gmx_parse_gpu_ids(&hw_opt->gpu_opt);
797 if (hw_opt->gpu_opt.n_dev_use > 0 && hw_opt->nthreads_tmpi == 0)
799 /* Set the number of MPI threads equal to the number of GPUs */
800 hw_opt->nthreads_tmpi = hw_opt->gpu_opt.n_dev_use;
802 if (hw_opt->nthreads_tot > 0 &&
803 hw_opt->nthreads_tmpi > hw_opt->nthreads_tot)
805 /* We have more GPUs than total threads requested.
806 * We choose to (later) generate a mismatch error,
807 * instead of launching more threads than requested.
809 hw_opt->nthreads_tmpi = hw_opt->nthreads_tot;
816 print_hw_opt(debug, hw_opt);
819 /* Asserting this simplifies the hardware resource division later
821 GMX_RELEASE_ASSERT(!(hw_opt->nthreads_omp_pme >= 1 && hw_opt->nthreads_omp <= 0),
822 "PME thread count should only be set when the normal thread count is also set");
825 /* Checks we can do when we know the cut-off scheme */
826 void check_and_update_hw_opt_2(gmx_hw_opt_t *hw_opt,
829 if (cutoff_scheme == ecutsGROUP)
831 /* We only have OpenMP support for PME only nodes */
832 if (hw_opt->nthreads_omp > 1)
834 gmx_fatal(FARGS, "OpenMP threads have been requested with cut-off scheme %s, but these are only supported with cut-off scheme %s",
835 ecutscheme_names[cutoff_scheme],
836 ecutscheme_names[ecutsVERLET]);
838 hw_opt->nthreads_omp = 1;
842 /* Checks we can do when we know the thread-MPI rank count */
843 void check_and_update_hw_opt_3(gmx_hw_opt_t gmx_unused *hw_opt)
846 GMX_RELEASE_ASSERT(hw_opt->nthreads_tmpi >= 1, "Must have at least one thread-MPI rank");
848 /* If the user set the total number of threads on the command line
849 * and did not specify the number of OpenMP threads, set the latter here.
851 if (hw_opt->nthreads_tot > 0 && hw_opt->nthreads_omp <= 0)
853 hw_opt->nthreads_omp = hw_opt->nthreads_tot/hw_opt->nthreads_tmpi;
855 if (!bHasOmpSupport && hw_opt->nthreads_omp > 1)
857 gmx_fatal(FARGS, "You (indirectly) asked for OpenMP threads by setting -nt > -ntmpi, but GROMACS was compiled without OpenMP support");
862 GMX_RELEASE_ASSERT(bHasOmpSupport || hw_opt->nthreads_omp == 1, "Without OpenMP support, only one thread per rank can be used");
864 /* We are done with updating nthreads_omp, we can set nthreads_omp_pme */
865 if (hw_opt->nthreads_omp_pme <= 0 && hw_opt->nthreads_omp > 0)
867 hw_opt->nthreads_omp_pme = hw_opt->nthreads_omp;
872 print_hw_opt(debug, hw_opt);