gasDynamicMaxwellDednerEqn

Defines the equations of inviscid fluid dynamics coupled to Maxwell’s equations in source term form with divergence cleaning:

\[\notag \begin{align} \frac{\partial \rho}{\partial t} + \nabla\cdot\left[ \rho\,\mathbf{u} \right] = 0 \\ \frac{\partial \left(\rho \mathbf{u} \right)}{\partial t} + \nabla\cdot\left[ \rho\,\mathbf{u}\,\mathbf{u}^{T} +\mathbb{I} P \right] = \sum_{\mathrm{species}} \left( q^{\mathrm{species}} \mathbf{E} + \mathbf{J}^{\mathrm{species}} \times \mathbf{B} \right) \\ \frac{\partial E}{\partial t} + \nabla\cdot\left[ \left(E + P \right) \mathbf{u} \right] = \sum_{\mathrm{species}} \mathbf{J}^{\mathrm{species}} \cdot \mathbf{E} \\ \frac{\partial \mathbf{B}^{\mathrm{plasma}}}{\partial t} + \nabla\times\mathbf{E} + \nabla \Psi= 0 \\ \frac{\partial \mathbf{E}^{\mathrm{plasma}}}{\partial t} - c^2 \nabla\times\mathbf{B} + \nabla \Phi= - \epsilon^{-1}_{0} \sum_{\mathrm{species}} \mathbf{J}^{\mathrm{species}} \\ \frac{\partial \Psi}{\partial t} + \nabla\cdot\left[ c^{2}_{\mathrm{fast}} \mathbf{B}^{\mathrm{plasma}} \right] = 0 \\ \frac{\partial \Phi}{\partial t} + \nabla\cdot\left[ c^{2}_{\mathrm{fast}} \mathbf{E}^{\mathrm{plasma}} \right] = \sum_{\mathrm{species}} q^{\mathrm{species}} \\ \end{align}\]

Here, \(q^{\mathrm{species}}\) is the species charge density, \(\mathbf{J}^{\mathrm{species}}\) is the species current density, \(\mathbb{I}\) is the identity matrix, \(P = \rho\epsilon(\gamma-1)\) is the pressure of an ideal gas, \(\epsilon\) is the specific internal energy and \(\gamma\) is the adiabatic index (ratio of specific heats). The quantity \(c_{\mathrm{fast}}\) corresponds to the fastest wave speed over the entire simulation domain; divergence errors are advected out of the domain with this speed.

In order to integrate these equations, USim casts them into flux-conservative form using the following standard identities (note that the use of these identities does not require an assumption of quasi-neutrality):

\[\notag \begin{align} \sum_{\mathrm{species}} \left( q^{\mathrm{species}} \mathbf{E} + \mathbf{J}^{\mathrm{species}} \times \mathbf{B} \right) = - \frac{\partial c^{-2} \mathbf{S}^\mathrm{EM}}{\partial t} + \nabla\cdot\mathcal{T}^{\mathrm{EM}} \\ \sum_{\mathrm{species}} \mathbf{J}^{\mathrm{species}} \cdot \mathbf{E} = -\frac{\partial E^{\mathrm{EM}}}{\partial t} - \nabla\cdot \mathbf{S}^{\mathrm{EM}} \end{align}\]

Here, \(\mathcal{T}^{\mathrm{EM}}\) is the electromagnetic stress tensor and \(\mathbf{S}^{\mathrm{EM}}\) is the electromagnetic energy (Poynting) flux vector, which are defined as:

\[\notag \begin{align} \mathcal{T}^{\mathrm{EM}} = \frac{1}{\mu_0} \left( \frac{\mathbf{E} \mathbf{E}^T}{c^2} + \mathbf{B} \mathbf{B}^T \right) + \mathbb{I} E_{\mathrm{EM}} = \frac{\mathbf{e} \mathbf{e}^T}{c^2} + \mathbf{b} \mathbf{b}^T + \mathbb{I} E_{\mathrm{EM}} \\ \mathbf{S}^{\mathrm{EM}} = \mu^{-1}_{0}\mathbf{E} \times \mathbf{B} = \mathbf{e} \times \mathbf{b} \\ E^{\mathrm{EM}} = \frac{1}{2\mu_0} \left( \frac{\left| \mathbf{E} \right|^2}{c^{2}} + \left| \mathbf{B} \right|^2 \right) = \frac{1}{2} \left( \frac{\left| \mathbf{e} \right|^2}{c^{2}} + \left| \mathbf{b} \right|^2 \right) \end{align}\]

