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

colvar Class Reference

A collective variable (main class); to be defined, it needs at least one object of a derived class of colvar::cvc; it calculates and returns a object. More...

#include <colvar.h>

Inheritance diagram for colvar:

colvarparse List of all members.

Public Types

enum  task {
  task_gradients, task_collect_gradients, task_fdiff_velocity, task_system_force,
  task_extended_lagrangian, task_langevin, task_output_energy, task_Jacobian_force,
  task_report_Jacobian_force, task_output_value, task_output_velocity, task_output_applied_force,
  task_output_system_force, task_lower_boundary, task_upper_boundary, task_grid,
  task_lower_wall, task_upper_wall, task_runave, task_corrfunc,
  task_ntot
}
 Options controlling the behaviour of the colvar during the simulation, which are set from outside the cvcs. More...

Public Member Functions

colvarvalue::Type type () const
 Type of value.
colvarvalue const & value () const
 Current value (previously obtained from calc() or read_traj()).
colvarvalue const & velocity () const
 Current velocity (previously obtained from calc() or read_traj()).
colvarvalue const & system_force () const
 Current system force (previously obtained from calc() or read_traj()). Note: this is calculated using the atomic forces from the last simulation step.
bool periodic_boundaries () const
 Is the interval defined by the two boundaries periodic?
bool periodic_boundaries (colvarvalue const &lb, colvarvalue const &ub) const
 Is the interval defined by the two boundaries periodic?
 colvar (std::string const &conf)
 Constructor.
void enable (colvar::task const &t)
 Enable the specified task.
void disable (colvar::task const &t)
 Disable the specified task.
 ~colvar ()
 Destructor.
void calc ()
 Calculate the colvar value and all the other requested quantities.
void calc_value (colvarvalue &xn)
 Calculate just the value, and store it in the argument.
void reset_bias_force ()
 Set the total biasing force to zero.
void add_bias_force (colvarvalue const &force)
 Add to the total force from biases.
cvm::real update ()
 Collect all forces on this colvar, integrate internal equations of motion of internal degrees of freedom; see also colvar::communicate_forces() return colvar energy if extended Lagrandian active.
void communicate_forces ()
 Communicate forces (previously calculated in colvar::update()) to the external degrees of freedom.
cvm::real dist2 (colvarvalue const &x1, colvarvalue const &x2) const
 Use the internal metrics (as from objects) to calculate square distances and gradients.
colvarvalue dist2_lgrad (colvarvalue const &x1, colvarvalue const &x2) const
 Use the internal metrics (as from objects) to calculate square distances and gradients.
colvarvalue dist2_rgrad (colvarvalue const &x1, colvarvalue const &x2) const
 Use the internal metrics (as from objects) to calculate square distances and gradients.
cvm::real compare (colvarvalue const &x1, colvarvalue const &x2) const
 Use the internal metrics (as from objects) to compare colvar values.
void wrap (colvarvalue &x) const
 Use the internal metrics (as from objects) to wrap a value into a standard interval.
void parse_analysis (std::string const &conf)
 Read the analysis tasks.
void analyse ()
 Perform analysis tasks.
std::istream & read_traj (std::istream &is)
 Read the value from a collective variable trajectory file.
std::ostream & write_traj (std::ostream &os)
 Output formatted values to the trajectory file.
std::ostream & write_traj_label (std::ostream &os)
 Write a label to the trajectory file (comment line).
std::istream & read_restart (std::istream &is)
 Read the collective variable from a restart file.
std::ostream & write_restart (std::ostream &os)
 Write the collective variable to a restart file.
size_t n_components () const

Public Attributes

std::string name
 Name.
cvm::real width
 Typical fluctuation amplitude for this collective variable (e.g. local width of a free energy basin).
bool b_linear
 True if this is a linear combination of elements.
bool b_inverse_gradients
 True if all objects are capable of calculating inverse gradients.
bool b_Jacobian_force
 True if all objects are capable of calculating Jacobian forces.
bool tasks [task_ntot]
 Tasks performed by this colvar.
colvarvalue fb
 Bias force; reset_bias_force() should be called before the biases are updated.
colvarvalue f
 Total applied force; fr (if task_extended_lagrangian is defined), fb (if biases are applied) and the walls' forces (if defined) contribute to it.
colvarvalue ft
 Total force, as derived from the atomic trajectory; should equal the total system force plus .
cvm::real period
 Period, if it is a constant.
bool b_periodic
 Same as above, but also takes into account components with a variable period, such as distanceZ.
bool expand_boundaries
 Expand the boundaries of multiples of width, to keep the value always within range.
colvarvalue lower_boundary
 Location of the lower boundary.
colvarvalue lower_wall
 Location of the lower wall.
cvm::real lower_wall_k
 Force constant for the lower boundary potential (|x-xb|^2).
colvarvalue upper_boundary
 Location of the upper boundary.
colvarvalue upper_wall
 Location of the upper wall.
cvm::real upper_wall_k
 Force constant for the upper boundary potential (|x-xb|^2).
std::vector< int > atom_ids
 Sorted array of (zero-based) IDs for all atoms involved.
std::vector< cvm::rvectoratomic_gradients
 Array of atomic gradients collected from all cvcs with appropriate components, rotations etc. For scalar variables only!

Protected Types

enum  acf_type_e { acf_notset, acf_vel, acf_coor, acf_p2coor }
 Type of autocorrelation function (ACF). More...

Protected Member Functions

colvarvalue fdiff_velocity (colvarvalue const &xold, colvarvalue const &xnew)
void calc_vel_acf (std::list< colvarvalue > &v_history, colvarvalue const &v)
 Velocity ACF, scalar product between v(0) and v(t).
void calc_coor_acf (std::list< colvarvalue > &x_history, colvarvalue const &x)
 Coordinate ACF, scalar product between x(0) and x(t) (does not work with scalar numbers).
void calc_p2coor_acf (std::list< colvarvalue > &x_history, colvarvalue const &x)
 Coordinate ACF, second order Legendre polynomial between x(0) and x(t) (does not work with scalar numbers).
void calc_acf ()
 Calculate the auto-correlation function (ACF).
void write_acf (std::ostream &os)
 Save the ACF to a file.
void calc_runave ()
 Calculate the running average and its standard deviation.
void build_atom_list (void)
 Initialize the sorted list of atom IDs for atoms involved in all cvcs (called when enabling task_collect_gradients).

Protected Attributes

colvarvalue x
 Value of the colvar.
colvarvalue x_reported
 Cached reported value (x may be manipulated).
colvarvalue v_fdiff
 Finite-difference velocity.
colvarvalue v_reported
 Cached reported velocity.
colvarvalue xr
 Restraint center.
colvarvalue vr
 Velocity of the restraint center.
cvm::real ext_mass
 Mass of the restraint center.
cvm::real ext_force_k
 Restraint force constant.
cvm::real ext_gamma
 Friction coefficient for Langevin extended dynamics.
cvm::real ext_sigma
 Amplitude of Gaussian white noise for Langevin extended dynamics.
colvarvalue fr
 Harmonic restraint force.
colvarvalue fj
 Jacobian force, when task_Jacobian_force is enabled.
colvarvalue ft_reported
 Cached reported system force.
colvarvalue x_old
 Previous value (to calculate velocities during analysis).
std::list< std::list< colvarvalue > > acf_x_history
std::list< std::list< colvarvalue > > acf_v_history
std::list< std::list< colvarvalue
> >::iterator 
acf_x_history_p
std::list< std::list< colvarvalue
> >::iterator 
acf_v_history_p
std::list< std::list< colvarvalue > > x_history
 Time series of values and velocities used in running averages.
std::list< std::list< colvarvalue
> >::iterator 
x_history_p
std::string acf_colvar_name
 Collective variable with which the correlation is calculated (default: itself).
size_t acf_length
 Length of autocorrelation function (ACF).
size_t acf_offset
 After how many steps the ACF starts.
size_t acf_stride
 How many timesteps separate two ACF values.
size_t acf_nframes
 Number of frames for each ACF point.
bool acf_normalize
 Normalize the ACF to a maximum value of 1?
std::vector< cvm::realacf
 ACF values.
std::string acf_outfile
 Name of the file to write the ACF.
acf_type_e acf_type
 Type of autocorrelation function (ACF).
size_t runave_length
 Length of running average series.
size_t runave_stride
 Timesteps to skip between two values in the running average series.
std::ofstream runave_os
 Name of the file to write the running average.
colvarvalue runave
 Current value of the running average.
cvm::real runave_variance
 Current value of the square deviation from the running average.
cvm::real kinetic_energy
 If extended Lagrangian active: colvar energies (kinetic and harmonic potential).
cvm::real potential_energy
std::vector< cvc * > cvcs
 Array of objects.

Detailed Description

A collective variable (main class); to be defined, it needs at least one object of a derived class of colvar::cvc; it calculates and returns a object.

This class parses the configuration, defines the behaviour and stores the value () and all related data of a collective variable. How the value is calculated is defined in and its derived classes. The object contains pointers to multiple derived objects, which can be combined together into one collective variable. This makes possible to implement new collective variables at runtime based on the existing ones. Currently, this possibility is limited to a polynomial, using the coefficients cvc::sup_coeff and the exponents cvc::sup_np. In case of non-scalar variables, only exponents equal to 1 are accepted.

Please note that most of its members are objects, i.e. they can handle different data types together, and must all be set to the same type of colvar::x by using the colvarvalue::type() member function before using them together in assignments or other operations; this is usually done automatically in the constructor. If you add a new member of type, you should also add its initialization line in the constructor.

Definition at line 41 of file colvar.h.


Member Enumeration Documentation

enum colvar::acf_type_e [protected]
 

Type of autocorrelation function (ACF).

Enumeration values:
acf_notset  Unset type.
acf_vel  Velocity ACF, scalar product between v(0) and v(t).
acf_coor  Coordinate ACF, scalar product between x(0) and x(t).
acf_p2coor  Coordinate ACF, second order Legendre polynomial between x(0) and x(t) (does not work with scalar numbers).

Definition at line 404 of file colvar.h.

00404                   {
00406     acf_notset,
00408     acf_vel,
00410     acf_coor,
00413     acf_p2coor
00414   };

enum colvar::task
 

Options controlling the behaviour of the colvar during the simulation, which are set from outside the cvcs.

Enumeration values:
task_gradients  Gradients are calculated and temporarily stored, so that external forces can be applied.
task_collect_gradients  Collect atomic gradient data from all cvcs into vector atomic_gradients.
task_fdiff_velocity  Calculate the velocity with finite differences.
task_system_force  The system force is calculated, projecting the atomic forces on the inverse gradients.
task_extended_lagrangian  The variable has a harmonic restraint around a moving center with fictitious mass; bias forces will be applied to the center.
task_langevin  The extended system coordinate undergoes Langevin dynamics.
task_output_energy  Output the potential and kinetic energies (for extended Lagrangian colvars only).
task_Jacobian_force  Compute analytically the "force" arising from the geometric entropy component (for example, that from the angular states orthogonal to a distance vector).
task_report_Jacobian_force  Report the Jacobian force as part of the system force if disabled, apply a correction internally to cancel it.
task_output_value  Output the value to the trajectory file (on by default).
task_output_velocity  Output the velocity to the trajectory file.
task_output_applied_force  Output the applied force to the trajectory file.
task_output_system_force  Output the system force to the trajectory file.
task_lower_boundary  A lower boundary is defined.
task_upper_boundary  An upper boundary is defined.
task_grid  Provide a discretization of the values of the colvar to be used by the biases or in analysis (needs lower and upper boundary).
task_lower_wall  Apply a restraining potential (|x-xb|^2) to the colvar when it goes below the lower wall.
task_upper_wall  Apply a restraining potential (|x-xb|^2) to the colvar when it goes above the upper wall.
task_runave  Compute running average.
task_corrfunc  Compute time correlation function.
task_ntot  Number of possible tasks.

