From 67f083a489abf1ff2fd0d201e2aa26975d4e81e5 Mon Sep 17 00:00:00 2001 From: "Viral B. Shah" Date: Sun, 25 May 2014 19:12:04 +0530 Subject: [PATCH] Reintroduce nnz and nonzeros as discussed in #6769 Deprecate nfilled for SparseMatrixCSC. Deprecate nonzeros for StridedArray and BitArray. Deprecate nnz for StridedArray. Update documentation for structural nonzeros. Update NEWS --- NEWS.md | 4 ++-- base/abstractarray.jl | 1 - base/array.jl | 18 --------------- base/bitarray.jl | 2 -- base/deprecated.jl | 17 +++++++++----- base/exports.jl | 2 +- base/linalg/cholmod.jl | 4 ++-- base/linalg/sparse.jl | 12 +++++----- base/sparse.jl | 2 +- base/sparse/csparse.jl | 14 +++++------ base/sparse/sparsematrix.jl | 46 +++++++++++++++++++------------------ doc/manual/arrays.rst | 15 +++++++----- doc/stdlib/base.rst | 16 ++++++------- doc/stdlib/sparse.rst | 2 ++ test/arrayops.jl | 2 +- test/sparse.jl | 4 ++-- 16 files changed, 75 insertions(+), 86 deletions(-) diff --git a/NEWS.md b/NEWS.md index 16a7ac7ce67bb..bcf59dc096ea3 100644 --- a/NEWS.md +++ b/NEWS.md @@ -320,7 +320,7 @@ Deprecated or removed * `factorize!` is deprecated in favor of `factorize`. ([#5526]) - * `nnz` is removed. Use `countnz` or `nfilled` instead ([#5538]) + * `nnz` counts the number of structural nonzeros in a sparse matrix. Use `countnz` for the actual number of nonzeros. ([#6769]) * `setfield` is renamed `setfield!` ([#5748]) @@ -383,7 +383,7 @@ Deprecated or removed [#4888]: https://github.com/JuliaLang/julia/pull/4888 [#5475]: https://github.com/JuliaLang/julia/pull/5475 [#5526]: https://github.com/JuliaLang/julia/pull/5526 -[#5538]: https://github.com/JuliaLang/julia/pull/5538 +[#6769]: https://github.com/JuliaLang/julia/pull/6769 [#5726]: https://github.com/JuliaLang/julia/pull/5726 [#5811]: https://github.com/JuliaLang/julia/pull/5811 [#5462]: https://github.com/JuliaLang/julia/pull/5462 diff --git a/base/abstractarray.jl b/base/abstractarray.jl index 9d2e2f8c81efd..010ad76efdbda 100644 --- a/base/abstractarray.jl +++ b/base/abstractarray.jl @@ -19,7 +19,6 @@ isreal{T<:Real,n}(x::AbstractArray{T,n}) = true ndims{T,n}(::AbstractArray{T,n}) = n ndims{T,n}(::Type{AbstractArray{T,n}}) = n ndims{T<:AbstractArray}(::Type{T}) = ndims(super(T)) -nfilled(t::AbstractArray) = length(t) length(t::AbstractArray) = prod(size(t))::Int endof(a::AbstractArray) = length(a) first(a::AbstractArray) = a[1] diff --git a/base/array.jl b/base/array.jl index 4c469d40229b1..f2f0f8daa6a18 100644 --- a/base/array.jl +++ b/base/array.jl @@ -1102,24 +1102,6 @@ function findnz{T}(A::StridedMatrix{T}) return (I, J, NZs) end -function nonzeros{T}(A::StridedArray{T}) - nnzA = countnz(A) - V = similar(A, T, nnzA) - count = 1 - if nnzA > 0 - for i=1:length(A) - Ai = A[i] - if Ai != 0 - V[count] = Ai - count += 1 - end - end - end - return V -end - -nonzeros(x::Number) = x == 0 ? Array(typeof(x),0) : [x] - function findmax(a) if isempty(a) error("array must be non-empty") diff --git a/base/bitarray.jl b/base/bitarray.jl index a83450b923abd..3253b08b9e46b 100644 --- a/base/bitarray.jl +++ b/base/bitarray.jl @@ -1385,8 +1385,6 @@ function findnz(B::BitMatrix) return I, J, trues(length(I)) end -nonzeros(B::BitArray) = trues(countnz(B)) - ## Reductions ## sum(A::BitArray, region) = reducedim(+, A, region, 0, Array(Int,reduced_dims(A,region))) diff --git a/base/deprecated.jl b/base/deprecated.jl index 4eea52cbea0c2..b6562d529d98f 100644 --- a/base/deprecated.jl +++ b/base/deprecated.jl @@ -380,6 +380,17 @@ end export mmread # 0.3 deprecations + +function nfilled(X) + depwarn("nfilled has been renamed to nnz", :nfilled) + nnz(X) +end +export nfilled + +@deprecate nonzeros(A::StridedArray) A[find(A)] +@deprecate nonzeros(B::BitArray) trues(countnz(B)) +@deprecate nnz(A::StridedArray) countnz(A) + @deprecate dense full export Stat @@ -448,12 +459,6 @@ Set{T<:Number}(xs::T...) = Set{T}(xs) # 0.3 discontinued functions -function nnz(X) - depwarn("nnz has been renamed to countnz and is no longer computed in constant time for sparse matrices. Instead, use nfilled() for the number of elements in a sparse matrix.", :nnz) - countnz(X) -end -export nnz - scale!{T<:Base.LinAlg.BlasReal}(X::Array{T}, s::Complex) = error("scale!: Cannot scale a real array by a complex value in-place. Use scale(X::Array{Real}, s::Complex) instead.") @deprecate which(f::Callable, args...) @which f(args...) diff --git a/base/exports.jl b/base/exports.jl index a0c47747b95cb..9fb5e4735adb0 100644 --- a/base/exports.jl +++ b/base/exports.jl @@ -536,7 +536,7 @@ export minimum, minmax, ndims, - nfilled, + nnz, nonzeros, nthperm!, nthperm, diff --git a/base/linalg/cholmod.jl b/base/linalg/cholmod.jl index 78bb732d48b12..fccc54c83c476 100644 --- a/base/linalg/cholmod.jl +++ b/base/linalg/cholmod.jl @@ -14,7 +14,7 @@ export # types using Base.LinAlg.UMFPACK # for decrement, increment, etc. import Base: (*), convert, copy, ctranspose, eltype, findnz, getindex, hcat, - isvalid, nfilled, show, size, sort!, transpose, vcat + isvalid, nnz, show, size, sort!, transpose, vcat import ..LinAlg: (\), A_mul_Bc, A_mul_Bt, Ac_ldiv_B, Ac_mul_B, At_ldiv_B, At_mul_B, cholfact, cholfact!, copy, det, diag, @@ -778,7 +778,7 @@ for Ti in (:Int32,:Int64) (Ptr{c_CholmodSparse{Tv,$Ti}},Ptr{c_CholmodSparse{Tv,$Ti}},Cint,Ptr{Uint8}), &A.c,&B.c,true,cmn($Ti)) end - function nfilled{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti}) + function nnz{Tv<:CHMVTypes}(A::CholmodSparse{Tv,$Ti}) ccall((@chm_nm "nnz" $Ti ,:libcholmod), Int, (Ptr{c_CholmodSparse{Tv,$Ti}},Ptr{Uint8}),&A.c,cmn($Ti)) end diff --git a/base/linalg/sparse.jl b/base/linalg/sparse.jl index afe2d7558c004..8b5119db88b30 100644 --- a/base/linalg/sparse.jl +++ b/base/linalg/sparse.jl @@ -347,7 +347,7 @@ function sparse_diff1{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}) m,n = size(S) m > 1 || return SparseMatrixCSC{Tv,Ti}(0, n, ones(n+1), Ti[], Tv[]) colptr = Array(Ti, n+1) - numnz = 2 * nfilled(S) # upper bound; will shrink later + numnz = 2 * nnz(S) # upper bound; will shrink later rowval = Array(Ti, numnz) nzval = Array(Tv, numnz) numnz = 0 @@ -387,7 +387,7 @@ function sparse_diff2{Tv,Ti}(a::SparseMatrixCSC{Tv,Ti}) m,n = size(a) colptr = Array(Ti, max(n,1)) - numnz = 2 * nfilled(a) # upper bound; will shrink later + numnz = 2 * nnz(a) # upper bound; will shrink later rowval = Array(Ti, numnz) nzval = Array(Tv, numnz) @@ -516,8 +516,8 @@ end # kron function kron{Tv,Ti}(a::SparseMatrixCSC{Tv,Ti}, b::SparseMatrixCSC{Tv,Ti}) - numnzA = nfilled(a) - numnzB = nfilled(b) + numnzA = nnz(a) + numnzB = nnz(b) numnz = numnzA * numnzB @@ -593,7 +593,7 @@ inv(A::SparseMatrixCSC) = error("The inverse of a sparse matrix can often be den function scale!{Tv,Ti}(C::SparseMatrixCSC{Tv,Ti}, A::SparseMatrixCSC, b::Vector) m, n = size(A) (n==length(b) && size(A)==size(C)) || throw(DimensionMismatch("")) - numnz = nfilled(A) + numnz = nnz(A) C.colptr = convert(Array{Ti}, A.colptr) C.rowval = convert(Array{Ti}, A.rowval) C.nzval = Array(Tv, numnz) @@ -606,7 +606,7 @@ end function scale!{Tv,Ti}(C::SparseMatrixCSC{Tv,Ti}, b::Vector, A::SparseMatrixCSC) m, n = size(A) (n==length(b) && size(A)==size(C)) || throw(DimensionMismatch("")) - numnz = nfilled(A) + numnz = nnz(A) C.colptr = convert(Array{Ti}, A.colptr) C.rowval = convert(Array{Ti}, A.rowval) C.nzval = Array(Tv, numnz) diff --git a/base/sparse.jl b/base/sparse.jl index 7a4f6df7fd363..aab1c16c64828 100644 --- a/base/sparse.jl +++ b/base/sparse.jl @@ -7,7 +7,7 @@ import Base.NonTupleType, Base.float export SparseMatrixCSC, blkdiag, dense, diag, diagm, droptol!, dropzeros!, etree, full, - getindex, ishermitian, issparse, issym, istril, istriu, + getindex, ishermitian, issparse, issym, istril, istriu, nnz, setindex!, sparse, sparsevec, spdiagm, speye, spones, sprand, sprandbool, sprandn, spzeros, symperm, trace, tril, tril!, triu, triu! diff --git a/base/sparse/csparse.jl b/base/sparse/csparse.jl index 94bc02dc5a650..e1410689debc5 100644 --- a/base/sparse/csparse.jl +++ b/base/sparse/csparse.jl @@ -128,7 +128,7 @@ function transpose!{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}, T::SparseMatrixCSC{Tv,Ti}) rowval_T = T.rowval nzval_T = T.nzval - nnzS = nfilled(S) + nnzS = nnz(S) colptr_S = S.colptr rowval_S = S.rowval nzval_S = S.nzval @@ -147,7 +147,7 @@ end function transpose{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}) (nT, mT) = size(S) - nnzS = nfilled(S) + nnzS = nnz(S) rowval_S = S.rowval rowval_T = Array(Ti, nnzS) @@ -155,7 +155,7 @@ function transpose{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}) colptr_T = zeros(Ti, nT+1) colptr_T[1] = 1 - @inbounds for i=1:nfilled(S) + @inbounds for i=1:nnz(S) colptr_T[rowval_S[i]+1] += 1 end colptr_T = cumsum(colptr_T) @@ -172,7 +172,7 @@ function ctranspose!{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}, T::SparseMatrixCSC{Tv,Ti} rowval_T = T.rowval nzval_T = T.nzval - nnzS = nfilled(S) + nnzS = nnz(S) colptr_S = S.colptr rowval_S = S.rowval nzval_S = S.nzval @@ -191,7 +191,7 @@ end function ctranspose{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}) (nT, mT) = size(S) - nnzS = nfilled(S) + nnzS = nnz(S) rowval_S = S.rowval rowval_T = Array(Ti, nnzS) @@ -199,7 +199,7 @@ function ctranspose{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}) colptr_T = zeros(Ti, nT+1) colptr_T[1] = 1 - @inbounds for i=1:nfilled(S) + @inbounds for i=1:nnz(S) colptr_T[rowval_S[i]+1] += 1 end colptr_T = cumsum(colptr_T) @@ -335,7 +335,7 @@ end # Section 2.7: Removing entries from a matrix # http://www.cise.ufl.edu/research/sparse/CSparse/ function fkeep!{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, f, other) - nzorig = nfilled(A) + nzorig = nnz(A) nz = 1 for j = 1:A.n p = A.colptr[j] # record current position diff --git a/base/sparse/sparsematrix.jl b/base/sparse/sparsematrix.jl index 0dfb97b0560ba..737db6832278a 100644 --- a/base/sparse/sparsematrix.jl +++ b/base/sparse/sparsematrix.jl @@ -15,7 +15,7 @@ SparseMatrixCSC{Tv,Ti}(m::Integer, n::Integer, colptr::Vector{Ti}, rowval::Vecto SparseMatrixCSC(int(m), int(n), colptr, rowval, nzval) size(S::SparseMatrixCSC) = (S.m, S.n) -nfilled(S::SparseMatrixCSC) = int(S.colptr[end]-1) +nnz(S::SparseMatrixCSC) = int(S.colptr[end]-1) countnz(S::SparseMatrixCSC) = countnz(S.nzval) function Base.showarray(io::IO, S::SparseMatrixCSC; @@ -24,7 +24,7 @@ function Base.showarray(io::IO, S::SparseMatrixCSC; # TODO: repr? if header - print(io, S.m, "x", S.n, " sparse matrix with ", nfilled(S), " ", eltype(S), " entries:") + print(io, S.m, "x", S.n, " sparse matrix with ", nnz(S), " ", eltype(S), " entries:") end if limit @@ -36,7 +36,7 @@ function Base.showarray(io::IO, S::SparseMatrixCSC; k = 0 sep = "\n\t" for col = 1:S.n, k = S.colptr[col] : (S.colptr[col+1]-1) - if k < half_screen_rows || k > nfilled(S)-half_screen_rows + if k < half_screen_rows || k > nnz(S)-half_screen_rows print(io, sep, '[', rpad(S.rowval[k], pad), ", ", lpad(col, pad), "] = ") showcompact(io, S.nzval[k]) elseif k == half_screen_rows @@ -96,7 +96,7 @@ function reinterpret{T,Tv,Ti,N}(::Type{T}, a::SparseMatrixCSC{Tv,Ti}, dims::NTup end mS,nS = dims mA,nA = size(a) - numnz = nfilled(a) + numnz = nnz(a) colptr = Array(Ti, nS+1) rowval = Array(Ti, numnz) nzval = reinterpret(T, a.nzval) @@ -112,7 +112,7 @@ function reshape{Tv,Ti}(a::SparseMatrixCSC{Tv,Ti}, dims::NTuple{2,Int}) end mS,nS = dims mA,nA = size(a) - numnz = nfilled(a) + numnz = nnz(a) colptr = Array(Ti, nS+1) rowval = Array(Ti, numnz) nzval = a.nzval @@ -203,7 +203,7 @@ function sparsevec(a::Vector) n = length(a) I = find(a) J = ones(Int, n) - V = nonzeros(a) + V = a[I] return sparse_IJ_sorted!(I,J,V,n,1,+) end @@ -287,7 +287,7 @@ function find(S::SparseMatrixCSC) end function findn{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}) - numnz = nfilled(S) + numnz = nnz(S) I = Array(Ti, numnz) J = Array(Ti, numnz) @@ -310,7 +310,7 @@ function findn{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}) end function findnz{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}) - numnz = nfilled(S) + numnz = nnz(S) I = Array(Ti, numnz) J = Array(Ti, numnz) V = Array(Tv, numnz) @@ -335,6 +335,8 @@ function findnz{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}) return (I, J, V) end +nonzeros(S::SparseMatrixCSC) = S.nzval + function sprand{T}(m::Integer, n::Integer, density::FloatingPoint, rng::Function,::Type{T}=eltype(rng(1))) 0 <= density <= 1 || throw(ArgumentError("$density not in [0,1]")) @@ -426,7 +428,7 @@ for (op, restype) in ((:iceil, Int), (:ceil, Nothing), @eval begin function ($op){Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}) - nfilledA = nfilled(A) + nfilledA = nnz(A) colptrB = Array(Ti, A.n+1) rowvalB = Array(Ti, nfilledA) nzvalB = Array($(restype==Nothing ? (:Tv) : restype), nfilledA) @@ -516,7 +518,7 @@ for (op, restype) in ( (:+, Nothing), (:-, Nothing), (:.*, Nothing), (m, n) = size(A) # TODO: Need better method to estimate result space - nnzS = nfilled(A) + nfilled(B) + nnzS = nnz(A) + nnz(B) colptrS = Array(Ti, A.n+1) rowvalS = Array(Ti, nnzS) nzvalS = Array($(restype==Nothing ? (:Tv) : restype), nnzS) @@ -671,7 +673,7 @@ function reducedim{Tv,Ti}(f::Function, A::SparseMatrixCSC{Tv,Ti}, region, v0) for i = 1 : A.n, j = A.colptr[i] : A.colptr[i+1]-1 S = f(S, A.nzval[j]) end - if nfilled(A) != A.m*A.n; S = f(S, zero(Tv)); end + if nnz(A) != A.m*A.n; S = f(S, zero(Tv)); end return fill(S, 1, 1) @@ -683,7 +685,7 @@ end function maximum{T}(A::SparseMatrixCSC{T}) isempty(A) && throw(ArgumentError("argument must not be empty")) m = maximum(A.nzval) - nfilled(A)!=length(A) ? max(m,zero(T)) : m + nnz(A)!=length(A) ? max(m,zero(T)) : m end maximum{T}(A::SparseMatrixCSC{T}, region) = @@ -692,7 +694,7 @@ maximum{T}(A::SparseMatrixCSC{T}, region) = function minimum{T}(A::SparseMatrixCSC{T}) isempty(A) && throw(ArgumentError("argument must not be empty")) m = minimum(A.nzval) - nfilled(A)!=length(A) ? min(m,zero(T)) : m + nnz(A)!=length(A) ? min(m,zero(T)) : m end minimum{T}(A::SparseMatrixCSC{T}, region) = @@ -701,7 +703,7 @@ minimum{T}(A::SparseMatrixCSC{T}, region) = sum{T}(A::SparseMatrixCSC{T}) = sum(A.nzval) sum{T}(A::SparseMatrixCSC{T}, region) = reducedim(+,A,region,zero(T)) -prod{T}(A::SparseMatrixCSC{T}) = nfilled(A)!=length(A) ? zero(T) : prod(A.nzval) +prod{T}(A::SparseMatrixCSC{T}) = nnz(A)!=length(A) ? zero(T) : prod(A.nzval) prod{T}(A::SparseMatrixCSC{T}, region) = reducedim(*,A,region,one(T)) #all(A::SparseMatrixCSC{Bool}, region) = reducedim(all,A,region,true) @@ -771,7 +773,7 @@ end # TODO: See if growing arrays is faster than pre-computing structure # and then populating nonzeros -# TODO: Use binary search in cases where nI >> nfilled(A[:,j]) or nI << nfilled(A[:,j]) +# TODO: Use binary search in cases where nI >> nnz(A[:,j]) or nI << nnz(A[:,j]) function getindex_I_sorted{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, I::AbstractVector, J::AbstractVector) (m, n) = size(A) @@ -1151,7 +1153,7 @@ function spset!{Tv,Ti<:Integer}(A::SparseMatrixCSC{Tv}, x::Tv, I::AbstractVector m, n = size(A) ((I[end] > m) || (J[end] > n)) && throw(DimensionMismatch("")) - nnzA = nfilled(A) + length(I) * length(J) + nnzA = nnz(A) + length(I) * length(J) colptrA = colptr = A.colptr rowvalA = rowval = A.rowval @@ -1254,7 +1256,7 @@ end function spdelete!{Tv,Ti<:Integer}(A::SparseMatrixCSC{Tv}, I::AbstractVector{Ti}, J::AbstractVector{Ti}) m, n = size(A) - nnzA = nfilled(A) + nnzA = nnz(A) (nnzA == 0) && (return A) !issorted(I) && (I = I[sortperm(I)]) @@ -1348,7 +1350,7 @@ function setindex!{Tv,Ti,T<:Integer}(A::SparseMatrixCSC{Tv,Ti}, B::SparseMatrixC colptrA = A.colptr; rowvalA = A.rowval; nzvalA = A.nzval colptrB = B.colptr; rowvalB = B.rowval; nzvalB = B.nzval - nnzS = nfilled(A) + nfilled(B) + nnzS = nnz(A) + nnz(B) colptrS = Array(Ti, n+1) rowvalS = Array(Ti, nnzS) nzvalS = Array(Tv, nnzS) @@ -1466,7 +1468,7 @@ function vcat(X::SparseMatrixCSC...) Ti = promote_type(map(x->eltype(x.rowval), X)...) colptr = Array(Ti, n + 1) - nnzX = [ nfilled(x) for x in X ] + nnzX = [ nnz(x) for x in X ] nnz_res = sum(nnzX) rowval = Array(Ti, nnz_res) nzval = Array(Tv, nnz_res) @@ -1504,7 +1506,7 @@ function hcat(X::SparseMatrixCSC...) Ti = promote_type(map(x->eltype(x.rowval), X)...) colptr = Array(Ti, n + 1) - nnzX = [ nfilled(x) for x in X ] + nnzX = [ nnz(x) for x in X ] nnz_res = sum(nnzX) rowval = Array(Ti, nnz_res) nzval = Array(Tv, nnz_res) @@ -1545,7 +1547,7 @@ function blkdiag(X::SparseMatrixCSC...) Ti = promote_type(map(x->eltype(x.rowval), X)...) colptr = Array(Ti, n + 1) - nnzX = [ nfilled(x) for x in X ] + nnzX = [ nnz(x) for x in X ] nnz_res = sum(nnzX) rowval = Array(Ti, nnz_res) nzval = Array(Tv, nnz_res) @@ -1718,7 +1720,7 @@ function diagm{Tv,Ti}(v::SparseMatrixCSC{Tv,Ti}) end n = length(v) - numnz = nfilled(v) + numnz = nnz(v) colptr = Array(Ti, n+1) rowval = Array(Ti, numnz) nzval = Array(Tv, numnz) diff --git a/doc/manual/arrays.rst b/doc/manual/arrays.rst index ee3bc98dc7e06..af6d8d50f5eca 100644 --- a/doc/manual/arrays.rst +++ b/doc/manual/arrays.rst @@ -528,12 +528,15 @@ The row indices in every column need to be sorted. If your `SparseMatrixCSC` object contains unsorted row indices, one quick way to sort them is by doing a double transpose. -In some applications, it is convenient to store explicit zero values in -a `SparseMatrixCSC`. These *are* accepted by functions in ``Base`` (but -there is no guarantee that they will be preserved in mutating operations). -Because of this, ``countnz`` is not a constant-time operation; instead, -``nfilled`` should be used to obtain the number of elements in a sparse -matrix. +In some applications, it is convenient to store explicit zero values +in a `SparseMatrixCSC`. These *are* accepted by functions in ``Base`` +(but there is no guarantee that they will be preserved in mutating +operations). Such explicitly stored zeros are treated as structural +nonzeros by many routines. The ``nnz`` function returns the number of +elements explicitly stored in the sparse data structure, +including structural nonzeros. In order to count the exact number of actual +values that are nonzero, use ``countnz``, which inspects every stored +element of a sparse matrix. Sparse matrix constructors -------------------------- diff --git a/doc/stdlib/base.rst b/doc/stdlib/base.rst index f6679fbe0e617..e24de4fadbc0f 100644 --- a/doc/stdlib/base.rst +++ b/doc/stdlib/base.rst @@ -3837,8 +3837,9 @@ Indexing, Assignment, and Concatenation .. function:: find(A) - Return a vector of the linear indexes of the non-zeros in ``A`` (determined by ``A[i]!=0``). - A common use of this is to convert a boolean array to an array of indexes of the ``true`` + Return a vector of the linear indexes of the non-zeros in ``A`` + (determined by ``A[i]!=0``). A common use of this is to convert a + boolean array to an array of indexes of the ``true`` elements. .. function:: find(f,A) @@ -3847,16 +3848,13 @@ Indexing, Assignment, and Concatenation .. function:: findn(A) - Return a vector of indexes for each dimension giving the locations of the non-zeros in ``A`` - (determined by ``A[i]!=0``). + Return a vector of indexes for each dimension giving the locations of the non-zeros in ``A`` (determined by ``A[i]!=0``). .. function:: findnz(A) - Return a tuple ``(I, J, V)`` where ``I`` and ``J`` are the row and column indexes of the non-zero values in matrix ``A``, and ``V`` is a vector of the non-zero values. - -.. function:: nonzeros(A) - - Return a vector of the non-zero values in array ``A`` (determined by ``A[i]!=0``). + Return a tuple ``(I, J, V)`` where ``I`` and ``J`` are the row and + column indexes of the non-zero values in matrix ``A``, and ``V`` is + a vector of the non-zero values. .. function:: findfirst(A) diff --git a/doc/stdlib/sparse.rst b/doc/stdlib/sparse.rst index cee6c6b8bc42d..0f6410eeeeae9 100644 --- a/doc/stdlib/sparse.rst +++ b/doc/stdlib/sparse.rst @@ -75,4 +75,6 @@ Sparse matrices support much of the same set of operations as dense matrices. Th Return the symmetric permutation of A, which is ``A[p,p]``. A should be symmetric and sparse, where only the upper triangular part of the matrix is stored. This algorithm ignores the lower triangular part of the matrix. Only the upper triangular part of the result is returned as well. +.. function:: nonzeros(A) + Return a vector of the structural nonzero values in sparse matrix ``A``. This includes zeros that are explicitly stored in the sparse matrix. The returned vector points directly to the internal nonzero storage of ``A``, and any modifications to the returned vector will mutate ``A`` as well. diff --git a/test/arrayops.jl b/test/arrayops.jl index ddcbb733d025c..e911538e34a40 100644 --- a/test/arrayops.jl +++ b/test/arrayops.jl @@ -3,7 +3,7 @@ ## basics @test length([1, 2, 3]) == 3 -@test nfilled([1, 2, 3]) == 3 +@test countnz([1, 2, 3]) == 3 a = ones(4) b = a+a diff --git a/test/sparse.jl b/test/sparse.jl index 63722d74d2d2a..be570f904c06e 100644 --- a/test/sparse.jl +++ b/test/sparse.jl @@ -184,14 +184,14 @@ mfe22 = eye(Float64, 2) @test_throws DimensionMismatch sparsevec([3,5,7],[0.1,0.0,3.2],4) # issue #5169 -@test nfilled(sparse([1,1],[1,2],[0.0,-0.0])) == 0 +@test nnz(sparse([1,1],[1,2],[0.0,-0.0])) == 0 # issue #5386 K,J,V = findnz(SparseMatrixCSC(2,1,[1,3],[1,2],[1.0,0.0])) @test length(K) == length(J) == length(V) == 1 # issue #5437 -@test nfilled(sparse([1,2,3],[1,2,3],[0.0,1.0,2.0])) == 2 +@test nnz(sparse([1,2,3],[1,2,3],[0.0,1.0,2.0])) == 2 # issue #5824 @test sprand(4,5,0.5).^0 == sparse(ones(4,5))