![]() |
ryujin 2.1.1 revision ef0fcd4010d109b860652ece4a7b8963fb7d46b1
|
#include <source/euler/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 |
template<typename DISPATCH , typename SPARSITY > | |
DEAL_II_ALWAYS_INLINE void | precomputation_loop (unsigned int cycle, const DISPATCH &dispatch_check, const SPARSITY &sparsity_simd, StateVector &state_vector, unsigned int left, unsigned int right, const bool skip_constrained_dofs) const |
template<int component> | |
DEAL_II_ALWAYS_INLINE auto | linearized_eigenvector (const state_type &U, const dealii::Tensor< 1, dim, Number > &normal) const -> std::array< state_type, 2 > |
template<int component> | |
DEAL_II_ALWAYS_INLINE auto | prescribe_riemann_characteristic (const state_type &U, const state_type &U_bar, const dealii::Tensor< 1, dim, Number > &normal) const -> state_type |
template<typename Lambda > | |
DEAL_II_ALWAYS_INLINE auto | apply_boundary_conditions (dealii::types::boundary_id id, const state_type &U, const dealii::Tensor< 1, dim, Number > &normal, const Lambda &get_dirichlet_data) const -> state_type |
template<typename ST > | |
auto | expand_state (const ST &state) const -> state_type |
template<typename ST > | |
DEAL_II_ALWAYS_INLINE auto | from_initial_state (const ST &initial_state) const -> state_type |
template<typename Lambda > | |
auto | apply_galilei_transform (const state_type &state, const Lambda &lambda) const -> state_type |
template<typename DISPATCH , typename SPARSITY > | |
DEAL_II_ALWAYS_INLINE void | precomputation_loop (unsigned int cycle, const DISPATCH &dispatch_check, const SPARSITY &sparsity_simd, StateVector &state_vector, unsigned int left, unsigned int right, const bool skip_constrained_dofs) const |
template<int component> | |
DEAL_II_ALWAYS_INLINE auto | 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 -> state_type |
template<typename Lambda > | |
DEAL_II_ALWAYS_INLINE auto | apply_boundary_conditions (dealii::types::boundary_id id, const state_type &U, const dealii::Tensor< 1, dim, Number > &normal, const Lambda &get_dirichlet_data) const -> state_type |
template<typename ST > | |
auto | expand_state (const ST &state) const -> state_type |
template<typename ST > | |
DEAL_II_ALWAYS_INLINE auto | from_initial_state (const ST &initial_state) const -> state_type |
template<typename Lambda > | |
auto | apply_galilei_transform (const state_type &state, const Lambda &lambda) const -> state_type |
template<typename DISPATCH , typename SPARSITY > | |
DEAL_II_ALWAYS_INLINE void | precomputation_loop (unsigned int cycle, const DISPATCH &dispatch_check, const SPARSITY &sparsity_simd, StateVector &state_vector, unsigned int left, unsigned int right, const bool skip_constrained_dofs) const |
template<typename Lambda > | |
DEAL_II_ALWAYS_INLINE auto | apply_boundary_conditions (dealii::types::boundary_id id, const state_type &U, const dealii::Tensor< 1, dim, Number > &, const Lambda &get_dirichlet_data) const -> state_type |
template<typename DISPATCH , typename SPARSITY > | |
DEAL_II_ALWAYS_INLINE void | precomputation_loop (unsigned int cycle, const DISPATCH &dispatch_check, const SPARSITY &sparsity_simd, StateVector &state_vector, unsigned int left, unsigned int right, const bool skip_constrained_dofs) const |
template<int component> | |
DEAL_II_ALWAYS_INLINE auto | prescribe_riemann_characteristic (const state_type &U, const state_type &U_bar, const dealii::Tensor< 1, dim, Number > &normal) const -> state_type |
template<typename Lambda > | |
DEAL_II_ALWAYS_INLINE auto | 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_type |
template<typename ST > | |
DEAL_II_ALWAYS_INLINE auto | expand_state (const ST &state) const -> state_type |
template<typename ST > | |
DEAL_II_ALWAYS_INLINE auto | from_initial_state (const ST &initial_state) const -> state_type |
template<typename Lambda > | |
DEAL_II_ALWAYS_INLINE auto | apply_galilei_transform (const state_type &state, const Lambda &lambda) const -> state_type |
Access to runtime parameters | |
DEAL_II_ALWAYS_INLINE ScalarNumber | gamma () 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 |
Special functions for boundary states | |
template<int component> | |
std::array< state_type, 2 > | linearized_eigenvector (const state_type &U, const dealii::Tensor< 1, dim, Number > &normal) const |
template<int component> | |
state_type | prescribe_riemann_characteristic (const state_type &U, const state_type &U_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 |
Access to cached inverses | |
A collection of commonly used expressions with gamma that would otherwise need to be recomputed many times putting unnecessary pressure on the div/sqrt ALU unit. | |
static constexpr bool | have_gamma = true |
static constexpr bool | have_eos_interpolation_b = false |
DEAL_II_ALWAYS_INLINE ScalarNumber | gamma_inverse () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | gamma_plus_one_inverse () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | gamma_minus_one_inverse () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | gamma_minus_one_over_gamma_plus_one () const |
Computing precomputed quantities | |
static constexpr unsigned int | n_precomputation_cycles = 1 |
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 bool skip_constrained_dofs=true) const |
Flux computations | |
static constexpr bool | have_high_order_flux = false |
flux_type | f (const state_type &U) const |
flux_contribution_type | flux_contribution (const PrecomputedVector &pv, const InitialPrecomputedVector &ipv, const unsigned int i, const state_type &U_i) const |
flux_contribution_type | flux_contribution (const PrecomputedVector &pv, const InitialPrecomputedVector &ipv, 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 |
Number | pressure (const state_type &U) const |
Number | speed_of_sound (const state_type &U) const |
Number | specific_entropy (const state_type &U) const |
Number | harten_entropy (const state_type &U) const |
state_type | harten_entropy_derivative (const state_type &U) const |
Number | mathematical_entropy (const state_type &U) const |
state_type | mathematical_entropy_derivative (const state_type &U) const |
bool | is_admissible (const state_type &U) 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 105 of file hyperbolic_system.h.
using ryujin::Euler::HyperbolicSystemView< dim, Number >::ScalarNumber = typename get_value_type<Number>::type |
The underlying scalar number type.
Definition at line 129 of file hyperbolic_system.h.
using ryujin::Euler::HyperbolicSystemView< dim, Number >::state_type = dealii::Tensor<1, problem_dimension, Number> |
Storage type for a (conserved) state vector \(\boldsymbol U\).
Definition at line 224 of file hyperbolic_system.h.
using ryujin::Euler::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 229 of file hyperbolic_system.h.
using ryujin::Euler::HyperbolicSystemView< dim, Number >::flux_contribution_type = flux_type |
The storage type used for flux contributions.
Definition at line 235 of file hyperbolic_system.h.
using ryujin::Euler::HyperbolicSystemView< dim, Number >::precomputed_type = std::array<Number, n_precomputed_values> |
Array type used for precomputed values.
Definition at line 275 of file hyperbolic_system.h.
using ryujin::Euler::HyperbolicSystemView< dim, Number >::initial_precomputed_type = std::array<Number, n_initial_precomputed_values> |
Array type used for precomputed initial values.
Definition at line 291 of file hyperbolic_system.h.
using ryujin::Euler::HyperbolicSystemView< dim, Number >::StateVector = Vectors:: StateVector<ScalarNumber, problem_dimension, n_precomputed_values> |
A compound state vector.
Definition at line 303 of file hyperbolic_system.h.
using ryujin::Euler::HyperbolicSystemView< dim, Number >::HyperbolicVector = Vectors::MultiComponentVector<ScalarNumber, problem_dimension> |
MulticomponentVector for storing the hyperbolic state vector:
Definition at line 309 of file hyperbolic_system.h.
using ryujin::Euler::HyperbolicSystemView< dim, Number >::PrecomputedVector = Vectors::MultiComponentVector<ScalarNumber, n_precomputed_values> |
MulticomponentVector for storing a vector of precomputed states:
Definition at line 315 of file hyperbolic_system.h.
using ryujin::Euler::HyperbolicSystemView< dim, Number >::InitialPrecomputedVector = Vectors::MultiComponentVector<ScalarNumber, n_initial_precomputed_values> |
MulticomponentVector for storing a vector of precomputed initial states:
Definition at line 322 of file hyperbolic_system.h.
|
inline |
Constructor taking a reference to the underlying HyperbolicSystem
Definition at line 112 of file hyperbolic_system.h.
|
inline |
Create a modified view from the current one:
Definition at line 121 of file hyperbolic_system.h.
|
inline |
Definition at line 136 of file hyperbolic_system.h.
|
inline |
Definition at line 141 of file hyperbolic_system.h.
|
inline |
Definition at line 147 of file hyperbolic_system.h.
|
inline |
Definition at line 153 of file hyperbolic_system.h.
|
inline |
Definition at line 168 of file hyperbolic_system.h.
|
inline |
Definition at line 173 of file hyperbolic_system.h.
|
inline |
Definition at line 178 of file hyperbolic_system.h.
|
inline |
Definition at line 184 of file hyperbolic_system.h.
void ryujin::Euler::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 bool | skip_constrained_dofs = true |
||
) | 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 745 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 753 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 767 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 778 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 786 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 800 of file hyperbolic_system.h.
|
inline |
For a given (2+dim dimensional) state vector U
, compute and return the pressure \(p\).
We assume that the pressure is given by a polytropic equation of state, i.e.,
\[ p = (\gamma - 1)\;(\rho e) \]
Definition at line 827 of file hyperbolic_system.h.
|
inline |
For a given (2+dim dimensional) state vector U
, compute the (physical) speed of sound:
\[ c^2 = \frac{\gamma\,p}{\rho} \]
Definition at line 836 of file hyperbolic_system.h.
|
inline |
For a given (2+dim dimensional) state vector U
, compute and return the (scaled) specific entropy
\[ e^{(\gamma-1)s} = \frac{\rho\,e}{\rho^\gamma}. \]
Definition at line 847 of file hyperbolic_system.h.
References ryujin::pow().
|
inline |
For a given (2+dim dimensional) state vector U
, compute and return the Harten-type entropy
\[ \eta = (\rho^2 e) ^ {1 / (\gamma + 1)}. \]
Definition at line 858 of file hyperbolic_system.h.
References ryujin::pow().
|
inline |
For a given (2+dim dimensional) state vector U
, compute and return the derivative \(\eta'\) of the Harten-type entropy
\[ \eta = (\rho^2 e) ^ {1 / (\gamma + 1)}. \]
Definition at line 873 of file hyperbolic_system.h.
References ryujin::pow().
|
inline |
For a given (2+dim dimensional) state vector U
, compute and return the entropy \(\eta = p^{1/\gamma}\).
Definition at line 910 of file hyperbolic_system.h.
References ryujin::pow().
|
inline |
For a given (2+dim dimensional) state vector U
, compute and return the derivative \(\eta'\) of the entropy \(\eta =
p^{1/\gamma}\).
Definition at line 921 of file hyperbolic_system.h.
References ryujin::pow().
|
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 958 of file hyperbolic_system.h.
std::array< state_type, 2 > ryujin::Euler::HyperbolicSystemView< dim, Number >::linearized_eigenvector | ( | const state_type & | U, |
const dealii::Tensor< 1, dim, Number > & | normal | ||
) | const |
For a given state U
and normal direction normal
returns the n-th pair of left and right eigenvectors of the linearized normal flux.
state_type ryujin::Euler::HyperbolicSystemView< dim, Number >::prescribe_riemann_characteristic | ( | const state_type & | U, |
const state_type & | U_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::Euler::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
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 1176 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 1198 of file hyperbolic_system.h.
|
inline |
Definition at line 1210 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 1222 of file hyperbolic_system.h.
References ryujin::add(), and ryujin::contract().
|
delete |
|
delete |
|
delete |
state_type ryujin::Euler::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::Euler::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].
This function simply calls from_primitive_state() and expand_state().
|
inline |
Given a primitive state [rho, u_1, ..., u_d, p] return a conserved state
Definition at line 1267 of file hyperbolic_system.h.
|
inline |
Given a conserved state return a primitive state [rho, u_1, ..., u_d, p]
Definition at line 1289 of file hyperbolic_system.h.
state_type ryujin::Euler::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.
|
inline |
Definition at line 707 of file hyperbolic_system.h.
References RYUJIN_OMP_FOR.
|
inline |
Definition at line 988 of file hyperbolic_system.h.
|
inline |
Definition at line 1043 of file hyperbolic_system.h.
References ryujin::pow().
|
inline |
Definition at line 1102 of file hyperbolic_system.h.
References ryujin::dirichlet, ryujin::dirichlet_momentum, ryujin::dynamic, ryujin::no_slip, and ryujin::slip.
auto ryujin::Euler::HyperbolicSystemView< dim, Number >::expand_state | ( | const ST & | state | ) | const -> state_type |
Definition at line 1233 of file hyperbolic_system.h.
|
inline |
Definition at line 1257 of file hyperbolic_system.h.
auto ryujin::Euler::HyperbolicSystemView< dim, Number >::apply_galilei_transform | ( | const state_type & | state, |
const Lambda & | lambda | ||
) | const -> state_type |
Definition at line 1309 of file hyperbolic_system.h.
|
inline |
Definition at line 878 of file hyperbolic_system.h.
References RYUJIN_OMP_FOR, and RYUJIN_OMP_SINGLE.
|
inline |
Definition at line 1294 of file hyperbolic_system.h.
References ryujin::pow().
|
inline |
Definition at line 1447 of file hyperbolic_system.h.
References ryujin::dirichlet, ryujin::dirichlet_momentum, ryujin::dynamic, ryujin::no_slip, and ryujin::slip.
auto ryujin::Euler::HyperbolicSystemView< dim, Number >::expand_state | ( | const ST & | state | ) | const -> state_type |
Definition at line 1602 of file hyperbolic_system.h.
|
inline |
Definition at line 1626 of file hyperbolic_system.h.
auto ryujin::Euler::HyperbolicSystemView< dim, Number >::apply_galilei_transform | ( | const state_type & | state, |
const Lambda & | lambda | ||
) | const -> state_type |
Definition at line 1686 of file hyperbolic_system.h.
|
inline |
Definition at line 603 of file hyperbolic_system.h.
References RYUJIN_OMP_FOR.
|
inline |
Definition at line 697 of file hyperbolic_system.h.
References ryujin::dirichlet, ryujin::dirichlet_momentum, ryujin::dynamic, ryujin::no_slip, and ryujin::slip.
|
inline |
Definition at line 670 of file hyperbolic_system.h.
References ryujin::pow(), and RYUJIN_OMP_FOR.
|
inline |
Definition at line 892 of file hyperbolic_system.h.
|
inline |
Definition at line 941 of file hyperbolic_system.h.
References ryujin::dirichlet, ryujin::dirichlet_momentum, ryujin::dynamic, ryujin::no_slip, and ryujin::slip.
|
inline |
Definition at line 1245 of file hyperbolic_system.h.
|
inline |
Definition at line 1267 of file hyperbolic_system.h.
|
inline |
Definition at line 1310 of file hyperbolic_system.h.
|
staticconstexpr |
constexpr boolean used in the EulerInitialStates namespace
Definition at line 193 of file hyperbolic_system.h.
|
staticconstexpr |
constexpr boolean used in the EulerInitialStates namespace
Definition at line 198 of file hyperbolic_system.h.
|
staticconstexpr |
The dimension of the state space.
Definition at line 219 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the conserved state as a string.
Definition at line 241 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the primitive state as a string.
Definition at line 256 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputed values.
Definition at line 270 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the precomputed values.
Definition at line 280 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputed initial values.
Definition at line 286 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the precomputed values.
Definition at line 297 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputation cycles.
Definition at line 335 of file hyperbolic_system.h.
|
staticconstexpr |
The low-order and high-order fluxes are the same
Definition at line 575 of file hyperbolic_system.h.
|
staticconstexpr |
We do not have source terms:
Definition at line 589 of file hyperbolic_system.h.