Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <limits>
7 : #include <ostream>
8 :
9 : #include "DataStructures/CachedTempBuffer.hpp"
10 : #include "DataStructures/DataBox/Prefixes.hpp"
11 : #include "DataStructures/DataVector.hpp"
12 : #include "DataStructures/TaggedTuple.hpp"
13 : #include "DataStructures/Tensor/Tensor.hpp"
14 : #include "Elliptic/Systems/Xcts/Tags.hpp"
15 : #include "NumericalAlgorithms/LinearOperators/PartialDerivatives.hpp"
16 : #include "Options/String.hpp"
17 : #include "PointwiseFunctions/AnalyticSolutions/RelativisticEuler/TovStar.hpp"
18 : #include "PointwiseFunctions/AnalyticSolutions/Xcts/CommonVariables.hpp"
19 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
20 : #include "PointwiseFunctions/GeneralRelativity/Tags/Conformal.hpp"
21 : #include "PointwiseFunctions/InitialDataUtilities/AnalyticSolution.hpp"
22 : #include "Utilities/Gsl.hpp"
23 : #include "Utilities/Serialization/CharmPupable.hpp"
24 : #include "Utilities/TMPL.hpp"
25 :
26 : /// \cond
27 : namespace PUP {
28 : class er;
29 : } // namespace PUP
30 : /// \endcond
31 :
32 : namespace Xcts::Solutions {
33 : namespace tov_detail {
34 :
35 : using TovCoordinates = RelativisticEuler::Solutions::TovCoordinates;
36 :
37 : template <typename DataType>
38 : using TovVariablesCache = cached_temp_buffer_from_typelist<tmpl::push_back<
39 : common_tags<DataType>,
40 : ::Tags::deriv<gr::Tags::Lapse<DataType>, tmpl::size_t<3>, Frame::Inertial>,
41 : gr::Tags::Conformal<gr::Tags::EnergyDensity<DataType>, 0>,
42 : gr::Tags::Conformal<gr::Tags::StressTrace<DataType>, 0>,
43 : gr::Tags::Conformal<gr::Tags::MomentumDensity<DataType, 3>, 0>>>;
44 :
45 : template <typename DataType>
46 : struct TovVariables : CommonVariables<DataType, TovVariablesCache<DataType>> {
47 : static constexpr size_t Dim = 3;
48 : static constexpr int ConformalMatterScale = 0;
49 : using Cache = TovVariablesCache<DataType>;
50 : using Base = CommonVariables<DataType, TovVariablesCache<DataType>>;
51 : using Base::operator();
52 :
53 : const tnsr::I<DataType, 3>& x;
54 : const DataType& radius;
55 : const RelativisticEuler::Solutions::TovStar& tov_star;
56 :
57 : TovVariables(
58 : std::optional<std::reference_wrapper<const Mesh<Dim>>> local_mesh,
59 : std::optional<std::reference_wrapper<const InverseJacobian<
60 : DataType, Dim, Frame::ElementLogical, Frame::Inertial>>>
61 : local_inv_jacobian,
62 : const tnsr::I<DataType, 3>& local_x, const DataType& local_radius,
63 : const RelativisticEuler::Solutions::TovStar& local_tov_star)
64 : : Base(std::move(local_mesh), std::move(local_inv_jacobian)),
65 : x(local_x),
66 : radius(local_radius),
67 : tov_star(local_tov_star) {}
68 :
69 : void operator()(gsl::not_null<tnsr::ii<DataType, 3>*> conformal_metric,
70 : gsl::not_null<Cache*> cache,
71 : Tags::ConformalMetric<DataType, 3, Frame::Inertial> /*meta*/)
72 : const override;
73 : void operator()(
74 : gsl::not_null<tnsr::II<DataType, 3>*> inv_conformal_metric,
75 : gsl::not_null<Cache*> cache,
76 : Tags::InverseConformalMetric<DataType, 3, Frame::Inertial> /*meta*/)
77 : const override;
78 : void operator()(
79 : gsl::not_null<tnsr::ijj<DataType, 3>*> deriv_conformal_metric,
80 : gsl::not_null<Cache*> cache,
81 : ::Tags::deriv<Tags::ConformalMetric<DataType, 3, Frame::Inertial>,
82 : tmpl::size_t<3>, Frame::Inertial> /*meta*/) const override;
83 : void operator()(
84 : gsl::not_null<tnsr::ii<DataType, 3>*> extrinsic_curvature,
85 : gsl::not_null<Cache*> cache,
86 : gr::Tags::ExtrinsicCurvature<DataType, 3> /*meta*/) const override;
87 : void operator()(
88 : gsl::not_null<Scalar<DataType>*> trace_extrinsic_curvature,
89 : gsl::not_null<Cache*> cache,
90 : gr::Tags::TraceExtrinsicCurvature<DataType> /*meta*/) const override;
91 : void operator()(
92 : gsl::not_null<tnsr::i<DataType, 3>*> deriv_trace_extrinsic_curvature,
93 : gsl::not_null<Cache*> cache,
94 : ::Tags::deriv<gr::Tags::TraceExtrinsicCurvature<DataType>,
95 : tmpl::size_t<3>, Frame::Inertial> /*meta*/) const override;
96 : void operator()(
97 : gsl::not_null<Scalar<DataType>*> dt_trace_extrinsic_curvature,
98 : gsl::not_null<Cache*> cache,
99 : ::Tags::dt<gr::Tags::TraceExtrinsicCurvature<DataType>> /*meta*/)
100 : const override;
101 : void operator()(gsl::not_null<Scalar<DataType>*> conformal_factor,
102 : gsl::not_null<Cache*> cache,
103 : Tags::ConformalFactor<DataType> /*meta*/) const override;
104 : void operator()(
105 : gsl::not_null<Scalar<DataType>*> conformal_factor_minus_one,
106 : gsl::not_null<Cache*> cache,
107 : Tags::ConformalFactorMinusOne<DataType> /*meta*/) const override;
108 : void operator()(
109 : gsl::not_null<tnsr::i<DataType, 3>*> deriv_conformal_factor,
110 : gsl::not_null<Cache*> cache,
111 : ::Tags::deriv<Xcts::Tags::ConformalFactorMinusOne<DataType>,
112 : tmpl::size_t<3>, Frame::Inertial> /*meta*/) const override;
113 : void operator()(gsl::not_null<Scalar<DataType>*> lapse,
114 : gsl::not_null<Cache*> cache,
115 : gr::Tags::Lapse<DataType> /*meta*/) const override;
116 : void operator()(gsl::not_null<tnsr::i<DataType, 3>*> deriv_lapse,
117 : gsl::not_null<Cache*> cache,
118 : ::Tags::deriv<gr::Tags::Lapse<DataType>, tmpl::size_t<3>,
119 : Frame::Inertial> /*meta*/) const;
120 : void operator()(
121 : gsl::not_null<Scalar<DataType>*> lapse_times_conformal_factor,
122 : gsl::not_null<Cache*> cache,
123 : Tags::LapseTimesConformalFactor<DataType> /*meta*/) const override;
124 : void operator()(
125 : gsl::not_null<Scalar<DataType>*> lapse_times_conformal_factor_minus_one,
126 : gsl::not_null<Cache*> cache,
127 : Tags::LapseTimesConformalFactorMinusOne<DataType> /*meta*/)
128 : const override;
129 : void operator()(
130 : gsl::not_null<tnsr::i<DataType, 3>*> deriv_lapse_times_conformal_factor,
131 : gsl::not_null<Cache*> cache,
132 : ::Tags::deriv<Tags::LapseTimesConformalFactorMinusOne<DataType>,
133 : tmpl::size_t<3>, Frame::Inertial> /*meta*/) const override;
134 : void operator()(gsl::not_null<tnsr::I<DataType, 3>*> shift_background,
135 : gsl::not_null<Cache*> cache,
136 : Tags::ShiftBackground<DataType, 3, Frame::Inertial> /*meta*/)
137 : const override;
138 : void operator()(
139 : gsl::not_null<tnsr::iJ<DataType, 3, Frame::Inertial>*>
140 : deriv_shift_background,
141 : gsl::not_null<Cache*> cache,
142 : ::Tags::deriv<Xcts::Tags::ShiftBackground<DataType, 3, Frame::Inertial>,
143 : tmpl::size_t<3>, Frame::Inertial> /*meta*/) const override;
144 : void operator()(gsl::not_null<tnsr::II<DataType, 3, Frame::Inertial>*>
145 : longitudinal_shift_background_minus_dt_conformal_metric,
146 : gsl::not_null<Cache*> cache,
147 : Tags::LongitudinalShiftBackgroundMinusDtConformalMetric<
148 : DataType, 3, Frame::Inertial> /*meta*/) const override;
149 : void operator()(
150 : gsl::not_null<tnsr::I<DataType, 3>*> shift_excess,
151 : gsl::not_null<Cache*> cache,
152 : Tags::ShiftExcess<DataType, 3, Frame::Inertial> /*meta*/) const override;
153 : void operator()(
154 : gsl::not_null<tnsr::iJ<DataType, 3>*> deriv_shift_excess,
155 : gsl::not_null<Cache*> cache,
156 : ::Tags::deriv<Tags::ShiftExcess<DataType, 3, Frame::Inertial>,
157 : tmpl::size_t<3>, Frame::Inertial> /*meta*/) const override;
158 : void operator()(gsl::not_null<Scalar<DataType>*> energy_density,
159 : gsl::not_null<Cache*> cache,
160 : gr::Tags::Conformal<gr::Tags::EnergyDensity<DataType>,
161 : ConformalMatterScale> /*meta*/) const;
162 : void operator()(gsl::not_null<Scalar<DataType>*> stress_trace,
163 : gsl::not_null<Cache*> cache,
164 : gr::Tags::Conformal<gr::Tags::StressTrace<DataType>,
165 : ConformalMatterScale> /*meta*/) const;
166 : void operator()(gsl::not_null<tnsr::I<DataType, 3>*> momentum_density,
167 : gsl::not_null<Cache*> cache,
168 : gr::Tags::Conformal<gr::Tags::MomentumDensity<DataType, 3>,
169 : ConformalMatterScale> /*meta*/) const;
170 :
171 : private:
172 : template <typename Tag>
173 : typename Tag::type get_tov_var(Tag /*meta*/) const {
174 : // Possible optimization: Access the cache of the RelEuler::TovStar solution
175 : // so its intermediate quantities don't have to be re-computed repeatedly
176 : return get<Tag>(tov_star.variables(
177 : x, std::numeric_limits<double>::signaling_NaN(), tmpl::list<Tag>{}));
178 : }
179 : };
180 :
181 : } // namespace tov_detail
182 :
183 : /*!
184 : * \brief TOV solution to the XCTS equations
185 : *
186 : * \see RelativisticEuler::Solutions::TovStar
187 : * \see gr::Solutions::TovSolution
188 : */
189 1 : class TovStar : public elliptic::analytic_data::AnalyticSolution {
190 : private:
191 0 : using RelEulerTovStar = RelativisticEuler::Solutions::TovStar;
192 :
193 : public:
194 0 : using options = RelEulerTovStar::options;
195 0 : static constexpr Options::String help = RelEulerTovStar::help;
196 :
197 0 : TovStar() = default;
198 0 : TovStar(const TovStar&) = default;
199 0 : TovStar& operator=(const TovStar&) = default;
200 0 : TovStar(TovStar&&) = default;
201 0 : TovStar& operator=(TovStar&&) = default;
202 0 : ~TovStar() = default;
203 :
204 0 : TovStar(double central_rest_mass_density,
205 : std::unique_ptr<EquationsOfState::EquationOfState<true, 1>>
206 : equation_of_state,
207 : const RelativisticEuler::Solutions::TovCoordinates coordinate_system)
208 : : tov_star(central_rest_mass_density, std::move(equation_of_state),
209 : coordinate_system) {}
210 :
211 0 : const EquationsOfState::EquationOfState<true, 1>& equation_of_state() const {
212 : return tov_star.equation_of_state();
213 : }
214 :
215 0 : const RelativisticEuler::Solutions::TovSolution& radial_solution() const {
216 : return tov_star.radial_solution();
217 : }
218 :
219 : /// \cond
220 : explicit TovStar(CkMigrateMessage* m)
221 : : elliptic::analytic_data::AnalyticSolution(m) {}
222 : using PUP::able::register_constructor;
223 : WRAPPED_PUPable_decl_template(TovStar);
224 : std::unique_ptr<elliptic::analytic_data::AnalyticSolution> get_clone()
225 : const override {
226 : return std::make_unique<TovStar>(*this);
227 : }
228 : /// \endcond
229 :
230 : template <typename DataType>
231 0 : using tags = typename tov_detail::TovVariablesCache<DataType>::tags_list;
232 :
233 : template <typename DataType, typename... RequestedTags>
234 0 : tuples::TaggedTuple<RequestedTags...> variables(
235 : const tnsr::I<DataType, 3, Frame::Inertial>& x,
236 : tmpl::list<RequestedTags...> /*meta*/) const {
237 : return variables_impl<DataType>(x, std::nullopt, std::nullopt,
238 : tmpl::list<RequestedTags...>{});
239 : }
240 :
241 : template <typename... RequestedTags>
242 0 : tuples::TaggedTuple<RequestedTags...> variables(
243 : const tnsr::I<DataVector, 3, Frame::Inertial>& x, const Mesh<3>& mesh,
244 : const InverseJacobian<DataVector, 3, Frame::ElementLogical,
245 : Frame::Inertial>& inv_jacobian,
246 : tmpl::list<RequestedTags...> /*meta*/) const {
247 : return variables_impl<DataVector>(x, mesh, inv_jacobian,
248 : tmpl::list<RequestedTags...>{});
249 : }
250 :
251 : // NOLINTNEXTLINE(google-runtime-references)
252 0 : void pup(PUP::er& p) override {
253 : elliptic::analytic_data::AnalyticSolution::pup(p);
254 : p | tov_star;
255 : }
256 :
257 : private:
258 : template <typename DataType, typename... RequestedTags>
259 0 : tuples::TaggedTuple<RequestedTags...> variables_impl(
260 : const tnsr::I<DataType, 3, Frame::Inertial>& x,
261 : std::optional<std::reference_wrapper<const Mesh<3>>> mesh,
262 : std::optional<std::reference_wrapper<const InverseJacobian<
263 : DataType, 3, Frame::ElementLogical, Frame::Inertial>>>
264 : inv_jacobian,
265 : tmpl::list<RequestedTags...> /*meta*/) const {
266 : using VarsComputer = tov_detail::TovVariables<DataType>;
267 : typename VarsComputer::Cache cache{get_size(*x.begin())};
268 : const DataType radius = get(magnitude(x));
269 : const VarsComputer computer{std::move(mesh), std::move(inv_jacobian), x,
270 : radius, tov_star};
271 : using unrequested_hydro_tags =
272 : tmpl::list_difference<hydro_tags<DataType>,
273 : tmpl::list<RequestedTags...>>;
274 : using requested_hydro_tags =
275 : tmpl::list_difference<hydro_tags<DataType>, unrequested_hydro_tags>;
276 : tuples::tagged_tuple_from_typelist<requested_hydro_tags> hydro_vars;
277 : if constexpr (not std::is_same_v<requested_hydro_tags, tmpl::list<>>) {
278 : hydro_vars =
279 : tov_star.variables(x, std::numeric_limits<double>::signaling_NaN(),
280 : requested_hydro_tags{});
281 : }
282 : const auto get_var = [&cache, &computer, &hydro_vars](auto tag_v) {
283 : using tag = std::decay_t<decltype(tag_v)>;
284 : if constexpr (tmpl::list_contains_v<hydro_tags<DataType>, tag>) {
285 : (void)cache;
286 : (void)computer;
287 : return get<tag>(hydro_vars);
288 : } else {
289 : (void)hydro_vars;
290 : return cache.get_var(computer, tag{});
291 : }
292 : };
293 : return {get_var(RequestedTags{})...};
294 : }
295 :
296 0 : friend bool operator==(const TovStar& lhs, const TovStar& rhs) {
297 : return lhs.tov_star == rhs.tov_star;
298 : }
299 :
300 : // Instead of inheriting from the RelEuler::TovStar we use an aggregate
301 : // pattern to avoid multiple-inheritance issues.
302 0 : RelativisticEuler::Solutions::TovStar tov_star{};
303 : };
304 :
305 0 : inline bool operator!=(const TovStar& lhs, const TovStar& rhs) {
306 : return not(lhs == rhs);
307 : }
308 :
309 : } // namespace Xcts::Solutions
|