Here, \(E^{\mathrm{EM}}\) is the electromagnetic energy density and the electromagnetic fields are defined as:

\[\notag \begin{align} \mathbf{b} = \mathbf{b}^{\mathrm{plasma}}+\mathbf{b}^{\mathrm{external}} = \mu^{-1/2}_{0} \left(\mathbf{B}^{\mathrm{plasma}}+\mathbf{B}^{\mathrm{external}} \right) = \mu^{-1/2}_{0} \mathbf{B} \\ \mathbf{e} = \mathbf{e}^{\mathrm{plasma}} + \mathbf{e}^{\mathrm{external}} = \mu^{-1/2}_{0} \left(\mathbf{E}^{\mathrm{plasma}} + \mathbf{E}^{\mathrm{external}} \right) = \mu^{-1/2}_{0} \mathbf{E} \end{align}\]

Here, \(\mathbf{b}^{\mathrm{plasma}}\) is the magnetic field induced in the plasma, \(\mathbf{e}^{\mathrm{plasma}}\) is the electric field associated with net charge in the plasma, while \(\mathbf{e}^{\mathrm{external}}\) and \(\mathbf{b}^{\mathrm{external}}\) are electromagnetic fields computed “externally” to Maxwell’s equations inside the plasma. With these identitifications, the fluid part of the gasDynamicMaxwellDednerEqn takes the form:

\[\notag \begin{align} \frac{\partial \rho}{\partial t} + \nabla\cdot\left[ \rho\,\mathbf{u} \right] = 0 \\ \frac{\partial \left(\rho \mathbf{u} + c^{-2} \mathbf{S}^\mathrm{EM} \right)}{\partial t} + \nabla\cdot\left[ \rho\,\mathbf{u}\,\mathbf{u}^{T} +\mathbb{I} P - \mathcal{T}^{\mathrm{EM}} \right] = 0 \\ \frac{\partial \left(E + E^{\mathrm{EM}} \right)}{\partial t} + \nabla\cdot\left[ \left(E + P \right) \mathbf{u} + \mathbf{S}^{\mathrm{EM}} \right] = 0 \\ \end{align}\]

The electromagnetic part of the system is solved in USim as:

\[\notag \begin{align} \frac{\partial \mathbf{b}^{\mathrm{plasma}}}{\partial t} - \nabla\times\mathbf{e} + \nabla \psi= 0 \\ \frac{\partial \mathbf{e}^{\mathrm{plasma}}}{\partial t} + f^2 c^2_{\mathrm{fast}} \nabla\times\mathbf{b} + \nabla \phi= - f^2 c^2_{\mathrm{fast}} \mu^{1/2}_0 \sum_{\mathrm{species}} \mathbf{J}^{\mathrm{species}} \\ \frac{\partial \psi}{\partial t} + \nabla\cdot\left[ c^{2}_{\mathrm{fast}} \mathbf{b}^{\mathrm{plasma}} \right] = 0 \\ \frac{\partial \phi}{\partial t} + \nabla\cdot\left[ c^{2}_{\mathrm{fast}} \mathbf{e}^{\mathrm{plasma}} \right] = \mu^{-1/2}_0 \sum_{\mathrm{species}} q^{\mathrm{species}} \\ \end{align}\]

Here, we have written \(c^2 = f^2 c^2_{\mathrm{fast}} = (\epsilon_0 \mu_0)^{-1}\), where \(f\) is a dimensionless number that defines the ratio of the speed of light to the fatest wave in the mesh and we have further defined \(\psi = \mu^{-1/2}_0 \Psi\) and \(\Phi = \mu^{-1/2}_0 \Phi\).

In order to close the electromagnetic part of the equations, a model for the current density and charge is required. An example of such a model that is provided with USim is mhdSrc. However, the user is also free to construct their own closure that returns:

\[\notag \begin{align} \mu^{-1/2}_0 \sum_{\mathrm{species}} q^{\mathrm{species}}; \;\; \mu^{1/2}_0 \sum_{\mathrm{species}} \mathbf{J}^{\mathrm{species}} \end{align}\]

Parameters

