C++ linear algebra library armadillo : how to use eig_pair to get the same result as eig function in Matlab?

StackOverflow https://stackoverflow.com/questions/22634638

Question

I try to use eig_pair to get Eigen decomposition for pair of general dense square matrices A and B of the same size, such that A*eigvec = B*eigvec*diagmat(eigval), but the result doesn't match the Matlab function eig. for example:

A= [1,2;3,4]     B=[2,4;5,8]

in Matlab:

[u,v] = eig(A,B)

result:

u =

 -1.0000   -0.0000

 0.5000   -1.0000

v =

    1.0000         0
         0    0.5000

in armadillo:

eig_pair(v,u,A,B)

result: u:

9.9301e-016     -1.0000e+000

1.0000e+000      5.0000e-001

v:

0.5000

1.0000

My question is: how to get the values of u and v that match the results in Matlab? Looking forwards to your reply!!!

Was it helpful?

Solution

Eigenvectors are not unique. If u is an eigenvector, so is m * u for all m != 0. Furthermore, the order that eig returns eigenvectors in Matlab is arbitrary. (I don't know what order Armadillo returns eigenvectors.) You could try and create a canonical order for the eigenvectors by sorting the eigenvalues, but that is problematic if you have complex eigenvalues. (Recall that real matrices can have complex eigenvalues.)

Thus, (-1.0000, 0.5000) (first column of u in Matlab) is the same eigenvector as ( -1.0000e+000, 5.0000e-001) (second column of u in Armadillo). Similarly, (-0.0000, -1.0000) is equivalent to (9.9301e-016, 1.0000e+000) when you scale by -1 and account for floating point errors. Note that there may be numerical precision errors which would cause the floating point values to compare not equal even if mathematically the numbers are equal.

If you want a canonical representation of eigenvectors, you could rescale them to have norm 1, and also multiply by -1 if the sign of the first element is negative. Of course, if the first element in the eigenvector is close to 0, this is again problematic since the value might have ended up just barely on the wrong side of zero due to numerical reasons. So come to think of it, it might be better to ensure that the largest element (after normalization)--rather than the first--is positive.

Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top