![]() |
ryujin 2.1.1 revision 15c5f3ea3ef91eaf08db04f6f4c18a1823a6f822
|
#include <source/shallow_water/riemann_solver.h>
Public Types | |
using | primitive_type = std::array< Number, riemann_data_size > |
using | state_type = HyperbolicSystem::state_type< dim, Number > |
using | precomputed_type = HyperbolicSystem::precomputed_type< dim, Number > |
using | ScalarNumber = typename get_value_type< Number >::type |
Public Member Functions | |
Compute wavespeed estimates | |
RiemannSolver (const HyperbolicSystem &hyperbolic_system, const MultiComponentVector< ScalarNumber, n_precomputed_values > &precomputed_values) | |
Number | compute (const primitive_type &riemann_data_i, const primitive_type &riemann_data_j) const |
Number | compute (const state_type &U_i, const state_type &U_j, const unsigned int i, const unsigned int *js, const dealii::Tensor< 1, dim, Number > &n_ij) const |
Static Public Attributes | |
static constexpr unsigned int | problem_dimension |
static constexpr unsigned int | riemann_data_size = 3 |
static constexpr unsigned int | n_precomputed_values |
Internal functions used in the Riemann solver | |
Number | f (const primitive_type &primitive_state, const Number &h_star) const |
Number | phi (const primitive_type &riemann_data_i, const primitive_type &riemann_data_j, const Number &h) const |
Number | lambda1_minus (const primitive_type &riemann_data, const Number h_star) const |
Number | lambda3_plus (const primitive_type &riemann_data, const Number h_star) const |
Number | compute_lambda (const primitive_type &riemann_data_i, const primitive_type &riemann_data_j, const Number h_star) const |
Number | h_star_two_rarefaction (const primitive_type &riemann_data_i, const primitive_type &riemann_data_j) const |
primitive_type | riemann_data_from_state (const state_type &U, const dealii::Tensor< 1, dim, Number > &n_ij) const |
A fast approximative solver for the 1D Riemann problem. The solver ensures that the estimate \(\lambda_{\text{max}}\) that is returned for the maximal wavespeed is a strict upper bound.
The solver is based on [4].
Definition at line 31 of file riemann_solver.h.
using ryujin::ShallowWater::RiemannSolver< dim, Number >::primitive_type = std::array<Number, riemann_data_size> |
Definition at line 41 of file riemann_solver.h.
using ryujin::ShallowWater::RiemannSolver< dim, Number >::state_type = HyperbolicSystem::state_type<dim, Number> |
The storage type used for a (conserved) state vector \(\boldsymbol U\).
Definition at line 46 of file riemann_solver.h.
using ryujin::ShallowWater::RiemannSolver< dim, Number >::precomputed_type = HyperbolicSystem::precomputed_type<dim, Number> |
Array type used for precomputed values.
Definition at line 51 of file riemann_solver.h.
using ryujin::ShallowWater::RiemannSolver< dim, Number >::ScalarNumber = typename get_value_type<Number>::type |
Definition at line 62 of file riemann_solver.h.
|
inline |
Constructor taking a HyperbolicSystem instance as argument
Definition at line 72 of file riemann_solver.h.
|
inline |
For two given 1D primitive states riemann_data_i and riemann_data_j, compute an estimation of an upper bound for the maximum wavespeed lambda.
Definition at line 211 of file riemann_solver.template.h.
Number ryujin::ShallowWater::RiemannSolver< dim, Number >::compute | ( | const state_type & | U_i, |
const state_type & | U_j, | ||
const unsigned int | i, | ||
const unsigned int * | js, | ||
const dealii::Tensor< 1, dim, Number > & | n_ij | ||
) | const |
For two given states U_i a U_j and a (normalized) "direction" n_ij compute an estimation of an upper bound for lambda.
Returns a tuple consisting of lambda max and the number of Newton iterations used in the solver to find it.
Definition at line 162 of file riemann_solver.h.
|
inline |
Definition at line 23 of file riemann_solver.template.h.
References ryujin::pow().
|
inline |
Definition at line 49 of file riemann_solver.template.h.
|
inline |
Definition at line 64 of file riemann_solver.template.h.
References ryujin::positive_part().
|
inline |
Definition at line 81 of file riemann_solver.template.h.
References ryujin::positive_part().
|
inline |
Definition at line 98 of file riemann_solver.template.h.
References ryujin::negative_part(), and ryujin::positive_part().
|
inline |
Definition at line 112 of file riemann_solver.template.h.
References ryujin::positive_part().
|
inline |
Definition at line 146 of file riemann_solver.h.
|
staticconstexpr |
The dimension of the state space.
Definition at line 37 of file riemann_solver.h.
|
staticconstexpr |
Definition at line 40 of file riemann_solver.h.
|
staticconstexpr |
The number of precomputed values.
Definition at line 56 of file riemann_solver.h.