Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #include "PointwiseFunctions/AnalyticSolutions/GeneralRelativity/KerrSchild.hpp"
5 :
6 : #include <cmath> // IWYU pragma: keep
7 : #include <numeric>
8 : #include <ostream>
9 : #include <utility>
10 :
11 : #include "DataStructures/DataBox/Prefixes.hpp"
12 : #include "DataStructures/DataVector.hpp" // IWYU pragma: keep
13 : #include "PointwiseFunctions/GeneralRelativity/ExtrinsicCurvature.hpp"
14 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
15 : #include "Utilities/ConstantExpressions.hpp"
16 : #include "Utilities/ContainerHelpers.hpp"
17 : #include "Utilities/Gsl.hpp"
18 : #include "Utilities/MakeWithValue.hpp"
19 : #include "Utilities/StdArrayHelpers.hpp"
20 : #include "Utilities/StdHelpers.hpp"
21 :
22 : /// \cond
23 : namespace gr::Solutions {
24 :
25 : KerrSchild::KerrSchild(const double mass,
26 : KerrSchild::Spin::type dimensionless_spin,
27 : KerrSchild::Center::type center,
28 : const Options::Context& context)
29 : : mass_(mass),
30 : // clang-tidy: do not std::move trivial types.
31 : dimensionless_spin_(std::move(dimensionless_spin)), // NOLINT
32 : // clang-tidy: do not std::move trivial types.
33 : center_(std::move(center)) // NOLINT
34 : {
35 : const double spin_magnitude = magnitude(dimensionless_spin_);
36 : if (spin_magnitude > 1.0) {
37 : PARSE_ERROR(context, "Spin magnitude must be < 1. Given spin: "
38 : << dimensionless_spin_ << " with magnitude "
39 : << spin_magnitude);
40 : }
41 : if (mass_ < 0.0) {
42 : PARSE_ERROR(context, "Mass must be non-negative. Given mass: " << mass_);
43 : }
44 : }
45 :
46 : void KerrSchild::pup(PUP::er& p) noexcept {
47 : p | mass_;
48 : p | dimensionless_spin_;
49 : p | center_;
50 : }
51 :
52 : template <typename DataType>
53 : KerrSchild::IntermediateComputer<DataType>::IntermediateComputer(
54 : const KerrSchild& solution, const tnsr::I<DataType, 3>& x,
55 : const double null_vector_0) noexcept
56 : : solution_(solution), x_(x), null_vector_0_(null_vector_0) {}
57 :
58 : template <typename DataType>
59 : void KerrSchild::IntermediateComputer<DataType>::operator()(
60 : const gsl::not_null<tnsr::I<DataType, 3>*> x_minus_center,
61 : const gsl::not_null<CachedBuffer*> /*cache*/,
62 : internal_tags::x_minus_center<DataType> /*meta*/) const noexcept {
63 : *x_minus_center = x_;
64 : for (size_t d = 0; d < 3; ++d) {
65 : x_minus_center->get(d) -= gsl::at(solution_.center(), d);
66 : }
67 : }
68 :
69 : template <typename DataType>
70 : void KerrSchild::IntermediateComputer<DataType>::operator()(
71 : const gsl::not_null<Scalar<DataType>*> a_dot_x,
72 : const gsl::not_null<CachedBuffer*> cache,
73 : internal_tags::a_dot_x<DataType> /*meta*/) const noexcept {
74 : const auto& x_minus_center =
75 : cache->get_var(internal_tags::x_minus_center<DataType>{});
76 :
77 : const auto spin_a = solution_.dimensionless_spin() * solution_.mass();
78 : get(*a_dot_x) = spin_a[0] * get<0>(x_minus_center) +
79 : spin_a[1] * get<1>(x_minus_center) +
80 : spin_a[2] * get<2>(x_minus_center);
81 : }
82 :
83 : template <typename DataType>
84 : void KerrSchild::IntermediateComputer<DataType>::operator()(
85 : const gsl::not_null<Scalar<DataType>*> a_dot_x_squared,
86 : const gsl::not_null<CachedBuffer*> cache,
87 : internal_tags::a_dot_x_squared<DataType> /*meta*/) const noexcept {
88 : const auto& a_dot_x = get(cache->get_var(internal_tags::a_dot_x<DataType>{}));
89 :
90 : get(*a_dot_x_squared) = square(a_dot_x);
91 : }
92 :
93 : template <typename DataType>
94 : void KerrSchild::IntermediateComputer<DataType>::operator()(
95 : const gsl::not_null<Scalar<DataType>*> half_xsq_minus_asq,
96 : const gsl::not_null<CachedBuffer*> cache,
97 : internal_tags::half_xsq_minus_asq<DataType> /*meta*/) const noexcept {
98 : const auto& x_minus_center =
99 : cache->get_var(internal_tags::x_minus_center<DataType>{});
100 :
101 : const auto spin_a = solution_.dimensionless_spin() * solution_.mass();
102 : const auto a_squared =
103 : std::inner_product(spin_a.begin(), spin_a.end(), spin_a.begin(), 0.);
104 :
105 : get(*half_xsq_minus_asq) =
106 : 0.5 * (square(get<0>(x_minus_center)) + square(get<1>(x_minus_center)) +
107 : square(get<2>(x_minus_center)) - a_squared);
108 : }
109 :
110 : template <typename DataType>
111 : void KerrSchild::IntermediateComputer<DataType>::operator()(
112 : const gsl::not_null<Scalar<DataType>*> r_squared,
113 : const gsl::not_null<CachedBuffer*> cache,
114 : internal_tags::r_squared<DataType> /*meta*/) const noexcept {
115 : const auto& half_xsq_minus_asq =
116 : get(cache->get_var(internal_tags::half_xsq_minus_asq<DataType>{}));
117 : const auto& a_dot_x_squared =
118 : get(cache->get_var(internal_tags::a_dot_x_squared<DataType>{}));
119 :
120 : get(*r_squared) =
121 : half_xsq_minus_asq + sqrt(square(half_xsq_minus_asq) + a_dot_x_squared);
122 : }
123 :
124 : template <typename DataType>
125 : void KerrSchild::IntermediateComputer<DataType>::operator()(
126 : const gsl::not_null<Scalar<DataType>*> r,
127 : const gsl::not_null<CachedBuffer*> cache,
128 : internal_tags::r<DataType> /*meta*/) const noexcept {
129 : const auto& r_squared =
130 : get(cache->get_var(internal_tags::r_squared<DataType>{}));
131 :
132 : get(*r) = sqrt(r_squared);
133 : }
134 :
135 : template <typename DataType>
136 : void KerrSchild::IntermediateComputer<DataType>::operator()(
137 : const gsl::not_null<Scalar<DataType>*> a_dot_x_over_rsquared,
138 : const gsl::not_null<CachedBuffer*> cache,
139 : internal_tags::a_dot_x_over_rsquared<DataType> /*meta*/) const noexcept {
140 : const auto& a_dot_x = get(cache->get_var(internal_tags::a_dot_x<DataType>{}));
141 : const auto& r_squared =
142 : get(cache->get_var(internal_tags::r_squared<DataType>{}));
143 :
144 : get(*a_dot_x_over_rsquared) = a_dot_x / r_squared;
145 : }
146 :
147 : template <typename DataType>
148 : void KerrSchild::IntermediateComputer<DataType>::operator()(
149 : const gsl::not_null<Scalar<DataType>*> deriv_log_r_denom,
150 : const gsl::not_null<CachedBuffer*> cache,
151 : internal_tags::deriv_log_r_denom<DataType> /*meta*/) const noexcept {
152 : const auto& r_squared =
153 : get(cache->get_var(internal_tags::r_squared<DataType>{}));
154 : const auto& half_xsq_minus_asq =
155 : get(cache->get_var(internal_tags::half_xsq_minus_asq<DataType>{}));
156 :
157 : get(*deriv_log_r_denom) = 0.5 / (r_squared - half_xsq_minus_asq);
158 : }
159 :
160 : template <typename DataType>
161 : void KerrSchild::IntermediateComputer<DataType>::operator()(
162 : const gsl::not_null<tnsr::i<DataType, 3>*> deriv_log_r,
163 : const gsl::not_null<CachedBuffer*> cache,
164 : internal_tags::deriv_log_r<DataType> /*meta*/) const noexcept {
165 : const auto& x_minus_center =
166 : cache->get_var(internal_tags::x_minus_center<DataType>{});
167 : const auto& a_dot_x_over_rsquared =
168 : get(cache->get_var(internal_tags::a_dot_x_over_rsquared<DataType>{}));
169 : const auto& deriv_log_r_denom =
170 : get(cache->get_var(internal_tags::deriv_log_r_denom<DataType>{}));
171 :
172 : const auto spin_a = solution_.dimensionless_spin() * solution_.mass();
173 : for (size_t i = 0; i < 3; ++i) {
174 : deriv_log_r->get(i) =
175 : deriv_log_r_denom *
176 : (x_minus_center.get(i) + gsl::at(spin_a, i) * a_dot_x_over_rsquared);
177 : }
178 : }
179 :
180 : template <typename DataType>
181 : void KerrSchild::IntermediateComputer<DataType>::operator()(
182 : const gsl::not_null<Scalar<DataType>*> H_denom,
183 : const gsl::not_null<CachedBuffer*> cache,
184 : internal_tags::H_denom<DataType> /*meta*/) const noexcept {
185 : const auto& r_squared =
186 : get(cache->get_var(internal_tags::r_squared<DataType>{}));
187 : const auto& a_dot_x_squared =
188 : get(cache->get_var(internal_tags::a_dot_x_squared<DataType>{}));
189 :
190 : get(*H_denom) = 1.0 / (square(r_squared) + a_dot_x_squared);
191 : }
192 :
193 : template <typename DataType>
194 : void KerrSchild::IntermediateComputer<DataType>::operator()(
195 : const gsl::not_null<Scalar<DataType>*> H,
196 : const gsl::not_null<CachedBuffer*> cache,
197 : internal_tags::H<DataType> /*meta*/) const noexcept {
198 : const auto& r = get(cache->get_var(internal_tags::r<DataType>{}));
199 : const auto& r_squared =
200 : get(cache->get_var(internal_tags::r_squared<DataType>{}));
201 : const auto& H_denom = get(cache->get_var(internal_tags::H_denom<DataType>{}));
202 :
203 : get(*H) = solution_.mass() * r * r_squared * H_denom;
204 : }
205 :
206 : template <typename DataType>
207 : void KerrSchild::IntermediateComputer<DataType>::operator()(
208 : const gsl::not_null<Scalar<DataType>*> deriv_H_temp1,
209 : const gsl::not_null<CachedBuffer*> cache,
210 : internal_tags::deriv_H_temp1<DataType> /*meta*/) const noexcept {
211 : const auto& H = get(cache->get_var(internal_tags::H<DataType>{}));
212 : const auto& r_squared =
213 : get(cache->get_var(internal_tags::r_squared<DataType>{}));
214 : const auto& H_denom = get(cache->get_var(internal_tags::H_denom<DataType>{}));
215 :
216 : get(*deriv_H_temp1) = H * (3.0 - 4.0 * square(r_squared) * H_denom);
217 : }
218 :
219 : template <typename DataType>
220 : void KerrSchild::IntermediateComputer<DataType>::operator()(
221 : const gsl::not_null<Scalar<DataType>*> deriv_H_temp2,
222 : const gsl::not_null<CachedBuffer*> cache,
223 : internal_tags::deriv_H_temp2<DataType> /*meta*/) const noexcept {
224 : const auto& H = get(cache->get_var(internal_tags::H<DataType>{}));
225 : const auto& H_denom = get(cache->get_var(internal_tags::H_denom<DataType>{}));
226 : const auto& a_dot_x = get(cache->get_var(internal_tags::a_dot_x<DataType>{}));
227 :
228 : get(*deriv_H_temp2) = H * (2.0 * H_denom * a_dot_x);
229 : }
230 :
231 : template <typename DataType>
232 : void KerrSchild::IntermediateComputer<DataType>::operator()(
233 : const gsl::not_null<tnsr::i<DataType, 3>*> deriv_H,
234 : const gsl::not_null<CachedBuffer*> cache,
235 : internal_tags::deriv_H<DataType> /*meta*/) const noexcept {
236 : const auto& deriv_log_r =
237 : cache->get_var(internal_tags::deriv_log_r<DataType>{});
238 : const auto& deriv_H_temp1 =
239 : get(cache->get_var(internal_tags::deriv_H_temp1<DataType>{}));
240 : const auto& deriv_H_temp2 =
241 : get(cache->get_var(internal_tags::deriv_H_temp2<DataType>{}));
242 :
243 : const auto spin_a = solution_.dimensionless_spin() * solution_.mass();
244 : for (size_t i = 0; i < 3; ++i) {
245 : deriv_H->get(i) =
246 : deriv_H_temp1 * deriv_log_r.get(i) - deriv_H_temp2 * gsl::at(spin_a, i);
247 : }
248 : }
249 :
250 : template <typename DataType>
251 : void KerrSchild::IntermediateComputer<DataType>::operator()(
252 : const gsl::not_null<Scalar<DataType>*> denom,
253 : const gsl::not_null<CachedBuffer*> cache,
254 : internal_tags::denom<DataType> /*meta*/) const noexcept {
255 : const auto& r_squared =
256 : get(cache->get_var(internal_tags::r_squared<DataType>{}));
257 :
258 : const auto spin_a = solution_.dimensionless_spin() * solution_.mass();
259 : const auto a_squared =
260 : std::inner_product(spin_a.begin(), spin_a.end(), spin_a.begin(), 0.);
261 :
262 : get(*denom) = 1.0 / (r_squared + a_squared);
263 : }
264 :
265 : template <typename DataType>
266 : void KerrSchild::IntermediateComputer<DataType>::operator()(
267 : const gsl::not_null<Scalar<DataType>*> a_dot_x_over_r,
268 : const gsl::not_null<CachedBuffer*> cache,
269 : internal_tags::a_dot_x_over_r<DataType> /*meta*/) const noexcept {
270 : const auto& a_dot_x = get(cache->get_var(internal_tags::a_dot_x<DataType>{}));
271 : const auto& r = get(cache->get_var(internal_tags::r<DataType>{}));
272 :
273 : get(*a_dot_x_over_r) = a_dot_x / r;
274 : }
275 :
276 : template <typename DataType>
277 : void KerrSchild::IntermediateComputer<DataType>::operator()(
278 : const gsl::not_null<tnsr::i<DataType, 3>*> null_form,
279 : const gsl::not_null<CachedBuffer*> cache,
280 : internal_tags::null_form<DataType> /*meta*/) const noexcept {
281 : const auto& a_dot_x_over_r =
282 : get(cache->get_var(internal_tags::a_dot_x_over_r<DataType>{}));
283 : const auto& r = get(cache->get_var(internal_tags::r<DataType>{}));
284 : const auto& x_minus_center =
285 : cache->get_var(internal_tags::x_minus_center<DataType>{});
286 : const auto& denom = get(cache->get_var(internal_tags::denom<DataType>{}));
287 :
288 : const auto spin_a = solution_.dimensionless_spin() * solution_.mass();
289 : for (size_t i = 0; i < 3; ++i) {
290 : const size_t cross_product_index_1 = (i + 1) % 3;
291 : const size_t cross_product_index_2 = (i + 2) % 3;
292 : null_form->get(i) = denom * (r * x_minus_center.get(i) -
293 : gsl::at(spin_a, cross_product_index_1) *
294 : x_minus_center.get(cross_product_index_2) +
295 : gsl::at(spin_a, cross_product_index_2) *
296 : x_minus_center.get(cross_product_index_1) +
297 : a_dot_x_over_r * gsl::at(spin_a, i));
298 : }
299 : }
300 :
301 : template <typename DataType>
302 : void KerrSchild::IntermediateComputer<DataType>::operator()(
303 : const gsl::not_null<tnsr::ij<DataType, 3>*> deriv_null_form,
304 : const gsl::not_null<CachedBuffer*> cache,
305 : internal_tags::deriv_null_form<DataType> /*meta*/) const noexcept {
306 : const auto spin_a = solution_.dimensionless_spin() * solution_.mass();
307 : const auto& denom = get(cache->get_var(internal_tags::denom<DataType>{}));
308 : const auto& r = get(cache->get_var(internal_tags::r<DataType>{}));
309 : const auto& x_minus_center =
310 : cache->get_var(internal_tags::x_minus_center<DataType>{});
311 : const auto& null_form = cache->get_var(internal_tags::null_form<DataType>{});
312 : const auto& a_dot_x_over_rsquared =
313 : get(cache->get_var(internal_tags::a_dot_x_over_rsquared<DataType>{}));
314 : const auto& deriv_log_r =
315 : cache->get_var(internal_tags::deriv_log_r<DataType>{});
316 :
317 : for (size_t i = 0; i < 3; i++) {
318 : for (size_t j = 0; j < 3; j++) {
319 : deriv_null_form->get(j, i) =
320 : denom * (gsl::at(spin_a, i) * gsl::at(spin_a, j) / r +
321 : (x_minus_center.get(i) - 2.0 * r * null_form.get(i) -
322 : a_dot_x_over_rsquared * gsl::at(spin_a, i)) *
323 : deriv_log_r.get(j) * r);
324 : if (i == j) {
325 : deriv_null_form->get(j, i) += denom * r;
326 : } else { // add denom*epsilon^ijk a_k
327 : size_t k = (j + 1) % 3;
328 : if (k == i) { // j+1 = i (cyclic), so choose minus sign
329 : k++;
330 : k = k % 3; // and set k to be neither i nor j
331 : deriv_null_form->get(j, i) -= denom * gsl::at(spin_a, k);
332 : } else { // i+1 = j (cyclic), so choose plus sign
333 : deriv_null_form->get(j, i) += denom * gsl::at(spin_a, k);
334 : }
335 : }
336 : }
337 : }
338 : }
339 :
340 : template <typename DataType>
341 : void KerrSchild::IntermediateComputer<DataType>::operator()(
342 : const gsl::not_null<Scalar<DataType>*> lapse_squared,
343 : const gsl::not_null<CachedBuffer*> cache,
344 : internal_tags::lapse_squared<DataType> /*meta*/) const noexcept {
345 : const auto& H = get(cache->get_var(internal_tags::H<DataType>{}));
346 : get(*lapse_squared) = 1.0 / (1.0 + 2.0 * square(null_vector_0_) * H);
347 : }
348 :
349 : template <typename DataType>
350 : void KerrSchild::IntermediateComputer<DataType>::operator()(
351 : const gsl::not_null<Scalar<DataType>*> lapse,
352 : const gsl::not_null<CachedBuffer*> cache,
353 : gr::Tags::Lapse<DataType> /*meta*/) const noexcept {
354 : const auto& lapse_squared =
355 : get(cache->get_var(internal_tags::lapse_squared<DataType>{}));
356 : get(*lapse) = sqrt(lapse_squared);
357 : }
358 :
359 : template <typename DataType>
360 : void KerrSchild::IntermediateComputer<DataType>::operator()(
361 : const gsl::not_null<Scalar<DataType>*> deriv_lapse_multiplier,
362 : const gsl::not_null<CachedBuffer*> cache,
363 : internal_tags::deriv_lapse_multiplier<DataType> /*meta*/) const noexcept {
364 : const auto& lapse = get(cache->get_var(gr::Tags::Lapse<DataType>{}));
365 : const auto& lapse_squared =
366 : get(cache->get_var(internal_tags::lapse_squared<DataType>{}));
367 : get(*deriv_lapse_multiplier) =
368 : -square(null_vector_0_) * lapse * lapse_squared;
369 : }
370 :
371 : template <typename DataType>
372 : void KerrSchild::IntermediateComputer<DataType>::operator()(
373 : const gsl::not_null<Scalar<DataType>*> shift_multiplier,
374 : const gsl::not_null<CachedBuffer*> cache,
375 : internal_tags::shift_multiplier<DataType> /*meta*/) const noexcept {
376 : const auto& H = get(cache->get_var(internal_tags::H<DataType>{}));
377 : const auto& lapse_squared =
378 : get(cache->get_var(internal_tags::lapse_squared<DataType>{}));
379 :
380 : get(*shift_multiplier) = -2.0 * null_vector_0_ * H * lapse_squared;
381 : }
382 :
383 : template <typename DataType>
384 : void KerrSchild::IntermediateComputer<DataType>::operator()(
385 : const gsl::not_null<tnsr::I<DataType, 3>*> shift,
386 : const gsl::not_null<CachedBuffer*> cache,
387 : gr::Tags::Shift<3, Frame::Inertial, DataType> /*meta*/) const noexcept {
388 : const auto& null_form = cache->get_var(internal_tags::null_form<DataType>{});
389 : const auto& shift_multiplier =
390 : get(cache->get_var(internal_tags::shift_multiplier<DataType>{}));
391 :
392 : for (size_t i = 0; i < 3; ++i) {
393 : shift->get(i) = shift_multiplier * null_form.get(i);
394 : }
395 : }
396 :
397 : template <typename DataType>
398 : void KerrSchild::IntermediateComputer<DataType>::operator()(
399 : const gsl::not_null<tnsr::iJ<DataType, 3>*> deriv_shift,
400 : const gsl::not_null<CachedBuffer*> cache,
401 : DerivShift<DataType> /*meta*/) const noexcept {
402 : const auto& H = get(cache->get_var(internal_tags::H<DataType>{}));
403 : const auto& null_form = cache->get_var(internal_tags::null_form<DataType>{});
404 : const auto& lapse_squared =
405 : get(cache->get_var(internal_tags::lapse_squared<DataType>{}));
406 : const auto& deriv_H = cache->get_var(internal_tags::deriv_H<DataType>{});
407 : const auto& deriv_null_form =
408 : cache->get_var(internal_tags::deriv_null_form<DataType>{});
409 :
410 : for (size_t m = 0; m < 3; ++m) {
411 : for (size_t i = 0; i < 3; ++i) {
412 : deriv_shift->get(m, i) = 4.0 * cube(null_vector_0_) * H *
413 : null_form.get(i) * square(lapse_squared) *
414 : deriv_H.get(m) -
415 : 2.0 * null_vector_0_ * lapse_squared *
416 : (null_form.get(i) * deriv_H.get(m) +
417 : H * deriv_null_form.get(m, i));
418 : }
419 : }
420 : }
421 :
422 : template <typename DataType>
423 : void KerrSchild::IntermediateComputer<DataType>::operator()(
424 : const gsl::not_null<tnsr::ii<DataType, 3>*> spatial_metric,
425 : const gsl::not_null<CachedBuffer*> cache,
426 : gr::Tags::SpatialMetric<3, Frame::Inertial, DataType> /*meta*/)
427 : const noexcept {
428 : const auto& H = get(cache->get_var(internal_tags::H<DataType>{}));
429 : const auto& null_form = cache->get_var(internal_tags::null_form<DataType>{});
430 :
431 : std::fill(spatial_metric->begin(), spatial_metric->end(), 0.);
432 : for (size_t i = 0; i < 3; ++i) {
433 : spatial_metric->get(i, i) = 1.;
434 : for (size_t j = i; j < 3; ++j) { // Symmetry
435 : spatial_metric->get(i, j) +=
436 : 2.0 * H * null_form.get(i) * null_form.get(j);
437 : }
438 : }
439 : }
440 :
441 : template <typename DataType>
442 : void KerrSchild::IntermediateComputer<DataType>::operator()(
443 : const gsl::not_null<tnsr::ijj<DataType, 3>*> deriv_spatial_metric,
444 : const gsl::not_null<CachedBuffer*> cache,
445 : DerivSpatialMetric<DataType> /*meta*/) const noexcept {
446 : const auto& null_form = cache->get_var(internal_tags::null_form<DataType>{});
447 : const auto& deriv_H = cache->get_var(internal_tags::deriv_H<DataType>{});
448 : const auto& H = get(cache->get_var(internal_tags::H<DataType>{}));
449 : const auto& deriv_null_form =
450 : cache->get_var(internal_tags::deriv_null_form<DataType>{});
451 :
452 : for (size_t i = 0; i < 3; ++i) {
453 : for (size_t j = i; j < 3; ++j) { // Symmetry
454 : for (size_t m = 0; m < 3; ++m) {
455 : deriv_spatial_metric->get(m, i, j) =
456 : 2.0 * null_form.get(i) * null_form.get(j) * deriv_H.get(m) +
457 : 2.0 * H *
458 : (null_form.get(i) * deriv_null_form.get(m, j) +
459 : null_form.get(j) * deriv_null_form.get(m, i));
460 : }
461 : }
462 : }
463 : }
464 :
465 : template <typename DataType>
466 : void KerrSchild::IntermediateComputer<DataType>::operator()(
467 : const gsl::not_null<tnsr::ii<DataType, 3>*> dt_spatial_metric,
468 : const gsl::not_null<CachedBuffer*> /*cache*/,
469 : ::Tags::dt<gr::Tags::SpatialMetric<3, Frame::Inertial, DataType>> /*meta*/)
470 : const noexcept {
471 : std::fill(dt_spatial_metric->begin(), dt_spatial_metric->end(), 0.);
472 : }
473 :
474 : template <typename DataType>
475 : KerrSchild::IntermediateVars<DataType>::IntermediateVars(
476 : const KerrSchild& solution, const tnsr::I<DataType, 3>& x) noexcept
477 : : CachedBuffer(get_size(::get<0>(x)), IntermediateComputer<DataType>(
478 : solution, x, null_vector_0_)) {}
479 :
480 : template <typename DataType>
481 : tnsr::i<DataType, 3> KerrSchild::IntermediateVars<DataType>::get_var(
482 : DerivLapse<DataType> /*meta*/) noexcept {
483 : tnsr::i<DataType, 3> result{};
484 : const auto& deriv_H = get_var(internal_tags::deriv_H<DataType>{});
485 : const auto& deriv_lapse_multiplier =
486 : get(get_var(internal_tags::deriv_lapse_multiplier<DataType>{}));
487 :
488 : for (size_t i = 0; i < 3; ++i) {
489 : result.get(i) = deriv_lapse_multiplier * deriv_H.get(i);
490 : }
491 : return result;
492 : }
493 :
494 : template <typename DataType>
495 : Scalar<DataType> KerrSchild::IntermediateVars<DataType>::get_var(
496 : ::Tags::dt<gr::Tags::Lapse<DataType>> /*meta*/) noexcept {
497 : const auto& H = get(get_var(internal_tags::H<DataType>{}));
498 : return make_with_value<Scalar<DataType>>(H, 0.);
499 : }
500 :
501 : template <typename DataType>
502 : tnsr::I<DataType, 3> KerrSchild::IntermediateVars<DataType>::get_var(
503 : ::Tags::dt<
504 : gr::Tags::Shift<3, Frame::Inertial, DataType>> /*meta*/) noexcept {
505 : const auto& H = get(get_var(internal_tags::H<DataType>()));
506 : return make_with_value<tnsr::I<DataType, 3>>(H, 0.);
507 : }
508 :
509 : template <typename DataType>
510 : Scalar<DataType> KerrSchild::IntermediateVars<DataType>::get_var(
511 : gr::Tags::SqrtDetSpatialMetric<DataType> /*meta*/) noexcept {
512 : return Scalar<DataType>(1.0 / get(get_var(gr::Tags::Lapse<DataType>{})));
513 : }
514 :
515 : template <typename DataType>
516 : tnsr::II<DataType, 3> KerrSchild::IntermediateVars<DataType>::get_var(
517 : gr::Tags::InverseSpatialMetric<3, Frame::Inertial,
518 : DataType> /*meta*/) noexcept {
519 : const auto& H = get(get_var(internal_tags::H<DataType>{}));
520 : const auto& lapse_squared =
521 : get(get_var(internal_tags::lapse_squared<DataType>{}));
522 : const auto& null_form = get_var(internal_tags::null_form<DataType>{});
523 :
524 : auto result = make_with_value<tnsr::II<DataType, 3>>(H, 0.);
525 : for (size_t i = 0; i < 3; ++i) {
526 : result.get(i, i) = 1.;
527 : for (size_t j = i; j < 3; ++j) { // Symmetry
528 : result.get(i, j) -=
529 : 2.0 * H * lapse_squared * null_form.get(i) * null_form.get(j);
530 : }
531 : }
532 :
533 : return result;
534 : }
535 :
536 : template <typename DataType>
537 : tnsr::ii<DataType, 3> KerrSchild::IntermediateVars<DataType>::get_var(
538 : gr::Tags::ExtrinsicCurvature<3, Frame::Inertial,
539 : DataType> /*meta*/) noexcept {
540 : return gr::extrinsic_curvature(
541 : get_var(gr::Tags::Lapse<DataType>{}),
542 : get_var(gr::Tags::Shift<3, Frame::Inertial, DataType>{}),
543 : get_var(DerivShift<DataType>{}),
544 : get_var(gr::Tags::SpatialMetric<3, Frame::Inertial, DataType>{}),
545 : get_var(
546 : ::Tags::dt<gr::Tags::SpatialMetric<3, Frame::Inertial, DataType>>{}),
547 : get_var(DerivSpatialMetric<DataType>{}));
548 : }
549 :
550 : template class KerrSchild::IntermediateVars<DataVector>;
551 : template class KerrSchild::IntermediateVars<double>;
552 : template class KerrSchild::IntermediateComputer<DataVector>;
553 : template class KerrSchild::IntermediateComputer<double>;
554 : } // namespace gr::Solutions
555 : /// \endcond
|
