ComputeTags.hpp
1 // 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"
16 #include "Domain/Tags.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"
24 #include "Utilities/TMPL.hpp"
25 
27 
28 /// The number of points a neighbor's subdomain extends into the element
29 template <size_t Dim, typename OptionsGroup>
31  IntrudingExtents<Dim, OptionsGroup> {
33  using return_type = typename base::type;
34  using argument_tags =
35  tmpl::list<domain::Tags::Mesh<Dim>, MaxOverlap<OptionsGroup>>;
36  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>
51  IntrudingOverlapWidths<Dim, OptionsGroup> {
53  using return_type = typename base::type;
54  using argument_tags =
55  tmpl::list<domain::Tags::Mesh<Dim>, IntrudingExtents<Dim, OptionsGroup>>;
56  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),
66  }
67  }
68 };
69 
70 /// Weighting field for data on the element
71 template <size_t Dim, typename OptionsGroup>
72 struct ElementWeightCompute : db::ComputeTag, Weight<OptionsGroup> {
73  using base = Weight<OptionsGroup>;
74  using return_type = typename base::type;
75  using argument_tags =
76  tmpl::list<domain::Tags::Element<Dim>,
80  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>
103  using base = Weight<OptionsGroup>;
104  using return_type = typename base::type;
105  using argument_tags =
106  tmpl::list<domain::Tags::Mesh<Dim>, domain::Tags::Element<Dim>,
111  using volume_tags =
112  tmpl::list<domain::Tags::Mesh<Dim>, domain::Tags::Element<Dim>,
116  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 =
127  logical_coords, mesh.extents(), gsl::at(intruding_extents, dim),
128  direction);
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>
144  : db::ComputeTag,
145  SummedIntrudingOverlapWeights<OptionsGroup> {
147  using return_type = typename base::type;
148  using argument_tags =
149  tmpl::list<domain::Tags::Interface<domain::Tags::InternalDirections<Dim>,
152  static void function(
153  const gsl::not_null<return_type*> summed_intruding_overlap_weights,
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;
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
LinearSolver::Schwarz::intruding_weight
DataVector intruding_weight(const DataVector &logical_coords, const double width, const Side &side) noexcept
Weights for the intruding solution of a neighboring element-centered subdomain, increasing from 0 to ...
Definition: Weighting.cpp:79
gsl::at
constexpr T & at(std::array< T, N > &arr, Size index)
Retrieve a entry from a container, with checks in Debug mode that the index being retrieved is valid.
Definition: Gsl.hpp:125
db::ComputeTag
Mark a struct as a compute tag by inheriting from this.
Definition: Tag.hpp:157
LinearSolver::Schwarz::Tags::SummedIntrudingOverlapWeightsCompute
A diagnostic quantity to check that weights are conserved.
Definition: ComputeTags.hpp:143
LinearSolver::Schwarz::overlap_extent
size_t overlap_extent(const size_t volume_extent, const size_t max_overlap) noexcept
The number of points that an overlap extends into the volume_extent
Definition: OverlapHelpers.cpp:21
get
constexpr Tag::type & get(Variables< TagList > &v) noexcept
Return Tag::type pointing into the contiguous array.
Definition: Variables.hpp:638
domain::Tags::Coordinates
Definition: Tags.hpp:130
Spectral::collocation_points
const DataVector & collocation_points(const size_t num_points) noexcept
Collocation points.
Definition: Spectral.cpp:282
LinearSolver::Schwarz::element_weight
void element_weight(gsl::not_null< Scalar< DataVector > * > element_weight, const tnsr::I< DataVector, Dim, Frame::Logical > &logical_coords, const std::array< double, Dim > &overlap_widths, const std::unordered_set< Direction< Dim >> &external_boundaries) noexcept
Weights for data on the central element of an element-centered subdomain.
Tags.hpp
domain::Tags::Element
Definition: Tags.hpp:97
LinearSolver::Schwarz::Tags::ElementWeightCompute
Weighting field for data on the element.
Definition: ComputeTags.hpp:72
LinearSolver::Schwarz::Tags::MaxOverlap
Number of points a subdomain can overlap with its neighbor.
Definition: Tags.hpp:48
LinearSolver::Schwarz::Tags::IntrudingOverlapWeightCompute
Weighting field for data on neighboring subdomains that overlap with the element. Intended to be used...
Definition: ComputeTags.hpp:102
domain::Tags::Mesh
The computational grid of the Element in the DataBox.
Definition: Tags.hpp:107
LinearSolver::Schwarz::Tags::IntrudingExtents
The number of points a neighbor's subdomain extends into the element.
Definition: Tags.hpp:97
Spectral.hpp
Direction
Definition: Direction.hpp:23
Element
Definition: Element.hpp:29
LinearSolver::Schwarz::Tags::IntrudingExtentsCompute
The number of points a neighbor's subdomain extends into the element.
Definition: ComputeTags.hpp:30
LinearSolver::Schwarz::overlap_width
double overlap_width(const size_t overlap_extent, const DataVector &collocation_points) noexcept
Width of an overlap extending overlap_extent points into the collocation_points from either side.
Definition: OverlapHelpers.cpp:49
cstddef
LinearSolver::Schwarz::data_on_overlap
void data_on_overlap(const gsl::not_null< Tensor< DataType, TensorStructure... > * > restricted_tensor, const Tensor< DataType, TensorStructure... > &tensor, const Index< Dim > &volume_extents, const size_t overlap_extent, const Direction< Dim > &direction) noexcept
The part of the tensor data that lies within the overlap region.
Definition: OverlapHelpers.hpp:139
MakeWithValue.hpp
array
LinearSolver::Schwarz::Tags::IntrudingOverlapWidthsCompute
The width in element-logical coordinates that a neighbor's subdomain extends into the element.
Definition: ComputeTags.hpp:49
LinearSolver::Schwarz::Tags::Weight
Weighting field for combining data from multiple overlapping subdomains.
Definition: Tags.hpp:116
LinearSolver::Schwarz::extended_overlap_data
void extended_overlap_data(const gsl::not_null< Variables< ExtendedTagsList > * > extended_data, const Variables< OverlapTagsList > &overlap_data, const Index< Dim > &volume_extents, const size_t overlap_extent, const Direction< Dim > &direction) noexcept
Extend the overlap data to the full mesh by filling it with zeros outside the overlap region.
Definition: OverlapHelpers.hpp:253
Spectral::collocation_points
const DataVector & collocation_points(const Mesh< 1 > &mesh) noexcept
Collocation points for a one-dimensional mesh.
Variables.hpp
Element.hpp
Mesh
Holds the number of grid points, basis, and quadrature in each direction of the computational grid.
Definition: Mesh.hpp:47
LinearSolver::Schwarz::Tags::SummedIntrudingOverlapWeights
A diagnostic quantity to check that weights are conserved.
Definition: Tags.hpp:131
Scalar
Tensor< T, Symmetry<>, index_list<> > Scalar
Definition: TypeAliases.hpp:21
Tags::TempTensor
Definition: TempTensor.hpp:21
Gsl.hpp
LinearSolver::Schwarz::Tags
Tags related to the Schwarz solver.
Definition: ComputeTags.hpp:26
domain::Tags::Direction
Definition: Tags.hpp:381
Tensor.hpp
LIKELY
#define LIKELY(x)
Definition: Gsl.hpp:67
Direction.hpp
Mesh.hpp
unordered_map
TMPL.hpp
gsl::not_null
Require a pointer to not be a nullptr
Definition: Gsl.hpp:183
destructive_resize_components
void destructive_resize_components(const gsl::not_null< Container * > container, const size_t new_size, DestructiveResizeFunction destructive_resize=ContainerDestructiveResize{}) noexcept
Checks the size of each component of the container, and resizes if necessary.
Definition: ContainerHelpers.hpp:177
LinearSolver::Schwarz::Tags::IntrudingOverlapWidths
The width in element-logical coordinates that a neighbor's subdomain extends into the element.
Definition: Tags.hpp:107