GiRaFFE: An Open-Source General Relativistic Force-Free Electrodynamics Code

Zachariah B. Etienne <zachetie *at* gmail *dot* com>
Mew-Bing Wan
Maria C. Babiuc
Sean T. McWilliams
Ashok Choudhary

November 19, 2018

Abstract

GiRaFFE solves the equations of General Relativistic Force-Free Electrodynamics (GRFFE) using a high-resolution shock capturing scheme. It is based on IllinoisGRMHD, which is a rewrite of the Illinois Numerical Relativity (ILNR) group’s GRMHD code.

GiRaFFE is particularly good at modeling black hole magnetospheres. Its conservative-to-primitive solver has also been modified to check the physicality of conservative variables prior to primitive inversion, and move them into the physical range if they become unphysical.

Currently GiRaFFE consists of

  1. the Piecewise Parabolic Method (PPM) for reconstruction,
  2. the Harten, Lax, van Leer (HLL/HLLE) approximate Riemann solver, and
  3. a modified HARM Conservative-to-Primitive solver.

GiRaFFEevolves the vector potential Aμ (on staggered grids) instead of the magnetic fields (Bi) directly, to guarantee that the magnetic fields will remain divergenceless even at AMR boundaries. On uniform resolution grids, this vector potential formulation produces results equivalent to those generated using the standard, staggered flux-CT scheme.

1 Paper Reference

The basic equations, algorithms, and code validation tests for GiRaFFE are described in the code announcement paper, which supplements the many code comments as an excellent resource for understanding this code. To find the code announcement paper, search the arXiv for “GiRaFFE: An Open-Source General Relativistic Force-Free Electrodynamics Code”.

The relevant sections of the paper regarding the GRFFE formalism we adopt as well as the basic algorithmic approach are pasted below.

2 Adopted GRFFE Formalism

All equations are written in geometrized units, such that G = c = 1, and Einstein notation is chosen for implied summation. Greek indices span all 4 spacetime dimensions, and Latin indices span only the 3 spatial dimensions.

The line element for our spacetime in standard 3+1 form is given by

ds2 = α2dt2 + γ ij(dxi + βidt)(dxj + βjdt), (1)

where αdt denotes the proper time interval between adjacent spatial hypersurfaces separated by coordinate times t = t0 and t = t0 + dt, βi the magnitude of the spatial coordinate shift between adjacent hypersurfaces, and γij is the three-metric within a given hypersurface at coordinate time t.

Our 3+1 GRFFE formalism is written in terms of electric and magnetic fields as measured by an observer co-moving and normal to the spatial hypersurface, with 4-velocity nμ = (1α,βiα). The stress-energy tensor of the electromagnetic field is defined as:

TEMμν = F σμFνσ 1 4gμνF σηFση. (2)

where Fμν is the electromagnetic (or Faraday) tensor:

Fμν = nμEν nνEμ 𝜖μνσηB σnη, (3)

In terms of the Faraday tensor, the electric Eν and magnetic Bν fields in this frame are given by

Eν = n μFμν (4) Bν = 1 2nμ𝜖νμσηF ση = nμFμν, (5)

where Fμν is the dual of the Faraday tensor, and 𝜖νμση = [νμση]|g| is the rank-4 Levi-Civita tensor, with [νμση] the regular permutation symbol. In ideal MHD, the electric field vanishes for observers moving parallel to the plasma with 4-velocity uν:

Fμνu ν = 4πE(u)μ = 0. (6)

I.e., when an observer moves in a direction parallel to the magnetic field lines, the electric field vanishes. This is simply a statement of Ohm’s law for the case of perfect conductivity. One implication of perfect conductivity is that “magnetic field lines remain attached to the fluid elements they connect”—the so-called “frozen-in condition of MHD”. In addition, so long as Fμν0, the ideal MHD condition implies that

FμνF μν = 2(B2 E2) > 0 B2 > E2,and (7) FμνF μν = 4EμBμ = 0. (8)

See, e.g., [??] for further discussion.

In addition to the ideal MHD condition E(u)μ = 0, force-free electrodynamics also assumes that the plasma dynamics are completely driven by the electromagnetic fields (as opposed to, e.g., hydrodynamic pressure gradients). This implies that the stress-energy of the plasma: Tμν = Tmatterμν + TEMμν, is completely dominated by the electromagnetic terms, which yields the conservation equation [??]:

νTμν νTEMμν = FμνJ ν = 0 (9) ρEi + n ν𝜖νijkJ j(3)B k = 0. (10)

where Ji(3) is the 3-current, and ρ the charge density. This formalism is valid in the tenuous plasma of the stellar magnetosphere, where the rest-mass density is vanishingly small and assumed to be zero. From these, the 4-current can be expressed Jν = ρnν + γjiJj = γμνJμ

The left-hand side of Eq. 10 is simply the general relativistic expression for the Lorentz force, indicating that indeed the Lorentz force is zero in GRFFE. Notice that if we assume we are not in electrovacuum (Jμ0), multiplying Eq. 10 by Bi yields the familiar BiEi = 0 constraint of force-free electrodynamics (see [??] for full derivation).

In summary, the force-free constraints can be written

BiEi = 0,and (11) B2 > E2. (12)

Under these constraints, the GRFFE evolution equations consist of the Cauchy momentum equation and the induction equation (see [???] for derivation):

  1. The Cauchy momentum equation follows directly from the spatial components of μTμν = μTEMμν = 0 (the time component yields the energy equation, which in GRFFE is redundant). We choose to write the momentum equation in conservative form and in terms of the densitized spatial Poynting flux one-form S̃i = γSi,
    tS̃i + j αγTEMj i = 1 2αγTEMμν igμν, (13)

    where Si can be derived from the expression of the Poynting one-form,

    Sμ = nνTEMμν. (14)
  2. The induction equation in the force-free limit is derived from the spatial components of μFμν = 0 (the time components yield the “no-monopoles” constraint), and can be written in terms of the densitized magnetic field Bĩ = γBi as
    tB̃i + j vjB̃i viB̃j = 0, (15)

    where vj = uju0. As detailed in Appendix A of[?], the force-free conditions do not uniquely constrain uμ, allowing for the freedom to choose from a one-parameter family. As in [??], we choose uμ to correspond to the minimum plasma 3-velocity that satisfies Fμνuν = 0. This choice of vj is often referred to as the drift velocity, which can be defined in terms of known variables as

    vi = 4παγijS̃ j γB2 βi. (16)

3 Numerical Algorithms

We briefly review the numerical algorithms employed in GiRaFFE to solve the equations of GRFFE as outlined in Sec. 2.

GiRaFFE fully supports Cartesian adaptive mesh refinement (AMR) grids via the Cactus/Carpet [?] infrastructure within the Einstein Toolkit [?].

As in IllinoisGRMHD, GiRaFFE guarantees that the magnetic fields remain divergenceless to roundoff error even on AMR grids by evolving the vector potential 𝒜μ = Φnμ + Aμ, where Aμ is spatial (Aμnμ = 0), instead of the magnetic fields directly.



Table 1: Storage location on grid of the magnetic field Bi and vector potential 𝒜μ. Note that P is the vector of primitive variables {vi}


Variable(s) storage location


Metric terms, P, S̃i(i,j,k)
Bx, B̃x (i + 1 2,j,k)
By, B̃y (i,j + 1 2,k)
Bz, B̃z (i,j,k + 1 2)
Ax (i,j + 1 2,k + 1 2)
Ay (i + 1 2,j,k + 1 2)
Az (i + 1 2,j + 1 2,k)
γΦ (i + 1 2,j + 1 2,k + 1 2)



