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