Main Page | Namespace List | Class Hierarchy | Alphabetical List | Class List | File List | Class Members | File Members

colvarbias_meta Class Reference

Metadynamics bias (implementation of ). More...

#include <colvarbias_meta.h>

Inheritance diagram for colvarbias_meta:

colvarbias colvarparse List of all members.

Public Types

typedef std::list< hill
>::iterator 
hill_iter
enum  Communication { single_replica, multiple_replicas }
 Communication between different replicas. More...

Public Member Functions

 colvarbias_meta (std::string const &conf, char const *key)
 Constructor.
virtual ~colvarbias_meta ()
 Destructor.
virtual void update ()
 Retrieve colvar values and calculate their biasing forces.
virtual void analyse ()
 Perform analysis.
virtual std::istream & read_restart (std::istream &is)
 Read the bias configuration from a restart file.
virtual std::ostream & write_restart (std::ostream &os)
 Write the bias configuration to a restart file.

Public Attributes

Communication comm
 Communication between different replicas.

Protected Member Functions

void recount_hills_off_grid (hill_iter h_first, hill_iter h_last, colvar_grid_scalar *ge)
 Regenerate the hills_off_grid list.
std::istream & read_hill (std::istream &is)
 Read a hill from a file.
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, communicate it to the others.
virtual std::list< hill
>::const_iterator 
delete_hill (std::list< hill >::iterator &h)
 Remove a previously saved hill (returns an iterator for the next hill in the list).
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 colvar_energy.
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 calc_hills each time the values of the colvars are changed.
void project_hills (hill_iter h_first, hill_iter h_last, colvar_grid_scalar *ge, colvar_grid_gradient *gf, cvm::real const scale_factor=1.0)
 Project the selected hills onto grids.
virtual void register_replica_file (std::string const &new_file)
 Add this replica to the registry (called only when != ).
virtual void update_replica_files_registry ()
 Read the names all replica output files which have been added, skipping the one which is currently being written to (called only when != ).
virtual void read_replica_files ()
 Read new data from replicas' files (called only when != ).

Protected Attributes

std::list< hillhills
 Parse analysis tasks and options List of hills used on this bias (total); if a grid is employed, these don't need to be updated at every time step.
