Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <cstddef>
7 :
8 : #include "DataStructures/VariablesTag.hpp"
9 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/BoundaryConditions/BoundaryCondition.hpp"
10 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/Characteristics.hpp"
11 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/ConservativeFromPrimitive.hpp"
12 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/NewmanHamlin.hpp"
13 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/PrimitiveFromConservative.hpp"
14 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/Tags.hpp"
15 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/TimeDerivativeTerms.hpp"
16 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
17 : #include "PointwiseFunctions/Hydro/Tags.hpp"
18 : #include "Utilities/TMPL.hpp"
19 :
20 : /// \ingroup EvolutionSystemsGroup
21 : /// \brief Items related to general relativistic magnetohydrodynamics (GRMHD)
22 : namespace grmhd {
23 : /// Items related to the Valencia formulation of ideal GRMHD with divergence
24 : /// cleaning coupled with electron fraction.
25 : ///
26 : /// References:
27 : /// - Numerical 3+1 General Relativistic Magnetohydrodynamics: A Local
28 : /// Characteristic Approach \cite Anton2006
29 : /// - GRHydro: a new open-source general-relativistic magnetohydrodynamics code
30 : /// for the Einstein toolkit \cite Moesta2014
31 : /// - Black hole-neutron star mergers with a hot nuclear equation of state:
32 : /// outflow and neutrino-cooled disk for a low-mass, high-spin case
33 : /// \cite Deaton2013
34 : ///
35 : namespace ValenciaDivClean {
36 :
37 : /*!
38 : * \brief Ideal general relativistic magnetohydrodynamics (GRMHD) system with
39 : * divergence cleaning coupled with electron fraction.
40 : *
41 : * We adopt the standard 3+1 form of metric
42 : * \f{align*}
43 : * ds^2 = -\alpha^2 dt^2 + \gamma_{ij}(dx^i + \beta^i dt)(dx^j + \beta^j dt)
44 : * \f}
45 : * where \f$\alpha\f$ is the lapse, \f$\beta^i\f$ is the shift vector, and
46 : * \f$\gamma_{ij}\f$ is the spatial metric.
47 : *
48 : * Primitive variables of the system are
49 : * - rest mass density \f$\rho\f$
50 : * - electron fraction \f$Y_e\f$
51 : * - the spatial velocity \f{align*} v^i =
52 : * \frac{1}{\alpha}\left(\frac{u^i}{u^0} + \beta^i\right) \f} measured by the
53 : * Eulerian observer
54 : * - specific internal energy density \f$\epsilon\f$
55 : * - magnetic field \f$B^i = -^*F^{ia}n_a = -\alpha (^*F^{0i})\f$
56 : * - the divergence cleaning field \f$\Phi\f$
57 : *
58 : * with corresponding derived physical quantities which frequently appear in
59 : * equations:
60 : * - The transport velocity \f$v^i_{tr} = \alpha v^i - \beta^i\f$
61 : * - The Lorentz factor \f{align*} W = -u^an_a = \alpha u^0 =
62 : * \frac{1}{\sqrt{1-\gamma_{ij}v^iv^j}} = \sqrt{1+\gamma^{ij}u_iu_j} =
63 : * \sqrt{1+\gamma^{ij}W^2v_iv_j} \f}
64 : * - The specific enthalpy \f$h = 1 + \epsilon + p/\rho\f$ where \f$p\f$ is the
65 : * pressure specified by a particular equation of state (EoS)
66 : * - The comoving magnetic field \f$b^b = -^*F^{ba} u_a\f$ in the ideal MHD
67 : * limit
68 : * \f{align*}
69 : * b^0 & = W B^i v_i / \alpha \\
70 : * b^i & = (B^i + \alpha b^0 u^i) / W \\
71 : * b^2 & = b^ab_a = B^2/W^2 + (B^iv_i)^2
72 : * \f}
73 : * - Augmented enthalpy density \f$(\rho h)^* = \rho h + b^2\f$ and augmented
74 : * pressure \f$p^* = p + b^2/2\f$ which include contributions from magnetic
75 : * field
76 : *
77 : * \note We are using the electromagnetic variables with the scaling convention
78 : * that the factor \f$4\pi\f$ does not appear in Maxwell equations and the
79 : * stress-energy tensor of the EM fields (geometrized Heaviside-Lorentz units).
80 : * To recover the physical value of magnetic field in the usual CGS Gaussian
81 : * unit, the conversion factor is
82 : * \f{align*}
83 : * \sqrt{4\pi}\frac{c^4}{G^{3/2}M_\odot} \approx 8.35 \times 10^{19}
84 : * \,\text{Gauss}
85 : * \f}
86 : * For example, magnetic field $10^{-5}$ with the code unit corresponds to the
87 : * $8.35\times 10^{14}\,\text{G}$ in the CGS Gaussian unit. See also
88 : * documentation of hydro::units::cgs::gauss_unit for details.
89 : *
90 : * The GRMHD equations can be written in a flux-balanced form
91 : * \f[
92 : * \partial_t U+ \partial_i F^i(U) = S(U).
93 : * \f]
94 : *
95 : * Evolved (conserved) variables \f$U\f$ are
96 : * \f{align*}{
97 : * U = \sqrt{\gamma}\left[\,\begin{matrix}
98 : * D \\
99 : * D Y_e \\
100 : * S_j \\
101 : * \tau \\
102 : * B^j \\
103 : * \Phi \\
104 : * \end{matrix}\,\right] \equiv \left[\,\,\begin{matrix}
105 : * \tilde{D} \\
106 : * \tilde{Y_e} \\
107 : * \tilde{S}_j \\
108 : * \tilde{\tau} \\
109 : * \tilde{B}^j \\
110 : * \tilde{\Phi} \\
111 : * \end{matrix}\,\,\right] = \sqrt{\gamma} \left[\,\,\begin{matrix}
112 : * \rho W \\
113 : * \rho W Y_e \\
114 : * (\rho h)^* W^2 v_j - \alpha b^0 b_j \\
115 : * (\rho h)^* W^2 - p^* - (\alpha b^0)^2 - \rho W \\
116 : * B^j \\
117 : * \Phi \\
118 : * \end{matrix}\,\,\right]
119 : * \f}
120 : *
121 : * where \f${\tilde D}\f$, \f${\tilde Y}_e\f$,\f${\tilde S}_i\f$, \f${\tilde
122 : * \tau}\f$, \f${\tilde B}^i\f$, and \f${\tilde \Phi}\f$ are a generalized
123 : * mass-energy density, electron fraction, momentum density, specific internal
124 : * energy density, magnetic field, and divergence cleaning field. Also,
125 : * \f$\gamma\f$ is the determinant of the spatial metric \f$\gamma_{ij}\f$.
126 : *
127 : * Corresponding fluxes \f$F^i(U)\f$ are
128 : *
129 : * \f{align*}
130 : * F^i({\tilde D}) &= {\tilde D} v^i_{tr} \\
131 : * F^i({\tilde Y}_e) &= {\tilde Y}_e v^i_{tr} \\
132 : * F^i({\tilde S}_j) &= {\tilde S}_j v^i_{tr} + \alpha \sqrt{\gamma} p^*
133 : * \delta^i_j - \alpha b_j \tilde{B}^i / W \\
134 : * F^i({\tilde \tau}) &= {\tilde \tau} v^i_{tr} + \alpha \sqrt{\gamma} p^* v^i
135 : * - \alpha^2 b^0 \tilde{B}^i / W \\
136 : * F^i({\tilde B}^j) &= {\tilde B}^j v^i_{tr} - \alpha v^j {\tilde B}^i
137 : * + \alpha \gamma^{ij} {\tilde \Phi} \\
138 : * F^i({\tilde \Phi}) &= \alpha {\tilde B^i} - {\tilde \Phi} \beta^i
139 : * \f}
140 : *
141 : * and source terms \f$S(U)\f$ are
142 : *
143 : * \f{align*}
144 : * S({\tilde D}) &= 0 \\
145 : * S({\tilde Y}_e) &= 0 \\
146 : * S({\tilde S}_j) &= \frac{\alpha}{2} {\tilde S}^{mn} \partial_j \gamma_{mn}
147 : * + {\tilde S}_k \partial_j \beta^k - ({\tilde D}
148 : * + {\tilde \tau}) \partial_j \alpha \\
149 : * S({\tilde \tau}) &= \alpha {\tilde S}^{mn} K_{mn}
150 : * - {\tilde S}^m \partial_m \alpha \\
151 : * S({\tilde B}^j) &= -{\tilde B}^m \partial_m \beta^j
152 : * + \Phi \partial_k (\alpha \sqrt{\gamma}\gamma^{jk}) \\
153 : * S({\tilde \Phi}) &= {\tilde B}^k \partial_k \alpha - \alpha K
154 : * {\tilde \Phi} - \alpha \kappa {\tilde \Phi}
155 : * \f}
156 : *
157 : * with
158 : * \f{align*}
159 : * {\tilde S}^{ij} = & \sqrt{\gamma} \left[ (\rho h)^* W^2 v^i v^j + p^*
160 : * \gamma^{ij} - \gamma^{ik}\gamma^{jl}b_kb_l \right]
161 : * \f}
162 : *
163 : * where \f$K\f$ is the trace of the extrinsic curvature \f$K_{ij}\f$, and
164 : * \f$\kappa\f$ is a damping parameter that damps violations of the
165 : * divergence-free (no-monopole) condition \f$\Phi = \partial_i {\tilde B}^i =
166 : * 0\f$.
167 : *
168 : * \note On the electron fraction side, the source term is currently set to
169 : * \f$S(\tilde{Y}_e) = 0\f$. Implementing the source term using neutrino scheme
170 : * is in progress (Last update : Oct 2022).
171 : *
172 : */
173 1 : struct System {
174 0 : static std::string name() { return "ValenciaDivClean"; }
175 :
176 0 : static constexpr bool is_in_flux_conservative_form = true;
177 0 : static constexpr bool has_primitive_and_conservative_vars = true;
178 0 : static constexpr size_t volume_dim = 3;
179 :
180 0 : using boundary_conditions_base = BoundaryConditions::BoundaryCondition;
181 :
182 0 : using variables_tag = ::Tags::Variables<
183 : tmpl::list<Tags::TildeD, Tags::TildeYe, Tags::TildeTau, Tags::TildeS<>,
184 : Tags::TildeB<>, Tags::TildePhi>>;
185 0 : using flux_variables =
186 : tmpl::list<Tags::TildeD, Tags::TildeYe, Tags::TildeTau, Tags::TildeS<>,
187 : Tags::TildeB<>, Tags::TildePhi>;
188 0 : using non_conservative_variables = tmpl::list<>;
189 0 : using gradient_variables = tmpl::list<>;
190 0 : using primitive_variables_tag =
191 : ::Tags::Variables<hydro::grmhd_tags<DataVector>>;
192 0 : using spacetime_variables_tag =
193 : ::Tags::Variables<gr::tags_for_hydro<volume_dim, DataVector>>;
194 0 : using flux_spacetime_variables_tag = ::Tags::Variables<
195 : tmpl::list<gr::Tags::Lapse<DataVector>, gr::Tags::Shift<DataVector, 3>,
196 : gr::Tags::SpatialMetric<DataVector, 3>,
197 : gr::Tags::SqrtDetSpatialMetric<DataVector>,
198 : gr::Tags::InverseSpatialMetric<DataVector, 3>>>;
199 :
200 0 : using compute_volume_time_derivative_terms = TimeDerivativeTerms;
201 :
202 0 : using conservative_from_primitive = ConservativeFromPrimitive;
203 : template <typename OrderedListOfPrimitiveRecoverySchemes>
204 0 : using primitive_from_conservative =
205 : PrimitiveFromConservative<OrderedListOfPrimitiveRecoverySchemes>;
206 :
207 0 : using compute_largest_characteristic_speed =
208 : Tags::ComputeLargestCharacteristicSpeed;
209 :
210 0 : using inverse_spatial_metric_tag =
211 : gr::Tags::InverseSpatialMetric<DataVector, volume_dim>;
212 : };
213 : } // namespace ValenciaDivClean
214 : } // namespace grmhd
|