ryujin 2.1.1 revision 2062b43aae59b7ee510dec4ae666e29d56c4f0be
hyperbolic_system.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
9#include <discretization.h>
11#include <patterns_conversion.h>
12#include <simd.h>
13
14#include <deal.II/base/parameter_acceptor.h>
15#include <deal.II/base/tensor.h>
16
17#include <array>
18#include <functional>
19
20namespace ryujin
21{
22 namespace Skeleton
23 {
24 template <int dim, typename Number>
25 class HyperbolicSystemView;
26
32 class HyperbolicSystem final : public dealii::ParameterAcceptor
33 {
34 public:
38 static inline const std::string problem_name =
39 "Skeleton Hyperbolic System";
40
44 HyperbolicSystem(const std::string &subsection = "/HyperbolicSystem")
45 : ParameterAcceptor(subsection)
46 {
47 }
48
55 template <int dim, typename Number>
56 auto view() const
57 {
59 }
60
61 private:
62 template <int dim, typename Number>
64 }; /* HyperbolicSystem */
65
66
72 template <int dim, typename Number>
74 {
75 public:
80 HyperbolicSystemView(const HyperbolicSystem &hyperbolic_system)
81 : hyperbolic_system_(hyperbolic_system)
82 {
83 }
84
88 template <int dim2, typename Number2>
89 auto view() const
90 {
91 return HyperbolicSystemView<dim2, Number2>{hyperbolic_system_};
92 }
93
98
99
100 private:
101 const HyperbolicSystem &hyperbolic_system_;
102
103
104 public:
109
113 static constexpr unsigned int problem_dimension = 1;
114
119 using state_type = dealii::Tensor<1, problem_dimension, Number>;
120
125
130 static inline const auto component_names =
131 []() -> std::array<std::string, problem_dimension> {
132 if constexpr (dim == 1)
133 return {"u"};
134 else if constexpr (dim == 2)
135 return {"u"};
136 else if constexpr (dim == 3)
137 return {"u"};
138 __builtin_trap();
139 }();
140
144 using primitive_state_type = dealii::Tensor<1, problem_dimension, Number>;
145
150 static inline const auto primitive_component_names =
151 std::array<std::string, problem_dimension>{"u"};
152
156 using flux_type =
157 dealii::Tensor<1, problem_dimension, dealii::Tensor<1, dim, Number>>;
158
163
165
169
173 static constexpr unsigned int n_precomputed_initial_values = 0;
174
179 std::array<Number, n_precomputed_initial_values>;
180
187
191 static inline const auto precomputed_initial_names =
192 std::array<std::string, n_precomputed_initial_values>{};
193
197 static constexpr unsigned int n_precomputed_values = 0;
198
202 using precomputed_state_type = std::array<Number, n_precomputed_values>;
203
209
213 static inline const auto precomputed_names =
214 std::array<std::string, n_precomputed_values>{};
215
219 static constexpr unsigned int n_precomputation_cycles = 0;
220
225 template <typename DISPATCH, typename SPARSITY>
226 void precomputation_loop(unsigned int /*cycle*/,
227 const DISPATCH &dispatch_check,
228 precomputed_vector_type & /*precomputed_values*/,
229 const SPARSITY & /*sparsity_simd*/,
230 const vector_type & /*U*/,
231 unsigned int /*left*/,
232 unsigned int /*right*/) const = delete;
233
235
239
245 bool is_admissible(const state_type & /*U*/) const
246 {
247 return true;
248 }
249
251
255
259 template <typename Lambda>
261 const dealii::types::boundary_id /*id*/,
262 const state_type &U,
263 const dealii::Tensor<1, dim, Number> & /*normal*/,
264 const Lambda & /*get_dirichlet_data*/) const
265 {
266 return U;
267 }
268
270
274
296 const precomputed_initial_vector_type & /*piv*/,
297 const unsigned int /*i*/,
298 const state_type & /*U_i*/) const
299 {
300 return flux_contribution_type{};
301 }
302
305 const precomputed_initial_vector_type & /*piv*/,
306 const unsigned int * /*js*/,
307 const state_type & /*U_j*/) const
308 {
309 return flux_contribution_type{};
310 }
311
317 const flux_contribution_type & /*flux_j*/) const
318 {
319 return flux_type{};
320 }
321
325 static constexpr bool have_high_order_flux = false;
326
328 const flux_contribution_type &) const = delete;
329
331
335
337 static constexpr bool have_source_terms = false;
338
340 const unsigned int i,
341 const state_type &U_i,
342 const ScalarNumber t,
343 const ScalarNumber tau) const = delete;
344
346 const unsigned int i,
347 const state_type &U_i,
348 const ScalarNumber t,
349 const ScalarNumber tau) const = delete;
350
352
356
367 template <typename ST>
368 state_type expand_state(const ST &state) const
369 {
370 return state;
371 }
372
378 from_primitive_state(const primitive_state_type &primitive_state) const
379 {
380 return primitive_state;
381 }
382
388 {
389 return state;
390 }
391
397 template <typename Lambda>
399 const Lambda & /*lambda*/) const
400 {
401 return state;
402 }
403
404 }; /* HyperbolicSystemView */
405 } // namespace Skeleton
406} // namespace ryujin
state_type apply_galilei_transform(const state_type &state, const Lambda &) const
bool is_admissible(const state_type &) const
primitive_state_type to_primitive_state(const state_type &state) const
static constexpr unsigned int n_precomputed_values
std::array< Number, n_precomputed_values > precomputed_state_type
dealii::Tensor< 1, problem_dimension, Number > state_type
flux_type high_order_flux(const flux_contribution_type &, const flux_contribution_type &) const =delete
HyperbolicSystemView(const HyperbolicSystem &hyperbolic_system)
flux_contribution_type flux_contribution(const precomputed_vector_type &, const precomputed_initial_vector_type &, const unsigned int *, const state_type &) const
dealii::Tensor< 1, problem_dimension, Number > primitive_state_type
state_type expand_state(const ST &state) const
state_type from_primitive_state(const primitive_state_type &primitive_state) const
typename get_value_type< Number >::type ScalarNumber
void precomputation_loop(unsigned int, const DISPATCH &dispatch_check, precomputed_vector_type &, const SPARSITY &, const vector_type &, unsigned int, unsigned int) const =delete
dealii::Tensor< 1, problem_dimension, dealii::Tensor< 1, dim, Number > > flux_type
state_type apply_boundary_conditions(const dealii::types::boundary_id, const state_type &U, const dealii::Tensor< 1, dim, Number > &, const Lambda &) const
flux_type flux(const flux_contribution_type &, const flux_contribution_type &) const
flux_contribution_type flux_contribution(const precomputed_vector_type &, const precomputed_initial_vector_type &, const unsigned int, const state_type &) const
static constexpr unsigned int n_precomputation_cycles
static constexpr unsigned int problem_dimension
static constexpr unsigned int n_precomputed_initial_values
std::array< Number, n_precomputed_initial_values > precomputed_initial_state_type
state_type high_order_source(const precomputed_vector_type &pv, const unsigned int i, const state_type &U_i, const ScalarNumber t, const ScalarNumber tau) const =delete
state_type low_order_source(const precomputed_vector_type &pv, const unsigned int i, const state_type &U_i, const ScalarNumber t, const ScalarNumber tau) const =delete
static const std::string problem_name
HyperbolicSystem(const std::string &subsection="/HyperbolicSystem")