ryujin 2.1.1 revision ef0fcd4010d109b860652ece4a7b8963fb7d46b1
mesh_adaptor.h
Go to the documentation of this file.
1//
2// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
3// Copyright (C) 2024 by the ryujin authors
4//
5
6#pragma once
7
8#include <compile_time_options.h>
9
10#include "mpi_ensemble.h"
11#include "offline_data.h"
12
13#include <deal.II/base/parameter_acceptor.h>
14#include <deal.II/base/smartpointer.h>
15
16#include <random>
17
18namespace ryujin
19{
25 enum class AdaptationStrategy {
30
37
42 };
43
51 enum class MarkingStrategy {
61 };
62
73
78 };
79} // namespace ryujin
80
81#ifndef DOXYGEN
82DECLARE_ENUM(
84 LIST({ryujin::AdaptationStrategy::global_refinement, "global refinement"},
86 {ryujin::AdaptationStrategy::kelly_estimator, "kelly estimator"}, ));
87
88DECLARE_ENUM(ryujin::MarkingStrategy,
89 LIST({ryujin::MarkingStrategy::fixed_number, "fixed number"},
91 "fixed fraction"}, ));
92
95 "fixed time points"},
97 "simulation cycle"}, ));
98#endif
99
100namespace ryujin
101{
108 template <typename Description, int dim, typename Number = double>
109 class MeshAdaptor final : public dealii::ParameterAcceptor
110 {
111 public:
116
119
120 using View =
121 typename Description::template HyperbolicSystemView<dim, Number>;
122
123 static constexpr auto problem_dimension = View::problem_dimension;
124
126 using InitialPrecomputedVector = typename View::InitialPrecomputedVector;
128
130
134
138 MeshAdaptor(const MPIEnsemble &mpi_ensemble,
139 const OfflineData<dim, Number> &offline_data,
140 const HyperbolicSystem &hyperbolic_system,
141 const ParabolicSystem &parabolic_system,
142 const InitialPrecomputedVector &initial_precomputed,
143 const ScalarVector &alpha,
144 const std::string &subsection = "/MeshAdaptor");
145
150 void prepare(const Number t);
151
157 void analyze(const StateVector &state_vector,
158 const Number t,
159 unsigned int cycle);
160
168
169
174 dealii::Triangulation<dim> &triangulation) const;
175
176 private:
181
182 AdaptationStrategy adaptation_strategy_;
183 std::uint_fast64_t random_adaptation_mersenne_twister_seed_;
184
185 MarkingStrategy marking_strategy_;
186 double refinement_fraction_;
187 double coarsening_fraction_;
188 unsigned int min_refinement_level_;
189 unsigned int max_refinement_level_;
190 unsigned int max_num_cells_;
191
192 TimePointSelectionStrategy time_point_selection_strategy_;
193 std::vector<Number> adaptation_time_points_;
194 unsigned int adaptation_cycle_interval_;
195
196 std::vector<std::string> kelly_quantities_;
197
199
203
204 const MPIEnsemble &mpi_ensemble_;
205
206 dealii::SmartPointer<const OfflineData<dim, Number>> offline_data_;
207 dealii::SmartPointer<const HyperbolicSystem> hyperbolic_system_;
208 dealii::SmartPointer<const ParabolicSystem> parabolic_system_;
209
210 bool need_mesh_adaptation_;
211
212 mutable dealii::Vector<float> indicators_;
213
214 /* random adaptation: */
215
216 void compute_random_indicators() const;
217
218 mutable std::mt19937_64 mersenne_twister_;
219
220 /* Kelly estimator: */
221
222 void populate_kelly_quantities(const StateVector &state_vector) const;
223 void compute_kelly_indicators() const;
224
225 const InitialPrecomputedVector &initial_precomputed_;
226 const ScalarVector &alpha_;
227
228 mutable std::vector<ScalarVector> kelly_components_;
230 };
231
232} // namespace ryujin
typename View::StateVector StateVector
Definition: mesh_adaptor.h:125
Vectors::ScalarVector< Number > ScalarVector
Definition: mesh_adaptor.h:127
MeshAdaptor(const MPIEnsemble &mpi_ensemble, const OfflineData< dim, Number > &offline_data, const HyperbolicSystem &hyperbolic_system, const ParabolicSystem &parabolic_system, const InitialPrecomputedVector &initial_precomputed, const ScalarVector &alpha, const std::string &subsection="/MeshAdaptor")
static constexpr auto problem_dimension
Definition: mesh_adaptor.h:123
void analyze(const StateVector &state_vector, const Number t, unsigned int cycle)
typename View::InitialPrecomputedVector InitialPrecomputedVector
Definition: mesh_adaptor.h:126
typename Description::ParabolicSystem ParabolicSystem
Definition: mesh_adaptor.h:118
typename Description::HyperbolicSystem HyperbolicSystem
Definition: mesh_adaptor.h:117
void mark_cells_for_coarsening_and_refinement(dealii::Triangulation< dim > &triangulation) const
typename Description::template HyperbolicSystemView< dim, Number > View
Definition: mesh_adaptor.h:121
const auto & need_mesh_adaptation() const
Definition: mesh_adaptor.h:167
void prepare(const Number t)
TimePointSelectionStrategy
Definition: mesh_adaptor.h:68
AdaptationStrategy
Definition: mesh_adaptor.h:25
MarkingStrategy
Definition: mesh_adaptor.h:51
#define ACCESSOR_READ_ONLY(member)
dealii::LinearAlgebra::distributed::Vector< Number > ScalarVector
Definition: state_vector.h:33
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