hill_iter new_hills_begin
 Iterator to the first of the "newest" hills (when using grids, those who haven't been mapped yet).
std::list< hillhills_off_grid
 List of hills used on this bias that are on the boundary edges; these are updated regardless of whether hills are used.
hill_iter new_hills_off_grid_begin
 Same as new_hills_begin, but for the off-grid ones.
cvm::real hill_weight
 Height of new hills.
bool use_grids
 Bin the hills on grids of energy and forces, and use them to force the colvars (as opposed to deriving the hills analytically).
bool rebin_grids
 Rebin the hills upon restarting.
bool expand_grids
 Should the grids be expanded if necessary?
size_t grids_freq
 How often the hills should be projected onto the grids.
bool keep_hills
 Whether to keep the hills in the restart file (e.g. to do meaningful accurate rebinning afterwards).
bool dump_fes
 Dump the free energy surface (.pmf file) every restartFrequency.
bool dump_fes_save
 Keep the free energy surface files at different iterations, appending a step number to each.
colvar_grid_scalarhills_energy
 Hill energy, cached on a grid.
colvar_grid_gradienthills_energy_gradients
 Hill forces, cached on a grid.
cvm::real hill_width
 width of a hill
size_t new_hill_freq
 Number of simulation steps between two hills.
bool b_hills_traj
 Write the hill logfile.
std::ofstream hills_traj_os
 Logfile of hill management (creation and deletion).
std::string replica
 Identifier for this replica.
size_t replica_update_freq
 Frequency at which output files from other replicas are checked.
std::string replica_files_registry
 List of hill files from all the replicas.
std::list< std::string > replica_files
 Replicas file names.
std::list< size_t > replica_files_pos
 Positions in replica files (files are reopened at these positions).
std::ofstream replica_out_file
 File to contain the hills created by this replica in this run.
std::string replica_out_file_name
 File to contain the hills created by this replica in this run (name).
size_t free_energy_begin
 Select only hills after this step.
size_t free_energy_end
 Select only hills before this step.
bool shift_fes
 Make the free energy surface be larger than or equal to zero.
cvm::real free_energy_offset
 Subtract this value from the free energy (after.
std::string free_energy_file
 Free energy surface output file.
std::string free_energy_gradients_file
 Free energy gradients output file.
std::string boltzmann_weights_file
 Boltzmann weights output file.
std::string boltzmann_counts_file
 Boltzmann counts output file.
cvm::real boltzmann_weights_scale
 Multiply the exponential weights by this constant.
cvm::real boltzmann_weights_temp
 Temperature (in K) for Boltzmann weights.
colvar_grid_scalarfree_energy
 Free energy values.
colvar_grid_gradientfree_energy_gradients
 Free energy gradients.
colvar_grid_scalarboltzmann_weights
 Boltzmann populations (real version).
colvar_grid_countboltzmann_counts
 Boltzmann populations (integer version).

Detailed Description

Metadynamics bias (implementation of ).

Definition at line 13 of file colvarbias_meta.h.


Member Typedef Documentation

typedef std::list<hill>::iterator colvarbias_meta::hill_iter
 

Definition at line 44 of file colvarbias_meta.h.

Referenced by calc_hills(), calc_hills_force(), create_hill(), read_hill(), and recount_hills_off_grid().


Member Enumeration Documentation

enum colvarbias_meta::Communication
 

Communication between different replicas.

Enumeration values:
single_replica  One replica (default).
multiple_replicas  Hills added concurrently by several replicas.

Definition at line 18 of file colvarbias_meta.h.

00018                      {
00020     single_replica,
00022     multiple_replicas
00023   };


Constructor & Destructor Documentation

colvarbias_meta::colvarbias_meta std::string const &  conf,
char const *  key
 

Constructor.

Definition at line 26 of file colvarbias_meta.C.

References b_hills_traj, comm, colvarmodule::debug(), dump_fes, dump_fes_save, expand_grids, colvarmodule::fatal_error(), GETCWD, grids_freq, hill_weight, hill_width, hills_energy, hills_energy_gradients, hills_traj_os, keep_hills, colvarmodule::log(), new_hill_freq, PATHSEP, PI, rebin_grids, register_replica_file(), replica, replica_files_registry, replica_out_file, replica_out_file_name, replica_update_freq, single_replica, and use_grids.

00027   : colvarbias (conf, key),
00028     comm (single_replica),
00029     new_hills_begin (hills.end())
00030 //     free_energy (NULL), free_energy_gradients (NULL),
00031 //     boltzmann_weights (NULL), boltzmann_counts (NULL)
00032 {
00033   if (cvm::n_abf_biases > 0)
00034     cvm::log ("Warning: ABF and metadynamics running at the "
00035               "same time can give inconsistent results.\n");
00036 
00037   get_keyval (conf, "hillWeight", hill_weight, 0.1);
00038   if (hill_weight == 0.0)
00039     cvm::log ("Warning: zero weight specified, "
00040               "this bias will have no effect.\n");
00041 
00042   get_keyval (conf, "newHillFrequency", new_hill_freq, 1000);
00043 
00044   get_keyval (conf, "hillWidth", hill_width, ::sqrt (2.0 * PI) / 2.0);
00045 
00046   {
00047     bool b_replicas = false;
00048     get_keyval (conf, "multipleReplicas", b_replicas, false);
00049     if (b_replicas)
00050       comm = multiple_replicas;
00051     else 
00052       comm = single_replica;
00053   }
00054 
00055   get_keyval (conf, "useGrids", use_grids, true);
00056   if (use_grids && (comm != single_replica)) {
00057     cvm::log ("Warning: calculations with a grid and "
00058               "multiple replicas is currently not supported: "
00059               "setting useGrids to \"no\".\n");
00060     use_grids = false;
00061   }
00062 
00063   if (use_grids) {
00064     get_keyval (conf, "gridsUpdateFrequency", grids_freq, new_hill_freq);
00065 
00066     get_keyval (conf, "rebinGrids", rebin_grids, false);
00067 
00068     expand_grids = false;
00069     for (size_t i = 0; i < colvars.size(); i++) {
00070       if (colvars[i]->expand_boundaries) {
00071         expand_grids = true;
00072         cvm::log ("Will expand the metadynamics grid when the colvar \""+
00073                   colvars[i]->name+"\" approaches the boundaries.\n");
00074       }
00075     }
00076 
00077     //    get_keyval (conf, "expandGrids", expand_grids, false);
00078 
00079     get_keyval (conf, "keepHills", keep_hills, false);
00080     get_keyval (conf, "dumpFreeEnergyFile", dump_fes, true);
00081     get_keyval (conf, "saveFreeEnergyFile", dump_fes_save, false);
00082 
00083     for (size_t i = 0; i < colvars.size(); i++) {
00084       colvars[i]->enable (colvar::task_grid);
00085     }
00086 
00087     hills_energy           = new colvar_grid_scalar   (colvars);
00088     hills_energy_gradients = new colvar_grid_gradient (colvars);
00089   }
00090 
00091   if (comm != single_replica) {
00092 
00093     get_keyval (conf, "replicaID", replica, std::string (""));
00094     if (!replica.size())
00095       cvm::fatal_error ("Error: you must define an id for this replica "
00096                         "when using more than one.\n");
00097 
00098     get_keyval (conf, "replicaFilesRegistry",
00099                 replica_files_registry, 
00100                 (this->name+".replica_files.txt"));
00101 
00102     get_keyval (conf, "replicaUpdateFrequency",
00103                 replica_update_freq, new_hill_freq);
00104 
00105     char *pwd = new char[321];
00106     if (GETCWD (pwd, 320) == NULL)
00107       cvm::fatal_error ("Error: cannot read the current working directory.\n");
00108     replica_out_file_name =
00109       (std::string (pwd)+std::string (PATHSEP)+
00110        cvm::output_prefix+".colvars."+this->name+
00111        "."+replica+".hills");
00112     delete pwd;
00113 
00114     replica_out_file.open (replica_out_file_name.c_str());
00115     if (!replica_out_file.good())
00116       cvm::fatal_error ("Error: in opening hills output file \""+
00117                         replica_out_file_name+"\".\n");
00118     register_replica_file (replica_out_file_name);
00119   }
00120 
00121   get_keyval (conf, "writeHillsTrajectory", b_hills_traj, false);
00122   if (b_hills_traj) {
00123 
00124     std::string const traj_file_name (cvm::output_prefix+
00125                                       ".colvars."+this->name+
00126                                       ( (comm != single_replica) ?
00127                                         ("."+replica) :
00128                                         ("") )+
00129                                       ".hills.traj");
00130     hills_traj_os.open (traj_file_name.c_str());
00131     if (!hills_traj_os.good())
00132       cvm::fatal_error ("Error: in opening hills output file \"" +
00133                         traj_file_name + "\".\n");
00134   }
00135               
00136   if (cvm::debug())
00137     cvm::log ("Done initializing the metadynamics bias \""+this->name+"\".\n");
00138 
00139   save_delimiters = false;
00140 }

colvarbias_meta::~colvarbias_meta  )  [virtual]
 

Destructor.

Definition at line 200 of file colvarbias_meta.C.

References hills_energy, hills_energy_gradients, hills_traj_os, and replica_out_file.

00201 {
00202   if (hills_energy) {
00203     delete hills_energy;
00204     hills_energy = NULL;
00205   }
00206 
00207   if (hills_energy_gradients) {
00208     delete hills_energy_gradients;
00209     hills_energy_gradients = NULL;
00210   }
00211 
00212 //   if (free_energy) {
00213 //     delete free_energy;
00214 //     free_energy = NULL;
00215 //   }
00216 
00217 //   if (free_energy_gradients) {
00218 //     delete free_energy_gradients;
00219 //     free_energy_gradients = NULL;
00220 //   }
00221 
00222 //   if (boltzmann_weights) {
00223 //     delete boltzmann_weights;
00224 //     boltzmann_weights = NULL;
00225 //   }
00226 
00227 //   if (boltzmann_counts) {
00228 //     delete boltzmann_counts;
00229 //     boltzmann_counts = NULL;
00230 //   }
00231 
00232   if (replica_out_file.good())
00233     replica_out_file.close();
00234 
00235   if (hills_traj_os.good())
00236     hills_traj_os.close();
00237 }


Member Function Documentation

void colvarbias_meta::analyse  )  [virtual]
 

Perform analysis.

Reimplemented from colvarbias.

Definition at line 1180 of file colvarbias_meta.C.

01181 {}

void colvarbias_meta::calc_hills hill_iter  h_first,
hill_iter  h_last,
cvm::real energy,
std::vector< colvarvalue > const &  values = std::vector< colvarvalue >(0)
[protected, virtual]
 

Calculate the values of the hills, incrementing colvar_energy.

Definition at line 525 of file colvarbias_meta.C.

References hill_iter, colvar_grid< T >::value(), and colvar_grid< T >::widths.

Referenced by project_hills(), and update().

00529 {
00530   std::vector<colvarvalue> curr_values (colvars.size());
00531   for (size_t i = 0; i < colvars.size(); i++) {
00532     curr_values[i].type (colvars[i]->type());
00533   }
00534 
00535   if (colvar_values.size()) {
00536     for (size_t i = 0; i < colvars.size(); i++) {
00537       curr_values[i] = colvar_values[i];
00538     }
00539   } else {
00540     for (size_t i = 0; i < colvars.size(); i++) {
00541       curr_values[i] = colvars[i]->value();
00542     }
00543   }
00544 
00545   for (hill_iter h = h_first; h != h_last; h++) {
00546       
00547     // compute the gaussian exponent
00548     cvm::real cv_sqdev = 0.0;
00549     for (size_t i = 0; i < colvars.size(); i++) {
00550       colvarvalue const &x  = curr_values[i];
00551       colvarvalue const &center = h->centers[i];
00552       cvm::real const    half_width = 0.5 * h->widths[i];
00553       cv_sqdev += (colvars[i]->dist2 (x, center)) / (half_width*half_width);
00554     }
00555 
00556     // compute the gaussian
00557     if (cv_sqdev > 23.0) {
00558       // set it to zero if the exponent is more negative than log(1.0E-05)
00559       h->value (0.0);
00560     } else {
00561       h->value (::exp (-0.5*cv_sqdev));
00562     }
00563     energy += h->energy();
00564   }
00565 }

