Rusanov.hpp
1 // Distributed under the MIT License.
2 // See LICENSE.txt for details.
3 
4 #pragma once
5 
6 #include <memory>
7 #include <optional>
8 
11 #include "Evolution/Systems/RelativisticEuler/Valencia/BoundaryCorrections/BoundaryCorrection.hpp"
12 #include "Evolution/Systems/RelativisticEuler/Valencia/Tags.hpp"
13 #include "NumericalAlgorithms/DiscontinuousGalerkin/Formulation.hpp"
14 #include "Options/Options.hpp"
16 #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
17 #include "PointwiseFunctions/Hydro/EquationsOfState/EquationOfState.hpp"
18 #include "PointwiseFunctions/Hydro/Tags.hpp"
19 #include "Utilities/Gsl.hpp"
20 #include "Utilities/TMPL.hpp"
21 
22 /// \cond
23 class DataVector;
24 namespace gsl {
25 template <typename T>
26 class not_null;
27 } // namespace gsl
28 namespace PUP {
29 class er;
30 } // namespace PUP
31 /// \endcond
32 
34 /*!
35  * \brief A Rusanov/local Lax-Friedrichs Riemann solver
36  *
37  * Let \f$U\f$ be the state vector of evolved variables, \f$F^i\f$ the
38  * corresponding fluxes, and \f$n_i\f$ be the outward directed unit normal to
39  * the interface. Denoting \f$F := n_i F^i\f$, the %Rusanov boundary correction
40  * is
41  *
42  * \f{align*}
43  * G_\text{Rusanov} = \frac{F_\text{int} - F_\text{ext}}{2} -
44  * \frac{\text{max}\left(\{|\lambda_\text{int}|\},
45  * \{|\lambda_\text{ext}|\}\right)}{2} \left(U_\text{ext} - U_\text{int}\right),
46  * \f}
47  *
48  * where "int" and "ext" stand for interior and exterior, and
49  * \f$\{|\lambda|\}\f$ is the set of characteristic/signal speeds. The minus
50  * sign in front of the \f$F_{\text{ext}}\f$ is necessary because the outward
51  * directed normal of the neighboring element has the opposite sign, i.e.
52  * \f$n_i^{\text{ext}}=-n_i^{\text{int}}\f$. The characteristic/signal speeds
53  * expressions are listed in the documentation of `characteristic_speeds`.
54  *
55  * \note In the strong form the `dg_boundary_terms` function returns
56  * \f$G - F_\text{int}\f$
57  */
58 template <size_t Dim>
59 class Rusanov final : public BoundaryCorrection<Dim> {
60  private:
61  struct AbsCharSpeed : db::SimpleTag {
62  using type = Scalar<DataVector>;
63  };
64 
65  public:
66  using options = tmpl::list<>;
67  static constexpr Options::String help = {
68  "Computes the Rusanov or local Lax-Friedrichs boundary correction term "
69  "for the Valencia formulation of the relativistic Euler/hydrodynamics "
70  "system."};
71 
72  Rusanov() = default;
73  Rusanov(const Rusanov&) = default;
74  Rusanov& operator=(const Rusanov&) = default;
75  Rusanov(Rusanov&&) = default;
76  Rusanov& operator=(Rusanov&&) = default;
77  ~Rusanov() override = default;
78 
79  /// \cond
80  explicit Rusanov(CkMigrateMessage* msg) noexcept;
81  using PUP::able::register_constructor;
83  /// \endcond
84  void pup(PUP::er& p) override; // NOLINT
85 
86  std::unique_ptr<BoundaryCorrection<Dim>> get_clone() const noexcept override;
87 
88  using dg_package_field_tags =
89  tmpl::list<Tags::TildeD, Tags::TildeTau, Tags::TildeS<Dim>,
93  using dg_package_data_temporary_tags =
94  tmpl::list<gr::Tags::Lapse<DataVector>, gr::Tags::Shift<Dim>,
96  using dg_package_data_primitive_tags =
97  tmpl::list<hydro::Tags::RestMassDensity<DataVector>,
101  using dg_package_data_volume_tags =
102  tmpl::list<hydro::Tags::EquationOfStateBase>;
103 
104  template <size_t ThermodynamicDim>
105  double dg_package_data(
106  gsl::not_null<Scalar<DataVector>*> packaged_tilde_d,
107  gsl::not_null<Scalar<DataVector>*> packaged_tilde_tau,
108  gsl::not_null<tnsr::i<DataVector, Dim, Frame::Inertial>*>
109  packaged_tilde_s,
110  gsl::not_null<Scalar<DataVector>*> packaged_normal_dot_flux_tilde_d,
111  gsl::not_null<Scalar<DataVector>*> packaged_normal_dot_flux_tilde_tau,
112  gsl::not_null<tnsr::i<DataVector, Dim, Frame::Inertial>*>
113  packaged_normal_dot_flux_tilde_s,
114  gsl::not_null<Scalar<DataVector>*> packaged_abs_char_speed,
115 
116  const Scalar<DataVector>& tilde_d, const Scalar<DataVector>& tilde_tau,
117  const tnsr::i<DataVector, Dim, Frame::Inertial>& tilde_s,
118 
119  const tnsr::I<DataVector, Dim, Frame::Inertial>& flux_tilde_d,
120  const tnsr::I<DataVector, Dim, Frame::Inertial>& flux_tilde_tau,
121  const tnsr::Ij<DataVector, Dim, Frame::Inertial>& flux_tilde_s,
122 
123  const Scalar<DataVector>& lapse,
124  const tnsr::I<DataVector, Dim, Frame::Inertial>& shift,
125  const tnsr::ii<DataVector, Dim, Frame::Inertial>& spatial_metric,
126  const Scalar<DataVector>& rest_mass_density,
127  const Scalar<DataVector>& specific_internal_energy,
128  const Scalar<DataVector>& specific_enthalpy,
129  const tnsr::I<DataVector, Dim, Frame::Inertial>& spatial_velocity,
130 
131  const tnsr::i<DataVector, Dim, Frame::Inertial>& normal_covector,
132  const tnsr::I<DataVector, Dim, Frame::Inertial>& normal_vector,
133  const std::optional<tnsr::I<DataVector, Dim, Frame::Inertial>>&
134  /*mesh_velocity*/,
135  const std::optional<Scalar<DataVector>>& normal_dot_mesh_velocity,
137  equation_of_state) const noexcept;
138 
139  void dg_boundary_terms(
140  gsl::not_null<Scalar<DataVector>*> boundary_correction_tilde_d,
141  gsl::not_null<Scalar<DataVector>*> boundary_correction_tilde_tau,
142  gsl::not_null<tnsr::i<DataVector, Dim, Frame::Inertial>*>
143  boundary_correction_tilde_s,
144  const Scalar<DataVector>& tilde_d_int,
145  const Scalar<DataVector>& tilde_tau_int,
146  const tnsr::i<DataVector, Dim, Frame::Inertial>& tilde_s_int,
147  const Scalar<DataVector>& normal_dot_flux_tilde_d_int,
148  const Scalar<DataVector>& normal_dot_flux_tilde_tau_int,
149  const tnsr::i<DataVector, Dim, Frame::Inertial>&
150  normal_dot_flux_tilde_s_int,
151  const Scalar<DataVector>& abs_char_speed_int,
152  const Scalar<DataVector>& tilde_d_ext,
153  const Scalar<DataVector>& tilde_tau_ext,
154  const tnsr::i<DataVector, Dim, Frame::Inertial>& tilde_s_ext,
155  const Scalar<DataVector>& normal_dot_flux_tilde_d_ext,
156  const Scalar<DataVector>& normal_dot_flux_tilde_tau_ext,
157  const tnsr::i<DataVector, Dim, Frame::Inertial>&
158  normal_dot_flux_tilde_s_ext,
159  const Scalar<DataVector>& abs_char_speed_ext,
160  dg::Formulation dg_formulation) const noexcept;
161 };
162 } // namespace RelativisticEuler::Valencia::BoundaryCorrections
CharmPupable.hpp
hydro::Tags::SpecificInternalEnergy
The specific internal energy .
Definition: Tags.hpp:173
gr::Tags::SpatialMetric
Definition: Tags.hpp:26
EquationsOfState::EquationOfState
Base class for equations of state depending on whether or not the system is relativistic,...
Definition: EquationOfState.hpp:63
Options.hpp
RelativisticEuler::Valencia::BoundaryCorrections
Boundary corrections/numerical fluxes.
Definition: BoundaryCorrection.hpp:13
db::SimpleTag
Mark a struct as a simple tag by inheriting from this.
Definition: Tag.hpp:36
Tags::NormalDotFlux
Prefix indicating a boundary unit normal vector dotted into the flux.
Definition: Prefixes.hpp:96
gr::lapse
Scalar< DataType > lapse(const tnsr::I< DataType, SpatialDim, Frame > &shift, const tnsr::aa< DataType, SpatialDim, Frame > &spacetime_metric) noexcept
Compute lapse from shift and spacetime metric.
RelativisticEuler::Valencia::BoundaryCorrections::Rusanov
A Rusanov/local Lax-Friedrichs Riemann solver.
Definition: Rusanov.hpp:59
WRAPPED_PUPable_decl_template
#define WRAPPED_PUPable_decl_template(className)
Mark derived classes as serializable.
Definition: CharmPupable.hpp:22
hydro::Tags::SpecificEnthalpy
The relativistic specific enthalpy .
Definition: Tags.hpp:167
DataVector
Stores a collection of function values.
Definition: DataVector.hpp:46
memory
gr::Tags::Shift
Definition: Tags.hpp:48
RelativisticEuler::Valencia::BoundaryCorrections::BoundaryCorrection
The base class used to create boundary corrections from input files and store them in the global cach...
Definition: BoundaryCorrection.hpp:24
gr::shift
tnsr::I< DataType, SpatialDim, Frame > shift(const tnsr::aa< DataType, SpatialDim, Frame > &spacetime_metric, const tnsr::II< DataType, SpatialDim, Frame > &inverse_spatial_metric) noexcept
Compute shift from spacetime metric and inverse spatial metric.
Scalar
Tensor< T, Symmetry<>, index_list<> > Scalar
Definition: TypeAliases.hpp:21
Gsl.hpp
TypeAliases.hpp
Options::String
const char *const String
The string used in option structs.
Definition: Options.hpp:32
hydro::Tags::SpatialVelocity
The spatial velocity of the fluid, where . Here is the spatial part of the 4-velocity of the fluid,...
Definition: Tags.hpp:142
gr::spatial_metric
tnsr::ii< DataType, SpatialDim, Frame > spatial_metric(const tnsr::aa< DataType, SpatialDim, Frame > &spacetime_metric) noexcept
Compute spatial metric from spacetime metric.
optional
dg::Formulation
Formulation
The DG formulation to use.
Definition: Formulation.hpp:34
std::unique_ptr
Prefixes.hpp
gsl
Implementations from the Guideline Support Library.
Definition: ReadSpecPiecewisePolynomial.hpp:11
TMPL.hpp
gsl::not_null
Require a pointer to not be a nullptr
Definition: ReadSpecPiecewisePolynomial.hpp:13