#!/usr/bin/env python
# coding: utf-8

# <script async src="https://www.googletagmanager.com/gtag/js?id=UA-59152712-8"></script>
# <script>
#   window.dataLayer = window.dataLayer || [];
#   function gtag(){dataLayer.push(arguments);}
#   gtag('js', new Date());
# 
#   gtag('config', 'UA-59152712-8');
# </script>
# 
# # Tutorial-IllinoisGRMHD: ID_converter_ILGRMHD ETKThorn
# 
# ## Authors: Leo Werneck & Zach Etienne
# 
# <font color='red'>**This module is currently under development**</font>
# 
# ## In this tutorial module we generate the ID_converter_ILGRMHD ETK thorn files, compatible with our latest implementation of IllinoisGRMHD
# 
# ### Required and recommended citations:
# 
# * **(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)).
# * **(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)).
# * **(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)).

# <a id='toc'></a>
# 
# # Table of Contents
# $$\label{toc}$$
# 
# This module is organized as follows
# 
# 0. [Step 0](#src_dir): **Source directory creation**
# 1. [Step 1](#introduction): **Introduction**
# 1. [Step 2](#convert_to_hydrobase__src): **`set_IllinoisGRMHD_metric_GRMHD_variables_based_on_HydroBase_and_ADMBase_variables.C`**
# 1. [Step 3](#convert_to_hydrobase__param): **`param.ccl`**
# 1. [Step 4](#convert_to_hydrobase__interface): **`interface.ccl`**
# 1. [Step 5](#convert_to_hydrobase__schedule): **`schedule.ccl`**
# 1. [Step 6](#convert_to_hydrobase__make): **`make.code.defn`**
# 1. [Step n-1](#code_validation): **Code validation**
# 1. [Step n](#latex_pdf_output): **Output this notebook to $\LaTeX$-formatted PDF file**

# <a id='src_dir'></a>
# 
# # 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.

# In[1]:


# Step 0: Creation of the IllinoisGRMHD source directory
# Step 0a: Load up cmdline_helper and create the directory
import os,sys
nrpy_dir_path = os.path.join("..","..")
if nrpy_dir_path not in sys.path:
    sys.path.append(nrpy_dir_path)

import cmdline_helper as cmd
IDcIGM_dir_path = os.path.join("..","ID_converter_ILGRMHD")
cmd.mkdir(IDcIGM_dir_path)
IDcIGM_src_dir_path = os.path.join(IDcIGM_dir_path,"src")
cmd.mkdir(IDcIGM_src_dir_path)

# Step 0b: Create the output file path
outfile_path__ID_converter_ILGRMHD__source    = os.path.join(IDcIGM_src_dir_path,"set_IllinoisGRMHD_metric_GRMHD_variables_based_on_HydroBase_and_ADMBase_variables.C")
outfile_path__ID_converter_ILGRMHD__make      = os.path.join(IDcIGM_src_dir_path,"make.code.defn")
outfile_path__ID_converter_ILGRMHD__param     = os.path.join(IDcIGM_dir_path,"param.ccl")
outfile_path__ID_converter_ILGRMHD__interface = os.path.join(IDcIGM_dir_path,"interface.ccl")
outfile_path__ID_converter_ILGRMHD__schedule  = os.path.join(IDcIGM_dir_path,"schedule.ccl")


# <a id='introduction'></a>
# 
# # Step 1: Introduction \[Back to [top](#toc)\]
# $$\label{introduction}$$

# <a id='convert_to_hydrobase__src'></a>
# 
# # Step 2: `set_IllinoisGRMHD_metric_GRMHD_variables _based_on_HydroBase_and_ADMBase_variables.C` \[Back to [top](#toc)\]
# $$\label{convert_to_hydrobase__src}$$

# In[2]:


