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 <functional>
9 : #include <optional>
10 : #include <unordered_map>
11 : #include <unordered_set>
12 : #include <vector>
13 :
14 : #include "Domain/Block.hpp"
15 : #include "Domain/ElementDistribution.hpp"
16 : #include "Domain/Structure/ElementId.hpp"
17 : #include "Domain/Structure/InitialElementIds.hpp"
18 : #include "Parallel/DomainDiagnosticInfo.hpp"
19 : #include "Parallel/GlobalCache.hpp"
20 : #include "Utilities/ErrorHandling/Assert.hpp"
21 : #include "Utilities/Numeric.hpp"
22 :
23 : /// \cond
24 : namespace Spectral {
25 : enum class Basis : uint8_t;
26 : enum class Quadrature : uint8_t;
27 : } // namespace Spectral
28 : /// \endcond
29 :
30 : namespace Parallel {
31 : /*!
32 : * \brief Creates elements using a chosen distribution.
33 : *
34 : * The `func` is called with `(element_id, target_proc, target_node)` allowing
35 : * the `func` to insert the element with `element_id` on the target processor
36 : * and node.
37 : *
38 : * The `weighting_extents_override`, if present, gives per-block extents to use
39 : * instead of `initial_extents` only when computing the cost of elements for
40 : * the distribution (e.g. to weight subcell-capable elements by their
41 : * finite-difference grid rather than their DG grid). Only blocks whose
42 : * weighting extents differ from `initial_extents` need an entry; blocks not
43 : * found in the map use `initial_extents`. The actual `initial_extents` are
44 : * always used for everything else, e.g. the diagnostics below.
45 : */
46 : template <typename F, size_t Dim, typename Metavariables>
47 1 : void create_elements_using_distribution(
48 : const F& func, const std::optional<domain::ElementWeight>& element_weight,
49 : const std::vector<Block<Dim>>& blocks,
50 : const std::vector<std::array<size_t, Dim>>& initial_extents,
51 : const std::vector<std::array<size_t, Dim>>& initial_refinement_levels,
52 : const Spectral::Basis i1_basis, const Spectral::Quadrature i1_quadrature,
53 :
54 : const std::unordered_set<size_t>& procs_to_ignore,
55 : const size_t number_of_procs, const size_t number_of_nodes,
56 : const size_t num_of_procs_to_use,
57 : const Parallel::GlobalCache<Metavariables>& local_cache,
58 : const bool print_diagnostics,
59 : const std::optional<std::unordered_map<size_t, std::array<size_t, Dim>>>&
60 : weighting_extents_override = std::nullopt) {
61 : // Only need the element distribution if the element weight has a value
62 : // because then we have to use the space filling curve and not just use round
63 : // robin.
64 : domain::BlockZCurveProcDistribution<Dim> element_distribution{};
65 : if (element_weight.has_value()) {
66 : std::optional<std::vector<std::array<size_t, Dim>>>
67 : weighting_extents_buffer{};
68 : if (weighting_extents_override.has_value()) {
69 : weighting_extents_buffer = initial_extents;
70 : const size_t num_blocks = weighting_extents_buffer->size();
71 : for (const auto& block_id_and_extents :
72 : weighting_extents_override.value()) {
73 : ASSERT(block_id_and_extents.first < num_blocks,
74 : "The weighting extents override contains block ID "
75 : << block_id_and_extents.first << " but there are only "
76 : << num_blocks << " blocks.");
77 : (*weighting_extents_buffer)[block_id_and_extents.first] =
78 : block_id_and_extents.second;
79 : }
80 : }
81 : const std::vector<std::array<size_t, Dim>>& weighting_extents =
82 : weighting_extents_buffer.has_value() ? weighting_extents_buffer.value()
83 : : initial_extents;
84 : const std::unordered_map<ElementId<Dim>, double> element_costs =
85 : domain::get_element_costs(blocks, initial_refinement_levels,
86 : weighting_extents, element_weight.value(),
87 : i1_basis, i1_quadrature);
88 : element_distribution = domain::BlockZCurveProcDistribution<Dim>{
89 : element_costs, num_of_procs_to_use, blocks, initial_refinement_levels,
90 : procs_to_ignore};
91 : }
92 :
93 : // Will be used to print domain diagnostic info
94 : std::vector<size_t> elements_per_core(number_of_procs, 0_st);
95 : std::vector<size_t> elements_per_node(number_of_nodes, 0_st);
96 : std::vector<size_t> grid_points_per_core(number_of_procs, 0_st);
97 : std::vector<size_t> grid_points_per_node(number_of_nodes, 0_st);
98 :
99 : size_t which_proc = 0;
100 : for (const auto& block : blocks) {
101 : const auto& initial_ref_levs = initial_refinement_levels[block.id()];
102 : const size_t grid_points_per_element =
103 : alg::accumulate(initial_extents[block.id()], 1_st, std::multiplies<>());
104 :
105 : const std::vector<ElementId<Dim>> element_ids =
106 : initial_element_ids(block.id(), initial_ref_levs);
107 :
108 : // Value means ZCurve. nullopt means round robin
109 : if (element_weight.has_value()) {
110 : for (const auto& element_id : element_ids) {
111 : const size_t target_proc =
112 : element_distribution.get_proc_for_element(element_id);
113 : const size_t target_node =
114 : Parallel::node_of<size_t>(target_proc, local_cache);
115 : func(element_id, target_proc, target_node);
116 :
117 : ++elements_per_core[target_proc];
118 : ++elements_per_node[target_node];
119 : grid_points_per_core[target_proc] += grid_points_per_element;
120 : grid_points_per_node[target_node] += grid_points_per_element;
121 : }
122 : } else {
123 : for (size_t i = 0; i < element_ids.size(); ++i) {
124 : while (procs_to_ignore.find(which_proc) != procs_to_ignore.end()) {
125 : which_proc = which_proc + 1 == number_of_procs ? 0 : which_proc + 1;
126 : }
127 : const size_t target_proc = which_proc;
128 : const size_t target_node =
129 : Parallel::node_of<size_t>(which_proc, local_cache);
130 : const ElementId<Dim> element_id(element_ids[i]);
131 : func(element_id, target_proc, target_node);
132 :
133 : ++elements_per_core[which_proc];
134 : ++elements_per_node[target_node];
135 : grid_points_per_core[which_proc] += grid_points_per_element;
136 : grid_points_per_node[target_node] += grid_points_per_element;
137 :
138 : which_proc = which_proc + 1 == number_of_procs ? 0 : which_proc + 1;
139 : }
140 : }
141 : }
142 :
143 : if (print_diagnostics) {
144 : Parallel::printf("\n%s\n", domain::diagnostic_info(
145 : blocks.size(), local_cache,
146 : elements_per_core, elements_per_node,
147 : grid_points_per_core, grid_points_per_node));
148 : }
149 : }
150 : } // namespace Parallel
|
