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 <unordered_map>
9 :
10 : #include "DataStructures/DataBox/DataBox.hpp"
11 : #include "DataStructures/DataBox/Prefixes.hpp"
12 : #include "DataStructures/DataVector.hpp"
13 : #include "DataStructures/ModalVector.hpp"
14 : #include "DataStructures/Variables.hpp"
15 : #include "Domain/Structure/ElementId.hpp"
16 : #include "Domain/Tags.hpp"
17 : #include "Evolution/Systems/CurvedScalarWave/Worldtube/Inboxes.hpp"
18 : #include "Evolution/Systems/CurvedScalarWave/Worldtube/SingletonActions/UpdateAcceleration.hpp"
19 : #include "Evolution/Systems/CurvedScalarWave/Worldtube/Tags.hpp"
20 : #include "NumericalAlgorithms/SphericalHarmonics/YlmToStf.hpp"
21 : #include "NumericalAlgorithms/Strahlkorper/Tags.hpp"
22 : #include "Parallel/AlgorithmExecution.hpp"
23 : #include "Parallel/GlobalCache.hpp"
24 : #include "ParallelAlgorithms/Actions/MutateApply.hpp"
25 : #include "ParallelAlgorithms/Initialization/MutateAssign.hpp"
26 : #include "Time/TimeStepId.hpp"
27 : #include "Utilities/ErrorHandling/Assert.hpp"
28 : #include "Utilities/Gsl.hpp"
29 : #include "Utilities/TMPL.hpp"
30 :
31 : /// \cond
32 : namespace Tags {
33 : struct TimeStepId;
34 : } // namespace Tags
35 : /// \endcond
36 :
37 : namespace CurvedScalarWave::Worldtube::Actions {
38 : /*!
39 : * \brief Adds up the spherical harmonic projections from the different elements
40 : * abutting the worldtube.
41 : *
42 : * \details This action currently assumes that there is no h-refinement
43 : * ocurring in the elements abutting the worldtubes. This could be accounted for
44 : * by checking that data from at least one element has been sent from each
45 : * abutting block and then using its `ElementId` to figure out the current
46 : * refinement level and therefore how many elements are expected to send data
47 : * for each block.
48 : *
49 : * DataBox:
50 : * - Uses:
51 : * - `Worldtube::Tags::ExpansionOrder`
52 : * - `Worldtube::Tags::ExcisionSphere`
53 : * - `Worldtube::Tags::ElementFacesGridCoordinates`
54 : * - `Tags::TimeStepId`
55 : * - Mutates:
56 : * - `Stf::Tags::StfTensor<Tags::PsiWorldtube, 0, Dim, Frame::Inertial>`
57 : * - `Stf::Tags::StfTensor<::Tags::dt<Tags::PsiWorldtube>, 0, Dim,
58 : * Frame::Inertial>`
59 : * - `Stf::Tags::StfTensor<Tags::PsiWorldtube, 1, Dim, Frame::Inertial>`
60 : * - `Stf::Tags::StfTensor<::Tags::dt<Tags::PsiWorldtube>, 1, Dim,
61 : * Frame::Inertial>`
62 : */
63 1 : struct ReceiveElementData {
64 0 : static constexpr size_t Dim = 3;
65 0 : using tags_list = tmpl::list<CurvedScalarWave::Tags::Psi,
66 : ::Tags::dt<CurvedScalarWave::Tags::Psi>>;
67 0 : using inbox_tags = tmpl::list<
68 : ::CurvedScalarWave::Worldtube::Tags::SphericalHarmonicsInbox<Dim>>;
69 0 : using simple_tags = tmpl::list<
70 : Stf::Tags::StfTensor<Tags::PsiWorldtube, 0, Dim, Frame::Inertial>,
71 : Stf::Tags::StfTensor<::Tags::dt<Tags::PsiWorldtube>, 0, Dim,
72 : Frame::Inertial>,
73 : Stf::Tags::StfTensor<Tags::PsiWorldtube, 1, Dim, Frame::Inertial>,
74 : Stf::Tags::StfTensor<::Tags::dt<Tags::PsiWorldtube>, 1, Dim,
75 : Frame::Inertial>>;
76 :
77 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
78 : typename ArrayIndex, typename ActionList,
79 : typename ParallelComponent>
80 0 : static Parallel::iterable_action_return_t apply(
81 : db::DataBox<DbTagsList>& box, tuples::TaggedTuple<InboxTags...>& inboxes,
82 : const Parallel::GlobalCache<Metavariables>& /*cache*/,
83 : const ArrayIndex& /*array_index*/, ActionList /*meta*/,
84 : const ParallelComponent* const /*meta*/) {
85 : const size_t expected_number_of_senders =
86 : db::get<Tags::ElementFacesGridCoordinates<Dim>>(box).size();
87 : const auto& time_step_id = db::get<::Tags::TimeStepId>(box);
88 : auto& inbox = tuples::get<Tags::SphericalHarmonicsInbox<Dim>>(inboxes);
89 : if (inbox.count(time_step_id) == 0 or
90 : inbox.at(time_step_id).size() < expected_number_of_senders) {
91 : return {Parallel::AlgorithmExecution::Retry, std::nullopt};
92 : }
93 : ASSERT(inbox.at(time_step_id).size() == expected_number_of_senders,
94 : "Expected data from "
95 : << expected_number_of_senders << " senders, but received "
96 : << inbox.at(time_step_id).size() << " for TimeStepId "
97 : << time_step_id);
98 : const size_t order = db::get<Tags::ExpansionOrder>(box);
99 : const size_t num_modes = (order + 1) * (order + 1);
100 :
101 : Variables<tags_list> external_ylm_coefs{num_modes, 0.};
102 : for (const auto& [_, element_ylm_coefs] : inbox.at(time_step_id)) {
103 : external_ylm_coefs += element_ylm_coefs;
104 : }
105 : const double wt_radius = db::get<Tags::WorldtubeRadius>(box);
106 : external_ylm_coefs /= wt_radius * wt_radius;
107 :
108 : DataVector& psi_ylm_coefs =
109 : get(get<CurvedScalarWave::Tags::Psi>(external_ylm_coefs));
110 : DataVector& dt_psi_ylm_coefs =
111 : get(get<::Tags::dt<CurvedScalarWave::Tags::Psi>>(external_ylm_coefs));
112 :
113 : const ModalVector psi_ylm_l0(&psi_ylm_coefs[0], 1);
114 : const ModalVector dt_psi_ylm_l0(&dt_psi_ylm_coefs[0], 1);
115 : tnsr::i<double, Dim, Frame::Inertial> psi_stf_l1{};
116 : tnsr::i<double, Dim, Frame::Inertial> dt_psi_stf_l1{};
117 : if (order > 0) {
118 : ModalVector psi_ylm_l1(&psi_ylm_coefs[1], 3);
119 : ModalVector dt_psi_ylm_l1(&dt_psi_ylm_coefs[1], 3);
120 : psi_ylm_l1 /= wt_radius;
121 : dt_psi_ylm_l1 /= wt_radius;
122 : psi_stf_l1 = ylm::ylm_to_stf_1<Frame::Inertial>(psi_ylm_l1);
123 : dt_psi_stf_l1 = ylm::ylm_to_stf_1<Frame::Inertial>(dt_psi_ylm_l1);
124 : }
125 :
126 : ::Initialization::mutate_assign<simple_tags>(
127 : make_not_null(&box), ylm::ylm_to_stf_0(psi_ylm_l0),
128 : ylm::ylm_to_stf_0(dt_psi_ylm_l0), std::move(psi_stf_l1),
129 : std::move(dt_psi_stf_l1));
130 : inbox.erase(time_step_id);
131 : if (db::get<Tags::CurrentIteration>(box) + 1 <
132 : db::get<Tags::MaxIterations>(box)) {
133 : db::mutate<Tags::CurrentIteration>(
134 : [](const gsl::not_null<size_t*> current_iteration) {
135 : *current_iteration += 1;
136 : },
137 : make_not_null(&box));
138 : // still iterating, go to `IterateAccelerationTerms`
139 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
140 :
141 : } else {
142 : db::mutate<Tags::CurrentIteration>(
143 : [](const gsl::not_null<size_t*> current_iteration) {
144 : *current_iteration = 0;
145 : },
146 : make_not_null(&box));
147 : // done iterating
148 : return {
149 : Parallel::AlgorithmExecution::Continue,
150 : tmpl::index_of<ActionList,
151 : ::Actions::MutateApply<UpdateAcceleration>>::value};
152 : }
153 : }
154 : };
155 : } // namespace CurvedScalarWave::Worldtube::Actions
|
