Defines the C++ interface to SPHEREPACK. More...

#include <YlmSpherepack.hpp>

Classes
struct	InterpolationInfo
	Struct to hold cached information at a set of target interpolation points. More...

Public Types
using	FirstDeriv = tnsr::i< DataVector, 2, Frame::Logical >
	Type returned by gradient function.

using	SecondDeriv = tnsr::ij< DataVector, 2, Frame::Logical >
	Type returned by second derivative function.

Public Member Functions
	YlmSpherepack (size_t l_max, size_t m_max) noexcept
	Here l_max and m_max are the largest fully-represented l and m in the Ylm expansion.

std::array< size_t, 2 >	physical_extents () const noexcept

void	gradient (const std::array< double , 2 > &df, gsl::not_null< const double > collocation_values, size_t physical_stride=1, size_t physical_offset=0) const noexcept
	Computes Pfaffian derivative (df/dtheta, csc(theta) df/dphi) at the collocation values. To act on a slice of the input and output arrays, specify stride and offset (assumed to be the same for input and output).

void	gradient_from_coefs (const std::array< double , 2 > &df, gsl::not_null< const double > spectral_coefs, size_t spectral_stride=1, size_t spectral_offset=0, size_t physical_stride=1, size_t physical_offset=0) const noexcept
	Same as `gradient`, but takes the spectral coefficients (rather than collocation values) of the function. This is more efficient if one happens to already have the spectral coefficients. To act on a slice of the input and output arrays, specify strides and offsets.

void	scalar_laplacian (gsl::not_null< double * > scalar_laplacian, gsl::not_null< const double * > collocation_values, size_t physical_stride=1, size_t physical_offset=0) const noexcept
	Computes Laplacian in physical space. To act on a slice of the input and output arrays, specify stride and offset (assumed to be the same for input and output).

void	scalar_laplacian_from_coefs (gsl::not_null< double * > scalar_laplacian, gsl::not_null< const double * > spectral_coefs, size_t spectral_stride=1, size_t spectral_offset=0, size_t physical_stride=1, size_t physical_offset=0) const noexcept
	Same as `scalar_laplacian` above, but the input is the spectral coefficients (rather than collocation values) of the function. This is more efficient if one happens to already have the spectral coefficients. To act on a slice of the input and output arrays, specify strides and offsets.

void	second_derivative (const std::array< double , 2 > &df, gsl::not_null< SecondDeriv > ddf, gsl::not_null< const double * > collocation_values, size_t physical_stride=1, size_t physical_offset=0) const noexcept
	Computes Pfaffian first and second derivative in physical space. The first derivative is \(df(i) = d_i f\), and the second derivative is \(ddf(i,j) = d_i (d_j f)\), where \(d_0 = d/d\theta\) and \(d_1 = csc(\theta) d/d\phi\). ddf is not symmetric. To act on a slice of the input and output arrays, specify stride and offset (assumed to be the same for input and output).

std::pair< FirstDeriv, SecondDeriv >	first_and_second_derivative (const DataVector &collocation_values) const noexcept
	Simpler, less general interface to second_derivative.

double	definite_integral (gsl::not_null< const double * > collocation_values, size_t physical_stride=1, size_t physical_offset=0) const noexcept
	Computes the integral over the sphere.

const std::vector< double > &	integration_weights () const noexcept
	Returns weights \(w_i\) such that \(sum_i (c_i w_i)\) is the definite integral, where \(c_i\) are collocation values at point i.

template<typename T >
InterpolationInfo< T >	set_up_interpolation_info (const std::array< T, 2 > &target_points) const noexcept
	Sets up the `InterpolationInfo` structure for interpolating onto a set of target \((\theta,\phi)\) points. Does not depend on the function being interpolated.

