Collective Variables Module - Developer Documentation
Loading...
Searching...
No Matches
colvarbias_meta.h
1// -*- c++ -*-
2
3// This file is part of the Collective Variables module (Colvars).
4// The original version of Colvars and its updates are located at:
5// https://github.com/Colvars/colvars
6// Please update all Colvars source files before making any changes.
7// If you wish to distribute your changes, please submit them to the
8// Colvars repository at GitHub.
9
10#ifndef COLVARBIAS_META_H
11#define COLVARBIAS_META_H
12
13#include <vector>
14#include <list>
15#include <iosfwd>
16
17#include "colvarbias.h"
18#include "colvargrid.h"
19
20
23 : public virtual colvarbias,
24 public virtual colvarbias_ti
25{
26
27public:
28
35 };
36
39
40 colvarbias_meta(char const *key);
41 virtual ~colvarbias_meta();
42
43 virtual int init(std::string const &conf);
44 virtual int init_replicas_params(std::string const &conf);
45 virtual int init_well_tempered_params(std::string const &conf);
46 virtual int init_ebmeta_params(std::string const &conf);
47
48 virtual int clear_state_data();
49
50 virtual int update();
51 virtual int update_grid_params();
52 virtual int update_bias();
53 virtual int update_grid_data();
54 virtual int replica_share();
55 virtual size_t replica_share_freq() const;
56
57 virtual int calc_energy(std::vector<colvarvalue> const *values);
58 virtual int calc_forces(std::vector<colvarvalue> const *values);
59
60 virtual std::string const get_state_params() const;
61 virtual int set_state_params(std::string const &state_conf);
62
63 virtual std::ostream &write_state_data(std::ostream &os);
65 virtual std::istream &read_state_data(std::istream &is);
67
68private:
69
70 template <typename IST, typename GT>
71 IST &read_grid_data_template_(IST &is, std::string const &key, GT *grid, GT *backup_grid);
72
73 template <typename IST> IST &read_state_data_template_(IST &is);
74
75 template <typename OST> OST &write_state_data_template_(OST &os);
76
77public:
78
81
82 virtual int setup_output();
83 virtual int write_output_files();
84 virtual void write_pmf();
85 virtual int write_state_to_replicas();
86
87 class hill;
88 typedef std::list<hill>::iterator hill_iter;
89
90protected:
91
97
99 std::vector<cvm::real> colvar_sigmas;
100
103
106
108 std::string const hills_traj_file_name() const;
109
112 std::list<hill> hills;
113
117
120 std::list<hill> hills_off_grid;
121
124
126 void recount_hills_off_grid(hill_iter h_first, hill_iter h_last,
128
129 template <typename OST> OST &write_hill_template_(OST &os, colvarbias_meta::hill const &h);
130
132 std::ostream &write_hill(std::ostream &os, hill const &h);
133
136
137 template <typename IST> IST &read_hill_template_(IST &is);
138
140 std::istream & read_hill(std::istream &is);
141
144
148 std::list<hill>::const_iterator add_hill(hill const &h);
149
152 std::list<hill>::const_iterator delete_hill(hill_iter &h);
153
156 virtual void calc_hills(hill_iter h_first,
157 hill_iter h_last,
158 cvm::real &energy,
159 std::vector<colvarvalue> const *values);
160
164 virtual void calc_hills_force(size_t const &i,
165 hill_iter h_first,
166 hill_iter h_last,
167 std::vector<colvarvalue> &forces,
168 std::vector<colvarvalue> const *values);
169
170
173
177
180
183
186
189
192
195
199
203
206
209
211 bool ebmeta;
212
215
218
219
224
227
230
232 void project_hills(hill_iter h_first, hill_iter h_last,
234 bool print_progress = false);
235
236
237 // Multiple Replicas variables and functions
238
240 std::string replica_id;
241
243 std::string replica_file_name;
244
246 virtual int update_replicas_registry();
247
249 virtual int read_replica_files();
250
252 virtual int write_replica_state_file();
253
255 virtual int reopen_replica_buffer_file();
256
262 std::vector<colvarbias_meta *> replicas;
263
266
271 std::string replicas_registry;
273 std::string replica_list_file;
274
280
284
289
292
294 std::ostringstream hills_traj_os_buf;
295};
296
297
298
299
302
303protected:
304
307
310
313
316
318 std::vector<colvarvalue> centers;
319
321 std::vector<cvm::real> sigmas;
322
324 std::string replica;
325
326public:
327
328 friend class colvarbias_meta;
329
337 std::vector<colvarvalue> const &cv_values,
338 std::vector<cvm::real> const &cv_sigmas,
339 std::string const &replica = "");
340
342 hill(colvarbias_meta::hill const &h);
343
345 ~hill();
346
349
352 {
353 return W * sW * hill_value;
354 }
355
357 inline cvm::real energy(cvm::real const &new_weight)
358 {
359 return new_weight * sW * hill_value;
360 }
361
363 inline cvm::real const &value()
364 {
365 return hill_value;
366 }
367
369 inline void value(cvm::real const &new_value)
370 {
371 hill_value = new_value;
372 }
373
376 {
377 return W * sW;
378 }
379
381 inline void scale(cvm::real const &new_scale_fac)
382 {
383 sW = new_scale_fac;
384 }
385
387 inline std::vector<colvarvalue> & center()
388 {
389 return centers;
390 }
391
393 inline colvarvalue & center(size_t const &i)
394 {
395 return centers[i];
396 }
397
399 inline friend bool operator < (hill const &h1, hill const &h2)
400 {
401 if (h1.it < h2.it) return true;
402 else return false;
403 }
404
406 inline friend bool operator <= (hill const &h1, hill const &h2)
407 {
408 if (h1.it <= h2.it) return true;
409 else return false;
410 }
411
413 inline friend bool operator > (hill const &h1, hill const &h2)
414 {
415 if (h1.it > h2.it) return true;
416 else return false;
417 }
418
420 inline friend bool operator >= (hill const &h1, hill const &h2)
421 {
422 if (h1.it >= h2.it) return true;
423 else return false;
424 }
425
427 inline friend bool operator == (hill const &h1, hill const &h2)
428 {
429 if ( (h1.it >= h2.it) && (h1.replica == h2.replica) ) return true;
430 else return false;
431 }
432
434 std::string output_traj();
435
436};
437
438
439#endif
Class for accumulating the gradient of a scalar function on a grid.
Definition: colvargrid.h:1555
Class for accumulating a scalar function on a grid.
Definition: colvargrid.h:1242
A hill for the metadynamics bias.
Definition: colvarbias_meta.h:301
colvarvalue & center(size_t const &i)
Get the i-th component of the center.
Definition: colvarbias_meta.h:393
friend bool operator>=(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:420
cvm::real energy(cvm::real const &new_weight)
Get the energy using another hill weight.
Definition: colvarbias_meta.h:357
cvm::real hill_value
Value of the hill function (ranges between 0 and 1)
Definition: colvarbias_meta.h:309
~hill()
Destructor.
Definition: colvarbias_meta.cpp:2125
cvm::step_number it
Time step at which this hill was added.
Definition: colvarbias_meta.h:306
void value(cvm::real const &new_value)
Set the hill value as specified.
Definition: colvarbias_meta.h:369
hill & operator=(colvarbias_meta::hill const &h)
Assignment operator.
Definition: colvarbias_meta.cpp:2111
cvm::real const & value()
Get the current hill value.
Definition: colvarbias_meta.h:363
friend bool operator<=(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:406
friend bool operator<(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:399
std::string replica
Identity of the replica who added this hill.
Definition: colvarbias_meta.h:324
std::vector< colvarvalue > centers
Centers of the hill in the collective variable space.
Definition: colvarbias_meta.h:318
cvm::real W
Maximum height in energy of the hill.
Definition: colvarbias_meta.h:315
cvm::real sW
Scale factor, which could be modified at runtime (default: 1)
Definition: colvarbias_meta.h:312
friend bool operator==(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:427
void scale(cvm::real const &new_scale_fac)
Scale the weight with this factor (by default 1.0 is used)
Definition: colvarbias_meta.h:381
cvm::real weight()
Get the weight.
Definition: colvarbias_meta.h:375
std::vector< colvarvalue > & center()
Get the center of the hill.
Definition: colvarbias_meta.h:387
std::string output_traj()
Represent the hill ina string suitable for a trajectory file.
Definition: colvarbias_meta.cpp:2038
friend bool operator>(hill const &h1, hill const &h2)
Comparison operator.
Definition: colvarbias_meta.h:413
std::vector< cvm::real > sigmas
Half-widths of the hill in the collective variable space.
Definition: colvarbias_meta.h:321
cvm::real energy()
Get the energy.
Definition: colvarbias_meta.h:351
Metadynamics bias (implementation of colvarbias)
Definition: colvarbias_meta.h:25
bool dump_replica_fes
Dump the free energy surface (.pmf file) every restartFrequency using only the hills from this replic...
Definition: colvarbias_meta.h:198
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:223
std::string replica_id
Identifier for this replica.
Definition: colvarbias_meta.h:240
colvar_grid_gradient * hills_energy_gradients
Hill forces, cached on a grid.
Definition: colvarbias_meta.h:229
virtual int reopen_replica_buffer_file()
Call this after write_replica_state_file()
Definition: colvarbias_meta.cpp:2020
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:120
size_t grids_freq
How often the hills should be projected onto the grids.
Definition: colvarbias_meta.h:185
virtual int update()
Definition: colvarbias_meta.cpp:422
cvm::real hill_weight
Height of new hills.
Definition: colvarbias_meta.h:172
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:176
size_t new_hill_freq
Number of simulation steps between two hills.
Definition: colvarbias_meta.h:102
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)
Calculate the forces acting on the i-th colvar, incrementing colvar_forces[i]; must be called after c...
Definition: colvarbias_meta.cpp:811
colvar_grid_scalar * target_dist
Target distribution for EBmeta.
Definition: colvarbias_meta.h:214
std::string replicas_registry_file
Definition: colvarbias_meta.h:269
virtual int calc_forces(std::vector< colvarvalue > const *values)
Definition: colvarbias_meta.cpp:711
bool ebmeta
Ensemble-biased metadynamics (EBmeta) flag.
Definition: colvarbias_meta.h:211
virtual int write_output_files()
Write any output files that this bias may have (e.g. PMF files)
Definition: colvarbias_meta.cpp:1891
std::string replica_list_file
List of files written by this replica.
Definition: colvarbias_meta.h:273
std::ostream & write_hill(std::ostream &os, hill const &h)
Write a hill to a formatted stream.
Definition: colvarbias_meta.cpp:1562
virtual std::ostream & write_state_data(std::ostream &os)
Write all mutable data not already written by get_state_params() to a formatted stream.
Definition: colvarbias_meta.cpp:1864
cvm::step_number ebmeta_equil_steps
Number of equilibration steps for EBmeta.
Definition: colvarbias_meta.h:217
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:886
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:964
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:112
virtual int init(std::string const &conf)
Parse config string and (re)initialize.
Definition: colvarbias_meta.cpp:51
bool keep_hills
Keep hills in the restart file (e.g. to accurately rebin later)
Definition: colvarbias_meta.h:188
std::string replica_hills_file
Definition: colvarbias_meta.h:288
virtual std::string const get_state_params() const
Write the values of specific mutable properties to a string.
Definition: colvarbias_meta.cpp:1817
virtual int clear_state_data()
Delete only the allocatable data (save memory)
Definition: colvarbias_meta.cpp:333
virtual std::istream & read_state_data(std::istream &is)
Read all mutable data not already set by set_state_params() from a formatted stream.
Definition: colvarbias_meta.cpp:1437
cvm::real bias_temperature
Biasing temperature in well-tempered metadynamics.
Definition: colvarbias_meta.h:208
bool rebin_grids
Rebin the hills upon restarting.
Definition: colvarbias_meta.h:179
virtual int read_replica_files()
Read new data from replicas' files.
Definition: colvarbias_meta.cpp:1141
std::list< hill >::const_iterator add_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:357
std::string replica_state_file
Definition: colvarbias_meta.h:277
std::istream & read_hill(std::istream &is)
Read a new hill from a formatted stream.
Definition: colvarbias_meta.cpp:1698
bool expand_grids
Should the grids be expanded if necessary?
Definition: colvarbias_meta.h:182
bool well_tempered
Whether to use well-tempered metadynamics.
Definition: colvarbias_meta.h:205
bool restart_keep_hills
value of keepHills saved in the most recent restart file
Definition: colvarbias_meta.h:191
std::ostringstream hills_traj_os_buf
Cache of the hills trajectory.
Definition: colvarbias_meta.h:294
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:202
virtual int setup_output()
(Re)initialize the output files (does not write them yet)
Definition: colvarbias_meta.cpp:1710
size_t update_status
Definition: colvarbias_meta.h:283
Communication
Communication between different replicas.
Definition: colvarbias_meta.h:30
@ multiple_replicas
Hills added concurrently by several replicas.
Definition: colvarbias_meta.h:34
@ single_replica
One replica (default)
Definition: colvarbias_meta.h:32
colvar_grid_scalar * hills_energy
Hill energy, cached on a grid.
Definition: colvarbias_meta.h:226
size_t replica_update_freq
Frequency at which data the "mirror" biases are updated.
Definition: colvarbias_meta.h:265
Communication comm
Communication between different replicas.
Definition: colvarbias_meta.h:38
virtual int write_replica_state_file()
Write full state information to be read by other replicas.
Definition: colvarbias_meta.cpp:1992
bool b_hills_traj
Write the hill logfile.
Definition: colvarbias_meta.h:105
virtual size_t replica_share_freq() const
Report the frequency at which this bias needs to communicate with replicas.
Definition: colvarbias_meta.cpp:1001
virtual void calc_hills(hill_iter h_first, hill_iter h_last, cvm::real &energy, std::vector< colvarvalue > const *values)
Calculate the values of the hills, incrementing bias_energy.
Definition: colvarbias_meta.cpp:781
cvm::real hill_width
Definition: colvarbias_meta.h:96
bool replica_state_file_in_sync
Whether a mirror bias has read the latest version of its state file.
Definition: colvarbias_meta.h:279
void rebin_grids_after_restart()
Function called by read_state_data() to execute rebinning (if requested)
Definition: colvarbias_meta.cpp:1449
std::streampos replica_hills_file_pos
Position within replica_hills_file (when reading it)
Definition: colvarbias_meta.h:291
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:390
virtual int write_state_to_replicas()
If this bias is communicating with other replicas through files, send it to them.
Definition: colvarbias_meta.cpp:1876
virtual int calc_energy(std::vector< colvarvalue > const *values)
Definition: colvarbias_meta.cpp:650
virtual int set_state_params(std::string const &state_conf)
Read the values of specific mutable properties from a string.
Definition: colvarbias_meta.cpp:1270
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:116
std::string replica_file_name
File containing the paths to the output files from this replica.
Definition: colvarbias_meta.h:243
std::string const hills_traj_file_name() const
Name of the hill logfile.
Definition: colvarbias_meta.cpp:1806
virtual int update_replicas_registry()
Read the existing replicas on registry.
Definition: colvarbias_meta.cpp:1007
std::string replicas_registry
List of replicas (and their output list files)
Definition: colvarbias_meta.h:271
bool dump_fes
Dump the free energy surface (.pmf file) every restartFrequency.
Definition: colvarbias_meta.h:194
hill_iter new_hills_off_grid_begin
Same as new_hills_begin, but for the off-grid ones.
Definition: colvarbias_meta.h:123
std::vector< colvarbias_meta * > replicas
Additional, "mirror" metadynamics biases, to collect info from the other replicas.
Definition: colvarbias_meta.h:262
std::vector< cvm::real > colvar_sigmas
The sigma parameters of the Gaussian hills.
Definition: colvarbias_meta.h:99
Base class for unconstrained thermodynamic-integration FE estimator.
Definition: colvarbias.h:333
Collective variable bias, base class.
Definition: colvarbias.h:23
double real
Defining an abstract real number allows to switch precision.
Definition: colvarmodule.h:95
long long step_number
Use a 64-bit integer to store the step number.
Definition: colvarmodule.h:92
Value of a collective variable: this is a metatype which can be set at runtime. By default it is set ...
Definition: colvarvalue.h:43
Definition: colvars_memstream.h:30