Fixed BlockDiagonal

Author: jagot

src/BlockBandedMatrices.jl

@@ -45,7 +45,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, BlockBidiagonal, isblockbanded
+        subblockbandwidth, subblockbandwidths, Ones, Zeros, Fill, Block, BlockDiagonal, BlockTridiagonal, BlockBidiagonal, isblockbanded
 
 
 const Block1 = Block{1,Int}

src/interfaceimpl.jl

@@ -25,15 +25,15 @@ 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
     (size.(A.diag, 1), size.(A.diag,2))
-end    
+end
 
 
 # Block Bi/Tridiagonal
@@ -43,7 +43,7 @@ const BlockBidiagonal{T,VT<:Matrix{T}} = BlockMatrix{T,<:Bidiagonal{VT}}
 BlockTridiagonal(A,B,C) = mortar(Tridiagonal(A,B,C))
 BlockBidiagonal(A, B, uplo) = mortar(Bidiagonal(A,B,uplo))
 
-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"))
@@ -53,7 +53,7 @@ function sizes_from_blocks(A::Tridiagonal, _)
     (size.(A.d, 1), size.(A.d,2))
 end
 
-function sizes_from_blocks(A::Bidiagonal, _) 
+function sizes_from_blocks(A::Bidiagonal, _)
     # 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"))
@@ -83,21 +83,44 @@ 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(broadcast($op,A.blocks.dl), broadcast($op, Ref(λ), A.blocks.d), broadcast($op,A.blocks.du)))
         end
-        function $op(A::BlockBidiagonal, λ::UniformScaling) 
+        function $op(A::BlockBidiagonal, λ::UniformScaling)
             checksquareblocks(A)
             mortar(Bidiagonal(broadcast($op, A.blocks.dv, Ref(λ)), A.blocks.ev, A.blocks.uplo))
         end
-        function $op(λ::UniformScaling, A::BlockBidiagonal) 
+        function $op(λ::UniformScaling, A::BlockBidiagonal)
             checksquareblocks(A)
             mortar(Bidiagonal(broadcast($op, Ref(λ), A.blocks.dv), broadcast($op,A.blocks.ev), A.blocks.uplo))
         end
     end
-end
+end
+
+function replace_in_print_matrix(A::BlockDiagonal, i::Integer, j::Integer, s::AbstractString)
+    bi = findblockindex.(axes(A), (i,j))
+    I,J = block.(bi)
+    i,j = blockindex.(bi)
+    Int(J-I) == 0 ? s : Base.replace_with_centered_mark(s)
+end
+
+function replace_in_print_matrix(A::BlockTridiagonal, i::Integer, j::Integer, s::AbstractString)
+    bi = findblockindex.(axes(A), (i,j))
+    I,J = block.(bi)
+    i,j = blockindex.(bi)
+    -1 ≤ Int(J-I) ≤ 1 ? s : Base.replace_with_centered_mark(s)
+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 Bidiagonal" begin
     Bu = BlockBidiagonal(fill([1 2],4), fill([3 4],3), :U)
     Bl = BlockBidiagonal(fill([1 2],4), fill([3 4],3), :L)
@@ -52,4 +63,4 @@ end
     @test A+I == I+A == mortar(Tridiagonal(fill([1 2; 1 2],3), fill([4 4; 3 5],4), fill([4 5; 4 5],3))) == Matrix(A) + I
     @test A-I == mortar(Tridiagonal(fill([1 2; 1 2],3), fill([2 4; 3 3],4), fill([4 5; 4 5],3))) == Matrix(A) - I
     @test I-A == mortar(Tridiagonal(fill(-[1 2; 1 2],3), fill([-2 -4; -3 -3],4), fill(-[4 5; 4 5],3))) == I - Matrix(A)
-end
+end