Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <cstddef>
       7             : 
       8             : #include "DataStructures/Tensor/Tensor.hpp"
       9             : #include "Evolution/Systems/NewtonianEuler/Sources/Source.hpp"
      10             : #include "Evolution/Systems/NewtonianEuler/TagsDeclarations.hpp"
      11             : #include "Options/String.hpp"
      12             : #include "Utilities/TMPL.hpp"
      13             : 
      14             : /// \cond
      15             : class DataVector;
      16             : 
      17             : namespace NewtonianEuler {
      18             : namespace Solutions {
      19             : struct LaneEmdenStar;
      20             : }  // namespace Solutions
      21             : }  // namespace NewtonianEuler
      22             : 
      23             : namespace PUP {
      24             : class er;
      25             : }  // namespace PUP
      26             : 
      27             : namespace Tags {
      28             : template <typename SolutionType>
      29             : struct AnalyticSolution;
      30             : }  // namespace Tags
      31             : 
      32             : namespace domain::Tags {
      33             : template <size_t Dim, typename Frame>
      34             : struct Coordinates;
      35             : }  // namespace domain::Tags
      36             : 
      37             : namespace gsl {
      38             : template <typename T>
      39             : class not_null;
      40             : }  // namespace gsl
      41             : /// \endcond
      42             : 
      43             : namespace NewtonianEuler::Sources {
      44             : 
      45             : /*!
      46             :  * \brief Source giving the acceleration due to gravity in the spherical,
      47             :  * Newtonian Lane-Emden star solution.
      48             :  *
      49             :  * The gravitational field \f$g^i\f$ enters the NewtonianEuler system as source
      50             :  * terms for the conserved momentum and energy:
      51             :  *
      52             :  * \f{align*}
      53             :  * \partial_t S^i + \partial_j F^{j}(S^i) &= S(S^i) = \rho g^i
      54             :  * \partial_t e + \partial_j F^{j}(e) &= S(e) = S_i g^i,
      55             :  * \f}
      56             :  *
      57             :  * where \f$S^i\f$ is the conserved momentum density, \f$e\f$ is the conserved
      58             :  * energy, \f$F^{j}(u)\f$ is the flux of the conserved quantity \f$u\f$, and
      59             :  * \f$\rho\f$ is the fluid mass density.
      60             :  *
      61             :  * \note This source is specialized to the Lane-Emden solution because it
      62             :  * queries a LaneEmdenStar analytic solution for the gravitational field that
      63             :  * generates the fluid acceleration. This source does *not* integrate the fluid
      64             :  * density to compute a self-consistent gravitational field (i.e., as if one
      65             :  * were solving a coupled Euler + Poisson system).
      66             :  */
      67           1 : class LaneEmdenGravitationalField : public Source<3> {
      68             :  public:
      69             :   /// The central mass density of the star.
      70           1 :   struct CentralMassDensity {
      71           0 :     using type = double;
      72           0 :     static constexpr Options::String help = {
      73             :         "The central mass density of the star."};
      74           0 :     static type lower_bound() { return 0.; }
      75             :   };
      76             : 
      77             :   /// The polytropic constant of the polytropic fluid.
      78           1 :   struct PolytropicConstant {
      79           0 :     using type = double;
      80           0 :     static constexpr Options::String help = {
      81             :         "The polytropic constant of the fluid."};
      82           0 :     static type lower_bound() { return 0.; }
      83             :   };
      84             : 
      85           0 :   using options = tmpl::list<CentralMassDensity, PolytropicConstant>;
      86             : 
      87           0 :   static constexpr Options::String help = {
      88             :       "The gravitational field corresponding to a static, "
      89             :       "spherically-symmetric star in Newtonian gravity, found by "
      90             :       "solving the Lane-Emden equations, with a given central density and "
      91             :       "polytropic fluid. The fluid has polytropic index 1, but the polytropic "
      92             :       "constant is specifiable"};
      93             : 
      94           0 :   LaneEmdenGravitationalField(double central_mass_density,
      95             :                               double polytropic_constant);
      96             : 
      97           0 :   LaneEmdenGravitationalField() = default;
      98           0 :   LaneEmdenGravitationalField(const LaneEmdenGravitationalField& /*rhs*/) =
      99             :       default;
     100           0 :   LaneEmdenGravitationalField& operator=(
     101             :       const LaneEmdenGravitationalField& /*rhs*/) = default;
     102           0 :   LaneEmdenGravitationalField(LaneEmdenGravitationalField&& /*rhs*/) = default;
     103           0 :   LaneEmdenGravitationalField& operator=(
     104             :       LaneEmdenGravitationalField&& /*rhs*/) = default;
     105           0 :   ~LaneEmdenGravitationalField() override = default;
     106             : 
     107             :   /// \cond
     108             :   explicit LaneEmdenGravitationalField(CkMigrateMessage* msg);
     109             :   using PUP::able::register_constructor;
     110             :   WRAPPED_PUPable_decl_template(LaneEmdenGravitationalField);
     111             :   /// \endcond
     112             : 
     113             :   // NOLINTNEXTLINE(google-runtime-references)
     114           0 :   void pup(PUP::er& p) override;
     115             : 
     116           0 :   auto get_clone() const -> std::unique_ptr<Source> override;
     117             : 
     118           0 :   void operator()(
     119             :       gsl::not_null<Scalar<DataVector>*> source_mass_density_cons,
     120             :       gsl::not_null<tnsr::I<DataVector, 3>*> source_momentum_density,
     121             :       gsl::not_null<Scalar<DataVector>*> source_energy_density,
     122             :       const Scalar<DataVector>& mass_density_cons,
     123             :       const tnsr::I<DataVector, 3>& momentum_density,
     124             :       const Scalar<DataVector>& energy_density,
     125             :       const tnsr::I<DataVector, 3>& velocity,
     126             :       const Scalar<DataVector>& pressure,
     127             :       const Scalar<DataVector>& specific_internal_energy,
     128             :       const EquationsOfState::EquationOfState<false, 2>& eos,
     129             :       const tnsr::I<DataVector, 3>& coords, double time) const override;
     130             : 
     131             :  private:
     132           0 :   tnsr::I<DataVector, 3> gravitational_field(
     133             :       const tnsr::I<DataVector, 3>& x) const;
     134             : 
     135           0 :   double central_mass_density_ = std::numeric_limits<double>::signaling_NaN();
     136           0 :   double polytropic_constant_ = std::numeric_limits<double>::signaling_NaN();
     137             : };
     138             : }  // namespace NewtonianEuler::Sources