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/TMPL.hpp"
      23             : 
      24             : /// \cond
      25             : namespace domain {
      26             : namespace CoordinateMaps {
      27             : class Interval;
      28             : template <typename Map1, typename Map2, typename Map3>
      29             : class ProductOf3Maps;
      30             : }  // namespace CoordinateMaps
      31             : 
      32             : template <typename SourceFrame, typename TargetFrame, typename... Maps>
      33             : class CoordinateMap;
      34             : }  // namespace domain
      35             : /// \endcond
      36             : 
      37             : namespace domain::creators {
      38             : /// Create a 3D Domain with its computational domain being the\f$x-y\f$
      39             : /// plane. The third dimension uses a Cartoon basis with Killing vector along
      40             : /// the \f$\phi\f$ direction.
      41           1 : class CartoonCylinder : public DomainCreator<3> {
      42             :  public:
      43           0 :   using maps_list = tmpl::list<domain::CoordinateMap<
      44             :       Frame::BlockLogical, Frame::Inertial,
      45             :       CoordinateMaps::ProductOf3Maps<CoordinateMaps::Interval,
      46             :                                      CoordinateMaps::Interval,
      47             :                                      CoordinateMaps::Identity<1>>>>;
      48             : 
      49           0 :   static std::string name() { return "CartoonCylinder"; }
      50             : 
      51           0 :   struct LowerBounds {
      52           0 :     using type = std::array<double, 2>;
      53           0 :     static constexpr Options::String help = {
      54             :         "Lower bound in the [x,y] dimensions."};
      55             :   };
      56             : 
      57           0 :   struct UpperBounds {
      58           0 :     using type = std::array<double, 2>;
      59           0 :     static constexpr Options::String help = {
      60             :         "Upper bound in the [x,y] dimensions."};
      61             :   };
      62             : 
      63           0 :   struct InitialRefinement {
      64           0 :     using type = std::array<size_t, 2>;
      65           0 :     static constexpr Options::String help = {
      66             :         "Initial refinement level for the [x,y] dimensions."};
      67             :   };
      68             : 
      69           0 :   struct InitialGridPoints {
      70           0 :     using type = std::array<size_t, 2>;
      71           0 :     static constexpr Options::String help = {
      72             :         "Initial number of grid points in the [x,y] dimensions."};
      73             :   };
      74             : 
      75           0 :   struct Distributions {
      76           0 :     using type = std::array<CoordinateMaps::Distribution, 2>;
      77           0 :     static constexpr Options::String help = {
      78             :         "Distribution of grid points in the [x,y] dimensions."};
      79             :   };
      80             : 
      81             :   template <typename BoundaryConditionsBase>
      82           0 :   struct LowerUpperBoundaryCondition {
      83           0 :     static constexpr Options::String help =
      84             :         "Lower and upper Boundary Conditions";
      85           0 :     struct LowerBC {
      86           0 :       using type = std::unique_ptr<BoundaryConditionsBase>;
      87           0 :       static constexpr Options::String help = "Lower Boundary Condition";
      88           0 :       static std::string name() { return "Lower"; };
      89             :     };
      90           0 :     struct UpperBC {
      91           0 :       using type = std::unique_ptr<BoundaryConditionsBase>;
      92           0 :       static constexpr Options::String help = "Upper Boundary Condition";
      93           0 :       static std::string name() { return "Upper"; };
      94             :     };
      95           0 :     LowerUpperBoundaryCondition(typename LowerBC::type lower_bc,
      96             :                                 typename UpperBC::type upper_bc)
      97             :         : lower(std::move(lower_bc)), upper(std::move(upper_bc)){};
      98           0 :     LowerUpperBoundaryCondition() = default;
      99           0 :     std::unique_ptr<BoundaryConditionsBase> lower;
     100           0 :     std::unique_ptr<BoundaryConditionsBase> upper;
     101           0 :     using options = tmpl::list<LowerBC, UpperBC>;
     102             :   };
     103             : 
     104             :   template <typename BoundaryConditionsBase>
     105           0 :   struct BoundaryConditions {
     106           0 :     static constexpr Options::String help = {
     107             :         "The boundary conditions to be imposed in the x & y dimensions. "
     108             :         "Either specify one B.C. to be imposed for both "
     109             :         "lower and upper boundary or a pair 'Lower:' and 'Upper:'."};
     110           0 :     using type = std::array<
     111             :         std::variant<std::unique_ptr<BoundaryConditionsBase>,
     112             :                      LowerUpperBoundaryCondition<BoundaryConditionsBase>>,
     113             :         2>;
     114             :   };
     115             : 
     116           0 :   struct TimeDependence {
     117           0 :     using type =
     118             :         std::unique_ptr<domain::creators::time_dependence::TimeDependence<3>>;
     119           0 :     static constexpr Options::String help = {
     120             :         "The time dependence of the moving mesh domain. Specify `None` for no "
     121             :         "time dependant maps."};
     122             :   };
     123             : 
     124           0 :   using basic_options =
     125             :       tmpl::list<LowerBounds, UpperBounds, InitialRefinement, InitialGridPoints,
     126             :                  Distributions, TimeDependence>;
     127             : 
     128             :   template <typename Metavariables>
     129           0 :   using options = tmpl::conditional_t<
     130             :       domain::BoundaryConditions::has_boundary_conditions_base_v<
     131             :           typename Metavariables::system>,
     132             :       tmpl::push_back<
     133             :           basic_options,
     134             :           BoundaryConditions<
     135             :               domain::BoundaryConditions::get_boundary_conditions_base<
     136             :                   typename Metavariables::system>>>,
     137             :       basic_options>;
     138             : 
     139           0 :   static constexpr Options::String help{
     140             :       "A cylinder domain that requires/enforces axial symmetry. The "
     141             :       "computational domain is the x-y plane, with Cartoon partial derivatives "
     142             :       "being used for the z direction."};
     143             : 
     144           0 :   CartoonCylinder(
     145             :       std::array<double, 2> lower_bounds, std::array<double, 2> upper_bounds,
     146             :       std::array<size_t, 2> initial_refinement_levels,
     147             :       std::array<size_t, 2> initial_num_points,
     148             :       std::array<CoordinateMaps::Distribution, 2> distributions = {},
     149             :       std::unique_ptr<domain::creators::time_dependence::TimeDependence<3>>
     150             :           time_dependence = nullptr,
     151             :       std::array<std::array<std::unique_ptr<
     152             :                                 domain::BoundaryConditions::BoundaryCondition>,
     153             :                             2>,
     154             :                  2>
     155             :           boundary_conditions = {},
     156             :       const Options::Context& context = {});
     157             : 
     158             :   template <typename BoundaryConditionsBase>
     159           0 :   CartoonCylinder(
     160             :       std::array<double, 2> lower_bounds, std::array<double, 2> upper_bounds,
     161             :       std::array<size_t, 2> initial_refinement_levels,
     162             :       std::array<size_t, 2> initial_num_points,
     163             :       std::array<CoordinateMaps::Distribution, 2> distributions = {},
     164             :       std::unique_ptr<domain::creators::time_dependence::TimeDependence<3>>
     165             :           time_dependence = nullptr,
     166             :       std::array<
     167             :           std::variant<std::unique_ptr<BoundaryConditionsBase>,
     168             :                        LowerUpperBoundaryCondition<BoundaryConditionsBase>>,
     169             :           2>
     170             :           boundary_conditions = {},
     171             :       const Options::Context& context = {})
     172             :       : CartoonCylinder(
     173             :             lower_bounds, upper_bounds, initial_refinement_levels,
     174             :             initial_num_points,
     175             :             distributions, std::move(time_dependence),
     176             :             transform_boundary_conditions(std::move(boundary_conditions)),
     177             :             context) {}
     178             : 
     179           0 :   CartoonCylinder() = default;
     180           0 :   CartoonCylinder(const CartoonCylinder&) = delete;
     181           0 :   CartoonCylinder(CartoonCylinder&&) = default;
     182           0 :   CartoonCylinder& operator=(const CartoonCylinder&) = delete;
     183           0 :   CartoonCylinder& operator=(CartoonCylinder&&) = default;
     184           0 :   ~CartoonCylinder() override = default;
     185             : 
     186           0 :   Domain<3> create_domain() const override;
     187             : 
     188             :   std::vector<DirectionMap<
     189             :       3, std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>>>
     190           1 :   external_boundary_conditions() const override;
     191             : 
     192           1 :   std::vector<std::array<size_t, 3>> initial_extents() const override;
     193             : 
     194           1 :   std::vector<std::array<size_t, 3>> initial_refinement_levels() const override;
     195             : 
     196           1 :   std::vector<std::string> block_names() const override { return block_names_; }
     197             : 
     198             :   std::unordered_map<std::string, std::unordered_set<std::string>>
     199           1 :   block_groups() const override {
     200             :     return {{name(), {name()}}};
     201             :   }
     202             : 
     203           1 :   auto functions_of_time(const std::unordered_map<std::string, double>&
     204             :                              initial_expiration_times = {}) const
     205             :       -> std::unordered_map<
     206             :           std::string,
     207             :           std::unique_ptr<domain::FunctionsOfTime::FunctionOfTime>> override;
     208             : 
     209             :   // Transforms from option-created boundary conditions to the type used in the
     210             :   // constructor
     211             :   template <typename BoundaryConditionsBase>
     212           0 :   static auto transform_boundary_conditions(
     213             :       std::array<
     214             :           std::variant<std::unique_ptr<BoundaryConditionsBase>,
     215             :                        LowerUpperBoundaryCondition<BoundaryConditionsBase>>,
     216             :           2>
     217             :           boundary_conditions)
     218             :       -> std::array<
     219             :           std::array<
     220             :               std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>,
     221             :               2>,
     222             :           2>;
     223             : 
     224             :  private:
     225           0 :   std::array<double, 2> lower_bounds_{};
     226           0 :   std::array<double, 2> upper_bounds_{};
     227           0 :   std::array<size_t, 2> initial_refinement_levels_{};
     228           0 :   std::array<size_t, 2> initial_num_points_{};
     229           0 :   bool is_periodic_in_y_{};
     230           0 :   std::array<CoordinateMaps::Distribution, 2> distributions_{};
     231             :   std::array<
     232             :       std::array<std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>,
     233             :                  2>,
     234             :       2>
     235           0 :       boundary_conditions_{};
     236             :   std::unique_ptr<domain::creators::time_dependence::TimeDependence<3>>
     237           0 :       time_dependence_;
     238           0 :   size_t num_blocks_{};
     239           0 :   inline static const std::vector<std::string> block_names_{name()};
     240             : };
     241             : 
     242             : template <typename BoundaryConditionsBase>
     243             : auto CartoonCylinder::transform_boundary_conditions(
     244             :     std::array<
     245             :         std::variant<std::unique_ptr<BoundaryConditionsBase>,
     246             :                      LowerUpperBoundaryCondition<BoundaryConditionsBase>>,
     247             :         2>
     248             :         boundary_conditions)
     249             :     -> std::array<
     250             :         std::array<
     251             :             std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>, 2>,
     252             :         2> {
     253             :   std::array<
     254             :       std::array<std::unique_ptr<domain::BoundaryConditions::BoundaryCondition>,
     255             :                  2>,
     256             :       2>
     257             :       result{};
     258             :   for (size_t d = 0; d < 2; ++d) {
     259             :     if (std::holds_alternative<std::unique_ptr<BoundaryConditionsBase>>(
     260             :             boundary_conditions[d])) {
     261             :       auto bc = std::move(std::get<std::unique_ptr<BoundaryConditionsBase>>(
     262             :           boundary_conditions[d]));
     263             :       gsl::at(gsl::at(result, d), 0) = bc->get_clone();
     264             :       gsl::at(gsl::at(result, d), 1) = std::move(bc);
     265             :     } else {
     266             :       auto& bc = std::get<LowerUpperBoundaryCondition<BoundaryConditionsBase>>(
     267             :           boundary_conditions[d]);
     268             :       gsl::at(gsl::at(result, d), 0) = std::move(bc.lower);
     269             :       gsl::at(gsl::at(result, d), 1) = std::move(bc.upper);
     270             :     }
     271             :   }
     272             :   return result;
     273             : }
     274             : 
     275             : }  // namespace domain::creators