In this part of the tutorial, we will construct an Einstein Toolkit (ETK) thorn (module) that will set up initial data for Maxwell's equations. In the Tutorial-MaxwellCartesian, tutorial notebook we used NRPy+ to contruct the SymPy expressions for plane-wave initial data.
We will construct this thorn in two steps.
After importing the core modules, since we are writing an ETK thorn, we'll need to set "grid::GridFuncMemAccess"
to "ETK"
. SymPy expressions for plane wave initial data are written inside Maxwell/MaxwellCartesian_ID.py, and we simply import them for use here.
# Step 1: Call on NRPy+ to convert the SymPy expression
# for the Maxwell's equations initial data into a C-code kernel
# Step 1a: Import needed NRPy+ core modules:
import NRPy_param_funcs as par # NRPy+: Parameter interface
import indexedexp as ixp # NRPy+: Symbolic indexed expression (e.g., tensors, vectors, etc.) support
import grid as gri # NRPy+: Functions having to do with numerical grids
import loop as lp # NRPy+: Generate C code loops
import sympy as sp # SymPy: The Python computer algebra package upon which NRPy+ depends
from outputC import lhrh # NRPy+: Core C code output module
import finite_difference as fin # NRPy+: Finite difference C code generation module
# Step 1b: This is an Einstein Toolkit (ETK) thorn. Here we
# tell NRPy+ that gridfunction memory access will
# therefore be in the "ETK" style.
par.set_parval_from_str("grid::GridFuncMemAccess","ETK")
#Set the spatial dimension parameter to 3.
par.set_parval_from_str("grid::DIM", 3)
DIM = par.parval_from_str("grid::DIM")
import Maxwell.MaxwellCartesian_Evol as mwrhs # This declare the parameter System_to_use needed for MaxwellCartesian_ID
# Step 1c: Call the MaxwellCartesian_ID() function from within the
# Maxwell/MaxwellCartesian_ID.py.py module.
import Maxwell.MaxwellCartesian_ID as mwid
# Step 3: Set up the plane wave initial data. This sets uu_ID and vv_ID.
mwid.MaxwellCartesian_ID()
# Step 4: Register gridfunctions so they can be written to by NRPy.
# System I:
AID = ixp.register_gridfunctions_for_single_rank1("EVOL","AID")
EID = ixp.register_gridfunctions_for_single_rank1("EVOL","EID")
psiI = gri.register_gridfunctions("EVOL","psiI")
# Step 5: FIXME: ADD DESCRIPTION
for i in range(DIM):
AID[i] = mwid.AidD[i]
EID[i] = mwid.EidD[i]
psiI = mwid.psi_ID
# Step 4: Register gridfunctions so they can be written to by NRPy.
# System II:
AIID = ixp.register_gridfunctions_for_single_rank1("EVOL","AIID")
EIID = ixp.register_gridfunctions_for_single_rank1("EVOL","EIID")
psiII = gri.register_gridfunctions("EVOL","psiII")
Gamma = gri.register_gridfunctions("EVOL","Gamma")
# Step 5: FIXME: ADD DESCRIPTION
for i in range(DIM):
AIID[i] = mwid.AidD[i]
EIID[i] = mwid.EidD[i]
psiII = mwid.psi_ID
Gamma = mwid.Gamma_ID
# Step 6: Create the C code output kernel.
Maxwell_ID_to_print = [\
lhrh(lhs=gri.gfaccess("out_gfs","AID0"),rhs=AID[0]),\
lhrh(lhs=gri.gfaccess("out_gfs","AID1"),rhs=AID[1]),\
lhrh(lhs=gri.gfaccess("out_gfs","AID2"),rhs=AID[2]),\
lhrh(lhs=gri.gfaccess("out_gfs","EID0"),rhs=EID[0]),\
lhrh(lhs=gri.gfaccess("out_gfs","EID1"),rhs=EID[1]),\
lhrh(lhs=gri.gfaccess("out_gfs","EID2"),rhs=EID[2]),\
lhrh(lhs=gri.gfaccess("out_gfs","psiI"),rhs=psiI),\
lhrh(lhs=gri.gfaccess("out_gfs","AIID0"),rhs=AIID[0]),\
lhrh(lhs=gri.gfaccess("out_gfs","AIID1"),rhs=AIID[1]),\
lhrh(lhs=gri.gfaccess("out_gfs","AIID2"),rhs=AIID[2]),\
lhrh(lhs=gri.gfaccess("out_gfs","EIID0"),rhs=EIID[0]),\
lhrh(lhs=gri.gfaccess("out_gfs","EIID1"),rhs=EIID[1]),\
lhrh(lhs=gri.gfaccess("out_gfs","EIID2"),rhs=EIID[2]),\
lhrh(lhs=gri.gfaccess("out_gfs","psiII"),rhs=psiII),\
lhrh(lhs=gri.gfaccess("out_gfs","Gamma"),rhs=Gamma),]
Maxwell_ID_CcodeKernel = fin.FD_outputC("returnstring",Maxwell_ID_to_print)
Maxwell_ID_looped = lp.loop(["i2","i1","i0"],["0","0","0"],["cctk_lsh[2]","cctk_lsh[1]","cctk_lsh[0]"],\
["1","1","1"],["#pragma omp parallel for","",""],"",\
Maxwell_ID_CcodeKernel.replace("time","cctk_time"))
# Step 7: Create directories for the thorn if they don't exist.
!mkdir MaxwellID 2>/dev/null # 2>/dev/null: Don't throw an error if the directory already exists.
!mkdir MaxwellID/src 2>/dev/null # 2>/dev/null: Don't throw an error if the directory already exists.
# Step 8: Write the C code kernel to file.
with open("MaxwellID/src/Maxwell_ID.h", "w") as file:
file.write(str(Maxwell_ID_looped))
%%writefile MaxwellID/src/InitialData.c
#include <math.h>
#include <stdio.h>
#include "cctk.h"
#include "cctk_Parameters.h"
#include "cctk_Arguments.h"
void Maxwell_InitialData(CCTK_ARGUMENTS)
{
DECLARE_CCTK_ARGUMENTS
DECLARE_CCTK_PARAMETERS
const CCTK_REAL *xGF = x;
const CCTK_REAL *yGF = y;
const CCTK_REAL *zGF = z;
const CCTK_REAL *gammaDD00GF = gxx;
const CCTK_REAL *gammaDD01GF = gxy;
const CCTK_REAL *gammaDD02GF = gxz;
const CCTK_REAL *gammaDD11GF = gyy;
const CCTK_REAL *gammaDD12GF = gyz;
const CCTK_REAL *gammaDD22GF = gzz;
#include "Maxwell_ID.h"
}
Writing MaxwellID/src/InitialData.c
Writing a module ("thorn") within the Einstein Toolkit requires that three "ccl" files be constructed, all in the root directory of the thorn:
interface.ccl
: defines the gridfunction groups needed, and provides keywords denoting what this thorn provides and what it should inherit from other thorns. Specifically, this file governs the interaction between this thorn and others; more information can be found in the official Einstein Toolkit documentation.With "implements", we give our thorn its unique name. By "inheriting" other thorns, we tell the Toolkit that we will rely on variables that exist and are declared "public" within those functions.
%%writefile MaxwellID/interface.ccl
implements: MaxwellID
inherits: admbase MaxwellEvol grid
Writing MaxwellID/interface.ccl
param.ccl
: specifies free parameters within the thorn, enabling them to be set at runtime. It is required to provide allowed ranges and default values for each parameter. More information on this file's syntax can be found in the official Einstein Toolkit documentation.%%writefile MaxwellID/param.ccl
shares: grid
USES KEYWORD type
restricted:
CCTK_KEYWORD initial_data "Type of initial data"
{
"toroid" :: "Toroidal Dipole field"
} "toroid"
restricted:
CCTK_REAL amp "The amplitude of the initial wave packet"
{
0.0:* :: "Should be positive"
} 1.0
restricted:
CCTK_REAL lam "A parameter describing the size of the wave package"
{
0.0:* :: "Should be positive"
} 1.0
Writing MaxwellID/param.ccl
schedule.ccl
: allocates storage for gridfunctions, defines how the thorn's functions should be scheduled in a broader simulation, and specifies the regions of memory written to or read from gridfunctions. schedule.ccl
's official documentation may be found here.We specify here the standardized ETK "scheduling bins" in which we want each of our thorn's functions to run.
%%writefile MaxwellID/schedule.ccl
schedule Maxwell_InitialData at CCTK_INITIAL as Maxwell_InitialData
{
STORAGE: MaxwellEvol::system_I[3]
STORAGE: MaxwellEvol::system_II[3]
LANG: C
READS: admbase::gxx(Everywhere)
READS: admbase::gxy(Everywhere)
READS: admbase::gxz(Everywhere)
READS: admbase::gyy(Everywhere)
READS: admbase::gyz(Everywhere)
READS: admbase::gzz(Everywhere)
READS: grid::x(Everywhere)
READS: grid::y(Everywhere)
READS: grid::y(Everywhere)
WRITES: AI0GF(Everywhere)
WRITES: AI1GF(Everywhere)
WRITES: AI2GF(Everywhere)
WRITES: EI0GF(Everywhere)
WRITES: EI1GF(Everywhere)
WRITES: EI2GF(Everywhere)
WRITES: psiIGF(Everywhere)
WRITES: AII0GF(Everywhere)
WRITES: AII1GF(Everywhere)
WRITES: AII2GF(Everywhere)
WRITES: EII0GF(Everywhere)
WRITES: EII1GF(Everywhere)
WRITES: EII2GF(Everywhere)
WRITES: psiIIGF(Everywhere)
WRITES: GammaGF(Everywhere)
} "Initial data for Maxwell's equations"
Writing MaxwellID/schedule.ccl
%%writefile MaxwellID/src/make.code.defn
SRCS = InitialData.c
Writing MaxwellID/src/make.code.defn
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-ETK_thorn-IDMaxwellNRPy.pdf (Note that clicking on this link may not work; you may need to open the PDF file through another means.)
import cmdline_helper as cmd # NRPy+: Multi-platform Python command-line interface
cmd.output_Jupyter_notebook_to_LaTeXed_PDF("Tutorial-ETK_thorn-IDMaxwellNRPy")
Created Tutorial-ETK_thorn-IDMaxwellNRPy.tex, and compiled LaTeX file to PDF file Tutorial-ETK_thorn-IDMaxwellNRPy.pdf