{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"\n",
"# Tutorial-IllinoisGRMHD: A_i_rhs_no_gauge_terms.C\n",
"\n",
"## Authors: Leo Werneck & Zach Etienne\n",
"\n",
"**This module is currently under development**\n",
"\n",
"## This notebook explains the construction of the right-hand side of the evolution equations of $A_{i}$.\n",
"\n",
"### Required and recommended citations:\n",
"\n",
"* **(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)).\n",
"* **(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)).\n",
"* **(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))."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Table of Contents\n",
"$$\\label{toc}$$\n",
"\n",
"This module is organized as follows\n",
"\n",
"0. [Step 0](#src_dir): **Source directory creation**\n",
"1. [Step 1](#introduction): **Introduction**\n",
"1. [Step 2](#a_i_rhs_no_gauge_terms__c): **`A_i_rhs_no_gauge_terms.C`**\n",
"1. [Step n-1](#code_validation): **Code validation**\n",
"1. [Step n](#latex_pdf_output): **Output this notebook to $\\LaTeX$-formatted PDF file**"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 0: Source directory creation \\[Back to [top](#toc)\\]\n",
"$$\\label{src_dir}$$\n",
"\n",
"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."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"# Step 0: Creation of the IllinoisGRMHD source directory\n",
"# Step 0a: Add NRPy's directory to the path\n",
"# https://stackoverflow.com/questions/16780014/import-file-from-parent-directory\n",
"import os,sys\n",
"nrpy_dir_path = os.path.join(\"..\",\"..\")\n",
"if nrpy_dir_path not in sys.path:\n",
" sys.path.append(nrpy_dir_path)\n",
"\n",
"# Step 0b: Load up cmdline_helper and create the directory\n",
"import cmdline_helper as cmd\n",
"IGM_src_dir_path = os.path.join(\"..\",\"src\")\n",
"cmd.mkdir(IGM_src_dir_path)\n",
"\n",
"# Step 0c: Create the output file path\n",
"outfile_path__A_i_rhs_no_gauge_terms__C = os.path.join(IGM_src_dir_path,\"A_i_rhs_no_gauge_terms.C\")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 1: Introduction \\[Back to [top](#toc)\\]\n",
"$$\\label{introduction}$$"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step 2: `A_i_rhs_no_gauge_terms.C` \\[Back to [top](#toc)\\]\n",
"$$\\label{a_i_rhs_no_gauge_terms__c}$$"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Writing ../src/A_i_rhs_no_gauge_terms.C\n"
]
}
],
"source": [
"%%writefile $outfile_path__A_i_rhs_no_gauge_terms__C\n",
"/* Compute the part of A_i_rhs that excludes the gauge terms. I.e., we set\n",
" * A_i_rhs = \\partial_t A_i = \\psi^{6} (v^z B^x - v^x B^z) here.\n",
" */\n",
"static void A_i_rhs_no_gauge_terms(const int A_dirn,const cGH *cctkGH,const int *cctk_lsh,const int *cctk_nghostzones,gf_and_gz_struct *out_prims_r,gf_and_gz_struct *out_prims_l,\n",
" CCTK_REAL *phi_interped,CCTK_REAL *cmax_1,CCTK_REAL *cmin_1,CCTK_REAL *cmax_2,CCTK_REAL *cmin_2, CCTK_REAL *A3_rhs) {\n",
" // If A_dirn=1, then v1_offset=1 (v1=VY) and v2_offset=2 (v2=VZ)\n",
" // If A_dirn=2, then v1_offset=2 (v1=VZ) and v2_offset=0 (v2=VX)\n",
" // If A_dirn=3, then v1_offset=0 (v1=VX) and v2_offset=1 (v2=VY)\n",
" int v1_offset = ((A_dirn-1)+1)%3, v2_offset = ((A_dirn-1)+2)%3;\n",
"\n",
" CCTK_REAL *v1rr=out_prims_r[VXR+v1_offset].gf, *v2rr=out_prims_r[VXR+v2_offset].gf;\n",
" CCTK_REAL *v1rl=out_prims_l[VXR+v1_offset].gf, *v2rl=out_prims_l[VXR+v2_offset].gf;\n",
" CCTK_REAL *v1lr=out_prims_r[VXL+v1_offset].gf, *v2lr=out_prims_r[VXL+v2_offset].gf;\n",
" CCTK_REAL *v1ll=out_prims_l[VXL+v1_offset].gf, *v2ll=out_prims_l[VXL+v2_offset].gf;\n",
"\n",
" CCTK_REAL *B1r=out_prims_r[BX_STAGGER+v1_offset].gf, *B1l=out_prims_l[BX_STAGGER+v1_offset].gf;\n",
" CCTK_REAL *B2r=out_prims_r[BX_STAGGER+v2_offset].gf, *B2l=out_prims_l[BX_STAGGER+v2_offset].gf;\n",
"\n",
" /**** V DEPENDENCIES ****/\n",
" /* In the case of Ax_rhs, we need v{y,z}{r,l} at (i,j+1/2,k+1/2).\n",
" * However, v{y,z}{r,l}{r,l} are defined at (i,j-1/2,k-1/2), so\n",
" * v{y,z}{r,l} at (i,j+1/2,k+1/2) is stored at v{y,z}{r,l}{r,l}(i,j+1,k+1).\n",
" * In the case of Ay_rhs, we need v{x,z}{r,l} at (i+1/2,j,k+1/2).\n",
" * However, v{x,z}{r,l}{r,l} are defined at (i-1/2,j,k-1/2), so\n",
" * v{x,z}{r,l} at (i+1/2,j,k+1/2) is stored at v{x,z}{r,l}{r,l}(i+1,j,k+1).\n",
" * In the case of Az_rhs, we need v{x,y}{r,l} at (i+1/2,j+1/2,k).\n",
" * However, v{x,y}{r,l}{r,l} are defined at (i-1/2,j-1/2,k), so\n",
" * v{x,y}{r,l} at (i+1/2,j+1/2,k) is stored at v{x,y}{r,l}{r,l}(i+1,j+1,k). */\n",
" static const int vs_ijk_offset[4][3] = { {0,0,0} , {0,1,1} , {1,0,1} , {1,1,0} };\n",
"\n",
" /**** B DEPENDENCIES ****/\n",
" /* In the case of Ax_rhs, we need B{y,z}{r,l} at (i,j+1/2,k+1/2).\n",
" * However, By_stagger{r,l} is defined at (i,j+1/2,k-1/2), and\n",
" * Bz_stagger{r,l} is defined at (i,j-1/2,k+1/2), so\n",
" * By_stagger{r,l} at (i,j+1/2,k+1/2) is stored at By_stagger{r,l}(i,j,k+1), and\n",
" * Bz_stagger{r,l} at (i,j+1/2,k+1/2) is stored at Bz_stagger{r,l}(i,j+1,k).\n",
" * In the case of Ay_rhs, we need B{z,x}_stagger{r,l} at (i+1/2,j,k+1/2).\n",
" * However, Bz_stagger{r,l} is defined at (i-1/2,j,k+1/2), and\n",
" * Bx_stagger{r,l} is defined at (i+1/2,j,k-1/2), so\n",
" * Bz_stagger{r,l} at (i+1/2,j,k+1/2) is stored at Bz_stagger{r,l}(i+1,j,k), and\n",
" * Bx_stagger{r,l} at (i+1/2,j,k+1/2) is stored at Bx_stagger{r,l}(i,j,k+1).\n",
" * In the case of Az_rhs, we need B{x,y}_stagger{r,l} at (i+1/2,j+1/2,k).\n",
" * However, Bx_stagger{r,l} is defined at (i+1/2,j-1/2,k), and\n",
" * By_stagger{r,l} is defined at (i-1/2,j+1/2,k), so\n",
" * Bx_stagger{r,l} at (i+1/2,j+1/2,k) is stored at Bx_stagger{r,l}(i,j+1,k), and\n",
" * By_stagger{r,l} at (i+1/2,j+1/2,k) is stored at By_stagger{r,l}(i+1,j,k).\n",
" */\n",
" static const int B1_ijk_offset[4][3] = { {0,0,0} , {0,0,1} , {1,0,0} , {0,1,0} };\n",
" static const int B2_ijk_offset[4][3] = { {0,0,0} , {0,1,0} , {0,0,1} , {1,0,0} };\n",
"\n",
"#pragma omp parallel for\n",
" for(int k=cctk_nghostzones[2];k\n",
"\n",
"# Step n-1: Code validation \\[Back to [top](#toc)\\]\n",
"$$\\label{code_validation}$$\n",
"\n",
"First, we download the original `IllinoisGRMHD` source code and then compare it to the source code generated by this tutorial notebook."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Validation test for A_i_rhs_no_gauge_terms.C: FAILED!\n",
"Diff:\n",
"20c20\n",
"< /* In the case of Ax_rhs, we need v{y,z}{r,l} at (i,j+1/2,k+1/2). \n",
"---\n",
"> /* In the case of Ax_rhs, we need v{y,z}{r,l} at (i,j+1/2,k+1/2).\n",
"23c23\n",
"< * In the case of Ay_rhs, we need v{x,z}{r,l} at (i+1/2,j,k+1/2). \n",
"---\n",
"> * In the case of Ay_rhs, we need v{x,z}{r,l} at (i+1/2,j,k+1/2).\n",
"26c26\n",
"< * In the case of Az_rhs, we need v{x,y}{r,l} at (i+1/2,j+1/2,k). \n",
"---\n",
"> * In the case of Az_rhs, we need v{x,y}{r,l} at (i+1/2,j+1/2,k).\n",
"33c33\n",
"< * However, By_stagger{r,l} is defined at (i,j+1/2,k-1/2), and \n",
"---\n",
"> * However, By_stagger{r,l} is defined at (i,j+1/2,k-1/2), and\n",
"43c43\n",
"< * However, Bx_stagger{r,l} is defined at (i+1/2,j-1/2,k), and \n",
"---\n",
"> * However, Bx_stagger{r,l} is defined at (i+1/2,j-1/2,k), and\n",
"86c86\n",
"< * Implement 2D HLL flux \n",
"---\n",
"> * Implement 2D HLL flux\n",
"89c89\n",
"< * Note that cmax/cmin (\\alpha^{\\pm} as defined in that paper) is at a slightly DIFFERENT \n",
"---\n",
"> * Note that cmax/cmin (\\alpha^{\\pm} as defined in that paper) is at a slightly DIFFERENT\n",
"96,97c96,97\n",
"< /( (cmax_1L+cmin_1L)*(cmax_2L+cmin_2L) ) \n",
"< - cmax_1L*cmin_1L*(B2tilder_minus_B2tildel)/(cmax_1L+cmin_1L) \n",
"---\n",
"> /( (cmax_1L+cmin_1L)*(cmax_2L+cmin_2L) )\n",
"> - cmax_1L*cmin_1L*(B2tilder_minus_B2tildel)/(cmax_1L+cmin_1L)\n",
"100a101\n",
"> \n"
]
}
],
"source": [
"# Verify if the code generated by this tutorial module\n",
"# matches the original IllinoisGRMHD source code\n",
"\n",
"# First download the original IllinoisGRMHD source code\n",
"import urllib\n",
"from os import path\n",
"\n",
"original_IGM_file_url = \"https://bitbucket.org/zach_etienne/wvuthorns/raw/5611b2f0b17135538c9d9d17c7da062abe0401b6/IllinoisGRMHD/src/A_i_rhs_no_gauge_terms.C\"\n",
"original_IGM_file_name = \"A_i_rhs_no_gauge_terms-original.C\"\n",
"original_IGM_file_path = os.path.join(IGM_src_dir_path,original_IGM_file_name)\n",
"\n",
"# Then download the original IllinoisGRMHD source code\n",
"# We try it here in a couple of ways in an attempt to keep\n",
"# the code more portable\n",
"try:\n",
" original_IGM_file_code = urllib.request.urlopen(original_IGM_file_url).read().decode(\"utf-8\")\n",
" # Write down the file the original IllinoisGRMHD source code\n",
" with open(original_IGM_file_path,\"w\") as file:\n",
" file.write(original_IGM_file_code)\n",
"except:\n",
" try:\n",
" original_IGM_file_code = urllib.urlopen(original_IGM_file_url).read().decode(\"utf-8\")\n",
" # Write down the file the original IllinoisGRMHD source code\n",
" with open(original_IGM_file_path,\"w\") as file:\n",
" file.write(original_IGM_file_code)\n",
" except:\n",
" # If all else fails, hope wget does the job\n",
" !wget -O $original_IGM_file_path $original_IGM_file_url\n",
"\n",
"# Perform validation\n",
"Validation__A_i_rhs_no_gauge_terms__C = !diff $original_IGM_file_path $outfile_path__A_i_rhs_no_gauge_terms__C\n",
"\n",
"if Validation__A_i_rhs_no_gauge_terms__C == []:\n",
" # If the validation passes, we do not need to store the original IGM source code file\n",
" !rm $original_IGM_file_path\n",
" print(\"Validation test for A_i_rhs_no_gauge_terms.C: PASSED!\")\n",
"else:\n",
" # If the validation fails, we keep the original IGM source code file\n",
" print(\"Validation test for A_i_rhs_no_gauge_terms.C: FAILED!\")\n",
" # We also print out the difference between the code generated\n",
" # in this tutorial module and the original IGM source code\n",
" print(\"Diff:\")\n",
" for diff_line in Validation__A_i_rhs_no_gauge_terms__C:\n",
" print(diff_line)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"\n",
"\n",
"# Step n: Output this notebook to $\\LaTeX$-formatted PDF file \\[Back to [top](#toc)\\]\n",
"$$\\label{latex_pdf_output}$$\n",
"\n",
"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\n",
"[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)."
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"latex_nrpy_style_path = os.path.join(nrpy_dir_path,\"latex_nrpy_style.tplx\")\n",
"#!jupyter nbconvert --to latex --template $latex_nrpy_style_path --log-level='WARN' Tutorial-IllinoisGRMHD__A_i_rhs_no_gauge_terms.ipynb\n",
"#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__A_i_rhs_no_gauge_terms.tex\n",
"#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__A_i_rhs_no_gauge_terms.tex\n",
"#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__A_i_rhs_no_gauge_terms.tex\n",
"!rm -f Tut*.out Tut*.aux Tut*.log"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.10.11"
}
},
"nbformat": 4,
"nbformat_minor": 4
}