search

JuliaApproximation / MultivariateOrthogonalPolynomials.jl

Supports approximating functions and solving differential equations on various higher dimensional domains such as disks and triangles
Other
17 stars 5 forks source link

JacobiTriangle() expansions failing in strange way #140

Open TSGut opened 1 year ago

TSGut commented 1 year ago

There is some strange behavior when expanding functions on triangles, apparently when it comes to globally defined variables? Here is a minimal example:

julia> using MultivariateOrthogonalPolynomials

julia> T = JacobiTriangle()
JacobiTriangle(0, 0, 0)

julia> xy = axes(T,1)
Inclusion(UnitSimplex(Val(2)))

julia> x, y = first.(xy), last.(xy)
(first.(Inclusion(UnitSimplex(Val(2)))), last.(Inclusion(UnitSimplex(Val(2)))))

julia> α = 2.
2.0

julia> K(x,y) = α*x*y
K (generic function with 1 method)

julia> K(α,x,y) = α*x*y
K (generic function with 2 methods)

julia> T \ (@. α*x*y)  # this works
ℵ₀-blocked ℵ₀-element LazyBandedMatrices.DiagTrav{Float64, 2, LazyArrays.ApplyArray{Float64, 2, typeof(Base.setindex), Tuple{FillArrays.Zeros{Float64, 2, Tuple{InfiniteArrays.OneToInf{Int64}, InfiniteArrays.OneToInf{Int64}}}, Matrix{Float64}, Base.OneTo{Int64}, Base.OneTo{Int64}}}}:
  0.1666666666666667    
 ───────────────────────
  0.06666666666666667   
  0.19999999999999998   
 ───────────────────────
 -0.09999999999999999   
  0.19999999999999993   
 -6.699811370772778e-18 
 ───────────────────────
  0.0                   
 -6.148984834775253e-18 
 -1.6452652161292338e-18
  ⋮

julia> T \ (@. K(α,x,y))  # this works
ℵ₀-blocked ℵ₀-element LazyBandedMatrices.DiagTrav{Float64, 2, LazyArrays.ApplyArray{Float64, 2, typeof(Base.setindex), Tuple{FillArrays.Zeros{Float64, 2, Tuple{InfiniteArrays.OneToInf{Int64}, InfiniteArrays.OneToInf{Int64}}}, Matrix{Float64}, Base.OneTo{Int64}, Base.OneTo{Int64}}}}:
  0.1666666666666667    
 ───────────────────────
  0.06666666666666667   
  0.19999999999999998   
 ───────────────────────
 -0.09999999999999999   
  0.19999999999999993   
 -6.699811370772778e-18 
 ───────────────────────
  0.0                   
 -6.148984834775253e-18 
 -1.6452652161292338e-18
  ⋮

julia> T \ (@. K(x,y))  # this fails
ERROR: MethodError: no method matching *(::FastTransforms.FTPlan{Float64, 2, FastTransforms.TRIANGLEANALYSIS}, ::Matrix{Any})
Closest candidates are:
  *(::Any, ::Any, ::Any, ::Any...) at operators.jl:591
  *(::Union{SparseArrays.AbstractSparseMatrixCSC{TA, Ti}, SubArray{TA, 2, <:SparseArrays.AbstractSparseMatrixCSC{TA, Ti}, Tuple{Base.Slice{Base.OneTo{Int64}}, I}} where I<:AbstractUnitRange} where Ti, ::Union{LinearAlgebra.Adjoint{<:Any, <:Union{LinearAlgebra.LowerTriangular, LinearAlgebra.UnitLowerTriangular, LinearAlgebra.UnitUpperTriangular, LinearAlgebra.UpperTriangular, StridedMatrix, BitMatrix}}, LinearAlgebra.LowerTriangular, LinearAlgebra.Transpose{<:Any, <:Union{LinearAlgebra.LowerTriangular, LinearAlgebra.UnitLowerTriangular, LinearAlgebra.UnitUpperTriangular, LinearAlgebra.UpperTriangular, StridedMatrix, BitMatrix}}, LinearAlgebra.UnitLowerTriangular, LinearAlgebra.UnitUpperTriangular, LinearAlgebra.UpperTriangular, StridedMatrix, BitMatrix}) where TA at ~/Documents/Backends/julia-1.8.5/share/julia/stdlib/v1.8/SparseArrays/src/linalg.jl:52
  *(::LinearAlgebra.Hermitian{<:Any, <:ArrayLayouts.LayoutMatrix}, ::AbstractMatrix) at ~/.julia/packages/ArrayLayouts/4IG3b/src/mul.jl:263
  ...