void colvarbias_meta::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)
[protected, virtual]
 

Calculate the forces acting on the i-th colvar, incrementing colvar_forces[i]; must be called after calc_hills each time the values of the colvars are changed.

Definition at line 568 of file colvarbias_meta.C.

References hill_iter.

Referenced by project_hills(), and update().

00573 {
00574   // Retrieve the value of the colvar
00575   colvarvalue x (values.size() ? values[i].type() : colvars[i]->type());
00576   x = (values.size() ? values[i] : colvars[i]->value());
00577 
00578   // do the type check only once (all colvarvalues in the hills series
00579   // were already saved with their types matching those in the
00580   // colvars)
00581 
00582   switch (colvars[i]->type()) {
00583 
00584   case colvarvalue::type_scalar:
00585     for (hill_iter h = h_first; h != h_last; h++) {
00586       if (h->value() == 0.0) continue;
00587       colvarvalue const &center = h->centers[i];
00588       cvm::real const    half_width = 0.5 * h->widths[i];
00589       forces[i].real_value += 
00590         ( h->weight() * h->value() * (0.5 / (half_width*half_width)) * 
00591           (colvars[i]->dist2_lgrad (x, center)).real_value );
00592     }
00593     break;
00594 
00595   case colvarvalue::type_vector:
00596   case colvarvalue::type_unitvector:
00597     for (hill_iter h = h_first; h != h_last; h++) {
00598       if (h->value() == 0.0) continue;
00599       colvarvalue const &center = h->centers[i];
00600       cvm::real const    half_width = 0.5 * h->widths[i];
00601       forces[i].rvector_value +=
00602         ( h->weight() * h->value() * (0.5 / (half_width*half_width)) *
00603           (colvars[i]->dist2_lgrad (x, center)).rvector_value );
00604     }
00605     break;
00606 
00607   case colvarvalue::type_quaternion:
00608     for (hill_iter h = h_first; h != h_last; h++) {
00609       if (h->value() == 0.0) continue;
00610       colvarvalue const &center = h->centers[i];
00611       cvm::real const    half_width = 0.5 * h->widths[i];
00612       forces[i].quaternion_value +=
00613         ( h->weight() * h->value() * (0.5 / (half_width*half_width)) *
00614           (colvars[i]->dist2_lgrad (x, center)).quaternion_value );
00615     }
00616     break;
00617 
00618   case colvarvalue::type_notset:
00619     break;
00620   }
00621 }

std::list< colvarbias_meta::hill >::const_iterator colvarbias_meta::create_hill hill const &  h  )  [protected, virtual]
 

Add a new hill; if a .hills trajectory is written, write it there; if there is more than one replica, communicate it to the others.

Definition at line 246 of file colvarbias_meta.C.

References colvar_grid< T >::bin_distance_from_boundaries(), colvarbias_meta::hill::centers, comm, colvarmodule::debug(), colvarmodule::fatal_error(), hill_iter, hills, hills_energy, hills_off_grid, hills_traj_os, new_hills_begin, replica, and replica_out_file.

Referenced by update().

00247 {
00248   hill_iter const hills_end = hills.end();
00249   hills.push_back (h);
00250   if (new_hills_begin == hills_end) {
00251     // set the beginning of the "new" hills
00252     new_hills_begin = hills.end();
00253     new_hills_begin--;
00254   }
00255 
00256   if (use_grids) {
00257 
00258     // also add it to the list of hills that are off-grid, which
00259     // receive special treatment (i.e. they are computed analytically)
00260 
00261     cvm::real const min_dist = hills_energy->bin_distance_from_boundaries (h.centers);
00262     if (min_dist < (3.0 * ::floor (hill_width)) + 1.0) {
00263       hills_off_grid.push_back (hills.back());
00264     }
00265   }
00266 
00267   // output to trajectory
00268   if (hills_traj_os.good()) {
00269     hills_traj_os << (hills.back()).output_traj();
00270     if (cvm::debug()) {
00271       hills_traj_os.flush();
00272     }
00273   }
00274 
00275   // output to the replica file
00276   if (comm != single_replica) {
00277     if (replica_out_file.good()) {
00278       replica_out_file << hills.back();
00279     } else {
00280       cvm::fatal_error ("Error in writing to the hills database for replica \""+
00281                         replica+"\".\n");
00282     }
00283   }
00284 
00285   return hills.end();
00286 }

virtual std::list<hill>::const_iterator colvarbias_meta::delete_hill std::list< hill >::iterator &  h  )  [protected, virtual]
 

Remove a previously saved hill (returns an iterator for the next hill in the list).

void colvarbias_meta::project_hills hill_iter  h_first,
hill_iter  h_last,
colvar_grid_scalar ge,
colvar_grid_gradient gf,
cvm::real const   scale_factor = 1.0
[protected]
 

Project the selected hills onto grids.

Definition at line 1115 of file colvarbias_meta.C.

References colvar_grid_gradient::acc_force(), colvar_grid_scalar::acc_value(), colvar_grid< T >::bin_to_value_scalar(), calc_hills(), calc_hills_force(), colvarmodule::debug(), hills, hills_energy, colvar_grid< T >::incr(), colvar_grid< T >::index_ok(), colvarmodule::log(), and colvar_grid< T >::new_index().

Referenced by read_restart(), and update().

01120 {
01121   if (cvm::debug())
01122     cvm::log ("Projecting hills.\n");
01123 
01124   std::vector<colvarvalue> colvar_values (colvars.size());
01125   std::vector<cvm::real> colvar_forces_scalar (colvars.size());
01126 
01127   std::vector<int> he_ix = he->new_index();
01128   std::vector<int> hg_ix = (hg != NULL) ? hg->new_index() : std::vector<int> (0);
01129   cvm::real hills_energy_here = 0.0;
01130   std::vector<colvarvalue> hills_forces_here (colvars.size(), 0.0);
01131   // loop over the points of the grid
01132   for ( ;
01133         he->index_ok (he_ix) && ((hg != NULL) ? hg->index_ok (hg_ix) : true);
01134         ) {
01135 
01136     for (size_t i = 0; i < colvars.size(); i++) {
01137       colvar_values[i] = hills_energy->bin_to_value_scalar (he_ix[i], i);
01138     }
01139 
01140     // loop over the hills and increment the energy grid locally
01141     hills_energy_here = 0.0;
01142     calc_hills (h_first, h_last, hills_energy_here, colvar_values);
01143     he->acc_value (he_ix, scale_factor * hills_energy_here);
01144 
01145     if (hg != NULL) {
01146       for (size_t i = 0; i < colvars.size(); i++) {
01147         hills_forces_here[i].reset();
01148         calc_hills_force (i, h_first, h_last, hills_forces_here, colvar_values);
01149         colvar_forces_scalar[i] = scale_factor * hills_forces_here[i].real_value;
01150       }
01151       hg->acc_force (hg_ix, &(colvar_forces_scalar.front()));
01152     }
01153 
01154     he->incr (he_ix);
01155     if (hg != NULL)
01156       hg->incr (hg_ix);
01157   }
01158 
01159   if (! keep_hills) {
01160     hills.erase (hills.begin(), hills.end());
01161   }
01162 }

