Gromacs
2025-dev-20240812-545ca5b
|
#include <gromacs/nbnxm/nbnxm.h>
Top-level non-bonded data structure for the Verlet-type cut-off scheme.
Public Member Functions | |
nonbonded_verlet_t (std::unique_ptr< PairlistSets > pairlistSets, std::unique_ptr< PairSearch > pairSearch, std::unique_ptr< nbnxn_atomdata_t > nbat, const Nbnxm::KernelSetup &kernelSetup, std::unique_ptr< ExclusionChecker > exclusionChecker, NbnxmGpu *gpu_nbv, gmx_wallcycle *wcycle) | |
Constructs an object from its components. More... | |
nonbonded_verlet_t (std::unique_ptr< PairlistSets > pairlistSets, std::unique_ptr< PairSearch > pairSearch, std::unique_ptr< nbnxn_atomdata_t > nbat, const Nbnxm::KernelSetup &kernelSetup, NbnxmGpu *gpu_nbv) | |
Constructs an object from its, minimal, components. More... | |
bool | useGpu () const |
Returns whether a GPU is use for the non-bonded calculations. | |
bool | emulateGpu () const |
Returns whether a GPU is emulated for the non-bonded calculations. | |
bool | pairlistIsSimple () const |
Return whether the pairlist is of simple, CPU type. | |
void | putAtomsOnGrid (const matrix box, int gridIndex, const gmx::RVec &lowerCorner, const gmx::RVec &upperCorner, const gmx::UpdateGroupsCog *updateGroupsCog, gmx::Range< int > atomRange, real atomDensity, gmx::ArrayRef< const int32_t > atomInfo, gmx::ArrayRef< const gmx::RVec > x, int numAtomsMoved, const int *move) |
Put the atoms on the pair search grid. More... | |
gmx::ArrayRef< const int > | getLocalAtomOrder () const |
Returns the order of the local atoms on the grid. | |
void | setLocalAtomOrder () const |
Sets the order of the local atoms to the order grid atom ordering. | |
gmx::ArrayRef< const int > | getGridIndices () const |
Returns the index position of the atoms on the search grid. | |
void | constructPairlist (gmx::InteractionLocality iLocality, const gmx::ListOfLists< int > &exclusions, int64_t step, t_nrnb *nrnb) const |
Constructs the pairlist for the given locality. More... | |
void | setAtomProperties (gmx::ArrayRef< const int > atomTypes, gmx::ArrayRef< const real > atomCharges, gmx::ArrayRef< const int32_t > atomInfo) const |
Updates all the atom properties in Nbnxm. | |
void | convertCoordinates (gmx::AtomLocality locality, gmx::ArrayRef< const gmx::RVec > coordinates) |
Convert the coordinates to NBNXM format for the given locality. More... | |
void | convertCoordinatesGpu (gmx::AtomLocality locality, DeviceBuffer< gmx::RVec > d_x, GpuEventSynchronizer *xReadyOnDevice) |
Convert the coordinates to NBNXM format on the GPU for the given locality. More... | |
void | atomdata_init_copy_x_to_nbat_x_gpu () const |
Init for GPU version of setup coordinates in Nbnxm. | |
const PairlistSets & | pairlistSets () const |
Returns a reference to the pairlist sets. | |
bool | isDynamicPruningStepCpu (int64_t step) const |
Returns whether step is a dynamic list pruning step, for CPU lists. | |
bool | isDynamicPruningStepGpu (int64_t step) const |
Returns whether step is a dynamic list pruning step, for GPU lists. | |
void | dispatchPruneKernelCpu (gmx::InteractionLocality iLocality, gmx::ArrayRef< const gmx::RVec > shift_vec) const |
Dispatches the dynamic pruning kernel for the given locality, for CPU lists. | |
void | dispatchPruneKernelGpu (int64_t step) |
Dispatches the dynamic pruning kernel for GPU lists. | |
void | dispatchNonbondedKernel (gmx::InteractionLocality iLocality, const interaction_const_t &ic, const gmx::StepWorkload &stepWork, int clearF, gmx::ArrayRef< const gmx::RVec > shiftvec, gmx::ArrayRef< real > repulsionDispersionSR, gmx::ArrayRef< real > CoulombSR, t_nrnb *nrnb) const |
Executes the non-bonded kernel of the GPU or launches it on the GPU. | |
void | dispatchFreeEnergyKernels (const gmx::ArrayRefWithPadding< const gmx::RVec > &coords, gmx::ForceWithShiftForces *forceWithShiftForces, bool useSimd, int ntype, const interaction_const_t &ic, gmx::ArrayRef< const gmx::RVec > shiftvec, gmx::ArrayRef< const real > nbfp, gmx::ArrayRef< const real > nbfp_grid, gmx::ArrayRef< const real > chargeA, gmx::ArrayRef< const real > chargeB, gmx::ArrayRef< const int > typeA, gmx::ArrayRef< const int > typeB, gmx::ArrayRef< const real > lambda, gmx_enerdata_t *enerd, const gmx::StepWorkload &stepWork, t_nrnb *nrnb) |
Executes the non-bonded free-energy kernels, local + non-local, always runs on the CPU. | |
void | atomdata_add_nbat_f_to_f (gmx::AtomLocality locality, gmx::ArrayRef< gmx::RVec > force) |
Add the forces stored in nbat to f, zeros the forces in nbat. More... | |
int | getNumAtoms (gmx::AtomLocality locality) const |
Get the number of atoms for a given locality. More... | |
const Nbnxm::KernelSetup & | kernelSetup () const |
Return the kernel setup. | |
real | pairlistInnerRadius () const |
Returns the outer radius for the pair list. | |
real | pairlistOuterRadius () const |
Returns the outer radius for the pair list. | |
void | changePairlistRadii (real rlistOuter, real rlistInner) const |
Changes the pair-list outer and inner radius. | |
void | setupGpuShortRangeWork (const gmx::ListedForcesGpu *listedForcesGpu, gmx::InteractionLocality iLocality) const |
Set up internal flags that indicate what type of short-range work there is. | |
void | setupFepThreadedForceBuffer (int numAtomsForce) |
nbnxn_atomdata_t & | nbat () |
Returns a reference to the nbnxn_atomdata_t object. | |
const NbnxmGpu * | gpuNbv () const |
Returns a pointer to the NbnxmGpu object, can return nullptr. | |
NbnxmGpu * | gpuNbv () |
Returns a pointer to the NbnxmGpu object, can return nullptr. | |
nonbonded_verlet_t::nonbonded_verlet_t | ( | std::unique_ptr< PairlistSets > | pairlistSets, |
std::unique_ptr< PairSearch > | pairSearch, | ||
std::unique_ptr< nbnxn_atomdata_t > | nbat, | ||
const Nbnxm::KernelSetup & | kernelSetup, | ||
std::unique_ptr< ExclusionChecker > | exclusionChecker, | ||
NbnxmGpu * | gpu_nbv, | ||
gmx_wallcycle * | wcycle | ||
) |
Constructs an object from its components.
[in] | pairlistSets | The pairlist sets, is consumed |
[in] | pairSearch | The pairsearch setup, is consumed |
[in] | nbat | The atom data, is consumed |
[in] | kernelSetup | The non-bonded kernel setup |
[in] | exclusionChecker | The FEP exclusion checker, is consumed, can be nullptr |
[in] | gpu_nbv | The GPU non-bonded setup, ownership is transferred, can be nullptr |
[in] | wcycle | Pointer to wallcycle counters, can be nullptr |
nonbonded_verlet_t::nonbonded_verlet_t | ( | std::unique_ptr< PairlistSets > | pairlistSets, |
std::unique_ptr< PairSearch > | pairSearch, | ||
std::unique_ptr< nbnxn_atomdata_t > | nbat, | ||
const Nbnxm::KernelSetup & | kernelSetup, | ||
NbnxmGpu * | gpu_nbv | ||
) |
Constructs an object from its, minimal, components.
[in] | pairlistSets | The pairlist sets, is consumed |
[in] | pairSearch | The pairsearch setup, is consumed |
[in] | nbat | The atom data, is consumed |
[in] | kernelSetup | The non-bonded kernel setup |
[in] | gpu_nbv | The GPU non-bonded setup, ownership is transferred, can be nullptr |
void nonbonded_verlet_t::atomdata_add_nbat_f_to_f | ( | gmx::AtomLocality | locality, |
gmx::ArrayRef< gmx::RVec > | force | ||
) |
Add the forces stored in nbat to f, zeros the forces in nbat.
[in] | locality | Local or non-local |
[in,out] | force | Force to be added to |
void nonbonded_verlet_t::constructPairlist | ( | gmx::InteractionLocality | iLocality, |
const gmx::ListOfLists< int > & | exclusions, | ||
int64_t | step, | ||
t_nrnb * | nrnb | ||
) | const |
Constructs the pairlist for the given locality.
When there are no non-self exclusions, exclusions
can be empty. Otherwise the number of lists in exclusions
should match the number of atoms when not using DD, or the total number of atoms in the i-zones when using DD.
[in] | iLocality | The interaction locality: local or non-local |
[in] | exclusions | Lists of exclusions for every atom. |
[in] | step | Used to set the list creation step |
[in,out] | nrnb | Flop accounting struct, can be nullptr |
void nonbonded_verlet_t::convertCoordinates | ( | gmx::AtomLocality | locality, |
gmx::ArrayRef< const gmx::RVec > | coordinates | ||
) |
Convert the coordinates to NBNXM format for the given locality.
The API function for the transformation of the coordinates from one layout to another.
[in] | locality | Whether coordinates for local or non-local atoms should be transformed. |
[in] | coordinates | Coordinates in plain rvec format to be transformed. |
void nonbonded_verlet_t::convertCoordinatesGpu | ( | gmx::AtomLocality | locality, |
DeviceBuffer< gmx::RVec > | d_x, | ||
GpuEventSynchronizer * | xReadyOnDevice | ||
) |
Convert the coordinates to NBNXM format on the GPU for the given locality.
The API function for the transformation of the coordinates from one layout to another in the GPU memory.
[in] | locality | Whether coordinates for local or non-local atoms should be transformed. |
[in] | d_x | GPU coordinates buffer in plain rvec format to be transformed. |
[in] | xReadyOnDevice | Event synchronizer indicating that the coordinates are ready in the device memory. |
int nonbonded_verlet_t::getNumAtoms | ( | gmx::AtomLocality | locality | ) | const |
Get the number of atoms for a given locality.
[in] | locality | Local or non-local |
void nonbonded_verlet_t::putAtomsOnGrid | ( | const matrix | box, |
int | gridIndex, | ||
const gmx::RVec & | lowerCorner, | ||
const gmx::RVec & | upperCorner, | ||
const gmx::UpdateGroupsCog * | updateGroupsCog, | ||
gmx::Range< int > | atomRange, | ||
real | atomDensity, | ||
gmx::ArrayRef< const int32_t > | atomInfo, | ||
gmx::ArrayRef< const gmx::RVec > | x, | ||
int | numAtomsMoved, | ||
const int * | move | ||
) |
Put the atoms on the pair search grid.
Only atoms with indices wihtin atomRange
in x are put on the grid. When updateGroupsCog
!= nullptr, atoms are put on the grid based on the center of geometry of the group they belong to. Atoms or COGs of groups should be within the bounding box provided, this is checked in debug builds when not using update groups. The atom density is used to determine the grid size when gridIndex
= 0. When atomDensity
<= 0, the density is determined from atomEnd-atomStart and the bounding box corners. With domain decomposition, part of the atoms might have migrated, but have not been removed yet. This count is given by numAtomsMoved
. When move
[i] < 0 particle i has migrated and will not be put on the grid.
[in] | box | Box used for periodic distance calculations |
[in] | gridIndex | The index of the grid to spread to, always 0 except with test particle insertion |
[in] | lowerCorner | Atom groups to be gridded should have coordinates >= this corner |
[in] | upperCorner | Atom groups to be gridded should have coordinates <= this corner |
[in] | updateGroupsCog | Centers of geometry for update groups, pass nullptr when not using update groups |
[in] | atomRange | Range of atoms to grid |
[in] | atomDensity | An estimate of the atom density, used for peformance optimization and only with gridIndex = 0 |
[in] | atomInfo | Atom information flags |
[in] | x | Coordinates for atoms to grid |
[in] | numAtomsMoved | The number of atoms that will move to another domain, pass 0 without DD |
[in] | move | Move flags for atoms, pass nullptr without DD |