# [Scipy-tickets] [SciPy] #1456: scipy.sparse.linalg.lobpcg fails for large sparse matrices

SciPy Trac scipy-tickets@scipy....
Wed Jun 8 09:52:24 CDT 2011

```#1456: scipy.sparse.linalg.lobpcg fails for large sparse matrices
---------------------------------+------------------------------------------
Reporter:  AndrewMacLean        |       Owner:  wnbell
Type:  defect               |      Status:  new
Priority:  normal               |   Milestone:
Component:  scipy.sparse.linalg  |     Version:  0.7.2
Keywords:                       |
---------------------------------+------------------------------------------
Hi,

I have been using ''lobpcg'' (scipy.sparse.linalg.lobpcg) to solve
symmetric generalized eigenvalue problems with large, sparse stiffness and
mass matrices, say '''A''' and '''B'''. The problem is of the form
'''Av''' =lambda'''BV'''.

They are both symmetric, and '''B''' is positive definite, and both are of
size 70 000 x 70 000. '''A''' has about 5 million non-zero values and
'''B''' has about 1.6 million. '''A''' has condition number on the order
of 10^9^ and '''B''' has a condition number on the order of 10^6^. I have
stored them both as ''csc'' type sparse matrices from the ''scipy.sparse''
library.

The part of my code using lobpcg to find the 20 smallest eigenvalues:

{{{
from scipy.sparse.linalg import lobpcg
from scipy import rand

X = rand(A.shape[0], 20)

W, V = lobpcg (A, X, B = B, largest = False, maxiter = 40)
}}}

I tested lobpcg on a "scaled down" version of my problem, with '''A''' and
'''B''' matrices of size 10 000 x 10 000, and got the correct values
(using maxiter = 20), as confirmed by using the ''eigs'' function in
Matlab. I used it here to find the smallest 10 eigenvalues, and here is a
table of my results, showing that the eigenvalues computed from lobpcg in
Python are very close to those using eigs in Matlab:

With full sized '''A''' and '''B''' matrices, I could not get the correct
values, and it became clear that increasing the maximum number of
iterations indefinitely would not work (see graph below). I made a graph
for the 20th smallest eigenvalue vs. the number of iterations. I compared
4 different guesses for X, 3 of which were just random matrices (as in the
code above), and a 4th orthonormalized one.

It appears that it will take a very large number of iterations to get the
correct eigenvalues.  As well, I tested lobpcg by using eigenvectors
generated by ''eigs'' in Matlab for X, and ''lobpcg'' returned the correct
values.

I am wondering if anyone else has encountered problems using matrices of
this size with ''lobpcg''.

Thanks for any suggestions on this,
Andrew

--
Ticket URL: <http://projects.scipy.org/scipy/ticket/1456>
SciPy <http://www.scipy.org>
SciPy is open-source software for mathematics, science, and engineering.
```