std::istream & colvarbias_meta::read_hill std::istream &  is  )  [protected]
 

Read a hill from a file.

Definition at line 904 of file colvarbias_meta.C.

References colvar_grid< T >::bin_distance_from_boundaries(), comm, hill_iter, hill_weight, hills, hills_energy, hills_off_grid, colvarparse::key_lookup(), new_hills_begin, and replica.

Referenced by read_replica_files(), and read_restart().

00905 {
00906   if (!is) return is; // do nothing if failbit is set
00907 
00908   size_t const start_pos = is.tellg();
00909 
00910   std::string data;
00911   if ( !(is >> read_block ("hill", data)) ) {
00912     is.clear();
00913     is.seekg (start_pos, std::ios::beg);
00914     is.setstate (std::ios::failbit);
00915     return is;
00916   }
00917 
00918   std::string h_replica = "";
00919 
00920   int h_it;
00921   get_keyval (data, "step", h_it, 0, parse_silent);
00922   cvm::real h_weight;
00923   get_keyval (data, "weight", h_weight, hill_weight, parse_silent);
00924 
00925   std::vector<colvarvalue> h_centers (colvars.size());
00926   for (size_t i = 0; i < colvars.size(); i++) {
00927     h_centers[i].type ((colvars[i]->value()).type()); 
00928   }
00929   {
00930     // it is safer to read colvarvalue objects one at a time; XXX
00931     // todo: check to remove it later
00932     std::string centers_input;
00933     key_lookup (data, "centers", centers_input);
00934     std::istringstream centers_is (centers_input);
00935     for (size_t i = 0; i < colvars.size(); i++) {
00936       centers_is >> h_centers[i];
00937     }
00938   }
00939 
00940   std::vector<cvm::real> h_widths (colvars.size());
00941   get_keyval (data, "widths", h_widths, std::vector<cvm::real> (colvars.size(), 1.0), parse_silent);
00942   
00943   if (comm != single_replica)
00944     get_keyval (data, "replica", h_replica, replica, parse_silent);
00945 
00946   hill_iter const hills_end = hills.end();
00947   hills.push_back (hill (h_it, h_weight, h_centers, h_widths, h_replica));
00948   if (new_hills_begin == hills_end) {
00949     // set the beginning of the "new" hills
00950     new_hills_begin = hills.end();
00951     new_hills_begin--;
00952   }
00953 
00954   if (use_grids) {
00955 
00956     // add this also to the list of hills that are off-grid, which will
00957     // be computed analytically
00958     cvm::real const min_dist =
00959       hills_energy->bin_distance_from_boundaries ((hills.back()).centers);
00960     if (min_dist < (3.0 * ::floor (hill_width)) + 1.0) {
00961       hills_off_grid.push_back (hills.back());
00962     }
00963   }
00964 
00965   return is;
00966 }

void colvarbias_meta::read_replica_files  )  [protected, virtual]
 

Read new data from replicas' files (called only when != ).

Definition at line 705 of file colvarbias_meta.C.

References colvarmodule::debug(), colvarmodule::fatal_error(), hills, hills_traj_os, colvarmodule::log(), read_hill(), replica_files, replica_files_pos, and colvarmodule::to_str().

Referenced by read_restart(), and update().

00706 {
00707 
00708   // read hills from the other replicas' files; for each file, resume
00709   // the position recorded last time
00710   std::list<std::string>::iterator rsi = replica_files.begin();
00711   std::list<size_t>::iterator      pi  = replica_files_pos.begin();
00712   for ( ; rsi != replica_files.end(); pi++, rsi++) {
00713 
00714     if (cvm::debug()) 
00715       cvm::log ("Checking for new hills in the replica file \""+
00716                 (*rsi)+"\".\n");
00717 
00718     std::ifstream is ((*rsi).c_str());
00719     if (is.good()) {
00720 
00721       is.seekg ((*pi), std::ios::beg);
00722       while (read_hill (is)) {
00723 
00724         if (cvm::debug())
00725           cvm::log ("Read a previously saved hill under the "
00726                     "metadynamics bias \""+
00727                     this->name+"\", created by replica "+
00728                     (hills.back()).replica+
00729                     " at step "+
00730                     cvm::to_str ((hills.back()).it)+".\n");
00731 
00732         // output to trajectory
00733         if (hills_traj_os.good()) {
00734           hills_traj_os << (hills.back()).output_traj();
00735           if (cvm::debug()) {
00736             hills_traj_os.flush();
00737           }
00738         }
00739       }
00740       is.clear();
00741       // store it for the next read
00742       *pi = is.tellg();
00743       if (cvm::debug())
00744         cvm::log ("Stopped reading at position "+
00745                   cvm::to_str (*pi)+".\n");
00746     } else {
00747       cvm::fatal_error ("Error: cannot read from file \""+(*rsi)+"\".\n");
00748     }
00749 
00750     is.close();
00751   }
00752 
00753 }

std::istream & colvarbias_meta::read_restart std::istream &  is  )  [virtual]
 

Read the bias configuration from a restart file.

Implements colvarbias.

Definition at line 760 of file colvarbias_meta.C.

References comm, colvarmodule::debug(), colvarmodule::fatal_error(), hills, hills_energy, hills_energy_gradients, keep_hills, colvarmodule::log(), colvar_grid< T >::map_grid(), new_hills_begin, project_hills(), read_hill(), read_replica_files(), colvar_grid_gradient::read_restart(), colvar_grid_scalar::read_restart(), recount_hills_off_grid(), and update_replica_files_registry().

