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 <memory>
9 :
10 : #include "DataStructures/TaggedTuple.hpp"
11 : #include "DataStructures/Tensor/TypeAliases.hpp"
12 : #include "Options/String.hpp"
13 : #include "PointwiseFunctions/AnalyticData/AnalyticData.hpp"
14 : #include "PointwiseFunctions/AnalyticData/GrMhd/AnalyticData.hpp"
15 : #include "PointwiseFunctions/AnalyticSolutions/GeneralRelativity/Minkowski.hpp"
16 : #include "PointwiseFunctions/Hydro/EquationsOfState/EquationOfState.hpp"
17 : #include "PointwiseFunctions/Hydro/EquationsOfState/IdealFluid.hpp"
18 : #include "PointwiseFunctions/Hydro/TagsDeclarations.hpp"
19 : #include "PointwiseFunctions/InitialDataUtilities/InitialData.hpp"
20 : #include "Utilities/Serialization/CharmPupable.hpp"
21 : #include "Utilities/TMPL.hpp"
22 :
23 : /// \cond
24 : namespace PUP {
25 : class er;
26 : } // namespace PUP
27 : /// \endcond
28 :
29 : namespace grmhd::AnalyticData {
30 :
31 : /*!
32 : * \brief Analytic initial data for a slab jet
33 : *
34 : * This test problem is described in \cite Komissarov1999, Section 7.4 and
35 : * Fig. 13. It involves a high Lorentz factor jet injected into an ambient fluid
36 : * with an initially uniform magnetic field.
37 : *
38 : * The parameters used in \cite Komissarov1999 are:
39 : *
40 : * ```yaml
41 : * AdiabaticIndex: 4. / 3.
42 : * AmbientDensity: 10.
43 : * AmbientPressure: 0.01
44 : * JetDensity: 0.1
45 : * JetPressure: 0.01
46 : * # u^i = [20, 0, 0] or W = sqrt(401)
47 : * JetVelocity: [0.9987523388778445, 0., 0.]
48 : * InletRadius: 1.
49 : * MagneticField: [1., 0., 0.]
50 : * ```
51 : *
52 : * In \cite Komissarov1999 an artificial dissipation of $\eta_u=0.2$ and
53 : * $\eta_b=0.15$ is used, which we don't use here (yet).
54 : *
55 : * \note The inlet is currently modeled as the region of the domain where
56 : * $x <= 0$ and $|y| <= \mathrm{inlet_radius}$. Since this class only sets
57 : * the initial data, no fluid will flow into the domain during the evolution.
58 : * This has to be modeled as a boundary condition and is not implemented yet.
59 : */
60 1 : class SlabJet : public evolution::initial_data::InitialData,
61 : public MarkAsAnalyticData,
62 : public AnalyticDataBase,
63 : public hydro::TemperatureInitialization<SlabJet> {
64 : public:
65 0 : using equation_of_state_type = EquationsOfState::IdealFluid<true>;
66 :
67 0 : struct AdiabaticIndex {
68 0 : using type = double;
69 0 : static constexpr Options::String help = {
70 : "The adiabatic index of the ideal fluid"};
71 0 : static double lower_bound() { return 1.; }
72 : };
73 0 : struct AmbientDensity {
74 0 : using type = double;
75 0 : static constexpr Options::String help = {
76 : "Fluid rest mass density outside the jet"};
77 0 : static double lower_bound() { return 0.; }
78 0 : static double suggested_value() { return 10.; }
79 : };
80 0 : struct AmbientPressure {
81 0 : using type = double;
82 0 : static constexpr Options::String help = {"Fluid pressure outside the jet"};
83 0 : static double lower_bound() { return 0.; }
84 0 : static double suggested_value() { return 0.01; }
85 : };
86 0 : struct AmbientElectronFraction {
87 0 : using type = double;
88 0 : static constexpr Options::String help = {
89 : "Electron fraction outside the jet"};
90 0 : static double lower_bound() { return 0.; }
91 0 : static double upper_bound() { return 1.; }
92 : };
93 0 : struct JetDensity {
94 0 : using type = double;
95 0 : static constexpr Options::String help = {
96 : "Fluid rest mass density of the jet inlet"};
97 0 : static double lower_bound() { return 0.; }
98 0 : static double suggested_value() { return 0.1; }
99 : };
100 0 : struct JetPressure {
101 0 : using type = double;
102 0 : static constexpr Options::String help = {"Fluid pressure of the jet inlet"};
103 0 : static double lower_bound() { return 0.; }
104 0 : static double suggested_value() { return 0.01; }
105 : };
106 0 : struct JetElectronFraction {
107 0 : using type = double;
108 0 : static constexpr Options::String help = {
109 : "Electron fraction of the jet inlet"};
110 0 : static double lower_bound() { return 0.; }
111 0 : static double upper_bound() { return 1.; }
112 : };
113 0 : struct JetVelocity {
114 0 : using type = std::array<double, 3>;
115 0 : static constexpr Options::String help = {
116 : "Fluid spatial velocity of the jet inlet"};
117 : };
118 0 : struct InletRadius {
119 0 : using type = double;
120 0 : static constexpr Options::String help = {
121 : "Radius of the jet inlet around y=0"};
122 0 : static double lower_bound() { return 0.; }
123 0 : static double suggested_value() { return 1.; }
124 : };
125 0 : struct MagneticField {
126 0 : using type = std::array<double, 3>;
127 0 : static constexpr Options::String help = {
128 : "Initially uniform magnetic field"};
129 0 : static std::array<double, 3> suggested_value() { return {{1., 0., 0.}}; }
130 : };
131 :
132 0 : using options =
133 : tmpl::list<AdiabaticIndex, AmbientDensity, AmbientPressure,
134 : AmbientElectronFraction, JetDensity, JetPressure,
135 : JetElectronFraction, JetVelocity, InletRadius, MagneticField>;
136 :
137 0 : static constexpr Options::String help = {
138 : "Analytic initial data for a jet test."};
139 :
140 0 : SlabJet() = default;
141 0 : SlabJet(const SlabJet& /*rhs*/) = default;
142 0 : SlabJet& operator=(const SlabJet& /*rhs*/) = default;
143 0 : SlabJet(SlabJet&& /*rhs*/) = default;
144 0 : SlabJet& operator=(SlabJet&& /*rhs*/) = default;
145 0 : ~SlabJet() = default;
146 :
147 0 : SlabJet(double adiabatic_index, double ambient_density,
148 : double ambient_pressure, double ambient_electron_fraction,
149 : double jet_density, double jet_pressure, double jet_electron_fraction,
150 : std::array<double, 3> jet_velocity, double inlet_radius,
151 : std::array<double, 3> magnetic_field);
152 :
153 0 : auto get_clone() const
154 : -> std::unique_ptr<evolution::initial_data::InitialData> override;
155 :
156 : /// \cond
157 : explicit SlabJet(CkMigrateMessage* msg);
158 : using PUP::able::register_constructor;
159 : WRAPPED_PUPable_decl_template(SlabJet);
160 : /// \endcond
161 :
162 : /// @{
163 : /// Retrieve the GRMHD variables at a given position.
164 : template <typename DataType>
165 1 : auto variables(const tnsr::I<DataType, 3>& x,
166 : tmpl::list<hydro::Tags::RestMassDensity<DataType>> /*meta*/)
167 : const -> tuples::TaggedTuple<hydro::Tags::RestMassDensity<DataType>>;
168 :
169 : template <typename DataType>
170 1 : auto variables(const tnsr::I<DataType, 3>& x,
171 : tmpl::list<hydro::Tags::ElectronFraction<DataType>> /*meta*/)
172 : const -> tuples::TaggedTuple<hydro::Tags::ElectronFraction<DataType>>;
173 :
174 : template <typename DataType>
175 1 : auto variables(
176 : const tnsr::I<DataType, 3>& x,
177 : tmpl::list<hydro::Tags::SpecificInternalEnergy<DataType>> /*meta*/) const
178 : -> tuples::TaggedTuple<hydro::Tags::SpecificInternalEnergy<DataType>>;
179 :
180 : template <typename DataType>
181 1 : auto variables(const tnsr::I<DataType, 3>& x,
182 : tmpl::list<hydro::Tags::Pressure<DataType>> /*meta*/) const
183 : -> tuples::TaggedTuple<hydro::Tags::Pressure<DataType>>;
184 :
185 : template <typename DataType>
186 1 : auto variables(const tnsr::I<DataType, 3>& x,
187 : tmpl::list<hydro::Tags::SpatialVelocity<DataType, 3>> /*meta*/)
188 : const -> tuples::TaggedTuple<hydro::Tags::SpatialVelocity<DataType, 3>>;
189 :
190 : template <typename DataType>
191 1 : auto variables(const tnsr::I<DataType, 3>& x,
192 : tmpl::list<hydro::Tags::MagneticField<DataType, 3>> /*meta*/)
193 : const -> tuples::TaggedTuple<hydro::Tags::MagneticField<DataType, 3>>;
194 :
195 : template <typename DataType>
196 1 : auto variables(
197 : const tnsr::I<DataType, 3>& x,
198 : tmpl::list<hydro::Tags::DivergenceCleaningField<DataType>> /*meta*/) const
199 : -> tuples::TaggedTuple<hydro::Tags::DivergenceCleaningField<DataType>>;
200 :
201 : template <typename DataType>
202 1 : auto variables(const tnsr::I<DataType, 3>& x,
203 : tmpl::list<hydro::Tags::LorentzFactor<DataType>> /*meta*/)
204 : const -> tuples::TaggedTuple<hydro::Tags::LorentzFactor<DataType>>;
205 :
206 : template <typename DataType>
207 1 : auto variables(const tnsr::I<DataType, 3>& x,
208 : tmpl::list<hydro::Tags::SpecificEnthalpy<DataType>> /*meta*/)
209 : const -> tuples::TaggedTuple<hydro::Tags::SpecificEnthalpy<DataType>>;
210 :
211 : template <typename DataType>
212 1 : auto variables(const tnsr::I<DataType, 3>& x,
213 : tmpl::list<hydro::Tags::Temperature<DataType>> /*meta*/) const
214 : -> tuples::TaggedTuple<hydro::Tags::Temperature<DataType>> {
215 : return TemperatureInitialization::variables(
216 : x, tmpl::list<hydro::Tags::Temperature<DataType>>{});
217 : }
218 : /// @}
219 :
220 : /// Retrieve a collection of hydrodynamic variables at position x
221 : template <typename DataType, typename Tag1, typename Tag2, typename... Tags>
222 1 : tuples::TaggedTuple<Tag1, Tag2, Tags...> variables(
223 : const tnsr::I<DataType, 3>& x,
224 : tmpl::list<Tag1, Tag2, Tags...> /*meta*/) const {
225 : return {tuples::get<Tag1>(variables(x, tmpl::list<Tag1>{})),
226 : tuples::get<Tag2>(variables(x, tmpl::list<Tag2>{})),
227 : tuples::get<Tags>(variables(x, tmpl::list<Tags>{}))...};
228 : }
229 :
230 : /// Retrieve the metric variables
231 : template <typename DataType, typename Tag,
232 : Requires<tmpl::list_contains_v<
233 : gr::analytic_solution_tags<3, DataType>, Tag>> = nullptr>
234 1 : tuples::TaggedTuple<Tag> variables(const tnsr::I<DataType, 3>& x,
235 : tmpl::list<Tag> /*meta*/) const {
236 : constexpr double dummy_time = 0.0;
237 : return background_spacetime_.variables(x, dummy_time, tmpl::list<Tag>{});
238 : }
239 :
240 0 : const EquationsOfState::IdealFluid<true>& equation_of_state() const {
241 : return equation_of_state_;
242 : }
243 :
244 : // NOLINTNEXTLINE(google-runtime-references)
245 0 : void pup(PUP::er& /*p*/) override;
246 :
247 : private:
248 0 : EquationsOfState::IdealFluid<true> equation_of_state_{};
249 0 : gr::Solutions::Minkowski<3> background_spacetime_{};
250 :
251 0 : double ambient_density_ = std::numeric_limits<double>::signaling_NaN();
252 0 : double ambient_pressure_ = std::numeric_limits<double>::signaling_NaN();
253 0 : double ambient_electron_fraction_ =
254 : std::numeric_limits<double>::signaling_NaN();
255 0 : double jet_density_ = std::numeric_limits<double>::signaling_NaN();
256 0 : double jet_pressure_ = std::numeric_limits<double>::signaling_NaN();
257 0 : double jet_electron_fraction_ = std::numeric_limits<double>::signaling_NaN();
258 0 : std::array<double, 3> jet_velocity_{
259 : {std::numeric_limits<double>::signaling_NaN(),
260 : std::numeric_limits<double>::signaling_NaN(),
261 : std::numeric_limits<double>::signaling_NaN()}};
262 0 : double inlet_radius_ = std::numeric_limits<double>::signaling_NaN();
263 0 : std::array<double, 3> magnetic_field_{
264 : {std::numeric_limits<double>::signaling_NaN(),
265 : std::numeric_limits<double>::signaling_NaN(),
266 : std::numeric_limits<double>::signaling_NaN()}};
267 :
268 0 : friend bool operator==(const SlabJet& lhs, const SlabJet& rhs);
269 :
270 0 : friend bool operator!=(const SlabJet& lhs, const SlabJet& rhs);
271 : };
272 :
273 : } // namespace grmhd::AnalyticData
|
