GiRaFFEfood Initial Data for GiRaFFE¶Notebook Status: Validated
Validation Notes: This tutorial notebook has been confirmed to be self-consistent with its corresponding NRPy+ module, as documented below. The initial data has validated against the original GiRaFFE, as documented here.
With the GiRaFFE evolution thorn constructed, we now need to "feed" our giraffe with initial data to evolve. There are several different choices of initial data we can use here; here, we will only be implementing the "Exact Wald" initial data, given by Table 3 in the original paper:
\begin{align}
A_{\phi} &= \frac{C_0}{2} r^2 \sin^2 \theta \\
E_{\phi} &= 2 M C_0 \left( 1+ \frac {2M}{r} \right)^{-1/2} \sin^2 \theta \\
\end{align}
(the unspecified components are set to 0). Here, $C_0$ is a constant that we will set to $1$ in our simulations. Now, to use this initial data scheme, we need to transform the above into the quantities actually tracked by GiRaFFE and HydroBase: $A_i$, $B^i$, $\tilde{S}_i$, $v^i$, and $\Phi$. Of these quantities, GiRaFFEfood will only set $A_i$, $v^i$, and $\Phi=0$, then call a separate function to calculate $\tilde{S}_i$; GiRaFFE itself will call a function to set $B^i$ before the time-evolution begins. This can be done with eqs. 16 and 18, here given in that same order:
\begin{align}
v^i &= \alpha \frac{\epsilon^{ijk} E_j B_k}{B^2} -\beta^i \\
B^i &= \frac{[ijk]}{\sqrt{\gamma}} \partial_j A_k \\
\end{align}
In the simulations, $B^i$ will be calculated numerically from $A_i$; however, it will be useful to analytically calculate $B^i$ to use calculating the initial $v^i$.
This notebook is organized as follows
Here, we will import the NRPy+ core modules, set the reference metric to Cartesian, and set commonly used NRPy+ parameters. We will also set up a parameter to determine what initial data is set up, although it won't do much yet.
# Step 0: Add NRPy's directory to the path
# https://stackoverflow.com/questions/16780014/import-file-from-parent-directory
import os,sys
nrpy_dir_path = os.path.join("..")
if nrpy_dir_path not in sys.path:
sys.path.append(nrpy_dir_path)
# Step 0.a: Import the NRPy+ core modules and set the reference metric to Cartesian
import sympy as sp # SymPy: The Python computer algebra package upon which NRPy+ depends
import NRPy_param_funcs as par # NRPy+: Parameter interface
import indexedexp as ixp # NRPy+: Symbolic indexed expression (e.g., tensors, vectors, etc.) support
import GiRaFFEfood_NRPy.GiRaFFEfood_NRPy_Common_Functions as gfcf # Some useful functions for GiRaFFE initial data.
import reference_metric as rfm
par.set_parval_from_str("reference_metric::CoordSystem","Cartesian")
rfm.reference_metric()
# Step 1a: Set commonly used parameters.
thismodule = "GiRaFFEfood_NRPy_Exact_Wald"
We will first build the fundamental vectors $A_i$ and $E_i$ in spherical coordinates (see Table 3). Note that we use reference_metric.py to set $r$ and $\theta$ in terms of Cartesian coordinates; this will save us a step later when we convert to Cartesian coordinates. Since $C_0 = 1$, $$A_{\phi} = \frac{1}{2} r^2 \sin^2 \theta.$$
# Step 2: Set the vectors A and E in Spherical coordinates
def Ar_EW(r,theta,phi, **params):
return sp.sympify(0)
def Ath_EW(r,theta,phi, **params):
return sp.sympify(0)
def Aph_EW(r,theta,phi, **params):
# 1/2 r^2 \sin^2 \theta
return sp.Rational(1,2) * (r * r * sp.sin(theta)**2)
Next, we will code up the Valencia velocity, which will require us to first code the electric and magnetic fields. The magnetic field is simply $$B^i = \frac{[ijk]}{\sqrt{\gamma}} \partial_j A_k.$$
The electric field is
$$E_{\phi} = 2 M \left( 1+ \frac {2M}{r} \right)^{-1/2} \sin^2 \theta,$$
Where the other components are $0$.
We can then calculate the the velocity as $$v_{(n)}^i = \frac{\epsilon^{ijk} E_j B_k}{B^2}.$$
#Step 3: Compute v^i from A_i and E_i
def ValenciavU_func_EW(**params):
M = params["M"]
gammaDD = params["gammaDD"] # Note that this must use a Cartesian basis!
sqrtgammaDET = params["sqrtgammaDET"]
KerrSchild_radial_shift = params["KerrSchild_radial_shift"]
r = rfm.xxSph[0] + KerrSchild_radial_shift # We are setting the data up in Shifted Kerr-Schild coordinates
theta = rfm.xxSph[1]
LeviCivitaTensorUUU = ixp.LeviCivitaTensorUUU_dim3_rank3(sqrtgammaDET)
AD = gfcf.Axyz_func_spherical(Ar_EW,Ath_EW,Aph_EW,False,**params)
# For the initial data, we can analytically take the derivatives of A_i
ADdD = ixp.zerorank2()
for i in range(3):
for j in range(3):
ADdD[i][j] = sp.simplify(sp.diff(AD[i],rfm.xx_to_Cart[j]))
BU = ixp.zerorank1()
for i in range(3):
for j in range(3):
for k in range(3):
BU[i] += LeviCivitaTensorUUU[i][j][k] * ADdD[k][j]
EsphD = ixp.zerorank1()
# 2 M ( 1+ 2M/r )^{-1/2} \sin^2 \theta
EsphD[2] = 2 * M * sp.sin(theta)**2 / sp.sqrt(1+2*M/r)
ED = rfm.basis_transform_vectorD_from_rfmbasis_to_Cartesian(gfcf.Jac_dUrfm_dDCartUD, EsphD)
return gfcf.compute_ValenciavU_from_ED_and_BU(ED, BU, gammaDD)
GiRaFFEfood_NRPy.GiRaFFEfood_NRPy NRPy+ module [Back to top]¶Here, as a code validation check, we verify agreement in the SymPy expressions for the GiRaFFE Exact Wald initial data equations we intend to use between
import GiRaFFEfood_NRPy.GiRaFFEfood_NRPy as gf
import BSSN.ShiftedKerrSchild as sks
sks.ShiftedKerrSchild(True)
import reference_metric as rfm # NRPy+: Reference metric support
par.set_parval_from_str("reference_metric::CoordSystem","Cartesian")
rfm.reference_metric()
gammaSphDD = ixp.zerorank2()
for i in range(3):
for j in range(3):
gammaSphDD[i][j] += sks.gammaSphDD[i][j].subs(sks.r,rfm.xxSph[0]).subs(sks.th,rfm.xxSph[1])
gammaDD = rfm.basis_transform_tensorDD_from_rfmbasis_to_Cartesian(gfcf.Jac_dUrfm_dDCartUD,gammaSphDD)
unused_gammaUU,gammaDET = ixp.symm_matrix_inverter3x3(gammaDD)
sqrtgammaDET = sp.sqrt(gammaDET)
A_ewD = gfcf.Axyz_func_spherical(Ar_EW,Ath_EW,Aph_EW,stagger_enable = True,KerrSchild_radial_shift=sks.r0)
Valenciav_ewD = ValenciavU_func_EW(M=sks.M,KerrSchild_radial_shift=sks.r0,gammaDD=gammaDD,sqrtgammaDET=sqrtgammaDET)
gf.GiRaFFEfood_NRPy_generate_initial_data(ID_type = "ExactWald", stagger_enable = True,M=sks.M,KerrSchild_radial_shift=sks.r0,gammaDD=gammaDD,sqrtgammaDET=sqrtgammaDET)
print("Consistency check between GiRaFFEfood_NRPy tutorial and NRPy+ module:")
def consistency_check(quantity1,quantity2,string):
if quantity1-quantity2==0:
print(string+" is in agreement!")
else:
print(string+" does not agree!")
sys.exit(1)
for i in range(3):
consistency_check(Valenciav_ewD[i],gf.ValenciavU[i],"ValenciavU"+str(i))
consistency_check(A_ewD[i],gf.AD[i],"AD"+str(i))
Consistency check between GiRaFFEfood_NRPy tutorial and NRPy+ module: ValenciavU0 is in agreement! AD0 is in agreement! ValenciavU1 is in agreement! AD1 is in agreement! ValenciavU2 is in agreement! AD2 is in agreement!
The following code cell converts this Jupyter notebook into a proper, clickable $\LaTeX$-formatted PDF file. After the cell is successfully run, the generated PDF may be found in the root NRPy+ tutorial directory, with filename Tutorial-GiRaFFEfood_NRPy.pdf (Note that clicking on this link may not work; you may need to open the PDF file through another means.)
import cmdline_helper as cmd # NRPy+: Multi-platform Python command-line interface
cmd.output_Jupyter_notebook_to_LaTeXed_PDF("Tutorial-GiRaFFEfood_NRPy-Exact_Wald",location_of_template_file=os.path.join(".."))
Created Tutorial-GiRaFFEfood_NRPy-Exact_Wald.tex, and compiled LaTeX file
to PDF file Tutorial-GiRaFFEfood_NRPy-Exact_Wald.pdf