![]() |
ryujin 2.1.1 revision 15c5f3ea3ef91eaf08db04f6f4c18a1823a6f822
|
#include <source/shallow_water/hyperbolic_system.h>
Public Types | |
template<int dim, typename Number > | |
using | state_type = dealii::Tensor< 1, problem_dimension< dim >, Number > |
template<int dim, typename Number > | |
using | primitive_state_type = dealii::Tensor< 1, problem_dimension< dim >, Number > |
template<int dim, typename Number > | |
using | flux_type = dealii::Tensor< 1, problem_dimension< dim >, dealii::Tensor< 1, dim, Number > > |
template<int dim, typename Number > | |
using | flux_contribution_type = std::tuple< state_type< dim, Number >, Number > |
Public Member Functions | |
HyperbolicSystem (const std::string &subsection="/HyperbolicSystem") | |
template<unsigned int cycle, typename Number , typename ScalarNumber , int problem_dim, typename MCV , typename SPARSITY > | |
DEAL_II_ALWAYS_INLINE void | precomputation (MCV &precomputed_values, const MultiComponentVector< ScalarNumber, problem_dim > &U, const SPARSITY &, unsigned int i) const |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE Number | water_depth (const dealii::Tensor< 1, problem_dim, Number > &U) |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE Number | inverse_water_depth_mollified (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE Number | water_depth_sharp (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE Number | inverse_water_depth_sharp (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<typename Number > | |
DEAL_II_ALWAYS_INLINE Number | filter_dry_water_depth (const Number &h) const |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim - 1, Number > | momentum (const dealii::Tensor< 1, problem_dim, Number > &U) |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE Number | kinetic_energy (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE Number | pressure (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE Number | speed_of_sound (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE Number | mathematical_entropy (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim, Number > | mathematical_entropy_derivative (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename T > | |
DEAL_II_ALWAYS_INLINE bool | is_admissible (const dealii::Tensor< 1, problem_dim, T > &U) const |
template<int component, int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim, Number > | prescribe_riemann_characteristic (const dealii::Tensor< 1, problem_dim, Number > &U, const dealii::Tensor< 1, problem_dim, Number > &U_bar, const dealii::Tensor< 1, problem_dim - 1, Number > &normal) const |
template<int problem_dim, typename Number , typename Lambda > | |
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim, Number > | apply_boundary_conditions (dealii::types::boundary_id id, dealii::Tensor< 1, problem_dim, Number > U, const dealii::Tensor< 1, problem_dim - 1, Number > &normal, Lambda get_dirichlet_data) const |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim, Number > | star_state (const dealii::Tensor< 1, problem_dim, Number > &U_left, const Number &Z_left, const Number &Z_right) const |
template<typename ST , int dim, typename T > | |
DEAL_II_ALWAYS_INLINE auto | f (const ST &U) const -> flux_type< dim, T > |
template<typename MCV , typename MICV , typename ST , int dim, typename T > | |
DEAL_II_ALWAYS_INLINE auto | flux_contribution (const MCV &, const MICV &precomputed_initial_values, const unsigned int i, const ST &U_i) const -> flux_contribution_type< dim, T > |
template<typename MCV , typename MICV , typename ST , int dim, typename T > | |
DEAL_II_ALWAYS_INLINE auto | flux_contribution (const MCV &, const MICV &precomputed_initial_values, const unsigned int *js, const ST &U_j) const -> flux_contribution_type< dim, T > |
template<typename ST , int dim, typename T > | |
DEAL_II_ALWAYS_INLINE auto | flux (const std::tuple< ST, T > &prec_i, const std::tuple< ST, T > &prec_j) const -> flux_type< dim, T > |
template<typename ST , int dim, typename T > | |
DEAL_II_ALWAYS_INLINE auto | high_order_flux (const std::tuple< ST, T > &prec_i, const std::tuple< ST, T > &prec_j) const -> flux_type< dim, T > |
template<typename ST , int dim, typename T > | |
DEAL_II_ALWAYS_INLINE std::array< ST, 2 > | equilibrated_states (const std::tuple< ST, T > &prec_i, const std::tuple< ST, T > &prec_j) const |
template<int dim, typename ST , typename T , int dim2, typename > | |
auto | expand_state (const ST &state) const -> state_type< dim, T > |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim, Number > | from_primitive_state (const dealii::Tensor< 1, problem_dim, Number > &primitive_state) const |
template<int problem_dim, typename Number > | |
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim, Number > | to_primitive_state (const dealii::Tensor< 1, problem_dim, Number > &state) const |
Special functions for boundary states | |
template<int component, int problem_dim, typename Number > | |
dealii::Tensor< 1, problem_dim, Number > | prescribe_riemann_characteristic (const dealii::Tensor< 1, problem_dim, Number > &U, const dealii::Tensor< 1, problem_dim, Number > &U_bar, const dealii::Tensor< 1, problem_dim - 1, Number > &normal) const |
template<int problem_dim, typename Number , typename Lambda > | |
dealii::Tensor< 1, problem_dim, Number > | apply_boundary_conditions (dealii::types::boundary_id id, dealii::Tensor< 1, problem_dim, Number > U, const dealii::Tensor< 1, problem_dim - 1, Number > &normal, Lambda get_dirichlet_data) const |
State transformations (primitive states, expanding | |
dimensionality, Galilei transform, etc.) | |
template<int dim, typename ST , typename T = typename ST::value_type, int dim2 = ST::dimension - 1, typename = typename std::enable_if<(dim >= dim2)>::type> | |
state_type< dim, T > | expand_state (const ST &state) const |
template<int problem_dim, typename Number > | |
dealii::Tensor< 1, problem_dim, Number > | from_primitive_state (const dealii::Tensor< 1, problem_dim, Number > &primitive_state) const |
template<int problem_dim, typename Number > | |
dealii::Tensor< 1, problem_dim, Number > | to_primitive_state (const dealii::Tensor< 1, problem_dim, Number > &state) const |
template<int problem_dim, typename Number , typename Lambda > | |
dealii::Tensor< 1, problem_dim, Number > | apply_galilei_transform (const dealii::Tensor< 1, problem_dim, Number > &state, const Lambda &lambda) const |
Static Public Attributes | |
static const std::string | problem_name |
template<int dim> | |
static constexpr unsigned int | problem_dimension = 1 + dim |
template<int dim> | |
static const std::array< std::string, problem_dimension< dim > > | component_names |
template<int dim> | |
static const std::array< std::string, problem_dimension< dim > > | primitive_component_names |
Precomputation of flux quantities | |
template<int dim, typename Number > | |
using | precomputed_initial_type = std::array< Number, n_precomputed_initial_values< dim > > |
template<int dim, typename Number > | |
using | precomputed_type = std::array< Number, n_precomputed_values< dim > > |
template<int dim> | |
static constexpr unsigned int | n_precomputed_initial_values = 1 |
template<int dim> | |
static const std::array< std::string, n_precomputed_initial_values< dim > > | precomputed_initial_names |
template<int dim> | |
static constexpr unsigned int | n_precomputed_values = 1 |
template<int dim> | |
static const std::array< std::string, n_precomputed_values< dim > > | precomputed_names |
static constexpr unsigned int | n_precomputation_cycles = 1 |
template<unsigned int cycle, typename Number , typename ScalarNumber = typename get_value_type<Number>::type, int problem_dim, typename MCV , typename SPARSITY > | |
void | precomputation (MCV &precomputed_values, const MultiComponentVector< ScalarNumber, problem_dim > &U, const SPARSITY &sparsity_simd, unsigned int i) const |
Computing fluxes. | |
static constexpr bool | have_high_order_flux = true |
static constexpr bool | have_equilibrated_states = true |
template<int problem_dim, typename Number > | |
dealii::Tensor< 1, problem_dim, Number > | star_state (const dealii::Tensor< 1, problem_dim, Number > &U, const Number &Z_left, const Number &Z_right) const |
template<typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type> | |
flux_type< dim, T > | f (const ST &U) const |
template<typename MultiComponentVector , typename MultiComponentVector2 , typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type> | |
flux_contribution_type< dim, T > | flux_contribution (const MultiComponentVector &precomputed_values, const MultiComponentVector2 &precomputed_initial_values, const unsigned int i, const ST &U_i) const |
template<typename MultiComponentVector , typename MultiComponentVector2 , typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type> | |
flux_contribution_type< dim, T > | flux_contribution (const MultiComponentVector &precomputed_values, const MultiComponentVector2 &precomputed_initial_values, const unsigned int *js, const ST &U_j) const |
template<typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type> | |
flux_type< dim, T > | flux (const std::tuple< ST, T > &prec_i, const std::tuple< ST, T > &prec_j) const |
template<typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type> | |
flux_type< dim, T > | high_order_flux (const std::tuple< ST, T > &prec_i, const std::tuple< ST, T > &prec_j) const |
template<typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type> | |
std::array< ST, 2 > | equilibrated_states (const std::tuple< ST, T > &prec_i, const std::tuple< ST, T > &prec_j) const |
Computing stencil source terms | |
static constexpr bool | have_source_terms = true |
template<typename MultiComponentVector , typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type> | |
ST | low_order_nodal_source (const MultiComponentVector &precomputed_values, const unsigned int i, const ST &U_i) const |
template<typename MultiComponentVector , typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type> | |
ST | high_order_nodal_source (const MultiComponentVector &precomputed_values, const unsigned int i, const ST &U_i) const |
template<typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type> | |
ST | low_order_stencil_source (const std::tuple< ST, T > &prec_i, const std::tuple< ST, T > &prec_j, const T &d_ij, const dealii::Tensor< 1, dim, T > &c_ij) const |
template<typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type> | |
ST | high_order_stencil_source (const std::tuple< ST, T > &prec_i, const std::tuple< ST, T > &prec_j, const T &d_ij, const dealii::Tensor< 1, dim, T > &c_ij) const |
template<typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type> | |
ST | affine_shift_stencil_source (const std::tuple< ST, T > &prec_i, const std::tuple< ST, T > &prec_j, const T &d_ij, const dealii::Tensor< 1, dim, T > &c_ij) const |
Computing derived physical quantities. | |
template<int problem_dim, typename Number > | |
Number | water_depth_sharp (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
Number | inverse_water_depth_mollified (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
Number | inverse_water_depth_sharp (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<typename Number > | |
Number | filter_dry_water_depth (const Number &h) const |
template<int problem_dim, typename Number > | |
Number | kinetic_energy (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
Number | pressure (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
Number | speed_of_sound (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
Number | mathematical_entropy (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
dealii::Tensor< 1, problem_dim, Number > | mathematical_entropy_derivative (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
bool | is_admissible (const dealii::Tensor< 1, problem_dim, Number > &U) const |
template<int problem_dim, typename Number > | |
static Number | water_depth (const dealii::Tensor< 1, problem_dim, Number > &U) |
template<int problem_dim, typename Number > | |
static dealii::Tensor< 1, problem_dim - 1, Number > | momentum (const dealii::Tensor< 1, problem_dim, Number > &U) |
Description of a dim
dimensional hyperbolic conservation law modeling the shallow water equations.
We have a (1 + dim) dimensional state space \([h, \textbf m]\), where \(h\) denotes the water depth, abd \(\textbf m\) is the momentum.
Definition at line 35 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystem::state_type = dealii::Tensor<1, problem_dimension<dim>, Number> |
The storage type used for a (conserved) state vector \(\boldsymbol U\).
Definition at line 54 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystem::primitive_state_type = dealii::Tensor<1, problem_dimension<dim>, Number> |
The storage type used for a primitive state vector.
Definition at line 67 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystem::flux_type = dealii:: Tensor<1, problem_dimension<dim>, dealii::Tensor<1, dim, Number> > |
The storage type used for the flux \(\mathbf{f}\).
Definition at line 81 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystem::flux_contribution_type = std::tuple<state_type<dim, Number>, Number> |
The storage type used for the flux precomputations.
Definition at line 88 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystem::precomputed_initial_type = std::array<Number, n_precomputed_initial_values<dim> > |
Array type used for precomputed initial values.
Definition at line 111 of file hyperbolic_system.h.
using ryujin::ShallowWater::HyperbolicSystem::precomputed_type = std::array<Number, n_precomputed_values<dim> > |
Array type used for precomputed values.
Definition at line 131 of file hyperbolic_system.h.
ryujin::ShallowWater::HyperbolicSystem::HyperbolicSystem | ( | const std::string & | subsection = "/HyperbolicSystem" | ) |
Constructor.
Definition at line 14 of file hyperbolic_system.template.h.
void ryujin::ShallowWater::HyperbolicSystem::precomputation | ( | MCV & | precomputed_values, |
const MultiComponentVector< ScalarNumber, problem_dim > & | U, | ||
const SPARSITY & | sparsity_simd, | ||
unsigned int | i | ||
) | const |
Precomputed values for a given state.
|
static |
For a given (1+dim dimensional) state vector U
, return the water depth U[0]
Referenced by ryujin::ShallowWater::Limiter< dim, Number >::accumulate(), high_order_stencil_source(), inverse_water_depth_mollified(), is_admissible(), kinetic_energy(), ryujin::ShallowWater::Limiter< dim, Number >::limit(), low_order_stencil_source(), star_state(), and water_depth_sharp().
Number ryujin::ShallowWater::HyperbolicSystem::water_depth_sharp | ( | const dealii::Tensor< 1, problem_dim, Number > & | U | ) | const |
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.
Referenced by inverse_water_depth_sharp().
Number ryujin::ShallowWater::HyperbolicSystem::inverse_water_depth_mollified | ( | const dealii::Tensor< 1, problem_dim, Number > & | U | ) | const |
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.
Referenced by star_state().
Number ryujin::ShallowWater::HyperbolicSystem::inverse_water_depth_sharp | ( | const dealii::Tensor< 1, problem_dim, Number > & | U | ) | const |
Number ryujin::ShallowWater::HyperbolicSystem::filter_dry_water_depth | ( | const Number & | h | ) | const |
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.
Referenced by is_admissible(), and ryujin::ShallowWater::Limiter< dim, Number >::limit().
|
static |
For a given (1+dim dimensional) state vector U
, return the momentum vector [U[1], ..., U[1+dim]]
.
Referenced by apply_boundary_conditions(), apply_galilei_transform(), kinetic_energy(), ryujin::ShallowWater::Limiter< dim, Number >::limit(), mathematical_entropy_derivative(), and prescribe_riemann_characteristic().
Number ryujin::ShallowWater::HyperbolicSystem::kinetic_energy | ( | const dealii::Tensor< 1, problem_dim, Number > & | U | ) | const |
For a given (1+dim dimensional) state vector U
, compute and return the kinetic energy.
\[ KE = 1/2 |m|^2 / h \]
Referenced by ryujin::ShallowWater::Limiter< dim, Number >::accumulate(), and mathematical_entropy().
Number ryujin::ShallowWater::HyperbolicSystem::pressure | ( | const dealii::Tensor< 1, problem_dim, Number > & | U | ) | const |
For a given (state dimensional) state vector U
, compute and return the hydrostatic pressure \(p\):
\[ p = 1/2 g h^2 \]
Referenced by mathematical_entropy().
Number ryujin::ShallowWater::HyperbolicSystem::speed_of_sound | ( | const dealii::Tensor< 1, problem_dim, Number > & | U | ) | const |
For a given (1+dim dimensional) state vector U
, compute the (physical) speed of sound:
\[ c^2 = g * h \]
Referenced by apply_boundary_conditions(), ryujin::ShallowWater::SlopingRampDamBreak< dim, Number, state_type >::compute(), ryujin::ShallowWater::ThreeBumpsDamBreak< dim, Number, state_type >::compute(), and prescribe_riemann_characteristic().
Number ryujin::ShallowWater::HyperbolicSystem::mathematical_entropy | ( | const dealii::Tensor< 1, problem_dim, Number > & | U | ) | const |
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\).
Referenced by precomputation().
dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::mathematical_entropy_derivative | ( | const dealii::Tensor< 1, problem_dim, Number > & | U | ) | const |
For a given (1+dim dimensional) state vector U
, compute and return the derivative \(\eta'\) of the entropy defined above.
bool ryujin::ShallowWater::HyperbolicSystem::is_admissible | ( | const dealii::Tensor< 1, problem_dim, Number > & | U | ) | const |
Returns whether the state U is admissible. If U is a vectorized state then U is admissible if all vectorized values are admissible.
dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::prescribe_riemann_characteristic | ( | const dealii::Tensor< 1, problem_dim, Number > & | U, |
const dealii::Tensor< 1, problem_dim, Number > & | U_bar, | ||
const dealii::Tensor< 1, problem_dim - 1, 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.
dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::apply_boundary_conditions | ( | dealii::types::boundary_id | id, |
dealii::Tensor< 1, problem_dim, Number > | U, | ||
const dealii::Tensor< 1, problem_dim - 1, Number > & | normal, | ||
Lambda | get_dirichlet_data | ||
) | const |
Apply boundary conditions.
For the compressible Euler equations we have:
dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::star_state | ( | const dealii::Tensor< 1, problem_dim, Number > & | U, |
const Number & | Z_left, | ||
const Number & | Z_right | ||
) | const |
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
Referenced by ryujin::ShallowWater::Limiter< dim, Number >::accumulate(), equilibrated_states(), and low_order_stencil_source().
flux_type< dim, T > ryujin::ShallowWater::HyperbolicSystem::f | ( | const ST & | U | ) | const |
Given a state U
compute the flux
\[ \begin{pmatrix} \textbf m \\ \textbf v\otimes \textbf m + p\mathbb{I}_d \\ \end{pmatrix}, \]
Referenced by ryujin::ShallowWater::Limiter< dim, Number >::accumulate(), and affine_shift_stencil_source().
flux_contribution_type< dim, T > ryujin::ShallowWater::HyperbolicSystem::flux_contribution | ( | const MultiComponentVector & | precomputed_values, |
const MultiComponentVector2 & | precomputed_initial_values, | ||
const unsigned int | i, | ||
const ST & | U_i | ||
) | const |
Given a state U_i
and an index i
compute flux contributions.
Intended usage:
flux_contribution_type< dim, T > ryujin::ShallowWater::HyperbolicSystem::flux_contribution | ( | const MultiComponentVector & | precomputed_values, |
const MultiComponentVector2 & | precomputed_initial_values, | ||
const unsigned int * | js, | ||
const ST & | U_j | ||
) | const |
flux_type< dim, T > ryujin::ShallowWater::HyperbolicSystem::flux | ( | const std::tuple< ST, T > & | prec_i, |
const std::tuple< ST, T > & | prec_j | ||
) | const |
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})\)
flux_type< dim, T > ryujin::ShallowWater::HyperbolicSystem::high_order_flux | ( | const std::tuple< ST, T > & | prec_i, |
const std::tuple< ST, T > & | prec_j | ||
) | const |
Given precomputed flux contributions prec_i
and prec_j
compute the high-order flux \((f(U_i}) + f(U_j\)
std::array< ST, 2 > ryujin::ShallowWater::HyperbolicSystem::equilibrated_states | ( | const std::tuple< ST, T > & | prec_i, |
const std::tuple< ST, T > & | prec_j | ||
) | const |
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})\)
Referenced by affine_shift_stencil_source().
ST ryujin::ShallowWater::HyperbolicSystem::low_order_nodal_source | ( | const MultiComponentVector & | precomputed_values, |
const unsigned int | i, | ||
const ST & | U_i | ||
) | const |
FIXME
Definition at line 1010 of file hyperbolic_system.h.
ST ryujin::ShallowWater::HyperbolicSystem::high_order_nodal_source | ( | const MultiComponentVector & | precomputed_values, |
const unsigned int | i, | ||
const ST & | U_i | ||
) | const |
FIXME
Definition at line 1021 of file hyperbolic_system.h.
ST ryujin::ShallowWater::HyperbolicSystem::low_order_stencil_source | ( | const std::tuple< ST, T > & | prec_i, |
const std::tuple< ST, T > & | prec_j, | ||
const T & | d_ij, | ||
const dealii::Tensor< 1, dim, T > & | c_ij | ||
) | const |
Given precomputed flux contributions prec_i
and prec_j
compute the equilibrated, low-order source term \(-g(H^{\ast,j}_i)^2c_ij\).
Definition at line 1032 of file hyperbolic_system.h.
References star_state(), and water_depth().
ST ryujin::ShallowWater::HyperbolicSystem::high_order_stencil_source | ( | const std::tuple< ST, T > & | prec_i, |
const std::tuple< ST, T > & | prec_j, | ||
const T & | d_ij, | ||
const dealii::Tensor< 1, dim, T > & | c_ij | ||
) | const |
Given precomputed flux contributions prec_i
and prec_j
compute the high-order source term \( g H_i Z_j c_ij\).
Definition at line 1056 of file hyperbolic_system.h.
References water_depth().
ST ryujin::ShallowWater::HyperbolicSystem::affine_shift_stencil_source | ( | const std::tuple< ST, T > & | prec_i, |
const std::tuple< ST, T > & | prec_j, | ||
const T & | d_ij, | ||
const dealii::Tensor< 1, dim, T > & | c_ij | ||
) | const |
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 1081 of file hyperbolic_system.h.
References ryujin::contract(), equilibrated_states(), and f().
state_type< dim, T > ryujin::ShallowWater::HyperbolicSystem::expand_state | ( | const ST & | state | ) | const |
dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::from_primitive_state | ( | const dealii::Tensor< 1, problem_dim, Number > & | primitive_state | ) | const |
dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::to_primitive_state | ( | const dealii::Tensor< 1, problem_dim, Number > & | state | ) | const |
dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::apply_galilei_transform | ( | const dealii::Tensor< 1, problem_dim, Number > & | state, |
const Lambda & | lambda | ||
) | const |
Definition at line 1144 of file hyperbolic_system.h.
References momentum().
|
inline |
Definition at line 580 of file hyperbolic_system.h.
References mathematical_entropy().
|
inline |
Definition at line 598 of file hyperbolic_system.h.
|
inline |
Definition at line 607 of file hyperbolic_system.h.
References ryujin::positive_part(), and water_depth().
|
inline |
Definition at line 624 of file hyperbolic_system.h.
References water_depth().
|
inline |
Definition at line 640 of file hyperbolic_system.h.
References water_depth_sharp().
|
inline |
Definition at line 650 of file hyperbolic_system.h.
|
inline |
Definition at line 665 of file hyperbolic_system.h.
|
inline |
Definition at line 678 of file hyperbolic_system.h.
References inverse_water_depth_sharp(), momentum(), and water_depth().
|
inline |
Definition at line 693 of file hyperbolic_system.h.
|
inline |
Definition at line 707 of file hyperbolic_system.h.
|
inline |
Definition at line 718 of file hyperbolic_system.h.
References kinetic_energy(), and pressure().
|
inline |
Definition at line 729 of file hyperbolic_system.h.
References inverse_water_depth_sharp(), and momentum().
|
inline |
Definition at line 765 of file hyperbolic_system.h.
References filter_dry_water_depth(), and water_depth().
|
inline |
Definition at line 788 of file hyperbolic_system.h.
References inverse_water_depth_sharp(), momentum(), and speed_of_sound().
|
inline |
Definition at line 842 of file hyperbolic_system.h.
References ryujin::dirichlet, ryujin::dynamic, inverse_water_depth_sharp(), momentum(), ryujin::no_slip, ryujin::slip, and speed_of_sound().
|
inline |
Definition at line 907 of file hyperbolic_system.h.
References inverse_water_depth_mollified(), and water_depth().
|
inline |
Definition at line 921 of file hyperbolic_system.h.
|
inline |
Definition at line 941 of file hyperbolic_system.h.
|
inline |
Definition at line 954 of file hyperbolic_system.h.
|
inline |
Definition at line 967 of file hyperbolic_system.h.
References ryujin::add().
|
inline |
Definition at line 982 of file hyperbolic_system.h.
References ryujin::add().
|
inline |
Definition at line 998 of file hyperbolic_system.h.
References star_state().
auto ryujin::ShallowWater::HyperbolicSystem::expand_state | ( | const ST & | state | ) | const -> state_type<dim, T> |
Definition at line 1098 of file hyperbolic_system.h.
|
inline |
Definition at line 1112 of file hyperbolic_system.h.
|
inline |
Definition at line 1128 of file hyperbolic_system.h.
References inverse_water_depth_sharp().
|
static |
The name of the hyperbolic system.
Definition at line 41 of file hyperbolic_system.h.
|
staticconstexpr |
The dimension of the state space.
Definition at line 47 of file hyperbolic_system.h.
|
static |
An array holding all component names of the conserved state.
Definition at line 61 of file hyperbolic_system.h.
|
static |
An array holding all component names of the primitive state.
Definition at line 75 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputed initial values.
Definition at line 105 of file hyperbolic_system.h.
|
static |
An array holding all component names of the precomputed values.
Definition at line 119 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputed values.
Definition at line 125 of file hyperbolic_system.h.
|
static |
An array holding all component names of the precomputed values.
Definition at line 138 of file hyperbolic_system.h.
|
staticconstexpr |
The number of precomputation cycles.
Definition at line 143 of file hyperbolic_system.h.
|
staticconstexpr |
Definition at line 394 of file hyperbolic_system.h.
|
staticconstexpr |
Definition at line 406 of file hyperbolic_system.h.
|
staticconstexpr |
Definition at line 426 of file hyperbolic_system.h.