Hello,

The class GeneralizedEigenSolver doesn't have the eigenvectors() method implemented. I noticed that the the class GeneralizedSelfAdjointEigenSolver has the method, but the matrices in my case are not selfadjoint.
Is there another way to use the eigen library to obtain the generalized eigenvectors if the generalized eigenvalues are available?

Hi,

You'll see my question is similar: viewtopic.php?f=74&t=117277

The relationship is A=PDP(-1) as stated at http://mathworld.wolfram.com/EigenDecomposition.html
Now I haven't come across how to get P(column eigenvectors) from known A and D where D is the the diagonal matrix of known eigenvalues. The solvers solve for both P and D and thus make no use of known eigenvalues.

The only way I know is taking the nullspace for each eigenvalue (A-lambda*I) and going to row reduced echelon form and pick the associated eigenvector colum; build up P one column at a time. Google eigenvectors from row reduced echelon form and eigenvector Gaussian elimination on the augmented matrix to get the algorithm. Now this might not be efficient for large matrices. Sparse may get some speed if there is an efficient way to do the reduction to echelon form. Though many applications are a mix of sparse and band matrices.

I'm using the unsupported DMRES with sparse and it makes sense it doesn't do eigenvectors along the way. Also I need only the first half of the eigenvalues and eigevectors.

The other aspect is dealing with multiplicity/degeneracy complicates the "Gaussian elimination on the augmented matrix" method but is important to involve since underneath krylov space solves for them. What complicates is floating point errors and many systems may appear to be full rank but are really gettting hit by the deterministic effects of binary computers of limited precision. Might need quadmath levels of precision for many applications. We need ternary machines!

Thank you very much for the tip, I'm trying to figure it out. I assume you haven't found any code that could do this? I hoped I wouldn't have to go through writing the code myself.
Dan

Haven't found specific code yet, no. Might have to spend some time with linear algebra journals to see if there are efficient ways to impliment. Or a mathematician around here knows what to do. Also find if there is a way to skip directly to P*P.transpose() from A and D which is really what is needed for rebuilding the next iteration.

The details of degeneracy complicate for me. When the GeneralizedSelfAdjointEigenSolver.info() says NoConvergencehttp://eigen.tuxfamily.org/dox-devel/group__enums.html#ga51bc1ac16f26ebe51eae1abb77bd037b then I go to DGMRES http://eigen.tuxfamily.org/dox/unsupported/classEigen_1_1DGMRES.html solver which gets the eigenvalues and now I'm learning what to do to handle the degeneracy/multiplicity [some call it repeated eigenvalues] for the eigenvectors per eigenvalue where they exist.

It may be effected by the precision limit for a given problem's order of magnitude. For double floating point b.bbbbbxxxxEbbb the computer's filling in of the xxxx with 'deterministic' http://en.wikipedia.org/wiki/Non-deterministic_Turing_machine digits can wreck havoc with lower precision eigenvalues and rank. I've got more test cases to run to see how my particular application is doing with this. Something to be mindful of.

I've been taking a look at "$eigenroot\Eigen\src\Eigenvalues\GeneralizedEigenSolver.h". It seems the generalized eigensolver comes-up with the eigenvectors as well.
It's doing some QZ decomposition, I haven't looked into the mathematical details, but it computes the eigenvecs implicitly.
In the class there is a variable called "m_eivec".
Inside the function ::Compute there is "m_eivec = m_realQZ.matrixZ().transpose();" so I'm guessing the eigenvectors should be around there somewhere.
I will print out the values to see if something is going on with them.
Maybe if you want you can give it some attention as well.

Dan

Ah just missing the public getter. In a unit test I used "#define protected public" above the includes and tested getting at the m_eivec you found and it seems populated with eigenvectors. If you have a known benchmark this can be verified

Code: Select all

Eigen::GeneralizedEigenSolver<Eigen::MatrixXd> ges;
Eigen::MatrixXd A = Eigen::MatrixXd::Random(4,4);
Eigen::MatrixXd B = Eigen::MatrixXd::Random(4,4);
ges.compute(A, B, true);
std::cout << "The (complex) numerators of the generalized eigenvalues are: " << ges.alphas().transpose() << std::endl;
std::cout << "The (real) denominators of the generalized eigenvalues are: " << ges.betas().transpose() << std::endl;
std::cout << "The (complex) generalized eigenvalues are (alphas./beta): " << ges.eigenvalues().transpose() << std::endl;
std::cout << "The (complex) generalized eigenvectors are (alphas./beta): " << ges.m_eivec << std::endl;


yields up:

Suite: gpu_tests
Test: testGeneral ...A:
0.680375 0.823295 -0.444451 -0.270431
-0.211234 -0.604897 0.10794 0.0268018
0.566198 -0.329554 -0.0452059 0.904459
0.59688 0.536459 0.257742 0.83239
B:
0.271423 -0.967399 -0.686642 0.997849
0.434594 -0.514226 -0.198111 -0.563486
-0.716795 -0.725537 -0.740419 0.0258648
0.213938 0.608354 -0.782382 0.678224
The (complex) numerators of the generalized eigenvalues are: (-0.650112,0.800296) (-0.650112,-0.800296) (-0.39761,0) (-1.1157,0)
The (real) denominators of the generalized eigenvalues are: 1.10971 1.10971 -1.24986 0.745646
The (complex) generalized eigenvalues are (alphas./beta): (-0.585838,0.721174) (-0.585838,-0.721174) (0.318124,-0) (-1.49628,0)
The (complex) generalized eigenvectors are (alphas./beta):
0.88967 0.362989 -0.268673 -0.0673911
0.0847702 0.268359 0.782561 -0.555334
0.438861 -0.66006 0.461525 0.398392
-0.0932864 0.600453 0.320012 0.726874

Then find eiegenvectors method and put it into GeneralizedEigenSolver(.h) and return the protected m_eivec. Should get the generalized eigenvectors going

Has anyone verified that the values in GeneralizedEigenSolver.m_eivec are indeed the correct eigenvectors corresponding to the calculated eigenvalues? I have been looking at these with a specific solution in mind, but am not getting results I expect. Maybe the results stored in m_eivec are not complete? I don't have the math background to know for sure, but seems like the method GeneralizedEigenSolver::eigenvectors() needs more code added to take what's stored in m_eivec and generate true eigenvectors. Any help or insights on this would be greatly appreciated.

Rus