Spectral

Namespaces
namespace	Spectral
	Functionality associated with a particular choice of basis functions and quadrature for spectral operations.
namespace	Spectral::Swsh
	Namespace for spin-weighted spherical harmonic utilities.
namespace	Spectral::filtering
	Matrices for filtering spectral coefficients.

Classes
class	Spectral::Chebyshev
	A collection of helper functions for Chebyshev polynomials. More...
class	Spectral::Fourier
	A collection of helper functions for using a Fourier series. More...
class	Spectral::HalfFourier
	A collection of helper functions for the half-Fourier spectral basis. More...
class	Spectral::Jacobi
	A collection of helper functions for Jacobi polyomials. More...
class	Spectral::Zernike< Dim >
	A collection of helper functions for the radial functions used in Zernike polyomials. More...
class	ylm::Spherepack
	Defines the C++ interface to SPHEREPACK. More...
class	ylm::SpherepackIterator
	Iterates over spectral coefficients stored in SPHEREPACK format. More...

Typedefs
template<typename DerivativeTagList, ComplexRepresentation Representation = ComplexRepresentation::Interleaved>
using	Spectral::Swsh::AngularDerivatives
	A DataBox mutate-compatible computational struct for computing a set of spin-weighted spherical harmonic derivatives by grouping and batch-computing spin-weighted spherical harmonic transforms.

Functions
ConvergenceInfo	PowerMonitors::convergence_rate_and_number_of_pile_up_modes (const DataVector &power_monitor, size_t number_of_filtered_modes=0)
	Returns the convergence rate and the number of pile up modes of a power monitor as a ConvergenceInfo.
void	Spectral::zernike_b2_disk_nodal_to_modal (gsl::not_null< DataVector * > modal, gsl::not_null< DataVector * > buf, const DataVector &u, size_t n_r, size_t n_phi, size_t n_r_max, size_t num_components=1)
	Transform nodal ZernikeB2 data on a disk to modal space.
template<int Spin>
auto	Spectral::Swsh::swsh_buffer (const size_t l_max, const size_t number_of_radial_points)
	Produces a SpinWeighted<ComplexModalVector, Spin> of the appropriate size to be used as a modal buffer for Spectral::Swsh::AngularDerivatives or Spectral::Swsh::angular_derivatives.
template<typename DerivativeKindList, ComplexRepresentation Representation = ComplexRepresentation::Interleaved, typename... ArgumentTypes>
void	Spectral::Swsh::angular_derivatives (const size_t l_max, const size_t number_of_radial_points, const ArgumentTypes &... arguments)
	Evaluate all of the spin-weighted derivatives in DerivKindList on input SpinWeighted<ComplexDataVector, Spin> collocation data, returning by pointer.
template<typename DerivKind, ComplexRepresentation Representation = ComplexRepresentation::Interleaved, int Spin>
SpinWeighted< ComplexDataVector, Tags::derivative_spin_weight< DerivKind >+Spin >	Spectral::Swsh::angular_derivative (size_t l_max, size_t number_of_radial_points, const SpinWeighted< ComplexDataVector, Spin > &to_differentiate)
	Evaluate the spin-weighted derivative DerivKind on the provided SpinWeighted<ComplexDataVector, Spin> collocation data, returning by value.
template<size_t Dim>
void	to_modal_coefficients (gsl::not_null< ComplexModalVector * > modal_coefficients, const ComplexDataVector &nodal_coefficients, const Mesh< Dim > &mesh)
	Compute the modal coefficients from the nodal coefficients.
template<size_t Dim>
void	to_modal_coefficients (gsl::not_null< ModalVector * > modal_coefficients, const DataVector &nodal_coefficients, const Mesh< Dim > &mesh)
	Compute the modal coefficients from the nodal coefficients.
template<size_t Dim>
ModalVector	to_modal_coefficients (const DataVector &nodal_coefficients, const Mesh< Dim > &mesh)
	Compute the modal coefficients from the nodal coefficients.
template<size_t Dim>
ComplexModalVector	to_modal_coefficients (const ComplexDataVector &nodal_coefficients, const Mesh< Dim > &mesh)
	Compute the modal coefficients from the nodal coefficients.
