2. 偏微分方程离散化回顾
迭代求解器的选型从来不是从算法本身出发的，而是从“离散之后得到的线性系统长什么样”出发的。同一个偏微分方程，用有限差分、有限元或有限体积离散，得到的矩阵在稀疏结构、对称性、条件数上可能差异巨大，而这些差异直接决定了应该用 Jacobi、Gauss-Seidel、共轭梯度、GMRES、多重网格还是区域分解。本节把三种主流离散方式产生的线性系统特征梳理清楚，为后续选型建立判断依据。
三种离散方式产生的矩阵从哪里来
考虑一个典型的椭圆型模型问题：在区域 Ω 上求解 −∇·(k∇u)=f，配以 Dirichlet 或 Neumann 边界条件。三种离散方式的出发点不同，但最终都归结为求解形如 Au = b 的线性系统，其中 A 是 N×N 矩阵，N 为未知量个数。
有限差分（FDM）直接在网格节点上近似导数。以二维 Poisson 方程在均匀网格（步长 h）上的五点格式为例，每个内部节点得到一个方程：
(4u[i,j] − u[i-1,j] − u[i+1,j] − u[i,j-1] − u[i,j+1]) / h² = f[i,j]每个方程只涉及自身和上下左右四个邻居，因此矩阵每行最多 5 个非零元。按自然行序编号时，矩阵呈现块三对角结构，带宽约为网格一行的节点数。FDM 矩阵的对称性取决于算子本身：常系数 Poisson 的五点格式矩阵是对称正定的；但若加入一阶对流项并用中心差分，对称性就被破坏。
有限元（FEM）从弱形式出发，在单元上做积分并用基函数展开。以线性三角形单元为例，刚度矩阵的元素为 A[i,j] = ∫∇φ_i·∇φ_j dΩ。由于基函数具有局部支撑，只有共享同一单元的节点之间才有非零耦合，因此矩阵同样是稀疏的，每行非零元个数与节点的邻接节点数相关。对于对称的微分算子（如扩散项），FEM 刚度矩阵天然对称；若采用质量集中或某些稳定化格式，对称性可能被改变。FEM 的一个重要特点是矩阵元素来自数值积分，系数可能不像 FDM 那样整齐，但结构上的稀疏模式与网格拓扑严格对应。
有限体积（FVM）在每个控制体积上对守恒型方程积分，通量在控制体界面上计算。以单元中心格式为例，每个控制体的方程涉及通过各界面的通量，因此矩阵每行的非零元对应本单元及其所有面邻居。对于结构网格上的二维问题，这与五点格式形式相似；但在非结构网格上，邻居数不固定，稀疏模式更不规则。FVM 的系数通常由界面通量系数组装，若采用迎风格式处理对流项，矩阵会失去对称性，且往往成为非对称的 M 矩阵。
稀疏性：不只是“零很多”
三种离散方式得到的矩阵都是稀疏的，但稀疏的“形状”不同，这直接影响存储格式和求解器的选择。FDM 在结构网格上产生规则带状矩阵，可以用带状存储或简单的 CSR 格式高效处理；FEM 和 FVM 在非结构网格上产生不规则稀疏模式，需要通用的 CSR/CSC 存储，且矩阵-向量乘法的间接寻址开销更大。
更关键的是，稀疏性决定了迭代求解器每步迭代的计算量。以 CSR 存储为例，一次 SpMV（稀疏矩阵向量乘）的浮点运算量约为 2·nnz（nnz 为非零元总数），而间接寻址次数与 nnz 同阶。对于三维问题，FDM 七点格式每行 7 个非零元，FEM 四面体单元每行非零元数约为邻接节点数加一，通常也在 10～15 之间。这意味着每步迭代的成本大致与未知量数 N 成正比，这是迭代法相对于直接法（O(N²) 或 O(N^1.5)）的核心优势。
实际实现时，一个常见的判断是：若矩阵每行非零元数不超过约 20，且 N 超过一万，迭代法通常比直接法更合适；若 N 较小或矩阵接近稠密，直接法可能更稳。这个阈值不是绝对的，但可以作为初步筛选的依据。
对称性：决定能不能用 CG
对称性是最“硬”的一个特征，因为它直接排除了共轭梯度法（CG）。CG 要求矩阵对称正定，若矩阵非对称，CG 的理论基础不成立，迭代可能不收敛或收敛到错误结果。
对于纯扩散问题，三种离散方式在适当格式下都能得到对称矩阵：FDM 的五点/七点格式对称，FEM 的刚度矩阵对称，FVM 在正交网格上纯扩散通量离散也对称。但一旦引入对流项，情况就分化了。FDM 用中心差分处理对流项时，矩阵仍可能保持对称（取决于离散细节），但可能产生非物理振荡；FVM 为保证守恒性和有界性，通常采用迎风格式，此时矩阵严格非对称。FEM 若采用 SUPG 等稳定化格式，也会引入非对称项。
因此，实际选型时的判断顺序是：先看矩阵是否对称。若对称正定，优先考虑 CG 配合预条件（如不完全 Cholesky、代数多重网格）；若非对称，则转向 GMRES、BiCGSTAB 等 Krylov 方法，并配合 ILU 类预条件。这个判断不需要等到矩阵组装完成，从离散格式和方程类型就能预先确定。
条件数：迭代收敛速度的“隐形天花板”
条件数 κ(A) 衡量的是矩阵对扰动的敏感程度，也直接控制 Krylov 方法的收敛速度。对于二阶椭圆问题，无论哪种离散方式，条件数都随网格加密而恶化，典型量级为 O(h⁻²)。这意味着 h 减半时，条件数约增大 4 倍，未预条件的 CG 迭代步数大致按 √κ 增长，即步数约翻倍。
这个规律对三种离散方式都成立，但常数因子不同。FEM 在畸形单元上条件数会额外恶化；FVM 在非结构网格上，界面系数差异大时条件数也会增大；FDM 在均匀网格上条件数最“干净”，但一旦网格拉伸或系数变化，同样会恶化。实际教学中可以做一个简单估算：二维 Poisson 问题，h=1/100，N≈10⁴，κ≈O(10⁴)，未预条件 CG 可能需要数百步；若用代数多重网格预条件，步数可降到十步量级，且几乎与 h 无关。
把三个特征串起来：一个选型判断表
把稀疏性、对称性、条件数放在一起，可以形成一个初步的选型判断路径。以下判断表以“离散方式 + 方程类型”为输入，给出矩阵特征和对应的求解器方向，可作为实际工程中的快速参考。
| 离散方式 | 方程类型 | 对称性 | 条件数量级 | 求解器方向 |
|---|---|---|---|---|
| FDM（结构网格） | 纯扩散 | 对称正定 | O(h⁻²) | CG + 几何/代数多重网格 |
| FDM（结构网格） | 对流-扩散（中心差分） | 可能对称，但振荡 | O(h⁻²)，Péclet 数大时恶化 | GMRES + ILU，或改用迎风 |
| FEM（非结构） | 纯扩散 | 对称正定 | O(h⁻²)，受单元质量影响 | CG + 代数多重网格 |
| FEM（非结构） | 对流-扩散（SUPG） | 非对称 | O(h⁻²) 以上 | GMRES / BiCGSTAB + ILU |
| FVM（非结构） | 纯扩散 | 正交网格对称，斜网格可能非对称 | O(h⁻²)，界面系数差异大时更大 | CG（对称时）/ GMRES + AMG |
| FVM（非结构） | 对流-扩散（迎风） | 非对称 | O(h⁻²) 以上，CFL 相关 | GMRES / BiCGSTAB + ILU |
需要强调的是，这张表给出的是“方向”而非“定论”。实际选型还要考虑实现成本、并行性、内存占用和问题规模。例如，代数多重网格在对称正定问题上表现优异，但实现复杂；ILU 预条件通用性强，但并行扩展性受限。下一节将进入具体迭代方法的讨论，届时会结合本节建立的矩阵特征，逐一分析各方法的适用条件和实际表现。