ryujin 2.1.1 revision 15c5f3ea3ef91eaf08db04f6f4c18a1823a6f822
Public Types | Public Member Functions | Static Public Attributes | List of all members
ryujin::ShallowWater::HyperbolicSystem Class Referencefinal

#include <source/shallow_water/hyperbolic_system.h>

Inheritance diagram for ryujin::ShallowWater::HyperbolicSystem:
Inheritance graph
[legend]
Collaboration diagram for ryujin::ShallowWater::HyperbolicSystem:
Collaboration graph
[legend]

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)
 

Detailed Description

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.

Member Typedef Documentation

◆ state_type

template<int dim, typename Number >
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.

◆ primitive_state_type

template<int dim, typename Number >
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.

◆ flux_type

template<int dim, typename Number >
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.

◆ flux_contribution_type

template<int dim, typename Number >
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.

◆ precomputed_initial_type

template<int dim, typename Number >
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.

◆ precomputed_type

template<int dim, typename Number >
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.

Constructor & Destructor Documentation

◆ HyperbolicSystem()

ryujin::ShallowWater::HyperbolicSystem::HyperbolicSystem ( const std::string &  subsection = "/HyperbolicSystem")

Constructor.

Definition at line 14 of file hyperbolic_system.template.h.

Member Function Documentation

◆ precomputation() [1/2]

template<unsigned int cycle, typename Number , typename ScalarNumber = typename get_value_type<Number>::type, int problem_dim, typename MCV , typename SPARSITY >
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.

◆ water_depth() [1/2]

template<int problem_dim, typename Number >
static Number ryujin::ShallowWater::HyperbolicSystem::water_depth ( const dealii::Tensor< 1, problem_dim, Number > &  U)
static

◆ water_depth_sharp() [1/2]

template<int problem_dim, typename Number >
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().

◆ inverse_water_depth_mollified() [1/2]

template<int problem_dim, typename Number >
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().

◆ inverse_water_depth_sharp() [1/2]

template<int problem_dim, typename Number >
Number ryujin::ShallowWater::HyperbolicSystem::inverse_water_depth_sharp ( const dealii::Tensor< 1, problem_dim, Number > &  U) const

◆ filter_dry_water_depth() [1/2]

template<typename Number >
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().

◆ momentum() [1/2]

template<int problem_dim, typename Number >
static dealii::Tensor< 1, problem_dim - 1, Number > ryujin::ShallowWater::HyperbolicSystem::momentum ( const dealii::Tensor< 1, problem_dim, Number > &  U)
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().

◆ kinetic_energy() [1/2]

template<int problem_dim, typename Number >
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().

◆ pressure() [1/2]

template<int problem_dim, typename Number >
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().

◆ speed_of_sound() [1/2]

template<int problem_dim, typename Number >
Number ryujin::ShallowWater::HyperbolicSystem::speed_of_sound ( const dealii::Tensor< 1, problem_dim, Number > &  U) const

◆ mathematical_entropy() [1/2]

template<int problem_dim, typename Number >
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().

◆ mathematical_entropy_derivative() [1/2]

template<int problem_dim, typename Number >
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.

◆ is_admissible() [1/2]

template<int problem_dim, typename Number >
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.

◆ prescribe_riemann_characteristic() [1/2]

template<int component, int problem_dim, typename Number >
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.

◆ apply_boundary_conditions() [1/2]

template<int problem_dim, typename Number , typename Lambda >
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:

  • Dirichlet boundary conditions by prescribing the return value of get_dirichlet_data() as is.
  • Slip boundary conditions where we remove the normal component of the momentum.
  • No slip boundary conditions where we set the momentum to 0.
  • "Dynamic boundary" conditions that prescribe different Riemann invariants from the return value of get_dirichlet_data() depending on the flow state (supersonic versus subsonic, outflow versus inflow).

◆ star_state() [1/2]

template<int problem_dim, typename Number >
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().

◆ f() [1/2]

template<typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type>
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() [1/4]

template<typename MultiComponentVector , typename MultiComponentVector2 , typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type>
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:

