diff --git a/LICENSE.md b/LICENSE.md index e560a3e359284..fbd0efad2ef9c 100644 --- a/LICENSE.md +++ b/LICENSE.md @@ -66,7 +66,6 @@ their own licenses: The following components of Julia's standard library have separate licenses: - base/fftw.jl (see [FFTW](http://fftw.org/doc/License-and-Copyright.html)) -- base/sparse/csparse.jl (LGPL-2.1+) - base/linalg/umfpack.jl (see [SUITESPARSE](http://faculty.cse.tamu.edu/davis/suitesparse.html)) - base/linalg/cholmod.jl (see [SUITESPARSE](http://faculty.cse.tamu.edu/davis/suitesparse.html)) diff --git a/base/deprecated.jl b/base/deprecated.jl index a0b22e31af6c3..b22c1a3cc6523 100644 --- a/base/deprecated.jl +++ b/base/deprecated.jl @@ -769,6 +769,25 @@ end @deprecate slice view @deprecate sub view +# Point users to SuiteSparse +function ereach{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, k::Integer, parent::Vector{Ti}) + error(string("ereach(A, k, parent) now lives in package SuiteSparse.jl. Run", + "Pkg.add(\"SuiteSparse\") to install SuiteSparse on Julia v0.5.")) +end +function etree{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, postorder::Bool) + error(string("etree(A[, post]) now lives in package SuiteSparse.jl. Run", + "Pkg.add(\"SuiteSparse\") to install SuiteSparse on Julia v0.5.")) +end +etree(A::SparseMatrixCSC) = etree(A, false) +function csc_permute{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, pinv::Vector{Ti}, q::Vector{Ti}) + error(string("csc_permute(A, pinv, q) now lives in package SuiteSparse.jl. Run", + "Pkg.add(\"SuiteSparse\") to install SuiteSparse on Julia v0.5.")) +end +function symperm{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, pinv::Vector{Ti}) + error(string("symperm(A, pinv) now lives in package SuiteSparse.jl. Run,", + "Pkg.add(\"SuiteSparse\") to install SuiteSparse on Julia v0.5.")) +end + # During the 0.5 development cycle, do not add any deprecations below this line # To be deprecated in 0.6 diff --git a/base/sparse/csparse.jl b/base/sparse/csparse.jl deleted file mode 100644 index ab037530dbcbe..0000000000000 --- a/base/sparse/csparse.jl +++ /dev/null @@ -1,178 +0,0 @@ -# These functions are based on C functions in the CSparse library by Tim Davis. -# These are pure Julia implementations, and do not link to the CSparse library. -# CSparse can be downloaded from http://www.cise.ufl.edu/research/sparse/CSparse/CSparse.tar.gz -# CSparse is Copyright (c) 2006-2007, Timothy A. Davis and released under -# Lesser GNU Public License, version 2.1 or later. A copy of the license can be -# downloaded from http://www.gnu.org/licenses/lgpl-2.1.html - -# Because these functions are based on code covered by LGPL-2.1+ the same license -# must apply to the code in this file which is -# Copyright (c) 2013-2014 Viral Shah, Douglas Bates and other contributors - -# Based on Direct Methods for Sparse Linear Systems, T. A. Davis, SIAM, Philadelphia, Sept. 2006. -# Section 2.4: Triplet form -# http://www.cise.ufl.edu/research/sparse/CSparse/ - - -# Compute the elimination tree of A using triu(A) returning the parent vector. -# A root node is indicated by 0. This tree may actually be a forest in that -# there may be more than one root, indicating complete separability. -# A trivial example is speye(n, n) in which every node is a root. - -""" - etree(A[, post]) - -Compute the elimination tree of a symmetric sparse matrix `A` from `triu(A)` -and, optionally, its post-ordering permutation. -""" -function etree{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, postorder::Bool) - m,n = size(A) - Ap = A.colptr - Ai = A.rowval - parent = zeros(Ti, n) - ancestor = zeros(Ti, n) - for k in 1:n, p in Ap[k]:(Ap[k+1] - 1) - i = Ai[p] - while i != 0 && i < k - inext = ancestor[i] # inext = ancestor of i - ancestor[i] = k # path compression - if (inext == 0) parent[i] = k end # no anc., parent is k - i = inext - end - end - if !postorder return parent end - head = zeros(Ti,n) # empty linked lists - next = zeros(Ti,n) - for j in n:-1:1 # traverse in reverse order - if (parent[j] == 0); continue; end # j is a root - next[j] = head[parent[j]] # add j to list of its parent - head[parent[j]] = j - end - stack = Ti[] - sizehint!(stack, n) - post = zeros(Ti,n) - k = 1 - for j in 1:n - if (parent[j] != 0) continue end # skip j if it is not a root - push!(stack, j) # place j on the stack - while (!isempty(stack)) # while (stack is not empty) - p = stack[end] # p = top of stack - i = head[p] # i = youngest child of p - if (i == 0) - pop!(stack) - post[k] = p # node p is the kth postordered node - k += 1 - else - head[p] = next[i] # remove i from children of p - push!(stack, i) - end - end - end - parent, post -end - -etree(A::SparseMatrixCSC) = etree(A, false) - -# find nonzero pattern of Cholesky L[k,1:k-1] using etree and triu(A[:,k]) -# based on cs_ereach p. 43, "Direct Methods for Sparse Linear Systems" -function ereach{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, k::Integer, parent::Vector{Ti}) - m,n = size(A); Ap = A.colptr; Ai = A.rowval - s = Ti[]; sizehint!(s, n) # to be used as a stack - visited = falses(n) - visited[k] = true - for p in Ap[k]:(Ap[k+1] - 1) - i = Ai[p] # A[i,k] is nonzero - if i > k continue end # only use upper triangular part of A - while !visited[i] # traverse up etree - push!(s,i) # L[k,i] is nonzero - visited[i] = true - i = parent[i] - end - end - s -end - -# based on cs_permute p. 21, "Direct Methods for Sparse Linear Systems" -function csc_permute{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, pinv::Vector{Ti}, q::Vector{Ti}) - m, n = size(A) - Ap = A.colptr - Ai = A.rowval - Ax = A.nzval - lpinv = length(pinv) - if m != lpinv - throw(DimensionMismatch( - "the number of rows of sparse matrix A must equal the length of pinv, $m != $lpinv")) - end - lq = length(q) - if n != lq - throw(DimensionMismatch( - "the number of columns of sparse matrix A must equal the length of q, $n != $lq")) - end - if !isperm(pinv) || !isperm(q) - throw(ArgumentError("both pinv and q must be permutations")) - end - C = copy(A); Cp = C.colptr; Ci = C.rowval; Cx = C.nzval - nz = one(Ti) - for k in 1:n - Cp[k] = nz - j = q[k] - for t = Ap[j]:(Ap[j+1]-1) - Cx[nz] = Ax[t] - Ci[nz] = pinv[Ai[t]] - nz += one(Ti) - end - end - Cp[n + 1] = nz - (C.').' # double transpose to order the columns -end - - -# based on cs_symperm p. 21, "Direct Methods for Sparse Linear Systems" -# form A[p,p] for a symmetric A stored in the upper triangle -""" - symperm(A, p) - -Return the symmetric permutation of `A`, which is `A[p,p]`. `A` should be -symmetric, sparse, and only contain nonzeros in 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. -""" -function symperm{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, pinv::Vector{Ti}) - m, n = size(A) - if m != n - throw(DimensionMismatch("sparse matrix A must be square")) - end - Ap = A.colptr - Ai = A.rowval - Ax = A.nzval - if !isperm(pinv) - throw(ArgumentError("pinv must be a permutation")) - end - lpinv = length(pinv) - if n != lpinv - throw(DimensionMismatch( - "dimensions of sparse matrix A must equal the length of pinv, $((m,n)) != $lpinv")) - end - C = copy(A); Cp = C.colptr; Ci = C.rowval; Cx = C.nzval - w = zeros(Ti,n) - for j in 1:n # count entries in each column of C - j2 = pinv[j] - for p in Ap[j]:(Ap[j+1]-1) - (i = Ai[p]) > j || (w[max(pinv[i],j2)] += one(Ti)) - end - end - Cp[:] = cumsum(vcat(one(Ti),w)) - copy!(w,Cp[1:n]) # needed to be consistent with cs_cumsum - for j in 1:n - j2 = pinv[j] - for p = Ap[j]:(Ap[j+1]-1) - (i = Ai[p]) > j && continue - i2 = pinv[i] - ind = max(i2,j2) - Ci[q = w[ind]] = min(i2,j2) - w[ind] += 1 - Cx[q] = Ax[p] - end - end - (C.').' # double transpose to order the columns -end diff --git a/base/sparse/sparse.jl b/base/sparse/sparse.jl index 9e262e840a5da..93340cbb467fc 100644 --- a/base/sparse/sparse.jl +++ b/base/sparse/sparse.jl @@ -36,7 +36,6 @@ export AbstractSparseArray, AbstractSparseMatrix, AbstractSparseVector, include("abstractsparse.jl") include("sparsematrix.jl") include("sparsevector.jl") -include("csparse.jl") include("linalg.jl") if Base.USE_GPL_LIBS diff --git a/contrib/add_license_to_files.jl b/contrib/add_license_to_files.jl index 897f5d52ea1fb..3ee803b132035 100644 --- a/contrib/add_license_to_files.jl +++ b/contrib/add_license_to_files.jl @@ -34,7 +34,6 @@ const skipfiles = [ # files to check - already copyright # see: https://github.com/JuliaLang/julia/pull/11073#issuecomment-98099389 "../base/special/trig.jl", - "../base/sparse/csparse.jl", "../base/linalg/givens.jl", # "../src/abi_llvm.cpp", diff --git a/doc/stdlib/arrays.rst b/doc/stdlib/arrays.rst index c258443fd465c..e8902dc4de511 100644 --- a/doc/stdlib/arrays.rst +++ b/doc/stdlib/arrays.rst @@ -986,18 +986,6 @@ dense counterparts. The following functions are specific to sparse arrays. Create a random sparse vector of length ``m`` or sparse matrix of size ``m`` by ``n`` with the specified (independent) probability ``p`` of any entry being nonzero, where nonzero values are sampled from the normal distribution. The optional ``rng`` argument specifies a random number generator, see :ref:`Random Numbers `\ . -.. function:: etree(A[, post]) - - .. Docstring generated from Julia source - - Compute the elimination tree of a symmetric sparse matrix ``A`` from ``triu(A)`` and, optionally, its post-ordering permutation. - -.. function:: symperm(A, p) - - .. Docstring generated from Julia source - - Return the symmetric permutation of ``A``\ , which is ``A[p,p]``\ . ``A`` should be symmetric, sparse, and only contain nonzeros in 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. - .. function:: nonzeros(A) .. Docstring generated from Julia source diff --git a/test/sparsedir/sparse.jl b/test/sparsedir/sparse.jl index 2f9b5755708d3..38d679956a452 100644 --- a/test/sparsedir/sparse.jl +++ b/test/sparsedir/sparse.jl @@ -343,20 +343,6 @@ end @test full(spdiagm((ones(2), ones(2)), (0, -1), 3, 3)) == [1.0 0.0 0.0; 1.0 1.0 0.0; 0.0 1.0 0.0] -# elimination tree -## upper triangle of the pattern test matrix from Figure 4.2 of -## "Direct Methods for Sparse Linear Systems" by Tim Davis, SIAM, 2006 -rowval = Int32[1,2,2,3,4,5,1,4,6,1,7,2,5,8,6,9,3,4,6,8,10,3,5,7,8,10,11] -colval = Int32[1,2,3,3,4,5,6,6,6,7,7,8,8,8,9,9,10,10,10,10,10,11,11,11,11,11,11] -A = sparse(rowval, colval, ones(length(rowval))) -p = etree(A) -P,post = etree(A, true) -@test P == p -@test P == Int32[6,3,8,6,8,7,9,10,10,11,0] -@test post == Int32[2,3,5,8,1,4,6,7,9,10,11] -@test isperm(post) - - # issue #4986, reinterpret sfe22 = speye(Float64, 2) mfe22 = eye(Float64, 2) @@ -1043,13 +1029,9 @@ A = sparse(tril(rand(5,5))) @test istril(A) @test !istril(sparse(ones(5,5))) -# symperm -srand(1234321) -A = triu(sprand(10,10,0.2)) # symperm operates on upper triangle -perm = randperm(10) -@test symperm(A,perm).colptr == [1,1,2,4,5,5,6,8,8,9,10] - # droptol +srand(1234321) +A = triu(sprand(10,10,0.2)) @test Base.droptol!(A,0.01).colptr == [1,1,1,2,2,3,4,6,6,7,9] @test isequal(Base.droptol!(sparse([1], [1], [1]), 1), SparseMatrixCSC(1,1,Int[1,1],Int[],Int[])) @@ -1291,12 +1273,6 @@ if Base.USE_GPL_LIBS end @test_throws DimensionMismatch Base.SparseArrays.normestinv(sprand(3,5,.9)) -# csc_permute -A = sprand(10,10,0.2) -p = randperm(10) -q = randperm(10) -@test Base.SparseArrays.csc_permute(A, invperm(p), q) == full(A)[p, q] - # issue #13008 @test_throws ArgumentError sparse(collect(1:100), collect(1:100), fill(5,100), 5, 5) @test_throws ArgumentError sparse(Int[], collect(1:5), collect(1:5))