Dear All,
Scott Noble and I wrote up some notes on MHD development in GRHydro. This is the basis of the discussion currently going on on the ET phone call.
- Christian
--------------------------
We propose the following strategy for implementing GRMHD in EinsteinEvolve/GRHydro:
* Follow the formalism laid out in Anton et al. 2006, i.e., the Valencia formulation of GRMHD.
Key components to implement / technical challenges will be: (a) Con2Prim (b) GRMHD Riemann problem (c) Constrained transport (d) Recovery procedure for handling unphysical/intractable updates (e) GRMHD with AMR, refluxing and all.
Logical Steps:
(I) Initially, focus on unigrid; no AMR. Do AMR in a second step. Get basic formalism working with a GRMHD Con2Prim scheme, TVD reconstruction, HLL Riemann solver, constrained transport.
(II) Extend to AMR.
Concrete steps in (I):
(1) Define comoving magnetic field b^a in GRHydro as one 3D scalar and one 3D vector GF (b0 and b^i).
(2) Implement (in its own routine) a Con2Prim routine for GRMHD based on routine(s) used in Scott Noble's HARM3D code. This will use the old EOS interface and will be tested with polytropes (no special routine for polytropes, though) and gamma-law EOS. Perform Con2Prim tests.
(3) Implement Prim2Con, reconstruction (add calls), GRMHD Riemann problem, and update of the magnetic field variables (via MoL). Implement divergence check; run code to see how div B grows in test problems. Initial data routines to setup test problems (e.g. those in Gammie++ 2003) will be written.
(4) Implement constrained transport on unigrid. Initial strategy will be to use the FluxCT algorithm (with parabolic EMF reconstruction) as in HARM3D and the lower-order version as in the original HARM code. Future plans may involve adding alternate divergence constraint methods.
(5) Once all the basic routines are tested and working together, we will add a more sophisticated recovery procedure for dealing with instances where either the floor is penetrated, Con2Prim fails to converge to a physical solution, or a sudden and significant loss of energy/momentum/mass conservation occurs due to runaway recovery operations. From prior experience, this procedure will likely require an evolutionary style of development which will have to be tailored to the problem at hand.
(I.1)-(I.5) should be done without touching GRHydro's GRHD routines/capabilities.
Concrete steps in (II): (Need Erik's input as it they will require expertise on Carpet's data structures/interfaces that are not readily accessible to Cactus-at-large methods).
(1) Support for "refluxing", aka child-to-parent flux communication along refinement boundaries. This procedure requires interface to Carpet's data structures. refs: Berger & Colella 1989.
(2) Injection/restriction routine for communicating children EMF data to parent for parent's constrained transport method. This will likely require breaking down the FluxCT procedure into two parts: one that interpolates the face-centered fluxes to calculate the EMFs at the cell vertices, and one that uses the EMFs to calculate the constraint-preserving numerical fluxes. The restriction step will come in between these two steps, as it will modify the parent's EMFs before they are used to calculate the parent's final numerical fluxes. refs: Toth 2000, Balsara & Spicer 1999
(3) Adding a mask function that indicates whether a particular cell is covered by finer cells. This will be useful for adjusting the rigor with which we recover the variables (see I.5). For example, if a cell's updated conserved variables fail to yield a physical set of primitive variables and that cell has "children" then we can immediately skip (I.5) and wait for the restriction operation.