template<size_t Dim>
void	to_nodal_coefficients (gsl::not_null< ComplexDataVector * > nodal_coefficients, const ComplexModalVector &modal_coefficients, const Mesh< Dim > &mesh)
	Compute the nodal coefficients from the modal coefficients.
template<size_t Dim>
void	to_nodal_coefficients (gsl::not_null< DataVector * > nodal_coefficients, const ModalVector &modal_coefficients, const Mesh< Dim > &mesh)
	Compute the nodal coefficients from the modal coefficients.
template<size_t Dim>
DataVector	to_nodal_coefficients (const ModalVector &modal_coefficients, const Mesh< Dim > &mesh)
	Compute the nodal coefficients from the modal coefficients.
template<size_t Dim>
ComplexDataVector	to_nodal_coefficients (const ComplexModalVector &modal_coefficients, const Mesh< Dim > &mesh)
	Compute the nodal coefficients from the modal coefficients.
template<typename VectorType, size_t Dim>
void	PowerMonitors::power_monitors (gsl::not_null< std::array< DataVector, Dim > * > result, const VectorType &u, const Mesh< Dim > &mesh)
	Returns array of power monitors in each spatial dimension.
template<typename VectorType, size_t Dim>
std::array< DataVector, Dim >	PowerMonitors::power_monitors (const VectorType &u, const Mesh< Dim > &mesh)
	Returns array of power monitors in each spatial dimension.
double	PowerMonitors::relative_truncation_error (const DataVector &power_monitor, size_t num_modes_to_use)
	Compute the relative truncation error.
template<typename VectorType, size_t Dim>
std::array< double, Dim >	PowerMonitors::absolute_truncation_error (const VectorType &tensor_component, const Mesh< Dim > &mesh)
	Returns an estimate of the absolute truncation error in each dimension.
void	Spectral::b2_power_monitor_radial (gsl::not_null< DataVector * > result, const DataVector &u, const Mesh< 2 > &mesh)
	Returns the radial B2 power monitor indexed by radial spectral level \(\ell = (n+1)/2\) for a function on a ZernikeB2 disk or cylinder mesh.
void	Spectral::b2_power_monitor_radial (gsl::not_null< DataVector * > result, const DataVector &u, const Mesh< 3 > &mesh)
	Returns the radial B2 power monitor indexed by radial spectral level \(\ell = (n+1)/2\) for a function on a ZernikeB2 disk or cylinder mesh.
void	Spectral::b2_power_monitor_angular (gsl::not_null< DataVector * > result, const DataVector &u, const Mesh< 2 > &mesh)
	Returns the B2 power monitor indexed by azimuthal wavenumber \(m\) for a function on a ZernikeB2 disk or cylinder mesh.
void	Spectral::b2_power_monitor_angular (gsl::not_null< DataVector * > result, const DataVector &u, const Mesh< 3 > &mesh)
	Returns the B2 power monitor indexed by azimuthal wavenumber \(m\) for a function on a ZernikeB2 disk or cylinder mesh.
void	ylm::real_spherical_harmonic (gsl::not_null< DataVector * > spherical_harmonic, const DataVector &theta, const DataVector &phi, size_t l, int m)
	Evaluates a real spherical harmonic of order l at the requested angles \(\theta\) and \(\phi\).
DataVector	ylm::real_spherical_harmonic (const DataVector &theta, const DataVector &phi, size_t l, int m)
	Evaluates a real spherical harmonic of order l at the requested angles \(\theta\) and \(\phi\).

Detailed Description

Things related to spectral transformations.

Typedef Documentation

AngularDerivatives

template<typename DerivativeTagList, ComplexRepresentation Representation = ComplexRepresentation::Interleaved>

using Spectral::Swsh::AngularDerivatives

Initial value:

 detail::AngularDerivativesImpl<
    DerivativeTagList,
    typename detail::unique_derived_from_list<DerivativeTagList>::type,
    Representation>

A DataBox mutate-compatible computational struct for computing a set of spin-weighted spherical harmonic derivatives by grouping and batch-computing spin-weighted spherical harmonic transforms.

Details

