Collective Variables Module - Developer Documentation
 All Classes Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
colvar.h
1 // -*- c++ -*-
2 
3 #ifndef COLVAR_H
4 #define COLVAR_H
5 
6 #include <iostream>
7 #include <iomanip>
8 #include <cmath>
9 
10 #include "colvarmodule.h"
11 #include "colvarvalue.h"
12 #include "colvarparse.h"
13 #include "colvardeps.h"
14 
15 
40 
41 class colvar : public colvarparse, public colvardeps {
42 
43 public:
44 
46  std::string name;
47 
49  colvarvalue const & value() const;
50 
52  colvarvalue const & actual_value() const;
53 
55  colvarvalue const & velocity() const;
56 
65  colvarvalue const & total_force() const;
66 
76 
78  static std::vector<feature *> cv_features;
79 
81  std::vector<feature *> &features() {
82  return cv_features;
83  }
84 
85  int refresh_deps();
86 
88  std::vector<colvarbias *> biases;
89 protected:
90 
91 
92  /*
93  extended:
94  H = H_{0} + \sum_{i} 1/2*\lambda*(S_i(x(t))-s_i(t))^2 \\
95  + \sum_{i} 1/2*m_i*(ds_i(t)/dt)^2 \\
96  + \sum_{t'<t} W * exp(-1/2*\sum_{i} (s_i(t')-s_i(t))^2/(\delta{}s_i)^2) \\
97  + \sum_{w} (\sum_{i}c_{w,i}s_i(t) - D_w)^M/(\Sigma_w)^M
98 
99  normal:
100  H = H_{0} + \sum_{t'<t} W * exp(-1/2*\sum_{i} (S_i(x(t'))-S_i(x(t)))^2/(\delta{}S_i)^2) \\
101  + \sum_{w} (\sum_{i}c_{w,i}S_i(t) - D_w)^M/(\Sigma_w)^M
102 
103  output:
104  H = H_{0} (only output S(x), no forces)
105 
106  Here:
107  S(x(t)) = x
108  s(t) = xr
109  DS = Ds = delta
110  */
111 
112 
115 
116  // TODO: implement functionality to treat these
117  // /// Vector of individual values from CVCs
118  // colvarvalue x_cvc;
119 
120  // /// Jacobian matrix of individual values from CVCs
121  // colvarvalue dx_cvc;
122 
125 
128 
129  inline colvarvalue fdiff_velocity(colvarvalue const &xold,
130  colvarvalue const &xnew)
131  {
132  // using the gradient of the square distance to calculate the
133  // velocity (non-scalar variables automatically taken into
134  // account)
135  cvm::real dt = cvm::dt();
136  return ( ( (dt > 0.0) ? (1.0/dt) : 1.0 ) *
137  0.5 * dist2_lgrad(xnew, xold) );
138  }
139 
142 
143  // Options for extended_lagrangian
156 
159 
162 
165 
166 public:
167 
168 
172 
175 
180 
183 
187 
188 
191 
192 
196 
205 
214 
216  bool periodic_boundaries() const;
217 
219  bool periodic_boundaries(colvarvalue const &lb, colvarvalue const &ub) const;
220 
221 
223  colvar();
224 
226  int init(std::string const &conf);
227 
229  int init_components(std::string const &conf);
230 
232  int init_grid_parameters(std::string const &conf);
233 
235  int init_extended_Lagrangian(std::string const &conf);
236 
238  int init_output_flags(std::string const &conf);
239 
240 private:
242  template<typename def_class_name> int init_components_type(std::string const &conf,
243  char const *def_desc,
244  char const *def_config_key);
245 
246 public:
247 
249  void setup();
250 
252  ~colvar();
253 
254 
256  int calc();
257 
261  int calc_cvcs(int first = 0, size_t num_cvcs = 0);
262 
264  int check_cvc_range(int first_cvc, size_t num_cvcs);
265 
267  int calc_cvc_values(int first, size_t num_cvcs);
269  int calc_cvc_gradients(int first, size_t num_cvcs);
271  int calc_cvc_total_force(int first, size_t num_cvcs);
273  int calc_cvc_Jacobians(int first, size_t num_cvcs);
274 
276  int collect_cvc_data();
277 
279  int collect_cvc_values();
281  int collect_cvc_gradients();
285  int collect_cvc_Jacobians();
288 
290  inline colvarvalue const applied_force() const
291  {
292  if (is_enabled(f_cv_extended_Lagrangian)) {
293  return fr;
294  }
295  return f;
296  }
297 
299  void reset_bias_force();
300 
302  void add_bias_force(colvarvalue const &force);
303 
305  void add_bias_force_actual_value(colvarvalue const &force);
306 
312 
315  void communicate_forces();
316 
318  int set_cvc_flags(std::vector<bool> const &flags);
319 
321  int update_cvc_flags();
322 
323 protected:
326 
329 
337 
340 
341 public:
343  inline size_t num_active_cvcs() const { return n_active_cvcs; }
344 
346  inline int get_time_step_factor() const {return time_step_factor;}
347 
352  cvm::real dist2(colvarvalue const &x1,
353  colvarvalue const &x2) const;
354 
360  colvarvalue const &x2) const;
361 
367  colvarvalue const &x2) const;
368 
373  void wrap(colvarvalue &x) const;
374 
375 
377  int parse_analysis(std::string const &conf);
379  void analyze();
380 
381 
383  std::istream & read_traj(std::istream &is);
385  std::ostream & write_traj(std::ostream &os);
387  std::ostream & write_traj_label(std::ostream &os);
388 
389 
391  std::istream & read_restart(std::istream &is);
393  std::ostream & write_restart(std::ostream &os);
394 
396  int write_output_files();
397 
398 
399 protected:
402 
405 
408 
411  std::list< std::list<colvarvalue> > acf_x_history, acf_v_history;
414  std::list< std::list<colvarvalue> >::iterator acf_x_history_p, acf_v_history_p;
415 
417  std::list< std::list<colvarvalue> > x_history;
420  std::list< std::list<colvarvalue> >::iterator x_history_p;
421 
424  std::string acf_colvar_name;
426  size_t acf_length;
428  size_t acf_offset;
430  size_t acf_stride;
432  size_t acf_nframes;
436  std::vector<cvm::real> acf;
438  std::string acf_outfile;
439 
441  enum acf_type_e {
451  };
452 
455 
457  int calc_vel_acf(std::list<colvarvalue> &v_history,
458  colvarvalue const &v);
459 
462  void calc_coor_acf(std::list<colvarvalue> &x_history,
463  colvarvalue const &x);
464 
467  void calc_p2coor_acf(std::list<colvarvalue> &x_history,
468  colvarvalue const &x);
469 
471  int calc_acf();
473  void write_acf(std::ostream &os);
474 
480  cvm::ofstream runave_os;
485 
487  void calc_runave();
488 
491  cvm::real potential_energy;
492 public:
493 
494 
495  // collective variable component base class
496  class cvc;
497 
498  // currently available collective variable components
499 
500  // scalar colvar components
501  class distance;
502  class distance_z;
503  class distance_xy;
504  class distance_inv;
505  class distance_pairs;
506  class angle;
507  class dipole_angle;
508  class dihedral;
509  class coordnum;
510  class selfcoordnum;
511  class groupcoordnum;
512  class h_bond;
513  class rmsd;
514  class orientation_angle;
515  class orientation_proj;
516  class tilt;
517  class spin_angle;
518  class gyration;
519  class inertia;
520  class inertia_z;
521  class eigenvector;
522  class alpha_dihedrals;
523  class alpha_angles;
524  class dihedPC;
525 
526  // non-scalar components
527  class distance_vec;
528  class distance_dir;
529  class cartesian;
530  class orientation;
531 
532 protected:
533 
535  std::vector<cvc *> cvcs;
536 
538  std::vector<bool> cvc_flags;
539 
542  void build_atom_list(void);
543 
544 private:
546  std::string scripted_function;
547 
550  std::vector<const colvarvalue *> sorted_cvc_values;
551 
552 public:
554  std::vector<int> atom_ids;
555 
559  std::vector<cvm::rvector> atomic_gradients;
560 
561  inline size_t n_components() const {
562  return cvcs.size();
563  }
564 };
565 
566 
567 inline colvarvalue const & colvar::value() const
568 {
569  return x_reported;
570 }
571 
572 
573 inline colvarvalue const & colvar::actual_value() const
574 {
575  return x;
576 }
577 
578 
579 inline colvarvalue const & colvar::velocity() const
580 {
581  return v_reported;
582 }
583 
584 
585 inline colvarvalue const & colvar::total_force() const
586 {
587  return ft_reported;
588 }
589 
590 
591 inline void colvar::add_bias_force(colvarvalue const &force)
592 {
593  if (cvm::debug()) {
594  cvm::log("Adding biasing force "+cvm::to_str(force)+" to colvar \""+name+"\".\n");
595  }
596  fb += force;
597 }
598 
599 
601 {
602  if (cvm::debug()) {
603  cvm::log("Adding biasing force "+cvm::to_str(force)+" to colvar \""+name+"\".\n");
604  }
605  fb_actual += force;
606 }
607 
608 
610  fb.type(value());
611  fb.reset();
612  fb_actual.type(value());
613  fb_actual.reset();
614 }
615 
616 #endif
617 
colvarvalue const & actual_value() const
Current actual value (not extended DOF)
Definition: colvar.h:573
cvm::real ext_sigma
Amplitude of Gaussian white noise for Langevin extended dynamics.
Definition: colvar.h:155
Colvar component: average distance between two groups of atoms, weighted as the sixth power...
Definition: colvarcomp.h:426
Colvar component: moment of inertia of an atom group (colvarvalue::type_scalar type, range [0:*))
Definition: colvarcomp.h:502
cvm::real period
Period, if this variable is periodic.
Definition: colvar.h:190
bool acf_normalize
Normalize the ACF to a maximum value of 1?
Definition: colvar.h:434
std::vector< int > atom_ids
Sorted array of (zero-based) IDs for all atoms involved.
Definition: colvar.h:554
Velocity ACF, scalar product between v(0) and v(t)
Definition: colvar.h:445
std::list< std::list< colvarvalue > > x_history
Time series of values and velocities used in running averages.
Definition: colvar.h:417
colvarvalue x_restart
Value read from the most recent state file (if any)
Definition: colvar.h:404
cvm::real dist2(colvarvalue const &x1, colvarvalue const &x2) const
Use the internal metrics (as from cvc objects) to calculate square distances and gradients.
Definition: colvar.cpp:1427
int time_step_factor
Definition: colvar.h:336
size_t acf_nframes
Number of frames for each ACF point.
Definition: colvar.h:432
colvarvalue upper_boundary
Location of the upper boundary.
Definition: colvar.h:207
cvm::real ext_mass
Mass of the restraint center.
Definition: colvar.h:149
Colvar component: hydrogen bond, defined as the product of a colvar::coordnum and 1/2*(1-cos((180-ang...
Definition: colvarcomp.h:897
Colvar component: orientation in space of an atom group, with respect to a set of reference coordinat...
Definition: colvarcomp.h:1047
void analyze()
Perform analysis tasks.
Definition: colvar.cpp:1756
size_t acf_offset
After how many steps the ACF starts.
Definition: colvar.h:428
Definition: colvarcomp.h:303
static real dt()
Time step of MD integrator (fs)
Definition: colvarmodule.h:669
colvarvalue runave
Current value of the running average.
Definition: colvar.h:482
std::string acf_colvar_name
Collective variable with which the correlation is calculated (default: itself)
Definition: colvar.h:424
colvarvalue const applied_force() const
Get the current applied force.
Definition: colvar.h:290
colvarvalue const & total_force() const
Current total force (previously obtained from calc() or read_traj()). Note: this is calculated using ...
Definition: colvar.h:585
acf_type_e acf_type
Type of autocorrelation function (ACF)
Definition: colvar.h:454
colvarvalue x_reported
Cached reported value (x may be manipulated)
Definition: colvar.h:124
Colvar component: coordination number between two groups (colvarvalue::type_scalar type...
Definition: colvarcomp.h:738
colvarvalue fb_actual
Bias force to the actual value (only useful with extended Lagrangian)
Definition: colvar.h:174
void add_bias_force_actual_value(colvarvalue const &force)
Apply a force to the actual value (only meaningful with extended Lagrangian)
Definition: colvar.h:600
cvm::real runave_variance
Current value of the square deviation from the running average.
Definition: colvar.h:484
int calc()
Calculate the colvar's value and related quantities.
Definition: colvar.cpp:793
int init_extended_Lagrangian(std::string const &conf)
Init extended Lagrangian parameters.
Definition: colvar.cpp:366
std::vector< cvm::real > acf
ACF values.
Definition: colvar.h:436
Colvar component: cosine of the angle of rotation with respect to a set of reference coordinates (col...
Definition: colvarcomp.h:1112
cvm::real active_cvc_square_norm
Sum of square coefficients for active cvcs.
Definition: colvar.h:328
colvarvalue v_fdiff
Finite-difference velocity.
Definition: colvar.h:127
colvarvalue f
Total applied force; fr (if extended_lagrangian is defined), fb (if biases are applied) and the walls...
Definition: colvar.h:179
cvm::ofstream runave_os
Name of the file to write the running average.
Definition: colvar.h:480
int parse_analysis(std::string const &conf)
Read the analysis tasks.
Definition: colvar.cpp:661
std::string acf_outfile
Name of the file to write the ACF.
Definition: colvar.h:438
size_t runave_stride
Timesteps to skip between two values in the running average series.
Definition: colvar.h:478
int init_components(std::string const &conf)
Parse the CVC configuration and allocate their data.
Definition: colvar.cpp:546
int init_grid_parameters(std::string const &conf)
Init defaults for grid options.
Definition: colvar.cpp:258
int write_output_files()
Write output files (if defined, e.g. in analysis mode)
Definition: colvar.cpp:1729
size_t acf_stride
How many timesteps separate two ACF values.
Definition: colvar.h:430
Colvar component: projection of the distance vector on a plane (colvarvalue::type_scalar type...
Definition: colvarcomp.h:397
cvm::real upper_wall_k
Force constant for the upper boundary potential (|x-xb|^2)
Definition: colvar.h:211
cvm::real lower_wall_k
Force constant for the lower boundary potential (|x-xb|^2)
Definition: colvar.h:202
int check_cvc_range(int first_cvc, size_t num_cvcs)
Ensure that the selected range of CVCs is consistent.
Definition: colvar.cpp:859
Colvar component: projection of the distance vector along an axis(colvarvalue::type_scalar type...
Definition: colvarcomp.h:354
void setup()
Get ready for a run and re-initialize internal data if needed.
Definition: colvar.cpp:746
colvar()
Constructor.
Definition: colvar.cpp:20
int collect_cvc_total_forces()
Same as colvar::collect_cvc_values but for total forces.
Definition: colvar.cpp:1080
size_t acf_length
Length of autocorrelation function (ACF)
Definition: colvar.h:426
int set_cvc_flags(std::vector< bool > const &flags)
Enables and disables individual CVCs based on flags.
Definition: colvar.cpp:1356
std::string name
Name.
Definition: colvar.h:46
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 num...
Definition: colvar.cpp:1934
void communicate_forces()
Communicate forces (previously calculated in colvar::update()) to the external degrees of freedom...
Definition: colvar.cpp:1294
Collective variables main module.
Colvar component: projection of the orientation vector onto a predefined axis (colvarvalue::type_scal...
Definition: colvarcomp.h:1135
int get_time_step_factor() const
returns time_step_factor
Definition: colvar.h:346
Definition: colvarcomp.h:1229
Coordinate ACF, scalar product between x(0) and x(t)
Definition: colvar.h:447
double real
Defining an abstract real number allows to switch precision.
Definition: colvarmodule.h:85
std::ostream & write_restart(std::ostream &os)
Write the collective variable to a restart file.
Definition: colvar.cpp:1591
std::ostream & write_traj_label(std::ostream &os)
Write a label to the trajectory file (comment line)
Definition: colvar.cpp:1624
void build_atom_list(void)
Initialize the sorted list of atom IDs for atoms involved in all cvcs (called when enabling f_cv_coll...
Definition: colvar.cpp:624
std::vector< colvarbias * > biases
List of biases that depend on this colvar.
Definition: colvar.h:88
int refresh_deps()
Definition: colvar.cpp:611
std::vector< cvc * > cvcs
Array of cvc objects.
Definition: colvar.h:530
int calc_cvcs(int first=0, size_t num_cvcs=0)
Calculate a subset of the colvar components (CVCs) currently active (default: all active CVCs) Note: ...
Definition: colvar.cpp:808
colvarvalue ft_reported
Cached reported total force.
Definition: colvar.h:164
int calc_colvar_properties()
Calculate the quantities associated to the colvar (but not to the CVCs)
Definition: colvar.cpp:1152
colvarvalue const & value() const
Current value (previously set by calc() or by read_traj())
Definition: colvar.h:567
Colvar component: moment of inertia of an atom group around a user-defined axis (colvarvalue::type_sc...
Definition: colvarcomp.h:525
Colvar component: alpha helix content of a contiguous segment of 5 or more residues, implemented as a sum of phi/psi dihedral angles and hydrogen bonds (colvarvalue::type_scalar type, range [0:1])
Definition: colvarcomp.h:975
size_t runave_length
Length of running average series.
Definition: colvar.h:476
colvarvalue v_reported
Cached reported velocity.
Definition: colvar.h:141
colvarvalue lower_boundary
Location of the lower boundary.
Definition: colvar.h:198
int collect_cvc_Jacobians()
Same as colvar::collect_cvc_values but for Jacobian derivatives/forces.
Definition: colvar.cpp:1131
The variable has a harmonic restraint around a moving center with fictitious mass; bias forces will b...
Definition: colvardeps.h:222
std::vector< feature * > & features()
Implementation of the feature list accessor for colvar.
Definition: colvar.h:81
bool hard_upper_boundary
Whether this colvar has a hard upper boundary.
Definition: colvar.h:213
colvarvalue lower_wall
Location of the lower wall.
Definition: colvar.h:200
Colvar component: projection of 3N coordinates onto an eigenvector(colvarvalue::type_scalar type...
Definition: colvarcomp.h:551
Value of a collective variable: this is a metatype which can be set at runtime. By default it is set ...
Definition: colvarvalue.h:34
colvarvalue dist2_rgrad(colvarvalue const &x1, colvarvalue const &x2) const
Use the internal metrics (as from cvc objects) to calculate square distances and gradients.
Definition: colvar.cpp:1447
int calc_cvc_Jacobians(int first, size_t num_cvcs)
Same as colvar::calc_cvc_values but for Jacobian derivatives/forces.
Definition: colvar.cpp:1110
int collect_cvc_gradients()
Same as colvar::collect_cvc_values but for gradients.
Definition: colvar.cpp:999
bool expand_boundaries
Expand the boundaries of multiples of width, to keep the value always within range.
Definition: colvar.h:195
Colvar component: dihedral between the centers of mass of four groups (colvarvalue::type_scalar type...
Definition: colvarcomp.h:687
colvarvalue fb
Bias force; reset_bias_force() should be called before the biases are updated.
Definition: colvar.h:171
colvarvalue ft
Total force, as derived from the atomic trajectory; should equal the system force plus f...
Definition: colvar.h:186
colvarvalue upper_wall
Location of the upper wall.
Definition: colvar.h:209
bool after_restart
True if a state file was just read.
Definition: colvar.h:407
colvarvalue fj
Jacobian force, when Jacobian_force is enabled.
Definition: colvar.h:161
std::list< std::list< colvarvalue > >::iterator acf_x_history_p
Definition: colvar.h:414
acf_type_e
Type of autocorrelation function (ACF)
Definition: colvar.h:441
void wrap(colvarvalue &x) const
Use the internal metrics (as from cvc objects) to wrap a value into a standard interval.
Definition: colvar.cpp:1457
colvarvalue vr
Velocity of the restraint center.
Definition: colvar.h:147
cvm::real ext_gamma
Friction coefficient for Langevin extended dynamics.
Definition: colvar.h:153
int calc_cvc_gradients(int first, size_t num_cvcs)
Same as colvar::calc_cvc_values but for gradients.
Definition: colvar.cpp:957
int calc_vel_acf(std::list< colvarvalue > &v_history, colvarvalue const &v)
Velocity ACF, scalar product between v(0) and v(t)
Definition: colvar.cpp:1888
cvm::real width
Typical fluctuation amplitude for this collective variable (e.g. local width of a free energy basin) ...
Definition: colvar.h:75
Colvar component: angle between the centers of mass of three groups (colvarvalue::type_scalar type...
Definition: colvarcomp.h:593
Colvar component: self-coordination number within a group (colvarvalue::type_scalar type...
Definition: colvarcomp.h:799
colvarvalue x_old
Previous value (to calculate velocities during analysis)
Definition: colvar.h:401
void reset()
Set to the null value for the data type currently defined.
Definition: colvarvalue.cpp:106
Colvar component: distance unit vector (direction) between centers of mass of two groups (colvarvalue...
Definition: colvarcomp.h:329
colvarvalue f_old
Applied force at the previous step (to be subtracted from total force if needed)
Definition: colvar.h:182
int calc_cvc_total_force(int first, size_t num_cvcs)
Same as colvar::calc_cvc_values but for total forces.
Definition: colvar.cpp:1051
Coordinate ACF, second order Legendre polynomial between x(0) and x(t) (does not work with scalar num...
Definition: colvar.h:450
std::vector< bool > cvc_flags
Flags to enable or disable cvcs at next colvar evaluation.
Definition: colvar.h:538
colvarvalue dist2_lgrad(colvarvalue const &x1, colvarvalue const &x2) const
Use the internal metrics (as from cvc objects) to calculate square distances and gradients.
Definition: colvar.cpp:1437
int collect_cvc_values()
Collect the values of the CVCs.
Definition: colvar.cpp:901
Unset type.
Definition: colvar.h:443
bool hard_lower_boundary
Whether this colvar has a hard lower boundary.
Definition: colvar.h:204
Colvar component: angle between the dipole of a molecule and an axis formed by two groups of atoms(co...
Definition: colvarcomp.h:641
void reset_bias_force()
Set the total biasing force to zero.
Definition: colvar.h:609
static std::string to_str(T const &x, size_t const &width=0, size_t const &prec=0)
Quick conversion of an object to a string.
Definition: colvarmodule.h:615
colvarvalue f_accumulated
Biasing force collected between updates, to be applied at next update for coarse-time-step colvars...
Definition: colvar.h:339
Colvar component: Radius of gyration of an atom group (colvarvalue::type_scalar type, range [0:*))
Definition: colvarcomp.h:474
Parent class for a member object of a bias, cv or cvc etc. containing features and their dependencies...
Definition: colvardeps.h:19
colvarvalue x
Value of the colvar.
Definition: colvar.h:114
A collective variable (main class); to be defined, it needs at least one object of a derived class of...
Definition: colvar.h:41
void add_bias_force(colvarvalue const &force)
Add to the total force from biases.
Definition: colvar.h:591
~colvar()
Destructor.
Definition: colvar.cpp:759
Colvar component: distance between the centers of mass of two groups (colvarvalue::type_scalar type...
Definition: colvarcomp.h:269
Colvar component: angle of rotation around a predefined axis (colvarvalue::type_scalar type...
Definition: colvarcomp.h:1162
static std::vector< feature * > cv_features
Implementation of the feature list for colvar.
Definition: colvar.h:78
std::ostream & write_traj(std::ostream &os)
Output formatted values to the trajectory file.
Definition: colvar.cpp:1675
colvarvalue const & velocity() const
Current velocity (previously set by calc() or by read_traj())
Definition: colvar.h:579
size_t n_active_cvcs
Number of CVC objects with an active flag.
Definition: colvar.h:325
bool periodic_boundaries() const
Is the interval defined by the two boundaries periodic?
Definition: colvar.cpp:1416
Colvar component (base class); most implementations of cvc utilize one or more colvarmodule::atom or ...
Definition: colvarcomp.h:63
colvarvalue xr
Restraint center.
Definition: colvar.h:145
Type type() const
Get the current type.
Definition: colvarvalue.h:159
size_t num_active_cvcs() const
Return the number of CVC objects with an active flag (as set by update_cvc_flags) ...
Definition: colvar.h:343
cvm::real ext_force_k
Restraint force constant.
Definition: colvar.h:151
int collect_cvc_data()
Collect quantities from CVCs and update aggregated data for the colvar.
Definition: colvar.cpp:836
int init_output_flags(std::string const &conf)
Init output flags.
Definition: colvar.cpp:431
std::list< std::list< colvarvalue > >::iterator x_history_p
Definition: colvar.h:420
int calc_cvc_values(int first, size_t num_cvcs)
Calculate the values of the given subset of CVCs.
Definition: colvar.cpp:870
Colvar component: dihedPC Projection of the config onto a dihedral principal component See e...
Definition: colvarcomp.h:1017
std::vector< cvm::rvector > atomic_gradients
Array of atomic gradients collected from all cvcs with appropriate components, rotations etc...
Definition: colvar.h:559
Colvar component: root mean square deviation (RMSD) of a group with respect to a set of reference coo...
Definition: colvarcomp.h:1196
Parsing functions for collective variables.
Base class containing parsing functions; all objects which need to parse input inherit from this...
Definition: colvarparse.h:18
void init()
Set the object ready to parse a new configuration string.
Definition: colvarparse.h:71
cvm::real update_forces_energy()
Collect all forces on this colvar, integrate internal equations of motion of internal degrees of free...
Definition: colvar.cpp:1199
static void log(std::string const &message)
Print a message to the main log.
Definition: colvarmodule.cpp:1467
Colvar component: coordination number between two groups (colvarvalue::type_scalar type...
Definition: colvarcomp.h:840
void write_acf(std::ostream &os)
Save the ACF to a file.
Definition: colvar.cpp:1956
Colvar component: angle of rotation with respect to a set of reference coordinates (colvarvalue::type...
Definition: colvarcomp.h:1088
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) ...
Definition: colvar.cpp:1914
std::istream & read_traj(std::istream &is)
Read the value from a collective variable trajectory file.
Definition: colvar.cpp:1541
void calc_runave()
Calculate the running average and its standard deviation.
Definition: colvar.cpp:1983
std::istream & read_restart(std::istream &is)
Read the collective variable from a restart file.
Definition: colvar.cpp:1468
colvarvalue fr
Harmonic restraint force.
Definition: colvar.h:158
cvm::real kinetic_energy
If extended Lagrangian active: colvar energies (kinetic and harmonic potential)
Definition: colvar.h:490
int calc_acf()
Calculate the auto-correlation function (ACF)
Definition: colvar.cpp:1785
int update_cvc_flags()
Updates the flags in the CVC objects, and their number.
Definition: colvar.cpp:1369
Colvar component: N1xN2 vector of pairwise distances (colvarvalue::type_vector type, range (0:*) for each component)
Definition: colvarcomp.h:451
static bool debug()
Whether debug output should be enabled (compile-time option)
Definition: colvarmodule.h:225
std::list< std::list< colvarvalue > > acf_x_history
Definition: colvar.h:411