Collective Variables Module - Developer Documentation
Searching...
No Matches
colvartypes.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 COLVARTYPES_H
11#define COLVARTYPES_H
12
13#include <sstream> // TODO specialize templates and replace this with iosfwd
14#include <vector>
15
16#ifdef COLVARS_LAMMPS
17// Use open-source Jacobi implementation
18#include "math_eigen_impl.h"
19#endif
20
21#include "colvarmodule.h"
22
23#ifndef PI
24#define PI 3.14159265358979323846
25#endif
26
27// ----------------------------------------------------------------------
30// ----------------------------------------------------------------------
31
32
36template <class T> class colvarmodule::vector1d
37{
38protected:
39
40 std::vector<T> data;
41
42public:
43
45 inline vector1d(size_t const n = 0)
46 {
47 data.resize(n);
48 reset();
49 }
50
52 inline vector1d(size_t const n, T const *t)
53 {
54 data.resize(n);
55 reset();
56 size_t i;
57 for (i = 0; i < size(); i++) {
58 data[i] = t[i];
59 }
60 }
61
63 inline vector1d(const vector1d&) = default;
64
66 inline vector1d& operator=(const vector1d&) = default;
67
69 inline T * c_array()
70 {
71 if (data.size() > 0) {
72 return &(data[0]);
73 } else {
74 return NULL;
75 }
76 }
77
79 inline std::vector<T> &data_array()
80 {
81 return data;
82 }
83
85 inline std::vector<T> const &data_array() const
86 {
87 return data;
88 }
89
90 inline ~vector1d()
91 {
92 data.clear();
93 }
94
96 inline void reset()
97 {
98 data.assign(data.size(), T(0.0));
99 }
100
101 inline size_t size() const
102 {
103 return data.size();
104 }
105
106 inline void resize(size_t const n)
107 {
108 data.resize(n);
109 }
110
111 inline void clear()
112 {
113 data.clear();
114 }
115
116 inline T & operator [] (size_t const i) {
117 return data[i];
118 }
119
120 inline T const & operator [] (size_t const i) const {
121 return data[i];
122 }
123
124 inline static void check_sizes(vector1d<T> const &v1, vector1d<T> const &v2)
125 {
126 if (v1.size() != v2.size()) {
127 cvm::error("Error: trying to perform an operation between vectors of different sizes, "+
128 cvm::to_str(v1.size())+" and "+cvm::to_str(v2.size())+".\n");
129 }
130 }
131
132 inline void operator += (vector1d<T> const &v)
133 {
134 check_sizes(*this, v);
135 size_t i;
136 for (i = 0; i < this->size(); i++) {
137 (*this)[i] += v[i];
138 }
139 }
140
141 inline void operator -= (vector1d<T> const &v)
142 {
143 check_sizes(*this, v);
144 size_t i;
145 for (i = 0; i < this->size(); i++) {
146 (*this)[i] -= v[i];
147 }
148 }
149
150 inline void operator *= (cvm::real a)
151 {
152 size_t i;
153 for (i = 0; i < this->size(); i++) {
154 (*this)[i] *= a;
155 }
156 }
157
158 inline void operator /= (cvm::real a)
159 {
160 size_t i;
161 for (i = 0; i < this->size(); i++) {
162 (*this)[i] /= a;
163 }
164 }
165
166 inline friend vector1d<T> operator + (vector1d<T> const &v1,
167 vector1d<T> const &v2)
168 {
169 check_sizes(v1.size(), v2.size());
170 vector1d<T> result(v1.size());
171 size_t i;
172 for (i = 0; i < v1.size(); i++) {
173 result[i] = v1[i] + v2[i];
174 }
175 return result;
176 }
177
178 inline friend vector1d<T> operator - (vector1d<T> const &v1,
179 vector1d<T> const &v2)
180 {
181 check_sizes(v1.size(), v2.size());
182 vector1d<T> result(v1.size());
183 size_t i;
184 for (i = 0; i < v1.size(); i++) {
185 result[i] = v1[i] - v2[i];
186 }
187 return result;
188 }
189
190 inline friend vector1d<T> operator * (vector1d<T> const &v, cvm::real a)
191 {
192 vector1d<T> result(v.size());
193 size_t i;
194 for (i = 0; i < v.size(); i++) {
195 result[i] = v[i] * a;
196 }
197 return result;
198 }
199
200 inline friend vector1d<T> operator * (cvm::real a, vector1d<T> const &v)
201 {
202 return v * a;
203 }
204
205 inline friend vector1d<T> operator / (vector1d<T> const &v, cvm::real a)
206 {
207 vector1d<T> result(v.size());
208 size_t i;
209 for (i = 0; i < v.size(); i++) {
210 result[i] = v[i] / a;
211 }
212 return result;
213 }
214
216 inline friend T operator * (vector1d<T> const &v1, vector1d<T> const &v2)
217 {
218 check_sizes(v1.size(), v2.size());
219 T prod(0.0);
220 size_t i;
221 for (i = 0; i < v1.size(); i++) {
222 prod += v1[i] * v2[i];
223 }
224 return prod;
225 }
226
228 inline cvm::real norm2() const
229 {
230 cvm::real result = 0.0;
231 size_t i;
232 for (i = 0; i < this->size(); i++) {
233 result += (*this)[i] * (*this)[i];
234 }
235 return result;
236 }
237
238 inline cvm::real norm() const
239 {
240 return cvm::sqrt(this->norm2());
241 }
242
243 inline cvm::real sum() const
244 {
245 cvm::real result = 0.0;
246 size_t i;
247 for (i = 0; i < this->size(); i++) {
248 result += (*this)[i];
249 }
250 return result;
251 }
252
254 inline vector1d<T> const slice(size_t const i1, size_t const i2) const
255 {
256 if ((i2 < i1) || (i2 >= this->size())) {
257 cvm::error("Error: trying to slice a vector using incorrect boundaries.\n");
258 }
259 vector1d<T> result(i2 - i1);
260 size_t i;
261 for (i = 0; i < (i2 - i1); i++) {
262 result[i] = (*this)[i1+i];
263 }
264 return result;
265 }
266
268 inline void sliceassign(size_t const i1, size_t const i2,
269 vector1d<T> const &v)
270 {
271 if ((i2 < i1) || (i2 >= this->size())) {
272 cvm::error("Error: trying to slice a vector using incorrect boundaries.\n");
273 }
274 size_t i;
275 for (i = 0; i < (i2 - i1); i++) {
276 (*this)[i1+i] = v[i];
277 }
278 }
279
281
282 inline size_t output_width(size_t real_width) const
283 {
284 return real_width*(this->size()) + 3*(this->size()-1) + 4;
285 }
286
287 inline friend std::istream & operator >> (std::istream &is,
289 {
290 if (v.size() == 0) return is;
291 std::streampos const start_pos = is.tellg();
292 char sep;
293 if ( !(is >> sep) || !(sep == '(') ) {
294 is.clear();
295 is.seekg(start_pos, std::ios::beg);
296 is.setstate(std::ios::failbit);
297 return is;
298 }
299 size_t count = 0;
300 while ( (is >> v[count]) &&
301 (count < (v.size()-1) ? ((is >> sep) && (sep == ',')) : true) ) {
302 if (++count == v.size()) break;
303 }
304 if (count < v.size()) {
305 is.clear();
306 is.seekg(start_pos, std::ios::beg);
307 is.setstate(std::ios::failbit);
308 }
309 return is;
310 }
311
312 inline friend std::ostream & operator << (std::ostream &os,
313 cvm::vector1d<T> const &v)
314 {
315 std::streamsize const w = os.width();
316 std::streamsize const p = os.precision();
317
318 os.width(2);
319 os << "( ";
320 size_t i;
321 for (i = 0; i < v.size()-1; i++) {
322 os.width(w); os.precision(p);
323 os << v[i] << " , ";
324 }
325 os.width(w); os.precision(p);
326 os << v[v.size()-1] << " )";
327 return os;
328 }
329
330 inline std::string to_simple_string() const
331 {
332 if (this->size() == 0) return std::string("");
333 std::ostringstream os;
334 os.setf(std::ios::scientific, std::ios::floatfield);
335 os.precision(cvm::cv_prec);
336 os << (*this)[0];
337 size_t i;
338 for (i = 1; i < this->size(); i++) {
339 os << " " << (*this)[i];
340 }
341 return os.str();
342 }
343
344 inline int from_simple_string(std::string const &s)
345 {
346 std::stringstream stream(s);
347 size_t i = 0;
348 if (this->size()) {
349 while ((stream >> (*this)[i]) && (i < this->size())) {
350 i++;
351 }
352 if (i < this->size()) {
353 return COLVARS_ERROR;
354 }
355 } else {
356 T input;
357 while (stream >> input) {
358 if ((i % 100) == 0) {
359 data.reserve(data.size()+100);
360 }
361 data.resize(data.size()+1);
362 data[i] = input;
363 i++;
364 }
365 }
366 return COLVARS_OK;
367 }
368
369};
370
371
375template <class T> class colvarmodule::matrix2d
376{
377public:
378
379 friend class row;
380 size_t outer_length;
381 size_t inner_length;
382
383protected:
384
385 class row {
386 public:
387 T * data;
388 size_t length;
389 inline row(T * const row_data, size_t const inner_length)
390 : data(row_data), length(inner_length)
391 {}
392 inline T & operator [] (size_t const j) {
393 return *(data+j);
394 }
395 inline T const & operator [] (size_t const j) const {
396 return *(data+j);
397 }
398 inline operator vector1d<T>() const
399 {
400 return vector1d<T>(length, data);
401 }
402 inline int set(cvm::vector1d<T> const &v) const
403 {
404 if (v.size() != length) {
405 return cvm::error("Error: setting a matrix row from a vector of "
406 "incompatible size.\n", COLVARS_BUG_ERROR);
407 }
408 for (size_t i = 0; i < length; i++) data[i] = v[i];
409 return COLVARS_OK;
410 }
411 };
412
413 std::vector<T> data;
414 std::vector<row> rows;
415 std::vector<T *> pointers;
416
417public:
418
420 inline void resize(size_t const ol, size_t const il)
421 {
422 if ((ol > 0) && (il > 0)) {
423
424 if (data.size() > 0) {
425 // copy previous data
426 size_t i, j;
427 std::vector<T> new_data(ol * il);
428 for (i = 0; i < outer_length; i++) {
429 for (j = 0; j < inner_length; j++) {
430 new_data[il*i+j] = data[inner_length*i+j];
431 }
432 }
433 data.resize(ol * il);
434 // copy them back
435 data = new_data;
436 } else {
437 data.resize(ol * il);
438 }
439
440 outer_length = ol;
441 inner_length = il;
442
443 if (data.size() > 0) {
444 // rebuild rows
445 size_t i;
446 rows.clear();
447 rows.reserve(outer_length);
448 pointers.clear();
449 pointers.reserve(outer_length);
450 for (i = 0; i < outer_length; i++) {
451 rows.push_back(row(&(data[0])+inner_length*i, inner_length));
452 pointers.push_back(&(data[0])+inner_length*i);
453 }
454 }
455 } else {
456 // zero size
457 data.clear();
458 rows.clear();
459 }
460 }
461
463 inline void clear() {
464 rows.clear();
465 data.clear();
466 }
467
469 inline void reset()
470 {
471 data.assign(data.size(), T(0.0));
472 }
473
474 inline size_t size() const
475 {
476 return data.size();
477 }
478
480 inline matrix2d()
481 : outer_length(0), inner_length(0)
482 {
483 this->resize(0, 0);
484 }
485
486 inline matrix2d(size_t const ol, size_t const il)
487 : outer_length(ol), inner_length(il)
488 {
489 this->resize(outer_length, inner_length);
490 reset();
491 }
492
494 inline matrix2d(matrix2d<T> const &m)
495 : outer_length(m.outer_length), inner_length(m.inner_length)
496 {
497 // reinitialize data and rows arrays
498 this->resize(outer_length, inner_length);
499 // copy data
500 data = m.data;
501 }
502
504 inline ~matrix2d() {
505 this->clear();
506 }
507
509 inline std::vector<T> &data_array()
510 {
511 return data;
512 }
513
515 inline std::vector<T> const &data_array() const
516 {
517 return data;
518 }
519
520 inline row & operator [] (size_t const i)
521 {
522 return rows[i];
523 }
524 inline row const & operator [] (size_t const i) const
525 {
526 return rows[i];
527 }
528
531 {
532 if ((outer_length != m.outer_length) || (inner_length != m.inner_length)){
533 this->clear();
534 outer_length = m.outer_length;
535 inner_length = m.inner_length;
536 this->resize(outer_length, inner_length);
537 }
538 data = m.data;
539 return *this;
540 }
541
543 inline T ** c_array() {
544 if (rows.size() > 0) {
545 return &(pointers[0]);
546 } else {
547 return NULL;
548 }
549 }
550
551 inline static void check_sizes(matrix2d<T> const &m1, matrix2d<T> const &m2)
552 {
553 if ((m1.outer_length != m2.outer_length) ||
554 (m1.inner_length != m2.inner_length)) {
555 cvm::error("Error: trying to perform an operation between "
556 "matrices of different sizes, "+
557 cvm::to_str(m1.outer_length)+"x"+
558 cvm::to_str(m1.inner_length)+" and "+
559 cvm::to_str(m2.outer_length)+"x"+
560 cvm::to_str(m2.inner_length)+".\n");
561 }
562 }
563
564 inline void operator += (matrix2d<T> const &m)
565 {
566 check_sizes(*this, m);
567 size_t i;
568 for (i = 0; i < data.size(); i++) {
569 data[i] += m.data[i];
570 }
571 }
572
573 inline void operator -= (matrix2d<T> const &m)
574 {
575 check_sizes(*this, m);
576 size_t i;
577 for (i = 0; i < data.size(); i++) {
578 data[i] -= m.data[i];
579 }
580 }
581
582 inline void operator *= (cvm::real a)
583 {
584 size_t i;
585 for (i = 0; i < data.size(); i++) {
586 data[i] *= a;
587 }
588 }
589
590 inline void operator /= (cvm::real a)
591 {
592 size_t i;
593 for (i = 0; i < data.size(); i++) {
594 data[i] /= a;
595 }
596 }
597
598 inline friend matrix2d<T> operator + (matrix2d<T> const &m1,
599 matrix2d<T> const &m2)
600 {
601 check_sizes(m1, m2);
602 matrix2d<T> result(m1.outer_length, m1.inner_length);
603 size_t i;
604 for (i = 0; i < m1.data.size(); i++) {
605 result.data[i] = m1.data[i] + m2.data[i];
606 }
607 return result;
608 }
609
610 inline friend matrix2d<T> operator - (matrix2d<T> const &m1,
611 matrix2d<T> const &m2)
612 {
613 check_sizes(m1, m2);
614 matrix2d<T> result(m1.outer_length, m1.inner_length);
615 size_t i;
616 for (i = 0; i < m1.data.size(); i++) {
617 result.data[i] = m1.data[i] - m1.data[i];
618 }
619 return result;
620 }
621
622 inline friend matrix2d<T> operator * (matrix2d<T> const &m, cvm::real a)
623 {
624 matrix2d<T> result(m.outer_length, m.inner_length);
625 size_t i;
626 for (i = 0; i < m.data.size(); i++) {
627 result.data[i] = m.data[i] * a;
628 }
629 return result;
630 }
631
632 inline friend matrix2d<T> operator * (cvm::real a, matrix2d<T> const &m)
633 {
634 return m * a;
635 }
636
637 inline friend matrix2d<T> operator / (matrix2d<T> const &m, cvm::real a)
638 {
639 matrix2d<T> result(m.outer_length, m.inner_length);
640 size_t i;
641 for (i = 0; i < m.data.size(); i++) {
642 result.data[i] = m.data[i] * a;
643 }
644 return result;
645 }
646
648 inline friend vector1d<T> operator * (vector1d<T> const &v,
649 matrix2d<T> const &m)
650 {
651 vector1d<T> result(m.inner_length);
652 if (m.outer_length != v.size()) {
653 cvm::error("Error: trying to multiply a vector and a matrix "
654 "of incompatible sizes, "+
655 cvm::to_str(v.size()) + " and " +
656 cvm::to_str(m.outer_length)+"x"+cvm::to_str(m.inner_length) +
657 ".\n");
658 } else {
659 size_t i, k;
660 for (i = 0; i < m.inner_length; i++) {
661 for (k = 0; k < m.outer_length; k++) {
662 result[i] += m[k][i] * v[k];
663 }
664 }
665 }
666 return result;
667 }
668
670 friend std::ostream & operator << (std::ostream &os,
671 matrix2d<T> const &m)
672 {
673 std::streamsize const w = os.width();
674 std::streamsize const p = os.precision();
675
676 os.width(2);
677 os << "( ";
678 size_t i;
679 for (i = 0; i < m.outer_length; i++) {
680 os << " ( ";
681 size_t j;
682 for (j = 0; j < m.inner_length-1; j++) {
683 os.width(w);
684 os.precision(p);
685 os << m[i][j] << " , ";
686 }
687 os.width(w);
688 os.precision(p);
689 os << m[i][m.inner_length-1] << " )";
690 }
691
692 os << " )";
693 return os;
694 }
695
696 inline std::string to_simple_string() const
697 {
698 if (this->size() == 0) return std::string("");
699 std::ostringstream os;
700 os.setf(std::ios::scientific, std::ios::floatfield);
701 os.precision(cvm::cv_prec);
702 os << (*this)[0];
703 size_t i;
704 for (i = 1; i < data.size(); i++) {
705 os << " " << data[i];
706 }
707 return os.str();
708 }
709
710 inline int from_simple_string(std::string const &s)
711 {
712 std::stringstream stream(s);
713 size_t i = 0;
714 while ((i < data.size()) && (stream >> data[i])) {
715 i++;
716 }
717 if (i < data.size()) {
718 return COLVARS_ERROR;
719 }
720 return COLVARS_OK;
721 }
722
723};
724
725
728
729public:
730
731 cvm::real x, y, z;
732
733 inline rvector()
734 {
735 reset();
736 }
737
739 inline void reset()
740 {
741 set(0.0);
742 }
743
744 inline rvector(cvm::real x_i, cvm::real y_i, cvm::real z_i)
745 {
746 set(x_i, y_i, z_i);
747 }
748
749 inline rvector(cvm::vector1d<cvm::real> const &v)
750 {
751 set(v[0], v[1], v[2]);
752 }
753
754 inline rvector(cvm::real t)
755 {
756 set(t);
757 }
758
760 inline void set(cvm::real value)
761 {
762 x = y = z = value;
763 }
764
766 inline void set(cvm::real x_i, cvm::real y_i, cvm::real z_i)
767 {
768 x = x_i;
769 y = y_i;
770 z = z_i;
771 }
772
774 inline cvm::real & operator [] (int i) {
775 return (i == 0) ? x : (i == 1) ? y : (i == 2) ? z : x;
776 }
777
779 inline cvm::real operator [] (int i) const {
780 return (i == 0) ? x : (i == 1) ? y : (i == 2) ? z : x;
781 }
782
783 inline cvm::vector1d<cvm::real> const as_vector() const
784 {
785 cvm::vector1d<cvm::real> result(3);
786 result[0] = this->x;
787 result[1] = this->y;
788 result[2] = this->z;
789 return result;
790 }
791
792 inline void operator += (cvm::rvector const &v)
793 {
794 x += v.x;
795 y += v.y;
796 z += v.z;
797 }
798
799 inline void operator -= (cvm::rvector const &v)
800 {
801 x -= v.x;
802 y -= v.y;
803 z -= v.z;
804 }
805
806 inline void operator *= (cvm::real v)
807 {
808 x *= v;
809 y *= v;
810 z *= v;
811 }
812
813 inline void operator /= (cvm::real const& v)
814 {
815 x /= v;
816 y /= v;
817 z /= v;
818 }
819
820 inline cvm::real norm2() const
821 {
822 return (x*x + y*y + z*z);
823 }
824
825 inline cvm::real norm() const
826 {
827 return cvm::sqrt(this->norm2());
828 }
829
830 inline cvm::rvector unit() const
831 {
832 const cvm::real n = this->norm();
833 return (n > 0. ? cvm::rvector(x, y, z)/n : cvm::rvector(1., 0., 0.));
834 }
835
836 static inline size_t output_width(size_t real_width)
837 {
838 return 3*real_width + 10;
839 }
840
841
842 static inline cvm::rvector outer(cvm::rvector const &v1,
843 cvm::rvector const &v2)
844 {
845 return cvm::rvector( v1.y*v2.z - v2.y*v1.z,
846 -v1.x*v2.z + v2.x*v1.z,
847 v1.x*v2.y - v2.x*v1.y);
848 }
849
850 friend inline cvm::rvector operator - (cvm::rvector const &v)
851 {
852 return cvm::rvector(-v.x, -v.y, -v.z);
853 }
854
855 friend inline cvm::rvector operator + (cvm::rvector const &v1,
856 cvm::rvector const &v2)
857 {
858 return cvm::rvector(v1.x + v2.x, v1.y + v2.y, v1.z + v2.z);
859 }
860 friend inline cvm::rvector operator - (cvm::rvector const &v1,
861 cvm::rvector const &v2)
862 {
863 return cvm::rvector(v1.x - v2.x, v1.y - v2.y, v1.z - v2.z);
864 }
865
867 friend inline cvm::real operator * (cvm::rvector const &v1,
868 cvm::rvector const &v2)
869 {
870 return v1.x * v2.x + v1.y * v2.y + v1.z * v2.z;
871 }
872
873 friend inline cvm::rvector operator * (cvm::real a, cvm::rvector const &v)
874 {
875 return cvm::rvector(a*v.x, a*v.y, a*v.z);
876 }
877
878 friend inline cvm::rvector operator * (cvm::rvector const &v, cvm::real a)
879 {
880 return cvm::rvector(a*v.x, a*v.y, a*v.z);
881 }
882
883 friend inline cvm::rvector operator / (cvm::rvector const &v, cvm::real a)
884 {
885 return cvm::rvector(v.x/a, v.y/a, v.z/a);
886 }
887
888 std::string to_simple_string() const;
889 int from_simple_string(std::string const &s);
890};
891
892
896
897public:
898
899 cvm::real xx, xy, xz, yx, yy, yz, zx, zy, zz;
900
902 inline rmatrix()
903 {
904 reset();
905 }
906
908 inline rmatrix(cvm::real xxi, cvm::real xyi, cvm::real xzi,
909 cvm::real yxi, cvm::real yyi, cvm::real yzi,
910 cvm::real zxi, cvm::real zyi, cvm::real zzi)
911 {
912 xx = xxi;
913 xy = xyi;
914 xz = xzi;
915 yx = yxi;
916 yy = yyi;
917 yz = yzi;
918 zx = zxi;
919 zy = zyi;
920 zz = zzi;
921 }
922
923
924 inline void reset()
925 {
926 xx = xy = xz = yx = yy = yz = zx = zy = zz = 0.0;
927 }
928
930 inline cvm::real determinant() const
931 {
932 return
933 ( xx * (yy*zz - zy*yz))
934 - (yx * (xy*zz - zy*xz))
935 + (zx * (xy*yz - yy*xz));
936 }
937
938 inline cvm::rmatrix transpose() const
939 {
940 return cvm::rmatrix(xx, yx, zx,
941 xy, yy, zy,
942 xz, yz, zz);
943 }
944
945 inline friend cvm::rvector operator * (cvm::rmatrix const &m,
946 cvm::rvector const &r)
947 {
948 return cvm::rvector(m.xx*r.x + m.xy*r.y + m.xz*r.z,
949 m.yx*r.x + m.yy*r.y + m.yz*r.z,
950 m.zx*r.x + m.zy*r.y + m.zz*r.z);
951 }
952};
953
954
955
959
960public:
961
962 cvm::real q0, q1, q2, q3;
963
966 : q0(0.0), q1(x), q2(y), q3(z)
967 {}
968
970 inline quaternion(cvm::real const qv[4])
971 : q0(qv[0]), q1(qv[1]), q2(qv[2]), q3(qv[3])
972 {}
973
976 cvm::real q1i,
977 cvm::real q2i,
978 cvm::real q3i)
979 : q0(q0i), q1(q1i), q2(q2i), q3(q3i)
980 {}
981
982 inline quaternion(cvm::vector1d<cvm::real> const &v)
983 : q0(v[0]), q1(v[1]), q2(v[2]), q3(v[3])
984 {}
985
990 cvm::real theta_in,
991 cvm::real psi_in)
992 {
993 q0 = ( (cvm::cos(phi_in/2.0)) * (cvm::cos(theta_in/2.0)) * (cvm::cos(psi_in/2.0)) +
994 (cvm::sin(phi_in/2.0)) * (cvm::sin(theta_in/2.0)) * (cvm::sin(psi_in/2.0)) );
995
996 q1 = ( (cvm::sin(phi_in/2.0)) * (cvm::cos(theta_in/2.0)) * (cvm::cos(psi_in/2.0)) -
997 (cvm::cos(phi_in/2.0)) * (cvm::sin(theta_in/2.0)) * (cvm::sin(psi_in/2.0)) );
998
999 q2 = ( (cvm::cos(phi_in/2.0)) * (cvm::sin(theta_in/2.0)) * (cvm::cos(psi_in/2.0)) +
1000 (cvm::sin(phi_in/2.0)) * (cvm::cos(theta_in/2.0)) * (cvm::sin(psi_in/2.0)) );
1001
1002 q3 = ( (cvm::cos(phi_in/2.0)) * (cvm::cos(theta_in/2.0)) * (cvm::sin(psi_in/2.0)) -
1003 (cvm::sin(phi_in/2.0)) * (cvm::sin(theta_in/2.0)) * (cvm::cos(psi_in/2.0)) );
1004 }
1005
1007 inline quaternion()
1008 {
1009 reset();
1010 }
1011
1013 inline void set(cvm::real value)
1014 {
1015 q0 = q1 = q2 = q3 = value;
1016 }
1017
1019 inline void reset()
1020 {
1021 set(0.0);
1022 }
1023
1026 inline void reset_rotation()
1027 {
1028 q0 = 1.0;
1029 q1 = q2 = q3 = 0.0;
1030 }
1031
1033 static inline size_t output_width(size_t real_width)
1034 {
1035 return 4*real_width + 13;
1036 }
1037
1038 std::string to_simple_string() const;
1039 int from_simple_string(std::string const &s);
1040
1042 friend std::ostream & operator << (std::ostream &os, cvm::quaternion const &q);
1044 friend std::istream & operator >> (std::istream &is, cvm::quaternion &q);
1045
1047 inline cvm::real & operator [] (int i) {
1048 switch (i) {
1049 case 0:
1050 return this->q0;
1051 case 1:
1052 return this->q1;
1053 case 2:
1054 return this->q2;
1055 case 3:
1056 return this->q3;
1057 default:
1058 cvm::error("Error: incorrect quaternion component.\n");
1059 return q0;
1060 }
1061 }
1062
1064 inline cvm::real operator [] (int i) const {
1065 switch (i) {
1066 case 0:
1067 return this->q0;
1068 case 1:
1069 return this->q1;
1070 case 2:
1071 return this->q2;
1072 case 3:
1073 return this->q3;
1074 default:
1075 cvm::error("Error: trying to access a quaternion "
1076 "component which is not between 0 and 3.\n");
1077 return 0.0;
1078 }
1079 }
1080
1081 inline cvm::vector1d<cvm::real> const as_vector() const
1082 {
1083 cvm::vector1d<cvm::real> result(4);
1084 result[0] = q0;
1085 result[1] = q1;
1086 result[2] = q2;
1087 result[3] = q3;
1088 return result;
1089 }
1090
1092 inline cvm::real norm2() const
1093 {
1094 return q0*q0 + q1*q1 + q2*q2 + q3*q3;
1095 }
1096
1098 inline cvm::real norm() const
1099 {
1100 return cvm::sqrt(this->norm2());
1101 }
1102
1105 {
1106 return cvm::quaternion(q0, -q1, -q2, -q3);
1107 }
1108
1109 inline void operator *= (cvm::real a)
1110 {
1111 q0 *= a; q1 *= a; q2 *= a; q3 *= a;
1112 }
1113
1114 inline void operator /= (cvm::real a)
1115 {
1116 q0 /= a; q1 /= a; q2 /= a; q3 /= a;
1117 }
1118
1119 inline void set_positive()
1120 {
1121 if (q0 > 0.0) return;
1122 q0 = -q0;
1123 q1 = -q1;
1124 q2 = -q2;
1125 q3 = -q3;
1126 }
1127
1128 inline void operator += (cvm::quaternion const &h)
1129 {
1130 q0+=h.q0; q1+=h.q1; q2+=h.q2; q3+=h.q3;
1131 }
1132 inline void operator -= (cvm::quaternion const &h)
1133 {
1134 q0-=h.q0; q1-=h.q1; q2-=h.q2; q3-=h.q3;
1135 }
1136
1139 {
1140 return cvm::rvector(q1, q2, q3);
1141 }
1142
1143
1144 friend inline cvm::quaternion operator + (cvm::quaternion const &h,
1145 cvm::quaternion const &q)
1146 {
1147 return cvm::quaternion(h.q0+q.q0, h.q1+q.q1, h.q2+q.q2, h.q3+q.q3);
1148 }
1149
1150 friend inline cvm::quaternion operator - (cvm::quaternion const &h,
1151 cvm::quaternion const &q)
1152 {
1153 return cvm::quaternion(h.q0-q.q0, h.q1-q.q1, h.q2-q.q2, h.q3-q.q3);
1154 }
1155
1159 cvm::quaternion const &q)
1160 {
1161 return cvm::quaternion(h.q0*q.q0 - h.q1*q.q1 - h.q2*q.q2 - h.q3*q.q3,
1162 h.q0*q.q1 + h.q1*q.q0 + h.q2*q.q3 - h.q3*q.q2,
1163 h.q0*q.q2 + h.q2*q.q0 + h.q3*q.q1 - h.q1*q.q3,
1164 h.q0*q.q3 + h.q3*q.q0 + h.q1*q.q2 - h.q2*q.q1);
1165 }
1166
1167 friend inline cvm::quaternion operator * (cvm::real c,
1168 cvm::quaternion const &q)
1169 {
1170 return cvm::quaternion(c*q.q0, c*q.q1, c*q.q2, c*q.q3);
1171 }
1172 friend inline cvm::quaternion operator * (cvm::quaternion const &q,
1173 cvm::real c)
1174 {
1175 return cvm::quaternion(q.q0*c, q.q1*c, q.q2*c, q.q3*c);
1176 }
1177 friend inline cvm::quaternion operator / (cvm::quaternion const &q,
1178 cvm::real c)
1179 {
1180 return cvm::quaternion(q.q0/c, q.q1/c, q.q2/c, q.q3/c);
1181 }
1182
1183
1186 inline cvm::rvector rotate(cvm::rvector const &v) const
1187 {
1188 return ( (*this) * cvm::quaternion(0.0, v.x, v.y, v.z) *
1189 this->conjugate() ).get_vector();
1190 }
1191
1195 {
1196 cvm::rvector const vq_rot = this->rotate(Q2.get_vector());
1197 return cvm::quaternion(Q2.q0, vq_rot.x, vq_rot.y, vq_rot.z);
1198 }
1199
1202 {
1203 cvm::rmatrix R;
1204
1205 R.xx = q0*q0 + q1*q1 - q2*q2 - q3*q3;
1206 R.yy = q0*q0 - q1*q1 + q2*q2 - q3*q3;
1207 R.zz = q0*q0 - q1*q1 - q2*q2 + q3*q3;
1208
1209 R.xy = 2.0 * (q1*q2 - q0*q3);
1210 R.xz = 2.0 * (q0*q2 + q1*q3);
1211
1212 R.yx = 2.0 * (q0*q3 + q1*q2);
1213 R.yz = 2.0 * (q2*q3 - q0*q1);
1214
1215 R.zx = 2.0 * (q1*q3 - q0*q2);
1216 R.zy = 2.0 * (q0*q1 + q2*q3);
1217
1218 return R;
1219 }
1220
1221
1225 cvm::rvector const &vec) const;
1226
1227
1230 inline cvm::real cosine(cvm::quaternion const &q) const
1231 {
1232 cvm::real const iprod = this->inner(q);
1233 return 2.0*iprod*iprod - 1.0;
1234 }
1235
1239 inline cvm::real dist2(cvm::quaternion const &Q2) const
1240 {
1241 cvm::real const cos_omega = this->q0*Q2.q0 + this->q1*Q2.q1 +
1242 this->q2*Q2.q2 + this->q3*Q2.q3;
1243
1244 cvm::real const omega = cvm::acos( (cos_omega > 1.0) ? 1.0 :
1245 ( (cos_omega < -1.0) ? -1.0 : cos_omega) );
1246
1247 // get the minimum distance: x and -x are the same quaternion
1248 if (cos_omega > 0.0)
1249 return omega * omega;
1250 else
1251 return (PI-omega) * (PI-omega);
1252 }
1253
1257 {
1258 cvm::real const cos_omega = this->q0*Q2.q0 + this->q1*Q2.q1 + this->q2*Q2.q2 + this->q3*Q2.q3;
1259 cvm::real const omega = cvm::acos( (cos_omega > 1.0) ? 1.0 :
1260 ( (cos_omega < -1.0) ? -1.0 : cos_omega) );
1261 cvm::real const sin_omega = cvm::sin(omega);
1262
1263 if (cvm::fabs(sin_omega) < 1.0E-14) {
1264 // return a null 4d vector
1265 return cvm::quaternion(0.0, 0.0, 0.0, 0.0);
1266 }
1267
1268 cvm::quaternion const
1270 (-1.0)*sin_omega*Q2.q1 + cos_omega*(this->q1-cos_omega*Q2.q1)/sin_omega,
1271 (-1.0)*sin_omega*Q2.q2 + cos_omega*(this->q2-cos_omega*Q2.q2)/sin_omega,
1272 (-1.0)*sin_omega*Q2.q3 + cos_omega*(this->q3-cos_omega*Q2.q3)/sin_omega);
1273
1274 if (cos_omega > 0.0) {
1276 } else {
1278 }
1279 }
1280
1283 inline void match(cvm::quaternion &Q2) const
1284 {
1285 cvm::real const cos_omega = this->q0*Q2.q0 + this->q1*Q2.q1 +
1286 this->q2*Q2.q2 + this->q3*Q2.q3;
1287 if (cos_omega < 0.0) Q2 *= -1.0;
1288 }
1289
1292 inline cvm::real inner(cvm::quaternion const &Q2) const
1293 {
1294 cvm::real const prod = this->q0*Q2.q0 + this->q1*Q2.q1 +
1295 this->q2*Q2.q2 + this->q3*Q2.q3;
1296 return prod;
1297 }
1298
1299
1300};
1301
1302#ifndef COLVARS_LAMMPS
1303namespace NR {
1304void diagonalize_matrix(cvm::real m[4][4],
1305 cvm::real eigval[4],
1306 cvm::real eigvec[4][4]);
1307}
1308#endif
1309
1310
1314{
1315private:
1318
1321
1324
1327
1330
1331public:
1334
1337
1338 template <typename T1, typename T2>
1339 friend struct rotation_derivative;
1340
1341 template<typename T1, typename T2>
1343 cvm::rotation &rot,
1344 const std::vector<T1> &pos1,
1345 const std::vector<T2> &pos2);
1346
1357 void calc_optimal_rotation(std::vector<atom_pos> const &pos1,
1358 std::vector<atom_pos> const &pos2);
1359 void calc_optimal_rotation(std::vector<cvm::atom> const &pos1,
1360 std::vector<atom_pos> const &pos2);
1361
1363 int init();
1364
1366 rotation();
1367
1369 rotation(cvm::quaternion const &qi);
1370
1372 rotation(cvm::real angle, cvm::rvector const &axis);
1373
1375 ~rotation();
1376
1378 inline cvm::rvector rotate(cvm::rvector const &v) const
1379 {
1380 return q.rotate(v);
1381 }
1382
1384 inline cvm::rotation inverse() const
1385 {
1386 return cvm::rotation(this->q.conjugate());
1387 }
1388
1390 inline cvm::rmatrix matrix() const
1391 {
1392 return q.rotation_matrix();
1393 }
1394
1397 inline cvm::real spin_angle(cvm::rvector const &axis) const
1398 {
1399 cvm::rvector const q_vec = q.get_vector();
1400 cvm::real alpha = (180.0/PI) * 2.0 * cvm::atan2(axis * q_vec, q.q0);
1401 while (alpha > 180.0) alpha -= 360;
1402 while (alpha < -180.0) alpha += 360;
1403 return alpha;
1404 }
1405
1409 {
1410 cvm::rvector const q_vec = q.get_vector();
1411 cvm::real const iprod = axis * q_vec;
1412
1413 if (q.q0 != 0.0) {
1414
1415 cvm::real const dspindx =
1416 (180.0/PI) * 2.0 * (1.0 / (1.0 + (iprod*iprod)/(q.q0*q.q0)));
1417
1418 return cvm::quaternion( dspindx * (iprod * (-1.0) / (q.q0*q.q0)),
1419 dspindx * ((1.0/q.q0) * axis.x),
1420 dspindx * ((1.0/q.q0) * axis.y),
1421 dspindx * ((1.0/q.q0) * axis.z));
1422 } else {
1423 // (1/(1+x^2)) ~ (1/x)^2
1424 // The documentation of spinAngle discourages its use when q_vec and
1425 // axis are not close
1426 return cvm::quaternion((180.0/PI) * 2.0 * ((-1.0)/iprod), 0.0, 0.0, 0.0);
1427 }
1428 }
1429
1432 inline cvm::real cos_theta(cvm::rvector const &axis) const
1433 {
1434 cvm::rvector const q_vec = q.get_vector();
1435 cvm::real const alpha =
1436 (180.0/PI) * 2.0 * cvm::atan2(axis * q_vec, q.q0);
1437
1438 cvm::real const cos_spin_2 = cvm::cos(alpha * (PI/180.0) * 0.5);
1439 cvm::real const cos_theta_2 = ( (cos_spin_2 != 0.0) ?
1440 (q.q0 / cos_spin_2) :
1441 (0.0) );
1442 // cos(2t) = 2*cos(t)^2 - 1
1443 return 2.0 * (cos_theta_2*cos_theta_2) - 1.0;
1444 }
1445
1448 {
1449 cvm::rvector const q_vec = q.get_vector();
1450 cvm::real const iprod = axis * q_vec;
1451
1452 cvm::real const cos_spin_2 = cvm::cos(cvm::atan2(iprod, q.q0));
1453
1454 if (q.q0 != 0.0) {
1455
1456 cvm::real const d_cos_theta_dq0 =
1457 (4.0 * q.q0 / (cos_spin_2*cos_spin_2)) *
1458 (1.0 - (iprod*iprod)/(q.q0*q.q0) / (1.0 + (iprod*iprod)/(q.q0*q.q0)));
1459
1460 cvm::real const d_cos_theta_dqn =
1461 (4.0 * q.q0 / (cos_spin_2*cos_spin_2) *
1462 (iprod/q.q0) / (1.0 + (iprod*iprod)/(q.q0*q.q0)));
1463
1464 return cvm::quaternion(d_cos_theta_dq0,
1465 d_cos_theta_dqn * axis.x,
1466 d_cos_theta_dqn * axis.y,
1467 d_cos_theta_dqn * axis.z);
1468 } else {
1469
1470 cvm::real const d_cos_theta_dqn =
1471 (4.0 / (cos_spin_2*cos_spin_2 * iprod));
1472
1473 return cvm::quaternion(0.0,
1474 d_cos_theta_dqn * axis.x,
1475 d_cos_theta_dqn * axis.y,
1476 d_cos_theta_dqn * axis.z);
1477 }
1478 }
1479
1484
1485protected:
1486
1491
1493 void build_correlation_matrix(std::vector<cvm::atom_pos> const &pos1,
1494 std::vector<cvm::atom_pos> const &pos2);
1495 void build_correlation_matrix(std::vector<cvm::atom> const &pos1,
1496 std::vector<cvm::atom_pos> const &pos2);
1497
1500
1503
1505 void *jacobi;
1506};
1507
1508
1509#endif
Definition: colvartypes.h:385
Arbitrary size array (two dimensions) suitable for linear algebra operations (i.e....
Definition: colvartypes.h:376
matrix2d(matrix2d< T > const &m)
Copy constructor.
Definition: colvartypes.h:494
T ** c_array()
Return the 2-d C array.
Definition: colvartypes.h:543
void clear()
Deallocation routine.
Definition: colvartypes.h:463
~matrix2d()
Destructor.
Definition: colvartypes.h:504
std::vector< T > const & data_array() const
Return a reference to the data.
Definition: colvartypes.h:515
matrix2d< T > & operator=(matrix2d< T > const &m)
Assignment.
Definition: colvartypes.h:530
matrix2d()
Default constructor.
Definition: colvartypes.h:480
std::vector< T > & data_array()
Return a reference to the data.
Definition: colvartypes.h:509
void resize(size_t const ol, size_t const il)
Allocation routine, used by all constructors.
Definition: colvartypes.h:420
friend std::ostream & operator<<(std::ostream &os, matrix2d< T > const &m)
Formatted output.
Definition: colvartypes.h:670
void reset()
Set all elements to zero.
Definition: colvartypes.h:469
1-dimensional vector of real numbers with four components and a quaternion algebra
Definition: colvartypes.h:958
cvm::real & operator[](int i)
Access the quaternion as a 4-d array (return a reference)
Definition: colvartypes.h:1047
cvm::quaternion position_derivative_inner(cvm::rvector const &pos, cvm::rvector const &vec) const
Multiply the given vector by the derivative of the given (rotated) position with respect to the quate...
Definition: colvartypes.cpp:141
void set_from_euler_angles(cvm::real phi_in, cvm::real theta_in, cvm::real psi_in)
Definition: colvartypes.h:989
cvm::real cosine(cvm::quaternion const &q) const
Return the cosine between the orientation frame associated to this quaternion and another.
Definition: colvartypes.h:1230
cvm::rvector get_vector() const
Return the vector component.
Definition: colvartypes.h:1138
friend std::istream & operator>>(std::istream &is, cvm::quaternion &q)
Formatted input operator.
Definition: colvartypes.cpp:123
cvm::quaternion rotate(cvm::quaternion const &Q2) const
Rotate Q2 through this quaternion (put it in the rotated reference frame)
Definition: colvartypes.h:1194
cvm::real norm() const
Norm of the quaternion.
Definition: colvartypes.h:1098
quaternion(cvm::real const qv[4])
Constructor component by component.
Definition: colvartypes.h:970
quaternion()
Default constructor.
Definition: colvartypes.h:1007
void reset()
Set all components to zero (null quaternion)
Definition: colvartypes.h:1019
friend std::ostream & operator<<(std::ostream &os, cvm::quaternion const &q)
Formatted output operator.
Definition: colvartypes.cpp:104
cvm::rvector rotate(cvm::rvector const &v) const
Rotate v through this quaternion (put it in the rotated reference frame)
Definition: colvartypes.h:1186
Definition: colvartypes.h:1256
void match(cvm::quaternion &Q2) const
Choose the closest between Q2 and -Q2 and save it back. Not required for dist2() and dist2_grad()
Definition: colvartypes.h:1283
static size_t output_width(size_t real_width)
Tell the number of characters required to print a quaternion, given that of a real number.
Definition: colvartypes.h:1033
cvm::real norm2() const
Square norm of the quaternion.
Definition: colvartypes.h:1092
cvm::rmatrix rotation_matrix() const
Return the 3x3 matrix associated to this quaternion.
Definition: colvartypes.h:1201
cvm::quaternion conjugate() const
Return the conjugate quaternion.
Definition: colvartypes.h:1104
void set(cvm::real value)
Set all components to a scalar.
Definition: colvartypes.h:1013
cvm::real inner(cvm::quaternion const &Q2) const
Inner product (as a 4-d vector) with Q2; requires match() if the largest overlap is looked for.
Definition: colvartypes.h:1292
quaternion(cvm::real q0i, cvm::real q1i, cvm::real q2i, cvm::real q3i)
Constructor component by component.
Definition: colvartypes.h:975
quaternion(cvm::real x, cvm::real y, cvm::real z)
Constructor from a 3-d vector.
Definition: colvartypes.h:965
friend cvm::quaternion operator*(cvm::quaternion const &h, cvm::quaternion const &q)
Provides the quaternion product. NOTE: for the inner product use: h.inner (q);
Definition: colvartypes.h:1158
cvm::real dist2(cvm::quaternion const &Q2) const
Square distance from another quaternion on the 4-dimensional unit sphere: returns the square of the a...
Definition: colvartypes.h:1239
void reset_rotation()
Set the q0 component to 1 and the others to 0 (quaternion representing no rotation)
Definition: colvartypes.h:1026
2-dimensional array of real numbers with three components along each dimension (works with colvarmodu...
Definition: colvartypes.h:895
rmatrix(cvm::real xxi, cvm::real xyi, cvm::real xzi, cvm::real yxi, cvm::real yyi, cvm::real yzi, cvm::real zxi, cvm::real zyi, cvm::real zzi)
Constructor component by component.
Definition: colvartypes.h:908
cvm::real determinant() const
Return the determinant.
Definition: colvartypes.h:930
rmatrix()
Default constructor.
Definition: colvartypes.h:902
A rotation between two sets of coordinates (for the moment a wrapper for colvarmodule::quaternion)
Definition: colvartypes.h:1314
~rotation()
Destructor.
Definition: colvartypes.cpp:263
cvm::rmatrix C
Correlation matrix C (3, 3)
Definition: colvartypes.h:1317
rotation()
Default constructor.
Definition: colvartypes.cpp:225
cvm::real cos_theta(cvm::rvector const &axis) const
Return the projection of the orientation vector onto a predefined axis.
Definition: colvartypes.h:1432
cvm::real S_eigval[4]
Eigenvalues of S.
Definition: colvartypes.h:1323
cvm::real S_backup[4][4]
Definition: colvartypes.h:1329
cvm::rvector rotate(cvm::rvector const &v) const
Return the rotated vector.
Definition: colvartypes.h:1378
int init()
Initialize member data.
Definition: colvartypes.cpp:216
cvm::real S[4][4]
Overlap matrix S (4, 4)
Definition: colvartypes.h:1320
cvm::quaternion dspin_angle_dq(cvm::rvector const &axis) const
Return the derivative of the spin angle with respect to the quaternion.
Definition: colvartypes.h:1408
friend void debug_gradients(cvm::rotation &rot, const std::vector< T1 > &pos1, const std::vector< T2 > &pos2)
Function for debugging gradients (allow using either std::vector<cvm::atom_pos> or std::vector<cvm::a...
Definition: colvar_rotation_derivative.h:538
cvm::quaternion dcos_theta_dq(cvm::rvector const &axis) const
Return the derivative of the tilt wrt the quaternion.
Definition: colvartypes.h:1447
cvm::real spin_angle(cvm::rvector const &axis) const
Return the spin angle (in degrees) with respect to the provided axis (which MUST be normalized)
Definition: colvartypes.h:1397
cvm::real S_eigvec[4][4]
Eigenvectors of S.
Definition: colvartypes.h:1326
void calc_optimal_rotation(std::vector< atom_pos > const &pos1, std::vector< atom_pos > const &pos2)
Calculate the optimal rotation and store the corresponding eigenvalue and eigenvector in the argument...
static cvm::real crossing_threshold
Threshold for the eigenvalue crossing test.
Definition: colvartypes.h:1483
cvm::quaternion q_old
Previous value of the rotation (used to warn the user when the structure changes too much,...
Definition: colvartypes.h:1490
static bool monitor_crossings
Whether to test for eigenvalue crossing.
Definition: colvartypes.h:1481
Definition: colvartypes.h:1333
void build_correlation_matrix(std::vector< cvm::atom_pos > const &pos1, std::vector< cvm::atom_pos > const &pos2)
Build the correlation matrix C (used by calc_optimal_rotation())
Definition: colvartypes.cpp:274
cvm::rmatrix matrix() const
Return the associated 3x3 matrix.
Definition: colvartypes.h:1390
void * jacobi
Pointer to instance of Jacobi solver.
Definition: colvartypes.h:1505
void compute_overlap_matrix()
Compute the overlap matrix S (used by calc_optimal_rotation())
Definition: colvartypes.cpp:313
void calc_optimal_rotation_impl()
Actual implementation of calc_optimal_rotation (and called by it)
Definition: colvartypes.cpp:406
cvm::rotation inverse() const
Return the inverse of this rotation.
Definition: colvartypes.h:1384
cvm::quaternion q
The rotation itself (implemented as a quaternion)
Definition: colvartypes.h:1336
vector of real numbers with three components
Definition: colvartypes.h:727
cvm::real & operator[](int i)
Access cartesian components by index.
Definition: colvartypes.h:774
void set(cvm::real x_i, cvm::real y_i, cvm::real z_i)
Assign all components.
Definition: colvartypes.h:766
void reset()
Set all components to zero.
Definition: colvartypes.h:739
void set(cvm::real value)
Set all components to a scalar.
Definition: colvartypes.h:760
friend cvm::real operator*(cvm::rvector const &v1, cvm::rvector const &v2)
Inner (dot) product.
Definition: colvartypes.h:867
Arbitrary size array (one dimensions) suitable for linear algebra operations (i.e....
Definition: colvartypes.h:37
void reset()
Set all elements to zero.
Definition: colvartypes.h:96
size_t output_width(size_t real_width) const
Formatted output.
Definition: colvartypes.h:282
cvm::real norm2() const
Squared norm.
Definition: colvartypes.h:228
void sliceassign(size_t const i1, size_t const i2, vector1d< T > const &v)
Assign a vector to a slice of this vector.
Definition: colvartypes.h:268
std::vector< T > & data_array()
Return a reference to the data.
Definition: colvartypes.h:79
vector1d & operator=(const vector1d &)=default
Explicit Copy assignement.
vector1d(const vector1d &)=default
Explicit Copy constructor.
vector1d(size_t const n, T const *t)
Constructor from C array.
Definition: colvartypes.h:52
vector1d(size_t const n=0)
Default constructor.
Definition: colvartypes.h:45
T * c_array()
Return a pointer to the data location.
Definition: colvartypes.h:69
vector1d< T > const slice(size_t const i1, size_t const i2) const
Slicing.
Definition: colvartypes.h:254
std::vector< T > const & data_array() const
Return a reference to the data.
Definition: colvartypes.h:85
double real
Defining an abstract real number allows to switch precision.
Definition: colvarmodule.h:95
static real cos(real const &x)
Reimplemented to work around MS compiler issues.
Definition: colvarmodule.h:145
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:2032
static real atan2(real const &x, real const &y)
Reimplemented to work around MS compiler issues.
Definition: colvarmodule.h:163
static size_t const cv_prec
Number of digits to represent a collective variables value(s)
Definition: colvarmodule.h:664
static real sqrt(real const &x)
Reimplemented to work around MS compiler issues.
Definition: colvarmodule.h:133
static real acos(real const &x)
Reimplemented to work around MS compiler issues.
Definition: colvarmodule.h:157
static real fabs(real const &x)
Reimplemented to work around MS compiler issues.
Definition: colvarmodule.h:127
static real sin(real const &x)
Reimplemented to work around MS compiler issues.
Definition: colvarmodule.h:139
static std::string to_str(char const *s)
Convert to string for output purposes.
Definition: colvarmodule.cpp:2378
Collective variables main module.
Helper class for calculating the derivative of rotation.
Definition: colvar_rotation_derivative.h:49