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 <limits>
8 : #include <pup.h>
9 :
10 : #include "DataStructures/CachedTempBuffer.hpp"
11 : #include "DataStructures/Tags/TempTensor.hpp"
12 : #include "DataStructures/Tensor/IndexType.hpp"
13 : #include "DataStructures/Tensor/TypeAliases.hpp"
14 : #include "Options/Context.hpp"
15 : #include "Options/String.hpp"
16 : #include "PointwiseFunctions/AnalyticSolutions/AnalyticSolution.hpp"
17 : #include "PointwiseFunctions/AnalyticSolutions/GeneralRelativity/Solutions.hpp"
18 : #include "PointwiseFunctions/GeneralRelativity/TagsDeclarations.hpp"
19 : #include "Utilities/ContainerHelpers.hpp"
20 : #include "Utilities/ForceInline.hpp"
21 : #include "Utilities/Literals.hpp"
22 : #include "Utilities/MakeArray.hpp"
23 : #include "Utilities/TMPL.hpp"
24 : #include "Utilities/TaggedTuple.hpp"
25 :
26 : /// \cond
27 : namespace Tags {
28 : template <typename Tag>
29 : struct dt;
30 : } // namespace Tags
31 : namespace gsl {
32 : template <class T>
33 : class not_null;
34 : } // namespace gsl
35 : /// \endcond
36 :
37 : namespace gr {
38 : namespace Solutions {
39 : /*!
40 : * \brief Schwarzschild black hole in Cartesian coordinates with harmonic gauge
41 : *
42 : * \details
43 : * This solution represents a Schwarzschild black hole in coordinates that are
44 : * harmonic in both space and time, as well as horizon-penetrating. Therefore,
45 : * this solution fulfills the harmonic coordinate conditions Eq. (4.42), (4.44),
46 : * and (4.45) in \cite BaumgarteShapiro :
47 : *
48 : * \f{align}
49 : * {}^{(4)}\Gamma^a &= 0
50 : * \f}
51 : *
52 : * \f{align}
53 : * (\partial_t - \beta^j \partial_j)\alpha &= -\alpha^2 K
54 : * \f}
55 : *
56 : * \f{align}
57 : * (\partial_t - \beta^j \partial_j)\beta^i &=
58 : * -\alpha^2 (\gamma^{ij} \partial_j \ln \alpha -
59 : * \gamma^{jk} \Gamma^i_{jk})
60 : * \f}
61 : *
62 : * (Note that Eq. 4.45 in \cite BaumgarteShapiro is missing a minus sign in
63 : * front of the \f$\Gamma\f$-contraction term)
64 : *
65 : * We implement Eqs. (45)--(50) in \cite Cook1997qc , which represent the
66 : * zero-spin limit of the time-harmonic and horizon-penetrating slices of Kerr
67 : * spacetime presented in the paper. We add the radial transformation
68 : * \f$r \to r + M\f$ to make the spatial coordinates harmonic as well (see
69 : * Eq. (43) in \cite Cook1997qc ), so the coordinates remain harmonic under
70 : * boosts.
71 : *
72 : * Consider a Schwarzschild black hole of mass \f$M\f$ and center \f$C^i\f$. The
73 : * spacetime will be specified using Cartesian coordinates \f$x^i\f$ that are
74 : * spatially and temporally harmonic. A radius centered on the black hole is
75 : *
76 : * \f{align}
77 : * r &= \sqrt{\delta_{ij} \left(x^i - C^i\right)\left(x^j - C^j\right)}
78 : * \f}
79 : *
80 : * For computing the spatial metric, we define the following quantities:
81 : *
82 : * \f{align}
83 : * \gamma_{rr} &= 1 + \frac{2M}{M+r} + \left(\frac{2M}{M+r}\right)^2
84 : * + \left(\frac{2M}{M+r}\right)^3,\\
85 : * \partial_r \gamma_{rr} &= -\frac{1}{2M}\left(\frac{2M}{M+r}\right)^2
86 : * -\frac{1}{M}\left(\frac{2M}{M+r}\right)^3
87 : * -\frac{3}{2M}\left(\frac{2M}{M+r}\right)^4,\\
88 : * \frac{X^i}{r} &= \frac{x^i - C^i}{r},\\
89 : * X_j &= X^i \delta_{ij}
90 : * \f}
91 : *
92 : * From these quantities, the spatial metric and its time derivative are
93 : * computed as
94 : *
95 : * \f{align}
96 : * \gamma_{ij} &=
97 : * \left(\gamma_{rr} - \left(1+\frac{M}{r}\right)^2\right)
98 : * \frac{X_i}{r} \frac{X_j}{r} +
99 : * \delta_{ij} \left(1+\frac{M}{r}\right)^2,\\
100 : * \partial_t \gamma_{ij} &= 0
101 : * \f}
102 : *
103 : * The spatial derivative is given in terms of the following quantities:
104 : *
105 : * \f{align}
106 : * f_0 &= \left(1+\frac{M}{r}\right)^2\\
107 : * \partial_r f_0 &=
108 : * 2 \left(1+\frac{M}{r}\right)\left(-\frac{M}{r^2}\right),\\
109 : * f_1 &= \frac{1}{r} \left(\gamma_{rr} - f_0\right),\\
110 : * f_2 &= \partial_r \gamma_{rr} - \partial_r f_0 - 2 f_1
111 : * \f}
112 : *
113 : * In terms of these, the spatial metric's spatial derivative is
114 : *
115 : * \f{align}
116 : * \partial_k \gamma_{ij} &=
117 : * f_2 \frac{X_i}{r} \frac{X_j}{r} \frac{X_k}{r} +
118 : * f_1 \frac{X_j}{r} \delta_{ik} + f_1 \frac{X_i}{r} \delta_{jk} +
119 : * \partial_r f_0 \frac{X_k}{r} \delta_{ij}
120 : * \f}
121 : *
122 : * The lapse and its derivatives are
123 : *
124 : * \f{align}
125 : * \alpha &= \gamma_{rr}^{-1/2},\\
126 : * \partial_t \alpha &= 0,\\
127 : * \partial_i \alpha &=
128 : * -\frac{1}{2} \gamma_{rr}^{-3/2} \partial_r \gamma_{rr} \frac{X_i}{r}
129 : * \f}
130 : *
131 : * The shift and its time derivative are
132 : *
133 : * \f{align}
134 : * \beta^i &= \left(\frac{2M}{M+r}\right)^2 \frac{X^i}{r}
135 : * \frac{1}{\gamma_{rr}},\\
136 : * \partial_t \beta^i &= 0
137 : * \f}
138 : *
139 : * The spatial derivative of the shift is computed in terms of the following
140 : * quantities:
141 : *
142 : * \f{align}
143 : * f_3 &= \frac{1}{r} \frac{1}{\gamma_{rr}} \left(\frac{2M}{M+r}\right)^2,\\
144 : * f_4 &=
145 : * -f_3 -
146 : * \frac{1}{M} \frac{1}{\gamma_{rr}}
147 : * \left(\frac{2M}{M+r}\right)^3 -
148 : * \partial_r \gamma_{rr} \left(\frac{2 M}{M+r}
149 : * \frac{1}{\gamma_{rr}}\right)^2
150 : * \f}
151 : *
152 : * In terms of these, the shift's spatial derivative is
153 : *
154 : * \f{align}
155 : * \partial_k \beta^i &= f_4 \frac{X^i}{r} \frac{X_k}{r} + \delta_k^i f_3
156 : * \f}
157 : */
158 1 : class HarmonicSchwarzschild : public AnalyticSolution<3_st>,
159 : public MarkAsAnalyticSolution {
160 : public:
161 0 : struct Mass {
162 0 : using type = double;
163 0 : static constexpr Options::String help = {"Mass of the black hole"};
164 0 : static type lower_bound() { return 0.; }
165 : };
166 0 : struct Center {
167 0 : using type = std::array<double, volume_dim>;
168 0 : static constexpr Options::String help = {
169 : "The [x,y,z] center of the black hole"};
170 : };
171 0 : using options = tmpl::list<Mass, Center>;
172 0 : static constexpr Options::String help{
173 : "Schwarzschild black hole in Cartesian coordinates with harmonic gauge"};
174 :
175 0 : HarmonicSchwarzschild(double mass,
176 : const std::array<double, volume_dim>& center,
177 : const Options::Context& context = {});
178 :
179 0 : HarmonicSchwarzschild() = default;
180 0 : HarmonicSchwarzschild(const HarmonicSchwarzschild& /*rhs*/) = default;
181 0 : HarmonicSchwarzschild& operator=(const HarmonicSchwarzschild& /*rhs*/) =
182 : default;
183 0 : HarmonicSchwarzschild(HarmonicSchwarzschild&& /*rhs*/) = default;
184 0 : HarmonicSchwarzschild& operator=(HarmonicSchwarzschild&& /*rhs*/) = default;
185 0 : ~HarmonicSchwarzschild() = default;
186 :
187 0 : explicit HarmonicSchwarzschild(CkMigrateMessage* /*msg*/);
188 :
189 : /*!
190 : * \brief Computes and returns spacetime quantities for a Schwarzschild black
191 : * hole with harmonic coordinates at a specific Cartesian position
192 : *
193 : * \param x Cartesian coordinates of the position at which to compute
194 : * spacetime quantities
195 : */
196 : template <typename DataType, typename Frame, typename... Tags>
197 1 : tuples::TaggedTuple<Tags...> variables(
198 : const tnsr::I<DataType, volume_dim, Frame>& x, double /*t*/,
199 : tmpl::list<Tags...> /*meta*/) const {
200 : static_assert(
201 : tmpl2::flat_all_v<
202 : tmpl::list_contains_v<tags<DataType, Frame>, Tags>...>,
203 : "At least one of the requested tags is not supported. The requested "
204 : "tags are listed as template parameters of the `variables` function.");
205 : IntermediateVars<DataType, Frame> cache(get_size(*x.begin()));
206 : IntermediateComputer<DataType, Frame> computer(*this, x);
207 : return {cache.get_var(computer, Tags{})...};
208 : }
209 :
210 : // NOLINTNEXTLINE(google-runtime-references)
211 0 : void pup(PUP::er& p);
212 :
213 : /*!
214 : * \brief Return the mass of the black hole
215 : */
216 1 : SPECTRE_ALWAYS_INLINE double mass() const { return mass_; }
217 : /*!
218 : * \brief Return the center of the black hole
219 : */
220 1 : SPECTRE_ALWAYS_INLINE const std::array<double, volume_dim>& center() const {
221 : return center_;
222 : }
223 :
224 : /*!
225 : * \brief Tags defined for intermediates specific to the harmonic
226 : * Schwarzschild solution
227 : */
228 1 : struct internal_tags {
229 : /*!
230 : * \brief Tag for the position of a point relative to the center of the
231 : * black hole
232 : *
233 : * \details Defined as \f$X^i = \left(x^i - C^i\right)\f$, where \f$C^i\f$
234 : * is the Cartesian coordinates of the center of the black hole
235 : * and \f$x^i\f$ is the Cartesian coordinates of the point where we're
236 : * wanting to compute spacetime quantities.
237 : */
238 : template <typename DataType, typename Frame = ::Frame::Inertial>
239 1 : using x_minus_center = ::Tags::TempI<0, 3, Frame, DataType>;
240 : /*!
241 : * \brief Tag for the radius corresponding to the position of a point
242 : * relative to the center of the black hole
243 : *
244 : * \details Defined as \f$r = \sqrt{\delta_{ij} X^i X^j}\f$, where \f$X^i\f$
245 : * is defined by `internal_tags::x_minus_center`.
246 : */
247 : template <typename DataType>
248 1 : using r = ::Tags::TempScalar<1, DataType>;
249 : /*!
250 : * \brief Tag for one over the radius corresponding to the position of a
251 : * point relative to the center of the black hole
252 : *
253 : * \details The quantity \f$r\f$ is the radius defined by
254 : * `internal_tags::r`.
255 : */
256 : template <typename DataType>
257 1 : using one_over_r = ::Tags::TempScalar<2, DataType>;
258 : /*!
259 : * \brief Tag for the intermediate \f$\frac{X^i}{r}\f$
260 : *
261 : * \details The quantity \f$X^i\f$ is defined by
262 : * `internal_tags::x_minus_center` and \f$r\f$ is the radius defined by
263 : * `internal_tags::r`.
264 : */
265 : template <typename DataType, typename Frame = ::Frame::Inertial>
266 1 : using x_over_r = ::Tags::TempI<3, 3, Frame, DataType>;
267 : /*!
268 : * \brief Tag for the intermediate \f$\frac{M}{r}\f$
269 : *
270 : * \details The quantity \f$M\f$ is the mass of the black hole and \f$r\f$
271 : * is the radius defined by `internal_tags::r`.
272 : */
273 : template <typename DataType>
274 1 : using m_over_r = ::Tags::TempScalar<4, DataType>;
275 : /*!
276 : * \brief Tag for the intermediate \f$\sqrt{f_0} = 1 + \frac{M}{r}\f$
277 : *
278 : * \details The quantity \f$M\f$ is the mass of the black hole, \f$r\f$ is
279 : * the radius defined by `internal_tags::r`, and \f$f_0\f$ is defined by
280 : * `internal_tags::f_0`.
281 : */
282 : template <typename DataType>
283 1 : using sqrt_f_0 = ::Tags::TempScalar<5, DataType>;
284 : /*!
285 : * \brief Tag for the intermediate
286 : * \f$f_0 = \left(1 + \frac{M}{r}\right)^2\f$
287 : *
288 : * \details The quantity \f$M\f$ is the mass of the black hole and \f$r\f$
289 : * is the radius defined by `internal_tags::r`.
290 : */
291 : template <typename DataType>
292 1 : using f_0 = ::Tags::TempScalar<6, DataType>;
293 : /*!
294 : * \brief Tag for the intermediate \f$\frac{2M}{M+r}\f$
295 : *
296 : * \details The quantity \f$M\f$ is the mass of the black hole and \f$r\f$
297 : * is the radius defined by `internal_tags::r`.
298 : */
299 : template <typename DataType>
300 1 : using two_m_over_m_plus_r = ::Tags::TempScalar<7, DataType>;
301 : /*!
302 : * \brief Tag for the intermediate \f$\left(\frac{2M}{M+r}\right)^2\f$
303 : *
304 : * \details The quantity \f$M\f$ is the mass of the black hole and \f$r\f$
305 : * is the radius defined by `internal_tags::r`.
306 : */
307 : template <typename DataType>
308 1 : using two_m_over_m_plus_r_squared = ::Tags::TempScalar<8, DataType>;
309 : /*!
310 : * \brief Tag for the intermediate \f$\left(\frac{2M}{M+r}\right)^3\f$
311 : *
312 : * \details The quantity \f$M\f$ is the mass of the black hole and \f$r\f$
313 : * is the radius defined by `internal_tags::r`.
314 : */
315 : template <typename DataType>
316 1 : using two_m_over_m_plus_r_cubed = ::Tags::TempScalar<9, DataType>;
317 : /*!
318 : * \brief Tag for the \f$\gamma_{rr}\f$ component of the spatial metric
319 : *
320 : * \details Defined as
321 : *
322 : * \f{align}
323 : * \gamma_{rr} &= 1 + \frac{2M}{M+r} + \left(\frac{2M}{M+r}\right)^2
324 : * + \left(\frac{2M}{M+r}\right)^3
325 : * \f}
326 : *
327 : * where \f$M\f$ is the mass of the black hole and \f$r\f$ is the radius
328 : * defined by `internal_tags::r`.
329 : */
330 : template <typename DataType>
331 1 : using spatial_metric_rr = ::Tags::TempScalar<10, DataType>;
332 : /*!
333 : * \brief Tag for the intermediate \f$\frac{1}{\gamma_{rr}}\f$
334 : *
335 : * \details The quantity \f$\gamma_{rr}\f$ is defined by
336 : * `internal_tags::spatial_metric_rr`.
337 : */
338 : template <typename DataType>
339 1 : using one_over_spatial_metric_rr = ::Tags::TempScalar<11, DataType>;
340 : /*!
341 : * \brief Tag for the intermediate \f$\gamma_{rr} - f_0\f$
342 : *
343 : * \details The quantity \f$\gamma_{rr}\f$ is defined by
344 : * `internal_tags::spatial_metric_rr` and \f$f_0\f$ is defined by
345 : * `internal_tags::f_0`.
346 : */
347 : template <typename DataType>
348 1 : using spatial_metric_rr_minus_f_0 = ::Tags::TempScalar<12, DataType>;
349 : /*!
350 : * \brief Tag for the intermediate \f$\partial_r \gamma_{rr}\f$
351 : *
352 : * \details Defined as
353 : *
354 : * \f{align}
355 : * \partial_r \gamma_{rr} &= -\frac{1}{2M}\left(\frac{2M}{M+r}\right)^2
356 : * -\frac{1}{M}\left(\frac{2M}{M+r}\right)^3
357 : * -\frac{3}{2M}\left(\frac{2M}{M+r}\right)^4
358 : * \f}
359 : *
360 : * where \f$M\f$ is the mass of the black hole and \f$r\f$ is the radius
361 : * defined by `internal_tags::r`.
362 : */
363 : template <typename DataType>
364 1 : using d_spatial_metric_rr = ::Tags::TempScalar<13, DataType>;
365 : /*!
366 : * \brief Tag for the intermediate \f$\partial_r f_0\f$
367 : *
368 : * \details Defined as
369 : *
370 : * \f{align}
371 : * \partial_r f_0 &=
372 : * 2 \left(1+\frac{M}{r}\right)\left(-\frac{M}{r^2}\right)
373 : * \f}
374 : *
375 : * where \f$M\f$ is the mass of the black hole and \f$r\f$ is the radius
376 : * defined by `internal_tags::r`.
377 : */
378 : template <typename DataType>
379 1 : using d_f_0 = ::Tags::TempScalar<14, DataType>;
380 : /*!
381 : * \brief Tag for the intermediate \f$\partial_r f_0 \frac{X_i}{r}\f$
382 : *
383 : * \details The quantity \f$r\f$ is the radius defined by
384 : * `internal_tags::r`, \f$\partial_r f_0\f$ is defined by
385 : * `internal_tags::d_f_0`, and \f$X_j = X^i \delta_{ij}\f$ where \f$X^i\f$
386 : * is defined by `internal_tags::x_minus_center`.
387 : */
388 : template <typename DataType, typename Frame = ::Frame::Inertial>
389 1 : using d_f_0_times_x_over_r = ::Tags::Tempi<15, 3, Frame, DataType>;
390 : /*!
391 : * \brief Tag for the intermediate
392 : * \f$f_1 = \frac{1}{r} \left(\gamma_{rr} - f_0\right)\f$
393 : *
394 : * \details The quantity \f$r\f$ is the radius defined by
395 : * `internal_tags::r`, \f$\gamma_{rr}\f$ is defined by
396 : * `internal_tags::spatial_metric_rr`, and \f$f_0\f$ is defined by
397 : * `internal_tags::f_0`.
398 : */
399 : template <typename DataType>
400 1 : using f_1 = ::Tags::TempScalar<16, DataType>;
401 : /*!
402 : * \brief Tag for the intermediate \f$f_1 \frac{X_i}{r}\f$
403 : *
404 : * \details The quantity \f$r\f$ is the radius defined by
405 : * `internal_tags::r`, \f$f_1\f$ is defined by `internal_tags::f_1`, and
406 : * \f$X_j = X^i \delta_{ij}\f$ where \f$X^i\f$ is defined by
407 : * `internal_tags::x_minus_center`.
408 : */
409 : template <typename DataType, typename Frame = ::Frame::Inertial>
410 1 : using f_1_times_x_over_r = ::Tags::Tempi<17, 3, Frame, DataType>;
411 : /*!
412 : * \brief Tag for the intermediate
413 : * \f$f_2 = \partial_r \gamma_{rr} - \partial_r f_0 - 2 f_1\f$
414 : *
415 : * \details The quantity \f$\partial_r \gamma_{rr}\f$ is defined by
416 : * `internal_tags::d_spatial_metric_rr`, \f$\partial_r f_0\f$ is defined by
417 : * `internal_tags::d_f_0`, and \f$f_1\f$ is defined by `internal_tags::f_1`.
418 : */
419 : template <typename DataType>
420 1 : using f_2 = ::Tags::TempScalar<18, DataType>;
421 : /*!
422 : * \brief Tag for the intermediate
423 : * \f$f_2 \frac{X_i}{r} \frac{X_j}{r} \frac{X_k}{r}\f$
424 : *
425 : * \details The quantity \f$r\f$ is the radius defined by
426 : * `internal_tags::r`, \f$f_2\f$ is defined by `internal_tags::f_2`, and
427 : * \f$X_j = X^i \delta_{ij}\f$ where \f$X^i\f$ is defined by
428 : * `internal_tags::x_minus_center`.
429 : */
430 : template <typename DataType, typename Frame = ::Frame::Inertial>
431 1 : using f_2_times_xxx_over_r_cubed = ::Tags::Tempiii<19, 3, Frame, DataType>;
432 : /*!
433 : * \brief Tag for the intermediate
434 : * \f$f_3 = \frac{1}{r}\frac{1}{\gamma_{rr}}\left(\frac{2M}{M+r}\right)^2\f$
435 : *
436 : * \details The quantity \f$M\f$ is the mass of the black hole, \f$r\f$ is
437 : * the radius defined by `internal_tags::r`, and \f$\gamma_{rr}\f$ is
438 : * defined by `internal_tags::spatial_metric_rr`.
439 : */
440 : template <typename DataType>
441 1 : using f_3 = ::Tags::TempScalar<20, DataType>;
442 : /*!
443 : * \brief Tag for the intermediate
444 : *
445 : * \f{align}
446 : * f_4 &=
447 : * -f_3 -
448 : * \frac{1}{M} \frac{1}{\gamma_{rr}}
449 : * \left(\frac{2M}{M+r}\right)^3 -
450 : * \partial_r \gamma_{rr} \left(\frac{2 M}{M+r}
451 : * \frac{1}{\gamma_{rr}}\right)^2
452 : * \f}
453 : *
454 : * \details The quantity \f$M\f$ is the mass of the black hole, \f$r\f$ is
455 : * the radius defined by `internal_tags::r`, \f$f_3\f$ is defined by
456 : * `internal_tags::f_3`, \f$\gamma_{rr}\f$ is defined by
457 : * `internal_tags::spatial_metric_rr`, and its derivative
458 : * \f$\partial_r \gamma_{rr}\f$ is defined by
459 : * `internal_tags::d_spatial_metric_rr`.
460 : */
461 : template <typename DataType>
462 1 : using f_4 = ::Tags::TempScalar<21, DataType>;
463 : /*!
464 : * \brief Tag for one over the determinant of the spatial metric
465 : */
466 : template <typename DataType>
467 1 : using one_over_det_spatial_metric = ::Tags::TempScalar<22, DataType>;
468 : /*!
469 : * \brief Tag for the intermediate
470 : * \f$-\frac{1}{2} \gamma_{rr}^{-3/2} \partial_r \gamma_{rr}\f$
471 : *
472 : * \details The lapse is defined as \f$\alpha = \gamma_{rr}^{-1/2}\f$,
473 : * \f$\gamma_{rr}\f$ is defined by `internal_tags::spatial_metric_rr` and
474 : * its derivative \f$\partial_r \gamma_{rr}\f$ is defined by
475 : * `internal_tags::d_spatial_metric_rr`.
476 : */
477 : template <typename DataType>
478 1 : using neg_half_lapse_cubed_times_d_spatial_metric_rr =
479 : ::Tags::TempScalar<23, DataType>;
480 : };
481 :
482 : /*!
483 : * \brief Buffer for caching computed intermediates and quantities that we do
484 : * not want to recompute across the solution's implementation
485 : *
486 : * \details See `internal_tags` documentation for details on what quantities
487 : * the internal tags represent
488 : */
489 : template <typename DataType, typename Frame = ::Frame::Inertial>
490 1 : using CachedBuffer = CachedTempBuffer<
491 : internal_tags::x_minus_center<DataType, Frame>,
492 : internal_tags::r<DataType>, internal_tags::one_over_r<DataType>,
493 : internal_tags::x_over_r<DataType, Frame>,
494 : internal_tags::m_over_r<DataType>, internal_tags::sqrt_f_0<DataType>,
495 : internal_tags::f_0<DataType>,
496 : internal_tags::two_m_over_m_plus_r<DataType>,
497 : internal_tags::two_m_over_m_plus_r_squared<DataType>,
498 : internal_tags::two_m_over_m_plus_r_cubed<DataType>,
499 : internal_tags::spatial_metric_rr<DataType>,
500 : internal_tags::one_over_spatial_metric_rr<DataType>,
501 : internal_tags::spatial_metric_rr_minus_f_0<DataType>,
502 : internal_tags::d_spatial_metric_rr<DataType>,
503 : internal_tags::d_f_0<DataType>,
504 : internal_tags::d_f_0_times_x_over_r<DataType, Frame>,
505 : internal_tags::f_1<DataType>,
506 : internal_tags::f_1_times_x_over_r<DataType, Frame>,
507 : internal_tags::f_2<DataType>,
508 : internal_tags::f_2_times_xxx_over_r_cubed<DataType, Frame>,
509 : internal_tags::f_3<DataType>, internal_tags::f_4<DataType>,
510 : gr::Tags::Lapse<DataType>,
511 : internal_tags::neg_half_lapse_cubed_times_d_spatial_metric_rr<DataType>,
512 : gr::Tags::Shift<DataType, 3, Frame>, DerivShift<DataType, Frame>,
513 : gr::Tags::SpatialMetric<DataType, 3, Frame>,
514 : DerivSpatialMetric<DataType, Frame>,
515 : ::Tags::dt<gr::Tags::SpatialMetric<DataType, 3, Frame>>,
516 : gr::Tags::DetSpatialMetric<DataType>,
517 : internal_tags::one_over_det_spatial_metric<DataType>>;
518 :
519 : /*!
520 : * \brief Computes the intermediates and quantities that we do not want to
521 : * recompute across the solution's implementation
522 : */
523 : template <typename DataType, typename Frame = ::Frame::Inertial>
524 1 : class IntermediateComputer {
525 : public:
526 0 : using CachedBuffer = HarmonicSchwarzschild::CachedBuffer<DataType, Frame>;
527 :
528 : /*!
529 : * \brief Constructs a computer for spacetime quantities of a given
530 : * `gr::Solutions::HarmonicSchwarzschild` solution at at a specific
531 : * Cartesian position
532 : *
533 : * \param solution the given `gr::Solutions::HarmonicSchwarzschild` solution
534 : * \param x Cartesian coordinates of the position at which to compute
535 : * spacetime quantities
536 : */
537 1 : IntermediateComputer(const HarmonicSchwarzschild& solution,
538 : const tnsr::I<DataType, 3, Frame>& x);
539 :
540 : /*!
541 : * \brief Computes the intermediate defined by
542 : * `internal_tags::x_minus_center`
543 : */
544 1 : void operator()(
545 : gsl::not_null<tnsr::I<DataType, 3, Frame>*> x_minus_center,
546 : gsl::not_null<CachedBuffer*> /*cache*/,
547 : internal_tags::x_minus_center<DataType, Frame> /*meta*/) const;
548 :
549 : /*!
550 : * \brief Computes the radius defined by `internal_tags::r`
551 : */
552 1 : void operator()(gsl::not_null<Scalar<DataType>*> r,
553 : gsl::not_null<CachedBuffer*> cache,
554 : internal_tags::r<DataType> /*meta*/) const;
555 :
556 : /*!
557 : * \brief Computes the intermediate defined by `internal_tags::one_over_r`
558 : */
559 1 : void operator()(gsl::not_null<Scalar<DataType>*> one_over_r,
560 : gsl::not_null<CachedBuffer*> cache,
561 : internal_tags::one_over_r<DataType> /*meta*/) const;
562 :
563 : /*!
564 : * \brief Computes the intermediate defined by `internal_tags::x_over_r`
565 : */
566 1 : void operator()(gsl::not_null<tnsr::I<DataType, 3, Frame>*> x_over_r,
567 : gsl::not_null<CachedBuffer*> cache,
568 : internal_tags::x_over_r<DataType, Frame> /*meta*/) const;
569 :
570 : /*!
571 : * \brief Computes the intermediate defined by `internal_tags::m_over_r`
572 : */
573 1 : void operator()(gsl::not_null<Scalar<DataType>*> m_over_r,
574 : gsl::not_null<CachedBuffer*> cache,
575 : internal_tags::m_over_r<DataType> /*meta*/) const;
576 :
577 : /*!
578 : * \brief Computes the intermediate defined by `internal_tags::sqrt_f_0`
579 : */
580 1 : void operator()(gsl::not_null<Scalar<DataType>*> sqrt_f_0,
581 : gsl::not_null<CachedBuffer*> cache,
582 : internal_tags::sqrt_f_0<DataType> /*meta*/) const;
583 :
584 : /*!
585 : * \brief Computes the intermediate defined by `internal_tags::f_0`
586 : */
587 1 : void operator()(gsl::not_null<Scalar<DataType>*> f_0,
588 : gsl::not_null<CachedBuffer*> cache,
589 : internal_tags::f_0<DataType> /*meta*/) const;
590 :
591 : /*!
592 : * \brief Computes the intermediate defined by
593 : * `internal_tags::two_m_over_m_plus_r`
594 : */
595 1 : void operator()(
596 : gsl::not_null<Scalar<DataType>*> two_m_over_m_plus_r,
597 : gsl::not_null<CachedBuffer*> cache,
598 : internal_tags::two_m_over_m_plus_r<DataType> /*meta*/) const;
599 :
600 : /*!
601 : * \brief Computes the intermediate defined by
602 : * `internal_tags::two_m_over_m_plus_r_squared`
603 : */
604 1 : void operator()(
605 : gsl::not_null<Scalar<DataType>*> two_m_over_m_plus_r_squared,
606 : gsl::not_null<CachedBuffer*> cache,
607 : internal_tags::two_m_over_m_plus_r_squared<DataType> /*meta*/) const;
608 :
609 : /*!
610 : * \brief Computes the intermediate defined by
611 : * `internal_tags::two_m_over_m_plus_r_cubed`
612 : */
613 1 : void operator()(
614 : gsl::not_null<Scalar<DataType>*> two_m_over_m_plus_r_cubed,
615 : gsl::not_null<CachedBuffer*> cache,
616 : internal_tags::two_m_over_m_plus_r_cubed<DataType> /*meta*/) const;
617 :
618 : /*!
619 : * \brief Computes the intermediate defined by
620 : * `internal_tags::spatial_metric_rr`
621 : */
622 1 : void operator()(gsl::not_null<Scalar<DataType>*> spatial_metric_rr,
623 : gsl::not_null<CachedBuffer*> cache,
624 : internal_tags::spatial_metric_rr<DataType> /*meta*/) const;
625 :
626 : /*!
627 : * \brief Computes the intermediate defined by
628 : * `internal_tags::one_over_spatial_metric_rr`
629 : */
630 1 : void operator()(
631 : gsl::not_null<Scalar<DataType>*> one_over_spatial_metric_rr,
632 : gsl::not_null<CachedBuffer*> cache,
633 : internal_tags::one_over_spatial_metric_rr<DataType> /*meta*/) const;
634 :
635 : /*!
636 : * \brief Computes the intermediate defined by
637 : * `internal_tags::spatial_metric_rr_minus_f_0`
638 : */
639 1 : void operator()(
640 : gsl::not_null<Scalar<DataType>*> spatial_metric_rr_minus_f_0,
641 : gsl::not_null<CachedBuffer*> cache,
642 : internal_tags::spatial_metric_rr_minus_f_0<DataType> /*meta*/) const;
643 :
644 : /*!
645 : * \brief Computes the intermediate defined by
646 : * `internal_tags::d_spatial_metric_rr`
647 : */
648 1 : void operator()(
649 : gsl::not_null<Scalar<DataType>*> d_spatial_metric_rr,
650 : gsl::not_null<CachedBuffer*> cache,
651 : internal_tags::d_spatial_metric_rr<DataType> /*meta*/) const;
652 :
653 : /*!
654 : * \brief Computes the intermediate defined by `internal_tags::d_f_0`
655 : */
656 1 : void operator()(gsl::not_null<Scalar<DataType>*> d_f_0,
657 : gsl::not_null<CachedBuffer*> cache,
658 : internal_tags::d_f_0<DataType> /*meta*/) const;
659 :
660 : /*!
661 : * \brief Computes the intermediate defined by
662 : * `internal_tags::d_f_0_times_x_over_r`
663 : */
664 1 : void operator()(
665 : gsl::not_null<tnsr::i<DataType, 3, Frame>*> d_f_0_times_x_over_r,
666 : gsl::not_null<CachedBuffer*> cache,
667 : internal_tags::d_f_0_times_x_over_r<DataType, Frame> /*meta*/) const;
668 :
669 : /*!
670 : * \brief Computes the intermediate defined by `internal_tags::f_1`
671 : */
672 1 : void operator()(gsl::not_null<Scalar<DataType>*> f_1,
673 : gsl::not_null<CachedBuffer*> cache,
674 : internal_tags::f_1<DataType> /*meta*/) const;
675 :
676 : /*!
677 : * \brief Computes the intermediate defined by
678 : * `internal_tags::f_1_times_x_over_r`
679 : */
680 1 : void operator()(
681 : gsl::not_null<tnsr::i<DataType, 3, Frame>*> f_1_times_x_over_r,
682 : gsl::not_null<CachedBuffer*> cache,
683 : internal_tags::f_1_times_x_over_r<DataType, Frame> /*meta*/) const;
684 :
685 : /*!
686 : * \brief Computes the intermediate defined by `internal_tags::f_2`
687 : */
688 1 : void operator()(gsl::not_null<Scalar<DataType>*> f_2,
689 : gsl::not_null<CachedBuffer*> cache,
690 : internal_tags::f_2<DataType> /*meta*/) const;
691 :
692 : /*!
693 : * \brief Computes the intermediate defined by
694 : * `internal_tags::f_2_times_xxx_over_r_cubed`
695 : */
696 1 : void operator()(
697 : gsl::not_null<tnsr::iii<DataType, 3, Frame>*>
698 : f_2_times_xxx_over_r_cubed,
699 : gsl::not_null<CachedBuffer*> cache,
700 : internal_tags::f_2_times_xxx_over_r_cubed<DataType, Frame> /*meta*/)
701 : const;
702 :
703 : /*!
704 : * \brief Computes the intermediate defined by `internal_tags::f_3`
705 : */
706 1 : void operator()(gsl::not_null<Scalar<DataType>*> f_3,
707 : gsl::not_null<CachedBuffer*> cache,
708 : internal_tags::f_3<DataType> /*meta*/) const;
709 :
710 : /*!
711 : * \brief Computes the intermediate defined by `internal_tags::f_4`
712 : */
713 1 : void operator()(gsl::not_null<Scalar<DataType>*> f_4,
714 : gsl::not_null<CachedBuffer*> cache,
715 : internal_tags::f_4<DataType> /*meta*/) const;
716 :
717 : /*!
718 : * \brief Computes the lapse
719 : *
720 : * \details Computed as
721 : *
722 : * \f{align}
723 : * \alpha &= \gamma_{rr}^{-1/2}
724 : * \f}
725 : *
726 : * where \f$\gamma_{rr}\f$ is a component of the spatial metric defined by
727 : * `internal_tags::spatial_metric_rr`.
728 : */
729 1 : void operator()(gsl::not_null<Scalar<DataType>*> lapse,
730 : gsl::not_null<CachedBuffer*> cache,
731 : gr::Tags::Lapse<DataType> /*meta*/) const;
732 :
733 : /*!
734 : * \brief Computes the intermediate defined by
735 : * `internal_tags::neg_half_lapse_cubed_times_d_spatial_metric_rr`
736 : */
737 1 : void operator()(
738 : gsl::not_null<Scalar<DataType>*>
739 : neg_half_lapse_cubed_times_d_spatial_metric_rr,
740 : gsl::not_null<CachedBuffer*> cache,
741 : internal_tags::neg_half_lapse_cubed_times_d_spatial_metric_rr<
742 : DataType> /*meta*/) const;
743 :
744 : /*!
745 : * \brief Computes the shift
746 : *
747 : * \details Computed as
748 : *
749 : * \f{align}
750 : * \beta^i &= \left(\frac{2M}{M+r}\right)^2 \frac{X^i}{r}
751 : * \frac{1}{\gamma_{rr}}
752 : * \f}
753 : *
754 : * where \f$M\f$ is the mass of the black hole, \f$r\f$ is the radius
755 : * defined by `internal_tags::r`, \f$X^i\f$ is defined by
756 : * `internal_tags::x_minus_center`, and \f$\gamma_{rr}\f$ is defined by
757 : * `internal_tags::spatial_metric_rr`.
758 : */
759 1 : void operator()(gsl::not_null<tnsr::I<DataType, 3, Frame>*> shift,
760 : gsl::not_null<CachedBuffer*> cache,
761 : gr::Tags::Shift<DataType, 3, Frame> /*meta*/) const;
762 :
763 : /*!
764 : * \brief Computes the spatial derivative of the shift
765 : *
766 : * \details Computed as
767 : *
768 : * \f{align}
769 : * \partial_k \beta^i &=
770 : * f_4 \frac{X^i}{r} \frac{X_k}{r} + \delta_k^i f_3
771 : * \f}
772 : *
773 : * where \f$r\f$ is the radius defined by `internal_tags::r`, \f$X^i\f$ is
774 : * defined by `internal_tags::x_minus_center`, \f$X_j = X^i \delta_{ij}\f$,
775 : * \f$f_3\f$ is defined by `internal_tags::f_3`, and \f$f_4\f$ is defined by
776 : * `internal_tags::f_4`.
777 : */
778 1 : void operator()(gsl::not_null<tnsr::iJ<DataType, 3, Frame>*> deriv_shift,
779 : gsl::not_null<CachedBuffer*> cache,
780 : DerivShift<DataType, Frame> /*meta*/) const;
781 :
782 : /*!
783 : * \brief Computes the spatial metric
784 : *
785 : * \details Computed as
786 : *
787 : * \f{align}
788 : * \gamma_{ij} &=
789 : * \left(\gamma_{rr} - \left(1+\frac{M}{r}\right)^2\right)
790 : * \frac{X_i}{r} \frac{X_j}{r} +
791 : * \delta_{ij} \left(1+\frac{M}{r}\right)^2
792 : * \f}
793 : *
794 : * where \f$M\f$ is the mass of the black hole, \f$r\f$ is the radius
795 : * defined by `internal_tags::r`, \f$\gamma_{rr}\f$ is defined by
796 : * `internal_tags::spatial_metric_rr`, and \f$X_j = X^i \delta_{ij}\f$ where
797 : * \f$X^i\f$ is defined by `internal_tags::x_minus_center`.
798 : */
799 1 : void operator()(gsl::not_null<tnsr::ii<DataType, 3, Frame>*> spatial_metric,
800 : gsl::not_null<CachedBuffer*> cache,
801 : gr::Tags::SpatialMetric<DataType, 3, Frame> /*meta*/) const;
802 :
803 : /*!
804 : * \brief Computes the spatial derivative of the spatial metric
805 : *
806 : * \details Computed as
807 : *
808 : * \f{align}
809 : * \partial_k \gamma_{ij} &=
810 : * f_2 \frac{X_i}{r} \frac{X_j}{r} \frac{X_k}{r} +
811 : * f_1 \frac{X_j}{r} \delta_{ik} + f_1 \frac{X_i}{r} \delta_{jk} +
812 : * \partial_r f_0 \frac{X_k}{r} \delta_{ij}
813 : * \f}
814 : *
815 : * where \f$r\f$ is the radius defined by `internal_tags::r`,
816 : * \f$\partial_r f_0\f$ is defined by `internal_tags::d_f_0`, \f$f_1\f$ is
817 : * defined by `internal_tags::f_1`, \f$f_2\f$ is defined by
818 : * `internal_tags::f_2`, and \f$X_j = X^i \delta_{ij}\f$ where \f$X^i\f$ is
819 : * defined by `internal_tags::x_minus_center`.
820 : */
821 1 : void operator()(
822 : gsl::not_null<tnsr::ijj<DataType, 3, Frame>*> deriv_spatial_metric,
823 : gsl::not_null<CachedBuffer*> cache,
824 : DerivSpatialMetric<DataType, Frame> /*meta*/) const;
825 :
826 : /*!
827 : * \brief Sets the time derivative of the spatial metric to 0
828 : */
829 1 : void operator()(
830 : gsl::not_null<tnsr::ii<DataType, 3, Frame>*> dt_spatial_metric,
831 : gsl::not_null<CachedBuffer*> cache,
832 : ::Tags::dt<gr::Tags::SpatialMetric<DataType, 3, Frame>> /*meta*/) const;
833 :
834 : /*!
835 : * \brief Computes the determinant of the spatial metric
836 : */
837 1 : void operator()(gsl::not_null<Scalar<DataType>*> det_spatial_metric,
838 : gsl::not_null<CachedBuffer*> cache,
839 : gr::Tags::DetSpatialMetric<DataType> /*meta*/) const;
840 :
841 : /*!
842 : * \brief Computes one over the determinant of the spatial metric
843 : */
844 1 : void operator()(
845 : gsl::not_null<Scalar<DataType>*> one_over_det_spatial_metric,
846 : gsl::not_null<CachedBuffer*> cache,
847 : internal_tags::one_over_det_spatial_metric<DataType> /*meta*/) const;
848 :
849 : private:
850 : /*!
851 : * \brief The harmonic Schwarzschild solution
852 : */
853 1 : const HarmonicSchwarzschild& solution_;
854 : /*!
855 : * \brief Cartesian coordinates of the position at which to compute
856 : * spacetime quantities
857 : */
858 1 : const tnsr::I<DataType, 3, Frame>& x_;
859 : };
860 :
861 : /*!
862 : * \brief Computes and returns spacetime quantities of interest
863 : */
864 : template <typename DataType, typename Frame = ::Frame::Inertial>
865 1 : class IntermediateVars : public CachedBuffer<DataType, Frame> {
866 : public:
867 0 : using CachedBuffer = HarmonicSchwarzschild::CachedBuffer<DataType, Frame>;
868 : using CachedBuffer::CachedBuffer;
869 1 : using CachedBuffer::get_var;
870 :
871 : /*!
872 : * \brief Computes and returns the spatial derivative of the lapse
873 : *
874 : * \details Computed as
875 : *
876 : * \f{align}
877 : * \partial_i \alpha &=
878 : * -\frac{1}{2} \gamma_{rr}^{-3/2}
879 : * \partial_r \gamma_{rr} \frac{X_i}{r}
880 : * \f}
881 : *
882 : * where \f$r\f$ is the radius defined by `internal_tags::r`,
883 : * \f$\gamma_{rr}\f$ is defined by `internal_tags::spatial_metric_rr`, its
884 : * derivative \f$\partial_r \gamma_{rr}\f$ is defined by
885 : * `internal_tags::d_spatial_metric_rr`, and \f$X_j = X^i \delta_{ij}\f$
886 : * where \f$X^i\f$ is defined by `internal_tags::x_minus_center`.
887 : */
888 1 : tnsr::i<DataType, 3, Frame> get_var(
889 : const IntermediateComputer<DataType, Frame>& computer,
890 : DerivLapse<DataType, Frame> /*meta*/);
891 :
892 : /*!
893 : * \brief Returns the time derivative of the lapse, which is 0
894 : */
895 1 : Scalar<DataType> get_var(
896 : const IntermediateComputer<DataType, Frame>& computer,
897 : ::Tags::dt<gr::Tags::Lapse<DataType>> /*meta*/);
898 :
899 : /*!
900 : * \brief Returns the time derivative of the shift, which is 0
901 : */
902 1 : tnsr::I<DataType, 3, Frame> get_var(
903 : const IntermediateComputer<DataType, Frame>& computer,
904 : ::Tags::dt<gr::Tags::Shift<DataType, 3, Frame>> /*meta*/);
905 :
906 : /*!
907 : * \brief Computes and returns the square root of the determinant of the
908 : * spatial metric
909 : */
910 1 : Scalar<DataType> get_var(
911 : const IntermediateComputer<DataType, Frame>& computer,
912 : gr::Tags::SqrtDetSpatialMetric<DataType> /*meta*/);
913 :
914 : /*!
915 : * \brief Computes and returns the inverse spatial metric
916 : */
917 1 : tnsr::II<DataType, 3, Frame> get_var(
918 : const IntermediateComputer<DataType, Frame>& computer,
919 : gr::Tags::InverseSpatialMetric<DataType, 3, Frame> /*meta*/);
920 :
921 : /*!
922 : * \brief Computes and returns the extrinsic curvature
923 : */
924 1 : tnsr::ii<DataType, 3, Frame> get_var(
925 : const IntermediateComputer<DataType, Frame>& computer,
926 : gr::Tags::ExtrinsicCurvature<DataType, 3, Frame> /*meta*/);
927 : };
928 :
929 : private:
930 : /*!
931 : * \brief Mass of the black hole
932 : */
933 1 : double mass_{std::numeric_limits<double>::signaling_NaN()};
934 : /*!
935 : * \brief Center of the black hole
936 : */
937 1 : std::array<double, volume_dim> center_ =
938 : make_array<volume_dim>(std::numeric_limits<double>::signaling_NaN());
939 : };
940 :
941 : /*!
942 : * \brief Return whether two harmonic Schwarzschild solutions are equivalent
943 : */
944 1 : SPECTRE_ALWAYS_INLINE bool operator==(const HarmonicSchwarzschild& lhs,
945 : const HarmonicSchwarzschild& rhs) {
946 : return lhs.mass() == rhs.mass() and lhs.center() == rhs.center();
947 : }
948 :
949 : /*!
950 : * \brief Return whether two harmonic Schwarzschild solutions are not equivalent
951 : */
952 1 : SPECTRE_ALWAYS_INLINE bool operator!=(const HarmonicSchwarzschild& lhs,
953 : const HarmonicSchwarzschild& rhs) {
954 : return not(lhs == rhs);
955 : }
956 : } // namespace Solutions
957 : } // namespace gr
|
