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 <unordered_map>
10 : #include <unordered_set>
11 :
12 : #include "DataStructures/DataBox/DataBox.hpp"
13 : #include "DataStructures/DataBox/PrefixHelpers.hpp"
14 : #include "DataStructures/FixedHashMap.hpp"
15 : #include "DataStructures/Matrix.hpp"
16 : #include "DataStructures/TaggedTuple.hpp"
17 : #include "Domain/Creators/Tags/InitialRefinementLevels.hpp"
18 : #include "Domain/Structure/ChildSize.hpp"
19 : #include "Domain/Structure/ElementId.hpp"
20 : #include "IO/Logging/Tags.hpp"
21 : #include "IO/Observer/Tags.hpp"
22 : #include "NumericalAlgorithms/Convergence/Tags.hpp"
23 : #include "NumericalAlgorithms/Spectral/Projection.hpp"
24 : #include "Parallel/AlgorithmExecution.hpp"
25 : #include "Parallel/GlobalCache.hpp"
26 : #include "Parallel/InboxInserters.hpp"
27 : #include "Parallel/Invoke.hpp"
28 : #include "Parallel/Printf/Printf.hpp"
29 : #include "ParallelAlgorithms/Actions/Goto.hpp"
30 : #include "ParallelAlgorithms/Amr/Protocols/Projector.hpp"
31 : #include "ParallelAlgorithms/Amr/Tags.hpp"
32 : #include "ParallelAlgorithms/LinearSolver/Multigrid/Actions/RestrictFields.hpp"
33 : #include "ParallelAlgorithms/LinearSolver/Multigrid/Hierarchy.hpp"
34 : #include "ParallelAlgorithms/LinearSolver/Multigrid/Tags.hpp"
35 : #include "ParallelAlgorithms/LinearSolver/Tags.hpp"
36 : #include "Utilities/ConstantExpressions.hpp"
37 : #include "Utilities/ErrorHandling/Assert.hpp"
38 : #include "Utilities/GetOutput.hpp"
39 : #include "Utilities/PrettyType.hpp"
40 : #include "Utilities/ProtocolHelpers.hpp"
41 : #include "Utilities/TypeTraits/IsA.hpp"
42 :
43 : namespace LinearSolver::multigrid::detail {
44 :
45 : /// \cond
46 : template <typename FieldsTag, typename OptionsGroup, typename SourceTag>
47 : struct SendCorrectionToFinerGrid;
48 : template <typename FieldsTag, typename OptionsGroup, typename SourceTag>
49 : struct SkipBottomSolver;
50 : template <typename FieldsTag, typename OptionsGroup, typename SourceTag>
51 : struct SkipPostSmoothingAtBottom;
52 : /// \endcond
53 :
54 : struct PostSmoothingBeginLabel {};
55 :
56 : template <size_t Dim, typename FieldsTag, typename OptionsGroup,
57 : typename SourceTag>
58 : struct InitializeElement : tt::ConformsTo<amr::protocols::Projector> {
59 : private:
60 : using VolumeDataVars =
61 : typename Tags::VolumeDataForOutput<OptionsGroup, FieldsTag>::type;
62 :
63 : public: // Iterable action
64 : using simple_tags_from_options =
65 : tmpl::list<Tags::ChildrenRefinementLevels<Dim>,
66 : Tags::ParentRefinementLevels<Dim>>;
67 : using simple_tags =
68 : tmpl::list<amr::Tags::ParentId<Dim>, amr::Tags::ChildIds<Dim>,
69 : amr::Tags::ParentMesh<Dim>,
70 : LinearSolver::Tags::ObservationId<OptionsGroup>,
71 : Tags::VolumeDataForOutput<OptionsGroup, FieldsTag>>;
72 : using compute_tags = tmpl::list<>;
73 : using const_global_cache_tags =
74 : tmpl::list<Tags::InitialCoarseLevels<OptionsGroup>,
75 : LinearSolver::Tags::OutputVolumeData<OptionsGroup>>;
76 :
77 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
78 : typename ActionList, typename ParallelComponent>
79 : static Parallel::iterable_action_return_t apply(
80 : db::DataBox<DbTagsList>& box,
81 : const tuples::TaggedTuple<InboxTags...>& /*inboxes*/,
82 : const Parallel::GlobalCache<Metavariables>& /*cache*/,
83 : const ElementId<Dim>& /*element_id*/, const ActionList /*meta*/,
84 : const ParallelComponent* const /*meta*/) {
85 : db::mutate_apply<InitializeElement>(make_not_null(&box));
86 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
87 : }
88 :
89 : public: // amr::protocols::Projector
90 : using argument_tags =
91 : tmpl::list<domain::Tags::Mesh<Dim>, domain::Tags::Element<Dim>,
92 : domain::Tags::InitialRefinementLevels<Dim>,
93 : LinearSolver::Tags::OutputVolumeData<OptionsGroup>>;
94 : using return_tags = tmpl::append<simple_tags, simple_tags_from_options>;
95 :
96 : template <typename... AmrData>
97 : static void apply(
98 : const gsl::not_null<std::optional<ElementId<Dim>>*> parent_id,
99 : const gsl::not_null<std::unordered_set<ElementId<Dim>>*> child_ids,
100 : const gsl::not_null<std::optional<Mesh<Dim>>*> parent_mesh,
101 : const gsl::not_null<size_t*> observation_id,
102 : const gsl::not_null<VolumeDataVars*> volume_data_for_output,
103 : const gsl::not_null<std::vector<std::array<size_t, Dim>>*>
104 : children_refinement_levels,
105 : const gsl::not_null<std::vector<std::array<size_t, Dim>>*>
106 : parent_refinement_levels,
107 : const Mesh<Dim>& mesh, const Element<Dim>& element,
108 : const std::vector<std::array<size_t, Dim>> initial_refinement_levels,
109 : const bool output_volume_data, const AmrData&... amr_data) {
110 : // Note: This initialization code runs on elements of the initial multigrid
111 : // hierarchy (created by LinearSolver::multigrid::ElementsAllocator) and on
112 : // elements created by AMR. Elements in a block of the initial domain are
113 : // assumed to have the same p-refinement.
114 :
115 : if constexpr (sizeof...(AmrData) == 0) {
116 : // Initialization: use initial domain to set up multigrid hierarchy
117 : const auto& element_id = element.id();
118 : const bool is_coarsest_grid =
119 : initial_refinement_levels == *parent_refinement_levels;
120 : const bool is_finest_grid =
121 : initial_refinement_levels == *children_refinement_levels;
122 : *parent_id = is_coarsest_grid
123 : ? std::nullopt
124 : : std::make_optional(multigrid::parent_id(element_id));
125 : *child_ids =
126 : is_finest_grid
127 : ? std::unordered_set<ElementId<Dim>>{}
128 : : multigrid::child_ids(
129 : element_id,
130 : (*children_refinement_levels)[element_id.block_id()]);
131 : *parent_mesh = is_coarsest_grid ? std::nullopt : std::make_optional(mesh);
132 : *observation_id = 0;
133 : } else {
134 : // These items are updated by AMR
135 : (void)parent_id;
136 : (void)child_ids;
137 : (void)parent_mesh;
138 : // These items are only needed during initialization
139 : (void)parent_refinement_levels;
140 : (void)children_refinement_levels;
141 : // Preserve state of observation ID
142 : if constexpr (tt::is_a_v<tuples::TaggedTuple, AmrData...>) {
143 : // h-refinement: copy from the parent
144 : *observation_id =
145 : get<LinearSolver::Tags::ObservationId<OptionsGroup>>(amr_data...);
146 : } else if constexpr (tt::is_a_v<std::unordered_map, AmrData...>) {
147 : // h-coarsening: copy from one of the children (doesn't matter which)
148 : *observation_id = get<LinearSolver::Tags::ObservationId<OptionsGroup>>(
149 : amr_data.begin()->second...);
150 : } else {
151 : (void)observation_id;
152 : }
153 : }
154 : // Initialize volume data output
155 : if (output_volume_data) {
156 : volume_data_for_output->initialize(mesh.number_of_grid_points());
157 : }
158 : }
159 : };
160 :
161 : // These two actions communicate and project the residual from the finer grid to
162 : // the coarser grid, storing it in the `SourceTag` on the coarser grid.
163 : template <typename FieldsTag, typename OptionsGroup,
164 : typename ResidualIsMassiveTag, typename SourceTag>
165 : using SendResidualToCoarserGrid = Actions::SendFieldsToCoarserGrid<
166 : tmpl::list<db::add_tag_prefix<LinearSolver::Tags::Residual, FieldsTag>>,
167 : OptionsGroup, ResidualIsMassiveTag, tmpl::list<SourceTag>>;
168 :
169 : template <size_t Dim, typename FieldsTag, typename OptionsGroup,
170 : typename SourceTag>
171 : using ReceiveResidualFromFinerGrid = Actions::ReceiveFieldsFromFinerGrid<
172 : Dim,
173 : tmpl::list<db::add_tag_prefix<LinearSolver::Tags::Residual, FieldsTag>>,
174 : OptionsGroup, tmpl::list<SourceTag>>;
175 :
176 : // Once the residual from the finer grid has been received and stored in the
177 : // `SourceTag`, this action prepares the pre-smoothing that will determine
178 : // an approximate solution on this grid. The pre-smoother is a separate
179 : // linear solver that runs independently after this action.
180 : template <typename FieldsTag, typename OptionsGroup, typename SourceTag,
181 : bool EnableBottomSolver>
182 : struct PreparePreSmoothing {
183 : private:
184 : using fields_tag = FieldsTag;
185 : using operator_applied_to_fields_tag =
186 : db::add_tag_prefix<LinearSolver::Tags::OperatorAppliedTo, fields_tag>;
187 : using source_tag = SourceTag;
188 :
189 : public:
190 : using const_global_cache_tags = tmpl::append<
191 : tmpl::list<
192 : LinearSolver::multigrid::Tags::EnablePreSmoothing<OptionsGroup>>,
193 : tmpl::conditional_t<EnableBottomSolver,
194 : tmpl::list<LinearSolver::multigrid::Tags::
195 : UseBottomSolver<OptionsGroup>>,
196 : tmpl::list<>>>;
197 :
198 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
199 : size_t Dim, typename ActionList, typename ParallelComponent>
200 : static Parallel::iterable_action_return_t apply(
201 : db::DataBox<DbTagsList>& box,
202 : tuples::TaggedTuple<InboxTags...>& /*inboxes*/,
203 : const Parallel::GlobalCache<Metavariables>& /*cache*/,
204 : const ElementId<Dim>& element_id, const ActionList /*meta*/,
205 : const ParallelComponent* const /*meta*/) {
206 : const size_t iteration_id =
207 : db::get<Convergence::Tags::IterationId<OptionsGroup>>(box);
208 : const bool is_coarsest_grid =
209 : not db::get<::amr::Tags::ParentId<Dim>>(box).has_value();
210 : if (UNLIKELY(db::get<logging::Tags::Verbosity<OptionsGroup>>(box) >=
211 : ::Verbosity::Debug)) {
212 : Parallel::printf("%s %s(%zu): Prepare %s\n", element_id,
213 : pretty_type::name<OptionsGroup>(), iteration_id,
214 : (EnableBottomSolver and is_coarsest_grid)
215 : ? "bottom-solver"
216 : : "pre-smoothing");
217 : }
218 :
219 : // On coarser grids the smoother solves for a correction to the finer-grid
220 : // fields, so we set its initial guess to zero. On the finest grid we smooth
221 : // the fields directly, so there's nothing to prepare.
222 : const bool is_finest_grid = db::get<amr::Tags::ChildIds<Dim>>(box).empty();
223 : if (not is_finest_grid) {
224 : db::mutate<fields_tag, operator_applied_to_fields_tag>(
225 : [](const auto fields, const auto operator_applied_to_fields,
226 : const auto& source) {
227 : *fields = make_with_value<typename fields_tag::type>(source, 0.);
228 : // We can set the linear operator applied to the initial fields to
229 : // zero as well, since it's linear
230 : *operator_applied_to_fields =
231 : make_with_value<typename operator_applied_to_fields_tag::type>(
232 : source, 0.);
233 : },
234 : make_not_null(&box), db::get<source_tag>(box));
235 : }
236 :
237 : // Record pre-smoothing initial fields and source
238 : if (db::get<LinearSolver::Tags::OutputVolumeData<OptionsGroup>>(box)) {
239 : db::mutate<Tags::VolumeDataForOutput<OptionsGroup, FieldsTag>>(
240 : [](const auto volume_data, const auto& initial_fields,
241 : const auto& source) {
242 : volume_data->assign_subset(
243 : Variables<db::wrap_tags_in<Tags::PreSmoothingInitial,
244 : typename fields_tag::tags_list>>(
245 : initial_fields));
246 : volume_data->assign_subset(
247 : Variables<db::wrap_tags_in<Tags::PreSmoothingSource,
248 : typename fields_tag::tags_list>>(
249 : source));
250 : },
251 : make_not_null(&box), db::get<fields_tag>(box),
252 : db::get<source_tag>(box));
253 : }
254 :
255 : // Skip pre-smoothing, if requested, or use bottom solver
256 : if constexpr (EnableBottomSolver) {
257 : if (db::get<LinearSolver::multigrid::Tags::UseBottomSolver<OptionsGroup>>(
258 : box) and
259 : is_coarsest_grid) {
260 : const size_t bottom_solver_index =
261 : tmpl::index_of<ActionList, SkipBottomSolver<FieldsTag, OptionsGroup,
262 : SourceTag>>::value +
263 : 1;
264 : return {Parallel::AlgorithmExecution::Continue, bottom_solver_index};
265 : }
266 : }
267 : const size_t first_action_after_pre_smoothing_index = tmpl::index_of<
268 : ActionList,
269 : SkipPostSmoothingAtBottom<FieldsTag, OptionsGroup, SourceTag>>::value;
270 : const size_t this_action_index =
271 : tmpl::index_of<ActionList, PreparePreSmoothing>::value;
272 : return {
273 : Parallel::AlgorithmExecution::Continue,
274 : db::get<
275 : LinearSolver::multigrid::Tags::EnablePreSmoothing<OptionsGroup>>(
276 : box)
277 : ? (this_action_index + 1)
278 : : first_action_after_pre_smoothing_index};
279 : }
280 : };
281 :
282 : // Once pre-smoothing is done, skip the bottom solver that comes next in the
283 : // action list. On the coarsest grid we directly jump to the bottom solver.
284 : template <typename FieldsTag, typename OptionsGroup, typename SourceTag>
285 : struct SkipBottomSolver {
286 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
287 : size_t Dim, typename ActionList, typename ParallelComponent>
288 : static Parallel::iterable_action_return_t apply(
289 : db::DataBox<DbTagsList>& /*box*/,
290 : const tuples::TaggedTuple<InboxTags...>& /*inboxes*/,
291 : const Parallel::GlobalCache<Metavariables>& /*cache*/,
292 : const ElementId<Dim>& /*element_id*/, const ActionList /*meta*/,
293 : const ParallelComponent* const /*meta*/) {
294 : const size_t first_action_after_bottom_solver_index = tmpl::index_of<
295 : ActionList,
296 : SkipPostSmoothingAtBottom<FieldsTag, OptionsGroup, SourceTag>>::value;
297 : return {Parallel::AlgorithmExecution::Continue,
298 : first_action_after_bottom_solver_index};
299 : }
300 : };
301 :
302 : // Once the pre-smoothing is done, we skip the second smoothing step on the
303 : // coarsest grid, i.e. at the "tip" of the V-cycle.
304 : template <typename FieldsTag, typename OptionsGroup, typename SourceTag>
305 : struct SkipPostSmoothingAtBottom {
306 : private:
307 : using fields_tag = FieldsTag;
308 : using residual_tag =
309 : db::add_tag_prefix<LinearSolver::Tags::Residual, fields_tag>;
310 :
311 : public:
312 : using const_global_cache_tags = tmpl::list<
313 : LinearSolver::multigrid::Tags::EnablePostSmoothingAtBottom<OptionsGroup>>;
314 :
315 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
316 : size_t Dim, typename ActionList, typename ParallelComponent>
317 : static Parallel::iterable_action_return_t apply(
318 : db::DataBox<DbTagsList>& box,
319 : const tuples::TaggedTuple<InboxTags...>& /*inboxes*/,
320 : const Parallel::GlobalCache<Metavariables>& /*cache*/,
321 : const ElementId<Dim>& /*element_id*/, const ActionList /*meta*/,
322 : const ParallelComponent* const /*meta*/) {
323 : const bool is_coarsest_grid =
324 : not db::get<amr::Tags::ParentId<Dim>>(box).has_value();
325 :
326 : // Record pre-smoothing result fields and residual
327 : if (db::get<LinearSolver::Tags::OutputVolumeData<OptionsGroup>>(box)) {
328 : db::mutate<Tags::VolumeDataForOutput<OptionsGroup, FieldsTag>>(
329 : [](const auto volume_data, const auto& result_fields,
330 : const auto& residuals) {
331 : volume_data->assign_subset(
332 : Variables<db::wrap_tags_in<Tags::PreSmoothingResult,
333 : typename fields_tag::tags_list>>(
334 : result_fields));
335 : volume_data->assign_subset(
336 : Variables<db::wrap_tags_in<Tags::PreSmoothingResidual,
337 : typename fields_tag::tags_list>>(
338 : residuals));
339 : },
340 : make_not_null(&box), db::get<fields_tag>(box),
341 : db::get<residual_tag>(box));
342 : }
343 :
344 : // Skip post-smoothing on the coarsest grid, if requested
345 : const size_t first_action_after_post_smoothing_index = tmpl::index_of<
346 : ActionList,
347 : SendCorrectionToFinerGrid<FieldsTag, OptionsGroup, SourceTag>>::value;
348 : const size_t post_smoothing_begin_index =
349 : tmpl::index_of<ActionList,
350 : ::Actions::Label<PostSmoothingBeginLabel>>::value +
351 : 1;
352 : const size_t this_action_index =
353 : tmpl::index_of<ActionList, SkipPostSmoothingAtBottom>::value;
354 : return {Parallel::AlgorithmExecution::Continue,
355 : is_coarsest_grid
356 : ? (db::get<LinearSolver::multigrid::Tags::
357 : EnablePostSmoothingAtBottom<OptionsGroup>>(box)
358 : ? post_smoothing_begin_index
359 : : first_action_after_post_smoothing_index)
360 : : (this_action_index + 1)};
361 : }
362 : };
363 :
364 : template <typename FieldsTag>
365 : struct CorrectionInboxTag
366 : : public Parallel::InboxInserters::Value<CorrectionInboxTag<FieldsTag>> {
367 : using temporal_id = size_t;
368 : using type = std::map<temporal_id, typename FieldsTag::type>;
369 : };
370 :
371 : // The next two actions communicate and project the coarse-grid correction, i.e.
372 : // the solution of the post-smoother, to the finer grid. The post-smoother on
373 : // finer grids runs after receiving this coarse-grid correction. Since the
374 : // post-smoother is skipped on the coarsest level, it directly sends the
375 : // solution of the pre-smoother to the finer grid, thus kicking off the
376 : // "ascending" branch of the V-cycle.
377 : template <typename FieldsTag, typename OptionsGroup, typename SourceTag>
378 : struct SendCorrectionToFinerGrid {
379 : private:
380 : using fields_tag = FieldsTag;
381 : using residual_tag =
382 : db::add_tag_prefix<LinearSolver::Tags::Residual, fields_tag>;
383 :
384 : public:
385 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
386 : size_t Dim, typename ActionList, typename ParallelComponent>
387 : static Parallel::iterable_action_return_t apply(
388 : db::DataBox<DbTagsList>& box,
389 : const tuples::TaggedTuple<InboxTags...>& /*inboxes*/,
390 : Parallel::GlobalCache<Metavariables>& cache,
391 : const ElementId<Dim>& element_id, const ActionList /*meta*/,
392 : const ParallelComponent* const /*meta*/) {
393 : const auto& child_ids = db::get<amr::Tags::ChildIds<Dim>>(box);
394 :
395 : // Record post-smoothing result fields and residual
396 : if (db::get<LinearSolver::Tags::OutputVolumeData<OptionsGroup>>(box)) {
397 : db::mutate<Tags::VolumeDataForOutput<OptionsGroup, FieldsTag>>(
398 : [](const auto volume_data, const auto& result_fields,
399 : const auto& residuals) {
400 : volume_data->assign_subset(
401 : Variables<db::wrap_tags_in<Tags::PostSmoothingResult,
402 : typename fields_tag::tags_list>>(
403 : result_fields));
404 : volume_data->assign_subset(
405 : Variables<db::wrap_tags_in<Tags::PostSmoothingResidual,
406 : typename fields_tag::tags_list>>(
407 : residuals));
408 : },
409 : make_not_null(&box), db::get<fields_tag>(box),
410 : db::get<residual_tag>(box));
411 : }
412 :
413 : if (child_ids.empty()) {
414 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
415 : }
416 :
417 : const size_t iteration_id =
418 : db::get<Convergence::Tags::IterationId<OptionsGroup>>(box);
419 : if (UNLIKELY(db::get<logging::Tags::Verbosity<OptionsGroup>>(box) >=
420 : ::Verbosity::Debug)) {
421 : Parallel::printf("%s %s(%zu): Send correction to children\n", element_id,
422 : pretty_type::name<OptionsGroup>(), iteration_id);
423 : }
424 :
425 : // Send a copy of the correction to all children
426 : auto& receiver_proxy =
427 : Parallel::get_parallel_component<ParallelComponent>(cache);
428 : for (const auto& child_id : child_ids) {
429 : auto coarse_grid_correction = db::get<fields_tag>(box);
430 : Parallel::receive_data<CorrectionInboxTag<FieldsTag>>(
431 : receiver_proxy[child_id], iteration_id,
432 : std::move(coarse_grid_correction));
433 : }
434 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
435 : }
436 : };
437 :
438 : template <size_t Dim, typename FieldsTag, typename OptionsGroup,
439 : typename SourceTag>
440 : struct ReceiveCorrectionFromCoarserGrid {
441 : private:
442 : using fields_tag = FieldsTag;
443 : using source_tag = SourceTag;
444 :
445 : public:
446 : using inbox_tags = tmpl::list<CorrectionInboxTag<FieldsTag>>;
447 :
448 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
449 : typename ActionList, typename ParallelComponent>
450 : static Parallel::iterable_action_return_t apply(
451 : db::DataBox<DbTagsList>& box, tuples::TaggedTuple<InboxTags...>& inboxes,
452 : const Parallel::GlobalCache<Metavariables>& /*cache*/,
453 : const ElementId<Dim>& element_id, const ActionList /*meta*/,
454 : const ParallelComponent* const /*meta*/) {
455 : const auto& parent_id = db::get<amr::Tags::ParentId<Dim>>(box);
456 : // We should always have a `parent_id` at this point because we skip this
457 : // part of the algorithm on the coarsest grid with the
458 : // `SkipPostSmoothingAtBottom` action
459 : ASSERT(parent_id.has_value(),
460 : "Trying to receive data from parent but no parent is set on element "
461 : << element_id << ".");
462 : const size_t iteration_id =
463 : db::get<Convergence::Tags::IterationId<OptionsGroup>>(box);
464 :
465 : // Wait for data from coarser grid
466 : auto& inbox = tuples::get<CorrectionInboxTag<FieldsTag>>(inboxes);
467 : if (inbox.find(iteration_id) == inbox.end()) {
468 : return {Parallel::AlgorithmExecution::Retry, std::nullopt};
469 : }
470 : auto parent_correction = std::move(inbox.extract(iteration_id).mapped());
471 :
472 : if (UNLIKELY(db::get<logging::Tags::Verbosity<OptionsGroup>>(box) >=
473 : ::Verbosity::Debug)) {
474 : Parallel::printf("%s %s(%zu): Prolongate correction from parent\n",
475 : element_id, pretty_type::name<OptionsGroup>(),
476 : iteration_id);
477 : }
478 :
479 : // Apply prolongation operator
480 : const auto& mesh = db::get<domain::Tags::Mesh<Dim>>(box);
481 : const auto& parent_mesh = db::get<amr::Tags::ParentMesh<Dim>>(box);
482 : ASSERT(
483 : parent_mesh.has_value(),
484 : "Should have a parent mesh, because a parent ID is set. This element: "
485 : << element_id << ", parent element: " << *parent_id);
486 : const auto child_size =
487 : domain::child_size(element_id.segment_ids(), parent_id->segment_ids());
488 : const auto prolongated_parent_correction =
489 : [&parent_correction, &parent_mesh, &mesh, &child_size]() {
490 : if (Spectral::needs_projection(*parent_mesh, mesh, child_size)) {
491 : return Spectral::project(
492 : parent_correction, *parent_mesh, mesh,
493 : make_array<Dim>(Spectral::SegmentSize::Full), child_size);
494 : } else {
495 : return std::move(parent_correction);
496 : }
497 : }();
498 :
499 : // Add correction to the solution on this grid
500 : db::mutate<fields_tag>(
501 : [&prolongated_parent_correction](const auto fields) {
502 : *fields += prolongated_parent_correction;
503 : },
504 : make_not_null(&box));
505 :
506 : // Record post-smoothing initial fields and source
507 : if (db::get<LinearSolver::Tags::OutputVolumeData<OptionsGroup>>(box)) {
508 : db::mutate<Tags::VolumeDataForOutput<OptionsGroup, FieldsTag>>(
509 : [](const auto volume_data, const auto& initial_fields,
510 : const auto& source) {
511 : volume_data->assign_subset(
512 : Variables<db::wrap_tags_in<Tags::PostSmoothingInitial,
513 : typename fields_tag::tags_list>>(
514 : initial_fields));
515 : volume_data->assign_subset(
516 : Variables<db::wrap_tags_in<Tags::PostSmoothingSource,
517 : typename fields_tag::tags_list>>(
518 : source));
519 : },
520 : make_not_null(&box), db::get<fields_tag>(box),
521 : db::get<source_tag>(box));
522 : }
523 :
524 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
525 : }
526 : };
527 :
528 : } // namespace LinearSolver::multigrid::detail
|
