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 <optional>
8 : #include <pup.h>
9 : #include <string>
10 : #include <unordered_map>
11 : #include <utility>
12 : #include <vector>
13 :
14 : #include "DataStructures/DataBox/DataBox.hpp"
15 : #include "DataStructures/DataBox/ObservationBox.hpp"
16 : #include "DataStructures/DataBox/TagName.hpp"
17 : #include "DataStructures/DataVector.hpp"
18 : #include "DataStructures/Tensor/Tensor.hpp"
19 : #include "Domain/Structure/ElementId.hpp"
20 : #include "Domain/Tags.hpp"
21 : #include "IO/Observer/GetSectionObservationKey.hpp"
22 : #include "IO/Observer/Helpers.hpp"
23 : #include "IO/Observer/ObservationId.hpp"
24 : #include "IO/Observer/ObserverComponent.hpp"
25 : #include "IO/Observer/ReductionActions.hpp"
26 : #include "IO/Observer/TypeOfObservation.hpp"
27 : #include "NumericalAlgorithms/LinearOperators/DefiniteIntegral.hpp"
28 : #include "NumericalAlgorithms/Spectral/Mesh.hpp"
29 : #include "NumericalAlgorithms/Spectral/Quadrature.hpp"
30 : #include "Options/String.hpp"
31 : #include "Parallel/ArrayIndex.hpp"
32 : #include "Parallel/GlobalCache.hpp"
33 : #include "Parallel/Invoke.hpp"
34 : #include "Parallel/Local.hpp"
35 : #include "Parallel/Reduction.hpp"
36 : #include "Parallel/TypeTraits.hpp"
37 : #include "ParallelAlgorithms/Events/Tags.hpp"
38 : #include "ParallelAlgorithms/EventsAndTriggers/Event.hpp"
39 : #include "Utilities/ErrorHandling/Assert.hpp"
40 : #include "Utilities/ErrorHandling/Error.hpp"
41 : #include "Utilities/Functional.hpp"
42 : #include "Utilities/OptionalHelpers.hpp"
43 : #include "Utilities/PrettyType.hpp"
44 : #include "Utilities/Serialization/CharmPupable.hpp"
45 : #include "Utilities/TMPL.hpp"
46 :
47 : namespace Events {
48 : /// @{
49 : /*!
50 : * \brief Compute norms of tensors in the DataBox and write them to disk.
51 : *
52 : * The L1 norm is computed as the mean absolute value, so
53 : *
54 : * \f{align*}{
55 : * L_1(u)=\frac{1}{N}\sum_{i=0}^{N-1} |u_i|
56 : * \f}
57 : *
58 : * where \f$N\f$ is the number of grid points.
59 : *
60 : * The L2 norm is computed as the RMS, so
61 : *
62 : * \f{align*}{
63 : * L_2(u)=\sqrt{\frac{1}{N}\sum_{i=0}^{N-1} u_i^2}
64 : * \f}
65 : *
66 : * The norm can be taken for each individual component, or summed over
67 : * components. For the max/min it is then the max/min over all components, while
68 : * for the L1 norm we have (for a 3d vector, 2d and 1d are similar)
69 : *
70 : * \f{align*}{
71 : * L_1(v^k)=\frac{1}{N}\sum_{i=0}^{N-1} \left[|v^x_i| + |v^y_i|
72 : * + |v^z_i|\right]
73 : * \f}
74 : *
75 : * and for the L2 norm
76 : *
77 : * \f{align*}{
78 : * L_2(v^k)=\sqrt{\frac{1}{N}\sum_{i=0}^{N-1} \left[(v^x_i)^2 + (v^y_i)^2
79 : * + (v^z_i)^2\right]}
80 : * \f}
81 : *
82 : * The L1 integral norm is:
83 : *
84 : * \begin{equation}
85 : * L_{1,\mathrm{int}}(v^k) = \frac{1}{V}\int_\Omega \left[
86 : * |v^x_i| + |v^y_i| + |v^z_i|\right] \mathrm{d}V
87 : * \end{equation}
88 : *
89 : * The L2 integral norm is:
90 : *
91 : * \begin{equation}
92 : * L_{2,\mathrm{int}}(v^k) = \sqrt{\frac{1}{V}\int_\Omega \left[
93 : * (v^x_i)^2 + (v^y_i)^2 + (v^z_i)^2\right] \mathrm{d}V}
94 : * \end{equation}
95 : *
96 : * where $V=\int_\Omega$ is the volume of the entire domain in inertial
97 : * coordinates.
98 : *
99 : * VolumeIntegral only computes the volume integral without any normalization.
100 : *
101 : * Here is an example of an input file:
102 : *
103 : * \snippet Test_ObserveNorms.cpp input_file_examples
104 : *
105 : * \note The `NonTensorComputeTags` are intended to be used for `Variables`
106 : * compute tags like `Tags::DerivCompute`
107 : *
108 : * \par Array sections
109 : * This event supports sections (see `Parallel::Section`). Set the
110 : * `ArraySectionIdTag` template parameter to split up observations into subsets
111 : * of elements. The `observers::Tags::ObservationKey<ArraySectionIdTag>` must be
112 : * available in the DataBox. It identifies the section and is used as a suffix
113 : * for the path in the output file.
114 : *
115 : * \par Option name
116 : * The `OptionName` template parameter is used to give the event a name in the
117 : * input file. If it is not specified, the name defaults to "ObserveNorms". If
118 : * you have multiple `ObserveNorms` events in the input file, you must specify a
119 : * unique name for each one. This can happen, for example, if you want to
120 : * observe norms the full domain and also over a section of the domain.
121 : */
122 : template <typename ObservableTensorTagsList,
123 : typename NonTensorComputeTagsList = tmpl::list<>,
124 : typename ArraySectionIdTag = void, typename OptionName = void>
125 1 : class ObserveNorms;
126 :
127 : template <typename... ObservableTensorTags, typename... NonTensorComputeTags,
128 : typename ArraySectionIdTag, typename OptionName>
129 0 : class ObserveNorms<tmpl::list<ObservableTensorTags...>,
130 : tmpl::list<NonTensorComputeTags...>, ArraySectionIdTag,
131 : OptionName> : public Event {
132 : private:
133 0 : struct ObserveTensor {
134 0 : static constexpr Options::String help = {
135 : "The tensor to reduce, and how to reduce it."};
136 :
137 0 : struct Name {
138 0 : using type = std::string;
139 0 : static constexpr Options::String help = {
140 : "The name of the tensor to observe."};
141 : };
142 0 : struct NormType {
143 0 : using type = std::string;
144 0 : static constexpr Options::String help = {
145 : "The type of norm to use. Must be one of Max, Min, L1Norm, "
146 : "L1IntegralNorm, L2Norm, L2IntegralNorm, or VolumeIntegral."};
147 : };
148 0 : struct Components {
149 0 : using type = std::string;
150 0 : static constexpr Options::String help = {
151 : "How to handle tensor components. Must be Individual or Sum."};
152 : };
153 :
154 0 : using options = tmpl::list<Name, NormType, Components>;
155 :
156 0 : ObserveTensor() = default;
157 :
158 0 : ObserveTensor(std::string in_tensor, std::string in_norm_type,
159 : std::string in_components,
160 : const Options::Context& context = {});
161 :
162 0 : std::string tensor{};
163 0 : std::string norm_type{};
164 0 : std::string components{};
165 : };
166 :
167 0 : using ReductionData = Parallel::ReductionData<
168 : // Observation value
169 : Parallel::ReductionDatum<double, funcl::AssertEqual<>>,
170 : // Number of grid points
171 : Parallel::ReductionDatum<size_t, funcl::Plus<>>,
172 : // Total volume
173 : Parallel::ReductionDatum<double, funcl::Plus<>>,
174 : // Max
175 : Parallel::ReductionDatum<std::vector<double>,
176 : funcl::ElementWise<funcl::Max<>>>,
177 : // Min
178 : Parallel::ReductionDatum<std::vector<double>,
179 : funcl::ElementWise<funcl::Min<>>>,
180 : // L1Norm
181 : Parallel::ReductionDatum<
182 : std::vector<double>, funcl::ElementWise<funcl::Plus<>>,
183 : funcl::ElementWise<funcl::Divides<>>, std::index_sequence<1>>,
184 : // L1IntegralNorm
185 : Parallel::ReductionDatum<
186 : std::vector<double>, funcl::ElementWise<funcl::Plus<>>,
187 : funcl::ElementWise<funcl::Divides<>>, std::index_sequence<2>>,
188 : // L2Norm
189 : Parallel::ReductionDatum<
190 : std::vector<double>, funcl::ElementWise<funcl::Plus<>>,
191 : funcl::ElementWise<funcl::Sqrt<funcl::Divides<>>>,
192 : std::index_sequence<1>>,
193 : // L2IntegralNorm
194 : Parallel::ReductionDatum<
195 : std::vector<double>, funcl::ElementWise<funcl::Plus<>>,
196 : funcl::ElementWise<funcl::Sqrt<funcl::Divides<>>>,
197 : std::index_sequence<2>>,
198 : // VolumeIntegral
199 : Parallel::ReductionDatum<std::vector<double>,
200 : funcl::ElementWise<funcl::Plus<>>>>;
201 :
202 : public:
203 0 : static std::string name() {
204 : if constexpr (std::is_same_v<OptionName, void>) {
205 : return "ObserveNorms";
206 : } else {
207 : return pretty_type::name<OptionName>();
208 : }
209 : }
210 :
211 : /// The name of the subfile inside the HDF5 file
212 1 : struct SubfileName {
213 0 : using type = std::string;
214 0 : static constexpr Options::String help = {
215 : "The name of the subfile inside the HDF5 file without an extension and "
216 : "without a preceding '/'."};
217 : };
218 : /// The tensor to observe and how to do the reduction
219 1 : struct TensorsToObserve {
220 0 : using type = std::vector<ObserveTensor>;
221 0 : static constexpr Options::String help = {
222 : "List specifying each tensor to observe and how it is reduced."};
223 : };
224 :
225 0 : explicit ObserveNorms(CkMigrateMessage* msg);
226 : using PUP::able::register_constructor;
227 0 : WRAPPED_PUPable_decl_template(ObserveNorms); // NOLINT
228 :
229 0 : using options = tmpl::list<SubfileName, TensorsToObserve>;
230 :
231 0 : static constexpr Options::String help =
232 : "Observe norms of tensors in the DataBox.\n"
233 : "\n"
234 : "You can choose the norm type for each observation. Note that 'L1Norm'\n"
235 : "(mean absolute value) and 'L2Norm' (root mean square) emphasize "
236 : "regions\n"
237 : "of the domain with many grid points, whereas 'L1IntegralNorm' and\n"
238 : "'L2IntegralNorm' emphasize regions of the domain with large volume.\n"
239 : "Choose wisely! When in doubt, try the 'L2Norm' first.\n"
240 : "\n"
241 : "Writes reduction quantities:\n"
242 : " * Observation value (e.g. Time or IterationId)\n"
243 : " * NumberOfPoints = total number of points in the domain\n"
244 : " * Volume = total volume of the domain in inertial coordinates\n"
245 : " * Max values\n"
246 : " * Min values\n"
247 : " * L1-norm values\n"
248 : " * L1 integral norm values\n"
249 : " * L2-norm values\n"
250 : " * L2 integral norm values\n"
251 : " * Volume integral values\n";
252 :
253 0 : ObserveNorms() = default;
254 :
255 0 : ObserveNorms(const std::string& subfile_name,
256 : const std::vector<ObserveTensor>& observe_tensors);
257 :
258 0 : using observed_reduction_data_tags =
259 : observers::make_reduction_data_tags<tmpl::list<ReductionData>>;
260 :
261 0 : using compute_tags_for_observation_box =
262 : tmpl::list<ObservableTensorTags..., NonTensorComputeTags...>;
263 :
264 0 : using return_tags = tmpl::list<>;
265 0 : using argument_tags = tmpl::list<::Tags::ObservationBox>;
266 :
267 : template <typename TensorToObserveTag, typename ComputeTagsList,
268 : typename DataBoxType, size_t Dim>
269 0 : void observe_norms_impl(
270 : gsl::not_null<
271 : std::unordered_map<std::string, std::pair<std::vector<double>,
272 : std::vector<std::string>>>*>
273 : norm_values_and_names,
274 : const ObservationBox<ComputeTagsList, DataBoxType>& box,
275 : const Mesh<Dim>& mesh, const DataVector& det_jacobian,
276 : size_t number_of_points) const;
277 :
278 : template <typename ComputeTagsList, typename DataBoxType,
279 : typename Metavariables, size_t VolumeDim,
280 : typename ParallelComponent>
281 0 : void operator()(const ObservationBox<ComputeTagsList, DataBoxType>& box,
282 : Parallel::GlobalCache<Metavariables>& cache,
283 : const ElementId<VolumeDim>& array_index,
284 : const ParallelComponent* const /*meta*/,
285 : const ObservationValue& observation_value) const;
286 :
287 0 : using observation_registration_tags = tmpl::list<::Tags::DataBox>;
288 :
289 : template <typename DbTagsList>
290 : std::optional<
291 : std::pair<observers::TypeOfObservation, observers::ObservationKey>>
292 0 : get_observation_type_and_key_for_registration(
293 : const db::DataBox<DbTagsList>& box) const {
294 : const std::optional<std::string> section_observation_key =
295 : observers::get_section_observation_key<ArraySectionIdTag>(box);
296 : if (not section_observation_key.has_value()) {
297 : return std::nullopt;
298 : }
299 : return {{observers::TypeOfObservation::Reduction,
300 : observers::ObservationKey(
301 : subfile_path_ + section_observation_key.value() + ".dat")}};
302 : }
303 :
304 0 : using is_ready_argument_tags = tmpl::list<>;
305 :
306 : template <typename Metavariables, typename ArrayIndex, typename Component>
307 0 : bool is_ready(Parallel::GlobalCache<Metavariables>& /*cache*/,
308 : const ArrayIndex& /*array_index*/,
309 : const Component* const /*meta*/) const {
310 : return true;
311 : }
312 :
313 1 : bool needs_evolved_variables() const override { return true; }
314 :
315 : // NOLINTNEXTLINE(google-runtime-references)
316 0 : void pup(PUP::er& p) override;
317 :
318 : private:
319 0 : std::string subfile_path_;
320 0 : std::vector<std::string> tensor_names_{};
321 0 : std::vector<std::string> tensor_norm_types_{};
322 0 : std::vector<std::string> tensor_components_{};
323 : };
324 : /// @}
325 :
326 : /// \cond
327 : template <typename... ObservableTensorTags, typename... NonTensorComputeTags,
328 : typename ArraySectionIdTag, typename OptionName>
329 : ObserveNorms<tmpl::list<ObservableTensorTags...>,
330 : tmpl::list<NonTensorComputeTags...>, ArraySectionIdTag,
331 : OptionName>::ObserveNorms(CkMigrateMessage* msg)
332 : : Event(msg) {}
333 :
334 : template <typename... ObservableTensorTags, typename... NonTensorComputeTags,
335 : typename ArraySectionIdTag, typename OptionName>
336 : ObserveNorms<tmpl::list<ObservableTensorTags...>,
337 : tmpl::list<NonTensorComputeTags...>, ArraySectionIdTag,
338 : OptionName>::ObserveNorms(const std::string& subfile_name,
339 : const std::vector<ObserveTensor>&
340 : observe_tensors)
341 : : subfile_path_("/" + subfile_name) {
342 : tensor_names_.reserve(observe_tensors.size());
343 : tensor_norm_types_.reserve(observe_tensors.size());
344 : tensor_components_.reserve(observe_tensors.size());
345 : for (const auto& observe_tensor : observe_tensors) {
346 : tensor_names_.push_back(observe_tensor.tensor);
347 : tensor_norm_types_.push_back(observe_tensor.norm_type);
348 : tensor_components_.push_back(observe_tensor.components);
349 : }
350 : }
351 :
352 : template <typename... ObservableTensorTags, typename... NonTensorComputeTags,
353 : typename ArraySectionIdTag, typename OptionName>
354 : ObserveNorms<
355 : tmpl::list<ObservableTensorTags...>, tmpl::list<NonTensorComputeTags...>,
356 : ArraySectionIdTag,
357 : OptionName>::ObserveTensor::ObserveTensor(std::string in_tensor,
358 : std::string in_norm_type,
359 : std::string in_components,
360 : const Options::Context& context)
361 : : tensor(std::move(in_tensor)),
362 : norm_type(std::move(in_norm_type)),
363 : components(std::move(in_components)) {
364 : if (((tensor != db::tag_name<ObservableTensorTags>()) and ...)) {
365 : PARSE_ERROR(
366 : context, "Tensor '"
367 : << tensor << "' is not known. Known tensors are: "
368 : << ((db::tag_name<ObservableTensorTags>() + ",") + ...));
369 : }
370 : if (norm_type != "Max" and norm_type != "Min" and norm_type != "L1Norm" and
371 : norm_type != "L1IntegralNorm" and norm_type != "L2Norm" and
372 : norm_type != "L2IntegralNorm" and norm_type != "VolumeIntegral") {
373 : PARSE_ERROR(
374 : context,
375 : "NormType must be one of Max, Min, L1Norm, L1IntegralNorm, L2Norm, "
376 : "L2IntegralNorm, or VolumeIntegral not "
377 : << norm_type);
378 : }
379 : if (components != "Individual" and components != "Sum") {
380 : PARSE_ERROR(context,
381 : "Components must be Individual or Sum, not " << components);
382 : }
383 : }
384 :
385 : // implementation of ObserveNorms::operator() factored out to save on compile
386 : // time and compile memory
387 : namespace ObserveNorms_impl {
388 : void check_norm_is_observable(const std::string& tensor_name,
389 : bool tag_has_value);
390 :
391 : template <size_t Dim>
392 : void fill_norm_values_and_names(
393 : gsl::not_null<std::unordered_map<
394 : std::string, std::pair<std::vector<double>, std::vector<std::string>>>*>
395 : norm_values_and_names,
396 : const std::pair<std::vector<std::string>, std::vector<DataVector>>&
397 : names_and_components,
398 : const Mesh<Dim>& mesh, const DataVector& det_jacobian,
399 : const std::string& tensor_name, const std::string& tensor_norm_type,
400 : const std::string& tensor_component, size_t number_of_points);
401 :
402 : // Expand complex data into real and imaginary parts, or just forward real data
403 : std::pair<std::vector<std::string>, std::vector<DataVector>>
404 : split_complex_vector_of_data(
405 : std::pair<std::vector<std::string>, std::vector<DataVector>>&&
406 : names_and_components);
407 : std::pair<std::vector<std::string>, std::vector<DataVector>>
408 : split_complex_vector_of_data(
409 : const std::pair<std::vector<std::string>, std::vector<ComplexDataVector>>&
410 : names_and_components);
411 : } // namespace ObserveNorms_impl
412 :
413 : template <typename... ObservableTensorTags, typename... NonTensorComputeTags,
414 : typename ArraySectionIdTag, typename OptionName>
415 : template <typename TensorToObserveTag, typename ComputeTagsList,
416 : typename DataBoxType, size_t Dim>
417 : void ObserveNorms<tmpl::list<ObservableTensorTags...>,
418 : tmpl::list<NonTensorComputeTags...>, ArraySectionIdTag,
419 : OptionName>::
420 : observe_norms_impl(
421 : const gsl::not_null<std::unordered_map<
422 : std::string,
423 : std::pair<std::vector<double>, std::vector<std::string>>>*>
424 : norm_values_and_names,
425 : const ObservationBox<ComputeTagsList, DataBoxType>& box,
426 : const Mesh<Dim>& mesh, const DataVector& det_jacobian,
427 : const size_t number_of_points) const {
428 : const std::string tensor_name = db::tag_name<TensorToObserveTag>();
429 : for (size_t i = 0; i < tensor_names_.size(); ++i) {
430 : if (tensor_name == tensor_names_[i]) {
431 : ObserveNorms_impl::check_norm_is_observable(
432 : tensor_name, has_value(get<TensorToObserveTag>(box)));
433 : ObserveNorms_impl::fill_norm_values_and_names(
434 : norm_values_and_names,
435 : ObserveNorms_impl::split_complex_vector_of_data(
436 : value(get<TensorToObserveTag>(box)).get_vector_of_data()),
437 : mesh, det_jacobian, tensor_name, tensor_norm_types_[i],
438 : tensor_components_[i], number_of_points);
439 : }
440 : }
441 : }
442 :
443 : template <typename... ObservableTensorTags, typename... NonTensorComputeTags,
444 : typename ArraySectionIdTag, typename OptionName>
445 : template <typename ComputeTagsList, typename DataBoxType,
446 : typename Metavariables, size_t VolumeDim, typename ParallelComponent>
447 : void ObserveNorms<tmpl::list<ObservableTensorTags...>,
448 : tmpl::list<NonTensorComputeTags...>, ArraySectionIdTag,
449 : OptionName>::
450 : operator()(const ObservationBox<ComputeTagsList, DataBoxType>& box,
451 : Parallel::GlobalCache<Metavariables>& cache,
452 : const ElementId<VolumeDim>& array_index,
453 : const ParallelComponent* const /*meta*/,
454 : const ObservationValue& observation_value) const {
455 : // Skip observation on elements that are not part of a section
456 : const std::optional<std::string> section_observation_key =
457 : observers::get_section_observation_key<ArraySectionIdTag>(box);
458 : if (not section_observation_key.has_value()) {
459 : return;
460 : }
461 :
462 : const auto& mesh = get<::Events::Tags::ObserverMesh<VolumeDim>>(box);
463 : const auto det_jacobian = [&box, &mesh]() -> DataVector {
464 : if constexpr (VolumeDim == 3 and
465 : db::tag_is_retrievable_v<::Events::Tags::ObserverCoordinates<
466 : VolumeDim, Frame::Inertial>,
467 : std::decay_t<decltype(box)>>) {
468 : if (mesh.basis(2) == Spectral::Basis::Cartoon) {
469 : if (mesh.quadrature(2) == Spectral::Quadrature::SphericalSymmetry) {
470 : // Spherical Symmetry, needs x^2 cartesian to spherical jacobian
471 : return square(get<0>(get<::Events::Tags::ObserverCoordinates<
472 : VolumeDim, Frame::Inertial>>(box))) /
473 : get(get<::Events::Tags::ObserverDetInvJacobian<
474 : Frame::ElementLogical, Frame::Inertial>>(box));
475 : } else {
476 : // Axial Symmetry, needs x cartesian to cylindrical jacobian
477 : ASSERT(mesh.quadrature(2) == Spectral::Quadrature::AxialSymmetry,
478 : "Unexpected quadrature " << mesh.quadrature(2)
479 : << " (expected AxialSymmetry)");
480 : return get<0>(get<::Events::Tags::ObserverCoordinates<
481 : VolumeDim, Frame::Inertial>>(box)) /
482 : get(get<::Events::Tags::ObserverDetInvJacobian<
483 : Frame::ElementLogical, Frame::Inertial>>(box));
484 : }
485 : } else {
486 : return 1. / get(get<::Events::Tags::ObserverDetInvJacobian<
487 : Frame::ElementLogical, Frame::Inertial>>(box));
488 : }
489 : } else {
490 : (void)mesh;
491 : return 1. / get(get<::Events::Tags::ObserverDetInvJacobian<
492 : Frame::ElementLogical, Frame::Inertial>>(box));
493 : }
494 : }();
495 : const size_t number_of_points = mesh.number_of_grid_points();
496 : const double local_volume = definite_integral(det_jacobian, mesh);
497 :
498 : std::unordered_map<std::string,
499 : std::pair<std::vector<double>, std::vector<std::string>>>
500 : norm_values_and_names{};
501 : // Loop over ObservableTensorTags and see if it was requested to be observed.
502 : // This approach allows us to delay evaluating any compute tags until they're
503 : // actually needed for observing.
504 : (observe_norms_impl<ObservableTensorTags>(
505 : make_not_null(&norm_values_and_names), box, mesh, det_jacobian,
506 : number_of_points),
507 : ...);
508 :
509 : // Concatenate the legend info together.
510 : std::vector<std::string> legend{observation_value.name, "NumberOfPoints",
511 : "Volume"};
512 : legend.insert(legend.end(), norm_values_and_names["Max"].second.begin(),
513 : norm_values_and_names["Max"].second.end());
514 : legend.insert(legend.end(), norm_values_and_names["Min"].second.begin(),
515 : norm_values_and_names["Min"].second.end());
516 : legend.insert(legend.end(), norm_values_and_names["L1Norm"].second.begin(),
517 : norm_values_and_names["L1Norm"].second.end());
518 : legend.insert(legend.end(),
519 : norm_values_and_names["L1IntegralNorm"].second.begin(),
520 : norm_values_and_names["L1IntegralNorm"].second.end());
521 : legend.insert(legend.end(), norm_values_and_names["L2Norm"].second.begin(),
522 : norm_values_and_names["L2Norm"].second.end());
523 : legend.insert(legend.end(),
524 : norm_values_and_names["L2IntegralNorm"].second.begin(),
525 : norm_values_and_names["L2IntegralNorm"].second.end());
526 : legend.insert(legend.end(),
527 : norm_values_and_names["VolumeIntegral"].second.begin(),
528 : norm_values_and_names["VolumeIntegral"].second.end());
529 :
530 : const std::string subfile_path_with_suffix =
531 : subfile_path_ + section_observation_key.value();
532 : // Send data to reduction observer
533 : auto& local_observer = *Parallel::local_branch(
534 : Parallel::get_parallel_component<
535 : tmpl::conditional_t<Parallel::is_nodegroup_v<ParallelComponent>,
536 : observers::ObserverWriter<Metavariables>,
537 : observers::Observer<Metavariables>>>(cache));
538 : observers::ObservationId observation_id{observation_value.value,
539 : subfile_path_with_suffix + ".dat"};
540 : Parallel::ArrayComponentId array_component_id{
541 : std::add_pointer_t<ParallelComponent>{nullptr},
542 : Parallel::ArrayIndex<ElementId<VolumeDim>>(array_index)};
543 : ReductionData reduction_data{
544 : observation_value.value,
545 : number_of_points,
546 : local_volume,
547 : std::move(norm_values_and_names["Max"].first),
548 : std::move(norm_values_and_names["Min"].first),
549 : std::move(norm_values_and_names["L1Norm"].first),
550 : std::move(norm_values_and_names["L1IntegralNorm"].first),
551 : std::move(norm_values_and_names["L2Norm"].first),
552 : std::move(norm_values_and_names["L2IntegralNorm"].first),
553 : std::move(norm_values_and_names["VolumeIntegral"].first)};
554 :
555 : if constexpr (Parallel::is_nodegroup_v<ParallelComponent>) {
556 : Parallel::threaded_action<
557 : observers::ThreadedActions::CollectReductionDataOnNode>(
558 : local_observer, std::move(observation_id),
559 : std::move(array_component_id), subfile_path_with_suffix,
560 : std::move(legend), std::move(reduction_data));
561 : } else {
562 : Parallel::simple_action<observers::Actions::ContributeReductionData>(
563 : local_observer, std::move(observation_id),
564 : std::move(array_component_id), subfile_path_with_suffix,
565 : std::move(legend), std::move(reduction_data));
566 : }
567 : }
568 :
569 : template <typename... ObservableTensorTags, typename... NonTensorComputeTags,
570 : typename ArraySectionIdTag, typename OptionName>
571 : void ObserveNorms<tmpl::list<ObservableTensorTags...>,
572 : tmpl::list<NonTensorComputeTags...>, ArraySectionIdTag,
573 : OptionName>::pup(PUP::er& p) {
574 : Event::pup(p);
575 : p | subfile_path_;
576 : p | tensor_names_;
577 : p | tensor_norm_types_;
578 : p | tensor_components_;
579 : }
580 :
581 : // NOLINTBEGIN(cppcoreguidelines-avoid-non-const-global-variables)
582 : template <typename... ObservableTensorTags, typename... NonTensorComputeTags,
583 : typename ArraySectionIdTag, typename OptionName>
584 : PUP::able::PUP_ID ObserveNorms<tmpl::list<ObservableTensorTags...>,
585 : tmpl::list<NonTensorComputeTags...>,
586 : ArraySectionIdTag, OptionName>::my_PUP_ID = 0;
587 : // NOLINTEND(cppcoreguidelines-avoid-non-const-global-variables)
588 : /// \endcond
589 : } // namespace Events
|