lightSpeed (float, optional)
The speed of light in m/s. Used to specify the speed of light in the fluid momentum and energy equations. Defaults to 2.99792458e8.
lightSpeedFactor (float, optional)
Dimensionless number, used to specify the ratio of the speed of light to the fastest wave speed in the grid. Defaults to 1.0e3.
basementPressure (float, optional)
The minimum pressure allowed. Pressures below this value will be replaced with this value for primitive state, eigensystem and flux computations. Defaults to zero.
basementDensity (float, optional)
The minimum density allowed. Densities below this value will be replaced with this value for primitive state, eigensystem and flux computations. Defaults to zero.
gasGamma (float, optional)
Specifies the adiabatic index (ratio of specific heats), \(\gamma\). Defaults to 5/3.
externalEfield (string, optional)
Specifies the name of the data structure containing the externally computed electric field, \(\mathbf{e}^{\mathrm{external}}\).
externalBfield (string, optional)
Specifies the name of the data structure containing the externally computed magnetic field, \(\mathbf{b}^{\mathrm{external}}\).

Parent Updater Data

in (string vector, required)
Vector of Conserved Quantities (nodalArray, 12-components, required)

The vector of conserved quantities, \(\mathbf{q}\) has 9 entries:

  1. \(\rho\): mass density
  2. \(\rho\,u_{\hat{\mathbf{i}}} + c^{-2} {S}^\mathrm{EM}_{\hat{\mathbf{i}}} = \left(\rho \mathbf{u} + c^{-2} \mathbf{S}^\mathrm{EM} \right) \cdot \hat{\mathbf{i}}\): total momentum density in the \(\hat{\mathbf{i}}\) direction
  3. \(\rho\,u_{\hat{\mathbf{j}}} + c^{-2} {S}^\mathrm{EM}_{\hat{\mathbf{j}}} = \left(\rho \mathbf{u} + c^{-2} \mathbf{S}^\mathrm{EM} \right) \cdot \hat{\mathbf{j}}\): total momentum density in the \(\hat{\mathbf{j}}\) direction
  4. \(\rho\,u_{\hat{\mathbf{k}}} + c^{-2} {S}^\mathrm{EM}_{\hat{\mathbf{k}}} = \left(\rho \mathbf{u} + c^{-2} \mathbf{S}^\mathrm{EM} \right) \cdot \hat{\mathbf{k}}\): total momentum density in the \(\hat{\mathbf{k}}\) direction
  5. \(E + E^{\mathrm{EM}} = \frac{P}{\gamma -1} + \tfrac{1}{2}\rho|\mathbf{u}|^2 + E^{\mathrm{EM}}\): total energy density
  6. \(b_{\hat{\mathbf{i}}} = \mathbf{b} \cdot \hat{\mathbf{i}} = \mu^{-1/2}_{0} \mathbf{B} \cdot \hat{\mathbf{i}}\): magnetic field normalized by permeability of free-space in the \(\hat{\mathbf{i}}\) direction
  7. \(b_{\hat{\mathbf{j}}} = \mathbf{b} \cdot \hat{\mathbf{j}} = \mu^{-1/2}_{0} \mathbf{B} \cdot \hat{\mathbf{j}}\): magnetic field normalized by permeability of free-space in the \(\hat{\mathbf{j}}\) direction
  8. \(b_{\hat{\mathbf{k}}} = \mathbf{b} \cdot \hat{\mathbf{k}} = \mu^{-1/2}_{0} \mathbf{B} \cdot \hat{\mathbf{k}}\): magnetic field normalized by permeability of free-space in the \(\hat{\mathbf{k}}\) direction
  9. \(e_{\hat{\mathbf{i}}} = \mathbf{e} \cdot \hat{\mathbf{i}} = \mu^{-1/2}_{0} \mathbf{E} \cdot \hat{\mathbf{i}}\): electric field normalized by permeability of free-space in the \(\hat{\mathbf{i}}\) direction
  10. \(e_{\hat{\mathbf{j}}} = \mathbf{e} \cdot \hat{\mathbf{j}} = \mu^{-1/2}_{0} \mathbf{E} \cdot \hat{\mathbf{j}}\): electric field normalized by permeability of free-space in the \(\hat{\mathbf{j}}\) direction
  11. \(e_{\hat{\mathbf{k}}} = \mathbf{e} \cdot \hat{\mathbf{k}} = \mu^{-1/2}_{0} \mathbf{E} \cdot \hat{\mathbf{k}}\): electric field normalized by permeability of free-space in the \(\hat{\mathbf{k}}\) direction
  12. \(\psi\): magnetic field correction potential
  13. \(\phi\): electric field correction potential
