Collective Variables Module - Developer Documentation
colvarbias_meta.h
1 // -*- c++ -*-
2 
3 #ifndef COLVARBIAS_META_H
4 #define COLVARBIAS_META_H
5 
6 #include <vector>
7 #include <list>
8 #include <sstream>
9 #include <fstream>
10 
11 #include "colvarbias.h"
12 #include "colvargrid.h"
13 
16  : public virtual colvarbias,
17  public virtual colvarbias_ti
18 {
19 
20 public:
21 
28  };
29 
32 
33  colvarbias_meta(char const *key);
34  virtual ~colvarbias_meta();
35 
36  virtual int init(std::string const &conf);
37  virtual int init_well_tempered_params(std::string const &conf);
38  virtual int init_ebmeta_params(std::string const &conf);
39 
40  virtual int clear_state_data();
41 
42  virtual int update();
43  virtual int update_grid_params();
44  virtual int update_bias();
45  virtual int update_grid_data();
46  virtual int replica_share();
47 
48  virtual int calc_energy(std::vector<colvarvalue> const &values =
49  std::vector<colvarvalue>(0));
50  virtual int calc_forces(std::vector<colvarvalue> const &values =
51  std::vector<colvarvalue>(0));
52 
53  virtual std::string const get_state_params() const;
54  virtual int set_state_params(std::string const &state_conf);
55  virtual std::ostream & write_state_data(std::ostream &os);
56  virtual std::istream & read_state_data(std::istream &os);
57 
58  virtual int setup_output();
59  virtual int write_output_files();
60  virtual void write_pmf();
61  virtual int write_state_to_replicas();
62 
63  class hill;
64  typedef std::list<hill>::iterator hill_iter;
65 
66 protected:
67 
73 
75  size_t new_hill_freq;
76 
80  std::ostream *hills_traj_os;
81 
83  std::string const hills_traj_file_name() const;
84 
87  std::list<hill> hills;
88 
91  hill_iter new_hills_begin;
92 
95  std::list<hill> hills_off_grid;
96 
99 
101  void recount_hills_off_grid(hill_iter h_first, hill_iter h_last,
102  colvar_grid_scalar *ge);
103 
105  std::istream & read_hill(std::istream &is);
106 
110  virtual std::list<hill>::const_iterator create_hill(hill const &h);
111 
114  virtual std::list<hill>::const_iterator delete_hill(hill_iter &h);
115 
118  virtual void calc_hills(hill_iter h_first,
119  hill_iter h_last,
120  cvm::real &energy,
121  std::vector<colvarvalue> const &values = std::vector<colvarvalue>(0));
122 
126  virtual void calc_hills_force(size_t const &i,
127  hill_iter h_first,
128  hill_iter h_last,
129  std::vector<colvarvalue> &forces,
130  std::vector<colvarvalue> const &values = std::vector<colvarvalue>(0));
131 
132 
135 
138  bool use_grids;
139 
142 
145 
147  size_t grids_freq;
148 
152 
154  bool dump_fes;
155 
159 
163 
166 
169 
170  // EBmeta parameters
171  bool ebmeta;
172  colvar_grid_scalar* target_dist;
173  std::string target_dist_file;
174  cvm::real target_dist_volume;
175  size_t ebmeta_equil_steps;
176 
177 
182 
185 
188 
190  void project_hills(hill_iter h_first, hill_iter h_last,
192  bool print_progress = false);
193 
194 
195  // Multiple Replicas variables and functions
196 
198  std::string replica_id;
199 
201  std::string replica_file_name;
202 
204  virtual void update_replicas_registry();
205 
207  virtual void read_replica_files();
208 
210  virtual int write_replica_state_file();
211 
217  std::vector<colvarbias_meta *> replicas;
218 
221 
226  std::string replicas_registry;
228  std::string replica_list_file;
229 
232  std::string replica_state_file;
235 
239 
243  std::string replica_hills_file;
244 
246  std::ostream *replica_hills_os;
247 
250 
251 };
252 
253 
254 
255 
258 
259 protected:
260 
263 
266 
269 
271  std::vector<colvarvalue> centers;
272 
274  std::vector<cvm::real> widths;
275 
276 public:
277 
278  friend class colvarbias_meta;
279 
281  size_t it;
282 
284  std::string replica;
285 
291  inline hill(cvm::real const &W_in,
292  std::vector<colvar *> &cv,
293  cvm::real const &hill_width,
294  std::string const &replica_in = "")
295  : sW(1.0),
296  W(W_in),
297  centers(cv.size()),
298  widths(cv.size()),
299  it(cvm::step_absolute()),
300  replica(replica_in)
301  {
302  for (size_t i = 0; i < cv.size(); i++) {
303  centers[i].type(cv[i]->value());
304  centers[i] = cv[i]->value();
305  widths[i] = cv[i]->width * hill_width;
306  }
307  if (cvm::debug())
308  cvm::log("New hill, applied to "+cvm::to_str(cv.size())+
309  " collective variables, with centers "+
310  cvm::to_str(centers)+", widths "+
311  cvm::to_str(widths)+" and weight "+
312  cvm::to_str(W)+".\n");
313  }
314 
321  inline hill(size_t const &it_in,
322  cvm::real const &W_in,
323  std::vector<colvarvalue> const &centers_in,
324  std::vector<cvm::real> const &widths_in,
325  std::string const &replica_in = "")
326  : sW(1.0),
327  W(W_in),
328  centers(centers_in),
329  widths(widths_in),
330  it(it_in),
331  replica(replica_in)
332  {}
333 
335  inline hill(colvarbias_meta::hill const &h)
336  : sW(1.0),
337  W(h.W),
338  centers(h.centers),
339  widths(h.widths),
340  it(h.it),
341  replica(h.replica)
342  {}
343 
345  inline ~hill()
346  {}
347 
349  inline cvm::real energy()
350  {
351  return W * sW * hill_value;
352  }
353 
355  inline cvm::real energy(cvm::real const &new_weight)
356  {
357  return new_weight * sW * hill_value;
358  }
359 
361  inline cvm::real const &value()
362  {
363  return hill_value;
364  }
365 
367  inline void value(cvm::real const &new_value)
368  {
369  hill_value = new_value;
370  }
371 
373  inline cvm::real weight()
374  {
375  return W * sW;
376  }
377 
379  inline void scale(cvm::real const &new_scale_fac)
380  {
381  sW = new_scale_fac;
382  }
383 
385  inline std::vector<colvarvalue> & center()
386  {
387  return centers;
388  }
389 
391  inline colvarvalue & center(size_t const &i)
392  {
393  return centers[i];
394  }
395 
397  inline friend bool operator < (hill const &h1, hill const &h2)
398  {
399  if (h1.it < h2.it) return true;
400  else return false;
401  }
402 
404  inline friend bool operator <= (hill const &h1, hill const &h2)
405  {
406  if (h1.it <= h2.it) return true;
407  else return false;
408  }
409 
411  inline friend bool operator > (hill const &h1, hill const &h2)
412  {
413  if (h1.it > h2.it) return true;
414  else return false;
415  }
416 
418  inline friend bool operator >= (hill const &h1, hill const &h2)
419  {
420  if (h1.it >= h2.it) return true;
421  else return false;
422  }
423 
425  inline friend bool operator == (hill const &h1, hill const &h2)
426  {
427  if ( (h1.it >= h2.it) && (h1.replica == h2.replica) ) return true;
428  else return false;
429  }
430 
432  std::string output_traj();
433 
435  inline friend std::ostream & operator << (std::ostream &os,
436  hill const &h);
437 
438 };
439 
440 
441 #endif
A hill for the metadynamics bias.
Definition: colvarbias_meta.h:257
std::string replicas_registry
List of replicas (and their output list files)
Definition: colvarbias_meta.h:226
bool well_tempered
Whether to use well-tempered metadynamics.
Definition: colvarbias_meta.h:165
size_t replica_update_freq
Frequency at which data the "mirror" biases are updated.
Definition: colvarbias_meta.h:220
friend bool operator<=(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:404
Metadynamics bias (implementation of colvarbias)
Definition: colvarbias_meta.h:15
friend bool operator==(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:425
bool dump_fes_save
Dump the free energy surface files at different time steps, appending the step number to each file...
Definition: colvarbias_meta.h:162
cvm::real hill_width
width of a hill
Definition: colvarbias_meta.h:72
virtual void read_replica_files()
Read new data from replicas&#39; files.
Definition: colvarbias_meta.cpp:1063
std::vector< cvm::real > widths
Widths of the hill in the collective variable space.
Definition: colvarbias_meta.h:274
colvarvalue & center(size_t const &i)
Get the i-th component of the center.
Definition: colvarbias_meta.h:391
std::string replica_hills_file
Definition: colvarbias_meta.h:243
void scale(cvm::real const &new_scale_fac)
Scale the weight with this factor (by default 1.0 is used)
Definition: colvarbias_meta.h:379
virtual void calc_hills_force(size_t const &i, hill_iter h_first, hill_iter h_last, std::vector< colvarvalue > &forces, std::vector< colvarvalue > const &values=std::vector< colvarvalue >(0))
Calculate the forces acting on the i-th colvar, incrementing colvar_forces[i]; must be called after c...
Definition: colvarbias_meta.cpp:713
virtual int update()
Definition: colvarbias_meta.cpp:337
std::string replica_id
Identifier for this replica.
Definition: colvarbias_meta.h:198
hill(size_t const &it_in, cvm::real const &W_in, std::vector< colvarvalue > const &centers_in, std::vector< cvm::real > const &widths_in, std::string const &replica_in="")
General constructor: all data are explicitly passed as arguments (used for instance when reading hill...
Definition: colvarbias_meta.h:321
bool keep_hills
Whether to keep the hills in the restart file (e.g. to do meaningful accurate rebinning afterwards) ...
Definition: colvarbias_meta.h:151
bool use_grids
Bin the hills on grids of energy and forces, and use them to force the colvars (as opposed to derivin...
Definition: colvarbias_meta.h:138
colvar_grid_gradient * hills_energy_gradients
Hill forces, cached on a grid.
Definition: colvarbias_meta.h:187
Collective variables module (main class)
Definition: colvarmodule.h:62
bool safely_read_restart
Try to read the restart information by allocating new grids before replacing the current ones (used e...
Definition: colvarbias_meta.h:181
cvm::real energy()
Get the energy.
Definition: colvarbias_meta.h:349
virtual int setup_output()
(Re)initialize the output files (does not write them yet)
Definition: colvarbias_meta.cpp:1488
std::list< hill > hills
List of hills used on this bias (total); if a grid is employed, these don&#39;t need to be updated at eve...
Definition: colvarbias_meta.h:87
hill(colvarbias_meta::hill const &h)
Copy constructor.
Definition: colvarbias_meta.h:335
int replica_hills_file_pos
Position within replica_hills_file (when reading it)
Definition: colvarbias_meta.h:249
virtual int set_state_params(std::string const &state_conf)
Read the values of specific mutable properties from a string.
Definition: colvarbias_meta.cpp:1179
size_t grids_freq
How often the hills should be projected onto the grids.
Definition: colvarbias_meta.h:147
virtual int write_output_files()
Write any output files that this bias may have (e.g. PMF files)
Definition: colvarbias_meta.cpp:1667
cvm::real energy(cvm::real const &new_weight)
Get the energy using another hill weight.
Definition: colvarbias_meta.h:355
One replica (default)
Definition: colvarbias_meta.h:25
double real
Defining an abstract real number allows to switch precision.
Definition: colvarmodule.h:76
cvm::real hill_weight
Height of new hills.
Definition: colvarbias_meta.h:134
cvm::real bias_temperature
Biasing temperature in well-tempered metadynamics.
Definition: colvarbias_meta.h:168
Collective variable bias, base class.
Definition: colvarbias.h:12
virtual std::ostream & write_state_data(std::ostream &os)
Write all mutable data not already written by get_state_params()
Definition: colvarbias_meta.cpp:1615
virtual void update_replicas_registry()
Read the existing replicas on registry.
Definition: colvarbias_meta.cpp:929
virtual void calc_hills(hill_iter h_first, hill_iter h_last, cvm::real &energy, std::vector< colvarvalue > const &values=std::vector< colvarvalue >(0))
Calculate the values of the hills, incrementing bias_energy.
Definition: colvarbias_meta.cpp:669
std::string replica_list_file
List of files written by this replica.
Definition: colvarbias_meta.h:228
std::ostream * replica_hills_os
Output stream corresponding to replica_hills_file.
Definition: colvarbias_meta.h:246
bool b_hills_traj
Write the hill logfile.
Definition: colvarbias_meta.h:78
Value of a collective variable: this is a metatype which can be set at runtime. By default it is set ...
Definition: colvarvalue.h:34
friend bool operator<(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:397
cvm::real hill_value
Value of the hill function (ranges between 0 and 1)
Definition: colvarbias_meta.h:262
void recount_hills_off_grid(hill_iter h_first, hill_iter h_last, colvar_grid_scalar *ge)
Regenerate the hills_off_grid list.
Definition: colvarbias_meta.cpp:892
cvm::real W
Maximum height in energy of the hill.
Definition: colvarbias_meta.h:268
size_t it
Time step at which this hill was added.
Definition: colvarbias_meta.h:281
std::string replica_file_name
File containing the paths to the output files from this replica.
Definition: colvarbias_meta.h:201
Communication comm
Communication between different replicas.
Definition: colvarbias_meta.h:31
std::string replica_state_file
Definition: colvarbias_meta.h:232
std::istream & read_hill(std::istream &is)
Read a hill from a file.
Definition: colvarbias_meta.cpp:1409
Hills added concurrently by several replicas.
Definition: colvarbias_meta.h:27
bool dump_replica_fes
Dump the free energy surface (.pmf file) every restartFrequency using only the hills from this replic...
Definition: colvarbias_meta.h:158
hill(cvm::real const &W_in, std::vector< colvar *> &cv, cvm::real const &hill_width, std::string const &replica_in="")
Runtime constructor: data are read directly from collective variables.
Definition: colvarbias_meta.h:291
std::ostream * hills_traj_os
Logfile of hill management (creation and deletion)
Definition: colvarbias_meta.h:80
virtual std::list< hill >::const_iterator delete_hill(hill_iter &h)
Remove a previously saved hill (returns an iterator for the next hill in the list) ...
Definition: colvarbias_meta.cpp:304
friend bool operator>(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:411
cvm::real weight()
Get the weight.
Definition: colvarbias_meta.h:373
size_t new_hill_freq
Number of simulation steps between two hills.
Definition: colvarbias_meta.h:75
static std::string to_str(T const &x, size_t const &width=0, size_t const &prec=0)
Quick conversion of an object to a string.
Definition: colvarmodule.h:607
hill_iter new_hills_off_grid_begin
Same as new_hills_begin, but for the off-grid ones.
Definition: colvarbias_meta.h:98
colvar_grid_scalar * hills_energy
Hill energy, cached on a grid.
Definition: colvarbias_meta.h:184
bool expand_grids
Should the grids be expanded if necessary?
Definition: colvarbias_meta.h:144
Class for accumulating a scalar function on a grid.
Definition: colvargrid.h:1249
size_t update_status
Definition: colvarbias_meta.h:238
std::string replica
Identity of the replica who added this hill (only in multiple replica simulations) ...
Definition: colvarbias_meta.h:284
Communication
Communication between different replicas.
Definition: colvarbias_meta.h:23
virtual int clear_state_data()
Delete only the allocatable data (save memory)
Definition: colvarbias_meta.cpp:247
friend bool operator>=(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:418
std::vector< colvarvalue > & center()
Get the center of the hill.
Definition: colvarbias_meta.h:385
virtual int write_replica_state_file()
Write data to other replicas.
Definition: colvarbias_meta.cpp:1736
virtual int calc_energy(std::vector< colvarvalue > const &values=std::vector< colvarvalue >(0))
Compute the energy of the bias with alternative values of the collective variables (suitable for bias...
Definition: colvarbias_meta.cpp:557
virtual std::istream & read_state_data(std::istream &os)
Read all mutable data not already set by set_state_params()
Definition: colvarbias_meta.cpp:1195
std::vector< colvarvalue > centers
Center of the hill in the collective variable space.
Definition: colvarbias_meta.h:271
std::list< hill > hills_off_grid
List of hills used on this bias that are on the boundary edges; these are updated regardless of wheth...
Definition: colvarbias_meta.h:95
friend std::ostream & operator<<(std::ostream &os, hill const &h)
Write the hill to an output stream.
Definition: colvarbias_meta.cpp:1840
~hill()
Destructor.
Definition: colvarbias_meta.h:345
virtual std::string const get_state_params() const
Write the values of specific mutable properties to a string.
Definition: colvarbias_meta.cpp:1606
void init()
Set the object ready to parse a new configuration string.
Definition: colvarparse.h:61
static void log(std::string const &message)
Print a message to the main log.
Definition: colvarmodule.cpp:1518
virtual std::list< hill >::const_iterator create_hill(hill const &h)
Add a new hill; if a .hills trajectory is written, write it there; if there is more than one replica...
Definition: colvarbias_meta.cpp:271
std::string const hills_traj_file_name() const
Name of the hill logfile.
Definition: colvarbias_meta.cpp:1595
std::string replicas_registry_file
Definition: colvarbias_meta.h:224
cvm::real sW
Scale factor, which could be modified at runtime (default: 1)
Definition: colvarbias_meta.h:265
virtual int write_state_to_replicas()
If this bias is communicating with other replicas through files, send it to them. ...
Definition: colvarbias_meta.cpp:1653
void project_hills(hill_iter h_first, hill_iter h_last, colvar_grid_scalar *ge, colvar_grid_gradient *gf, bool print_progress=false)
Project the selected hills onto grids.
Definition: colvarbias_meta.cpp:788
cvm::real const & value()
Get the current hill value.
Definition: colvarbias_meta.h:361
Class for accumulating the gradient of a scalar function on a grid.
Definition: colvargrid.h:1378
std::vector< colvarbias_meta * > replicas
Additional, "mirror" metadynamics biases, to collect info from the other replicas.
Definition: colvarbias_meta.h:217
bool replica_state_file_in_sync
Whether a mirror bias has read the latest version of its state file.
Definition: colvarbias_meta.h:234
bool dump_fes
Dump the free energy surface (.pmf file) every restartFrequency.
Definition: colvarbias_meta.h:154
void value(cvm::real const &new_value)
Set the hill value as specified.
Definition: colvarbias_meta.h:367
std::string output_traj()
Represent the hill ina string suitable for a trajectory file.
Definition: colvarbias_meta.cpp:1809
bool rebin_grids
Rebin the hills upon restarting.
Definition: colvarbias_meta.h:141
hill_iter new_hills_begin
Iterator to the first of the "newest" hills (when using grids, those who haven&#39;t been mapped yet) ...
Definition: colvarbias_meta.h:91
static bool debug()
Whether debug output should be enabled (compile-time option)
Definition: colvarmodule.h:237
Base class for unconstrained thermodynamic-integration FE estimator.
Definition: colvarbias.h:222