ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

PINN求解微分方程的Python实战:损失函数与自动微分解析

PINN求解微分方程的Python实战:损失函数与自动微分解析 简介面向Python开发者与科学计算研究者这份基于PINN的微分方程求解方法资源包系统展示了物理学信息神经网络在常微分方程、偏微分方程及PDAE系统等场景中的落地实现。内容覆盖欧拉梁方程、扩散方程、拉普拉斯方程、泊松方程、Lorenz系统等典型问题并给出DeepXDE框架的Poisson方程Dirichlet边界、Neumann边界、Robin边界及周期边界等多种案例适合希望将深度学习与物理建模结合的入门及进阶学习者。包内共26个文件以17个ipynb教学笔记为核心辅以3个py源码脚本、4个zbak备份文件以及1个png效果图和1个md说明文档整体约891KB。通过自动微分与损失函数设计读者可直接运行并理解数据驱动与物理约束结合的完整求解流程也可参考Jacobian-Hessian方法对比实验用于改进自身模型训练与数值精度。已有217人学习适合作为PINN方向快速上手的实操参考。1. PINN求解微分方程的Python实战资源一个能直接跑起来的算例包做力学仿真或者热传导分析的朋友应该都有体会碰上一个没有解析解的偏微分方程传统有限元要画网格、定边界、调求解器一套流程下来半天就没了。我拆完这套基于PINNPhysics-Informed Neural Networks物理信息神经网络的Python资源包之后最大的感受是它把用神经网络求解微分方程这件事变成了十几个可以直接运行、直接改参数的notebook算例。里面既有ysin(πx)·cos(πx)这种简单的一阶ODE也有Poisson方程在Dirichlet、Neumann、Robin、Periodic四种边界条件下的变体还有Lorenz系统、扩散方程、欧拉梁方程以及DeepXDE封装的对比版本。适合谁用两类人最合适一类是刚接触PINN、想看看损失函数和边界条件到底怎么写的新手另一类是已经在用有限元做数值计算、想对比一下PINN与传统解法差异的工程师。这套资源不需要你从零搭框架改方程、改边界、改网络宽度跑通一个算例基本就能迁移到自己的问题上。2. PINN代码骨架损失函数四项、自动微分与Adam/LBFGS训练节奏2.1 PINN为什么能解微分方程从物理约束到损失函数四项先抛开复杂的数学把PINN的思路压缩成一句话用一个全连接神经网络 u(x;θ) 去逼近微分方程的解然后让这个网络既满足方程本身又满足边界条件和初始条件。训练的时候不是拿真实解当标签而是把方程左边减右边等于零当成目标这就是物理约束进入神经网络的方式。资源包里PDE.py、model.py、geometry.py这几个文件就是把这件事拆成了三块PDE.py负责定义方程和边界条件model.py负责定义网络结构geometry.py负责在求解域里采样。对于一维问题x是网络输入输出是u(x)对于扩散方程这种时间依赖问题输入变成(x,t)输出还是u(x,t)时间维度和空间维度一样作为坐标处理。损失函数也可以拆成几项来看。对于常微分方程一般形式残差项ODE residualR du/dx - f(x)加上边界条件项如果有观测数据再加一个数据项。资源包里的算例绝大多数只用前两项少数带解析解的会额外加数据匹配。实现的时候关键是把这些约束写成Python里可微的损失函数下面给一个PyTorch版本的简化骨架import torch import torch.nn as nn # 1. 网络结构3层全连接每层50个神经元tanh激活 class FNN(nn.Module): def __init__(self, layers[1, 50, 50, 50, 1]): super().__init__() self.net nn.Sequential() # 中间层线性 tanh交替堆叠 for i in range(len(layers)-2): self.net.add_module(flinear_{i}, nn.Linear(layers[i], layers[i1])) self.net.add_module(fact_{i}, nn.Tanh()) # 输出层不加激活 self.net.add_module(linear_out, nn.Linear(layers[-2], layers[-1])) def forward(self, x): return self.net(x) model FNN([1, 50, 50, 50, 1]) # 2. 残差采样点在[-1, 1]区间均匀取100个点 x_pde torch.linspace(-1, 1, 100, requires_gradTrue).reshape(-1, 1) # 边界点 x_bc torch.tensor([[-1.0], [1.0]]) # 3. 用自动微分计算一阶导 u model(x_pde) u_x torch.autograd.grad(u, x_pde, grad_outputstorch.ones_like(u), create_graphTrue)[0] # 4. 以 y sin(pi*x)*cos(pi*x) 为例 f torch.sin(torch.pi * x_pde) * torch.cos(torch.pi * x_pde) residual u_x - f loss_pde torch.mean(residual**2) # 5. 边界损失y(-1)0 u_bc model(x_bc) loss_bc torch.mean((u_bc[0] - 0.0)**2) total_loss loss_pde loss_bc这里最关键的一行是torch.autograd.grad。PyTorch会自动从 u 到 x 做反向传播求梯度所以不需要手推导数表达式。注意create_graphTrue这个参数它表示把求导过程本身也记录进计算图这样后续如果还需要对一阶导再求导比如二阶ODE、扩散方程就能继续用grad叠加。如果你去掉这个参数第一次求导后再对u_x求导就会直接报错。网络结构方面资源包里这些算例用的都是全连接网络加tanh激活。为什么普遍选tanh而不是ReLU因为ReLU的导数在零点不连续而PINN的损失函数里全是导数项激活函数不够光滑计算高阶导数时误差会变大。tanh的导数是连续的二阶导也存在配合自动微分不会有ReLU那种断裂的问题。层数一般3到5层每层50到100个神经元这个配置在资源包里跑得比较稳。2.2 从PDE.py到model.py这套资源的代码结构与关键参数打开资源包先别急着跑notebook把根目录下的文件扫一遍。README.md里面写了每个算例对应的方程和边界条件PDE.py、model.py、geometry.py这三个是公共模块会被notebook引用带.zbak后缀的其实是备份文件比如PDE.py.zbak、model.py.zbak说明作者在调整过程中存了旧版本。你如果新写的代码报错可以跟备份文件对比一下差异排查改动点。文件用途我整理了一个表方便你快速定位文件作用里面有什么PDE.py定义方程类型、残差计算、边界条件方程名、PDE残差函数、初边值定义model.py定义神经网络结构与初始化网络层数、每层宽度、激活函数、权重初始化方式geometry.py求解域采样区间端点、采样点数、随机种子README.md算例索引与运行说明每个notebook对应方程和边界条件notebook完整训练流程数据生成、损失构建、优化器、训练循环、绘图一个比较常见的调整思路是改方程只动PDE.py改网络只动model.py改采样只动geometry.py。我拆的时候试过把某个算例的采样点数从100加到500损失下降曲线明显更平滑但训练时间也翻了将近五倍。PINN这里有个本质矛盾残差点太少方程约束学不充分残差点太多每次迭代算损失都慢。常见做法是先300个点左右起步看损失能不能压到目标量级再决定要不要加密。2.3 训练循环与优化器选择Adam与LBFGS的配合节奏optimizer torch.optim.Adam(model.parameters(), lr1e-3) for epoch in range(5000): optimizer.zero_grad() loss compute_total_loss(model, x_pde, x_bc) # 由PDE.py组装 loss.backward() optimizer.step() if epoch % 500 0: print(epoch, loss.item())这段代码几乎在每个notebook里都能看到但真正决定能不能收敛的往往是后段操作。资源包里好几个算例是先Adam后LBFGS的组合Adam先跑两三千步把损失粗降到1e-3左右然后换成LBFGS做精调一般能再往下压两三个数量级。原因在于Adam对学习率不敏感适合前期从乱序的损失面里找到大致方向而LBFGS利用二阶信息在损失已经比较低的时候收敛更稳。切优化器时有一个容易踩的坑LBFGS在PyTorch里的用法和Adam不一样它不是每次手动zero_grad而是把损失函数和epoch数都传进去由内部循环控制。常见写法是def closure(): optimizer.zero_grad() loss compute_total_loss(model, x_pde, x_bc) loss.backward() return loss optimizer torch.optim.LBFGS(model.parameters(), lr0.5, max_iter20) for inner in range(100): optimizer.step(closure)这里的lr0.5在LBFGS里对应的是初始步长缩放跟Adam的1e-3完全不是一个量级别照搬。切换的时候建议先把Adam的训练轮数固定下来记录一下切LBFGS之前的loss如果LBFGS跑20步后不降反升说明Adam阶段还没收敛到合适的碗底把Adam轮数再加长。提示切换优化器前先把当前loss打印出来存档。如果LBFGS阶段loss没有低于这个存档值说明切换时机不对回退到Adam继续跑。3. 手把手复现三类算例ODE、PDE与特征值问题的参数设置3.1 常微分方程从ysin(πx)到二阶ODE的边界条件处理先从一个最简单的一阶ODE开始。文件ysin(πx)_cos(πx), x∈[-1,1], y(-1)0.ipynb方程是 y sin(πx)·cos(πx)边界条件是 y(-1)0。这个方程的解析解可以手推sin(πx)cos(πx) (1/2)sin(2πx)积分后 y sin²(πx)/(2π) C代入 y(-1)0 得 C0所以 y sin²(πx)/(2π)。这个解析解可以直接用来画对比图验证网络预测是不是对的。跑这个算例的时候我一般先把问题简化到极致网络就一层50个神经元2000步Adam不切LBFGS看基本趋势能不能对上。如果这都跑不对先去查梯度计算和边界条件代码而不是急着调网络宽度。注意这里有个细节边界点只有一个约束很弱如果残差损失和边界损失不加权重网络经常会学出一个整体偏移的解就是形状对但整体上下平移。这时候把边界权重调大10倍偏移问题基本就消失了。第二个ODE是y-y-y2x, x∈[-5,5], y(-5)1, y(5)5.ipynb。这里面最大的变化是从一阶变二阶损失函数里同时出现 u 和 u。在PyTorch里算二阶导的标准套路是嵌套grad先对u求一阶导得到u_x再对u_x求一次grad得到u_xx。代码段如下u model(x) u_x torch.autograd.grad(u, x, grad_outputstorch.ones_like(u), create_graphTrue)[0] u_xx torch.autograd.grad(u_x, x, grad_outputstorch.ones_like(u_x), create_graphTrue)[0] # 残差u - u - u - 2x 0 residual u_xx - u_x - u - 2.0 * x loss_pde torch.mean(residual**2)注意这里有个细节第一次grad必须设create_graphTrue否则u_x就是普通张量而不是计算图的一部分第二次grad会直接报 element 0 of tensors does not require grad。资源包里的文件没有特别注释这一点我第一次复现时就在这里栽了跟头。二阶ODE的边界条件在x∈[-5,5]的两端采样时记得把x_pde的范围扩到[-5,5]不然方程约束根本覆盖不到边界附近。第三个-yπ²sin(πx), -1≤x≤1, y(-1)0, y(1)0.ipynb是带系数和齐次Dirichlet边界条件的二阶问题它的解析解就是 u sin(πx)。这个算例的一个技巧是把网络输出直接构造成在边界上自动等于零的形式比如令 u(x) (x1)(x-1)·N(x)这样无论N输出什么u(±1) 都等于0边界损失就可以直接从损失函数里删掉。这种hard constraint做法能显著降低边界附近误差但注意它改变了网络输出与解之间的映射关系网络要学习的函数从 u(x) 变成了 N(x)如果N(x)在边界附近变化剧烈不一定比soft constraint更稳。资源包里还有一个Euler Beam.ipynb是四阶梁方程要对u连续求四阶导写法就是嵌套四层grad每层都要create_graphTrue。虽然算得慢但能帮你确认自动微分链路是完整的。3.2 偏微分方程Poisson方程四种边界条件的损失函数差异资源包里有一组Poisson方程算例非常值得逐个跑Dirichlet、Neumann、Robin、Periodic。这四种边界条件写进损失函数的方式完全不同也是PINN最容易翻车的地方。以二维Poisson方程 Δu f 为例先给一个Dirichlet和Neumann混合的损失函数框架def poisson_loss(model, x_interior, x_dirichlet, x_neumann): u model(x_interior) # 手动算拉普拉斯先求一阶梯度的每个分量再叠加 u_x torch.autograd.grad(u, x_interior, grad_outputstorch.ones_like(u), create_graphTrue)[0] u_xx torch.autograd.grad(u_x[:, 0], x_interior, grad_outputstorch.ones_like(u_x[:, 0]), create_graphTrue)[0][:, 0] u_yy torch.autograd.grad(u_x[:, 1], x_interior, grad_outputstorch.ones_like(u_x[:, 1]), create_graphTrue)[0][:, 1] lap_u u_xx u_yy f compute_source(x_interior) # 由PDE.py提供 loss_pde torch.mean((lap_u - f)**2) # Dirichlet约束函数值 u_d model(x_dirichlet) loss_dirichlet torch.mean((u_d - g_dirichlet(x_dirichlet))**2) # Neumann约束法向导数 u_n model(x_neumann) u_nx torch.autograd.grad(u_n, x_neumann, grad_outputstorch.ones_like(u_n), create_graphTrue)[0] loss_neumann torch.mean((u_nx - g_neumann(x_neumann))**2) return loss_pde loss_dirichlet loss_neumann这里有几个关键点。Neumann边界条件的本质是约束解在边界上的法向导数所以损失里出现的是u_nx而不是u_n。Robin条件更特殊它是 u α·u_n g 的组合形式损失函数里要同时包含函数值和导数值的线性组合。Periodic条件则要额外处理它要求 u 在两端以及各阶导数都相等通常做法是在损失里加 u(x_left) - u(x_right) 和 u_x(x_left) - u_x(x_right) 的均方误差。资源包里的robin plot.png是Robin算例的边界误差可视化图跑完之后可以对照着看边界点上的误差分布。Δu 2, x∈[-1,1], u(-1)0, u(1)4.ipynb这个算例在一维情形下就是 u2手推解析解 u x²2x1。跑通之后可以拿预测曲线和这个解析解直接对比。我实测中发现Neumann边界那一端总是比Dirichlet端误差大后来把边界点的权重单独调大了10倍情况才好转。这是PINN的常见现象Dirichlet边界约束的是函数值本身网络容易学Neumann约束的是导数值信息要经由梯度传回去等效权重会低需要人为放大。3.3 特征值问题与Laplace圆域几何域采样与极坐标转换Laplace_equation_on_a_disk.ipynb 是另一个容易让人困惑的算例因为求解域是圆盘而不是矩形。直接在直角坐标里对圆盘采样边界点很难均匀分布在圆周上而且靠近边界的采样点会被网格化的坐标遗漏。常见做法是转极坐标x r·cos(θ), y r·sin(θ)在 r∈[0,R], θ∈[0,2π] 上采样再映射回直角坐标作为网络输入。R 1.0 r torch.rand(500, 1) * R # 半径 [0, R] theta torch.rand(500, 1) * 2 * torch.pi x r * torch.cos(theta) y r * torch.sin(theta) x_interior torch.cat([x, y], dim1) # 网络输入是 (x, y)圆域Laplace方程的难点还在于产生网络输入的方法是极坐标采样但损失函数里算梯度用的仍然是直角坐标(x,y)所以映射关系要保持一致。采样点如果过于不均匀比如部分区域点密集、部分区域稀疏会直接导致残差损失在不同区域权重失衡。这时候可以固定随机种子并且在每次迭代中重新采样让网络在训练过程中见到更多不同的点。资源包里data_domain.ipynb就是做数据域可视化与采样检查用的跑完采样后把散点图画出来看一眼是判断几何域对不对的最快方法。4. 方程组与封装框架PDAE、Lorenz系统、DeepXDE与Jacobian-Hessian4.1 方程组与系统A_simple_ODE_system和Lorenz系统的损失叠加单方程跑通之后下一个台阶是方程组。A_simple_ODE_system.ipynb和PDAE_system.ipynb讲的就是多分量输出怎么组织损失。以Lorenz系统为例三个方程耦合在一起网络输出不再是单个u而是三个分量 u,v,wdef lorenz_residual(model, t): t t.requires_grad_(True) # 确保时间维度可求导 out model(t) # out[:, 0]x, out[:, 1]y, out[:, 2]z x, y, z out[:, 0:1], out[:, 1:2], out[:, 2:3] x_t torch.autograd.grad(x, t, grad_outputstorch.ones_like(x), create_graphTrue)[0] y_t torch.autograd.grad(y, t, grad_outputstorch.ones_like(y), create_graphTrue)[0] z_t torch.autograd.grad(z, t, grad_outputstorch.ones_like(z), create_graphTrue)[0] sigma, rho, beta 10.0, 28.0, 8.0 / 3.0 # 标准Lorenz方程组残差 r1 x_t - sigma * (y - x) r2 y_t - x * (rho - z) y r3 z_t - x * y beta * z return r1, r2, r3三个残差各自算MSE再相加就是总损失。训练这种系统比起单方程更要注意损失项之间的量级平衡。Lorenz系统的典型参数 ρ28 会让 y 的量级远大于 x、z如果三个残差不做任何归一化网络会把大部分学习能力都花在量级最大的那个分量上另外两个分量收敛很慢。我一般会先跑一小段看一下三个残差各自的初始loss然后给每个残差除以它的初始值让它们从同一个起点出发。时间序列类的方程比如扩散方程和Lorenz系统还有一个容易被忽略的细节初始条件与时间域采样的比例。如果你把整个时间域均匀采样但只在t0处放一批初始条件点这些点在整个训练集里占比太小约束会被残差损失淹没。常见做法是把初始条件点单独作为一个batch每轮训练都从初始点集合里重新采样一批而不是一次性生成后用到结束。4.2 DeepXDE封装框架几行代码定义方程与边界DeepXDE Poisson_equation_Dirichlet.ipynb属于资源包里的锦上添花版本它用DeepXDE库把PINN的流程进一步封装。手写PyTorch的好处是一切可控但代码量大DeepXDE则把方程定义、边界条件、采样、模型训练都收敛成高层API。同一个Poisson方程在DeepXDE里大概长这样import deepxde as dde # 几何域二维矩形 geom dde.geometry.Rectangle([-1.0, -1.0], [1.0, 1.0]) # 方程Δu f def pde(x, y): dy_xx dde.grad.hessian(y, x, component0, i0, j0) dy_yy dde.grad.hessian(y, x, component0, i1, j1) return dy_xx dy_yy - f(x) # 边界条件四周Dirichlet值为0 bc dde.icbc.DirichletBC(geom, lambda x: 0.0, lambda _, on_boundary: on_boundary) data dde.data.PDE(geom, pde, bc, num_domain1000, num_boundary100) net dde.nn.FNN([2] [50]*3 [1], tanh, Glorot normal) model dde.Model(data, net) model.compile(adam, lr1e-3) model.train(iterations3000)注意DeepXDE库现在版本迭代比较快老版本里dde.icbc.DirichletBC曾被写作dde.DirichletBC如果导入报错先检查版本。另外DeepXDE的dde.grad.hessian返回的已经是二阶偏导不用再手写嵌套grad。这种封装方式的优点是写起来快适合作对比实验缺点是出了问题不好debug如果网络输出了NaN你很难判断是方程定义错了还是采样问题最终还是得回到手写版去定位。我的习惯是手写版作为基准DeepXDE版作为交叉验证工具两边结果一致了才放心换问题。4.3 Jacobian-Hessian方法二阶导数计算的两种数值路线Method Testing Jacobian-Hessian methods for Laplace equation.ipynb和Jacobian-Hessian methods for ODE system.ipynb这两个文件名字里的Jacobian和Hessian值得展开讲。PINN计算一阶导本质是求输出u对输入x的梯度对一个多输入多输出的网络来说完整梯度是一个Jacobian矩阵计算二阶导得到的是Hessian矩阵。资源包把这两种方法单独拎出来做测试核心是让你比较两种计算路线的精度和开销。第一种路线是用torch.autograd.grad嵌套这也是前面代码里一直用的方式先求一阶梯度保留计算图再对梯度分量重复求导。精度高但两次反向传播会累积计算开销和内存。第二种路线是解析地组织Jacobian和Hessian的乘积计算。对于Laplace方程 Δu u_xx u_yy如果你先把网络输出u对x、y的一阶导分别算出来再按坐标组合出二阶导PyTorch内部会复用一部分中间结果。在输入维度低、网络宽度大时第二种方式的内存占用会更友好。对于二维Laplace问题两种路线的最终loss差异通常不大但在训练中段会出现微小的数值分歧这是因为浮点运算顺序不同会引入舍入误差。我个人的建议是先跑通嵌套grad版本因为它的逻辑和数学表达式一一对应不容易出错等需要压内存或者算高阶动力学方程时再换成Jacobian-Hessian的显式版本并用小算例对比两者的loss曲线确认一致性。5. PINN避坑指南损失不降、边界不满足、高阶导数算错5.1 损失降到1e-4但解不对权重分配与残差点采样现象训练完总loss已经压到1e-4把预测曲线画出来和解析解完全对不上形状明显偏离。原因损失函数里各项没有加权。边界条件只有两个点约束极强loss很容易降到很低但PDE残差有几百个点如果残差项权重相对较小网络会牺牲方程约束来满足边界总loss看似收敛实际上方程根本没学会。解决给各项设置可调权重比如 loss λ_pde * loss_pde λ_bc * loss_bc。我一般从λ_pde1.0、λ_bc10.0起步然后观察单项loss如果loss_pde长期比loss_bc高一两个量级就把λ_pde加大或者把残差采样点加密。权重调节的目标不是总loss最小而是每项都降到同一量级。判断标准很简单把loss_pde和loss_bc分别打印出来两者在同一数量级内解才是同时满足物理和边界的平衡点。5.2 边界条件总是被破坏hard constraint与soft constraint的取舍现象边界点上的误差始终比其他区域大一个量级把边界点数量从10加到500也没改善多少。原因soft constraint是用损失项软性约束边界值网络在边界附近拟合的是一个平滑过渡不是精确强制。边界点再多梯度传播也会被PDE残差项稀释等效边界权重始终偏低。解决改用hard constraint把解构造成自动满足边界条件的形式。一维问题最常用的是 u(x) (x-a)(x-b)·N(x) 线性插值项构造后边界条件天然成立损失函数里直接删掉边界项。但要提醒一句hard constraint会改变网络需要学习的函数形状如果N(x)在边界附近出现数值突增反而需要更宽的网络才能稳住。最稳妥的做法是两种约束都试一遍比较验证集误差再选。资源包里-yπ²sin(πx)那个算例就是hard constraint的典型示例可以拿它做参照。5.3 训练到一半loss突然上涨学习率与优化器切换现象Adam跑了1000步loss稳步下降切到LBFGS后第一个epoch loss直接涨到初始值附近再也没下来。原因Adam和LBFGS对参数空间的理解不同。Adam会把参数推到它认为的平坦谷底但LBFGS在迭代时对曲率更敏感如果切换时参数还处于狭长谷底LBFGS一步就跨出了有效区域loss瞬间弹回高位。解决切换前先把Adam的学习率降到1e-4再多跑200步做平滑过渡切LBFGS时把max_iter从默认20调小到5观察第一个closure的loss变化确认下降趋势再逐步加大。另一个思路是干脆换回Adam把总训练步数直接加长虽然慢一点但稳定。这个资源包里凡是跑到Lorenz系统的notebook都能复现这个现象所以训练循环里最好加一个切换后loss高于切换前就回滚的判断逻辑。5.4 高阶导数计算缓慢或报错create_graph与grad设置现象算二阶导时报 requires_grad 相关错误或者一阶导能算但二阶导的loss曲线抖动剧烈。原因第一次torch.autograd.grad忘了设create_graphTrue导致后续导数无法继续求或者对分量的梯度计算方式写错比如二维问题里对u_x张量整体求导交叉梯度的分量会互相干扰。解决把所有中间梯度都加上create_graphTrue设计时就假设后面还要继续求导不要等到写二阶导代码再回头补。对二维问题的二阶导建议先分别取出u_x的第0列和第1列再对它们各自求一次grad不要直接对整个u_x张量求导。Euler Beam那个四阶算例连续四层grad嵌套每一步都要检查create_graph有没有被不小心关掉这是资源包里最容易触发这个报错的地方。5.5 不同算例结果差异大初始化与网络宽度的影响现象同一个方程换一台机器或重跑一次结果忽好忽坏没有一个稳定的复现结果。原因PyTorch默认的权重初始化是随机的PINN的训练又是非凸优化不同初始化可能落到不同的局部最优。另外网络宽度如果太小表达能力不足也会造成结果抖动。解决固定随机种子是第一步torch.manual_seed(42) np.random.seed(42)但这只是保证可复现不代表结果最优。更实际的做法是跑3到5个不同种子选验证误差最小的模型同时把网络宽度从50加到100观察loss下降曲线的稳定性变化。如果3个种子的结果差异很大基本可以判断网络结构或损失权重需要调整而不是运气问题。注意固定随机种子只保证可复现不保证最优。PINN的调参顺序是先确认损失函数各项结构正确再平衡权重最后才动网络宽度和种子。6. 验证三板斧解析解对比、残差检查与收敛曲线判断6.1 用解析解做数值对比两个可直接套用的指标对于资源包里那些有解析解的算例验证最简单直接。对ODE算例把预测曲线的均方误差和最大误差都算出来对Poisson算例除了整体误差还要单独看边界上的最大偏差。u_pred model(x_test).detach().numpy() u_true exact_solution(x_test).detach().numpy() mse np.mean((u_pred - u_true)**2) max_err np.max(np.abs(u_pred - u_true))严格来说PINN训练得到的模型在测试点上的误差同时包含网络逼近误差和优化误差。如果mse停留在1e-4左右下不去先别急着怀疑代码回去检查残差点的数量和分布是不是覆盖了全域其次是边界权重是否过强压住了内部方程。ysin(πx)·cos(πx)和Δu2这两个算例都有解析解是最容易做对比的两个。6.2 在测试点重新计算物理残差有些算例没有解析解这时候不能只看loss要回到物理本身把测试点代入网络预测再重新算一遍方程残差观察残差的空间分布。二维问题可以把残差的绝对值做成热力图看它是不是集中在某个角落一维问题直接画残差曲线。残差存在明显的区域性峰值说明那一带的采样点密度不够或者边界条件的权重压得太重。这个检查和训练过程中的loss_pde不是一回事——训练loss是优化目标测试残差是模型在未见过点上的真实物理偏差。6.3 收敛曲线的三个判据与参数调整方向最后看训练过程的loss曲线我一般用三个判据第一曲线是否单调下降如果出现周期性抖动说明学习率偏大第二总loss和各分量loss之间是否在同一量级如果某项长期比其他项大十倍以上权重就该调第三切LBFGS后是否还能继续下降如果切了之后纹丝不动说明Adam阶段已经收敛得差不多了再跑Adam意义也不大。这个资源包我拆完后最大的收获不是代码本身而是明白了一个道理PINN的调试顺序永远是先检查损失函数各项的构成对不对再谈优化器、网络宽度这些偏方。现在我做任何新算例都强制自己走一遍解析解对比、物理残差分布、各分量loss量级这三板斧跑通之后再往下做。资源包里的十几个notebook可以按ODE到PDE、单方程到方程组、手写到DeepXDE的顺序逐个复现每一步都有对照结果非常适合当PINN的入门调试基准。希望帮到你。本文还有配套的精品资源点击获取
返回列表