~ubuntu-branches/debian/sid/lammps/sid

« back to all changes in this revision

Viewing changes to lib/colvars/colvar.h

  • Committer: Package Import Robot
  • Author(s): Anton Gladky
  • Date: 2015-04-29 23:44:49 UTC
  • mfrom: (5.1.3 experimental)
  • Revision ID: package-import@ubuntu.com-20150429234449-mbhy9utku6hp6oq8
Tags: 0~20150313.gitfa668e1-1
Upload into unstable.

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
// -*- c++ -*-
 
1
/// -*- c++ -*-
2
2
 
3
3
#ifndef COLVAR_H
4
4
#define COLVAR_H
31
31
///
32
32
/// Please note that most of its members are \link colvarvalue
33
33
/// \endlink objects, i.e. they can handle different data types
34
 
/// together, and must all be set to the same type of colvar::x by
35
 
/// using the colvarvalue::type() member function before using them
36
 
/// together in assignments or other operations; this is usually done
 
34
/// together, and must all be set to the same type of colvar::value()
 
35
/// before using them together in assignments or other operations; this is usually done
37
36
/// automatically in the constructor.  If you add a new member of
38
37
/// \link colvarvalue \endlink type, you should also add its
39
38
/// initialization line in the \link colvar \endlink constructor.
45
44
  /// Name
46
45
  std::string name;
47
46
 
48
 
  /// Type of value
49
 
  colvarvalue::Type type() const;
50
 
 
51
 
  /// \brief Current value (previously obtained from calc() or read_traj())
 
47
  /// \brief Current value (previously set by calc() or by read_traj())
52
48
  colvarvalue const & value() const;
53
49
 
54
50
  /// \brief Current actual value (not extended DOF)
55
51
  colvarvalue const & actual_value() const;
56
52
 
57
 
  /// \brief Current velocity (previously obtained from calc() or read_traj())
 
53
  /// \brief Current velocity (previously set by calc() or by read_traj())
58
54
  colvarvalue const & velocity() const;
59
55
 
60
56
  /// \brief Current system force (previously obtained from calc() or
83
79
  /// combination of \link cvc \endlink elements
84
80
  bool b_linear;
85
81
 
 
82
  /// \brief True if this \link colvar \endlink is a linear
 
83
  /// combination of cvcs with coefficients 1 or -1
 
84
  bool b_homogeneous;
 
85
 
86
86
  /// \brief True if all \link cvc \endlink objects are capable
87
87
  /// of calculating inverse gradients
88
88
  bool b_inverse_gradients;
148
148
    task_runave,
149
149
    /// \brief Compute time correlation function
150
150
    task_corrfunc,
 
151
    /// \brief Value and gradients computed by user script
 
152
    task_scripted,
151
153
    /// \brief Number of possible tasks
152
154
    task_ntot
153
155
  };
155
157
  /// Tasks performed by this colvar
156
158
  bool tasks[task_ntot];
157
159
 
 
160
  /// List of biases that depend on this colvar
 
161
  std::vector<colvarbias *> biases;
158
162
protected:
159
163
 
160
164
 
162
166
    extended:
163
167
    H = H_{0} + \sum_{i} 1/2*\lambda*(S_i(x(t))-s_i(t))^2 \\
164
168
    + \sum_{i} 1/2*m_i*(ds_i(t)/dt)^2 \\
165
 
    + \sum_{t'<t} W * exp (-1/2*\sum_{i} (s_i(t')-s_i(t))^2/(\delta{}s_i)^2) \\
 
169
    + \sum_{t'<t} W * exp(-1/2*\sum_{i} (s_i(t')-s_i(t))^2/(\delta{}s_i)^2) \\
166
170
    + \sum_{w} (\sum_{i}c_{w,i}s_i(t) - D_w)^M/(\Sigma_w)^M
167
171
 
168
172
    normal:
169
 
    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) \\
 
173
    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) \\
170
174
    + \sum_{w} (\sum_{i}c_{w,i}S_i(t) - D_w)^M/(\Sigma_w)^M
171
175
 
172
176
    output:
182
186
  /// Value of the colvar
183
187
  colvarvalue x;
184
188
 
 
189
  // TODO: implement functionality to treat these
 
190
  // /// Vector of individual values from CVCs
 
191
  // colvarvalue x_cvc;
 
192
 
 
193
  // /// Jacobian matrix of individual values from CVCs
 
194
  // colvarvalue dx_cvc;
 
195
 
185
196
  /// Cached reported value (x may be manipulated)
186
197
  colvarvalue x_reported;
187
198
 
188
199
  /// Finite-difference velocity
189
200
  colvarvalue v_fdiff;
190
201
 