get_ipython().run_cell_magic('writefile', '$outfile_path__ID_converter_ILGRMHD__source', '/********************************\n * CONVERT ET ID TO IllinoisGRMHD\n *\n * Written in 2014 by Zachariah B. Etienne\n *\n * Sets metric & MHD variables needed\n * by IllinoisGRMHD, converting from\n * HydroBase and ADMBase.\n ********************************/\n\n#include <stdio.h>\n#include <stdlib.h>\n#include <math.h>\n#include <sys/time.h>\n#include "cctk.h"\n#include "cctk_Parameters.h"\n#include "cctk_Arguments.h"\n#include "IllinoisGRMHD_headers.h"\n\nextern "C" void set_IllinoisGRMHD_metric_GRMHD_variables_based_on_HydroBase_and_ADMBase_variables(CCTK_ARGUMENTS) {\n\n  DECLARE_CCTK_ARGUMENTS;\n  DECLARE_CCTK_PARAMETERS;\n\n  if(rho_b_atm > 1e199) {\n    CCTK_VError(VERR_DEF_PARAMS, "You MUST set rho_b_atm to some reasonable value in your param.ccl file.\\n");\n  }\n\n  // Convert ADM variables (from ADMBase) to the BSSN-based variables expected by this routine.\n  IllinoisGRMHD_convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij(cctkGH,cctk_lsh,  gxx,gxy,gxz,gyy,gyz,gzz,alp,\n                                                                              gtxx,gtxy,gtxz,gtyy,gtyz,gtzz,\n                                                                              gtupxx,gtupxy,gtupxz,gtupyy,gtupyz,gtupzz,\n                                                                              phi_bssn,psi_bssn,lapm1);\n\n  /***************\n   * PPEOS Patch *\n   ***************/\n  eos_struct eos;\n  initialize_EOS_struct_from_input(eos);\n\n  if(pure_hydro_run) {\n#pragma omp parallel for\n    for(int k=0;k<cctk_lsh[2];k++) for(int j=0;j<cctk_lsh[1];j++) for(int i=0;i<cctk_lsh[0];i++) {\n          int index=CCTK_GFINDEX3D(cctkGH,i,j,k);\n          Avec[CCTK_GFINDEX4D(cctkGH,i,j,k,0)]=0;\n          Avec[CCTK_GFINDEX4D(cctkGH,i,j,k,1)]=0;\n          Avec[CCTK_GFINDEX4D(cctkGH,i,j,k,2)]=0;\n          Aphi[index]=0;\n        }\n  }\n\n#pragma omp parallel for\n  for(int k=0;k<cctk_lsh[2];k++) for(int j=0;j<cctk_lsh[1];j++) for(int i=0;i<cctk_lsh[0];i++) {\n        int index=CCTK_GFINDEX3D(cctkGH,i,j,k);\n\n        rho_b[index] = rho[index];\n        P[index] = press[index];\n\n        /***************\n         * PPEOS Patch *\n         ***************\n         * We now verify that the initial data\n         * provided by the user is indeed "cold",\n         * i.e. it contains no Thermal part and\n         * P = P_cold.\n         */\n        /* Compute P_cold */\n        const int polytropic_index = find_polytropic_K_and_Gamma_index(eos, rho_b[index]);\n        const double K_poly     = eos.K_ppoly_tab[polytropic_index];\n        const double Gamma_poly = eos.Gamma_ppoly_tab[polytropic_index];\n        const double P_cold     = K_poly*pow(rho_b[index],Gamma_poly);\n\n        /* Compare P and P_cold */\n        double P_rel_error = fabs(P[index] - P_cold)/P[index];\n        if( rho_b[index] > rho_b_atm && P_rel_error > 1e-2 ) {\n          const double Gamma_poly_local = log(P[index]/K_poly) / log(rho_b[index]);\n          /* Determine the value of Gamma_poly_local associated with P[index] */\n          CCTK_VWarn(CCTK_WARN_ALERT, __LINE__, __FILE__, CCTK_THORNSTRING,\n"Expected a PP EOS with local Gamma_poly = %.15e, but found a point such that Gamma_poly_local = %.15e.\\n",\n                     Gamma_poly, Gamma_poly_local);\n          CCTK_VWarn(CCTK_WARN_ALERT, __LINE__, __FILE__, CCTK_THORNSTRING,\n"{rho_b; rho_b_atm; P; P_cold; P_rel_Error} = %.10e %e %.10e %.10e %e\\n",\n                      rho_b[index], rho_b_atm, P[index],P_cold,P_rel_error);\n        }\n\n        Ax[index] = Avec[CCTK_GFINDEX4D(cctkGH,i,j,k,0)];\n        Ay[index] = Avec[CCTK_GFINDEX4D(cctkGH,i,j,k,1)];\n        Az[index] = Avec[CCTK_GFINDEX4D(cctkGH,i,j,k,2)];\n        psi6phi[index] = Aphi[index];\n\n        double ETvx = vel[CCTK_GFINDEX4D(cctkGH,i,j,k,0)];\n        double ETvy = vel[CCTK_GFINDEX4D(cctkGH,i,j,k,1)];\n        double ETvz = vel[CCTK_GFINDEX4D(cctkGH,i,j,k,2)];\n\n        // IllinoisGRMHD defines v^i = u^i/u^0.\n\n        // Meanwhile, the ET/HydroBase formalism, called the Valencia\n        // formalism, splits the 4 velocity into a purely spatial part\n        // and a part that is normal to the spatial hypersurface:\n        // u^a = G (n^a + U^a), (Eq. 14 of arXiv:1304.5544; G=W, U^a=v^a)\n        // where n^a is the unit normal vector to the spatial hypersurface,\n        // n_a = {-\\alpha,0,0,0}, and U^a is the purely spatial part, which\n        // is defined in HydroBase as the vel[] vector gridfunction.\n        // Then u^a n_a = - \\alpha u^0 = G n^a n_a = -G, and\n        // of course \\alpha u^0 = 1/sqrt(1+γ^ij u_i u_j) = \\Gamma,\n        // the standard Lorentz factor.\n\n        // Note that n^i = - \\beta^i / \\alpha, so\n        // u^a = \\Gamma (n^a + U^a)\n        // -> u^i = \\Gamma ( U^i - \\beta^i / \\alpha )\n        // which implies\n        // v^i = u^i/u^0\n        //     = \\Gamma/u^0 ( U^i - \\beta^i / \\alpha ) <- \\Gamma = \\alpha u^0\n        //     = \\alpha ( U^i - \\beta^i / \\alpha )\n        //     = \\alpha U^i - \\beta^i\n\n        vx[index] = alp[index]*ETvx - betax[index];\n        vy[index] = alp[index]*ETvy - betay[index];\n        vz[index] = alp[index]*ETvz - betaz[index];\n\n      }\n\n  // Neat feature for debugging: Add a roundoff-error perturbation\n  //    to the initial data.\n  // Set random_pert variable to ~1e-14 for a random 15th digit\n  //    perturbation.\n  srand(random_seed); // Use srand() as rand() is thread-safe.\n  for(int k=0;k<cctk_lsh[2];k++)\n    for(int j=0;j<cctk_lsh[1];j++)\n      for(int i=0;i<cctk_lsh[0];i++) {\n        int index=CCTK_GFINDEX3D(cctkGH,i,j,k);\n        double pert = (random_pert*(double)rand() / RAND_MAX);\n        double one_plus_pert=(1.0+pert);\n        rho[index]*=one_plus_pert;\n        vx[index]*=one_plus_pert;\n        vy[index]*=one_plus_pert;\n        vz[index]*=one_plus_pert;\n\n        psi6phi[index]*=one_plus_pert;\n        Ax[index]*=one_plus_pert;\n        Ay[index]*=one_plus_pert;\n        Az[index]*=one_plus_pert;\n      }\n\n  // Next compute B & B_stagger from A_i. Note that this routine also depends on\n  //   the psi_bssn[] gridfunction being set to exp(phi).\n\n  double dxi = 1.0/CCTK_DELTA_SPACE(0);\n  double dyi = 1.0/CCTK_DELTA_SPACE(1);\n  double dzi = 1.0/CCTK_DELTA_SPACE(2);\n\n#pragma omp parallel for\n  for(int k=0;k<cctk_lsh[2];k++)\n    for(int j=0;j<cctk_lsh[1];j++)\n      for(int i=0;i<cctk_lsh[0];i++) {\n        // Look Mom, no if() statements!\n        int shiftedim1 = (i-1)*(i!=0); // This way, i=0 yields shiftedim1=0 and shiftedi=1, used below for our COPY boundary condition.\n        int shiftedi   = shiftedim1+1;\n\n        int shiftedjm1 = (j-1)*(j!=0);\n        int shiftedj   = shiftedjm1+1;\n\n        int shiftedkm1 = (k-1)*(k!=0);\n        int shiftedk   = shiftedkm1+1;\n\n        int index,indexim1,indexjm1,indexkm1;\n\n        int actual_index = CCTK_GFINDEX3D(cctkGH,i,j,k);\n\n        double Psi = psi_bssn[actual_index];\n        double Psim3 = 1.0/(Psi*Psi*Psi);\n\n        // For the lower boundaries, the following applies a "copy"\n        //    boundary condition on Bi_stagger where needed.\n        //    E.g., Bx_stagger(i,jmin,k) = Bx_stagger(i,jmin+1,k)\n        //    We find the copy BC works better than extrapolation.\n        // For the upper boundaries, we do the following copy:\n        //    E.g., Psi(imax+1,j,k)=Psi(imax,j,k)\n        /**************/\n        /* Bx_stagger */\n        /**************/\n\n        index    = CCTK_GFINDEX3D(cctkGH,i,shiftedj,shiftedk);\n        indexjm1 = CCTK_GFINDEX3D(cctkGH,i,shiftedjm1,shiftedk);\n        indexkm1 = CCTK_GFINDEX3D(cctkGH,i,shiftedj,shiftedkm1);\n        // Set Bx_stagger = \\partial_y A_z - partial_z A_y\n        // "Grid" Ax(i,j,k) is actually Ax(i,j+1/2,k+1/2)\n        // "Grid" Ay(i,j,k) is actually Ay(i+1/2,j,k+1/2)\n        // "Grid" Az(i,j,k) is actually Ay(i+1/2,j+1/2,k)\n        // Therefore, the 2nd order derivative \\partial_z A_y at (i+1/2,j,k) is:\n        //          ["Grid" Ay(i,j,k) - "Grid" Ay(i,j,k-1)]/dZ\n        Bx_stagger[actual_index] = (Az[index]-Az[indexjm1])*dyi - (Ay[index]-Ay[indexkm1])*dzi;\n\n        // Now multiply Bx and Bx_stagger by 1/sqrt(gamma(i+1/2,j,k)]) = 1/sqrt(1/2 [gamma + gamma_ip1]) = exp(-6 x 1/2 [phi + phi_ip1] )\n        int imax_minus_i = (cctk_lsh[0]-1)-i;\n        int indexip1jk = CCTK_GFINDEX3D(cctkGH,i + ( (imax_minus_i > 0) - (0 > imax_minus_i) ),j,k);\n        double Psi_ip1 = psi_bssn[indexip1jk];\n        Bx_stagger[actual_index] *= Psim3/(Psi_ip1*Psi_ip1*Psi_ip1);\n\n        /**************/\n        /* By_stagger */\n        /**************/\n\n        index    = CCTK_GFINDEX3D(cctkGH,shiftedi,j,shiftedk);\n        indexim1 = CCTK_GFINDEX3D(cctkGH,shiftedim1,j,shiftedk);\n        indexkm1 = CCTK_GFINDEX3D(cctkGH,shiftedi,j,shiftedkm1);\n        // Set By_stagger = \\partial_z A_x - \\partial_x A_z\n        By_stagger[actual_index] = (Ax[index]-Ax[indexkm1])*dzi - (Az[index]-Az[indexim1])*dxi;\n\n        // Now multiply By and By_stagger by 1/sqrt(gamma(i,j+1/2,k)]) = 1/sqrt(1/2 [gamma + gamma_jp1]) = exp(-6 x 1/2 [phi + phi_jp1] )\n        int jmax_minus_j = (cctk_lsh[1]-1)-j;\n        int indexijp1k = CCTK_GFINDEX3D(cctkGH,i,j + ( (jmax_minus_j > 0) - (0 > jmax_minus_j) ),k);\n        double Psi_jp1 = psi_bssn[indexijp1k];\n        By_stagger[actual_index] *= Psim3/(Psi_jp1*Psi_jp1*Psi_jp1);\n\n\n        /**************/\n        /* Bz_stagger */\n        /**************/\n\n        index    = CCTK_GFINDEX3D(cctkGH,shiftedi,shiftedj,k);\n        indexim1 = CCTK_GFINDEX3D(cctkGH,shiftedim1,shiftedj,k);\n        indexjm1 = CCTK_GFINDEX3D(cctkGH,shiftedi,shiftedjm1,k);\n        // Set Bz_stagger = \\partial_x A_y - \\partial_y A_x\n        Bz_stagger[actual_index] = (Ay[index]-Ay[indexim1])*dxi - (Ax[index]-Ax[indexjm1])*dyi;\n\n        // Now multiply Bz_stagger by 1/sqrt(gamma(i,j,k+1/2)]) = 1/sqrt(1/2 [gamma + gamma_kp1]) = exp(-6 x 1/2 [phi + phi_kp1] )\n        int kmax_minus_k = (cctk_lsh[2]-1)-k;\n        int indexijkp1 = CCTK_GFINDEX3D(cctkGH,i,j,k + ( (kmax_minus_k > 0) - (0 > kmax_minus_k) ));\n        double Psi_kp1 = psi_bssn[indexijkp1];\n        Bz_stagger[actual_index] *= Psim3/(Psi_kp1*Psi_kp1*Psi_kp1);\n\n      }\n\n#pragma omp parallel for\n  for(int k=0;k<cctk_lsh[2];k++)\n    for(int j=0;j<cctk_lsh[1];j++)\n      for(int i=0;i<cctk_lsh[0];i++) {\n        // Look Mom, no if() statements!\n        int shiftedim1 = (i-1)*(i!=0); // This way, i=0 yields shiftedim1=0 and shiftedi=1, used below for our COPY boundary condition.\n        int shiftedi   = shiftedim1+1;\n\n        int shiftedjm1 = (j-1)*(j!=0);\n        int shiftedj   = shiftedjm1+1;\n\n        int shiftedkm1 = (k-1)*(k!=0);\n        int shiftedk   = shiftedkm1+1;\n\n        int index,indexim1,indexjm1,indexkm1;\n\n        int actual_index = CCTK_GFINDEX3D(cctkGH,i,j,k);\n\n        // For the lower boundaries, the following applies a "copy"\n        //    boundary condition on Bi and Bi_stagger where needed.\n        //    E.g., Bx(imin,j,k) = Bx(imin+1,j,k)\n        //    We find the copy BC works better than extrapolation.\n        /******/\n        /* Bx */\n        /******/\n        index = CCTK_GFINDEX3D(cctkGH,shiftedi,j,k);\n        indexim1 = CCTK_GFINDEX3D(cctkGH,shiftedim1,j,k);\n        // Set Bx = 0.5 ( Bx_stagger + Bx_stagger_im1 )\n        // "Grid" Bx_stagger(i,j,k) is actually Bx_stagger(i+1/2,j,k)\n        Bx[actual_index] = 0.5 * ( Bx_stagger[index] + Bx_stagger[indexim1] );\n\n        /******/\n        /* By */\n        /******/\n        index = CCTK_GFINDEX3D(cctkGH,i,shiftedj,k);\n        indexjm1 = CCTK_GFINDEX3D(cctkGH,i,shiftedjm1,k);\n        // Set By = 0.5 ( By_stagger + By_stagger_im1 )\n        // "Grid" By_stagger(i,j,k) is actually By_stagger(i,j+1/2,k)\n        By[actual_index] = 0.5 * ( By_stagger[index] + By_stagger[indexjm1] );\n\n        /******/\n        /* Bz */\n        /******/\n        index = CCTK_GFINDEX3D(cctkGH,i,j,shiftedk);\n        indexkm1 = CCTK_GFINDEX3D(cctkGH,i,j,shiftedkm1);\n        // Set Bz = 0.5 ( Bz_stagger + Bz_stagger_im1 )\n        // "Grid" Bz_stagger(i,j,k) is actually Bz_stagger(i,j+1/2,k)\n        Bz[actual_index] = 0.5 * ( Bz_stagger[index] + Bz_stagger[indexkm1] );\n      }\n\n  // Finally, enforce limits on primitives & compute conservative variables.\n#pragma omp parallel for\n  for(int k=0;k<cctk_lsh[2];k++)\n    for(int j=0;j<cctk_lsh[1];j++)\n      for(int i=0;i<cctk_lsh[0];i++) {\n        static const int zero_int=0;\n        int index = CCTK_GFINDEX3D(cctkGH,i,j,k);\n\n        int ww;\n\n        double PRIMS[MAXNUMVARS];\n        ww=0;\n        PRIMS[ww] = rho_b[index]; ww++;\n        PRIMS[ww] = P[index];     ww++;\n        PRIMS[ww] = vx[index];    ww++;\n        PRIMS[ww] = vy[index];    ww++;\n        PRIMS[ww] = vz[index];    ww++;\n        PRIMS[ww] = Bx[index];    ww++;\n        PRIMS[ww] = By[index];    ww++;\n        PRIMS[ww] = Bz[index];    ww++;\n\n        double METRIC[NUMVARS_FOR_METRIC],dummy=0;\n        ww=0;\n        // FIXME: NECESSARY?\n        //psi_bssn[index] = exp(phi[index]);\n        METRIC[ww] = phi_bssn[index];ww++;\n        METRIC[ww] = dummy;          ww++; // Don\'t need to set psi.\n        METRIC[ww] = gtxx[index];    ww++;\n        METRIC[ww] = gtxy[index];    ww++;\n        METRIC[ww] = gtxz[index];    ww++;\n        METRIC[ww] = gtyy[index];    ww++;\n        METRIC[ww] = gtyz[index];    ww++;\n        METRIC[ww] = gtzz[index];    ww++;\n        METRIC[ww] = lapm1[index];   ww++;\n        METRIC[ww] = betax[index];   ww++;\n        METRIC[ww] = betay[index];   ww++;\n        METRIC[ww] = betaz[index];   ww++;\n        METRIC[ww] = gtupxx[index];  ww++;\n        METRIC[ww] = gtupyy[index];  ww++;\n        METRIC[ww] = gtupzz[index];  ww++;\n        METRIC[ww] = gtupxy[index];  ww++;\n        METRIC[ww] = gtupxz[index];  ww++;\n        METRIC[ww] = gtupyz[index];  ww++;\n\n        double CONSERVS[NUM_CONSERVS] = {0,0,0,0,0};\n        double g4dn[4][4];\n        double g4up[4][4];\n        double TUPMUNU[10],TDNMUNU[10];\n\n        struct output_stats stats; stats.failure_checker=0;\n        IllinoisGRMHD_enforce_limits_on_primitives_and_recompute_conservs(zero_int,PRIMS,stats,eos,\n                                                                          METRIC,g4dn,g4up,TUPMUNU,TDNMUNU,CONSERVS);\n        rho_b[index] = PRIMS[RHOB];\n        P[index]     = PRIMS[PRESSURE];\n        vx[index]    = PRIMS[VX];\n        vy[index]    = PRIMS[VY];\n        vz[index]    = PRIMS[VZ];\n\n        rho_star[index] = CONSERVS[RHOSTAR];\n        mhd_st_x[index] = CONSERVS[STILDEX];\n        mhd_st_y[index] = CONSERVS[STILDEY];\n        mhd_st_z[index] = CONSERVS[STILDEZ];\n        tau[index]      = CONSERVS[TAUENERGY];\n\n        if(update_Tmunu) {\n          ww=0;\n          eTtt[index] = TDNMUNU[ww]; ww++;\n          eTtx[index] = TDNMUNU[ww]; ww++;\n          eTty[index] = TDNMUNU[ww]; ww++;\n          eTtz[index] = TDNMUNU[ww]; ww++;\n          eTxx[index] = TDNMUNU[ww]; ww++;\n          eTxy[index] = TDNMUNU[ww]; ww++;\n          eTxz[index] = TDNMUNU[ww]; ww++;\n          eTyy[index] = TDNMUNU[ww]; ww++;\n          eTyz[index] = TDNMUNU[ww]; ww++;\n          eTzz[index] = TDNMUNU[ww];\n        }\n      }\n}\n\n')


# <a id='convert_to_hydrobase__param'></a>
# 
# # Step 3: `param.ccl` \[Back to [top](#toc)\]
# $$\label{convert_to_hydrobase__param}$$

# In[3]:


get_ipython().run_cell_magic('writefile', '$outfile_path__ID_converter_ILGRMHD__param', '# Parameter definitions for thorn ID_converter_ILGRMHD\n\nshares: IllinoisGRMHD\nUSES KEYWORD rho_b_max\nUSES KEYWORD rho_b_atm\nUSES KEYWORD tau_atm\nUSES KEYWORD neos\nUSES KEYWORD K_ppoly_tab0\nUSES KEYWORD rho_ppoly_tab_in[10]\nUSES KEYWORD Gamma_ppoly_tab_in[10]\nUSES KEYWORD Sym_Bz\nUSES KEYWORD GAMMA_SPEED_LIMIT\nUSES KEYWORD Psi6threshold\nUSES KEYWORD update_Tmunu\n\nprivate:\n\nINT random_seed "Random seed for random, generally roundoff-level perturbation on initial data. Seeds srand(), and rand() is used for the RNG."\n{\n 0:99999999 :: "Anything unsigned goes."\n} 0\n\nREAL random_pert "Random perturbation atop data"\n{\n *:* :: "Anything goes."\n} 0\n\nBOOLEAN pure_hydro_run "Set the vector potential and corresponding EM gauge quantity to zero"\n{\n} "no"\n\n')


