PeTar
N-body code for collisional gravitational systems
tidal_tensor.hpp
Go to the documentation of this file.
1 #pragma once
2 
4 class TidalTensor{
5 private:
6  PS::F64 T1[3]; // 0 constant
7  PS::F64 T2[9]; // 1st (9) general tensor
8  //PS::F64 T2[6]; // 1st (6) symmetry tensor
9 #ifdef TIDAL_TENSOR_3RD
10  PS::F64 T3[10]; // 2nd Tensor (10)
11 #endif
12 public:
13  PS::F64vec pos; // position of c.m.
14  PS::F64 group_id; // indicate which group use the tensor
15 
16  TidalTensor(): T1{0.0, 0.0, 0.0},
17  T2{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
18 #ifdef TIDAL_TENSOR_3RD
19  T3{0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
20 #endif
21  pos(0.0), group_id(0.0) {}
22 
23  void dump(FILE *fp){
24  fwrite(this, sizeof(TidalTensor),1,fp);
25  }
26 
27  void read(FILE *fp){
28  size_t rcount = fread(this, sizeof(TidalTensor),1,fp);
29  if (rcount<1) {
30  std::cerr<<"Error: Data reading fails! requiring data number is 1, only obtain "<<rcount<<".\n";
31  abort();
32  }
33  }
34 
35  void clear(){
36  T1[0] = T1[1] = T1[2] = 0;
37  for(PS::S32 i=0; i<9; i++) T2[i] = 0;
38 #ifdef TIDAL_TENSOR_3RD
39  for(PS::S32 i=0; i<10; i++) T3[i] = 0;
40 #endif
41  pos = 0.0;
42  group_id = 0.0;
43  }
44 
46 
50  template<class Tptcl>
51  static void createTidalTensorMeasureParticles(Tptcl* _ptcl_tt, const Tptcl& _ptcl_cm, const PS::F64 _size) {
52 #ifdef TIDAL_TENSOR_3RD
53  // Then the edge length=_size/(2*sqrt(2))
55  // */
56 
57  PS::F64 lscale = 0.16*_size; // half edge size
58 
59  // set box
60  _ptcl_tt[0].pos = PS::F64vec(lscale, 0, -lscale) + _ptcl_cm.pos;
61  _ptcl_tt[1].pos = PS::F64vec(0, lscale, -lscale) + _ptcl_cm.pos;
62  _ptcl_tt[2].pos = PS::F64vec(-lscale, 0, -lscale) + _ptcl_cm.pos;
63  _ptcl_tt[3].pos = PS::F64vec(0, -lscale, -lscale) + _ptcl_cm.pos;
64  _ptcl_tt[4].pos = PS::F64vec(lscale, 0, lscale) + _ptcl_cm.pos;
65  _ptcl_tt[5].pos = PS::F64vec(0, lscale, lscale) + _ptcl_cm.pos;
66  _ptcl_tt[6].pos = PS::F64vec(-lscale, 0, lscale) + _ptcl_cm.pos;
67  _ptcl_tt[7].pos = PS::F64vec(0, -lscale, lscale) + _ptcl_cm.pos;
68 
69 #else
70  // Then the edge length=_size
72  // */
73 
74  PS::F64 lscale = 0.5*_size; // half edge size
75 
76  // set tetrahedron
77  _ptcl_tt[0].pos = PS::F64vec( lscale, 0, -lscale*0.707106781186548) + _ptcl_cm.pos;
78  _ptcl_tt[1].pos = PS::F64vec(-lscale, 0, -lscale*0.707106781186548) + _ptcl_cm.pos;
79  _ptcl_tt[2].pos = PS::F64vec(0, lscale, lscale*0.707106781186548) + _ptcl_cm.pos;
80  _ptcl_tt[3].pos = PS::F64vec(0, -lscale, lscale*0.707106781186548) + _ptcl_cm.pos;
81 
82 #endif
83 
84  for (int i=0; i<getParticleN(); i++) {
85  // co-moving velocity
86  _ptcl_tt[i].vel = _ptcl_cm.vel;
87  // no mass
88  _ptcl_tt[i].mass = 0.0;
89  }
90  }
91 
93  template<class Tptcl>
94  static void subtractCMForce(Tptcl* _ptcl_tt, const Tptcl& _ptcl_cm) {
95  for (int k=0; k<getParticleN(); k++) _ptcl_tt[k].acc -= _ptcl_cm.acc;
96  }
97 
99 
127  template<class Tptcl>
128  void fit(Tptcl* _ptcl_tt, Tptcl& _ptcl_cm, const PS::F64 _size) {
129  // get c.m. position
130  pos = _ptcl_cm.pos;
131 
132  PS::S32 n_point = getParticleN();
133 
134  PS::F64vec fi[n_point];
135 
136  // get acceleration
137  for (PS::S32 i=0; i<n_point; i++) fi[i] = _ptcl_tt[i].acc;
138 
139  // get cofficients
140  // T1, assume input force already remove the c.m.
141  T1[0] = T1[1] = T1[2] = 0.0;
142 
143 #ifdef TIDAL_TENSOR_3RD
144  // T2, general form
145  // 0 1 2
146  T2[0] = 0.250000000000000*fi[0][0] + -0.250000000000000*fi[2][0] + 0.250000000000000*fi[4][0] + -0.250000000000000*fi[6][0];
147  T2[1] = 0.125000000000000*fi[0][1] + 0.125000000000000*fi[1][0] + -0.125000000000000*fi[2][1] + -0.125000000000000*fi[3][0]
148  + 0.125000000000000*fi[4][1] + 0.125000000000000*fi[5][0] + -0.125000000000000*fi[6][1] + -0.125000000000000*fi[7][0];
149  T2[2] = -0.083333333333333*fi[0][0] + 0.083333333333333*fi[0][2] + -0.083333333333333*fi[1][0] + -0.083333333333333*fi[2][0]
150  + -0.083333333333333*fi[2][2] + -0.083333333333333*fi[3][0] + 0.083333333333333*fi[4][0] + 0.083333333333333*fi[4][2]
151  + 0.083333333333333*fi[5][0] + 0.083333333333333*fi[6][0] + -0.083333333333333*fi[6][2] + 0.083333333333333*fi[7][0];
152 
153  // 3 4 5
154  T2[3] = T2[1];
155  T2[4] = 0.250000000000000*fi[1][1] + -0.250000000000000*fi[3][1] + 0.250000000000000*fi[5][1] + -0.250000000000000*fi[7][1];
156  T2[5] = -0.083333333333333*fi[0][1] + -0.083333333333333*fi[1][1] + 0.083333333333333*fi[1][2] + -0.083333333333333*fi[2][1]
157  + -0.083333333333333*fi[3][1] + -0.083333333333333*fi[3][2] + 0.083333333333333*fi[4][1] + 0.083333333333333*fi[5][1]
158  + 0.083333333333333*fi[5][2] + 0.083333333333333*fi[6][1] + 0.083333333333333*fi[7][1] + -0.083333333333333*fi[7][2];
159 
160  // 6 7 8
161  T2[6] = T2[2];
162  T2[7] = T2[5];
163  T2[8] = -0.125000000000000*fi[0][2] + -0.125000000000000*fi[1][2] + -0.125000000000000*fi[2][2] + -0.125000000000000*fi[3][2]
164  + 0.125000000000000*fi[4][2] + 0.125000000000000*fi[5][2] + 0.125000000000000*fi[6][2] + 0.125000000000000*fi[7][2];
165 
166  // T3, symmetry form
167  T3[0] = 0.250000000000000*fi[0][0] + 0.125000000000000*fi[0][2] + 0.250000000000000*fi[2][0] + -0.125000000000000*fi[2][2]
168  + 0.250000000000000*fi[4][0] + -0.125000000000000*fi[4][2] + 0.250000000000000*fi[6][0] + 0.125000000000000*fi[6][2];
169  T3[1] = 0.250000000000000*fi[0][1] + 0.125000000000000*fi[1][2] + 0.250000000000000*fi[2][1] + -0.125000000000000*fi[3][2]
170  + 0.250000000000000*fi[4][1] + -0.125000000000000*fi[5][2] + 0.250000000000000*fi[6][1] + 0.125000000000000*fi[7][2];
171  T3[2] = -0.112500000000000*fi[0][0] + 0.025000000000000*fi[0][2] + -0.012500000000000*fi[1][1] + -0.025000000000000*fi[1][2]
172  + 0.112500000000000*fi[2][0] + 0.025000000000000*fi[2][2] + 0.012500000000000*fi[3][1] + -0.025000000000000*fi[3][2]
173  + 0.112500000000000*fi[4][0] + 0.025000000000000*fi[4][2] + 0.012500000000000*fi[5][1] + -0.025000000000000*fi[5][2]
174  + -0.112500000000000*fi[6][0] + 0.025000000000000*fi[6][2] + -0.012500000000000*fi[7][1] + -0.025000000000000*fi[7][2];
175  T3[3] = 0.125000000000000*fi[0][2] + 0.250000000000000*fi[1][0] + -0.125000000000000*fi[2][2] + 0.250000000000000*fi[3][0]
176  + -0.125000000000000*fi[4][2] + 0.250000000000000*fi[5][0] + 0.125000000000000*fi[6][2] + 0.250000000000000*fi[7][0];
177  T3[4] = -0.062500000000000*fi[0][1] + -0.062500000000000*fi[1][0] + 0.062500000000000*fi[2][1] + 0.062500000000000*fi[3][0]
178  + 0.062500000000000*fi[4][1] + 0.062500000000000*fi[5][0] + -0.062500000000000*fi[6][1] + -0.062500000000000*fi[7][0];
179  T3[5] = -0.125000000000000*fi[0][2] + 0.125000000000000*fi[2][2] + 0.125000000000000*fi[4][2] + -0.125000000000000*fi[6][2];
180  T3[6] = 0.250000000000000*fi[1][1] + 0.125000000000000*fi[1][2] + 0.250000000000000*fi[3][1] + -0.125000000000000*fi[3][2]
181  + 0.250000000000000*fi[5][1] + -0.125000000000000*fi[5][2] + 0.250000000000000*fi[7][1] + 0.125000000000000*fi[7][2];
182  T3[7] = -0.012500000000000*fi[0][0] + -0.025000000000000*fi[0][2] + -0.112500000000000*fi[1][1] + 0.025000000000000*fi[1][2]
183  + 0.012500000000000*fi[2][0] + -0.025000000000000*fi[2][2] + 0.112500000000000*fi[3][1] + 0.025000000000000*fi[3][2]
184  + 0.012500000000000*fi[4][0] + -0.025000000000000*fi[4][2] + 0.112500000000000*fi[5][1] + 0.025000000000000*fi[5][2]
185  + -0.012500000000000*fi[6][0] + -0.025000000000000*fi[6][2] + -0.112500000000000*fi[7][1] + 0.025000000000000*fi[7][2];
186  T3[8] = -0.125000000000000*fi[1][2] + 0.125000000000000*fi[3][2] + 0.125000000000000*fi[5][2] + -0.125000000000000*fi[7][2];
187  T3[9] = 0.062500000000000*fi[0][0] + 0.125000000000000*fi[0][2] + 0.062500000000000*fi[1][1] + 0.125000000000000*fi[1][2]
188  + -0.062500000000000*fi[2][0] + 0.125000000000000*fi[2][2] + -0.062500000000000*fi[3][1] + 0.125000000000000*fi[3][2]
189  + -0.062500000000000*fi[4][0] + 0.125000000000000*fi[4][2] + -0.062500000000000*fi[5][1] + 0.125000000000000*fi[5][2]
190  + 0.062500000000000*fi[6][0] + 0.125000000000000*fi[6][2] + 0.062500000000000*fi[7][1] + 0.125000000000000*fi[7][2];
191 
192  // Rescale
193  //PS::F64 T2S = 1.0/(_bin.semi*(1+_bin.ecc)*0.35);
194  PS::F64 T2S = 1.0/(_size*0.16);
195  PS::F64 T3S = T2S*T2S;
196  for (PS::S32 i=0; i<9; i++) T2[i] *= T2S;
197  for (PS::S32 i=0; i<10; i++) T3[i] *= T3S;
198 
199 #else
200 
201  // T2, general form
202  // 0 1 2
203  T2[0] = 0.500000000000000*fi[0][0]+ -0.500000000000000*fi[1][0];
204  T2[1] = 0.250000000000000*fi[0][1]+ -0.250000000000000*fi[1][1]+ 0.250000000000000*fi[2][0]+ -0.250000000000000*fi[3][0];
205  T2[2] = -0.176776695296637*fi[0][0]+ 0.250000000000000*fi[0][2]+ -0.176776695296637*fi[1][0]+ -0.250000000000000*fi[1][2]+ 0.176776695296637*fi[2][0]+ 0.176776695296637*fi[3][0];
206 
207  // 3 4 5
208  T2[3] = T2[1];
209  T2[4] = 0.500000000000000*fi[2][1]+ -0.500000000000000*fi[3][1];
210  T2[5] = -0.176776695296637*fi[0][1]+ -0.176776695296637*fi[1][1]+ 0.176776695296637*fi[2][1]+ 0.250000000000000*fi[2][2]+ 0.176776695296637*fi[3][1]+ -0.250000000000000*fi[3][2];
211 
212  // 6 7 8
213  T2[6] = T2[2];
214  T2[7] = T2[5];
215  T2[8] = -0.353553390593274*fi[0][2]+ -0.353553390593274*fi[1][2]+ 0.353553390593274*fi[2][2]+ 0.353553390593274*fi[3][2];
216 
217  // Rescale
218  PS::F64 T2S = 2.0/_size;
219  for (PS::S32 i=0; i<9; i++) T2[i] *= T2S;
220 #endif
221  }
222 
224 
247  void shiftCM(const PS::F64vec & _pos) {
248 #ifdef TIDAL_TENSOR_3RD
249  PS::F64vec dr = _pos-pos;
250 
251  PS::F64 x = dr.x;
252  PS::F64 y = dr.y;
253  PS::F64 z = dr.z;
254  //PS::F64 x2 = x*x;
255  //PS::F64 xy = x*y;
256  //PS::F64 xz = x*z;
257  //PS::F64 y2 = y*y;
258  //PS::F64 yz = y*z;
259  //PS::F64 z2 = z*z;
260 
261  // T1 += T2^dr + dr^T3^dr
271  // T2 += 2*dr^T3
272  T2[0] += 2.0*(T3[0]*x + T3[1]*y + T3[2]*z); // xx: xxx*x + xyx*y + xzx*z
273  T2[1] += 2.0*(T3[1]*x + T3[3]*y + T3[4]*z); // xy: xxy*x + xyy*y + xzy*z
274  T2[2] += 2.0*(T3[2]*x + T3[4]*y + T3[5]*z); // xy: xxz*x + xyz*y + xzz*z
275 
276  T2[3] += 2.0*(T3[1]*x + T3[3]*y + T3[4]*z); // yx: yxx*x + yyx*y + yzx*z
277  T2[4] += 2.0*(T3[3]*x + T3[6]*y + T3[7]*z); // yy: yxy*x + yyy*y + yzy*z
278  T2[5] += 2.0*(T3[4]*x + T3[7]*y + T3[8]*z); // yy: yxz*x + yyz*y + yzz*z
279 
280  T2[6] += 2.0*(T3[2]*x + T3[4]*y + T3[5]*z); // zx: zxx*x + zyx*y + zzx*z
281  T2[7] += 2.0*(T3[4]*x + T3[7]*y + T3[8]*z); // zy: zxy*x + zyy*y + zzy*z
282  T2[8] += 2.0*(T3[5]*x + T3[8]*y + T3[9]*z); // zy: zxz*x + zyz*y + zzz*z
283 
284 #endif
285  // update c.m.
286  pos = _pos;
287  }
288 
289  void eval(PS::F64* acc, const PS::F64vec &pos) const {
290  PS::F64 acc0=acc[0];
291  PS::F64 acc1=acc[1];
292  PS::F64 acc2=acc[2];
293 
294  PS::F64 x = pos.x;
295  PS::F64 y = pos.y;
296  PS::F64 z = pos.z;
297 
298 #ifdef TIDAL_TENSOR_3RD
299  /*
300  T2:
301  [[0 1 2]
302  [3 4 5]
303  [6 7 8]]
304 
305  T3:
306  [[[6 7 8 ]
307  [7 9 10]
308  [8 10 11]]
309 
310  [[7 9 10]
311  [9 12 13]
312  [10 13 14]]
313 
314  [[8 10 11]
315  [10 13 14]
316  [11 14 15]]]
317 
318  */
319  PS::F64 x2 = x*x;
320  PS::F64 xy = x*y;
321  PS::F64 xz = x*z;
322  PS::F64 y2 = y*y;
323  PS::F64 yz = y*z;
324  PS::F64 z2 = z*z;
325 
326  acc0 += T1[0] + T2[0]*x + T2[1]*y + T2[2]*z
327  + T3[0]*x2 + 2*T3[1]*xy + 2*T3[2]*xz + T3[3]*y2 + 2*T3[4]*yz + T3[5]*z2;
328  acc1 += T1[1] + T2[3]*x + T2[4]*y + T2[5]*z
329  + T3[1]*x2 + 2*T3[3]*xy + 2*T3[4]*xz + T3[6]*y2 + 2*T3[7]*yz + T3[8]*z2;
330  acc2 += T1[2] + T2[6]*x + T2[7]*y + T2[8]*z
331  + T3[2]*x2 + 2*T3[4]*xy + 2*T3[5]*xz + T3[7]*y2 + 2*T3[8]*yz + T3[9]*z2;
332 
333 #else
334  acc0 += T1[0] + T2[0]*x + T2[1]*y + T2[2]*z;
335  acc1 += T1[1] + T2[3]*x + T2[4]*y + T2[5]*z;
336  acc2 += T1[2] + T2[6]*x + T2[7]*y + T2[8]*z;
337 #endif
338 
339  acc[0] = acc0;
340  acc[1] = acc1;
341  acc[2] = acc2;
342  }
343 
344  PS::F64 evalPot(const PS::F64vec &pos) const {
345  PS::F64 x = pos.x;
346  PS::F64 y = pos.y;
347  PS::F64 z = pos.z;
348 #ifdef TIDAL_TENSOR_3RD
349  PS::F64 x2 = x*x;
350  PS::F64 xy = x*y;
351  PS::F64 xz = x*z;
352  PS::F64 y2 = y*y;
353  PS::F64 yz = y*z;
354  PS::F64 z2 = z*z;
355 
356  PS::F64 acc0 = T1[0] + 0.5*(T2[0]*x + T2[1]*y + T2[2]*z)
357  + (T3[0]*x2 + 2*T3[1]*xy + 2*T3[2]*xz + T3[3]*y2 + 2*T3[4]*yz + T3[5]*z2)/3.0;
358  PS::F64 acc1 = T1[1] + 0.5*(T2[3]*x + T2[4]*y + T2[5]*z)
359  + (T3[1]*x2 + 2*T3[3]*xy + 2*T3[4]*xz + T3[6]*y2 + 2*T3[7]*yz + T3[8]*z2)/3.0;
360  PS::F64 acc2 = T1[2] + 0.5*(T2[6]*x + T2[7]*y + T2[8]*z)
361  + (T3[2]*x2 + 2*T3[4]*xy + 2*T3[5]*xz + T3[7]*y2 + 2*T3[8]*yz + T3[9]*z2)/3.0;
362 #else
363  PS::F64 acc0 = T1[0] + 0.5*(T2[0]*x + T2[1]*y + T2[2]*z);
364  PS::F64 acc1 = T1[1] + 0.5*(T2[3]*x + T2[4]*y + T2[5]*z);
365  PS::F64 acc2 = T1[2] + 0.5*(T2[6]*x + T2[7]*y + T2[8]*z);
366 #endif
367 
368  return - x*acc0 - y*acc1 - z*acc2;
369  }
370 
371  void print(std::ostream & _fout, const int _width) const{
372  _fout<<"T1: \n"
373  <<std::setw(_width)<<T1[0]
374  <<std::setw(_width)<<T1[1]
375  <<std::setw(_width)<<T1[2]
376  <<std::endl
377  <<"T2: \n"
378  <<std::setw(_width)<<T2[0]
379  <<std::setw(_width)<<T2[1]
380  <<std::setw(_width)<<T2[2]
381  <<std::endl
382  <<std::setw(_width)<<T2[3]
383  <<std::setw(_width)<<T2[4]
384  <<std::setw(_width)<<T2[5]
385  <<std::endl
386  <<std::setw(_width)<<T2[6]
387  <<std::setw(_width)<<T2[7]
388  <<std::setw(_width)<<T2[8]
389 #ifdef TIDAL_TENSOR_3RD
390  <<std::endl
391  <<"T3: \n"
392  <<"x: \n"
393  <<std::setw(_width)<<T3[0]
394  <<std::setw(_width)<<T3[1]
395  <<std::setw(_width)<<T3[2]
396  <<std::endl
397  <<std::setw(_width)<<T3[1]
398  <<std::setw(_width)<<T3[3]
399  <<std::setw(_width)<<T3[4]
400  <<std::endl
401  <<std::setw(_width)<<T3[2]
402  <<std::setw(_width)<<T3[4]
403  <<std::setw(_width)<<T3[5]
404  <<std::endl
405  <<"y: \n"
406  <<std::setw(_width)<<T3[1]
407  <<std::setw(_width)<<T3[3]
408  <<std::setw(_width)<<T3[4]
409  <<std::endl
410  <<std::setw(_width)<<T3[3]
411  <<std::setw(_width)<<T3[6]
412  <<std::setw(_width)<<T3[7]
413  <<std::endl
414  <<std::setw(_width)<<T3[4]
415  <<std::setw(_width)<<T3[7]
416  <<std::setw(_width)<<T3[8]
417  <<std::endl
418  <<"x: \n"
419  <<std::setw(_width)<<T3[2]
420  <<std::setw(_width)<<T3[4]
421  <<std::setw(_width)<<T3[5]
422  <<std::endl
423  <<std::setw(_width)<<T3[4]
424  <<std::setw(_width)<<T3[7]
425  <<std::setw(_width)<<T3[8]
426  <<std::endl
427  <<std::setw(_width)<<T3[5]
428  <<std::setw(_width)<<T3[8]
429  <<std::setw(_width)<<T3[9]
430 #endif
431  <<std::endl;
432  }
433 
436 #ifdef TIDAL_TENSOR_3RD
437  return 8;
438 #else
439  return 4;
440 #endif
441  }
442 };
443 
TidalTensor::getParticleN
static PS::S32 getParticleN()
get particle number
Definition: tidal_tensor.hpp:435
TidalTensor::createTidalTensorMeasureParticles
static void createTidalTensorMeasureParticles(Tptcl *_ptcl_tt, const Tptcl &_ptcl_cm, const PS::F64 _size)
create tidal tensor measurement particles
Definition: tidal_tensor.hpp:51
TidalTensor::read
void read(FILE *fp)
Definition: tidal_tensor.hpp:27
PIKG::S32
int32_t S32
Definition: pikg_vector.hpp:24
galpy_pot_movie.fp
fp
Definition: galpy_pot_movie.py:131
PIKG::F64vec
Vector3< F64 > F64vec
Definition: pikg_vector.hpp:167
TidalTensor::dump
void dump(FILE *fp)
Definition: tidal_tensor.hpp:23
TidalTensor::group_id
PS::F64 group_id
Definition: tidal_tensor.hpp:14
PIKG::F64
double F64
Definition: pikg_vector.hpp:17
TidalTensor::shiftCM
void shiftCM(const PS::F64vec &_pos)
Shift c.m. to new reference position.
Definition: tidal_tensor.hpp:247
TidalTensor::fit
void fit(Tptcl *_ptcl_tt, Tptcl &_ptcl_cm, const PS::F64 _size)
tidal tensor fitting function,
Definition: tidal_tensor.hpp:128
TidalTensor::print
void print(std::ostream &_fout, const int _width) const
Definition: tidal_tensor.hpp:371
TidalTensor::TidalTensor
TidalTensor()
Definition: tidal_tensor.hpp:16
TidalTensor::subtractCMForce
static void subtractCMForce(Tptcl *_ptcl_tt, const Tptcl &_ptcl_cm)
subtract c.m. force from measure points
Definition: tidal_tensor.hpp:94
TidalTensor
Tidal tensor perterbation for AR.
Definition: tidal_tensor.hpp:4
TidalTensor::eval
void eval(PS::F64 *acc, const PS::F64vec &pos) const
Definition: tidal_tensor.hpp:289
TidalTensor::evalPot
PS::F64 evalPot(const PS::F64vec &pos) const
Definition: tidal_tensor.hpp:344
TidalTensor::clear
void clear()
Definition: tidal_tensor.hpp:35
TidalTensor::pos
PS::F64vec pos
Definition: tidal_tensor.hpp:13