191
 
  inline colvarvalue fdiff_velocity (colvarvalue const &xold,
 
202
  inline colvarvalue fdiff_velocity(colvarvalue const &xold,
192
203
                                     colvarvalue const &xnew)
193
204
  {
194
205
    // using the gradient of the square distance to calculate the
196
207
    // account)
197
208
    cvm::real dt = cvm::dt();
198
209
    return ( ( (dt > 0.0) ? (1.0/dt) : 1.0 ) *
199
 
             0.5 * dist2_lgrad (xnew, xold) );
 
210
             0.5 * dist2_lgrad(xnew, xold) );
200
211
  }
201
212
 
202
213
  /// Cached reported velocity
276
287
  bool periodic_boundaries() const;
277
288
 
278
289
  /// \brief Is the interval defined by the two boundaries periodic?
279
 
  bool periodic_boundaries (colvarvalue const &lb, colvarvalue const &ub) const;
 
290
  bool periodic_boundaries(colvarvalue const &lb, colvarvalue const &ub) const;
280
291
 
281
292
 
282
293
  /// Constructor
283
 
  colvar (std::string const &conf);
 
294
  colvar(std::string const &conf);
284
295
 
285
296
  /// Enable the specified task
286
 
  void enable (colvar::task const &t);
 
297
  int enable(colvar::task const &t);
287
298
 
288
299
  /// Disable the specified task
289
 
  void disable (colvar::task const &t);
 
300
  void disable(colvar::task const &t);
290
301
 
291
302
  /// Get ready for a run and possibly re-initialize internal data
292
303
  void setup();
300
311
  void calc();
301
312
 
302
313
  /// \brief Calculate just the value, and store it in the argument
303
 
  void calc_value (colvarvalue &xn);
 
314
  void calc_value(colvarvalue &xn);
304
315
 
305
316
  /// Set the total biasing force to zero
306
317
  void reset_bias_force();
307
318
 
308
319
  /// Add to the total force from biases
309
 
  void add_bias_force (colvarvalue const &force);
 
320
  void add_bias_force(colvarvalue const &force);
310
321
 
311
322
  /// \brief Collect all forces on this colvar, integrate internal
312
323
  /// equations of motion of internal degrees of freedom; see also
323
334
  /// \endlink objects) to calculate square distances and gradients
324
335
  ///
325
336
  /// Handles correctly symmetries and periodic boundary conditions
