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 <memory>
8 : #include <unordered_map>
9 : #include <utility>
10 :
11 : #include "DataStructures/TaggedTuple.hpp"
12 : #include "DataStructures/Tags.hpp"
13 : #include "Domain/Structure/Direction.hpp"
14 : #include "Domain/Structure/DirectionMap.hpp"
15 : #include "Domain/Structure/DirectionalId.hpp"
16 : #include "Utilities/Gsl.hpp"
17 :
18 : /// \cond
19 : class DataVector;
20 : template <size_t VolumeDim>
21 : class Element;
22 : template <size_t VolumeDim>
23 : class ElementId;
24 : template <size_t VolumeDim>
25 : class Mesh;
26 : enum class Side : uint8_t;
27 :
28 : namespace boost {
29 : template <class T>
30 : struct hash;
31 : } // namespace boost
32 : /// \endcond
33 :
34 : namespace Limiters::Minmod_detail {
35 :
36 : // Encodes the return status of the tvb_corrected_minmod function.
37 : struct MinmodResult {
38 : const double value;
39 : const bool activated;
40 : };
41 :
42 : // The TVB-corrected minmod function, see e.g. Cockburn reference Eq. 2.26.
43 : MinmodResult tvb_corrected_minmod(double a, double b, double c,
44 : double tvb_scale);
45 :
46 : // Holds various optimization-related allocations for the Minmod TCI.
47 : // There is no pup::er, because these allocations should be short-lived (i.e.,
48 : // scoped within a single limiter invocation).
49 : template <size_t VolumeDim>
50 : class BufferWrapper {
51 : public:
52 : BufferWrapper() = delete;
53 : explicit BufferWrapper(const Mesh<VolumeDim>& mesh);
54 :
55 : private:
56 : std::unique_ptr<double[]> contiguous_boundary_buffer_;
57 : const std::pair<std::unique_ptr<std::pair<size_t, size_t>[]>,
58 : std::array<std::pair<gsl::span<std::pair<size_t, size_t>>,
59 : gsl::span<std::pair<size_t, size_t>>>,
60 : VolumeDim>>
61 : volume_and_slice_buffer_and_indices_ {};
62 :
63 : public:
64 : std::array<DataVector, VolumeDim> boundary_buffers{};
65 : const std::array<std::pair<gsl::span<std::pair<size_t, size_t>>,
66 : gsl::span<std::pair<size_t, size_t>>>,
67 : VolumeDim>& volume_and_slice_indices{};
68 : };
69 :
70 : // In each direction, average the size of all different neighbors in that
71 : // direction. Note that only the component of neighor_size that is normal
72 : // to the face is needed (and, therefore, computed).
73 : //
74 : // Expects type `PackagedData` to contain a variable `element_size` that is a
75 : // `std::array<double, VolumeDim>`.
76 : template <size_t VolumeDim, typename PackagedData>
77 : DirectionMap<VolumeDim, double> compute_effective_neighbor_sizes(
78 : const Element<VolumeDim>& element,
79 : const std::unordered_map<DirectionalId<VolumeDim>, PackagedData,
80 : boost::hash<DirectionalId<VolumeDim>>>&
81 : neighbor_data) {
82 : DirectionMap<VolumeDim, double> result;
83 : for (const auto& dir : Direction<VolumeDim>::all_directions()) {
84 : const bool neighbors_in_this_dir = element.neighbors().contains(dir);
85 : if (neighbors_in_this_dir) {
86 : const double effective_neighbor_size = [&dir, &element,
87 : &neighbor_data]() {
88 : const size_t dim = dir.dimension();
89 : const auto& neighbor_ids = element.neighbors().at(dir).ids();
90 : double size_accumulate = 0.;
91 : for (const auto& id : neighbor_ids) {
92 : size_accumulate += gsl::at(
93 : neighbor_data.at(DirectionalId<VolumeDim>{dir, id}).element_size,
94 : dim);
95 : }
96 : return size_accumulate / neighbor_ids.size();
97 : }();
98 : result.insert(std::make_pair(dir, effective_neighbor_size));
99 : }
100 : }
101 : return result;
102 : }
103 :
104 : // In each direction, average the mean of the specified tensor component over
105 : // all different neighbors that direction. This produces one effective neighbor
106 : // per direction.
107 : //
108 : // Expects type `PackagedData` to contain a variable `means` that is a
109 : // `TaggedTuple<Tags::Mean<Tags>...>`. Tags... must contain Tag, the tag
110 : // specifying the tensor to work with.
111 : template <typename Tag, size_t VolumeDim, typename PackagedData>
112 : DirectionMap<VolumeDim, double> compute_effective_neighbor_means(
113 : const size_t tensor_storage_index, const Element<VolumeDim>& element,
114 : const std::unordered_map<DirectionalId<VolumeDim>, PackagedData,
115 : boost::hash<DirectionalId<VolumeDim>>>&
116 : neighbor_data) {
117 : DirectionMap<VolumeDim, double> result;
118 : for (const auto& dir : Direction<VolumeDim>::all_directions()) {
119 : const bool neighbors_in_this_dir = element.neighbors().contains(dir);
120 : if (neighbors_in_this_dir) {
121 : const double effective_neighbor_mean =
122 : [&dir, &element, &neighbor_data, &tensor_storage_index]() {
123 : const auto& neighbor_ids = element.neighbors().at(dir).ids();
124 : double mean_accumulate = 0.0;
125 : for (const auto& id : neighbor_ids) {
126 : mean_accumulate += tuples::get<::Tags::Mean<Tag>>(
127 : neighbor_data.at(DirectionalId<VolumeDim>{dir, id})
128 : .means)[tensor_storage_index];
129 : }
130 : return mean_accumulate / neighbor_ids.size();
131 : }();
132 : result.insert(std::make_pair(dir, effective_neighbor_mean));
133 : }
134 : }
135 : return result;
136 : }
137 :
138 : // Compute an effective element-center-to-neighbor-center distance that accounts
139 : // for the possibility of different refinement levels or discontinuous maps
140 : // (e.g., at Block boundaries). Treated naively, these domain features can make
141 : // a smooth solution appear to be non-smooth in the logical coordinates, which
142 : // could potentially lead to the limiter triggering erroneously. This effective
143 : // distance is used to scale the difference in the means, so that a linear
144 : // function at a refinement or Block boundary will still appear smooth to the
145 : // limiter. The factor is normalized to be 1.0 on a uniform grid.
146 : //
147 : // Note that this is not "by the book" Minmod, but an attempt to
148 : // generalize Minmod to work on non-uniform grids.
149 : template <size_t VolumeDim>
150 : double effective_difference_to_neighbor(
151 : double u_mean, const Element<VolumeDim>& element,
152 : const std::array<double, VolumeDim>& element_size, size_t dim,
153 : const Side& side,
154 : const DirectionMap<VolumeDim, double>& effective_neighbor_means,
155 : const DirectionMap<VolumeDim, double>& effective_neighbor_sizes);
156 :
157 : } // namespace Limiters::Minmod_detail
|
