Merge pull request #17685 from Sacha0/specialcat

Fix #17675 (concatenations involving special matrices) and #17738 (concatenations involving sparse vectors)

base/sparse/sparsematrix.jl

@@ -3264,37 +3264,6 @@ function hcat(X::SparseMatrixCSC...)
     SparseMatrixCSC(m, n, colptr, rowval, nzval)
 end
 
-
-# Sparse/dense concatenation
-
-function hcat(Xin::Union{Vector, Matrix, SparseMatrixCSC}...)
-    X = SparseMatrixCSC[issparse(x) ? x : sparse(x) for x in Xin]
-    hcat(X...)
-end
-
-function vcat(Xin::Union{Vector, Matrix, SparseMatrixCSC}...)
-    X = SparseMatrixCSC[issparse(x) ? x : sparse(x) for x in Xin]
-    vcat(X...)
-end
-
-function hvcat(rows::Tuple{Vararg{Int}}, X::Union{Vector, Matrix, SparseMatrixCSC}...)
-    nbr = length(rows)  # number of block rows
-
-    tmp_rows = Array{SparseMatrixCSC}(nbr)
-    k = 0
-    @inbounds for i = 1 : nbr
-        tmp_rows[i] = hcat(X[(1 : rows[i]) + k]...)
-        k += rows[i]
-    end
-    vcat(tmp_rows...)
-end
-
-function cat(catdims, Xin::Union{Vector, Matrix, SparseMatrixCSC}...)
-    X = SparseMatrixCSC[issparse(x) ? x : sparse(x) for x in Xin]
-    T = promote_eltype(Xin...)
-    Base.cat_t(catdims, T, X...)
-end
-
 """
     blkdiag(A...)
 

base/sparse/sparsevector.jl

@@ -720,7 +720,13 @@ complex(x::AbstractSparseVector) =
 
 ### Concatenation
 
-function hcat{Tv,Ti}(X::AbstractSparseVector{Tv,Ti}...)
+# Without the first of these methods, horizontal concatenations of SparseVectors fall
+# back to the horizontal concatenation method that ensures that combinations of
+# sparse/special/dense matrix/vector types concatenate to SparseMatrixCSCs, instead
+# of _absspvec_hcat below. The <:Integer qualifications are necessary for correct dispatch.
+hcat{Tv,Ti<:Integer}(X::SparseVector{Tv,Ti}...) = _absspvec_hcat(X...)
+hcat{Tv,Ti<:Integer}(X::AbstractSparseVector{Tv,Ti}...) = _absspvec_hcat(X...)
+function _absspvec_hcat{Tv,Ti}(X::AbstractSparseVector{Tv,Ti}...)
     # check sizes
     n = length(X)
     m = length(X[1])
@@ -749,7 +755,13 @@ function hcat{Tv,Ti}(X::AbstractSparseVector{Tv,Ti}...)
     SparseMatrixCSC{Tv,Ti}(m, n, colptr, nzrow, nzval)
 end
 
-function vcat{Tv,Ti}(X::AbstractSparseVector{Tv,Ti}...)
+# Without the first of these methods, vertical concatenations of SparseVectors fall
+# back to the vertical concatenation method that ensures that combinations of
+# sparse/special/dense matrix/vector types concatenate to SparseMatrixCSCs, instead
+# of _absspvec_vcat below. The <:Integer qualifications are necessary for correct dispatch.
+vcat{Tv,Ti<:Integer}(X::SparseVector{Tv,Ti}...) = _absspvec_vcat(X...)
+vcat{Tv,Ti<:Integer}(X::AbstractSparseVector{Tv,Ti}...) = _absspvec_vcat(X...)
+function _absspvec_vcat{Tv,Ti}(X::AbstractSparseVector{Tv,Ti}...)
     # check sizes
     n = length(X)
     tnnz = 0
@@ -777,11 +789,49 @@ function vcat{Tv,Ti}(X::AbstractSparseVector{Tv,Ti}...)
     SparseVector(len, rnzind, rnzval)
 end
 
-hcat(Xin::Union{AbstractSparseVector, SparseMatrixCSC}...) = hcat(map(SparseMatrixCSC, Xin)...)
-vcat(Xin::Union{AbstractSparseVector, SparseMatrixCSC}...) = vcat(map(SparseMatrixCSC, Xin)...)
 hcat(Xin::Union{Vector, AbstractSparseVector}...) = hcat(map(sparse, Xin)...)
 vcat(Xin::Union{Vector, AbstractSparseVector}...) = vcat(map(sparse, Xin)...)
 
+
+### Sparse/special/dense vector/matrix concatenation
+
+# TODO: These methods should be moved to a more appropriate location, particularly some
+# future equivalent of base/linalg/special.jl dedicated to interactions between a broader
+# set of matrix types.
+
+# TODO: A similar definition also exists in base/linalg/bidiag.jl. These definitions should
+# be consolidated in a more appropriate location, for example base/linalg/special.jl.
+SpecialArrays = Union{Diagonal, Bidiagonal, Tridiagonal, SymTridiagonal}
+
+function hcat(Xin::Union{Vector, Matrix, SparseVector, SparseMatrixCSC, SpecialArrays}...)
+    X = SparseMatrixCSC[issparse(x) ? x : sparse(x) for x in Xin]
+    hcat(X...)
+end
+
+function vcat(Xin::Union{Vector, Matrix, SparseVector, SparseMatrixCSC, SpecialArrays}...)
+    X = SparseMatrixCSC[issparse(x) ? x : sparse(x) for x in Xin]
+    vcat(X...)
+end
+
+function hvcat(rows::Tuple{Vararg{Int}}, X::Union{Vector, Matrix, SparseVector, SparseMatrixCSC, SpecialArrays}...)
+    nbr = length(rows)  # number of block rows
+
+    tmp_rows = Array{SparseMatrixCSC}(nbr)
+    k = 0
+    @inbounds for i = 1 : nbr
+        tmp_rows[i] = hcat(X[(1 : rows[i]) + k]...)
+        k += rows[i]
+    end
+    vcat(tmp_rows...)
+end
+
+function cat(catdims, Xin::Union{Vector, Matrix, SparseVector, SparseMatrixCSC, SpecialArrays}...)
+    X = SparseMatrixCSC[issparse(x) ? x : sparse(x) for x in Xin]
+    T = promote_eltype(Xin...)
+    Base.cat_t(catdims, T, X...)
+end
+
+
 ### math functions
 
 ### Unary Map

test/linalg/special.jl

@@ -128,3 +128,41 @@ for typ in [UpperTriangular,LowerTriangular,Base.LinAlg.UnitUpperTriangular,Base
     @test Base.LinAlg.A_mul_Bc(atri,qrb[:Q]) ≈ full(atri) * qrb[:Q]'
     @test Base.LinAlg.A_mul_Bc!(copy(atri),qrb[:Q]) ≈ full(atri) * qrb[:Q]'
 end
+
+# Test that concatenations of combinations of special and other matrix types yield sparse arrays
+let
+    N = 4
+    # Test concatenating pairwise combinations of special matrices
+    diagmat = Diagonal(ones(N))
+    bidiagmat = Bidiagonal(ones(N), ones(N-1), true)
+    tridiagmat = Tridiagonal(ones(N-1), ones(N), ones(N-1))
+    symtridiagmat = SymTridiagonal(ones(N), ones(N-1))
+    specialmats = (diagmat, bidiagmat, tridiagmat, symtridiagmat)
+    for specialmata in specialmats, specialmatb in specialmats
+        @test issparse(hcat(specialmata, specialmatb))
+        @test issparse(vcat(specialmata, specialmatb))
+        @test issparse(hvcat((1,1), specialmata, specialmatb))
+        @test issparse(cat((1,2), specialmata, specialmatb))
+    end
+    # Test concatenating pairwise combinations of special matrices with sparse matrices,
+    # dense matrices, or dense vectors
+    densevec = ones(N)
+    densemat = diagm(ones(N))
+    spmat = spdiagm(ones(N))
+    for specialmat in specialmats
+        # --> Tests applicable only to pairs of matrices
+        for othermat in (spmat, densemat)
+            @test issparse(vcat(specialmat, othermat))
+            @test issparse(vcat(othermat, specialmat))
+        end
+        # --> Tests applicable also to pairs including vectors
+        for specialmat in specialmats, othermatorvec in (spmat, densemat, densevec)
+            @test issparse(hcat(specialmat, othermatorvec))
+            @test issparse(hcat(othermatorvec, specialmat))
+            @test issparse(hvcat((2,), specialmat, othermatorvec))
+            @test issparse(hvcat((2,), othermatorvec, specialmat))
+            @test issparse(cat((1,2), specialmat, othermatorvec))
+            @test issparse(cat((1,2), othermatorvec, specialmat))
+        end
+    end
+end

test/sparsedir/sparsevector.jl

@@ -404,6 +404,41 @@ let m = 80, n = 100
     @test full(V) == Vr
 end
 
+# Test that concatenations of combinations of sparse vectors with various other
+# matrix/vector types yield sparse arrays
+let
+    N = 4
+    spvec = spzeros(N)
+    spmat = spzeros(N, 1)
+    densevec = ones(N)
+    densemat = ones(N, 1)
+    diagmat = Diagonal(ones(4))
+    # Test that concatenations of pairwise combinations of sparse vectors with dense
+    # vectors/matrices, sparse matrices, or special matrices yield sparse arrays
+    for othervecormat in (densevec, densemat, spmat)
+        @test issparse(vcat(spvec, othervecormat))
+        @test issparse(vcat(othervecormat, spvec))
+    end
+    for othervecormat in (densevec, densemat, spmat, diagmat)
+        @test issparse(hcat(spvec, othervecormat))
+        @test issparse(hcat(othervecormat, spvec))
+        @test issparse(hvcat((2,), spvec, othervecormat))
+        @test issparse(hvcat((2,), othervecormat, spvec))
+        @test issparse(cat((1,2), spvec, othervecormat))
+        @test issparse(cat((1,2), othervecormat, spvec))
+    end
+    # The preceding tests should cover multi-way combinations of those types, but for good
+    # measure test a few multi-way combinations involving those types
+    @test issparse(vcat(spvec, densevec, spmat, densemat))
+    @test issparse(vcat(densevec, spvec, densemat, spmat))
+    @test issparse(hcat(spvec, densemat, spmat, densevec, diagmat))
+    @test issparse(hcat(densemat, spmat, spvec, densevec, diagmat))
+    @test issparse(hvcat((5,), diagmat, densevec, spvec, densemat, spmat))
+    @test issparse(hvcat((5,), spvec, densemat, diagmat, densevec, spmat))
+    @test issparse(cat((1,2), densemat, diagmat, spmat, densevec, spvec))
+    @test issparse(cat((1,2), spvec, diagmat, densevec, spmat, densemat))
+end
+
 
 ## sparsemat: combinations with sparse matrix