00761 {
00762   size_t const start_pos = is.tellg();
00763 
00764   cvm::log ("Restarting metadynamics bias \""+
00765             this->name+"\".\n");
00766   std::string key, brace, conf;
00767   if ( !(is >> key)   || !(key == "metadynamics") ||
00768        !(is >> brace) || !(brace == "{") ||
00769        !(is >> colvarparse::read_block ("configuration", conf)) ) {
00770 
00771     cvm::log ("Error: in reading restart configuration for metadynamics bias \""+
00772               this->name+"\" at position "+
00773               cvm::to_str (start_pos)+" in stream.\n");
00774     is.clear();
00775     is.seekg (start_pos, std::ios::beg);
00776     is.setstate (std::ios::failbit);
00777     return is;
00778   }
00779 
00780   std::string name = "";
00781   if ( colvarparse::get_keyval (conf, "name", name, std::string (""), colvarparse::parse_silent) &&
00782        (name != this->name) )
00783     cvm::fatal_error ("Error: in the restart file, the "
00784                       "\"metadynamics\" block has a wrong name: different system?\n");
00785 
00786   if (name.size() == 0) {
00787     cvm::fatal_error ("Error: \"metadynamics\" block in the restart file "
00788                       "has no identifiers.\n");
00789   }
00790 
00791   if (this->comm == single_replica) {
00792 
00793     if (use_grids) {
00794 
00795       if (expand_grids) {
00796         // the boundaries of the colvars may have been changed (note:
00797         // this second reallocation may be deleted when the new
00798         // restart format for the grids has kicked in)
00799         delete hills_energy;
00800         delete hills_energy_gradients;
00801         hills_energy = new colvar_grid_scalar (colvars);
00802         hills_energy_gradients = new colvar_grid_gradient (colvars);
00803       }
00804 
00805       if ( !(is >> key) ||
00806            !(key == std::string ("hills_energy")) ||
00807            !(hills_energy->read_restart (is)) ) {
00808         cvm::log ("Error: in reading restart information for metadynamics bias \""+
00809                   this->name+"\".\n");
00810         is.clear();
00811         is.seekg (start_pos, std::ios::beg);
00812         is.setstate (std::ios::failbit);
00813         return is;
00814       }
00815 
00816       if ( !(is >> key) ||
00817            !(key == std::string ("hills_energy_gradients")) || 
00818            !(hills_energy_gradients->read_restart (is)) ) {
00819         cvm::log ("Error: in reading restart information for metadynamics bias \""+
00820                   this->name+"\".\n");
00821         is.clear();
00822         is.seekg (start_pos, std::ios::beg);
00823         is.setstate (std::ios::failbit);
00824         return is;
00825       }
00826 
00827     }
00828 
00829     // read the hills explicitly (if there are any)
00830     while (this->read_hill (is)) {
00831       if (cvm::debug()) 
00832         cvm::log ("Read a previously saved hill under the "
00833                   "metadynamics bias \""+
00834                   this->name+"\", created at step "+
00835                   cvm::to_str ((hills.back()).it)+".\n");
00836     }
00837     is.clear();
00838     new_hills_begin = hills.end();
00839 
00840     if (rebin_grids) {
00841 
00842       // allocate new grids, based on the CURRENT boundaries and
00843       // widths, as per the restarts, and project the grids from the
00844       // restart file onto them
00845       colvar_grid_scalar   *new_hills_energy =
00846         new colvar_grid_scalar (colvars);
00847       colvar_grid_gradient *new_hills_energy_gradients =
00848         new colvar_grid_gradient (colvars);
00849       
00850       if (keep_hills && (hills.size() > 0)) {
00851 
00852         // if there are hills, recompute the new grids from them
00853         cvm::log ("rebinGrids and keepHills defined, recomputing "
00854                   "analytically the energy and force grids; this may take a while...\n");
00855         project_hills (hills.begin(), hills.end(),
00856                        new_hills_energy, new_hills_energy_gradients);
00857         cvm::log ("rebinning done.\n");
00858 
00859       } else {
00860 
00861         // otherwise, use the grids in the restart file
00862 
00863         cvm::log ("rebinGrids defined, mapping energy and forces grid to new grids.\n");
00864         new_hills_energy->map_grid (*hills_energy);
00865         new_hills_energy_gradients->map_grid (*hills_energy_gradients);
00866       }
00867 
00868       delete hills_energy;
00869       delete hills_energy_gradients;
00870       hills_energy = new_hills_energy;
00871       hills_energy_gradients = new_hills_energy_gradients;
00872 
00873       if (hills.size())
00874         recount_hills_off_grid (hills.begin(), hills.end(), hills_energy);
00875     }
00876 
00877   } else {
00878 
00879     // read all hills from the replica files
00880     cvm::log ("Multiple replicas have been defined for the "
00881               "metadynamics bias \""+this->name+"\", skipping "
00882               "the restart file and reading hills from the replica "
00883               "files, whose list is contained in \""+
00884               this->replica_files_registry+"\".\n");
00885     this->update_replica_files_registry();
00886     this->read_replica_files();
00887   }
00888 
00889   is >> brace;
00890   if (brace != "}") {
00891     cvm::fatal_error ("Error: corrupted restart information for metadynamics bias \""+
00892                       this->name+"\": no matching brace at position "+
00893                       cvm::to_str (is.tellg())+" in the restart file.\n");
00894     is.setstate (std::ios::failbit);
00895   }
00896 
00897   if (cvm::debug())
00898     cvm::log ("colvarbias_meta::read_restart() done\n");
00899 
00900   return is;
00901 }  

void colvarbias_meta::recount_hills_off_grid hill_iter  h_first,
hill_iter  h_last,
colvar_grid_scalar ge
[protected]
 

Regenerate the hills_off_grid list.

Definition at line 1165 of file colvarbias_meta.C.

References colvar_grid< T >::bin_distance_from_boundaries(), hill_iter, hills_energy, and hills_off_grid.

Referenced by read_restart().

01168 {
01169   hills_off_grid.clear();
01170 
01171   for (hill_iter h = h_first; h != h_last; h++) {
01172     cvm::real const min_dist = hills_energy->bin_distance_from_boundaries (h->centers);
01173     if (min_dist < (3.0 * ::floor (hill_width)) + 1.0) {
01174       hills_off_grid.push_back (*h);
01175     }
01176   }
01177 }

void colvarbias_meta::register_replica_file std::string const &  new_file  )  [protected, virtual]
 

Add this replica to the registry (called only when != ).

Definition at line 629 of file colvarbias_meta.C.

References replica_files_registry.

Referenced by colvarbias_meta().

00630 {
00631   std::ifstream reg_is (replica_files_registry.c_str());
00632   if (reg_is.good()) {
00633     std::string s;
00634     while (reg_is >> s)
00635       if (new_file == s) {
00636         // already found this replica in registry
00637         return;
00638       }
00639   }
00640   reg_is.close();
00641 
00642   std::ofstream reg_os (replica_files_registry.c_str(), std::ios::app);
00643   reg_os << new_file << "\n";
00644   reg_os.close();
00645 }

void colvarbias_meta::update  )  [virtual]
 

Retrieve colvar values and calculate their biasing forces.

Implements colvarbias.

Definition at line 323 of file colvarbias_meta.C.

References calc_hills(), calc_hills_force(), comm, colvar_grid< T >::create(), create_hill(), colvarmodule::debug(), colvar_grid< T >::get_colvars_index(), hill_weight, hill_width, hills, hills_energy, hills_energy_gradients, hills_off_grid, colvar_grid< T >::index_ok(), colvarmodule::log(), colvar_grid< T >::lower_boundaries, colvar_grid< T >::map_grid(), multiple_replicas, new_hills_begin, project_hills(), read_replica_files(), replica, replica_out_file, single_replica, colvar_grid< T >::sizes(), colvarmodule::step_relative(), update_replica_files_registry(), colvar_grid< T >::upper_boundaries, and colvar_grid< T >::value().

