diff --git a/src/symmetriceigen.jl b/src/symmetriceigen.jl index 1d92e4f8..7a01cd8c 100644 --- a/src/symmetriceigen.jl +++ b/src/symmetriceigen.jl @@ -331,6 +331,12 @@ function eigvals!(A::Hermitian{T,S}, B::Hermitian{T,S}; sortby::Union{Function,N end eigvecs(A::HermOrSym) = eigvecs(eigen(A)) +function eigvecs(A::RealHermSymComplexHerm, eigvals::AbstractVector{<:Real}) + F = hessenberg(A) # transform to SymTridiagonal form + X = eigvecs(F.H, eigvals) + return F.Q * X # transform eigvecs of F.H back to eigvecs of A +end + function eigvals(A::AbstractMatrix, C::Cholesky; sortby::Union{Function,Nothing}=nothing) if ishermitian(A) eigvals!(eigencopy_oftype(Hermitian(A), eigtype(eltype(A))), C; sortby) diff --git a/src/tridiag.jl b/src/tridiag.jl index b45e6bc2..31209310 100644 --- a/src/tridiag.jl +++ b/src/tridiag.jl @@ -311,7 +311,7 @@ eigmin(A::SymTridiagonal) = eigvals(A, 1:1)[1] #Compute selected eigenvectors only corresponding to particular eigenvalues """ - eigvecs(A::SymTridiagonal[, eigvals]) -> Matrix + eigvecs(A::Union{Symmetric, Hermitian, SymTridiagonal}[, eigvals])::Matrix Return a matrix `M` whose columns are the eigenvectors of `A`. (The `k`th eigenvector can be obtained from the slice `M[:, k]`.) @@ -346,7 +346,7 @@ julia> eigvecs(A, [1.]) -0.5547001962252291 ``` """ -eigvecs(A::SymTridiagonal{<:BlasFloat,<:StridedVector}, eigvals::Vector{<:Real}) = LAPACK.stein!(A.dv, A.ev, eigvals) +eigvecs(A::SymTridiagonal{<:BlasFloat,<:StridedVector}, eigvals::StridedVector{<:Real}) = LAPACK.stein!(A.dv, A.ev, eigvals) function svdvals!(A::SymTridiagonal) vals = eigvals!(A) diff --git a/test/symmetriceigen.jl b/test/symmetriceigen.jl index 7a2b88a5..2150457a 100644 --- a/test/symmetriceigen.jl +++ b/test/symmetriceigen.jl @@ -199,4 +199,23 @@ end end end +@testset "eigvecs for specific eigvals" begin + function testeigvecs(S, vals) + V = eigvecs(S, vals) + @test S * V ≈ V * Diagonal(vals) + end + for T in (Symmetric, Hermitian) + S = T(rand(3,3)) + vals = eigvals(S) + testeigvecs(S, vals) + testeigvecs(S, vals[1:2]) + testeigvecs(S, @view vals[2:2]) + end + H = Hermitian(rand(ComplexF64,3,3)) + vals = eigvals(H) + testeigvecs(H, vals) + testeigvecs(H, vals[1:2]) + testeigvecs(H, @view vals[2:2]) +end + end # module TestSymmetricEigen