From a21a0b83ba7882287db376dd7b748fc5539bee5d Mon Sep 17 00:00:00 2001 From: Daniel VandenHeuvel <95613936+DanielVandH@users.noreply.github.com> Date: Wed, 3 Sep 2025 15:11:04 +0100 Subject: [PATCH 1/4] Support Jacobi(-1,0)\Ultraspherical(-1/2) --- Project.toml | 2 +- src/classical/ultraspherical.jl | 12 ++++++++++-- test/test_ultraspherical.jl | 10 ++++++++++ 3 files changed, 21 insertions(+), 3 deletions(-) diff --git a/Project.toml b/Project.toml index 3fa63d9d..0f6ebfe2 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.3" +version = "0.15.4" [deps] diff --git a/src/classical/ultraspherical.jl b/src/classical/ultraspherical.jl index de525a53..b73693f9 100644 --- a/src/classical/ultraspherical.jl +++ b/src/classical/ultraspherical.jl @@ -199,8 +199,16 @@ function \(A::Ultraspherical, B::Jacobi) Diagonal(Ã[1,:]./A[1,:]) * (Ã\B) end function \(A::Jacobi, B::Ultraspherical) - B̃ = Jacobi(B) - (A\B̃)*Diagonal(B[1,:]./B̃[1,:]) + if isone(-2B.λ) && (A.a, A.b) == (-1, 0) + T = promote_type(eltype(A), eltype(B)) + n = -2one(T) ./ (2 .* (2:∞) .- one(T)) + dv = Vcat(one(T), -2one(T), n) + ev = Vcat(-one(T), n) + LazyBandedMatrices.Bidiagonal(dv, ev, :U) + else + B̃ = Jacobi(B) + (A\B̃)*Diagonal(B[1,:]./B̃[1,:]) + end end function \(wA::Weighted{<:Any,<:Ultraspherical}, wB::Weighted{<:Any,<:Jacobi}) diff --git a/test/test_ultraspherical.jl b/test/test_ultraspherical.jl index 526eaed4..95c77627 100644 --- a/test/test_ultraspherical.jl +++ b/test/test_ultraspherical.jl @@ -207,4 +207,14 @@ using ClassicalOrthogonalPolynomials: grammatrix @testset "ChebyshevInterval constructior" begin @test ultraspherical(2,ChebyshevInterval()) ≡ Ultraspherical(2) end +end + +@testset "Jacobi(-1, 0) \ Ultraspherical(-1/2)" begin + A = Jacobi(-1, 0) + B = Ultraspherical(-1/2) + R = A \ B + AR = A * R + lhs = AR[0.3918, 1:100] + rhs = B[0.3918, 1:100] + @test lhs ≈ rhs end \ No newline at end of file From 2f91512dd5ad97861a5c51bab87e58fa898d1fe5 Mon Sep 17 00:00:00 2001 From: Daniel VandenHeuvel <95613936+DanielVandH@users.noreply.github.com> Date: Wed, 3 Sep 2025 15:13:08 +0100 Subject: [PATCH 2/4] Fix slash --- test/test_ultraspherical.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/test_ultraspherical.jl b/test/test_ultraspherical.jl index 95c77627..6976a2f0 100644 --- a/test/test_ultraspherical.jl +++ b/test/test_ultraspherical.jl @@ -209,7 +209,7 @@ using ClassicalOrthogonalPolynomials: grammatrix end end -@testset "Jacobi(-1, 0) \ Ultraspherical(-1/2)" begin +@testset "Jacobi(-1, 0) / Ultraspherical(-1/2)" begin A = Jacobi(-1, 0) B = Ultraspherical(-1/2) R = A \ B From 00fb966ccfd8e90e2fb8260bd1f725f6a38f203a Mon Sep 17 00:00:00 2001 From: Daniel VandenHeuvel <95613936+DanielVandH@users.noreply.github.com> Date: Wed, 3 Sep 2025 16:51:42 +0100 Subject: [PATCH 3/4] Add Jacobi(0,-1) --- src/classical/ultraspherical.jl | 5 +++-- test/test_ultraspherical.jl | 17 +++++++++-------- 2 files changed, 12 insertions(+), 10 deletions(-) diff --git a/src/classical/ultraspherical.jl b/src/classical/ultraspherical.jl index b73693f9..76ec69ec 100644 --- a/src/classical/ultraspherical.jl +++ b/src/classical/ultraspherical.jl @@ -199,11 +199,12 @@ function \(A::Ultraspherical, B::Jacobi) Diagonal(Ã[1,:]./A[1,:]) * (Ã\B) end function \(A::Jacobi, B::Ultraspherical) - if isone(-2B.λ) && (A.a, A.b) == (-1, 0) + if B == Ultraspherical(-1/2) && (A == Jacobi(-1, 0) || A == Jacobi(0, -1)) T = promote_type(eltype(A), eltype(B)) n = -2one(T) ./ (2 .* (2:∞) .- one(T)) + sgn = A == Jacobi(-1, 0) ? one(T) : -one(T) dv = Vcat(one(T), -2one(T), n) - ev = Vcat(-one(T), n) + ev = Vcat(-sgn, sgn .* n) LazyBandedMatrices.Bidiagonal(dv, ev, :U) else B̃ = Jacobi(B) diff --git a/test/test_ultraspherical.jl b/test/test_ultraspherical.jl index 6976a2f0..89ca9617 100644 --- a/test/test_ultraspherical.jl +++ b/test/test_ultraspherical.jl @@ -209,12 +209,13 @@ using ClassicalOrthogonalPolynomials: grammatrix end end -@testset "Jacobi(-1, 0) / Ultraspherical(-1/2)" begin - A = Jacobi(-1, 0) - B = Ultraspherical(-1/2) - R = A \ B - AR = A * R - lhs = AR[0.3918, 1:100] - rhs = B[0.3918, 1:100] - @test lhs ≈ rhs +@testset "Jacobi(-1, 0) \\ Ultraspherical(-1/2) and Jacobi(0, -1) \\ Ultraspherical(-1/2)" begin + for A in (Jacobi(-1,0),Jacobi(0,-1)) + B = Ultraspherical(-1/2) + R = A \ B + AR = A * R + lhs = AR[0.3918, 1:100] + rhs = B[0.3918, 1:100] + @test lhs ≈ rhs + end end \ No newline at end of file From cbc74241dedc64d792ca4c4372e42b1e9d13bf85 Mon Sep 17 00:00:00 2001 From: Daniel VandenHeuvel <95613936+DanielVandH@users.noreply.github.com> Date: Wed, 3 Sep 2025 20:10:40 +0100 Subject: [PATCH 4/4] Clarify --- src/classical/ultraspherical.jl | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/classical/ultraspherical.jl b/src/classical/ultraspherical.jl index 76ec69ec..a2a857e4 100644 --- a/src/classical/ultraspherical.jl +++ b/src/classical/ultraspherical.jl @@ -200,6 +200,8 @@ function \(A::Ultraspherical, B::Jacobi) end function \(A::Jacobi, B::Ultraspherical) if B == Ultraspherical(-1/2) && (A == Jacobi(-1, 0) || A == Jacobi(0, -1)) + # In this case, Jacobi(-1, -1) is (currently) undefined, so the conversion via B̃ = Jacobi(B) leads to NaNs + # from evaluating in B̃[1, :] T = promote_type(eltype(A), eltype(B)) n = -2one(T) ./ (2 .* (2:∞) .- one(T)) sgn = A == Jacobi(-1, 0) ? one(T) : -one(T)