# A misunderstanding or a bug in LAPACK's solver for generalized eigenvalue problems?

In my application, I have two general real matrices $$A$$,$$B$$ defined as follows, $$A=\begin{bmatrix} -s I_3 & A_0 & 0 & 0 \\ A_0^T & -s I_3 & 0 & 0 \\ 0 & A_1 & -s I_3 & A_0 \\ A_1^T & 0 & A_0^T & -s I_3 \end{bmatrix}, B=\begin{bmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ I_3 & 0 & 0 & 0 \\ 0 & I_3& 0 & 0 \end{bmatrix},$$ where $$A_0,A_1 \in R^{3\times3}$$ are general matrices, $$s$$ is a real scalar, and $$I_3 \in R^{3\times3}$$ is the identity matrix. I would like to solve the eigenvalue problem $$AX= B X D$$ numerically, where $$X$$ is the eigenvector matrix and $$D$$ is the eigenvalue matrix. I only want the right (not left) eigenvalues and eigenvectors. Consider the following as a concrete example, where $$A_0=\begin{bmatrix} 1 & & \\ & 2 & \\ & & 3 \end{bmatrix}, A_1=\begin{bmatrix} 0 & 1 & 0\\ 0 & 0 & 0\\ 0 & 0 & 0 \end{bmatrix}, s=1,$$ From my application, I know there exists at least one eigenvalue $$d_1$$ and one eigenvector $$x_1$$, which are $$x_1=\left(0.625543,0,0,0.625543, 0, 0, 0, -0.417029, 0, 0, -0.208514, 0\right), d_1=0,$$ because $$A x_1 -d_1 B x_1 = \textbf{0},\left\|x_1\right\|=1$$. However, no matter which software I used, I couldn't get the correct eigenvalues, even for $$d_1$$. The software I have tried are Matlab, Eigen(c++ library) and LAPACK.

1. The eigenvalues from Matlab are infinity.
2. The right eigenvalues from Eigen and LAPACK are represented by $$d_i=\alpha_i / \beta_i$$, where $$\alpha_i$$ is a complex number and $$\beta_i$$ is a real number. The output $$\beta_i$$ is zero for all eigenvalues. In other words, all right eigenvalues are infinity.

I'm not sure which part is wrong. Is it that the software cannot compute a correct eigenvalue or that $$d_1$$ cannot be considered as an eigenvalue?

Crosspost my own question: A misunderstanding or a bug in LAPACK's solver for generalized eigenvalue problems?, as it may be a better fit here.

• The $x_1$ you provide does not appear to be a valid eigenvector. For instance, $A x_1 = \left(0,0,0,0,0,0,0,10^{-6},0,0,-10^{-6},0\right)^T$. Therefore $Ax_1 - d_1 B x_1$ is not 0.
– vibe
Oct 7, 2020 at 17:04
• @vibe In my opinion, your result of $Ax_1$ can be considered as zero because it is sufficiently small given the fact that I only provide 6 digits for each entry of the eigenvector. I can give you $x_1$ with 15 digits, $x_1=(0.625543242171224,0,0,0.625543242171224,0,0,0,-0.417028828114149,0,0,-0.208514414057075,0)^T$. Oct 7, 2020 at 19:15
• You can check in Matlab that $\left\|A x_1\right\|$\approx10^{-15}$. Oct 7, 2020 at 19:16 • Moreover, you can also find that det(A) is zero. There are two zero singular values of$A$. In other words, there are non-zero nullspace vectors for$A$. Any nullspace vectors$v$can make$Av=0$, so (v, 0) is a valid solution to the eigenvalue problem. Oct 7, 2020 at 19:19 • I computed the generalized Schur decomposition (using GSL):$A = QSZ^T, B=QTZ^T$. The$\alpha$values come from the diagonal elements of$S$, while the$\beta$values are the diagonal elements of$T$. I found that all of the diagonal elements of$T$are zero, while a few (but not all) of the diagonal elements of$S$are zero. According to Golub and Van Loan (4th ed), theorem 7.7.1, if$T_{kk}$and$S_{kk}$are both zero for any$k$, then the set of eigenvalues$\lambda(A,B) = \mathbb{C}\$.
– vibe
Oct 8, 2020 at 0:38

I suspect the root of your trouble is what has been detected in the comments by Vibe: For any number $$\omega\in \mathbb{K}$$ (with $$\mathbb{K}= \mathbb{R}$$ or $$\mathbb{C}$$) you can find $$\boldsymbol{X}$$ such that $$AX = \omega BX$$ (with $$A$$ and $$B$$ taken in your concrete example).

You have already decomposed the problem in 4 blocks of 3 variables. Then let us define $$X = (X_1, X_2, X_3, X_4)^T$$ where $$X_i\in \mathbb{K}^3$$.

The first 6 lines of your system yield $$X_1 = X_2 = (a,0,0)$$, with $$a\in\mathbb{K}$$ an arbitrary number.

With this in mind, it becomes obvious upon inspection that the last 6 lines of the system are of the form $$C (X_3, X_4)^T = R(a,\omega)$$ where $$C$$ is a $$6\times6$$ non singular matrix (essentially the bottom right blocks in $$A$$), and $$R\in\mathbb{K}^6$$ is a right-hand side that depends on $$a,\omega\in\mathbb{K^2}$$. Since $$C$$ is invertible, you can always find solutions to this non homogeneous problem.

In conclusion, for all $$\omega$$, you can invert this system and compute $$X\ne 0$$ such that $$AX=\omega BX$$

The other answers already tell you what went wrong, but I will add a terminology note: the term for what is happening is that the pencil $$A - \lambda B$$ is a singular matrix pencil, i.e., $$\det (A - \lambda B)$$ is identically equal to zero. So there are no generalized eigenvalues (or, at least, they cannot be defined as usual as the roots of the generalized characteristic polynomial).

An useful tool to analyze singular pencils is the Kronecker canonical form, which is a generalization of the Jordan canonical form to pencils and can include singular (rectangular) blocks.

Note that small perturbations of a singular pencil can give regular (=non-singular) pencils with any eigenvalue in $$\mathbb{C} \cup \{\infty\}$$, so numerically computing those eigenvalues is an ill-posed problem, and in practice you can expect just about any value to show up among the computed eigenvalues.