网格收敛确认应该用多少个网格水平进行比较？

最少3个水平。粗（coarse）、中（medium）、细（fine）的3个水平来计算GCI（Grid Convergence Index），估计离散化误差带。ASME V&V 20将这种3网格法定为标准步骤。

网格细密化比（refinement ratio）应该设为多少？

r ≥ 1.3是最低要求。ASME V&V 20推荐√2 ≈ 1.41。如果r太小，将导致舍入误差和离散化误差的变化微小，观测次数p的估计精度严重下降。

什么是GCI（Grid Convergence Index）？

GCI是Roache提出的离散化不确定性指标，用GCI = Fs·|ε|/(r^p − 1)计算。以Fs=1.25（3网格法）给出95%置信区间相当的误差带。

什么是Richardson外推？

Richardson外推是用于估计网格大小→0极限解（exact solution）的方法。从3个网格的解f1、f2、f3推断观测次数p，用f_exact ≈ f1 + (f1−f2)/(r^p−1)进行外推值计算。

网格细密化比与GCI收敛确认方法

分类: V&V·网格收敛 | 整合版 2026-04-12

Mesh convergence study diagram showing three mesh levels with GCI error bands and Richardson extrapolation

网格细密化比 r 与 GCI 的离散化误差带概念图

网格细密化比与GCI收敛确认方法的理论基础

概述 — 什么是网格收敛确认

🧑‍🎓

老师，网格收敛确认应该用多少个网格进行比较呢？我听说要"试3个"，为什么2个不行呢？

🎓

最少需要3个网格水平。仅用2个网格比较只能说"结果很接近，应该没问题"。有3个的话就能推断出收敛速度（次数 p），还能外推预测网格大小→0时的极限解。这就是Richardson外推法，其误差带用GCI（Grid Convergence Index）定量化。

🧑‍🎓

我明白了，2个点只能划直线，3个点能看出曲线的弧度。

🎓

正是这样。ASME V&V 20（计算流体力学和固体力学的验证、妥当性确认指南）规定这种3网格法为标准步骤。在论文和项目报告中声称"没有网格依赖性"时，展示GCI数值已成为国际通用语言。

网格收敛确认（Mesh Convergence Study）是定量评估离散化误差、保证解的可信度的V&V过程的核心。其基本要素可归纳为以下3点：

细密化比 r — 粗网格与细网格代表大小的比值。$r \geq 1.3$（推荐 $r = \sqrt{2} \approx 1.41$）
Richardson外推 — 从3个或以上网格的结果推估网格大小→0时的极限解
GCI — 将离散化不确定性定量为95%置信区间

🧑‍🎓

细密化比就是网格变细的"倍数"吗？$r = 2$ 的话把单元大小减半，应该很简单呀。

🎓

结构网格的话 $r = 2$ 确实简单，但单元数会在3D中增加 8倍。计算成本爆炸。所以ASME V&V 20推荐 $r = \sqrt{2} \approx 1.41$。这样单元数只增加约2.8倍。反过来，$r < 1.3$ 的话网格间差异太小，会被舍入误差淹没，观测次数 $p$ 的估计就变得不稳定。

网格细密化比（grid refinement ratio）定义如下：

$$ r = \frac{h_{\text{coarse}}}{h_{\text{fine}}} $$

其中 $h$ 是代表网格大小（representative grid size）。总是令 $r > 1$，分子为粗网格、分母为细网格。

非结构网格的代表大小

🧑‍🎓

结构网格单元大小统一，$h$ 很明确。但四面体网格这样的非结构网格怎样定义 $h$ 呢？

🎓

最标准的方法是，从所有单元体积的平均值反推等效长度。$d$ 是空间维数，2D用面积，3D用体积。

$$ h = \left( \frac{1}{N} \sum_{i=1}^{N} V_i \right)^{1/d} $$

$N$ 是单元数，$V_i$ 是各单元体积（2D情况下为面积），$d$ 是空间维数（2或3）。这样即使是非结构网格，也能得到标量形式的代表大小。

🧑‍🎓

那么非结构网格的细密化比 $r$ 能从单元数的比计算出来吗？

🎓

可以的。等方向细密化时，单元数比与 $r$ 的关系是这样的：

$$ r = \left( \frac{N_{\text{fine}}}{N_{\text{coarse}}} \right)^{1/d} $$

比如3D模型，$N_{\text{fine}} = 800{,}000$、$N_{\text{coarse}} = 200{,}000$，则 $r = (800{,}000 / 200{,}000)^{1/3} = 4^{1/3} \approx 1.587$，细密化比充分了。

Richardson外推

🧑‍🎓

Richardson外推就是"网格无限细的话答案是多少"的预测方法，对吗？具体怎么算？

🎓

离散化的解 $f_h$ 对真解 $f_{\text{exact}}$ 可展开为Taylor级数：

$$ f_h = f_{\text{exact}} + g_p \, h^p + g_{p+1} \, h^{p+1} + \cdots $$

其中 $p$ 是离散化格式的理论次数，$g_p$ 是未知系数。忽略高阶项，从3个网格（$h_1 < h_2 < h_3$, $f_1, f_2, f_3$）就能得到外推值。

$$ f_{\text{exact}} \approx f_1 + \frac{f_1 - f_2}{r_{21}^{\hat{p}} - 1} $$

这里 $r_{21} = h_2 / h_1$，$\hat{p}$ 是从3个解推估的观测次数（后述）。

🧑‍🎓

也就是把最细网格的结果 $f_1$ 进一步校正，使其更接近真解。

🎓

对的。但高阶项被舍弃了，外推值是"估计"而不是真值。所以要用GCI来定量化不确定性。

观测次数 p 的计算

在细密化比相等的情况下（$r_{21} = r_{32} = r$），观测次数 $\hat{p}$ 可直接计算：

$$ \hat{p} = \frac{\ln \left| \dfrac{f_3 - f_2}{f_2 - f_1} \right|}{\ln(r)} $$

$f_1$：细网格解，$f_2$：中网格解，$f_3$：粗网格解。若 $\hat{p}$ 接近理论次数 $p_{\text{formal}}$（如2次单元为2），则解已进入渐近范围（asymptotic range）。

🧑‍🎓

如果 $r_{21} \neq r_{32}$（细密化比不相等）呢？实际工作中网格往往不会完全等比。

🎓

那时解不出来，得用迭代法（不动点迭代或牛顿法）求解以下非线性方程：

$$ \hat{p} = \frac{1}{\ln(r_{21})} \left| \ln \left| \frac{\varepsilon_{32}}{\varepsilon_{21}} \right| + \ln \left( \frac{r_{21}^{\hat{p}} - s}{r_{32}^{\hat{p}} - s} \right) \right| $$

其中 $\varepsilon_{32} = f_3 - f_2$、$\varepsilon_{21} = f_2 - f_1$、$s = \text{sign}(\varepsilon_{32}/\varepsilon_{21})$。通常10次以内迭代就收敛。

GCI（Grid Convergence Index）

🧑‍🎓

GCI公式里有安全系数 $F_s$。这是什么意思？

🎓

$F_s$ 是提高置信度的"保险系数"。Roache的原论文推荐2网格法用 $F_s = 3$（保守），3网格法用 $F_s = 1.25$。3网格法已经验证了观测次数 $\hat{p}$ 的一致性，外推精度更高，所以能用更小的安全系数。

细网格侧的GCI（$\text{GCI}_{21}^{\text{fine}}$）定义为：

$$ \text{GCI}_{21}^{\text{fine}} = \frac{F_s \left| \varepsilon \right|}{r_{21}^{\hat{p}} - 1} $$

其中 $\varepsilon$ 是相对误差：

$$ \varepsilon = \frac{f_2 - f_1}{f_1} $$

方法	安全系数 $F_s$	含义
2网格法	3.0	$\hat{p}$ 未知，保守
3网格法	1.25	$\hat{p}$ 已估计，相当于95%置信区间

🧑‍🎓

想用具体例子算一下。假如3个网格的最大应力分别是 248.3 MPa、251.6 MPa、260.1 MPa？

🎓

设 $f_1 = 248.3$（细）、$f_2 = 251.6$（中）、$f_3 = 260.1$（粗）、$r = \sqrt{2}$：

$\varepsilon_{21} = 251.6 - 248.3 = 3.3$、$\varepsilon_{32} = 260.1 - 251.6 = 8.5$

$\hat{p} = \ln(8.5/3.3) / \ln(\sqrt{2}) = \ln(2.576)/0.3466 \approx 2.73$

$\varepsilon = (251.6 - 248.3)/248.3 = 0.01329$

$\text{GCI}_{21}^{\text{fine}} = 1.25 \times 0.01329 / (1.414^{2.73} - 1) = 0.01661 / 1.634 \approx 0.0102 = 1.02\%$

即细网格结果 248.3 MPa 有 $\pm 1.0\%$ 的离散化不确定性。

渐近范围检查

判断是否进入渐近范围，用以下指标：

$$ \text{AR} = \frac{\text{GCI}_{32}}{r^{\hat{p}} \cdot \text{GCI}_{21}} \approx 1 $$

该值接近1（目标：$0.95 \leq \text{AR} \leq 1.05$）则解已进入渐近范围，GCI和Richardson外推的可信度高。若显著偏离1，说明网格还不够细，或模型含有奇异点。

🧑‍🎓

如果没进入渐近范围怎么办？

🎓

有3个选择：(1) 添加更细的网格重新计算，(2) 把评估点远离奇异点（如改用截面平均应力而非应力集中点），(3) 用更高阶的单元以更快到达渐近范围。实际上(2)成本效益最好。

网格细密化比与GCI收敛确认方法的数值计算方法

3网格法的步骤

基于ASME V&V 20的3网格法标准步骤如下：

生成3水平网格：细密化比 $r \geq 1.3$（推荐 $\sqrt{2}$）的粗、中、细3级。网格拓扑保持一致
统一条件分析：边界条件、材料、求解器设置完全相同，分别运行3个算例
提取评估量：记录关注的物理量（最大应力、平均温度、压力损失等）的值 $f_1, f_2, f_3$
计算代表网格大小 $h_1, h_2, h_3$：结构网格用单元边长，非结构网格用前述体积平均法
计算观测次数 $\hat{p}$
计算Richardson外推值 $f_{\text{ext}}$
计算GCI**：$\text{GCI}_{21}^{\text{fine}}$ 和 $\text{GCI}_{32}^{\text{coarse}}$
渐近范围检查：计算AR，接近1为佳

🧑‍🎓

步骤明白了。但结构网格和非结构网格都要"保持相同拓扑"，非结构网格这样做困难吗？

🎓

非结构网格确实不能完全保持拓扑。实际做法是"全局网格大小统一乘以1/r"，局部加密区域在3个水平也采用相同区域和相同比例加密。关键是让网格生成"配方"统一。

非均匀细密化比的扩展

当 $r_{21} \neq r_{32}$ 时，观测次数 $\hat{p}$ 需用迭代法求解前述非线性方程。实际中固定点迭代即可：

$$ \hat{p}_{k+1} = \frac{1}{\ln(r_{21})} \left| \ln \left| \frac{\varepsilon_{32}}{\varepsilon_{21}} \right| + \ln \left( \frac{r_{21}^{\hat{p}_k} - s}{r_{32}^{\hat{p}_k} - s} \right) \right| $$

初值 $\hat{p}_0$ 用等比公式计算，迭代至 $|\hat{p}_{k+1} - \hat{p}_k| < 10^{-6}$ 停止。

多个评估量的处理

🧑‍🎓

最大应力的GCI是1%，但挠度的GCI是5%。哪个结果该报告？

🎓

全部报告才对。GCI是逐个评估量独立计算的，不能说"整个模型GCI 1%"。局部量（应力集中点的应力）和全局量（总反力、最大挠度）的收敛性完全不同，局部量收敛慢得多。影响设计判断的量至少评估3个，逐个报告GCI。

网格细密化比与GCI收敛确认方法的工程应用

网格收敛确认的工作流程

步骤	作业内容	判定标准
1. 生成网格	粗、中、细3水平（$r \geq 1.3$）	确认 $r$，网格质量指标（长宽比 < 5等）
2. 运行分析	同一边界条件和求解器设置运行3例	全部算例正常收敛
3. 计算 $\hat{p}$	观测次数的推算	$0 < \hat{p} \leq p_{\text{formal}} + 1$
4. 计算GCI	$\text{GCI}_{21}^{\text{fine}}$ 的推算	与设计余量相符
5. AR检查	渐近范围的确认	$0.95 \leq \text{AR} \leq 1.05$
6. 报告	GCI值和网格信息记录	第三者可复现的细节度

常见失败与对策

🧑‍🎓

前辈说"做网格收敛，结果GCI是50%"。怎么会这样？

🎓

GCI过大的原因前3名如下：

在奇异点评估 — 应力集中点（缺口根部、荷载点）的理论解趋于无穷，再细网格也不收敛。对策：在离奇异点的位置评估，或用截面平均应力
细密化比太小 — $r < 1.1$ 时舍入误差和离散化误差分不开。对策：保证 $r \geq 1.3$
分析未收敛 — 非线性分析求解器收敛判定太松，各网格的解本身不准。对策：严格设置余数

🧑‍🎓

奇异点那个太关键了。角部应力是无穷大，网格越细数值越大…

🎓

