Added equations for the field tensors.

examples/theory/multipoles.tex

@@ -9,7 +9,10 @@
 \newcommand{\p}[1]{\mathbf{p}_{#1}}
 \newcommand{\acc}[1]{\mathbf{a}_{#1}}
 \newcommand{\muu}{\boldsymbol{\mu}}
-
+\newcommand{\qq}[1]{\mathbf{Q}_{#1}}
+\newcommand{\jj}[1]{\underline{J_{#1}}}
+\newcommand{\ii}[1]{\underline{I_{#1}}}
+\newcommand{\identity}{\rm{Id_3}}
 \title{B-H and FMM equations up to quadrupole terms}
 \author{Matthieu Schaller}
 
@@ -18,7 +21,8 @@
 \maketitle
 
 Bold quantities are vectors, underlined quantities are matrices. The indices 
-$\alpha,\beta$ run over the directions $x,y,z$.
+$\alpha,\beta$ run over the directions $x,y,z$. $\identity$ is the identity 
+matrix in 3D.
 
 \section{Construction of multipoles}
 
@@ -37,33 +41,33 @@ centre of mass $\muu_A$. The first LHS term uses Dehnen's notation.\\
 
 Monopole:
 \begin{equation}
-M_{(0,0,0)} = M_{{\rm tot}, A} 
+M_{(0,0,0)}(\muu_A) = M_{{\rm tot}, A} 
 \end{equation}
 
 Dipole:
 \begin{eqnarray}
-M_{(1,0,0)} &=& P_{A,x}~= 0\\
-M_{(0,1,0)} &=& P_{A,y}~=0\\
-M_{(0,0,1)} &=& P_{A,z}~=0
+M_{(1,0,0)}(\muu_A) &=& P_{A,x}~= 0\\
+M_{(0,1,0)}(\muu_A) &=& P_{A,y}~=0\\
+M_{(0,0,1)}(\muu_A) &=& P_{A,z}~=0
 \end{eqnarray}
 
 Quadrupole:
 \begin{eqnarray}
-M_{(2,0,0)} &=& \frac{1}{2}I_{A,xx}~= \frac{1}{2}\sum_{i\in A}m_i ( 
+M_{(2,0,0)}(\muu_A) &=& \frac{1}{2}I_{A,xx}~= \frac{1}{2}\sum_{i\in A}m_i ( 
 p_{i,x}-\mu_{A,x})^2\\
-M_{(0,2,0)} &=& \frac{1}{2}I_{A,yy}~= \frac{1}{2}\sum_{i\in A}m_i ( 
+M_{(0,2,0)}(\muu_A) &=& \frac{1}{2}I_{A,yy}~= \frac{1}{2}\sum_{i\in A}m_i ( 
 p_{i,y}-\mu_{A,y})^2\\
-M_{(0,0,2)} &=& \frac{1}{2}I_{A,zz}~= \frac{1}{2}\sum_{i\in A}m_i ( 
+M_{(0,0,2)}(\muu_A) &=& \frac{1}{2}I_{A,zz}~= \frac{1}{2}\sum_{i\in A}m_i ( 
 p_{i,z}-\mu_{A,z})^2\\
-M_{(1,1,0)} &=& \frac{1}{2}I_{A,xy}~= \frac{1}{2}\sum_{i\in A}m_i ( 
+M_{(1,1,0)}(\muu_A) &=& \frac{1}{2}I_{A,xy}~= \frac{1}{2}\sum_{i\in A}m_i ( 
 p_{i,x}-\mu_{A,x})( p_{i,y}-\mu_{A,y})\\
-M_{(0,1,1)} &=& \frac{1}{2}I_{A,yz}~= \frac{1}{2}\sum_{i\in A}m_i ( 
+M_{(0,1,1)}(\muu_A) &=& \frac{1}{2}I_{A,yz}~= \frac{1}{2}\sum_{i\in A}m_i ( 
 p_{i,y}-\mu_{A,y})( p_{i,z}-\mu_{A,z})\\
-M_{(1,0,1)} &=& \frac{1}{2}I_{A,xz}~= \frac{1}{2}\sum_{i\in A}m_i ( 
+M_{(1,0,1)}(\muu_A) &=& \frac{1}{2}I_{A,xz}~= \frac{1}{2}\sum_{i\in A}m_i ( 
 p_{i,x}-\mu_{A,x})( p_{i,z}-\mu_{A,z})
 \end{eqnarray}
 
-\section{Recursive construction of the quadrupoles}
+\section{Recursive construction of the multipoles}
 
 Given a set of multipoles $B$ expressed around their centre of masses 
 $\muu_{B}$, we can construct the total multipoles around the centre of mass
@@ -127,7 +131,7 @@ In the B-H approximation, the potential at position $\p{i}$ due to a set of
 particles $A$ is given by
 
 \begin{equation}
- \phi(\p{i}) = -\sum_{\bf n} M_{A, \bf n} D_{\bf n}(\muu_A - \p{i})
+ \phi(\p{i}) = -\sum_{\bf n} M_{\bf n}(\muu_A) D_{\bf n}(\muu_A - \p{i})
 \end{equation}
 
 Keeping only the terms up to second order (i.e. letting the sum run over 
@@ -146,20 +150,101 @@ I_{A,\alpha} \left(\frac{3Gr_\alpha^2}{|\rr_{i,A}|^5} -
  & & - \frac{1}{2}\sum_{\alpha,\beta} I_{A,\alpha\beta}  \frac{3Gr_\alpha 
 r_\beta}{|\rr_{i,A}|^5} \\
 &=& -G \left[ \frac{M_{{\rm tot},A}}{|\rr_{i,A}|} - \frac{1}{2} \frac{{\rm 
-tr}(\underline{I_A})}{|\rr_{i,A}|^3} + \frac{3}{2}\frac{\rr_{i,A}^T \cdot 
-\underline{I_A} \cdot \rr_{i,A}}{|\rr_{i,A}|^5}\right]
+tr}(\ii{A})}{|\rr_{i,A}|^3} + \frac{3}{2}\frac{\rr_{i,A}^T \cdot 
+\ii{A} \cdot \rr_{i,A}}{|\rr_{i,A}|^5}\right]
 \end{eqnarray}
 
 The accelerations $\acc{i} = -\nabla_{\rr_{i,A}}\phi(\p{i})$ are then given by:
 
 \begin{eqnarray}
 \acc{i} = G\left[\frac{M_{{\rm tot},A} \rr_{i,A}}{|\rr_{i,A}|} 
--\frac{3}{2}\frac{{\rm tr}(\underline{I_A})\rr_{i,A}}{|\rr_{i,A}|^5}
--\frac{3\underline{I_A}\cdot\rr_{i,A}}{|\rr_{i,A}|^5}
-+\frac{15}{2}\frac{(\rr_{i,A}^T \cdot \underline{I_A} \cdot 
+-\frac{3}{2}\frac{{\rm tr}(\ii{A})\rr_{i,A}}{|\rr_{i,A}|^5}
+-\frac{3\ii{A}\cdot\rr_{i,A}}{|\rr_{i,A}|^5}
++\frac{15}{2}\frac{(\rr_{i,A}^T \cdot \ii{A} \cdot 
 \rr_{i,A})\rr_{i,A}}{|\rr_{i,A}|^7}\right]
 \end{eqnarray}
 
+The last two expressions are used in the Quickshed example code as well as in 
+Bonsai. Gadget only uses the first term in each expression, avoiding the 
+construction of the matrices $\ii{A}$. In practice, the 2-nd order 
+accurate B-H method requires the storage of 10 variables per cell ($M_{{\rm 
+tot},A}, \muu_A, \ii{A}$).
+
+\section{FMM Field tensors}
+
+Instead of computing the potential and accelerations of each particle, field 
+tensors (generated by the set of particles A) at position $\muu_B$ can 
+be derived. The accelerations of each particle can then be obtained from 
+these field tensors. These tensors are given by:
+
+\begin{equation}
+ F_{\bf n}(\muu_B) = \sum_{|{\bf m}| + |{\bf n}|\leq2} M_{\bf m}(\muu_A)D_{{\bf 
+n}+{\bf m}}(\muu_B - \muu_A)
+\end{equation}
+
+Writing these out explicitly with $\rr_{BA} = \muu_B - \muu_A$, we obtain the 
+following set of expressions.\\
+
+Monopole:
+\begin{eqnarray}
+ F_{(0,0,0)}(\muu_B) ~=~ N_B &=& \sum_{|{\bf m}| \leq 2} M_{\bf 
+m}(\muu_A)D_{{\bf m}}(\rr_{BA}) \\
+&=&  M_{(0,0,0)}(\muu_A)D_{(0,0,0)}(\rr_{BA})\\ 
+& & +M_{(2,0,0)}(\muu_A)D_{(2,0,0)}(\rr_{BA}) \\
+& & +M_{(1,1,0)}(\muu_A)D_{(1,1,0)}(\rr_{BA})\\
+& & + \dots \\
+&=&M_{{\rm tot}, A}\phi(\rr_{BA}) + 
+\frac{1}{2}\sum_{\alpha,\beta}I_{A,\alpha\beta}\partial_{\alpha\beta}\phi(\rr_{
+BA }) \\
+&=& \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|} -\frac{1}{2}\frac{{\rm 
+tr}(\ii{A})}{|\rr_{BA}|^3} + 
+\frac{3}{2}\frac{\rr_{BA}^T \cdot 
+\ii{A} \cdot \rr_{BA}}{|\rr_{BA}|^5}
+\end{eqnarray}
+
+Dipole:
+\begin{eqnarray}
+  F_{(1,0,0)}(\muu_B) ~=~ Q_{B,x} &=& \sum_{|{\bf m}| \leq 1} M_{\bf 
+m}(\muu_A)D_{{\bf 
+m}+(1,0,0)}(\rr_{BA}) \\
+&=&  M_{(0,0,0)}(\muu_A)D_{(1,0,0)}(\rr_{BA})\\ 
+&=& M_{{\rm tot}, A} \partial_x \phi(\rr_{BA}) \\
+&=& M_{{\rm tot}, A} \frac{G}{|\rr_{BA}|^3}r_{BA,x}
+\end{eqnarray}
+
+Quadrupole (diagonal term):
+\begin{eqnarray}
+ F_{(2,0,0)}(\muu_B)~=~ J_{B,xx} &=& \sum_{|{\bf m}| \leq 0} M_{\bf 
+m}(\muu_A)D_{{\bf 
+m}+(2,0,0)}(\rr_{BA}) \\
+&=& M_{(0,0,0)}(\muu_A)D_{(2,0,0)}(\rr_{BA}) \\
+&=& M_{{\rm tot}, A} \partial_{xx} \phi(\rr_{BA}) \\
+&=& M_{{\rm tot}, A} \left(\frac{3Gr_{BA,x}^2}{|\rr_{BA}|^5} 
+- \frac{G}{|\rr_{BA}|^3}\right)
+\end{eqnarray}
+
+
+Quadrupole (off-diagonal term):
+\begin{eqnarray}
+ F_{(1,1,0)}(\muu_B)~=~J_{B,xy} &=& \sum_{|{\bf m}| \leq 0} M_{\bf 
+m}(\muu_A)D_{{\bf 
+m}+(1,1,0)}(\rr_{BA}) \\
+&=& M_{(0,0,0)}(\muu_A)D_{(1,1,0)}(\rr_{BA}) \\
+&=& M_{{\rm tot}, A} \partial_{xy} \phi(\rr_{BA}) \\
+&=& M_{{\rm tot}, A}\frac{3Gr_{BA,x}r_{BA,y}}{|\rr|^5}
+\end{eqnarray}
+
+All these terms can be written using a more compact notation.
+\begin{eqnarray}
+N_B &=&  \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|} -\frac{1}{2}\frac{{\rm 
+tr}(\ii{A})}{|\rr_{BA}|^3} + 
+\frac{3}{2}\frac{\rr_{BA}^T \cdot 
+\ii{A} \cdot \rr_{BA}}{|\rr_{BA}|^5} \\
+\qq{B} &=&  \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|^3} \rr_{BA}\\
+\jj{B} &=& \frac{3GM_{{\rm tot}, A}}{|\rr_{BA}|^5} \rr_{BA} 
+\otimes\rr_{BA} - \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|^3} \identity
+\end{eqnarray}
+