2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2019, 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 Common functions for the different NBNXN GPU implementations.
38 * \author Berk Hess <hess@kth.se>
40 * \ingroup module_nbnxm
45 #include "gromacs/domdec/domdec.h"
46 #include "gromacs/domdec/domdec_struct.h"
47 #include "gromacs/hardware/hw_info.h"
48 #include "gromacs/mdlib/gmx_omp_nthreads.h"
49 #include "gromacs/mdtypes/commrec.h"
50 #include "gromacs/mdtypes/forcerec.h"
51 #include "gromacs/mdtypes/inputrec.h"
52 #include "gromacs/nbnxm/atomdata.h"
53 #include "gromacs/nbnxm/gpu_data_mgmt.h"
54 #include "gromacs/nbnxm/nbnxm.h"
55 #include "gromacs/nbnxm/nbnxm_geometry.h"
56 #include "gromacs/nbnxm/nbnxm_simd.h"
57 #include "gromacs/nbnxm/pairlist_tuning.h"
58 #include "gromacs/simd/simd.h"
59 #include "gromacs/utility/fatalerror.h"
60 #include "gromacs/utility/logger.h"
64 /*! \brief Returns whether CPU SIMD support exists for the given inputrec
66 * If the return value is FALSE and fplog/cr != NULL, prints a fallback
67 * message to fplog/stderr.
69 static gmx_bool nbnxn_simd_supported(const gmx::MDLogger &mdlog,
72 if (ir->vdwtype == evdwPME && ir->ljpme_combination_rule == eljpmeLB)
74 /* LJ PME with LB combination rule does 7 mesh operations.
75 * This so slow that we don't compile SIMD non-bonded kernels
77 GMX_LOG(mdlog.warning).asParagraph().appendText("LJ-PME with Lorentz-Berthelot is not supported with SIMD kernels, falling back to plain C kernels");
84 /*! \brief Returns the most suitable CPU kernel type and Ewald handling */
85 static void pick_nbnxn_kernel_cpu(const t_inputrec gmx_unused *ir,
88 const gmx_hw_info_t gmx_unused &hardwareInfo)
90 *kernel_type = nbnxnk4x4_PlainC;
91 *ewald_excl = ewaldexclTable;
95 #ifdef GMX_NBNXN_SIMD_4XN
96 *kernel_type = nbnxnk4xN_SIMD_4xN;
98 #ifdef GMX_NBNXN_SIMD_2XNN
99 *kernel_type = nbnxnk4xN_SIMD_2xNN;
102 #if defined GMX_NBNXN_SIMD_2XNN && defined GMX_NBNXN_SIMD_4XN
103 /* We need to choose if we want 2x(N+N) or 4xN kernels.
104 * This is based on the SIMD acceleration choice and CPU information
105 * detected at runtime.
107 * 4xN calculates more (zero) interactions, but has less pair-search
108 * work and much better kernel instruction scheduling.
110 * Up till now we have only seen that on Intel Sandy/Ivy Bridge,
111 * which doesn't have FMA, both the analytical and tabulated Ewald
112 * kernels have similar pair rates for 4x8 and 2x(4+4), so we choose
113 * 2x(4+4) because it results in significantly fewer pairs.
114 * For RF, the raw pair rate of the 4x8 kernel is higher than 2x(4+4),
115 * 10% with HT, 50% without HT. As we currently don't detect the actual
116 * use of HT, use 4x8 to avoid a potential performance hit.
117 * On Intel Haswell 4x8 is always faster.
119 *kernel_type = nbnxnk4xN_SIMD_4xN;
121 #if !GMX_SIMD_HAVE_FMA
122 if (EEL_PME_EWALD(ir->coulombtype) ||
123 EVDW_PME(ir->vdwtype))
125 /* We have Ewald kernels without FMA (Intel Sandy/Ivy Bridge).
126 * There are enough instructions to make 2x(4+4) efficient.
128 *kernel_type = nbnxnk4xN_SIMD_2xNN;
131 if (hardwareInfo.haveAmdZenCpu)
133 /* One 256-bit FMA per cycle makes 2xNN faster */
134 *kernel_type = nbnxnk4xN_SIMD_2xNN;
136 #endif /* GMX_NBNXN_SIMD_2XNN && GMX_NBNXN_SIMD_4XN */
139 if (getenv("GMX_NBNXN_SIMD_4XN") != nullptr)
141 #ifdef GMX_NBNXN_SIMD_4XN
142 *kernel_type = nbnxnk4xN_SIMD_4xN;
144 gmx_fatal(FARGS, "SIMD 4xN kernels requested, but GROMACS has been compiled without support for these kernels");
147 if (getenv("GMX_NBNXN_SIMD_2XNN") != nullptr)
149 #ifdef GMX_NBNXN_SIMD_2XNN
150 *kernel_type = nbnxnk4xN_SIMD_2xNN;
152 gmx_fatal(FARGS, "SIMD 2x(N+N) kernels requested, but GROMACS has been compiled without support for these kernels");
156 /* Analytical Ewald exclusion correction is only an option in
158 * Since table lookup's don't parallelize with SIMD, analytical
159 * will probably always be faster for a SIMD width of 8 or more.
160 * With FMA analytical is sometimes faster for a width if 4 as well.
161 * In single precision, this is faster on Bulldozer.
163 #if GMX_SIMD_REAL_WIDTH >= 8 || \
164 (GMX_SIMD_REAL_WIDTH >= 4 && GMX_SIMD_HAVE_FMA && !GMX_DOUBLE)
165 /* On AMD Zen, tabulated Ewald kernels are faster on all 4 combinations
166 * of single or double precision and 128 or 256-bit AVX2.
168 if (!hardwareInfo.haveAmdZenCpu)
170 *ewald_excl = ewaldexclAnalytical;
173 if (getenv("GMX_NBNXN_EWALD_TABLE") != nullptr)
175 *ewald_excl = ewaldexclTable;
177 if (getenv("GMX_NBNXN_EWALD_ANALYTICAL") != nullptr)
179 *ewald_excl = ewaldexclAnalytical;
186 const char *lookup_nbnxn_kernel_name(int kernel_type)
188 const char *returnvalue = nullptr;
192 returnvalue = "not set";
194 case nbnxnk4x4_PlainC:
195 returnvalue = "plain C";
197 case nbnxnk4xN_SIMD_4xN:
198 case nbnxnk4xN_SIMD_2xNN:
200 returnvalue = "SIMD";
202 returnvalue = "not available";
205 case nbnxnk8x8x8_GPU: returnvalue = "GPU"; break;
206 case nbnxnk8x8x8_PlainC: returnvalue = "plain C"; break;
210 gmx_fatal(FARGS, "Illegal kernel type selected");
215 /*! \brief Returns the most suitable kernel type and Ewald handling */
216 static void pick_nbnxn_kernel(const gmx::MDLogger &mdlog,
217 gmx_bool use_simd_kernels,
218 const gmx_hw_info_t &hardwareInfo,
220 EmulateGpuNonbonded emulateGpu,
221 const t_inputrec *ir,
224 gmx_bool bDoNonbonded)
226 GMX_RELEASE_ASSERT(kernel_type, "Need a valid kernel_type pointer");
228 *kernel_type = nbnxnkNotSet;
229 *ewald_excl = ewaldexclTable;
231 if (emulateGpu == EmulateGpuNonbonded::Yes)
233 *kernel_type = nbnxnk8x8x8_PlainC;
237 GMX_LOG(mdlog.warning).asParagraph().appendText("Emulating a GPU run on the CPU (slow)");
242 *kernel_type = nbnxnk8x8x8_GPU;
245 if (*kernel_type == nbnxnkNotSet)
247 if (use_simd_kernels &&
248 nbnxn_simd_supported(mdlog, ir))
250 pick_nbnxn_kernel_cpu(ir, kernel_type, ewald_excl, hardwareInfo);
254 *kernel_type = nbnxnk4x4_PlainC;
260 GMX_LOG(mdlog.info).asParagraph().appendTextFormatted(
261 "Using %s %dx%d nonbonded short-range kernels",
262 lookup_nbnxn_kernel_name(*kernel_type),
263 nbnxn_kernel_to_cluster_i_size(*kernel_type),
264 nbnxn_kernel_to_cluster_j_size(*kernel_type));
266 if (nbnxnk4x4_PlainC == *kernel_type ||
267 nbnxnk8x8x8_PlainC == *kernel_type)
269 GMX_LOG(mdlog.warning).asParagraph().appendTextFormatted(
270 "WARNING: Using the slow %s kernels. This should\n"
271 "not happen during routine usage on supported platforms.",
272 lookup_nbnxn_kernel_name(*kernel_type));
277 void init_nb_verlet(const gmx::MDLogger &mdlog,
278 nonbonded_verlet_t **nb_verlet,
279 gmx_bool bFEP_NonBonded,
280 const t_inputrec *ir,
281 const t_forcerec *fr,
283 const gmx_hw_info_t &hardwareInfo,
284 const gmx_device_info_t *deviceInfo,
285 const gmx_mtop_t *mtop,
288 nonbonded_verlet_t *nbv;
291 nbv = new nonbonded_verlet_t();
293 nbv->emulateGpu = ((getenv("GMX_EMULATE_GPU") != nullptr) ? EmulateGpuNonbonded::Yes : EmulateGpuNonbonded::No);
294 nbv->bUseGPU = deviceInfo != nullptr;
296 GMX_RELEASE_ASSERT(!(nbv->emulateGpu == EmulateGpuNonbonded::Yes && nbv->bUseGPU), "When GPU emulation is active, there cannot be a GPU assignment");
299 nbv->min_ci_balanced = 0;
301 nbv->ngrp = (DOMAINDECOMP(cr) ? 2 : 1);
302 for (int i = 0; i < nbv->ngrp; i++)
304 nbv->grp[i].nbl_lists.nnbl = 0;
305 nbv->grp[i].kernel_type = nbnxnkNotSet;
307 if (i == 0) /* local */
309 pick_nbnxn_kernel(mdlog, fr->use_simd_kernels, hardwareInfo,
310 nbv->bUseGPU, nbv->emulateGpu, ir,
311 &nbv->grp[i].kernel_type,
312 &nbv->grp[i].ewald_excl,
317 /* Use the same kernel for local and non-local interactions */
318 nbv->grp[i].kernel_type = nbv->grp[0].kernel_type;
319 nbv->grp[i].ewald_excl = nbv->grp[0].ewald_excl;
323 nbv->listParams = std::make_unique<NbnxnListParameters>(ir->rlist);
324 setupDynamicPairlistPruning(mdlog, ir, mtop, box, nbv->grp[0].kernel_type, fr->ic,
325 nbv->listParams.get());
327 nbv->nbs = std::make_unique<nbnxn_search>(DOMAINDECOMP(cr) ? &cr->dd->nc : nullptr,
328 DOMAINDECOMP(cr) ? domdec_zones(cr->dd) : nullptr,
330 gmx_omp_nthreads_get(emntPairsearch));
332 for (int i = 0; i < nbv->ngrp; i++)
334 nbnxn_init_pairlist_set(&nbv->grp[i].nbl_lists,
335 nbnxn_kernel_pairlist_simple(nbv->grp[i].kernel_type),
336 /* 8x8x8 "non-simple" lists are ATM always combined */
337 !nbnxn_kernel_pairlist_simple(nbv->grp[i].kernel_type));
340 int enbnxninitcombrule;
341 if (fr->ic->vdwtype == evdwCUT &&
342 (fr->ic->vdw_modifier == eintmodNONE ||
343 fr->ic->vdw_modifier == eintmodPOTSHIFT) &&
344 getenv("GMX_NO_LJ_COMB_RULE") == nullptr)
346 /* Plain LJ cut-off: we can optimize with combination rules */
347 enbnxninitcombrule = enbnxninitcombruleDETECT;
349 else if (fr->ic->vdwtype == evdwPME)
351 /* LJ-PME: we need to use a combination rule for the grid */
352 if (fr->ljpme_combination_rule == eljpmeGEOM)
354 enbnxninitcombrule = enbnxninitcombruleGEOM;
358 enbnxninitcombrule = enbnxninitcombruleLB;
363 /* We use a full combination matrix: no rule required */
364 enbnxninitcombrule = enbnxninitcombruleNONE;
367 nbv->nbat = new nbnxn_atomdata_t(nbv->bUseGPU ? gmx::PinningPolicy::PinnedIfSupported : gmx::PinningPolicy::CannotBePinned);
368 int mimimumNumEnergyGroupNonbonded = ir->opts.ngener;
369 if (ir->opts.ngener - ir->nwall == 1)
371 /* We have only one non-wall energy group, we do not need energy group
372 * support in the non-bondeds kernels, since all non-bonded energy
373 * contributions go to the first element of the energy group matrix.
375 mimimumNumEnergyGroupNonbonded = 1;
377 bool bSimpleList = nbnxn_kernel_pairlist_simple(nbv->grp[0].kernel_type);
378 nbnxn_atomdata_init(mdlog,
380 nbv->grp[0].kernel_type,
383 mimimumNumEnergyGroupNonbonded,
384 bSimpleList ? gmx_omp_nthreads_get(emntNonbonded) : 1);
388 /* init the NxN GPU data; the last argument tells whether we'll have
389 * both local and non-local NB calculation on GPU */
390 nbnxn_gpu_init(&nbv->gpu_nbv,
393 nbv->listParams.get(),
398 if ((env = getenv("GMX_NB_MIN_CI")) != nullptr)
402 nbv->min_ci_balanced = strtol(env, &end, 10);
403 if (!end || (*end != 0) || nbv->min_ci_balanced < 0)
405 gmx_fatal(FARGS, "Invalid value passed in GMX_NB_MIN_CI=%s, non-negative integer required", env);
410 fprintf(debug, "Neighbor-list balancing parameter: %d (passed as env. var.)\n",
411 nbv->min_ci_balanced);
416 nbv->min_ci_balanced = nbnxn_gpu_min_ci_balanced(nbv->gpu_nbv);
419 fprintf(debug, "Neighbor-list balancing parameter: %d (auto-adjusted to the number of GPU multi-processors)\n",
420 nbv->min_ci_balanced);