Gromacs  2025-dev-20241003-bd59e46
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
List of all members | Public Types | Public Member Functions | Static Public Member Functions
ListedForces Class Reference

#include <gromacs/listed_forces/listed_forces.h>

Description

Class for calculating listed interactions, uses OpenMP parallelization.

Listed interactions can be divided over multiple instances of ListedForces using the selection flags passed to the constructor.

Public Types

enum  InteractionGroup : int {
  InteractionGroup::Pairs, InteractionGroup::Dihedrals, InteractionGroup::Angles, InteractionGroup::Rest,
  InteractionGroup::Count
}
 Enum for selecting groups of listed interaction types. More...
 
using InteractionSelection = std::bitset< static_cast< int >(InteractionGroup::Count)>
 Type for specifying selections of groups of interaction types.
 

Public Member Functions

 ListedForces (const gmx_ffparams_t &ffparams, int numEnergyGroups, int numThreads, InteractionSelection interactionSelection, FILE *fplog)
 Constructor. More...
 
 ListedForces (ListedForces &&o) noexcept
 Move constructor, default, but in the source file to hide implementation classes.
 
 ~ListedForces ()
 Destructor which is actually default but in the source file to hide implementation classes.
 
void setup (const InteractionDefinitions &domainIdef, int numAtomsForce, bool useGpu)
 Copy the listed interactions from idef and set up the thread parallelization. More...
 
void calculate (struct gmx_wallcycle *wcycle, const matrix box, const t_commrec *cr, const gmx_multisim_t *ms, gmx::ArrayRefWithPadding< const gmx::RVec > coordinates, gmx::ArrayRef< const gmx::RVec > xWholeMolecules, t_fcdata *fcdata, const history_t *hist, gmx::ForceOutputs *forceOutputs, const t_forcerec *fr, const struct t_pbc *pbc, gmx_enerdata_t *enerd, t_nrnb *nrnb, gmx::ArrayRef< const real > lambda, gmx::ArrayRef< const real > chargeA, gmx::ArrayRef< const real > chargeB, gmx::ArrayRef< const bool > atomIsPerturbed, gmx::ArrayRef< const unsigned short > cENER, int nPerturbed, int *global_atom_index, const gmx::StepWorkload &stepWork)
 Do all aspects of energy and force calculations for mdrun on the set of listed interactions. More...
 
bool haveCpuBondeds () const
 Returns whether bonded interactions are assigned to the CPU.
 
bool haveCpuListedForces (const t_fcdata &fcdata) const
 Returns whether listed forces are computed on the CPU. More...
 
bool haveRestraints (const t_fcdata &fcdata) const
 Returns true if there are position, distance or orientation restraints.
 

Static Public Member Functions

static InteractionSelection interactionSelectionAll ()
 Returns a selection with all listed interaction types selected.
 

Member Enumeration Documentation

enum ListedForces::InteractionGroup : int
strong

Enum for selecting groups of listed interaction types.

Enumerator
Pairs 

Pair interactions.

Dihedrals 

Dihedrals, including cmap.

Angles 

Angles.

Rest 

All listed interactions that are not any of the above.

Count 

The number of items above.

Constructor & Destructor Documentation

ListedForces::ListedForces ( const gmx_ffparams_t ffparams,
int  numEnergyGroups,
int  numThreads,
InteractionSelection  interactionSelection,
FILE *  fplog 
)

Constructor.

Parameters
[in]ffparamsThe force field parameters
[in]numEnergyGroupsThe number of energy groups, used for storage of pair energies
[in]numThreadsThe number of threads used for computed listed interactions
[in]interactionSelectionSelect of interaction groups through bits set
[in]fplogLog file for printing env.var. override, can be nullptr

Member Function Documentation

void ListedForces::calculate ( struct gmx_wallcycle *  wcycle,
const matrix  box,
const t_commrec *  cr,
const gmx_multisim_t ms,
gmx::ArrayRefWithPadding< const gmx::RVec coordinates,
gmx::ArrayRef< const gmx::RVec xWholeMolecules,
t_fcdata *  fcdata,
const history_t hist,
gmx::ForceOutputs forceOutputs,
const t_forcerec *  fr,
const struct t_pbc pbc,
gmx_enerdata_t enerd,
t_nrnb *  nrnb,
gmx::ArrayRef< const real lambda,
gmx::ArrayRef< const real chargeA,
gmx::ArrayRef< const real chargeB,
gmx::ArrayRef< const bool >  atomIsPerturbed,
gmx::ArrayRef< const unsigned short >  cENER,
int  nPerturbed,
int *  global_atom_index,
const gmx::StepWorkload stepWork 
)

Do all aspects of energy and force calculations for mdrun on the set of listed interactions.

xWholeMolecules only needs to contain whole molecules when orientation restraints need to be computed and can be empty otherwise.

bool ListedForces::haveCpuListedForces ( const t_fcdata &  fcdata) const

Returns whether listed forces are computed on the CPU.

NOTE: the current implementation returns true if there are position restraints or any bonded interactions computed on the CPU.

void ListedForces::setup ( const InteractionDefinitions domainIdef,
int  numAtomsForce,
bool  useGpu 
)

Copy the listed interactions from idef and set up the thread parallelization.

Parameters
[in]domainIdefInteraction definitions for all listed interactions to be computed on this domain/rank
[in]numAtomsForceForce are, potentially, computed for atoms 0 to numAtomsForce
[in]useGpuWhether a GPU is used to compute (part of) the listed interactions

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