Collective Variables Module - Developer Documentation
 All Classes Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
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 
193  {
194  save_delimiters = false;
195  nd = nt = 0;
196  mult = 1;
197  this->setup();
198  }
199 
201  virtual ~colvar_grid()
202  {}
203 
207  colvar_grid(colvar_grid<T> const &g) : nd(g.nd),
208  nx(g.nx),
209  mult(g.mult),
210  data(),
211  cv(g.cv),
215  periodic(g.periodic),
218  widths(g.widths),
219  has_data(false)
220  {
221  save_delimiters = false;
222  }
223 
228  colvar_grid(std::vector<int> const &nx_i,
229  T const &t = T(),
230  size_t mult_i = 1)
231  : has_data(false)
232  {
233  save_delimiters = false;
234  this->setup(nx_i, t, mult_i);
235  }
236 
238  colvar_grid(std::vector<colvar *> const &colvars,
239  T const &t = T(),
240  size_t mult_i = 1,
241  bool margin = false)
242  : has_data(false)
243  {
244  save_delimiters = false;
245  this->init_from_colvars(colvars, t, mult_i, margin);
246  }
247 
248  int init_from_colvars(std::vector<colvar *> const &colvars,
249  T const &t = T(),
250  size_t mult_i = 1,
251  bool margin = false)
252  {
253  if (cvm::debug()) {
254  cvm::log("Reading grid configuration from collective variables.\n");
255  }
256 
257  cv = colvars;
258  nd = colvars.size();
259  mult = mult_i;
260 
261  size_t i;
262 
263  if (cvm::debug()) {
264  cvm::log("Allocating a grid for "+cvm::to_str(colvars.size())+
265  " collective variables, multiplicity = "+cvm::to_str(mult_i)+".\n");
266  }
267 
268  for (i = 0; i < cv.size(); i++) {
269 
270  if (cv[i]->value().type() != colvarvalue::type_scalar) {
271  cvm::error("Colvar grids can only be automatically "
272  "constructed for scalar variables. "
273  "ABF and histogram can not be used; metadynamics "
274  "can be used with useGrids disabled.\n", INPUT_ERROR);
275  return COLVARS_ERROR;
276  }
277 
278  if (cv[i]->width <= 0.0) {
279  cvm::error("Tried to initialize a grid on a "
280  "variable with negative or zero width.\n", INPUT_ERROR);
281  return COLVARS_ERROR;
282  }
283 
284  widths.push_back(cv[i]->width);
285  hard_lower_boundaries.push_back(cv[i]->hard_lower_boundary);
286  hard_upper_boundaries.push_back(cv[i]->hard_upper_boundary);
287  periodic.push_back(cv[i]->periodic_boundaries());
288 
289  // By default, get reported colvar value (for extended Lagrangian colvars)
290  actual_value.push_back(false);
291 
292  // except if a colvar is specified twice in a row
293  // then the first instance is the actual value
294  // For histograms of extended-system coordinates
295  if (i > 0 && cv[i-1] == cv[i]) {
296  actual_value[i-1] = true;
297  }
298 
299  if (margin) {
300  if (periodic[i]) {
301  // Shift the grid by half the bin width (values at edges instead of center of bins)
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  } else {
305  // Make this grid larger by one bin width
306  lower_boundaries.push_back(cv[i]->lower_boundary.real_value - 0.5 * widths[i]);
307  upper_boundaries.push_back(cv[i]->upper_boundary.real_value + 0.5 * widths[i]);
308  }
309  } else {
310  lower_boundaries.push_back(cv[i]->lower_boundary);
311  upper_boundaries.push_back(cv[i]->upper_boundary);
312  }
313  }
314 
315  this->init_from_boundaries();
316  return this->setup();
317  }
318 
319  int init_from_boundaries(T const &t = T(),
320  size_t const &mult_i = 1)
321  {
322  if (cvm::debug()) {
323  cvm::log("Configuring grid dimensions from colvars boundaries.\n");
324  }
325 
326  // these will have to be updated
327  nx.clear();
328  nxc.clear();
329  nt = 0;
330 
331  for (size_t i = 0; i < lower_boundaries.size(); i++) {
332 
333  cvm::real nbins = ( upper_boundaries[i].real_value -
334  lower_boundaries[i].real_value ) / widths[i];
335  int nbins_round = (int)(nbins+0.5);
336 
337  if (std::fabs(nbins - cvm::real(nbins_round)) > 1.0E-10) {
338  cvm::log("Warning: grid interval("+
341  ") is not commensurate to its bin width("+
343  upper_boundaries[i].real_value = lower_boundaries[i].real_value +
344  (nbins_round * widths[i]);
345  }
346 
347  if (cvm::debug())
348  cvm::log("Number of points is "+cvm::to_str((int) nbins_round)+
349  " for the colvar no. "+cvm::to_str(i+1)+".\n");
350 
351  nx.push_back(nbins_round);
352  }
353 
354  return COLVARS_OK;
355  }
356 
359  inline void wrap(std::vector<int> & ix) const
360  {
361  for (size_t i = 0; i < nd; i++) {
362  if (periodic[i]) {
363  ix[i] = (ix[i] + nx[i]) % nx[i]; //to ensure non-negative result
364  } else {
365  if (ix[i] < 0 || ix[i] >= nx[i]) {
366  cvm::error("Trying to wrap illegal index vector (non-PBC) for a grid point: "
367  + cvm::to_str(ix), BUG_ERROR);
368  return;
369  }
370  }
371  }
372  }
373 
374 
376  inline int current_bin_scalar(int const i) const
377  {
378  return value_to_bin_scalar(actual_value[i] ? cv[i]->actual_value() : cv[i]->value(), i);
379  }
380 
383  inline int current_bin_scalar_bound(int const i) const
384  {
385  return value_to_bin_scalar_bound(actual_value[i] ? cv[i]->actual_value() : cv[i]->value(), i);
386  }
387 
389  inline int current_bin_scalar(int const i, int const iv) const
390  {
392  cv[i]->actual_value().vector1d_value[iv] :
393  cv[i]->value().vector1d_value[iv], i);
394  }
395 
398  inline int value_to_bin_scalar(colvarvalue const &value, const int i) const
399  {
400  return (int) std::floor( (value.real_value - lower_boundaries[i].real_value) / widths[i] );
401  }
402 
405  inline int value_to_bin_scalar_bound(colvarvalue const &value, const int i) const
406  {
407  int bin_index = std::floor( (value.real_value - lower_boundaries[i].real_value) / widths[i] );
408  if (bin_index < 0) bin_index=0;
409  if (bin_index >=int(nx[i])) bin_index=int(nx[i])-1;
410  return (int) bin_index;
411  }
412 
415  colvarvalue const &new_offset,
416  cvm::real const &new_width) const
417  {
418  return (int) std::floor( (value.real_value - new_offset.real_value) / new_width );
419  }
420 
423  inline colvarvalue bin_to_value_scalar(int const &i_bin, int const i) const
424  {
425  return lower_boundaries[i].real_value + widths[i] * (0.5 + i_bin);
426  }
427 
429  inline colvarvalue bin_to_value_scalar(int const &i_bin,
430  colvarvalue const &new_offset,
431  cvm::real const &new_width) const
432  {
433  return new_offset.real_value + new_width * (0.5 + i_bin);
434  }
435 
437  inline void set_value(std::vector<int> const &ix,
438  T const &t,
439  size_t const &imult = 0)
440  {
441  data[this->address(ix)+imult] = t;
442  has_data = true;
443  }
444 
449  void delta_grid(colvar_grid<T> const &other_grid)
450  {
451 
452  if (other_grid.multiplicity() != this->multiplicity()) {
453  cvm::error("Error: trying to subtract two grids with "
454  "different multiplicity.\n");
455  return;
456  }
457 
458  if (other_grid.data.size() != this->data.size()) {
459  cvm::error("Error: trying to subtract two grids with "
460  "different size.\n");
461  return;
462  }
463 
464  for (size_t i = 0; i < data.size(); i++) {
465  data[i] = other_grid.data[i] - data[i];
466  }
467  has_data = true;
468  }
469 
473  void copy_grid(colvar_grid<T> const &other_grid)
474  {
475  if (other_grid.multiplicity() != this->multiplicity()) {
476  cvm::error("Error: trying to copy two grids with "
477  "different multiplicity.\n");
478  return;
479  }
480 
481  if (other_grid.data.size() != this->data.size()) {
482  cvm::error("Error: trying to copy two grids with "
483  "different size.\n");
484  return;
485  }
486 
487 
488  for (size_t i = 0; i < data.size(); i++) {
489  data[i] = other_grid.data[i];
490  }
491  has_data = true;
492  }
493 
496  void raw_data_out(T* out_data) const
497  {
498  for (size_t i = 0; i < data.size(); i++) out_data[i] = data[i];
499  }
501  void raw_data_in(const T* in_data)
502  {
503  for (size_t i = 0; i < data.size(); i++) data[i] = in_data[i];
504  has_data = true;
505  }
507  size_t raw_data_num() const { return data.size(); }
508 
509 
512  inline T const & value(std::vector<int> const &ix,
513  size_t const &imult = 0) const
514  {
515  return data[this->address(ix) + imult];
516  }
517 
518 
520  inline void add_constant(T const &t)
521  {
522  for (size_t i = 0; i < nt; i++)
523  data[i] += t;
524  has_data = true;
525  }
526 
528  inline void multiply_constant(cvm::real const &a)
529  {
530  for (size_t i = 0; i < nt; i++)
531  data[i] *= a;
532  }
533 
535  inline void remove_zeros(cvm::real const &a)
536  {
537  for (size_t i = 0; i < nt; i++)
538  if(data[i]==0) data[i] = a;
539  }
540 
541 
544  inline std::vector<int> const get_colvars_index(std::vector<colvarvalue> const &values) const
545  {
546  std::vector<int> index = new_index();
547  for (size_t i = 0; i < nd; i++) {
548  index[i] = value_to_bin_scalar(values[i], i);
549  }
550  return index;
551  }
552 
555  inline std::vector<int> const get_colvars_index() const
556  {
557  std::vector<int> index = new_index();
558  for (size_t i = 0; i < nd; i++) {
559  index[i] = current_bin_scalar(i);
560  }
561  return index;
562  }
563 
566  inline std::vector<int> const get_colvars_index_bound() const
567  {
568  std::vector<int> index = new_index();
569  for (size_t i = 0; i < nd; i++) {
570  index[i] = current_bin_scalar_bound(i);
571  }
572  return index;
573  }
574 
578  inline cvm::real bin_distance_from_boundaries(std::vector<colvarvalue> const &values,
579  bool skip_hard_boundaries = false)
580  {
581  cvm::real minimum = 1.0E+16;
582  for (size_t i = 0; i < nd; i++) {
583 
584  if (periodic[i]) continue;
585 
586  cvm::real dl = std::sqrt(cv[i]->dist2(values[i], lower_boundaries[i])) / widths[i];
587  cvm::real du = std::sqrt(cv[i]->dist2(values[i], upper_boundaries[i])) / widths[i];
588 
589  if (values[i].real_value < lower_boundaries[i])
590  dl *= -1.0;
591  if (values[i].real_value > upper_boundaries[i])
592  du *= -1.0;
593 
594  if ( ((!skip_hard_boundaries) || (!hard_lower_boundaries[i])) && (dl < minimum))
595  minimum = dl;
596  if ( ((!skip_hard_boundaries) || (!hard_upper_boundaries[i])) && (du < minimum))
597  minimum = du;
598  }
599 
600  return minimum;
601  }
602 
603 
608  void map_grid(colvar_grid<T> const &other_grid)
609  {
610  if (other_grid.multiplicity() != this->multiplicity()) {
611  cvm::error("Error: trying to merge two grids with values of "
612  "different multiplicity.\n");
613  return;
614  }
615 
616  std::vector<colvarvalue> const &gb = this->lower_boundaries;
617  std::vector<cvm::real> const &gw = this->widths;
618  std::vector<colvarvalue> const &ogb = other_grid.lower_boundaries;
619  std::vector<cvm::real> const &ogw = other_grid.widths;
620 
621  std::vector<int> ix = this->new_index();
622  std::vector<int> oix = other_grid.new_index();
623 
624  if (cvm::debug())
625  cvm::log("Remapping grid...\n");
626  for ( ; this->index_ok(ix); this->incr(ix)) {
627 
628  for (size_t i = 0; i < nd; i++) {
629  oix[i] =
630  value_to_bin_scalar(bin_to_value_scalar(ix[i], gb[i], gw[i]),
631  ogb[i],
632  ogw[i]);
633  }
634 
635  if (! other_grid.index_ok(oix)) {
636  continue;
637  }
638 
639  for (size_t im = 0; im < mult; im++) {
640  this->set_value(ix, other_grid.value(oix, im), im);
641  }
642  }
643 
644  has_data = true;
645  if (cvm::debug())
646  cvm::log("Remapping done.\n");
647  }
648 
651  void add_grid(colvar_grid<T> const &other_grid,
652  cvm::real scale_factor = 1.0)
653  {
654  if (other_grid.multiplicity() != this->multiplicity()) {
655  cvm::error("Error: trying to sum togetehr two grids with values of "
656  "different multiplicity.\n");
657  return;
658  }
659  if (scale_factor != 1.0)
660  for (size_t i = 0; i < data.size(); i++) {
661  data[i] += scale_factor * other_grid.data[i];
662  }
663  else
664  // skip multiplication if possible
665  for (size_t i = 0; i < data.size(); i++) {
666  data[i] += other_grid.data[i];
667  }
668  has_data = true;
669  }
670 
673  virtual inline T value_output(std::vector<int> const &ix,
674  size_t const &imult = 0)
675  {
676  return value(ix, imult);
677  }
678 
682  virtual inline void value_input(std::vector<int> const &ix,
683  T const &t,
684  size_t const &imult = 0,
685  bool add = false)
686  {
687  if ( add )
688  data[address(ix) + imult] += t;
689  else
690  data[address(ix) + imult] = t;
691  has_data = true;
692  }
693 
694  // /// Get the pointer to the binned value indexed by ix
695  // inline T const *value_p (std::vector<int> const &ix)
696  // {
697  // return &(data[address (ix)]);
698  // }
699 
702  inline std::vector<int> const new_index() const
703  {
704  return std::vector<int> (nd, 0);
705  }
706 
709  inline bool index_ok(std::vector<int> const &ix) const
710  {
711  for (size_t i = 0; i < nd; i++) {
712  if ( (ix[i] < 0) || (ix[i] >= int(nx[i])) )
713  return false;
714  }
715  return true;
716  }
717 
720  inline void incr(std::vector<int> &ix) const
721  {
722  for (int i = ix.size()-1; i >= 0; i--) {
723 
724  ix[i]++;
725 
726  if (ix[i] >= nx[i]) {
727 
728  if (i > 0) {
729  ix[i] = 0;
730  continue;
731  } else {
732  // this is the last iteration, a non-valid index is being
733  // set for the outer index, which will be caught by
734  // index_ok()
735  ix[0] = nx[0];
736  return;
737  }
738  } else {
739  return;
740  }
741  }
742  }
743 
745  std::ostream & write_params(std::ostream &os)
746  {
747  size_t i;
748  os << "grid_parameters {\n n_colvars " << nd << "\n";
749 
750  os << " lower_boundaries ";
751  for (i = 0; i < nd; i++)
752  os << " " << lower_boundaries[i];
753  os << "\n";
754 
755  os << " upper_boundaries ";
756  for (i = 0; i < nd; i++)
757  os << " " << upper_boundaries[i];
758  os << "\n";
759 
760  os << " widths ";
761  for (i = 0; i < nd; i++)
762  os << " " << widths[i];
763  os << "\n";
764 
765  os << " sizes ";
766  for (i = 0; i < nd; i++)
767  os << " " << nx[i];
768  os << "\n";
769 
770  os << "}\n";
771  return os;
772  }
773 
775  int parse_params(std::string const &conf,
777  {
778  if (cvm::debug()) cvm::log("Reading grid configuration from string.\n");
779 
780  std::vector<int> old_nx = nx;
781  std::vector<colvarvalue> old_lb = lower_boundaries;
782 
783  {
784  size_t nd_in = 0;
785  // this is only used in state files
786  colvarparse::get_keyval(conf, "n_colvars", nd_in, nd, colvarparse::parse_silent);
787  if (nd_in != nd) {
788  cvm::error("Error: trying to read data for a grid "
789  "that contains a different number of colvars ("+
790  cvm::to_str(nd_in)+") than the grid defined "
791  "in the configuration file("+cvm::to_str(nd)+
792  ").\n");
793  return COLVARS_ERROR;
794  }
795  }
796 
797  // underscore keywords are used in state file
798  colvarparse::get_keyval(conf, "lower_boundaries",
800  colvarparse::get_keyval(conf, "upper_boundaries",
802 
803  // camel case keywords are used in config file
804  colvarparse::get_keyval(conf, "lowerBoundaries",
805  lower_boundaries, lower_boundaries, parse_mode);
806  colvarparse::get_keyval(conf, "upperBoundaries",
807  upper_boundaries, upper_boundaries, parse_mode);
808 
809  colvarparse::get_keyval(conf, "widths", widths, widths, parse_mode);
810 
811  // only used in state file
812  colvarparse::get_keyval(conf, "sizes", nx, nx, colvarparse::parse_silent);
813 
814  if (nd < lower_boundaries.size()) nd = lower_boundaries.size();
815 
816  if (! actual_value.size()) actual_value.assign(nd, false);
817  if (! periodic.size()) periodic.assign(nd, false);
818  if (! widths.size()) widths.assign(nd, 1.0);
819 
820  bool new_params = false;
821  if (old_nx.size()) {
822  for (size_t i = 0; i < nd; i++) {
823  if ( (old_nx[i] != nx[i]) ||
824  (std::sqrt(cv[i]->dist2(old_lb[i],
825  lower_boundaries[i])) > 1.0E-10) ) {
826  new_params = true;
827  }
828  }
829  } else {
830  new_params = true;
831  }
832 
833  // reallocate the array in case the grid params have just changed
834  if (new_params) {
835  init_from_boundaries();
836  // data.resize(0); // no longer needed: setup calls clear()
837  return this->setup(nx, T(), mult);
838  }
839 
840  return COLVARS_OK;
841  }
842 
847  {
848  for (size_t i = 0; i < nd; i++) {
849  if ( (std::sqrt(cv[i]->dist2(cv[i]->lower_boundary,
850  lower_boundaries[i])) > 1.0E-10) ||
851  (std::sqrt(cv[i]->dist2(cv[i]->upper_boundary,
852  upper_boundaries[i])) > 1.0E-10) ||
853  (std::sqrt(cv[i]->dist2(cv[i]->width,
854  widths[i])) > 1.0E-10) ) {
855  cvm::error("Error: restart information for a grid is "
856  "inconsistent with that of its colvars.\n");
857  return;
858  }
859  }
860  }
861 
862 
865  void check_consistency(colvar_grid<T> const &other_grid)
866  {
867  for (size_t i = 0; i < nd; i++) {
868  // we skip dist2(), because periodicities and the like should
869  // matter: boundaries should be EXACTLY the same (otherwise,
870  // map_grid() should be used)
871  if ( (std::fabs(other_grid.lower_boundaries[i] -
872  lower_boundaries[i]) > 1.0E-10) ||
873  (std::fabs(other_grid.upper_boundaries[i] -
874  upper_boundaries[i]) > 1.0E-10) ||
875  (std::fabs(other_grid.widths[i] -
876  widths[i]) > 1.0E-10) ||
877  (data.size() != other_grid.data.size()) ) {
878  cvm::error("Error: inconsistency between "
879  "two grids that are supposed to be equal, "
880  "aside from the data stored.\n");
881  return;
882  }
883  }
884  }
885 
886 
888  std::istream & read_restart(std::istream &is)
889  {
890  size_t const start_pos = is.tellg();
891  std::string key, conf;
892  if ((is >> key) && (key == std::string("grid_parameters"))) {
893  is.seekg(start_pos, std::ios::beg);
894  is >> colvarparse::read_block("grid_parameters", conf);
896  } else {
897  cvm::log("Grid parameters are missing in the restart file, using those from the configuration.\n");
898  is.seekg(start_pos, std::ios::beg);
899  }
900  read_raw(is);
901  return is;
902  }
903 
905  std::ostream & write_restart(std::ostream &os)
906  {
907  write_params(os);
908  write_raw(os);
909  return os;
910  }
911 
912 
916  std::ostream & write_raw(std::ostream &os,
917  size_t const buf_size = 3)
918  {
919  std::streamsize const w = os.width();
920  std::streamsize const p = os.precision();
921 
922  std::vector<int> ix = new_index();
923  size_t count = 0;
924  for ( ; index_ok(ix); incr(ix)) {
925  for (size_t imult = 0; imult < mult; imult++) {
926  os << " "
927  << std::setw(w) << std::setprecision(p)
928  << value_output(ix, imult);
929  if (((++count) % buf_size) == 0)
930  os << "\n";
931  }
932  }
933  // write a final newline only if buffer is not empty
934  if ((count % buf_size) != 0)
935  os << "\n";
936 
937  return os;
938  }
939 
941  std::istream & read_raw(std::istream &is)
942  {
943  size_t const start_pos = is.tellg();
944 
945  for (std::vector<int> ix = new_index(); index_ok(ix); incr(ix)) {
946  for (size_t imult = 0; imult < mult; imult++) {
947  T new_value;
948  if (is >> new_value) {
949  value_input(ix, new_value, imult);
950  } else {
951  is.clear();
952  is.seekg(start_pos, std::ios::beg);
953  is.setstate(std::ios::failbit);
954  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");
955  return is;
956  }
957  }
958  }
959 
960  has_data = true;
961  return is;
962  }
963 
966  void write_multicol(std::ostream &os)
967  {
968  std::streamsize const w = os.width();
969  std::streamsize const p = os.precision();
970 
971  // Data in the header: nColvars, then for each
972  // xiMin, dXi, nPoints, periodic
973 
974  os << std::setw(2) << "# " << nd << "\n";
975  for (size_t i = 0; i < nd; i++) {
976  os << "# "
977  << std::setw(10) << lower_boundaries[i]
978  << std::setw(10) << widths[i]
979  << std::setw(10) << nx[i] << " "
980  << periodic[i] << "\n";
981  }
982 
983 
984  for (std::vector<int> ix = new_index(); index_ok(ix); incr(ix) ) {
985 
986  if (ix.back() == 0) {
987  // if the last index is 0, add a new line to mark the new record
988  os << "\n";
989  }
990 
991  for (size_t i = 0; i < nd; i++) {
992  os << " "
993  << std::setw(w) << std::setprecision(p)
994  << bin_to_value_scalar(ix[i], i);
995  }
996  os << " ";
997  for (size_t imult = 0; imult < mult; imult++) {
998  os << " "
999  << std::setw(w) << std::setprecision(p)
1000  << value_output(ix, imult);
1001  }
1002  os << "\n";
1003  }
1004  }
1005 
1008  std::istream & read_multicol(std::istream &is, bool add = false)
1009  {
1010  // Data in the header: nColvars, then for each
1011  // xiMin, dXi, nPoints, periodic
1012 
1013  std::string hash;
1014  cvm::real lower, width, x;
1015  size_t n, periodic;
1016  bool remap;
1017  std::vector<T> new_value;
1018  std::vector<int> nx_read;
1019  std::vector<int> bin;
1020 
1021  if ( cv.size() != nd ) {
1022  cvm::error("Cannot read grid file: missing reference to colvars.");
1023  return is;
1024  }
1025 
1026  if ( !(is >> hash) || (hash != "#") ) {
1027  cvm::error("Error reading grid at position "+
1028  cvm::to_str(is.tellg())+" in stream(read \"" + hash + "\")\n");
1029  return is;
1030  }
1031 
1032  is >> n;
1033  if ( n != nd ) {
1034  cvm::error("Error reading grid: wrong number of collective variables.\n");
1035  return is;
1036  }
1037 
1038  nx_read.resize(n);
1039  bin.resize(n);
1040  new_value.resize(mult);
1041 
1042  if (this->has_parent_data && add) {
1043  new_data.resize(data.size());
1044  }
1045 
1046  remap = false;
1047  for (size_t i = 0; i < nd; i++ ) {
1048  if ( !(is >> hash) || (hash != "#") ) {
1049  cvm::error("Error reading grid at position "+
1050  cvm::to_str(is.tellg())+" in stream(read \"" + hash + "\")\n");
1051  return is;
1052  }
1053 
1054  is >> lower >> width >> nx_read[i] >> periodic;
1055 
1056 
1057  if ( (std::fabs(lower - lower_boundaries[i].real_value) > 1.0e-10) ||
1058  (std::fabs(width - widths[i] ) > 1.0e-10) ||
1059  (nx_read[i] != nx[i]) ) {
1060  cvm::log("Warning: reading from different grid definition (colvar "
1061  + cvm::to_str(i+1) + "); remapping data on new grid.\n");
1062  remap = true;
1063  }
1064  }
1065 
1066  if ( remap ) {
1067  // re-grid data
1068  while (is.good()) {
1069  bool end_of_file = false;
1070 
1071  for (size_t i = 0; i < nd; i++ ) {
1072  if ( !(is >> x) ) end_of_file = true;
1073  bin[i] = value_to_bin_scalar(x, i);
1074  }
1075  if (end_of_file) break;
1076 
1077  for (size_t imult = 0; imult < mult; imult++) {
1078  is >> new_value[imult];
1079  }
1080 
1081  if ( index_ok(bin) ) {
1082  for (size_t imult = 0; imult < mult; imult++) {
1083  value_input(bin, new_value[imult], imult, add);
1084  }
1085  }
1086  }
1087  } else {
1088  // do not re-grid the data but assume the same grid is used
1089  for (std::vector<int> ix = new_index(); index_ok(ix); incr(ix) ) {
1090  for (size_t i = 0; i < nd; i++ ) {
1091  is >> x;
1092  }
1093  for (size_t imult = 0; imult < mult; imult++) {
1094  is >> new_value[imult];
1095  value_input(ix, new_value[imult], imult, add);
1096  }
1097  }
1098  }
1099  has_data = true;
1100  return is;
1101  }
1102 
1105  std::ostream & write_opendx(std::ostream &os)
1106  {
1107  // write the header
1108  os << "object 1 class gridpositions counts";
1109  size_t icv;
1110  for (icv = 0; icv < number_of_colvars(); icv++) {
1111  os << " " << number_of_points(icv);
1112  }
1113  os << "\n";
1114 
1115  os << "origin";
1116  for (icv = 0; icv < number_of_colvars(); icv++) {
1117  os << " " << (lower_boundaries[icv].real_value + 0.5 * widths[icv]);
1118  }
1119  os << "\n";
1120 
1121  for (icv = 0; icv < number_of_colvars(); icv++) {
1122  os << "delta";
1123  for (size_t icv2 = 0; icv2 < number_of_colvars(); icv2++) {
1124  if (icv == icv2) os << " " << widths[icv];
1125  else os << " " << 0.0;
1126  }
1127  os << "\n";
1128  }
1129 
1130  os << "object 2 class gridconnections counts";
1131  for (icv = 0; icv < number_of_colvars(); icv++) {
1132  os << " " << number_of_points(icv);
1133  }
1134  os << "\n";
1135 
1136  os << "object 3 class array type double rank 0 items "
1137  << number_of_points() << " data follows\n";
1138 
1139  write_raw(os);
1140 
1141  os << "object \"collective variables scalar field\" class field\n";
1142  return os;
1143  }
1144 };
1145 
1146 
1147 
1150 class colvar_grid_count : public colvar_grid<size_t>
1151 {
1152 public:
1153 
1156 
1158  virtual inline ~colvar_grid_count()
1159  {}
1160 
1162  colvar_grid_count(std::vector<int> const &nx_i,
1163  size_t const &def_count = 0);
1164 
1166  colvar_grid_count(std::vector<colvar *> &colvars,
1167  size_t const &def_count = 0);
1168 
1170  inline void incr_count(std::vector<int> const &ix)
1171  {
1172  ++(data[this->address(ix)]);
1173  }
1174 
1176  inline size_t const & new_count(std::vector<int> const &ix,
1177  size_t const &imult = 0)
1178  {
1179  return new_data[address(ix) + imult];
1180  }
1181 
1185  virtual inline void value_input(std::vector<int> const &ix,
1186  size_t const &t,
1187  size_t const &imult = 0,
1188  bool add = false)
1189  {
1190  if (add) {
1191  data[address(ix)] += t;
1192  if (this->has_parent_data) {
1193  // save newly read data for inputting parent grid
1194  new_data[address(ix)] = t;
1195  }
1196  } else {
1197  data[address(ix)] = t;
1198  }
1199  has_data = true;
1200  }
1201 
1204  inline cvm::real log_gradient_finite_diff(const std::vector<int> &ix0,
1205  int n = 0)
1206  {
1207  int A0, A1, A2;
1208  std::vector<int> ix = ix0;
1209 
1210  if (periodic[n]) {
1211  ix[n]--; wrap(ix);
1212  A0 = data[address(ix)];
1213  ix = ix0;
1214  ix[n]++; wrap(ix);
1215  A1 = data[address(ix)];
1216  if (A0 * A1 == 0) {
1217  return 0.; // can't handle empty bins
1218  } else {
1219  return (std::log((cvm::real)A1) - std::log((cvm::real)A0))
1220  / (widths[n] * 2.);
1221  }
1222  } else if (ix[n] > 0 && ix[n] < nx[n]-1) { // not an edge
1223  ix[n]--;
1224  A0 = data[address(ix)];
1225  ix = ix0;
1226  ix[n]++;
1227  A1 = data[address(ix)];
1228  if (A0 * A1 == 0) {
1229  return 0.; // can't handle empty bins
1230  } else {
1231  return (std::log((cvm::real)A1) - std::log((cvm::real)A0))
1232  / (widths[n] * 2.);
1233  }
1234  } else {
1235  // edge: use 2nd order derivative
1236  int increment = (ix[n] == 0 ? 1 : -1);
1237  // move right from left edge, or the other way around
1238  A0 = data[address(ix)];
1239  ix[n] += increment; A1 = data[address(ix)];
1240  ix[n] += increment; A2 = data[address(ix)];
1241  if (A0 * A1 * A2 == 0) {
1242  return 0.; // can't handle empty bins
1243  } else {
1244  return (-1.5 * std::log((cvm::real)A0) + 2. * std::log((cvm::real)A1)
1245  - 0.5 * std::log((cvm::real)A2)) * increment / widths[n];
1246  }
1247  }
1248  }
1249 };
1250 
1251 
1253 class colvar_grid_scalar : public colvar_grid<cvm::real>
1254 {
1255 public:
1256 
1260 
1263 
1266 
1269 
1271  colvar_grid_scalar(std::vector<int> const &nx_i);
1272 
1274  colvar_grid_scalar(std::vector<colvar *> &colvars,
1275  bool margin = 0);
1276 
1278  inline void acc_value(std::vector<int> const &ix,
1279  cvm::real const &new_value,
1280  size_t const &imult = 0)
1281  {
1282  // only legal value of imult here is 0
1283  data[address(ix)] += new_value;
1284  if (samples)
1285  samples->incr_count(ix);
1286  has_data = true;
1287  }
1288 
1290  inline const cvm::real * gradient_finite_diff( const std::vector<int> &ix0 )
1291  {
1292  cvm::real A0, A1;
1293  std::vector<int> ix;
1294  if (nd != 2) {
1295  cvm::error("Finite differences available in dimension 2 only.");
1296  return grad;
1297  }
1298  for (unsigned int n = 0; n < nd; n++) {
1299  ix = ix0;
1300  A0 = data[address(ix)];
1301  ix[n]++; wrap(ix);
1302  A1 = data[address(ix)];
1303  ix[1-n]++; wrap(ix);
1304  A1 += data[address(ix)];
1305  ix[n]--; wrap(ix);
1306  A0 += data[address(ix)];
1307  grad[n] = 0.5 * (A1 - A0) / widths[n];
1308  }
1309  return grad;
1310  }
1311 
1314  virtual cvm::real value_output(std::vector<int> const &ix,
1315  size_t const &imult = 0)
1316  {
1317  if (imult > 0) {
1318  cvm::error("Error: trying to access a component "
1319  "larger than 1 in a scalar data grid.\n");
1320  return 0.;
1321  }
1322  if (samples) {
1323  return (samples->value(ix) > 0) ?
1324  (data[address(ix)] / cvm::real(samples->value(ix))) :
1325  0.0;
1326  } else {
1327  return data[address(ix)];
1328  }
1329  }
1330 
1334  virtual void value_input(std::vector<int> const &ix,
1335  cvm::real const &new_value,
1336  size_t const &imult = 0,
1337  bool add = false)
1338  {
1339  if (imult > 0) {
1340  cvm::error("Error: trying to access a component "
1341  "larger than 1 in a scalar data grid.\n");
1342  return;
1343  }
1344  if (add) {
1345  if (samples)
1346  data[address(ix)] += new_value * samples->new_count(ix);
1347  else
1348  data[address(ix)] += new_value;
1349  } else {
1350  if (samples)
1351  data[address(ix)] = new_value * samples->value(ix);
1352  else
1353  data[address(ix)] = new_value;
1354  }
1355  has_data = true;
1356  }
1357 
1359  cvm::real maximum_value() const;
1360 
1362  cvm::real minimum_value() const;
1363 
1365  cvm::real minimum_pos_value() const;
1366 
1368  cvm::real integral() const;
1369 
1372  cvm::real entropy() const;
1373 
1374 private:
1375  // gradient
1376  cvm::real * grad;
1377 };
1378 
1379 
1380 
1382 class colvar_grid_gradient : public colvar_grid<cvm::real>
1383 {
1384 public:
1385 
1389 
1392 
1394  virtual inline ~colvar_grid_gradient()
1395  {}
1396 
1398  colvar_grid_gradient(std::vector<int> const &nx_i);
1399 
1401  colvar_grid_gradient(std::vector<colvar *> &colvars);
1402 
1404  inline void acc_grad(std::vector<int> const &ix, cvm::real const *grads) {
1405  for (size_t imult = 0; imult < mult; imult++) {
1406  data[address(ix) + imult] += grads[imult];
1407  }
1408  if (samples)
1409  samples->incr_count(ix);
1410  }
1411 
1414  inline void acc_force(std::vector<int> const &ix, cvm::real const *forces) {
1415  for (size_t imult = 0; imult < mult; imult++) {
1416  data[address(ix) + imult] -= forces[imult];
1417  }
1418  if (samples)
1419  samples->incr_count(ix);
1420  }
1421 
1424  virtual inline cvm::real value_output(std::vector<int> const &ix,
1425  size_t const &imult = 0)
1426  {
1427  if (samples)
1428  return (samples->value(ix) > 0) ?
1429  (data[address(ix) + imult] / cvm::real(samples->value(ix))) :
1430  0.0;
1431  else
1432  return data[address(ix) + imult];
1433  }
1434 
1438  virtual inline void value_input(std::vector<int> const &ix,
1439  cvm::real const &new_value,
1440  size_t const &imult = 0,
1441  bool add = false)
1442  {
1443  if (add) {
1444  if (samples)
1445  data[address(ix) + imult] += new_value * samples->new_count(ix);
1446  else
1447  data[address(ix) + imult] += new_value;
1448  } else {
1449  if (samples)
1450  data[address(ix) + imult] = new_value * samples->value(ix);
1451  else
1452  data[address(ix) + imult] = new_value;
1453  }
1454  has_data = true;
1455  }
1456 
1457 
1460  {
1461  size_t n = 0;
1462 
1463  if (nd != 1 || nx[0] == 0) {
1464  return 0.0;
1465  }
1466 
1467  cvm::real sum = 0.0;
1468  std::vector<int> ix = new_index();
1469  if (samples) {
1470  for ( ; index_ok(ix); incr(ix)) {
1471  if ( (n = samples->value(ix)) )
1472  sum += value(ix) / n;
1473  }
1474  } else {
1475  for ( ; index_ok(ix); incr(ix)) {
1476  sum += value(ix);
1477  }
1478  }
1479  return (sum / cvm::real(nx[0]));
1480  }
1481 
1484  void write_1D_integral(std::ostream &os);
1485 
1486 };
1487 
1488 
1489 #endif
1490 
std::ostream & write_params(std::ostream &os)
Write the grid parameters (number of colvars, boundaries, width and number of points) ...
Definition: colvargrid.h:745
Scalar number, implemented as colvarmodule::real (default)
Definition: colvarvalue.h:47
Like parse_normal, but don't send any message to the log (useful e.g. in restart files when such mess...
Definition: colvarparse.h:100
std::istream & read_restart(std::istream &is)
Read grid entry in restart file.
Definition: colvargrid.h:888
std::ostream & write_opendx(std::ostream &os)
Write the grid data without labels, as they are represented in memory.
Definition: colvargrid.h:1105
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:228
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:1314
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:429
cvm::real minimum_value() const
Return the lowest value.
Definition: colvargrid.cpp:67
const cvm::real * gradient_finite_diff(const std::vector< int > &ix0)
Return the gradient of the scalar field from finite differences.
Definition: colvargrid.h:1290
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:682
cvm::real average()
Compute and return average value for a 1D gradient grid.
Definition: colvargrid.h:1459
std::vector< int > const & sizes() const
Get the sizes in each direction.
Definition: colvargrid.h:111
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:496
size_t address(std::vector< int > const &ix) const
Get the low-level index corresponding to an index.
Definition: colvargrid.h:52
size_t mult
Multiplicity of each datum (allow the binning of non-scalar types such as atomic gradients) ...
Definition: colvargrid.h:34
size_t number_of_points(int const icv=-1) const
Definition: colvargrid.h:101
colvar_grid_count * samples
Provide the associated sample count by which each binned value should be divided. ...
Definition: colvargrid.h:1259
cvm::real real_value
Real data member.
Definition: colvarvalue.h:68
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:398
std::vector< bool > hard_upper_boundaries
Whether some colvars have hard upper boundaries.
Definition: colvargrid.h:82
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:389
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:1334
virtual ~colvar_grid()
Destructor.
Definition: colvargrid.h:201
virtual ~colvar_grid_count()
Destructor.
Definition: colvargrid.h:1158
bool index_ok(std::vector< int > const &ix) const
Check that the index is within range in each of the dimensions.
Definition: colvargrid.h:709
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:1204
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:865
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:473
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:1008
~colvar_grid_scalar()
Destructor.
Definition: colvargrid.cpp:49
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:238
Grid of values of a function of several collective variables.
Definition: colvargrid.h:19
colvar_grid(colvar_grid< T > const &g)
"Almost copy-constructor": only copies configuration parameters from another grid, but doesn't reallocate stuff; setup() must be called after that;
Definition: colvargrid.h:207
colvar_grid_count()
Default constructor.
Definition: colvargrid.cpp:11
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:966
cvm::real minimum_pos_value() const
Return the lowest positive value.
Definition: colvargrid.cpp:76
std::vector< bool > actual_value
Do we request actual value (for extended-system colvars)?
Definition: colvargrid.h:49
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:702
colvar_grid_scalar()
Default constructor.
Definition: colvargrid.cpp:27
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:85
cvm::real entropy() const
Assuming that the map is a normalized probability density, calculates the entropy (uses widths if the...
Definition: colvargrid.cpp:106
void raw_data_in(const T *in_data)
Input the data as they are represented in memory.
Definition: colvargrid.h:501
std::vector< int > const get_colvars_index() const
Get the bin indices corresponding to the current values of the colvars.
Definition: colvargrid.h:555
void acc_grad(std::vector< int > const &ix, cvm::real const *grads)
Accumulate the gradient.
Definition: colvargrid.h:1404
bool save_delimiters
Whether or not to accumulate data_begin_pos and data_end_pos in key_lookup(); it may be useful to dis...
Definition: colvarparse.h:43
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:720
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:916
int current_bin_scalar(int const i) const
Report the bin corresponding to the current value of variable i.
Definition: colvargrid.h:376
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:512
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:775
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
cvm::real integral() const
Calculates the integral of the map (uses widths if they are defined)
Definition: colvargrid.cpp:92
void check_consistency()
Check that the grid information inside (boundaries, widths, ...) is consistent with the current setti...
Definition: colvargrid.h:846
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:423
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:1424
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:383
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:414
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:437
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:1176
int setup()
Allocate data (allow initialization also after construction)
Definition: colvargrid.h:179
cvm::real maximum_value() const
Return the highest value.
Definition: colvargrid.cpp:57
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:1438
colvar_grid_gradient()
Default constructor.
Definition: colvargrid.cpp:120
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:263
size_t nt
Total number of grid points.
Definition: colvargrid.h:37
virtual ~colvar_grid_gradient()
Destructor.
Definition: colvargrid.h:1394
std::ostream & write_restart(std::ostream &os)
Write grid entry in restart file.
Definition: colvargrid.h:905
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:566
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:673
void map_grid(colvar_grid< T > const &other_grid)
Add data from another grid of the same type.
Definition: colvargrid.h:608
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:1388
(default) Read the first instance of a keyword (if any), report its value, and print a warning when t...
Definition: colvarparse.h:96
size_t raw_data_num() const
Size of the data as they are represented in memory.
Definition: colvargrid.h:507
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:1414
static void error(std::string const &message, int code=COLVARS_ERROR)
Print a message to the main log and set global error code.
Definition: colvarmodule.cpp:1536
void acc_value(std::vector< int > const &ix, cvm::real const &new_value, size_t const &imult=0)
Accumulate the value.
Definition: colvargrid.h:1278
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:651
std::istream & read_raw(std::istream &is)
Read data written by colvar_grid::write_raw()
Definition: colvargrid.h:941
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:615
Class for accumulating a scalar function on a grid.
Definition: colvargrid.h:1253
void add_constant(T const &t)
Add a constant to all elements (fast loop)
Definition: colvargrid.h:520
void wrap(std::vector< int > &ix) const
Definition: colvargrid.h:359
Colvar_grid derived class to hold counters in discrete n-dim colvar space.
Definition: colvargrid.h:1150
static size_t const cv_prec
Number of digits to represent a collective variables value(s)
Definition: colvarmodule.h:415
Parse_Mode
How a keyword is parsed in a string.
Definition: colvarparse.h:92
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:1185
size_t multiplicity() const
Return the multiplicity of the type used.
Definition: colvargrid.h:123
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
size_t number_of_colvars() const
Return the number of colvars.
Definition: colvargrid.h:94
static void log(std::string const &message)
Print a message to the main log.
Definition: colvarmodule.cpp:1467
void write_1D_integral(std::ostream &os)
If the grid is 1-dimensional, integrate it and write the integral to a file.
Definition: colvargrid.cpp:132
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:1170
void multiply_constant(cvm::real const &a)
Multiply all elements by a scalar constant (fast loop)
Definition: colvargrid.h:528
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:405
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:449
std::vector< cvm::real > widths
Widths of the colvars in this grid.
Definition: colvargrid.h:85
Class for accumulating the gradient of a scalar function on a grid.
Definition: colvargrid.h:1382
void remove_zeros(cvm::real const &a)
Assign all zero elements a scalar constant (fast loop)
Definition: colvargrid.h:535
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:417
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:544
static bool debug()
Whether debug output should be enabled (compile-time option)
Definition: colvarmodule.h:225
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:578