389 lines
17 KiB
TeX
389 lines
17 KiB
TeX
\documentclass[10pt]{beamer}
|
|
|
|
\input{mystyle_beamer.tex}
|
|
|
|
\usetikzlibrary{arrows.meta}
|
|
\usetikzlibrary{decorations.pathmorphing}
|
|
\usetikzlibrary{calc}
|
|
|
|
\usetikzlibrary{external}
|
|
\tikzexternalize[prefix=tikzcache/]
|
|
|
|
\usepackage{appendixnumberbeamer}
|
|
|
|
\tikzset{zigzag/.style={decorate, decoration=zigzag}}
|
|
\def \L {2.}
|
|
|
|
|
|
\title{Towards Moving Picture Simulations in a Discontinuous Galerkin Framework}
|
|
\subtitle{A Progress Report}
|
|
\author{Yingjie Wang}
|
|
\institute{FAU}
|
|
|
|
\addbibresource{ref.bib}
|
|
|
|
\begin{document}
|
|
\maketitle
|
|
|
|
\begin{frame}{Contents}
|
|
\tableofcontents
|
|
\end{frame}
|
|
|
|
\begin{frame}{Outline of the Project}
|
|
\section{Introduction}
|
|
\begin{enumerate}
|
|
\item A first order evolution system, to be named FOZ4c
|
|
\item An implementation in our discontinuous Galerkin code \texttt{nmesh}
|
|
\item Simulations of binary black holes with large mass ratio with the moving puncture gauge
|
|
\end{enumerate}
|
|
\end{frame}
|
|
|
|
\begin{frame}{Final Goal: binary black holes with large mass ratio}
|
|
\subsection{binary black holes with large mass ratio}
|
|
|
|
Why binary black holes?
|
|
|
|
\begin{itemize}
|
|
\item Main sources of gravitational waves
|
|
\item Simplest two-body problem in general relativity
|
|
\item Connecting strong field GR with post-Newtonian approximations
|
|
\end{itemize}
|
|
|
|
Why large mass ratio?
|
|
\begin{itemize}
|
|
\item LISA is able to detect extreme mass ratio inspirals (EMRIs)
|
|
\item Compare with perturbation theory
|
|
\end{itemize}
|
|
\end{frame}
|
|
|
|
\begin{frame}{Methods on evolving black holes}
|
|
\subsection{Methods on evolving black holes}
|
|
There are two main methods to evolve black holes in numerical relativity:
|
|
\begin{itemize}
|
|
\item Excision method: excise the black hole interior from the computational domain, and impose boundary conditions on the excision surface
|
|
\item Moving puncture method: evolve the black hole as a puncture, and use a suitable gauge condition to avoid the singularity.
|
|
\end{itemize}
|
|
\end{frame}
|
|
|
|
\begin{frame}{Methods on evolving black holes: excision}
|
|
\begin{figure}
|
|
\centering
|
|
\begin{tikzpicture}[>=Latex, line cap=round, line join=round]
|
|
% causal diamond
|
|
\draw[thick,red,zigzag] (-\L,\L) coordinate(stl) -- (\L,\L) coordinate (str);
|
|
\draw[thick,black] (\L,-\L) coordinate (sbr)
|
|
-- (0,0) coordinate (bif) -- (stl);
|
|
\draw[thick,black,fill=blue, fill opacity=0.2,text opacity=1]
|
|
(bif) -- (str) -- (2*\L,0) node[right] (io) {$i^0$} -- (sbr);
|
|
|
|
% null labels
|
|
\draw[black] (1.4*\L,0.7*\L) node[right] (scrip) {$\mathcal{I}^+$}
|
|
(1.5*\L,-0.6*\L) node[right] (scrip) {$\mathcal{I}^-$}
|
|
(0.2*\L,-0.6*\L) node[right] (scrip) {$\mathcal{H}^-$}
|
|
(0.5*\L,0.85*\L) node[right] (scrip) {$\mathcal{H}^+$};
|
|
|
|
% singularity label
|
|
\draw[thick,red,<-] (0,1.05*\L)
|
|
-- (0,1.2*\L) node[above] {\color{red} singularity};
|
|
% % Scwharzschild surface
|
|
% \draw[thick,blue] (bif) .. controls (1.*\L,-0.35*\L) .. (2*\L,0);
|
|
% \draw[thick,blue,<-] (1.75*\L,-0.1*\L) -- (1.9*\L,-0.5*\L)
|
|
% -- (2*\L,-0.5*\L) node[right,align=left]
|
|
% {$t=$ constant\\in Schwarzschild\\coordinates};
|
|
% excision surface
|
|
\draw[thick,dashed,red] (-0.3*\L,0.3*\L) -- (0.4*\L,\L);
|
|
\draw[thick,red,<-] (-0.33*\L,0.3*\L)
|
|
-- (-0.5*\L,0.26*\L) node[left,align=right] {excision\\surface};
|
|
% Kerr-Schild surface
|
|
\draw[green,thick] (0.325*\L,0.325*\L) .. controls (\L,0) .. (2*\L,0);
|
|
\draw[green,dashed,thick] (0.325*\L,0.325*\L) -- (-0.051*\L,0.5*\L);
|
|
% Kerr-Schild label
|
|
\draw[green,thick,<-] (0.95*\L,0.15*\L) -- (1.2*\L,0.5*\L)
|
|
-- (2*\L,0.5*\L) node[right,align=left]
|
|
{time slice};
|
|
\end{tikzpicture}
|
|
\caption{excision on a single Schwarchild black hole}
|
|
\end{figure}
|
|
\end{frame}
|
|
|
|
\begin{frame}{Methods on evolving black holes: excision}
|
|
\begin{figure}
|
|
\centering
|
|
\includegraphics[width=0.7\textwidth]{imgs/black_hole_excision_mesh.png}
|
|
\caption{An example mesh using excision method from \cite{Hemberger_2013}.}
|
|
\end{figure}
|
|
\end{frame}
|
|
|
|
|
|
\begin{frame}{Methods on evolving black holes: moving puncture}
|
|
In the isotropic coordinate for Schwarzschild black hole, the spical metric
|
|
\begin{equation*}
|
|
\dd{l^2} = \left( 1+\frac{M}{2r} \right)^4 (\dd{r^2} + r^2 \dd{\Omega^2})
|
|
\end{equation*}
|
|
is conformally flat, with the conformal factor
|
|
\begin{equation*}
|
|
\psi = 1 + \frac{M}{2r}.
|
|
\end{equation*}
|
|
|
|
General $N$ black hole puncture initial data:
|
|
\begin{equation*}
|
|
\psi = \sum_{i=1}^N \frac{M_i}{2|\myvec{r}-\myvec{r}_i|} + (\text{regular part}),
|
|
\end{equation*}
|
|
again, if we factor out the conformal factor $\psi$, the remaining conformal metric is regular everywhere.
|
|
\end{frame}
|
|
|
|
\begin{frame}{Methods on evolving black holes: moving puncture}
|
|
|
|
In the moving puncture mothod, we won't cut out the black hole singularity, but choose a gauge condition that makes the singularity invisible to the numerical evolution.
|
|
\begin{figure}
|
|
\centering
|
|
\includegraphics[width=0.7\textwidth, trim=4bp 4bp 4bp 4bp, clip]{imgs/moving_puncture_penrose_1.png}
|
|
\caption{The numerical time slices in the moving puncture method in the Penrose diagram, from \cite{Hannam_2008}.}
|
|
\end{figure}
|
|
\end{frame}
|
|
|
|
\begin{frame}{Methods on evolving black holes: moving puncture}
|
|
\begin{figure}
|
|
\centering
|
|
\includegraphics[width=0.42\textwidth, trim=5bp 4bp 4bp 4bp, clip]{imgs/moving_puncture_init.png}
|
|
\includegraphics[width=0.45\textwidth, trim=4bp 4bp 4bp 4bp, clip]{imgs/moving_puncture_stable.png}
|
|
\caption{The embedding of the initial slice (left) and the final slice (right) in the moving puncture method, from \cite{Hannam_2008}.}
|
|
\end{figure}
|
|
\end{frame}
|
|
|
|
\begin{frame}{Discontinuous Galerkin method and \texttt{nmesh}}
|
|
\subsection{Discontinuous Galerkin method and \texttt{nmesh}}
|
|
There are two main numerical methods in numerical relativity:
|
|
\begin{itemize}
|
|
\item Finite Difference / Finite Volume
|
|
\begin{itemize}
|
|
\item Easy to implement, and mature codes available
|
|
\end{itemize}
|
|
\item Spectral / Discontinuous Galerkin
|
|
\begin{itemize}
|
|
\item High accuracy for smooth solutions
|
|
\end{itemize}
|
|
\end{itemize}
|
|
\end{frame}
|
|
|
|
\begin{frame}{Discontinuous Galerkin method and \texttt{nmesh}}
|
|
Core idea of DG method: if we have a first-order PDE system:
|
|
\begin{equation*}
|
|
\pdvt{u} + A^i \pdv{u}{x^i} = S,
|
|
\end{equation*}
|
|
then for some test functions $\{v_a\}$, we have
|
|
\begin{equation*}
|
|
\left(v_a, \pdvt{u}\right) + \left(v_a, A^i \pdv{u}{x^i}\right) = (v_a, S),
|
|
\end{equation*}
|
|
which reduce to a linear system of $\{ u_a = (v_a, u)\}$ after intergrating by parts. Boundary terms are replaced by numerical fluxes.
|
|
|
|
In \texttt{nmesh}, we use Lagrange polynomials over Gauss-Legendre points on each element.
|
|
\begin{center}
|
|
\begin{tikzpicture}[scale=4]
|
|
\draw[thick, blue] (0,0) -- (1,0);
|
|
\draw[thick, red] (1,0) -- (2,0);
|
|
|
|
\pgfmathsetmacro{\a}{0.5 - 0.5*sqrt(3/7)}
|
|
\pgfmathsetmacro{\b}{0.5 + 0.5*sqrt(3/7)}
|
|
\pgfmathsetmacro{\c}{1.5 - 0.5*sqrt(3/7)}
|
|
\pgfmathsetmacro{\d}{1.5 + 0.5*sqrt(3/7)}
|
|
|
|
\foreach \x in {0,\a,0.5,\b,1}
|
|
\fill[blue] (\x,0) circle (0.018);
|
|
|
|
\foreach \x in {1,\c,1.5,\d,2}
|
|
\fill[red] (\x,0) circle (0.018);
|
|
|
|
\fill[blue] (1,0) circle (0.022);
|
|
\fill[red] (1,0) circle (0.013);
|
|
\end{tikzpicture}
|
|
\end{center}
|
|
Information exchange only happens at the boundaries of elements through numerical fluxes, which makes it efficient for parallel computing.
|
|
\end{frame}
|
|
|
|
\begin{frame}{Discontinuous Galerkin method and \texttt{nmesh}}
|
|
\texttt{nmesh} has two features that are useful for our project:
|
|
\begin{itemize}
|
|
\item Adaptive mesh refinement (AMR)
|
|
\item DG + FD/FV dynamically switching
|
|
\end{itemize}
|
|
\begin{figure}
|
|
\centering
|
|
\includegraphics[width=0.65\textwidth]{imgs/MPA1_W-9sn12l5_GRHD_D_t0400.pdf}
|
|
\caption{The mesh in a binary star emulation with \texttt{nmesh}}
|
|
\end{figure}
|
|
\end{frame}
|
|
|
|
\begin{frame}{Moving puncture evolution in \texttt{nmesh}?}
|
|
\subsection{Moving puncture evolution in \texttt{nmesh}?}
|
|
There are many $3+1$ formulations of Einstein's equations,
|
|
\begin{itemize}
|
|
\item We want DG method:
|
|
\begin{itemize}
|
|
\item a {\color{blue}first order} formulation is needed, like the generalized harmonic (GH) formulation.
|
|
\end{itemize}
|
|
\item We want moving puncture evolution:
|
|
\begin{itemize}
|
|
\item the {\color{blue}moving puncture gauge condition} is needed, and not compatible with the GH formulation.
|
|
\item People usually use BSSN, Z4c or CCZ4 formulation for moving puncture evolutions, but they are second order in space, and not directly suitable for DG method.
|
|
\end{itemize}
|
|
\item We want stable evolution:
|
|
\begin{itemize}
|
|
\item a {\color{blue}strongly hyperbolic} formulation is needed.
|
|
\item {\color{blue} constraint damping} is needed.
|
|
\end{itemize}
|
|
\end{itemize}
|
|
$\implies$ we need to find a first order, strongly hyperbolic formulation of Einstein's equations with constraint damping and compatible with the moving puncture gauge condition.
|
|
\end{frame}
|
|
|
|
\begin{frame}{Searching for a first order Z4c}{Hyperbolicity of first order systems}
|
|
\section{Searching for a first order Z4c}
|
|
\subsection{Hyperbolicity of first order systems}
|
|
Say we have a list of variables $u^I$, and a first order PDE system
|
|
\begin{equation}
|
|
\pdvt{u^I} + \tensor{A}{^i^I_J} \pdv{u^J}{x^i} = S^I, \label{first-order-PDE}
|
|
\end{equation}
|
|
where both $A$ and $S$ can depend on $u^I$ but not their derivatives. For any unit covector $\xi_i$, the \emph{Principal symbol} matrix is defined as $\tensor{P}{^I_J}(\xi) := \tensor{A}{^i^I_J} \xi_i$.
|
|
The system is said to be:
|
|
\begin{itemize}
|
|
\item \emph{weakly hyperbolic} if for any $\xi_i$, $\tensor{P}{^I_J}(\xi)$ has only real eigenvalues.
|
|
\item \emph{strongly hyperbolic} if for any $\xi_i$, $\tensor{P}{^I_J}(\xi)$ has only real eigenvalues and a complete set of eigenvectors, i.e., it is diagonalizable.
|
|
\item \emph{symmetrically hyperbolic} if there exists a positive definite symmetrizer matrix $H_{IJ}$ such that $H_{IK} \tensor{P}{^K_J}(\xi)$ is always symmetric.
|
|
\end{itemize}
|
|
|
|
\begin{theorem}
|
|
System~\eqref{first-order-PDE} is well-posed in the $L^2$ sense if and only if it is strongly hyperbolic.
|
|
\end{theorem}
|
|
\end{frame}
|
|
|
|
\begin{frame}{Searching for a first order Z4c}{First order reduction}
|
|
\subsection{First order reduction: wave equation as an example}
|
|
We can introduce auxiliary variables to reduce a second order PDE system to a first order one. For example, for the wave equation on flat spacetime:
|
|
\begin{equation}
|
|
\Partial{t}^2 {\phi} - \tensor{\delta}{^i^j} \Partial{i} \Partial{j} \phi = 0,
|
|
\end{equation}
|
|
we can introduce $\pi := \Partial{t} \phi$ and $\tensor{\psi}{_i} := \Partial{_i} \phi$, and rewrite the wave equation as
|
|
\begin{equation}
|
|
\begin{cases}
|
|
\Partial{t}{\phi} = \pi, \\
|
|
\Partial{t}{\pi} = \tensor{\delta}{^i^j} \Partial{i} \tensor{\psi}{_j},\\
|
|
\Partial{t}{\tensor{\psi}{_i}} = \Partial{i} \pi.
|
|
\end{cases}
|
|
\end{equation}
|
|
write $u^I = (\phi, \pi, \tensor{\psi}{_1}, \tensor{\psi}{_2}, \tensor{\psi}{_3})$, then the system can be written in the form of~\eqref{first-order-PDE} with
|
|
\begin{equation}
|
|
\tensor{A}{^i^I_J} = \begin{pmatrix}
|
|
0 & 0 & 0 \\
|
|
0 & 0 & -\tensor{e}{_i} \\
|
|
0 & -\tensor{e}{^i} & 0
|
|
\end{pmatrix} \qc
|
|
S=0,
|
|
\end{equation}
|
|
where $\tensor{e}{^i}$ is the standard column vector basis of $\setR^3$, while $\tensor{e}{_i}$ is the standard row vector basis of $(\setR^3)^*$.
|
|
\end{frame}
|
|
|
|
\begin{frame}{Searching for a first order Z4c}{First order reduction}
|
|
The principal symbol matrix is
|
|
\begin{equation}
|
|
\tensor{P}{^I_J}(\xi) = \begin{pmatrix}
|
|
0 & 0 & 0 \\
|
|
0 & 0 & -\myvec{\xi}^T \\
|
|
0 & -\myvec{\xi} & 0
|
|
\end{pmatrix},
|
|
\end{equation}
|
|
it's very clear that the principal symbol matrix is already symmetric, and thus the system is symmetrically hyperbolic, and well-posed in the $L^2$ sense.
|
|
|
|
The eigenvalues are $\lambda_1 = \lambda_2 = \lambda_3 = 0$, $\lambda_4 = -1$ and $\lambda_5 = 1$. The eigenvectors are
|
|
\begin{equation}
|
|
\begin{gathered}
|
|
\myvec{e}_1 = \mqty(1,0,0,0,0)^T \qc
|
|
\myvec{e}_2 = \mqty(0,0,-\xi_3,0,\xi_1)^T \qc
|
|
\myvec{e}_3 = \mqty(0,0,-\xi_2,\xi_1,0)^T,\\
|
|
\myvec{e}_4 = \mqty(0,1,\xi_1,\xi_2,\xi_3)^T \qc
|
|
\myvec{e}_5 = \mqty(0,-1,\xi_1,\xi_2,\xi_3)^T.
|
|
\end{gathered}
|
|
\end{equation}
|
|
\end{frame}
|
|
|
|
\begin{frame}{Searching for a first order Z4c}{constraints during the reduction}
|
|
During the first order reduction, we introduce new variables. A solution to the new system is a solution to the original system if and only if the new variables satisfy some constraints. For example, in the wave equation case, we have the constraints
|
|
\begin{equation}
|
|
\mathcal{C}_i \definedby \tensor{\psi}{_i} - \Partial{i} \phi = 0.
|
|
\end{equation}
|
|
The evolution of the constraints is given by
|
|
\begin{equation}
|
|
\Partial{t} \mathcal{C}_i = \Partial{t} \psi_i - \Partial{i} \Partial{t} \phi = \Partial{i} \pi - \Partial{i} \pi = 0,
|
|
\end{equation}
|
|
which means that if the constraints are satisfied initially, they will be satisfied for all time. However, in numerical simulations, there will always be some constraint violation due to numerical errors, and thus it's better to add some constraint damping terms to the evolution equations to suppress the growth of constraint violation.
|
|
\end{frame}
|
|
|
|
\begin{frame}{Searching for a first order Z4c}{constraints during the reduction}
|
|
For example, we can add a constraint damping term $-\gamma \mathcal{C}_i$ to the evolution equation of $\tensor{\psi}{_i}$,
|
|
\begin{equation}
|
|
\begin{cases}
|
|
\Partial{t}{\phi} = \pi, \\
|
|
\Partial{t}{\pi} = \tensor{\delta}{^i^j} \Partial{i} \tensor{\psi}{_j},\\
|
|
\Partial{t}{\tensor{\psi}{_i}} = \Partial{i} \pi - \gamma \left( \tensor{\psi}{_i} - \Partial{i} \phi \right),
|
|
\end{cases}
|
|
\end{equation}
|
|
where $\gamma > 0$ is a constant. Then the evolution of the constraints becomes
|
|
\begin{equation}
|
|
\Partial{t} \mathcal{C}_i = -\gamma \mathcal{C}_i,
|
|
\end{equation}
|
|
which means that the constraint violation will decay exponentially with time, and thus the system is more stable for numerical simulations.
|
|
\end{frame}
|
|
|
|
\begin{frame}{Searching for a first order Z4c}{constraints during the reduction}
|
|
We have to check if the constraint damping term will change the hyperbolicity of the system. The new principal symbol matrix is
|
|
\begin{equation}
|
|
\tensor{P}{^I_J}(\xi) = \begin{pmatrix}
|
|
0 & 0 & 0 \\
|
|
0 & 0 & -\myvec{\xi}^T \\
|
|
-\gamma \myvec{\xi} & -\myvec{\xi} & 0
|
|
\end{pmatrix},
|
|
\end{equation}
|
|
which is no longer symmetric, but luckily it's still diagonalizable with exactly the same eigenvalues, and the eigenvectors are also almost the same, except that
|
|
\begin{equation}
|
|
\myvec{e}_1 = \mqty(1,-\gamma,0,0,0)^T,
|
|
\end{equation}
|
|
and thus the system is still strongly hyperbolic, and well-posed in the $L^2$ sense.
|
|
\end{frame}
|
|
|
|
\begin{frame}{Existing first order formulations}
|
|
\subsection{Existing first order formulations}
|
|
\begin{itemize}
|
|
\item GH
|
|
\begin{itemize}
|
|
\item Not compatible with moving puncture mothod
|
|
\end{itemize}
|
|
\item FOCCZ4
|
|
\begin{itemize}
|
|
\item No constraint damping
|
|
\end{itemize}
|
|
\item FOZ4
|
|
\begin{itemize}
|
|
\item No constraint damping
|
|
\end{itemize}
|
|
\item FOBSSN
|
|
\begin{itemize}
|
|
\item Probably works; not well tested in the community; WIP in \texttt{nmesh}
|
|
\end{itemize}
|
|
\end{itemize}
|
|
\end{frame}
|
|
|
|
\begin{frame}{References}
|
|
\printbibliography
|
|
\end{frame}
|
|
|
|
\appendix
|
|
|
|
\begin{frame}{Appendix}
|
|
\section{Appendix}
|
|
\begin{figure}
|
|
\centering
|
|
\includegraphics[width=0.7\textwidth, trim=4bp 4bp 4bp 4bp, clip]{imgs/moving_puncture_penrose.png}
|
|
\caption{The numerical time slices in the moving puncture method in the Penrose diagram, from \cite{Hannam_2008}.}
|
|
\end{figure}
|
|
\end{frame}
|
|
|
|
\end{document} |