{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "\n", "# Tutorial-IllinoisGRMHD: Convert_to_HydroBase ETKThorn\n", "\n", "## Authors: Leo Werneck & Zach Etienne\n", "\n", "**This module is currently under development**\n", "\n", "## This tutorial module generates the `Convert_to_HydroBase` ETK thorn files, compatible with the latest implementation of IllinoisGRMHD.\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](#convert_to_hydrobase__src): **`set_IllinoisGRMHD_metric_GRMHD_variables_based_on_HydroBase_and_ADMBase_variables.C`**\n", "1. [Step 3](#convert_to_hydrobase__param): **`param.ccl`**\n", "1. [Step 4](#convert_to_hydrobase__interface): **`interface.ccl`**\n", "1. [Step 5](#convert_to_hydrobase__schedule): **`schedule.ccl`**\n", "1. [Step 6](#convert_to_hydrobase__make): **`make.code.defn`**\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: Load up cmdline_helper and create the 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", "import cmdline_helper as cmd\n", "CtH_dir_path = os.path.join(\"..\",\"Convert_to_HydroBase\")\n", "cmd.mkdir(CtH_dir_path)\n", "CtH_src_dir_path = os.path.join(CtH_dir_path,\"src\")\n", "cmd.mkdir(CtH_src_dir_path)\n", "\n", "# Step 0b: Create the output file path\n", "outfile_path__Convert_to_HydroBase__source = os.path.join(CtH_src_dir_path,\"Convert_to_HydroBase.C\")\n", "outfile_path__Convert_to_HydroBase__make = os.path.join(CtH_src_dir_path,\"make.code.defn\")\n", "outfile_path__Convert_to_HydroBase__param = os.path.join(CtH_dir_path,\"param.ccl\")\n", "outfile_path__Convert_to_HydroBase__interface = os.path.join(CtH_dir_path,\"interface.ccl\")\n", "outfile_path__Convert_to_HydroBase__schedule = os.path.join(CtH_dir_path,\"schedule.ccl\")" ] }, { "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: `set_IllinoisGRMHD_metric_GRMHD_variables _based_on_HydroBase_and_ADMBase_variables.C` \\[Back to [top](#toc)\\]\n", "$$\\label{convert_to_hydrobase__src}$$" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Overwriting ../Convert_to_HydroBase/src/Convert_to_HydroBase.C\n" ] } ], "source": [ "%%writefile $outfile_path__Convert_to_HydroBase__source\n", "#include \"cctk.h\"\n", "#include \n", "#include \n", "#include \n", "#include \n", "#include \"cctk_Arguments.h\"\n", "#include \"cctk_Parameters.h\"\n", "\n", "#include \"IllinoisGRMHD_headers.h\"\n", "\n", "void Convert_to_HydroBase(CCTK_ARGUMENTS) {\n", "\n", " DECLARE_CCTK_ARGUMENTS;\n", " DECLARE_CCTK_PARAMETERS;\n", "\n", " // Generally, we only need the HydroBase variables for diagnostic purposes, so we run the below loop only at iterations in which diagnostics are run.\n", " if(Convert_to_HydroBase_every==0 || cctk_iteration%Convert_to_HydroBase_every!=0) return;\n", "\n", " /***************\n", " * PPEOS Patch *\n", " ***************\n", " * We will need to set up our EOS in\n", " * order to be able to compute eps below\n", " */\n", " eos_struct eos;\n", " initialize_EOS_struct_from_input(eos);\n", "\n", "#pragma omp parallel for\n", " for(int k=0;k u^i = \\Gamma ( U^i - \\beta^i / \\alpha )\n", " // which implies\n", " // v^i = u^i/u^0\n", " // = \\Gamma/u^0 ( U^i - \\beta^i / \\alpha ) <- \\Gamma = \\alpha u^0\n", " // = \\alpha ( U^i - \\beta^i / \\alpha )\n", " // = \\alpha U^i - \\beta^i\n", " CCTK_REAL lapseL=alp[index];\n", " CCTK_REAL lapseL_inv=1.0/lapseL;\n", " vel[CCTK_GFINDEX4D(cctkGH,i,j,k,0)] = (PRIMS[VX] + betax[index])*lapseL_inv;\n", " vel[CCTK_GFINDEX4D(cctkGH,i,j,k,1)] = (PRIMS[VY] + betay[index])*lapseL_inv;\n", " vel[CCTK_GFINDEX4D(cctkGH,i,j,k,2)] = (PRIMS[VZ] + betaz[index])*lapseL_inv;\n", "\n", " // \\alpha u^0 = 1/sqrt(1+γ^ij u_i u_j) = \\Gamma = w_lorentz\n", " // First compute u^0:\n", " // Derivation of first equation:\n", " // \\gamma_{ij} (v^i + \\beta^i)(v^j + \\beta^j)/(\\alpha)^2\n", " // = \\gamma_{ij} 1/(u^0)^2 ( \\gamma^{ik} u_k \\gamma^{jl} u_l /(\\alpha)^2 <- Using Eq. 53 of arXiv:astro-ph/0503420\n", " // = 1/(u^0 \\alpha)^2 u_j u_l \\gamma^{jl} <- Since \\gamma_{ij} \\gamma^{ik} = \\delta^k_j\n", " // = 1/(u^0 \\alpha)^2 ( (u^0 \\alpha)^2 - 1 ) <- Using Eq. 56 of arXiv:astro-ph/0503420\n", " // = 1 - 1/(u^0 \\alpha)^2 <= 1\n", " CCTK_REAL shiftxL = betax[index];\n", " CCTK_REAL shiftyL = betay[index];\n", " CCTK_REAL shiftzL = betaz[index];\n", "\n", " CCTK_REAL gxxL = gxx[index];\n", " CCTK_REAL gxyL = gxy[index];\n", " CCTK_REAL gxzL = gxz[index];\n", " CCTK_REAL gyyL = gyy[index];\n", " CCTK_REAL gyzL = gyz[index];\n", " CCTK_REAL gzzL = gzz[index];\n", "\n", " CCTK_REAL one_minus_one_over_alpha_u0_squared = (gxxL* SQR(PRIMS[VX] + shiftxL) +\n", " 2.0*gxyL*(PRIMS[VX] + shiftxL)*(PRIMS[VY] + shiftyL) +\n", " 2.0*gxzL*(PRIMS[VX] + shiftxL)*(PRIMS[VZ] + shiftzL) +\n", " gyyL* SQR(PRIMS[VY] + shiftyL) +\n", " 2.0*gyzL*(PRIMS[VY] + shiftyL)*(PRIMS[VZ] + shiftzL) +\n", " gzzL* SQR(PRIMS[VZ] + shiftzL) )*SQR(lapseL_inv);\n", " /*** Check for superluminal velocity ***/\n", " //FIXME: Instead of >1.0, should be one_minus_one_over_alpha_u0_squared > ONE_MINUS_ONE_OVER_GAMMA_SPEED_LIMIT_SQUARED, for consistency with conserv_to_prims routines\n", "\n", " if(one_minus_one_over_alpha_u0_squared > 1.0) {\n", " CCTK_VInfo(CCTK_THORNSTRING,\"Convert_to_HydroBase WARNING: Found superluminal velocity. This should have been caught by IllinoisGRMHD.\");\n", " }\n", "\n", " // A = 1.0-one_minus_one_over_alpha_u0_squared = 1-(1-1/(al u0)^2) = 1/(al u0)^2\n", " // 1/sqrt(A) = al u0\n", " CCTK_REAL alpha_u0 = 1.0/sqrt(1.0-one_minus_one_over_alpha_u0_squared);\n", " if(std::isnan(alpha_u0*lapseL_inv)) printf(\"BAD FOUND NAN ALPHAU0 CALC: %.15e %.15e %.15e\\n\",alpha_u0,lapseL_inv,one_minus_one_over_alpha_u0_squared);\n", "\n", " w_lorentz[index] = alpha_u0;\n", "\n", " Bvec[CCTK_GFINDEX4D(cctkGH,i,j,k,0)] = PRIMS[BX_CENTER];\n", " Bvec[CCTK_GFINDEX4D(cctkGH,i,j,k,1)] = PRIMS[BY_CENTER];\n", " Bvec[CCTK_GFINDEX4D(cctkGH,i,j,k,2)] = PRIMS[BZ_CENTER];\n", "\n", " }\n", "}\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "# Step 3: `param.ccl` \\[Back to [top](#toc)\\]\n", "$$\\label{convert_to_hydrobase__param}$$" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Overwriting ../Convert_to_HydroBase/param.ccl\n" ] } ], "source": [ "%%writefile $outfile_path__Convert_to_HydroBase__param\n", "# Parameter definitions for thorn convert_to_HydroBase\n", "# $Header:$\n", "\n", "#############################################################################\n", "### import HydroBase & ADMBase parameters\n", "\n", "shares: HydroBase\n", "USES CCTK_INT timelevels\n", "\n", "shares: ADMBase\n", "USES CCTK_INT lapse_timelevels\n", "USES CCTK_INT shift_timelevels\n", "USES CCTK_INT metric_timelevels\n", "\n", "shares: IllinoisGRMHD\n", "USES CCTK_INT neos\n", "USES CCTK_REAL Gamma_th\n", "USES CCTK_REAL K_ppoly_tab0\n", "USES CCTK_REAL rho_ppoly_tab_in[10]\n", "USES CCTK_REAL Gamma_ppoly_tab_in[10]\n", "#############################################################################\n", "\n", "private:\n", "INT Convert_to_HydroBase_every \"How often to convert IllinoisGRMHD primitive variables to HydroBase (Valencia formulation) primitive variables? Needed for some ET-based diagnostics. NOT needed for pure IllinoisGRMHD runs.\"\n", "{\n", " 0:* :: \"zero (disable) or positive (every N iterations)\"\n", "} 0\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "# Step 4: `interface.ccl` \\[Back to [top](#toc)\\]\n", "$$\\label{convert_to_hydrobase__interface}$$" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Overwriting ../Convert_to_HydroBase/interface.ccl\n" ] } ], "source": [ "%%writefile $outfile_path__Convert_to_HydroBase__interface\n", "# Interface definition for thorn Convert_to_HydroBase\n", "# $Header:$\n", "\n", "implements: Convert_to_HydroBase\n", "inherits: grid HydroBase ADMBase IllinoisGRMHD\n", "\n", "uses include header: IllinoisGRMHD_headers.h\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "# Step 5: `schedule.ccl` \\[Back to [top](#toc)\\]\n", "$$\\label{convert_to_hydrobase__schedule}$$" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Overwriting ../Convert_to_HydroBase/schedule.ccl\n" ] } ], "source": [ "%%writefile $outfile_path__Convert_to_HydroBase__schedule\n", "# Schedule definitions for thorn Convert_to_HydroBase\n", "# $Header:$\n", "\n", "SCHEDULE Convert_to_HydroBase AT CCTK_INITIAL AFTER SetTmunu\n", "{\n", " LANG: C\n", "} \"Convert IllinoisGRMHD-native variables to HydroBase\"\n", "\n", "SCHEDULE Convert_to_HydroBase AT CCTK_ANALYSIS BEFORE compute_bi_b2_Poyn_fluxET BEFORE particle_tracerET BEFORE VolumeIntegralGroup BEFORE convert_to_MHD_3velocity AFTER ML_BSSN_evolCalcGroup\n", "{\n", " OPTIONS: GLOBAL-EARLY,LOOP-LOCAL\n", " LANG: C\n", "} \"Convert IllinoisGRMHD-native variables to HydroBase\"\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "# Step 6: `make.code.defn` \\[Back to [top](#toc)\\]\n", "$$\\label{convert_to_hydrobase__make}$$" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Overwriting ../Convert_to_HydroBase/src/make.code.defn\n" ] } ], "source": [ "%%writefile $outfile_path__Convert_to_HydroBase__make\n", "# Main make.code.defn file for thorn Convert_to_HydroBase\n", "# $Header:$\n", "\n", "# Source files in this directory\n", "SRCS = Convert_to_HydroBase.C\n", "\n", "# Subdirectories containing source files\n", "SUBDIRS =\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "\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": 7, "metadata": {}, "outputs": [], "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": 8, "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.11.1" } }, "nbformat": 4, "nbformat_minor": 4 }