Merge pull request #12236 from mfasi/logm_test_redirect

Fix warnings due to logm in tests/linalg3.jl

Author: tkelman

base/linalg/dense.jl

@@ -286,9 +286,14 @@ function logm(A::StridedMatrix)
 
     # Use Schur decomposition
     n = chksquare(A)
-    S,Q,d = schur(complex(A))
-    R = logm(UpperTriangular(S))
-    retmat = Q * R * Q'
+    if istriu(A)
+        retmat = full(logm(UpperTriangular(complex(A))))
+        d = diag(A)
+    else
+        S,Q,d = schur(complex(A))
+        R = logm(UpperTriangular(S))
+        retmat = Q * R * Q'
+    end
 
     # Check whether the matrix has nonpositive real eigs
     np_real_eigs = false

base/linalg/diagonal.jl

@@ -119,6 +119,7 @@ end
 eye{T}(::Type{Diagonal{T}}, n::Int) = Diagonal(ones(T,n))
 
 expm(D::Diagonal) = Diagonal(exp(D.diag))
+logm(D::Diagonal) = Diagonal(log(D.diag))
 sqrtm(D::Diagonal) = Diagonal(sqrt(D.diag))
 
 #Linear solver

base/linalg/triangular.jl

@@ -891,11 +891,14 @@ end
 
 # Complex matrix logarithm for the upper triangular factor, see:
 #   Al-Mohy and Higham, "Improved inverse  scaling and squaring algorithms for
-# the matrix logarithm", SIAM J. Sci. Comput., 34(4), (2012), pp. C153-C169.
+#     the matrix logarithm", SIAM J. Sci. Comput., 34(4), (2012), pp. C153-C169.
 #   Al-Mohy, Higham and Relton, "Computing the Frechet derivative of the matrix
-# logarithm and estimating the condition number", SIAM J. Sci. Comput., 35(4),
-# (2013), C394-C410.
-#   http://eprints.ma.man.ac.uk/1687/02/logm.zip
+#     logarithm and estimating the condition number", SIAM J. Sci. Comput., 35(4),
+#     (2013), C394-C410.
+#
+# Based on the code available at http://eprints.ma.man.ac.uk/1851/02/logm.zip,
+# Copyright (c) 2011, Awad H. Al-Mohy and Nicholas J. Higham
+# Julia version relicensed with permission from original authors
 function logm{T<:Union{Float64,Complex{Float64}}}(A0::UpperTriangular{T})
     maxsqrt = 100
     theta = [1.586970738772063e-005,
@@ -907,7 +910,7 @@ function logm{T<:Union{Float64,Complex{Float64}}}(A0::UpperTriangular{T})
          2.879093714241194e-001]
     tmax = size(theta, 1)
     n = size(A0, 1)
-    A = A0
+    A = copy(A0)
     p = 0
     m = 0
 
@@ -916,7 +919,7 @@ function logm{T<:Union{Float64,Complex{Float64}}}(A0::UpperTriangular{T})
     dm1 = Array(T, n)
     s = 0
     for i = 1:n
-        dm1[i] = dm1[i] - 1.
+        dm1[i] = d[i] - 1.
     end
     while norm(dm1, Inf) > theta[tmax]
         for i = 1:n
@@ -1075,6 +1078,7 @@ function logm{T<:Union{Float64,Complex{Float64}}}(A0::UpperTriangular{T})
     return UpperTriangular(Y)
 
 end
+logm(A::LowerTriangular) = logm(A.').'
 
 function sqrtm{T}(A::UpperTriangular{T})
     n = size(A, 1)

test/linalg3.jl

@@ -162,20 +162,29 @@ for elty in (Float32, Float64, Complex64, Complex128)
 end
 
 for elty in (Float64, Complex{Float64})
+
     A4  = convert(Matrix{elty}, [1/2 1/3 1/4 1/5+eps();
                                  1/3 1/4 1/5 1/6;
                                  1/4 1/5 1/6 1/7;
                                  1/5 1/6 1/7 1/8])
     @test_approx_eq expm(logm(A4)) A4
 
+    OLD_STDERR = STDERR
+    rd,wr = redirect_stderr()
     A5  = convert(Matrix{elty}, [1 1 0 1; 0 1 1 0; 0 0 1 1; 1 0 0 1])
     @test_approx_eq expm(logm(A5)) A5
+    s = readline(rd)
+    @test contains(s, "WARNING: Matrix with nonpositive real eigenvalues, a nonprimary matrix logarithm will be returned.")
 
     A6  = convert(Matrix{elty}, [-5 2 0 0 ; 1/2 -7 3 0; 0 1/3 -9 4; 0 0 1/4 -11])
     @test_approx_eq expm(logm(A6)) A6
+    s = readline(rd)
+    @test contains(s, "WARNING: Matrix with nonpositive real eigenvalues, a nonprimary matrix logarithm will be returned.")
+    redirect_stderr(OLD_STDERR)
 
     A7  = convert(Matrix{elty}, [1 0 0 1e-8; 0 1 0 0; 0 0 1 0; 0 0 0 1])
     @test_approx_eq expm(logm(A7)) A7
+
 end
 
 A8 = 100 * [-1+1im 0 0 1e-8; 0 1 0 0; 0 0 1 0; 0 0 0 1];