A derivative is evaluated for each tag in DerivativeTagList. All entries in DerivativeTagList must be the tag Spectral::Swsh::Tags::Derivative<Tag, DerivativeKind> prefixing the Tag to be differentiated, and indicating the spin-weighted derivative DerivativeKind to be taken. A DataBox on which this struct is invoked must contain:

• each of the tags in DerivativeTagList (the results of the computation)
• each of the tags Tag prefixed by Spectral::Swsh::Tags::Derivative in DerivativeTagList (the inputs of the computation).
• each of the tags Spectral::Swsh::Tags::SwshTransform<DerivativeTag> for DerivativeTagin DerivativeTagList (the buffers for the derivative applied to the modes)
• each of the tags Spectral::Swsh::Tags::SwshTransform<Tag> for Tag prefixed by any DerivativeTag in DerivativeTagList (the buffers for the transforms of the input data).

This function optimizes the derivative taking process by clustering like spins of tags, forward-transforming each spin cluster together, applying the factor for the derivative to each modal vector, re-clustering according to the new spin weights (the derivatives alter the spin weights), and finally inverse-transforming in clusters.

Function Documentation

absolute_truncation_error()

template<typename VectorType, size_t Dim>

std::array< double, Dim > PowerMonitors::absolute_truncation_error	(	const VectorType &	tensor_component,
		const Mesh< Dim > &	mesh )

Returns an estimate of the absolute truncation error in each dimension.

The estimate of the numerical error is given by

\begin{align*} \mathcal{E}\left[P_k\right] = u_\mathrm{max} \times 10^{- \mathcal{T}[P_k]}, \end{align*}

where \( u_\mathrm{max} = \mathrm{max} |u|\) in the corresponding element and \( \mathcal{T}[P_k] \) is the relative error estimate computed from the power monitors \( P_k \).

Warning

This estimate is intended for visualization purposes only.

angular_derivatives()

template<typename DerivativeKindList, ComplexRepresentation Representation = ComplexRepresentation::Interleaved, typename... ArgumentTypes>

void Spectral::Swsh::angular_derivatives	(	const size_t	l_max,
		const size_t	number_of_radial_points,
		const ArgumentTypes &...	arguments )

Evaluate all of the spin-weighted derivatives in DerivKindList on input SpinWeighted<ComplexDataVector, Spin> collocation data, returning by pointer.

Details

This function provides two interfaces, one in which the caller provides the intermediate coefficient buffers needed during the computation of the derivatives, and one in which those buffers are temporarily allocated during the derivative function calls.

For the interface in which the caller does not provide buffers, the arguments must take the following structure (enforced by internal function calls):

• size_t l_max : angular resolution for the spherical representation
• size_t number_of_radial_points : radial resolution (number of consecutive blocks to evaluate derivatives, for each input vector )
• for each DerivKind in DerivKindList, a gsl::not_null<SpinWeighted<ComplexDataVector, Spin + Tags::derivative_spin_weight<DerivKind>>> : the output of the derivative evaluation
• for each DerivKind in DerivKindList, a const SpinWeighted<ComplexDataVector, Spin>& (where the Spin for these arguments matches the corresponding vector from the previous set) : the input to the derivative evaluation.

For the interface in which the caller does provide buffers, the arguments must take the following structure (enforced by internal function calls):

• size_t l_max : angular resolution for the spherical representation
• size_t number_of_radial_points : radial resolution (number of consecutive blocks to evaluate derivatives, for each input vector )
• for each DerivKind in DerivKindList, a gsl::not_null<SpinWeighted<ComplexModalVector, Spin + Tags::derivative_spin_weight<DerivKind>>> : the buffer for the spin-weighted spherical harmonic modes of the derivative quantities.
• for each DerivKind in DerivKindList, a const SpinWeighted<ComplexModalVector, Spin> (where the Spin for these arguments matches the corresponding vector from the previous set) : the buffer for the spin-weighted spherical harmonic modes of the input quantities.
• for each DerivKind in DerivKindList, a gsl::not_null<SpinWeighted<ComplexDataVector, Spin + Tags::derivative_spin_weight<DerivKind>>> : the output of the derivative evaluation
• for each DerivKind in DerivKindList, a const SpinWeighted<ComplexDataVector, Spin> (where the Spin for these arguments matches the corresponding vector from the previous set) : the input to the derivative evaluation.

