Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <limits>
       7             : 
       8             : #include "DataStructures/Tensor/Tensor.hpp"
       9             : #include "NumericalAlgorithms/Interpolation/BarycentricRational.hpp"
      10             : #include "PointwiseFunctions/AnalyticSolutions/RelativisticEuler/TovStar.hpp"
      11             : #include "PointwiseFunctions/Hydro/EquationsOfState/EquationOfState.hpp"  // IWYU pragma: keep
      12             : 
      13             : /// \cond
      14             : namespace PUP {
      15             : class er;
      16             : }  // namespace PUP
      17             : /// \endcond
      18             : 
      19             : // IWYU pragma: no_forward_declare EquationsOfState::EquationOfState
      20             : // IWYU pragma: no_forward_declare Tensor
      21             : 
      22             : namespace gr {
      23             : namespace Solutions {
      24             : 
      25             : /*!
      26             :  * \brief TOV solver based on Lindblom's method
      27             :  *
      28             :  * Uses Lindblom's method of integrating the TOV equations from
      29             :  * \cite Lindblom1998dp
      30             :  *
      31             :  * Instead of integrating \f$m(r)\f$ and \f$p(r)\f$
      32             :  * (\f$r\f$=radius, \f$m\f$=mass, \f$p\f$=pressure)
      33             :  * Lindblom introduces the variables \f$u\f$ and \f$v\f$, with \f$u=r^{2}\f$ and
      34             :  * \f$v=m/r\f$.
      35             :  * The integration is then done with the log of the specific enthalpy
      36             :  * (\f$\mathrm{log}(h)\f$) as the independent variable.
      37             :  *
      38             :  * Lindblom's paper simply labels the independent variable as \f$h\f$.
      39             :  * The \f$h\f$ in Lindblom's paper is NOT the specific enthalpy.
      40             :  * Rather, Lindblom's \f$h\f$ is in fact \f$\mathrm{log}(h)\f$.
      41             :  */
      42           1 : class TovSolution {
      43             :  public:
      44           0 :   TovSolution(
      45             :       const EquationsOfState::EquationOfState<true, 1>& equation_of_state,
      46             :       double central_mass_density, double log_enthalpy_at_outer_radius = 0.0,
      47             :       double absolute_tolerance = 1.0e-14, double relative_tolerance = 1.0e-14);
      48             : 
      49           0 :   TovSolution() = default;
      50           0 :   TovSolution(const TovSolution& /*rhs*/) = delete;
      51           0 :   TovSolution& operator=(const TovSolution& /*rhs*/) = delete;
      52           0 :   TovSolution(TovSolution&& /*rhs*/) noexcept = default;
      53           0 :   TovSolution& operator=(TovSolution&& /*rhs*/) noexcept = default;
      54           0 :   ~TovSolution() = default;
      55             : 
      56             :   /// \brief The outer radius of the solution.
      57             :   ///
      58             :   /// \note This is the radius at which `log_specific_enthalpy` is equal
      59             :   /// to the value of `log_enthalpy_at_outer_radius` that was given when
      60             :   /// constructing this TovSolution
      61           1 :   double outer_radius() const noexcept;
      62             : 
      63             :   /// \brief The mass inside the given radius over the radius
      64             :   /// \f$\frac{m(r)}{r}\f$
      65             :   ///
      66             :   /// \note `r` should be non-negative and not greater than `outer_radius()`.
      67           1 :   double mass_over_radius(double r) const noexcept;
      68             : 
      69             :   /// \brief The mass inside the given radius \f$m(r)\f$
      70             :   ///
      71             :   /// \warning When computing \f$\frac{m(r)}{r}\f$, use the `mass_over_radius`
      72             :   /// function instead for greater accuracy.
      73             :   ///
      74             :   /// \note `r` should be non-negative and not greater than `outer_radius()`
      75           1 :   double mass(double r) const noexcept;
      76             : 
      77             :   /// \brief The log of the specific enthalpy at the given radius
      78             :   ///
      79             :   /// \note `r` should be non-negative and not greater than `outer_radius()`
      80           1 :   double log_specific_enthalpy(double r) const noexcept;
      81             : 
      82             :   /// \brief The radial variables from which the hydrodynamic quantities and
      83             :   /// spacetime metric can be computed.
      84             :   ///
      85             :   /// For radii greater than the `outer_radius()`, this returns the appropriate
      86             :   /// vacuum spacetime.
      87             :   ///
      88             :   /// \note This solution of the TOV equations is a function of areal radius.
      89             :   template <typename DataType>
      90             :   RelativisticEuler::Solutions::TovStar<TovSolution>::RadialVariables<DataType>
      91           1 :   radial_variables(
      92             :       const EquationsOfState::EquationOfState<true, 1>& equation_of_state,
      93             :       const tnsr::I<DataType, 3>& x) const noexcept;
      94             : 
      95             :   // NOLINTNEXTLINE(google-runtime-references)
      96           0 :   void pup(PUP::er& p) noexcept;
      97             : 
      98             :  private:
      99           0 :   double outer_radius_{std::numeric_limits<double>::signaling_NaN()};
     100           0 :   double total_mass_{std::numeric_limits<double>::signaling_NaN()};
     101           0 :   double log_lapse_at_outer_radius_{
     102             :       std::numeric_limits<double>::signaling_NaN()};
     103           0 :   intrp::BarycentricRational mass_over_radius_interpolant_;
     104           0 :   intrp::BarycentricRational log_enthalpy_interpolant_;
     105             : };
     106             : 
     107             : }  // namespace Solutions
     108             : }  // namespace gr