Namespaces | Enumerations | Functions | Variables
Spectral Namespace Reference

Functionality associated with a particular choice of basis functions and quadrature for spectral operations. More...

Namespaces

 Swsh
 Namespace for spin-weighted spherical harmonic utilities.
 

Enumerations

enum  MortarSize { Full, UpperHalf, LowerHalf }
 The portion of an element covered by a mortar.
 
enum  Basis { Chebyshev, Legendre }
 The choice of basis functions for computing collocation points and weights. More...
 
enum  Quadrature { Gauss, GaussLobatto }
 The choice of quadrature method to compute integration weights. More...
 

Functions

const Matrixprojection_matrix_mortar_to_element (MortarSize size, const Mesh< 1 > &element_mesh, const Mesh< 1 > &mortar_mesh) noexcept
 The projection matrix from a 1D mortar to an element. More...
 
const Matrixprojection_matrix_element_to_mortar (MortarSize size, const Mesh< 1 > &mortar_mesh, const Mesh< 1 > &element_mesh) noexcept
 The projection matrix from a 1D element to a mortar. See projection_matrix_mortar_to_element() for details.
 
std::ostreamoperator<< (std::ostream &os, const Basis &basis) noexcept
 
std::ostreamoperator<< (std::ostream &os, const Quadrature &quadrature) noexcept
 
template<Basis BasisType, Quadrature QuadratureType>
const DataVectorcollocation_points (size_t num_points) noexcept
 Collocation points. More...
 
template<Basis BasisType, Quadrature QuadratureType>
const DataVectorquadrature_weights (size_t num_points) noexcept
 Weights to compute definite integrals. More...
 
template<Basis BasisType, Quadrature QuadratureType, typename T >
Matrix interpolation_matrix (size_t num_points, const T &target_points) noexcept
 Matrix used to interpolate to the target_points. More...
 
template<typename T >
Matrix interpolation_matrix (const Mesh< 1 > &mesh, const T &target_points) noexcept
 Interpolation matrix to the target_points for a one-dimensional mesh. More...
 
template<Basis BasisType, typename T >
compute_basis_function_value (size_t k, const T &x) noexcept
 Compute the function values of the basis function \(\Phi_k(x)\) (zero-indexed).
 
template<Basis >
DataVector compute_inverse_weight_function_values (const DataVector &) noexcept
 Compute the inverse of the weight function \(w(x)\) w.r.t. which the basis functions are orthogonal. See the description of quadrature_weights(size_t) for details.
 
template<Basis BasisType>
double compute_basis_function_normalization_square (size_t k) noexcept
 Compute the normalization square of the basis function \(\Phi_k\) (zero-indexed), i.e. the weighted definite integral over its square.
 
template<Basis BasisType, Quadrature QuadratureType>
std::pair< DataVector, DataVectorcompute_collocation_points_and_weights (size_t num_points) noexcept
 Compute the collocation points and weights associated to the basis and quadrature. More...
 
const DataVectorcollocation_points (const Mesh< 1 > &mesh) noexcept
 Collocation points for a one-dimensional mesh. More...
 
const DataVectorquadrature_weights (const Mesh< 1 > &mesh) noexcept
 Quadrature weights for a one-dimensional mesh. More...
 
template<Basis BasisType, Quadrature QuadratureType>
const Matrixdifferentiation_matrix (size_t num_points) noexcept
 Matrix used to compute the derivative of a function. More...
 
const Matrixdifferentiation_matrix (const Mesh< 1 > &mesh) noexcept
 Differentiation matrix for a one-dimensional mesh. More...
 
template<Basis BasisType, Quadrature QuadratureType>
const Matrixspectral_to_grid_points_matrix (size_t num_points) noexcept
 Matrix used to transform from the spectral coefficients (modes) of a function to its nodal coefficients. Also referred to as the Vandermonde matrix. More...
 