The function swsh_buffer assists in generating the modal buffers of appropriate size.

b2_power_monitor_angular() [1/2]

void Spectral::b2_power_monitor_angular	(	gsl::not_null< DataVector * >	result,
		const DataVector &	u,
		const Mesh< 2 > &	mesh )

Returns the B2 power monitor indexed by azimuthal wavenumber \(m\) for a function on a ZernikeB2 disk or cylinder mesh.

Details

For functions represented on a filled disk by ZernikeB2 basis functions, the radial and angular spectral spaces are coupled. This function transforms to the combined ZernikeB2 spectral space and returns the root-mean-square of the spectral coefficients at each azimuthal wavenumber \(m = 0, 1, \ldots, M\), where \(M = N_\phi / 2\) and \(N_\phi\) is the number of azimuthal grid points.

The returned DataVector has \(M + 1\) entries. The \(m\)-th entry is

\begin{align*} P_m[\psi] = \sqrt{ \frac{1}{S_m} \sum_n \left| c_{n,m} \right|^2 }, \end{align*}

where \(c_{n,m}\) are the ZernikeB2 spectral coefficients of \(\psi\) at angular wavenumber \(m\) (summing over both cosine and sine components for \(m \geq 1\)), and \(S_m\) is the number of terms in the sum.

For the 3D (cylinder) overload, coefficients are pooled across all \(z\)-slices: \(S_m\) is multiplied by the number of \(z\) points.

b2_power_monitor_angular() [2/2]

void Spectral::b2_power_monitor_angular	(	gsl::not_null< DataVector * >	result,
		const DataVector &	u,
		const Mesh< 3 > &	mesh )

Returns the B2 power monitor indexed by azimuthal wavenumber \(m\) for a function on a ZernikeB2 disk or cylinder mesh.

Details

For functions represented on a filled disk by ZernikeB2 basis functions, the radial and angular spectral spaces are coupled. This function transforms to the combined ZernikeB2 spectral space and returns the root-mean-square of the spectral coefficients at each azimuthal wavenumber \(m = 0, 1, \ldots, M\), where \(M = N_\phi / 2\) and \(N_\phi\) is the number of azimuthal grid points.

The returned DataVector has \(M + 1\) entries. The \(m\)-th entry is

\begin{align*} P_m[\psi] = \sqrt{ \frac{1}{S_m} \sum_n \left| c_{n,m} \right|^2 }, \end{align*}

where \(c_{n,m}\) are the ZernikeB2 spectral coefficients of \(\psi\) at angular wavenumber \(m\) (summing over both cosine and sine components for \(m \geq 1\)), and \(S_m\) is the number of terms in the sum.

For the 3D (cylinder) overload, coefficients are pooled across all \(z\)-slices: \(S_m\) is multiplied by the number of \(z\) points.

b2_power_monitor_radial() [1/2]

void Spectral::b2_power_monitor_radial	(	gsl::not_null< DataVector * >	result,
		const DataVector &	u,
		const Mesh< 2 > &	mesh )

Returns the radial B2 power monitor indexed by radial spectral level \(\ell = (n+1)/2\) for a function on a ZernikeB2 disk or cylinder mesh.

Details

For functions represented on a filled disk by ZernikeB2 basis functions, the radial and angular spectral spaces are coupled. This function transforms to the combined ZernikeB2 spectral space and groups based on spectral level, meaning all spectral modes \((n, m)\) satisfying \((n+1)/2 = \ell\) (using integer division), across all \(m\) and both cosine and sine components, are pooled together.

The returned DataVector has \(N_r\) entries \((\ell = 0, 1, \ldots, N_r - 1)\). The \(\ell\)-th entry is

\begin{align*} P_\ell[\psi] = \sqrt{ \frac{1}{S_\ell} \sum_{\substack{n,m: \\ (n+1)/2 = \ell}} \left| c_{n,m} \right|^2 }, \end{align*}