Fastest Wave Speed (dynVector, 1-component, required)
The fastest wave speed across the entire simulation domain, \(c_{\mathrm{fast}}\). Can be computed using hyperbolic (1d, 2d, 3d) (see below).
Externally Computed Electric Field (nodalArray, 3-components, optional)

Additional contribution to the electric field, \(\mathbf{e}^{\mathrm{external}}\), which is not evolved by Ampere’s equation, but does contribution to the induction equation, the Lorentz force and the work done on the plasma. The data structure containing \(\mathbf{e}^{\mathrm{external}}\) is specified by the “externalEField” option described below.

  1. \({e}^{\mathrm{external}}_{\hat{\mathbf{i}}} = \mathbf{e}^{\mathrm{external}} \cdot \hat{\mathbf{i}} = \mu^{-1/2}_{0} \mathbf{E}^{\mathrm{external}} \cdot \hat{\mathbf{i}}\): “externally” computed electric field normalized by permeability of free-space in the \(\hat{\mathbf{i}}\) direction.
  2. \({e}^{\mathrm{external}}_{\hat{\mathbf{j}}} =\mathbf{e}^{\mathrm{external}} \cdot \hat{\mathbf{j}} = \mu^{-1/2}_{0} \mathbf{E}^{\mathrm{external}} \cdot \hat{\mathbf{j}}\):”externally” computed electric field normalized by permeability of free-space in the \(\hat{\mathbf{j}}\) direction
  3. \({e}^{\mathrm{external}}_{\hat{\mathbf{k}}} = \mathbf{e}^{\mathrm{external}} \cdot \hat{\mathbf{k}} = \mu^{-1/2}_{0} \mathbf{E}^{\mathrm{external}} \cdot \hat{\mathbf{k}}\): “externally” computed electric field normalized by permeability of free-space in the \(\hat{\mathbf{k}}\) direction
Externally Computed Magnetic Field (nodalArray, 3-components, optional)

Additional contribution to the magnetic field, \(\mathbf{b}^{\mathrm{external}}\), which is not evolved by the induction equation, but does contribute to Ampere’s equation, the Lorentz force and the work done on the plasma. The data structure containing \(\mathbf{b}^{\mathrm{external}}\) is specified by the “externalBField” option described below.

  1. \({b}^{\mathrm{external}}_{\hat{\mathbf{i}}} = \mathbf{b}^{\mathrm{external}} \cdot \hat{\mathbf{i}} = \mu^{-1/2}_{0} \mathbf{B}^{\mathrm{external}} \cdot \hat{\mathbf{i}}\): magnetic field normalized by permeability of free-space in the \(\hat{\mathbf{i}}\) direction
  2. \({b}^{\mathrm{external}}_{\hat{\mathbf{j}}} =\mathbf{b}^{\mathrm{external}} \cdot \hat{\mathbf{j}} = \mu^{-1/2}_{0} \mathbf{B}^{\mathrm{external}} \cdot \hat{\mathbf{j}}\): magnetic field normalized by permeability of free-space in the \(\hat{\mathbf{j}}\) direction
  3. \({b}^{\mathrm{external}}_{\hat{\mathbf{k}}} = \mathbf{b}^{\mathrm{external}} \cdot \hat{\mathbf{k}} = \mu^{-1/2}_{0} \mathbf{B}^{\mathrm{external}} \cdot \hat{\mathbf{k}}\): magnetic field normalized by permeability of free-space in the \(\hat{\mathbf{k}}\) direction
out (string vector, required)

For the gasDynamicMhdDednerEqn, one of four output variables are computed, depending on whether the equation is combined with an updater capable of computing fluxes (classicMusclUpdater (1d, 2d, 3d)), primitive variables (computePrimitiveState(1d, 2d, 3d)), the time step associated with the CFL condition (timeStepRestrictionUpdater (1d, 2d, 3d)) or the fastest wave speed in the grid (hyperbolic (1d, 2d, 3d)).

Vector of Fluxes (nodalArray, 12-components)

