The Driver

The Driver object/package is responsible for constructing the list of tasks that are to be performed each (sub-)step. This means that if there is an operation that is to be executed every step, it has to be called from the driver in some way or another; either it is explicitly added to a TaskCollection made by the driver, or it can be a ``callback'' registered by a package, e.g. pkg->BlockUtoP, which will be called whenever KHARMA needs the primitive variables as a function of the conserved variables. The `driver` class contains a public function `Execute` that is called from `main`. In KHARMA, we currently have a single `KHARMADriver` class, which implements three different algorithms - a default `KHARMA` step, an `ImEx` step, and a `Simple` step for testing. One type of step is chosen in any given simulation, based on the parameter `driver/type` provided at runtime. Colloquially, the different steps are still referred to as "drivers" -- this is a holdover from when they were each different `Driver` subclasses.

TODO package+object combo, weird

The KHARMADriver class inherits from Parthenon's MultiStageDriver class, with each step creating a different TaskCollection, corresponding to a single sub-step (i.e., the task list is called twice in a second-order scheme like vl2 or rk2). We now discuss the algorithms available in KHARMA:

KHARMA Step

This is the default mode, used for most ideal GRMHD simulations. Unlike MHD, GRMHD must keep two forms of the variables: the conserved variables, and a set of "primitive" variables representing more tangible things like the density and velocity (primitive variables are also smoother and ). To evolve the fluid, the code must:

1. Reconstruct the right- and left-going components at zone faces, given the primitive variables
2. Calculate the fluxes of conserved quantities through the faces 2a. Apply any fixes to fluxes (e.g., for the magnetic field)
3. Update conserved variables using their prior values the divergence of conserved fluxes 3a. Apply any source terms (e.g., the geometric term in GRMHD) (Exchange ghost zones)
4. Recover primitive variables 4a. Apply any stability limits (floors) 4b. Fix any errors in recovering the primitives, re-apply floors to fixed zones
5. Apply any source terms (KEL), or calculate outputs (jcon) which use the primitive variables

This is before any synchronization between different blocks, etc, etc.

ImEx Step

Much like the KHARMA step, the Imex step starts off with counting the packages that are loaded and allocates the required number of fluid states each sub-step. There are however three key differences,

• The Imex step needs the implicit package to be loaded. This is because there are certain flags that the driver needs to know during execution in order to update variables appropriately (implicitly vs explicitly).
• It also needs two additional fluid states - one for the implicit solver, and one for the linesearch operation performed in the solver.
• The Imex step handles primitives as fundamental variables. This means that it is the primitives that are synced during boundary exchanges. This is specified by declaring the primitives with the Metadata::FillGhost flag.

It is probably worth noting that while the driver can handle ImEx schemes, it can also perform a fully explicit update. The advantage of using the Imex step over the KHARMA step for a fully explicit update has to do with the fact that this driver treats primitives as fundamental variables. This is favorable in the case of problems where one needs to define the primitives at physical boundaries manually eg., the conducting atmosphere test.

This is followed by a large synchronous TaskRegion where we perform several operations over the entire mesh,

1. Calculating the fluxes through all faces and the local wavespeeds.
2. Performing Flux-CT to ensure divergence-free magnetic fields at zone corners.
3. Computing the flux divergence.
4. Evaluating all the explicit sources based on the packages (physics) invoked.
5. Updating all explicit variables and computing the primitives so that they are in lockstep with the conserved variables.
6. Copying an initial guess for the implicit solver to the md_solver fluid state ONLY for the variables that are tagged isImplicit, and performing an implicit update.
7. Following the implicit update, md_solver contains the updated primitives. These are copied into the final state of the sub-step.
8. If fluid variables are evolved explicitly, copying a guess for the UtoP solver.

The next TaskRegion is an asynchronous one that is performed over all meshblocks concurrently. It executes the PreFillDerived, FillDerived and PostFillDerived calls for all the packages by calling Parthenon's Update::FillDerived function. This call rather surreptitiously carries out a lot of vital operations, applies floors (and instability limits for EMHD problems), performs UtoP for fluid variables and syncs pflag to name a few.

The driver then performs a boundary sync. This needs to be a synchronous region.

This is followed by an asynchronous region that,

1. Fixes zones where UtoP failed (this is only done if we updated fluid variables explicitly) over the entire domain.
2. Applies physical boundary conditions.
3. Performs any operations that need to be done post-sub-step, for eg, heating electrons.
4. Has a global PtoU call to ensure conserved vars and primitives are in sync.
5. Computes time step.

The final TaskRegion is an additional boundary sync that is executed only if two_sync is enabled in the GRMHD package. This is a safeguard to ensure the ghost zone primitives are identical to their physical counterpart, a strategy similar to the one employed in iharm3d.

Simple Step

The simple step is