Gromacs  2020.4
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Functions
gmx::anonymous_namespace{biasstate.cpp} Namespace Reference

Functions

void sumOverSimulations (gmx::ArrayRef< int > arrayRef, const gmx_multisim_t *multiSimComm)
 Sum an array over all simulations on the master rank of each simulation. More...
 
void sumOverSimulations (gmx::ArrayRef< double > arrayRef, const gmx_multisim_t *multiSimComm)
 Sum an array over all simulations on the master rank of each simulation. More...
 
template<typename T >
void sumOverSimulations (gmx::ArrayRef< T > arrayRef, const t_commrec *commRecord, const gmx_multisim_t *multiSimComm)
 Sum an array over all simulations on all ranks of each simulation. More...
 
void sumPmf (gmx::ArrayRef< PointState > pointState, int numSharedUpdate, const t_commrec *commRecord, const gmx_multisim_t *multiSimComm)
 Sum PMF over multiple simulations, when requested. More...
 
double freeEnergyMinimumValue (gmx::ArrayRef< const PointState > pointState)
 Find the minimum free energy value. More...
 
double biasedLogWeightFromPoint (const std::vector< DimParams > &dimParams, const std::vector< PointState > &points, const Grid &grid, int pointIndex, double pointBias, const awh_dvec value)
 Find and return the log of the probability weight of a point given a coordinate value. More...
 
void updateTargetDistribution (gmx::ArrayRef< PointState > pointState, const BiasParams &params)
 Updates the target distribution for all points. More...
 
std::string gridPointValueString (const Grid &grid, int point)
 Puts together a string describing a grid point. More...
 
void setHistogramUpdateScaleFactors (const BiasParams &params, double newHistogramSize, double oldHistogramSize, double *weightHistScaling, double *logPmfSumScaling)
 Sets the histogram rescaling factors needed to control the histogram size. More...
 
void mergeSharedUpdateLists (std::vector< int > *updateList, int numPoints, const t_commrec *commRecord, const gmx_multisim_t *multiSimComm)
 Merge update lists from multiple sharing simulations. More...
 
void makeLocalUpdateList (const Grid &grid, const std::vector< PointState > &points, const awh_ivec originUpdatelist, const awh_ivec endUpdatelist, std::vector< int > *updateList)
 Generate an update list of points sampled since the last update. More...
 
void sumHistograms (gmx::ArrayRef< PointState > pointState, gmx::ArrayRef< double > weightSumCovering, int numSharedUpdate, const t_commrec *commRecord, const gmx_multisim_t *multiSimComm, const std::vector< int > &localUpdateList)
 Add partial histograms (accumulating between updates) to accumulating histograms. More...
 
void labelCoveredPoints (const std::vector< bool > &visited, const std::vector< bool > &checkCovering, int numPoints, int period, int coverRadius, gmx::ArrayRef< int > covered)
 Label points along an axis as covered or not. More...
 

Function Documentation

double gmx::anonymous_namespace{biasstate.cpp}::biasedLogWeightFromPoint ( const std::vector< DimParams > &  dimParams,
const std::vector< PointState > &  points,
const Grid grid,
int  pointIndex,
double  pointBias,
const awh_dvec  value 
)

Find and return the log of the probability weight of a point given a coordinate value.

The unnormalized weight is given by w(point|value) = exp(bias(point) - U(value,point)), where U is a harmonic umbrella potential.

Parameters
[in]dimParamsThe bias dimensions parameters
[in]pointsThe point state.
[in]gridThe grid.
[in]pointIndexPoint to evaluate probability weight for.
[in]pointBiasBias for the point (as a log weight).
[in]valueCoordinate value.
Returns
the log of the biased probability weight.
double gmx::anonymous_namespace{biasstate.cpp}::freeEnergyMinimumValue ( gmx::ArrayRef< const PointState >  pointState)

Find the minimum free energy value.

Parameters
[in]pointStateThe state of the points.
Returns
the minimum free energy value.
std::string gmx::anonymous_namespace{biasstate.cpp}::gridPointValueString ( const Grid grid,
int  point 
)

Puts together a string describing a grid point.

Parameters
[in]gridThe grid.
[in]pointGrid point index.
Returns
a string for the point.
void gmx::anonymous_namespace{biasstate.cpp}::labelCoveredPoints ( const std::vector< bool > &  visited,
const std::vector< bool > &  checkCovering,
int  numPoints,
int  period,
int  coverRadius,
gmx::ArrayRef< int >  covered 
)

Label points along an axis as covered or not.

A point is covered if it is surrounded by visited points up to a radius = coverRadius.

Parameters
[in]visitedVisited? For each point.
[in]checkCoveringCheck for covering? For each point.
[in]numPointsThe number of grid points along this dimension.
[in]periodPeriod in number of points.
[in]coverRadiusCover radius, in points, needed for defining a point as covered.
[in,out]coveredIn this array elements are 1 for covered points and 0 for non-covered points, this routine assumes that covered has at least size numPoints.
void gmx::anonymous_namespace{biasstate.cpp}::makeLocalUpdateList ( const Grid grid,
const std::vector< PointState > &  points,
const awh_ivec  originUpdatelist,
const awh_ivec  endUpdatelist,
std::vector< int > *  updateList 
)

