{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"\n",
"# Tutorial-IllinoisGRMHD: `convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C`\n",
"\n",
"## Authors: Leo Werneck & Zach Etienne\n",
"\n",
"**This module is currently under development**\n",
"\n",
"## This tutorial module explains the algorithm used to get the BSSN variables from the ADM variables. Further, it outlines how to impose the constraint that the conformal metric has a unit determinant.\n",
"\n",
"### Required and recommended citations:\n",
"\n",
"* **(Required)** Etienne, Z. B., Paschalidis, V., Haas R., Mösta P., and Shapiro, S. L. IllinoisGRMHD: an open-source, user-friendly GRMHD code for dynamical spacetimes. Class. Quantum Grav. 32 (2015) 175009. ([arxiv:1501.07276](http://arxiv.org/abs/1501.07276)).\n",
"* **(Required)** Noble, S. C., Gammie, C. F., McKinney, J. C., Del Zanna, L. Primitive Variable Solvers for Conservative General Relativistic Magnetohydrodynamics. Astrophysical Journal, 641, 626 (2006) ([astro-ph/0512420](https://arxiv.org/abs/astro-ph/0512420)).\n",
"* **(Recommended)** Del Zanna, L., Bucciantini N., Londrillo, P. An efficient shock-capturing central-type scheme for multidimensional relativistic flows - II. Magnetohydrodynamics. A&A 400 (2) 397-413 (2003). DOI: 10.1051/0004-6361:20021641 ([astro-ph/0210618](https://arxiv.org/abs/astro-ph/0210618))."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Table of Contents\n",
"$$\\label{toc}$$\n",
"\n",
"This module is organized as follows\n",
"\n",
"0. [Step 0](#src_dir): **Source directory creation**\n",
"1. [Step 1](#introduction): **Introduction**\n",
"1. [Step 2](#convert_adm_to_bssn__det_gammabar_eq_1): **`convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C`**\n",
" 1. [Step 2.a](#physical_metric_and_its_determinant): *Loading the physical metric $\\gamma_{ij}$ and computing $\\gamma = \\det\\left(\\gamma_{ij}\\right)$*\n",
" 1. [Step 2.b](#phi_and_psi): *Computing $\\phi$ and $\\psi$*\n",
" 1. [Step 2.c](#enforce_det_gammabar_eq_1): *Enforcing the constraint $\\bar\\gamma = 1$*\n",
" 1. [Step 2.d](#update_gfs): *Updating the gridfunctions*\n",
"1. [Step 3](#code_validation): **Code validation**\n",
"1. [Step 4](#latex_pdf_output): **Output this notebook to $\\LaTeX$-formatted PDF file**"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 0: Source directory creation \\[Back to [top](#toc)\\]\n",
"$$\\label{src_dir}$$\n",
"\n",
"We will now use the [cmdline_helper.py NRPy+ module](Tutorial-Tutorial-cmdline_helper.ipynb) to create the source directory within the `IllinoisGRMHD` NRPy+ directory if it does not exist yet."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"# Step 0: Creation of the IllinoisGRMHD source directory\n",
"# Step 0a: Add NRPy's directory to the path\n",
"# https://stackoverflow.com/questions/16780014/import-file-from-parent-directory\n",
"import os,sys\n",
"nrpy_dir_path = os.path.join(\"..\",\"..\")\n",
"if nrpy_dir_path not in sys.path:\n",
" sys.path.append(nrpy_dir_path)\n",
"\n",
"# Step 0b: Load up cmdline_helper and create the directory\n",
"import cmdline_helper as cmd\n",
"IGM_src_dir_path = os.path.join(\"..\",\"src\")\n",
"cmd.mkdir(IGM_src_dir_path)\n",
"\n",
"# Step 0c: Create the output file path\n",
"outfile_path__ADM_to_BSSN__det_gammabar_eq_1__C = os.path.join(IGM_src_dir_path,\"convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 1: Introduction \\[Back to [top](#toc)\\]\n",
"$$\\label{introduction}$$\n",
"\n",
"In this module, we will explain the procedure used within `IllinoisGRMHD` to compute the conformal metric, $\\bar\\gamma_{ij}$, and its inverse, $\\bar\\gamma^{ij}$, from the physical metric $\\gamma_{ij}$. We will also describe how to compute $\\phi$, the conformal factor, and $\\psi\\equiv e^{\\phi}$. Finally, we also explain the procedure used to enforce the constraint $\\bar\\gamma = \\det\\left(\\bar\\gamma_{ij}\\right) = 1$.\n",
"\n",
"**A note on notation**: The notation used throughout the NRPy tutorial notebooks for the conformal metric is $\\bar\\gamma_{ij}$. However, in the literature, the notation $\\tilde\\gamma_{ij}$ is also used to represent the conformal metric. It is important to note that in `IllinoisGRMHD` we refer to the conformal metric using the *latter* notation, i.e. ${\\rm gtij} := \\tilde\\gamma_{ij}$ and ${\\rm gtupij} := \\tilde\\gamma^{ij}$. To keep the discussion consistent with the other notebooks, however, we will still present the discussion using the notation $\\tilde\\gamma_{ij} \\to \\bar\\gamma_{ij}$. Bottom line, here $\\tilde\\gamma_{ij}$ and $\\bar\\gamma_{ij}$ represent exactly the same quantity, the conformal metric."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 2: `convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C` \\[Back to [top](#toc)\\]\n",
"$$\\label{convert_adm_to_bssn__det_gammabar_eq_1}$$\n",
"\n",
"\n",
"\n",
"## Step 2.a: Loading the physical metric $\\gamma_{ij}$ and computing $\\gamma = \\det\\left(\\gamma_{ij}\\right)$ \\[Back to [top](#toc)\\]\n",
"$$\\label{physical_metric_and_its_determinant}$$\n",
"\n",
"We start by reading in the physical metric ${\\rm gij\\_physL} := \\gamma_{ij}$ and computing its determinant\n",
"\n",
"$$\n",
"\\boxed{\n",
"\\begin{align}\n",
"\\gamma = \\det\\left(\\gamma_{ij}\\right)\n",
"&= \\gamma_{xx}\\gamma_{yy}\\gamma_{zz} + \\gamma_{xy}\\gamma_{yz}\\gamma_{xz} + \\gamma_{xz}\\gamma_{xy}\\gamma_{yz}\\\\\n",
"&- \\gamma_{xz}\\gamma_{yy}\\gamma_{xz} - \\gamma_{xy}\\gamma_{xy}\\gamma_{zz} - \\gamma_{xx}\\gamma_{yz}\\gamma_{yz}\n",
"\\end{align}\n",
"}\\ .\n",
"$$\n",
"\n",
"Notice that we have used the fact that $\\gamma_{ij}$ is symmetric above, i.e. $\\gamma_{ij}=\\gamma_{ji}$."
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Writing ../src/convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C\n"
]
}
],
"source": [
"%%writefile $outfile_path__ADM_to_BSSN__det_gammabar_eq_1__C\n",
"#include \n",
"#include \n",
"#include \n",
"#ifndef ENABLE_STANDALONE_IGM_C2P_SOLVER\n",
"#include \"cctk.h\"\n",
"#include \"cctk_Parameters.h\"\n",
"#endif\n",
"\n",
"void IllinoisGRMHD_convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij\n",
"(const cGH *cctkGH,const int *cctk_lsh,\n",
" CCTK_REAL *gxx,CCTK_REAL *gxy,CCTK_REAL *gxz,CCTK_REAL *gyy,CCTK_REAL *gyz,CCTK_REAL *gzz,CCTK_REAL *alp,\n",
" CCTK_REAL *gtxx,CCTK_REAL *gtxy,CCTK_REAL *gtxz,CCTK_REAL *gtyy,CCTK_REAL *gtyz,CCTK_REAL *gtzz,\n",
" CCTK_REAL *gtupxx,CCTK_REAL *gtupxy,CCTK_REAL *gtupxz,CCTK_REAL *gtupyy,CCTK_REAL *gtupyz,CCTK_REAL *gtupzz,\n",
" CCTK_REAL *phi,CCTK_REAL *psi,CCTK_REAL *lapm1) {\n",
"\n",
"#ifndef ENABLE_STANDALONE_IGM_C2P_SOLVER\n",
" DECLARE_CCTK_PARAMETERS;\n",
"#endif\n",
"\n",
"\n",
"#pragma omp parallel for\n",
" for(int k=0;k\n",
"\n",
"## Step 2.b: Computing $\\phi$ and $\\psi$ \\[Back to [top](#toc)\\]\n",
"$$\\label{phi_and_psi}$$\n",
"\n",
"The BSSN formalism relates the physical metric, $\\gamma_{ij}$, to a conformal metric, $\\bar\\gamma{ij}$, via the relation\n",
"\n",
"$$\n",
"\\boxed{\\bar\\gamma_{ij} = e^{-4\\phi}\\gamma_{ij}}\\ ,\n",
"$$\n",
"\n",
"with the constraint that $\\bar\\gamma \\equiv \\det\\left(\\bar\\gamma_{ij}\\right) = 1$. This immediately implies that\n",
"\n",
"$$\n",
"e^{-12\\phi} \\gamma = \\bar\\gamma = 1 \\implies -12\\phi = \\log\\left(\\frac{1}{\\gamma}\\right) = - \\log\\gamma\n",
"\\implies \\boxed{\\phi = \\frac{1}{12}\\log\\gamma}\\ .\n",
"$$\n",
"\n",
"Useful quantities to compute are\n",
"\n",
"$$\n",
"\\boxed{\n",
"\\begin{align}\n",
"\\psi &\\equiv e^{\\phi}\\\\\n",
"\\psi^{-4} &\\equiv e^{-4\\phi}\n",
"\\end{align}\n",
"}\\ .\n",
"$$\n",
"\n",
"All the boxed quantities above are evaluated below."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Appending to ../src/convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C\n"
]
}
],
"source": [
"%%writefile -a $outfile_path__ADM_to_BSSN__det_gammabar_eq_1__C\n",
"\n",
"\n",
" CCTK_REAL phiL = (1.0/12.0) * log(gijdet);\n",
" CCTK_REAL psiL = exp(phiL);\n",
"\n",
" CCTK_REAL Psim4 = 1.0/(psiL*psiL*psiL*psiL);\n",
" CCTK_REAL gtxxL = gxx_physL*Psim4;\n",
" CCTK_REAL gtxyL = gxy_physL*Psim4;\n",
" CCTK_REAL gtxzL = gxz_physL*Psim4;\n",
" CCTK_REAL gtyyL = gyy_physL*Psim4;\n",
" CCTK_REAL gtyzL = gyz_physL*Psim4;\n",
" CCTK_REAL gtzzL = gzz_physL*Psim4;"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"## Step 2.c: Enforcing the constraint $\\bar\\gamma = 1$ \\[Back to [top](#toc)\\]\n",
"$$\\label{enforce_det_gammabar_eq_1}$$\n",
"\n",
"The BSSN formalism evolves the determinant of the conformal metric through $\\partial_{t}\\bar\\gamma = 0$. Since initially $\\bar\\gamma=1$, one would expect that $\\bar\\gamma=1$ at all times. However, numerical errors can cause the value of the determinant of the conformal metric to drift away from unity. To remedy this, we enforce the constraint $\\bar\\gamma=1$ by transforming the conformal metric according to the algebraic condition (cf. eq. 53 of [Ruchlin, Etienne, and Baumgarte (2018)](https://arxiv.org/pdf/1712.07658.pdf) with $\\hat\\gamma = 1$).\n",
"\n",
"$$\n",
"\\boxed{\\bar{\\gamma}_{ij} \\to \\left(\\frac{1}{\\bar{\\gamma}}\\right)^{1/3} \\bar{\\gamma}_{ij}}\\ .\n",
"$$\n",
"\n",
"We then start by evaluating\n",
"\n",
"$$\n",
"\\boxed{\n",
"\\begin{align}\n",
"\\bar\\gamma = \\det\\left(\\bar\\gamma_{ij}\\right)\n",
"&= \\bar\\gamma_{xx}\\bar\\gamma_{yy}\\bar\\gamma_{zz} + \\bar\\gamma_{xy}\\bar\\gamma_{yz}\\bar\\gamma_{xz} + \\bar\\gamma_{xz}\\bar\\gamma_{xy}\\bar\\gamma_{yz}\\\\\n",
"&- \\bar\\gamma_{xz}\\bar\\gamma_{yy}\\bar\\gamma_{xz} - \\bar\\gamma_{xy}\\bar\\gamma_{xy}\\bar\\gamma_{zz} - \\bar\\gamma_{xx}\\bar\\gamma_{yz}\\bar\\gamma_{yz}\n",
"\\end{align}\n",
"}\\ .\n",
"$$\n",
"\n",
"In the code below we also define\n",
"\n",
"$$\n",
"\\boxed{{\\rm gtijdet\\_Fm1o3} = \\left(\\frac{1}{\\bar{\\gamma}}\\right)^{1/3}}\\ .\n",
"$$"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Appending to ../src/convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C\n"
]
}
],
"source": [
"%%writefile -a $outfile_path__ADM_to_BSSN__det_gammabar_eq_1__C\n",
"\n",
"\n",
" /*********************************\n",
" * Apply det gtij = 1 constraint *\n",
" *********************************/\n",
" CCTK_REAL gtijdet = gtxxL * gtyyL * gtzzL + gtxyL * gtyzL * gtxzL + gtxzL * gtxyL * gtyzL -\n",
" gtxzL * gtyyL * gtxzL - gtxyL * gtxyL * gtzzL - gtxxL * gtyzL * gtyzL;\n",
"\n",
" CCTK_REAL gtijdet_Fm1o3 = fabs(1.0/cbrt(gtijdet));\n",
"\n",
" gtxxL = gtxxL * gtijdet_Fm1o3;\n",
" gtxyL = gtxyL * gtijdet_Fm1o3;\n",
" gtxzL = gtxzL * gtijdet_Fm1o3;\n",
" gtyyL = gtyyL * gtijdet_Fm1o3;\n",
" gtyzL = gtyzL * gtijdet_Fm1o3;\n",
" gtzzL = gtzzL * gtijdet_Fm1o3;\n",
"\n",
" if(gtijdet<0.0) {\n",
"#ifndef ENABLE_STANDALONE_IGM_C2P_SOLVER\n",
" CCTK_VWarn(CCTK_WARN_ALERT,__LINE__, __FILE__, CCTK_THORNSTRING,\n",
"#else\n",
" fprintf(stderr,\n",
"#endif\n",
"\"WARNING: det[3-metric]<0.0 at point %d %d %d | cctk_lsh: %d %d %d. Hopefully this is occurring in gz's! gtij_phys = %.2e %.2e %.2e %.2e %.2e %.2e gtij_new = %.2e %.2e %.2e %.2e %.2e %.2e | gijdet = %.2e | gtijdet = %.2e\\n\",\n",
"i,j,k,cctk_lsh[0],cctk_lsh[1],cctk_lsh[2],gxx_physL,gxy_physL,gxz_physL,gyy_physL,gyz_physL,gzz_physL,gtxxL,gtxyL,gtxzL,gtyyL,gtyzL,gtzzL,-gijdet,gtijdet);\n",
"}\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"## Step 2.d: Updating the gridfunctions \\[Back to [top](#toc)\\]\n",
"$$\\label{update_gfs}$$\n",
"\n",
"Finally, we update the gridfunctions for $\\phi$, $\\psi$, $\\bar\\gamma_{ij}$, $\\gamma_{ij}$, and $\\bar\\gamma^{ij}$, respectively, in the code below.\n",
"\n",
"If $\\mathcal{M}$ is a $3\\times3$ symmetric matrix and $\\det(\\mathcal{M})\\neq0$, [recall that we have](https://en.wikipedia.org/wiki/Invertible_matrix#Inversion_of_3_%C3%97_3_matrices)\n",
"\n",
"$$\n",
"\\mathcal{M}^{-1} = \n",
"\\begin{pmatrix}\n",
"a & b & c\\\\\n",
"b & d & e\\\\\n",
"c & e & f\n",
"\\end{pmatrix}^{-1}\n",
"=\n",
"\\frac{1}{\\det(\\mathcal{M})}\n",
"\\begin{pmatrix}\n",
"A & B & C\\\\\n",
"B & D & E\\\\\n",
"C & E & F\n",
"\\end{pmatrix}\\ ,\n",
"$$\n",
"\n",
"where\n",
"\n",
"$$\n",
"\\boxed{\n",
"\\begin{align}\n",
"A &= \\ \\ \\left(df - ee\\right)\\\\\n",
"B &= -\\left(bf - ce\\right)\\\\\n",
"C &= \\ \\ \\left(be - cd\\right)\\\\\n",
"D &= \\ \\ \\left(af - cc\\right)\\\\\n",
"E &= -\\left(ae - bc\\right)\\\\\n",
"F &= \\ \\ \\left(ad - bb\\right)\n",
"\\end{align}\n",
"}\\ .\n",
"$$\n",
"\n",
"Notice that if we replace $\\mathcal{M}\\to\\bar\\gamma_{ij}$, then we have also $\\det(\\mathcal{M})\\to1$."
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Appending to ../src/convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C\n"
]
}
],
"source": [
"%%writefile -a $outfile_path__ADM_to_BSSN__det_gammabar_eq_1__C\n",
"\n",
"\n",
" CCTK_REAL Psi4 = psiL*psiL*psiL*psiL;\n",
" /*****************************************\n",
" * Set all the needed BSSN gridfunctions *\n",
" *****************************************/\n",
" phi[index] = phiL;\n",
" psi[index] = psiL;\n",
"\n",
" lapm1[index] = alp[index] - 1.0;\n",
"\n",
" gtxx[index] = gtxxL;\n",
" gtxy[index] = gtxyL;\n",
" gtxz[index] = gtxzL;\n",
" gtyy[index] = gtyyL;\n",
" gtyz[index] = gtyzL;\n",
" gtzz[index] = gtzzL;\n",
"\n",
" gxx[index] = gtxxL*Psi4;\n",
" gxy[index] = gtxyL*Psi4;\n",
" gxz[index] = gtxzL*Psi4;\n",
" gyy[index] = gtyyL*Psi4;\n",
" gyz[index] = gtyzL*Psi4;\n",
" gzz[index] = gtzzL*Psi4;\n",
"\n",
" gtupxx[index] = ( gtyyL * gtzzL - gtyzL * gtyzL );\n",
" gtupxy[index] = - ( gtxyL * gtzzL - gtyzL * gtxzL );\n",
" gtupxz[index] = ( gtxyL * gtyzL - gtyyL * gtxzL );\n",
" gtupyy[index] = ( gtxxL * gtzzL - gtxzL * gtxzL );\n",
" gtupyz[index] = - ( gtxxL * gtyzL - gtxyL * gtxzL );\n",
" gtupzz[index] = ( gtxxL * gtyyL - gtxyL * gtxyL );\n",
" }\n",
"}\n",
"\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 3: Code validation \\[Back to [top](#toc)\\]\n",
"$$\\label{code_validation}$$\n",
"\n",
"First, we download the original `IllinoisGRMHD` source code and then compare it to the source code generated by this tutorial notebook."
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Validation test for convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C: FAILED!\n",
"Diff:\n",
"1d0\n",
"< #include \"cctk.h\"\n",
"4a4,5\n",
"> #ifndef ENABLE_STANDALONE_IGM_C2P_SOLVER\n",
"> #include \"cctk.h\"\n",
"5a7\n",
"> #endif\n",
"13c15,20\n",
"< DECLARE_CCTK_PARAMETERS; \n",
"---\n",
"> \n",
"> #ifndef ENABLE_STANDALONE_IGM_C2P_SOLVER\n",
"> DECLARE_CCTK_PARAMETERS;\n",
"> #endif\n",
"> \n",
"> \n",
"31a39\n",
"> \n",
"42c50,51\n",
"< \n",
"---\n",
"> \n",
"> \n",
"46c55\n",
"< CCTK_REAL gtijdet = gtxxL * gtyyL * gtzzL + gtxyL * gtyzL * gtxzL + gtxzL * gtxyL * gtyzL - \n",
"---\n",
"> CCTK_REAL gtijdet = gtxxL * gtyyL * gtzzL + gtxyL * gtyzL * gtxzL + gtxzL * gtxyL * gtyzL -\n",
"52,55c61,64\n",
"< gtxyL = gtxyL * gtijdet_Fm1o3; \n",
"< gtxzL = gtxzL * gtijdet_Fm1o3; \n",
"< gtyyL = gtyyL * gtijdet_Fm1o3; \n",
"< gtyzL = gtyzL * gtijdet_Fm1o3; \n",
"---\n",
"> gtxyL = gtxyL * gtijdet_Fm1o3;\n",
"> gtxzL = gtxzL * gtijdet_Fm1o3;\n",
"> gtyyL = gtyyL * gtijdet_Fm1o3;\n",
"> gtyzL = gtyzL * gtijdet_Fm1o3;\n",
"58,60c67,76\n",
"< if(gtijdet<0.0) { CCTK_VWarn(CCTK_WARN_ALERT,__LINE__, __FILE__, CCTK_THORNSTRING,\n",
"< \"WARNING: det[3-metric]<0.0 at point %d %d %d | cctk_lsh: %d %d %d. Hopefully this is occurring in gz's! gtij_phys = %.2e %.2e %.2e %.2e %.2e %.2e gtij_new = %.2e %.2e %.2e %.2e %.2e %.2e | gijdet = %.2e | gtijdet = %.2e\",\n",
"< \t\t\t\t i,j,k,cctk_lsh[0],cctk_lsh[1],cctk_lsh[2],gxx_physL,gxy_physL,gxz_physL,gyy_physL,gyz_physL,gzz_physL,gtxxL,gtxyL,gtxzL,gtyyL,gtyzL,gtzzL,-gijdet,gtijdet); }\n",
"---\n",
"> if(gtijdet<0.0) {\n",
"> #ifndef ENABLE_STANDALONE_IGM_C2P_SOLVER\n",
"> CCTK_VWarn(CCTK_WARN_ALERT,__LINE__, __FILE__, CCTK_THORNSTRING,\n",
"> #else\n",
"> fprintf(stderr,\n",
"> #endif\n",
"> \"WARNING: det[3-metric]<0.0 at point %d %d %d | cctk_lsh: %d %d %d. Hopefully this is occurring in gz's! gtij_phys = %.2e %.2e %.2e %.2e %.2e %.2e gtij_new = %.2e %.2e %.2e %.2e %.2e %.2e | gijdet = %.2e | gtijdet = %.2e\\n\",\n",
"> i,j,k,cctk_lsh[0],cctk_lsh[1],cctk_lsh[2],gxx_physL,gxy_physL,gxz_physL,gyy_physL,gyz_physL,gzz_physL,gtxxL,gtxyL,gtxzL,gtyyL,gtyzL,gtzzL,-gijdet,gtijdet);\n",
"> }\n",
"> \n",
"92a109\n",
"> \n"
]
}
],
"source": [
"# Verify if the code generated by this tutorial module\n",
"# matches the original IllinoisGRMHD source code\n",
"\n",
"# First download the original IllinoisGRMHD source code\n",
"import urllib\n",
"from os import path\n",
"\n",
"original_IGM_file_url = \"https://bitbucket.org/zach_etienne/wvuthorns/raw/5611b2f0b17135538c9d9d17c7da062abe0401b6/IllinoisGRMHD/src/convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C\"\n",
"original_IGM_file_name = \"convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij-original.C\"\n",
"original_IGM_file_path = os.path.join(IGM_src_dir_path,original_IGM_file_name)\n",
"\n",
"# Then download the original IllinoisGRMHD source code\n",
"# We try it here in a couple of ways in an attempt to keep\n",
"# the code more portable\n",
"try:\n",
" original_IGM_file_code = urllib.request.urlopen(original_IGM_file_url).read().decode(\"utf-8\")\n",
" # Write down the file the original IllinoisGRMHD source code\n",
" with open(original_IGM_file_path,\"w\") as file:\n",
" file.write(original_IGM_file_code)\n",
"except:\n",
" try:\n",
" original_IGM_file_code = urllib.urlopen(original_IGM_file_url).read().decode(\"utf-8\")\n",
" # Write down the file the original IllinoisGRMHD source code\n",
" with open(original_IGM_file_path,\"w\") as file:\n",
" file.write(original_IGM_file_code)\n",
" except:\n",
" # If all else fails, hope wget does the job\n",
" !wget -O $original_IGM_file_path $original_IGM_file_url\n",
"\n",
"# Perform validation\n",
"Validation__ADM_to_BSSN__det_gammabar_eq_1__C = !diff $original_IGM_file_path $outfile_path__ADM_to_BSSN__det_gammabar_eq_1__C\n",
"\n",
"if Validation__ADM_to_BSSN__det_gammabar_eq_1__C == []:\n",
" # If the validation passes, we do not need to store the original IGM source code file\n",
" !rm $original_IGM_file_path\n",
" print(\"Validation test for convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C: PASSED!\")\n",
"else:\n",
" # If the validation fails, we keep the original IGM source code file\n",
" print(\"Validation test for convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.C: FAILED!\")\n",
" # We also print out the difference between the code generated\n",
" # in this tutorial module and the original IGM source code\n",
" print(\"Diff:\")\n",
" for diff_line in Validation__ADM_to_BSSN__det_gammabar_eq_1__C:\n",
" print(diff_line)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 4: 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-IllinoisGRMHD__convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.pdf](Tutorial-IllinoisGRMHD__convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.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": 7,
"metadata": {},
"outputs": [],
"source": [
"latex_nrpy_style_path = os.path.join(nrpy_dir_path,\"latex_nrpy_style.tplx\")\n",
"#!jupyter nbconvert --to latex --template $latex_nrpy_style_path --log-level='WARN' Tutorial-IllinoisGRMHD__convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.ipynb\n",
"#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.tex\n",
"#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.tex\n",
"#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij.tex\n",
"!rm -f Tut*.out Tut*.aux Tut*.log"
]
}
],
"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
}