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/RadiationTransport/M1Grey/BoundaryCorrections/BoundaryCorrection.hpp"
12 #include "Evolution/Systems/RadiationTransport/M1Grey/Tags.hpp"
13 #include "NumericalAlgorithms/DiscontinuousGalerkin/Formulation.hpp"
14 #include "Options/Options.hpp"
16 #include "PointwiseFunctions/Hydro/EquationsOfState/EquationOfState.hpp"
17 #include "PointwiseFunctions/Hydro/Tags.hpp"
18 #include "Utilities/Gsl.hpp"
19 #include "Utilities/TMPL.hpp"
20 
21 /// \cond
22 class DataVector;
23 namespace gsl {
24 template <typename T>
25 class not_null;
26 } // namespace gsl
27 namespace PUP {
28 class er;
29 } // namespace PUP
30 /// \endcond
31 
33 
34 namespace Rusanov_detail {
35 void dg_package_data_impl(
36  gsl::not_null<Scalar<DataVector>*> packaged_tilde_e,
37  gsl::not_null<tnsr::i<DataVector, 3, Frame::Inertial>*> packaged_tilde_s,
38  gsl::not_null<Scalar<DataVector>*> packaged_normal_dot_flux_tilde_e,
39  gsl::not_null<tnsr::i<DataVector, 3, Frame::Inertial>*>
40  packaged_normal_dot_flux_tilde_s,
41  const Scalar<DataVector>& tilde_e,
42  const tnsr::i<DataVector, 3, Frame::Inertial>& tilde_s,
43  const tnsr::I<DataVector, 3, Frame::Inertial>& flux_tilde_e,
44  const tnsr::Ij<DataVector, 3, Frame::Inertial>& flux_tilde_s,
45  const tnsr::i<DataVector, 3, Frame::Inertial>& normal_covector,
46  const tnsr::I<DataVector, 3, Frame::Inertial>& normal_vector,
47  const std::optional<tnsr::I<DataVector, 3, Frame::Inertial>>& mesh_velocity,
48  const std::optional<Scalar<DataVector>>& normal_dot_mesh_velocity) noexcept;
49 
50 void dg_boundary_terms_impl(
51  gsl::not_null<Scalar<DataVector>*> boundary_correction_tilde_e,
52  gsl::not_null<tnsr::i<DataVector, 3, Frame::Inertial>*>
53  boundary_correction_tilde_s,
54  const Scalar<DataVector>& tilde_e_int,
55  const tnsr::i<DataVector, 3, Frame::Inertial>& tilde_s_int,
56  const Scalar<DataVector>& normal_dot_flux_tilde_e_int,
57  const tnsr::i<DataVector, 3, Frame::Inertial>& normal_dot_flux_tilde_s_int,
58  const Scalar<DataVector>& tilde_e_ext,
59  const tnsr::i<DataVector, 3, Frame::Inertial>& tilde_s_ext,
60  const Scalar<DataVector>& normal_dot_flux_tilde_e_ext,
61  const tnsr::i<DataVector, 3, Frame::Inertial>& normal_dot_flux_tilde_s_ext,
62  dg::Formulation dg_formulation) noexcept;
63 } // namespace Rusanov_detail
64 
65 /// \cond
66 template <typename NeutrinoSpeciesList>
67 class Rusanov;
68 /// \endcond
69 
70 /*!
71  * \brief A Rusanov/local Lax-Friedrichs Riemann solver
72  *
73  * Let \f$U\f$ be the state vector of evolved variables, \f$F^i\f$ the
74  * corresponding fluxes, and \f$n_i\f$ be the outward directed unit normal to
75  * the interface. Denoting \f$F := n_i F^i\f$, the %Rusanov boundary correction
76  * is
77  *
78  * \f{align*}
79  * G_\text{Rusanov} = \frac{F_\text{int} - F_\text{ext}}{2} -
80  * \frac{\text{max}\left(\{|\lambda_\text{int}|\},
81  * \{|\lambda_\text{ext}|\}\right)}{2} \left(U_\text{ext} - U_\text{int}\right),
82  * \f}
83  *
84  * where "int" and "ext" stand for interior and exterior, and
85  * \f$\{|\lambda|\}\f$ is the set of characteristic/signal speeds. The minus
86  * sign in front of the \f$F_{\text{ext}}\f$ is necessary because the outward
87  * directed normal of the neighboring element has the opposite sign, i.e.
88  * \f$n_i^{\text{ext}}=-n_i^{\text{int}}\f$.
89  *
90  * For radiation, the max characteristic/signal speed is 1 (the speed of light).
91  * Note that the characteristic speeds of this system are *not* yet fully
92  * implemented, see the warning in the documentation for the characteristics.
93  *
94  * \note In the strong form the `dg_boundary_terms` function returns
95  * \f$G - F_\text{int}\f$
96  */
97 template <typename... NeutrinoSpecies>
98 class Rusanov<tmpl::list<NeutrinoSpecies...>> final
99  : public BoundaryCorrection<tmpl::list<NeutrinoSpecies...>> {
100  public:
101  using options = tmpl::list<>;
102  static constexpr Options::String help = {
103  "Computes the Rusanov or local Lax-Friedrichs boundary correction term "
104  "for the M1Grey radiation transport system."};
105 
106  Rusanov() = default;
107  Rusanov(const Rusanov&) = default;
108  Rusanov& operator=(const Rusanov&) = default;
109  Rusanov(Rusanov&&) = default;
110  Rusanov& operator=(Rusanov&&) = default;
111  ~Rusanov() override = default;
112 
113  /// \cond
114  explicit Rusanov(CkMigrateMessage* msg) noexcept
115  : BoundaryCorrection<tmpl::list<NeutrinoSpecies...>>(msg) {}
116  using PUP::able::register_constructor;
117  WRAPPED_PUPable_decl_template(Rusanov); // NOLINT
118  /// \endcond
119  void pup(PUP::er& p) override { // NOLINT
120  BoundaryCorrection<tmpl::list<NeutrinoSpecies...>>::pup(p);
121  }
122 
123  std::unique_ptr<BoundaryCorrection<tmpl::list<NeutrinoSpecies...>>>
124  get_clone() const noexcept override {
125  return std::make_unique<Rusanov>(*this);
126  }
127 
128  using dg_package_field_tags = tmpl::list<
133  using dg_package_data_temporary_tags = tmpl::list<>;
134  using dg_package_data_primitive_tags = tmpl::list<>;
135  using dg_package_data_volume_tags = tmpl::list<>;
136 
137  double dg_package_data(
138  const gsl::not_null<typename Tags::TildeE<
139  Frame::Inertial, NeutrinoSpecies>::type*>... packaged_tilde_e,
140  const gsl::not_null<typename Tags::TildeS<
141  Frame::Inertial, NeutrinoSpecies>::type*>... packaged_tilde_s,
142  const gsl::not_null<typename Tags::TildeE<
144  NeutrinoSpecies>::type*>... packaged_normal_dot_flux_tilde_e,
145  const gsl::not_null<typename Tags::TildeS<
147  NeutrinoSpecies>::type*>... packaged_normal_dot_flux_tilde_s,
148 
149  const typename Tags::TildeE<Frame::Inertial,
150  NeutrinoSpecies>::type&... tilde_e,
151  const typename Tags::TildeS<Frame::Inertial,
152  NeutrinoSpecies>::type&... tilde_s,
153  const typename ::Tags::Flux<
155  Frame::Inertial>::type&... flux_tilde_e,
156  const typename ::Tags::Flux<
158  Frame::Inertial>::type&... flux_tilde_s,
159 
160  const tnsr::i<DataVector, 3, Frame::Inertial>& normal_covector,
161  const tnsr::I<DataVector, 3, Frame::Inertial>& normal_vector,
162  const std::optional<tnsr::I<DataVector, 3, Frame::Inertial>>&
163  mesh_velocity,
164  const std::optional<Scalar<DataVector>>& normal_dot_mesh_velocity)
165  const noexcept {
166  EXPAND_PACK_LEFT_TO_RIGHT(Rusanov_detail::dg_package_data_impl(
167  packaged_tilde_e, packaged_tilde_s, packaged_normal_dot_flux_tilde_e,
168  packaged_normal_dot_flux_tilde_s, tilde_e, tilde_s, flux_tilde_e,
169  flux_tilde_s, normal_covector, normal_vector, mesh_velocity,
170  normal_dot_mesh_velocity));
171  // max speed = speed of light
172  return 1.0;
173  }
174 
175  void dg_boundary_terms(
176  const gsl::not_null<typename Tags::TildeE<
178  NeutrinoSpecies>::type*>... boundary_correction_tilde_e,
179  const gsl::not_null<typename Tags::TildeS<
181  NeutrinoSpecies>::type*>... boundary_correction_tilde_s,
182  const typename Tags::TildeE<Frame::Inertial,
183  NeutrinoSpecies>::type&... tilde_e_int,
184  const typename Tags::TildeS<Frame::Inertial,
185  NeutrinoSpecies>::type&... tilde_s_int,
187  type&... normal_dot_flux_tilde_e_int,
189  type&... normal_dot_flux_tilde_s_int,
190  const typename Tags::TildeE<Frame::Inertial,
191  NeutrinoSpecies>::type&... tilde_e_ext,
192  const typename Tags::TildeS<Frame::Inertial,
193  NeutrinoSpecies>::type&... tilde_s_ext,
195  type&... normal_dot_flux_tilde_e_ext,
197  type&... normal_dot_flux_tilde_s_ext,
198  dg::Formulation dg_formulation) const noexcept {
199  EXPAND_PACK_LEFT_TO_RIGHT(Rusanov_detail::dg_boundary_terms_impl(
200  boundary_correction_tilde_e, boundary_correction_tilde_s, tilde_e_int,
201  tilde_s_int, normal_dot_flux_tilde_e_int, normal_dot_flux_tilde_s_int,
202  tilde_e_ext, tilde_s_ext, normal_dot_flux_tilde_e_ext,
203  normal_dot_flux_tilde_s_ext, dg_formulation));
204  }
205 };
206 
207 /// \cond
208 template <typename... NeutrinoSpecies>
209 // NOLINTNEXTLINE
210 PUP::able::PUP_ID Rusanov<tmpl::list<NeutrinoSpecies...>>::my_PUP_ID = 0;
211 /// \endcond
212 
213 } // namespace RadiationTransport::M1Grey::BoundaryCorrections
CharmPupable.hpp
EXPAND_PACK_LEFT_TO_RIGHT
#define EXPAND_PACK_LEFT_TO_RIGHT(...)
Expand a parameter pack evaluating the terms from left to right.
Definition: TMPL.hpp:601
Frame::Inertial
Definition: IndexType.hpp:44
Options.hpp
Tags::NormalDotFlux
Prefix indicating a boundary unit normal vector dotted into the flux.
Definition: Prefixes.hpp:96
RadiationTransport::M1Grey::Tags::TildeE
The densitized energy density of neutrinos of a given species .
Definition: Tags.hpp:31
RadiationTransport::M1Grey::Tags::TildeS
The densitized momentum density of neutrinos of a given species .
Definition: Tags.hpp:41
WRAPPED_PUPable_decl_template
#define WRAPPED_PUPable_decl_template(className)
Mark derived classes as serializable.
Definition: CharmPupable.hpp:22
RadiationTransport::M1Grey::BoundaryCorrections
Boundary corrections/numerical fluxes.
Definition: BoundaryCorrection.hpp:13
RadiationTransport::M1Grey::BoundaryCorrections::BoundaryCorrection
The base class used to create boundary corrections from input files and store them in the global cach...
Definition: BoundaryCorrection.hpp:24
DataVector
Stores a collection of function values.
Definition: DataVector.hpp:46
memory
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
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