Converting Exact ADM Initial Data in the Spherical or Cartesian Basis to BSSN Initial Data in the Desired Curvilinear Basis

Author: Zach Etienne

Formatting improvements courtesy Brandon Clark

This module is meant for use only with initial data that can be represented exactly in ADM form, either in the Spherical or Cartesian basis. I.e., the ADM variables are given $\left\{\gamma_{ij}, K_{ij}, \alpha, \beta^i\right\}$ exactly as functions of $(r,\theta,\phi)$ or $(x,y,z)$, respectively. If instead the initial data are given only numerically (e.g., through an initial data solver), then the Numerical-ADM-Spherical/Cartesian-to-BSSNCurvilinear module will need to be used instead.

Notebook Status: Self-Validated

Validation Notes: This tutorial notebook has been confirmed to be self-consistent with its corresponding NRPy+ module, as documented below. Additional validation tests may have been performed, but are as yet, undocumented. (TODO)

NRPy+ Source Code for this module: BSSN/ADM_Exact_Spherical_or_Cartesian_to_BSSNCurvilinear.py

Introduction:

Given the ADM variables:

$$\left\{\gamma_{ij}, K_{ij}, \alpha, \beta^i\right\}$$

in the Spherical or Cartesian basis, and as functions of $(r,\theta,\phi)$ or $(x,y,z)$, respectively, this module documents their conversion to the BSSN variables

$$\left\{\bar{\gamma}_{i j},\bar{A}_{i j},\phi, K, \bar{\Lambda}^{i}, \alpha, \beta^i, B^i\right\},$$

in the desired curvilinear basis (given by reference_metric::CoordSystem). Then it rescales the resulting BSSNCurvilinear variables (as defined in the covariant BSSN formulation tutorial) into the form needed for solving Einstein's equations with the BSSN formulation:

$$\left\{h_{i j},a_{i j},\phi, K, \lambda^{i}, \alpha, \mathcal{V}^i, \mathcal{B}^i\right\}.$$

We will use as our core example in this module UIUC initial data, which are (as documented in their NRPy+ initial data module) given in terms of ADM variables in Spherical coordinates.

Table of Contents

$$\label{toc}$$

This notebook is organized as follows

  1. Step 1: Initialize core Python/NRPy+ modules
  2. Step 2: Desired output BSSN Curvilinear coordinate system set to Cylindrical, as a proof-of-principle
  3. Step 3: Converting ADM variables to functions of (xx0,xx1,xx2)
  4. Step 4: Applying Jacobian transformations to get in the correct xx0,xx1,xx2 basis
  5. Step 5: Call functions within BSSN.BSSN_in_terms_of_ADM (tutorial) to perform the ADM-to-BSSN conversion
  6. Step 6: Code Validation against BSSN.ADM_Exact_Spherical_or_Cartesian_to_BSSNCurvilinear NRPy+ module
  7. Step 7: Output this notebook to $\LaTeX$-formatted PDF file

Step 1: Initialize core Python/NRPy+ modules [Back to top]

$$\label{initializenrpy}$$
In [1]:
# Step P1: Import needed NRPy+ core modules:
import sympy as sp                # SymPy: The Python computer algebra package upon which NRPy+ depends
import indexedexp as ixp          # NRPy+: Symbolic indexed expression (e.g., tensors, vectors, etc.) support
import reference_metric as rfm    # NRPy+: Reference metric support
import NRPy_param_funcs as par    # NRPy+: Parameter interface
import sys                        # Standard Python module for multiplatform OS-level functions

Step 2: Desired output BSSN Curvilinear coordinate system set to Cylindrical, as a proof-of-principle [Back to top]

$$\label{cylindrical}$$
In [2]:
# The ADM & BSSN formalisms only work in 3D; they are 3+1 decompositions of Einstein's equations.
#    To implement axisymmetry or spherical symmetry, simply set all spatial derivatives in
#    the relevant angular directions to zero; DO NOT SET DIM TO ANYTHING BUT 3.

