Line data Source code
1 1 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : /// \file
5 : /// Defines the class CylindricalFlatSide.
6 :
7 : #pragma once
8 :
9 : #include <array>
10 : #include <cstddef>
11 : #include <limits>
12 : #include <optional>
13 :
14 : #include "DataStructures/Tensor/TypeAliases.hpp"
15 : #include "Domain/CoordinateMaps/FocallyLiftedFlatSide.hpp"
16 : #include "Domain/CoordinateMaps/FocallyLiftedMap.hpp"
17 :
18 : /// \cond
19 : namespace PUP {
20 : class er;
21 : } // namespace PUP
22 : /// \endcond
23 :
24 : namespace domain::CoordinateMaps {
25 :
26 : /*!
27 : * \ingroup CoordinateMapsGroup
28 : *
29 : * \brief Map from 3D unit right cylindrical shell to a volume that connects
30 : * a portion of an annulus to a portion of a spherical surface.
31 : *
32 : * \image html CylindricalFlatSide.svg "A cylinder maps to the shaded region."
33 : *
34 : * \details Consider a 2D annulus in 3D space that is normal to the
35 : * \f$z\f$ axis and has (3D) center \f$C_1\f$, inner radius
36 : * \f$R_\mathrm{in}\f$ and outer radius \f$R_\mathrm{out}\f$
37 : * Also consider a sphere with center \f$C_2\f$, and radius \f$R_2\f$.
38 : * Also let there be a projection point \f$P\f$.
39 : *
40 : * CylindricalFlatSide maps a 3D unit right cylindrical shell (with
41 : * coordinates \f$(\bar{x},\bar{y},\bar{z})\f$ such that
42 : * \f$-1\leq\bar{z}\leq 1\f$ and \f$1 \leq \bar{x}^2+\bar{y}^2 \leq
43 : * 4\f$) to the shaded area in the figure above (with coordinates
44 : * \f$(x,y,z)\f$). The "bottom" of the cylinder \f$\bar{z}=-1\f$ is
45 : * mapped to the interior of the annulus with radii
46 : * \f$R_\mathrm{in}\f$ and \f$R_\mathrm{out}\f$. Curves of constant
47 : * \f$(\bar{x},\bar{y})\f$ are mapped to portions of lines that pass
48 : * through \f$P\f$. Along each of these curves, \f$\bar{z}=-1\f$ is
49 : * mapped to a point inside the annulus and \f$\bar{z}=+1\f$ is mapped to a
50 : * point on the sphere.
51 : *
52 : * CylindricalFlatSide is intended to be composed with Wedge2D maps to
53 : * construct a portion of a cylindrical domain for a binary system.
54 : *
55 : * CylindricalFlatSide is described briefly in the Appendix of
56 : * \cite Buchman:2012dw.
57 : * CylindricalFlatSide is used to construct the blocks labeled 'ME
58 : * cylinder' in Figure 20 of that paper.
59 : *
60 : * CylindricalFlatSide is implemented using `FocallyLiftedMap`
61 : * and `FocallyLiftedInnerMaps::FlatSide`; see those classes for
62 : * details.
63 : *
64 : * ### Restrictions on map parameters.
65 : *
66 : * We demand that:
67 : * - The sphere is at a larger value of \f$z\f$ (plus 5
68 : * percent of the sphere radius) than the plane containing the
69 : * annulus.
70 : * - The projection point \f$z_\mathrm{P}\f$ is
71 : * inside the sphere and more than 15 percent away from the boundary
72 : * of the sphere.
73 : * - The center of the annulus is contained in the circle that results from
74 : * projecting the sphere into the \f$xy\f$ plane.
75 : * - The outer radius of the annulus is larger than 5 percent of the distance
76 : * between the center of the annulus and the projection point.
77 : * - The inner radius of the annulus is less than 95 percent of the outer
78 : * radius, larger than 5 percent of the outer radius, and larger than one
79 : * percent of the distance between the center of the annulus and the
80 : * projection point. The last condition means that the angle subtended by
81 : * the inner radius with respect to the projection point is not too small.
82 : *
83 : * It is possible to construct a valid map without these assumptions,
84 : * but some of these assumptions simplify the code and others eliminate
85 : * edge cases where Jacobians become large or small.
86 : *
87 : */
88 1 : class CylindricalFlatSide {
89 : public:
90 0 : static constexpr size_t dim = 3;
91 0 : CylindricalFlatSide(const std::array<double, 3>& center_one,
92 : const std::array<double, 3>& center_two,
93 : const std::array<double, 3>& proj_center,
94 : const double inner_radius, const double outer_radius,
95 : const double radius_two);
96 :
97 0 : CylindricalFlatSide() = default;
98 0 : ~CylindricalFlatSide() = default;
99 0 : CylindricalFlatSide(CylindricalFlatSide&&) = default;
100 0 : CylindricalFlatSide(const CylindricalFlatSide&) = default;
101 0 : CylindricalFlatSide& operator=(const CylindricalFlatSide&) = default;
102 0 : CylindricalFlatSide& operator=(CylindricalFlatSide&&) = default;
103 :
104 : template <typename T>
105 0 : std::array<T, 3> operator()(const std::array<T, 3>& source_coords) const;
106 :
107 0 : std::optional<std::array<double, 3>> inverse(
108 : const std::array<double, 3>& target_coords) const;
109 :
110 : template <typename T>
111 0 : tnsr::Ij<T, 3, Frame::NoFrame> jacobian(
112 : const std::array<T, 3>& source_coords) const;
113 :
114 : template <typename T>
115 0 : tnsr::Ij<T, 3, Frame::NoFrame> inv_jacobian(
116 : const std::array<T, 3>& source_coords) const;
117 :
118 : // NOLINTNEXTLINE(google-runtime-references)
119 0 : void pup(PUP::er& p);
120 :
121 0 : static bool is_identity() { return false; }
122 :
123 0 : static constexpr bool supports_hessian{false};
124 :
125 : private:
126 0 : friend bool operator==(const CylindricalFlatSide& lhs,
127 : const CylindricalFlatSide& rhs);
128 0 : FocallyLiftedMap<FocallyLiftedInnerMaps::FlatSide> impl_;
129 : };
130 0 : bool operator!=(const CylindricalFlatSide& lhs, const CylindricalFlatSide& rhs);
131 :
132 : } // namespace domain::CoordinateMaps
|