const Matrixspectral_to_grid_points_matrix (const Mesh< 1 > &mesh) noexcept
 Transformation matrix from modal to nodal coefficients for a one-dimensional mesh. More...
 
template<Basis BasisType, Quadrature QuadratureType>
const Matrixgrid_points_to_spectral_matrix (size_t num_points) noexcept
 Matrix used to transform from the nodal coefficients of a function to its spectral coefficients (modes). Also referred to as the inverse Vandermonde matrix. More...
 
const Matrixgrid_points_to_spectral_matrix (const Mesh< 1 > &mesh) noexcept
 Transformation matrix from nodal to modal coefficients for a one-dimensional mesh. More...
 
template<Basis BasisType, Quadrature QuadratureType>
const Matrixlinear_filter_matrix (size_t num_points) noexcept
 Matrix used to linearize a function. More...
 
const Matrixlinear_filter_matrix (const Mesh< 1 > &mesh) noexcept
 Linear filter matrix for a one-dimensional mesh. More...
 

Variables

template<Basis , Quadrature >
constexpr size_t minimum_number_of_points {std::numeric_limits<size_t>::max()}
 Minimum number of possible collocation points for a quadrature type. More...
 
template<Basis >
constexpr size_t maximum_number_of_points = 12
 Maximum number of allowed collocation points.
 

Detailed Description

Functionality associated with a particular choice of basis functions and quadrature for spectral operations.

Details

The functions in this namespace provide low-level access to collocation points, quadrature weights and associated matrices, such as for differentiation and interpolation. They are available in two versions: either templated directly on the enum cases of the Spectral::Basis and Spectral::Quadrature types, or taking a one-dimensional Mesh as their argument.

Note
Generally you should prefer working with a Mesh. Use its Mesh::slice_through() method to retrieve the mesh in a particular dimension:
const Mesh<2> mesh2d{
{{3, 4}},
{{Spectral::Basis::Legendre, Spectral::Basis::Legendre}},
{{Spectral::Quadrature::Gauss, Spectral::Quadrature::GaussLobatto}}};
const auto collocation_points_in_first_dim =
Spectral::collocation_points(mesh2d.slice_through(0));

Most algorithms in this namespace are adapted from [16].

Function Documentation

◆ collocation_points() [1/2]

const DataVector& Spectral::collocation_points ( const Mesh< 1 > &  mesh)
noexcept

Collocation points for a one-dimensional mesh.

See also
collocation_points(size_t)

◆ collocation_points() [2/2]

template<Basis BasisType, Quadrature QuadratureType>
const DataVector & Spectral::collocation_points ( size_t  num_points)
noexcept

Collocation points.

Parameters
num_pointsThe number of collocation points

◆ compute_collocation_points_and_weights()

template<Basis BasisType, Quadrature QuadratureType>
std::pair<DataVector, DataVector> Spectral::compute_collocation_points_and_weights ( size_t  num_points)
noexcept

Compute the collocation points and weights associated to the basis and quadrature.

Details

This function is expected to return the tuple \((\xi_k,w_k)\) where the \(\xi_k\) are the collocation points and the \(w_k\) are defined in the description of quadrature_weights(size_t).

◆ differentiation_matrix() [1/2]

template<Basis BasisType, Quadrature QuadratureType>
const Matrix& Spectral::differentiation_matrix ( size_t  num_points)
noexcept

Matrix used to compute the derivative of a function.

Details

For a function represented by the nodal coefficients \(u_j\) a matrix multiplication with the differentiation matrix \(D_{ij}\) gives the coefficients of the function's derivative. Since \(u(x)\) is expanded in Lagrange polynomials \(u(x)=\sum_j u_j l_j(x)\) the differentiation matrix is computed as \(D_{ij}=l_j^\prime(\xi_i)\) where the \(\xi_i\) are the collocation points.

Parameters
num_pointsThe number of collocation points

◆ differentiation_matrix() [2/2]

const Matrix& Spectral::differentiation_matrix ( const Mesh< 1 > &  mesh)
noexcept