When combined with an updater that computes \(\nabla \cdot \mathcal{F}\left(\mathbf{w} \right)\) (e.g. classicMusclUpdater (1d, 2d, 3d)), the equation system returns:

  1. \(\nabla \cdot \mathcal{F}\left( \rho \right)\): mass flux
  2. \(\nabla \cdot \mathcal{F}\left( \rho\,u_{\hat{\mathbf{i}}} + c^{-2} {S}^\mathrm{EM}_{\hat{\mathbf{i}}} \right)\): \(\hat{\mathbf{i}}\) momentum flux
  3. \(\nabla \cdot \mathcal{F}\left( \rho\,u_{\hat{\mathbf{j}}} + c^{-2} {S}^\mathrm{EM}_{\hat{\mathbf{j}}} \right)\): \(\hat{\mathbf{j}}\) momentum flux
  4. \(\nabla \cdot \mathcal{F}\left( \rho\,u_{\hat{\mathbf{k}}} + c^{-2} {S}^\mathrm{EM}_{\hat{\mathbf{k}}} \right)\): \(\hat{\mathbf{k}}\) momentum flux
  5. \(\nabla \cdot \mathcal{F}\left( E \right)\): total energy flux
  6. \(\nabla \cdot \mathcal{F}\left( b_{\hat{\mathbf{i}}} \right)\): \(\hat{\mathbf{i}}\) magnetic field flux
  7. \(\nabla \cdot \mathcal{F}\left( b_{\hat{\mathbf{j}}} \right)\): \(\hat{\mathbf{j}}\) magnetic field flux
  8. \(\nabla \cdot \mathcal{F}\left( b_{\hat{\mathbf{k}}} \right)\): \(\hat{\mathbf{k}}\) magnetic field flux
  9. \(\nabla \cdot \mathcal{F}\left( e_{\hat{\mathbf{i}}} \right)\): \(\hat{\mathbf{i}}\) electric field flux
  10. \(\nabla \cdot \mathcal{F}\left( e_{\hat{\mathbf{j}}} \right)\): \(\hat{\mathbf{j}}\) electric field flux
  11. \(\nabla \cdot \mathcal{F}\left( e_{\hat{\mathbf{k}}} \right)\): \(\hat{\mathbf{k}}\) electric field flux
  12. \(\nabla \cdot \mathcal{F}\left(\psi \right)\): magnetic correction potential flux
  13. \(\nabla \cdot \mathcal{F}\left(\phi \right)\): electric correction potential flux
Vector of Primitive States (nodalArray, 9-components)

When combined with an updater that computes \(\mathbf{w} = \mathbf{w}(\mathbf{q})\) (e.g. computePrimitiveState(1d, 2d, 3d)), the equation systen returns:

  1. \(\rho\): mass density
  2. \(u_{\hat{\mathbf{i}}} = \mathbf{u} \cdot \hat{\mathbf{i}}\): velocity in the \(\hat{\mathbf{i}}\) direction
  3. \(u_{\hat{\mathbf{j}}} = \mathbf{u} \cdot \hat{\mathbf{j}}\): velocity in the \(\hat{\mathbf{j}}\) direction
  4. \(u_{\hat{\mathbf{k}}} = \mathbf{u} \cdot \hat{\mathbf{k}}\): velocity in the \(\hat{\mathbf{k}}\) direction
  5. \(P = \rho\epsilon(\gamma-1)\): ideal gas pressure
  6. \(b_{\hat{\mathbf{i}}} = \mathbf{b} \cdot \hat{\mathbf{i}} = \mu^{-1/2}_{0} \mathbf{B} \cdot \hat{\mathbf{i}}\): magnetic field normalized by permeability of free-space in the \(\hat{\mathbf{i}}\) direction
  7. \(b_{\hat{\mathbf{j}}} = \mathbf{b} \cdot \hat{\mathbf{j}} = \mu^{-1/2}_{0} \mathbf{B} \cdot \hat{\mathbf{j}}\): magnetic field normalized by permeability of free-space in the \(\hat{\mathbf{j}}\) direction
  8. \(b_{\hat{\mathbf{k}}} = \mathbf{b} \cdot \hat{\mathbf{k}} = \mu^{-1/2}_{0} \mathbf{B} \cdot \hat{\mathbf{k}}\): magnetic field normalized by permeability of free-space in the \(\hat{\mathbf{k}}\) direction
  9. \(e_{\hat{\mathbf{i}}} = \mathbf{e} \cdot \hat{\mathbf{i}} = \mu^{-1/2}_{0} \mathbf{E} \cdot \hat{\mathbf{i}}\): electric field normalized by permeability of free-space in the \(\hat{\mathbf{i}}\) direction
  10. \(e_{\hat{\mathbf{j}}} = \mathbf{e} \cdot \hat{\mathbf{j}} = \mu^{-1/2}_{0} \mathbf{E} \cdot \hat{\mathbf{j}}\): electric field normalized by permeability of free-space in the \(\hat{\mathbf{j}}\) direction
  11. \(e_{\hat{\mathbf{k}}} = \mathbf{e} \cdot \hat{\mathbf{k}} = \mu^{-1/2}_{0} \mathbf{E} \cdot \hat{\mathbf{k}}\): electric field normalized by permeability of free-space in the \(\hat{\mathbf{k}}\) direction
  12. \(\psi\): magnetic field correction potential
  13. \(\phi\): electric field correction potential
