Gromacs  2023.2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Classes | Enumerations | Functions | Variables
#include <memory>
#include "gromacs/gpu_utils/devicebuffer_datatype.h"
#include "gromacs/math/vectypes.h"
#include "gromacs/mdtypes/locality.h"
#include "gromacs/utility/arrayref.h"
#include "gromacs/utility/enumerationhelpers.h"
#include "gromacs/utility/real.h"
+ Include dependency graph for nbnxm.h:
+ This graph shows which files directly or indirectly include this file:

Description

This file contains the public interface of the nbnxm module that implements the NxM atom cluster non-bonded algorithm to efficiently compute pair forces.

Author
Berk Hess hess@.nosp@m.kth..nosp@m.se
Szilárd Páll pall..nosp@m.szil.nosp@m.ard@g.nosp@m.mail.nosp@m..com

Classes

class  gmx::ArrayRefWithPadding< typename >
 Interface to a C array of T (or part of a std container of T), that includes padding that is suitable for the kinds of SIMD operations GROMACS uses. More...
 
class  gmx::ListOfLists< typename >
 A list of lists, optimized for performance. More...
 
class  gmx::Range< typename >
 Defines a range of integer numbers and accompanying operations. More...
 
struct  nonbonded_verlet_t
 Top-level non-bonded data structure for the Verlet-type cut-off scheme. More...
 

Enumerations

enum  Nbnxm::ElecType : int {
  Nbnxm::ElecType::Cut, Nbnxm::ElecType::RF, Nbnxm::ElecType::EwaldTab, Nbnxm::ElecType::EwaldTabTwin,
  Nbnxm::ElecType::EwaldAna, Nbnxm::ElecType::EwaldAnaTwin, Nbnxm::ElecType::Count
}
 Nbnxm electrostatic GPU kernel flavors. More...
 
enum  Nbnxm::VdwType : int {
  Nbnxm::VdwType::Cut, Nbnxm::VdwType::CutCombGeom, Nbnxm::VdwType::CutCombLB, Nbnxm::VdwType::FSwitch,
  Nbnxm::VdwType::PSwitch, Nbnxm::VdwType::EwaldGeom, Nbnxm::VdwType::EwaldLB, Nbnxm::VdwType::Count
}
 Nbnxm VdW GPU kernel flavors. More...
 
enum  Nbnxm::KernelType : int {
  NotSet = 0, Cpu4x4_PlainC, Cpu4xN_Simd_4xN, Cpu4xN_Simd_2xNN,
  Gpu8x8x8, Cpu8x8x8_PlainC, Count
}
 Nonbonded NxN kernel types: plain C, CPU SIMD, GPU, GPU emulation.
 
enum  Nbnxm::EwaldExclusionType : int { NotSet = 0, Table, Analytical, DecidedByGpuModule }
 Ewald exclusion types.
 
enum  { enbvClearFNo, enbvClearFYes }
 Flag to tell the nonbonded kernels whether to clear the force output buffers.
 

Functions

const char * Nbnxm::lookup_kernel_name (Nbnxm::KernelType kernelType)
 Return a string identifying the kernel type. More...
 
std::unique_ptr
< nonbonded_verlet_t
Nbnxm::init_nb_verlet (const gmx::MDLogger &mdlog, const t_inputrec &inputrec, const t_forcerec &forcerec, const t_commrec *commrec, const gmx_hw_info_t &hardwareInfo, bool useGpuForNonbonded, const gmx::DeviceStreamManager *deviceStreamManager, const gmx_mtop_t &mtop, gmx::ObservablesReducerBuilder *observablesReducerBuilder, gmx::ArrayRef< const gmx::RVec > coordinates, matrix box, gmx_wallcycle *wcycle)
 Creates an Nbnxm object.
 
void nbnxn_put_on_grid (nonbonded_verlet_t *nb_verlet, const matrix box, int gridIndex, const rvec lowerCorner, const rvec upperCorner, const gmx::UpdateGroupsCog *updateGroupsCog, gmx::Range< int > atomRange, real atomDensity, gmx::ArrayRef< const int64_t > atomInfo, gmx::ArrayRef< const gmx::RVec > x, int numAtomsMoved, const int *move)
 Put the atoms on the pair search grid. More...
 
void nbnxn_put_on_grid_nonlocal (nonbonded_verlet_t *nb_verlet, const struct gmx_domdec_zones_t *zones, gmx::ArrayRef< const int64_t > atomInfo, gmx::ArrayRef< const gmx::RVec > x)
 As nbnxn_put_on_grid, but for the non-local atoms. More...
 
bool buildSupportsNonbondedOnGpu (std::string *error)
 Check if GROMACS has been built with GPU support. More...
 

Variables

constexpr int Nbnxm::c_numElecTypes = static_cast<int>(ElecType::Count)
 Number of possible ElecType values.
 
constexpr int Nbnxm::c_numVdwTypes = static_cast<int>(VdwType::Count)
 Number of possible VdwType values.
 

Function Documentation

bool buildSupportsNonbondedOnGpu ( std::string *  error)

Check if GROMACS has been built with GPU support.

Parameters
[in]errorPointer to error string or nullptr.
Todo:
Move this to NB module once it exists.
void nbnxn_put_on_grid ( nonbonded_verlet_t nb_verlet,
const matrix  box,
int  gridIndex,
const rvec  lowerCorner,
const rvec  upperCorner,
const gmx::UpdateGroupsCog updateGroupsCog,
gmx::Range< int >  atomRange,
real  atomDensity,
gmx::ArrayRef< const int64_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.

Parameters
[in,out]nb_verletThe non-bonded object
[in]boxBox used for periodic distance calculations
[in]gridIndexThe index of the grid to spread to, always 0 except with test particle insertion
[in]lowerCornerAtom groups to be gridded should have coordinates >= this corner
[in]upperCornerAtom groups to be gridded should have coordinates <= this corner
[in]updateGroupsCogCenters of geometry for update groups, pass nullptr when not using update groups
[in]atomRangeRange of atoms to grid
[in]atomDensityAn estimate of the atom density, used for peformance optimization and only with gridIndex = 0
[in]atomInfoAtom information flags
[in]xCoordinates for atoms to grid
[in]numAtomsMovedThe number of atoms that will move to another domain, pass 0 without DD
[in]moveMove flags for atoms, pass nullptr without DD
void nbnxn_put_on_grid_nonlocal ( nonbonded_verlet_t nb_verlet,
const struct gmx_domdec_zones_t *  zones,
gmx::ArrayRef< const int64_t >  atomInfo,
gmx::ArrayRef< const gmx::RVec x 
)

As nbnxn_put_on_grid, but for the non-local atoms.

with domain decomposition. Should be called after calling nbnxn_search_put_on_grid for the local atoms / home zone.