Definition at line 93 of file colvar.h.

00093             {
00096     task_gradients,
00099     task_collect_gradients,
00101     task_fdiff_velocity,
00104     task_system_force,
00108     task_extended_lagrangian,
00111     task_langevin,
00114     task_output_energy,
00118     task_Jacobian_force,
00121     task_report_Jacobian_force,
00123     task_output_value,
00125     task_output_velocity,
00127     task_output_applied_force,
00129     task_output_system_force,
00131     task_lower_boundary,
00133     task_upper_boundary,
00137     task_grid,
00140     task_lower_wall,
00143     task_upper_wall,
00145     task_runave,
00147     task_corrfunc,
00149     task_ntot
00150   };


Constructor & Destructor Documentation

colvar::colvar std::string const &  conf  ) 
 

Constructor.

Definition at line 14 of file colvar.C.

References b_inverse_gradients, b_Jacobian_force, b_linear, b_periodic, colvarmodule::boltzmann(), cvcs, colvarmodule::debug(), disable(), dist2(), enable(), expand_boundaries, ext_force_k, ext_gamma, ext_mass, ext_sigma, colvarmodule::fatal_error(), fr, initialize_components, j, kinetic_energy, colvarmodule::log(), lower_boundary, lower_wall, lower_wall_k, parse_analysis(), period, periodic_boundaries(), PI, potential_energy, task_extended_lagrangian, task_langevin, task_lower_boundary, task_lower_wall, task_output_applied_force, task_output_energy, task_output_system_force, task_output_value, task_output_velocity, task_upper_boundary, task_upper_wall, tasks, colvarvalue::type(), type(), upper_boundary, upper_wall, upper_wall_k, vr, width, x, x_reported, and xr.

00015 {
00016   cvm::log ("Initializing a new collective variable.\n");
00017   
00018   get_keyval (conf, "name", this->name,
00019               (std::string ("colvar")+cvm::to_str (cvm::colvars.size()+1)));
00020 
00021   for (std::vector<colvar *>::iterator cvi = cvm::colvars.begin();
00022        cvi < cvm::colvars.end();
00023        cvi++) {
00024     if ((*cvi)->name == this->name)
00025       cvm::fatal_error ("Error: this colvar has the same name, \""+this->name+
00026                         "\", as another colvar.\n");
00027   }
00028 
00029   // all tasks disabled by default 
00030   for (size_t i = 0; i < task_ntot; i++) {
00031     tasks[i] = false;
00032   }
00033 
00034   kinetic_energy = 0.0;
00035   potential_energy = 0.0;
00036 
00037   // read the configuration and set up corresponding instances, for
00038   // each type of component implemented
00039 #define initialize_components(def_desc,def_config_key,def_class_name)   \
00040   {                                                                     \
00041     size_t def_count = 0;                                               \
00042     std::string def_conf = "";                                          \
00043     size_t pos = 0;                                                     \
00044     while ( this->key_lookup (conf,                                     \
00045                               def_config_key,                           \
00046                               def_conf,                                 \
00047                               pos) ) {                                  \
00048       if (!def_conf.size()) continue;                                   \
00049       cvm::log ("Initializing "                                         \
00050                 "a new \""+std::string (def_config_key)+"\" component"+ \
00051                 (cvm::debug() ? ", with configuration:\n"+def_conf      \
00052                  : ".\n"));                                             \
00053       cvm::increase_depth();                                            \
00054       cvc *cvcp = new colvar::def_class_name (def_conf);                \
00055       if (cvcp != NULL) {                                               \
00056         cvcs.push_back (cvcp);                                          \
00057         cvcp->check_keywords (def_conf, def_config_key);                \
00058         cvm::decrease_depth();                                          \
00059       } else {                                                          \
00060         cvm::fatal_error ("Error: in allocating component \""           \
00061                           def_config_key"\".\n");                       \
00062       }                                                                 \
00063       if ( (cvcp->period != 0.0) || (cvcp->wrap_center != 0.0) ) {      \
00064         if ( (cvcp->function_type != std::string ("distance_z")) &&     \
00065              (cvcp->function_type != std::string ("dihedral")) &&       \
00066              (cvcp->function_type != std::string ("spin_angle")) ) {    \
00067           cvm::fatal_error ("Error: invalid use of period and/or "      \
00068                             "wrapAround in a \""+                       \
00069                             std::string (def_config_key)+               \
00070                             "\" component.\n"+                          \
00071                             "Period: "+cvm::to_str(cvcp->period) +      \
00072                         " wrapAround: "+cvm::to_str(cvcp->wrap_center));\
00073         }                                                               \
00074       }                                                                 \
00075       if ( ! cvcs.back()->name.size())                                  \
00076         cvcs.back()->name = std::string (def_config_key)+               \
00077           (cvm::to_str (++def_count));                                  \
00078       if (cvm::debug())                                                 \
00079         cvm::log ("Done initializing a \""+                             \
00080                   std::string (def_config_key)+                         \
00081                   "\" component"+                                       \
00082                   (cvm::debug() ?                                       \
00083                    ", named \""+cvcs.back()->name+"\""                  \
00084                    : "")+".\n");                                        \
00085       def_conf = "";                                                    \
00086       if (cvm::debug())                                                 \
00087         cvm::log ("Parsed "+cvm::to_str (cvcs.size())+                  \
00088                   " components at this time.\n");                       \
00089     }                                                                   \
00090   }
00091 
00092 
00093   initialize_components ("distance",         "distance",       distance);
00094   initialize_components ("distance vector",  "distanceVec",    distance_vec);
00095   initialize_components ("distance vector "
00096                          "direction",        "distanceDir",    distance_dir);
00097   initialize_components ("distance projection "
00098                          "on an axis",       "distanceZ",      distance_z);
00099   initialize_components ("distance projection "
00100                          "on a plane",       "distanceXY",     distance_xy);
00101   initialize_components ("minimum distance", "minDistance",    min_distance);
00102 
00103   initialize_components ("coordination "
00104                          "number",           "coordNum",       coordnum);
00105   initialize_components ("self-coordination "
00106                          "number",           "selfCoordNum",   selfcoordnum);
00107 
00108   initialize_components ("angle",            "angle",          angle);
00109   initialize_components ("dihedral",         "dihedral",       dihedral);
00110 
00111   initialize_components ("hydrogen bond",    "hBond",          h_bond);
00112 
00113   //  initialize_components ("alpha helix",      "alphaDihedrals", alpha_dihedrals);
00114   initialize_components ("alpha helix",      "alpha",          alpha_angles);
00115 
00116   initialize_components ("dihedral principal "
00117                          "component",        "dihedralPC",     dihedPC);
00118 
00119   initialize_components ("orientation",      "orientation",    orientation);
00120   initialize_components ("orientation "
00121                          "angle",            "orientationAngle",orientation_angle);
00122   initialize_components ("tilt",             "tilt",           tilt);
00123   initialize_components ("spin angle",       "spinAngle",      spin_angle);
00124 
00125   initialize_components ("RMSD",             "rmsd",           rmsd);
00126 
00127   //  initialize_components ("logarithm of MSD", "logmsd",         logmsd);
00128 
00129   initialize_components ("radius of "
00130                          "gyration",         "gyration",       gyration);
00131   initialize_components ("eigenvector",      "eigenvector",    eigenvector);
00132 
00133   if (!cvcs.size())
00134     cvm::fatal_error ("Error: no valid components were provided "
00135                       "for this collective variable.\n");
00136 
00137   cvm::log ("All components initialized.\n");
00138 
00139 
00140   // this is set false if any of the components has an exponent
00141   // different from 1 in the polynomial
00142   b_linear = true;
00143 
00144   // these will be set to false if any of the cvcs has them false
00145   b_inverse_gradients = true;
00146   b_Jacobian_force    = true;
00147 
00148   // Decide whether the colvar is periodic
00149   // Used to wrap extended DOF if extendedLagrangian is on
00150   if (cvcs.size() == 1 && (cvcs[0])->b_periodic && (cvcs[0])->sup_np == 1 
00151                                                 && (cvcs[0])->sup_coeff == 1.0 ) {
00152     this->b_periodic = true;
00153     this->period = (cvcs[0])->period;
00154     // TODO write explicit wrap() function for colvars to allow for
00155     // sup_coeff different from 1
00156     // this->period = (cvcs[0])->period * (cvcs[0])->sup_coeff;
00157   } else {
00158     this->b_periodic = false;
00159     this->period = 0.0;
00160   }
00161 
00162   // check the available features of each cvc
00163   for (size_t i = 0; i < cvcs.size(); i++) {
00164 
00165     if ((cvcs[i])->sup_np != 1) {
00166       if (cvm::debug() && b_linear)
00167         cvm::log ("Warning: You are using a non-linear polynomial "
00168                   "combination to define this collective variable, "
00169                   "some biasing methods may be unavailable.\n");
00170       b_linear = false;
00171 
00172       if ((cvcs[i])->sup_np < 0) {
00173         cvm::log ("Warning: you chose a negative exponent in the combination; "
00174                   "if you apply forces, the simulation may become unstable "
00175                   "when the component \""+
00176                   (cvcs[i])->function_type+"\" approaches zero.\n");
00177       }
00178     }
00179 
00180     if ((cvcs[i])->b_periodic && !b_periodic) {
00181         cvm::log ("Warning: although this component is periodic, the colvar will "
00182                   "not be treated as periodic, either because the exponent is not "
00183                   "1, or because multiple components are present. Make sure that "
00184                   "you know what you are doing!");
00185     }
00186 
00187     if (! (cvcs[i])->b_inverse_gradients)
00188       b_inverse_gradients = false;
00189 
00190     if (! (cvcs[i])->b_Jacobian_derivative)
00191       b_Jacobian_force = false;
00192 
00193     for (size_t j = i; j < cvcs.size(); j++) {
00194       if ( (cvcs[i])->type() != (cvcs[j])->type() ) {
00195         cvm::fatal_error ("ERROR: you are definining this collective variable "
00196                           "by using components of different types, \""+
00197                           colvarvalue::type_desc[(cvcs[i])->type()]+
00198                           "\" and \""+
00199                           colvarvalue::type_desc[(cvcs[j])->type()]+
00200                           "\". "
00201                           "You must use the same type in order to "
00202                           " sum them together.\n");
00203       }
00204     }
00205   }
00206 
00207 
00208   {
00209     colvarvalue::Type const value_type = (cvcs[0])->type();
00210     if (cvm::debug())
00211       cvm::log ("This collective variable is a "+
00212                 colvarvalue::type_desc[value_type]+", corresponding to "+
00213                 cvm::to_str (colvarvalue::dof_num[value_type])+
00214                 " internal degrees of freedom.\n");
00215     x.type (value_type);
00216     x_reported.type (value_type);
00217   }
00218 
00219   get_keyval (conf, "width", width, 1.0);
00220   if (width <= 0.0)
00221     cvm::fatal_error ("Error: \"width\" must be positive.\n");
00222 
00223   lower_boundary.type (this->type());
00224   lower_wall.type     (this->type());
00225 
00226   upper_boundary.type (this->type());
00227   upper_wall.type     (this->type());
00228 
00229   if (this->type() == colvarvalue::type_scalar) {
00230 
00231     if (get_keyval (conf, "lowerBoundary", lower_boundary, lower_boundary)) {
00232       enable (task_lower_boundary);
00233     }
00234 
00235     get_keyval (conf, "lowerWallConstant", lower_wall_k, 0.0);
00236     if (lower_wall_k > 0.0) {
00237       get_keyval (conf, "lowerWall", lower_wall, lower_boundary);
00238       enable (task_lower_wall);
00239     }
00240 
00241     if (get_keyval (conf, "upperBoundary", upper_boundary, upper_boundary)) {
00242       enable (task_upper_boundary);
00243     }
00244 
00245     get_keyval (conf, "upperWallConstant", upper_wall_k, 0.0);
00246     if (upper_wall_k > 0.0) {
00247       get_keyval (conf, "upperWall", upper_wall, upper_boundary);
00248       enable (task_upper_wall);
00249     }
00250   }
00251 
00252   // consistency checks for boundaries and walls
00253   if (tasks[task_lower_boundary] && tasks[task_upper_boundary]) {
00254     if (lower_boundary >= upper_boundary) {
00255       cvm::fatal_error ("Error: the upper boundary, "+
00256                         cvm::to_str (upper_boundary)+
00257                         ", is not higher than the lower boundary, "+
00258                         cvm::to_str (lower_boundary)+".\n");
00259     }
00260   }
00261 
00262   if (tasks[task_lower_wall] && tasks[task_upper_wall]) {
00263     if (lower_wall >= upper_wall) {
00264       cvm::fatal_error ("Error: the upper wall, "+
00265                         cvm::to_str (upper_wall)+
00266                         ", is not higher than the lower wall, "+
00267                         cvm::to_str (lower_wall)+".\n");
00268     }
00269 
00270     if (dist2 (lower_wall, upper_wall) < 1.0E-12) {
00271       cvm::log ("Lower wall and upper wall are equal "
00272                 "in the periodic domain of the colvar: disabling walls.\n");
00273       disable (task_lower_wall);
00274       disable (task_upper_wall);
00275     }
00276   }
00277 
00278   get_keyval (conf, "expandBoundaries", expand_boundaries, false);
00279   if (expand_boundaries && periodic_boundaries()) {
00280     cvm::fatal_error ("Error: trying to expand boundaries that already "
00281                       "cover a whole period of a periodic colvar.\n");
00282   }
00283 
00284   {
00285     bool b_extended_lagrangian;
00286     get_keyval (conf, "extendedLagrangian", b_extended_lagrangian, false);
00287 
00288     if (b_extended_lagrangian) {
00289       cvm::real temp, tolerance, period;
00290 
00291       cvm::log ("Enabling the extended Lagrangian term for colvar \""+
00292                 this->name+"\".\n");
00293     
00294       enable (task_extended_lagrangian);
00295 
00296       xr.type (this->type());
00297       vr.type (this->type());
00298       fr.type (this->type());
00299 
00300       const bool found = get_keyval (conf, "extendedTemp", temp, cvm::temperature());
00301       if (temp <= 0.0) {
00302         if (found)
00303           cvm::fatal_error ("Error: \"extendedTemp\" must be positive.\n");
00304         else
00305           cvm::fatal_error ("Error: a positive temperature must be provided, either "
00306                             "by enabling a thermostat, or through \"extendedTemp\".\n");
00307       }
00308 
00309       get_keyval (conf, "extendedFluctuation", tolerance, 0.2*width);
00310       if (tolerance <= 0.0)
00311         cvm::fatal_error ("Error: \"extendedFluctuation\" must be positive.\n");
00312       ext_force_k = cvm::boltzmann() * temp / (tolerance * tolerance);
00313       cvm::log ("Computed extended system force constant: " + cvm::to_str(ext_force_k) + " kcal/mol/U^2");
00314 
00315       get_keyval (conf, "extendedTimeConstant", period, 40.0 * cvm::dt());
00316       if (period <= 0.0)
00317         cvm::fatal_error ("Error: \"extendedTimeConstant\" must be positive.\n");
00318       ext_mass = (cvm::boltzmann() * temp * period * period)
00319                  / (4.0 * PI * PI * tolerance * tolerance);
00320       cvm::log ("Computed fictitious mass: " + cvm::to_str(ext_mass) + " kcal/mol/(U/fs)^2   (U: colvar unit)");
00321 
00322       {
00323         bool b_output_energy;
00324         get_keyval (conf, "outputEnergy", b_output_energy, false);
00325         if (b_output_energy) {
00326           enable (task_output_energy);
00327         }
00328       }
00329 
00330       get_keyval (conf, "extendedLangevinDamping", ext_gamma, 0.0);
00331       if (ext_gamma < 0.0)
00332         cvm::fatal_error ("Error: \"extendedLangevinDamping\" may not be negative.\n");
00333       if (ext_gamma != 0.0) {
00334         enable (task_langevin);
00335         ext_gamma *= 1.0e-3; // convert from ps-1 to fs-1
00336         ext_sigma = std::sqrt(2.0 * cvm::boltzmann() * temp * ext_gamma * ext_mass / cvm::dt());
00337       }
00338     }
00339   }
00340 
00341   {
00342     bool b_output_value;
00343     get_keyval (conf, "outputValue", b_output_value, true);
00344     if (b_output_value) {
00345       enable (task_output_value);
00346     }
00347   }
00348 
00349   {
00350     bool b_output_velocity;
00351     get_keyval (conf, "outputVelocity", b_output_velocity, false);
00352     if (b_output_velocity) {
00353       enable (task_output_velocity);
00354     }
00355   }
00356 
00357   {
00358     bool b_output_system_force;
00359     get_keyval (conf, "outputSystemForce", b_output_system_force, false);
00360     if (b_output_system_force) {
00361       enable (task_output_system_force);
00362     }
00363   }
00364 
00365   {
00366     bool b_output_applied_force;
00367     get_keyval (conf, "outputAppliedForce", b_output_applied_force, false);
00368     if (b_output_applied_force) {
00369       enable (task_output_applied_force);
00370     }
00371   }
00372 
00373   if (cvm::b_analysis)
00374     parse_analysis (conf);
00375 
00376   if (cvm::debug())
00377     cvm::log ("Done initializing collective variable \""+this->name+"\".\n");
00378 }

