{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"\n",
"# ADM Quantities in terms of BSSN Quantities\n",
"\n",
"## Author: Zach Etienne\n",
"### Formatting improvements courtesy Brandon Clark\n",
"\n",
"## This notebook demonstrates the conversion between the 4-metric $g_{\\mu\\nu}$ and the ADM or BSSN variables. Initially, ADM quantities are defined directly or in terms of BSSN quantities, followed by writing the 4-metric and its inverse using these quantities. The notebook then details the conversion of the ADM/BSSN metric quantities back into the 4-metric and validates the result against the NRPy+ BSSN.ADMBSSN_tofrom_4metric module.\n",
"\n",
"**Notebook Status:** Self-Validated \n",
"\n",
"**Validation Notes:** This tutorial notebook has been confirmed to be self-consistent with its corresponding NRPy+ module, as documented [below](#code_validation). **Additional validation tests may have been performed, but are as yet, undocumented. (TODO)**\n",
"\n",
"### NRPy+ Source Code for this module: [ADM_in_terms_of_BSSN.py](../edit/BSSN/ADM_in_terms_of_BSSN.py)\n",
"\n",
"## Introduction:\n",
"This tutorial notebook constructs all quantities in the [ADM formalism](https://en.wikipedia.org/wiki/ADM_formalism) (see also Chapter 2 in Baumgarte & Shapiro's book *Numerical Relativity*) in terms of quantities in our adopted (covariant, tensor-rescaled) BSSN formalism. That is to say, we will write the ADM quantities $\\left\\{\\gamma_{ij},K_{ij},\\alpha,\\beta^i\\right\\}$ and their derivatives in terms of the BSSN quantities $\\left\\{\\bar{\\gamma}_{ij},\\text{cf},\\bar{A}_{ij},\\text{tr}K,\\alpha,\\beta^i\\right\\}$ and their derivatives.\n",
"\n",
"### A Note on Notation:\n",
"\n",
"As is standard in NRPy+, \n",
"\n",
"* Greek indices refer to four-dimensional quantities where the zeroth component indicates the temporal (time) component.\n",
"* Latin indices refer to three-dimensional quantities. This is somewhat counterintuitive since Python always indexes its lists starting from 0. As a result, the zeroth component of three-dimensional quantities will necessarily indicate the first *spatial* direction.\n",
"\n",
"As a corollary, any expressions in NRPy+ involving mixed Greek and Latin indices will need to offset one set of indices by one; a Latin index in a four-vector will be incremented and a Greek index in a three-vector will be decremented (however, the latter case does not occur in this tutorial notebook)."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Table of Contents\n",
"$$\\label{toc}$$\n",
"\n",
"This notebook is organized as follows\n",
"\n",
"1. [Step 1](#initializenrpy): Initialize core Python/NRPy+ modules\n",
"1. [Step 2](#threemetric): The ADM three-metric $\\gamma_{ij}$ and its derivatives in terms of rescaled BSSN quantities\n",
" 1. [Step 2.a](#derivatives_e4phi): Derivatives of $e^{4\\phi}$\n",
" 1. [Step 2.b](#derivatives_adm_3metric): Derivatives of the ADM three-metric: $\\gamma_{ij,k}$ and $\\gamma_{ij,kl}$\n",
" 1. [Step 2.c](#christoffel): Christoffel symbols $\\Gamma^i_{jk}$ associated with the ADM 3-metric $\\gamma_{ij}$\n",
"1. [Step 3](#extrinsiccurvature): The ADM extrinsic curvature $K_{ij}$ and its derivatives in terms of rescaled BSSN quantities\n",
"1. [Step 4](#code_validation): Code Validation against `BSSN.ADM_in_terms_of_BSSN` NRPy+ module\n",
"1. [Step 5](#latex_pdf_output): Output this notebook to $\\LaTeX$-formatted PDF file"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 1: Initialize core Python/NRPy+ modules \\[Back to [top](#toc)\\]\n",
"$$\\label{initializenrpy}$$\n",
"\n",
"Let's start by importing all the needed modules from Python/NRPy+:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:40.170776Z",
"iopub.status.busy": "2021-03-07T17:06:40.169649Z",
"iopub.status.idle": "2021-03-07T17:06:41.077742Z",
"shell.execute_reply": "2021-03-07T17:06:41.077043Z"
}
},
"outputs": [],
"source": [
"# Step 1.a: Import all needed modules from NRPy+\n",
"import NRPy_param_funcs as par # NRPy+: parameter interface\n",
"import sympy as sp # SymPy: The Python computer algebra package upon which NRPy+ depends\n",
"import indexedexp as ixp # NRPy+: Symbolic indexed expression (e.g., tensors, vectors, etc.) support\n",
"import reference_metric as rfm # NRPy+: Reference metric support\n",
"import sys # Standard Python module for multiplatform OS-level functions\n",
"\n",
"# Step 1.b: Set the coordinate system for the numerical grid\n",
"par.set_parval_from_str(\"reference_metric::CoordSystem\",\"Spherical\")\n",
"\n",
"# Step 1.c: Given the chosen coordinate system, set up\n",
"# corresponding reference metric and needed\n",
"# reference metric quantities\n",
"# The following function call sets up the reference metric\n",
"# and related quantities, including rescaling matrices ReDD,\n",
"# ReU, and hatted quantities.\n",
"rfm.reference_metric()\n",
"\n",
"# Step 1.d: Set spatial dimension (must be 3 for BSSN, as BSSN is\n",
"# a 3+1-dimensional decomposition of the general\n",
"# relativistic field equations)\n",
"DIM = 3\n",
"\n",
"# Step 1.e: Import all basic (unrescaled) BSSN scalars & tensors\n",
"import BSSN.BSSN_quantities as Bq\n",
"Bq.BSSN_basic_tensors()\n",
"gammabarDD = Bq.gammabarDD\n",
"cf = Bq.cf\n",
"AbarDD = Bq.AbarDD\n",
"trK = Bq.trK\n",
"\n",
"Bq.gammabar__inverse_and_derivs()\n",
"gammabarDD_dD = Bq.gammabarDD_dD\n",
"gammabarDD_dDD = Bq.gammabarDD_dDD\n",
"\n",
"Bq.AbarUU_AbarUD_trAbar_AbarDD_dD()\n",
"AbarDD_dD = Bq.AbarDD_dD"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 2: The ADM three-metric $\\gamma_{ij}$ and its derivatives in terms of rescaled BSSN quantities. \\[Back to [top](#toc)\\]\n",
"$$\\label{threemetric}$$\n",
"\n",
"The ADM three-metric is written in terms of the covariant BSSN three-metric tensor as (Eqs. 2 and 3 of [Ruchlin *et al.*](https://arxiv.org/pdf/1712.07658.pdf)):\n",
"$$\n",
"\\gamma_{ij} = \\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right)^{1/3} \\bar{\\gamma}_{i j},\n",
"$$\n",
"where $\\gamma=\\det{\\gamma_{ij}}$ and $\\bar{\\gamma}=\\det{\\bar{\\gamma}_{ij}}$. \n",
"\n",
"The \"standard\" BSSN conformal factor $\\phi$ is given by (Eq. 3 of [Ruchlin *et al.*](https://arxiv.org/pdf/1712.07658.pdf)):\n",
"\n",
"\\begin{align}\n",
"\\phi &= \\frac{1}{12} \\log\\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right) \\\\\n",
"\\implies e^{\\phi} &= \\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right)^{1/12} \\\\\n",
"\\implies e^{4 \\phi} &= \\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right)^{1/3}\n",
"\\end{align}\n",
"\n",
"Thus the ADM three-metric may be written in terms of the BSSN three-metric and conformal factor $\\phi$ as\n",
"\n",
"$$\n",
"\\gamma_{ij} = e^{4 \\phi} \\bar{\\gamma}_{i j}.\n",
"$$\n",
"\n",
"NRPy+'s implementation of BSSN allows for $\\phi$ and two other alternative conformal factors to be defined:\n",
"\n",
"\\begin{align}\n",
"\\chi &= e^{-4\\phi} \\\\\n",
"W &= e^{-2\\phi},\n",
"\\end{align}\n",
"\n",
"Thus if `\"BSSN_quantities::EvolvedConformalFactor_cf\"` is set to `\"chi\"`, then\n",
"\n",
"\\begin{align}\n",
"\\gamma_{ij} &= \\frac{1}{\\chi} \\bar{\\gamma}_{i j} \\\\\n",
"&= \\frac{1}{\\text{cf}} \\bar{\\gamma}_{i j},\n",
"\\end{align}\n",
"\n",
"and if `\"BSSN_quantities::EvolvedConformalFactor_cf\"` is set to `\"W\"`, then\n",
"\\begin{align}\n",
"\\gamma_{ij} &= \\frac{1}{W^2} \\bar{\\gamma}_{i j} \\\\\n",
"&= \\frac{1}{\\text{cf}^2} \\bar{\\gamma}_{i j}.\n",
"\\end{align}"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:41.086856Z",
"iopub.status.busy": "2021-03-07T17:06:41.085949Z",
"iopub.status.idle": "2021-03-07T17:06:41.088286Z",
"shell.execute_reply": "2021-03-07T17:06:41.088841Z"
}
},
"outputs": [],
"source": [
"# Step 2: The ADM three-metric gammaDD and its\n",
"# derivatives in terms of BSSN quantities.\n",
"gammaDD = ixp.zerorank2()\n",
"\n",
"exp4phi = sp.sympify(0)\n",
"if par.parval_from_str(\"EvolvedConformalFactor_cf\") == \"phi\":\n",
" exp4phi = sp.exp(4*cf)\n",
"elif par.parval_from_str(\"EvolvedConformalFactor_cf\") == \"chi\":\n",
" exp4phi = (1 / cf)\n",
"elif par.parval_from_str(\"EvolvedConformalFactor_cf\") == \"W\":\n",
" exp4phi = (1 / cf**2)\n",
"else:\n",
" print(\"Error EvolvedConformalFactor_cf type = \\\"\"+par.parval_from_str(\"EvolvedConformalFactor_cf\")+\"\\\" unknown.\")\n",
" sys.exit(1)\n",
"\n",
"for i in range(DIM):\n",
" for j in range(DIM):\n",
" gammaDD[i][j] = exp4phi*gammabarDD[i][j]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"## Step 2.a: Derivatives of $e^{4\\phi}$ \\[Back to [top](#toc)\\]\n",
"$$\\label{derivatives_e4phi}$$\n",
"\n",
"To compute derivatives of $\\gamma_{ij}$ in terms of BSSN variables and their derivatives, we will first need derivatives of $e^{4\\phi}$ in terms of the conformal BSSN variable `cf`.\n",
"\n",
"\\begin{align}\n",
"\\frac{\\partial}{\\partial x^i} e^{4\\phi} &= 4 e^{4\\phi} \\phi_{,i} \\\\\n",
"\\implies \\frac{\\partial}{\\partial x^j} \\frac{\\partial}{\\partial x^i} e^{4\\phi} &= \\frac{\\partial}{\\partial x^j} \\left(4 e^{4\\phi} \\phi_{,i}\\right) \\\\\n",
"&= 16 e^{4\\phi} \\phi_{,i} \\phi_{,j} + 4 e^{4\\phi} \\phi_{,ij}\n",
"\\end{align}\n",
"\n",
"Thus computing first and second derivatives of $e^{4\\phi}$ in terms of the BSSN quantity `cf` requires only that we evaluate $\\phi_{,i}$ and $\\phi_{,ij}$ in terms of $e^{4\\phi}$ (computed above in terms of `cf`) and derivatives of `cf`:\n",
"\n",
"If `\"BSSN_quantities::EvolvedConformalFactor_cf\"` is set to `\"phi\"`, then\n",
"\\begin{align}\n",
"\\phi_{,i} &= \\text{cf}_{,i} \\\\\n",
"\\phi_{,ij} &= \\text{cf}_{,ij}\n",
"\\end{align}\n",
"\n",
"If `\"BSSN_quantities::EvolvedConformalFactor_cf\"` is set to `\"chi\"`, then\n",
"\\begin{align}\n",
"\\text{cf} = e^{-4\\phi} \\implies \\text{cf}_{,i} &= -4 e^{-4\\phi} \\phi_{,i} \\\\\n",
"\\implies \\phi_{,i} &= -\\frac{e^{4\\phi}}{4} \\text{cf}_{,i} \\\\\n",
"\\implies \\phi_{,ij} &= -e^{4\\phi} \\phi_{,j} \\text{cf}_{,i} -\\frac{e^{4\\phi}}{4} \\text{cf}_{,ij}\\\\\n",
"&= -e^{4\\phi} \\left(-\\frac{e^{4\\phi}}{4} \\text{cf}_{,j}\\right) \\text{cf}_{,i} -\\frac{e^{4\\phi}}{4} \\text{cf}_{,ij} \\\\\n",
"&= \\frac{1}{4} \\left[\\left(e^{4\\phi}\\right)^2 \\text{cf}_{,i} \\text{cf}_{,j} -e^{4\\phi} \\text{cf}_{,ij}\\right] \\\\\n",
"\\end{align}\n",
"\n",
"If `\"BSSN_quantities::EvolvedConformalFactor_cf\"` is set to `\"W\"`, then\n",
"\\begin{align}\n",
"\\text{cf} = e^{-2\\phi} \\implies \\text{cf}_{,i} &= -2 e^{-2\\phi} \\phi_{,i} \\\\\n",
"\\implies \\phi_{,i} &= -\\frac{e^{2\\phi}}{2} \\text{cf}_{,i} \\\\\n",
"\\implies \\phi_{,ij} &= -e^{2\\phi} \\phi_{,j} \\text{cf}_{,i} -\\frac{e^{2\\phi}}{2} \\text{cf}_{,ij}\\\\\n",
"&= -e^{2\\phi} \\left(-\\frac{e^{2\\phi}}{2} \\text{cf}_{,j}\\right) \\text{cf}_{,i} -\\frac{e^{2\\phi}}{2} \\text{cf}_{,ij} \\\\\n",
"&= \\frac{1}{2} \\left[e^{4\\phi} \\text{cf}_{,i} \\text{cf}_{,j} -e^{2\\phi} \\text{cf}_{,ij}\\right] \\\\\n",
"\\end{align}"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:41.113499Z",
"iopub.status.busy": "2021-03-07T17:06:41.112528Z",
"iopub.status.idle": "2021-03-07T17:06:41.114886Z",
"shell.execute_reply": "2021-03-07T17:06:41.115417Z"
}
},
"outputs": [],
"source": [
"# Step 2.a: Derivatives of $e^{4\\phi}$\n",
"phidD = ixp.zerorank1()\n",
"phidDD = ixp.zerorank2()\n",
"cf_dD = ixp.declarerank1(\"cf_dD\")\n",
"cf_dDD = ixp.declarerank2(\"cf_dDD\",\"sym01\")\n",
"if par.parval_from_str(\"EvolvedConformalFactor_cf\") == \"phi\":\n",
" for i in range(DIM):\n",
" phidD[i] = cf_dD[i]\n",
" for j in range(DIM):\n",
" phidDD[i][j] = cf_dDD[i][j]\n",
"elif par.parval_from_str(\"EvolvedConformalFactor_cf\") == \"chi\":\n",
" for i in range(DIM):\n",
" phidD[i] = -sp.Rational(1,4)*exp4phi*cf_dD[i]\n",
" for j in range(DIM):\n",
" phidDD[i][j] = sp.Rational(1,4)*( exp4phi**2*cf_dD[i]*cf_dD[j] - exp4phi*cf_dDD[i][j] )\n",
"elif par.parval_from_str(\"EvolvedConformalFactor_cf\") == \"W\":\n",
" exp2phi = (1 / cf)\n",
" for i in range(DIM):\n",
" phidD[i] = -sp.Rational(1,2)*exp2phi*cf_dD[i]\n",
" for j in range(DIM):\n",
" phidDD[i][j] = sp.Rational(1,2)*( exp4phi*cf_dD[i]*cf_dD[j] - exp2phi*cf_dDD[i][j] )\n",
"else:\n",
" print(\"Error EvolvedConformalFactor_cf type = \\\"\"+par.parval_from_str(\"EvolvedConformalFactor_cf\")+\"\\\" unknown.\")\n",
" sys.exit(1)\n",
"\n",
"exp4phidD = ixp.zerorank1()\n",
"exp4phidDD = ixp.zerorank2()\n",
"for i in range(DIM):\n",
" exp4phidD[i] = 4*exp4phi*phidD[i]\n",
" for j in range(DIM):\n",
" exp4phidDD[i][j] = 16*exp4phi*phidD[i]*phidD[j] + 4*exp4phi*phidDD[i][j]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"## Step 2.b: Derivatives of the ADM three-metric: $\\gamma_{ij,k}$ and $\\gamma_{ij,kl}$ \\[Back to [top](#toc)\\]\n",
"$$\\label{derivatives_adm_3metric}$$\n",
"\n",
"Recall the relation between the ADM three-metric $\\gamma_{ij}$, the BSSN conformal three-metric $\\bar{\\gamma}_{i j}$, and the BSSN conformal factor $\\phi$:\n",
"\n",
"$$\n",
"\\gamma_{ij} = e^{4 \\phi} \\bar{\\gamma}_{i j}.\n",
"$$\n",
"\n",
"Now that we have constructed derivatives of $e^{4 \\phi}$ in terms of the chosen BSSN conformal factor `cf`, and the [BSSN.BSSN_quantities module](../edit/BSSN/BSSN_quantities.py) ([**tutorial**](Tutorial-BSSN_quantities.ipynb)) defines derivatives of $\\bar{\\gamma}_{ij}$ in terms of rescaled BSSN variables, derivatives of $\\gamma_{ij}$ can be immediately constructed using the product rule:\n",
"\n",
"\\begin{align}\n",
"\\gamma_{ij,k} &= \\left(e^{4 \\phi}\\right)_{,k} \\bar{\\gamma}_{i j} + e^{4 \\phi} \\bar{\\gamma}_{ij,k} \\\\\n",
"\\gamma_{ij,kl} &= \\left(e^{4 \\phi}\\right)_{,kl} \\bar{\\gamma}_{i j} + \\left(e^{4 \\phi}\\right)_{,k} \\bar{\\gamma}_{i j,l} + \\left(e^{4 \\phi}\\right)_{,l} \\bar{\\gamma}_{ij,k} + e^{4 \\phi} \\bar{\\gamma}_{ij,kl}\n",
"\\end{align}"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:41.177860Z",
"iopub.status.busy": "2021-03-07T17:06:41.154025Z",
"iopub.status.idle": "2021-03-07T17:06:41.179989Z",
"shell.execute_reply": "2021-03-07T17:06:41.180544Z"
}
},
"outputs": [],
"source": [
"# Step 2.b: Derivatives of gammaDD, the ADM three-metric\n",
"gammaDDdD = ixp.zerorank3()\n",
"gammaDDdDD = ixp.zerorank4()\n",
"\n",
"for i in range(DIM):\n",
" for j in range(DIM):\n",
" for k in range(DIM):\n",
" gammaDDdD[i][j][k] = exp4phidD[k]*gammabarDD[i][j] + exp4phi*gammabarDD_dD[i][j][k]\n",
" for l in range(DIM):\n",
" gammaDDdDD[i][j][k][l] = exp4phidDD[k][l]*gammabarDD[i][j] + \\\n",
" exp4phidD[k]*gammabarDD_dD[i][j][l] + \\\n",
" exp4phidD[l]*gammabarDD_dD[i][j][k] + \\\n",
" exp4phi*gammabarDD_dDD[i][j][k][l]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"## Step 2.c: Christoffel symbols $\\Gamma^i_{jk}$ associated with the ADM 3-metric $\\gamma_{ij}$ \\[Back to [top](#toc)\\]\n",
"$$\\label{christoffel}$$\n",
"\n",
"The 3-metric analog to the definition of Christoffel symbol (Eq. 1.18) in Baumgarte & Shapiro's *Numerical Relativity* is given by\n",
"$$\n",
"\\Gamma^i_{jk} = \\frac{1}{2} \\gamma^{il} \\left(\\gamma_{lj,k} + \\gamma_{lk,j} - \\gamma_{jk,l} \\right),\n",
"$$\n",
"which we implement here:"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:41.235764Z",
"iopub.status.busy": "2021-03-07T17:06:41.214121Z",
"iopub.status.idle": "2021-03-07T17:06:41.238563Z",
"shell.execute_reply": "2021-03-07T17:06:41.237912Z"
}
},
"outputs": [],
"source": [
"# Step 2.c: 3-Christoffel symbols associated with ADM 3-metric gammaDD\n",
"# Step 2.c.i: First compute the inverse 3-metric gammaUU:\n",
"gammaUU, detgamma = ixp.symm_matrix_inverter3x3(gammaDD)\n",
"\n",
"GammaUDD = ixp.zerorank3()\n",
"\n",
"for i in range(DIM):\n",
" for j in range(DIM):\n",
" for k in range(DIM):\n",
" for l in range(DIM):\n",
" GammaUDD[i][j][k] += sp.Rational(1,2)*gammaUU[i][l]* \\\n",
" (gammaDDdD[l][j][k] + gammaDDdD[l][k][j] - gammaDDdD[j][k][l])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 3: The ADM extrinsic curvature $K_{ij}$ and its derivatives in terms of rescaled BSSN quantities. \\[Back to [top](#toc)\\]\n",
"$$\\label{extrinsiccurvature}$$\n",
"\n",
"The ADM extrinsic curvature may be written in terms of the BSSN trace-free extrinsic curvature tensor $\\bar{A}_{ij}$ and the trace of the ADM extrinsic curvature $K$:\n",
"\n",
"\\begin{align}\n",
"K_{ij} &= \\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right)^{1/3} \\bar{A}_{ij} + \\frac{1}{3} \\gamma_{ij} K \\\\\n",
"&= e^{4\\phi} \\bar{A}_{ij} + \\frac{1}{3} \\gamma_{ij} K \\\\\n",
"\\end{align}\n",
"\n",
"We only compute first spatial derivatives of $K_{ij}$, as higher-derivatives are generally not needed:\n",
"$$\n",
"K_{ij,k} = \\left(e^{4\\phi}\\right)_{,k} \\bar{A}_{ij} + e^{4\\phi} \\bar{A}_{ij,k} + \\frac{1}{3} \\left(\\gamma_{ij,k} K + \\gamma_{ij} K_{,k}\\right)\n",
"$$\n",
"which is expressed in terms of quantities already defined."
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:41.256999Z",
"iopub.status.busy": "2021-03-07T17:06:41.246585Z",
"iopub.status.idle": "2021-03-07T17:06:41.306372Z",
"shell.execute_reply": "2021-03-07T17:06:41.305791Z"
}
},
"outputs": [],
"source": [
"# Step 3: Define ADM extrinsic curvature KDD and\n",
"# its first spatial derivatives KDDdD\n",
"# in terms of BSSN quantities\n",
"KDD = ixp.zerorank2()\n",
"\n",
"for i in range(DIM):\n",
" for j in range(DIM):\n",
" KDD[i][j] = exp4phi*AbarDD[i][j] + sp.Rational(1,3)*gammaDD[i][j]*trK\n",
"\n",
"KDDdD = ixp.zerorank3()\n",
"trK_dD = ixp.declarerank1(\"trK_dD\")\n",
"for i in range(DIM):\n",
" for j in range(DIM):\n",
" for k in range(DIM):\n",
" KDDdD[i][j][k] = exp4phidD[k]*AbarDD[i][j] + exp4phi*AbarDD_dD[i][j][k] + \\\n",
" sp.Rational(1,3)*(gammaDDdD[i][j][k]*trK + gammaDD[i][j]*trK_dD[k])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 4: Code Validation against `BSSN.ADM_in_terms_of_BSSN` NRPy+ module \\[Back to [top](#toc)\\]\n",
"$$\\label{code_validation}$$\n",
"\n",
"Here, as a code validation check, we verify agreement in the SymPy expressions between\n",
"\n",
"1. this tutorial and \n",
"2. the NRPy+ [BSSN.ADM_in_terms_of_BSSN](../edit/BSSN/ADM_in_terms_of_BSSN.py) module.\n"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:41.322765Z",
"iopub.status.busy": "2021-03-07T17:06:41.322057Z",
"iopub.status.idle": "2021-03-07T17:06:41.713780Z",
"shell.execute_reply": "2021-03-07T17:06:41.713253Z"
},
"scrolled": true
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"ALL TESTS PASSED!\n"
]
}
],
"source": [
"def comp_func(expr1,expr2,basename,prefixname2=\"Bq.\"):\n",
" if str(expr1-expr2)!=\"0\":\n",
" print(basename+\" - \"+prefixname2+basename+\" = \"+ str(expr1-expr2))\n",
" return 1\n",
" return 0\n",
"\n",
"def gfnm(basename,idx1,idx2=None,idx3=None,idx4=None):\n",
" if idx2 is None:\n",
" return basename+\"[\"+str(idx1)+\"]\"\n",
" if idx3 is None:\n",
" return basename+\"[\"+str(idx1)+\"][\"+str(idx2)+\"]\"\n",
" if idx4 is None:\n",
" return basename+\"[\"+str(idx1)+\"][\"+str(idx2)+\"][\"+str(idx3)+\"]\"\n",
" return basename+\"[\"+str(idx1)+\"][\"+str(idx2)+\"][\"+str(idx3)+\"][\"+str(idx4)+\"]\"\n",
"\n",
"expr_list = []\n",
"exprcheck_list = []\n",
"namecheck_list = []\n",
"\n",
"import BSSN.ADM_in_terms_of_BSSN as AB\n",
"AB.ADM_in_terms_of_BSSN()\n",
"\n",
"namecheck_list.extend([\"detgamma\"])\n",
"exprcheck_list.extend([AB.detgamma])\n",
"expr_list.extend([detgamma])\n",
"for i in range(DIM):\n",
" for j in range(DIM):\n",
" namecheck_list.extend([gfnm(\"gammaDD\",i,j),gfnm(\"gammaUU\",i,j),gfnm(\"KDD\",i,j)])\n",
" exprcheck_list.extend([AB.gammaDD[i][j],AB.gammaUU[i][j],AB.KDD[i][j]])\n",
" expr_list.extend([gammaDD[i][j],gammaUU[i][j],KDD[i][j]])\n",
" for k in range(DIM):\n",
" namecheck_list.extend([gfnm(\"gammaDDdD\",i,j,k),gfnm(\"GammaUDD\",i,j,k),gfnm(\"KDDdD\",i,j,k)])\n",
" exprcheck_list.extend([AB.gammaDDdD[i][j][k],AB.GammaUDD[i][j][k],AB.KDDdD[i][j][k]])\n",
" expr_list.extend([gammaDDdD[i][j][k],GammaUDD[i][j][k],KDDdD[i][j][k]])\n",
" for l in range(DIM):\n",
" namecheck_list.extend([gfnm(\"gammaDDdDD\",i,j,k,l)])\n",
" exprcheck_list.extend([AB.gammaDDdDD[i][j][k][l]])\n",
" expr_list.extend([gammaDDdDD[i][j][k][l]])\n",
"\n",
"num_failures = 0\n",
"for i in range(len(expr_list)):\n",
" num_failures += comp_func(expr_list[i],exprcheck_list[i],namecheck_list[i])\n",
"\n",
"\n",
"if num_failures == 0:\n",
" print(\"ALL TESTS PASSED!\")\n",
"else:\n",
" print(\"ERROR. \" + str(num_failures) + \" TESTS FAILED\")\n",
" sys.exit(1)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 5: Output this notebook to $\\LaTeX$-formatted PDF file \\[Back to [top](#toc)\\]\n",
"$$\\label{latex_pdf_output}$$\n",
"\n",
"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\n",
"[Tutorial-ADM_in_terms_of_BSSN.pdf](Tutorial-ADM_in_terms_of_BSSN.pdf) (Note that clicking on this link may not work; you may need to open the PDF file through another means.)"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"execution": {
"iopub.execute_input": "2021-03-07T17:06:41.718813Z",
"iopub.status.busy": "2021-03-07T17:06:41.717758Z",
"iopub.status.idle": "2021-03-07T17:06:45.208343Z",
"shell.execute_reply": "2021-03-07T17:06:45.208956Z"
},
"scrolled": true
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Created Tutorial-ADM_in_terms_of_BSSN.tex, and compiled LaTeX file to PDF\n",
" file Tutorial-ADM_in_terms_of_BSSN.pdf\n"
]
}
],
"source": [
"import cmdline_helper as cmd # NRPy+: Multi-platform Python command-line interface\n",
"cmd.output_Jupyter_notebook_to_LaTeXed_PDF(\"Tutorial-ADM_in_terms_of_BSSN\")"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.1"
}
},
"nbformat": 4,
"nbformat_minor": 4
}