shenfun.legendre package

Submodules

shenfun.legendre.bases module

Module for defining function spaces in the Legendre family.

A function is approximated in the Legendre basis as

\[u(x) = \sum_{i=0}^{N-1} \hat{u}_i L_i(x)\]

where \(L_i(x)\) is the i’th Legendre polynomial of the first kind. The Legendre polynomials are orthogonal with weight \(\omega=1\)

\[\int_{-1}^1 L_i L_k dx = \frac{2}{2k+1} \delta_{ki}.\]

All other bases defined in this module are combinations of \(L_i\)’s. For example, a Dirichlet basis is

\[\phi_i = L_i - L_{i+2}\]

The basis is implemented using a stencil matrix \(K \in \mathbb{R}^{N-2 \times N}\), such that

\[\boldsymbol{\phi} = K \boldsymbol{L},\]

where \(\boldsymbol{\phi}=(\phi_0, \phi_1, \ldots, \phi_{N-3})\) and \(\boldsymbol{L}=(L_0, L_1, \ldots, L_{N-1})\). For the Dirichlet basis \(K = (\delta_{i, j} - \delta_{i+2, j})_{i,j=0}^{N-2, N}\).

The stencil matrix is used to transfer any composite basis back and forth to the orthogonal basis.

class shenfun.legendre.bases.BCGeneric(N, bc=None, domain=None, alpha=0, beta=0, **kw)

Bases: CompositeBase

Function space for setting inhomogeneous boundary conditions

Parameters:

• N (int) – Number of quadrature points in the homogeneous space.

• bc (dict) – The boundary conditions in dictionary form, see BoundaryConditions.

• domain (2-tuple of numbers, optional) – The domain of the homogeneous space.

• alpha (number, optional) – Parameter of the Jacobi polynomial

• beta (number, optional) – Parameter of the Jacobi polynomial

basis_function(i=0, x=x)

Return basis function i

Parameters:

• i (int, optional) – The degree of freedom of the basis function

• x (sympy Symbol, optional)

static boundary_condition()

property dim_ortho

eval(x, u, output_array=None)

Evaluate Function u at position x

Parameters:

• x (float or array of floats)

• u (array) – Expansion coefficients or instance of Function

• output_array (array, optional) – Function values at points

Returns:

output_array

Return type:

array

evaluate_basis(x, i=0, output_array=None)

Evaluate basis function i at points x

Parameters:

• x (float or array of floats)

• i (int, optional) – Basis function number

• output_array (array, optional) – Return result in output_array if provided

Returns:

output_array

Return type:

array

evaluate_basis_derivative(x=None, i=0, k=0, output_array=None)

Evaluate k’th derivative of basis function i at x or all quadrature points

Parameters:

• x (float or array of floats, optional) – If not provided use quadrature points of self

• i (int, optional) – Basis function number

• k (int, optional) – k’th derivative

• output_array (array, optional) – return array

Returns:

output_array

Return type:

array

get_orthogonal(**kwargs)

Return orthogonal space (otherwise as self)

Returns:

The orthogonal space in the same family, and otherwise as self.

Return type:

SpectralBase

property is_boundary_basis

shape(forward_output=True)

Return the allocated shape of arrays used for self

Parameters:

forward_output (bool, optional) – If True then return allocated shape of spectral space (the result of a forward transform). If False then return allocated shape of physical space (the input to a forward transform).

static short_name()

slice()

Return index set of current space

stencil_matrix(N=None)

Return stencil matrix in SparseMatrix format

Parameters:

N (int or None, optional) – Shape (N, N) of the stencil matrix. Using self.N if None

to_ortho(input_array, output_array=None)

Project to orthogonal basis

Parameters:

• input_array (array) – Expansion coefficients of input basis

• output_array (array, optional) – Expansion coefficients in orthogonal basis

Returns:

output_array

Return type:

array

vandermonde(x)

Return Vandermonde matrix based on the primary (orthogonal) basis of the family.

Evaluates basis function \(\psi_k(x)\) for all wavenumbers, and all x. Returned Vandermonde matrix is an N x M matrix with N the length of x and M the number of bases.

\[\begin{split}\begin{bmatrix} \psi_0(x_0) & \psi_1(x_0) & \ldots & \psi_{M-1}(x_0)\\ \psi_0(x_1) & \psi_1(x_1) & \ldots & \psi_{M-1}(x_1)\\ \vdots & \ldots \\ \psi_{0}(x_{N-1}) & \psi_1(x_{N-1}) & \ldots & \psi_{M-1}(x_{N-1}) \end{bmatrix}\end{split}\]

Parameters:

x (array of floats) – points for evaluation

Note

This function returns a matrix of evaluated orthogonal basis functions for either family. The true Vandermonde matrix of a basis is obtained through SpectralBase.evaluate_basis_all().

class shenfun.legendre.bases.BeamFixedFree(N, quad='LG', bc=(0, 0, 0, 0), domain=(-1, 1), padding_factor=1, dealias_direct=False, dtype=<class 'float'>, coordinates=None, **kw)

Bases: CompositeBase

Function space for fixed free beams

The basis \(\{\phi_k\}_{k=0}^{N-1}\) is

\[\begin{split}\phi_k &= L_k + a_kL_{k+1} + b_k L_{k+2} + c_k L_{k+3} + d_k L_{k+4} , \, k=0, 1, \ldots, N-5, \\ \phi_{N-4} &= \tfrac{1}{2}L_0-\tfrac{3}{5}L_1+\tfrac{1}{10}L_3, \\ \phi_{N-3} &= \tfrac{1}{6}L_0-\tfrac{1}{10}L_1-\tfrac{1}{6}L_2+\tfrac{1}{10}L_3, \\ \phi_{N-2} &= \tfrac{1}{2}L_0+\tfrac{3}{5}L_1-\tfrac{1}{10}L_3), \\ \phi_{N-1} &= -\tfrac{1}{6}L_0-\tfrac{1}{10}L_1+\tfrac{1}{6}L_2+\tfrac{1}{10}L_3,\end{split}\]

where

\[\begin{split}a_k &= \frac{4 \left(2 n + 3\right)}{\left(n + 3\right)^{2}}, \\ b_k &= -\frac{2 \left(n - 1\right) \left(n + 1\right) \left(n + 6\right) \left(2 n + 5\right)}{\left(n + 3\right)^{2} \left(n + 4\right) \left(2 n + 7\right)}, \\ c_k &= -\frac{4 \left(n + 1\right)^{2} \left(2 n + 3\right)}{\left(n + 3\right)^{2} \left(n + 4\right)^{2}}, \\ d_k &= \frac{\left(n + 1\right)^{2} \left(n + 2\right)^{2} \left(2 n + 3\right)}{\left(n + 3\right)^{2} \left(n + 4\right)^{2} \left(2 n + 7\right)}.\end{split}\]

