Collective Variables Module - Developer Documentation
 All Classes Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
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 
15 class colvarbias_meta : public colvarbias {
16 
17 public:
18 
25  };
26 
29 
30  colvarbias_meta(char const *key);
31  virtual int init(std::string const &conf);
32  virtual int init_well_tempered_params(std::string const &conf);
33  virtual int init_ebmeta_params(std::string const &conf);
34  virtual ~colvarbias_meta();
35 
36  virtual int update();
37  virtual int update_grid_params();
38  virtual int update_bias();
39  virtual int update_grid_data();
40  virtual int replica_share();
41 
42  virtual int calc_energy(std::vector<colvarvalue> const &values =
43  std::vector<colvarvalue>(0));
44  virtual int calc_forces(std::vector<colvarvalue> const &values =
45  std::vector<colvarvalue>(0));
46 
47  virtual std::string const get_state_params() const;
48  virtual int set_state_params(std::string const &state_conf);
49  virtual std::ostream & write_state_data(std::ostream &os);
50  virtual std::istream & read_state_data(std::istream &os);
51 
52  virtual int setup_output();
53  virtual int write_output_files();
54  virtual void write_pmf();
55  virtual int write_state_to_replicas();
56 
57  class hill;
58  typedef std::list<hill>::iterator hill_iter;
59 
60 protected:
61 
67 
69  size_t new_hill_freq;
70 
74  cvm::ofstream hills_traj_os;
75 
78  std::list<hill> hills;
79 
82  hill_iter new_hills_begin;
83 
86  std::list<hill> hills_off_grid;
87 
90 
92  void recount_hills_off_grid(hill_iter h_first, hill_iter h_last,
93  colvar_grid_scalar *ge);
94 
96  std::istream & read_hill(std::istream &is);
97 
101  virtual std::list<hill>::const_iterator create_hill(hill const &h);
102 
105  virtual std::list<hill>::const_iterator delete_hill(hill_iter &h);
106 
109  virtual void calc_hills(hill_iter h_first,
110  hill_iter h_last,
111  cvm::real &energy,
112  std::vector<colvarvalue> const &values = std::vector<colvarvalue>(0));
113 
117  virtual void calc_hills_force(size_t const &i,
118  hill_iter h_first,
119  hill_iter h_last,
120  std::vector<colvarvalue> &forces,
121  std::vector<colvarvalue> const &values = std::vector<colvarvalue>(0));
122 
123 
126 
129  bool use_grids;
130 
133 
136 
138  size_t grids_freq;
139 
143 
145  bool dump_fes;
146 
150 
154 
157 
160 
161  // EBmeta parameters
162  bool ebmeta;
163  colvar_grid_scalar* target_dist;
164  std::string target_dist_file;
165  cvm::real target_dist_volume;
166  size_t ebmeta_equil_steps;
167 
168 
173 
176 
179 
181  void project_hills(hill_iter h_first, hill_iter h_last,
183  bool print_progress = false);
184 
185 
186  // Multiple Replicas variables and functions
187 
189  std::string replica_id;
190 
192  std::string replica_file_name;
193 
195  virtual void update_replicas_registry();
196 
198  virtual void read_replica_files();
199 
201  virtual int write_replica_state_file();
202 
208  std::vector<colvarbias_meta *> replicas;
209 
212 
217  std::string replicas_registry;
219  std::string replica_list_file;
220 
223  std::string replica_state_file;
226 
230 
234  std::string replica_hills_file;
235 
237  cvm::ofstream replica_hills_os;
238 
241 
242 };
243 
244 
245 
246 
249 
250 protected:
251 
254 
257 
260 
262  std::vector<colvarvalue> centers;
263 
265  std::vector<cvm::real> widths;
266 
267 public:
268 
269  friend class colvarbias_meta;
270 
272  size_t it;
273 
275  std::string replica;
276 
282  inline hill(cvm::real const &W_in,
283  std::vector<colvar *> &cv,
284  cvm::real const &hill_width,
285  std::string const &replica_in = "")
286  : sW(1.0),
287  W(W_in),
288  centers(cv.size()),
289  widths(cv.size()),
290  it(cvm::step_absolute()),
291  replica(replica_in)
292  {
293  for (size_t i = 0; i < cv.size(); i++) {
294  centers[i].type(cv[i]->value());
295  centers[i] = cv[i]->value();
296  widths[i] = cv[i]->width * hill_width;
297  }
298  if (cvm::debug())
299  cvm::log("New hill, applied to "+cvm::to_str(cv.size())+
300  " collective variables, with centers "+
301  cvm::to_str(centers)+", widths "+
302  cvm::to_str(widths)+" and weight "+
303  cvm::to_str(W)+".\n");
304  }
305 
312  inline hill(size_t const &it_in,
313  cvm::real const &W_in,
314  std::vector<colvarvalue> const &centers_in,
315  std::vector<cvm::real> const &widths_in,
316  std::string const &replica_in = "")
317  : sW(1.0),
318  W(W_in),
319  centers(centers_in),
320  widths(widths_in),
321  it(it_in),
322  replica(replica_in)
323  {}
324 
326  inline hill(colvarbias_meta::hill const &h)
327  : sW(1.0),
328  W(h.W),
329  centers(h.centers),
330  widths(h.widths),
331  it(h.it),
332  replica(h.replica)
333  {}
334 
336  inline ~hill()
337  {}
338 
340  inline cvm::real energy()
341  {
342  return W * sW * hill_value;
343  }
344 
346  inline cvm::real energy(cvm::real const &new_weight)
347  {
348  return new_weight * sW * hill_value;
349  }
350 
352  inline cvm::real const &value()
353  {
354  return hill_value;
355  }
356 
358  inline void value(cvm::real const &new_value)
359  {
360  hill_value = new_value;
361  }
362 
364  inline cvm::real weight()
365  {
366  return W * sW;
367  }
368 
370  inline void scale(cvm::real const &new_scale_fac)
371  {
372  sW = new_scale_fac;
373  }
374 
376  inline std::vector<colvarvalue> & center()
377  {
378  return centers;
379  }
380 
382  inline colvarvalue & center(size_t const &i)
383  {
384  return centers[i];
385  }
386 
388  inline friend bool operator < (hill const &h1, hill const &h2)
389  {
390  if (h1.it < h2.it) return true;
391  else return false;
392  }
393 
395  inline friend bool operator <= (hill const &h1, hill const &h2)
396  {
397  if (h1.it <= h2.it) return true;
398  else return false;
399  }
400 
402  inline friend bool operator > (hill const &h1, hill const &h2)
403  {
404  if (h1.it > h2.it) return true;
405  else return false;
406  }
407 
409  inline friend bool operator >= (hill const &h1, hill const &h2)
410  {
411  if (h1.it >= h2.it) return true;
412  else return false;
413  }
414 
416  inline friend bool operator == (hill const &h1, hill const &h2)
417  {
418  if ( (h1.it >= h2.it) && (h1.replica == h2.replica) ) return true;
419  else return false;
420  }
421 
423  std::string output_traj();
424 
426  inline friend std::ostream & operator << (std::ostream &os,
427  hill const &h);
428 
429 };
430 
431 
432 #endif
A hill for the metadynamics bias.
Definition: colvarbias_meta.h:248
std::string replicas_registry
List of replicas (and their output list files)
Definition: colvarbias_meta.h:217
bool well_tempered
Whether to use well-tempered metadynamics.
Definition: colvarbias_meta.h:156
size_t replica_update_freq
Frequency at which data the "mirror" biases are updated.
Definition: colvarbias_meta.h:211
friend bool operator<=(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:395
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:416
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:153
cvm::real hill_width
width of a hill
Definition: colvarbias_meta.h:66
virtual void read_replica_files()
Read new data from replicas' files.
Definition: colvarbias_meta.cpp:1038
std::vector< cvm::real > widths
Widths of the hill in the collective variable space.
Definition: colvarbias_meta.h:265
colvarvalue & center(size_t const &i)
Get the i-th component of the center.
Definition: colvarbias_meta.h:382
std::string replica_hills_file
Definition: colvarbias_meta.h:234
void scale(cvm::real const &new_scale_fac)
Scale the weight with this factor (by default 1.0 is used)
Definition: colvarbias_meta.h:370
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:693
virtual int update()
Definition: colvarbias_meta.cpp:320
std::string replica_id
Identifier for this replica.
Definition: colvarbias_meta.h:189
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:312
bool keep_hills
Whether to keep the hills in the restart file (e.g. to do meaningful accurate rebinning afterwards) ...
Definition: colvarbias_meta.h:142
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:129
virtual std::string const get_state_params() const
Write the values of specific mutable properties to a string.
Definition: colvarbias_meta.cpp:1575
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:282
colvar_grid_gradient * hills_energy_gradients
Hill forces, cached on a grid.
Definition: colvarbias_meta.h:178
Collective variables module (main class)
Definition: colvarmodule.h:71
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:172
cvm::real energy()
Get the energy.
Definition: colvarbias_meta.h:340
virtual int setup_output()
(Re)initialize the output files (does not write them yet)
Definition: colvarbias_meta.cpp:1461
std::list< hill > hills
List of hills used on this bias (total); if a grid is employed, these don't need to be updated at eve...
Definition: colvarbias_meta.h:78
hill(colvarbias_meta::hill const &h)
Copy constructor.
Definition: colvarbias_meta.h:326
int replica_hills_file_pos
Position within replica_hills_file (when reading it)
Definition: colvarbias_meta.h:240
virtual int set_state_params(std::string const &state_conf)
Read the values of specific mutable properties from a string.
Definition: colvarbias_meta.cpp:1154
size_t grids_freq
How often the hills should be projected onto the grids.
Definition: colvarbias_meta.h:138
virtual int write_output_files()
Write any output files that this bias may have (e.g. PMF files)
Definition: colvarbias_meta.cpp:1635
cvm::real energy(cvm::real const &new_weight)
Get the energy using another hill weight.
Definition: colvarbias_meta.h:346
One replica (default)
Definition: colvarbias_meta.h:22
double real
Defining an abstract real number allows to switch precision.
Definition: colvarmodule.h:85
cvm::real hill_weight
Height of new hills.
Definition: colvarbias_meta.h:125
cvm::real bias_temperature
Biasing temperature in well-tempered metadynamics.
Definition: colvarbias_meta.h:159
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:1584
virtual void update_replicas_registry()
Read the existing replicas on registry.
Definition: colvarbias_meta.cpp:908
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:649
std::string replica_list_file
List of files written by this replica.
Definition: colvarbias_meta.h:219
bool b_hills_traj
Write the hill logfile.
Definition: colvarbias_meta.h:72
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:388
cvm::real hill_value
Value of the hill function (ranges between 0 and 1)
Definition: colvarbias_meta.h:253
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:872
cvm::real W
Maximum height in energy of the hill.
Definition: colvarbias_meta.h:259
size_t it
Time step at which this hill was added.
Definition: colvarbias_meta.h:272
std::string replica_file_name
File containing the paths to the output files from this replica.
Definition: colvarbias_meta.h:192
Communication comm
Communication between different replicas.
Definition: colvarbias_meta.h:28
std::string replica_state_file
Definition: colvarbias_meta.h:223
std::istream & read_hill(std::istream &is)
Read a hill from a file.
Definition: colvarbias_meta.cpp:1382
Hills added concurrently by several replicas.
Definition: colvarbias_meta.h:24
bool dump_replica_fes
Dump the free energy surface (.pmf file) every restartFrequency using only the hills from this replic...
Definition: colvarbias_meta.h:149
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:287
friend bool operator>(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:402
cvm::ofstream hills_traj_os
Logfile of hill management (creation and deletion)
Definition: colvarbias_meta.h:74
cvm::real weight()
Get the weight.
Definition: colvarbias_meta.h:364
cvm::ofstream replica_hills_os
Output stream corresponding to replica_hills_file.
Definition: colvarbias_meta.h:237
size_t new_hill_freq
Number of simulation steps between two hills.
Definition: colvarbias_meta.h:69
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:615
hill_iter new_hills_off_grid_begin
Same as new_hills_begin, but for the off-grid ones.
Definition: colvarbias_meta.h:89
colvar_grid_scalar * hills_energy
Hill energy, cached on a grid.
Definition: colvarbias_meta.h:175
bool expand_grids
Should the grids be expanded if necessary?
Definition: colvarbias_meta.h:135
Class for accumulating a scalar function on a grid.
Definition: colvargrid.h:1253
size_t update_status
Definition: colvarbias_meta.h:229
std::string replica
Identity of the replica who added this hill (only in multiple replica simulations) ...
Definition: colvarbias_meta.h:275
Communication
Communication between different replicas.
Definition: colvarbias_meta.h:20
friend bool operator>=(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:409
std::vector< colvarvalue > & center()
Get the center of the hill.
Definition: colvarbias_meta.h:376
virtual int write_replica_state_file()
Write data to other replicas.
Definition: colvarbias_meta.cpp:1702
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:537
virtual std::istream & read_state_data(std::istream &os)
Read all mutable data not already set by set_state_params()
Definition: colvarbias_meta.cpp:1170
std::vector< colvarvalue > centers
Center of the hill in the collective variable space.
Definition: colvarbias_meta.h:262
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:86
friend std::ostream & operator<<(std::ostream &os, hill const &h)
Write the hill to an output stream.
Definition: colvarbias_meta.cpp:1808
~hill()
Destructor.
Definition: colvarbias_meta.h:336
void init()
Set the object ready to parse a new configuration string.
Definition: colvarparse.h:71
static void log(std::string const &message)
Print a message to the main log.
Definition: colvarmodule.cpp:1467
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:254
std::string replicas_registry_file
Definition: colvarbias_meta.h:215
cvm::real sW
Scale factor, which could be modified at runtime (default: 1)
Definition: colvarbias_meta.h:256
virtual int write_state_to_replicas()
If this bias is communicating with other replicas through files, send it to them. ...
Definition: colvarbias_meta.cpp:1621
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:768
cvm::real const & value()
Get the current hill value.
Definition: colvarbias_meta.h:352
Class for accumulating the gradient of a scalar function on a grid.
Definition: colvargrid.h:1382
std::vector< colvarbias_meta * > replicas
Additional, "mirror" metadynamics biases, to collect info from the other replicas.
Definition: colvarbias_meta.h:208
bool replica_state_file_in_sync
Whether a mirror bias has read the latest version of its state file.
Definition: colvarbias_meta.h:225
bool dump_fes
Dump the free energy surface (.pmf file) every restartFrequency.
Definition: colvarbias_meta.h:145
void value(cvm::real const &new_value)
Set the hill value as specified.
Definition: colvarbias_meta.h:358
std::string output_traj()
Represent the hill ina string suitable for a trajectory file.
Definition: colvarbias_meta.cpp:1777
bool rebin_grids
Rebin the hills upon restarting.
Definition: colvarbias_meta.h:132
hill_iter new_hills_begin
Iterator to the first of the "newest" hills (when using grids, those who haven't been mapped yet) ...
Definition: colvarbias_meta.h:82
static bool debug()
Whether debug output should be enabled (compile-time option)
Definition: colvarmodule.h:225