Python Interface¶
-
Primme.
eigsh
(A, k=6, M=None, sigma=None, which='LM', v0=None, ncv=None, maxiter=None, tol=0, return_eigenvectors=True, Minv=None, OPinv=None, mode='normal', ortho=None, return_stats=False, maxBlockSize=0, minRestartSize=0, maxPrevRetain=0, method=None, return_history=False, **kargs)¶ Find k eigenvalues and eigenvectors of the real symmetric square matrix or complex Hermitian matrix A.
Solves
A * x[i] = w[i] * x[i]
, the standard eigenvalue problem for w[i] eigenvalues with corresponding eigenvectors x[i].If M is specified, solves
A * x[i] = w[i] * M * x[i]
, the generalized eigenvalue problem for w[i] eigenvalues with corresponding eigenvectors x[i]Parameters: - A (An N x N matrix, array, sparse matrix, or LinearOperator) -- the operation A * x, where A is a real symmetric matrix or complex Hermitian.
- k (int, optional) -- The number of eigenvalues and eigenvectors to be computed. Must be 1 <= k < min(A.shape).
- M (An N x N matrix, array, sparse matrix, or LinearOperator) --
(not supported yet) 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. For best results, the data type of M should be the same as that of A.
- sigma (real, optional) -- Find eigenvalues near sigma.
- v0 (N x i, ndarray, optional) -- Initial guesses to the eigenvectors.
- ncv (int, optional) -- The maximum size of the basis
- which (str ['LM' | 'SM' | 'LA' | 'SA']) --
Which k eigenvectors and eigenvalues to find:
'LM' : Largest in magnitude eigenvalues; the farthest from sigma'SM' : Smallest in magnitude eigenvalues; the closest to sigma
'LA' : Largest algebraic eigenvalues
'SA' : Smallest algebraic eigenvalues
'CLT' : closest but left to sigma
'CGT' : closest but greater than sigma
When sigma == None, 'LM', 'SM', 'CLT', and 'CGT' treat sigma as zero.
- maxiter (int, optional) -- Maximum number of iterations.
- tol (float) -- Required accuracy for eigenpairs (stopping criterion). The default value is sqrt of machine precision.
- Minv ((not supported yet)) -- The inverse of M in the generalized eigenproblem.
- OPinv (N x N matrix, array, sparse matrix, or LinearOperator, optional) -- Preconditioner to accelerate the convergence. Usually it is an approximation of the inverse of (A - sigma*M).
- return_eigenvectors (bool, optional) -- Return eigenvectors (True) in addition to eigenvalues
- mode (string ['normal' | 'buckling' | 'cayley']) -- Only 'normal' mode is supported.
- ortho (N x i, ndarray, optional) -- Seek the eigenvectors orthogonal to these ones. The provided vectors should be orthonormal. Useful to avoid converging to previously computed solutions.
- maxBlockSize (int, optional) -- Maximum number of vectors added at every iteration.
- minRestartSize (int, optional) -- Number of approximate eigenvectors kept during restart.
- maxPrevRetain (int, optional) -- Number of approximate eigenvectors kept from previous iteration in restart. Also referred as +k vectors in GD+k.
- method (int, optional) --
Preset method, one of:
- DEFAULT_MIN_TIME : a variant of JDQMR,
- DEFAULT_MIN_MATVECS : GD+k
- DYNAMIC : choose dynamically between these previous methods.
See a detailed description of the methods and other possible values in [2].
- return_stats (bool, optional) -- If True, the function returns extra information (see stats in Returns).
- return_history (bool, optional) -- If True, the function returns performance information at every iteration (see hist in Returns).
Returns: - w (array) -- Array of k eigenvalues ordered to best satisfy "which".
- v (array) -- An array representing the k eigenvectors. The column
v[:, i]
is the eigenvector corresponding to the eigenvaluew[i]
. - stats (dict, optional (if return_stats)) -- Extra information reported by PRIMME:
- "numOuterIterations": number of outer iterations
- "numRestarts": number of restarts
- "numMatvecs": number of A*v
- "numPreconds": number of OPinv*v
- "elapsedTime": time that took
- "estimateMinEVal": the leftmost Ritz value seen
- "estimateMaxEVal": the rightmost Ritz value seen
- "estimateLargestSVal": the largest singular value seen
- "rnorms" : ||A*x[i] - x[i]*w[i]||
- "hist" : (if return_history) report at every outer iteration of:
- "elapsedTime": time spent up to now
- "numMatvecs": number of A*v spent up to now
- "nconv": number of converged pair
- "eval": eigenvalue of the first unconverged pair
- "resNorm": residual norm of the first unconverged pair
Raises: PrimmeError -- When the requested convergence is not obtained. The PRIMME error code can be found as
err
attribute of the exception object.See also
scipy.sparse.linalg.eigs()
- eigenvalues and eigenvectors for a general (nonsymmetric) matrix A
Primme.svds()
- singular value decomposition for a matrix A
Notes
This function is a wrapper to PRIMME functions to find the eigenvalues and eigenvectors [1].
References
[1] PRIMME Software, https://github.com/primme/primme [2] Preset Methods, http://www.cs.wm.edu/~andreas/software/doc/readme.html#preset-methods [3] A. Stathopoulos and J. R. McCombs PRIMME: PReconditioned Iterative MultiMethod Eigensolver: Methods and software description, ACM Transaction on Mathematical Software Vol. 37, No. 2, (2010), 21:1-21:30. Examples
>>> import Primme, scipy.sparse >>> A = scipy.sparse.spdiags(range(100), [0], 100, 100) # sparse diag. matrix >>> evals, evecs = Primme.eigsh(A, 3, tol=1e-6, which='LA') >>> evals # the three largest eigenvalues of A array([ 99., 98., 97.]) >>> new_evals, new_evecs = Primme.eigsh(A, 3, tol=1e-6, which='LA', ortho=evecs) >>> new_evals # the next three largest eigenvalues array([ 96., 95., 94.])