JOVANA
Library Glossary Getting Started Three Levels Fields How it works Mission
Join the mission
All guides

向解迭代逼近:從分裂法到 Krylov

當矩陣大到無法分解時,你只能猜測再改進。本篇從簡單的定常迭代——Jacobi 與 Gauss-Seidel——建立直覺,然後揭示每一種現代求解器背後的結構性思想:Krylov 子空間。

迭代的交易

直接法變換 A 直到答案落出;迭代法從一個猜測 x_0 出發並精化它,產生收斂到真實 x 的 x_1, x_2, ...。對於巨型稀疏系統這筆交易很誘人:每次迭代的主要成本是一次矩陣-向量乘積 A v,對稀疏 A 而言僅需 O(nnz)——與非零元個數成正比,而非 n^2。你從不構造 L 或 U,從不承受填充。你只透過乘法接觸 A。

定常迭代:Jacobi 與 Gauss-Seidel

最古老的思想是矩陣分裂:寫 A = M - N,其中 M 易求逆,把 A x = b 化為 M x = N x + b,並迭代 x_{k+1} = M^-1 (N x_k + b)。在 Jacobi 迭代中,M 就是 A 的對角部分,於是每個新分量僅用舊值求解。Gauss-Seidel取 M 為下三角部分,在同一遍掃描內立即複用剛更新的分量——通常以約一半的步數收斂。

Jacobi sweep for  A x = b  (a_ii on the diagonal):

  for i = 1..n:
      x_new[i] = ( b[i] - sum_{j != i} a[i][j] * x_old[j] ) / a[i][i]
  x_old <- x_new      # all components use OLD values

Gauss-Seidel sweep (in place):

  for i = 1..n:
      x[i] = ( b[i] - sum_{j < i} a[i][j]*x[j]      # already updated
                    - sum_{j > i} a[i][j]*x[j]       # still old
             ) / a[i][i]

Convergence guaranteed if A is strictly diagonally dominant
(|a_ii| > sum_{j != i} |a_ij| for every row).
Jacobi 每遍只用舊值;Gauss-Seidel 立即複用新值,通常使迭代次數減半。

飛躍:Krylov 子空間

定常方法浪費資訊:在第 k 步它只記住 x_k,丟棄了之前的每一個向量。現代洞見是把它們全部保留。從殘差 r_0 = b - A x_0 出發,反覆乘以 A 生成 r_0, A r_0, A^2 r_0, ...,它們的張成就是 Krylov 子空間 K_m = span{r_0, A r_0, ..., A^{m-1} r_0}。

為何這是正確的空間?因為你想要的答案 A^-1 b 可以寫成 A 的多項式作用於 b(Cayley-Hamilton 保證了這一點),而該多項式中 A 的各次冪恰好落在 Krylov 子空間內。於是 K_m 包含著對真解越來越好的逼近,且它僅用你已經負擔得起的稀疏矩陣-向量乘積建構。一個 Krylov 方法按某種準則在 K_m 內挑出最佳逼近——而準則的選擇正是下一篇中區分共軛梯度GMRES 的關鍵。