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 <cstddef>
8 : #include <optional>
9 : #include <string>
10 :
11 : #include "DataStructures/Tensor/TypeAliases.hpp"
12 : #include "NumericalAlgorithms/LinearOperators/PartialDerivatives.hpp"
13 : #include "Options/Options.hpp"
14 : #include "Options/String.hpp"
15 : #include "PointwiseFunctions/AnalyticSolutions/AnalyticSolution.hpp"
16 : #include "PointwiseFunctions/AnalyticSolutions/RelativisticEuler/Solutions.hpp"
17 : #include "PointwiseFunctions/Hydro/EquationsOfState/Factory.hpp"
18 : #include "PointwiseFunctions/Hydro/Tags.hpp"
19 : #include "PointwiseFunctions/Hydro/Temperature.hpp"
20 : #include "PointwiseFunctions/InitialDataUtilities/InitialData.hpp"
21 : #include "Utilities/TMPL.hpp"
22 : #include "Utilities/TaggedTuple.hpp"
23 :
24 : namespace RelativisticEuler::Solutions {
25 : namespace detail {
26 : /*!
27 : * \brief Read the CST rotating neutron star solution from file, rescaling the
28 : * solution. If the EoS in RotNS is set to PT (polytrope), then the user
29 : * should use `is_polytrope=true`, so that the normalization is done using the
30 : * method in \cite Cook1992. In this case, RotNS already uses geometric units,
31 : * letting \f$K\f$ (i.e. \f$P=K\rho^\gamma\f$) be a free parameter. RotNS
32 : * assumes \f$K=1\f$ for its calculations. If the user doesn't want to rescale,
33 : * they can just set `polytropic_constant=1.` and nothing will change.
34 : *
35 : * If any other EoS is used, then it's done according to the method in
36 : * \cite Cook1994. In this case, RotNS instead normalizes with CGS, so
37 : * `is_polytrope=false` converts from CGS into geometric units. The user would
38 : * not want to use this for a PT EoS, even if they don't want to rescale.
39 : */
40 : class CstSolution {
41 : public:
42 : CstSolution() = default;
43 : CstSolution(const std::string& filename, bool is_polytrope,
44 : double polytropic_constant);
45 :
46 : std::array<double, 6> interpolate(double target_radius,
47 : double target_cos_theta,
48 : bool interpolate_hydro_vars) const;
49 :
50 : double polytropic_index() const { return polytropic_index_; }
51 :
52 : double equatorial_radius() const { return equatorial_radius_; }
53 :
54 : void pup(PUP::er& p);
55 :
56 : private:
57 : double maximum_radius_{std::numeric_limits<double>::signaling_NaN()};
58 : double max_density_ratio_for_linear_interpolation_ = 1.0e2;
59 :
60 : double equatorial_radius_{std::numeric_limits<double>::signaling_NaN()};
61 : double polytropic_index_{std::numeric_limits<double>::signaling_NaN()};
62 : double central_angular_speed_{std::numeric_limits<double>::signaling_NaN()};
63 : double rotation_profile_{std::numeric_limits<double>::signaling_NaN()};
64 : size_t num_radial_points_{std::numeric_limits<size_t>::max()};
65 : size_t num_angular_points_{std::numeric_limits<size_t>::max()};
66 : size_t num_grid_points_{std::numeric_limits<size_t>::max()};
67 :
68 : DataVector radius_;
69 : DataVector cos_theta_; // Note that cos(theta) is between 0 and 1.
70 : DataVector rest_mass_density_;
71 : DataVector fluid_velocity_;
72 : DataVector alpha_;
73 : DataVector rho_;
74 : DataVector gamma_;
75 : DataVector omega_;
76 : };
77 : } // namespace detail
78 :
79 : /*!
80 : * \brief A solution obtained by reading in rotating neutron star initial data
81 : * from the RotNS code based on \cite Cook1992 and \cite Cook1994.
82 : *
83 : * The code that generates the initial data is part of a private SXS repository
84 : * called `RotNS`.
85 : *
86 : * The metric in spherical coordinates is given by \cite Cook1992
87 : *
88 : * \f{align}{
89 : * ds^2=-e^{\gamma+\rho}dt^2+e^{2\alpha}(dr^2+r^2d\theta^2)
90 : * +e^{\gamma-\rho}r^2\sin^2(\theta)(d\phi-\omega dt)^2.
91 : * \f}
92 : *
93 : * We use rotation about the \f$z\f$-axis. That is,
94 : *
95 : * \f{align}{
96 : * g_{tt}
97 : * &=-e^{\gamma+\rho} + e^{\gamma-\rho}r^2\sin^2(\theta)\omega^2 \\
98 : * g_{rr}
99 : * &=e^{2\alpha} \\
100 : * g_{\theta\theta}
101 : * &=e^{2\alpha}r^2 \\
102 : * g_{\phi\phi}
103 : * &=e^{\gamma-\rho}r^2\sin^2(\theta) \\
104 : * g_{t\phi}
105 : * &=-e^{\gamma-\rho}r^2\sin^2(\theta)\omega.
106 : * \f}
107 : *
108 : * We can transform from spherical to Cartesian coordinates using
109 : *
110 : * \f{align}{
111 : * \label{eq:Jacobian}
112 : * \frac{\partial (r,\theta,\phi)}{\partial (x,y,z)}=
113 : * \begin{pmatrix}
114 : * \cos(\phi) \sin(\theta) & \sin(\theta)\sin(\phi) & \cos(\theta) \\
115 : * \tfrac{\cos(\phi)\cos(\theta)}{r} & \tfrac{\sin(\phi)\cos(\theta)}{r}
116 : * & -\tfrac{\sin(\theta)}{r} \\
117 : * -\tfrac{\sin(\phi)}{r\sin(\theta)} & \tfrac{\cos(\phi)}{r\sin(\theta)}
118 : * & 0
119 : * \end{pmatrix}
120 : * \f}
121 : *
122 : * and
123 : *
124 : * \f{align}{
125 : * \frac{\partial (x,y,z)}{\partial (r,\theta,\phi)}=
126 : * \begin{pmatrix}
127 : * \cos(\phi) \sin(\theta) & r\cos(\phi)\cos(\theta) &
128 : * -r\sin(\theta)\sin(\phi)
129 : * \\
130 : * \sin(\theta)\sin(\phi) & r\sin(\phi)\cos(\theta) &
131 : * r\cos(\phi)\sin(\theta)
132 : * \\
133 : * \cos(\theta) & -r\sin(\theta) & 0
134 : * \end{pmatrix}
135 : * \f}
136 : *
137 : *
138 : * We denote the lapse as \f$N\f$ since \f$\alpha\f$ is already being used,
139 : *
140 : * \f{align}{
141 : * N = e^{(\gamma+\rho)/2}.
142 : * \f}
143 : *
144 : * The shift is
145 : *
146 : * \f{align}{
147 : * \beta^\phi =& -\omega
148 : * \f}
149 : *
150 : * so
151 : *
152 : * \f{align}{
153 : * \beta^x &= \partial_\phi x \beta^\phi = \sin(\phi)\omega r \sin(\theta), \\
154 : * \beta^y &= \partial_\phi y \beta^\phi = -\cos(\phi)\omega r \sin(\theta),
155 : * \\ \beta^z &= 0. \f}
156 : *
157 : * The spatial metric is
158 : *
159 : * \f{align}{
160 : * \gamma_{xx}
161 : * &= \sin^2(\theta)\cos^2(\phi) e^{(2 \alpha)}
162 : * + \cos^2(\theta)\cos^2(\phi) e^{(2 \alpha)}
163 : * + \sin^2(\phi) e^{(\gamma-\rho)} \notag \\
164 : * &=\cos^2(\phi)e^{2\alpha} + \sin^2(\phi) e^{(\gamma-\rho)} \\
165 : * \gamma_{xy}
166 : * &= \sin^2(\theta)\cos(\phi)\sin(\phi) e^{(2 \alpha)}
167 : * + \cos^2(\theta)\cos(\phi)\sin(\phi) e^{(2 \alpha)}
168 : * - \sin(\phi)\cos(\phi) e^{(\gamma-\rho)} \notag \\
169 : * &=\cos(\phi)\sin(\phi)e^{2\alpha} - \sin(\phi)\cos(\phi) e^{(\gamma-\rho)}
170 : * \\ \gamma_{yy}
171 : * &= \sin^2(\theta)\sin^2(\phi) e^{(2 \alpha)}
172 : * + \cos^2(\theta)\sin^2(\phi) e^{(2 \alpha)}
173 : * + \cos^2(\phi) e^{(\gamma-\rho)} \notag \\
174 : * &=\sin^2(\phi)e^{2\alpha} + \cos^2(\phi) e^{(\gamma-\rho)} \\
175 : * \gamma_{xz}
176 : * &= \sin(\theta)\cos(\phi)\cos(\theta) e^{2\alpha}
177 : * -\cos(\theta)\cos(\phi)\sin(\theta)e^{2\alpha} = 0 \\
178 : * \gamma_{yz}
179 : * &= \sin(\theta)\sin(\phi)\cos(\theta) e^{2\alpha}
180 : * -\cos(\theta)\sin(\phi)\sin(\theta)e^{2\alpha} = 0 \\
181 : * \gamma_{zz}
182 : * &=\cos^2(\theta)e^{2\alpha} + \sin^2(\theta) e^{2\alpha}
183 : * = e^{2\alpha}
184 : * \f}
185 : *
186 : * and its determinant is
187 : *
188 : * \f{align}{
189 : * \gamma = e^{4\alpha + (\gamma-\rho)} = e^{4\alpha}e^{(\gamma-\rho)}.
190 : * \f}
191 : *
192 : * At \f$r=0\f$ we have \f$2\alpha=\gamma-\rho\f$ and so the
193 : * \f$\gamma_{xx}=\gamma_{yy}=\gamma_{zz} = e^{2\alpha}\f$ and all other
194 : * components are zero. The inverse spatial metric is given by
195 : *
196 : * \f{align}{
197 : * \gamma^{xx}
198 : * &= \frac{\gamma_{yy}}{e^{2\alpha}e^{(\gamma-\rho)}} =
199 : * \left[\sin^2(\phi)e^{2\alpha} + \cos^2(\phi) e^{(\gamma-\rho)}\right]
200 : * e^{-2\alpha} e^{-(\gamma-\rho)} \notag \\
201 : * &=\sin^2(\phi)e^{-(\gamma-\rho)} + \cos^2(\phi) e^{-2\alpha} \\
202 : * \gamma^{yy}
203 : * &= \frac{\gamma_{xx}}{e^{2\alpha}e^{(\gamma-\rho)}} =
204 : * \left[\cos^2(\phi)e^{2\alpha} + \sin^2(\phi) e^{(\gamma-\rho)}\right]
205 : * e^{-2\alpha} e^{-(\gamma-\rho)} \notag \\
206 : * &=\cos^2(\phi) e^{-(\gamma-\rho)} + \sin^2(\phi) e^{-2\alpha} \\
207 : * \gamma^{xy}
208 : * &=\frac{-\gamma_{xy}}{e^{2\alpha}e^{(\gamma-\rho)}} =
209 : * -\left[\cos(\phi)\sin(\phi)e^{2\alpha} - \sin(\phi)\cos(\phi)
210 : * e^{(\gamma-\rho)}\right] e^{-2\alpha} e^{-(\gamma-\rho)} \notag \\
211 : * &=-\cos(\phi)\sin(\phi)e^{-(\gamma-\rho)} -
212 : * \sin(\phi)\cos(\phi)e^{-2\alpha} \notag \\
213 : * &=\cos(\phi)\sin(\phi) \left[e^{-2\alpha} -
214 : * e^{-(\gamma-\rho)}\right] \\
215 : * \gamma^{xz}
216 : * &= 0 \\
217 : * \gamma^{yz}
218 : * &= 0 \\
219 : * \gamma^{zz} &= e^{-2\alpha}.
220 : * \f}
221 : *
222 : * The 4-velocity in spherical coordinates is given by
223 : *
224 : * \f{align}{
225 : * u^{\bar{a}}=\frac{e^{-(\rho+\gamma)/2}}{\sqrt{1-v^2}}
226 : * \left[1,0,0,\Omega\right],
227 : * \f}
228 : *
229 : * where
230 : *
231 : * \f{align}{
232 : * v=(\Omega-\omega)r\sin(\theta)e^{-\rho}.
233 : * \f}
234 : *
235 : * Transforming to Cartesian coordinates we have
236 : *
237 : * \f{align}{
238 : * u^t
239 : * &=\frac{e^{-(\rho+\gamma)/2}}{\sqrt{1-v^2}} \\
240 : * u^x
241 : * &=\partial_\phi x u^\phi = -r\sin(\theta)\sin(\phi) u^t\Omega \\
242 : * u^y
243 : * &=\partial_\phi y u^\phi = r\sin(\theta)\cos(\phi) u^t\Omega \\
244 : * u^z &= 0.
245 : * \f}
246 : *
247 : * The Lorentz factor is given by
248 : *
249 : * \f{align}{
250 : * W
251 : * &=Nu^t=e^{(\gamma+\rho)/2}\frac{e^{-(\rho+\gamma)/2}}{\sqrt{1-v^2}} \notag
252 : * \\
253 : * &=\frac{1}{\sqrt{1-v^2}}.
254 : * \f}
255 : *
256 : * Using
257 : *
258 : * \f{align}{
259 : * v^i = \frac{1}{N}\left(\frac{u^i}{u^t} + \beta^i\right)
260 : * \f}
261 : *
262 : * we get
263 : *
264 : * \f{align}{
265 : * v^x
266 : * &= -e^{-(\gamma+\rho)/2}r\sin(\theta)\sin(\phi)(\Omega-\omega)
267 : * =-e^{-(\gamma-\rho)/2}\sin(\phi) v\\
268 : * v^y
269 : * &= e^{-(\gamma+\rho)/2}r\sin(\theta)\cos(\phi)(\Omega-\omega)
270 : * = e^{-(\gamma-\rho)/2}\cos(\phi)v \\
271 : * v^z&=0.
272 : * \f}
273 : *
274 : * Lowering with the spatial metric we get
275 : *
276 : * \f{align}{
277 : * v_x
278 : * &=\gamma_{xx} v^x + \gamma_{xy} v^y \notag \\
279 : * &=-\left[\cos^2(\phi)e^{2\alpha}+\sin^2(\phi)e^{\gamma-\rho}\right]
280 : * e^{-(\gamma-\rho)/2}\sin(\phi) v \notag \\
281 : * &+\left[\cos(\phi)\sin(\phi)e^{2\alpha} -
282 : * \sin(\phi)\cos(\phi)e^{\gamma-\rho}\right]
283 : * e^{-(\gamma-\rho)/2}\cos(\phi)v \notag \\
284 : * &=-e^{-(\gamma-\rho)/2}v\sin(\phi)e^{\gamma-\rho}
285 : * \left[\sin^2(\phi)+\cos^2(\phi)\right] \notag \\
286 : * &=-e^{(\gamma-\rho)/2}v\sin(\phi) \\
287 : * v_y
288 : * &=\gamma_{yx} v^x + \gamma_{yy} v^y \notag \\
289 : * &=-\left[\cos(\phi)\sin(\phi)e^{2\alpha} - \sin(\phi)\cos(\phi)
290 : * e^{(\gamma-\rho)}\right] e^{-(\gamma-\rho)/2}\sin(\phi) v \notag \\
291 : * &+\left[\sin^2(\phi)e^{2\alpha} + \cos^2(\phi) e^{(\gamma-\rho)}\right]
292 : * e^{-(\gamma-\rho)/2}\cos(\phi)v \notag \\
293 : * &=e^{(\gamma-\rho)/2}v\cos(\phi) \\
294 : * v_z &= 0.
295 : * \f}
296 : *
297 : * This is consistent with the Lorentz factor read off from \f$u^t\f$ since
298 : * \f$v^iv_i=v^2\f$. For completeness, \f$u_i=Wv_i\f$ so
299 : *
300 : * \f{align}{
301 : * u_x
302 : * &=-\frac{e^{(\gamma-\rho)/2}v\sin(\phi)}{\sqrt{1-v^2}} \\
303 : * u_y
304 : * &=\frac{e^{(\gamma-\rho)/2}v\cos(\phi)}{\sqrt{1-v^2}} \\
305 : * u_z&=0.
306 : * \f}
307 : *
308 : * \warning Near (within `1e-2`) \f$r=0\f$ the numerical errors from
309 : * interpolation and computing the metric derivatives by finite difference no
310 : * longer cancel out and so the `tilde_s` time derivative only vanishes to
311 : * roughly `1e-8` rather than machine precision. Computing the Cartesian
312 : * derivatives from analytic differentiation of the radial and angular
313 : * polynomial fits might improve the situation but is a decent about of work to
314 : * implement.
315 : */
316 1 : class RotatingStar : public virtual evolution::initial_data::InitialData,
317 : public MarkAsAnalyticSolution,
318 : public AnalyticSolution<3>,
319 : public hydro::TemperatureInitialization<RotatingStar> {
320 : template <typename DataType>
321 0 : struct IntermediateVariables {
322 0 : IntermediateVariables(
323 : const tnsr::I<DataType, 3, Frame::Inertial>& in_coords,
324 : double in_delta_r);
325 :
326 0 : struct MetricData {
327 0 : MetricData() = default;
328 0 : MetricData(size_t num_points);
329 0 : DataType alpha;
330 0 : DataType rho;
331 0 : DataType gamma;
332 0 : DataType omega;
333 : };
334 :
335 0 : const tnsr::I<DataType, 3, Frame::Inertial>& coords;
336 0 : DataType radius;
337 0 : DataType phi;
338 0 : DataType cos_theta;
339 0 : DataType sin_theta;
340 : // Data at coords
341 0 : std::optional<DataType> rest_mass_density;
342 0 : std::optional<DataType> fluid_velocity;
343 0 : std::optional<MetricData> metric_data;
344 : // Data for 2nd-order FD derivatives.
345 : //
346 : // Technically we could do full-order accurate derivatives by using
347 : // Jacobians to transform quantities, but since we really want the initial
348 : // data to be solved natively in SpECTRE and satisfy the Einstein equations
349 : // with neutrinos and magnetic fields, doing the 2nd order FD like SpEC
350 : // should be fine.
351 : //
352 : // Note: We guard all the upper/lower values by checking if
353 : // metric_data_upper is computed. While not perfect, this simplifies the
354 : // code a lot.
355 0 : double delta_r;
356 0 : std::optional<std::array<MetricData, 3>> metric_data_upper{};
357 0 : std::optional<std::array<MetricData, 3>> metric_data_lower{};
358 0 : std::array<DataType, 3> radius_upper{};
359 0 : std::array<DataType, 3> radius_lower{};
360 0 : std::array<DataType, 3> cos_theta_upper{};
361 0 : std::array<DataType, 3> cos_theta_lower{};
362 0 : std::array<DataType, 3> sin_theta_upper{};
363 0 : std::array<DataType, 3> sin_theta_lower{};
364 0 : std::array<DataType, 3> phi_upper{};
365 0 : std::array<DataType, 3> phi_lower{};
366 : };
367 :
368 : public:
369 : /// The path to the RotNS data file.
370 1 : struct RotNsFilename {
371 0 : using type = std::string;
372 0 : static constexpr Options::String help = {
373 : "The path to the RotNS data file."};
374 : };
375 :
376 : /// The polytropic constant of the fluid.
377 : ///
378 : /// The data in the RotNS file will be rescaled.
379 1 : struct PolytropicConstant {
380 0 : using type = double;
381 0 : static constexpr Options::String help = {
382 : "The polytropic constant of the fluid."};
383 0 : static type lower_bound() { return 0.; }
384 : };
385 :
386 0 : using options = tmpl::list<Options::Alternatives<
387 : tmpl::list<RotNsFilename, PolytropicConstant>,
388 : tmpl::list<RotNsFilename,
389 : hydro::OptionTags::InitialDataEquationOfState<true, 1>>>>;
390 :
391 0 : static constexpr Options::String help = {
392 : "Rotating neutron star initial data solved by the RotNS solver. The data "
393 : "is read in from disk. If RotNS used a polytropic equation of state, use "
394 : "the PolytropicConstant option, so that the polytropic index may be read "
395 : "directly from the RotNS output. Otherwise, set the equation of state to "
396 : "whichever one you would like to use to load the initial data. Note that "
397 : "if a polytrope is put into the EquationOfState option rather than using "
398 : "the PolytropicConstant option, the generic unit conversion will be used "
399 : "instead of the specific polytrope one."};
400 :
401 0 : RotatingStar() = default;
402 0 : RotatingStar(const RotatingStar& /*rhs*/);
403 0 : RotatingStar& operator=(const RotatingStar& /*rhs*/);
404 0 : RotatingStar(RotatingStar&& /*rhs*/) = default;
405 0 : RotatingStar& operator=(RotatingStar&& /*rhs*/) = default;
406 0 : ~RotatingStar() override = default;
407 :
408 0 : RotatingStar(std::string rot_ns_filename,
409 : std::unique_ptr<EquationsOfState::EquationOfState<true, 1>>
410 : equation_of_state);
411 0 : RotatingStar(std::string rot_ns_filename, double polytropic_constant);
412 :
413 0 : auto get_clone() const
414 : -> std::unique_ptr<evolution::initial_data::InitialData> override;
415 :
416 : /// \cond
417 : explicit RotatingStar(CkMigrateMessage* msg);
418 : using PUP::able::register_constructor;
419 : WRAPPED_PUPable_decl_template(RotatingStar);
420 : /// \endcond
421 :
422 : /// Retrieve a collection of variables at `(x, t)`
423 : template <typename DataType, typename... Tags>
424 1 : tuples::TaggedTuple<Tags...> variables(const tnsr::I<DataType, 3>& x,
425 : const double /*t*/,
426 : tmpl::list<Tags...> /*meta*/) const {
427 : IntermediateVariables<DataType> intermediate_vars{
428 : x, 1.0e-4 * cst_solution_.equatorial_radius()};
429 : return {get<Tags>(variables(make_not_null(&intermediate_vars), x,
430 : tmpl::list<Tags>{}))...};
431 : }
432 :
433 : // NOLINTNEXTLINE(google-runtime-references)
434 0 : void pup(PUP::er& p) override;
435 :
436 0 : const EquationsOfState::EquationOfState<true, 1>& equation_of_state() const {
437 : return *equation_of_state_;
438 : }
439 :
440 0 : double equatorial_radius() const {
441 : return cst_solution_.equatorial_radius();
442 : }
443 :
444 : protected:
445 : template <typename DataType>
446 0 : using DerivLapse = ::Tags::deriv<gr::Tags::Lapse<DataType>, tmpl::size_t<3>,
447 : Frame::Inertial>;
448 :
449 : template <typename DataType>
450 0 : using DerivShift = ::Tags::deriv<gr::Tags::Shift<DataType, 3>,
451 : tmpl::size_t<3>, Frame::Inertial>;
452 :
453 : template <typename DataType>
454 0 : using DerivSpatialMetric = ::Tags::deriv<gr::Tags::SpatialMetric<DataType, 3>,
455 : tmpl::size_t<3>, Frame::Inertial>;
456 :
457 : template <typename DataType>
458 0 : void interpolate_vars_if_necessary(
459 : gsl::not_null<IntermediateVariables<DataType>*> vars) const;
460 :
461 : template <typename DataType>
462 0 : void interpolate_deriv_vars_if_necessary(
463 : gsl::not_null<IntermediateVariables<DataType>*> vars) const;
464 :
465 : template <typename DataType>
466 0 : Scalar<DataType> lapse(const DataType& gamma, const DataType& rho) const;
467 :
468 : template <typename DataType>
469 0 : tnsr::I<DataType, 3, Frame::Inertial> shift(const DataType& omega,
470 : const DataType& phi,
471 : const DataType& radius,
472 : const DataType& sin_theta) const;
473 :
474 : template <typename DataType>
475 0 : tnsr::ii<DataType, 3, Frame::Inertial> spatial_metric(
476 : const DataType& gamma, const DataType& rho, const DataType& alpha,
477 : const DataType& phi) const;
478 :
479 : template <typename DataType>
480 0 : tnsr::II<DataType, 3, Frame::Inertial> inverse_spatial_metric(
481 : const DataType& gamma, const DataType& rho, const DataType& alpha,
482 : const DataType& phi) const;
483 :
484 : template <typename DataType>
485 0 : Scalar<DataType> sqrt_det_spatial_metric(const DataType& gamma,
486 : const DataType& rho,
487 : const DataType& alpha) const;
488 :
489 : template <typename DataType>
490 0 : auto make_metric_data(size_t num_points) const ->
491 : typename IntermediateVariables<DataType>::MetricData;
492 :
493 : template <typename DataType>
494 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
495 : const tnsr::I<DataType, 3>& x,
496 : tmpl::list<hydro::Tags::RestMassDensity<DataType>> /*meta*/)
497 : const -> tuples::TaggedTuple<hydro::Tags::RestMassDensity<DataType>>;
498 :
499 : template <typename DataType>
500 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
501 : const tnsr::I<DataType, 3>& x,
502 : tmpl::list<hydro::Tags::ElectronFraction<DataType>> /*meta*/)
503 : const -> tuples::TaggedTuple<hydro::Tags::ElectronFraction<DataType>>;
504 :
505 : template <typename DataType>
506 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
507 : const tnsr::I<DataType, 3>& x,
508 : tmpl::list<hydro::Tags::SpecificEnthalpy<DataType>> /*meta*/)
509 : const -> tuples::TaggedTuple<hydro::Tags::SpecificEnthalpy<DataType>>;
510 :
511 : template <typename DataType>
512 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
513 : const tnsr::I<DataType, 3>& x,
514 : tmpl::list<hydro::Tags::Temperature<DataType>> /*meta*/) const
515 : -> tuples::TaggedTuple<hydro::Tags::Temperature<DataType>>;
516 :
517 : template <typename DataType>
518 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
519 : const tnsr::I<DataType, 3>& x,
520 : tmpl::list<hydro::Tags::Pressure<DataType>> /*meta*/) const
521 : -> tuples::TaggedTuple<hydro::Tags::Pressure<DataType>>;
522 :
523 : template <typename DataType>
524 0 : auto variables(
525 : gsl::not_null<IntermediateVariables<DataType>*> vars,
526 : const tnsr::I<DataType, 3>& x,
527 : tmpl::list<hydro::Tags::SpecificInternalEnergy<DataType>> /*meta*/) const
528 : -> tuples::TaggedTuple<hydro::Tags::SpecificInternalEnergy<DataType>>;
529 :
530 : template <typename DataType>
531 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
532 : const tnsr::I<DataType, 3>& x,
533 : tmpl::list<hydro::Tags::SpatialVelocity<DataType, 3>> /*meta*/)
534 : const -> tuples::TaggedTuple<hydro::Tags::SpatialVelocity<DataType, 3>>;
535 :
536 : template <typename DataType>
537 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
538 : const tnsr::I<DataType, 3>& x,
539 : tmpl::list<hydro::Tags::LorentzFactor<DataType>> /*meta*/)
540 : const -> tuples::TaggedTuple<hydro::Tags::LorentzFactor<DataType>>;
541 :
542 : template <typename DataType>
543 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
544 : const tnsr::I<DataType, 3>& x,
545 : tmpl::list<hydro::Tags::MagneticField<DataType, 3>> /*meta*/)
546 : const -> tuples::TaggedTuple<hydro::Tags::MagneticField<DataType, 3>>;
547 :
548 : template <typename DataType>
549 0 : auto variables(
550 : gsl::not_null<IntermediateVariables<DataType>*> vars,
551 : const tnsr::I<DataType, 3>& x,
552 : tmpl::list<hydro::Tags::DivergenceCleaningField<DataType>> /*meta*/) const
553 : -> tuples::TaggedTuple<hydro::Tags::DivergenceCleaningField<DataType>>;
554 :
555 : template <typename DataType>
556 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
557 : const tnsr::I<DataType, 3>& x,
558 : tmpl::list<gr::Tags::Lapse<DataType>> /*meta*/) const
559 : -> tuples::TaggedTuple<gr::Tags::Lapse<DataType>>;
560 :
561 : template <typename DataType>
562 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
563 : const tnsr::I<DataType, 3>& x,
564 : tmpl::list<gr::Tags::Shift<DataType, 3>> /*meta*/) const
565 : -> tuples::TaggedTuple<gr::Tags::Shift<DataType, 3>>;
566 :
567 : template <typename DataType>
568 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
569 : const tnsr::I<DataType, 3>& x,
570 : tmpl::list<gr::Tags::SpatialMetric<DataType, 3>> /*meta*/)
571 : const -> tuples::TaggedTuple<gr::Tags::SpatialMetric<DataType, 3>>;
572 :
573 : template <typename DataType>
574 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
575 : const tnsr::I<DataType, 3>& x,
576 : tmpl::list<gr::Tags::SqrtDetSpatialMetric<DataType>> /*meta*/)
577 : const -> tuples::TaggedTuple<gr::Tags::SqrtDetSpatialMetric<DataType>>;
578 :
579 : template <typename DataType>
580 0 : auto variables(
581 : gsl::not_null<IntermediateVariables<DataType>*> vars,
582 : const tnsr::I<DataType, 3>& x,
583 : tmpl::list<gr::Tags::InverseSpatialMetric<DataType, 3>> /*meta*/) const
584 : -> tuples::TaggedTuple<gr::Tags::InverseSpatialMetric<DataType, 3>>;
585 :
586 : template <typename DataType>
587 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
588 : const tnsr::I<DataType, 3>& x,
589 : tmpl::list<::Tags::dt<gr::Tags::Lapse<DataType>>> /*meta*/)
590 : const -> tuples::TaggedTuple<::Tags::dt<gr::Tags::Lapse<DataType>>>;
591 :
592 : template <typename DataType>
593 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
594 : const tnsr::I<DataType, 3>& x,
595 : tmpl::list<DerivLapse<DataType>> /*meta*/) const
596 : -> tuples::TaggedTuple<DerivLapse<DataType>>;
597 :
598 : template <typename DataType>
599 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
600 : const tnsr::I<DataType, 3>& x,
601 : tmpl::list<::Tags::dt<gr::Tags::Shift<DataType, 3>>> /*meta*/)
602 : const -> tuples::TaggedTuple<::Tags::dt<gr::Tags::Shift<DataType, 3>>>;
603 :
604 : template <typename DataType>
605 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
606 : const tnsr::I<DataType, 3>& x,
607 : tmpl::list<DerivShift<DataType>> /*meta*/) const
608 : -> tuples::TaggedTuple<DerivShift<DataType>>;
609 :
610 : template <typename DataType>
611 0 : auto variables(
612 : gsl::not_null<IntermediateVariables<DataType>*> vars,
613 : const tnsr::I<DataType, 3>& x,
614 : tmpl::list<::Tags::dt<gr::Tags::SpatialMetric<DataType, 3>>> /*meta*/)
615 : const
616 : -> tuples::TaggedTuple<::Tags::dt<gr::Tags::SpatialMetric<DataType, 3>>>;
617 :
618 : template <typename DataType>
619 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
620 : const tnsr::I<DataType, 3>& x,
621 : tmpl::list<DerivSpatialMetric<DataType>> /*meta*/) const
622 : -> tuples::TaggedTuple<DerivSpatialMetric<DataType>>;
623 :
624 : template <typename DataType>
625 0 : auto variables(gsl::not_null<IntermediateVariables<DataType>*> vars,
626 : const tnsr::I<DataType, 3>& x,
627 : tmpl::list<gr::Tags::ExtrinsicCurvature<DataType, 3>> /*meta*/)
628 : const -> tuples::TaggedTuple<gr::Tags::ExtrinsicCurvature<DataType, 3>>;
629 :
630 0 : friend bool operator==(const RotatingStar& lhs, const RotatingStar& rhs);
631 :
632 0 : std::string rot_ns_filename_{};
633 0 : detail::CstSolution cst_solution_{};
634 0 : double polytropic_constant_ = std::numeric_limits<double>::signaling_NaN();
635 0 : double polytropic_exponent_ = std::numeric_limits<double>::signaling_NaN();
636 0 : bool is_polytrope_{};
637 : std::unique_ptr<EquationsOfState::EquationOfState<true, 1>>
638 0 : equation_of_state_;
639 : // Floor value to protect EoS from encountering FPEs when computing state
640 : // variables in the atmosphere
641 0 : static constexpr double atmosphere_floor_ = 1.e-50;
642 : };
643 :
644 0 : bool operator!=(const RotatingStar& lhs, const RotatingStar& rhs);
645 : } // namespace RelativisticEuler::Solutions
|
