[SciPy-User] scipy.sparse.linalg.eigen error code?
Pauli Virtanen
pav@iki...
Wed Oct 13 13:37:09 CDT 2010
Wed, 13 Oct 2010 15:55:42 +0200, Nico Schlömer wrote:
> I'm using scipy.sparse.linalg.eigen to compute the lowest magnitude
> eigenvalue of a matrix, and I noticed that sometimes the code would
> return 0+0j where I really didn't expect it. Turns out that tweaking the
> number of iterations results in some more meaningful value here,
> suggesting that the underlying Arnoldi (ARPACK?) iteration failed with
> the maxiter value previously given. As far as I can see, there's no way
> to tell that the iteration actually failed.
>
> Is that correct?
Yep. It seems to try to raise a warning, though,
elif self.info == -1:
warnings.warn("Maximum number of iterations taken: %s" %
self.iparam[2])
But this is bad practice. You'll only see the warning *once*, and there
is indeed no way to catch it.
The Arpack interface seems to need some fixing, too. Exceeding the number
of iterations is really an error condition. The interface also apparently
doesn't check that ncv is sensible, so in some cases (small matrices) it
can give "Error -3".
As a work-around, you can probably use
def eigen_fixed(A, k=6, M=None, sigma=None, which='LM', v0=None,
ncv=None, maxiter=None, tol=0,
return_eigenvectors=True):
import numpy as np
from scipy.sparse.linalg.interface import aslinearoperator
from scipy.sparse.linalg.eigen.arpack.arpack import _UnsymmetricArpackParams
A = aslinearoperator(A)
if A.shape[0] != A.shape[1]:
raise ValueError('expected square matrix (shape=%s)' % (A.shape,))
n = A.shape[0]
matvec = lambda x : A.matvec(x)
params = _UnsymmetricArpackParams(n, k, A.dtype.char, matvec, sigma,
ncv, v0, maxiter, which, tol)
if M is not None:
raise NotImplementedError("generalized eigenproblem not supported yet")
while not params.converged:
params.iterate()
if params.info == -1:
raise RuntimeError("Did not converge")
return params.extract(return_eigenvectors)
--
Pauli Virtanen
More information about the SciPy-User
mailing list