10#include <compile_time_options.h>
18#include <deal.II/base/parameter_acceptor.h>
19#include <deal.II/base/tensor.h>
27 template <
int dim,
typename Number>
28 class HyperbolicSystemView;
48 "Compressible Euler equations (arbitrary EOS)";
61 template <
int dim,
typename Number>
73 std::string equation_of_state_;
74 double reference_density_;
75 double vacuum_state_relaxation_small_;
76 double vacuum_state_relaxation_large_;
77 bool compute_strict_bounds_;
80 equation_of_state_list_;
83 std::shared_ptr<EquationOfState> selected_equation_of_state_;
85 template <
int dim,
typename Number>
107 template <
int dim,
typename Number>
116 : hyperbolic_system_(hyperbolic_system)
123 template <
int dim2,
typename Number2>
141 return hyperbolic_system_.equation_of_state_;
146 return hyperbolic_system_.reference_density_;
152 return hyperbolic_system_.vacuum_state_relaxation_small_;
158 return hyperbolic_system_.vacuum_state_relaxation_large_;
163 return hyperbolic_system_.compute_strict_bounds_;
177 const Number &e)
const
179 const auto &eos = hyperbolic_system_.selected_equation_of_state_;
181 if constexpr (std::is_same_v<ScalarNumber, Number>) {
185 for (
unsigned int k = 0; k < Number::size(); ++k) {
196 DEAL_II_ALWAYS_INLINE
inline Number
199 const auto &eos = hyperbolic_system_.selected_equation_of_state_;
201 if constexpr (std::is_same_v<ScalarNumber, Number>) {
202 return ScalarNumber(eos->specific_internal_energy(rho, p));
205 for (
unsigned int k = 0; k < Number::size(); ++k) {
206 e[k] =
ScalarNumber(eos->specific_internal_energy(rho[k], p[k]));
217 const Number &e)
const
219 const auto &eos = hyperbolic_system_.selected_equation_of_state_;
221 if constexpr (std::is_same_v<ScalarNumber, Number>) {
225 for (
unsigned int k = 0; k < Number::size(); ++k) {
236 DEAL_II_ALWAYS_INLINE
inline Number
239 const auto &eos = hyperbolic_system_.selected_equation_of_state_;
241 if constexpr (std::is_same_v<ScalarNumber, Number>) {
245 for (
unsigned int k = 0; k < Number::size(); ++k) {
257 const auto &eos = hyperbolic_system_.selected_equation_of_state_;
295 using state_type = dealii::Tensor<1, problem_dimension, Number>;
307 []() -> std::array<std::string, problem_dimension> {
308 if constexpr (dim == 1)
309 return {
"rho",
"m",
"E"};
310 else if constexpr (dim == 2)
311 return {
"rho",
"m_1",
"m_2",
"E"};
312 else if constexpr (dim == 3)
313 return {
"rho",
"m_1",
"m_2",
"m_3",
"E"};
324 []() -> std::array<std::string, problem_dimension> {
325 if constexpr (dim == 1)
326 return {
"rho",
"v",
"e"};
327 else if constexpr (dim == 2)
328 return {
"rho",
"v_1",
"v_2",
"e"};
329 else if constexpr (dim == 3)
330 return {
"rho",
"v_1",
"v_2",
"v_3",
"e"};
338 dealii::Tensor<1, problem_dimension, dealii::Tensor<1, dim, Number>>;
360 std::array<Number, n_precomputed_initial_values>;
373 std::array<std::string, n_precomputed_initial_values>{};
395 std::array<std::string, n_precomputed_values>{
398 "surrogate_specific_entropy",
399 "surrogate_harten_entropy"}};
410 template <
typename DISPATCH,
typename SPARSITY>
412 const DISPATCH &dispatch_check,
414 const SPARSITY &sparsity_simd,
417 unsigned int right)
const;
482 const Number &gamma_min)
const;
493 const Number &gamma_min)
const;
506 const Number &gamma_min)
const;
555 template <
int component>
559 const dealii::Tensor<1, dim, Number> &normal)
const;
579 template <
typename Lambda>
583 const dealii::Tensor<1, dim, Number> &normal,
584 const Lambda &get_dirichlet_data)
const;
626 const unsigned int i,
632 const unsigned int *js,
642 const dealii::Tensor<1, dim, Number> &c_ij)
const;
652 const dealii::Tensor<1, dim, Number> &c_ij)
const =
delete;
663 const unsigned int i,
668 const unsigned int *js,
688 template <
typename ST>
703 template <
typename ST>
723 template <
typename Lambda>
725 const Lambda &lambda)
const;
738 const std::string &subsection )
739 : ParameterAcceptor(subsection)
741 equation_of_state_ =
"polytropic gas";
745 "The equation of state. Valid names are given by any of the "
746 "subsections defined below");
748 compute_strict_bounds_ =
true;
750 "compute strict bounds",
751 compute_strict_bounds_,
752 "Compute strict, but significantly more expensive bounds at various "
753 "places: (a) an expensive, but better upper wavespeed estimate in "
754 "the approximate RiemannSolver; (b) entropy viscosity-commutator "
755 "with correct gamma_min over the stencil; (c) mathematically correct "
756 "surrogate specific entropy minimum with gamma_min over the "
759 reference_density_ = 1.;
760 add_parameter(
"reference density",
762 "Problem specific density reference");
764 vacuum_state_relaxation_small_ = 1.e2;
765 add_parameter(
"vacuum state relaxation small",
766 vacuum_state_relaxation_small_,
767 "Problem specific vacuum relaxation parameter");
769 vacuum_state_relaxation_large_ = 1.e4;
770 add_parameter(
"vacuum state relaxation large",
771 vacuum_state_relaxation_large_,
772 "Problem specific vacuum relaxation parameter");
779 equation_of_state_list_, subsection);
781 const auto populate_functions = [
this]() {
782 bool initialized =
false;
783 for (
auto &it : equation_of_state_list_)
786 if (it->name() == equation_of_state_) {
787 selected_equation_of_state_ = it;
789 "Compressible Euler equations (" + it->name() +
" EOS)";
797 "Could not find an equation of state description with name \"" +
798 equation_of_state_ +
"\""));
801 ParameterAcceptor::parse_parameters_call_back.connect(populate_functions);
802 populate_functions();
806 template <
int dim,
typename Number>
807 template <
typename DISPATCH,
typename SPARSITY>
808 DEAL_II_ALWAYS_INLINE
inline void
810 unsigned int cycle [[maybe_unused]],
811 const DISPATCH &dispatch_check,
813 const SPARSITY &sparsity_simd,
816 unsigned int right)
const
818 Assert(cycle == 0 || cycle == 1, dealii::ExcInternalError());
822 const auto &eos = hyperbolic_system_.selected_equation_of_state_;
823 unsigned int stride_size = get_stride_size<Number>;
826 if (eos->prefer_vector_interface()) {
831 const auto offset = left;
832 const auto size = right - left;
834 static std::vector<double> p;
835 static std::vector<double> rho;
836 static std::vector<double> e;
845 for (
unsigned int i = 0; i < size; i += stride_size) {
846 const auto U_i = U.template get_tensor<Number>(offset + i);
847 const auto rho_i = density(U_i);
848 const auto e_i = internal_energy(U_i) / rho_i;
854 write_entry<Number>(rho, rho_i, i);
855 write_entry<Number>(e, e_i, i);
861 eos->pressure(p, rho, e);
865 for (
unsigned int i = 0; i < size; i += stride_size) {
867 const unsigned int row_length = sparsity_simd.row_length(i);
874 const auto U_i = U.template get_tensor<Number>(offset + i);
875 const auto p_i = get_entry<Number>(p, i);
876 const auto gamma_i = surrogate_gamma(U_i, p_i);
877 const PT prec_i{p_i, gamma_i, Number(0.), Number(0.)};
878 precomputed_values.template write_tensor<Number>(prec_i,
888 for (
unsigned int i = left; i < right; i += stride_size) {
890 const unsigned int row_length = sparsity_simd.row_length(i);
896 const auto U_i = U.template get_tensor<Number>(i);
897 const auto rho_i = density(U_i);
898 const auto e_i = internal_energy(U_i) / rho_i;
899 const auto p_i = eos_pressure(rho_i, e_i);
901 const auto gamma_i = surrogate_gamma(U_i, p_i);
903 const PT prec_i{p_i, gamma_i, Number(0.), Number(0.)};
904 precomputed_values.template write_tensor<Number>(prec_i, i);
911 for (
unsigned int i = left; i < right; i += stride_size) {
915 const unsigned int row_length = sparsity_simd.row_length(i);
921 const auto U_i = U.template get_tensor<Number>(i);
922 auto prec_i = precomputed_values.template get_tensor<Number, PT>(i);
923 auto &[p_i, gamma_min_i, s_i, eta_i] = prec_i;
925 const unsigned int *js = sparsity_simd.columns(i) + stride_size;
926 for (
unsigned int col_idx = 1; col_idx < row_length;
927 ++col_idx, js += stride_size) {
929 const auto U_j = U.template get_tensor<Number>(js);
931 precomputed_values.template get_tensor<Number, PT>(js);
932 auto &[p_j, gamma_min_j, s_j, eta_j] = prec_j;
933 const auto gamma_j = surrogate_gamma(U_j, p_j);
934 gamma_min_i = std::min(gamma_min_i, gamma_j);
937 s_i = surrogate_specific_entropy(U_i, gamma_min_i);
938 eta_i = surrogate_harten_entropy(U_i, gamma_min_i);
939 precomputed_values.template write_tensor<Number>(prec_i, i);
945 template <
int dim,
typename Number>
946 DEAL_II_ALWAYS_INLINE
inline Number
953 template <
int dim,
typename Number>
954 DEAL_II_ALWAYS_INLINE
inline Number
956 const Number &rho)
const
958 constexpr ScalarNumber eps = std::numeric_limits<ScalarNumber>::epsilon();
959 const Number rho_cutoff_large =
960 reference_density() * vacuum_state_relaxation_large() * eps;
962 return dealii::compare_and_apply_mask<dealii::SIMDComparison::less_than>(
963 std::abs(rho), rho_cutoff_large, Number(0.), rho);
967 template <
int dim,
typename Number>
968 DEAL_II_ALWAYS_INLINE
inline dealii::Tensor<1, dim, Number>
971 dealii::Tensor<1, dim, Number> result;
972 for (
unsigned int i = 0; i < dim; ++i)
973 result[i] = U[1 + i];
978 template <
int dim,
typename Number>
979 DEAL_II_ALWAYS_INLINE
inline Number
986 template <
int dim,
typename Number>
987 DEAL_II_ALWAYS_INLINE
inline Number
993 const Number rho_inverse =
ScalarNumber(1.) / density(U);
994 const auto m = momentum(U);
995 const Number E = total_energy(U);
996 return E -
ScalarNumber(0.5) * m.norm_square() * rho_inverse;
1000 template <
int dim,
typename Number>
1001 DEAL_II_ALWAYS_INLINE
inline auto
1012 const Number rho_inverse =
ScalarNumber(1.) / density(U);
1013 const auto u = momentum(U) * rho_inverse;
1018 for (
unsigned int i = 0; i < dim; ++i) {
1019 result[1 + i] = -u[i];
1027 template <
int dim,
typename Number>
1028 DEAL_II_ALWAYS_INLINE
inline Number
1030 const state_type &U,
const Number &gamma_min)
const
1034 const auto rho = density(U);
1036 const auto interpolation_b = Number(eos_interpolation_b());
1037 const auto covolume = Number(1.) - interpolation_b * rho;
1039 return internal_energy(U) *
1040 ryujin::pow(rho_inverse - interpolation_b, gamma_min) / covolume;
1044 template <
int dim,
typename Number>
1045 DEAL_II_ALWAYS_INLINE
inline Number
1047 const state_type &U,
const Number &gamma_min)
const
1049 const auto rho = density(U);
1050 const auto m = momentum(U);
1051 const auto E = total_energy(U);
1052 const auto rho_rho_e = rho * E -
ScalarNumber(0.5) * m.norm_square();
1054 const auto exponent =
ScalarNumber(1.) / (gamma_min + Number(1.));
1056 const auto interpolation_b = Number(eos_interpolation_b());
1057 const auto covolume = Number(1.) - interpolation_b * rho;
1058 const auto covolume_term =
ryujin::pow(covolume, gamma_min - Number(1.));
1060 return ryujin::pow(rho_rho_e * covolume_term, exponent);
1064 template <
int dim,
typename Number>
1065 DEAL_II_ALWAYS_INLINE
inline auto
1067 const state_type &U,
const Number &eta,
const Number &gamma_min)
const
1084 const auto rho = density(U);
1085 const auto m = momentum(U);
1086 const auto E = total_energy(U);
1087 const auto rho_rho_e = rho * E -
ScalarNumber(0.5) * m.norm_square();
1089 const auto interpolation_b = Number(eos_interpolation_b());
1090 const auto covolume = Number(1.) - interpolation_b * rho;
1091 const auto covolume_inverse = Number(1.) / covolume;
1093 const auto factor =
ryujin::pow(eta * covolume_inverse, -gamma_min) *
1094 fixed_power<2>(covolume_inverse) /
1095 (gamma_min + Number(1.));
1099 result[0] = factor * (covolume * E - (gamma_min - Number(1.)) *
1100 rho_rho_e * interpolation_b);
1101 for (
unsigned int i = 0; i < dim; ++i)
1102 result[1 + i] = -factor * covolume * m[i];
1103 result[dim + 1] = factor * covolume * rho;
1109 template <
int dim,
typename Number>
1110 DEAL_II_ALWAYS_INLINE
inline Number
1112 const Number &p)
const
1114 const auto rho = density(U);
1115 const auto rho_e = internal_energy(U);
1116 const auto interpolation_b = Number(eos_interpolation_b());
1117 const auto covolume = Number(1.) - interpolation_b * rho;
1119 return Number(1.) + p * covolume / rho_e;
1123 template <
int dim,
typename Number>
1124 DEAL_II_ALWAYS_INLINE
inline Number
1126 const state_type &U,
const Number &gamma)
const
1128 const auto rho = density(U);
1129 const auto rho_e = internal_energy(U);
1130 const auto interpolation_b = Number(eos_interpolation_b());
1131 const auto covolume = Number(1.) - interpolation_b * rho;
1133 return (gamma - Number(1.)) * rho_e / covolume;
1137 template <
int dim,
typename Number>
1138 DEAL_II_ALWAYS_INLINE
inline bool
1141 const auto rho = density(U);
1142 const auto e = internal_energy(U);
1144 constexpr auto gt = dealii::SIMDComparison::greater_than;
1147 dealii::compare_and_apply_mask<gt>(rho, T(0.), T(0.), T(-1.)) +
1148 dealii::compare_and_apply_mask<gt>(e, T(0.), T(0.), T(-1.));
1151 if (!(test == Number(0.))) {
1152 std::cout << std::fixed << std::setprecision(16);
1153 std::cout <<
"Bounds violation: Negative state [rho, e] detected!\n";
1154 std::cout <<
"\t\trho: " << rho <<
"\n";
1155 std::cout <<
"\t\tint: " << e <<
"\n";
1159 return (test == Number(0.));
1163 template <
int dim,
typename Number>
1164 template <
int component>
1165 DEAL_II_ALWAYS_INLINE
inline auto
1169 const dealii::Tensor<1, dim, Number> &normal)
const ->
state_type
1173 static_assert(component == 1 || component == 2,
1174 "component has to be 1 or 2");
1176 const auto m = momentum(U);
1177 const auto rho = density(U);
1178 const auto vn = m * normal / rho;
1180 const auto p = surrogate_pressure(U);
1181 const auto gamma = surrogate_gamma(U, p);
1182 const auto interpolation_b =
ScalarNumber(eos_interpolation_b());
1183 const auto x = Number(1.) - interpolation_b * rho;
1184 const auto a = std::sqrt(gamma * p / (rho * x));
1187 const auto m_bar = momentum(U_bar);
1188 const auto rho_bar = density(U_bar);
1189 const auto vn_bar = m_bar * normal / rho_bar;
1191 const auto p_bar = surrogate_pressure(U_bar);
1192 const auto gamma_bar = surrogate_gamma(U_bar, p_bar);
1193 const auto x_bar = Number(1.) - interpolation_b * rho_bar;
1194 const auto a_bar = std::sqrt(gamma_bar * p_bar / (rho_bar * x_bar));
1198 const auto R_1 = component == 1
1199 ? vn_bar - 2. * a_bar / (gamma_bar - 1.) * x_bar
1200 : vn - 2. * a / (gamma - 1.) * x;
1204 const auto R_2 = component == 2
1205 ? vn_bar + 2. * a_bar / (gamma_bar - 1.) * x_bar
1206 : vn + 2. * a / (gamma - 1.) * x;
1210 const auto vperp = m / rho - vn * normal;
1212 const auto vn_new = 0.5 * (R_1 + R_2);
1215 1. / (gamma * s) *
ryujin::pow((gamma - 1.) / 4. * (R_2 - R_1), 2.);
1217 1. / (
ryujin::pow(rho_new, 1. / (1. - gamma)) + interpolation_b);
1219 const auto x_new = 1. - interpolation_b * rho_new;
1226 for (
unsigned int d = 0; d < dim; ++d) {
1227 U_new[1 + d] = rho_new * (vn_new * normal + vperp)[d];
1229 U_new[1 + dim] = p_new / (gamma - 1.) +
1230 0.5 * rho_new * (vn_new * vn_new + vperp.norm_square());
1236 template <
int dim,
typename Number>
1237 template <
typename Lambda>
1238 DEAL_II_ALWAYS_INLINE
inline auto
1240 dealii::types::boundary_id
id,
1242 const dealii::Tensor<1, dim, Number> &normal,
1243 const Lambda &get_dirichlet_data)
const ->
state_type
1248 result = get_dirichlet_data();
1251 auto m = momentum(U);
1252 m -= 1. * (m * normal) * normal;
1253 for (
unsigned int k = 0; k < dim; ++k)
1254 result[k + 1] = m[k];
1257 for (
unsigned int k = 0; k < dim; ++k)
1258 result[k + 1] = Number(0.);
1275 const auto m = momentum(U);
1276 const auto rho = density(U);
1277 const auto a = speed_of_sound(U);
1278 const auto vn = m * normal / rho;
1282 result = get_dirichlet_data();
1286 if (vn >= -a && vn <= 0.) {
1287 const auto U_dirichlet = get_dirichlet_data();
1288 result = prescribe_riemann_characteristic<2>(U_dirichlet, U, normal);
1292 if (vn > 0. && vn <= a) {
1293 const auto U_dirichlet = get_dirichlet_data();
1294 result = prescribe_riemann_characteristic<1>(U, U_dirichlet, normal);
1305 template <
int dim,
typename Number>
1306 DEAL_II_ALWAYS_INLINE
inline auto
1310 const auto rho_inverse =
ScalarNumber(1.) / density(U);
1311 const auto m = momentum(U);
1312 const auto E = total_energy(U);
1317 for (
unsigned int i = 0; i < dim; ++i) {
1318 result[1 + i] = m * (m[i] * rho_inverse);
1319 result[1 + i][i] += p;
1321 result[dim + 1] = m * (rho_inverse * (E + p));
1327 template <
int dim,
typename Number>
1328 DEAL_II_ALWAYS_INLINE
inline auto
1332 const unsigned int i,
1335 const auto &[p_i, gamma_min_i, s_i, eta_i] =
1336 pv.template get_tensor<Number, precomputed_state_type>(i);
1341 template <
int dim,
typename Number>
1342 DEAL_II_ALWAYS_INLINE
inline auto
1346 const unsigned int *js,
1349 const auto &[p_j, gamma_min_j, s_j, eta_j] =
1350 pv.template get_tensor<Number, precomputed_state_type>(js);
1355 template <
int dim,
typename Number>
1356 DEAL_II_ALWAYS_INLINE
inline auto
1360 const dealii::Tensor<1, dim, Number> &c_ij)
const ->
state_type
1366 template <
int dim,
typename Number>
1367 template <
typename ST>
1371 using T =
typename ST::value_type;
1372 static_assert(std::is_same_v<Number, T>,
"template mismatch");
1374 constexpr auto dim2 = ST::dimension - 2;
1375 static_assert(dim >= dim2,
1376 "the space dimension of the argument state must not be "
1377 "larger than the one of the target state");
1380 result[0] = state[0];
1381 result[dim + 1] = state[dim2 + 1];
1382 for (
unsigned int i = 1; i < dim2 + 1; ++i)
1383 result[i] = state[i];
1389 template <
int dim,
typename Number>
1390 template <
typename ST>
1391 DEAL_II_ALWAYS_INLINE
inline auto
1395 auto primitive_state = expand_state(initial_state);
1398 const auto rho = density(primitive_state);
1399 const auto p = total_energy(primitive_state);
1400 const auto e = eos_specific_internal_energy(rho, p);
1401 primitive_state[dim + 1] = e;
1403 return from_primitive_state(primitive_state);
1407 template <
int dim,
typename Number>
1408 DEAL_II_ALWAYS_INLINE
inline auto
1412 const auto rho = density(primitive_state);
1414 const auto u = momentum(primitive_state);
1416 const auto &e = total_energy(primitive_state);
1418 auto state = primitive_state;
1420 for (
unsigned int i = 1; i < dim + 1; ++i)
1424 state[dim + 1] = rho * e + Number(0.5) * rho * u * u;
1430 template <
int dim,
typename Number>
1431 DEAL_II_ALWAYS_INLINE
inline auto
1435 const auto rho = density(state);
1436 const auto rho_inverse = Number(1.) / rho;
1437 const auto rho_e = internal_energy(state);
1439 auto primitive_state = state;
1441 for (
unsigned int i = 1; i < dim + 1; ++i)
1442 primitive_state[i] *= rho_inverse;
1444 primitive_state[dim + 1] = rho_e * rho_inverse;
1446 return primitive_state;
1450 template <
int dim,
typename Number>
1451 template <
typename Lambda>
1455 auto result = state;
1456 const auto M = lambda(momentum(state));
1457 for (
unsigned int d = 0; d < dim; ++d)
1458 result[1 + d] = M[d];
dealii::Tensor< 1, problem_dimension, Number > state_type
state_type from_initial_state(const ST &initial_state) const
Number surrogate_harten_entropy(const state_type &U, const Number &gamma_min) const
std::array< Number, n_precomputed_values > precomputed_state_type
static Number internal_energy(const state_type &U)
state_type to_primitive_state(const state_type &state) 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 =delete
static constexpr bool have_high_order_flux
std::array< Number, n_precomputed_initial_values > precomputed_initial_state_type
typename get_value_type< Number >::type ScalarNumber
state_type nodal_source(const precomputed_vector_type &pv, const unsigned int i, const state_type &U_i, const ScalarNumber tau) const =delete
static dealii::Tensor< 1, dim, Number > momentum(const state_type &U)
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
Number surrogate_pressure(const state_type &U, const Number &gamma) const
DEAL_II_ALWAYS_INLINE const std::string & equation_of_state() const
DEAL_II_ALWAYS_INLINE ScalarNumber vacuum_state_relaxation_large() const
static constexpr unsigned int n_precomputed_values
DEAL_II_ALWAYS_INLINE ScalarNumber reference_density() const
static const auto component_names
state_type from_primitive_state(const state_type &primitive_state) const
static const auto precomputed_initial_names
static constexpr bool have_eos_interpolation_b
static state_type internal_energy_derivative(const state_type &U)
DEAL_II_ALWAYS_INLINE Number eos_temperature(const Number &rho, const Number &e) const
DEAL_II_ALWAYS_INLINE ScalarNumber eos_interpolation_b() const
static constexpr bool have_source_terms
DEAL_II_ALWAYS_INLINE Number eos_pressure(const Number &rho, const Number &e) const
static constexpr unsigned int problem_dimension
flux_type flux_contribution_type
Number surrogate_specific_entropy(const state_type &U, const Number &gamma_min) const
static const auto precomputed_names
static const auto primitive_component_names
static Number total_energy(const state_type &U)
bool is_admissible(const state_type &U) const
DEAL_II_ALWAYS_INLINE bool compute_strict_bounds() 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
DEAL_II_ALWAYS_INLINE Number eos_specific_internal_energy(const Number &rho, const Number &p) const
state_type nodal_source(const precomputed_vector_type &pv, const unsigned int *js, const state_type &U_j, const ScalarNumber tau) const =delete
static constexpr unsigned int n_precomputation_cycles
HyperbolicSystemView(const HyperbolicSystem &hyperbolic_system)
Number filter_vacuum_density(const Number &rho) const
flux_type f(const state_type &U, const Number &p) 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
state_type expand_state(const ST &state) const
state_type surrogate_harten_entropy_derivative(const state_type &U, const Number &eta, const Number &gamma_min) const
DEAL_II_ALWAYS_INLINE ScalarNumber vacuum_state_relaxation_small() const
Number surrogate_gamma(const state_type &U, const Number &p) const
state_type prescribe_riemann_characteristic(const state_type &U, const state_type &U_bar, const dealii::Tensor< 1, dim, Number > &normal) const
DEAL_II_ALWAYS_INLINE Number eos_speed_of_sound(const Number &rho, const Number &e) const
static constexpr unsigned int n_precomputed_initial_values
dealii::Tensor< 1, problem_dimension, dealii::Tensor< 1, dim, Number > > flux_type
static Number density(const state_type &U)
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
static constexpr bool have_gamma
state_type apply_galilei_transform(const state_type &state, const Lambda &lambda) const
HyperbolicSystem(const std::string &subsection="/HyperbolicSystem")
static std::string problem_name
typename get_value_type< Number >::type ScalarNumber
std::array< Number, n_precomputed_values > precomputed_state_type
dealii::Tensor< 1, problem_dimension, Number > state_type
void populate_equation_of_state_list(equation_of_state_list_type &equation_of_state_list, const std::string &subsection)
#define RYUJIN_OMP_SINGLE
T pow(const T x, const T b)
std::set< std::shared_ptr< EquationOfState > > equation_of_state_list_type
DEAL_II_ALWAYS_INLINE FT add(const FT &flux_left_ij, const FT &flux_right_ij)
DEAL_II_ALWAYS_INLINE dealii::Tensor< 1, problem_dim, T > contract(const FT &flux_ij, const TT &c_ij)