Hllc.hpp
1 // Distributed under the MIT License.
2 // See LICENSE.txt for details.
3 
4 #pragma once
5 
6 #include <cstddef>
7 #include <string>
8 
10 #include "Domain/FaceNormal.hpp"
11 #include "Evolution/Systems/NewtonianEuler/Characteristics.hpp"
12 #include "Evolution/Systems/NewtonianEuler/TagsDeclarations.hpp"
13 #include "Options/Options.hpp"
14 #include "Utilities/Gsl.hpp"
15 #include "Utilities/TMPL.hpp"
16 
17 /// \cond
18 namespace PUP {
19 class er;
20 } // namespace PUP
21 
22 namespace Tags {
23 template <typename>
24 struct NormalDotFlux;
25 template <typename>
26 struct Normalized;
27 } // namespace Tags
28 
29 class DataVector;
30 template <typename>
31 class Variables;
32 /// \endcond
33 
34 namespace NewtonianEuler {
35 namespace NumericalFluxes {
36 
37 /*!
38  * \brief Compute the HLLC numerical flux
39  *
40  * Class implementing the HLLC flux for the Newtonian Euler equations,
41  * originally introduced by E. F. Toro, M. Spruce and W. Speares \cite Toro1994.
42  * Let \f$F^k = [F^k(\rho), F^k(S^i), F^k(e)]^T\f$ be the fluxes of
43  * the conservative quantities \f$U = [\rho, S^i, e]^T\f$
44  * (see NewtonianEuler::ComputeFluxes). Denoting \f$v_n = n_k v^k\f$ and
45  * \f$F = n_kF^k\f$, where \f$v^k\f$ is the velocity and
46  * \f$n_k\f$ is the interface unit normal, the HLLC flux is
47  *
48  * \f{align*}
49  * G_\text{HLLC} =
50  * \begin{cases}
51  * F_\text{int}, & S_\text{min} > 0, \\
52  * G_{*\text{int}}, & S_\text{min} \leq 0\quad\text{and}\quad S_* > 0, \\
53  * G_{*\text{ext}}, & S_* \leq 0\quad\text{and}\quad S_\text{max} > 0, \\
54  * F_\text{ext}, & S_\text{max} \leq 0,
55  * \end{cases}
56  * \f}
57  *
58  * where
59  *
60  * \f{align*}
61  * G_{*\text{int}} &= \frac{S_*\left(F_\text{int}
62  * - S_\text{min}U_\text{int}\right)
63  * - S_\text{min} p_{*\text{int}}D_*}{S_* - S_\text{min}},\\
64  * G_{*\text{ext}} &= \frac{S_*\left(F_\text{ext}
65  * - S_\text{max}U_\text{ext}\right)
66  * - S_\text{max} p_{*\text{ext}}D_*}{S_* - S_\text{max}},
67  * \f}
68  *
69  * with
70  *
71  * \f{align*}
72  * p_{*\text{int}} &\equiv p_\text{int} + \rho_\text{int}
73  * \left[(v_n)_\text{int} - S_\text{min}\right]
74  * \left[(v_n)_\text{int} - S_*\right], \\
75  * p_{*\text{ext}} &\equiv p_\text{ext} + \rho_\text{ext}
76  * \left[(v_n)_\text{ext} - S_\text{max}\right]
77  * \left[(v_n)_\text{ext} - S_*\right], \\
78  * S_* &\equiv \frac{p_\text{int} + F(\rho)_\text{int}\left[(v_n)_\text{int}
79  * - S_\text{min}\right] - p_\text{ext}
80  * - F(\rho)_\text{ext}\left[(v_n)_\text{ext} -
81  * S_\text{max}\right]}{F(\rho)_\text{int} - \rho_\text{int}S_\text{min}
82  * - F(\rho)_\text{ext} + \rho_\text{ext}S_\text{max}},\\
83  * D_* &\equiv \left[\begin{array}{c}
84  * 0\\ n^i\\ S_*
85  * \end{array}\right],
86  * \f}
87  *
88  * and \f$S_\text{min}\f$ and \f$S_\text{max}\f$ are estimates of the minimum
89  * and maximum signal speeds bounding the interior-moving and exterior-moving
90  * wavespeeds that arise when solving the Riemann problem. One requires
91  * \f$S_\text{min} \leq S_* \leq S_\text{max}\f$. As estimates, we use
92  *
93  * \f{align*}
94  * S_\text{min} &=
95  * \text{min}\left(\{\lambda_\text{int}\},\{\lambda_\text{ext}\}, 0\right)\\
96  * S_\text{max} &=
97  * \text{max}\left(\{\lambda_\text{int}\},\{\lambda_\text{ext}\}, 0\right),
98  * \f}
99  *
100  * where \f$\{\lambda\}\f$ is the set of all the characteristic speeds along a
101  * given normal. This way, the definition of \f$G_\text{HLLC}\f$ simplifies to
102  *
103  * \f{align*}
104  * G_\text{HLLC} =
105  * \begin{cases}
106  * \dfrac{S_*\left(F_\text{int} - S_\text{min}U_\text{int}\right)
107  * - S_\text{min} p_{*\text{int}}D_*}{S_* - S_\text{min}}, & S_* > 0, \\
108  * \dfrac{S_*\left(F_\text{ext} - S_\text{max}U_\text{ext}\right)
109  * - S_\text{max} p_{*\text{ext}}D_*}{S_* - S_\text{max}}, & S_* \leq 0. \\
110  * \end{cases}
111  * \f}
112  *
113  * \warning The HLLC flux implemented here does not incorporate any cure for
114  * the Carbuncle phenomenon or other shock instabilities reported in the
115  * literature. Prefer using another numerical flux in more than 1-d.
116  */
117 template <size_t Dim, typename Frame>
118 struct Hllc {
120 
121  /// Estimate for one of the signal speeds
123  using type = Scalar<DataVector>;
124  static std::string name() noexcept { return "FirstSpeedEstimate"; }
125  };
126 
127  /// Estimate for the other signal speed
129  using type = Scalar<DataVector>;
130  static std::string name() noexcept { return "SecondSpeedEstimate"; }
131  };
132 
133  /// The normal component of the velocity
135  using type = Scalar<DataVector>;
136  static std::string name() noexcept { return "NormalVelocity"; }
137  };
138 
139  using options = tmpl::list<>;
140  static constexpr OptionString help = {
141  "Compute the HLLC flux for the Newtonian Euler system."};
142 
143  // clang-tidy: non-const reference
144  void pup(PUP::er& /*p*/) noexcept {} // NOLINT
145 
146  using package_tags = tmpl::list<
155 
156  using argument_tags = tmpl::list<
157  ::Tags::NormalDotFlux<Tags::MassDensityCons<DataVector>>,
158  ::Tags::NormalDotFlux<Tags::MomentumDensity<DataVector, Dim, Frame>>,
159  ::Tags::NormalDotFlux<Tags::EnergyDensity<DataVector>>,
160  Tags::MassDensityCons<DataVector>,
161  Tags::MomentumDensity<DataVector, Dim, Frame>,
162  Tags::EnergyDensity<DataVector>, Tags::Velocity<DataVector, Dim, Frame>,
163  Tags::Pressure<DataVector>, char_speeds_tag,
164  ::Tags::Normalized<::Tags::UnnormalizedFaceNormal<Dim, Frame>>>;
165 
166  void package_data(
167  gsl::not_null<Variables<package_tags>*> packaged_data,
168  const Scalar<DataVector>& normal_dot_flux_mass_density,
169  const tnsr::I<DataVector, Dim, Frame>& normal_dot_flux_momentum_density,
170  const Scalar<DataVector>& normal_dot_flux_energy_density,
171  const Scalar<DataVector>& mass_density,
172  const tnsr::I<DataVector, Dim, Frame>& momentum_density,
173  const Scalar<DataVector>& energy_density,
174  const tnsr::I<DataVector, Dim, Frame>& velocity,
175  const Scalar<DataVector>& pressure,
176  const db::const_item_type<char_speeds_tag>& characteristic_speeds,
177  const tnsr::i<DataVector, Dim, Frame>& interface_unit_normal) const
178  noexcept;
179 
180  void operator()(
181  gsl::not_null<Scalar<DataVector>*> normal_dot_numerical_flux_mass_density,
182  gsl::not_null<tnsr::I<DataVector, Dim, Frame>*>
183  normal_dot_numerical_flux_momentum_density,
185  normal_dot_numerical_flux_energy_density,
186  const Scalar<DataVector>& normal_dot_flux_mass_density_int,
187  const tnsr::I<DataVector, Dim, Frame>&
188  normal_dot_flux_momentum_density_int,
189  const Scalar<DataVector>& normal_dot_flux_energy_density_int,
190  const Scalar<DataVector>& mass_density_int,
191  const tnsr::I<DataVector, Dim, Frame>& momentum_density_int,
192  const Scalar<DataVector>& energy_density_int,
193  const Scalar<DataVector>& pressure_int,
194  const tnsr::i<DataVector, Dim, Frame>& interface_unit_normal,
195  const Scalar<DataVector>& normal_velocity_int,
196  const Scalar<DataVector>& min_signal_speed_int,
197  const Scalar<DataVector>& max_signal_speed_int,
198  const Scalar<DataVector>& minus_normal_dot_flux_mass_density_ext,
199  const tnsr::I<DataVector, Dim, Frame>&
200  minus_normal_dot_flux_momentum_density_ext,
201  const Scalar<DataVector>& minus_normal_dot_flux_energy_density_ext,
202  const Scalar<DataVector>& mass_density_ext,
203  const tnsr::I<DataVector, Dim, Frame>& momentum_density_ext,
204  const Scalar<DataVector>& energy_density_ext,
205  const Scalar<DataVector>& pressure_ext,
206  const tnsr::i<DataVector, Dim, Frame>& minus_interface_unit_normal,
207  const Scalar<DataVector>& minus_normal_velocity_ext,
208  // names are inverted w.r.t interior data. See package_data()
209  const Scalar<DataVector>& minus_max_signal_speed_ext,
210  const Scalar<DataVector>& minus_min_signal_speed_ext) const noexcept;
211 };
212 
213 } // namespace NumericalFluxes
214 } // namespace NewtonianEuler
Prefix indicating a boundary unit normal vector dotted into the flux.
Definition: Prefixes.hpp:122
The momentum density of the fluid.
Definition: Tags.hpp:43
Definition: Strahlkorper.hpp:14
The normalized (co)vector represented by Tag.
Definition: Magnitude.hpp:119
The macroscopic or flow velocity of the fluid.
Definition: Tags.hpp:59
Defines classes and functions for making classes creatable from input files.
Tags for the DataBox inherit from this type.
Definition: DataBoxTag.hpp:64
const char *const OptionString
The string used in option structs.
Definition: Options.hpp:29
The energy density of the fluid.
Definition: Tags.hpp:52
Estimate for one of the signal speeds.
Definition: Hllc.hpp:122
Estimate for the other signal speed.
Definition: Hllc.hpp:128
Definition: DataBoxTag.hpp:29
The fluid pressure.
Definition: Tags.hpp:75
Defines a list of useful type aliases for tensors.
Declares function unnormalized_face_normal.
Definition: Characteristics.hpp:120
Compute the HLLC numerical flux.
Definition: Hllc.hpp:118
Stores a collection of function values.
Definition: DataVector.hpp:42
Wraps the template metaprogramming library used (brigand)
Defines functions and classes from the GSL.
Tensor< T, Symmetry<>, index_list<> > Scalar
Scalar type.
Definition: TypeAliases.hpp:21
The mass density of the fluid (as a conservative variable).
Definition: Tags.hpp:36
The normal component of the velocity.
Definition: Hllc.hpp:134
Require a pointer to not be a nullptr
Definition: Gsl.hpp:182
Items related to evolving the Newtonian Euler system.
Definition: EvolveNewtonianEulerFwd.hpp:8