This module is currently under development
IllinoisGRMHD
handles symmetry options and the ghostzones of staggered gridfunctions. This module will likely be absorbed by another one once we finish documenting the code.¶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__symmetry__set_gzs_staggered_gfs__C = os.path.join(IGM_src_dir_path,"symmetry__set_gzs_staggered_gfs.C")
%%writefile $outfile_path__symmetry__set_gzs_staggered_gfs__C
#include "cctk.h"
#include "cctk_Parameters.h"
#include <cstdio>
#include <cstdlib>
#include "IllinoisGRMHD_headers.h"
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) {
DECLARE_CCTK_PARAMETERS;
if(CCTK_EQUALS(Symmetry, "equatorial"))
CCTK_VError(VERR_DEF_PARAMS,"Warning: Symmetry==equatorial not supported! USE AT YOUR OWN RISK. You will need to comment this error message out.");
// No symmetries -> return.
if(CCTK_EQUALS(Symmetry, "none")) return;
CCTK_REAL dz = Z[CCTK_GFINDEX3D(cctkGH,0,0,1)] - Z[CCTK_GFINDEX3D(cctkGH,0,0,0)];
CCTK_REAL z_offset = dz*0.5*stagger_z;
int num_gzs=0;
//FIXME: Might want to use cctk_nghostzones instead...
while( (Z[CCTK_GFINDEX3D(cctkGH,0,0,num_gzs)]+z_offset) < -dz*0.1 && num_gzs<cctk_lsh[2]) num_gzs++;
if(num_gzs*2>=cctk_lsh[2]) CCTK_VError(VERR_DEF_PARAMS,"ERROR in symmetry__set_gzs_staggered_gfs.C");
#pragma omp parallel for
for(int k=0;k<num_gzs;k++) for(int j=0;j<cctk_lsh[1];j++) for(int i=0;i<cctk_lsh[0];i++) {
int index_inside__sym_gz = CCTK_GFINDEX3D(cctkGH,i,j,k);
/* This loop sets symmetry ghostzones, regardless of how the gridfunction is staggered.
*
* STAGGERED PATTERN:
* if num_gzs==1 && stagger_z==1:
* z[] = {-dz/2,dz/2,3dz/2, etc} -> gridfunc[index 0] = gridfunc_syms[2]*gridfunc[index 1]
*
* if num_gzs==2 && stagger_z==1:
* z[] = {-3dz/2,-dz/2,dz/2,3dz/2 etc}
* -> gridfunc[index 0] = gridfunc_syms[2]*gridfunc[index 3]
* -> gridfunc[index 1] = gridfunc_syms[2]*gridfunc[index 2]
* .
* .
* .
* -> gridfunc[i] = gridfunc_syms[2]*gridfunc[(num_gz*2-1)-i]
*
* UNSTAGGERED PATTERN:
* if num_gzs==1 && stagger_z==0:
* z[] = {-dz,0,dz, etc} -> gridfunc[index 0] = gridfunc_syms[2]*gridfunc[index 2]
*
* if num_gzs==2 && stagger_z==0:
* z[] = {-2dz,-dz,0,dz,2dz, etc} -> gridfunc[index 0] = gridfunc_syms[2]*gridfunc[index 4]
* z[] = {-2dz,-dz,0,dz,2dz, etc} -> gridfunc[index 1] = gridfunc_syms[2]*gridfunc[index 3]
* .
* .
* .
* -> gridfunc[i] = gridfunc_syms[2]*gridfunc[(num_gz*2)-i]
*
* OVERALL PATTERN: gridfunc[i] = gridfunc_syms[2]*gridfunc[(num_gz*2-stagger_z)-i] */
int matching_index_outside_sym_gz = CCTK_GFINDEX3D(cctkGH,i,j,(num_gzs*2-stagger_z)-k);
gridfunc[index_inside__sym_gz] = gridfunc_syms[2]*gridfunc[matching_index_outside_sym_gz];
}
}
Writing ../src/symmetry__set_gzs_staggered_gfs.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/symmetry__set_gzs_staggered_gfs.C"
original_IGM_file_name = "symmetry__set_gzs_staggered_gfs-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__symmetry__set_gzs_staggered_gfs__C = !diff $original_IGM_file_path $outfile_path__symmetry__set_gzs_staggered_gfs__C
if Validation__symmetry__set_gzs_staggered_gfs__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 symmetry__set_gzs_staggered_gfs.C: PASSED!")
else:
# If the validation fails, we keep the original IGM source code file
print("Validation test for symmetry__set_gzs_staggered_gfs.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__symmetry__set_gzs_staggered_gfs__C:
print(diff_line)
Validation test for symmetry__set_gzs_staggered_gfs.C: FAILED! Diff: 12c12 < if(CCTK_EQUALS(Symmetry, "equatorial")) --- > if(CCTK_EQUALS(Symmetry, "equatorial")) 34c34 < * if num_gzs==1 && stagger_z==1: --- > * if num_gzs==1 && stagger_z==1: 37,38c37,38 < * if num_gzs==2 && stagger_z==1: < * z[] = {-3dz/2,-dz/2,dz/2,3dz/2 etc} --- > * if num_gzs==2 && stagger_z==1: > * z[] = {-3dz/2,-dz/2,dz/2,3dz/2 etc} 64a65 >
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__symmetry__set_gzs_staggered_gfs.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__symmetry__set_gzs_staggered_gfs.ipynb
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__symmetry__set_gzs_staggered_gfs.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__symmetry__set_gzs_staggered_gfs.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__symmetry__set_gzs_staggered_gfs.tex
!rm -f Tut*.out Tut*.aux Tut*.log