\begin{frame}\frametitle{Appendix I: Gravitational Softening}

\textbf{The Motivation:} Most of the computational effort in N-body
simulations is spent on binary stars.\vspace{2mm}


\textbf{Solution:} Prevent the formation of binary star systems.

\vspace{2mm}

\begin{block}
<1->{}
\begin{itemize}
\item Replace point-particles with extended {}``fuzzy'' particles.
\item Particles now represent small \emph{groups} of stars rather than single
stars.
\item For short distances, replace Newtonian force with another expression
corresponding to extended objects.
\end{itemize}
\end{block}

\end{frame}

%----------------------------


\begin{frame}\frametitle{Appendix I: Gravitational Softening}

\textbf{More Details: }\vspace{2mm}


Point particles are ``smeared out'' into spheres of radius $\epsilon$.
The parameter $\epsilon$ is known as the \underbar{softening length}.

\vspace{4mm}

\underbar{Examples of Softening Methods}:\vspace{5mm}

\begin{block}
<1->{Plummer Softening}
\begin{center}Replace $\;\frac{Gm_{i}m_{j}\mathbf{r}_{ij}}{r_{ij}^{3}}\;$
with $\;\frac{Gm_{i}m_{j}\mathbf{r}_{ij}}{\left(r_{ij}^{2}+\epsilon^{2}\right)^{3/2}}$\end{center}
\end{block}

\end{frame}

%----------------------------


\begin{frame}\frametitle{Appendix I:  Gravitational Softening}
\begin{block}
<1->{Plummer Softening Spline}
\begin{center} $\mathbf{F}_{ij}=\begin{cases}
-\frac{Gm_{i}m_{j}\mathbf{r}_{ij}}{r_{ij}^{3}} & \;\;\mbox{for }r>2\epsilon\\
-\frac{Gm_{i}m_{j}\mathbf{r}_{ij}}{\left(r_{ij}^{2}+\epsilon^{2}\right)^{3/2}} & \;\;\mbox{for }r\leq2\epsilon\end{cases}$ \end{center}
\end{block}
\vspace{6mm}

\begin{itemize}
\item Other methods are also utilized, often incorporating Gaussians and
splines.
\end{itemize}

\end{frame}

%----------------------------


\begin{frame}\frametitle{Appendix II: Leap-Frog Method}

Using the centered difference, compute acceleration at time $t$,\[
\mathbf{a}\left(\mathbf{r}\left(t\right)\right)=\frac{\mathbf{v}\left(t+\delta t\right)-\mathbf{v}\left(t-\delta t\right)}{2\delta t}+\mathcal{O}\left(\delta t^{2}\right)\]
and the velocity at time $t+\delta t$, \[
\mathbf{v}\left(t+\delta t\right)=\frac{\mathbf{r}\left(t+2\delta t\right)-\mathbf{r}\left(t\right)}{2\delta t}+\mathcal{O}\left(\delta t^{2}\right)\]
Then rearrange terms so that future values are on the left\[
\mathbf{v}\left(t+\delta t\right)=\mathbf{v}\left(t-\delta t\right)+2\mathbf{a}\left(\mathbf{r}\left(t\right)\right)\delta t+\mathcal{O}\left(\delta t^{3}\right)\]
\[
\mathbf{r}\left(t+2\delta t\right)=\mathbf{r}\left(t\right)+2\mathbf{v}\left(t+\delta t\right)\delta t+\mathcal{O}\left(\delta t^{3}\right)\]



\end{frame}

%----------------------------


\begin{frame}\frametitle{Appendix II: Leap-Frog Method}

Using index notation, \[
\mathbf{v}_{i+1}=\mathbf{v}_{i-1}+2\mathbf{a}_{i}\delta t+\mathcal{O}\left(\delta t^{3}\right)\]
\[
\mathbf{r}_{i+2}=\mathbf{r}_{i}+2\mathbf{v}_{i+1}\delta t+\mathcal{O}\left(\delta t^{3}\right)\]
This is sometimes written as\[
\mathbf{v}_{i+\frac{1}{2}}=\mathbf{v}_{i-\frac{1}{2}}+\mathbf{a}_{i}\delta t+\mathcal{O}\left(\delta t^{3}\right)\]
\[
\mathbf{r}_{i+1}=\mathbf{r}_{i}+\mathbf{v}_{i+\frac{1}{2}}\delta t+\mathcal{O}\left(\delta t^{3}\right)\]
and other ways...

