Gromacs  2024.4
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
List of all members | Public Member Functions | Static Public Member Functions
gmx::UpdateConstrainGpu::Impl Class Reference

#include <gromacs/mdlib/update_constrain_gpu_impl.h>

Description

Class with interfaces and data for GPU version of Update-Constraint.

Public Member Functions

 Impl (const t_inputrec &ir, const gmx_mtop_t &mtop, int numTempScaleValues, const DeviceContext &deviceContext, const DeviceStream &deviceStream, gmx_wallcycle *wcycle)
 Create Update-Constrain object. More...
 
void integrate (GpuEventSynchronizer *fReadyOnDevice, real dt, bool updateVelocities, bool computeVirial, tensor virial, bool doTemperatureScaling, gmx::ArrayRef< const t_grp_tcstat > tcstat, bool doParrinelloRahman, float dtPressureCouple, const gmx::Matrix3x3 &prVelocityScalingMatrix)
 Integrate. More...
 
void scaleCoordinates (const Matrix3x3 &scalingMatrix)
 Scale coordinates on the GPU for the pressure coupling. More...
 
void scaleVelocities (const Matrix3x3 &scalingMatrix)
 Scale velocities on the GPU for the pressure coupling. More...
 
void set (DeviceBuffer< Float3 > d_x, DeviceBuffer< Float3 > d_v, DeviceBuffer< Float3 > d_f, const InteractionDefinitions &idef, const t_mdatoms &md)
 Set the pointers and update data-structures (e.g. after NB search step). More...
 
void setPbc (PbcType pbcType, const matrix box)
 Update PBC data. More...
 
GpuEventSynchronizer * xUpdatedOnDeviceEvent ()
 Return the synchronizer associated with the event indicated that the coordinates are ready on the device.
 

Static Public Member Functions

static bool isNumCoupledConstraintsSupported (const gmx_mtop_t &mtop)
 Returns whether the maximum number of coupled constraints is supported by the GPU LINCS code. More...
 

Constructor & Destructor Documentation

gmx::UpdateConstrainGpu::Impl::Impl ( const t_inputrec &  ir,
const gmx_mtop_t &  mtop,
int  numTempScaleValues,
const DeviceContext &  deviceContext,
const DeviceStream deviceStream,
gmx_wallcycle *  wcycle 
)

Create Update-Constrain object.

The constructor is given a non-nullptr deviceStream, in which all the update and constrain routines are executed.

Parameters
[in]irInput record data: LINCS takes number of iterations and order of projection from it.
[in]mtopTopology of the system: SETTLE gets the masses for O and H atoms and target O-H and H-H distances from this object.
[in]numTempScaleValuesNumber of temperature scaling groups. Set zero for no temperature coupling.
[in]deviceContextGPU device context.
[in]deviceStreamGPU stream to use.
[in]wcycleThe wallclock counter

Member Function Documentation

void gmx::UpdateConstrainGpu::Impl::integrate ( GpuEventSynchronizer *  fReadyOnDevice,
real  dt,
bool  updateVelocities,
bool  computeVirial,
tensor  virial,
bool  doTemperatureScaling,
gmx::ArrayRef< const t_grp_tcstat >  tcstat,
bool  doParrinelloRahman,
float  dtPressureCouple,
const gmx::Matrix3x3 prVelocityScalingMatrix 
)

Integrate.

Integrates the equation of motion using Leap-Frog algorithm and applies LINCS and SETTLE constraints. If computeVirial is true, constraints virial is written at the provided pointer. doTempCouple should be true if:

  1. The temperature coupling is enabled.
  2. This is the temperature coupling step. Parameters virial/lambdas can be nullptr if computeVirial/doTempCouple are false.
Parameters
[in]fReadyOnDeviceEvent synchronizer indicating that the forces are ready in the device memory.
[in]dtTimestep.
[in]updateVelocitiesIf the velocities should be constrained.
[in]computeVirialIf virial should be updated.
[out]virialPlace to save virial tensor.
[in]doTemperatureScalingIf velocities should be scaled for temperature coupling.
[in]tcstatTemperature coupling data.
[in]doParrinelloRahmanIf current step is a Parrinello-Rahman pressure coupling step.
[in]dtPressureCouplePeriod between pressure coupling steps.
[in]prVelocityScalingMatrixParrinello-Rahman velocity scaling matrix.
static bool gmx::UpdateConstrainGpu::Impl::isNumCoupledConstraintsSupported ( const gmx_mtop_t &  mtop)
static

Returns whether the maximum number of coupled constraints is supported by the GPU LINCS code.

Parameters
[in]mtopThe molecular topology
void gmx::UpdateConstrainGpu::Impl::scaleCoordinates ( const Matrix3x3 scalingMatrix)

Scale coordinates on the GPU for the pressure coupling.

After pressure coupling step, the box size may change. Hence, the coordinates should be scaled so that all the particles fit in the new box.

Parameters
[in]scalingMatrixCoordinates scaling matrix.
void gmx::UpdateConstrainGpu::Impl::scaleVelocities ( const Matrix3x3 scalingMatrix)

Scale velocities on the GPU for the pressure coupling.

After pressure coupling step, the box size may change. In the C-Rescale algorithm, velocities should be scaled.

Parameters
[in]scalingMatrixVelocities scaling matrix.
void gmx::UpdateConstrainGpu::Impl::set ( DeviceBuffer< Float3 d_x,
DeviceBuffer< Float3 d_v,
DeviceBuffer< Float3 d_f,
const InteractionDefinitions &  idef,
const t_mdatoms md 
)

Set the pointers and update data-structures (e.g. after NB search step).

Parameters
[in,out]d_xDevice buffer with coordinates.
[in,out]d_vDevice buffer with velocities.
[in]d_fDevice buffer with forces.
[in]idefSystem topology
[in]mdAtoms data.
void gmx::UpdateConstrainGpu::Impl::setPbc ( PbcType  pbcType,
const matrix  box 
)

Update PBC data.

Converts PBC data from t_pbc into the PbcAiuc format and stores the latter.

Parameters
[in]pbcTypeThe type of the periodic boundary.
[in]boxThe periodic boundary box matrix.

The documentation for this class was generated from the following files: