3. 迭代求解器家族概览

求解大规模稀疏线性方程组 Ax = b 时，直接法（如 LU 分解）在三维问题中往往因填充元爆炸而不可行，迭代法成为主流选择。迭代求解器可粗分为两大阵营：定常迭代（stationary iteration）与非定常 Krylov 子空间方法。理解两者的定位差异，是选型的第一步。

定常迭代：以分裂矩阵为核心的简单迭代

定常迭代的基本思想是把系数矩阵分裂为 A = M − N，其中 M 易于求逆，然后构造迭代格式 x^(k+1) = M⁻¹N x^(k) + M⁻¹b。迭代矩阵 G = M⁻¹N 在整个迭代过程中保持不变，这正是“定常”的含义。三种经典方法对应不同的分裂方式：

• Jacobi：取 M = D（对角部分），每次用上一轮的全部旧值更新所有分量，天然可并行，但收敛最慢，且仅对对角占优等较窄矩阵类保证收敛。
• Gauss-Seidel：取 M = D + L（下三角部分），新算出的分量立即参与后续计算，收敛速度通常约为 Jacobi 的两倍，但存在数据依赖，并行性差。
• SOR：在 Gauss-Seidel 基础上引入松弛因子 ω，格式为 x_i^(k+1) = (1−ω)x_i^(k) + ω·(Gauss-Seidel 更新值)。当 0 < ω < 2 时可能收敛，ω = 1 退化为 Gauss-Seidel。对 Poisson 型矩阵，最优 ω 可使谱半径从 1−O(h²) 改善到 1−O(h)，迭代次数从 O(h⁻²) 降到 O(h⁻¹)。

定常迭代的优点是每步计算量小、内存占用低、实现简单，且每次迭代都在单调地“磨平”误差的高频分量。但其收敛速度强烈依赖矩阵谱性质，对对流占优、强非对称或病态问题可能极其缓慢甚至发散。因此它们在现代求解器中更多扮演平滑器（smoother）角色，嵌套在多重网格或作为 Krylov 方法的预条件子，而非独立求解器。

非定常 Krylov 方法：在子空间中寻找最优近似

Krylov 子空间方法从初始残差 r₀ = b − Ax₀ 出发，构造子空间 K_m(A, r₀) = span{r₀, Ar₀, A²r₀, …, A^(m−1)r₀}，并在其中寻找满足某种最优性条件的近似解。由于子空间随迭代步数增长，迭代系数每步都在变化，故称“非定常”。不同方法的核心区别在于适用的矩阵性质与最优性准则：

• CG（共轭梯度）：要求 A 对称正定。在 Krylov 子空间中最小化误差的 A-范数，等价于残差正交于子空间。短递推，每步仅需一次矩阵向量乘和少量向量更新，内存为常数级，是对称正定问题的首选。
• GMRES（广义最小残差）：适用于一般非奇异矩阵，在子空间中最小化残差 2-范数 ‖b − Ax_m‖₂。为保证正交基，需存储全部 Arnoldi 向量，内存和每步计算量随迭代数线性增长，通常配合重启（GMRES(m)）使用，但重启会损失最优性、可能停滞。
• BiCGSTAB（双共轭梯度稳定化）：面向非对称问题，通过双正交化与局部残差最小化结合，保持短递推、常数内存，每步约两次矩阵向量乘。收敛行为不如 GMRES 单调，可能出现残差震荡，但对许多非对称问题效率高于重启 GMRES。

定位与适用边界：一张图看清家族关系

下面的信号链式示意图按“矩阵性质 → 候选方法 → 关键约束”组织，帮助快速定位：

如何在实际问题中做初步判断

选型时先看矩阵是否对称正定：若是，CG 几乎总是第一选择，配合不完全 Cholesky 或代数多重网格预条件可处理条件数达 10⁶ 量级的 Poisson 型问题。若矩阵非对称，则需在 GMRES 与 BiCGSTAB 之间权衡：GMRES 残差单调下降、行为可预测，但内存随迭代步线性增长，重启后可能停滞；BiCGSTAB 内存恒定、每步两次矩阵向量乘，适合内存敏感场景，但残差曲线可能剧烈震荡，需设置合理的相对残差容差（如 1e−8）与最大迭代数（如 1000）防止无效循环。

一个常见的教学示例是二维对流扩散方程离散后得到的非对称矩阵。若对流项较弱（Péclet 数小），矩阵接近对称，BiCGSTAB 通常几十步收敛；若对流占优，非对称性增强，GMRES(30) 可能需要数百步，此时应优先考虑加入 ILU(0) 预条件，而非单纯更换 Krylov 方法。这说明方法选择与预条件设计必须联合考虑，单靠换求解器往往收效有限。