#!/usr/bin/env python
# coding: utf-8
#
#
#
# # Tutorial-IllinoisGRMHD: MoL_registration.C
#
# ## Authors: Leo Werneck & Zach Etienne
#
# **This module is currently under development**
#
# ## This tutorial module explains the Method of Lines (MoL) options used by the `IllinoisGRMHD` codes.
#
# ### Note: This module will likely be absorbed by another one once we finish documenting the code.
#
# ### 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)).
#
#
# # 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](#mol_registration__c): **`MoL_registration.C`**
# 1. [Step n-1](#code_validation): **Code validation**
# 1. [Step n](#latex_pdf_output): **Output this notebook to $\LaTeX$-formatted PDF file**
#
#
# # 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: 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__MoL_registration__C = os.path.join(IGM_src_dir_path,"MoL_registration.C")
#
#
# # Step 1: Introduction \[Back to [top](#toc)\]
# $$\label{introduction}$$
#
#
# # Step 2: `MoL_registration.C` \[Back to [top](#toc)\]
# $$\label{mol_registration__c}$$
# In[2]:
get_ipython().run_cell_magic('writefile', '$outfile_path__MoL_registration__C', '//--------------------------------------------------------------------------\n// Register with the time stepper\n// (MoL thorn, found in arrangements/CactusBase/MoL)\n// To understand this, read documentation in arrangements/CactusBase/MoL/doc\n//--------------------------------------------------------------------------\n\n#include "cctk.h"\n#include \n#include \n#include \n#include "cctk_Parameters.h"\n#include "cctk_Arguments.h"\n\n#include "Symmetry.h"\n\nextern "C" void IllinoisGRMHD_RegisterVars(CCTK_ARGUMENTS)\n{\n DECLARE_CCTK_ARGUMENTS;\n DECLARE_CCTK_PARAMETERS;\n\n CCTK_INT ierr = 0, group, rhs;\n\n // Register evolution & RHS gridfunction groups\n\n /* Ax and Ax_rhs */\n group = CCTK_GroupIndex("IllinoisGRMHD::em_Ax");\n rhs = CCTK_GroupIndex("IllinoisGRMHD::em_Ax_rhs");\n ierr += MoLRegisterEvolvedGroup(group, rhs);\n\n /* Ay and Ay_rhs */\n group = CCTK_GroupIndex("IllinoisGRMHD::em_Ay");\n rhs = CCTK_GroupIndex("IllinoisGRMHD::em_Ay_rhs");\n ierr += MoLRegisterEvolvedGroup(group, rhs);\n\n /* Az and Az_rhs */\n group = CCTK_GroupIndex("IllinoisGRMHD::em_Az");\n rhs = CCTK_GroupIndex("IllinoisGRMHD::em_Az_rhs");\n ierr += MoLRegisterEvolvedGroup(group, rhs);\n\n /* psi6phi and psi6phi_rhs */\n group = CCTK_GroupIndex("IllinoisGRMHD::em_psi6phi");\n rhs = CCTK_GroupIndex("IllinoisGRMHD::em_psi6phi_rhs");\n ierr += MoLRegisterEvolvedGroup(group, rhs);\n\n /* ALL OTHER EVOLVED VARIABLES (rho_star,tau,mhd_st_x,mhd_st_y,mhd_st_z) */\n group = CCTK_GroupIndex("IllinoisGRMHD::grmhd_conservatives");\n rhs = CCTK_GroupIndex("IllinoisGRMHD::grmhd_conservatives_rhs");\n ierr += MoLRegisterEvolvedGroup(group, rhs);\n\n if (ierr) CCTK_ERROR("Problems registering with MoL");\n //***********************************************\n\n //***********************************************\n // Next register ADMBase variables needed by\n // IllinoisGRMHD as SaveAndRestore, so that\n // they are not set to NaN at the start of\n // each timestep (requiring that they be\n // e.g., recomputed from BSSN variables\n // in the BSSN solver, like Baikal or\n // ML_BSSN)\n ierr += MoLRegisterSaveAndRestoreGroup(CCTK_GroupIndex("admbase::lapse"));\n ierr += MoLRegisterSaveAndRestoreGroup(CCTK_GroupIndex("admbase::shift"));\n ierr += MoLRegisterSaveAndRestoreGroup(CCTK_GroupIndex("admbase::metric"));\n ierr += MoLRegisterSaveAndRestoreGroup(CCTK_GroupIndex("admbase::curv"));\n if (ierr) CCTK_ERROR("Problems registering with MoLRegisterSaveAndRestoreGroup");\n}\n')
#
#
# # 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[3]:
# 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/MoL_registration.C"
original_IGM_file_name = "MoL_registration-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
get_ipython().system('wget -O $original_IGM_file_path $original_IGM_file_url')
# Perform validation
Validation__MoL_registration__C = get_ipython().getoutput('diff $original_IGM_file_path $outfile_path__MoL_registration__C')
if Validation__MoL_registration__C == []:
# If the validation passes, we do not need to store the original IGM source code file
get_ipython().system('rm $original_IGM_file_path')
print("Validation test for MoL_registration.C: PASSED!")
else:
# If the validation fails, we keep the original IGM source code file
print("Validation test for MoL_registration.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__MoL_registration__C:
print(diff_line)
#
#
# # 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__MoL_registration.pdf](Tutorial-IllinoisGRMHD__MoL_registration.pdf) (Note that clicking on this link may not work; you may need to open the PDF file through another means).
# In[4]:
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__MoL_registration.ipynb
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__MoL_registration.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__MoL_registration.tex
#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__MoL_registration.tex
get_ipython().system('rm -f Tut*.out Tut*.aux Tut*.log')