Differentiation matrix for a one-dimensional mesh.

See also
differentiation_matrix(size_t)

◆ grid_points_to_spectral_matrix() [1/2]

template<Basis BasisType, Quadrature QuadratureType>
const Matrix& Spectral::grid_points_to_spectral_matrix ( size_t  num_points)
noexcept

Matrix used to transform from the nodal coefficients of a function to its spectral coefficients (modes). Also referred to as the inverse Vandermonde matrix.

Details

This is the inverse to the Vandermonde matrix \(\mathcal{V}\) computed in spectral_to_grid_points_matrix(size_t). It can be computed analytically for Gauss quadrature by evaluating \(\sum_j\mathcal{V}^{-1}_{ij}u_j=\widetilde{u}_i= \frac{(u,\Phi_i)}{\gamma_i}\) for a Lagrange basis function \(u(x)=l_k(x)\) to find \(\mathcal{V}^{-1}_{ij}=\mathcal{V}_{ji}\frac{w_j}{\gamma_i}\) where the \(w_j\) are the Gauss quadrature weights and \(\gamma_i\) is the norm square of the spectral basis function \(\Phi_i\).

Parameters
num_pointsThe number of collocation points
See also
spectral_to_grid_points_matrix(size_t)

◆ grid_points_to_spectral_matrix() [2/2]

const Matrix& Spectral::grid_points_to_spectral_matrix ( const Mesh< 1 > &  mesh)
noexcept

Transformation matrix from nodal to modal coefficients for a one-dimensional mesh.

See also
grid_points_to_spectral_matrix(size_t)

◆ interpolation_matrix() [1/2]

template<Basis BasisType, Quadrature QuadratureType, typename T >
Matrix Spectral::interpolation_matrix ( size_t  num_points,
const T &  target_points 
)
noexcept

Matrix used to interpolate to the target_points.

Warning
For each target point located outside of the logical coordinate interval covered by BasisType (often \([-1,1]\)), the resulting matrix performs polynomial extrapolation instead of interpolation. The extapolation will be correct but may suffer from reduced accuracy, especially for higher-order polynomials (i.e., larger values of num_points).
Parameters
num_pointsThe number of collocation points
target_pointsThe points to interpolate to

◆ interpolation_matrix() [2/2]

template<typename T >
Matrix Spectral::interpolation_matrix ( const Mesh< 1 > &  mesh,
const T &  target_points 
)
noexcept

Interpolation matrix to the target_points for a one-dimensional mesh.

See also
interpolation_matrix(size_t, const T&)

◆ linear_filter_matrix() [1/2]

template<Basis BasisType, Quadrature QuadratureType>
const Matrix& Spectral::linear_filter_matrix ( size_t  num_points)
noexcept

Matrix used to linearize a function.

Details

Filters out all except the lowest two modes by applying \(\mathcal{V}^{-1}\cdot\mathrm{diag}(1,1,0,0,...)\cdot\mathcal{V}\) to the nodal coefficients, where \(\mathcal{V}\) is the Vandermonde matrix computed in spectral_to_grid_points_matrix(size_t).

Parameters
num_pointsThe number of collocation points
See also
spectral_to_grid_points_matrix(size_t)
grid_points_to_spectral_matrix(size_t)

◆ linear_filter_matrix() [2/2]

const Matrix& Spectral::linear_filter_matrix ( const Mesh< 1 > &  mesh)
noexcept

Linear filter matrix for a one-dimensional mesh.

See also
linear_filter_matrix(size_t)

◆ projection_matrix_mortar_to_element()

const Matrix & Spectral::projection_matrix_mortar_to_element ( MortarSize  size,
const Mesh< 1 > &  element_mesh,
const Mesh< 1 > &  mortar_mesh 
)
noexcept

The projection matrix from a 1D mortar to an element.

Details

