SSAGES  0.8.3
Software Suite for Advanced General Ensemble Simulations
BasisFunc.h
1 
21 #pragma once
22 
23 #include "Method.h"
24 #include "Utility/Basis.h"
25 #include "Grids/Grid.h"
26 #include <fstream>
27 #include <vector>
28 
29 namespace SSAGES
30 {
32 
38  class BFS : public Method
39  {
40  private:
41 
43 
47 
49 
53 
55 
59 
61 
66  std::vector<double> unbias_;
67 
69  std::vector<double> coeffArr_;
70 
72 
77 
79 
84  std::vector<double> restraint_;
85 
87  std::vector<double> boundUp_;
88 
90  std::vector<double> boundLow_;
91 
93  unsigned int cyclefreq_;
94 
96  unsigned int cyclefreqTmp_;
97 
99  unsigned int step_;
100 
102  unsigned int mpiid_;
103 
105 
110  double weight_;
111 
113 
118  double temperature_;
119 
121  double tol_;
122 
124  bool bounds_;
125 
128 
130  void ProjectBias(const CVList& cvs, const double beta);
131 
133  void InBounds(const CVList& cvs);
134 
136 
140  void PrintBias(const CVList& cvs, const double beta);
141 
144  std::string bnme_;
145 
147  unsigned int iteration_;
148 
149  public:
151 
174  BFS(const MPI_Comm& world,
175  const MPI_Comm& comm,
177  Grid<std::vector<double>> *f,
178  Grid<double> *b,
179  const std::vector<BasisFunction*>& functions,
180  const std::vector<double>& restraint,
181  const std::vector<double>& boundUp,
182  const std::vector<double>& boundLow,
183  unsigned int cyclefreq,
184  unsigned int frequency,
185  const std::string bnme,
186  const double temperature,
187  const double tol,
188  const double weight,
189  bool converge) :
190  Method(frequency, world, comm),
191  h_(h), b_(b), f_(f), unbias_(), coeffArr_(), evaluator_(functions),
192  restraint_(restraint), boundUp_(boundUp), boundLow_(boundLow),
193  cyclefreq_(cyclefreq), mpiid_(0), weight_(weight),
194  temperature_(temperature), tol_(tol),
195  convergeExit_(converge), bnme_(bnme), iteration_(0)
196  {
197  }
198 
200 
204  void PreSimulation(Snapshot* snapshot, const class CVManager& cvmanager) override;
205 
207 
211  void PostIntegration(Snapshot* snapshot, const class CVManager& cvmanager) override;
212 
214 
218  void PostSimulation(Snapshot* snapshot, const class CVManager& cvmanager) override;
219 
221 
227  void SetIteration(const int iter)
228  {
229  iteration_ = iter;
230  }
231 
233 
240  void SetBasis(const std::vector<double>&coeff, std::vector<double>&unbias)
241  {
242  evaluator_.SetCoeff(coeff);
243  unbias_ = unbias;
244  }
245 
247  static BFS* Build(const Json::Value& json,
248  const MPI_Comm& world,
249  const MPI_Comm& comm,
250  const std::string& path);
251 
254  {
255  //delete h_;
256  //delete f_;
257  //delete b_;
258  }
259  };
260 }
261 
Collective variable manager.
Definition: CVManager.h:40
Grid< std::vector< double > > * f_
Stored gradients.
Definition: BasisFunc.h:58
std::vector< CollectiveVariable * > CVList
List of Collective Variables.
Definition: types.h:51
std::vector< double > boundLow_
Lower position of the spring restraint.
Definition: BasisFunc.h:90
unsigned int step_
Step counter for the cyclefrequency.
Definition: BasisFunc.h:99
Class containing a snapshot of the current simulation in time.
Definition: Snapshot.h:43
void SetBasis(const std::vector< double > &coeff, std::vector< double > &unbias)
Set the values for the basis.
Definition: BasisFunc.h:240
void ProjectBias(const CVList &cvs, const double beta)
The functions which calculates the updated bias and coefficients and then stores them.
Definition: BasisFunc.cpp:166
void InBounds(const CVList &cvs)
A function which checks to see if the CVs are still in bounds.
Definition: BasisFunc.cpp:272
unsigned int iteration_
Iteration counter.
Definition: BasisFunc.h:147
double weight_
Weighting for potentially faster sampling.
Definition: BasisFunc.h:110
Interface for Method implementations.
Definition: Method.h:43
BFS(const MPI_Comm &world, const MPI_Comm &comm, Grid< unsigned int > *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:174
BasisEvaluator evaluator_
The basis evaluator class.
Definition: BasisFunc.h:76
void PreSimulation(Snapshot *snapshot, const class CVManager &cvmanager) override
Pre-simulation hook.
Definition: BasisFunc.cpp:36
std::vector< double > coeffArr_
The coefficient array for restart runs.
Definition: BasisFunc.h:69
unsigned int mpiid_
The node that the current system belongs to, primarily for printing and debugging.
Definition: BasisFunc.h:102
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:127
Grid< double > * b_
Stored bias potential.
Definition: BasisFunc.h:52
std::vector< double > unbias_
The biased histogram of states.
Definition: BasisFunc.h:66
Calculates the inner product of all the basis functions and the histogram.
Definition: Basis.h:332
void SetIteration(const int iter)
Set the current iteration.
Definition: BasisFunc.h:227
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:121
unsigned int cyclefreqTmp_
Temporary cyclefrequency that ensures sweeps are of proper length for restrained simulations.
Definition: BasisFunc.h:96
void PostSimulation(Snapshot *snapshot, const class CVManager &cvmanager) override
Post-simulation hook.
Definition: BasisFunc.cpp:160
void PostIntegration(Snapshot *snapshot, const class CVManager &cvmanager) override
Post-integration hook.
Definition: BasisFunc.cpp:81
bool bounds_
A variable to check to see if the simulation is in bounds or not.
Definition: BasisFunc.h:124
static BFS * Build(const Json::Value &json, const MPI_Comm &world, const MPI_Comm &comm, const std::string &path)
Build a derived method from JSON node.
Definition: BasisFunc.cpp:308
std::vector< double > restraint_
Spring constants for restrained system.
Definition: BasisFunc.h:84
Basis Function Sampling Algorithm.
Definition: BasisFunc.h:38
void SetCoeff(const std::vector< double > &coeff)
Set the coefficient array in the event of a restart run.
Definition: Basis.h:383
~BFS()
Destructor.
Definition: BasisFunc.h:253
void PrintBias(const CVList &cvs, const double beta)
Prints the current bias to a defined file from the JSON.
Definition: BasisFunc.cpp:227
std::string bnme_
Definition: BasisFunc.h:144
Grid< unsigned int > * h_
Grid of visited states.
Definition: BasisFunc.h:46
double temperature_
Self-defined temperature.
Definition: BasisFunc.h:118