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 struct gmx_domdec_zones_t;
125 struct gmx_enerdata_t;
126 struct gmx_hw_info_t;
129 struct gmx_wallcycle;
130 struct interaction_const_t;
131 enum class LJCombinationRule;
132 struct nbnxn_atomdata_t;
133 struct nonbonded_verlet_t;
141 struct gmx_grppairener_t;
143 class GpuEventSynchronizer;
147 class DeviceStreamManager;
148 class ForceWithShiftForces;
156 class UpdateGroupsCog;
159 //! Namespace for non-bonded kernels
162 enum class KernelType;
164 /*! \brief Nbnxm electrostatic GPU kernel flavors.
166 * Types of electrostatics implementations available in the GPU non-bonded
167 * force kernels. These represent both the electrostatics types implemented
168 * by the kernels (cut-off, RF, and Ewald - a subset of what's defined in
169 * enums.h) as well as encode implementation details analytical/tabulated
170 * and single or twin cut-off (for Ewald kernels).
171 * Note that the cut-off and RF kernels have only analytical flavor and unlike
172 * in the CPU kernels, the tabulated kernels are ATM Ewald-only.
174 * The row-order of pointers to different electrostatic kernels defined in
175 * nbnxn_cuda.cu by the nb_*_kfunc_ptr function pointer table
176 * should match the order of enumerated types below.
178 enum class ElecType : int
180 Cut, //!< Plain cut-off
181 RF, //!< Reaction field
182 EwaldTab, //!< Tabulated Ewald with single cut-off
183 EwaldTabTwin, //!< Tabulated Ewald with twin cut-off
184 EwaldAna, //!< Analytical Ewald with single cut-off
185 EwaldAnaTwin, //!< Analytical Ewald with twin cut-off
186 Count //!< Number of valid values
189 //! Number of possible \ref ElecType values.
190 constexpr int c_numElecTypes = static_cast<int>(ElecType::Count);
192 /*! \brief Nbnxm VdW GPU kernel flavors.
194 * The enumerates values correspond to the LJ implementations in the GPU non-bonded
197 * The column-order of pointers to different electrostatic kernels defined in
198 * nbnxn_cuda_ocl.cpp/.cu by the nb_*_kfunc_ptr function pointer table
199 * should match the order of enumerated types below.
201 enum class VdwType : int
203 Cut, //!< Plain cut-off
204 CutCombGeom, //!< Cut-off with geometric combination rules
205 CutCombLB, //!< Cut-off with Lorentz-Berthelot combination rules
206 FSwitch, //!< Smooth force switch
207 PSwitch, //!< Smooth potential switch
208 EwaldGeom, //!< Ewald with geometric combination rules
209 EwaldLB, //!< Ewald with Lorentz-Berthelot combination rules
210 Count //!< Number of valid values
213 //! Number of possible \ref VdwType values.
214 constexpr int c_numVdwTypes = static_cast<int>(VdwType::Count);
216 /*! \brief Nonbonded NxN kernel types: plain C, CPU SIMD, GPU, GPU emulation */
217 enum class KernelType : int
228 /*! \brief Ewald exclusion types */
229 enum class EwaldExclusionType : int
237 /* \brief The non-bonded setup, also affects the pairlist construction kernel */
240 //! The non-bonded type, also affects the pairlist construction kernel
241 KernelType kernelType = KernelType::NotSet;
242 //! Ewald exclusion computation handling type, currently only used for CPU
243 EwaldExclusionType ewaldExclusionType = EwaldExclusionType::NotSet;
246 /*! \brief Return a string identifying the kernel type.
248 * \param [in] kernelType nonbonded kernel type, takes values from the nbnxn_kernel_type enum
249 * \returns a string identifying the kernel corresponding to the type passed as argument
251 const char* lookup_kernel_name(Nbnxm::KernelType kernelType);
255 /*! \brief Flag to tell the nonbonded kernels whether to clear the force output buffers */
263 * \brief Top-level non-bonded data structure for the Verlet-type cut-off scheme. */
264 struct nonbonded_verlet_t
267 //! Constructs an object from its components
268 nonbonded_verlet_t(std::unique_ptr<PairlistSets> pairlistSets,
269 std::unique_ptr<PairSearch> pairSearch,
270 std::unique_ptr<nbnxn_atomdata_t> nbat,
271 const Nbnxm::KernelSetup& kernelSetup,
273 gmx_wallcycle* wcycle);
275 ~nonbonded_verlet_t();
277 //! Returns whether a GPU is use for the non-bonded calculations
278 bool useGpu() const { return kernelSetup_.kernelType == Nbnxm::KernelType::Gpu8x8x8; }
280 //! Returns whether a GPU is emulated for the non-bonded calculations
281 bool emulateGpu() const
283 return kernelSetup_.kernelType == Nbnxm::KernelType::Cpu8x8x8_PlainC;
286 //! Return whether the pairlist is of simple, CPU type
287 bool pairlistIsSimple() const { return !useGpu() && !emulateGpu(); }
290 //! Returns the order of the local atoms on the grid
291 gmx::ArrayRef<const int> getLocalAtomOrder() const;
293 //! Sets the order of the local atoms to the order grid atom ordering
294 void setLocalAtomOrder() const;
296 //! Returns the index position of the atoms on the search grid
297 gmx::ArrayRef<const int> getGridIndices() const;
299 /*! \brief Constructs the pairlist for the given locality
301 * When there are no non-self exclusions, \p exclusions can be empty.
302 * Otherwise the number of lists in \p exclusions should match the number
303 * of atoms when not using DD, or the total number of atoms in the i-zones
306 * \param[in] iLocality The interaction locality: local or non-local
307 * \param[in] exclusions Lists of exclusions for every atom.
308 * \param[in] step Used to set the list creation step
309 * \param[in,out] nrnb Flop accounting struct, can be nullptr
311 void constructPairlist(gmx::InteractionLocality iLocality,
312 const gmx::ListOfLists<int>& exclusions,
316 //! Updates all the atom properties in Nbnxm
317 void setAtomProperties(gmx::ArrayRef<const int> atomTypes,
318 gmx::ArrayRef<const real> atomCharges,
319 gmx::ArrayRef<const int> atomInfo) const;
321 /*!\brief Convert the coordinates to NBNXM format for the given locality.
323 * The API function for the transformation of the coordinates from one layout to another.
325 * \param[in] locality Whether coordinates for local or non-local atoms should be
326 * transformed. \param[in] coordinates Coordinates in plain rvec format to be transformed.
328 void convertCoordinates(gmx::AtomLocality locality, gmx::ArrayRef<const gmx::RVec> coordinates);
330 /*!\brief Convert the coordinates to NBNXM format on the GPU for the given locality
332 * The API function for the transformation of the coordinates from one layout to another in the GPU memory.
334 * \param[in] locality Whether coordinates for local or non-local atoms should be transformed.
335 * \param[in] d_x GPU coordinates buffer in plain rvec format to be transformed.
336 * \param[in] xReadyOnDevice Event synchronizer indicating that the coordinates are ready in the device memory.
338 void convertCoordinatesGpu(gmx::AtomLocality locality,
339 DeviceBuffer<gmx::RVec> d_x,
340 GpuEventSynchronizer* xReadyOnDevice);
342 //! Init for GPU version of setup coordinates in Nbnxm
343 void atomdata_init_copy_x_to_nbat_x_gpu() const;
345 //! Returns a reference to the pairlist sets
346 const PairlistSets& pairlistSets() const { return *pairlistSets_; }
348 //! Returns whether step is a dynamic list pruning step, for CPU lists
349 bool isDynamicPruningStepCpu(int64_t step) const;
351 //! Returns whether step is a dynamic list pruning step, for GPU lists
352 bool isDynamicPruningStepGpu(int64_t step) const;
354 //! Dispatches the dynamic pruning kernel for the given locality, for CPU lists
355 void dispatchPruneKernelCpu(gmx::InteractionLocality iLocality,
356 gmx::ArrayRef<const gmx::RVec> shift_vec) const;
358 //! Dispatches the dynamic pruning kernel for GPU lists
359 void dispatchPruneKernelGpu(int64_t step);
361 //! \brief Executes the non-bonded kernel of the GPU or launches it on the GPU
362 void dispatchNonbondedKernel(gmx::InteractionLocality iLocality,
363 const interaction_const_t& ic,
364 const gmx::StepWorkload& stepWork,
366 gmx::ArrayRef<const gmx::RVec> shiftvec,
367 gmx::ArrayRef<real> repulsionDispersionSR,
368 gmx::ArrayRef<real> CoulombSR,
371 //! Executes the non-bonded free-energy kernel, always runs on the CPU
372 void dispatchFreeEnergyKernel(gmx::InteractionLocality iLocality,
373 gmx::ArrayRef<const gmx::RVec> coords,
374 gmx::ForceWithShiftForces* forceWithShiftForces,
378 const interaction_const_t& ic,
379 gmx::ArrayRef<const gmx::RVec> shiftvec,
380 gmx::ArrayRef<const real> nbfp,
381 gmx::ArrayRef<const real> nbfp_grid,
382 gmx::ArrayRef<const real> chargeA,
383 gmx::ArrayRef<const real> chargeB,
384 gmx::ArrayRef<const int> typeA,
385 gmx::ArrayRef<const int> typeB,
387 gmx::ArrayRef<const real> lambda,
388 gmx_enerdata_t* enerd,
389 const gmx::StepWorkload& stepWork,
392 /*! \brief Add the forces stored in nbat to f, zeros the forces in nbat
393 * \param [in] locality Local or non-local
394 * \param [inout] force Force to be added to
396 void atomdata_add_nbat_f_to_f(gmx::AtomLocality locality, gmx::ArrayRef<gmx::RVec> force);
398 /*! \brief Get the number of atoms for a given locality
400 * \param [in] locality Local or non-local
401 * \returns The number of atoms for given locality
403 int getNumAtoms(gmx::AtomLocality locality) const;
405 //! Return the kernel setup
406 const Nbnxm::KernelSetup& kernelSetup() const { return kernelSetup_; }
408 //! Returns the outer radius for the pair list
409 real pairlistInnerRadius() const;
411 //! Returns the outer radius for the pair list
412 real pairlistOuterRadius() const;
414 //! Changes the pair-list outer and inner radius
415 void changePairlistRadii(real rlistOuter, real rlistInner) const;
417 //! Set up internal flags that indicate what type of short-range work there is.
418 void setupGpuShortRangeWork(const gmx::GpuBonded* gpuBonded, gmx::InteractionLocality iLocality) const;
420 // TODO: Make all data members private
421 //! All data related to the pair lists
422 std::unique_ptr<PairlistSets> pairlistSets_;
423 //! Working data for constructing the pairlists
424 std::unique_ptr<PairSearch> pairSearch_;
426 std::unique_ptr<nbnxn_atomdata_t> nbat;
429 //! The non-bonded setup, also affects the pairlist construction kernel
430 Nbnxm::KernelSetup kernelSetup_;
431 //! \brief Pointer to wallcycle structure.
432 gmx_wallcycle* wcycle_;
433 //! Temporary array for storing foreign lambda group pair energies
434 std::unique_ptr<gmx_grppairener_t> foreignEnergyGroups_;
437 //! GPU Nbnxm data, only used with a physical GPU (TODO: use unique_ptr)
444 /*! \brief Creates an Nbnxm object */
445 std::unique_ptr<nonbonded_verlet_t> init_nb_verlet(const gmx::MDLogger& mdlog,
446 const t_inputrec& inputrec,
447 const t_forcerec& forcerec,
448 const t_commrec* commrec,
449 const gmx_hw_info_t& hardwareInfo,
450 bool useGpuForNonbonded,
451 const gmx::DeviceStreamManager* deviceStreamManager,
452 const gmx_mtop_t& mtop,
454 gmx_wallcycle* wcycle);
458 /*! \brief Put the atoms on the pair search grid.
460 * Only atoms with indices wihtin \p atomRange in x are put on the grid.
461 * When \p updateGroupsCog != nullptr, atoms are put on the grid
462 * based on the center of geometry of the group they belong to.
463 * Atoms or COGs of groups should be within the bounding box provided,
464 * this is checked in debug builds when not using update groups.
465 * The atom density is used to determine the grid size when \p gridIndex = 0.
466 * When \p atomDensity <= 0, the density is determined from atomEnd-atomStart
467 * and the bounding box corners.
468 * With domain decomposition, part of the atoms might have migrated,
469 * but have not been removed yet. This count is given by \p numAtomsMoved.
470 * When \p move[i] < 0 particle i has migrated and will not be put on the grid.
472 * \param[in,out] nb_verlet The non-bonded object
473 * \param[in] box Box used for periodic distance calculations
474 * \param[in] gridIndex The index of the grid to spread to, always 0 except with test particle insertion
475 * \param[in] lowerCorner Atom groups to be gridded should have coordinates >= this corner
476 * \param[in] upperCorner Atom groups to be gridded should have coordinates <= this corner
477 * \param[in] updateGroupsCog Centers of geometry for update groups, pass nullptr when not using update groups
478 * \param[in] atomRange Range of atoms to grid
479 * \param[in] atomDensity An estimate of the atom density, used for peformance optimization and only with \p gridIndex = 0
480 * \param[in] atomInfo Atom information flags
481 * \param[in] x Coordinates for atoms to grid
482 * \param[in] numAtomsMoved The number of atoms that will move to another domain, pass 0 without DD
483 * \param[in] move Move flags for atoms, pass nullptr without DD
485 void nbnxn_put_on_grid(nonbonded_verlet_t* nb_verlet,
488 const rvec lowerCorner,
489 const rvec upperCorner,
490 const gmx::UpdateGroupsCog* updateGroupsCog,
491 gmx::Range<int> atomRange,
493 gmx::ArrayRef<const int> atomInfo,
494 gmx::ArrayRef<const gmx::RVec> x,
498 /*! \brief As nbnxn_put_on_grid, but for the non-local atoms
500 * with domain decomposition. Should be called after calling
501 * nbnxn_search_put_on_grid for the local atoms / home zone.
503 void nbnxn_put_on_grid_nonlocal(nonbonded_verlet_t* nb_verlet,
504 const struct gmx_domdec_zones_t* zones,
505 gmx::ArrayRef<const int> atomInfo,
506 gmx::ArrayRef<const gmx::RVec> x);
508 /*! \brief Check if GROMACS has been built with GPU support.
510 * \param[in] error Pointer to error string or nullptr.
511 * \todo Move this to NB module once it exists.
513 bool buildSupportsNonbondedOnGpu(std::string* error);
515 #endif // GMX_NBNXN_NBNXM_H