Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <array>
7 : #include <cstddef>
8 : #include <memory>
9 : #include <string>
10 : #include <unordered_map>
11 : #include <vector>
12 :
13 : #include "Domain/BoundaryConditions/BoundaryCondition.hpp"
14 : #include "Domain/BoundaryConditions/GetBoundaryConditionsBase.hpp"
15 : #include "Domain/CoordinateMaps/Distribution.hpp"
16 : #include "Domain/Creators/DomainCreator.hpp"
17 : #include "Domain/Creators/TimeDependence/TimeDependence.hpp"
18 : #include "Domain/Domain.hpp"
19 : #include "Domain/Structure/DirectionMap.hpp"
20 : #include "Options/Context.hpp"
21 : #include "Options/String.hpp"
22 : #include "Utilities/GetOutput.hpp"
23 : #include "Utilities/MakeArray.hpp"
24 : #include "Utilities/TMPL.hpp"
25 :
26 : /// \cond
27 : namespace domain {
28 : namespace CoordinateMaps {
29 : class Affine;
30 : class Interval;
31 : template <typename Map1, typename Map2>
32 : class ProductOf2Maps;
33 : template <typename Map1, typename Map2, typename Map3>
34 : class ProductOf3Maps;
35 : } // namespace CoordinateMaps
36 :
37 : template <typename SourceFrame, typename TargetFrame, typename... Maps>
38 : class CoordinateMap;
39 : } // namespace domain
40 : /// \endcond
41 :
42 : namespace domain::creators {
43 :
44 : /// Create a domain consisting of a single Block in `Dim` dimensions.
45 : template <size_t Dim>
46 1 : class Rectilinear : public DomainCreator<Dim> {
47 : private:
48 : static_assert(Dim == 1 or Dim == 2 or Dim == 3,
49 : "Rectilinear domain is only implemented in 1, 2, or 3 "
50 : "dimensions.");
51 :
52 0 : using Interval = CoordinateMaps::Interval;
53 0 : using Interval2D = CoordinateMaps::ProductOf2Maps<Interval, Interval>;
54 0 : using Interval3D =
55 : CoordinateMaps::ProductOf3Maps<Interval, Interval, Interval>;
56 0 : using Affine = CoordinateMaps::Affine;
57 0 : using Affine2D = CoordinateMaps::ProductOf2Maps<Affine, Affine>;
58 0 : using Affine3D = CoordinateMaps::ProductOf3Maps<Affine, Affine, Affine>;
59 :
60 : public:
61 0 : using maps_list = tmpl::list<
62 : domain::CoordinateMap<
63 : Frame::BlockLogical, Frame::Inertial,
64 : tmpl::conditional_t<
65 : Dim == 1, Interval,
66 : tmpl::conditional_t<Dim == 2, Interval2D, Interval3D>>>,
67 : // Previous version of the domain used the `Affine` map, so we
68 : // need to keep it in this list for backwards compatibility.
69 : domain::CoordinateMap<
70 : Frame::BlockLogical, Frame::Inertial,
71 : tmpl::conditional_t<
72 : Dim == 1, Affine,
73 : tmpl::conditional_t<Dim == 2, Affine2D, Affine3D>>>>;
74 :
75 0 : static std::string name() {
76 : if constexpr (Dim == 1) {
77 : return "Interval";
78 : } else if constexpr (Dim == 2) {
79 : return "Rectangle";
80 : } else {
81 : return "Brick";
82 : }
83 : }
84 :
85 0 : struct LowerBound {
86 0 : using type = std::array<double, Dim>;
87 0 : static constexpr Options::String help = {"Lower bound in each dimension."};
88 : };
89 :
90 0 : struct UpperBound {
91 0 : using type = std::array<double, Dim>;
92 0 : static constexpr Options::String help = {"Upper bound in each dimension."};
93 : };
94 :
95 0 : struct Distribution {
96 0 : using type =
97 : std::array<CoordinateMaps::DistributionAndSingularityPosition, Dim>;
98 0 : static constexpr Options::String help = {
99 : "Distribution of grid points in each dimension"};
100 : };
101 :
102 0 : struct IsPeriodicIn {
103 0 : using type = std::array<bool, Dim>;
104 0 : static constexpr Options::String help = {"Periodicity in each dimension."};
105 : };
106 :
107 0 : struct InitialRefinement {
108 0 : using type = std::array<size_t, Dim>;
109 0 : static constexpr Options::String help = {
110 : "Initial refinement level in each dimension."};
111 : };
112 :
113 0 : struct InitialGridPoints {
114 0 : using type = std::array<size_t, Dim>;
115 0 : static constexpr Options::String help = {
116 : "Initial number of grid points in each dimension."};
117 : };
118 :
119 0 : struct TimeDependence {
120 0 : using type =
121 : std::unique_ptr<domain::creators::time_dependence::TimeDependence<Dim>>;
122 0 : static constexpr Options::String help = {
123 : "The time dependence of the moving mesh domain."};
124 : };
125 :
126 : template <typename BoundaryConditionsBase>
127 0 : struct LowerUpperBoundaryCondition {
128 0 : static constexpr Options::String help =
129 : "Lower and upper Boundary Conditions";
130 0 : struct LowerBC {
131 0 : using type = std::unique_ptr<BoundaryConditionsBase>;
132 0 : static constexpr Options::String help = "Lower Boundary Condition";
133 0 : static std::string name() { return "Lower"; };
134 : };
135 0 : struct UpperBC {
136 0 : using type = std::unique_ptr<BoundaryConditionsBase>;
137 0 : static constexpr Options::String help = "Upper Boundary Condition";
138 0 : static std::string name() { return "Upper"; };
139 : };
140 0 : LowerUpperBoundaryCondition(typename LowerBC::type lower_bc,
141 : typename UpperBC::type upper_bc)
142 : : lower(std::move(lower_bc)), upper(std::move(upper_bc)){};
143 0 : LowerUpperBoundaryCondition() = default;
144 0 : std::unique_ptr<BoundaryConditionsBase> lower;
145 0 : std::unique_ptr<BoundaryConditionsBase> upper;
146 0 : using options = tmpl::list<LowerBC, UpperBC>;
147 : };
148 :
149 : template <typename BoundaryConditionsBase>
150 0 : struct BoundaryConditions {
151 0 : static constexpr Options::String help = {
152 : "The boundary conditions to be imposed in each dimension. "
153 : "Either specify one B.C. to be imposed for both "
154 : "lower and upper boundary or a pair 'Lower:' and 'Upper:'."};
155 0 : using type = std::array<
156 : std::variant<std::unique_ptr<BoundaryConditionsBase>,
157 : LowerUpperBoundaryCondition<BoundaryConditionsBase>>,
158 : Dim>;
159 : };
160 :
161 : template <typename Metavariables>
162 0 : using options = tmpl::list<
163 : LowerBound, UpperBound, InitialRefinement, InitialGridPoints,
164 : tmpl::conditional_t<
165 : domain::BoundaryConditions::has_boundary_conditions_base_v<
166 : typename Metavariables::system>,
167 : BoundaryConditions<
168 : domain::BoundaryConditions::get_boundary_conditions_base<
169 : typename Metavariables::system>>,
170 : IsPeriodicIn>,
171 : Distribution, TimeDependence>;
172 :
173 0 : static constexpr Options::String help{"A rectilinear domain."};
174 :
175 0 : Rectilinear(
176 : std::array<double, Dim> lower_bounds,
177 : std::array<double, Dim> upper_bounds,
178 : std::array<size_t, Dim> initial_refinement_levels,
179 : std::array<size_t, Dim> initial_num_points,
180 : std::array<bool, Dim> is_periodic = make_array<Dim>(false),
181 : std::array<CoordinateMaps::DistributionAndSingularityPosition, Dim>
182 : distributions = {},
183 : std::unique_ptr<domain::creators::time_dependence::TimeDependence<Dim>>
184 : time_dependence = nullptr,
185 : const Options::Context& context = {});
186 :
187 0 : Rectilinear(
188 : std::array<double, Dim> lower_bounds,
189 : std::array<double, Dim> upper_bounds,
190 : std::array<size_t, Dim> initial_refinement_levels,
191 : std::array<size_t, Dim> initial_num_points,
192 : std::array<std::array<std::unique_ptr<
193 : domain::BoundaryConditions::BoundaryCondition>,
194 : 2>,
195 : Dim>
196 : boundary_conditions,
197 : std::array<CoordinateMaps::DistributionAndSingularityPosition, Dim>
198 : distributions = {},
199 : std::unique_ptr<domain::creators::time_dependence::TimeDependence<Dim>>
200 : time_dependence = nullptr,
201 : const Options::Context& context = {});
202 :
203 : template <typename BoundaryConditionsBase>
204 0 : Rectilinear(
205 : std::array<double, Dim> lower_bounds,
206 : std::array<double, Dim> upper_bounds,
207 : std::array<size_t, Dim> initial_refinement_levels,
208 : std::array<size_t, Dim> initial_num_points,
209 : std::array<
210 : std::variant<std::unique_ptr<BoundaryConditionsBase>,
211 : LowerUpperBoundaryCondition<BoundaryConditionsBase>>,
212 : Dim>
213 : boundary_conditions,
214 : std::array<CoordinateMaps::DistributionAndSingularityPosition, Dim>
215 : distributions = {},
216 : std::unique_ptr<domain::creators::time_dependence::TimeDependence<Dim>>
217 : time_dependence = nullptr,
218 : const Options::Context& context = {})
219 : : Rectilinear(
220 : lower_bounds, upper_bounds, initial_refinement_levels,
221 : initial_num_points,
222 : transform_boundary_conditions(std::move(boundary_conditions)),
223 : distributions, std::move(time_dependence), context) {}
224 :
225 0 : Rectilinear() = default;
226 0 : Rectilinear(const Rectilinear&) = delete;
227 0 : Rectilinear(Rectilinear&&) = default;
228 0 : Rectilinear& operator=(const Rectilinear&) = delete;
229 0 : Rectilinear& operator=(Rectilinear&&) = default;
230 0 : ~Rectilinear() override = default;
231 :
232 0 : Domain<Dim> create_domain() const override;
233 :
234 : std::vector<DirectionMap<
235 : Dim, std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>>>
236 1 : external_boundary_conditions() const override;
237 :
238 1 : std::vector<std::array<size_t, Dim>> initial_extents() const override;
239 :
240 1 : std::vector<std::array<size_t, Dim>> initial_refinement_levels()
241 : const override;
242 :
243 1 : auto functions_of_time(const std::unordered_map<std::string, double>&
244 : initial_expiration_times = {}) const
245 : -> std::unordered_map<
246 : std::string,
247 : std::unique_ptr<domain::FunctionsOfTime::FunctionOfTime>> override;
248 :
249 1 : std::vector<std::string> block_names() const override { return block_names_; }
250 :
251 : std::unordered_map<std::string, std::unordered_set<std::string>>
252 1 : block_groups() const override {
253 : return {{name(), {name()}}};
254 : }
255 :
256 : // Transforms from option-created boundary conditions to the type used in the
257 : // constructor
258 : template <typename BoundaryConditionsBase>
259 0 : static auto transform_boundary_conditions(
260 : std::array<
261 : std::variant<std::unique_ptr<BoundaryConditionsBase>,
262 : LowerUpperBoundaryCondition<BoundaryConditionsBase>>,
263 : Dim>
264 : boundary_conditions)
265 : -> std::array<
266 : std::array<
267 : std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>,
268 : 2>,
269 : Dim>;
270 :
271 : private:
272 0 : std::array<double, Dim> lower_bounds_{};
273 0 : std::array<double, Dim> upper_bounds_{};
274 : std::array<CoordinateMaps::DistributionAndSingularityPosition, Dim>
275 0 : distributions_{};
276 0 : std::array<bool, Dim> is_periodic_{};
277 0 : std::array<size_t, Dim> initial_refinement_levels_{};
278 0 : std::array<size_t, Dim> initial_num_points_{};
279 : std::unique_ptr<domain::creators::time_dependence::TimeDependence<Dim>>
280 0 : time_dependence_;
281 : std::array<
282 : std::array<std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>,
283 : 2>,
284 : Dim>
285 0 : boundary_conditions_{};
286 0 : inline static const std::vector<std::string> block_names_{name()};
287 : };
288 :
289 : template <size_t Dim>
290 : template <typename BoundaryConditionsBase>
291 : auto Rectilinear<Dim>::transform_boundary_conditions(
292 : std::array<
293 : std::variant<std::unique_ptr<BoundaryConditionsBase>,
294 : LowerUpperBoundaryCondition<BoundaryConditionsBase>>,
295 : Dim>
296 : boundary_conditions)
297 : -> std::array<
298 : std::array<
299 : std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>, 2>,
300 : Dim> {
301 : std::array<
302 : std::array<std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>,
303 : 2>,
304 : Dim>
305 : result{};
306 : for (size_t d = 0; d < Dim; ++d) {
307 : if (std::holds_alternative<std::unique_ptr<BoundaryConditionsBase>>(
308 : boundary_conditions[d])) {
309 : auto bc = std::move(std::get<std::unique_ptr<BoundaryConditionsBase>>(
310 : boundary_conditions[d]));
311 : result[d][0] = bc->get_clone();
312 : result[d][1] = std::move(bc);
313 : } else {
314 : auto& bc = std::get<LowerUpperBoundaryCondition<BoundaryConditionsBase>>(
315 : boundary_conditions[d]);
316 : result[d][0] = std::move(bc.lower);
317 : result[d][1] = std::move(bc.upper);
318 : }
319 : }
320 : return result;
321 : }
322 :
323 0 : using Interval = Rectilinear<1>;
324 0 : using Rectangle = Rectilinear<2>;
325 0 : using Brick = Rectilinear<3>;
326 :
327 : } // namespace domain::creators
|