for (unsigned int i = n_internal; i < n_owned; ++i) {
// ...
const auto flux_i = flux_contribution(precomputed..., i, U_i);
for (unsigned int col_idx = 1; col_idx < row_length; ++col_idx) {
// ...
const auto flux_j = flux_contribution(precomputed..., js, U_j);
const auto flux_ij = flux(flux_i, flux_j);
}
}
flux_type< dim, T > flux(const std::tuple< ST, T > &prec_i, const std::tuple< ST, T > &prec_j) const
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

◆ flux_contribution() [2/4]

template<typename MultiComponentVector , typename MultiComponentVector2 , typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type>
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() [1/2]

template<typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type>
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})\)

◆ high_order_flux() [1/2]

template<typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type>
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\)

◆ equilibrated_states() [1/2]

template<typename ST , int dim = ST::dimension - 1, typename T = typename ST::value_type>
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().

◆ low_order_nodal_source()

template<typename MultiComponentVector , typename ST , int dim, typename T >
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.

◆ high_order_nodal_source()

template<typename MultiComponentVector , typename ST , int dim, typename T >
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.

◆ low_order_stencil_source()

template<typename ST , int dim, typename T >
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().

◆ high_order_stencil_source()

template<typename ST , int dim, typename T >
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().

◆ affine_shift_stencil_source()

template<typename ST , int dim, typename T >
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().

◆ expand_state() [1/2]

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 > ryujin::ShallowWater::HyperbolicSystem::expand_state ( const ST &  state) const

◆ from_primitive_state() [1/2]

template<int problem_dim, typename Number >
dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::from_primitive_state ( const dealii::Tensor< 1, problem_dim, Number > &  primitive_state) const

◆ to_primitive_state() [1/2]

template<int problem_dim, typename Number >
dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::to_primitive_state ( const dealii::Tensor< 1, problem_dim, Number > &  state) const

◆ apply_galilei_transform()

template<int problem_dim, typename Number , typename Lambda >
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().

◆ precomputation() [2/2]

template<unsigned int cycle, typename Number , typename ScalarNumber , int problem_dim, typename MCV , typename SPARSITY >
DEAL_II_ALWAYS_INLINE void ryujin::ShallowWater::HyperbolicSystem::precomputation ( MCV &  precomputed_values,
const MultiComponentVector< ScalarNumber, problem_dim > &  U,
const SPARSITY &  ,
unsigned int  i 
) const
inline

Definition at line 580 of file hyperbolic_system.h.

References mathematical_entropy().

◆ water_depth() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE Number ryujin::ShallowWater::HyperbolicSystem::water_depth ( const dealii::Tensor< 1, problem_dim, Number > &  U)
inline

Definition at line 598 of file hyperbolic_system.h.

◆ inverse_water_depth_mollified() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE Number ryujin::ShallowWater::HyperbolicSystem::inverse_water_depth_mollified ( const dealii::Tensor< 1, problem_dim, Number > &  U) const
inline

Definition at line 607 of file hyperbolic_system.h.

References ryujin::positive_part(), and water_depth().

◆ water_depth_sharp() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE Number ryujin::ShallowWater::HyperbolicSystem::water_depth_sharp ( const dealii::Tensor< 1, problem_dim, Number > &  U) const
inline

Definition at line 624 of file hyperbolic_system.h.

References water_depth().

◆ inverse_water_depth_sharp() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE Number ryujin::ShallowWater::HyperbolicSystem::inverse_water_depth_sharp ( const dealii::Tensor< 1, problem_dim, Number > &  U) const
inline

Definition at line 640 of file hyperbolic_system.h.

References water_depth_sharp().

◆ filter_dry_water_depth() [2/2]

template<typename Number >
DEAL_II_ALWAYS_INLINE Number ryujin::ShallowWater::HyperbolicSystem::filter_dry_water_depth ( const Number &  h) const
inline

Definition at line 650 of file hyperbolic_system.h.

◆ momentum() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim - 1, Number > ryujin::ShallowWater::HyperbolicSystem::momentum ( const dealii::Tensor< 1, problem_dim, Number > &  U)
inline

Definition at line 665 of file hyperbolic_system.h.

◆ kinetic_energy() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE Number ryujin::ShallowWater::HyperbolicSystem::kinetic_energy ( const dealii::Tensor< 1, problem_dim, Number > &  U) const
inline

Definition at line 678 of file hyperbolic_system.h.

