Numerical Relativity

Running Binary Black Hole

AthenaK requires external modules to provide Cauchy initial date for the binary black hole problems. Currently, we can read in initial data from TWOPUNCTURE and SpECTRE code. The TWOPUNCTURE initial data is calculated on the fly at the start of the evolution, while SpECTRE initial data needs to be pre-computated. Below is a tutorial for evolving the TWOPUNCTURE initial data.

TWOPUNCTURE Initial Data

Two external libraries are required for building the z4c_two_puncture problem, i.e. gsl and twopuncturesc initial data solver (in this order).

cd $HOME && mkdir -p usr/gsl && mkdir codes && cd codes
# grab source
wget ftp://ftp.gnu.org/gnu/gsl/gsl-2.5.tar.gz

# extract and configure for local install
tar -zxvf gsl-2.5.tar.gz

cd gsl-2.5

./configure --prefix=$HOME/usr/gsl

make -j8
# make check
make install

# link gsl into athenak
ln -s $HOME/usr/gsl ${path_to_athenak}

twopuncturesc:

cd $HOME && cd usr
git clone git@bitbucket.org:bernuzzi/twopuncturesc.git
cd twopuncturesc
make -j8

# link twopuncturesc into athenak
ln -s $HOME/usr/twopuncturesc ${path_to_athenak}

Now create build directory and configure with cmake

mkdir build_z4c_twopunc && cd build_z4c_twopunc

Finally, we can configure and build Athenak

$ cmake ../ -DPROBLEM=z4c_two_puncture

(TO DO (HZ): add tutorial for SpECTRE ID)

Running GRMHD problems in Dynamical Spacetimes

The Valencia solver can be enabled by starting with a standard GRMHD parfile and adding an <adm> or <z4c> block. There are few new <mhd> parameters available for the Valencia solver:

mhd/dyn_eos: the PrimitiveSolver EOS (currently supports ideal for an ideal gas, piecewise_poly for piecewise polytropes, and eos_compose for a tabulated EOS in the PyCompOSE format). Note that mhd/eos must still be defined as ideal no matter what dyn_eos is; this inconvenience is for historical reasons and will hopefully be patched in a future update.
mhd/dyn_error: the error policy used by PrimitiveSolver. Currently reset_floor is the only supported option.
mhd/tfloor: set an atmosphere for temperature (mhd/pfloor is unused by the Valencia solver).
mhd/s<n>_atmosphere: set the atmosphere for scalar #<n>. Scalars are not directly floored (though a limiter defined by the EOS is applied during the conservative-to-primitive inversion), so this is a true atmosphere value and not a floor.
mhd/dthreshold: set the floor thresholding coefficient $f_\mathrm{thr}$ such that the density is reset to atmosphere when $\rho < f_\mathrm{thr}\rho_\mathrm{atmo}$. Default is 1.0.
mhd/enforce_maximum: use FOFC to enforce a relaxed discreted maximum principle (see arXiv:2409.10384). Not usually necessary, but can improve stability in some circumstances.
mhd/dmp_M: the constant $M \ge 1$ for the relaxed DMP above; this sets the permitted range without enabling FOFC to $A\in [A_\mathrm{min}/M,M A_\mathrm{max}]$, where $A\in\{\tilde{D},\tilde{\tau}\}$.

If only the <adm> block is present, the Valencia GRMHD solver is enabled but the spacetime is not evolved, making it possible to do GRMHD in a fixed but otherwise generic spacetime as defined by the problem generator. Note that AMR does not currently work for the Valencia solver unless the spacetime is also evolved. This is because the ADM variables, which are defined in the problem generator, are not currently updated during remeshing (see #600).

TOV problems

A standard TOV neutron star (in either Cartesian Schwarzschild or isotropic coordinates) with an optional poloidal magnetic field is available via the problem generator -DPROBLEM=dyngr_tov. The solver uses the shooting method with an RK4 integrator to solve the TOV equations assuming a barotropic EOS based on the EOS defined by mhd/dyn_eos. The parameters under <problem> are as follows:

Integration parameters

rhoc: central density of the TOV
npoints: buffer size for TOV integration
dr: radial step for the RK4 integrator
dfloor: density floor for determining the stopping condition during integration. This is particularly useful for tabulated equations of state.

EOS parameters

kappa: if mhd/dyn_eos = ideal, the corresponding polytropic constant for an EOS of the form $P = K\rho^\Gamma$, where $\Gamma$ is defined by mhd/gamma.
table: if using mhd/dyn_eos = eos_compose, an appropriate 1D slice in the PyCompOSE format stored in a .athtab binary table.

No additional parameters are necessary for piecewise polytropes, which use the same parameters as those specified in mhd/dyn_eos.

Other parameters and physics options

isotropic: whether to use isotropic coordinates (true) or Cartesian Schwarzschild (false, default).
minkowski: solve the correct TOV equations but initialize the metric with Minkowski space instead (false). This can be useful for debugging.
v_pert: magnitude $U$ of a velocity perturbation of the form $Wv^r = \frac{1}{2} U (3(r/R_\ast) - (r/R_\ast)^3)$, where $R_\ast$ is the TOV radius in the appropriate coordinates.
p_pert: magnitude $P_\mathrm{pert}$ of a random perturbation to the pressure.
b_norm: magnitude of a poloidal magnetic field initialized with the vector potential $A_\phi = \max\{P-P_\mathrm{cut},0\}\times(1-\rho/\rho_c)^m$
pcut: pressure cutoff $P_\mathrm{cut}$ for the magnetic field.
use_pcut_rel: whether $P_\mathrm{cut}$ is an absolute cutoff (false, default) or relative to the central pressure $P_c$.
magindex: exponent $m$ for magnetic field.

The TOV solver is independent of the problem generator and available via utils/tov/tov.hpp, so it is also possible to write your own problem generator using the TOV solver. See pgen/dyngr_tov.cpp for an example on how to use it. (Note: this is currently only true in the branch tov_refactor. main currently stores the TOV solver inside pgen/dyngr_tov.cpp and may struggle with tabulated EOSs.)

(TO DO(JF): add information on BNS problems)