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

#include <gromacs/applied_forces/awh/bias.h>

Description

A bias acting on a multidimensional coordinate.

At each step AWH should provide its biases with updated values of their coordinates. Each bias provides AWH with an updated bias forces and the corresponding potential.

See the user manual for details on the algorithm and equations.

The bias is responsible for keeping and updating a free energy estimate along the coordinate. The bias potential is basically a function of the free energy estimate and so also changes by the update. The free energy update is based on information from coordinate samples collected at a constant bias potential, between updates.

The bias keeps a grid with coordinate points that organizes spatial information about the coordinate. The grid has the same geometry as the coordinate, i.e. they have the same dimensionality and periodicity (if any). The number of points in the grid sets the resolution of the collected data and its extent defines the sampling region of interest.

Each coordinate point has further statistical properties and function values which a grid point does not know about. E.g., for the bias each coordinate point is associated with values of the bias, free energy and target distribution, accumulated sampling weight, etc. For this the bias attaches to each grid point a state. The grid + vector of point states are the bias coordinate points.

The bias has a fairly complex global state keeping track of where the system (coordinate) currently is (CoordState), where it has sampled since the last update (BiasState) and controlling the free energy convergence rate (HistogramSize).

Partly, the complexity comes from the bias having two convergence stages: an initial stage which in an heuristic, non-deterministic way restricts the early convergence rate for sake of robustness; and a final stage where the convergence rate is constant. The length of the initial stage depends on the sampling and is unknown beforehand.

Another complexity comes from the fact that coordinate points, for sake of efficiency in the case of many grid points, are typically only accessed in recently sampled regions even though the free energy update is inherently global and affects all points. The bias allows points thay are non-local at the time the update was issued to postpone ("skip", as it is called in the code) the update. A non-local point is defined as a point which has not been sampled since the last update. Local points are points that have been sampled since the last update. The (current) set of local points are kept track of by the bias state and reset after every update. An update is called local if it only updates local points. Non-local points will temporarily "skip" the update until next time they are local (or when a global update is issued). For this to work, the bias keeps a global "clock" (in HistogramSize) of the number of issued updates. Each PointState also has its own local "clock" with the counting the number of updates it has pulled through. When a point updates its state it asserts that its local clock is synchronized with the global clock.

Public Types

enum  ThisRankWillDoIO { ThisRankWillDoIO::No, ThisRankWillDoIO::Yes }
 Enum for requesting Bias set up with(out) I/O on this rank. More...
 

Public Member Functions

 Bias (int biasIndexInCollection, const AwhParams &awhParams, const AwhBiasParams &awhBiasParams, ArrayRef< const DimParams > dimParams, double beta, double mdTimeStep, const BiasSharing *biasSharing, const std::string &biasInitFilename, ThisRankWillDoIO thisRankWillDoIO, BiasParams::DisableUpdateSkips disableUpdateSkips=BiasParams::DisableUpdateSkips::no)
 Constructor. More...
 
void printInitializationToLog (FILE *fplog) const
 Print information about initialization to log file. More...
 
gmx::ArrayRef< const double > calcForceAndUpdateBias (const awh_dvec coordValue, ArrayRef< const double > neighborLambdaEnergies, ArrayRef< const double > neighborLambdaDhdl, double *awhPotential, double *potentialJump, double t, int64_t step, int64_t seed, FILE *fplog)
 Evolves the bias at every step. More...
 
double calcConvolvedBias (const awh_dvec &coordValue) const
 Calculates the convolved bias for a given coordinate value. More...
 
void restoreStateFromHistory (const AwhBiasHistory *biasHistory, const t_commrec *cr)
 Restore the bias state from history on the main rank and broadcast it. More...
 
void initHistoryFromState (AwhBiasHistory *biasHistory) const
 Allocate and initialize a bias history with the given bias state. More...
 
void updateHistory (AwhBiasHistory *biasHistory) const
 Update the bias history with the current state. More...
 
void doSkippedUpdatesForAllPoints ()
 Do all previously skipped updates. Public for use by tests.
 
int ndim () const
 Returns the dimensionality of the bias.
 
ArrayRef< const DimParamsdimParams () const
 Returns the dimension parameters.
 
const BiasParamsparams () const
 Returns the bias parameters.
 
const BiasStatestate () const
 Returns the global state of the bias.
 
int biasIndex () const
 Returns the index of the bias.
 
const awh_dvecgetGridCoordValue (size_t gridPointIndex) const
 Return the coordinate value for a grid point. More...
 
void updateBiasStateSharedCorrelationTensorTimeIntegral ()
 Update the correlation tensor time integral shared across multiple AWH walkers.
 
const CorrelationGridforceCorrelationGrid () const
 Return a const reference to the force correlation grid.
 
int numEnergySubblocksToWrite () const
 Return the number of data blocks that have been prepared for writing.
 
int writeToEnergySubblocks (t_enxsubblock *subblock) const
 Write bias data blocks to energy subblocks. More...
 
bool hasFepLambdaDimension () const
 Returns true if the bias has a free energy lambda state dimension at all.
 
bool isSampleCoordStep (int64_t step) const
 Returns whether we should sample the coordinate. More...
 

