Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <array>
7 : #include <cstddef>
8 : #include <map>
9 : #include <optional>
10 : #include <tuple>
11 : #include <type_traits>
12 : #include <utility>
13 :
14 : #include "DataStructures/DataBox/DataBox.hpp"
15 : #include "DataStructures/DataBox/PrefixHelpers.hpp"
16 : #include "DataStructures/FixedHashMap.hpp"
17 : #include "DataStructures/TaggedTuple.hpp"
18 : #include "Domain/Structure/ChildSize.hpp"
19 : #include "Domain/Structure/ElementId.hpp"
20 : #include "Domain/Tags.hpp"
21 : #include "IO/Logging/Tags.hpp"
22 : #include "IO/Logging/Verbosity.hpp"
23 : #include "IO/Observer/Tags.hpp"
24 : #include "NumericalAlgorithms/Convergence/Tags.hpp"
25 : #include "NumericalAlgorithms/Spectral/Projection.hpp"
26 : #include "NumericalAlgorithms/Spectral/SegmentSize.hpp"
27 : #include "Parallel/AlgorithmExecution.hpp"
28 : #include "Parallel/GlobalCache.hpp"
29 : #include "Parallel/InboxInserters.hpp"
30 : #include "Parallel/Invoke.hpp"
31 : #include "Parallel/Printf/Printf.hpp"
32 : #include "ParallelAlgorithms/Amr/Tags.hpp"
33 : #include "ParallelAlgorithms/LinearSolver/Multigrid/Tags.hpp"
34 : #include "ParallelAlgorithms/LinearSolver/Tags.hpp"
35 : #include "Utilities/ConstantExpressions.hpp"
36 : #include "Utilities/ErrorHandling/Assert.hpp"
37 : #include "Utilities/GetOutput.hpp"
38 : #include "Utilities/Gsl.hpp"
39 : #include "Utilities/MakeArray.hpp"
40 : #include "Utilities/PrettyType.hpp"
41 : #include "Utilities/TMPL.hpp"
42 :
43 : /// \cond
44 : template <size_t Dim>
45 : struct ElementId;
46 : /// \endcond
47 :
48 1 : namespace LinearSolver::multigrid {
49 :
50 : template <size_t Dim, typename ReceiveTags>
51 0 : struct DataFromChildrenInboxTag
52 : : public Parallel::InboxInserters::Map<
53 : DataFromChildrenInboxTag<Dim, ReceiveTags>> {
54 0 : using temporal_id = size_t;
55 0 : using type =
56 : std::map<temporal_id,
57 : FixedHashMap<two_to_the(Dim), ElementId<Dim>,
58 : tuples::tagged_tuple_from_typelist<ReceiveTags>,
59 : boost::hash<ElementId<Dim>>>>;
60 : };
61 :
62 : /// Actions related to the Multigrid linear solver
63 1 : namespace Actions {
64 : /*!
65 : * \brief Communicate and project the `FieldsTags` to the next-coarser grid
66 : * in the multigrid hierarchy
67 : *
68 : * \tparam FieldsTags These tags will be communicated and projected. They can
69 : * hold any type that works with `::apply_matrices` and supports addition, e.g.
70 : * `Variables`.
71 : * \tparam OptionsGroup The option group identifying the multigrid solver
72 : * \tparam FieldsAreMassiveTag A boolean tag in the DataBox that indicates
73 : * whether or not the `FieldsTags` have already been multiplied by the mass
74 : * matrix. This setting influences the way the fields are projected. In
75 : * particular, the mass matrix already includes a Jacobian factor, so the
76 : * difference in size between the parent and the child element is already
77 : * accounted for.
78 : * \tparam ReceiveTags The projected fields will be stored in these tags
79 : * (default: `FieldsTags`).
80 : */
81 : template <typename FieldsTags, typename OptionsGroup,
82 : typename FieldsAreMassiveTag, typename ReceiveTags = FieldsTags>
83 1 : struct SendFieldsToCoarserGrid;
84 :
85 : /// \cond
86 : template <typename... FieldsTags, typename OptionsGroup,
87 : typename FieldsAreMassiveTag, typename... ReceiveTags>
88 : struct SendFieldsToCoarserGrid<tmpl::list<FieldsTags...>, OptionsGroup,
89 : FieldsAreMassiveTag,
90 : tmpl::list<ReceiveTags...>> {
91 : using const_global_cache_tags =
92 : tmpl::list<logging::Tags::Verbosity<OptionsGroup>>;
93 :
94 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
95 : size_t Dim, typename ActionList, typename ParallelComponent>
96 : static Parallel::iterable_action_return_t apply(
97 : db::DataBox<DbTagsList>& box,
98 : const tuples::TaggedTuple<InboxTags...>& /*inboxes*/,
99 : Parallel::GlobalCache<Metavariables>& cache,
100 : const ElementId<Dim>& element_id, const ActionList /*meta*/,
101 : const ParallelComponent* const /*meta*/) {
102 : // Skip restriction on coarsest level
103 : const auto& parent_id = db::get<amr::Tags::ParentId<Dim>>(box);
104 : if (not parent_id.has_value()) {
105 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
106 : }
107 :
108 : const size_t iteration_id =
109 : db::get<Convergence::Tags::IterationId<OptionsGroup>>(box);
110 : if (UNLIKELY(db::get<logging::Tags::Verbosity<OptionsGroup>>(box) >=
111 : ::Verbosity::Debug)) {
112 : Parallel::printf("%s %s(%zu): Send fields to coarser grid %s\n",
113 : element_id, pretty_type::name<OptionsGroup>(),
114 : iteration_id, parent_id);
115 : }
116 :
117 : // Restrict the fields to the coarser (parent) grid.
118 : // We restrict before sending the data so the restriction operation is
119 : // parellelized. The parent only needs to sum up all child contributions.
120 : const auto& mesh = db::get<domain::Tags::Mesh<Dim>>(box);
121 : const auto& parent_mesh = db::get<amr::Tags::ParentMesh<Dim>>(box);
122 : ASSERT(
123 : parent_mesh.has_value(),
124 : "Should have a parent mesh, because a parent ID is set. This element: "
125 : << element_id << ", parent element: " << *parent_id);
126 : const auto child_size =
127 : domain::child_size(element_id.segment_ids(), parent_id->segment_ids());
128 : bool massive = false;
129 : if constexpr (not std::is_same_v<FieldsAreMassiveTag, void>) {
130 : massive = db::get<FieldsAreMassiveTag>(box);
131 : }
132 : tuples::TaggedTuple<ReceiveTags...> restricted_fields{};
133 : if (Spectral::needs_projection(mesh, *parent_mesh, child_size)) {
134 : const auto restrict_fields =
135 : [&restricted_fields, &mesh, &parent_mesh, &child_size, massive](
136 : const auto receive_tag_v, const auto& fields) {
137 : using receive_tag = std::decay_t<decltype(receive_tag_v)>;
138 : get<receive_tag>(restricted_fields) =
139 : typename receive_tag::type(Spectral::project(
140 : fields, mesh, *parent_mesh, child_size,
141 : make_array<Dim>(Spectral::SegmentSize::Full), massive));
142 : return '0';
143 : };
144 : expand_pack(restrict_fields(ReceiveTags{}, db::get<FieldsTags>(box))...);
145 : } else {
146 : expand_pack(
147 : (get<ReceiveTags>(restricted_fields) =
148 : typename ReceiveTags::type(db::get<FieldsTags>(box)))...);
149 : }
150 :
151 : // Send restricted fields to the parent
152 : auto& receiver_proxy =
153 : Parallel::get_parallel_component<ParallelComponent>(cache);
154 : Parallel::receive_data<
155 : DataFromChildrenInboxTag<Dim, tmpl::list<ReceiveTags...>>>(
156 : receiver_proxy[*parent_id], iteration_id,
157 : std::make_pair(element_id, std::move(restricted_fields)));
158 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
159 : }
160 : };
161 : /// \endcond
162 :
163 : /// Receive the `FieldsTags` communicated from the finer grid in the multigrid
164 : /// hierarchy.
165 : ///
166 : /// \see LinearSolver::multigrid::Actions::SendFieldsToCoarserGrid
167 : template <size_t Dim, typename FieldsTags, typename OptionsGroup,
168 : typename ReceiveTags = FieldsTags>
169 1 : struct ReceiveFieldsFromFinerGrid;
170 :
171 : /// \cond
172 : template <size_t Dim, typename FieldsTags, typename OptionsGroup,
173 : typename... ReceiveTags>
174 : struct ReceiveFieldsFromFinerGrid<Dim, FieldsTags, OptionsGroup,
175 : tmpl::list<ReceiveTags...>> {
176 : using inbox_tags =
177 : tmpl::list<DataFromChildrenInboxTag<Dim, tmpl::list<ReceiveTags...>>>;
178 : using const_global_cache_tags =
179 : tmpl::list<logging::Tags::Verbosity<OptionsGroup>>;
180 :
181 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
182 : typename ActionList, typename ParallelComponent>
183 : static Parallel::iterable_action_return_t apply(
184 : db::DataBox<DbTagsList>& box, tuples::TaggedTuple<InboxTags...>& inboxes,
185 : const Parallel::GlobalCache<Metavariables>& /*cache*/,
186 : const ElementId<Dim>& element_id, const ActionList /*meta*/,
187 : const ParallelComponent* const /*meta*/) {
188 : // Skip on finest grid
189 : const auto& child_ids = db::get<amr::Tags::ChildIds<Dim>>(box);
190 : const bool is_finest_grid = child_ids.empty();
191 : if (is_finest_grid) {
192 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
193 : }
194 :
195 : // Wait for data from finer grid
196 : const size_t iteration_id =
197 : db::get<Convergence::Tags::IterationId<OptionsGroup>>(box);
198 : auto& inbox =
199 : tuples::get<DataFromChildrenInboxTag<Dim, tmpl::list<ReceiveTags...>>>(
200 : inboxes);
201 : const auto received_this_iteration = inbox.find(iteration_id);
202 : if (received_this_iteration == inbox.end()) {
203 : if (UNLIKELY(db::get<logging::Tags::Verbosity<OptionsGroup>>(box) >=
204 : ::Verbosity::Debug)) {
205 : Parallel::printf(
206 : "%s %s(%zu): Waiting for fine-grid data (still empty)\n",
207 : element_id, pretty_type::name<OptionsGroup>(), iteration_id);
208 : }
209 : return {Parallel::AlgorithmExecution::Retry, std::nullopt};
210 : }
211 : const auto& received_children_data = received_this_iteration->second;
212 : for (const auto& child_id : child_ids) {
213 : if (received_children_data.find(child_id) ==
214 : received_children_data.end()) {
215 : if (UNLIKELY(db::get<logging::Tags::Verbosity<OptionsGroup>>(box) >=
216 : ::Verbosity::Debug)) {
217 : Parallel::printf(
218 : "%s %s(%zu): Waiting for fine-grid data from child %s\n",
219 : element_id, pretty_type::name<OptionsGroup>(), iteration_id,
220 : child_id);
221 : }
222 : return {Parallel::AlgorithmExecution::Retry, std::nullopt};
223 : }
224 : }
225 : auto children_data = std::move(inbox.extract(iteration_id).mapped());
226 :
227 : if (UNLIKELY(db::get<logging::Tags::Verbosity<OptionsGroup>>(box) >=
228 : ::Verbosity::Debug)) {
229 : Parallel::printf("%s %s(%zu): Receive fields from finer grid\n",
230 : element_id, pretty_type::name<OptionsGroup>(),
231 : iteration_id);
232 : }
233 :
234 : // Assemble restricted data from children
235 : const auto assemble_children_data =
236 : [&children_data](const auto source, const auto receive_tag_v) {
237 : using receive_tag = std::decay_t<decltype(receive_tag_v)>;
238 : // Move the first child data directly into the buffer, then add the
239 : // data from the remaining children.
240 : auto child_id_and_data = children_data.begin();
241 : *source = std::move(get<receive_tag>(child_id_and_data->second));
242 : ++child_id_and_data;
243 : while (child_id_and_data != children_data.end()) {
244 : *source += get<receive_tag>(child_id_and_data->second);
245 : ++child_id_and_data;
246 : }
247 : return '0';
248 : };
249 : expand_pack(db::mutate<ReceiveTags>(assemble_children_data,
250 : make_not_null(&box), ReceiveTags{})...);
251 :
252 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
253 : }
254 : };
255 : /// \endcond
256 :
257 : } // namespace Actions
258 : } // namespace LinearSolver::multigrid
|