gasDynamicMhdDednerEqn

Defines the equations of inviscid fluid dynamics coupled to pre-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}^{\mathrm{species}} \\ \frac{\partial \mathbf{B^{\mathrm{plasma}}}}{\partial t} + \nabla\times\mathbf{E} + \nabla \psi= 0 \\ \frac{\partial \psi}{\partial t} + \nabla\cdot\left[ c^{2}_{\mathrm{fast}} \mathbf{b} \right] = 0 \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) \\ \mathbf{e} = - \mathbf{u} \times \mathbf{b} + \mathbf{e}^{\mathrm{external}} = \mu^{-1/2}_{0} \left( - \mathbf{u} \times \mathbf{B} + \mathbf{E}^{\mathrm{external}} \right) \end{align}\]

Here, \(\mathbf{b}^{\mathrm{plasma}}\) is the magnetic field induced in the plasma by the inductive electric field, \(\mathbf{e}\), while \(\mathbf{e}^{\mathrm{external}}\) and \(\mathbf{b}^{\mathrm{external}}\) are electromagnetic fields computed “externally” to the pre-Maxwell equations.

With these identitifications, the gasDynamicMhdDednerEqn 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 \\ \frac{\partial \mathbf{b^{\mathrm{plasma}}}}{\partial t} + \nabla\times\mathbf{e} + \nabla \psi= 0 \\ \frac{\partial \psi}{\partial t} + \nabla\cdot\left[ c^{2}_{\mathrm{fast}} \mathbf{b} \right] = 0 \end{align}\]

This flux-conservative formulation is implemented in USim.

Parameters

lightSpeed (float, optional)
The speed of light in m/s. Defaults to 2.99792458e8.
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, 9-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. \(\psi\): 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 terms in the generalized Ohm’s law, \(\mathbf{E}^{\mathrm{external}}\), computed “externally” to the ideal magnetohydrodynamic system. 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 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, 9-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(\psi \right)\): 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. \(\psi\): 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 mhdDednerEqn 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=gasDynamicMhdDednerEqn
    gasGamma=1.4
    externalBfield="backgroundB"
  </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=gasDynamicMhdDednerEqn
    gasGamma= 1.4
    externalBfield=True
    includeInTimeStep=False
  </TimeStepRestriction>
</Updater>