Add the issym/ishermitian based on norm(A - A', Inf)<tol*norm(A,Inf). fix #10298.

Author: dhoegh

base/docs/helpdb.jl

@@ -7211,9 +7211,9 @@ Given `@osx? a : b`, do `a` on OS X and `b` elsewhere. See documentation for Han
 :@osx
 
 doc"""
-    ishermitian(A) -> Bool
+    ishermitian(A; tol=eps()*norm(A,1)) -> Bool
 
-Test whether a matrix is Hermitian.
+Test whether a matrix is Hermitian. The `tol` keyword is only used for matrices of floating point or `Complex{T<:AbstractFloat}` numbers, where it tests `norm(A - A', Inf) < tol*norm(A,Inf)`.
 """
 ishermitian
 
@@ -7345,9 +7345,9 @@ The arguments to a function or constructor are outside the valid domain.
 DomainError
 
 doc"""
-    issym(A) -> Bool
+    issym(A; tol=eps()*norm(A,1)) -> Bool
 
-Test whether a matrix is symmetric.
+Test whether a matrix is symmetric. The `tol` keyword is only used for matrices of floating point or `Complex{T<:AbstractFloat}` numbers, where it tests `norm(A - A', Inf) < tol*norm(A,Inf)`.
 """
 issym
 

base/linalg/generic.jl

@@ -350,7 +350,7 @@ condskeel{T<:Integer}(A::AbstractMatrix{T}, p::Real=Inf) = norm(abs(inv(float(A)
 condskeel(A::AbstractMatrix, x::AbstractVector, p::Real=Inf) = norm(abs(inv(A))*abs(A)*abs(x), p)
 condskeel{T<:Integer}(A::AbstractMatrix{T}, x::AbstractVector, p::Real=Inf) = norm(abs(inv(float(A)))*abs(A)*abs(x), p)
 
-function issym(A::AbstractMatrix)
+function _issym(A::AbstractMatrix)
     m, n = size(A)
     if m != n
         return false
@@ -365,7 +365,7 @@ end
 
 issym(x::Number) = true
 
-function ishermitian(A::AbstractMatrix)
+function _ishermitian(A::AbstractMatrix)
     m, n = size(A)
     if m != n
         return false
@@ -378,6 +378,37 @@ function ishermitian(A::AbstractMatrix)
     return true
 end
 
+ishermitian(A::AbstractMatrix) = _ishermitian(A::AbstractMatrix)
+issym(A::AbstractMatrix) = _issym(A::AbstractMatrix)
+
+for (f,t) in ((:ishermitian, :ctranspose),(:issym, :transpose))
+    eval(quote
+        # ishermitian and issym functions test whether norm(A - A', Inf)<tol*norm(A,Inf) for FloatingPoint matrices
+        function $f{T<:AbstractFloat}(A::Union{AbstractMatrix{Complex{T}}, AbstractMatrix{T}},tol=eps()*norm(A,1))
+            # This is performed due to if tol==0. _issym/_ishermitian is a lot faster for exact equality.
+            tol == 0. && return $(symbol("_", f))(A)
+            m,n = size(A)
+            Tnorm = typeof(float(real(zero(T))))
+            Tsum = promote_type(Float64,Tnorm)
+            nrm_diff::Tsum = 0
+            nrm_A::Tsum = 0
+            @inbounds begin
+                for i = 1:m
+                    nrmi_diff::Tsum = 0
+                    nrmi_A::Tsum = 0
+                    for j = 1:n
+                        nrmi_diff += abs(A[i,j]-$t(A[j,i]))
+                        nrmi_A += abs(A[i,j])
+                    end
+                    nrm_diff = max(nrm_diff,nrmi_diff)
+                    nrm_A = max(nrm_A,nrmi_A)
+                end
+            end
+            return nrm_diff < tol*nrm_A
+        end
+    end)
+end
+
 ishermitian(x::Number) = (x == conj(x))
 
 function istriu(A::AbstractMatrix)

base/linalg/symmetric.jl

@@ -105,6 +105,10 @@ function triu(A::Symmetric, k::Integer=0)
     end
 end
 
+#Test whether a matrix is positive-definite
+isposdef!(A::HermOrSym) = isposdef!(A.data, symbol(A.uplo))
+isposdef{T}(A::HermOrSym{T}) = (S = typeof(sqrt(one(T))); isposdef!(S == T ? copy(A.data) : convert(AbstractMatrix{S}, A.data), symbol(A.uplo)))
+
 ## Matvec
 A_mul_B!{T<:BlasFloat,S<:StridedMatrix}(y::StridedVector{T}, A::Symmetric{T,S}, x::StridedVector{T}) = BLAS.symv!(A.uplo, one(T), A.data, x, zero(T), y)
 A_mul_B!{T<:BlasComplex,S<:StridedMatrix}(y::StridedVector{T}, A::Hermitian{T,S}, x::StridedVector{T}) = BLAS.hemv!(A.uplo, one(T), A.data, x, zero(T), y)

test/linalg/symmetric.jl

@@ -209,3 +209,25 @@ let A = [1.0+im 2.0; 2.0 0.0]
     @test !ishermitian(A)
     @test_throws ArgumentError Hermitian(A)
 end
+
+# Tests for #10298
+for elT in [Float32, Float64]
+    A1 = ones(elT,10,10)
+    A = convert(Array{Complex{elT}}, A1)
+    A[end,1] = Complex{elT}(nextfloat((A[end,1].re)), A[end,1].im)
+    @test issym(A)
+    A = A + triu(A1)im - tril(A1)im
+    @test ishermitian(A)
+
+    D=diagm(rand(100))*100
+    D[end,1] = nextfloat((D[end,1]))
+
+    for HOrS in [Hermitian, Symmetric]
+        @test isposdef(HOrS(D))
+        @test isposdef!(copy(HOrS(D)))
+    end
+    @test ishermitian(D)
+    @test issym(D)
+    @test isposdef(D)
+end
+