SSAGES  0.1
A MetaDynamics Package
Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes | List of all members
SSAGES::BFS Class Reference

Basis Function Sampling Algorithm. More...

#include <BasisFunc.h>

Inheritance diagram for SSAGES::BFS:
Inheritance graph
[legend]

Public Member Functions

 BFS (const MPI_Comm &world, const MPI_Comm &comm, Grid< uint > *h, Grid< std::vector< double >> *f, Grid< double > *b, const std::vector< BasisFunction * > &functions, const std::vector< double > &restraint, const std::vector< double > &boundUp, const std::vector< double > &boundLow, unsigned int cyclefreq, unsigned int frequency, const std::string bnme, const double temperature, const double tol, const double weight, bool converge)
 Constructor. More...
 
void PreSimulation (Snapshot *snapshot, const class CVManager &cvmanager) override
 Pre-simulation hook. More...
 
void PostIntegration (Snapshot *snapshot, const class CVManager &cvmanager) override
 Post-integration hook. More...
 
void PostSimulation (Snapshot *snapshot, const class CVManager &cvmanager) override
 Post-simulation hook. More...
 
void SetIteration (const int iter)
 Set the current iteration. More...
 
void SetBasis (const std::vector< double > &coeff, std::vector< double > &unbias)
 Set the values for the basis. More...
 
 ~BFS ()
 Destructor.
 
- Public Member Functions inherited from SSAGES::Method
 Method (uint frequency, const MPI_Comm &world, const MPI_Comm &comm)
 Constructor. More...
 
void SetCVMask (const std::vector< uint > &mask)
 Sets the collective variable mask.
 
virtual ~Method ()
 Destructor.
 
- Public Member Functions inherited from SSAGES::EventListener
 EventListener (uint frequency)
 Constructor. More...
 
uint GetFrequency () const
 Get frequency of event listener. More...
 
virtual ~EventListener ()
 Destructor.
 

Static Public Member Functions

static BFSBuild (const Json::Value &json, const MPI_Comm &world, const MPI_Comm &comm, const std::string &path)
 
- Static Public Member Functions inherited from SSAGES::Method
static MethodBuildMethod (const Json::Value &json, const MPI_Comm &world, const MPI_Comm &comm, const std::string &path)
 Build a derived method from JSON node. More...
 

Private Member Functions

void ProjectBias (const CVList &cvs, const double beta)
 The functions which calculates the updated bias and coefficients and then stores them.
 
void InBounds (const CVList &cvs)
 A function which checks to see if the CVs are still in bounds.
 
void PrintBias (const CVList &cvs, const double beta)
 Prints the current bias to a defined file from the JSON. More...
 

Private Attributes

Grid< uint > * h_
 Grid of visited states. More...
 
Grid< double > * b_
 Stored bias potential. More...
 
Grid< std::vector< double > > * f_
 Stored gradients. More...
 
BasisEvaluator evaluator_
 The basis evaluator class. More...
 
std::vector< double > unbias_
 The biased histogram of states. More...
 
std::vector< double > coeffArr_
 The coefficient array for restart runs.
 
std::vector< double > restraint_
 Spring constants for restrained system. More...
 
std::vector< double > boundUp_
 Upper position of the spring restraint.
 
std::vector< double > boundLow_
 Lower position of the spring restraint.
 
unsigned int cyclefreq_
 Frequency of coefficient updates.
 
unsigned int mpiid_
 The node that the current system belongs to, primarily for printing and debugging.
 
double weight_
 Weighting for potentially faster sampling. More...
 
double temperature_
 Self-defined temperature. More...
 
double tol_
 The tolerance criteria for the system .
 
bool bounds_
 A variable to check to see if the simulation is in bounds or not.
 
bool convergeExit_
 A check to see if you want the system to end when it reaches the convergence criteria.
 
std::string bnme_
 
uint iteration_
 Iteration counter.
 

Additional Inherited Members

- Protected Attributes inherited from SSAGES::Method
mxx::comm world_
 Global MPI communicator.
 
mxx::comm comm_
 Local MPI communicator.
 
std::vector< uint > cvmask_
 Mask which identifies which CVs to act on.
 

Detailed Description

Basis Function Sampling Algorithm.

Implementation of the Basis Function Sampling Method based on [4].

Definition at line 38 of file BasisFunc.h.

Constructor & Destructor Documentation

