Gromacs  2022.2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Macros | Typedefs | Enumerations | Functions | Variables
atomdata.h File Reference
#include <cstdio>
#include "gromacs/gpu_utils/devicebuffer_datatype.h"
#include "gromacs/gpu_utils/hostallocator.h"
#include "gromacs/math/vectypes.h"
#include "gromacs/mdtypes/locality.h"
#include "gromacs/utility/bitmask.h"
#include "gromacs/utility/real.h"
+ Include dependency graph for atomdata.h:
+ This graph shows which files directly or indirectly include this file:

Description

Functionality for per-atom data in the nbnxm module.

Author
Berk Hess hess@.nosp@m.kth..nosp@m.se

Macros

#define NBNXN_BUFFERFLAG_SIZE   16
 Block size in atoms for the non-bonded thread force-buffer reduction. More...
 
#define NBNXN_BUFFERFLAG_MAX_THREADS   (BITMASK_SIZE)
 We store the reduction flags as gmx_bitmask_t. This limits the number of flags to BITMASK_SIZE.
 

Typedefs

template<typename T >
using AlignedVector = std::vector< T, gmx::AlignedAllocator< T >>
 Convenience type for vector with aligned memory.
 

Enumerations

enum  { nbatXYZ, nbatXYZQ, nbatX4, nbatX8 }
 
enum  LJCombinationRule : int { LJCombinationRule::Geometric, LJCombinationRule::LorentzBerthelot, LJCombinationRule::None, LJCombinationRule::Count }
 LJ combination rules. More...
 
enum  { enbnxninitcombruleDETECT, enbnxninitcombruleGEOM, enbnxninitcombruleLB, enbnxninitcombruleNONE }
 Describes the combination rule in use by this force field.
 

Functions

template<int packSize>
static int atom_to_x_index (int a)
 Returns the index in a coordinate array corresponding to atom a.
 
const char * enumValueToString (LJCombinationRule enumValue)
 String corresponding to LJ combination rule.
 
void copy_rvec_to_nbat_real (const int *a, int na, int na_round, const rvec *x, int nbatFormat, real *xnb, int a0)
 Copy na rvec elements from x to xnb using nbatFormat, start dest a0, and fills up to na_round with coordinates that are far away.
 
void nbnxn_atomdata_set (nbnxn_atomdata_t *nbat, const Nbnxm::GridSet &gridSet, gmx::ArrayRef< const int > atomTypes, gmx::ArrayRef< const real > atomCharges, gmx::ArrayRef< const int64_t > atomInfo)
 Sets the atomdata after pair search.
 
void nbnxn_atomdata_copy_shiftvec (bool dynamic_box, gmx::ArrayRef< gmx::RVec > shift_vec, nbnxn_atomdata_t *nbat)
 Copy the shift vectors to nbat.
 
void nbnxn_atomdata_copy_x_to_nbat_x (const Nbnxm::GridSet &gridSet, gmx::AtomLocality locality, const rvec *coordinates, nbnxn_atomdata_t *nbat)
 Transform coordinates to xbat layout. More...
 
void nbnxn_atomdata_x_to_nbat_x_gpu (const Nbnxm::GridSet &gridSet, gmx::AtomLocality locality, NbnxmGpu *gpu_nbv, DeviceBuffer< gmx::RVec > d_x, GpuEventSynchronizer *xReadyOnDevice)
 Transform coordinates to xbat layout on GPU. More...
 
void reduceForces (nbnxn_atomdata_t *nbat, gmx::AtomLocality locality, const Nbnxm::GridSet &gridSet, rvec *totalForce)
 Add the computed forces to f, an internal reduction might be performed as well. More...
 
void nbnxn_atomdata_add_nbat_fshift_to_fshift (const nbnxn_atomdata_t &nbat, gmx::ArrayRef< gmx::RVec > fshift)
 Add the fshift force stored in nbat to fshift.
 
void nbnxn_get_atom_range (gmx::AtomLocality atomLocality, const Nbnxm::GridSet &gridSet, int *atomStart, int *nAtoms)
 Get the atom start index and number of atoms for a given locality.
 

Variables

static constexpr int STRIDE_XYZ = 3
 Stride for coordinate/force arrays with xyz coordinate storage.
 
static constexpr int STRIDE_XYZQ = 4
 Stride for coordinate/force arrays with xyzq coordinate storage.
 
static constexpr int c_packX4 = 4
 Size of packs of x, y or z with SIMD 4-grouped packed coordinates/forces.
 
static constexpr int c_packX8 = 8
 Size of packs of x, y or z with SIMD 8-grouped packed coordinates/forces.
 
static constexpr int STRIDE_P4 = DIM * c_packX4
 Stridefor a pack of 4 coordinates/forces.
 
static constexpr int STRIDE_P8 = DIM * c_packX8
 Stridefor a pack of 8 coordinates/forces.
 

Macro Definition Documentation

#define NBNXN_BUFFERFLAG_SIZE   16

Block size in atoms for the non-bonded thread force-buffer reduction.

Should be a multiple of all cell and x86 SIMD sizes (i.e. 2, 4 and 8). Should be small to reduce the reduction and zeroing cost, but too small will result in overhead. Currently the block size is NBNXN_BUFFERFLAG_SIZE*3*sizeof(real)=192 bytes.

Enumeration Type Documentation

enum LJCombinationRule : int
strong

LJ combination rules.

Enumerator
Geometric 

Geometric.

LorentzBerthelot 

Lorentz-Berthelot.

None 

No rule.

Count 

Size of the enum.

Function Documentation

void nbnxn_atomdata_copy_x_to_nbat_x ( const Nbnxm::GridSet &  gridSet,
gmx::AtomLocality  locality,
const rvec *  coordinates,
nbnxn_atomdata_t *  nbat 
)

Transform coordinates to xbat layout.

Creates a copy of the coordinates buffer using short-range ordering.

Parameters
[in]gridSetThe grids data.
[in]localityIf the transformation should be applied to local or non local coordinates.
[in]coordinatesCoordinates in plain rvec format.
[in,out]nbatData in NBNXM format, used for mapping formats and to locate the output buffer.
void nbnxn_atomdata_x_to_nbat_x_gpu ( const Nbnxm::GridSet &  gridSet,
gmx::AtomLocality  locality,
NbnxmGpu *  gpu_nbv,
DeviceBuffer< gmx::RVec d_x,
GpuEventSynchronizer xReadyOnDevice 
)

Transform coordinates to xbat layout on GPU.

Creates a GPU copy of the coordinates buffer using short-range ordering. As input, uses coordinates in plain rvec format in GPU memory.

Parameters
[in]gridSetThe grids data.
[in]localityIf the transformation should be applied to local or non local coordinates.
[in,out]gpu_nbvThe NBNXM GPU data structure.
[in]d_xCoordinates to be copied (in plain rvec format).
[in]xReadyOnDeviceEvent synchronizer indicating that the coordinates are ready in the device memory.
void reduceForces ( nbnxn_atomdata_t *  nbat,
gmx::AtomLocality  locality,
const Nbnxm::GridSet &  gridSet,
rvec *  totalForce 
)

Add the computed forces to f, an internal reduction might be performed as well.

Parameters
[in]nbatAtom data in NBNXM format.
[in]localityIf the reduction should be performed on local or non-local atoms.
[in]gridSetThe grids data.
[out]totalForceBuffer to accumulate resulting force