Schwarzschild.hpp
1 // 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"
11 #include "DataStructures/DataVector.hpp"
13 #include "Elliptic/Systems/Xcts/Tags.hpp"
15 #include "Options/Options.hpp"
16 #include "Options/ParseOptions.hpp"
18 #include "PointwiseFunctions/AnalyticSolutions/Xcts/AnalyticSolution.hpp"
19 #include "PointwiseFunctions/AnalyticSolutions/Xcts/CommonVariables.hpp"
20 #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
21 #include "Utilities/ContainerHelpers.hpp"
22 #include "Utilities/Gsl.hpp"
23 #include "Utilities/TMPL.hpp"
24 #include "Utilities/TaggedTuple.hpp"
25 
26 /// \cond
27 namespace PUP {
28 class er;
29 } // namespace PUP
30 /// \endcond
31 
32 namespace Xcts::Solutions {
33 
34 /// Various coordinate systems in which to express the Schwarzschild solution
36  /*!
37  * \brief Isotropic Schwarzschild coordinates
38  *
39  * These arise from the canonical Schwarzschild coordinates by the radial
40  * transformation
41  *
42  * \f{equation}
43  * r = \bar{r}\left(1+\frac{M}{2\bar{r}}\right)^2
44  * \f}
45  *
46  * (Eq. (1.61) in \cite BaumgarteShapiro) where \f$r\f$ is the canonical
47  * Schwarzschild radius, also referred to as "areal" radius because it is
48  * defined such that spheres with constant \f$r\f$ have the area \f$4\pi
49  * r^2\f$, and \f$\bar{r}\f$ is the "isotropic" radius. In the isotropic
50  * radius the Schwarzschild spatial metric is conformally flat:
51  *
52  * \f{equation}
53  * \gamma_{ij}=\psi^4\eta_{ij} \quad \text{with conformal factor} \quad
54  * \psi=1+\frac{M}{2\bar{r}}
55  * \f}
56  *
57  * (Table 2.1 in \cite BaumgarteShapiro). Its lapse transforms to
58  *
59  * \f{equation}
60  * \alpha=\frac{1-M/(2\bar{r})}{1+M/(2\bar{r})}
61  * \f}
62  *
63  * and the shift vanishes (\f$\beta^i=0\f$) as it does in areal Schwarzschild
64  * coordinates. The solution also remains maximally sliced, i.e. \f$K=0\f$.
65  *
66  * The Schwarzschild horizon in these coordinates is at
67  * \f$\bar{r}=\frac{M}{2}\f$ due to the radial transformation from \f$r=2M\f$.
68  */
69  Isotropic,
70 };
71 
72 std::ostream& operator<<(std::ostream& os,
73  SchwarzschildCoordinates coords) noexcept;
74 
75 } // namespace Xcts::Solutions
76 
77 template <>
78 struct Options::create_from_yaml<Xcts::Solutions::SchwarzschildCoordinates> {
79  template <typename Metavariables>
81  const Options::Option& options) {
82  return create<void>(options);
83  }
84 };
85 
86 template <>
89  void>(const Options::Option& options);
90 
91 namespace Xcts::Solutions {
92 
93 namespace detail {
94 
95 struct SchwarzschildImpl {
96  struct Mass {
97  using type = double;
98  static constexpr Options::String help = "Mass parameter M";
99  };
100 
101  struct CoordinateSystem {
102  static std::string name() noexcept { return "Coordinates"; }
103  using type = SchwarzschildCoordinates;
104  static constexpr Options::String help =
105  "The coordinate system used to describe the solution";
106  };
107 
108  using options = tmpl::list<Mass, CoordinateSystem>;
109  static constexpr Options::String help{
110  "Schwarzschild spacetime in general relativity"};
111 
112  SchwarzschildImpl() = default;
113  SchwarzschildImpl(const SchwarzschildImpl&) noexcept = default;
114  SchwarzschildImpl& operator=(const SchwarzschildImpl&) noexcept = default;
115  SchwarzschildImpl(SchwarzschildImpl&&) noexcept = default;
116  SchwarzschildImpl& operator=(SchwarzschildImpl&&) noexcept = default;
117  ~SchwarzschildImpl() noexcept = default;
118 
119  explicit SchwarzschildImpl(
120  double mass, SchwarzschildCoordinates coordinate_system) noexcept;
121 
122  /// The mass parameter \f$M\f$.
123  double mass() const noexcept;
124 
125  SchwarzschildCoordinates coordinate_system() const noexcept;
126 
127  /// The radius of the Schwarzschild horizon in the given coordinates.
128  double radius_at_horizon() const noexcept;
129 
130  // NOLINTNEXTLINE(google-runtime-references)
131  void pup(PUP::er& p) noexcept;
132 
133  protected:
135  SchwarzschildCoordinates coordinate_system_{};
136 };
137 
138 bool operator==(const SchwarzschildImpl& lhs,
139  const SchwarzschildImpl& rhs) noexcept;
140 
141 bool operator!=(const SchwarzschildImpl& lhs,
142  const SchwarzschildImpl& rhs) noexcept;
143 
144 template <typename DataType>
145 struct SchwarzschildVariables;
146 
147 template <typename DataType>
148 using SchwarzschildVariablesCache = cached_temp_buffer_from_typelist<
149  SchwarzschildVariables<DataType>,
151  common_tags<DataType>,
152  Tags::ConformalMetric<DataType, 3, Frame::Inertial>,
153  Tags::InverseConformalMetric<DataType, 3, Frame::Inertial>,
155  tmpl::size_t<3>, Frame::Inertial>,
158  Tags::ConformalFactor<DataType>,
161  gr::Tags::Lapse<DataType>, Tags::LapseTimesConformalFactor<DataType>,
163  tmpl::size_t<3>, Frame::Inertial>,
164  Tags::ShiftBackground<DataType, 3, Frame::Inertial>,
165  Tags::LongitudinalShiftBackgroundMinusDtConformalMetric<
166  DataType, 3, Frame::Inertial>,
167  Tags::ShiftExcess<DataType, 3, Frame::Inertial>,
168  Tags::ShiftStrain<DataType, 3, Frame::Inertial>,
171 
172 template <typename DataType>
173 struct SchwarzschildVariables
174  : CommonVariables<DataType, SchwarzschildVariablesCache<DataType>> {
175  static constexpr size_t Dim = 3;
176  using Cache = SchwarzschildVariablesCache<DataType>;
177  using CommonVariables<DataType,
178  SchwarzschildVariablesCache<DataType>>::operator();
179 
180  const tnsr::I<DataType, 3>& x;
181  double mass;
182  SchwarzschildCoordinates coordinate_system;
183 
184  void operator()(gsl::not_null<tnsr::ii<DataType, 3>*> conformal_metric,
186  Tags::ConformalMetric<DataType, 3, Frame::Inertial> /*meta*/)
187  const noexcept;
188  void operator()(
189  gsl::not_null<tnsr::II<DataType, 3>*> inv_conformal_metric,
191  Tags::InverseConformalMetric<DataType, 3, Frame::Inertial> /*meta*/)
192  const noexcept;
193  void operator()(
194  gsl::not_null<tnsr::ijj<DataType, 3>*> deriv_conformal_metric,
196  ::Tags::deriv<Tags::ConformalMetric<DataType, 3, Frame::Inertial>,
197  tmpl::size_t<3>, Frame::Inertial> /*meta*/) const noexcept;
198  void operator()(
199  gsl::not_null<Scalar<DataType>*> trace_extrinsic_curvature,
200  gsl::not_null<Cache*> cache,
201  gr::Tags::TraceExtrinsicCurvature<DataType> /*meta*/) const noexcept;
202  void operator()(
203  gsl::not_null<tnsr::i<DataType, 3>*> trace_extrinsic_curvature_gradient,
204  gsl::not_null<Cache*> cache,
205  ::Tags::deriv<gr::Tags::TraceExtrinsicCurvature<DataType>,
206  tmpl::size_t<3>, Frame::Inertial> /*meta*/) const noexcept;
207  void operator()(
208  gsl::not_null<Scalar<DataType>*> dt_trace_extrinsic_curvature,
209  gsl::not_null<Cache*> cache,
210  ::Tags::dt<gr::Tags::TraceExtrinsicCurvature<DataType>> /*meta*/)
211  const noexcept;
212  void operator()(gsl::not_null<Scalar<DataType>*> conformal_factor,
213  gsl::not_null<Cache*> cache,
214  Tags::ConformalFactor<DataType> /*meta*/) const noexcept;
215  void operator()(
216  gsl::not_null<tnsr::i<DataType, 3>*> conformal_factor_gradient,
217  gsl::not_null<Cache*> cache,
218  ::Tags::deriv<Xcts::Tags::ConformalFactor<DataType>, tmpl::size_t<3>,
219  Frame::Inertial> /*meta*/) const noexcept;
220  void operator()(gsl::not_null<Scalar<DataType>*> lapse,
221  gsl::not_null<Cache*> cache,
222  gr::Tags::Lapse<DataType> /*meta*/) const noexcept;
223  void operator()(
224  gsl::not_null<Scalar<DataType>*> lapse_times_conformal_factor,
225  gsl::not_null<Cache*> cache,
226  Tags::LapseTimesConformalFactor<DataType> /*meta*/) const noexcept;
227  void operator()(
228  gsl::not_null<tnsr::i<DataType, 3>*>
229  lapse_times_conformal_factor_gradient,
230  gsl::not_null<Cache*> cache,
231  ::Tags::deriv<Tags::LapseTimesConformalFactor<DataType>, tmpl::size_t<3>,
232  Frame::Inertial> /*meta*/) const noexcept;
233  void operator()(gsl::not_null<tnsr::I<DataType, 3>*> shift_background,
234  gsl::not_null<Cache*> cache,
235  Tags::ShiftBackground<DataType, 3, Frame::Inertial> /*meta*/)
236  const noexcept;
237  void operator()(gsl::not_null<tnsr::II<DataType, 3, Frame::Inertial>*>
238  longitudinal_shift_background_minus_dt_conformal_metric,
239  gsl::not_null<Cache*> cache,
240  Tags::LongitudinalShiftBackgroundMinusDtConformalMetric<
241  DataType, 3, Frame::Inertial> /*meta*/) const noexcept;
242  void operator()(
243  gsl::not_null<tnsr::I<DataType, 3>*> shift_excess,
244  gsl::not_null<Cache*> cache,
245  Tags::ShiftExcess<DataType, 3, Frame::Inertial> /*meta*/) const noexcept;
246  void operator()(
247  gsl::not_null<tnsr::ii<DataType, 3>*> shift_strain,
248  gsl::not_null<Cache*> cache,
249  Tags::ShiftStrain<DataType, 3, Frame::Inertial> /*meta*/) const noexcept;
250  void operator()(gsl::not_null<Scalar<DataType>*> energy_density,
251  gsl::not_null<Cache*> cache,
252  gr::Tags::EnergyDensity<DataType> /*meta*/) const noexcept;
253  void operator()(gsl::not_null<Scalar<DataType>*> stress_trace,
254  gsl::not_null<Cache*> cache,
255  gr::Tags::StressTrace<DataType> /*meta*/) const noexcept;
256  void operator()(gsl::not_null<tnsr::I<DataType, 3>*> momentum_density,
257  gsl::not_null<Cache*> cache,
258  gr::Tags::MomentumDensity<3, Frame::Inertial,
259  DataType> /*meta*/) const noexcept;
260 };
261 
262 } // namespace detail
263 
264 // The following implements the registration and factory-creation mechanism
265 
266 /// \cond
267 template <typename Registrars>
268 struct Schwarzschild;
269 
270 namespace Registrars {
271 struct Schwarzschild {
272  template <typename Registrars>
273  using f = Solutions::Schwarzschild<Registrars>;
274 };
275 } // namespace Registrars
276 /// \endcond
277 
278 /*!
279  * \brief Schwarzschild spacetime in general relativity
280  *
281  * This class implements the Schwarzschild solution with mass parameter
282  * \f$M\f$ in various coordinate systems. See the entries of the
283  * `Xcts::Solutions::SchwarzschildCoordinates` enum for the available coordinate
284  * systems and for the solution variables in the respective coordinates.
285  */
286 template <typename Registrars =
287  tmpl::list<Solutions::Registrars::Schwarzschild>>
288 class Schwarzschild : public AnalyticSolution<Registrars>,
289  public detail::SchwarzschildImpl {
290  private:
292 
293  public:
294  Schwarzschild() = default;
295  Schwarzschild(const Schwarzschild&) noexcept = default;
296  Schwarzschild& operator=(const Schwarzschild&) noexcept = default;
297  Schwarzschild(Schwarzschild&&) noexcept = default;
298  Schwarzschild& operator=(Schwarzschild&&) noexcept = default;
299  ~Schwarzschild() noexcept = default;
300 
301  using SchwarzschildImpl::SchwarzschildImpl;
302 
303  /// \cond
304  explicit Schwarzschild(CkMigrateMessage* m) noexcept : Base(m) {}
305  using PUP::able::register_constructor;
307  /// \endcond
308 
309  template <typename DataType, typename... RequestedTags>
310  tuples::TaggedTuple<RequestedTags...> variables(
311  const tnsr::I<DataType, 3, Frame::Inertial>& x,
312  tmpl::list<RequestedTags...> /*meta*/) const noexcept {
313  using VarsComputer = detail::SchwarzschildVariables<DataType>;
314  typename VarsComputer::Cache cache{
315  get_size(*x.begin()),
316  VarsComputer{
317  {{std::nullopt, std::nullopt}}, x, mass_, coordinate_system_}};
318  return {cache.get_var(RequestedTags{})...};
319  }
320 
321  template <typename... RequestedTags>
322  tuples::TaggedTuple<RequestedTags...> variables(
323  const tnsr::I<DataVector, 3, Frame::Inertial>& x, const Mesh<3>& mesh,
324  const InverseJacobian<DataVector, 3, Frame::Logical, Frame::Inertial>&
325  inv_jacobian,
326  tmpl::list<RequestedTags...> /*meta*/) const noexcept {
327  using VarsComputer = detail::SchwarzschildVariables<DataVector>;
328  typename VarsComputer::Cache cache{
329  get_size(*x.begin()),
330  VarsComputer{{{mesh, inv_jacobian}}, x, mass_, coordinate_system_}};
331  return {cache.get_var(RequestedTags{})...};
332  }
333 
334  void pup(PUP::er& p) noexcept override {
335  Base::pup(p);
336  detail::SchwarzschildImpl::pup(p);
337  }
338 };
339 
340 /// \cond
341 template <typename Registrars>
342 PUP::able::PUP_ID Schwarzschild<Registrars>::my_PUP_ID = 0; // NOLINT
343 /// \endcond
344 
345 } // namespace Xcts::Solutions
domain::push_front
CoordinateMap< SourceFrame, TargetFrame, NewMap, Maps... > push_front(CoordinateMap< SourceFrame, TargetFrame, Maps... > old_map, NewMap new_map) noexcept
Creates a CoordinateMap by prepending the new map to the beginning of the old maps.
gr::Tags::TraceExtrinsicCurvature
Definition: Tags.hpp:120
Xcts
Items related to solving the Extended Conformal Thin Sandwich (XCTS) decomposition of the Einstein co...
Definition: Flatness.hpp:22
std::string
CharmPupable.hpp
Frame::Inertial
Definition: IndexType.hpp:44
std::rel_ops::operator!=
T operator!=(T... args)
Options.hpp
gr::Tags::EnergyDensity
The energy density , where denotes the normal to the spatial hypersurface.
Definition: Tags.hpp:138
ParseOptions.hpp
get_size
decltype(auto) get_size(const T &t, SizeFunction size=GetContainerSize{}) noexcept
Retrieve the size of t if t.size() is a valid expression, otherwise if T is fundamental or a std::com...
Definition: ContainerHelpers.hpp:145
Options::Option
Definition: Options.hpp:108
gr::lapse
Scalar< DataType > lapse(const tnsr::I< DataType, SpatialDim, Frame > &shift, const tnsr::aa< DataType, SpatialDim, Frame > &spacetime_metric) noexcept
Compute lapse from shift and spacetime metric.
WRAPPED_PUPable_decl_template
#define WRAPPED_PUPable_decl_template(className)
Mark derived classes as serializable.
Definition: CharmPupable.hpp:22
std::ostream
Options::create_from_yaml
Definition: MinmodType.hpp:11
Xcts::Solutions::SchwarzschildCoordinates
SchwarzschildCoordinates
Various coordinate systems in which to express the Schwarzschild solution.
Definition: Schwarzschild.hpp:35
Xcts::Solutions::AnalyticSolution
Base class for analytic solutions of the XCTS equations.
Definition: AnalyticSolution.hpp:34
tuples::TaggedTuple
An associative container that is indexed by structs.
Definition: TaggedTuple.hpp:271
std::numeric_limits::signaling_NaN
T signaling_NaN(T... args)
Xcts::Solutions::SchwarzschildCoordinates::Isotropic
@ Isotropic
Isotropic Schwarzschild coordinates.
Mesh
Holds the number of grid points, basis, and quadrature in each direction of the computational grid.
Definition: Mesh.hpp:48
tnsr
Type aliases to construct common Tensors.
Definition: TypeAliases.hpp:31
Tags::dt
Prefix indicating a time derivative.
Definition: Prefixes.hpp:29
ActionTesting::cache
Parallel::GlobalCache< Metavariables > & cache(MockRuntimeSystem< Metavariables > &runner, const ArrayIndex &array_index) noexcept
Returns the GlobalCache of Component with index array_index.
Definition: MockRuntimeSystemFreeFunctions.hpp:380
Scalar
Tensor< T, Symmetry<>, index_list<> > Scalar
Definition: TypeAliases.hpp:21
limits
gr::Tags::MomentumDensity
The spatial momentum density , where denotes the normal to the spatial hypersurface.
Definition: Tags.hpp:156
Gsl.hpp
Options::String
const char *const String
The string used in option structs.
Definition: Options.hpp:32
Frame
Definition: IndexType.hpp:36
gr
Definition: GaugeWave.hpp:27
Tensor.hpp
Tags::deriv
Prefix indicating spatial derivatives.
Definition: PartialDerivatives.hpp:52
PartialDerivatives.hpp
Xcts::Solutions
Analytic solutions of the XCTS equations.
Definition: AnalyticSolution.hpp:18
ostream
Prefixes.hpp
gr::Tags::StressTrace
The trace of the spatial stress-energy tensor .
Definition: Tags.hpp:147
gsl
Implementations from the Guideline Support Library.
Definition: ReadSpecThirdOrderPiecewisePolynomial.hpp:11
gr::Tags::Lapse
Definition: Tags.hpp:52
TMPL.hpp
Xcts::Solutions::Schwarzschild
Schwarzschild spacetime in general relativity.
Definition: Schwarzschild.hpp:288
gsl::not_null
Require a pointer to not be a nullptr
Definition: ReadSpecThirdOrderPiecewisePolynomial.hpp:13