# <a id='convert_to_hydrobase__interface'></a>
# 
# # Step 4: `interface.ccl` \[Back to [top](#toc)\]
# $$\label{convert_to_hydrobase__interface}$$

# In[4]:


get_ipython().run_cell_magic('writefile', '$outfile_path__ID_converter_ILGRMHD__interface', '# Interface definition for thorn ID_converter_ILGRMHD\n\nimplements: ID_converter_ILGRMHD\n\ninherits: ADMBase, Boundary, SpaceMask, Tmunubase, HydroBase, grid, IllinoisGRMHD\n\nuses include header: IllinoisGRMHD_headers.h\nUSES INCLUDE: Symmetry.h\n\n')


# <a id='convert_to_hydrobase__schedule'></a>
# 
# # Step 5: `schedule.ccl` \[Back to [top](#toc)\]
# $$\label{convert_to_hydrobase__schedule}$$

# In[5]:


get_ipython().run_cell_magic('writefile', '$outfile_path__ID_converter_ILGRMHD__schedule', '# Schedule definitions for thorn ID_converter_ILGRMHD\n\nschedule group IllinoisGRMHD_ID_Converter at CCTK_INITIAL after HydroBase_Initial before Convert_to_HydroBase\n{\n} "Translate ET-generated, HydroBase-compatible initial data and convert into variables used by IllinoisGRMHD"\n\nschedule set_IllinoisGRMHD_metric_GRMHD_variables_based_on_HydroBase_and_ADMBase_variables IN IllinoisGRMHD_ID_Converter as first_initialdata before TOV_Initial_Data\n{\n   LANG:       C\n   OPTIONS:    LOCAL\n   # What the heck, let\'s synchronize everything!\n   SYNC: IllinoisGRMHD::grmhd_primitives_Bi, IllinoisGRMHD::grmhd_primitives_Bi_stagger, IllinoisGRMHD::grmhd_primitives_allbutBi, IllinoisGRMHD::em_Ax,IllinoisGRMHD::em_Ay,IllinoisGRMHD::em_Az,IllinoisGRMHD::em_psi6phi,IllinoisGRMHD::grmhd_conservatives,IllinoisGRMHD::BSSN_quantities,ADMBase::metric,ADMBase::lapse,ADMBase::shift,ADMBase::curv\n} "Convert HydroBase initial data (ID) to ID that IllinoisGRMHD can read."\n\nschedule IllinoisGRMHD_InitSymBound IN IllinoisGRMHD_ID_Converter as third_initialdata after second_initialdata\n{\n  SYNC: IllinoisGRMHD::grmhd_conservatives,IllinoisGRMHD::em_Ax,IllinoisGRMHD::em_Ay,IllinoisGRMHD::em_Az,IllinoisGRMHD::em_psi6phi\n  LANG: C\n} "Schedule symmetries -- Actually just a placeholder function to ensure prolongation / processor syncs are done BEFORE the primitives solver."\n\nschedule IllinoisGRMHD_compute_B_and_Bstagger_from_A IN IllinoisGRMHD_ID_Converter as fourth_initialdata after third_initialdata\n{\n  SYNC: IllinoisGRMHD::grmhd_primitives_Bi, IllinoisGRMHD::grmhd_primitives_Bi_stagger\n  LANG: C\n} "Compute B and B_stagger from A"\n\nschedule IllinoisGRMHD_conserv_to_prims IN IllinoisGRMHD_ID_Converter as fifth_initialdata after fourth_initialdata\n{\n  LANG: C\n} "Compute primitive variables from conservatives. This is non-trivial, requiring a Newton-Raphson root-finder."\n\n')