The vector potential fields exist on a staggered grid (as defined in Table 1 such that our magnetic fields are evolved according to the flux constrained transport (FluxCT) algorithm of Refs. [??].

As in IllinoisGRMHD, this choice of staggering in GiRaFFE means that standard symmetry thorns within the Einstein Toolkit cannot be used. To add this, one option is to modify src/symmetry__set_gzs_staggered_gfs.C.

Our choice of vector potential requires that we solve the vector potential version of the induction equation

tAi = 𝜖ijkvjBk i(αΦ βjA j), (17)

where 𝜖ijk = [ijk]γ is the anti-symmetric Levi-Civita tensor and γ is the 3-metric determinant, which in a flat spacetime in Cartesian coordinates reduces to 1. Bk in Eq. 17 is computed from the vector potential via

Bi = 𝜖ijk jAk = [ijk] γ jAk. (18)

Φ is evolved via an additional electromagnetic gauge evolution equation, which was devised specifically to avoid the buildup of numerical errors due to zero-speed characteristic modes [?] on AMR grids. Our electromagnetic gauge is identical to the Lorenz gauge, but with an exponential damping term with damping constant ξ [?]:

t γΦ + j αγAj βj γΦ = ξα γΦ. (19)

Step 0: Initial data: In addition to 3+1 metric quantities in the Arnowitt-Deser-Misner (ADM) formalism [?], GiRaFFE requires that the “Valencia” 3-velocity v̄i and vector potential one-form Aμ be set initially. Regarding the former, the “Valencia” 3-velocity v̄i is related to the 3-velocity appearing in the induction equation vi via

vi = ui u0 = αv̄i βi. (20)

As for Aμ, for all cases in this paper, we set the evolution variable γΦ = 0 initially, and Ai is set based on our initial physical scenario.

After vi and Aμ are set, Bi is computed via Eq. 18, and then the evolution variable S̃i is given by

S̃i = γij(vj + βj)γB2 4πα . (21)

Step 1: Evaluation of evolution equations: In tandem with the high-resolution shock-capturing (HRSC) scheme within GiRaFFE, the Runge-Kutta fourth-order (RK4) scheme is chosen to march our evolution variables Ai and S̃i forward in time from their initial values, adopting precisely the same reconstruction and Riemann solver algorithms as in IllinoisGRMHD (see Ref. [?] for more details). In short, Ai and S̃i are evolved forward in time using the Piecewise Parabolic Method (PPM) [?] for reconstruction and a Harten-Lax-van Leer (HLL)-based algorithm [??] for approximately solving the Riemann problem. Meanwhile, spatial derivatives within γΦ’s evolution equation (Eq. 19) are evaluated via finite difference techniques (as in IllinoisGRMHD).

Step 2: Boundary conditions on Aμ: At the end of each RK4 substep, the evolved variables Ai and S̃i have been updated at all points except the outer boundaries. So next the outer boundary conditions on Ai and γΦ are applied. As no exact outer boundary conditions typically exist for systems of interest to GiRaFFE, we typically take advantage of AMR and push our outer boundary out of causal contact from the physical system of interest. However, to retain good numerical stability, we apply “reasonable” outer boundary conditions. Specifically, values of Ai and γΦ in the interior grid are linearly extrapolated to the outer boundary.

Step 3: Computing Bi: Bi is next computed from Ai via Eq. 18.

Step 4: Applying GRFFE constraints & computing vi: Truncation, roundoff, and under-sampling errors will at times push physical quantities into regions that violate the GRFFE constraints. To nudge the variables back into a physically realistic domain, we apply the same strategy as was devised in Ref. [?] to guarantee that the GRFFE constraints remain satisfied:

First, we adjust S̃i via

S̃i S̃i (S̃jB̃j)B̃i B̃2 (22)

to enforce BiSi = 0, which as shown by Ref. [?], is equivalent to the GRFFE constraint Eq. 11.

Next, we limit the Lorentz factor of the plasma, typically to be 2,000, by rescaling S̃i according to Eq. 92 in Ref. [?]. After S̃i is rescaled the 3-velocity vi is recomputed via Eq. 16.

Finally, errors within our numerical scheme dissipate sharp features, so when current sheets appear, they are quickly and unphysically dissipated. This is unfortunate because current sheets lie at the heart of many GRFFE phenomena. So to remedy the situation, we apply the basic strategy of McKinney [?] (that was also adopted by Paschalidis and Shapiro [?]) and set the velocity perpendicular to the current sheet v to zero. For example, if the current sheet exists on the z = 0 plane, then v = vz, which we set to zero via nivi = 0, where ni = γijnj is a unit normal one-form with nj = δjz. Specifically, in the case of a current sheet on the z = 0 plane, we set

vz = γxzvx + γ yzvy γzz (23)

at all gridpoints that lie within |z| 4Δz of the current sheet.

At present the code addresses numerical dissipation of current sheets only if they appear on the equatorial plane. For cases in which current sheets appear off of the equatorial plane or spontaneously, [?] suggest the development of algorithms akin to reconnection-capturing strategies [?]. We intend to explore such approaches in future work.

Step 5: Boundary conditions on vi: vi is set to zero at a given face of our outermost AMR grid cube unless the velocity is outgoing. Otherwise the value for the velocity is simply copied from the interior grid to the nearest neighbor on a face-by-face basis.

After boundary conditions on vi are updated, all data needed for the next RK4 substep have been generated, so we return to Step 1.

4 Parameters




damp_lorenz
Scope: private  REAL



Description: Damping factor for the generalized Lorenz gauge. Has units of 1/length = 1/M. Typically set this parameter to 1.5/(maximum Delta t on AMR grids).



Range  Default: 0.0
*:*
any real






current_sheet_null_v
Scope: restricted  BOOLEAN



Description: Shall we null the velocity normal to the current sheet?



 Default: no






em_bc
Scope: restricted  KEYWORD



Description: EM field boundary condition



Range  Default: copy
copy
Copy data from nearest boundary point
frozen
Frozen boundaries






gamma_speed_limit
Scope: restricted  REAL



Description: Maximum relativistic gamma factor. Note the default is much higher than IllinoisGRMHD. (GRFFE can handle higher Lorentz factors)



Range  Default: 2000.0
1:*
Positive > 1, though you’ll likely have troubles far above 2000.






min_radius_inside_of_which_conserv_to_prims_ffe_and_ffe_evolution_is_disabled
Scope: restricted  REAL



Description: As parameter suggests, this is the minimum radius inside of which the conservatives-to-primitives solver is disabled. In the Aligned Rotator test, this should be set equal to R_NS_aligned_rotator.



Range  Default: -1.
-1.
”disable the conservative-to-prim itive solver modification”
(0:*
any positive value






sym_bz
Scope: restricted  REAL



Description: In-progress equatorial symmetry support: Symmetry parameter across z axis for magnetic fields = +/- 1



Range  Default: 1.0
-1.0:1.0
Set to +1 or -1.






symmetry
Scope: restricted  KEYWORD



Description: Currently only no symmetry supported, though work has begun in adding equatorial-symmetry support. FIXME: Extend ET symmetry interface to support symmetries on staggered gridfunctions



Range  Default: none
none
no symmetry, full 3d domain






velocity_bc
Scope: restricted  KEYWORD



Description: Chosen fluid velocity boundary condition



Range  Default: outflow
outflow
Outflow boundary conditions
copy
Copy data from nearest boundary point
frozen
Frozen boundaries






lapse_timelevels
Scope: shared from ADMBASE INT






metric_timelevels
Scope: shared from ADMBASE INT






shift_timelevels
Scope: shared from ADMBASE INT



5 Interfaces

General

Implements:

giraffe

Inherits:

admbase

boundary

spacemask

tmunubase

hydrobase

grid

Grid Variables

5.0.1 PRIVATE GROUPS




  Group Names     Variable Names    Details   




diagnostic_gfs   compact0
failure_checker   descriptionGridfunction to track conservative-to-primitives solver fixes. Beware that this gridfunction is overwritten at each RK substep.
  dimensions3
  distributionDEFAULT
  group typeGF
  tagsprolongation=”none” Checkpoint=”no”
  timelevels1
 variable typeREAL




grmhd_primitives_reconstructed_temps   compact0
temporary   descriptionTemporary variables used for primitives reconstruction
vxr   dimensions3
vyr   distributionDEFAULT
vzr   group typeGF
Bxr   tagsprolongation=”none” Checkpoint=”no”
Byr   timelevels1
Bzr  variable typeREAL




grmhd_conservatives_rhs   compact0
st_x_rhs   descriptionStorage for the right-hand side of the partial_t tilde{S}_i equations. Needed for MoL timestepping.
st_y_rhs   dimensions3
st_z_rhs   distributionDEFAULT
  group typeGF
  tagsprolongation=”none” Checkpoint=”no”
  timelevels1
 variable typeREAL




em_ax_rhs   compact0
Ax_rhs   descriptionStorage for the right-hand side of the partial_t A_x equation. Needed for MoL timestepping.
  dimensions3
  distributionDEFAULT
  group typeGF
  tagsprolongation=”none” Checkpoint=”no”
  timelevels1
 variable typeREAL




em_ay_rhs   compact0
Ay_rhs   descriptionStorage for the right-hand side of the partial_t A_y equation. Needed for MoL timestepping.
  dimensions3
  distributionDEFAULT
  group typeGF
  tagsprolongation=”none” Checkpoint=”no”
  timelevels1
 variable typeREAL




em_az_rhs   compact0
Az_rhs   descriptionStorage for the right-hand side of the partial_t A_z equation. Needed for MoL timestepping.
  dimensions3
  distributionDEFAULT
  group typeGF
  tagsprolongation=”none” Checkpoint=”no”
  timelevels1
 variable typeREAL








  Group Names     Variable Names    Details   




em_psi6phi_rhs   compact0
psi6phi_rhs   descriptionStorage for the right-hand side of the partial_t (psi6̂ Phi) equation. Needed for MoL timestepping.
  dimensions3
  distributionDEFAULT
  group typeGF
  tagsprolongation=”none” Checkpoint=”no”
  timelevels1
 variable typeREAL




grmhd_cmin_cmax_temps   compact0
cmin_x   descriptionStore min and max characteristic speeds in all three directions.
cmax_x   dimensions3
cmin_y   distributionDEFAULT
cmax_y   group typeGF
cmin_z   tagsprolongation=”none” Checkpoint=”no”
cmax_z   timelevels1
 variable typeREAL




grmhd_flux_temps   compact0
st_x_flux   descriptionTemporary variables for storing the flux terms of tilde{S}_i.
st_y_flux   dimensions3
st_z_flux   distributionDEFAULT
  group typeGF
  tagsprolongation=”none” Checkpoint=”no”
  timelevels1
 variable typeREAL




tupmunu   compact0
TUPtt   descriptionT{̂mu nu}
    descriptionstored to avoid expensive recomputation
TUPtx   dimensions3
TUPty   distributionDEFAULT
TUPtz   group typeGF
TUPxx   tagsprolongation=”none” Checkpoint=”no”
TUPxy   timelevels1
TUPxz  variable typeREAL




5.0.2 PUBLIC GROUPS




  Group Names     Variable Names    Details   




grmhd_conservatives   compact0
mhd_st_x   descriptionEvolved mhd variables
mhd_st_y   dimensions3
mhd_st_z   distributionDEFAULT
  group typeGF
  timelevels3
 variable typeREAL




em_ax   compact0
Ax   descriptionx-component of the vector potential
    descriptionevolved when constrained_transport_scheme==3
  dimensions3
  distributionDEFAULT
  group typeGF
  tagsProlongation=”STAGGER011”
  timelevels3
 variable typeREAL




em_ay   compact0
Ay   descriptiony-component of the vector potential
    descriptionevolved when constrained_transport_scheme==3
  dimensions3
  distributionDEFAULT
  group typeGF
  tagsProlongation=”STAGGER101”
  timelevels3
 variable typeREAL




em_az   compact0
Az   descriptionz-component of the vector potential
    descriptionevolved when constrained_transport_scheme==3
  dimensions3
  distributionDEFAULT
  group typeGF
  tagsProlongation=”STAGGER110”
  timelevels3
 variable typeREAL




em_psi6phi   compact0
psi6phi   descriptionsqrt{gamma} Phi
    descriptionwhere Phi is the em scalar potential
  dimensions3
  distributionDEFAULT
  group typeGF
  tagsProlongation=”STAGGER111”
  timelevels3
 variable typeREAL




grmhd_primitives_allbutbi   compact0
vx   descriptionComponents of three velocity
    descriptiondefined in terms of 4-velocity as: vî = uî/u0̂. Note that this definition differs from the Valencia formalism.
vy   dimensions3
vz   distributionDEFAULT
  group typeGF
  tagsInterpNumTimelevels=1 prolongation=”none”
  timelevels1
 variable typeREAL








  Group Names    Variable Names    Details   




grmhd_primitives_bi   compact0
Bx   descriptionB-field components defined both at vertices (Bx
    descriptionBy
Bx   descriptionBz) as well as at staggered points [Bx_stagger at (i+1/2
Bx   descriptionj
Bx   descriptionk)
Bx   descriptionBy_stagger at (i
Bx   descriptionj+1/2
Bx   descriptionk)
Bx   descriptionBz_stagger at (i
Bx   descriptionj
Bx   descriptionk+1/2)].
By   dimensions3
Bz   distributionDEFAULT
Bx_stagger   group typeGF
By_stagger   tagsInterpNumTimelevels=1 prolongation=”none”
Bz_stagger   timelevels1
 variable typeREAL




bssn_quantities   compact0
gtxx   descriptionBSSN quantities
    descriptioncomputed from ADM quantities
gtxy   dimensions3
gtxz   distributionDEFAULT
gtyy   group typeGF
gtyz   tagsprolongation=”none” Checkpoint=”no”
gtzz   timelevels1
gtupxx  variable typeREAL




Adds header:

GiRaFFE_headers.h

Uses header:

Symmetry.h

6 Schedule

This section lists all the variables which are assigned storage by thorn WVUThorns/GiRaFFE. Storage can either last for the duration of the run (Always means that if this thorn is activated storage will be assigned, Conditional means that if this thorn is activated storage will be assigned for the duration of the run if some condition is met), or can be turned on for the duration of a schedule function.

Storage

 

Always:  
ADMBase::metric[metric_timelevels] ADMBase::curv[metric_timelevels] ADMBase::lapse[lapse_timelevels] ADMBase::shift[shift_timelevels]  
GiRaFFE::BSSN_quantities  
grmhd_conservatives[3] em_Ax[3] em_Ay[3] em_Az[3] em_psi6phi[3]  
grmhd_primitives_allbutBi grmhd_primitives_Bi grmhd_primitives_reconstructed_temps grmhd_conservatives_rhs em_Ax_rhs em_Ay_rhs em_Az_rhs em_psi6phi_rhs grmhd_cmin_cmax_temps grmhd_flux_temps TUPmunu diagnostic_gfs 
   

Scheduled Functions

MoL_Register

  giraffe_registervars

  register evolved, rhs variables in giraffe for mol

 

 After: bssn_registervars
   lapse_registervars
 Language:c
 Options: meta
 Type: function

CCTK_BASEGRID

  giraffe_initsymbound

  schedule symmetries

 

 After: lapse_initsymbound
 Language:c
 Type: function

HydroBase_Boundaries

  giraffe_compute_b_and_bstagger_from_a

  compute b and b_stagger from a

 

 After: giraffe_outer_boundaries_on_a_mu
 Language:c
 Sync: grmhd_primitives_bi
 Type: function

AddToTmunu

  giraffe_conserv_to_prims_ffe

  applies the ffe condition b2̂>e2̂ and recomputes the velocities

 

 After: giraffe_compute_b_and_bstagger_from_a
 Language:c
 Type: function

AddToTmunu

  giraffe_outer_boundaries_on_vx_vy_vz

  apply outflow-only, flat bcs on {vx,vy,vz}. outflow only == velocities pointed inward from outer boundary are set to zero.

 

 After: conserv_to_prims
 Language:c
 Sync: grmhd_primitives_allbutbi
   grmhd_conservatives
 Type: function

CCTK_POSTPOSTINITIAL

  giraffe_postpostinitial

  hydrobase_con2prim in cctk_postpostinitial sets conserv to prim then outer boundaries (obs, which are technically disabled). the post ob syncs actually reprolongate the conservative variables, making cons and prims inconsistent. so here we redo the con2prim, avoiding the sync afterward, then copy the result to other timelevels

 

 After: hydrobase_con2prim
 Before:mol_poststep
 Type: group

GiRaFFE_PostPostInitial

  giraffe_initsymbound

  schedule symmetries – actually just a placeholder function to ensure prolongations / processor syncs are done before outer boundaries are updated.

 

 Before: compute_b
 Language:c
 Sync: grmhd_conservatives
   em_ax
   em_ay
   em_az
   em_psi6phi
 Type: function

GiRaFFE_PostPostInitial

  giraffe_compute_b_and_bstagger_from_a

  compute b and b_stagger from a sync: grmhd_primitives_bi

 

 After: postid
   empostid
   lapsepostid
 Language:c
 Sync: grmhd_primitives_bi
 Type: function

GiRaFFE_PostPostInitial

  giraffe_conserv_to_prims_ffe

  applies the ffe condition b2̂>e2̂ and recomputes the velocities

 

 After: compute_b
 Language:c
 Type: function

GiRaFFE_PostPostInitial

  giraffe_postpostinitial_set_symmetries__copy_timelevels

  compute post-initialdata quantities

 

 After: compute_b
   p2c
 Language:c
 Type: function

MoL_CalcRHS

  giraffe_driver_evaluate_ffe_rhs

  evaluate rhss of grffe equations

 

 After: bssn_rhs
    shift_rhs
 Language:c
 Type: function

HydroBase_Boundaries

  giraffe_initsymbound

  schedule symmetries – actually just a placeholder function to ensure prolongations / processor syncs are done before outer boundaries are updated.

 

 After: hydrobase_applybcs
 Before: giraffe_outer_boundaries_on_a_mu
 Language:c
 Sync: grmhd_conservatives
   em_ax
   em_ay
   em_az
   em_psi6phi
 Type: function

HydroBase_Boundaries

  giraffe_outer_boundaries_on_a_mu

  apply linear extrapolation bcs on a_mu, so that bcs are flat on bî.

 

 After: giraffe_initsymbound
 Before: giraffe_compute_b_and_bstagger_from_a
 Language:c
 Sync: em_ax
   em_ay
   em_az
   em_psi6phi
 Type: function

Aliased Functions

 

Alias Name:         Function Name:
GiRaFFE_PostPostInitial_Set_Symmetries__Copy_Timelevelsmhdpostid
GiRaFFE_compute_B_and_Bstagger_from_A compute_b
GiRaFFE_driver_evaluate_FFE_rhs GiRaFFE_RHS_eval