{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "\n", "\n", "\n", "# Tutorial-IllinoisGRMHD: postpostinitial__set_symmetries__copy_timelevels.C\n", "\n", "## Authors: Leo Werneck & Zach Etienne\n", "\n", "**This module is currently under development**\n", "\n", "## In this tutorial module we explain tasks that are performed just after the initial data has been set up. This module will likely be absorbed by another one once we finish documenting the code.\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](#postpostinitial__set_symmetries__copy_timelevels__c): **`postpostinitial__set_symmetries__copy_timelevels.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__postpostinitial__set_symmetries__copy_timelevels__C = os.path.join(IGM_src_dir_path,\"postpostinitial__set_symmetries__copy_timelevels.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: `postpostinitial__set_symmetries__copy_timelevels.C` \\[Back to [top](#toc)\\]\n", "$$\\label{postpostinitial__set_symmetries__copy_timelevels__c}$$" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Writing ../src/postpostinitial__set_symmetries__copy_timelevels.C\n" ] } ], "source": [ "%%writefile $outfile_path__postpostinitial__set_symmetries__copy_timelevels__C\n", "//-------------------------------------------------\n", "// Stuff to run right after initial data is set up\n", "//-------------------------------------------------\n", "\n", "#include \"cctk.h\"\n", "#include \n", "#include \n", "#include \"cctk_Arguments.h\"\n", "#include \"cctk_Functions.h\"\n", "#include \"cctk_Parameters.h\"\n", "#include \"Symmetry.h\"\n", "#include \"IllinoisGRMHD_headers.h\"\n", "#include \"IllinoisGRMHD_EoS_lowlevel_functs.C\"\n", "\n", "extern \"C\" void IllinoisGRMHD_PostPostInitial_Set_Symmetries__Copy_Timelevels(CCTK_ARGUMENTS) {\n", " DECLARE_CCTK_ARGUMENTS;\n", " DECLARE_CCTK_PARAMETERS;\n", "\n", " /**********************************\n", " * Piecewise Polytropic EOS Patch *\n", " * Printing the EOS table *\n", " **********************************/\n", " /*\n", " * The short piece of code below takes care\n", " * of initializing the EOS parameters.\n", " * Please refer to the \"inlined_functions.C\"\n", " * source file for the documentation on the\n", " * function.\n", " */\n", " eos_struct eos;\n", " initialize_EOS_struct_from_input(eos);\n", "\n", " /* For diagnostic and user convenience purposes, we print\n", " * out the EOS parameters (rho_ppoly_tab, K_ppoly_tab,\n", " * and Gamma_ppoly_tab) at t=0.\n", " */\n", " if(cctk_iteration==0 && (int)GetRefinementLevel(cctkGH)==0) { print_EOS_table(eos); }\n", "\n", " if(Gamma_th<0)\n", " 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\");\n", "\n", " //For emfields, we assume that you've set Bx, By, Bz (the UN-tilded B^i's)\n", " // or Ax, Ay, Az (if using constrained transport scheme of Del Zanna)\n", "\n", " if(CCTK_EQUALS(Symmetry,\"equatorial\")) {\n", " // SET SYMMETRY GHOSTZONES ON ALL CONSERVATIVE AND PRIMIIVE VARIABLES!\n", " int ierr;\n", " ierr=CartSymGN(cctkGH,\"IllinoisGRMHD::grmhd_conservatives\"); if(ierr!=0) CCTK_VError(VERR_DEF_PARAMS,\"Microsoft error code #1874109358120048. Grep it in the source code\");\n", " 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\");\n", "\n", " // Finish up by setting symmetry ghostzones on Bx, By, Bz, and their staggered variants.\n", " CCTK_REAL gridfunc_syms_Bx[3] = {-1, 1,-Sym_Bz};\n", " IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Bx , gridfunc_syms_Bx,0,0,0);\n", " IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Bx_stagger, gridfunc_syms_Bx,1,0,0);\n", " CCTK_REAL gridfunc_syms_By[3] = { 1,-1,-Sym_Bz};\n", " IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, By , gridfunc_syms_Bx,0,0,0);\n", " IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, By_stagger, gridfunc_syms_By,0,1,0);\n", " CCTK_REAL gridfunc_syms_Bz[3] = { 1, 1, Sym_Bz};\n", " IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Bz , gridfunc_syms_Bz,0,0,0);\n", " IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Bz_stagger, gridfunc_syms_Bz,0,0,1);\n", "\n", " CCTK_REAL gridfunc_syms_psi6phi[3] = { 1, 1, 1};\n", " IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z,psi6phi , gridfunc_syms_psi6phi,1,1,1);\n", " CCTK_REAL gridfunc_syms_Ax[3] = {-1, 1, Sym_Bz};\n", " IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Ax , gridfunc_syms_Ax,0,1,1);\n", " CCTK_REAL gridfunc_syms_Ay[3] = { 1,-1, Sym_Bz};\n", " IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Ay , gridfunc_syms_Ay,1,0,1);\n", " CCTK_REAL gridfunc_syms_Az[3] = { 1, 1,-Sym_Bz};\n", " IllinoisGRMHD_set_symmetry_gzs_staggered(cctkGH,cctk_lsh,x,y,z, Az , gridfunc_syms_Az,1,1,0);\n", " }\n", "\n", "\n", " //------------------------------------------------------------------\n", " // FILL _p AND _p_p TIMELEVELS. Probably don't need to do this if\n", " // Carpet::init_fill_timelevels=yes and\n", " // MoL::initial_data_is_crap = yes\n", " // NOTE: We don't fill metric data here.\n", " // FIXME: Do we really need this?\n", "\n", "#pragma omp parallel for\n", " for(int k=0;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 postpostinitial__set_symmetries__copy_timelevels.C: FAILED!\n", "Diff:\n", "12a13\n", "> #include \"IllinoisGRMHD_EoS_lowlevel_functs.C\"\n", "18,19c19,40\n", "< if(gamma_th<0)\n", "< 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\");\n", "---\n", "> /**********************************\n", "> * Piecewise Polytropic EOS Patch *\n", "> * Printing the EOS table *\n", "> **********************************/\n", "> /*\n", "> * The short piece of code below takes care\n", "> * of initializing the EOS parameters.\n", "> * Please refer to the \"inlined_functions.C\"\n", "> * source file for the documentation on the\n", "> * function.\n", "> */\n", "> eos_struct eos;\n", "> initialize_EOS_struct_from_input(eos);\n", "> \n", "> /* For diagnostic and user convenience purposes, we print\n", "> * out the EOS parameters (rho_ppoly_tab, K_ppoly_tab,\n", "> * and Gamma_ppoly_tab) at t=0.\n", "> */\n", "> if(cctk_iteration==0 && (int)GetRefinementLevel(cctkGH)==0) { print_EOS_table(eos); }\n", "> \n", "> if(Gamma_th<0)\n", "> 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\");\n", "53c74\n", "< // FILL _p AND _p_p TIMELEVELS. Probably don't need to do this if \n", "---\n", "> // FILL _p AND _p_p TIMELEVELS. Probably don't need to do this if\n", "65c86\n", "< mhd_st_x_p[index] = mhd_st_x[index]; \n", "---\n", "> mhd_st_x_p[index] = mhd_st_x[index];\n", "70c91\n", "< Ax_p[index] = Ax[index]; \n", "---\n", "> Ax_p[index] = Ax[index];\n", "76c97\n", "< mhd_st_x_p_p[index] = mhd_st_x[index]; \n", "---\n", "> mhd_st_x_p_p[index] = mhd_st_x[index];\n", "81c102\n", "< Ax_p_p[index] = Ax[index]; \n", "---\n", "> Ax_p_p[index] = Ax[index];\n", "85a107\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/postpostinitial__set_symmetries__copy_timelevels.C\"\n", "original_IGM_file_name = \"postpostinitial__set_symmetries__copy_timelevels-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__postpostinitial__set_symmetries__copy_timelevels__C = !diff $original_IGM_file_path $outfile_path__postpostinitial__set_symmetries__copy_timelevels__C\n", "\n", "if Validation__postpostinitial__set_symmetries__copy_timelevels__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 postpostinitial__set_symmetries__copy_timelevels.C: PASSED!\")\n", "else:\n", " # If the validation fails, we keep the original IGM source code file\n", " print(\"Validation test for postpostinitial__set_symmetries__copy_timelevels.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__postpostinitial__set_symmetries__copy_timelevels__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__postpostinitial__set_symmetries__copy_timelevels.pdf](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)." ] }, { "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__postpostinitial__set_symmetries__copy_timelevels.ipynb\n", "#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__postpostinitial__set_symmetries__copy_timelevels.tex\n", "#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__postpostinitial__set_symmetries__copy_timelevels.tex\n", "#!pdflatex -interaction=batchmode Tutorial-IllinoisGRMHD__postpostinitial__set_symmetries__copy_timelevels.tex\n", "!rm -f Tut*.out Tut*.aux Tut*.log" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "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.8.3" } }, "nbformat": 4, "nbformat_minor": 2 }