scipy.sparse.linalg.

lsmr#

scipy.sparse.linalg.lsmr(A, b, damp=0.0, atol=1e-06, btol=1e-06, conlim=100000000.0, maxiter=None, show=False, x0=None)[source]#

用於最小平方問題的迭代求解器。

lsmr 求解線性方程式系統 Ax = b。如果系統不一致,它會求解最小平方問題 min ||b - Ax||_2A 是一個 m 乘 n 維度的矩陣,其中允許所有情況:m = n、m > n 或 m < n。b 是一個長度為 m 的向量。矩陣 A 可以是稠密的或稀疏的(通常是稀疏的)。

參數:
A{稀疏陣列, ndarray, LinearOperator}

線性系統中的矩陣 A。或者,A 可以是一個線性運算子,它可以使用例如 scipy.sparse.linalg.LinearOperator 來產生 AxA^H x

barray_like,形狀 (m,)

線性系統中的向量 b

dampfloat

正則化最小平方的阻尼因子。lsmr 求解正則化最小平方問題

min ||(b) - (  A   )x||
    ||(0)   (damp*I) ||_2

其中 damp 是一個純量。如果 damp 為 None 或 0,則系統在沒有正則化的情況下求解。預設值為 0。

atol, btolfloat,可選

停止容忍度。lsmr 會持續迭代,直到某個後向誤差估計值小於取決於 atol 和 btol 的某個量。令 r = b - Ax 為目前近似解 x 的殘差向量。如果 Ax = b 似乎是一致的,則當 norm(r) <= atol * norm(A) * norm(x) + btol * norm(b) 時,lsmr 終止。否則,當 norm(A^H r) <= atol * norm(A) * norm(r) 時,lsmr 終止。如果兩個容忍度均為 1.0e-6(預設值),則最終的 norm(r) 應精確到約 6 位數。(最終的 x 通常會具有較少的正確位數,具體取決於 cond(A) 和 LAMBDA 的大小。)如果 atolbtol 為 None,則將使用預設值 1.0e-6。理想情況下,它們應該分別是 Ab 條目中相對誤差的估計值。例如,如果 A 的條目具有 7 位正確數字,請設定 atol = 1e-7。這可以防止演算法執行超出輸入資料不確定性的不必要工作。

conlimfloat,可選

lsmr 如果 cond(A) 的估計值超過 conlim,則終止。對於相容系統 Ax = b,conlim 可能高達 1.0e+12(例如)。對於最小平方問題,conlim 應小於 1.0e+8。如果 conlim 為 None,則預設值為 1e+8。可以通過設定 atol = btol = conlim = 0 獲得最大精度,但迭代次數可能會過多。預設值為 1e8。

maxiterint,可選

lsmr 如果迭代次數達到 maxiter,則終止。預設值為 maxiter = min(m, n)。對於病態系統,可能需要更大的 maxiter 值。預設值為 False。

showbool,可選

如果 show=True,則印出迭代日誌。預設值為 False。

x0array_like,形狀 (n,),可選

x 的初始猜測值,如果為 None,則使用零。預設值為 None。

在版本 1.0.0 中新增。

回傳值:
xfloat 的 ndarray

回傳最小平方解。

istopint

istop 給出停止的原因

istop   = 0 means x=0 is a solution.  If x0 was given, then x=x0 is a
            solution.
        = 1 means x is an approximate solution to A@x = B,
            according to atol and btol.
        = 2 means x approximately solves the least-squares problem
            according to atol.
        = 3 means COND(A) seems to be greater than CONLIM.
        = 4 is the same as 1 with atol = btol = eps (machine
            precision)
        = 5 is the same as 2 with atol = eps.
        = 6 is the same as 3 with CONLIM = 1/eps.
        = 7 means ITN reached maxiter before the other stopping
            conditions were satisfied.
itnint

使用的迭代次數。

normrfloat

norm(b-Ax)

normarfloat

norm(A^H (b - Ax))

normafloat

norm(A)

condafloat

A 的條件數。

normxfloat

norm(x)

註解

在版本 0.11.0 中新增。

參考文獻

[1]

D. C.-L. Fong and M. A. Saunders, “LSMR: An iterative algorithm for sparse least-squares problems”, SIAM J. Sci. Comput., vol. 33, pp. 2950-2971, 2011. arXiv:1006.0758

範例

>>> import numpy as np
>>> from scipy.sparse import csc_array
>>> from scipy.sparse.linalg import lsmr
>>> A = csc_array([[1., 0.], [1., 1.], [0., 1.]], dtype=float)

第一個範例具有平凡解 [0, 0]

>>> b = np.array([0., 0., 0.], dtype=float)
>>> x, istop, itn, normr = lsmr(A, b)[:4]
>>> istop
0
>>> x
array([0., 0.])

回傳的停止代碼 istop=0 表示找到一個零向量作為解。回傳的解 x 確實包含 [0., 0.]。下一個範例具有非平凡解

>>> b = np.array([1., 0., -1.], dtype=float)
>>> x, istop, itn, normr = lsmr(A, b)[:4]
>>> istop
1
>>> x
array([ 1., -1.])
>>> itn
1
>>> normr
4.440892098500627e-16

istop=1 所示,lsmr 找到一個符合容忍度限制的解。給定的解 [1., -1.] 顯然求解了方程式。剩餘的回傳值包含關於迭代次數 (itn=1) 以及已求解方程式左右兩側剩餘差異的資訊。最後一個範例示範了在方程式無解的情況下的行為

>>> b = np.array([1., 0.01, -1.], dtype=float)
>>> x, istop, itn, normr = lsmr(A, b)[:4]
>>> istop
2
>>> x
array([ 1.00333333, -0.99666667])
>>> A.dot(x)-b
array([ 0.00333333, -0.00333333,  0.00333333])
>>> normr
0.005773502691896255

istop 表示系統不一致,因此 x 更像是對應最小平方問題的近似解。normr 包含找到的最小距離。