00001
00002
00003 #ifndef COLVAR_H
00004 #define COLVAR_H
00005
00006 #include <iostream>
00007 #include <iomanip>
00008 #include <cmath>
00009
00010 #include "colvarmodule.h"
00011 #include "colvarvalue.h"
00012 #include "colvarparse.h"
00013
00014
00040
00041 class colvar : public colvarparse {
00042
00043 public:
00044
00046 std::string name;
00047
00049 colvarvalue::Type type() const;
00050
00052 colvarvalue const & value() const;
00053
00055 colvarvalue const & velocity() const;
00056
00065 colvarvalue const & system_force() const;
00066
00068
00077 cvm::real width;
00078
00081 bool b_linear;
00082
00085 bool b_inverse_gradients;
00086
00089 bool b_Jacobian_force;
00090
00093 enum task {
00096 task_gradients,
00098 task_fdiff_velocity,
00101 task_system_force,
00105 task_extended_lagrangian,
00109 task_Jacobian_force,
00112 task_report_Jacobian_force,
00114 task_output_value,
00116 task_output_velocity,
00118 task_output_applied_force,
00120 task_output_system_force,
00122 task_lower_boundary,
00124 task_upper_boundary,
00128 task_grid,
00131 task_lower_wall,
00134 task_upper_wall,
00136 task_runave,
00138 task_corrfunc,
00140 task_ntot
00141 };
00142
00144 bool tasks[task_ntot];
00145
00146 protected:
00147
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00171 colvarvalue x;
00172
00174 colvarvalue x_reported;
00175
00177 colvarvalue v_fdiff;
00178
00179 inline colvarvalue fdiff_velocity (colvarvalue const &xold,
00180 colvarvalue const &xnew)
00181 {
00182
00183
00184
00185 return ( ( (cvm::dt > 0.0) ? (1.0/cvm::dt) : 1.0 ) *
00186 0.5 * dist2_lgrad (xnew, xold) );
00187 }
00188
00190 colvarvalue v_reported;
00191
00192
00194 colvarvalue xr;
00196 colvarvalue vr;
00198 cvm::real ext_mass;
00200 cvm::real ext_force_k;
00202 colvarvalue fr;
00203
00205 colvarvalue fj;
00206
00208 colvarvalue ft_reported;
00209
00210 public:
00211
00212
00215 colvarvalue fb;
00216
00220 colvarvalue f;
00221
00224 colvarvalue ft;
00225
00226
00228 cvm::real period;
00229
00232 bool b_periodic;
00233
00234
00237 bool expand_boundaries;
00238
00240 colvarvalue lower_boundary;
00242 colvarvalue lower_wall;
00244 cvm::real lower_wall_k;
00245
00247 colvarvalue upper_boundary;
00249 colvarvalue upper_wall;
00251 cvm::real upper_wall_k;
00252
00253
00254
00255
00258 int current_bin_scalar() const;
00259
00262 int value_to_bin_scalar (colvarvalue const &value) const;
00263
00265 int value_to_bin_scalar (colvarvalue const &value,
00266 colvarvalue const &new_offset,
00267 cvm::real const &new_width) const;
00268
00271 colvarvalue bin_to_value_scalar (int const &i_bin) const;
00272
00274 colvarvalue bin_to_value_scalar (int const &i_bin,
00275 colvarvalue const &new_offset,
00276 cvm::real const &new_width) const;
00277
00278
00280 bool periodic_boundaries() const;
00281
00283 bool periodic_boundaries (colvarvalue const &lb, colvarvalue const &ub) const;
00284
00285
00287 colvar (std::string const &conf);
00288
00290 void enable (colvar::task const &t);
00291
00293 void disable (colvar::task const &t);
00294
00296 ~colvar();
00297
00298
00301 void calc();
00302
00304 void calc_value (colvarvalue &xn);
00305
00307 void reset_bias_force();
00308
00310 void add_bias_force (colvarvalue const &force);
00311
00315 void update();
00316
00319 void communicate_forces();
00320
00321
00326 cvm::real dist2 (colvarvalue const &x1,
00327 colvarvalue const &x2) const;
00328
00333 colvarvalue dist2_lgrad (colvarvalue const &x1,
00334 colvarvalue const &x2) const;
00335
00340 colvarvalue dist2_rgrad (colvarvalue const &x1,
00341 colvarvalue const &x2) const;
00342
00347 cvm::real compare (colvarvalue const &x1,
00348 colvarvalue const &x2) const;
00349
00350
00351
00352
00354 void parse_analysis (std::string const &conf);
00356 void analyse();
00357
00358
00360 std::istream & read_traj (std::istream &is);
00362 std::ostream & write_traj (std::ostream &os);
00364 std::ostream & write_traj_label (std::ostream &os);
00365
00366
00368 std::istream & read_restart (std::istream &is);
00370 std::ostream & write_restart (std::ostream &os);
00371
00372
00373 protected:
00374
00376 colvarvalue x_old;
00377
00380 std::list< std::list<colvarvalue> > acf_x_history, acf_v_history;
00383 std::list< std::list<colvarvalue> >::iterator acf_x_history_p, acf_v_history_p;
00384
00386 std::list< std::list<colvarvalue> > x_history;
00389 std::list< std::list<colvarvalue> >::iterator x_history_p;
00390
00393 std::string acf_colvar_name;
00395 size_t acf_length;
00397 size_t acf_offset;
00399 size_t acf_stride;
00401 size_t acf_nframes;
00403 bool acf_normalize;
00405 std::vector<cvm::real> acf;
00407 std::string acf_outfile;
00408
00410 enum acf_type_e {
00412 acf_notset,
00414 acf_vel,
00416 acf_coor,
00419 acf_p2coor
00420 };
00421
00423 acf_type_e acf_type;
00424
00426 void calc_vel_acf (std::list<colvarvalue> &v_history,
00427 colvarvalue const &v);
00428
00431 void calc_coor_acf (std::list<colvarvalue> &x_history,
00432 colvarvalue const &x);
00433
00436 void calc_p2coor_acf (std::list<colvarvalue> &x_history,
00437 colvarvalue const &x);
00438
00440 void calc_acf();
00442 void write_acf (std::ostream &os);
00443
00445 size_t runave_length;
00447 size_t runave_stride;
00449 std::ofstream runave_os;
00451 colvarvalue runave;
00453 cvm::real runave_variance;
00454
00456 void calc_runave();
00457
00458
00459 public:
00460
00461
00462
00463 class cvc;
00464
00465
00466
00467
00468 class distance;
00469 class distance_z;
00470 class distance_xy;
00471 class min_distance;
00472 class angle;
00473 class dihedral;
00474 class coordnum;
00475 class h_bond;
00476 class rmsd;
00477 class logmsd;
00478 class orientation_angle;
00479 class gyration;
00480 class eigenvector;
00481 class alpha_dihedrals;
00482 class alpha_angles;
00483
00484
00485
00486 class distance_dir;
00487 class orientation;
00488
00489 protected:
00490
00492 std::vector<cvc *> cvcs;
00493
00494 public:
00495
00496 inline size_t n_components () const {
00497 return cvcs.size();
00498 }
00499 };
00500
00501
00502 inline colvar * cvm::colvar_p (std::string const &name)
00503 {
00504 for (std::vector<colvar *>::iterator cvi = cvm::colvars.begin();
00505 cvi != cvm::colvars.end();
00506 cvi++) {
00507 if ((*cvi)->name == name) {
00508 return (*cvi);
00509 }
00510 }
00511 return NULL;
00512 }
00513
00514
00515 inline colvarvalue::Type colvar::type() const
00516 {
00517 return x.type();
00518 }
00519
00520
00521 inline colvarvalue const & colvar::value() const
00522 {
00523 return x_reported;
00524 }
00525
00526
00527 inline colvarvalue const & colvar::velocity() const
00528 {
00529 return v_reported;
00530 }
00531
00532
00533 inline colvarvalue const & colvar::system_force() const
00534 {
00535 return ft_reported;
00536 }
00537
00538
00539
00540 inline int colvar::current_bin_scalar() const
00541 {
00542 return this->value_to_bin_scalar (this->value());
00543 }
00544
00545 inline int colvar::value_to_bin_scalar (colvarvalue const &value) const
00546 {
00547 return (int) ::floor ( (value.real_value - lower_boundary.real_value) / width );
00548 }
00549
00550 inline int colvar::value_to_bin_scalar (colvarvalue const &value,
00551 colvarvalue const &new_offset,
00552 cvm::real const &new_width) const
00553 {
00554 return (int) ::floor ( (value.real_value - new_offset.real_value) / new_width );
00555 }
00556
00557
00558 inline colvarvalue colvar::bin_to_value_scalar (int const &i_bin) const
00559 {
00560 return lower_boundary.real_value + width * (0.5 + i_bin);
00561 }
00562
00563 inline colvarvalue colvar::bin_to_value_scalar (int const &i_bin,
00564 colvarvalue const &new_offset,
00565 cvm::real const &new_width) const
00566 {
00567 return new_offset.real_value + new_width * (0.5 + i_bin);
00568 }
00569
00570
00571 inline void colvar::add_bias_force (colvarvalue const &force)
00572 {
00573 fb += force;
00574 }
00575
00576 inline void colvar::reset_bias_force() {
00577 fb.reset();
00578 }
00579
00580 #endif
00581