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)