WARPXM v1.10.0
Loading...
Searching...
No Matches
wxm::apps::five_moment Namespace Reference

Namespaces

namespace  bc
 
namespace  limiters
 

Classes

class  AppSetupHelper
 Helper class that automates common setup tasks for 5-moment fluid physics applications. More...
 
class  braginskii_flux
 
class  ConservedPrimitiveConversion
 
class  DensityDiffusionCylSource
 
class  DensityDiffusionFlux
 
class  Euler
 Euler flux. More...
 
class  Euler1D_Riemann_Problem_Analytic_Solution
 Class that calculates the analytic solution to the Riemann Problem for the Euler Equations. More...
 
class  EulerCylSource
 Implements the source terms in the Euler equations arising from a cylindrical geometry. More...
 
class  field_source_t
 
class  five_moment_to_kinetic
 5-Moment to Kinetic Boundary Condition App More...
 
class  general_source_t
 
class  InterspeciesCollisions
 
class  IntraspeciesCollisions
 Implements Navier-Stokes style collisions. More...
 
class  IntraspeciesCyl
 Implements the source terms for intraspecies collisions arising from a cylindrical geometry. More...
 
class  navier_stokes_t
 
class  reaction_source_mol_t
 
class  reaction_source_t
 Implements source terms for atomic reactions. More...
 

Functions

real getRhoFloored (const real *q, const real rho_min=std::numeric_limits< real >::epsilon())
 Get the floored mass density.
 
