<jupyter_start><jupyter_text>window.dataLayer = window.dataLayer || [];  function gtag(){dataLayer.push(arguments);}  gtag('js', new Date());  gtag('config', 'UA-59152712-8'); Tutorial-IllinoisGRMHD: compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.C Authors: <NAME> & <NAME>**This module is currently under development** In this tutorial module we explain how we evaluate $\partial_{t}\tilde{\tau}$, terms that depend on the extrinsic curvature, $K_{ij}$, and the energy-momentum tensor $T^{\mu\nu}$ Required and recommended citations:* **(Required)** <NAME>., <NAME>., <NAME>., <NAME>., and <NAME>. 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)).* **(Required)** <NAME>., <NAME>., <NAME>., <NAME>. Primitive Variable Solvers for Conservative General Relativistic Magnetohydrodynamics. Astrophysical Journal, 641, 626 (2006) ([astro-ph/0512420](https://arxiv.org/abs/astro-ph/0512420)).* **(Recommended)** <NAME>., <NAME>., <NAME>. 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)).If using the version of `IllinoisGRMHD` with piecewise polytropic *or* tabulated (coming soon!) EOS support, then the following citation is also required:* **(Required)** <NAME>., <NAME>., <NAME>., <NAME>., <NAME>., and <NAME>., *IllinoisGRMHD github repository* (2019). Source Code URL: https://github.com/zachetienne/nrpytutorial/tree/master/IllinoisGRMHD/.  DependenciesThe files generated in this tutorial notebook depend on the following files:*  [GRHD/equations NRPy+ module](/edit/NRPyIGM/GRHD/equations.py) \[[**tutorial**](/notebooks/NRPyIGM/Tutorial-GRHD_Equations-Cartesian.ipynb)\]*  [GRFFE/equations NRPy+ module](/edit/NRPyIGM/GRFFE/equations.py) \[[**tutorial**](/notebooks/NRPyIGM/Tutorial-GRFFE_Equations-Cartesian.ipynb)\]*  [GRMHD/equations NRPy+ module](/edit/NRPyIGM/GRMHD/equations.py) \[[**tutorial**](/notebooks/NRPyIGM/Tutorial-GRMHD_Equations-Cartesian.ipynb)\]* `IllinoisGRMHD_EoS_lowlevel_functs.C` \[[**tutorial**](Tutorial-IllinoisGRMHD__EoS_lowlevel_functs.ipynb)\]  Table of Contents$$\label{toc}$$This module is organized as follows0. [Step 0](src_dir): **Source directory creation**1. [Step 1](compute_extrinsic_curvature_terms_using_nrpy): **Computing the extrinsic curvature terms in $\partial_{t}\tilde{\tau}$ RHS using NRPy+**    1. [Step 1.a](load_python_nrpy_modules): *Load necessary Python/NRPy+ modules*    1. [Step 1.b](adm_3metric_from_conf_metric): *Compute ADM 3-metric quantities from conformal metric*    1. [Step 1.c](alpha_beta_gammadn_from_conf_metric_header): *The `ADM_3METRIC_in_terms_of_CONF_METRIC.h` file*    1. [Step 1.d](setting_primitives): *Setting primitive variables*    1. [Step 1.e](setting_eos_quantities): *Setting EOS quantities*    1. [Step 1.f](compute_smallb4_and_smallbsquared): *The `compute_smallb4U_and_smallbsquared_from_ADM_3METRIC_and_PRIMS.h` file*    1. [Step 1.g](grmhd_tupmunu): *Computing $T^{\mu\nu}_{\rm GRMHD}$*    1. [Step 1.h](compute_tau_rhs_extrinsic_curvature_header): *The `compute_tau_rhs_extrinsic_curvature_terms_and_TUPmunu.h` file*1. [Step 2](compute_tau_rhs_extrinsic_curvature_terms_and_tupmunu__c): **`compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.C`**1. [Step 3](code_validation): **Code validation**1. [Step 4](latex_pdf_output): **Output this notebook to $\LaTeX$-formatted PDF file**  Step 0: Source directory creation \[Back to [top](toc)\]$$\label{src_dir}$$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.<jupyter_code># Step 0: Creation of the IllinoisGRMHD source directory
# Step 0a: 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 0b: Load up cmdline_helper and create the directory
import cmdline_helper as cmd
outdir = os.path.join("..","src")
cmd.mkdir(outdir)
NRPy_headers_dir_path = os.path.join(outdir,"NRPy_generated_headers")
cmd.mkdir(NRPy_headers_dir_path)<jupyter_output><empty_output><jupyter_text>Step 1: Computing the extrinsic curvature terms in $\partial_{t}\tilde{\tau}$ RHS using NRPy+ \[Back to [top](toc)\]$$\label{compute_extrinsic_curvature_terms_using_nrpy}$$$$\partial_{t}\tilde\tau = \left[\rm Flux\ term\right] + s\ ,$$where$$s = \underbrace{\alpha\sqrt{\gamma}\left(T^{00}\beta^{i}\beta^{j} + 2T^{0i}\beta^{j} + T^{ij}\right)K_{ij}}_{\text{Extrinsic curvature terms}} - \underbrace{\left(T^{00}\beta^{i} + T^{0i}\right)\partial_{i}\alpha}_{\text{Will be implemented later}}\ .$$To compute the extrinsic curvature terms, we will be using the following NRPy+ modules:* [GRHD.equations NRPy+ module](/edit/NRPyIGM/GRHD/equations.py): used to compute $u^{\mu}$ and the extrinsic curvature terms in $s$.* [GRFFE.equations NRPy+ module](/edit/NRPyIGM/GRFFE/equations.py): used to compute $b^{\mu}$ and $b^{2}\equiv b^{\mu}b_{\mu}$.* [GRMHD.equations NRPy+ module](/edit/NRPyIGM/GRMHD/equations.py): used to compute $T^{\mu\nu}_{\rm GRMHD} = T^{\mu\nu}_{\rm GRHD} + T^{\mu\nu}_{\rm GRFFE}$.  Step 1.a: Load necessary Python/NRPy+ modules \[Back to [top](toc)\]$$\label{load_python_nrpy_modules}$$We now load all necessary Python/NRPy+ modules needed to perform the calculation of $s$.<jupyter_code># Imported needed Python modules
import time        # Python module: we will use the time() function to estimate code execution time
import sympy as sp # Python module: symbolic expressions capabilities

# Imported needed NRPy+ modules
from outputC import *                           # NRPy+ module: C code output routines
import indexedexp as ixp                        # NRPy+ module: funcionality to work with vectors/tensors
import IllinoisGRMHD_output_functions as IGMout # NRPy+ module: IllinoisGRMHD output file functions
import GRHD.equations  as GRHD                  # NRPy+ module: Collection GRHD related functions
import GRFFE.equations as GRFFE                 # NRPy+ module: Collection GRFFE related functions
import GRMHD.equations as GRMHD                 # NRPy+ module: Collection GRMHD related functions<jupyter_output><empty_output><jupyter_text>Step 1.b: Compute ADM 3-metric quantities from conformal metric \[Back to [top](toc)\]$$\label{adm_3metric_from_conf_metric}$$We now compute the ADM 3 metric quantities $\left(\alpha,\beta^{i},\gamma_{ij},\gamma^{ij},\sqrt{\gamma}\right)$ from the conformal metric quantities $\left(\alpha-1,\beta^{i},\bar\gamma_{ij},\bar\gamma^{ij},\phi,\psi\right)$.Setting $\alpha$ and $\beta^{i}$ is straightfoward. For $\gamma_{ij}$, $\gamma^{ij}$, and $\sqrt{\gamma}$ we have\begin{align}\gamma_{ij}   &= \psi^{4}\bar\gamma_{ij} \ ,\\\gamma^{ij}   &= \psi^{-4}\bar\gamma^{ij}\ ,\\\sqrt{\gamma} &= \psi^{6}\ ,\end{align}where the last equality holds true because $\sqrt{\bar\gamma}=1$.<jupyter_code># Set DIM = 3
DIM = 3

# Set up alpha in terms of conformal metric
alpha_ito_CONF_METRIC = sp.Symbol("CONF_METRIC[CM_LAPM1]",real=True)+1.0

# Set up beta^{i} in terms of conformal metric
beta_ito_CONF_METRIC_U = ixp.zerorank1()
for i in range(DIM):
    beta_ito_CONF_METRIC_U[i] = sp.Symbol("CONF_METRIC[CM_SHIFT"+chr(ord('X')+i)+"]",real=True)

# Set up psi^{4}, psi^{-4}, and psi^{6}
psi4  = sp.Symbol("CONF_METRIC[CM_PSI]",real=True)**4
psim4 = 1/psi4
psi6  = sp.Symbol("CONF_METRIC[CM_PSI]",real=True)**6

# Set up gamma_{ij} and gamma^{ij} in terms of conformal metric
gamma_ito_CONF_METRIC_DD = ixp.zerorank2()
gamma_ito_CONF_METRIC_UU = ixp.zerorank2()
for i in range(DIM):
    for j in range(i,DIM):
        gamma_ito_CONF_METRIC_DD[i][j] = gamma_ito_CONF_METRIC_DD[j][i] = psi4  * sp.Symbol("CONF_METRIC[CM_GAMMATILDE"+chr(ord('X')+i)+chr(ord('X')+j)+"]",real=True)
        gamma_ito_CONF_METRIC_UU[i][j] = gamma_ito_CONF_METRIC_UU[j][i] = psim4 * sp.Symbol("CONF_METRIC[CM_GAMMATILDEUP"+chr(ord('X')+i)+chr(ord('X')+j)+"]",real=True)
        
# Set up sqrtgamma
sqrtgamma_ito_CONF_METRIC = psi6<jupyter_output><empty_output><jupyter_text>Step 1.c: The `ADM_3METRIC_in_terms_of_CONF_METRIC.h` file \[Back to [top](toc)\]$$\label{alpha_beta_gammadn_from_conf_metric_header}$$We now generate the `ADM_3METRIC_in_terms_of_CONF_METRIC.h` file, which allows `IllinoisGRMHD` to populate the `ADM_3METRIC` array in terms of the `CONF_METRIC` array. The structure of the following code cell is the following:1. We start by setting up variables that point to the `ADM_3METRIC` array.1. We then set up a list of expressions and variables, to be converted into C code.1. We then convert these lists to C code using the `outputC()` function from the [outputC NRPy+ module](/edit/NRPyIGM/outputC.py).1. Finally, we output the resulting string to a `IllinoisGRMHD` header file using the `NRPy_IGM_write_to_file()` function from the [IllinoisGRMHD_output_functions NRPy+ module](/edit/NRPyIGM/IllinoisGRMHD/doc/IllinoisGRMHD_output_functions.py).<jupyter_code># Declare basic ADM variables to be used by IllinoisGRMHD
# Set spatial dimension to 3

# Set up alpha
ADMalpha  = sp.Symbol("ADM_3METRIC[ALPHA]",real=True)
# Set up beta^{i}
ADMbetaU   = ixp.zerorank1()
for i in range(DIM):
    ADMbetaU[i] = sp.Symbol('ADM_3METRIC[BETA'+chr(ord('X')+i)+"]",real=True)

# Set up gamma_{ij}
ADMgammaDD = ixp.zerorank2()
for i in range(DIM):
    for j in range(i,DIM):
        ADMgammaDD[i][j] = ADMgammaDD[j][i] = sp.Symbol('ADM_3METRIC[GAMMA'+chr(ord('X')+i)+chr(ord('X')+j)+"]",real=True)

# Set up gamma^{ij}
ADMgammaUU = ixp.zerorank2()
for i in range(3):
    for j in range(i,3):
        ADMgammaUU[i][j] = ADMgammaUU[j][i] = sp.Symbol('ADM_3METRIC[GAMMAUP'+chr(ord('X')+i)+chr(ord('X')+j)+"]",real=True)
        
# Set up sqrt(gamma)
ADMsqrtgamma = sp.Symbol("ADM_3METRIC[SQRTGAMMA]",real=True)

# Set up expressions and variables lists
exprlist = [alpha_ito_CONF_METRIC]
varslist = [ADMalpha]
for i in range(DIM):
    exprlist.append(beta_ito_CONF_METRIC_U[i])
    varslist.append(ADMbetaU[i])
for i in range(DIM):
    for j in range(i,DIM):
        exprlist.append(gamma_ito_CONF_METRIC_DD[i][j])
        varslist.append(ADMgammaDD[i][j])
for i in range(DIM):
    for j in range(i,DIM):
        exprlist.append(gamma_ito_CONF_METRIC_UU[i][j])
        varslist.append(ADMgammaUU[i][j])
exprlist.append(sqrtgamma_ito_CONF_METRIC)
varslist.append(ADMsqrtgamma)
        
string = outputC(exprlist,varslist,"returnstring", params="outCverbose=False")

# Output to file
filename = "ADM_3METRIC_in_terms_of_CONF_METRIC.h"
filepath = os.path.join(NRPy_headers_dir_path,filename)
IGMout.NRPy_IGM_write_to_file(filepath,filename,string)<jupyter_output>Just generated the file: ../src/NRPy_generated_headers/ADM_3METRIC_in_terms_of_CONF_METRIC.h<jupyter_text>Step 1.d: Setting primitive variables \[Back to [top](toc)\]$$\label{setting_primitives}$$At this stage in the calculation, we will have access to the `PRIMS` array, which is an array of pointers to the gridfunctions `prims`. Thus, we will set up our primitive variables $\left(\rho_{b},P,v^{i}B^{i}\right)$ using `PRIMS`. We will also set up the constant `sqrt4pi` which, as the name suggests, will be set to $\sqrt{4\pi}$ in the C code.<jupyter_code># Set up primitives
# Set rho and P
rhob, P = sp.symbols("PRIMS[RHOB] PRIMS[PRESSURE]",real=True)

# Set v^{i} and B^{i}
vU         = ixp.zerorank1()
B_notildeU = ixp.zerorank1()
for i in range(DIM):
    vU[i]         = sp.Symbol("PRIMS[V"+chr(ord('X')+i)+"]",real=True)
    B_notildeU[i] = sp.Symbol("PRIMS[B"+chr(ord('X')+i)+"_CENTER]",real=True)

# Set sqrt(4*pi)
sqrt4pi = sp.symbols('sqrt4pi', real=True)<jupyter_output><empty_output><jupyter_text>Step 1.e: Setting EOS quantities \[Back to [top](toc)\]$$\label{setting_eos_quantities}$$We now set up needed equation of state (EOS) quantities, namely $\left(\epsilon,P_{\rm cold},\epsilon_{\rm cold},\Gamma_{\rm th}\right)$. $\epsilon$ is not defined as a symbol, but instead it is computed using the hybrid EOS$$P = P_{\rm cold} + \left(\Gamma_{\rm th}-1\right)\rho_{b}\epsilon_{\rm th} = P_{\rm cold} + \left(\Gamma_{\rm th}-1\right)\rho_{b}\left(\epsilon - \epsilon_{\rm cold}\right)\implies\boxed{\epsilon = \epsilon_{\rm cold} + \frac{P-P_{\rm cold}}{\rho_{b}\left(\Gamma_{\rm th}-1\right)}}\ .$$<jupyter_code># Set up EOS quantities
# Set P_cold, eps_cold, Gamma_th
P_cold, eps_cold, Gamma_th = sp.symbols("P_cold eps_cold Gamma_th",real=True)

# Compute epsilon
epsilon = eps_cold + (P - P_cold)/(rhob * (Gamma_th - 1))<jupyter_output><empty_output><jupyter_text>Step 1.f: The `compute_smallb4U_and_smallbsquared_from_ADM_3METRIC_and_PRIMS.h` file \[Back to [top](toc)\]$$\label{compute_smallb4_and_smallbsquared}$$We now compute $b^{\mu}$ and $b^{2}\equiv b^{\mu}b_{\mu}$ from the ADM quantities and primitives. This is done below in three separate steps:1. Compute $u^{\mu}$ using the `u4U_in_terms_of_vU__rescale_vU_by_applying_speed_limit()` function from the [GRHD.equations NRPy+ module](/edit/NRPyIGM/GRHD/equations.py).1. Compute $b^{\mu}$ using the `compute_smallb4U()` function from the [GRFFE.equations NRPy+ module](/edit/NRPyIGM/GRFFE/equations.py).1. Compute $b^{2}$ using the `compute_smallbsquared()` function from the [GRFFE.equations NRPy+ module](/edit/NRPyIGM/GRFFE/equations.py).Then, we generate the `compute_smallb4U_and_smallbsquared_from_ADM_3METRIC_and_PRIMS.h` file, which allow the quantities above to be used by `IllinoisGRMHD`.<jupyter_code># Define smallb4U and smallbsquared to be used later when
# computing T^{\mu\nu}, thus simplifying the expression
# generation
smallbsquared = sp.Symbol("smallbsquared",real=True)
smallb4U = ixp.zerorank1(DIM=4)
for mu in range(4):
    smallb4U[mu] = sp.Symbol("smallb4U"+str(mu),real=True)

# Then compute u^{\mu}
GRHD.u4U_in_terms_of_vU__rescale_vU_by_applying_speed_limit(ADMalpha,ADMbetaU,ADMgammaDD, vU)
u4U = GRHD.u4U_ito_vU

# Compute b^{\mu}
GRFFE.compute_smallb4U(ADMgammaDD,ADMbetaU,ADMalpha, u4U,B_notildeU, sqrt4pi)

# Compute b^{\mu}b_{\mu}
GRFFE.compute_smallbsquared(ADMgammaDD,ADMbetaU,ADMalpha, GRFFE.smallb4U)

# Output b^{\mu} and b^{2} to string
string = outputC([GRFFE.smallb4U[0],GRFFE.smallb4U[1],GRFFE.smallb4U[2],GRFFE.smallb4U[3],GRFFE.smallbsquared],
                 [smallb4U[0],      smallb4U[1],      smallb4U[2],      smallb4U[3],      smallbsquared],
                 filename="returnstring", params="outCverbose=False")

# Output string to file
filename="compute_smallb4U_and_smallbsquared_from_ADM_3METRIC_and_PRIMS.h"
filepath = os.path.join(NRPy_headers_dir_path,filename)
IGMout.NRPy_IGM_write_to_file(filepath,filename,string)<jupyter_output>Just generated the file: ../src/NRPy_generated_headers/compute_smallb4U_and_smallbsquared_from_ADM_3METRIC_and_PRIMS.h<jupyter_text>Step 1.g: Computing $T^{\mu\nu}_{\rm GRMHD}$ \[Back to [top](toc)\]$$\label{grmhd_tupmunu}$$Now that we have already defined and have access to the ADM quantities $\left(\alpha,\beta^{i},\gamma_{ij},\gamma^{ij},\sqrt{\gamma}\right)$, primitive variables $\left(P,\rho_{b},v^{i},B^{i}\right)$, EOS quantities $\left(P_{\rm cold},\epsilon_{\rm cold},\Gamma_{\rm th},\epsilon\right)$, fluid 4-velocity $u^{\mu}$, and the rescaled magnetic field measured by an observer comoving with the fluid, $\left(b^{\mu},b^{2}\right)$, we are ready to compute the GRMHD energy-momentum tensor:$$T^{\mu\nu}_{\rm GRMHD}=T^{\mu\nu}_{\rm GRHD} + T^{\mu\nu}_{\rm GRFFE}=\left(\rho_{b}h + b^{2}\right)u^{\mu}u^{\nu} + \left(P + \frac{b^{2}}{2}\right)g^{\mu\nu} - b^{\mu}b^{\nu}\ .$$This is done in a single step below, by calling the `compute_GRMHD_T4UU()` function from the [GRMHD.equations NRPy+ module](/edit/NRPyIGM/GRMHD/equations.py).<jupyter_code># Compute the GRMHD energy-momentum tensor
GRMHD.compute_GRMHD_T4UU(ADMgammaDD, ADMbetaU, ADMalpha, rhob, P, epsilon, u4U, smallb4U, smallbsquared)<jupyter_output><empty_output><jupyter_text>Step 1.h: The `compute_tau_rhs_extrinsic_curvature_terms_and_TUPmunu.h` file \[Back to [top](toc)\]$$\label{compute_tau_rhs_extrinsic_curvature_header}$$We now generate the `compute_tau_rhs_extrinsic_curvature_terms_and_TUPmunu.h` file, to be used by `IllinoisGRMHD` to compute the extrinsic curvature terms on the RHSs of the evolution equation of $\tilde{\tau}$. This is done in the following steps:1. Define $K_{ij}$ according to `IllinoisGRMHD`'s notation.1. Set $\partial_{i}\alpha=0$, since we are interested only in the extrinsic curvature terms in $s$.1. Compute the curvature terms in $s$ using the `compute_s_source_term()` function from the [GRHD.equations NRPy+ module](/edit/NRPyIGM/GRHD/equations.py).1. Add $T^{\mu\nu}$ and $s$ to expressions and variables lists.1. We then convert these lists to C code using the `outputC()` function from the [outputC NRPy+ module](/edit/NRPyIGM/outputC.py).1. Finally, we output the resulting string to a `IllinoisGRMHD` header file using the `NRPy_IGM_write_to_file()` function from the [IllinoisGRMHD_output_functions NRPy+ module](/edit/NRPyIGM/IllinoisGRMHD/doc/IllinoisGRMHD_output_functions.py).<jupyter_code># Set K_{ij}
KLDD = ixp.zerorank2()
for i in range(DIM):
    for j in range(i,DIM):
        KLDD[i][j] = KLDD[j][i] = sp.Symbol('K'+chr(ord('x')+i)+chr(ord('x')+j)+"L",real=True)

# Set up the derivative of alpha to zero, as we wish
# only to compute the extrinsic curvature terms of
# \partial_{t}\tilde{\tau}
alpha_zero_derivD = ixp.zerorank1()
GRHD.compute_s_source_term(KLDD,ADMbetaU,ADMalpha, ADMsqrtgamma,alpha_zero_derivD, GRMHD.T4UU)

# Add T^{\mu\nu} to the variables/expressions lists
exprlist = []
varslist = []
counter = 0
for mu in range(4):
    for nu in range(mu,4):
        exprlist.append(GRMHD.T4UU[mu][nu])
        varslist.append("TUPMUNU["+str(counter)+"][index]")
        counter += 1

# Add tau source term to variables/expressions lists
exprlist.append(GRHD.s_source_term)
varslist.append("""
if(k<cctk_lsh[2]-cctk_nghostzones[2] && j<cctk_lsh[1]-cctk_nghostzones[1] && i<cctk_lsh[0]-cctk_nghostzones[0]) {

tau_rhs[index]""")

# Output expressions to string
string = outputC(exprlist,varslist,filename="returnstring", params="outCverbose=False", poststring="\n   }\n\n")

# Output string to file
filename = "compute_tau_rhs_extrinsic_curvature_terms_and_TUPmunu.h"
filepath = os.path.join(NRPy_headers_dir_path,filename)
IGMout.NRPy_IGM_write_to_file(filepath,filename,string)<jupyter_output>Just generated the file: ../src/NRPy_generated_headers/compute_tau_rhs_extrinsic_curvature_terms_and_TUPmunu.h<jupyter_text>Step 2: The `compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.C` file \[Back to [top](toc)\]$$\label{compute_tau_rhs_extrinsic_curvature_terms_and_tupmunu__c}$$We now write the `compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.C` file, which puts together all the header files we set up this far. The overall structure of the file is the following:* Looping over all gridpoints:    1. Read in $K_{ij}$ from gridfunctions    1. Read in `CONF_METRIC` quantities from input    1. Compute `ADM_3METRIC` using `ADM_3METRIC_in_terms_of_CONF_METRIC.h`, generated in [Step 1.c](alpha_beta_gammadn_from_conf_metric_header)    1. Read in `PRIMS` from input    1. Compute $b^{\mu}$ and $b^{2}$ using `compute_smallb4U_and_smallbsquared_from_ADM_3METRIC_and_PRIMS.h`, generated in [Step 1.f](compute_smallb4_and_smallbsquared)    1. Compute EOS quantities    1. Compute $T^{\mu\nu}$ and $s$, using `compute_tau_rhs_extrinsic_curvature_terms_and_TUPmunu.h`, generated in [Step 1.h](compute_tau_rhs_extrinsic_curvature_header)<jupyter_code>%%writefile $outdir/compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.C
static void compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU
(const cGH *cctkGH,const int *cctk_lsh,const int *cctk_nghostzones,CCTK_REAL *dX, CCTK_REAL **CONF_METRIC_IN,gf_and_gz_struct *prims,
 CCTK_REAL **TUPMUNU,eos_struct &eos,
 CCTK_REAL *gupxy,CCTK_REAL *gupxz,CCTK_REAL *gupyz,
 CCTK_REAL *kxx,CCTK_REAL *kxy,CCTK_REAL *kxz,CCTK_REAL *kyy,CCTK_REAL *kyz,CCTK_REAL *kzz,
 CCTK_REAL *tau_rhs) {

  // Set sqrt4pi
  const CCTK_REAL sqrt4pi = sqrt(4.0*M_PI);

  // These loop extents must be consistent with add_fluxes_and_source_terms_to_hydro_rhss(), since we use TUPMUNU there as well.
#pragma omp parallel for
  for(int k=cctk_nghostzones[2];k<cctk_lsh[2]-(cctk_nghostzones[2]-1);k++) for(int j=cctk_nghostzones[1];j<cctk_lsh[1]-(cctk_nghostzones[1]-1);j++) for(int i=cctk_nghostzones[0];i<cctk_lsh[0]-(cctk_nghostzones[0]-1);i++) {
        int index = CCTK_GFINDEX3D(cctkGH,i,j,k);

        DECLARE_CCTK_PARAMETERS;

        // First we pull in needed hydrodynamic and metric variables from memory: PART 1.
        // Reading from main memory is a SLOW operation, usually resulting in
        //   cache misses, which
        //   will waste precious runtime. Note that cache misses will often not
        //   show up when using, e.g., gprof. The slowdown due to cache misses
        //   can more than double the amount of time in this routine, so instead
        //   of reading in variables from main memory multiple times, the below
        //   forces us to read in only once, storing in local variables 
        //   U{p,m}{1,2,3}, METRIC{p,m}{1,2}, U, and METRIC.

        // Read in K_{ij} from gridfunctions
        CCTK_REAL KxxL = kxx[index];
        CCTK_REAL KxyL = kxy[index];
        CCTK_REAL KxzL = kxz[index];
        CCTK_REAL KyyL = kyy[index];
        CCTK_REAL KyzL = kyz[index];
        CCTK_REAL KzzL = kzz[index];

        // Read in conformal metric quantitites from gridfunctions
        CCTK_REAL CONF_METRIC[NUMVARS_FOR_CONF_METRIC_FACEVALS];
        for(int ii=0;ii<NUMVARS_FOR_CONF_METRIC_FACEVALS;ii++) CONF_METRIC[ii] = CONF_METRIC_IN[ii][index];

        // Compute ADM_3METRIC quantities from CONF_METRIC
        CCTK_REAL ADM_3METRIC[NUMVARS_FOR_ADM_3METRIC];
#include "NRPy_generated_headers/ADM_3METRIC_in_terms_of_CONF_METRIC.h"

        // The "vector" PRIMS represents the primitive variables: rho, P, vx, vy, vz, Bx, By, and Bz.
        CCTK_REAL PRIMS[8]; // 8 primitives in the set: {rho_b,P,vx,vy,vz,Bx,By,Bz}
        for(int ii=0;ii<8;ii++) PRIMS[ii] = prims[ii].gf[index];

        // Define b^{\mu} variables and compute b^{\mu} and b^{2} and compute them
        CCTK_REAL smallb4U0, smallb4U1, smallb4U2, smallb4U3, smallbsquared;
#include "NRPy_generated_headers/compute_smallb4U_and_smallbsquared_from_ADM_3METRIC_and_PRIMS.h"

        struct output_stats stats; stats.failure_checker=0;

        // Set EOS quantities for the computation of T^{\mu\nu}
        CCTK_REAL h=0,  P_cold,eps_cold,dPcold_drho,eps_th,Gamma_cold; /* <- Note that in setting h, we need to define several 
                                                                        *    other variables, even though they will be unused later
                                                                        *    in this function. */
        compute_P_cold__eps_cold__dPcold_drho__eps_th__h__Gamma_cold(PRIMS,eos,Gamma_th,P_cold,eps_cold,dPcold_drho,eps_th,h,Gamma_cold);
      
        CCTK_REAL eps = eps_cold + eps_th;

        if((i==14) && (j==14) && (k==14)) {
          printf("eeee %e %e %e %e %e %e %e %e %e %e\n",PRIMS[PRESSURE],PRIMS[RHOB],P_cold,eps_cold,dPcold_drho,eps_th,Gamma_cold,h, eps, eps_cold + eps_th);
          printf("ffff ");
          for(int ll=0; ll<8; ll++) printf("%e ",PRIMS[ll]);
          printf("\n");
          printf("gggg ");
          for(int ll=0; ll<NUMVARS_FOR_CONF_METRIC_FACEVALS; ll++) printf("%e ",CONF_METRIC[ll]);
          printf("\n");
        }

        // Compute T^{\mu\nu} and the extrinsic curvatures in tau_rhs
#include "NRPy_generated_headers/compute_tau_rhs_extrinsic_curvature_terms_and_TUPmunu.h"

      }
}<jupyter_output>Overwriting ../src/compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.C<jupyter_text>Step 3: Code validation \[Back to [top](toc)\]$$\label{code_validation}$$First we download the original `IllinoisGRMHD` source code and then compare it to the source code generated by this tutorial notebook.<jupyter_code># # Verify if the code generated by this tutorial module
# # matches the original IllinoisGRMHD source code

# # First download the original IllinoisGRMHD source code
# import urllib
# from os import path

# original_IGM_file_url  = "https://bitbucket.org/zach_etienne/wvuthorns/raw/5611b2f0b17135538c9d9d17c7da062abe0401b6/IllinoisGRMHD/src/compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.C"
# original_IGM_file_name = "compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU-original.C"
# original_IGM_file_path = os.path.join(IGM_src_dir_path,original_IGM_file_name)

# # Then download the original IllinoisGRMHD source code
# # We try it here in a couple of ways in an attempt to keep
# # the code more portable
# try:
#     original_IGM_file_code = urllib.request.urlopen(original_IGM_file_url).read().decode("utf-8")
#     # Write down the file the original IllinoisGRMHD source code
#     with open(original_IGM_file_path,"w") as file:
#         file.write(original_IGM_file_code)
# except:
#     try:
#         original_IGM_file_code = urllib.urlopen(original_IGM_file_url).read().decode("utf-8")
#         # Write down the file the original IllinoisGRMHD source code
#         with open(original_IGM_file_path,"w") as file:
#             file.write(original_IGM_file_code)
#     except:
#         # If all else fails, hope wget does the job
#         !wget -O $original_IGM_file_path $original_IGM_file_url

# # Perform validation
# Validation__compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU__C  = !diff $original_IGM_file_path $outfile_path__compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU__C

# if Validation__compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU__C == []:
#     # If the validation passes, we do not need to store the original IGM source code file
#     !rm $original_IGM_file_path
#     print("Validation test for compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.C: PASSED!")
# else:
#     # If the validation fails, we keep the original IGM source code file
#     print("Validation test for compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.C: FAILED!")
#     # We also print out the difference between the code generated
#     # in this tutorial module and the original IGM source code
#     print("Diff:")
#     for diff_line in Validation__compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU__C:
#         print(diff_line)<jupyter_output><empty_output><jupyter_text>Step 4: Output this notebook to $\LaTeX$-formatted PDF file \[Back to [top](toc)\]$$\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-IllinoisGRMHD__compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.pdf](Tutorial-IllinoisGRMHD__compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.pdf) (Note that clicking on this link may not work; you may need to open the PDF file through another means).<jupyter_code>latex_nrpy_style_path = os.path.join(nrpy_dir_path,"latex_nrpy_style.tplx")
#!jupyter nbconvert --to latex --template $latex_nrpy_style_path --log-level='WARN' Tutorial-IllinoisGRMHD__compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.ipynb
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.tex
!rm -f Tut*.out Tut*.aux Tut*.log<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1315114373263023, 'fraction_numerical': 0.014893465402978694, 'mean_word_length': 5.2440489432703, 'pattern_counts': {'":': 0, '<': 76, '<?xml version=': 0, '>': 65, 'https://': 5, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 10, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '12', 'dir': 'jupyter-structured-clean-dedup', 'id': '279797', 'license': 'BSD-2-Clause', 'n_tokens_mistral': 10758, 'n_tokens_neox': 10250, 'n_words': 2237, 'path': 'IllinoisGRMHD/doc/Tutorial-IllinoisGRMHD__compute_tau_rhs_extrinsic_curvature_terms_and_TUPMUNU.ipynb', 'repo_name': 'leowerneck/NRPyIGM', 'size': '28171'}