00324 {
00325   if (cvm::debug())
00326     cvm::log ("Updating the metadynamics bias \""+this->name+"\".\n");
00327 
00328   if (use_grids) {
00329 
00330     std::vector<int> curr_bin = hills_energy->get_colvars_index();
00331 
00332     if (expand_grids) {
00333 
00334       // first of all, expand the grids, if specified
00335       if (cvm::debug())
00336         cvm::log ("Current coordinates on the grid: "+
00337                   cvm::to_str (curr_bin)+".\n");
00338 
00339       if (cvm::debug())
00340         cvm::log ("Checking if the grid is big enough.\n");
00341 
00342       bool changed_grids = false;
00343       int const min_buffer =
00344         (3 * (size_t) ::floor (hill_width)) + 1;
00345 
00346       std::vector<int>         new_sizes (hills_energy->sizes());
00347       std::vector<colvarvalue> new_lower_boundaries (hills_energy->lower_boundaries);
00348       std::vector<colvarvalue> new_upper_boundaries (hills_energy->upper_boundaries);
00349 
00350       for (size_t i = 0; i < colvars.size(); i++) {
00351 
00352         if (! colvars[i]->expand_boundaries)
00353           continue;
00354 
00355         cvm::real &new_lb   = new_lower_boundaries[i].real_value;
00356         cvm::real &new_ub   = new_upper_boundaries[i].real_value;
00357         int       &new_size = new_sizes[i];
00358         bool changed_lb = false, changed_ub = false;
00359 
00360         if (curr_bin[i] < min_buffer) {
00361           int const extra_points = (min_buffer - curr_bin[i]);
00362           new_lb -= extra_points * colvars[i]->width;
00363           new_size += extra_points;
00364           // changed offset in this direction => the pointer needs to
00365           // be changed, too
00366           curr_bin[i] += extra_points;
00367 
00368           changed_lb = true;
00369           cvm::log ("Metadynamics"+((cvm::n_meta_biases > 1) ? " \""+name+"\"" : "")+
00370                     ": new lower boundary for colvar \""+
00371                     colvars[i]->name+"\", at "+
00372                     cvm::to_str (new_lower_boundaries[i])+".\n");
00373         }
00374 
00375         if (curr_bin[i] > new_size - min_buffer - 1) {
00376           int const extra_points = (curr_bin[i] - (new_size - 1) + min_buffer);
00377           new_ub += extra_points * colvars[i]->width;
00378           new_size += extra_points;
00379 
00380           changed_ub = true;
00381           cvm::log ("Metadynamics"+((cvm::n_meta_biases > 1) ? " \""+name+"\"" : "")+
00382                     ": new upper boundary for colvar \""+
00383                     colvars[i]->name+"\", at "+
00384                     cvm::to_str (new_upper_boundaries[i])+".\n");
00385         }
00386 
00387         if (changed_lb || changed_ub)
00388           changed_grids = true;
00389       }
00390 
00391       if (changed_grids) {
00392 
00393         // map everything into new grids
00394         if (cvm::debug())
00395           cvm::log ("Expanding the grids for the metadynamics bias \""+
00396                     name+"\".\n");
00397 
00398         colvar_grid_scalar *new_hills_energy =
00399           new colvar_grid_scalar (*hills_energy);
00400         colvar_grid_gradient *new_hills_energy_gradients =
00401           new colvar_grid_gradient (*hills_energy_gradients);
00402 
00403         // supply new boundaries to the new grids
00404 
00405         new_hills_energy->lower_boundaries = new_lower_boundaries;
00406         new_hills_energy->upper_boundaries = new_upper_boundaries;
00407         new_hills_energy->create (new_sizes, 0.0, 1);
00408 
00409         new_hills_energy_gradients->lower_boundaries = new_lower_boundaries;
00410         new_hills_energy_gradients->upper_boundaries = new_upper_boundaries;
00411         new_hills_energy_gradients->create (new_sizes, 0.0, colvars.size());
00412 
00413         new_hills_energy->map_grid (*hills_energy);
00414         new_hills_energy_gradients->map_grid (*hills_energy_gradients);
00415 
00416         delete hills_energy;
00417         delete hills_energy_gradients;
00418         hills_energy = new_hills_energy;
00419         hills_energy_gradients = new_hills_energy_gradients;
00420 
00421         curr_bin = hills_energy->get_colvars_index();
00422         if (cvm::debug())
00423           cvm::log ("Coordinates on the new grid: "+
00424                     cvm::to_str (curr_bin)+".\n");
00425       }
00426     }
00427   }
00428 
00429   // add a new hill if the required time interval has passed
00430   if ((cvm::it % new_hill_freq) == 0) {
00431 
00432     if (cvm::debug())
00433       cvm::log ("Adding a new hill under the bias \""+
00434                 this->name+"\", at step "+cvm::to_str (cvm::it)+".\n");
00435 
00436     switch (comm) {
00437 
00438     case single_replica:
00439       create_hill (hill (hill_weight, colvars, hill_width));
00440       break;
00441 
00442     case multiple_replicas:
00443       create_hill (hill (hill_weight, colvars, hill_width, replica));
00444       break;
00445     }
00446   }
00447 
00448   // syncronise with the other replicas if specified
00449   if (comm != single_replica) {
00450     if ((cvm::it % replica_update_freq) == 0) {
00451       // the buffer should be always emptied to keep the other
00452       // replicas as much in sync as possible
00453       replica_out_file.flush();
00454       // read all replica files, except those for this replica; hills
00455       // from previous run were already loaded in memory by restarts
00456       update_replica_files_registry();
00457       read_replica_files();
00458     }
00459   }
00460 
00461   // calculate the biasing energy and forces
00462   colvar_energy = 0.0;
00463   for (size_t i = 0; i < colvars.size(); i++) {
00464     colvar_forces[i].reset();
00465   }
00466 
00467 
00468   if (use_grids) {
00469 
00470     // get the forces from the grid
00471 
00472     if ((cvm::step_relative() % grids_freq) == 0) {
00473       // map the most recent gaussians to the grids
00474       project_hills (new_hills_begin, hills.end(),
00475                      hills_energy,    hills_energy_gradients);
00476       new_hills_begin = hills.end();
00477     }
00478 
00479     std::vector<int> curr_bin = hills_energy->get_colvars_index();
00480     if (cvm::debug())
00481       cvm::log ("Current coordinates on the grid: "+
00482                 cvm::to_str (curr_bin)+".\n");
00483 
00484     if (hills_energy->index_ok (curr_bin)) {
00485 
00486       // within the grid: add the energy and the forces from there
00487 
00488       colvar_energy += hills_energy->value (curr_bin);
00489       for (size_t i = 0; i < colvars.size(); i++) {
00490         cvm::real const *f = &(hills_energy_gradients->value (curr_bin));
00491         colvar_forces[i].real_value += -1.0 * f[i]; // the gradients
00492                                                     // are stored, not
00493                                                     // the forces
00494       }
00495 
00496     } else {
00497 
00498       // we're off the grid, computing analytically only the hills
00499       // within range
00500 
00501       calc_hills (hills_off_grid.begin(), hills_off_grid.end(), colvar_energy);
00502       for (size_t i = 0; i < colvars.size(); i++) {
00503         calc_hills_force (i, hills_off_grid.begin(), hills_off_grid.end(), colvar_forces);
00504       }
00505     }
00506 
00507   } else {
00508 
00509     // calculate the current value of each hill and add it to colvar_energy
00510     calc_hills (hills.begin(), hills.end(), colvar_energy);
00511   
00512     // calculate the current derivatives of each hill and add them to
00513     // colvar_forces
00514     for (size_t i = 0; i < colvars.size(); i++) {
00515       calc_hills_force (i, hills.begin(), hills.end(), colvar_forces);
00516     }
00517   }
00518 
00519   if (cvm::debug()) 
00520     cvm::log ("Hills energy = "+cvm::to_str (colvar_energy)+
00521               ", hills forces = "+cvm::to_str (colvar_forces)+".\n");
00522 }

