Reorganize code and remove unused _BandedBlockMatrix fn (#189)

* Reorganize code and remove unused _BandedBlockMatrix fn

* Add tests

* Update BlockSkyLine docstring example

* Adjust comma

* fix type namespace in docstring

Author: jishnub

src/BlockSkylineMatrix.jl

@@ -103,75 +103,40 @@ function bb_numentries(B::BlockSkylineSizes)
     numentries
 end
 
-
-
-
-function _BandedBlockMatrix end
-
-
-#  A block matrix where only the bands are nonzero
-#   isomorphic to BandedMatrix{Matrix{T}}
-struct BlockSkylineMatrix{T, DATA<:AbstractVector{T}, BS<:BlockSkylineSizes} <: AbstractBlockBandedMatrix{T}
-    data::DATA
-    block_sizes::BS
-
-    global function _BlockSkylineMatrix(data::DATA, block_sizes::BS) where {T,DATA<:AbstractVector{T}, BS<:BlockSkylineSizes}
-        new{T,DATA,BS}(data, block_sizes)
-    end
-end
-
-"""
-    BlockBandedMatrix
-
-A `BlockBandedMatrix` is a subtype of `BlockMatrix` of `BlockArrays.jl` whose
-layout of non-zero blocks is banded.
-"""
-const BlockBandedMatrix{T} = BlockSkylineMatrix{T, Vector{T}, BlockBandedSizes}
-
-# Auxiliary outer constructors
-@inline _BlockBandedMatrix(data::AbstractVector, bs::BlockBandedSizes) =
-    _BlockSkylineMatrix(data, bs)
-
-@inline _BlockBandedMatrix(data::AbstractVector, kr::AbstractVector{Int}, jr::AbstractVector{Int}, (l,u)::NTuple{2, Int}) =
-    _BlockBandedMatrix(data, BlockBandedSizes(kr,jr, l,u))
-
-@inline BlockSkylineMatrix{T}(::UndefInitializer, block_sizes::BlockSkylineSizes) where T =
-    _BlockSkylineMatrix(Vector{T}(undef, bb_numentries(block_sizes)), block_sizes)
-
-@inline BlockBandedMatrix{T}(::UndefInitializer, block_sizes::BlockBandedSizes) where T =
-    _BlockSkylineMatrix(Vector{T}(undef, bb_numentries(block_sizes)), block_sizes)
-
 """
     BlockSkylineMatrix{T,LL,UU}(M::Union{UndefInitializer,UniformScaling,AbstractMatrix},
                                 rows, cols, (l::LL, u::UU))
 
 returns a `sum(rows)`×`sum(cols)` block-banded matrix `A` having elements of type `T`,
-with block-bandwidths `(l,u)` and where `A[Block(K,J)]` is a
+with block-bandwidths `(l,u)`, and where `A[Block(K,J)]` is a
 `Matrix{T}` of size `rows[K]`×`cols[J]`.
 
-`(l,u)` may be integers for constant bandwidths or integer vectors of
-lengths `rows` and `cols`, respectively, for ragged bands.
+`(l,u)` may be integers for constant bandwidths, or integer vectors of length
+`length(cols)` for ragged bands. In the latter case, `l` and `u` represent the
+number of sub and super-block-bands in each column.
 
 # Examples
 
 ```jldoctest
 julia> using LinearAlgebra, FillArrays
 
-julia> BlockSkylineMatrix(I, [2,3,4], [1,2,3], ([2,0,1],[1,1,1]))
-3×3-blocked 9×6 BlockSkylineMatrix{Bool, Vector{Bool}, BlockBandedMatrices.BlockSkylineSizes{Tuple{BlockArrays.BlockedUnitRange{Vector{Int64}}, BlockArrays.BlockedUnitRange{Vector{Int64}}}, Vector{Int64}, Vector{Int64}, BandedMatrices.BandedMatrix{Int64, Matrix{Int64}, Base.OneTo{Int64}}, Vector{Int64}}}:
+julia> BlockSkylineMatrix(I, [2,2,2,4], [1,2,3], ([2,0,1],[0,1,1]))
+4×3-blocked 10×6 BlockSkylineMatrix{Bool, Vector{Bool}, BlockBandedMatrices.BlockSkylineSizes{Tuple{BlockArrays.BlockedUnitRange{Vector{Int64}}, BlockArrays.BlockedUnitRange{Vector{Int64}}}, Vector{Int64}, Vector{Int64}, BandedMatrices.BandedMatrix{Int64, Matrix{Int64}, Base.OneTo{Int64}}, Vector{Int64}}}:
  1  │  0  0  │  ⋅  ⋅  ⋅
  0  │  1  0  │  ⋅  ⋅  ⋅
  ───┼────────┼─────────
  0  │  0  1  │  0  0  0
  0  │  0  0  │  1  0  0
- 0  │  0  0  │  0  1  0
  ───┼────────┼─────────
+ 0  │  ⋅  ⋅  │  0  1  0
  0  │  ⋅  ⋅  │  0  0  1
- 0  │  ⋅  ⋅  │  0  0  0
- 0  │  ⋅  ⋅  │  0  0  0
- 0  │  ⋅  ⋅  │  0  0  0
+ ───┼────────┼─────────
+ ⋅  │  ⋅  ⋅  │  0  0  0
+ ⋅  │  ⋅  ⋅  │  0  0  0
+ ⋅  │  ⋅  ⋅  │  0  0  0
+ ⋅  │  ⋅  ⋅  │  0  0  0
 
-julia> BlockSkylineMatrix(Ones(9,6), [2,3,4], [1,2,3], ([2,0,1],[1,1,1]))
+julia> BlockSkylineMatrix(Ones(9,6), [2,3,4], [1,2,3], ([2,0,0],[0,1,1]))
 3×3-blocked 9×6 BlockSkylineMatrix{Float64, Vector{Float64}, BlockBandedMatrices.BlockSkylineSizes{Tuple{BlockArrays.BlockedUnitRange{Vector{Int64}}, BlockArrays.BlockedUnitRange{Vector{Int64}}}, Vector{Int64}, Vector{Int64}, BandedMatrices.BandedMatrix{Int64, Matrix{Int64}, Base.OneTo{Int64}}, Vector{Int64}}}:
  1.0  │  1.0  1.0  │   ⋅    ⋅    ⋅
  1.0  │  1.0  1.0  │   ⋅    ⋅    ⋅
@@ -188,6 +153,38 @@ julia> BlockSkylineMatrix(Ones(9,6), [2,3,4], [1,2,3], ([2,0,1],[1,1,1]))
 """
 BlockSkylineMatrix
 
+#  A block matrix where only the bands are nonzero
+#   isomorphic to BandedMatrix{Matrix{T}}
+struct BlockSkylineMatrix{T, DATA<:AbstractVector{T}, BS<:BlockSkylineSizes} <: AbstractBlockBandedMatrix{T}
+    data::DATA
+    block_sizes::BS
+
+    global function _BlockSkylineMatrix(data::DATA, block_sizes::BS) where {T,DATA<:AbstractVector{T}, BS<:BlockSkylineSizes}
+        new{T,DATA,BS}(data, block_sizes)
+    end
+end
+
+"""
+    BlockBandedMatrix
+
+A `BlockBandedMatrix` is a subtype of `BlockMatrix` of `BlockArrays.jl` whose
+layout of non-zero blocks is banded.
+"""
+const BlockBandedMatrix{T} = BlockSkylineMatrix{T, Vector{T}, BlockBandedSizes}
+
+# Auxiliary outer constructors
+@inline _BlockBandedMatrix(data::AbstractVector, bs::BlockBandedSizes) =
+    _BlockSkylineMatrix(data, bs)
+
+@inline _BlockBandedMatrix(data::AbstractVector, kr::AbstractVector{Int}, jr::AbstractVector{Int}, (l,u)::NTuple{2, Int}) =
+    _BlockBandedMatrix(data, BlockBandedSizes(kr,jr, l,u))
+
+@inline BlockSkylineMatrix{T}(::UndefInitializer, block_sizes::BlockSkylineSizes) where T =
+    _BlockSkylineMatrix(Vector{T}(undef, bb_numentries(block_sizes)), block_sizes)
+
+@inline BlockBandedMatrix{T}(::UndefInitializer, block_sizes::BlockBandedSizes) where T =
+    _BlockSkylineMatrix(Vector{T}(undef, bb_numentries(block_sizes)), block_sizes)
+
 @inline BlockBandedMatrix{T}(::UndefInitializer, axes::NTuple{2,AbstractUnitRange{Int}}, lu::NTuple{2, Int}) where T =
     BlockSkylineMatrix{T}(undef, BlockBandedSizes(axes, lu...))
 

test/test_blockbanded.jl

@@ -35,6 +35,19 @@ import BlockBandedMatrices: MemoryLayout, ColumnMajor, BroadcastStyle,
 
         @test Matrix(BlockBandedMatrix{Int}(I, rows,cols, (l,u))) ==
             Matrix{Int}(I, 10, 10)
+
+        BS = BlockBandedMatrices.BlockBandedSizes(1:3, 1:3, 1, 1)
+        B = BlockBandedMatrix{Float64}(undef, BS)
+        @test blockbandwidths(B) == (1,1)
+        blks = blocks(B)
+        for (ind, blk) in zip(CartesianIndices(blks), blks)
+            @test size(blk) == Tuple(ind)
+        end
+
+        B2 = BlockBandedMatrices._BlockBandedMatrix(1:30, 1:3, 1:3, (1, 1))
+        @test @view(B2[Block(1,1)]) == reshape(1:1, 1, 1)
+        @test @view(B2[Block(2,2)]) == hcat(5:6, 11:12)
+        @test @view(B2[Block(3,3)]) == hcat(18:20, 23:25, 28:30)
     end
 
     @testset "BlockBandedMatrix block indexing" begin