ryujin 2.1.1 revision d0a94ad2ccc0c4c2e8c2485c52b06b90e2fc9853
|
#include <source/shallow_water/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 ScalarNumber | gravity () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | manning_friction_coefficient () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | reference_water_depth () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | dry_state_relaxation_factor () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | dry_state_relaxation_small () const |
DEAL_II_ALWAYS_INLINE ScalarNumber | dry_state_relaxation_large () const |
Special functions for boundary states | |
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 (primitive states, expanding | |
dimensionality, Galilei transform, etc.) | |
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 |
Types and compile time constants | |
using | state_type = dealii::Tensor< 1, problem_dimension, Number > |
using | vector_type = MultiComponentVector< ScalarNumber, problem_dimension > |
using | flux_type = dealii::Tensor< 1, problem_dimension, dealii::Tensor< 1, dim, Number > > |
using | flux_contribution_type = std::tuple< state_type, Number > |
static constexpr unsigned int | problem_dimension = dim + 1 |
static const auto | component_names |
static const auto | primitive_component_names |
Precomputed quantities | |
using | precomputed_initial_state_type = std::array< Number, n_precomputed_initial_values > |
using | precomputed_initial_vector_type = MultiComponentVector< ScalarNumber, n_precomputed_initial_values > |
using | precomputed_state_type = std::array< Number, n_precomputed_values > |
using | precomputed_vector_type = MultiComponentVector< ScalarNumber, n_precomputed_values > |
static constexpr unsigned int | n_precomputed_initial_values = 1 |
static const auto | precomputed_initial_names |
static constexpr unsigned int | n_precomputed_values = 2 |
static const auto | precomputed_names |
static constexpr unsigned int | n_precomputation_cycles = 1 |
template<typename DISPATCH , typename SPARSITY > | |
void | precomputation_loop (unsigned int cycle, const DISPATCH &dispatch_check, precomputed_vector_type &precomputed_values, const SPARSITY &sparsity_simd, const vector_type &U, unsigned int left, unsigned int right) const |
Flux computations | |
static constexpr bool | have_high_order_flux = true |
flux_type | f (const state_type &U) const |
flux_type | g (const state_type &U) const |
state_type | star_state (const state_type &U, const Number &Z_left, const Number &Z_right) const |
std::array< state_type, 2 > | equilibrated_states (const flux_contribution_type &, const flux_contribution_type &) const |
flux_contribution_type | flux_contribution (const precomputed_vector_type &pv, const precomputed_initial_vector_type &piv, const unsigned int i, const state_type &U_i) const |
flux_contribution_type | flux_contribution (const precomputed_vector_type &pv, const precomputed_initial_vector_type &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 |
state_type | affine_shift (const flux_contribution_type &flux_i, const flux_contribution_type &flux_j, const dealii::Tensor< 1, dim, Number > &c_ij, const Number &d_ij) const |
Computing source terms | |
static constexpr bool | have_source_terms = true |
state_type | manning_friction (const state_type &U, const Number &h_star, const ScalarNumber tau) const |
state_type | nodal_source (const precomputed_vector_type &pv, const unsigned int i, const state_type &U_i, const ScalarNumber tau) const |
state_type | nodal_source (const precomputed_vector_type &pv, const unsigned int *js, const state_type &U_j, const ScalarNumber tau) const |
Computing derived physical quantities | |
Number | inverse_water_depth_mollified (const state_type &U) const |
Number | water_depth_sharp (const state_type &U) const |
Number | inverse_water_depth_sharp (const state_type &U) const |
Number | filter_dry_water_depth (const Number &h) const |
Number | kinetic_energy (const state_type &U) const |
Number | pressure (const state_type &U) const |
Number | speed_of_sound (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 | water_depth (const state_type &U) |
static dealii::Tensor< 1, dim, Number > | momentum (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 100 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystemView< dim, Number >::ScalarNumber = typename get_value_type<Number>::type |
The underlying scalar number type.
Definition at line 124 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystemView< dim, Number >::state_type = dealii::Tensor<1, problem_dimension, Number> |
The storage type used for a (conserved) state vector \(\boldsymbol U\).
Definition at line 190 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystemView< dim, Number >::vector_type = MultiComponentVector<ScalarNumber, problem_dimension> |
MulticomponentVector for storing a vector of conserved states:
Definition at line 195 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystemView< dim, Number >::flux_type = dealii::Tensor<1, problem_dimension, dealii::Tensor<1, dim, Number> > |
The storage type used for the flux \(\mathbf{f}\).
Definition at line 230 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystemView< dim, Number >::flux_contribution_type = std::tuple<state_type, Number> |
The storage type used for flux contributions.
Definition at line 236 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystemView< dim, Number >::precomputed_initial_state_type = std::array<Number, n_precomputed_initial_values> |
Array type used for precomputed initial values.
Definition at line 252 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystemView< dim, Number >::precomputed_initial_vector_type = MultiComponentVector<ScalarNumber, n_precomputed_initial_values> |
MulticomponentVector for storing a vector of precomputed initial states:
Definition at line 259 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystemView< dim, Number >::precomputed_state_type = std::array<Number, n_precomputed_values> |
Array type used for precomputed values.
Definition at line 276 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystemView< dim, Number >::precomputed_vector_type = MultiComponentVector<ScalarNumber, n_precomputed_values> |
MulticomponentVector for storing a vector of precomputed states:
Definition at line 281 of file hyperbolic_system.h.
|
inline |
Constructor taking a reference to the underlying HyperbolicSystem
Definition at line 107 of file hyperbolic_system.h.
|
inline |
Create a modified view from the current one:
Definition at line 116 of file hyperbolic_system.h.
|
inline |
Definition at line 131 of file hyperbolic_system.h.
|
inline |
Definition at line 137 of file hyperbolic_system.h.
|
inline |
Definition at line 142 of file hyperbolic_system.h.
|
inline |
Definition at line 148 of file hyperbolic_system.h.
|
inline |
Definition at line 154 of file hyperbolic_system.h.
|
inline |
Definition at line 160 of file hyperbolic_system.h.
void ryujin::ShallowWater::HyperbolicSystemView< dim, Number >::precomputation_loop | ( | unsigned int | cycle, |
const DISPATCH & | dispatch_check, | ||
precomputed_vector_type & | precomputed_values, | ||
const SPARSITY & | sparsity_simd, | ||
const vector_type & | U, | ||
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 (1+dim dimensional) state vector U
, return the water depth U[0]
Definition at line 712 of file hyperbolic_system.h.
|
inline |
For a given (1+dim dimensional) state vector U
, return a regularized inverse of the water depth. This function returns 2h / (h^2+max(h, h_cutoff)^2), where h_cutoff is the reference water depth multiplied by eps.
Definition at line 720 of file hyperbolic_system.h.
References ryujin::positive_part().
|
inline |
For a given (1+dim dimensional) state vector U
, return the regularized water depth U[0]
This function returns max(h, h_cutoff), where h_cutoff is the reference water depth multiplied by eps.
Definition at line 738 of file hyperbolic_system.h.
|
inline |
For a given (1+dim dimensional) state vector U
, return a regularized inverse of the water depth. This function returns 1 / max(h, h_cutoff), where h_cutoff is the reference water depth multiplied by eps.
Definition at line 754 of file hyperbolic_system.h.
|
inline |
Given a water depth h this function returns 0 if h is in the interval [-relaxation * h_cutoff, relaxation * h_cutoff], otherwise h is returned unmodified. Here, h_cutoff is the reference water depth multiplied by eps.
Definition at line 763 of file hyperbolic_system.h.
|
inlinestatic |
For a given (1+dim dimensional) state vector U
, return the momentum vector [U[1], ..., U[1+dim]]
.
Definition at line 779 of file hyperbolic_system.h.
|
inline |
For a given (1+dim dimensional) state vector U
, compute and return the kinetic energy.
\[ KE = 1/2 |m|^2 / h \]
Definition at line 791 of file hyperbolic_system.h.
|
inline |
For a given (state dimensional) state vector U
, compute and return the hydrostatic pressure \(p\):
\[ p = 1/2 g h^2 \]
Definition at line 803 of file hyperbolic_system.h.
|
inline |
For a given (1+dim dimensional) state vector U
, compute the (physical) speed of sound:
\[ c^2 = g * h \]
Definition at line 814 of file hyperbolic_system.h.
|
inline |
For a given (1+dim dimensional) state vector U
, compute and return the entropy \(\eta = 1/2 g h^2 + 1/2 |m|^2 / h\).
Definition at line 823 of file hyperbolic_system.h.
|
inline |
For a given (1+dim dimensional) state vector U
, compute and return the derivative \(\eta'\) of the entropy defined above.
Definition at line 834 of file hyperbolic_system.h.
|
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 867 of file hyperbolic_system.h.
state_type ryujin::ShallowWater::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.
state_type ryujin::ShallowWater::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.
|
inline |
Given a state U
compute the flux
\[ \begin{pmatrix} \textbf m \\ \textbf v\otimes \textbf m + p\mathbb{I}_d \\ \end{pmatrix}, \]
Definition at line 1015 of file hyperbolic_system.h.
|
inline |
Given a state U
compute the flux \g[ \begin{pmatrix} \textbf m \ \textbf v\otimes \textbf m \ \end{pmatrix}, \f]
Definition at line 1034 of file hyperbolic_system.h.
|
inline |
For a given (1+dim dimensional) state vector U
and left/right topography states Z_left
and Z_right
, return the star_state U_star
Definition at line 1051 of file hyperbolic_system.h.
|
inline |
Given precomputed flux contributions prec_i
and prec_j
compute the equilibrated states \(U_i^{\ast,j}\) and \(U_j^{\ast,i}\).
Definition at line 1066 of file hyperbolic_system.h.
|
inline |
Given a state U_i
and an index i
compute flux contributions.
Intended usage:
For the Shallow water equations we simply retrieve the bathymetry and return, both, state and bathymetry.
Definition at line 1082 of file hyperbolic_system.h.
|
inline |
Definition at line 1095 of file hyperbolic_system.h.
|
inline |
Given precomputed flux contributions prec_i
and prec_j
compute the equilibrated, low-order flux \((f(U_i^{\ast,j}) +
f(U_j^{\ast,i})\)
Definition at line 1108 of file hyperbolic_system.h.
References ryujin::add(), and ryujin::contract().
|
inline |
Given precomputed flux contributions prec_i
and prec_j
compute the high-order flux \((f(U_i^{\ast,j}) +
f(U_j^{\ast,i})\)
Definition at line 1140 of file hyperbolic_system.h.
References ryujin::add(), and ryujin::contract().
|
inline |
Given precomputed flux contributions prec_i
and prec_j
compute the equilibrated, low-order affine shift \( B_{ij} = -2d_ij(U^{\ast,j}_i)-2f((U^{\ast,j}_i))c_ij\).
Definition at line 1167 of file hyperbolic_system.h.
|
inline |
A regularized Gauckler-Manning friction coefficient.
FIXME: details
Definition at line 1187 of file hyperbolic_system.h.
|
inline |
Definition at line 1214 of file hyperbolic_system.h.
|
inline |
Definition at line 1229 of file hyperbolic_system.h.
state_type ryujin::ShallowWater::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::ShallowWater::HyperbolicSystemView< dim, Number >::from_initial_state | ( | const ST & | initial_state | ) | const |
Given an initial state [h, u_1, ..., u_d] return a conserved state [h, m_1, ..., m_d].
This function simply calls from_primitive_state() and expand_state().
|
inline |
Given a primitive state [h, u_1, ..., u_d] return a conserved state
Definition at line 1277 of file hyperbolic_system.h.
|
inline |
Given a conserved state return a primitive state [h, u_1, ..., u_d]
Definition at line 1293 of file hyperbolic_system.h.
state_type ryujin::ShallowWater::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 |
The dimension of the state space.
Definition at line 184 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the conserved state as a string.
Definition at line 201 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the primitive state as a string.
Definition at line 216 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputed initial values.
Definition at line 247 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the precomputed values.
Definition at line 265 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputed values.
Definition at line 271 of file hyperbolic_system.h.
|
inlinestatic |
An array holding all component names of the precomputed values.
Definition at line 287 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputation cycles.
Definition at line 293 of file hyperbolic_system.h.
|
staticconstexpr |
The low-order and high-order fluxes differ:
Definition at line 522 of file hyperbolic_system.h.
|
staticconstexpr |
We do have source terms
Definition at line 551 of file hyperbolic_system.h.