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 <optional>
9 : #include <tuple>
10 : #include <unordered_map>
11 : #include <utility>
12 : #include <vector>
13 :
14 : #include "DataStructures/ApplyMatrices.hpp"
15 : #include "DataStructures/DataBox/DataBox.hpp"
16 : #include "DataStructures/DataBox/DataBoxTag.hpp"
17 : #include "DataStructures/DataBox/PrefixHelpers.hpp"
18 : #include "DataStructures/DataBox/Prefixes.hpp"
19 : #include "DataStructures/Variables.hpp"
20 : #include "Domain/Creators/Tags/Domain.hpp"
21 : #include "Domain/Structure/ChildSize.hpp"
22 : #include "Domain/Structure/Direction.hpp"
23 : #include "Domain/Structure/Element.hpp"
24 : #include "Domain/Structure/Neighbors.hpp"
25 : #include "Domain/Structure/OrientationMap.hpp"
26 : #include "Domain/Structure/TrimMap.hpp"
27 : #include "Domain/Tags.hpp"
28 : #include "Domain/Tags/NeighborMesh.hpp"
29 : #include "Evolution/DiscontinuousGalerkin/InboxTags.hpp"
30 : #include "Evolution/DiscontinuousGalerkin/Initialization/QuadratureTag.hpp"
31 : #include "Evolution/DiscontinuousGalerkin/MortarData.hpp"
32 : #include "Evolution/DiscontinuousGalerkin/MortarDataHolder.hpp"
33 : #include "Evolution/DiscontinuousGalerkin/MortarInfo.hpp"
34 : #include "Evolution/DiscontinuousGalerkin/MortarTags.hpp"
35 : #include "Evolution/DiscontinuousGalerkin/NormalVectorTags.hpp"
36 : #include "NumericalAlgorithms/DiscontinuousGalerkin/MortarHelpers.hpp"
37 : #include "NumericalAlgorithms/Spectral/Mesh.hpp"
38 : #include "NumericalAlgorithms/Spectral/SegmentSize.hpp"
39 : #include "Parallel/AlgorithmExecution.hpp"
40 : #include "ParallelAlgorithms/Amr/Protocols/Projector.hpp"
41 : #include "ParallelAlgorithms/Initialization/MutateAssign.hpp"
42 : #include "Time/BoundaryHistory.hpp"
43 : #include "Time/TimeStepId.hpp"
44 : #include "Utilities/ErrorHandling/Assert.hpp"
45 : #include "Utilities/Gsl.hpp"
46 : #include "Utilities/MakeArray.hpp"
47 : #include "Utilities/TMPL.hpp"
48 :
49 : /// \cond
50 : template <size_t Dim>
51 : class Domain;
52 : namespace Parallel {
53 : template <typename Metavariables>
54 : class GlobalCache;
55 : } // namespace Parallel
56 : namespace Spectral {
57 : enum class Quadrature : uint8_t;
58 : } // namespace Spectral
59 : namespace Tags {
60 : struct TimeStepId;
61 : } // namespace Tags
62 : namespace tuples {
63 : template <class... Tags>
64 : class TaggedTuple;
65 : } // namespace tuples
66 : /// \endcond
67 :
68 1 : namespace evolution::dg::Initialization {
69 : namespace detail {
70 : template <size_t Dim>
71 : ::dg::MortarMap<Dim, evolution::dg::MortarDataHolder<Dim>> empty_mortar_data(
72 : const Element<Dim>& element);
73 :
74 : template <size_t Dim>
75 : ::dg::MortarMap<Dim, MortarInfo<Dim>> mortar_infos(
76 : const Domain<Dim>& domain, const Element<Dim>& element,
77 : const Mesh<Dim>& volume_mesh,
78 : const ::dg::MortarMap<Dim, Mesh<Dim>>& neighbor_mesh);
79 :
80 : template <size_t Dim>
81 : std::tuple<::dg::MortarMap<Dim, Mesh<Dim - 1>>,
82 : ::dg::MortarMap<Dim, TimeStepId>,
83 : DirectionMap<Dim, std::optional<Variables<tmpl::list<
84 : evolution::dg::Tags::MagnitudeOfNormal,
85 : evolution::dg::Tags::NormalCovector<Dim>>>>>>
86 : mortars_apply_impl(const Element<Dim>& element,
87 : const TimeStepId& next_temporal_id,
88 : const Mesh<Dim>& volume_mesh,
89 : const ::dg::MortarMap<Dim, Mesh<Dim>>& neighbor_mesh);
90 :
91 : template <size_t Dim>
92 : void h_refine_structure(
93 : gsl::not_null<::dg::MortarMap<Dim, evolution::dg::MortarDataHolder<Dim>>*>
94 : mortar_data,
95 : gsl::not_null<::dg::MortarMap<Dim, Mesh<Dim - 1>>*> mortar_mesh,
96 : gsl::not_null<::dg::MortarMap<Dim, MortarInfo<Dim>>*> mortar_infos,
97 : gsl::not_null<::dg::MortarMap<Dim, TimeStepId>*> mortar_next_temporal_id,
98 : gsl::not_null<
99 : DirectionMap<Dim, std::optional<::Variables<tmpl::list<
100 : ::evolution::dg::Tags::MagnitudeOfNormal,
101 : ::evolution::dg::Tags::NormalCovector<Dim>>>>>*>
102 : normal_covector_and_magnitude,
103 : const Domain<Dim>& domain, const Mesh<Dim>& new_mesh,
104 : const Element<Dim>& new_element,
105 : const ::dg::MortarMap<Dim, Mesh<Dim>>& neighbor_mesh,
106 : const TimeStepId& current_temporal_id);
107 : } // namespace detail
108 :
109 : /*!
110 : * \brief Initialize mortars between elements for exchanging boundary correction
111 : * terms.
112 : *
113 : * Uses:
114 : * - DataBox:
115 : * - `Tags::Element<Dim>`
116 : * - `Tags::Mesh<Dim>`
117 : * - `BoundaryScheme::receive_temporal_id`
118 : *
119 : * DataBox changes:
120 : * - Adds:
121 : * - `Tags::MortarData<Dim>`
122 : * - `Tags::MortarMesh<Dim>`
123 : * - `Tags::MortarInfo<Dim>`
124 : * - `Tags::MortarNextTemporalId<Dim>`
125 : * - `evolution::dg::Tags::NormalCovectorAndMagnitude<Dim>`
126 : * - Removes: nothing
127 : * - Modifies: nothing
128 : */
129 : template <size_t Dim, typename System>
130 1 : struct Mortars {
131 : public:
132 0 : using const_global_cache_tags = tmpl::list<domain::Tags::Domain<Dim>>;
133 0 : using simple_tags_from_options = tmpl::list<>;
134 :
135 0 : using simple_tags = tmpl::list<
136 : Tags::MortarData<Dim>, Tags::MortarMesh<Dim>, Tags::MortarInfo<Dim>,
137 : Tags::MortarNextTemporalId<Dim>,
138 : evolution::dg::Tags::NormalCovectorAndMagnitude<Dim>,
139 : Tags::MortarDataHistory<
140 : Dim, typename db::add_tag_prefix<
141 : ::Tags::dt, typename System::variables_tag>::type>>;
142 0 : using compute_tags = tmpl::list<>;
143 :
144 : template <typename DbTagsList, typename... InboxTags, typename Metavariables,
145 : typename ArrayIndex, typename ActionList,
146 : typename ParallelComponent>
147 0 : static Parallel::iterable_action_return_t apply(
148 : db::DataBox<DbTagsList>& box,
149 : const tuples::TaggedTuple<InboxTags...>& /*inboxes*/,
150 : const Parallel::GlobalCache<Metavariables>& /*cache*/,
151 : const ArrayIndex& /*array_index*/, ActionList /*meta*/,
152 : const ParallelComponent* const /*meta*/) {
153 : const auto& domain = db::get<domain::Tags::Domain<Dim>>(box);
154 : const auto& element = db::get<::domain::Tags::Element<Dim>>(box);
155 : const auto& volume_mesh = db::get<domain::Tags::Mesh<Dim>>(box);
156 : const auto& neighbor_mesh = db::get<domain::Tags::NeighborMesh<Dim>>(box);
157 : auto mortar_data = detail::empty_mortar_data(element);
158 : auto mortar_infos =
159 : detail::mortar_infos(domain, element, volume_mesh, neighbor_mesh);
160 : auto [mortar_meshes, mortar_next_temporal_ids, normal_covector_quantities] =
161 : detail::mortars_apply_impl(
162 : element, db::get<::Tags::Next<::Tags::TimeStepId>>(box),
163 : db::get<::domain::Tags::Mesh<Dim>>(box),
164 : db::get<::domain::Tags::NeighborMesh<Dim>>(box));
165 : typename Tags::MortarDataHistory<
166 : Dim, typename db::add_tag_prefix<
167 : ::Tags::dt, typename System::variables_tag>::type>::type
168 : boundary_data_history{};
169 : if (Metavariables::local_time_stepping) {
170 : for (const auto& mortar_id_and_data : mortar_data) {
171 : // default initialize data
172 : boundary_data_history[mortar_id_and_data.first];
173 : }
174 : }
175 : ::Initialization::mutate_assign<simple_tags>(
176 : make_not_null(&box), std::move(mortar_data), std::move(mortar_meshes),
177 : std::move(mortar_infos), std::move(mortar_next_temporal_ids),
178 : std::move(normal_covector_quantities),
179 : std::move(boundary_data_history));
180 : return {Parallel::AlgorithmExecution::Continue, std::nullopt};
181 : }
182 : };
183 :
184 : /// \brief Initialize/update items related to mortars after an AMR change
185 : ///
186 : /// Mutates:
187 : /// - Tags::MortarData<dim>
188 : /// - Tags::MortarMesh<dim>
189 : /// - Tags::MortarInfo<dim>
190 : /// - Tags::MortarNextTemporalId<dim>
191 : /// - evolution::dg::Tags::NormalCovectorAndMagnitude<dim>
192 : /// - Tags::MortarDataHistory<dim, typename dt_variables_tag::type>>
193 : ///
194 : /// For p-refined interfaces:
195 : /// - Regenerates MortarData and MortarInfo (should have no effect)
196 : /// - Sets the NormalCovectorAndMagnitude to std::nullopt if the face mesh
197 : /// changed
198 : /// - Projects the local geometric data (but not the data on the mortar-mesh)
199 : /// in the MortarDataHistory, if present
200 : /// - Does nothing to MortarMesh and MortarNextTemporalId
201 : ///
202 : /// For h-refined interfaces:
203 : /// - Regenerates MortarData and MortarInfo
204 : /// - Sets the NormalCovectorAndMagnitude to std::nullopt
205 : /// - Calculates MortarMesh
206 : /// - Sets MortarNextTemporalId to the current temporal id
207 : /// - For local time-stepping:
208 : /// - Removes MortarDataHistory data corresponding to split or joined
209 : /// elements
210 : /// - Projects MortarDataHistory data corresponding to non-h-refined
211 : /// elements onto refined mortars (both geometric and mortar-mesh data)
212 : /// - Creates empty histories for new mortars between two h-refined
213 : /// elements
214 : template <typename Metavariables>
215 1 : struct ProjectMortars : tt::ConformsTo<amr::protocols::Projector> {
216 : private:
217 0 : using magnitude_and_normal_type = ::Variables<tmpl::list<
218 : ::evolution::dg::Tags::MagnitudeOfNormal,
219 : ::evolution::dg::Tags::NormalCovector<Metavariables::volume_dim>>>;
220 :
221 : public:
222 0 : static constexpr size_t dim = Metavariables::volume_dim;
223 0 : using dt_variables_tag = typename db::add_tag_prefix<
224 : ::Tags::dt, typename Metavariables::system::variables_tag>;
225 0 : using mortar_data_history_tag =
226 : Tags::MortarDataHistory<dim, typename dt_variables_tag::type>;
227 0 : using mortar_data_history_type = typename mortar_data_history_tag::type;
228 :
229 0 : using return_tags =
230 : tmpl::list<Tags::MortarData<dim>, Tags::MortarMesh<dim>,
231 : Tags::MortarInfo<dim>, Tags::MortarNextTemporalId<dim>,
232 : evolution::dg::Tags::NormalCovectorAndMagnitude<dim>,
233 : Tags::MortarDataHistory<dim, typename dt_variables_tag::type>>;
234 0 : using argument_tags =
235 : tmpl::list<domain::Tags::Domain<dim>, domain::Tags::Mesh<dim>,
236 : domain::Tags::Element<dim>, domain::Tags::NeighborMesh<dim>,
237 : ::Tags::TimeStepId>;
238 :
239 0 : static void apply(
240 : const gsl::not_null<
241 : ::dg::MortarMap<dim, evolution::dg::MortarDataHolder<dim>>*>
242 : mortar_data,
243 : const gsl::not_null<::dg::MortarMap<dim, Mesh<dim - 1>>*> mortar_mesh,
244 : const gsl::not_null<::dg::MortarMap<dim, MortarInfo<dim>>*> mortar_infos,
245 : const gsl::not_null<::dg::MortarMap<dim, TimeStepId>*>
246 : mortar_next_temporal_id,
247 : const gsl::not_null<
248 : DirectionMap<dim, std::optional<magnitude_and_normal_type>>*>
249 : normal_covector_and_magnitude,
250 : const gsl::not_null<mortar_data_history_type*> mortar_data_history,
251 : const Domain<dim>& domain, const Mesh<dim>& new_mesh,
252 : const Element<dim>& new_element,
253 : const ::dg::MortarMap<dim, Mesh<dim>>& neighbor_mesh,
254 : const TimeStepId& current_temporal_id,
255 : const std::pair<Mesh<dim>, Element<dim>>& old_mesh_and_element) {
256 : const auto& [old_mesh, old_element] = old_mesh_and_element;
257 : ASSERT(old_element.id() == new_element.id(),
258 : "p-refinement should not have changed the element id");
259 :
260 : const bool mesh_changed = old_mesh != new_mesh;
261 :
262 : auto new_mortar_infos =
263 : detail::mortar_infos(domain, new_element, new_mesh, neighbor_mesh);
264 :
265 : // The old mortars must be removed from the MortarMaps before the
266 : // new ones can be added to avoid potentially exceeding the
267 : // MortarMap capacity, so we collect the new data and add it at
268 : // the end.
269 : using NewMortarEntry = std::tuple<
270 : DirectionalId<dim>, Mesh<dim - 1>,
271 : std::optional<typename mortar_data_history_type::mapped_type>>;
272 : std::vector<NewMortarEntry> new_mortars{};
273 :
274 : for (const auto& [direction, neighbors] : new_element.neighbors()) {
275 : const auto sliced_away_dimension = direction.dimension();
276 : const auto old_face_mesh = old_mesh.slice_away(sliced_away_dimension);
277 : const auto new_face_mesh = new_mesh.slice_away(sliced_away_dimension);
278 : const bool face_mesh_changed = old_face_mesh != new_face_mesh;
279 : if (face_mesh_changed) {
280 : (*normal_covector_and_magnitude)[direction] = std::nullopt;
281 : }
282 : for (const auto& neighbor : neighbors) {
283 : const DirectionalId<dim> mortar_id{direction, neighbor};
284 : const auto& new_neighbor_mesh = neighbor_mesh.at(mortar_id);
285 : const auto new_mortar_mesh = ::dg::mortar_mesh(
286 : new_face_mesh, new_neighbor_mesh.slice_away(sliced_away_dimension));
287 : if (mortar_mesh->contains(mortar_id)) {
288 : // Set the mortar mesh, but do not project any existing mesh
289 : // data. The mesh needs to have a valid value in order to
290 : // project data when we send to our neighbors. If the mesh
291 : // resolution needs to be changed afterwards because of a
292 : // change in the neighbor, that's fine, because
293 : // element-to-mortar projections are lossless. Projecting
294 : // existing data would be bad, because we might erroneously
295 : // decrease the mesh resolution, losing the high-order
296 : // modes.
297 : mortar_mesh->at(mortar_id) = new_mortar_mesh;
298 : // mortar_data does not need projecting as it has already been used
299 : // and will be resized automatically
300 : if (mesh_changed and not mortar_data_history->empty()) {
301 : auto& boundary_history = mortar_data_history->at(mortar_id);
302 : auto local_history = boundary_history.local();
303 : const auto project_local_boundary_data =
304 : [&new_face_mesh, &new_mesh](
305 : const TimeStepId& /* id */,
306 : const gsl::not_null<::evolution::dg::MortarData<dim>*>
307 : data) {
308 : return p_project_geometric_data(data, new_face_mesh,
309 : new_mesh);
310 : };
311 : local_history.for_each(project_local_boundary_data);
312 : }
313 : } else {
314 : const auto& new_mortar_size =
315 : new_mortar_infos.at(mortar_id).mortar_size();
316 :
317 : std::optional<typename mortar_data_history_type::mapped_type>
318 : new_history{};
319 : for (const auto& [old_mortar_id, old_history] :
320 : *mortar_data_history) {
321 : if (old_mortar_id.direction() != direction or
322 : not overlapping(old_mortar_id.id(), neighbor)) {
323 : continue;
324 : }
325 : const auto& old_mortar_size =
326 : mortar_infos->at(old_mortar_id).mortar_size();
327 : if (not new_history.has_value()) {
328 : new_history.emplace(old_history);
329 : new_history->remote().clear();
330 : auto local_history = new_history->local();
331 : if (mesh_changed) {
332 : const auto project_face_data =
333 : [&new_face_mesh, &new_mesh](
334 : const TimeStepId& /* id */,
335 : const gsl::not_null<::evolution::dg::MortarData<dim>*>
336 : data) {
337 : return p_project_geometric_data(data, new_face_mesh,
338 : new_mesh);
339 : };
340 : local_history.for_each(project_face_data);
341 : }
342 : const auto project_mortar_data =
343 : [&new_mortar_mesh, &new_mortar_size, &old_mortar_size](
344 : const TimeStepId& /* id */,
345 : const gsl::not_null<::evolution::dg::MortarData<dim>*>
346 : data) {
347 : const auto& old_mortar_mesh = data->mortar_mesh.value();
348 : const auto mortar_projection_matrices =
349 : Spectral::projection_matrices(
350 : old_mortar_mesh, new_mortar_mesh, old_mortar_size,
351 : new_mortar_size);
352 : DataVector& vars = data->mortar_data.value();
353 : vars = apply_matrices(mortar_projection_matrices, vars,
354 : old_mortar_mesh.extents());
355 : data->mortar_mesh = new_mortar_mesh;
356 : return true;
357 : };
358 : local_history.for_each(project_mortar_data);
359 : } else {
360 : auto local_history = new_history->local();
361 : const auto old_local_history = old_history.local();
362 : const auto project_local_mortar_data =
363 : [&new_mortar_mesh, &new_mortar_size, &old_local_history,
364 : &old_mortar_size](
365 : const TimeStepId& id,
366 : const gsl::not_null<::evolution::dg::MortarData<dim>*>
367 : data) {
368 : const auto& old_data = old_local_history.data(id);
369 : const auto& old_mortar_mesh = old_data.mortar_mesh.value();
370 : const auto mortar_projection_matrices =
371 : Spectral::projection_matrices(
372 : old_mortar_mesh, new_mortar_mesh, old_mortar_size,
373 : new_mortar_size);
374 : data->mortar_data.value() +=
375 : apply_matrices(mortar_projection_matrices,
376 : old_data.mortar_data.value(),
377 : old_mortar_mesh.extents());
378 : return true;
379 : };
380 : local_history.for_each(project_local_mortar_data);
381 : }
382 : }
383 :
384 : new_mortars.emplace_back(mortar_id, new_mortar_mesh, new_history);
385 : }
386 : }
387 : }
388 :
389 : if (not new_mortars.empty()) {
390 : domain::remove_nonexistent_neighbors(mortar_mesh, new_element);
391 : domain::remove_nonexistent_neighbors(mortar_next_temporal_id,
392 : new_element);
393 : domain::remove_nonexistent_neighbors(mortar_data_history, new_element);
394 : for (auto& [mortar_id, new_mortar_mesh, new_history] : new_mortars) {
395 : mortar_mesh->emplace(mortar_id, std::move(new_mortar_mesh));
396 : // We only do h refinement at a slab boundary, so we know all
397 : // the neighbors are aligned with us temporally.
398 : mortar_next_temporal_id->emplace(mortar_id, current_temporal_id);
399 : if (new_history.has_value()) {
400 : mortar_data_history->emplace(mortar_id, *new_history);
401 : }
402 : }
403 : }
404 :
405 : *mortar_data = detail::empty_mortar_data(new_element);
406 : *mortar_infos = std::move(new_mortar_infos);
407 :
408 : for (const auto& direction : new_element.external_boundaries()) {
409 : const auto sliced_away_dimension = direction.dimension();
410 : const auto old_face_mesh = old_mesh.slice_away(sliced_away_dimension);
411 : const auto new_face_mesh = new_mesh.slice_away(sliced_away_dimension);
412 : const bool face_mesh_changed = old_face_mesh != new_face_mesh;
413 : if (face_mesh_changed) {
414 : (*normal_covector_and_magnitude)[direction] = std::nullopt;
415 : }
416 : }
417 : }
418 :
419 : template <typename... ParentTags>
420 0 : static void apply(
421 : const gsl::not_null<
422 : ::dg::MortarMap<dim, evolution::dg::MortarDataHolder<dim>>*>
423 : mortar_data,
424 : const gsl::not_null<::dg::MortarMap<dim, Mesh<dim - 1>>*> mortar_mesh,
425 : const gsl::not_null<::dg::MortarMap<dim, MortarInfo<dim>>*> mortar_infos,
426 : const gsl::not_null<::dg::MortarMap<dim, TimeStepId>*>
427 : mortar_next_temporal_id,
428 : const gsl::not_null<
429 : DirectionMap<dim, std::optional<magnitude_and_normal_type>>*>
430 : normal_covector_and_magnitude,
431 : const gsl::not_null<mortar_data_history_type*> mortar_data_history,
432 : const Domain<dim>& domain, const Mesh<dim>& new_mesh,
433 : const Element<dim>& new_element,
434 : const ::dg::MortarMap<dim, Mesh<dim>>& neighbor_mesh,
435 : const TimeStepId& /*possibly_unset*/,
436 : const tuples::TaggedTuple<ParentTags...>& parent_items) {
437 : detail::h_refine_structure(
438 : mortar_data, mortar_mesh, mortar_infos, mortar_next_temporal_id,
439 : normal_covector_and_magnitude, domain, new_mesh, new_element,
440 : neighbor_mesh, get<::Tags::TimeStepId>(parent_items));
441 :
442 : if (Metavariables::local_time_stepping) {
443 : const auto& old_element = get<domain::Tags::Element<dim>>(parent_items);
444 : const auto& old_histories = get<mortar_data_history_tag>(parent_items);
445 : for (const auto& [direction, neighbors] : new_element.neighbors()) {
446 : for (const auto& neighbor : neighbors) {
447 : const DirectionalId<dim> mortar_id{direction, neighbor};
448 : if (const auto old_history = old_histories.find(mortar_id);
449 : old_history != old_histories.end()) {
450 : // The neighbor did not h-refine, so we have to project
451 : // its mortar data from our parent.
452 : auto& new_history =
453 : mortar_data_history->emplace(mortar_id, old_history->second)
454 : .first->second;
455 : new_history.local().clear();
456 : auto remote_history = new_history.remote();
457 : const auto& new_mortar_mesh = mortar_mesh->at(mortar_id);
458 : const auto& orientation = neighbors.orientation(neighbor);
459 : const auto new_mortar_size = domain::child_size(
460 : ::dg::mortar_segments(new_element.id(), neighbor,
461 : direction.dimension(), orientation),
462 : ::dg::mortar_segments(old_element.id(), neighbor,
463 : direction.dimension(), orientation));
464 : const auto project_mortar_data =
465 : [&new_mortar_mesh, &new_mortar_size](
466 : const TimeStepId& /* id */,
467 : const gsl::not_null<::evolution::dg::MortarData<dim>*>
468 : data) {
469 : const auto& old_mortar_mesh = data->mortar_mesh.value();
470 : const auto mortar_projection_matrices =
471 : Spectral::projection_matrices(
472 : old_mortar_mesh, new_mortar_mesh,
473 : make_array<dim - 1>(Spectral::SegmentSize::Full),
474 : new_mortar_size);
475 : DataVector& vars = data->mortar_data.value();
476 : vars = apply_matrices(mortar_projection_matrices, vars,
477 : old_mortar_mesh.extents());
478 : data->mortar_mesh = new_mortar_mesh;
479 : return true;
480 : };
481 : remote_history.for_each(project_mortar_data);
482 : } else {
483 : // Neither this element nor the neighbor existed before
484 : // refinement.
485 : mortar_data_history->emplace(
486 : mortar_id, typename mortar_data_history_type::mapped_type{});
487 : }
488 : }
489 : }
490 : }
491 : }
492 :
493 : template <typename... ChildTags>
494 0 : static void apply(
495 : const gsl::not_null<
496 : ::dg::MortarMap<dim, evolution::dg::MortarDataHolder<dim>>*>
497 : mortar_data,
498 : const gsl::not_null<::dg::MortarMap<dim, Mesh<dim - 1>>*> mortar_mesh,
499 : const gsl::not_null<::dg::MortarMap<dim, MortarInfo<dim>>*> mortar_infos,
500 : const gsl::not_null<::dg::MortarMap<dim, TimeStepId>*>
501 : mortar_next_temporal_id,
502 : const gsl::not_null<
503 : DirectionMap<dim, std::optional<magnitude_and_normal_type>>*>
504 : normal_covector_and_magnitude,
505 : const gsl::not_null<mortar_data_history_type*> mortar_data_history,
506 : const Domain<dim>& domain, const Mesh<dim>& new_mesh,
507 : const Element<dim>& new_element,
508 : const ::dg::MortarMap<dim, Mesh<dim>>& neighbor_mesh,
509 : const TimeStepId& /*possibly_unset*/,
510 : const std::unordered_map<
511 : ElementId<dim>, tuples::TaggedTuple<ChildTags...>>& children_items) {
512 : detail::h_refine_structure(
513 : mortar_data, mortar_mesh, mortar_infos, mortar_next_temporal_id,
514 : normal_covector_and_magnitude, domain, new_mesh, new_element,
515 : neighbor_mesh, get<::Tags::TimeStepId>(children_items.begin()->second));
516 :
517 : if (Metavariables::local_time_stepping) {
518 : for (const auto& [direction, neighbors] : new_element.neighbors()) {
519 : for (const auto& neighbor : neighbors) {
520 : const DirectionalId<dim> mortar_id{direction, neighbor};
521 : std::optional<typename mortar_data_history_type::mapped_type>
522 : new_history{};
523 : for (const auto& [child, child_items] : children_items) {
524 : const auto& old_histories =
525 : get<mortar_data_history_tag>(child_items);
526 : if (const auto old_history = old_histories.find(mortar_id);
527 : old_history != old_histories.end()) {
528 : // The neighbor did not h-refine, so we have to project
529 : // its mortar data from our children.
530 : const auto& new_mortar_mesh = mortar_mesh->at(mortar_id);
531 : const auto& orientation = neighbors.orientation(neighbor);
532 : const auto old_mortar_size = domain::child_size(
533 : ::dg::mortar_segments(child, neighbor, direction.dimension(),
534 : orientation),
535 : ::dg::mortar_segments(new_element.id(), neighbor,
536 : direction.dimension(), orientation));
537 : if (not new_history.has_value()) {
538 : new_history.emplace(old_history->second);
539 : new_history->local().clear();
540 : auto remote_history = new_history->remote();
541 : const auto project_mortar_data =
542 : [&new_mortar_mesh, &old_mortar_size](
543 : const TimeStepId& /* id */,
544 : const gsl::not_null<::evolution::dg::MortarData<dim>*>
545 : data) {
546 : const auto& old_mortar_mesh = data->mortar_mesh.value();
547 : const auto mortar_projection_matrices =
548 : Spectral::projection_matrices(
549 : old_mortar_mesh, new_mortar_mesh, old_mortar_size,
550 : make_array<dim - 1>(Spectral::SegmentSize::Full));
551 : DataVector& vars = data->mortar_data.value();
552 : vars = apply_matrices(mortar_projection_matrices, vars,
553 : old_mortar_mesh.extents());
554 : data->mortar_mesh = new_mortar_mesh;
555 : return true;
556 : };
557 : remote_history.for_each(project_mortar_data);
558 : } else {
559 : auto remote_history = new_history->remote();
560 : const auto old_remote_history = old_history->second.remote();
561 : const auto project_mortar_data =
562 : [&new_mortar_mesh, &old_mortar_size, &old_remote_history](
563 : const TimeStepId& id,
564 : const gsl::not_null<::evolution::dg::MortarData<dim>*>
565 : data) {
566 : const auto& old_data = old_remote_history.data(id);
567 : const auto& old_mortar_mesh =
568 : old_data.mortar_mesh.value();
569 : const auto mortar_projection_matrices =
570 : Spectral::projection_matrices(
571 : old_mortar_mesh, new_mortar_mesh, old_mortar_size,
572 : make_array<dim - 1>(Spectral::SegmentSize::Full));
573 : data->mortar_data.value() +=
574 : apply_matrices(mortar_projection_matrices,
575 : old_data.mortar_data.value(),
576 : old_mortar_mesh.extents());
577 : return true;
578 : };
579 : remote_history.for_each(project_mortar_data);
580 : }
581 : }
582 : }
583 :
584 : if (new_history.has_value()) {
585 : mortar_data_history->emplace(mortar_id, std::move(*new_history));
586 : } else {
587 : // Neither this element nor the neighbor existed before
588 : // refinement.
589 : mortar_data_history->emplace(
590 : mortar_id, typename mortar_data_history_type::mapped_type{});
591 : }
592 : }
593 : }
594 : }
595 : }
596 : };
597 : } // namespace evolution::dg::Initialization
|
