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 <optional>
9 : #include <string>
10 : #include <utility>
11 :
12 : #include "DataStructures/ComplexModalVector.hpp"
13 : #include "DataStructures/DataBox/PrefixHelpers.hpp"
14 : #include "DataStructures/DataBox/Prefixes.hpp"
15 : #include "DataStructures/DataBox/Tag.hpp"
16 : #include "DataStructures/DataBox/TagName.hpp"
17 : #include "DataStructures/DataVector.hpp"
18 : #include "DataStructures/Matrix.hpp"
19 : #include "Evolution/Systems/Cce/Tags.hpp"
20 : #include "IO/H5/Dat.hpp"
21 : #include "IO/H5/File.hpp"
22 : #include "IO/H5/Version.hpp"
23 : #include "NumericalAlgorithms/Spectral/SwshTags.hpp"
24 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
25 : #include "Utilities/ErrorHandling/Assert.hpp"
26 : #include "Utilities/Gsl.hpp"
27 : #include "Utilities/Serialization/CharmPupable.hpp"
28 : #include "Utilities/Serialization/PupStlCpp17.hpp"
29 : #include "Utilities/TMPL.hpp"
30 : #include "Utilities/TaggedTuple.hpp"
31 :
32 : namespace Cce {
33 : namespace Tags {
34 : namespace detail {
35 : // tags for use in the buffers for the modal input worldtube data management
36 : // classes
37 : using SpatialMetric = gr::Tags::SpatialMetric<ComplexModalVector, 3>;
38 : using Shift = gr::Tags::Shift<ComplexModalVector, 3>;
39 : using Lapse = gr::Tags::Lapse<ComplexModalVector>;
40 :
41 : // radial derivative prefix tag to be used with the modal input worldtube data
42 : template <typename Tag>
43 : struct Dr : db::SimpleTag, db::PrefixTag {
44 : using type = typename Tag::type;
45 : using tag = Tag;
46 : };
47 :
48 : // tag for the string for accessing the quantity associated with `Tag` in
49 : // worldtube h5 file
50 : template <typename Tag>
51 : struct InputDataSet : db::SimpleTag, db::PrefixTag {
52 : using type = std::string;
53 : using tag = Tag;
54 : };
55 : } // namespace detail
56 : } // namespace Tags
57 :
58 : namespace detail {
59 : // generates the component dataset name in the worldtube file based on the
60 : // tensor indices requested. For instance, if called with arguments ("/g", 0,1),
61 : // it returns the dataset name "/gxy".
62 : template <typename... T>
63 : std::string dataset_name_for_component(std::string base_name,
64 : const T... indices) { // NOLINT
65 : const auto add_index = [&base_name](size_t index) {
66 : ASSERT(index < 3, "The character-arithmetic index must be less than 3.");
67 : base_name += static_cast<char>('x' + index);
68 : };
69 : EXPAND_PACK_LEFT_TO_RIGHT(add_index(indices));
70 : // void cast so that compilers can tell it's used.
71 : (void)add_index;
72 : return base_name;
73 : }
74 :
75 : // creates a pair of indices such that the difference is `2 *
76 : // interpolator_length + pad`, centered around `time`, and bounded by
77 : // `lower_bound` and `upper_bound`. If it cannot be centered, it gives a span
78 : // that is appropriately sized and bounded by the supplied bounds. If the bounds
79 : // are too constraining for the necessary size, it gives a span that is the
80 : // correct size starting at `lower bound`, but not constrained by `upper_bound`
81 : std::pair<size_t, size_t> create_span_for_time_value(
82 : double time, size_t pad, size_t interpolator_length, size_t lower_bound,
83 : size_t upper_bound, const DataVector& time_buffer);
84 :
85 : // retrieves the extraction radius from the specified file.
86 : // We assume that the filename has the extraction radius encoded as an
87 : // integer between the last occurrence of 'R' and the last occurrence of
88 : // '.'. This is the format provided by SpEC.
89 : std::string get_text_radius(const std::string& cce_data_filename);
90 :
91 : // retrieves time stamps and lmax the from the specified file.
92 : void set_time_buffer_and_lmax(gsl::not_null<DataVector*> time_buffer,
93 : size_t& l_max, const h5::Dat& data);
94 :
95 : // retrieves modal data from Bondi or Klein-Gordon worldtube H5 file.
96 : void update_buffer_with_modal_data(
97 : gsl::not_null<ComplexModalVector*> buffer_to_update,
98 : const h5::Dat& read_data, size_t computation_l_max, size_t l_max,
99 : size_t time_span_start, size_t time_span_end, bool is_real);
100 :
101 : // updates `time_span_start` and `time_span_end` based on the provided `time`,
102 : // and inserts the cooresponding modal data (for `InputTags`) from worldtube H5
103 : // file into `buffers`. The function is used by Bondi and Klein-Gordon systems.
104 : template <typename InputTags>
105 : double update_buffers_for_time(
106 : gsl::not_null<Variables<InputTags>*> buffers,
107 : gsl::not_null<size_t*> time_span_start,
108 : gsl::not_null<size_t*> time_span_end, double time, size_t computation_l_max,
109 : size_t l_max, size_t interpolator_length, size_t buffer_depth,
110 : const DataVector& time_buffer,
111 : const tuples::tagged_tuple_from_typelist<
112 : db::wrap_tags_in<Tags::detail::InputDataSet, InputTags>>& dataset_names,
113 : const h5::H5File<h5::AccessType::ReadOnly>& cce_data_file);
114 : } // namespace detail
115 :
116 : /// the full set of tensors to be extracted from the worldtube h5 file
117 1 : using cce_metric_input_tags = tmpl::list<
118 : Tags::detail::SpatialMetric, Tags::detail::Dr<Tags::detail::SpatialMetric>,
119 : ::Tags::dt<Tags::detail::SpatialMetric>, Tags::detail::Shift,
120 : Tags::detail::Dr<Tags::detail::Shift>, ::Tags::dt<Tags::detail::Shift>,
121 : Tags::detail::Lapse, Tags::detail::Dr<Tags::detail::Lapse>,
122 : ::Tags::dt<Tags::detail::Lapse>>;
123 :
124 : /// the full set of tensors to be extracted from the reduced form of the
125 : /// worldtube h5 file
126 1 : using cce_bondi_input_tags =
127 : tmpl::list<Spectral::Swsh::Tags::SwshTransform<Tags::BondiBeta>,
128 : Spectral::Swsh::Tags::SwshTransform<Tags::BondiU>,
129 : Spectral::Swsh::Tags::SwshTransform<Tags::BondiQ>,
130 : Spectral::Swsh::Tags::SwshTransform<Tags::BondiW>,
131 : Spectral::Swsh::Tags::SwshTransform<Tags::BondiJ>,
132 : Spectral::Swsh::Tags::SwshTransform<Tags::Dr<Tags::BondiJ>>,
133 : Spectral::Swsh::Tags::SwshTransform<Tags::Du<Tags::BondiJ>>,
134 : Spectral::Swsh::Tags::SwshTransform<Tags::BondiR>,
135 : Spectral::Swsh::Tags::SwshTransform<Tags::Du<Tags::BondiR>>>;
136 :
137 0 : using klein_gordon_input_tags =
138 : tmpl::list<Spectral::Swsh::Tags::SwshTransform<Tags::KleinGordonPsi>,
139 : Spectral::Swsh::Tags::SwshTransform<Tags::KleinGordonPi>>;
140 :
141 : /// \cond
142 : class MetricWorldtubeH5BufferUpdater;
143 : class BondiWorldtubeH5BufferUpdater;
144 : class KleinGordonWorldtubeH5BufferUpdater;
145 : /// \endcond
146 :
147 : /*!
148 : * \brief Abstract base class for utilities that are able to perform the buffer
149 : * updating procedure needed by the `WorldtubeDataManager`.
150 : *
151 : * \details The methods that are required to be overridden in the derived
152 : * classes are:
153 : * - `WorldtubeBufferUpdater::update_buffers_for_time()`:
154 : * updates the buffers passed by pointer and the `time_span_start` and
155 : * `time_span_end` to be appropriate for the requested `time`,
156 : * `interpolator_length`, and `buffer_depth`.
157 : * - `WorldtubeBufferUpdater::get_clone()`
158 : * clone function to obtain a `std::unique_ptr` of the base
159 : * `WorldtubeBufferUpdater`, needed to pass around the factory-created
160 : * object.
161 : * - `WorldtubeBufferUpdater::time_is_outside_range()`
162 : * the override should return `true` if the `time` could be used in a
163 : * `update_buffers_for_time` call given the data available to the derived
164 : * class, and `false` otherwise
165 : * - `WorldtubeBufferUpdater::get_l_max()`
166 : * The override should return the `l_max` it uses in the
167 : * Goldberg modal data placed in the buffers.
168 : * - `WorldtubeBufferUpdater::get_extraction_radius()`
169 : * The override should return the coordinate radius associated with the modal
170 : * worldtube data that it supplies in the buffer update function. This is
171 : * currently assumed to be a single double, but may be generalized in future
172 : * to be time-dependent.
173 : * - `WorldtubeBufferUpdater::get_time_buffer`
174 : * The override should return the vector of times that it can produce modal
175 : * data at. For instance, if associated with a file input, this will be the
176 : * times at each of the rows of the time-series data.
177 : */
178 : template <typename BufferTags>
179 1 : class WorldtubeBufferUpdater : public PUP::able {
180 : public:
181 0 : using creatable_classes =
182 : tmpl::list<MetricWorldtubeH5BufferUpdater, BondiWorldtubeH5BufferUpdater,
183 : KleinGordonWorldtubeH5BufferUpdater>;
184 :
185 0 : WRAPPED_PUPable_abstract(WorldtubeBufferUpdater); // NOLINT
186 :
187 0 : virtual double update_buffers_for_time(
188 : gsl::not_null<Variables<BufferTags>*> buffers,
189 : gsl::not_null<size_t*> time_span_start,
190 : gsl::not_null<size_t*> time_span_end, double time,
191 : size_t computation_l_max, size_t interpolator_length,
192 : size_t buffer_depth) const = 0;
193 :
194 0 : virtual std::unique_ptr<WorldtubeBufferUpdater> get_clone() const = 0;
195 :
196 0 : virtual bool time_is_outside_range(double time) const = 0;
197 :
198 0 : virtual size_t get_l_max() const = 0;
199 :
200 0 : virtual double get_extraction_radius() const = 0;
201 :
202 0 : virtual bool has_version_history() const = 0;
203 :
204 0 : virtual DataVector& get_time_buffer() = 0;
205 : };
206 :
207 : /// A `WorldtubeBufferUpdater` specialized to the CCE input worldtube H5 file
208 : /// produced by SpEC.
209 1 : class MetricWorldtubeH5BufferUpdater
210 : : public WorldtubeBufferUpdater<cce_metric_input_tags> {
211 : public:
212 : // charm needs the empty constructor
213 0 : MetricWorldtubeH5BufferUpdater() = default;
214 :
215 : /// The constructor takes the filename of the SpEC h5 file that will be used
216 : /// for boundary data. The extraction radius can either be passed in directly,
217 : /// or if it takes the value `std::nullopt`, then the extraction radius is
218 : /// retrieved as an integer in the filename.
219 1 : explicit MetricWorldtubeH5BufferUpdater(
220 : const std::string& cce_data_filename,
221 : std::optional<double> extraction_radius = std::nullopt);
222 :
223 0 : WRAPPED_PUPable_decl_template(MetricWorldtubeH5BufferUpdater); // NOLINT
224 :
225 0 : explicit MetricWorldtubeH5BufferUpdater(CkMigrateMessage* /*unused*/) {}
226 :
227 : /// update the `buffers`, `time_span_start`, and `time_span_end` with
228 : /// time-varies-fastest, Goldberg modal data and the start and end index in
229 : /// the member `time_buffer_` covered by the newly updated `buffers`. The
230 : /// function returns the next time at which a full update will occur. If
231 : /// called again at times earlier than the next full update time, it will
232 : /// leave the `buffers` unchanged and again return the next needed time.
233 1 : double update_buffers_for_time(
234 : gsl::not_null<Variables<cce_metric_input_tags>*> buffers,
235 : gsl::not_null<size_t*> time_span_start,
236 : gsl::not_null<size_t*> time_span_end, double time,
237 : size_t computation_l_max, size_t interpolator_length,
238 : size_t buffer_depth) const override;
239 :
240 0 : std::unique_ptr<WorldtubeBufferUpdater<cce_metric_input_tags>> get_clone()
241 : const override;
242 :
243 : /// The time can only be supported in the buffer update if it is between the
244 : /// first and last time of the input file.
245 1 : bool time_is_outside_range(double time) const override;
246 :
247 : /// retrieves the l_max of the input file
248 1 : size_t get_l_max() const override { return l_max_; }
249 :
250 : /// retrieves the extraction radius
251 1 : double get_extraction_radius() const override { return extraction_radius_; }
252 :
253 : /// The time buffer is supplied by non-const reference to allow views to
254 : /// easily point into the buffer.
255 : ///
256 : /// \warning Altering this buffer outside of the constructor of this class
257 : /// results in undefined behavior! This should be supplied by const reference
258 : /// once there is a convenient method of producing a const view of a vector
259 : /// type.
260 1 : DataVector& get_time_buffer() override { return time_buffer_; }
261 :
262 0 : bool has_version_history() const override { return has_version_history_; }
263 :
264 : /// Serialization for Charm++.
265 1 : void pup(PUP::er& p) override;
266 :
267 : private:
268 0 : void update_buffer(gsl::not_null<ComplexModalVector*> buffer_to_update,
269 : const h5::Dat& read_data, size_t computation_l_max,
270 : size_t time_span_start, size_t time_span_end) const;
271 :
272 0 : bool has_version_history_ = true;
273 0 : double extraction_radius_ = std::numeric_limits<double>::signaling_NaN();
274 0 : size_t l_max_ = 0;
275 :
276 0 : h5::H5File<h5::AccessType::ReadOnly> cce_data_file_;
277 0 : std::string filename_;
278 :
279 : tuples::tagged_tuple_from_typelist<
280 : db::wrap_tags_in<Tags::detail::InputDataSet, cce_metric_input_tags>>
281 0 : dataset_names_;
282 :
283 : // stores all the times in the input file
284 0 : DataVector time_buffer_;
285 : };
286 :
287 : /// A `WorldtubeBufferUpdater` specialized to the CCE input worldtube H5 file
288 : /// produced by the reduced SpEC format.
289 1 : class BondiWorldtubeH5BufferUpdater
290 : : public WorldtubeBufferUpdater<cce_bondi_input_tags> {
291 : public:
292 : // charm needs the empty constructor
293 0 : BondiWorldtubeH5BufferUpdater() = default;
294 :
295 : /// The constructor takes the filename of the SpEC h5 file that will be used
296 : /// for boundary data. The extraction radius can either be passed in directly,
297 : /// or if it takes the value `std::nullopt`, then the extraction radius is
298 : /// retrieved as an integer in the filename.
299 1 : explicit BondiWorldtubeH5BufferUpdater(
300 : const std::string& cce_data_filename,
301 : std::optional<double> extraction_radius = std::nullopt);
302 :
303 0 : WRAPPED_PUPable_decl_template(BondiWorldtubeH5BufferUpdater); // NOLINT
304 :
305 0 : explicit BondiWorldtubeH5BufferUpdater(CkMigrateMessage* /*unused*/) {}
306 :
307 : /// update the `buffers`, `time_span_start`, and `time_span_end` with
308 : /// time-varies-fastest, Goldberg modal data and the start and end index in
309 : /// the member `time_buffer_` covered by the newly updated `buffers`.
310 1 : double update_buffers_for_time(
311 : gsl::not_null<Variables<cce_bondi_input_tags>*> buffers,
312 : gsl::not_null<size_t*> time_span_start,
313 : gsl::not_null<size_t*> time_span_end, double time,
314 : size_t computation_l_max, size_t interpolator_length,
315 : size_t buffer_depth) const override;
316 :
317 0 : std::unique_ptr<WorldtubeBufferUpdater<cce_bondi_input_tags>> get_clone()
318 : const override {
319 : return std::make_unique<BondiWorldtubeH5BufferUpdater>(filename_);
320 : }
321 :
322 : /// The time can only be supported in the buffer update if it is between the
323 : /// first and last time of the input file.
324 1 : bool time_is_outside_range(const double time) const override {
325 : return time < time_buffer_[0] or
326 : time > time_buffer_[time_buffer_.size() - 1];
327 : }
328 :
329 : /// retrieves the l_max of the input file
330 1 : size_t get_l_max() const override { return l_max_; }
331 :
332 : /// retrieves the extraction radius. In most normal circumstances, this will
333 : /// not be needed for Bondi data.
334 1 : double get_extraction_radius() const override {
335 : if (not extraction_radius_.has_value()) {
336 : ERROR(
337 : "Extraction radius has not been set, and was not successfully parsed "
338 : "from the filename. The extraction radius has been used, so must be "
339 : "set either by the input file or via the filename.");
340 : }
341 : return *extraction_radius_;
342 : }
343 :
344 : /// The time buffer is supplied by non-const reference to allow views to
345 : /// easily point into the buffer.
346 : ///
347 : /// \warning Altering this buffer outside of the constructor of this class
348 : /// results in undefined behavior! This should be supplied by const reference
349 : /// once there is a convenient method of producing a const view of a vector
350 : /// type.
351 1 : DataVector& get_time_buffer() override { return time_buffer_; }
352 :
353 0 : bool has_version_history() const override { return true; }
354 :
355 : /// Serialization for Charm++.
356 1 : void pup(PUP::er& p) override;
357 :
358 : private:
359 0 : void update_buffer(gsl::not_null<ComplexModalVector*> buffer_to_update,
360 : const h5::Dat& read_data, size_t computation_l_max,
361 : size_t time_span_start, size_t time_span_end,
362 : bool is_real) const;
363 :
364 0 : std::optional<double> extraction_radius_ = std::nullopt;
365 0 : size_t l_max_ = 0;
366 :
367 0 : h5::H5File<h5::AccessType::ReadOnly> cce_data_file_;
368 0 : std::string filename_;
369 :
370 : tuples::tagged_tuple_from_typelist<
371 : db::wrap_tags_in<Tags::detail::InputDataSet, cce_bondi_input_tags>>
372 0 : dataset_names_;
373 :
374 : // stores all the times in the input file
375 0 : DataVector time_buffer_;
376 : };
377 :
378 : /// A `WorldtubeBufferUpdater` specialized to the Klein-Gordon input worldtube
379 : /// H5 file produced by the SpEC format. We assume the scalar field is
380 : /// real-valued.
381 1 : class KleinGordonWorldtubeH5BufferUpdater
382 : : public WorldtubeBufferUpdater<klein_gordon_input_tags> {
383 : public:
384 : // charm needs the empty constructor
385 0 : KleinGordonWorldtubeH5BufferUpdater() = default;
386 :
387 : /// The constructor takes the filename of the SpEC h5 file that will be used
388 : /// for boundary data. The extraction radius can either be passed in directly,
389 : /// or if it takes the value `std::nullopt`, then the extraction radius is
390 : /// retrieved as an integer in the filename.
391 1 : explicit KleinGordonWorldtubeH5BufferUpdater(
392 : const std::string& cce_data_filename,
393 : std::optional<double> extraction_radius = std::nullopt);
394 :
395 0 : WRAPPED_PUPable_decl_template(KleinGordonWorldtubeH5BufferUpdater); // NOLINT
396 :
397 0 : explicit KleinGordonWorldtubeH5BufferUpdater(CkMigrateMessage* /*unused*/) {}
398 :
399 : /// update the `buffers`, `time_span_start`, and `time_span_end` with
400 : /// time-varies-fastest, Goldberg modal data and the start and end index in
401 : /// the member `time_buffer_` covered by the newly updated `buffers`.
402 1 : double update_buffers_for_time(
403 : gsl::not_null<Variables<klein_gordon_input_tags>*> buffers,
404 : gsl::not_null<size_t*> time_span_start,
405 : gsl::not_null<size_t*> time_span_end, double time,
406 : size_t computation_l_max, size_t interpolator_length,
407 : size_t buffer_depth) const override;
408 :
409 0 : std::unique_ptr<WorldtubeBufferUpdater<klein_gordon_input_tags>> get_clone()
410 : const override {
411 : return std::make_unique<KleinGordonWorldtubeH5BufferUpdater>(filename_);
412 : }
413 :
414 : /// The time can only be supported in the buffer update if it is between the
415 : /// first and last time of the input file.
416 1 : bool time_is_outside_range(const double time) const override {
417 : return time < time_buffer_[0] or
418 : time > time_buffer_[time_buffer_.size() - 1];
419 : }
420 :
421 : /// retrieves the l_max of the input file
422 1 : size_t get_l_max() const override { return l_max_; }
423 :
424 : /// retrieves the extraction radius. In most normal circumstances, this will
425 : /// not be needed for Klein-Gordon data.
426 1 : double get_extraction_radius() const override {
427 : if (not static_cast<bool>(extraction_radius_)) {
428 : ERROR(
429 : "Extraction radius has not been set, and was not successfully parsed "
430 : "from the filename. The extraction radius has been used, so must be "
431 : "set either by the input file or via the filename.");
432 : }
433 : return *extraction_radius_;
434 : }
435 :
436 : /// The time buffer is supplied by non-const reference to allow views to
437 : /// easily point into the buffer.
438 : ///
439 : /// \warning Altering this buffer outside of the constructor of this class
440 : /// results in undefined behavior! This should be supplied by const reference
441 : /// once there is a convenient method of producing a const view of a vector
442 : /// type.
443 1 : DataVector& get_time_buffer() override { return time_buffer_; }
444 :
445 0 : bool has_version_history() const override { return true; }
446 :
447 : /// Serialization for Charm++.
448 1 : void pup(PUP::er& p) override;
449 :
450 : private:
451 : // The scalar field is assumed to be real-valued.
452 0 : void update_buffer(gsl::not_null<ComplexModalVector*> buffer_to_update,
453 : const h5::Dat& read_data, size_t computation_l_max,
454 : size_t time_span_start, size_t time_span_end) const;
455 :
456 0 : std::optional<double> extraction_radius_ = std::nullopt;
457 0 : size_t l_max_ = 0;
458 :
459 0 : h5::H5File<h5::AccessType::ReadOnly> cce_data_file_;
460 0 : std::string filename_;
461 :
462 : tuples::tagged_tuple_from_typelist<
463 : db::wrap_tags_in<Tags::detail::InputDataSet, klein_gordon_input_tags>>
464 0 : dataset_names_;
465 :
466 : // stores all the times in the input file
467 0 : DataVector time_buffer_;
468 : };
469 : } // namespace Cce
|