SSAGES::BFS::BFS ( const MPI_Comm &  world,
const MPI_Comm &  comm,
Grid< uint > *  h,
Grid< std::vector< double >> *  f,
Grid< double > *  b,
const std::vector< BasisFunction * > &  functions,
const std::vector< double > &  restraint,
const std::vector< double > &  boundUp,
const std::vector< double > &  boundLow,
unsigned int  cyclefreq,
unsigned int  frequency,
const std::string  bnme,
const double  temperature,
const double  tol,
const double  weight,
bool  converge 
)
inline

Constructor.

Parameters
worldMPI global communicator.
commMPI local communicator.
polyordOrder of Legendre polynomials.
restraintRestraint spring constants.
boundUpUpper bounds of restraint springs.
boundLowLower bounds of restraint springs.
cyclefreqCycle frequency.
frequencyFrequency with which this Method is applied.
bnmeBasis file name.
cnmeCoefficient file name.
temperatureAutomatically set temperature.
tolThreshold for tolerance criterion.
weightWeight for improved sampling.
convergeIf True quit on convergence.

Constructs an instance of the Basis function sampling method. The coefficients describes the basis projection of the system. This is updated once every cyclefreq_. For now, only the Legendre polynomial is implemented. Others will be added later.

Definition at line 166 of file BasisFunc.h.

References PostIntegration(), PostSimulation(), and PreSimulation().

181  :
182  Method(frequency, world, comm),
183  h_(h), b_(b), f_(f), unbias_(), coeffArr_(), evaluator_(functions),
184  restraint_(restraint), boundUp_(boundUp), boundLow_(boundLow),
185  cyclefreq_(cyclefreq), mpiid_(0), weight_(weight),
186  temperature_(temperature), tol_(tol),
187  convergeExit_(converge), bnme_(bnme), iteration_(0)
188  {
189  }
Grid< std::vector< double > > * f_
Stored gradients.
Definition: BasisFunc.h:58
std::vector< double > boundLow_
Lower position of the spring restraint.
Definition: BasisFunc.h:90
double weight_
Weighting for potentially faster sampling.
Definition: BasisFunc.h:104
BasisEvaluator evaluator_
The basis evaluator class.
Definition: BasisFunc.h:65
std::vector< double > coeffArr_
The coefficient array for restart runs.
Definition: BasisFunc.h:76
unsigned int mpiid_
The node that the current system belongs to, primarily for printing and debugging.
Definition: BasisFunc.h:96
unsigned int cyclefreq_
Frequency of coefficient updates.
Definition: BasisFunc.h:93
bool convergeExit_
A check to see if you want the system to end when it reaches the convergence criteria.
Definition: BasisFunc.h:121
Grid< uint > * h_
Grid of visited states.
Definition: BasisFunc.h:46
uint iteration_
Iteration counter.
Definition: BasisFunc.h:141
Grid< double > * b_
Stored bias potential.
Definition: BasisFunc.h:52
Method(uint frequency, const MPI_Comm &world, const MPI_Comm &comm)
Constructor.
Definition: Method.h:61
std::vector< double > unbias_
The biased histogram of states.
Definition: BasisFunc.h:73
std::vector< double > boundUp_
Upper position of the spring restraint.
Definition: BasisFunc.h:87
double tol_
The tolerance criteria for the system .
Definition: BasisFunc.h:115
std::vector< double > restraint_
Spring constants for restrained system.
Definition: BasisFunc.h:84
std::string bnme_
Definition: BasisFunc.h:138
double temperature_
Self-defined temperature.
Definition: BasisFunc.h:112

Here is the call graph for this function:

Member Function Documentation

BFS * SSAGES::BFS::Build ( const Json::Value &  json,
const MPI_Comm &  world,
const MPI_Comm &  comm,
const std::string &  path 
)
static

Definition at line 302 of file BasisFunc.cpp.

References Json::Requirement::GetErrors(), SSAGES::GridBase< T >::GetNumPoints(), Json::Requirement::HasErrors(), Json::ObjectRequirement::Parse(), and Json::ObjectRequirement::Validate().

Referenced by SetBasis().