# Step P1: Set spatial dimension (must be 3 for BSSN)
DIM = 3

# Set the desired *output* coordinate system to Cylindrical:
par.set_parval_from_str("reference_metric::CoordSystem","Cylindrical")
rfm.reference_metric()

# Import UIUC Black Hole initial data
import BSSN.UIUCBlackHole as uibh
uibh.UIUCBlackHole(ComputeADMGlobalsOnly=True)
Sph_r_th_ph_or_Cart_xyz = [uibh.r,uibh.th,uibh.ph]
alphaSphorCart    = uibh.alphaSph
betaSphorCartU    = uibh.betaSphU
BSphorCartU       = uibh.BSphU
gammaSphorCartDD  = uibh.gammaSphDD
KSphorCartDD      = uibh.KSphDD

Step 3: Converting ADM variables to functions of ${\rm xx0},{\rm xx1},{\rm xx2}$ [Back to top]

$$\label{admxx0xx1xx2}$$

ADM variables are given as functions of $(r,\theta,\phi)$ or $(x,y,z)$. We convert them to functions of (xx0,xx1,xx2) using SymPy's subs() function.

In [3]:
# Step 3: All input quantities are in terms of r,th,ph or x,y,z. We want them in terms
#         of xx0,xx1,xx2, so here we call sympify_integers__replace_rthph() to replace
#         r,th,ph or x,y,z, respectively, with the appropriate functions of xx0,xx1,xx2
#         as defined for this particular reference metric in reference_metric.py's
#         xxSph[] or xxCart[], respectively:

# UIUC Black Hole initial data are given in Spherical coordinates.
CoordType_in = "Spherical"

# Make sure that rfm.reference_metric() has been called.
#    We'll need the variables it defines throughout this module.
if rfm.have_already_called_reference_metric_function == False:
    print("Error. Called Convert_Spherical_ADM_to_BSSN_curvilinear() without")
    print("       first setting up reference metric, by calling rfm.reference_metric().")
    sys.exit(1)

# Note that substitution only works when the variable is not an integer. Hence the
#         if isinstance(...,...) stuff:
def sympify_integers__replace_rthph_or_Cartxyz(obj, rthph_or_xyz, rthph_or_xyz_of_xx):
    if isinstance(obj, int):
        return sp.sympify(obj)
    else:
        return obj.subs(rthph_or_xyz[0], rthph_or_xyz_of_xx[0]).\
            subs(rthph_or_xyz[1], rthph_or_xyz_of_xx[1]).\
            subs(rthph_or_xyz[2], rthph_or_xyz_of_xx[2])

r_th_ph_or_Cart_xyz_of_xx = []
if CoordType_in == "Spherical":
    r_th_ph_or_Cart_xyz_of_xx = rfm.xxSph
elif CoordType_in == "Cartesian":
    r_th_ph_or_Cart_xyz_of_xx = rfm.xxCart
else:
    print("Error: Can only convert ADM Cartesian or Spherical initial data to BSSN Curvilinear coords.")
    sys.exit(1)

alphaSphorCart = sympify_integers__replace_rthph_or_Cartxyz(
    alphaSphorCart, Sph_r_th_ph_or_Cart_xyz, r_th_ph_or_Cart_xyz_of_xx)
for i in range(DIM):
    betaSphorCartU[i] = sympify_integers__replace_rthph_or_Cartxyz(
        betaSphorCartU[i], Sph_r_th_ph_or_Cart_xyz, r_th_ph_or_Cart_xyz_of_xx)
    BSphorCartU[i]    = sympify_integers__replace_rthph_or_Cartxyz(
        BSphorCartU[i],    Sph_r_th_ph_or_Cart_xyz, r_th_ph_or_Cart_xyz_of_xx)
    for j in range(DIM):
        gammaSphorCartDD[i][j] = sympify_integers__replace_rthph_or_Cartxyz(
            gammaSphorCartDD[i][j], Sph_r_th_ph_or_Cart_xyz, r_th_ph_or_Cart_xyz_of_xx)
        KSphorCartDD[i][j]     = sympify_integers__replace_rthph_or_Cartxyz(
            KSphorCartDD[i][j],     Sph_r_th_ph_or_Cart_xyz, r_th_ph_or_Cart_xyz_of_xx)

