colvarbias_meta.h

Go to the documentation of this file.

#ifndef COLVARBIAS_META_H
#define COLVARBIAS_META_H

#include <vector>
#include <list>
#include <sstream>
#include <fstream>

#include "colvarbias.h"
#include "colvargrid.h"

class colvarbias_meta : public colvarbias {

public:

  enum Communication {
    single_replica,
    multiple_replicas
  };

  Communication comm;

  colvarbias_meta (std::string const &conf, char const *key);

  virtual ~colvarbias_meta();
  
  virtual void update();

  virtual void analyse();

  virtual std::istream & read_restart (std::istream &is);

  virtual std::ostream & write_restart (std::ostream &os);

  class hill;
  typedef std::list<hill>::iterator hill_iter;

protected:

//   virtual void parse_analysis (std::string const &conf);

  std::list<hill> hills;

  hill_iter new_hills_begin;

  std::list<hill> hills_off_grid;

  hill_iter new_hills_off_grid_begin;

  void recount_hills_off_grid (hill_iter h_first, hill_iter h_last,
                               colvar_grid_scalar *ge);

  std::istream & read_hill (std::istream &is);

  virtual std::list<hill>::const_iterator create_hill (hill const &h);

  virtual std::list<hill>::const_iterator delete_hill (std::list<hill>::iterator &h);

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

  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));

  cvm::real  hill_weight;

  bool       use_grids;

  bool       rebin_grids;

  bool       expand_grids;

  size_t     grids_freq;

  bool       keep_hills;

  bool       dump_fes;

  bool       dump_fes_save;

  colvar_grid_scalar    *hills_energy;

  colvar_grid_gradient  *hills_energy_gradients;

  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);



  cvm::real  hill_width;

  size_t     new_hill_freq;

  bool           b_hills_traj;
  std::ofstream  hills_traj_os;


  std::string    replica;

  virtual void register_replica_file (std::string const &new_file);

  virtual void update_replica_files_registry();

  virtual void read_replica_files();

  size_t                 replica_update_freq;

  std::string            replica_files_registry;

  std::list<std::string> replica_files;

  std::list<size_t>      replica_files_pos;

  std::ofstream          replica_out_file;

  std::string            replica_out_file_name;


  // Analysis

  size_t                 free_energy_begin;

  size_t                 free_energy_end;

  bool                   shift_fes;

  cvm::real              free_energy_offset;

  std::string            free_energy_file;

  std::string            free_energy_gradients_file;

  std::string            boltzmann_weights_file;

  std::string            boltzmann_counts_file;

  // weight = scale * exp (-(fe-offset)/(kB*temp))

  cvm::real              boltzmann_weights_scale;

  cvm::real              boltzmann_weights_temp;

  colvar_grid_scalar    *free_energy;

  colvar_grid_gradient  *free_energy_gradients;

  colvar_grid_scalar    *boltzmann_weights;

  colvar_grid_count     *boltzmann_counts;

};




class colvarbias_meta::hill {

protected:

  cvm::real hill_value;

  cvm::real sW;

  cvm::real W;

  std::vector<colvarvalue>  centers;

  std::vector<cvm::real>    widths;

public:

  friend class colvarbias_meta;

  size_t      it;

  std::string replica;

  inline hill (cvm::real             const &W_in,
               std::vector<colvar *>       &cv,
               cvm::real             const &hill_width,
               std::string           const &replica_in = "")
    : sW (1.0),
      W (W_in),
      centers (cv.size()),
      widths (cv.size()),
      it (cvm::it),
      replica (replica_in)
  {
    for (size_t i = 0; i < cv.size(); i++) {
      centers[i].type (cv[i]->type());
      centers[i] = cv[i]->value();
      widths[i] = cv[i]->width * hill_width;
    }
    if (cvm::debug()) 
      cvm::log ("New hill, applied to "+cvm::to_str (cv.size())+
                " collective variables, with reference values "+
                cvm::to_str (centers)+", widths "+
                cvm::to_str (widths)+" and weight "+
                cvm::to_str (W)+".\n");
  }

  inline 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 = "")
    : sW (1.0),
      W (W_in),
      centers (centers_in),
      widths (widths_in),
      it (it_in),
      replica (replica_in)
  {}

  inline hill (colvarbias_meta::hill const &h)
    : sW (1.0),
      W (h.W),
      centers (h.centers),
      widths (h.widths),
      it (h.it),
      replica (h.replica)
  {}

  inline ~hill()
  {}
  
  inline cvm::real energy()
  {
    return W * sW * hill_value;
  }

  inline cvm::real energy (cvm::real const &new_weight)
  {
    return new_weight * sW * hill_value;
  }

  inline cvm::real const &value()
  {
    return hill_value;
  }

  inline void value (cvm::real const &new_value)
  {
    hill_value = new_value;
  }

  inline cvm::real weight()
  {
    return W * sW;
  }

  inline void scale (cvm::real const &new_scale_fac)
  {
    sW = new_scale_fac;
  }

  inline std::vector<colvarvalue> & center()
  {
    return centers;
  }

  inline colvarvalue & center (size_t const &i)
  {
    return centers[i];
  }

  inline friend bool operator < (hill const &h1, hill const &h2)
  {
    if (h1.it < h2.it) return true;
    else return false;
  }

  inline friend bool operator <= (hill const &h1, hill const &h2)
  {
    if (h1.it <= h2.it) return true;
    else return false;
  }

  inline friend bool operator > (hill const &h1, hill const &h2)
  {
    if (h1.it > h2.it) return true;
    else return false;
  }

  inline friend bool operator >= (hill const &h1, hill const &h2)
  {
    if (h1.it >= h2.it) return true;
    else return false;
  }

  inline friend bool operator == (hill const &h1, hill const &h2)
  {
    if ( (h1.it >= h2.it) && (h1.replica == h2.replica) ) return true;
    else return false;
  }

  std::string output_traj();

  inline friend std::ostream & operator << (std::ostream &os,
                                            hill const &h);

};


#endif


// Emacs
// Local Variables:
// mode: C++
// End:

---