306  {
307  ObjectRequirement validator;
308  Value schema;
309  Reader reader;
310 
311  reader.parse(JsonSchema::BFSMethod, schema);
312  validator.Parse(schema, path);
313 
314  //Validate Inputs
315  validator.Validate(json, path);
316  if(validator.HasErrors())
317  throw BuildException(validator.GetErrors());
318 
319  std::vector<double> restrCV(0);
320  for(auto& restr : json["CV_restraint_spring_constants"])
321  restrCV.push_back(restr.asDouble());
322 
323  std::vector<double> boundLow(0);
324  for(auto& bndl : json["CV_restraint_minimums"])
325  boundLow.push_back(bndl.asDouble());
326 
327  std::vector<double> boundUp(0);
328  for(auto& bndu : json["CV_restraint_maximums"])
329  boundUp.push_back(bndu.asDouble());
330 
331  auto cyclefreq = json.get("cycle_frequency", 100000).asInt();
332  auto freq = json.get("frequency", 1).asInt();
333  auto wght = json.get("weight", 1.0).asDouble();
334  auto bnme = json.get("basis_filename", "").asString();
335  auto temp = json.get("temperature", 0.0).asDouble();
336  auto tol = json.get("tolerance", 1e-6).asDouble();
337  auto conv = json.get("convergence_exit", false).asBool();
338 
339  Grid<uint> *h = Grid<uint>::BuildGrid(
340  json.get("grid", Json::Value()) );
341 
342  Grid<std::vector<double>> *f = Grid<std::vector<double>>::BuildGrid(
343  json.get("grid", Json::Value()) );
344 
345  Grid<double> *b = Grid<double>::BuildGrid(
346  json.get("grid", Json::Value()) );
347 
348  size_t ii = 0;
349  std::vector<BasisFunction*> functions;
350  for(auto& m : json["basis_functions"])
351  {
352  auto *bf = BasisFunction::Build(m, path, b->GetNumPoints(ii));
353  functions.push_back(bf);
354  ii++;
355  }
356 
357 
358  auto* m = new BFS(world, comm, h, f, b, functions, restrCV, boundUp, boundLow,
359  cyclefreq, freq, bnme, temp, tol, wght,
360  conv);
361 
362  if(json.isMember("iteration"))
363  m->SetIteration(json.get("iteration",0).asInt());
364 
365  // Check if previously saved grids exist. If so, check that data match and load grids.
366  std::ifstream restrFile;
367  restrFile.open("restart"+bnme+".out");
368  if(restrFile.is_open())
369  {
370  std::string line;
371  std::vector<double> coeff;
372  std::vector<double> unbias;
373  bool coeffFlag = false;
374  bool basisFlag = false;
375  std::cout << "Attempting to load data from a previous run of BFS." << std::endl;
376 
377  while(!std::getline(restrFile,line).eof())
378  {
379  if(line.find("COEFFICIENTS") != std::string::npos) {
380  std::getline(restrFile,line);
381  coeffFlag = true;
382  }
383  else if (line.find("HISTOGRAM") != std::string::npos) {
384  std::getline(restrFile,line);
385  basisFlag = true;
386  coeffFlag = false;
387  }
388  if(coeffFlag) {
389  coeff.push_back(std::stod(line));
390  }
391  else if (basisFlag) {
392  unbias.push_back(std::stod(line));
393  }
394  }
395  m->SetBasis(coeff, unbias);
396  restrFile.close();
397  }
398 
399  return m;
400  }
bool HasErrors()
Check if errors have occured.
Definition: Requirement.h:86
virtual void Parse(Value json, const std::string &path) override
Parse JSON value to generate Requirement(s).
std::vector< std::string > GetErrors()
Get list of error messages.
Definition: Requirement.h:92
BFS(const MPI_Comm &world, const MPI_Comm &comm, Grid< uint > *h, Grid< std::vector< double >> *f, Grid< double > *b, const std::vector< BasisFunction * > &functions, const std::vector< double > &restraint, const std::vector< double > &boundUp, const std::vector< double > &boundLow, unsigned int cyclefreq, unsigned int frequency, const std::string bnme, const double temperature, const double tol, const double weight, bool converge)
Constructor.
Definition: BasisFunc.h:166
Requirements on an object.
static Grid< T > * BuildGrid(const Json::Value &json)
Set up the grid.
Definition: Grid.h:127
virtual void Validate(const Value &json, const std::string &path) override
Validate JSON value.

Here is the call graph for this function:

Here is the caller graph for this function:

void SSAGES::BFS::PostIntegration ( Snapshot snapshot,
const class CVManager cvmanager 
)
overridevirtual

Post-integration hook.

Parameters
snapshotSimulation snapshot.
cvmanagerCollective variable manager.

Implements SSAGES::Method.

Definition at line 77 of file BasisFunc.cpp.

References SSAGES::CVManager::GetCVs(), SSAGES::Snapshot::GetForces(), SSAGES::Snapshot::GetIteration(), SSAGES::Snapshot::GetKb(), and SSAGES::Snapshot::GetVirial().

