From c31a0a3f0e5b8680de85f8dfcfbf4bb3b61f7815 Mon Sep 17 00:00:00 2001 From: Dahua Lin Date: Sun, 24 May 2015 14:01:02 +0800 Subject: [PATCH 1/4] add several missing functors: SubFun, DivFun, and PowFun --- base/functors.jl | 12 ++++++++++++ 1 file changed, 12 insertions(+) diff --git a/base/functors.jl b/base/functors.jl index 646facee9d749..7bff9b87f4e30 100644 --- a/base/functors.jl +++ b/base/functors.jl @@ -36,9 +36,18 @@ call(::OrFun, x, y) = x | y immutable AddFun <: Func{2} end call(::AddFun, x, y) = x + y +immutable SubFun <: Func{2} end +call(::SubFun, x, y) = x - y + immutable MulFun <: Func{2} end call(::MulFun, x, y) = x * y +immutable DivFun <: Func{2} end +call(::DivFun, x, y) = x / y + +immutable PowFun <: Func{2} end +call(::PowFun, x, y) = x ^ y + immutable MaxFun <: Func{2} end call(::MaxFun, x, y) = scalarmax(x,y) @@ -119,7 +128,10 @@ function specialized_unary(f::Function) end function specialized_binary(f::Function) is(f, +) ? AddFun() : + is(f, -) ? SubFun() : is(f, *) ? MulFun() : + is(f, /) ? DivFun() : + is(f, ^) ? PowFun() : is(f, &) ? AndFun() : is(f, |) ? OrFun() : UnspecializedFun{2}(f) From d0c2443d11c47569c861e06232ea9bd531498628 Mon Sep 17 00:00:00 2001 From: Dahua Lin Date: Sun, 24 May 2015 14:33:05 +0800 Subject: [PATCH 2/4] Add sparse vector type and methods --- base/sparse.jl | 6 +- base/sparse/sparsevector.jl | 359 ++++++++++++++++++++++++++++++++++++ 2 files changed, 363 insertions(+), 2 deletions(-) create mode 100644 base/sparse/sparsevector.jl diff --git a/base/sparse.jl b/base/sparse.jl index d931e43082e4e..d25452c6a678e 100644 --- a/base/sparse.jl +++ b/base/sparse.jl @@ -2,7 +2,7 @@ module SparseMatrix -using Base: Func, AddFun, OrFun, ConjFun, IdFun +using Base: Func, AddFun, SubFun, OrFun, ConjFun, IdFun using Base.Sort: Forward using Base.LinAlg: AbstractTriangular @@ -12,12 +12,14 @@ import Base.promote_eltype import Base.@get! import Base.Broadcast.eltype_plus, Base.Broadcast.broadcast_shape -export AbstractSparseArray, AbstractSparseMatrix, AbstractSparseVector, SparseMatrixCSC, +export AbstractSparseArray, AbstractSparseMatrix, AbstractSparseVector, + SparseVector, SparseMatrixCSC, blkdiag, dense, droptol!, dropzeros!, etree, issparse, nnz, nonzeros, nzrange, rowvals, sparse, sparsevec, spdiagm, speye, spones, sprand, sprandbool, sprandn, spzeros, symperm include("sparse/abstractsparse.jl") +include("sparse/sparsevector.jl") include("sparse/sparsematrix.jl") include("sparse/csparse.jl") diff --git a/base/sparse/sparsevector.jl b/base/sparse/sparsevector.jl new file mode 100644 index 0000000000000..0f05fa11df8f2 --- /dev/null +++ b/base/sparse/sparsevector.jl @@ -0,0 +1,359 @@ + +# Sparse vector with compressed storage both indices and values. +# +# nzind[i] is the index of the i-th stored element, whose value is nzval[i] +# +# It assumes that nzind is in ascending order. +# +type SparseVector{Tv,Ti<:Integer} <: AbstractSparseVector{Tv,Ti} + n::Int # the number of elements + nzind::Vector{Ti} # the indices of nonzeros + nzval::Vector{Tv} # the values of nonzeros + + function SparseVector(n::Integer, nzind::Vector{Ti}, nzval::Vector{Tv}) + n >= 0 || throw(ArgumentError("The number of elements must be non-negative.")) + length(nzind) == length(nzval) || + throw(DimensionMismatch("The lengths of nzind and nzval are inconsistent.")) + new(Int(n), nzind, nzval) + end +end + +## Construction + +SparseVector{Tv,Ti<:Integer}(n::Integer, nzind::Vector{Ti}, nzval::Vector{Tv}) = + SparseVector{Tv,Ti}(n, nzind, nzval) + +SparseVector(n::Integer) = SparseVector(n, Int[], Float64[]) + +SparseVector{Tv}(::Type{Tv}, n::Integer) = SparseVector(n, Int[], Tv[]) + +SparseVector{Tv,Ti<:Integer}(::Type{Tv}, ::Type{Ti}, n::Integer) = SparseVector(n, Ti[], Tv[]) + +# construct from an associative map (Ti => Tv) +function SparseVector{Tv,Ti<:Integer}(n::Integer, s::Associative{Ti,Tv}) + m = length(s) + nzind = Array(Ti, 0) + nzval = Array(Tv, 0) + sizehint!(nzind, m) + sizehint!(nzval, m) + + for (k, v) in s + if v != zero(v) + push!(nzind, k) + push!(nzval, v) + end + end + + p = sortperm(nzind) + permute!(nzind, p) + permute!(nzval, p) + + return SparseVector{Tv,Ti}(Int(n), nzind, nzval) +end + + +## Basic properties + +length(x::SparseVector) = x.n +size(x::SparseVector) = (x.n,) + +nnz(x::SparseVector) = length(x.nzval) +countnz(x::SparseVector) = countnz(x.nzval) +nonzeros(x::SparseVector) = x.nzval + +## Element access + +function getindex{Tv}(x::SparseVector{Tv}, i::Int) + m = length(x.nzind) + ii = searchsortedfirst(x.nzind, i) + (ii <= m && x.nzind[ii] == i) ? x.nzval[ii] : zero(Tv) +end + +getindex(x::SparseVector, i::Integer) = x[Int(i)] + +function setindex!{Tv,Ti<:Integer}(x::SparseVector{Tv,Ti}, v::Tv, i::Ti) + nzind = x.nzind + nzval = x.nzval + + m = length(nzind) + k = searchsortedfirst(nzind, i) + if 1 <= k <= m && nzind[k] == i # i found + if v == zero(v) + deleteat!(nzind, k) + deleteat!(nzval, k) + else + nzval[k] = v + end + else # i not found + if v != zero(v) + insert!(nzind, k, i) + insert!(nzval, k, v) + end + end + x +end + +setindex!{Tv, Ti<:Integer}(x::SparseVector{Tv,Ti}, v, i::Integer) = + setindex!(x, convert(Tv, v), Ti(i)) + + +## Conversion + +# convert Vector to SparseVector +function convert{Tv,Ti<:Integer}(::Type{SparseVector{Tv,Ti}}, s::Vector) + n = length(s) + nzind = Array(Ti, 0) + nzval = Array(Tv, 0) + for i = 1:n + @inbounds v = s[i] + if v != zero(v) + push!(nzind, convert(Ti, i)) + push!(nzval, convert(Tv, v)) + end + end + return SparseVector(n, nzind, nzval) +end + +convert{Tv}(::Type{SparseVector{Tv}}, s::Vector{Tv}) = + convert(SparseVector{Tv,Int}, s) + +convert{Tv}(::Type{SparseVector}, s::Vector{Tv}) = + convert(SparseVector{Tv,Int}, s) + + +# convert between different types of SparseVector +convert{Tv,Ti,TvS,TiS}(::Type{SparseVector{Tv,Ti}}, s::SparseVector{TvS,TiS}) = + SparseVector{Tv,Ti}(s.n, convert(Vector{Ti}, s.nzind), convert(Vector{Tv}, s.nzval)) + +convert{Tv,TvS,TiS}(::Type{SparseVector{Tv}}, s::SparseVector{TvS,TiS}) = + SparseVector{Tv,TiS}(s.n, s.nzind, convert(Vector{Tv}, s.nzval)) + + +## Random sparse vector + +function sprand{T}(n::Integer, p::FloatingPoint, rfn::Function, ::Type{T}) + I = randsubseq(1:convert(Int, n), p) + V = rfn(T, length(I)) + SparseVector(n, I, V) +end + +function sprand(n::Integer, p::FloatingPoint, rfn::Function) + I = randsubseq(1:convert(Int, n), p) + V = rfn(length(I)) + SparseVector(n, I, V) +end + +sprand{T}(n::Integer, p::FloatingPoint, ::Type{T}) = sprand(n, p, rand, T) +sprand(n::Integer, p::FloatingPoint) = sprand(n, p, rand) + +sprandn(n::Integer, p::FloatingPoint) = sprand(n, p, randn) + + +## Array operations + +function full{Tv}(x::SparseVector{Tv}) + n = x.n + nzind = x.nzind + nzval = x.nzval + r = zeros(Tv, n) + for i = 1:length(nzind) + r[nzind[i]] = nzval[i] + end + return r +end + +vec(x::SparseVector) = x + +copy(x::SparseVector) = SparseVector(x.n, copy(x.nzind), copy(x.nzval)) + + +## Show + +function showarray(io::IO, x::SparseVector; + header::Bool=true, limit::Bool=Base._limit_output, + rows = Base.tty_size()[1], repr=false) + if header + print(io, "Sparse vector, length = ", x.n, + ", with ", nnz(x), " ", eltype(x), " entries:", "\n") + end + + if limit + half_screen_rows = div(rows - 8, 2) + else + half_screen_rows = typemax(Int) + end + pad = ndigits(x.n) + k = 0 + + for k = 1:length(x.nzind) + if k < half_screen_rows || k > nnz(x)-half_screen_rows + print(io, "\t", '[', rpad(x.nzind[k], pad), "] = ") + showcompact(io, x.nzval[k]) + print(io, "\n") + elseif k == half_screen_rows + print(io, sep, '\u22ee') + end + k += 1 + end +end + +show(io::IO, x::SparseVector) = Base.showarray(io, x) +writemime(io::IO, ::MIME"text/plain", x::SparseVector) = show(io, x) + + +## Computation + +Base.abs(x::SparseVector) = SparseVector(x.n, copy(x.nzind), abs(x.nzval)) +Base.abs2(x::SparseVector) = SparseVector(x.n, copy(x.nzind), abs2(x.nzval)) + +# zero-preserved: f(0, 0) -> 0 +function zero_preserve_map{Tx,Ty}(f, x::SparseVector{Tx}, y::SparseVector{Ty}) + R = typeof(_eval(op, zero(Tx), zero(Ty))) + n = length(x) + length(y) == n || throw(DimensionMismatch()) + + xnzind = x.nzind + xnzval = x.nzval + ynzind = y.nzind + ynzval = y.nzval + mx = length(xnzind) + my = length(ynzind) + + ix = 1 + iy = 1 + rind = Array(Int, 0) + rval = Array(R, 0) + sizehint!(rind, mx + my) + sizehint!(rval, mx + my) + + @inbounds while ix <= mx && iy <= my + jx = xnzind[ix] + jy = ynzind[iy] + + if jx == jy + v = f(xnzval[ix], ynzval[iy]) + if v != zero(v) + push!(rind, jx) + push!(rval, v) + end + ix += 1 + iy += 1 + elseif jx < jy + v = f(xnzval[i], zero(Ty)) + if v != zero(v) + push!(rind, jx) + push!(rval, v) + end + ix += 1 + else + v = f(zero(Tx), ynzval[iy]) + if v != zero(v) + push!(rind, jy) + push!(rval, v) + end + iy += 1 + end + end + + @inbounds while ix <= mx + v = f(xnzval[ix], zero(Ty)) + if v != zero(v) + push!(rind, xnzind[ix]) + push!(rval, v) + end + ix += 1 + end + + @inbounds while iy <= my + v = f(zero(Tx), ynzval[iy]) + if v != zero(v) + push!(rind, ynzind[iy]) + push!(rval, v) + end + iy += 1 + end + + return SparseVector(n, rind, rval) +end + +function map{Tx,Ty}(f, x::StridedVector{Tx}, y::SparseVector{Ty}) + R = typeof(_eval(op, zero(Tx), zero(Ty))) + n = length(x) + length(y) == n || throw(DimensionMismatch()) + + ynzind = y.nzind + ynzval = y.nzval + m = length(ynzind) + + dst = Array(R, n) + ii = 1 + @inbounds for i = 1:m + j = ynzind[i] + while ii < j + dst[ii] = f(x[ii], zero(Ty)) + ii += 1 + end + # at this point: ii == j + dst[j] = f(x[j], ynzval[i]) + ii += 1 + end + + @inbounds while ii <= n + dst[ii] = f(x[ii], zero(Ty)) + ii += 1 + end + return dst +end + +function map{Tx,Ty}(f, x::SparseVector{Tx}, y::StridedVector{Ty}) + R = typeof(_eval(op, zero(Tx), zero(Ty))) + n = length(x) + length(y) == n || throw(DimensionMismatch()) + + xnzind = x.nzind + xnzval = x.nzval + m = length(xnzind) + + dst = Array(R, n) + ii = 1 + @inbounds for i = 1:m + j = xnzind[i] + while ii < j + dst[ii] = f(zero(Tx), y[ii]) + ii += 1 + end + # at this point: ii == j + dst[j] = f(xnzval[i], y[j]) + ii += 1 + end + + @inbounds while ii <= n + dst[ii] = f(zero(Tx), y[ii]) + ii += 1 + end + return dst +end + + ++ (x::SparseVector, y::SparseVector) = zero_preserve_map(AddFun(), x, y) +- (x::SparseVector, y::SparseVector) = zero_preserve_map(SubFun(), x, y) ++ (x::StridedVector, y::SparseVector) = map(AddFun(), x, y) +- (x::StridedVector, y::SparseVector) = map(SubFun(), x, y) ++ (x::SparseVector, y::StridedVector) = map(AddFun(), x, y) +- (x::SparseVector, y::StridedVector) = map(SubFun(), x, y) + +.+ (x::SparseVector, y::SparseVector) = (x + y) +.- (x::SparseVector, y::SparseVector) = (x - y) + +.+ (x::StridedVector, y::SparseVector) = (x + y) +.- (x::StridedVector, y::SparseVector) = (x - y) + +.+ (x::SparseVector, y::StridedVector) = (x + y) +.- (x::SparseVector, y::StridedVector) = (x - y) + + +sum(x::SparseVector) = sum(x.nzval) +sumabs(x::SparseVector) = sumabs(x.nzval) +sumabs2(x::SparseVector) = sumabs2(x.nzval) + +vecnorm(x::SparseVector, p::Real=2) = vecnorm(x.nzval, p) From 81b79616d4571b4592aa9c67bdb7d5dda4eb5281 Mon Sep 17 00:00:00 2001 From: Dahua Lin Date: Sun, 24 May 2015 15:08:49 +0800 Subject: [PATCH 3/4] rename module Base.SparseMatrix to Base.Sparse, and export SparseVector --- base/exports.jl | 3 ++- base/sparse.jl | 2 +- base/sparse/cholmod.jl | 4 ++-- base/sparse/spqr.jl | 8 ++++---- base/sparse/umfpack.jl | 4 ++-- base/sysimg.jl | 2 +- 6 files changed, 12 insertions(+), 11 deletions(-) diff --git a/base/exports.jl b/base/exports.jl index 285134b6e5ac4..3d8b726b1d996 100644 --- a/base/exports.jl +++ b/base/exports.jl @@ -18,7 +18,7 @@ export BLAS, LAPACK, Serializer, - SparseMatrix, + Sparse, Docs, Markdown, @@ -103,6 +103,7 @@ export SharedArray, SharedMatrix, SharedVector, + SparseVector, SparseMatrixCSC, StatStruct, StepRange, diff --git a/base/sparse.jl b/base/sparse.jl index d25452c6a678e..cd0f0d29e9dbb 100644 --- a/base/sparse.jl +++ b/base/sparse.jl @@ -1,6 +1,6 @@ # This file is a part of Julia. License is MIT: http://julialang.org/license -module SparseMatrix +module Sparse using Base: Func, AddFun, SubFun, OrFun, ConjFun, IdFun using Base.Sort: Forward diff --git a/base/sparse/cholmod.jl b/base/sparse/cholmod.jl index eaff5b67f7f37..e0c607d130ab2 100644 --- a/base/sparse/cholmod.jl +++ b/base/sparse/cholmod.jl @@ -9,14 +9,14 @@ import Base.LinAlg: (\), A_mul_Bc, A_mul_Bt, Ac_ldiv_B, Ac_mul_B, At_ldiv_B, At_ cholfact, det, diag, ishermitian, isposdef, issym, ldltfact, logdet -import Base.SparseMatrix: sparse, nnz +import Base.Sparse: sparse, nnz export Dense, Factor, Sparse -using Base.SparseMatrix: AbstractSparseMatrix, SparseMatrixCSC, increment, indtype +using Base.Sparse: AbstractSparseMatrix, SparseMatrixCSC, increment, indtype ######### # Setup # diff --git a/base/sparse/spqr.jl b/base/sparse/spqr.jl index 39d27975aca39..eff5ab09259c6 100644 --- a/base/sparse/spqr.jl +++ b/base/sparse/spqr.jl @@ -36,12 +36,12 @@ const RTX_EQUALS_B = Int32(2) # solve R'*X=B or X = R'\B const RTX_EQUALS_ETB = Int32(3) # solve R'*X=E'*B or X = R'\(E'*B) -using Base.SparseMatrix: SparseMatrixCSC -using Base.SparseMatrix.CHOLMOD: C_Dense, C_Sparse, Dense, ITypes, Sparse, VTypes, common +using Base.Sparse: SparseMatrixCSC +using Base.Sparse.CHOLMOD: C_Dense, C_Sparse, Dense, ITypes, Sparse, VTypes, common import Base: size import Base.LinAlg: qrfact -import Base.SparseMatrix.CHOLMOD: convert, free! +import Base.Sparse.CHOLMOD: convert, free! @@ -143,4 +143,4 @@ function (\){T}(F::Factorization{T}, B::StridedVecOrMat{T}) convert(typeof(B), solve(RETX_EQUALS_B, F, QtB)) end -end # module \ No newline at end of file +end # module diff --git a/base/sparse/umfpack.jl b/base/sparse/umfpack.jl index 2a1e17c3628e8..e9ee7bd5d15c0 100644 --- a/base/sparse/umfpack.jl +++ b/base/sparse/umfpack.jl @@ -7,8 +7,8 @@ export UmfpackLU import Base: (\), Ac_ldiv_B, At_ldiv_B, findnz, getindex, show, size import Base.LinAlg: A_ldiv_B!, Ac_ldiv_B!, At_ldiv_B!, Factorization, det, lufact -importall Base.SparseMatrix -import Base.SparseMatrix: increment, increment!, decrement, decrement!, nnz +importall Base.Sparse +import Base.Sparse: increment, increment!, decrement, decrement!, nnz include("umfpack_h.jl") type MatrixIllConditionedException <: Exception diff --git a/base/sysimg.jl b/base/sysimg.jl index db3bc2f1ccff0..ef1afc69e0dc9 100644 --- a/base/sysimg.jl +++ b/base/sysimg.jl @@ -267,7 +267,7 @@ include("statistics.jl") # sparse matrices and sparse linear algebra include("sparse.jl") -importall .SparseMatrix +importall .Sparse # signal processing if USE_GPL_LIBS From 73cfb40cee2331cd34e80e6345feec9b7c508d0f Mon Sep 17 00:00:00 2001 From: Dahua Lin Date: Sun, 24 May 2015 15:36:51 +0800 Subject: [PATCH 4/4] add a test suite to SparseVector --- base/sparse/linalg.jl | 2 +- base/sparse/sparsevector.jl | 8 +- test/sparse.jl | 1 + test/sparsedir/cholmod.jl | 20 ++-- test/sparsedir/sparse.jl | 10 +- test/sparsedir/sparsevec.jl | 176 ++++++++++++++++++++++++++++++++++++ test/sparsedir/spqr.jl | 4 +- test/sparsedir/umfpack.jl | 6 +- 8 files changed, 202 insertions(+), 25 deletions(-) create mode 100644 test/sparsedir/sparsevec.jl diff --git a/base/sparse/linalg.jl b/base/sparse/linalg.jl index dfc53df672715..2c06d5f0094f8 100644 --- a/base/sparse/linalg.jl +++ b/base/sparse/linalg.jl @@ -173,7 +173,7 @@ function spmatmul{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, B::SparseMatrixCSC{Tv,Ti}; # The Gustavson algorithm does not guarantee the product to have sorted row indices. Cunsorted = SparseMatrixCSC(mA, nB, colptrC, rowvalC, nzvalC) - C = Base.SparseMatrix.sortSparseMatrixCSC!(Cunsorted, sortindices=sortindices) + C = Base.Sparse.sortSparseMatrixCSC!(Cunsorted, sortindices=sortindices) return C end diff --git a/base/sparse/sparsevector.jl b/base/sparse/sparsevector.jl index 0f05fa11df8f2..da650fa82e35e 100644 --- a/base/sparse/sparsevector.jl +++ b/base/sparse/sparsevector.jl @@ -208,7 +208,7 @@ Base.abs2(x::SparseVector) = SparseVector(x.n, copy(x.nzind), abs2(x.nzval)) # zero-preserved: f(0, 0) -> 0 function zero_preserve_map{Tx,Ty}(f, x::SparseVector{Tx}, y::SparseVector{Ty}) - R = typeof(_eval(op, zero(Tx), zero(Ty))) + R = typeof(f(zero(Tx), zero(Ty))) n = length(x) length(y) == n || throw(DimensionMismatch()) @@ -239,7 +239,7 @@ function zero_preserve_map{Tx,Ty}(f, x::SparseVector{Tx}, y::SparseVector{Ty}) ix += 1 iy += 1 elseif jx < jy - v = f(xnzval[i], zero(Ty)) + v = f(xnzval[ix], zero(Ty)) if v != zero(v) push!(rind, jx) push!(rval, v) @@ -277,7 +277,7 @@ function zero_preserve_map{Tx,Ty}(f, x::SparseVector{Tx}, y::SparseVector{Ty}) end function map{Tx,Ty}(f, x::StridedVector{Tx}, y::SparseVector{Ty}) - R = typeof(_eval(op, zero(Tx), zero(Ty))) + R = typeof(f(zero(Tx), zero(Ty))) n = length(x) length(y) == n || throw(DimensionMismatch()) @@ -306,7 +306,7 @@ function map{Tx,Ty}(f, x::StridedVector{Tx}, y::SparseVector{Ty}) end function map{Tx,Ty}(f, x::SparseVector{Tx}, y::StridedVector{Ty}) - R = typeof(_eval(op, zero(Tx), zero(Ty))) + R = typeof(f(zero(Tx), zero(Ty))) n = length(x) length(y) == n || throw(DimensionMismatch()) diff --git a/test/sparse.jl b/test/sparse.jl index 79dcd9c27f413..e091c94da3b68 100644 --- a/test/sparse.jl +++ b/test/sparse.jl @@ -1,5 +1,6 @@ # This file is a part of Julia. License is MIT: http://julialang.org/license +include("sparsedir/sparsevec.jl") include("sparsedir/sparse.jl") if Base.USE_GPL_LIBS include("sparsedir/umfpack.jl") diff --git a/test/sparsedir/cholmod.jl b/test/sparsedir/cholmod.jl index dce3b312528ce..5b8141b590304 100644 --- a/test/sparsedir/cholmod.jl +++ b/test/sparsedir/cholmod.jl @@ -3,7 +3,7 @@ srand(123) using Base.Test -using Base.SparseMatrix.CHOLMOD +using Base.Sparse.CHOLMOD # based on deps/SuiteSparse-4.0.2/CHOLMOD/Demo/ @@ -185,10 +185,10 @@ end # test Sparse constructor Symmetric and Hermitian input (and issym and ishermitian) ACSC = sprandn(10, 10, 0.3) + I -@test issym(Sparse(Symmetric(ACSC, :L))) -@test issym(Sparse(Symmetric(ACSC, :U))) -@test ishermitian(Sparse(Hermitian(complex(ACSC), :L))) -@test ishermitian(Sparse(Hermitian(complex(ACSC), :U))) +@test issym(CHOLMOD.Sparse(Symmetric(ACSC, :L))) +@test issym(CHOLMOD.Sparse(Symmetric(ACSC, :U))) +@test ishermitian(CHOLMOD.Sparse(Hermitian(complex(ACSC), :L))) +@test ishermitian(CHOLMOD.Sparse(Hermitian(complex(ACSC), :U))) # test Sparse constructor for c_SparseVoid (and read_sparse) let testfile = joinpath(tempdir(), "tmp.mtx") @@ -318,8 +318,8 @@ for elty in (Float64, Complex{Float64}) A1pdSparse = CHOLMOD.Sparse( A1pd.m, A1pd.n, - Base.SparseMatrix.decrement(A1pd.colptr), - Base.SparseMatrix.decrement(A1pd.rowval), + Base.Sparse.decrement(A1pd.colptr), + Base.Sparse.decrement(A1pd.rowval), A1pd.nzval) ## High level interface @@ -388,7 +388,7 @@ for elty in (Float64, Complex{Float64}) @test !isposdef(A1 + A1' |> t -> t - 2eigmax(full(t))*I) if elty <: Real - @test CHOLMOD.issym(Sparse(A1pd, 0)) + @test CHOLMOD.issym(CHOLMOD.Sparse(A1pd, 0)) @test CHOLMOD.Sparse(cholfact(Symmetric(A1pd, :L))) == CHOLMOD.Sparse(cholfact(A1pd)) F1 = CHOLMOD.Sparse(cholfact(Symmetric(A1pd, :L), shift=2)) F2 = CHOLMOD.Sparse(cholfact(A1pd, shift=2)) @@ -398,8 +398,8 @@ for elty in (Float64, Complex{Float64}) F2 = CHOLMOD.Sparse(ldltfact(A1pd, shift=2)) @test F1 == F2 else - @test !CHOLMOD.issym(Sparse(A1pd, 0)) - @test CHOLMOD.ishermitian(Sparse(A1pd, 0)) + @test !CHOLMOD.issym(CHOLMOD.Sparse(A1pd, 0)) + @test CHOLMOD.ishermitian(CHOLMOD.Sparse(A1pd, 0)) @test CHOLMOD.Sparse(cholfact(Hermitian(A1pd, :L))) == CHOLMOD.Sparse(cholfact(A1pd)) F1 = CHOLMOD.Sparse(cholfact(Hermitian(A1pd, :L), shift=2)) F2 = CHOLMOD.Sparse(cholfact(A1pd, shift=2)) diff --git a/test/sparsedir/sparse.jl b/test/sparsedir/sparse.jl index 23dade78deb39..82ac60364749d 100644 --- a/test/sparsedir/sparse.jl +++ b/test/sparsedir/sparse.jl @@ -188,8 +188,8 @@ for i = 1:5 a = sprand(10, 5, 0.7) b = sprand(5, 15, 0.3) @test maximum(abs(a*b - full(a)*full(b))) < 100*eps() - @test maximum(abs(Base.SparseMatrix.spmatmul(a,b,sortindices=:sortcols) - full(a)*full(b))) < 100*eps() - @test maximum(abs(Base.SparseMatrix.spmatmul(a,b,sortindices=:doubletranspose) - full(a)*full(b))) < 100*eps() + @test maximum(abs(Base.Sparse.spmatmul(a,b,sortindices=:sortcols) - full(a)*full(b))) < 100*eps() + @test maximum(abs(Base.Sparse.spmatmul(a,b,sortindices=:doubletranspose) - full(a)*full(b))) < 100*eps() @test full(kron(a,b)) == kron(full(a), full(b)) @test full(kron(full(a),b)) == kron(full(a), full(b)) @test full(kron(a,full(b))) == kron(full(a), full(b)) @@ -635,9 +635,9 @@ function test_getindex_algs{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, I::AbstractVector, ((minj < 1) || (maxj > n)) && BoundsError() end - (alg == 0) ? Base.SparseMatrix.getindex_I_sorted_bsearch_A(A, I, J) : - (alg == 1) ? Base.SparseMatrix.getindex_I_sorted_bsearch_I(A, I, J) : - Base.SparseMatrix.getindex_I_sorted_linear(A, I, J) + (alg == 0) ? Base.Sparse.getindex_I_sorted_bsearch_A(A, I, J) : + (alg == 1) ? Base.Sparse.getindex_I_sorted_bsearch_I(A, I, J) : + Base.Sparse.getindex_I_sorted_linear(A, I, J) end let M=2^14, N=2^4 diff --git a/test/sparsedir/sparsevec.jl b/test/sparsedir/sparsevec.jl new file mode 100644 index 0000000000000..72b4809302415 --- /dev/null +++ b/test/sparsedir/sparsevec.jl @@ -0,0 +1,176 @@ +# helpers + +function exact_equal(x::SparseVector, y::SparseVector) + x.n == y.n && x.nzind == y.nzind && x.nzval == y.nzval +end + +# empty sparse vectors + +x0 = SparseVector(8) +@test isa(x0, SparseVector{Float64,Int}) +@test length(x0) == 8 +@test nnz(x0) == 0 + +x0 = SparseVector(Float32, 8) +@test isa(x0, SparseVector{Float32,Int}) +@test length(x0) == 8 +@test nnz(x0) == 0 + +x0 = SparseVector(Float32, Int32, 8) +@test isa(x0, SparseVector{Float32, Int32}) +@test length(x0) == 8 +@test nnz(x0) == 0 + + +# construction + +x = SparseVector(8, [2, 5, 6], [1.25, -0.75, 3.5]) +x2 = SparseVector(8, [1, 2, 6, 7], [3.25, 4.0, -5.5, -6.0]) + +@test eltype(x) == Float64 +@test ndims(x) == 1 +@test length(x) == 8 +@test size(x) == (8,) +@test size(x,1) == 8 +@test size(x,2) == 1 +@test !isempty(x) + +@test countnz(x) == 3 +@test nnz(x) == 3 +@test nonzeros(x) == [1.25, -0.75, 3.5] + +dct = Dict{Int,Float64}() +dct[2] = 1.25 +dct[5] = -0.75 +dct[6] = 3.5 +xc = SparseVector(8, dct) +@test isa(xc, SparseVector{Float64,Int}) +@test exact_equal(x, xc) + +# full + +xf = zeros(8) +xf[2] = 1.25 +xf[5] = -0.75 +xf[6] = 3.5 +@test isa(full(x), Vector{Float64}) +@test full(x) == xf + +xf2 = zeros(8) +xf2[1] = 3.25 +xf2[2] = 4.0 +xf2[6] = -5.5 +xf2[7] = -6.0 +@test isa(full(x2), Vector{Float64}) +@test full(x2) == xf2 + + +# conversion + +xc = convert(SparseVector, xf) +@test isa(xc, SparseVector{Float64,Int}) +@test exact_equal(x, xc) + +xc = convert(SparseVector{Float32,Int}, x) +@test isa(xc, SparseVector{Float32,Int}) +@test exact_equal(x, xc) + +xc = convert(SparseVector{Float32}, x) +@test isa(xc, SparseVector{Float32,Int}) +@test exact_equal(x, xc) + +# copy + +xc = copy(x) +@test isa(xc, SparseVector{Float64,Int}) +@test !is(x.nzind, xc.nzval) +@test !is(x.nzval, xc.nzval) +@test exact_equal(x, xc) + +# getindex + +for i = 1:length(x) + @test x[i] == xf[i] +end + +# setindex + +xc = SparseVector(8) +xc[3] = 2.0 +@test exact_equal(xc, SparseVector(8, [3], [2.0])) + +xc = copy(x) +xc[5] = 2.0 +@test exact_equal(xc, SparseVector(8, [2, 5, 6], [1.25, 2.0, 3.5])) + +xc = copy(x) +xc[3] = 4.0 +@test exact_equal(xc, SparseVector(8, [2, 3, 5, 6], [1.25, 4.0, -0.75, 3.5])) + +xc[1] = 6.0 +@test exact_equal(xc, SparseVector(8, [1, 2, 3, 5, 6], [6.0, 1.25, 4.0, -0.75, 3.5])) + +xc[8] = -1.5 +@test exact_equal(xc, SparseVector(8, [1, 2, 3, 5, 6, 8], [6.0, 1.25, 4.0, -0.75, 3.5, -1.5])) + +xc = copy(x) +xc[5] = 0.0 +@test exact_equal(xc, SparseVector(8, [2, 6], [1.25, 3.5])) + +xc[6] = 0.0 +@test exact_equal(xc, SparseVector(8, [2], [1.25])) + +xc[2] = 0.0 +@test exact_equal(xc, SparseVector(8, Int[], Float64[])) + + +# sprand + +xr = sprand(1000, 0.3) +@test isa(xr, SparseVector{Float64,Int}) +@test length(xr) == 1000 +@test all(nonzeros(xr) .>= 0.0) + +xr = sprand(1000, 0.3, Float32) +@test isa(xr, SparseVector{Float32,Int}) +@test length(xr) == 1000 +@test all(nonzeros(xr) .>= 0.0) + +xr = sprandn(1000, 0.3) +@test isa(xr, SparseVector{Float64,Int}) +@test length(xr) == 1000 +@test any(nonzeros(xr) .> 0.0) && any(nonzeros(xr) .< 0.0) + +# abs and abs2 + +@test exact_equal(abs(x), SparseVector(8, [2, 5, 6], abs([1.25, -0.75, 3.5]))) +@test exact_equal(abs2(x), SparseVector(8, [2, 5, 6], abs2([1.25, -0.75, 3.5]))) + +# plus and minus + +xa = SparseVector(8, [1,2,5,6,7], [3.25,5.25,-0.75,-2.0,-6.0]) + +@test exact_equal(x + x, SparseVector(8, [2,5,6], [2.5,-1.5,7.0])) +@test exact_equal(x + x2, xa) + +xb = SparseVector(8, [1,2,5,6,7], [-3.25,-2.75,-0.75,9.0,6.0]) + +@test exact_equal(x - x, SparseVector(8, Int[], Float64[])) +@test exact_equal(x - x2, xb) + +@test full(x) + x2 == full(xa) +@test full(x) - x2 == full(xb) +@test x + full(x2) == full(xa) +@test x - full(x2) == full(xb) + + +# reduction + +@test sum(x) == 4.0 +@test sumabs(x) == 5.5 +@test sumabs2(x) == 14.375 + +@test vecnorm(x) == sqrt(14.375) +@test vecnorm(x, 1) == 5.5 +@test vecnorm(x, 2) == sqrt(14.375) +@test vecnorm(x, Inf) == 3.5 diff --git a/test/sparsedir/spqr.jl b/test/sparsedir/spqr.jl index 2a9559ee79cd4..8f97e80cb7626 100644 --- a/test/sparsedir/spqr.jl +++ b/test/sparsedir/spqr.jl @@ -2,8 +2,8 @@ using Base.Test -using Base.SparseMatrix.SPQR -using Base.SparseMatrix.CHOLMOD +using Base.Sparse.SPQR +using Base.Sparse.CHOLMOD m, n = 100, 10 nn = 100 diff --git a/test/sparsedir/umfpack.jl b/test/sparsedir/umfpack.jl index 69df2a54a5c53..f53982f51d160 100644 --- a/test/sparsedir/umfpack.jl +++ b/test/sparsedir/umfpack.jl @@ -6,14 +6,14 @@ do33 = ones(3) # based on deps/Suitesparse-4.0.2/UMFPACK/Demo/umfpack_di_demo.c -using Base.SparseMatrix.UMFPACK.increment! +using Base.Sparse.UMFPACK.increment! A0 = sparse(increment!([0,4,1,1,2,2,0,1,2,3,4,4]), increment!([0,4,0,2,1,2,1,4,3,2,1,2]), [2.,1.,3.,4.,-1.,-3.,3.,6.,2.,1.,4.,2.], 5, 5) for Tv in (Float64, Complex128) - for Ti in Base.SparseMatrix.UMFPACK.UMFITypes.types + for Ti in Base.Sparse.UMFPACK.UMFITypes.types A = convert(SparseMatrixCSC{Tv,Ti}, A0) lua = lufact(A) @test nnz(lua) == 18 @@ -37,7 +37,7 @@ for Tv in (Float64, Complex128) end Ac0 = complex(A0,A0) -for Ti in Base.SparseMatrix.UMFPACK.UMFITypes.types +for Ti in Base.Sparse.UMFPACK.UMFITypes.types Ac = convert(SparseMatrixCSC{Complex128,Ti}, Ac0) lua = lufact(Ac) L,U,p,q,Rs = lua[:(:)]