Eigen returns non-eigenvectors for certain matrices

[moble]

I have a test that I'm trying to get working for BigFloats. It somewhat randomly constructs 5,000 symmetric real 4x4 matrices of BigFloats and extracts the dominant eigenvector of each. For 3 of these instances, I get results for which the answer is not an eigenvector at all (never mind with the given eigenvalue). There doesn't seem to be anything special about these matrices. I've extracted one to provide a concrete MWE. The eigenvalues of this matrix are expected to be (and reported as) [-1, -1, -1, 3], so I don't know that there would be any numerical problem here.

Here's the MWE:

using LinearAlgebra
using GenericLinearAlgebra
using StaticArrays

M = Symmetric(@SMatrix[
    big"-0.4080898675832881399369478084264191594976530854542904557798567397269356887436951";
        big"-0.1032324294981949906363774065395184125237581835226155628209100984396171211818558";
        big"-1.0795157507124452910839896877334667387210301781514938067860918240771876343947";
        big"0.9172086645212876240254394768180975107502376572771647296150618931226550446699544";;
    big"-0.1032324294981949906363774065395184125237581835226155628209100984396171211818558";
        big"-0.9819956883377066621250198846550622559246996804965712336465013506629992739010227";
        big"0.1882735697944729855991976669864503854920622386133987141371224931350749728226066";
        big"-0.1599663084136352437739757607131301560774255778371317602542426234968564801904052";;
    big"-1.0795157507124452910839896877334667387210301781514938067860918240771876343947";
        big"0.1882735697944729855991976669864503854920622386133987141371224931350749728226066";
        big"0.9688026817149176598146701814747478080649943014810992426739997593840858865727305";
        big"-1.672789745967021000172452940954243617442140494364475046869527486458478435262502";;
    big"0.9172086645212876240254394768180975107502376572771647296150618931226550446699544";
        big"-0.1599663084136352437739757607131301560774255778371317602542426234968564801904052";
        big"-1.672789745967021000172452940954243617442140494364475046869527486458478435262502";
        big"0.4212828742060771422472975116067336073573584644697624467523583310058490760719874"
])

eigensystem = eigen(M)

Here, the supposed dominant eigenvector is

julia> eigensystem.vectors[:, 4]
4-element SVector{4, BigFloat} with indices SOneTo(4):
  0.5293186034906468513583044045578991187176767353991885285160593337162140639833333
  0.6076028636690888706511331534102457699010874314100936245656011216542549302322659
 -0.5375888975041810576235384352977728036722978161994470931619050808663662918353576
  0.2482715314732502633344123259979948815899264503590173204924838651280614259806134

and its reported eigenvalue (3) is nonzero, so

(M * eigensystem.vectors[:, 4]) ./ eigensystem.vectors[:, 4]

should be a vector with entries that are at least approximately equal. Instead, I get

4-element SVector{4, BigFloat} with indices SOneTo(4):
  0.9999999999999999999999999999999999999999999999999999999999999999999999999999309
 -1.303869922808737255865818137316681872433951050681740949789081294567406034061292
  2.591452043646991636425616621765959618553871235916851252573571248344097579507887
  5.607430516593620629517859036040424284010328224975389025670154530756814862123736

It's interesting to note eigensystem.vectors[:, 1] is also not an eigenvector.

The dominant eigenvector (with eigenvalue 3) that I expect in this case is

x4 = [
    big"-0.3846784801677603137399569700744195708135504357934165621816525583943309062098062",
    big"0.06709007315224313025838456104648289261272108291672656267831069781563681483410404",
    big"0.7015701464776914343825844884840609580397271879944290132648084132383274894856047",
    big"-0.5960878446600964581372146249105424317644837852918997676276656289252753627271997",
]

and indeed I get the expected eigenvector behavior with the given matrix:

julia> (M * x4) ./ x4
4-element Vector{BigFloat}:
 3.0
 2.999999999999999999999999999999999999999999999999999999999999999999999999999862
 2.999999999999999999999999999999999999999999999999999999999999999999999999999965
 3.000000000000000000000000000000000000000000000000000000000000000000000000000035

---

julia> versioninfo()
Julia Version 1.9.0
Commit 8e630552924 (2023-05-07 11:25 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin22.4.0)
  CPU: 8 × Intel(R) Core(TM) i7-4770HQ CPU @ 2.20GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-14.0.6 (ORCJIT, haswell)
  Threads: 1 on 8 virtual cores

(jl_6G5ebb) pkg> status
Status `/private/var/folders/76/v4_xw_mx1ml49pbk_yd5pnrh0000gn/T/jl_6G5ebb/Project.toml`
  [14197337] GenericLinearAlgebra v0.3.10
  [90137ffa] StaticArrays v1.5.25
  [37e2e46d] LinearAlgebra

[andreasnoack]

Thanks for reporting this. In #124 I'm getting

julia> (M * eigensystem.vectors[:, 4]) ./ eigensystem.vectors[:, 4]
4-element SVector{4, BigFloat} with indices SOneTo(4):
 2.999999999999999999999999999999999999999999999999999999999999999999999999999896
 3.000000000000000000000000000000000000000000000000000000000000000000000000000691
 2.999999999999999999999999999999999999999999999999999999999999999999999999999931
 3.000000000000000000000000000000000000000000000000000000000000000000000000000069

[moble]

Awesome! Thanks very much!