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