Gromacs
2024.4
|
#include <gromacs/mdlib/update_constrain_gpu_impl.h>
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... | |
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.
[in] | ir | Input record data: LINCS takes number of iterations and order of projection from it. |
[in] | mtop | Topology of the system: SETTLE gets the masses for O and H atoms and target O-H and H-H distances from this object. |
[in] | numTempScaleValues | Number of temperature scaling groups. Set zero for no temperature coupling. |
[in] | deviceContext | GPU device context. |
[in] | deviceStream | GPU stream to use. |
[in] | wcycle | The wallclock counter |
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:
[in] | fReadyOnDevice | Event synchronizer indicating that the forces are ready in the device memory. |
[in] | dt | Timestep. |
[in] | updateVelocities | If the velocities should be constrained. |
[in] | computeVirial | If virial should be updated. |
[out] | virial | Place to save virial tensor. |
[in] | doTemperatureScaling | If velocities should be scaled for temperature coupling. |
[in] | tcstat | Temperature coupling data. |
[in] | doParrinelloRahman | If current step is a Parrinello-Rahman pressure coupling step. |
[in] | dtPressureCouple | Period between pressure coupling steps. |
[in] | prVelocityScalingMatrix | Parrinello-Rahman velocity scaling matrix. |
|
static |
Returns whether the maximum number of coupled constraints is supported by the GPU LINCS code.
[in] | mtop | The 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.
[in] | scalingMatrix | Coordinates 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.
[in] | scalingMatrix | Velocities 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).
[in,out] | d_x | Device buffer with coordinates. |
[in,out] | d_v | Device buffer with velocities. |
[in] | d_f | Device buffer with forces. |
[in] | idef | System topology |
[in] | md | Atoms data. |
void gmx::UpdateConstrainGpu::Impl::setPbc | ( | PbcType | pbcType, |
const matrix | box | ||
) |