SwshDerivatives.hpp
1 // Distributed under the MIT License.
2 // See LICENSE.txt for details.
3 
4 #pragma once
5 
6 #include <cmath>
7 #include <cstddef>
8 #include <type_traits>
9 
10 #include "DataStructures/ComplexDataVector.hpp"
12 #include "DataStructures/Tags.hpp"
15 #include "Domain/Mesh.hpp"
16 #include "Evolution/Systems/Cce/IntegrandInputSteps.hpp"
17 #include "Evolution/Systems/Cce/Tags.hpp"
19 #include "NumericalAlgorithms/Spectral/SwshDerivatives.hpp"
20 #include "NumericalAlgorithms/Spectral/SwshTags.hpp"
22 #include "Utilities/Gsl.hpp"
23 #include "Utilities/TypeTraits.hpp"
24 #include "Utilities/TypeTraits/IsA.hpp"
25 
26 namespace Cce {
27 namespace detail {
28 // Precomputation routines for supplying the additional quantities necessary
29 // to correct the output of the angular derivative routines from the angular
30 // derivatives evaluated at constant numerical coordinates (which is what is
31 // returned after the libsharp evaluation) to the angular derivatives at
32 // constant Bondi radius (which is what appears in the literature equations
33 // and is simple to combine to obtain the hypersurface integrands).
34 //
35 // Warning: this 'on demand' template is a way of taking advantage of the blaze
36 // expression templates in a generic, modular way. However, this can be
37 // dangerous. The returned value MUST be a blaze expression template directly,
38 // and not a wrapper type (like `SpinWeighted`). Otherwise, some information is
39 // lost on the stack and the expression template is corrupted. So, in these 'on
40 // demand' returns, the arguments must be fully unpacked to vector types.
41 //
42 // The `Select` operand is a `cpp17::bool_constant` that ensures mutual
43 // exclusivity of the template specializations.
44 template <typename Tag, typename SpinConstant, typename Select>
45 struct OnDemandInputsForSwshJacobianImpl;
46 
47 // default to just retrieving it from the box if not providing an expression
48 // template shortcut
49 template <typename Tag>
50 struct OnDemandInputsForSwshJacobianImpl<
51  Tags::Dy<Tag>, std::integral_constant<int, Tag::type::type::spin>,
52  cpp17::bool_constant<not tt::is_a_v<::Tags::Multiplies, Tag> and
53  not tt::is_a_v<Tags::Dy, Tag>>> {
54  template <typename DataBoxTagList>
55  SPECTRE_ALWAYS_INLINE decltype(auto) operator()(
56  const db::DataBox<DataBoxTagList>& box) noexcept {
57  return get(db::get<Tags::Dy<Tag>>(box)).data();
58  }
59 };
60 
61 // default to retrieving from the box if the requested tag is a second
62 // derivative without an additional evaluation channel
63 template <typename Tag>
64 struct OnDemandInputsForSwshJacobianImpl<
65  Tags::Dy<Tags::Dy<Tag>>, std::integral_constant<int, Tag::type::type::spin>,
66  cpp17::bool_constant<not tt::is_a_v<::Tags::Multiplies, Tag>>> {
67  template <typename DataBoxTagList>
68  SPECTRE_ALWAYS_INLINE decltype(auto) operator()(
69  const db::DataBox<DataBoxTagList>& box) noexcept {
70  return get(db::get<Tags::Dy<Tags::Dy<Tag>>>(box)).data();
71  }
72 };
73 
74 // use the product rule to provide an expression template for derivatives of
75 // products
76 template <typename LhsTag, typename RhsTag>
77 struct OnDemandInputsForSwshJacobianImpl<
78  Tags::Dy<::Tags::Multiplies<LhsTag, RhsTag>>,
80  LhsTag::type::type::spin + RhsTag::type::type::spin>,
81  cpp17::bool_constant<not cpp17::is_same_v<LhsTag, Tags::BondiJbar> and
82  not cpp17::is_same_v<LhsTag, Tags::BondiUbar> and
83  not cpp17::is_same_v<RhsTag, Tags::BondiJbar>>> {
84  template <typename DataBoxTagList>
85  SPECTRE_ALWAYS_INLINE decltype(auto) operator()(
86  const db::DataBox<DataBoxTagList>& box) noexcept {
87  decltype(auto) lhs = get(db::get<LhsTag>(box)).data();
88  decltype(auto) dy_lhs = get(db::get<Tags::Dy<LhsTag>>(box)).data();
89  decltype(auto) rhs = get(db::get<RhsTag>(box)).data();
90  decltype(auto) dy_rhs = get(db::get<Tags::Dy<RhsTag>>(box)).data();
91  return lhs * dy_rhs + dy_lhs * rhs;
92  }
93 };
94 
95 // use the product rule and an explicit conjugate to provide an expression
96 // template for derivatives of products with `Tags::BondiJbar` as the right-hand
97 // operand
98 template <typename LhsTag>
99 struct OnDemandInputsForSwshJacobianImpl<
100  Tags::Dy<::Tags::Multiplies<LhsTag, Tags::BondiJbar>>,
101  std::integral_constant<int, LhsTag::type::type::spin - 2>, std::true_type> {
102  template <typename DataBoxTagList>
103  SPECTRE_ALWAYS_INLINE decltype(auto) operator()(
104  const db::DataBox<DataBoxTagList>& box) noexcept {
105  decltype(auto) lhs = get(get<LhsTag>(box)).data();
106  decltype(auto) dy_lhs = get(get<Tags::Dy<LhsTag>>(box)).data();
107  decltype(auto) jbar = conj(get(get<Tags::BondiJ>(box)).data());
108  decltype(auto) dy_jbar = conj(get(get<Tags::Dy<Tags::BondiJ>>(box)).data());
109  return lhs * dy_jbar + dy_lhs * jbar;
110  }
111 };
112 
113 // use the product rule and an explicit conjugate to provide an expression
114 // template for derivatives of products with `Tags::BondiJbar` as the left-hand
115 // operand.
116 template <typename RhsTag>
117 struct OnDemandInputsForSwshJacobianImpl<
118  Tags::Dy<::Tags::Multiplies<Tags::BondiJbar, RhsTag>>,
119  std::integral_constant<int, RhsTag::type::type::spin - 2>, std::true_type> {
120  template <typename DataBoxTagList>
121  SPECTRE_ALWAYS_INLINE decltype(auto) operator()(
122  const db::DataBox<DataBoxTagList>& box) noexcept {
123  decltype(auto) rhs = get(get<RhsTag>(box)).data();
124  decltype(auto) dy_rhs = get(get<Tags::Dy<RhsTag>>(box)).data();
125  decltype(auto) jbar = conj(get(get<Tags::BondiJ>(box)).data());
126  decltype(auto) dy_jbar = conj(get(get<Tags::Dy<Tags::BondiJ>>(box)).data());
127  return dy_jbar * rhs + jbar * dy_rhs;
128  }
129 };
130 
131 // use the product rule and an explicit conjugate to provide an expression
132 // template for derivatives of products with `Tags::BondiUbar` as the left-hand
133 // operand.
134 template <typename RhsTag>
135 struct OnDemandInputsForSwshJacobianImpl<
136  Tags::Dy<::Tags::Multiplies<Tags::BondiUbar, RhsTag>>,
137  std::integral_constant<int, RhsTag::type::type::spin - 1>, std::true_type> {
138  template <typename DataBoxTagList>
139  SPECTRE_ALWAYS_INLINE decltype(auto) operator()(
140  const db::DataBox<DataBoxTagList>& box) noexcept {
141  decltype(auto) ubar = conj(get(get<Tags::BondiU>(box)).data());
142  decltype(auto) dy_ubar = conj(get(get<Tags::Dy<Tags::BondiU>>(box)).data());
143  decltype(auto) rhs = get(get<RhsTag>(box)).data();
144  decltype(auto) dy_rhs = get(get<Tags::Dy<RhsTag>>(box)).data();
145  return ubar * dy_rhs + dy_ubar * rhs;
146  }
147 };
148 
149 // use the product rule and an explicit conjugate to provide an expression
150 // template for second derivatives of products with `Tags::BondiJbar` as the
151 // right-hand operand.
152 template <typename LhsTag>
153 struct OnDemandInputsForSwshJacobianImpl<
154  Tags::Dy<Tags::Dy<::Tags::Multiplies<LhsTag, Tags::BondiJbar>>>,
155  std::integral_constant<int, LhsTag::type::type::spin - 2>, std::true_type> {
156  template <typename DataBoxTagList>
157  SPECTRE_ALWAYS_INLINE decltype(auto) operator()(
158  const db::DataBox<DataBoxTagList>& box) noexcept {
159  decltype(auto) lhs = get(get<LhsTag>(box)).data();
160  decltype(auto) dy_lhs = get(get<Tags::Dy<LhsTag>>(box)).data();
161  decltype(auto) dy_dy_lhs = get(get<Tags::Dy<Tags::Dy<LhsTag>>>(box)).data();
162  decltype(auto) jbar = conj(get(get<Tags::BondiJ>(box)).data());
163  decltype(auto) dy_jbar = conj(get(get<Tags::Dy<Tags::BondiJ>>(box)).data());
164  decltype(auto) dy_dy_jbar =
165  conj(get(get<Tags::Dy<Tags::Dy<Tags::BondiJ>>>(box)).data());
166  return lhs * dy_dy_jbar + 2.0 * dy_lhs * dy_jbar + dy_dy_lhs * jbar;
167  }
168 };
169 
170 // default to extracting directly from the box for spin-weighted derivatives of
171 // radial (y) derivatives
172 template <typename Tag, typename DerivKind>
173 struct OnDemandInputsForSwshJacobianImpl<
174  Spectral::Swsh::Tags::Derivative<Tags::Dy<Tag>, DerivKind>,
175  std::integral_constant<
176  int, Spectral::Swsh::Tags::Derivative<Tags::Dy<Tag>, DerivKind>::spin>,
177  std::true_type> {
178  template <typename DataBoxTagList>
179  SPECTRE_ALWAYS_INLINE decltype(auto) operator()(
180  const db::DataBox<DataBoxTagList>& box) noexcept {
181  return get(get<Spectral::Swsh::Tags::Derivative<Tags::Dy<Tag>, DerivKind>>(
182  box))
183  .data();
184  }
185 };
186 
187 // compute the derivative of the `jbar * (q - 2 eth_beta)` using the product
188 // rule and the commutation rule for the partials with respect to y and the
189 // spin_weighted derivatives
190 template <>
191 struct OnDemandInputsForSwshJacobianImpl<Tags::Dy<Tags::JbarQMinus2EthBeta>,
192  std::integral_constant<int, -1>,
193  std::true_type> {
194  template <typename DataBoxTagList>
195  SPECTRE_ALWAYS_INLINE decltype(auto) operator()(
196  const db::DataBox<DataBoxTagList>& box) noexcept {
197  decltype(auto) dy_beta = get(get<Tags::Dy<Tags::BondiBeta>>(box)).data();
198  decltype(auto) dy_j = get(get<Tags::Dy<Tags::BondiJ>>(box)).data();
199  decltype(auto) dy_q = get(get<Tags::Dy<Tags::BondiQ>>(box)).data();
200  decltype(auto) eth_beta =
201  get(get<Spectral::Swsh::Tags::Derivative<Tags::BondiBeta,
202  Spectral::Swsh::Tags::Eth>>(
203  box))
204  .data();
205  decltype(auto) eth_dy_beta =
206  get(get<Spectral::Swsh::Tags::Derivative<Tags::Dy<Tags::BondiBeta>,
207  Spectral::Swsh::Tags::Eth>>(
208  box))
209  .data();
210  decltype(auto) eth_r_divided_by_r =
211  get(get<Tags::EthRDividedByR>(box)).data();
212  decltype(auto) j = get(get<Tags::BondiJ>(box)).data();
213  decltype(auto) q = get(get<Tags::BondiQ>(box)).data();
214  return conj(j) * dy_q + conj(dy_j) * q - 2.0 * conj(j) * eth_dy_beta -
215  2.0 * eth_beta * conj(dy_j) -
216  2.0 * conj(j) * eth_r_divided_by_r * dy_beta;
217  }
218 };
219 } // namespace detail
220 
221 /// Provide an expression template or reference to `Tag`, intended for
222 /// situations for which a repeated computation is more
223 /// desirable than storing a value in the \ref DataBoxGroup (e.g. for
224 /// conjugation and simple product rule expansion).
225 template <typename Tag>
226 using OnDemandInputsForSwshJacobian = detail::OnDemandInputsForSwshJacobianImpl<
227  Tag, std::integral_constant<int, Tag::type::type::spin>, std::true_type>;
228 
229 /*!
230  * \brief Performs a mutation to a spin-weighted spherical harmonic derivative
231  * value from the numerical coordinate (the spin-weighted derivative at
232  * fixed \f$y\f$) to the Bondi coordinates (the spin-weighted derivative at
233  * fixed \f$r\f$), inplace to the requested tag.
234  *
235  * \details This should be performed only once for each derivative evaluation
236  * for each tag, as a repeated inplace evaluation will compound and result in
237  * incorrect values in the \ref DataBoxGroup. This is compatible with acting as
238  * a mutation in `db::mutate_apply`.
239  * \note In each specialization, there is an additional type alias
240  * `on_demand_argument_tags` that contains tags that represent additional
241  * quantities to be passed as arguments that need not be in the \ref
242  * DataBoxGroup. These quantities are suggested to be evaluated by the 'on
243  * demand' mechanism provided by `Cce::OnDemandInputsForSwshJacobian`, which
244  * provides the additional quantities as blaze expression templates rather than
245  * unnecessarily caching intermediate results that aren't re-used.
246  */
247 template <typename DerivativeTag>
249 
250 /*!
251  * \brief Specialization for the spin-weighted derivative \f$\eth\f$.
252  *
253  * \details The implemented equation is:
254  *
255  * \f[ \eth F = \eth^\prime F - (1 - y) \frac{\eth R}{R} \partial_y F,
256  * \f]
257  *
258  * where \f$\eth\f$ is the derivative at constant Bondi radius \f$r\f$ and
259  * \f$\eth^\prime\f$ is the derivative at constant numerical radius \f$y\f$.
260  */
261 template <typename ArgumentTag>
262 struct ApplySwshJacobianInplace<
263  Spectral::Swsh::Tags::Derivative<ArgumentTag, Spectral::Swsh::Tags::Eth>> {
264  using pre_swsh_derivative_tags = tmpl::list<>;
265  using swsh_derivative_tags = tmpl::list<>;
266  using integration_independent_tags =
267  tmpl::list<Tags::OneMinusY, Tags::EthRDividedByR>;
268 
269  using return_tags = tmpl::list<
271  using argument_tags = tmpl::append<integration_independent_tags>;
272  using on_demand_argument_tags = tmpl::list<Tags::Dy<ArgumentTag>>;
273 
274  static constexpr int spin =
277  template <typename DyArgumentType>
278  static void apply(
280  eth_argument,
281  const Scalar<SpinWeighted<ComplexDataVector, 0>>& one_minus_y,
282  const Scalar<SpinWeighted<ComplexDataVector, 1>>& eth_r_divided_by_r,
283  const DyArgumentType& dy_argument) {
284  get(*eth_argument) -= get(one_minus_y) * get(eth_r_divided_by_r) *
285  SpinWeighted<DyArgumentType, spin - 1>{dy_argument};
286  }
287 };
288 
289 /*!
290  * \brief Specialization for the spin-weighted derivative \f$\bar{\eth}\f$.
291  *
292  * \details The implemented equation is:
293  *
294  * \f[
295  * \bar{\eth} F = \bar{\eth}^\prime F
296  * - (1 - y) \frac{\bar{\eth} R}{R} \partial_y F,
297  *\f]
298  *
299  * where \f$\bar{\eth}\f$ is the derivative at constant Bondi radius \f$r\f$ and
300  * \f$\bar{\eth}^\prime\f$ is the derivative at constant numerical radius
301  * \f$y\f$.
302  */
303 template <typename ArgumentTag>
304 struct ApplySwshJacobianInplace<Spectral::Swsh::Tags::Derivative<
305  ArgumentTag, Spectral::Swsh::Tags::Ethbar>> {
306  using pre_swsh_derivative_tags = tmpl::list<>;
307  using swsh_derivative_tags = tmpl::list<>;
308  using integration_independent_tags =
309  tmpl::list<Tags::OneMinusY, Tags::EthRDividedByR>;
310 
311  using return_tags = tmpl::list<Spectral::Swsh::Tags::Derivative<
312  ArgumentTag, Spectral::Swsh::Tags::Ethbar>>;
313  using argument_tags = tmpl::append<integration_independent_tags>;
314  using on_demand_argument_tags = tmpl::list<Tags::Dy<ArgumentTag>>;
315 
316  static constexpr int spin =
318  Spectral::Swsh::Tags::Ethbar>::spin;
319  template <typename DyArgumentType>
320  static void apply(
322  ethbar_argument,
323  const Scalar<SpinWeighted<ComplexDataVector, 0>>& one_minus_y,
324  const Scalar<SpinWeighted<ComplexDataVector, 1>>& eth_r_divided_by_r,
325  const DyArgumentType& dy_argument) {
326  get(*ethbar_argument) -=
327  get(one_minus_y) * conj(get(eth_r_divided_by_r)) *
329  }
330 };
331 
332 /*!
333  * \brief Specialization for the spin-weighted derivative \f$\eth \bar{\eth}\f$.
334  *
335  * \details The implemented equation is:
336  *
337  * \f[
338  * \eth \bar{\eth} F = \eth^\prime \bar{\eth}^\prime F
339  * - \frac{\eth R \bar{\eth} R}{R^2} (1 - y)^2 \partial_y^2 F
340  * - (1 - y)\left(\frac{\eth R}{R} \bar{\eth} \partial_y F
341  * + \frac{\bar{\eth} R}{R} \eth \partial_y F
342  * + \frac{\eth \bar\eth R}{R} \partial_y F\right),
343  * \f]
344  *
345  * where \f$\eth \bar{\eth}\f$ is the derivative at constant Bondi radius
346  * \f$r\f$ and \f$\eth^\prime \bar{\eth}^\prime\f$ is the derivative at constant
347  * numerical radius \f$y\f$.
348  */
349 template <typename ArgumentTag>
350 struct ApplySwshJacobianInplace<Spectral::Swsh::Tags::Derivative<
351  ArgumentTag, Spectral::Swsh::Tags::EthEthbar>> {
352  using pre_swsh_derivative_tags = tmpl::list<>;
353  using swsh_derivative_tags = tmpl::list<>;
354  using integration_independent_tags =
357 
358  using return_tags = tmpl::list<Spectral::Swsh::Tags::Derivative<
359  ArgumentTag, Spectral::Swsh::Tags::EthEthbar>>;
360  using argument_tags = tmpl::append<integration_independent_tags>;
361  using on_demand_argument_tags =
362  tmpl::list<Tags::Dy<ArgumentTag>, Tags::Dy<Tags::Dy<ArgumentTag>>,
367 
368  static constexpr int spin =
370  Spectral::Swsh::Tags::EthEthbar>::spin;
371  template <typename DyArgumentType, typename DyDyArgumentType,
372  typename EthDyArgumentType, typename EthbarDyArgumentType>
373  static void apply(
375  eth_ethbar_argument,
376  const Scalar<SpinWeighted<ComplexDataVector, 0>>& one_minus_y,
377  const Scalar<SpinWeighted<ComplexDataVector, 1>>& eth_r_divided_by_r,
379  eth_ethbar_r_divided_by_r,
380  const DyArgumentType& dy_argument, const DyDyArgumentType& dy_dy_argument,
381  const EthDyArgumentType& eth_dy_argument,
382  const EthbarDyArgumentType ethbar_dy_argument) {
383  get(*eth_ethbar_argument) -=
384  get(eth_r_divided_by_r) * conj(get(eth_r_divided_by_r)) *
385  (square(get(one_minus_y)) *
386  SpinWeighted<DyDyArgumentType, spin>{dy_dy_argument}) +
387  get(one_minus_y) *
388  (get(eth_r_divided_by_r) *
390  ethbar_dy_argument} +
391  conj(get(eth_r_divided_by_r)) *
393  get(eth_ethbar_r_divided_by_r) *
395  }
396 };
397 
398 /*!
399  * \brief Specialization for the spin-weighted derivative \f$\bar{\eth} \eth\f$.
400  *
401  * \details The implemented equation is:
402  *
403  * \f[
404  * \bar{\eth} \eth F = \bar{\eth}^\prime \eth^\prime F
405  * - \frac{\eth R \bar{\eth} R}{R^2} (1 - y)^2 \partial_y^2 F
406  * - (1 - y)\left(\frac{\eth R}{R} \bar{\eth} \partial_y F
407  * + \frac{\bar{\eth} R}{R} \eth \partial_y F
408  * + \frac{\eth \bar\eth R}{R} \partial_y F\right),
409  * \f]
410  *
411  * where \f$\bar{\eth} \eth\f$ is the derivative at constant Bondi radius
412  * \f$r\f$ and \f$\bar{\eth}^\prime \eth^\prime\f$ is the derivative at constant
413  * numerical radius \f$y\f$.
414  */
415 template <typename ArgumentTag>
416 struct ApplySwshJacobianInplace<Spectral::Swsh::Tags::Derivative<
417  ArgumentTag, Spectral::Swsh::Tags::EthbarEth>> {
418  using pre_swsh_derivative_tags = tmpl::list<>;
419  using swsh_derivative_tags = tmpl::list<>;
420  using integration_independent_tags =
423 
424  using return_tags = tmpl::list<Spectral::Swsh::Tags::Derivative<
425  ArgumentTag, Spectral::Swsh::Tags::EthbarEth>>;
426  using argument_tags = tmpl::append<integration_independent_tags>;
427  using on_demand_argument_tags =
428  tmpl::list<Tags::Dy<ArgumentTag>, Tags::Dy<Tags::Dy<ArgumentTag>>,
433 
434  static constexpr int spin =
436  Spectral::Swsh::Tags::EthbarEth>::spin;
437  template <typename DyArgumentType, typename DyDyArgumentType,
438  typename EthDyArgumentType, typename EthbarDyArgumentType>
439  static void apply(
441  ethbar_eth_argument,
442  const Scalar<SpinWeighted<ComplexDataVector, 0>>& one_minus_y,
443  const Scalar<SpinWeighted<ComplexDataVector, 1>>& eth_r_divided_by_r,
445  eth_ethbar_r_divided_by_r,
446  const DyArgumentType& dy_argument, const DyDyArgumentType& dy_dy_argument,
447  const EthDyArgumentType& eth_dy_argument,
448  const EthbarDyArgumentType ethbar_dy_argument) {
449  get(*ethbar_eth_argument) -=
450  get(eth_r_divided_by_r) * conj(get(eth_r_divided_by_r)) *
451  (square(get(one_minus_y)) *
452  SpinWeighted<DyDyArgumentType, spin>{dy_dy_argument}) +
453  get(one_minus_y) *
454  (get(eth_r_divided_by_r) *
456  ethbar_dy_argument} +
457  conj(get(eth_r_divided_by_r)) *
459  get(eth_ethbar_r_divided_by_r) *
461  }
462 };
463 
464 /*!
465  * \brief Specialization for the spin-weighted derivative \f$\eth \eth\f$.
466  *
467  * \details The implemented equation is:
468  *
469  * \f[
470  * \eth \eth F = \eth^\prime \eth^\prime F
471  * - (1 - y)^2 \frac{(\eth R)^2}{R^2} \partial_y^2 F
472  * - (1 - y) \left( 2 \frac{\eth R}{R} \eth \partial_y F
473  * + \frac{\eth \eth R}{R} \partial_y F\right),
474  * \f]
475  *
476  * where \f$\eth \eth\f$ is the derivative at constant Bondi radius \f$r\f$ and
477  * \f$\eth^\prime \eth^\prime\f$ is the derivative at constant numerical radius
478  * \f$y\f$.
479  */
480 template <typename ArgumentTag>
481 struct ApplySwshJacobianInplace<Spectral::Swsh::Tags::Derivative<
482  ArgumentTag, Spectral::Swsh::Tags::EthEth>> {
483  using pre_swsh_derivative_tags = tmpl::list<>;
484  using swsh_derivative_tags = tmpl::list<>;
485  using integration_independent_tags =
488 
489  using return_tags = tmpl::list<Spectral::Swsh::Tags::Derivative<
490  ArgumentTag, Spectral::Swsh::Tags::EthEth>>;
491  using argument_tags = tmpl::append<integration_independent_tags>;
492  using on_demand_argument_tags =
493  tmpl::list<Tags::Dy<ArgumentTag>, Tags::Dy<Tags::Dy<ArgumentTag>>,
496 
497  static constexpr int spin =
499  Spectral::Swsh::Tags::EthEth>::spin;
500  template <typename DyArgumentType, typename DyDyArgumentType,
501  typename EthDyArgumentType>
502  static void apply(
504  eth_eth_argument,
505  const Scalar<SpinWeighted<ComplexDataVector, 0>>& one_minus_y,
506  const Scalar<SpinWeighted<ComplexDataVector, 1>>& eth_r_divided_by_r,
507  const Scalar<SpinWeighted<ComplexDataVector, 2>>& eth_eth_r_divided_by_r,
508  const DyArgumentType& dy_argument, const DyDyArgumentType& dy_dy_argument,
509  const EthDyArgumentType& eth_dy_argument) {
510  get(*eth_eth_argument) -=
511  square(get(eth_r_divided_by_r)) *
512  (square(get(one_minus_y)) *
513  SpinWeighted<DyDyArgumentType, spin - 2>{dy_dy_argument}) +
514  get(one_minus_y) *
515  (2.0 * get(eth_r_divided_by_r) *
517  get(eth_eth_r_divided_by_r) *
519  }
520 };
521 
522 /*!
523  * \brief Specialization for the spin-weighted derivative \f$\bar{\eth}
524  * \bar{\eth}\f$.
525  *
526  * \details The implemented equation is:
527  *
528  * \f[
529  * \bar{\eth} \bar{\eth} F = \bar{\eth}^\prime \bar{\eth}^\prime F
530  * - (1 - y)^2 \frac{(\bar{\eth} R)^2}{R^2} \partial_y^2 F
531  * - (1 - y) \left( 2 \frac{\bar{\eth} R}{R} \bar{\eth} \partial_y F
532  * + \frac{\bar{\eth} \bar{\eth} R}{R} \partial_y F\right),
533  * \f]
534  *
535  * where \f$\bar{\eth} \bar{\eth}\f$ is the derivative at constant Bondi radius
536  * \f$r\f$ and \f$\bar{\eth}^\prime \bar{\eth}^\prime\f$ is the derivative at
537  * constant numerical radius \f$y\f$.
538  */
539 template <typename ArgumentTag>
540 struct ApplySwshJacobianInplace<Spectral::Swsh::Tags::Derivative<
541  ArgumentTag, Spectral::Swsh::Tags::EthbarEthbar>> {
542  using pre_swsh_derivative_tags = tmpl::list<>;
543  using swsh_derivative_tags = tmpl::list<>;
544  using integration_independent_tags =
547 
548  using return_tags = tmpl::list<Spectral::Swsh::Tags::Derivative<
549  ArgumentTag, Spectral::Swsh::Tags::EthbarEthbar>>;
550  using argument_tags = tmpl::append<integration_independent_tags>;
551  using on_demand_argument_tags =
552  tmpl::list<Tags::Dy<ArgumentTag>, Tags::Dy<Tags::Dy<ArgumentTag>>,
555 
556  static constexpr int spin = Spectral::Swsh::Tags::Derivative<
557  ArgumentTag, Spectral::Swsh::Tags::EthbarEthbar>::spin;
558  template <typename DyArgumentType, typename DyDyArgumentType,
559  typename EthbarDyArgumentType>
560  static void apply(
562  ethbar_ethbar_argument,
563  const Scalar<SpinWeighted<ComplexDataVector, 0>>& one_minus_y,
564  const Scalar<SpinWeighted<ComplexDataVector, 1>>& eth_r_divided_by_r,
565  const Scalar<SpinWeighted<ComplexDataVector, 2>>& eth_eth_r_divided_by_r,
566  const DyArgumentType& dy_argument, const DyDyArgumentType& dy_dy_argument,
567  const EthbarDyArgumentType& ethbar_dy_argument) {
568  get(*ethbar_ethbar_argument) -=
569  square(conj(get(eth_r_divided_by_r))) *
570  (square(get(one_minus_y)) *
572  get(one_minus_y) *
573  (2.0 * conj(get(eth_r_divided_by_r)) *
575  ethbar_dy_argument} +
576  conj(get(eth_eth_r_divided_by_r)) *
578  }
579 };
580 
581 namespace detail {
582 // A helper to forward to the `ApplySwshJacobianInplace` mutators that takes
583 // advantage of the `OnDemandInputsForSwshJacobian`, which computes blaze
584 // template expressions for those quantities for which it is anticipated to be
585 // sufficiently cheap to repeatedly compute that it is worth saving the cost of
586 // additional storage (e.g. conjugates and derivatives for which the product
587 // rule applies)
588 template <typename DerivativeTag, typename DataBoxTagList,
589  typename... OnDemandTags>
590 void apply_swsh_jacobian_helper(
592  tmpl::list<OnDemandTags...> /*meta*/) noexcept {
593  db::mutate_apply<ApplySwshJacobianInplace<DerivativeTag>>(
595 }
596 } // namespace detail
597 
598 /*!
599  * \brief This routine evaluates the set of inputs to the CCE integrand for
600  * `BondiValueTag` which are spin-weighted angular derivatives.
601 
602  * \details This function is called on the \ref DataBoxGroup holding the
603  * relevant CCE data during each hypersurface integration step, after evaluating
604  * `mutate_all_pre_swsh_derivatives_for_tag()` with template argument
605  * `BondiValueTag` and before evaluating `ComputeBondiIntegrand<BondiValueTag>`.
606  * Provided a \ref DataBoxGroup with the appropriate tags (including
607  * `Cce::all_pre_swsh_derivative_tags`, `Cce::all_swsh_derivative_tags`,
608  * `Cce::all_transform_buffer_tags`, `Cce::pre_computation_tags`, and
609  * `Cce::Tags::LMax`), this function will apply all of the necessary
610  * mutations to update
611  * `Cce::single_swsh_derivative_tags_to_compute_for<BondiValueTag>` and
612  * `Cce::second_swsh_derivative_tags_to_compute_for<BondiValueTag>` to their
613  * correct values for the current values of the remaining (input) tags.
614  */
615 template <typename BondiValueTag, typename DataBoxTagList>
617  const gsl::not_null<db::DataBox<DataBoxTagList>*> box) noexcept {
618  // The collection of spin-weighted derivatives cannot be applied as individual
619  // compute items, because it is better to aggregate similar spins and dispatch
620  // to libsharp in groups. So, we supply a bulk mutate operation which takes in
621  // multiple Variables from the presumed DataBox, and alters their values as
622  // necessary.
624  single_swsh_derivative_tags_to_compute_for_t<BondiValueTag>>>(box);
625  tmpl::for_each<single_swsh_derivative_tags_to_compute_for_t<BondiValueTag>>(
626  [&box](auto derivative_tag_v) noexcept {
627  using derivative_tag = typename decltype(derivative_tag_v)::type;
628  detail::apply_swsh_jacobian_helper<derivative_tag>(
629  box, typename ApplySwshJacobianInplace<
630  derivative_tag>::on_demand_argument_tags{});
631  });
632 
634  second_swsh_derivative_tags_to_compute_for_t<BondiValueTag>>>(box);
635  tmpl::for_each<second_swsh_derivative_tags_to_compute_for_t<BondiValueTag>>(
636  [&box](auto derivative_tag_v) noexcept {
637  using derivative_tag = typename decltype(derivative_tag_v)::type;
638  detail::apply_swsh_jacobian_helper<derivative_tag>(
639  box, typename ApplySwshJacobianInplace<
640  derivative_tag>::on_demand_argument_tags{});
641  });
642 }
643 } // namespace Cce
detail::OnDemandInputsForSwshJacobianImpl< Tag, std::integral_constant< int, Tag::type::type::spin >, std::true_type > OnDemandInputsForSwshJacobian
Provide an expression template or reference to Tag, intended for situations for which a repeated comp...
Definition: SwshDerivatives.hpp:227
The set of utilities for performing Cauchy characteristic evolution and Cauchy characteristic matchin...
Definition: BoundaryComputeAndSendToEvolution.hpp:28
The value , where is Bondi radius of the worldtube.
Definition: Tags.hpp:312
Struct for labeling the spin-weighted derivative in tags.
Definition: SwshTags.hpp:49
constexpr Tag::type & get(Variables< TagList > &v) noexcept
Return Tag::type pointing into the contiguous array.
Definition: Variables.hpp:688
void mutate_all_swsh_derivatives_for_tag(const gsl::not_null< db::DataBox< DataBoxTagList > *> box) noexcept
This routine evaluates the set of inputs to the CCE integrand for BondiValueTag which are spin-weight...
Definition: SwshDerivatives.hpp:616
Struct for labeling the spin-weighted derivative in tags.
Definition: SwshTags.hpp:41
Prefix tag representing the spin-weighted derivative of a spin-weighted scalar.
Definition: SwshTags.hpp:152
detail::AngularDerivativesImpl< DerivativeTagList, typename detail::unique_derived_from_list< DerivativeTagList >::type, Representation > AngularDerivatives
A DataBox mutate-compatible computational struct for computing a set of spin-weighted spherical harmo...
Definition: SwshDerivatives.hpp:340
Coordinate value , which will be cached and sent to the implementing functions.
Definition: Tags.hpp:277
Performs a mutation to a spin-weighted spherical harmonic derivative value from the numerical coordin...
Definition: SwshDerivatives.hpp:248
Definition: Determinant.hpp:11
Make a spin-weighted type T with spin-weight Spin. Mathematical operators are restricted to addition...
Definition: SpinWeighted.hpp:24
constexpr auto mutate_apply(F &&f, const gsl::not_null< DataBox< BoxTags > *> box, Args &&... args) noexcept(DataBox_detail::check_mutate_apply_mutate_tags(BoxTags{}, MutateTags{}) and DataBox_detail::check_mutate_apply_argument_tags(BoxTags{}, ArgumentTags{}) and noexcept(DataBox_detail::mutate_apply(f, box, MutateTags{}, ArgumentTags{}, std::forward< Args >(args)...)))
Apply the invokable f mutating items MutateTags and taking as additional arguments ArgumentTags and a...
Definition: DataBox.hpp:1662
#define SPECTRE_ALWAYS_INLINE
Always inline a function. Only use this if you benchmarked the code.
Definition: ForceInline.hpp:16
Defines classes and functions used for manipulating DataBox&#39;s.
The value , where is Bondi radius of the worldtube.
Definition: Tags.hpp:296
The value , where is Bondi radius of the worldtube.
Definition: Tags.hpp:304
Definition: InterpolationTargetWedgeSectionTorus.hpp:25
Defines macro to always inline a function.
Definition: DataBoxTag.hpp:27
Defines the class template Mesh.
The derivative with respect to the numerical coordinate , where is Bondi radius of the worldtube...
Definition: Tags.hpp:115
Defines classes for Tensor.
Defines a list of useful type aliases for tensors.
constexpr auto apply(F &&f, const DataBox< BoxTags > &box, Args &&... args) noexcept
Apply the invokable f with argument Tags TagsList from DataBox box
Definition: DataBox.hpp:1444
Struct for labeling the spin-weighted derivative in tags.
Definition: SwshTags.hpp:30
Namespace for DataBox related things.
Definition: DataBox.hpp:42
Defines functions and classes from the GSL.
Struct for labeling the spin-weighted derivative in tags.
Definition: SwshTags.hpp:37
Functionality associated with a particular choice of basis functions and quadrature for spectral oper...
Definition: Chebyshev.cpp:19
Struct for labeling the spin-weighted derivative in tags.
Definition: SwshTags.hpp:45
decltype(auto) constexpr square(const T &x)
Compute the square of x
Definition: ConstantExpressions.hpp:54
Tensor< T, Symmetry<>, index_list<> > Scalar
Scalar type.
Definition: TypeAliases.hpp:21
Struct for labeling the spin-weighted derivative in tags.
Definition: SwshTags.hpp:34
Require a pointer to not be a nullptr
Definition: Gsl.hpp:182