References inverse_water_depth_sharp(), momentum(), and water_depth().

◆ pressure() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE Number ryujin::ShallowWater::HyperbolicSystem::pressure ( const dealii::Tensor< 1, problem_dim, Number > &  U) const
inline

Definition at line 693 of file hyperbolic_system.h.

◆ speed_of_sound() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE Number ryujin::ShallowWater::HyperbolicSystem::speed_of_sound ( const dealii::Tensor< 1, problem_dim, Number > &  U) const
inline

Definition at line 707 of file hyperbolic_system.h.

◆ mathematical_entropy() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE Number ryujin::ShallowWater::HyperbolicSystem::mathematical_entropy ( const dealii::Tensor< 1, problem_dim, Number > &  U) const
inline

Definition at line 718 of file hyperbolic_system.h.

References kinetic_energy(), and pressure().

◆ mathematical_entropy_derivative() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::mathematical_entropy_derivative ( const dealii::Tensor< 1, problem_dim, Number > &  U) const
inline

Definition at line 729 of file hyperbolic_system.h.

References inverse_water_depth_sharp(), and momentum().

◆ is_admissible() [2/2]

template<int problem_dim, typename T >
DEAL_II_ALWAYS_INLINE bool ryujin::ShallowWater::HyperbolicSystem::is_admissible ( const dealii::Tensor< 1, problem_dim, T > &  U) const
inline

Definition at line 765 of file hyperbolic_system.h.

References filter_dry_water_depth(), and water_depth().

◆ prescribe_riemann_characteristic() [2/2]

template<int component, int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE 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
inline

Definition at line 788 of file hyperbolic_system.h.

References inverse_water_depth_sharp(), momentum(), and speed_of_sound().

◆ apply_boundary_conditions() [2/2]

template<int problem_dim, typename Number , typename Lambda >
DEAL_II_ALWAYS_INLINE 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
inline

◆ star_state() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::star_state ( const dealii::Tensor< 1, problem_dim, Number > &  U_left,
const Number &  Z_left,
const Number &  Z_right 
) const
inline

Definition at line 907 of file hyperbolic_system.h.

References inverse_water_depth_mollified(), and water_depth().

◆ f() [2/2]

template<typename ST , int dim, typename T >
DEAL_II_ALWAYS_INLINE auto ryujin::ShallowWater::HyperbolicSystem::f ( const ST &  U) const -> flux_type<dim, T>
inline

Definition at line 921 of file hyperbolic_system.h.

◆ flux_contribution() [3/4]

template<typename MCV , typename MICV , typename ST , int dim, typename T >
DEAL_II_ALWAYS_INLINE auto ryujin::ShallowWater::HyperbolicSystem::flux_contribution ( const MCV &  ,
const MICV &  precomputed_initial_values,
const unsigned int  i,
const ST &  U_i 
) const -> flux_contribution_type<dim, T>
inline

Definition at line 941 of file hyperbolic_system.h.

◆ flux_contribution() [4/4]

template<typename MCV , typename MICV , typename ST , int dim, typename T >
DEAL_II_ALWAYS_INLINE auto ryujin::ShallowWater::HyperbolicSystem::flux_contribution ( const MCV &  ,
const MICV &  precomputed_initial_values,
const unsigned int *  js,
const ST &  U_j 
) const -> flux_contribution_type<dim, T>
inline

Definition at line 954 of file hyperbolic_system.h.

◆ flux() [2/2]

template<typename ST , int dim, typename T >
DEAL_II_ALWAYS_INLINE auto ryujin::ShallowWater::HyperbolicSystem::flux ( const std::tuple< ST, T > &  prec_i,
const std::tuple< ST, T > &  prec_j 
) const -> flux_type<dim, T>
inline

Definition at line 967 of file hyperbolic_system.h.

References ryujin::add().

◆ high_order_flux() [2/2]

template<typename ST , int dim, typename T >
DEAL_II_ALWAYS_INLINE auto ryujin::ShallowWater::HyperbolicSystem::high_order_flux ( const std::tuple< ST, T > &  prec_i,
const std::tuple< ST, T > &  prec_j 
) const -> flux_type<dim, T>
inline

Definition at line 982 of file hyperbolic_system.h.

References ryujin::add().

