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 <limits>
9 : #include <optional>
10 : #include <string>
11 : #include <variant>
12 :
13 : #include "DataStructures/DataVector.hpp"
14 : #include "Domain/Creators/TimeDependentOptions/FromVolumeFile.hpp"
15 : #include "Domain/Structure/ObjectLabel.hpp"
16 : #include "NumericalAlgorithms/SphericalHarmonics/IO/ReadSurfaceYlm.hpp"
17 : #include "NumericalAlgorithms/SphericalHarmonics/Strahlkorper.hpp"
18 : #include "Options/Auto.hpp"
19 : #include "Options/Context.hpp"
20 : #include "Options/String.hpp"
21 : #include "Utilities/Gsl.hpp"
22 : #include "Utilities/TMPL.hpp"
23 :
24 : namespace domain::creators::time_dependent_options {
25 : /*!
26 : * \brief Mass and spin necessary for calculating the \f$ Y_{lm} \f$
27 : * coefficients of a Kerr horizon of certain Boyer-Lindquist radius for the
28 : * shape map of the Sphere domain creator.
29 : */
30 1 : struct KerrSchildFromBoyerLindquist {
31 : /// \brief The mass of the Kerr black hole.
32 1 : struct Mass {
33 0 : using type = double;
34 0 : static constexpr Options::String help = {"The mass of the Kerr BH."};
35 : };
36 : /// \brief The dimensionless spin of the Kerr black hole.
37 1 : struct Spin {
38 0 : using type = std::array<double, 3>;
39 0 : static constexpr Options::String help = {
40 : "The dim'less spin of the Kerr BH."};
41 : };
42 :
43 0 : using options = tmpl::list<Mass, Spin>;
44 :
45 0 : static constexpr Options::String help = {
46 : "Conform to an ellipsoid of constant Boyer-Lindquist radius in "
47 : "Kerr-Schild coordinates. This Boyer-Lindquist radius is chosen as the "
48 : "value of the 'InnerRadius'. To conform to the outer Kerr horizon, "
49 : "choose an 'InnerRadius' of r_+ = M + sqrt(M^2-a^2)."};
50 :
51 0 : KerrSchildFromBoyerLindquist();
52 0 : KerrSchildFromBoyerLindquist(double mass_in, std::array<double, 3> spin_in);
53 :
54 0 : double mass{std::numeric_limits<double>::signaling_NaN()};
55 0 : std::array<double, 3> spin{std::numeric_limits<double>::signaling_NaN(),
56 : std::numeric_limits<double>::signaling_NaN(),
57 : std::numeric_limits<double>::signaling_NaN()};
58 : };
59 :
60 : /// Label for shape map options
61 1 : struct Spherical {};
62 :
63 0 : struct YlmsFromFile {
64 0 : struct H5Filename {
65 0 : using type = std::string;
66 0 : static constexpr Options::String help =
67 : "Path to the data file containing the ylm coefficients and their "
68 : "derivatives.";
69 : };
70 :
71 0 : struct SubfileNames {
72 0 : using type = std::vector<std::string>;
73 0 : static constexpr Options::String help =
74 : "Subfile names for the different order derivatives of the ylm "
75 : "coefficients. You must specify the subfile name for the ylm "
76 : "coefficients themselves, and can optionally specify the subfile name "
77 : "for the first and second time derivatives as well, in that order. If "
78 : "you don't specify a derivative subfile, those coefficients will be "
79 : "defaulted to zero.";
80 0 : static size_t lower_bound_on_size() { return 1; }
81 0 : static size_t upper_bound_on_size() { return 3; }
82 : };
83 :
84 0 : struct MatchTime {
85 0 : using type = double;
86 0 : static constexpr Options::String help =
87 : "Time in the H5File to get the coefficients at. Will likely be the "
88 : "same as the initial time";
89 : };
90 :
91 0 : struct MatchTimeEpsilon {
92 0 : using type = Options::Auto<double>;
93 0 : static constexpr Options::String help =
94 : "Look for times in the H5File within this epsilon of the match time. "
95 : "This is to avoid having to know the exact time to all digits. Default "
96 : "is 1e-12.";
97 : };
98 :
99 0 : struct SetL1CoefsToZero {
100 0 : using type = bool;
101 0 : static constexpr Options::String help =
102 : "Whether to set the L=1 coefs to zero or not. This may be desirable "
103 : "because L=1 is degenerate with a translation of the BH.";
104 : };
105 :
106 0 : struct CheckFrame {
107 0 : using type = bool;
108 0 : static constexpr Options::String help =
109 : "Whether to check if the frame of the Strahlkorper in the file matches "
110 : "the Distorted frame.";
111 : };
112 :
113 0 : using options = tmpl::list<H5Filename, SubfileNames, MatchTime,
114 : MatchTimeEpsilon, SetL1CoefsToZero, CheckFrame>;
115 :
116 0 : static constexpr Options::String help = {
117 : "Read the Y_lm coefficients of a Strahlkorper from file and use those to "
118 : "initialize the coefficients of a shape map."};
119 0 : YlmsFromFile();
120 0 : YlmsFromFile(std::string h5_filename_in,
121 : std::vector<std::string> subfile_names_in, double match_time_in,
122 : std::optional<double> match_time_epsilon_in,
123 : bool set_l1_coefs_to_zero_in, bool check_frame_in = true);
124 :
125 0 : std::string h5_filename{};
126 0 : std::vector<std::string> subfile_names{};
127 0 : double match_time{};
128 0 : std::optional<double> match_time_epsilon{};
129 0 : bool set_l1_coefs_to_zero{};
130 0 : bool check_frame{true};
131 : };
132 :
133 0 : struct YlmsFromSpEC {
134 0 : struct DatFilename {
135 0 : using type = std::string;
136 0 : static constexpr Options::String help =
137 : "Name of the SpEC dat file holding the coefficients. Note that this "
138 : "isn't a Dat file within an H5 file. This must be an actual '.dat' "
139 : "file on disk.";
140 : };
141 :
142 0 : struct MatchTime {
143 0 : using type = double;
144 0 : static constexpr Options::String help =
145 : "Time in the H5File to get the coefficients at. Will likely be the "
146 : "same as the initial time";
147 : };
148 :
149 0 : struct MatchTimeEpsilon {
150 0 : using type = Options::Auto<double>;
151 0 : static constexpr Options::String help =
152 : "Look for times in the H5File within this epsilon of the match time. "
153 : "This is to avoid having to know the exact time to all digits. Default "
154 : "is 1e-12.";
155 : };
156 :
157 0 : struct SetL1CoefsToZero {
158 0 : using type = bool;
159 0 : static constexpr Options::String help =
160 : "Whether to set the L=1 coefs to zero or not. This may be desirable "
161 : "because L=1 is degenerate with a translation of the BH.";
162 : };
163 :
164 0 : using options =
165 : tmpl::list<DatFilename, MatchTime, MatchTimeEpsilon, SetL1CoefsToZero>;
166 :
167 0 : static constexpr Options::String help = {
168 : "Read the Y_lm coefficients of a Strahlkorper from file and use those to "
169 : "initialize the coefficients of a shape map."};
170 0 : YlmsFromSpEC();
171 0 : YlmsFromSpEC(std::string dat_filename_in, double match_time_in,
172 : std::optional<double> match_time_epsilon_in,
173 : bool set_l1_coefs_to_zero_in);
174 :
175 0 : std::string dat_filename{};
176 0 : double match_time{};
177 0 : std::optional<double> match_time_epsilon{};
178 0 : bool set_l1_coefs_to_zero{};
179 : };
180 :
181 : /*!
182 : * \brief Class to be used as an option for initializing shape map coefficients.
183 : *
184 : * \tparam IncludeTransitionEndsAtCube This is mainly added for the
185 : * `domain::creators::BinaryCompactObject` domain.
186 : * \tparam Object Which object that this shape map represents. Use
187 : * `domain::ObjectLabel::None` if there is only a single object in your
188 : * simulation.
189 : */
190 : template <bool IncludeTransitionEndsAtCube, domain::ObjectLabel Object>
191 1 : struct ShapeMapOptions {
192 0 : using type = Options::Auto<ShapeMapOptions, Options::AutoLabel::None>;
193 0 : static std::string name() { return "ShapeMap" + get_output(Object); }
194 0 : static constexpr Options::String help = {
195 : "Options for a time-dependent distortion (shape) map about the "
196 : "specified object. Specify 'None' to not use this map."};
197 :
198 0 : struct LMax {
199 0 : using type = size_t;
200 0 : static constexpr Options::String help = {
201 : "LMax used for the number of spherical harmonic coefficients of the "
202 : "distortion map."};
203 : };
204 :
205 0 : struct InitialValues {
206 0 : using type = Options::Auto<
207 : std::variant<KerrSchildFromBoyerLindquist, YlmsFromFile, YlmsFromSpEC,
208 : FromVolumeFile<names::ShapeSize<Object>>>,
209 : Spherical>;
210 0 : static constexpr Options::String help = {
211 : "Initial Ylm coefficients for the shape map. Specify 'Spherical' for "
212 : "all coefficients to be initialized to zero."};
213 : };
214 :
215 0 : struct SizeInitialValues {
216 0 : using type = Options::Auto<std::array<double, 3>>;
217 0 : static constexpr Options::String help = {
218 : "Initial value and two derivatives of the 00 coefficient. Specify "
219 : "'Auto' to use the 00 coefficient specified in the 'InitialValues' "
220 : "option."};
221 : };
222 :
223 0 : struct TransitionEndsAtCube {
224 0 : using type = bool;
225 0 : static constexpr Options::String help = {
226 : "If 'true', the shape map transition function will be 0 at the cubical "
227 : "boundary around the object. If 'false' the transition function will "
228 : "be 0 at the outer radius of the inner sphere around the object"};
229 : };
230 :
231 0 : using common_options = tmpl::list<LMax, InitialValues, SizeInitialValues>;
232 :
233 0 : using options =
234 : tmpl::conditional_t<IncludeTransitionEndsAtCube,
235 : tmpl::push_back<common_options, TransitionEndsAtCube>,
236 : common_options>;
237 0 : ShapeMapOptions() = default;
238 0 : ShapeMapOptions(
239 : size_t l_max_in,
240 : std::optional<
241 : std::variant<KerrSchildFromBoyerLindquist, YlmsFromFile, YlmsFromSpEC,
242 : FromVolumeFile<names::ShapeSize<Object>>>>
243 : initial_values_in,
244 : std::optional<std::array<double, 3>> initial_size_values_in =
245 : std::nullopt,
246 : bool transition_ends_at_cube_in = false)
247 : : l_max(l_max_in),
248 : initial_values(std::move(initial_values_in)),
249 : initial_size_values(initial_size_values_in),
250 : transition_ends_at_cube(transition_ends_at_cube_in) {}
251 :
252 0 : size_t l_max{};
253 : std::optional<
254 : std::variant<KerrSchildFromBoyerLindquist, YlmsFromFile, YlmsFromSpEC,
255 : FromVolumeFile<names::ShapeSize<Object>>>>
256 0 : initial_values{};
257 0 : std::optional<std::array<double, 3>> initial_size_values{};
258 0 : bool transition_ends_at_cube{false};
259 : };
260 :
261 : template <bool IncludeTransitionEndsAtCube, domain::ObjectLabel Object>
262 : std::pair<std::array<DataVector, 3>, std::array<DataVector, 4>>
263 0 : initial_shape_and_size_funcs(
264 : const ShapeMapOptions<IncludeTransitionEndsAtCube, Object>& shape_options,
265 : double inner_radius);
266 : } // namespace domain::creators::time_dependent_options
|