Clenshaw2D

Project.toml

@@ -1,11 +1,12 @@
 name = "RecurrenceRelationshipArrays"
 uuid = "b889d2dc-af3c-4820-88a8-238fa91d3518"
 authors = ["Sheehan Olver <[email protected]>"]
-version = "0.1.2"
+version = "0.1.3"
 
 [deps]
 ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"
 BandedMatrices = "aae01518-5342-5314-be14-df237901396f"
+BlockBandedMatrices = "ffab5731-97b5-5995-9138-79e8c1846df0"
 FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b"
 InfiniteArrays = "4858937d-0d70-526a-a4dd-2d5cb5dd786c"
 LazyArrays = "5078a376-72f3-5289-bfd5-ec5146d43c02"
@@ -15,11 +16,12 @@ RecurrenceRelationships = "807425ed-42ea-44d6-a357-6771516d7b2c"
 [compat]
 ArrayLayouts = "1"
 BandedMatrices = "1"
+BlockBandedMatrices = "0.13.3"
 DynamicPolynomials = "0.6.1"
 FillArrays = "1"
 InfiniteArrays = "0.14, 0.15"
 LazyArrays = "2"
-RecurrenceRelationships = "0.1"
+RecurrenceRelationships = "0.2"
 SpecialFunctions = "1, 2"
 julia = "1.10"
 

src/RecurrenceRelationshipArrays.jl

@@ -1,13 +1,13 @@
 module RecurrenceRelationshipArrays
-using LinearAlgebra, InfiniteArrays, LazyArrays, FillArrays, ArrayLayouts, BandedMatrices
+using LinearAlgebra, InfiniteArrays, LazyArrays, FillArrays, ArrayLayouts, BandedMatrices, BlockBandedMatrices
 import LazyArrays: AbstractCachedArray, LazyArrayStyle, AbstractPaddedLayout, resizedata!, zero!, paddeddata
 import Base: axes, size, getindex, broadcasted, copy, view, +, -, Slice, tail
 import FillArrays: AbstractFill, getindex_value
 import BandedMatrices: AbstractBandedMatrix, bandwidths, _BandedMatrix
 import InfiniteArrays: OneToInf, AbstractInfUnitRange
 import ArrayLayouts: MatMulVecAdd, sublayout, MemoryLayout, sub_materialize, transposelayout, materialize!, _fill_lmul!
 import RecurrenceRelationships: forwardrecurrence_next, forwardrecurrence_partial!, check_clenshaw_recurrences, clenshaw, polynomialtype
-export RecurrenceArray, Clenshaw
+export RecurrenceArray, Clenshaw, ClenshawKron
 
 const LazyArraysBandedMatricesExt = Base.get_extension(LazyArrays, :LazyArraysBandedMatricesExt)
 const AbstractLazyBandedLayout = LazyArraysBandedMatricesExt.AbstractLazyBandedLayout

src/clenshaw.jl

@@ -233,4 +233,5 @@ end
 # Banded dot is slow
 ###
 
-LinearAlgebra.dot(x::AbstractVector, A::Clenshaw, y::AbstractVector) = dot(x, mul(A, y))
+LinearAlgebra.dot(x::AbstractVector, A::Clenshaw, y::AbstractVector) = dot(x, mul(A, y))
+

src/recurrence.jl

@@ -24,7 +24,7 @@ const RecurrenceMatrix{T, Z<:AbstractVector, A<:AbstractVector, B<:AbstractVecto
 RecurrenceArray(z, A, B, C, data::Array{T,N}, datasize, p0, p1) where {T,N} = RecurrenceArray{T,N,typeof(z),typeof(A),typeof(B),typeof(C)}(z, A, B, C, data, datasize, p0, p1, T[])
 
 function initiateforwardrecurrence(N, A, B, C, x, μ)
-    T = polynomialtype(eltype(A), eltype(B), eltype(C), typeof(x))
+    T = polynomialtype(promote_type(eltype(A), eltype(B), eltype(C)), typeof(x))
     p0 = convert(T, μ)
     N == 0 && return zero(T), p0
     p1 = convert(T, muladd(A[1],x,B[1])*p0)

test/runtests.jl

@@ -50,6 +50,13 @@ end
     end
 end
 
+@testset "ClenshawKron" begin
+    C = [1 2 3; 4 5 0; 6 0 0]
+    X_U = SymTridiagonal(Fill(0.0,∞), Fill(1/2,∞))
+    ClenshawKron(C, (rec_U, rec_U), (X_U,X_U))
+end
+
+
 @testset "DynamicPolynomials" begin
     @polyvar x
     @test @inferred(RecurrenceRelationshipArrays.initiateforwardrecurrence(5, rec_T..., x, 1))[2] == 5x - 20x^3 + 16x^5