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_set>
       9             : 
      10             : #include "DataStructures/DataVector.hpp"
      11             : #include "DataStructures/Tensor/Tensor.hpp"
      12             : #include "Domain/Structure/Direction.hpp"
      13             : #include "Utilities/Gsl.hpp"
      14             : 
      15             : namespace LinearSolver::Schwarz {
      16             : 
      17             : /*!
      18             :  * \brief Weights for the solution on an element-centered subdomain, decreasing
      19             :  * from 1 to 0.5 towards the `side` over the logical distance `width`, and
      20             :  * further to 0 over the same distance outside the element.
      21             :  *
      22             :  * The weighting function over a full element-centered subdomain is
      23             :  *
      24             :  * \f{equation}
      25             :  * w(\xi) = \frac{1}{2}\left( \phi\left( \frac{\xi + 1}{\delta} \right) -
      26             :  * \phi\left( \frac{\xi - 1}{\delta} \right) \right) \f}
      27             :  *
      28             :  * where \f$\phi(\xi)\f$ is a second-order `::smoothstep`, i.e. the quintic
      29             :  * polynomial
      30             :  *
      31             :  * \f{align*}
      32             :  * \phi(\xi) = \begin{cases} \mathrm{sign}(\xi) \quad \text{for}
      33             :  * \quad |\xi| > 1 \\
      34             :  * \frac{1}{8}\left(15\xi - 10\xi^3 + 3\xi^5\right) \end{cases}
      35             :  * \f}
      36             :  *
      37             :  * (see Eq. (39) in \cite Vincent2019qpd).
      38             :  *
      39             :  * The `LinearSolver::Schwarz::extruding_weight` and
      40             :  * `LinearSolver::Schwarz::intruding_weight` functions each compute one of the
      41             :  * two terms in \f$w(\xi)\f$. For example, consider an element-centered
      42             :  * subdomain `A` that overlaps with a neighboring element-centered subdomain
      43             :  * `B`. To combine solutions on `A` and `B` to a weighted solution on `A`,
      44             :  * multiply the solution on `A` with the `extruding_weight` and the solution on
      45             :  * `B` with the `intruding_weight`, both evaluated at the logical coordinates in
      46             :  * `A` and at the `side` of `A` that faces `B`.
      47             :  */
      48           1 : DataVector extruding_weight(const DataVector& logical_coords, double width,
      49             :                             const Side& side);
      50             : 
      51             : /// @{
      52             : /*!
      53             :  * \brief Weights for data on the central element of an element-centered
      54             :  * subdomain
      55             :  *
      56             :  * Constructs the weighting field
      57             :  *
      58             :  * \f{equation}
      59             :  * W(\boldsymbol{\xi}) = \prod^d_{i=0} w(\xi^i)
      60             :  * \f}
      61             :  *
      62             :  * where \f$w(\xi^i)\f$ is the one-dimensional weighting function described in
      63             :  * `LinearSolver::Schwarz::extruding_weight` and \f$\xi^i\f$ are the
      64             :  * element-logical coordinates (see Eq. (41) in \cite Vincent2019qpd).
      65             :  */
      66             : template <size_t Dim>
      67           1 : void element_weight(
      68             :     gsl::not_null<Scalar<DataVector>*> element_weight,
      69             :     const tnsr::I<DataVector, Dim, Frame::ElementLogical>& logical_coords,
      70             :     const std::array<double, Dim>& overlap_widths,
      71             :     const std::unordered_set<Direction<Dim>>& external_boundaries);
      72             : 
      73             : template <size_t Dim>
      74           1 : Scalar<DataVector> element_weight(
      75             :     const tnsr::I<DataVector, Dim, Frame::ElementLogical>& logical_coords,
      76             :     const std::array<double, Dim>& overlap_widths,
      77             :     const std::unordered_set<Direction<Dim>>& external_boundaries);
      78             : /// @}
      79             : 
      80             : /*!
      81             :  * \brief Weights for the intruding solution of a neighboring element-centered
      82             :  * subdomain, increasing from 0 to 0.5 towards the `side` over the logical
      83             :  * distance `width`, and further to 1 over the same distance outside the
      84             :  * element.
      85             :  *
      86             :  * \see `LinearSolver::Schwarz::extruding_weight`
      87             :  */
      88           1 : DataVector intruding_weight(const DataVector& logical_coords, double width,
      89             :                             const Side& side);
      90             : 
      91             : /// @{
      92             : /*!
      93             :  * \brief Weights for data on overlap regions intruding into an element-centered
      94             :  * subdomain
      95             :  *
      96             :  * Constructs the weighting field \f$W(\xi)\f$ as described in
      97             :  * `LinearSolver::Schwarz::element_weight` for the data that overlaps with the
      98             :  * central element of an element-centered subdomain. The weights are constructed
      99             :  * in such a way that all weights at a grid point sum to one, i.e. the weight is
     100             :  * conserved. The `logical_coords` are the element-logical coordinates of the
     101             :  * central element.
     102             :  *
     103             :  * This function assumes that corner- and edge-neighbors of the central element
     104             :  * are not part of the subdomain, which means that no contributions from those
     105             :  * neighbors are expected although the weighting field is non-zero in overlap
     106             :  * regions with those neighbors. Therefore, to retain conservation we must
     107             :  * account for this missing weight by adding it to the central element, to the
     108             :  * intruding overlaps from face-neighbors, or split it between the two. We
     109             :  * choose to add the weight to the intruding overlaps, since that's where
     110             :  * information from the corner- and edge-regions propagates through in a DG
     111             :  * context.
     112             :  */
     113             : template <size_t Dim>
     114           1 : void intruding_weight(
     115             :     gsl::not_null<Scalar<DataVector>*> weight,
     116             :     const tnsr::I<DataVector, Dim, Frame::ElementLogical>& logical_coords,
     117             :     const Direction<Dim>& direction,
     118             :     const std::array<double, Dim>& overlap_widths,
     119             :     size_t num_intruding_overlaps,
     120             :     const std::unordered_set<Direction<Dim>>& external_boundaries);
     121             : 
     122             : template <size_t Dim>
     123           1 : Scalar<DataVector> intruding_weight(
     124             :     const tnsr::I<DataVector, Dim, Frame::ElementLogical>& logical_coords,
     125             :     const Direction<Dim>& direction,
     126             :     const std::array<double, Dim>& overlap_widths,
     127             :     size_t num_intruding_overlaps,
     128             :     const std::unordered_set<Direction<Dim>>& external_boundaries);
     129             : /// @}
     130             : 
     131             : }  // namespace LinearSolver::Schwarz