Skip to content

Commit

Permalink
improve performance of issym and ishermitian
Browse files Browse the repository at this point in the history
  • Loading branch information
KristofferC committed May 20, 2015
1 parent 8a1508d commit 230c8f1
Show file tree
Hide file tree
Showing 4 changed files with 73 additions and 8 deletions.
3 changes: 3 additions & 0 deletions base/functors.jl
Original file line number Diff line number Diff line change
Expand Up @@ -24,6 +24,9 @@ call(::ExpFun, x) = exp(x)
immutable LogFun <: Func{1} end
call(::LogFun, x) = log(x)

immutable ConjFun <: Func{1} end
call(::ConjFun, x) = conj(x)

immutable AndFun <: Func{2} end
call(::AndFun, x, y) = x & y

Expand Down
2 changes: 1 addition & 1 deletion base/sparse.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@

module SparseMatrix

using Base: Func, AddFun, OrFun
using Base: Func, AddFun, OrFun, ConjFun, IdFun
using Base.Sort: Forward
using Base.LinAlg: AbstractTriangular

Expand Down
38 changes: 31 additions & 7 deletions base/sparse/sparsematrix.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2376,17 +2376,41 @@ end
squeeze(S::SparseMatrixCSC, dims::Dims) = throw(ArgumentError("squeeze is not available for sparse matrices"))

## Structure query functions
issym(A::SparseMatrixCSC) = is_hermsym(A, IdFun())

function issym(A::SparseMatrixCSC)
m, n = size(A)
if m != n; return false; end
return countnz(A - A.') == 0
end
ishermitian(A::SparseMatrixCSC) = is_hermsym(A, ConjFun())

function ishermitian(A::SparseMatrixCSC)
function is_hermsym(A::SparseMatrixCSC, check::Func)
m, n = size(A)
if m != n; return false; end
return countnz(A - A') == 0

colptr = A.colptr
rowval = A.rowval
nzval = A.nzval
tracker = copy(A.colptr)
@inbounds for col = 1:A.n
for p = tracker[col]:colptr[col+1]-1
val = nzval[p]
val == 0 && continue # In case of explicit zeros
row = rowval[p]
row < col && return false
if row == col
if val != check(val)
return false
end
else
row2 = rowval[tracker[row]]
row2 > col && return false
if row2 == col
if val != check(nzval[tracker[row]])
return false
end
tracker[row] += 1
end
end
end
end
return true
end

function istriu{Tv}(A::SparseMatrixCSC{Tv})
Expand Down
38 changes: 38 additions & 0 deletions test/sparsedir/sparse.jl
Original file line number Diff line number Diff line change
Expand Up @@ -930,3 +930,41 @@ A = sparse([1.0])
@test norm(A) == 1.0
@test_throws ArgumentError norm(sprand(5,5,0.2),3)
@test_throws ArgumentError norm(sprand(5,5,0.2),2)

# test ishermitian and issym real matrices
A = speye(5,5)
@test ishermitian(A) == true
@test issym(A) == true
A[1,3] = 1.0
@test ishermitian(A) == false
@test issym(A) == false
A[3,1] = 1.0
@test ishermitian(A) == true
@test issym(A) == true

# test ishermitian and issym complex matrices
A = speye(5,5) + im*speye(5,5)
@test ishermitian(A) == false
@test issym(A) == true
A[1,4] = 1.0 + im
@test ishermitian(A) == false
@test issym(A) == false

A = speye(Complex128, 5,5)
A[3,2] = 1.0 + im
@test ishermitian(A) == false
@test issym(A) == false
A[2,3] = 1.0 - im
@test ishermitian(A) == true
@test issym(A) == false

A = sparse(zeros(5,5))
@test ishermitian(A) == true
@test issym(A) == true

# Test with explicit zeros
A = speye(Complex128, 5,5)
A[3,1] = 2
A.nzval[2] = 0.0
@test ishermitian(A) == true
@test issym(A) == true

0 comments on commit 230c8f1

Please sign in to comment.