Referenced by BFS().

78  {
79  auto cvs = cvmanager.GetCVs(cvmask_);
80  std::vector<double> x(cvs.size(),0);
81 
82  /*The binned cv space is updated at every step
83  *After a certain number of steps has been passed, the system updates a
84  *bias projection based on the visited histogram states
85  */
86  for(size_t i = 0; i < cvs.size(); ++i)
87  {
88  x[i] = cvs[i]->GetValue();
89  }
90 
91  // The histogram is updated based on the index
92  InBounds(cvs);
93  if (bounds_)
94  h_->at(x) += 1;
95 
96  // Update the basis projection after a predefined number of steps
97  if(snapshot->GetIteration() % cyclefreq_ == 0) {
98  double beta;
99  beta = 1.0 / (temperature_ * snapshot->GetKb());
100 
101  // For systems with poorly defined temperature (ie: 1 particle) the
102  // user needs to define their own temperature. This is a hack that
103  // will be removed in future versions.
104 
105  iteration_+= 1;
106  ProjectBias(cvs,beta);
107  std::cout<<"Node: ["<<mpiid_<<"]"<<std::setw(10)<<"\tSweep: "<<iteration_<<std::endl;
108  }
109 
110  // This gets the bias force from the grid
111  std::vector<double> bias_grad(cvs.size(),0);
112 
113  if (bounds_)
114  {
115  auto& bias = f_->at(x);
116  for (size_t ii = 0; ii<bias.size(); ii++)
117  bias_grad[ii] = bias[ii];
118  }
119  else
120  {
121  // This is where the wall potentials are going to be thrown into the method if the system is not a periodic CV
122  for(size_t j = 0; j < cvs.size(); ++j)
123  {
124  if(!h_->GetPeriodic(j))
125  {
126  if(x[j] > boundUp_[j])
127  bias_grad[j] = -restraint_[j] * (x[j] - boundUp_[j]);
128  else if(x[j] < boundLow_[j])
129  bias_grad[j] = restraint_[j] * (x[j] - boundLow_[j]);
130  }
131  }
132  }
133 
134  // Take each CV and add its biased forces to the atoms using the chain rule
135  auto& forces = snapshot->GetForces();
136  auto& virial = snapshot->GetVirial();
137  for(size_t i = 0; i < cvs.size(); ++i)
138  {
139  auto& grad = cvs[i]->GetGradient();
140  auto& boxgrad = cvs[i]->GetBoxGradient();
141 
142  /* Update the forces in snapshot by adding in the force bias from each
143  *CV to each atom based on the gradient of the CV.
144  */
145  for (size_t j = 0; j < forces.size(); ++j)
146  forces[j] += bias_grad[i]*grad[j];
147 
148  virial -= bias_grad[i]*boxgrad;
149  }
150  }
Grid< std::vector< double > > * f_
Stored gradients.
Definition: BasisFunc.h:58
std::vector< double > boundLow_
Lower position of the spring restraint.
Definition: BasisFunc.h:90
void ProjectBias(const CVList &cvs, const double beta)
The functions which calculates the updated bias and coefficients and then stores them.
Definition: BasisFunc.cpp:159
void InBounds(const CVList &cvs)
A function which checks to see if the CVs are still in bounds.
Definition: BasisFunc.cpp:265
unsigned int mpiid_
The node that the current system belongs to, primarily for printing and debugging.
Definition: BasisFunc.h:96
unsigned int cyclefreq_
Frequency of coefficient updates.
Definition: BasisFunc.h:93
Grid< uint > * h_
Grid of visited states.
Definition: BasisFunc.h:46
uint iteration_
Iteration counter.
Definition: BasisFunc.h:141
const std::vector< bool > & GetPeriodic() const
Return the periodicity of the Grid.
Definition: GridBase.h:291
std::vector< double > boundUp_
Upper position of the spring restraint.
Definition: BasisFunc.h:87
bool bounds_
A variable to check to see if the simulation is in bounds or not.
Definition: BasisFunc.h:118
std::vector< uint > cvmask_
Mask which identifies which CVs to act on.
Definition: Method.h:50
std::vector< double > restraint_
Spring constants for restrained system.
Definition: BasisFunc.h:84
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:536
double temperature_
Self-defined temperature.
Definition: BasisFunc.h:112

Here is the call graph for this function:

Here is the caller graph for this function:

