This module is currently under development
IllinoisGRMHD
header file.¶We will now use the cmdline_helper.py NRPy+ module to create the source directory within the IllinoisGRMHD
NRPy+ directory if it does not exist yet.
# 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
IGM_src_dir_path = os.path.join("..","src")
cmd.mkdir(IGM_src_dir_path)
# Step 0c: Create the output file path
outfile_path__IllinoisGRMHD_headers__h = os.path.join(IGM_src_dir_path,"IllinoisGRMHD_headers.h")
%%writefile $outfile_path__IllinoisGRMHD_headers__h
// To safeguard against double-including this header file:
#ifndef ILLINOISGRMHD_HEADERS_H_
#define ILLINOISGRMHD_HEADERS_H_
#define MIN(a,b) ( ((a) < (b)) ? (a) : (b) )
#define MAX(a,b) ( ((a) > (b)) ? (a) : (b) )
#define SQR(x) ((x) * (x))
#define ONE_OVER_SQRT_4PI 0.282094791773878143474039725780
#define VERR_DEF_PARAMS __LINE__, __FILE__, CCTK_THORNSTRING
// The order here MATTERS, as we assume that GUPXX+1=GUPYY, etc.
static const int PHI=0,PSI=1,GXX=2,GXY=3,GXZ=4,GYY=5,GYZ=6,GZZ=7,
LAPM1=8,SHIFTX=9,SHIFTY=10,SHIFTZ=11,GUPXX=12,GUPYY=13,GUPZZ=14,
NUMVARS_FOR_METRIC_FACEVALS=15; //<-- Be _sure_ to set this correctly, or you'll have memory access bugs!
// These are not used for facevals in the reconstruction step, but boy are they useful anyway.
static const int GUPXY=15,GUPXZ=16,GUPYZ=17,
NUMVARS_FOR_METRIC=18; //<-- Be _sure_ to set this correctly, or you'll have memory access bugs!
// The order here MATTERS, and must be consistent with the order in the in_prims[] array in driver_evaluate_MHD_rhs.C.
static const int RHOB=0,PRESSURE=1,VX=2,VY=3,VZ=4,
BX_CENTER=5,BY_CENTER=6,BZ_CENTER=7,BX_STAGGER=8,BY_STAGGER=9,BZ_STAGGER=10,
VXR=11,VYR=12,VZR=13,VXL=14,VYL=15,VZL=16,MAXNUMVARS=17; //<-- Be _sure_ to define MAXNUMVARS appropriately!
static const int UT=0,UX=1,UY=2,UZ=3;
// The "I" suffix denotes interpolation. In other words, these
// definitions are used for interpolation ONLY. The order here
// matters as well!
static const int SHIFTXI=0,SHIFTYI=1,SHIFTZI=2,GUPXXI=3,GUPXYI=4,GUPXZI=5,GUPYYI=6,GUPYZI=7,GUPZZI=8,
PSII=9,LAPM1I=10,A_XI=11,A_YI=12,A_ZI=13,LAPSE_PSI2I=14,LAPSE_OVER_PSI6I=15,MAXNUMINTERP=16;
// Again, the order here MATTERS, since we assume in the code that, e.g., smallb[0]=b^t, smallb[3]=b^z, etc.
static const int SMALLBT=0,SMALLBX=1,SMALLBY=2,SMALLBZ=3,SMALLB2=4,NUMVARS_SMALLB=5;
// Again, the order here MATTERS, since we assume in the code that, CONSERV[STILDEX+1] = \tilde{S}_y
static const int RHOSTAR=0,STILDEX=1,STILDEY=2,STILDEZ=3,TAUENERGY=4,NUM_CONSERVS=5;
static const int LAPSE=0,PSI2=1,PSI4=2,PSI6=3,PSIM4=4,LAPSEINV=5,NUMVARS_METRIC_AUX=6;
#define SET_LAPSE_PSI4(array_name,METRIC) { \
array_name[LAPSE] = METRIC[LAPM1]+1.0; \
array_name[PSI2] = exp(2.0*METRIC[PHI]); \
array_name[PSI4] = SQR(array_name[PSI2]); \
array_name[PSI6] = array_name[PSI4]*array_name[PSI2]; \
array_name[PSIM4] = 1.0/array_name[PSI4]; \
array_name[LAPSEINV] = 1.0/array_name[LAPSE]; \
}
// Keeping track of ghostzones between routines is a nightmare, so
// we instead attach ghostzone info to each gridfunction and set
// the ghostzone information correctly within each routine.
struct gf_and_gz_struct {
CCTK_REAL *gf;
int gz_lo[4],gz_hi[4];
};
#define MAX_EOS_PARAMS 10
struct eos_struct {
int neos;
CCTK_REAL rho_ppoly_tab[MAX_EOS_PARAMS-1];
CCTK_REAL eps_integ_const[MAX_EOS_PARAMS],K_ppoly_tab[MAX_EOS_PARAMS],Gamma_ppoly_tab[MAX_EOS_PARAMS];
};
struct output_stats {
int font_fixed,vel_limited,failure_checker;
long n_iter;
};
// FIXME: For cosmetic purposes, we might want to make everything either zero-offset or one-offset, instead of a mixture.
const int kronecker_delta[4][3] = { { 0,0,0 },
{ 1,0,0 },
{ 0,1,0 },
{ 0,0,1 } };
/* PUBLIC FUNCTIONS, USED OUTSIDE IllinoisGRMHD AS WELL */
void IllinoisGRMHD_enforce_limits_on_primitives_and_recompute_conservs(const int already_computed_physical_metric_and_inverse,CCTK_REAL *U,struct output_stats &stats,eos_struct &eos,
CCTK_REAL *METRIC,CCTK_REAL g4dn[4][4],CCTK_REAL g4up[4][4], CCTK_REAL *TUPMUNU,CCTK_REAL *TDNMUNU,CCTK_REAL *CONSERVS);
void IllinoisGRMHD_convert_ADM_to_BSSN__enforce_detgtij_eq_1__and_compute_gtupij
(const cGH *cctkGH,const int *cctk_lsh,
CCTK_REAL *gxx,CCTK_REAL *gxy,CCTK_REAL *gxz,CCTK_REAL *gyy,CCTK_REAL *gyz,CCTK_REAL *gzz,CCTK_REAL *alp,
CCTK_REAL *gtxx,CCTK_REAL *gtxy,CCTK_REAL *gtxz,CCTK_REAL *gtyy,CCTK_REAL *gtyz,CCTK_REAL *gtzz,
CCTK_REAL *gtupxx,CCTK_REAL *gtupxy,CCTK_REAL *gtupxz,CCTK_REAL *gtupyy,CCTK_REAL *gtupyz,CCTK_REAL *gtupzz,
CCTK_REAL *phi,CCTK_REAL *psi,CCTK_REAL *lapm1);
void IllinoisGRMHD_set_symmetry_gzs_staggered(const cGH *cctkGH, const int *cctk_lsh,CCTK_REAL *X,CCTK_REAL *Y,CCTK_REAL *Z, CCTK_REAL *gridfunc,
CCTK_REAL *gridfunc_syms,int stagger_x,int stagger_y,int stagger_z);
#include "IllinoisGRMHD_EoS_lowlevel_functs.C"
#endif // ILLINOISGRMHD_HEADERS_H
Writing ../src/IllinoisGRMHD_headers.h
# 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/IllinoisGRMHD_headers.h"
original_IGM_file_name = "IllinoisGRMHD_headers-original.h"
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__IllinoisGRMHD_headers__h = !diff $original_IGM_file_path $outfile_path__IllinoisGRMHD_headers__h
if Validation__IllinoisGRMHD_headers__h == []:
# 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 IllinoisGRMHD_headers.h: PASSED!")
else:
# If the validation fails, we keep the original IGM source code file
print("Validation test for IllinoisGRMHD_headers.h: 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__IllinoisGRMHD_headers__h:
print(diff_line)
Validation test for IllinoisGRMHD_headers.h: FAILED! Diff: 17c17 < // These are not used for facevals in the reconstruction step, but boy are they useful anyway. --- > // These are not used for facevals in the reconstruction step, but boy are they useful anyway. 61,64c61,63 < CCTK_REAL rho_tab[MAX_EOS_PARAMS],P_tab[MAX_EOS_PARAMS],eps_tab[MAX_EOS_PARAMS],k_tab[MAX_EOS_PARAMS],gamma_tab[MAX_EOS_PARAMS]; < CCTK_REAL gamma_th; < CCTK_REAL K_poly; < }; --- > CCTK_REAL rho_ppoly_tab[MAX_EOS_PARAMS-1]; > CCTK_REAL eps_integ_const[MAX_EOS_PARAMS],K_ppoly_tab[MAX_EOS_PARAMS],Gamma_ppoly_tab[MAX_EOS_PARAMS]; > }; 91a91 > #include "IllinoisGRMHD_EoS_lowlevel_functs.C" 92a93 >
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__IllinoisGRMHD_headers.pdf (Note that clicking on this link may not work; you may need to open the PDF file through another means).
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__IllinoisGRMHD_headers.ipynb
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__IllinoisGRMHD_headers.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__IllinoisGRMHD_headers.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__IllinoisGRMHD_headers.tex
!rm -f Tut*.out Tut*.aux Tut*.log