10 #include "DataStructures/DataBox/Tag.hpp"
11 #include "DataStructures/Tags/TempTensor.hpp"
19 #include "ParallelAlgorithms/LinearSolver/Schwarz/OverlapHelpers.hpp"
20 #include "ParallelAlgorithms/LinearSolver/Schwarz/Tags.hpp"
21 #include "ParallelAlgorithms/LinearSolver/Schwarz/Weighting.hpp"
29 template <
size_t Dim,
typename OptionsGroup>
36 static constexpr
void function(
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) =
48 template <
size_t Dim,
typename OptionsGroup>
56 static constexpr
void function(
60 for (
size_t d = 0; d < Dim; ++d) {
63 gsl::at(*intruding_overlap_widths, d) =
71 template <
size_t Dim,
typename OptionsGroup>
74 using return_type =
typename base::type;
76 tmpl::list<domain::Tags::Element<Dim>,
80 static constexpr
void function(
83 const tnsr::I<DataVector, Dim, Frame::Logical>& logical_coords,
85 const size_t max_overlap) noexcept {
89 if (
LIKELY(max_overlap > 0)) {
91 intruding_overlap_widths,
92 element.external_boundaries());
94 *element_weight = make_with_value<return_type>(logical_coords, 1.);
101 template <
size_t Dim,
typename OptionsGroup>
104 using return_type =
typename base::type;
105 using argument_tags =
116 static constexpr
void function(
121 const tnsr::I<DataVector, Dim, Frame::Logical>& logical_coords,
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),
130 intruding_overlap_weight, intruding_logical_coords, direction,
131 intruding_overlap_widths, element.neighbors().at(direction).size(),
132 element.external_boundaries());
142 template <
size_t Dim,
typename 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(
155 all_intruding_weights,
159 mesh.number_of_grid_points());
160 get(*summed_intruding_overlap_weights) = 0.;
172 temp_vars, mesh.extents(),
173 gsl::at(all_intruding_extents, direction.dimension()), direction);
175 get(*summed_intruding_overlap_weights) +=
get(get<temp_tag>(temp_vars));
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
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
Mark a struct as a compute tag by inheriting from this.
Definition: Tag.hpp:157
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
constexpr Tag::type & get(Variables< TagList > &v) noexcept
Return Tag::type pointing into the contiguous array.
Definition: Variables.hpp:660
const DataVector & collocation_points(const size_t num_points) noexcept
Collocation points.
Definition: Spectral.cpp:283
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.
The computational grid of the Element in the DataBox.
Definition: Tags.hpp:107
Definition: Element.hpp:29
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
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
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
const DataVector & collocation_points(const Mesh< 1 > &mesh) noexcept
Collocation points for a one-dimensional mesh.
Holds the number of grid points, basis, and quadrature in each direction of the computational grid.
Definition: Mesh.hpp:48
Tensor< T, Symmetry<>, index_list<> > Scalar
Definition: TypeAliases.hpp:21
#define LIKELY(x)
Definition: Gsl.hpp:67
Require a pointer to not be a nullptr
Definition: ReadSpecThirdOrderPiecewisePolynomial.hpp:13
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