ryujin 2.1.1 revision ef0fcd4010d109b860652ece4a7b8963fb7d46b1
parabolic_module.h
Go to the documentation of this file.
1//
2// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
3// Copyright (C) 2023 - 2024 by the ryujin authors
4//
5
6#pragma once
7
8#include <compile_time_options.h>
9
10#include "convenience_macros.h"
11#include "hyperbolic_module.h"
12#include "initial_values.h"
13#include "mpi_ensemble.h"
14#include "offline_data.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{
31 template <typename Description, int dim, typename Number = double>
32 class ParabolicModule final : public dealii::ParameterAcceptor
33 {
34 public:
39
41
42 using View =
43 typename Description::template HyperbolicSystemView<dim, Number>;
44
46
48 typename Description::template ParabolicSolver<dim, Number>;
49
51
53
57
62 const MPIEnsemble &mpi_ensemble,
63 std::map<std::string, dealii::Timer> &computing_timer,
64 const OfflineData<dim, Number> &offline_data,
65 const HyperbolicSystem &hyperbolic_system,
66 const ParabolicSystem &parabolic_system,
67 const InitialValues<Description, dim, Number> &initial_values,
68 const std::string &subsection = "/ParabolicModule");
69
75 void prepare();
76
78
82
90 void prepare_state_vector(StateVector &state_vector, Number t) const;
91
101 template <int stages>
102 void
103 backward_euler_step(const StateVector &old_state_vector,
104 const Number old_t,
105 std::array<std::reference_wrapper<const StateVector>,
106 stages> stage_state_vectors,
107 const std::array<Number, stages> stage_weights,
108 StateVector &new_state_vector,
109 Number tau) const;
110
119 void crank_nicolson_step(const StateVector &old_state_vector,
120 const Number old_t,
121 StateVector &new_state_vector,
122 Number tau) const;
123
129 void print_solver_statistics(std::ostream &output) const;
130
132
140
141
142
148
154
155 // FIXME: refactor to function
157
158 private:
160
164
165 // FIXME: refactor contents into this class.
166 ParabolicSolver parabolic_solver_;
167 mutable unsigned int cycle_;
168
169 mutable unsigned int n_restarts_;
170 mutable unsigned int n_corrections_;
171 mutable unsigned int n_warnings_;
172
174 };
175
176} /* namespace ryujin */
void prepare_state_vector(StateVector &state_vector, Number t) const
typename View::StateVector StateVector
void crank_nicolson_step(const StateVector &old_state_vector, const Number old_t, StateVector &new_state_vector, Number tau) const
const auto & n_warnings() const
const auto & n_corrections() const
IDViolationStrategy id_violation_strategy_
typename Description::ParabolicSystem ParabolicSystem
const auto & n_restarts() const
typename Description::template ParabolicSolver< dim, Number > ParabolicSolver
typename Description::template HyperbolicSystemView< dim, Number > View
ParabolicModule(const MPIEnsemble &mpi_ensemble, std::map< std::string, dealii::Timer > &computing_timer, const OfflineData< dim, Number > &offline_data, const HyperbolicSystem &hyperbolic_system, const ParabolicSystem &parabolic_system, const InitialValues< Description, dim, Number > &initial_values, const std::string &subsection="/ParabolicModule")
typename Description::HyperbolicSystem HyperbolicSystem
void print_solver_statistics(std::ostream &output) const
void backward_euler_step(const StateVector &old_state_vector, const Number old_t, std::array< std::reference_wrapper< const StateVector >, stages > stage_state_vectors, const std::array< Number, stages > stage_weights, StateVector &new_state_vector, Number tau) const
#define ACCESSOR_READ_ONLY(member)
std::tuple< MultiComponentVector< Number, problem_dim >, MultiComponentVector< Number, prec_dim >, BlockVector< Number > > StateVector
Definition: state_vector.h:53
Euler::HyperbolicSystem HyperbolicSystem
Definition: description.h:34
ryujin::StubParabolicSystem ParabolicSystem
Definition: description.h:39