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/Tensor/Tensor.hpp"
9 : #include "Elliptic/Systems/Xcts/Geometry.hpp"
10 : #include "Elliptic/Systems/Xcts/Tags.hpp"
11 : #include "NumericalAlgorithms/LinearOperators/Divergence.hpp"
12 : #include "NumericalAlgorithms/LinearOperators/PartialDerivatives.hpp"
13 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
14 : #include "PointwiseFunctions/GeneralRelativity/Tags/Conformal.hpp"
15 : #include "Utilities/Gsl.hpp"
16 : #include "Utilities/TMPL.hpp"
17 :
18 : /// \cond
19 : class DataVector;
20 : /// \endcond
21 :
22 : namespace Xcts {
23 :
24 : /// Indicates a subset of the XCTS equations
25 1 : enum class Equations {
26 : /// Only the Hamiltonian constraint, solved for \f$\psi\f$
27 : Hamiltonian,
28 : /// Both the Hamiltonian constraint and the lapse equation, solved for
29 : /// \f$\psi\f$ and \f$\alpha\psi\f$
30 : HamiltonianAndLapse,
31 : /// The full XCTS equations, solved for \f$\psi\f$, \f$\alpha\psi\f$ and
32 : /// \f$\beta_\mathrm{excess}\f$
33 : HamiltonianLapseAndShift
34 : };
35 :
36 : /// The fluxes \f$F^i\f$ for the first-order formulation of the XCTS equations.
37 : ///
38 : /// \see Xcts::FirstOrderSystem
39 : template <Equations EnabledEquations, Geometry ConformalGeometry>
40 1 : struct Fluxes;
41 :
42 : /// \cond
43 : template <>
44 : struct Fluxes<Equations::Hamiltonian, Geometry::FlatCartesian> {
45 : using argument_tags = tmpl::list<>;
46 : using volume_tags = tmpl::list<>;
47 : using const_global_cache_tags = tmpl::list<>;
48 : static constexpr bool is_trivial = true;
49 : static constexpr bool is_discontinuous = false;
50 : static void apply(
51 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
52 : const Scalar<DataVector>& conformal_factor_minus_one,
53 : const tnsr::i<DataVector, 3>& conformal_factor_gradient);
54 : static void apply(
55 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
56 : const tnsr::i<DataVector, 3>& face_normal,
57 : const tnsr::I<DataVector, 3>& face_normal_vector,
58 : const Scalar<DataVector>& conformal_factor_minus_one);
59 : };
60 :
61 : template <>
62 : struct Fluxes<Equations::Hamiltonian, Geometry::Curved> {
63 : using argument_tags =
64 : tmpl::list<Tags::InverseConformalMetric<DataVector, 3, Frame::Inertial>>;
65 : using volume_tags = tmpl::list<>;
66 : using const_global_cache_tags = tmpl::list<>;
67 : static constexpr bool is_trivial = true;
68 : static constexpr bool is_discontinuous = false;
69 : static void apply(
70 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
71 : const tnsr::II<DataVector, 3>& inv_conformal_metric,
72 : const Scalar<DataVector>& conformal_factor_minus_one,
73 : const tnsr::i<DataVector, 3>& conformal_factor_gradient);
74 : static void apply(
75 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
76 : const tnsr::II<DataVector, 3>& inv_conformal_metric,
77 : const tnsr::i<DataVector, 3>& face_normal,
78 : const tnsr::I<DataVector, 3>& face_normal_vector,
79 : const Scalar<DataVector>& conformal_factor_minus_one);
80 : };
81 :
82 : template <>
83 : struct Fluxes<Equations::HamiltonianAndLapse, Geometry::FlatCartesian> {
84 : using argument_tags = tmpl::list<>;
85 : using volume_tags = tmpl::list<>;
86 : using const_global_cache_tags = tmpl::list<>;
87 : static constexpr bool is_trivial = true;
88 : static constexpr bool is_discontinuous = false;
89 : static void apply(
90 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
91 : gsl::not_null<tnsr::I<DataVector, 3>*>
92 : flux_for_lapse_times_conformal_factor,
93 : const Scalar<DataVector>& conformal_factor_minus_one,
94 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
95 : const tnsr::i<DataVector, 3>& conformal_factor_gradient,
96 : const tnsr::i<DataVector, 3>& lapse_times_conformal_factor_gradient);
97 : static void apply(
98 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
99 : gsl::not_null<tnsr::I<DataVector, 3>*>
100 : flux_for_lapse_times_conformal_factor,
101 : const tnsr::i<DataVector, 3>& face_normal,
102 : const tnsr::I<DataVector, 3>& face_normal_vector,
103 : const Scalar<DataVector>& conformal_factor_minus_one,
104 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one);
105 : };
106 :
107 : template <>
108 : struct Fluxes<Equations::HamiltonianAndLapse, Geometry::Curved> {
109 : using argument_tags =
110 : tmpl::list<Tags::InverseConformalMetric<DataVector, 3, Frame::Inertial>>;
111 : using volume_tags = tmpl::list<>;
112 : using const_global_cache_tags = tmpl::list<>;
113 : static constexpr bool is_trivial = true;
114 : static constexpr bool is_discontinuous = false;
115 : static void apply(
116 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
117 : gsl::not_null<tnsr::I<DataVector, 3>*>
118 : flux_for_lapse_times_conformal_factor,
119 : const tnsr::II<DataVector, 3>& inv_conformal_metric,
120 : const Scalar<DataVector>& conformal_factor_minus_one,
121 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
122 : const tnsr::i<DataVector, 3>& conformal_factor_gradient,
123 : const tnsr::i<DataVector, 3>& lapse_times_conformal_factor_gradient);
124 : static void apply(
125 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
126 : gsl::not_null<tnsr::I<DataVector, 3>*>
127 : flux_for_lapse_times_conformal_factor,
128 : const tnsr::II<DataVector, 3>& inv_conformal_metric,
129 : const tnsr::i<DataVector, 3>& face_normal,
130 : const tnsr::I<DataVector, 3>& face_normal_vector,
131 : const Scalar<DataVector>& conformal_factor_minus_one,
132 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one);
133 : };
134 :
135 : template <>
136 : struct Fluxes<Equations::HamiltonianLapseAndShift, Geometry::FlatCartesian> {
137 : using argument_tags = tmpl::list<>;
138 : using volume_tags = tmpl::list<>;
139 : using const_global_cache_tags = tmpl::list<>;
140 : static constexpr bool is_trivial = false;
141 : static constexpr bool is_discontinuous = false;
142 : static void apply(
143 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
144 : gsl::not_null<tnsr::I<DataVector, 3>*>
145 : flux_for_lapse_times_conformal_factor,
146 : gsl::not_null<tnsr::II<DataVector, 3>*> longitudinal_shift_excess,
147 : const Scalar<DataVector>& conformal_factor_minus_one,
148 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
149 : const tnsr::I<DataVector, 3>& shift_excess,
150 : const tnsr::i<DataVector, 3>& conformal_factor_gradient,
151 : const tnsr::i<DataVector, 3>& lapse_times_conformal_factor_gradient,
152 : const tnsr::iJ<DataVector, 3>& deriv_shift_excess);
153 : static void apply(
154 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
155 : gsl::not_null<tnsr::I<DataVector, 3>*>
156 : flux_for_lapse_times_conformal_factor,
157 : gsl::not_null<tnsr::II<DataVector, 3>*> longitudinal_shift_excess,
158 : const tnsr::i<DataVector, 3>& face_normal,
159 : const tnsr::I<DataVector, 3>& face_normal_vector,
160 : const Scalar<DataVector>& conformal_factor_minus_one,
161 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
162 : const tnsr::I<DataVector, 3>& shift_excess);
163 : };
164 :
165 : template <>
166 : struct Fluxes<Equations::HamiltonianLapseAndShift, Geometry::Curved> {
167 : using argument_tags = tmpl::list<
168 : Tags::ConformalMetric<DataVector, 3, Frame::Inertial>,
169 : Tags::InverseConformalMetric<DataVector, 3, Frame::Inertial>,
170 : Tags::ConformalChristoffelSecondKind<DataVector, 3, Frame::Inertial>>;
171 : using volume_tags = tmpl::list<>;
172 : using const_global_cache_tags = tmpl::list<>;
173 : static constexpr bool is_trivial = false;
174 : static constexpr bool is_discontinuous = false;
175 : static void apply(
176 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
177 : gsl::not_null<tnsr::I<DataVector, 3>*>
178 : flux_for_lapse_times_conformal_factor,
179 : gsl::not_null<tnsr::II<DataVector, 3>*> longitudinal_shift_excess,
180 : const tnsr::ii<DataVector, 3>& /*conformal_metric*/,
181 : const tnsr::II<DataVector, 3>& inv_conformal_metric,
182 : const tnsr::Ijj<DataVector, 3>& christoffel_second_kind,
183 : const Scalar<DataVector>& conformal_factor_minus_one,
184 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
185 : const tnsr::I<DataVector, 3>& shift_excess,
186 : const tnsr::i<DataVector, 3>& conformal_factor_gradient,
187 : const tnsr::i<DataVector, 3>& lapse_times_conformal_factor_gradient,
188 : const tnsr::iJ<DataVector, 3>& deriv_shift_excess);
189 : static void apply(
190 : gsl::not_null<tnsr::I<DataVector, 3>*> flux_for_conformal_factor,
191 : gsl::not_null<tnsr::I<DataVector, 3>*>
192 : flux_for_lapse_times_conformal_factor,
193 : gsl::not_null<tnsr::II<DataVector, 3>*> longitudinal_shift_excess,
194 : const tnsr::ii<DataVector, 3>& /*conformal_metric*/,
195 : const tnsr::II<DataVector, 3>& inv_conformal_metric,
196 : const tnsr::Ijj<DataVector, 3>& christoffel_second_kind,
197 : const tnsr::i<DataVector, 3>& face_normal,
198 : const tnsr::I<DataVector, 3>& face_normal_vector,
199 : const Scalar<DataVector>& conformal_factor_minus_one,
200 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
201 : const tnsr::I<DataVector, 3>& shift_excess);
202 : };
203 : /// \endcond
204 :
205 : /// The sources \f$S\f$ for the first-order formulation of the XCTS equations.
206 : ///
207 : /// \see Xcts::FirstOrderSystem
208 : template <Equations EnabledEquations, Geometry ConformalGeometry,
209 : int ConformalMatterScale>
210 1 : struct Sources;
211 :
212 : /// \cond
213 : template <int ConformalMatterScale>
214 : struct Sources<Equations::Hamiltonian, Geometry::FlatCartesian,
215 : ConformalMatterScale> {
216 : using argument_tags = tmpl::list<
217 : gr::Tags::Conformal<gr::Tags::EnergyDensity<DataVector>,
218 : ConformalMatterScale>,
219 : gr::Tags::TraceExtrinsicCurvature<DataVector>,
220 : Tags::LongitudinalShiftMinusDtConformalMetricOverLapseSquare<DataVector>>;
221 : static void apply(
222 : gsl::not_null<Scalar<DataVector>*> hamiltonian_constraint,
223 : const Scalar<DataVector>& conformal_energy_density,
224 : const Scalar<DataVector>& extrinsic_curvature_trace,
225 : const Scalar<DataVector>&
226 : longitudinal_shift_minus_dt_conformal_metric_over_lapse_square,
227 : const Scalar<DataVector>& conformal_factor_minus_one,
228 : const tnsr::I<DataVector, 3>& conformal_factor_flux);
229 : };
230 :
231 : template <int ConformalMatterScale>
232 : struct Sources<Equations::Hamiltonian, Geometry::Curved, ConformalMatterScale> {
233 : using argument_tags = tmpl::push_back<
234 : typename Sources<Equations::Hamiltonian, Geometry::FlatCartesian,
235 : ConformalMatterScale>::argument_tags,
236 : Tags::ConformalChristoffelContracted<DataVector, 3, Frame::Inertial>,
237 : Tags::ConformalRicciScalar<DataVector>>;
238 : static void apply(
239 : gsl::not_null<Scalar<DataVector>*> hamiltonian_constraint,
240 : const Scalar<DataVector>& conformal_energy_density,
241 : const Scalar<DataVector>& extrinsic_curvature_trace,
242 : const Scalar<DataVector>&
243 : longitudinal_shift_minus_dt_conformal_metric_over_lapse_square,
244 : const tnsr::i<DataVector, 3>& conformal_christoffel_contracted,
245 : const Scalar<DataVector>& conformal_ricci_scalar,
246 : const Scalar<DataVector>& conformal_factor_minus_one,
247 : const tnsr::I<DataVector, 3>& conformal_factor_flux);
248 : };
249 :
250 : template <int ConformalMatterScale>
251 : struct Sources<Equations::HamiltonianAndLapse, Geometry::FlatCartesian,
252 : ConformalMatterScale> {
253 : using argument_tags = tmpl::list<
254 : gr::Tags::Conformal<gr::Tags::EnergyDensity<DataVector>,
255 : ConformalMatterScale>,
256 : gr::Tags::Conformal<gr::Tags::StressTrace<DataVector>,
257 : ConformalMatterScale>,
258 : gr::Tags::TraceExtrinsicCurvature<DataVector>,
259 : ::Tags::dt<gr::Tags::TraceExtrinsicCurvature<DataVector>>,
260 : Tags::LongitudinalShiftMinusDtConformalMetricSquare<DataVector>,
261 : Tags::ShiftDotDerivExtrinsicCurvatureTrace<DataVector>>;
262 : static void apply(
263 : gsl::not_null<Scalar<DataVector>*> hamiltonian_constraint,
264 : gsl::not_null<Scalar<DataVector>*> lapse_equation,
265 : const Scalar<DataVector>& conformal_energy_density,
266 : const Scalar<DataVector>& conformal_stress_trace,
267 : const Scalar<DataVector>& extrinsic_curvature_trace,
268 : const Scalar<DataVector>& dt_extrinsic_curvature_trace,
269 : const Scalar<DataVector>&
270 : longitudinal_shift_minus_dt_conformal_metric_square,
271 : const Scalar<DataVector>& shift_dot_deriv_extrinsic_curvature_trace,
272 : const Scalar<DataVector>& conformal_factor_minus_one,
273 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
274 : const tnsr::I<DataVector, 3>& conformal_factor_flux,
275 : const tnsr::I<DataVector, 3>& lapse_times_conformal_factor_flux);
276 : };
277 :
278 : template <int ConformalMatterScale>
279 : struct Sources<Equations::HamiltonianAndLapse, Geometry::Curved,
280 : ConformalMatterScale> {
281 : using argument_tags = tmpl::push_back<
282 : typename Sources<Equations::HamiltonianAndLapse, Geometry::FlatCartesian,
283 : ConformalMatterScale>::argument_tags,
284 : Tags::ConformalChristoffelContracted<DataVector, 3, Frame::Inertial>,
285 : Tags::ConformalRicciScalar<DataVector>>;
286 : static void apply(
287 : gsl::not_null<Scalar<DataVector>*> hamiltonian_constraint,
288 : gsl::not_null<Scalar<DataVector>*> lapse_equation,
289 : const Scalar<DataVector>& conformal_energy_density,
290 : const Scalar<DataVector>& conformal_stress_trace,
291 : const Scalar<DataVector>& extrinsic_curvature_trace,
292 : const Scalar<DataVector>& dt_extrinsic_curvature_trace,
293 : const Scalar<DataVector>&
294 : longitudinal_shift_minus_dt_conformal_metric_square,
295 : const Scalar<DataVector>& shift_dot_deriv_extrinsic_curvature_trace,
296 : const tnsr::i<DataVector, 3>& conformal_christoffel_contracted,
297 : const Scalar<DataVector>& conformal_ricci_scalar,
298 : const Scalar<DataVector>& conformal_factor_minus_one,
299 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
300 : const tnsr::I<DataVector, 3>& conformal_factor_flux,
301 : const tnsr::I<DataVector, 3>& lapse_times_conformal_factor_flux);
302 : };
303 :
304 : template <int ConformalMatterScale>
305 : struct Sources<Equations::HamiltonianLapseAndShift, Geometry::FlatCartesian,
306 : ConformalMatterScale> {
307 : using argument_tags = tmpl::list<
308 : gr::Tags::Conformal<gr::Tags::EnergyDensity<DataVector>,
309 : ConformalMatterScale>,
310 : gr::Tags::Conformal<gr::Tags::StressTrace<DataVector>,
311 : ConformalMatterScale>,
312 : gr::Tags::Conformal<gr::Tags::MomentumDensity<DataVector, 3>,
313 : ConformalMatterScale>,
314 : gr::Tags::TraceExtrinsicCurvature<DataVector>,
315 : ::Tags::dt<gr::Tags::TraceExtrinsicCurvature<DataVector>>,
316 : ::Tags::deriv<gr::Tags::TraceExtrinsicCurvature<DataVector>,
317 : tmpl::size_t<3>, Frame::Inertial>,
318 : Tags::ShiftBackground<DataVector, 3, Frame::Inertial>,
319 : Tags::LongitudinalShiftBackgroundMinusDtConformalMetric<DataVector, 3,
320 : Frame::Inertial>,
321 : ::Tags::div<Tags::LongitudinalShiftBackgroundMinusDtConformalMetric<
322 : DataVector, 3, Frame::Inertial>>>;
323 : static void apply(
324 : gsl::not_null<Scalar<DataVector>*> hamiltonian_constraint,
325 : gsl::not_null<Scalar<DataVector>*> lapse_equation,
326 : gsl::not_null<tnsr::I<DataVector, 3>*> momentum_constraint,
327 : const Scalar<DataVector>& conformal_energy_density,
328 : const Scalar<DataVector>& conformal_stress_trace,
329 : const tnsr::I<DataVector, 3>& conformal_momentum_density,
330 : const Scalar<DataVector>& extrinsic_curvature_trace,
331 : const Scalar<DataVector>& dt_extrinsic_curvature_trace,
332 : const tnsr::i<DataVector, 3>& extrinsic_curvature_trace_gradient,
333 : const tnsr::I<DataVector, 3>& shift_background,
334 : const tnsr::II<DataVector, 3>&
335 : longitudinal_shift_background_minus_dt_conformal_metric,
336 : const tnsr::I<DataVector, 3>&
337 : div_longitudinal_shift_background_minus_dt_conformal_metric,
338 : const Scalar<DataVector>& conformal_factor_minus_one,
339 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
340 : const tnsr::I<DataVector, 3>& shift_excess,
341 : const tnsr::I<DataVector, 3>& conformal_factor_flux,
342 : const tnsr::I<DataVector, 3>& lapse_times_conformal_factor_flux,
343 : const tnsr::II<DataVector, 3>& longitudinal_shift_excess);
344 : };
345 :
346 : template <int ConformalMatterScale>
347 : struct Sources<Equations::HamiltonianLapseAndShift, Geometry::Curved,
348 : ConformalMatterScale> {
349 : using argument_tags = tmpl::push_back<
350 : typename Sources<Equations::HamiltonianLapseAndShift,
351 : Geometry::FlatCartesian,
352 : ConformalMatterScale>::argument_tags,
353 : Tags::ConformalMetric<DataVector, 3, Frame::Inertial>,
354 : Tags::InverseConformalMetric<DataVector, 3, Frame::Inertial>,
355 : Tags::ConformalChristoffelFirstKind<DataVector, 3, Frame::Inertial>,
356 : Tags::ConformalChristoffelSecondKind<DataVector, 3, Frame::Inertial>,
357 : Tags::ConformalChristoffelContracted<DataVector, 3, Frame::Inertial>,
358 : Tags::ConformalRicciScalar<DataVector>>;
359 : static void apply(
360 : gsl::not_null<Scalar<DataVector>*> hamiltonian_constraint,
361 : gsl::not_null<Scalar<DataVector>*> lapse_equation,
362 : gsl::not_null<tnsr::I<DataVector, 3>*> momentum_constraint,
363 : const Scalar<DataVector>& conformal_energy_density,
364 : const Scalar<DataVector>& conformal_stress_trace,
365 : const tnsr::I<DataVector, 3>& conformal_momentum_density,
366 : const Scalar<DataVector>& extrinsic_curvature_trace,
367 : const Scalar<DataVector>& dt_extrinsic_curvature_trace,
368 : const tnsr::i<DataVector, 3>& extrinsic_curvature_trace_gradient,
369 : const tnsr::I<DataVector, 3>& shift_background,
370 : const tnsr::II<DataVector, 3>&
371 : longitudinal_shift_background_minus_dt_conformal_metric,
372 : const tnsr::I<DataVector, 3>&
373 : div_longitudinal_shift_background_minus_dt_conformal_metric,
374 : const tnsr::ii<DataVector, 3>& conformal_metric,
375 : const tnsr::II<DataVector, 3>& inv_conformal_metric,
376 : const tnsr::ijj<DataVector, 3>& /*conformal_christoffel_first_kind*/,
377 : const tnsr::Ijj<DataVector, 3>& conformal_christoffel_second_kind,
378 : const tnsr::i<DataVector, 3>& conformal_christoffel_contracted,
379 : const Scalar<DataVector>& conformal_ricci_scalar,
380 : const Scalar<DataVector>& conformal_factor_minus_one,
381 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
382 : const tnsr::I<DataVector, 3>& shift_excess,
383 : const tnsr::I<DataVector, 3>& conformal_factor_flux,
384 : const tnsr::I<DataVector, 3>& lapse_times_conformal_factor_flux,
385 : const tnsr::II<DataVector, 3>& longitudinal_shift_excess);
386 : };
387 : /// \endcond
388 :
389 : /// The linearization of the sources \f$S\f$ for the first-order formulation of
390 : /// the XCTS equations.
391 : ///
392 : /// \see Xcts::FirstOrderSystem
393 : template <Equations EnabledEquations, Geometry ConformalGeometry,
394 : int ConformalMatterScale>
395 1 : struct LinearizedSources;
396 :
397 : /// \cond
398 : template <int ConformalMatterScale>
399 : struct LinearizedSources<Equations::Hamiltonian, Geometry::FlatCartesian,
400 : ConformalMatterScale> {
401 : using argument_tags = tmpl::push_back<
402 : typename Sources<Equations::Hamiltonian, Geometry::FlatCartesian,
403 : ConformalMatterScale>::argument_tags,
404 : Tags::ConformalFactorMinusOne<DataVector>>;
405 : static void apply(
406 : gsl::not_null<Scalar<DataVector>*> linearized_hamiltonian_constraint,
407 : const Scalar<DataVector>& conformal_energy_density,
408 : const Scalar<DataVector>& extrinsic_curvature_trace,
409 : const Scalar<DataVector>&
410 : longitudinal_shift_minus_dt_conformal_metric_over_lapse_square,
411 : const Scalar<DataVector>& conformal_factor_minus_one,
412 : const Scalar<DataVector>& conformal_factor_correction,
413 : const tnsr::I<DataVector, 3>&
414 : /*conformal_factor_flux_correction*/);
415 : };
416 :
417 : template <int ConformalMatterScale>
418 : struct LinearizedSources<Equations::Hamiltonian, Geometry::Curved,
419 : ConformalMatterScale> {
420 : using argument_tags =
421 : tmpl::push_back<typename Sources<Equations::Hamiltonian, Geometry::Curved,
422 : ConformalMatterScale>::argument_tags,
423 : Tags::ConformalFactorMinusOne<DataVector>>;
424 : static void apply(
425 : gsl::not_null<Scalar<DataVector>*> linearized_hamiltonian_constraint,
426 : const Scalar<DataVector>& conformal_energy_density,
427 : const Scalar<DataVector>& extrinsic_curvature_trace,
428 : const Scalar<DataVector>&
429 : longitudinal_shift_minus_dt_conformal_metric_over_lapse_square,
430 : const tnsr::i<DataVector, 3>& conformal_christoffel_contracted,
431 : const Scalar<DataVector>& conformal_ricci_scalar,
432 : const Scalar<DataVector>& conformal_factor_minus_one,
433 : const Scalar<DataVector>& conformal_factor_correction,
434 : const tnsr::I<DataVector, 3>& conformal_factor_flux_correction);
435 : };
436 :
437 : template <int ConformalMatterScale>
438 : struct LinearizedSources<Equations::HamiltonianAndLapse,
439 : Geometry::FlatCartesian, ConformalMatterScale> {
440 : using argument_tags = tmpl::push_back<
441 : typename Sources<Equations::HamiltonianAndLapse, Geometry::FlatCartesian,
442 : ConformalMatterScale>::argument_tags,
443 : Tags::ConformalFactorMinusOne<DataVector>,
444 : Tags::LapseTimesConformalFactorMinusOne<DataVector>>;
445 : static void apply(
446 : gsl::not_null<Scalar<DataVector>*> linearized_hamiltonian_constraint,
447 : gsl::not_null<Scalar<DataVector>*> linearized_lapse_equation,
448 : const Scalar<DataVector>& conformal_energy_density,
449 : const Scalar<DataVector>& conformal_stress_trace,
450 : const Scalar<DataVector>& extrinsic_curvature_trace,
451 : const Scalar<DataVector>& dt_extrinsic_curvature_trace,
452 : const Scalar<DataVector>&
453 : longitudinal_shift_minus_dt_conformal_metric_square,
454 : const Scalar<DataVector>& shift_dot_deriv_extrinsic_curvature_trace,
455 : const Scalar<DataVector>& conformal_factor_minus_one,
456 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
457 : const Scalar<DataVector>& conformal_factor_correction,
458 : const Scalar<DataVector>& lapse_times_conformal_factor_correction,
459 : const tnsr::I<DataVector, 3>& /*conformal_factor_flux_correction*/,
460 : const tnsr::I<DataVector, 3>&
461 : lapse_times_conformal_factor_flux_correction);
462 : };
463 :
464 : template <int ConformalMatterScale>
465 : struct LinearizedSources<Equations::HamiltonianAndLapse, Geometry::Curved,
466 : ConformalMatterScale> {
467 : using argument_tags = tmpl::push_back<
468 : typename Sources<Equations::HamiltonianAndLapse, Geometry::Curved,
469 : ConformalMatterScale>::argument_tags,
470 : Tags::ConformalFactorMinusOne<DataVector>,
471 : Tags::LapseTimesConformalFactorMinusOne<DataVector>>;
472 : static void apply(
473 : gsl::not_null<Scalar<DataVector>*> linearized_hamiltonian_constraint,
474 : gsl::not_null<Scalar<DataVector>*> linearized_lapse_equation,
475 : const Scalar<DataVector>& conformal_energy_density,
476 : const Scalar<DataVector>& conformal_stress_trace,
477 : const Scalar<DataVector>& extrinsic_curvature_trace,
478 : const Scalar<DataVector>& dt_extrinsic_curvature_trace,
479 : const Scalar<DataVector>&
480 : longitudinal_shift_minus_dt_conformal_metric_square,
481 : const Scalar<DataVector>& shift_dot_deriv_extrinsic_curvature_trace,
482 : const tnsr::i<DataVector, 3>& conformal_christoffel_contracted,
483 : const Scalar<DataVector>& conformal_ricci_scalar,
484 : const Scalar<DataVector>& conformal_factor_minus_one,
485 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
486 : const Scalar<DataVector>& conformal_factor_correction,
487 : const Scalar<DataVector>& lapse_times_conformal_factor_correction,
488 : const tnsr::I<DataVector, 3>& conformal_factor_flux_correction,
489 : const tnsr::I<DataVector, 3>&
490 : lapse_times_conformal_factor_flux_correction);
491 : };
492 :
493 : template <int ConformalMatterScale>
494 : struct LinearizedSources<Equations::HamiltonianLapseAndShift,
495 : Geometry::FlatCartesian, ConformalMatterScale> {
496 : using argument_tags = tmpl::push_back<
497 : typename Sources<Equations::HamiltonianLapseAndShift,
498 : Geometry::FlatCartesian,
499 : ConformalMatterScale>::argument_tags,
500 : Tags::ConformalFactorMinusOne<DataVector>,
501 : Tags::LapseTimesConformalFactorMinusOne<DataVector>,
502 : Tags::ShiftExcess<DataVector, 3, Frame::Inertial>,
503 : ::Tags::Flux<Tags::ConformalFactorMinusOne<DataVector>, tmpl::size_t<3>,
504 : Frame::Inertial>,
505 : ::Tags::Flux<Tags::LapseTimesConformalFactorMinusOne<DataVector>,
506 : tmpl::size_t<3>, Frame::Inertial>,
507 : Tags::LongitudinalShiftExcess<DataVector, 3, Frame::Inertial>>;
508 : static void apply(
509 : gsl::not_null<Scalar<DataVector>*> linearized_hamiltonian_constraint,
510 : gsl::not_null<Scalar<DataVector>*> linearized_lapse_equation,
511 : gsl::not_null<tnsr::I<DataVector, 3>*> linearized_momentum_constraint,
512 : const Scalar<DataVector>& conformal_energy_density,
513 : const Scalar<DataVector>& conformal_stress_trace,
514 : const tnsr::I<DataVector, 3>& conformal_momentum_density,
515 : const Scalar<DataVector>& extrinsic_curvature_trace,
516 : const Scalar<DataVector>& dt_extrinsic_curvature_trace,
517 : const tnsr::i<DataVector, 3>& extrinsic_curvature_trace_gradient,
518 : const tnsr::I<DataVector, 3>& shift_background,
519 : const tnsr::II<DataVector, 3>&
520 : longitudinal_shift_background_minus_dt_conformal_metric,
521 : const tnsr::I<DataVector, 3>&
522 : div_longitudinal_shift_background_minus_dt_conformal_metric,
523 : const Scalar<DataVector>& conformal_factor_minus_one,
524 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
525 : const tnsr::I<DataVector, 3>& shift_excess,
526 : const tnsr::I<DataVector, 3>& conformal_factor_flux,
527 : const tnsr::I<DataVector, 3>& lapse_times_conformal_factor_flux,
528 : const tnsr::II<DataVector, 3>& longitudinal_shift_excess,
529 : const Scalar<DataVector>& conformal_factor_correction,
530 : const Scalar<DataVector>& lapse_times_conformal_factor_correction,
531 : const tnsr::I<DataVector, 3>& shift_excess_correction,
532 : const tnsr::I<DataVector, 3>& conformal_factor_flux_correction,
533 : const tnsr::I<DataVector, 3>&
534 : lapse_times_conformal_factor_flux_correction,
535 : const tnsr::II<DataVector, 3>& longitudinal_shift_excess_correction);
536 : };
537 :
538 : template <int ConformalMatterScale>
539 : struct LinearizedSources<Equations::HamiltonianLapseAndShift, Geometry::Curved,
540 : ConformalMatterScale> {
541 : using argument_tags = tmpl::push_back<
542 : typename Sources<Equations::HamiltonianLapseAndShift, Geometry::Curved,
543 : ConformalMatterScale>::argument_tags,
544 : Tags::ConformalFactorMinusOne<DataVector>,
545 : Tags::LapseTimesConformalFactorMinusOne<DataVector>,
546 : Tags::ShiftExcess<DataVector, 3, Frame::Inertial>,
547 : ::Tags::Flux<Tags::ConformalFactorMinusOne<DataVector>, tmpl::size_t<3>,
548 : Frame::Inertial>,
549 : ::Tags::Flux<Tags::LapseTimesConformalFactorMinusOne<DataVector>,
550 : tmpl::size_t<3>, Frame::Inertial>,
551 : Tags::LongitudinalShiftExcess<DataVector, 3, Frame::Inertial>>;
552 : static void apply(
553 : gsl::not_null<Scalar<DataVector>*> linearized_hamiltonian_constraint,
554 : gsl::not_null<Scalar<DataVector>*> linearized_lapse_equation,
555 : gsl::not_null<tnsr::I<DataVector, 3>*> linearized_momentum_constraint,
556 : const Scalar<DataVector>& conformal_energy_density,
557 : const Scalar<DataVector>& conformal_stress_trace,
558 : const tnsr::I<DataVector, 3>& conformal_momentum_density,
559 : const Scalar<DataVector>& extrinsic_curvature_trace,
560 : const Scalar<DataVector>& dt_extrinsic_curvature_trace,
561 : const tnsr::i<DataVector, 3>& extrinsic_curvature_trace_gradient,
562 : const tnsr::I<DataVector, 3>& shift_background,
563 : const tnsr::II<DataVector, 3>&
564 : longitudinal_shift_background_minus_dt_conformal_metric,
565 : const tnsr::I<DataVector, 3>&
566 : /*div_longitudinal_shift_background_minus_dt_conformal_metric*/,
567 : const tnsr::ii<DataVector, 3>& conformal_metric,
568 : const tnsr::II<DataVector, 3>& inv_conformal_metric,
569 : const tnsr::ijj<DataVector, 3>& /*conformal_christoffel_first_kind*/,
570 : const tnsr::Ijj<DataVector, 3>& conformal_christoffel_second_kind,
571 : const tnsr::i<DataVector, 3>& conformal_christoffel_contracted,
572 : const Scalar<DataVector>& conformal_ricci_scalar,
573 : const Scalar<DataVector>& conformal_factor_minus_one,
574 : const Scalar<DataVector>& lapse_times_conformal_factor_minus_one,
575 : const tnsr::I<DataVector, 3>& shift_excess,
576 : const tnsr::I<DataVector, 3>& conformal_factor_flux,
577 : const tnsr::I<DataVector, 3>& lapse_times_conformal_factor_flux,
578 : const tnsr::II<DataVector, 3>& longitudinal_shift_excess,
579 : const Scalar<DataVector>& conformal_factor_correction,
580 : const Scalar<DataVector>& lapse_times_conformal_factor_correction,
581 : const tnsr::I<DataVector, 3>& shift_excess_correction,
582 : const tnsr::I<DataVector, 3>& conformal_factor_flux_correction,
583 : const tnsr::I<DataVector, 3>&
584 : lapse_times_conformal_factor_flux_correction,
585 : const tnsr::II<DataVector, 3>& longitudinal_shift_excess_correction);
586 : };
587 : /// \endcond
588 :
589 : } // namespace Xcts
|