对。初学CAE的常见错误就是"用最大应力做收敛确认"，特别是模型有锐角的时候，最大应力会随网格细密化单调上升，永远不收敛。评估量的选择决定了网格收敛确认的成败。

报告书的记载事项

网格收敛确认结果报告时，必须包含以下信息：

各网格的单元数、代表大小 $h$
细密化比 $r_{21}$, $r_{32}$
评估量（物理量名、位置、方向）
各网格评估量的值 $f_1, f_2, f_3$
观测次数 $\hat{p}$
Richardson外推值 $f_{\text{ext}}$
$\text{GCI}_{21}^{\text{fine}}$（百分比表示）
渐近范围检查 AR
网格质量指标（最小长宽比、最大歪斜度等）

软件对应

主要CAE工具的网格收敛功能

工具	自动网格收敛	GCI计算	备注
Ansys Mechanical	Convergence Tool（自适应细密化）	手动或ACT Extension	Workbench中可参数化扫描网格大小
Abaqus	Python脚本 + 自适应重网格	手动（Python脚本）	通过abaqus.rpy参数化网格大小并自动化
COMSOL Multiphysics	参数化扫描 + 自适应网格	手动	以网格大小为参数的扫描很便利
OpenFOAM	snappyHexMesh的细密度指定	手动	pyFoam等工具的公开脚本可自动化GCI
Ansys Fluent	自适应网格细密化(AMR)	手动或日志文件	可与求解适应型细密化结合
Simcenter STAR-CCM+	网格流程 + 设计管理器	手动	设计管理器可自动扫描网格大小

🧑‍🎓

所有软件都是"GCI手动计算"？没有自动的工具吗？

🎓

GCI计算本身只是代入公式，Excel或脚本可轻松搞定。真正重要的是"统一网格生成配方、高效运行3个算例"的机制。这方面Ansys Workbench的Design Points或COMSOL的参数化扫描有优势。而GCI计算可以自己写简单脚本。

网格细密化比与GCI收敛确认方法的先端研究

自适应网格与GCI的融合

求解适应型网格细密化（Adaptive Mesh Refinement: AMR）是根据误差估计量局部细密化网格的方法，与GCI本质不同。但近来有融合两者的研究。

🧑‍🎓

AMR时局部网格不同，"代表大小 $h$ 定义"难吗？

🎓

确实。AMR时GCI应用有约束。一个思路是"固定AMR最大细密化级，仅改变基础网格大小乘以 $r$"，这样能直接套用GCI框架。另一个思路是用AMR加事后误差估计量（ZZ估计或SPR法），不用GCI直接估计离散化误差。

最小二乘GCI（多网格法）

当网格数≥4时，用最小二乘法估计 $\hat{p}$ 和 $g_p$，可大幅提高外推稳定性。特别是观测次数不稳定（未进入渐近范围）时有效。Eqa等(2014)已体系化：

$$ \min_{\hat{p}, \, g_p} \sum_{k=1}^{M} \left( f_k - f_{\text{ext}} - g_p \, h_k^{\hat{p}} \right)^2 $$

用 $M \geq 4$ 个网格水平求解上式，能获得比3网格法更稳健的 $\hat{p}$ 和 $f_{\text{ext}}$。成本高但可信度也高，用于认证分析（certification analysis）。

网格细密化比与GCI收敛确认方法的故障排除

观测次数发散

🧑‍🎓

$\hat{p}$ 跳到8甚至15，远超理论次数，咋回事？

🎓

$\hat{p}$ 异常大说明 $\varepsilon_{21}$ 与 $\varepsilon_{32}$ 比值接近1。即中网格和细网格差别很小，但粗网格差别大。原因多是进入舍入误差区或中细网格实际分辨率接近。对策：增大 $r$ 或提高精度（双精度→四精度）。

p 为负值

$\hat{p} < 0$ 表示 $\varepsilon_{32}$ 和 $\varepsilon_{21}$ 符号相反（振荡收敛），此时单调收敛假设破缺，GCI不适用。

🎓

振荡收敛常见于CFD中格式次数不足或应力评估点靠近网格边界。改变评估点位置或升级离散化格式通常能解决。

非单调收敛

🧑‍🎓

细→中→粗的结果不单调变化，Richardson外推就用不了吧？

🎓

确实用不了。非单调收敛时要检查：

网格质量：细网格是否产生了高长宽比或薄片单元
边界条件一致性：网格改变时荷载面或约束面定义有没有偏移
物理模型：接触、湍流模型有没有网格大小依赖
舍入误差：解的差 $10^{-12}$ 以下说明双精度已用尽

全部检查完还非单调，改用4网格以上的最小二乘法更稳妥。