Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Inconsistent behavior for t<1 #37

Open
TSGut opened this issue Mar 29, 2021 · 1 comment
Open

Inconsistent behavior for t<1 #37

TSGut opened this issue Mar 29, 2021 · 1 comment

Comments

@TSGut
Copy link
Member

TSGut commented Mar 29, 2021

I don't think t<1 is meant to work to begin with but we should ensure the package behavior with such input is still consistent.

Right now this is definitely not the case. For example, this throws a Domain Error:

julia> testcase1(x) = SemiclassicalJacobi(0.51,0,1,1)[x,2]
testcase1 (generic function with 1 method)

julia> testcase1(0.4)
ERROR: LoadError: DomainError with -0.08526165895336571:
sqrt will only return a complex result if called with a complex argument. Try sqrt(Complex(x)).
Stacktrace:
 [1] throw_complex_domainerror(::Symbol, ::Float64) at ./math.jl:33
 [2] sqrt(::Float64) at ./math.jl:573
 [3] lanczos!(::UnitRange{Int64}, ::LazyBandedMatrices.SymTridiagonal{Float64,ApplyArray{Float64,1,typeof(vcat),Tuple{Float64,BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(-),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{Float64,BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Fill{Float64,1,Tuple{InfiniteArrays.OneToInf{Int64}}}}},Float64}}}},BroadcastArray{Float64,1,typeof(sqrt),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}},ApplyArray{Float64,1,typeof(vcat),Tuple{Float64,BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Int64,Int64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}}}}}}}}}, ::ClassicalOrthogonalPolynomials.Clenshaw{Float64,SubArray{Float64,1,Array{Float64,1},Tuple{Base.OneTo{Int64}},true},BroadcastArray{Float64,1,typeof(inv),Tuple{BroadcastArray{Float64,1,typeof(sqrt),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}},ApplyArray{Float64,1,typeof(vcat),Tuple{Float64,BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Int64,Int64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}}}}}}}}}},BroadcastArray{Float64,1,typeof(-),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{ApplyArray{Float64,1,typeof(vcat),Tuple{Float64,BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(-),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{Float64,BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Fill{Float64,1,Tuple{InfiniteArrays.OneToInf{Int64}}}}},Float64}}}},BroadcastArray{Float64,1,typeof(sqrt),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}},ApplyArray{Float64,1,typeof(vcat),Tuple{Float64,BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Int64,Int64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}}}}}}}}}}}},BroadcastArray{Float64,1,typeof(/),Tuple{ApplyArray{Float64,1,typeof(vcat),Tuple{Float64,BroadcastArray{Float64,1,typeof(sqrt),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}},ApplyArray{Float64,1,typeof(vcat),Tuple{Float64,BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Int64,Int64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}}}}}}}}}},BroadcastArray{Float64,1,typeof(sqrt),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}},ApplyArray{Float64,1,typeof(vcat),Tuple{Float64,BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Int64,Int64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}}}}}}}}}},LazyBandedMatrices.SymTridiagonal{Float64,ApplyArray{Float64,1,typeof(vcat),Tuple{Float64,BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(-),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{Float64,BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Fill{Float64,1,Tuple{InfiniteArrays.OneToInf{Int64}}}}},Float64}}}},BroadcastArray{Float64,1,typeof(sqrt),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}},ApplyArray{Float64,1,typeof(vcat),Tuple{Float64,BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Int64,Int64},InfiniteArrays.InfStepRange{Float64,Float64}}},BroadcastArray{Float64,1,typeof(*),Tuple{InfiniteArrays.InfStepRange{Float64,Float64},InfiniteArrays.InfStepRange{Float64,Float64}}}}},Float64}}}}}}}}}}, ::LazyArrays.CachedArray{Float64,1,Array{Float64,1},Zeros{Float64,1,Tuple{InfiniteArrays.OneToInf{Int64}}}}, ::LazyArrays.CachedArray{Float64,1,Array{Float64,1},Zeros{Float64,1,Tuple{InfiniteArrays.OneToInf{Int64}}}}, ::UpperTriangular{Float64,LazyArrays.CachedArray{Float64,2,Array{Float64,2},Zeros{Float64,2,Tuple{InfiniteArrays.OneToInf{Int64},InfiniteArrays.OneToInf{Int64}}}}}) at /home/timon/.julia/packages/ClassicalOrthogonalPolynomials/vc0Uq/src/lanczos.jl:19
 [4] resizedata!(::ClassicalOrthogonalPolynomials.LanczosData{Float64}, ::Int64) at /home/timon/.julia/packages/ClassicalOrthogonalPolynomials/vc0Uq/src/lanczos.jl:61
 [5] resizedata! at /home/timon/.julia/packages/ClassicalOrthogonalPolynomials/vc0Uq/src/lanczos.jl:132 [inlined]
 [6] _lanczos_getindex at /home/timon/.julia/packages/ClassicalOrthogonalPolynomials/vc0Uq/src/lanczos.jl:137 [inlined]
 [7] getindex at /home/timon/.julia/packages/ClassicalOrthogonalPolynomials/vc0Uq/src/lanczos.jl:145 [inlined]
 [8] _getindex at ./abstractarray.jl:1100 [inlined]
 [9] getindex at ./abstractarray.jl:1060 [inlined]
 [10] _broadcast_getindex at ./broadcast.jl:601 [inlined]
 [11] _getindex at ./broadcast.jl:645 [inlined]
 [12] _broadcast_getindex at ./broadcast.jl:620 [inlined]
 [13] getindex at ./broadcast.jl:575 [inlined]
 [14] getindex at /home/timon/.julia/packages/LazyArrays/lAt5J/src/lazybroadcasting.jl:88 [inlined]
 [15] forwardrecurrence!(::Array{Float64,1}, ::BroadcastArray{Float64,1,typeof(inv),Tuple{ClassicalOrthogonalPolynomials.LanczosJacobiBand{Float64}}}, ::BroadcastArray{Float64,1,typeof(-),Tuple{BroadcastArray{Float64,1,typeof(/),Tuple{ClassicalOrthogonalPolynomials.LanczosJacobiBand{Float64},ClassicalOrthogonalPolynomials.LanczosJacobiBand{Float64}}}}}, ::BroadcastArray{Float64,1,typeof(/),Tuple{ApplyArray{Float64,1,typeof(vcat),Tuple{Float64,ClassicalOrthogonalPolynomials.LanczosJacobiBand{Float64}}},ClassicalOrthogonalPolynomials.LanczosJacobiBand{Float64}}}, ::Float64, ::Float64) at /home/timon/.julia/packages/FastTransforms/DH8k1/src/clenshaw.jl:13
 [16] copyto!(::Array{Float64,1}, ::SubArray{Float64,1,SemiclassicalJacobi{Float64},Tuple{Float64,Base.OneTo{Int64}},false}) at /home/timon/.julia/packages/ClassicalOrthogonalPolynomials/vc0Uq/src/clenshaw.jl:38
 [17] getindex at /home/timon/.julia/packages/ClassicalOrthogonalPolynomials/vc0Uq/src/clenshaw.jl:26 [inlined]
 [18] getindex(::SemiclassicalJacobi{Float64}, ::Float64, ::Int64) at /home/timon/.julia/packages/ClassicalOrthogonalPolynomials/vc0Uq/src/clenshaw.jl:70
 [19] testcase1(::Float64) at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/test/runtests.jl:435
 [20] top-level scope at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/test/runtests.jl:437