Generate an update list of points sampled since the last update.

Parameters
[in]gridThe AWH bias.
[in]pointsThe point state.
[in]originUpdatelistThe origin of the rectangular region that has been sampled since last update.
[in]endUpdatelistThe end of the rectangular that has been sampled since last update.
[in,out]updateListLocal update list to set (assumed >= npoints long).
void gmx::anonymous_namespace{biasstate.cpp}::mergeSharedUpdateLists ( std::vector< int > *  updateList,
int  numPoints,
const t_commrec *  commRecord,
const gmx_multisim_t multiSimComm 
)

Merge update lists from multiple sharing simulations.

Parameters
[in,out]updateListUpdate list for this simulation (assumed >= npoints long).
[in]numPointsTotal number of points.
[in]commRecordStruct for intra-simulation communication.
[in]multiSimCommStruct for multi-simulation communication.
void gmx::anonymous_namespace{biasstate.cpp}::setHistogramUpdateScaleFactors ( const BiasParams &  params,
double  newHistogramSize,
double  oldHistogramSize,
double *  weightHistScaling,
double *  logPmfSumScaling 
)

Sets the histogram rescaling factors needed to control the histogram size.

For sake of robustness, the reference weight histogram can grow at a rate different from the actual sampling rate. Typically this happens for a limited initial time, alternatively growth is scaled down by a constant factor for all times. Since the size of the reference histogram sets the size of the free energy update this should be reflected also in the PMF. Thus the PMF histogram needs to be rescaled too.

This function should only be called by the bias update function or wrapped by a function that knows what scale factors should be applied when, e.g, getSkippedUpdateHistogramScaleFactors().

Parameters
[in]paramsThe bias parameters.
[in]newHistogramSizeNew reference weight histogram size.
[in]oldHistogramSizePrevious reference weight histogram size (before adding new samples).
[out]weightHistScalingScaling factor for the reference weight histogram.
[out]logPmfSumScalingLog of the scaling factor for the PMF histogram.
void gmx::anonymous_namespace{biasstate.cpp}::sumHistograms ( gmx::ArrayRef< PointState >  pointState,
gmx::ArrayRef< double >  weightSumCovering,
int  numSharedUpdate,
const t_commrec *  commRecord,
const gmx_multisim_t multiSimComm,
const std::vector< int > &  localUpdateList 
)

Add partial histograms (accumulating between updates) to accumulating histograms.

Parameters
[in,out]pointStateThe state of the points in the bias.
[in,out]weightSumCoveringThe weights for checking covering.
[in]numSharedUpdateThe number of biases sharing the histrogram.
[in]commRecordStruct for intra-simulation communication.
[in]multiSimCommStruct for multi-simulation communication.
[in]localUpdateListList of points with data.
void gmx::anonymous_namespace{biasstate.cpp}::sumOverSimulations ( gmx::ArrayRef< int >  arrayRef,
const gmx_multisim_t multiSimComm 
)

Sum an array over all simulations on the master rank of each simulation.

Parameters
[in,out]arrayRefThe data to sum.
[in]multiSimCommStruct for multi-simulation communication.
void gmx::anonymous_namespace{biasstate.cpp}::sumOverSimulations ( gmx::ArrayRef< double >  arrayRef,
const gmx_multisim_t multiSimComm 
)

Sum an array over all simulations on the master rank of each simulation.

Parameters
[in,out]arrayRefThe data to sum.
[in]multiSimCommStruct for multi-simulation communication.
template<typename T >
void gmx::anonymous_namespace{biasstate.cpp}::sumOverSimulations ( gmx::ArrayRef< T >  arrayRef,
const t_commrec *  commRecord,
const gmx_multisim_t multiSimComm 
)

Sum an array over all simulations on all ranks of each simulation.

This assumes the data is identical on all ranks within each simulation.

Parameters
[in,out]arrayRefThe data to sum.
[in]commRecordStruct for intra-simulation communication.
[in]multiSimCommStruct for multi-simulation communication.
void gmx::anonymous_namespace{biasstate.cpp}::sumPmf ( gmx::ArrayRef< PointState >  pointState,
int  numSharedUpdate,
const t_commrec *  commRecord,
const gmx_multisim_t multiSimComm 
)

Sum PMF over multiple simulations, when requested.

Parameters
[in,out]pointStateThe state of the points in the bias.
[in]numSharedUpdateThe number of biases sharing the histogram.
[in]commRecordStruct for intra-simulation communication.
[in]multiSimCommStruct for multi-simulation communication.
void gmx::anonymous_namespace{biasstate.cpp}::updateTargetDistribution ( gmx::ArrayRef< PointState >  pointState,
const BiasParams &  params 
)

Updates the target distribution for all points.

The target distribution is always updated for all points at the same time.

Parameters
[in,out]pointStateThe state of all points.
[in]paramsThe bias parameters.