Step 4: Applying Jacobian transformations to get in the correct xx0,xx1,xx2 basis [Back to top]

$$\label{adm_jacobian}$$

All ADM initial data quantities are now functions of xx0,xx1,xx2, but they are still in the Spherical or Cartesian basis. We can now directly apply Jacobian transformations to get them in the correct xx0,xx1,xx2 basis. The following discussion holds for either Spherical or Cartesian input data, so for simplicity let's just assume the data are given in Spherical coordinates.

All ADM tensors and vectors are in the Spherical coordinate basis $x^i_{\rm Sph} = (r,\theta,\phi)$, but we need them in the curvilinear coordinate basis $x^i_{\rm rfm}=$(xx0,xx1,xx2) set by the "reference_metric::CoordSystem" variable. Empirically speaking, it is far easier to write (x(xx0,xx1,xx2),y(xx0,xx1, xx2),z(xx0,xx1,xx2)) than the inverse, so we will compute the Jacobian matrix

$$ {\rm Jac\_dUSph\_dDrfmUD[i][j]} = \frac{\partial x^i_{\rm Sph}}{\partial x^j_{\rm rfm}}, $$

via exact differentiation (courtesy SymPy), and the inverse Jacobian $$ {\rm Jac\_dUrfm\_dDSphUD[i][j]} = \frac{\partial x^i_{\rm rfm}}{\partial x^j_{\rm Sph}}, $$

using NRPy+'s generic_matrix_inverter3x3() function. In terms of these, the transformation of BSSN tensors from Spherical to "reference_metric::CoordSystem" coordinates may be written:

\begin{align} \beta^i_{\rm rfm} &= \frac{\partial x^i_{\rm rfm}}{\partial x^\ell_{\rm Sph}} \beta^\ell_{\rm Sph}\\ B^i_{\rm rfm} &= \frac{\partial x^i_{\rm rfm}}{\partial x^\ell_{\rm Sph}} B^\ell_{\rm Sph}\\ \gamma^{\rm rfm}_{ij} &= \frac{\partial x^\ell_{\rm Sph}}{\partial x^i_{\rm rfm}} \frac{\partial x^m_{\rm Sph}}{\partial x^j_{\rm rfm}} \gamma^{\rm Sph}_{\ell m}\\ K^{\rm rfm}_{ij} &= \frac{\partial x^\ell_{\rm Sph}}{\partial x^i_{\rm rfm}} \frac{\partial x^m_{\rm Sph}}{\partial x^j_{\rm rfm}} K^{\rm Sph}_{\ell m} \end{align}
In [4]:
# Step 2: All ADM initial data quantities are now functions of xx0,xx1,xx2, but
#         they are still in the Spherical or Cartesian basis. We can now directly apply
#         Jacobian transformations to get them in the correct xx0,xx1,xx2 basis:

# alpha is a scalar, so no Jacobian transformation is necessary.
alpha = alphaSphorCart

Jac_dUSphorCart_dDrfmUD = ixp.zerorank2()
for i in range(DIM):
    for j in range(DIM):
        Jac_dUSphorCart_dDrfmUD[i][j] = sp.diff(r_th_ph_or_Cart_xyz_of_xx[i],rfm.xx[j])

Jac_dUrfm_dDSphorCartUD, dummyDET = ixp.generic_matrix_inverter3x3(Jac_dUSphorCart_dDrfmUD)

