4. 代码验证的方法与实践

代码验证（Code Verification）回答的是一个非常具体的问题：我写的求解器，是否正确地求解了我所设定的方程？它不关心物理模型是否选对，只关心离散方程、迭代过程、边界条件实现和编程实现是否与数学方程一致。因此，代码验证的核心手段是构造“已知答案”的问题，再观察数值解是否按预期收敛到该答案。本节围绕五类常用方法展开：网格收敛指数（GCI）、制造解方法（MMS）、解析解对比、代码间基准测试，以及单位与量纲一致性检查。

一、网格收敛指数：用三套网格估计离散误差

GCI 的理论基础是 Richardson 外推。设网格特征尺寸为 h，数值解为 f，则离散误差可写成 f_exact ≈ f_h + C h^p，其中 p 为观测精度阶。实际操作时，取三套按固定比例 r 加密的网格（常用 r = 2 或 1.5），得到 f₁（细）、f₂（中）、f₃（粗），先计算收敛比：

R = (f2 - f1) / (f3 - f2)

若 0 < R < 1，说明单调收敛，可求观测阶：
p = ln|(f3 - f2)/(f2 - f1)| / ln(r)

Richardson 外推值：
f_exact ≈ f1 + (f1 - f2) / (r^p - 1)

相对离散误差：
e_a = |(f1 - f2) / f1|

GCI_fine = Fs * e_a / (r^p - 1)，安全因子 Fs 通常取 1.25

举例说明：某二维通道流动中，监测出口截面平均速度，三套网格结果分别为 f₁ = 1.982 m/s、f₂ = 1.965 m/s、f₃ = 1.921 m/s，r = 2。则 R = (1.965−1.982)/(1.921−1.965) = 0.386，落在 (0,1) 内，单调收敛。观测阶 p = ln|(1.921−1.965)/(1.965−1.982)| / ln2 = ln(2.588)/0.693 ≈ 1.37。外推值约为 1.982 + (1.982−1.965)/(2^1.37−1) ≈ 1.992 m/s，GCI_fine ≈ 1.25×|(1.982−1.965)/1.982|/(2^1.37−1) ≈ 1.6%。这说明细网格解约有 1.6% 的离散误差带，若工程允许误差为 5%，该网格可接受。需要强调：以上数值为教学示例，用于演示计算流程，并非某次真实仿真的实测结果。

GCI 的实践要点有三：一是必须使用同一套离散格式和同一求解设置，只改变网格；二是监测量的选择要有代表性，优先选对目标量敏感的位置（如驻点压力、分离点位置、壁面热流）；三是若 R 为负或大于 1，说明网格未进入渐近收敛区，需继续加密或检查是否存在振荡、边界层分辨率不足等问题。

二、制造解方法：给代码“出题”

MMS 的思路是反向构造：先人为指定一个解析形式的解 u*(x,y,t)，代入控制方程，求出残余项作为源项 S，再把 S 和对应的边界条件一起加入代码，最后比较数值解与 u*。这样即使原方程没有解析解，也能获得精确参照。

以二维不可压 Navier–Stokes 方程为例，设制造解为：

u*(x,y) = sin(πx) cos(πy)
v*(x,y) = -cos(πx) sin(πy)
p*(x,y) = cos(πx) cos(πy)

将其代入动量方程，得到的残余项即为应加入代码的源项。随后在 N = 16、32、64、128 四套网格上运行，计算 L₂ 误差范数，若误差随 h 的下降斜率接近格式理论阶（如二阶格式斜率约 −2），则说明离散实现正确。MMS 的最大价值在于：它能同时检验对流项、扩散项、源项、时间推进和边界条件实现，且不依赖任何实验数据。

三、解析解对比与代码间基准测试

解析解对比适用于有经典精确解的场景，如 Poiseuille 流、Couette 流、Stokes 第一问题、一维热传导等。做法是设置与解析解一致的边界条件和物性，比较速度剖面、温度分布或壁面梯度的误差。它比 MMS 更贴近物理，但覆盖的方程项有限。

代码间基准测试（Code-to-Code Benchmark）则用于没有解析解、又需要交叉确认的复杂问题。选择同一几何、同一工况、同一物性，由两套独立代码分别计算，比较关键监测量。判断标准不是“完全相等”，而是差异是否落在各自 GCI 估计的误差带内。若差异显著超出双方 GCI，则至少有一套代码存在问题，需回到网格、边界条件或离散格式逐项排查。

四、单位与量纲一致性检查

量纲错误是代码验证中最容易被忽视、却最常导致“结果看起来合理但整体偏移”的原因。检查应贯穿输入、方程和输出三个环节：输入物性是否统一到 SI（如动力粘度 Pa·s 与运动粘度 m²/s 混用）；源项量纲是否与方程各项一致；输出量是否与监测量定义匹配。一个实用做法是建立量纲检查表，对每个变量标注符号、单位、量纲和换算系数，并在代码中显式写出单位注释。

五、实践中的组合策略与常见陷阱

五类方法并非互斥，成熟的做法是组合使用：先用 MMS 验证代码实现，再用解析解验证物理项，接着用 GCI 量化离散误差，最后用代码间基准测试覆盖复杂工况，量纲检查则贯穿始终。常见陷阱包括：只做一套网格就宣称“网格无关”；MMS 源项符号写反导致误差不降反升；代码间比较时两套代码用了不同的边界条件却未察觉；以及单位混用造成结果整体偏移一个常数因子。识别这些陷阱的关键，是让每一步都有明确的“已知参照”和“可量化判据”，而不是凭经验判断结果“看起来对”。