4. 深度集成PINN用于反应扩散方程

任务目标

训练多个随机初始化的PINN（物理信息神经网络）组成深度集成，求解一维反应扩散方程，输出解的均值预测与预测方差，并分析集成大小（成员数）对不确定性估计的影响。

问题定义与必要依赖

考虑一维反应扩散方程（Fisher方程）：

∂u/∂t = D·∂²u/∂x² + r·u·(1-u)，x∈[0,1]，t∈[0,1]

初值：u(x,0) = exp(-50(x-0.5)²)；边界：u(0,t)=u(1,t)=0。取 D=0.01，r=1.0。

依赖：torch、numpy、matplotlib（可选，用于绘图）。

核心实现

每个集成成员使用相同的网络架构（4层全连接，每层50个神经元，tanh激活），但采用不同的随机种子初始化权重。每个成员独立训练，损失函数为PINN残差损失（反应扩散方程残差）+ 初值损失 + 边界损失。训练完成后，对同一查询点集合，所有成员各自输出预测，计算均值与方差。

import torch
import numpy as np

torch.manual_seed(42)
np.random.seed(42)

# 问题参数
D = 0.01
r = 1.0
N_ENSEMBLE = 10  # 集成大小，可修改为 1,3,5,10 进行对比

# 网络定义
class PINN(torch.nn.Module):
    def __init__(self):
        super().__init__()
        self.net = torch.nn.Sequential(
            torch.nn.Linear(2, 50), torch.nn.Tanh(),
            torch.nn.Linear(50, 50), torch.nn.Tanh(),
            torch.nn.Linear(50, 50), torch.nn.Tanh(),
            torch.nn.Linear(50, 1)
        )
    def forward(self, x, t):
        xt = torch.cat([x, t], dim=1)
        return self.net(xt)

# 采样点生成
N_f = 2000   # 内部残差点
N_b = 200    # 边界点
N_0 = 300    # 初始点

x_f = torch.rand(N_f, 1)
t_f = torch.rand(N_f, 1)
x_b0 = torch.zeros(N_b, 1); t_b0 = torch.rand(N_b, 1)
x_b1 = torch.ones(N_b, 1);  t_b1 = torch.rand(N_b, 1)
x_0 = torch.rand(N_0, 1);   t_0 = torch.zeros(N_0, 1)
u_0 = torch.exp(-50*(x_0-0.5)**2)

# 自动微分求导工具
def grad(u, x, t, order=1):
    g = torch.autograd.grad(u, [x, t], grad_outputs=torch.ones_like(u), create_graph=True)[0]
    if order == 1:
        return g
    u_x = g[:, 0:1]
    u_xx = torch.autograd.grad(u_x, x, grad_outputs=torch.ones_like(u_x), create_graph=True)[0]
    return u_x, u_xx

# 训练单个成员
def train_member(seed):
    torch.manual_seed(seed)
    model = PINN()
    optimizer = torch.optim.Adam(model.parameters(), lr=1e-3)
    # 残差点需要梯度
    x_f.requires_grad_(True); t_f.requires_grad_(True)
    for epoch in range(3000):
        optimizer.zero_grad()
        # 残差损失
        u_f = model(x_f, t_f)
        u_x, u_xx = grad(u_f, x_f, t_f)
        u_t = torch.autograd.grad(u_f, t_f, grad_outputs=torch.ones_like(u_f), create_graph=True)[0]
        f = u_t - D*u_xx - r*u_f*(1-u_f)
        loss_f = torch.mean(f**2)
        # 边界损失
        u_b0 = model(x_b0, t_b0); u_b1 = model(x_b1, t_b1)
        loss_b = torch.mean(u_b0**2) + torch.mean(u_b1**2)
        # 初始损失
        u_0_pred = model(x_0, t_0)
        loss_0 = torch.mean((u_0_pred - u_0)**2)
        loss = loss_f + loss_b + loss_0
        loss.backward()
        optimizer.step()
    return model

# 训练集成
ensemble = [train_member(seed) for seed in range(N_ENSEMBLE)]

# 在网格上预测
x_grid = torch.linspace(0, 1, 100).repeat(100, 1).reshape(-1, 1)
t_grid = torch.linspace(0, 1, 100).repeat_interleave(100).reshape(-1, 1)

with torch.no_grad():
    preds = torch.stack([m(x_grid, t_grid).squeeze() for m in ensemble])  # shape: [N_ENSEMBLE, 10000]
mean_pred = preds.mean(dim=0).reshape(100, 100)
var_pred = preds.var(dim=0).reshape(100, 100)

# 输出统计摘要（确定性输出）
print(f"集成大小: {N_ENSEMBLE}")
print(f"预测均值范围: [{mean_pred.min().item():.4f}, {mean_pred.max().item():.4f}]")
print(f"预测方差范围: [{var_pred.min().item():.6f}, {var_pred.max().item():.6f}]")
print(f"平均预测方差: {var_pred.mean().item():.6f}")

运行与验证

运行方式：直接执行上述脚本。修改 N_ENSEMBLE 为 1、3、5、10 分别运行，记录每次输出的平均预测方差。

验证方法：

  1. 解的正确性：观察 mean_pred 是否呈现从初始高斯峰向两侧扩散并趋于均匀分布 u≈1 的趋势（反应扩散方程稳态解为1）。在 t=0 处，预测均值应接近初始条件 exp(-50(x-0.5)²)。
  2. 不确定性行为：方差在解变化剧烈的区域（初始陡峭梯度处）应较大，在解趋于平滑后减小。
  3. 集成大小影响：随着 N_ENSEMBLE 从1增加到10，平均预测方差应逐渐减小并趋于稳定，表明集成降低了预测不确定性。

通过条件：当 N_ENSEMBLE=1 时方差为0（单模型无集成方差）；N_ENSEMBLE=10 时平均方差显著小于 N_ENSEMBLE=3 的情况，且预测均值在物理上合理（不出现负值或大于1.2的异常值）。

集成大小对不确定性的影响分析

深度集成的预测方差由两部分构成：数据不确定性（噪声）与模型不确定性（参数随机性）。由于本问题为确定性方程，方差主要反映模型不确定性——即不同随机初始化导致的预测分歧。集成大小增加时，极端预测被平均稀释，方差下降；但下降速率递减，说明存在收益递减点。实际应用中，可通过绘制“集成大小 vs 平均方差”曲线，选择方差曲线拐点处的集成大小作为计算资源与不确定性精度的平衡点。

集成大小（成员数） 平均预测方差 1 3 5 8 10 方差快速下降区 收益递减区 拐点≈5

上图展示了集成大小与平均预测方差的典型关系：当成员数从1增至5时，方差显著下降；超过5后下降趋于平缓。实际应用中建议至少使用5个成员以平衡计算成本与不确定性估计质量。

← 上一章 📖 返回目录 下一章 →