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 <unordered_map>
9 :
10 : #include "DataStructures/DataBox/Tag.hpp"
11 : #include "DataStructures/Tags/TempTensor.hpp"
12 : #include "DataStructures/Tensor/Tensor.hpp"
13 : #include "DataStructures/Variables.hpp"
14 : #include "Domain/Structure/Direction.hpp"
15 : #include "Domain/Structure/Element.hpp"
16 : #include "Domain/Tags.hpp"
17 : #include "NumericalAlgorithms/Spectral/Mesh.hpp"
18 : #include "NumericalAlgorithms/Spectral/Spectral.hpp"
19 : #include "ParallelAlgorithms/LinearSolver/Schwarz/OverlapHelpers.hpp"
20 : #include "ParallelAlgorithms/LinearSolver/Schwarz/Tags.hpp"
21 : #include "ParallelAlgorithms/LinearSolver/Schwarz/Weighting.hpp"
22 : #include "Utilities/Gsl.hpp"
23 : #include "Utilities/MakeWithValue.hpp"
24 : #include "Utilities/TMPL.hpp"
25 :
26 1 : namespace LinearSolver::Schwarz::Tags {
27 :
28 : /// The number of points a neighbor's subdomain extends into the element
29 : template <size_t Dim, typename OptionsGroup>
30 : struct IntrudingExtentsCompute : db::ComputeTag,
31 1 : IntrudingExtents<Dim, OptionsGroup> {
32 0 : using base = IntrudingExtents<Dim, OptionsGroup>;
33 0 : using return_type = typename base::type;
34 0 : using argument_tags =
35 : tmpl::list<domain::Tags::Mesh<Dim>, MaxOverlap<OptionsGroup>>;
36 0 : static constexpr void function(
37 : const gsl::not_null<return_type*> intruding_extents,
38 : const Mesh<Dim>& mesh, const size_t max_overlap) noexcept {
39 : for (size_t d = 0; d < Dim; ++d) {
40 : gsl::at(*intruding_extents, d) =
41 : LinearSolver::Schwarz::overlap_extent(mesh.extents(d), max_overlap);
42 : }
43 : }
44 : };
45 :
46 : /// The width in element-logical coordinates that a neighbor's subdomain extends
47 : /// into the element
48 : template <size_t Dim, typename OptionsGroup>
49 : struct IntrudingOverlapWidthsCompute
50 : : db::ComputeTag,
51 1 : IntrudingOverlapWidths<Dim, OptionsGroup> {
52 0 : using base = IntrudingOverlapWidths<Dim, OptionsGroup>;
53 0 : using return_type = typename base::type;
54 0 : using argument_tags =
55 : tmpl::list<domain::Tags::Mesh<Dim>, IntrudingExtents<Dim, OptionsGroup>>;
56 0 : static constexpr void function(
57 : const gsl::not_null<return_type*> intruding_overlap_widths,
58 : const Mesh<Dim>& mesh,
59 : const std::array<size_t, Dim>& intruding_extents) noexcept {
60 : for (size_t d = 0; d < Dim; ++d) {
61 : const auto& collocation_points =
62 : Spectral::collocation_points(mesh.slice_through(d));
63 : gsl::at(*intruding_overlap_widths, d) =
64 : LinearSolver::Schwarz::overlap_width(gsl::at(intruding_extents, d),
65 : collocation_points);
66 : }
67 : }
68 : };
69 :
70 : /// Weighting field for data on the element
71 : template <size_t Dim, typename OptionsGroup>
72 1 : struct ElementWeightCompute : db::ComputeTag, Weight<OptionsGroup> {
73 0 : using base = Weight<OptionsGroup>;
74 0 : using return_type = typename base::type;
75 0 : using argument_tags =
76 : tmpl::list<domain::Tags::Element<Dim>,
77 : domain::Tags::Coordinates<Dim, Frame::Logical>,
78 : IntrudingOverlapWidths<Dim, OptionsGroup>,
79 : MaxOverlap<OptionsGroup>>;
80 0 : static constexpr void function(
81 : const gsl::not_null<return_type*> element_weight,
82 : const Element<Dim>& element,
83 : const tnsr::I<DataVector, Dim, Frame::Logical>& logical_coords,
84 : const std::array<double, Dim>& intruding_overlap_widths,
85 : const size_t max_overlap) noexcept {
86 : // For max_overlap > 0 all overlaps will have non-zero extents on an LGL
87 : // mesh (because it has at least 2 points per dimension), so we don't need
88 : // to check their extents are non-zero individually
89 : if (LIKELY(max_overlap > 0)) {
90 : LinearSolver::Schwarz::element_weight(element_weight, logical_coords,
91 : intruding_overlap_widths,
92 : element.external_boundaries());
93 : } else {
94 : *element_weight = make_with_value<return_type>(logical_coords, 1.);
95 : }
96 : }
97 : };
98 :
99 : /// Weighting field for data on neighboring subdomains that overlap with the
100 : /// element. Intended to be used with `domain::Tags::InterfaceCompute`.
101 : template <size_t Dim, typename OptionsGroup>
102 1 : struct IntrudingOverlapWeightCompute : db::ComputeTag, Weight<OptionsGroup> {
103 0 : using base = Weight<OptionsGroup>;
104 0 : using return_type = typename base::type;
105 0 : using argument_tags =
106 : tmpl::list<domain::Tags::Mesh<Dim>, domain::Tags::Element<Dim>,
107 : IntrudingExtents<Dim, OptionsGroup>,
108 : domain::Tags::Direction<Dim>,
109 : domain::Tags::Coordinates<Dim, Frame::Logical>,
110 : IntrudingOverlapWidths<Dim, OptionsGroup>>;
111 0 : using volume_tags =
112 : tmpl::list<domain::Tags::Mesh<Dim>, domain::Tags::Element<Dim>,
113 : IntrudingExtents<Dim, OptionsGroup>,
114 : domain::Tags::Coordinates<Dim, Frame::Logical>,
115 : IntrudingOverlapWidths<Dim, OptionsGroup>>;
116 0 : static constexpr void function(
117 : const gsl::not_null<return_type*> intruding_overlap_weight,
118 : const Mesh<Dim>& mesh, const Element<Dim>& element,
119 : const std::array<size_t, Dim>& intruding_extents,
120 : const Direction<Dim>& direction,
121 : const tnsr::I<DataVector, Dim, Frame::Logical>& logical_coords,
122 : const std::array<double, Dim> intruding_overlap_widths) noexcept {
123 : const size_t dim = direction.dimension();
124 : if (gsl::at(intruding_extents, dim) > 0) {
125 : const auto intruding_logical_coords =
126 : LinearSolver::Schwarz::data_on_overlap(
127 : logical_coords, mesh.extents(), gsl::at(intruding_extents, dim),
128 : direction);
129 : LinearSolver::Schwarz::intruding_weight(
130 : intruding_overlap_weight, intruding_logical_coords, direction,
131 : intruding_overlap_widths, element.neighbors().at(direction).size(),
132 : element.external_boundaries());
133 : }
134 : }
135 : };
136 :
137 : /*!
138 : * \brief A diagnostic quantity to check that weights are conserved
139 : *
140 : * \see `LinearSolver::Schwarz::Tags::SummedIntrudingOverlapWeights`
141 : */
142 : template <size_t Dim, typename OptionsGroup>
143 : struct SummedIntrudingOverlapWeightsCompute
144 : : db::ComputeTag,
145 1 : SummedIntrudingOverlapWeights<OptionsGroup> {
146 0 : using base = SummedIntrudingOverlapWeights<OptionsGroup>;
147 0 : using return_type = typename base::type;
148 0 : using argument_tags =
149 : tmpl::list<domain::Tags::Interface<domain::Tags::InternalDirections<Dim>,
150 : Weight<OptionsGroup>>,
151 : domain::Tags::Mesh<Dim>, IntrudingExtents<Dim, OptionsGroup>>;
152 0 : static void function(
153 : const gsl::not_null<return_type*> summed_intruding_overlap_weights,
154 : const std::unordered_map<Direction<Dim>, Scalar<DataVector>>&
155 : all_intruding_weights,
156 : const Mesh<Dim>& mesh,
157 : const std::array<size_t, Dim>& all_intruding_extents) noexcept {
158 : destructive_resize_components(summed_intruding_overlap_weights,
159 : mesh.number_of_grid_points());
160 : get(*summed_intruding_overlap_weights) = 0.;
161 : for (const auto& [direction, intruding_weight] : all_intruding_weights) {
162 : // Extend intruding weight to full extents
163 : // There's not function to extend a single tensor, so we create a
164 : // temporary Variables. This is only a diagnostic quantity that won't be
165 : // used in production code, so it doesn't need to be particularly
166 : // efficient and we don't need to add an overload of
167 : // `LinearSolver::Schwarz::extended_overlap_data` just for this purpose.
168 : using temp_tag = ::Tags::TempScalar<0>;
169 : Variables<tmpl::list<temp_tag>> temp_vars{get(intruding_weight).size()};
170 : get<temp_tag>(temp_vars) = intruding_weight;
171 : temp_vars = LinearSolver::Schwarz::extended_overlap_data(
172 : temp_vars, mesh.extents(),
173 : gsl::at(all_intruding_extents, direction.dimension()), direction);
174 : // Contribute to conserved weight
175 : get(*summed_intruding_overlap_weights) += get(get<temp_tag>(temp_vars));
176 : }
177 : }
178 : };
179 :
180 : } // namespace LinearSolver::Schwarz::Tags
|