Time Step (dynVector, 1-component)
When combined with timeStepRestrictionUpdater (1d, 2d, 3d), the equation system returns the time step consisten with the CFL condition across the entire simulation domain.
Fastest Wave Speed (dynVector, 1-component)
When combined with hyperbolic (1d, 2d, 3d), the equation system returns the fastest wave speed across the entire simulation domain, \(c_{\mathrm{fast}}\).

Examples

The following block demonstrates the gasDynamicMaxwellDednerEqn used in combination with classicMusclUpdater (1d, 2d, 3d) to compute \(\nabla \cdot \mathcal{F}\left(\mathbf{w} \right)\) with an externally supplied magnetic field:

<Updater hyper>
  kind=classicMuscl1d
  onGrid=domain

  # input nodal component arrays
  in=[q   backgroundB]

  # output nodal component array
  out=[qnew]

  # input dynVector containing fastest wave speed
  waveSpeeds=[waveSpeed]

  # the numerical flux to use
  numericalFlux= hlldFlux

  # CFL number to use
  cfl=0.3
  # Form of variables to limit
  variableForm= primitive

  # Limiter; one per input nodal component array
  limiter=[minmod   minmod]

  # list of equations to solve
  equations=[mhd]

  <Equation mhd>
    kind = gasDynamicMaxwellDednerEqn
    gasGamma = GAS_GAMMA
    lightSpeedFactor  =LIGHT_SPEED_FACTOR
    externalBfield = EXTERNAL_FIELD
    basementPressure = BASEMENT_PRESSURE
    basementDensity = BASEMENT_DENSITY
  </Equation>

  <Source mhdSrc>
    kind = gasDynamicMhdDednerSrc
    scalarConductivity = $1.0/OHMIC_RESISTIVITY$
  </Source>

  <Source mhdClean>
    kind = mhdSrc
    model = mhdDednerEqn
    momentumEnergySource = 1
    inputVariables = [q,divB,gradPsi]
  </Source>

</Updater>

The following block demonstrates the gasDynamicMaxwellDednerEqn used in combination with computePrimitiveState(1d, 2d, 3d) to compute \(\mathbf{w}\):

<Updater computePrimitiveState>
  kind = computePrimitiveState$NDIM$d

  onGrid = domain
  # input array
  in = [q]

  # ouput data-structures
  out = [w]

  <Equation fluid>
    kind = gasDynamicMaxwellDednerEqn
    gasGamma = GAS_GAMMA
    lightSpeedFactor  =LIGHT_SPEED_FACTOR
    externalBfield = EXTERNAL_FIELD
    basementPressure = BASEMENT_PRESSURE
    basementDensity = BASEMENT_DENSITY
  </Equation>

</Updater>

The following block demonstrates the gasDynamicMhdDednerEqn used in combination with timeStepRestrictionUpdater (1d, 2d, 3d) and hyperbolic (1d, 2d, 3d) to compute \(c_{\mathrm{fast}}\) with an externally supplied magnetic field:

<Updater getWaveSpeed>
  kind=timeStepRestrictionUpdater1d
  onGrid=domain

  # input nodal component arrays
  in=[q   backgroundB]

  # output dynVector containing fastest wave speed
  waveSpeeds=[waveSpeed]

  # list of equations to compute fastest wave speed for
  restrictions=[idealMhd]

  # courant condition to apply to the timestep
  courantCondition=1.0

 <TimeStepRestriction idealMhd>
   kind = hyperbolic1d
   model = gasDynamicMaxwellDednerEqn
    gasGamma = GAS_GAMMA
    lightSpeedFactor  =LIGHT_SPEED_FACTOR
    externalBfield = EXTERNAL_FIELD
    basementPressure = BASEMENT_PRESSURE
    basementDensity = BASEMENT_DENSITY
 </TimeStepRestriction>
</Updater>