Line data Source code
1 0 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : #pragma once
5 :
6 : #include <cstddef>
7 : #include <string>
8 : #include <vector>
9 :
10 : #include "DataStructures/ComplexDataVector.hpp"
11 : #include "DataStructures/Tensor/TypeAliases.hpp"
12 : #include "Domain/Tags.hpp"
13 : #include "Elliptic/BoundaryConditions/BoundaryCondition.hpp"
14 : #include "Elliptic/BoundaryConditions/BoundaryConditionType.hpp"
15 : #include "Elliptic/Systems/SelfForce/GeneralRelativity/Tags.hpp"
16 : #include "Options/String.hpp"
17 : #include "Utilities/Gsl.hpp"
18 : #include "Utilities/TMPL.hpp"
19 :
20 : namespace GrSelfForce::BoundaryConditions {
21 :
22 : /*!
23 : * \brief Radial Sommerfeld boundary conditions for the m-mode metric
24 : * perturbations.
25 : *
26 : * These radial boundary conditions enforce an outgoing direction of radiation
27 : * propagation. They apply both near the Kerr horizon (inner radial
28 : * boundary) and at large distance (outer radial boundary):
29 : *
30 : * \begin{equation}
31 : * n_i F^i = i m \Omega \Psi_m
32 : * \end{equation}
33 : *
34 : * With hyperboloidal slicing (see `HyperboloidalSlicing`), the outgoing
35 : * radiation condition above is replaced by a Neumann condition
36 : * ($F^i = 0$) at first order (`Order = 1`). At second order (`Order = 2`),
37 : * it instead imposes
38 : *
39 : * \begin{equation}
40 : * \partial_{x^i}\left(\alpha^{i}\partial_{x^i}\Psi\right) = 0,
41 : * \end{equation}
42 : *
43 : * i.e. $\gamma^{i}\cdot\partial_{x^i}\Psi = -\beta\cdot\Psi$, where $x^i$ is
44 : * the radial coordinate ($r_\star$ or $r$, depending on
45 : * `PenetratingHorizon`; see `CircularOrbit`). Since $\gamma^{i}$ and $\beta$
46 : * couple all tensor components of $\Psi$ together, we invert this dense
47 : * matrix in tensor-component space to solve for $\partial_{x^i}\Psi$, then
48 : * set the flux $F^{i} = \alpha^{i}\partial_{x^i}\Psi$.
49 : *
50 : * These boundary conditions currently assume a circular equatorial orbit.
51 : */
52 1 : class Sommerfeld : public elliptic::BoundaryConditions::BoundaryCondition<2> {
53 : private:
54 0 : using Base = elliptic::BoundaryConditions::BoundaryCondition<2>;
55 0 : using GradTensorType =
56 : TensorMetafunctions::prepend_spatial_index<tnsr::aa<ComplexDataVector, 3>,
57 : 2, UpLo::Lo, Frame::Inertial>;
58 :
59 : public:
60 0 : struct BlackHoleMass {
61 0 : static constexpr Options::String help =
62 : "Kerr mass parameter 'M' of the black hole";
63 0 : using type = double;
64 : };
65 0 : struct BlackHoleSpin {
66 0 : static constexpr Options::String help =
67 : "Kerr dimensionless spin parameter 'chi' of the black hole";
68 0 : using type = double;
69 : };
70 0 : struct OrbitalRadius {
71 0 : static constexpr Options::String help =
72 : "Radius 'r_0' of the circular orbit";
73 0 : using type = double;
74 : };
75 0 : struct MModeNumber {
76 0 : static constexpr Options::String help =
77 : "Mode number 'm' of the scalar field";
78 0 : using type = int;
79 : };
80 0 : struct HyperboloidalSlicing {
81 0 : static constexpr Options::String help =
82 : "Whether hyperboloidal slicing is applied. If true, a simple Neumann "
83 : "boundary condition is applied.";
84 0 : using type = bool;
85 : };
86 0 : struct Order {
87 0 : static constexpr Options::String help =
88 : "Order of the boundary condition. First order (Order=1) implements "
89 : "just the leading 'i m Omega' term. Second order (Order=2) includes "
90 : "the next-to-leading '1/r' term as well (Robin-type).";
91 0 : using type = int;
92 : };
93 :
94 0 : static constexpr Options::String help =
95 : "Radial Sommerfeld boundary condition";
96 0 : using options = tmpl::list<BlackHoleMass, BlackHoleSpin, OrbitalRadius,
97 : MModeNumber, HyperboloidalSlicing, Order>;
98 0 : Sommerfeld() = default;
99 0 : Sommerfeld(const Sommerfeld&) = default;
100 0 : Sommerfeld& operator=(const Sommerfeld&) = default;
101 0 : Sommerfeld(Sommerfeld&&) = default;
102 0 : Sommerfeld& operator=(Sommerfeld&&) = default;
103 0 : ~Sommerfeld() override = default;
104 :
105 0 : explicit Sommerfeld(double black_hole_mass, double black_hole_spin,
106 : double orbital_radius, int m_mode_number,
107 : bool hyperboloidal_slicing, int order);
108 :
109 0 : double black_hole_mass() const { return black_hole_mass_; }
110 0 : double black_hole_spin() const { return black_hole_spin_; }
111 0 : double orbital_radius() const { return orbital_radius_; }
112 0 : int m_mode_number() const { return m_mode_number_; }
113 0 : bool hyperboloidal_slicing() const { return hyperboloidal_slicing_; }
114 0 : int order() const { return order_; }
115 :
116 : /// \cond
117 : explicit Sommerfeld(CkMigrateMessage* m);
118 : using PUP::able::register_constructor;
119 : WRAPPED_PUPable_decl_template(Sommerfeld);
120 : /// \endcond
121 :
122 0 : std::unique_ptr<domain::BoundaryConditions::BoundaryCondition> get_clone()
123 : const override;
124 :
125 0 : std::vector<elliptic::BoundaryConditionType> boundary_condition_types()
126 : const override {
127 : return {elliptic::BoundaryConditionType::Neumann};
128 : }
129 :
130 0 : using argument_tags = tmpl::list<Tags::Alpha, Tags::Beta, Tags::GammaRstar>;
131 0 : using volume_tags = tmpl::list<>;
132 :
133 0 : void apply(gsl::not_null<tnsr::aa<ComplexDataVector, 3>*> field,
134 : gsl::not_null<tnsr::aa<ComplexDataVector, 3>*> n_dot_flux,
135 : const GradTensorType& deriv_field,
136 : const tnsr::I<ComplexDataVector, 2>& alpha,
137 : const tnsr::aaBB<ComplexDataVector, 3>& beta,
138 : const tnsr::aaBB<ComplexDataVector, 3>& gamma_rstar) const;
139 :
140 0 : using argument_tags_linearized =
141 : tmpl::list<Tags::Alpha, Tags::Beta, Tags::GammaRstar>;
142 0 : using volume_tags_linearized = tmpl::list<>;
143 :
144 0 : void apply_linearized(
145 : gsl::not_null<tnsr::aa<ComplexDataVector, 3>*> field_correction,
146 : gsl::not_null<tnsr::aa<ComplexDataVector, 3>*> n_dot_flux_correction,
147 : const GradTensorType& deriv_field_correction,
148 : const tnsr::I<ComplexDataVector, 2>& alpha,
149 : const tnsr::aaBB<ComplexDataVector, 3>& beta,
150 : const tnsr::aaBB<ComplexDataVector, 3>& gamma_rstar) const;
151 :
152 : // NOLINTNEXTLINE
153 0 : void pup(PUP::er& p) override;
154 :
155 : private:
156 0 : friend bool operator==(const Sommerfeld& lhs, const Sommerfeld& rhs);
157 :
158 0 : double black_hole_mass_{std::numeric_limits<double>::signaling_NaN()};
159 0 : double black_hole_spin_{std::numeric_limits<double>::signaling_NaN()};
160 0 : double orbital_radius_{std::numeric_limits<double>::signaling_NaN()};
161 0 : int m_mode_number_{};
162 0 : bool hyperboloidal_slicing_{false};
163 0 : int order_{1};
164 : };
165 :
166 0 : bool operator!=(const Sommerfeld& lhs, const Sommerfeld& rhs);
167 :
168 : } // namespace GrSelfForce::BoundaryConditions
|