ryujin 2.1.1 revision 0348cbb53a3e4b1da2a4c037e81f88f2d21ce219
|
#include <source/euler_aeos/hyperbolic_system.h>
Public Types | |
using | ScalarNumber = typename get_value_type< Number >::type |
Public Member Functions | |
HyperbolicSystemView (const HyperbolicSystem &hyperbolic_system) | |
template<int dim2, typename Number2 > | |
auto | view () const |
Access to runtime parameters | |
DEAL_II_ALWAYS_INLINE const std::string & | equation_of_state () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | reference_density () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | vacuum_state_relaxation_small () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | vacuum_state_relaxation_large () const |
DEAL_II_ALWAYS_INLINE bool | compute_strict_bounds () const |
Surrogate functions for computing various interpolatory | |
physical quantities that are needed for Riemann solver, indicator and limiter. | |
Number | surrogate_specific_entropy (const state_type &U, const Number &gamma_min) const |
Number | surrogate_harten_entropy (const state_type &U, const Number &gamma_min) const |
state_type | surrogate_harten_entropy_derivative (const state_type &U, const Number &eta, const Number &gamma_min) const |
Number | surrogate_gamma (const state_type &U, const Number &p) const |
Number | surrogate_pressure (const state_type &U, const Number &gamma) const |
Number | surrogate_speed_of_sound (const state_type &U, const Number &gamma) const |
bool | is_admissible (const state_type &U) const |
Special functions for boundary states | |
template<int component> | |
state_type | prescribe_riemann_characteristic (const state_type &U, const Number &p, const state_type &U_bar, const Number &p_bar, const dealii::Tensor< 1, dim, Number > &normal) const |
template<typename Lambda > | |
state_type | apply_boundary_conditions (const dealii::types::boundary_id id, const state_type &U, const dealii::Tensor< 1, dim, Number > &normal, const Lambda &get_dirichlet_data) const |
State transformations | |
template<typename ST > | |
state_type | expand_state (const ST &state) const |
template<typename ST > | |
state_type | from_initial_state (const ST &initial_state) const |
state_type | from_primitive_state (const state_type &primitive_state) const |
state_type | to_primitive_state (const state_type &state) const |
template<typename Lambda > | |
state_type | apply_galilei_transform (const state_type &state, const Lambda &lambda) const |
Low-level access to the selected equation of state. | |
static constexpr bool | have_gamma = false |
static constexpr bool | have_eos_interpolation_b = true |
DEAL_II_ALWAYS_INLINE Number | eos_pressure (const Number &rho, const Number &e) const |
DEAL_II_ALWAYS_INLINE Number | eos_specific_internal_energy (const Number &rho, const Number &p) const |
DEAL_II_ALWAYS_INLINE Number | eos_temperature (const Number &rho, const Number &e) const |
DEAL_II_ALWAYS_INLINE Number | eos_speed_of_sound (const Number &rho, const Number &e) const |
DEAL_II_ALWAYS_INLINE ScalarNumber | eos_interpolation_b () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | eos_interpolation_pinfty () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | eos_interpolation_q () const |
Computing precomputed quantities | |
static constexpr unsigned int | n_precomputation_cycles = 2 |
template<typename DISPATCH , typename SPARSITY > | |
void | precomputation_loop (unsigned int cycle, const DISPATCH &dispatch_check, const SPARSITY &sparsity_simd, StateVector &state_vector, unsigned int left, unsigned int right) const |
Flux computations | |
static constexpr bool | have_high_order_flux = false |
flux_type | f (const state_type &U, const Number &p) const |
flux_contribution_type | flux_contribution (const PrecomputedVector &pv, const InitialPrecomputedVector &piv, const unsigned int i, const state_type &U_i) const |
flux_contribution_type | flux_contribution (const PrecomputedVector &pv, const InitialPrecomputedVector &piv, const unsigned int *js, const state_type &U_j) const |
state_type | flux_divergence (const flux_contribution_type &flux_i, const flux_contribution_type &flux_j, const dealii::Tensor< 1, dim, Number > &c_ij) const |
state_type | high_order_flux_divergence (const flux_contribution_type &flux_i, const flux_contribution_type &flux_j, const dealii::Tensor< 1, dim, Number > &c_ij) const =delete |
Computing stencil source terms | |
static constexpr bool | have_source_terms = false |
state_type | nodal_source (const PrecomputedVector &pv, const unsigned int i, const state_type &U_i, const ScalarNumber tau) const =delete |
state_type | nodal_source (const PrecomputedVector &pv, const unsigned int *js, const state_type &U_j, const ScalarNumber tau) const =delete |
Computing derived physical quantities | |
Number | filter_vacuum_density (const Number &rho) const |
static Number | density (const state_type &U) |
static dealii::Tensor< 1, dim, Number > | momentum (const state_type &U) |
static Number | total_energy (const state_type &U) |
static Number | internal_energy (const state_type &U) |
static state_type | internal_energy_derivative (const state_type &U) |
A view of the HyperbolicSystem that makes methods available for a given dimension dim
and choice of number type Number
(which can be a scalar float, or double, as well as a VectorizedArray holding packed scalars.
Intended usage:
Definition at line 130 of file hyperbolic_system.h.
using ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::ScalarNumber = typename get_value_type<Number>::type |
The underlying scalar number type.
Definition at line 154 of file hyperbolic_system.h.
using ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::state_type = dealii::Tensor<1, problem_dimension, Number> |
Storage type for a (conserved) state vector \(\boldsymbol U\).
Definition at line 338 of file hyperbolic_system.h.
using ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::flux_type = dealii::Tensor<1, problem_dimension, dealii::Tensor<1, dim, Number> > |
Storage type for the flux \(\mathbf{f}\).
Definition at line 343 of file hyperbolic_system.h.
using ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::flux_contribution_type = flux_type |
The storage type used for flux contributions.
Definition at line 349 of file hyperbolic_system.h.
using ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::precomputed_type = std::array<Number, n_precomputed_values> |
Array type used for precomputed values.
Definition at line 389 of file hyperbolic_system.h.
using ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::initial_precomputed_type = std::array<Number, n_initial_precomputed_values> |
Array type used for precomputed initial values.
Definition at line 409 of file hyperbolic_system.h.
using ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::StateVector = Vectors:: StateVector<ScalarNumber, problem_dimension, n_precomputed_values> |
A compound state vector.
Definition at line 421 of file hyperbolic_system.h.
using ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::HyperbolicVector = Vectors::MultiComponentVector<ScalarNumber, problem_dimension> |
MulticomponentVector for storing the hyperbolic state vector:
Definition at line 427 of file hyperbolic_system.h.
using ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::PrecomputedVector = Vectors::MultiComponentVector<ScalarNumber, n_precomputed_values> |
MulticomponentVector for storing a vector of precomputed states:
Definition at line 433 of file hyperbolic_system.h.
using ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::InitialPrecomputedVector = Vectors::MultiComponentVector<ScalarNumber, n_initial_precomputed_values> |
MulticomponentVector for storing a vector of precomputed initial states:
Definition at line 440 of file hyperbolic_system.h.
|
inline |
Constructor taking a reference to the underlying HyperbolicSystem
Definition at line 137 of file hyperbolic_system.h.
|
inline |
Create a modified view from the current one:
Definition at line 146 of file hyperbolic_system.h.
|
inline |
Definition at line 161 of file hyperbolic_system.h.
|
inline |
Definition at line 166 of file hyperbolic_system.h.
|
inline |
Definition at line 172 of file hyperbolic_system.h.
|
inline |
Definition at line 178 of file hyperbolic_system.h.
|
inline |
Definition at line 183 of file hyperbolic_system.h.
|
inline |
For a given density \(\rho\) and specific internal energy \(e\) return the pressure \(p\).
Definition at line 198 of file hyperbolic_system.h.
|
inline |
For a given density \(\rho\) and pressure \(p\) return the specific internal energy \(e\).
Definition at line 219 of file hyperbolic_system.h.
|
inline |
For a given density \(\rho\) and specific internal energy \(e\) return the temperature \(T\).
Definition at line 238 of file hyperbolic_system.h.
|
inline |
For a given density \(\rho\) and specific internal energy \(e\) return the sound speed \(a\).
Definition at line 259 of file hyperbolic_system.h.
|
inline |
Return the interpolatory covolume \(b_{\text{interp}}\).
Definition at line 277 of file hyperbolic_system.h.
|
inline |
Return the interpolatory reference pressure \(p_{\infty}\).
Definition at line 286 of file hyperbolic_system.h.
|
inline |
Return the interpolatory reference specific internal energy \(q\).
Definition at line 296 of file hyperbolic_system.h.
void ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::precomputation_loop | ( | unsigned int | cycle, |
const DISPATCH & | dispatch_check, | ||
const SPARSITY & | sparsity_simd, | ||
StateVector & | state_vector, | ||
unsigned int | left, | ||
unsigned int | right | ||
) | const |
Step 0: precompute values for hyperbolic update. This routine is called within our usual loop() idiom in HyperbolicModule
|
inlinestatic |
For a given (2+dim dimensional) state vector U
, return the density U[0]
Definition at line 1014 of file hyperbolic_system.h.
|
inline |
Given a density rho
this function returns 0 if the magniatude of rho is smaller or equal than relaxation_large * rho_cutoff. Otherwise rho is returned unmodified. Here, rho_cutoff is the reference density multiplied by eps.
Definition at line 1022 of file hyperbolic_system.h.
|
inlinestatic |
For a given (2+dim dimensional) state vector U
, return the momentum vector [U[1], ..., U[1+dim]]
.
Definition at line 1036 of file hyperbolic_system.h.
|
inlinestatic |
For a given (2+dim dimensional) state vector U
, return the total energy U[1+dim]
Definition at line 1047 of file hyperbolic_system.h.
|
inlinestatic |
For a given (2+dim dimensional) state vector U
, compute and return the internal energy \(\varepsilon = (\rho e)\).
Definition at line 1055 of file hyperbolic_system.h.
|
inlinestatic |
For a given (2+dim dimensional) state vector U
, compute and return the derivative of the internal energy \(\varepsilon = (\rho e)\).
Definition at line 1069 of file hyperbolic_system.h.
|
inline |
For a given (2+dim dimensional) state vector U
, compute and return a (scaled) surrogate specific entropy
\[ e^{(\gamma_{\text{min}} - 1)s} = \frac{\rho\,(e-q)-p_{\infty}(1-b\rho)}{\rho^\gamma_{\text{min}}} (1 - b\,\rho)^{\gamma_{\text{min}} -1}. \]
Definition at line 1096 of file hyperbolic_system.h.
References ryujin::pow().
|
inline |
For a given (2+dim dimensional) state vector U
, compute and return a surrogate Harten-type entropy
\[ \eta = (1-b\,\rho)^{\frac{\gamma_{\text{min}-1}}{\gamma_{\text{min}}+1}} \big(\rho^2 (e-q) - \rho p_{\infty}(1-b\,\rho)\big) ^{1/(\gamma_{\text{min}}+1)} \]
Definition at line 1116 of file hyperbolic_system.h.
References ryujin::positive_part(), and ryujin::pow().
|
inline |
For a given (2+dim dimensional) state vector U
, compute and return the derivative \(\eta'\) of the Harten-type entropy
\[ \eta = (1-b\,\rho)^{\frac{\gamma_{\text{min}-1}}{\gamma_{\text{min}}+1}} \big(\rho^2 (e-q) - \rho p_{\infty}(1-b\,\rho)\big) ^{1/(\gamma_{\text{min}}+1)} \]
Definition at line 1143 of file hyperbolic_system.h.
References ryujin::pow().
|
inline |
For a given (2+dim dimensional) state vector U
and pressure p
, compute a surrogate gamma:
\[ \gamma(\rho, e, p) = 1 + \frac{(p + p_{\infty})(1 - b \rho)} {\rho (e-q) - p_{\infty}(1-b \rho)} \]
This function is used in various places to create interpolations of the pressure.
Definition at line 1199 of file hyperbolic_system.h.
References ryujin::EulerAEOS::safe_division().
|
inline |
For a given (2+dim dimensional) state vector U
and gamma gamma
, compute a surrogate pressure:
\[ p(\rho, e, \gamma) = (\gamma - 1) \frac{\rho (e - q)}{1 - b \rho} -\gamma\,p_{\infty} \]
This function is the complementary function to surrogate_gamma(), meaning the following property holds true:
Definition at line 1218 of file hyperbolic_system.h.
References ryujin::positive_part(), and ryujin::EulerAEOS::safe_division().
|
inline |
For a given (2+dim dimensional) state vector U
and gamma gamma
, compute a surrogate speed of sound:
\begin{align} c^2(\rho, e, \gamma) = \frac{\gamma (p + p_\infty)}{\rho X} = \frac{\gamma (\gamma -1)[\rho (e - q) - p_\infty X]}{\rho X^2} \end{align}
Definition at line 1237 of file hyperbolic_system.h.
References ryujin::positive_part().
|
inline |
Returns whether the state U
is admissible. If U
is a vectorized state then U
is admissible if all vectorized values are admissible.
Definition at line 1257 of file hyperbolic_system.h.
state_type ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::prescribe_riemann_characteristic | ( | const state_type & | U, |
const Number & | p, | ||
const state_type & | U_bar, | ||
const Number & | p_bar, | ||
const dealii::Tensor< 1, dim, Number > & | normal | ||
) | const |
Decomposes a given state U
into Riemann invariants and then replaces the first or second Riemann characteristic from the one taken from U_bar
state. Note that the U_bar
state is just the prescribed dirichlet values.
state_type ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::apply_boundary_conditions | ( | const dealii::types::boundary_id | id, |
const state_type & | U, | ||
const dealii::Tensor< 1, dim, Number > & | normal, | ||
const Lambda & | get_dirichlet_data | ||
) | const |
Apply boundary conditions.
For the compressible Euler equations we have:
|
inline |
Given a state U
and a pressure p
compute the flux
\[ \begin{pmatrix} \textbf m \\ \textbf v\otimes \textbf m + p\mathbb{I}_d \\ \textbf v(E+p) \end{pmatrix}, \]
Definition at line 1538 of file hyperbolic_system.h.
|
inline |
Given a state U_i
and an index i
compute flux contributions.
Intended usage:
For the Euler equations we simply compute f(U_i)
.
Definition at line 1560 of file hyperbolic_system.h.
|
inline |
Definition at line 1574 of file hyperbolic_system.h.
|
inline |
Given flux contributions flux_i
and flux_j
compute the flux (-f(U_i) - f(U_j)
Definition at line 1588 of file hyperbolic_system.h.
References ryujin::add(), and ryujin::contract().
|
delete |
|
delete |
|
delete |
state_type ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::expand_state | ( | const ST & | state | ) | const |
Given a state vector associated with a different spatial dimensions than the current one, return an "expanded" version of the state vector associated with dim spatial dimensions where the momentum vector of the conserved state state
is expaned with zeros to a total length of dim entries.
state_type ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::from_initial_state | ( | const ST & | initial_state | ) | const |
Given an initial state [rho, u_1, ..., u_d, p] return a conserved state [rho, m_1, ..., m_d, E]. Most notably, the specific equation of state oracle is queried to convert the pressure value into a specific internal energy.
|
inline |
Given a primitive state [rho, u_1, ..., u_d, e] return a conserved state.
Definition at line 1640 of file hyperbolic_system.h.
|
inline |
Given a conserved state return a primitive state [rho, u_1, ..., u_d, e]
Definition at line 1663 of file hyperbolic_system.h.
state_type ryujin::EulerAEOS::HyperbolicSystemView< dim, Number >::apply_galilei_transform | ( | const state_type & | state, |
const Lambda & | lambda | ||
) | const |
Transform the current state according to a given operator lambda
acting on a dim dimensional momentum (or velocity) vector.
|
staticconstexpr |
constexpr boolean used in the EulerInitialStates namespace
Definition at line 306 of file hyperbolic_system.h.
|
staticconstexpr |
constexpr boolean used in the EulerInitialStates namespace
Definition at line 311 of file hyperbolic_system.h.
|
staticconstexpr |
The dimension of the state space.
Definition at line 333 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the conserved state as a string.
Definition at line 355 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the primitive state as a string.
Definition at line 370 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputed values.
Definition at line 384 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the precomputed values.
Definition at line 394 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputed initial values.
Definition at line 404 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the precomputed values.
Definition at line 415 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputation cycles.
Definition at line 453 of file hyperbolic_system.h.
|
staticconstexpr |
The low-order and high-order fluxes are the same:
Definition at line 714 of file hyperbolic_system.h.
|
staticconstexpr |
We do not have source terms
Definition at line 727 of file hyperbolic_system.h.