0、前言

博主略懂一些微分方程的求解方法，主要是用有限元方法。但是最近刷B站老给我推送一些PINN神经网络求解微分方程的方法，因此简单研究下是什么高端的技术。

1、PINN 基础理论

在物理信息神经网络（Physics-Informed Neural Networks, PINN）中，我们希望利用神经网络来求解偏微分方程（PDE）。其核心思想是：
让神经网络同时满足 PDE、边界条件、初始条件，以及（可选的）观测数据。

---

1.1、PINN 的基本 PDE 形式

假设 PDE 为：

$$\mathcal{N}[u](x,t)=0,(x,t)\in \Omega$$

其中：

• $u(x,t)$：PDE 的未知解
• $\mathcal{N}[\cdot ]$：微分算子（如 $u_t-u_{xx}$、$u_{tt}-c^2{u}_{xx}$ 等）

神经网络输出记为 $u_{\theta }(x,t)$。

---

1.2、PINN 常见的 Loss 构成

PINN 的总损失一般由以下部分组成。

---

1.2.1、PDE 残差损失 $L_{\text{PDE}}$

$$L_{\text{PDE}}=\frac{1}{N_r}\sum _{i=1}^{N_r}{\left|\mathcal{N}(u_{\theta }(x_i,t_i))\right|}^2$$

其中 $(x_i,t_i)$ 是 PDE 内部点（collocation points）。

---

1.2.2、 边界条件损失 $L_{\text{BC}}$

对于 Dirichlet 边界：

$$L_{\text{BC}}=\frac{1}{N_b}\sum _{i=1}^{N_b}{\left|u_{\theta }(x_i^{(b)},t_i^{(b)})-g(x_i^{(b)},t_i^{(b)})\right|}^2$$

在边界条件中，Neumann 条件不是规定函数值 $u$，而是规定解在边界方向的导数（通量）。例如：

• 一维情况下：
  $$u_x(x,t)=h(x,t)$$

• 多维情况下（法向导数）：
  $$\frac{\partial u}{\partial n}=\nabla u\cdot n=u_x{n}_x+u_y{n}_y(+u_z{n}_z)$$

因此，PINN 中的 Neumann 边界损失为：

$$L_{\text{BC}}=\frac{1}{N_b}\sum _{i=1}^{N_b}{\left|\frac{\partial u_{\theta }}{\partial n}(x_i^{(b)},t_i^{(b)})-h(x_i^{(b)},t_i^{(b)})\right|}^2$$

---

1.2.3、初始条件损失 $L_{\text{IC}}$

用于随时间演化的 PDE（如热方程）：

$$L_{\text{IC}}=\frac{1}{N_0}\sum _{i=1}^{N_0}{\left|u_{\theta }(x_i,0)-u_0(x_i)\right|}^2$$

---

1.2.4、数据点损失（可选）$L_{\text{data}}$

当有实验数据时，可加入监督项：

$$L_{\text{data}}=\frac{1}{N_d}\sum _{i=1}^{N_d}{\left|u_{\theta }(x_i^{(d)},t_i^{(d)})-y_i^{(d)}\right|}^2$$

这一项使 PINN 能执行半监督学习。

---

1.3、总损失函数

总 loss 为加权求和：

$$L={\lambda }_{\text{PDE}}{L}_{\text{PDE}}+{\lambda }_{\text{BC}}{L}_{\text{BC}}+{\lambda }_{\text{IC}}{L}_{\text{IC}}+{\lambda }_{\text{data}}{L}_{\text{data}}$$

不同 $\lambda$ 用于平衡各个损失项，常见调法包括：

• 边界条件权重可以更高
• PDE 残差通常是主导项
• 若有数据点，也需调 ${\lambda }_{\text{data}}$

---

1.4、Loss 实现方法（概念层面）

PINN 的实现依赖深度学习框架的自动微分（PyTorch / TensorFlow）。

① 自动微分求导数

例如通过 autograd 求 $u_x,u_t,u_{xx}$ 等。

② 计算 PDE 残差

将网络输出和求得导数代入 PDE 公式。

③ 构造 MSE 损失

把 PDE、边界、初始条件全部转换成 MSE。

④ 反向传播优化参数

用 Adam 或 LBFGS 等优化器训练神经网络。

---

✔ 总结

PINN 的核心思想就是：把 PDE、边界条件、初始条件变成损失函数，通过自动微分和优化算法逼近 PDE 的解。

2、PINN 代码与解释（完整示例）

本节我们结合一个简单的偏微分方程（PDE）示例，完整展示 PINN 的代码实现过程，并在每一步给出详细解释。

我们要求解的 PDE 是：

