How "chol" and "qz" MATLAB algorithms are utilised in the "eig" MATLAB function?

The MATLAB function "eig" when used this way: [V,D] = eig(A,B), gives us the eigenvalues and eigenvectors of an eigenvalue problem A*V = B*V*D, where D is the diagonal matrix of eigenvalues and matrix V columns are the corresponding right eigenvectors.

In the documentation, it is menioned that the function "eig" uses the algorithms "chol" (Cholesky factorization), and "qz" (QZ factorization for generalized eigenvalues / generalized Schur decomposition).

However, when reading about the methods of Cholesky factorization and QZ factorization for generalized eigenvalues, the first method only decomposes a matrix into a product of where L is a lower triangular matrix, and the second method computes the eigenvalues and what about the eigenvectors?

So, for this reason, I don't understand how using only these algorithms could return us both eigenvectos and eigenvalues of an eigenvalue problem.

Answer 1

[Bruno Luong]

See stage 3 described in https://netlib.org/lapack/lug/node56.html or somewhere in this file

[Zahraa]

Thank you! It shows how LAPACK fortran library uses generalized Schur decomposition method.

But LAPACK isn't in MATLAB, right? My concern was how the algorithms are used inside eig MATLAB function.

[Bruno Luong]

MATLAB uses Lapack in many linear algebra for dense matrix

https://www.mathworks.com/help/matlab/math/lapack-in-matlab.html

I beleive they develop more for sparse matrix.

[Zahraa]

So, LAPACK is included in the "eig" function?

[Bruno Luong]

Yes

[Zahraa]

Thank you. Thats right. MATLAB uses LAPACK in some linear algebra functions such as eig . In an answer here, it is mentioned that MATLAB function eig checks for symmetry, and then uses an appropriate call to LAPACK.

In a documentation, it is mentioned that: The primary tool for the solution of generalized matrix eigenvalue problems, Ax = λBx, is the function qz(A,B). MATLAB 5.3 has only one version of the QZ algorithm. It always uses complex arithmetic and produces complex results, even if the input matrices are real. LAPACK provides a real QZ algorithm for real matrices, so we can now have qz(A,B,'real') and qz(A,B,'complex').

Back to eig documentation, when it is said that eig uses the algorithms of the functions "chol" or "qz", why it doesn't also mention that it uses other functions in LAPACK, because in the stage 3 of the link you sent, the eigenvectors are computed using xTGEVC and xTGEVC solvers of LAPACK. So, "qz" does give us the eigenvalues, but alone doesn't give eigenvectors. And "chol" only decomposes a matrix as I wrote in my question.

There must be information missing in the section of algorithms of eig function of its documentaion?

Then, also:

In a documentation, it is mentioned that: The extensive list of functions now available with LAPACK means that MATLAB's space saving general-purpose codes can be replaced by faster, more focused routines. There are now 16 different code paths underlying the eig function, depending on whether there are one or two arguments, whether the arguments are real or complex, whether the problem is symmetric and whether the eigenvectors are requested.

This made me think, what do they mean with eig function having "16 different code paths".

And just a little more, here it is mentioned that "chol" uses some of the LAPACK routines.

[Bruno Luong]

MATLAB documentation mentions only qz and chol because it has the user switch of eig. It doesn't mention other stage because it judges most users do not care about knowing such thing. ou can consider it as missing if you want.

I don't know exactly how MATLAB uses chol.

All I know is it can transform generalized eigen decomposition of sdp matrix (which has specific algorithm) with Cholesky decomposition:

% Generate random sdp matrices A and B
A=rand(5);
A=A'*A;
B=rand(5);
B=B'*B;

eig(A,B)

ans = 5×1

0.0622 0.1110 1.2352 3.8498 7.6685

BB=chol(B); C=BB'\(A/BB); 
C = (C + C')/2; % spd C
eig(C)

ans = 5×1

0.0622 0.1110 1.2352 3.8498 7.6685

Answer 2

[zhangbaby]

Mathematically, for any symmetric matrix, after taking the Cholesky decomposition , we can further write the triangular matrix as where D is a diagonal matrix and its diagonal is that of R. A simple proof provides the proof that we can get the eigenvalues and eigenvectors from this decomposition.

Considering coding, I should say it is far more complicated. The early version MATLAB document directly provides the one-to-one mapping from the LAPACK function to MATLAB function, which you can see from Which algorithm do DGGEV or DSYGV Eigen solvers in LAPACK implement. For symmetric matrix they use ?syev and for non-symmetric matrix they use ?ggev (they are not using ?tgevc).