void SSAGES::BFS::PostSimulation ( Snapshot snapshot,
const class CVManager cvmanager 
)
overridevirtual

Post-simulation hook.

Parameters
snapshotSimulation snapshot.
cvmanagerCollective variable manager.

Implements SSAGES::Method.

Definition at line 153 of file BasisFunc.cpp.

Referenced by BFS().

154  {
155  std::cout<<"Run has finished"<<std::endl;
156  }

Here is the caller graph for this function:

void SSAGES::BFS::PreSimulation ( Snapshot snapshot,
const class CVManager cvmanager 
)
overridevirtual

Pre-simulation hook.

Parameters
snapshotSimulation snapshot.
cvmanagerCollective variable manager.

Implements SSAGES::Method.

Definition at line 36 of file BasisFunc.cpp.

References SSAGES::CVManager::GetCVs(), and SSAGES::Snapshot::GetWalkerID().

Referenced by BFS().

37  {
38  auto cvs = cvmanager.GetCVs(cvmask_);
39  // For print statements and file I/O, the walker IDs are used
40  mpiid_ = snapshot->GetWalkerID();
41 
42  // Make sure the iteration index is set correctly
43  iteration_ = 0;
44 
45  // There are a few error messages / checks that are in place with
46  // defining CVs and grids
47  if(h_->GetDimension() != cvs.size())
48  {
49  std::cerr<<"ERROR: Grid dimensions doesn't match number of CVS."<<std::endl;
50  MPI_Abort(world_, EXIT_FAILURE);
51  }
52 
53  // This is to check for non-periodic bounds. It comes into play in the update bias function
54  bounds_ = true;
55  unbias_.resize(h_->size(),0);
56  std::fill(unbias_.begin(),unbias_.end(),1.0);
57 
58  // Initialize the mapping for the hist function
59  for (auto &val : *h_) {
60  val = 0;
61  }
62 
63  // Reset the values in the force grid
64  for (auto &val : *f_) {
65  for (size_t i = 0; i <val.size(); i++) {
66  val[i] = 0;
67  }
68  }
69 
70  // Reset the values in the bias potential grid
71  for (auto &val : *b_) {
72  val = 0;
73  }
74  }
Grid< std::vector< double > > * f_
Stored gradients.
Definition: BasisFunc.h:58
mxx::comm world_
Global MPI communicator.
Definition: Method.h:46
unsigned int mpiid_
The node that the current system belongs to, primarily for printing and debugging.
Definition: BasisFunc.h:96
Grid< uint > * h_
Grid of visited states.
Definition: BasisFunc.h:46
uint iteration_
Iteration counter.
Definition: BasisFunc.h:141
size_t size() const
Get the size of the internal storage vector.
Definition: GridBase.h:322
Grid< double > * b_
Stored bias potential.
Definition: BasisFunc.h:52
std::vector< double > unbias_
The biased histogram of states.
Definition: BasisFunc.h:73
bool bounds_
A variable to check to see if the simulation is in bounds or not.
Definition: BasisFunc.h:118
size_t GetDimension() const
Get the dimension.
Definition: GridBase.h:186
std::vector< uint > cvmask_
Mask which identifies which CVs to act on.
Definition: Method.h:50

Here is the call graph for this function:

Here is the caller graph for this function:

void SSAGES::BFS::PrintBias ( const CVList cvs,
const double  beta 
)
private

Prints the current bias to a defined file from the JSON.

Parameters
cvsList of collective variables.
betaScale parameter.

Definition at line 220 of file BasisFunc.cpp.

References SSAGES::Grid< T >::begin().

