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 <boost/functional/hash.hpp>
8 : #include <cstddef>
9 : #include <iomanip>
10 : #include <map>
11 : #include <memory>
12 : #include <optional>
13 : #include <sstream>
14 : #include <string>
15 : #include <type_traits>
16 : #include <utility>
17 :
18 : #include "DataStructures/DataBox/Tag.hpp"
19 : #include "DataStructures/FixedHashMap.hpp"
20 : #include "Domain/Structure/Direction.hpp"
21 : #include "Domain/Structure/DirectionalId.hpp"
22 : #include "Domain/Structure/DirectionalIdMap.hpp"
23 : #include "Domain/Structure/ElementId.hpp"
24 : #include "Evolution/DiscontinuousGalerkin/AtomicInboxBoundaryData.hpp"
25 : #include "Evolution/DiscontinuousGalerkin/BoundaryData.hpp"
26 : #include "Evolution/DiscontinuousGalerkin/Messages/BoundaryMessage.hpp"
27 : #include "NumericalAlgorithms/Spectral/Mesh.hpp"
28 : #include "Parallel/InboxInserters.hpp"
29 : #include "Parallel/StaticSpscQueue.hpp"
30 : #include "Time/TimeStepId.hpp"
31 : #include "Utilities/TMPL.hpp"
32 :
33 : /// \cond
34 : class DataVector;
35 : /// \endcond
36 :
37 1 : namespace evolution::dg::Tags {
38 : /*!
39 : * \brief The inbox tag for boundary correction communication and DG-subcell
40 : * ghost zone cells.
41 : *
42 : * The stored data consists of the following:
43 : *
44 : * 1. the volume mesh of the element.
45 : * 2. the volume mesh corresponding to the ghost cell data. This allows eliding
46 : * projection when all neighboring elements are doing DG.
47 : * 3. the mortar mesh of the data on the mortar
48 : * 4. the variables at the ghost zone cells for finite difference/volume
49 : * reconstruction
50 : * 5. the data on the mortar needed for computing the boundary corrections (e.g.
51 : * fluxes, characteristic speeds, conserved variables)
52 : * 6. the TimeStepId beyond which the boundary terms are no longer valid, when
53 : * using local time stepping.
54 : * 7. the troublade cell indicator status using for determining halos around
55 : * troubled cells.
56 : * 8. the integration order of the time-stepper.
57 : * 9. the InterpolatedBoundaryData sent by a non-conforming Element that
58 : * interpolates its data to a subset of the points of the Element receiving
59 : * this BoundaryData
60 : *
61 : * The TimeStepId is the neighboring element's next time step. When using local
62 : * time stepping, the neighbor's boundary data is valid up until this time,
63 : * which may include multiple local time steps. By receiving and storing the
64 : * neighbor time step, the local element knows whether or not it should remove
65 : * boundary data and expect new data to be sent from the neighbor.
66 : *
67 : * The ghost cell data will be valid whenever a DG-subcell scheme is being used.
68 : * Whenever a DG-subcell scheme is being used, elements using DG and not FD/FV
69 : * always send both the ghost cells and boundary correction data together.
70 : * Elements using FD/FV send the ghost cells first followed by the boundary
71 : * correction data once the element has received all neighbor ghost cell data.
72 : * Note that the second send/receive only modifies the flux and the TimeStepId
73 : * used for the flux validity range.
74 : *
75 : * When only a DG scheme (not a DG-subcell scheme) is used the ghost cell data
76 : * will never be valid.
77 : *
78 : * In the DG-subcell scheme this tag is used both for communicating the ghost
79 : * cell data needed for the FD/FV reconstruction step and the data needed for
80 : * the boundary corrections.
81 : * - For an element using DG, both ghost cells and boundary corrections are
82 : * sent using a single communication. After receiving all neighbor
83 : * boundary corrections the element can finish computing the time step.
84 : * The ghost cell data from neighbors is unused.
85 : * - For an element using FD/FV, first the ghost cells are sent. Once all
86 : * neighboring ghost cell data is received, reconstruction is done and the
87 : * boundary terms are computed and sent to the neighbors. After receiving all
88 : * neighbor boundary corrections the element can finish computing the time
89 : * step.
90 : * - Whether or not an extra communication is needed when an element switches
91 : * from DG to FD/FV depends on how exactly the decision to switch is
92 : * implemented. If the volume terms are integrated and verified to be
93 : * valid before a DG element sends ghost cell and boundary data then no
94 : * additional communication is needed when switching from DG to FD/FV. In this
95 : * case a second check of the data that includes the boundary correction needs
96 : * to be done. If the second check determines a switch from DG to FD/FV is
97 : * needed, we can continue to use the DG fluxes since the evolution in the
98 : * small was valid, thus avoiding an additional communication. However, to
99 : * fully guarantee physical realizability a second communication or evolving
100 : * the neighboring ghost cells needs to be done. We have not yet decided how
101 : * to deal with the possible need for an additional communication since it
102 : * also depends on whether or not we decide to couple to Voronoi instead of
103 : * just Neumann neighbors.
104 : * - The data for the inbox tags is erased after the boundary correction is
105 : * complete and the solution has been verified to be valid at the new time
106 : * step. The ghost cells could be invalidated immediately after
107 : * reconstruction, thus using the ghost cell data after reconstruction is
108 : * complete is considered undefined behavior. That is, we make no guarantee as
109 : * to whether or not it will work.
110 : * - The reason for minimizing the number of communications rather than having a
111 : * more uniform implementation between DG and FD/FV is that it is the number
112 : * of communications that adds the most overhead, not the size of each
113 : * communication. Thus, one large communication is cheaper than several small
114 : * communications.
115 : *
116 : * #### Return Values:
117 : * - In the case that the type is `type_map` the `insert_into_inbox` function
118 : * returns the size of the inbox.
119 : * - In the case that the type is `type_spsc` the `insert_into_inbox` function
120 : * returns the number of neighbor data contributions made that allow the
121 : * element to take its next time step/needs a message called on it. When
122 : * this number is equal to the number of neighbors, a Charm++ message must
123 : * be sent to the element to have it continue the algorithm. This is so as
124 : * to minimize the number of communications made through Charm++ and instead
125 : * to move data atomically between neighbors whenever possible.
126 : *
127 : * #### DG Element Nodegroup Support
128 : * If you are using the `DgElementCollection` then you must set
129 : * `UseNodegroupDgElements` to `true`. The actions that use this tag check
130 : * that the parallel component and the `UseNodegroupDgElements` is consistent.
131 : */
132 : template <size_t Dim, bool UseNodegroupDgElements>
133 1 : struct BoundaryCorrectionAndGhostCellsInbox {
134 0 : using stored_type = evolution::dg::BoundaryData<Dim>;
135 :
136 : public:
137 0 : using temporal_id = TimeStepId;
138 : // Used by array implementation
139 0 : using type_map = std::map<TimeStepId, DirectionalIdMap<Dim, stored_type>>;
140 :
141 : // Used by nodegroup implementation
142 0 : using type_spsc = evolution::dg::AtomicInboxBoundaryData<Dim>;
143 :
144 : // The actual type being used.
145 0 : using type = tmpl::conditional_t<UseNodegroupDgElements, type_spsc, type_map>;
146 0 : using value_type = type;
147 :
148 : template <typename ReceiveDataType>
149 0 : static size_t insert_into_inbox(const gsl::not_null<type_spsc*> inbox,
150 : const temporal_id& time_step_id,
151 : ReceiveDataType&& data) {
152 : const DirectionalId<Dim>& neighbor_id = data.first;
153 : // Note: This assumes the neighbor_id is oriented into our (the element
154 : // whose inbox this is) frame.
155 : const size_t neighbor_index = inbox->index(neighbor_id);
156 : if (UNLIKELY(not gsl::at(inbox->boundary_data_in_directions, neighbor_index)
157 : .try_emplace(time_step_id, std::move(data.second),
158 : std::move(data.first)))) {
159 : ERROR(
160 : "Failed to emplace data into inbox. neighbor_id: ("
161 : << neighbor_id.direction() << ',' << neighbor_id.id()
162 : << ") at TimeStepID: " << time_step_id << " the size of the inbox is "
163 : << gsl::at(inbox->boundary_data_in_directions, neighbor_index).size()
164 : << " the message count is " << inbox->message_count.load()
165 : << " and the number of neighbors is "
166 : << inbox->number_of_neighbors.load());
167 : }
168 : // Notes:
169 : // 1. fetch_add does a post-increment.
170 : // 2. We need thread synchronization here, so doing relaxed_order would be a
171 : // bug.
172 : return inbox->message_count.fetch_add(1, std::memory_order_acq_rel) + 1;
173 : }
174 :
175 : template <typename ReceiveDataType>
176 0 : static size_t insert_into_inbox(const gsl::not_null<type_map*> inbox,
177 : const temporal_id& time_step_id,
178 : ReceiveDataType&& data) {
179 : auto& current_inbox = (*inbox)[time_step_id];
180 : if (auto it = current_inbox.find(data.first); it != current_inbox.end()) {
181 : auto& [volume_mesh, volume_mesh_ghost_cell_data, boundary_correction_mesh,
182 : ghost_cell_data, boundary_correction_data, validity_range,
183 : tci_status, integration_order, interpolated_boundary_data] =
184 : data.second;
185 : (void)ghost_cell_data;
186 : auto& [current_volume_mesh, current_volume_mesh_ghost_cell_data,
187 : current_boundary_correction_mesh, current_ghost_cell_data,
188 : current_boundary_correction_data, current_validity_range,
189 : current_tci_status, current_integration_order,
190 : current_interpolated_boundary_data] = it->second;
191 : (void)current_volume_mesh_ghost_cell_data; // Need to use when
192 : // optimizing subcell
193 : // We have already received some data at this time. Receiving data twice
194 : // at the same time should only occur when receiving fluxes after having
195 : // previously received ghost cells. We sanity check that the data we
196 : // already have is the ghost cells and that we have not yet received flux
197 : // data.
198 : //
199 : // This is used if a 2-send implementation is used (which we don't right
200 : // now!). We generally find that the number of communications is more
201 : // important than the size of each communication, and so a single
202 : // communication per time/sub step is preferred.
203 : ASSERT(current_ghost_cell_data.has_value(),
204 : "Have not yet received ghost cells at time step "
205 : << time_step_id
206 : << " but the inbox entry already exists. This is a bug in the "
207 : "ordering of the actions.");
208 : ASSERT(not current_boundary_correction_data.has_value() and
209 : not current_boundary_correction_mesh.has_value(),
210 : "The fluxes have already been received at time step "
211 : << time_step_id
212 : << ". They are either being received for a second time, there "
213 : "is a bug in the ordering of the actions (though a "
214 : "different ASSERT should've caught that), or the incorrect "
215 : "temporal ID is being sent.");
216 :
217 : ASSERT(current_volume_mesh == volume_mesh,
218 : "The mesh being received for the fluxes is different than the "
219 : "mesh received for the ghost cells. Mesh for fluxes: "
220 : << volume_mesh << " mesh for ghost cells "
221 : << current_volume_mesh);
222 : ASSERT(current_volume_mesh_ghost_cell_data == volume_mesh_ghost_cell_data,
223 : "The mesh being received for the ghost cell data is different "
224 : "than the mesh received previously. Mesh for received when we got "
225 : "fluxes: "
226 : << volume_mesh_ghost_cell_data
227 : << " mesh received when we got ghost cells "
228 : << current_volume_mesh_ghost_cell_data);
229 :
230 : // We always move here since we take ownership of the data and moves
231 : // implicitly decay to copies
232 : current_boundary_correction_data = std::move(boundary_correction_data);
233 : current_validity_range = validity_range;
234 : current_tci_status = tci_status;
235 : current_integration_order = integration_order;
236 : current_interpolated_boundary_data = interpolated_boundary_data;
237 : } else {
238 : // We have not received ghost cells or fluxes at this time.
239 : if (not current_inbox.insert(std::forward<ReceiveDataType>(data))
240 : .second) {
241 : ERROR("Failed to insert data to receive at instance '"
242 : << time_step_id
243 : << "' with tag 'BoundaryCorrectionAndGhostCellsInbox'.\n");
244 : }
245 : }
246 : return current_inbox.size();
247 : }
248 :
249 0 : static std::string output_inbox(const type_spsc& inbox,
250 : const size_t padding_size) {
251 : std::stringstream ss{};
252 : const std::string pad(padding_size, ' ');
253 : ss << std::scientific << std::setprecision(16);
254 : ss << pad << "BoundaryCorrectionAndGhostCellInbox:\n";
255 : ss << pad
256 : << "Warning: Printing atomic state is not possible in general so data "
257 : "printed is limited.\n";
258 : for (size_t i = 0; i < inbox.boundary_data_in_directions.size(); ++i) {
259 : const auto& data_in_direction =
260 : gsl::at(inbox.boundary_data_in_directions, i);
261 : ss << pad << "Id: "
262 : << "Approximate size: " << data_in_direction.size() << "\n";
263 : }
264 :
265 : return ss.str();
266 : }
267 :
268 0 : static std::string output_inbox(const type_map& inbox,
269 : const size_t padding_size) {
270 : std::stringstream ss{};
271 : const std::string pad(padding_size, ' ');
272 : ss << std::scientific << std::setprecision(16);
273 : ss << pad << "BoundaryCorrectionAndGhostCellInbox:\n";
274 :
275 : for (const auto& [current_time_step_id, hash_map] : inbox) {
276 : ss << pad << " Current time: " << current_time_step_id << "\n";
277 : // We only care about the next time because that's important for deadlock
278 : // detection. The data itself isn't super important
279 : for (const auto& [key, boundary_data] : hash_map) {
280 : ss << pad << " Key: " << key
281 : << ", next time: " << boundary_data.validity_range << "\n";
282 : }
283 : }
284 :
285 : return ss.str();
286 : }
287 :
288 0 : void pup(PUP::er& /*p*/) {}
289 : };
290 :
291 : /*!
292 : * \brief Simple tag used to store inbox data in the DataBox.
293 : *
294 : * Since the inbox data can be received asynchronously and lockfree ordered
295 : * containers are slow and challenging to implement, we instead use an
296 : * unordered container for the inbox, then transfer the data into an ordered map
297 : * in the DataBox.
298 : */
299 : template <size_t Dim>
300 1 : struct BoundaryData : db::SimpleTag {
301 0 : using type =
302 : typename BoundaryCorrectionAndGhostCellsInbox<Dim, false>::type_map;
303 : };
304 :
305 : /*!
306 : * \brief The inbox tag for boundary correction communication and DG-subcell
307 : * ghost zone cells using a `BoundaryMessage` object
308 : *
309 : * To see what is stored within a `BoundaryMessage`, see its documentation.
310 : *
311 : * This inbox tag is very similar to `BoundaryCorrectionAndGhostCellsInbox` in
312 : * that it stores subcell/DG data sent from neighboring elements. To see exactly
313 : * when data is stored and how it's used, see the docs for
314 : * `BoundaryCorrectionAndGhostCellsInbox`. This inbox tag is different than
315 : * `BoundaryCorrectionAndGhostCellsInbox` in that it only takes a pointer to a
316 : * `BoundaryMessage` as an argument to `insert_into_inbox` and stores a
317 : * `std::unique_ptr<BoundaryMessage>` inside the inbox.
318 : *
319 : * This inbox tag is meant to be used to avoid unnecessary copies between
320 : * elements on the same node which share a block of memory. If two elements
321 : * aren't on the same node, then a copy/send is done regardless.
322 : *
323 : * \warning The `boundary_message` argument to `insert_into_inbox()` will be
324 : * invalid after the function is called because a `std::unique_ptr` now controls
325 : * the memory. Calling a method on the `boundary_message` pointer after the
326 : * `insert_into_inbox()` function is called can result in undefined behaviour.
327 : */
328 : template <size_t Dim>
329 1 : struct BoundaryMessageInbox {
330 0 : using stored_type = std::unique_ptr<BoundaryMessage<Dim>>;
331 :
332 : public:
333 0 : using temporal_id = TimeStepId;
334 0 : using type = std::map<TimeStepId, DirectionalIdMap<Dim, stored_type>>;
335 0 : using message_type = BoundaryMessage<Dim>;
336 :
337 : template <typename Inbox>
338 0 : static void insert_into_inbox(const gsl::not_null<Inbox*> inbox,
339 : BoundaryMessage<Dim>* boundary_message) {
340 : const auto& time_step_id = boundary_message->current_time_step_id;
341 : auto& current_inbox = (*inbox)[time_step_id];
342 :
343 : const auto key = DirectionalId<Dim>{boundary_message->neighbor_direction,
344 : boundary_message->element_id};
345 :
346 : if (auto it = current_inbox.find(key); it != current_inbox.end()) {
347 : auto& current_boundary_data = it->second;
348 : // We have already received some data at this time. Receiving data twice
349 : // at the same time should only occur when receiving fluxes after having
350 : // previously received ghost cells. We sanity check that the data we
351 : // already have is the ghost cells and that we have not yet received flux
352 : // data.
353 : //
354 : // This is used if a 2-send implementation is used (which we don't right
355 : // now!). We generally find that the number of communications is more
356 : // important than the size of each communication, and so a single
357 : // communication per time/sub step is preferred.
358 : ASSERT(current_boundary_data->subcell_ghost_data != nullptr,
359 : "Have not yet received ghost cells at time step "
360 : << time_step_id
361 : << " but the inbox entry already exists. This is a bug in the "
362 : "ordering of the actions.");
363 : ASSERT(current_boundary_data->dg_flux_data == nullptr,
364 : "The fluxes have already been received at time step "
365 : << time_step_id
366 : << ". They are either being received for a second time, there "
367 : "is a bug in the ordering of the actions (though a "
368 : "different ASSERT should've caught that), or the incorrect "
369 : "temporal ID is being sent.");
370 : ASSERT(boundary_message->subcell_ghost_data == nullptr,
371 : "Have already received ghost cells at time step "
372 : << time_step_id
373 : << ", but we are being sent ghost cells again. This data will "
374 : "not be cleaned up while the dg flux data is being "
375 : "inserted into the inbox which will cause a memory leak.");
376 :
377 : ASSERT(current_boundary_data->interface_mesh ==
378 : boundary_message->interface_mesh,
379 : "The mesh being received for the fluxes is different than the "
380 : "mesh received for the ghost cells. Mesh for fluxes: "
381 : << boundary_message->interface_mesh << " mesh for ghost cells "
382 : << current_boundary_data->interface_mesh);
383 : ASSERT(current_boundary_data->volume_or_ghost_mesh ==
384 : boundary_message->volume_or_ghost_mesh,
385 : "The mesh being received for the ghost cell data is different "
386 : "than the mesh received previously. Mesh for received when we got "
387 : "fluxes: "
388 : << boundary_message->volume_or_ghost_mesh
389 : << " mesh received when we got ghost cells "
390 : << current_boundary_data->volume_or_ghost_mesh);
391 :
392 : // Just need to set the pointer. No need to worry about the data being
393 : // deleted upon destruction of boundary_message because the destructor of
394 : // a BoundaryMessage will only delete the pointer stored inside the
395 : // BoundaryMessage. It won't delete the actual data it points to. Now the
396 : // unique_ptr owns the data that was previously pointed to by
397 : // boundary_message
398 : current_boundary_data->dg_flux_data_size =
399 : boundary_message->dg_flux_data_size;
400 : current_boundary_data->dg_flux_data = boundary_message->dg_flux_data;
401 : current_boundary_data->next_time_step_id =
402 : boundary_message->next_time_step_id;
403 : current_boundary_data->tci_status = boundary_message->tci_status;
404 : current_boundary_data->integration_order =
405 : boundary_message->integration_order;
406 : } else {
407 : // We have not received ghost cells or fluxes at this time.
408 : // Once we insert boundary_message into the unique_ptr we cannot use
409 : // boundary_message anymore because it is invalidated. The unique_ptr now
410 : // owns the memory.
411 : if (not current_inbox
412 : .insert(std::pair{key, std::unique_ptr<BoundaryMessage<Dim>>(
413 : boundary_message)})
414 : .second) {
415 : ERROR("Failed to insert data to receive at instance '"
416 : << time_step_id << "' with tag 'BoundaryMessageInbox'.\n");
417 : }
418 : }
419 : }
420 : };
421 : } // namespace evolution::dg::Tags
|
