Fix constructing from b=-1 (#117)

* Support construction from b=-1

* Fix some broken tests

* Inefficient

* Add comment explaining one(Q.b), give testset a title

* and this one

---------

Co-authored-by: Sheehan Olver <[email protected]>

Author: DanielVandH

Project.toml

@@ -1,7 +1,7 @@
 name = "SemiclassicalOrthogonalPolynomials"
 uuid = "291c01f3-23f6-4eb6-aeb0-063a639b53f2"
 authors = ["Sheehan Olver <[email protected]>"]
-version = "0.5.7"
+version = "0.5.8"
 
 [deps]
 ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"

src/SemiclassicalOrthogonalPolynomials.jl

@@ -180,7 +180,10 @@ function semiclassical_jacobimatrix(Q::SemiclassicalJacobi, a, b, c)
         return Q.X
     elseif b == -one(eltype(Q.t))
         return semiclassical_jacobimatrix(Q.t, a, b, c)
-    end 
+    elseif Q.b == -one(eltype(Q.t))
+        newQ = SemiclassicalJacobi(Q.t, Q.a, one(Q.b), Q.c, Q) # Constructing first b = 1 from b = -1 is efficient since b = -1 is defined directly in terms of b = 1.
+        return semiclassical_jacobimatrix(newQ, a, b, c)
+    end
 
     if isone(Δa/2) && iszero(Δb) && iszero(Δc)  # raising by 2
         qr_jacobimatrix(Q.X,Q)

test/runtests.jl

@@ -610,7 +610,7 @@ end
     end
 end
 
-@testset "Issue #103" begin 
+@testset "Issue #103: Missing show method for SemiclassicalJacobi" begin 
     w = SemiclassicalJacobiWeight(2.0, 2.3, 5.3, 0.4)
     P = SemiclassicalJacobi(2.0, 2.3, 5.3, 0.4)
     W = Weighted(P)

test/test_neg1b.jl

@@ -234,6 +234,18 @@ end
     @test_throws ArgumentError expand(HalfWeighted{:b}(P), g)
 end
 
+@testset "Issue #115: Constructing from b = -1" begin
+    P = SemiclassicalJacobi(2.0, -1/2, -1.0, -1/2)
+    Q1 = SemiclassicalJacobi(2.0, -1/2, 0.0, -1/2, P)
+    R1 = SemiclassicalJacobi(2.0, -1/2, 0.0, -1/2)
+    Q2 = SemiclassicalJacobi(2.0, 3/2, 2.0, 3/2, P)
+    R2 = SemiclassicalJacobi(2.0, 3/2, 2.0, 3/2)
+    Q3 = SemiclassicalJacobi(2.0, 5/2, 3.0, 0.0, P)
+    R3 = SemiclassicalJacobi(2.0, 5/2, 3.0, 0.0)
+    @test Q1.X[1:100, 1:100] ≈ R1.X[1:100, 1:100]
+    @test Q2.X[1:100, 1:100] ≈ R2.X[1:100, 1:100]
+    @test Q3.X[1:100, 1:100] ≈ R3.X[1:100, 1:100]
+end
 @testset "Weighted conversion between b=-1" begin
     for (t, a, b, c) in ((2.0, 1 / 2, -1.0, 1 / 2), (2.5, 3 / 2, -1.0, 1 / 2), (2.5, 1.0, -1.0, 2.0))
         Q = SemiclassicalJacobi(t, a, b, c)