Skip to content

Instantly share code, notes, and snippets.

/-

Created November 18, 2014 08:43
Show Gist options
  • Save anonymous/aeb8cc230b506e921a84 to your computer and use it in GitHub Desktop.
Save anonymous/aeb8cc230b506e921a84 to your computer and use it in GitHub Desktop.
diff --cc base/sparse/csparse.jl
index e59f988,2cd1816..0000000
--- a/base/sparse/csparse.jl
+++ b/base/sparse/csparse.jl
@@@ -123,10 -123,13 +123,19 @@@ en
# Based on Direct Methods for Sparse Linear Systems, T. A. Davis, SIAM, Philadelphia, Sept. 2006.
# Section 2.5: Transpose
# http://www.cise.ufl.edu/research/sparse/CSparse/
++<<<<<<< HEAD
++=======
+ function transpose!{Tv,Ti}(T::SparseMatrixCSC{Tv,Ti}, S::SparseMatrixCSC{Tv,Ti})
+ (mS, nS) = size(S)
+ nnzS = nnz(S)
+ colptr_S = S.colptr
+ rowval_S = S.rowval
+ nzval_S = S.nzval
++>>>>>>> ba1f731... fix transpose! for sparse
+# NOTE: When calling transpose!(S,T), the colptr in the result matrix T must be set up
+# by counting the nonzeros in every row of S.
+function transpose!{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}, T::SparseMatrixCSC{Tv,Ti})
(mT, nT) = size(T)
colptr_T = T.colptr
rowval_T = T.rowval
@@@ -157,20 -161,17 +166,30 @@@ function transpose{Tv,Ti}(S::SparseMatr
rowval_T = Array(Ti, nnzS)
nzval_T = Array(Tv, nnzS)
+ colptr_T = zeros(Ti, nT+1)
+ colptr_T[1] = 1
+ @simd for i=1:nnz(S)
+ @inbounds colptr_T[rowval_S[i]+1] += 1
+ end
+ colptr_T = cumsum(colptr_T)
+
T = SparseMatrixCSC(mT, nT, colptr_T, rowval_T, nzval_T)
- return transpose!(S, T)
+ return transpose!(T, S)
end
++<<<<<<< HEAD
+# NOTE: When calling ctranspose!(S,T), the colptr in the result matrix T must be set up
+# by counting the nonzeros in every row of S.
+function ctranspose!{Tv,Ti}(S::SparseMatrixCSC{Tv,Ti}, T::SparseMatrixCSC{Tv,Ti})
++=======
+ function ctranspose!{Tv,Ti}(T::SparseMatrixCSC{Tv,Ti}, S::SparseMatrixCSC{Tv,Ti})
+ (mS, nS) = size(S)
+ nnzS = nnz(S)
+ colptr_S = S.colptr
+ rowval_S = S.rowval
+ nzval_S = S.nzval
+
++>>>>>>> ba1f731... fix transpose! for sparse
(mT, nT) = size(T)
colptr_T = T.colptr
rowval_T = T.rowval
@@@ -201,15 -203,8 +220,15 @@@ function ctranspose{Tv,Ti}(S::SparseMat
rowval_T = Array(Ti, nnzS)
nzval_T = Array(Tv, nnzS)
+ colptr_T = zeros(Ti, nT+1)
+ colptr_T[1] = 1
+ @inbounds for i=1:nnz(S)
+ colptr_T[rowval_S[i]+1] += 1
+ end
+ colptr_T = cumsum(colptr_T)
+
T = SparseMatrixCSC(mT, nT, colptr_T, rowval_T, nzval_T)
- return ctranspose!(S, T)
+ return ctranspose!(T, S)
end
# Compute the elimination tree of A using triu(A) returning the parent vector.
diff --cc base/sparse/sparsematrix.jl
index 2837a66,0613faf..0000000
--- a/base/sparse/sparsematrix.jl
+++ b/base/sparse/sparsematrix.jl
@@@ -1946,3 -2059,60 +1946,63 @@@ function diagm{Tv,Ti}(v::SparseMatrixCS
return SparseMatrixCSC{Tv,Ti}(n, n, colptr, rowval, nzval)
end
++<<<<<<< HEAD
++=======
+
+ # Sort all the indices in each column of a CSC sparse matrix
+ # sortSparseMatrixCSC!(A, sortindices = :sortcols) # Sort each column with sort()
+ # sortSparseMatrixCSC!(A, sortindices = :doubletranspose) # Sort with a double transpose
+ function sortSparseMatrixCSC!{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}; sortindices::Symbol = :sortcols)
+ if sortindices == :doubletranspose
+ nB, mB = size(A)
+ B = SparseMatrixCSC(mB, nB, Array(Ti, nB+1), similar(A.rowval), similar(A.nzval))
+ transpose!(B, A)
+ transpose!(A, B)
+ return A
+ end
+
+ m, n = size(A)
+ colptr = A.colptr; rowval = A.rowval; nzval = A.nzval
+
+ index = zeros(Ti, m)
+ row = zeros(Ti, m)
+ val = zeros(Tv, m)
+
+ for i = 1:n
+ @inbounds col_start = colptr[i]
+ @inbounds col_end = (colptr[i+1] - 1)
+
+ numrows = col_end - col_start + 1
+ if numrows <= 1
+ continue
+ elseif numrows == 2
+ f = col_start
+ s = f+1
+ if rowval[f] > rowval[s]
+ @inbounds rowval[f], rowval[s] = rowval[s], rowval[f]
+ @inbounds nzval[f], nzval[s] = nzval[s], nzval[f]
+ end
+ continue
+ end
+
+ jj = 1
+ @simd for j = col_start:col_end
+ @inbounds row[jj] = rowval[j]
+ @inbounds val[jj] = nzval[j]
+ jj += 1
+ end
+
+ sortperm!(pointer_to_array(pointer(index), numrows),
+ pointer_to_array(pointer(row), numrows))
+
+ jj = 1;
+ @simd for j = col_start:col_end
+ @inbounds rowval[j] = row[index[jj]]
+ @inbounds nzval[j] = val[index[jj]]
+ jj += 1
+ end
+ end
+
+ return A
+ end
++>>>>>>> ba1f731... fix transpose! for sparse
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment