scipy.linalg.eigvalsh¶
-
scipy.linalg.
eigvalsh
(a, b=None, lower=True, overwrite_a=False, overwrite_b=False, turbo=True, eigvals=None, type=1, check_finite=True)[source]¶ Solve an ordinary or generalized eigenvalue problem for a complex Hermitian or real symmetric matrix.
Find eigenvalues w of matrix a, where b is positive definite:
a v[:,i] = w[i] b v[:,i] v[i,:].conj() a v[:,i] = w[i] v[i,:].conj() b v[:,i] = 1
- Parameters
- a(M, M) array_like
A complex Hermitian or real symmetric matrix whose eigenvalues and eigenvectors will be computed.
- b(M, M) array_like, optional
A complex Hermitian or real symmetric definite positive matrix in. If omitted, identity matrix is assumed.
- lowerbool, optional
Whether the pertinent array data is taken from the lower or upper triangle of a. (Default: lower)
- turbobool, optional
Use divide and conquer algorithm (faster but expensive in memory, only for generalized eigenvalue problem and if eigvals=None)
- eigvalstuple (lo, hi), optional
Indexes of the smallest and largest (in ascending order) eigenvalues and corresponding eigenvectors to be returned: 0 <= lo < hi <= M-1. If omitted, all eigenvalues and eigenvectors are returned.
- typeint, optional
Specifies the problem type to be solved:
type = 1: a v[:,i] = w[i] b v[:,i]
type = 2: a b v[:,i] = w[i] v[:,i]
type = 3: b a v[:,i] = w[i] v[:,i]
- overwrite_abool, optional
Whether to overwrite data in a (may improve performance)
- overwrite_bbool, optional
Whether to overwrite data in b (may improve performance)
- check_finitebool, optional
Whether to check that the input matrices contain only finite numbers. Disabling may give a performance gain, but may result in problems (crashes, non-termination) if the inputs do contain infinities or NaNs.
- Returns
- w(N,) float ndarray
The N (1<=N<=M) selected eigenvalues, in ascending order, each repeated according to its multiplicity.
- Raises
- LinAlgError
If eigenvalue computation does not converge, an error occurred, or b matrix is not definite positive. Note that if input matrices are not symmetric or hermitian, no error is reported but results will be wrong.
See also
eigh
eigenvalues and right eigenvectors for symmetric/Hermitian arrays
eigvals
eigenvalues of general arrays
eigvals_banded
eigenvalues for symmetric/Hermitian band matrices
eigvalsh_tridiagonal
eigenvalues of symmetric/Hermitian tridiagonal matrices
Notes
This function does not check the input array for being hermitian/symmetric in order to allow for representing arrays with only their upper/lower triangular parts.
Examples
>>> from scipy.linalg import eigvalsh >>> A = np.array([[6, 3, 1, 5], [3, 0, 5, 1], [1, 5, 6, 2], [5, 1, 2, 2]]) >>> w = eigvalsh(A) >>> w array([-3.74637491, -0.76263923, 6.08502336, 12.42399079])