-
Notifications
You must be signed in to change notification settings - Fork 184
/
PunctureField.hpp
59 lines (54 loc) · 2.59 KB
/
PunctureField.hpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
// Distributed under the MIT License.
// See LICENSE.txt for details.
#pragma once
#include <cstddef>
#include "DataStructures/DataBox/Prefixes.hpp"
#include "DataStructures/DataVector.hpp"
#include "DataStructures/Tensor/Tensor.hpp"
#include "DataStructures/Variables.hpp"
#include "Evolution/Systems/CurvedScalarWave/Tags.hpp"
#include "NumericalAlgorithms/LinearOperators/PartialDerivatives.hpp"
#include "Utilities/Gsl.hpp"
namespace CurvedScalarWave::Worldtube {
/*!
* \brief Computes the puncture/singular field \f$\Psi^\mathcal{P}\f$ of a
* scalar charge in circular orbit around a Schwarzschild black hole as
* described in \ref Detweiler2003.
*
* \details The field is computed using a Detweiler-Whiting singular
* Green's function and perturbatively expanded in the geodesic distance from
* the particle. It solves the inhomogeneous wave equation
* \f{align*}{
* \box \Psi^\mathcal{P} = -4 \pi q \int \sqrt{-g} \delta^4(x^i, z(\tau)) d \tau
* \f}
* where \f$q\f$ is the scalar charge and \f$z(\tau)\f$ is the worldline of the
* particle. The expression is expanded up to order n in geodesic distance and
* transformed to Kerr-Schild coordinates.
*
* The function given here assumes that the particle has scalar charge \f$q=1$\f
* and is on a fixed geodesic orbit \f$z(t) = (R \cos{\omega t}, R \sin{\omega
* t}, 0)\f$ around a Schwarzschild black hole, where R is the orbital radius
* and \f$ \omega = R^(-3/2) \f$. It returns the singular field at the requested
* coordinates as well as its time and spatial derivative.
*
* \note The expressions were computed with Mathematica and optimized by
* applying manual common subexpression elimination with sympy.
* The calculation is currently not vectorized to avoid a large amount of
* temporary allocations but this may be added as a future optimization. In its
* current form it does not noticeably impact the simulation time.
*/
void compute_puncture_field(
gsl::not_null<Variables<tmpl::list<
Tags::Psi, ::Tags::dt<Tags::Psi>,
::Tags::deriv<Tags::Psi, tmpl::size_t<3>, Frame::Inertial>>>*>
evolved_vars,
const tnsr::I<DataVector, 3, Frame::Inertial>& coords, const double time,
const double orbital_radius, const double bh_mass, size_t order);
void compute_puncture_field_0(
gsl::not_null<Variables<tmpl::list<
Tags::Psi, ::Tags::dt<Tags::Psi>,
::Tags::deriv<Tags::Psi, tmpl::size_t<3>, Frame::Inertial>>>*>
evolved_vars,
const tnsr::I<DataVector, 3, Frame::Inertial>& coords, const double time,
const double orbital_radius, const double bh_mass);
} // namespace CurvedScalarWave::Worldtube