Gromacs
2025-dev-20241003-bd59e46
|
#include <cstdio>
#include <optional>
#include "gromacs/math/vectypes.h"
#include "gromacs/mdtypes/pull_params.h"
#include "gromacs/utility/arrayref.h"
#include "gromacs/utility/basedefinitions.h"
#include "gromacs/utility/real.h"
This file contains datatypes and function declarations necessary for mdrun to interface with the pull code.
Classes | |
class | gmx::ArrayRef< typename > |
STL-like interface to a C array of T (or part of a std container of T). More... | |
Functions | |
const char * | pull_coordinate_units (const t_pull_coord &pcrd) |
Returns the units of the pull coordinate. More... | |
double | pull_conversion_factor_userinput2internal (const t_pull_coord &pcrd) |
Returns the conversion factor from the pull coord init/rate unit to internal value unit. More... | |
double | pull_conversion_factor_internal2userinput (const t_pull_coord &pcrd) |
Returns the conversion factor from the pull coord internal value unit to the init/rate unit. More... | |
double | get_pull_coord_value (pull_t *pull, int coordIndex, const t_pbc &pbc, double t) |
Get the value for pull coord coord_ind. More... | |
double | get_pull_coord_value (pull_t *pull, int coordIndex, const t_pbc &pbc) |
Get the value for pull coord coord_ind. More... | |
void | register_external_pull_potential (struct pull_t *pull, int coord_index, const char *provider) |
Registers the provider of an external potential for a coordinate. More... | |
void | apply_external_pull_coord_force (pull_t *pull, int coord_index, double coord_force) |
Apply forces of an external potential to a pull coordinate. More... | |
void | clear_pull_forces (pull_t *pull) |
Set the all the pull forces to zero. More... | |
real | pull_potential (pull_t *pull, gmx::ArrayRef< const real > masses, const t_pbc &pbc, const t_commrec *cr, double t, real lambda, gmx::ArrayRef< const gmx::RVec > x, real *dvdlambda) |
Computes the COM pull forces, returns the potential. More... | |
void | pull_apply_forces (struct pull_t *pull, gmx::ArrayRef< const real > masses, const t_commrec *cr, gmx::ForceWithVirial *force) |
Applies the computed COM pull forces to the atoms and accumulates the virial. More... | |
void | pull_constraint (struct pull_t *pull, gmx::ArrayRef< const real > masses, const t_pbc &pbc, const t_commrec *cr, double dt, double t, gmx::ArrayRef< gmx::RVec > x, gmx::ArrayRef< gmx::RVec > xp, gmx::ArrayRef< gmx::RVec > v, tensor vir) |
Constrain the coordinates xp in the directions in x and also constrain v when v != NULL. More... | |
void | dd_make_local_pull_groups (const t_commrec *cr, pull_t *pull) |
Make a selection of the home atoms for all pull groups. Should be called at every domain decomposition. More... | |
struct pull_t * | init_pull (FILE *fplog, const pull_params_t *pull_params, const t_inputrec *ir, const gmx_mtop_t &mtop, const t_commrec *cr, gmx::LocalAtomSetManager *atomSets, real lambda) |
Allocate, initialize and return a pull work struct. More... | |
void | finish_pull (struct pull_t *pull) |
Close the pull output files and delete pull. More... | |
void | pull_calc_coms (const t_commrec *cr, pull_t *pull, gmx::ArrayRef< const real > masses, const t_pbc &pbc, double t, gmx::ArrayRef< const gmx::RVec > x, gmx::ArrayRef< gmx::RVec > xp) |
Calculates centers of mass all pull groups. More... | |
int | pullCheckPbcWithinGroups (const pull_t &pull, gmx::ArrayRef< const gmx::RVec > x, const t_pbc &pbc, real pbcMargin) |
Checks whether all groups that use a reference atom are within PBC restrictions. More... | |
bool | pullCheckPbcWithinGroup (const pull_t &pull, gmx::ArrayRef< const gmx::RVec > x, const t_pbc &pbc, int groupNr, real pbcMargin) |
Checks whether a specific group that uses a reference atom is within PBC restrictions. More... | |
bool | pull_have_potential (const pull_t &pull) |
Returns if we have pull coordinates with potential pulling. More... | |
bool | pull_have_constraint (const pull_t &pull) |
Returns if we have pull coordinates with constraint pulling. More... | |
bool | pull_have_constraint (const pull_params_t &pullParameters) |
Returns if inputrec has pull coordinates with constraint pulling. More... | |
real | max_pull_distance2 (const pull_coord_work_t &pcrd, const t_pbc &pbc) |
Returns the maxing distance for pulling. More... | |
void | updatePrevStepPullCom (pull_t *pull, std::optional< gmx::ArrayRef< double >> comPreviousStep) |
Sets the previous step COM in pull to the current COM, and optionally updates it in the provided ArrayRef. More... | |
std::vector< double > | prevStepPullCom (const pull_t *pull) |
Returns a copy of the previous step pull COM as flat vector. More... | |
void | setPrevStepPullCom (pull_t *pull, gmx::ArrayRef< const double > prevStepPullCom) |
Set the previous step pull COM from a flat vector. More... | |
void | preparePrevStepPullCom (const t_inputrec *ir, pull_t *pull_work, gmx::ArrayRef< const real > masses, t_state *state, const t_state *state_global, const t_commrec *cr, bool startingFromCheckpoint) |
Allocates, initializes and communicates the previous step pull COM (if that option is set to true). More... | |
void | initPullComFromPrevStep (const t_commrec *cr, pull_t *pull, gmx::ArrayRef< const real > masses, const t_pbc &pbc, gmx::ArrayRef< const gmx::RVec > x) |
Initializes the COM of the previous step (set to initial COM) More... | |
void | preparePrevStepPullComNewSimulation (const t_commrec *cr, pull_t *pull_work, gmx::ArrayRef< const real > masses, gmx::ArrayRef< const gmx::RVec > x, const matrix box, PbcType pbcType, std::optional< gmx::ArrayRef< double >> &&comPreviousStep) |
Initializes the previous step pull COM for new simulations (no reading from checkpoint). More... | |
Variables | |
static constexpr real | c_pullGroupPbcMargin = 0.9 |
Margin for checking pull group PBC distances compared to half the box size. | |
static constexpr real | c_pullGroupSmallGroupThreshold = 0.5 |
Threshold (as a factor of half the box size) for accepting pull groups without explicitly set refatom. | |
void apply_external_pull_coord_force | ( | pull_t * | pull, |
int | coord_index, | ||
double | coord_force | ||
) |
Apply forces of an external potential to a pull coordinate.
This function applies the external scalar force coord_force
to the pull coordinate. The corresponding potential energy value should be added to the pull or the module's potential energy term separately by the module itself. This function should be called after pull_potential() has been called and, before calling pull_apply_forces().
[in,out] | pull | The pull struct. |
[in] | coord_index | The pull coordinate index to set the force for. |
[in] | coord_force | The scalar force for the pull coordinate. |
void clear_pull_forces | ( | pull_t * | pull | ) |
Set the all the pull forces to zero.
pull | The pull group. |
void dd_make_local_pull_groups | ( | const t_commrec * | cr, |
pull_t * | pull | ||
) |
Make a selection of the home atoms for all pull groups. Should be called at every domain decomposition.
cr | Structure for communication info. |
pull | The pull group. |
void finish_pull | ( | struct pull_t * | pull | ) |
Close the pull output files and delete pull.
pull | The pull data structure. |
double get_pull_coord_value | ( | pull_t * | pull, |
int | coordIndex, | ||
const t_pbc & | pbc, | ||
double | t | ||
) |
Get the value for pull coord coord_ind.
[in,out] | pull | The pull struct. |
[in] | coordIndex | Index of the pull coordinate in the list of coordinates |
[in] | pbc | Information structure about periodicity. |
[in] | t | The time |
double get_pull_coord_value | ( | pull_t * | pull, |
int | coordIndex, | ||
const t_pbc & | pbc | ||
) |
Get the value for pull coord coord_ind.
Should only be called when time is not allowed to be used as a transformation coordinate variable. This condition is (release) asserted upon.
[in,out] | pull | The pull struct. |
[in] | coordIndex | Index of the pull coordinate in the list of coordinates |
[in] | pbc | Information structure about periodicity. |
struct pull_t* init_pull | ( | FILE * | fplog, |
const pull_params_t * | pull_params, | ||
const t_inputrec * | ir, | ||
const gmx_mtop_t & | mtop, | ||
const t_commrec * | cr, | ||
gmx::LocalAtomSetManager * | atomSets, | ||
real | lambda | ||
) |
Allocate, initialize and return a pull work struct.
fplog | General output file, normally md.log. |
pull_params | The pull input parameters containing all pull settings. |
ir | The inputrec. |
mtop | The topology of the whole system. |
cr | Struct for communication info. |
atomSets | The manager that handles the pull atom sets |
lambda | FEP lambda. |
void initPullComFromPrevStep | ( | const t_commrec * | cr, |
pull_t * | pull, | ||
gmx::ArrayRef< const real > | masses, | ||
const t_pbc & | pbc, | ||
gmx::ArrayRef< const gmx::RVec > | x | ||
) |
Initializes the COM of the previous step (set to initial COM)
[in] | cr | Struct for communication info. |
[in] | pull | The pull data structure. |
[in] | masses | Atoms masses. |
[in] | pbc | Information struct about periodicity. |
[in] | x | The local positions. |
real max_pull_distance2 | ( | const pull_coord_work_t & | pcrd, |
const t_pbc & | pbc | ||
) |
Returns the maxing distance for pulling.
For distance geometries, only dimensions with pcrd->params[dim]=1 are included in the distance calculation. For directional geometries, only dimensions with pcrd->vec[dim]!=0 are included in the distance calculation.
[in] | pcrd | Pulling data structure |
[in] | pbc | Information on periodic boundary conditions |
void preparePrevStepPullCom | ( | const t_inputrec * | ir, |
pull_t * | pull_work, | ||
gmx::ArrayRef< const real > | masses, | ||
t_state * | state, | ||
const t_state * | state_global, | ||
const t_commrec * | cr, | ||
bool | startingFromCheckpoint | ||
) |
Allocates, initializes and communicates the previous step pull COM (if that option is set to true).
If ir->pull->bSetPbcRefToPrevStepCOM is not true nothing is done.
[in] | ir | The input options/settings of the simulation. |
[in] | pull_work | The COM pull force calculation data structure |
[in] | masses | Atoms masses. |
[in] | state | The local (to this rank) state. |
[in] | state_global | The global state. |
[in] | cr | Struct for communication info. |
[in] | startingFromCheckpoint | Is the simulation starting from a checkpoint? |
void preparePrevStepPullComNewSimulation | ( | const t_commrec * | cr, |
pull_t * | pull_work, | ||
gmx::ArrayRef< const real > | masses, | ||
gmx::ArrayRef< const gmx::RVec > | x, | ||
const matrix | box, | ||
PbcType | pbcType, | ||
std::optional< gmx::ArrayRef< double >> && | comPreviousStep | ||
) |
Initializes the previous step pull COM for new simulations (no reading from checkpoint).
[in] | cr | Struct for communication info. |
[in] | pull_work | The COM pull force calculation data structure. |
[in] | masses | Atoms masses. |
[in] | x | The local positions. |
[in] | box | The current box matrix. |
[in] | pbcType | The type of periodic boundary conditions. |
[in] | comPreviousStep | The COM of the previous step of each pull group. |
std::vector<double> prevStepPullCom | ( | const pull_t * | pull | ) |
Returns a copy of the previous step pull COM as flat vector.
Used for modular simulator checkpointing. Allows to keep the implementation details of pull_t hidden from its users.
[in] | pull | The COM pull force calculation data structure |
void pull_apply_forces | ( | struct pull_t * | pull, |
gmx::ArrayRef< const real > | masses, | ||
const t_commrec * | cr, | ||
gmx::ForceWithVirial * | force | ||
) |
Applies the computed COM pull forces to the atoms and accumulates the virial.
When force!=nullptr
, distributes the pull force on the COM of each normal pull group to the atoms in the group (using mass weighting).
Also performs the recursion for transformation pull coordinates, when present, distributing the force on transformation coordinates to the COM of groups involved.
This function should be called after calling pull_potential()
and also after other modules, e.g. AWH, have called apply_external_pull_coord_force()
.
Note: this function is fully local and does not perform MPI communication.
[in,out] | pull | The pull struct. |
[in] | masses | Atoms masses. |
[in] | cr | Struct for communication info. |
[in,out] | force | Forces and virial. |
void pull_calc_coms | ( | const t_commrec * | cr, |
pull_t * | pull, | ||
gmx::ArrayRef< const real > | masses, | ||
const t_pbc & | pbc, | ||
double | t, | ||
gmx::ArrayRef< const gmx::RVec > | x, | ||
gmx::ArrayRef< gmx::RVec > | xp | ||
) |
Calculates centers of mass all pull groups.
[in] | cr | Struct for communication info. |
[in] | pull | The pull data structure. |
[in] | masses | Atoms masses. |
[in] | pbc | Information struct about periodicity. |
[in] | t | Time, only used for cylinder ref. |
[in] | x | The local positions. |
[in,out] | xp | Updated x, can be NULL. |
void pull_constraint | ( | struct pull_t * | pull, |
gmx::ArrayRef< const real > | masses, | ||
const t_pbc & | pbc, | ||
const t_commrec * | cr, | ||
double | dt, | ||
double | t, | ||
gmx::ArrayRef< gmx::RVec > | x, | ||
gmx::ArrayRef< gmx::RVec > | xp, | ||
gmx::ArrayRef< gmx::RVec > | v, | ||
tensor | vir | ||
) |
Constrain the coordinates xp in the directions in x and also constrain v when v != NULL.
[in,out] | pull | The pull data. |
[in] | masses | Atoms masses. |
[in] | pbc | Information struct about periodicity. |
[in] | cr | Struct for communication info. |
[in] | dt | The time step length. |
[in] | t | The time. |
[in] | x | Positions. |
[in,out] | xp | Updated x, can be NULL. |
[in,out] | v | Velocities, which may get a pull correction. |
[in,out] | vir | The virial, which, if != NULL, gets a pull correction. |
double pull_conversion_factor_internal2userinput | ( | const t_pull_coord & | pcrd | ) |
Returns the conversion factor from the pull coord internal value unit to the init/rate unit.
[in] | pcrd | The pull coordinate to get the conversion factor for. |
double pull_conversion_factor_userinput2internal | ( | const t_pull_coord & | pcrd | ) |
Returns the conversion factor from the pull coord init/rate unit to internal value unit.
[in] | pcrd | The pull coordinate to get the conversion factor for. |
const char* pull_coordinate_units | ( | const t_pull_coord & | pcrd | ) |
Returns the units of the pull coordinate.
[in] | pcrd | The pull coordinate to query the units for. |
bool pull_have_constraint | ( | const pull_t & | pull | ) |
Returns if we have pull coordinates with constraint pulling.
[in] | pull | The pull data structure. |
bool pull_have_constraint | ( | const pull_params_t & | pullParameters | ) |
Returns if inputrec has pull coordinates with constraint pulling.
[in] | pullParameters | Pulling input parameters from input record. |
bool pull_have_potential | ( | const pull_t & | pull | ) |
Returns if we have pull coordinates with potential pulling.
[in] | pull | The pull data structure. |
real pull_potential | ( | pull_t * | pull, |
gmx::ArrayRef< const real > | masses, | ||
const t_pbc & | pbc, | ||
const t_commrec * | cr, | ||
double | t, | ||
real | lambda, | ||
gmx::ArrayRef< const gmx::RVec > | x, | ||
real * | dvdlambda | ||
) |
Computes the COM pull forces, returns the potential.
The function computes the COMs of the pull groups and the potentials and forces acting on the pull groups, except for external potential coordinates, which forces are set by calls to apply_external_pull_coord_force()
after calling this function. To finalize the pull application, a call to pull_apply_forces()
is required to distribute the forces on the COMs to the atoms.
Note: performance global MPI communication, potentially on a subset of the MPI ranks.
[in,out] | pull | The pull struct. |
[in] | masses | Atoms masses. |
[in] | pbc | Information struct about periodicity. |
[in] | cr | Struct for communication info. |
[in] | t | Time. |
[in] | lambda | The value of lambda in FEP calculations. |
[in] | x | Positions. |
[out] | dvdlambda | Pull contribution to dV/d(lambda). |
bool pullCheckPbcWithinGroup | ( | const pull_t & | pull, |
gmx::ArrayRef< const gmx::RVec > | x, | ||
const t_pbc & | pbc, | ||
int | groupNr, | ||
real | pbcMargin | ||
) |
Checks whether a specific group that uses a reference atom is within PBC restrictions.
Groups that use a reference atom for determining PBC should have all their atoms within half the box size from the PBC atom. The box size is used per dimension for rectangular boxes, but can be a combination of dimensions for triclinic boxes, depending on which dimensions are involved in the pull coordinates a group is involved in. A margin is specified to ensure that atoms are not too close to the maximum distance. Only one group is checked.
Should be called without MPI parallelization and after pull_calc_coms() has been called at least once.
[in] | pull | The pull data structure |
[in] | x | The coordinates |
[in] | pbc | Information struct about periodicity |
[in] | groupNr | The index of the group (in pull.group[]) to check |
[in] | pbcMargin | The minimum margin (as a fraction) to half the box size |
int pullCheckPbcWithinGroups | ( | const pull_t & | pull, |
gmx::ArrayRef< const gmx::RVec > | x, | ||
const t_pbc & | pbc, | ||
real | pbcMargin | ||
) |
Checks whether all groups that use a reference atom are within PBC restrictions.
Groups that use a reference atom for determining PBC should have all their atoms within half the box size from the PBC atom. The box size is used per dimension for rectangular boxes, but can be a combination of dimensions for triclinic boxes, depending on which dimensions are involved in the pull coordinates a group is involved in. A margin is specified to ensure that atoms are not too close to the maximum distance.
Should be called without MPI parallelization and after pull_calc_coms() has been called at least once.
[in] | pull | The pull data structure |
[in] | x | The coordinates |
[in] | pbc | Information struct about periodicity |
[in] | pbcMargin | The minimum margin (as a fraction) to half the box size |
void register_external_pull_potential | ( | struct pull_t * | pull, |
int | coord_index, | ||
const char * | provider | ||
) |
Registers the provider of an external potential for a coordinate.
This function is only used for checking the consistency of the pull setup. For each pull coordinate of type external-potential, selected by the user in the mdp file, there has to be a module that provides this potential. The module registers itself as the provider by calling this function. The passed provider
string has to match the string that the user passed with the potential-provider pull coordinate mdp option. This function should be called after init_pull has been called and before pull_potential is called for the first time. This function does many consistency checks and when it returns and the first call to do_potential passes, the pull setup is guaranteed to be correct (unless the module doesn't call apply_external_pull_coord_force every step or calls it with incorrect forces). This registering function will exit with a (release) assertion failure when used incorrely or with a fatal error when the user (mdp) input in inconsistent.
Thread-safe for simultaneous registration from multiple threads.
[in,out] | pull | The pull struct. |
[in] | coord_index | The pull coordinate index to register the external potential for. |
[in] | provider | Provider string, should match the potential-provider pull coordinate mdp option. |
void setPrevStepPullCom | ( | pull_t * | pull, |
gmx::ArrayRef< const double > | prevStepPullCom | ||
) |
Set the previous step pull COM from a flat vector.
Used to restore modular simulator checkpoints. Allows to keep the implementation details of pull_t hidden from its users.
[in] | pull | The COM pull force calculation data structure |
[in] | prevStepPullCom | The previous step COM to set |
void updatePrevStepPullCom | ( | pull_t * | pull, |
std::optional< gmx::ArrayRef< double >> | comPreviousStep | ||
) |
Sets the previous step COM in pull to the current COM, and optionally updates it in the provided ArrayRef.
[in] | pull | The COM pull force calculation data structure |
[in] | comPreviousStep | The COM of the previous step of each pull group |