Stacktrace:
  [1] *(T::MultivariateOrthogonalPolynomials.TriPlan{Float64}, F::Matrix{Any})
    @ MultivariateOrthogonalPolynomials ~/.julia/packages/MultivariateOrthogonalPolynomials/tNdGo/src/triangle.jl:661
  [2] \(a::ContinuumArrays.TransformFactorization{Any, Matrix{StaticArraysCore.SVector{2, Float64}}, MultivariateOrthogonalPolynomials.TriPlan{Float64}}, b::QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}})
    @ ContinuumArrays ~/.julia/packages/ContinuumArrays/pWRm2/src/plans.jl:26
  [3] transform_ldiv(V::QuasiArrays.SubQuasiArray{Float64, 2, JacobiTriangle{Float64, Int64}, Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}, BlockArrays.BlockSlice{BlockArrays.BlockRange{1, Tuple{Base.OneTo{Int64}}}, BlockArrays.BlockedUnitRange{ArrayLayouts.RangeCumsum{Int64, Base.OneTo{Int64}}}}}, false}, B::QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}}, #unused#::Tuple{Infinities.InfiniteCardinal{1}, Int64})
    @ HarmonicOrthogonalPolynomials ~/.julia/packages/HarmonicOrthogonalPolynomials/51vXS/src/multivariateops.jl:91
  [4] transform_ldiv(A::QuasiArrays.SubQuasiArray{Float64, 2, JacobiTriangle{Float64, Int64}, Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}, BlockArrays.BlockSlice{BlockArrays.BlockRange{1, Tuple{Base.OneTo{Int64}}}, BlockArrays.BlockedUnitRange{ArrayLayouts.RangeCumsum{Int64, Base.OneTo{Int64}}}}}, false}, B::QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}})
    @ ContinuumArrays ~/.julia/packages/ContinuumArrays/pWRm2/src/bases/bases.jl:258
  [5] transform_ldiv_if_columns(L::ArrayLayouts.Ldiv{ContinuumArrays.SubBasisLayout, LazyArrays.BroadcastLayout{typeof(K)}, QuasiArrays.SubQuasiArray{Float64, 2, JacobiTriangle{Float64, Int64}, Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}, BlockArrays.BlockSlice{BlockArrays.BlockRange{1, Tuple{Base.OneTo{Int64}}}, BlockArrays.BlockedUnitRange{ArrayLayouts.RangeCumsum{Int64, Base.OneTo{Int64}}}}}, false}, QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}}}, #unused#::Base.OneTo{Int64})
    @ ContinuumArrays ~/.julia/packages/ContinuumArrays/pWRm2/src/bases/bases.jl:289
  [6] transform_ldiv_if_columns(L::ArrayLayouts.Ldiv{ContinuumArrays.SubBasisLayout, LazyArrays.BroadcastLayout{typeof(K)}, QuasiArrays.SubQuasiArray{Float64, 2, JacobiTriangle{Float64, Int64}, Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}, BlockArrays.BlockSlice{BlockArrays.BlockRange{1, Tuple{Base.OneTo{Int64}}}, BlockArrays.BlockedUnitRange{ArrayLayouts.RangeCumsum{Int64, Base.OneTo{Int64}}}}}, false}, QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}}})
    @ ContinuumArrays ~/.julia/packages/ContinuumArrays/pWRm2/src/bases/bases.jl:290
  [7] copy(L::ArrayLayouts.Ldiv{ContinuumArrays.SubBasisLayout, LazyArrays.BroadcastLayout{typeof(K)}, QuasiArrays.SubQuasiArray{Float64, 2, JacobiTriangle{Float64, Int64}, Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}, BlockArrays.BlockSlice{BlockArrays.BlockRange{1, Tuple{Base.OneTo{Int64}}}, BlockArrays.BlockedUnitRange{ArrayLayouts.RangeCumsum{Int64, Base.OneTo{Int64}}}}}, false}, QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}}})
    @ ContinuumArrays ~/.julia/packages/ContinuumArrays/pWRm2/src/bases/bases.jl:291
  [8] materialize(M::ArrayLayouts.Ldiv{ContinuumArrays.SubBasisLayout, LazyArrays.BroadcastLayout{typeof(K)}, QuasiArrays.SubQuasiArray{Float64, 2, JacobiTriangle{Float64, Int64}, Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}, BlockArrays.BlockSlice{BlockArrays.BlockRange{1, Tuple{Base.OneTo{Int64}}}, BlockArrays.BlockedUnitRange{ArrayLayouts.RangeCumsum{Int64, Base.OneTo{Int64}}}}}, false}, QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}}}; kwds::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ ArrayLayouts ~/.julia/packages/ArrayLayouts/4IG3b/src/ldiv.jl:22
  [9] materialize
    @ ~/.julia/packages/ArrayLayouts/4IG3b/src/ldiv.jl:22 [inlined]
 [10] #ldiv#18
    @ ~/.julia/packages/ArrayLayouts/4IG3b/src/ldiv.jl:86 [inlined]
 [11] ldiv
    @ ~/.julia/packages/ArrayLayouts/4IG3b/src/ldiv.jl:86 [inlined]
 [12] \
    @ ~/.julia/packages/QuasiArrays/tF3vb/src/matmul.jl:34 [inlined]
 [13] adaptivetransform_ldiv(A::JacobiTriangle{Float64, Int64}, f::QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}})
    @ ClassicalOrthogonalPolynomials ~/.julia/packages/ClassicalOrthogonalPolynomials/ix48P/src/adaptivetransform.jl:35
 [14] transform_ldiv(A::JacobiTriangle{Float64, Int64}, f::QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}}, #unused#::Tuple{Infinities.InfiniteCardinal{1}, Infinities.InfiniteCardinal{0}})
    @ ClassicalOrthogonalPolynomials ~/.julia/packages/ClassicalOrthogonalPolynomials/ix48P/src/adaptivetransform.jl:2
 [15] transform_ldiv(A::JacobiTriangle{Float64, Int64}, B::QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}})
    @ ContinuumArrays ~/.julia/packages/ContinuumArrays/pWRm2/src/bases/bases.jl:258
 [16] transform_ldiv_if_columns(L::ArrayLayouts.Ldiv{HarmonicOrthogonalPolynomials.MultivariateOPLayout{2}, LazyArrays.BroadcastLayout{typeof(K)}, JacobiTriangle{Float64, Int64}, QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}}}, #unused#::Base.OneTo{Int64})
    @ ContinuumArrays ~/.julia/packages/ContinuumArrays/pWRm2/src/bases/bases.jl:289
 [17] transform_ldiv_if_columns(L::ArrayLayouts.Ldiv{HarmonicOrthogonalPolynomials.MultivariateOPLayout{2}, LazyArrays.BroadcastLayout{typeof(K)}, JacobiTriangle{Float64, Int64}, QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}}})
    @ ContinuumArrays ~/.julia/packages/ContinuumArrays/pWRm2/src/bases/bases.jl:290
 [18] copy(L::ArrayLayouts.Ldiv{HarmonicOrthogonalPolynomials.MultivariateOPLayout{2}, LazyArrays.BroadcastLayout{typeof(K)}, JacobiTriangle{Float64, Int64}, QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}}})
    @ ContinuumArrays ~/.julia/packages/ContinuumArrays/pWRm2/src/bases/bases.jl:291
 [19] #materialize#11
    @ ~/.julia/packages/ArrayLayouts/4IG3b/src/ldiv.jl:22 [inlined]
 [20] materialize(M::ArrayLayouts.Ldiv{HarmonicOrthogonalPolynomials.MultivariateOPLayout{2}, LazyArrays.BroadcastLayout{typeof(K)}, JacobiTriangle{Float64, Int64}, QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}}})
    @ ArrayLayouts ~/.julia/packages/ArrayLayouts/4IG3b/src/ldiv.jl:22
 [21] ldiv(A::JacobiTriangle{Float64, Int64}, B::QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}}; kwds::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ ArrayLayouts ~/.julia/packages/ArrayLayouts/4IG3b/src/ldiv.jl:86
 [22] ldiv
    @ ~/.julia/packages/ArrayLayouts/4IG3b/src/ldiv.jl:86 [inlined]
 [23] \(A::JacobiTriangle{Float64, Int64}, B::QuasiArrays.BroadcastQuasiVector{Any, typeof(K), Tuple{QuasiArrays.BroadcastQuasiVector{Float64, typeof(first), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}, QuasiArrays.BroadcastQuasiVector{Float64, typeof(last), Tuple{QuasiArrays.Inclusion{StaticArraysCore.SVector{2, Float64}, DomainSets.EuclideanUnitSimplex{2, Float64, :closed}}}}}})
    @ QuasiArrays ~/.julia/packages/QuasiArrays/tF3vb/src/matmul.jl:34
 [24] top-level scope
    @ ~/Documents/Projects/SparseVolterraExamples.jl/examples/Example 0 - Some linear Volterra integral equation examples.jl:82

