Collective Variables Module - Developer Documentation
colvar.h
1 // -*- c++ -*-
2 
3 #ifndef COLVAR_H
4 #define COLVAR_H
5 
6 #include <iostream>
7 
8 #include "colvarmodule.h"
9 #include "colvarvalue.h"
10 #include "colvarparse.h"
11 #include "colvardeps.h"
12 
13 #ifdef LEPTON
14 #include "Lepton.h" // for runtime custom expressions
15 #endif
16 
41 
42 class colvar : public colvarparse, public colvardeps {
43 
44 public:
45 
47  std::string name;
48 
50  colvarvalue const & value() const;
51 
53  colvarvalue const & actual_value() const;
54 
56  cvm::real const & force_constant() const;
57 
59  colvarvalue const & velocity() const;
60 
69  colvarvalue const & total_force() const;
70 
80 
82  static std::vector<feature *> cv_features;
83 
85  virtual const std::vector<feature *> &features()
86  {
87  return cv_features;
88  }
89  virtual std::vector<feature *> &modify_features()
90  {
91  return cv_features;
92  }
93  static void delete_features() {
94  for (size_t i=0; i < cv_features.size(); i++) {
95  delete cv_features[i];
96  }
97  cv_features.clear();
98  }
99 
103  void do_feature_side_effects(int id);
104 
106  std::vector<colvarbias *> biases;
107 protected:
108 
109 
110  /*
111  extended:
112  H = H_{0} + \sum_{i} 1/2*\lambda*(S_i(x(t))-s_i(t))^2 \\
113  + \sum_{i} 1/2*m_i*(ds_i(t)/dt)^2 \\
114  + \sum_{t'<t} W * exp(-1/2*\sum_{i} (s_i(t')-s_i(t))^2/(\delta{}s_i)^2) \\
115  + \sum_{w} (\sum_{i}c_{w,i}s_i(t) - D_w)^M/(\Sigma_w)^M
116 
117  normal:
118  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) \\
119  + \sum_{w} (\sum_{i}c_{w,i}S_i(t) - D_w)^M/(\Sigma_w)^M
120 
121  output:
122  H = H_{0} (only output S(x), no forces)
123 
124  Here:
125  S(x(t)) = x
126  s(t) = xr
127  DS = Ds = delta
128  */
129 
130 
133 
134  // TODO: implement functionality to treat these
135  // /// Vector of individual values from CVCs
136  // colvarvalue x_cvc;
137 
138  // /// Jacobian matrix of individual values from CVCs
139  // colvarvalue dx_cvc;
140 
143 
146 
147  inline colvarvalue fdiff_velocity(colvarvalue const &xold,
148  colvarvalue const &xnew)
149  {
150  // using the gradient of the square distance to calculate the
151  // velocity (non-scalar variables automatically taken into
152  // account)
153  cvm::real dt = cvm::dt();
154  return ( ( (dt > 0.0) ? (1.0/dt) : 1.0 ) *
155  0.5 * dist2_lgrad(xnew, xold) );
156  }
157 
160 
161  // Options for extended_lagrangian
174 
177 
180 
183 
184 public:
185 
186 
190 
193 
198 
201 
205 
206 
209 
210 
214 
223 
232 
234  bool periodic_boundaries() const;
235 
237  bool periodic_boundaries(colvarvalue const &lb, colvarvalue const &ub) const;
238 
239 
241  colvar();
242 
244  int init(std::string const &conf);
245 
247  int init_components(std::string const &conf);
248 
250  int init_custom_function(std::string const &conf);
251 
253  int init_grid_parameters(std::string const &conf);
254 
256  int init_extended_Lagrangian(std::string const &conf);
257 
259  int init_output_flags(std::string const &conf);
260 
261 private:
263  template<typename def_class_name> int init_components_type(std::string const &conf,
264  char const *def_desc,
265  char const *def_config_key);
266 
267 public:
268 
270  void setup();
271 
273  ~colvar();
274 
275 
277  int calc();
278 
282  int calc_cvcs(int first = 0, size_t num_cvcs = 0);
283 
285  int check_cvc_range(int first_cvc, size_t num_cvcs);
286 
288  int calc_cvc_values(int first, size_t num_cvcs);
290  int calc_cvc_gradients(int first, size_t num_cvcs);
292  int calc_cvc_total_force(int first, size_t num_cvcs);
294  int calc_cvc_Jacobians(int first, size_t num_cvcs);
295 
297  int collect_cvc_data();
298 
300  int collect_cvc_values();
302  int collect_cvc_gradients();
306  int collect_cvc_Jacobians();
309 
311  inline colvarvalue const applied_force() const
312  {
313  if (is_enabled(f_cv_extended_Lagrangian)) {
314  return fr;
315  }
316  return f;
317  }
318 
320  void reset_bias_force();
321 
323  void add_bias_force(colvarvalue const &force);
324 
326  void add_bias_force_actual_value(colvarvalue const &force);
327 
333 
336  void communicate_forces();
337 
339  int set_cvc_flags(std::vector<bool> const &flags);
340 
342  int update_cvc_flags();
343 
344 protected:
347 
350 
353 
354 public:
356  inline size_t num_active_cvcs() const { return n_active_cvcs; }
357 
362  cvm::real dist2(colvarvalue const &x1,
363  colvarvalue const &x2) const;
364 
370  colvarvalue const &x2) const;
371 
377  colvarvalue const &x2) const;
378 
383  void wrap(colvarvalue &x) const;
384 
385 
387  int parse_analysis(std::string const &conf);
389  void analyze();
390 
391 
393  std::istream & read_traj(std::istream &is);
395  std::ostream & write_traj(std::ostream &os);
397  std::ostream & write_traj_label(std::ostream &os);
398 
399 
401  std::istream & read_restart(std::istream &is);
403  std::ostream & write_restart(std::ostream &os);
404 
406  int write_output_files();
407 
408 
409 protected:
412 
415 
418 
421  std::list< std::list<colvarvalue> > acf_x_history, acf_v_history;
424  std::list< std::list<colvarvalue> >::iterator acf_x_history_p, acf_v_history_p;
425 
427  std::list< std::list<colvarvalue> > x_history;
430  std::list< std::list<colvarvalue> >::iterator x_history_p;
431 
434  std::string acf_colvar_name;
436  size_t acf_length;
438  size_t acf_offset;
440  size_t acf_stride;
442  size_t acf_nframes;
446  std::vector<cvm::real> acf;
448  std::string acf_outfile;
449 
451  enum acf_type_e {
461  };
462 
465 
467  int calc_vel_acf(std::list<colvarvalue> &v_history,
468  colvarvalue const &v);
469 
472  void calc_coor_acf(std::list<colvarvalue> &x_history,
473  colvarvalue const &x);
474 
477  void calc_p2coor_acf(std::list<colvarvalue> &x_history,
478  colvarvalue const &x);
479 
481  int calc_acf();
483  void write_acf(std::ostream &os);
484 
490  std::string runave_outfile;
492  std::ostream *runave_os;
497 
499  void calc_runave();
500 
503  cvm::real potential_energy;
504 public:
505 
506 
507  // collective variable component base class
508  class cvc;
509 
510  // list of available collective variable components
511 
512  // scalar colvar components
513  class distance;
514  class distance_z;
515  class distance_xy;
516  class polar_theta;
517  class polar_phi;
518  class distance_inv;
519  class distance_pairs;
520  class angle;
521  class dipole_angle;
522  class dihedral;
523  class coordnum;
524  class selfcoordnum;
525  class groupcoordnum;
526  class h_bond;
527  class rmsd;
528  class orientation_angle;
529  class orientation_proj;
530  class tilt;
531  class spin_angle;
532  class gyration;
533  class inertia;
534  class inertia_z;
535  class eigenvector;
536  class alpha_dihedrals;
537  class alpha_angles;
538  class dihedPC;
539 
540  // non-scalar components
541  class distance_vec;
542  class distance_dir;
543  class cartesian;
544  class orientation;
545 
546 protected:
547 
549  std::vector<cvc *> cvcs;
550 
552  std::vector<bool> cvc_flags;
553 
556  void build_atom_list(void);
557 
558 private:
560  std::string scripted_function;
561 
564  std::vector<const colvarvalue *> sorted_cvc_values;
565 
566 #ifdef LEPTON
567  std::vector<Lepton::CompiledExpression *> value_evaluators;
569 
571  std::vector<Lepton::CompiledExpression *> gradient_evaluators;
572 
574  std::vector<double *> value_eval_var_refs;
575  std::vector<double *> grad_eval_var_refs;
576 
578  double dev_null;
579 #endif
580 
581 public:
583  std::vector<int> atom_ids;
584 
588  std::vector<cvm::rvector> atomic_gradients;
589 
590  inline size_t n_components() const {
591  return cvcs.size();
592  }
593 };
594 
595 inline cvm::real const & colvar::force_constant() const
596 {
597  return ext_force_k;
598 }
599 
600 inline colvarvalue const & colvar::value() const
601 {
602  return x_reported;
603 }
604 
605 
606 inline colvarvalue const & colvar::actual_value() const
607 {
608  return x;
609 }
610 
611 
612 inline colvarvalue const & colvar::velocity() const
613 {
614  return v_reported;
615 }
616 
617 
618 inline colvarvalue const & colvar::total_force() const
619 {
620  return ft_reported;
621 }
622 
623 
624 inline void colvar::add_bias_force(colvarvalue const &force)
625 {
626  if (cvm::debug()) {
627  cvm::log("Adding biasing force "+cvm::to_str(force)+" to colvar \""+name+"\".\n");
628  }
629  fb += force;
630 }
631 
632 
634 {
635  if (cvm::debug()) {
636  cvm::log("Adding biasing force "+cvm::to_str(force)+" to colvar \""+name+"\".\n");
637  }
638  fb_actual += force;
639 }
640 
641 
643  fb.type(value());
644  fb.reset();
645  fb_actual.type(value());
646  fb_actual.reset();
647 }
648 
649 #endif
650 
cvm::real ext_sigma
Amplitude of Gaussian white noise for Langevin extended dynamics.
Definition: colvar.h:173
Colvar component: average distance between two groups of atoms, weighted as the sixth power...
Definition: colvarcomp.h:498
Colvar component: moment of inertia of an atom group (colvarvalue::type_scalar type, range [0:*))
Definition: colvarcomp.h:580
cvm::real period
Period, if this variable is periodic.
Definition: colvar.h:208
bool acf_normalize
Normalize the ACF to a maximum value of 1?
Definition: colvar.h:444
std::vector< int > atom_ids
Sorted array of (zero-based) IDs for all atoms involved.
Definition: colvar.h:583
Velocity ACF, scalar product between v(0) and v(t)
Definition: colvar.h:455
std::list< std::list< colvarvalue > > x_history
Time series of values and velocities used in running averages.
Definition: colvar.h:427
colvarvalue x_restart
Value read from the most recent state file (if any)
Definition: colvar.h:414
virtual const std::vector< feature * > & features()
Implementation of the feature list accessor for colvar.
Definition: colvar.h:85
size_t acf_nframes
Number of frames for each ACF point.
Definition: colvar.h:442
colvarvalue upper_boundary
Location of the upper boundary.
Definition: colvar.h:225
cvm::real ext_mass
Mass of the restraint center.
Definition: colvar.h:167
Colvar component: hydrogen bond, defined as the product of a colvar::coordnum and 1/2*(1-cos((180-ang...
Definition: colvarcomp.h:975
colvarvalue const & total_force() const
Current total force (previously obtained from calc() or read_traj()). Note: this is calculated using ...
Definition: colvar.h:618
Colvar component: orientation in space of an atom group, with respect to a set of reference coordinat...
Definition: colvarcomp.h:1125
void analyze()
Perform analysis tasks.
Definition: colvar.cpp:2006
size_t acf_offset
After how many steps the ACF starts.
Definition: colvar.h:438
Definition: colvarcomp.h:319
static real dt()
Time step of MD integrator (fs)
Definition: colvarmodule.cpp:1756
colvarvalue runave
Current value of the running average.
Definition: colvar.h:494
std::string acf_colvar_name
Collective variable with which the correlation is calculated (default: itself)
Definition: colvar.h:434
Type type() const
Get the current type.
Definition: colvarvalue.h:159
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:1678
Colvar component: polar coordinate theta of a group (colvarvalue::type_scalar type, range [0:180])
Definition: colvarcomp.h:471
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:1698
acf_type_e acf_type
Type of autocorrelation function (ACF)
Definition: colvar.h:464
colvarvalue x_reported
Cached reported value (x may be manipulated)
Definition: colvar.h:142
cvm::real const & force_constant() const
Force constant of the spring.
Definition: colvar.h:595
Colvar component: coordination number between two groups (colvarvalue::type_scalar type...
Definition: colvarcomp.h:816
colvarvalue fb_actual
Bias force to the actual value (only useful with extended Lagrangian)
Definition: colvar.h:192
void add_bias_force_actual_value(colvarvalue const &force)
Apply a force to the actual value (only meaningful with extended Lagrangian)
Definition: colvar.h:633
void do_feature_side_effects(int id)
Definition: colvar.cpp:784
cvm::real runave_variance
Current value of the square deviation from the running average.
Definition: colvar.h:496
int calc()
Calculate the colvar&#39;s value and related quantities.
Definition: colvar.cpp:987
int init_extended_Lagrangian(std::string const &conf)
Init extended Lagrangian parameters.
Definition: colvar.cpp:533
std::vector< cvm::real > acf
ACF values.
Definition: colvar.h:446
Colvar component: cosine of the angle of rotation with respect to a set of reference coordinates (col...
Definition: colvarcomp.h:1190
cvm::real active_cvc_square_norm
Sum of square coefficients for active cvcs.
Definition: colvar.h:349
colvarvalue v_fdiff
Finite-difference velocity.
Definition: colvar.h:145
colvarvalue f
Total applied force; fr (if extended_lagrangian is defined), fb (if biases are applied) and the walls...
Definition: colvar.h:197
int parse_analysis(std::string const &conf)
Read the analysis tasks.
Definition: colvar.cpp:836
std::string acf_outfile
Name of the file to write the ACF.
Definition: colvar.h:448
size_t runave_stride
Timesteps to skip between two values in the running average series.
Definition: colvar.h:488
int init_components(std::string const &conf)
Parse the CVC configuration and allocate their data.
Definition: colvar.cpp:715
int init_grid_parameters(std::string const &conf)
Init defaults for grid options.
Definition: colvar.cpp:424
int write_output_files()
Write output files (if defined, e.g. in analysis mode)
Definition: colvar.cpp:1980
size_t acf_stride
How many timesteps separate two ACF values.
Definition: colvar.h:440
Colvar component: projection of the distance vector on a plane (colvarvalue::type_scalar type...
Definition: colvarcomp.h:413
cvm::real upper_wall_k
Force constant for the upper boundary potential (|x-xb|^2)
Definition: colvar.h:229
cvm::real lower_wall_k
Force constant for the lower boundary potential (|x-xb|^2)
Definition: colvar.h:220
int check_cvc_range(int first_cvc, size_t num_cvcs)
Ensure that the selected range of CVCs is consistent.
Definition: colvar.cpp:1060
Colvar component: projection of the distance vector along an axis(colvarvalue::type_scalar type...
Definition: colvarcomp.h:370
void setup()
Get ready for a run and re-initialize internal data if needed.
Definition: colvar.cpp:920
colvar()
Constructor.
Definition: colvar.cpp:21
int collect_cvc_total_forces()
Same as colvar::collect_cvc_values but for total forces.
Definition: colvar.cpp:1283
size_t acf_length
Length of autocorrelation function (ACF)
Definition: colvar.h:436
int set_cvc_flags(std::vector< bool > const &flags)
Enables and disables individual CVCs based on flags.
Definition: colvar.cpp:1607
std::string name
Name.
Definition: colvar.h:47
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:2184
void communicate_forces()
Communicate forces (previously calculated in colvar::update()) to the external degrees of freedom...
Definition: colvar.cpp:1521
Collective variables main module.
Colvar component: projection of the orientation vector onto a predefined axis (colvarvalue::type_scal...
Definition: colvarcomp.h:1213
Definition: colvarcomp.h:1307
Coordinate ACF, scalar product between x(0) and x(t)
Definition: colvar.h:457
double real
Defining an abstract real number allows to switch precision.
Definition: colvarmodule.h:76
std::ostream & write_restart(std::ostream &os)
Write the collective variable to a restart file.
Definition: colvar.cpp:1842
std::ostream & write_traj_label(std::ostream &os)
Write a label to the trajectory file (comment line)
Definition: colvar.cpp:1875
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:799
std::vector< colvarbias * > biases
List of biases that depend on this colvar.
Definition: colvar.h:106
std::vector< cvc * > cvcs
Array of cvc objects.
Definition: colvar.h:544
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:1002
colvarvalue ft_reported
Cached reported total force.
Definition: colvar.h:182
int calc_colvar_properties()
Calculate the quantities associated to the colvar (but not to the CVCs)
Definition: colvar.cpp:1355
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:356
Colvar component: moment of inertia of an atom group around a user-defined axis (colvarvalue::type_sc...
Definition: colvarcomp.h:603
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:1053
size_t runave_length
Length of running average series.
Definition: colvar.h:486
colvarvalue v_reported
Cached reported velocity.
Definition: colvar.h:159
void wrap(colvarvalue &x) const
Use the internal metrics (as from cvc objects) to wrap a value into a standard interval.
Definition: colvar.cpp:1708
colvarvalue lower_boundary
Location of the lower boundary.
Definition: colvar.h:216
int collect_cvc_Jacobians()
Same as colvar::collect_cvc_values but for Jacobian derivatives/forces.
Definition: colvar.cpp:1334
The variable has a harmonic restraint around a moving center with fictitious mass; bias forces will b...
Definition: colvardeps.h:266
bool hard_upper_boundary
Whether this colvar has a hard upper boundary.
Definition: colvar.h:231
colvarvalue lower_wall
Location of the lower wall.
Definition: colvar.h:218
Colvar component: projection of 3N coordinates onto an eigenvector(colvarvalue::type_scalar type...
Definition: colvarcomp.h:629
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 const & value() const
Current value (previously set by calc() or by read_traj())
Definition: colvar.h:600
int calc_cvc_Jacobians(int first, size_t num_cvcs)
Same as colvar::calc_cvc_values but for Jacobian derivatives/forces.
Definition: colvar.cpp:1313
int collect_cvc_gradients()
Same as colvar::collect_cvc_values but for gradients.
Definition: colvar.cpp:1209
bool expand_boundaries
Expand the boundaries of multiples of width, to keep the value always within range.
Definition: colvar.h:213
Colvar component: dihedral between the centers of mass of four groups (colvarvalue::type_scalar type...
Definition: colvarcomp.h:765
colvarvalue fb
Bias force; reset_bias_force() should be called before the biases are updated.
Definition: colvar.h:189
colvarvalue ft
Total force, as derived from the atomic trajectory; should equal the system force plus f...
Definition: colvar.h:204
Colvar component: polar coordinate phi of a group (colvarvalue::type_scalar type, range [-180:180]) ...
Definition: colvarcomp.h:441
colvarvalue upper_wall
Location of the upper wall.
Definition: colvar.h:227
bool after_restart
True if a state file was just read.
Definition: colvar.h:417
colvarvalue fj
Jacobian force, when Jacobian_force is enabled.
Definition: colvar.h:179
std::string runave_outfile
Name of the file to write the running average.
Definition: colvar.h:490
std::list< std::list< colvarvalue > >::iterator acf_x_history_p
Definition: colvar.h:424
acf_type_e
Type of autocorrelation function (ACF)
Definition: colvar.h:451
colvarvalue const & actual_value() const
Current actual value (not extended DOF)
Definition: colvar.h:606
colvarvalue vr
Velocity of the restraint center.
Definition: colvar.h:165
cvm::real ext_gamma
Friction coefficient for Langevin extended dynamics.
Definition: colvar.h:171
int calc_cvc_gradients(int first, size_t num_cvcs)
Same as colvar::calc_cvc_values but for gradients.
Definition: colvar.cpp:1174
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:2138
cvm::real width
Typical fluctuation amplitude for this collective variable (e.g. local width of a free energy basin) ...
Definition: colvar.h:79
Colvar component: angle between the centers of mass of three groups (colvarvalue::type_scalar type...
Definition: colvarcomp.h:671
Colvar component: self-coordination number within a group (colvarvalue::type_scalar type...
Definition: colvarcomp.h:877
colvarvalue x_old
Previous value (to calculate velocities during analysis)
Definition: colvar.h:411
int prev_timestep
Absolute timestep number when this colvar was last updated.
Definition: colvar.h:352
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:345
colvarvalue const applied_force() const
Get the current applied force.
Definition: colvar.h:311
colvarvalue f_old
Applied force at the previous step (to be subtracted from total force if needed)
Definition: colvar.h:200
int calc_cvc_total_force(int first, size_t num_cvcs)
Same as colvar::calc_cvc_values but for total forces.
Definition: colvar.cpp:1254
Coordinate ACF, second order Legendre polynomial between x(0) and x(t) (does not work with scalar num...
Definition: colvar.h:460
std::vector< bool > cvc_flags
Flags to enable or disable cvcs at next colvar evaluation.
Definition: colvar.h:552
int collect_cvc_values()
Collect the values of the CVCs.
Definition: colvar.cpp:1102
Unset type.
Definition: colvar.h:453
bool periodic_boundaries() const
Is the interval defined by the two boundaries periodic?
Definition: colvar.cpp:1667
int init_custom_function(std::string const &conf)
Parse parameters for custom function with Lepton.
Definition: colvar.cpp:416
bool hard_lower_boundary
Whether this colvar has a hard lower boundary.
Definition: colvar.h:222
Colvar component: angle between the dipole of a molecule and an axis formed by two groups of atoms(co...
Definition: colvarcomp.h:719
void reset_bias_force()
Set the total biasing force to zero.
Definition: colvar.h:642
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:607
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:1688
Colvar component: Radius of gyration of an atom group (colvarvalue::type_scalar type, range [0:*))
Definition: colvarcomp.h:552
Parent class for a member object of a bias, cv or cvc etc. containing features and their dependencies...
Definition: colvardeps.h:23
colvarvalue x
Value of the colvar.
Definition: colvar.h:132
A collective variable (main class); to be defined, it needs at least one object of a derived class of...
Definition: colvar.h:42
void add_bias_force(colvarvalue const &force)
Add to the total force from biases.
Definition: colvar.h:624
~colvar()
Destructor.
Definition: colvar.cpp:933
Colvar component: distance between the centers of mass of two groups (colvarvalue::type_scalar type...
Definition: colvarcomp.h:285
Colvar component: angle of rotation around a predefined axis (colvarvalue::type_scalar type...
Definition: colvarcomp.h:1240
static std::vector< feature * > cv_features
Implementation of the feature list for colvar.
Definition: colvar.h:82
std::ostream & write_traj(std::ostream &os)
Output formatted values to the trajectory file.
Definition: colvar.cpp:1926
size_t n_active_cvcs
Number of CVC objects with an active flag.
Definition: colvar.h:346
Colvar component (base class for collective variables)
Definition: colvarcomp.h:59
colvarvalue xr
Restraint center.
Definition: colvar.h:163
cvm::real ext_force_k
Restraint force constant.
Definition: colvar.h:169
int collect_cvc_data()
Collect quantities from CVCs and update aggregated data for the colvar.
Definition: colvar.cpp:1036
int init_output_flags(std::string const &conf)
Init output flags.
Definition: colvar.cpp:599
std::list< std::list< colvarvalue > >::iterator x_history_p
Definition: colvar.h:430
int calc_cvc_values(int first, size_t num_cvcs)
Calculate the values of the given subset of CVCs.
Definition: colvar.cpp:1071
Colvar component: dihedPC Projection of the config onto a dihedral principal component See e...
Definition: colvarcomp.h:1095
std::vector< cvm::rvector > atomic_gradients
Array of atomic gradients collected from all cvcs with appropriate components, rotations etc...
Definition: colvar.h:588
Colvar component: root mean square deviation (RMSD) of a group with respect to a set of reference coo...
Definition: colvarcomp.h:1274
std::ostream * runave_os
File to write the running average.
Definition: colvar.h:492
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:61
cvm::real update_forces_energy()
Collect all forces on this colvar, integrate internal equations of motion of internal degrees of free...
Definition: colvar.cpp:1402
static void log(std::string const &message)
Print a message to the main log.
Definition: colvarmodule.cpp:1518
Colvar component: coordination number between two groups (colvarvalue::type_scalar type...
Definition: colvarcomp.h:918
void write_acf(std::ostream &os)
Save the ACF to a file.
Definition: colvar.cpp:2206
Colvar component: angle of rotation with respect to a set of reference coordinates (colvarvalue::type...
Definition: colvarcomp.h:1166
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:2164
std::istream & read_traj(std::istream &is)
Read the value from a collective variable trajectory file.
Definition: colvar.cpp:1792
void calc_runave()
Calculate the running average and its standard deviation.
Definition: colvar.cpp:2233
std::istream & read_restart(std::istream &is)
Read the collective variable from a restart file.
Definition: colvar.cpp:1719
colvarvalue fr
Harmonic restraint force.
Definition: colvar.h:176
colvarvalue const & velocity() const
Current velocity (previously set by calc() or by read_traj())
Definition: colvar.h:612
cvm::real kinetic_energy
If extended Lagrangian active: colvar energies (kinetic and harmonic potential)
Definition: colvar.h:502
int calc_acf()
Calculate the auto-correlation function (ACF)
Definition: colvar.cpp:2035
int update_cvc_flags()
Updates the flags in the CVC objects, and their number.
Definition: colvar.cpp:1620
Colvar component: N1xN2 vector of pairwise distances (colvarvalue::type_vector type, range (0:*) for each component)
Definition: colvarcomp.h:529
static bool debug()
Whether debug output should be enabled (compile-time option)
Definition: colvarmodule.h:237
std::list< std::list< colvarvalue > > acf_x_history
Definition: colvar.h:421