1. 一维热传导方程的热扩散系数辨识

任务目标

已知一维热传导方程 u_t = D·u_xx 在区间 [0, 1] 上的温度分布观测数据，利用物理信息神经网络（PINN）反演未知热扩散系数 D，并评估辨识精度。本案例采用合成数据验证：先给定真值 D_true = 0.4 生成观测数据，再用 PINN 反演并与真值对比。

问题设定与数据生成

考虑如下初边值问题：

u_t = D·u_xx,  x ∈ [0, 1],  t ∈ [0, 1]
u(x, 0) = sin(πx)
u(0, t) = u(1, t) = 0

该问题有解析解 u(x,t) = sin(πx)·exp(-Dπ²t)。在 D_true = 0.4 下，在 x 方向取 21 个等距点、t 方向取 21 个等距点，共 441 个时空点，计算解析解作为“观测数据”（不额外加噪，以便清晰评估辨识精度）。

PINN 结构与损失函数

神经网络 N(x,t; θ) 输出近似解 u_θ(x,t)。损失函数由三部分组成：

• 数据损失：MSE_data = mean((u_θ(x_data, t_data) - u_obs)²)，约束网络输出匹配观测温度；
• 方程残差损失：MSE_pde = mean((∂u_θ/∂t - D·∂²u_θ/∂x²)²)，在配点（本案例使用与数据点相同的网格点）上强制满足热传导方程；
• 初边值损失：MSE_ic = mean((u_θ(x, 0) - sin(πx))²)，MSE_bc = mean(u_θ(0, t)²) + mean(u_θ(1, t)²)。

总损失为 L = MSE_data + MSE_pde + MSE_ic + MSE_bc。其中 D 作为可训练参数，与网络权重一起通过梯度下降优化。

关键完整代码

import numpy as np
import torch
import torch.nn as nn

# ========== 1. 生成合成观测数据（解析解） ==========
D_true = 0.4
x = np.linspace(0, 1, 21)
t = np.linspace(0, 1, 21)
X, T = np.meshgrid(x, t)
x_flat = X.flatten()[:, None]
t_flat = T.flatten()[:, None]
u_exact = np.sin(np.pi * x_flat) * np.exp(-D_true * np.pi**2 * t_flat)

# 转为 tensor
x_tensor = torch.tensor(x_flat, dtype=torch.float32, requires_grad=False)
t_tensor = torch.tensor(t_flat, dtype=torch.float32, requires_grad=False)
u_tensor = torch.tensor(u_exact, dtype=torch.float32)

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

# ========== 3. 初始化模型与可训练参数 D ==========
model = PINN()
D_log = torch.tensor([0.0], dtype=torch.float32, requires_grad=True)  # D = exp(D_log) > 0

optimizer = torch.optim.Adam(list(model.parameters()) + [D_log], lr=1e-3)

# ========== 4. 训练 ==========
def pde_residual(x, t):
    x = x.clone().requires_grad_(True)
    t = t.clone().requires_grad_(True)
    u = model(x, t)
    u_t = torch.autograd.grad(u, t, grad_outputs=torch.ones_like(u), create_graph=True)[0]
    u_x = torch.autograd.grad(u, x, grad_outputs=torch.ones_like(u), create_graph=True)[0]
    u_xx = torch.autograd.grad(u_x, x, grad_outputs=torch.ones_like(u), create_graph=True)[0]
    D = torch.exp(D_log)
    return u_t - D * u_xx

# 初边值配点
x_ic = torch.tensor(x_flat, dtype=torch.float32)
t_ic = torch.zeros_like(x_ic)
u_ic = torch.sin(np.pi * x_ic)

x_bc0 = torch.zeros_like(t_tensor)
x_bc1 = torch.ones_like(t_tensor)

# 训练循环
epochs = 3000
for epoch in range(epochs):
    optimizer.zero_grad()

    # 数据损失
    u_pred = model(x_tensor, t_tensor)
    loss_data = torch.mean((u_pred - u_tensor) ** 2)

    # PDE 残差损失
    res = pde_residual(x_tensor, t_tensor)
    loss_pde = torch.mean(res ** 2)

    # 初值损失
    u_ic_pred = model(x_ic, t_ic)
    loss_ic = torch.mean((u_ic_pred - u_ic) ** 2)

    # 边值损失
    u_bc0 = model(x_bc0, t_tensor)
    u_bc1 = model(x_bc1, t_tensor)
    loss_bc = torch.mean(u_bc0 ** 2) + torch.mean(u_bc1 ** 2)

    loss = loss_data + loss_pde + loss_ic + loss_bc
    loss.backward()
    optimizer.step()

    if (epoch + 1) % 500 == 0:
        D_cur = torch.exp(D_log).item()
        print(f"Epoch {epoch+1}, Loss: {loss.item():.6e}, D = {D_cur:.4f}")

# ========== 5. 输出辨识结果与精度评估 ==========
D_estimated = torch.exp(D_log).item()
print(f"\n真值 D_true = {D_true}")
print(f"辨识值 D_estimated = {D_estimated:.4f}")
print(f"相对误差 = {abs(D_estimated - D_true) / D_true * 100:.2f}%")

运行方式与依赖

需要安装 numpy 和 torch（CPU 版本即可）。将上述代码保存为 pinn_diffusion.py，在终端运行：

python pinn_diffusion.py

验证方法与通过标准

训练约 3000 轮后，观察终端输出的 D 值。判定通过的标准：

• 损失函数下降至约 1e-5 以下（无噪声数据下可达到）；
• 辨识出的 D_estimated 与真值 0.4 的相对误差在 5% 以内。

由于本案例使用无噪声解析解数据且配点覆盖整个时空域，PINN 能够以较高精度反演扩散系数。若增加噪声或减少配点，误差会相应增大，属正常现象。

常见问题

• D 不收敛或为负：通过对数参数化 D = exp(D_log) 保证正性，若仍不收敛可尝试增大训练轮数或调整学习率。
• 训练不稳定：可尝试增加网络宽度/深度，或使用学习率调度器（如 torch.optim.lr_scheduler.StepLR）。
• 初值损失权重过大：若初边值条件与数据点重叠较多，可适当降低初边值损失权重（如乘以 0.1），但本案例中各项权重均为 1 即可收敛。