Basic Grids
CompositeGrids.SimpleGrid
— ModuleBasic grids including common grids like arbitrary grids, uniform grids, log grids, and optimized grids like barycheb for interpolation and gausslegendre for integration.
CompositeGrids.SimpleGrid.AbstractGrid
— TypeAll Grids are derived from AbstractGrid; ClosedGrid has bound[1], bound[2] == grid[1], grid[end], while OpenGrid has bound[1]<grid[1]<grid[end]<bound[2]
CompositeGrids.SimpleGrid.Arbitrary
— Typestruct Arbitrary{T<:AbstractFloat} <: ClosedGrid
Arbitrary grid generated from given sorted grid.
#Members:
bound
: boundary of the gridsize
: number of grid pointsgrid
: grid pointsweight
: integration weight
CompositeGrids.SimpleGrid.BaryCheb
— Typestruct BaryCheb{T<:AbstractFloat} <: OpenGrid
BaryCheb grid generated on [bound[1], bound[2]] with order N.
#Members:
bound
: boundary of the gridsize
: number of grid pointsgrid
: grid pointsweight
: interpolation weight
CompositeGrids.SimpleGrid.GaussLegendre
— Typestruct GaussLegendre{T<:AbstractFloat} <: OpenGrid
GaussLegendre grid generated on [bound[1], bound[2]] with order N.
#Members:
bound
: boundary of the gridsize
: number of grid pointsgrid
: grid pointsweight
: integration weight
CompositeGrids.SimpleGrid.Log
— Typestruct Log{T<:AbstractFloat} <: ClosedGrid
Log grid generated on [bound[1], bound[2]] with N grid points.
Minimal interval is set to be minterval. Dense to sparse if d2s, vice versa.
On [0, 1], a typical d2s Log grid looks like [0, λ^(N-1), ..., λ^2, λ, 1].
#Members:
bound
: boundary of the gridsize
: number of grid pointsgrid
: grid pointsweight
: integration weightλ
: scale parameterd2s
: dense to sparse or not
CompositeGrids.SimpleGrid.Uniform
— Typestruct Uniform{T<:AbstractFloat} <: ClosedGrid
Uniform grid generated on [bound[1], bound[2]] with N points
#Members:
bound
: boundary of the gridsize
: number of grid pointsgrid
: grid pointsweight
: integration weight
Base.floor
— Methodfunction Base.floor(grid::AbstractGrid, x) #where {T}
use basic searchsorted function to find the index of largest
grid point smaller than x.
return 1 for x<grid[1] and grid.size-1 for x>grid[end].
Base.floor
— Methodfunction Base.floor(grid::Log{T}, x) where {T}
find the index of largest
grid point smaller than x.
return 1 for x<grid[1] and grid.size-1 for x>grid[end].
Base.floor
— Methodfunction Base.floor(grid::Uniform{T}, x) where {T}
find the index of largest
grid point smaller than x.
return 1 for x<grid[1] and grid.size-1 for x>grid[end].
CompositeGrids.SimpleGrid.barycheb
— Methodfunction barycheb(n, x, f, wc, xc)
Barycentric Lagrange interpolation at Chebyshev nodes Reference: Berrut, J.P. and Trefethen, L.N., 2004. Barycentric lagrange interpolation. SIAM review, 46(3), pp.501-517.
Arguments
n
: order of the Chebyshev interpolationx
: coordinate to interpolatef
: array of size n, function at the Chebyshev nodeswc
: array of size n, Barycentric Lagrange interpolation weightsxc
: array of size n, coordinates of Chebyshev nodes
Returns
- Interpolation result
CompositeGrids.SimpleGrid.barychebinit
— Methodbarychebinit(n)
Get Chebyshev nodes of first kind and corresponding barycentric Lagrange interpolation weights. Reference: Berrut, J.P. and Trefethen, L.N., 2004. Barycentric lagrange interpolation. SIAM review, 46(3), pp.501-517.
Arguments
n
: order of the Chebyshev interpolation
Returns
- Chebyshev nodes
- Barycentric Lagrange interpolation weights