Minmod.hpp
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"
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 {
28  const double* const q, const int stride) noexcept {
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() noexcept {
41  return 3;
42  }
43 };
44 } // namespace detail
45
46 /*!
47  * \ingroup FiniteDifferenceGroup
48  * \brief Performs minmod reconstruction on the volume_vars in each direction.
49  *
50  * On a 1d mesh with spacing \f$\Delta x\f$ we denote the solution at the
51  * \f$j\f$th point by \f$u_j\f$. The minmod reconstruction \cite RezzollaBook
52  * first computes:
53  *
54  * \f{align}
55  * \sigma_j\equiv \textrm{minmod}\left(\frac{u_j-u_{j-1}}{\Delta x},
56  * \frac{u_{j+1}-u_j}{\Delta x}\right)
57  * \f}
58  *
59  * where
60  *
61  * \f{align}
62  * \mathrm{minmod}(a,b)=
63  * \left\{
64  * \begin{array}{ll}
65  * \mathrm{sgn}(a)\min(\lvert a\rvert, \lvert b\rvert)
66  * & \mathrm{if} \; \mathrm{sgn}(a)=\mathrm{sgn}(b) \\
67  * 0 & \mathrm{otherwise}
68  * \end{array}\right.
69  * \f}
70  *
71  * The reconstructed solution \f$u_j(x)\f$ in the \f$j\f$th cell is given by
72  *
73  * \f{align}
74  * u_j(x)=u_j + \sigma_j (x-x_j),
75  * \f}
76  *
77  * where \f$x_j\f$ is the coordinate \f$x\f$ at the center of the cell.
78  */
79 template <size_t Dim>
81  reconstructed_upper_side_of_face_vars,
83  reconstructed_lower_side_of_face_vars,
84  const gsl::span<const double>& volume_vars,
85  const DirectionMap<Dim, gsl::span<const double>>& ghost_cell_vars,
86  const Index<Dim>& volume_extents,
87  const size_t number_of_variables) noexcept {
88  detail::reconstruct<detail::MinmodReconstructor>(
89  reconstructed_upper_side_of_face_vars,
90  reconstructed_lower_side_of_face_vars, volume_vars, ghost_cell_vars,
91  volume_extents, number_of_variables);
92 }
93 } // namespace fd::reconstruction
utility
fd::reconstruction
Variable and flux vector splitting reconstruction schemes for finite difference methods.
Definition: AoWeno.hpp:25
Index
Definition: Index.hpp:31
cmath
sgn
constexpr T sgn(const T &val) noexcept
Compute the sign function of val defined as 1 if val > 0, 0 if val == 0, and -1 if val < 0.
Definition: Math.hpp:163
Direction
Definition: Direction.hpp:23
SPECTRE_ALWAYS_INLINE
#define SPECTRE_ALWAYS_INLINE
Definition: ForceInline.hpp:16
cstddef
array
std::min
T min(T... args)
gsl::span
Create a span/view on a range, which is cheap to copy (one pointer).
Definition: Gsl.hpp:292
DirectionMap
Definition: DirectionMap.hpp:15
Gsl.hpp
ForceInline.hpp
fd::reconstruction::minmod
void minmod(const gsl::not_null< std::array< gsl::span< double >, Dim > * > reconstructed_upper_side_of_face_vars, const gsl::not_null< std::array< gsl::span< double >, Dim > * > reconstructed_lower_side_of_face_vars, const gsl::span< const double > &volume_vars, const DirectionMap< Dim, gsl::span< const double >> &ghost_cell_vars, const Index< Dim > &volume_extents, const size_t number_of_variables) noexcept
Performs minmod reconstruction on the volume_vars in each direction.
Definition: Minmod.hpp:80
TMPL.hpp
gsl::not_null
Require a pointer to not be a nullptr