We have

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(-1)&=a, u'(-1) = b, u''(1)=c, u'''(1) = d.\end{split}\]

The last four basis functions are for boundary conditions and only used if a, b, c or d are different from 0. In one dimension \(\hat{u}_{N-4}=a\), \(\hat{u}_{N-3}=b\), \(\hat{u}_{N-2}=c\) and \(\hat{u}_{N-1}=d\).

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• 4-tuple of numbers, optional – The values of the 4 boundary conditions u(-1) = a, u’(-1) = b, u’’(1) = c, u’’’(1) = d

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.DirichletNeumann(N, quad='LG', bc=(0, 0), domain=(-1, 1), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space for mixed Dirichlet/Neumann boundary conditions

The basis \(\{\phi_k\}_{k=0}^{N-1}\) is

\[\begin{split}\phi_k &= L_{k} + \frac{2n+3}{\left(n+2\right)^{2}}L_{k+1} - \frac{\left(n+1\right)^{2}}{\left(n+2\right)^{2}} L_{k+2}, \, k=0, 1, \ldots, N-3, \\ \phi_{N-2} &= L_0, \\ \phi_{N-1} &= L_0+L_1,\end{split}\]

such that

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(-1) &= a, u'(1)=b.\end{split}\]

The last two bases are for boundary conditions and only used if a or b are different from 0. In one dimension \(\hat{u}_{N-2}=a\) and \(\hat{u}_{N-1}=b\).

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (tuple of numbers) – Boundary conditions at edges of domain. Dirichlet first.

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.Generic(N, quad='LG', bc={}, domain=(-1, 1), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space for space with any boundary conditions

Any combination of Dirichlet and Neumann is possible.

Parameters:

• N (int, optional) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (dict, optional) – The dictionary must have keys ‘left’ and ‘right’, to describe boundary conditions on the left and right boundaries. Specify Dirichlet on both ends with

  {‘left’: {‘D’: a}, ‘right’: {‘D’: b}}

  for some values a and b, that will be neglected in the current function. Specify mixed Neumann and Dirichlet as

  {‘left’: {‘N’: a}, ‘right’: {‘N’: b}}

  For both conditions on the right do

  {‘right’: {‘N’: a, ‘D’: b}}

  Any combination should be possible, and it should also be possible to use second derivatives N2. See BoundaryConditions.

• domain (2-tuple of numbers, optional) – The computational domain

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

Note

A test function is always using homogeneous boundary conditions.

static boundary_condition()

static short_name()

class shenfun.legendre.bases.LowerDirichlet(N, quad='LG', bc=(0, None), domain=(-1, 1), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space with single Dirichlet boundary condition

The basis \(\{\phi_k\}_{k=0}^{N-1}\) is

\[\begin{split}\phi_k &= L_{k} + L_{k+1}, \, k=0, 1, \ldots, N-2, \\ \phi_{N-1} &= L_0,\end{split}\]

such that

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(-1) &= a.\end{split}\]

The last basis function is for boundary condition and only used if a is different from 0. In one dimension \(\hat{u}_{N-1}=a\).

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (tuple of (number, None)) – Boundary conditions at edges of domain.

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.NeumannDirichlet(N, quad='LG', bc=(0, 0), domain=(-1, 1), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space for mixed Neumann/Dirichlet boundary conditions

The basis \(\{\phi_k\}_{k=0}^{N-1}\) is

\[\begin{split}\phi_k &= L_{k} - \frac{2n+3}{\left(n+2\right)^{2}}L_{k+1} - \frac{\left(n+1\right)^{2}}{\left(n+2\right)^{2}}L_{k+2}, \, k=0, 1, \ldots, N-3, \\ \phi_{N-2} &= -L_0+L_1, \\ \phi_{N-1} &= L_0,\end{split}\]

such that

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u'(-1) &= a, u(1)=b.\end{split}\]

The last two bases are for boundary conditions and only used if a or b are different from 0. In one dimension \(\hat{u}_{N-2}=a\) and \(\hat{u}_{N-1}=b\).

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (tuple of numbers) – Boundary conditions at edges of domain. Neumann first.

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.Orthogonal(N, quad='LG', domain=(-1, 1), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: JacobiBase

Function space for a regular Legendre series

The orthogonal basis is

\[L_k, \quad k = 0, 1, \ldots, N-1,\]

where \(L_k\) is the \(k\)’th Legendre polynomial.

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

L2_norm_sq(i)

Return square of L2-norm

\[\| \phi_i \|^2_{\omega} = (\phi_i, \phi_i)_{\omega} = \int_{I} \phi_i \overline{\phi}_i \omega dx\]

where \(\phi_i\) is the i’th orthogonal basis function for the orthogonal basis in the given family.

Parameters:

i (int) – The number of the orthogonal basis function

eval(x, u, output_array=None)

Evaluate Function u at position x

Parameters:

• x (float or array of floats)

• u (array) – Expansion coefficients or instance of Function

• output_array (array, optional) – Function values at points

Returns:

output_array

Return type:

array

evaluate_basis(x, i=0, output_array=None)

Evaluate basis function i at points x

Parameters:

• x (float or array of floats)

• i (int, optional) – Basis function number

• output_array (array, optional) – Return result in output_array if provided

Returns:

output_array

Return type:

array

evaluate_basis_derivative(x=None, i=0, k=0, output_array=None)

Evaluate k’th derivative of basis function i at x or all quadrature points

Parameters:

• x (float or array of floats, optional) – If not provided use quadrature points of self

• i (int, optional) – Basis function number

• k (int, optional) – k’th derivative

• output_array (array, optional) – return array

Returns:

output_array

Return type:

array

evaluate_basis_derivative_all(x=None, k=0, argument=0)

Return k’th derivative of basis evaluated at x or all quadrature points as a Vandermonde matrix.

Parameters:

• x (float or array of floats, optional) – If not provided use quadrature points of self

• k (int, optional) – k’th derivative

• argument (int) – Zero for test and 1 for trialfunction

Returns:

Vandermonde matrix

Return type:

array

static family()

get_bc_space()

get_orthogonal(**kwargs)

Return orthogonal space (otherwise as self)

Returns:

The orthogonal space in the same family, and otherwise as self.

Return type:

SpectralBase

get_recursion_matrix(M, N)

property is_orthogonal

l2_norm_sq(i=None)

Return square of l2-norm

\[\| u \|^2_{N,\omega} = (u, u)_{N,\omega} = \sun_{j=0}^{N-1} u(x_j) \overline{u}(x_j) \omega_j\]

where \(u=\{\phi_i\}_{i=0}^{N-1}\) and \(\phi_i\) is the i’th orthogonal basis function in the given family.

Parameters:

i (None or int) – If None then return the square of the l2-norm for all i=0, 1, …, N-1. Else, return for given i.

orthogonal_basis_function(i=0, x=x)

Return the orthogonal basis function i

Parameters:

• i (int, optional) – The degree of freedom of the basis function

• x (sympy Symbol, optional)

plan(shape, axis, dtype, options)

Plan transform

Allocate work arrays for transforms and set up methods forward, backward and scalar_product with or without padding

Parameters:

• shape (array) – Local shape of global array

• axis (int) – This base’s axis in global TensorProductSpace

• dtype (numpy.dtype) – Type of array

• options (dict) – Options for planning transforms

points_and_weights(N=None, map_true_domain=False, weighted=True, **kw)

Return points and weights of quadrature for weighted integral

\[\int_{\Omega} f(x) w(x) dx \approx \sum_{i} f(x_i) w_i\]

Parameters:

• N (int, optional) – Number of quadrature points

• map_true_domain (bool, optional) – Whether or not to map points to true domain

• weighted (bool, optional) – Whether to use quadrature weights for a weighted inner product (default), or a regular, non-weighted inner product.

Note

The weight of the space is included in the returned quadrature weights.

static short_name()

stencil_matrix(N=None)

sympy_stencil(i=i, j=j)

to_chebyshev(input_array, output_array=None)

to_ortho(input_array, output_array=None)

Project to orthogonal basis

Parameters:

• input_array (array) – Expansion coefficients of input basis

• output_array (array, optional) – Expansion coefficients in orthogonal basis

Returns:

output_array

Return type:

array

vandermonde(x)

Return Vandermonde matrix based on the primary (orthogonal) basis of the family.

Evaluates basis function \(\psi_k(x)\) for all wavenumbers, and all x. Returned Vandermonde matrix is an N x M matrix with N the length of x and M the number of bases.

\[\begin{split}\begin{bmatrix} \psi_0(x_0) & \psi_1(x_0) & \ldots & \psi_{M-1}(x_0)\\ \psi_0(x_1) & \psi_1(x_1) & \ldots & \psi_{M-1}(x_1)\\ \vdots & \ldots \\ \psi_{0}(x_{N-1}) & \psi_1(x_{N-1}) & \ldots & \psi_{M-1}(x_{N-1}) \end{bmatrix}\end{split}\]

Parameters:

x (array of floats) – points for evaluation

Note

This function returns a matrix of evaluated orthogonal basis functions for either family. The true Vandermonde matrix of a basis is obtained through SpectralBase.evaluate_basis_all().

class shenfun.legendre.bases.Phi1(N, quad='LG', bc=(0, 0), domain=(-1, 1), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space for Dirichlet boundary conditions

The basis \(\{\phi_k\}_{k=0}^{N-1}\) is

\[\begin{split}\phi_k &= \frac{1}{2}(L_k - L_{k+2}) = \frac{(2k+3)(1-x^2)}{2(k+1)(k+2)} L'_{k+1}, \, k=0, 1, \ldots, N-3, \\ \phi_{N-2} &= \frac{1}{2}(L_0-L_1), \\ \phi_{N-1} &= \frac{1}{2}(L_0+L_1),\end{split}\]

such that

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(-1) &= a \text{ and } u(1) = b.\end{split}\]

The last two basis functions are for boundary conditions and only used if a or b are different from 0. In one dimension \(\hat{u}_{N-2}=a\) and \(\hat{u}_{N-1}=b\).

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (tuple of numbers) – Boundary conditions at edges of domain

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.Phi2(N, quad='LG', bc=(0, 0, 0, 0), domain=(-1, 1), padding_factor=1, dealias_direct=False, dtype=<class 'float'>, coordinates=None, **kw)

Bases: CompositeBase

Function space for biharmonic equation

The basis functions \(\phi_k\) for \(k=0, 1, \ldots, N-5\) are

\[\begin{split}\phi_k &= \frac{(1-x^2)^2 L''_{k+2}}{h^{(2)}_{k+2}}, \\ h^{(2)}_{k+2} &= \int_{-1}^1 L''_{k+2} L''_{k+2} (1-x^2)^2 dx, \\ &= \frac{2 (k+1)(k+2)(k+3)(k+4)}{2k+5},\end{split}\]

which (along with boundary functions) becomes the basis

\[\begin{split}\phi_k &= \frac{1}{2(2k+3)}\left(L_k - \frac{2(2k+5)}{2k+7}L_{k+2} + \frac{2k+3}{2k+7}L_{k+4}\right), \, k=0, 1, \ldots, N-5, \\ \phi_{N-4} &= \tfrac{1}{2}L_0-\tfrac{3}{5}L_1+\tfrac{1}{10}L_3, \\ \phi_{N-3} &= \tfrac{1}{6}L_0-\tfrac{1}{10}L_1-\tfrac{1}{6}L_2+\tfrac{1}{10}L_3, \\ \phi_{N-2} &= \tfrac{1}{2}L_0+\tfrac{3}{5}L_1-\tfrac{1}{10}L_3, \\ \phi_{N-1} &= -\tfrac{1}{6}L_0-\tfrac{1}{10}L_1+\tfrac{1}{6}L_2+\tfrac{1}{10}L_3,\end{split}\]

such that

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(-1)&=a, u'(-1) = b, u(1)=c, u'(1) = d.\end{split}\]

The last four basis functions are for boundary conditions and only used if a, b, c or d are different from 0. In one dimension \(\hat{u}_{N-4}=a\), \(\hat{u}_{N-3}=b\), \(\hat{u}_{N-2}=c\) and \(\hat{u}_{N-1}=d\).

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (4-tuple of numbers, optional) – The values of the 4 boundary conditions at x=(-1, 1). The two on x=-1 first and then x=1. (a, b, c, d)

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.Phi3(N, quad='LG', bc=(0, 0, 0, 0, 0, 0), domain=(-1, 1), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space for 6th order equations

The basis functions \(\phi_k\) for \(k=0, 1, \ldots, N-7\) are

\[\begin{split}\phi_k &= \frac{(1-x^2)^3}{h^{(3)}_{k+3}} L^{(3)}_{k+3}}, \, k=0, 1, \ldots, N-7, \\ h^{(3)}_{k+3} &= \frac{2\Gamma(k+7)}{\Gamma(k+1)(2k+7)} = \int_{-1}^1 L^{(3)}_{k+3} L^{(3)}_{k+3}(1-x^2)^3 dx,\end{split}\]

where \(L^{(3)}_k\) is the 3’rd derivative of \(L_k\). The 6 boundary basis functions are computed using jacobi.findbasis.get_bc_basis(), but they are too messy to print here. We have

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(-1) &= a, u'(-1)=b, u''(-1)=c, u(1)=d u'(1)=e, u''(1)=f.\end{split}\]

The last 6 basis functions are for boundary conditions and only used if there are nonzero boundary conditions.

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature - LG - Legendre-Gauss - GL - Legendre-Gauss-Lobatto

• bc (6-tuple of numbers, optional) – Boundary conditions.

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.Phi4(N, quad='LG', bc=(0, 0, 0, 0, 0, 0, 0, 0), domain=(-1, 1), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space with 2 Dirichlet and 6 Neumann boundary conditions

The basis functions \(\phi_k\) for \(k=0, 1, \ldots, N-9\) are

\[\begin{split}\phi_k &= \frac{(1-x^2)^4}{h^{(4)}_{k+4}} L^{(4)}_{k+4}, \\ h^{(4)}_{k+4} &= \frac{2\Gamma(k+9)}{\Gamma(k+1)(2k+9)} = \int_{-1}^1 L^{(4)}_{k+4} L^{(4)}_{k+4} (1-x^2)^4 dx,\end{split}\]

where \(L^{(4)}_k\) is the 4’th derivative of \(L_k\). The boundary basis for inhomogeneous boundary conditions is too messy to print, but can be obtained using get_bc_basis(). We have

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(-1) &= a, u'(-1)=b, u''(-1)=c, u'''(-1)=d, u(1)=e u'(1)=f, u''(1)=g, u'''(1)=h.\end{split}\]

The last 8 basis functions are for boundary conditions and only used if there are nonzero boundary conditions.

Parameters:

• N (int, optional) – Number of quadrature points

• quad (str, optional) – Type of quadrature - LG - Legendre-Gauss - GL - Legendre-Gauss-Lobatto

• bc (8-tuple of numbers)

• domain (2-tuple of numbers, optional) – The computational domain

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.Phi6(N, quad='LG', bc=(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), domain=(-1, 1), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space for 12th order equation

The basis functions \(\phi_k\) for \(k=0, 1, \ldots, N-9\) are

\[\begin{split}\phi_k &= \frac{(1-x^2)^6}{h^{(6)}_{k+6}} L^{(6)}_{k+6}, \\ h^{(6)}_{k+6} &= \int_{-1}^1 L^{(6)}_{k+6} L^{(6)}_{k+6} (1-x^2)^6 dx,\end{split}\]

where \(L^{(6)}_k\) is the 6’th derivative of \(L_k\). The boundary basis for inhomogeneous boundary conditions is too messy to print, but can be obtained using get_bc_basis().

Parameters:

• N (int, optional) – Number of quadrature points

• quad (str, optional) – Type of quadrature - LG - Legendre-Gauss - GL - Legendre-Gauss-Lobatto

• bc (12-tuple of numbers)

• domain (2-tuple of numbers, optional) – The computational domain

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.ShenBiPolar(N, quad='LG', domain=(-1, 1), bc=(0, 0, 0, 0), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space for the Biharmonic equation

The basis \(\{\phi_k\}_{k=0}^{N-1}\) is

\[\begin{split}\phi_k &= (1-x^2)^2 L'_{k+1}, \quad k=0,1, \ldots, N-5, \\ \phi_{N-4} &= \tfrac{1}{2}L_0-\tfrac{3}{5}L_1+\tfrac{1}{10}L_3, \\ \phi_{N-3} &= \tfrac{1}{6}L_0-\tfrac{1}{10}L_1-\tfrac{1}{6}L_2+\tfrac{1}{10}L_3, \\ \phi_{N-2} &= \tfrac{1}{2}L_0+\tfrac{3}{5}L_1-\tfrac{1}{10}L_3), \\ \phi_{N-1} &= -\tfrac{1}{6}L_0-\tfrac{1}{10}L_1+\tfrac{1}{6}L_2+\tfrac{1}{10}L_3,\end{split}\]

such that

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(-1)&=a, u'(-1) = b, u(1)=c, u'(1) = d.\end{split}\]

The last four bases are for boundary conditions and only used if a, b, c or d are different from 0. In one dimension \(\hat{u}_{N-4}=a\), \(\hat{u}_{N-3}=b\), \(\hat{u}_{N-2}=c\) and \(\hat{u}_{N-1}=d\).

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (4-tuple of numbers, optional) – The values of the 4 boundary conditions at x=(-1, 1). The two on x=-1 first and then x=1. (a, b, c, d)

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

stencil_matrix(N=None)

Return stencil matrix in SparseMatrix format

Parameters:

N (int or None, optional) – Shape (N, N) of the stencil matrix. Using self.N if None

sympy_stencil(i=i, j=j)

Return stencil matrix as a Sympy matrix

Parameters:

• i, j (Sympy symbols) – indices for row and column

• implicit (bool or str, optional) – Whether to use an unevaluated Sympy function instead of the actual value of the stencil. This makes the matrix prettier, and it can still be evaluated. If implicit is not False, then it must be a string of length one. This string represents the diagonals of the matrix. If implicit=’a’ and the stencil matrix has two diagonals in the main diagonal and the first upper diagonal, then these are called ‘a0’ and ‘a1’.

Example

>>> from shenfun import FunctionSpace
>>> import sympy as sp
>>> i, j = sp.symbols('i,j', integer=True)
>>> D = FunctionSpace(8, 'L', bc=(0, 0), scaled=True)
>>> D._stencil
{0: 1/sqrt(4*n + 6), 2: -1/sqrt(4*n + 6)}
>>> D.sympy_stencil()
KroneckerDelta(i, j)/sqrt(4*i + 6) - KroneckerDelta(j, i + 2)/sqrt(4*i + 6)
>>> D.sympy_stencil(implicit='a')
KroneckerDelta(i, j)*a0(i) + KroneckerDelta(j, i + 2)*a2(i)

Get the main diagonal

>>> D.sympy_stencil(implicit=False).subs(j, i)
1/sqrt(4*i + 6)

class shenfun.legendre.bases.ShenBiharmonic(N, quad='LG', bc=(0, 0, 0, 0), domain=(-1, 1), padding_factor=1, dealias_direct=False, dtype=<class 'float'>, coordinates=None, **kw)

Bases: CompositeBase

Function space for biharmonic equation

The basis \(\{\phi_k\}_{k=0}^{N-1}\) is

\[\begin{split}\phi_k &= L_n - \frac{4n+10}{2n+7}L_{n+2}+\frac{2 n + 3}{2 n + 7}L_{n+4}, \, k=0, 1, \ldots, N-5, \\ \phi_{N-4} &= \tfrac{1}{2}L_0-\tfrac{3}{5}L_1+\tfrac{1}{10}L_3, \\ \phi_{N-3} &= \tfrac{1}{6}L_0-\tfrac{1}{10}L_1-\tfrac{1}{6}L_2+\tfrac{1}{10}L_3, \\ \phi_{N-2} &= \tfrac{1}{2}L_0+\tfrac{3}{5}L_1-\tfrac{1}{10}L_3), \\ \phi_{N-1} &= -\tfrac{1}{6}L_0-\tfrac{1}{10}L_1+\tfrac{1}{6}L_2+\tfrac{1}{10}L_3,\end{split}\]

such that

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(-1)&=a, u'(-1) = b, u(1)=c, u'(1) = d.\end{split}\]

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (4-tuple of numbers, optional) – The values of the 4 boundary conditions at x=(-1, 1). The two conditions on x=-1 first, and then x=1. With (a, b, c, d) corresponding to bc = {‘left’: [(‘D’, a), (‘N’, b)], ‘right’: [(‘D’, c), (‘N’, d)]}

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.ShenDirichlet(N, quad='LG', bc=(0, 0), domain=(-1, 1), dtype=<class 'float'>, scaled=False, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space for Dirichlet boundary conditions

The basis \(\{\phi_k\}_{k=0}^{N-1}\) is

\[\begin{split}\phi_k &= L_k - L_{k+2}, \, k=0, 1, \ldots, N-3, \\ \phi_{N-2} &= \frac{1}{2}(L_0-L_1), \\ \phi_{N-1} &= \frac{1}{2}(L_0+L_1),\end{split}\]

such that

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(-1) &= a \text{ and } u(1) = b.\end{split}\]

The last two basis functions are for boundary conditions and only used if a or b are different from 0. In one dimension \(\hat{u}_{N-2}=a\) and \(\hat{u}_{N-1}=b\).

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (tuple of numbers) – Boundary conditions at edges of domain

• domain (2-tuple of numbers, optional) – The computational domain

• scaled (bool, optional) – Whether or not to scale test functions with 1/sqrt(4k+6). Scaled test functions give a stiffness matrix equal to the identity matrix.

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.ShenNeumann(N, quad='LG', bc=(0, 0), domain=(-1, 1), padding_factor=1, dealias_direct=False, dtype=<class 'float'>, coordinates=None, **kw)

Bases: CompositeBase

Function space for Neumann boundary conditions

The basis \(\{\phi_k\}_{k=0}^{N-1}\) is

\[\begin{split}\phi_k &= L_{k} - \frac{k(k+1)}{(k+2)(k+3)}L_{k+2}, \, k=0, 1, \ldots, N-3, \\ \phi_{N-2} &= \frac{1}{6}(3L_1-L_2), \\ \phi_{N-1} &= \frac{1}{6}(3L_1+L_2),\end{split}\]

such that

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u'(-1) &= a \text{ and } u'(1) = b.\end{split}\]

The last two bases are for boundary conditions and only used if a or b are different from 0. In one dimension \(\hat{u}_{N-2}=a\) and \(\hat{u}_{N-1}=b\).

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (2-tuple of numbers) – Boundary conditions at edges of domain

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.UpperDirichlet(N, quad='LG', bc=(None, 0), domain=(-1, 1), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space with single Dirichlet on upper edge

The basis \(\{\phi_k\}_{k=0}^{N-1}\) is

\[\begin{split}\phi_k &= L_{k} - L_{k+1}, \, k=0, 1, \ldots, N-2, \\ \phi_{N-1} &= L_0,\end{split}\]

such that

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(1) &= a.\end{split}\]

The last basis function is for boundary condition and only used if a is different from 0. In one dimension \(\hat{u}_{N-1}=a\).

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (2-tuple of (None, number), optional) – The number is the boundary condition value

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

static boundary_condition()

static short_name()

class shenfun.legendre.bases.UpperDirichletNeumann(N, quad='LG', bc=(0, 0), domain=(-1, 1), dtype=<class 'float'>, padding_factor=1, dealias_direct=False, coordinates=None, **kw)

Bases: CompositeBase

Function space for both Dirichlet and Neumann boundary conditions on the right hand side.

The basis \(\{\phi_k\}_{k=0}^{N-1}\) is

\[\begin{split}\phi_k &= L_{k} - \frac{2k+3}{k+2}L_{k+1} + \frac{k+1}{k+2}L_{k+2}, \, k=0, 1, \ldots, N-3, \\ \phi_{N-2} &= L_0, \\ \phi_{N-1} &= -L_0+L_1,\end{split}\]

such that

\[\begin{split}u(x) &= \sum_{k=0}^{N-1} \hat{u}_k \phi_k(x), \\ u(1) &= a, u'(1)=b.\end{split}\]

The last two bases are for boundary conditions and only used if a or b are different from 0. In one dimension \(\hat{u}_{N-2}=a\) and \(\hat{u}_{N-1}=b\).

Parameters:

• N (int) – Number of quadrature points

• quad (str, optional) – Type of quadrature

  • LG - Legendre-Gauss

  • GL - Legendre-Gauss-Lobatto

• bc (tuple of numbers) – Boundary conditions at edges of domain, Dirichlet first.

• domain (2-tuple of numbers, optional) – The computational domain

• padding_factor (float, optional) – Factor for padding backward transforms.

• dealias_direct (bool, optional) – Set upper 1/3 of coefficients to zero before backward transform

• dtype (data-type, optional) – Type of input data in real physical space. Will be overloaded when basis is part of a TensorProductSpace.

• coordinates (2- or 3-tuple (coordinate, position vector (, sympy assumptions)), optional) – Map for curvilinear coordinatesystem, and parameters to Coordinates

Note

This basis is not recommended as it leads to a poorly conditioned stiffness matrix.

static boundary_condition()

static short_name()

shenfun.legendre.la module

class shenfun.legendre.la.Helmholtz_2dirichlet(matrices)

Bases: object

Helmholtz solver for 2-dimensional problems with 2 Dirichlet bases.

\[a_0 BUB + a_1 AUB + a_2 BUA^T = F\]

Somewhat experimental.

__call__(b, u, solver=1)

Call self as a function.

solve_eigen_problem(A, B, solver)

Solve the eigen problem

shenfun.legendre.dlt module

class shenfun.legendre.dlt.Cheb2Leg

Bases: FMMLevel

__call__(*args, **kwargs)

Call self as a function.

apply(input_array, output_array)

class shenfun.legendre.dlt.DLT(input_array, s=None, axes=(-1,), threads=1, kind='backward', flags=(0, 16), output_array=None)

Bases: object

Discrete Legendre Transform

A class for performing fast FFT-based discrete Legendre transforms, both forwards and backwards. Based on:

Nicholas Hale and Alex Townsend "A fast FFT-based discrete Legendre
transform", IMA Journal of Numerical Analysis (2015)
(https://arxiv.org/abs/1505.00354)

Parameters:

• input_array (real or complex array)

• s (sequence of ints, optional) – Not used - included for compatibility with Numpy

• axes (integer or 1-tuple of int, optional) – Axis over which to compute the DLT. Named axes for compatibility.

• threads (int, optional) – Number of threads used in computing DLT.

• kind (str, optional) –

  Either one of

  • ‘forward’

  • ‘backward’

  • ‘scalar product’

  The scalar product is exactly like the forward transform, except that the Legendre mass matrix is not applied to the output.

• flags (sequence of ints, optional) –

  Flags from

  • FFTW_MEASURE

  • FFTW_EXHAUSTIVE

  • FFTW_PATIENT

  • FFTW_DESTROY_INPUT

  • FFTW_PRESERVE_INPUT

  • FFTW_UNALIGNED

  • FFTW_CONSERVE_MEMORY

  • FFTW_ESTIMATE

• output_array (complex array, optional) – Array to be used as output array. Must be of correct shape, type, strides and alignment

Note

The Legendre series expansion of \(f(x)\) is

\[f(x) = \sum_{k=0}^{N-1} \hat{f}_k L_k(x)\]

The series evaluated at quadrature points \(\{x_j\}_{j=0}^{N-1}\) gives the vector \(\{f(x_j)\}_{j=0}^{N-1}\). We define a forward transform as computing \(\{\hat{f}_k\}_{k=0}^{N-1}\) from \(\{f(x_j)\}_{j=0}^{N-1}\). The other way around, a backward transform, computes \(\{f(x_j)\}_{j=0}^{N-1}\) given \(\{\hat{f}_k\}_{k=0}^{N-1}\). This is in agreement with shenfun’s definition of forward/backward directions, but it disagrees with the definitions used by Hale and Townsend.

Also note that this fast transform is actually slower than the default recursive version for approximately \(N<1000\).

Example

>>> import numpy as np
>>> from shenfun import legendre, FunctionSpace, Function, Array
>>> N = 4
>>> L = FunctionSpace(N, 'L')
>>> u = Function(L)
>>> u[1] = 1
>>> c = Array(L)
>>> c = L.backward(u, c, kind='fast')
>>> print(c)
[-0.86113631 -0.33998104  0.33998104  0.86113631]
>>> u[:] = 1
>>> np.alltrue(np.abs(u - u.backward().forward()) < 1e-12)
True

__call__(input_array=None, output_array=None, **kw)

Call self as a function.

property input_array

property output_array

plan(U, V, kind, threads, flags)

class shenfun.legendre.dlt.FMMCheb2Leg(input, output_array=None, diagonals=16, domains=None, levels=None, maxs=100, axis=0, use_direct=-1)

Bases: FMMLevel

Transform Chebyshev coefficients to Legendre coefficients

Parameters:

• input (int or input array) – The length of the array to transform or the array itself

• diagonals (int) – The number of neglected diagonals, that are treated using a direct approach

• domains (None, int or sequence of ints) – The domain sizes for all levels

  If domains=None then an appropriate domain size is computed according to the given number of levels and maxs. If the number of levels is also None, then domains is set to 2 and the number of levels is computed from maxs.

  If domains is an integer, then this integer is used for each level, and the number of levels is either given or computed according to maxs

  If domains is a sequence of integers, then these are the domain sizes for all the levels and the length of this sequence is the number of levels.

• levels (None or int) – The number of levels in the hierarchical matrix

  If levels is None, then it is computed according to domains and maxs

• l2c (bool) – If True, the transform goes from Legendre to Chebyshev, and vice versa if False

• maxs (int) – The maximum size of the smallest submatrices (on the highest level). This number is used if the number of levels or domain sizes need to be computed.

• axis (int) – The axis over which to apply the transform if the input array is a multidimensional array.

• use_direct (int) – Use direct method if N is smaller than this number

__call__(input_array, output_array=None)

Execute transform

Parameters:

• input_array (array) – The array to be transformed

• output_array (None or array) – The return array. Will be created if None.

apply(input_array, output_array)

class shenfun.legendre.dlt.FMMLeg2Cheb(input, output_array=None, domains=None, levels=None, maxs=100, axis=0, use_direct=-1)

Bases: FMMLevel

Transform Legendre coefficients to Chebyshev coefficients

Parameters:

• input (int or input array) – The length of the array to transform or the array itself

• diagonals (int) – The number of neglected diagonals, that are treated using a direct approach

• domains (None, int or sequence of ints) – The domain sizes for all levels

  If domains=None then an appropriate domain size is computed according to the given number of levels and maxs. If the number of levels is also None, then domains is set to 2 and the number of levels is computed from maxs.

  If domains is an integer, then this integer is used for each level, and the number of levels is either given or computed according to maxs

  If domains is a sequence of integers, then these are the domain sizes for all the levels and the length of this sequence is the number of levels.

• levels (None or int) – The number of levels in the hierarchical matrix

  If levels is None, then it is computed according to domains and maxs

• l2c (bool) – If True, the transform goes from Legendre to Chebyshev, and vice versa if False

• maxs (int) – The maximum size of the smallest submatrices (on the highest level). This number is used if the number of levels or domain sizes need to be computed.

• axis (int) – The axis over which to apply the transform if the input array is a multidimensional array.

• use_direct (int) – Use direct method if N is smaller than this number

__call__(input_array, output_array=None, transpose=False)

Execute transform

Parameters:

• input_array (array) – The array to be transformed

• output_array (None or array) – The return array. Will be created if None.

• transpose (bool) – If True, then apply the transpose operation \(A^Tu\) instead of the default \(Au\), where \(A\) is the matrix and \(u\) is the array of Legendre coefficients.

apply(input_array, output_array, transpose)

class shenfun.legendre.dlt.Leg2Cheb

Bases: FMMLevel

__call__(*args, **kwargs)

Call self as a function.

class shenfun.legendre.dlt.Leg2chebHaleTownsend(input_array, output_array=None, axis=0, nM=50, Nmin=40000)

Bases: object

Class for computing Chebyshev coefficients from Legendre coefficients

Algorithm from:

Nicholas Hale and Alex Townsend "A fast, simple and stable Chebyshev-
Legendre transform using an asymptotic formula", SIAM J Sci Comput (2014)
(https://epubs.siam.org/doi/pdf/10.1137/130932223)

Parameters:

• input_array (array) – Legendre coefficients

• output_array (array) – The returned array

• axis (int) – The axis over which to perform the computation in case the input_array is multidimensional.

• nM (int) – Parameter, see Hale and Townsend (2014). Note that one must have N >> nM.

• Nmin (int) – Parameter. Choose direct matvec approach for N < Nmin

__call__(input_array=None, output_array=None, transpose=False)

Compute Chebyshev coefficients from Legendre. That is, compute

\[\hat{c}^{cheb} = M \hat{c}^{leg}\]

where \(\hat{c}^{cheb} \in \mathbb{R}^N\) are the Chebyshev coefficients, \(\hat{c}^{leg} \in \mathbb{R}^N\) the Legendre coefficients and \(M\in\mathbb{R}^{N \times N}\) the matrix for the conversion. Note that if keyword ‘transpose’ is true, then we compute

\[\hat{a} = M^T \hat{b}\]

for some vectors \(\hat{a}\) and \(\hat{b}\).

The Chebyshev and Legendre coefficients are the regular coefficients to the series

\[\begin{split}p_l(x) = \sum_{k=0}^{N} \hat{c}_k^{leg}L_k(x) \\ p_c(x) = \sum_{k=0}^{N} \hat{c}_k^{cheb}T_k(x)\end{split}\]

and we get \(\{\hat{c}_k^{cheb}\}_{k=0}^N\) by setting \(p_l(x)=p_c(x)\) for \(x=\{x_i\}_{i=0}^N\), where \(x_i=\cos(i+0.5)\pi/N\).

Parameters:

• input_array (array)

• output_array (array)

• transpose (bool) – Whether to compute the transpose of the regular transform

Note

For small N we use a direct method that costs approximately \(0.25 N^2\) operations. For larger N (see ‘Nmin’ parameter) we use the fast routine of

Hale and Townsend ‘A fast, simple and stable Chebyshev-Legendre transform using an asymptotic formula’, SIAM J Sci Comput (2014)

property input_array

property output_array

plan(input_array)

shenfun.legendre.lobatto module

Module contains some useful methods for Legendre quadrature.

shenfun.legendre.lobatto.legendre_lobatto_nodes_and_weights(N, tol=1e-16)

Return points and weights for Legendre-Lobatto quadrature

Parameters:

N (int) – Number of quadrature points

shenfun.legendre.matrices module

This module contains specific inner product matrices for the different bases in the Legendre family.

A naming convention is used for the first three capital letters for all matrices. The first letter refers to type of matrix.

• Mass matrices start with B

• One derivative start with C

• Stiffness - One derivative for test and trial - start with A

• Biharmonic - Two derivatives for test and trial - start with S

A matrix may consist of different types of test and trialfunctions. The next letters in the matrix name uses the ‘short_name’ method for all these different bases, see legendre.bases.py.

So a mass matrix using ShenDirichlet test and ShenNeumann trial is named BSDSNmat.

All matrices in this module may be looked up using the ‘mat’ dictionary, which takes test and trialfunctions along with the number of derivatives to be applied to each. As such the mass matrix BSDSDmat may be looked up as

>>> import numpy as np
>>> from shenfun.legendre.matrices import mat
>>> from shenfun.legendre.bases import ShenDirichlet as SD
>>> B = mat[((SD, 0), (SD, 0))]

and an instance of the matrix can be created as

>>> B0 = SD(10)
>>> BM = B((B0, 0), (B0, 0))
>>> d = {-2: np.array([-0.4, -0.28571429, -0.22222222, -0.18181818, -0.15384615, -0.13333333]),
...       0: np.array([2.4, 0.95238095, 0.62222222, 0.46753247, 0.37606838, 0.31515152, 0.27149321, 0.23859649]),
...       2: np.array([-0.4, -0.28571429, -0.22222222, -0.18181818, -0.15384615, -0.13333333])}
>>> [np.all(abs(BM[k]-v) < 1e-7) for k, v in d.items()]
[True, True, True]

However, this way of creating matrices is not reccommended use. It is far more elegant to use the TrialFunction/TestFunction interface, and to generate the matrix as an inner product:

>>> from shenfun import TrialFunction, TestFunction, inner
>>> u = TrialFunction(B0)
>>> v = TestFunction(B0)
>>> BM = inner(u, v)
>>> [np.all(abs(BM[k]-v) < 1e-7) for k, v in d.items()]
[True, True, True]

To see that this is in fact the BSDSDmat:

>>> print(BM.__class__)
<class 'shenfun.legendre.matrices.BSDSDmat'>

class shenfun.legendre.matrices.ADNDNmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Stiffness matrix \(A=(a_{kj}) \in \mathbb{R}^{M \times N}\), where

\[a_{kj} = (\phi'_j, \phi'_k)\]

where \(\phi_k \in\) legendre.bases.DirichletNeumann, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.ALLmat(test, trial, scale=1, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Stiffness matrix \(A=(a_{kj}) \in \mathbb{R}^{M \times N}\), where

\[\begin{split}a_{kj} &= (L''_j, L_k), \text{ or } \\ a_{kj} &= (L_j, L''_k)\end{split}\]

where \(L_k \in\) legendre.bases.Orthogonal, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.ASBSBmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Stiffness matrix \(A=(a_{kj}) \in \mathbb{R}^{M \times N}\), where

\[a_{kj} = (\phi'_j, \phi'_k)\]

where \(\phi_k \in\) legendre.bases.ShenBiharmonic, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.ASDSD2Trp1mat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Stiffness matrix \(A=(a_{kj}) \in \mathbb{R}^{M \times N}\), where

\[a_{kj} = (\phi''_j, (1+x)\phi_k)\]

where \(\phi_k \in\) legendre.bases.ShenDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.ASDSD2rp1mat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Stiffness matrix \(A=(a_{kj}) \in \mathbb{R}^{M \times N}\), where

\[a_{kj} = (\phi_j, (1+x)\phi''_k)\]

where \(\phi_k \in\) legendre.bases.ShenDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.ASDSDmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Stiffness matrix \(A=(a_{kj}) \in \mathbb{R}^{M \times N}\), where

\[a_{kj} = (\phi'_j, \phi'_k)\]

where \(\phi_k \in\) legendre.bases.ShenDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.ASDSDrp1mat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Stiffness matrix \(A=(a_{kj}) \in \mathbb{R}^{M \times N}\), where

\[a_{kj} = (\phi'_j, (1+x)\phi'_k)\]

where \(\phi_k \in\) legendre.bases.ShenDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.ASNSNmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Stiffness matrix \(A=(a_{kj}) \in \mathbb{R}^{M \times N}\), where

\[a_{kj} = (\phi'_j, \phi'_k)\]

where \(\phi_k \in\) legendre.bases.ShenNeumann, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.AUDUDrp1mat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Stiffness matrix \(A=(a_{kj}) \in \mathbb{R}^{M \times N}\), where

\[a_{kj} = (\phi'_j, (1+x)\phi'_k)\]

where \(\phi_k \in\) legendre.bases.UpperDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.AUDUDrp1smat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Stiffness matrix \(A=(a_{kj}) \in \mathbb{R}^{M \times N}\), where

\[a_{kj} = (\phi'_j, (1+x)^2\phi'_k)\]

where \(\phi_k \in\) legendre.bases.UpperDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.BDNDNmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj}=(\phi_j, \phi_k),\]

where \(\phi_k \in\) legendre.bases.DirichletNeumann, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.BLDLDmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj} = (\phi_j, \phi_k)\]

where \(\phi_k \in\) legendre.bases.LowerDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.BLLmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj}=(L_j, L_k),\]

