eig#
- scipy.linalg.eig(a, b=None, left=False, right=True, overwrite_a=False, overwrite_b=False, check_finite=True, homogeneous_eigvals=False)[source]#
Solve an ordinary or generalized eigenvalue problem of a square matrix.
Find eigenvalues w and right or left eigenvectors of a general matrix:
a vr[:,i] = w[i] b vr[:,i] a.H vl[:,i] = w[i].conj() b.H vl[:,i]
where
.H
is the Hermitian conjugation.The documentation is written assuming array arguments are of specified “core” shapes. However, array argument(s) of this function may have additional “batch” dimensions prepended to the core shape. In this case, the array is treated as a batch of lower-dimensional slices; see Batched Linear Operations for details.
- Parameters:
- a(M, M) array_like
A complex or real matrix whose eigenvalues and eigenvectors will be computed.
- b(M, M) array_like, optional
Right-hand side matrix in a generalized eigenvalue problem. Default is None, identity matrix is assumed.
- leftbool, optional
Whether to calculate and return left eigenvectors. Default is False.
- rightbool, optional
Whether to calculate and return right eigenvectors. Default is True.
- overwrite_abool, optional
Whether to overwrite a; may improve performance. Default is False.
- overwrite_bbool, optional
Whether to overwrite b; may improve performance. Default is False.
- 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.
- homogeneous_eigvalsbool, optional
If True, return the eigenvalues in homogeneous coordinates. In this case
w
is a (2, M) array so that:w[1,i] a vr[:,i] = w[0,i] b vr[:,i]
Default is False.
- Returns:
- w(M,) or (2, M) double or complex ndarray
The eigenvalues, each repeated according to its multiplicity. The shape is (M,) unless
homogeneous_eigvals=True
.- vl(M, M) double or complex ndarray
The left eigenvector corresponding to the eigenvalue
w[i]
is the columnvl[:,i]
. Only returned ifleft=True
. The left eigenvector is not normalized.- vr(M, M) double or complex ndarray
The normalized right eigenvector corresponding to the eigenvalue
w[i]
is the columnvr[:,i]
. Only returned ifright=True
.
- Raises:
- LinAlgError
If eigenvalue computation does not converge.
See also
eigvals
eigenvalues of general arrays
eigh
Eigenvalues and right eigenvectors for symmetric/Hermitian arrays.
eig_banded
eigenvalues and right eigenvectors for symmetric/Hermitian band matrices
eigh_tridiagonal
eigenvalues and right eiegenvectors for symmetric/Hermitian tridiagonal matrices
Examples
>>> import numpy as np >>> from scipy import linalg >>> a = np.array([[0., -1.], [1., 0.]]) >>> linalg.eigvals(a) array([0.+1.j, 0.-1.j])
>>> b = np.array([[0., 1.], [1., 1.]]) >>> linalg.eigvals(a, b) array([ 1.+0.j, -1.+0.j])
>>> a = np.array([[3., 0., 0.], [0., 8., 0.], [0., 0., 7.]]) >>> linalg.eigvals(a, homogeneous_eigvals=True) array([[3.+0.j, 8.+0.j, 7.+0.j], [1.+0.j, 1.+0.j, 1.+0.j]])
>>> a = np.array([[0., -1.], [1., 0.]]) >>> linalg.eigvals(a) == linalg.eig(a)[0] array([ True, True]) >>> linalg.eig(a, left=True, right=False)[1] # normalized left eigenvector array([[-0.70710678+0.j , -0.70710678-0.j ], [-0. +0.70710678j, -0. -0.70710678j]]) >>> linalg.eig(a, left=False, right=True)[1] # normalized right eigenvector array([[0.70710678+0.j , 0.70710678-0.j ], [0. -0.70710678j, 0. +0.70710678j]])