From 55ee7d92db6ec829a9fbb23e2c0d223d403b4bce Mon Sep 17 00:00:00 2001 From: Sheehan Olver Date: Sun, 11 Oct 2026 10:52:41 +0100 Subject: [PATCH 1/7] Support FillArrays.OffDiagonal and remove sublayout piracy for DiagonalLayout Require ArrayLayouts v1.15, which adds OffDiagonalLayout and returns it for views of a DiagonalLayout matrix by unit ranges. Remove the pirated sublayout(::DiagonalLayout, ...) methods, which would otherwise overwrite ArrayLayouts' method, and add bandwidths, inbands_getindex/inbands_setindex! and BandedMatrix conversion for OffDiagonalLayout. Co-Authored-By: Claude Opus 5.5 --- Project.toml | 6 +++--- src/BandedMatrices.jl | 5 +++-- src/banded/BandedMatrix.jl | 9 +++++++++ src/interfaceimpl.jl | 9 ++++----- test/test_interface.jl | 35 ++++++++++++++++++++++++++++++++--- 5 files changed, 51 insertions(+), 13 deletions(-) diff --git a/Project.toml b/Project.toml index a7d5c29c..c6fc0e14 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "BandedMatrices" uuid = "aae01518-5342-5314-be14-df237901396f" -version = "1.13.1" +version = "1.14.0" [deps] ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a" @@ -18,11 +18,11 @@ CliqueTreesExt = "CliqueTrees" [compat] Aqua = "0.8" -ArrayLayouts = "1.6.0" +ArrayLayouts = "1.15" BlockArrays = "1" CliqueTrees = "0.5, 1" Documenter = "1" -FillArrays = "1.3" +FillArrays = "1.18" GenericLinearAlgebra = "0.3, 0.4" InfiniteArrays = "0.12, 0.13, 0.14, 0.15" LinearAlgebra = "1" diff --git a/src/BandedMatrices.jl b/src/BandedMatrices.jl index f0a20e29..65779f20 100644 --- a/src/BandedMatrices.jl +++ b/src/BandedMatrices.jl @@ -32,9 +32,10 @@ import ArrayLayouts: AbstractTridiagonalLayout, BidiagonalLayout, BlasMatLdivVec hermitianlayout, materialize, materialize!, reflector!, reflectorApply!, rowsupport, sub_materialize, subdiagonaldata, sublayout, supdiagonaldata, symmetricdata, symmetriclayout, symmetricuplo, transposelayout, triangulardata, triangularlayout, zero!, - QRPackedQLayout, AdjQRPackedQLayout, LayoutVecOrMats + QRPackedQLayout, AdjQRPackedQLayout, LayoutVecOrMats, OffDiagonalLayout, offdiagonaldata, + offdiagonalindex -import FillArrays: AbstractFill, getindex_value, _broadcasted_zeros, unique_value, OneElement, RectDiagonal, OneElementMatrix, OneElementVector, ZerosMatrix, ZerosVector +import FillArrays: AbstractFill, getindex_value, _broadcasted_zeros, unique_value, OneElement, RectDiagonal, OffDiagonal, OneElementMatrix, OneElementVector, ZerosMatrix, ZerosVector const libblas = LinearAlgebra.BLAS.libblas const liblapack = LinearAlgebra.BLAS.liblapack diff --git a/src/banded/BandedMatrix.jl b/src/banded/BandedMatrix.jl index bb0bb0bc..e3adbbbf 100644 --- a/src/banded/BandedMatrix.jl +++ b/src/banded/BandedMatrix.jl @@ -295,6 +295,15 @@ function _BandedMatrix(::DiagonalLayout, A::AbstractMatrix{T}) where T copyto!(vec(dat), diagonaldata(A)) _BandedMatrix(dat, m, 0, 0) end +function _BandedMatrix(::OffDiagonalLayout, A::AbstractMatrix{T}) where T + m,n = size(A) + k = offdiagonalindex(A) + dat = zeros(T, 1, n) + # the entry in column j lies in row j-k + band = offdiagonaldata(A) + copyto!(view(dat, 1, max(0,k) .+ axes(band,1)), band) + _BandedMatrix(dat, m, -k, k) +end function _BandedMatrix(::BidiagonalLayout, A::AbstractMatrix{T}) where T m,n = size(A) dat = Matrix{T}(undef, 2, n) diff --git a/src/interfaceimpl.jl b/src/interfaceimpl.jl index 3a3a7b15..7ab474b7 100644 --- a/src/interfaceimpl.jl +++ b/src/interfaceimpl.jl @@ -35,11 +35,10 @@ bandedrowsdata(Ac::Transpose{<:Number}) = Transpose(bandeddata(parent(Ac))) bandedrowsdata(Ac::Adjoint{<:Real}) = Transpose(bandeddata(parent(Ac))) bandedrowsdata(V::SubArray{<:Number}) = Transpose(bandeddata(view(parent(parent(V)), reverse(parentindices(V))...))) -# treat subinds as banded -sublayout(::DiagonalLayout{L}, inds::Type) where L = sublayout(bandedcolumns(L()), inds) -sublayout(::DiagonalLayout{L}, inds::Type{<:NTuple{2,AbstractUnitRange{Int}}}) where L = sublayout(bandedcolumns(L()), inds) - -# bandeddata(V::SubArray{<:Any,2,<:Diagonal}) = view(bandeddata(parent(V)), :, parentindices(V)[2]) +# OffDiagonal, and views of Diagonal/OffDiagonal by unit ranges, are supported on a single band +bandwidths(::OffDiagonalLayout, A) = (-offdiagonalindex(A), offdiagonalindex(A)) +inbands_getindex(A::OffDiagonal, k::Integer, j::Integer) = A.band[min(k,j)] +inbands_setindex!(A::OffDiagonal, v, k::Integer, j::Integer) = (A.band[min(k,j)] = v) isbanded(::Bidiagonal) = true bandwidths(A::Bidiagonal) = A.uplo == 'U' ? (0,1) : (1,0) diff --git a/test/test_interface.jl b/test/test_interface.jl index 21d58538..30a0d460 100644 --- a/test/test_interface.jl +++ b/test/test_interface.jl @@ -53,7 +53,7 @@ LinearAlgebra.fill!(A::PseudoBandedMatrix, v) = fill!(A.data,v) V = view(Eye(5),2:3,:) @test bandwidths(V) == (-1,1) - @test MemoryLayout(typeof(V)) isa BandedColumns{OnesLayout} + @test MemoryLayout(typeof(V)) isa OffDiagonalLayout{OnesLayout} @test BandedMatrix(V) == [0 1 0 0 0; 0 0 1 0 0] B = brand(5,5,1,2) @@ -79,7 +79,9 @@ LinearAlgebra.fill!(A::PseudoBandedMatrix, v) = fill!(A.data,v) @test bandeddata(A) == bandeddata(Adjoint(A)) == bandeddata(Transpose(A)) == diag(A)' - @test MemoryLayout(view(A, 1:3,2:4)) isa BandedColumns{DenseColumnMajor} + @test MemoryLayout(view(A, 1:3,2:4)) isa OffDiagonalLayout{DenseColumnMajor} + @test bandwidths(view(A, 1:3,2:4)) == (1,-1) + @test BandedMatrix(view(A, 1:3,2:4)) == A[1:3,2:4] @test MemoryLayout(view(A, [1,2,3],2:4)) isa UnknownLayout A[band(0)][1] = 3 @@ -99,6 +101,32 @@ LinearAlgebra.fill!(A::PseudoBandedMatrix, v) = fill!(A.data,v) @test @inferred(B[Band(0)]) == 1:∞ end + @testset "OffDiagonal" begin + for A in (FillArrays.OffDiagonal([1.,2,3], 1), FillArrays.OffDiagonal([1.,2,3], -2), FillArrays.OffDiagonal(Fill(2.0,3), 1, 3, 5), + FillArrays.OffDiagonal([1.,2,3], 1)', view(Diagonal([1.,2,3,4,5]), 2:5, 1:3)) + k = ArrayLayouts.offdiagonalindex(A) + @test isbanded(A) + @test bandwidths(A) == (-k,k) + @test BandedMatrix(A) isa BandedMatrix + @test BandedMatrix(A) == A + @test bandwidths(BandedMatrix(A)) == (-k,k) + B = brand(size(A,2), 4, 1, 1) + C = brand(4, size(A,1), 1, 1) + @test A*B ≈ Matrix(A)*B + @test C*A ≈ C*Matrix(A) + end + + A = FillArrays.OffDiagonal([1.,2,3], 1) + @test MemoryLayout(A) isa OffDiagonalLayout{DenseColumnMajor} + @test BandedMatrices.inbands_getindex(A, 2, 3) == 2 + BandedMatrices.inbands_setindex!(A, 4, 2, 3) + @test A[2,3] == 4 + @test A[band(1)] == [1,4,3] + @test A[band(0)] == zeros(4) + @test MemoryLayout(view(A, 1:3, 2:4)) isa OffDiagonalLayout{DenseColumnMajor} + @test bandwidths(view(A, 1:3, 2:4)) == (0,0) + end + @testset "SymTridiagonal" begin A = SymTridiagonal([1,2,3],[4,5]) @test @inferred(A[Band(0)]) == [1,2,3] @@ -303,7 +331,8 @@ end @test bandwidths(A) == bandwidths(Ã) == (0,2) @test à == A == Matrix(B) + Matrix(Bu) - @test layout_getindex(D,1:10,1:10) isa BandedMatrix + @test layout_getindex(D,1:10,1:10) isa FillArrays.OffDiagonal + @test layout_getindex(D,1:10,1:10) == D @test layout_getindex(B,1:10,1:10) isa BandedMatrix @test layout_getindex(Bu,1:10,1:10) isa BandedMatrix @test layout_getindex(J,1:10,1:10) isa BandedMatrix From 74210996a7ac93f8afa848bb62daffe0a777e5af Mon Sep 17 00:00:00 2001 From: Sheehan Olver Date: Mon, 12 Oct 2026 22:04:55 +0100 Subject: [PATCH 2/7] Update downstream.yml --- .github/workflows/downstream.yml | 1 + 1 file changed, 1 insertion(+) diff --git a/.github/workflows/downstream.yml b/.github/workflows/downstream.yml index 9b71a64d..cd30f7d2 100644 --- a/.github/workflows/downstream.yml +++ b/.github/workflows/downstream.yml @@ -40,6 +40,7 @@ jobs: os: [ubuntu-latest] package: - {repo: BlockBandedMatrices.jl, group: JuliaLinearAlgebra} + - {repo: LazyArrays.jl, group: JuliaArrays} - {repo: LazyBandedMatrices.jl, group: JuliaLinearAlgebra} - {repo: ApproxFunBase.jl, group: JuliaApproximation} - {repo: ApproxFun.jl, group: JuliaApproximation} From 875bf0cf5a5a39f9013e0744f9e3165db300fb7f Mon Sep 17 00:00:00 2001 From: Sheehan Olver Date: Tue, 13 Oct 2026 08:12:34 +0100 Subject: [PATCH 3/7] Update interfaceimpl.jl --- src/interfaceimpl.jl | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/src/interfaceimpl.jl b/src/interfaceimpl.jl index 7ab474b7..bd357039 100644 --- a/src/interfaceimpl.jl +++ b/src/interfaceimpl.jl @@ -36,7 +36,15 @@ bandedrowsdata(Ac::Adjoint{<:Real}) = Transpose(bandeddata(parent(Ac))) bandedrowsdata(V::SubArray{<:Number}) = Transpose(bandeddata(view(parent(parent(V)), reverse(parentindices(V))...))) # OffDiagonal, and views of Diagonal/OffDiagonal by unit ranges, are supported on a single band -bandwidths(::OffDiagonalLayout, A) = (-offdiagonalindex(A), offdiagonalindex(A)) +function bandwidths(::OffDiagonalLayout, A) + b = offdiagonalindex(A) + m,n = size(A) + if -b ≥ m || b ≥ n + (-720, -720) # empty + else + (-b, b) + end +end inbands_getindex(A::OffDiagonal, k::Integer, j::Integer) = A.band[min(k,j)] inbands_setindex!(A::OffDiagonal, v, k::Integer, j::Integer) = (A.band[min(k,j)] = v) From 9c51978c3cc607414550587e6040d2f7a4c6200c Mon Sep 17 00:00:00 2001 From: Sheehan Olver Date: Tue, 13 Oct 2026 08:50:20 +0100 Subject: [PATCH 4/7] Update BandedMatrix.jl --- src/banded/BandedMatrix.jl | 7 +------ 1 file changed, 1 insertion(+), 6 deletions(-) diff --git a/src/banded/BandedMatrix.jl b/src/banded/BandedMatrix.jl index e3adbbbf..0d7e6bc7 100644 --- a/src/banded/BandedMatrix.jl +++ b/src/banded/BandedMatrix.jl @@ -296,13 +296,8 @@ function _BandedMatrix(::DiagonalLayout, A::AbstractMatrix{T}) where T _BandedMatrix(dat, m, 0, 0) end function _BandedMatrix(::OffDiagonalLayout, A::AbstractMatrix{T}) where T - m,n = size(A) k = offdiagonalindex(A) - dat = zeros(T, 1, n) - # the entry in column j lies in row j-k - band = offdiagonaldata(A) - copyto!(view(dat, 1, max(0,k) .+ axes(band,1)), band) - _BandedMatrix(dat, m, -k, k) + _BandedMatrix(copy(permutedims(offdiagonaldata(A))), axes(A,1), -k, k) end function _BandedMatrix(::BidiagonalLayout, A::AbstractMatrix{T}) where T m,n = size(A) From ef938dce482564118e3226d1ee89a76a33392b30 Mon Sep 17 00:00:00 2001 From: Sheehan Olver Date: Tue, 13 Oct 2026 09:05:27 +0100 Subject: [PATCH 5/7] Update BandedMatrix.jl --- src/banded/BandedMatrix.jl | 39 +++++++++++++++++++------------------- 1 file changed, 20 insertions(+), 19 deletions(-) diff --git a/src/banded/BandedMatrix.jl b/src/banded/BandedMatrix.jl index 0d7e6bc7..9db33589 100644 --- a/src/banded/BandedMatrix.jl +++ b/src/banded/BandedMatrix.jl @@ -290,37 +290,38 @@ function _BandedMatrix(::BandedRows, A::AbstractMatrix) return _BandedMatrix(data_new, axes(A, 1), u, ℓ) end function _BandedMatrix(::DiagonalLayout, A::AbstractMatrix{T}) where T - m,n = size(A) - dat = Matrix{T}(undef, 1, n) - copyto!(vec(dat), diagonaldata(A)) - _BandedMatrix(dat, m, 0, 0) + ret = BandedMatrix{T}(undef, axes(A), (0,0)) + copyto!(ret.data, diagonaldata(A)) + ret end function _BandedMatrix(::OffDiagonalLayout, A::AbstractMatrix{T}) where T k = offdiagonalindex(A) - _BandedMatrix(copy(permutedims(offdiagonaldata(A))), axes(A,1), -k, k) + band = offdiagonaldata(A) + ret = BandedMatrix{T}(undef, axes(A), (-k,k)) + copyto!(view(ret.data, range(max(1,k+1); length=length(band))), band) + ret end function _BandedMatrix(::BidiagonalLayout, A::AbstractMatrix{T}) where T m,n = size(A) - dat = Matrix{T}(undef, 2, n) if bidiagonaluplo(A) == 'L' - copyto!(view(dat, 1, :), diagonaldata(A)) - copyto!(view(dat, 2, 1:n-1), subdiagonaldata(A)) - _BandedMatrix(dat, m, 1, 0) + ret = BandedMatrix{T}(undef, axes(A), (1,0)) + copyto!(view(ret.data, 1, :), diagonaldata(A)) + copyto!(view(ret.data, 2, 1:n-1), subdiagonaldata(A)) + ret else # u - copyto!(view(dat, 1, 2:n), supdiagonaldata(A)) - copyto!(view(dat, 2, :), diagonaldata(A)) - _BandedMatrix(dat, m, 0, 1) + ret = BandedMatrix{T}(undef, axes(A), (0,1)) + copyto!(view(ret.data, 1, 2:n), supdiagonaldata(A)) + copyto!(view(ret.data, 2, :), diagonaldata(A)) + ret end end function _BandedMatrix(::AbstractTridiagonalLayout, A::AbstractMatrix{T}) where T - m,n = size(A) - dat = Matrix{T}(undef, 3, n) - - copyto!(view(dat, 1, 2:n), supdiagonaldata(A)) - copyto!(view(dat, 2, :), diagonaldata(A)) - copyto!(view(dat, 3, 1:n-1), subdiagonaldata(A)) - _BandedMatrix(dat, m, 1, 1) + ret = BandedMatrix{T}(undef, axes(A), (1,1)) + copyto!(view(ret.data, 1, 2:n), supdiagonaldata(A)) + copyto!(view(ret.data, 2, :), diagonaldata(A)) + copyto!(view(ret.data, 3, 1:n-1), subdiagonaldata(A)) + ret end From 56a217460fc76a7d62e2f68f902db71592af5a30 Mon Sep 17 00:00:00 2001 From: Sheehan Olver Date: Tue, 13 Oct 2026 09:56:22 +0100 Subject: [PATCH 6/7] Update BandedMatrix.jl --- src/banded/BandedMatrix.jl | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/src/banded/BandedMatrix.jl b/src/banded/BandedMatrix.jl index 9db33589..fae81976 100644 --- a/src/banded/BandedMatrix.jl +++ b/src/banded/BandedMatrix.jl @@ -318,8 +318,10 @@ end function _BandedMatrix(::AbstractTridiagonalLayout, A::AbstractMatrix{T}) where T ret = BandedMatrix{T}(undef, axes(A), (1,1)) + d = diagonaldata(A) + n = length(d) copyto!(view(ret.data, 1, 2:n), supdiagonaldata(A)) - copyto!(view(ret.data, 2, :), diagonaldata(A)) + copyto!(view(ret.data, 2, :), d) copyto!(view(ret.data, 3, 1:n-1), subdiagonaldata(A)) ret end From 39b996891bf072eb430eb6d84e3f83170d4cd8f0 Mon Sep 17 00:00:00 2001 From: Sheehan Olver Date: Tue, 13 Oct 2026 11:52:35 +0100 Subject: [PATCH 7/7] BandedStyle wins over StructuredMatrixStyle{OffDiagonal} Keeps BandedMatrix .+ OffDiagonal banded once OffDiagonal has its own broadcast style (JuliaArrays/FillArrays.jl#458). Co-Authored-By: Claude Opus 5.5 --- src/interfaceimpl.jl | 2 +- test/test_interface.jl | 12 ++++++++++++ 2 files changed, 13 insertions(+), 1 deletion(-) diff --git a/src/interfaceimpl.jl b/src/interfaceimpl.jl index bd357039..4d4ae066 100644 --- a/src/interfaceimpl.jl +++ b/src/interfaceimpl.jl @@ -4,7 +4,7 @@ # Sparse BroadcastStyle ## -for Typ in (:Diagonal, :SymTridiagonal, :Tridiagonal, :Bidiagonal) +for Typ in (:Diagonal, :SymTridiagonal, :Tridiagonal, :Bidiagonal, :OffDiagonal) @eval begin BroadcastStyle(::StructuredMatrixStyle{<:$Typ}, ::BandedStyle) = BandedStyle() diff --git a/test/test_interface.jl b/test/test_interface.jl index 30a0d460..6e71e685 100644 --- a/test/test_interface.jl +++ b/test/test_interface.jl @@ -125,6 +125,18 @@ LinearAlgebra.fill!(A::PseudoBandedMatrix, v) = fill!(A.data,v) @test A[band(0)] == zeros(4) @test MemoryLayout(view(A, 1:3, 2:4)) isa OffDiagonalLayout{DenseColumnMajor} @test bandwidths(view(A, 1:3, 2:4)) == (0,0) + + @testset "broadcast with BandedMatrix" begin + A = FillArrays.OffDiagonal([1.,2,3], 1) + B = brand(4, 4, 1, 1) + for C in (A .+ B, B .+ A, A .* B, sin.(A) .+ 2 .* B) + @test C isa BandedMatrix + end + @test A .+ B == Matrix(A) .+ Matrix(B) + @test bandwidths(A .+ B) == (1,1) + @test A .* B == Matrix(A) .* Matrix(B) + @test sin.(A) .+ 2 .* B ≈ sin.(Matrix(A)) .+ 2 .* Matrix(B) + end end @testset "SymTridiagonal" begin