Any idea what is happening here?

HadrienNU commented 1 year ago

This follow from the type instability of global variable. From your example

julia> typeof(K.(α,x,y))
BroadcastQuasiVector{Float64, typeof(K), Tuple{Float64, BroadcastQuasiVector{Float64, typeof(first), Tuple{Inclusion{SVector{2, Float64}, EuclideanUnitSimplex{2, Float64, :closed}}}}, BroadcastQuasiVector{Float64, typeof(last), Tuple{Inclusion{SVector{2, Float64}, EuclideanUnitSimplex{2, Float64, :closed}}}}}}

whereas

julia> typeof(K.(x,y))
BroadcastQuasiVector{Any, typeof(K), Tuple{BroadcastQuasiVector{Float64, typeof(first), Tuple{Inclusion{SVector{2, Float64}, EuclideanUnitSimplex{2, Float64, :closed}}}}, BroadcastQuasiVector{Float64, typeof(last), Tuple{Inclusion{SVector{2, Float64}, EuclideanUnitSimplex{2, Float64, :closed}}}}}}

Note the Any in the second line of code. However if α is declare as a constant 

const α =2.0

it works

julia> T \ (@. K(x,y))
ℵ₀-blocked ℵ₀-element LazyBandedMatrices.DiagTrav{Float64, 2, LazyArrays.ApplyArray{Float64, 2, typeof(Base.setindex), Tuple{FillArrays.Zeros{Float64, 2, Tuple{InfiniteArrays.OneToInf{Int64}, InfiniteArrays.OneToInf{Int64}}}, Matrix{Float64}, Base.OneTo{Int64}, Base.OneTo{Int64}}}}:
  0.1666666666666667    
 ───────────────────────
  0.06666666666666667   
  0.19999999999999998   
 ───────────────────────
 -0.09999999999999999   
  0.19999999999999993   
 -6.699811370772778e-18 
 ───────────────────────
  0.0                   
 -6.148984834775253e-18 
 -1.6452652161292338e-18
  0.0                   
 ───────────────────────
  ⋮
dlfivefifty commented 1 year ago

Note I fixed essentially the same bug back in August: https://github.com/JuliaApproximation/ClassicalOrthogonalPolynomials.jl/commit/2e33f50b7326706775597c75f284eb57e3dc0597

So I suspect something similar needs to be done here, just need to step through the transforms and see where the promotion to Any happens.