where \(L_k \in\) legendre.bases.Orthogonal, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.BLSDmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj}=(\psi_j, L_k),\]

where the test function \(L_k \in\) legendre.bases.Orthogonal, the trial function \(\psi_j \in\) legendre.bases.ShenDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.BSBSBmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj}=(\phi_j, \phi_k),\]

where \(\phi_k \in\) legendre.bases.ShenBiharmonic, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

get_solver()

Return appropriate solver for self

Note

Fall back on generic Solve, which is using Scipy sparse matrices with splu/spsolve. This is still pretty fast.

class shenfun.legendre.matrices.BSDLmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj}=(L_j, \phi_k),\]

where the test function \(\phi_k \in\) legendre.bases.ShenDirichlet, the trial function \(L_j \in\) legendre.bases.Orthogonal, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.BSDSD1orp1mat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj} = (\phi_j, \frac{1}{1+x}\phi_k)\]

where \(\phi_k \in\) legendre.bases.ShenDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.BSDSDmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj}=(\phi_j, \phi_k),\]

where \(\phi_k \in\) legendre.bases.ShenDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

get_solver()

Return appropriate solver for self

Note

Fall back on generic Solve, which is using Scipy sparse matrices with splu/spsolve. This is still pretty fast.

class shenfun.legendre.matrices.BSDSDrp1mat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj} = (\phi_j, (1+x)\phi_k)\]