betaU   = ixp.zerorank1()
BU      = ixp.zerorank1()
gammaDD = ixp.zerorank2()
KDD     = ixp.zerorank2()
for i in range(DIM):
    for j in range(DIM):
        betaU[i] += Jac_dUrfm_dDSphorCartUD[i][j] * betaSphorCartU[j]
        BU[i]    += Jac_dUrfm_dDSphorCartUD[i][j] * BSphorCartU[j]
        for k in range(DIM):
            for l in range(DIM):
                gammaDD[i][j] += Jac_dUSphorCart_dDrfmUD[k][i]*Jac_dUSphorCart_dDrfmUD[l][j] * gammaSphorCartDD[k][l]
                KDD[i][j]     += Jac_dUSphorCart_dDrfmUD[k][i]*Jac_dUSphorCart_dDrfmUD[l][j] *     KSphorCartDD[k][l]

Step 5: Call functions within BSSN.BSSN_in_terms_of_ADM (tutorial) to perform the ADM-to-BSSN conversion [Back to top]

$$\label{adm2bssn}$$

All ADM quantities were input into this function in the Spherical or Cartesian basis, as functions of $r,\theta,\phi$ or $x,y,z$, respectively. In Step 3 and Step 4 above, we converted them to the xx0,xx1,xx2 basis, and as functions of xx0,xx1,xx2. Here we convert ADM 3-metric, extrinsic curvature, and gauge quantities in the xx0,xx1,xx2 (a.k.a. "rfm") basis to their BSSN Curvilinear counterparts, in the same basis.

In [5]:
import BSSN.BSSN_in_terms_of_ADM as BitoA
BitoA.gammabarDD_hDD(                   gammaDD)
BitoA.trK_AbarDD_aDD(                   gammaDD,KDD)
BitoA.LambdabarU_lambdaU__exact_gammaDD(gammaDD)
BitoA.cf_from_gammaDD(                  gammaDD)
BitoA.betU_vetU(                        betaU,BU)
hDD     = BitoA.hDD
trK     = BitoA.trK
aDD     = BitoA.aDD
lambdaU = BitoA.lambdaU
cf      = BitoA.cf
vetU    = BitoA.vetU
betU    = BitoA.betU

Step 6: Code Validation against BSSN.ADM_Exact_Spherical_or_Cartesian_to_BSSNCurvilinear module [Back to top]

$$\label{code_validation}$$

Here, as a code validation check, we verify agreement in the SymPy expressions for BrillLindquist initial data between

  1. this tutorial and
  2. the NRPy+ BSSN.ADM_Exact_Spherical_or_Cartesian_to_BSSNCurvilinear module.

By default, we analyze these expressions in Spherical coordinates, though other coordinate systems may be chosen.

In [6]:
import BSSN.UIUCBlackHole as uibh
import BSSN.ADM_Exact_Spherical_or_Cartesian_to_BSSNCurvilinear as ADMtoBSSN
returnfunction = uibh.UIUCBlackHole()
mod_cf,mod_hDD,mod_lambdaU,mod_aDD,mod_trK,mod_alpha,mod_vetU,mod_betU = \
    ADMtoBSSN.Convert_Spherical_or_Cartesian_ADM_to_BSSN_curvilinear("Spherical",uibh.Sph_r_th_ph,
                                            uibh.gammaSphDD, uibh.KSphDD, uibh.alphaSph, uibh.betaSphU, uibh.BSphU)

print("Consistency check between this tutorial notebook and BSSN.ADM_Exact_Spherical_or_Cartesian_to_BSSNCurvilinear NRPy+ module: ALL SHOULD BE ZERO.")

print("cf - mod_cf = " + str(cf - mod_cf))
print("trK - mod_trK = " + str(trK - mod_trK))
print("alpha - mod_alpha = " + str(alpha - mod_alpha))

for i in range(DIM):
    print("vetU["+str(i)+"] - mod_vetU["+str(i)+"] = " + str(vetU[i] - mod_vetU[i]))
    print("betU["+str(i)+"] - mod_betU["+str(i)+"] = " + str(betU[i] - mod_betU[i]))
    print("lambdaU["+str(i)+"] - mod_lambdaU["+str(i)+"] = " + str(lambdaU[i] - mod_lambdaU[i]))
    for j in range(DIM):
        print("hDD["+str(i)+"]["+str(j)+"] - mod_hDD["+str(i)+"]["+str(j)+"] = "
              + str(hDD[i][j] - mod_hDD[i][j]))
        print("aDD["+str(i)+"]["+str(j)+"] - mod_aDD["+str(i)+"]["+str(j)+"] = "
              + str(aDD[i][j] - mod_aDD[i][j]))