The projection matrices returned by this function (and by projection_matrix_element_to_mortar()) define orthogonal projection operators between the spaces of functions on the boundary of the element and the mortar. These projections are usually the correct way to transfer data onto and off of the mortars.

The half-interval projections are based on an equation derived by Saul. This shows that the projection from the spectral basis for the entire interval to the spectral basis for the upper half interval is

\begin{equation*} T_{jk} = \frac{2 j + 1}{2} 2^j \sum_{n=0}^{j-k} \binom{j}{k+n} \binom{(j + k + n - 1)/2}{j} \frac{(k + n)!^2}{(2 k + n + 1)! n!} \end{equation*}

◆ quadrature_weights() [1/2]

const DataVector& Spectral::quadrature_weights ( const Mesh< 1 > &  mesh)
noexcept

Quadrature weights for a one-dimensional mesh.

See also
quadrature_weights(size_t)

◆ quadrature_weights() [2/2]

template<Basis BasisType, Quadrature QuadratureType>
const DataVector & Spectral::quadrature_weights ( size_t  num_points)
noexcept

Weights to compute definite integrals.

Details

These are the coefficients to contract with the nodal function values \(f_k\) to approximate the definite integral \(I[f]=\int f(x)\mathrm{d}x\).

Note that the term quadrature also often refers to the quantity \(Q[f]=\int f(x)w(x)\mathrm{d}x\approx \sum_k f_k w_k\). Here, \(w(x)\) denotes the basis-specific weight function w.r.t. to which the basis functions \(\Phi_k\) are orthogonal, i.e \(\int\Phi_i(x)\Phi_j(x)w(x)=0\) for \(i\neq j\). The weights \(w_k\) approximate this inner product. To approximate the definite integral \(I[f]\) we must employ the coefficients \(\frac{w_k}{w(\xi_k)}\) instead, where the \(\xi_k\) are the collocation points. These are the coefficients this function returns. Only for a unit weight function \(w(x)=1\), i.e. a Legendre basis, is \(I[f]=Q[f]\) so this function returns the \(w_k\) identically.

Parameters
num_pointsThe number of collocation points

◆ spectral_to_grid_points_matrix() [1/2]

template<Basis BasisType, Quadrature QuadratureType>
const Matrix& Spectral::spectral_to_grid_points_matrix ( size_t  num_points)
noexcept

Matrix used to transform from the spectral coefficients (modes) of a function to its nodal coefficients. Also referred to as the Vandermonde matrix.

Details

The Vandermonde matrix is computed as \(\mathcal{V}_{ij}=\Phi_j(\xi_i)\) where the \(\Phi_j(x)\) are the spectral basis functions used in the modal expansion \(u(x)=\sum_j \widetilde{u}_j\Phi_j(x)\), e.g. normalized Legendre polynomials, and the \(\xi_i\) are the collocation points used to construct the interpolating Lagrange polynomials in the nodal expansion \(u(x)=\sum_j u_j l_j(x)\). Then the Vandermonde matrix arises as \(u(\xi_i)=u_i=\sum_j \widetilde{u}_j\Phi_j(\xi_i)=\sum_j \mathcal{V}_{ij}\widetilde{u}_j\).

Parameters
num_pointsThe number of collocation points
See also
grid_points_to_spectral_matrix(size_t)

◆ spectral_to_grid_points_matrix() [2/2]

const Matrix& Spectral::spectral_to_grid_points_matrix ( const Mesh< 1 > &  mesh)
noexcept

Transformation matrix from modal to nodal coefficients for a one-dimensional mesh.

See also
spectral_to_grid_points_matrix(size_t)

Variable Documentation

◆ minimum_number_of_points

template<Basis , Quadrature >
constexpr size_t Spectral::minimum_number_of_points {std::numeric_limits<size_t>::max()}

Minimum number of possible collocation points for a quadrature type.

Details

Since Gauss-Lobatto quadrature has points on the domain boundaries it must have at least two collocation points. Gauss quadrature can have only one collocation point.