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