Skip to content

KHARMA Specific Conventions

Ben Prather edited this page Sep 3, 2024 · 15 revisions

There are a few conventions used throughout the KHARMA source. Knowing them will make reading the code much easier. If you know how to edit other Parthenon-based codes, this is the quick summary to get you started on KHARMA.

As with most Parthenon-based codes, an overview of the algorithm can be gleaned from the Driver class, in this case HARMDriver. The primary package is "GRMHD", in namespace and folder of the same name, and the reconstruction and Riemann solver are in fluxes.hpp (see GetFlux section).

Host and Device functions:

Device-side (technically, host/device) functions are named in snake_case, except where doing so conflicts with Parthenon's naming scheme (so far, only in reconstruction implementations). They are defined almost exclusively in headers, except when they are particular to a single file, and are rarely if ever used host-side. However, they are still generally namespaced, under the package most relevant to their functionality.

Host-side-only functions in KHARMA are named with InitialCaps. They generally operate on:

  1. Block-wise fluid states of type MeshBlockData<Real> *rc or
  2. Mesh-wide packed states of type MeshData<Real> *rc

Most functions take MeshData as this is faster -- the exceptions are the problem initialization functions, which are run separately over each block.

Macros and Types

KHARMA has several macros and universal types used throughout the code. These are defined in the universal header files decs.hpp and types.hpp. A few of these mirror definitions in iharm3D, some of which were themselves inherited from previous implementations. Both universal headers are short and worth reading in entirety.

VarMap

One struct in particular deserves attention here. To avoid the overhead of dealing with strings on the device side, KHARMA defines a struct with all possible variables as int8 members. It then defines a constructor which takes the Parthenon PackIndexMap, finding each variable by name and assigning it to its struct member.

Using the VarMap struct allows for code that looks similar to iharm3D or any other array-based code. In place of macro definitions like RHO, UU there are map members m_p.RHO, m_p.UU. Eventually this struct will be replaced by a more flexible system which maps variables to types in order to assign indices on the fly, but the result should look and behave nearly identically.

Flag("Name") and EndFlag()

KHARMA also defines a function representing every time the code enters or exits a major operation (most host-side functions). This enables the tracing functionality mentioned in the pages on compiling & running KHARMA, as well as helping with profiling.

Generally, package callbacks do not need to call Flag for themselves, as the code-wide functions will do it for them (e.g. your UtoP function in a new package will be flagged as UtoP_Packagename automatically). However, if writing new functionality called in some other way (e.g., directly by the Driver), it's good practice to call Flag with the function name at the beginning of the function, and EndFlag before returning.

Good flagging also helps in profiling KHARMA! If you suspect a region is expensive, you can Flag it even if it is not a separate function, and it will show up in kernel timing stats when you run ./run.sh prof [ARGS].

Real vs GReal vs double

Currently, there is very little distinction between Real, GReal and straight double variables in KHARMA -- all map to double type when it is compiled. However, in the future, when running on consumer GPUs for example, it may become advantageous to define Real as float instead, in order to perform most fluid operations in single-precision. This would almost certainly be with the exception of the magnetic field, and certainly with the exception of the geometry variables and timestep -- hence using separate GReal and double declarations for these. The distinction is not tested now, and needn't be kept perfectly yet.

Variables and Flags

KHARMA splits the usual set of 8 primitive variables for GRMHD into four Parthenon fields: scalars "rho" and "u", and 3-vectors "uvec" and "B". Primitive versions are prefixed with "prims.", and conserved with "cons.". Choosing fields with the flag "HD" will select the first three, and "MHD" will select all four. Examples of packing any given subset of variables are available by looking around in the code.

Variables which must be exchanged during startup, but not at each step are given a label "StartupOnly". This flag is applied to all variables associated with cleaning up the magnetic field when resizing a simulation, when they are not needed after the initial cleanup operation.

Each variable declared within a package will also be flagged with the package name, e.g. "Flux" for everything within the Flux package. A few additional user flags are added in specific packages, and each is generally only used within its own package. Feel free to add more user flags.

Adding Fluid Variables

New variables which should be advected by KHARMA should be marked with the Metadata::Conserved flag. Primitive versions should be marked with the "Primitive" user flag. Generally, new packages should use the vectors of flags "prims_flags" and "cons_flags" stored in the "Driver" package. This is because, confusingly, KHARMA may exchange either the primitive variables or conserved variables at boundaries, depending on which driver is selected.

Note that when adding variables, you'll also need to add conversion functions UtoP and prim_to_flux; the former a host-side function registered by your package as a callback, the latter a device-side function added manually to Flux::prim_to_flux so it can be called in GetFlux on the reconstructed primitive values at zone faces.

GetFlux

GetFlux, in the file get_flux.hpp, performs reconstruction on the primitive variables, converts the values at faces to conserved variables, and applies the LLF or HLLE Riemann solver to calculate the fluxes through faces of all conserved variables. It does all this in just a few kernels, which account together for about 60% of KHARMA's runtime (when evolving explicit ideal MHD in the absence of other performance issues).

It is defined in a header, despite being a host function, because it is a template on the direction to reconstruct (1,2,3) and reconstruction scheme (e.g. WENO5, linear). This makes the reconstruction and direction options constexpr within the function, which avoids large kernels full of branching if statements -- instead, different kernels are generated for each direction and reconstruction combo.

Unless you want to go spelunking for cycles (which you are welcome to do!) the only interaction you should need to have with this file is adding a new prim_to_flux implementation in flux_functions.hpp as described above.

Other function templates

Several other host-side functions are defined in headers, usually those ending in *_impl.hpp. In each instance, the function is templated on a single enum, and calls through to a device-side function templated on the same enum. This switches out different device-side implementations (e.g., initializing different magnetic fields), without the slowness of a large if statement within the kernel, or the significant boilerplate of a large if statement full of copies of the same loop, each calling a different function inside.

Templating host-side functions is not the simplest or most flexible solution, but it works and it is fast. To add an implementation to one of these template monsters, see the checklist at the top of seed_B.hpp.

Coordinates

KHARMA replaces Parthenon's UniformCoordinates object with its own subclass, GRCoordinates. GRCoordinates objects define several extra functions, in addition to those present in UniformCartesian:

  1. A translator for any zone number (and location within, e.g. the corner or X1 face) into coordinates x1, x2, x3 in the native (regular) coordinate system
  2. A cache of the local metric at several locations in each zone, and connection coefficients at zone centers. These cached values are carried around as ParArray2D objects (2D since all supported metrics are uniform in X3).
  3. A coordinate system object, which can translate between the native coordinate system (for example, "Modified" Kerr-Schild Coordinates, with a transformation applied for code speed and accuracy) and the embedding coordinate system (Kerr-Schild Coordinates). The coordinate system object also keeps track of coordinate properties (e.g., spherical or not) and spacetime properties (e.g., event horizon).

There is a compile-time option, FAST_CARTESIAN, which prevents geometry caching and simply hard-codes values for Cartesian Minkowski spacetime. This is not often used (KHARMA can of course also handle Cartesian Minkowski spacetime when compiled normally), but it can provide a modest speedup for SRMHD simulations.