Contents
Moments App: Multifluid-moment-Maxwell model¶
Summary of model equations¶
The Moment
app solves high-moment multifluid equations on a Cartesian grid, coupled to Maxwell’s equations through the Lorentz force.
Each fluid could be either five-moment, where the plasma pressure is assumed to a scalar (see [Hakim+2006]), or ten-moment, where the full anisotropic and nongyrotropic plasma pressure tensors (see [Hakim2008]) are evolved.
The hyperbolic system is solved in converative forms of the five-moment (Euler) equations, and/or the ten-moment equations, and the perfectly hyperbolic Maxwell’s equations.
The sources that couples the plasma species momenta and electromagnetic fields are described here and more comprehensively in [Wang+2020].
This App solves the hyperbolic and source parts parts of the coupled system separately and apply high accuracy schemes on both.
Ignoring sources, the homogeneous equations of the fluid moments and electromagneic fields are solved separately using a high-resolution wave-propagation finite-volume method described in [Hakim+2006]. The main time-step size constraint comes from the speed of light.
The sources are evolved using a locally implicit, time-centered solver to step over the constraining scales like the Debye length and plasma frequency, etc. See Handling two-fluid five-moment and ten-moment source terms or [Wang+2020] for more details.
Additonial sources can be added if needed, for example, for the ten-moment model, we may apply the closure to relax the pressure tensor towards a scalar pressure (see [Wang+2015]).
We then apply an Strang-type operator-splitting sequence to combine the hyperbolic and source parts to achive second-order accuracy in time:
\[\exp\left(\mathcal{L}_{S}\Delta t/2\right)\exp\left(\mathcal{L}_{H}\Delta t\right)\exp\left(\mathcal{L}_{S}\Delta t/2\right).\]Here, we represent the homogeneous update schematically as the operator \(\exp\left(\mathcal{L}_{H}\Delta t\right)\) and the source update as \(\exp\left(\mathcal{L}_{S}\Delta t\right)\).
Overall structure of the Moments app¶
--------------
-- Preamble --
--------------
-- The Moments app wraps fluid and field objects, and tells the
-- the program how to evolve and couple them
local Moments = require("App.PlasmaOnCartGrid").Moments()
local TenMoment = require "Eq.TenMoment" -- TenMoment or Euler
-- Create the app
momentApp = Moments.App {
------------
-- COMMON --
------------
-- basic parameters, e.g., time step, grid, domain decomposition
-- Description of each species: names are arbitrary
electron = Moments.Species {
-- species parameters, equations, and boundary conditions
},
-- Repeat to add more species
hydrogen = Moments.Species { ... },
oxygen = Moments.Species { ... },
-- EM fields (optional, can be omitted for neutral fluids)
field = Moments.Field {
-- EM field parameters, equations, and boundary conditions
},
-- Basic source that couple the fluids and EM fields
emSource = Moments.CollisionlessEmSource {
-- names of the species to be coupled
species = {"electron", "hydorgen", "oxygen"},
-- other specifications
},
-- Additional sources if needed
elc10mSource = Moments.TenMomentRelaxSource {
species = {"elctron"},
-- other specifications
},
}
-- run the app
momentApp:run()
Examples¶
Five-moment modeling of the GEM challenge magnetic reconnection problem.
Ten-moment modeling of the GEM challenge magnetic reconnection problem. This simulation uses a simplified closure appropriate for this problem.
Basic parameters¶
Parameter |
Description |
Default |
---|---|---|
|
If set to true, log messages are written to log file |
|
|
End time of simulation |
|
|
Initial suggested time-step. Adjusted as simulation progresses. |
|
|
Number of frames of data to write. Initial conditions are always written. For more fine-grained control over species and field output, see below. |
|
|
CDIM length table with lower-left configuration space coordinates |
|
|
CDIM length table with upper-right configuration space coordinates |
|
|
CDIM length table with number of configuration space cells |
|
|
CFL number to use. This parameter should be avoided and ``cflFrac`` used instead. |
Determined from |
|
Fraction (usually 1.0) to multiply CFL determined time-step. |
Determined from |
|
Hard limit of time step size. |
|
|
The multifluid-Maxwell model currently only supports the dimensional-
splitting finite-volume method, i.e., |
|
|
CDIM length table with number of processors to use in each configuration space direction. |
|
|
Set to |
|
|
Periodic directions. Note: X is 1, Y is 2 and Z is 3. E.g., |
|
Note
In general, you should not specify
cfl
orcflFrac
, unless either doing tests or explicitly controlling the time-step. The app will determine the time-step automatically.When
useShared=true
thedecompCuts
must specify the number of nodes and not number of processors. That is, the total number of processors will be determined fromdecompCuts
and the number of threads per node.
Species parameters¶
Parameter |
Description |
Default |
---|---|---|
|
Species charge (ignored for neutral particles) |
|
|
Species mass (ignored for neutral particles) |
|
|
The type of default moment equation for this species, e.g.,
|
|
|
Backup equation that guarantees positivity in case it is violated when
the default |
|
|
Function with signature |
|
|
Length two table with BCs in X direction. See details on BCs below. |
|
|
Length two table with BCs in Y direction. Only needed if CDIM>1 |
|
|
Length two table with BCs in Z direction. Only needed if CDIM>2 |
|
|
If set to |
|
|
If set to |
|
Electromagnetic field parameters¶
Parameter |
Description |
Default |
---|---|---|
|
These many field outputs will be written during simulation. If
not specified, top-level |
|
|
Vacuum permittivity (\(\epsilon_0\)) |
|
|
Vacuum permeability (\(\mu_0\)) |
|
|
Factor specifying speed for magnetic field divergence error correction |
|
|
Factor specifying speed for electric field divergence error correction |
|
|
Function with signature |
|
|
Length two table with BCs in X direction. See details on BCs below. |
|
|
Length two table with BCs in Y direction. Only needed if CDIM>1 |
|
|
Length two table with BCs in Z direction. Only needed if CDIM>2 |
|
|
If set to |
|
|
If set to |
|
App output¶
The app will write snapshots of moments for each species and the EM fields at specified time intervals. Diagnostics like integrated fluid moments and field energy are recorded for each time-step and written in one file for each species/field object.
The output format is ADIOS BP files. Say your input file is called “5m.lua” and your species are called “elc” and “ion”. Then, over specified time invertals the app will write out the following files:
5m_elc_N.bp
5m_ion_N.bp
5m_field_N.bp
Where N
is the frame number (frame 0 is the initial
conditions). Note that if a species or the field is not evolved, then
only initial conditions will be written unless the forceWrite
option
is set to true
.
In addition, integrated moments for each species are written:
vlasov_elc_intMom_N.bp
For the field, the electromagnetic energy components \(E_x^2\), \(E_y^2\), \(E_z^2\), \(B_x^2\), \(B_y^2\), and \(B_z^2\) (integrated over configuration space) are stored in the file:
vlasov_fieldEnergy_N.bp
These can be plotted using postgkyl in the usual way.
References¶
Hakim, A., Loverich, J., & Shumlak, U. (2006). A high resolution wave propagation scheme for ideal Two-Fluid plasma equations. Journal of Computational Physics, 219, 418–442. https://doi.org/10.1016/j.jcp.2006.03.036
Hakim, A. H. (2008). Extended MHD modeling with the ten-moment equations. Journal of Fusion Energy, 27, 36–43. https://doi.org/10.1007/s10894-007-9116-z
Wang, L., Hakim, A. H., Ng, J., Dong, C., & Germaschewski, K. (2020). Exact and locally implicit source term solvers for multifluid-Maxwell systems. Journal of Computational Physics. https://doi.org/10.1016/j.jcp.2020.109510
Wang, L., Hakim, A. H. A. H., Bhattacharjee, A., & Germaschewski, K. (2015). Comparison of multi-fluid moment models with particle-in-cell simulations of collisionless magnetic reconnection. Physics of Plasmas, 22(1), 012108. https://doi.org/10.1063/1.4906063