void colvarbias_meta::update_replica_files_registry  )  [protected, virtual]
 

Read the names all replica output files which have been added, skipping the one which is currently being written to (called only when != ).

Definition at line 648 of file colvarbias_meta.C.

References colvarmodule::debug(), colvarmodule::fatal_error(), colvarmodule::log(), replica_files, replica_files_pos, and replica_files_registry.

Referenced by read_restart(), and update().

00649 {
00650   if (cvm::debug())
00651     cvm::log ("Updating the list of replica files, currently containing "+
00652               cvm::to_str (replica_files.size())+" elements.\n");
00653 
00654   // update the list of replica files by reading those created after restarting
00655   std::ifstream reg_is (replica_files_registry.c_str());
00656 
00657   if (reg_is.good()) {
00658 
00659     std::string s ("");
00660     while ((reg_is >> s) && s.size()) {
00661 
00662       if (s == replica_out_file_name) {
00663         // this is the file from this same job, skip it
00664         if (cvm::debug())
00665           cvm::log ("Skipping this job's replica file, \""+s+"\".\n");
00666         s = "";
00667         continue;
00668       }
00669 
00670       bool already_loaded = false;
00671       std::list<std::string>::iterator rsi = replica_files.begin();
00672       for ( ; rsi != replica_files.end(); rsi++) {
00673         if (s == *rsi) {
00674           // this file was already added to the list
00675           if (cvm::debug())
00676             cvm::log ("Skipping a replica file already loaded, \""+(*rsi)+"\".\n");
00677           already_loaded = true;;
00678         }
00679       }
00680 
00681       if ( (replica_files.size() == 0) || (!already_loaded) ) {
00682         // add this file to the registry
00683         if (cvm::debug())
00684           cvm::log ("New replica file found: \""+s+"\".\n");
00685         replica_files.push_back (s);
00686         replica_files_pos.push_back (0);
00687       }
00688 
00689     }
00690 
00691     s = "";
00692   } else {
00693     cvm::fatal_error ("Error: cannot read the replica files registry, \""+
00694                       replica_files_registry+"\".\n");
00695   }
00696 
00697   reg_is.close();
00698 
00699   if (cvm::debug())
00700     cvm::log ("The list of replica files now contains "+
00701               cvm::to_str (replica_files.size())+" elements.\n");
00702 }

std::ostream & colvarbias_meta::write_restart std::ostream &  os  )  [virtual]
 

Write the bias configuration to a restart file.

Implements colvarbias.

Definition at line 973 of file colvarbias_meta.C.

References colvar_grid< T >::add_constant(), comm, dump_fes_save, hills, hills_energy, hills_energy_gradients, hills_off_grid, colvarmodule::log(), colvar_grid_scalar::maximum_value(), colvar_grid< T >::multiply_constant(), replica, replica_out_file_name, colvar_grid< T >::write_multicol(), colvar_grid_gradient::write_restart(), and colvar_grid_scalar::write_restart().

00974 {
00975   os << "metadynamics {\n"
00976      << "  configuration {\n"
00977     //      << "    id " << this->id << "\n"
00978      << "    name " << this->name << "\n";
00979   if (this->comm != single_replica)
00980     os << "    replica " << this->replica << "\n";
00981   os << "  }\n";
00982 
00983   if (this->comm == single_replica) {
00984 
00985     if (use_grids) {
00986 
00987       os << "  hills_energy\n";
00988       hills_energy->write_restart (os);
00989       os << "  hills_energy_gradients\n";
00990       hills_energy_gradients->write_restart (os);
00991 
00992       if (dump_fes) {
00993         cvm::real const max = hills_energy->maximum_value();
00994         hills_energy->add_constant (-1.0 * max);
00995         hills_energy->multiply_constant (-1.0);
00996         // if this is the only free energy integrator, the pmf file
00997         // name is general, otherwise there is a label with the bias
00998         // name
00999         std::string const fes_file_name =
01000           ((cvm::n_meta_biases == 1) && (cvm::n_abf_biases == 0)) ?
01001           std::string (cvm::output_prefix+
01002                        (dump_fes_save ? "."+cvm::to_str (cvm::step_absolute()) : "")+
01003                        ".pmf") :
01004           std::string (cvm::output_prefix+"."+this->name+
01005                        (dump_fes_save ? "."+cvm::to_str (cvm::step_absolute()) : "")+
01006                        ".pmf");
01007         std::ofstream fes_os (fes_file_name.c_str());
01008         hills_energy->write_multicol (fes_os);
01009 
01010         // restore the grid to original values
01011         hills_energy->multiply_constant (-1.0);
01012         hills_energy->add_constant (max);
01013       }
01014 
01015     }
01016 
01017     if ( (!use_grids) || keep_hills ) {
01018 
01019       // write all the hills in memory
01020       for (std::list<hill>::const_iterator h = this->hills.begin();
01021            h != this->hills.end();
01022            h++) {
01023         os << *h;
01024       }
01025 
01026     } else {
01027 
01028       // write just those that need to be computed analytically
01029       for (std::list<hill>::const_iterator h = this->hills_off_grid.begin();
01030            h != this->hills_off_grid.end();
01031            h++) {
01032         os << *h;
01033       }
01034     }
01035 
01036   } else {
01037     cvm::log ("Hills from the metadynamics bias \""+
01038               this->name+"\" have been previously saved in the file \""+
01039               replica_out_file_name+"\".\n");
01040   }
01041 
01042   os << "}\n\n";
01043   return os;
01044 }  


Member Data Documentation

bool colvarbias_meta::b_hills_traj [protected]
 

Write the hill logfile.

Definition at line 148 of file colvarbias_meta.h.

Referenced by colvarbias_meta().

colvar_grid_count* colvarbias_meta::boltzmann_counts [protected]
 

Boltzmann populations (integer version).

Definition at line 240 of file colvarbias_meta.h.

std::string colvarbias_meta::boltzmann_counts_file [protected]
 

Boltzmann counts output file.

These are discretized Boltzmann weights, useful e.g. to initialize the samples of an ABF calculation. boltzmann_weights_scale must be much larger than 1 for a proper discretization

Definition at line 220 of file colvarbias_meta.h.

colvar_grid_scalar* colvarbias_meta::boltzmann_weights [protected]
 

Boltzmann populations (real version).

Definition at line 237 of file colvarbias_meta.h.

std::string colvarbias_meta::boltzmann_weights_file [protected]
 