colvar::~colvar  ) 
 

Destructor.

Definition at line 650 of file colvar.C.

References acf, acf_outfile, cvcs, colvarmodule::fatal_error(), colvarmodule::log(), runave_os, and write_acf().

00651 {
00652   if (cvm::b_analysis) {
00653 
00654     if (acf.size()) {
00655       cvm::log ("Writing acf to file \""+acf_outfile+"\".\n");
00656 
00657       std::ofstream acf_os (acf_outfile.c_str());
00658       if (! acf_os.good())
00659         cvm::fatal_error ("Cannot open file \""+acf_outfile+"\".\n");
00660       write_acf (acf_os);
00661       acf_os.close();
00662     }
00663 
00664     if (runave_os.good()) {
00665       runave_os.close();
00666     }
00667   }
00668 
00669   for (size_t i = 0; i < cvcs.size(); i++) {
00670     delete cvcs[i];
00671   }
00672 }  


Member Function Documentation

void colvar::add_bias_force colvarvalue const &  force  )  [inline]
 

Add to the total force from biases.

Definition at line 550 of file colvar.h.

References fb.

Referenced by colvarbias::communicate_forces().

00551 {
00552   fb += force;
00553 }

void colvar::analyse  ) 
 

Perform analysis tasks.

Definition at line 1325 of file colvar.C.

References calc_acf(), calc_runave(), and tasks.

01326 {
01327   if (tasks[task_runave]) {
01328     calc_runave();
01329   }
01330 
01331   if (tasks[task_corrfunc]) {
01332     calc_acf();
01333   }
01334 }

void colvar::build_atom_list void   )  [protected]
 

Initialize the sorted list of atom IDs for atoms involved in all cvcs (called when enabling task_collect_gradients).

Definition at line 382 of file colvar.C.

References atom_ids, atomic_gradients, cvcs, colvarmodule::debug(), j, and colvarmodule::log().

Referenced by enable().

00383 {
00384   // If atomic gradients are requested, build full list of atom ids from all cvcs
00385   std::list<int> temp_id_list;
00386 
00387   for (size_t i = 0; i < cvcs.size(); i++) {
00388     for (size_t j = 0; j < cvcs[i]->atom_groups.size(); j++) {
00389       for (size_t k = 0; k < cvcs[i]->atom_groups[j]->size(); k++) {
00390         temp_id_list.push_back (cvcs[i]->atom_groups[j]->at(k).id);
00391       }
00392     }
00393   }
00394 
00395   temp_id_list.sort();
00396   temp_id_list.unique();
00397   atom_ids = std::vector<int> (temp_id_list.begin(), temp_id_list.end());
00398   temp_id_list.clear();
00399 
00400   atomic_gradients.resize (atom_ids.size());
00401   if (atom_ids.size()) {
00402     if (cvm::debug())
00403       cvm::log ("Colvar: created atom list with " + cvm::to_str(atom_ids.size()) + " atoms.\n");
00404   } else {
00405     cvm::log ("Warning: colvar components communicated no atom IDs.\n");
00406   }
00407 }

void colvar::calc  ) 
 

Calculate the colvar value and all the other requested quantities.

Definition at line 679 of file colvar.C.

References atom_ids, atomic_gradients, cvcs, colvarmodule::debug(), colvarmodule::decrease_depth(), dist2_lgrad(), fdiff_velocity(), ft, ft_reported, colvarmodule::increase_depth(), colvarvalue::inverse(), j, colvarmodule::log(), colvarmodule::real, colvarvalue::real_value, colvarvalue::reset(), colvarmodule::rotation::rotate(), colvarmodule::step_relative(), tasks, colvarvalue::type(), v_fdiff, v_reported, value(), vr, x, x_old, x_reported, and xr.