template<typename T >
void	interpolate (gsl::not_null< T * > result, gsl::not_null< const double * > collocation_values, const InterpolationInfo< T > &interpolation_info, size_t physical_stride=1, size_t physical_offset=0) const noexcept
	Interpolates from `collocation_values` onto the points that have been passed into the `set_up_interpolation_info` function. To interpolate a different function on the same spectral grid, there is no need to recompute `interpolation_info`. If you specify stride and offset, acts on a slice of the input values. The output has unit stride.

template<typename T , typename R >
void	interpolate_from_coefs (gsl::not_null< T * > result, const R &spectral_coefs, const InterpolationInfo< T > &interpolation_info, size_t spectral_stride=1, size_t spectral_offset=0) const noexcept
	Same as `interpolate`, but assumes you have spectral coefficients. This is more efficient if you already have the spectral coefficients available. If you specify stride and offset, acts on a slice of the input coefs. The output has unit stride.

template<typename T >
T	interpolate (const DataVector &collocation_values, const std::array< T, 2 > &target_points) const noexcept
	Simpler interface to `interpolate`. If you need to call this repeatedly on different `spectral_coefs` or `collocation_values` for the same target points, this is inefficient; instead use `set_up_interpolation_info` and the functions that use `InterpolationInfo`.

template<typename T >
T	interpolate_from_coefs (const DataVector &spectral_coefs, const std::array< T, 2 > &target_points) const noexcept

DataVector	prolong_or_restrict (const DataVector &spectral_coefs, const YlmSpherepack &target) const noexcept
	Takes spectral coefficients compatible with `*this`, and either prolongs them or restricts them to be compatible with `target`. This is done by truncation (restriction) or padding with zeros (prolongation).


size_t	l_max () const noexcept
	Sizes in physical and spectral space for this instance.

size_t	m_max () const noexcept
	Sizes in physical and spectral space for this instance.

size_t	physical_size () const noexcept
	Sizes in physical and spectral space for this instance.

size_t	spectral_size () const noexcept
	Sizes in physical and spectral space for this instance.


const std::vector< double > &	theta_points () const noexcept
	Collocation points theta and phi. More...

const std::vector< double > &	phi_points () const noexcept
	Collocation points theta and phi. More...

std::array< DataVector, 2 >	theta_phi_points () const noexcept
	Collocation points theta and phi. More...


void	phys_to_spec (gsl::not_null< double * > spectral_coefs, gsl::not_null< const double * > collocation_values, size_t physical_stride=1, size_t physical_offset=0, size_t spectral_stride=1, size_t spectral_offset=0) const noexcept
	Spectral transformations. To act on a slice of the input and output arrays, specify strides and offsets.

void	spec_to_phys (gsl::not_null< double * > collocation_values, gsl::not_null< const double * > spectral_coefs, size_t spectral_stride=1, size_t spectral_offset=0, size_t physical_stride=1, size_t physical_offset=0) const noexcept
	Spectral transformations. To act on a slice of the input and output arrays, specify strides and offsets.


void	phys_to_spec_all_offsets (gsl::not_null< double * > spectral_coefs, gsl::not_null< const double * > collocation_values, size_t stride) const noexcept
	Spectral transformations where `collocation_values` and `spectral_coefs` are assumed to point to 3-dimensional arrays (I1 x S2 topology), and the transformations are done for all 'radial' points at once by internally looping over all values of the offset from zero to `stride`-1 (the physical and spectral strides are equal and are called `stride`).

void	spec_to_phys_all_offsets (gsl::not_null< double * > collocation_values, gsl::not_null< const double * > spectral_coefs, size_t stride) const noexcept
	Spectral transformations where `collocation_values` and `spectral_coefs` are assumed to point to 3-dimensional arrays (I1 x S2 topology), and the transformations are done for all 'radial' points at once by internally looping over all values of the offset from zero to `stride`-1 (the physical and spectral strides are equal and are called `stride`).


DataVector	phys_to_spec (const DataVector &collocation_values, size_t physical_stride=1, size_t physical_offset=0) const noexcept
	Simpler, less general interfaces to `phys_to_spec` and `spec_to_phys`. Acts on a slice of the input and returns a unit-stride result.

