Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <array>
7 : #include <cstddef>
8 :
9 : #include "DataStructures/DataBox/Tag.hpp"
10 : #include "DataStructures/Tensor/EagerMath/Magnitude.hpp"
11 : #include "DataStructures/Tensor/TypeAliases.hpp"
12 : #include "Domain/FaceNormal.hpp"
13 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/Tags.hpp"
14 : #include "PointwiseFunctions/GeneralRelativity/TagsDeclarations.hpp"
15 : #include "PointwiseFunctions/Hydro/EquationsOfState/EquationOfState.hpp"
16 : #include "PointwiseFunctions/Hydro/TagsDeclarations.hpp"
17 : #include "Utilities/Gsl.hpp"
18 : #include "Utilities/TMPL.hpp"
19 :
20 : /// \cond
21 : class DataVector;
22 : namespace Tags {
23 : template <typename Tag>
24 : struct Normalized;
25 : } // namespace Tags
26 : /// \endcond
27 :
28 : namespace grmhd {
29 : namespace ValenciaDivClean {
30 : /// @{
31 : /*!
32 : * \brief Compute the characteristic speeds for the Valencia formulation of
33 : * GRMHD with divergence cleaning.
34 : *
35 : * Obtaining the exact form of the characteristic speeds involves the solution
36 : * of a nontrivial quartic equation for the fast and slow modes. Here we make
37 : * use of a common approximation in the literature (e.g. \cite Gammie2003)
38 : * where the resulting characteristic speeds are analogous to those of the
39 : * Valencia formulation of the 3-D relativistic Euler system
40 : * (see RelativisticEuler::Valencia::characteristic_speeds),
41 : *
42 : * \f{align*}
43 : * \lambda_2 &= \alpha \Lambda^- - \beta_n,\\
44 : * \lambda_{3, 4, 5, 6, 7} &= \alpha v_n - \beta_n,\\
45 : * \lambda_{8} &= \alpha \Lambda^+ - \beta_n,
46 : * \f}
47 : *
48 : * with the substitution
49 : *
50 : * \f{align*}
51 : * c_s^2 \longrightarrow c_s^2 + v_A^2(1 - c_s^2)
52 : * \f}
53 : *
54 : * in the definition of \f$\Lambda^\pm\f$. Here \f$v_A\f$ is the Alfvén
55 : * speed. In addition, two more speeds corresponding to the divergence cleaning
56 : * mode and the longitudinal magnetic field are added,
57 : *
58 : * \f{align*}
59 : * \lambda_1 = -\alpha - \beta_n,\\
60 : * \lambda_9 = \alpha - \beta_n.
61 : * \f}
62 : *
63 : * \note The ordering assumed here is such that, in the Newtonian limit,
64 : * the exact expressions for \f$\lambda_{2, 8}\f$, \f$\lambda_{3, 7}\f$,
65 : * and \f$\lambda_{4, 6}\f$ should reduce to the
66 : * corresponding fast modes, Alfvén modes, and slow modes, respectively.
67 : * See \cite Dedner2002 for a detailed description of the hyperbolic
68 : * characterization of Newtonian MHD. In terms of the primitive variables:
69 : *
70 : * \f{align*}
71 : * v^2 &= \gamma_{mn} v^m v^n \\
72 : * c_s^2 &= \frac{1}{h} \left[ \left( \frac{\partial p}{\partial \rho}
73 : * \right)_\epsilon +
74 : * \frac{p}{\rho^2} \left(\frac{\partial p}{\partial \epsilon}
75 : * \right)_\rho \right] \\
76 : * v_A^2 &= \frac{b^2}{b^2 + \rho h} \\
77 : * b^2 &= \frac{1}{W^2} \gamma_{mn} B^m B^n + \left( \gamma_{mn} B^m v^n
78 : * \right)^2
79 : * \f}
80 : *
81 : * where \f$\gamma_{mn}\f$ is the spatial metric, \f$\rho\f$ is the rest
82 : * mass density, \f$W = 1/\sqrt{1-v_i v^i}\f$ is the Lorentz factor, \f$h = 1 +
83 : * \epsilon + \frac{p}{\rho}\f$ is the specific enthalpy, \f$v^i\f$ is the
84 : * spatial velocity, \f$\epsilon\f$ is the specific internal energy, \f$p\f$ is
85 : * the pressure, and \f$B^i\f$ is the spatial magnetic field measured by an
86 : * Eulerian observer.
87 : */
88 :
89 : template <size_t ThermodynamicDim>
90 1 : std::array<DataVector, 9> characteristic_speeds_approximate_mhd(
91 : const Scalar<DataVector>& rest_mass_density,
92 : const Scalar<DataVector>& electron_fraction,
93 : const Scalar<DataVector>& specific_internal_energy,
94 : const Scalar<DataVector>& specific_enthalpy,
95 : const tnsr::I<DataVector, 3, Frame::Inertial>& spatial_velocity,
96 : const Scalar<DataVector>& lorentz_factor,
97 : const tnsr::I<DataVector, 3, Frame::Inertial>& magnetic_field,
98 : const Scalar<DataVector>& lapse, const tnsr::I<DataVector, 3>& shift,
99 : const tnsr::ii<DataVector, 3, Frame::Inertial>& spatial_metric,
100 : const tnsr::i<DataVector, 3>& unit_normal,
101 : const EquationsOfState::EquationOfState<true, ThermodynamicDim>&
102 : equation_of_state);
103 :
104 : template <size_t ThermodynamicDim>
105 1 : void characteristic_speeds_approximate_mhd(
106 : gsl::not_null<std::array<DataVector, 9>*> char_speeds,
107 : const Scalar<DataVector>& rest_mass_density,
108 : const Scalar<DataVector>& electron_fraction,
109 : const Scalar<DataVector>& specific_internal_energy,
110 : const Scalar<DataVector>& specific_enthalpy,
111 : const tnsr::I<DataVector, 3, Frame::Inertial>& spatial_velocity,
112 : const Scalar<DataVector>& lorentz_factor,
113 : const tnsr::I<DataVector, 3, Frame::Inertial>& magnetic_field,
114 : const Scalar<DataVector>& lapse, const tnsr::I<DataVector, 3>& shift,
115 : const tnsr::ii<DataVector, 3, Frame::Inertial>& spatial_metric,
116 : const tnsr::i<DataVector, 3>& unit_normal,
117 : const EquationsOfState::EquationOfState<true, ThermodynamicDim>&
118 : equation_of_state);
119 : /// @}
120 :
121 : /// @{
122 : /*!
123 : * \brief Compute the characteristic speeds for the relativistic hydrodynamics
124 : * system in the Eulerian frame.
125 : *
126 : * These are the eigenvalues of the flux Jacobian projected along a spatial
127 : * direction with unit normal \f$ s_i\f$, measured by an Eulerian observer.
128 : * They consist of four degenerate modes and two non-degenerate modes:
129 : *
130 : * \f{align*}
131 : * \lambda_\mathrm{deg} &= v_n ,\\
132 : * \lambda_{\pm} &=
133 : * \frac{ (1 - c_s^2)\,v_n \pm c_s\,d / W^2 }
134 : * { 1 - v^2 c_s^2 }.
135 : * \f}
136 : *
137 : * where
138 : *
139 : * \f{align*}
140 : * d &\equiv
141 : * W\,\sqrt{ 1 - v^2 c_s^2 - v_n^2(1 - c_s^2) } .
142 : * \f}
143 : *
144 : * The variables are defined as:
145 : *
146 : * \f{align*}
147 : * v^2 &= \gamma_{mn} v^m v^n ,\\
148 : * v_n &= v^a s_a ,\\
149 : * W &= \frac{1}{\sqrt{1 - v^a v_a}} ,\\
150 : * h &= 1 + \epsilon + \frac{p}{\rho} ,\\
151 : * c_s^2 &= \frac{1}{h}\left(\chi + \kappa \frac{p}{\rho^2}\right) ,
152 : * \f}
153 : *
154 : * The returned array has size 3 and contains the unique characteristic speeds:
155 : * - `char_speeds[0]`: \f$v_n\f$ (degenerate eigenvalue, multiplicity 4 if the
156 : * electron fraction is included),
157 : * - `char_speeds[1]`: \f$\lambda_+\f$ (right-propagating acoustic mode),
158 : * - `char_speeds[2]`: \f$\lambda_-\f$ (left-propagating acoustic mode).
159 : */
160 :
161 : template <size_t ThermodynamicDim>
162 1 : void characteristic_speeds_hydro(
163 : gsl::not_null<std::array<DataVector, 3>*> pchar_speeds,
164 : const tnsr::I<DataVector, 3, Frame::Inertial>& spatial_velocity,
165 : const Scalar<DataVector>& rest_mass_density,
166 : const Scalar<DataVector>& specific_internal_energy,
167 : const Scalar<DataVector>& specific_enthalpy,
168 : const Scalar<DataVector>& electron_fraction,
169 : const Scalar<DataVector>& lorentz_factor,
170 : const tnsr::i<DataVector, 3>& unit_normal,
171 : const tnsr::ii<DataVector, 3, Frame::Inertial>& spatial_metric,
172 : const EquationsOfState::EquationOfState<true, ThermodynamicDim>&
173 : equation_of_state);
174 :
175 : template <size_t ThermodynamicDim>
176 1 : std::array<DataVector, 3> characteristic_speeds_hydro(
177 : const tnsr::I<DataVector, 3, Frame::Inertial>& spatial_velocity,
178 : const Scalar<DataVector>& rest_mass_density,
179 : const Scalar<DataVector>& specific_internal_energy,
180 : const Scalar<DataVector>& specific_enthalpy,
181 : const Scalar<DataVector>& electron_fraction,
182 : const Scalar<DataVector>& lorentz_factor,
183 : const tnsr::i<DataVector, 3>& unit_normal,
184 : const tnsr::ii<DataVector, 3, Frame::Inertial>& spatial_metric,
185 : const EquationsOfState::EquationOfState<true, ThermodynamicDim>&
186 : equation_of_state);
187 : /// @}
188 :
189 1 : namespace Tags {
190 : /// \brief Compute the characteristic speeds for the Valencia formulation of
191 : /// GRMHD with divergence cleaning.
192 : ///
193 : /// \details see grmhd::ValenciaDivClean::characteristic_speeds
194 1 : struct CharacteristicSpeedsCompute : Tags::CharacteristicSpeeds,
195 : db::ComputeTag {
196 0 : using base = Tags::CharacteristicSpeeds;
197 0 : using argument_tags =
198 : tmpl::list<hydro::Tags::RestMassDensity<DataVector>,
199 : hydro::Tags::ElectronFraction<DataVector>,
200 : hydro::Tags::SpecificInternalEnergy<DataVector>,
201 : hydro::Tags::SpecificEnthalpy<DataVector>,
202 : hydro::Tags::SpatialVelocity<DataVector, 3>,
203 : hydro::Tags::LorentzFactor<DataVector>,
204 : hydro::Tags::MagneticField<DataVector, 3>,
205 : gr::Tags::Lapse<DataVector>, gr::Tags::Shift<DataVector, 3>,
206 : gr::Tags::SpatialMetric<DataVector, 3>,
207 : ::Tags::Normalized<domain::Tags::UnnormalizedFaceNormal<3>>,
208 : hydro::Tags::GrmhdEquationOfState>;
209 :
210 0 : using volume_tags = tmpl::list<hydro::Tags::GrmhdEquationOfState>;
211 :
212 0 : using return_type = std::array<DataVector, 9>;
213 :
214 : template <size_t ThermodynamicDim>
215 0 : void function(gsl::not_null<return_type*> result,
216 : const Scalar<DataVector>& rest_mass_density,
217 : const Scalar<DataVector>& /* electron_fraction */,
218 : const Scalar<DataVector>& specific_internal_energy,
219 : const Scalar<DataVector>& specific_enthalpy,
220 : const tnsr::I<DataVector, 3, Frame::Inertial>& spatial_velocity,
221 : const Scalar<DataVector>& lorentz_factor,
222 : const tnsr::I<DataVector, 3, Frame::Inertial>& magnetic_field,
223 : const Scalar<DataVector>& lapse,
224 : const tnsr::I<DataVector, 3>& shift,
225 : const tnsr::ii<DataVector, 3, Frame::Inertial>& spatial_metric,
226 : const tnsr::i<DataVector, 3>& unit_normal,
227 : const EquationsOfState::EquationOfState<true, ThermodynamicDim>&
228 : equation_of_state);
229 : };
230 :
231 0 : struct LargestCharacteristicSpeed : db::SimpleTag {
232 0 : using type = double;
233 : };
234 :
235 0 : struct ComputeLargestCharacteristicSpeed : db::ComputeTag,
236 : LargestCharacteristicSpeed {
237 0 : using argument_tags =
238 : tmpl::list<gr::Tags::Lapse<DataVector>, gr::Tags::Shift<DataVector, 3>,
239 : gr::Tags::SpatialMetric<DataVector, 3>>;
240 0 : using return_type = double;
241 0 : using base = LargestCharacteristicSpeed;
242 0 : static void function(gsl::not_null<double*> speed,
243 : const Scalar<DataVector>& lapse,
244 : const tnsr::I<DataVector, 3>& shift,
245 : const tnsr::ii<DataVector, 3>& spatial_metric) {
246 : const auto shift_magnitude = magnitude(shift, spatial_metric);
247 : *speed = max(get(shift_magnitude) + get(lapse));
248 : }
249 : };
250 : } // namespace Tags
251 : } // namespace ValenciaDivClean
252 : } // namespace grmhd
|