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 <cstddef>
8 :
9 : #include "DataStructures/Tensor/TypeAliases.hpp"
10 : #include "Utilities/ContainerHelpers.hpp"
11 :
12 : /// \cond
13 : class DataVector;
14 : namespace gsl {
15 : template <class T>
16 : class not_null;
17 : } // namespace gsl
18 : /// \endcond
19 :
20 : namespace gh::BoundaryConditions {
21 : /// \brief Detailed implementation of Bjorhus-type boundary corrections
22 1 : namespace Bjorhus {
23 : /*!
24 : * \brief Computes the expression needed to set boundary conditions on the time
25 : * derivative of the characteristic field \f$v^{\psi}_{ab}\f$
26 : *
27 : * \details In the Bjorhus scheme, the time derivatives of evolved variables are
28 : * characteristic projected. A constraint-preserving correction term is added
29 : * here to the resulting characteristic (time-derivative) field:
30 : *
31 : * \f{align}
32 : * \Delta \partial_t v^{\psi}_{ab} = \lambda_{\psi} n^i C_{iab}
33 : * \f}
34 : *
35 : * where \f$n^i\f$ is the local unit normal to the external boundary,
36 : * \f$C_{iab} = \partial_i \psi_{ab} - \Phi_{iab}\f$ is the three-index
37 : * constraint, and \f$\lambda_{\psi}\f$ is the characteristic speed of the field
38 : * \f$v^{\psi}_{ab}\f$.
39 : */
40 : template <size_t VolumeDim, typename DataType>
41 1 : void constraint_preserving_bjorhus_corrections_dt_v_psi(
42 : gsl::not_null<tnsr::aa<DataType, VolumeDim, Frame::Inertial>*> bc_dt_v_psi,
43 : const tnsr::I<DataType, VolumeDim, Frame::Inertial>&
44 : unit_interface_normal_vector,
45 : const tnsr::iaa<DataType, VolumeDim, Frame::Inertial>&
46 : three_index_constraint,
47 : const std::array<DataType, 4>& char_speeds);
48 :
49 : /*!
50 : * \brief Computes the expression needed to set boundary conditions on the time
51 : * derivative of the characteristic field \f$v^{0}_{iab}\f$
52 : *
53 : * \details In the Bjorhus scheme, the time derivatives of evolved variables are
54 : * characteristic projected. A constraint-preserving correction term is added
55 : * here to the resulting characteristic (time-derivative) field \f$v^0_{iab}\f$:
56 : *
57 : * \f{align}
58 : * \Delta \partial_t v^{0}_{iab} = \lambda_{0} n^j C_{jiab}
59 : * \f}
60 : *
61 : * where \f$n^i\f$ is the local unit normal to the external boundary,
62 : * \f$C_{ijab} = \partial_i\Phi_{jab} - \partial_j\Phi_{iab}\f$ is the
63 : * four-index constraint, and \f$\lambda_{0}\f$ is the characteristic speed of
64 : * the field \f$v^0_{iab}\f$.
65 : *
66 : * \note In 3D, only the non-zero and unique components of the four-index
67 : * constraint are stored, as \f$\hat{C}_{iab} = \epsilon_i^{jk} C_{jkab}\f$. In
68 : * 2D the input expected here is \f$\hat{C}_{0ab} = C_{01ab}, \hat{C}_{1ab} =
69 : * C_{10ab}\f$, and in 1D \f$\hat{C}_{0ab} = 0\f$.
70 : */
71 : template <size_t VolumeDim, typename DataType>
72 1 : void constraint_preserving_bjorhus_corrections_dt_v_zero(
73 : gsl::not_null<tnsr::iaa<DataType, VolumeDim, Frame::Inertial>*>
74 : bc_dt_v_zero,
75 : const tnsr::I<DataType, VolumeDim, Frame::Inertial>&
76 : unit_interface_normal_vector,
77 : const tnsr::iaa<DataType, VolumeDim, Frame::Inertial>&
78 : four_index_constraint,
79 : const std::array<DataType, 4>& char_speeds);
80 :
81 : /// @{
82 : /*!
83 : * \brief Computes the expression needed to set boundary conditions on the time
84 : * derivative of the characteristic field \f$v^{-}_{ab}\f$
85 : *
86 : * \details In the Bjorhus scheme, the time derivatives of evolved variables are
87 : * characteristic projected. Constraint-preserving correction terms
88 : * \f$T^{\mathrm C}_{ab}\f$, Sommerfeld condition terms for gauge degrees of
89 : * freedom \f$T^{\mathrm G}_{ab}\f$, and terms constraining physical degrees of
90 : * freedom \f$T^{\mathrm P}_{ab}\f$ are added here to the resulting
91 : * characteristic (time-derivative) field:
92 : *
93 : * \f{align}
94 : * \Delta \partial_t v^{-}_{ab} = T^{\mathrm C}_{ab} + T^{\mathrm G}_{ab}
95 : * + T^{\mathrm P}_{ab}.
96 : * \f}
97 : *
98 : * These terms are given by Eq. (64) of \cite Lindblom2005qh :
99 : *
100 : * \f{eqnarray}
101 : * T^{\mathrm C}_{ab} &=& \frac{1}{2}
102 : * \left(2 k^c k^d l_a l_b - k^c l_b P^d_a - k^c l_a P^d_b - k^d l_b P^c_a
103 : * - k^d l_a P^c_b + P^{cd} P_{ab}\right) \partial_t v^{-}_{cd}\\
104 : * &&+ \frac{1}{\sqrt{2}} \lambda_{-}c^{\hat{0}-}_c \left(l_a l_b k^c +
105 : * P_{ab} l^c - P^c_b l_a - P^c_a l_b \right) \nonumber
106 : * \f}
107 : *
108 : * where \f$l^a (l_a)\f$ is the outgoing null vector (one-form), \f$k^a (k_a)\f$
109 : * is the incoming null vector (one-form), \f$P_{ab}, P^{ab}, P^a_b\f$ are
110 : * spacetime projection operators defined in `transverse_projection_operator()`,
111 : * \f$\partial_t v^{-}_{ab}\f$ is the characteristic projected time derivative
112 : * of evolved variables (corresponding to the \f$v^{-}\f$ field), and
113 : * \f$c^{\hat{0}\pm}_a\f$ are characteristic modes of the constraint evolution
114 : * system:
115 : *
116 : * \f{align}\nonumber
117 : * c^{\hat{0}\pm}_a = F_a \mp n^k C_{ka},
118 : * \f}
119 : *
120 : * where \f$F_a\f$ is the generalized-harmonic (GH) F constraint [Eq. (43) of
121 : * \cite Lindblom2005qh], \f$C_{ka}\f$ is the GH 2-index constraint [Eq. (44)
122 : * of \cite Lindblom2005qh], and \f$n^k\f$ is the unit spatial normal to the
123 : * outer boundary. Boundary correction terms that prevent strong reflections of
124 : * gauge perturbations are given by Eq. (25) of \cite Rinne2007ui :
125 : *
126 : * \f{align}
127 : * T^{\mathrm G}_{ab} =
128 : * \left(k_a P^c_b l^d + k_b P^c_a l^d -
129 : * \left(k_a l_b k^c l^d + k_b l_a k_c l^d + k_a k_b l^c l^d
130 : * \right) \right) \left(\gamma_2 - \frac{1}{r}
131 : * \right) \partial_t v^{\psi}_{cd}
132 : * \f}
133 : *
134 : * where \f$r\f$ is the radial coordinate at the outer boundary, which is
135 : * assumed to be spherical, \f$\gamma_2\f$ is a GH constraint damping parameter,
136 : * and \f$\partial_t v^{\psi}_{ab}\f$ is the characteristic projected time
137 : * derivative of evolved variables (corresponding to the \f$v^{\psi}\f$ field).
138 : * Finally, we constrain physical degrees of freedom using corrections from
139 : * Eq. (68) of \cite Lindblom2005qh :
140 : *
141 : * \f{align}
142 : * T^{\mathrm P}_{ab} = \left( P^c_a P^d_b - \frac{1}{2} P_{ab} P^{cd} \right)
143 : * \left(\partial_t v^{-}_{cd} + \lambda_{-}
144 : * \left(U^{3-}_{cd} - \gamma_2 n^i C_{icd}\right)\right),
145 : * \f}
146 : *
147 : * where \f$C_{icd}\f$ is the GH 3-index constraint [c.f. Eq. (26) of
148 : * \cite Lindblom2005qh, see also `three_index_constraint()`], and
149 : *
150 : * \f{align}\nonumber
151 : * U^{3-}_{ab} = 2 P^i_a P^j_b U^{8-}_{ij}
152 : * \f}
153 : *
154 : * is the inward propagating characteristic mode of the Weyl tensor evolution
155 : * \f$U^{8-}_{ab}\f$ [c.f. Eq. 75 of \cite Kidder2004rw ] projected onto the
156 : * outer boundary using the spatial-spacetime projection operators \f$P^i_a\f$.
157 : * Note that (A) the covariant derivative of extrinsic curvature needed to
158 : * get the Weyl propagating modes is calculated substituting the evolved
159 : * variable \f$\Phi_{iab}\f$ for spatial derivatives of the spacetime metric;
160 : * and (B) the spatial Ricci tensor used in the same calculation is also
161 : * calculated using the same substituion, and also includes corrections
162 : * proportional to the GH 4-index constraint \f$C_{ijab}\f$ [c.f. Eq. (45) of
163 : * \cite Lindblom2005qh , see also `four_index_constraint()`]:
164 : *
165 : * \f{align}\nonumber
166 : * R_{ij} \rightarrow R_{ij} + C_{(iklj)} + \frac{1}{2} n^k q^a C_{(ikj)a},
167 : * \f}
168 : *
169 : * where \f$q^a\f$ is the future-directed spacetime normal vector.
170 : */
171 : template <size_t VolumeDim, typename DataType>
172 1 : void constraint_preserving_bjorhus_corrections_dt_v_minus(
173 : gsl::not_null<tnsr::aa<DataType, VolumeDim, Frame::Inertial>*>
174 : bc_dt_v_minus,
175 : const Scalar<DataType>& gamma2,
176 : const tnsr::I<DataType, VolumeDim, Frame::Inertial>& inertial_coords,
177 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>& incoming_null_one_form,
178 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>& outgoing_null_one_form,
179 : const tnsr::A<DataType, VolumeDim, Frame::Inertial>& incoming_null_vector,
180 : const tnsr::A<DataType, VolumeDim, Frame::Inertial>& outgoing_null_vector,
181 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>& projection_ab,
182 : const tnsr::Ab<DataType, VolumeDim, Frame::Inertial>& projection_Ab,
183 : const tnsr::AA<DataType, VolumeDim, Frame::Inertial>& projection_AB,
184 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>&
185 : char_projected_rhs_dt_v_psi,
186 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>&
187 : char_projected_rhs_dt_v_minus,
188 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>&
189 : constraint_char_zero_plus,
190 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>&
191 : constraint_char_zero_minus,
192 : const std::array<DataType, 4>& char_speeds);
193 :
194 : template <size_t VolumeDim, typename DataType>
195 1 : void constraint_preserving_physical_bjorhus_corrections_dt_v_minus(
196 : gsl::not_null<tnsr::aa<DataType, VolumeDim, Frame::Inertial>*>
197 : bc_dt_v_minus,
198 : const Scalar<DataType>& gamma2,
199 : const tnsr::I<DataType, VolumeDim, Frame::Inertial>& inertial_coords,
200 : const tnsr::i<DataType, VolumeDim, Frame::Inertial>&
201 : unit_interface_normal_one_form,
202 : const tnsr::I<DataType, VolumeDim, Frame::Inertial>&
203 : unit_interface_normal_vector,
204 : const tnsr::A<DataType, VolumeDim, Frame::Inertial>&
205 : spacetime_unit_normal_vector,
206 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>& incoming_null_one_form,
207 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>& outgoing_null_one_form,
208 : const tnsr::A<DataType, VolumeDim, Frame::Inertial>& incoming_null_vector,
209 : const tnsr::A<DataType, VolumeDim, Frame::Inertial>& outgoing_null_vector,
210 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>& projection_ab,
211 : const tnsr::Ab<DataType, VolumeDim, Frame::Inertial>& projection_Ab,
212 : const tnsr::AA<DataType, VolumeDim, Frame::Inertial>& projection_AB,
213 : const tnsr::II<DataType, VolumeDim, Frame::Inertial>&
214 : inverse_spatial_metric,
215 : const tnsr::ii<DataType, VolumeDim, Frame::Inertial>& extrinsic_curvature,
216 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>& spacetime_metric,
217 : const tnsr::AA<DataType, VolumeDim, Frame::Inertial>&
218 : inverse_spacetime_metric,
219 : const tnsr::iaa<DataType, VolumeDim, Frame::Inertial>&
220 : three_index_constraint,
221 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>&
222 : char_projected_rhs_dt_v_psi,
223 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>&
224 : char_projected_rhs_dt_v_minus,
225 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>&
226 : constraint_char_zero_plus,
227 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>&
228 : constraint_char_zero_minus,
229 : const tnsr::iaa<DataType, VolumeDim, Frame::Inertial>& phi,
230 : const tnsr::ijaa<DataType, VolumeDim, Frame::Inertial>& d_phi,
231 : const tnsr::iaa<DataType, VolumeDim, Frame::Inertial>& d_pi,
232 : const std::array<DataType, 4>& char_speeds);
233 : /// @}
234 :
235 : namespace detail {
236 : template <size_t VolumeDim, typename DataType>
237 : void add_gauge_sommerfeld_terms_to_dt_v_minus(
238 : const gsl::not_null<tnsr::aa<DataType, VolumeDim, Frame::Inertial>*>
239 : bc_dt_v_minus,
240 : const Scalar<DataType>& gamma2,
241 : const tnsr::I<DataType, VolumeDim, Frame::Inertial>& inertial_coords,
242 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>& incoming_null_one_form,
243 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>& outgoing_null_one_form,
244 : const tnsr::A<DataType, VolumeDim, Frame::Inertial>& incoming_null_vector,
245 : const tnsr::A<DataType, VolumeDim, Frame::Inertial>& outgoing_null_vector,
246 : const tnsr::Ab<DataType, VolumeDim, Frame::Inertial>& projection_Ab,
247 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>&
248 : char_projected_rhs_dt_v_psi);
249 :
250 : template <size_t VolumeDim, typename DataType>
251 : void add_constraint_dependent_terms_to_dt_v_minus(
252 : const gsl::not_null<tnsr::aa<DataType, VolumeDim, Frame::Inertial>*>
253 : bc_dt_v_minus,
254 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>& outgoing_null_one_form,
255 : const tnsr::A<DataType, VolumeDim, Frame::Inertial>& incoming_null_vector,
256 : const tnsr::A<DataType, VolumeDim, Frame::Inertial>& outgoing_null_vector,
257 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>& projection_ab,
258 : const tnsr::Ab<DataType, VolumeDim, Frame::Inertial>& projection_Ab,
259 : const tnsr::AA<DataType, VolumeDim, Frame::Inertial>& projection_AB,
260 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>&
261 : constraint_char_zero_plus,
262 : const tnsr::a<DataType, VolumeDim, Frame::Inertial>&
263 : constraint_char_zero_minus,
264 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>&
265 : char_projected_rhs_dt_v_minus,
266 : const std::array<DataType, 4>& char_speeds);
267 :
268 : template <size_t VolumeDim, typename DataType>
269 : void add_physical_terms_to_dt_v_minus(
270 : const gsl::not_null<tnsr::aa<DataType, VolumeDim, Frame::Inertial>*>
271 : bc_dt_v_minus,
272 : const Scalar<DataType>& gamma2,
273 : const tnsr::i<DataType, VolumeDim, Frame::Inertial>&
274 : unit_interface_normal_one_form,
275 : const tnsr::I<DataType, VolumeDim, Frame::Inertial>&
276 : unit_interface_normal_vector,
277 : const tnsr::A<DataType, VolumeDim, Frame::Inertial>&
278 : spacetime_unit_normal_vector,
279 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>& projection_ab,
280 : const tnsr::Ab<DataType, VolumeDim, Frame::Inertial>& projection_Ab,
281 : const tnsr::AA<DataType, VolumeDim, Frame::Inertial>& projection_AB,
282 : const tnsr::II<DataType, VolumeDim, Frame::Inertial>&
283 : inverse_spatial_metric,
284 : const tnsr::ii<DataType, VolumeDim, Frame::Inertial>& extrinsic_curvature,
285 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>& spacetime_metric,
286 : const tnsr::AA<DataType, VolumeDim, Frame::Inertial>&
287 : inverse_spacetime_metric,
288 : const tnsr::iaa<DataType, VolumeDim, Frame::Inertial>&
289 : three_index_constraint,
290 : const tnsr::aa<DataType, VolumeDim, Frame::Inertial>&
291 : char_projected_rhs_dt_v_minus,
292 : const tnsr::iaa<DataType, VolumeDim, Frame::Inertial>& phi,
293 : const tnsr::ijaa<DataType, VolumeDim, Frame::Inertial>& d_phi,
294 : const tnsr::iaa<DataType, VolumeDim, Frame::Inertial>& d_pi,
295 : const std::array<DataType, 4>& char_speeds);
296 : } // namespace detail
297 : } // namespace Bjorhus
298 : } // namespace gh::BoundaryConditions
|