# <a id='convert_to_hydrobase__make'></a>
# 
# # Step 6: `make.code.defn` \[Back to [top](#toc)\]
# $$\label{convert_to_hydrobase__make}$$

# In[6]:


get_ipython().run_cell_magic('writefile', '$outfile_path__ID_converter_ILGRMHD__make', '# Main make.code.defn file for thorn ID_converter_ILGRMHD\n\n# Source files in this directory\nSRCS =  set_IllinoisGRMHD_metric_GRMHD_variables_based_on_HydroBase_and_ADMBase_variables.C\n\n')


# <a id='code_validation'></a>
# 
# # Step n-1: 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.

# In[7]:


# # 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/A_i_rhs_no_gauge_terms.C"
# original_IGM_file_name = "A_i_rhs_no_gauge_terms-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__A_i_rhs_no_gauge_terms__C  = !diff $original_IGM_file_path $outfile_path__A_i_rhs_no_gauge_terms__C

# if Validation__A_i_rhs_no_gauge_terms__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 A_i_rhs_no_gauge_terms.C: PASSED!")
# else:
#     # If the validation fails, we keep the original IGM source code file
#     print("Validation test for A_i_rhs_no_gauge_terms.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__A_i_rhs_no_gauge_terms__C:
#         print(diff_line)


# <a id='latex_pdf_output'></a>
# 
# # Step n: 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__A_i_rhs_no_gauge_terms.pdf](Tutorial-IllinoisGRMHD__A_i_rhs_no_gauge_terms.pdf) (Note that clicking on this link may not work; you may need to open the PDF file through another means).

# In[8]:


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__A_i_rhs_no_gauge_terms.ipynb
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__A_i_rhs_no_gauge_terms.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__A_i_rhs_no_gauge_terms.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__A_i_rhs_no_gauge_terms.tex
get_ipython().system('rm -f Tut*.out Tut*.aux Tut*.log')