Collective Variables Module - Developer Documentation
Loading...
Searching...
No Matches
colvargrid.h
1// -*- c++ -*-
2
3// This file is part of the Collective Variables module (Colvars).
4// The original version of Colvars and its updates are located at:
5// https://github.com/Colvars/colvars
6// Please update all Colvars source files before making any changes.
7// If you wish to distribute your changes, please submit them to the
8// Colvars repository at GitHub.
9
10#ifndef COLVARGRID_H
11#define COLVARGRID_H
12
13#include <iosfwd>
14#include <memory>
15
16#include "colvar.h"
17#include "colvarmodule.h"
18#include "colvarvalue.h"
19#include "colvarparse.h"
20
21
25
26public:
28 size_t nd = 0;
29
31 std::vector<int> nx;
32
34 std::vector<int> nxc;
35
37 std::vector<colvarvalue> lower_boundaries;
38
40 std::vector<colvarvalue> upper_boundaries;
41
43 std::vector<cvm::real> widths;
44};
45
46
52template <class T> class colvar_grid : public colvar_grid_params, public colvarparse {
53
54 //protected:
55public: // TODO create accessors for these after all instantiations work
56
59 size_t mult;
60
62 size_t nt;
63
65 std::vector<T> data;
66
68 std::vector<size_t> new_data;
69
71 std::vector<colvar *> cv;
72
74 std::vector<bool> use_actual_value;
75
77 inline size_t address(std::vector<int> const &ix) const
78 {
79 size_t addr = 0;
80 for (size_t i = 0; i < nd; i++) {
81 addr += ix[i]*static_cast<size_t>(nxc[i]);
82 if (cvm::debug()) {
83 if (ix[i] >= nx[i]) {
84 cvm::error("Error: exceeding bounds in colvar_grid.\n", COLVARS_BUG_ERROR);
85 return 0;
86 }
87 }
88 }
89 return addr;
90 }
91
92public:
94 std::vector<bool> periodic;
95
97 std::vector<bool> hard_lower_boundaries;
98
100 std::vector<bool> hard_upper_boundaries;
101
104
107
109 inline size_t num_variables() const
110 {
111 return nd;
112 }
113
115 inline std::vector<int> const &number_of_points_vec() const
116 {
117 return nx;
118 }
119
122 inline size_t number_of_points(int const icv = -1) const
123 {
124 if (icv < 0) {
125 return nt;
126 } else {
127 return nx[icv];
128 }
129 }
130
132 inline std::vector<int> const & sizes() const
133 {
134 return nx;
135 }
136
138 inline void set_sizes(std::vector<int> const &new_sizes)
139 {
140 nx = new_sizes;
141 }
142
144 inline size_t multiplicity() const
145 {
146 return mult;
147 }
148
150 inline void request_actual_value(bool b = true)
151 {
152 size_t i;
153 for (i = 0; i < use_actual_value.size(); i++) {
154 use_actual_value[i] = b;
155 }
156 }
157
159 int setup(std::vector<int> const &nx_i,
160 T const &t = T(),
161 size_t const &mult_i = 1)
162 {
163 if (cvm::debug()) {
164 cvm::log("Allocating grid: multiplicity = "+cvm::to_str(mult_i)+
165 ", dimensions = "+cvm::to_str(nx_i)+".\n");
166 }
167
168 mult = mult_i;
169
170 data.clear();
171
172 nx = nx_i;
173 nd = nx.size();
174
175 nxc.resize(nd);
176
177 // setup dimensions
178 nt = mult;
179 for (int i = nd-1; i >= 0; i--) {
180 if (nx[i] <= 0) {
181 cvm::error("Error: providing an invalid number of grid points, "+
182 cvm::to_str(nx[i])+".\n", COLVARS_BUG_ERROR);
183 return COLVARS_ERROR;
184 }
185 nxc[i] = nt;
186 nt *= nx[i];
187 }
188
189 if (cvm::debug()) {
190 cvm::log("Total number of grid elements = "+cvm::to_str(nt)+".\n");
191 }
192
193 data.reserve(nt);
194 data.assign(nt, t);
195
196 return COLVARS_OK;
197 }
198
200 int setup()
201 {
202 return setup(this->nx, T(), this->mult);
203 }
204
206 void reset(T const &t = T())
207 {
208 data.assign(nt, t);
209 }
210
211
214 {
215 nd = nt = 0;
216 mult = 1;
217 has_parent_data = false;
218 this->setup();
219 }
220
222 virtual ~colvar_grid()
223 {}
224
229 colvarparse(),
230 mult(g.mult),
231 data(),
232 cv(g.cv),
237 has_parent_data(false),
238 has_data(false)
239 {}
240
245 colvar_grid(std::vector<int> const &nx_i,
246 T const &t = T(),
247 size_t mult_i = 1)
248 : has_parent_data(false), has_data(false)
249 {
250 this->setup(nx_i, t, mult_i);
251 }
252
256 colvar_grid(std::vector<colvar *> const &colvars,
257 T const &t = T(),
258 size_t mult_i = 1,
259 bool add_extra_bin = false,
260 std::shared_ptr<const colvar_grid_params> params = nullptr,
261 std::string config = std::string())
262 : has_parent_data(false), has_data(false)
263 {
264 (void) t;
265 this->init_from_colvars(colvars, mult_i, add_extra_bin, params, config);
266 }
267
271 colvar_grid(std::string const &filename, size_t mult_i = 1);
272
273 int init_from_colvars(std::vector<colvar *> const &colvars,
274 size_t mult_i = 1,
275 bool add_extra_bin = false,
276 std::shared_ptr<const colvar_grid_params> params = nullptr,
277 std::string config = std::string())
278 {
279 if (cvm::debug()) {
280 cvm::log("Reading grid configuration from collective variables.\n");
281 }
282
283 cv = colvars;
284 nd = colvars.size();
285 mult = mult_i;
286
287 size_t i;
288
289 if (cvm::debug()) {
290 cvm::log("Allocating a grid for "+cvm::to_str(colvars.size())+
291 " collective variables, multiplicity = "+cvm::to_str(mult_i)+".\n");
292 }
293
294 for (i = 0; i < nd; i++) {
295 if (cv[i]->value().type() != colvarvalue::type_scalar) {
296 cvm::error("Colvar grids can only be automatically "
297 "constructed for scalar variables. "
298 "ABF and histogram can not be used; metadynamics "
299 "can be used with useGrids disabled.\n", COLVARS_INPUT_ERROR);
300 return COLVARS_ERROR;
301 }
302
303 if (cv[i]->width <= 0.0) {
304 cvm::error("Tried to initialize a grid on a "
305 "variable with negative or zero width.\n", COLVARS_INPUT_ERROR);
306 return COLVARS_ERROR;
307 }
308
309 widths.push_back(cv[i]->width);
312
313 // By default, get reported colvar value (for extended Lagrangian colvars)
314 use_actual_value.push_back(false);
315
316 // except if a colvar is specified twice in a row
317 // then the first instance is the actual value
318 // For histograms of extended-system coordinates
319 if (i > 0 && cv[i-1] == cv[i]) {
320 use_actual_value[i-1] = true;
321 }
322
323 // This needs to work if the boundaries are undefined in the colvars
324 lower_boundaries.push_back(cv[i]->lower_boundary);
325 upper_boundaries.push_back(cv[i]->upper_boundary);
326 }
327
328 // Replace widths and boundaries with optional custom configuration
329 if (!config.empty()) {
330 this->parse_params(config);
331 this->check_keywords(config, "grid");
332
333 if (params) {
334 cvm::error("Error: init_from_colvars was passed both a grid config and a template grid.", COLVARS_BUG_ERROR);
335 return COLVARS_BUG_ERROR;
336 }
337 } else if (params) {
338 // Match grid sizes with template
339
340 if (params->nd != nd) {
341 cvm::error("Trying to initialize grid from template with wrong dimension (" +
342 cvm::to_str(params->nd) + " instead of " +
343 cvm::to_str(this->nd) + ").");
344 return COLVARS_ERROR;
345 }
346
347 widths =params->widths;
348 lower_boundaries =params->lower_boundaries;
349 upper_boundaries =params->upper_boundaries;
350 }
351
352 // Only now can we determine periodicity
353 for (i = 0; i < nd; i++) {
354 periodic.push_back(cv[i]->periodic_boundaries(lower_boundaries[i].real_value,
355 upper_boundaries[i].real_value));
356
357 if (add_extra_bin) {
358 // Shift the grid by half the bin width (values at edges instead of center of bins)
359 lower_boundaries[i] -= 0.5 * widths[i];
360
361 if (periodic[i]) {
362 // Just shift
363 upper_boundaries[i] -= 0.5 * widths[i];
364 } else {
365 // Widen grid by one bin width
366 upper_boundaries[i] += 0.5 * widths[i];
367 }
368 }
369 }
370
371 // Reset grid sizes based on widths and boundaries
372 this->init_from_boundaries();
373 return this->setup();
374 }
375
376 int init_from_boundaries()
377 {
378 if (cvm::debug()) {
379 cvm::log("Configuring grid dimensions from colvars boundaries.\n");
380 }
381
382 // these will have to be updated
383 nx.clear();
384 nxc.clear();
385 nt = 0;
386
387 for (size_t i = 0; i < lower_boundaries.size(); i++) {
388 // Re-compute periodicity using current grid boundaries
389 periodic[i] = cv[i]->periodic_boundaries(lower_boundaries[i].real_value,
390 upper_boundaries[i].real_value);
391
392 cvm::real nbins = ( upper_boundaries[i].real_value -
393 lower_boundaries[i].real_value ) / widths[i];
394 int nbins_round = (int)(nbins+0.5);
395
396 if (cvm::fabs(nbins - cvm::real(nbins_round)) > 1.0E-10) {
397 cvm::log("Warning: grid interval("+
400 ") is not commensurate to its bin width("+
402 upper_boundaries[i].real_value = lower_boundaries[i].real_value +
403 (nbins_round * widths[i]);
404 }
405
406 if (cvm::debug())
407 cvm::log("Number of points is "+cvm::to_str((int) nbins_round)+
408 " for the colvar no. "+cvm::to_str(i+1)+".\n");
409
410 nx.push_back(nbins_round);
411 }
412
413 return COLVARS_OK;
414 }
415
418 inline void wrap(std::vector<int> & ix) const
419 {
420 for (size_t i = 0; i < nd; i++) {
421 if (periodic[i]) {
422 ix[i] = (ix[i] + nx[i]) % nx[i]; // Avoid modulo with negative operands (implementation-defined)
423 } else {
424 if (ix[i] < 0 || ix[i] >= nx[i]) {
425 cvm::error("Trying to wrap illegal index vector (non-PBC) for a grid point: "
426 + cvm::to_str(ix), COLVARS_BUG_ERROR);
427 return;
428 }
429 }
430 }
431 }
432
435 inline bool wrap_detect_edge(std::vector<int> & ix) const
436 {
437 bool edge = false;
438 for (size_t i = 0; i < nd; i++) {
439 if (periodic[i]) {
440 ix[i] = (ix[i] + nx[i]) % nx[i]; // Avoid modulo with negative operands (implementation-defined)
441 } else if (ix[i] < 0 || ix[i] >= nx[i]) {
442 edge = true;
443 }
444 }
445 return edge;
446 }
447
450 inline bool wrap_to_edge(std::vector<int> & ix, std::vector<int> & edge_bin) const
451 {
452 bool edge = false;
453 edge_bin = ix;
454 for (size_t i = 0; i < nd; i++) {
455 if (periodic[i]) {
456 ix[i] = (ix[i] + nx[i]) % nx[i]; // Avoid modulo with negative operands (implementation-defined)
457 edge_bin[i] = ix[i];
458 } else if (ix[i] < 0) {
459 edge = true;
460 edge_bin[i] = 0;
461 } else if (ix[i] >= nx[i]) {
462 edge = true;
463 edge_bin[i] = nx[i] - 1;
464 }
465 }
466 return edge;
467 }
468
469
471 inline int current_bin_scalar(int const i) const
472 {
473 return value_to_bin_scalar(use_actual_value[i] ? cv[i]->actual_value() : cv[i]->value(), i);
474 }
475
478 inline int current_bin_flat_bound() const
479 {
480 std::vector<int> index = new_index();
481 for (size_t i = 0; i < nd; i++) {
482 index[i] = current_bin_scalar_bound(i);
483 }
484 return address(index);
485 }
486
489 inline int current_bin_scalar_bound(int const i) const
490 {
491 return value_to_bin_scalar_bound(use_actual_value[i] ? cv[i]->actual_value() : cv[i]->value(), i);
492 }
493
495 inline int current_bin_scalar(int const i, int const iv) const
496 {
498 cv[i]->actual_value().vector1d_value[iv] :
499 cv[i]->value().vector1d_value[iv], i);
500 }
501
504 inline int value_to_bin_scalar(colvarvalue const &value, const int i) const
505 {
506 return (int) cvm::floor( (value.real_value - lower_boundaries[i].real_value) / widths[i] );
507 }
508
510 inline cvm::real current_bin_scalar_fraction(int const i) const
511 {
512 return value_to_bin_scalar_fraction(use_actual_value[i] ? cv[i]->actual_value() : cv[i]->value(), i);
513 }
514
517 inline cvm::real value_to_bin_scalar_fraction(colvarvalue const &value, const int i) const
518 {
519 cvm::real x = (value.real_value - lower_boundaries[i].real_value) / widths[i];
520 return x - cvm::floor(x);
521 }
522
525 inline int value_to_bin_scalar_bound(colvarvalue const &value, const int i) const
526 {
527 int bin_index = cvm::floor( (value.real_value - lower_boundaries[i].real_value) / widths[i] );
528
529 // Wrap bins for periodic dimensions before truncating
530 if (periodic[i]) bin_index %= nx[i];
531 if (bin_index < 0) bin_index=0;
532 if (bin_index >=int(nx[i])) bin_index=int(nx[i])-1;
533 return (int) bin_index;
534 }
535
537 inline int value_to_bin_scalar(colvarvalue const &value,
538 colvarvalue const &new_offset,
539 cvm::real const &new_width) const
540 {
541 return (int) cvm::floor( (value.real_value - new_offset.real_value) / new_width );
542 }
543
546 inline colvarvalue bin_to_value_scalar(int const &i_bin, int const i) const
547 {
548 return lower_boundaries[i].real_value + widths[i] * (0.5 + i_bin);
549 }
550
552 inline colvarvalue bin_to_value_scalar(int const &i_bin,
553 colvarvalue const &new_offset,
554 cvm::real const &new_width) const
555 {
556 return new_offset.real_value + new_width * (0.5 + i_bin);
557 }
558
560 inline void set_value(std::vector<int> const &ix,
561 T const &t,
562 size_t const &imult = 0)
563 {
564 data[this->address(ix)+imult] = t;
565 has_data = true;
566 }
567
569 inline void set_value(size_t i, T const &t)
570 {
571 data[i] = t;
572 }
573
575 inline T get_value(size_t i) const
576 {
577 return data[i];
578 }
579
580
585 void delta_grid(colvar_grid<T> const &other_grid)
586 {
587
588 if (other_grid.multiplicity() != this->multiplicity()) {
589 cvm::error("Error: trying to subtract two grids with "
590 "different multiplicity.\n");
591 return;
592 }
593
594 if (other_grid.data.size() != this->data.size()) {
595 cvm::error("Error: trying to subtract two grids with "
596 "different size.\n");
597 return;
598 }
599
600 for (size_t i = 0; i < data.size(); i++) {
601 data[i] = other_grid.data[i] - data[i];
602 }
603 has_data = true;
604 }
605
606
610 void copy_grid(colvar_grid<T> const &other_grid)
611 {
612 if (other_grid.multiplicity() != this->multiplicity()) {
613 cvm::error("Error: trying to copy two grids with "
614 "different multiplicity.\n");
615 return;
616 }
617
618 if (other_grid.data.size() != this->data.size()) {
619 cvm::error("Error: trying to copy two grids with "
620 "different size.\n");
621 return;
622 }
623
624
625 for (size_t i = 0; i < data.size(); i++) {
626 data[i] = other_grid.data[i];
627 }
628 has_data = true;
629 }
630
633 void raw_data_out(T* out_data) const
634 {
635 for (size_t i = 0; i < data.size(); i++) out_data[i] = data[i];
636 }
637 void raw_data_out(std::vector<T>& out_data) const
638 {
639 out_data = data;
640 }
642 void raw_data_in(const T* in_data)
643 {
644 for (size_t i = 0; i < data.size(); i++) data[i] = in_data[i];
645 has_data = true;
646 }
647 void raw_data_in(const std::vector<T>& in_data)
648 {
649 data = in_data;
650 has_data = true;
651 }
653 size_t raw_data_num() const { return data.size(); }
654
655
658 inline T const & value(std::vector<int> const &ix,
659 size_t const &imult = 0) const
660 {
661 return data[this->address(ix) + imult];
662 }
663
665 inline T const & value(size_t i) const
666 {
667 return data[i];
668 }
669
671 inline void add_constant(T const &t)
672 {
673 for (size_t i = 0; i < nt; i++)
674 data[i] += t;
675 has_data = true;
676 }
677
679 inline void multiply_constant(cvm::real const &a)
680 {
681 for (size_t i = 0; i < nt; i++)
682 data[i] *= a;
683 }
684
686 inline void remove_small_values(cvm::real const &a)
687 {
688 for (size_t i = 0; i < nt; i++)
689 if(data[i]<a) data[i] = a;
690 }
691
692
695 inline std::vector<int> const get_colvars_index(std::vector<colvarvalue> const &values) const
696 {
697 std::vector<int> index = new_index();
698 for (size_t i = 0; i < nd; i++) {
699 index[i] = value_to_bin_scalar(values[i], i);
700 }
701 return index;
702 }
703
706 inline std::vector<int> const get_colvars_index() const
707 {
708 std::vector<int> index = new_index();
709 for (size_t i = 0; i < nd; i++) {
710 index[i] = current_bin_scalar(i);
711 }
712 return index;
713 }
714
717 inline std::vector<int> const get_colvars_index_bound() const
718 {
719 std::vector<int> index = new_index();
720 for (size_t i = 0; i < nd; i++) {
721 index[i] = current_bin_scalar_bound(i);
722 }
723 return index;
724 }
725
729 inline cvm::real bin_distance_from_boundaries(std::vector<colvarvalue> const &values,
730 bool skip_hard_boundaries = false)
731 {
732 cvm::real minimum = 1.0E+16;
733 for (size_t i = 0; i < nd; i++) {
734
735 if (periodic[i]) continue;
736
737 cvm::real dl = cvm::sqrt(cv[i]->dist2(values[i], lower_boundaries[i])) / widths[i];
738 cvm::real du = cvm::sqrt(cv[i]->dist2(values[i], upper_boundaries[i])) / widths[i];
739
740 if (values[i].real_value < lower_boundaries[i])
741 dl *= -1.0;
742 if (values[i].real_value > upper_boundaries[i])
743 du *= -1.0;
744
745 if ( ((!skip_hard_boundaries) || (!hard_lower_boundaries[i])) && (dl < minimum))
746 minimum = dl;
747 if ( ((!skip_hard_boundaries) || (!hard_upper_boundaries[i])) && (du < minimum))
748 minimum = du;
749 }
750
751 return minimum;
752 }
753
754
759 void map_grid(colvar_grid<T> const &other_grid)
760 {
761 if (other_grid.multiplicity() != this->multiplicity()) {
762 cvm::error("Error: trying to merge two grids with values of "
763 "different multiplicity.\n");
764 return;
765 }
766
767 std::vector<colvarvalue> const &gb = this->lower_boundaries;
768 std::vector<cvm::real> const &gw = this->widths;
769 std::vector<colvarvalue> const &ogb = other_grid.lower_boundaries;
770 std::vector<cvm::real> const &ogw = other_grid.widths;
771
772 std::vector<int> ix = this->new_index();
773 std::vector<int> oix = other_grid.new_index();
774
775 if (cvm::debug())
776 cvm::log("Remapping grid...\n");
777 for ( ; this->index_ok(ix); this->incr(ix)) {
778
779 for (size_t i = 0; i < nd; i++) {
780 oix[i] =
781 value_to_bin_scalar(bin_to_value_scalar(ix[i], gb[i], gw[i]),
782 ogb[i],
783 ogw[i]);
784 }
785
786 if (! other_grid.index_ok(oix)) {
787 continue;
788 }
789
790 for (size_t im = 0; im < mult; im++) {
791 this->set_value(ix, other_grid.value(oix, im), im);
792 }
793 }
794
795 has_data = true;
796 if (cvm::debug())
797 cvm::log("Remapping done.\n");
798 }
799
802 void add_grid(colvar_grid<T> const &other_grid,
803 cvm::real scale_factor = 1.0)
804 {
805 if (other_grid.multiplicity() != this->multiplicity()) {
806 cvm::error("Error: trying to sum togetehr two grids with values of "
807 "different multiplicity.\n");
808 return;
809 }
810 if (scale_factor != 1.0)
811 for (size_t i = 0; i < data.size(); i++) {
812 data[i] += static_cast<T>(scale_factor * other_grid.data[i]);
813 }
814 else
815 // skip multiplication if possible
816 for (size_t i = 0; i < data.size(); i++) {
817 data[i] += other_grid.data[i];
818 }
819 has_data = true;
820 }
821
824 virtual T value_output(std::vector<int> const &ix,
825 size_t const &imult = 0) const
826 {
827 return value(ix, imult);
828 }
829
833 virtual void value_input(std::vector<int> const &ix,
834 T const &t,
835 size_t const &imult = 0,
836 bool add = false)
837 {
838 if ( add )
839 data[address(ix) + imult] += t;
840 else
841 data[address(ix) + imult] = t;
842 has_data = true;
843 }
844
845
846 // /// Get the pointer to the binned value indexed by ix
847 // inline T const *value_p (std::vector<int> const &ix)
848 // {
849 // return &(data[address (ix)]);
850 // }
851
854 inline std::vector<int> const new_index() const
855 {
856 return std::vector<int> (nd, 0);
857 }
858
861 inline bool index_ok(std::vector<int> const &ix) const
862 {
863 for (size_t i = 0; i < nd; i++) {
864 if ( (ix[i] < 0) || (ix[i] >= int(nx[i])) )
865 return false;
866 }
867 return true;
868 }
869
872 inline void incr(std::vector<int> &ix) const
873 {
874 for (int i = ix.size()-1; i >= 0; i--) {
875
876 ix[i]++;
877
878 if (ix[i] >= nx[i]) {
879
880 if (i > 0) {
881 ix[i] = 0;
882 continue;
883 } else {
884 // this is the last iteration, a non-valid index is being
885 // set for the outer index, which will be caught by
886 // index_ok()
887 ix[0] = nx[0];
888 return;
889 }
890 } else {
891 return;
892 }
893 }
894 }
895
897 std::string get_state_params() const;
898
900 int parse_params(std::string const &conf,
902
907 {
908 for (size_t i = 0; i < nd; i++) {
909 if ( (cvm::sqrt(cv[i]->dist2(cv[i]->lower_boundary,
910 lower_boundaries[i])) > 1.0E-10) ||
911 (cvm::sqrt(cv[i]->dist2(cv[i]->upper_boundary,
912 upper_boundaries[i])) > 1.0E-10) ||
913 (cvm::sqrt(cv[i]->dist2(cv[i]->width,
914 widths[i])) > 1.0E-10) ) {
915 cvm::error("Error: restart information for a grid is "
916 "inconsistent with that of its colvars.\n");
917 return;
918 }
919 }
920 }
921
922
925 void check_consistency(colvar_grid<T> const &other_grid)
926 {
927 for (size_t i = 0; i < nd; i++) {
928 // we skip dist2(), because periodicities and the like should
929 // matter: boundaries should be EXACTLY the same (otherwise,
930 // map_grid() should be used)
931 if ( (cvm::fabs(other_grid.lower_boundaries[i] -
932 lower_boundaries[i]) > 1.0E-10) ||
933 (cvm::fabs(other_grid.upper_boundaries[i] -
934 upper_boundaries[i]) > 1.0E-10) ||
935 (cvm::fabs(other_grid.widths[i] -
936 widths[i]) > 1.0E-10) ||
937 (data.size() != other_grid.data.size()) ) {
938 cvm::error("Error: inconsistency between "
939 "two grids that are supposed to be equal, "
940 "aside from the data stored.\n");
941 return;
942 }
943 }
944 }
945
947 std::istream & read_restart(std::istream &is);
948
951
953 std::ostream & write_restart(std::ostream &os);
954
957
959 std::istream &read_raw(std::istream &is);
960
963
967 std::ostream &write_raw(std::ostream &os, size_t const buf_size = 3) const;
968
972 cvm::memory_stream &write_raw(cvm::memory_stream &os, size_t const buf_size = 3) const;
973
975 std::istream & read_multicol(std::istream &is, bool add = false);
976
978 int read_multicol(std::string const &filename,
979 std::string description = "grid file",
980 bool add = false);
981
983 std::ostream & write_multicol(std::ostream &os) const;
984
986 int write_multicol(std::string const &filename,
987 std::string description = "grid file") const;
988
990 std::ostream & write_opendx(std::ostream &os) const;
991
993 int write_opendx(std::string const &filename,
994 std::string description = "grid file") const;
995};
996
997
998
1001class colvar_grid_count : public colvar_grid<size_t>
1002{
1003public:
1004
1007
1010 {}
1011
1013 colvar_grid_count(std::vector<colvar *> &colvars,
1014 std::shared_ptr<const colvar_grid_params> params = nullptr);
1015
1016 colvar_grid_count(std::vector<colvar *> &colvars,
1017 std::string config);
1018
1020 inline void incr_count(std::vector<int> const &ix)
1021 {
1022 ++(data[this->address(ix)]);
1023 }
1024
1026 inline size_t const & new_value(std::vector<int> const &ix)
1027 {
1028 return new_data[address(ix)];
1029 }
1030
1032 std::string get_state_params() const;
1033
1035 int parse_params(std::string const &conf,
1037
1039 std::istream & read_restart(std::istream &is);
1040
1043
1045 std::ostream & write_restart(std::ostream &os);
1046
1049
1051 std::istream &read_raw(std::istream &is);
1052
1055
1059 std::ostream &write_raw(std::ostream &os, size_t const buf_size = 3) const;
1060
1064 cvm::memory_stream &write_raw(cvm::memory_stream &os, size_t const buf_size = 3) const;
1065
1067 std::istream & read_multicol(std::istream &is, bool add = false);
1068
1070 int read_multicol(std::string const &filename,
1071 std::string description = "grid file",
1072 bool add = false);
1073
1075 std::ostream & write_multicol(std::ostream &os) const;
1076
1078 int write_multicol(std::string const &filename,
1079 std::string description = "grid file") const;
1080
1082 std::ostream & write_opendx(std::ostream &os) const;
1083
1085 int write_opendx(std::string const &filename,
1086 std::string description = "grid file") const;
1087
1089 virtual void value_input(std::vector<int> const &ix,
1090 size_t const &t,
1091 size_t const &imult = 0,
1092 bool add = false)
1093 {
1094 (void) imult;
1095 if (add) {
1096 data[address(ix)] += t;
1097 if (this->has_parent_data) {
1098 // save newly read data for inputting parent grid
1099 new_data[address(ix)] = t;
1100 }
1101 } else {
1102 data[address(ix)] = t;
1103 }
1104 has_data = true;
1105 }
1106
1109 inline int local_sample_count(int radius)
1110 {
1111 std::vector<int> ix0 = new_index();
1112 std::vector<int> ix = new_index();
1113
1114 for (size_t i = 0; i < nd; i++) {
1115 ix0[i] = current_bin_scalar_bound(i);
1116 }
1117 if (radius < 1) {
1118 // Simple case: no averaging
1119 if (index_ok(ix0))
1120 return value(ix0);
1121 else
1122 return 0;
1123 }
1124 size_t count = 0;
1125 size_t nbins = 0;
1126 int i, j, k;
1127 bool edge;
1128 ix = ix0;
1129 // Treat each dimension separately to simplify code
1130 switch (nd)
1131 {
1132 case 1:
1133 for (i = -radius; i <= radius; i++) {
1134 ix[0] = ix0[0] + i;
1135 edge = wrap_detect_edge(ix);
1136 if (!edge) {
1137 nbins++;
1138 count += value(ix);
1139 }
1140 }
1141 break;
1142 case 2:
1143 for (i = -radius; i <= radius; i++) {
1144 ix[0] = ix0[0] + i;
1145 for (j = -radius; j <= radius; j++) {
1146 ix[1] = ix0[1] + j;
1147 edge = wrap_detect_edge(ix);
1148 if (!edge) {
1149 nbins++;
1150 count += value(ix);
1151 }
1152 }
1153 }
1154 break;
1155 case 3:
1156 for (i = -radius; i <= radius; i++) {
1157 ix[0] = ix0[0] + i;
1158 for (j = -radius; j <= radius; j++) {
1159 ix[1] = ix0[1] + j;
1160 for (k = -radius; k <= radius; k++) {
1161 ix[2] = ix0[2] + k;
1162 edge = wrap_detect_edge(ix);
1163 if (!edge) {
1164 nbins++;
1165 count += value(ix);
1166 }
1167 }
1168 }
1169 }
1170 break;
1171 default:
1172 cvm::error("Error: local_sample_count is not implemented for grids of dimension > 3", COLVARS_NOT_IMPLEMENTED);
1173 break;
1174 }
1175
1176 if (nbins)
1177 // Integer division - an error on the order of 1 doesn't matter
1178 return count / nbins;
1179 else
1180 return 0.0;
1181 }
1182
1183
1187 inline cvm::real log_gradient_finite_diff(const std::vector<int> &ix0,
1188 int n = 0, int offset = 0)
1189 {
1190 cvm::real A0, A1, A2;
1191 std::vector<int> ix = ix0;
1192
1193 // TODO this can be rewritten more concisely with wrap_edge()
1194 if (periodic[n]) {
1195 ix[n]--; wrap(ix);
1196 A0 = value(ix) + offset;
1197 ix = ix0;
1198 ix[n]++; wrap(ix);
1199 A1 = value(ix) + offset;
1200 if (A0 * A1 == 0) {
1201 return 0.; // can't handle empty bins
1202 } else {
1203 return (cvm::logn(A1) - cvm::logn(A0))
1204 / (widths[n] * 2.);
1205 }
1206 } else if (ix[n] > 0 && ix[n] < nx[n]-1) { // not an edge
1207 ix[n]--;
1208 A0 = value(ix) + offset;
1209 ix = ix0;
1210 ix[n]++;
1211 A1 = value(ix) + offset;
1212 if (A0 * A1 == 0) {
1213 return 0.; // can't handle empty bins
1214 } else {
1215 return (cvm::logn(A1) - cvm::logn(A0))
1216 / (widths[n] * 2.);
1217 }
1218 } else {
1219 // edge: use 2nd order derivative
1220 int increment = (ix[n] == 0 ? 1 : -1);
1221 // move right from left edge, or the other way around
1222 A0 = value(ix) + offset;
1223 ix[n] += increment; A1 = value(ix) + offset;
1224 ix[n] += increment; A2 = value(ix) + offset;
1225 if (A0 * A1 * A2 == 0) {
1226 return 0.; // can't handle empty bins
1227 } else {
1228 return (-1.5 * cvm::logn(A0) + 2. * cvm::logn(A1)
1229 - 0.5 * cvm::logn(A2)) * increment / widths[n];
1230 }
1231 }
1232 }
1233
1234
1238 inline cvm::real gradient_finite_diff(const std::vector<int> &ix0,
1239 int n = 0)
1240 {
1241 cvm::real A0, A1, A2;
1242 std::vector<int> ix = ix0;
1243
1244 // FIXME this can be rewritten more concisely with wrap_edge()
1245 if (periodic[n]) {
1246 ix[n]--; wrap(ix);
1247 A0 = value(ix);
1248 ix = ix0;
1249 ix[n]++; wrap(ix);
1250 A1 = value(ix);
1251 if (A0 * A1 == 0) {
1252 return 0.; // can't handle empty bins
1253 } else {
1254 return (A1 - A0) / (widths[n] * 2.);
1255 }
1256 } else if (ix[n] > 0 && ix[n] < nx[n]-1) { // not an edge
1257 ix[n]--;
1258 A0 = value(ix);
1259 ix = ix0;
1260 ix[n]++;
1261 A1 = value(ix);
1262 if (A0 * A1 == 0) {
1263 return 0.; // can't handle empty bins
1264 } else {
1265 return (A1 - A0) / (widths[n] * 2.);
1266 }
1267 } else {
1268 // edge: use 2nd order derivative
1269 int increment = (ix[n] == 0 ? 1 : -1);
1270 // move right from left edge, or the other way around
1271 A0 = value(ix);
1272 ix[n] += increment; A1 = value(ix);
1273 ix[n] += increment; A2 = value(ix);
1274 return (-1.5 * A0 + 2. * A1
1275 - 0.5 * A2) * increment / widths[n];
1276 }
1277 }
1278};
1279
1280
1282class colvar_grid_scalar : public colvar_grid<cvm::real>
1283{
1284public:
1285
1289
1292
1295
1297 virtual ~colvar_grid_scalar();
1298
1300 colvar_grid_scalar(std::vector<colvar *> &colvars,
1301 std::shared_ptr<const colvar_grid_params> params = nullptr,
1302 bool add_extra_bin = false,
1303 std::string config = std::string());
1304
1306 colvar_grid_scalar(std::string const &filename);
1307
1309 inline void acc_value(std::vector<int> const &ix,
1310 cvm::real const &new_value,
1311 size_t const &imult = 0)
1312 {
1313 (void) imult;
1314 // only legal value of imult here is 0
1315 data[address(ix)] += new_value;
1316 if (samples)
1317 samples->incr_count(ix);
1318 has_data = true;
1319 }
1320
1322 std::string get_state_params() const;
1323
1325 int parse_params(std::string const &conf,
1327
1329 std::istream & read_restart(std::istream &is);
1330
1333
1335 std::ostream & write_restart(std::ostream &os);
1336
1339
1341 std::istream &read_raw(std::istream &is);
1342
1345
1349 std::ostream &write_raw(std::ostream &os, size_t const buf_size = 3) const;
1350
1354 cvm::memory_stream &write_raw(cvm::memory_stream &os, size_t const buf_size = 3) const;
1355
1357 std::istream & read_multicol(std::istream &is, bool add = false);
1358
1360 int read_multicol(std::string const &filename,
1361 std::string description = "grid file",
1362 bool add = false);
1363
1365 std::ostream & write_multicol(std::ostream &os) const;
1366
1368 int write_multicol(std::string const &filename,
1369 std::string description = "grid file") const;
1370
1372 std::ostream & write_opendx(std::ostream &os) const;
1373
1375 int write_opendx(std::string const &filename,
1376 std::string description = "grid file") const;
1377
1382 inline void vector_gradient_finite_diff( const std::vector<int> &ix0, std::vector<cvm::real> &grad)
1383 {
1384 cvm::real A0, A1;
1385 std::vector<int> ix;
1386 size_t i, j, k, n;
1387
1388 if (nd == 2) {
1389 for (n = 0; n < 2; n++) {
1390 ix = ix0;
1391 A0 = value(ix);
1392 ix[n]++; wrap(ix);
1393 A1 = value(ix);
1394 ix[1-n]++; wrap(ix);
1395 A1 += value(ix);
1396 ix[n]--; wrap(ix);
1397 A0 += value(ix);
1398 grad[n] = 0.5 * (A1 - A0) / widths[n];
1399 }
1400 } else if (nd == 3) {
1401
1402 cvm::real p[8]; // potential values within cube, indexed in binary (4 i + 2 j + k)
1403 ix = ix0;
1404 int index = 0;
1405 for (i = 0; i<2; i++) {
1406 ix[1] = ix0[1];
1407 for (j = 0; j<2; j++) {
1408 ix[2] = ix0[2];
1409 for (k = 0; k<2; k++) {
1410 wrap(ix);
1411 p[index++] = value(ix);
1412 ix[2]++;
1413 }
1414 ix[1]++;
1415 }
1416 ix[0]++;
1417 }
1418
1419 // The following would be easier to read using binary literals
1420 // 100 101 110 111 000 001 010 011
1421 grad[0] = 0.25 * ((p[4] + p[5] + p[6] + p[7]) - (p[0] + p[1] + p[2] + p[3])) / widths[0];
1422 // 010 011 110 111 000 001 100 101
1423 grad[1] = 0.25 * ((p[2] + p[3] + p[6] + p[7]) - (p[0] + p[1] + p[4] + p[5])) / widths[1];
1424 // 001 011 101 111 000 010 100 110
1425 grad[2] = 0.25 * ((p[1] + p[3] + p[5] + p[7]) - (p[0] + p[2] + p[4] + p[6])) / widths[2];
1426 } else {
1427 cvm::error("Finite differences available in dimension 2 and 3 only.");
1428 }
1429 }
1430
1431
1435 inline cvm::real log_gradient_finite_diff(const std::vector<int> &ix0,
1436 int n = 0, int offset = 0)
1437 {
1438 cvm::real A0, A1, A2;
1439 std::vector<int> ix = ix0;
1440
1441 // TODO this can be rewritten more concisely with wrap_edge()
1442 if (periodic[n]) {
1443 ix[n]--; wrap(ix);
1444 A0 = value(ix) + offset;
1445 ix = ix0;
1446 ix[n]++; wrap(ix);
1447 A1 = value(ix) + offset;
1448 if (A0 * A1 == 0) {
1449 return 0.; // can't handle empty bins
1450 } else {
1451 return (cvm::logn(A1) - cvm::logn(A0))
1452 / (widths[n] * 2.);
1453 }
1454 } else if (ix[n] > 0 && ix[n] < nx[n]-1) { // not an edge
1455 ix[n]--;
1456 A0 = value(ix) + offset;
1457 ix = ix0;
1458 ix[n]++;
1459 A1 = value(ix) + offset;
1460 if (A0 * A1 == 0) {
1461 return 0.; // can't handle empty bins
1462 } else {
1463 return (cvm::logn(A1) - cvm::logn(A0))
1464 / (widths[n] * 2.);
1465 }
1466 } else {
1467 // edge: use 2nd order derivative
1468 int increment = (ix[n] == 0 ? 1 : -1);
1469 // move right from left edge, or the other way around
1470 A0 = value(ix) + offset;
1471 ix[n] += increment; A1 = value(ix) + offset;
1472 ix[n] += increment; A2 = value(ix) + offset;
1473 if (A0 * A1 * A2 == 0) {
1474 return 0.; // can't handle empty bins
1475 } else {
1476 return (-1.5 * cvm::logn(A0) + 2. * cvm::logn(A1)
1477 - 0.5 * cvm::logn(A2)) * increment / widths[n];
1478 }
1479 }
1480 }
1481
1482
1486 inline cvm::real gradient_finite_diff(const std::vector<int> &ix0,
1487 int n = 0)
1488 {
1489 cvm::real A0, A1, A2;
1490 std::vector<int> ix = ix0;
1491
1492 // FIXME this can be rewritten more concisely with wrap_edge()
1493 if (periodic[n]) {
1494 ix[n]--; wrap(ix);
1495 A0 = value(ix);
1496 ix = ix0;
1497 ix[n]++; wrap(ix);
1498 A1 = value(ix);
1499 if (A0 * A1 == 0) {
1500 return 0.; // can't handle empty bins
1501 } else {
1502 return (A1 - A0) / (widths[n] * 2.);
1503 }
1504 } else if (ix[n] > 0 && ix[n] < nx[n]-1) { // not an edge
1505 ix[n]--;
1506 A0 = value(ix);
1507 ix = ix0;
1508 ix[n]++;
1509 A1 = value(ix);
1510 if (A0 * A1 == 0) {
1511 return 0.; // can't handle empty bins
1512 } else {
1513 return cvm::real(A1 - A0) / (widths[n] * 2.);
1514 }
1515 } else {
1516 // edge: use 2nd order derivative
1517 int increment = (ix[n] == 0 ? 1 : -1);
1518 // move right from left edge, or the other way around
1519 A0 = value(ix);
1520 ix[n] += increment; A1 = value(ix);
1521 ix[n] += increment; A2 = value(ix);
1522 return (-1.5 * cvm::real(A0) + 2. * cvm::real(A1)
1523 - 0.5 * cvm::real(A2)) * increment / widths[n];
1524 }
1525 }
1526
1527
1530 virtual inline cvm::real value_output(std::vector<int> const &ix,
1531 size_t const &imult = 0) const override
1532 {
1533 int s;
1534 if (imult > 0) {
1535 cvm::error("Error: trying to access a component "
1536 "larger than 1 in a scalar data grid.\n");
1537 return 0.;
1538 }
1539 if (samples) {
1540 return ( (s = samples->value(ix)) > 0) ?
1541 (data[address(ix) + imult] / cvm::real(s)) :
1542 0.0;
1543 } else {
1544 return data[address(ix) + imult];
1545 }
1546 }
1547
1549 virtual void value_input(std::vector<int> const &ix,
1550 cvm::real const &new_value,
1551 size_t const &imult = 0,
1552 bool add = false) override
1553 {
1554 if (imult > 0) {
1555 cvm::error("Error: trying to access a component "
1556 "larger than 1 in a scalar data grid.\n");
1557 return;
1558 }
1559 if (add) {
1560 if (samples)
1561 data[address(ix)] += new_value * samples->new_value(ix);
1562 else
1563 data[address(ix)] += new_value;
1564 } else {
1565 if (samples)
1566 data[address(ix)] = new_value * samples->value(ix);
1567 else
1568 data[address(ix)] = new_value;
1569 }
1570 has_data = true;
1571 }
1572
1574 cvm::real maximum_value() const;
1575
1577 cvm::real minimum_value() const;
1578
1581
1583 cvm::real integral() const;
1584
1587 cvm::real entropy() const;
1588
1591 cvm::real grid_rmsd(colvar_grid_scalar const &other_grid) const;
1592};
1593
1594
1595
1597class colvar_grid_gradient : public colvar_grid<cvm::real>
1598{
1599public:
1600
1603 std::shared_ptr<colvar_grid_count> samples;
1604
1607
1610 {}
1611
1612 // /// Constructor from specific sizes arrays
1613 // colvar_grid_gradient(std::vector<int> const &nx_i);
1614
1615 // /// Constructor from a vector of colvars
1616 // colvar_grid_gradient(std::vector<colvar *> &colvars,
1617 // std::string config = std::string());
1618
1620 colvar_grid_gradient(std::string const &filename);
1621
1623 colvar_grid_gradient(std::vector<colvar *> &colvars,
1624 std::shared_ptr<colvar_grid_count> samples_in = nullptr,
1625 std::shared_ptr<const colvar_grid_params> params = nullptr,
1626 std::string config = std::string());
1627
1630 int min_samples;
1631
1633 std::string get_state_params() const;
1634
1636 int parse_params(std::string const &conf,
1638
1640 std::istream & read_restart(std::istream &is);
1641
1644
1646 std::ostream & write_restart(std::ostream &os);
1647
1650
1652 std::istream &read_raw(std::istream &is);
1653
1656
1660 std::ostream &write_raw(std::ostream &os, size_t const buf_size = 3) const;
1661
1665 cvm::memory_stream &write_raw(cvm::memory_stream &os, size_t const buf_size = 3) const;
1666
1668 std::istream & read_multicol(std::istream &is, bool add = false);
1669
1671 int read_multicol(std::string const &filename,
1672 std::string description = "grid file",
1673 bool add = false);
1674
1676 std::ostream & write_multicol(std::ostream &os) const;
1677
1679 int write_multicol(std::string const &filename,
1680 std::string description = "grid file") const;
1681
1683 std::ostream & write_opendx(std::ostream &os) const;
1684
1686 int write_opendx(std::string const &filename,
1687 std::string description = "grid file") const;
1688
1690 inline void vector_value(std::vector<int> const &ix, std::vector<cvm::real> &v) const
1691 {
1692 cvm::real const * p = &value(ix);
1693 if (samples) {
1694 int count = samples->value(ix);
1695 if (count) {
1696 cvm::real invcount = 1.0 / count;
1697 for (size_t i = 0; i < mult; i++) {
1698 v[i] = invcount * p[i];
1699 }
1700 } else {
1701 for (size_t i = 0; i < mult; i++) {
1702 v[i] = 0.0;
1703 }
1704 }
1705 } else {
1706 for (size_t i = 0; i < mult; i++) {
1707 v[i] = p[i];
1708 }
1709 }
1710 }
1711
1712
1714 inline void acc_value(std::vector<int> const &ix, std::vector<colvarvalue> const &values) {
1715 for (size_t imult = 0; imult < mult; imult++) {
1716 data[address(ix) + imult] += values[imult].real_value;
1717 }
1718 if (samples)
1719 samples->incr_count(ix);
1720 }
1721
1724 inline void acc_force(std::vector<int> const &ix, cvm::real const *forces) {
1725 for (size_t imult = 0; imult < mult; imult++) {
1726 data[address(ix) + imult] -= forces[imult];
1727 }
1728 if (samples)
1729 samples->incr_count(ix);
1730 }
1731
1734 virtual cvm::real value_output(std::vector<int> const &ix,
1735 size_t const &imult = 0) const override
1736 {
1737 int s;
1738 if (samples) {
1739 return ( (s = samples->value(ix)) > 0) ?
1740 (data[address(ix) + imult] / cvm::real(s)) :
1741 0.0;
1742 } else {
1743 return data[address(ix) + imult];
1744 }
1745 }
1746
1750 {
1751 cvm::real fact;
1752 if ( weight <= min_samples ) {
1753 fact = 0.0;
1754 } else if ( weight < full_samples ) {
1755 fact = (weight - min_samples) / (weight * cvm::real(full_samples - min_samples));
1756 } else {
1757 fact = 1.0 / weight;
1758 }
1759 return fact;
1760 }
1761
1762
1767 virtual inline cvm::real value_output_smoothed(std::vector<int> const &ix, bool smoothed = true)
1768 {
1769 cvm::real weight, fact;
1770
1771 if (samples) {
1772 weight = cvm::real(samples->value(ix));
1773 } else {
1774 weight = 1.;
1775 }
1776
1777 if (smoothed) {
1778 fact = smooth_inverse_weight(weight);
1779 } else {
1780 fact = weight > 0. ? 1. / weight : 0.;
1781 }
1782
1783 return fact * data[address(ix)];
1784 }
1785
1789 inline void vector_value_smoothed(std::vector<int> const &ix, cvm::real *grad, bool smoothed = true)
1790 {
1791 cvm::real weight, fact;
1792
1793 if (samples) {
1794 weight = cvm::real(samples->value(ix));
1795 } else {
1796 weight = 1.;
1797 }
1798
1799 if (smoothed) {
1800 fact = smooth_inverse_weight(weight);
1801 } else {
1802 fact = weight > 0. ? 1. / weight : 0.;
1803 }
1804
1805 cvm::real *p = &(data[address(ix)]);
1806
1807 // Appease Clang analyzer, which likes to assume that mult is zero
1808 #ifdef __clang_analyzer__
1809 assert(mult > 0);
1810 #endif
1811
1812 for (size_t imult = 0; imult < mult; imult++) {
1813 grad[imult] = fact * p[imult];
1814 }
1815 }
1816
1820 virtual void value_input(std::vector<int> const &ix,
1821 cvm::real const &new_value,
1822 size_t const &imult = 0,
1823 bool add = false) override
1824 {
1825 if (add) {
1826 if (samples)
1827 data[address(ix) + imult] += new_value * samples->new_value(ix);
1828 else
1829 data[address(ix) + imult] += new_value;
1830 } else {
1831 if (samples)
1832 data[address(ix) + imult] = new_value * samples->value(ix);
1833 else
1834 data[address(ix) + imult] = new_value;
1835 }
1836 has_data = true;
1837 }
1838
1839
1841 inline cvm::real average(bool smoothed = false)
1842 {
1843 if (nd != 1 || nx[0] == 0) {
1844 return 0.0;
1845 }
1846
1847 cvm::real sum = 0.0;
1848 for (std::vector<int> ix = new_index(); index_ok(ix); incr(ix)) {
1849 sum += value_output_smoothed(ix, smoothed);
1850 }
1851
1852 return (sum / cvm::real(nx[0]));
1853 }
1854
1857 cvm::real grid_rmsd(colvar_grid_gradient const &other_grid) const;
1858
1861 void write_1D_integral(std::ostream &os);
1862
1863};
1864
1865
1866
1868
1870{
1871 public:
1872
1874
1875 virtual ~integrate_potential()
1876 {}
1877
1879 integrate_potential(std::vector<colvar *> &colvars,
1880 std::shared_ptr<colvar_grid_gradient> gradients);
1881
1883 integrate_potential(std::shared_ptr<colvar_grid_gradient> gradients);
1884
1886 int integrate(const int itmax, const cvm::real & tol, cvm::real & err, bool verbose = true);
1887
1890 void update_div_neighbors(const std::vector<int> &ix);
1891
1894 void update_div_local(const std::vector<int> &ix);
1895
1898 void set_div();
1899
1902 inline void set_zero_minimum() {
1903 add_constant(-1.0 * minimum_value());
1904 }
1905
1908
1909
1910 protected:
1911
1912 // Reference to gradient grid
1913 std::shared_ptr<colvar_grid_gradient> gradients;
1914
1916 std::vector<cvm::real> divergence;
1917
1918// std::vector<cvm::real> inv_lap_diag; // Inverse of the diagonal of the Laplacian; for conditioning
1919
1923 void get_grad(cvm::real * g, std::vector<int> &ix);
1924
1926 void nr_linbcg_sym(const std::vector<cvm::real> &b, std::vector<cvm::real> &x,
1927 const cvm::real &tol, const int itmax, int &iter, cvm::real &err);
1928
1930 cvm::real l2norm(const std::vector<cvm::real> &x);
1931
1933 void atimes(const std::vector<cvm::real> &x, std::vector<cvm::real> &r);
1934
1935// /// Inversion of preconditioner matrix
1936// void asolve(const std::vector<cvm::real> &b, std::vector<cvm::real> &x);
1937};
1938
1939#endif
1940
Colvar_grid derived class to hold counters in discrete n-dim colvar space.
Definition: colvargrid.h:1002
colvar_grid_count()
Default constructor.
Definition: colvargrid.cpp:21
int local_sample_count(int radius)
Return the average number of samples in a given "radius" around current bin Really a hypercube of len...
Definition: colvargrid.h:1109
std::istream & read_raw(std::istream &is)
Read all grid parameters and data from a formatted stream.
Definition: colvargrid.cpp:68
std::string get_state_params() const
Write the current grid parameters to a string.
Definition: colvargrid.cpp:37
cvm::real gradient_finite_diff(const std::vector< int > &ix0, int n=0)
Return the gradient of discrete count from finite differences on the same grid for dimension n (colva...
Definition: colvargrid.h:1238
std::ostream & write_restart(std::ostream &os)
Write all grid parameters and data to a formatted stream.
Definition: colvargrid.cpp:58
std::ostream & write_opendx(std::ostream &os) const
Write the grid data without labels, as they are represented in memory.
Definition: colvargrid.cpp:112
std::ostream & write_raw(std::ostream &os, size_t const buf_size=3) const
Definition: colvargrid.cpp:78
int parse_params(std::string const &conf, colvarparse::Parse_Mode const parse_mode=colvarparse::parse_normal)
Read new grid parameters from a string.
Definition: colvargrid.cpp:42
virtual void value_input(std::vector< int > const &ix, size_t const &t, size_t const &imult=0, bool add=false)
Enter or add a value, but also handle parent grid.
Definition: colvargrid.h:1089
std::istream & read_restart(std::istream &is)
Read all grid parameters and data from a formatted stream.
Definition: colvargrid.cpp:48
std::istream & read_multicol(std::istream &is, bool add=false)
Read a grid written by write_multicol(), incrementin if data is true.
Definition: colvargrid.cpp:89
std::ostream & write_multicol(std::ostream &os) const
Write grid in a format which is both human-readable and gnuplot-friendly.
Definition: colvargrid.cpp:101
cvm::real log_gradient_finite_diff(const std::vector< int > &ix0, int n=0, int offset=0)
Return the log-gradient from finite differences on the same grid for dimension n (colvar_grid_count)
Definition: colvargrid.h:1187
size_t const & new_value(std::vector< int > const &ix)
Get the binned count indexed by ix from the newly read data.
Definition: colvargrid.h:1026
virtual ~colvar_grid_count()
Destructor.
Definition: colvargrid.h:1009
void incr_count(std::vector< int > const &ix)
Increment the counter at given position.
Definition: colvargrid.h:1020
Class for accumulating the gradient of a scalar function on a grid.
Definition: colvargrid.h:1598
std::string get_state_params() const
Write the current grid parameters to a string.
Definition: colvargrid.cpp:369
cvm::real smooth_inverse_weight(cvm::real weight)
Definition: colvargrid.h:1749
void vector_value_smoothed(std::vector< int > const &ix, cvm::real *grad, bool smoothed=true)
Obtain the vector value of the function at ix divided by its number of samples (if the count grid is ...
Definition: colvargrid.h:1789
virtual cvm::real value_output_smoothed(std::vector< int > const &ix, bool smoothed=true)
Return the scalar value of the function at ix divided by its number of samples (if the count grid is ...
Definition: colvargrid.h:1767
void acc_value(std::vector< int > const &ix, std::vector< colvarvalue > const &values)
Accumulate the value.
Definition: colvargrid.h:1714
cvm::real average(bool smoothed=false)
Compute and return average value for a 1D gradient grid.
Definition: colvargrid.h:1841
std::shared_ptr< colvar_grid_count > samples
Provide the sample count by which each binned value should be divided.
Definition: colvargrid.h:1603
std::istream & read_multicol(std::istream &is, bool add=false)
Read a grid written by write_multicol(), incrementin if data is true.
Definition: colvargrid.cpp:421
colvar_grid_gradient()
Default constructor.
Definition: colvargrid.cpp:335
std::ostream & write_restart(std::ostream &os)
Write all grid parameters and data to a formatted stream.
Definition: colvargrid.cpp:390
void write_1D_integral(std::ostream &os)
If the grid is 1-dimensional, integrate it and write the integral to a file (DEPRECATED by the integr...
Definition: colvargrid.cpp:456
int parse_params(std::string const &conf, colvarparse::Parse_Mode const parse_mode=colvarparse::parse_normal)
Read new grid parameters from a string.
Definition: colvargrid.cpp:374
virtual cvm::real value_output(std::vector< int > const &ix, size_t const &imult=0) const override
Return the value of the function at ix divided by its number of samples (if the count grid is defined...
Definition: colvargrid.h:1734
std::istream & read_raw(std::istream &is)
Read all grid parameters and data from a formatted stream.
Definition: colvargrid.cpp:400
virtual ~colvar_grid_gradient()
Destructor.
Definition: colvargrid.h:1609
virtual void value_input(std::vector< int > const &ix, cvm::real const &new_value, size_t const &imult=0, bool add=false) override
Get the value from a formatted output and transform it into the internal representation (it may have ...
Definition: colvargrid.h:1820
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:1724
int full_samples
Parameters for smoothing data with low sampling.
Definition: colvargrid.h:1629
std::ostream & write_raw(std::ostream &os, size_t const buf_size=3) const
Definition: colvargrid.cpp:410
std::ostream & write_multicol(std::ostream &os) const
Write grid in a format which is both human-readable and gnuplot-friendly.
Definition: colvargrid.cpp:433
std::istream & read_restart(std::istream &is)
Read all grid parameters and data from a formatted stream.
Definition: colvargrid.cpp:380
cvm::real grid_rmsd(colvar_grid_gradient const &other_grid) const
Return the RMSD between this grid's data and another one Grids must have the same dimensions.
Definition: colvargrid.cpp:513
void vector_value(std::vector< int > const &ix, std::vector< cvm::real > &v) const
Get a vector with the binned value(s) indexed by ix, normalized if applicable.
Definition: colvargrid.h:1690
std::ostream & write_opendx(std::ostream &os) const
Write the grid data without labels, as they are represented in memory.
Definition: colvargrid.cpp:444
Unified base class for grid of values of a function of several collective variables.
Definition: colvargrid.h:24
size_t nd
Number of dimensions.
Definition: colvargrid.h:28
std::vector< int > nxc
Cumulative number of points along each dimension.
Definition: colvargrid.h:34
std::vector< colvarvalue > lower_boundaries
Lower boundaries of the colvars in this grid.
Definition: colvargrid.h:37
std::vector< cvm::real > widths
Widths of the colvars in this grid.
Definition: colvargrid.h:43
std::vector< colvarvalue > upper_boundaries
Upper boundaries of the colvars in this grid.
Definition: colvargrid.h:40
std::vector< int > nx
Number of points along each dimension.
Definition: colvargrid.h:31
Class for accumulating a scalar function on a grid.
Definition: colvargrid.h:1283
void acc_value(std::vector< int > const &ix, cvm::real const &new_value, size_t const &imult=0)
Accumulate the value.
Definition: colvargrid.h:1309
cvm::real minimum_value() const
Return the lowest value.
Definition: colvargrid.cpp:249
cvm::real gradient_finite_diff(const std::vector< int > &ix0, int n=0)
Return the gradient of discrete count from finite differences on the same grid for dimension n (colva...
Definition: colvargrid.h:1486
std::istream & read_restart(std::istream &is)
Read all grid parameters and data from a formatted stream.
Definition: colvargrid.cpp:163
std::string get_state_params() const
Write the current grid parameters to a string.
Definition: colvargrid.cpp:152
virtual ~colvar_grid_scalar()
Destructor.
Definition: colvargrid.cpp:148
cvm::real minimum_pos_value() const
Return the lowest positive value.
Definition: colvargrid.cpp:258
virtual cvm::real value_output(std::vector< int > const &ix, size_t const &imult=0) const override
Return the value of the function at ix divided by its number of samples (if the count grid is defined...
Definition: colvargrid.h:1530
cvm::real grid_rmsd(colvar_grid_scalar const &other_grid) const
Return the RMSD between this grid's data and another one Grids must have the same dimensions.
Definition: colvargrid.cpp:305
std::ostream & write_opendx(std::ostream &os) const
Write the grid data without labels, as they are represented in memory.
Definition: colvargrid.cpp:227
std::ostream & write_restart(std::ostream &os)
Write all grid parameters and data to a formatted stream.
Definition: colvargrid.cpp:173
std::ostream & write_multicol(std::ostream &os) const
Write grid in a format which is both human-readable and gnuplot-friendly.
Definition: colvargrid.cpp:216
std::istream & read_multicol(std::istream &is, bool add=false)
Read a grid written by write_multicol(), incrementin if data is true.
Definition: colvargrid.cpp:204
cvm::real integral() const
Calculates the integral of the map (uses widths if they are defined)
Definition: colvargrid.cpp:274
int parse_params(std::string const &conf, colvarparse::Parse_Mode const parse_mode=colvarparse::parse_normal)
Read new grid parameters from a string.
Definition: colvargrid.cpp:157
std::ostream & write_raw(std::ostream &os, size_t const buf_size=3) const
Definition: colvargrid.cpp:193
virtual void value_input(std::vector< int > const &ix, cvm::real const &new_value, size_t const &imult=0, bool add=false) override
Enter or add value but also deal with count grid.
Definition: colvargrid.h:1549
cvm::real entropy() const
Assuming that the map is a normalized probability density, calculates the entropy (uses widths if the...
Definition: colvargrid.cpp:288
colvar_grid_scalar()
Default constructor.
Definition: colvargrid.cpp:125
cvm::real maximum_value() const
Return the highest value.
Definition: colvargrid.cpp:239
cvm::real log_gradient_finite_diff(const std::vector< int > &ix0, int n=0, int offset=0)
Return the log-gradient from finite differences on the same grid for dimension n (colvar_grid_scalar)
Definition: colvargrid.h:1435
colvar_grid_count * samples
Provide the associated sample count by which each binned value should be divided.
Definition: colvargrid.h:1288
void vector_gradient_finite_diff(const std::vector< int > &ix0, std::vector< cvm::real > &grad)
Return the gradient of the scalar field from finite differences Input coordinates are those of gradie...
Definition: colvargrid.h:1382
std::istream & read_raw(std::istream &is)
Read all grid parameters and data from a formatted stream.
Definition: colvargrid.cpp:183
Grid of values of a function of several collective variables.
Definition: colvargrid.h:52
int setup()
Allocate data (allow initialization also after construction)
Definition: colvargrid.h:200
std::istream & read_restart(std::istream &is)
Read all grid parameters and data from a formatted stream.
Definition: colvargrid_def.h:100
void set_value(size_t i, T const &t)
Set the value at the point with linear address i (for speed)
Definition: colvargrid.h:569
void check_consistency()
Check that the grid information inside (boundaries, widths, ...) is consistent with the current setti...
Definition: colvargrid.h:906
std::string get_state_params() const
Write the current grid parameters to a string.
Definition: colvargrid_def.h:202
std::ostream & write_raw(std::ostream &os, size_t const buf_size=3) const
Definition: colvargrid_def.h:169
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:489
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:854
void add_constant(T const &t)
Add a constant to all elements (fast loop)
Definition: colvargrid.h:671
std::vector< size_t > new_data
Newly read data (used for count grids, when adding several grids read from disk)
Definition: colvargrid.h:68
cvm::real value_to_bin_scalar_fraction(colvarvalue const &value, const int i) const
Use the lower boundary and the width to report the fraction of bin beyond value_to_bin_scalar() that ...
Definition: colvargrid.h:517
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:872
size_t number_of_points(int const icv=-1) const
Definition: colvargrid.h:122
std::vector< bool > hard_upper_boundaries
Whether some colvars have hard upper boundaries.
Definition: colvargrid.h:100
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:833
T const & value(size_t i) const
Get the binned value indexed by linear address i.
Definition: colvargrid.h:665
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:717
int current_bin_flat_bound() const
Report the flattened bin address corresponding to the current value of all variables and assign first...
Definition: colvargrid.h:478
size_t mult
Multiplicity of each datum (allow the binning of non-scalar types such as atomic gradients)
Definition: colvargrid.h:59
size_t multiplicity() const
Return the multiplicity of the type used.
Definition: colvargrid.h:144
std::vector< int > const & sizes() const
Get the sizes in each direction.
Definition: colvargrid.h:132
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:560
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:546
std::vector< T > data
Low-level array of values.
Definition: colvargrid.h:65
void multiply_constant(cvm::real const &a)
Multiply all elements by a scalar constant (fast loop)
Definition: colvargrid.h:679
std::vector< int > const & number_of_points_vec() const
Return the numbers of points in all dimensions.
Definition: colvargrid.h:115
bool index_ok(std::vector< int > const &ix) const
Check that the index is within range in each of the dimensions.
Definition: colvargrid.h:861
cvm::real current_bin_scalar_fraction(int const i) const
Report the fraction of bin beyond current_bin_scalar()
Definition: colvargrid.h:510
bool has_data
Whether this grid has been filled with data or is still empty.
Definition: colvargrid.h:106
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,...
Definition: colvargrid.h:802
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:695
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:925
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:245
int parse_params(std::string const &conf, colvarparse::Parse_Mode const parse_mode=colvarparse::parse_normal)
Read new grid parameters from a string.
Definition: colvargrid_def.h:232
size_t nt
Total number of grid points.
Definition: colvargrid.h:62
T get_value(size_t i) const
Get the value at the point with linear address i (for speed)
Definition: colvargrid.h:575
void reset(T const &t=T())
Reset data (in case the grid is being reused)
Definition: colvargrid.h:206
size_t num_variables() const
Return the number of colvar objects.
Definition: colvargrid.h:109
int setup(std::vector< int > const &nx_i, T const &t=T(), size_t const &mult_i=1)
Allocate data.
Definition: colvargrid.h:159
colvar_grid(std::vector< colvar * > const &colvars, T const &t=T(), size_t mult_i=1, bool add_extra_bin=false, std::shared_ptr< const colvar_grid_params > params=nullptr, std::string config=std::string())
Constructor from a vector of colvars or an optional grid config string.
Definition: colvargrid.h:256
colvar_grid(colvar_grid< T > const &g)
"Almost copy-constructor": only copies configuration parameters from another grid,...
Definition: colvargrid.h:228
bool has_parent_data
True if this is a count grid related to another grid of data.
Definition: colvargrid.h:103
std::ostream & write_opendx(std::ostream &os) const
Write the grid data without labels, as they are represented in memory.
Definition: colvargrid_def.h:493
bool wrap_to_edge(std::vector< int > &ix, std::vector< int > &edge_bin) const
Definition: colvargrid.h:450
std::vector< int > const get_colvars_index() const
Get the bin indices corresponding to the current values of the colvars.
Definition: colvargrid.h:706
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:495
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:585
void request_actual_value(bool b=true)
Request grid to use actual values of extended coords.
Definition: colvargrid.h:150
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:729
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:525
int current_bin_scalar(int const i) const
Report the bin corresponding to the current value of variable i.
Definition: colvargrid.h:471
std::istream & read_raw(std::istream &is)
Read all grid parameters and data from a formatted stream.
Definition: colvargrid_def.h:156
std::vector< bool > hard_lower_boundaries
Whether some colvars have hard lower boundaries.
Definition: colvargrid.h:97
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:658
std::vector< bool > use_actual_value
Do we request actual value (for extended-system colvars)?
Definition: colvargrid.h:74
size_t raw_data_num() const
Size of the data as they are represented in memory.
Definition: colvargrid.h:653
size_t address(std::vector< int > const &ix) const
Get the low-level index corresponding to an index.
Definition: colvargrid.h:77
void remove_small_values(cvm::real const &a)
Assign values that are smaller than scalar constant the latter value (fast loop)
Definition: colvargrid.h:686
std::ostream & write_multicol(std::ostream &os) const
Write grid in a format which is both human-readable and gnuplot-friendly.
Definition: colvargrid_def.h:429
virtual ~colvar_grid()
Destructor.
Definition: colvargrid.h:222
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:504
virtual T value_output(std::vector< int > const &ix, size_t const &imult=0) const
Return the value suitable for output purposes (so that it may be rescaled or manipulated without chan...
Definition: colvargrid.h:824
void set_sizes(std::vector< int > const &new_sizes)
Set the sizes in each direction.
Definition: colvargrid.h:138
void copy_grid(colvar_grid< T > const &other_grid)
Copy data from another grid of the same type, AND identical definition (boundaries,...
Definition: colvargrid.h:610
bool wrap_detect_edge(std::vector< int > &ix) const
Definition: colvargrid.h:435
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:633
void map_grid(colvar_grid< T > const &other_grid)
Add data from another grid of the same type.
Definition: colvargrid.h:759
std::istream & read_multicol(std::istream &is, bool add=false)
Read a grid written by write_multicol(), incrementing if add is true.
Definition: colvargrid_def.h:310
std::vector< bool > periodic
Whether some colvars are periodic.
Definition: colvargrid.h:94
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:537
void wrap(std::vector< int > &ix) const
Definition: colvargrid.h:418
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:552
colvar_grid()
Default constructor.
Definition: colvargrid.h:213
std::vector< colvar * > cv
Colvars collected in this grid.
Definition: colvargrid.h:71
std::ostream & write_restart(std::ostream &os)
Write all grid parameters and data to a formatted stream.
Definition: colvargrid_def.h:112
void raw_data_in(const T *in_data)
Input the data as they are represented in memory.
Definition: colvargrid.h:642
@ f_cv_hard_lower_boundary
The lower boundary is not defined from user's choice.
Definition: colvardeps.h:326
@ f_cv_hard_upper_boundary
The upper boundary is not defined from user's choice.
Definition: colvardeps.h:328
static real logn(real const &x)
Definition: colvarmodule.h:209
double real
Defining an abstract real number allows to switch precision.
Definition: colvarmodule.h:100
static real floor(real const &x)
Reimplemented to work around MS compiler issues.
Definition: colvarmodule.h:126
static void log(std::string const &message, int min_log_level=10)
Definition: colvarmodule.cpp:1993
static int error(std::string const &message, int code=-1)
Print a message to the main log and set global error code.
Definition: colvarmodule.cpp:2070
static size_t const cv_prec
Number of digits to represent a collective variables value(s)
Definition: colvarmodule.h:699
static real sqrt(real const &x)
Reimplemented to work around MS compiler issues.
Definition: colvarmodule.h:138
static size_t const cv_width
Number of characters to represent a collective variables value(s)
Definition: colvarmodule.h:701
static bool debug()
Whether debug output should be enabled (compile-time option)
Definition: colvarmodule.h:365
static real fabs(real const &x)
Reimplemented to work around MS compiler issues.
Definition: colvarmodule.h:132
static std::string to_str(char const *s)
Convert to string for output purposes.
Definition: colvarmodule.cpp:2416
Base class containing parsing functions; all objects which need to parse input inherit from this.
Definition: colvarparse.h:27
int check_keywords(std::string &conf, char const *key)
Check that all the keywords within "conf" are in the list of allowed keywords; this will invoke strip...
Definition: colvarparse.cpp:591
Parse_Mode
How a keyword is parsed in a string.
Definition: colvarparse.h:53
@ parse_normal
Alias for old default behavior (should be phased out)
Definition: colvarparse.h:72
Value of a collective variable: this is a metatype which can be set at runtime. By default it is set ...
Definition: colvarvalue.h:43
@ type_scalar
Scalar number, implemented as colvarmodule::real (default)
Definition: colvarvalue.h:56
cvm::real real_value
Real data member.
Definition: colvarvalue.h:77
Definition: colvars_memstream.h:30
Integrate (1D, 2D or 3D) gradients.
Definition: colvargrid.h:1870
void set_zero_minimum()
Add constant to potential so that its minimum value is zero Useful e.g. for output.
Definition: colvargrid.h:1902
void update_div_local(const std::vector< int > &ix)
Update matrix containing divergence and boundary conditions called by update_div_neighbors and by col...
Definition: colvargrid.cpp:699
std::vector< cvm::real > divergence
Array holding divergence + boundary terms (modified Neumann) if not periodic.
Definition: colvargrid.h:1916
void set_div()
Set matrix containing divergence and boundary conditions based on complete gradient grid.
Definition: colvargrid.cpp:636
void update_div_neighbors(const std::vector< int > &ix)
Update matrix containing divergence and boundary conditions based on new gradient point value,...
Definition: colvargrid.cpp:645
cvm::real l2norm(const std::vector< cvm::real > &x)
l2 norm of a vector
Definition: colvargrid.cpp:1220
bool b_smoothed
Flag requesting the use of a smoothed version of the gradient (default: false)
Definition: colvargrid.h:1907
int integrate(const int itmax, const cvm::real &tol, cvm::real &err, bool verbose=true)
Calculate potential from divergence (in 2D); return number of steps.
Definition: colvargrid.cpp:596
void get_grad(cvm::real *g, std::vector< int > &ix)
Definition: colvargrid.cpp:683
void nr_linbcg_sym(const std::vector< cvm::real > &b, std::vector< cvm::real > &x, const cvm::real &tol, const int itmax, int &iter, cvm::real &err)
Solve linear system based on CG, valid for symmetric matrices only.
Definition: colvargrid.cpp:1167
void atimes(const std::vector< cvm::real > &x, std::vector< cvm::real > &r)
Multiplication by sparse matrix representing Lagrangian (or its transpose)
Definition: colvargrid.cpp:751
Collective variables main module.
Parsing functions for collective variables.