Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <memory>
7 : #include <optional>
8 :
9 : #include "DataStructures/DataBox/Prefixes.hpp"
10 : #include "DataStructures/Tensor/TypeAliases.hpp"
11 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/BoundaryCorrections/BoundaryCorrection.hpp"
12 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/Tags.hpp"
13 : #include "NumericalAlgorithms/DiscontinuousGalerkin/Formulation.hpp"
14 : #include "Options/String.hpp"
15 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
16 : #include "Utilities/Gsl.hpp"
17 : #include "Utilities/Serialization/CharmPupable.hpp"
18 : #include "Utilities/TMPL.hpp"
19 :
20 : /// \cond
21 : class DataVector;
22 : namespace gsl {
23 : template <typename T>
24 : class not_null;
25 : } // namespace gsl
26 : namespace PUP {
27 : class er;
28 : } // namespace PUP
29 : /// \endcond
30 :
31 : namespace grmhd::ValenciaDivClean::BoundaryCorrections {
32 : /*!
33 : * \brief A Rusanov/local Lax-Friedrichs Riemann solver
34 : *
35 : * Let \f$U\f$ be the state vector of evolved variables, \f$F^i\f$ the
36 : * corresponding fluxes, and \f$n_i\f$ be the outward directed unit normal to
37 : * the interface. Denoting \f$F := n_i F^i\f$, the %Rusanov boundary correction
38 : * is
39 : *
40 : * \f{align*}
41 : * G_\text{Rusanov} = \frac{F_\text{int} - F_\text{ext}}{2} -
42 : * \frac{\text{max}\left(\{|\lambda_\text{int}|\},
43 : * \{|\lambda_\text{ext}|\}\right)}{2} \left(U_\text{ext} - U_\text{int}\right),
44 : * \f}
45 : *
46 : * where "int" and "ext" stand for interior and exterior, and
47 : * \f$\{|\lambda|\}\f$ is the set of characteristic/signal speeds. The minus
48 : * sign in front of the \f$F_{\text{ext}}\f$ is necessary because the outward
49 : * directed normal of the neighboring element has the opposite sign, i.e.
50 : * \f$n_i^{\text{ext}}=-n_i^{\text{int}}\f$. The characteristic/signal speeds
51 : * are given in the documentation for
52 : * grmhd::ValenciaDivClean::characteristic_speeds(). Since the fluid is
53 : * travelling slower than the speed of light, the speeds we are interested in
54 : * are
55 : *
56 : * \f{align*}{
57 : * \lambda_{\pm}&=\pm\alpha-\beta^i n_i,
58 : * \f}
59 : *
60 : * which correspond to the divergence cleaning field.
61 : *
62 : * \note
63 : * - In the strong form the `dg_boundary_terms` function returns
64 : * \f$G - F_\text{int}\f$
65 : * - It may be possible to use the slower speeds for the magnetic field and
66 : * fluid part of the system in order to make the flux less dissipative for
67 : * those variables.
68 : */
69 1 : class Rusanov final : public BoundaryCorrection {
70 : private:
71 0 : struct AbsCharSpeed : db::SimpleTag {
72 0 : using type = Scalar<DataVector>;
73 : };
74 :
75 : public:
76 0 : using options = tmpl::list<>;
77 0 : static constexpr Options::String help = {
78 : "Computes the Rusanov or local Lax-Friedrichs boundary correction term "
79 : "for the GRMHD system."};
80 :
81 0 : Rusanov() = default;
82 0 : Rusanov(const Rusanov&) = default;
83 0 : Rusanov& operator=(const Rusanov&) = default;
84 0 : Rusanov(Rusanov&&) = default;
85 0 : Rusanov& operator=(Rusanov&&) = default;
86 0 : ~Rusanov() override = default;
87 :
88 : /// \cond
89 : explicit Rusanov(CkMigrateMessage* /*unused*/);
90 : using PUP::able::register_constructor;
91 : WRAPPED_PUPable_decl_template(Rusanov); // NOLINT
92 : /// \endcond
93 0 : void pup(PUP::er& p) override; // NOLINT
94 :
95 0 : std::unique_ptr<BoundaryCorrection> get_clone() const override;
96 :
97 0 : using dg_package_field_tags =
98 : tmpl::list<Tags::TildeD, Tags::TildeYe, Tags::TildeTau,
99 : Tags::TildeS<Frame::Inertial>, Tags::TildeB<Frame::Inertial>,
100 : Tags::TildePhi, ::Tags::NormalDotFlux<Tags::TildeD>,
101 : ::Tags::NormalDotFlux<Tags::TildeYe>,
102 : ::Tags::NormalDotFlux<Tags::TildeTau>,
103 : ::Tags::NormalDotFlux<Tags::TildeS<Frame::Inertial>>,
104 : ::Tags::NormalDotFlux<Tags::TildeB<Frame::Inertial>>,
105 : ::Tags::NormalDotFlux<Tags::TildePhi>, AbsCharSpeed>;
106 0 : using dg_package_data_temporary_tags =
107 : tmpl::list<gr::Tags::Lapse<DataVector>, gr::Tags::Shift<DataVector, 3>>;
108 0 : using dg_package_data_primitive_tags = tmpl::list<>;
109 0 : using dg_package_data_volume_tags = tmpl::list<>;
110 0 : using dg_boundary_terms_volume_tags = tmpl::list<>;
111 :
112 0 : static double dg_package_data(
113 : gsl::not_null<Scalar<DataVector>*> packaged_tilde_d,
114 : gsl::not_null<Scalar<DataVector>*> packaged_tilde_ye,
115 : gsl::not_null<Scalar<DataVector>*> packaged_tilde_tau,
116 : gsl::not_null<tnsr::i<DataVector, 3, Frame::Inertial>*> packaged_tilde_s,
117 : gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*> packaged_tilde_b,
118 : gsl::not_null<Scalar<DataVector>*> packaged_tilde_phi,
119 : gsl::not_null<Scalar<DataVector>*> packaged_normal_dot_flux_tilde_d,
120 : gsl::not_null<Scalar<DataVector>*> packaged_normal_dot_flux_tilde_ye,
121 : gsl::not_null<Scalar<DataVector>*> packaged_normal_dot_flux_tilde_tau,
122 : gsl::not_null<tnsr::i<DataVector, 3, Frame::Inertial>*>
123 : packaged_normal_dot_flux_tilde_s,
124 : gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*>
125 : packaged_normal_dot_flux_tilde_b,
126 : gsl::not_null<Scalar<DataVector>*> packaged_normal_dot_flux_tilde_phi,
127 : gsl::not_null<Scalar<DataVector>*> packaged_abs_char_speed,
128 :
129 : const Scalar<DataVector>& tilde_d, const Scalar<DataVector>& tilde_ye,
130 : const Scalar<DataVector>& tilde_tau,
131 : const tnsr::i<DataVector, 3, Frame::Inertial>& tilde_s,
132 : const tnsr::I<DataVector, 3, Frame::Inertial>& tilde_b,
133 : const Scalar<DataVector>& tilde_phi,
134 :
135 : const tnsr::I<DataVector, 3, Frame::Inertial>& flux_tilde_d,
136 : const tnsr::I<DataVector, 3, Frame::Inertial>& flux_tilde_ye,
137 : const tnsr::I<DataVector, 3, Frame::Inertial>& flux_tilde_tau,
138 : const tnsr::Ij<DataVector, 3, Frame::Inertial>& flux_tilde_s,
139 : const tnsr::IJ<DataVector, 3, Frame::Inertial>& flux_tilde_b,
140 : const tnsr::I<DataVector, 3, Frame::Inertial>& flux_tilde_phi,
141 :
142 : const Scalar<DataVector>& lapse,
143 : const tnsr::I<DataVector, 3, Frame::Inertial>& shift,
144 :
145 : const tnsr::i<DataVector, 3, Frame::Inertial>& normal_covector,
146 : const tnsr::I<DataVector, 3, Frame::Inertial>& normal_vector,
147 : const std::optional<tnsr::I<DataVector, 3, Frame::Inertial>>&
148 : /*mesh_velocity*/,
149 : const std::optional<Scalar<DataVector>>& normal_dot_mesh_velocity);
150 :
151 0 : static void dg_boundary_terms(
152 : gsl::not_null<Scalar<DataVector>*> boundary_correction_tilde_d,
153 : gsl::not_null<Scalar<DataVector>*> boundary_correction_tilde_ye,
154 : gsl::not_null<Scalar<DataVector>*> boundary_correction_tilde_tau,
155 : gsl::not_null<tnsr::i<DataVector, 3, Frame::Inertial>*>
156 : boundary_correction_tilde_s,
157 : gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*>
158 : boundary_correction_tilde_b,
159 : gsl::not_null<Scalar<DataVector>*> boundary_correction_tilde_phi,
160 : const Scalar<DataVector>& tilde_d_int,
161 : const Scalar<DataVector>& tilde_ye_int,
162 : const Scalar<DataVector>& tilde_tau_int,
163 : const tnsr::i<DataVector, 3, Frame::Inertial>& tilde_s_int,
164 : const tnsr::I<DataVector, 3, Frame::Inertial>& tilde_b_int,
165 : const Scalar<DataVector>& tilde_phi_int,
166 : const Scalar<DataVector>& normal_dot_flux_tilde_d_int,
167 : const Scalar<DataVector>& normal_dot_flux_tilde_ye_int,
168 : const Scalar<DataVector>& normal_dot_flux_tilde_tau_int,
169 : const tnsr::i<DataVector, 3, Frame::Inertial>&
170 : normal_dot_flux_tilde_s_int,
171 : const tnsr::I<DataVector, 3, Frame::Inertial>&
172 : normal_dot_flux_tilde_b_int,
173 : const Scalar<DataVector>& normal_dot_flux_tilde_phi_int,
174 : const Scalar<DataVector>& abs_char_speed_int,
175 : const Scalar<DataVector>& tilde_d_ext,
176 : const Scalar<DataVector>& tilde_ye_ext,
177 : const Scalar<DataVector>& tilde_tau_ext,
178 : const tnsr::i<DataVector, 3, Frame::Inertial>& tilde_s_ext,
179 : const tnsr::I<DataVector, 3, Frame::Inertial>& tilde_b_ext,
180 : const Scalar<DataVector>& tilde_phi_ext,
181 : const Scalar<DataVector>& normal_dot_flux_tilde_d_ext,
182 : const Scalar<DataVector>& normal_dot_flux_tilde_ye_ext,
183 : const Scalar<DataVector>& normal_dot_flux_tilde_tau_ext,
184 : const tnsr::i<DataVector, 3, Frame::Inertial>&
185 : normal_dot_flux_tilde_s_ext,
186 : const tnsr::I<DataVector, 3, Frame::Inertial>&
187 : normal_dot_flux_tilde_b_ext,
188 : const Scalar<DataVector>& normal_dot_flux_tilde_phi_ext,
189 : const Scalar<DataVector>& abs_char_speed_ext,
190 : dg::Formulation dg_formulation);
191 : };
192 :
193 0 : bool operator==(const Rusanov& lhs, const Rusanov& rhs);
194 0 : bool operator!=(const Rusanov& lhs, const Rusanov& rhs);
195 : } // namespace grmhd::ValenciaDivClean::BoundaryCorrections
|