From 757fb9ae29518abf761971fe48593530268fc634 Mon Sep 17 00:00:00 2001
From: Daniel VandenHeuvel <95613936+DanielVandH@users.noreply.github.com>
Date: Fri, 31 Oct 2025 10:54:17 +0000
Subject: [PATCH] Function varriers

---
 src/SemiclassicalOrthogonalPolynomials.jl |  2 +-
 src/derivatives.jl                        | 18 ++++++++++++------
 test/test_derivative.jl                   | 10 +++++++++-
 3 files changed, 22 insertions(+), 8 deletions(-)

diff --git a/src/SemiclassicalOrthogonalPolynomials.jl b/src/SemiclassicalOrthogonalPolynomials.jl
index 9dcb571..3920723 100644
--- a/src/SemiclassicalOrthogonalPolynomials.jl
+++ b/src/SemiclassicalOrthogonalPolynomials.jl
@@ -7,7 +7,7 @@ import Base: getindex, axes, size, \, /, *, +, -, summary, show, ==, copy, sum,
 
 import ArrayLayouts: MemoryLayout, ldiv, diagonaldata, subdiagonaldata, supdiagonaldata
 import BandedMatrices: bandwidths, AbstractBandedMatrix, BandedLayout, _BandedMatrix
-import LazyArrays: resizedata!, paddeddata, CachedVector, CachedMatrix, CachedAbstractVector, LazyMatrix, LazyVector, arguments, ApplyLayout, colsupport, AbstractCachedVector, ApplyArray,
+import LazyArrays: resizedata!, paddeddata, maybe_cacheddata, CachedVector, CachedMatrix, CachedAbstractVector, LazyMatrix, LazyVector, arguments, ApplyLayout, colsupport, AbstractCachedVector, ApplyArray,
                     AccumulateAbstractVector, LazyVector, AbstractCachedMatrix, BroadcastLayout, simplifiable
 import ClassicalOrthogonalPolynomials: OrthogonalPolynomial, recurrencecoefficients, jacobimatrix, normalize, _p0, UnitInterval, orthogonalityweight, NormalizedOPLayout, MappedOPLayout,
                                         Bidiagonal, Tridiagonal, SymTridiagonal, symtridiagonalize, normalizationconstant, LanczosPolynomial,
diff --git a/src/derivatives.jl b/src/derivatives.jl
index 6cfd25b..81d503b 100644
--- a/src/derivatives.jl
+++ b/src/derivatives.jl
@@ -1,7 +1,7 @@
 """
-MulAddAccumulate(μ, A, B)
+    MulAddAccumulate(μ, A, B)
 
-represents the vector satisfying v[k+1] == A[k]*v[k]+B[k] with v[1] == μ
+represents the vector satisfying v[k+1] == A[k+1]*v[k]+B[k+1] with v[1] == μ
 """
 mutable struct MulAddAccumulate{T} <: AbstractCachedVector{T}
     data::Vector{T}
@@ -21,11 +21,17 @@ MulAddAccumulate(data::Vector{T}, A::AbstractVector, B::AbstractVector, datasize
 MulAddAccumulate(μ, A, B) = MulAddAccumulate([μ], A, B, (1,))
 MulAddAccumulate(A, B) = MulAddAccumulate(A[1]+B[1], A, B)
 
-function LazyArrays.cache_filldata!(K::MulAddAccumulate, inds)
-    A,B = K.A,K.B
+LazyArrays.cache_filldata!(K::MulAddAccumulate, inds) = _muladdaccumulate_filldata!(K.data, K.A, K.B, inds)
+function _muladdaccumulate_filldata!(K, A, B, inds)
+    n = maximum(inds)
+    resizedata!(A, n)
+    resizedata!(B, n)
+    __muladdaccumulate_filldata!(K, maybe_cacheddata(A), maybe_cacheddata(B), inds)
+end
+function __muladdaccumulate_filldata!(data, A, B, inds)
     @inbounds for k in inds
-        K.data[k] = muladd(A[k], K.data[k-1], B[k])
-    end
+        data[k] = muladd(A[k], data[k-1], B[k])
+    end 
 end
 
 """
diff --git a/test/test_derivative.jl b/test/test_derivative.jl
index 2ea3d8b..ae84ef7 100644
--- a/test/test_derivative.jl
+++ b/test/test_derivative.jl
@@ -1,7 +1,15 @@
 using SemiclassicalOrthogonalPolynomials, ClassicalOrthogonalPolynomials, LazyArrays, Test
 import ClassicalOrthogonalPolynomials: recurrencecoefficients, _BandedMatrix, _p0, Weighted
 import LazyArrays: Accumulate, AccumulateAbstractVector
-import SemiclassicalOrthogonalPolynomials: MulAddAccumulate, HalfWeighted, toclassical
+import SemiclassicalOrthogonalPolynomials: MulAddAccumulate, HalfWeighted, toclassical, MulAddAccumulate
+
+@testset "MulAddAccumulate" begin 
+    A = Vcat(1, cache(2:3), cache(4:∞) ./ cache(3:∞))
+    B = Vcat(cache(1:5), [1, 2], cache(3:∞))
+    v = MulAddAccumulate(A, B)
+    @test v[1] ≈ 2
+    @test all(k -> v[k+1] ≈ A[k+1] * v[k] + B[k+1], 2:1000)
+end
 
 @testset "Derivative" begin
     @testset "basics" begin