00680 {
00681   if (cvm::debug())
00682     cvm::log ("Calculating colvar \""+this->name+"\".\n");
00683 
00684   // calculate the value of the colvar
00685 
00686   x.reset();
00687   if (x.type() == colvarvalue::type_scalar) {
00688 
00689     // scalar variable, polynomial combination allowed
00690     for (size_t i = 0; i < cvcs.size(); i++) {
00691       cvm::increase_depth();
00692       (cvcs[i])->calc_value();
00693       cvm::decrease_depth();
00694       if (cvm::debug())
00695         cvm::log ("Colvar component no. "+cvm::to_str (i+1)+
00696                   " within colvar \""+this->name+"\" has value "+
00697                   cvm::to_str ((cvcs[i])->value(),
00698                                cvm::cv_width, cvm::cv_prec)+".\n");
00699       x += (cvcs[i])->sup_coeff *
00700         ( ((cvcs[i])->sup_np != 1) ?
00701           std::pow ((cvcs[i])->value().real_value, (cvcs[i])->sup_np) :
00702           (cvcs[i])->value().real_value );
00703     } 
00704   } else {
00705 
00706     for (size_t i = 0; i < cvcs.size(); i++) {
00707       cvm::increase_depth();
00708       (cvcs[i])->calc_value();
00709       cvm::decrease_depth();
00710       if (cvm::debug())
00711         cvm::log ("Colvar component no. "+cvm::to_str (i+1)+
00712                   " within colvar \""+this->name+"\" has value "+
00713                   cvm::to_str ((cvcs[i])->value(),
00714                                cvm::cv_width, cvm::cv_prec)+".\n");
00715       x += (cvcs[i])->sup_coeff * (cvcs[i])->value();
00716     } 
00717   }
00718 
00719   if (cvm::debug())
00720     cvm::log ("Colvar \""+this->name+"\" has value "+
00721               cvm::to_str (x, cvm::cv_width, cvm::cv_prec)+".\n");
00722 
00723   if (tasks[task_gradients]) {
00724     // calculate the gradients
00725     for (size_t i = 0; i < cvcs.size(); i++) {
00726       cvm::increase_depth();
00727       (cvcs[i])->calc_gradients();
00728       cvm::decrease_depth();
00729     }
00730     if (tasks[task_collect_gradients]) {
00731       // Collect the atomic gradients inside colvar object
00732       for (int a = 0; a < atomic_gradients.size(); a++) {
00733         atomic_gradients[a].reset();
00734       }
00735       for (size_t i = 0; i < cvcs.size(); i++) {
00736         // Coefficient: d(a * x^n) = a * n * x^(n-1) * dx
00737         cvm::real coeff = (cvcs[i])->sup_coeff * cvm::real ((cvcs[i])->sup_np) *
00738           std::pow ((cvcs[i])->value().real_value, (cvcs[i])->sup_np-1);
00739 
00740         for (size_t j = 0; j < cvcs[i]->atom_groups.size(); j++) {
00741 
00742           // If necessary, apply inverse rotation to get atomic
00743           // gradient in the laboratory frame
00744           if (cvcs[i]->atom_groups[j]->b_rotate) {
00745             cvm::rotation const rot_inv = cvcs[i]->atom_groups[j]->rot.inverse();
00746 
00747             for (size_t k = 0; k < cvcs[i]->atom_groups[j]->size(); k++) {
00748               int a = std::lower_bound (atom_ids.begin(), atom_ids.end(),
00749                   cvcs[i]->atom_groups[j]->at(k).id) - atom_ids.begin();
00750               atomic_gradients[a] += coeff *
00751                 rot_inv.rotate (cvcs[i]->atom_groups[j]->at(k).grad);
00752             }
00753 
00754           } else {
00755 
00756             for (size_t k = 0; k < cvcs[i]->atom_groups[j]->size(); k++) {
00757               int a = std::lower_bound (atom_ids.begin(), atom_ids.end(),
00758                   cvcs[i]->atom_groups[j]->at(k).id) - atom_ids.begin();
00759               atomic_gradients[a] += coeff * cvcs[i]->atom_groups[j]->at(k).grad; 
00760             }
00761           }
00762         }
00763       }
00764     }
00765   }
00766 
00767   if (tasks[task_system_force]) {
00768     ft.reset();
00769 
00770     if(!tasks[task_extended_lagrangian] && (cvm::step_relative() > 0)) {
00771       // get from the cvcs the system forces from the PREVIOUS step
00772       for (size_t i = 0; i < cvcs.size(); i++) {
00773         (cvcs[i])->calc_force_invgrads();
00774         // linear combination is assumed
00775         cvm::increase_depth();
00776         ft += (cvcs[i])->system_force() / ((cvcs[i])->sup_coeff * cvm::real (cvcs.size()));
00777         cvm::decrease_depth();
00778       }
00779     }
00780 
00781     if (tasks[task_report_Jacobian_force]) {
00782       // add the Jacobian force to the system force, and don't apply any silent
00783       // correction internally: biases such as colvarbias_abf will handle it
00784       ft += fj;
00785     }
00786   }
00787 
00788   if (tasks[task_fdiff_velocity]) {
00789     // calculate the velocity by finite differences
00790     if (cvm::step_relative() == 0)
00791       x_old = x;
00792     else {
00793       v_fdiff = fdiff_velocity (x_old, x);
00794       v_reported = v_fdiff;
00795     }
00796   }
00797 
00798   if (tasks[task_extended_lagrangian]) {
00799 
00800     // initialize the restraint center in the first step to the value
00801     // just calculated from the cvcs
00802     // TODO: put it in the restart information
00803     if (cvm::step_relative() == 0) {
00804       xr = x;
00805       vr = 0.0; // (already 0; added for clarity)
00806     }
00807 
00808     // report the restraint center as "value"
00809     x_reported = xr;
00810     v_reported = vr;
00811     // the "system force" with the extended Lagrangian is just the
00812     // harmonic term acting on the extended coordinate
00813     // Note: this is the force for current timestep
00814     ft_reported = (-0.5 * ext_force_k) * this->dist2_lgrad (xr, x);
00815 
00816   } else {
00817 
00818     x_reported = x;
00819     ft_reported = ft;
00820   }
00821 
00822   if (cvm::debug())
00823     cvm::log ("Done calculating colvar \""+this->name+"\".\n");
00824 }

void colvar::calc_acf  )  [protected]
 

Calculate the auto-correlation function (ACF).

Definition at line 1354 of file colvar.C.

References acf, acf_colvar_name, acf_coor, acf_length, acf_nframes, acf_offset, acf_p2coor, acf_v_history, acf_v_history_p, acf_vel, acf_x_history, acf_x_history_p, calc_coor_acf(), calc_p2coor_acf(), calc_vel_acf(), colvarmodule::colvar_p(), colvarmodule::fatal_error(), fdiff_velocity(), history_add_value(), history_incr(), colvarmodule::log(), name, tasks, type(), v_fdiff, v_reported, value(), velocity(), and x_old.

Referenced by analyse().

01355 {
01356   // using here an acf_stride-long list of vectors for either
01357   // coordinates (acf_x_history) or velocities (acf_v_history); each vector can
01358   // contain up to acf_length values, which are contiguous in memory
01359   // representation but separated by acf_stride in the time series;
01360   // the pointer to each vector is changed at every step
01361 
01362   if (! (acf_x_history.size() || acf_v_history.size()) ) {
01363 
01364     // first-step operations
01365 
01366     colvar *cfcv = (acf_colvar_name.size() ? 
01367                     cvm::colvar_p (acf_colvar_name) :
01368                     this);
01369     if (cfcv->type() != this->type())
01370       cvm::fatal_error ("Error: correlation function between \""+cfcv->name+
01371                         "\" and \""+this->name+"\" cannot be calculated, "
01372                         "because their value types are different.\n");
01373     acf_nframes = 0;
01374 
01375     cvm::log ("Colvar \""+this->name+"\": initializing ACF calculation.\n");
01376 
01377     if (acf.size() < acf_length+1)
01378       acf.resize (acf_length+1, 0.0);
01379 
01380     switch (acf_type) {
01381 
01382     case acf_vel:
01383       // allocate space for the velocities history
01384       for (size_t i = 0; i < acf_stride; i++) {
01385         acf_v_history.push_back (std::list<colvarvalue>());
01386       }
01387       acf_v_history_p = acf_v_history.begin();
01388       break;
01389 
01390     case acf_coor:
01391     case acf_p2coor:
01392       // allocate space for the coordinates history
01393       for (size_t i = 0; i < acf_stride; i++) {
01394         acf_x_history.push_back (std::list<colvarvalue>());
01395       }
01396       acf_x_history_p = acf_x_history.begin();
01397       break;
01398 
01399     default:
01400       break;
01401     }
01402 
01403   } else {
01404 
01405     colvar *cfcv = (acf_colvar_name.size() ? 
01406                     cvm::colvar_p (acf_colvar_name) :
01407                     this);
01408 
01409     switch (acf_type) {
01410 
01411     case acf_vel:
01412 
01413       if (tasks[task_fdiff_velocity]) {
01414         // calc() should do this already, but this only happens in a
01415         // simulation; better do it again in case a trajectory is
01416         // being read
01417         v_reported = v_fdiff = fdiff_velocity (x_old, cfcv->value());
01418       }
01419 
01420       calc_vel_acf ((*acf_v_history_p), cfcv->velocity());
01421       // store this value in the history
01422       history_add_value (acf_length+acf_offset, *acf_v_history_p, cfcv->velocity());
01423       // if stride is larger than one, cycle among different histories
01424       history_incr (acf_v_history, acf_v_history_p);
01425       break;
01426 
01427     case acf_coor:
01428 
01429       calc_coor_acf ((*acf_x_history_p), cfcv->value());
01430       history_add_value (acf_length+acf_offset, *acf_x_history_p, cfcv->value());
01431       history_incr (acf_x_history, acf_x_history_p);
01432       break;
01433 
01434     case acf_p2coor:
01435 
01436       calc_p2coor_acf ((*acf_x_history_p), cfcv->value());
01437       history_add_value (acf_length+acf_offset, *acf_x_history_p, cfcv->value());
01438       history_incr (acf_x_history, acf_x_history_p);
01439       break;
01440 
01441     default:
01442       break;
01443     }
01444   }
01445 
01446   if (tasks[task_fdiff_velocity]) {
01447     // set it for the next step
01448     x_old = x;
01449   }
01450 }

void colvar::calc_coor_acf std::list< colvarvalue > &  x_history,
colvarvalue const &  x
[protected]
 

Coordinate ACF, scalar product between x(0) and x(t) (does not work with scalar numbers).

Definition at line 1477 of file colvar.C.

References acf, acf_length, acf_nframes, colvarvalue::inner_opt(), and colvarvalue::norm2().

Referenced by calc_acf().

01479 {
01480   // same as above but for coordinates
01481   if (x_list.size() >= acf_length+acf_offset) {
01482     std::list<colvarvalue>::iterator  xs_i = x_list.begin();
01483     std::vector<cvm::real>::iterator acf_i = acf.begin();
01484 
01485     for (size_t i = 0; i < acf_offset; i++)
01486       xs_i++;
01487 
01488     *(acf_i++) += x.norm2();
01489 
01490     colvarvalue::inner_opt (x, xs_i, x_list.end(), acf_i); 
01491 
01492     acf_nframes++;
01493   }
01494 }

void colvar::calc_p2coor_acf std::list< colvarvalue > &  x_history,
colvarvalue const &  x
[protected]
 

Coordinate ACF, second order Legendre polynomial between x(0) and x(t) (does not work with scalar numbers).

Definition at line 1497 of file colvar.C.

References acf, acf_length, acf_nframes, and colvarvalue::p2leg_opt().

Referenced by calc_acf().

01499 {
01500   // same as above but with second order Legendre polynomial instead
01501   // of just the scalar product
01502   if (x_list.size() >= acf_length+acf_offset) {
01503     std::list<colvarvalue>::iterator  xs_i = x_list.begin();
01504     std::vector<cvm::real>::iterator acf_i = acf.begin();
01505 
01506     for (size_t i = 0; i < acf_offset; i++)
01507       xs_i++;
01508 
01509     // value of P2(0) = 1
01510     *(acf_i++) += 1.0;
01511 
01512     colvarvalue::p2leg_opt (x, xs_i, x_list.end(), acf_i); 
01513 
01514     acf_nframes++;
01515   }
01516 }

void colvar::calc_runave  )  [protected]
 

Calculate the running average and its standard deviation.

Definition at line 1546 of file colvar.C.

References acf_nframes, colvarvalue::apply_constraints(), colvarmodule::debug(), dist2(), history_add_value(), colvarmodule::log(), colvarmodule::real, colvarvalue::reset(), runave, runave_length, runave_os, runave_variance, colvarmodule::step_relative(), colvarvalue::type(), x, x_history, and x_history_p.

Referenced by analyse().

01547 {
01548   if (!x_history.size()) {
01549 
01550     runave.type (x.type());
01551     runave.reset();
01552 
01553     // first-step operations
01554 
01555     if (cvm::debug())
01556       cvm::log ("Colvar \""+this->name+
01557                 "\": initializing running average calculation.\n");
01558 
01559     acf_nframes = 0;
01560 
01561     x_history.push_back (std::list<colvarvalue>());
01562     x_history_p = x_history.begin();
01563 
01564   } else {
01565 
01566     if ( (cvm::step_relative() % runave_stride) == 0) {
01567 
01568       if ((*x_history_p).size() >= runave_length-1) {
01569 
01570         runave = x;
01571         for (std::list<colvarvalue>::iterator xs_i = (*x_history_p).begin();
01572              xs_i != (*x_history_p).end(); xs_i++) {
01573           runave += (*xs_i);
01574         }
01575         runave *= 1.0 / cvm::real (runave_length);
01576         runave.apply_constraints();
01577 
01578         runave_variance = 0.0;
01579         runave_variance += this->dist2 (x, runave);
01580         for (std::list<colvarvalue>::iterator xs_i = (*x_history_p).begin();
01581              xs_i != (*x_history_p).end(); xs_i++) {
01582           runave_variance += this->dist2 (x, (*xs_i));
01583         }
01584         runave_variance *= 1.0 / cvm::real (runave_length-1);
01585 
01586         runave_os << std::setw (cvm::it_width) << cvm::step_relative()
01587                   << "  "
01588                   << std::setprecision (cvm::cv_prec) << std::setw (cvm::cv_width)
01589                   << runave << " "
01590                   << std::setprecision (cvm::cv_prec) << std::setw (cvm::cv_width)
01591                   << std::sqrt (runave_variance) << "\n";
01592       }
01593 
01594       history_add_value (runave_length, *x_history_p, x);
01595     }
01596   }
01597 
01598 }

