Skip to content
Snippets Groups Projects
Commit 6e5380ba authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Corrected acceleration equations in the BH case.

parent be7a8787
No related branches found
No related tags found
No related merge requests found
...@@ -106,9 +106,9 @@ notation.\\ ...@@ -106,9 +106,9 @@ notation.\\
1-st order: 1-st order:
\begin{eqnarray} \begin{eqnarray}
D_{(1,0,0)}(\rr) &=& \partial_x \phi(\rr)~=\frac{G}{|\rr|^3}r_x\\ D_{(1,0,0)}(\rr) &=& \partial_x \phi(\rr)~=-\frac{G}{|\rr|^3}r_x\\
D_{(0,1,0)}(\rr) &=& \partial_y \phi(\rr)~=\frac{G}{|\rr|^3}r_y\\ D_{(0,1,0)}(\rr) &=& \partial_y \phi(\rr)~=-\frac{G}{|\rr|^3}r_y\\
D_{(0,0,1)}(\rr) &=& \partial_z \phi(\rr)~=\frac{G}{|\rr|^3}r_z D_{(0,0,1)}(\rr) &=& \partial_z \phi(\rr)~=-\frac{G}{|\rr|^3}r_z
\end{eqnarray} \end{eqnarray}
2-nd order: 2-nd order:
...@@ -124,6 +124,23 @@ D_{(0,1,1)} &=& \partial_{yz} \phi(\rr)~= \frac{3Gr_yr_z}{|\rr|^5}\\ ...@@ -124,6 +124,23 @@ D_{(0,1,1)} &=& \partial_{yz} \phi(\rr)~= \frac{3Gr_yr_z}{|\rr|^5}\\
D_{(1,0,1)} &=& \partial_{xz} \phi(\rr)~= \frac{3Gr_xr_z}{|\rr|^5} D_{(1,0,1)} &=& \partial_{xz} \phi(\rr)~= \frac{3Gr_xr_z}{|\rr|^5}
\end{eqnarray} \end{eqnarray}
3-rd order:
\begin{eqnarray}
D_{(3,0,0)} &=& \partial_{xxx} \phi(\rr)~= -\frac{15Gr_x^3}{|\rr|^7} +
\frac{9Gr_x}{|\rr|^5}\\
D_{(1,2,0)} &=& \partial_{xyy} \phi(\rr)~= -\frac{15Gr_xr_y^2}{|\rr|^7} +
\frac{3Gr_x}{|\rr|^5}\\
D_{(1,0,2)} &=& \partial_{xzz} \phi(\rr)~= -\frac{15Gr_xr_z^2}{|\rr|^7} +
\frac{3Gr_x}{|\rr|^5}\\
D_{(2,1,0)} &=& \partial_{xxy} \phi(\rr)~= -\frac{15Gr_x^2r_y}{|\rr|^7} +
\frac{3Gr_y}{|\rr|^5}\\
D_{(2,0,1)} &=& \partial_{xxz} \phi(\rr)~= -\frac{15Gr_x^2r_z}{|\rr|^7} +
\frac{3Gr_z}{|\rr|^5}\\
D_{(1,1,1)} &=& \partial_{xyz} \phi(\rr)~= -\frac{15Gr_xr_yr_z}{|\rr|^7}
\end{eqnarray}
\section{B-H potential and accelerations} \section{B-H potential and accelerations}
...@@ -168,9 +185,21 @@ The last two expressions are used in the Quickshed example code as well as in ...@@ -168,9 +185,21 @@ 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 Bonsai. Gadget only uses the first term in each expression, avoiding the
construction of the matrices $\ii{A}$. In practice, the 2-nd order 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 accurate B-H method requires the storage of 10 variables per cell ($M_{{\rm
tot},A}, \muu_A, \ii{A}$). tot},A}, \muu_A, \ii{A}$).\\
\section{FMM Field tensors} The accelerations can also written using Dehnen's short notation:
\begin{eqnarray}
a_{i,x} &=& \sum_{\bf n} M_{\bf n}(\muu_A) D_{{\bf n}+(1,0,0)}(\muu_A -
\p{i}) \\
a_{i,y} &=& \sum_{\bf n} M_{\bf n}(\muu_A) D_{{\bf n}+(0,1,0)}(\muu_A -
\p{i}) \\
a_{i,z} &=& \sum_{\bf n} M_{\bf n}(\muu_A) D_{{\bf n}+(0,0,1)}(\muu_A -
\p{i})
\end{eqnarray}
\section{FMM Field tensors for potential}
Instead of computing the potential and accelerations of each particle, field Instead of computing the potential and accelerations of each particle, field
tensors (generated by the set of particles A) at position $\muu_B$ can tensors (generated by the set of particles A) at position $\muu_B$ can
...@@ -187,7 +216,7 @@ following set of expressions.\\ ...@@ -187,7 +216,7 @@ following set of expressions.\\
Monopole: Monopole:
\begin{eqnarray} \begin{eqnarray}
F_{(0,0,0)}(\muu_B) ~=~ N_B &=& \sum_{|{\bf m}| \leq 2} M_{\bf F_{(0,0,0)}(\muu_B) ~=~ N_{BA} &=& \sum_{|{\bf m}| \leq 2} M_{\bf
m}(\muu_A)D_{{\bf m}}(\rr_{BA}) \\ m}(\muu_A)D_{{\bf m}}(\rr_{BA}) \\
&=& M_{(0,0,0)}(\muu_A)D_{(0,0,0)}(\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_{(2,0,0)}(\muu_A)D_{(2,0,0)}(\rr_{BA}) \\
...@@ -196,15 +225,15 @@ m}(\muu_A)D_{{\bf m}}(\rr_{BA}) \\ ...@@ -196,15 +225,15 @@ m}(\muu_A)D_{{\bf m}}(\rr_{BA}) \\
&=&M_{{\rm tot}, A}\phi(\rr_{BA}) + &=&M_{{\rm tot}, A}\phi(\rr_{BA}) +
\frac{1}{2}\sum_{\alpha,\beta}I_{A,\alpha\beta}\partial_{\alpha\beta}\phi(\rr_{ \frac{1}{2}\sum_{\alpha,\beta}I_{A,\alpha\beta}\partial_{\alpha\beta}\phi(\rr_{
BA }) \\ BA }) \\
&=& \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|} -\frac{1}{2}\frac{{\rm &=& \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|} -\frac{G}{2}\frac{{\rm
tr}(\ii{A})}{|\rr_{BA}|^3} + tr}(\ii{A})}{|\rr_{BA}|^3} +
\frac{3}{2}\frac{\rr_{BA}^T \cdot \frac{3G}{2}\frac{\rr_{BA}^T \cdot
\ii{A} \cdot \rr_{BA}}{|\rr_{BA}|^5} \ii{A} \cdot \rr_{BA}}{|\rr_{BA}|^5}
\end{eqnarray} \end{eqnarray}
Dipole: Dipole:
\begin{eqnarray} \begin{eqnarray}
F_{(1,0,0)}(\muu_B) ~=~ Q_{B,x} &=& \sum_{|{\bf m}| \leq 1} M_{\bf F_{(1,0,0)}(\muu_B) ~=~ Q_{BA,x} &=& \sum_{|{\bf m}| \leq 1} M_{\bf
m}(\muu_A)D_{{\bf m}(\muu_A)D_{{\bf
m}+(1,0,0)}(\rr_{BA}) \\ m}+(1,0,0)}(\rr_{BA}) \\
&=& M_{(0,0,0)}(\muu_A)D_{(1,0,0)}(\rr_{BA})\\ &=& M_{(0,0,0)}(\muu_A)D_{(1,0,0)}(\rr_{BA})\\
...@@ -214,7 +243,7 @@ m}+(1,0,0)}(\rr_{BA}) \\ ...@@ -214,7 +243,7 @@ m}+(1,0,0)}(\rr_{BA}) \\
Quadrupole (diagonal term): Quadrupole (diagonal term):
\begin{eqnarray} \begin{eqnarray}
F_{(2,0,0)}(\muu_B)~=~ J_{B,xx} &=& \sum_{|{\bf m}| \leq 0} M_{\bf F_{(2,0,0)}(\muu_B)~=~ J_{BA,xx} &=& \sum_{|{\bf m}| \leq 0} M_{\bf
m}(\muu_A)D_{{\bf m}(\muu_A)D_{{\bf
m}+(2,0,0)}(\rr_{BA}) \\ m}+(2,0,0)}(\rr_{BA}) \\
&=& M_{(0,0,0)}(\muu_A)D_{(2,0,0)}(\rr_{BA}) \\ &=& M_{(0,0,0)}(\muu_A)D_{(2,0,0)}(\rr_{BA}) \\
...@@ -226,7 +255,7 @@ m}+(2,0,0)}(\rr_{BA}) \\ ...@@ -226,7 +255,7 @@ m}+(2,0,0)}(\rr_{BA}) \\
Quadrupole (off-diagonal term): Quadrupole (off-diagonal term):
\begin{eqnarray} \begin{eqnarray}
F_{(1,1,0)}(\muu_B)~=~J_{B,xy} &=& \sum_{|{\bf m}| \leq 0} M_{\bf F_{(1,1,0)}(\muu_B)~=~J_{BA,xy} &=& \sum_{|{\bf m}| \leq 0} M_{\bf
m}(\muu_A)D_{{\bf m}(\muu_A)D_{{\bf
m}+(1,1,0)}(\rr_{BA}) \\ m}+(1,1,0)}(\rr_{BA}) \\
&=& M_{(0,0,0)}(\muu_A)D_{(1,1,0)}(\rr_{BA}) \\ &=& M_{(0,0,0)}(\muu_A)D_{(1,1,0)}(\rr_{BA}) \\
...@@ -236,16 +265,48 @@ m}+(1,1,0)}(\rr_{BA}) \\ ...@@ -236,16 +265,48 @@ m}+(1,1,0)}(\rr_{BA}) \\
All these terms can be written using a more compact notation. All these terms can be written using a more compact notation.
\begin{eqnarray} \begin{eqnarray}
N_B &=& \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|} -\frac{1}{2}\frac{{\rm N_{BA} &=& \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|} -\frac{G}{2}\frac{{\rm
tr}(\ii{A})}{|\rr_{BA}|^3} + tr}(\ii{A})}{|\rr_{BA}|^3} +
\frac{3}{2}\frac{\rr_{BA}^T \cdot \frac{3G}{2}\frac{\rr_{BA}^T \cdot
\ii{A} \cdot \rr_{BA}}{|\rr_{BA}|^5} \\ \ii{A} \cdot \rr_{BA}}{|\rr_{BA}|^5} \\
\qq{B} &=& \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|^3} \rr_{BA}\\ \qq{BA} &=& \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|^3} \rr_{BA}\\
\jj{B} &=& \frac{3GM_{{\rm tot}, A}}{|\rr_{BA}|^5} \rr_{BA} \jj{BA} &=& \frac{3GM_{{\rm tot}, A}}{|\rr_{BA}|^5} \rr_{BA}
\otimes\rr_{BA} - \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|^3} \identity \otimes\rr_{BA} - \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|^3} \identity
\end{eqnarray} \end{eqnarray}
Note that $\jj{B}$ is symmetric, meaning that $10$ additional numbers per cell
only have to be stored. That is a total of $10$ for the multipoles and $10$ for
the field tensors.
\section{Gravitational potential in FMM}
The potential at the position $\p{i}$ in the vicinity of $\muu_B$ due to the
mass clustered around $\muu_A$ can be computed using the field tensors created
by the mutipoles at position $\muu_A$ on the position $\muu_B$.
\begin{equation}
\phi_A(\p{i}) = \sum_{|{\bf n}| \leq2} \frac{1}{ {\bf n}!} (\p{i} -
\muu_B)^{\bf n} F_{\bf n}(\muu_B)
\end{equation}
Setting $\dd_{i,B}= \p{i} - \muu_B$ and $\rr_{BA} = \muu_B - \muu_A$, this sum
expands into
\begin{eqnarray}
\phi_A(\p{i}) &=& N_{BA} + \dd_{i,B}\cdot\qq{BA} + \dd_{i,B}^T \cdot \jj{BA}
\cdot \dd_{i,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} \\
& &+ \frac{GM_{{\rm tot}, A}}{|\rr_{BA}|^3} \rr_{BA}\cdot\dd_{i,B} \\
& &+\frac{3}{2}\frac{GM_{{\rm tot}, A}}{|\rr_{BA}|^5} \left(\dd_{i,B}^T \cdot
\rr_{BA}\otimes\rr_{BA}\cdot \dd_{i,B}\right) \\
& & - \frac{1}{2}\frac{GM_{{\rm tot}, A}}{|\rr_{BA}|^3} \left(\dd_{i,B}^T
\cdot\identity \cdot \dd_{i,B}\right)
\end{eqnarray}
Note that if $\p{i}=\muu_B$ (i.e. $d_{i,B}={\bf 0}$), the expression for the
potential in the B-H formalism is recovered.
\end{document} \end{document}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment