ZGGEV fails but Matlab works

[CyLith]

For a pair of tridiagonal matrices obtained from a finite element discretization, computing the generalized eigenvalues with ZGGEV fails with info > 0, but using eig in Matlab succeeds. I would like to know what Matlab is doing differently. I have attached the matrices in text form output from C++.

Edit: I'm also attaching a Matlab script to create the matrices for convenience.

[Julien Langou]

Hi Victor,

Thanks for the bug report.

I can reproduce a nonconvergence of ZGGEV for the pencil you provided. This confirms what you are saying. Nice pencil, thanks . . .

Your pencil is made of two complex symmetric tridiagonal matrix. (I do mean complex symmetric, I do not mean Hermitian.) I guess you knew that.

Results:

1) Matlab is working (as you mentioned).

2) The error INFO is 481 on my computer using the svn-up-to-date LAPACK and associated reference BLAS,

3) The error INFO is 454 (so a little better) on my computer using the svn-up-to-date LAPACK and vecLib BLAS,

4) The error INFO is 454 (so a little better) on my computer using vecLib LAPACK+BLAS

=> bug 116 is confirmed.

Victor, what is your error INFO? (Just for fun . . . )

MathWorks is using a Nonsymmetric Generalized Eigensolver which is based on LAPACK, but is improved. This is Cleve Moler's legacy I guess there. We are in contact with MathWorks, and we are trying to synch our version with theirs. We did this a few times ago but I guess we did not fully do it since there seem to still be some problems in LAPACK which are not in MathWorks version. Shoo . . . Well, we will need to look in to this.

Cheers,
Julien.

[jean-christophe]

Hi everyone,

I have developed a small python code to perform some linear stability analysis of the Navier-Stokes equations. For that purpose, I have a generalized eigenvalue problem to solve, which I do using scipy.linalg.eig(A,B). This python function relies on the zggev subroutine from lapack. Unfortunately, for some set of parameters (which are surprisingly machine-dependent), the generalized eigensolver actually fails. In order to make sure this was not an issue with my matrices (or python), I have saved them and processed them with matlab which runs normally. I have stored the matrices (256x256) in a matlab .mat file. Please let me know if someone wants to take a look at them and I will send them to whoever is interested.

Any chance someone would have an idea on how to work this around?

Thanks a lot,
JC