◆ equilibrated_states() [2/2]

template<typename ST , int dim, typename T >
DEAL_II_ALWAYS_INLINE std::array< ST, 2 > ryujin::ShallowWater::HyperbolicSystem::equilibrated_states ( const std::tuple< ST, T > &  prec_i,
const std::tuple< ST, T > &  prec_j 
) const
inline

Definition at line 998 of file hyperbolic_system.h.

References star_state().

◆ expand_state() [2/2]

template<int dim, typename ST , typename T , int dim2, typename >
auto ryujin::ShallowWater::HyperbolicSystem::expand_state ( const ST &  state) const -> state_type<dim, T>

Definition at line 1098 of file hyperbolic_system.h.

◆ from_primitive_state() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::from_primitive_state ( const dealii::Tensor< 1, problem_dim, Number > &  primitive_state) const
inline

Definition at line 1112 of file hyperbolic_system.h.

◆ to_primitive_state() [2/2]

template<int problem_dim, typename Number >
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim, Number > ryujin::ShallowWater::HyperbolicSystem::to_primitive_state ( const dealii::Tensor< 1, problem_dim, Number > &  state) const
inline

Definition at line 1128 of file hyperbolic_system.h.

References inverse_water_depth_sharp().

Member Data Documentation

◆ problem_name

const std::string ryujin::ShallowWater::HyperbolicSystem::problem_name
static

The name of the hyperbolic system.

Definition at line 41 of file hyperbolic_system.h.

◆ problem_dimension

template<int dim>
constexpr unsigned int ryujin::ShallowWater::HyperbolicSystem::problem_dimension = 1 + dim
staticconstexpr

The dimension of the state space.

Definition at line 47 of file hyperbolic_system.h.

◆ component_names

template<int dim>
const std::array<std::string, problem_dimension<dim> > ryujin::ShallowWater::HyperbolicSystem::component_names
static

An array holding all component names of the conserved state.

Definition at line 61 of file hyperbolic_system.h.

◆ primitive_component_names

template<int dim>
const std::array<std::string, problem_dimension<dim> > ryujin::ShallowWater::HyperbolicSystem::primitive_component_names
static

An array holding all component names of the primitive state.

Definition at line 75 of file hyperbolic_system.h.

◆ n_precomputed_initial_values

template<int dim>
constexpr unsigned int ryujin::ShallowWater::HyperbolicSystem::n_precomputed_initial_values = 1
staticconstexpr

The number of precomputed initial values.

Definition at line 105 of file hyperbolic_system.h.

◆ precomputed_initial_names

template<int dim>
const std::array<std::string, n_precomputed_initial_values<dim> > ryujin::ShallowWater::HyperbolicSystem::precomputed_initial_names
static

An array holding all component names of the precomputed values.

Definition at line 119 of file hyperbolic_system.h.

◆ n_precomputed_values

template<int dim>
constexpr unsigned int ryujin::ShallowWater::HyperbolicSystem::n_precomputed_values = 1
staticconstexpr

The number of precomputed values.

Definition at line 125 of file hyperbolic_system.h.

◆ precomputed_names

template<int dim>
const std::array<std::string, n_precomputed_values<dim> > ryujin::ShallowWater::HyperbolicSystem::precomputed_names
static

An array holding all component names of the precomputed values.

Definition at line 138 of file hyperbolic_system.h.

◆ n_precomputation_cycles

constexpr unsigned int ryujin::ShallowWater::HyperbolicSystem::n_precomputation_cycles = 1
staticconstexpr

The number of precomputation cycles.

Definition at line 143 of file hyperbolic_system.h.

◆ have_high_order_flux

constexpr bool ryujin::ShallowWater::HyperbolicSystem::have_high_order_flux = true
staticconstexpr

Definition at line 394 of file hyperbolic_system.h.

◆ have_equilibrated_states

constexpr bool ryujin::ShallowWater::HyperbolicSystem::have_equilibrated_states = true
staticconstexpr

Definition at line 406 of file hyperbolic_system.h.

◆ have_source_terms

constexpr bool ryujin::ShallowWater::HyperbolicSystem::have_source_terms = true
staticconstexpr

Definition at line 426 of file hyperbolic_system.h.


The documentation for this class was generated from the following files: