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