where \(c_{n,m}\) are the ZernikeB2 spectral coefficients summing over both cosine and sine components for \(m \geq 1\), and \(S_\ell\) is the total number of spectral coefficient slots at level \(\ell\) (including slots that are zero).

For the 3D (cylinder) overload, coefficients are pooled across all \(z\)-slices: \(S_\ell\) is multiplied by the number of \(z\) points.

b2_power_monitor_radial() [2/2]

void Spectral::b2_power_monitor_radial	(	gsl::not_null< DataVector * >	result,
		const DataVector &	u,
		const Mesh< 3 > &	mesh )

Returns the radial B2 power monitor indexed by radial spectral level \(\ell = (n+1)/2\) for a function on a ZernikeB2 disk or cylinder mesh.

Details

For functions represented on a filled disk by ZernikeB2 basis functions, the radial and angular spectral spaces are coupled. This function transforms to the combined ZernikeB2 spectral space and groups based on spectral level, meaning all spectral modes \((n, m)\) satisfying \((n+1)/2 = \ell\) (using integer division), across all \(m\) and both cosine and sine components, are pooled together.

The returned DataVector has \(N_r\) entries \((\ell = 0, 1, \ldots, N_r - 1)\). The \(\ell\)-th entry is

\begin{align*} P_\ell[\psi] = \sqrt{ \frac{1}{S_\ell} \sum_{\substack{n,m: \\ (n+1)/2 = \ell}} \left| c_{n,m} \right|^2 }, \end{align*}

where \(c_{n,m}\) are the ZernikeB2 spectral coefficients summing over both cosine and sine components for \(m \geq 1\), and \(S_\ell\) is the total number of spectral coefficient slots at level \(\ell\) (including slots that are zero).

For the 3D (cylinder) overload, coefficients are pooled across all \(z\)-slices: \(S_\ell\) is multiplied by the number of \(z\) points.

convergence_rate_and_number_of_pile_up_modes()

ConvergenceInfo PowerMonitors::convergence_rate_and_number_of_pile_up_modes	(	const DataVector &	power_monitor,
		size_t	number_of_filtered_modes = 0 )

Returns the convergence rate and the number of pile up modes of a power monitor as a ConvergenceInfo.

Details

Computes the convergence rate of a power monitor as a weighted average of slopes measured using different subsets of spectral modes in the power monitor. Equation (53) of [202] gives the convergence rate \(\mathcal{C}\) in terms of a power monitor \(P_k\) as

\begin{equation}\mathcal{C}(P_k) = -\frac{\sum_{k_1=0}^2\sum_{k_2=\tilde{k}_1}^{\tilde{N}-1} \frac{\mathcal{S}(k_1,k_2)}{\epsilon + \mathcal{E}(k_1,k_2)}}{ \sum_{k_1=0}^2\sum_{k_2=\tilde{k}_1}^{\tilde{N}-1} \frac{1}{\epsilon + \mathcal{E}(k_1,k_2)}}. \end{equation}

Here, \(\mathcal{S}(k_1,k_2)\) is the slope of a linear regression fit of \(\log_{10}(P_k)\) with \(k\) satisfying \(k_1\leq k \leq k_2\), \(\mathcal{E}(k_1,k_2)\) is the error of the slope in that fit, \(\epsilon=\max\left(10^{-3}\max\left(\mathcal{E}(k_1,k_2)\right), 10^{-15}\right)\) is a small number to avoid dividing by zero in the event the fit errors vanish, \(\max\left(\mathcal{E}(k_1,k_2)\right)\) is the maximum fit error of each fit whose slope is included in the summation, \(\tilde{k}_1 = \min\left(k_1+4,\tilde{N}-1\right)\), \(\tilde{N} = N-N_f\), \(N\) is the number of modes in the power monitor, and the highest \(N_f\) modes are filtered. Note that the way \(\epsilon\) is defined is so that it matches SpEC's definition, while also ensuring that it is nonzero even if the error in the slope fit is exactly zero.

Also computes the number of pile up modes in a power monitor. Pile up modes are modes where the power is no longer converging at the overall convergence rate. Following Eq. (56) of [202], the number of pile up modes \(\mathcal{P}\) is defined as

