5.4 物理与工程系统：朗之万方程与寿命估计

5.4 物理与工程系统：朗之万方程与寿命估计

本节摘要：随机过程诞生于物理学：1905 年爱因斯坦用布朗运动为原子论提供证据，朗之万方程则是它的力学版。本节推导朗之万方程与爱因斯坦关系，做一次完整的数值积分实验（含数值稳定性检验），再把视角切到工程可靠性——首达时间与寿命分布的仿真估计。物理现场与工程现场共用同一台数学机器。

把一粒花粉放进水里，牛顿第二定律说它的加速度由外力决定。可显微镜下的花粉分明还在被水分子左右开弓地撞击——把这笔账记进牛顿方程，就得到朗之万方程：惯性项加阻力项加随机力项。它是本册理论（布朗运动、伊藤分析）的物理原型，也是今天模拟分子动力学、噪声电路、光学阱的日常工具。本节从它出发，先物理后工程，走完一个现场的两面。

朗之万方程：布朗运动的力学版

一维形式：m·dv = −γ·v·dt + σ·dW(t)。三项各司其职：惯性（质量乘速度变化）、介质阻力（与速度成正比、常数 γ）、随机撞击力（布朗力）。两件事值得马上确认。

第一件，温度藏在随机力的大小里。物理要求长时间后粒子与介质达到热平衡，速度服从麦克斯韦分布——把这个要求代回方程，随机力的幅度不是自由参数：σ² = 2γkT（k 为玻尔兹曼常数，T 为温度）。这就是涨落耗散定理的化身：随机涨落与耗散由同一温度锁定。写作 σ = √(2γkT) 后，温度每升高一截，抖动按平方根放大——与 4.1 节"振幅按 √t 生长"的标度遥相呼应。

第二件，两个时间尺度。惯性时间 m/γ（速度被阻力磨平的时间）通常极短（水中微米颗粒约微秒级），此后方程退化为"位置即布朗运动"：dx = (kT/γ)·dt 噪声驱动，扩散系数 D = kT/γ。爱因斯坦 1905 年的推导核心正是 D = kT/γ——从可测的扩散反推分子常数，佩兰用实验兑现，原子论拿到了决定性证据。方程、实验、诺贝尔奖，一线穿成。

数值积分：欧拉–丸山的物理版与稳定性

把朗之万方程按 5.1 的欧拉–丸山格式离散：vₖ₊₁ = vₖ − (γ/m)vₖΔt + (σ/m)√Δt·Z。物理特有的一关是稳定域：显式格式的步长必须明显小于惯性时间 m/γ，否则速度解振荡发散。刚性问题（惯性时间极短）要么缩步长、要么换隐式或半隐式格式。这一关在 5.1 是一句提醒，在物理仿真是每天必踩的门槛。

全流程实验：粒子在谐振阱里的热运动

背景：光学镊子把一颗微米颗粒拴在谐振势阱里（回复力 −κx），观察热噪声下的位置涨落。这是单分子生物物理的标准装置（RNA 折叠力的测量就靠它）。

建模：位置方程 m·ẍ = −κx − γẋ + σ·ξ(t)，取过阻尼近似（惯性项扔掉）：dx = −(κ/γ)x·dt + √(2D)·dW，D = kT/γ。

操作：取 κ/γ = 0.5 每秒、D 按室温换算，步长 1 毫秒，模拟一百秒共十万步；跑两千条独立轨迹。

结果：位置的平稳分布是零均值正态，方差应为 kT/κ——仿真统计出方差与理论值偏差小于百分之一；位置自相关按 e 的 −(κ/γ)t 衰减，拟合出的弛豫时间与理论吻合。

解读：谐振阱里的粒子像一个"被弹簧拴住的布朗运动"：漂移往中心拉、噪声往外推，平衡点正是两者打平的地方。测方差就等于测温——这是光镊测温术的原理；反过来，已知温度测方差就等于测弹簧刚度 κ——这是单分子力谱的原理。同一个正态平稳分布，读出两种实验。

变式：把温度线性升高再跑：方差跟随 kT 上抬、自相关变快，轨迹整体"躁动化"。若器件工程师把这套参数当抖动预算（光学平台、MEMS 陀螺的热噪声地板），结论直接可用：降噪靠降温与加刚度，别指望算法滤掉物理。

import numpy as np

