前段时间我在整理结构参数快速评估方案时又一次把悬臂梁挠度问题翻了出来。这根一端固定、一端自由的小梁在均布载荷下的公式解其实十分钟就能算完但当我需要批量扫描不同载荷和刚度还希望后续能扩展到更复杂的载荷分布时解析法就不够灵活了。于是我用物理信息神经网络PINN重写了一遍一维悬臂梁挠度的求解过程让网络自己“学会”满足梁方程而不是靠人工构造经验公式。整篇分享将从控制方程推导、PINN 损失设计、Python 实现到踩坑记录逐一展开适合想入门 PINN 的读者也适合做结构计算但还没碰过深度学习的同学。1. 悬臂梁挠度问题的控制方程与解析解1.1 梁的挠度方程为什么是四阶导悬臂梁建模里最经典的是欧拉-伯努利梁理论。假设梁沿 x 方向的截面转角很小挠度用 w(x) 表示那么梁截面的弯矩 M 和挠度曲率近似成正比M EIw″其中 EI 是抗弯刚度。对这段梁取微元做力平衡弯矩沿 x 方向的变化率就是剪力 VV M′ EIw‴剪力沿 x 方向的变化率又等于分布载荷 q(x)于是 V′ q就得到 EI w″″ q。这里连续用了两次“对导数求导”所以最后出现四阶导。把物理图像翻译成直白的话就是载荷分布 q 经过材料刚度 EI 的“积分平滑”最终变成位移场 w而位置越靠里的内力平衡条件越严格约束阶数就越高。很多不熟悉梁理论的人第一次看到 w″″ 会懵其实只要记住“弯矩、剪力、分布载荷”这三个量与位移导数之间的层层对应关系就够了。符号约定上要特别注意我这里取载荷向下为正、挠度向下为正方程右侧就是正的 q。有人习惯取向上为正那方程右边就会变成 -q解析解整体反号部署到代码里很容易出现结果方向反了的问题。我建议先固定符号约定再开始推导和编码。针对一维悬臂梁典型的边界条件如下x 0固定端w(0)0w′(0)0位移和转角都被约束住。x L自由端w″(L)0w‴(L)0表示自由端弯矩和剪力为零。这四条边界条件加上四阶常微分方程正好组成一个适定的边界值问题。解析上对求解域没有特殊限制但 PINN 训练时这些边界值会作为约束条件塞到损失函数里处理得好不好直接决定训练能否收敛。1.2 无量纲化与解析解推导如果直接带着实际量纲训练神经网络EI、q、L 的数值范围可能相差好几个数量级梯度训练非常容易被大数主导。所以我在所有实验里先把方程无量纲化。设 x̂ x/Lŵ w / (qL⁴/EI)则原方程 EI w″″ q 简化为ŵ″″ 1这里的导数都是对 x̂ 求导。无量纲化之后控制方程与具体材料参数彻底解耦训练出来的网络形状是通用的想要恢复真实挠度只需要乘以 qL⁴/EI。这个操作看似简单实际上是 PINN 工程里最容易忽略的一环。很多刚接触的人把 EI100、q5、L4 直接代入损失函数结果损失值在不同物理量纲之间来回横跳训练极不稳定。这个方程可以连续积分得到封闭解。积分四次并代入四条边界条件后ŵ(x̂) x̂⁴/24 - x̂³/6 x̂²/4这个结果可以看作 PINN 的“标准答案”。自由端 x1 处ŵ(1)1/24 - 1/6 1/4 0.125也就是无量纲自由端挠度是 1/8记住这个数后面验证网络学得对不对非常方便。换成有量纲形式真实挠度为w(x) ŵ(x̂) · qL⁴/EI1.3 为什么用 PINN 而不是解析法或有限元对均布载荷这种简单情况解析法显然效率最高。但工程里常常会遇到分段载荷、局部集中力甚至随位置变化的刚度解析解要么分段积分要么根本推不出显式表达式。传统有限元需要划分网格、组装刚度矩阵每改一次载荷或边界条件都要重新处理前处理流程在快速原型阶段显得笨重。PINN 的思路完全不同用一个可微的神经网络 w_θ(x) 近似真解把方程残差和边界条件残差变成损失函数直接用梯度下降优化网络参数。它不需要网格换载荷只要改损失函数里的 q(x) 项换边界只要改边界损失代码骨架完全不动。一维悬臂梁挠度问题正好是验证这套思路的最小实验有解析解可以对比方程形式清晰训练开销也不大非常适合作为 PINN 的入门练习。另外一个潜在优势在逆问题场景如果手里只有部分测量点想反推荷载或刚度分布解析法会非常吃力而 PINN 可以把未知参数直接定义为变量和网络一起优化。一维悬臂梁看似简单却是理解这套扩展思路的最好载体。2. PINN 求解悬臂梁的核心思路损失函数、自动微分与边界约束2.1 把神经网络当成一个可微的挠度场在 PINN 里网络输入是坐标 x输出是挠度的预测值 w_net(x)。网络的前向传播本质上就是定义一个带大量参数的函数我们希望这个函数尽量满足控制方程和边界条件。为了让网络满足“高阶可导”的要求激活函数我选了 tanh。至于为什么不用 ReLUReLU 的一阶导是阶跃函数二阶导就变成零四阶导更是无从谈起对光滑性要求很高的梁问题完全不适用。tanh 在实数域上无穷阶可导特别适合配合自动微分做多次求导。求导这一步是 PINN 与普通机器学习最大的区别。手写有限差分会有截断误差网格太粗容易把导数算歪。自动微分在计算图上精确求导只要网络函数够光滑二阶导、四阶导都能直接用链式法则得到。这里还牵扯到一个网络设计选择把挠度和 θ 放在同一个网络的两个输出分支还是拆成两个独立网络我推荐前者。虽然两个物理量可以用两套独立网络分别拟合但它们背后描述的是同一根梁共享隐藏层可以让低阶特征互相复用训练时的信息传递更充分。后面代码里我会直接在一个网络里输出两个分量。2.2 损失函数到底在惩罚什么把控制方程 ŵ″″1 写成一个残差形式L_pde (1/N) Σ (ŵ″″(xᵢ) - 1)²网络预测结果如果完全满足方程这个损失项应该接近 0。边界条件也要写进去L_bc ŵ(0)² ŵ′(0)² ŵ″(1)² ŵ‴(1)²总损失通常写成L_total λ_pde · L_pde λ_bc · L_bcλ_pde 和 λ_bc 是权重。如果两个损失的量级差得太大梯度更新就会被大的那项主导小的一项虽然重要却始终“学不到位”。这个问题在 PINN 里非常常见后面第 4 章我会专门聊权重怎么调。如果手头还有实验测量数据还可以在损失函数里加一项 L_data用来拟合观测点。但纯物理已知的情况下不需要这是 PINN 区别于普通监督学习的地方没有标签数据也能训练因为物理方程本身就是标签。一维悬臂梁的解析解在这里反而是验证工具不是训练必需输入。2.3 四阶导是常见坑用降阶规避我最早写的版本是直接用自动微分求 ŵ 的直到四阶导再把 ŵ″″-1 的残差塞进损失。理论上没问题但实践里明显感觉到训练不稳。原因不难理解网络每多求一次导数反向传播路径就会叠加更多乘性梯度高阶导数值的方差被放大小扰动在梯度更新时会被成倍放大网络参数稍微一动损失就可能剧烈震荡。更稳的做法是降阶。令 θ ŵ″那么原来四阶方程被拆成两个二阶方程ŵ″ θθ″ 1对应的边界条件也改写固定端ŵ(0)0ŵ′(0)0自由端θ(1)ŵ″(1)0θ′(1)ŵ‴(1)0这样网络要输出两个物理量一个是挠度 ŵ另一个是 θ正比于弯矩。自动微分只需要求到二阶导数值稳定性明显好了一大截。我在实验里对比过直接四阶方案在同样网络结构下 loss 很难稳定下降降阶方案则可以很快压到 10⁻⁵ 量级以下。能降阶就降阶是求解高阶 PINN 问题的一个基本原则。这个思路跟有限元里的“降阶积分”有一点神似本质都是在降低高阶导数带来的强梯度变化。对于工程人员来说θ 本身也有物理意义它正比于弯矩所以降阶不是单纯的数值技巧还顺便给出了弯矩分布。2.4 硬约束边界条件的构造技巧边界条件可以用软约束损失项也可以硬约束。软约束的问题在于权重系数怎么选都可能有偏差边界条件永远只是“尽量满足”而不是“必然满足”。硬约束的思路更优雅把网络的输出经过一个显式构造的数学映射让边界条件在数学上自动成立从而直接把边界损失从总损失里删掉。对固定端构造 ŵ(x) x² · U_net(x)。因为当 x0 时无论 U_net 输出多少x² 都能保证 ŵ(0)0同时 ŵ′(0) 2x·U x²·U′在 x0 处也自动为 0。所以固定端的两个边界条件被牢牢绑死。对自由端由于 θ ŵ″我需要 θ(1)0 和 θ′(1)0于是构造 θ(x) (1-x)² · T_net(x)。当 x1 时θ(1)0θ′(x) 包含 (1-x) 因子所以 θ′(1) 也自动等于 0。这样一来四条边界条件都不需要进损失函数训练只剩内部残差项优化目标干净许多。这种“用零点因子去强制边界”的技巧在 PINN 里非常实用。你可以把它理解成在梁的端点处“预先固定了答案”网络只需要在内部把方程配平训练负担大幅下降。后面第 3 章的代码就是基于这套硬约束写的。3. 基于 Python 的一维悬臂梁挠度实现与复现3.1 环境准备代码依赖四个基础库PyTorch、NumPy、Matplotlib。PyTorch 主要负责自动微分和网络训练NumPy 负责数据比较Matplotlib 画图验证。实际上这个一维问题非常轻量CPU 上跑完全足够不需要 GPU。安装指令如下pip install torch numpy matplotlib代码里的模型非常简单三层隐藏层、每层 64 个神经元对于一维问题已经属于“豪华配置”。如果只想快速验证概念把宽度改成 32 也可以训练时间更短精度不会差太多。3.2 完整可运行代码下面这段代码就是完整实现。为了方便直接复现我做了三件事一是无量纲化二是在网络输出上直接构造硬约束边界条件三是每次训练迭代都重新随机采样坐标点。import torch import torch.nn as nn import numpy as np import matplotlib.pyplot as plt torch.manual_seed(42) class PINN(nn.Module): def __init__(self, width64): super().__init__() self.net nn.Sequential( nn.Linear(1, width), nn.Tanh(), nn.Linear(width, width), nn.Tanh(), nn.Linear(width, width), nn.Tanh(), nn.Linear(width, 2) ) def forward(self, x): raw self.net(x) return raw[:, 0:1], raw[:, 1:2] model PINN(width64) optimizer torch.optim.Adam(model.parameters(), lr1e-3) N 128 loss_history [] for it in range(20000): x torch.rand(N, 1) x.requires_grad_(True) U, T model(x) u x ** 2 * U theta (1 - x) ** 2 * T u_d1 torch.autograd.grad( u, x, grad_outputstorch.ones_like(u), create_graphTrue)[0] u_d2 torch.autograd.grad( u_d1, x, grad_outputstorch.ones_like(u_d1), create_graphTrue)[0] theta_d1 torch.autograd.grad( theta, x, grad_outputstorch.ones_like(theta), create_graphTrue)[0] theta_d2 torch.autograd.grad( theta_d1, x, grad_outputstorch.ones_like(theta_d1), create_graphTrue)[0] loss_pde1 torch.mean((u_d2 - theta) ** 2) loss_pde2 torch.mean((theta_d2 - 1.0) ** 2) loss loss_pde1 loss_pde2 optimizer.zero_grad() loss.backward() optimizer.step() loss_history.append(loss.item()) if it % 2000 0: print(fiter {it:6d}, pde1{loss_pde1.item():.3e}, pde2{loss_pde2.item():.3e}) x_test torch.linspace(0, 1, 101).reshape(-1, 1) with torch.no_grad(): U, T model(x_test) u_pred (x_test ** 2 * U).numpy() theta_pred ((1 - x_test) ** 2 * T).numpy() x_vals x_test.numpy() u_exact x_vals ** 4 / 24 - x_vals ** 3 / 6 x_vals ** 2 / 4 theta_exact x_vals ** 2 / 2 - x_vals 0.5 err np.max(np.abs(u_pred - u_exact)) print(f最大挠度误差: {err:.3e}) plt.figure(figsize(12, 4)) plt.subplot(1, 3, 1) plt.plot(x_vals, u_exact, k--, label解析解) plt.plot(x_vals, u_pred, r-, labelPINN) plt.xlabel(x/L) plt.ylabel(w/(qL^4/EI)) plt.legend() plt.title(挠度曲线) plt.subplot(1, 3, 2) plt.plot(x_vals, theta_exact, k--, label解析解) plt.plot(x_vals, theta_pred, b-, labelPINN) plt.xlabel(x/L) plt.ylabel(theta) plt.legend() plt.title(theta 曲线正比于弯矩) plt.subplot(1, 3, 3) plt.semilogy(loss_history) plt.xlabel(iteration) plt.ylabel(loss) plt.title(训练损失曲线) plt.tight_layout() plt.show()3.3 代码里的关键细节先看网络输出。forward 返回两个分量第一个经过 x² 修饰输出挠度 u第二个经过 (1-x)² 修饰输出 θ。这里有一个容易忽略的点U 和 T 来自同一个共享网络的前四层后面的训练会同时更新共享特征让挠度残差与弯矩残差互相渗透整体收敛更快。如果拆成两个独立网络每一套特征都要从头学实验里往往需要更多迭代才能达到同等精度。再看求导部分。PyTorch 的 autograd.grad 默认针对标量求导当 u 是一个向量时需要传 grad_outputstorch.ones_like(u)。它的实际效果等价于先对 u 的所有采样点求和再对 x 求梯度。这样每个点都能得到自己的导数并且 create_graphTrue 保留二阶导数图才能继续算 u_d2 和 theta_d2。最后说采样策略。代码里 x 每次都从 [0,1] 均匀随机采样而不是固定点集。这样做的好处是每一轮训练网络都在接触不同的“约束点”泛化能力更好不容易在固定点集上过拟合。代价是损失曲线会有明显振荡看整体趋势而不是单点数值。如果你更在意损失曲线平稳可以改成 torch.linspace(0, 1, 129) 这样固定点集但注意点太少时训练容易欠拟合。训练完成后最大挠度误差通常能控制得很小。误差量级取决于随机种子和训练时间我实测用上面参数跑到 20000 步最大误差在 10⁻⁴ 量级附近。如果追求更高精度继续减小学习率或者用 L-BFGS 微调就可以。3.4 如何恢复有量纲结果代码里训练的是无量纲挠度 ŵ工程上还是要换成真实物理量。举一个具体例子假设悬臂梁长度 L2 米抗弯刚度 EI10000 N·m²均布载荷 q500 N/m。先算系数qL⁴/EI 500 × 16 / 10000 0.8 m所以自由端真实挠度 ŵ(1) × 0.8 0.125 × 0.8 0.1 m。如果 x0.5 处网络预测 ŵ(0.5) 0.052083那么真实挠度就是 0.052083 × 0.8 0.041666 m。这个恢复过程很简单但别忘了方向约定我这里取向下为正输出正值就表示梁在向下挠曲。有一种常见的错误是恢复前把符号弄反或者忘记乘系数。建议在部署代码里单独写一个函数 complete(x_physical, w_hat, L, q, EI)用来统一处理坐标缩放和挠度缩放别把乘系数逻辑散落在各处。3.5 可选优化L-BFGS 微调与固定采样点很多时候 Adam 跑到后期收敛变慢因为一阶优化器在接近最优解时更新步长不够精细。这时可以切换到 L-BFGS 做精调。L-BFGS 是一种拟牛顿法利用历史梯度信息近似曲率在中小规模问题上收敛更锐利。一个常见的做法是先用 Adam 跑几千步做“热身”然后用固定点集切到 L-BFGSx_fixed torch.linspace(0, 1, 129).reshape(-1, 1) x_fixed.requires_grad_(True) opt_bfgs torch.optim.LBFGS(model.parameters(), lr0.5, max_iter300) def closure(): opt_bfgs.zero_grad() U, T model(x_fixed) u x_fixed ** 2 * U theta (1 - x_fixed) ** 2 * T # 按上文相同的 autograd.grad 方式求 u_d2 和 theta_d2 loss torch.mean((u_d2 - theta) ** 2) torch.mean((theta_d2 - 1.0) ** 2) loss.backward() return loss opt_bfgs.step(closure)需要注意的是 L-BFGS 每一步都要传入 closure 函数重新计算损失和梯度不能像 Adam 那样直接在循环里 forward 一次。因为拟牛顿法需要梯度信息来做线搜索这个机制是必须的。Adam 热身加 L-BFGS 精调是我在一维物理场问题里最推荐的标准组合。4. 常见问题与排查技巧实录4.1 训练不收敛怎么办一维悬臂梁问题其实不算难训练不收敛多半可以从几个方向排查。先看是直接四阶还是降阶。如果还在用 ŵ″″1 直接求残差我建议直接换降阶方案能少掉 80% 的麻烦。再看学习率Adam 默认的 1e-3 在大多数情况下可用但如果 loss 频繁出现 NaN把学习率降到 1e-4 试试。另外注意网络初始化不要手动给线性层设太大的初始权重PINN 默认的初始化通常够用。排查顺序我一般这样来先确认解析解没错避免“标准答案”本身就是错的。打印网络前向输出形状确认 U、T 不是标量而是二维列向量。检查硬约束映射固定端 x0 时 u 是否严格为 0自由端 x1 时 theta 是否严格为 0。看 pde2 单项损失θ″-1 残差如果一直很大说明网络还没有学会基本的均布载荷关系。如果损失卡在某个量级迟迟不下去检查采样点数是不是太少。N128 对这个例子足够但有人为了追求“快”改成 N16会明显发现曲线形状能看出来但精度上不去。反过来说一味加到 N1024 也不会带来质的提升因为点太多后随机采样的信息冗余变高训练变慢。4.2 边界条件满足但曲线形状不对如果你的预测曲线在边界处看起来合理但整体形态偏离解析解问题大概率出在控制方程残差和边界残差的比例失衡或者方程无量纲化尺度不对。我之前遇到过一次把所有物理量直接塞进损失函数EI 是 300、q 是 5、L 是 4结果 pde2 的损失量级在 10⁻³ 附近边界损失量级却在 10⁻⁶训练结果边界条件完美但位移完全不对。这就是数量级失衡的典型表现。解决方式就是无量纲化把 ŵ、x̂ 都变成 O(1) 的量所有损失在同一个尺度下竞争。还有一种情况是曲线形状有“波浪感”。这通常意味着网络容量不足或激活函数选择不当。先用 tanh 和三层隐藏层如果还是震荡可以适当加宽到 128但不要盲目加深到五六层。一维问题的真实复杂度很低网络容量超过需求后反而容易在小样本约束下产生过拟合式的波动。4.3 损失权重到底怎么定如果用软约束λ_pde 和 λ_bc 的取值会直接影响结果。理论上两个损失应该平衡但不同问题的真实数值尺度不同没有万能常数。常用做法是先让网络跑一个很短的前向记录两类损失的初始值然后令 λ_pde 1λ_bc 初始 pde 损失 / 初始 bc 损失让它们初始量级一致。之后如果发现边界被牺牲就把 λ_bc 适当调大。在本文的硬约束方案里其实根本不需要考虑这个权重这也是我推荐硬约束的重要原因。边界被零点因子固定后优化器只需要集中精力把内部方程配平不再需要在两类目标之间做折中。软约束不是不能用只是调试成本更高新手建议先用硬约束把流程跑通再去尝试软约束对比。4.4 网络结构、激活函数和采样点的影响网络宽度从 20 加到 64精度通常有明显提升但继续加到 256 对一维问题收益不大反而更容易过拟合。深度方面三层隐藏层足够太深会导致训练初期梯度信号传递困难。激活函数首选 tanh也可以试试 sin。sin 激活在部分高频问题上表现更好但训练对初始化敏感不建议新手一上来就用。采样点数方面每次随机采样 128 个点是我常用的配置。如果你的问题载荷变化更剧烈比如高梯度区域应该在局部加密或者采用自适应采样把模型误差大的区域加更多点。一维问题用均匀随机采样已经够用。如果你发现某个区域残差特别大也可以刻意在这个区间加采样点这在工程上是性价比很高的手动干预手段。4.5 还能怎么扩展这个例子这次做的是均布载荷下的等截面悬臂梁但实现框架可以非常自然地扩展。换成集中载荷时可以把载荷 q(x) 写成一个高斯形式的窄带近似近似成为分布载荷。或者更稳妥的做法是分段处理集中力作用点两侧分别设内部残差再用连续的位移和转角匹配条件作为额外边界约束。换变截面梁时把 EI 从常数改成 EI(x)在 θ″q(x)/EI(x) 里体现 x 依赖即可。换二维问题时把网络输入变成 x 和 y 两个坐标输出挠度曲面再把偏微分方程残差写成二维形式整体思路完全一致。也就是说这篇一维悬臂梁的 Python 实现本质上是 PINN 的一种最小可运行模板。控制方程怎么改、边界条件怎么加逻辑链条是不变的。我在实际调试这个例子的过程中最大的体会是“物理结构的约束是训练最好的正则化器”。同样一组网络参数如果只给数据点拟合可能过拟合到完全无法外推但一旦把梁方程和边界条件写进损失模型就被牢牢锁定在物理事实附近哪怕采样点少一点预测结果也不会离谱。这也是 PINN 相比普通监督学习的魅力所在。最后再分享一个小技巧如果你在训练过程中发现不管怎么调都达不到预期先把数据可视化做出来把预测挠度、解析解、损失曲线三张图放在一起看大多数问题一眼就能定位。一维问题花不了多少时间但它能把 PINN 的核心流程练得很熟。希望这篇拆解能帮你把这个入门例子真正跑通再遇到更高维的物理场问题时不至于手足无措。