DataVector	spec_to_phys (const DataVector &spectral_coefs, size_t spectral_stride=1, size_t spectral_offset=0) const noexcept
	Simpler, less general interfaces to `phys_to_spec` and `spec_to_phys`. Acts on a slice of the input and returns a unit-stride result.


DataVector	phys_to_spec_all_offsets (const DataVector &collocation_values, size_t stride) const noexcept
	Simpler, less general interfaces to `phys_to_spec_all_offsets` and `spec_to_phys_all_offsets`. Result has the same stride as the input.

DataVector	spec_to_phys_all_offsets (const DataVector &spectral_coefs, size_t stride) const noexcept
	Simpler, less general interfaces to `phys_to_spec_all_offsets` and `spec_to_phys_all_offsets`. Result has the same stride as the input.


void	gradient_all_offsets (const std::array< double , 2 > &df, gsl::not_null< const double > collocation_values, size_t stride=1) const noexcept
	Same as `gradient` but pointers are assumed to point to 3-dimensional arrays (I1 x S2 topology), and the gradient is done for all 'radial' points at once by internally looping over all values of the offset from zero to `stride`-1.

void	gradient_from_coefs_all_offsets (const std::array< double , 2 > &df, gsl::not_null< const double > spectral_coefs, size_t stride=1) const noexcept
	Same as `gradient` but pointers are assumed to point to 3-dimensional arrays (I1 x S2 topology), and the gradient is done for all 'radial' points at once by internally looping over all values of the offset from zero to `stride`-1.


FirstDeriv	gradient (const DataVector &collocation_values, size_t physical_stride=1, size_t physical_offset=0) const noexcept
	Simpler, less general interfaces to `gradient`. Acts on a slice of the input and returns a unit-stride result.

FirstDeriv	gradient_from_coefs (const DataVector &spectral_coefs, size_t spectral_stride=1, size_t spectral_offset=0) const noexcept
	Simpler, less general interfaces to `gradient`. Acts on a slice of the input and returns a unit-stride result.


FirstDeriv	gradient_all_offsets (const DataVector &collocation_values, size_t stride=1) const noexcept
	Simpler, less general interfaces to `gradient_all_offsets`. Result has the same stride as the input.

FirstDeriv	gradient_from_coefs_all_offsets (const DataVector &spectral_coefs, size_t stride=1) const noexcept
	Simpler, less general interfaces to `gradient_all_offsets`. Result has the same stride as the input.


DataVector	scalar_laplacian (const DataVector &collocation_values, size_t physical_stride=1, size_t physical_offset=0) const noexcept
	Simpler, less general interfaces to `scalar_laplacian`. Acts on a slice of the input and returns a unit-stride result.

DataVector	scalar_laplacian_from_coefs (const DataVector &spectral_coefs, size_t spectral_stride=1, size_t spectral_offset=0) const noexcept
	Simpler, less general interfaces to `scalar_laplacian`. Acts on a slice of the input and returns a unit-stride result.

Static Public Member Functions
static void	add_constant (const gsl::not_null< DataVector * > spectral_coefs, const double c) noexcept
	Adds a constant (i.e. \(f(\theta,\phi)\) += \(c\)) to the function given by the spectral coefficients, by modifying the coefficients.

static double	average (const DataVector &spectral_coefs) noexcept
	Returns the average of \(f(\theta,\phi)\) over \((\theta,\phi)\).


static constexpr size_t	physical_size (const size_t l_max, const size_t m_max) noexcept
	Static functions to return the correct sizes of vectors of collocation points and spectral coefficients for a given l_max and m_max. Useful for allocating space without having to create a YlmSpherepack.

static constexpr size_t	spectral_size (const size_t l_max, const size_t m_max) noexcept

Detailed Description

Defines the C++ interface to SPHEREPACK.

Details

The class YlmSpherepack defines the C++ interface to the fortran library SPHEREPACK used for computations on the surface of a sphere.

