PolynomialBases.ClenshawCurtis
— TypeClenshawCurtis{T}
The nodal basis corresponding to the Clenshaw Curtis quadrature in [-1,1] with scalar type T
.
PolynomialBases.ClenshawCurtis
— TypeClenshawCurtis(p::Int, T=Float64)
Generate the ClenshawCurtis
basis of degree p
with scalar type T
.
PolynomialBases.ClosedNewtonCotes
— TypeClosedNewtonCotes{T}
The nodal basis corresponding to the closed Newton Cotes quadrature in [-1,1] with scalar type T
.
PolynomialBases.ClosedNewtonCotes
— TypeClosedNewtonCotes(p::Int, T=Float64)
Generate the ClosedNewtonCotes
basis of degree p
with scalar type T
.
PolynomialBases.GaussJacobi
— TypeGaussJacobi{T<:Real}
The nodal basis corresponding to Jacobi Gauss quadrature in [-1,1] with parameters α
, β
and scalar type T
.
PolynomialBases.GaussJacobi
— TypeGaussJacobi(p::Int, α, β, T=Float64)
Generate the JacobiLegendre
basis of degree p
with parameters α
, β
and scalar type T
.
PolynomialBases.GaussLegendre
— TypeGaussLegendre{T}
The nodal basis corresponding to Legendre Gauss quadrature in [-1,1] with scalar type T
.
PolynomialBases.GaussLegendre
— TypeGaussLegendre(p::Int, T=Float64)
Generate the GaussLegendre
basis of degree p
with scalar type T
.
PolynomialBases.LobattoLegendre
— TypeLobattoLegendre{T}
The nodal basis corresponding to Legendre Gauss Lobatto quadrature in [-1,1] with scalar type T
.
PolynomialBases.LobattoLegendre
— TypeLobattoLegendre(p::Int, T=Float64)
Generate the LobattoLegendre
basis of degree p
with scalar type T
.
PolynomialBases.barycentric_weights
— Methodbarycentric_weights{T<:Real}(x::AbstractVector{T})
Compute the barycentric weights corresponding to the nodes x
. [Kopriva, Implementing Spectral Methods for PDEs, Algorithm 30].
PolynomialBases.change_basis!
— Methodchange_basis!(ret, dest_basis::NodalBasis{Domain},
values, src_basis::NodalBasis{Domain}) where {Domain<:AbstractDomain}
Perform the change of basis for the coefficients values
from src_basis
to dest_basis
and store the resulting coefficients in ret
.
PolynomialBases.change_basis
— Methodfunction change_basis(dest_basis::NodalBasis{Domain},
values, src_basis::NodalBasis{Domain}) where {Domain<:AbstractDomain}
Perform the change of basis for the coefficients values
from src_basis
to dest_basis
.
PolynomialBases.compute_coefficients!
— Methodcompute_coefficients!(uval::AbstractVector, u, basis::NodalBasis{Line})
Compute the nodal values of the function u
at the nodes corresponding to the nodal basis basis
and store the result in uval
.
PolynomialBases.compute_coefficients
— Methodcompute_coefficients(u, basis::NodalBasis{Line})
Compute the nodal values of the function u
at the nodes corresponding to the nodal basis basis
.
PolynomialBases.degree
— Methoddegree(basis::NodalBasis{Line})
Return the polynomial degree used by basis
.
PolynomialBases.derivative_at
— Methodderivative_at(x::Real, values, nodes, baryweights)
Compute the derivative of the function represented by values
on the nodes
at x
using the corresponding barycentric weights baryweights
. [Kopriva, Implementing Spectral Methods for PDEs, Algorithm 36].
PolynomialBases.derivative_matrix
— Methodderivative_matrix(nodes, baryweights)
Compute the derivative matrix corresponding to nodes
and baryweights
. [Kopriva, Implementing Spectral Methods for PDEs, Algorithm 37].
PolynomialBases.evaluate_coefficients
— Functionevaluate_coefficients(u, basis::NodalBasis{Line}, npoints=2*length(basis.nodes))
Evaluate the coefficients u
of the nodal basis basis
at npoints
equally spaced nodes. Returns xplot, uplot
, where xplot
contains the equally spaced nodes and uplot
the corresponding values of u
.
PolynomialBases.evaluate_coefficients!
— Methodevaluate_coefficients!(xplot, uplot, u, basis::NodalBasis{Line})
Evaluate the coefficients u
of the nodal basis basis
at npoints
equally spaced nodes and store the result in xplot, uplot
. Returns xplot, uplot
, where xplot
contains the equally spaced nodes and uplot
the corresponding values of u
.
PolynomialBases.gauss_jacobi_nodes_and_weights
— Functiongauss_jacobi_nodes_and_weights(p, α, β, T=Float64::Type, tol=4*eps(T), maxit=100)
Compute the Gauss-Jacobi nodes and weights for polynomials of degree p
with parameters α
, β
using the scalar type T
, tolerance tol
and maximal number of Newton iterations maxit
[Karniadakis and Sherwin, Spectral/hp Element Methods for CFD, Appendix B].
PolynomialBases.gauss_legendre_nodes_and_weights
— Functiongauss_legendre_nodes_and_weights(p, T=Float64::Type, tol=4*eps(T), maxit=100)
Compute the Gauss-Legendre nodes and weights for polynomials of degree p
using the scalar type T
, tolerance tol
and maximal number of Newton iterations maxit
[Kopriva, Implementing Spectral Methods for PDEs, Algorithm 23].
PolynomialBases.gegenbauer
— Methodgegenbauer(x, p::Integer)
Evaluate the Gegenbauer polynomial with parameter α
of degree p
at x
using the three term recursion.
PolynomialBases.hahn
— Methodhahn(x, p::Integer, α, β, N::Integer)
Evaluate the Hahn polynomial with parameters α
, β
, N
of degree p
at x
using the three term recursion [Öffner, Zweidimensionale klassische und diskrete orthogonale Polynome, Chapter 5].
PolynomialBases.hermite
— Methodhermite(x, p::Integer)
Evaluate the Hermite polynomial of degree p
at x
using the three term recursion.
PolynomialBases.integrate
— Methodintegrate(func, u, weights)
Map the function func
to the coefficients u
and integrate with respect to the quadrature rule given by weights
.
PolynomialBases.interpolate
— Methodinterpolate(x::Real, nodes, values, baryweights)
Interpolate the function represented by values
on the nodes
using the corresponding barycentric weights baryweights
. [Kopriva, Implementing Spectral Methods for PDEs, Algorithm 31].
PolynomialBases.interpolation_matrix
— Methodinterpolation_matrix(destination, source, baryweights)
Compute the matrix performing interpolation from src
to dest
, where baryweights
are the barycentric weights corresponding to src
. [Kopriva, Implementing Spectral Methods for PDEs, Algorithm 32].
PolynomialBases.jacobi
— Methodjacobi(x, p::Integer, α, β)
Evaluate the Legendre polynomial with parameters α
, β
of degree p
at x
using the three term recursion [Karniadakis and Sherwin, Spectral/hp Element Methods for CFD, Appendix A].
PolynomialBases.jacobi_and_derivative
— Methodjacobi_and_derivative(x, p::Integer, α, β)
Evaluate the Jacobi polynomial with parameters α
, β
of degree p
and its derivative at x
using the three term recursion [Karniadakis and Sherwin, Spectral/hp Element Methods for CFD, Appendix A].
PolynomialBases.jacobi_vandermonde
— Methodjacobi_vandermonde(nodes, α, β)
Computes the Vandermonde matrix with respect to the Jacobi polynomials with parameters α
, β
and the nodal basis on nodes
. The Vandermonde matrix V
is the transformation matrix from the modal Jacobi basis to the nodal Lagrange basis associated with nodes
.
PolynomialBases.laguerre
— Methodlaguerre(x, p::Integer, α)
Evaluate the generalised Laguerre polynomial with parameter α
of degree p
at x
using the three term recursion.
PolynomialBases.laguerre
— Methodlaguerre(x, p::Integer)
Evaluate the Laguerre polynomial of degree p
at x
using the three term recursion.
PolynomialBases.legendre
— Methodlegendre(x, p::Integer)
Evaluate the Legendre polynomial of degree p
at x
using the three term recursion [Kopriva, Implementing Spectral Methods for PDEs, Algorithm 20].
PolynomialBases.legendre_D
— Functionlegendre_D(p, T=Float64)
Computes the derivative matrix in the modal Legendre basis up to degree p
using the scalar type T
.
PolynomialBases.legendre_M
— Functionlegendre_M(p, T=Float64)
Computes the diagonal mass matrix in the modal Legendre basis up to degree p
using the scalar type T
.
PolynomialBases.legendre_and_derivative
— Methodlegendre_and_derivative(x, p::Integer)
Evaluate the Legendre polynomial of degree p
and its derivative at x
using the three term recursion [Kopriva, Implementing Spectral Methods for PDEs, Algorithm 22].
PolynomialBases.legendre_vandermonde
— Methodlegendre_vandermonde(nodes)
Computes the Vandermonde matrix with respect to the Legendre polynomials and the nodal basis on nodes
. The Vandermonde matrix V
is the transformation matrix from the modal Legendre basis to the nodal Lagrange basis associated with nodes
.
PolynomialBases.lobatto_legendre_nodes_and_weights
— Functionlobatto_legendre_nodes_and_weights(p, T=Float64::Type, tol=4*eps(T), maxit=100)
Compute the Lobatto-Legendre nodes and weights for polynomials of degree p
using the scalar type T
, tolerance tol
and maximal number of Newton iterations maxit
[Kopriva, Implementing Spectral Methods for PDEs, Algorithm 25].
PolynomialBases.map_from_canonical!
— Methodmap_from_canonical!(x, ξ, xmin, xmax, basis::NodalBasis{Line})
Map the given canonical coordinate ξ
in the interval [-1, 1] to the corresponding coordinate x
in the interval [xmin
, xmax
], updating x
.
PolynomialBases.map_from_canonical
— Methodmap_from_canonical(ξ, xmin, xmax, basis::NodalBasis{Line})
Map the given canonical coordinate ξ
in the interval [-1, 1] to the corresponding coordinate x
in the interval [xmin
, xmax
].
PolynomialBases.map_to_canonical!
— Methodmap_to_canonical!(ξ, x, xmin, xmax, basis::NodalBasis{Line})
Map the given coordinate x
in the interval [xmin
, xmax
] to the corresponding canonical coordinate ξ
in the interval [-1, 1], updating ξ
.
PolynomialBases.map_to_canonical
— Methodmap_to_canonical(x, xmin, xmax, basis::NodalBasis{Line})
Map the given coordinate x
in the interval [xmin
, xmax
] to the corresponding canonical coordinate ξ
in the interval [-1, 1].
PolynomialBases.utility_matrices
— Methodutility_matrices(basis::NodalBasis{Line})
Return the matrices
D
, derivative matrixM
, mass matrixR
, restriction matrix (interpolation to the boundaries)B
, boundary normal matrixMinvRtB = M \ R' * B
used in the formulation of flux reconstruction / correction procedure via reconstruction using summation-by-parts operators, cf. Ranocha, Öffner, Sonar (2016) Summation-by-parts operators for correction procedure via reconstruction, Journal of Computational Physics 311, 299-328.