$$u^{\prime \prime }(x)+{\pi }^2\sin (\pi x)=0,x\in (0,1)$$

边界条件：

$$u(0)=0,u(1)=0$$

该 PDE 的真实解是：

$$u(x)=\sin (\pi x)$$

下面用 PINN 来逼近这个解析解。

---

2.1、导入库与固定随机种子

import torch
import torch.nn as nn
import torch.optim as optim
import math

from matplotlib import pyplot as plt

# 固定随机种子（方便复现）
torch.manual_seed(0)

• torch.autograd 将用于自动微分计算 $u_x,u_{xx}$
• matplotlib 用来可视化最终的结果

---

2.2、定义 PINN 网络结构

class PINN(nn.Module):
    def __init__(self, hidden_dim=3, num_layers=1):
        super().__init__()
        layers = []
        input_dim = 1
        output_dim = 1

        layers.append(nn.Linear(input_dim, hidden_dim))
        layers.append(nn.Tanh())
        for _ in range(num_layers - 1):
            layers.append(nn.Linear(hidden_dim, hidden_dim))
            layers.append(nn.Tanh())
        layers.append(nn.Linear(hidden_dim, output_dim))

        self.net = nn.Sequential(*layers)

    def forward(self, x):
        return self.net(x)

说明：

• 输入是 $x$，输出是 $u_{\theta }(x)$
• 激活函数使用 tanh，因为 PINN 中要求二阶导数连续
• 网络深度非常浅（只有 1 隐藏层）用于演示 PINN 在简单网络下的表现

---

2.3、构造模型与优化器

model = PINN()
optimizer = optim.Adam(model.parameters(), lr=1e-3)

• 使用 Adam 作为优化器
• 学习率设为 $1\times 10^{-3}$

---

2.4、训练参数设置

num_epochs = 5000
n_residual = 100  # PDE 内部点
n_bc = 2  # 边界点：0 和 1

说明：

• 内部点（collocation points）用于 enforcing PDE
• 边界点只有两个：$x=0$ 和 $x=1$

---

2.5、PINN 训练循环（核心部分）

下面代码逐行解释 PINN 如何构建 PDE 和 BC 的损失。

---

① 构造 PDE 残差点（内部点）

x_r = torch.rand(n_residual, 1, requires_grad=True)
u_r = model(x_r)

• 在 $(0,1)$ 采样 100 个点
• 网络输出 $u_{\theta }(x_r)$

---

② 自动微分求一阶导 $u_x$

u_x = torch.autograd.grad(
    outputs=u_r,
    inputs=x_r,
    grad_outputs=torch.ones_like(u_r),
    create_graph=True,
    retain_graph=True
)[0]

PINN 用自动微分（autograd）求导：

• create_graph=True 允许继续求二阶导
• 得到 $u_x$

---

③ 自动微分求二阶导 $u_{xx}$

u_xx = torch.autograd.grad(
    outputs=u_x,
    inputs=x_r,
    grad_outputs=torch.ones_like(u_x),
    create_graph=True,
    retain_graph=True
)[0]

二阶导数直接使用 autograd 再求一次梯度。

---

④ 构造 PDE 残差

residual = u_xx + (math.pi ** 2) * torch.sin(math.pi * x_r)
loss_pde = torch.mean(residual ** 2)

这是 PDE：

$$u_{xx}(x)+{\pi }^2\sin (\pi x)=0$$

PINN 的做法：

• 把 PDE 当成一个“误差”
• 用 MSE 让残差尽量为 0

---

2.6、构造边界条件损失

x_bc = torch.tensor([[0.0], [1.0]])
u_bc = model(x_bc)
target_bc = torch.zeros_like(u_bc)

loss_bc = torch.mean((u_bc - target_bc) ** 2)

对应边界条件：

$$u(0)=0,u(1)=0$$

使用 MSE 强制网络满足该约束。

---

2.7、总 Loss 与反向传播

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

即：

$$L=L_{\text{PDE}}+L_{\text{BC}}$$

PINN 训练本质上就是最小化这个 loss。

---

2.8、打印训练日志

if (epoch + 1) % 500 == 0:
    print(f"...")

用于监控：

• 总损失
• PDE 残差
• 边界条件残差

---

2.9、模型测试：与真实解对比

x_test = torch.linspace(0, 1, 5).view(-1, 1)
u_pred = model(x_test)
u_true = torch.sin(math.pi * x_test)

打印：

• $u_{\theta }(x)$ 预测值
• 真实解 $\sin (\pi x)$

用于验证 PINN 是否学到了 PDE 结构。

---

2.10、绘制对比图（PINN vs 真解）

