Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include "DataStructures/Tensor/TypeAliases.hpp"
7 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/TagsDeclarations.hpp"
8 : #include "NumericalAlgorithms/LinearOperators/PartialDerivatives.hpp"
9 : #include "PointwiseFunctions/GeneralRelativity/TagsDeclarations.hpp"
10 : #include "PointwiseFunctions/Hydro/Tags.hpp"
11 : #include "PointwiseFunctions/Hydro/TagsDeclarations.hpp"
12 : #include "Utilities/TMPL.hpp"
13 :
14 : /// \cond
15 : class DataVector;
16 :
17 : namespace Tags {
18 : template <typename>
19 : struct Source;
20 : } // namespace Tags
21 : /// \endcond
22 :
23 : namespace grmhd {
24 : namespace ValenciaDivClean {
25 : namespace detail {
26 : void sources_impl(
27 : gsl::not_null<Scalar<DataVector>*> source_tilde_tau,
28 : gsl::not_null<tnsr::i<DataVector, 3, Frame::Inertial>*> source_tilde_s,
29 : gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*> source_tilde_b,
30 : gsl::not_null<Scalar<DataVector>*> source_tilde_phi,
31 :
32 : gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*> tilde_s_up,
33 : gsl::not_null<tnsr::II<DataVector, 3, Frame::Inertial>*> densitized_stress,
34 : gsl::not_null<Scalar<DataVector>*> h_rho_w_squared_plus_b_squared,
35 :
36 : const Scalar<DataVector>& magnetic_field_dot_spatial_velocity,
37 : const Scalar<DataVector>& magnetic_field_squared,
38 : const Scalar<DataVector>& one_over_w_squared,
39 : const Scalar<DataVector>& pressure_star,
40 : const tnsr::I<DataVector, 3, Frame::Inertial>&
41 : trace_spatial_christoffel_second,
42 :
43 : const Scalar<DataVector>& tilde_d, const Scalar<DataVector>& tilde_ye,
44 : const Scalar<DataVector>& tilde_tau,
45 : const tnsr::i<DataVector, 3, Frame::Inertial>& tilde_s,
46 : const tnsr::I<DataVector, 3, Frame::Inertial>& tilde_b,
47 : const Scalar<DataVector>& tilde_phi, const Scalar<DataVector>& lapse,
48 : const Scalar<DataVector>& sqrt_det_spatial_metric,
49 : const tnsr::II<DataVector, 3, Frame::Inertial>& inv_spatial_metric,
50 : const tnsr::i<DataVector, 3, Frame::Inertial>& d_lapse,
51 : const tnsr::iJ<DataVector, 3, Frame::Inertial>& d_shift,
52 : const tnsr::ijj<DataVector, 3, Frame::Inertial>& d_spatial_metric,
53 : const tnsr::I<DataVector, 3, Frame::Inertial>& spatial_velocity,
54 : const Scalar<DataVector>& lorentz_factor,
55 : const tnsr::I<DataVector, 3, Frame::Inertial>& magnetic_field,
56 :
57 : const Scalar<DataVector>& rest_mass_density,
58 : const Scalar<DataVector>& electron_fraction,
59 : const Scalar<DataVector>& pressure,
60 : const Scalar<DataVector>& specific_internal_energy,
61 : const tnsr::ii<DataVector, 3, Frame::Inertial>& extrinsic_curvature,
62 : double constraint_damping_parameter);
63 : } // namespace detail
64 :
65 : /*!
66 : * \brief Compute the source terms for the flux-conservative Valencia
67 : * formulation of GRMHD with divergence cleaning, coupled with electron
68 : * fraction.
69 : *
70 : * A flux-conservative system has the generic form:
71 : * \f[
72 : * \partial_t U_i + \partial_m F^m(U_i) = S(U_i)
73 : * \f]
74 : *
75 : * where \f$F^a()\f$ denotes the flux of a conserved variable \f$U_i\f$ and
76 : * \f$S()\f$ denotes the source term for the conserved variable.
77 : *
78 : * For the Valencia formulation:
79 : * \f{align*}
80 : * S({\tilde D}) = & 0\\
81 : * S({\tilde S}_i) = & \frac{1}{2} \alpha {\tilde S}^{mn} \partial_i \gamma_{mn}
82 : * + {\tilde S}_m \partial_i \beta^m - ({\tilde D} + {\tilde \tau}) \partial_i
83 : * \alpha \\
84 : * S({\tilde \tau}) = & \alpha {\tilde S}^{mn} K_{mn} - {\tilde S}^m \partial_m
85 : * \alpha \\
86 : * S({\tilde B}^i) = & {\tilde \Phi}
87 : * \gamma^{im} \partial_m \alpha + \alpha {\tilde \Phi} \left( \frac{1}{2}
88 : * \gamma^{il} \gamma^{jk} - \gamma^{ij} \gamma^{lk} \right) \partial_l
89 : * \gamma_{jk} \\ S({\tilde \Phi}) = & {\tilde B}^k \partial_k \alpha - \alpha K
90 : * {\tilde \Phi}
91 : * - \alpha \kappa {\tilde \Phi}
92 : * \f}
93 : *
94 : * where
95 : * \f{align*}
96 : * {\tilde S}^i = & {\tilde S}_m \gamma^{im} \\
97 : * {\tilde S}^{ij} = & \sqrt{\gamma} \left[ \left(h \rho W^2 + B^n B_n \right)
98 : * v^i v^j + \left(p + p_m \right) \gamma^{ij} - B^n v_n \left( B^i v^j + B^j
99 : * v^i \right) - \frac{B^i B^j}{W^2} \right]
100 : * \f}
101 : *
102 : * where \f${\tilde D}\f$, \f${\tilde S}_i\f$, \f${\tilde \tau}\f$, \f${\tilde
103 : * B}^i\f$, and \f${\tilde \Phi}\f$ are a generalized mass-energy density,
104 : * momentum density, specific internal energy density, magnetic field, and
105 : * divergence cleaning field. Furthermore, \f$\gamma\f$ is the determinant of
106 : * the spatial metric \f$\gamma_{ij}\f$, \f$\rho\f$ is the rest mass density,
107 : * \f$W\f$ is the Lorentz factor, \f$h\f$ is the specific enthalpy, \f$v^i\f$ is
108 : * the spatial velocity, \f$B^k\f$ is the magnetic field, \f$p\f$ is the
109 : * pressure, \f$p_m = \frac{1}{2} \left[ \left( B^n v_n \right)^2 + B^n B_n /
110 : * W^2 \right]\f$ is the magnetic pressure, \f$\alpha\f$ is the lapse,
111 : * \f$\beta^i\f$ is the shift, \f$K\f$ is the trace of the extrinsic curvature
112 : * \f$K_{ij}\f$, and \f$\kappa\f$ is a damping parameter that damps violations
113 : * of the divergence-free (no-monopole) condition \f$\Phi = \partial_i {\tilde
114 : * B}^i = 0\f$ .
115 : *
116 : * \note For the electron fraction side, the source term is currently set to
117 : * \f$S(\tilde{Y}_e) = 0\f$ where the conserved variable \f$\tilde{Y}_e\f$ is a
118 : * generalized electron fraction. Implementing the source term using neutrino
119 : * scheme is in progress (Last update : Oct 2022).
120 : *
121 : */
122 1 : struct ComputeSources {
123 0 : using return_tags =
124 : tmpl::list<::Tags::Source<grmhd::ValenciaDivClean::Tags::TildeTau>,
125 : ::Tags::Source<grmhd::ValenciaDivClean::Tags::TildeS<>>,
126 : ::Tags::Source<grmhd::ValenciaDivClean::Tags::TildeB<>>,
127 : ::Tags::Source<grmhd::ValenciaDivClean::Tags::TildePhi>>;
128 :
129 0 : using argument_tags =
130 : tmpl::list<grmhd::ValenciaDivClean::Tags::TildeD,
131 : grmhd::ValenciaDivClean::Tags::TildeYe,
132 : grmhd::ValenciaDivClean::Tags::TildeTau,
133 : grmhd::ValenciaDivClean::Tags::TildeS<>,
134 : grmhd::ValenciaDivClean::Tags::TildeB<>,
135 : grmhd::ValenciaDivClean::Tags::TildePhi,
136 : hydro::Tags::SpatialVelocity<DataVector, 3>,
137 : hydro::Tags::MagneticField<DataVector, 3>,
138 : hydro::Tags::RestMassDensity<DataVector>,
139 : hydro::Tags::ElectronFraction<DataVector>,
140 : hydro::Tags::SpecificInternalEnergy<DataVector>,
141 : hydro::Tags::LorentzFactor<DataVector>,
142 : hydro::Tags::Pressure<DataVector>, gr::Tags::Lapse<DataVector>,
143 : ::Tags::deriv<gr::Tags::Lapse<DataVector>, tmpl::size_t<3>,
144 : Frame::Inertial>,
145 : ::Tags::deriv<gr::Tags::Shift<DataVector, 3>, tmpl::size_t<3>,
146 : Frame::Inertial>,
147 : gr::Tags::SpatialMetric<DataVector, 3>,
148 : ::Tags::deriv<gr::Tags::SpatialMetric<DataVector, 3>,
149 : tmpl::size_t<3>, Frame::Inertial>,
150 : gr::Tags::InverseSpatialMetric<DataVector, 3>,
151 : gr::Tags::SqrtDetSpatialMetric<DataVector>,
152 : gr::Tags::ExtrinsicCurvature<DataVector, 3>,
153 : grmhd::ValenciaDivClean::Tags::ConstraintDampingParameter>;
154 :
155 0 : static void apply(
156 : gsl::not_null<Scalar<DataVector>*> source_tilde_tau,
157 : gsl::not_null<tnsr::i<DataVector, 3, Frame::Inertial>*> source_tilde_s,
158 : gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*> source_tilde_b,
159 : gsl::not_null<Scalar<DataVector>*> source_tilde_phi,
160 : const Scalar<DataVector>& tilde_d, const Scalar<DataVector>& tilde_ye,
161 : const Scalar<DataVector>& tilde_tau,
162 : const tnsr::i<DataVector, 3, Frame::Inertial>& tilde_s,
163 : const tnsr::I<DataVector, 3, Frame::Inertial>& tilde_b,
164 : const Scalar<DataVector>& tilde_phi,
165 : const tnsr::I<DataVector, 3, Frame::Inertial>& spatial_velocity,
166 : const tnsr::I<DataVector, 3, Frame::Inertial>& magnetic_field,
167 : const Scalar<DataVector>& rest_mass_density,
168 : const Scalar<DataVector>& electron_fraction,
169 : const Scalar<DataVector>& specific_internal_energy,
170 : const Scalar<DataVector>& lorentz_factor,
171 : const Scalar<DataVector>& pressure, const Scalar<DataVector>& lapse,
172 : const tnsr::i<DataVector, 3, Frame::Inertial>& d_lapse,
173 : const tnsr::iJ<DataVector, 3, Frame::Inertial>& d_shift,
174 : const tnsr::ii<DataVector, 3, Frame::Inertial>& spatial_metric,
175 : const tnsr::ijj<DataVector, 3, Frame::Inertial>& d_spatial_metric,
176 : const tnsr::II<DataVector, 3, Frame::Inertial>& inv_spatial_metric,
177 : const Scalar<DataVector>& sqrt_det_spatial_metric,
178 : const tnsr::ii<DataVector, 3, Frame::Inertial>& extrinsic_curvature,
179 : double constraint_damping_parameter);
180 : };
181 : } // namespace ValenciaDivClean
182 : } // namespace grmhd
|
