CylindricalBinaryCompactObject.hpp
1 // 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 <variant>
12 #include <vector>
13 
14 #include "Domain/BoundaryConditions/BoundaryCondition.hpp"
15 #include "Domain/BoundaryConditions/GetBoundaryConditionsBase.hpp"
18 #include "Domain/Creators/TimeDependence/TimeDependence.hpp"
19 #include "Domain/Domain.hpp"
20 #include "Options/Options.hpp"
21 #include "Utilities/TMPL.hpp"
22 
23 /// \cond
24 namespace domain {
25 namespace CoordinateMaps {
26 class Affine;
27 template <typename Map1, typename Map2>
28 class ProductOf2Maps;
29 template <typename Map1, typename Map2, typename Map3>
30 class ProductOf3Maps;
31 template <size_t VolumeDim>
32 class Wedge;
33 template <size_t VolumeDim>
34 class DiscreteRotation;
35 class CylindricalEndcap;
36 class CylindricalFlatEndcap;
37 class CylindricalFlatSide;
38 class CylindricalSide;
39 } // namespace CoordinateMaps
40 
41 template <typename SourceFrame, typename TargetFrame, typename... Maps>
42 class CoordinateMap;
43 
44 namespace FunctionsOfTime {
45 class FunctionOfTime;
46 } // namespace FunctionsOfTime
47 } // namespace domain
48 
49 namespace Frame {
50 struct Inertial;
51 struct Logical;
52 } // namespace Frame
53 /// \endcond
54 
55 namespace domain::creators {
56 
57 /*!
58  * \ingroup ComputationalDomainGroup
59  *
60  * \brief A general domain for two compact objects based on cylinders.
61  *
62  * Creates a 3D Domain that represents a binary compact object
63  * solution. This domain is described briefly in the Appendix of
64  * \cite Buchman:2012dw, and is illustrated in Figure 20 of that
65  * paper.
66  *
67  * In the code and options below, `ObjectA` and `ObjectB` refer to the
68  * two compact objects. In the grid frame, `ObjectA` is located to the
69  * right of (i.e. a more positive value of the x-coordinate than)
70  * `ObjectB`. The inner edge of the Blocks surrounding each of
71  * `ObjectA` and `ObjectB` is spherical in grid coordinates; the
72  * user must specify the center and radius of this surface for both
73  * `ObjectA` and `ObjectB`, and the user must specify the outer boundary
74  * radius. The outer boundary is a sphere centered at the origin.
75  *
76  * Note that Figure 20 of \cite Buchman:2012dw illustrates additional
77  * spherical shells inside the "EA" and "EB" blocks, and the caption
78  * of Figure 20 indicates that there are additional spherical shells
79  * outside the "CA" and "CB" blocks; `CylindricalBinaryCompactObject`
80  * does not have any of these spherical shells: the "EA" and "EB"
81  * blocks extend to the excision boundary and the "CA" and "CB" blocks
82  * extend to the outer boundary.
83  *
84  * The Blocks are named as follows:
85  * - Each of CAFilledCylinder EAFilledCylinder, EBFilledCylinder,
86  * MAFilledCylinder, MBFilledCylinder, and CBFilledCylinder consists
87  * of 5 blocks, named 'Center', 'East', 'North', 'West', and
88  * 'South', so an example of a valid block name is
89  * 'CAFilledCylinderCenter'.
90  * - Each of CACylinder, EACylinder, EBCylinder, and CBCylinder
91  * consists of 4 blocks, named 'East', 'North', 'West', and 'South',
92  * so an example of a valid block name is 'CACylinderEast'.
93  * - The Block group called "Outer" consists of all the CA and CB blocks. They
94  * all border the outer boundary.
95  * - The Block group called "InnerA" consists of all the EA, and MA
96  * blocks. They all border the inner boundary "A".
97  * - The Block group called "InnerB" consists of all the EB, and MB
98  * blocks. They all border the inner boundary "B".
99  *
100  * If \f$c_A\f$ and \f$c_B\f$ are the input parameters center_A and
101  * center_B, \f$r_A\f$ and \f$r_B\f$ are the input parameters radius_A and
102  * radius_B, and \f$R\f$ is the outer boundary radius, we demand the
103  * following restrictions on parameters:
104  * - \f$c_A^0>0\f$; this is a convention to simplify the code.
105  * - \f$c_B^0<0\f$; this is a convention to simplify the code.
106  * - \f$|c_A^0|\le|c_B^0|\f$. We should roughly have \f$r_A c_A^0 + r_B c_B^0\f$
107  * close to zero; that is, for BBHs (where \f$r_A\f$ is roughly twice the
108  * mass of the heavier object A, and \f$r_B\f$ is roughly twice the mass
109  * of the lighter object B) the center of mass should be roughly
110  * at the origin.
111  * - \f$0 < r_B < r_A\f$
112  * - \f$R \ge 3(|c_A^0|-|c_B^0|)\f$; otherwise the blocks will be too compressed
113  * near the outer boundary.
114  *
115  */
117  public:
118  using maps_list = tmpl::list<
155 
156  struct CenterA {
157  using type = std::array<double, 3>;
158  static constexpr Options::String help = {
159  "Grid coordinates of center for Object A, which is at x>0."};
160  };
161  struct CenterB {
162  using type = std::array<double, 3>;
163  static constexpr Options::String help = {
164  "Grid coordinates of center for Object B, which is at x<0."};
165  };
166  struct RadiusA {
167  using type = double;
168  static constexpr Options::String help = {
169  "Grid-coordinate radius of grid boundary around Object A."};
170  };
171  struct RadiusB {
172  using type = double;
173  static constexpr Options::String help = {
174  "Grid-coordinate radius of grid boundary around Object B."};
175  };
176  struct OuterRadius {
177  using type = double;
178  static constexpr Options::String help = {
179  "Grid-coordinate radius of outer boundary."};
180  };
181 
183  using type =
187  static constexpr Options::String help = {
188  "Initial refinement level. Specify one of: a single number, a list "
189  "representing [r, theta, perp], or such a list for every block in the "
190  "domain. Here 'r' is the radial direction normal to the inner and "
191  "outer boundaries, 'theta' is the periodic direction, and 'perp' is "
192  "the third direction."};
193  };
195  using type =
199  static constexpr Options::String help = {
200  "Initial number of grid points. Specify one of: a single number, a "
201  "list representing [r, theta, perp], or such a list for every block in "
202  "the domain. Here 'r' is the radial direction normal to the inner and "
203  "outer boundaries, 'theta' is the periodic direction, and 'perp' is "
204  "the third direction."};
205  };
206 
207  struct TimeDependence {
208  using type =
210  static constexpr Options::String help = {
211  "The time dependence of the moving mesh domain."};
212  };
213 
215  static constexpr Options::String help = "The boundary conditions to apply.";
216  };
217  template <typename BoundaryConditionsBase>
219  static std::string name() noexcept { return "InnerBoundary"; }
220  static constexpr Options::String help =
221  "Options for the inner boundary conditions.";
223  using group = BoundaryConditions;
224  };
225 
226  template <typename BoundaryConditionsBase>
228  static std::string name() noexcept { return "OuterBoundary"; }
229  static constexpr Options::String help =
230  "Options for the outer boundary conditions.";
232  using group = BoundaryConditions;
233  };
234 
235  using basic_options =
236  tmpl::list<CenterA, CenterB, RadiusA, RadiusB, OuterRadius,
238 
239  template <typename Metavariables>
240  using options = tmpl::conditional_t<
241  domain::BoundaryConditions::has_boundary_conditions_base_v<
242  typename Metavariables::system>,
244  basic_options,
247  typename Metavariables::system>>,
250  typename Metavariables::system>>>,
251  basic_options>;
252 
253  static constexpr Options::String help{
254  "The CylindricalBinaryCompactObject domain is a general domain for "
255  "two compact objects. The user must provide the (grid-frame) "
256  "centers and radii of the spherical inner edge of the grid surrounding "
257  "each of the two compact objects A and B."};
258 
260  typename CenterA::type center_A, typename CenterB::type center_B,
261  typename RadiusA::type radius_A, typename RadiusB::type radius_B,
262  typename OuterRadius::type outer_radius,
263  const typename InitialRefinement::type& initial_refinement,
264  const typename InitialGridPoints::type& initial_grid_points,
266  time_dependence = nullptr,
268  inner_boundary_condition = nullptr,
270  outer_boundary_condition = nullptr,
271  const Options::Context& context = {});
272 
273  CylindricalBinaryCompactObject() = default;
274  CylindricalBinaryCompactObject(const CylindricalBinaryCompactObject&) =
275  delete;
276  CylindricalBinaryCompactObject(CylindricalBinaryCompactObject&&) noexcept =
277  default;
278  CylindricalBinaryCompactObject& operator=(
279  const CylindricalBinaryCompactObject&) = delete;
280  CylindricalBinaryCompactObject& operator=(
281  CylindricalBinaryCompactObject&&) noexcept = default;
282  ~CylindricalBinaryCompactObject() noexcept override = default;
283 
284  Domain<3> create_domain() const noexcept override;
285 
286  std::vector<std::array<size_t, 3>> initial_extents() const noexcept override;
287 
288  std::vector<std::array<size_t, 3>> initial_refinement_levels()
289  const noexcept override;
290 
291  auto functions_of_time() const noexcept -> std::unordered_map<
292  std::string,
293  std::unique_ptr<domain::FunctionsOfTime::FunctionOfTime>> override;
294 
295  private:
296  // Note that center_A_ and center_B_ are rotated with respect to the
297  // input centers (which are in the grid frame), so that we can
298  // construct the map in a frame where the centers are offset in the
299  // z direction. At the end, there will be another rotation back to
300  // the grid frame (where the centers are offset in the x direction).
301  typename CenterA::type center_A_{};
302  typename CenterB::type center_B_{};
303  typename RadiusA::type radius_A_{};
304  typename RadiusB::type radius_B_{};
305  typename OuterRadius::type outer_radius_{};
306  typename std::vector<std::array<size_t, 3>> initial_refinement_{};
307  typename std::vector<std::array<size_t, 3>> initial_grid_points_{};
308  // cut_spheres_offset_factor_ is eta in Eq. (A.9) of
309  // https://arxiv.org/abs/1206.3015. cut_spheres_offset_factor_
310  // could be set to unity to simplify the equations. Here we fix it
311  // to the value 0.99 used in SpEC, so that we reproduce SpEC's
312  // domain decomposition.
313  double cut_spheres_offset_factor_{0.99};
314  // z_cutting_plane_ is x_C in Eq. (A.9) of
315  // https://arxiv.org/abs/1206.3015 (but rotated to the z-axis).
316  double z_cutting_plane_{};
317  // number_of_blocks_ could be eliminated or just set to its
318  // constant value of 46. But this value will change with the
319  // next PR that adds support for domains with unequal-sized objects.
320  size_t number_of_blocks_{};
322  time_dependence_;
324  inner_boundary_condition_;
326  outer_boundary_condition_;
327 };
328 } // namespace domain::creators
domain::CoordinateMaps::ProductOf3Maps
Product of three one-dimensional CoordinateMaps.
Definition: ProductMaps.hpp:90
domain::creators::CylindricalBinaryCompactObject
A general domain for two compact objects based on cylinders.
Definition: CylindricalBinaryCompactObject.hpp:116
domain::BoundaryConditions::get_boundary_conditions_base
typename detail::get_boundary_conditions_base< T >::type get_boundary_conditions_base
Returns T::boundary_condition_base or a placeholder class.
Definition: GetBoundaryConditionsBase.hpp:33
std::string
Frame::Inertial
Definition: IndexType.hpp:44
domain::CoordinateMaps::CylindricalSide
Map from a 3D unit right cylindrical shell to a volume that connects portions of two spherical surfac...
Definition: CylindricalSide.hpp:91
Options.hpp
vector
Domain.hpp
domain::CoordinateMaps::ProductOf2Maps
Product of two codimension=0 CoordinateMaps.
Definition: ProductMaps.hpp:35
domain::creators
Defines classes that create Domains.
Definition: AlignedLattice.hpp:38
domain::creators::CylindricalBinaryCompactObject::InitialGridPoints
Definition: CylindricalBinaryCompactObject.hpp:194
CoordinateMap.hpp
Options::Context
Definition: Options.hpp:41
domain::creators::CylindricalBinaryCompactObject::BoundaryConditions
Definition: CylindricalBinaryCompactObject.hpp:214
domain::CoordinateMaps::DiscreteRotation
A CoordinateMap that swaps/negates the coordinate axes.
Definition: DiscreteRotation.hpp:32
domain::CoordinateMaps::Affine
Affine map from .
Definition: Affine.hpp:37
domain::push_back
CoordinateMap< SourceFrame, TargetFrame, Maps..., NewMap > push_back(CoordinateMap< SourceFrame, TargetFrame, Maps... > old_map, NewMap new_map) noexcept
Creates a CoordinateMap by appending the new map to the end of the old maps.
domain::creators::CylindricalBinaryCompactObject::OuterRadius
Definition: CylindricalBinaryCompactObject.hpp:176
domain::creators::CylindricalBinaryCompactObject::RadiusB
Definition: CylindricalBinaryCompactObject.hpp:171
domain::creators::CylindricalBinaryCompactObject::OuterBoundaryCondition
Definition: CylindricalBinaryCompactObject.hpp:227
cstddef
domain::creators::CylindricalBinaryCompactObject::InnerBoundaryCondition
Definition: CylindricalBinaryCompactObject.hpp:218
domain::CoordinateMaps::CylindricalFlatSide
Map from 3D unit right cylindrical shell to a volume that connects a portion of an annulus to a porti...
Definition: CylindricalFlatSide.hpp:89
domain::creators::time_dependence::TimeDependence< 3 >
array
domain::creators::CylindricalBinaryCompactObject::CenterA
Definition: CylindricalBinaryCompactObject.hpp:156
Domain
A wrapper around a vector of Blocks that represent the computational domain.
Definition: Domain.hpp:40
memory
DomainCreator
Base class for creating Domains from an option string.
Definition: DomainCreator.hpp:32
domain::creators::CylindricalBinaryCompactObject::RadiusA
Definition: CylindricalBinaryCompactObject.hpp:166
domain::creators::CylindricalBinaryCompactObject::TimeDependence
Definition: CylindricalBinaryCompactObject.hpp:207
domain::creators::CylindricalBinaryCompactObject::CenterB
Definition: CylindricalBinaryCompactObject.hpp:161
Options::String
const char *const String
The string used in option structs.
Definition: Options.hpp:32
Frame
Definition: IndexType.hpp:36
Frame::Logical
Definition: IndexType.hpp:42
DomainCreator.hpp
domain::CoordinateMaps::CylindricalEndcap
Map from 3D unit right cylinder to a volume that connects portions of two spherical surfaces.
Definition: CylindricalEndcap.hpp:135
domain::creators::CylindricalBinaryCompactObject::InitialRefinement
Definition: CylindricalBinaryCompactObject.hpp:182
std::unique_ptr< domain::creators::time_dependence::TimeDependence< 3 > >
unordered_map
TMPL.hpp
domain::CoordinateMaps::CylindricalFlatEndcap
Map from 3D unit right cylinder to a volume that connects a portion of a circle to a portion of a sph...
Definition: CylindricalFlatEndcap.hpp:82
variant
domain::CoordinateMap
A coordinate map or composition of coordinate maps.
Definition: CoordinateMap.hpp:240
string