SpECTRE Documentation Coverage Report
Current view: top level - PointwiseFunctions/AnalyticSolutions/GeneralRelativity - KerrSchild.cpp Hit Total Coverage
Commit: 2db722c93a8e9b106e406b439b79c8e05c2057fb Lines: 0 1 0.0 %
Date: 2021-03-03 22:01:00
Legend: Lines: hit not hit

          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

Generated by: LCOV version 1.14