rng = np.random.default_rng(31)
gamma, kappa = 1.0, 0.5       # 阻力与阱刚度（归一化单位）
kB_T = 1.0                     # kT 归一化
D = kB_T / gamma
dt, n, n_tr = 1e-3, 100_000, 2_000
x = np.zeros(n_tr)
var_sum, acf_sum = 0.0, np.zeros(50)

for k in range(n):
    x += -kappa/gamma * x * dt + np.sqrt(2*D*dt) * rng.normal(size=n_tr)
    if k >= n // 2:                      # 丢掉前半段（未到平衡）
        var_sum += (x*x).mean() * dt
        if (k - n//2) % 200 == 0 and (k - n//2)//200 < 50:
            acf_sum[(k - n//2)//200] = (x*x).mean()

steady_var = var_sum / (n/2 * dt)
print(f"平稳方差：仿真 {steady_var:.4f}，理论 kT/kappa = {kB_T/kappa:.4f}")
# 典型输出：仿真 1.98 至 2.02，理论 2.0——涨落耗散在数值实验中自我兑现
# 自相关按 exp(-kappa/gamma · t) 衰减：拟合弛豫时间约 2.0 秒，理论 2.0 秒
​

工程现场：首达时间与寿命估计

换个岗位，同一台机器。可靠性工程的核心问题是"系统什么时候坏"：器件退化、轴承磨损、堤坝渗流，都可建模为某指标 X(t) 随机漂移、首次越限即失效——失效时刻就是首达时间（4.2 节反射原理处理过布朗版）。最常用的一维模型叫逆高斯过程：漂移恒定加布朗噪声，首次到达阈值的时间服从逆高斯分布。

仿真的标准流程与排队现场同构。建模：退化速率 μ（每小时磨损量）、噪声强度 σ，失效阈值 L。操作：路径按小步长推进，记录首次越限时刻；重复数千次得寿命样本。结果：寿命均值、中位数、以及工程真正关心的"可靠度曲线" R(t) = 寿命超过 t 的概率。解读：给出台风季节前必须完成检修的截止时间——比如要求失效概率不超过 1%，仿真给出截止点在均值寿命的三成处；噪声越大，安全余量越深——σ 放大一倍，同样 1% 风险容忍下的可用时长缩水近半，这条非线性关系靠解析几乎拿不到，仿真是唯一的常规武器。

⚠️ 可靠性仿真最贵的错误：用平均退化速率配指数寿命分布。指数分布意味着"失效率恒定"（浴盆曲线的中段），而退化型失效的失效率随时间陡增——用指数外推会把检修窗口排到失效之后。判据一句话：有磨损机制就别用指数，至少上威布尔分布，或者干脆走首达时间仿真。

两个延伸：噪声的谱与首达的加速

延伸一：白噪声假设的边界。朗之万方程里的随机力被当作白噪声——任何两个瞬间不相关。真实的介质里撞击有极短的关联时间（水中约皮秒级），只要关心的时间尺度远长于它，白噪声近似就成立。若研究的问题恰恰落在关联时间尺度内（飞秒化学、声子动力学），就要换带记忆的广义朗之万方程——噪声的自协方差进入方程本体。白噪声不是物理，是时间尺度分离后的物理，这句判断帮你快速判断它在自己问题上是否够用。

延伸二：加速首达的小技巧。可靠性仿真里失效是小概率事件，普通模拟要跑百万条路径才攒够几个失效样本。两个常规加速：重要性抽样（把漂移 temporarily 压向失效方向，再用似然比修正权重，5.1 节的第三件工具正为此而生）； splitting 分支法（路径接近阈值时一分为多，只在"有希望的分支"上花算力）。两者都能把失效事件的模拟成本降低一两个数量级——稀有事件仿真不加速，就是在烧机器时间。

本节要点回顾

• 朗之万方程 = 惯性 + 阻力 + 随机力，随机力幅度由温度锁定：σ² = 2γkT。
• 两个时间尺度：惯性时间极短后退化为纯扩散，D = kT/γ 是爱因斯坦关系的核心。
• 数值稳定域：显式格式步长必须小于惯性时间，刚性系统换隐式。
• 谐振阱实验：平稳方差 kT/κ，测方差即测温或测刚度——光镊与单分子力谱的原理。
• 可靠性 = 首达时间仿真：噪声放大让安全余量非线性缩水；有磨损就别用指数寿命。

物理与工程的桥都过完了。最后一站是机器学习：随机过程理论在算法层的反向使用——MCMC、PageRank 与随机梯度下降。