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 <optional>
8 :
9 : #include "DataStructures/DataBox/Tag.hpp"
10 : #include "DataStructures/Tensor/TypeAliases.hpp"
11 : #include "Domain/FaceNormal.hpp"
12 : #include "Domain/TagsTimeDependent.hpp"
13 : #include "Evolution/Systems/GeneralizedHarmonic/ConstraintDamping/Tags.hpp"
14 : #include "Evolution/Systems/GeneralizedHarmonic/Tags.hpp"
15 : #include "ParallelAlgorithms/ApparentHorizonFinder/HorizonAliases.hpp"
16 : #include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
17 : #include "Utilities/TMPL.hpp"
18 :
19 : /// \cond
20 : class DataVector;
21 : namespace gsl {
22 : template <class T>
23 : class not_null;
24 : } // namespace gsl
25 : namespace Tags {
26 : template <typename Tag>
27 : struct Normalized;
28 : } // namespace Tags
29 : /// \endcond
30 :
31 : // IWYU pragma: no_forward_declare Tensor
32 :
33 : namespace gh {
34 : /// @{
35 : /*!
36 : * \brief Compute the characteristic speeds for the generalized harmonic system.
37 : *
38 : * Computes the speeds as described in "A New Generalized Harmonic
39 : * Evolution System" by Lindblom et. al \cite Lindblom2005qh
40 : * [see text following Eq.(34)]. The characteristic fields' names used here
41 : * differ from this paper:
42 : *
43 : * \f{align*}
44 : * \mathrm{SpECTRE} && \mathrm{Lindblom} \\
45 : * u^{\psi}_{ab} && u^\hat{0}_{ab} \\
46 : * u^0_{iab} && u^\hat{2}_{iab} \\
47 : * u^{\pm}_{ab} && u^{\hat{1}\pm}_{ab}
48 : * \f}
49 : *
50 : * The corresponding characteristic speeds \f$v\f$ are given in the text between
51 : * Eq.(34) and Eq.(35) of \cite Lindblom2005qh :
52 : *
53 : * \f{align*}
54 : * v_{\psi} =& -(1 + \gamma_1) n_k \beta^k - \gamma_1 n_k v^k_g \\
55 : * v_{0} =& -n_k \beta^k \\
56 : * v_{\pm} =& -n_k \beta^k \pm \alpha
57 : * \f}
58 : *
59 : * where \f$\alpha, \beta^k\f$ are the lapse and shift respectively,
60 : * \f$\gamma_1\f$ is a constraint damping parameter, \f$n_k\f$ is the unit
61 : * normal to the surface, and $v^k_g$ is the (optional) mesh velocity.
62 : */
63 : template <size_t Dim, typename Frame>
64 1 : std::array<DataVector, 4> characteristic_speeds(
65 : const Scalar<DataVector>& gamma_1, const Scalar<DataVector>& lapse,
66 : const tnsr::I<DataVector, Dim, Frame>& shift,
67 : const tnsr::i<DataVector, Dim, Frame>& unit_normal_one_form,
68 : const std::optional<tnsr::I<DataVector, Dim, Frame>>& mesh_velocity);
69 :
70 : template <size_t Dim, typename Frame>
71 1 : void characteristic_speeds(
72 : gsl::not_null<std::array<DataVector, 4>*> char_speeds,
73 : const Scalar<DataVector>& gamma_1, const Scalar<DataVector>& lapse,
74 : const tnsr::I<DataVector, Dim, Frame>& shift,
75 : const tnsr::i<DataVector, Dim, Frame>& unit_normal_one_form,
76 : const std::optional<tnsr::I<DataVector, Dim, Frame>>& mesh_velocity);
77 :
78 : template <size_t Dim, typename Frame>
79 0 : struct CharacteristicSpeedsCompute
80 : : Tags::CharacteristicSpeeds<DataVector, Dim, Frame>,
81 : db::ComputeTag {
82 0 : using base = Tags::CharacteristicSpeeds<DataVector, Dim, Frame>;
83 0 : using type = typename base::type;
84 0 : using argument_tags = tmpl::list<
85 : ::gh::ConstraintDamping::Tags::ConstraintGamma1,
86 : gr::Tags::Lapse<DataVector>, gr::Tags::Shift<DataVector, Dim, Frame>,
87 : ::Tags::Normalized<domain::Tags::UnnormalizedFaceNormal<Dim, Frame>>,
88 : domain::Tags::MeshVelocity<Dim, Frame>>;
89 :
90 0 : using return_type = typename base::type;
91 :
92 0 : static void function(
93 : const gsl::not_null<return_type*> result,
94 : const Scalar<DataVector>& gamma_1, const Scalar<DataVector>& lapse,
95 : const tnsr::I<DataVector, Dim, Frame>& shift,
96 : const tnsr::i<DataVector, Dim, Frame>& unit_normal_one_form,
97 : const std::optional<tnsr::I<DataVector, Dim, Frame>>& mesh_velocity) {
98 : characteristic_speeds(result, gamma_1, lapse, shift, unit_normal_one_form,
99 : mesh_velocity);
100 : };
101 : };
102 :
103 : // Simple tag used when observing the characteristic speeds on a Strahlkorper
104 : template <typename Frame>
105 0 : struct CharacteristicSpeedsOnStrahlkorper : db::SimpleTag {
106 0 : using type = tnsr::a<DataVector, 3, Frame>;
107 : };
108 :
109 : // Compute tag used when computing the characteristic speeds on a Strahlkorper
110 : // from gamma1, the lapse, shift, and the unit normal one form.
111 : template <size_t Dim, typename Frame>
112 0 : struct CharacteristicSpeedsOnStrahlkorperCompute
113 : : CharacteristicSpeedsOnStrahlkorper<Frame>,
114 : db::ComputeTag {
115 0 : using base = CharacteristicSpeedsOnStrahlkorper<Frame>;
116 0 : using type = typename base::type;
117 0 : using argument_tags =
118 : tmpl::list<::gh::ConstraintDamping::Tags::ConstraintGamma1,
119 : gr::Tags::Lapse<DataVector>,
120 : gr::Tags::Shift<DataVector, Dim, Frame>,
121 : ::ylm::Tags::UnitNormalOneForm<Frame>>;
122 :
123 0 : using return_type = typename base::type;
124 :
125 0 : static void function(
126 : const gsl::not_null<return_type*> result,
127 : const Scalar<DataVector>& gamma_1, const Scalar<DataVector>& lapse,
128 : const tnsr::I<DataVector, Dim, Frame>& shift,
129 : const tnsr::i<DataVector, Dim, Frame>& unit_normal_one_form) {
130 : auto array_char_speeds = make_array<4>(DataVector(get(lapse).size(), 0.0));
131 : characteristic_speeds(make_not_null(&array_char_speeds), gamma_1, lapse,
132 : shift, unit_normal_one_form, {});
133 : for (size_t i = 0; i < 4; ++i) {
134 : (*result)[i] = array_char_speeds[i];
135 : }
136 : }
137 : };
138 :
139 : /// @}
140 :
141 : /// @{
142 : /*!
143 : * \brief Computes characteristic fields from evolved fields
144 : *
145 : * \ref CharacteristicFieldsCompute and
146 : * \ref EvolvedFieldsFromCharacteristicFieldsCompute convert between
147 : * characteristic and evolved fields for the generalized harmonic system.
148 : *
149 : * \ref CharacteristicFieldsCompute computes
150 : * characteristic fields as described in "A New Generalized Harmonic
151 : * Evolution System" by Lindblom et. al \cite Lindblom2005qh .
152 : * Their names used here differ from this paper:
153 : *
154 : * \f{align*}
155 : * \mathrm{SpECTRE} && \mathrm{Lindblom} \\
156 : * u^{\psi}_{ab} && u^\hat{0}_{ab} \\
157 : * u^0_{iab} && u^\hat{2}_{iab} \\
158 : * u^{\pm}_{ab} && u^{\hat{1}\pm}_{ab}
159 : * \f}
160 : *
161 : * The characteristic fields \f$u\f$ are given in terms of the evolved fields by
162 : * Eq.(32) - (34) of \cite Lindblom2005qh, respectively:
163 : * \f{align*}
164 : * u^{\psi}_{ab} =& \psi_{ab} \\
165 : * u^0_{iab} =& (\delta^k_i - n_i n^k) \Phi_{kab} := P^k_i \Phi_{kab} \\
166 : * u^{\pm}_{ab} =& \Pi_{ab} \pm n^i \Phi_{iab} - \gamma_2\psi_{ab}
167 : * \f}
168 : *
169 : * where \f$\psi_{ab}\f$ is the spacetime metric, \f$\Pi_{ab}\f$ and
170 : * \f$\Phi_{iab}\f$ are evolved generalized harmonic fields introduced by first
171 : * derivatives of \f$\psi_{ab}\f$, \f$\gamma_2\f$ is a constraint damping
172 : * parameter, and \f$n_k\f$ is the unit normal to the surface.
173 : *
174 : * \ref EvolvedFieldsFromCharacteristicFieldsCompute computes evolved fields
175 : * \f$w\f$ in terms of the characteristic fields. This uses the inverse of
176 : * above relations:
177 : *
178 : * \f{align*}
179 : * \psi_{ab} =& u^{\psi}_{ab}, \\
180 : * \Pi_{ab} =& \frac{1}{2}(u^{+}_{ab} + u^{-}_{ab}) + \gamma_2 u^{\psi}_{ab}, \\
181 : * \Phi_{iab} =& \frac{1}{2}(u^{+}_{ab} - u^{-}_{ab}) n_i + u^0_{iab}.
182 : * \f}
183 : *
184 : * The corresponding characteristic speeds \f$v\f$ are computed by
185 : * \ref CharacteristicSpeedsCompute .
186 : */
187 : template <size_t Dim, typename Frame>
188 : typename Tags::CharacteristicFields<DataVector, Dim, Frame>::type
189 1 : characteristic_fields(
190 : const Scalar<DataVector>& gamma_2,
191 : const tnsr::II<DataVector, Dim, Frame>& inverse_spatial_metric,
192 : const tnsr::aa<DataVector, Dim, Frame>& spacetime_metric,
193 : const tnsr::aa<DataVector, Dim, Frame>& pi,
194 : const tnsr::iaa<DataVector, Dim, Frame>& phi,
195 : const tnsr::i<DataVector, Dim, Frame>& unit_normal_one_form);
196 :
197 : template <size_t Dim, typename Frame>
198 1 : void characteristic_fields(
199 : gsl::not_null<
200 : typename Tags::CharacteristicFields<DataVector, Dim, Frame>::type*>
201 : char_fields,
202 : const Scalar<DataVector>& gamma_2,
203 : const tnsr::II<DataVector, Dim, Frame>& inverse_spatial_metric,
204 : const tnsr::aa<DataVector, Dim, Frame>& spacetime_metric,
205 : const tnsr::aa<DataVector, Dim, Frame>& pi,
206 : const tnsr::iaa<DataVector, Dim, Frame>& phi,
207 : const tnsr::i<DataVector, Dim, Frame>& unit_normal_one_form);
208 :
209 : template <size_t Dim, typename Frame>
210 0 : struct CharacteristicFieldsCompute
211 : : Tags::CharacteristicFields<DataVector, Dim, Frame>,
212 : db::ComputeTag {
213 0 : using base = Tags::CharacteristicFields<DataVector, Dim, Frame>;
214 0 : using return_type = typename base::type;
215 0 : using argument_tags = tmpl::list<
216 : ::gh::ConstraintDamping::Tags::ConstraintGamma2,
217 : gr::Tags::InverseSpatialMetric<DataVector, Dim, Frame>,
218 : gr::Tags::SpacetimeMetric<DataVector, Dim, Frame>,
219 : Tags::Pi<DataVector, Dim, Frame>, Tags::Phi<DataVector, Dim, Frame>,
220 : ::Tags::Normalized<domain::Tags::UnnormalizedFaceNormal<Dim, Frame>>>;
221 :
222 0 : static constexpr auto function = static_cast<void (*)(
223 : const gsl::not_null<return_type*>, const Scalar<DataVector>&,
224 : const tnsr::II<DataVector, Dim, Frame>&,
225 : const tnsr::aa<DataVector, Dim, Frame>&,
226 : const tnsr::aa<DataVector, Dim, Frame>&,
227 : const tnsr::iaa<DataVector, Dim, Frame>&,
228 : const tnsr::i<DataVector, Dim, Frame>&)>(&characteristic_fields);
229 : };
230 : /// @}
231 :
232 : /// @{
233 : /*!
234 : * \brief For expressions used here to compute evolved fields from
235 : * characteristic ones, see \ref CharacteristicFieldsCompute.
236 : */
237 : template <size_t Dim, typename Frame>
238 : typename Tags::EvolvedFieldsFromCharacteristicFields<DataVector, Dim,
239 : Frame>::type
240 1 : evolved_fields_from_characteristic_fields(
241 : const Scalar<DataVector>& gamma_2,
242 : const tnsr::aa<DataVector, Dim, Frame>& u_psi,
243 : const tnsr::iaa<DataVector, Dim, Frame>& u_zero,
244 : const tnsr::aa<DataVector, Dim, Frame>& u_plus,
245 : const tnsr::aa<DataVector, Dim, Frame>& u_minus,
246 : const tnsr::i<DataVector, Dim, Frame>& unit_normal_one_form);
247 :
248 : template <size_t Dim, typename Frame>
249 1 : void evolved_fields_from_characteristic_fields(
250 : gsl::not_null<typename Tags::EvolvedFieldsFromCharacteristicFields<
251 : DataVector, Dim, Frame>::type*>
252 : evolved_fields,
253 : const Scalar<DataVector>& gamma_2,
254 : const tnsr::aa<DataVector, Dim, Frame>& u_psi,
255 : const tnsr::iaa<DataVector, Dim, Frame>& u_zero,
256 : const tnsr::aa<DataVector, Dim, Frame>& u_plus,
257 : const tnsr::aa<DataVector, Dim, Frame>& u_minus,
258 : const tnsr::i<DataVector, Dim, Frame>& unit_normal_one_form);
259 :
260 : template <size_t Dim, typename Frame>
261 0 : struct EvolvedFieldsFromCharacteristicFieldsCompute
262 : : Tags::EvolvedFieldsFromCharacteristicFields<DataVector, Dim, Frame>,
263 : db::ComputeTag {
264 0 : using base =
265 : Tags::EvolvedFieldsFromCharacteristicFields<DataVector, Dim, Frame>;
266 0 : using return_type = typename base::type;
267 0 : using argument_tags = tmpl::list<
268 : gh::ConstraintDamping::Tags::ConstraintGamma2,
269 : Tags::VSpacetimeMetric<DataVector, Dim, Frame>,
270 : Tags::VZero<DataVector, Dim, Frame>, Tags::VPlus<DataVector, Dim, Frame>,
271 : Tags::VMinus<DataVector, Dim, Frame>,
272 : ::Tags::Normalized<domain::Tags::UnnormalizedFaceNormal<Dim, Frame>>>;
273 :
274 0 : static constexpr auto function = static_cast<void (*)(
275 : const gsl::not_null<return_type*>, const Scalar<DataVector>& gamma_2,
276 : const tnsr::aa<DataVector, Dim, Frame>& u_psi,
277 : const tnsr::iaa<DataVector, Dim, Frame>& u_zero,
278 : const tnsr::aa<DataVector, Dim, Frame>& u_plus,
279 : const tnsr::aa<DataVector, Dim, Frame>& u_minus,
280 : const tnsr::i<DataVector, Dim, Frame>& unit_normal_one_form)>(
281 : &evolved_fields_from_characteristic_fields);
282 : };
283 : /// @}
284 :
285 1 : namespace Tags{
286 0 : struct LargestCharacteristicSpeed : db::SimpleTag {
287 0 : using type = double;
288 : };
289 : /*!
290 : * \brief Computes the largest magnitude of the characteristic speeds.
291 : */
292 : template <size_t Dim, typename Frame>
293 1 : struct ComputeLargestCharacteristicSpeed : db::ComputeTag,
294 : LargestCharacteristicSpeed {
295 0 : using argument_tags =
296 : tmpl::list<::gh::ConstraintDamping::Tags::ConstraintGamma1,
297 : gr::Tags::Lapse<DataVector>,
298 : gr::Tags::Shift<DataVector, Dim, Frame>,
299 : gr::Tags::SpatialMetric<DataVector, Dim, Frame>>;
300 0 : using return_type = double;
301 0 : using base = LargestCharacteristicSpeed;
302 0 : static void function(const gsl::not_null<double*> speed,
303 : const Scalar<DataVector>& gamma_1,
304 : const Scalar<DataVector>& lapse,
305 : const tnsr::I<DataVector, Dim, Frame>& shift,
306 : const tnsr::ii<DataVector, Dim, Frame>& spatial_metric);
307 : };
308 : } // namespace Tags
309 : } // namespace gh
|
