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 <optional>
9 :
10 : #include "DataStructures/CachedTempBuffer.hpp"
11 : #include "DataStructures/DataBox/Prefixes.hpp"
12 : #include "DataStructures/TempBuffer.hpp"
13 : #include "DataStructures/Tensor/EagerMath/Magnitude.hpp"
14 : #include "DataStructures/Tensor/Tensor.hpp"
15 : #include "Elliptic/Systems/Xcts/Tags.hpp"
16 : #include "NumericalAlgorithms/LinearOperators/PartialDerivatives.hpp"
17 : #include "Options/Auto.hpp"
18 : #include "Options/Context.hpp"
19 : #include "Options/ParseError.hpp"
20 : #include "Options/String.hpp"
21 : #include "PointwiseFunctions/AnalyticData/Xcts/CommonVariables.hpp"
22 : #include "PointwiseFunctions/AnalyticSolutions/Xcts/Flatness.hpp"
23 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
24 : #include "PointwiseFunctions/GeneralRelativity/Tags/Conformal.hpp"
25 : #include "PointwiseFunctions/InitialDataUtilities/Background.hpp"
26 : #include "PointwiseFunctions/InitialDataUtilities/InitialGuess.hpp"
27 : #include "Utilities/CallWithDynamicType.hpp"
28 : #include "Utilities/Requires.hpp"
29 : #include "Utilities/Serialization/CharmPupable.hpp"
30 : #include "Utilities/Serialization/PupStlCpp17.hpp"
31 : #include "Utilities/TMPL.hpp"
32 : #include "Utilities/TaggedTuple.hpp"
33 :
34 : /// \cond
35 : namespace PUP {
36 : class er;
37 : } // namespace PUP
38 : /// \endcond
39 :
40 : namespace Xcts::AnalyticData {
41 :
42 : namespace detail {
43 :
44 : template <typename DataType>
45 : using BinaryVariablesCache = cached_temp_buffer_from_typelist<tmpl::append<
46 : common_tags<DataType>,
47 : tmpl::list<
48 : ::Tags::deriv<Tags::ShiftBackground<DataType, 3, Frame::Inertial>,
49 : tmpl::size_t<3>, Frame::Inertial>,
50 : gr::Tags::Conformal<gr::Tags::EnergyDensity<DataType>, 0>,
51 : gr::Tags::Conformal<gr::Tags::StressTrace<DataType>, 0>,
52 : gr::Tags::Conformal<gr::Tags::MomentumDensity<DataType, 3>, 0>,
53 : // For initial guesses
54 : Tags::ConformalFactor<DataType>,
55 : Tags::LapseTimesConformalFactor<DataType>,
56 : Tags::ShiftExcess<DataType, 3, Frame::Inertial>>,
57 : hydro_tags<DataType>>>;
58 :
59 : template <typename DataType>
60 : struct BinaryVariables
61 : : CommonVariables<DataType, BinaryVariablesCache<DataType>> {
62 : static constexpr size_t Dim = 3;
63 : using Cache = BinaryVariablesCache<DataType>;
64 : using Base = CommonVariables<DataType, BinaryVariablesCache<DataType>>;
65 : using Base::operator();
66 :
67 : using superposed_tags = tmpl::append<
68 : tmpl::list<
69 : Tags::ConformalMetric<DataType, Dim, Frame::Inertial>,
70 : ::Tags::deriv<Tags::ConformalMetric<DataType, Dim, Frame::Inertial>,
71 : tmpl::size_t<Dim>, Frame::Inertial>,
72 : gr::Tags::TraceExtrinsicCurvature<DataType>,
73 : ::Tags::dt<gr::Tags::TraceExtrinsicCurvature<DataType>>,
74 : gr::Tags::Conformal<gr::Tags::EnergyDensity<DataType>, 0>,
75 : gr::Tags::Conformal<gr::Tags::StressTrace<DataType>, 0>,
76 : gr::Tags::Conformal<gr::Tags::MomentumDensity<DataType, Dim>, 0>,
77 : Tags::ConformalFactor<DataType>,
78 : Tags::LapseTimesConformalFactor<DataType>,
79 : Tags::ShiftExcess<DataType, Dim, Frame::Inertial>>,
80 : hydro_tags<DataType>>;
81 :
82 : BinaryVariables(
83 : std::optional<std::reference_wrapper<const Mesh<Dim>>> local_mesh,
84 : std::optional<std::reference_wrapper<const InverseJacobian<
85 : DataType, Dim, Frame::ElementLogical, Frame::Inertial>>>
86 : local_inv_jacobian,
87 : const tnsr::I<DataVector, Dim>& local_x,
88 : const double local_angular_velocity, const double local_expansion,
89 : const std::array<double, 3> local_linear_velocity,
90 : std::optional<std::array<double, 2>> local_falloff_widths,
91 : std::array<tnsr::I<DataVector, Dim>, 2> local_x_isolated,
92 : std::array<DataVector, 2> local_windows,
93 : tuples::tagged_tuple_from_typelist<superposed_tags> local_flat_vars,
94 : std::array<tuples::tagged_tuple_from_typelist<superposed_tags>, 2>
95 : local_isolated_vars)
96 : : Base(std::move(local_mesh), std::move(local_inv_jacobian)),
97 : x(local_x),
98 : angular_velocity(local_angular_velocity),
99 : expansion(local_expansion),
100 : linear_velocity(local_linear_velocity),
101 : falloff_widths(std::move(local_falloff_widths)),
102 : x_isolated(std::move(local_x_isolated)),
103 : windows(std::move(local_windows)),
104 : flat_vars(std::move(local_flat_vars)),
105 : isolated_vars(std::move(local_isolated_vars)) {}
106 :
107 : const tnsr::I<DataVector, Dim>& x;
108 : const double angular_velocity;
109 : const double expansion;
110 : const std::array<double, 3> linear_velocity;
111 : const std::optional<std::array<double, 2>> falloff_widths;
112 : const std::array<tnsr::I<DataVector, Dim>, 2> x_isolated;
113 : const std::array<DataVector, 2> windows;
114 : const tuples::tagged_tuple_from_typelist<superposed_tags> flat_vars;
115 : const std::array<tuples::tagged_tuple_from_typelist<superposed_tags>, 2>
116 : isolated_vars;
117 :
118 : template <bool ApplyWindow = true, typename Tag,
119 : Requires<tmpl::list_contains_v<superposed_tags, Tag>> = nullptr>
120 : void superposition(gsl::not_null<typename Tag::type*> superposed_var,
121 : gsl::not_null<Cache*> /*cache*/, Tag /*meta*/) const {
122 : for (size_t i = 0; i < superposed_var->size(); ++i) {
123 : if constexpr (ApplyWindow) {
124 : (*superposed_var)[i] =
125 : get<Tag>(flat_vars)[i] +
126 : windows[0] *
127 : (get<Tag>(isolated_vars[0])[i] - get<Tag>(flat_vars)[i]) +
128 : windows[1] *
129 : (get<Tag>(isolated_vars[1])[i] - get<Tag>(flat_vars)[i]);
130 : } else {
131 : (*superposed_var)[i] = get<Tag>(isolated_vars[0])[i] +
132 : get<Tag>(isolated_vars[1])[i] -
133 : get<Tag>(flat_vars)[i];
134 : }
135 : }
136 : }
137 :
138 : void operator()(
139 : const gsl::not_null<tnsr::ii<DataType, Dim>*> conformal_metric,
140 : const gsl::not_null<Cache*> cache,
141 : Tags::ConformalMetric<DataType, Dim, Frame::Inertial> meta)
142 : const override {
143 : superposition(conformal_metric, cache, meta);
144 : }
145 : void operator()(
146 : const gsl::not_null<tnsr::ijj<DataType, Dim>*> deriv_conformal_metric,
147 : const gsl::not_null<Cache*> cache,
148 : ::Tags::deriv<Tags::ConformalMetric<DataType, Dim, Frame::Inertial>,
149 : tmpl::size_t<Dim>, Frame::Inertial>
150 : meta) const override {
151 : superposition(deriv_conformal_metric, cache, meta);
152 : add_deriv_of_window_function(deriv_conformal_metric);
153 : }
154 : void operator()(
155 : const gsl::not_null<Scalar<DataType>*> extrinsic_curvature_trace,
156 : const gsl::not_null<Cache*> cache,
157 : gr::Tags::TraceExtrinsicCurvature<DataType> meta) const override {
158 : superposition(extrinsic_curvature_trace, cache, meta);
159 : }
160 : void operator()(
161 : const gsl::not_null<Scalar<DataType>*> dt_extrinsic_curvature_trace,
162 : const gsl::not_null<Cache*> cache,
163 : ::Tags::dt<gr::Tags::TraceExtrinsicCurvature<DataType>> meta)
164 : const override {
165 : superposition(dt_extrinsic_curvature_trace, cache, meta);
166 : }
167 : void operator()(
168 : gsl::not_null<tnsr::I<DataType, Dim>*> shift_background,
169 : gsl::not_null<Cache*> cache,
170 : Tags::ShiftBackground<DataType, Dim, Frame::Inertial> /*meta*/)
171 : const override;
172 : void operator()(
173 : gsl::not_null<tnsr::iJ<DataType, Dim>*> deriv_shift_background,
174 : gsl::not_null<Cache*> cache,
175 : ::Tags::deriv<Tags::ShiftBackground<DataType, Dim, Frame::Inertial>,
176 : tmpl::size_t<Dim>, Frame::Inertial> /*meta*/) const;
177 : void operator()(gsl::not_null<tnsr::II<DataType, Dim, Frame::Inertial>*>
178 : longitudinal_shift_background_minus_dt_conformal_metric,
179 : gsl::not_null<Cache*> cache,
180 : Tags::LongitudinalShiftBackgroundMinusDtConformalMetric<
181 : DataType, Dim, Frame::Inertial> /*meta*/) const override;
182 : void operator()(
183 : const gsl::not_null<Scalar<DataType>*> conformal_energy_density,
184 : const gsl::not_null<Cache*> cache,
185 : gr::Tags::Conformal<gr::Tags::EnergyDensity<DataType>, 0> meta) const {
186 : superposition<false>(conformal_energy_density, cache, meta);
187 : }
188 : void operator()(
189 : const gsl::not_null<Scalar<DataType>*> conformal_stress_trace,
190 : const gsl::not_null<Cache*> cache,
191 : gr::Tags::Conformal<gr::Tags::StressTrace<DataType>, 0> meta) const {
192 : superposition<false>(conformal_stress_trace, cache, meta);
193 : }
194 : void operator()(
195 : const gsl::not_null<tnsr::I<DataType, Dim>*> conformal_momentum_density,
196 : const gsl::not_null<Cache*> cache,
197 : gr::Tags::Conformal<gr::Tags::MomentumDensity<DataType, Dim>, 0> meta)
198 : const {
199 : superposition<false>(conformal_momentum_density, cache, meta);
200 : }
201 : void operator()(const gsl::not_null<Scalar<DataType>*> conformal_factor,
202 : const gsl::not_null<Cache*> cache,
203 : Tags::ConformalFactor<DataType> meta) const {
204 : superposition(conformal_factor, cache, meta);
205 : }
206 : void operator()(
207 : const gsl::not_null<Scalar<DataType>*> lapse_times_conformal_factor,
208 : const gsl::not_null<Cache*> cache,
209 : Tags::LapseTimesConformalFactor<DataType> meta) const {
210 : superposition(lapse_times_conformal_factor, cache, meta);
211 : }
212 : void operator()(
213 : const gsl::not_null<tnsr::I<DataType, Dim>*> shift_excess,
214 : const gsl::not_null<Cache*> cache,
215 : Tags::ShiftExcess<DataType, Dim, Frame::Inertial> meta) const {
216 : superposition(shift_excess, cache, meta);
217 : }
218 : void operator()(const gsl::not_null<Scalar<DataType>*> rest_mass_density,
219 : const gsl::not_null<Cache*> cache,
220 : hydro::Tags::RestMassDensity<DataType> meta) const {
221 : superposition<false>(rest_mass_density, cache, meta);
222 : }
223 : void operator()(const gsl::not_null<Scalar<DataType>*> specific_enthalpy,
224 : const gsl::not_null<Cache*> cache,
225 : hydro::Tags::SpecificEnthalpy<DataType> meta) const {
226 : superposition<false>(specific_enthalpy, cache, meta);
227 : }
228 : void operator()(const gsl::not_null<Scalar<DataType>*> pressure,
229 : const gsl::not_null<Cache*> cache,
230 : hydro::Tags::Pressure<DataType> meta) const {
231 : superposition<false>(pressure, cache, meta);
232 : }
233 : void operator()(const gsl::not_null<tnsr::I<DataType, 3>*> spatial_velocity,
234 : const gsl::not_null<Cache*> cache,
235 : hydro::Tags::SpatialVelocity<DataType, 3> meta) const {
236 : superposition<false>(spatial_velocity, cache, meta);
237 : }
238 : void operator()(const gsl::not_null<Scalar<DataType>*> lorentz_factor,
239 : const gsl::not_null<Cache*> cache,
240 : hydro::Tags::LorentzFactor<DataType> meta) const {
241 : superposition<false>(lorentz_factor, cache, meta);
242 : }
243 : void operator()(const gsl::not_null<tnsr::I<DataType, 3>*> magnetic_field,
244 : const gsl::not_null<Cache*> cache,
245 : hydro::Tags::MagneticField<DataType, 3> meta) const {
246 : superposition<false>(magnetic_field, cache, meta);
247 : }
248 :
249 : private:
250 : void add_deriv_of_window_function(
251 : gsl::not_null<tnsr::ijj<DataType, Dim>*> deriv_conformal_metric) const;
252 : };
253 : } // namespace detail
254 :
255 : /*!
256 : * \brief Binary compact-object data in general relativity, constructed from
257 : * superpositions of two isolated objects.
258 : *
259 : * This class implements background data for the XCTS equations describing two
260 : * objects in a quasi-equilibrium orbit, i.e. with \f$\bar{u}=0\f$ and
261 : * \f$\partial_t K=0\f$. Both objects can be chosen from the list of
262 : * `IsolatedObjectRegistrars`, e.g. they can be black-hole or neutron-star
263 : * solutions in different coordinates. Most quantities are constructed by
264 : * superposing the two isolated solutions (see e.g. Eq. (8-9) in
265 : * \cite Varma2018sqd or Eq. (45-46) in \cite Lovelace2008tw):
266 : *
267 : * \f{align}
268 : * \bar{\gamma}_{ij} &= f_{ij} + \sum_{\alpha=1}^2
269 : * e^{-r_\alpha^2 / w_\alpha^2}\left(\gamma^\alpha_{ij} - f_{ij}\right) \\
270 : * K &= \sum_{\alpha=1}^2 e^{-r_\alpha^2 / w_\alpha^2}K^\alpha
271 : * \f}
272 : *
273 : * where \f$\gamma^\alpha_{ij}\f$ and \f$K^\alpha\f$ denote the spatial metric
274 : * and extrinsic-curvature trace of the two individual solutions, \f$r_\alpha\f$
275 : * is the Euclidean coordinate-distance from the center of each object and
276 : * \f$w_\alpha\f$ are parameters describing the falloff widths of Gaussian
277 : * window functions. The window functions
278 : * facilitate that the influence of either of the two objects
279 : * at the position of the other is strongly damped, and they also avoid
280 : * logarithmic scaling of the solution at large distances where we would
281 : * typically employ an inverse-radial coordinate map and asymptotically-flat
282 : * boundary conditions. The falloff-widths are chosen in terms of the Newtonian
283 : * Lagrange points of the two objects in \cite Varma2018sqd and
284 : * \cite Lovelace2008tw, and they are input parameters in this implementation.
285 : * The falloff can be disabled by passing `std::nullopt` to the constructor, or
286 : * `None` in the input file.
287 : *
288 : * \par Matter sources
289 : * Matter sources are superposed without the window functions. The analytic
290 : * matter sources are of
291 : * limited use anyway, because in a binary setting they don't take the
292 : * gravitational influence of the other body into account. Therefore, the matter
293 : * sources should typically be solved-for alongside the gravity sector to impose
294 : * conditions such as hydrostatic equilibrium. For scenarios where we just want
295 : * to superpose the isolated matter solutions and compute the resulting gravity,
296 : * the matter sources are simply added.
297 : *
298 : * \par Orbital motion
299 : * The remaining quantities that this class implements relate to the orbital
300 : * motion of the two objects. To obtain initial data in "co-rotating"
301 : * coordinates where the two objects are initially at rest we prescribe the
302 : * background shift
303 : *
304 : * \f{equation} \beta^i_\mathrm{background} = (-\Omega y, \Omega x, 0) +
305 : * \dot{a}_0 x^i + v^i_0 \f}
306 : *
307 : * where \f$\Omega\f$ is the angular-velocity parameter and \f$\dot{a}_0\f$
308 : * is an expansion parameter. Both control the eccentricity of the orbit.
309 : * The parameter \f$v^i_0\f$ is a constant velocity that can be used to
310 : * control the linear momentum of the system (see Eq. (28) in
311 : * \cite Ossokine2015yla).
312 : */
313 : template <typename IsolatedObjectBase, typename IsolatedObjectClasses>
314 1 : class Binary : public elliptic::analytic_data::Background,
315 : public elliptic::analytic_data::InitialGuess {
316 : public:
317 0 : struct XCoords {
318 0 : static constexpr Options::String help =
319 : "The coordinates on the x-axis where the two objects are placed";
320 0 : using type = std::array<double, 2>;
321 : };
322 0 : struct ObjectLeft {
323 0 : static constexpr Options::String help =
324 : "The object placed on the negative x-axis";
325 0 : using type = std::unique_ptr<IsolatedObjectBase>;
326 : };
327 0 : struct ObjectRight {
328 0 : static constexpr Options::String help =
329 : "The object placed on the positive x-axis";
330 0 : using type = std::unique_ptr<IsolatedObjectBase>;
331 : };
332 0 : struct AngularVelocity {
333 0 : static constexpr Options::String help =
334 : "Orbital angular velocity 'Omega0' about the z-axis. Added to the "
335 : "background shift as a term 'Omega0 x r'.";
336 0 : using type = double;
337 : };
338 0 : struct Expansion {
339 0 : static constexpr Options::String help =
340 : "The expansion parameter 'adot0', which is a radial velocity over "
341 : "radius. Added to the background shift as a term 'adot0 r^i'";
342 0 : using type = double;
343 : };
344 0 : struct LinearVelocity {
345 0 : static constexpr Options::String help =
346 : "Constant velocity 'v0' added to the background shift to control the "
347 : "linear momentum of the system.";
348 0 : using type = std::array<double, 3>;
349 : };
350 0 : struct FalloffWidths {
351 0 : static constexpr Options::String help =
352 : "The widths for the window functions around the two objects, or 'None' "
353 : "to disable the Gaussian falloff.";
354 0 : using type = Options::Auto<std::array<double, 2>, Options::AutoLabel::None>;
355 : };
356 0 : using options = tmpl::list<XCoords, ObjectLeft, ObjectRight, AngularVelocity,
357 : Expansion, LinearVelocity, FalloffWidths>;
358 0 : static constexpr Options::String help =
359 : "Binary compact-object data in general relativity, constructed from "
360 : "superpositions of two isolated objects.";
361 :
362 0 : Binary() = default;
363 0 : Binary(const Binary&) = delete;
364 0 : Binary& operator=(const Binary&) = delete;
365 0 : Binary(Binary&&) = default;
366 0 : Binary& operator=(Binary&&) = default;
367 0 : ~Binary() = default;
368 :
369 0 : Binary(const std::array<double, 2> xcoords,
370 : std::unique_ptr<IsolatedObjectBase> object_left,
371 : std::unique_ptr<IsolatedObjectBase> object_right,
372 : const double angular_velocity, const double expansion,
373 : const std::array<double, 3> linear_velocity,
374 : const std::optional<std::array<double, 2>> falloff_widths,
375 : const Options::Context& context = {})
376 : : xcoords_(xcoords),
377 : superposed_objects_({std::move(object_left), std::move(object_right)}),
378 : angular_velocity_(angular_velocity),
379 : expansion_(expansion),
380 : linear_velocity_(linear_velocity),
381 : falloff_widths_(falloff_widths) {
382 : if (xcoords_[0] >= xcoords_[1]) {
383 : PARSE_ERROR(context, "Specify 'XCoords' ascending from left to right.");
384 : }
385 : }
386 :
387 0 : explicit Binary(CkMigrateMessage* m)
388 : : elliptic::analytic_data::Background(m),
389 : elliptic::analytic_data::InitialGuess(m) {}
390 : using PUP::able::register_constructor;
391 0 : WRAPPED_PUPable_decl_template(Binary);
392 :
393 : template <typename DataType, typename... RequestedTags>
394 0 : tuples::TaggedTuple<RequestedTags...> variables(
395 : const tnsr::I<DataType, 3, Frame::Inertial>& x,
396 : tmpl::list<RequestedTags...> /*meta*/) const {
397 : return variables_impl<DataType>(x, std::nullopt, std::nullopt,
398 : tmpl::list<RequestedTags...>{});
399 : }
400 : template <typename... RequestedTags>
401 0 : tuples::TaggedTuple<RequestedTags...> variables(
402 : const tnsr::I<DataVector, 3, Frame::Inertial>& x, const Mesh<3>& mesh,
403 : const InverseJacobian<DataVector, 3, Frame::ElementLogical,
404 : Frame::Inertial>& inv_jacobian,
405 : tmpl::list<RequestedTags...> /*meta*/) const {
406 : return variables_impl<DataVector>(x, mesh, inv_jacobian,
407 : tmpl::list<RequestedTags...>{});
408 : }
409 :
410 : // NOLINTNEXTLINE
411 0 : void pup(PUP::er& p) override {
412 : elliptic::analytic_data::Background::pup(p);
413 : elliptic::analytic_data::InitialGuess::pup(p);
414 : p | xcoords_;
415 : p | superposed_objects_;
416 : p | angular_velocity_;
417 : p | expansion_;
418 : p | linear_velocity_;
419 : p | falloff_widths_;
420 : }
421 :
422 : /// Coordinates of the objects, ascending left to right
423 1 : const std::array<double, 2>& x_coords() const { return xcoords_; }
424 : /// The two objects. First entry is the left object, second entry is the right
425 : /// object.
426 1 : const std::array<std::unique_ptr<IsolatedObjectBase>, 2>& superposed_objects()
427 : const {
428 : return superposed_objects_;
429 : }
430 0 : double angular_velocity() const { return angular_velocity_; }
431 0 : double expansion() const { return expansion_; }
432 0 : const std::array<double, 3>& linear_velocity() const {
433 : return linear_velocity_;
434 : }
435 0 : const std::optional<std::array<double, 2>>& falloff_widths() const {
436 : return falloff_widths_;
437 : }
438 :
439 : private:
440 0 : std::array<double, 2> xcoords_{};
441 0 : std::array<std::unique_ptr<IsolatedObjectBase>, 2> superposed_objects_{};
442 0 : Xcts::Solutions::Flatness flatness_{};
443 0 : double angular_velocity_ = std::numeric_limits<double>::signaling_NaN();
444 0 : double expansion_ = std::numeric_limits<double>::signaling_NaN();
445 0 : std::array<double, 3> linear_velocity_{};
446 0 : std::optional<std::array<double, 2>> falloff_widths_{};
447 :
448 : template <typename DataType, typename... RequestedTags>
449 0 : tuples::TaggedTuple<RequestedTags...> variables_impl(
450 : const tnsr::I<DataType, 3, Frame::Inertial>& x,
451 : std::optional<std::reference_wrapper<const Mesh<3>>> mesh,
452 : std::optional<std::reference_wrapper<const InverseJacobian<
453 : DataType, 3, Frame::ElementLogical, Frame::Inertial>>>
454 : inv_jacobian,
455 : tmpl::list<RequestedTags...> /*meta*/) const {
456 : std::array<tnsr::I<DataVector, 3>, 2> x_isolated{{x, x}};
457 : std::array<DataVector, 2> euclidean_distance{};
458 : std::array<DataVector, 2> windows{};
459 : // Possible optimization: Only retrieve those superposed tags from the
460 : // isolated solutions that are actually needed. This needs some dependency
461 : // logic, because some of the non-superposed tags depend on superposed tags.
462 : using VarsComputer = detail::BinaryVariables<DataType>;
463 : using requested_superposed_tags = typename VarsComputer::superposed_tags;
464 : std::array<tuples::tagged_tuple_from_typelist<requested_superposed_tags>, 2>
465 : isolated_vars;
466 : for (size_t i = 0; i < 2; ++i) {
467 : get<0>(gsl::at(x_isolated, i)) -= gsl::at(xcoords_, i);
468 : gsl::at(euclidean_distance, i) = get(magnitude(gsl::at(x_isolated, i)));
469 : if (falloff_widths_.has_value()) {
470 : gsl::at(windows, i) = exp(-square(gsl::at(euclidean_distance, i)) /
471 : square(gsl::at(*falloff_widths_, i)));
472 : } else {
473 : gsl::at(windows, i) = make_with_value<DataVector>(x, 1.);
474 : }
475 : gsl::at(isolated_vars, i) = get_isolated_vars<requested_superposed_tags>(
476 : *gsl::at(superposed_objects_, i), gsl::at(x_isolated, i));
477 : }
478 : auto flat_vars = flatness_.variables(x, requested_superposed_tags{});
479 : typename VarsComputer::Cache cache{get_size(*x.begin())};
480 : const VarsComputer computer{std::move(mesh),
481 : std::move(inv_jacobian),
482 : x,
483 : angular_velocity_,
484 : expansion_,
485 : linear_velocity_,
486 : falloff_widths_,
487 : std::move(x_isolated),
488 : std::move(windows),
489 : std::move(flat_vars),
490 : std::move(isolated_vars)};
491 : return {cache.get_var(computer, RequestedTags{})...};
492 : }
493 :
494 : template <typename TagsList, typename... Args>
495 0 : tuples::tagged_tuple_from_typelist<TagsList> get_isolated_vars(
496 : const IsolatedObjectBase& isolated_object, const Args&... args) const {
497 : return call_with_dynamic_type<tuples::tagged_tuple_from_typelist<TagsList>,
498 : IsolatedObjectClasses>(
499 : &isolated_object, [&args...](const auto* const derived) {
500 : return derived->variables(args..., TagsList{});
501 : });
502 : }
503 : };
504 :
505 : /// \cond
506 : template <typename IsolatedObjectBase, typename IsolatedObjectClasses>
507 : PUP::able::PUP_ID Binary<IsolatedObjectBase, IsolatedObjectClasses>::my_PUP_ID =
508 : 0; // NOLINT
509 : /// \endcond
510 :
511 : } // namespace Xcts::AnalyticData
|