real getPressureFloored (const real *q, const real gas_gamma, const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 Get the pressure from Euler (5-moment) variables using floors for pressure and density.
 
real getTemperatureFloored (const real *q, const real gas_gamma, const real A=1.0, const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 Get the temperature from Euler (5-moment) variables using floors for pressure and density.
 
void NCtoC (const constants_5moment_t *pC, const real QNC[5], real QC[5])
 Converts nonconserved (primitive) variables to conserved variables.
 
void CtoNC (const constants_5moment_t *pC, const real QC[5], real QNC[5], const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 Converts conserved variables to nonconserved (primitive) variables.
 
void NCtoCDifferential (const constants_5moment_t *pC, const real QNC_center[5], const real dQNC[5], const real QC_center[5], real dQC[5])
 This currently isn't used anywhere.
 
void rotateVector (const real v[5], const solverVariables_t *pSV, real Rv[5])
 Rotates a 5-moment variable set from world frame to local frame.
 
void antirotateVector (const real Rv[5], const solverVariables_t *pSV, real v[5])
 Rotates a 5-moment variable set from local frame to world frame.
 
void eulerFlux (const constants_5moment_t *pC, const real QC[5], real flux[5], const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 Flux along x-axis for Euler fluid equations using conserved variables.
 
void eulerRoeFluxSolver (const constants_5moment_t *pC, const real QC_l[5], const real QC_r[5], real *s, real amdq[5], real apdq[5], const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 Solves the full Roe decomposition for the euler fluid problem using conserved variables - assumes an entropy fix - function is stolen from WARPX.
 
void eulerRoeAverages (const constants_5moment_t *pC, const real ql[5], const real qr[5], const real &pl, const real &pr, real &rho, real &u, real &v, real &w, real &enth, real &a, const real rho_min=std::numeric_limits< real >::epsilon())
 
void hartenHymanEntropyFix (const constants_5moment_t *pC, const real ql[5], const real qr[5], const real s[3], const real wave[][3], real amdq[5], real apdq[5], const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 
real eulerNumericalFlux_x (const constants_5moment_t *pC, const real QC_l[5], const real QC_r[5], void(*numerical_flux_scheme)(const constants_5moment_t *pC, const real QC_l[5], const real QC_r[5], real *s, real amdq[5], real apdq[5], const real rho_min, const real press_min), real flux[5], const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 Calculates the Euler Numerical flux given a scheme along the x-axis which is to be used as a numerical flux.
 
real eulerNumericalFlux (const constants_5moment_t *pC, const solverVariables_t *pSV, const real QC_l[5], const real QC_r[5], void(*numerical_flux_scheme)(const constants_5moment_t *pC, const real QC_l[5], const real QC_r[5], real *s, real amdq[5], real apdq[5], const real rho_min, const real press_min), real numericalFlux[5], const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 Calculates the Numerical Flux for the Euler fluid equations using Roe's Method as given by Leveque in "Finite Volume Methods for Hyperbolic Problems" equation 15.40 for a given numerical Flux scheme (Roe, HLLE, etc.)
 
void eulerHLLEFluxSolver (const constants_5moment_t *pC, const real QC_l[5], const real QC_r[5], real *s, real amdq[5], real apdq[5], const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 Calculates the HLLE flux for the Euler fluid equations.
 
void eulerHLLECFluxSolver (const constants_5moment_t *pC, const real QC_l[5], const real QC_r[5], real *s, real amdq[5], real apdq[5], const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 
real eulerRusanovFlux (const constants_5moment_t *pC, const solverVariables_t *pSV, const real QC_l[5], const real QC_r[5], real numericalFlux[5], const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 Calculates the Rusanov flux for the Euler fluid equations directly.
 
real eulerInternalFlux (const constants_5moment_t *pC, const solverVariables_t *pSV, const real QC[5], std::vector< std::vector< real > > &internalFlux, const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 Calculates the internal euler flux at a given point.
 
void navier_stokes_anisotropic_pressure_tensor (const real viscosity, const real *q, const real *grad_rho, const real *grad_p[3], real pi[3][3])
 
void navier_stokes_viscous_flux_tensor (const real *q, const solverVariables_t *pSV, const real pi[][3], const real h[3], std::vector< std::vector< real > > &flux_tensor)
 
void navier_stokes_internal_flux (const constants_5moment_t *pC, const solverVariables_t *pSV, const real mu, const real k, const real *q, const real u[3], const real *grad_u[3], const real div_u, const real grad_T[3], std::vector< std::vector< real > > &internal_flux, const real rho_min=std::numeric_limits< real >::epsilon())
 
real navier_stokes_numerical_flux (constants_5moment_t *pC, const solverVariables_t *pSV, const real viscosity, const real thermal_conductivity, const real *q_l, const real *q_r, const real *aux_l, const real *aux_r, real *numerical_flux, const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 
real interspeciesCollisions (const real nu_ab, const constants_5moment_t *pC_a, const constants_5moment_t *pC_b, const real *q_a, const real *q_b, real *source, const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 
void calcGradT (const constants_5moment_t *pC, const real *q, const real *aux, const real *grad_u, real *grad_T, const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 Computes temperature gradients from 5 moment conserved variable gradients/values.
 
real calcDivU (const real *q, const real *aux, const real rho_min=std::numeric_limits< real >::epsilon())
 
void calcGradU (const real *q, const real *aux, real *grad_u[3], const real rho_min=std::numeric_limits< real >::epsilon())
 
void get_mu_and_k (const constants_5moment_t *pC, const real *q, const bool is_neutral, const real collision_coefficient, const real minimum_frequency, const real diffusivity_min, const real diffusivity_max, const bool is_const_mu, const real const_mu, const bool is_const_k, const real const_k, real &mu, real &k, const real rho_min=std::numeric_limits< real >::epsilon(), const real press_min=std::numeric_limits< real >::epsilon())
 
real positivity_preserving_source_dt (const real *q, const real *dqdt, real lambda, const real rho_min=std::numeric_limits< real >::epsilon())
 Computes the maximum allowable dt for a source term to preserve positivity for the positivity-preserving DG scheme.
 

Function Documentation

◆ antirotateVector()

void wxm::apps::five_moment::antirotateVector ( const real  Rv[5],
const solverVariables_t pSV,
real  v[5] 
)

Rotates a 5-moment variable set from local frame to world frame.

Parameters
RvVariables in local frame
pSVSolver variables which define the rotation
vResultant variables in world frame

◆ calcDivU()

real wxm::apps::five_moment::calcDivU ( const real q,
const real aux,
const real  rho_min = std::numeric_limits< real >::epsilon() 
)

◆ calcGradT()

void wxm::apps::five_moment::calcGradT ( const constants_5moment_t pC,
const real q,
const real aux,
const real grad_u,
real grad_T,
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

Computes temperature gradients from 5 moment conserved variable gradients/values.

Parameters
qconserved variables
auxgradients of conserved variables
grad_uvelocity gradients (pre-computed using grad_u)
grad_Toutput temperature gradients
massspecies normalized mass
gammaspecies gas gamma

◆ calcGradU()

void wxm::apps::five_moment::calcGradU ( const real q,
const real aux,
real grad_u[3],
const real  rho_min = std::numeric_limits< real >::epsilon() 
)

◆ CtoNC()

void wxm::apps::five_moment::CtoNC ( const constants_5moment_t pC,
const real  QC[5],
real  QNC[5],
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

Converts conserved variables to nonconserved (primitive) variables.

Parameters
QCResultant conserved variables
QNCNonconserved variables
rho_minNormalized minimum mass density

◆ eulerFlux()

void wxm::apps::five_moment::eulerFlux ( const constants_5moment_t pC,
const real  QC[5],
real  flux[5],
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

Flux along x-axis for Euler fluid equations using conserved variables.

Parameters
pCConstants that define this specific gas
QCConserved variables
fluxResultant flux
rho_minNormalized minimum mass density.
press_minNormalized minimum mass density.

◆ eulerHLLECFluxSolver()

void wxm::apps::five_moment::eulerHLLECFluxSolver ( const constants_5moment_t pC,
const real  QC_l[5],
const real  QC_r[5],
real s,
real  amdq[5],
real  apdq[5],
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

◆ eulerHLLEFluxSolver()

void wxm::apps::five_moment::eulerHLLEFluxSolver ( const constants_5moment_t pC,
const real  QC_l[5],
const real  QC_r[5],
real s,
real  amdq[5],
real  apdq[5],
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

Calculates the HLLE flux for the Euler fluid equations.

Parameters
pCConstants that define this specific gas
pSVSolver variables
QC_lConserved variables of left state
QC_rConserved variables of right state
numericalFluxResultant numerical flux out of the face
Returns
Suggested dt for this time step

◆ eulerInternalFlux()

real wxm::apps::five_moment::eulerInternalFlux ( const constants_5moment_t pC,
const solverVariables_t pSV,
const real  QC[5],
std::vector< std::vector< real > > &  internalFlux,
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

Calculates the internal euler flux at a given point.

Parameters
pCConstants that define this specific gas
pSVSolver variables
QCConserved variables
internalFluxResultant internal flux
Returns
Suggested dt for this time step

◆ eulerNumericalFlux()

real wxm::apps::five_moment::eulerNumericalFlux ( const constants_5moment_t pC,
const solverVariables_t pSV,
const real  QC_l[5],
const real  QC_r[5],
void(*)(const constants_5moment_t *pC, const real QC_l[5], const real QC_r[5], real *s, real amdq[5], real apdq[5], const real rho_min, const real press_min)  numerical_flux_scheme,
real  numericalFlux[5],
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

Calculates the Numerical Flux for the Euler fluid equations using Roe's Method as given by Leveque in "Finite Volume Methods for Hyperbolic Problems" equation 15.40 for a given numerical Flux scheme (Roe, HLLE, etc.)

Parameters
pCConstants that define this specific gas
pSVSolver variables
QC_lConserved variables of left state
QC_rConserved variables of right state
numerical_flux_schemeNumerical Flux Scheme
numericalFluxResultant numerical flux out of the face
Returns
Suggested dt for this time step

◆ eulerNumericalFlux_x()

real wxm::apps::five_moment::eulerNumericalFlux_x ( const constants_5moment_t pC,
const real  QC_l[5],
const real  QC_r[5],
void(*)(const constants_5moment_t *pC, const real QC_l[5], const real QC_r[5], real *s, real amdq[5], real apdq[5], const real rho_min, const real press_min)  numerical_flux_scheme,
real  flux[5],
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

Calculates the Euler Numerical flux given a scheme along the x-axis which is to be used as a numerical flux.

Parameters
pCConstants that define this specific gas
QC_lConserved variables for left state
QC_rConserved variables for right state
fluxResultant numerical flux traveling to the right
numerical_flux_schemeNumerical Flux Scheme
Returns
Fastest absolute wave speed

◆ eulerRoeAverages()

void wxm::apps::five_moment::eulerRoeAverages ( const constants_5moment_t pC,
const real  ql[5],
const real  qr[5],
const real pl,
const real pr,
real rho,
real u,
real v,
real w,
real enth,
real a,
const real  rho_min = std::numeric_limits< real >::epsilon() 
)

◆ eulerRoeFluxSolver()

void wxm::apps::five_moment::eulerRoeFluxSolver ( const constants_5moment_t pC,
const real  QC_l[5],
const real  QC_r[5],
real s,
real  amdq[5],
real  apdq[5],
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

Solves the full Roe decomposition for the euler fluid problem using conserved variables - assumes an entropy fix - function is stolen from WARPX.

Parameters
pCConstants that define this specific gas
QC_lConserved variables describing left state
QC_rConserved variables describing right state
sResultant wave speeds
amdqLeft traveling fluctuations
apdqRight traveling fluctuation

◆ eulerRusanovFlux()

real wxm::apps::five_moment::eulerRusanovFlux ( const constants_5moment_t pC,
const solverVariables_t pSV,
const real  QC_l[5],
const real  QC_r[5],
real  numericalFlux[5],
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

Calculates the Rusanov flux for the Euler fluid equations directly.

Parameters
pCConstants that define this specific gas
pSVSolver variables
QC_lConserved variables of left state
QC_rConserved variables of right state
numericalFluxResultant numerical flux out of the face
Returns
Suggested dt for this time step

◆ get_mu_and_k()

void wxm::apps::five_moment::get_mu_and_k ( const constants_5moment_t pC,
const real q,
const bool  is_neutral,
const real  collision_coefficient,
const real  minimum_frequency,
const real  diffusivity_min,
const real  diffusivity_max,
const bool  is_const_mu,
const real  const_mu,
const bool  is_const_k,
const real  const_k,
real mu,
real k,
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

◆ getPressureFloored()

real wxm::apps::five_moment::getPressureFloored ( const real q,
const real  gas_gamma,
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

Get the pressure from Euler (5-moment) variables using floors for pressure and density.

The floored pressure is given as

\[ p_{fl} = \max(p, p_{min}), \]

where

\[ p = (\gamma - 1) \left( e - \frac{1}{2} \max(\rho, \rho_{min}) v^2) \]

Assumes \( \rho_{min}, \; p_{min} > 0 \).

Parameters
qThe raw 5-moment variables.
gas_gammaRatio of specific heats.
rho_minNormalized minimum mass density.
press_minNormalized minimum mass density.
Returns
Pressure calculated assuming pressure and density floors.

◆ getRhoFloored()

real wxm::apps::five_moment::getRhoFloored ( const real q,
const real  rho_min = std::numeric_limits< real >::epsilon() 
)

Get the floored mass density.

The floored mass density is given as

\[ \rho_{fl} = \max(\rho, \rho_{min}), \]

Assumes \( \rho_{min} > 0 \).

Parameters
qThe raw 5-moment variables.
rho_minNormalized minimum mass density.
Returns
Floored mass density.

◆ getTemperatureFloored()

real wxm::apps::five_moment::getTemperatureFloored ( const real q,
const real  gas_gamma,
const real  A = 1.0,
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

Get the temperature from Euler (5-moment) variables using floors for pressure and density.

Temperature is calculated assuming ideal gas law: \( p = n T \)

Therefore

\[ T = \frac{p_{fl}}{\rho_{fl}} A \]

where

\[ p_{fl} = \max(p, p_{min}), \quad \text{and} \quad \rho_{fl} = \max(\rho, \rho_{min}). \]

Assumes \( \rho_{min}, \; p_{min} > 0 \).

Parameters
qThe raw 5-moment variables.
gas_gammaRatio of specific heats.
AProton normalized ion mass.
rho_minNormalized minimum mass density.
press_minNormalized minimum mass density.
Returns
Temperature calculated assuming pressure and density floors.

◆ hartenHymanEntropyFix()

void wxm::apps::five_moment::hartenHymanEntropyFix ( const constants_5moment_t pC,
const real  ql[5],
const real  qr[5],
const real  s[3],
const real  wave[][3],
real  amdq[5],
real  apdq[5],
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

◆ interspeciesCollisions()

real wxm::apps::five_moment::interspeciesCollisions ( const real  nu_ab,
const constants_5moment_t pC_a,
const constants_5moment_t pC_b,
const real q_a,
const real q_b,
real source,
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

◆ navier_stokes_anisotropic_pressure_tensor()

void wxm::apps::five_moment::navier_stokes_anisotropic_pressure_tensor ( const real  viscosity,
const real q,
const real grad_rho,
const real grad_p[3],
real  pi[3][3] 
)

◆ navier_stokes_internal_flux()

void wxm::apps::five_moment::navier_stokes_internal_flux ( const constants_5moment_t pC,
const solverVariables_t pSV,
const real  mu,
const real  k,
const real q,
const real  u[3],
const real grad_u[3],
const real  div_u,
const real  grad_T[3],
std::vector< std::vector< real > > &  internal_flux,
const real  rho_min = std::numeric_limits< real >::epsilon() 
)

◆ navier_stokes_numerical_flux()

real wxm::apps::five_moment::navier_stokes_numerical_flux ( constants_5moment_t pC,
const solverVariables_t pSV,
const real  viscosity,
const real  thermal_conductivity,
const real q_l,
const real q_r,
const real aux_l,
const real aux_r,
real numerical_flux,
const real  rho_min = std::numeric_limits< real >::epsilon(),
const real  press_min = std::numeric_limits< real >::epsilon() 
)

◆ navier_stokes_viscous_flux_tensor()

void wxm::apps::five_moment::navier_stokes_viscous_flux_tensor ( const real q,
const solverVariables_t pSV,
const real  pi[][3],
const real  h[3],
std::vector< std::vector< real > > &  flux_tensor 
)

◆ NCtoC()

void wxm::apps::five_moment::NCtoC ( const constants_5moment_t pC,
const real  QNC[5],
real  QC[5] 
)

Converts nonconserved (primitive) variables to conserved variables.

Parameters
QNCNonconserved variables
QCResultant conserved variables
rho_minNormalized minimum mass density

◆ NCtoCDifferential()

void wxm::apps::five_moment::NCtoCDifferential ( const constants_5moment_t pC,
const real  QNC_center[5],
const real  dQNC[5],
const real  QC_center[5],
real  dQC[5] 
)

This currently isn't used anywhere.

Can it be deleted?

◆ positivity_preserving_source_dt()

real wxm::apps::five_moment::positivity_preserving_source_dt ( const real q,
const real dqdt,
real  lambda,
const real  rho_min = std::numeric_limits< real >::epsilon() 
)

Computes the maximum allowable dt for a source term to preserve positivity for the positivity-preserving DG scheme.

Parameters
qthe conserved variables: density, momentum and total (thermal + kinetic) energy, in an array of length 5.
dqdtthe time derivatives of q
lambdathe proportion of the positive densities which should remain after the timestep.

The timestep is computed to be the maximum which satisfies several constraints:

  • That \(\rho^n + 2 * \Delta t \dot{\rho^n} > 0\). The factor of 2 is required for the positivity-preserving DG method.
  • That \(\rho^n + \Delta t \dot{\rho^n} > \lambda \rho^n\).
  • That \(p(\mathbf{q}^n + 2\Delta t \mathbf{q}^n) > 0\).
  • That \(p(\mathbf{q}^n + 2\Delta t \mathbf{q}^n) > (1 - 2(1 - \lambda)) p(\mathbf{q}^n)\).

For details on these constraints and how we find a dt that satisfies them, see the writeup.

For MHD< variables, one should subtract off the magnetic energy \(|B|^2/2\) from the fifth element of the variable, q[4], before passing it in to this function.

See also
Writeup document, "Timestep limit for source terms".

◆ rotateVector()

void wxm::apps::five_moment::rotateVector ( const real  v[5],
const solverVariables_t pSV,
real  Rv[5] 
)

Rotates a 5-moment variable set from world frame to local frame.

Parameters
vVariables in world frame
pSVSolver variables which define the rotation
RvResultant variables in local frame