Given a real-valued, scalar function \(g(\theta, \phi)\), SPHEREPACK expands it as:

\begin{align*} g(\theta, \phi) &=\frac{1}{2}\sum_{l=0}^{l_{\max}}\bar P_l^0(\cos\theta) a_{l0} +\sum_{l=1}^{l_{\max}}\sum_{m=1}^{\min(l, m_{\max})}\bar P_l^m(\cos\theta)\{ a_{lm}\cos m\phi -b_{lm}\sin m\phi\}\nonumber \end{align*}

where \(a_{lm}\) and \(b_{lm}\) are spectral coefficient arrays used by SPHEREPACK, \(P_l^m(x)\) are defined as

\begin{align*} \bar P_l^m(x)&=\sqrt{\frac{(2l+1)(l-m)!}{2(l+m)!}}\;P_{lm}(x) \end{align*}

and \(P_{nm}(x)\) are the associated Legendre polynomials.

Internally, SPHEREPACK can represent such an expansion in two ways which we will refer to as modal and nodal representations:

modal: The spectral coefficient arrays \(a_{lm}\) and \(b_{lm}\), referred to as spectral_coefs in the methods below. For this C++ interface, they are saved in a single DataVector. To help you index the coefficients as expected by this interface, use the class SpherepackIterator.
nodal: The values at certain collocation points, referred to as collocation_values in the methods below. This is an array of the expanded function \(g(\theta,\phi)\) evaluated at collocation values \((\theta_i,\phi_j)\), where \(\theta_i\) are Gauss-Legendre quadrature nodes in the interval \((0, \pi)\) with \(i = 0, ..., l_{\max}\), and \(\phi_j\) is distributed uniformly in \((0, 2\pi)\) with \(i = 0, ..., 2m_{\max}\). The angles of the collocation points can be computed with the method theta_phi_points.

To convert between the two representations the methods spec_to_phys and phys_to_spec can be used. For internal calculations SPHEREPACK will usually convert to spectral coefficients first, so it is in general more efficient to use these directly.

Most methods of SPHEREPACK will compute the requested values of e.g. gradient or scalar_laplacian at the collocation points, effectively returning an expansion in nodal form as defined above. To evaluate the function at arbitrary angles \(\theta\), \(\phi\), these values have to be "interpolated" (i.e. the new expansion evaluated) using interpolate.

Member Function Documentation

◆ phi_points()

const std::vector< double > & YlmSpherepack::phi_points ( ) const

noexcept

Collocation points theta and phi.

The phi points are uniform in phi, with the first point at phi=0.

The theta points are Gauss-Legendre in \(\cos(\theta)\), so there are no points at the poles.

◆ spectral_size()

static constexpr size_t YlmSpherepack::spectral_size	(	const size_t	l_max,
		const size_t	m_max
	)

inlinestaticconstexprnoexcept

Note: spectral_size is the size of the buffer that holds the coefficients; it is not the number of coefficients (which is \(m_{\rm max}^2+(l_{\rm max}-m_{\rm max})(2m_{\rm max}+1)\)). To simplify its internal indexing, SPHEREPACK uses a buffer with more space than necessary. See SpherepackIterator for how to index the coefficients in the buffer.

◆ theta_phi_points()

std::array< DataVector, 2 > YlmSpherepack::theta_phi_points ( ) const

noexcept

Collocation points theta and phi.

The phi points are uniform in phi, with the first point at phi=0.

The theta points are Gauss-Legendre in \(\cos(\theta)\), so there are no points at the poles.

◆ theta_points()

const std::vector< double > & YlmSpherepack::theta_points ( ) const

noexcept

Collocation points theta and phi.

The phi points are uniform in phi, with the first point at phi=0.

The theta points are Gauss-Legendre in \(\cos(\theta)\), so there are no points at the poles.

Classes

Public Types

Public Member Functions

Static Public Member Functions

Detailed Description

Details

Member Function Documentation

◆ phi_points()

◆ spectral_size()

◆ theta_phi_points()

◆ theta_points()