Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <atomic>
7 : #include <cstddef>
8 : #include <limits>
9 : #include <map>
10 : #include <mutex>
11 : #include <optional>
12 : #include <tuple>
13 : #include <type_traits>
14 : #include <utility>
15 : #include <vector>
16 :
17 : #include "DataStructures/DataBox/AsAccess.hpp"
18 : #include "DataStructures/DataBox/DataBox.hpp"
19 : #include "DataStructures/DataBox/PrefixHelpers.hpp"
20 : #include "DataStructures/DataBox/Prefixes.hpp"
21 : #include "DataStructures/Tensor/EagerMath/Magnitude.hpp"
22 : #include "Domain/FaceNormal.hpp"
23 : #include "Domain/Structure/DirectionalIdMap.hpp"
24 : #include "Domain/Structure/Element.hpp"
25 : #include "Domain/Structure/ElementId.hpp"
26 : #include "Domain/Structure/Topology.hpp"
27 : #include "Domain/Structure/TrimMap.hpp"
28 : #include "Domain/Tags.hpp"
29 : #include "Domain/Tags/NeighborMesh.hpp"
30 : #include "Evolution/BoundaryCorrection.hpp"
31 : #include "Evolution/BoundaryCorrectionTags.hpp"
32 : #include "Evolution/DiscontinuousGalerkin/BoundaryData.hpp"
33 : #include "Evolution/DiscontinuousGalerkin/InboxTags.hpp"
34 : #include "Evolution/DiscontinuousGalerkin/MortarData.hpp"
35 : #include "Evolution/DiscontinuousGalerkin/MortarDataHolder.hpp"
36 : #include "Evolution/DiscontinuousGalerkin/MortarTags.hpp"
37 : #include "Evolution/DiscontinuousGalerkin/NormalVectorTags.hpp"
38 : #include "Evolution/DiscontinuousGalerkin/UsingSubcell.hpp"
39 : #include "NumericalAlgorithms/DiscontinuousGalerkin/Formulation.hpp"
40 : #include "NumericalAlgorithms/DiscontinuousGalerkin/LiftFlux.hpp"
41 : #include "NumericalAlgorithms/DiscontinuousGalerkin/LiftFromBoundary.hpp"
42 : #include "NumericalAlgorithms/DiscontinuousGalerkin/MortarHelpers.hpp"
43 : #include "NumericalAlgorithms/DiscontinuousGalerkin/Tags/Formulation.hpp"
44 : #include "NumericalAlgorithms/Spectral/BoundaryInterpolationMatrices.hpp"
45 : #include "NumericalAlgorithms/Spectral/Mesh.hpp"
46 : #include "NumericalAlgorithms/Spectral/Quadrature.hpp"
47 : #include "NumericalAlgorithms/Spectral/SegmentSize.hpp"
48 : #include "Parallel/AlgorithmExecution.hpp"
49 : #include "Parallel/ArrayCollection/IsDgElementCollection.hpp"
50 : #include "Parallel/GlobalCache.hpp"
51 : #include "Time/BoundaryHistory.hpp"
52 : #include "Time/EvolutionOrdering.hpp"
53 : #include "Time/SelfStart.hpp"
54 : #include "Time/Time.hpp"
55 : #include "Time/TimeStepId.hpp"
56 : #include "Time/TimeSteppers/LtsTimeStepper.hpp"
57 : #include "Time/TimeSteppers/TimeStepper.hpp"
58 : #include "Utilities/Algorithm.hpp"
59 : #include "Utilities/CallWithDynamicType.hpp"
60 : #include "Utilities/ErrorHandling/Error.hpp"
61 : #include "Utilities/Gsl.hpp"
62 : #include "Utilities/MakeArray.hpp"
63 : #include "Utilities/TMPL.hpp"
64 : #include "Utilities/TaggedTuple.hpp"
65 :
66 : /// \cond
67 : namespace Tags {
68 : struct Time;
69 : struct TimeStep;
70 : struct TimeStepId;
71 : template <typename StepperInterface>
72 : struct TimeStepper;
73 : } // namespace Tags
74 :
75 : namespace evolution::dg::subcell {
76 : // We use a forward declaration instead of including a header file to avoid
77 : // coupling to the DG-subcell libraries for executables that don't use subcell.
78 : template <size_t VolumeDim, typename DgComputeSubcellNeighborPackagedData>
79 : void neighbor_reconstructed_face_solution(
80 : gsl::not_null<db::Access*> box,
81 : gsl::not_null<std::pair<
82 : TimeStepId,
83 : DirectionalIdMap<VolumeDim, evolution::dg::BoundaryData<VolumeDim>>>*>
84 : received_temporal_id_and_data);
85 : template <size_t Dim>
86 : void neighbor_tci_decision(
87 : gsl::not_null<db::Access*> box,
88 : const std::pair<TimeStepId,
89 : DirectionalIdMap<Dim, evolution::dg::BoundaryData<Dim>>>&
90 : received_temporal_id_and_data);
91 : } // namespace evolution::dg::subcell
92 : /// \endcond
93 :
94 : namespace evolution::dg {
95 : namespace detail {
96 : template <typename BoundaryCorrectionClass>
97 : struct get_dg_boundary_terms {
98 : using type = typename BoundaryCorrectionClass::dg_boundary_terms_volume_tags;
99 : };
100 :
101 : template <typename Tag, typename Type = db::const_item_type<Tag, tmpl::list<>>>
102 : struct TemporaryReference {
103 : using tag = Tag;
104 : using type = const Type&;
105 : };
106 : } // namespace detail
107 :
108 : /// Receive boundary data for global time-stepping. Returns true if
109 : /// all necessary data has been received.
110 : template <bool UseNodegroupDgElements, typename Metavariables,
111 : typename DbTagsList, typename... InboxTags>
112 1 : bool receive_boundary_data_global_time_stepping(
113 : const gsl::not_null<db::DataBox<DbTagsList>*> box,
114 : const gsl::not_null<tuples::TaggedTuple<InboxTags...>*> inboxes) {
115 : constexpr size_t volume_dim = Metavariables::system::volume_dim;
116 :
117 : const TimeStepId& temporal_id = get<::Tags::TimeStepId>(*box);
118 :
119 : const auto number_of_neighbors =
120 : db::get<domain::Tags::Element<volume_dim>>(*box).number_of_neighbors();
121 :
122 : auto& inbox =
123 : tuples::get<evolution::dg::Tags::BoundaryCorrectionAndGhostCellsInbox<
124 : volume_dim, UseNodegroupDgElements>>(*inboxes);
125 : collect_messages:
126 : inbox.collect_messages();
127 : const auto received_record = inbox.messages.find(temporal_id);
128 : if (received_record == inbox.messages.end()) {
129 : if (inbox.set_missing_messages(number_of_neighbors)) {
130 : // We've received new messages while this function was running,
131 : // so try again.
132 : goto collect_messages; // NOLINT(cppcoreguidelines-avoid-goto)
133 : }
134 : return false;
135 : }
136 : auto& received_neighbor_data = received_record->second;
137 : if (received_neighbor_data.size() != number_of_neighbors) {
138 : ASSERT(received_neighbor_data.size() < number_of_neighbors,
139 : "Received too many messages: " << received_neighbor_data);
140 : if (inbox.set_missing_messages(number_of_neighbors -
141 : received_neighbor_data.size())) {
142 : // We've received new messages while this function was running,
143 : // so try again.
144 : goto collect_messages; // NOLINT(cppcoreguidelines-avoid-goto)
145 : }
146 : return false;
147 : }
148 :
149 : std::pair received_temporal_id_and_data{temporal_id,
150 : std::move(received_neighbor_data)};
151 : inbox.messages.erase(received_record);
152 :
153 : // Move inbox contents into the DataBox
154 : if constexpr (using_subcell_v<Metavariables>) {
155 : evolution::dg::subcell::neighbor_reconstructed_face_solution<
156 : volume_dim, typename Metavariables::SubcellOptions::
157 : DgComputeSubcellNeighborPackagedData>(
158 : &db::as_access(*box), make_not_null(&received_temporal_id_and_data));
159 : evolution::dg::subcell::neighbor_tci_decision<volume_dim>(
160 : make_not_null(&db::as_access(*box)), received_temporal_id_and_data);
161 : }
162 :
163 : db::mutate<evolution::dg::Tags::MortarMesh<volume_dim>,
164 : evolution::dg::Tags::MortarData<volume_dim>,
165 : evolution::dg::Tags::MortarNextTemporalId<volume_dim>,
166 : domain::Tags::NeighborMesh<volume_dim>>(
167 : [&received_temporal_id_and_data](
168 : const gsl::not_null<
169 : DirectionalIdMap<volume_dim, Mesh<volume_dim - 1>>*>
170 : mortar_meshes,
171 : const gsl::not_null<DirectionalIdMap<
172 : volume_dim, evolution::dg::MortarDataHolder<volume_dim>>*>
173 : mortar_data,
174 : const gsl::not_null<DirectionalIdMap<volume_dim, TimeStepId>*>
175 : mortar_next_time_step_id,
176 : const gsl::not_null<DirectionalIdMap<volume_dim, Mesh<volume_dim>>*>
177 : neighbor_mesh,
178 : const Mesh<volume_dim>& volume_mesh) {
179 : neighbor_mesh->clear();
180 : for (auto& received_mortar_data :
181 : received_temporal_id_and_data.second) {
182 : const auto& mortar_id = received_mortar_data.first;
183 : const size_t sliced_away_dim = mortar_id.direction().dimension();
184 : const Mesh<volume_dim - 1> face_mesh =
185 : volume_mesh.slice_away(sliced_away_dim);
186 : const Mesh<volume_dim - 1> neighbor_face_mesh =
187 : received_mortar_data.second.volume_mesh.slice_away(
188 : sliced_away_dim);
189 : const Mesh<volume_dim - 1> mortar_mesh =
190 : ::dg::mortar_mesh(face_mesh, neighbor_face_mesh);
191 : mortar_meshes->at(mortar_id) = mortar_mesh;
192 : p_project_mortar_data(
193 : make_not_null(&mortar_data->at(mortar_id).local()), mortar_mesh);
194 : neighbor_mesh->insert_or_assign(
195 : mortar_id, received_mortar_data.second.volume_mesh);
196 : mortar_next_time_step_id->at(mortar_id) =
197 : received_mortar_data.second.validity_range;
198 : ASSERT(using_subcell_v<Metavariables> or
199 : received_mortar_data.second.boundary_correction_data
200 : .has_value(),
201 : "Must receive number boundary correction data when not using "
202 : "DG-subcell. Mortar ID is: ("
203 : << mortar_id.direction() << "," << mortar_id.id()
204 : << ") and TimeStepId is "
205 : << received_temporal_id_and_data.first);
206 : if (received_mortar_data.second.boundary_correction_data
207 : .has_value()) {
208 : mortar_data->at(mortar_id).neighbor().face_mesh =
209 : neighbor_face_mesh;
210 : mortar_data->at(mortar_id).neighbor().mortar_mesh =
211 : received_mortar_data.second.boundary_correction_mesh.value();
212 : mortar_data->at(mortar_id).neighbor().mortar_data = std::move(
213 : received_mortar_data.second.boundary_correction_data.value());
214 : p_project_mortar_data(
215 : make_not_null(&mortar_data->at(mortar_id).neighbor()),
216 : mortar_mesh);
217 : }
218 : }
219 : },
220 : box, db::get<domain::Tags::Mesh<volume_dim>>(*box));
221 : return true;
222 : }
223 :
224 : /// Receive boundary data for local time-stepping. Returns true if
225 : /// all necessary data has been received.
226 : ///
227 : /// Setting \p DenseOutput to true receives data required for output
228 : /// at `::Tags::Time` instead of `::Tags::Next<::Tags::TimeStepId>`.
229 : template <bool UseNodegroupDgElements, typename System, size_t Dim,
230 : bool DenseOutput, typename DbTagsList, typename... InboxTags>
231 1 : bool receive_boundary_data_local_time_stepping(
232 : const gsl::not_null<db::DataBox<DbTagsList>*> box,
233 : const gsl::not_null<tuples::TaggedTuple<InboxTags...>*> inboxes) {
234 : using variables_tag = typename System::variables_tag;
235 : using dt_variables_tag = db::add_tag_prefix<::Tags::dt, variables_tag>;
236 :
237 : const auto needed_time = [&box]() {
238 : const LtsTimeStepper& time_stepper =
239 : db::get<::Tags::TimeStepper<LtsTimeStepper>>(*box);
240 : if constexpr (DenseOutput) {
241 : const auto& dense_output_time = db::get<::Tags::Time>(*box);
242 : return [&dense_output_time, &time_stepper](const TimeStepId& id) {
243 : return time_stepper.neighbor_data_required(dense_output_time, id);
244 : };
245 : } else {
246 : const auto& next_temporal_id =
247 : db::get<::Tags::Next<::Tags::TimeStepId>>(*box);
248 : return [&next_temporal_id, &time_stepper](const TimeStepId& id) {
249 : return time_stepper.neighbor_data_required(next_temporal_id, id);
250 : };
251 : }
252 : }();
253 :
254 : auto& inbox =
255 : tuples::get<evolution::dg::Tags::BoundaryCorrectionAndGhostCellsInbox<
256 : Dim, UseNodegroupDgElements>>(*inboxes);
257 :
258 : size_t missing_messages{};
259 : do {
260 : // The boundary history coupling computation (which computes the _lifted_
261 : // boundary correction) returns a Variables<dt<EvolvedVars>> instead of
262 : // using the `NormalDotNumericalFlux` prefix tag. This is because the
263 : // returned quantity is more a `dt` quantity than a
264 : // `NormalDotNormalDotFlux` since it's been lifted to the volume.
265 : using InboxMap =
266 : std::map<TimeStepId,
267 : DirectionalIdMap<Dim, evolution::dg::BoundaryData<Dim>>>;
268 : inbox.collect_messages();
269 : InboxMap& inbox_data = inbox.messages;
270 :
271 : missing_messages = 0;
272 :
273 : db::mutate<evolution::dg::Tags::MortarMesh<Dim>,
274 : evolution::dg::Tags::MortarDataHistory<
275 : Dim, typename dt_variables_tag::type>,
276 : evolution::dg::Tags::MortarNextTemporalId<Dim>,
277 : domain::Tags::NeighborMesh<Dim>>(
278 : [&inbox_data, &missing_messages, &needed_time](
279 : const gsl::not_null<DirectionalIdMap<Dim, Mesh<Dim - 1>>*>
280 : mortar_meshes,
281 : const gsl::not_null<
282 : DirectionalIdMap<Dim, TimeSteppers::BoundaryHistory<
283 : evolution::dg::MortarData<Dim>,
284 : evolution::dg::MortarData<Dim>,
285 : typename dt_variables_tag::type>>*>
286 : boundary_data_history,
287 : const gsl::not_null<DirectionalIdMap<Dim, TimeStepId>*>
288 : mortar_next_time_step_ids,
289 : const gsl::not_null<DirectionalIdMap<Dim, Mesh<Dim>>*>
290 : neighbor_mesh,
291 : const Element<Dim>& element, const Mesh<Dim>& volume_mesh) {
292 : // Remove neighbor meshes for neighbors that don't exist anymore
293 : domain::remove_nonexistent_neighbors(neighbor_mesh, element);
294 :
295 : // Move received boundary data into boundary history.
296 : for (auto& [mortar_id, mortar_next_time_step_id] :
297 : *mortar_next_time_step_ids) {
298 : if (mortar_id.id() == ElementId<Dim>::external_boundary_id()) {
299 : continue;
300 : }
301 : const size_t sliced_away_dim = mortar_id.direction().dimension();
302 : const Mesh<Dim - 1> face_mesh =
303 : volume_mesh.slice_away(sliced_away_dim);
304 : while (needed_time(mortar_next_time_step_id)) {
305 : const auto time_entry = inbox_data.find(mortar_next_time_step_id);
306 : if (time_entry == inbox_data.end()) {
307 : ++missing_messages;
308 : break;
309 : }
310 : const auto received_mortar_data =
311 : time_entry->second.find(mortar_id);
312 : if (received_mortar_data == time_entry->second.end()) {
313 : ++missing_messages;
314 : break;
315 : }
316 :
317 : const Mesh<Dim - 1> neighbor_face_mesh =
318 : received_mortar_data->second.volume_mesh.slice_away(
319 : sliced_away_dim);
320 : const Mesh<Dim - 1> mortar_mesh =
321 : ::dg::mortar_mesh(face_mesh, neighbor_face_mesh);
322 :
323 : const auto project_boundary_mortar_data =
324 : [&mortar_mesh](
325 : const TimeStepId& /*id*/,
326 : const gsl::not_null<::evolution::dg::MortarData<Dim>*>
327 : mortar_data) {
328 : return p_project_mortar_data(mortar_data, mortar_mesh);
329 : };
330 :
331 : mortar_meshes->at(mortar_id) = mortar_mesh;
332 : boundary_data_history->at(mortar_id).local().for_each(
333 : project_boundary_mortar_data);
334 :
335 : MortarData<Dim> neighbor_mortar_data{};
336 : // Insert:
337 : // - the current TimeStepId of the neighbor
338 : // - the current face mesh of the neighbor
339 : // - the current boundary correction data of the neighbor
340 : ASSERT(received_mortar_data->second.boundary_correction_data
341 : .has_value(),
342 : "Did not receive boundary correction data from the "
343 : "neighbor\nMortarId: "
344 : << mortar_id
345 : << "\nTimeStepId: " << mortar_next_time_step_id);
346 : neighbor_mesh->insert_or_assign(
347 : mortar_id, received_mortar_data->second.volume_mesh);
348 : neighbor_mortar_data.mortar_mesh =
349 : received_mortar_data->second.boundary_correction_mesh.value();
350 : neighbor_mortar_data.mortar_data =
351 : std::move(received_mortar_data->second
352 : .boundary_correction_data.value());
353 : boundary_data_history->at(mortar_id).remote().insert(
354 : time_entry->first,
355 : received_mortar_data->second.integration_order,
356 : std::move(neighbor_mortar_data));
357 : boundary_data_history->at(mortar_id).remote().for_each(
358 : project_boundary_mortar_data);
359 : mortar_next_time_step_id =
360 : received_mortar_data->second.validity_range;
361 : time_entry->second.erase(received_mortar_data);
362 : if (time_entry->second.empty()) {
363 : inbox_data.erase(time_entry);
364 : }
365 : }
366 : }
367 : },
368 : box, db::get<::domain::Tags::Element<Dim>>(*box),
369 : db::get<domain::Tags::Mesh<Dim>>(*box));
370 :
371 : if (missing_messages == 0) {
372 : return true;
373 : }
374 : } while (inbox.set_missing_messages(missing_messages));
375 : return false;
376 : }
377 :
378 : /// Apply corrections from boundary communication.
379 : ///
380 : /// If `LocalTimeStepping` is false, updates the derivative of the variables,
381 : /// which should be done before taking a time step. If
382 : /// `LocalTimeStepping` is true, updates the variables themselves, which should
383 : /// be done after the volume update.
384 : ///
385 : /// Setting \p DenseOutput to true receives data required for output
386 : /// at ::Tags::Time instead of performing a full step. This is only
387 : /// used for local time-stepping.
388 : template <bool LocalTimeStepping, typename Metavariables, size_t VolumeDim,
389 : bool DenseOutput>
390 1 : struct ApplyBoundaryCorrections {
391 0 : static constexpr bool local_time_stepping = LocalTimeStepping;
392 : static_assert(local_time_stepping or not DenseOutput,
393 : "GTS does not use ApplyBoundaryCorrections for dense output.");
394 :
395 0 : using system = typename Metavariables::system;
396 0 : static constexpr size_t volume_dim = VolumeDim;
397 0 : using variables_tag = typename system::variables_tag;
398 0 : using dt_variables_tag = db::add_tag_prefix<::Tags::dt, variables_tag>;
399 0 : using DtVariables = typename dt_variables_tag::type;
400 0 : using derived_boundary_corrections =
401 : tmpl::at<typename Metavariables::factory_creation::factory_classes,
402 : evolution::BoundaryCorrection>;
403 0 : using volume_tags_for_dg_boundary_terms = tmpl::remove_duplicates<
404 : tmpl::flatten<tmpl::transform<derived_boundary_corrections,
405 : detail::get_dg_boundary_terms<tmpl::_1>>>>;
406 :
407 0 : using TimeStepperType =
408 : tmpl::conditional_t<local_time_stepping, LtsTimeStepper, TimeStepper>;
409 :
410 0 : using tag_to_update =
411 : tmpl::conditional_t<local_time_stepping, variables_tag, dt_variables_tag>;
412 0 : using mortar_data_tag = tmpl::conditional_t<
413 : local_time_stepping,
414 : evolution::dg::Tags::MortarDataHistory<volume_dim, DtVariables>,
415 : evolution::dg::Tags::MortarData<volume_dim>>;
416 0 : using MortarDataType =
417 : tmpl::conditional_t<DenseOutput, const typename mortar_data_tag::type,
418 : typename mortar_data_tag::type>;
419 :
420 0 : using return_tags =
421 : tmpl::conditional_t<DenseOutput, tmpl::list<tag_to_update>,
422 : tmpl::list<tag_to_update, mortar_data_tag>>;
423 0 : using argument_tags = tmpl::append<
424 : tmpl::flatten<tmpl::list<
425 : tmpl::conditional_t<DenseOutput, mortar_data_tag, tmpl::list<>>,
426 : domain::Tags::Mesh<volume_dim>, domain::Tags::Element<volume_dim>,
427 : Tags::MortarMesh<volume_dim>, Tags::MortarInfo<volume_dim>,
428 : ::dg::Tags::Formulation,
429 : evolution::dg::Tags::NormalCovectorAndMagnitude<volume_dim>,
430 : ::Tags::TimeStepper<TimeStepperType>,
431 : evolution::Tags::BoundaryCorrection,
432 : tmpl::conditional_t<DenseOutput, ::Tags::Time, ::Tags::TimeStep>,
433 : tmpl::conditional_t<local_time_stepping, tmpl::list<>,
434 : domain::Tags::DetInvJacobian<
435 : Frame::ElementLogical, Frame::Inertial>>>>,
436 : volume_tags_for_dg_boundary_terms>;
437 :
438 : // full step
439 : template <typename... VolumeArgs>
440 0 : static void apply(
441 : const gsl::not_null<typename tag_to_update::type*> vars_to_update,
442 : const gsl::not_null<MortarDataType*> mortar_data,
443 : const Mesh<volume_dim>& volume_mesh, const Element<volume_dim>& element,
444 : const typename Tags::MortarMesh<volume_dim>::type& mortar_meshes,
445 : const typename Tags::MortarInfo<volume_dim>::type& mortar_infos,
446 : const ::dg::Formulation dg_formulation,
447 : const DirectionMap<
448 : volume_dim, std::optional<Variables<tmpl::list<
449 : evolution::dg::Tags::MagnitudeOfNormal,
450 : evolution::dg::Tags::NormalCovector<volume_dim>>>>>&
451 : face_normal_covector_and_magnitude,
452 : const TimeStepperType& time_stepper,
453 : const evolution::BoundaryCorrection& boundary_correction,
454 : const TimeDelta& time_step,
455 : const Scalar<DataVector>& gts_det_inv_jacobian,
456 : const VolumeArgs&... volume_args) {
457 : apply_impl(vars_to_update, mortar_data, volume_mesh, element, mortar_meshes,
458 : mortar_infos, dg_formulation, face_normal_covector_and_magnitude,
459 : time_stepper, boundary_correction, time_step,
460 : std::numeric_limits<double>::signaling_NaN(),
461 : gts_det_inv_jacobian, volume_args...);
462 : }
463 :
464 : template <typename... VolumeArgs>
465 0 : static void apply(
466 : const gsl::not_null<typename tag_to_update::type*> vars_to_update,
467 : const gsl::not_null<MortarDataType*> mortar_data,
468 : const Mesh<volume_dim>& volume_mesh, const Element<volume_dim>& element,
469 : const typename Tags::MortarMesh<volume_dim>::type& mortar_meshes,
470 : const typename Tags::MortarInfo<volume_dim>::type& mortar_infos,
471 : const ::dg::Formulation dg_formulation,
472 : const DirectionMap<
473 : volume_dim, std::optional<Variables<tmpl::list<
474 : evolution::dg::Tags::MagnitudeOfNormal,
475 : evolution::dg::Tags::NormalCovector<volume_dim>>>>>&
476 : face_normal_covector_and_magnitude,
477 : const TimeStepperType& time_stepper,
478 : const evolution::BoundaryCorrection& boundary_correction,
479 : const TimeDelta& time_step, const VolumeArgs&... volume_args) {
480 : apply_impl(vars_to_update, mortar_data, volume_mesh, element, mortar_meshes,
481 : mortar_infos, dg_formulation, face_normal_covector_and_magnitude,
482 : time_stepper, boundary_correction, time_step,
483 : std::numeric_limits<double>::signaling_NaN(), {},
484 : volume_args...);
485 : }
486 :
487 : // dense output (LTS only)
488 : template <typename... VolumeArgs>
489 0 : static void apply(
490 : const gsl::not_null<typename variables_tag::type*> vars_to_update,
491 : const MortarDataType& mortar_data, const Mesh<volume_dim>& volume_mesh,
492 : const Element<volume_dim>& element,
493 : const typename Tags::MortarMesh<volume_dim>::type& mortar_meshes,
494 : const typename Tags::MortarInfo<volume_dim>::type& mortar_infos,
495 : const ::dg::Formulation dg_formulation,
496 : const DirectionMap<
497 : volume_dim, std::optional<Variables<tmpl::list<
498 : evolution::dg::Tags::MagnitudeOfNormal,
499 : evolution::dg::Tags::NormalCovector<volume_dim>>>>>&
500 : face_normal_covector_and_magnitude,
501 : const LtsTimeStepper& time_stepper,
502 : const evolution::BoundaryCorrection& boundary_correction,
503 : const double dense_output_time, const VolumeArgs&... volume_args) {
504 : apply_impl(vars_to_update, &mortar_data, volume_mesh, element,
505 : mortar_meshes, mortar_infos, dg_formulation,
506 : face_normal_covector_and_magnitude, time_stepper,
507 : boundary_correction, TimeDelta{}, dense_output_time, {},
508 : volume_args...);
509 : }
510 :
511 : template <typename DbTagsList, typename... InboxTags, typename ArrayIndex,
512 : typename ParallelComponent>
513 0 : static bool is_ready(
514 : const gsl::not_null<db::DataBox<DbTagsList>*> box,
515 : const gsl::not_null<tuples::TaggedTuple<InboxTags...>*> inboxes,
516 : Parallel::GlobalCache<Metavariables>& /*cache*/,
517 : const ArrayIndex& /*array_index*/,
518 : const ParallelComponent* const /*component*/) {
519 : if constexpr (local_time_stepping) {
520 : return receive_boundary_data_local_time_stepping<
521 : Parallel::is_dg_element_collection_v<ParallelComponent>, system,
522 : VolumeDim, DenseOutput>(box, inboxes);
523 : } else {
524 : return receive_boundary_data_global_time_stepping<
525 : Parallel::is_dg_element_collection_v<ParallelComponent>,
526 : Metavariables>(box, inboxes);
527 : }
528 : }
529 :
530 : private:
531 : template <typename... VolumeArgs>
532 0 : static void apply_impl(
533 : const gsl::not_null<typename tag_to_update::type*> vars_to_update,
534 : const gsl::not_null<MortarDataType*> mortar_data,
535 : const Mesh<volume_dim>& volume_mesh, const Element<volume_dim>& element,
536 : const typename Tags::MortarMesh<volume_dim>::type& mortar_meshes,
537 : const typename Tags::MortarInfo<volume_dim>::type& mortar_infos,
538 : const ::dg::Formulation dg_formulation,
539 : const DirectionMap<
540 : volume_dim, std::optional<Variables<tmpl::list<
541 : evolution::dg::Tags::MagnitudeOfNormal,
542 : evolution::dg::Tags::NormalCovector<volume_dim>>>>>&
543 : face_normal_covector_and_magnitude,
544 : const TimeStepperType& time_stepper,
545 : const evolution::BoundaryCorrection& boundary_correction,
546 : const TimeDelta& time_step, const double dense_output_time,
547 : const Scalar<DataVector>& gts_det_inv_jacobian,
548 : const VolumeArgs&... volume_args) {
549 : tuples::tagged_tuple_from_typelist<db::wrap_tags_in<
550 : detail::TemporaryReference, volume_tags_for_dg_boundary_terms>>
551 : volume_args_tuple{volume_args...};
552 :
553 : // Set up helper lambda that will compute and lift the boundary corrections
554 : ASSERT(
555 : volume_mesh.quadrature() ==
556 : make_array<volume_dim>(volume_mesh.quadrature(0)) or
557 : element.topologies() != domain::topologies::hypercube<volume_dim>,
558 : "Must have isotropic quadrature, but got volume mesh: " << volume_mesh);
559 : const bool using_gauss_lobatto_points =
560 : volume_mesh.quadrature(0) == Spectral::Quadrature::GaussLobatto;
561 :
562 : Scalar<DataVector> volume_det_inv_jacobian{};
563 : Scalar<DataVector> volume_det_jacobian{};
564 : if constexpr (not local_time_stepping) {
565 : if (not using_gauss_lobatto_points) {
566 : get(volume_det_inv_jacobian)
567 : .set_data_ref(make_not_null(
568 : // NOLINTNEXTLINE(cppcoreguidelines-pro-type-const-cast)
569 : &const_cast<DataVector&>(get(gts_det_inv_jacobian))));
570 : get(volume_det_jacobian) = 1.0 / get(volume_det_inv_jacobian);
571 : }
572 : }
573 :
574 : static_assert(
575 : tmpl::all<derived_boundary_corrections, std::is_final<tmpl::_1>>::value,
576 : "All createable classes for boundary corrections must be marked "
577 : "final.");
578 : call_with_dynamic_type<void, derived_boundary_corrections>(
579 : &boundary_correction,
580 : [&dense_output_time, &dg_formulation,
581 : &face_normal_covector_and_magnitude, &mortar_data, &mortar_meshes,
582 : &mortar_infos, &time_step, &time_stepper, using_gauss_lobatto_points,
583 : &vars_to_update, &volume_args_tuple, &volume_det_jacobian,
584 : &volume_det_inv_jacobian,
585 : &volume_mesh](auto* typed_boundary_correction) {
586 : using BcType = std::decay_t<decltype(*typed_boundary_correction)>;
587 : // Compute internal boundary quantities on the mortar for sides of
588 : // the element that have neighbors, i.e. they are not an external
589 : // side.
590 : using mortar_tags_list = typename BcType::dg_package_field_tags;
591 :
592 : // Variables for reusing allocations. The actual values are
593 : // not reused.
594 : DtVariables dt_boundary_correction_on_mortar{};
595 : DtVariables volume_dt_correction{};
596 : // These variables may change size for each mortar and require
597 : // a new memory allocation, but they may also happen to need
598 : // to be the same size twice in a row, in which case holding
599 : // on to the allocation is a win.
600 : Scalar<DataVector> face_det_jacobian{};
601 : Variables<mortar_tags_list> local_data_on_mortar{};
602 : Variables<mortar_tags_list> neighbor_data_on_mortar{};
603 :
604 : for (auto& mortar_id_and_data : *mortar_data) {
605 : const auto& mortar_id = mortar_id_and_data.first;
606 : const auto& direction = mortar_id.direction();
607 : if (UNLIKELY(mortar_id.id() ==
608 : ElementId<volume_dim>::external_boundary_id())) {
609 : ERROR(
610 : "Cannot impose boundary conditions on external boundary in "
611 : "direction "
612 : << direction
613 : << " in the ApplyBoundaryCorrections action. Boundary "
614 : "conditions are applied in the ComputeTimeDerivative "
615 : "action "
616 : "instead. You may have unintentionally added external "
617 : "mortars in one of the initialization actions.");
618 : }
619 :
620 : const Mesh<volume_dim - 1> face_mesh =
621 : volume_mesh.slice_away(direction.dimension());
622 :
623 : const auto compute_correction_coupling =
624 : [&typed_boundary_correction, &direction, dg_formulation,
625 : &dt_boundary_correction_on_mortar, &face_det_jacobian,
626 : &face_mesh, &face_normal_covector_and_magnitude,
627 : &local_data_on_mortar, &mortar_id, &mortar_meshes,
628 : &mortar_infos, &neighbor_data_on_mortar,
629 : using_gauss_lobatto_points, &volume_args_tuple,
630 : &volume_det_jacobian, &volume_det_inv_jacobian,
631 : &volume_dt_correction, &volume_mesh](
632 : const MortarData<volume_dim>& local_mortar_data,
633 : const MortarData<volume_dim>& neighbor_mortar_data)
634 : -> DtVariables {
635 : if (local_time_stepping and not using_gauss_lobatto_points) {
636 : // This needs to be updated every call because the Jacobian
637 : // may be time-dependent. In the case of time-independent maps
638 : // and local time stepping we could first perform the integral
639 : // on the boundaries, and then lift to the volume. This is
640 : // left as a future optimization.
641 : volume_det_inv_jacobian =
642 : local_mortar_data.volume_det_inv_jacobian.value();
643 : get(volume_det_jacobian) = 1.0 / get(volume_det_inv_jacobian);
644 : }
645 : const auto& mortar_mesh = mortar_meshes.at(mortar_id);
646 :
647 : // Extract local and neighbor data, copy into Variables because
648 : // we store them in a std::vector for type erasure.
649 : ASSERT(*local_mortar_data.mortar_mesh ==
650 : *neighbor_mortar_data.mortar_mesh and
651 : *local_mortar_data.mortar_mesh == mortar_mesh,
652 : "local mortar mesh: " << *local_mortar_data.mortar_mesh
653 : << "\nneighbor mortar mesh: "
654 : << *neighbor_mortar_data.mortar_mesh
655 : << "\nmortar mesh: " << mortar_mesh
656 : << "\n");
657 : const DataVector& local_data = *local_mortar_data.mortar_data;
658 : const DataVector& neighbor_data =
659 : *neighbor_mortar_data.mortar_data;
660 : ASSERT(local_data.size() == neighbor_data.size(),
661 : "local data size: "
662 : << local_data.size()
663 : << "\nneighbor_data: " << neighbor_data.size()
664 : << "\n mortar_mesh: " << mortar_mesh << "\n");
665 : ASSERT(local_data_on_mortar.number_of_grid_points() ==
666 : neighbor_data_on_mortar.number_of_grid_points(),
667 : "Local data size = "
668 : << local_data_on_mortar.number_of_grid_points()
669 : << ", but neighbor size = "
670 : << neighbor_data_on_mortar.number_of_grid_points());
671 : local_data_on_mortar.set_data_ref(
672 : // NOLINTNEXTLINE(cppcoreguidelines-pro-type-const-cast)
673 : const_cast<double*>(local_data.data()), local_data.size());
674 : neighbor_data_on_mortar.set_data_ref(
675 : // NOLINTNEXTLINE(cppcoreguidelines-pro-type-const-cast)
676 : const_cast<double*>(neighbor_data.data()),
677 : neighbor_data.size());
678 :
679 : // The boundary computations and lifting can be further
680 : // optimized by in the h-refinement case having only one
681 : // allocation for the face and having the projection from the
682 : // mortar to the face be done in place. E.g.
683 : // local_data_on_mortar and neighbor_data_on_mortar could be
684 : // allocated fewer times, as well as `needs_projection` section
685 : // below could do an in-place projection.
686 : dt_boundary_correction_on_mortar.initialize(
687 : mortar_mesh.number_of_grid_points());
688 :
689 : call_boundary_correction(
690 : make_not_null(&dt_boundary_correction_on_mortar),
691 : local_data_on_mortar, neighbor_data_on_mortar,
692 : *typed_boundary_correction, dg_formulation, volume_args_tuple,
693 : typename BcType::dg_boundary_terms_volume_tags{});
694 :
695 : const std::array<Spectral::SegmentSize, volume_dim - 1>&
696 : mortar_size = mortar_infos.at(mortar_id).mortar_size();
697 :
698 : // This cannot reuse an allocation because it is initialized
699 : // via move-assignment. (If it is used at all.)
700 : DtVariables dt_boundary_correction_projected_onto_face{};
701 : auto& dt_boundary_correction =
702 : [&dt_boundary_correction_on_mortar,
703 : &dt_boundary_correction_projected_onto_face, &face_mesh,
704 : &mortar_mesh, &mortar_size]() -> DtVariables& {
705 : if (Spectral::needs_projection(face_mesh, mortar_mesh,
706 : mortar_size)) {
707 : dt_boundary_correction_projected_onto_face =
708 : ::dg::project_from_mortar(
709 : dt_boundary_correction_on_mortar, face_mesh,
710 : mortar_mesh, mortar_size);
711 : return dt_boundary_correction_projected_onto_face;
712 : }
713 : return dt_boundary_correction_on_mortar;
714 : }();
715 :
716 : // Both paths initialize this to be non-owning.
717 : Scalar<DataVector> magnitude_of_face_normal{};
718 : if constexpr (local_time_stepping) {
719 : (void)face_normal_covector_and_magnitude;
720 : get(magnitude_of_face_normal)
721 : .set_data_ref(make_not_null(&const_cast<DataVector&>(
722 : get(local_mortar_data.face_normal_magnitude.value()))));
723 : } else {
724 : ASSERT(
725 : face_normal_covector_and_magnitude.count(direction) == 1 and
726 : face_normal_covector_and_magnitude.at(direction)
727 : .has_value(),
728 : "Face normal covector and magnitude not set in "
729 : "direction: "
730 : << direction);
731 : get(magnitude_of_face_normal)
732 : .set_data_ref(make_not_null(&const_cast<DataVector&>(
733 : get(get<evolution::dg::Tags::MagnitudeOfNormal>(
734 : *face_normal_covector_and_magnitude.at(
735 : direction))))));
736 : }
737 :
738 : if (using_gauss_lobatto_points) {
739 : // The lift_flux function lifts only on the slice, it does not
740 : // add the contribution to the volume.
741 : ::dg::lift_flux(make_not_null(&dt_boundary_correction),
742 : volume_mesh.extents(direction.dimension()),
743 : magnitude_of_face_normal);
744 : return std::move(dt_boundary_correction);
745 : } else {
746 : // We are using Gauss points.
747 : //
748 : // Notes:
749 : // - We should really lift both sides simultaneously since this
750 : // reduces memory accesses. Lifting all sides at the same
751 : // time is unlikely to improve performance since we lift by
752 : // jumping through slices. There may also be compatibility
753 : // issues with local time stepping.
754 : // - If we lift both sides at the same time we first need to
755 : // deal with projecting from mortars to the face, then lift
756 : // off the faces. With non-owning Variables memory
757 : // allocations could be significantly reduced in this code.
758 : if constexpr (local_time_stepping) {
759 : ASSERT(get(volume_det_inv_jacobian).size() > 0,
760 : "For local time stepping the volume determinant of "
761 : "the inverse Jacobian has not been set.");
762 :
763 : get(face_det_jacobian)
764 : .set_data_ref(make_not_null(&const_cast<DataVector&>(
765 : get(local_mortar_data.face_det_jacobian.value()))));
766 : } else {
767 : // Project the determinant of the Jacobian to the face. This
768 : // could be optimized by caching in the time-independent case.
769 : get(face_det_jacobian)
770 : .destructive_resize(face_mesh.number_of_grid_points());
771 : const Matrix identity{};
772 : auto interpolation_matrices =
773 : make_array<volume_dim>(std::cref(identity));
774 : const std::pair<Matrix, Matrix>& matrices =
775 : Spectral::boundary_interpolation_matrices(
776 : volume_mesh.slice_through(direction.dimension()));
777 : gsl::at(interpolation_matrices, direction.dimension()) =
778 : direction.side() == Side::Upper ? matrices.second
779 : : matrices.first;
780 : apply_matrices(make_not_null(&get(face_det_jacobian)),
781 : interpolation_matrices,
782 : get(volume_det_jacobian),
783 : volume_mesh.extents());
784 : }
785 :
786 : volume_dt_correction.initialize(
787 : volume_mesh.number_of_grid_points(), 0.0);
788 : ::dg::lift_boundary_terms_gauss_points(
789 : make_not_null(&volume_dt_correction),
790 : volume_det_inv_jacobian, volume_mesh, direction,
791 : dt_boundary_correction, magnitude_of_face_normal,
792 : face_det_jacobian);
793 : return std::move(volume_dt_correction);
794 : }
795 : };
796 :
797 : if constexpr (local_time_stepping) {
798 : typename variables_tag::type lgl_lifted_data{};
799 : auto& lifted_data = using_gauss_lobatto_points ? lgl_lifted_data
800 : : *vars_to_update;
801 : if (using_gauss_lobatto_points) {
802 : lifted_data.initialize(face_mesh.number_of_grid_points(), 0.0);
803 : }
804 :
805 : auto& mortar_data_history = mortar_id_and_data.second;
806 : if constexpr (DenseOutput) {
807 : (void)time_step;
808 : time_stepper.boundary_dense_output(
809 : &lifted_data, mortar_data_history, dense_output_time,
810 : compute_correction_coupling);
811 : } else {
812 : (void)dense_output_time;
813 : time_stepper.add_boundary_delta(&lifted_data,
814 : mortar_data_history, time_step,
815 : compute_correction_coupling);
816 : }
817 :
818 : if (using_gauss_lobatto_points) {
819 : // Add the flux contribution to the volume data
820 : add_slice_to_data(
821 : vars_to_update, lifted_data, volume_mesh.extents(),
822 : direction.dimension(),
823 : index_to_slice_at(volume_mesh.extents(), direction));
824 : }
825 : } else {
826 : (void)time_step;
827 : (void)time_stepper;
828 : (void)dense_output_time;
829 :
830 : // Choose an allocation cache that may be empty, so we
831 : // might be able to reuse the allocation obtained for the
832 : // lifted data. This may result in a self assignment,
833 : // depending on the code paths taken, but handling the
834 : // results this way makes the GTS and LTS paths more
835 : // similar because the LTS code always stores the result
836 : // in the history and so sometimes benefits from moving
837 : // into the return value of compute_correction_coupling.
838 : auto& lifted_data = using_gauss_lobatto_points
839 : ? dt_boundary_correction_on_mortar
840 : : volume_dt_correction;
841 : lifted_data = compute_correction_coupling(
842 : mortar_id_and_data.second.local(),
843 : mortar_id_and_data.second.neighbor());
844 :
845 : if (using_gauss_lobatto_points) {
846 : // Add the flux contribution to the volume data
847 : add_slice_to_data(
848 : vars_to_update, lifted_data, volume_mesh.extents(),
849 : direction.dimension(),
850 : index_to_slice_at(volume_mesh.extents(), direction));
851 : } else {
852 : *vars_to_update += lifted_data;
853 : }
854 : }
855 : }
856 : });
857 : }
858 :
859 : template <typename... BoundaryCorrectionTags, typename... Tags,
860 : typename BoundaryCorrection, typename... AllVolumeArgs,
861 : typename... VolumeTagsForCorrection>
862 0 : static void call_boundary_correction(
863 : const gsl::not_null<Variables<tmpl::list<BoundaryCorrectionTags...>>*>
864 : boundary_corrections_on_mortar,
865 : const Variables<tmpl::list<Tags...>>& local_boundary_data,
866 : const Variables<tmpl::list<Tags...>>& neighbor_boundary_data,
867 : const BoundaryCorrection& boundary_correction,
868 : const ::dg::Formulation dg_formulation,
869 : const tuples::TaggedTuple<detail::TemporaryReference<AllVolumeArgs>...>&
870 : volume_args_tuple,
871 : tmpl::list<VolumeTagsForCorrection...> /*meta*/) {
872 : boundary_correction.dg_boundary_terms(
873 : make_not_null(
874 : &get<BoundaryCorrectionTags>(*boundary_corrections_on_mortar))...,
875 : get<Tags>(local_boundary_data)..., get<Tags>(neighbor_boundary_data)...,
876 : dg_formulation,
877 : tuples::get<detail::TemporaryReference<VolumeTagsForCorrection>>(
878 : volume_args_tuple)...);
879 : }
880 : };
881 :
882 1 : namespace Actions {
883 : /*!
884 : * \brief Computes the boundary corrections for global time-stepping
885 : * and adds them to the time derivative.
886 : */
887 : template <size_t VolumeDim, bool DenseOutput, bool UseNodegroupDgElements>
888 1 : struct ApplyBoundaryCorrectionsToTimeDerivative {
889 0 : using inbox_tags =
890 : tmpl::list<evolution::dg::Tags::BoundaryCorrectionAndGhostCellsInbox<
891 : VolumeDim, UseNodegroupDgElements>>;
892 0 : using const_global_cache_tags =
893 : tmpl::list<evolution::Tags::BoundaryCorrection, ::dg::Tags::Formulation>;
894 :
895 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
896 : typename ArrayIndex, typename ActionList,
897 : typename ParallelComponent>
898 0 : static Parallel::iterable_action_return_t apply(
899 : db::DataBox<DbTagsList>& box, tuples::TaggedTuple<InboxTags...>& inboxes,
900 : const Parallel::GlobalCache<Metavariables>& /*cache*/,
901 : const ArrayIndex& /*array_index*/, ActionList /*meta*/,
902 : const ParallelComponent* const /*meta*/) {
903 : static_assert(
904 : UseNodegroupDgElements ==
905 : Parallel::is_dg_element_collection_v<ParallelComponent>,
906 : "The action ApplyBoundaryCorrectionsToTimeDerivative is told by the "
907 : "template parameter UseNodegroupDgElements that it is being "
908 : "used with a DgElementCollection, but the ParallelComponent "
909 : "is not a DgElementCollection. You need to change the template "
910 : "parameter on the ApplyBoundaryCorrectionsToTimeDerivative action "
911 : "in your action list.");
912 : constexpr size_t volume_dim = Metavariables::system::volume_dim;
913 : const Element<volume_dim>& element =
914 : db::get<domain::Tags::Element<volume_dim>>(box);
915 :
916 : if (UNLIKELY(element.number_of_neighbors() == 0)) {
917 : // We have no neighbors, yay!
918 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
919 : }
920 :
921 : if (not receive_boundary_data_global_time_stepping<
922 : Parallel::is_dg_element_collection_v<ParallelComponent>,
923 : Metavariables>(make_not_null(&box), make_not_null(&inboxes))) {
924 : return {Parallel::AlgorithmExecution::Retry, std::nullopt};
925 : }
926 :
927 : db::mutate_apply<
928 : ApplyBoundaryCorrections<false, Metavariables, VolumeDim, DenseOutput>>(
929 : make_not_null(&box));
930 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
931 : }
932 : };
933 :
934 : /*!
935 : * \brief Computes the boundary corrections for local time-stepping
936 : * and adds them to the variables.
937 : *
938 : * When using local time stepping the neighbor sends data at the neighbor's
939 : * current temporal id. Along with the boundary data, the next temporal id at
940 : * which the neighbor will send data is also sent. This is equal to the
941 : * neighbor's `::Tags::Next<::Tags::TimeStepId>`. When inserting into the mortar
942 : * data history, we insert the received temporal id, that is, the current time
943 : * of the neighbor, along with the boundary correction data.
944 : */
945 : template <size_t VolumeDim, bool DenseOutput, bool UseNodegroupDgElements>
946 1 : struct ApplyLtsBoundaryCorrections {
947 0 : using inbox_tags =
948 : tmpl::list<evolution::dg::Tags::BoundaryCorrectionAndGhostCellsInbox<
949 : VolumeDim, UseNodegroupDgElements>>;
950 0 : using const_global_cache_tags =
951 : tmpl::list<evolution::Tags::BoundaryCorrection, ::dg::Tags::Formulation>;
952 :
953 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
954 : typename ArrayIndex, typename ActionList,
955 : typename ParallelComponent>
956 0 : static Parallel::iterable_action_return_t apply(
957 : db::DataBox<DbTagsList>& box, tuples::TaggedTuple<InboxTags...>& inboxes,
958 : const Parallel::GlobalCache<Metavariables>& /*cache*/,
959 : const ArrayIndex& /*array_index*/, ActionList /*meta*/,
960 : const ParallelComponent* const /*meta*/) {
961 : static_assert(
962 : UseNodegroupDgElements ==
963 : Parallel::is_dg_element_collection_v<ParallelComponent>,
964 : "The action ApplyLtsBoundaryCorrections is told by the "
965 : "template parameter UseNodegroupDgElements that it is being "
966 : "used with a DgElementCollection, but the ParallelComponent "
967 : "is not a DgElementCollection. You need to change the "
968 : "template parameter on the ApplyLtsBoundaryCorrections action "
969 : "in your action list.");
970 : constexpr size_t volume_dim = Metavariables::system::volume_dim;
971 : const Element<volume_dim>& element =
972 : db::get<domain::Tags::Element<volume_dim>>(box);
973 :
974 : if (UNLIKELY(element.number_of_neighbors() == 0)) {
975 : // We have no neighbors, yay!
976 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
977 : }
978 :
979 : if (not receive_boundary_data_local_time_stepping<
980 : Parallel::is_dg_element_collection_v<ParallelComponent>,
981 : typename Metavariables::system, VolumeDim, false>(
982 : make_not_null(&box), make_not_null(&inboxes))) {
983 : return {Parallel::AlgorithmExecution::Retry, std::nullopt};
984 : }
985 :
986 : if (::SelfStart::step_unused(
987 : db::get<::Tags::TimeStepId>(box),
988 : db::get<::Tags::Next<::Tags::TimeStepId>>(box))) {
989 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
990 : }
991 :
992 : db::mutate_apply<
993 : ApplyBoundaryCorrections<true, Metavariables, VolumeDim, DenseOutput>>(
994 : make_not_null(&box));
995 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
996 : }
997 : };
998 : } // namespace Actions
999 : } // namespace evolution::dg
|
