|
SpECTRE
v2024.04.12
|
Performs interpolation for spin-weighted spherical harmonics by taking advantage of the Clenshaw method of expanding recurrence relations. More...
#include <SwshInterpolation.hpp>
Public Member Functions | |
| SwshInterpolator (const SwshInterpolator &)=default | |
| SwshInterpolator (SwshInterpolator &&)=default | |
| SwshInterpolator & | operator= (const SwshInterpolator &)=default |
| SwshInterpolator & | operator= (SwshInterpolator &&)=default |
| SwshInterpolator (const DataVector &theta, const DataVector &phi, size_t l_max) | |
| template<int Spin> | |
| void | interpolate (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > interpolated, const SpinWeighted< ComplexModalVector, Spin > &goldberg_modes) const |
Perform the Clenshaw recurrence sum, returning by pointer interpolated of interpolating the goldberg_modes at the collocation points passed to the constructor. More... | |
| template<int Spin> | |
| void | interpolate (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > interpolated, const SpinWeighted< ComplexDataVector, Spin > &libsharp_collocation) const |
Perform the Clenshaw recurrence sum, returning by pointer interpolated of interpolating function represented by libsharp_collocation at the target points passed to the constructor. More... | |
| template<int Spin> | |
| void | direct_evaluation_swsh_at_l_min (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > harmonic, int m) const |
| Evaluate the SWSH function at the lowest \(l\) value for a given \(m\) at the target interpolation points. More... | |
| template<int Spin> | |
| void | evaluate_swsh_at_l_min_plus_one (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > harmonic, const SpinWeighted< ComplexDataVector, Spin > &harmonic_at_l_min, int m) const |
| Evaluate the SWSH function at the next-to-lowest \(l\) value for a given \(m\) at the target interpolation points, given input harmonic values for the lowest \(l\) value. More... | |
| template<int Spin> | |
| void | evaluate_swsh_m_recurrence_at_l_min (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > harmonic, int m) const |
| Evaluate the SWSH function at the lowest \(l\) value for a given \(m\) at the target interpolation points, given harmonic data at the next lower \(m\) (by magnitude), passed in by the same pointer used for the return. More... | |
| template<int Spin> | |
| void | clenshaw_sum (gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > interpolation, const SpinWeighted< ComplexDataVector, Spin > &l_min_harmonic, const SpinWeighted< ComplexDataVector, Spin > &l_min_plus_one_harmonic, const SpinWeighted< ComplexModalVector, Spin > &goldberg_modes, int m) const |
Perform the core Clenshaw interpolation at fixed \(m\), accumulating the result in interpolation. More... | |
| void | pup (PUP::er &p) |
| Serialization for Charm++. | |
Performs interpolation for spin-weighted spherical harmonics by taking advantage of the Clenshaw method of expanding recurrence relations.
During construction, we cache several functions of the target interpolation points that will be used during the Clenshaw evaluation. A new SwshInterpolator object must be created for each new set of target points, but the member function SwshInterpolator::interpolate() may be called on several different coefficients or collocation sets, and of different spin-weights.
This utility obtains the Clenshaw interpolation constants from a StaticCache, so that 'universal' quantities can be calculated only once per execution and re-used on each interpolation.
We evaluate the recurrence coefficients \(\alpha_l^{(a,b)}\) and \(\beta_l^{(a,b)}\), where \(a = |s + m|\), \(b = |s - m|\), and
\[ {}_sY_{l, m}(\theta, \phi) = \alpha_l^{(a,b)}(\theta) {}_s Y_{l-1, m} + \beta_l^{(a, b)} {}_s Y_{l - 2, m}(\theta, \phi) \]
The core Clenshaw recurrence is in the \(l\) modes of the spin-weighted spherical harmonic. For a set of modes \(a_l^{(a,b)}\) at fixed \(m\), the function value is evaluated by the recurrence:
\begin{align*} y^{(a, b)}_{l_\text{max} + 2} &= y^{(a, b)}_{l_\text{max} + 1} = 0 \\ y^{(a, b)}_{l}(\theta) &= \alpha_{l + 1}^{(a, b)} y^{(a, b)}_{l + 1}(\theta) + \beta_{l + 2}^{(a,b)} y^{(a, b)}_{l + 2}(\theta) + a_l^{(a, b)} \\ f_m(\theta, \phi) &= \beta_{l_{\text{min}} + 2} {}_s Y_{l_{\text{min}}, m}(\theta, \phi) y^{(a, b)}_2(\theta) + {}_s Y_{l_{\text{min}} + 1, m}(\theta, \phi) y^{(a, b)}_1(\theta) + a_0^{(a, b)} {}_s Y_{l_{\text{min}}, m}(\theta, \phi), \end{align*}
where \(l_{\text{min}} = \max(|m|, |s|)\) is the lowest nonvanishing \(l\) mode for a given \(m\) and \(s\).
The coefficients to cache are inferred from a mechanical but lengthy calculation involving the recurrence relation for the Jacobi polynomials. The result is:
\begin{align*} \alpha_l^{(a, b)} &= \frac{\sqrt{(2l + 1)(2l - 1)}} {2\sqrt{(l + k)(l + k + a + b)(l + k + a)(l + k + b)}} \left[2 l \cos(\theta) + \frac{a^2 - b^2}{2 l - 2}\right] \\ \beta_l^{(a, b)} & = - \sqrt{\frac{(2 l + 1)(l + k + a - 1)(l + k + b - 1) (l + k - 1)(l + k + a + b - 1)} {(2l - 3)(l + k)(l + k + a + b)(l + k + a)(l + k + b)}} \frac{2 l}{2 l - 2}, \end{align*}
where \(k = - (a + b)/2\) (which is always integral due to the properties of \(a\) and \(b\)). Note that because the values of \(\alpha\) and \(\beta\) in the recurrence relation are not needed for any value below \(l_{\text{min}} + 2\), so none of the values in the square-roots or denominators take pathological values for any of the coefficients we require.
The \(\beta\) constants are filled in member variable beta_constant. The \(\alpha\) values are separately stored as the prefactor for the \(\cos(\theta)\) term and the constant term in alpha_prefactor and alpha_constant member variables.
In addition, it is efficient to cache recurrence coefficients necessary for generating the first couple of spin-weighted spherical harmonic functions for each \(m\) used in the Clenshaw sum.
The member variable harmonic_at_l_min_prefactors holds the prefactors for directly evaluating the harmonics at \(s >= m\),
\begin{align*} {}_s Y_{|s| m} = {}_s C_{m} e^{i m \phi} \sin^a(\theta/2) \cos^b(\theta/2), \end{align*}
where \({}_s C_m\) are the cached prefactors that take the values
\begin{align*} {}_s C_m = (-1)^{m + \lambda(m)} \sqrt{\frac{(2 |s| + 1) (|s| + k)!} {4 \pi (|s| + k + a)! (|s| + k + b)!}} \end{align*}
and
\begin{align*} \lambda(m) = \left\{ \begin{array}{ll} 0 &\text{ for } s \ge -m\\ s + m &\text{ for } s < -m \end{array} \right.. \end{align*}
The member variable harmonic_m_recurrence_prefactors holds the prefactors necessary to evaluate the lowest harmonics for each \(m\) from the next lowest-in-magnitude \(m\), allowing most leading harmonics to recursively derived.
\begin{align*} {}_s Y_{|m|, m} = {}_s D_m \sin(\theta / 2) \cos(\theta / 2) \left\{ \begin{array}{ll} e^{i \phi} {}_s Y_{|m| - 1, m - 1} &\text{ for } m > 0\\ e^{-i \phi} {}_s Y_{|m| - 1, m + 1} &\text{ for } m < 0 \end{array}\right., \end{align*}
where \({}_s D_m\) are the cached prefactors that take the values
\begin{align*} {}_s D_m = -(-1)^{\Delta \lambda(m)} \sqrt{ \frac{(2 |m| + 1)(k + |m| + a + b - 1)(k + |m| + a + b)} {(2 |m| - 1)(k + |m| + a)(k + |m| + b)}}, \end{align*}
and \(\Delta \lambda(m)\) is the difference in \(\lambda(m)\) between the harmonic on the left-hand side and right-hand side of the above recurrence relation, that is \(\lambda(m) - \lambda(m - 1)\) for positive \(m\) and \(\lambda(m) - \lambda(m + 1)\) for negative \(m\).
Finally, the member variable harmonic_at_l_min_plus_one_recurrence_prefactors holds the prefactors necessary to evaluate parts of the recurrence relations from the lowest \(l\) for a given \(m\) to the next-to-lowest \(l\).
\begin{align*} {}_s Y_{l_{\min} + 1, m} = {}_s E_m \left[(a + 1) + (a + b + 2) \frac{(\cos(\theta) - 1)}{2}\right] {}_s Y_{l_{\min}, m}. \end{align*}
where \({}_s E_m\) are the cached prefactors that take the values
\begin{align*} {}_s E_{m} = \sqrt{\frac{2 l_{\min} + 3}{ 2 l_{\min} + 1}} \sqrt{\frac{(l_{\min} + k + 1)(l_{\min} + k + a + b + 1)} {(l_{\min} + k + a + 1)(l_{\min} + k + b + 1)}} \end{align*}
| void Spectral::Swsh::SwshInterpolator::clenshaw_sum | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | interpolation, |
| const SpinWeighted< ComplexDataVector, Spin > & | l_min_harmonic, | ||
| const SpinWeighted< ComplexDataVector, Spin > & | l_min_plus_one_harmonic, | ||
| const SpinWeighted< ComplexModalVector, Spin > & | goldberg_modes, | ||
| int | m | ||
| ) | const |
Perform the core Clenshaw interpolation at fixed \(m\), accumulating the result in interpolation.
Included in the public interface for thorough testing, most use cases should just use the interpolate member function.
| void Spectral::Swsh::SwshInterpolator::direct_evaluation_swsh_at_l_min | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | harmonic, |
| int | m | ||
| ) | const |
Evaluate the SWSH function at the lowest \(l\) value for a given \(m\) at the target interpolation points.
Included in the public interface for thorough testing, most use cases should just use the SwshInterpolator::interpolate() member function.
| void Spectral::Swsh::SwshInterpolator::evaluate_swsh_at_l_min_plus_one | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | harmonic, |
| const SpinWeighted< ComplexDataVector, Spin > & | harmonic_at_l_min, | ||
| int | m | ||
| ) | const |
Evaluate the SWSH function at the next-to-lowest \(l\) value for a given \(m\) at the target interpolation points, given input harmonic values for the lowest \(l\) value.
Included in the public interface for thorough testing, most use cases should just use the SwshInterpolator::interpolate() member function.
| void Spectral::Swsh::SwshInterpolator::evaluate_swsh_m_recurrence_at_l_min | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | harmonic, |
| int | m | ||
| ) | const |
Evaluate the SWSH function at the lowest \(l\) value for a given \(m\) at the target interpolation points, given harmonic data at the next lower \(m\) (by magnitude), passed in by the same pointer used for the return.
Included in the public interface for thorough testing, most use cases should just use the SwshInterpolator::interpolate() member function.
| void Spectral::Swsh::SwshInterpolator::interpolate | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | interpolated, |
| const SpinWeighted< ComplexDataVector, Spin > & | libsharp_collocation | ||
| ) | const |
Perform the Clenshaw recurrence sum, returning by pointer interpolated of interpolating function represented by libsharp_collocation at the target points passed to the constructor.
The core Clenshaw recurrence is in the \(l\) modes of the spin-weighted spherical harmonic. For a set of modes \(a_l^{(a,b)}\) at fixed \(m\), the function value is evaluated by the recurrence:
\begin{align*} y^{(a, b)}_{l_\text{max} + 2} &= y^{(a, b)}_{l_\text{max} + 1} = 0 \\ y^{(a, b)}_{l}(\theta) &= \alpha_{l + 1}^{(a, b)} y^{(a, b)}_{l + 1}(\theta) + \beta_{l + 2}^{(a,b)} y^{(a, b)}_{l + 2}(\theta) + a_l^{(a, b)} \\ f_m(\theta, \phi) &= \beta_{l_{\text{min}} + 2} {}_s Y_{l_{\text{min}}, m}(\theta, \phi) y^{(a, b)}_2(\theta) + {}_s Y_{l_{\text{min}} + 1, m}(\theta, \phi) y^{(a, b)}_1(\theta) + a_0^{(a, b)} {}_s Y_{l_{\text{min}}, m}(\theta, \phi) \end{align*}
The recurrence in \(l\) accomplishes much of the work, but for full efficiency, we also recursively evaluate the lowest handful of \(l\)s for each \(m\). The details of those additional recurrence tricks can be found in the documentation for ClenshawRecurrenceConstants.
| void Spectral::Swsh::SwshInterpolator::interpolate | ( | gsl::not_null< SpinWeighted< ComplexDataVector, Spin > * > | interpolated, |
| const SpinWeighted< ComplexModalVector, Spin > & | goldberg_modes | ||
| ) | const |
Perform the Clenshaw recurrence sum, returning by pointer interpolated of interpolating the goldberg_modes at the collocation points passed to the constructor.
The core Clenshaw recurrence is in the \(l\) modes of the spin-weighted spherical harmonic. For a set of modes \(a_l^{(a,b)}\) at fixed \(m\), the function value is evaluated by the recurrence:
\begin{align*} y^{(a, b)}_{l_\text{max} + 2} &= y^{(a, b)}_{l_\text{max} + 1} = 0 \\ y^{(a, b)}_{l}(\theta) &= \alpha_{l + 1}^{(a, b)} y^{(a, b)}_{l + 1}(\theta) + \beta_{l + 2}^{(a,b)} y^{(a, b)}_{l + 2}(\theta) + a_l^{(a, b)} \\ f_m(\theta, \phi) &= \beta_{l_{\text{min}} + 2} {}_s Y_{l_{\text{min}}, m}(\theta, \phi) y^{(a, b)}_2(\theta) + {}_s Y_{l_{\text{min}} + 1, m}(\theta, \phi) y^{(a, b)}_1(\theta) + a_0^{(a, b)} {}_s Y_{l_{\text{min}}, m}(\theta, \phi) \end{align*}
The recurrence in \(l\) accomplishes much of the work, but for full efficiency, we also recursively evaluate the lowest handful of \(l\)s for each \(m\). The details of those additional recurrence tricks can be found in the documentation for ClenshawRecurrenceConstants.