Fixed BlockDiagonal

Author: jagot

src/BlockBandedMatrices.jl

@@ -42,7 +42,7 @@ import BandedMatrices: isbanded, bandwidths, bandwidth, banded_getindex, colrang
                         _banded_colval, _banded_rowval, _banded_nzval # for sparse
 
 export BandedBlockBandedMatrix, BlockBandedMatrix, BlockSkylineMatrix, blockbandwidth, blockbandwidths,
-        subblockbandwidth, subblockbandwidths, Ones, Zeros, Fill, Block, BlockTridiagonal, isblockbanded
+        subblockbandwidth, subblockbandwidths, Ones, Zeros, Fill, Block, BlockDiagonal, BlockTridiagonal, isblockbanded
 
 
 const Block1 = Block{1,Int}

src/interfaceimpl.jl

@@ -25,23 +25,23 @@ const BlockDiagonal{T,VT<:Matrix{T}} = BlockMatrix{T,<:Diagonal{VT}}
 
 BlockDiagonal(A) = mortar(Diagonal(A))
 
-function sizes_from_blocks(A::Diagonal, _) 
+function sizes_from_blocks(A::Diagonal, _)
     # for k = 1:length(A.du)
     #     size(A.du[k],1) == sz[1][k] || throw(ArgumentError("block sizes of upper diagonal inconsisent with diagonal"))
     #     size(A.du[k],2) == sz[2][k+1] || throw(ArgumentError("block sizes of upper diagonal inconsisent with diagonal"))
     #     size(A.dl[k],1) == sz[1][k+1] || throw(ArgumentError("block sizes of lower diagonal inconsisent with diagonal"))
     #     size(A.dl[k],2) == sz[2][k] || throw(ArgumentError("block sizes of lower diagonal inconsisent with diagonal"))
     # end
-    BlockSizes(size.(A.diag, 1), size.(A.diag,2))
-end    
+    (size.(A.diag, 1), size.(A.diag,2))
+end
 
 
 # Block Tridiagonal
 const BlockTridiagonal{T,VT<:Matrix{T}} = BlockMatrix{T,<:Tridiagonal{VT}}
 
 BlockTridiagonal(A,B,C) = mortar(Tridiagonal(A,B,C))
 
-function sizes_from_blocks(A::Tridiagonal, _) 
+function sizes_from_blocks(A::Tridiagonal, _)
     # for k = 1:length(A.du)
     #     size(A.du[k],1) == sz[1][k] || throw(ArgumentError("block sizes of upper diagonal inconsisent with diagonal"))
     #     size(A.du[k],2) == sz[2][k+1] || throw(ArgumentError("block sizes of upper diagonal inconsisent with diagonal"))
@@ -64,11 +64,20 @@ checksquareblocks(A) = blockisequal(axes(A)...) || throw(DimensionMismatch("bloc
 
 for op in (:-, :+)
     @eval begin
-        function $op(A::BlockTridiagonal, λ::UniformScaling) 
+        function $op(A::BlockDiagonal, λ::UniformScaling)
+            checksquareblocks(A)
+            mortar(Diagonal(broadcast($op, A.blocks.diag, Ref(λ))))
+        end
+        function $op(λ::UniformScaling, A::BlockDiagonal)
+            checksquareblocks(A)
+            mortar(Diagonal(broadcast($op, Ref(λ), A.blocks.diag)))
+        end
+
+        function $op(A::BlockTridiagonal, λ::UniformScaling)
             checksquareblocks(A)
             mortar(Tridiagonal(A.blocks.dl, broadcast($op, A.blocks.d, Ref(λ)), A.blocks.du))
         end
-        function $op(λ::UniformScaling, A::BlockTridiagonal) 
+        function $op(λ::UniformScaling, A::BlockTridiagonal)
             checksquareblocks(A)
             mortar(Tridiagonal(A.blocks.dl, broadcast($op, Ref(λ), A.blocks.d), A.blocks.du))
         end
@@ -87,4 +96,4 @@ function replace_in_print_matrix(A::BlockTridiagonal, i::Integer, j::Integer, s:
     I,J = block.(bi)
     i,j = blockindex.(bi)
     -1 ≤ Int(J-I) ≤ 1 ? s : Base.replace_with_centered_mark(s)
-end
+end

test/test_misc.jl

@@ -21,6 +21,17 @@ import BandedMatrices: bandwidths, AbstractBandedMatrix, BandedStyle, bandeddata
     @test BandedMatrix(V) == V
 end
 
+@testset "Block Diagonal" begin
+    A = BlockDiagonal(fill([1 2],3))
+    @test blockbandwidths(A) == (0,0)
+    @test isblockbanded(A)
+    @test A[Block(1,1)] == [1 2]
+    @test @inferred(getblock(A,1,2)) == @inferred(A[Block(1,2)]) == [0 0]
+    @test_throws DimensionMismatch A+I
+    A = BlockDiagonal(fill([1 2; 1 2],3))
+    @test A+I == I+A == mortar(Diagonal(fill([2 2; 1 3],3))) == Matrix(A) + I
+end
+
 @testset "Block Tridiagonal" begin
     A = BlockTridiagonal(fill([1 2],3), fill([3 4],4), fill([4 5],3))
     @test blockbandwidths(A) == (1,1)