Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <array>
7 : #include <limits>
8 : #include <pup.h>
9 :
10 : #include "DataStructures/Tensor/TypeAliases.hpp"
11 : #include "Options/String.hpp"
12 : #include "PointwiseFunctions/AnalyticSolutions/AnalyticSolution.hpp"
13 : #include "PointwiseFunctions/Hydro/EquationsOfState/IdealFluid.hpp"
14 : #include "PointwiseFunctions/Hydro/TagsDeclarations.hpp"
15 : #include "Utilities/MakeArray.hpp"
16 : #include "Utilities/TMPL.hpp"
17 : #include "Utilities/TaggedTuple.hpp"
18 :
19 1 : namespace hydro::Solutions {
20 : /*!
21 : * \brief Smooth sinusoidal density wave.
22 : *
23 : * This is the generic infrastructure for a smooth flow solution that can be
24 : * used by the hydro systems to avoid code duplication. The solution has a
25 : * constant pressure and uniform spatial velocity provided that the rest mass
26 : * density satisfies the advection equation
27 : *
28 : * \f{align*}{
29 : * \partial_t\rho + v^i\partial_i\rho = 0,
30 : * \f}
31 : *
32 : * and the specific internal energy is a function of the rest mass density only,
33 : * \f$\epsilon = \epsilon(\rho)\f$. For testing purposes, this class implements
34 : * this solution for the case where \f$\rho\f$ is a sine wave. The user
35 : * specifies the mean flow velocity of the fluid, the wavevector of the density
36 : * profile, and the amplitude \f$A\f$ of the density profile. In Cartesian
37 : * coordinates \f$(x, y, z)\f$, and using dimensionless units, the primitive
38 : * variables at a given time \f$t\f$ are then
39 : *
40 : * \f{align*}{
41 : * \rho(\vec{x},t) &= 1 + A \sin(\vec{k}\cdot(\vec{x} - \vec{v}t)) \\
42 : * \vec{v}(\vec{x},t) &= [v_x, v_y, v_z]^{T},\\
43 : * P(\vec{x},t) &= P, \\
44 : * \epsilon(\vec{x}, t) &= \frac{P}{(\gamma - 1)\rho}\\
45 : * \f}
46 : *
47 : * where we have assumed \f$\epsilon\f$ and \f$\rho\f$ to be related through an
48 : * equation mathematically equivalent to the equation of state of an ideal gas,
49 : * where the pressure is held constant.
50 : */
51 : template <size_t Dim, bool IsRelativistic>
52 1 : class SmoothFlow : virtual public MarkAsAnalyticSolution {
53 : public:
54 0 : SmoothFlow() = default;
55 0 : SmoothFlow(const SmoothFlow& /*rhs*/) = default;
56 0 : SmoothFlow& operator=(const SmoothFlow& /*rhs*/) = default;
57 0 : SmoothFlow(SmoothFlow&& /*rhs*/) = default;
58 0 : SmoothFlow& operator=(SmoothFlow&& /*rhs*/) = default;
59 0 : ~SmoothFlow() = default;
60 :
61 0 : explicit SmoothFlow(CkMigrateMessage* /*unused*/);
62 :
63 : // clang-tidy: no runtime references
64 0 : void pup(PUP::er& /*p*/); // NOLINT
65 :
66 : protected:
67 0 : using equation_of_state_type = EquationsOfState::IdealFluid<IsRelativistic>;
68 :
69 : /// The mean flow velocity.
70 1 : struct MeanVelocity {
71 0 : using type = std::array<double, Dim>;
72 0 : static constexpr Options::String help = {"The mean flow velocity."};
73 : };
74 :
75 : /// The wave vector of the profile.
76 1 : struct WaveVector {
77 0 : using type = std::array<double, Dim>;
78 0 : static constexpr Options::String help = {"The wave vector of the profile."};
79 : };
80 :
81 : /// The constant pressure throughout the fluid.
82 1 : struct Pressure {
83 0 : using type = double;
84 0 : static constexpr Options::String help = {
85 : "The constant pressure throughout the fluid."};
86 0 : static type lower_bound() { return 0.0; }
87 : };
88 :
89 : /// The adiabatic index for the ideal fluid.
90 1 : struct AdiabaticIndex {
91 0 : using type = double;
92 0 : static constexpr Options::String help = {
93 : "The adiabatic index for the ideal fluid."};
94 0 : static type lower_bound() { return 1.0; }
95 : };
96 :
97 : /// The perturbation amplitude of the rest mass density of the fluid.
98 1 : struct PerturbationSize {
99 0 : using type = double;
100 0 : static constexpr Options::String help = {
101 : "The perturbation size of the rest mass density."};
102 0 : static type lower_bound() { return -1.0; }
103 0 : static type upper_bound() { return 1.0; }
104 : };
105 :
106 0 : using options = tmpl::list<MeanVelocity, WaveVector, Pressure, AdiabaticIndex,
107 : PerturbationSize>;
108 :
109 0 : SmoothFlow(const std::array<double, Dim>& mean_velocity,
110 : const std::array<double, Dim>& wavevector, double pressure,
111 : double adiabatic_index, double perturbation_size);
112 :
113 : /// @{
114 : /// Retrieve hydro variable at `(x, t)`
115 : template <typename DataType>
116 1 : auto variables(const tnsr::I<DataType, Dim>& x, double t,
117 : tmpl::list<hydro::Tags::RestMassDensity<DataType>> /*meta*/)
118 : const -> tuples::TaggedTuple<hydro::Tags::RestMassDensity<DataType>>;
119 :
120 : template <typename DataType>
121 1 : auto variables(
122 : const tnsr::I<DataType, Dim>& x, double t,
123 : tmpl::list<hydro::Tags::SpecificInternalEnergy<DataType>> /*meta*/) const
124 : -> tuples::TaggedTuple<hydro::Tags::SpecificInternalEnergy<DataType>>;
125 :
126 : template <typename DataType>
127 1 : auto variables(const tnsr::I<DataType, Dim>& x, double /*t*/,
128 : tmpl::list<hydro::Tags::Pressure<DataType>> /*meta*/) const
129 : -> tuples::TaggedTuple<hydro::Tags::Pressure<DataType>>;
130 :
131 : template <typename DataType>
132 1 : auto variables(
133 : const tnsr::I<DataType, Dim>& x, double /*t*/,
134 : tmpl::list<hydro::Tags::SpatialVelocity<DataType, Dim>> /*meta*/) const
135 : -> tuples::TaggedTuple<hydro::Tags::SpatialVelocity<DataType, Dim>>;
136 :
137 : template <typename DataType, bool LocalIsRelativistic = IsRelativistic,
138 : Requires<IsRelativistic and IsRelativistic == LocalIsRelativistic> =
139 : nullptr>
140 1 : auto variables(const tnsr::I<DataType, Dim>& x, double /*t*/,
141 : tmpl::list<hydro::Tags::LorentzFactor<DataType>> /*meta*/)
142 : const -> tuples::TaggedTuple<hydro::Tags::LorentzFactor<DataType>>;
143 :
144 : template <typename DataType>
145 1 : auto variables(const tnsr::I<DataType, Dim>& x, double t,
146 : tmpl::list<hydro::Tags::SpecificEnthalpy<DataType>> /*meta*/)
147 : const -> tuples::TaggedTuple<hydro::Tags::SpecificEnthalpy<DataType>>;
148 : /// @}
149 :
150 0 : const EquationsOfState::IdealFluid<IsRelativistic>& equation_of_state()
151 : const {
152 : return equation_of_state_;
153 : }
154 :
155 : private:
156 : template <size_t LocalDim, bool LocalIsRelativistic>
157 : friend bool
158 0 : operator==( // NOLINT (clang-tidy: readability-redundant-declaration)
159 : const SmoothFlow<LocalDim, LocalIsRelativistic>& lhs,
160 : const SmoothFlow<LocalDim, LocalIsRelativistic>& rhs);
161 :
162 : // Computes the phase.
163 : template <typename DataType>
164 0 : DataType k_dot_x_minus_vt(const tnsr::I<DataType, Dim>& x, double t) const;
165 :
166 0 : std::array<double, Dim> mean_velocity_ =
167 : make_array<Dim>(std::numeric_limits<double>::signaling_NaN());
168 0 : std::array<double, Dim> wavevector_ =
169 : make_array<Dim>(std::numeric_limits<double>::signaling_NaN());
170 0 : double pressure_ = std::numeric_limits<double>::signaling_NaN();
171 0 : double adiabatic_index_ = std::numeric_limits<double>::signaling_NaN();
172 0 : double perturbation_size_ = std::numeric_limits<double>::signaling_NaN();
173 : // The angular frequency.
174 0 : double k_dot_v_ = std::numeric_limits<double>::signaling_NaN();
175 0 : EquationsOfState::IdealFluid<IsRelativistic> equation_of_state_{};
176 : };
177 :
178 : template <size_t Dim, bool IsRelativistic>
179 0 : bool operator!=(const SmoothFlow<Dim, IsRelativistic>& lhs,
180 : const SmoothFlow<Dim, IsRelativistic>& rhs);
181 : } // namespace hydro::Solutions
|