\end{frame}

%----------------------------


\begin{frame}\frametitle{Appendix II: Leap-Frog Method}

\vspace{4mm}
\begin{columns}[T]   
 \begin{column}{5cm}\textbf{drift-kick-drift version:}\vspace{1mm}
\textrm{$\mathbf{r}_{i+\frac{1}{2}}=\mathbf{r}_{i}+\mathbf{a}_{i}\frac{\delta t}{2}$}\vspace{1mm}


\textrm{$\mathbf{v}_{i+1}=\mathbf{v}_{i}+\mathbf{a}_{i+\frac{1}{2}}\delta t$}\vspace{1mm}


\textrm{$\mathbf{r}_{i+1}=\mathbf{r}_{i+\frac{1}{2}}+\mathbf{v}_{i+1}\frac{\delta t}{2}$}\vspace{1mm}


where

\textrm{$\mathbf{\;\;\;\;\; a}_{i+\frac{1}{2}}=\mathbf{a}\left(\mathbf{r}_{i+\frac{1}{2}}\right)$}

\end{column}   
 \begin{column}{5cm}     

\textbf{kick-drift-kick version:}

\vspace{3mm}


\textrm{$\mathbf{v}_{i+\frac{1}{2}}=\mathbf{v}_{i}+\mathbf{a}_{i}\frac{\delta t}{2}$}\vspace{1mm}


\textrm{$\mathbf{r}_{i+1}=\mathbf{r}_{i}+\mathbf{v}_{i+\frac{1}{2}}\frac{\delta t}{2}$}\vspace{1mm}


\textrm{$\mathbf{v}_{i+1}=\mathbf{v}_{i+\frac{1}{2}}+\mathbf{a}_{i+1}\frac{\delta t}{2}$}\vspace{1mm}


where 

\textrm{$\;\;\;\;\;\mathbf{a}_{i+1}=\mathbf{a}\left(\mathbf{r}_{i+1}\right)$}

\end{column} 
\end{columns} 


\end{frame}

%----------------------------


\begin{frame}\frametitle{Appendix II: Leap-Frog Method}

\textbf{Comments:}
\begin{itemize}
\item Even though it is only second order, Leap-Frog yields better results
than RK4
\item Leap-frog method conserves energy (it is symplectic\ldots but only for a fixed step size)
\item Kick-drift-kick is the most stable 
\end{itemize}

\end{frame}

%----------------------------


\begin{frame}\frametitle{Appendix III: Smoothed-Particle Hydrodynamics (SPH)}

\begin{columns}[T]   
 \begin{column}{5cm}

\textbf{Basic principle:}
\begin{itemize}
\item Fluid quantities are sampled at discrete points with fixed mass.
\item Particles move with the bulk flow of the fluid.
\item Particles are \emph{smoothed} over spherical regions of radius $h$,
called {}``smoothing length''.
\item Particles are required to overlap. 
\end{itemize}
\end{column}   
 \begin{column}{5cm}     \includegraphics[width=5cm]{images/sph-overlap.jpg}
\begin{itemize}
\item Specify minimum number of required overlaps, then scale $h$.
\end{itemize}
\end{column} 
\end{columns} 


\end{frame}

%----------------------------

\begin{frame}\frametitle{Appendix III: SPH, Basic Formulation}

