4.1 Jacobi 與 Gauss-Seidel 法#

概念說明#

把矩陣 \(A\) 拆成 \(A = D + R\) ，其中 \(D\) 是對角部分、 \(R\) 是其餘。則 \(A\mathbf{x} = \mathbf{b}\) 等價於

這正是不動點形式 \(\mathbf{x} = g(\mathbf{x})\) ，迭代之即得 Jacobi 法：

Gauss-Seidel 法只改一件事：算出第 \(i\) 個新分量後，立刻用它更新同一輪後續的計算（用 \(\mathbf{x}^{(k+1)}\) 的已算部分取代舊值）。這個小改動通常讓收斂快約一倍，且省去另一份儲存。

收斂條件由壓縮映像原理給出。Jacobi 的迭代矩陣為 \(T_J = -D^{-1}R\) ，Gauss-Seidel 為 \(T_{GS} = -(D + L)^{-1}U\) （ \(A = D + L + U\) ）。定理：**迭代收斂 \(\iff\) 迭代矩陣的譜半徑 \(\rho(T) < 1\) **。兩個常用充分條件：

• \(A\) 嚴格對角佔優（每列對角元素絕對值大於其餘之和）→ Jacobi 與 Gauss-Seidel 皆收斂，且 \(\rho(T_{GS}) \le \rho(T_J) < 1\) 。
• \(A\) 對稱正定 → Gauss-Seidel 收斂（Jacobi 不一定）。

直觀理解：每輪迭代是「用鄰居的現況猜自己」，對角佔優表示「自己受自己的影響最大」，猜錯不會被鄰居放大——系統自然平息。這也是為什麼偏微分方程離散化（五點差分）產生的系統特別適合這類方法。

生活類比#

一群人各自猜一個數，規則是「我的猜測 = 大家平均值的某種修正」。每輪大家同時公布新猜測是 Jacobi（同步更新）；誰算出來就馬上被別人使用是 Gauss-Seidel（非同步更新）。只要規則「以自己為主」，謠言會逐步平息、猜測趨於一致。

演算法步驟#

Gauss-Seidel 法：

1. 給初值 \(\mathbf{x}^{(0)}\) （常用 0 或前次解的熱啟動）。
2. 對 \(i = 1, \dots, n\) ： \(x_i \leftarrow \dfrac{1}{a_{ii}}\left(b_i - \sum_{j<i} a_{ij}x_j^{(\text{new})} - \sum_{j>i} a_{ij}x_j^{(\text{old})}\right)\) 。
3. 計算殘差 \(\|\mathbf{b} - A\mathbf{x}\|\) ；小於容許值則停止。
4. 回到步驟 2。

Python 範例#

import numpy as np

def jacobi(A, b, x0=None, tol=1e-10, max_iter=1000):
    x = np.zeros_like(b, dtype=float) if x0 is None else x0.copy()
    D = np.diag(A).astype(float)
    R = A - np.diag(D)
    for k in range(1, max_iter + 1):
        x_new = (b - R @ x) / D
        if np.linalg.norm(x_new - x, np.inf) < tol:
            return x_new, k
        x = x_new
    return x, max_iter

def gauss_seidel(A, b, tol=1e-10, max_iter=1000):
    n = len(b)
    x = np.zeros(n)
    for k in range(1, max_iter + 1):
        x_old = x.copy()
        for i in range(n):
            s = A[i] @ x - A[i, i] * x[i]   # 用最新 x
            x[i] = (b[i] - s) / A[i, i]
        if np.linalg.norm(x - x_old, np.inf) < tol:
            return x, k
    return x, k

A = np.array([[10., -1, 2], [-1, 11, -1], [2, -1, 10]])
b = np.array([6., 25, -11])
xj, nj = jacobi(A, b)
xg, ng = gauss_seidel(A, b)
print(f"Jacobi      : {xj}, 迭代 {nj} 次")
print(f"Gauss-Seidel: {xg}, 迭代 {ng} 次")
print(f"驗證誤差: {np.linalg.norm(A @ xg - b):.2e}")

此矩陣嚴格對角佔優，兩法都收斂；實測 Gauss-Seidel 迭代次數約為 Jacobi 的一半。

複雜度與收斂性分析#

每輪 Jacobi 成本為一次矩陣-向量乘， \(O(\text{nnz})\) （非零元素數）；Gauss-Seidel 相同但順序存取快取較友善。記憶體：Jacobi 需兩份解向量，Gauss-Seidel 一份即可。收斂速度由譜半徑決定，誤差 \(\|\mathbf{e}_k\| \lesssim \rho^k \|\mathbf{e}_0\|\) （線性收斂）。弱點：對「接近奇異」或「對角不佔優」的系統 \(\rho\) 接近甚至超過 1，迭代極慢或發散。改善手段：矩陣重排、對角縮放（前處理）、或直接升級到 SOR 與共軛梯度。

練習#

1. 對矩陣 \(\begin{pmatrix} 2 & 1 \\ 5 & 3 \end{pmatrix}\) 分別寫出 Jacobi 與 Gauss-Seidel 的迭代矩陣，算譜半徑，判斷誰收斂誰發散。
2. 把五點差分拉普拉斯系統（ \(n = 50\) 的網格）用 Gauss-Seidel 解到殘差 \(10^{-8}\) ，記錄迭代次數；改變 \(n\) 觀察迭代數如何增長。
3. 實作「逐列檢查對角佔優」的函數，若不佔優則對列做重排（若可能）後再解，並用 4.1 的例子測試。