From a1a0a324a546f1eb42e6ec42bb65695acb80f06d Mon Sep 17 00:00:00 2001 From: "John F. Gibson" Date: Tue, 7 Sep 2021 09:48:25 -0400 Subject: [PATCH] Add wrapper function eigen(Matrix A) where A is hermitian by value but not by type --- src/eigenSelfAdjoint.jl | 7 ++++++- test/eigengeneral.jl | 10 +++++++++- 2 files changed, 15 insertions(+), 2 deletions(-) diff --git a/src/eigenSelfAdjoint.jl b/src/eigenSelfAdjoint.jl index 0b6dd3c..52f3970 100644 --- a/src/eigenSelfAdjoint.jl +++ b/src/eigenSelfAdjoint.jl @@ -593,7 +593,6 @@ LinearAlgebra.eigen!(A::Hermitian; tol = eps(real(eltype(A))), debug = false) = _eigen!(A, tol = tol, debug = debug) - function eigen2!(A::SymmetricTridiagonalFactorization; tol = eps(real(float(one(eltype(A))))), debug = false) @@ -649,6 +648,12 @@ function LinearAlgebra.eigen(A::Hermitian) T = typeof(sqrt(zero(eltype(A)))) return eigen!(LinearAlgebra.copy_oftype(A, T)) end +function LinearAlgebra.eigen(A::Matrix) + if !ishermitian(A) + throw(ArgumentError("eigen not implement for non-Hermitian matrices of generic type (e.g. Matrix{BigFloat})")) + end + return eigen(Hermitian(A)) +end # Aux (should go somewhere else at some point) function LinearAlgebra.givensAlgorithm(f::Real, g::Real) diff --git a/test/eigengeneral.jl b/test/eigengeneral.jl index 91ae90b..f8d0613 100644 --- a/test/eigengeneral.jl +++ b/test/eigengeneral.jl @@ -150,4 +150,12 @@ end @test sort(imag(vals)) ≈ sort(imag(λs)) atol=1e-25 end -end \ No newline at end of file +@testset "Wrapper function for eigen(Matrix A) where A is hermitian by value but not by type" begin + # just need to verify that the wrapper function correctly calls eigen(Hermitian(A)) + # so no need for multiple numeric tests + a = big(2.0); + A = [1 a; a 0] # A == A^*, but type is Matrix{BigFloat}, not Hermitian{BigFloat, Matrix{BigFloat}} + λs = [(1 - √(1+4a^2))/2 ; (1 + √(1+4a^2))/2] + vals = eigvals(A) + @test sort(vals) ≈ sort(λs) atol=1e-25 +end