SpECTRE Documentation Coverage Report
Current view: top level - Domain/CoordinateMaps - KerrHorizonConforming.hpp Hit Total Coverage
Commit: d0fc80462417e83e5cddfa1b9901bb4a9b6af4d6 Lines: 2 16 12.5 %
Date: 2024-03-29 00:33:31
Legend: Lines: hit not hit

          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 <cstddef>
       8             : #include <optional>
       9             : 
      10             : #include "DataStructures/Tensor/Tensor.hpp"
      11             : #include "Utilities/Gsl.hpp"
      12             : #include "Utilities/TypeTraits/RemoveReferenceWrapper.hpp"
      13             : 
      14             : namespace domain::CoordinateMaps {
      15             : 
      16             : /*!
      17             :  * \brief Distorts cartesian coordinates \f$x^i\f$ such that a coordinate sphere
      18             :  * \f$\delta_{ij}x^ix^j=C^2\f$ is mapped to an ellipsoid of constant
      19             :  * Kerr-Schild radius \f$r=C\f$.
      20             :  *
      21             :  * The Kerr-Schild radius \f$r\f$ is defined as the largest positive
      22             :  * root of
      23             :  *
      24             :  *  \f{equation*}{
      25             :  *  r^4 - r^2 (x^2 - a^2) - (\vec{a}\cdot \vec{x})^2 = 0.
      26             :  *  \f}
      27             :  *
      28             :  * In this equation, \f$\vec{x}\f$ are the coordinates of the (distorted)
      29             :  * surface, and \f$\vec{a}=\vec{S}/M\f$ is the spin-parameter of the black hole.
      30             :  * This is equivalent to the implicit definition of \f$r\f$ in Eq. (47) of
      31             :  * \cite Lovelace2008tw. The black hole is assumed to be at the origin.
      32             :  *
      33             :  * \details Given a spin vector \f$\vec{a}\f$, we define the Kerr-Schild radius:
      34             :  *
      35             :  * \f{equation}{r(\vec{x}) = \sqrt{\frac{x^2 - a^2 + \sqrt{(x^2 -
      36             :  * a^2)^2 + 4 (\vec{x} \cdot \vec{a})^2}}{2}} \f}
      37             :  *
      38             :  * We also define the auxiliary variable:
      39             :  *
      40             :  * \f{equation}{s(\vec{x}) = \frac{x^2(x^2 + a^2)}{x^4 + (\vec{x}
      41             :  * \cdot \vec{a})^2} \f}
      42             :  *
      43             :  * The map is then given by:
      44             :  *
      45             :  * \f{equation}{\vec{x}(\vec{\xi}) = \vec{\xi} \sqrt{s(\vec{\xi})} \f}
      46             :  *
      47             :  * The inverse map is given by:
      48             :  * \f{equation}{\vec{\xi}(\vec{x}) = \vec{x} \frac{r(\vec{x})}{|\vec{x}|}
      49             :  * \f}
      50             :  *
      51             :  * The jacobian is:
      52             :  * \f{equation}{ \frac{\partial \xi^i}{\partial x^j} =
      53             :  * \delta_{ij} \sqrt{s(\vec{\xi})} +
      54             :  * \frac{\xi_j}{2 \sqrt{s(\vec{\xi})}} \frac{4\xi^2(1 -
      55             :  * s(\vec{\xi})) \xi_j + 2a^2\xi_i + 2 s(\vec{\xi})
      56             :  * (\vec{\xi} \cdot \vec{a}) a_i}{\xi^4 + (\vec{\xi} \cdot
      57             :  * \vec{a})^2} \f}
      58             :  *
      59             :  * The inverse jacobian is:
      60             :  * \f{equation}{\frac{\partial x^i}{\partial \xi^j} =
      61             :  * \delta_{ij}\frac{r(\vec{x})}{|\vec{x}|}
      62             :  * + \frac{r(\vec{x})^2 x_i + (\vec{x} \cdot \vec{a})a_i}{2r(\vec{x})^3
      63             :  * - r(\vec{x}) (x^2 - a^2)}\frac{x_j}{|\vec{x}|} - \frac{r(\vec{x})
      64             :  * x_i x_j}{x^3} \f}
      65             :  */
      66           1 : class KerrHorizonConforming {
      67             :  public:
      68           0 :   KerrHorizonConforming() = default;
      69           0 :   static constexpr size_t dim = 3;
      70             :   /*!
      71             :    * \brief Constructs a Kerr horizon conforming map.
      72             :    *
      73             :    * \param mass The Kerr mass parameter $M$
      74             :    * \param dimensionless_spin The dimensionless spin $\vec{\chi} = \vec{a} / M
      75             :    * = \vec{S} / M^2$, where $M$ is the Kerr mass parameter, $\vec{S}$ is the
      76             :    * angular momentum or quasilocal spin, and $\vec{a}$ is the Kerr spin
      77             :    * parameter.
      78             :    *
      79             :    * \note The horizon depends only on the dimensionful spin parameter $\vec{a}
      80             :    * = M \vec{\chi}$. This constructor takes $M$ and $\vec{\chi}$ separately for
      81             :    * consistency with other code such as gr::Solutions::KerrSchild, and hence to
      82             :    * avoid bugs where the wrong spin quantity is used accidentally.
      83             :    */
      84           1 :   explicit KerrHorizonConforming(const double mass,
      85             :                                  std::array<double, 3> dimensionless_spin);
      86             : 
      87             :   template <typename T>
      88           0 :   std::array<tt::remove_cvref_wrap_t<T>, 3> operator()(
      89             :       const std::array<T, 3>& source_coords) const;
      90             : 
      91           0 :   std::optional<std::array<double, 3>> inverse(
      92             :       const std::array<double, 3>& target_coords) const;
      93             : 
      94             :   template <typename T>
      95           0 :   tnsr::Ij<tt::remove_cvref_wrap_t<T>, 3, Frame::NoFrame> jacobian(
      96             :       const std::array<T, 3>& source_coords) const;
      97             : 
      98             :   template <typename T>
      99           0 :   tnsr::Ij<tt::remove_cvref_wrap_t<T>, 3, Frame::NoFrame> inv_jacobian(
     100             :       const std::array<T, 3>& source_coords) const;
     101             : 
     102           0 :   bool is_identity() const {
     103             :     return spin_parameter_ == std::array<double, 3>{0., 0., 0.};
     104             :   }
     105             : 
     106           0 :   friend bool operator==(const KerrHorizonConforming& lhs,
     107             :                          const KerrHorizonConforming& rhs);
     108             : 
     109           0 :   void pup(PUP::er& p);
     110             : 
     111             :  private:
     112             :   template <typename T>
     113           0 :   void stretch_factor_square(
     114             :       const gsl::not_null<tt::remove_cvref_wrap_t<T>*> result,
     115             :       const std::array<T, 3>& source_coords) const;
     116             : 
     117           0 :   std::array<double, 3> spin_parameter_;
     118           0 :   double spin_mag_sq_;
     119             : };
     120           0 : bool operator!=(const KerrHorizonConforming& lhs,
     121             :                 const KerrHorizonConforming& rhs);
     122             : 
     123             : }  // namespace domain::CoordinateMaps

Generated by: LCOV version 1.14