SpECTRE Documentation Coverage Report
Current view: top level - PointwiseFunctions/Hydro/EquationsOfState - EquationOfState.hpp Hit Total Coverage
Commit: c3e43f8d41800b0ecefb9d1393f1de1d5a280c8f Lines: 70 107 65.4 %
Date: 2026-07-24 22:09:25
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 <boost/preprocessor/arithmetic/sub.hpp>
       7             : #include <boost/preprocessor/list/for_each.hpp>
       8             : #include <boost/preprocessor/punctuation/comma_if.hpp>
       9             : #include <boost/preprocessor/repetition/repeat.hpp>
      10             : #include <boost/preprocessor/tuple/enum.hpp>
      11             : #include <boost/preprocessor/tuple/to_list.hpp>
      12             : 
      13             : #include "DataStructures/DataVector.hpp"
      14             : #include "DataStructures/Tensor/Tensor.hpp"
      15             : #include "PointwiseFunctions/Hydro/Units.hpp"
      16             : #include "Utilities/CallWithDynamicType.hpp"
      17             : #include "Utilities/Serialization/CharmPupable.hpp"
      18             : #include "Utilities/TMPL.hpp"
      19             : #include "Utilities/TypeTraits.hpp"
      20             : 
      21             : /// \cond
      22             : namespace EquationsOfState {
      23             : template <typename ColdEos>
      24             : class Barotropic2D;
      25             : template <typename ColdEquilEos>
      26             : class Barotropic3D;
      27             : template <typename EquilEos>
      28             : class Equilibrium3D;
      29             : template <typename ColdEquationOfState>
      30             : class HybridEos;
      31             : template <bool IsRelativistic>
      32             : class IdealFluid;
      33             : template <bool IsRelativistic>
      34             : class PolytropicFluid;
      35             : template <bool IsRelativistic>
      36             : class PiecewisePolytropicFluid;
      37             : class Spectral;
      38             : template <typename LowDensityEoS>
      39             : class Enthalpy;
      40             : template <bool IsRelativistic>
      41             : class Tabulated3D;
      42             : }  // namespace EquationsOfState
      43             : /// \endcond
      44             : 
      45             : /// Contains all equations of state, including base class
      46             : namespace EquationsOfState {
      47             : 
      48             : namespace detail {
      49             : template <bool IsRelativistic, size_t ThermodynamicDim>
      50             : struct DerivedClasses {};
      51             : 
      52             : template <>
      53             : struct DerivedClasses<true, 1> {
      54             :   using type = tmpl::list<PolytropicFluid<true>, PiecewisePolytropicFluid<true>,
      55             :                           Spectral, Enthalpy<PolytropicFluid<true>>,
      56             :                           Enthalpy<Enthalpy<PolytropicFluid<true>>>,
      57             :                           Enthalpy<Enthalpy<Enthalpy<PolytropicFluid<true>>>>,
      58             :                           Enthalpy<Spectral>, Enthalpy<Enthalpy<Spectral>>,
      59             :                           Enthalpy<Enthalpy<Enthalpy<Spectral>>>>;
      60             : };
      61             : 
      62             : template <>
      63             : struct DerivedClasses<false, 1> {
      64             :   using type =
      65             :       tmpl::list<PolytropicFluid<false>, PiecewisePolytropicFluid<false>>;
      66             : };
      67             : 
      68             : template <>
      69             : struct DerivedClasses<true, 2> {
      70             :   using type = tmpl::list<
      71             :       IdealFluid<true>, Barotropic2D<PolytropicFluid<true>>,
      72             :       Barotropic2D<PiecewisePolytropicFluid<true>>, Barotropic2D<Spectral>,
      73             :       Barotropic2D<Enthalpy<PolytropicFluid<true>>>,
      74             :       Barotropic2D<Enthalpy<Enthalpy<Enthalpy<PolytropicFluid<true>>>>>,
      75             :       Barotropic2D<Enthalpy<Spectral>>,
      76             :       Barotropic2D<Enthalpy<Enthalpy<Spectral>>>,
      77             :       Barotropic2D<Enthalpy<Enthalpy<Enthalpy<Spectral>>>>,
      78             :       HybridEos<PolytropicFluid<true>>, HybridEos<Spectral>,
      79             :       HybridEos<Enthalpy<PolytropicFluid<true>>>,
      80             :       HybridEos<Enthalpy<Enthalpy<Enthalpy<PolytropicFluid<true>>>>>,
      81             :       HybridEos<Enthalpy<Spectral>>, HybridEos<Enthalpy<Enthalpy<Spectral>>>,
      82             :       HybridEos<Enthalpy<Enthalpy<Enthalpy<Spectral>>>>>;
      83             : };
      84             : 
      85             : template <>
      86             : struct DerivedClasses<false, 2> {
      87             :   using type =
      88             :       tmpl::list<IdealFluid<false>, Barotropic2D<PolytropicFluid<false>>,
      89             :                  Barotropic2D<PiecewisePolytropicFluid<false>>,
      90             :                  HybridEos<PolytropicFluid<false>>>;
      91             : };
      92             : 
      93             : template <>
      94             : struct DerivedClasses<true, 3> {
      95             :   using type = tmpl::list<
      96             :       Tabulated3D<true>, Barotropic3D<PolytropicFluid<true>>,
      97             :       Barotropic3D<PiecewisePolytropicFluid<true>>, Barotropic3D<Spectral>,
      98             :       Barotropic3D<Enthalpy<PolytropicFluid<true>>>,
      99             :       Barotropic3D<Enthalpy<Enthalpy<Enthalpy<PolytropicFluid<true>>>>>,
     100             :       Barotropic3D<Enthalpy<Spectral>>,
     101             :       Barotropic3D<Enthalpy<Enthalpy<Spectral>>>,
     102             :       Barotropic3D<Enthalpy<Enthalpy<Enthalpy<Spectral>>>>,
     103             :       Equilibrium3D<IdealFluid<true>>,
     104             :       Equilibrium3D<HybridEos<PolytropicFluid<true>>>,
     105             :       Equilibrium3D<HybridEos<Enthalpy<PolytropicFluid<true>>>>,
     106             :       Equilibrium3D<
     107             :           HybridEos<Enthalpy<Enthalpy<Enthalpy<PolytropicFluid<true>>>>>>,
     108             :       Equilibrium3D<HybridEos<Enthalpy<Spectral>>>,
     109             :       Equilibrium3D<HybridEos<Enthalpy<Enthalpy<Spectral>>>>,
     110             :       Equilibrium3D<HybridEos<Enthalpy<Enthalpy<Enthalpy<Spectral>>>>>,
     111             :       Equilibrium3D<HybridEos<Spectral>>>;
     112             : };
     113             : 
     114             : template <>
     115             : struct DerivedClasses<false, 3> {
     116             :   using type =
     117             :       tmpl::list<Tabulated3D<false>, Barotropic3D<PolytropicFluid<false>>,
     118             :                  Barotropic3D<PiecewisePolytropicFluid<false>>,
     119             :                  Equilibrium3D<IdealFluid<false>>,
     120             :                  Equilibrium3D<HybridEos<PolytropicFluid<false>>>>;
     121             : };
     122             : 
     123             : }  // namespace detail
     124             : 
     125             : /*!
     126             :  * \ingroup EquationsOfStateGroup
     127             :  * \brief Base class for equations of state depending on whether or not the
     128             :  * system is relativistic, and the number of independent thermodynamic variables
     129             :  * (`ThermodynamicDim`) needed to determine the pressure.
     130             :  *
     131             :  * The template parameter `IsRelativistic` is `true` for relativistic equations
     132             :  * of state and `false` for non-relativistic equations of state.
     133             :  */
     134             : template <bool IsRelativistic, size_t ThermodynamicDim>
     135           1 : class EquationOfState;
     136             : 
     137             : template <typename T>
     138           0 : struct get_eos_base_impl {
     139           0 :   using type = EquationsOfState::EquationOfState<T::is_relativistic,
     140             :                                                  T::thermodynamic_dim>;
     141             : };
     142             : 
     143             : template <bool IsRelativistic, size_t ThermodynamicDim>
     144           0 : struct get_eos_base_impl<
     145             :     EquationsOfState::EquationOfState<IsRelativistic, ThermodynamicDim>> {
     146           0 :   using type =
     147             :       EquationsOfState::EquationOfState<IsRelativistic, ThermodynamicDim>;
     148             : };
     149             : 
     150             : template <typename T>
     151           0 : using get_eos_base = typename get_eos_base_impl<T>::type;
     152             : 
     153             : /*!
     154             :  * \ingroup EquationsOfStateGroup
     155             :  * \brief Base class for equations of state which need one thermodynamic
     156             :  * variable in order to determine the pressure.
     157             :  *
     158             :  * The template parameter `IsRelativistic` is `true` for relativistic equations
     159             :  * of state and `false` for non-relativistic equations of state.
     160             :  */
     161             : template <bool IsRelativistic>
     162           1 : class EquationOfState<IsRelativistic, 1> : public PUP::able {
     163             :  public:
     164           0 :   static constexpr bool is_relativistic = IsRelativistic;
     165           0 :   static constexpr size_t thermodynamic_dim = 1;
     166           0 :   using creatable_classes =
     167             :       typename detail::DerivedClasses<IsRelativistic, 1>::type;
     168             : 
     169           0 :   EquationOfState() = default;
     170           0 :   EquationOfState(const EquationOfState&) = default;
     171           0 :   EquationOfState& operator=(const EquationOfState&) = default;
     172           0 :   EquationOfState(EquationOfState&&) = default;
     173           0 :   EquationOfState& operator=(EquationOfState&&) = default;
     174           0 :   ~EquationOfState() override = default;
     175             : 
     176           0 :   explicit EquationOfState(CkMigrateMessage* msg) : PUP::able(msg) {}
     177             : 
     178           0 :   WRAPPED_PUPable_abstract(EquationOfState);  // NOLINT
     179             : 
     180           0 :   virtual std::unique_ptr<EquationOfState<IsRelativistic, 1>> get_clone()
     181             :       const = 0;
     182             : 
     183           0 :   virtual bool is_equal(
     184             :       const EquationOfState<IsRelativistic, 1>& rhs) const = 0;
     185             : 
     186             :   /// Create a 3D EOS from the 1D EOS
     187             :   virtual std::unique_ptr<EquationOfState<IsRelativistic, 3>>
     188           1 :   promote_to_3d_eos() const = 0;
     189             : 
     190             :   /// Create a 2D EOS from the 1D EOS
     191             :   virtual std::unique_ptr<EquationOfState<IsRelativistic, 2>>
     192           1 :   promote_to_2d_eos() const = 0;
     193             : 
     194             :   /// \brief Returns `true` if the EOS is barotropic
     195           1 :   bool is_barotropic() const { return true; }
     196             :   /// @{
     197             :   /*!
     198             :    * Computes the electron fraction in beta-equilibrium \f$Y_e^{\rm eq}\f$ from
     199             :    * the rest mass density \f$\rho\f$.
     200             :    */
     201           1 :   virtual Scalar<double> equilibrium_electron_fraction_from_density_temperature(
     202             :       const Scalar<double>& rest_mass_density,
     203             :       const Scalar<double>& /*temperature*/) const {
     204             :     return make_with_value<Scalar<double>>(rest_mass_density, 0.1);
     205             :   }
     206             : 
     207             :   virtual Scalar<DataVector>
     208           1 :   equilibrium_electron_fraction_from_density_temperature(
     209             :       const Scalar<DataVector>& rest_mass_density,
     210             :       const Scalar<DataVector>& /*temperature*/) const {
     211             :     return make_with_value<Scalar<DataVector>>(rest_mass_density, 0.1);
     212             :   }
     213             :   /// @}
     214             : 
     215             :   /// @{
     216             :   /*!
     217             :    * Computes the pressure \f$p\f$ from the rest mass density \f$\rho\f$.
     218             :    */
     219           1 :   virtual Scalar<double> pressure_from_density(
     220             :       const Scalar<double>& /*rest_mass_density*/) const = 0;
     221           1 :   virtual Scalar<DataVector> pressure_from_density(
     222             :       const Scalar<DataVector>& /*rest_mass_density*/) const = 0;
     223             :   /// @}
     224             : 
     225             :   /// @{
     226             :   /*!
     227             :    * Computes the rest mass density \f$\rho\f$ from the specific enthalpy
     228             :    * \f$h\f$.
     229             :    */
     230           1 :   virtual Scalar<double> rest_mass_density_from_enthalpy(
     231             :       const Scalar<double>& /*specific_enthalpy*/) const = 0;
     232           1 :   virtual Scalar<DataVector> rest_mass_density_from_enthalpy(
     233             :       const Scalar<DataVector>& /*specific_enthalpy*/) const = 0;
     234             :   /// @}
     235             : 
     236             :   /// @{
     237             :   /*!
     238             :    * Computes the specific entropy \f$s\f$ from the rest mass density
     239             :    * \f$\rho\f$.
     240             :    */
     241           1 :   virtual Scalar<double> specific_entropy_from_density(
     242             :       const Scalar<double>& /*rest_mass_density*/) const {
     243             :     return Scalar<double>{0.0};
     244             :   }
     245           1 :   virtual Scalar<DataVector> specific_entropy_from_density(
     246             :       const Scalar<DataVector>& rest_mass_density) const {
     247             :     return make_with_value<Scalar<DataVector>>(rest_mass_density, 0.0);
     248             :   }
     249             :   /// @}
     250             : 
     251             :   /// @{
     252             :   /*!
     253             :    * Computes the specific entropy \f$s\f$ from the specific internal energy
     254             :    * \f$\epsilon\f$.
     255             :    */
     256           1 :   virtual Scalar<double> specific_entropy_from_specific_internal_energy(
     257             :       const Scalar<double>& /*specific_internal_energy*/) const {
     258             :     return Scalar<double>{0.0};
     259             :   }
     260           1 :   virtual Scalar<DataVector> specific_entropy_from_specific_internal_energy(
     261             :       const Scalar<DataVector>& specific_internal_energy) const {
     262             :     return make_with_value<Scalar<DataVector>>(specific_internal_energy, 0.0);
     263             :   }
     264             :   /// @}
     265             : 
     266             :   /// @{
     267             :   /*!
     268             :    * Computes the specific internal energy \f$\epsilon\f$ from the rest mass
     269             :    * density \f$\rho\f$.
     270             :    */
     271           1 :   virtual Scalar<double> specific_internal_energy_from_density(
     272             :       const Scalar<double>& /*rest_mass_density*/) const = 0;
     273           1 :   virtual Scalar<DataVector> specific_internal_energy_from_density(
     274             :       const Scalar<DataVector>& /*rest_mass_density*/) const = 0;
     275             :   /// @}
     276             : 
     277             :   /// @{
     278             :   /*!
     279             :    * Computes the temperature \f$T\f$ from the rest mass
     280             :    * density \f$\rho\f$.
     281             :    */
     282           1 :   virtual Scalar<double> temperature_from_density(
     283             :       const Scalar<double>& /*rest_mass_density*/) const {
     284             :     return Scalar<double>{0.0};
     285             :   }
     286           1 :   virtual Scalar<DataVector> temperature_from_density(
     287             :       const Scalar<DataVector>& rest_mass_density) const {
     288             :     return make_with_value<Scalar<DataVector>>(rest_mass_density, 0.0);
     289             :   }
     290             :   /// @}
     291             : 
     292             :   /// @{
     293             :   /*!
     294             :    * Computes the temperature \f$\T\f$ from the specific internal energy
     295             :    * \f$\epsilon\f$.
     296             :    */
     297           1 :   virtual Scalar<double> temperature_from_specific_internal_energy(
     298             :       const Scalar<double>& /*specific_internal_energy*/) const {
     299             :     return Scalar<double>{0.0};
     300             :   }
     301           1 :   virtual Scalar<DataVector> temperature_from_specific_internal_energy(
     302             :       const Scalar<DataVector>& specific_internal_energy) const {
     303             :     return make_with_value<Scalar<DataVector>>(specific_internal_energy, 0.0);
     304             :   }
     305             :   /// @}
     306             : 
     307             :   /// @{
     308             :   /*!
     309             :    * Computes \f$\chi=\partial p / \partial \rho\f$ from \f$\rho\f$, where
     310             :    * \f$p\f$ is the pressure and \f$\rho\f$ is the rest mass density.
     311             :    */
     312           1 :   virtual Scalar<double> chi_from_density(
     313             :       const Scalar<double>& /*rest_mass_density*/) const = 0;
     314           1 :   virtual Scalar<DataVector> chi_from_density(
     315             :       const Scalar<DataVector>& /*rest_mass_density*/) const = 0;
     316             :   /// @}
     317             : 
     318             :   /// @{
     319             :   /*!
     320             :    * Computes \f$\kappa p/\rho^2=(p/\rho^2)\partial p / \partial \epsilon\f$
     321             :    * from \f$\rho\f$, where \f$p\f$ is the pressure, \f$\rho\f$ is the rest mass
     322             :    * density, and \f$\epsilon\f$ is the specific internal energy.
     323             :    *
     324             :    * The reason for not returning just
     325             :    * \f$\kappa=\partial p / \partial \epsilon\f$ is to avoid division by zero
     326             :    * for small values of \f$\rho\f$ when assembling the speed of sound with
     327             :    * some equations of state.
     328             :    */
     329           1 :   virtual Scalar<double> kappa_times_p_over_rho_squared_from_density(
     330             :       const Scalar<double>& /*rest_mass_density*/) const = 0;
     331           1 :   virtual Scalar<DataVector> kappa_times_p_over_rho_squared_from_density(
     332             :       const Scalar<DataVector>& /*rest_mass_density*/) const = 0;
     333             : 
     334             :   /// The lower bound of the electron fraction that is valid for this EOS
     335           1 :   virtual double electron_fraction_lower_bound() const { return 0.0; }
     336             : 
     337             :   /// The upper bound of the electron fraction that is valid for this EOS
     338           1 :   virtual double electron_fraction_upper_bound() const { return 1.0; }
     339             : 
     340             :   /// The lower bound of the rest mass density that is valid for this EOS
     341           1 :   virtual double rest_mass_density_lower_bound() const = 0;
     342             : 
     343             :   /// The upper bound of the rest mass density that is valid for this EOS
     344           1 :   virtual double rest_mass_density_upper_bound() const = 0;
     345             : 
     346             :   /// The lower bound of the temperature that is valid for this EOS
     347           1 :   virtual double temperature_lower_bound() const { return 0.0; };
     348             : 
     349             :   /// The upper bound of the temperature that is valid for this EOS
     350           1 :   virtual double temperature_upper_bound() const {
     351             :     return std::numeric_limits<double>::max();
     352             :   };
     353             : 
     354             :   /// The lower bound of the specific internal energy that is valid for this EOS
     355           1 :   virtual double specific_internal_energy_lower_bound() const { return 0.0; };
     356             : 
     357             :   /// The upper bound of the specific internal energy that is valid for this EOS
     358           1 :   virtual double specific_internal_energy_upper_bound() const {
     359             :     return std::numeric_limits<double>::max();
     360             :   };
     361             : 
     362             :   /// The lower bound of the specific enthalpy that is valid for this EOS
     363           1 :   virtual double specific_enthalpy_lower_bound() const = 0;
     364             : 
     365             :   /// The vacuum mass of a baryon for this EOS
     366           1 :   virtual double baryon_mass() const {
     367             :     return hydro::units::geometric::default_baryon_mass;
     368             :   }
     369             : };
     370             : 
     371             : /*!
     372             :  * \ingroup EquationsOfStateGroup
     373             :  * \brief Base class for equations of state which need two independent
     374             :  * thermodynamic variables in order to determine the pressure.
     375             :  *
     376             :  * The template parameter `IsRelativistic` is `true` for relativistic equations
     377             :  * of state and `false` for non-relativistic equations of state.
     378             :  */
     379             : template <bool IsRelativistic>
     380           1 : class EquationOfState<IsRelativistic, 2> : public PUP::able {
     381             :  public:
     382           0 :   static constexpr bool is_relativistic = IsRelativistic;
     383           0 :   static constexpr size_t thermodynamic_dim = 2;
     384           0 :   using creatable_classes =
     385             :       typename detail::DerivedClasses<IsRelativistic, 2>::type;
     386             : 
     387           0 :   EquationOfState() = default;
     388           0 :   EquationOfState(const EquationOfState&) = default;
     389           0 :   EquationOfState& operator=(const EquationOfState&) = default;
     390           0 :   EquationOfState(EquationOfState&&) = default;
     391           0 :   EquationOfState& operator=(EquationOfState&&) = default;
     392           0 :   ~EquationOfState() override = default;
     393             : 
     394           0 :   explicit EquationOfState(CkMigrateMessage* msg) : PUP::able(msg) {}
     395             : 
     396           0 :   WRAPPED_PUPable_abstract(EquationOfState);  // NOLINT
     397             : 
     398           0 :   virtual inline std::unique_ptr<EquationOfState<IsRelativistic, 2>> get_clone()
     399             :       const = 0;
     400             : 
     401           0 :   virtual bool is_equal(
     402             :       const EquationOfState<IsRelativistic, 2>& rhs) const = 0;
     403             : 
     404             :   virtual std::unique_ptr<EquationOfState<IsRelativistic, 2>>
     405           0 :   promote_to_2d_eos() const {
     406             :     return this->get_clone();
     407             :   }
     408             : 
     409             :   virtual std::unique_ptr<EquationOfState<IsRelativistic, 3>>
     410           0 :   promote_to_3d_eos() const = 0;
     411             : 
     412             :   /// \brief Returns `true` if the EOS is barotropic
     413           1 :   virtual bool is_barotropic() const = 0;
     414             : 
     415             :   /// \brief Returns `true` if the EOS is in beta-equilibrium
     416           1 :   virtual bool is_equilibrium() const { return true; }
     417             : 
     418             :   /// @{
     419             :   /*!
     420             :    * Computes the electron fraction in beta-equilibrium \f$Y_e^{\rm eq}\f$ from
     421             :    * the rest mass density \f$\rho\f$ and the temperature \f$T\f$.
     422             :    */
     423           1 :   virtual Scalar<double> equilibrium_electron_fraction_from_density_temperature(
     424             :       const Scalar<double>& rest_mass_density,
     425             :       const Scalar<double>& /*temperature*/) const {
     426             :     return make_with_value<Scalar<double>>(rest_mass_density, 0.1);
     427             :   }
     428             : 
     429             :   virtual Scalar<DataVector>
     430           1 :   equilibrium_electron_fraction_from_density_temperature(
     431             :       const Scalar<DataVector>& rest_mass_density,
     432             :       const Scalar<DataVector>& /*temperature*/) const {
     433             :     return make_with_value<Scalar<DataVector>>(rest_mass_density, 0.1);
     434             :   }
     435             :   /// @}
     436             : 
     437             :   /// @{
     438             :   /*!
     439             :    * Computes the pressure \f$p\f$ from the rest mass density \f$\rho\f$ and the
     440             :    * specific internal energy \f$\epsilon\f$.
     441             :    */
     442           1 :   virtual Scalar<double> pressure_from_density_and_energy(
     443             :       const Scalar<double>& /*rest_mass_density*/,
     444             :       const Scalar<double>& /*specific_internal_energy*/) const = 0;
     445           1 :   virtual Scalar<DataVector> pressure_from_density_and_energy(
     446             :       const Scalar<DataVector>& /*rest_mass_density*/,
     447             :       const Scalar<DataVector>& /*specific_internal_energy*/) const = 0;
     448             :   /// @}
     449             : 
     450             :   /// @{
     451             :   /*!
     452             :    * Computes the pressure \f$p\f$ from the rest mass density \f$\rho\f$ and the
     453             :    * specific enthalpy \f$h\f$.
     454             :    */
     455           1 :   virtual Scalar<double> pressure_from_density_and_enthalpy(
     456             :       const Scalar<double>& /*rest_mass_density*/,
     457             :       const Scalar<double>& /*specific_enthalpy*/) const = 0;
     458           1 :   virtual Scalar<DataVector> pressure_from_density_and_enthalpy(
     459             :       const Scalar<DataVector>& /*rest_mass_density*/,
     460             :       const Scalar<DataVector>& /*specific_enthalpy*/) const = 0;
     461             :   /// @}
     462             : 
     463             :   /// @{
     464             :   /*!
     465             :    * Computes the specific entropy \f$s\f$ from the rest mass density \f$\rho\f$
     466             :    * and the specific internal energy \f$\epsilon\f$.
     467             :    */
     468           1 :   virtual Scalar<double> specific_entropy_from_density_and_energy(
     469             :       const Scalar<double>& /*rest_mass_density*/,
     470             :       const Scalar<double>& /*specific_internal_energy*/) const = 0;
     471             : 
     472           1 :   virtual Scalar<DataVector> specific_entropy_from_density_and_energy(
     473             :       const Scalar<DataVector>& /*rest_mass_density*/,
     474             :       const Scalar<DataVector>& /*specific_internal_energy*/) const = 0;
     475             :   /// @}
     476             : 
     477             :   /// @{
     478             :   /*!
     479             :    * Computes the specific entropy \f$s\f$ from the rest mass density \f$\rho\f$
     480             :    * and the temperature \f$T\f$.
     481             :    */
     482           1 :   virtual Scalar<double> specific_entropy_from_density_and_temperature(
     483             :       const Scalar<double>& /*rest_mass_density*/,
     484             :       const Scalar<double>& /*temperature*/) const = 0;
     485             : 
     486           1 :   virtual Scalar<DataVector> specific_entropy_from_density_and_temperature(
     487             :       const Scalar<DataVector>& /*rest_mass_density*/,
     488             :       const Scalar<DataVector>& /*temperature*/) const = 0;
     489             :   /// @}
     490             : 
     491             :   /// @{
     492             :   /*!
     493             :    * Computes the specific internal energy \f$\epsilon\f$ from the rest mass
     494             :    * density \f$\rho\f$ and the pressure \f$p\f$.
     495             :    */
     496           1 :   virtual Scalar<double> specific_internal_energy_from_density_and_pressure(
     497             :       const Scalar<double>& /*rest_mass_density*/,
     498             :       const Scalar<double>& /*pressure*/) const = 0;
     499           1 :   virtual Scalar<DataVector> specific_internal_energy_from_density_and_pressure(
     500             :       const Scalar<DataVector>& /*rest_mass_density*/,
     501             :       const Scalar<DataVector>& /*pressure*/) const = 0;
     502             :   /// @}
     503             : 
     504             :   /// @{
     505             :   /*!
     506             :    * Computes the temperature \f$T\f$ from the rest mass
     507             :    * density \f$\rho\f$ and the specific internal energy \f$\epsilon\f$.
     508             :    */
     509           1 :   virtual Scalar<double> temperature_from_density_and_energy(
     510             :       const Scalar<double>& /*rest_mass_density*/,
     511             :       const Scalar<double>& /*specific_internal_energy*/) const = 0;
     512           1 :   virtual Scalar<DataVector> temperature_from_density_and_energy(
     513             :       const Scalar<DataVector>& /*rest_mass_density*/,
     514             :       const Scalar<DataVector>& /*specific_internal_energy*/) const = 0;
     515             :   /// @}
     516             : 
     517             :   /// @{
     518             :   /*!
     519             :    * Computes the specific internal energy \f$\epsilon\f$ from the rest mass
     520             :    * density \f$\rho\f$ and the temperature \f$T\f$.
     521             :    */
     522           1 :   virtual Scalar<double> specific_internal_energy_from_density_and_temperature(
     523             :       const Scalar<double>& /*rest_mass_density*/,
     524             :       const Scalar<double>& /*temperature*/) const = 0;
     525             :   virtual Scalar<DataVector>
     526           1 :   specific_internal_energy_from_density_and_temperature(
     527             :       const Scalar<DataVector>& /*rest_mass_density*/,
     528             :       const Scalar<DataVector>& /*temperature*/) const = 0;
     529             :   /// @}
     530             : 
     531             :   /// @{
     532             :   /*!
     533             :    * Computes \f$\chi=\partial p / \partial \rho |_{\epsilon}\f$ from the
     534             :    * \f$\rho\f$ and \f$\epsilon\f$, where \f$p\f$ is the pressure, \f$\rho\f$ is
     535             :    * the rest mass density, and \f$\epsilon\f$ is the specific internal energy.
     536             :    */
     537           1 :   virtual Scalar<double> chi_from_density_and_energy(
     538             :       const Scalar<double>& /*rest_mass_density*/,
     539             :       const Scalar<double>& /*specific_internal_energy*/) const = 0;
     540           1 :   virtual Scalar<DataVector> chi_from_density_and_energy(
     541             :       const Scalar<DataVector>& /*rest_mass_density*/,
     542             :       const Scalar<DataVector>& /*specific_internal_energy*/) const = 0;
     543             :   /// @}
     544             : 
     545             :   /// @{
     546             :   /*!
     547             :    * Computes \f$\kappa p/\rho^2=(p/\rho^2)\partial p / \partial \epsilon
     548             :    * |_{\rho}\f$ from \f$\rho\f$ and \f$\epsilon\f$, where \f$p\f$ is the
     549             :    * pressure, \f$\rho\f$ is the rest mass density, and \f$\epsilon\f$ is the
     550             :    * specific internal energy.
     551             :    *
     552             :    * The reason for not returning just
     553             :    * \f$\kappa=\partial p / \partial \epsilon\f$ is to avoid division by zero
     554             :    * for small values of \f$\rho\f$ when assembling the speed of sound with
     555             :    * some equations of state.
     556             :    */
     557           1 :   virtual Scalar<double> kappa_times_p_over_rho_squared_from_density_and_energy(
     558             :       const Scalar<double>& /*rest_mass_density*/,
     559             :       const Scalar<double>& /*specific_internal_energy*/) const = 0;
     560             :   virtual Scalar<DataVector>
     561           1 :   kappa_times_p_over_rho_squared_from_density_and_energy(
     562             :       const Scalar<DataVector>& /*rest_mass_density*/,
     563             :       const Scalar<DataVector>& /*specific_internal_energy*/) const = 0;
     564             :   /// @}
     565             : 
     566             :   /// The lower bound of the electron fraction that is valid for this EOS
     567           1 :   virtual double electron_fraction_lower_bound() const { return 0.0; }
     568             : 
     569             :   /// The upper bound of the electron fraction that is valid for this EOS
     570           1 :   virtual double electron_fraction_upper_bound() const { return 1.0; }
     571             : 
     572             :   /// The lower bound of the rest mass density that is valid for this EOS
     573           1 :   virtual double rest_mass_density_lower_bound() const = 0;
     574             : 
     575             :   /// The upper bound of the rest mass density that is valid for this EOS
     576           1 :   virtual double rest_mass_density_upper_bound() const = 0;
     577             : 
     578             :   /// The lower bound of the specific internal energy that is valid for this EOS
     579             :   /// at the given rest mass density \f$\rho\f$
     580           1 :   virtual double specific_internal_energy_lower_bound(
     581             :       const double rest_mass_density) const = 0;
     582             : 
     583             :   /// The upper bound of the specific internal energy that is valid for this EOS
     584             :   /// at the given rest mass density \f$\rho\f$
     585           1 :   virtual double specific_internal_energy_upper_bound(
     586             :       const double rest_mass_density) const = 0;
     587             : 
     588             :   /// The lower bound of the temperature that is valid for this EOS
     589           1 :   virtual double temperature_lower_bound() const { return 0.0; };
     590             : 
     591             :   /// The upper bound of the temperature that is valid for this EOS
     592           1 :   virtual double temperature_upper_bound() const {
     593             :     return std::numeric_limits<double>::max();
     594             :   };
     595             : 
     596             :   /// The lower bound of the specific enthalpy that is valid for this EOS
     597           1 :   virtual double specific_enthalpy_lower_bound() const = 0;
     598             : 
     599             :   /// The vacuum mass of a baryon for this EOS
     600           1 :   virtual double baryon_mass() const {
     601             :     return hydro::units::geometric::default_baryon_mass;
     602             :   }
     603             : };
     604             : 
     605             : /*!
     606             :  * \ingroup EquationsOfStateGroup
     607             :  * \brief Base class for equations of state which need three independent
     608             :  * thermodynamic variables in order to determine the pressure.
     609             :  *
     610             :  * The template parameter `IsRelativistic` is `true` for relativistic equations
     611             :  * of state and `false` for non-relativistic equations of state.
     612             :  */
     613             : template <bool IsRelativistic>
     614             : class EquationOfState<IsRelativistic, 3> : public PUP::able {
     615             :  public:
     616             :   static constexpr bool is_relativistic = IsRelativistic;
     617             :   static constexpr size_t thermodynamic_dim = 3;
     618             :   using creatable_classes =
     619             :       typename detail::DerivedClasses<IsRelativistic, 3>::type;
     620             : 
     621             :   EquationOfState() = default;
     622             :   EquationOfState(const EquationOfState&) = default;
     623             :   EquationOfState& operator=(const EquationOfState&) = default;
     624             :   EquationOfState(EquationOfState&&) = default;
     625             :   EquationOfState& operator=(EquationOfState&&) = default;
     626             :   ~EquationOfState() override = default;
     627             : 
     628             :   explicit EquationOfState(CkMigrateMessage* msg) : PUP::able(msg) {}
     629             : 
     630             :   WRAPPED_PUPable_abstract(EquationOfState);  // NOLINT
     631             : 
     632             :   virtual inline std::unique_ptr<EquationOfState<IsRelativistic, 3>> get_clone()
     633             :       const = 0;
     634             : 
     635             :   virtual bool is_equal(
     636             :       const EquationOfState<IsRelativistic, 3>& rhs) const = 0;
     637             : 
     638             :   virtual std::unique_ptr<EquationOfState<IsRelativistic, 3>>
     639             :   promote_to_3d_eos() const {
     640             :     return this->get_clone();
     641             :   }
     642             : 
     643             :   /// \brief Returns `true` if the EOS is barotropic
     644             :   virtual bool is_barotropic() const = 0;
     645             : 
     646             :   /// \brief Returns `true` if the EOS is in beta-equilibrium
     647             :   virtual bool is_equilibrium() const = 0;
     648             : 
     649             :   /// @{
     650             :   /*!
     651             :    * Computes the electron fraction in beta-equilibrium \f$Y_e^{\rm eq}\f$ from
     652             :    * the rest mass density \f$\rho\f$ and the temperature \f$T\f$.
     653             :    */
     654             :   virtual Scalar<double> equilibrium_electron_fraction_from_density_temperature(
     655             :       const Scalar<double>& /*rest_mass_density*/,
     656             :       const Scalar<double>& /*temperature*/) const = 0;
     657             : 
     658             :   virtual Scalar<DataVector>
     659             :   equilibrium_electron_fraction_from_density_temperature(
     660             :       const Scalar<DataVector>& /*rest_mass_density*/,
     661             :       const Scalar<DataVector>& /*temperature*/) const = 0;
     662             :   /// @}
     663             : 
     664             :   /// @{
     665             :   /*!
     666             :    * Computes the pressure \f$p\f$ from the rest mass density \f$\rho\f$, the
     667             :    * specific internal energy \f$\epsilon\f$ and electron fraction \f$Y_e\f$.
     668             :    */
     669             :   virtual Scalar<double> pressure_from_density_and_energy(
     670             :       const Scalar<double>& /*rest_mass_density*/,
     671             :       const Scalar<double>& /*specific_internal_energy*/,
     672             :       const Scalar<double>& /*electron_fraction*/) const = 0;
     673             :   virtual Scalar<DataVector> pressure_from_density_and_energy(
     674             :       const Scalar<DataVector>& /*rest_mass_density*/,
     675             :       const Scalar<DataVector>& /*specific_internal_energy*/,
     676             :       const Scalar<DataVector>& /*electron_fraction*/) const = 0;
     677             :   /// @}
     678             : 
     679             :   /// @{
     680             :   /*!
     681             :    * Computes the pressure \f$p\f$ from the rest mass density \f$\rho\f$, the
     682             :    * temperature \f$T\f$, and electron fraction \f$Y_e\f$.
     683             :    */
     684             :   virtual Scalar<double> pressure_from_density_and_temperature(
     685             :       const Scalar<double>& /*rest_mass_density*/,
     686             :       const Scalar<double>& /*temperature*/,
     687             :       const Scalar<double>& /*electron_fraction*/) const = 0;
     688             :   virtual Scalar<DataVector> pressure_from_density_and_temperature(
     689             :       const Scalar<DataVector>& /*rest_mass_density*/,
     690             :       const Scalar<DataVector>& /*temperature*/,
     691             :       const Scalar<DataVector>& /*electron_fraction*/) const = 0;
     692             :   /// @}
     693             : 
     694             :   /// @{
     695             :   /*!
     696             :    * Computes the specific entropy \f$s\f$ from the rest mass density \f$\rho\f$
     697             :    * and the specific internal energy \f$\epsilon\f$ and electron fraction
     698             :    * \f$Y_e\f$.
     699             :    */
     700             :   virtual Scalar<double> specific_entropy_from_density_and_energy(
     701             :       const Scalar<double>& /*rest_mass_density*/,
     702             :       const Scalar<double>& /*specific_internal_energy*/,
     703             :       const Scalar<double>& /*electron_fraction*/) const = 0;
     704             : 
     705             :   virtual Scalar<DataVector> specific_entropy_from_density_and_energy(
     706             :       const Scalar<DataVector>& /*rest_mass_density*/,
     707             :       const Scalar<DataVector>& /*specific_internal_energy*/,
     708             :       const Scalar<DataVector>& /*electron_fraction*/) const = 0;
     709             :   /// @}
     710             : 
     711             :   /// @{
     712             :   /*!
     713             :    * Computes the specific entropy \f$s\f$ from the rest mass density \f$\rho\f$
     714             :    * and the temperature \f$T\f$ and electron fraction \f$Y_e\f$.
     715             :    */
     716             :   virtual Scalar<double> specific_entropy_from_density_and_temperature(
     717             :       const Scalar<double>& /*rest_mass_density*/,
     718             :       const Scalar<double>& /*temperature*/,
     719             :       const Scalar<double>& /*electron_fraction*/) const = 0;
     720             : 
     721             :   virtual Scalar<DataVector> specific_entropy_from_density_and_temperature(
     722             :       const Scalar<DataVector>& /*rest_mass_density*/,
     723             :       const Scalar<DataVector>& /*temperature*/,
     724             :       const Scalar<DataVector>& /*electron_fraction*/) const = 0;
     725             :   /// @}
     726             : 
     727             :   /// @{
     728             :   /*!
     729             :    * Computes the temperature \f$T\f$ from the rest mass
     730             :    * density \f$\rho\f$, the specific internal energy \f$\epsilon\f$,
     731             :    * and electron fraction \f$Y_e\f$.
     732             :    */
     733             :   virtual Scalar<double> temperature_from_density_and_energy(
     734             :       const Scalar<double>& /*rest_mass_density*/,
     735             :       const Scalar<double>& /*specific_internal_energy*/,
     736             :       const Scalar<double>& /*electron_fraction*/) const = 0;
     737             :   virtual Scalar<DataVector> temperature_from_density_and_energy(
     738             :       const Scalar<DataVector>& /*rest_mass_density*/,
     739             :       const Scalar<DataVector>& /*specific_internal_energy*/,
     740             :       const Scalar<DataVector>& /*electron_fraction*/) const = 0;
     741             :   /// @}
     742             : 
     743             :   /// @{
     744             :   /*!
     745             :    * Computes the specific internal energy \f$\epsilon\f$ from the rest mass
     746             :    * density \f$\rho\f$, the temperature \f$T\f$, and electron fraction
     747             :    * \f$Y_e\f$.
     748             :    */
     749             :   virtual Scalar<double> specific_internal_energy_from_density_and_temperature(
     750             :       const Scalar<double>& /*rest_mass_density*/,
     751             :       const Scalar<double>& /*temperature*/,
     752             :       const Scalar<double>& /*electron_fraction*/
     753             :   ) const = 0;
     754             :   virtual Scalar<DataVector>
     755             :   specific_internal_energy_from_density_and_temperature(
     756             :       const Scalar<DataVector>& /*rest_mass_density*/,
     757             :       const Scalar<DataVector>& /*temperature*/,
     758             :       const Scalar<DataVector>& /*electron_fraction*/
     759             :   ) const = 0;
     760             :   /// @}
     761             : 
     762             :   /// @{
     763             :   /*!
     764             :    * Computes adiabatic sound speed squared
     765             :    * \f[
     766             :    * c_s^2  \equiv \frac{\partial p}{\partial e} |_{s, Y_e} =
     767             :    * \frac{\rho}{h}\frac{\partial p}{\partial \rho} |_{e, Y_e} +
     768             :    * \frac{\partial p}{\partial e}|_{\rho, Y_e}
     769             :    * \f].
     770             :    * With \f$p, e\f$ the pressure and energy density respectively,
     771             :    * \f$s\f$ the entropy density, \f$Y_e\f$ the electron fraction
     772             :    * \f$\rho\f$ the rest-mass density, and \f$h\f$ the enthalpy density.
     773             :    * Note that \f$e\f$ is the total energy density and not the internal energy,
     774             :    * therefore
     775             :    * \f[
     776             :    * \frac{\partial p}{\partial \rho} |_{e, Y_e} \neq \chi \equiv \frac{\partial
     777             :    * p}{\partial \rho} |_{\epsilon, Y_e}
     778             :    * \f]
     779             :    * as defined in the 2-d EoS above. By definition
     780             :    * \f$ e = (1+\epsilon) \rho \f$ so holding \f$e\f$ constant
     781             :    * \f[
     782             :    *  0 = \frac{d e}{d \rho} = \frac{\partial e}{\partial \rho} +
     783             :    *     \frac{\partial e}{\partial \epsilon} \frac{\partial \epsilon}{\rho}.
     784             :    * \f]
     785             :    * (where we have suppressed \f$ Y_e\f$ dependence)
     786             :    * So \f$ \partial \epsilon /  \partial \rho |_{e} = (1 + \epsilon)/\rho \f$
     787             :    * and we can expand
     788             :    * \f[
     789             :    * \frac{\partial p}{\partial \rho} |_{e, Y_e} = \frac{\partial e}{\partial
     790             :    * \rho}_{\epsilon, Y_e} + \frac{(1 + \epsilon)}{\rho} \frac{\partial
     791             :    * e}{\partial \epsilon}|_{\rho, Y_e}
     792             :    * \f]
     793             :    *  Finally, we can rewrite the entire sound speed using only the rest-mass
     794             :    * density, specific internal energy, and electron fraction as variables,
     795             :    * by using \f$ \frac{\partial e}{\partial \epsilon}|_{\rho, Y_e} = 1 \f$
     796             :    * \f[
     797             :    * c_s^2   =
     798             :    * \frac{\rho}{h}\frac{\partial p}{\partial \rho} |_{\epsilon, Y_e} +
     799             :    * \frac{1}{\rho} \frac{\partial p}{\partial \epsilon}|_{\rho, Y_e} \left(
     800             :    * 1 - \frac{(1 + \epsilon)\rho}{h}\right)
     801             :    * \f]
     802             :    * Which reduces to our preferred form
     803             :    * \f[
     804             :    * c_s^2 =
     805             :    * \frac{\rho}{h}(\chi + \kappa)
     806             :    * \f]
     807             :    *
     808             :    * Computed as a function of temperature, rest-mass density and electron
     809             :    * fraction. Note that this will break thermodynamic consistency if the
     810             :    * pressure and internal energy interpolated separately. The precise impact of
     811             :    * this will depend on the EoS and numerical scheme used for the evolution.
     812             :    */
     813             :   virtual Scalar<double> sound_speed_squared_from_density_and_temperature(
     814             :       const Scalar<double>& /*rest_mass_density*/,
     815             :       const Scalar<double>& /*temperature*/,
     816             :       const Scalar<double>& /*electron_fraction*/) const = 0;
     817             :   virtual Scalar<DataVector> sound_speed_squared_from_density_and_temperature(
     818             :       const Scalar<DataVector>& /*rest_mass_density*/,
     819             :       const Scalar<DataVector>& /*temperature*/,
     820             :       const Scalar<DataVector>& /*electron_fraction*/) const = 0;
     821             :   /// @}
     822             : 
     823             :   /// @{
     824             :   /*!
     825             :    * Computes
     826             :    * \f[
     827             :    * \kappa \equiv \left.\frac{\partial p}{\partial \epsilon}\right|_{\rho,Y_e},
     828             :    * \f].
     829             :    * Unlike the 1D and 2D interfaces which return
     830             :    * \f$\kappa p/\rho^2\f$, the 3D interface returns \f$\kappa\f$ directly:
     831             :    * 3D tables typically store \f$\kappa\f$ as an independent quantity, and
     832             :    * callers that need \f$\kappa p/\rho^2\f$ can form it themselves.
     833             :    */
     834             :   virtual Scalar<double> kappa_from_density_and_temperature(
     835             :       const Scalar<double>& /*rest_mass_density*/,
     836             :       const Scalar<double>& /*temperature*/,
     837             :       const Scalar<double>& /*electron_fraction*/) const = 0;
     838             : 
     839             :   virtual Scalar<DataVector> kappa_from_density_and_temperature(
     840             :       const Scalar<DataVector>& /*rest_mass_density*/,
     841             :       const Scalar<DataVector>& /*temperature*/,
     842             :       const Scalar<DataVector>& /*electron_fraction*/) const = 0;
     843             :   /// @}
     844             : 
     845             :   /// @{
     846             :   /*!
     847             :    * Computes
     848             :    * \f[
     849             :    * \zeta \equiv \left.\frac{\partial p}{\partial Y_e}\right|_{\rho,\epsilon},
     850             :    * \f].
     851             :    */
     852             :   virtual Scalar<double> zeta_from_density_and_temperature(
     853             :       const Scalar<double>& /*rest_mass_density*/,
     854             :       const Scalar<double>& /*temperature*/,
     855             :       const Scalar<double>& /*electron_fraction*/) const = 0;
     856             : 
     857             :   virtual Scalar<DataVector> zeta_from_density_and_temperature(
     858             :       const Scalar<DataVector>& /*rest_mass_density*/,
     859             :       const Scalar<DataVector>& /*temperature*/,
     860             :       const Scalar<DataVector>& /*electron_fraction*/) const = 0;
     861             :   /// @}
     862             : 
     863             :   /// The lower bound of the electron fraction that is valid for this EOS
     864             :   virtual double electron_fraction_lower_bound() const = 0;
     865             : 
     866             :   /// The upper bound of the electron fraction that is valid for this EOS
     867             :   virtual double electron_fraction_upper_bound() const = 0;
     868             : 
     869             :   /// The lower bound of the rest mass density that is valid for this EOS
     870             :   virtual double rest_mass_density_lower_bound() const = 0;
     871             : 
     872             :   /// The upper bound of the rest mass density that is valid for this EOS
     873             :   virtual double rest_mass_density_upper_bound() const = 0;
     874             : 
     875             :   /// The lower bound of the specific internal energy that is valid for this EOS
     876             :   /// at the given rest mass density \f$\rho\f$ and electron fraction \f$Y_e\f$.
     877             :   virtual double specific_internal_energy_lower_bound(
     878             :       const double rest_mass_density, const double electron_fraction) const = 0;
     879             : 
     880             :   /// The upper bound of the specific internal energy that is valid for this EOS
     881             :   /// at the given rest mass density \f$\rho\f$ and electron fraction \f$Y_e\f$.
     882             :   virtual double specific_internal_energy_upper_bound(
     883             :       const double rest_mass_density, const double electron_fraction) const = 0;
     884             : 
     885             :   /// The lower bound of the specific enthalpy that is valid for this EOS
     886             :   virtual double specific_enthalpy_lower_bound() const = 0;
     887             : 
     888             :   /// The lower bound of the temperature that is valid for this EOS
     889             :   virtual double temperature_lower_bound() const = 0;
     890             : 
     891             :   /// The upper bound of the temperature that is valid for this EOS
     892             :   virtual double temperature_upper_bound() const = 0;
     893             : 
     894             :   /// The vacuum mass of a baryon for this EOS
     895             :   virtual double baryon_mass() const {
     896             :     return hydro::units::geometric::default_baryon_mass;
     897             :   }
     898             : };
     899             : 
     900             : /// Compare two equations of state for equality
     901             : template <bool IsRelLhs, bool IsRelRhs, size_t ThermoDimLhs,
     902             :           size_t ThermoDimRhs>
     903           1 : bool operator==(const EquationOfState<IsRelLhs, ThermoDimLhs>& lhs,
     904             :                 const EquationOfState<IsRelRhs, ThermoDimRhs>& rhs) {
     905             :   if constexpr (IsRelLhs == IsRelRhs and ThermoDimLhs == ThermoDimRhs) {
     906             :     return typeid(lhs) == typeid(rhs) and lhs.is_equal(rhs);
     907             :   } else {
     908             :     return false;
     909             :   }
     910             : }
     911             : template <bool IsRelLhs, bool IsRelRhs, size_t ThermoDimLhs,
     912             :           size_t ThermoDimRhs>
     913           0 : bool operator!=(const EquationOfState<IsRelLhs, ThermoDimLhs>& lhs,
     914             :                 const EquationOfState<IsRelRhs, ThermoDimRhs>& rhs) {
     915             :   return not(lhs == rhs);
     916             : }
     917             : }  // namespace EquationsOfState
     918             : 
     919             : /// \cond
     920             : #define EQUATION_OF_STATE_FUNCTIONS_1D                      \
     921             :   (pressure_from_density, rest_mass_density_from_enthalpy,  \
     922             :    specific_internal_energy_from_density, chi_from_density, \
     923             :    kappa_times_p_over_rho_squared_from_density)
     924             : 
     925             : #define EQUATION_OF_STATE_FUNCTIONS_2D                                   \
     926             :   (pressure_from_density_and_energy, pressure_from_density_and_enthalpy, \
     927             :    specific_entropy_from_density_and_energy,                             \
     928             :    specific_entropy_from_density_and_temperature,                        \
     929             :    specific_internal_energy_from_density_and_pressure,                   \
     930             :    temperature_from_density_and_energy,                                  \
     931             :    specific_internal_energy_from_density_and_temperature,                \
     932             :    chi_from_density_and_energy,                                          \
     933             :    kappa_times_p_over_rho_squared_from_density_and_energy)
     934             : 
     935             : #define EQUATION_OF_STATE_FUNCTIONS_3D                                      \
     936             :   (pressure_from_density_and_energy, pressure_from_density_and_temperature, \
     937             :    specific_entropy_from_density_and_energy,                                \
     938             :    specific_entropy_from_density_and_temperature,                           \
     939             :    temperature_from_density_and_energy,                                     \
     940             :    specific_internal_energy_from_density_and_temperature,                   \
     941             :    sound_speed_squared_from_density_and_temperature,                        \
     942             :    kappa_from_density_and_temperature, zeta_from_density_and_temperature)
     943             : 
     944             : #define EQUATION_OF_STATE_ARGUMENTS_EXPAND(z, n, type) \
     945             :   BOOST_PP_COMMA_IF(n) const Scalar<type>&
     946             : 
     947             : #define EQUATION_OF_STATE_FORWARD_DECLARE_MEMBERS_HELPER(r, DIM,        \
     948             :                                                          FUNCTION_NAME) \
     949             :   Scalar<double> FUNCTION_NAME(BOOST_PP_REPEAT(                         \
     950             :       DIM, EQUATION_OF_STATE_ARGUMENTS_EXPAND, double)) const override; \
     951             :   Scalar<DataVector> FUNCTION_NAME(BOOST_PP_REPEAT(                     \
     952             :       DIM, EQUATION_OF_STATE_ARGUMENTS_EXPAND, DataVector)) const override;
     953             : 
     954             : /// \endcond
     955             : 
     956             : /*!
     957             :  * \ingroup EquationsOfStateGroup
     958             :  * \brief Macro used to generate forward declarations of member functions in
     959             :  * derived classes
     960             :  */
     961           1 : #define EQUATION_OF_STATE_FORWARD_DECLARE_MEMBERS(DERIVED, DIM)            \
     962             :   BOOST_PP_LIST_FOR_EACH(                                                  \
     963             :       EQUATION_OF_STATE_FORWARD_DECLARE_MEMBERS_HELPER, DIM,               \
     964             :       BOOST_PP_TUPLE_TO_LIST(BOOST_PP_TUPLE_ELEM(                          \
     965             :           BOOST_PP_SUB(DIM, 1),                                            \
     966             :           (EQUATION_OF_STATE_FUNCTIONS_1D, EQUATION_OF_STATE_FUNCTIONS_2D, \
     967             :            EQUATION_OF_STATE_FUNCTIONS_3D))))                              \
     968             :                                                                            \
     969             :   /* clang-tidy: do not use non-const references */                        \
     970             :   void pup(PUP::er& p) override; /* NOLINT */                              \
     971             :                                                                            \
     972             :   explicit DERIVED(CkMigrateMessage* msg);
     973             : 
     974             : /// \cond
     975             : #define EQUATION_OF_STATE_FORWARD_ARGUMENTS(z, n, unused) \
     976             :   BOOST_PP_COMMA_IF(n) arg##n
     977             : 
     978             : #define EQUATION_OF_STATE_ARGUMENTS_EXPAND_NAMED(z, n, type) \
     979             :   BOOST_PP_COMMA_IF(n) const Scalar<type>& arg##n
     980             : 
     981             : #define EQUATION_OF_STATE_MEMBER_DEFINITIONS_HELPER(                        \
     982             :     TEMPLATE, DERIVED, DATA_TYPE, DIM, FUNCTION_NAME)                       \
     983             :   TEMPLATE                                                                  \
     984             :   Scalar<DATA_TYPE> DERIVED::FUNCTION_NAME(BOOST_PP_REPEAT(                 \
     985             :       DIM, EQUATION_OF_STATE_ARGUMENTS_EXPAND_NAMED, DATA_TYPE)) const {    \
     986             :     return FUNCTION_NAME##_impl(                                            \
     987             :         BOOST_PP_REPEAT(DIM, EQUATION_OF_STATE_FORWARD_ARGUMENTS, UNUSED)); \
     988             :   }
     989             : 
     990             : #define EQUATION_OF_STATE_MEMBER_DEFINITIONS_HELPER_2(r, ARGS, FUNCTION_NAME) \
     991             :   EQUATION_OF_STATE_MEMBER_DEFINITIONS_HELPER(                                \
     992             :       BOOST_PP_TUPLE_ELEM(0, ARGS), BOOST_PP_TUPLE_ELEM(1, ARGS),             \
     993             :       BOOST_PP_TUPLE_ELEM(2, ARGS), BOOST_PP_TUPLE_ELEM(3, ARGS),             \
     994             :       FUNCTION_NAME)
     995             : /// \endcond
     996             : 
     997             : #define EQUATION_OF_STATE_MEMBER_DEFINITIONS(TEMPLATE, DERIVED, DATA_TYPE, \
     998           0 :                                              DIM)                          \
     999             :   BOOST_PP_LIST_FOR_EACH(                                                  \
    1000             :       EQUATION_OF_STATE_MEMBER_DEFINITIONS_HELPER_2,                       \
    1001             :       (TEMPLATE, DERIVED, DATA_TYPE, DIM),                                 \
    1002             :       BOOST_PP_TUPLE_TO_LIST(BOOST_PP_TUPLE_ELEM(                          \
    1003             :           BOOST_PP_SUB(DIM, 1),                                            \
    1004             :           (EQUATION_OF_STATE_FUNCTIONS_1D, EQUATION_OF_STATE_FUNCTIONS_2D, \
    1005             :            EQUATION_OF_STATE_FUNCTIONS_3D))))
    1006             : 
    1007             : /// \cond
    1008             : #define EQUATION_OF_STATE_FORWARD_DECLARE_MEMBER_IMPLS_HELPER(r, DIM,        \
    1009             :                                                               FUNCTION_NAME) \
    1010             :   template <class DataType>                                                  \
    1011             :   Scalar<DataType> FUNCTION_NAME##_impl(BOOST_PP_REPEAT(                     \
    1012             :       DIM, EQUATION_OF_STATE_ARGUMENTS_EXPAND, DataType)) const;
    1013             : /// \endcond
    1014             : 
    1015           0 : #define EQUATION_OF_STATE_FORWARD_DECLARE_MEMBER_IMPLS(DIM)                \
    1016             :   BOOST_PP_LIST_FOR_EACH(                                                  \
    1017             :       EQUATION_OF_STATE_FORWARD_DECLARE_MEMBER_IMPLS_HELPER, DIM,          \
    1018             :       BOOST_PP_TUPLE_TO_LIST(BOOST_PP_TUPLE_ELEM(                          \
    1019             :           BOOST_PP_SUB(DIM, 1),                                            \
    1020             :           (EQUATION_OF_STATE_FUNCTIONS_1D, EQUATION_OF_STATE_FUNCTIONS_2D, \
    1021             :            EQUATION_OF_STATE_FUNCTIONS_3D))))

Generated by: LCOV version 1.14