NAMD
Tensor.h
Go to the documentation of this file.
1 
7 #ifndef TENSOR_H
8 #define TENSOR_H
9 
10 #include <math.h>
11 #include <stdio.h>
12 #include "common.h"
13 #include "Vector.h"
14 
15 class Tensor {
16  public:
20 
21  inline Tensor(void) {
22  xx=xy=xz=yx=yy=yz=zx=zy=zz=0.0;
23  }
24 
25  inline Tensor(const Tensor &t2) {
26  xx = t2.xx; xy = t2.xy; xz = t2.xz;
27  yx = t2.yx; yy = t2.yy; yz = t2.yz;
28  zx = t2.zx; zy = t2.zy; zz = t2.zz;
29  }
30 
31  static inline Tensor identity(BigReal v1 = 1.0) {
32  Tensor tmp;
33  tmp.xx = tmp.yy = tmp.zz = v1;
34  return tmp;
35  }
36 
37  static inline Tensor diagonal(const Vector &v1) {
38  Tensor tmp;
39  tmp.xx = v1.x; tmp.xy = 0; tmp.xz = 0;
40  tmp.yx = 0; tmp.yy = v1.y; tmp.yz = 0;
41  tmp.zx = 0; tmp.zy = 0; tmp.zz = v1.z;
42  return tmp;
43  }
44 
45  static inline Tensor symmetric(const Vector &v1, const Vector &v2) {
46  Tensor tmp;
47  tmp.xx = v1.x; tmp.xy = v2.x; tmp.xz = v2.y;
48  tmp.yx = v2.x; tmp.yy = v1.y; tmp.yz = v2.z;
49  tmp.zx = v2.y; tmp.zy = v2.z; tmp.zz = v1.z;
50  return tmp;
51  }
52 
53  static inline Tensor triangular(const Vector &v1, const Vector &v2) {
54  Tensor tmp;
55  tmp.xx = v1.x; tmp.xy = v2.x; tmp.xz = v2.y;
56  tmp.yx = 0; tmp.yy = v1.y; tmp.yz = v2.z;
57  tmp.zx = 0; tmp.zy = 0; tmp.zz = v1.z;
58  return tmp;
59  }
60 
61  ~Tensor(void) { }
62 
63  inline Tensor& operator=(const Tensor &t2) {
64  xx = t2.xx; xy = t2.xy; xz = t2.xz;
65  yx = t2.yx; yy = t2.yy; yz = t2.yz;
66  zx = t2.zx; zy = t2.zy; zz = t2.zz;
67  return *this;
68  }
69 
70  inline Tensor& operator=(const BigReal &r2) {
71  xx=xy=xz=yx=yy=yz=zx=zy=zz=r2;
72  return *this;
73  }
74 
75  inline Tensor& operator+=(const Tensor &t2) {
76  xx += t2.xx; xy += t2.xy; xz += t2.xz;
77  yx += t2.yx; yy += t2.yy; yz += t2.yz;
78  zx += t2.zx; zy += t2.zy; zz += t2.zz;
79  return *this;
80  }
81 
82  inline Tensor& operator-=(const Tensor &t2) {
83  xx -= t2.xx; xy -= t2.xy; xz -= t2.xz;
84  yx -= t2.yx; yy -= t2.yy; yz -= t2.yz;
85  zx -= t2.zx; zy -= t2.zy; zz -= t2.zz;
86  return *this;
87  }
88 
89  inline Tensor& operator*=(const BigReal &r2) {
90  xx *= r2; xy *= r2; xz *= r2;
91  yx *= r2; yy *= r2; yz *= r2;
92  zx *= r2; zy *= r2; zz *= r2;
93  return *this;
94  }
95 
96  inline Tensor& operator/=(const BigReal &r2) {
97  xx /= r2; xy /= r2; xz /= r2;
98  yx /= r2; yy /= r2; yz /= r2;
99  zx /= r2; zy /= r2; zz /= r2;
100  return *this;
101  }
102 
103  inline void outerAdd (BigReal scale, const Vector &v1, const Vector &v2);
104 
105  inline friend int operator==(const Tensor &t1, const Tensor &t2) {
106  return (
107  t1.xx == t2.xx && t1.xy == t2.xy && t1.xz == t2.xz &&
108  t1.yx == t2.yx && t1.yy == t2.yy && t1.yz == t2.yz &&
109  t1.zx == t2.zx && t1.zy == t2.zy && t1.zz == t2.zz );
110  }
111 
112  inline friend int operator!=(const Tensor &t1, const Tensor &t2) {
113  return ( ! ( t1 == t2 ) );
114  }
115 
116  inline friend Tensor operator+(const Tensor& t1, const Tensor& t2) {
117  Tensor tmp(t1);
118  tmp += t2;
119  return tmp;
120  }
121 
122  inline friend Tensor operator-(const Tensor& t1, const Tensor& t2) {
123  Tensor tmp(t1);
124  tmp -= t2;
125  return tmp;
126  }
127 
128  inline friend Tensor operator-(const Tensor &t1) {
129  Tensor tmp(t1);
130  tmp *= -1.0;
131  return tmp;
132  }
133 
134  inline friend Tensor operator*(const BigReal &r1, const Tensor &t2) {
135  Tensor tmp(t2);
136  tmp *= r1;
137  return tmp;
138  }
139 
140  inline friend Tensor operator*(const Tensor &t1, const BigReal &r2) {
141  Tensor tmp(t1);
142  tmp *= r2;
143  return tmp;
144  }
145 
146  inline friend Tensor operator/(const Tensor &t1, const BigReal &r2) {
147  Tensor tmp(t1);
148  tmp /= r2;
149  return tmp;
150  }
151 
152  inline friend Vector operator*(const Tensor &t1, const Vector &v2) {
153  Vector tmp;
154  tmp.x = t1.xx * v2.x + t1.xy * v2.y + t1.xz * v2.z;
155  tmp.y = t1.yx * v2.x + t1.yy * v2.y + t1.yz * v2.z;
156  tmp.z = t1.zx * v2.x + t1.zy * v2.y + t1.zz * v2.z;
157  return tmp;
158  }
159 
160  inline friend Vector operator*(const Vector &v1, const Tensor &t2) {
161  Vector tmp;
162  tmp.x = t2.xx * v1.x + t2.yx * v1.y + t2.zx * v1.z;
163  tmp.y = t2.xy * v1.x + t2.yy * v1.y + t2.zy * v1.z;
164  tmp.z = t2.xz * v1.x + t2.yz * v1.y + t2.zz * v1.z;
165  return tmp;
166  }
167 
168  inline friend Tensor outer(const Vector &v1, const Vector &v2);
169 
170  inline friend Tensor transpose(const Tensor &t1) {
171  Tensor tmp;
172  tmp.xx = t1.xx; tmp.yx = t1.xy; tmp.zx = t1.xz;
173  tmp.xy = t1.yx; tmp.yy = t1.yy; tmp.zy = t1.yz;
174  tmp.xz = t1.zx; tmp.yz = t1.zy; tmp.zz = t1.zz;
175  return tmp;
176  }
177 
178  inline friend Tensor symmetric(const Tensor &t1) {
179  Tensor tmp;
180  tmp.xx = t1.xx; tmp.xy = 0.5*(t1.xy+t1.yx); tmp.xz = 0.5*(t1.xz+t1.zx);
181  tmp.yx = tmp.xy; tmp.yy = t1.yy; tmp.yz = 0.5*(t1.yz+t1.zy);
182  tmp.zx = tmp.xz; tmp.zy = tmp.yz; tmp.zz = t1.zz;
183  return tmp;
184  }
185 
186  inline friend Tensor triangular(const Tensor &t1) {
187  Tensor tmp;
188  tmp.xx = t1.xx; tmp.xy = 0.5*(t1.xy+t1.yx); tmp.xz = 0.5*(t1.xz+t1.zx);
189  tmp.yx = 0; tmp.yy = t1.yy; tmp.yz = 0.5*(t1.yz+t1.zy);
190  tmp.zx = 0; tmp.zy = 0; tmp.zz = t1.zz;
191  return tmp;
192  }
193 
194  inline friend Vector diagonal(const Tensor &t1) {
195  return Vector(t1.xx,t1.yy,t1.zz);
196  }
197 
198  inline friend Vector off_diagonal(const Tensor &t1) {
199  return Vector(t1.xy,t1.xz,t1.yz);
200  }
201 
202  inline friend BigReal trace(const Tensor &t1) {
203  return (t1.xx + t1.yy + t1.zz);
204  }
205 
206 /*
207  // set the vector based on a string. If bad, return FALSE
208  // the string can be in the form "x y z" or "x, y, z"
209  Bool set(const char *s) {
210  double a[3]; // read into doubles, since I don't know what
211  char tmp[100]; // a "BigReal" is in real life
212  // cheap way to get commas, etc. a poor regex
213  int i=sscanf(s, "%lf%99[ \t,]%lf%99[ \t,]%lf%99s",
214  a, tmp, a+1, tmp, a+2, tmp);
215  if (i != 5) return FALSE;
216  const char *t = s; // now count commas (for "1,,,,2, , 3")
217  int flg = 0; // and check for "1 2,,3"
218  i = 0;
219  for (;*t;t++) {
220  if (*t == ',') {
221  if (flg == 0) { // expecting non-whitespace
222  return FALSE; // so error
223  }
224  flg = 0; // now expect non-whitespace
225  i++; // and increment comma counter
226  }
227  else if (*t != ' ' && *t != '\t') { // got non-whitespace
228  flg = 1; // so the next can be whitespace or commas
229  }
230  }
231  if (i == 0 || i == 2) { // allow "1 2 3" or "1, 2,3" forms
232  x = a[0]; y = a[1]; z = a[2];
233  return TRUE;
234  }
235  return FALSE;
236  }
237 */
238 
239 };
240 
241 inline Tensor outer(const Vector &v1, const Vector &v2) {
242  Tensor tmp;
243  tmp.xx = v1.x * v2.x;
244  tmp.xy = v1.x * v2.y;
245  tmp.xz = v1.x * v2.z;
246  tmp.yx = v1.y * v2.x;
247  tmp.yy = v1.y * v2.y;
248  tmp.yz = v1.y * v2.z;
249  tmp.zx = v1.z * v2.x;
250  tmp.zy = v1.z * v2.y;
251  tmp.zz = v1.z * v2.z;
252  return tmp;
253 }
254 
255 inline void Tensor::outerAdd (BigReal scale, const Vector &v1, const Vector &v2) {
256  xx += v1.x * v2.x * scale;
257  xy += v1.x * v2.y * scale;
258  xz += v1.x * v2.z * scale;
259  yx += v1.y * v2.x * scale;
260  yy += v1.y * v2.y * scale;
261  yz += v1.y * v2.z * scale;
262  zx += v1.z * v2.x * scale;
263  zy += v1.z * v2.y * scale;
264  zz += v1.z * v2.z * scale;
265 }
266 
267 #endif
268 
friend Vector operator*(const Tensor &t1, const Vector &v2)
Definition: Tensor.h:152
BigReal zy
Definition: Tensor.h:19
BigReal xz
Definition: Tensor.h:17
static Tensor diagonal(const Vector &v1)
Definition: Tensor.h:37
Tensor(void)
Definition: Tensor.h:21
Definition: Vector.h:64
Tensor & operator=(const Tensor &t2)
Definition: Tensor.h:63
BigReal z
Definition: Vector.h:66
BigReal yz
Definition: Tensor.h:18
void outerAdd(BigReal scale, const Vector &v1, const Vector &v2)
Definition: Tensor.h:255
Tensor & operator*=(const BigReal &r2)
Definition: Tensor.h:89
friend Tensor operator+(const Tensor &t1, const Tensor &t2)
Definition: Tensor.h:116
friend Tensor operator-(const Tensor &t1, const Tensor &t2)
Definition: Tensor.h:122
static Tensor triangular(const Vector &v1, const Vector &v2)
Definition: Tensor.h:53
friend Tensor operator/(const Tensor &t1, const BigReal &r2)
Definition: Tensor.h:146
~Tensor(void)
Definition: Tensor.h:61
friend Tensor transpose(const Tensor &t1)
Definition: Tensor.h:170
BigReal yx
Definition: Tensor.h:18
Tensor & operator/=(const BigReal &r2)
Definition: Tensor.h:96
Tensor & operator=(const BigReal &r2)
Definition: Tensor.h:70
friend Tensor operator*(const Tensor &t1, const BigReal &r2)
Definition: Tensor.h:140
BigReal x
Definition: Vector.h:66
friend Tensor operator*(const BigReal &r1, const Tensor &t2)
Definition: Tensor.h:134
friend Tensor operator-(const Tensor &t1)
Definition: Tensor.h:128
static Tensor identity(BigReal v1=1.0)
Definition: Tensor.h:31
friend Tensor triangular(const Tensor &t1)
Definition: Tensor.h:186
Tensor outer(const Vector &v1, const Vector &v2)
Definition: Tensor.h:241
BigReal xx
Definition: Tensor.h:17
Tensor & operator+=(const Tensor &t2)
Definition: Tensor.h:75
Tensor(const Tensor &t2)
Definition: Tensor.h:25
friend Vector operator*(const Vector &v1, const Tensor &t2)
Definition: Tensor.h:160
BigReal zz
Definition: Tensor.h:19
friend Tensor symmetric(const Tensor &t1)
Definition: Tensor.h:178
friend Tensor outer(const Vector &v1, const Vector &v2)
Definition: Tensor.h:241
Definition: Tensor.h:15
BigReal xy
Definition: Tensor.h:17
friend Vector diagonal(const Tensor &t1)
Definition: Tensor.h:194
friend Vector off_diagonal(const Tensor &t1)
Definition: Tensor.h:198
BigReal y
Definition: Vector.h:66
Tensor & operator-=(const Tensor &t2)
Definition: Tensor.h:82
BigReal yy
Definition: Tensor.h:18
friend int operator!=(const Tensor &t1, const Tensor &t2)
Definition: Tensor.h:112
static Tensor symmetric(const Vector &v1, const Vector &v2)
Definition: Tensor.h:45
friend int operator==(const Tensor &t1, const Tensor &t2)
Definition: Tensor.h:105
BigReal zx
Definition: Tensor.h:19
double BigReal
Definition: common.h:112
friend BigReal trace(const Tensor &t1)
Definition: Tensor.h:202