where \(\phi_k \in\) legendre.bases.ShenDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.BSNSNmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj}=(\phi_j, \phi_k),\]

where \(\phi_k \in\) legendre.bases.ShenNeumann, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.BUDUDmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj} = (\phi_j, \phi_k)\]

where \(\phi_k \in\) legendre.bases.UpperDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.BUDUDrp1mat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Mass matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj} = (\phi_j, (1+x)\phi_k)\]

where \(\phi_k \in\) legendre.bases.UpperDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.BUDUDrp1smat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Matrix \(B=(b_{kj}) \in \mathbb{R}^{M \times N}\), where

\[b_{kj} = (\phi_j, (1+x)^2\phi_k)\]

where \(\phi_k \in\) legendre.bases.UpperDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.CLLmat(test, trial, scale=1, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Derivative matrix \(C=(c_{kj}) \in \mathbb{R}^{M \times N}\), where

\[c_{kj}=(L'_j, L_k),\]

where \(L_k \in\) legendre.bases.Orthogonal, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

matvec(v, c, format='self', axis=0)

Matrix vector product

Returns c = dot(self, v)

Parameters:

• v (array) – Numpy input array of ndim>=1

• c (array) – Numpy output array of same shape as v

• format (str, optional) – Choice for computation

  • csr - Compressed sparse row format

  • dia - Sparse matrix with DIAgonal storage

  • python - Use numpy and vectorization

  • self - To be implemented in subclass

  • cython - Cython implementation that may be implemented in subclass

  • numba - Numba implementation that may be implemented in subclass

  Using config['matrix']['sparse']['matvec'] setting if format is None

• axis (int, optional) – The axis over which to take the matrix vector product

class shenfun.legendre.matrices.CLLmatT(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Derivative matrix \(C=(c_{kj}) \in \mathbb{R}^{M \times N}\), where

\[c_{kj}=(L_j, L'_k),\]

where \(L_k \in\) legendre.bases.Orthogonal, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.CLSDmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Derivative matrix \(C=(c_{kj}) \in \mathbb{R}^{M \times N}\), where

\[c_{kj}=(\psi'_j, L_k),\]

where the test function \(L_k \in\) legendre.bases.Orthogonal, the trial function \(\psi_j \in\) legendre.bases.ShenDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.CSDLmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Derivative matrix \(C=(c_{kj}) \in \mathbb{R}^{M \times N}\), where

\[c_{kj}=(L'_j, \phi_k),\]

where the test function \(\phi_k \in\) legendre.bases.ShenDirichlet, the trial function \(L_j \in\) legendre.bases.Orthogonal, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.CSDSDTmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Derivative matrix \(C=(c_{kj}) \in \mathbb{R}^{M \times N}\), where

\[c_{kj}=(\phi_j, \phi'_k),\]

where \(\phi_k \in\) legendre.bases.ShenDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.CSDSDmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Derivative matrix \(C=(c_{kj}) \in \mathbb{R}^{M \times N}\), where

\[c_{kj}=(\phi'_j, \phi_k),\]

where \(\phi_k \in\) legendre.bases.ShenDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.CUDUDrp1mat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Derivative matrix \(C=(c_{kj}) \in \mathbb{R}^{M \times N}\), where

\[c_{kj} = (\phi_j, (1+x)\phi'_k)\]

where \(\phi_k \in\) legendre.bases.UpperDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.GUDUDrp1smat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Stiffness matrix \(A=(a_{kj}) \in \mathbb{R}^{M \times N}\), where

\[a_{kj} = (\phi_j, (1+x)^2\phi''_k)\]

where \(\phi_k \in\) legendre.bases.UpperDirichlet, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.SBFBFmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Biharmonic matrix \(S=(s_{kj}) \in \mathbb{R}^{M \times N}\), where

\[s_{kj} = (\phi''_j, \phi''_k)\]

where \(\phi_k \in\) legendre.bases.BeamFixedFree, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

class shenfun.legendre.matrices.SSBSBmat(test, trial, scale=1.0, measure=1, assemble=None, kind=None, fixed_resolution=None)

Bases: SpectralMatrix

Biharmonic matrix \(S=(s_{kj}) \in \mathbb{R}^{M \times N}\), where

\[s_{kj} = (\phi''_j, \phi''_k)\]

where \(\phi_k \in\) legendre.bases.ShenBiharmonic, and test and trial spaces have dimensions of M and N, respectively.

assemble(method)

Return diagonals of SpectralMatrix

Parameters:

method (str) – Type of integration

• ‘exact’

• ‘quadrature’

Note

Subclass SpectralMatrix and overload this method in order to provide a fast and accurate implementation of the matrix representing an inner product. See the matrix modules in either one of

• legendre.matrix

• chebyshev.matrix

• chebyshevu.matrix

• ultraspherical.matrix

• fourier.matrix

• laguerre.matrix

• hermite.matrix

• jacobi.matrix

Example

The mass matrix for Chebyshev polynomials is

\[(T_j, T_i)_{\omega} = \frac{c_i \pi}{2}\delta_{ij},\]

where \(c_0=2\) and \(c_i=1\) for integer \(i>0\). We can implement this as

>>> from shenfun import SpectralMatrix
>>> class Bmat(SpectralMatrix):
...     def assemble(self, method):
...         test, trial = self.testfunction, self.trialfunction
...         ci = np.ones(test[0].N)
...         ci[0] = 2
...         if test[0].quad == 'GL' and method != 'exact':
...             # Gauss-Lobatto quadrature inexact at highest polynomial order
...             ci[-1] = 2
...         return {0: ci*np.pi/2}

Here {0: ci*np.pi/2} is the 0’th diagonal of the matrix. Note that test and trial are two-tuples of (instance of :class:.SpectralBase`, number)`, where the number represents the number of derivatives. For the mass matrix the number will be 0. Also note that the length of the diagonal must be correct.

Module contents

Functionality for working with Legendre bases