2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012,2013,2014,2015,2017 by the GROMACS development team.
5 * Copyright (c) 2018,2019,2020,2021, by the GROMACS development team, led by
6 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
7 * and including many others, as listed in the AUTHORS file in the
8 * top-level source directory and at http://www.gromacs.org.
10 * GROMACS is free software; you can redistribute it and/or
11 * modify it under the terms of the GNU Lesser General Public License
12 * as published by the Free Software Foundation; either version 2.1
13 * of the License, or (at your option) any later version.
15 * GROMACS is distributed in the hope that it will be useful,
16 * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18 * Lesser General Public License for more details.
20 * You should have received a copy of the GNU Lesser General Public
21 * License along with GROMACS; if not, see
22 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
23 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
25 * If you want to redistribute modifications to GROMACS, please
26 * consider that scientific software is very special. Version
27 * control is crucial - bugs must be traceable. We will be happy to
28 * consider code for inclusion in the official distribution, but
29 * derived work must not be called official GROMACS. Details are found
30 * in the README & COPYING files - if they are missing, get the
31 * official version at http://www.gromacs.org.
33 * To help us fund GROMACS development, we humbly ask that you cite
34 * the research papers on the package. Check out http://www.gromacs.org.
37 // FIXME: remove the "__" prefix in front of the group def when we move the
38 // nonbonded code into separate dir.
40 /*! \libinternal \defgroup __module_nbnxm Short-range non-bonded interaction module
41 * \ingroup group_mdrun
43 * \brief Computes forces and energies for short-range pair-interactions
44 * based on the Verlet algorithm. The algorithm uses pair-lists generated
45 * at fixed intervals as well as various flavors of pair interaction kernels
46 * implemented for a wide range of CPU and GPU architectures.
48 * The module includes support for flavors of Coulomb and Lennard-Jones interaction
49 * treatment implemented for a large range of SIMD instruction sets for CPU
50 * architectures as well as in CUDA and OpenCL for GPU architectures.
51 * Additionally there is a reference CPU non-SIMD and a reference CPU
52 * for GPU pair-list setup interaction kernel.
54 * The implementation of the kernels is based on the cluster non-bonded algorithm
55 * which in the code is referred to as the NxM algorithms ("nbnxm_" prefix);
56 * for details of the algorithm see DOI:10.1016/j.cpc.2013.06.003.
58 * Algorithmically, the non-bonded computation has two different modes:
59 * A "classical" mode: generate a list every nstlist steps containing at least
60 * all atom pairs up to a distance of rlistOuter and compute pair interactions
61 * for all pairs that are within the interaction cut-off.
62 * A "dynamic pruning" mode: generate an "outer-list" up to cut-off rlistOuter
63 * every nstlist steps and prune the outer-list using a cut-off of rlistInner
64 * every nstlistPrune steps to obtain a, smaller, "inner-list". This
65 * results in fewer interaction computations and allows for a larger nstlist.
66 * On a GPU, this dynamic pruning is performed in a rolling fashion, pruning
67 * only a sub-part of the list each (second) step. This way it can often
68 * overlap with integration and constraints on the CPU.
69 * Currently a simple heuristic determines which mode will be used.
71 * TODO: add a summary list and brief descriptions of the different submodules:
72 * search, CPU kernels, GPU glue code + kernels.
74 * \author Berk Hess <hess@kth.se>
75 * \author Szilárd Páll <pall.szilard@gmail.com>
76 * \author Mark Abraham <mark.j.abraham@gmail.com>
77 * \author Anca Hamuraru <anca@streamcomputing.eu>
78 * \author Teemu Virolainen <teemu@streamcomputing.eu>
79 * \author Dimitrios Karkoulis <dimitris.karkoulis@gmail.com>
81 * TODO: add more authors!
85 * \defgroup module_nbnxm Non-bonded pair interactions
86 * \ingroup group_mdrun
88 * Implements non-bonded pair interaction functionality for NxM atom clusters.
90 * This module provides methods to, very efficiently, compute non-bonded
91 * pair interactions on CPUs as well as accelerators. It also provides
92 * a method to construct the NxM atom-cluster pair-list required for
93 * computing these non-bonded iteractions.
96 /*! \libinternal \file
98 * \brief This file contains the public interface of the nbnxm module
99 * that implements the NxM atom cluster non-bonded algorithm to efficiently
100 * compute pair forces.
103 * \author Berk Hess <hess@kth.se>
104 * \author Szilárd Páll <pall.szilard@gmail.com>
107 * \ingroup module_nbnxm
111 #ifndef GMX_NBNXM_NBNXM_H
112 #define GMX_NBNXM_NBNXM_H
116 #include "gromacs/gpu_utils/devicebuffer_datatype.h"
117 #include "gromacs/math/vectypes.h"
118 #include "gromacs/mdtypes/locality.h"
119 #include "gromacs/utility/arrayref.h"
120 #include "gromacs/utility/enumerationhelpers.h"
121 #include "gromacs/utility/real.h"
123 struct DeviceInformation;
124 class FreeEnergyDispatch;
125 struct gmx_domdec_zones_t;
126 struct gmx_enerdata_t;
127 struct gmx_hw_info_t;
130 struct gmx_wallcycle;
131 struct interaction_const_t;
132 enum class LJCombinationRule;
133 struct nbnxn_atomdata_t;
134 struct nonbonded_verlet_t;
142 struct gmx_grppairener_t;
144 class GpuEventSynchronizer;
149 class ArrayRefWithPadding;
150 class DeviceStreamManager;
151 class ForceWithShiftForces;
152 class ListedForcesGpu;
159 class UpdateGroupsCog;
162 //! Namespace for non-bonded kernels
165 enum class KernelType;
167 /*! \brief Nbnxm electrostatic GPU kernel flavors.
169 * Types of electrostatics implementations available in the GPU non-bonded
170 * force kernels. These represent both the electrostatics types implemented
171 * by the kernels (cut-off, RF, and Ewald - a subset of what's defined in
172 * enums.h) as well as encode implementation details analytical/tabulated
173 * and single or twin cut-off (for Ewald kernels).
174 * Note that the cut-off and RF kernels have only analytical flavor and unlike
175 * in the CPU kernels, the tabulated kernels are ATM Ewald-only.
177 * The row-order of pointers to different electrostatic kernels defined in
178 * nbnxn_cuda.cu by the nb_*_kfunc_ptr function pointer table
179 * should match the order of enumerated types below.
181 enum class ElecType : int
183 Cut, //!< Plain cut-off
184 RF, //!< Reaction field
185 EwaldTab, //!< Tabulated Ewald with single cut-off
186 EwaldTabTwin, //!< Tabulated Ewald with twin cut-off
187 EwaldAna, //!< Analytical Ewald with single cut-off
188 EwaldAnaTwin, //!< Analytical Ewald with twin cut-off
189 Count //!< Number of valid values
192 //! Number of possible \ref ElecType values.
193 constexpr int c_numElecTypes = static_cast<int>(ElecType::Count);
195 /*! \brief Nbnxm VdW GPU kernel flavors.
197 * The enumerates values correspond to the LJ implementations in the GPU non-bonded
200 * The column-order of pointers to different electrostatic kernels defined in
201 * nbnxn_cuda_ocl.cpp/.cu by the nb_*_kfunc_ptr function pointer table
202 * should match the order of enumerated types below.
204 enum class VdwType : int
206 Cut, //!< Plain cut-off
207 CutCombGeom, //!< Cut-off with geometric combination rules
208 CutCombLB, //!< Cut-off with Lorentz-Berthelot combination rules
209 FSwitch, //!< Smooth force switch
210 PSwitch, //!< Smooth potential switch
211 EwaldGeom, //!< Ewald with geometric combination rules
212 EwaldLB, //!< Ewald with Lorentz-Berthelot combination rules
213 Count //!< Number of valid values
216 //! Number of possible \ref VdwType values.
217 constexpr int c_numVdwTypes = static_cast<int>(VdwType::Count);
219 /*! \brief Nonbonded NxN kernel types: plain C, CPU SIMD, GPU, GPU emulation */
220 enum class KernelType : int
231 /*! \brief Ewald exclusion types */
232 enum class EwaldExclusionType : int
240 /* \brief The non-bonded setup, also affects the pairlist construction kernel */
243 //! The non-bonded type, also affects the pairlist construction kernel
244 KernelType kernelType = KernelType::NotSet;
245 //! Ewald exclusion computation handling type, currently only used for CPU
246 EwaldExclusionType ewaldExclusionType = EwaldExclusionType::NotSet;
249 /*! \brief Return a string identifying the kernel type.
251 * \param [in] kernelType nonbonded kernel type, takes values from the nbnxn_kernel_type enum
252 * \returns a string identifying the kernel corresponding to the type passed as argument
254 const char* lookup_kernel_name(Nbnxm::KernelType kernelType);
258 /*! \brief Flag to tell the nonbonded kernels whether to clear the force output buffers */
266 * \brief Top-level non-bonded data structure for the Verlet-type cut-off scheme. */
267 struct nonbonded_verlet_t
270 //! Constructs an object from its components
271 nonbonded_verlet_t(std::unique_ptr<PairlistSets> pairlistSets,
272 std::unique_ptr<PairSearch> pairSearch,
273 std::unique_ptr<nbnxn_atomdata_t> nbat,
274 const Nbnxm::KernelSetup& kernelSetup,
276 gmx_wallcycle* wcycle);
278 ~nonbonded_verlet_t();
280 //! Returns whether a GPU is use for the non-bonded calculations
281 bool useGpu() const { return kernelSetup_.kernelType == Nbnxm::KernelType::Gpu8x8x8; }
283 //! Returns whether a GPU is emulated for the non-bonded calculations
284 bool emulateGpu() const
286 return kernelSetup_.kernelType == Nbnxm::KernelType::Cpu8x8x8_PlainC;
289 //! Return whether the pairlist is of simple, CPU type
290 bool pairlistIsSimple() const { return !useGpu() && !emulateGpu(); }
293 //! Returns the order of the local atoms on the grid
294 gmx::ArrayRef<const int> getLocalAtomOrder() const;
296 //! Sets the order of the local atoms to the order grid atom ordering
297 void setLocalAtomOrder() const;
299 //! Returns the index position of the atoms on the search grid
300 gmx::ArrayRef<const int> getGridIndices() const;
302 /*! \brief Constructs the pairlist for the given locality
304 * When there are no non-self exclusions, \p exclusions can be empty.
305 * Otherwise the number of lists in \p exclusions should match the number
306 * of atoms when not using DD, or the total number of atoms in the i-zones
309 * \param[in] iLocality The interaction locality: local or non-local
310 * \param[in] exclusions Lists of exclusions for every atom.
311 * \param[in] step Used to set the list creation step
312 * \param[in,out] nrnb Flop accounting struct, can be nullptr
314 void constructPairlist(gmx::InteractionLocality iLocality,
315 const gmx::ListOfLists<int>& exclusions,
319 //! Updates all the atom properties in Nbnxm
320 void setAtomProperties(gmx::ArrayRef<const int> atomTypes,
321 gmx::ArrayRef<const real> atomCharges,
322 gmx::ArrayRef<const int64_t> atomInfo) const;
324 /*!\brief Convert the coordinates to NBNXM format for the given locality.
326 * The API function for the transformation of the coordinates from one layout to another.
328 * \param[in] locality Whether coordinates for local or non-local atoms should be
329 * transformed. \param[in] coordinates Coordinates in plain rvec format to be transformed.
331 void convertCoordinates(gmx::AtomLocality locality, gmx::ArrayRef<const gmx::RVec> coordinates);
333 /*!\brief Convert the coordinates to NBNXM format on the GPU for the given locality
335 * The API function for the transformation of the coordinates from one layout to another in the GPU memory.
337 * \param[in] locality Whether coordinates for local or non-local atoms should be transformed.
338 * \param[in] d_x GPU coordinates buffer in plain rvec format to be transformed.
339 * \param[in] xReadyOnDevice Event synchronizer indicating that the coordinates are ready in the device memory.
341 void convertCoordinatesGpu(gmx::AtomLocality locality,
342 DeviceBuffer<gmx::RVec> d_x,
343 GpuEventSynchronizer* xReadyOnDevice);
345 //! Init for GPU version of setup coordinates in Nbnxm
346 void atomdata_init_copy_x_to_nbat_x_gpu() const;
348 //! Returns a reference to the pairlist sets
349 const PairlistSets& pairlistSets() const { return *pairlistSets_; }
351 //! Returns whether step is a dynamic list pruning step, for CPU lists
352 bool isDynamicPruningStepCpu(int64_t step) const;
354 //! Returns whether step is a dynamic list pruning step, for GPU lists
355 bool isDynamicPruningStepGpu(int64_t step) const;
357 //! Dispatches the dynamic pruning kernel for the given locality, for CPU lists
358 void dispatchPruneKernelCpu(gmx::InteractionLocality iLocality,
359 gmx::ArrayRef<const gmx::RVec> shift_vec) const;
361 //! Dispatches the dynamic pruning kernel for GPU lists
362 void dispatchPruneKernelGpu(int64_t step);
364 //! \brief Executes the non-bonded kernel of the GPU or launches it on the GPU
365 void dispatchNonbondedKernel(gmx::InteractionLocality iLocality,
366 const interaction_const_t& ic,
367 const gmx::StepWorkload& stepWork,
369 gmx::ArrayRef<const gmx::RVec> shiftvec,
370 gmx::ArrayRef<real> repulsionDispersionSR,
371 gmx::ArrayRef<real> CoulombSR,
374 //! Executes the non-bonded free-energy kernels, local + non-local, always runs on the CPU
375 void dispatchFreeEnergyKernels(const gmx::ArrayRefWithPadding<const gmx::RVec>& coords,
376 gmx::ForceWithShiftForces* forceWithShiftForces,
380 const interaction_const_t& ic,
381 gmx::ArrayRef<const gmx::RVec> shiftvec,
382 gmx::ArrayRef<const real> nbfp,
383 gmx::ArrayRef<const real> nbfp_grid,
384 gmx::ArrayRef<const real> chargeA,
385 gmx::ArrayRef<const real> chargeB,
386 gmx::ArrayRef<const int> typeA,
387 gmx::ArrayRef<const int> typeB,
389 gmx::ArrayRef<const real> lambda,
390 gmx_enerdata_t* enerd,
391 const gmx::StepWorkload& stepWork,
394 /*! \brief Add the forces stored in nbat to f, zeros the forces in nbat
395 * \param [in] locality Local or non-local
396 * \param [inout] force Force to be added to
398 void atomdata_add_nbat_f_to_f(gmx::AtomLocality locality, gmx::ArrayRef<gmx::RVec> force);
400 /*! \brief Get the number of atoms for a given locality
402 * \param [in] locality Local or non-local
403 * \returns The number of atoms for given locality
405 int getNumAtoms(gmx::AtomLocality locality) const;
407 //! Return the kernel setup
408 const Nbnxm::KernelSetup& kernelSetup() const { return kernelSetup_; }
410 //! Returns the outer radius for the pair list
411 real pairlistInnerRadius() const;
413 //! Returns the outer radius for the pair list
414 real pairlistOuterRadius() const;
416 //! Changes the pair-list outer and inner radius
417 void changePairlistRadii(real rlistOuter, real rlistInner) const;
419 //! Set up internal flags that indicate what type of short-range work there is.
420 void setupGpuShortRangeWork(const gmx::ListedForcesGpu* listedForcesGpu,
421 gmx::InteractionLocality iLocality) const;
423 void setupFepThreadedForceBuffer(int numAtomsForce);
425 // TODO: Make all data members private
426 //! All data related to the pair lists
427 std::unique_ptr<PairlistSets> pairlistSets_;
428 //! Working data for constructing the pairlists
429 std::unique_ptr<PairSearch> pairSearch_;
431 std::unique_ptr<nbnxn_atomdata_t> nbat;
434 //! The non-bonded setup, also affects the pairlist construction kernel
435 Nbnxm::KernelSetup kernelSetup_;
437 //! \brief Pointer to wallcycle structure.
438 gmx_wallcycle* wcycle_;
440 //! \brief The non-bonded free-energy kernel dispatcher
441 std::unique_ptr<FreeEnergyDispatch> freeEnergyDispatch_;
444 //! GPU Nbnxm data, only used with a physical GPU (TODO: use unique_ptr)
451 /*! \brief Creates an Nbnxm object */
452 std::unique_ptr<nonbonded_verlet_t> init_nb_verlet(const gmx::MDLogger& mdlog,
453 const t_inputrec& inputrec,
454 const t_forcerec& forcerec,
455 const t_commrec* commrec,
456 const gmx_hw_info_t& hardwareInfo,
457 bool useGpuForNonbonded,
458 const gmx::DeviceStreamManager* deviceStreamManager,
459 const gmx_mtop_t& mtop,
461 gmx_wallcycle* wcycle);
465 /*! \brief Put the atoms on the pair search grid.
467 * Only atoms with indices wihtin \p atomRange in x are put on the grid.
468 * When \p updateGroupsCog != nullptr, atoms are put on the grid
469 * based on the center of geometry of the group they belong to.
470 * Atoms or COGs of groups should be within the bounding box provided,
471 * this is checked in debug builds when not using update groups.
472 * The atom density is used to determine the grid size when \p gridIndex = 0.
473 * When \p atomDensity <= 0, the density is determined from atomEnd-atomStart
474 * and the bounding box corners.
475 * With domain decomposition, part of the atoms might have migrated,
476 * but have not been removed yet. This count is given by \p numAtomsMoved.
477 * When \p move[i] < 0 particle i has migrated and will not be put on the grid.
479 * \param[in,out] nb_verlet The non-bonded object
480 * \param[in] box Box used for periodic distance calculations
481 * \param[in] gridIndex The index of the grid to spread to, always 0 except with test particle insertion
482 * \param[in] lowerCorner Atom groups to be gridded should have coordinates >= this corner
483 * \param[in] upperCorner Atom groups to be gridded should have coordinates <= this corner
484 * \param[in] updateGroupsCog Centers of geometry for update groups, pass nullptr when not using update groups
485 * \param[in] atomRange Range of atoms to grid
486 * \param[in] atomDensity An estimate of the atom density, used for peformance optimization and only with \p gridIndex = 0
487 * \param[in] atomInfo Atom information flags
488 * \param[in] x Coordinates for atoms to grid
489 * \param[in] numAtomsMoved The number of atoms that will move to another domain, pass 0 without DD
490 * \param[in] move Move flags for atoms, pass nullptr without DD
492 void nbnxn_put_on_grid(nonbonded_verlet_t* nb_verlet,
495 const rvec lowerCorner,
496 const rvec upperCorner,
497 const gmx::UpdateGroupsCog* updateGroupsCog,
498 gmx::Range<int> atomRange,
500 gmx::ArrayRef<const int64_t> atomInfo,
501 gmx::ArrayRef<const gmx::RVec> x,
505 /*! \brief As nbnxn_put_on_grid, but for the non-local atoms
507 * with domain decomposition. Should be called after calling
508 * nbnxn_search_put_on_grid for the local atoms / home zone.
510 void nbnxn_put_on_grid_nonlocal(nonbonded_verlet_t* nb_verlet,
511 const struct gmx_domdec_zones_t* zones,
512 gmx::ArrayRef<const int64_t> atomInfo,
513 gmx::ArrayRef<const gmx::RVec> x);
515 /*! \brief Check if GROMACS has been built with GPU support.
517 * \param[in] error Pointer to error string or nullptr.
518 * \todo Move this to NB module once it exists.
520 bool buildSupportsNonbondedOnGpu(std::string* error);
522 #endif // GMX_NBNXN_NBNXM_H