\begin{equation}\mathcal{P}(P_k) = \sum_{j=2}^{\tilde{N}-2} \exp\left[-32\left(\frac{\tilde{\mathcal{C}}_j} {\mathcal{C}(P_k)}\right)^2\right], \end{equation}

where \(\mathcal{C}(P_k)\) is the convergence rate of the power monitor \(P_k\), the local convergence rate \(\tilde{\mathcal{C}}_j\) of mode \(j\) is

\begin{equation}\tilde{\mathcal{C}}_j = -\mathcal{S}(j,\min(\tilde{N}-1,j+4)), \end{equation}

\(\mathcal{S}(k_1,k_2)\) is the slope of a linear regression fit of \(\log_{10}(P_k)\) with \(k\) satisfying \(k_1\leq k \leq k_2\), \(\tilde{N} = N-N_f\), \(N\) is the number of modes in the power monitor, and the highest \(N_f\) modes are filtered. The motivation of this definition is the following: if the local convergence rate \(\tilde{\mathcal{C}}_j\) is comparable to the overall convergence rate \(\mathcal{C}(P_k)\), then the \(j^{\rm th}\) term in the summation becomes \(\approx \exp(-32) \approx 10^{-14}\), while if \(\tilde{\mathcal{C}}_j \ll \mathcal{C}(P_k)\), then the \(j^{\rm th}\) term in the summation is \(\approx \exp(0) = 1\). Note that the coefficient value 32 is chosen to agree with SpEC.

Note

The summation goes up to \(\tilde{N}-2\) so that there is it least one larger unfiltered mode for use in computing the slope. The highest mode used when computing the slope is the highest unfiltered mode, \(\tilde{N}-1\). Mode numbers in the power monitor are zero based. These choices are off by one vs. Eqs. (55) and (56) of [202], because those formulas apparently assume one-based indexing.

Parameters

power_monitor	The power monitor.
number_of_filtered_modes	How many of the highest modes of the power monitor are filtered (default 0).

power_monitors() [1/2]

template<typename VectorType, size_t Dim>

std::array< DataVector, Dim > PowerMonitors::power_monitors	(	const VectorType &	u,
		const Mesh< Dim > &	mesh )

Returns array of power monitors in each spatial dimension.

Computed following Sec. 5.1 of Ref. [202]. For example, in the x dimension (indexed by \( k_0 \)), we compute

\begin{align*} P_{k_0}[\psi] = \sqrt{ \frac{1}{N_1 N_2} \sum_{k_1,k_2} \left| C_{k_0,k_1,k_2} \right|^2} , \end{align*}

where \( C_{k_0,k_1,k_2}\) are the modal coefficients of variable \( \psi \).

For dimensions using a Fourier basis, the cosine and sine power for each wavenumber \( k \) are combined via \( P_k = \sqrt{P_{\cos,k}^2 + P_{\sin,k}^2} \), so the output size for a Fourier dimension of \( N \) points is \( N/2 + 1 \) (integer division) rather than \( N \).

power_monitors() [2/2]

template<typename VectorType, size_t Dim>

void PowerMonitors::power_monitors	(	gsl::not_null< std::array< DataVector, Dim > * >	result,
		const VectorType &	u,
		const Mesh< Dim > &	mesh )

Returns array of power monitors in each spatial dimension.

Computed following Sec. 5.1 of Ref. [202]. For example, in the x dimension (indexed by \( k_0 \)), we compute

\begin{align*} P_{k_0}[\psi] = \sqrt{ \frac{1}{N_1 N_2} \sum_{k_1,k_2} \left| C_{k_0,k_1,k_2} \right|^2} , \end{align*}

where \( C_{k_0,k_1,k_2}\) are the modal coefficients of variable \( \psi \).

For dimensions using a Fourier basis, the cosine and sine power for each wavenumber \( k \) are combined via \( P_k = \sqrt{P_{\cos,k}^2 + P_{\sin,k}^2} \), so the output size for a Fourier dimension of \( N \) points is \( N/2 + 1 \) (integer division) rather than \( N \).

real_spherical_harmonic() [1/2]

