Collective Variables Module - Developer Documentation
Loading...
Searching...
No Matches
colvar.h
1// -*- c++ -*-
2
3// This file is part of the Collective Variables module (Colvars).
4// The original version of Colvars and its updates are located at:
5// https://github.com/Colvars/colvars
6// Please update all Colvars source files before making any changes.
7// If you wish to distribute your changes, please submit them to the
8// Colvars repository at GitHub.
9
10#ifndef COLVAR_H
11#define COLVAR_H
12
13#include <functional>
14#include <list>
15#include <iosfwd>
16#include <map>
17#include <memory>
18
19#include "colvarmodule.h"
20#include "colvarvalue.h"
21#include "colvarparse.h"
22#include "colvardeps.h"
23
24#ifdef LEPTON
25#include "Lepton.h" // for runtime custom expressions
26#endif
27
52
53class colvar : public colvarparse, public colvardeps {
54
55public:
56
58 std::string name;
59
61 colvarvalue const & value() const;
62
64 colvarvalue const & actual_value() const;
65
67 colvarvalue const & run_ave() const;
68
70 cvm::real const & force_constant() const;
71
73 colvarvalue const & velocity() const;
74
83 colvarvalue const & total_force() const;
84
94
96 static std::vector<feature *> cv_features;
97
99 virtual const std::vector<feature *> &features() const
100 {
101 return cv_features;
102 }
103 virtual std::vector<feature *> &modify_features()
104 {
105 return cv_features;
106 }
107 static void delete_features() {
108 for (size_t i=0; i < cv_features.size(); i++) {
109 delete cv_features[i];
110 }
111 cv_features.clear();
112 }
113
117 void do_feature_side_effects(int id);
118
120 std::vector<colvarbias *> biases;
121
122protected:
123
124
125 /*
126 extended:
127 H = H_{0} + \sum_{i} 1/2*\lambda*(S_i(x(t))-s_i(t))^2 \\
128 + \sum_{i} 1/2*m_i*(ds_i(t)/dt)^2 \\
129 + \sum_{t'<t} W * exp(-1/2*\sum_{i} (s_i(t')-s_i(t))^2/(\delta{}s_i)^2) \\
130 + \sum_{w} (\sum_{i}c_{w,i}s_i(t) - D_w)^M/(\Sigma_w)^M
131
132 normal:
133 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) \\
134 + \sum_{w} (\sum_{i}c_{w,i}S_i(t) - D_w)^M/(\Sigma_w)^M
135
136 output:
137 H = H_{0} (only output S(x), no forces)
138
139 Here:
140 S(x(t)) = x
141 s(t) = x_ext
142 DS = Ds = delta
143 */
144
145
148
149 // TODO: implement functionality to treat these
150 // /// Vector of individual values from CVCs
151 // colvarvalue x_cvc;
152
153 // /// Jacobian matrix of individual values from CVCs
154 // colvarvalue dx_cvc;
155
158
161
162 inline colvarvalue fdiff_velocity(colvarvalue const &xold,
163 colvarvalue const &xnew)
164 {
165 // using the gradient of the square distance to calculate the
166 // velocity (non-scalar variables automatically taken into
167 // account)
168 cvm::real dt = cvm::dt();
169 return ( ( (dt > 0.0) ? (1.0/dt) : 1.0 ) *
170 0.5 * dist2_lgrad(xnew, xold) );
171 }
172
175
176 // Options for extended_lagrangian
193
196
199
202
203public:
204
205
209
212
217
220
224
227
230
233 bool expand_boundaries = false;
234
237
240
242 bool periodic_boundaries() const;
243
245 bool periodic_boundaries(colvarvalue const &lb, colvarvalue const &ub) const;
246
247
249 colvar();
250
252 int init(std::string const &conf);
253
256
258 int init_components(std::string const &conf);
259
261 int init_custom_function(std::string const &conf);
262
264 int init_grid_parameters(std::string const &conf);
265
267 int init_extended_Lagrangian(std::string const &conf);
268
270 int init_output_flags(std::string const &conf);
271
273 virtual int init_dependencies();
274
275private:
276
278 template <typename def_class_name>
279 void add_component_type(char const *description, char const *config_key);
280
282 int init_components_type(const std::string &conf, const char *config_key);
283
284public:
285
287 void setup();
288
290 ~colvar();
291
292
294 int calc();
295
297 int end_of_step();
298
302 int calc_cvcs(int first = 0, size_t num_cvcs = 0);
303
305 int check_cvc_range(int first_cvc, size_t num_cvcs);
306
308 int calc_cvc_values(int first, size_t num_cvcs);
310 int calc_cvc_gradients(int first, size_t num_cvcs);
312 int calc_cvc_total_force(int first, size_t num_cvcs);
314 int calc_cvc_Jacobians(int first, size_t num_cvcs);
315
317 int collect_cvc_data();
318
320 int collect_cvc_values();
329
331 inline colvarvalue const applied_force() const
332 {
333 if (is_enabled(f_cv_extended_Lagrangian)) {
334 return fr;
335 }
336 return f;
337 }
338
340 void reset_bias_force();
341
343 void add_bias_force(colvarvalue const &force);
344
346 void add_bias_force_actual_value(colvarvalue const &force);
347
353
356
359 void communicate_forces();
360
362 int set_cvc_flags(std::vector<bool> const &flags);
363
365 int update_cvc_flags();
366
368 int update_cvc_config(std::vector<std::string> const &confs);
369
371 int cvc_param_exists(std::string const &param_name);
372
374 cvm::real get_cvc_param(std::string const &param_name);
375
377 void const *get_cvc_param_ptr(std::string const &param_name);
378
380 colvarvalue const *get_cvc_param_grad(std::string const &param_name);
381
383 int set_cvc_param(std::string const &param_name, void const *new_value);
384
385protected:
386
388 size_t n_active_cvcs = 0;
389
392
395
398
399public:
400
402 inline size_t num_dimensions() const
403 {
404 return value().size();
405 }
406
408 inline size_t num_cvcs() const
409 {
410 return cvcs.size();
411 }
412
415 inline size_t num_active_cvcs() const
416 {
417 return n_active_cvcs;
418 }
419
424 cvm::real dist2(colvarvalue const &x1,
425 colvarvalue const &x2) const;
426
432 colvarvalue const &x2) const;
433
439 colvarvalue const &x2) const;
440
445 void wrap(colvarvalue &x_unwrapped) const;
446
448 int parse_analysis(std::string const &conf);
449
451 int analyze();
452
454 std::istream & read_traj(std::istream &is);
455
457 std::ostream & write_traj(std::ostream &os);
459 std::ostream & write_traj_label(std::ostream &os);
460
462 std::istream & read_state(std::istream &is);
463
466
468 int check_matching_state(std::string const &state_conf);
469
471 int set_state_params(std::string const &state_conf);
472
474 std::string const get_state_params() const;
475
477 std::ostream & write_state(std::ostream &os) const;
478
481
483 int write_output_files();
484
485protected:
486
489
492
495
498
501 std::list< std::list<colvarvalue> > acf_x_history, acf_v_history;
504 std::list< std::list<colvarvalue> >::iterator acf_x_history_p, acf_v_history_p;
505
507 std::list< std::list<colvarvalue> > x_history;
510 std::list< std::list<colvarvalue> >::iterator x_history_p;
511
514 std::string acf_colvar_name;
526 std::vector<cvm::real> acf;
528 std::string acf_outfile;
529
541 };
542
545
547 void calc_vel_acf(std::list<colvarvalue> &v_history,
548 colvarvalue const &v);
549
552 void calc_coor_acf(std::list<colvarvalue> &x_history,
553 colvarvalue const &x);
554
557 void calc_p2coor_acf(std::list<colvarvalue> &x_history,
558 colvarvalue const &x);
559
561 int calc_acf();
563 int write_acf(std::ostream &os);
564
570 std::string runave_outfile;
575
577 int calc_runave();
578
583
584public:
585
586 // collective variable component base class
587 class cvc;
588
589 // list of available collective variable components
590
591 // scalar colvar components
592 class distance;
593 class distance_z;
594 class distance_xy;
595 class polar_theta;
596 class polar_phi;
597 class distance_inv;
598 class distance_pairs;
599 class dipole_magnitude;
600 class angle;
601 class dipole_angle;
602 class dihedral;
603 class coordnum;
604 class selfcoordnum;
605 class groupcoordnum;
606 class h_bond;
607 class rmsd;
608 class orientation_angle;
609 class orientation_proj;
610 class tilt;
611 class spin_angle;
612 class gyration;
613 class inertia;
614 class inertia_z;
615 class eigenvector;
616 class alpha_dihedrals;
617 class alpha_angles;
618 class dihedPC;
619 class alch_lambda;
620 class alch_Flambda;
621 class CartesianBasedPath;
622 class aspath;
623 class azpath;
624 class gspath;
625 class gzpath;
626 class linearCombination;
627 class CVBasedPath;
628 class gspathCV;
629 class gzpathCV;
630 class aspathCV;
631 class azpathCV;
632 class euler_phi;
633 class euler_psi;
634 class euler_theta;
635 class neuralNetwork;
636 class customColvar;
637
638 // non-scalar components
639 class distance_vec;
640 class distance_dir;
641 class cartesian;
642 class orientation;
643
644 // components that do not handle any atoms directly
645 class map_total;
646
648 static const std::map<std::string, std::function<colvar::cvc *()>> &get_global_cvc_map()
649 {
650 return global_cvc_map;
651 }
652
654 static bool compare_cvc(const colvar::cvc* const i, const colvar::cvc* const j);
655
656protected:
657
659 std::vector<std::shared_ptr<colvar::cvc>> cvcs;
660
662 std::vector<bool> cvc_flags;
663
666 void build_atom_list(void);
667
669 std::string scripted_function;
670
673 std::vector<const colvarvalue *> sorted_cvc_values;
674
675#ifdef LEPTON
677 std::vector<Lepton::CompiledExpression *> value_evaluators;
678
680 std::vector<Lepton::CompiledExpression *> gradient_evaluators;
681
683 std::vector<double *> value_eval_var_refs;
684 std::vector<double *> grad_eval_var_refs;
685
687 double dev_null;
688#endif
689
691 static std::map<std::string, std::function<colvar::cvc *()>> global_cvc_map;
692
694 static std::map<std::string, std::string> global_cvc_desc_map;
695
697 std::vector<int> volmap_ids_;
698
699public:
700
702 std::vector<int> atom_ids;
703
707 std::vector<cvm::rvector> atomic_gradients;
708
710 virtual std::vector<std::vector<int> > get_atom_lists();
711
713 std::vector<int> const &get_volmap_ids();
714
715};
716
717
718inline cvm::real const & colvar::force_constant() const
719{
720 return ext_force_k;
721}
722
723
724inline colvarvalue const & colvar::value() const
725{
726 return x_reported;
727}
728
729
730inline colvarvalue const & colvar::actual_value() const
731{
732 return x;
733}
734
735
736inline colvarvalue const & colvar::run_ave() const
737{
738 return runave;
739}
740
741
742inline colvarvalue const & colvar::velocity() const
743{
744 return v_reported;
745}
746
747
748inline colvarvalue const & colvar::total_force() const
749{
750 return ft_reported;
751}
752
753
754inline void colvar::add_bias_force(colvarvalue const &force)
755{
757 std::string("applying a force to the variable \""+name+"\""));
758 if (cvm::debug()) {
759 cvm::log("Adding biasing force "+cvm::to_str(force)+" to colvar \""+name+"\".\n");
760 }
761 fb += force;
762}
763
764
766{
767 if (cvm::debug()) {
768 cvm::log("Adding biasing force "+cvm::to_str(force)+" to colvar \""+name+"\".\n");
769 }
770 fb_actual += force;
771}
772
773
775 fb.type(value());
776 fb.reset();
779}
780
781#endif
Definition: colvarcomp.h:1369
Definition: colvarcomp.h:1240
Definition: colvarcomp.h:1226
Definition: colvarcomp.h:1210
Colvar component: alpha helix content of a contiguous segment of 5 or more residues,...
Definition: colvarcomp.h:913
Colvar component: angle between the centers of mass of three groups (colvarvalue::type_scalar type,...
Definition: colvarcomp.h:635
Definition: colvarcomp.h:1472
Definition: colvarcomp.h:1442
Definition: colvarcomp.h:1488
Definition: colvarcomp.h:1457
Definition: colvarcomp.h:1179
Colvar component: coordination number between two groups (colvarvalue::type_scalar type,...
Definition: colvarcomp.h:748
custom expression of colvars
Definition: colvarcomp.h:1343
Colvar component (base class for collective variables)
Definition: colvarcomp.h:70
Colvar component: dihedPC Projection of the config onto a dihedral principal component See e....
Definition: colvarcomp.h:960
Colvar component: dihedral between the centers of mass of four groups (colvarvalue::type_scalar type,...
Definition: colvarcomp.h:711
Colvar component: angle between the dipole of a molecule and an axis formed by two groups of atoms(co...
Definition: colvarcomp.h:675
Colvar component: dipole magnitude of a molecule.
Definition: colvarcomp.h:532
Colvar component: distance unit vector (direction) between centers of mass of two groups (colvarvalue...
Definition: colvarcomp.h:372
Colvar component: average distance between two groups of atoms, weighted as the sixth power,...
Definition: colvarcomp.h:483
Colvar component: N1xN2 vector of pairwise distances (colvarvalue::type_vector type,...
Definition: colvarcomp.h:505
Definition: colvarcomp.h:349
Colvar component: projection of the distance vector on a plane (colvarvalue::type_scalar type,...
Definition: colvarcomp.h:427
Colvar component: projection of the distance vector along an axis(colvarvalue::type_scalar type,...
Definition: colvarcomp.h:395
Colvar component: distance between the centers of mass of two groups (colvarvalue::type_scalar type,...
Definition: colvarcomp.h:324
Colvar component: projection of 3N coordinates onto an eigenvector(colvarvalue::type_scalar type,...
Definition: colvarcomp.h:603
Definition: colvarcomp.h:1107
Definition: colvarcomp.h:1118
Definition: colvarcomp.h:1129
Colvar component: coordination number between two groups (colvarvalue::type_scalar type,...
Definition: colvarcomp.h:858
Colvar component: alternative path collective variable using geometry, variable s Allow any combinati...
Definition: colvarcomp.h:1408
Colvar component: alternative path collective variable using geometry, variable s For more informatio...
Definition: colvarcomp.h:1269
Colvar component: Radius of gyration of an atom group (colvarvalue::type_scalar type,...
Definition: colvarcomp.h:551
Definition: colvarcomp.h:1425
Colvar component: alternative path collective variable using geometry, variable z This should be merg...
Definition: colvarcomp.h:1291
Colvar component: hydrogen bond, defined as the product of a colvar::coordnum and 1/2*(1-cos((180-ang...
Definition: colvarcomp.h:886
Colvar component: moment of inertia of an atom group around a user-defined axis (colvarvalue::type_sc...
Definition: colvarcomp.h:585
Colvar component: moment of inertia of an atom group (colvarvalue::type_scalar type,...
Definition: colvarcomp.h:571
Current only linear combination of sub-CVCs is available.
Definition: colvarcomp.h:1311
Definition: colvarcomp.h:1540
Definition: colvarcomp.h:1509
Definition: colvarcomp.h:1035
Colvar component: cosine of the angle of rotation with respect to a set of reference coordinates (col...
Definition: colvarcomp.h:1060
Colvar component: orientation in space of an atom group, with respect to a set of reference coordinat...
Definition: colvarcomp.h:986
Colvar component: polar coordinate phi of a group (colvarvalue::type_scalar type, range [-180:180])
Definition: colvarcomp.h:447
Colvar component: polar coordinate theta of a group (colvarvalue::type_scalar type,...
Definition: colvarcomp.h:465
Colvar component: root mean square deviation (RMSD) of a group with respect to a set of reference coo...
Definition: colvarcomp.h:1144
Colvar component: self-coordination number within a group (colvarvalue::type_scalar type,...
Definition: colvarcomp.h:824
Colvar component: angle of rotation around a predefined axis (colvarvalue::type_scalar type,...
Definition: colvarcomp.h:1095
Colvar component: projection of the orientation vector onto a predefined axis (colvarvalue::type_scal...
Definition: colvarcomp.h:1075
A collective variable (main class); to be defined, it needs at least one object of a derived class of...
Definition: colvar.h:53
bool periodic_boundaries() const
Is the interval defined by the two boundaries periodic?
Definition: colvar.cpp:2202
void update_active_cvc_square_norm()
Update the sum of square coefficients for active cvcs.
Definition: colvar.cpp:2063
void communicate_forces()
Communicate forces (previously calculated in colvar::update()) to the external degrees of freedom.
Definition: colvar.cpp:1964
size_t num_cvcs() const
Number of CVC objects defined.
Definition: colvar.h:408
int end_of_step()
Carry out operations needed before next step is run.
Definition: colvar.cpp:1946
colvarvalue const * get_cvc_param_grad(std::string const &param_name)
Pointer to the derivative of the variable with respect to param_name.
Definition: colvar.cpp:2165
colvarvalue const & actual_value() const
Current actual value (not extended DOF)
Definition: colvar.h:730
colvarvalue const & run_ave() const
Current running average (if calculated as set by analysis flag)
Definition: colvar.h:736
colvarvalue const & velocity() const
Current velocity (previously set by calc() or by read_traj())
Definition: colvar.h:742
cvm::real active_cvc_square_norm
Sum of square coefficients for active cvcs.
Definition: colvar.h:391
colvarvalue x_restart
Value read from the most recent state file (if any)
Definition: colvar.h:494
size_t num_dimensions() const
Number of dimensions of the value of this colvar.
Definition: colvar.h:402
cvm::real runave_variance
Current value of the square deviation from the running average.
Definition: colvar.h:574
int calc_cvc_total_force(int first, size_t num_cvcs)
Same as colvar::calc_cvc_values but for total forces.
Definition: colvar.cpp:1584
colvarvalue f
Total applied force; fr (if extended_lagrangian is defined), fb (if biases are applied) and the walls...
Definition: colvar.h:216
std::list< std::list< colvarvalue > >::iterator x_history_p
Definition: colvar.h:510
std::vector< colvarbias * > biases
List of biases that depend on this colvar.
Definition: colvar.h:120
std::vector< const colvarvalue * > sorted_cvc_values
Definition: colvar.h:673
colvarvalue const & value() const
Current value (previously set by calc() or by read_traj())
Definition: colvar.h:724
int write_output_files()
Write output files (if defined, e.g. in analysis mode)
Definition: colvar.cpp:2640
std::ostream & write_state(std::ostream &os) const
Write the colvar's state to a formatted output stream.
Definition: colvar.cpp:2479
int parse_analysis(std::string const &conf)
Read the analysis tasks.
Definition: colvar.cpp:1016
int init_custom_function(std::string const &conf)
Parse parameters for custom function with Lepton.
Definition: colvar.cpp:484
cvm::real ext_gamma
Friction coefficient for Langevin extended dynamics.
Definition: colvar.h:190
std::string acf_outfile
Name of the file to write the ACF.
Definition: colvar.h:528
cvm::real ext_mass
Mass of the restraint center.
Definition: colvar.h:186
int init_components_type(const std::string &conf, const char *config_key)
Initialize any CVC objects matching the given key.
Definition: colvar.cpp:778
std::ostream & write_traj(std::ostream &os)
Output formatted values to the trajectory file.
Definition: colvar.cpp:2586
cvm::real dist2(colvarvalue const &x1, colvarvalue const &x2) const
Use the internal metrics (as from colvar::cvc objects) to calculate square distances and gradients.
Definition: colvar.cpp:2213
bool after_restart
True if a state file was just read.
Definition: colvar.h:497
std::string name
Name.
Definition: colvar.h:58
colvarvalue v_ext
Velocity of the restraint center.
Definition: colvar.h:182
virtual std::vector< std::vector< int > > get_atom_lists()
Get vector of vectors of atom IDs for all atom groups.
Definition: colvar.cpp:1252
size_t acf_nframes
Number of frames for each ACF point.
Definition: colvar.h:522
int set_state_params(std::string const &state_conf)
Read the values of colvar mutable data from a string (used by both versions of read_state())
Definition: colvar.cpp:2344
bool matching_state
Flag used to tell if the state string being read is for this colvar.
Definition: colvar.h:488
int calc()
Calculate the colvar's value and related quantities.
Definition: colvar.cpp:1339
cvm::real period
Period, if this variable is periodic.
Definition: colvar.h:226
cvm::real kinetic_energy
If extended Lagrangian active: colvar kinetic energy.
Definition: colvar.h:580
void calc_vel_acf(std::list< colvarvalue > &v_history, colvarvalue const &v)
Velocity ACF, scalar product between v(0) and v(t)
Definition: colvar.cpp:2799
colvarvalue lower_boundary
Location of the lower boundary.
Definition: colvar.h:236
int update_cvc_config(std::vector< std::string > const &confs)
Modify the configuration of CVCs (currently, only base class data)
Definition: colvar.cpp:2098
int init_output_flags(std::string const &conf)
Init output flags.
Definition: colvar.cpp:731
int collect_cvc_values()
Collect the values of the CVCs.
Definition: colvar.cpp:1459
std::istream & read_traj(std::istream &is)
Read the value from a collective variable trajectory file.
Definition: colvar.cpp:2429
int check_cvc_range(int first_cvc, size_t num_cvcs)
Ensure that the selected range of CVCs is consistent.
Definition: colvar.cpp:1417
colvarvalue fb
Bias force; reset_bias_force() should be called before the biases are updated.
Definition: colvar.h:208
size_t n_active_cvcs
Number of CVC objects with an active flag.
Definition: colvar.h:388
std::vector< cvm::real > acf
ACF values.
Definition: colvar.h:526
int cvc_param_exists(std::string const &param_name)
Whether this named parameter exists (in the first and only component)
Definition: colvar.cpp:2133
std::list< std::list< colvarvalue > >::iterator acf_x_history_p
Definition: colvar.h:504
colvarvalue dist2_rgrad(colvarvalue const &x1, colvarvalue const &x2) const
Use the internal metrics (as from colvar::cvc objects) to calculate square distances and gradients.
Definition: colvar.cpp:2251
size_t runave_length
Length of running average series.
Definition: colvar.h:566
acf_type_e
Type of autocorrelation function (ACF)
Definition: colvar.h:531
@ acf_coor
Coordinate ACF, scalar product between x(0) and x(t)
Definition: colvar.h:537
@ acf_notset
Unset type.
Definition: colvar.h:533
@ acf_vel
Velocity ACF, scalar product between v(0) and v(t)
Definition: colvar.h:535
@ acf_p2coor
Coordinate ACF, second order Legendre polynomial between x(0) and x(t) (does not work with scalar num...
Definition: colvar.h:540
colvarvalue fj
Jacobian force, when Jacobian_force is enabled.
Definition: colvar.h:198
colvarvalue const applied_force() const
Get the current applied force.
Definition: colvar.h:331
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:2844
int init_grid_parameters(std::string const &conf)
Init defaults for grid options.
Definition: colvar.cpp:500
int collect_cvc_Jacobians()
Same as colvar::collect_cvc_values but for Jacobian derivatives/forces.
Definition: colvar.cpp:1666
void wrap(colvarvalue &x_unwrapped) const
Use the internal metrics (as from colvar::cvc objects) to wrap a value into a standard interval.
Definition: colvar.cpp:2271
void do_feature_side_effects(int id)
Definition: colvar.cpp:959
int calc_runave()
Calculate the running average and its standard deviation.
Definition: colvar.cpp:2924
cvm::step_number prev_timestep
Absolute timestep number when this colvar was last updated.
Definition: colvar.h:397
colvarvalue prev_x_ext
Previous value of the restraint center;.
Definition: colvar.h:180
int check_matching_state(std::string const &state_conf)
Check the name of the bias vs. the given string, set the matching_state flag accordingly.
Definition: colvar.cpp:2320
std::string scripted_function
Name of scripted function to be used.
Definition: colvar.h:669
void add_bias_force(colvarvalue const &force)
Add to the total force from biases.
Definition: colvar.h:754
virtual int init_dependencies()
Initialize dependency tree.
Definition: colvar.cpp:1089
std::ostream & write_traj_label(std::ostream &os)
Write a label to the trajectory file (comment line)
Definition: colvar.cpp:2535
int calc_colvar_properties()
Calculate the quantities associated to the colvar (but not to the CVCs)
Definition: colvar.cpp:1688
colvar()
Constructor.
Definition: colvar.cpp:32
int calc_cvc_Jacobians(int first, size_t num_cvcs)
Same as colvar::calc_cvc_values but for Jacobian derivatives/forces.
Definition: colvar.cpp:1645
int collect_cvc_total_forces()
Same as colvar::collect_cvc_values but for total forces.
Definition: colvar.cpp:1613
bool expand_boundaries
Expand the boundaries of multiples of width, to keep the value always within range.
Definition: colvar.h:233
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:2824
static bool compare_cvc(const colvar::cvc *const i, const colvar::cvc *const j)
function for sorting cvcs by their names
Definition: colvar.cpp:55
static std::map< std::string, std::function< colvar::cvc *()> > global_cvc_map
A global mapping of cvc names to the cvc constructors.
Definition: colvar.h:691
std::list< std::list< colvarvalue > > x_history
Time series of values and velocities used in running averages.
Definition: colvar.h:507
colvarvalue v_reported
Cached reported velocity.
Definition: colvar.h:174
colvarvalue fr
Applied force on extended DOF, for output (unscaled if using MTS)
Definition: colvar.h:195
colvarvalue upper_boundary
Location of the upper boundary.
Definition: colvar.h:239
std::list< std::list< colvarvalue > > acf_x_history
Definition: colvar.h:501
std::string const get_state_params() const
Write the state information of this colvar in a block of text, suitable for later parsing.
Definition: colvar.cpp:2491
bool acf_normalize
Normalize the ACF to a maximum value of 1?
Definition: colvar.h:524
colvarvalue v_fdiff
Finite-difference velocity.
Definition: colvar.h:160
colvarvalue prev_v_ext
Previous velocity of the restraint center.
Definition: colvar.h:184
int init_extended_Lagrangian(std::string const &conf)
Init extended Lagrangian parameters.
Definition: colvar.cpp:648
static std::vector< feature * > cv_features
Implementation of the feature list for colvar.
Definition: colvar.h:96
cvm::real ext_sigma
Amplitude of Gaussian white noise for Langevin extended dynamics.
Definition: colvar.h:192
std::vector< cvm::rvector > atomic_gradients
Array of atomic gradients collected from all cvcs with appropriate components, rotations etc....
Definition: colvar.h:707
cvm::real width
Typical fluctuation amplitude for this collective variable (e.g. local width of a free energy basin)
Definition: colvar.h:93
colvarvalue runave
Current value of the running average.
Definition: colvar.h:572
void const * get_cvc_param_ptr(std::string const &param_name)
Get a pointer to the named parameter (from the first and only component)
Definition: colvar.cpp:2154
size_t acf_offset
After how many steps the ACF starts.
Definition: colvar.h:518
colvarvalue ft_reported
Cached reported total force.
Definition: colvar.h:201
std::vector< std::shared_ptr< colvar::cvc > > cvcs
Array of components objects.
Definition: colvar.h:659
static std::map< std::string, std::string > global_cvc_desc_map
A global mapping of cvc names to the corresponding descriptions.
Definition: colvar.h:694
acf_type_e acf_type
Type of autocorrelation function (ACF)
Definition: colvar.h:544
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:976
std::istream & read_state(std::istream &is)
Read the colvar's state from a formatted input stream.
Definition: colvar.cpp:2290
cvm::real get_cvc_param(std::string const &param_name)
Get the value of the named parameter (from the first and only component)
Definition: colvar.cpp:2143
cvm::real wrap_center
Center of wrapping, if this variable is periodic.
Definition: colvar.h:229
void add_bias_force_actual_value(colvarvalue const &force)
Apply a force to the actual value (only meaningful with extended Lagrangian)
Definition: colvar.h:765
colvarvalue f_old
Applied force at the previous step (to be subtracted from total force if needed)
Definition: colvar.h:219
int update_cvc_flags()
Updates the flags in the CVC objects, and their number.
Definition: colvar.cpp:2074
cvm::real const & force_constant() const
Force constant of the spring.
Definition: colvar.h:718
cvm::real ext_force_k
Restraint force constant.
Definition: colvar.h:188
colvarvalue dist2_lgrad(colvarvalue const &x1, colvarvalue const &x2) const
Use the internal metrics (as from colvar::cvc objects) to calculate square distances and gradients.
Definition: colvar.cpp:2232
void define_component_types()
Populate the map of available CVC types.
Definition: colvar.cpp:838
~colvar()
Destructor.
Definition: colvar.cpp:1278
int calc_acf()
Calculate the auto-correlation function (ACF)
Definition: colvar.cpp:2704
void reset_bias_force()
Set the total biasing force to zero.
Definition: colvar.h:774
colvarvalue ft
Total force, as derived from the atomic trajectory; should equal the system force plus f.
Definition: colvar.h:223
std::string acf_colvar_name
Collective variable with which the correlation is calculated (default: itself)
Definition: colvar.h:514
colvarvalue const & total_force() const
Current total force (previously obtained from calc() or read_traj()). Note: this is calculated using ...
Definition: colvar.h:748
int init_components(std::string const &conf)
Parse the CVC configuration and allocate their data.
Definition: colvar.cpp:903
int collect_cvc_gradients()
Same as colvar::collect_cvc_values but for gradients.
Definition: colvar.cpp:1567
colvarvalue fb_actual
Bias force to the actual value (only useful with extended Lagrangian)
Definition: colvar.h:211
int collect_cvc_data()
Collect quantities from CVCs and update aggregated data for the colvar.
Definition: colvar.cpp:1388
size_t runave_stride
Timesteps to skip between two values in the running average series.
Definition: colvar.h:568
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:1354
colvarvalue x_reported
Cached reported value (x may be manipulated)
Definition: colvar.h:157
colvarvalue x_old
Previous value (to calculate velocities during analysis)
Definition: colvar.h:491
int init(std::string const &conf)
Main init function.
Definition: colvar.cpp: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:1765
std::vector< int > const & get_volmap_ids()
Volmap numeric IDs, one for each CVC (-1 if not available)
Definition: colvar.cpp:1263
cvm::real potential_energy
If extended Lagrangian active: colvar harmonic potential.
Definition: colvar.h:582
colvarvalue x
Value of the colvar.
Definition: colvar.h:147
static const std::map< std::string, std::function< colvar::cvc *()> > & get_global_cvc_map()
A global mapping of cvc names to the cvc constructors.
Definition: colvar.h:648
size_t acf_stride
How many timesteps separate two ACF values.
Definition: colvar.h:520
size_t num_active_cvcs() const
number of CVC objects with an active flag (as set by update_cvc_flags)
Definition: colvar.h:415
void update_extended_Lagrangian()
Integrate equations of motion of extended Lagrangian coordinate if needed.
Definition: colvar.cpp:1811
int calc_cvc_values(int first, size_t num_cvcs)
Calculate the values of the given subset of CVCs.
Definition: colvar.cpp:1428
void add_component_type(char const *description, char const *config_key)
Declare an available CVC type and its description, register them in the global map.
Definition: colvar.cpp:767
std::vector< bool > cvc_flags
Flags to enable or disable cvcs at next colvar evaluation.
Definition: colvar.h:662
std::vector< int > volmap_ids_
Volmap numeric IDs, one for each CVC (-1 if not available)
Definition: colvar.h:697
size_t acf_length
Length of autocorrelation function (ACF)
Definition: colvar.h:516
std::string runave_outfile
Name of the file to write the running average.
Definition: colvar.h:570
int set_cvc_param(std::string const &param_name, void const *new_value)
Set the named parameter in the first and only component to the given value.
Definition: colvar.cpp:2176
int write_acf(std::ostream &os)
Save the ACF to a file.
Definition: colvar.cpp:2866
int calc_cvc_gradients(int first, size_t num_cvcs)
Same as colvar::calc_cvc_values but for gradients.
Definition: colvar.cpp:1532
virtual const std::vector< feature * > & features() const
Implementation of the feature list accessor for colvar.
Definition: colvar.h:99
colvarvalue x_ext
Restraint center.
Definition: colvar.h:178
int set_cvc_flags(std::vector< bool > const &flags)
Enables and disables individual CVCs based on the given array.
Definition: colvar.cpp:2050
int analyze()
Perform analysis tasks.
Definition: colvar.cpp:2670
std::vector< int > atom_ids
Sorted array of (zero-based) IDs for all atoms involved.
Definition: colvar.h:702
void setup()
Get ready for a run and re-initialize internal data if needed.
Definition: colvar.cpp:1238
Parent class for a member object of a bias, cv or cvc etc. containing features and their dependencies...
Definition: colvardeps.h:34
void check_enabled(int f, std::string const &reason) const
Check that a feature is enabled, raising COLVARS_BUG_ERROR if not.
Definition: colvardeps.h:434
@ f_cv_extended_Lagrangian
The variable has a harmonic restraint around a moving center with fictitious mass; bias forces will b...
Definition: colvardeps.h:291
@ f_cv_gradient
Gradients are calculated and temporarily stored, so that external forces can be applied.
Definition: colvardeps.h:268
double real
Defining an abstract real number allows to switch precision.
Definition: colvarmodule.h:95
static real dt()
Time step of MD integrator (fs)
static void log(std::string const &message, int min_log_level=10)
Definition: colvarmodule.cpp:1970
static bool debug()
Whether debug output should be enabled (compile-time option)
Definition: colvarmodule.h:330
static std::string to_str(char const *s)
Convert to string for output purposes.
Definition: colvarmodule.cpp:2393
long long step_number
Use a 64-bit integer to store the step number.
Definition: colvarmodule.h:92
Base class containing parsing functions; all objects which need to parse input inherit from this.
Definition: colvarparse.h:27
Value of a collective variable: this is a metatype which can be set at runtime. By default it is set ...
Definition: colvarvalue.h:43
Type type() const
Get the current type.
Definition: colvarvalue.h:154
void reset()
Set to the null value for the data type currently defined.
Definition: colvarvalue.cpp:178
size_t size() const
Number of dimensions of this variable.
Definition: colvarvalue.h:373
Definition: colvars_memstream.h:30
Collective variables main module.
Parsing functions for collective variables.