Why EIGS is not able to reproduce the same result as EIG for a generalised eigenvalue problem?

Question

I have an generalised eigenvalue problem:

A**q = alpha**B*q

where A and B are complex square matrices. Using eig I can easily find the correct answer for a given A and B matrices:

from scipy import linalg
alpha, q = linalg.eig(A,B)

And if I plot the result eigenvalues alpha in a complex plan where alpha_r is the real part of the eigenvalue and alpha_i is the imaginary part of the eigenvalue:

Now I would like to solve the same problem but using `eigs` sparse algorithms from Scipy and something weird happens:

from scipy.linalg.sparse import eigs
alpha, q = eigs(A,k=40,B,sigma=1.0+0.0J)

I know that's probably something related with the characteristics of the matrices A and B but I don't know exactly what. Using eigs is suppose to be more efficient specially if the matrices are sparse.

Here you can find the A and B data, eigsolution and eigssolution:

https://drive.google.com/drive/folders/1p_3AZzRtFfF64P7Sh5tDEA-I-voRtfAi?usp=share_link

To load the data:

import numpy as np

data  = np.load('YourPath/matrix.npz')
A     = data['A']
B     = data['B']
data  = np.load('YourPath/eig_result.npz')
alpha = data['alpha']
q     = data['q']
data  = np.load('YourPath/eigs_result.npz')
alpha = data['alpha']
q     = data['q']

Answer 1

At first, I didn't think there was anything strange here. After correcting the syntax and import errors:

from scipy.sparse.linalg import eigs
alpha, q = eigs(A, 40, B, sigma=1 + 0j)

alpha contains 40 eigenvalues near the specified sigma. It's not immediately obvious from your plot, but upon closer inspection, these don't appear to be eigenvalues returned by scipy.linalg.eig. In fact, when I specify a different sigma, say 4+4j, alpha is returned with 40 values near sigma - but these are obviously not present in your plot.

I agree that this is not desirable behavior, but it is probably because your M matrix (B) does not satisfy the requirements specified in the documentation.

M: An array, sparse matrix, or LinearOperator representing the operation M@x for the generalized eigenvalue problem

A @ x = w * M @ x.

M must represent a real symmetric matrix if A is real, and must represent a complex Hermitian matrix if A is complex.

A is complex, but your B matrix is not Hermitian; all the nonzero diagonal components are 0 - 1j.