changa 3.5
Loading...
Searching...
No Matches
MultipoleMoments.h
Go to the documentation of this file.
1
7
8#ifndef MULTIPOLEMOMENTS_H
9#define MULTIPOLEMOMENTS_H
10
11#include <cmath>
12#include <assert.h>
13#include <pup.h>
14
15#include <OrientedBox.h>
16#include <Vector3D.h>
17#ifdef HEXADECAPOLE
18#include "moments.h"
19#endif
20
21#include "SSEdefs.h"
22
23#if CMK_SSE && defined(HEXADECAPOLE)
24/*
25 ** This is a new fast version of QEVAL which evaluates
26 ** the interaction due to the reduced moment 'm'.
27 ** This version is nearly two times as fast as a naive
28 ** implementation.
29 **
30 ** OpCount = (*,+) = (103,72) = 175 - 8 = 167
31 */
32inline
33void momEvalMomr(MOMR *m,SSEcosmoType dir0,SSEcosmoType x,SSEcosmoType y,
34 SSEcosmoType z,SSEcosmoType *fPot,SSEcosmoType *ax,
35 SSEcosmoType *ay,SSEcosmoType *az)
36{
37 const SSEcosmoType onethird = 1.0/3.0;
38 SSEcosmoType xx,xy,xz,yy,yz,zz;
39 SSEcosmoType xxx,xxy,xxz,xyy,yyy,yyz,xyz;
40 SSEcosmoType tx,ty,tz,dir2,g2,g3,g4;
41 SSEcosmoType dir;
42
43 dir = -dir0;
44 dir2 = dir*dir;
45 g2 = 3.0*dir*dir2*dir2;
46 g3 = -5.0*g2*dir2;
47 g4 = -7.0*g3*dir2;
48 /*
49 ** Calculate the funky distance terms.
50 */
51 xx = 0.5*x*x;
52 xy = x*y;
53 xz = x*z;
54 yy = 0.5*y*y;
55 yz = y*z;
56 zz = 0.5*z*z;
57 xxx = x*(onethird*xx - zz);
58 xxz = z*(xx - onethird*zz);
59 yyy = y*(onethird*yy - zz);
60 yyz = z*(yy - onethird*zz);
61 xx -= zz;
62 yy -= zz;
63 xxy = y*xx;
64 xyy = x*yy;
65 xyz = xy*z;
66 /*
67 ** Now calculate the interaction up to Hexadecapole order.
68 */
69 tx = g4*(m->xxxx*xxx + m->xyyy*yyy + m->xxxy*xxy + m->xxxz*xxz + m->xxyy*xyy + m->xxyz*xyz + m->xyyz*yyz);
70 ty = g4*(m->xyyy*xyy + m->xxxy*xxx + m->yyyy*yyy + m->yyyz*yyz + m->xxyy*xxy + m->xxyz*xxz + m->xyyz*xyz);
71 tz = g4*(-m->xxxx*xxz - (m->xyyy + m->xxxy)*xyz - m->yyyy*yyz + m->xxxz*xxx + m->yyyz*yyy - m->xxyy*(xxz + yyz) + m->xxyz*xxy + m->xyyz*xyy);
72 g4 = 0.25*(tx*x + ty*y + tz*z);
73 xxx = g3*(m->xxx*xx + m->xyy*yy + m->xxy*xy + m->xxz*xz + m->xyz*yz);
74 xxy = g3*(m->xyy*xy + m->xxy*xx + m->yyy*yy + m->yyz*yz + m->xyz*xz);
75 xxz = g3*(-(m->xxx + m->xyy)*xz - (m->xxy + m->yyy)*yz + m->xxz*xx + m->yyz*yy + m->xyz*xy);
76 g3 = onethird*(xxx*x + xxy*y + xxz*z);
77 xx = g2*(m->xx*x + m->xy*y + m->xz*z);
78 xy = g2*(m->yy*y + m->xy*x + m->yz*z);
79 xz = g2*(-(m->xx + m->yy)*z + m->xz*x + m->yz*y);
80 g2 = 0.5*(xx*x + xy*y + xz*z);
81 dir *= m->m;
82 dir2 *= -(dir + 5.0*g2 + 7.0*g3 + 9.0*g4);
83 *fPot += dir + g2 + g3 + g4;
84 *ax += xx + xxx + tx + x*dir2;
85 *ay += xy + xxy + ty + y*dir2;
86 *az += xz + xxz + tz + z*dir2;
87 }
88/*
89 ** This is a new fast version of QEVAL which evaluates
90 ** the interaction due to the reduced moment 'm'.
91 ** This version is nearly two times as fast as a naive
92 ** implementation.
93 **
94 ** March 23, 2007: This function now uses unit vectors
95 ** which reduces the required precision in the exponent
96 ** since the highest power of r is now 5 (g4 ~ r^(-5)).
97 **
98 ** OpCount = (*,+) = (106,72) = 178 - 8 = 170
99 **
100 */
101inline
102void momEvalFmomrcm(FMOMR *m,SSEcosmoType u,SSEcosmoType dir,SSEcosmoType x,
103 SSEcosmoType y,SSEcosmoType z,
104 SSEcosmoType *fPot,SSEcosmoType *ax,SSEcosmoType *ay,
105 SSEcosmoType *az,SSEcosmoType *magai) {
106 const SSEcosmoType onethird = 1.0f/3.0f;
107 SSEcosmoType xx,xy,xz,yy,yz,zz;
108 SSEcosmoType xxx,xxy,xxz,xyy,yyy,yyz,xyz;
109 SSEcosmoType tx,ty,tz,g0,g2,g3,g4;
110
111 u *= dir;
112 g0 = dir;
113 g2 = 3*dir*u*u;
114 g3 = 5*g2*u;
115 g4 = 7*g3*u;
116 /*
117 ** Calculate the trace-free distance terms.
118 */
119 x *= dir;
120 y *= dir;
121 z *= dir;
122 xx = 0.5*x*x;
123 xy = x*y;
124 xz = x*z;
125 yy = 0.5*y*y;
126 yz = y*z;
127 zz = 0.5*z*z;
128 xxx = x*(onethird*xx - zz);
129 xxz = z*(xx - onethird*zz);
130 yyy = y*(onethird*yy - zz);
131 yyz = z*(yy - onethird*zz);
132 xx -= zz;
133 yy -= zz;
134 xxy = y*xx;
135 xyy = x*yy;
136 xyz = xy*z;
137 /*
138 ** Now calculate the interaction up to Hexadecapole order.
139 */
140 tx = g4*(m->xxxx*xxx + m->xyyy*yyy + m->xxxy*xxy + m->xxxz*xxz + m->xxyy*xyy + m->xxyz*xyz + m->xyyz*yyz);
141 ty = g4*(m->xyyy*xyy + m->xxxy*xxx + m->yyyy*yyy + m->yyyz*yyz + m->xxyy*xxy + m->xxyz*xxz + m->xyyz*xyz);
142 tz = g4*(-m->xxxx*xxz - (m->xyyy + m->xxxy)*xyz - m->yyyy*yyz + m->xxxz*xxx + m->yyyz*yyy - m->xxyy*(xxz + yyz) + m->xxyz*xxy + m->xyyz*xyy);
143 g4 = 0.25*(tx*x + ty*y + tz*z);
144 xxx = g3*(m->xxx*xx + m->xyy*yy + m->xxy*xy + m->xxz*xz + m->xyz*yz);
145 xxy = g3*(m->xyy*xy + m->xxy*xx + m->yyy*yy + m->yyz*yz + m->xyz*xz);
146 xxz = g3*(-(m->xxx + m->xyy)*xz - (m->xxy + m->yyy)*yz + m->xxz*xx + m->yyz*yy + m->xyz*xy);
147 g3 = onethird*(xxx*x + xxy*y + xxz*z);
148 xx = g2*(m->xx*x + m->xy*y + m->xz*z);
149 xy = g2*(m->yy*y + m->xy*x + m->yz*z);
150 xz = g2*(-(m->xx + m->yy)*z + m->xz*x + m->yz*y);
151 g2 = 0.5*(xx*x + xy*y + xz*z);
152 g0 *= m->m;
153 *fPot += -(g0 + g2 + g3 + g4);
154 g0 += 5*g2 + 7*g3 + 9*g4;
155 *ax += dir*(xx + xxx + tx - x*g0);
156 *ay += dir*(xy + xxy + ty - y*g0);
157 *az += dir*(xz + xxz + tz - z*g0);
158 *magai = g0*dir;
159 }
160#endif
161
163class MultipoleMoments {
167 cosmoType radius;
168public:
169 cosmoType soft; /* Effective softening */
170
172 cosmoType totalMass;
174 Vector3D<cosmoType> cm;
175
176#ifdef COOLING_MOLECULARH
177 double totalLW, totalgas;
178 Vector3D<double> cLW, cgas;
179 double xxgas, yygas, zzgas;
180#endif /*COOLING_MOLECULARH*/
181
182#ifdef HEXADECAPOLE
183 FMOMR mom;
184#else \
185 //Tensor for higher order moments goes here
186 double xx, xy, xz, yy, yz, zz;
187#endif
188
189 MultipoleMoments() : radius(0), totalMass(0) {
190 soft = 0;
191 cm.x = cm.y = cm.z = 0;
192#ifdef COOLING_MOLECULARH
193 totalLW = 0, totalgas = 0;
194 cLW.x = cLW.y = cLW.z = 0;
195 cgas.x = cgas.y = cgas.z = 0;
196 xxgas = yygas = zzgas = 0;
197#endif /*COOLING_MOLECULARH*/
198 //clear higher order components here
199#ifdef HEXADECAPOLE
200 momClearFmomr(&mom);
201#else
202 xx = xy = xz = yy = yz = zz = 0;
203#endif
204 }
205
207 MultipoleMoments& operator+=(const MultipoleMoments& m) {
208 //radius gets set by external function
209 cosmoType m1 = totalMass;
210 totalMass += m.totalMass;
211 if(totalMass == 0.0) {
212 soft = 0.5*(soft + m.soft);
213 cm = 0.5*(cm + m.cm);
214 return *this;
215 }
216 if(m1 == 0.0) { /* just copy over argument */
217 *this = m;
218 return *this;
219 }
220 if(m.totalMass == 0.0)
221 return *this;
222 soft = (m1*soft + m.totalMass*m.soft)/totalMass;
223 Vector3D<cosmoType> cm1 = cm;
224 cm = (m1*cm + m.totalMass*m.cm)/totalMass;
225#ifdef COOLING_MOLECULARH
226 double totalLW1 = totalLW;
227 if (m.totalLW != 0) { /*except this is log, so I'll have to figure out another way to mark this, CC*/
228 totalLW = log10(pow(10,totalLW - 30.0) + pow(10,m.totalLW - 30.0)) + 30.0; /*The thirty is just there to keep numbers in bounds*/
229 cLW = (pow(10,totalLW1 - 30.0)*cLW + pow(10,m.totalLW - 30.0)*m.cLW)/pow(10,totalLW - 30.0);
230 }
231
232 double totalgas1 = totalgas;
233 totalgas += m.totalgas;
234 Vector3D<double> cgas1 = cgas;
235 if (totalgas != 0) {
236 cgas = (totalgas1*cgas1 + m.totalgas*m.cgas)/totalgas;
237 }
238 Vector3D<double> drgas = cgas1 - cgas; /*Distance between previous center and new center*/
239 xxgas += totalgas1*drgas[0]*drgas[0];
240 yygas += totalgas1*drgas[1]*drgas[1];
241 zzgas += totalgas1*drgas[2]*drgas[2];
242 drgas = m.cgas - cgas; /*Distance between added node center and new center*/
243 xxgas += m.xxgas + m.totalgas*drgas[0]*drgas[0];
244 yygas += m.yygas + m.totalgas*drgas[1]*drgas[1];
245 zzgas += m.zzgas + m.totalgas*drgas[2]*drgas[2];
246#endif /*COOLING_MOLECULARH*/
247#ifdef HEXADECAPOLE
248 Vector3D<cosmoType> dr = cm1 - cm;
249 momShiftFmomr(&mom, radius, dr.x, dr.y, dr.z);
250 FMOMR mom2 = m.mom;
251 dr = m.cm - cm;
252 momShiftFmomr(&mom2, m.radius, dr.x, dr.y, dr.z);
253 momScaledAddFmomr(&mom, radius, &mom2, m.radius);
254#else
255 //add higher order components here
256 Vector3D<double> dr = cm1 - cm;
257 xx += m1*dr[0]*dr[0];
258 yy += m1*dr[1]*dr[1];
259 zz += m1*dr[2]*dr[2];
260 xy += m1*dr[0]*dr[1];
261 xz += m1*dr[0]*dr[2];
262 yz += m1*dr[1]*dr[2];
263 dr = m.cm - cm;
264 xx += m.xx + m.totalMass*dr[0]*dr[0];
265 yy += m.yy + m.totalMass*dr[1]*dr[1];
266 zz += m.zz + m.totalMass*dr[2]*dr[2];
267 xy += m.xy + m.totalMass*dr[0]*dr[1];
268 xz += m.xz + m.totalMass*dr[0]*dr[2];
269 yz += m.yz + m.totalMass*dr[1]*dr[2];
270#endif
271 return *this;
272 }
273
275 template <typename ParticleType>
276 MultipoleMoments& operator+=(const ParticleType& p) {
277 cosmoType m1 = totalMass;
278 totalMass += p.mass;
279 if(totalMass == 0.0) {
280 soft = 0.5*(soft + p.soft);
281 cm = 0.5*(cm + p.position);
282 return *this;
283 }
284 soft = (m1*soft + p.mass*p.soft)/totalMass;
285 Vector3D<cosmoType> cm1 = cm;
286 cm = (m1*cm + p.mass * p.position)/totalMass;
287#ifdef COOLING_MOLECULARH
288 double totalLW1 = totalLW;
289 if (p.isStar()) {
290 if (p.dStarLymanWerner() != 0) { /*except this is log, so I'll have to figure out another way to mark this*/
291 totalLW = log10(pow(10,totalLW - 30.0) + pow(10,p.dStarLymanWerner() - 30.0)) + 30.0; /*The thirty is just there to keep numbers in bounds*/
292 cLW = (pow(10,totalLW1 - 30.0)*cLW + pow(10,p.dStarLymanWerner() - 30.0)*p.position)/pow(10,totalLW - 30.0);
293 }
294 }
295
296 double totalgas1 = totalgas;
297 Vector3D<double> cgas1 = cgas;
298 if (p.isGas()) {
299 totalgas += p.mass;
300 if (totalgas!= 0) {
301 cgas = (totalgas1*cgas1 + p.mass * p.position)/totalgas;
302 }
303
304 Vector3D<double> drgas = cgas1 - cgas;
305
306 xxgas += totalgas1*drgas[0]*drgas[0];
307 yygas += totalgas1*drgas[1]*drgas[1];
308 zzgas += totalgas1*drgas[2]*drgas[2];
309
310 drgas = p.position - cgas;
311
312 xxgas += p.mass*drgas[0]*drgas[0];
313 yygas += p.mass*drgas[1]*drgas[1];
314 zzgas += p.mass*drgas[2]*drgas[2];
315 }
316
317#endif /*COOLING_MOLECULARH*/
318#ifdef HEXADECAPOLE
319 // XXX this isn't the most efficient way, but it
320 // retains the semantics of this function. It would
321 // be better to do this many particles at a time, then
322 // you could first determine the center of mass, then
323 // do a momMakeMomr(); momAddMomr() for each particle.
324 Vector3D<cosmoType> dr = cm1 - cm;
325 momShiftFmomr(&mom, radius, dr.x, dr.y, dr.z);
326 dr = p.position - cm;
327 FMOMR momPart;
328 momMakeFmomr(&momPart, p.mass, radius, dr.x, dr.y, dr.z);
329 momAddFmomr(&mom, &momPart);
330#else
331 //add higher order components here
332 Vector3D<double> dr = cm1 - cm;
333
334 xx += m1*dr[0]*dr[0];
335 yy += m1*dr[1]*dr[1];
336 zz += m1*dr[2]*dr[2];
337 xy += m1*dr[0]*dr[1];
338 xz += m1*dr[0]*dr[2];
339 yz += m1*dr[1]*dr[2];
340
341 dr = p.position - cm;
342
343 xx += p.mass*dr[0]*dr[0];
344 yy += p.mass*dr[1]*dr[1];
345 zz += p.mass*dr[2]*dr[2];
346 xy += p.mass*dr[0]*dr[1];
347 xz += p.mass*dr[0]*dr[2];
348 yz += p.mass*dr[1]*dr[2];
349#endif
350
351 return *this;
352 }
353
355 MultipoleMoments operator-(const MultipoleMoments& m) {
356 MultipoleMoments newMoments;
357 newMoments.totalMass = totalMass - m.totalMass;
358 if(newMoments.totalMass == 0.0) {
359 soft = 0.5*(soft - m.soft);
360 cm = 0.5*(cm - m.cm);
361 return *this;
362 }
363 newMoments.soft = (totalMass*soft - m.totalMass*m.soft)
364 /newMoments.totalMass;
365 newMoments.cm = (totalMass*cm - m.totalMass*m.cm)
366 /newMoments.totalMass;
367#ifdef HEXADECAPOLE
368 Vector3D<cosmoType> dr = cm - newMoments.cm;
369 newMoments.mom = mom;
370 momShiftFmomr(&mom, radius, dr.x, dr.y, dr.z);
371 FMOMR mom2 = m.mom;
372 dr = m.cm - newMoments.cm;
373 momShiftFmomr(&mom2, m.radius, dr.x, dr.y, dr.z);
374 momScaledSubFmomr(&newMoments.mom, radius, &mom2, m.radius);
375#else
376 //subtract off higher order components here
377 Vector3D<double> dr = cm - newMoments.cm;
378 newMoments.xx = xx + totalMass*dr[0]*dr[0];
379 newMoments.yy = yy + totalMass*dr[1]*dr[1];
380 newMoments.zz = zz + totalMass*dr[2]*dr[2];
381 newMoments.xy = xy + totalMass*dr[0]*dr[1];
382 newMoments.xz = xz + totalMass*dr[0]*dr[2];
383 newMoments.yz = yz + totalMass*dr[1]*dr[2];
384 dr = m.cm - newMoments.cm;
385 newMoments.xx -= m.xx + m.totalMass*dr[0]*dr[0];
386 newMoments.yy -= m.yy + m.totalMass*dr[1]*dr[1];
387 newMoments.zz -= m.zz + m.totalMass*dr[2]*dr[2];
388 newMoments.xy -= m.xy + m.totalMass*dr[0]*dr[1];
389 newMoments.xz -= m.xz + m.totalMass*dr[0]*dr[2];
390 newMoments.yz -= m.yz + m.totalMass*dr[1]*dr[2];
391#endif
392 return newMoments;
393 }
394
396 void clear() {
397 soft = 0;
398 radius = 0;
399 totalMass = 0;
400 cm.x = cm.y = cm.z = 0;
401 //clear higher order components here
402#ifdef HEXADECAPOLE
403 momClearFmomr(&mom);
404#else
405 xx = xy = xz = yy = yz = zz = 0;
406#endif
407 }
408 inline cosmoType getRadius() const {return radius;}
409 inline void setRadius(cosmoType r) {radius = r;}
410 friend void operator|(PUP::er& p, MultipoleMoments& m);
411 friend void calculateRadiusFarthestCorner(MultipoleMoments& m,
412 const OrientedBox<double>& box);
413 template<typename ParticleType>
414 friend void calculateRadiusFarthestParticle(MultipoleMoments& m,
415 const ParticleType * begin,
416 const ParticleType * end);
417 friend void calculateRadiusBox(MultipoleMoments& m,
418 const OrientedBox<double>& box);
419 template<typename ParticleType>
420 friend void calculateRadiusFirstParticle(MultipoleMoments& m,
421 const ParticleType* begin,
422 const ParticleType* end);
423};
424
425#ifdef __CHARMC__
426
427#include "pup.h"
428
429inline void operator|(PUP::er& p, MultipoleMoments& m) {
430 p | m.radius;
431 p | m.totalMass;
432 p | m.soft;
433 p | m.cm;
434#ifdef COOLING_MOLECULARH
435 p | m.totalLW;
436 p | m.totalgas;
437 p | m.cLW;
438 p | m.cgas;
439 p | m.xxgas;
440 p | m.yygas;
441 p | m.zzgas;
442#endif /*COOLING_MOLECULARH*/
443#ifdef HEXADECAPOLE
444 p((char *) &m.mom, sizeof(m.mom)); /* PUPs as bytes */
445#else
446 p | m.xx;
447 p | m.xy;
448 p | m.xz;
449 p | m.yy;
450 p | m.yz;
451 p | m.zz;
452#endif
453}
454
455#endif //__CHARMC__
456
457//What follows are criteria for deciding the size of a multipole
458
460inline void calculateRadiusFarthestCorner(MultipoleMoments& m, const OrientedBox<double>& box) {
461 Vector3D<cosmoType> delta1 = m.cm - box.lesser_corner;
462 Vector3D<cosmoType> delta2 = box.greater_corner - m.cm;
463 delta1.x = (delta1.x > delta2.x ? delta1.x : delta2.x);
464 delta1.y = (delta1.y > delta2.y ? delta1.y : delta2.y);
465 delta1.z = (delta1.z > delta2.z ? delta1.z : delta2.z);
466 cosmoType newradius = delta1.length();
467#ifdef HEXADECAPOLE
468 momRescaleFmomr(&m.mom, newradius, m.radius);
469#endif
470 m.radius = newradius;
471}
472
475inline void calculateRadiusBox(MultipoleMoments& m,
476 const OrientedBox<double>& box) {
477 Vector3D<cosmoType> delta = box.greater_corner - box.lesser_corner;
478 cosmoType newradius = 0.5*delta.length();
479#ifdef HEXADECAPOLE
480 if(m.totalMass > 0.0)
481 momRescaleFmomr(&m.mom, newradius, m.radius);
482#endif
483 m.radius = newradius;
484}
485
487template <typename ParticleType>
488inline void calculateRadiusFarthestParticle(MultipoleMoments& m, const ParticleType* begin, const ParticleType* end) {
489 Vector3D<cosmoType> cm = m.cm;
490 cosmoType d;
491 cosmoType newradius = 0;
492 for(const ParticleType* iter = begin; iter != end; ++iter) {
493 d = (cm - iter->position).lengthSquared();
494 if(d > newradius)
495 newradius = d;
496 }
497 if(newradius > 0.0) {
498 newradius = sqrt(newradius);
499#ifdef HEXADECAPOLE
500 momRescaleFmomr(&m.mom, newradius, m.radius);
501#endif
502 m.radius = newradius;
503 }
504}
505
509template<typename ParticleType>
510inline void calculateRadiusFirstParticle(MultipoleMoments& m,
511 const ParticleType* begin,
512 const ParticleType* end) {
513 const ParticleType* iter = begin;
514 Vector3D<cosmoType> pos1 = iter->position;
515 cosmoType newradius = 0.0;
516 iter++;
517 for(; iter != end; ++iter) {
518 cosmoType d = (pos1 - iter->position).lengthSquared();
519 if(d > newradius)
520 newradius = d;
521 }
522 if(newradius > 0.0) {
523 newradius = sqrt(newradius);
524 m.radius = newradius;
525 }
526}
527#endif //MULTIPOLEMOMENTS_H
void calculateRadiusFarthestCorner(MultipoleMoments &m, const OrientedBox< double > &box)
Given an enclosing box, set the multipole expansion size to the distance from the center of mass to t...
Definition MultipoleMoments.h:460
void calculateRadiusBox(MultipoleMoments &m, const OrientedBox< double > &box)
Definition MultipoleMoments.h:475
void calculateRadiusFirstParticle(MultipoleMoments &m, const ParticleType *begin, const ParticleType *end)
Definition MultipoleMoments.h:510
void calculateRadiusFarthestParticle(MultipoleMoments &m, const ParticleType *begin, const ParticleType *end)
Given the positions that make up a multipole expansion, set the distance to the farthest particle fro...
Definition MultipoleMoments.h:488
A representation of a multipole expansion.
Definition MultipoleMoments.h:163
Vector3D< cosmoType > cm
The center of mass (zeroth order multipole)
Definition MultipoleMoments.h:174
friend void calculateRadiusFarthestCorner(MultipoleMoments &m, const OrientedBox< double > &box)
Given an enclosing box, set the multipole expansion size to the distance from the center of mass to t...
Definition MultipoleMoments.h:460
MultipoleMoments operator-(const MultipoleMoments &m)
Subtract an expansion from this larger one, yielding the leftover.
Definition MultipoleMoments.h:355
friend class CudaMultipoleMoments
Version of MultipoleMoments using cudatype.
Definition MultipoleMoments.h:164
MultipoleMoments & operator+=(const ParticleType &p)
Add the contribution of a particle to this multipole expansion.
Definition MultipoleMoments.h:276
friend void calculateRadiusBox(MultipoleMoments &m, const OrientedBox< double > &box)
Definition MultipoleMoments.h:475
cosmoType totalMass
The total mass represented by this expansion.
Definition MultipoleMoments.h:172
void clear()
Reset this expansion to nothing.
Definition MultipoleMoments.h:396
friend void calculateRadiusFirstParticle(MultipoleMoments &m, const ParticleType *begin, const ParticleType *end)
Definition MultipoleMoments.h:510
MultipoleMoments & operator+=(const MultipoleMoments &m)
Add two expansions together, using parallel axis theorem.
Definition MultipoleMoments.h:207
friend void calculateRadiusFarthestParticle(MultipoleMoments &m, const ParticleType *begin, const ParticleType *end)
Given the positions that make up a multipole expansion, set the distance to the farthest particle fro...
Definition MultipoleMoments.h:488