Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <cstddef>
7 : #include <memory>
8 : #include <string>
9 :
10 : #include "DataStructures/SpinWeighted.hpp"
11 : #include "DataStructures/Tensor/TypeAliases.hpp"
12 : #include "Evolution/Systems/Cce/Initialize/InitializeJ.hpp"
13 : #include "Options/Options.hpp"
14 : #include "Options/String.hpp"
15 : #include "Utilities/Gsl.hpp"
16 : #include "Utilities/Serialization/CharmPupable.hpp"
17 : #include "Utilities/TMPL.hpp"
18 :
19 : /// \cond
20 : class ComplexDataVector;
21 : /// \endcond
22 :
23 : namespace Cce {
24 1 : namespace InitializeJ {
25 :
26 : /// Possible iteration heuristics to use for optimizing the value of the
27 : /// conformal factor \f$\omega\f$ to fix the initial data.
28 1 : enum class ConformalFactorIterationHeuristic {
29 : /// Assumes that the spin-weighted Jacobian perturbations obey
30 : /// \f$c = \hat \eth f\f$,\f$d = \hat{\bar\eth} f\f$, for some spin-weight-1
31 : /// value\f$f\f$.
32 : SpinWeight1CoordPerturbation,
33 : /// Varies only the \f$d\f$ spin-weighted Jacobian when constructing the
34 : /// itertion heuristic, leaving \f$c\f$ fixed
35 : OnlyVaryGaugeD
36 : };
37 :
38 : /*!
39 : * \brief Generate initial data that has a conformal factor \f$\omega\f$ chosen
40 : * to compensate for the boundary value of \f$\beta\f$ so that the initial time
41 : * coordinate is approximately inertial at \f$I^+\f$.
42 : *
43 : * \details The core calculation for this initial data choice is the iterative
44 : * optimization of the angular conformal factor \f$\omega\f$ such that it
45 : * cancels some portion of the value of \f$e^{2\beta}\f$ that contributes to the
46 : * definition of the asymptotically inertial time.
47 : * The initial data generation process proceeds slightly differently depending
48 : * on the set of input options that are used:
49 : * - If `UseBetaIntegralEstimate` is false, the conformal factor will be
50 : * optimized to minimize the transformed value of \f$\beta\f$ on the worldtube
51 : * boundary.
52 : * - If `UseBetaIntegralEstimate` is true, the conformal factor will be
53 : * optimized to minimize an estimate of the asymptotic value of \f$\beta\f$ in
54 : * the evolved coordinates.
55 : * - `OptimizeL0Mode` indicates whether the \f$l=0\f$ mode of the conformal
56 : * factor shoudl be included in the optimization. This option is useful because
57 : * the optimization can usually find a better solution when the \f$l=0\f$ mode
58 : * is ignored, and the \f$l=0\f$ should not contribute significantly to the
59 : * resulting waveform.
60 : *
61 : * - If `UseInputModes` is false, the \f$J\f$ value on the initial hypersurface
62 : * will be set by an \f$A/r + B/r^3\f$ ansatz, chosen to match the worldtube
63 : * boundary value of \f$J\f$ and \f$\partial_r J\f$ in the new coordinates.
64 : * In this case, the alternative arguments `InputModes` or `InputModesFromFile`
65 : * are ignored.
66 : * - If `UseInputModes` is true, the \f$1/r\f$ part of \f$J\f$ will be set to
67 : * spin-weighted spherical harmonic modes specified by either an input h5 file
68 : * (in the case of using the input option `InputModesFromFile`) or from a list
69 : * of complex values specified in the input file (in the case of using the input
70 : * option `InputModes`). Then, the \f$1/r^3\f$ and \f$1/r^4\f$ parts of \f$J\f$
71 : * are chosen to match the boundary value of \f$J\f$ and \f$\partial_r J\f$ on
72 : * the worldtube boundary in the new coordinates.
73 : */
74 1 : struct ConformalFactor : InitializeJ<false> {
75 0 : struct AngularCoordinateTolerance {
76 0 : using type = double;
77 0 : static std::string name() { return "AngularCoordTolerance"; }
78 0 : static constexpr Options::String help = {
79 : "Tolerance of initial angular coordinates for CCE"};
80 0 : static type lower_bound() { return 1.0e-14; }
81 0 : static type upper_bound() { return 1.0e-3; }
82 : };
83 0 : struct MaxIterations {
84 0 : using type = size_t;
85 0 : static constexpr Options::String help = {
86 : "Number of linearized inversion iterations."};
87 0 : static type lower_bound() { return 10; }
88 0 : static type upper_bound() { return 1000; }
89 0 : static type suggested_value() { return 300; }
90 : };
91 0 : struct RequireConvergence {
92 0 : using type = bool;
93 0 : static constexpr Options::String help = {
94 : "If true, initialization will error if it hits MaxIterations"};
95 0 : static type suggested_value() { return true; }
96 : };
97 0 : struct OptimizeL0Mode {
98 0 : using type = bool;
99 0 : static constexpr Options::String help = {
100 : "If true, the average value of the conformal factor will be included "
101 : "during optimization; otherwise it will be omitted (filtered)."};
102 0 : static type suggested_value() { return false; }
103 : };
104 0 : struct UseBetaIntegralEstimate {
105 0 : using type = bool;
106 0 : static constexpr Options::String help = {
107 : "If true, the iterative algorithm will calculate an estimate of the "
108 : "asymptotic beta value using the 1/r part of the initial J."};
109 0 : static type suggested_value() { return true; }
110 : };
111 0 : struct ConformalFactorIterationHeuristic {
112 0 : using type = ::Cce::InitializeJ::ConformalFactorIterationHeuristic;
113 0 : static constexpr Options::String help = {
114 : "The heuristic method used to set the spin-weighted Jacobian factors "
115 : "when iterating to minimize the asymptotic conformal factor."};
116 0 : static type suggested_value() {
117 : return ::Cce::InitializeJ::ConformalFactorIterationHeuristic::
118 : SpinWeight1CoordPerturbation;
119 : }
120 : };
121 0 : struct UseInputModes {
122 0 : using type = bool;
123 0 : static constexpr Options::String help = {
124 : "If true, the 1/r part of J will be set using modes read from the "
125 : "input file, or from a specified h5 file. If false, the inverse cubic "
126 : "scheme will determine the 1/r part of J."};
127 : };
128 0 : struct InputModesFromFile {
129 0 : using type = std::string;
130 0 : static constexpr Options::String help = {
131 : "A filename from which to retrieve a set of modes (from InitialJ.dat) "
132 : "to use to determine the 1/r part of J on the initial hypersurface. "
133 : "The modes are parsed in l-ascending, m-ascending, m-varies-fastest, "
134 : "real then imaginary part order."};
135 : };
136 0 : struct InputModes {
137 0 : using type = std::vector<std::complex<double>>;
138 0 : static constexpr Options::String help = {
139 : "An explicit list of modes to use to set the 1/r part of J on the "
140 : "initial hypersurface. They are parsed in l-ascending, m-ascending, "
141 : "m-varies-fastest order."};
142 : };
143 :
144 0 : using options =
145 : tmpl::list<AngularCoordinateTolerance, MaxIterations, RequireConvergence,
146 : OptimizeL0Mode, UseBetaIntegralEstimate,
147 : ConformalFactorIterationHeuristic, UseInputModes,
148 : Options::Alternatives<tmpl::list<InputModesFromFile>,
149 : tmpl::list<InputModes>>>;
150 0 : static constexpr Options::String help = {
151 : "Generate CCE initial data based on choosing an angular conformal factor "
152 : "based on the value of the CCE scalar beta in an attempt to make the "
153 : "time variable approximately asymptotically inertial"};
154 :
155 0 : WRAPPED_PUPable_decl_template(ConformalFactor); // NOLINT
156 0 : explicit ConformalFactor(CkMigrateMessage* msg);
157 :
158 0 : ConformalFactor() = default;
159 0 : ConformalFactor(
160 : double angular_coordinate_tolerance, size_t max_iterations,
161 : bool require_convergence, bool optimize_l_0_mode,
162 : bool use_beta_integral_estimate,
163 : ::Cce::InitializeJ::ConformalFactorIterationHeuristic iteration_heuristic,
164 : bool use_input_modes, std::string input_mode_filename);
165 :
166 0 : ConformalFactor(
167 : double angular_coordinate_tolerance, size_t max_iterations,
168 : bool require_convergence, bool optimize_l_0_mode,
169 : bool use_beta_integral_estimate,
170 : ::Cce::InitializeJ::ConformalFactorIterationHeuristic iteration_heuristic,
171 : bool use_input_modes, std::vector<std::complex<double>> input_modes);
172 :
173 0 : std::unique_ptr<InitializeJ> get_clone() const override;
174 :
175 0 : void operator()(
176 : gsl::not_null<Scalar<SpinWeighted<ComplexDataVector, 2>>*> j,
177 : gsl::not_null<tnsr::i<DataVector, 3>*> cartesian_cauchy_coordinates,
178 : gsl::not_null<
179 : tnsr::i<DataVector, 2, ::Frame::Spherical<::Frame::Inertial>>*>
180 : angular_cauchy_coordinates,
181 : const Scalar<SpinWeighted<ComplexDataVector, 2>>& boundary_j,
182 : const Scalar<SpinWeighted<ComplexDataVector, 2>>& boundary_dr_j,
183 : const Scalar<SpinWeighted<ComplexDataVector, 0>>& r,
184 : const Scalar<SpinWeighted<ComplexDataVector, 0>>& beta, size_t l_max,
185 : size_t number_of_radial_points,
186 : gsl::not_null<Parallel::NodeLock*> hdf5_lock) const override;
187 :
188 0 : void pup(PUP::er& p) override;
189 :
190 : private:
191 0 : double angular_coordinate_tolerance_ = 1.0e-11;
192 0 : size_t max_iterations_ = 300;
193 0 : bool require_convergence_ = true;
194 0 : bool optimize_l_0_mode_ = false;
195 0 : bool use_beta_integral_estimate_ = true;
196 0 : ::Cce::InitializeJ::ConformalFactorIterationHeuristic iteration_heuristic_ =
197 : ::Cce::InitializeJ::ConformalFactorIterationHeuristic::
198 : SpinWeight1CoordPerturbation;
199 0 : bool use_input_modes_ = false;
200 0 : std::vector<std::complex<double>> input_modes_;
201 0 : std::optional<std::string> input_mode_filename_;
202 : };
203 :
204 0 : std::ostream& operator<<(
205 : std::ostream& os,
206 : const Cce::InitializeJ::ConformalFactorIterationHeuristic& heuristic_type);
207 :
208 : } // namespace InitializeJ
209 : } // namespace Cce
210 :
211 : template <>
212 0 : struct Options::create_from_yaml<
213 : Cce::InitializeJ::ConformalFactorIterationHeuristic> {
214 : template <typename Metavariables>
215 0 : static Cce::InitializeJ::ConformalFactorIterationHeuristic create(
216 : const Options::Option& options) {
217 : return create<void>(options);
218 : }
219 : };
220 : template <>
221 0 : Cce::InitializeJ::ConformalFactorIterationHeuristic
222 : Options::create_from_yaml<Cce::InitializeJ::ConformalFactorIterationHeuristic>::
223 : create<void>(const Options::Option& options);
|