in expression starting at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/test/runtests.jl:437

But this actually gives us a value:

julia> testcase2(x) = SemiclassicalJacobi(0.91,0,1,1)[x,2]
testcase2 (generic function with 1 method)

julia> testcase2(0.4)
0.8915394431603408

But of course it's not like expansions actually work:

julia> P = SemiclassicalJacobi(0.91,0,1,1)
SemiclassicalJacobi with weight x^0.0 * (1-x)^1.0 * (0.91-x)^1.0

julia> x = axes(P,1)
Inclusion(0.0..1.0 (Unit))

julia> P \ x
ERROR: LoadError: DomainError with -1.0989010989010988:
log1p will only return a complex result if called with a complex argument. Try log1p(Complex(x)).
Stacktrace:
 [1] throw_complex_domainerror(::Symbol, ::Float64) at ./math.jl:33
 [2] log1p(::Float64) at ./special/log.jl:356
 [3] expm1nlog1p at /home/timon/.julia/packages/HypergeometricFunctions/g0fPc/src/specialfunctions.jl:99 [inlined]
 [4] _₂F₁(::Float64, ::Float64, ::Float64, ::Float64) at /home/timon/.julia/packages/HypergeometricFunctions/g0fPc/src/gauss.jl:39
 [5] _₂F₁(::Float64, ::Float64, ::Float64, ::Float64) at /home/timon/.julia/packages/HypergeometricFunctions/g0fPc/src/gauss.jl:8
 [6] sum(::SemiclassicalJacobiWeight{Float64}) at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:49
 [7] massmatrix at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:343 [inlined]
 [8] copy(::Ldiv{SemiclassicalOrthogonalPolynomials.SemiclassicalJacobiLayout,SemiclassicalOrthogonalPolynomials.SemiclassicalJacobiLayout,SemiclassicalJacobi{Float64},SemiclassicalJacobi{Float64}}) at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:300
 [9] materialize(::Ldiv{SemiclassicalOrthogonalPolynomials.SemiclassicalJacobiLayout,SemiclassicalOrthogonalPolynomials.SemiclassicalJacobiLayout,SemiclassicalJacobi{Float64},SemiclassicalJacobi{Float64}}) at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:22
 [10] ldiv(::SemiclassicalJacobi{Float64}, ::SemiclassicalJacobi{Float64}) at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:86
 [11] \ at /home/timon/.julia/packages/QuasiArrays/BZtPi/src/matmul.jl:34 [inlined]
 [12] semijacobi_ldiv(::SemiclassicalJacobi{Float64}, ::SubQuasiArray{Float64,2,Jacobi{Float64},Tuple{ContinuumArrays.AffineMap{Float64,Inclusion{Float64,DomainSets.UnitInterval{Float64}},Inclusion{Float64,ChebyshevInterval{Float64}}},Base.Slice{InfiniteArrays.OneToInf{Int64}}},false}) at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:276
 [13] copy at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:290 [inlined]
 [14] materialize at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:22 [inlined]
 [15] ldiv at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:86 [inlined]
 [16] \ at /home/timon/.julia/packages/QuasiArrays/BZtPi/src/matmul.jl:34 [inlined]
 [17] _copy_ldiv_mul at /home/timon/.julia/packages/LazyArrays/lAt5J/src/linalg/inv.jl:98 [inlined]
 [18] copy at /home/timon/.julia/packages/LazyArrays/lAt5J/src/linalg/inv.jl:100 [inlined]
 [19] copy(::Ldiv{ContinuumArrays.BasisLayout,LazyArrays.ApplyLayout{typeof(*)},SemiclassicalJacobi{Float64},Normalized{Float64,SubQuasiArray{Float64,2,Jacobi{Float64},Tuple{ContinuumArrays.AffineMap{Float64,Inclusion{Float64,DomainSets.UnitInterval{Float64}},Inclusion{Float64,ChebyshevInterval{Float64}}},Base.Slice{InfiniteArrays.OneToInf{Int64}}},false},Accumulate{Float64,1,typeof(*),Array{Float64,1},AbstractArray{Float64,1}}}}) at /home/timon/.julia/packages/ContinuumArrays/cNtkq/src/bases/bases.jl:200
 [20] copy at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:286 [inlined]
 [21] copy at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:283 [inlined]
 [22] materialize(::Ldiv{SemiclassicalOrthogonalPolynomials.SemiclassicalJacobiLayout,ClassicalOrthogonalPolynomials.NormalizedBasisLayout{ContinuumArrays.MappedBasisLayout},SemiclassicalJacobi{Float64},Normalized{Float64,SubQuasiArray{Float64,2,Jacobi{Float64},Tuple{ContinuumArrays.AffineMap{Float64,Inclusion{Float64,DomainSets.UnitInterval{Float64}},Inclusion{Float64,ChebyshevInterval{Float64}}},Base.Slice{InfiniteArrays.OneToInf{Int64}}},false},Accumulate{Float64,1,typeof(*),Array{Float64,1},AbstractArray{Float64,1}}}}) at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:22
 [23] ldiv at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:86 [inlined]
 [24] \(::SemiclassicalJacobi{Float64}, ::Normalized{Float64,SubQuasiArray{Float64,2,Jacobi{Float64},Tuple{ContinuumArrays.AffineMap{Float64,Inclusion{Float64,DomainSets.UnitInterval{Float64}},Inclusion{Float64,ChebyshevInterval{Float64}}},Base.Slice{InfiniteArrays.OneToInf{Int64}}},false},Accumulate{Float64,1,typeof(*),Array{Float64,1},AbstractArray{Float64,1}}}) at /home/timon/.julia/packages/QuasiArrays/BZtPi/src/matmul.jl:34
 [25] ldiv(::SemiclassicalJacobi{Float64}, ::Inclusion{Float64,DomainSets.UnitInterval{Float64}}) at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:362
 [26] \(::SemiclassicalJacobi{Float64}, ::Inclusion{Float64,DomainSets.UnitInterval{Float64}}) at /home/timon/.julia/packages/QuasiArrays/BZtPi/src/matmul.jl:34
 [27] top-level scope at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/test/runtests.jl:438
