ryujin 2.1.1 revision 9dcb748690310d6a540ebb8b066d1a0834fc7604
riemann_solver.h
Go to the documentation of this file.
1//
2// SPDX-License-Identifier: MIT
3// Copyright (C) 2020 - 2023 by the ryujin authors
4//
5
6#pragma once
7
8#include "hyperbolic_system.h"
9
10#include <simd.h>
11
12#include <deal.II/base/point.h>
13#include <deal.II/base/tensor.h>
14
15#include <functional>
16
17namespace ryujin
18{
19 namespace ShallowWater
20 {
30 template <int dim, typename Number = double>
32 {
33 public:
37 static constexpr unsigned int problem_dimension =
38 HyperbolicSystem::problem_dimension<dim>;
39
40 static constexpr unsigned int riemann_data_size = 3;
41 using primitive_type = std::array<Number, riemann_data_size>;
42
47
52
56 static constexpr unsigned int n_precomputed_values =
57 HyperbolicSystem::n_precomputed_values<dim>;
58
63
68
73 const HyperbolicSystem &hyperbolic_system,
75 &precomputed_values)
76 : hyperbolic_system(hyperbolic_system)
77 , precomputed_values(precomputed_values)
78 {
79 }
80
86 Number compute(const primitive_type &riemann_data_i,
87 const primitive_type &riemann_data_j) const;
88
96 Number compute(const state_type &U_i,
97 const state_type &U_j,
98 const unsigned int i,
99 const unsigned int *js,
100 const dealii::Tensor<1, dim, Number> &n_ij) const;
101
103
107
108 Number f(const primitive_type &primitive_state,
109 const Number &h_star) const;
110
111 Number phi(const primitive_type &riemann_data_i,
112 const primitive_type &riemann_data_j,
113 const Number &h) const;
114
115 Number lambda1_minus(const primitive_type &riemann_data,
116 const Number h_star) const;
117
118 Number lambda3_plus(const primitive_type &riemann_data,
119 const Number h_star) const;
120
121 Number compute_lambda(const primitive_type &riemann_data_i,
122 const primitive_type &riemann_data_j,
123 const Number h_star) const;
124
125 Number h_star_two_rarefaction(const primitive_type &riemann_data_i,
126 const primitive_type &riemann_data_j) const;
127
130 const dealii::Tensor<1, dim, Number> &n_ij) const;
131
132 private:
133 const HyperbolicSystem &hyperbolic_system;
134
136 &precomputed_values;
138 };
139
140
141 /* Inline definitions */
142
143
144 template <int dim, typename Number>
145 DEAL_II_ALWAYS_INLINE inline auto
147 const state_type &U, const dealii::Tensor<1, dim, Number> &n_ij) const
149 {
150 const Number h = hyperbolic_system.water_depth_sharp(U);
151 const Number gravity = hyperbolic_system.gravity();
152
153 const auto vel = hyperbolic_system.momentum(U) / h;
154 const auto proj_vel = n_ij * vel;
155 const auto a = std::sqrt(h * gravity);
156
157 return {{h, proj_vel, a}};
158 }
159
160
161 template <int dim, typename Number>
163 const state_type &U_i,
164 const state_type &U_j,
165 const unsigned int /*i*/,
166 const unsigned int * /*js*/,
167 const dealii::Tensor<1, dim, Number> &n_ij) const
168 {
169 const auto riemann_data_i = riemann_data_from_state(U_i, n_ij);
170 const auto riemann_data_j = riemann_data_from_state(U_j, n_ij);
171 return compute(riemann_data_i, riemann_data_j);
172 }
173
174
175 } // namespace ShallowWater
176} // namespace ryujin
std::array< Number, n_precomputed_values< dim > > precomputed_type
dealii::Tensor< 1, problem_dimension< dim >, Number > state_type
Number phi(const primitive_type &riemann_data_i, const primitive_type &riemann_data_j, const Number &h) const
static constexpr unsigned int problem_dimension
primitive_type riemann_data_from_state(const state_type &U, const dealii::Tensor< 1, dim, Number > &n_ij) const
Number lambda3_plus(const primitive_type &riemann_data, const Number h_star) const
Number lambda1_minus(const primitive_type &riemann_data, const Number h_star) const
static constexpr unsigned int n_precomputed_values
static constexpr unsigned int riemann_data_size
Number h_star_two_rarefaction(const primitive_type &riemann_data_i, const primitive_type &riemann_data_j) const
Number compute_lambda(const primitive_type &riemann_data_i, const primitive_type &riemann_data_j, const Number h_star) const
std::array< Number, riemann_data_size > primitive_type
typename get_value_type< Number >::type ScalarNumber
HyperbolicSystem::precomputed_type< dim, Number > precomputed_type
RiemannSolver(const HyperbolicSystem &hyperbolic_system, const MultiComponentVector< ScalarNumber, n_precomputed_values > &precomputed_values)
HyperbolicSystem::state_type< dim, Number > state_type
Number compute(const primitive_type &riemann_data_i, const primitive_type &riemann_data_j) const
Number f(const primitive_type &primitive_state, const Number &h_star) const