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 : #include <utility>
8 :
9 : #include "DataStructures/DataBox/PrefixHelpers.hpp"
10 : #include "DataStructures/DataBox/Prefixes.hpp"
11 : #include "DataStructures/DataVector.hpp"
12 : #include "DataStructures/TaggedContainers.hpp"
13 : #include "DataStructures/Tensor/Tensor.hpp"
14 : #include "Evolution/PassVariables.hpp"
15 : #include "Evolution/Systems/GeneralizedHarmonic/GaugeSourceFunctions/Harmonic.hpp"
16 : #include "Evolution/Systems/GeneralizedHarmonic/System.hpp"
17 : #include "Evolution/Systems/GeneralizedHarmonic/TimeDerivative.hpp"
18 : #include "Evolution/Systems/GrMhd/GhValenciaDivClean/StressEnergy.hpp"
19 : #include "Evolution/Systems/GrMhd/GhValenciaDivClean/System.hpp"
20 : #include "Evolution/Systems/GrMhd/GhValenciaDivClean/Tags.hpp"
21 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/System.hpp"
22 : #include "Evolution/Systems/GrMhd/ValenciaDivClean/TimeDerivativeTerms.hpp"
23 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
24 : // Tag obtained from gh::TimeDerivative needs to be complete here to
25 : // be used in TemporaryReference.
26 : #include "PointwiseFunctions/Hydro/Tags.hpp"
27 : #include "Time/Tags/Time.hpp"
28 : #include "Utilities/Gsl.hpp"
29 : #include "Utilities/Literals.hpp"
30 : #include "Utilities/TMPL.hpp"
31 : #include "Utilities/TaggedTuple.hpp"
32 :
33 : namespace grmhd::GhValenciaDivClean {
34 : namespace detail {
35 : // Some temporary tags appear both in the GRMHD temporary list and the
36 : // GeneralizedHarmonic temporary list, so we wrap the GRMHD temporaries in this
37 : // prefix tag to avoid collisions in data structures used to store the
38 : // temporaries.
39 : template <typename Tag>
40 : struct ValenciaTempTag : db::SimpleTag, db::PrefixTag {
41 : using tag = Tag;
42 : using type = typename Tag::type;
43 : };
44 :
45 : template <typename GhDtTagList, typename ValenciaDtTagList,
46 : typename ValenciaFluxTagList, typename GhTempTagList,
47 : typename ValenciaTempTagList, typename GhGradientTagList,
48 : typename GhArgTagList, typename ValenciaArgTagList,
49 : typename ValenciaTimeDerivativeArgTagList,
50 : typename TraceReversedStressResultTagsList,
51 : typename TraceReversedStressArgumentTagsList>
52 : struct TimeDerivativeTermsImpl;
53 :
54 : template <typename... GhDtTags, typename... ValenciaDtTags,
55 : typename... ValenciaFluxTags, typename... GhTempTags,
56 : typename... ValenciaTempTags, typename... GhGradientTags,
57 : typename... GhArgTags, typename... ValenciaArgTags,
58 : typename... ValenciaTimeDerivativeArgTags,
59 : typename... TraceReversedStressResultTags,
60 : typename... TraceReversedStressArgumentTags>
61 : struct TimeDerivativeTermsImpl<
62 : tmpl::list<GhDtTags...>, tmpl::list<ValenciaDtTags...>,
63 : tmpl::list<ValenciaFluxTags...>, tmpl::list<GhTempTags...>,
64 : tmpl::list<ValenciaTempTags...>, tmpl::list<GhGradientTags...>,
65 : tmpl::list<GhArgTags...>, tmpl::list<ValenciaArgTags...>,
66 : tmpl::list<ValenciaTimeDerivativeArgTags...>,
67 : tmpl::list<TraceReversedStressResultTags...>,
68 : tmpl::list<TraceReversedStressArgumentTags...>> {
69 : template <typename TemporaryTagsList, typename... ExtraTags>
70 : static void apply(
71 : const gsl::not_null<
72 : Variables<tmpl::list<GhDtTags..., ValenciaDtTags...>>*>
73 : dt_vars_ptr,
74 : const gsl::not_null<Variables<db::wrap_tags_in<
75 : ::Tags::Flux, typename ValenciaDivClean::System::flux_variables,
76 : tmpl::size_t<3>, Frame::Inertial>>*>
77 : fluxes_ptr,
78 : const gsl::not_null<Variables<TemporaryTagsList>*> temps_ptr,
79 :
80 : const tnsr::iaa<DataVector, 3>& d_spacetime_metric,
81 : const tnsr::iaa<DataVector, 3>& d_pi,
82 : const tnsr::ijaa<DataVector, 3>& d_phi,
83 :
84 : const tuples::TaggedTuple<ExtraTags...>& arguments) {
85 : gh::TimeDerivative<3_st>::apply(
86 : get<GhDtTags>(dt_vars_ptr)..., get<GhTempTags>(temps_ptr)...,
87 : d_spacetime_metric, d_pi, d_phi,
88 : get<Tags::detail::TemporaryReference<GhArgTags>>(arguments)...);
89 :
90 : if (get<Tags::detail::TemporaryReference<gh::gauges::Tags::GaugeCondition>>(
91 : arguments)
92 : .is_harmonic()) {
93 : get(get<gr::Tags::SqrtDetSpatialMetric<DataVector>>(*temps_ptr)) =
94 : sqrt(get(get<gr::Tags::DetSpatialMetric<DataVector>>(*temps_ptr)));
95 : }
96 :
97 : for (size_t i = 0; i < 3; ++i) {
98 : get<::Tags::deriv<gr::Tags::Lapse<DataVector>, tmpl::size_t<3>,
99 : Frame::Inertial>>(*temps_ptr)
100 : .get(i) = -get(get<gr::Tags::Lapse<DataVector>>(*temps_ptr)) *
101 : get<gh::Tags::HalfPhiTwoNormals<3>>(*temps_ptr).get(i);
102 : }
103 : const auto& phi =
104 : get<Tags::detail::TemporaryReference<gh::Tags::Phi<DataVector, 3>>>(
105 : arguments);
106 : const auto& inv_spatial_metric =
107 : get<gr::Tags::InverseSpatialMetric<DataVector, 3>>(*temps_ptr);
108 : const auto& shift = get<gr::Tags::Shift<DataVector, 3>>(*temps_ptr);
109 : for (size_t i = 0; i < 3; ++i) {
110 : for (size_t j = 0; j < 3; ++j) {
111 : get<::Tags::deriv<gr::Tags::Shift<DataVector, 3>, tmpl::size_t<3>,
112 : Frame::Inertial>>(*temps_ptr)
113 : .get(i, j) = inv_spatial_metric.get(j, 0) * phi.get(i, 0, 1);
114 : for (size_t k = 1; k < 3; ++k) {
115 : get<::Tags::deriv<gr::Tags::Shift<DataVector, 3>, tmpl::size_t<3>,
116 : Frame::Inertial>>(*temps_ptr)
117 : .get(i, j) += inv_spatial_metric.get(j, k) * phi.get(i, 0, k + 1);
118 : }
119 : for (size_t k = 0; k < 3; ++k) {
120 : for (size_t l = 0; l < 3; ++l) {
121 : get<::Tags::deriv<gr::Tags::Shift<DataVector, 3>, tmpl::size_t<3>,
122 : Frame::Inertial>>(*temps_ptr)
123 : .get(i, j) -= shift.get(k) * inv_spatial_metric.get(j, l) *
124 : phi.get(i, l + 1, k + 1);
125 : }
126 : }
127 : }
128 : }
129 :
130 : const auto& pi =
131 : get<Tags::detail::TemporaryReference<gh::Tags::Pi<DataVector, 3>>>(
132 : arguments);
133 : for (size_t i = 0; i < 3; ++i) {
134 : for (size_t j = i; j < 3; ++j) {
135 : get<gr::Tags::ExtrinsicCurvature<DataVector, 3>>(*temps_ptr).get(i, j) =
136 : 0.5 * (pi.get(i + 1, j + 1) +
137 : get<gh::Tags::PhiOneNormal<3>>(*temps_ptr).get(i, j + 1) +
138 : get<gh::Tags::PhiOneNormal<3>>(*temps_ptr).get(j, i + 1));
139 : }
140 : }
141 :
142 : using extra_tags_list = tmpl::list<ExtraTags...>;
143 :
144 : grmhd::ValenciaDivClean::TimeDerivativeTerms::apply(
145 : get<ValenciaDtTags>(dt_vars_ptr)...,
146 : get<ValenciaFluxTags>(fluxes_ptr)...,
147 : get<ValenciaTempTags>(temps_ptr)...,
148 :
149 : get<tmpl::conditional_t<
150 : tmpl::list_contains_v<extra_tags_list,
151 : Tags::detail::TemporaryReference<
152 : ValenciaTimeDerivativeArgTags>>,
153 : Tags::detail::TemporaryReference<ValenciaTimeDerivativeArgTags>,
154 : ValenciaTimeDerivativeArgTags>>(arguments, *temps_ptr)...);
155 :
156 : trace_reversed_stress_energy(
157 : get<TraceReversedStressResultTags>(temps_ptr)...,
158 : get<tmpl::conditional_t<
159 : tmpl::list_contains_v<extra_tags_list,
160 : Tags::detail::TemporaryReference<
161 : TraceReversedStressArgumentTags>>,
162 : Tags::detail::TemporaryReference<TraceReversedStressArgumentTags>,
163 : TraceReversedStressArgumentTags>>(*temps_ptr, arguments)...);
164 :
165 : // The addition to dt Pi is independent of the specific form of the stress
166 : // tensor.
167 : add_stress_energy_term_to_dt_pi(
168 : get<::Tags::dt<gh::Tags::Pi<DataVector, 3>>>(dt_vars_ptr),
169 : get<grmhd::GhValenciaDivClean::Tags::TraceReversedStressEnergy>(
170 : *temps_ptr),
171 : get<gr::Tags::Lapse<DataVector>>(*temps_ptr));
172 : }
173 : };
174 : } // namespace detail
175 :
176 : /*!
177 : * \brief Compute the RHS terms and flux values for both the Generalized
178 : * Harmonic formulation of Einstein's equations and the Valencia formulation of
179 : * the GRMHD equations with divergence cleaning.
180 : *
181 : * \details The bulk of the computations in this class dispatch to
182 : * `gh::TimeDerivative` and
183 : * `grmhd::ValenciaDivClean::TimeDerivativeTerms` as a 'product system' -- each
184 : * independently operating on its own subset of the supplied variable
185 : * collections.
186 : * The additional step is taken to compute the trace-reversed stress energy
187 : * tensor associated with the GRMHD part of the system and add its contribution
188 : * to the \f$\partial_t \Pi_{a b}\f$ variable in the Generalized Harmonic
189 : * system, which is the only explicit coupling required to back-react the effect
190 : * of matter on the spacetime solution.
191 : *
192 : * \note The MHD calculation reuses any spacetime quantities in its
193 : * argument_tags that are computed by the GH time derivative. However, other
194 : * quantities that aren't computed by the GH time derivative like the extrinsic
195 : * curvature are currently still retrieved from the DataBox. Those calculations
196 : * can be explicitly inlined here to reduce memory pressure and the number of
197 : * compute tags.
198 : */
199 1 : struct TimeDerivativeTerms : evolution::PassVariables {
200 0 : using gh_dt_tags =
201 : db::wrap_tags_in<::Tags::dt,
202 : typename gh::System<3_st>::variables_tag::tags_list>;
203 0 : using valencia_dt_tags = db::wrap_tags_in<
204 : ::Tags::dt,
205 : typename grmhd::ValenciaDivClean::System::variables_tag::tags_list>;
206 :
207 0 : using dt_tags = tmpl::append<gh_dt_tags, valencia_dt_tags>;
208 :
209 0 : using d_spatial_metric = ::Tags::deriv<gr::Tags::SpatialMetric<DataVector, 3>,
210 : tmpl::size_t<3>, Frame::Inertial>;
211 :
212 0 : using valencia_flux_tags = tmpl::transform<
213 : typename grmhd::ValenciaDivClean::System::flux_variables,
214 : tmpl::bind<::Tags::Flux, tmpl::_1, tmpl::pin<tmpl::size_t<3_st>>,
215 : tmpl::pin<Frame::Inertial>>>;
216 :
217 0 : using gh_temp_tags = typename gh::TimeDerivative<3_st>::temporary_tags;
218 0 : using gh_gradient_tags = typename gh::System<3_st>::gradients_tags;
219 0 : using gh_arg_tags = typename gh::TimeDerivative<3_st>::argument_tags;
220 :
221 0 : using valencia_temp_tags =
222 : typename grmhd::ValenciaDivClean::TimeDerivativeTerms::temporary_tags;
223 : // Additional temp tags are the derivatives of the metric since GH doesn't
224 : // explicitly calculate those.
225 0 : using valencia_extra_temp_tags =
226 : tmpl::list<::Tags::deriv<gr::Tags::Lapse<DataVector>, tmpl::size_t<3>,
227 : Frame::Inertial>,
228 : ::Tags::deriv<gr::Tags::Shift<DataVector, 3>, tmpl::size_t<3>,
229 : Frame::Inertial>,
230 : gr::Tags::ExtrinsicCurvature<DataVector, 3>>;
231 0 : using valencia_arg_tags = tmpl::list_difference<
232 : typename grmhd::ValenciaDivClean::TimeDerivativeTerms::argument_tags,
233 : tmpl::append<gh_temp_tags, valencia_extra_temp_tags>>;
234 :
235 0 : using trace_reversed_stress_result_tags =
236 : tmpl::list<Tags::TraceReversedStressEnergy, Tags::FourVelocityOneForm,
237 : Tags::ComovingMagneticFieldOneForm>;
238 0 : using trace_reversed_stress_argument_tags = tmpl::list<
239 : hydro::Tags::RestMassDensity<DataVector>,
240 : hydro::Tags::SpatialVelocityOneForm<DataVector, 3_st, Frame::Inertial>,
241 : hydro::Tags::MagneticFieldOneForm<DataVector, 3_st, Frame::Inertial>,
242 : hydro::Tags::MagneticFieldSquared<DataVector>,
243 : hydro::Tags::MagneticFieldDotSpatialVelocity<DataVector>,
244 : hydro::Tags::LorentzFactor<DataVector>,
245 : grmhd::ValenciaDivClean::TimeDerivativeTerms::OneOverLorentzFactorSquared,
246 : hydro::Tags::Pressure<DataVector>,
247 : hydro::Tags::SpecificInternalEnergy<DataVector>,
248 : gr::Tags::SpacetimeMetric<DataVector, 3>,
249 : gr::Tags::Shift<DataVector, 3_st>, gr::Tags::Lapse<DataVector>>;
250 0 : using extra_temp_tags = tmpl::list<gr::Tags::SpatialMetric<DataVector, 3>>;
251 :
252 0 : using temporary_tags = tmpl::remove<
253 : tmpl::remove_duplicates<tmpl::append<
254 : gh_temp_tags, valencia_temp_tags, valencia_extra_temp_tags,
255 : trace_reversed_stress_result_tags, extra_temp_tags>>,
256 : gr::Tags::SpatialMetric<DataVector, 3>>;
257 0 : using argument_tags =
258 : tmpl::remove<tmpl::remove<tmpl::append<gh_arg_tags,
259 :
260 : valencia_arg_tags>,
261 : gr::Tags::SpatialMetric<DataVector, 3>>,
262 : d_spatial_metric>;
263 :
264 : template <typename... Args>
265 0 : static void apply(
266 : const gsl::not_null<Variables<dt_tags>*> dt_vars_ptr,
267 : const gsl::not_null<Variables<db::wrap_tags_in<
268 : ::Tags::Flux, typename ValenciaDivClean::System::flux_variables,
269 : tmpl::size_t<3>, Frame::Inertial>>*>
270 : fluxes_ptr,
271 : const gsl::not_null<Variables<temporary_tags>*> temps_ptr,
272 : const tnsr::iaa<DataVector, 3>& d_spacetime_metric,
273 : const tnsr::iaa<DataVector, 3>& d_pi,
274 : const tnsr::ijaa<DataVector, 3>& d_phi, const Args&... args) {
275 : using args_list = tmpl::push_back<
276 : db::wrap_tags_in<Tags::detail::TemporaryReference, argument_tags>,
277 : gr::Tags::SpatialMetric<DataVector, 3>, d_spatial_metric>;
278 : tuples::tagged_tuple_from_typelist<args_list> arguments{
279 : args..., typename gr::Tags::SpatialMetric<DataVector, 3>::type{},
280 : typename d_spatial_metric::type{}};
281 : const size_t number_of_points = get<Tags::detail::TemporaryReference<
282 : gr::Tags::SpacetimeMetric<DataVector, 3>>>(arguments)[0]
283 : .size();
284 : for (size_t i = 0; i < 3; ++i) {
285 : for (size_t j = i; j < 3; ++j) {
286 : make_const_view(
287 : make_not_null(
288 : &std::as_const(
289 : get<gr::Tags::SpatialMetric<DataVector, 3>>(arguments))
290 : .get(i, j)),
291 : get<Tags::detail::TemporaryReference<
292 : gr::Tags::SpacetimeMetric<DataVector, 3>>>(arguments)
293 : .get(i + 1, j + 1),
294 : 0, number_of_points);
295 : }
296 : }
297 : for (size_t i = 0; i < 3; ++i) {
298 : for (size_t j = 0; j < 3; ++j) {
299 : for (size_t k = j; k < 3; ++k) {
300 : make_const_view(
301 : make_not_null(&std::as_const(get<d_spatial_metric>(arguments))
302 : .get(i, j, k)),
303 : get<Tags::detail::TemporaryReference<
304 : gh::Tags::Phi<DataVector, 3>>>(arguments)
305 : .get(i, j + 1, k + 1),
306 : 0, number_of_points);
307 : }
308 : }
309 : }
310 :
311 : detail::TimeDerivativeTermsImpl<
312 : gh_dt_tags, valencia_dt_tags, valencia_flux_tags, gh_temp_tags,
313 : valencia_temp_tags, gh_gradient_tags, gh_arg_tags, valencia_arg_tags,
314 : typename grmhd::ValenciaDivClean::TimeDerivativeTerms::argument_tags,
315 : trace_reversed_stress_result_tags,
316 : trace_reversed_stress_argument_tags>::apply(dt_vars_ptr, fluxes_ptr,
317 : temps_ptr,
318 : d_spacetime_metric, d_pi,
319 : d_phi, arguments);
320 : }
321 : };
322 : } // namespace grmhd::GhValenciaDivClean
|