void colvar::calc_value colvarvalue xn  ) 
 

Calculate just the value, and store it in the argument.

void colvar::calc_vel_acf std::list< colvarvalue > &  v_history,
colvarvalue const &  v
[protected]
 

Velocity ACF, scalar product between v(0) and v(t).

Definition at line 1453 of file colvar.C.

References acf, acf_length, acf_nframes, colvarvalue::inner_opt(), and colvarvalue::norm2().

Referenced by calc_acf().

01455 {
01456   // loop over stored velocities and add to the ACF, but only the
01457   // length is sufficient to hold an entire row of ACF values
01458   if (v_list.size() >= acf_length+acf_offset) {
01459     std::list<colvarvalue>::iterator  vs_i = v_list.begin();
01460     std::vector<cvm::real>::iterator acf_i = acf.begin();
01461 
01462     for (size_t i = 0; i < acf_offset; i++)
01463       vs_i++;
01464 
01465     // current vel with itself
01466     *(acf_i++) += v.norm2();
01467 
01468     // inner products of previous velocities with current (acf_i and
01469     // vs_i are updated)
01470     colvarvalue::inner_opt (v, vs_i, v_list.end(), acf_i); 
01471 
01472     acf_nframes++;
01473   }
01474 }

void colvar::communicate_forces  ) 
 

Communicate forces (previously calculated in colvar::update()) to the external degrees of freedom.

Definition at line 937 of file colvar.C.

References cvcs, colvarmodule::debug(), colvarmodule::decrease_depth(), f, colvarmodule::increase_depth(), colvarmodule::log(), colvarvalue::type(), value(), and x.

00938 {
00939   if (cvm::debug())
00940     cvm::log ("Communicating forces from colvar \""+this->name+"\".\n"); 
00941 
00942   if (x.type() == colvarvalue::type_scalar) {
00943 
00944     for (size_t i = 0; i < cvcs.size(); i++) {
00945       cvm::increase_depth();
00946       (cvcs[i])->apply_force (f * (cvcs[i])->sup_coeff * 
00947                               cvm::real ((cvcs[i])->sup_np) *
00948                               (std::pow ((cvcs[i])->value().real_value,
00949                                       (cvcs[i])->sup_np-1)) );
00950       cvm::decrease_depth();
00951     }
00952 
00953   } else {
00954 
00955     for (size_t i = 0; i < cvcs.size(); i++) {
00956       cvm::increase_depth();
00957       (cvcs[i])->apply_force (f * (cvcs[i])->sup_coeff);
00958       cvm::decrease_depth();
00959     }
00960   }
00961 
00962   if (cvm::debug())
00963     cvm::log ("Done communicating forces from colvar \""+this->name+"\".\n"); 
00964 }

cvm::real colvar::compare colvarvalue const &  x1,
colvarvalue const &  x2
const
 

Use the internal metrics (as from objects) to compare colvar values.

Handles correctly symmetries and periodic boundary conditions

Definition at line 1022 of file colvar.C.

References cvcs.

01024 {
01025   return (cvcs[0])->compare (x1, x2);
01026 }

void colvar::disable colvar::task const &  t  ) 
 

Disable the specified task.

Definition at line 602 of file colvar.C.

References task_corrfunc, task_extended_lagrangian, task_fdiff_velocity, task_gradients, task_grid, task_Jacobian_force, task_lower_boundary, task_lower_wall, task_ntot, task_output_applied_force, task_output_system_force, task_output_value, task_output_velocity, task_report_Jacobian_force, task_runave, task_system_force, task_upper_boundary, task_upper_wall, and tasks.

Referenced by colvar().

00603 {
00604   // check dependencies
00605   switch (t) {
00606   case task_gradients:
00607     disable (task_upper_wall);
00608     disable (task_lower_wall);
00609     disable (task_output_applied_force);
00610     disable (task_system_force);
00611     disable (task_Jacobian_force);
00612     break;
00613 
00614   case task_system_force:
00615     disable (task_output_system_force);
00616     break;
00617 
00618   case task_Jacobian_force:
00619     disable (task_report_Jacobian_force);
00620     break;
00621 
00622   case task_fdiff_velocity:
00623     disable (task_output_velocity);
00624     break;
00625 
00626   case task_lower_boundary:
00627   case task_upper_boundary:
00628     disable (task_grid);
00629     break;
00630 
00631   case task_extended_lagrangian:
00632   case task_report_Jacobian_force:
00633   case task_output_value:
00634   case task_output_velocity:
00635   case task_output_applied_force:
00636   case task_output_system_force:
00637   case task_runave:
00638   case task_corrfunc:
00639   case task_grid:
00640   case task_lower_wall:
00641   case task_upper_wall:
00642   case task_ntot:
00643     break;
00644   }
00645 
00646   tasks[t] = false;
00647 }

cvm::real colvar::dist2 colvarvalue const &  x1,
colvarvalue const &  x2
const
 

Use the internal metrics (as from objects) to calculate square distances and gradients.

Handles correctly symmetries and periodic boundary conditions

Definition at line 1004 of file colvar.C.

References cvcs, and colvarvalue::dist2().

Referenced by calc_runave(), colvar(), colvarbias_harmonic::energy_difference(), colvarbias_harmonic::update(), and update().

01006 {
01007   return (cvcs[0])->dist2 (x1, x2);
01008 }

colvarvalue colvar::dist2_lgrad colvarvalue const &  x1,
colvarvalue const &  x2
const
 

Use the internal metrics (as from objects) to calculate square distances and gradients.

Handles correctly symmetries and periodic boundary conditions

Definition at line 1010 of file colvar.C.

References cvcs.

Referenced by calc(), fdiff_velocity(), and update().

01012 {
01013   return (cvcs[0])->dist2_lgrad (x1, x2);
01014 }

colvarvalue colvar::dist2_rgrad colvarvalue const &  x1,
colvarvalue const &  x2
const
 

Use the internal metrics (as from objects) to calculate square distances and gradients.

Handles correctly symmetries and periodic boundary conditions

Definition at line 1016 of file colvar.C.

References cvcs.

Referenced by update().

01018 {
01019   return (cvcs[0])->dist2_rgrad (x1, x2);
01020 }

void colvar::enable colvar::task const &  t  ) 
 

Enable the specified task.

Definition at line 490 of file colvar.C.

References atom_ids, build_atom_list(), colvarmodule::debug(), f, colvarmodule::fatal_error(), fb, fj, ft, ft_reported, colvarmodule::log(), colvarmodule::request_system_force(), task_collect_gradients, task_corrfunc, task_extended_lagrangian, task_fdiff_velocity, task_gradients, task_grid, task_Jacobian_force, task_lower_boundary, task_lower_wall, task_ntot, task_output_applied_force, task_output_system_force, task_output_value, task_output_velocity, task_report_Jacobian_force, task_runave, task_system_force, task_upper_boundary, task_upper_wall, tasks, type(), colvarvalue::type(), v_fdiff, v_reported, and x_old.

Referenced by colvarbias::add_colvar(), colvar(), and parse_analysis().

00491 {
00492   switch (t) {
00493 
00494   case task_output_system_force:
00495     enable (task_system_force);
00496     break;
00497 
00498   case task_report_Jacobian_force:
00499     enable (task_Jacobian_force);
00500     enable (task_system_force);
00501     if (cvm::debug())
00502       cvm::log ("Adding the Jacobian force to the system force, "
00503                 "rather than applying its opposite silently.\n");
00504     break;
00505 
00506   case task_Jacobian_force:
00507     // checks below do not apply to extended-system colvars
00508     if ( !tasks[task_extended_lagrangian] ) {
00509       enable (task_gradients);
00510 
00511       if (!b_Jacobian_force) 
00512         cvm::fatal_error ("Error: colvar \""+this->name+
00513                           "\" does not have Jacobian forces implemented.\n");
00514       if (!b_linear) 
00515         cvm::fatal_error ("Error: colvar \""+this->name+
00516                           "\" must be defined as a linear combination "
00517                           "to calculate the Jacobian force.\n");
00518       if (cvm::debug())
00519         cvm::log ("Enabling calculation of the Jacobian force "
00520                   "on this colvar.\n");
00521     }
00522     fj.type (this->type());
00523     break;
00524 
00525   case task_system_force:
00526     if (!tasks[task_extended_lagrangian]) {
00527       if (!b_inverse_gradients) {
00528         cvm::fatal_error ("Error: one or more of the components of "
00529                           "colvar \""+this->name+
00530                           "\" does not support system force calculation.\n");
00531       }
00532       cvm::request_system_force();
00533     }
00534     ft.type (this->type());
00535     ft_reported.type (this->type());
00536     break;
00537 
00538   case task_output_applied_force:
00539   case task_lower_wall:
00540   case task_upper_wall:
00541     // all of the above require gradients
00542     enable (task_gradients);
00543     break;
00544 
00545   case task_fdiff_velocity:
00546     x_old.type (this->type());
00547     v_fdiff.type (this->type());
00548     v_reported.type (this->type());
00549     break;
00550 
00551   case task_output_velocity:
00552     enable (task_fdiff_velocity);
00553     break;
00554 
00555   case task_grid:
00556     if (this->type() != colvarvalue::type_scalar)
00557       cvm::fatal_error ("Cannot calculate a grid for collective variable, \""+
00558                         this->name+"\", because its value is not a scalar number.\n");
00559     break;
00560 
00561   case task_extended_lagrangian:
00562     enable (task_gradients);
00563     v_reported.type (this->type());
00564     break;
00565 
00566   case task_lower_boundary:
00567   case task_upper_boundary:
00568     if (this->type() != colvarvalue::type_scalar) {
00569       cvm::fatal_error ("Error: this colvar is not a scalar value "
00570                         "and cannot produce a grid.\n");
00571     }
00572     break;
00573 
00574 
00575   case task_output_value:
00576   case task_runave:
00577   case task_corrfunc:
00578   case task_ntot:
00579     break;
00580 
00581   case task_gradients:
00582     f.type  (this->type());
00583     fb.type (this->type());
00584     break;
00585 
00586   case task_collect_gradients:
00587     if (this->type() != colvarvalue::type_scalar)
00588       cvm::fatal_error ("Collecting atomic gradients for non-scalar collective variable \""+
00589                         this->name+"\" is not yet implemented.\n");
00590     enable (task_gradients);
00591     if (atom_ids.size() == 0) {
00592       build_atom_list();
00593     }
00594     break;
00595   }
00596 
00597 
00598   tasks[t] = true;
00599 }

colvarvalue colvar::fdiff_velocity colvarvalue const &  xold,
colvarvalue const &  xnew
[inline, protected]
 

Definition at line 188 of file colvar.h.

References dist2_lgrad(), and colvarmodule::dt().

Referenced by calc(), and calc_acf().

00190   {
00191     // using the gradient of the square distance to calculate the
00192     // velocity (non-scalar variables automatically taken into
00193     // account)
00194     cvm::real dt = cvm::dt();
00195     return ( ( (dt > 0.0) ? (1.0/dt) : 1.0 ) *
00196              0.5 * dist2_lgrad (xnew, xold) );
00197   }

size_t colvar::n_components  )  const [inline]
 

