Overload _BandedMatrix(::BandedRows, A::AbstactMatrix) (#486)

* BandedRows

* EOL

* Add new interface function

* Use permutedims; cleanup implementation a bit

* Combine iterator

* don't test on pre anymore

* Test transpose as well

---------

Co-authored-by: Sheehan Olver <solver@mac.com>

Author: DanielVandH

.github/workflows/ci.yml

@@ -45,7 +45,6 @@ jobs:
         version:
           - 'lts'
           - '1'
-          - 'pre'
         os:
           - ubuntu-latest
           - macOS-latest

Project.toml

@@ -1,6 +1,6 @@
 name = "BandedMatrices"
 uuid = "aae01518-5342-5314-be14-df237901396f"
-version = "1.9.5"
+version = "1.10.0"
 
 [deps]
 ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"

docs/src/index.md

@@ -218,6 +218,7 @@ interface consists of the following:
 | `inbands_setindex!(A, v, k, j)` | Unsafe: set `A[k,j] = v`, without the need to check if we are inside the bands |
 | `BandedMatrices.MemoryLayout(A)` | Override to get banded lazy linear algebra, e.g. `y .= Mul(A,x)` |
 | `BandedMatrices.bandeddata(A)` | Override to return a matrix of the entries in BLAS format. Required if `MemoryLayout(A)` returns `BandedColumnMajor` |
+| `BandedMatrices.bandedrowsdata(A)` | Override to return a matrix of the entries reshaped to a row-major format. Required if `MemoryLayout(A)` returns `BandedRowsMajor` |
 
 Note that certain `SubArray`s of `BandedMatrix` are also banded matrices.
 The banded matrix interface is implemented for such `SubArray`s to take advantage of this.

src/banded/BandedMatrix.jl

@@ -264,6 +264,16 @@ 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)
+    bdata = bandedrowsdata(A)
+    data = similar(bdata, eltype(bdata), reverse(size(bdata)))
+    u, ℓ = bandwidths(A)
+    n = size(A, 2)
+    for j in axes(A, 1), i in max(1, j - u):min(n, j + ℓ)
+        data[ℓ + 1 + j - i, i] = bdata[j, u+1+i-j]
+    end
+    return _BandedMatrix(data, axes(A, 1), bandwidths(A)...)
+end
 function _BandedMatrix(::DiagonalLayout, A::AbstractMatrix{T}) where T
     m,n = size(A)
     dat = Matrix{T}(undef, 1, n)

src/interfaceimpl.jl

@@ -30,6 +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)))
+
 # 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)

test/test_banded.jl

@@ -600,6 +600,16 @@ include("mymatrix.jl")
         @test B2 == B
         @test typeof(B2) == typeof(B)
     end
+
+    @testset "_BandedMatrix(::BandedRows)" begin
+        A = brand(5, 5, 2, 1)
+        for B in (A', transpose(A))
+            @test BandedMatrices._BandedMatrix(MemoryLayout(B), B) == B
+            @test BandedMatrix(B) == B
+            @test BandedMatrix(B) !== B
+            @test BandedMatrices.bandeddata(BandedMatrix(B)) !== BandedMatrices.bandeddata(A)
+        end
+    end
 end
 
 end # module