一条数值结论有多层来源
计算物理不是把公式输入程序后得到一个“精确答案”。实验对象先被理想化为物理模型，再写成连续或随机方程，随后经过空间与时间离散、代数或非线性迭代、有限精度运算及统计汇总。每层回答不同问题：
建模误差来自忽略的物理、闭合关系和适用范围。
参数与初边值不确定度来自测量、标定或估计。
截断或离散误差来自用有限网格、时间步和有限基逼近连续问题。
迭代误差来自在线性、非线性或优化求解器未完全收敛时停止。
舍入误差来自有限浮点表示、消去和归约次序。
统计误差来自有限随机样本或有限相关轨迹。
实现错误不是应当写进预算后接受的第七类误差。数组越界、边界符号写反或单位转换错误应通过测试和 verification 找出并修复。误差预算描述在已知实现下仍不可避免或有意保留的近似。
若目标量 Q Q Q W m − 2 \mathrm{W\,m^{-2}} W m − 2 Q h − Q Q_h-Q Q h − Q W m − 2 \mathrm{W\,m^{-2}} W m − 2 ∣ Q h − Q ∣ / Q r e f |Q_h-Q|/Q_{\mathrm{ref}} ∣ Q h − Q ∣/ Q ref
Taylor 截断误差与一致性
以光滑函数的一阶导数为例，中心差分
D h f ( x ) = f ( x + h ) − f ( x − h ) 2 h D_hf(x)=\frac{f(x+h)-f(x-h)}{2h} D h f ( x ) = 2 h f ( x + h ) − f ( x − h ) 由 Taylor 展开得到
D h f = f ′ ( x ) + h 2 6 f ( 3 ) ( x ) + O ( h 4 ) . D_hf=f'(x)+\frac{h^2}{6}f^{(3)}(x)+O(h^4). D h f = f ′ ( x ) + 6 h 2 f ( 3 ) ( x ) + O ( h 4 ) . 若 x , h x,h x , h f f f D h f D_hf D h f K m − 1 \mathrm{K\,m^{-1}} K m − 1 h , Δ t → 0 h,\Delta t\to0 h , Δ t → 0
例 1：用两次加密观测中心差分阶数
取无量纲变量 f ( x ) = sin x f(x)=\sin x f ( x ) = sin x x = 0.5 x=0.5 x = 0.5 f ′ ( 0.5 ) = cos 0.5 ≈ 0.877583 f'(0.5)=\cos0.5\approx0.877583 f ′ ( 0.5 ) = cos 0.5 ≈ 0.877583 h = 0.10 h=0.10 h = 0.10
D 0.10 f = sin 0.60 − sin 0.40 0.20 ≈ 0.876121 , D_{0.10}f
=\frac{\sin0.60-\sin0.40}{0.20}
\approx0.876121, D 0.10 f = 0.20 sin 0.60 − sin 0.40 ≈ 0.876121 , 绝对误差约 1.46 × 10 − 3 1.46\times10^{-3} 1.46 × 1 0 − 3 h = 0.05 h=0.05 h = 0.05
D 0.05 f = sin 0.55 − sin 0.45 0.10 ≈ 0.877217 , D_{0.05}f
=\frac{\sin0.55-\sin0.45}{0.10}
\approx0.877217, D 0.05 f = 0.10 sin 0.55 − sin 0.45 ≈ 0.877217 , 误差约 3.66 × 10 − 4 3.66\times10^{-4} 3.66 × 1 0 − 4 3.99 3.99 3.99 p o b s = log 2 ( 3.99 ) ≈ 2.00 p_{\mathrm{obs}}=\log_2(3.99)\approx2.00 p obs = log 2 ( 3.99 ) ≈ 2.00 h h h
稳定性、收敛与适用前提
稳定性研究离散推进对初值误差、边界误差和每步舍入扰动的放大。对线性、适定的初值问题，Lax 等价框架说明一致且稳定的线性差分格式收敛，反之亦然。这个结论有明确前提，不能直接替代非线性方程、移动间断、随机算法或非适定逆问题的专门分析。
考虑一维热方程
∂ u ∂ t = α ∂ 2 u ∂ x 2 , [ α ] = m 2 s − 1 . \frac{\partial u}{\partial t}
=\alpha\frac{\partial^2u}{\partial x^2},
\qquad [\alpha]=\mathrm{m^2\,s^{-1}}. ∂ t ∂ u = α ∂ x 2 ∂ 2 u , [ α ] = m 2 s − 1 . 显式中心格式
u i n + 1 = u i n + r ( u i + 1 n − 2 u i n + u i − 1 n ) , r = α Δ t Δ x 2 u_i^{n+1}
=u_i^n+r(u_{i+1}^n-2u_i^n+u_{i-1}^n),
\qquad
r=\frac{\alpha\Delta t}{\Delta x^2} u i n + 1 = u i n + r ( u i + 1 n − 2 u i n + u i − 1 n ) , r = Δ x 2 α Δ t 在一维周期或相容边界下要求 0 ≤ r ≤ 1 / 2 0\le r\le1/2 0 ≤ r ≤ 1/2 r r r
例 2：热扩散显式步长上限
取 α = 1.00 × 10 − 5 m 2 s − 1 \alpha=1.00\times10^{-5}\,\mathrm{m^2\,s^{-1}} α = 1.00 × 1 0 − 5 m 2 s − 1 Δ x = 1.00 m m = 10 − 3 m \Delta x=1.00\,\mathrm{mm}=10^{-3}\,\mathrm m Δ x = 1.00 mm = 1 0 − 3 m
Δ t ≤ Δ x 2 2 α = 10 − 6 2.00 × 10 − 5 = 0.0500 s . \Delta t\le
\frac{\Delta x^2}{2\alpha}
=\frac{10^{-6}}{2.00\times10^{-5}}
=0.0500\,\mathrm s. Δ t ≤ 2 α Δ x 2 = 2.00 × 1 0 − 5 1 0 − 6 = 0.0500 s . 选 Δ t = 0.0400 s \Delta t=0.0400\,\mathrm s Δ t = 0.0400 s r = 0.400 r=0.400 r = 0.400 Δ x / 2 \Delta x/2 Δ x /2 r r r
三层加密与观测收敛阶
设某标量输出在渐近区满足
Q h = Q + C h p + O ( h p + 1 ) . Q_h=Q+C h^p+O(h^{p+1}). Q h = Q + C h p + O ( h p + 1 ) . 对比例为 2 的三层网格 h , h / 2 , h / 4 h,h/2,h/4 h , h /2 , h /4
p o b s = ln ∣ ( Q h − Q h / 2 ) / ( Q h / 2 − Q h / 4 ) ∣ ln 2 . p_{\mathrm{obs}}
=\frac{\ln\left|
(Q_h-Q_{h/2})/(Q_{h/2}-Q_{h/4})
\right|}{\ln2}. p obs = ln 2 ln ( Q h − Q h /2 ) / ( Q h /2 − Q h /4 ) . 至少三层才能同时观察阶数与细化趋势。若差值变号、阶数剧烈变化或局部自适应网格不能用单一 h h h
场误差还必须声明范数。对控制体误差 e i e_i e i V i V_i V i
∥ e ∥ 2 , V = ( ∑ i e i 2 V i ∑ i V i ) 1 / 2 , ∥ e ∥ ∞ = max i ∣ e i ∣ . \|e\|_{2,V}
=\left(
\frac{\sum_i e_i^2V_i}{\sum_iV_i}
\right)^{1/2},
\qquad
\|e\|_\infty=\max_i|e_i|. ∥ e ∥ 2 , V = ( ∑ i V i ∑ i e i 2 V i ) 1/2 , ∥ e ∥ ∞ = i max ∣ e i ∣. 二者单位都与被比较场相同。直接对网格向量用未加权 Euclidean 范数，会随节点数改变尺度，跨网格不一定可比。L 2 L^2 L 2 L ∞ L^\infty L ∞
一个积分输出偶然达到高阶，也不能证明整个场同阶收敛。反之，局部尖角附近点值阶数降低，某个守恒通量仍可能准确。验证计划应同时选与决策相关的目标量、至少一个全场范数和守恒残差，并预先写明接受阈值。
例 3：热流的观测阶与 Richardson 外推
某壁面平均热流在三层网格上为
Q h = 1.200 , Q h / 2 = 1.050 , Q h / 4 = 1.0125 k W m − 2 . Q_h=1.200,\quad
Q_{h/2}=1.050,\quad
Q_{h/4}=1.0125
\quad \mathrm{kW\,m^{-2}}. Q h = 1.200 , Q h /2 = 1.050 , Q h /4 = 1.0125 kW m − 2 . 相邻差值为 0.150 0.150 0.150 0.0375 k W m − 2 0.0375\,\mathrm{kW\,m^{-2}} 0.0375 kW m − 2 p o b s = 2 p_{\mathrm{obs}}=2 p obs = 2
Q e x t = Q h / 4 + Q h / 4 − Q h / 2 2 p − 1 = 1.0000 k W m − 2 . Q_{\mathrm{ext}}
=Q_{h/4}
+\frac{Q_{h/4}-Q_{h/2}}{2^p-1}
=1.0000\,\mathrm{kW\,m^{-2}}. Q ext = Q h /4 + 2 p − 1 Q h /4 − Q h /2 = 1.0000 kW m − 2 . 最细网格相对外推值高 1.25 % 1.25\% 1.25%
迭代误差、条件数与舍入下限
离散后常得到 A x = b A\boldsymbol x=\boldsymbol b A x = b r = b − A x ^ \boldsymbol r=\boldsymbol b-A\hat{\boldsymbol x} r = b − A x ^ κ ( A ) \kappa(A) κ ( A ) 10 − 3 10^{-3} 1 0 − 3 10 − 12 10^{-12} 1 0 − 12 10 − 2 10^{-2} 1 0 − 2
减小 h h h
守恒残差是带单位的独立检查
对有限体积或粒子系统，全局质量、动量、能量与电荷平衡提供不依赖局部点值的检查。固定控制体质量缺陷可写为
δ m = M ( t 2 ) − M ( t 1 ) − ∫ t 1 t 2 ( m ˙ i n − m ˙ o u t ) d t , \delta_m
=M(t_2)-M(t_1)
-\int_{t_1}^{t_2}
(\dot m_{\mathrm{in}}-\dot m_{\mathrm{out}})\,\mathrm dt, δ m = M ( t 2 ) − M ( t 1 ) − ∫ t 1 t 2 ( m ˙ in − m ˙ out ) d t , 单位为千克。若按时间报告平衡残差
δ m / ( t 2 − t 1 ) \delta_m/(t_2-t_1) δ m / ( t 2 − t 1 ) k g s − 1 \mathrm{kg\,s^{-1}} kg s − 1
例 4：储罐计算的质量平衡缺陷
初始质量 10.0 k g 10.0\,\mathrm{kg} 10.0 kg 2.00 k g s − 1 2.00\,\mathrm{kg\,s^{-1}} 2.00 kg s − 1 0.500 k g s − 1 0.500\,\mathrm{kg\,s^{-1}} 0.500 kg s − 1
M ⋆ = 10.0 + ( 2.00 − 0.500 ) × 4.00 = 16.0 k g . M_\star=10.0+(2.00-0.500)\times4.00
=16.0\,\mathrm{kg}. M ⋆ = 10.0 + ( 2.00 − 0.500 ) × 4.00 = 16.0 kg . 某计算给 15.94 k g 15.94\,\mathrm{kg} 15.94 kg δ m = − 0.060 k g \delta_m=-0.060\,\mathrm{kg} δ m = − 0.060 kg − 0.015 k g s − 1 -0.015\,\mathrm{kg\,s^{-1}} − 0.015 kg s − 1 − 0.375 % -0.375\% − 0.375%
守恒到机器精度也不等于解正确：一个符号写错但内部成对抵消的程序仍可能守恒。守恒、收敛阶和独立基准必须联合使用。
制造解、解析基准与 verification
制造解方法先任取足够光滑且满足所需边界类型的函数
u m ( x , t ) u_m(\boldsymbol x,t) u m ( x , t )
f m = L ( u m ) . f_m=\mathcal L(u_m). f m = L ( u m ) . 数值代码用 f m f_m f m u m u_m u m
解析解、独立高精度求解器、对称极限和已知守恒解也是基准。Code verification 关注实现与离散阶；solution verification 估计某次生产计算的网格、时间步和迭代不确定度。只在一个网格与一个参考数值吻合，无法区分偶然误差抵消。
Validation 的证据边界
Validation 比较模型预测与实验或观察数据，问的是“模型在指定用途和条件下是否足够”。它必须记录测量不确定度、参数校准数据与验证数据是否独立、比较量和单位、空间时间对齐以及接受标准。用同一数据既拟合参数又宣称独立验证，会高估证据。
verification 通过不能保证 validation 通过：方程可以被准确求解但遗漏关键物理。反过来，单个实验点吻合也不能证明程序正确，离散误差和模型偏差可能抵消。Validation 支持的是限定工况、输出和误差容限内的使用主张，不是对模型“永远正确”的认证。
误差预算与报告
每个误差项至少记录来源、估计方法、单位、符号或区间、是否随机以及相关性。只有当各项可视为独立零均值标准不确定度且单位一致时，才可用平方和开根号：
u c = u 1 2 + ⋯ + u n 2 . u_c=\sqrt{u_1^2+\cdots+u_n^2}. u c = u 1 2 + ⋯ + u n 2 . 已知偏差应校正或单列，有界模型差应作为情景区间，不能都伪装成 Gaussian 标准差。
例 5：温度预测的分层误差预算
某温度输出为 T = 350.0 K T=350.0\,\mathrm K T = 350.0 K 0.40 K 0.40\,\mathrm K 0.40 K 0.05 K 0.05\,\mathrm K 0.05 K 0.30 K 0.30\,\mathrm K 0.30 K 0.01 K 0.01\,\mathrm K 0.01 K
u n u m ≈ 0.40 2 + 0.05 2 + 0.30 2 + 0.01 2 = 0.502 K . u_{\mathrm{num}}
\approx\sqrt{0.40^2+0.05^2+0.30^2+0.01^2}
=0.502\,\mathrm K. u num ≈ 0.4 0 2 + 0.0 5 2 + 0.3 0 2 + 0.0 1 2 = 0.502 K . 另有模型简化引起的有向偏差区间约
[ − 1.0 , 0.2 ] K [-1.0,0.2]\,\mathrm K [ − 1.0 , 0.2 ] K 0.50 K 0.50\,\mathrm K 0.50 K ± 1.1 K \pm1.1\,\mathrm K ± 1.1 K
常见误区
常见误区
“残差降到机器精度就表示物理解精确。”残差只说明离散代数方程被解得很紧，模型和离散误差仍可更大。
常见误区
“两个网格给相近结果已经证明收敛阶。”至少需要三层观察阶数，并确认处于渐近区且其他误差更小。
常见误区
“实验吻合就是 verification。”实验比较属于 validation；代码与离散实现需要制造解、解析基准和加密研究独立验证。
练习：从截断误差到验证主张
练习
所属知识 误差分类
难度 2/5 分别给“忽略辐射、有限网格、求解器提前停止、相近浮点数相减、有限样本均值”分类。
查看提示 判断误差发生在物理方程、网格、代数求解、浮点或采样哪一层。
查看解答 忽略辐射是建模误差；有限
Δ x \Delta x Δ x 是离散误差；线性迭代提前停止是迭代误差；相消损失是舍入误差；有限随机样本均值波动是统计误差。
练习
所属知识 观测阶
难度 3/5 三层输出 1.120 , 1.040 , 1.020 P a 1.120,1.040,1.020\,\mathrm{Pa} 1.120 , 1.040 , 1.020 Pa
查看提示 用
p = ln ( ∣ Q h − Q h / 2 ∣ / ∣ Q h / 2 − Q h / 4 ∣ ) / ln 2 p=\ln(|Q_h-Q_h/2|/|Q_h/2-Q_h/4|)/\ln 2 p = ln ( ∣ Q h − Q h /2∣/∣ Q h /2 − Q h /4∣ ) / ln 2 。
查看解答 差值为 0.080 和 0.020，比例 4，所以 p=2。外推值为
1.020 + ( 1.020 − 1.040 ) / 3 = 1.01333 1.020+(1.020-1.040)/3=1.01333 1.020 + ( 1.020 − 1.040 ) /3 = 1.01333 ，单位与 Q 相同。
练习
所属知识 扩散稳定性
难度 3/5 扩散率 2.0 × 10 − 6 m 2 s − 1 2.0\times10^{-6}\,\mathrm{m^2\,s^{-1}} 2.0 × 1 0 − 6 m 2 s − 1 0.50 m m 0.50\,\mathrm{mm} 0.50 mm
查看提示 一维显式中心格式要求
α Δ t / Δ x 2 ≤ 1 / 2 \alpha \Delta t/\Delta x^{2}\le 1/2 α Δ t /Δ x 2 ≤ 1/2 。
查看解答 α = 2 × 10 − 6 m 2 / s \alpha=2\times 10^{-6} m^{2}/s α = 2 × 1 0 − 6 m 2 / s 、
Δ x = 0.5 m m \Delta x=0.5 mm Δ x = 0.5 mm 时
Δ t ≤ ( 2.5 × 10 − 7 ) / ( 4 × 10 − 6 ) = 0.0625 s \Delta t\le(2.5\times 10^{-7})/(4\times 10^{-6})=0.0625\,\mathrm{s} Δ t ≤ ( 2.5 × 1 0 − 7 ) / ( 4 × 1 0 − 6 ) = 0.0625 s 。空间加密一半时上限变为 0.015625 s。
练习
所属知识 守恒单位
难度 3/5 写出质量平衡缺陷及缺陷率的单位，并说明净通量接近零时如何报告。
查看提示 存量变化和时间积分通量都是 kg；除以时间后才是 kg/s。
查看解答 δ m = M 2 − M 1 − ∫ ( m ˙ i n − m ˙ o u t ) d t \delta_m=M_2-M_1-\int(\dot{m}_{in}-\dot{m}_{out})dt δ m = M 2 − M 1 − ∫ ( m ˙ in − m ˙ o u t ) d t ，单位 kg；
δ m / Δ t \delta_m/\Delta t δ m /Δ t 单位 kg/s。若累计净通量近零，不宜用它归一，应报告绝对缺陷或选非零特征质量。
练习
所属知识 制造解
难度 4/5 说明制造解的构造步骤，以及它为何不能替代 validation。
查看提示 先选
u m u_m u m ，再由目标算子计算源项与边界值；不要声称它验证自然模型。
查看解答 制造解检查算子、源项、边界和预期收敛阶，即 code verification；因为源项是为所选函数人工构造，它不证明原物理模型与实验一致。
练习
所属知识 误差预算
难度 4/5 给出把离散、迭代、测量和模型误差合并前必须检查的条件，并指出哪些项不宜直接平方和。
查看提示 先统一单位并区分独立标准不确定度、有向偏差和有界区间。
查看解答 独立零均值标准项可平方和开根号；已知偏差应校正或单列；模型区间和强相关项不能当作独立 Gaussian 项。报告需保留来源、单位、相关性和适用条件。
知识连接与资源
课程 · 2015 Numerical Fluid Mechanics Pierre Lermusiaux
用于核对 P11 的离散化误差、场方程网格算法、时间步稳定性、边界条件和数值验证流程。
打开官方来源 MIT OpenCourseWare 2.29 覆盖数值误差、稳定性、网格、时间推进与验证，可用于核对本章离散分析和加密流程。本章示例使用明示教学参数；未给实验来源的数字不构成真实系统 validation 数据。
课程 · 2024 Computational Methods of Scientific Programming Thomas Herring, Chris Hill
用于核对 P11 的程序验证、并行性能测量、环境记录、结果传播和可复现计算实验要求。
打开官方来源 MIT OpenCourseWare 12.010《Computational Methods of Scientific Programming》支撑本章把离散化验证落实为可执行实验：同一物理量要在成组网格或时间步上自动重算，保存误差范数、收敛斜率与输入配置，并用回归测试防止代码修改破坏已知阶数。该课程在此支持的是验证流程和程序证据链；与实验观测的一致性仍需单独的 validation 数据，不能由网格收敛代替。