LAPACK fails calculating eigenvectors

Question

I wrote some code to solve the general eigenvalue problem and now I am comparing my results against LAPACK's DSPGVX function. I just worked with this example.

So I obtained the 4 auto vectors

{
 {-0.0319133, -0.265466, -0.713483,  0.64765},
 {-0.425628,  -0.520961, -0.714215,  0.193227},
 { 0.32702,    0.565845, -0.37129,  -0.659561},
 {-0.682699,  -0.056645,  0.0771025, 0.724409}
}

and auto values

{-2.22545, 1.12704, -0.454756, 0.100076}

both with my code and with Mathematica and results agree.

But in the previous link, auto vectors reported from LAPACK are completely different.

 Eigenvalues
    -0.4548  0.1001
 Selected eigenvectors
          1       2
 1   0.3080  0.4469
 2   0.5329  0.0371
 3  -0.3496 -0.0505
 4  -0.6211 -0.4743

Whom should I trust?

P.S. I also checked that my auto values/autovectors are correct since they yield A*x-lambda*B*x=0, while the values from LAPACK do not.

score 2 · Answer 1 · answered Mar 13 '12 at 10:25

2

I don't know why you think that LAPACK is giving incorrect answers, they seem fine to me. Using the four figure decimals that you quote I get residuals (r = A*x - lambda*B*x) such that

norm(r1) = 1.5921e-04, norm(r2) = 6.0842e-05.

Since norm(A) = 1.2994 and norm(B) = 7.9874, these residuals seem very satisfactory.

The eigenvectors produced by DSPGVX are normalized so that

norm(x'*B*x) = 1.

answered Mar 13 '12 at 10:25

Sven Hammarling

21
1

I get residues like 1e-15 with my results. If theirs are normalized, they should just notice it. – Open the way Mar 13 '12 at 10:39
thanks for your answer. My main concern is; since my solution gives a smalle residue, why yours does not, since this is the one reported in the website Otherwise I assume I can take any of them as good, do you agree? – Open the way Mar 14 '12 at 09:23

score 1 · Answer 2 · answered Mar 12 '12 at 18:14

It looks like Lapack's results do actually correspond to the last two eigenvalues generated by your code and Mathematica, though with lower-order bits coming out quite different. The corresponding vectors are quite close, just scaled differently.

Clearly, if your/Mathematica's values check out and Lapack's do not, you should trust the one that produces the correct answers. It might be valuable to investigate what it is about your problem and and Lapack's algorithms that make it provide very imprecise answers.

score 1 · Accepted Answer · answered Mar 12 '12 at 18:29

It looks like DSGPVX is solving A*lambda = B*x*lambda; Matlab gives the DSGPVX solution to your problem using "eig", though Matlab's documentation is correct. My guess is this a bug in the DSGPVX documentation.

>> a=[0.24 0.39 0.42 -0.16;0.39 -0.11 0.79 0.63;0.42 0.79 -0.25 0.48;-0.16 0.63 0.48 -0.03];
>> b=[4.16 -3.12 0.56 -0.1;-3.12 5.03 -0.83 1.09;0.56 -0.83 0.76 0.34;-0.1 1.09 0.34 1.18];
>> [v,d]=eig(a,b)

v =

   -0.0690    0.3080   -0.4469   -0.5528
   -0.5740    0.5329   -0.0371   -0.6766
   -1.5428   -0.3496    0.0505   -0.9276
    1.4004   -0.6211    0.4743    0.2510


d =

   -2.2254         0         0         0
         0   -0.4548         0         0
         0         0    0.1001         0
         0         0         0    1.1270

>> norm(a*v-b*v*d)

ans =

   1.5001e-15

LAPACK fails calculating eigenvectors

3 Answers3