# If you wish to generate & analyze C code output, uncomment the following:
# import os, shutil            # Standard Python modules for multiplatform OS-level functions
# import cmdline_helper as cmd # NRPy+: Multi-platform Python command-line interface
# # Step P2: Create C code output directory:
# Ccodesdir = os.path.join("BSSN_Exact_ADM_validation/")
# # First remove C code output directory if it exists
# # Courtesy https://stackoverflow.com/questions/303200/how-do-i-remove-delete-a-folder-that-is-not-empty
# # !rm -r ScalarWaveCurvilinear_Playground_Ccodes
# shutil.rmtree(Ccodesdir, ignore_errors=True)
# # Then create a fresh directory
# cmd.mkdir(Ccodesdir)
# with open(os.path.join(Ccodedir,"UIUCBlackHole-CylindricalTest.h"),"w") as file:
#     file.write(uibh.returnfunction)
Consistency check between this tutorial notebook and BSSN.ADM_Exact_Spherical_or_Cartesian_to_BSSNCurvilinear NRPy+ module: ALL SHOULD BE ZERO.
cf - mod_cf = 0
trK - mod_trK = 0
alpha - mod_alpha = 0
vetU[0] - mod_vetU[0] = 0
betU[0] - mod_betU[0] = 0
lambdaU[0] - mod_lambdaU[0] = 0
hDD[0][0] - mod_hDD[0][0] = 0
aDD[0][0] - mod_aDD[0][0] = 0
hDD[0][1] - mod_hDD[0][1] = 0
aDD[0][1] - mod_aDD[0][1] = 0
hDD[0][2] - mod_hDD[0][2] = 0
aDD[0][2] - mod_aDD[0][2] = 0
vetU[1] - mod_vetU[1] = 0
betU[1] - mod_betU[1] = 0
lambdaU[1] - mod_lambdaU[1] = 0
hDD[1][0] - mod_hDD[1][0] = 0
aDD[1][0] - mod_aDD[1][0] = 0
hDD[1][1] - mod_hDD[1][1] = 0
aDD[1][1] - mod_aDD[1][1] = 0
hDD[1][2] - mod_hDD[1][2] = 0
aDD[1][2] - mod_aDD[1][2] = 0
vetU[2] - mod_vetU[2] = 0
betU[2] - mod_betU[2] = 0
lambdaU[2] - mod_lambdaU[2] = 0
hDD[2][0] - mod_hDD[2][0] = 0
aDD[2][0] - mod_aDD[2][0] = 0
hDD[2][1] - mod_hDD[2][1] = 0
aDD[2][1] - mod_aDD[2][1] = 0
hDD[2][2] - mod_hDD[2][2] = 0
aDD[2][2] - mod_aDD[2][2] = 0

Step 7: Output this notebook to $\LaTeX$-formatted PDF file [Back to top]

$$\label{latex_pdf_output}$$

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-ADM_Initial_Data-Converting_Exact_ADM_Spherical_or_Cartesian_to_BSSNCurvilinear.pdf (Note that clicking on this link may not work; you may need to open the PDF file through another means.)

In [7]:
import cmdline_helper as cmd    # NRPy+: Multi-platform Python command-line interface
cmd.output_Jupyter_notebook_to_LaTeXed_PDF("Tutorial-ADM_Initial_Data-Converting_Exact_ADM_Spherical_or_Cartesian_to_BSSNCurvilinear")
Created Tutorial-ADM_Initial_Data-
    Converting_Exact_ADM_Spherical_or_Cartesian_to_BSSNCurvilinear.tex, and
    compiled LaTeX file to PDF file Tutorial-ADM_Initial_Data-
    Converting_Exact_ADM_Spherical_or_Cartesian_to_BSSNCurvilinear.pdf