3. 不确定性传播的基本概念

在工程实践中，我们很少直接面对“输出”的不确定性，而是先面对一堆“输入”的不确定性：材料强度的离散、几何尺寸的加工公差、载荷的随机波动、边界条件的模糊。不确定性传播（Uncertainty Propagation）研究的正是这样一个问题——当输入被看作随机变量或随机场时，它们的不确定性如何经过一个确定性的物理模型或计算模型，传递、放大、衰减，最终体现在输出量的分布上。理解这条“从输入到输出”的通道，是后续所有量化方法（抽样、谱方法、代理模型）的共同出发点。

一、输入不确定性如何影响输出

设一个工程模型可以抽象为映射 y = g(x)，其中输入向量 x = (x₁, x₂, …, x_d) 带有不确定性，输出 y 是我们关心的响应量（如最大应力、固有频率、疲劳寿命）。输入的不确定性对输出的影响，通常通过三种机制体现：

• 敏感度差异：不同输入对输出的影响权重不同。若 g 对某个输入偏导很大，该输入的微小波动就会被显著放大；反之则被“稀释”。
• 非线性放大：当模型存在阈值、屈曲、接触等非线性环节时，输入的对称分布可能产生严重偏斜的输出分布，均值也可能偏离名义值。
• 相关性耦合：输入之间若存在相关（如材料强度与弹性模量同向变化），其联合不确定性对输出的影响不能简单按独立情形叠加。

一个直观的示例：某简支梁跨中挠度 δ = 5qL⁴/(384EI)。若长度 L 的相对标准差为 1%，则挠度对其的敏感度约为 4 倍，即挠度相对标准差约 4%；而弹性模量 E 的相对标准差 5% 只带来约 5% 的挠度波动。这说明“输入波动大小”与“输出波动大小”之间并非等比例，敏感度指数才是关键。

二、传播问题的数学表述

把上面的直观认识形式化。设输入随机向量 X = (X₁, …, X_d) 具有联合概率密度 f_X(x)，模型 Y = g(X) 为确定性映射。不确定性传播要回答的核心问题是：求输出 Y 的分布 f_Y(y)，或至少求其若干统计矩与分位数。

最一般的表述是分布变换。当 g 可逆且单调时，可由变量替换公式得到

f_Y(y) = f_X(g⁻¹(y)) · |d g⁻¹(y) / dy|

但工程模型几乎不可能满足可逆单调，因此实际中更常用的是矩传播与概率积分。一阶二阶矩法（泰勒展开）在名义点 x₀ 处线性化：

μ_Y ≈ g(x₀) + ½ Σᵢ Σⱼ (∂²g/∂xᵢ∂xⱼ) Cov(Xᵢ, Xⱼ)
σ_Y² ≈ Σᵢ Σⱼ (∂g/∂xᵢ)(∂g/∂xⱼ) Cov(Xᵢ, Xⱼ)

若输入相互独立，协方差项为零，方差退化为各偏导平方与方差乘积之和。这正是“敏感度平方加权”的数学来源。需要提醒的是，一阶矩法在强非线性或大变异系数下误差明显，它给出的是近似，不能替代完整分布。

更完整的表述是失效概率积分：设极限状态函数 Z = g(X)，失效事件为 Z ≤ 0，则

P_f = ∫_{g(x) ≤ 0} f_X(x) dx

这个 d 维积分通常没有解析解，正是各种数值传播方法（蒙特卡洛、一次二阶矩、响应面、多项式混沌）要解决的共同对象。理解这一点，就能明白不同方法只是对同一个积分采取不同的逼近策略。

三、蒙特卡洛模拟的基本思想与简单实现

蒙特卡洛方法的思想极其朴素：既然不知道输出的分布，那就“造”大量符合输入分布的样本，逐个代入模型算出输出，再用输出样本的经验分布去逼近真实分布。大数定律保证，当样本量 N 足够大时，经验均值依概率收敛到真实均值，经验分布函数逐点收敛到真实分布函数。

其误差具有统计性质：均值估计的标准误约为 σ_Y/√N，因此样本量每增加 4 倍，误差才减半。这也是蒙特卡洛“稳健但昂贵”的根源——它对维度不敏感，但收敛速度慢，对高可靠度问题（P_f 很小）需要极大样本量。

下面给出一个可直接运行的简单实现，用 Python 对简支梁挠度做传播分析。示例参数为教学假设值，非实测数据：

import numpy as np

rng = np.random.default_rng(2024)
N = 200000

# 输入不确定性（示例假设）
L = rng.normal(6.0, 0.06, N)          # 跨度 m，变异系数 1%
E = rng.normal(2.1e11, 1.05e10, N)    # 弹性模量 Pa，变异系数 5%
q = rng.lognormal(np.log(12e3), 0.15, N)  # 均布荷载 N/m
I = rng.normal(8.5e-5, 1.7e-6, N)     # 惯性矩 m^4，变异系数 2%

# 确定性模型：简支梁跨中挠度
delta = 5 * q * L**4 / (384 * E * I)

print("均值 μ =", delta.mean())
print("标准差 σ =", delta.std())
print("变异系数 =", delta.std() / delta.mean())
print("5% 分位 =", np.percentile(delta, 5))
print("95% 分位 =", np.percentile(delta, 95))
print("P(delta > L/250) =", (delta > L/250).mean())

这段代码的关键点有三处。第一，输入抽样必须与设定的分布一致：正态用 normal，对数正态用 lognormal 并注意其参数是底层正态的均值与标准差，而非直接给定期望。第二，模型计算是逐样本向量化的，避免循环以提升效率。第三，结果判读不只看均值，还要看分位数与超限概率——工程上“挠度超过 L/250”的概率往往比均值更有决策意义。

运行后（示例结果，仅用于说明判读方式）可能得到挠度均值约 14.6 mm、标准差约 1.1 mm、变异系数约 7.5%，95% 分位约 16.4 mm，超限概率约 0.3%。可以看到，尽管各输入变异系数只有 1%～5%，输出变异系数却达到 7.5%，这正是长度四次方敏感度放大的结果。若把 L 的变异系数从 1% 降到 0.5%，输出变异系数会明显下降，这为“控制哪个输入最有效”提供了定量依据。

至此，不确定性传播的链条已经清晰：输入分布经敏感度与非线性映射，形成输出分布；数学上它归结为分布变换、矩近似或失效概率积分；蒙特卡洛则用“抽样—计算—统计”三步给出通用而直观的数值答案。掌握这条主线，后续的一次二阶矩、响应面与多项式混沌，都可以理解为在“精度”与“计算成本”之间对同一问题做出的不同权衡。