Member Enumeration Documentation

Enum for requesting Bias set up with(out) I/O on this rank.

Enumerator
No 

This rank will not do I/O.

Yes 

This rank will do I/O.

Constructor & Destructor Documentation

gmx::Bias::Bias ( int  biasIndexInCollection,
const AwhParams awhParams,
const AwhBiasParams &  awhBiasParams,
ArrayRef< const DimParams dimParams,
double  beta,
double  mdTimeStep,
const BiasSharing *  biasSharing,
const std::string &  biasInitFilename,
ThisRankWillDoIO  thisRankWillDoIO,
BiasParams::DisableUpdateSkips  disableUpdateSkips = BiasParams::DisableUpdateSkips::no 
)

Constructor.

Parameters
[in]biasIndexInCollectionIndex of the bias in collection.
[in]awhParamsAWH parameters.
[in]awhBiasParamsBias parameters.
[in]dimParamsBias dimension parameters.
[in]beta1/(k_B T).
[in]mdTimeStepThe MD time step.
[in]biasSharingPointer to object for sharing bias over simulations, can be nullptr
[in]biasInitFilenameName of file to read PMF and target from.
[in]thisRankWillDoIOTells whether this MPI rank will do I/O (checkpointing, AWH output), normally (only) the main rank does I/O.
[in]disableUpdateSkipsIf to disable update skips, useful for testing.

Member Function Documentation

double gmx::Bias::calcConvolvedBias ( const awh_dvec coordValue) const
inline

Calculates the convolved bias for a given coordinate value.

The convolved bias is the effective bias acting on the coordinate. Since the bias here has arbitrary normalization, this only makes sense as a relative, to other coordinate values, measure of the bias.

Parameters
[in]coordValueThe coordinate value.
Returns
the convolved bias >= -GMX_FLOAT_MAX.
gmx::ArrayRef< const double > gmx::Bias::calcForceAndUpdateBias ( const awh_dvec  coordValue,
ArrayRef< const double >  neighborLambdaEnergies,
ArrayRef< const double >  neighborLambdaDhdl,
double *  awhPotential,
double *  potentialJump,
double  t,
int64_t  step,
int64_t  seed,
FILE *  fplog 
)

Evolves the bias at every step.

At each step the bias step needs to:

  • set the bias force and potential;
  • update the free energy and bias if needed;
  • reweight samples to extract the PMF.
Parameters
[in]coordValueThe current coordinate value(s).
[in]neighborLambdaEnergiesAn array containing the energy of the system in neighboring lambdas. The array is of length numLambdas+1, where numLambdas is the number of free energy lambda states. Element 0 in the array is the energy of the current state and elements 1..numLambdas contain the energy of the system in the neighboring lambda states (also including the current state). When there are no free energy lambda state dimensions this can be empty.
[in]neighborLambdaDhdlAn array containing the dHdL at the neighboring lambda points. The array is of length numLambdas+1, where numLambdas is the number of free energy lambda states. Element 0 in the array is the dHdL of the current state and elements 1..numLambdas contain the dHdL of the system in the neighboring lambda states (also including the current state). When there are no free energy lambda state dimensions this can be empty.
[out]awhPotentialBias potential.
[out]potentialJumpChange in bias potential for this bias.
[in]tTime.
[in]stepTime step.
[in]seedRandom seed.
[in,out]fplogLog file.
Returns
a reference to the bias force, size ndim(), valid until the next call of this method or destruction of Bias, whichever comes first.
const awh_dvec& gmx::Bias::getGridCoordValue ( size_t  gridPointIndex) const
inline

Return the coordinate value for a grid point.

Parameters
[in]gridPointIndexThe index of the grid point.
void gmx::Bias::initHistoryFromState ( AwhBiasHistory *  biasHistory) const

Allocate and initialize a bias history with the given bias state.

This function will be called at the start of a new simulation. Note that only constant data will be initialized here. History data is set by updateHistory.

Parameters
[in,out]biasHistoryAWH history to initialize.
bool gmx::Bias::isSampleCoordStep ( int64_t  step) const

Returns whether we should sample the coordinate.

Parameters
[in]stepThe MD step number.
void gmx::Bias::printInitializationToLog ( FILE *  fplog) const

Print information about initialization to log file.

Prints information about AWH variables that are set internally but might be of interest to the user.

Parameters
[in,out]fplogLog file, can be nullptr.
void gmx::Bias::restoreStateFromHistory ( const AwhBiasHistory *  biasHistory,
const t_commrec *  cr 
)

Restore the bias state from history on the main rank and broadcast it.

Parameters
[in]biasHistoryBias history struct, only allowed to be nullptr on worker ranks.
[in]crThe communication record.
void gmx::Bias::updateHistory ( AwhBiasHistory *  biasHistory) const

Update the bias history with the current state.

Parameters
[out]biasHistoryBias history struct.
int gmx::Bias::writeToEnergySubblocks ( t_enxsubblock *  subblock) const

Write bias data blocks to energy subblocks.

Parameters
[in,out]subblockEnergy subblocks to write to.
Returns
the number of subblocks written.

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