ryujin 2.1.1 revision 9dcb748690310d6a540ebb8b066d1a0834fc7604
hyperbolic_module.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 <compile_time_options.h>
9
10#include "convenience_macros.h"
11#include "initial_values.h"
12#include "offline_data.h"
13#include "simd.h"
14#include "sparse_matrix_simd.h"
15
16#include <deal.II/base/parameter_acceptor.h>
17#include <deal.II/base/timer.h>
18#include <deal.II/lac/sparse_matrix.templates.h>
19#include <deal.II/lac/vector.h>
20
21#include <functional>
22
23namespace ryujin
24{
44 warn,
45
52 };
53
54
61 class Restart final
62 {
63 };
64
65
74 template <typename Description, int dim, typename Number = double>
75 class HyperbolicModule final : public dealii::ParameterAcceptor
76 {
77 public:
82
87 typename HyperbolicSystem::template View<dim, Number>;
88
92 static constexpr unsigned int problem_dimension =
93 HyperbolicSystemView::problem_dimension;
94
98 using state_type = typename HyperbolicSystemView::state_type;
99
103 using flux_type = typename HyperbolicSystemView::flux_type;
104
109
113 using vector_type = typename HyperbolicSystemView::vector_type;
114
118 static constexpr unsigned int n_precomputed_values =
119 HyperbolicSystemView::n_precomputed_values;
120
124 static constexpr unsigned int n_precomputation_cycles =
125 HyperbolicSystemView::n_precomputation_cycles;
126
131 typename HyperbolicSystemView::precomputed_vector_type;
132
136 static constexpr unsigned int n_precomputed_initial_values =
137 HyperbolicSystemView::n_precomputed_initial_values;
138
143 typename HyperbolicSystemView::precomputed_initial_vector_type;
144
145
150 const MPI_Comm &mpi_communicator,
151 std::map<std::string, dealii::Timer> &computing_timer,
154 const InitialValues<Description, dim, Number> &initial_values,
155 const std::string &subsection = "/HyperbolicModule");
156
162 void prepare();
163
168
227 template <int stages>
228 Number
229 step(const vector_type &old_U,
230 std::array<std::reference_wrapper<const vector_type>, stages> stage_U,
231 std::array<std::reference_wrapper<const precomputed_vector_type>,
232 stages> stage_precomputed,
233 const std::array<Number, stages> stage_weights,
234 vector_type &new_U,
235 precomputed_vector_type &new_precomputed,
236 Number tau = Number(0.)) const;
237
251 void apply_boundary_conditions(vector_type &U, Number t) const;
252
260 void cfl(Number new_cfl) const
261 {
262 Assert(cfl_ > Number(0.), dealii::ExcInternalError());
263 cfl_ = new_cfl;
264 }
265
271
272
276
281
286
294
299
305
306 // FIXME: refactor to function
307 mutable bool precompute_only_;
308
309 // FIXME: refactor to function
311
312 private:
314
318
319 unsigned int limiter_iter_;
320 Number limiter_newton_tolerance_;
321 unsigned int limiter_newton_max_iter_;
322 Number limiter_relaxation_factor_;
323
324 bool cfl_with_boundary_dofs_;
325
327
329
333
334 const MPI_Comm &mpi_communicator_;
335 std::map<std::string, dealii::Timer> &computing_timer_;
336
337 dealii::SmartPointer<const OfflineData<dim, Number>> offline_data_;
338 dealii::SmartPointer<const HyperbolicSystem> hyperbolic_system_;
339 dealii::SmartPointer<const InitialValues<Description, dim, Number>>
340 initial_values_;
341
342 mutable Number cfl_;
343
344 mutable unsigned int n_restarts_;
345
346 mutable unsigned int n_warnings_;
347
348 precomputed_initial_vector_type precomputed_initial_;
349
350 mutable scalar_type alpha_;
351
352 static constexpr auto n_bounds =
353 Description::template Limiter<dim, Number>::n_bounds;
354 mutable MultiComponentVector<Number, n_bounds> bounds_;
355
356 mutable vector_type source_;
357
358 mutable vector_type r_;
359 mutable vector_type source_r_;
360
361 mutable SparseMatrixSIMD<Number> dij_matrix_;
362 mutable SparseMatrixSIMD<Number> lij_matrix_;
363 mutable SparseMatrixSIMD<Number> lij_matrix_next_;
364 mutable SparseMatrixSIMD<Number, problem_dimension> pij_matrix_;
365 mutable SparseMatrixSIMD<Number, problem_dimension> qij_matrix_;
366
368 };
369
370} /* namespace ryujin */
auto & hyperbolic_system() const
static constexpr unsigned int problem_dimension
static constexpr unsigned int n_precomputation_cycles
typename Description::HyperbolicSystem HyperbolicSystem
HyperbolicModule(const MPI_Comm &mpi_communicator, std::map< std::string, dealii::Timer > &computing_timer, const OfflineData< dim, Number > &offline_data, const HyperbolicSystem &hyperbolic_system, const InitialValues< Description, dim, Number > &initial_values, const std::string &subsection="/HyperbolicModule")
typename HyperbolicSystemView::precomputed_vector_type precomputed_vector_type
typename HyperbolicSystem::template View< dim, Number > HyperbolicSystemView
typename HyperbolicSystemView::precomputed_initial_vector_type precomputed_initial_vector_type
void apply_boundary_conditions(vector_type &U, Number t) const
auto & precomputed_initial() const
typename HyperbolicSystemView::state_type state_type
static constexpr unsigned int n_precomputed_values
typename OfflineData< dim, Number >::scalar_type scalar_type
void cfl(Number new_cfl) const
IDViolationStrategy id_violation_strategy_
typename HyperbolicSystemView::vector_type vector_type
Number step(const vector_type &old_U, std::array< std::reference_wrapper< const vector_type >, stages > stage_U, std::array< std::reference_wrapper< const precomputed_vector_type >, stages > stage_precomputed, const std::array< Number, stages > stage_weights, vector_type &new_U, precomputed_vector_type &new_precomputed, Number tau=Number(0.)) const
typename HyperbolicSystemView::flux_type flux_type
static constexpr unsigned int n_precomputed_initial_values
#define ACCESSOR_READ_ONLY(member)
Euler::HyperbolicSystem HyperbolicSystem
Definition: description.h:32