DataVector ylm::real_spherical_harmonic	(	const DataVector &	theta,
		const DataVector &	phi,
		size_t	l,
		int	m )

Evaluates a real spherical harmonic of order l at the requested angles \(\theta\) and \(\phi\).

The real spherical harmonics are defined as:

\begin{equation}Y_{lm}(\theta, \phi) = \begin{cases} \sqrt{2} (-1)^m \mathcal{I}[Y_l^{|m|}] & \text{if } m < 0 \\ Y_l^{0} & \text{if } m = 0 \\ \sqrt{2} (-1)^m \mathcal{R}[Y_l^{m}] & \text{if } m > 0 \end{cases} \end{equation}

where \(Y_l^m\) are the complex spherical harmonics and \(\mathcal{R}\) denotes the real part and \(\mathcal{I}\) denotes the imaginary part.

Note

This implementation uses the boost implementation of spherical harmonics. The calculation is not vectorized and may be slow. Stability has not been tested for large l.

real_spherical_harmonic() [2/2]

void ylm::real_spherical_harmonic	(	gsl::not_null< DataVector * >	spherical_harmonic,
		const DataVector &	theta,
		const DataVector &	phi,
		size_t	l,
		int	m )

Evaluates a real spherical harmonic of order l at the requested angles \(\theta\) and \(\phi\).

The real spherical harmonics are defined as:

\begin{equation}Y_{lm}(\theta, \phi) = \begin{cases} \sqrt{2} (-1)^m \mathcal{I}[Y_l^{|m|}] & \text{if } m < 0 \\ Y_l^{0} & \text{if } m = 0 \\ \sqrt{2} (-1)^m \mathcal{R}[Y_l^{m}] & \text{if } m > 0 \end{cases} \end{equation}

where \(Y_l^m\) are the complex spherical harmonics and \(\mathcal{R}\) denotes the real part and \(\mathcal{I}\) denotes the imaginary part.

Note

This implementation uses the boost implementation of spherical harmonics. The calculation is not vectorized and may be slow. Stability has not been tested for large l.

relative_truncation_error()

double PowerMonitors::relative_truncation_error	(	const DataVector &	power_monitor,
		size_t	num_modes_to_use )

Compute the relative truncation error.

The negative logarithm of this quantity is defined by Eqs. (57) and (58) of Ref. [202], i.e.,

\begin{align*} \mathcal{T}\left[P_k\right] = \log_{10} \max \left(P_0, P_1\right) - \dfrac{\sum_{j=0}^{j_{\text{max}, k}} \log_{10} \left(P_j\right) w_j} {\sum_{j=0}^{j_{\text{max}, k}} w_j} , \end{align*}

with weights

\begin{align*} w_j = \exp\left[ - \left(j - j_{\text{max}, k} + \dfrac{1}{2}\right)^2 \right] . \end{align*}

where \( j_{\text{max}, k} = N_k - 1 \) and \( N_k \) is the number of modes or gridpoints in dimension k. Here the second term is a weighted average with larger weights toward the highest modes.

Note

Modes below a cutoff of \(100 \epsilon \mathrm{max}_k(P_k)\) are ignored in the weighted average, where \(\epsilon\) is the machine epsilon. This ensures that we don't underestimate the truncation error if some modes are zero (e.g. by symmetry). Furthermore, if the last two or more modes are zero, we assume that the function is represented exactly and return a relative truncation error of zero.

Details

The number of modes (num_modes_to_use) argument needs to be less or equal than the total number of power monitors (power_monitor.size()). In contrast with Ref. [202], here we index the modes starting from zero.

swsh_buffer()

template<int Spin>

auto Spectral::Swsh::swsh_buffer	(	const size_t	l_max,
		const size_t	number_of_radial_points )

Produces a SpinWeighted<ComplexModalVector, Spin> of the appropriate size to be used as a modal buffer for Spectral::Swsh::AngularDerivatives or Spectral::Swsh::angular_derivatives.

Details

The Spectral::Swsh::angular_derivatives and Spectral::Swsh::AngularDerivatives interfaces require that calling code provides a buffer for the intermediate transform results, to ensure that callers are aware of the allocations and can suitably reuse buffers if possible. This utility eases the creation of those buffers.