Definition at line 507 of file colvar.h.

References cvcs.

00507                                       {
00508     return cvcs.size();
00509   }

void colvar::parse_analysis std::string const &  conf  ) 
 

Read the analysis tasks.

Definition at line 410 of file colvar.C.

References acf_length, acf_normalize, acf_offset, acf_outfile, acf_stride, acf_type, colvarmodule::colvar_p(), enable(), colvarmodule::fatal_error(), colvarmodule::log(), colvarvalue::output_width(), runave_length, runave_os, runave_stride, task_corrfunc, task_fdiff_velocity, task_runave, colvarparse::to_lower_cppstr(), colvarmodule::wrap_string(), and x.

Referenced by colvar().

00411 {
00412 
00413   //   if (cvm::debug())
00414   //     cvm::log ("Parsing analysis flags for collective variable \""+
00415   //               this->name+"\".\n");
00416 
00417   runave_length = 0;
00418   bool b_runave = false;
00419   if (get_keyval (conf, "runAve", b_runave) && b_runave) {
00420 
00421     enable (task_runave);
00422 
00423     get_keyval (conf, "runAveLength", runave_length, 1000);
00424     get_keyval (conf, "runAveStride", runave_stride, 1);
00425 
00426     if ((cvm::restart_out_freq % runave_stride) != 0)
00427       cvm::fatal_error ("Error: runAveStride must be commensurate with the restart frequency.\n");
00428 
00429     std::string runave_outfile;
00430     get_keyval (conf, "runAveOutputFile", runave_outfile,
00431                 std::string (cvm::output_prefix+"."+
00432                              this->name+".runave.traj"));
00433 
00434     size_t const this_cv_width = x.output_width (cvm::cv_width);
00435     runave_os.open (runave_outfile.c_str());
00436     runave_os << "# " << cvm::wrap_string ("step", cvm::it_width-2)
00437               << "  "
00438               << cvm::wrap_string ("running average", this_cv_width)
00439               << " "
00440               << cvm::wrap_string ("running stddev", this_cv_width)
00441               << "\n";
00442   }
00443 
00444   acf_length = 0;
00445   bool b_acf = false;
00446   if (get_keyval (conf, "corrFunc", b_acf) && b_acf) {
00447 
00448     enable (task_corrfunc);
00449 
00450     std::string acf_colvar_name;
00451     get_keyval (conf, "corrFuncWithColvar", acf_colvar_name, this->name);
00452     if (acf_colvar_name == this->name) {
00453       cvm::log ("Calculating auto-correlation function.\n");
00454     } else {
00455       cvm::log ("Calculating correlation function with \""+
00456                 this->name+"\".\n");
00457     }
00458 
00459     std::string acf_type_str;
00460     get_keyval (conf, "corrFuncType", acf_type_str, to_lower_cppstr (std::string ("velocity")));
00461     if (acf_type_str == to_lower_cppstr (std::string ("coordinate"))) {
00462       acf_type = acf_coor;
00463     } else if (acf_type_str == to_lower_cppstr (std::string ("velocity"))) {
00464       acf_type = acf_vel;
00465       enable (task_fdiff_velocity);
00466       if (acf_colvar_name.size())
00467         (cvm::colvar_p (acf_colvar_name))->enable (task_fdiff_velocity);
00468     } else if (acf_type_str == to_lower_cppstr (std::string ("coordinate_p2"))) {
00469       acf_type = acf_p2coor;
00470     } else {
00471       cvm::fatal_error ("Unknown type of correlation function, \""+
00472                         acf_type_str+"\".\n");
00473     }
00474 
00475     get_keyval (conf, "corrFuncOffset", acf_offset, 0);
00476     get_keyval (conf, "corrFuncLength", acf_length, 1000);
00477     get_keyval (conf, "corrFuncStride", acf_stride, 1);
00478 
00479     if ((cvm::restart_out_freq % acf_stride) != 0)
00480       cvm::fatal_error ("Error: corrFuncStride must be commensurate with the restart frequency.\n");
00481 
00482     get_keyval (conf, "corrFuncNormalize", acf_normalize, true);
00483     get_keyval (conf, "corrFuncOutputFile", acf_outfile,
00484                 std::string (cvm::output_prefix+"."+this->name+
00485                              ".corrfunc.dat"));
00486   }
00487 }

bool colvar::periodic_boundaries colvarvalue const &  lb,
colvarvalue const &  ub
const
 

Is the interval defined by the two boundaries periodic?

Definition at line 972 of file colvar.C.

References colvarmodule::fatal_error(), period, and tasks.

00973 {
00974   if ( (!tasks[task_lower_boundary]) || (!tasks[task_upper_boundary]) ) {
00975     cvm::fatal_error ("Error: requesting to histogram the "
00976                       "collective variable \""+this->name+"\", but a "
00977                       "pair of lower and upper boundaries must be "
00978                       "defined.\n");
00979   }
00980 
00981   if (period > 0.0) {
00982     if ( ((std::sqrt (this->dist2 (lb, ub))) / this->width)
00983          < 1.0E-10 ) {
00984       return true;
00985     }
00986   }
00987 
00988   return false;
00989 }

bool colvar::periodic_boundaries  )  const
 

Is the interval defined by the two boundaries periodic?

Definition at line 991 of file colvar.C.

References colvarmodule::fatal_error(), lower_boundary, tasks, and upper_boundary.

Referenced by colvar().

00992 {
00993   if ( (!tasks[task_lower_boundary]) || (!tasks[task_upper_boundary]) ) {
00994     cvm::fatal_error ("Error: requesting to histogram the "
00995                       "collective variable \""+this->name+"\", but a "
00996                       "pair of lower and upper boundaries must be "
00997                       "defined.\n");
00998   }
00999 
01000   return periodic_boundaries (lower_boundary, upper_boundary);
01001 }

std::istream & colvar::read_restart std::istream &  is  ) 
 

Read the collective variable from a restart file.

Definition at line 1037 of file colvar.C.

References colvarmodule::fatal_error(), colvarmodule::log(), name, tasks, colvarvalue::type(), v_fdiff, v_reported, vr, x, x_reported, and xr.

01038 {
01039   size_t const start_pos = is.tellg();
01040 
01041   std::string conf;
01042   if ( !(is >> colvarparse::read_block ("colvar", conf)) ) {
01043     // this is not a colvar block
01044     is.clear();
01045     is.seekg (start_pos, std::ios::beg);
01046     is.setstate (std::ios::failbit);
01047     return is;
01048   }
01049 
01050   {
01051     std::string check_name = "";
01052     if ( (get_keyval (conf, "name", check_name,
01053                       std::string (""), colvarparse::parse_silent)) &&
01054          (check_name != name) )  {
01055       cvm::fatal_error ("Error: the state file does not match the "
01056                         "configuration file, at colvar \""+name+"\".\n");
01057     }
01058     if (check_name.size() == 0) {
01059       cvm::fatal_error ("Error: Collective variable in the "
01060                         "restart file without any identifier.\n");
01061     }
01062   }
01063 
01064   if ( !(get_keyval (conf, "x", x,
01065                      colvarvalue (x.type()), colvarparse::parse_silent)) ) {
01066     cvm::log ("Error: restart file does not contain "
01067               "the value of the colvar \""+
01068               name+"\" .\n");
01069   } else {
01070     cvm::log ("Restarting collective variable \""+name+"\" from value: "+
01071               cvm::to_str (x)+"\n");
01072   }
01073 
01074   if (tasks[colvar::task_extended_lagrangian]) {
01075 
01076     if ( !(get_keyval (conf, "extended_x", xr,
01077                        colvarvalue (x.type()), colvarparse::parse_silent)) &&
01078          !(get_keyval (conf, "extended_v", vr,
01079                        colvarvalue (x.type()), colvarparse::parse_silent)) ) {
01080       cvm::log ("Error: restart file does not contain "
01081                 "\"extended_x\" or \"extended_v\" for the colvar \""+
01082                 name+"\", but you requested \"extendedLagrangian\".\n");
01083     }
01084   }
01085 
01086   if (tasks[task_extended_lagrangian]) {
01087     x_reported = xr;
01088   } else {
01089     x_reported = x;
01090   }
01091 
01092   if (tasks[task_output_velocity]) {
01093 
01094     if ( !(get_keyval (conf, "v", v_fdiff,
01095                        colvarvalue (x.type()), colvarparse::parse_silent)) ) {
01096       cvm::log ("Error: restart file does not contain "
01097                 "the velocity for the colvar \""+
01098                 name+"\", but you requested \"outputVelocity\".\n");
01099     }
01100 
01101     if (tasks[task_extended_lagrangian]) {
01102       v_reported = vr;
01103     } else {
01104       v_reported = v_fdiff;
01105     }
01106   }
01107 
01108   return is;
01109 }

std::istream & colvar::read_traj std::istream &  is  ) 
 

Read the value from a collective variable trajectory file.

Definition at line 1112 of file colvar.C.

References ft_reported, colvarmodule::log(), tasks, v_reported, and x_reported.

01113 {
01114   size_t const start_pos = is.tellg();
01115 
01116   if (tasks[task_output_value]) {
01117 
01118     if (!(is >> x)) {
01119       cvm::log ("Error: in reading the value of colvar \""+
01120                 this->name+"\" from trajectory.\n");
01121       is.clear();
01122       is.seekg (start_pos, std::ios::beg);
01123       is.setstate (std::ios::failbit);
01124       return is;
01125     }
01126 
01127     if (tasks[task_extended_lagrangian]) {
01128       is >> xr;
01129       x_reported = xr;
01130     } else {
01131       x_reported = x;
01132     }
01133   }
01134 
01135   if (tasks[task_output_velocity]) {
01136 
01137     is >> v_fdiff;
01138 
01139     if (tasks[task_extended_lagrangian]) {
01140       is >> vr;
01141       v_reported = vr;
01142     } else {
01143       v_reported = v_fdiff;
01144     }
01145   }
01146 
01147   if (tasks[task_output_system_force]) {
01148 
01149     is >> ft;
01150 
01151     if (tasks[task_extended_lagrangian]) { 
01152       is >> fr;
01153       ft_reported = fr;
01154     } else {
01155       ft_reported = ft;
01156     }
01157   }
01158 
01159   if (tasks[task_output_applied_force]) {
01160     is >> f;
01161   }
01162 
01163   return is;
01164 }

void colvar::reset_bias_force  )  [inline]
 

Set the total biasing force to zero.

Definition at line 556 of file colvar.h.

References fb, and colvarvalue::reset().

00556                                      {
00557   fb.reset();
00558 }

colvarvalue const & colvar::system_force  )  const [inline]
 

Current system force (previously obtained from calc() or read_traj()). Note: this is calculated using the atomic forces from the last simulation step.

Total atom forces are read from the MD program, the total force acting on the collective variable is calculated summing those from all colvar components, the bias and walls forces are subtracted.

Definition at line 544 of file colvar.h.

00545 {
00546   return ft_reported;
00547 }

colvarvalue::Type colvar::type  )  const [inline]
 

Type of value.

Definition at line 526 of file colvar.h.

References colvarvalue::type(), and x.

Referenced by colvarbias::add_colvar(), calc_acf(), colvar(), colvarbias_harmonic::colvarbias_harmonic(), enable(), colvarbias_harmonic::energy_difference(), and colvarbias_harmonic::update().

00527 {
00528   return x.type();
00529 }

cvm::real colvar::update  ) 
 

Collect all forces on this colvar, integrate internal equations of motion of internal degrees of freedom; see also colvar::communicate_forces() return colvar energy if extended Lagrandian active.

Definition at line 827 of file colvar.C.

