Gromacs  2025-dev-20240913-b871546
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
List of all members | Public Member Functions
gmx::WholeMoleculeTransform Class Reference

#include <gromacs/mdlib/wholemoleculetransform.h>

Description

This class manages a coordinate buffer with molecules not split over periodic boundary conditions for use in force calculations which require whole molecules.

Note: This class should not be used for computation of forces which have virial contributions through shift forces.

Public Member Functions

 WholeMoleculeTransform (const gmx_mtop_t &mtop, PbcType pbcType, bool useAtomReordering)
 Constructor. More...
 
void updateAtomOrder (ArrayRef< const int > globalAtomIndices, const gmx_ga2la_t &ga2la)
 Changes the atom order to the one provided. More...
 
void updateForAtomPbcJumps (ArrayRef< const RVec > x, const matrix box)
 Updates the graph when atoms have been shifted by periodic vectors.
 
ArrayRef< const RVecwholeMoleculeCoordinates (ArrayRef< const RVec > x, const matrix box)
 Create and return coordinates with whole molecules for input coordinates x. More...
 

Constructor & Destructor Documentation

gmx::WholeMoleculeTransform::WholeMoleculeTransform ( const gmx_mtop_t &  mtop,
PbcType  pbcType,
bool  useAtomReordering 
)

Constructor.

Parameters
[in]mtopThe global topology use for getting the connections between atoms
[in]pbcTypeThe type of PBC
[in]useAtomReorderingWhether we will use atom reordering

Member Function Documentation

void gmx::WholeMoleculeTransform::updateAtomOrder ( ArrayRef< const int >  globalAtomIndices,
const gmx_ga2la_t ga2la 
)

Changes the atom order to the one provided.

This method is called after domain repartitioning. The object should have been constructed with useAtomReordering set to true.

Parameters
[in]globalAtomIndicesThe global atom indices for the local atoms, size should be the system size
[in]ga2laGlobal to local atom index lookup (the inverse of globalAtomIndices)
ArrayRef< const RVec > gmx::WholeMoleculeTransform::wholeMoleculeCoordinates ( ArrayRef< const RVec x,
const matrix  box 
)

Create and return coordinates with whole molecules for input coordinates x.

Parameters
[in]xInput coordinates, should not have periodic displacement compared with the coordinates passed in the last call to updateForAtomPbcJumps().
[in]boxThe current periodic image vectors

Note: this operation is not free. If you need whole molecules coordinates more than once during the force calculation, store the result and reuse it.


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