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 "Utilities/Serialization/CharmPupable.hpp"
16 : #include "Utilities/TMPL.hpp"
17 : #include "Utilities/TypeTraits.hpp"
18 :
19 : /// \cond
20 : namespace EquationsOfState {
21 : template <bool IsRelativistic>
22 : class DarkEnergyFluid;
23 : template <typename ColdEquationOfState>
24 : class HybridEos;
25 : template <bool IsRelativistic>
26 : class IdealFluid;
27 : template <bool IsRelativistic>
28 : class PolytropicFluid;
29 : template <bool IsRelativistic>
30 : class PiecewisePolytropicFluid;
31 : class Spectral;
32 : template <typename LowDensityEoS>
33 : class Enthalpy;
34 : template <bool IsRelativistic>
35 : class Tabulated3D;
36 : } // namespace EquationsOfState
37 : /// \endcond
38 :
39 : /// Contains all equations of state, including base class
40 : namespace EquationsOfState {
41 :
42 : namespace detail {
43 : template <bool IsRelativistic, size_t ThermodynamicDim>
44 : struct DerivedClasses {};
45 :
46 : template <>
47 : struct DerivedClasses<true, 1> {
48 : using type = tmpl::list<Enthalpy<Enthalpy<Enthalpy<Spectral>>>,
49 : Enthalpy<Enthalpy<Spectral>>, Enthalpy<Spectral>,
50 : Enthalpy<PolytropicFluid<true>>,
51 : PiecewisePolytropicFluid<true>, PolytropicFluid<true>,
52 : Spectral>;
53 : };
54 :
55 : template <>
56 : struct DerivedClasses<false, 1> {
57 : using type =
58 : tmpl::list<PiecewisePolytropicFluid<false>, PolytropicFluid<false>>;
59 : };
60 :
61 : template <>
62 : struct DerivedClasses<true, 2> {
63 : using type = tmpl::list<DarkEnergyFluid<true>, IdealFluid<true>,
64 : HybridEos<PolytropicFluid<true>>>;
65 : };
66 :
67 : template <>
68 : struct DerivedClasses<false, 2> {
69 : using type = tmpl::list<IdealFluid<false>, HybridEos<PolytropicFluid<false>>>;
70 : };
71 :
72 : template <>
73 : struct DerivedClasses<true, 3> {
74 : using type = tmpl::list<Tabulated3D<true>>;
75 : };
76 :
77 : template <>
78 : struct DerivedClasses<false, 3> {
79 : using type = tmpl::list<Tabulated3D<false>>;
80 : };
81 :
82 : } // namespace detail
83 :
84 : /*!
85 : * \ingroup EquationsOfStateGroup
86 : * \brief Base class for equations of state depending on whether or not the
87 : * system is relativistic, and the number of independent thermodynamic variables
88 : * (`ThermodynamicDim`) needed to determine the pressure.
89 : *
90 : * The template parameter `IsRelativistic` is `true` for relativistic equations
91 : * of state and `false` for non-relativistic equations of state.
92 : */
93 : template <bool IsRelativistic, size_t ThermodynamicDim>
94 1 : class EquationOfState;
95 :
96 : template <typename T>
97 0 : struct get_eos_base_impl {
98 0 : using type = EquationsOfState::EquationOfState<T::is_relativistic,
99 : T::thermodynamic_dim>;
100 : };
101 :
102 : template <bool IsRelativistic, size_t ThermodynamicDim>
103 0 : struct get_eos_base_impl<
104 : EquationsOfState::EquationOfState<IsRelativistic, ThermodynamicDim>> {
105 0 : using type =
106 : EquationsOfState::EquationOfState<IsRelativistic, ThermodynamicDim>;
107 : };
108 :
109 : template <typename T>
110 0 : using get_eos_base = typename get_eos_base_impl<T>::type;
111 :
112 : /*!
113 : * \ingroup EquationsOfStateGroup
114 : * \brief Base class for equations of state which need one thermodynamic
115 : * variable in order to determine the pressure.
116 : *
117 : * The template parameter `IsRelativistic` is `true` for relativistic equations
118 : * of state and `false` for non-relativistic equations of state.
119 : */
120 : template <bool IsRelativistic>
121 1 : class EquationOfState<IsRelativistic, 1> : public PUP::able {
122 : public:
123 0 : static constexpr bool is_relativistic = IsRelativistic;
124 0 : static constexpr size_t thermodynamic_dim = 1;
125 0 : using creatable_classes =
126 : typename detail::DerivedClasses<IsRelativistic, 1>::type;
127 :
128 0 : EquationOfState() = default;
129 0 : EquationOfState(const EquationOfState&) = default;
130 0 : EquationOfState& operator=(const EquationOfState&) = default;
131 0 : EquationOfState(EquationOfState&&) = default;
132 0 : EquationOfState& operator=(EquationOfState&&) = default;
133 0 : ~EquationOfState() override = default;
134 :
135 0 : explicit EquationOfState(CkMigrateMessage* msg) : PUP::able(msg) {}
136 :
137 0 : WRAPPED_PUPable_abstract(EquationOfState); // NOLINT
138 :
139 0 : virtual std::unique_ptr<EquationOfState<IsRelativistic, 1>> get_clone()
140 : const = 0;
141 :
142 0 : virtual bool is_equal(
143 : const EquationOfState<IsRelativistic, 1>& rhs) const = 0;
144 :
145 : /// @{
146 : /*!
147 : * Computes the electron fraction in beta-equilibrium \f$Y_e^{\rm eq}\f$ from
148 : * the rest mass density \f$\rho\f$.
149 : */
150 1 : virtual Scalar<double> equilibrium_electron_fraction_from_density_temperature(
151 : const Scalar<double>& rest_mass_density,
152 : const Scalar<double>& /*temperature*/) const {
153 : return make_with_value<Scalar<double>>(rest_mass_density, 0.1);
154 : }
155 :
156 : virtual Scalar<DataVector>
157 1 : equilibrium_electron_fraction_from_density_temperature(
158 : const Scalar<DataVector>& rest_mass_density,
159 : const Scalar<DataVector>& /*temperature*/) const {
160 : return make_with_value<Scalar<DataVector>>(rest_mass_density, 0.1);
161 : }
162 : /// @}
163 :
164 : /// @{
165 : /*!
166 : * Computes the pressure \f$p\f$ from the rest mass density \f$\rho\f$.
167 : */
168 1 : virtual Scalar<double> pressure_from_density(
169 : const Scalar<double>& /*rest_mass_density*/) const = 0;
170 1 : virtual Scalar<DataVector> pressure_from_density(
171 : const Scalar<DataVector>& /*rest_mass_density*/) const = 0;
172 : /// @}
173 :
174 : /// @{
175 : /*!
176 : * Computes the rest mass density \f$\rho\f$ from the specific enthalpy
177 : * \f$h\f$.
178 : */
179 1 : virtual Scalar<double> rest_mass_density_from_enthalpy(
180 : const Scalar<double>& /*specific_enthalpy*/) const = 0;
181 1 : virtual Scalar<DataVector> rest_mass_density_from_enthalpy(
182 : const Scalar<DataVector>& /*specific_enthalpy*/) const = 0;
183 : /// @}
184 :
185 : /// @{
186 : /*!
187 : * Computes the specific internal energy \f$\epsilon\f$ from the rest mass
188 : * density \f$\rho\f$.
189 : */
190 1 : virtual Scalar<double> specific_internal_energy_from_density(
191 : const Scalar<double>& /*rest_mass_density*/) const = 0;
192 1 : virtual Scalar<DataVector> specific_internal_energy_from_density(
193 : const Scalar<DataVector>& /*rest_mass_density*/) const = 0;
194 : /// @}
195 :
196 : /// @{
197 : /*!
198 : * Computes the temperature \f$T\f$ from the rest mass
199 : * density \f$\rho\f$.
200 : */
201 1 : virtual Scalar<double> temperature_from_density(
202 : const Scalar<double>& /*rest_mass_density*/) const {
203 : return Scalar<double>{0.0};
204 : }
205 1 : virtual Scalar<DataVector> temperature_from_density(
206 : const Scalar<DataVector>& rest_mass_density) const {
207 : return Scalar<DataVector>{DataVector{get(rest_mass_density).size(), 0.0}};
208 : }
209 : /// @}
210 :
211 : /// @{
212 : /*!
213 : * Computes the temperature \f$\T\f$ from the specific internal energy
214 : * \f$\epsilon\f$.
215 : */
216 1 : virtual Scalar<double> temperature_from_specific_internal_energy(
217 : const Scalar<double>& /*specific_internal_energy*/) const {
218 : return Scalar<double>{0.0};
219 : }
220 1 : virtual Scalar<DataVector> temperature_from_specific_internal_energy(
221 : const Scalar<DataVector>& specific_internal_energy) const {
222 : return Scalar<DataVector>{
223 : DataVector{get(specific_internal_energy).size(), 0.0}};
224 : }
225 : /// @}
226 :
227 : /// @{
228 : /*!
229 : * Computes \f$\chi=\partial p / \partial \rho\f$ from \f$\rho\f$, where
230 : * \f$p\f$ is the pressure and \f$\rho\f$ is the rest mass density.
231 : */
232 1 : virtual Scalar<double> chi_from_density(
233 : const Scalar<double>& /*rest_mass_density*/) const = 0;
234 1 : virtual Scalar<DataVector> chi_from_density(
235 : const Scalar<DataVector>& /*rest_mass_density*/) const = 0;
236 : /// @}
237 :
238 : /// @{
239 : /*!
240 : * Computes \f$\kappa p/\rho^2=(p/\rho^2)\partial p / \partial \epsilon\f$
241 : * from \f$\rho\f$, where \f$p\f$ is the pressure, \f$\rho\f$ is the rest mass
242 : * density, and \f$\epsilon\f$ is the specific internal energy.
243 : *
244 : * The reason for not returning just
245 : * \f$\kappa=\partial p / \partial \epsilon\f$ is to avoid division by zero
246 : * for small values of \f$\rho\f$ when assembling the speed of sound with
247 : * some equations of state.
248 : */
249 1 : virtual Scalar<double> kappa_times_p_over_rho_squared_from_density(
250 : const Scalar<double>& /*rest_mass_density*/) const = 0;
251 1 : virtual Scalar<DataVector> kappa_times_p_over_rho_squared_from_density(
252 : const Scalar<DataVector>& /*rest_mass_density*/) const = 0;
253 :
254 : /// The lower bound of the electron fraction that is valid for this EOS
255 1 : virtual double electron_fraction_lower_bound() const { return 0.0; }
256 :
257 : /// The upper bound of the electron fraction that is valid for this EOS
258 1 : virtual double electron_fraction_upper_bound() const { return 0.5; }
259 :
260 : /// The lower bound of the rest mass density that is valid for this EOS
261 1 : virtual double rest_mass_density_lower_bound() const = 0;
262 :
263 : /// The upper bound of the rest mass density that is valid for this EOS
264 1 : virtual double rest_mass_density_upper_bound() const = 0;
265 :
266 : /// The lower bound of the specific internal energy that is valid for this EOS
267 : /// at the given rest mass density \f$\rho\f$
268 1 : virtual double specific_internal_energy_lower_bound(
269 : const double rest_mass_density) const = 0;
270 :
271 : /// The upper bound of the specific internal energy that is valid for this EOS
272 : /// at the given rest mass density \f$\rho\f$
273 1 : virtual double specific_internal_energy_upper_bound(
274 : const double rest_mass_density) const = 0;
275 :
276 : /// The lower bound of the specific enthalpy that is valid for this EOS
277 1 : virtual double specific_enthalpy_lower_bound() const = 0;
278 : };
279 :
280 : /*!
281 : * \ingroup EquationsOfStateGroup
282 : * \brief Base class for equations of state which need two independent
283 : * thermodynamic variables in order to determine the pressure.
284 : *
285 : * The template parameter `IsRelativistic` is `true` for relativistic equations
286 : * of state and `false` for non-relativistic equations of state.
287 : */
288 : template <bool IsRelativistic>
289 1 : class EquationOfState<IsRelativistic, 2> : public PUP::able {
290 : public:
291 0 : static constexpr bool is_relativistic = IsRelativistic;
292 0 : static constexpr size_t thermodynamic_dim = 2;
293 0 : using creatable_classes =
294 : typename detail::DerivedClasses<IsRelativistic, 2>::type;
295 :
296 0 : EquationOfState() = default;
297 0 : EquationOfState(const EquationOfState&) = default;
298 0 : EquationOfState& operator=(const EquationOfState&) = default;
299 0 : EquationOfState(EquationOfState&&) = default;
300 0 : EquationOfState& operator=(EquationOfState&&) = default;
301 0 : ~EquationOfState() override = default;
302 :
303 0 : explicit EquationOfState(CkMigrateMessage* msg) : PUP::able(msg) {}
304 :
305 0 : WRAPPED_PUPable_abstract(EquationOfState); // NOLINT
306 :
307 0 : virtual inline std::unique_ptr<EquationOfState<IsRelativistic, 2>> get_clone()
308 : const = 0;
309 :
310 0 : virtual bool is_equal(
311 : const EquationOfState<IsRelativistic, 2>& rhs) const = 0;
312 :
313 : /// @{
314 : /*!
315 : * Computes the electron fraction in beta-equilibrium \f$Y_e^{\rm eq}\f$ from
316 : * the rest mass density \f$\rho\f$ and the temperature \f$T\f$.
317 : */
318 1 : virtual Scalar<double> equilibrium_electron_fraction_from_density_temperature(
319 : const Scalar<double>& rest_mass_density,
320 : const Scalar<double>& /*temperature*/) const {
321 : return make_with_value<Scalar<double>>(rest_mass_density, 0.1);
322 : }
323 :
324 : virtual Scalar<DataVector>
325 1 : equilibrium_electron_fraction_from_density_temperature(
326 : const Scalar<DataVector>& rest_mass_density,
327 : const Scalar<DataVector>& /*temperature*/) const {
328 : return make_with_value<Scalar<DataVector>>(rest_mass_density, 0.1);
329 : }
330 : /// @}
331 :
332 : /// @{
333 : /*!
334 : * Computes the pressure \f$p\f$ from the rest mass density \f$\rho\f$ and the
335 : * specific internal energy \f$\epsilon\f$.
336 : */
337 1 : virtual Scalar<double> pressure_from_density_and_energy(
338 : const Scalar<double>& /*rest_mass_density*/,
339 : const Scalar<double>& /*specific_internal_energy*/) const = 0;
340 1 : virtual Scalar<DataVector> pressure_from_density_and_energy(
341 : const Scalar<DataVector>& /*rest_mass_density*/,
342 : const Scalar<DataVector>& /*specific_internal_energy*/) const = 0;
343 : /// @}
344 :
345 : /// @{
346 : /*!
347 : * Computes the pressure \f$p\f$ from the rest mass density \f$\rho\f$ and the
348 : * specific enthalpy \f$h\f$.
349 : */
350 1 : virtual Scalar<double> pressure_from_density_and_enthalpy(
351 : const Scalar<double>& /*rest_mass_density*/,
352 : const Scalar<double>& /*specific_enthalpy*/) const = 0;
353 1 : virtual Scalar<DataVector> pressure_from_density_and_enthalpy(
354 : const Scalar<DataVector>& /*rest_mass_density*/,
355 : const Scalar<DataVector>& /*specific_enthalpy*/) const = 0;
356 : /// @}
357 :
358 : /// @{
359 : /*!
360 : * Computes the specific internal energy \f$\epsilon\f$ from the rest mass
361 : * density \f$\rho\f$ and the pressure \f$p\f$.
362 : */
363 1 : virtual Scalar<double> specific_internal_energy_from_density_and_pressure(
364 : const Scalar<double>& /*rest_mass_density*/,
365 : const Scalar<double>& /*pressure*/) const = 0;
366 1 : virtual Scalar<DataVector> specific_internal_energy_from_density_and_pressure(
367 : const Scalar<DataVector>& /*rest_mass_density*/,
368 : const Scalar<DataVector>& /*pressure*/) const = 0;
369 : /// @}
370 :
371 : /// @{
372 : /*!
373 : * Computes the temperature \f$T\f$ from the rest mass
374 : * density \f$\rho\f$ and the specific internal energy \f$\epsilon\f$.
375 : */
376 1 : virtual Scalar<double> temperature_from_density_and_energy(
377 : const Scalar<double>& /*rest_mass_density*/,
378 : const Scalar<double>& /*specific_internal_energy*/) const = 0;
379 1 : virtual Scalar<DataVector> temperature_from_density_and_energy(
380 : const Scalar<DataVector>& /*rest_mass_density*/,
381 : const Scalar<DataVector>& /*specific_internal_energy*/) const = 0;
382 : /// @}
383 :
384 : /// @{
385 : /*!
386 : * Computes the specific internal energy \f$\epsilon\f$ from the rest mass
387 : * density \f$\rho\f$ and the temperature \f$T\f$.
388 : */
389 1 : virtual Scalar<double> specific_internal_energy_from_density_and_temperature(
390 : const Scalar<double>& /*rest_mass_density*/,
391 : const Scalar<double>& /*temperature*/) const = 0;
392 : virtual Scalar<DataVector>
393 1 : specific_internal_energy_from_density_and_temperature(
394 : const Scalar<DataVector>& /*rest_mass_density*/,
395 : const Scalar<DataVector>& /*temperature*/) const = 0;
396 : /// @}
397 :
398 : /// @{
399 : /*!
400 : * Computes \f$\chi=\partial p / \partial \rho\f$ from the \f$\rho\f$ and
401 : * \f$\epsilon\f$, where \f$p\f$ is the pressure, \f$\rho\f$ is the rest mass
402 : * density, and \f$\epsilon\f$ is the specific internal energy.
403 : */
404 1 : virtual Scalar<double> chi_from_density_and_energy(
405 : const Scalar<double>& /*rest_mass_density*/,
406 : const Scalar<double>& /*specific_internal_energy*/) const = 0;
407 1 : virtual Scalar<DataVector> chi_from_density_and_energy(
408 : const Scalar<DataVector>& /*rest_mass_density*/,
409 : const Scalar<DataVector>& /*specific_internal_energy*/) const = 0;
410 : /// @}
411 :
412 : /// @{
413 : /*!
414 : * Computes \f$\kappa p/\rho^2=(p/\rho^2)\partial p / \partial \epsilon\f$
415 : * from \f$\rho\f$ and \f$\epsilon\f$, where \f$p\f$ is the pressure,
416 : * \f$\rho\f$ is the rest mass density, and \f$\epsilon\f$ is the specific
417 : * internal energy.
418 : *
419 : * The reason for not returning just
420 : * \f$\kappa=\partial p / \partial \epsilon\f$ is to avoid division by zero
421 : * for small values of \f$\rho\f$ when assembling the speed of sound with
422 : * some equations of state.
423 : */
424 1 : virtual Scalar<double> kappa_times_p_over_rho_squared_from_density_and_energy(
425 : const Scalar<double>& /*rest_mass_density*/,
426 : const Scalar<double>& /*specific_internal_energy*/) const = 0;
427 : virtual Scalar<DataVector>
428 1 : kappa_times_p_over_rho_squared_from_density_and_energy(
429 : const Scalar<DataVector>& /*rest_mass_density*/,
430 : const Scalar<DataVector>& /*specific_internal_energy*/) const = 0;
431 : /// @}
432 :
433 : /// The lower bound of the electron fraction that is valid for this EOS
434 1 : virtual double electron_fraction_lower_bound() const { return 0.0; }
435 :
436 : /// The upper bound of the electron fraction that is valid for this EOS
437 1 : virtual double electron_fraction_upper_bound() const { return 0.5; }
438 :
439 : /// The lower bound of the rest mass density that is valid for this EOS
440 1 : virtual double rest_mass_density_lower_bound() const = 0;
441 :
442 : /// The upper bound of the rest mass density that is valid for this EOS
443 1 : virtual double rest_mass_density_upper_bound() const = 0;
444 :
445 : /// The lower bound of the specific internal energy that is valid for this EOS
446 : /// at the given rest mass density \f$\rho\f$
447 1 : virtual double specific_internal_energy_lower_bound(
448 : const double rest_mass_density) const = 0;
449 :
450 : /// The upper bound of the specific internal energy that is valid for this EOS
451 : /// at the given rest mass density \f$\rho\f$
452 1 : virtual double specific_internal_energy_upper_bound(
453 : const double rest_mass_density) const = 0;
454 :
455 : /// The lower bound of the specific enthalpy that is valid for this EOS
456 1 : virtual double specific_enthalpy_lower_bound() const = 0;
457 : };
458 :
459 : /*!
460 : * \ingroup EquationsOfStateGroup
461 : * \brief Base class for equations of state which need three independent
462 : * thermodynamic variables in order to determine the pressure.
463 : *
464 : * The template parameter `IsRelativistic` is `true` for relativistic equations
465 : * of state and `false` for non-relativistic equations of state.
466 : */
467 : template <bool IsRelativistic>
468 1 : class EquationOfState<IsRelativistic, 3> : public PUP::able {
469 : public:
470 0 : static constexpr bool is_relativistic = IsRelativistic;
471 0 : static constexpr size_t thermodynamic_dim = 3;
472 0 : using creatable_classes =
473 : typename detail::DerivedClasses<IsRelativistic, 3>::type;
474 :
475 0 : EquationOfState() = default;
476 0 : EquationOfState(const EquationOfState&) = default;
477 0 : EquationOfState& operator=(const EquationOfState&) = default;
478 0 : EquationOfState(EquationOfState&&) = default;
479 0 : EquationOfState& operator=(EquationOfState&&) = default;
480 0 : ~EquationOfState() override = default;
481 :
482 0 : explicit EquationOfState(CkMigrateMessage* msg) : PUP::able(msg) {}
483 :
484 0 : WRAPPED_PUPable_abstract(EquationOfState); // NOLINT
485 :
486 0 : virtual inline std::unique_ptr<EquationOfState<IsRelativistic, 3>> get_clone()
487 : const = 0;
488 :
489 0 : virtual bool is_equal(
490 : const EquationOfState<IsRelativistic, 3>& rhs) const = 0;
491 : /// @{
492 : /*!
493 : * Computes the pressure \f$p\f$ from the rest mass density \f$\rho\f$, the
494 : * specific internal energy \f$\epsilon\f$ and electron fraction \f$Y_e\f$.
495 : */
496 1 : virtual Scalar<double> pressure_from_density_and_energy(
497 : const Scalar<double>& /*rest_mass_density*/,
498 : const Scalar<double>& /*specific_internal_energy*/,
499 : const Scalar<double>& /*electron_fraction*/) const = 0;
500 1 : virtual Scalar<DataVector> pressure_from_density_and_energy(
501 : const Scalar<DataVector>& /*rest_mass_density*/,
502 : const Scalar<DataVector>& /*specific_internal_energy*/,
503 : const Scalar<DataVector>& /*electron_fraction*/) const = 0;
504 : /// @}
505 :
506 : /// @{
507 : /*!
508 : * Computes the pressure \f$p\f$ from the rest mass density \f$\rho\f$, the
509 : * temperature \f$T\f$, and electron fraction \f$Y_e\f$.
510 : */
511 1 : virtual Scalar<double> pressure_from_density_and_temperature(
512 : const Scalar<double>& /*rest_mass_density*/,
513 : const Scalar<double>& /*temperature*/,
514 : const Scalar<double>& /*electron_fraction*/) const = 0;
515 1 : virtual Scalar<DataVector> pressure_from_density_and_temperature(
516 : const Scalar<DataVector>& /*rest_mass_density*/,
517 : const Scalar<DataVector>& /*temperature*/,
518 : const Scalar<DataVector>& /*electron_fraction*/) const = 0;
519 : /// @}
520 :
521 : /// @{
522 : /*!
523 : * Computes the temperature \f$T\f$ from the rest mass
524 : * density \f$\rho\f$, the specific internal energy \f$\epsilon\f$,
525 : * and electron fraction \f$Y_e\f$.
526 : */
527 1 : virtual Scalar<double> temperature_from_density_and_energy(
528 : const Scalar<double>& /*rest_mass_density*/,
529 : const Scalar<double>& /*specific_internal_energy*/,
530 : const Scalar<double>& /*electron_fraction*/) const = 0;
531 1 : virtual Scalar<DataVector> temperature_from_density_and_energy(
532 : const Scalar<DataVector>& /*rest_mass_density*/,
533 : const Scalar<DataVector>& /*specific_internal_energy*/,
534 : const Scalar<DataVector>& /*electron_fraction*/) const = 0;
535 : /// @}
536 :
537 : /// @{
538 : /*!
539 : * Computes the specific internal energy \f$\epsilon\f$ from the rest mass
540 : * density \f$\rho\f$, the temperature \f$T\f$, and electron fraction
541 : * \f$Y_e\f$.
542 : */
543 1 : virtual Scalar<double> specific_internal_energy_from_density_and_temperature(
544 : const Scalar<double>& /*rest_mass_density*/,
545 : const Scalar<double>& /*temperature*/,
546 : const Scalar<double>& /*electron_fraction*/
547 : ) const = 0;
548 : virtual Scalar<DataVector>
549 1 : specific_internal_energy_from_density_and_temperature(
550 : const Scalar<DataVector>& /*rest_mass_density*/,
551 : const Scalar<DataVector>& /*temperature*/,
552 : const Scalar<DataVector>& /*electron_fraction*/
553 : ) const = 0;
554 : /// @}
555 :
556 : /// @{
557 : /*!
558 : * Computes sound speed squared \f$ c_s^2 \f$.
559 : *
560 : * Interpolate as a function of temperature, rest-mass density and electron
561 : * fraction. Note that this will break thermodynamic consistency with the
562 : * pressure and internal energy interpolated separately. The precise impact of
563 : * this will depend on the EoS and numerical scheme used for the evolution.
564 : */
565 1 : virtual Scalar<double> sound_speed_squared_from_density_and_temperature(
566 : const Scalar<double>& /*rest_mass_density*/,
567 : const Scalar<double>& /*temperature*/,
568 : const Scalar<double>& /*electron_fraction*/) const = 0;
569 1 : virtual Scalar<DataVector> sound_speed_squared_from_density_and_temperature(
570 : const Scalar<DataVector>& /*rest_mass_density*/,
571 : const Scalar<DataVector>& /*temperature*/,
572 : const Scalar<DataVector>& /*electron_fraction*/) const = 0;
573 : /// @}
574 :
575 : /// The lower bound of the electron fraction that is valid for this EOS
576 1 : virtual double electron_fraction_lower_bound() const = 0;
577 :
578 : /// The upper bound of the electron fraction that is valid for this EOS
579 1 : virtual double electron_fraction_upper_bound() const = 0;
580 :
581 : /// The lower bound of the rest mass density that is valid for this EOS
582 1 : virtual double rest_mass_density_lower_bound() const = 0;
583 :
584 : /// The upper bound of the rest mass density that is valid for this EOS
585 1 : virtual double rest_mass_density_upper_bound() const = 0;
586 :
587 : /// The lower bound of the specific internal energy that is valid for this EOS
588 : /// at the given rest mass density \f$\rho\f$ and electron fraction \f$Y_e\f$.
589 1 : virtual double specific_internal_energy_lower_bound(
590 : const double rest_mass_density, const double electron_fraction) const = 0;
591 :
592 : /// The upper bound of the specific internal energy that is valid for this EOS
593 : /// at the given rest mass density \f$\rho\f$ and electron fraction \f$Y_e\f$.
594 1 : virtual double specific_internal_energy_upper_bound(
595 : const double rest_mass_density, const double electron_fraction) const = 0;
596 :
597 : /// The lower bound of the specific enthalpy that is valid for this EOS
598 1 : virtual double specific_enthalpy_lower_bound() const = 0;
599 :
600 : /// The lower bound of the temperature that is valid for this EOS
601 1 : virtual double temperature_lower_bound() const = 0;
602 :
603 : /// The upper bound of the temperature that is valid for this EOS
604 1 : virtual double temperature_upper_bound() const = 0;
605 : };
606 :
607 : /// Compare two equations of state for equality
608 : template <bool IsRelLhs, bool IsRelRhs, size_t ThermoDimLhs,
609 : size_t ThermoDimRhs>
610 1 : bool operator==(const EquationOfState<IsRelLhs, ThermoDimLhs>& lhs,
611 : const EquationOfState<IsRelRhs, ThermoDimRhs>& rhs) {
612 : if constexpr (IsRelLhs == IsRelRhs and ThermoDimLhs == ThermoDimRhs) {
613 : return typeid(lhs) == typeid(rhs) and lhs.is_equal(rhs);
614 : } else {
615 : return false;
616 : }
617 : }
618 : template <bool IsRelLhs, bool IsRelRhs, size_t ThermoDimLhs,
619 : size_t ThermoDimRhs>
620 0 : bool operator!=(const EquationOfState<IsRelLhs, ThermoDimLhs>& lhs,
621 : const EquationOfState<IsRelRhs, ThermoDimRhs>& rhs) {
622 : return not(lhs == rhs);
623 : }
624 : } // namespace EquationsOfState
625 :
626 : /// \cond
627 : #define EQUATION_OF_STATE_FUNCTIONS_1D \
628 : (pressure_from_density, rest_mass_density_from_enthalpy, \
629 : specific_internal_energy_from_density, chi_from_density, \
630 : kappa_times_p_over_rho_squared_from_density)
631 :
632 : #define EQUATION_OF_STATE_FUNCTIONS_2D \
633 : (pressure_from_density_and_energy, pressure_from_density_and_enthalpy, \
634 : specific_internal_energy_from_density_and_pressure, \
635 : temperature_from_density_and_energy, \
636 : specific_internal_energy_from_density_and_temperature, \
637 : chi_from_density_and_energy, \
638 : kappa_times_p_over_rho_squared_from_density_and_energy)
639 :
640 : #define EQUATION_OF_STATE_FUNCTIONS_3D \
641 : (pressure_from_density_and_energy, pressure_from_density_and_temperature, \
642 : temperature_from_density_and_energy, \
643 : specific_internal_energy_from_density_and_temperature, \
644 : sound_speed_squared_from_density_and_temperature)
645 :
646 : #define EQUATION_OF_STATE_ARGUMENTS_EXPAND(z, n, type) \
647 : BOOST_PP_COMMA_IF(n) const Scalar<type>&
648 :
649 : #define EQUATION_OF_STATE_FORWARD_DECLARE_MEMBERS_HELPER(r, DIM, \
650 : FUNCTION_NAME) \
651 : Scalar<double> FUNCTION_NAME(BOOST_PP_REPEAT( \
652 : DIM, EQUATION_OF_STATE_ARGUMENTS_EXPAND, double)) const override; \
653 : Scalar<DataVector> FUNCTION_NAME(BOOST_PP_REPEAT( \
654 : DIM, EQUATION_OF_STATE_ARGUMENTS_EXPAND, DataVector)) const override;
655 :
656 : /// \endcond
657 :
658 : /*!
659 : * \ingroup EquationsOfStateGroup
660 : * \brief Macro used to generate forward declarations of member functions in
661 : * derived classes
662 : */
663 1 : #define EQUATION_OF_STATE_FORWARD_DECLARE_MEMBERS(DERIVED, DIM) \
664 : BOOST_PP_LIST_FOR_EACH( \
665 : EQUATION_OF_STATE_FORWARD_DECLARE_MEMBERS_HELPER, DIM, \
666 : BOOST_PP_TUPLE_TO_LIST(BOOST_PP_TUPLE_ELEM( \
667 : BOOST_PP_SUB(DIM, 1), \
668 : (EQUATION_OF_STATE_FUNCTIONS_1D, EQUATION_OF_STATE_FUNCTIONS_2D, \
669 : EQUATION_OF_STATE_FUNCTIONS_3D)))) \
670 : \
671 : /* clang-tidy: do not use non-const references */ \
672 : void pup(PUP::er& p) override; /* NOLINT */ \
673 : \
674 : explicit DERIVED(CkMigrateMessage* msg);
675 :
676 : /// \cond
677 : #define EQUATION_OF_STATE_FORWARD_ARGUMENTS(z, n, unused) \
678 : BOOST_PP_COMMA_IF(n) arg##n
679 :
680 : #define EQUATION_OF_STATE_ARGUMENTS_EXPAND_NAMED(z, n, type) \
681 : BOOST_PP_COMMA_IF(n) const Scalar<type>& arg##n
682 :
683 : #define EQUATION_OF_STATE_MEMBER_DEFINITIONS_HELPER( \
684 : TEMPLATE, DERIVED, DATA_TYPE, DIM, FUNCTION_NAME) \
685 : TEMPLATE \
686 : Scalar<DATA_TYPE> DERIVED::FUNCTION_NAME(BOOST_PP_REPEAT( \
687 : DIM, EQUATION_OF_STATE_ARGUMENTS_EXPAND_NAMED, DATA_TYPE)) const { \
688 : return FUNCTION_NAME##_impl( \
689 : BOOST_PP_REPEAT(DIM, EQUATION_OF_STATE_FORWARD_ARGUMENTS, UNUSED)); \
690 : }
691 :
692 : #define EQUATION_OF_STATE_MEMBER_DEFINITIONS_HELPER_2(r, ARGS, FUNCTION_NAME) \
693 : EQUATION_OF_STATE_MEMBER_DEFINITIONS_HELPER( \
694 : BOOST_PP_TUPLE_ELEM(0, ARGS), BOOST_PP_TUPLE_ELEM(1, ARGS), \
695 : BOOST_PP_TUPLE_ELEM(2, ARGS), BOOST_PP_TUPLE_ELEM(3, ARGS), \
696 : FUNCTION_NAME)
697 : /// \endcond
698 :
699 : #define EQUATION_OF_STATE_MEMBER_DEFINITIONS(TEMPLATE, DERIVED, DATA_TYPE, \
700 0 : DIM) \
701 : BOOST_PP_LIST_FOR_EACH( \
702 : EQUATION_OF_STATE_MEMBER_DEFINITIONS_HELPER_2, \
703 : (TEMPLATE, DERIVED, DATA_TYPE, DIM), \
704 : BOOST_PP_TUPLE_TO_LIST(BOOST_PP_TUPLE_ELEM( \
705 : BOOST_PP_SUB(DIM, 1), \
706 : (EQUATION_OF_STATE_FUNCTIONS_1D, EQUATION_OF_STATE_FUNCTIONS_2D, \
707 : EQUATION_OF_STATE_FUNCTIONS_3D))))
708 :
709 : /// \cond
710 : #define EQUATION_OF_STATE_FORWARD_DECLARE_MEMBER_IMPLS_HELPER(r, DIM, \
711 : FUNCTION_NAME) \
712 : template <class DataType> \
713 : Scalar<DataType> FUNCTION_NAME##_impl(BOOST_PP_REPEAT( \
714 : DIM, EQUATION_OF_STATE_ARGUMENTS_EXPAND, DataType)) const;
715 : /// \endcond
716 :
717 0 : #define EQUATION_OF_STATE_FORWARD_DECLARE_MEMBER_IMPLS(DIM) \
718 : BOOST_PP_LIST_FOR_EACH( \
719 : EQUATION_OF_STATE_FORWARD_DECLARE_MEMBER_IMPLS_HELPER, DIM, \
720 : BOOST_PP_TUPLE_TO_LIST(BOOST_PP_TUPLE_ELEM( \
721 : BOOST_PP_SUB(DIM, 1), \
722 : (EQUATION_OF_STATE_FUNCTIONS_1D, EQUATION_OF_STATE_FUNCTIONS_2D, \
723 : EQUATION_OF_STATE_FUNCTIONS_3D))))
|