Weno.hpp
1 // Distributed under the MIT License.
2 // See LICENSE.txt for details.
3 
4 #pragma once
5 
6 #include <array>
7 #include <boost/functional/hash.hpp> // IWYU pragma: keep
8 #include <cstddef>
9 #include <optional>
10 #include <string>
11 #include <unordered_map>
12 #include <utility>
13 
15 #include "Domain/SizeOfElement.hpp"
16 #include "Domain/Tags.hpp" // IWYU pragma: keep
17 #include "Evolution/DiscontinuousGalerkin/Limiters/Tags.hpp"
18 #include "Evolution/DiscontinuousGalerkin/Limiters/Weno.hpp"
19 #include "Evolution/DiscontinuousGalerkin/Limiters/WenoType.hpp"
20 #include "Evolution/DiscontinuousGalerkin/NormalVectorTags.hpp"
21 #include "Evolution/Systems/NewtonianEuler/Limiters/VariablesToLimit.hpp"
22 #include "Evolution/Systems/NewtonianEuler/Tags.hpp"
23 #include "Options/Auto.hpp"
24 #include "Options/Options.hpp"
25 #include "PointwiseFunctions/Hydro/TagsDeclarations.hpp"
26 #include "Utilities/TMPL.hpp"
27 
28 /// \cond
29 class DataVector;
30 template <size_t VolumeDim>
31 class Direction;
32 template <size_t VolumeDim>
33 class Element;
34 template <size_t VolumeDim>
35 class ElementId;
36 template <size_t VolumeDim>
37 class Mesh;
38 template <size_t VolumeDim>
39 class OrientationMap;
40 
41 namespace EquationsOfState {
42 template <bool IsRelativistic, size_t ThermodynamicDim>
43 class EquationOfState;
44 } // namespace EquationsOfState
45 
46 namespace gsl {
47 template <typename T>
48 class not_null;
49 } // namespace gsl
50 
51 namespace PUP {
52 class er;
53 } // namespace PUP
54 /// \endcond
55 
56 namespace NewtonianEuler {
57 namespace Limiters {
58 
59 /// \ingroup LimitersGroup
60 /// \brief A compact-stencil WENO limiter for the NewtonianEuler system.
61 ///
62 /// Implements the simple WENO limiter of \cite Zhong2013 and the Hermite WENO
63 /// (HWENO) limiter of \cite Zhu2016. See the documentation of the
64 /// system-agnostic ::Limiters::Weno limiter for a general discussion of the
65 /// algorithm and the various options that control the action of the limiter.
66 //
67 /// This implemention is specialized to the NewtonianEuler evolution system.
68 /// By specializing the limiter to the system, we can add a few features that
69 /// improve its robustness:
70 /// - the troubled-cell indicator (TCI) can be specialized to the features of
71 /// the evolution system.
72 /// - the limiter can be applied to the system's characteristic variables. This
73 /// is the recommendation of the reference, because it reduces spurious
74 /// oscillations in the post-limiter solution.
75 /// - after limiting, the solution can be processed to remove any remaining
76 /// unphysical values like negative densities and pressures. We do this by
77 /// scaling the solution around its mean (a "flattener" or "bounds-preserving"
78 /// filter). Note: the flattener is applied to all elements, including those
79 /// where the limiter did not act to reduce the solution's slopes.
80 ///
81 /// The matrix of TCI, variables to limit, post-processing, etc. choices can
82 /// rapidly grow large. Here we reduce the possibilities by tying the TCI to the
83 /// limiter in keeping with each limiter's main reference: HWENO uses the KXRCF
84 /// TCI and simple WENO uses the TVB TCI. To fully explore the matrix of
85 /// possibilities, the source code could be generalized --- however, experience
86 /// suggests it is unlikely that there exists one combination that will perform
87 /// remarkably better than the others.
88 template <size_t VolumeDim>
89 class Weno {
90  public:
92  VolumeDim, tmpl::list<NewtonianEuler::Tags::MassDensityCons,
95 
98  static type suggested_value() noexcept {
99  return NewtonianEuler::Limiters::VariablesToLimit::Characteristic;
100  }
101  static constexpr Options::String help = {
102  "Variable representation on which to apply the limiter"};
103  };
104  // Future design improvement: attach the TvbConstant/KxrcfConstant to the
105  // limiter type, so that it isn't necessary to specify both (but with one
106  // required to be 'None') in each input file.
107  struct TvbConstant {
109  static constexpr Options::String help = {
110  "Constant in RHS of the TVB minmod TCI, used when Type = SimpleWeno"};
111  };
112  struct KxrcfConstant {
114  static constexpr Options::String help = {
115  "Constant in RHS of KXRCF TCI, used when Type = Hweno"};
116  };
117  struct ApplyFlattener {
118  using type = bool;
119  static constexpr Options::String help = {
120  "Flatten after limiting to restore pointwise positivity"};
121  };
122  using options =
123  tmpl::list<typename ConservativeVarsWeno::Type, VariablesToLimit,
124  typename ConservativeVarsWeno::NeighborWeight, TvbConstant,
126  typename ConservativeVarsWeno::DisableForDebugging>;
127  static constexpr Options::String help = {
128  "A WENO limiter specialized to the NewtonianEuler system"};
129  static std::string name() noexcept { return "NewtonianEulerWeno"; };
130 
131  Weno(::Limiters::WenoType weno_type,
133  double neighbor_linear_weight, std::optional<double> tvb_constant,
134  std::optional<double> kxrcf_constant, bool apply_flattener,
135  bool disable_for_debugging = false,
136  const Options::Context& context = {});
137 
138  Weno() noexcept = default;
139  Weno(const Weno& /*rhs*/) = default;
140  Weno& operator=(const Weno& /*rhs*/) = default;
141  Weno(Weno&& /*rhs*/) noexcept = default;
142  Weno& operator=(Weno&& /*rhs*/) noexcept = default;
143  ~Weno() = default;
144 
145  // NOLINTNEXTLINE(google-runtime-references)
146  void pup(PUP::er& p) noexcept;
147 
148  using PackagedData = typename ConservativeVarsWeno::PackagedData;
149  using package_argument_tags =
150  typename ConservativeVarsWeno::package_argument_tags;
151 
152  /// \brief Package data for sending to neighbor elements
153  void package_data(
154  gsl::not_null<PackagedData*> packaged_data,
155  const Scalar<DataVector>& mass_density_cons,
156  const tnsr::I<DataVector, VolumeDim>& momentum_density,
157  const Scalar<DataVector>& energy_density, const Mesh<VolumeDim>& mesh,
158  const std::array<double, VolumeDim>& element_size,
159  const OrientationMap<VolumeDim>& orientation_map) const noexcept;
160 
161  using limit_tags =
162  tmpl::list<NewtonianEuler::Tags::MassDensityCons,
163  NewtonianEuler::Tags::MomentumDensity<VolumeDim>,
164  NewtonianEuler::Tags::EnergyDensity>;
165  using limit_argument_tags =
166  tmpl::list<domain::Tags::Mesh<VolumeDim>,
167  domain::Tags::Element<VolumeDim>,
168  domain::Tags::SizeOfElement<VolumeDim>,
169  domain::Tags::DetInvJacobian<Frame::Logical, Frame::Inertial>,
170  evolution::dg::Tags::NormalCovectorAndMagnitude<VolumeDim>,
171  ::hydro::Tags::EquationOfStateBase>;
172 
173  /// \brief Limit the solution on the element
174  template <size_t ThermodynamicDim>
175  bool operator()(
176  gsl::not_null<Scalar<DataVector>*> mass_density_cons,
177  gsl::not_null<tnsr::I<DataVector, VolumeDim>*> momentum_density,
178  gsl::not_null<Scalar<DataVector>*> energy_density,
179  const Mesh<VolumeDim>& mesh, const Element<VolumeDim>& element,
180  const std::array<double, VolumeDim>& element_size,
181  const Scalar<DataVector>& det_inv_logical_to_inertial_jacobian,
182  const typename evolution::dg::Tags::NormalCovectorAndMagnitude<
183  VolumeDim>::type& normals_and_magnitudes,
184  const EquationsOfState::EquationOfState<false, ThermodynamicDim>&
185  equation_of_state,
186  const std::unordered_map<
187  std::pair<Direction<VolumeDim>, ElementId<VolumeDim>>, PackagedData,
188  boost::hash<std::pair<Direction<VolumeDim>, ElementId<VolumeDim>>>>&
189  neighbor_data) const noexcept;
190 
191  private:
192  template <size_t LocalDim>
193  // NOLINTNEXTLINE(readability-redundant-declaration) false positive
194  friend bool operator==(const Weno<LocalDim>& lhs,
195  const Weno<LocalDim>& rhs) noexcept;
196 
197  ::Limiters::WenoType weno_type_;
198  NewtonianEuler::Limiters::VariablesToLimit vars_to_limit_;
199  double neighbor_linear_weight_;
200  std::optional<double> tvb_constant_;
201  std::optional<double> kxrcf_constant_;
202  bool apply_flattener_;
203  bool disable_for_debugging_;
204  // Note: conservative_vars_weno_ is always used for calls to package_data, and
205  // is also used when limiting cons vars with the simple WENO algorithm. So we
206  // construct conservative_vars_weno_ with the correct TVB constant when it is
207  // used for limiting (precisely, when weno_type_ == SimpleWeno), and with a
208  // dummy TVB constant value of NaN otherwise (weno_type_ == Hweno). This lets
209  // us construct the conservative_vars_weno_ variable and use it for delegating
210  // the package_data work even when the specialized limiter has no TVB constant
211  // and is possible because package_data doesn't depend on the TCI.
212  ConservativeVarsWeno conservative_vars_weno_;
213 };
214 
215 template <size_t VolumeDim>
216 bool operator!=(const Weno<VolumeDim>& lhs,
217  const Weno<VolumeDim>& rhs) noexcept;
218 
219 } // namespace Limiters
220 } // namespace NewtonianEuler
EquationsOfState
Contains all equations of state, including base class.
Definition: DarkEnergyFluid.hpp:26
Limiters
Things relating to limiting.
Definition: HwenoImpl.hpp:42
std::string
utility
NewtonianEuler::Limiters::Weno::VariablesToLimit
Definition: Weno.hpp:96
evolution
Functionality for evolving hyperbolic partial differential equations.
Definition: RunEventsAndDenseTriggers.hpp:30
NewtonianEuler::Limiters::Weno::package_data
void package_data(gsl::not_null< PackagedData * > packaged_data, const Scalar< DataVector > &mass_density_cons, const tnsr::I< DataVector, VolumeDim > &momentum_density, const Scalar< DataVector > &energy_density, const Mesh< VolumeDim > &mesh, const std::array< double, VolumeDim > &element_size, const OrientationMap< VolumeDim > &orientation_map) const noexcept
Package data for sending to neighbor elements.
Options.hpp
Tags.hpp
NewtonianEuler::Limiters::Weno::TvbConstant
Definition: Weno.hpp:107
NewtonianEuler
Items related to evolving the Newtonian Euler system.
Definition: EvolveNewtonianEulerFwd.hpp:8
OrientationMap
A mapping of the logical coordinate axes of a host to the logical coordinate axes of a neighbor of th...
Definition: OrientationMap.hpp:34
NewtonianEuler::Tags::MassDensityCons
The mass density of the fluid (as a conservative variable).
Definition: Tags.hpp:31
Options::Context
Definition: Options.hpp:41
Direction
Definition: Direction.hpp:23
NewtonianEuler::Limiters::Weno::ApplyFlattener
Definition: Weno.hpp:117
dg
Functionality related to discontinuous Galerkin schemes.
Definition: ApplyMassMatrix.hpp:13
Element
Definition: Element.hpp:29
ElementId
An ElementId uniquely labels an Element.
Definition: ElementId.hpp:51
NewtonianEuler::Limiters::VariablesToLimit
VariablesToLimit
Type of NewtonianEuler variables to apply limiter to.
Definition: VariablesToLimit.hpp:20
NewtonianEuler::Limiters::Weno
A compact-stencil WENO limiter for the NewtonianEuler system.
Definition: Weno.hpp:89
cstddef
array
DataVector
Stores a collection of function values.
Definition: DataVector.hpp:46
hydro
Items related to hydrodynamic systems.
Definition: SmoothFlow.hpp:19
NewtonianEuler::Tags::EnergyDensity
The energy density of the fluid.
Definition: Tags.hpp:45
Options::Auto
A class indicating that a parsed value can be automatically computed instead of specified.
Definition: Auto.hpp:36
Mesh
Holds the number of grid points, basis, and quadrature in each direction of the computational grid.
Definition: Mesh.hpp:49
tnsr
Type aliases to construct common Tensors.
Definition: TypeAliases.hpp:31
NewtonianEuler::Tags::MomentumDensity
The momentum density of the fluid.
Definition: Tags.hpp:37
Scalar
Tensor< T, Symmetry<>, index_list<> > Scalar
Definition: TypeAliases.hpp:21
TypeAliases.hpp
Options::String
const char *const String
The string used in option structs.
Definition: Options.hpp:32
Limiters::WenoType
WenoType
Possible types of the WENO limiter.
Definition: WenoType.hpp:21
Frame
Definition: IndexType.hpp:36
NewtonianEuler::Limiters::Weno::KxrcfConstant
Definition: Weno.hpp:112
optional
unordered_map
gsl
Implementations from the Guideline Support Library.
Definition: ReadSpecPiecewisePolynomial.hpp:11
TMPL.hpp
string