to_modal_coefficients() [1/4]

template<size_t Dim>

ComplexModalVector to_modal_coefficients	(	const ComplexDataVector &	nodal_coefficients,
		const Mesh< Dim > &	mesh )

Compute the modal coefficients from the nodal coefficients.

See also

Spectral::nodal_to_modal_matrix

to_modal_coefficients() [2/4]

template<size_t Dim>

ModalVector to_modal_coefficients	(	const DataVector &	nodal_coefficients,
		const Mesh< Dim > &	mesh )

Compute the modal coefficients from the nodal coefficients.

See also

Spectral::nodal_to_modal_matrix

to_modal_coefficients() [3/4]

template<size_t Dim>

void to_modal_coefficients	(	gsl::not_null< ComplexModalVector * >	modal_coefficients,
		const ComplexDataVector &	nodal_coefficients,
		const Mesh< Dim > &	mesh )

Compute the modal coefficients from the nodal coefficients.

See also

Spectral::nodal_to_modal_matrix

to_modal_coefficients() [4/4]

template<size_t Dim>

void to_modal_coefficients	(	gsl::not_null< ModalVector * >	modal_coefficients,
		const DataVector &	nodal_coefficients,
		const Mesh< Dim > &	mesh )

Compute the modal coefficients from the nodal coefficients.

See also

Spectral::nodal_to_modal_matrix

to_nodal_coefficients() [1/4]

template<size_t Dim>

ComplexDataVector to_nodal_coefficients	(	const ComplexModalVector &	modal_coefficients,
		const Mesh< Dim > &	mesh )

Compute the nodal coefficients from the modal coefficients.

See also

Spectral::modal_to_nodal_matrix

to_nodal_coefficients() [2/4]

template<size_t Dim>

DataVector to_nodal_coefficients	(	const ModalVector &	modal_coefficients,
		const Mesh< Dim > &	mesh )

Compute the nodal coefficients from the modal coefficients.

See also

Spectral::modal_to_nodal_matrix

to_nodal_coefficients() [3/4]

template<size_t Dim>

void to_nodal_coefficients	(	gsl::not_null< ComplexDataVector * >	nodal_coefficients,
		const ComplexModalVector &	modal_coefficients,
		const Mesh< Dim > &	mesh )

Compute the nodal coefficients from the modal coefficients.

See also

Spectral::modal_to_nodal_matrix

to_nodal_coefficients() [4/4]

template<size_t Dim>

void to_nodal_coefficients	(	gsl::not_null< DataVector * >	nodal_coefficients,
		const ModalVector &	modal_coefficients,
		const Mesh< Dim > &	mesh )

Compute the nodal coefficients from the modal coefficients.

See also

Spectral::modal_to_nodal_matrix

zernike_b2_disk_nodal_to_modal()

void Spectral::zernike_b2_disk_nodal_to_modal	(	gsl::not_null< DataVector * >	modal,
		gsl::not_null< DataVector * >	buf,
		const DataVector &	u,
		size_t	n_r,
		size_t	n_phi,
		size_t	n_r_max,
		size_t	num_components = 1 )

Transform nodal ZernikeB2 data on a disk to modal space.

Transforms num_components independent disk slices simultaneously. This is the same operation as done in the disk filtering, while cylinder filtering batches more optimally so we choose not to refactor them to use this looped version. The Fourier nodal-to-modal step is applied to all components in a single batched DGEMM while the Zernike nodal-to-modal step loops over components.

Parameters

modal	Output spectral coefficients. Size: num_components * zernike_b2_disk_spectral_size(n_r_max, n_phi/2). Layout: modal[comp * spectral_size + spec_index].
buf	Scratch buffer. Size: >= 2 * num_components * n_r * n_phi.
u	Input nodal values. Size: num_components * n_r * n_phi. Layout: u[comp * n_r * n_phi + i_r + n_r * j_phi].
n_r	Number of radial grid points.
n_phi	Number of azimuthal grid points (must be odd).
n_r_max	Maximum Zernike degree.
num_components	Number of independent disk slices to transform.