From ad9809561c100346e7045aa25d812d648bbf1cac Mon Sep 17 00:00:00 2001 From: Daniel VandenHeuvel <95613936+DanielVandH@users.noreply.github.com> Date: Wed, 19 Nov 2025 22:24:41 +0000 Subject: [PATCH 1/7] Use a copyto! approach, and switch to PermutedDimsArray --- Project.toml | 2 +- src/banded/BandedMatrix.jl | 21 ++++++++++++++------- src/interfaceimpl.jl | 6 +++--- test/test_banded.jl | 3 +++ 4 files changed, 21 insertions(+), 11 deletions(-) diff --git a/Project.toml b/Project.toml index 3a125161..d3dc1960 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "BandedMatrices" uuid = "aae01518-5342-5314-be14-df237901396f" -version = "1.10.1" +version = "1.10.2" [deps] ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a" diff --git a/src/banded/BandedMatrix.jl b/src/banded/BandedMatrix.jl index 25b72a05..8ac5894e 100644 --- a/src/banded/BandedMatrix.jl +++ b/src/banded/BandedMatrix.jl @@ -266,14 +266,21 @@ BandedMatrix(A::AbstractMatrix) = _BandedMatrix(MemoryLayout(A), A) _BandedMatrix(::BandedColumns, A::AbstractMatrix) = _BandedMatrix(copy(bandeddata(A)), axes(A,1), bandwidths(A)...) function _BandedMatrix(::BandedRows, A::AbstractMatrix) bdata = bandedrowsdata(A) - m,ν = size(bdata) - ℓ,u = bandwidths(A) - data = similar(bdata, eltype(bdata), ν, size(A,2)) - n = size(A, 2) - for j in axes(A, 1), i in max(1, j - ℓ):min(n, j + u) - data[u + 1 + j - i, i] = bdata[j, ℓ+1+i-j] + u, ℓ = bandwidths(A) + m, n = size(A) + data_new = similar(bdata, eltype(bdata), size(bdata, 2), size(A, 2)) + for i in 1:(ℓ+u+1) + d = u + 1 - i + row_new = ℓ + 1 + d + k_start = max(1, 1 - d) + k_end = min(m, n - d) + if k_start ≤ k_end + col_start_A = k_start + d + col_end_A = k_end + d + copyto!(view(data_new, row_new, k_start:k_end), view(bdata, col_start_A:col_end_A, i)) + end end - return _BandedMatrix(data, axes(A, 1), ℓ, u) + return _BandedMatrix(data_new, axes(A, 2), u, ℓ) end function _BandedMatrix(::DiagonalLayout, A::AbstractMatrix{T}) where T m,n = size(A) diff --git a/src/interfaceimpl.jl b/src/interfaceimpl.jl index 6f54d624..efdfc04a 100644 --- a/src/interfaceimpl.jl +++ b/src/interfaceimpl.jl @@ -30,9 +30,9 @@ function bandeddata(D::Union{Diagonal, Transpose{<:Any, <:Diagonal}, Adjoint{<:A permutedims(diagonaldata(D)) end -bandedrowsdata(Ac::Transpose) = permutedims(bandeddata(parent(Ac))) -bandedrowsdata(Ac::Adjoint{<:Real}) = permutedims(bandeddata(parent(Ac))) -bandedrowsdata(V::SubArray) = permutedims(bandeddata(view(parent(parent(V)), reverse(parentindices(V))...))) +bandedrowsdata(Ac::Transpose) = PermutedDimsArray(bandeddata(parent(Ac)), (2, 1)) +bandedrowsdata(Ac::Adjoint{<:Real}) = PermutedDimsArray(bandeddata(parent(Ac)), (2, 1)) +bandedrowsdata(V::SubArray) = PermutedDimsArray(bandeddata(view(parent(parent(V)), reverse(parentindices(V))...)), (2, 1)) # treat subinds as banded sublayout(::DiagonalLayout{L}, inds::Type) where L = sublayout(bandedcolumns(L()), inds) diff --git a/test/test_banded.jl b/test/test_banded.jl index 84b259b2..e98596a3 100644 --- a/test/test_banded.jl +++ b/test/test_banded.jl @@ -620,6 +620,9 @@ include("mymatrix.jl") @test BandedMatrix(V) == V end + @test BandedMatrices.bandedrowsdata(A') isa PermutedDimsArray + @test BandedMatrices.bandedrowsdata(transpose(A)) isa PermutedDimsArray + @test BandedMatrices.bandedrowsdata(view(A', 1:2, 1:2)) isa PermutedDimsArray end end From 874ddb3f81f50d27ead818f4743f1bd9a443ea8f Mon Sep 17 00:00:00 2001 From: Daniel VandenHeuvel <95613936+DanielVandH@users.noreply.github.com> Date: Wed, 19 Nov 2025 22:43:43 +0000 Subject: [PATCH 2/7] Fix for rectangular matrices --- src/banded/BandedMatrix.jl | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/src/banded/BandedMatrix.jl b/src/banded/BandedMatrix.jl index 8ac5894e..ac778892 100644 --- a/src/banded/BandedMatrix.jl +++ b/src/banded/BandedMatrix.jl @@ -264,23 +264,23 @@ BandedMatrix(A::AbstractMatrix) = _BandedMatrix(MemoryLayout(A), A) ## specialised # use bandeddata if possible _BandedMatrix(::BandedColumns, A::AbstractMatrix) = _BandedMatrix(copy(bandeddata(A)), axes(A,1), bandwidths(A)...) -function _BandedMatrix(::BandedRows, A::AbstractMatrix) +function _BandedMatrix(::BandedMatrices.BandedRows, A::AbstractMatrix) bdata = bandedrowsdata(A) u, ℓ = bandwidths(A) m, n = size(A) - data_new = similar(bdata, eltype(bdata), size(bdata, 2), size(A, 2)) + data_new = similar(bdata, eltype(bdata), ℓ + u + 1, n) for i in 1:(ℓ+u+1) - d = u + 1 - i + d = u + 1 - i row_new = ℓ + 1 + d - k_start = max(1, 1 - d) - k_end = min(m, n - d) + k_start = max(1, 1 - d) + k_end = min(n, m - d) if k_start ≤ k_end col_start_A = k_start + d col_end_A = k_end + d copyto!(view(data_new, row_new, k_start:k_end), view(bdata, col_start_A:col_end_A, i)) end end - return _BandedMatrix(data_new, axes(A, 2), u, ℓ) + return _BandedMatrix(data_new, axes(A, 1), u, ℓ) end function _BandedMatrix(::DiagonalLayout, A::AbstractMatrix{T}) where T m,n = size(A) From 561d510d3eb015218c1859649d8bff6cc628fa12 Mon Sep 17 00:00:00 2001 From: Daniel VandenHeuvel <95613936+DanielVandH@users.noreply.github.com> Date: Wed, 19 Nov 2025 22:44:15 +0000 Subject: [PATCH 3/7] Remove Module prefix --- src/banded/BandedMatrix.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/banded/BandedMatrix.jl b/src/banded/BandedMatrix.jl index ac778892..7e503d09 100644 --- a/src/banded/BandedMatrix.jl +++ b/src/banded/BandedMatrix.jl @@ -264,7 +264,7 @@ BandedMatrix(A::AbstractMatrix) = _BandedMatrix(MemoryLayout(A), A) ## specialised # use bandeddata if possible _BandedMatrix(::BandedColumns, A::AbstractMatrix) = _BandedMatrix(copy(bandeddata(A)), axes(A,1), bandwidths(A)...) -function _BandedMatrix(::BandedMatrices.BandedRows, A::AbstractMatrix) +function _BandedMatrix(::BandedRows, A::AbstractMatrix) bdata = bandedrowsdata(A) u, ℓ = bandwidths(A) m, n = size(A) From 28337912069844b54eccff7ffa2ef09591c45962 Mon Sep 17 00:00:00 2001 From: Daniel VandenHeuvel <95613936+DanielVandH@users.noreply.github.com> Date: Wed, 19 Nov 2025 23:38:23 +0000 Subject: [PATCH 4/7] Use size instead of sum --- src/banded/BandedMatrix.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/banded/BandedMatrix.jl b/src/banded/BandedMatrix.jl index 7e503d09..58059b38 100644 --- a/src/banded/BandedMatrix.jl +++ b/src/banded/BandedMatrix.jl @@ -268,7 +268,7 @@ function _BandedMatrix(::BandedRows, A::AbstractMatrix) bdata = bandedrowsdata(A) u, ℓ = bandwidths(A) m, n = size(A) - data_new = similar(bdata, eltype(bdata), ℓ + u + 1, n) + data_new = similar(bdata, eltype(bdata), size(bdata, 2), n) for i in 1:(ℓ+u+1) d = u + 1 - i row_new = ℓ + 1 + d From acbad0057d799f9799a760c7f05b872fd5fa4f98 Mon Sep 17 00:00:00 2001 From: Daniel VandenHeuvel <95613936+DanielVandH@users.noreply.github.com> Date: Thu, 20 Nov 2025 00:09:51 +0000 Subject: [PATCH 5/7] Use axes --- src/banded/BandedMatrix.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/banded/BandedMatrix.jl b/src/banded/BandedMatrix.jl index 58059b38..85f80a0b 100644 --- a/src/banded/BandedMatrix.jl +++ b/src/banded/BandedMatrix.jl @@ -269,7 +269,7 @@ function _BandedMatrix(::BandedRows, A::AbstractMatrix) u, ℓ = bandwidths(A) m, n = size(A) data_new = similar(bdata, eltype(bdata), size(bdata, 2), n) - for i in 1:(ℓ+u+1) + for i in axes(bdata, 2) d = u + 1 - i row_new = ℓ + 1 + d k_start = max(1, 1 - d) From ec8d28af04b7d07143319ff65d99ffb28b41ad7f Mon Sep 17 00:00:00 2001 From: Sheehan Olver Date: Thu, 20 Nov 2025 10:58:03 +0000 Subject: [PATCH 6/7] Use Transpose for now --- src/interfaceimpl.jl | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/src/interfaceimpl.jl b/src/interfaceimpl.jl index efdfc04a..3a3a7b15 100644 --- a/src/interfaceimpl.jl +++ b/src/interfaceimpl.jl @@ -30,9 +30,10 @@ function bandeddata(D::Union{Diagonal, Transpose{<:Any, <:Diagonal}, Adjoint{<:A permutedims(diagonaldata(D)) end -bandedrowsdata(Ac::Transpose) = PermutedDimsArray(bandeddata(parent(Ac)), (2, 1)) -bandedrowsdata(Ac::Adjoint{<:Real}) = PermutedDimsArray(bandeddata(parent(Ac)), (2, 1)) -bandedrowsdata(V::SubArray) = PermutedDimsArray(bandeddata(view(parent(parent(V)), reverse(parentindices(V))...)), (2, 1)) +#TODO: PermutedDimsArray(bandeddata(parent(Ac)), (2, 1)) +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) From 57d63c94d9ae7566afa25af82fd5210680fabba3 Mon Sep 17 00:00:00 2001 From: Daniel VandenHeuvel <95613936+DanielVandH@users.noreply.github.com> Date: Thu, 20 Nov 2025 11:54:05 +0000 Subject: [PATCH 7/7] Change test --- test/test_banded.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/test/test_banded.jl b/test/test_banded.jl index e98596a3..b8789d10 100644 --- a/test/test_banded.jl +++ b/test/test_banded.jl @@ -620,9 +620,9 @@ include("mymatrix.jl") @test BandedMatrix(V) == V end - @test BandedMatrices.bandedrowsdata(A') isa PermutedDimsArray - @test BandedMatrices.bandedrowsdata(transpose(A)) isa PermutedDimsArray - @test BandedMatrices.bandedrowsdata(view(A', 1:2, 1:2)) isa PermutedDimsArray + @test BandedMatrices.bandedrowsdata(A') isa Transpose + @test BandedMatrices.bandedrowsdata(transpose(A)) isa Transpose + @test BandedMatrices.bandedrowsdata(view(A', 1:2, 1:2)) isa Transpose end end