References colvarvalue::apply_constraints(), colvarmodule::boltzmann(), cvcs, colvarmodule::debug(), colvarmodule::decrease_depth(), dist2(), dist2_lgrad(), dist2_rgrad(), colvarmodule::dt(), ext_force_k, ext_gamma, ext_mass, ext_sigma, f, fj, fr, colvarmodule::increase_depth(), kinetic_energy, colvarmodule::log(), lower_wall, lower_wall_k, potential_energy, colvarmodule::rand_gaussian(), colvarmodule::real, colvarvalue::real_value, colvarvalue::reset(), tasks, colvarmodule::temperature(), upper_wall, upper_wall_k, vr, wrap(), x, x_old, and xr.

00828 {
00829   if (cvm::debug())
00830     cvm::log ("Updating colvar \""+this->name+"\".\n");
00831 
00832 
00833   // set to zero the applied force
00834   f.reset();
00835 
00836   // add the biases' force, which at this point should already have
00837   // been summed over each bias using this colvar
00838   f += fb;
00839 
00840 
00841   if (tasks[task_lower_wall] || tasks[task_upper_wall]) {
00842 
00843     // wall force
00844     colvarvalue fw (this->type());
00845 
00846     // if the two walls are applied concurrently, decide which is the
00847     // closer one (on a periodic colvar, both walls may be applicable
00848     // at the same time)
00849     if ( (!tasks[task_upper_wall]) ||
00850          (this->dist2 (x, lower_wall) < this->dist2 (x, upper_wall)) ) {
00851 
00852       cvm::real const grad = this->dist2_lgrad (x, lower_wall);
00853       if (grad < 0.0) {
00854         fw = -0.5 * lower_wall_k * grad;
00855         if (cvm::debug())
00856           cvm::log ("Applying a lower wall force ("+
00857                     cvm::to_str (fw)+") to \""+this->name+"\".\n");
00858         f += fw;
00859 
00860       }
00861 
00862     } else {
00863 
00864       cvm::real const grad = this->dist2_lgrad (x, upper_wall);
00865       if (grad > 0.0) {
00866         fw = -0.5 * upper_wall_k * grad;
00867         if (cvm::debug())
00868           cvm::log ("Applying an upper wall force ("+
00869                     cvm::to_str (fw)+") to \""+this->name+"\".\n");
00870         f += fw;
00871       }
00872     }
00873   }
00874 
00875   if (tasks[task_Jacobian_force]) {
00876     
00877     cvm::increase_depth();
00878     for (size_t i = 0; i < cvcs.size(); i++) {
00879       (cvcs[i])->calc_Jacobian_derivative();
00880     }
00881     cvm::decrease_depth();
00882 
00883     fj.reset();
00884     for (size_t i = 0; i < cvcs.size(); i++) {
00885       // linear combination is assumed
00886       fj += 1.0 / ( cvm::real (cvcs.size()) *  cvm::real ((cvcs[i])->sup_coeff) ) *
00887         (cvcs[i])->Jacobian_derivative();
00888     }
00889     fj *= cvm::boltzmann() * cvm::temperature();
00890 
00891     // the instantaneous Jacobian force was not included in the reported system force;
00892     // instead, it is subtracted from the applied force (silent Jacobian correction) 
00893     if (! tasks[task_report_Jacobian_force]) 
00894       f -= fj;
00895   }
00896 
00897 
00898   if (tasks[task_extended_lagrangian]) {
00899 
00900     cvm::real dt = cvm::dt();
00901 
00902     // the total force is applied to the fictitious mass, while the
00903     // atoms only feel the harmonic force
00904     // fr: extended coordinate force; f: colvar force applied to atomic coordinates
00905     fr   = f;
00906     fr  += (-0.5 * ext_force_k) * this->dist2_lgrad (xr, x);
00907     f    = (-0.5 * ext_force_k) * this->dist2_rgrad (xr, x);
00908 
00909     // leap frog: starting from x_i, f_i, v_(i-1/2)
00910     vr  += (0.5 * dt) * fr / ext_mass;
00911     // Because of leapfrog, kinetic energy at time i is approximate
00912     kinetic_energy = 0.5 * ext_mass * vr * vr;
00913     potential_energy = 0.5 * ext_force_k * this->dist2(xr, x);
00914     // leap to v_(i+1/2)
00915     if (tasks[task_langevin]) {
00916       vr -= dt * ext_gamma * vr.real_value;
00917       vr += dt * ext_sigma * cvm::rand_gaussian() / ext_mass; 
00918     } 
00919     vr  += (0.5 * dt) * fr / ext_mass; 
00920     xr  += dt * vr;
00921     xr.apply_constraints();
00922     if (this->b_periodic) this->wrap (xr);
00923   }
00924 
00925 
00926   if (tasks[task_fdiff_velocity]) {
00927     // set it for the next step
00928     x_old = x;
00929   }
00930 
00931   if (cvm::debug())
00932     cvm::log ("Done updating colvar \""+this->name+"\".\n");
00933   return (potential_energy + kinetic_energy);
00934 }

colvarvalue const & colvar::value  )  const [inline]
 

Current value (previously obtained from calc() or read_traj()).

Definition at line 532 of file colvar.h.

Referenced by calc(), calc_acf(), and communicate_forces().

00533 {
00534   return x_reported;
00535 }

colvarvalue const & colvar::velocity  )  const [inline]
 

Current velocity (previously obtained from calc() or read_traj()).

Definition at line 538 of file colvar.h.

Referenced by calc_acf().

00539 {
00540   return v_reported;
00541 }

void colvar::wrap colvarvalue x  )  const
 

Use the internal metrics (as from objects) to wrap a value into a standard interval.

Handles correctly symmetries and periodic boundary conditions

Definition at line 1028 of file colvar.C.

References cvcs.

Referenced by colvarbias_harmonic::update(), and update().

01029 {
01030   (cvcs[0])->wrap (x);
01031   return;
01032 }

void colvar::write_acf std::ostream &  os  )  [protected]
 

Save the ACF to a file.

Definition at line 1519 of file colvar.C.

References acf, acf_nframes, acf_normalize, acf_stride, colvarmodule::log(), name, and colvarmodule::real.

Referenced by ~colvar().

01520 {
01521   if (!acf_nframes)
01522     cvm::log ("Warning: ACF was not calculated (insufficient frames).\n");
01523   os.setf (std::ios::scientific, std::ios::floatfield);
01524   os << "# Autocorrelation function for collective variable \""
01525      << this->name << "\"\n";
01526   // one frame is used for normalization, the statistical sample is
01527   // hence decreased
01528   os << "# nframes = " << (acf_normalize ?
01529                            acf_nframes - 1 :
01530                            acf_nframes) << "\n";
01531 
01532   cvm::real const acf_norm = acf.front() / cvm::real (acf_nframes);
01533   std::vector<cvm::real>::iterator acf_i;
01534   size_t it = acf_offset;
01535   for (acf_i = acf.begin(); acf_i != acf.end(); acf_i++) {
01536     os << std::setw (cvm::it_width) << acf_stride * (it++) << " "
01537        << std::setprecision (cvm::cv_prec)
01538        << std::setw (cvm::cv_width)
01539        << ( acf_normalize ?
01540             (*acf_i)/(acf_norm * cvm::real (acf_nframes)) :
01541             (*acf_i)/(cvm::real (acf_nframes)) ) << "\n";
01542   }
01543 }

std::ostream & colvar::write_restart std::ostream &  os  ) 
 

Write the collective variable to a restart file.

Definition at line 1169 of file colvar.C.

References name, tasks, v_reported, vr, x, and xr.

01169                                                   {
01170 
01171   os << "colvar {\n"
01172      << "  name " << name << "\n"
01173      << "  x "
01174      << std::setprecision (cvm::cv_prec)
01175      << std::setw (cvm::cv_width)
01176      << x << "\n";
01177 
01178   if (tasks[task_output_velocity]) {
01179     os << "  v "
01180        << std::setprecision (cvm::cv_prec)
01181        << std::setw (cvm::cv_width)
01182        << v_reported << "\n";
01183   }
01184 
01185   if (tasks[task_extended_lagrangian]) {
01186     os << "  extended_x "
01187        << std::setprecision (cvm::cv_prec)
01188        << std::setw (cvm::cv_width)
01189        << xr << "\n"
01190        << "  extended_v "
01191        << std::setprecision (cvm::cv_prec)
01192        << std::setw (cvm::cv_width)
01193        << vr << "\n";
01194   }
01195 
01196   os << "}\n\n";
01197 
01198   return os;
01199 }  

std::ostream & colvar::write_traj std::ostream &  os  ) 
 

Output formatted values to the trajectory file.

Definition at line 1260 of file colvar.C.

References potential_energy, and tasks.

01261 {
01262   os << " ";
01263 
01264   if (tasks[task_output_value]) {
01265 
01266     if (tasks[task_extended_lagrangian]) {
01267       os << " "
01268          << std::setprecision (cvm::cv_prec) << std::setw (cvm::cv_width)
01269          << x;
01270     }
01271 
01272     os << " "
01273        << std::setprecision (cvm::cv_prec) << std::setw (cvm::cv_width)
01274        << x_reported;
01275   }
01276 
01277   if (tasks[task_output_velocity]) {
01278 
01279     if (tasks[task_extended_lagrangian]) {
01280       os << " "
01281          << std::setprecision (cvm::cv_prec) << std::setw (cvm::cv_width)
01282          << v_fdiff;
01283     }
01284 
01285     os << " "
01286        << std::setprecision (cvm::cv_prec) << std::setw (cvm::cv_width)
01287        << v_reported;
01288   }
01289 
01290   if (tasks[task_output_energy]) {
01291       os << " "
01292          << std::setprecision (cvm::cv_prec) << std::setw (cvm::cv_width)
01293          << potential_energy
01294          << " "
01295          << kinetic_energy;
01296   }
01297 
01298 
01299   if (tasks[task_output_system_force]) {
01300 
01301     if (tasks[task_extended_lagrangian]) {
01302       os << " "
01303          << std::setprecision (cvm::cv_prec) << std::setw (cvm::cv_width)
01304          << ft;
01305     }
01306 
01307     os << " "
01308        << std::setprecision (cvm::cv_prec) << std::setw (cvm::cv_width)
01309        << ft_reported;
01310   }
01311 
01312   if (tasks[task_output_applied_force]) {
01313     os << " "
01314        << std::setprecision (cvm::cv_prec) << std::setw (cvm::cv_width)
01315        << f;
01316   }
01317 
01318   return os;
01319 }

std::ostream & colvar::write_traj_label std::ostream &  os  ) 
 

Write a label to the trajectory file (comment line).

Definition at line 1202 of file colvar.C.

References colvarvalue::output_width(), tasks, colvarmodule::wrap_string(), and x.

01203 {
01204   size_t const this_cv_width = x.output_width (cvm::cv_width);
01205 
01206   os << " ";
01207 
01208   if (tasks[task_output_value]) {
01209 
01210     os << " "
01211        << cvm::wrap_string (this->name, this_cv_width);
01212 
01213     if (tasks[task_extended_lagrangian]) {
01214       // restraint center
01215       os << " r_"
01216          << cvm::wrap_string (this->name, this_cv_width-2);
01217     }
01218   }
01219 
01220   if (tasks[task_output_velocity]) {
01221 
01222     os << " v_"
01223        << cvm::wrap_string (this->name, this_cv_width-2);
01224 
01225     if (tasks[task_extended_lagrangian]) {
01226       // restraint center
01227       os << " vr_"
01228          << cvm::wrap_string (this->name, this_cv_width-3);
01229     }
01230   }
01231 
01232   if (tasks[task_output_energy]) {
01233       os << " Ep_"
01234          << cvm::wrap_string (this->name, this_cv_width-3)
01235          << " Ek_"
01236          << cvm::wrap_string (this->name, this_cv_width-3);
01237   }
01238 
01239   if (tasks[task_output_system_force]) {
01240 
01241     os << " fs_"
01242        << cvm::wrap_string (this->name, this_cv_width-2);
01243 
01244     if (tasks[task_extended_lagrangian]) {
01245       // restraint center
01246       os << " fr_"
01247          << cvm::wrap_string (this->name, this_cv_width-3);
01248     }
01249   }
01250 
01251   if (tasks[task_output_applied_force]) {
01252     os << " fa_"
01253        << cvm::wrap_string (this->name, this_cv_width-2);
01254   }
01255 
01256   return os;
01257 }