in expression starting at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/test/runtests.jl:438

julia> P \ x.^2
ERROR: LoadError: DomainError with -1.0989010989010988:
log1p will only return a complex result if called with a complex argument. Try log1p(Complex(x)).
Stacktrace:
 [1] throw_complex_domainerror(::Symbol, ::Float64) at ./math.jl:33
 [2] log1p(::Float64) at ./special/log.jl:356
 [3] expm1nlog1p at /home/timon/.julia/packages/HypergeometricFunctions/g0fPc/src/specialfunctions.jl:99 [inlined]
 [4] _₂F₁(::Float64, ::Float64, ::Float64, ::Float64) at /home/timon/.julia/packages/HypergeometricFunctions/g0fPc/src/gauss.jl:39
 [5] _₂F₁(::Float64, ::Float64, ::Float64, ::Float64) at /home/timon/.julia/packages/HypergeometricFunctions/g0fPc/src/gauss.jl:8
 [6] sum(::SemiclassicalJacobiWeight{Float64}) at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:49
 [7] massmatrix at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:343 [inlined]
 [8] copy(::Ldiv{SemiclassicalOrthogonalPolynomials.SemiclassicalJacobiLayout,SemiclassicalOrthogonalPolynomials.SemiclassicalJacobiLayout,SemiclassicalJacobi{Float64},SemiclassicalJacobi{Float64}}) at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:300
 [9] materialize(::Ldiv{SemiclassicalOrthogonalPolynomials.SemiclassicalJacobiLayout,SemiclassicalOrthogonalPolynomials.SemiclassicalJacobiLayout,SemiclassicalJacobi{Float64},SemiclassicalJacobi{Float64}}) at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:22
 [10] ldiv(::SemiclassicalJacobi{Float64}, ::SemiclassicalJacobi{Float64}) at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:86
 [11] \ at /home/timon/.julia/packages/QuasiArrays/BZtPi/src/matmul.jl:34 [inlined]
 [12] semijacobi_ldiv(::SemiclassicalJacobi{Float64}, ::SubQuasiArray{Float64,2,Jacobi{Float64},Tuple{ContinuumArrays.AffineMap{Float64,Inclusion{Float64,DomainSets.UnitInterval{Float64}},Inclusion{Float64,ChebyshevInterval{Float64}}},Base.Slice{InfiniteArrays.OneToInf{Int64}}},false}) at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:276
 [13] copy at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:290 [inlined]
 [14] materialize at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:22 [inlined]
 [15] ldiv at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:86 [inlined]
 [16] \ at /home/timon/.julia/packages/QuasiArrays/BZtPi/src/matmul.jl:34 [inlined]
 [17] _copy_ldiv_mul at /home/timon/.julia/packages/LazyArrays/lAt5J/src/linalg/inv.jl:98 [inlined]
 [18] copy at /home/timon/.julia/packages/LazyArrays/lAt5J/src/linalg/inv.jl:100 [inlined]
 [19] copy(::Ldiv{ContinuumArrays.BasisLayout,LazyArrays.ApplyLayout{typeof(*)},SemiclassicalJacobi{Float64},Normalized{Float64,SubQuasiArray{Float64,2,Jacobi{Float64},Tuple{ContinuumArrays.AffineMap{Float64,Inclusion{Float64,DomainSets.UnitInterval{Float64}},Inclusion{Float64,ChebyshevInterval{Float64}}},Base.Slice{InfiniteArrays.OneToInf{Int64}}},false},Accumulate{Float64,1,typeof(*),Array{Float64,1},AbstractArray{Float64,1}}}}) at /home/timon/.julia/packages/ContinuumArrays/cNtkq/src/bases/bases.jl:200
 [20] copy at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:286 [inlined]
 [21] copy at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:283 [inlined]
 [22] materialize(::Ldiv{SemiclassicalOrthogonalPolynomials.SemiclassicalJacobiLayout,ClassicalOrthogonalPolynomials.NormalizedBasisLayout{ContinuumArrays.MappedBasisLayout},SemiclassicalJacobi{Float64},Normalized{Float64,SubQuasiArray{Float64,2,Jacobi{Float64},Tuple{ContinuumArrays.AffineMap{Float64,Inclusion{Float64,DomainSets.UnitInterval{Float64}},Inclusion{Float64,ChebyshevInterval{Float64}}},Base.Slice{InfiniteArrays.OneToInf{Int64}}},false},Accumulate{Float64,1,typeof(*),Array{Float64,1},AbstractArray{Float64,1}}}}) at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:22
 [23] ldiv at /home/timon/.julia/packages/ArrayLayouts/wBWqq/src/ldiv.jl:86 [inlined]
 [24] \(::SemiclassicalJacobi{Float64}, ::Normalized{Float64,SubQuasiArray{Float64,2,Jacobi{Float64},Tuple{ContinuumArrays.AffineMap{Float64,Inclusion{Float64,DomainSets.UnitInterval{Float64}},Inclusion{Float64,ChebyshevInterval{Float64}}},Base.Slice{InfiniteArrays.OneToInf{Int64}}},false},Accumulate{Float64,1,typeof(*),Array{Float64,1},AbstractArray{Float64,1}}}) at /home/timon/.julia/packages/QuasiArrays/BZtPi/src/matmul.jl:34
 [25] ldiv(::SemiclassicalJacobi{Float64}, ::BroadcastQuasiArray{Float64,1,typeof(Base.literal_pow),Tuple{Base.RefValue{typeof(^)},Inclusion{Float64,DomainSets.UnitInterval{Float64}},Base.RefValue{Val{2}}}}) at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/src/SemiclassicalOrthogonalPolynomials.jl:362
 [26] \(::SemiclassicalJacobi{Float64}, ::BroadcastQuasiArray{Float64,1,typeof(Base.literal_pow),Tuple{Base.RefValue{typeof(^)},Inclusion{Float64,DomainSets.UnitInterval{Float64}},Base.RefValue{Val{2}}}}) at /home/timon/.julia/packages/QuasiArrays/BZtPi/src/matmul.jl:34
 [27] top-level scope at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/test/runtests.jl:438
in expression starting at /home/timon/Documents/Projects/SemiclassicalOrthogonalPolynomials.jl/test/runtests.jl:438

Some other random examples for a value being returned:

julia> P = SemiclassicalJacobi(0.89,1,3,1)
SemiclassicalJacobi with weight x^1.0 * (1-x)^3.0 * (0.89-x)^1.0

julia> P[0.8,4]
10.435052411380479

julia> P = SemiclassicalJacobi(0.89,1,3,2)
SemiclassicalJacobi with weight x^1.0 * (1-x)^3.0 * (0.89-x)^2.0

julia> P[0.8,4]
18.75150830259011

Easy fix by just making SemiclassicalJacobi check t>1 when not being built from a pre-existing basis?

@dlfivefifty
Copy link
Member

👍

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants