{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "\n", "# BSSN Quantities in terms of ADM Quantities\n", "## Author: Zach Etienne\n", "\n", "## This notebook showcases the conversion of ADM to BSSN variables using NRPy+, preparing them for solving Einstein's equations with the BSSN formulation. It covers the rescaling of 3-metric, extrinsic curvature, and gauge quantities, with an outline of the rescaling processes from covariant BSSN formulation and the underlying references.\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: [BSSN_in_terms_of_ADM.py](../edit/BSSN/BSSN_in_terms_of_ADM.py)\n", "\n", "## Introduction:\n", "This module documents the conversion of ADM variables:\n", "\n", "$$\\left\\{\\gamma_{ij}, K_{ij}, \\alpha, \\beta^i\\right\\}$$\n", "\n", "into BSSN variables\n", "\n", "$$\\left\\{\\bar{\\gamma}_{i j},\\bar{A}_{i j},\\phi, K, \\bar{\\Lambda}^{i}, \\alpha, \\beta^i, B^i\\right\\},$$ \n", "\n", "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](Tutorial-BSSN_formulation.ipynb)) into the form needed for solving Einstein's equations with the BSSN formulation:\n", "\n", "$$\\left\\{h_{i j},a_{i j},\\phi, K, \\lambda^{i}, \\alpha, \\mathcal{V}^i, \\mathcal{B}^i\\right\\}.$$" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "# 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; set desired output BSSN Curvilinear coordinate system set to Spherical\n", "1. [Step 2](#adm2bssn): Perform the ADM-to-BSSN conversion for 3-metric, extrinsic curvature, and gauge quantities\n", " 1. [Step 2.a](#adm2bssn_gamma): Convert ADM $\\gamma_{ij}$ to BSSN $\\bar{\\gamma}_{ij}$; rescale to get $h_{ij}$\n", " 1. [Step 2.b](#admexcurv_convert): Convert the ADM extrinsic curvature $K_{ij}$ to BSSN $\\bar{A}_{ij}$ and $K$; rescale to get $a_{ij}$, $K$.\n", " 1. [Step 2.c](#lambda): Define $\\bar{\\Lambda}^i$\n", " 1. [Step 2.d](#conformal): Define the conformal factor variable `cf`\n", "1. [Step 3](#code_validation): Code Validation against `BSSN.BSSN_in_terms_of_ADM` NRPy+ module\n", "1. [Step 4](#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" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:21:27.103925Z", "iopub.status.busy": "2021-03-07T17:21:27.103154Z", "iopub.status.idle": "2021-03-07T17:21:27.427720Z", "shell.execute_reply": "2021-03-07T17:21:27.428217Z" } }, "outputs": [], "source": [ "# Step 1: Import needed core NRPy+ modules\n", "import sympy as sp # SymPy: The Python computer algebra package upon which NRPy+ depends\n", "import NRPy_param_funcs as par # NRPy+: Parameter interface\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 modules for multiplatform OS-level functions\n", "import BSSN.BSSN_quantities as Bq # NRPy+: This module depends on the parameter EvolvedConformalFactor_cf,\n", " # which is defined in BSSN.BSSN_quantities\n", "\n", "# Step 1.a: Set DIM=3, as we're using a 3+1 decomposition of Einstein's equations\n", "DIM=3" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "# Step 2: Perform the ADM-to-BSSN conversion for 3-metric, extrinsic curvature, and gauge quantities \\[Back to [top](#toc)\\]\n", "$$\\label{adm2bssn}$$\n", "\n", "Here we convert ADM quantities to their BSSN Curvilinear counterparts." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "## Step 2.a: Convert ADM $\\gamma_{ij}$ to BSSN $\\bar{\\gamma}_{ij}$; rescale to get $h_{ij}$ \\[Back to [top](#toc)\\]\n", "$$\\label{adm2bssn_gamma}$$\n", "\n", "We have (Eqs. 2 and 3 of [Ruchlin *et al.*](https://arxiv.org/pdf/1712.07658.pdf)):\n", "$$\n", "\\bar{\\gamma}_{i j} = \\left(\\frac{\\bar{\\gamma}}{\\gamma}\\right)^{1/3} \\gamma_{ij},\n", "$$\n", "where we always make the choice $\\bar{\\gamma} = \\hat{\\gamma}$.\n", "\n", "After constructing $\\bar{\\gamma}_{ij}$, we rescale to get $h_{ij}$ according to the prescription described in the [the covariant BSSN formulation tutorial](Tutorial-BSSN_formulation.ipynb) (also [Ruchlin *et al.*](https://arxiv.org/pdf/1712.07658.pdf)):\n", "\n", "$$\n", "h_{ij} = (\\bar{\\gamma}_{ij} - \\hat{\\gamma}_{ij})/\\text{ReDD[i][j]}.\n", "$$" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:21:27.436089Z", "iopub.status.busy": "2021-03-07T17:21:27.435407Z", "iopub.status.idle": "2021-03-07T17:21:27.437595Z", "shell.execute_reply": "2021-03-07T17:21:27.438163Z" } }, "outputs": [], "source": [ "# Step 2: All ADM quantities were input into this function in the Spherical or Cartesian\n", "# basis, as functions of r,th,ph or x,y,z, respectively. In Steps 1 and 2 above,\n", "# we converted them to the xx0,xx1,xx2 basis, and as functions of xx0,xx1,xx2.\n", "# Here we convert ADM quantities to their BSSN Curvilinear counterparts:\n", "\n", "# Step 2.a: Convert ADM $\\gamma_{ij}$ to BSSN $\\bar{gamma}_{ij}$:\n", "# We have (Eqs. 2 and 3 of [Ruchlin *et al.*](https://arxiv.org/pdf/1712.07658.pdf)):\n", "def gammabarDD_hDD(gammaDD):\n", " global gammabarDD,hDD\n", " if rfm.have_already_called_reference_metric_function == False:\n", " print(\"BSSN.BSSN_in_terms_of_ADM.hDD_given_ADM(): Must call reference_metric() first!\")\n", " sys.exit(1)\n", " # \\bar{gamma}_{ij} = (\\frac{\\bar{gamma}}{gamma})^{1/3}*gamma_{ij}.\n", " _gammaUU, gammaDET = ixp.symm_matrix_inverter3x3(gammaDD)\n", " gammabarDD = ixp.zerorank2()\n", " hDD = ixp.zerorank2()\n", " for i in range(DIM):\n", " for j in range(DIM):\n", " gammabarDD[i][j] = (rfm.detgammahat/gammaDET)**(sp.Rational(1,3))*gammaDD[i][j]\n", " hDD[i][j] = (gammabarDD[i][j] - rfm.ghatDD[i][j]) / rfm.ReDD[i][j]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "## Step 2.b: Convert the ADM extrinsic curvature $K_{ij}$ to BSSN quantities $\\bar{A}_{ij}$ and $K={\\rm tr}(K_{ij})$; rescale $\\bar{A}_{ij}$ to get $a_{ij}$ \\[Back to [top](#toc)\\]\n", "$$\\label{admexcurv_convert}$$\n", "\n", "Convert the ADM extrinsic curvature $K_{ij}$ to the trace-free extrinsic curvature $\\bar{A}_{ij}$, plus the trace of the extrinsic curvature $K$, where (Eq. 3 of [Baumgarte *et al.*](https://arxiv.org/pdf/1211.6632.pdf)):\n", "\\begin{align}\n", "K &= \\gamma^{ij} K_{ij} \\\\\n", "\\bar{A}_{ij} &= \\left(\\frac{\\bar{\\gamma}}{\\gamma}\\right)^{1/3} \\left(K_{ij} - \\frac{1}{3} \\gamma_{ij} K \\right)\n", "\\end{align}\n", "\n", "After constructing $\\bar{A}_{ij}$, we rescale to get $a_{ij}$ according to the prescription described in the [the covariant BSSN formulation tutorial](Tutorial-BSSN_formulation.ipynb) (also [Ruchlin *et al.*](https://arxiv.org/pdf/1712.07658.pdf)):\n", "\n", "$$\n", "a_{ij} = \\bar{A}_{ij}/\\text{ReDD[i][j]}.\n", "$$" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:21:27.447400Z", "iopub.status.busy": "2021-03-07T17:21:27.446610Z", "iopub.status.idle": "2021-03-07T17:21:27.448700Z", "shell.execute_reply": "2021-03-07T17:21:27.449265Z" } }, "outputs": [], "source": [ "# Step 2.b: Convert the extrinsic curvature K_{ij} to the trace-free extrinsic\n", "# curvature \\bar{A}_{ij}, plus the trace of the extrinsic curvature K,\n", "# where (Eq. 3 of [Baumgarte *et al.*](https://arxiv.org/pdf/1211.6632.pdf)):\n", "def trK_AbarDD_aDD(gammaDD,KDD):\n", " global trK,AbarDD,aDD\n", " if rfm.have_already_called_reference_metric_function == False:\n", " print(\"BSSN.BSSN_in_terms_of_ADM.trK_AbarDD(): Must call reference_metric() first!\")\n", " sys.exit(1)\n", " # \\bar{gamma}_{ij} = (\\frac{\\bar{gamma}}{gamma})^{1/3}*gamma_{ij}.\n", " gammaUU, gammaDET = ixp.symm_matrix_inverter3x3(gammaDD)\n", " # K = gamma^{ij} K_{ij}, and\n", " # \\bar{A}_{ij} &= (\\frac{\\bar{gamma}}{gamma})^{1/3}*(K_{ij} - \\frac{1}{3}*gamma_{ij}*K)\n", " trK = sp.sympify(0)\n", " for i in range(DIM):\n", " for j in range(DIM):\n", " trK += gammaUU[i][j]*KDD[i][j]\n", "\n", " AbarDD = ixp.zerorank2()\n", " aDD = ixp.zerorank2()\n", " for i in range(DIM):\n", " for j in range(DIM):\n", " AbarDD[i][j] = (rfm.detgammahat/gammaDET)**(sp.Rational(1,3))*(KDD[i][j] - sp.Rational(1,3)*gammaDD[i][j]*trK)\n", " aDD[i][j] = AbarDD[i][j] / rfm.ReDD[i][j]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "## Step 2.c: Assuming the ADM 3-metric $\\gamma_{ij}$ is given as an explicit function of `(xx0,xx1,xx2)`, convert to BSSN $\\bar{\\Lambda}^i$; rescale to compute $\\lambda^i$ \\[Back to [top](#toc)\\]\n", "$$\\label{lambda}$$\n", "\n", "To define $\\bar{\\Lambda}^i$ we implement Eqs. 4 and 5 of [Baumgarte *et al.*](https://arxiv.org/pdf/1211.6632.pdf):\n", "$$\n", "\\bar{\\Lambda}^i = \\bar{\\gamma}^{jk}\\left(\\bar{\\Gamma}^i_{jk} - \\hat{\\Gamma}^i_{jk}\\right).\n", "$$\n", "\n", "The [reference_metric.py](../edit/reference_metric.py) module provides us with exact, closed-form expressions for $\\hat{\\Gamma}^i_{jk}$, so here we need only compute exact expressions for $\\bar{\\Gamma}^i_{jk}$, based on $\\gamma_{ij}$ given as an explicit function of `(xx0,xx1,xx2)`. This is particularly useful when setting up initial data.\n", "\n", "After constructing $\\bar{\\Lambda}^i$, we rescale to get $\\lambda^i$ according to the prescription described in the [the covariant BSSN formulation tutorial](Tutorial-BSSN_formulation.ipynb) (also [Ruchlin *et al.*](https://arxiv.org/pdf/1712.07658.pdf)):\n", "\n", "$$\n", "\\lambda^i = \\bar{\\Lambda}^i/\\text{ReU[i]}.\n", "$$" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:21:27.459613Z", "iopub.status.busy": "2021-03-07T17:21:27.458929Z", "iopub.status.idle": "2021-03-07T17:21:27.461197Z", "shell.execute_reply": "2021-03-07T17:21:27.461662Z" } }, "outputs": [], "source": [ "# Step 2.c: Define \\bar{Lambda}^i (Eqs. 4 and 5 of [Baumgarte *et al.*](https://arxiv.org/pdf/1211.6632.pdf)):\n", "def LambdabarU_lambdaU__exact_gammaDD(gammaDD):\n", " global LambdabarU,lambdaU\n", "\n", " # \\bar{Lambda}^i = \\bar{gamma}^{jk}(\\bar{Gamma}^i_{jk} - \\hat{Gamma}^i_{jk}).\n", " gammabarDD_hDD(gammaDD)\n", " gammabarUU, _gammabarDET = ixp.symm_matrix_inverter3x3(gammabarDD)\n", "\n", " # First compute Christoffel symbols \\bar{Gamma}^i_{jk}, with respect to barred metric:\n", " GammabarUDD = ixp.zerorank3()\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", " GammabarUDD[i][j][k] += sp.Rational(1,2)*gammabarUU[i][l]*( sp.diff(gammabarDD[l][j],rfm.xx[k]) +\n", " sp.diff(gammabarDD[l][k],rfm.xx[j]) -\n", " sp.diff(gammabarDD[j][k],rfm.xx[l]) )\n", " # Next evaluate \\bar{Lambda}^i, based on GammabarUDD above and GammahatUDD\n", " # (from the reference metric):\n", " LambdabarU = ixp.zerorank1()\n", " for i in range(DIM):\n", " for j in range(DIM):\n", " for k in range(DIM):\n", " LambdabarU[i] += gammabarUU[j][k] * (GammabarUDD[i][j][k] - rfm.GammahatUDD[i][j][k])\n", " for i in range(DIM):\n", " # We evaluate LambdabarU[i] here to ensure proper cancellations. If these cancellations\n", " # are not applied, certain expressions (e.g., lambdaU[0] in StaticTrumpet) will\n", " # cause SymPy's (v1.5+) CSE algorithm to hang\n", " LambdabarU[i] = LambdabarU[i].doit()\n", " lambdaU = ixp.zerorank1()\n", " for i in range(DIM):\n", " lambdaU[i] = LambdabarU[i] / rfm.ReU[i]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "## Step 2.d: Define the conformal factor variable `cf` \\[Back to [top](#toc)\\]\n", "$$\\label{conformal}$$\n", "\n", "We define the conformal factor variable `cf` based on the setting of the `\"BSSN_quantities::EvolvedConformalFactor_cf\"` parameter.\n", "\n", "For example if `\"BSSN_quantities::EvolvedConformalFactor_cf\"` is set to `\"phi\"`, we can use Eq. 3 of [Ruchlin *et al.*](https://arxiv.org/pdf/1712.07658.pdf), which in arbitrary coordinates is written:\n", "\n", "$$\n", "\\phi = \\frac{1}{12} \\log\\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right).\n", "$$\n", "\n", "Alternatively if `\"BSSN_quantities::EvolvedConformalFactor_cf\"` is set to `\"chi\"`, then\n", "$$\n", "\\chi = e^{-4 \\phi} = \\exp\\left(-4 \\frac{1}{12} \\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right)\\right) \n", "= \\exp\\left(-\\frac{1}{3} \\log\\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right)\\right) = \\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right)^{-1/3}.\n", "$$\n", "\n", "Finally if `\"BSSN_quantities::EvolvedConformalFactor_cf\"` is set to `\"W\"`, then\n", "$$\n", "W = e^{-2 \\phi} = \\exp\\left(-2 \\frac{1}{12} \\log\\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right)\\right) = \n", "\\exp\\left(-\\frac{1}{6} \\log\\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right)\\right) = \n", "\\left(\\frac{\\gamma}{\\bar{\\gamma}}\\right)^{-1/6}.\n", "$$" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:21:27.469865Z", "iopub.status.busy": "2021-03-07T17:21:27.468755Z", "iopub.status.idle": "2021-03-07T17:21:27.471869Z", "shell.execute_reply": "2021-03-07T17:21:27.472581Z" } }, "outputs": [], "source": [ "# Step 2.d: Set the conformal factor variable cf, which is set\n", "# by the \"BSSN_quantities::EvolvedConformalFactor_cf\" parameter. For example if\n", "# \"EvolvedConformalFactor_cf\" is set to \"phi\", we can use Eq. 3 of\n", "# [Ruchlin *et al.*](https://arxiv.org/pdf/1712.07658.pdf),\n", "# which in arbitrary coordinates is written:\n", "def cf_from_gammaDD(gammaDD):\n", " global cf\n", "\n", " # \\bar{Lambda}^i = \\bar{gamma}^{jk}(\\bar{Gamma}^i_{jk} - \\hat{Gamma}^i_{jk}).\n", " gammabarDD_hDD(gammaDD)\n", " _gammabarUU, gammabarDET = ixp.symm_matrix_inverter3x3(gammabarDD)\n", " _gammaUU, gammaDET = ixp.symm_matrix_inverter3x3(gammaDD)\n", "\n", " cf = sp.sympify(0)\n", "\n", " if par.parval_from_str(\"EvolvedConformalFactor_cf\") == \"phi\":\n", " # phi = \\frac{1}{12} log(\\frac{gamma}{\\bar{gamma}}).\n", " cf = sp.Rational(1,12)*sp.log(gammaDET/gammabarDET)\n", " elif par.parval_from_str(\"EvolvedConformalFactor_cf\") == \"chi\":\n", " # chi = exp(-4*phi) = exp(-4*\\frac{1}{12}*(\\frac{gamma}{\\bar{gamma}}))\n", " # = exp(-\\frac{1}{3}*log(\\frac{gamma}{\\bar{gamma}})) = (\\frac{gamma}{\\bar{gamma}})^{-1/3}.\n", " #\n", " cf = (gammaDET/gammabarDET)**(-sp.Rational(1,3))\n", " elif par.parval_from_str(\"EvolvedConformalFactor_cf\") == \"W\":\n", " # W = exp(-2*phi) = exp(-2*\\frac{1}{12}*log(\\frac{gamma}{\\bar{gamma}}))\n", " # = exp(-\\frac{1}{6}*log(\\frac{gamma}{\\bar{gamma}})) = (\\frac{gamma}{bar{gamma}})^{-1/6}.\n", " cf = (gammaDET/gammabarDET)**(-sp.Rational(1,6))\n", " else:\n", " print(\"Error EvolvedConformalFactor_cf type = \\\"\"+par.parval_from_str(\"EvolvedConformalFactor_cf\")+\"\\\" unknown.\")\n", " sys.exit(1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "## Step 2.e: Rescale $\\beta^i$ and $B^i$ to compute $\\mathcal{V}^i={\\rm vet}^i$ and $\\mathcal{B}^i={\\rm bet}^i$, respectively \\[Back to [top](#toc)\\]\n", "$$\\label{betvet}$$\n", "\n", "We rescale $\\beta^i$ and $B^i$ according to the prescription described in the [the covariant BSSN formulation tutorial](Tutorial-BSSN_formulation.ipynb) (also [Ruchlin *et al.*](https://arxiv.org/pdf/1712.07658.pdf)):\n", "\\begin{align}\n", "\\mathcal{V}^i &= \\beta^i/\\text{ReU[i]}\\\\\n", "\\mathcal{B}^i &= B^i/\\text{ReU[i]}.\n", "\\end{align}" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:21:27.480259Z", "iopub.status.busy": "2021-03-07T17:21:27.479612Z", "iopub.status.idle": "2021-03-07T17:21:27.481888Z", "shell.execute_reply": "2021-03-07T17:21:27.482389Z" } }, "outputs": [], "source": [ "# Step 2.e: Rescale beta^i and B^i according to the prescription described in\n", "# the [BSSN in curvilinear coordinates tutorial notebook](Tutorial-BSSNCurvilinear.ipynb)\n", "# (also [Ruchlin *et al.*](https://arxiv.org/pdf/1712.07658.pdf)):\n", "#\n", "# \\mathcal{V}^i &= beta^i/(ReU[i])\n", "# \\mathcal{B}^i &= B^i/(ReU[i])\n", "def betU_vetU(betaU,BU):\n", " global vetU,betU\n", "\n", " if rfm.have_already_called_reference_metric_function == False:\n", " print(\"BSSN.BSSN_in_terms_of_ADM.bet_vet(): Must call reference_metric() first!\")\n", " sys.exit(1)\n", " vetU = ixp.zerorank1()\n", " betU = ixp.zerorank1()\n", " for i in range(DIM):\n", " vetU[i] = betaU[i] / rfm.ReU[i]\n", " betU[i] = BU[i] / rfm.ReU[i]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "# Step 3: Code Validation against `BSSN.BSSN_in_terms_of_ADM` module \\[Back to [top](#toc)\\] \n", "$$\\label{code_validation}$$\n", "\n", "Here, as a code validation check, we verify agreement in the SymPy expressions for [UIUC initial data](Tutorial-ADM_Initial_Data-UIUC_BlackHole.ipynb) between\n", "1. this tutorial and \n", "2. the NRPy+ [BSSN.BSSN_in_terms_of_ADM](../edit/BSSN/BSSN_in_terms_of_ADM.py) module.\n", "\n", "As no basis transformation is performed, we analyze these expressions in their native, Spherical coordinates." ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:21:27.497813Z", "iopub.status.busy": "2021-03-07T17:21:27.497089Z", "iopub.status.idle": "2021-03-07T17:21:28.240964Z", "shell.execute_reply": "2021-03-07T17:21:28.241439Z" } }, "outputs": [], "source": [ "# Step 3.a: Set the desired *output* coordinate system to Spherical:\n", "par.set_parval_from_str(\"reference_metric::CoordSystem\",\"Spherical\")\n", "rfm.reference_metric()\n", "\n", "# Step 3.b: Set up initial data; assume UIUC spinning black hole initial data\n", "import BSSN.UIUCBlackHole as uibh\n", "uibh.UIUCBlackHole()\n", "\n", "# Step 3.c: Call above functions to convert ADM to BSSN curvilinear\n", "gammabarDD_hDD( uibh.gammaDD)\n", "trK_AbarDD_aDD( uibh.gammaDD,uibh.KDD)\n", "LambdabarU_lambdaU__exact_gammaDD(uibh.gammaDD)\n", "cf_from_gammaDD( uibh.gammaDD)\n", "betU_vetU( uibh.betaU,uibh.BU)\n", "\n", "# Step 3.d: Now load the BSSN_in_terms_of_ADM module and perform the same conversion\n", "import BSSN.BSSN_in_terms_of_ADM as BitoA\n", "BitoA.gammabarDD_hDD( uibh.gammaDD)\n", "BitoA.trK_AbarDD_aDD( uibh.gammaDD,uibh.KDD)\n", "BitoA.LambdabarU_lambdaU__exact_gammaDD(uibh.gammaDD)\n", "BitoA.cf_from_gammaDD( uibh.gammaDD)\n", "BitoA.betU_vetU( uibh.betaU,uibh.BU)" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:21:27.497813Z", "iopub.status.busy": "2021-03-07T17:21:27.497089Z", "iopub.status.idle": "2021-03-07T17:21:28.240964Z", "shell.execute_reply": "2021-03-07T17:21:28.241439Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Consistency check between this tutorial notebook and BSSN.BSSN_in_terms_of_ADM NRPy+ module:\n", "ALL TESTS PASS\n" ] } ], "source": [ "# Step 3.e: Perform the consistency check\n", "def compare(q1, q2, q1name, q2name):\n", " if sp.simplify(q1 - q2) != 0:\n", " print(\"Error: \"+q1name+\" - \"+q2name+\" = \"+str(sp.simplify(q1 - q2)))\n", " sys.exit(1) \n", "\n", "print(\"Consistency check between this tutorial notebook and BSSN.BSSN_in_terms_of_ADM NRPy+ module:\")\n", "\n", "compare(cf, BitoA.cf, \"cf\", \"BitoA.cf\")\n", "compare(trK, BitoA.trK, \"trK\", \"BitoA.trK\")\n", "# alpha is the only variable that remains unchanged:\n", "# print(\"alpha - BitoA.alpha = \" + str(alpha - BitoA.alpha))\n", "for i in range(DIM):\n", " compare(vetU[i], BitoA.vetU[i], \"vetU\"+str(i), \"BitoA.vetU\"+str(i))\n", " compare(betU[i], BitoA.betU[i], \"betU\"+str(i), \"BitoA.betU\"+str(i))\n", " compare(lambdaU[i], BitoA.lambdaU[i], \"lambdaU\"+str(i), \"BitoA.lambdaU\"+str(i))\n", " for j in range(DIM):\n", " compare(hDD[i][j], BitoA.hDD[i][j], \"hDD\"+str(i)+str(j), \"BitoA.hDD\"+str(i)+str(j))\n", " compare(aDD[i][j], BitoA.aDD[i][j], \"aDD\"+str(i)+str(j), \"BitoA.aDD\"+str(i)+str(j))\n", "\n", "print(\"ALL TESTS PASS\")" ] }, { "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 [Tutorial-BSSN_in_terms_of_ADM.pdf](Tutorial-BSSN_in_terms_of_ADM.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": 9, "metadata": { "execution": { "iopub.execute_input": "2021-03-07T17:21:28.246581Z", "iopub.status.busy": "2021-03-07T17:21:28.245547Z", "iopub.status.idle": "2021-03-07T17:21:31.738562Z", "shell.execute_reply": "2021-03-07T17:21:31.739462Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Created Tutorial-BSSN_in_terms_of_ADM.tex, and compiled LaTeX file to PDF\n", " file Tutorial-BSSN_in_terms_of_ADM.pdf\n" ] } ], "source": [ "import cmdline_helper as cmd # NRPy+: Multi-platform Python command-line interface\n", "cmd.output_Jupyter_notebook_to_LaTeXed_PDF(\"Tutorial-BSSN_in_terms_of_ADM\")" ] } ], "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 }