This module is currently under development
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__postpostinitial__set_symmetries__copy_timelevels__C = os.path.join(IGM_src_dir_path,"postpostinitial__set_symmetries__copy_timelevels.C")
%%writefile $outfile_path__postpostinitial__set_symmetries__copy_timelevels__C
//-------------------------------------------------
// Stuff to run right after initial data is set up
//-------------------------------------------------
#include "cctk.h"
#include <cstdio>
#include <cstdlib>
#include "cctk_Arguments.h"
#include "cctk_Functions.h"
#include "cctk_Parameters.h"
#include "Symmetry.h"
#include "IllinoisGRMHD_headers.h"
#include "IllinoisGRMHD_EoS_lowlevel_functs.C"
extern "C" void IllinoisGRMHD_PostPostInitial_Set_Symmetries__Copy_Timelevels(CCTK_ARGUMENTS) {
DECLARE_CCTK_ARGUMENTS;
DECLARE_CCTK_PARAMETERS;
/**********************************
* Piecewise Polytropic EOS Patch *
* Printing the EOS table *
**********************************/
/*
* The short piece of code below takes care
* of initializing the EOS parameters.
* Please refer to the "inlined_functions.C"
* source file for the documentation on the
* function.
*/
eos_struct eos;
initialize_EOS_struct_from_input(eos);
/* For diagnostic and user convenience purposes, we print
* out the EOS parameters (rho_ppoly_tab, K_ppoly_tab,
* and Gamma_ppoly_tab) at t=0.
*/
if(cctk_iteration==0 && (int)GetRefinementLevel(cctkGH)==0) { print_EOS_table(eos); }
if(Gamma_th<0)
CCTK_VError(VERR_DEF_PARAMS,"ERROR. Default Gamma_th (=-1) detected. You must set Gamma_th to the appropriate value in your initial data thorn, or your .par file!\n");
//For emfields, we assume that you've set Bx, By, Bz (the UN-tilded B^i's)
// or Ax, Ay, Az (if using constrained transport scheme of Del Zanna)
if(CCTK_EQUALS(Symmetry,"equatorial")) {
// SET SYMMETRY GHOSTZONES ON ALL CONSERVATIVE AND PRIMIIVE VARIABLES!
int ierr;
ierr=CartSymGN(cctkGH,"IllinoisGRMHD::grmhd_conservatives"); if(ierr!=0) CCTK_VError(VERR_DEF_PARAMS,"Microsoft error code #1874109358120048. Grep it in the source code");
ierr=CartSymGN(cctkGH,"IllinoisGRMHD::grmhd_primitives_allbutBi"); if(ierr!=0) CCTK_VError(VERR_DEF_PARAMS,"Microsoft error code #1874109358120049. Grep it in the source code");
// Finish up by setting symmetry ghostzones on Bx, By, Bz, and their staggered variants.
CCTK_REAL gridfunc_syms_Bx[3] = {-1, 1,-Sym_Bz};
IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Bx , gridfunc_syms_Bx,0,0,0);
IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Bx_stagger, gridfunc_syms_Bx,1,0,0);
CCTK_REAL gridfunc_syms_By[3] = { 1,-1,-Sym_Bz};
IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, By , gridfunc_syms_Bx,0,0,0);
IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, By_stagger, gridfunc_syms_By,0,1,0);
CCTK_REAL gridfunc_syms_Bz[3] = { 1, 1, Sym_Bz};
IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Bz , gridfunc_syms_Bz,0,0,0);
IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Bz_stagger, gridfunc_syms_Bz,0,0,1);
CCTK_REAL gridfunc_syms_psi6phi[3] = { 1, 1, 1};
IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z,psi6phi , gridfunc_syms_psi6phi,1,1,1);
CCTK_REAL gridfunc_syms_Ax[3] = {-1, 1, Sym_Bz};
IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Ax , gridfunc_syms_Ax,0,1,1);
CCTK_REAL gridfunc_syms_Ay[3] = { 1,-1, Sym_Bz};
IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Ay , gridfunc_syms_Ay,1,0,1);
CCTK_REAL gridfunc_syms_Az[3] = { 1, 1,-Sym_Bz};
IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Az , gridfunc_syms_Az,1,1,0);
}
//------------------------------------------------------------------
// FILL _p AND _p_p TIMELEVELS. Probably don't need to do this if
// Carpet::init_fill_timelevels=yes and
// MoL::initial_data_is_crap = yes
// NOTE: We don't fill metric data here.
// FIXME: Do we really need this?
#pragma omp parallel for
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++) {
int index = CCTK_GFINDEX3D(cctkGH,i,j,k);
rho_star_p[index] = rho_star[index];
tau_p[index] = tau[index];
mhd_st_x_p[index] = mhd_st_x[index];
mhd_st_y_p[index] = mhd_st_y[index];
mhd_st_z_p[index] = mhd_st_z[index];
psi6phi_p[index] = psi6phi[index];
Ax_p[index] = Ax[index];
Ay_p[index] = Ay[index];
Az_p[index] = Az[index];
rho_star_p_p[index] = rho_star[index];
tau_p_p[index] = tau[index];
mhd_st_x_p_p[index] = mhd_st_x[index];
mhd_st_y_p_p[index] = mhd_st_y[index];
mhd_st_z_p_p[index] = mhd_st_z[index];
psi6phi_p_p[index] = psi6phi[index];
Ax_p_p[index] = Ax[index];
Ay_p_p[index] = Ay[index];
Az_p_p[index] = Az[index];
}
}
Writing ../src/postpostinitial__set_symmetries__copy_timelevels.C
# 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/postpostinitial__set_symmetries__copy_timelevels.C"
original_IGM_file_name = "postpostinitial__set_symmetries__copy_timelevels-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__postpostinitial__set_symmetries__copy_timelevels__C = !diff $original_IGM_file_path $outfile_path__postpostinitial__set_symmetries__copy_timelevels__C
if Validation__postpostinitial__set_symmetries__copy_timelevels__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 postpostinitial__set_symmetries__copy_timelevels.C: PASSED!")
else:
# If the validation fails, we keep the original IGM source code file
print("Validation test for postpostinitial__set_symmetries__copy_timelevels.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__postpostinitial__set_symmetries__copy_timelevels__C:
print(diff_line)
Validation test for postpostinitial__set_symmetries__copy_timelevels.C: FAILED! Diff: 12a13 > #include "IllinoisGRMHD_EoS_lowlevel_functs.C" 18,19c19,40 < if(gamma_th<0) < CCTK_VError(VERR_DEF_PARAMS,"ERROR. Default gamma_th (=-1) detected. You must set gamma_th to the appropriate value in your initial data thorn, or your .par file!\n"); --- > /********************************** > * Piecewise Polytropic EOS Patch * > * Printing the EOS table * > **********************************/ > /* > * The short piece of code below takes care > * of initializing the EOS parameters. > * Please refer to the "inlined_functions.C" > * source file for the documentation on the > * function. > */ > eos_struct eos; > initialize_EOS_struct_from_input(eos); > > /* For diagnostic and user convenience purposes, we print > * out the EOS parameters (rho_ppoly_tab, K_ppoly_tab, > * and Gamma_ppoly_tab) at t=0. > */ > if(cctk_iteration==0 && (int)GetRefinementLevel(cctkGH)==0) { print_EOS_table(eos); } > > if(Gamma_th<0) > CCTK_VError(VERR_DEF_PARAMS,"ERROR. Default Gamma_th (=-1) detected. You must set Gamma_th to the appropriate value in your initial data thorn, or your .par file!\n"); 53c74 < // FILL _p AND _p_p TIMELEVELS. Probably don't need to do this if --- > // FILL _p AND _p_p TIMELEVELS. Probably don't need to do this if 65c86 < mhd_st_x_p[index] = mhd_st_x[index]; --- > mhd_st_x_p[index] = mhd_st_x[index]; 70c91 < Ax_p[index] = Ax[index]; --- > Ax_p[index] = Ax[index]; 76c97 < mhd_st_x_p_p[index] = mhd_st_x[index]; --- > mhd_st_x_p_p[index] = mhd_st_x[index]; 81c102 < Ax_p_p[index] = Ax[index]; --- > Ax_p_p[index] = Ax[index]; 85a107 >
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__postpostinitial__set_symmetries__copy_timelevels.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__postpostinitial__set_symmetries__copy_timelevels.ipynb
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__postpostinitial__set_symmetries__copy_timelevels.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__postpostinitial__set_symmetries__copy_timelevels.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__postpostinitial__set_symmetries__copy_timelevels.tex
!rm -f Tut*.out Tut*.aux Tut*.log