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 <utility>
8 :
9 : #include "DataStructures/DataBox/DataBox.hpp"
10 : #include "DataStructures/DataBox/Tag.hpp"
11 : #include "DataStructures/Tags/TempTensor.hpp"
12 : #include "DataStructures/Variables.hpp"
13 : #include "IO/Logging/Tags.hpp"
14 : #include "IO/Logging/Verbosity.hpp"
15 : #include "NumericalAlgorithms/SphericalHarmonics/Strahlkorper.hpp"
16 : #include "NumericalAlgorithms/SphericalHarmonics/Tags.hpp"
17 : #include "Options/String.hpp"
18 : #include "Parallel/GlobalCache.hpp"
19 : #include "ParallelAlgorithms/ApparentHorizonFinder/FastFlow.hpp"
20 : #include "ParallelAlgorithms/ApparentHorizonFinder/Tags.hpp"
21 : #include "ParallelAlgorithms/Initialization/MutateAssign.hpp"
22 : #include "ParallelAlgorithms/Interpolation/Protocols/ComputeTargetPoints.hpp"
23 : #include "ParallelAlgorithms/Interpolation/Tags.hpp"
24 : #include "PointwiseFunctions/GeneralRelativity/Surfaces/Tags.hpp"
25 : #include "Utilities/PrettyType.hpp"
26 : #include "Utilities/ProtocolHelpers.hpp"
27 : #include "Utilities/Requires.hpp"
28 : #include "Utilities/TMPL.hpp"
29 : #include "Utilities/TaggedTuple.hpp"
30 :
31 : /// \cond
32 : class DataVector;
33 :
34 : namespace PUP {
35 : class er;
36 : } // namespace PUP
37 : namespace db {
38 : template <typename TagsList>
39 : class DataBox;
40 : } // namespace db
41 : namespace intrp {
42 : namespace Tags {
43 : template <typename TemporalId>
44 : struct TemporalIds;
45 : } // namespace Tags
46 : } // namespace intrp
47 : namespace Tags {
48 : struct Verbosity;
49 : } // namespace Tags
50 : /// \endcond
51 :
52 : namespace intrp {
53 :
54 0 : namespace OptionHolders {
55 : /// Options for finding an apparent horizon.
56 : template <typename Frame>
57 1 : struct ApparentHorizon {
58 : /// See Strahlkorper for suboptions.
59 1 : struct InitialGuess {
60 0 : static constexpr Options::String help = {"Initial guess"};
61 0 : using type = ylm::Strahlkorper<Frame>;
62 : };
63 : /// See ::FastFlow for suboptions.
64 1 : struct FastFlow {
65 0 : static constexpr Options::String help = {"FastFlow options"};
66 0 : using type = ::FastFlow;
67 : };
68 0 : struct Verbosity {
69 0 : static constexpr Options::String help = {"Verbosity"};
70 0 : using type = ::Verbosity;
71 : };
72 0 : using options = tmpl::list<InitialGuess, FastFlow, Verbosity>;
73 0 : static constexpr Options::String help = {
74 : "Provide an initial guess for the apparent horizon surface\n"
75 : "(Strahlkorper) and apparent-horizon-finding-algorithm (FastFlow)\n"
76 : "options."};
77 :
78 0 : ApparentHorizon(ylm::Strahlkorper<Frame> initial_guess_in,
79 : ::FastFlow fast_flow_in, ::Verbosity verbosity_in);
80 :
81 0 : ApparentHorizon() = default;
82 0 : ApparentHorizon(const ApparentHorizon& /*rhs*/) = default;
83 0 : ApparentHorizon& operator=(const ApparentHorizon& /*rhs*/) = delete;
84 0 : ApparentHorizon(ApparentHorizon&& /*rhs*/) = default;
85 0 : ApparentHorizon& operator=(ApparentHorizon&& /*rhs*/) = default;
86 0 : ~ApparentHorizon() = default;
87 :
88 : // NOLINTNEXTLINE(google-runtime-references)
89 0 : void pup(PUP::er& p);
90 :
91 0 : ylm::Strahlkorper<Frame> initial_guess{};
92 0 : ::FastFlow fast_flow{};
93 0 : ::Verbosity verbosity{::Verbosity::Quiet};
94 : };
95 :
96 : template <typename Frame>
97 0 : bool operator==(const ApparentHorizon<Frame>& lhs,
98 : const ApparentHorizon<Frame>& rhs);
99 : template <typename Frame>
100 0 : bool operator!=(const ApparentHorizon<Frame>& lhs,
101 : const ApparentHorizon<Frame>& rhs);
102 :
103 : } // namespace OptionHolders
104 :
105 0 : namespace OptionTags {
106 0 : struct ApparentHorizons {
107 0 : static constexpr Options::String help{"Options for apparent horizon finders"};
108 : };
109 :
110 : template <typename InterpolationTargetTag, typename Frame>
111 0 : struct ApparentHorizon {
112 0 : using type = OptionHolders::ApparentHorizon<Frame>;
113 0 : static constexpr Options::String help{
114 : "Options for interpolation onto apparent horizon."};
115 0 : static std::string name() {
116 : return pretty_type::name<InterpolationTargetTag>();
117 : }
118 0 : using group = ApparentHorizons;
119 : };
120 : } // namespace OptionTags
121 :
122 1 : namespace Tags {
123 : template <typename InterpolationTargetTag, typename Frame>
124 0 : struct ApparentHorizon : db::SimpleTag {
125 0 : using type = OptionHolders::ApparentHorizon<Frame>;
126 0 : using option_tags =
127 : tmpl::list<OptionTags::ApparentHorizon<InterpolationTargetTag, Frame>>;
128 :
129 0 : static constexpr bool pass_metavariables = false;
130 0 : static type create_from_options(const type& option) { return option; }
131 : };
132 : } // namespace Tags
133 :
134 0 : namespace TargetPoints {
135 : /// \brief Computes points on a trial apparent horizon`.
136 : ///
137 : /// This differs from `KerrHorizon` in the following ways:
138 : /// - It supplies points on a prolonged Strahlkorper, at a higher resolution
139 : /// than the Strahlkorper in the DataBox, as needed for horizon finding.
140 : /// - It uses a `FastFlow` in the DataBox.
141 : /// - It has different options (including those for `FastFlow`).
142 : ///
143 : /// Conforms to the intrp::protocols::ComputeTargetPoints protocol
144 : ///
145 : /// For requirements on InterpolationTargetTag, see
146 : /// intrp::protocols::InterpolationTargetTag
147 : template <typename InterpolationTargetTag, typename Frame>
148 1 : struct ApparentHorizon : tt::ConformsTo<intrp::protocols::ComputeTargetPoints> {
149 0 : using const_global_cache_tags =
150 : tmpl::list<Tags::ApparentHorizon<InterpolationTargetTag, Frame>>;
151 0 : using is_sequential = std::true_type;
152 0 : using frame = Frame;
153 :
154 0 : using common_tags =
155 : tmpl::push_back<ylm::Tags::items_tags<Frame>, ::ah::Tags::FastFlow,
156 : logging::Tags::Verbosity<InterpolationTargetTag>,
157 : ylm::Tags::PreviousStrahlkorpers<Frame>>;
158 0 : using simple_tags = tmpl::append<
159 : common_tags,
160 : tmpl::conditional_t<
161 : std::is_same_v<Frame, ::Frame::Inertial>, tmpl::list<>,
162 : tmpl::list<ylm::Tags::CartesianCoords<::Frame::Inertial>>>>;
163 0 : using compute_tags =
164 : tmpl::append<typename ylm::Tags::compute_items_tags<Frame>,
165 : ylm::Tags::TimeDerivStrahlkorperCompute<Frame>>;
166 :
167 : template <typename DbTags, typename Metavariables>
168 0 : static void initialize(const gsl::not_null<db::DataBox<DbTags>*> box,
169 : const Parallel::GlobalCache<Metavariables>& cache) {
170 : const auto& options =
171 : Parallel::get<Tags::ApparentHorizon<InterpolationTargetTag, Frame>>(
172 : cache);
173 :
174 : // Put Strahlkorper and its ComputeItems, FastFlow, and verbosity
175 : // into a new DataBox. The first element of PreviousStrahlkorpers
176 : // is initialized to (time=NaN, strahlkorper=options.initial_guess).
177 : // The NaN is a sentinel value which indicates that the
178 : // PreviousStrahlkorper has not been computed but is instead the
179 : // supplied initial guess. Note that the NaN must be quiet_NaN,
180 : // so we can test for it later without generating an FPE.
181 : //
182 : // Note that if frame is not inertial,
183 : // ylm::Tags::Strahlkorper<::Frame::Inertial> is already
184 : // default initialized so there is no need to do anything special
185 : // here for ylm::Tags::Strahlkorper<::Frame::Inertial>.
186 : Initialization::mutate_assign<common_tags>(
187 : box, options.initial_guess, options.fast_flow, options.verbosity,
188 : std::deque<std::pair<double, ylm::Strahlkorper<Frame>>>{std::make_pair(
189 : std::numeric_limits<double>::quiet_NaN(), options.initial_guess)});
190 : }
191 :
192 : template <typename Metavariables, typename DbTags, typename TemporalId>
193 0 : static tnsr::I<DataVector, 3, Frame> points(
194 : const db::DataBox<DbTags>& box,
195 : const tmpl::type_<Metavariables>& /*meta*/,
196 : const TemporalId& /*temporal_id*/) {
197 : const auto& fast_flow = db::get<::ah::Tags::FastFlow>(box);
198 : const auto& strahlkorper = db::get<ylm::Tags::Strahlkorper<Frame>>(box);
199 :
200 : const size_t L_mesh = fast_flow.current_l_mesh(strahlkorper);
201 : const auto prolonged_strahlkorper =
202 : ylm::Strahlkorper<Frame>(L_mesh, L_mesh, strahlkorper);
203 :
204 : Variables<tmpl::list<::Tags::Tempi<0, 2, ::Frame::Spherical<Frame>>,
205 : ::Tags::Tempi<1, 3, Frame>, ::Tags::TempScalar<2>>>
206 : temp_buffer(prolonged_strahlkorper.ylm_spherepack().physical_size());
207 :
208 : auto& theta_phi =
209 : get<::Tags::Tempi<0, 2, ::Frame::Spherical<Frame>>>(temp_buffer);
210 : auto& r_hat = get<::Tags::Tempi<1, 3, Frame>>(temp_buffer);
211 : auto& radius = get<::Tags::TempScalar<2>>(temp_buffer);
212 : ylm::Tags::ThetaPhiCompute<Frame>::function(make_not_null(&theta_phi),
213 : prolonged_strahlkorper);
214 : ylm::Tags::RhatCompute<Frame>::function(make_not_null(&r_hat), theta_phi);
215 : ylm::Tags::RadiusCompute<Frame>::function(make_not_null(&radius),
216 : prolonged_strahlkorper);
217 :
218 : tnsr::I<DataVector, 3, Frame> prolonged_coords{};
219 : ylm::Tags::CartesianCoordsCompute<Frame>::function(
220 : make_not_null(&prolonged_coords), prolonged_strahlkorper, radius,
221 : r_hat);
222 :
223 : return prolonged_coords;
224 : }
225 : };
226 :
227 : } // namespace TargetPoints
228 : } // namespace intrp
|