Start with the trivial identity:\[
A(\mathbf{r})=\int A(\mathbf{r}')\delta\left(\mathbf{r}-\mathbf{r}'\right)d\mathbf{r}'\]
This can be approximated\[
\left\langle A(\mathbf{r})\right\rangle \approx\int A(\mathbf{r}')w\left(\mathbf{r}-\mathbf{r}',h\right)d\mathbf{r}'\]
where the smoothing kernel \textrm{$w(r,h)$} obeys \[
\int w(r,h)dV=1\]
and\[
\lim_{h\rightarrow0}w(r,h)\;\longrightarrow\;\delta\left(r\right)\]


\end{frame}

%----------------------------


\begin{frame}\frametitle{Appendix III: SPH, Basic Formulation}


\[
\left\langle A(\mathbf{r})\right\rangle \approx\int A(\mathbf{r}')w\left(\mathbf{r}-\mathbf{r}',h\right)d\mathbf{r}'\]


\[
=\int\frac{A(\mathbf{r}')}{\rho\left(\mathbf{r}'\right)}\rho\left(\mathbf{r}'\right)w\left(\mathbf{r}-\mathbf{r}',h\right)d\mathbf{r}'\]
This is approximated as\[
\left\langle A(\mathbf{r}_{i})\right\rangle \approx\sum_{j}\left[\rho\left(\mathbf{r}{}_{j}\right)V_{j}\right]\frac{A\left(\mathbf{r}{}_{j}\right)}{\rho\left(\mathbf{r}{}_{j}\right)}w\left(\mathbf{r}_{i}-\mathbf{r}{}_{j},h\right)\]
\[
=\sum_{j}m_{j}\frac{A_{j}}{\rho_{j}}w\left(\mathbf{r}_{i}-\mathbf{r}_{j},h\right)\]


\end{frame}

%----------------------------


\begin{frame}\frametitle{Appendix III: SPH, Basic Formulation}

We arrive at the approximation\[
\left\langle A(\mathbf{r}_{i})\right\rangle \approx\sum_{j}m_{j}\frac{A_{j}}{\rho_{j}}w\left(\mathbf{r}_{i}-\mathbf{r}_{j},h\right)\]
this is analogous to FEM expansion in terms of shape functions.\vspace{3mm}

\begin{block}
<1->{Example: calculate $\left\langle \rho(\mathbf{r}_{i})\right\rangle $}\textrm{\[
\left\langle \rho(\mathbf{r}_{i})\right\rangle \approx\sum_{j}m_{j}\frac{\rho_{j}}{\rho_{j}}w\left(\mathbf{r}_{i}-\mathbf{r}_{j},h\right)=\sum_{j}m_{j}w\left(\mathbf{r}_{i}-\mathbf{r}_{j},h\right)\]
}
\end{block}

\end{frame}

%----------------------------



\begin{frame}\frametitle{Appendix III: SPH, Basic Formulation}

\begin{block}
<1->{Example: calculate $\left\langle \mathbf{v}_{i}\right\rangle  $}\textrm{\[
\left\langle \mathbf{v}_{i}\right\rangle \approx\sum_{j}\frac{m_{j}}{\rho_{j}}\mathbf{v}_{j}w\left(\mathbf{r}_{i}-\mathbf{r}_{j},h\right)\]
}
\end{block}
Gradient and divergence are approximated by returning to integral
formulation, using integration by parts and vector identities to manipulate
the expressions. The results are:

\[
\nabla A(\mathbf{r}_{i})\approx\left\langle \nabla A(\mathbf{r}_{i})\right\rangle \approx\sum_{j}m_{j}\frac{A_{j}}{\rho_{j}}\nabla_{i}w\left(\mathbf{r}_{i}-\mathbf{r}_{j},h\right)
\]
\[
\nabla\cdot\mathbf{v}_{i}\approx\frac{1}{\rho_{i}}\sum_{j}\frac{m_{j}}{\rho_{j}}\left(\mathbf{v}_{i}-\mathbf{v}_{j}\right)\cdot\nabla_{i}w\left(\mathbf{r}_{i}-\mathbf{r}_{j},h\right)
\]


\end{frame}

%----------------------------

\begin{frame}\frametitle{Appendix IV: Bondi-Hoyle-Lyttleton Accretion}


The Bondi-Hoyle-Lyttleton Accretion rate:

\begin{center}

\[\dot{M}=\frac{4\pi\alpha G^{2}M^{2}\rho}{\left(c_{s}^{2}+v_{\infty}^{2}\right)^{3/2}}\]

\begin{itemize}
 \item $\alpha$ is a resolution mis-match correction term.  If the Bondi radius is resolved, $\alpha=1$
 \item $c_s$ is the local speed of sound.
\item $v_{\infty}$ is the relative speed of the surrounding gas and the black hole. 
\item $\rho$ is the density of the local gas.
\end{itemize}

\end{center}

\end{frame}

%----------------------

\begin{frame}\frametitle{Appendix V: The Jeans Equation} 

For completely collisionless, dissipationless systems, the evolution is described by the \textbf{collisionless Boltzmann equation} 

$$\partfrac{f}{t}+\mathbf{v}\cdot\nabla f -\nabla \Phi \cdot \partfrac{f}{\mathbf{v}}=0$$

Where $f = f(\mathbf{x},\mathbf{v},t)$ is the phase space density.  

\pause

\vspace{0.3cm}

In component form,

\begin{equation}
\partfrac{f}{t}+ \sum_{i=1}^3\left( v_i \partfrac{f}{x_i} - \partfrac{\Phi}{x_i}\partfrac{f}{v_i}\right)=0\label{eq:vlasov}
\end{equation}

\end{frame}

%----------------------

\begin{frame}\frametitle{Appendix V: The Jeans Equation}

$$\partfrac{f}{t}+ \sum_{i=1}^3\left( v_i \partfrac{f}{x_i} - \partfrac{\Phi}{x_i}\partfrac{f}{v_i}\right)=0$$

\vspace{0.25cm}

Integrating over all velocities and using the summation convention,

$$\int \partfrac{f}{t}d^3\mathbf{v} + \int v_i \partfrac{f}{x_i} d^3\mathbf{v}-\partfrac{\Phi}{x_i}\int\partfrac{f}{v_i}d^3\mathbf{v}=0$$

\pause

\vspace{0.25cm}

Note that $v_i$ does not depend on $x_i$ and the integration domain (all velocities) is independent of $t$. Then\ldots 

\pause

$$\partfrac{}{t} \int f d^3\mathbf{v} + \partfrac{}{x_i} \int v_i f d^3 \mathbf{v}-\partfrac{\Phi}{x_i}\int\partfrac{f}{v_i}d^3\mathbf{v}=0$$

\vspace{0.25cm}

 
\end{frame}




%----------------------

\begin{frame}\frametitle{Appendix V: The Jeans Equation}

$$\partfrac{}{t} \int f d^3\mathbf{v} + \partfrac{}{x_i} \int v_i f d^3 \mathbf{v}-\partfrac{\Phi}{x_i}{\color{red}{\int\partfrac{f}{v_i}d^3\mathbf{v}}}=0$$

\vspace{0.25cm}

The red term vanishes upon using the divergence theorem and noting that stars have zero probability of having infinite speed. 

\pause

\vspace{0.25cm}

Thus $$\partfrac{}{t} \int f d^3\mathbf{v} + \partfrac{}{x_i} \int v_i f d^3 \mathbf{v}=0$$


\end{frame}

%----------------------

\begin{frame}\frametitle{Appendix V: The Jeans Equation}


$$\partfrac{}{t} \int f d^3\mathbf{v} + \partfrac{}{x_i} \int v_i f d^3 \mathbf{v}=0$$

\vspace{0.25cm}

Defining the stellar number density $\nu(\mathbf{x})$ and the mean stellar velocity $\bar{\mathbf{v}}(\mathbf{x})$: 

$$\nu\equiv\int f d^3\mathbf{v}\qquad\mbox{and}\qquad \bar{v}_i\equiv\frac{1}{\nu}\int fv_i d^3\mathbf{v}$$

\pause

The equation becomes

\begin{equation}
 \partfrac{\nu}{t} + \partfrac{(\nu\bar{v}_i)}{x_i}=0\label{eq:conserve}
\end{equation}



This is a statement of the conservation of particle number. 


\end{frame}


%------------------------


\begin{frame}\frametitle{Appendix V: The Jeans Equation}


Multiplying the collisionless Boltzmann equation (\ref{eq:vlasov}) by $v_i$ and repeating the process of integrating over all velocities, yields

\vspace{0.25cm}

$$\partfrac{(\nu\bar{v}_j)}{t} + \partfrac{(\nu\overline{v_i v_j})}{x_i} + \partfrac{\Phi}{x_i}\nu\delta_{ij}=0$$

where
$$\overline{v_i v_j}\equiv \frac{1}{\nu}\int v_i v_j f  d^3\mathbf{v}\mbox{.}$$ \pause

Simplifying,

\begin{equation}
 \partfrac{(\nu\bar{v}_j)}{t} + \partfrac{(\nu\overline{v_i v_j})}{x_i} + \nu \partfrac{\Phi}{x_j}=0\label{eq:vj}
\end{equation}

\end{frame}


%----------------------

\begin{frame}\frametitle{Appendix V: The Jeans Equation}

In summary:
 
\begin{eqnarray*}
\partfrac{\nu}{t} + \partfrac{(\nu\bar{v}_i)}{x_i} &=& 0\qquad\qquad\qquad\mbox{(\ref{eq:conserve})}\\
\partfrac{(\nu\bar{v}_j)}{t} + \partfrac{(\nu\overline{v_i v_j})}{x_i} + \nu \partfrac{\Phi}{x_j} &=& 0 \qquad\qquad\qquad\mbox{(\ref{eq:vj})}
\end{eqnarray*}

\pause
Multiplying eq.(\ref{eq:conserve}) by $\bar{v}_j$, subtracting this from eq.(\ref{eq:vj}), and simplifying 
\pause
\begin{equation}
 \nu\partfrac{\bar{v}_j}{t} + \partfrac{(\nu\overline{v_i v_j})}{x_i} + \nu \partfrac{\Phi}{x_j} -  \bar{v}_j\partfrac{(\nu\bar{v}_i)}{x_i}=0\label{eq:jeans1}
\end{equation}

\end{frame}

%----------------------

\begin{frame}\frametitle{Appendix V: The Jeans Equation}

\begin{equation*}
 \nu\partfrac{\bar{v}_j}{t} + \partfrac{(\nu\overline{v_i v_j})}{x_i} + \nu \partfrac{\Phi}{x_j} -  \bar{v}_j\partfrac{(\nu\bar{v}_i)}{x_i}=0
\end{equation*}

\vspace{0.2cm}

The term $\overline{v_i v_j}$ can be written as the sum of streaming and random parts: 
\begin{equation}
  \overline{v_i v_j}=\bar{v}_i\bar{v}_j+\sigma^2_{ij}\label{eq:barij}
\end{equation}

where $$\sigma^2_{ij}\equiv\overline{(v_i-\bar{v}_i)(v_j-\bar{v}_j)}$$


\vspace{0.25cm}

$\sigma^2_{ij}$ is the covariance of $(v_i,v_j)$, AKA the velocity dispersion tensor.

\pause
\vspace{0.25cm}

Now we substitute Eq.(\ref{eq:barij}) into our equation\ldots

\end{frame}


%----------------------

\begin{frame}\frametitle{Appendix V: The Jeans Equation}

\begin{equation*}
\nu\partfrac{\bar{v}_{j}}{t}+\partfrac{(\nu\bar{v}_i\bar{v}_{j})}{x_i}+\partfrac{(\nu\sigma_{ij}^{2})}{x_i}+\nu\partfrac{\Phi}{x_{j}}-\bar{v}_{j}\partfrac{(\nu\bar{v}_i)}{x_i}=0
\end{equation*}

\pause

Simplifying,

\[
\nu\partfrac{\bar{v}_{j}}{t}+\nu\bar{v}_i\partfrac{\bar{v}_{j}}{x_i}+\partfrac{(\nu\sigma_{ij}^{2})}{x_i}+\nu\partfrac{\Phi}{x_{j}}=0
\]

\pause

Upon rearranging, this becomes

\begin{equation}
\nu\partfrac{\bar{v}_{j}}{t}+\nu\bar{v}_i\partfrac{\bar{v}_{j}}{x_i}=-\nu\partfrac{\Phi}{x_{j}}-\partfrac{(\nu\sigma_{ij}^{2})}{x_i}\label{eq:jeans}
\end{equation}

This is known as the \textbf{Jeans equation}.

\end{frame}