221  {
222  // The filenames will have a standard name, with a user-defined suffix
223  std::string filename1 = "basis"+bnme_+".out";
224  std::string filename2 = "restart"+bnme_+".out";
225  std::ofstream basisout;
226  std::ofstream coeffout;
227  basisout.precision(5);
228  coeffout.precision(5);
229  basisout.open(filename1.c_str());
230  coeffout.open(filename2.c_str());
231 
232  // The CV values, PMF projection, PMF, and biased histogram are output for the user
233  basisout << "The information stored in this file is the output of a Basis Function Sampling run" << std::endl;
234  basisout << "CV Values" << std::setw(15*cvs.size()) << "Basis Set Bias" << std::setw(15) << "PMF Estimate" << std::endl;
235  coeffout << "The information stored in this file is for the purpose of restarting simulations with BFS" << std::endl;
236  coeffout << "***COEFFICIENTS***" << std::endl;
237 
238  size_t j = 0;
239  for(Grid<double>::iterator it = b_->begin(); it != b_->end(); ++it, ++j)
240  {
241  for(size_t k = 0; k < cvs.size(); ++k)
242  {
243  // Evaluate the CV values for printing purposes
244  basisout << it.coordinate(k) << std::setw(15);
245  }
246  basisout << -(*it) << std::setw(15) <<
247  -1.0/beta*log(unbias_[j]) <<
248  std::endl;
249  }
250 
251  std::vector<double> coeff = evaluator_.GetCoeff();
252  for (auto& val : coeff) {
253  coeffout << val << std::endl;
254  }
255  coeffout << "***HISTOGRAM***" << std::endl;
256  for (auto& val : unbias_) {
257  coeffout << val << std::endl;
258  }
259 
260  basisout.close();
261  coeffout.close();
262  }
BasisEvaluator evaluator_
The basis evaluator class.
Definition: BasisFunc.h:65
GridIterator< double > iterator
Custom iterator over a grid.
Definition: Grid.h:512
iterator begin()
Return iterator at first grid point.
Definition: Grid.h:524
Grid< double > * b_
Stored bias potential.
Definition: BasisFunc.h:52
std::vector< double > unbias_
The biased histogram of states.
Definition: BasisFunc.h:73
std::string bnme_
Definition: BasisFunc.h:138
iterator end()
Return iterator after last valid grid point.
Definition: Grid.h:537

Here is the call graph for this function:

void SSAGES::BFS::SetBasis ( const std::vector< double > &  coeff,
std::vector< double > &  unbias 
)
inline

Set the values for the basis.

Parameters
coeffList of coefficients.
unbiasList of unbiased values.

This function is used to set starting values at the beginning of a run. For example when continuing from a restart value.

Definition at line 232 of file BasisFunc.h.

References Build().

233  {
234  evaluator_.SetCoeff(coeff);
235  unbias_ = unbias;
236  }
BasisEvaluator evaluator_
The basis evaluator class.
Definition: BasisFunc.h:65
std::vector< double > unbias_
The biased histogram of states.
Definition: BasisFunc.h:73

Here is the call graph for this function:

void SSAGES::BFS::SetIteration ( const int  iter)
inline

Set the current iteration.

Parameters
iterNew value for the current iteration.

This function is used to set the current iteration, for example when continuing from a restart.

Definition at line 219 of file BasisFunc.h.

220  {
221  iteration_ = iter;
222  }
uint iteration_
Iteration counter.
Definition: BasisFunc.h:141

Member Data Documentation

Grid<double>* SSAGES::BFS::b_
private

Stored bias potential.

The sum of basis functions that adds up to the bias potential of the surface

Definition at line 52 of file BasisFunc.h.

std::string SSAGES::BFS::bnme_
private

The option to name both the basis and coefficient files will be given Basis filename

Definition at line 138 of file BasisFunc.h.

BasisEvaluator SSAGES::BFS::evaluator_
private

The basis evaluator class.

A class which holds the basis functions and updates the coefficients at each cycle

Definition at line 65 of file BasisFunc.h.

Grid<std::vector<double> >* SSAGES::BFS::f_
private

Stored gradients.

The gradients corresponding to each point on the grid

Definition at line 58 of file BasisFunc.h.

Grid<uint>* SSAGES::BFS::h_
private

Grid of visited states.

Grid is stored locally.

Definition at line 46 of file BasisFunc.h.

std::vector<double> SSAGES::BFS::restraint_
private

Spring constants for restrained system.

The system uses this to determine if the system is to be restrained on the defined interval. The user inputs the spring constants if the system is not periodic.

Definition at line 84 of file BasisFunc.h.

double SSAGES::BFS::temperature_
private

Self-defined temperature.

In the case of the MD engine using a poorly defined temperature, the option to throw it into the method is available. If not provided it takes the value from the engine.

Definition at line 112 of file BasisFunc.h.

std::vector<double> SSAGES::BFS::unbias_
private

The biased histogram of states.

The biased histogram of states has the form hist_*exp(phi*beta), where phi is the bias potential and beta is the inverse of the temperature. It is defined globally.

Definition at line 73 of file BasisFunc.h.

double SSAGES::BFS::weight_
private

Weighting for potentially faster sampling.

Weighting can be used to potentially sample faster, however it can cause the simulation to explode. By default this value will be set to 1.0

Definition at line 104 of file BasisFunc.h.


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