Boltzmann weights output file.

Definition at line 212 of file colvarbias_meta.h.

cvm::real colvarbias_meta::boltzmann_weights_scale [protected]
 

Multiply the exponential weights by this constant.

Definition at line 225 of file colvarbias_meta.h.

cvm::real colvarbias_meta::boltzmann_weights_temp [protected]
 

Temperature (in K) for Boltzmann weights.

Definition at line 228 of file colvarbias_meta.h.

Communication colvarbias_meta::comm
 

Communication between different replicas.

Definition at line 26 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), create_hill(), read_hill(), read_restart(), update(), and write_restart().

bool colvarbias_meta::dump_fes [protected]
 

Dump the free energy surface (.pmf file) every restartFrequency.

Definition at line 119 of file colvarbias_meta.h.

Referenced by colvarbias_meta().

bool colvarbias_meta::dump_fes_save [protected]
 

Keep the free energy surface files at different iterations, appending a step number to each.

Definition at line 123 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), and write_restart().

bool colvarbias_meta::expand_grids [protected]
 

Should the grids be expanded if necessary?

Definition at line 109 of file colvarbias_meta.h.

Referenced by colvarbias_meta().

colvar_grid_scalar* colvarbias_meta::free_energy [protected]
 

Free energy values.

Definition at line 231 of file colvarbias_meta.h.

size_t colvarbias_meta::free_energy_begin [protected]
 

Select only hills after this step.

Definition at line 194 of file colvarbias_meta.h.

size_t colvarbias_meta::free_energy_end [protected]
 

Select only hills before this step.

Definition at line 197 of file colvarbias_meta.h.

std::string colvarbias_meta::free_energy_file [protected]
 

Free energy surface output file.

Definition at line 206 of file colvarbias_meta.h.

colvar_grid_gradient* colvarbias_meta::free_energy_gradients [protected]
 

Free energy gradients.

Definition at line 234 of file colvarbias_meta.h.

std::string colvarbias_meta::free_energy_gradients_file [protected]
 

Free energy gradients output file.

Definition at line 209 of file colvarbias_meta.h.

cvm::real colvarbias_meta::free_energy_offset [protected]
 

Subtract this value from the free energy (after.

Definition at line 203 of file colvarbias_meta.h.

size_t colvarbias_meta::grids_freq [protected]
 

How often the hills should be projected onto the grids.

Definition at line 112 of file colvarbias_meta.h.

Referenced by colvarbias_meta().

cvm::real colvarbias_meta::hill_weight [protected]
 

Height of new hills.

Definition at line 99 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), read_hill(), and update().

cvm::real colvarbias_meta::hill_width [protected]
 

width of a hill

The local width of each collective variable, multiplied by this number, provides the hill width along that direction

Definition at line 142 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), and update().

std::list<hill> colvarbias_meta::hills [protected]
 

Parse analysis tasks and options List of hills used on this bias (total); if a grid is employed, these don't need to be updated at every time step.

Definition at line 53 of file colvarbias_meta.h.

Referenced by create_hill(), project_hills(), read_hill(), read_replica_files(), read_restart(), update(), and write_restart().

colvar_grid_scalar* colvarbias_meta::hills_energy [protected]
 

Hill energy, cached on a grid.

Definition at line 126 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), create_hill(), project_hills(), read_hill(), read_restart(), recount_hills_off_grid(), update(), write_restart(), and ~colvarbias_meta().

colvar_grid_gradient* colvarbias_meta::hills_energy_gradients [protected]
 

Hill forces, cached on a grid.

Definition at line 129 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), read_restart(), update(), write_restart(), and ~colvarbias_meta().

std::list<hill> colvarbias_meta::hills_off_grid [protected]
 

List of hills used on this bias that are on the boundary edges; these are updated regardless of whether hills are used.

Definition at line 61 of file colvarbias_meta.h.

Referenced by create_hill(), read_hill(), recount_hills_off_grid(), update(), and write_restart().

std::ofstream colvarbias_meta::hills_traj_os [protected]
 

Logfile of hill management (creation and deletion).

Definition at line 150 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), create_hill(), read_replica_files(), and ~colvarbias_meta().

bool colvarbias_meta::keep_hills [protected]
 

Whether to keep the hills in the restart file (e.g. to do meaningful accurate rebinning afterwards).

Definition at line 116 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), and read_restart().

size_t colvarbias_meta::new_hill_freq [protected]
 

Number of simulation steps between two hills.

Definition at line 145 of file colvarbias_meta.h.

Referenced by colvarbias_meta().

hill_iter colvarbias_meta::new_hills_begin [protected]
 

Iterator to the first of the "newest" hills (when using grids, those who haven't been mapped yet).

Definition at line 57 of file colvarbias_meta.h.

Referenced by create_hill(), read_hill(), read_restart(), and update().

hill_iter colvarbias_meta::new_hills_off_grid_begin [protected]
 

Same as new_hills_begin, but for the off-grid ones.

Definition at line 64 of file colvarbias_meta.h.

bool colvarbias_meta::rebin_grids [protected]
 

Rebin the hills upon restarting.

Definition at line 106 of file colvarbias_meta.h.

Referenced by colvarbias_meta().

std::string colvarbias_meta::replica [protected]
 

Identifier for this replica.

Definition at line 154 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), create_hill(), read_hill(), update(), and write_restart().

std::list<std::string> colvarbias_meta::replica_files [protected]
 

Replicas file names.

Definition at line 179 of file colvarbias_meta.h.

Referenced by read_replica_files(), and update_replica_files_registry().

std::list<size_t> colvarbias_meta::replica_files_pos [protected]
 

Positions in replica files (files are reopened at these positions).

Definition at line 182 of file colvarbias_meta.h.

Referenced by read_replica_files(), and update_replica_files_registry().

std::string colvarbias_meta::replica_files_registry [protected]
 

List of hill files from all the replicas.

Definition at line 176 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), register_replica_file(), and update_replica_files_registry().

std::ofstream colvarbias_meta::replica_out_file [protected]
 

File to contain the hills created by this replica in this run.

Definition at line 185 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), create_hill(), update(), and ~colvarbias_meta().

std::string colvarbias_meta::replica_out_file_name [protected]
 

File to contain the hills created by this replica in this run (name).

Definition at line 188 of file colvarbias_meta.h.

Referenced by colvarbias_meta(), and write_restart().

size_t colvarbias_meta::replica_update_freq [protected]
 

Frequency at which output files from other replicas are checked.

Definition at line 173 of file colvarbias_meta.h.

Referenced by colvarbias_meta().

bool colvarbias_meta::shift_fes [protected]
 

Make the free energy surface be larger than or equal to zero.

Definition at line 200 of file colvarbias_meta.h.

bool colvarbias_meta::use_grids [protected]
 

Bin the hills on grids of energy and forces, and use them to force the colvars (as opposed to deriving the hills analytically).

Definition at line 103 of file colvarbias_meta.h.

Referenced by colvarbias_meta().


The documentation for this class was generated from the following files:
Generated on Mon Nov 23 04:59:33 2009 for NAMD by  doxygen 1.3.9.1