From 9e9f05e855dbb6f32381096e4e6cc71fdeea4ba3 Mon Sep 17 00:00:00 2001 From: Sheehan Olver Date: Wed, 5 Feb 2025 14:29:20 +0000 Subject: [PATCH] Add normalized OP functions --- Project.toml | 2 +- src/ClassicalOrthogonalPolynomials.jl | 3 ++- src/classical/hermite.jl | 8 ++++++++ src/classical/laguerre.jl | 16 ++++++++++++++++ src/classical/legendre.jl | 9 +++++++++ src/classical/ultraspherical.jl | 1 + 6 files changed, 37 insertions(+), 2 deletions(-) diff --git a/Project.toml b/Project.toml index 4a65c44..e71be4d 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "ClassicalOrthogonalPolynomials" uuid = "b30e2e7b-c4ee-47da-9d5f-2c5c27239acd" authors = ["Sheehan Olver "] -version = "0.15" +version = "0.15.1" [deps] diff --git a/src/ClassicalOrthogonalPolynomials.jl b/src/ClassicalOrthogonalPolynomials.jl index f0b688b..a368e3f 100644 --- a/src/ClassicalOrthogonalPolynomials.jl +++ b/src/ClassicalOrthogonalPolynomials.jl @@ -54,7 +54,8 @@ export OrthogonalPolynomial, Normalized, LanczosPolynomial, WeightedUltraspherical, WeightedChebyshev, WeightedChebyshevT, WeightedChebyshevU, WeightedJacobi, ∞, Derivative, .., Inclusion, chebyshevt, chebyshevu, legendre, jacobi, ultraspherical, - legendrep, jacobip, ultrasphericalc, laguerrel,hermiteh, normalizedjacobip, + legendrep, jacobip, ultrasphericalc, laguerrel,hermiteh, + normalizedlegendrep, normalizedjacobip, normalizedultrasphericalc, normalizedlaguerrel, normalizedhermiteh, jacobimatrix, jacobiweight, legendreweight, chebyshevtweight, chebyshevuweight, Weighted, PiecewiseInterlace, plan_transform, expand, transform diff --git a/src/classical/hermite.jl b/src/classical/hermite.jl index ecb5c7b..18ec489 100644 --- a/src/classical/hermite.jl +++ b/src/classical/hermite.jl @@ -48,6 +48,14 @@ respec to `exp(-x^2)`, at `z`. """ hermiteh(n::Integer, z) = Base.unsafe_getindex(Hermite{polynomialtype(typeof(z))}(), z, n+1) +""" + normalizedhermiteh(n, z) + +computes the normalized `n`-th Hermite polynomial, orthogonal with +respec to `exp(-x^2)`, at `z`. +""" +normalizedhermiteh(n::Integer, z) = Base.unsafe_getindex(Normalized(Hermite{polynomialtype(typeof(z))}()), z, n+1) + broadcasted(::LazyQuasiArrayStyle{2}, ::typeof(*), ::HermiteWeight{T}, ::Hermite{V}) where {T,V} = Weighted(Hermite{promote_type(T,V)}()) # H_{n+1} = 2x H_n - 2n H_{n-1} diff --git a/src/classical/laguerre.jl b/src/classical/laguerre.jl index b6ac179..42d1db5 100644 --- a/src/classical/laguerre.jl +++ b/src/classical/laguerre.jl @@ -55,6 +55,22 @@ respec to `exp(-x)`, at `z`. """ laguerrel(n::Integer, z::Number) = laguerrel(n, 0, z) +""" + normalizedlaguerrel(n, α, z) + +computes the normalized `n`-th generalized Laguerre polynomial, orthogonal with +respec to `x^α * exp(-x)`, at `z`. +""" +normalizedlaguerrel(n::Integer, α, z::Number) = Base.unsafe_getindex(Normalized(Laguerre{polynomialtype(typeof(α), typeof(z))}(α)), z, n+1) + +""" + normalizedlaguerrel(n, z) + +computes the normalized `n`-th Laguerre polynomial, orthogonal with +respec to `exp(-x)`, at `z`. +""" +normalizedlaguerrel(n::Integer, z::Number) = normalizedlaguerrel(n, 0, z) + # L_{n+1} = (-1/(n+1) x + (2n+α+1)/(n+1)) L_n - (n+α)/(n+1) L_{n-1} # - (n+α) L_{n-1} + (2n+α+1)* L_n -(n+1) L_{n+1} = x L_n diff --git a/src/classical/legendre.jl b/src/classical/legendre.jl index a11267d..cc1b165 100644 --- a/src/classical/legendre.jl +++ b/src/classical/legendre.jl @@ -167,6 +167,15 @@ computes the `n`-th Legendre polynomial at `z`. """ legendrep(n::Integer, z) = Base.unsafe_getindex(Legendre{typeof(z)}(), z, n+1) +""" + normalizedlegendrep(n, z) + +computes the normalized `n`-th Legendre polynomial at `z`. +""" +normalizedlegendrep(n::Integer, z) = Base.unsafe_getindex(Normalized(Legendre{typeof(z)}()), z, n+1) + + + show(io::IO, w::Legendre{Float64}) = summary(io, w) summary(io::IO, ::Legendre{Float64}) = print(io, "Legendre()") diff --git a/src/classical/ultraspherical.jl b/src/classical/ultraspherical.jl index 7b75f0e..1bafd3c 100644 --- a/src/classical/ultraspherical.jl +++ b/src/classical/ultraspherical.jl @@ -54,6 +54,7 @@ const WeightedUltraspherical{T} = WeightedBasis{T,<:UltrasphericalWeight,<:Ultra orthogonalityweight(C::Ultraspherical) = UltrasphericalWeight(C.λ) +normalizedultrasphericalc(n::Integer, λ, z) = Base.unsafe_getindex(Normalized(Ultraspherical{polynomialtype(typeof(λ),typeof(z))}(λ)), z, n+1) ultrasphericalc(n::Integer, λ, z) = Base.unsafe_getindex(Ultraspherical{polynomialtype(typeof(λ),typeof(z))}(λ), z, n+1) ultraspherical(λ, d::AbstractInterval{T}) where T = Ultraspherical{float(promote_type(eltype(λ),T))}(λ)[affine(d,ChebyshevInterval{T}()), :]