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