Member Data Documentation

std::vector<cvm::real> colvar::acf [protected]
 

ACF values.

Definition at line 399 of file colvar.h.

Referenced by calc_acf(), calc_coor_acf(), calc_p2coor_acf(), calc_vel_acf(), write_acf(), and ~colvar().

std::string colvar::acf_colvar_name [protected]
 

Collective variable with which the correlation is calculated (default: itself).

Definition at line 387 of file colvar.h.

Referenced by calc_acf().

size_t colvar::acf_length [protected]
 

Length of autocorrelation function (ACF).

Definition at line 389 of file colvar.h.

Referenced by calc_acf(), calc_coor_acf(), calc_p2coor_acf(), calc_vel_acf(), and parse_analysis().

size_t colvar::acf_nframes [protected]
 

Number of frames for each ACF point.

Definition at line 395 of file colvar.h.

Referenced by calc_acf(), calc_coor_acf(), calc_p2coor_acf(), calc_runave(), calc_vel_acf(), and write_acf().

bool colvar::acf_normalize [protected]
 

Normalize the ACF to a maximum value of 1?

Definition at line 397 of file colvar.h.

Referenced by parse_analysis(), and write_acf().

size_t colvar::acf_offset [protected]
 

After how many steps the ACF starts.

Definition at line 391 of file colvar.h.

Referenced by calc_acf(), and parse_analysis().

std::string colvar::acf_outfile [protected]
 

Name of the file to write the ACF.

Definition at line 401 of file colvar.h.

Referenced by parse_analysis(), and ~colvar().

size_t colvar::acf_stride [protected]
 

How many timesteps separate two ACF values.

Definition at line 393 of file colvar.h.

Referenced by parse_analysis(), and write_acf().

acf_type_e colvar::acf_type [protected]
 

Type of autocorrelation function (ACF).

Definition at line 417 of file colvar.h.

Referenced by parse_analysis().

std::list< std::list<colvarvalue> > colvar::acf_v_history [protected]
 

Time series of values and velocities used in correlation functions

Definition at line 374 of file colvar.h.

Referenced by calc_acf().

std::list< std::list<colvarvalue> >::iterator colvar::acf_v_history_p [protected]
 

Time series of values and velocities used in correlation functions (pointers)x

Definition at line 377 of file colvar.h.

Referenced by calc_acf().

std::list< std::list<colvarvalue> > colvar::acf_x_history [protected]
 

Time series of values and velocities used in correlation functions

Definition at line 374 of file colvar.h.

Referenced by calc_acf().

std::list< std::list<colvarvalue> >::iterator colvar::acf_x_history_p [protected]
 

Time series of values and velocities used in correlation functions (pointers)x

Definition at line 377 of file colvar.h.

Referenced by calc_acf().

std::vector<int> colvar::atom_ids
 

Sorted array of (zero-based) IDs for all atoms involved.

Definition at line 500 of file colvar.h.

Referenced by build_atom_list(), calc(), and enable().

std::vector<cvm::rvector> colvar::atomic_gradients
 

Array of atomic gradients collected from all cvcs with appropriate components, rotations etc. For scalar variables only!

Definition at line 505 of file colvar.h.

Referenced by build_atom_list(), and calc().

bool colvar::b_inverse_gradients
 

True if all objects are capable of calculating inverse gradients.

Definition at line 85 of file colvar.h.

Referenced by colvar(), and colvar::cvc::cvc().

bool colvar::b_Jacobian_force
 

True if all objects are capable of calculating Jacobian forces.

Definition at line 89 of file colvar.h.

Referenced by colvar().

bool colvar::b_linear
 

True if this is a linear combination of elements.

Definition at line 81 of file colvar.h.

Referenced by colvar().

bool colvar::b_periodic
 

Same as above, but also takes into account components with a variable period, such as distanceZ.

Definition at line 247 of file colvar.h.

Referenced by colvar(), and colvar::cvc::cvc().

std::vector<cvc *> colvar::cvcs [protected]
 

Array of objects.

Definition at line 492 of file colvar.h.

Referenced by build_atom_list(), calc(), colvar(), communicate_forces(), compare(), dist2(), dist2_lgrad(), dist2_rgrad(), n_components(), update(), wrap(), and ~colvar().

bool colvar::expand_boundaries
 

Expand the boundaries of multiples of width, to keep the value always within range.

Definition at line 252 of file colvar.h.

Referenced by colvar().

cvm::real colvar::ext_force_k [protected]
 

Restraint force constant.

Definition at line 210 of file colvar.h.

Referenced by colvar(), and update().

cvm::real colvar::ext_gamma [protected]
 

Friction coefficient for Langevin extended dynamics.

Definition at line 212 of file colvar.h.

Referenced by colvar(), and update().

cvm::real colvar::ext_mass [protected]
 

Mass of the restraint center.

Definition at line 208 of file colvar.h.

Referenced by colvar(), and update().

cvm::real colvar::ext_sigma [protected]
 

Amplitude of Gaussian white noise for Langevin extended dynamics.

Definition at line 214 of file colvar.h.

Referenced by colvar(), and update().

colvarvalue colvar::f
 

Total applied force; fr (if task_extended_lagrangian is defined), fb (if biases are applied) and the walls' forces (if defined) contribute to it.

Definition at line 235 of file colvar.h.

Referenced by communicate_forces(), enable(), and update().

colvarvalue colvar::fb
 

Bias force; reset_bias_force() should be called before the biases are updated.

Definition at line 230 of file colvar.h.

Referenced by add_bias_force(), enable(), and reset_bias_force().

colvarvalue colvar::fj [protected]
 

Jacobian force, when task_Jacobian_force is enabled.

Definition at line 220 of file colvar.h.

Referenced by enable(), and update().

colvarvalue colvar::fr [protected]
 

Harmonic restraint force.

Definition at line 217 of file colvar.h.

Referenced by colvar(), and update().

colvarvalue colvar::ft
 

Total force, as derived from the atomic trajectory; should equal the total system force plus .

Definition at line 239 of file colvar.h.

Referenced by calc(), and enable().

colvarvalue colvar::ft_reported [protected]
 

Cached reported system force.

Definition at line 223 of file colvar.h.

Referenced by calc(), enable(), and read_traj().

cvm::real colvar::kinetic_energy [protected]
 

If extended Lagrangian active: colvar energies (kinetic and harmonic potential).

Definition at line 453 of file colvar.h.

Referenced by colvar(), and update().

colvarvalue colvar::lower_boundary
 

Location of the lower boundary.

Definition at line 255 of file colvar.h.

Referenced by colvar(), and periodic_boundaries().

colvarvalue colvar::lower_wall
 

Location of the lower wall.

Definition at line 257 of file colvar.h.

Referenced by colvar(), and update().

cvm::real colvar::lower_wall_k
 

Force constant for the lower boundary potential (|x-xb|^2).

Definition at line 259 of file colvar.h.

Referenced by colvar(), and update().

std::string colvar::name
 

Name.

Definition at line 46 of file colvar.h.

Referenced by colvarbias::add_colvar(), calc_acf(), read_restart(), write_acf(), and write_restart().

cvm::real colvar::period
 

Period, if it is a constant.

Definition at line 243 of file colvar.h.

Referenced by colvar(), and periodic_boundaries().

cvm::real colvar::potential_energy [protected]
 

Definition at line 454 of file colvar.h.

Referenced by colvar(), update(), and write_traj().

colvarvalue colvar::runave [protected]
 

Current value of the running average.

Definition at line 445 of file colvar.h.

Referenced by calc_runave().

size_t colvar::runave_length [protected]
 

Length of running average series.

Definition at line 439 of file colvar.h.

Referenced by calc_runave(), and parse_analysis().

std::ofstream colvar::runave_os [protected]
 

Name of the file to write the running average.

Definition at line 443 of file colvar.h.

Referenced by calc_runave(), parse_analysis(), and ~colvar().

size_t colvar::runave_stride [protected]
 

Timesteps to skip between two values in the running average series.

Definition at line 441 of file colvar.h.

Referenced by parse_analysis().

cvm::real colvar::runave_variance [protected]
 

Current value of the square deviation from the running average.

Definition at line 447 of file colvar.h.

Referenced by calc_runave().

bool colvar::tasks[task_ntot]
 

Tasks performed by this colvar.

Definition at line 153 of file colvar.h.

Referenced by analyse(), calc(), calc_acf(), colvar(), disable(), enable(), periodic_boundaries(), read_restart(), read_traj(), update(), write_restart(), write_traj(), and write_traj_label().

colvarvalue colvar::upper_boundary
 

Location of the upper boundary.

Definition at line 262 of file colvar.h.

Referenced by colvar(), and periodic_boundaries().

colvarvalue colvar::upper_wall
 

Location of the upper wall.

Definition at line 264 of file colvar.h.

Referenced by colvar(), and update().

cvm::real colvar::upper_wall_k
 

Force constant for the upper boundary potential (|x-xb|^2).

Definition at line 266 of file colvar.h.

Referenced by colvar(), and update().

colvarvalue colvar::v_fdiff [protected]
 

Finite-difference velocity.

Definition at line 186 of file colvar.h.

Referenced by calc(), calc_acf(), enable(), and read_restart().

colvarvalue colvar::v_reported [protected]
 

Cached reported velocity.

Definition at line 200 of file colvar.h.

Referenced by calc(), calc_acf(), enable(), read_restart(), read_traj(), and write_restart().

colvarvalue colvar::vr [protected]
 

Velocity of the restraint center.

Definition at line 206 of file colvar.h.

Referenced by calc(), colvar(), read_restart(), update(), and write_restart().

cvm::real colvar::width
 

Typical fluctuation amplitude for this collective variable (e.g. local width of a free energy basin).

In metadynamics calculations, , this value is used to calculate the width of a gaussian. In ABF calculations, , it is used to calculate the grid spacing in the direction of this collective variable.

Definition at line 77 of file colvar.h.

Referenced by colvar().

colvarvalue colvar::x [protected]
 

Value of the colvar.

Definition at line 180 of file colvar.h.

Referenced by calc(), calc_runave(), colvar(), communicate_forces(), parse_analysis(), read_restart(), type(), update(), write_restart(), and write_traj_label().

std::list< std::list<colvarvalue> > colvar::x_history [protected]
 

Time series of values and velocities used in running averages.

Definition at line 380 of file colvar.h.

Referenced by calc_runave().

std::list< std::list<colvarvalue> >::iterator colvar::x_history_p [protected]
 

Time series of values and velocities used in correlation functions (pointers)x

Definition at line 383 of file colvar.h.

Referenced by calc_runave().

colvarvalue colvar::x_old [protected]
 

Previous value (to calculate velocities during analysis).

Definition at line 370 of file colvar.h.

Referenced by calc(), calc_acf(), enable(), and update().

colvarvalue colvar::x_reported [protected]
 

Cached reported value (x may be manipulated).

Definition at line 183 of file colvar.h.

Referenced by calc(), colvar(), read_restart(), and read_traj().

colvarvalue colvar::xr [protected]
 

Restraint center.

Definition at line 204 of file colvar.h.

Referenced by calc(), colvar(), read_restart(), update(), and write_restart().


The documentation for this class was generated from the following files:
Generated on Fri May 25 04:07:19 2012 for NAMD by  doxygen 1.3.9.1