前段时间一直在折腾物理信息神经网络PINN在传热问题里的实际落地。手头有个场景是两块平行平板之间的二维稳态对流传热要预测温度场。传统做法是画网格跑CFD但临时搭个求解器实在费劲于是我从零用Python和PyTorch写了一个软PINNsoft PINN来解这个二维稳态对流传热方程。所谓软PINN简单说就是把控制方程残差和边界条件当作软约束加权塞进损失函数里让神经网络在拟合数据的同时也遵守物理规律。整个过程踩了不少坑损失项互相打架、自动微分得到nan、边界角点残差降不下去……这篇文章把我最终跑通的实现思路、完整代码和调参经验全部分享出来适合那些不想被有限元网格折磨、想快速拿一个温度场参考解的工程师也适合刚接触PINN的研究者。文中的代码基本在CPU上就能跑通不依赖高级显卡重点在于理解软PINN训练的逻辑。1. 为什么是软PINN平板对流传热问题里的约束哲学1.1 二维对流传热模拟中硬约束的局限性传统PINN大致分成两种角色硬约束和软约束。硬PINN会通过修改网络输出结构把边界条件直接嵌入到预测值中比如用一个距离函数乘以网络输出再加上边界温度从而保证边界条件自动满足。听起来很完美但在平板间对流传热这样的问题里边界不会是单纯的一条直线入口、出口、上下壁面有Dirichlet边界也有Neumann边界。特别是出口的∂T/∂x0要构造距离函数并保证法向导数恒成立非常麻烦。更致命的是多层网络的梯度在乘以距离函数后会被严重削弱导致内部区域更新缓慢整体收敛变差。我在早期尝试硬PINN时对这一点体会极深。1.2 软约束的本质把物理方程变成一项损失软PINN的思路则更简单粗暴把控制方程残差r(x,y)变成损失项把边界残差也变成损失项训练过程就是让这两组残差同时下降。损失函数写成L L_pde λ_bc * L_bc其中L_pde是计算域内部配点上控制方程残差的均方误差L_bc是边界条件残差。你可以随时在损失里添加数据项、界面连续性项而不用去改网络结构。对于本文这种边界不算极端复杂、但又要处理混合边界条件的场景软PINN是更省事的选择。当然省事的代价是损失中各项权重需要调否则很容易出现“PDE学得很好但边界一塌糊涂”的尴尬局面。后面我会单独用一个小节来聊权重问题。1.3 什么情况下我会建议你优先用软PINN根据我的实践以下情况优先考虑软PINN。第一边界条件比较复杂混合了Dirichlet、Neumann甚至Robin条件第二你需要在现有模型上快速扩展参数反演比如反推热流密度或换热系数第三作为教学或验证手段想快速看到PINN的训练过程。反过来如果计算域非常简单比如矩形域且边界条件全是Dirichlet那硬PINN会更稳毕竟它从结构上保证边界精确满足。本文最终用软PINN还有一个教学上的理由软约束的训练过程更容易观察损失分量变化帮助你理解PINN为什么能收敛又为什么可能不收敛。2. 二维稳态对流传热方程控制方程与无量纲化2.1 从能量守恒到对流-扩散方程在二维稳态不可压缩流体中温度场T(x,y)满足能量守恒方程。忽略粘性耗散、体积力和辐射温度由对流与扩散共同决定ρ c_p (u ∂T/∂x v ∂T/∂y) k (∂²T/∂x² ∂²T/∂y²)其中ρ是流体密度c_p是比热容k是导热系数u和v是速度场分量。引入热扩散率α k/(ρ c_p)方程简化为u ∂T/∂x v ∂T/∂y α (∂²T/∂x² ∂²T/∂y²)如果速度场本身也需要求解那就要耦合动量方程。本文直接给定流速分布只关心温度场这是典型的“已知流场求温度”的简化。工程上很多换热器初步估算都走这条路。2.2 无量纲化PINN训练的第一步也是最容易被跳过的一步很多初学者直接把物理量塞进网络结果训练老不稳定。原因很简单长度用米温度用摄氏度速度用米每秒这些数值量级从10⁻⁵到10²不等神经网络的权重初始化和梯度更新很难同时适配这么宽的取值范围。解决办法就是无量纲化。取参考长度L_ref、参考速度U_ref、参考温差ΔT定义x x / L_ref, y y / L_ref, θ (T - T_cold) / (T_hot - T_cold), u u / U_ref, v v / U_ref代入方程后得到无量纲形式u ∂θ/∂x v ∂θ/∂y (1/Pe) (∂²θ/∂x² ∂²θ/∂y²)其中Péclet数 Pe U_ref * L_ref / α它表征对流输运与扩散输运的相对大小。Pe越大对流越强方程越接近一阶双曲型对PINN来说也越难训练这就是为什么很多对流传热PINN案例喜欢从低Pe开始。本文算例取L_ref0.1mU_ref0.01m/sα2.0×10⁻⁵ m²/s则Pe 0.01×0.1 / 2e-5 50属于对流主导但还没到极度刚硬的程度适合用来展示软PINN的能力边界。2.3 边界条件的无量纲化与软约束展开算例几何是两块平行平板之间的矩形区域入口在左侧出口在右侧上下为恒温壁面。无量纲边界条件为入口 x0θ1Dirichlet出口 x1∂θ/∂x 0Neumann上下壁面 y0 或 y0.5θ0Dirichlet这里有个细节上下壁面的无量纲y范围取决于几何。假设平板间距0.05m则y_max 0.05/0.1 0.5。所有坐标都落到[0,1]×[0,0.5]这个矩形里。如果不做归一化网络输入x范围是[0,1]y范围是[0,0.5]虽然差异不算太大但为了后续同时处理多种几何统一归一化到[0,1]会更好。我在代码里会把y坐标除以最大y让域变成单位正方形两条边的长度尺度就一致了。这一点看似简单实际对训练稳定性影响很大。3. 软PINN的PyTorch网络设计与损失函数构造3.1 网络结构一个带Tanh隐藏层的MLP选择网络结构时不要一上来就堆大模型。针对二维二阶PDE一个3到5层的全连接网络足够。激活函数我第一个选的是Tanh而不是ReLU。原因是ReLU的一阶导数是阶跃函数二阶导数处处为零而损失函数里需要求温度的拉普拉斯算子∂²θ/∂x²如果激活函数是ReLU网络的二阶导数会被直接截断成零PDE残差项的信息全部丢失训练会变成“只有对流项在起作用”。用Tanh或Swish能保证一阶导、二阶导都有连续的非零梯度。以下是我用的网络结构import torch import torch.nn as nn class SoftPINN(nn.Module): def __init__(self, in_dim2, hidden_dim64, out_dim1): super().__init__() self.net nn.Sequential( nn.Linear(in_dim, hidden_dim), nn.Tanh(), nn.Linear(hidden_dim, hidden_dim), nn.Tanh(), nn.Linear(hidden_dim, hidden_dim), nn.Tanh(), nn.Linear(hidden_dim, hidden_dim), nn.Tanh(), nn.Linear(hidden_dim, out_dim) ) # 用Xavier初始化减轻梯度消失 for layer in self.net: if isinstance(layer, nn.Linear): nn.init.xavier_uniform_(layer.weight) nn.init.zeros_(layer.bias) def forward(self, x, y): coords torch.cat([x, y], dim-1) return self.net(coords)这里每层64个神经元、4个隐藏层是为了平衡二阶导计算带来的算力开销。你可以试试更宽的层但训练时间会显著上升。初始化用Xavier是因为它适配Tanh这类sigmoid型激活函数让激活值在训练初期保持合理范围。3.2 自动微分用torch.autograd.grad求偏导数PINN之所以强大是因为神经网络的导数可以用自动微分精确求出而不是像有限差分那样有截断误差。在PyTorch里我们对网络输出θ分别求∂θ/∂x、∂θ/∂y再求二阶导需要用到create_graphTruedef compute_pde_residual(model, x, y, u, v, alpha): x,y: 内部配点形状 [N,1] u,v: 速度场在对应点的值这里假设恒定速度 alpha: 热扩散系数 x.requires_grad_(True) y.requires_grad_(True) theta model(x, y) # 输出温度 theta_x torch.autograd.grad(theta, x, grad_outputstorch.ones_like(theta), create_graphTrue, retain_graphTrue)[0] theta_y torch.autograd.grad(theta, y, grad_outputstorch.ones_like(theta), create_graphTrue, retain_graphTrue)[0] theta_xx torch.autograd.grad(theta_x, x, grad_outputstorch.ones_like(theta_x), create_graphTrue, retain_graphTrue)[0] theta_yy torch.autograd.grad(theta_y, y, grad_outputstorch.ones_like(theta_y), create_graphTrue, retain_graphTrue)[0] # 无量纲方程u * T_x v * T_y (1/Pe) (T_xx T_yy) residual u * theta_x v * theta_y - alpha * (theta_xx theta_yy) return residualretain_graphTrue必须加上因为在同一张计算图里连续做了四次grad如果不保留中间图的梯度第二次grad可能会报错。这是初学者最容易踩的坑。3.3 损失函数三路并进权重要有策略软PINN的损失函数由多个部分组成PDE残差损失L_pde mean(‖r(x_i,y_i)‖²)边界损失Dirichlet边界用L_bc_dir mean(‖θ - g‖²)Neumann边界用mean(‖∂θ/∂n - g_n‖²)可选的数据损失L_data mean(‖θ - θ_obs‖²)组合时建议先不要一开始就用自适应权重而是固定简单权重比如loss pde_loss 50.0 * bc_loss这里边界权重取50原因是在平板传热问题中边界点数量通常远小于内部配点且边界条件需要被相对严格地满足。实际训练中你可以每隔500步打印一次各分量的量级再酌情调整。后续的“调参坑”里我会展开一个更系统的统计方法。4. 从零实现软PINN求解平板间温度场的完整代码4.1 环境准备与复现性设置我的环境是Python 3.10 PyTorch 2.1 CPU训练完全能跑通这个规模的问题。如果你用GPU训练会快不少但二阶自动微分带来的显存开销也不小小显存的显卡未必比CPU舒服。建议先用CPU验证逻辑再上GPU。import torch import torch.nn as nn import numpy as np import matplotlib.pyplot as plt torch.manual_seed(42) np.random.seed(42)4.2 数据生成内部配点与边界采样训练软PINN不需要“训练集标签”只需要坐标点。域内随机采样产生内部点边界用等间距或随机采样。这里的关键是角点处理角点同时属于两条边界如果你把角点同时算入两边的边界损失会要求网络同时满足两个不同的值比如入口θ1和上壁θ0这会造成不可调和的冲突。我的做法是把角点从边界集合里剔除让网络在角点附近用PDE残差自然过渡。Nf 10000 # 内部配点数量 Nb 200 # 每条边采样点数量 x_min, x_max 0.0, 1.0 y_min, y_max 0.0, 0.5 # 内部点 x_f torch.rand(Nf, 1) * (x_max - x_min) x_min y_f torch.rand(Nf, 1) * (y_max - y_min) y_min # 边界点分别采样入口、出口、上壁、下壁 x_in torch.zeros(Nb, 1) y_in torch.rand(Nb, 1) * (y_max - y_min) y_min x_out torch.ones(Nb, 1) y_out torch.rand(Nb, 1) * (y_max - y_min) y_min x_top torch.rand(Nb, 1) * (x_max - x_min) x_min y_top torch.ones(Nb, 1) * y_max x_bottom torch.rand(Nb, 1) * (x_max - x_min) x_min y_bottom torch.zeros(Nb, 1) x_bc torch.cat([x_in, x_out, x_top, x_bottom], dim0) y_bc torch.cat([y_in, y_out, y_top, y_bottom], dim0)因为每条边是分开采样的角点不会同时出现在同一边界集合里所以边界点不会被入口和上壁的约束同时覆盖这个冲突在数据层面就避开了。如果你把四条边直接拼接成一个集合则必须做去重。一个更省事的做法是先构造一个只包含20个配点的粗算例快速验证损失函数能降确认无误后再放大到1万个配点。这能帮你省下大量调试时间。4.3 训练主循环Adam热启动L-BFGS精调我最终跑的循环是“Adam先走到一个比较平滑的loss平台再切到L-BFGS精调”。Adam对不平衡的损失项容忍度高适合前期把残差从O(1)压到O(0.01)L-BFGS作为拟牛顿法在接近最优解时能更精确地减少PDE残差但对损失项量级更敏感。完整训练循环如下model SoftPINN() optimizer_adam torch.optim.Adam(model.parameters(), lr1e-3) # 速度场恒定u0.3, v0.0无量纲化后 u 0.3 v 0.0 alpha 1.0 / 50.0 pde_weight 1.0 bc_weight 50.0 def closure(): optimizer_adam.zero_grad() # PDE residual residual compute_pde_residual(model, x_f, y_f, u, v, alpha) pde_loss torch.mean(residual ** 2) # 边界损失 theta_in model(x_in, y_in) bc_in torch.mean((theta_in - 1.0) ** 2) theta_out model(x_out, y_out) # 出口 Neumann: dT/dx0需要自动微分 theta_out_x torch.autograd.grad(theta_out, x_out, grad_outputstorch.ones_like(theta_out), create_graphTrue)[0] bc_out torch.mean(theta_out_x ** 2) theta_top model(x_top, y_top) theta_bottom model(x_bottom, y_bottom) bc_wall torch.mean(theta_top ** 2) torch.mean(theta_bottom ** 2) bc_loss bc_in bc_out bc_wall loss pde_weight * pde_loss bc_weight * bc_loss loss.backward() return loss, pde_loss, bc_loss # Adam 阶段 for step in range(10000): loss, pde_loss, bc_loss closure() optimizer_adam.step() if step % 1000 0: print(fAdam step {step}: loss{loss.item():.4e}, pde{pde_loss.item():.4e}, bc{bc_loss.item():.4e}) # L-BFGS 精调阶段 optimizer_lbfgs torch.optim.LBFGS(model.parameters(), lr0.1, max_iter500) def lbfgs_closure(): optimizer_lbfgs.zero_grad() loss, _, _ closure() loss.backward() return loss optimizer_lbfgs.step(lbfgs_closure)需要注意两点一是Neumann边界上的自动微分必须带上create_graphTrue否则无法在二阶损失里继续回传二是如果Adam训练到后期边界损失仍然降得慢不要急着切L-BFGS先回头调边界权重。L-BFGS对平坦区域改善明显但无法解决权重失衡。4.4 预测与可视化训练完以后在网格上生成预测温度场并和边界条件一起画出来。可视化既是为了验收也能帮你快速发现训练是否在某个角落“跑飞”了。nx, ny 100, 50 x_plot torch.linspace(0, 1, nx) y_plot torch.linspace(0, 0.5, ny) X, Y torch.meshgrid(x_plot, y_plot, indexingij) X_flat X.reshape(-1, 1) Y_flat Y.reshape(-1, 1) with torch.no_grad(): T_flat model(X_flat, Y_flat) T T_flat.reshape(nx, ny).numpy() plt.figure(figsize(8, 4)) plt.contourf(X.numpy(), Y.numpy(), T, levels50, cmaphot) plt.colorbar(labeltheta) plt.xlabel(x) plt.ylabel(y) plt.title(Soft PINN predicted temperature field) plt.savefig(pinn_temperature_field.png, dpi150)如果这个温度场看起来平滑、入口扩散均匀、出口近似等温线垂直说明训练基本成功。如果云图出现波纹或局部凸起多半是配点太少或边界权重不足。5. 我踩过的软PINN调参坑逐条排查全记录5.1 所有loss变成nan先查输入范围再查学习率第一次跑训练我遇到的是loss变成nan。第一反应是梯度爆炸但我把学习率降到1e-5还是nan。后来查了输入数据发现y坐标用的是物理尺寸0.05而x坐标用归一化后的1两者差20倍。网络一开始面对这种尺度失衡梯度方向被y方向主导权重更新很快就溢出。这个问题在把y也归一化到0~1后就消失了。所以第一优先级永远是把输入坐标归一化而不是去调学习率。5.2 PDE loss降得到底但边界凸不起来权重比例失衡典型现象是训练日志里pde_loss已经掉到1e-5但bc_loss还在0.1量级。我看输出温度场时发现入口边界附近温度居然从1变成了0.8明显是边界条件没学好。为什么会这样因为内部配点有10000个边界点一共只有800个如果统一权重PDE残差对总loss的贡献天然比边界的贡献大网络自然优先满足内部方程。解决办法不是直接疯狂加大边界权重而是先统计各项的量级。我在每次迭代时打印print(fpde_loss{pde_loss.item():.3e}, bc_loss{bc_loss.item():.3e})然后把边界权重调到能让两项在训练后期处于同一数量级。我的经验是在本文算例里bc_weight设成50~100比较合适。边界权重过大会导致边界点学习了但内部区域被边缘点的强约束“拉出”不符合PDE的解所以也不要一味地调大。5.3 出口Neumann边界不消失梯度成了“摆设”出口∂θ/∂x0这个边界条件我一开始直接用网络对x的导数当损失但发现训练结束后出口附近仍然有明显梯度。原因在于出口边界上的点只参与边界损失但出口附近内部区域的配点没有显式地要求温度“平直”网络在高阶导上找到的解依然是光滑的并不会天然满足∂θ/∂x0。后来我在出口区域增加了“缓冲区”配点即离出口很近的内部点也额外加权PDE残差相当于把Neumann条件沿上游方向传播了一点。这个方法很土但有效。你也可以在边界损失里用一阶导匹配代替零梯度或者提高出口边界点的数量。5.4 角点处残差永远降不下去接受它别较劲这是软PINN的固有问题。入口和上壁的交点处一边要求θ1另一边要求θ0无论怎么调权重网络只能折中成0.5左右。这是数学上不可调和的不连续点。我一开始拼命加角点附近的配点结果只是让角点产生尖锐但陡峭的过渡整体误差反而变大。后来我接受这个事实训练时把角点从边界集合剔除让角点附近的解由内点PDE残差和附近边界点共同决定整体误差反而更小。如果你想严格满足角点连续性需要先对问题做角点奇异性处理那就不是软PINN一个模型能解决的了。5.5 训练时间过长合理裁剪算力二阶自动微分是算力大头。同样10000个配点一次backward要计算四个一阶导数、两个二阶导数计算量比普通分类网络大好几个量级。我实测用CPU训练10000步Adam加500步L-BFGS大概花了一夜。如果想快速验证我建议把配点降到3000隐藏层宽度降到32Adam训5000步十几分钟就能看到整体趋势。等趋势对了再跑正式版本。这个“先用小模型验证、再上大模型”的习惯能帮你节省大量等待时间。5.6 用解析解或高精度数值解做交叉验证软PINN调参是不是成功了不能只看loss曲线。我会额外构造一个简单场景令Pe很小例如Pe0.1方程退化为纯扩散问题此时平板间稳态温度分布有解析解。把软PINN预测的等温线和解析解叠在一起计算整体相对L2误差如果小于5%说明训练和代码实现基本正确。这个交叉验证步骤非常重要它能把“loss在下降”和“物理对了”两件事分开。为了更快定位问题我把常见的现象和对策整理成了一张速查表现象可能原因处理办法lossnan输入坐标未归一化先归一化输入再降学习率pde_loss低但bc_loss高边界权重不足调整λ_bc使两项量级接近出口法向梯度不为0出口缓冲区点太少增加出口附近配点权重角点残差大边界条件冲突从边界集合中剔除角点训练极慢二阶导计算过大减小网络宽度和配点数6. 下一步扩展软PINN还能怎么玩出新东西6.1 往瞬态传热扩展把时间当成第三个坐标平板传热如果做非稳态只需把网络输入从(x,y)变成(x,y,t)控制方程加上∂θ/∂t项损失函数里多一个初始条件约束。采样时内部配点从二维平面变成三维时空区域[0,1]×[0,0.5]×[0,t_max]其他逻辑完全不变。我试过类似问题瞬态任务的难度主要在初始条件与边界条件在t0时刻的兼容。如果初始温度跳变剧烈也会出现类似角点残差的问题处理思路同样是先让初始条件软约束主导再让PDE项逐步接管。6.2 反问题根据稀疏温度测量反推边界热流软PINN的一个天然优势是容易处理反问题。假如你在实验室测了少量点的温度想把入口温度或壁面热流反推出来那么只需要把未知参数例如入口温度T_in定义成nn.Parameter在损失函数里加一项数据拟合损失T_in nn.Parameter(torch.tensor(1.0))然后训练时同时更新网络权重和T_in反演就完成了。这个思路能扩展到反推换热系数、热源强度等多个参数是很多工程检测场景的刚需。我之前在一个小项目里用反向PINN同时反推两个边界热流效果比传统的优化迭代好很多几乎不用手动设置初值。6.3 区域分解和自适应配点进阶的“大力出奇迹”当计算域变大、对流变强单一网络会很吃力。可以考虑把计算域切成几个子区域每个子区域用一个子网络子网络之间通过内部界面上的温度连续性和热流连续性耦合。这有点像“神经有限元”但实现复杂度高出不少。对我来说最实用的技巧还是自适应配点训练到一半把当前PDE残差大的地方多采点样重新构成配点集合继续训练。这相当于把网格加密的精神引入PINN效果好实现也不难with torch.no_grad(): residual compute_pde_residual(model, x_f, y_f, u, v, alpha) _, indices torch.topk(residual.abs().flatten(), kNf // 5) new_points torch.stack([x_f[indices], y_f[indices]], dim-1)把这些高残差点的坐标补充进原有配点集然后按更小的学习率继续训练几轮通常能把局部尖峰压下去。最后再分享一点个人体会。软PINN这个东西表面上是“改一下损失函数”实际操作起来却有大量细节坐标归一化、边界权重、角点冲突、二阶导计算、训练策略……每一步都可能把训练带偏。我的建议是先跑通小规模算例打印各项损失分量的趋势再用解析解或高精度数值解做交叉验证。当你看到自己的软PINN在torch里一点点把温度场从噪声中“雕刻”出来时那种感觉比单纯调出一个更低loss爽得多。如果你也是刚接触物理信息神经网络希望这篇记录能让你少走几个弯路。
