Collective Variables Module - Developer Documentation
colvargrid.h
1 // -*- c++ -*-
2 
3 #ifndef COLVARGRID_H
4 #define COLVARGRID_H
5 
6 #include <iostream>
7 #include <iomanip>
8 #include <cmath>
9 
10 #include "colvar.h"
11 #include "colvarmodule.h"
12 #include "colvarvalue.h"
13 #include "colvarparse.h"
14 
19 template <class T> class colvar_grid : public colvarparse {
20 
21 protected:
22 
24  size_t nd;
25 
27  std::vector<int> nx;
28 
30  std::vector<int> nxc;
31 
34  size_t mult;
35 
37  size_t nt;
38 
40  std::vector<T> data;
41 
43  std::vector<size_t> new_data;
44 
46  std::vector<colvar *> cv;
47 
49  std::vector<bool> actual_value;
50 
52  inline size_t address(std::vector<int> const &ix) const
53  {
54  size_t addr = 0;
55  for (size_t i = 0; i < nd; i++) {
56  addr += ix[i]*nxc[i];
57  if (cvm::debug()) {
58  if (ix[i] >= nx[i]) {
59  cvm::error("Error: exceeding bounds in colvar_grid.\n", BUG_ERROR);
60  return 0;
61  }
62  }
63  }
64  return addr;
65  }
66 
67 public:
68 
70  std::vector<colvarvalue> lower_boundaries;
71 
73  std::vector<colvarvalue> upper_boundaries;
74 
76  std::vector<bool> periodic;
77 
79  std::vector<bool> hard_lower_boundaries;
80 
82  std::vector<bool> hard_upper_boundaries;
83 
85  std::vector<cvm::real> widths;
86 
89 
91  bool has_data;
92 
94  inline size_t number_of_colvars() const
95  {
96  return nd;
97  }
98 
101  inline size_t number_of_points(int const icv = -1) const
102  {
103  if (icv < 0) {
104  return nt;
105  } else {
106  return nx[icv];
107  }
108  }
109 
111  inline std::vector<int> const & sizes() const
112  {
113  return nx;
114  }
115 
117  inline void set_sizes(std::vector<int> const &new_sizes)
118  {
119  nx = new_sizes;
120  }
121 
123  inline size_t multiplicity() const
124  {
125  return mult;
126  }
127 
129  inline void request_actual_value(bool b = true)
130  {
131  size_t i;
132  for (i = 0; i < actual_value.size(); i++) {
133  actual_value[i] = b;
134  }
135  }
136 
138  int setup(std::vector<int> const &nx_i,
139  T const &t = T(),
140  size_t const &mult_i = 1)
141  {
142  if (cvm::debug()) {
143  cvm::log("Allocating grid: multiplicity = "+cvm::to_str(mult_i)+
144  ", dimensions = "+cvm::to_str(nx_i)+".\n");
145  }
146 
147  mult = mult_i;
148 
149  data.clear();
150 
151  nx = nx_i;
152  nd = nx.size();
153 
154  nxc.resize(nd);
155 
156  // setup dimensions
157  nt = mult;
158  for (int i = nd-1; i >= 0; i--) {
159  if (nx[i] <= 0) {
160  cvm::error("Error: providing an invalid number of grid points, "+
161  cvm::to_str(nx[i])+".\n", BUG_ERROR);
162  return COLVARS_ERROR;
163  }
164  nxc[i] = nt;
165  nt *= nx[i];
166  }
167 
168  if (cvm::debug()) {
169  cvm::log("Total number of grid elements = "+cvm::to_str(nt)+".\n");
170  }
171 
172  data.reserve(nt);
173  data.assign(nt, t);
174 
175  return COLVARS_OK;
176  }
177 
179  int setup()
180  {
181  return setup(this->nx, T(), this->mult);
182  }
183 
185  void reset(T const &t = T())
186  {
187  data.assign(nt, t);
188  }
189 
190 
192  colvar_grid() : has_data(false)
193  {
194  nd = nt = 0;
195  mult = 1;
196  this->setup();
197  }
198 
200  virtual ~colvar_grid()
201  {}
202 
206  colvar_grid(colvar_grid<T> const &g) : nd(g.nd),
207  nx(g.nx),
208  mult(g.mult),
209  data(),
210  cv(g.cv),
211  actual_value(g.actual_value),
212  lower_boundaries(g.lower_boundaries),
213  upper_boundaries(g.upper_boundaries),
214  periodic(g.periodic),
215  hard_lower_boundaries(g.hard_lower_boundaries),
216  hard_upper_boundaries(g.hard_upper_boundaries),
217  widths(g.widths),
218  has_data(false)
219  {
220  }
221 
226  colvar_grid(std::vector<int> const &nx_i,
227  T const &t = T(),
228  size_t mult_i = 1)
229  : has_data(false)
230  {
231  this->setup(nx_i, t, mult_i);
232  }
233 
235  colvar_grid(std::vector<colvar *> const &colvars,
236  T const &t = T(),
237  size_t mult_i = 1,
238  bool margin = false)
239  : has_data(false)
240  {
241  this->init_from_colvars(colvars, t, mult_i, margin);
242  }
243 
244  int init_from_colvars(std::vector<colvar *> const &colvars,
245  T const &t = T(),
246  size_t mult_i = 1,
247  bool margin = false)
248  {
249  if (cvm::debug()) {
250  cvm::log("Reading grid configuration from collective variables.\n");
251  }
252 
253  cv = colvars;
254  nd = colvars.size();
255  mult = mult_i;
256 
257  size_t i;
258 
259  if (cvm::debug()) {
260  cvm::log("Allocating a grid for "+cvm::to_str(colvars.size())+
261  " collective variables, multiplicity = "+cvm::to_str(mult_i)+".\n");
262  }
263 
264  for (i = 0; i < cv.size(); i++) {
265 
266  if (cv[i]->value().type() != colvarvalue::type_scalar) {
267  cvm::error("Colvar grids can only be automatically "
268  "constructed for scalar variables. "
269  "ABF and histogram can not be used; metadynamics "
270  "can be used with useGrids disabled.\n", INPUT_ERROR);
271  return COLVARS_ERROR;
272  }
273 
274  if (cv[i]->width <= 0.0) {
275  cvm::error("Tried to initialize a grid on a "
276  "variable with negative or zero width.\n", INPUT_ERROR);
277  return COLVARS_ERROR;
278  }
279 
280  widths.push_back(cv[i]->width);
281  hard_lower_boundaries.push_back(cv[i]->hard_lower_boundary);
282  hard_upper_boundaries.push_back(cv[i]->hard_upper_boundary);
283  periodic.push_back(cv[i]->periodic_boundaries());
284 
285  // By default, get reported colvar value (for extended Lagrangian colvars)
286  actual_value.push_back(false);
287 
288  // except if a colvar is specified twice in a row
289  // then the first instance is the actual value
290  // For histograms of extended-system coordinates
291  if (i > 0 && cv[i-1] == cv[i]) {
292  actual_value[i-1] = true;
293  }
294 
295  if (margin) {
296  if (periodic[i]) {
297  // Shift the grid by half the bin width (values at edges instead of center of bins)
298  lower_boundaries.push_back(cv[i]->lower_boundary.real_value - 0.5 * widths[i]);
299  upper_boundaries.push_back(cv[i]->upper_boundary.real_value - 0.5 * widths[i]);
300  } else {
301  // Make this grid larger by one bin width
302  lower_boundaries.push_back(cv[i]->lower_boundary.real_value - 0.5 * widths[i]);
303  upper_boundaries.push_back(cv[i]->upper_boundary.real_value + 0.5 * widths[i]);
304  }
305  } else {
306  lower_boundaries.push_back(cv[i]->lower_boundary);
307  upper_boundaries.push_back(cv[i]->upper_boundary);
308  }
309  }
310 
311  this->init_from_boundaries();
312  return this->setup();
313  }
314 
315  int init_from_boundaries(T const &t = T(),
316  size_t const &mult_i = 1)
317  {
318  if (cvm::debug()) {
319  cvm::log("Configuring grid dimensions from colvars boundaries.\n");
320  }
321 
322  // these will have to be updated
323  nx.clear();
324  nxc.clear();
325  nt = 0;
326 
327  for (size_t i = 0; i < lower_boundaries.size(); i++) {
328 
329  cvm::real nbins = ( upper_boundaries[i].real_value -
330  lower_boundaries[i].real_value ) / widths[i];
331  int nbins_round = (int)(nbins+0.5);
332 
333  if (std::fabs(nbins - cvm::real(nbins_round)) > 1.0E-10) {
334  cvm::log("Warning: grid interval("+
335  cvm::to_str(lower_boundaries[i], cvm::cv_width, cvm::cv_prec)+" - "+
336  cvm::to_str(upper_boundaries[i], cvm::cv_width, cvm::cv_prec)+
337  ") is not commensurate to its bin width("+
338  cvm::to_str(widths[i], cvm::cv_width, cvm::cv_prec)+").\n");
339  upper_boundaries[i].real_value = lower_boundaries[i].real_value +
340  (nbins_round * widths[i]);
341  }
342 
343  if (cvm::debug())
344  cvm::log("Number of points is "+cvm::to_str((int) nbins_round)+
345  " for the colvar no. "+cvm::to_str(i+1)+".\n");
346 
347  nx.push_back(nbins_round);
348  }
349 
350  return COLVARS_OK;
351  }
352 
355  inline void wrap(std::vector<int> & ix) const
356  {
357  for (size_t i = 0; i < nd; i++) {
358  if (periodic[i]) {
359  ix[i] = (ix[i] + nx[i]) % nx[i]; //to ensure non-negative result
360  } else {
361  if (ix[i] < 0 || ix[i] >= nx[i]) {
362  cvm::error("Trying to wrap illegal index vector (non-PBC) for a grid point: "
363  + cvm::to_str(ix), BUG_ERROR);
364  return;
365  }
366  }
367  }
368  }
369 
370 
372  inline int current_bin_scalar(int const i) const
373  {
374  return value_to_bin_scalar(actual_value[i] ? cv[i]->actual_value() : cv[i]->value(), i);
375  }
376 
379  inline int current_bin_scalar_bound(int const i) const
380  {
381  return value_to_bin_scalar_bound(actual_value[i] ? cv[i]->actual_value() : cv[i]->value(), i);
382  }
383 
385  inline int current_bin_scalar(int const i, int const iv) const
386  {
387  return value_to_bin_scalar(actual_value[i] ?
388  cv[i]->actual_value().vector1d_value[iv] :
389  cv[i]->value().vector1d_value[iv], i);
390  }
391 
394  inline int value_to_bin_scalar(colvarvalue const &value, const int i) const
395  {
396  return (int) std::floor( (value.real_value - lower_boundaries[i].real_value) / widths[i] );
397  }
398 
401  inline int value_to_bin_scalar_bound(colvarvalue const &value, const int i) const
402  {
403  int bin_index = std::floor( (value.real_value - lower_boundaries[i].real_value) / widths[i] );
404  if (bin_index < 0) bin_index=0;
405  if (bin_index >=int(nx[i])) bin_index=int(nx[i])-1;
406  return (int) bin_index;
407  }
408 
411  colvarvalue const &new_offset,
412  cvm::real const &new_width) const
413  {
414  return (int) std::floor( (value.real_value - new_offset.real_value) / new_width );
415  }
416 
419  inline colvarvalue bin_to_value_scalar(int const &i_bin, int const i) const
420  {
421  return lower_boundaries[i].real_value + widths[i] * (0.5 + i_bin);
422  }
423 
425  inline colvarvalue bin_to_value_scalar(int const &i_bin,
426  colvarvalue const &new_offset,
427  cvm::real const &new_width) const
428  {
429  return new_offset.real_value + new_width * (0.5 + i_bin);
430  }
431 
433  inline void set_value(std::vector<int> const &ix,
434  T const &t,
435  size_t const &imult = 0)
436  {
437  data[this->address(ix)+imult] = t;
438  has_data = true;
439  }
440 
445  void delta_grid(colvar_grid<T> const &other_grid)
446  {
447 
448  if (other_grid.multiplicity() != this->multiplicity()) {
449  cvm::error("Error: trying to subtract two grids with "
450  "different multiplicity.\n");
451  return;
452  }
453 
454  if (other_grid.data.size() != this->data.size()) {
455  cvm::error("Error: trying to subtract two grids with "
456  "different size.\n");
457  return;
458  }
459 
460  for (size_t i = 0; i < data.size(); i++) {
461  data[i] = other_grid.data[i] - data[i];
462  }
463  has_data = true;
464  }
465 
469  void copy_grid(colvar_grid<T> const &other_grid)
470  {
471  if (other_grid.multiplicity() != this->multiplicity()) {
472  cvm::error("Error: trying to copy two grids with "
473  "different multiplicity.\n");
474  return;
475  }
476 
477  if (other_grid.data.size() != this->data.size()) {
478  cvm::error("Error: trying to copy two grids with "
479  "different size.\n");
480  return;
481  }
482 
483 
484  for (size_t i = 0; i < data.size(); i++) {
485  data[i] = other_grid.data[i];
486  }
487  has_data = true;
488  }
489 
492  void raw_data_out(T* out_data) const
493  {
494  for (size_t i = 0; i < data.size(); i++) out_data[i] = data[i];
495  }
497  void raw_data_in(const T* in_data)
498  {
499  for (size_t i = 0; i < data.size(); i++) data[i] = in_data[i];
500  has_data = true;
501  }
503  size_t raw_data_num() const { return data.size(); }
504 
505 
508  inline T const & value(std::vector<int> const &ix,
509  size_t const &imult = 0) const
510  {
511  return data[this->address(ix) + imult];
512  }
513 
514 
516  inline void add_constant(T const &t)
517  {
518  for (size_t i = 0; i < nt; i++)
519  data[i] += t;
520  has_data = true;
521  }
522 
524  inline void multiply_constant(cvm::real const &a)
525  {
526  for (size_t i = 0; i < nt; i++)
527  data[i] *= a;
528  }
529 
531  inline void remove_zeros(cvm::real const &a)
532  {
533  for (size_t i = 0; i < nt; i++)
534  if(data[i]==0) data[i] = a;
535  }
536 
537 
540  inline std::vector<int> const get_colvars_index(std::vector<colvarvalue> const &values) const
541  {
542  std::vector<int> index = new_index();
543  for (size_t i = 0; i < nd; i++) {
544  index[i] = value_to_bin_scalar(values[i], i);
545  }
546  return index;
547  }
548 
551  inline std::vector<int> const get_colvars_index() const
552  {
553  std::vector<int> index = new_index();
554  for (size_t i = 0; i < nd; i++) {
555  index[i] = current_bin_scalar(i);
556  }
557  return index;
558  }
559 
562  inline std::vector<int> const get_colvars_index_bound() const
563  {
564  std::vector<int> index = new_index();
565  for (size_t i = 0; i < nd; i++) {
566  index[i] = current_bin_scalar_bound(i);
567  }
568  return index;
569  }
570 
574  inline cvm::real bin_distance_from_boundaries(std::vector<colvarvalue> const &values,
575  bool skip_hard_boundaries = false)
576  {
577  cvm::real minimum = 1.0E+16;
578  for (size_t i = 0; i < nd; i++) {
579 
580  if (periodic[i]) continue;
581 
582  cvm::real dl = std::sqrt(cv[i]->dist2(values[i], lower_boundaries[i])) / widths[i];
583  cvm::real du = std::sqrt(cv[i]->dist2(values[i], upper_boundaries[i])) / widths[i];
584 
585  if (values[i].real_value < lower_boundaries[i])
586  dl *= -1.0;
587  if (values[i].real_value > upper_boundaries[i])
588  du *= -1.0;
589 
590  if ( ((!skip_hard_boundaries) || (!hard_lower_boundaries[i])) && (dl < minimum))
591  minimum = dl;
592  if ( ((!skip_hard_boundaries) || (!hard_upper_boundaries[i])) && (du < minimum))
593  minimum = du;
594  }
595 
596  return minimum;
597  }
598 
599 
604  void map_grid(colvar_grid<T> const &other_grid)
605  {
606  if (other_grid.multiplicity() != this->multiplicity()) {
607  cvm::error("Error: trying to merge two grids with values of "
608  "different multiplicity.\n");
609  return;
610  }
611 
612  std::vector<colvarvalue> const &gb = this->lower_boundaries;
613  std::vector<cvm::real> const &gw = this->widths;
614  std::vector<colvarvalue> const &ogb = other_grid.lower_boundaries;
615  std::vector<cvm::real> const &ogw = other_grid.widths;
616 
617  std::vector<int> ix = this->new_index();
618  std::vector<int> oix = other_grid.new_index();
619 
620  if (cvm::debug())
621  cvm::log("Remapping grid...\n");
622  for ( ; this->index_ok(ix); this->incr(ix)) {
623 
624  for (size_t i = 0; i < nd; i++) {
625  oix[i] =
626  value_to_bin_scalar(bin_to_value_scalar(ix[i], gb[i], gw[i]),
627  ogb[i],
628  ogw[i]);
629  }
630 
631  if (! other_grid.index_ok(oix)) {
632  continue;
633  }
634 
635  for (size_t im = 0; im < mult; im++) {
636  this->set_value(ix, other_grid.value(oix, im), im);
637  }
638  }
639 
640  has_data = true;
641  if (cvm::debug())
642  cvm::log("Remapping done.\n");
643  }
644 
647  void add_grid(colvar_grid<T> const &other_grid,
648  cvm::real scale_factor = 1.0)
649  {
650  if (other_grid.multiplicity() != this->multiplicity()) {
651  cvm::error("Error: trying to sum togetehr two grids with values of "
652  "different multiplicity.\n");
653  return;
654  }
655  if (scale_factor != 1.0)
656  for (size_t i = 0; i < data.size(); i++) {
657  data[i] += scale_factor * other_grid.data[i];
658  }
659  else
660  // skip multiplication if possible
661  for (size_t i = 0; i < data.size(); i++) {
662  data[i] += other_grid.data[i];
663  }
664  has_data = true;
665  }
666 
669  virtual inline T value_output(std::vector<int> const &ix,
670  size_t const &imult = 0)
671  {
672  return value(ix, imult);
673  }
674 
678  virtual inline void value_input(std::vector<int> const &ix,
679  T const &t,
680  size_t const &imult = 0,
681  bool add = false)
682  {
683  if ( add )
684  data[address(ix) + imult] += t;
685  else
686  data[address(ix) + imult] = t;
687  has_data = true;
688  }
689 
690  // /// Get the pointer to the binned value indexed by ix
691  // inline T const *value_p (std::vector<int> const &ix)
692  // {
693  // return &(data[address (ix)]);
694  // }
695 
698  inline std::vector<int> const new_index() const
699  {
700  return std::vector<int> (nd, 0);
701  }
702 
705  inline bool index_ok(std::vector<int> const &ix) const
706  {
707  for (size_t i = 0; i < nd; i++) {
708  if ( (ix[i] < 0) || (ix[i] >= int(nx[i])) )
709  return false;
710  }
711  return true;
712  }
713 
716  inline void incr(std::vector<int> &ix) const
717  {
718  for (int i = ix.size()-1; i >= 0; i--) {
719 
720  ix[i]++;
721 
722  if (ix[i] >= nx[i]) {
723 
724  if (i > 0) {
725  ix[i] = 0;
726  continue;
727  } else {
728  // this is the last iteration, a non-valid index is being
729  // set for the outer index, which will be caught by
730  // index_ok()
731  ix[0] = nx[0];
732  return;
733  }
734  } else {
735  return;
736  }
737  }
738  }
739 
741  std::ostream & write_params(std::ostream &os)
742  {
743  size_t i;
744  os << "grid_parameters {\n n_colvars " << nd << "\n";
745 
746  os << " lower_boundaries ";
747  for (i = 0; i < nd; i++)
748  os << " " << lower_boundaries[i];
749  os << "\n";
750 
751  os << " upper_boundaries ";
752  for (i = 0; i < nd; i++)
753  os << " " << upper_boundaries[i];
754  os << "\n";
755 
756  os << " widths ";
757  for (i = 0; i < nd; i++)
758  os << " " << widths[i];
759  os << "\n";
760 
761  os << " sizes ";
762  for (i = 0; i < nd; i++)
763  os << " " << nx[i];
764  os << "\n";
765 
766  os << "}\n";
767  return os;
768  }
769 
771  int parse_params(std::string const &conf,
773  {
774  if (cvm::debug()) cvm::log("Reading grid configuration from string.\n");
775 
776  std::vector<int> old_nx = nx;
777  std::vector<colvarvalue> old_lb = lower_boundaries;
778 
779  {
780  size_t nd_in = 0;
781  // this is only used in state files
782  colvarparse::get_keyval(conf, "n_colvars", nd_in, nd, colvarparse::parse_silent);
783  if (nd_in != nd) {
784  cvm::error("Error: trying to read data for a grid "
785  "that contains a different number of colvars ("+
786  cvm::to_str(nd_in)+") than the grid defined "
787  "in the configuration file("+cvm::to_str(nd)+
788  ").\n");
789  return COLVARS_ERROR;
790  }
791  }
792 
793  // underscore keywords are used in state file
794  colvarparse::get_keyval(conf, "lower_boundaries",
795  lower_boundaries, lower_boundaries, colvarparse::parse_silent);
796  colvarparse::get_keyval(conf, "upper_boundaries",
797  upper_boundaries, upper_boundaries, colvarparse::parse_silent);
798 
799  // camel case keywords are used in config file
800  colvarparse::get_keyval(conf, "lowerBoundaries",
801  lower_boundaries, lower_boundaries, parse_mode);
802  colvarparse::get_keyval(conf, "upperBoundaries",
803  upper_boundaries, upper_boundaries, parse_mode);
804 
805  colvarparse::get_keyval(conf, "widths", widths, widths, parse_mode);
806 
807  // only used in state file
808  colvarparse::get_keyval(conf, "sizes", nx, nx, colvarparse::parse_silent);
809 
810  if (nd < lower_boundaries.size()) nd = lower_boundaries.size();
811 
812  if (! actual_value.size()) actual_value.assign(nd, false);
813  if (! periodic.size()) periodic.assign(nd, false);
814  if (! widths.size()) widths.assign(nd, 1.0);
815 
816  bool new_params = false;
817  if (old_nx.size()) {
818  for (size_t i = 0; i < nd; i++) {
819  if ( (old_nx[i] != nx[i]) ||
820  (std::sqrt(cv[i]->dist2(old_lb[i],
821  lower_boundaries[i])) > 1.0E-10) ) {
822  new_params = true;
823  }
824  }
825  } else {
826  new_params = true;
827  }
828 
829  // reallocate the array in case the grid params have just changed
830  if (new_params) {
831  init_from_boundaries();
832  // data.clear(); // no longer needed: setup calls clear()
833  return this->setup(nx, T(), mult);
834  }
835 
836  return COLVARS_OK;
837  }
838 
843  {
844  for (size_t i = 0; i < nd; i++) {
845  if ( (std::sqrt(cv[i]->dist2(cv[i]->lower_boundary,
846  lower_boundaries[i])) > 1.0E-10) ||
847  (std::sqrt(cv[i]->dist2(cv[i]->upper_boundary,
848  upper_boundaries[i])) > 1.0E-10) ||
849  (std::sqrt(cv[i]->dist2(cv[i]->width,
850  widths[i])) > 1.0E-10) ) {
851  cvm::error("Error: restart information for a grid is "
852  "inconsistent with that of its colvars.\n");
853  return;
854  }
855  }
856  }
857 
858 
861  void check_consistency(colvar_grid<T> const &other_grid)
862  {
863  for (size_t i = 0; i < nd; i++) {
864  // we skip dist2(), because periodicities and the like should
865  // matter: boundaries should be EXACTLY the same (otherwise,
866  // map_grid() should be used)
867  if ( (std::fabs(other_grid.lower_boundaries[i] -
868  lower_boundaries[i]) > 1.0E-10) ||
869  (std::fabs(other_grid.upper_boundaries[i] -
870  upper_boundaries[i]) > 1.0E-10) ||
871  (std::fabs(other_grid.widths[i] -
872  widths[i]) > 1.0E-10) ||
873  (data.size() != other_grid.data.size()) ) {
874  cvm::error("Error: inconsistency between "
875  "two grids that are supposed to be equal, "
876  "aside from the data stored.\n");
877  return;
878  }
879  }
880  }
881 
882 
884  std::istream & read_restart(std::istream &is)
885  {
886  size_t const start_pos = is.tellg();
887  std::string key, conf;
888  if ((is >> key) && (key == std::string("grid_parameters"))) {
889  is.seekg(start_pos, std::ios::beg);
890  is >> colvarparse::read_block("grid_parameters", conf);
892  } else {
893  cvm::log("Grid parameters are missing in the restart file, using those from the configuration.\n");
894  is.seekg(start_pos, std::ios::beg);
895  }
896  read_raw(is);
897  return is;
898  }
899 
901  std::ostream & write_restart(std::ostream &os)
902  {
903  write_params(os);
904  write_raw(os);
905  return os;
906  }
907 
908 
912  std::ostream & write_raw(std::ostream &os,
913  size_t const buf_size = 3)
914  {
915  std::streamsize const w = os.width();
916  std::streamsize const p = os.precision();
917 
918  std::vector<int> ix = new_index();
919  size_t count = 0;
920  for ( ; index_ok(ix); incr(ix)) {
921  for (size_t imult = 0; imult < mult; imult++) {
922  os << " "
923  << std::setw(w) << std::setprecision(p)
924  << value_output(ix, imult);
925  if (((++count) % buf_size) == 0)
926  os << "\n";
927  }
928  }
929  // write a final newline only if buffer is not empty
930  if ((count % buf_size) != 0)
931  os << "\n";
932 
933  return os;
934  }
935 
937  std::istream & read_raw(std::istream &is)
938  {
939  size_t const start_pos = is.tellg();
940 
941  for (std::vector<int> ix = new_index(); index_ok(ix); incr(ix)) {
942  for (size_t imult = 0; imult < mult; imult++) {
943  T new_value;
944  if (is >> new_value) {
945  value_input(ix, new_value, imult);
946  } else {
947  is.clear();
948  is.seekg(start_pos, std::ios::beg);
949  is.setstate(std::ios::failbit);
950  cvm::error("Error: failed to read all of the grid points from file. Possible explanations: grid parameters in the configuration (lowerBoundary, upperBoundary, width) are different from those in the file, or the file is corrupt/incomplete.\n");
951  return is;
952  }
953  }
954  }
955 
956  has_data = true;
957  return is;
958  }
959 
962  void write_multicol(std::ostream &os)
963  {
964  std::streamsize const w = os.width();
965  std::streamsize const p = os.precision();
966 
967  // Data in the header: nColvars, then for each
968  // xiMin, dXi, nPoints, periodic
969 
970  os << std::setw(2) << "# " << nd << "\n";
971  for (size_t i = 0; i < nd; i++) {
972  os << "# "
973  << std::setw(10) << lower_boundaries[i]
974  << std::setw(10) << widths[i]
975  << std::setw(10) << nx[i] << " "
976  << periodic[i] << "\n";
977  }
978 
979 
980  for (std::vector<int> ix = new_index(); index_ok(ix); incr(ix) ) {
981 
982  if (ix.back() == 0) {
983  // if the last index is 0, add a new line to mark the new record
984  os << "\n";
985  }
986 
987  for (size_t i = 0; i < nd; i++) {
988  os << " "
989  << std::setw(w) << std::setprecision(p)
990  << bin_to_value_scalar(ix[i], i);
991  }
992  os << " ";
993  for (size_t imult = 0; imult < mult; imult++) {
994  os << " "
995  << std::setw(w) << std::setprecision(p)
996  << value_output(ix, imult);
997  }
998  os << "\n";
999  }
1000  }
1001 
1004  std::istream & read_multicol(std::istream &is, bool add = false)
1005  {
1006  // Data in the header: nColvars, then for each
1007  // xiMin, dXi, nPoints, periodic
1008 
1009  std::string hash;
1010  cvm::real lower, width, x;
1011  size_t n, periodic;
1012  bool remap;
1013  std::vector<T> new_value;
1014  std::vector<int> nx_read;
1015  std::vector<int> bin;
1016 
1017  if ( cv.size() != nd ) {
1018  cvm::error("Cannot read grid file: missing reference to colvars.");
1019  return is;
1020  }
1021 
1022  if ( !(is >> hash) || (hash != "#") ) {
1023  cvm::error("Error reading grid at position "+
1024  cvm::to_str(is.tellg())+" in stream(read \"" + hash + "\")\n");
1025  return is;
1026  }
1027 
1028  is >> n;
1029  if ( n != nd ) {
1030  cvm::error("Error reading grid: wrong number of collective variables.\n");
1031  return is;
1032  }
1033 
1034  nx_read.resize(n);
1035  bin.resize(n);
1036  new_value.resize(mult);
1037 
1038  if (this->has_parent_data && add) {
1039  new_data.resize(data.size());
1040  }
1041 
1042  remap = false;
1043  for (size_t i = 0; i < nd; i++ ) {
1044  if ( !(is >> hash) || (hash != "#") ) {
1045  cvm::error("Error reading grid at position "+
1046  cvm::to_str(is.tellg())+" in stream(read \"" + hash + "\")\n");
1047  return is;
1048  }
1049 
1050  is >> lower >> width >> nx_read[i] >> periodic;
1051 
1052 
1053  if ( (std::fabs(lower - lower_boundaries[i].real_value) > 1.0e-10) ||
1054  (std::fabs(width - widths[i] ) > 1.0e-10) ||
1055  (nx_read[i] != nx[i]) ) {
1056  cvm::log("Warning: reading from different grid definition (colvar "
1057  + cvm::to_str(i+1) + "); remapping data on new grid.\n");
1058  remap = true;
1059  }
1060  }
1061 
1062  if ( remap ) {
1063  // re-grid data
1064  while (is.good()) {
1065  bool end_of_file = false;
1066 
1067  for (size_t i = 0; i < nd; i++ ) {
1068  if ( !(is >> x) ) end_of_file = true;
1069  bin[i] = value_to_bin_scalar(x, i);
1070  }
1071  if (end_of_file) break;
1072 
1073  for (size_t imult = 0; imult < mult; imult++) {
1074  is >> new_value[imult];
1075  }
1076 
1077  if ( index_ok(bin) ) {
1078  for (size_t imult = 0; imult < mult; imult++) {
1079  value_input(bin, new_value[imult], imult, add);
1080  }
1081  }
1082  }
1083  } else {
1084  // do not re-grid the data but assume the same grid is used
1085  for (std::vector<int> ix = new_index(); index_ok(ix); incr(ix) ) {
1086  for (size_t i = 0; i < nd; i++ ) {
1087  is >> x;
1088  }
1089  for (size_t imult = 0; imult < mult; imult++) {
1090  is >> new_value[imult];
1091  value_input(ix, new_value[imult], imult, add);
1092  }
1093  }
1094  }
1095  has_data = true;
1096  return is;
1097  }
1098 
1101  std::ostream & write_opendx(std::ostream &os)
1102  {
1103  // write the header
1104  os << "object 1 class gridpositions counts";
1105  size_t icv;
1106  for (icv = 0; icv < number_of_colvars(); icv++) {
1107  os << " " << number_of_points(icv);
1108  }
1109  os << "\n";
1110 
1111  os << "origin";
1112  for (icv = 0; icv < number_of_colvars(); icv++) {
1113  os << " " << (lower_boundaries[icv].real_value + 0.5 * widths[icv]);
1114  }
1115  os << "\n";
1116 
1117  for (icv = 0; icv < number_of_colvars(); icv++) {
1118  os << "delta";
1119  for (size_t icv2 = 0; icv2 < number_of_colvars(); icv2++) {
1120  if (icv == icv2) os << " " << widths[icv];
1121  else os << " " << 0.0;
1122  }
1123  os << "\n";
1124  }
1125 
1126  os << "object 2 class gridconnections counts";
1127  for (icv = 0; icv < number_of_colvars(); icv++) {
1128  os << " " << number_of_points(icv);
1129  }
1130  os << "\n";
1131 
1132  os << "object 3 class array type double rank 0 items "
1133  << number_of_points() << " data follows\n";
1134 
1135  write_raw(os);
1136 
1137  os << "object \"collective variables scalar field\" class field\n";
1138  return os;
1139  }
1140 };
1141 
1142 
1143 
1146 class colvar_grid_count : public colvar_grid<size_t>
1147 {
1148 public:
1149 
1152 
1154  virtual inline ~colvar_grid_count()
1155  {}
1156 
1158  colvar_grid_count(std::vector<int> const &nx_i,
1159  size_t const &def_count = 0);
1160 
1162  colvar_grid_count(std::vector<colvar *> &colvars,
1163  size_t const &def_count = 0);
1164 
1166  inline void incr_count(std::vector<int> const &ix)
1167  {
1168  ++(data[this->address(ix)]);
1169  }
1170 
1172  inline size_t const & new_count(std::vector<int> const &ix,
1173  size_t const &imult = 0)
1174  {
1175  return new_data[address(ix) + imult];
1176  }
1177 
1181  virtual inline void value_input(std::vector<int> const &ix,
1182  size_t const &t,
1183  size_t const &imult = 0,
1184  bool add = false)
1185  {
1186  if (add) {
1187  data[address(ix)] += t;
1188  if (this->has_parent_data) {
1189  // save newly read data for inputting parent grid
1190  new_data[address(ix)] = t;
1191  }
1192  } else {
1193  data[address(ix)] = t;
1194  }
1195  has_data = true;
1196  }
1197 
1200  inline cvm::real log_gradient_finite_diff(const std::vector<int> &ix0,
1201  int n = 0)
1202  {
1203  int A0, A1, A2;
1204  std::vector<int> ix = ix0;
1205 
1206  if (periodic[n]) {
1207  ix[n]--; wrap(ix);
1208  A0 = data[address(ix)];
1209  ix = ix0;
1210  ix[n]++; wrap(ix);
1211  A1 = data[address(ix)];
1212  if (A0 * A1 == 0) {
1213  return 0.; // can't handle empty bins
1214  } else {
1215  return (std::log((cvm::real)A1) - std::log((cvm::real)A0))
1216  / (widths[n] * 2.);
1217  }
1218  } else if (ix[n] > 0 && ix[n] < nx[n]-1) { // not an edge
1219  ix[n]--;
1220  A0 = data[address(ix)];
1221  ix = ix0;
1222  ix[n]++;
1223  A1 = data[address(ix)];
1224  if (A0 * A1 == 0) {
1225  return 0.; // can't handle empty bins
1226  } else {
1227  return (std::log((cvm::real)A1) - std::log((cvm::real)A0))
1228  / (widths[n] * 2.);
1229  }
1230  } else {
1231  // edge: use 2nd order derivative
1232  int increment = (ix[n] == 0 ? 1 : -1);
1233  // move right from left edge, or the other way around
1234  A0 = data[address(ix)];
1235  ix[n] += increment; A1 = data[address(ix)];
1236  ix[n] += increment; A2 = data[address(ix)];
1237  if (A0 * A1 * A2 == 0) {
1238  return 0.; // can't handle empty bins
1239  } else {
1240  return (-1.5 * std::log((cvm::real)A0) + 2. * std::log((cvm::real)A1)
1241  - 0.5 * std::log((cvm::real)A2)) * increment / widths[n];
1242  }
1243  }
1244  }
1245 };
1246 
1247 
1249 class colvar_grid_scalar : public colvar_grid<cvm::real>
1250 {
1251 public:
1252 
1256 
1259 
1262 
1264  ~colvar_grid_scalar();
1265 
1267  colvar_grid_scalar(std::vector<int> const &nx_i);
1268 
1270  colvar_grid_scalar(std::vector<colvar *> &colvars,
1271  bool margin = 0);
1272 
1274  inline void acc_value(std::vector<int> const &ix,
1275  cvm::real const &new_value,
1276  size_t const &imult = 0)
1277  {
1278  // only legal value of imult here is 0
1279  data[address(ix)] += new_value;
1280  if (samples)
1281  samples->incr_count(ix);
1282  has_data = true;
1283  }
1284 
1286  inline const cvm::real * gradient_finite_diff( const std::vector<int> &ix0 )
1287  {
1288  cvm::real A0, A1;
1289  std::vector<int> ix;
1290  if (nd != 2) {
1291  cvm::error("Finite differences available in dimension 2 only.");
1292  return grad;
1293  }
1294  for (unsigned int n = 0; n < nd; n++) {
1295  ix = ix0;
1296  A0 = data[address(ix)];
1297  ix[n]++; wrap(ix);
1298  A1 = data[address(ix)];
1299  ix[1-n]++; wrap(ix);
1300  A1 += data[address(ix)];
1301  ix[n]--; wrap(ix);
1302  A0 += data[address(ix)];
1303  grad[n] = 0.5 * (A1 - A0) / widths[n];
1304  }
1305  return grad;
1306  }
1307 
1310  virtual cvm::real value_output(std::vector<int> const &ix,
1311  size_t const &imult = 0)
1312  {
1313  if (imult > 0) {
1314  cvm::error("Error: trying to access a component "
1315  "larger than 1 in a scalar data grid.\n");
1316  return 0.;
1317  }
1318  if (samples) {
1319  return (samples->value(ix) > 0) ?
1320  (data[address(ix)] / cvm::real(samples->value(ix))) :
1321  0.0;
1322  } else {
1323  return data[address(ix)];
1324  }
1325  }
1326 
1330  virtual void value_input(std::vector<int> const &ix,
1331  cvm::real const &new_value,
1332  size_t const &imult = 0,
1333  bool add = false)
1334  {
1335  if (imult > 0) {
1336  cvm::error("Error: trying to access a component "
1337  "larger than 1 in a scalar data grid.\n");
1338  return;
1339  }
1340  if (add) {
1341  if (samples)
1342  data[address(ix)] += new_value * samples->new_count(ix);
1343  else
1344  data[address(ix)] += new_value;
1345  } else {
1346  if (samples)
1347  data[address(ix)] = new_value * samples->value(ix);
1348  else
1349  data[address(ix)] = new_value;
1350  }
1351  has_data = true;
1352  }
1353 
1355  cvm::real maximum_value() const;
1356 
1358  cvm::real minimum_value() const;
1359 
1361  cvm::real minimum_pos_value() const;
1362 
1364  cvm::real integral() const;
1365 
1368  cvm::real entropy() const;
1369 
1370 private:
1371  // gradient
1372  cvm::real * grad;
1373 };
1374 
1375 
1376 
1378 class colvar_grid_gradient : public colvar_grid<cvm::real>
1379 {
1380 public:
1381 
1385 
1388 
1390  virtual inline ~colvar_grid_gradient()
1391  {}
1392 
1394  colvar_grid_gradient(std::vector<int> const &nx_i);
1395 
1397  colvar_grid_gradient(std::vector<colvar *> &colvars);
1398 
1400  inline void acc_value(std::vector<int> const &ix, std::vector<colvarvalue> const &values) {
1401  for (size_t imult = 0; imult < mult; imult++) {
1402  data[address(ix) + imult] += values[imult].real_value;
1403  }
1404  if (samples)
1405  samples->incr_count(ix);
1406  }
1407 
1409  inline void acc_grad(std::vector<int> const &ix, cvm::real const *grads) {
1410  for (size_t imult = 0; imult < mult; imult++) {
1411  data[address(ix) + imult] += grads[imult];
1412  }
1413  if (samples)
1414  samples->incr_count(ix);
1415  }
1416 
1419  inline void acc_force(std::vector<int> const &ix, cvm::real const *forces) {
1420  for (size_t imult = 0; imult < mult; imult++) {
1421  data[address(ix) + imult] -= forces[imult];
1422  }
1423  if (samples)
1424  samples->incr_count(ix);
1425  }
1426 
1429  virtual inline cvm::real value_output(std::vector<int> const &ix,
1430  size_t const &imult = 0)
1431  {
1432  if (samples)
1433  return (samples->value(ix) > 0) ?
1434  (data[address(ix) + imult] / cvm::real(samples->value(ix))) :
1435  0.0;
1436  else
1437  return data[address(ix) + imult];
1438  }
1439 
1443  virtual inline void value_input(std::vector<int> const &ix,
1444  cvm::real const &new_value,
1445  size_t const &imult = 0,
1446  bool add = false)
1447  {
1448  if (add) {
1449  if (samples)
1450  data[address(ix) + imult] += new_value * samples->new_count(ix);
1451  else
1452  data[address(ix) + imult] += new_value;
1453  } else {
1454  if (samples)
1455  data[address(ix) + imult] = new_value * samples->value(ix);
1456  else
1457  data[address(ix) + imult] = new_value;
1458  }
1459  has_data = true;
1460  }
1461 
1462 
1465  {
1466  size_t n = 0;
1467 
1468  if (nd != 1 || nx[0] == 0) {
1469  return 0.0;
1470  }
1471 
1472  cvm::real sum = 0.0;
1473  std::vector<int> ix = new_index();
1474  if (samples) {
1475  for ( ; index_ok(ix); incr(ix)) {
1476  if ( (n = samples->value(ix)) )
1477  sum += value(ix) / n;
1478  }
1479  } else {
1480  for ( ; index_ok(ix); incr(ix)) {
1481  sum += value(ix);
1482  }
1483  }
1484  return (sum / cvm::real(nx[0]));
1485  }
1486 
1489  void write_1D_integral(std::ostream &os);
1490 
1491 };
1492 
1493 
1494 #endif
1495 
std::ostream & write_params(std::ostream &os)
Write the grid parameters (number of colvars, boundaries, width and number of points) ...
Definition: colvargrid.h:741
Scalar number, implemented as colvarmodule::real (default)
Definition: colvarvalue.h:47
Like parse_normal, but don&#39;t send any message to the log (useful e.g. in restart files when such mess...
Definition: colvarparse.h:90
std::istream & read_restart(std::istream &is)
Read grid entry in restart file.
Definition: colvargrid.h:884
std::ostream & write_opendx(std::ostream &os)
Write the grid data without labels, as they are represented in memory.
Definition: colvargrid.h:1101
std::vector< bool > periodic
Whether some colvars are periodic.
Definition: colvargrid.h:76
colvar_grid(std::vector< int > const &nx_i, T const &t=T(), size_t mult_i=1)
Constructor from explicit grid sizes.
Definition: colvargrid.h:226
virtual cvm::real value_output(std::vector< int > const &ix, size_t const &imult=0)
Return the value of the function at ix divided by its number of samples (if the count grid is defined...
Definition: colvargrid.h:1310
const cvm::real * gradient_finite_diff(const std::vector< int > &ix0)
Return the gradient of the scalar field from finite differences.
Definition: colvargrid.h:1286
virtual void value_input(std::vector< int > const &ix, T const &t, size_t const &imult=0, bool add=false)
Get the value from a formatted output and transform it into the internal representation (the two may ...
Definition: colvargrid.h:678
cvm::real average()
Compute and return average value for a 1D gradient grid.
Definition: colvargrid.h:1464
std::vector< int > const new_index() const
Get the index corresponding to the "first" bin, to be used as the initial value for an index in loopi...
Definition: colvargrid.h:698
size_t mult
Multiplicity of each datum (allow the binning of non-scalar types such as atomic gradients) ...
Definition: colvargrid.h:34
colvar_grid_count * samples
Provide the associated sample count by which each binned value should be divided. ...
Definition: colvargrid.h:1255
cvm::real real_value
Real data member.
Definition: colvarvalue.h:68
std::vector< bool > hard_upper_boundaries
Whether some colvars have hard upper boundaries.
Definition: colvargrid.h:82
virtual void value_input(std::vector< int > const &ix, cvm::real const &new_value, size_t const &imult=0, bool add=false)
Get the value from a formatted output and transform it into the internal representation (it may have ...
Definition: colvargrid.h:1330
virtual ~colvar_grid()
Destructor.
Definition: colvargrid.h:200
virtual ~colvar_grid_count()
Destructor.
Definition: colvargrid.h:1154
std::vector< int > const get_colvars_index(std::vector< colvarvalue > const &values) const
Get the bin indices corresponding to the provided values of the colvars.
Definition: colvargrid.h:540
colvarvalue bin_to_value_scalar(int const &i_bin, int const i) const
Use the two boundaries and the width to report the central value corresponding to a bin index...
Definition: colvargrid.h:419
void acc_value(std::vector< int > const &ix, std::vector< colvarvalue > const &values)
Accumulate the value.
Definition: colvargrid.h:1400
size_t multiplicity() const
Return the multiplicity of the type used.
Definition: colvargrid.h:123
void set_sizes(std::vector< int > const &new_sizes)
Set the sizes in each direction.
Definition: colvargrid.h:117
cvm::real log_gradient_finite_diff(const std::vector< int > &ix0, int n=0)
Return the log-gradient from finite differences on the same grid for dimension n. ...
Definition: colvargrid.h:1200
int current_bin_scalar(int const i) const
Report the bin corresponding to the current value of variable i.
Definition: colvargrid.h:372
std::vector< int > nx
Number of points along each dimension.
Definition: colvargrid.h:27
void check_consistency(colvar_grid< T > const &other_grid)
Check that the grid information inside (boundaries, widths, ...) is consistent with that of another g...
Definition: colvargrid.h:861
bool has_parent_data
True if this is a count grid related to another grid of data.
Definition: colvargrid.h:88
void copy_grid(colvar_grid< T > const &other_grid)
Copy data from another grid of the same type, AND identical definition (boundaries, widths) Added for shared ABF.
Definition: colvargrid.h:469
std::istream & read_multicol(std::istream &is, bool add=false)
Read a grid written by colvar_grid::write_multicol() Adding data if add is true, replacing if false...
Definition: colvargrid.h:1004
int value_to_bin_scalar(colvarvalue const &value, const int i) const
Use the lower boundary and the width to report which bin the provided value is in.
Definition: colvargrid.h:394
Grid of values of a function of several collective variables.
Definition: colvargrid.h:19
size_t address(std::vector< int > const &ix) const
Get the low-level index corresponding to an index.
Definition: colvargrid.h:52
colvar_grid(colvar_grid< T > const &g)
"Almost copy-constructor": only copies configuration parameters from another grid, but doesn&#39;t reallocate stuff; setup() must be called after that;
Definition: colvargrid.h:206
void write_multicol(std::ostream &os)
Write the grid in a format which is both human readable and suitable for visualization e...
Definition: colvargrid.h:962
size_t number_of_points(int const icv=-1) const
Definition: colvargrid.h:101
std::vector< bool > actual_value
Do we request actual value (for extended-system colvars)?
Definition: colvargrid.h:49
std::vector< int > const get_colvars_index() const
Get the bin indices corresponding to the current values of the colvars.
Definition: colvargrid.h:551
int setup(std::vector< int > const &nx_i, T const &t=T(), size_t const &mult_i=1)
Allocate data.
Definition: colvargrid.h:138
Collective variables main module.
double real
Defining an abstract real number allows to switch precision.
Definition: colvarmodule.h:76
static int error(std::string const &message, int code=COLVARS_ERROR)
Print a message to the main log and set global error code.
Definition: colvarmodule.cpp:1587
void raw_data_in(const T *in_data)
Input the data as they are represented in memory.
Definition: colvargrid.h:497
bool index_ok(std::vector< int > const &ix) const
Check that the index is within range in each of the dimensions.
Definition: colvargrid.h:705
void acc_grad(std::vector< int > const &ix, cvm::real const *grads)
Accumulate the gradient.
Definition: colvargrid.h:1409
std::ostream & write_raw(std::ostream &os, size_t const buf_size=3)
Write the grid data without labels, as they are represented in memory.
Definition: colvargrid.h:912
colvarvalue bin_to_value_scalar(int const &i_bin, colvarvalue const &new_offset, cvm::real const &new_width) const
Same as the standard version, but uses different parameters.
Definition: colvargrid.h:425
int parse_params(std::string const &conf, colvarparse::Parse_Mode const parse_mode=colvarparse::parse_normal)
Read a grid definition from a config string.
Definition: colvargrid.h:771
colvar_grid()
Default constructor.
Definition: colvargrid.h:192
bool has_data
Whether this grid has been filled with data or is still empty.
Definition: colvargrid.h:91
Value of a collective variable: this is a metatype which can be set at runtime. By default it is set ...
Definition: colvarvalue.h:34
void check_consistency()
Check that the grid information inside (boundaries, widths, ...) is consistent with the current setti...
Definition: colvargrid.h:842
virtual cvm::real value_output(std::vector< int > const &ix, size_t const &imult=0)
Return the value of the function at ix divided by its number of samples (if the count grid is defined...
Definition: colvargrid.h:1429
void set_value(std::vector< int > const &ix, T const &t, size_t const &imult=0)
Set the value at the point with index ix.
Definition: colvargrid.h:433
size_t const & new_count(std::vector< int > const &ix, size_t const &imult=0)
Get the binned count indexed by ix from the newly read data.
Definition: colvargrid.h:1172
int setup()
Allocate data (allow initialization also after construction)
Definition: colvargrid.h:179
int value_to_bin_scalar(colvarvalue const &value, colvarvalue const &new_offset, cvm::real const &new_width) const
Same as the standard version, but uses another grid definition.
Definition: colvargrid.h:410
virtual void value_input(std::vector< int > const &ix, cvm::real const &new_value, size_t const &imult=0, bool add=false)
Get the value from a formatted output and transform it into the internal representation (it may have ...
Definition: colvargrid.h:1443
size_t number_of_colvars() const
Return the number of colvars.
Definition: colvargrid.h:94
std::vector< colvarvalue > upper_boundaries
Upper boundaries of the colvars in this grid.
Definition: colvargrid.h:73
Helper class to read a block of the type "key { ... }" from a stream and store it in a string...
Definition: colvarparse.h:251
size_t nt
Total number of grid points.
Definition: colvargrid.h:37
virtual ~colvar_grid_gradient()
Destructor.
Definition: colvargrid.h:1390
std::ostream & write_restart(std::ostream &os)
Write grid entry in restart file.
Definition: colvargrid.h:901
void request_actual_value(bool b=true)
Request grid to use actual values of extended coords.
Definition: colvargrid.h:129
virtual T value_output(std::vector< int > const &ix, size_t const &imult=0)
Return the value suitable for output purposes (so that it may be rescaled or manipulated without chan...
Definition: colvargrid.h:669
void map_grid(colvar_grid< T > const &other_grid)
Add data from another grid of the same type.
Definition: colvargrid.h:604
size_t nd
Number of dimensions.
Definition: colvargrid.h:24
void reset(T const &t=T())
Reset data (in case the grid is being reused)
Definition: colvargrid.h:185
std::vector< int > nxc
Cumulative number of points along each dimension.
Definition: colvargrid.h:30
colvar_grid_count * samples
Provide the sample count by which each binned value should be divided.
Definition: colvargrid.h:1384
(default) Read the first instance of a keyword (if any), report its value, and print a warning when t...
Definition: colvarparse.h:86
size_t raw_data_num() const
Size of the data as they are represented in memory.
Definition: colvargrid.h:503
void acc_force(std::vector< int > const &ix, cvm::real const *forces)
Accumulate the gradient based on the force (i.e. sums the opposite of the force)
Definition: colvargrid.h:1419
void acc_value(std::vector< int > const &ix, cvm::real const &new_value, size_t const &imult=0)
Accumulate the value.
Definition: colvargrid.h:1274
void add_grid(colvar_grid< T > const &other_grid, cvm::real scale_factor=1.0)
Add data from another grid of the same type, AND identical definition (boundaries, widths)
Definition: colvargrid.h:647
void incr(std::vector< int > &ix) const
Increment the index, in a way that will make it loop over the whole nd-dimensional array...
Definition: colvargrid.h:716
std::istream & read_raw(std::istream &is)
Read data written by colvar_grid::write_raw()
Definition: colvargrid.h:937
std::vector< bool > hard_lower_boundaries
Whether some colvars have hard lower boundaries.
Definition: colvargrid.h:79
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
int value_to_bin_scalar_bound(colvarvalue const &value, const int i) const
Use the lower boundary and the width to report which bin the provided value is in and assign first or...
Definition: colvargrid.h:401
Class for accumulating a scalar function on a grid.
Definition: colvargrid.h:1249
int current_bin_scalar(int const i, int const iv) const
Report the bin corresponding to the current value of item iv in variable i.
Definition: colvargrid.h:385
void add_constant(T const &t)
Add a constant to all elements (fast loop)
Definition: colvargrid.h:516
Colvar_grid derived class to hold counters in discrete n-dim colvar space.
Definition: colvargrid.h:1146
static size_t const cv_prec
Number of digits to represent a collective variables value(s)
Definition: colvarmodule.h:427
Parse_Mode
How a keyword is parsed in a string.
Definition: colvarparse.h:82
virtual void value_input(std::vector< int > const &ix, size_t const &t, size_t const &imult=0, bool add=false)
Get the value from a formatted output and transform it into the internal representation (it may have ...
Definition: colvargrid.h:1181
std::vector< int > const get_colvars_index_bound() const
Get the bin indices corresponding to the provided values of the colvars and assign first or last bin ...
Definition: colvargrid.h:562
Parsing functions for collective variables.
Base class containing parsing functions; all objects which need to parse input inherit from this...
Definition: colvarparse.h:18
std::vector< size_t > new_data
Newly read data (used for count grids, when adding several grids read from disk)
Definition: colvargrid.h:43
static void log(std::string const &message)
Print a message to the main log.
Definition: colvarmodule.cpp:1518
colvar_grid(std::vector< colvar *> const &colvars, T const &t=T(), size_t mult_i=1, bool margin=false)
Constructor from a vector of colvars.
Definition: colvargrid.h:235
std::vector< colvarvalue > lower_boundaries
Lower boundaries of the colvars in this grid.
Definition: colvargrid.h:70
void incr_count(std::vector< int > const &ix)
Increment the counter at given position.
Definition: colvargrid.h:1166
void multiply_constant(cvm::real const &a)
Multiply all elements by a scalar constant (fast loop)
Definition: colvargrid.h:524
std::vector< int > const & sizes() const
Get the sizes in each direction.
Definition: colvargrid.h:111
std::vector< colvar * > cv
Colvars collected in this grid.
Definition: colvargrid.h:46
void delta_grid(colvar_grid< T > const &other_grid)
Get the change from this to other_grid and store the result in this. this_grid := other_grid - this_g...
Definition: colvargrid.h:445
std::vector< cvm::real > widths
Widths of the colvars in this grid.
Definition: colvargrid.h:85
int current_bin_scalar_bound(int const i) const
Report the bin corresponding to the current value of variable i and assign first or last bin if out o...
Definition: colvargrid.h:379
Class for accumulating the gradient of a scalar function on a grid.
Definition: colvargrid.h:1378
T const & value(std::vector< int > const &ix, size_t const &imult=0) const
Get the binned value indexed by ix, or the first of them if the multiplicity is larger than 1...
Definition: colvargrid.h:508
void remove_zeros(cvm::real const &a)
Assign all zero elements a scalar constant (fast loop)
Definition: colvargrid.h:531
void raw_data_out(T *out_data) const
Extract the grid data as they are represented in memory. Put the results in "out_data".
Definition: colvargrid.h:492
std::vector< T > data
Low-level array of values.
Definition: colvargrid.h:40
static size_t const cv_width
Number of characters to represent a collective variables value(s)
Definition: colvarmodule.h:429
void wrap(std::vector< int > &ix) const
Definition: colvargrid.h:355
static bool debug()
Whether debug output should be enabled (compile-time option)
Definition: colvarmodule.h:237
cvm::real bin_distance_from_boundaries(std::vector< colvarvalue > const &values, bool skip_hard_boundaries=false)
Get the minimal distance (in number of bins) from the boundaries; a negative number is returned if th...
Definition: colvargrid.h:574