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 <cmath>
       8             : #include <cstddef>
       9             : #include <utility>
      10             : 
      11             : #include "NumericalAlgorithms/FiniteDifference/Reconstruct.hpp"
      12             : #include "Utilities/ForceInline.hpp"
      13             : #include "Utilities/Gsl.hpp"
      14             : #include "Utilities/Math.hpp"
      15             : #include "Utilities/TMPL.hpp"
      16             : 
      17             : /// \cond
      18             : template <size_t Dim>
      19             : class Direction;
      20             : template <size_t Dim>
      21             : class Index;
      22             : /// \endcond
      23             : 
      24             : namespace fd::reconstruction {
      25             : namespace detail {
      26             : struct MinmodReconstructor {
      27             :   SPECTRE_ALWAYS_INLINE static std::array<double, 2> pointwise(
      28             :       const double* const q, const int stride) {
      29             :     using std::min;
      30             :     using std::abs;
      31             :     // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic)
      32             :     const double a = q[stride] - q[0];
      33             :     // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic)
      34             :     const double b = q[0] - q[-stride];
      35             :     const double slope = 0.5 * (sgn(a) + sgn(b)) * min(abs(a), abs(b));
      36             :     // NOLINTNEXTLINE(cppcoreguidelines-pro-bounds-pointer-arithmetic)
      37             :     return {{q[0] - 0.5 * slope, q[0] + 0.5 * slope}};
      38             :   }
      39             : 
      40             :   SPECTRE_ALWAYS_INLINE static constexpr size_t stencil_width() { return 3; }
      41             : };
      42             : }  // namespace detail
      43             : 
      44             : /*!
      45             :  * \ingroup FiniteDifferenceGroup
      46             :  * \brief Performs minmod reconstruction on the `volume_vars` in each direction.
      47             :  *
      48             :  * On a 1d mesh with spacing \f$\Delta x\f$ we denote the solution at the
      49             :  * \f$j\f$th point by \f$u_j\f$. The minmod reconstruction \cite RezzollaBook
      50             :  * first computes:
      51             :  *
      52             :  * \f{align}
      53             :  * \sigma_j\equiv \textrm{minmod}\left(\frac{u_j-u_{j-1}}{\Delta x},
      54             :  *           \frac{u_{j+1}-u_j}{\Delta x}\right)
      55             :  * \f}
      56             :  *
      57             :  * where
      58             :  *
      59             :  * \f{align}
      60             :  *  \mathrm{minmod}(a,b)=
      61             :  * \left\{
      62             :  *   \begin{array}{ll}
      63             :  *    \mathrm{sgn}(a)\min(\lvert a\rvert, \lvert b\rvert)
      64             :  *      & \mathrm{if} \; \mathrm{sgn}(a)=\mathrm{sgn}(b) \\
      65             :  *    0 & \mathrm{otherwise}
      66             :  *   \end{array}\right.
      67             :  * \f}
      68             :  *
      69             :  * The reconstructed solution \f$u_j(x)\f$ in the \f$j\f$th cell is given by
      70             :  *
      71             :  * \f{align}
      72             :  * u_j(x)=u_j + \sigma_j (x-x_j),
      73             :  * \f}
      74             :  *
      75             :  * where \f$x_j\f$ is the coordinate \f$x\f$ at the center of the cell.
      76             :  */
      77             : template <size_t Dim>
      78           1 : void minmod(const gsl::not_null<std::array<gsl::span<double>, Dim>*>
      79             :                 reconstructed_upper_side_of_face_vars,
      80             :             const gsl::not_null<std::array<gsl::span<double>, Dim>*>
      81             :                 reconstructed_lower_side_of_face_vars,
      82             :             const gsl::span<const double>& volume_vars,
      83             :             const DirectionMap<Dim, gsl::span<const double>>& ghost_cell_vars,
      84             :             const Index<Dim>& volume_extents,
      85             :             const size_t number_of_variables) {
      86             :   detail::reconstruct<detail::MinmodReconstructor>(
      87             :       reconstructed_upper_side_of_face_vars,
      88             :       reconstructed_lower_side_of_face_vars, volume_vars, ghost_cell_vars,
      89             :       volume_extents, number_of_variables);
      90             : }
      91             : }  // namespace fd::reconstruction