Python SciPy linalg.eigsh用法及代码示例

本文简要介绍 python 语言中 scipy.sparse.linalg.eigsh 的用法。

用法:

scipy.sparse.linalg.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')#

求实对称方阵或复 Hermitian 矩阵 A 的 k 个特征值和特征向量。

解决 A @ x[i] = w[i] * x[i] ，即 w[i] 特征值与相应特征向量 x[i] 的标准特征值问题。

如果指定了 M，则解决 A @ x[i] = w[i] * M @ x[i] ，即 w[i] 个具有相应特征向量 x[i] 的特征值的广义特征值问题。

请注意，当 A 是复 Hermitian 矩阵时，没有专门的例程。在这种情况下，eigsh() 将调用eigs() 并返回由此获得的特征值的实部。

参数 ：：

A： ndarray，稀疏矩阵或LinearOperator

表示操作 A @ x 的平方运算符，其中 A 是实对称或复 Hermitian。对于屈曲模式(见下文)A 必须另外是正定的。

k： 整数，可选

所需的特征值和特征向量的数量。 k 必须小于 N。不可能计算矩阵的所有特征向量。

返回 ：：

w： 数组

k 个特征值的数组。

v： 数组

一个数组表示k特征向量。专栏v[:, i]是特征值对应的特征向量w[i].

其他参数 ：：

M： 一个 N x N 矩阵、数组、稀疏矩阵或线性运算符表示

广义特征值问题的操作M @ x

A @ x = w * M @ x.

如果 A 是实数，则 M 必须表示一个实对称矩阵，如果 A 是复数，则 M 必须表示一个复数 Hermitian 矩阵。为获得最佳结果，M 的数据类型应与 A 的数据类型相同。另外：

If sigma is None, M is symmetric positive definite.

If sigma is specified, M is symmetric positive semi-definite.

In buckling mode, M is symmetric indefinite.

如果 sigma 为无，eigsh 需要一个运算符来计算线性方程 M @ x = b 的解。这是通过对显式矩阵 M 的(稀疏)LU 分解或通过对一般线性算子的迭代求解器在内部完成的。或者，用户可以提供矩阵或运算符 Minv，它给出 x = Minv @ b = M^-1 @ b 。

sigma： 真实的

使用 shift-invert 模式查找 sigma 附近的特征值。这需要运算符计算线性系统 [A - sigma * M] x = b 的解，其中 M 是单位矩阵(如果未指定)。这是通过显式矩阵 A 和 M 的(稀疏)LU 分解在内部计算的，或者如果 A 或 M 是一般线性算子，则通过迭代求解器计算。或者，用户可以提供矩阵或运算符 OPinv，它给出 x = OPinv @ b = [A - sigma * M]^-1 @ b 。请注意，当指定 sigma 时，关键字 ‘which’ 指的是移位的特征值 w'[i]，其中：

if mode == ‘normal’, w'[i] = 1 / (w[i] - sigma).

if mode == ‘cayley’, w'[i] = (w[i] + sigma) / (w[i] - sigma).

if mode == ‘buckling’, w'[i] = w[i] / (w[i] - sigma).

(见下文‘mode’ 中的进一步讨论)

v0： ndarray，可选

迭代的起始向量。默认值：随机

ncv： 整数，可选

生成的Lanczos向量个数ncv必须大于k且小于n；建议使用 ncv > 2*k 。默认值：min(n, max(2*k + 1, 20))

which： str ['LM' | 'SM' | '洛杉矶' | 'SA' | '是']

如果 A 是复 Hermitian 矩阵，则“BE”无效。要找到哪些 k 个特征向量和特征值：

‘LM’ : Largest (in magnitude) eigenvalues.

‘SM’ : Smallest (in magnitude) eigenvalues.

‘LA’ : Largest (algebraic) eigenvalues.

‘SA’ : Smallest (algebraic) eigenvalues.

‘BE’ : Half (k/2) from each end of the spectrum.

当 k 为奇数时，从高端返回一个 (k/2+1)。当 sigma != None 时，‘which’ 指的是移位的特征值 w'[i](参见上文 ‘sigma’ 中的讨论)。 ARPACK 通常比小值更擅长查找大值。如果需要小特征值，请考虑使用shift-invert 模式以获得更好的性能。

maxiter： 整数，可选

允许的最大 Arnoldi 更新迭代次数。默认值：n*10

tol： 浮点数

特征值的相对精度(停止标准)。默认值 0 表示机器精度。

Minv： N x N 矩阵、数组、稀疏矩阵或LinearOperator

见上文 M 中的注释。

OPinv： N x N 矩阵、数组、稀疏矩阵或LinearOperator

请参见上文 sigma 中的注释。

return_eigenvectors： bool

除了特征值之外，还返回特征向量 (True)。此值确定特征值的排序顺序。排序顺序也取决于哪一个多变的。

For which = ‘LM’ or ‘SA’:

If return_eigenvectors is True, eigenvalues are sorted by algebraic value.

If return_eigenvectors is False, eigenvalues are sorted by absolute value.

For which = ‘BE’ or ‘LA’:

eigenvalues are always sorted by algebraic value.

For which = ‘SM’:

If return_eigenvectors is True, eigenvalues are sorted by algebraic value.

If return_eigenvectors is False, eigenvalues are sorted by decreasing absolute value.

mode： 字符串 [‘normal’ | ‘buckling’ | ‘cayley’]

指定用于shift-invert 模式的策略。此参数仅适用于实值 A 和 sigma != None。对于 shift-invert 模式，ARPACK 在内部解决了特征值问题 OP @ x'[i] = w'[i] * B @ x'[i] 并将生成的 Ritz 向量 x'[i] 和 Ritz 值 w'[i] 转换为问题 A @ x[i] = w[i] * M @ x[i] 的所需特征向量和特征值。模式如下：

‘normal’ :

OP = [A - sigma * M]^-1 @ M, B = M, w’[i] = 1 / (w[i] - sigma)

‘buckling’ :

OP = [A - sigma * M]^-1 @ A, B = A, w’[i] = w[i] / (w[i] - sigma)

‘cayley’ :

OP = [A - sigma * M]^-1 @ [A + sigma * M], B = M, w’[i] = (w[i] + sigma) / (w[i] - sigma)

模式的选择会影响关键字‘which’选择的特征值，也会影响收敛的稳定性(参见[2]的讨论)。

抛出 ：：

ArpackNoConvergence

当未获得要求的收敛时。

当前收敛的特征值和特征向量可以在异常对象的eigenvalues和eigenvectors属性中找到。

注意：

该函数是 ARPACK [1] SSEUPD 和 DSEUPD 函数的包装，它们使用隐式重启 Lanczos 方法来查找特征值和特征向量 [2]。

参考：

[1]

ARPACK 软件，https://github.com/opencollab/arpack-ng

[2]

R. B. Lehoucq、D. C. Sorensen 和 C. Yang，ARPACK 用户指南：通过隐式重启 Arnoldi 方法解决大规模特征值问题。暹罗，宾夕法尼亚州费城，1998 年。

例子：

>>> import numpy as np
>>> from scipy.sparse.linalg import eigsh
>>> identity = np.eye(13)
>>> eigenvalues, eigenvectors = eigsh(identity, k=6)
>>> eigenvalues
array([1., 1., 1., 1., 1., 1.])
>>> eigenvectors.shape
(13, 6)