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 <unordered_set>
12 : #include <variant>
13 : #include <vector>
14 :
15 : #include "DataStructures/Tensor/IndexType.hpp"
16 : #include "Domain/BoundaryConditions/BoundaryCondition.hpp"
17 : #include "Domain/BoundaryConditions/GetBoundaryConditionsBase.hpp"
18 : #include "Domain/CoordinateMaps/Distribution.hpp"
19 : #include "Domain/Creators/DomainCreator.hpp"
20 : #include "Domain/Creators/TimeDependence/TimeDependence.hpp"
21 : #include "Domain/Domain.hpp"
22 : #include "Domain/Structure/DirectionMap.hpp"
23 : #include "Options/Context.hpp"
24 : #include "Options/String.hpp"
25 : #include "Utilities/TMPL.hpp"
26 :
27 : /// \cond
28 : namespace domain {
29 : namespace CoordinateMaps {
30 : class Affine;
31 : template <size_t Dim>
32 : class Identity;
33 : class Interval;
34 : class PolarToCartesian;
35 : template <typename Map1, typename Map2>
36 : class ProductOf2Maps;
37 : template <typename Map1, typename Map2, typename Map3>
38 : class ProductOf3Maps;
39 : } // namespace CoordinateMaps
40 :
41 : template <typename SourceFrame, typename TargetFrame, typename... Maps>
42 : class CoordinateMap;
43 : } // namespace domain
44 : /// \endcond
45 :
46 : namespace domain::creators {
47 : /*!
48 : * \brief Create a 3D filled cylinder domain with radial partitioning using a
49 : * B2/I1 filled cylinder at the center and Fourier hollow cylinders surrounding
50 : * it.
51 : */
52 1 : class AngularCylinder : public DomainCreator<3> {
53 : public:
54 0 : using maps_list = tmpl::list<domain::CoordinateMap<
55 : Frame::BlockLogical, Frame::Inertial,
56 : domain::CoordinateMaps::ProductOf3Maps<
57 : domain::CoordinateMaps::Affine, domain::CoordinateMaps::Identity<1>,
58 : domain::CoordinateMaps::Interval>,
59 : domain::CoordinateMaps::ProductOf2Maps<
60 : domain::CoordinateMaps::PolarToCartesian,
61 : domain::CoordinateMaps::Identity<1>>>>;
62 :
63 : /*!
64 : * \brief Radius of the cylinder's outer edge
65 : */
66 1 : struct OuterRadius {
67 0 : using type = double;
68 0 : static constexpr Options::String help = {"Radius of the cylinder."};
69 : };
70 :
71 : /*!
72 : * \brief Lower z-coordinate of the cylinder's base
73 : */
74 1 : struct LowerZBound {
75 0 : using type = double;
76 0 : static constexpr Options::String help = {
77 : "z-coordinate of the base of the cylinder."};
78 : };
79 :
80 : /*!
81 : * \brief Upper z-coordinate of the cylinder's top
82 : */
83 1 : struct UpperZBound {
84 0 : using type = double;
85 0 : static constexpr Options::String help = {
86 : "z-coordinate of the top of the cylinder."};
87 : };
88 :
89 : /*!
90 : * \brief Whether the cylinder is periodic in the z direction
91 : */
92 1 : struct IsPeriodicInZ {
93 0 : using type = bool;
94 0 : static constexpr Options::String help = {
95 : "True if periodic in the cylindrical z direction."};
96 : };
97 :
98 : /*!
99 : * \brief Radial coordinates of the boundaries splitting elements
100 : */
101 1 : struct RadialPartitioning {
102 0 : using type = std::vector<double>;
103 0 : static constexpr Options::String help = {
104 : "Radial coordinates of the boundaries splitting the inner cylinder and "
105 : "potential hollow cylinder shells"};
106 : };
107 :
108 : /*!
109 : * \brief Z-coordinates of the boundaries splitting the domain into layers
110 : */
111 1 : struct PartitioningInZ {
112 0 : using type = std::vector<double>;
113 0 : static constexpr Options::String help = {
114 : "z-coordinates of the boundaries splitting the domain into layers "
115 : "between LowerZBound and UpperZBound."};
116 : };
117 :
118 : /*!
119 : * \brief Initial number of \f$\theta\f$ gridpoints for filled cylinder. The
120 : * number for \f$r\f$ is accordingly set to match spectral space sizes. This
121 : * is enforced to be odd for numerical stability.
122 : */
123 1 : struct InitialCylinderThetaGridPoints {
124 0 : using type = size_t;
125 0 : static constexpr Options::String help = {
126 : "Initial number of grid points in theta (r is accordingly set). It "
127 : "must be an odd number for stability."};
128 : };
129 :
130 : /*!
131 : * \brief Initial number of z gridpoints for the cylinder
132 : */
133 1 : struct InitialCylinderZGridPoints {
134 0 : using type = size_t;
135 0 : static constexpr Options::String help = {
136 : "Initial number of grid points in z"};
137 : };
138 :
139 : /*!
140 : * \brief Initial number of \f$[r, \theta, z]\f$ gridpoints for hollow
141 : * cylinders. Can be one triplet which is applied to all shells, or each can
142 : * be specified. The theta component must be odd for numerical stability.
143 : */
144 1 : struct InitialHollowCylinderGridPoints {
145 0 : using type =
146 : std::variant<std::array<size_t, 3>, std::vector<std::array<size_t, 3>>>;
147 0 : static constexpr Options::String help = {
148 : "Initial number of grid points in [r, theta, z] for each hollow "
149 : "cylinder. If one triplet is given, it will be applied to all blocks, "
150 : "otherwise each hollow cylinder can be specified. The theta component "
151 : "must be odd for better numerical stability."};
152 : };
153 :
154 : /*!
155 : * \brief Grid point distribution along the z-axis for each layer
156 : */
157 1 : struct DistributionInZ {
158 0 : using type = std::vector<domain::CoordinateMaps::Distribution>;
159 0 : static constexpr Options::String help = {
160 : "Select the distribution of grid points along the z-axis in each "
161 : "layer. The lowermost layer must have a 'Linear' distribution, because "
162 : "both a 'Logarithmic' and 'Inverse' distribution places its "
163 : "singularity at 'LowerZBound'. The 'PartitioningInZ' determines the "
164 : "number of layers."};
165 0 : static size_t lower_bound_on_size() { return 1; }
166 : };
167 :
168 : /*!
169 : * \brief Initial refinement levels in the z direction
170 : */
171 1 : struct InitialRefinementInZ {
172 0 : using type = std::variant<size_t, std::vector<size_t>>;
173 0 : static constexpr Options::String help = {
174 : "Initial refinement level in the z-direction. Can be either a single "
175 : "value applied to all layers, or a vector with one value per layer. "
176 : "The number of layers is determined by 'PartitioningInZ'."};
177 : };
178 :
179 : /*!
180 : * \brief Time dependence of the domain
181 : */
182 1 : struct TimeDependence {
183 0 : using type =
184 : std::unique_ptr<domain::creators::time_dependence::TimeDependence<3>>;
185 0 : static constexpr Options::String help = {
186 : "The time dependence of the moving mesh domain. Specify `None` for no "
187 : "time dependent maps."};
188 : };
189 :
190 : /*!
191 : * \brief Boundary conditions group
192 : */
193 1 : struct BoundaryConditions {
194 0 : static constexpr Options::String help =
195 : "Options for the boundary conditions";
196 : };
197 :
198 : /*!
199 : * \brief Boundary condition on the lower base
200 : */
201 : template <typename BoundaryConditionsBase>
202 1 : struct LowerZBoundaryCondition {
203 0 : using group = BoundaryConditions;
204 0 : static std::string name() { return "LowerZ"; }
205 0 : static constexpr Options::String help =
206 : "The boundary condition to be imposed on the lower base of the "
207 : "cylinder, i.e. at the `LowerZBound` in the z-direction.";
208 0 : using type = std::unique_ptr<BoundaryConditionsBase>;
209 : };
210 :
211 : /*!
212 : * \brief Boundary condition on the upper base
213 : */
214 : template <typename BoundaryConditionsBase>
215 1 : struct UpperZBoundaryCondition {
216 0 : using group = BoundaryConditions;
217 0 : static std::string name() { return "UpperZ"; }
218 0 : static constexpr Options::String help =
219 : "The boundary condition to be imposed on the upper base of the "
220 : "cylinder, i.e. at the `UpperZBound` in the z-direction.";
221 0 : using type = std::unique_ptr<BoundaryConditionsBase>;
222 : };
223 :
224 : /*!
225 : * \brief Boundary condition on the radial boundary
226 : */
227 : template <typename BoundaryConditionsBase>
228 1 : struct MantleBoundaryCondition {
229 0 : using group = BoundaryConditions;
230 0 : static std::string name() { return "Mantle"; }
231 0 : static constexpr Options::String help =
232 : "The boundary condition to be imposed on the mantle of the "
233 : "cylinder, i.e. at the `OuterRadius` in the radial direction.";
234 0 : using type = std::unique_ptr<BoundaryConditionsBase>;
235 : };
236 :
237 0 : using basic_options =
238 : tmpl::list<OuterRadius, LowerZBound, UpperZBound, RadialPartitioning,
239 : PartitioningInZ, InitialCylinderThetaGridPoints,
240 : InitialCylinderZGridPoints, InitialHollowCylinderGridPoints,
241 : DistributionInZ, InitialRefinementInZ, TimeDependence>;
242 :
243 : template <typename Metavariables>
244 0 : using options = tmpl::append<
245 : basic_options,
246 : tmpl::conditional_t<
247 : domain::BoundaryConditions::has_boundary_conditions_base_v<
248 : typename Metavariables::system>,
249 : tmpl::list<
250 : LowerZBoundaryCondition<
251 : domain::BoundaryConditions::get_boundary_conditions_base<
252 : typename Metavariables::system>>,
253 : UpperZBoundaryCondition<
254 : domain::BoundaryConditions::get_boundary_conditions_base<
255 : typename Metavariables::system>>,
256 : MantleBoundaryCondition<
257 : domain::BoundaryConditions::get_boundary_conditions_base<
258 : typename Metavariables::system>>>,
259 : tmpl::list<IsPeriodicInZ>>>;
260 :
261 0 : static constexpr Options::String help{
262 : "Creates a cylinder using a Zernike basis radially, Fourier in the "
263 : "angular direction, and I1 in the z direction"};
264 :
265 0 : AngularCylinder(
266 : typename OuterRadius::type outer_radius,
267 : typename LowerZBound::type lower_z_bound,
268 : typename UpperZBound::type upper_z_bound,
269 : typename RadialPartitioning::type radial_partitioning,
270 : typename PartitioningInZ::type partitioning_in_z,
271 : typename InitialCylinderThetaGridPoints::type
272 : initial_cylinder_theta_grid_points,
273 : typename InitialCylinderZGridPoints::type initial_cylinder_z_grid_points,
274 : typename InitialHollowCylinderGridPoints::type
275 : initial_hollow_cylinder_grid_points,
276 : typename DistributionInZ::type distribution_in_z,
277 : typename InitialRefinementInZ::type initial_refinement_in_z,
278 : std::unique_ptr<domain::creators::time_dependence::TimeDependence<3>>
279 : time_dependence = nullptr,
280 : bool is_periodic_in_z = false,
281 : const Options::Context& context = {});
282 :
283 0 : AngularCylinder(
284 : typename OuterRadius::type outer_radius,
285 : typename LowerZBound::type lower_z_bound,
286 : typename UpperZBound::type upper_z_bound,
287 : typename RadialPartitioning::type radial_partitioning,
288 : typename PartitioningInZ::type partitioning_in_z,
289 : typename InitialCylinderThetaGridPoints::type
290 : initial_cylinder_theta_grid_points,
291 : typename InitialCylinderZGridPoints::type initial_cylinder_z_grid_points,
292 : typename InitialHollowCylinderGridPoints::type
293 : initial_hollow_cylinder_grid_points,
294 : typename DistributionInZ::type distribution_in_z,
295 : typename InitialRefinementInZ::type initial_refinement_in_z,
296 : std::unique_ptr<domain::creators::time_dependence::TimeDependence<3>>
297 : time_dependence = nullptr,
298 : std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>
299 : lower_z_boundary_condition = nullptr,
300 : std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>
301 : upper_z_boundary_condition = nullptr,
302 : std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>
303 : mantle_boundary_condition = nullptr,
304 : const Options::Context& context = {});
305 :
306 0 : AngularCylinder() = default;
307 0 : AngularCylinder(const AngularCylinder&) = delete;
308 0 : AngularCylinder(AngularCylinder&&) = default;
309 0 : AngularCylinder& operator=(const AngularCylinder&) = delete;
310 0 : AngularCylinder& operator=(AngularCylinder&&) = default;
311 0 : ~AngularCylinder() override = default;
312 :
313 0 : Domain<3> create_domain() const override;
314 :
315 : std::vector<DirectionMap<
316 : 3, std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>>>
317 1 : external_boundary_conditions() const override;
318 :
319 1 : std::vector<std::array<size_t, 3>> initial_extents() const override;
320 :
321 1 : std::vector<std::array<size_t, 3>> initial_refinement_levels() const override;
322 :
323 1 : std::vector<std::string> block_names() const override { return block_names_; }
324 :
325 : std::unordered_map<std::string, std::unordered_set<std::string>>
326 1 : block_groups() const override {
327 : return block_groups_;
328 : }
329 :
330 : private:
331 0 : typename OuterRadius::type outer_radius_{};
332 0 : typename LowerZBound::type lower_z_bound_{};
333 0 : typename UpperZBound::type upper_z_bound_{};
334 0 : typename IsPeriodicInZ::type is_periodic_in_z_{};
335 0 : typename RadialPartitioning::type radial_partitioning_{};
336 0 : typename PartitioningInZ::type partitioning_in_z_{};
337 : typename InitialCylinderThetaGridPoints::type
338 0 : initial_cylinder_theta_grid_points_{};
339 0 : typename InitialCylinderZGridPoints::type initial_cylinder_z_grid_points_{};
340 0 : std::vector<std::array<size_t, 3>> initial_hollow_cylinder_grid_points_{};
341 0 : DistributionInZ::type distribution_in_z_{};
342 0 : std::vector<size_t> initial_refinement_in_z_{};
343 : std::unique_ptr<domain::creators::time_dependence::TimeDependence<3>>
344 0 : time_dependence_ = nullptr;
345 : std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>
346 0 : lower_z_boundary_condition_{};
347 : std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>
348 0 : upper_z_boundary_condition_{};
349 : std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>
350 0 : mantle_boundary_condition_{};
351 0 : std::vector<std::string> block_names_;
352 : std::unordered_map<std::string, std::unordered_set<std::string>>
353 0 : block_groups_;
354 0 : size_t num_blocks_{};
355 : };
356 : } // namespace domain::creators
|