x_plot = torch.linspace(0, 1, 200).view(-1, 1)
u_pred = model(x_plot)
u_true = torch.sin(math.pi * x_plot)

plt.plot(x_np, u_true_np, label='True u(x) = sin(pi x)')
plt.plot(x_np, u_pred_np, '--', label='PINN prediction')
plt.legend()
plt.xlabel('x')
plt.ylabel('u(x)')
plt.title('PINN vs True Solution')
plt.show()

结果：

• 蓝色曲线是真解
• 虚线曲线是 PINN 拟合结果
  

尽管网络非常小，PINN 仍能准确复现 PDE 的解析解形状。

3、全部代码

import torch
import torch.nn as nn
import torch.optim as optim
import math

from matplotlib import pyplot as plt

# 固定随机种子（方便复现）
torch.manual_seed(0)


# 简单的全连接网络，输入 x，输出 u(x)
class PINN(nn.Module):
    def __init__(self, hidden_dim=3, num_layers=1):
        super().__init__()
        layers = []
        input_dim = 1
        output_dim = 1

        layers.append(nn.Linear(input_dim, hidden_dim))
        layers.append(nn.Tanh())
        for _ in range(num_layers - 1):
            layers.append(nn.Linear(hidden_dim, hidden_dim))
            layers.append(nn.Tanh())
        layers.append(nn.Linear(hidden_dim, output_dim))

        self.net = nn.Sequential(*layers)

    def forward(self, x):
        # x: [N, 1]
        return self.net(x)


# 构造网络和优化器
model = PINN()
optimizer = optim.Adam(model.parameters(), lr=1e-3)

# 训练参数
num_epochs = 5000
n_residual = 100  # PDE 内部点
n_bc = 2  # 边界点：0 和 1

for epoch in range(num_epochs):
    optimizer.zero_grad()

    # ===== 1. PDE 残差点（内部点） =====
    # 在 (0,1) 随机采样
    x_r = torch.rand(n_residual, 1, requires_grad=True)  # [N, 1]

    u_r = model(x_r)  # u(x)

    # 一阶导 u_x
    u_x = torch.autograd.grad(
        outputs=u_r,
        inputs=x_r,
        grad_outputs=torch.ones_like(u_r),
        create_graph=True,
        retain_graph=True
    )[0]

    # 二阶导 u_xx
    u_xx = torch.autograd.grad(
        outputs=u_x,
        inputs=x_r,
        grad_outputs=torch.ones_like(u_x),
        create_graph=True,
        retain_graph=True
    )[0]

    # PDE 残差：u_xx + pi^2 * sin(pi x) = 0
    residual = u_xx + (math.pi ** 2) * torch.sin(math.pi * x_r)

    # PDE loss
    loss_pde = torch.mean(residual ** 2)

    # ===== 2. 边界条件点（x=0,1） =====
    x_bc = torch.tensor([[0.0], [1.0]])  # [2,1]
    u_bc = model(x_bc)  # 应该接近 0
    target_bc = torch.zeros_like(u_bc)

    loss_bc = torch.mean((u_bc - target_bc) ** 2)

    # ===== 3. 总 loss =====
    loss = loss_pde + loss_bc

    # 反向传播和更新
    loss.backward()
    optimizer.step()

    # 每隔一段时间打印一下
    if (epoch + 1) % 500 == 0:
        print(f"Epoch {epoch + 1}/{num_epochs} | "
              f"Loss: {loss.item():.4e} | "
              f"PDE: {loss_pde.item():.4e} | "
              f"BC: {loss_bc.item():.4e}")

# 训练结束后，简单测试一下某些点的误差（与真解 sin(pi x) 比较）
with torch.no_grad():
    x_test = torch.linspace(0, 1, 5).view(-1, 1)
    u_pred = model(x_test)
    u_true = torch.sin(math.pi * x_test)
    print("\n测试点对比：")
    for x_val, up, ut in zip(x_test, u_pred, u_true):
        print(f"x={x_val.item():.2f}, u_pred={up.item():.4f}, u_true={ut.item():.4f}")

with torch.no_grad():
    x_plot = torch.linspace(0, 1, 200).view(-1, 1)
    u_pred = model(x_plot)
    u_true = torch.sin(math.pi * x_plot)

x_np = x_plot.numpy().flatten()
u_pred_np = u_pred.numpy().flatten()
u_true_np = u_true.numpy().flatten()

plt.plot(x_np, u_true_np, label='True u(x) = sin(pi x)')
plt.plot(x_np, u_pred_np, '--', label='PINN prediction')
plt.legend()
plt.xlabel('x')
plt.ylabel('u(x)')
plt.title('PINN vs True Solution')
plt.show()