326
 
  cvm::real dist2 (colvarvalue const &x1,
 
337
  cvm::real dist2(colvarvalue const &x1,
327
338
                   colvarvalue const &x2) const;
328
339
 
329
340
  /// \brief Use the internal metrics (as from \link cvc
330
341
  /// \endlink objects) to calculate square distances and gradients
331
342
  ///
332
343
  /// Handles correctly symmetries and periodic boundary conditions
333
 
  colvarvalue dist2_lgrad (colvarvalue const &x1,
 
344
  colvarvalue dist2_lgrad(colvarvalue const &x1,
334
345
                           colvarvalue const &x2) const;
335
346
 
336
347
  /// \brief Use the internal metrics (as from \link cvc
337
348
  /// \endlink objects) to calculate square distances and gradients
338
349
  ///
339
350
  /// Handles correctly symmetries and periodic boundary conditions
340
 
  colvarvalue dist2_rgrad (colvarvalue const &x1,
 
351
  colvarvalue dist2_rgrad(colvarvalue const &x1,
341
352
                           colvarvalue const &x2) const;
342
353
 
343
354
  /// \brief Use the internal metrics (as from \link cvc
344
 
  /// \endlink objects) to compare colvar values
345
 
  ///
346
 
  /// Handles correctly symmetries and periodic boundary conditions
347
 
  cvm::real compare (colvarvalue const &x1,
348
 
                     colvarvalue const &x2) const;
349
 
 
350
 
  /// \brief Use the internal metrics (as from \link cvc
351
355
  /// \endlink objects) to wrap a value into a standard interval
352
356
  ///
353
357
  /// Handles correctly symmetries and periodic boundary conditions
354
 
  void wrap (colvarvalue &x) const;
 
358
  void wrap(colvarvalue &x) const;
355
359
 
356
360
 
357
361
  /// Read the analysis tasks
358
 
  void parse_analysis (std::string const &conf);
 
362
  int parse_analysis(std::string const &conf);
359
363
  /// Perform analysis tasks
360
364
  void analyse();
361
365
 
362
366
 
363
367
  /// Read the value from a collective variable trajectory file
364
 
  std::istream & read_traj (std::istream &is);
 
368
  std::istream & read_traj(std::istream &is);
365
369
  /// Output formatted values to the trajectory file
366
 
  std::ostream & write_traj (std::ostream &os);
 
370
  std::ostream & write_traj(std::ostream &os);
367
371
  /// Write a label to the trajectory file (comment line)
368
 
  std::ostream & write_traj_label (std::ostream &os);
 
372
  std::ostream & write_traj_label(std::ostream &os);
369
373
 
370
374
 
371
375
  /// Read the collective variable from a restart file
372
 
  std::istream & read_restart (std::istream &is);
 
376
  std::istream & read_restart(std::istream &is);
373
377
  /// Write the collective variable to a restart file
374
 
  std::ostream & write_restart (std::ostream &os);
 
378
  std::ostream & write_restart(std::ostream &os);
375
379
 
376
380
  /// Write output files (if defined, e.g. in analysis mode)
377
 
  void write_output_files();
 
381
  int write_output_files();
378
382
 
379
383
 
380
384
protected:
430
434
  acf_type_e             acf_type;
431
435
 
432
436
  /// \brief Velocity ACF, scalar product between v(0) and v(t)
433
 
  void calc_vel_acf (std::list<colvarvalue> &v_history,
 
437
  int calc_vel_acf(std::list<colvarvalue> &v_history,
434
438
                     colvarvalue const      &v);
435
439
 
436
440
  /// \brief Coordinate ACF, scalar product between x(0) and x(t)
437
441
  /// (does not work with scalar numbers)
438
 
  void calc_coor_acf (std::list<colvarvalue> &x_history,
 
442
  void calc_coor_acf(std::list<colvarvalue> &x_history,
439
443
                      colvarvalue const      &x);
440
444
 
441
445
  /// \brief Coordinate ACF, second order Legendre polynomial between
442
446
  /// x(0) and x(t) (does not work with scalar numbers)
443
 
  void calc_p2coor_acf (std::list<colvarvalue> &x_history,
 
447
  void calc_p2coor_acf(std::list<colvarvalue> &x_history,
444
448
                        colvarvalue const      &x);
445
449
 
446
450
  /// Calculate the auto-correlation function (ACF)
447
 
  void calc_acf();
 
451
  int calc_acf();
448
452
  /// Save the ACF to a file
449
 
  void write_acf (std::ostream &os);
 
453
  void write_acf(std::ostream &os);
450
454
 
451
455
  /// Length of running average series
452
456
  size_t         runave_length;
453
457
  /// Timesteps to skip between two values in the running average series
454
458
  size_t         runave_stride;
455
459
  /// Name of the file to write the running average
456
 
  std::ofstream  runave_os;
 
460
  cvm::ofstream  runave_os;
457
461
  /// Current value of the running average
458
462
  colvarvalue    runave;
459
463
  /// Current value of the square deviation from the running average
478
482
  class distance_z;
479
483
  class distance_xy;
480
484
  class distance_inv;
 
485
  class distance_pairs;
481
486
  class angle;
482
487
  class dihedral;
483
488
  class coordnum;
485
490
  class h_bond;
486
491
  class rmsd;
487
492
  class orientation_angle;
 
493
  class orientation_proj;
488
494
  class tilt;
489
495
  class spin_angle;
490
496
  class gyration;
498
504
  // non-scalar components
499
505
  class distance_vec;
500
506
  class distance_dir;
 
507
  class cartesian;
501
508
  class orientation;
502
509
 
503
510
protected:
507
514
 
508
515
  /// \brief Initialize the sorted list of atom IDs for atoms involved
509
516
  /// in all cvcs (called when enabling task_collect_gradients)
510
 
  void build_atom_list (void);
 
517
  void build_atom_list(void);
 
518
 
 
519
private:
 
520
  /// Name of scripted function to be used
 
521
  std::string scripted_function;
 
522
 
 
523
  /// Current cvc values in the order requested by script
 
524
  /// when using scriptedFunction
 
525
  std::vector<const colvarvalue *> sorted_cvc_values;
511
526
 
512
527
public:
513
528
  /// \brief Sorted array of (zero-based) IDs for all atoms involved
518
533
  /// For scalar variables only!
519
534
  std::vector<cvm::rvector> atomic_gradients;
520
535
 
521
 
  inline size_t n_components () const {
 
536
  inline size_t n_components() const {
522
537
    return cvcs.size();
523
538
  }
524
539
};
525
540
 
526
541
 
527
 
inline colvar * cvm::colvar_p (std::string const &name)
528
 
{
529
 
  for (std::vector<colvar *>::iterator cvi = cvm::colvars.begin();
530
 
       cvi != cvm::colvars.end();
531
 
       cvi++) {
532
 
    if ((*cvi)->name == name) {
533
 
      return (*cvi);
534
 
    }
535
 
  }
536
 
  return NULL;
537
 
}
538
 
 
539
 
 
540
 
inline colvarvalue::Type colvar::type() const
541
 
{
542
 
  return x.type();
543
 
}
544
 
 
545
 
 
546
542
inline colvarvalue const & colvar::value() const
547
543
{
548
544
  return x_reported;
567
563
}
568
564
 
569
565
 
570
 
inline void colvar::add_bias_force (colvarvalue const &force)
 
566
inline void colvar::add_bias_force(colvarvalue const &force)
571
567
{
572
568
  fb += force;
573
569
}
574
570
 
575
571
 
576
572
inline void colvar::reset_bias_force() {
 
573
  fb.type(value());
577
574
  fb.reset();
578
575
}
579
576