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/dec.hpp>
7 : #include <boost/preprocessor/arithmetic/inc.hpp>
8 : #include <boost/preprocessor/control/expr_iif.hpp>
9 : #include <boost/preprocessor/list/adt.hpp>
10 : #include <boost/preprocessor/repetition/for.hpp>
11 : #include <boost/preprocessor/repetition/repeat.hpp>
12 : #include <boost/preprocessor/tuple/to_list.hpp>
13 : #include <limits>
14 : #include <pup.h>
15 :
16 : #include "DataStructures/Tensor/TypeAliases.hpp"
17 : #include "Options/String.hpp"
18 : #include "PointwiseFunctions/Hydro/EquationsOfState/EquationOfState.hpp"
19 : #include "PointwiseFunctions/Hydro/Units.hpp"
20 : #include "Utilities/Serialization/CharmPupable.hpp"
21 : #include "Utilities/TMPL.hpp"
22 :
23 : /// \cond
24 : class DataVector;
25 : /// \endcond
26 :
27 : namespace EquationsOfState {
28 : /*!
29 : * \ingroup EquationsOfStateGroup
30 : * \brief Equation of state for an ideal fluid
31 : *
32 : * An ideal fluid equation of state:
33 : *
34 : * \f[
35 : * p = \rho \epsilon(\gamma-1)
36 : * \f]
37 : *
38 : * where \f$\rho\f$ is the rest mass density, \f$\epsilon\f$ is the specific
39 : * internal energy, and \f$\gamma\f$ is the adiabatic index.
40 : *
41 : * The temperature \f$T\f$ is defined as
42 : *
43 : * \f[
44 : * T = (\gamma - 1) \epsilon
45 : * \f]
46 : */
47 : template <bool IsRelativistic>
48 1 : class IdealFluid : public EquationOfState<IsRelativistic, 2> {
49 : public:
50 0 : static constexpr size_t thermodynamic_dim = 2;
51 0 : static constexpr bool is_relativistic = IsRelativistic;
52 :
53 0 : struct AdiabaticIndex {
54 0 : using type = double;
55 0 : static constexpr Options::String help = {"Adiabatic index gamma"};
56 : };
57 :
58 0 : struct MinTemperature {
59 0 : using type = double;
60 0 : static type lower_bound() { return 0.0; }
61 0 : static constexpr Options::String help = {
62 : "Minimum temperature. "
63 : "This value must be non-negative."};
64 : };
65 :
66 0 : static constexpr Options::String help = {
67 : "An ideal fluid equation of state.\n"
68 : "The pressure is related to the rest mass density by p = rho * epsilon * "
69 : "(gamma - 1), where p is the pressure, rho is the rest mass density, "
70 : "epsilon is the specific internal energy, and gamma is the adiabatic "
71 : "index.\n"
72 : "The temperature T is defined as T=(gamma-1) * epsilon."};
73 :
74 0 : using options = tmpl::list<AdiabaticIndex, MinTemperature>;
75 :
76 0 : IdealFluid() = default;
77 0 : IdealFluid(const IdealFluid&) = default;
78 0 : IdealFluid& operator=(const IdealFluid&) = default;
79 0 : IdealFluid(IdealFluid&&) = default;
80 0 : IdealFluid& operator=(IdealFluid&&) = default;
81 0 : ~IdealFluid() override = default;
82 :
83 0 : explicit IdealFluid(double adiabatic_index, double min_temperature = 0.0);
84 :
85 : EQUATION_OF_STATE_FORWARD_DECLARE_MEMBERS(IdealFluid, 2)
86 :
87 0 : std::unique_ptr<EquationOfState<IsRelativistic, 2>> get_clone()
88 : const override;
89 :
90 0 : std::unique_ptr<EquationOfState<IsRelativistic, 3>> promote_to_3d_eos()
91 : const override;
92 :
93 : /// \brief Returns `true` if the EOS is barotropic
94 1 : bool is_barotropic() const override { return false; }
95 :
96 0 : bool operator==(const IdealFluid<IsRelativistic>& rhs) const;
97 :
98 0 : bool operator!=(const IdealFluid<IsRelativistic>& rhs) const;
99 :
100 0 : bool is_equal(const EquationOfState<IsRelativistic, 2>& rhs) const override;
101 :
102 0 : WRAPPED_PUPable_decl_base_template( // NOLINT
103 : SINGLE_ARG(EquationOfState<IsRelativistic, 2>), IdealFluid);
104 :
105 : /// The lower bound of the rest mass density that is valid for this EOS
106 1 : double rest_mass_density_lower_bound() const override { return 0.0; }
107 :
108 : /// The upper bound of the rest mass density that is valid for this EOS
109 1 : double rest_mass_density_upper_bound() const override {
110 : return std::numeric_limits<double>::max();
111 : }
112 :
113 : /// The lower bound of the specific internal energy that is valid for this EOS
114 : /// at the given rest mass density \f$\rho\f$
115 : /// If non-zero lower bound for temperature is provided, then the lower bound
116 : /// for specific internal energy is also non-zero accordingly.
117 1 : double specific_internal_energy_lower_bound(
118 : const double /* rest_mass_density */) const override {
119 : return (min_temperature_) / (adiabatic_index_ - 1.0);
120 : }
121 :
122 : /// The upper bound of the specific internal energy that is valid for this EOS
123 : /// at the given rest mass density \f$\rho\f$
124 1 : double specific_internal_energy_upper_bound(
125 : double rest_mass_density) const override;
126 :
127 : /// The lower bound of the specific enthalpy that is valid for this EOS
128 : /// If non-zero lower bound for temperature is provided, then the lower bound
129 : /// for specific internal enthalpy is also non-zero accordingly.
130 1 : double specific_enthalpy_lower_bound() const override {
131 : return IsRelativistic ? 1.0 + (adiabatic_index_ * min_temperature_) /
132 : (adiabatic_index_ - 1.0)
133 : : (adiabatic_index_ * min_temperature_) /
134 : (adiabatic_index_ - 1.0);
135 : }
136 :
137 : /// The lower bound of the temperature that is valid for this EOS.
138 : /// Non-zero lower bound could be set to impose floor on the specific
139 : /// internal energy.
140 1 : double temperature_lower_bound() const override { return min_temperature_; }
141 :
142 : /// The vacuum baryon mass for this EoS
143 1 : double baryon_mass() const override {
144 : return hydro::units::geometric::default_baryon_mass;
145 : }
146 :
147 : private:
148 : EQUATION_OF_STATE_FORWARD_DECLARE_MEMBER_IMPLS(2)
149 :
150 0 : double specific_entropy_from_density_and_energy(double density,
151 : double energy) const;
152 :
153 0 : double adiabatic_index_ = std::numeric_limits<double>::signaling_NaN();
154 0 : double min_temperature_ = std::numeric_limits<double>::signaling_NaN();
155 : };
156 :
157 : /// \cond
158 : template <bool IsRelativistic>
159 : PUP::able::PUP_ID EquationsOfState::IdealFluid<IsRelativistic>::my_PUP_ID = 0;
160 : /// \endcond
161 : } // namespace EquationsOfState
|
