Skip to content

Commit

Permalink
add sparse_*cat
Browse files Browse the repository at this point in the history
  • Loading branch information
@dkarrasch @LilithHafner
dkarrasch authored and LilithHafner committed Feb 22, 2022
1 parent 11758c3 commit 3a2e006
Show file tree
Hide file tree
Showing 6 changed files with 102 additions and 28 deletions.
6 changes: 6 additions & 0 deletions NEWS.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -141,6 +141,12 @@ Standard library changes

#### SparseArrays

* New sparse concatenation functions `sparse_hcat`, `sparse_vcat`, and `sparse_hvcat` return
`SparseMatrixCSC` output independent from the types of the input arguments. They make
concatenation behavior available, in which the presence of some special "sparse" matrix
argument resulted in sparse output by multiple dispatch. This is no longer possible after
making `LinearAlgebra.jl` independent from `SparseArrays.jl` ([#43127]).

#### Dates

#### Downloads
Expand Down
11 changes: 5 additions & 6 deletions stdlib/LinearAlgebra/src/uniformscaling.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -407,7 +407,7 @@ for (f, _f, dim, name) in ((:hcat, :_hcat, 1, "rows"), (:vcat, :_vcat, 2, "cols"
@eval begin
@inline $f(A::Union{AbstractVecOrMat,UniformScaling}...) = $_f(A...)
@inline $f(A::Union{AbstractVecOrMat,UniformScaling,Number}...) = $_f(A...)
function $_f(A::Union{AbstractVecOrMat,UniformScaling,Number}...)
function $_f(A::Union{AbstractVecOrMat,UniformScaling,Number}...; array_type = promote_to_array_type(A))
n = -1
for a in A
if !isa(a, UniformScaling)
Expand All @@ -420,14 +420,14 @@ for (f, _f, dim, name) in ((:hcat, :_hcat, 1, "rows"), (:vcat, :_vcat, 2, "cols"
end
end
n == -1 && throw(ArgumentError($("$f of only UniformScaling objects cannot determine the matrix size")))
return cat(promote_to_arrays(fill(n, length(A)), 1, promote_to_array_type(A), A...)..., dims=Val(3-$dim))
return cat(promote_to_arrays(fill(n, length(A)), 1, array_type, A...)..., dims=Val(3-$dim))
end
end
end

hvcat(rows::Tuple{Vararg{Int}}, A::Union{AbstractVecOrMat,UniformScaling}...) = _hvcat(rows, A...)
hvcat(rows::Tuple{Vararg{Int}}, A::Union{AbstractVecOrMat,UniformScaling,Number}...) = _hvcat(rows, A...)
function _hvcat(rows::Tuple{Vararg{Int}}, A::Union{AbstractVecOrMat,UniformScaling,Number}...)
function _hvcat(rows::Tuple{Vararg{Int}}, A::Union{AbstractVecOrMat,UniformScaling,Number}...; array_type = promote_to_array_type(A))
require_one_based_indexing(A...)
nr = length(rows)
sum(rows) == length(A) || throw(ArgumentError("mismatch between row sizes and number of arguments"))
Expand Down Expand Up @@ -479,14 +479,13 @@ function _hvcat(rows::Tuple{Vararg{Int}}, A::Union{AbstractVecOrMat,UniformScali
j += rows[i]
end
end
Atyp = promote_to_array_type(A)
Amat = promote_to_arrays(n, 1, Atyp, A...)
Amat = promote_to_arrays(n, 1, array_type, A...)
# We have two methods for promote_to_array_type, one returning Matrix and
# another one returning SparseMatrixCSC (in SparseArrays.jl). In the dense
# case, we cannot call hvcat for the promoted UniformScalings because this
# causes a stack overflow. In the sparse case, however, we cannot call
# typed_hvcat because we need a sparse output.
if Atyp == Matrix
if array_type == Matrix
return typed_hvcat(promote_eltype(Amat...), rows, Amat...)
else
return hvcat(rows, Amat...)
Expand Down
3 changes: 3 additions & 0 deletions stdlib/SparseArrays/docs/src/index.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -214,6 +214,9 @@ SparseArrays.nnz
SparseArrays.findnz
SparseArrays.spzeros
SparseArrays.spdiagm
SparseArrays.sparse_hcat
SparseArrays.sparse_vcat
SparseArrays.sparse_hvcat
SparseArrays.blockdiag
SparseArrays.sprand
SparseArrays.sprandn
Expand Down
3 changes: 2 additions & 1 deletion stdlib/SparseArrays/src/SparseArrays.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,8 @@ using Random: default_rng, AbstractRNG, randsubseq, randsubseq!
export AbstractSparseArray, AbstractSparseMatrix, AbstractSparseVector,
SparseMatrixCSC, SparseVector, blockdiag, droptol!, dropzeros!, dropzeros,
issparse, nonzeros, nzrange, rowvals, sparse, sparsevec, spdiagm,
sprand, sprandn, spzeros, nnz, permute, findnz
sprand, sprandn, spzeros, nnz, permute, findnz,
sparse_hcat, sparse_vcat, sparse_hvcat

include("abstractsparse.jl")
include("sparsematrix.jl")
Expand Down
49 changes: 49 additions & 0 deletions stdlib/SparseArrays/src/sparsevector.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1112,6 +1112,55 @@ _hvcat_rows(::Tuple{}, X::_SparseConcatGroup...) = ()
promote_to_array_type(A::Tuple{Vararg{Union{_SparseConcatGroup,UniformScaling}}}) = SparseMatrixCSC
promote_to_arrays_(n::Int, ::Type{SparseMatrixCSC}, J::UniformScaling) = sparse(J, n, n)

"""
sparse_hcat(A...)
Concatenate along dimension 2. Return a SparseMatrixCSC object.
!!! compat "Julia 1.8"
This method was added in Julia 1.8. It mimicks previous concatenation behavior, where
the concatenation with specialized "sparse" matrix types from LinearAlgebra.jl
automatically yielded sparse output even in the absence of any SparseArray argument.
"""
sparse_hcat(Xin::Union{AbstractVecOrMat,Number}...) = cat(map(_makesparse, Xin)..., dims=Val(2))
function sparse_hcat(X::Union{AbstractVecOrMat,UniformScaling,Number}...)
LinearAlgebra._hcat(X...; array_type = SparseMatrixCSC)
end

"""
sparse_vcat(A...)
Concatenate along dimension 1. Return a SparseMatrixCSC object.
!!! compat "Julia 1.8"
This method was added in Julia 1.8. It mimicks previous concatenation behavior, where
the concatenation with specialized "sparse" matrix types from LinearAlgebra.jl
automatically yielded sparse output even in the absence of any SparseArray argument.
"""
sparse_vcat(Xin::Union{AbstractVecOrMat,Number}...) = cat(map(_makesparse, Xin)..., dims=Val(1))
function sparse_vcat(X::Union{AbstractVecOrMat,UniformScaling,Number}...)
LinearAlgebra._vcat(X...; array_type = SparseMatrixCSC)
end

"""
sparse_hvcat(rows::Tuple{Vararg{Int}}, values...)
Sparse horizontal and vertical concatenation in one call. This function is called
for block matrix syntax. The first argument specifies the number of
arguments to concatenate in each block row.
!!! compat "Julia 1.8"
This method was added in Julia 1.8. It mimicks previous concatenation behavior, where
the concatenation with specialized "sparse" matrix types from LinearAlgebra.jl
automatically yielded sparse output even in the absence of any SparseArray argument.
"""
function sparse_hvcat(rows::Tuple{Vararg{Int}}, Xin::Union{AbstractVecOrMat,Number}...)
hvcat(rows, map(_makesparse, Xin)...)
end
function sparse_hvcat(rows::Tuple{Vararg{Int}}, X::Union{AbstractVecOrMat,UniformScaling,Number}...)
LinearAlgebra._hvcat(rows, X...; array_type = SparseMatrixCSC)
end

### math functions

### Unary Map
Expand Down
58 changes: 37 additions & 21 deletions stdlib/SparseArrays/test/sparse.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -262,19 +262,27 @@ end
# Test concatenating pairwise combinations of special matrices with sparse matrices,
# dense matrices, or dense vectors
spmat = spdiagm(0 => fill(1., N))
dmat = Array(spmat)
spvec = sparse(fill(1., N))
dvec = Array(spvec)
for specialmat in specialmats
# --> Tests applicable only to pairs of matrices
@test issparse(vcat(specialmat, spmat))
@test issparse(vcat(spmat, specialmat))
@test sparse_vcat(specialmat, dmat)::SparseMatrixCSC == vcat(specialmat, spmat)
@test sparse_vcat(dmat, specialmat)::SparseMatrixCSC == vcat(spmat, specialmat)
# --> Tests applicable also to pairs including vectors
for specialmat in specialmats, othermatorvec in (spmat, spvec)
@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(specialmat, othermatorvec; dims=(1,2)))
@test issparse(cat(othermatorvec, specialmat; dims=(1,2)))
for specialmat in specialmats, (smatorvec, dmatorvec) in ((spmat, dmat), (spvec, dvec))
@test issparse(hcat(specialmat, smatorvec))
@test sparse_hcat(specialmat, dmatorvec)::SparseMatrixCSC == hcat(specialmat, smatorvec)
@test issparse(hcat(smatorvec, specialmat))
@test sparse_hcat(dmatorvec, specialmat)::SparseMatrixCSC == hcat(smatorvec, specialmat)
@test issparse(hvcat((2,), specialmat, smatorvec))
@test sparse_hvcat((2,), specialmat, dmatorvec)::SparseMatrixCSC == hvcat((2,), specialmat, smatorvec)
@test issparse(hvcat((2,), smatorvec, specialmat))
@test sparse_hvcat((2,), dmatorvec, specialmat)::SparseMatrixCSC == hvcat((2,), smatorvec, specialmat)
@test issparse(cat(specialmat, smatorvec; dims=(1,2)))
@test issparse(cat(smatorvec, specialmat; dims=(1,2)))
end
end
end
Expand All @@ -300,8 +308,8 @@ end
symtridiagmat = SymTridiagonal(1:N, 1:(N-1))
sparseconcatmats = testfull ? (spmat, diagmat, bidiagmat, tridiagmat, symtridiagmat) : (spmat, diagmat)
# Concatenations involving strictly these types, un/annotated, should yield dense arrays
densevec = fill(1., N)
densemat = fill(1., N, N)
densevec = Array(spvec)
densemat = Array(spmat)
# Annotated collections
annodmats = [annot(densemat) for annot in annotations]
annospcmats = [annot(spmat) for annot in annotations]
Expand All @@ -321,6 +329,14 @@ end
@test issparse(vcat(annospcmat, othermat))
@test issparse(vcat(othermat, annospcmat))
end
for (smat, dmat) in zip(annospcmats, annodmats), specialmat in sparseconcatmats
@test sparse_hcat(dmat, specialmat)::SparseMatrixCSC == hcat(smat, specialmat)
@test sparse_hcat(specialmat, dmat)::SparseMatrixCSC == hcat(specialmat, smat)
@test sparse_vcat(dmat, specialmat)::SparseMatrixCSC == vcat(smat, specialmat)
@test sparse_vcat(specialmat, dmat)::SparseMatrixCSC == vcat(specialmat, smat)
@test sparse_hvcat((2,), dmat, specialmat)::SparseMatrixCSC == hvcat((2,), smat, specialmat)
@test sparse_hvcat((2,), specialmat, dmat)::SparseMatrixCSC == hvcat((2,), specialmat, smat)
end
# --> Tests applicable to pairs including other vectors or matrices
for other in (spvec, densevec, densemat, annodmats..., sparseconcatmats...)
@test issparse(hcat(annospcmat, other))
Expand Down Expand Up @@ -357,30 +373,30 @@ end
E = SparseMatrixCSC(rand(1,3))
F = SparseMatrixCSC(rand(3,1))
α = rand()
@test (hcat(A, 2I))::SparseMatrixCSC == hcat(A, Matrix(2I, 3, 3))
@test (hcat(A, 2I, I(3)))::SparseMatrixCSC == hcat(A, Matrix(2I, 3, 3), Matrix(I, 3, 3))
@test (hcat(E, α))::SparseMatrixCSC == hcat(E, [α])
@test (hcat(E, α, 2I))::SparseMatrixCSC == hcat(E, [α], fill(2, 1, 1))
@test (vcat(A, 2I))::SparseMatrixCSC == vcat(A, Matrix(2I, 4, 4))
@test (vcat(A, 2I))::SparseMatrixCSC == (vcat(A, 2I(4)))::SparseMatrixCSC == vcat(A, Matrix(2I, 4, 4))
@test (vcat(F, α))::SparseMatrixCSC == vcat(F, [α])
@test (vcat(F, α, 2I))::SparseMatrixCSC == vcat(F, [α], fill(2, 1, 1))
@test (vcat(F, α, 2I))::SparseMatrixCSC == (vcat(F, α, 2I(1)))::SparseMatrixCSC == vcat(F, [α], fill(2, 1, 1))
@test (hcat(C, 2I))::SparseMatrixCSC == C
@test_throws DimensionMismatch hcat(C, α)
@test (vcat(D, 2I))::SparseMatrixCSC == D
@test_throws DimensionMismatch vcat(D, α)
@test (hcat(I, 3I, A, 2I))::SparseMatrixCSC == hcat(Matrix(I, 3, 3), Matrix(3I, 3, 3), A, Matrix(2I, 3, 3))
@test (vcat(I, 3I, A, 2I))::SparseMatrixCSC == vcat(Matrix(I, 4, 4), Matrix(3I, 4, 4), A, Matrix(2I, 4, 4))
@test (hvcat((2,1,2), B, 2I, I, 3I, 4I))::SparseMatrixCSC ==
@test (hvcat((2,1,2), B, 2I, I(6), 3I, 4I))::SparseMatrixCSC ==
hvcat((2,1,2), B, Matrix(2I, 3, 3), Matrix(I, 6, 6), Matrix(3I, 3, 3), Matrix(4I, 3, 3))
@test hvcat((3,1), C, C, I, 3I)::SparseMatrixCSC == hvcat((2,1), C, C, Matrix(3I, 6,6))
@test hvcat((3,1), C, C, I, 3I)::SparseMatrixCSC == hvcat((2,1), C, C, Matrix(3I, 6, 6))
@test hvcat((2,2,2), I, 2I, 3I, 4I, C, C)::SparseMatrixCSC ==
hvcat((2,2,2), Matrix(I, 3, 3), Matrix(2I, 3,3 ), Matrix(3I, 3,3), Matrix(4I, 3,3), C, C)
@test hvcat((2,2,4), C, C, I, 2I, 3I, 4I, 5I, D)::SparseMatrixCSC ==
hvcat((2,2,4), C, C, Matrix(I, 3, 3), Matrix(2I,3,3),
Matrix(3I, 2, 2), Matrix(4I, 2, 2), Matrix(5I,2,2), D)
@test (hvcat((2,3,2), B, 2I, C, C, I, 3I, 4I))::SparseMatrixCSC ==
hvcat((2,2,2), Matrix(I, 3, 3), Matrix(2I, 3, 3), Matrix(3I, 3, 3), Matrix(4I, 3, 3), C, C)
@test hvcat((2,2,4), C, C, I(3), 2I, 3I, 4I, 5I, D)::SparseMatrixCSC ==
hvcat((2,2,4), C, C, Matrix(I, 3, 3), Matrix(2I, 3, 3),
Matrix(3I, 2, 2), Matrix(4I, 2, 2), Matrix(5I, 2, 2), D)
@test (hvcat((2,3,2), B, 2I(3), C, C, I, 3I, 4I))::SparseMatrixCSC ==
hvcat((2,2,2), B, Matrix(2I, 3, 3), C, C, Matrix(3I, 3, 3), Matrix(4I, 3, 3))
@test hvcat((3,2,1), C, C, I, B ,3I, 2I)::SparseMatrixCSC ==
hvcat((2,2,1), C, C, B, Matrix(3I,3,3), Matrix(2I,6,6))
@test hvcat((3,2,1), C, C, I, B, 3I(3), 2I)::SparseMatrixCSC ==
hvcat((2,2,1), C, C, B, Matrix(3I, 3, 3), Matrix(2I, 6, 6))
@test (hvcat((1,2), A, E, α))::SparseMatrixCSC == hvcat((1,2), A, E, [α]) == hvcat((1,2), A, E, α*I)
@test (hvcat((2,2), α, E, F, 3I))::SparseMatrixCSC == hvcat((2,2), [α], E, F, Matrix(3I, 3, 3))
@test (hvcat((2,2), 3I, F, E, α))::SparseMatrixCSC == hvcat((2,2), Matrix(3I, 3, 3), F, E, [α])
Expand Down

0 comments on commit 3a2e006

Please sign in to comment.