☰
物理信息神经网络求解一维悬臂梁挠度问题
2026/10/22 13:37:50 网站建设 项目流程

最近在做结构力学相关研究时,随手拿了一个经典的一维悬臂梁挠度问题来测试物理信息神经网络(PINN)。这个题目看起来简单,但作为PINN入门的练手案例非常合适——它既有明确的控制方程,又有典型的边界条件,而且解析解是现成的,可以很方便地验证神经网络算得对不对。这篇文章把我从搭建模型到调通代码的完整过程记录下来,包括网络结构怎么设计、边界条件怎么处理、损失函数每一项是怎么算的、以及我在训练中踩过的几个坑。如果你刚接触PINN,想找一个能快速跑通又看得懂原理的案例,这篇应该能帮到你。

1. 问题定义与整体思路

1.1 为什么选一维悬臂梁

悬臂梁是结构力学里最基础也最经典的构件之一。一根梁一端固定、一端自由,在自由端施加一个集中力,梁会发生弯曲变形,我们需要求解的是梁在每一个位置处的挠度。

这个问题之所以适合做PINN的入门案例,有几个原因:

  • 控制方程是一个四阶常微分方程,形式固定,不像二维三维问题那样需要处理复杂的几何边界。
  • 边界条件非常典型:固定端位移为零、转角为零,自由端弯矩为零、剪力等于施加的外力。
  • 解析解可以直接由材料力学公式给出,用来验证神经网络的计算精度非常方便。

我当时的目标很简单:用神经网络代替传统数值方法,输入是梁上的位置坐标x,输出是对应的挠度w(x),通过让网络输出满足控制方程和边界条件来训练模型。

1.2 为什么用PINN而不是有限元

传统方法里,求解这种问题最常用的是有限元法。有限元需要划分网格,需要对单元进行离散和组装,对于复杂几何和复杂边界,网格生成本身就是一个不小的工作量。PINN的思路完全不同:它把控制方程和边界条件直接编码到损失函数里,通过优化神经网络的参数来满足这些方程,整个过程不需要任何网格划分,也不需要单元矩阵组装。

PINN有一个很有吸引力的特点:只要损失函数下降到足够小,网络给出的是一个连续函数解,而不是只在离散节点上有值。这意味着你可以在任意位置查询挠度,不需要插值。

当然,PINN也有它的短板。四阶方程本身对网络的二阶甚至三阶求导要求很高,训练过程中梯度传播路径比较长,收敛速度和稳定性往往不如成熟的有限元求解器。不过作为方法验证和学习理解,这个案例的价值在于把PINN的整个工作流程完整走通。

2. 控制方程与网络结构设计

2.1 悬臂梁挠度的控制方程

一维悬臂梁在自由端受集中力F作用时,其挠度控制方程为:

$$EI \frac{d^4 w}{dx^4} = 0$$

其中E是弹性模量,I是截面惯性矩,w是挠度。

这个方程看起来很简单,但关键在于边界条件的组合。两端边界条件分别为:

  • 固定端x=0处:w=0,dw/dx=0
  • 自由端x=L处:弯矩为零,即d²w/dx²=0;剪力等于外力,即EI·d³w/dx³=F

为了计算方便,我对方程做了无量纲化处理。设梁长L=1m,EI=1,自由端剪力归一化为1。这样控制方程就简化为:

$$ \frac{d^4 w}{dx^4} = 0 $$

边界条件变为:

  • w(0)=0
  • dw/dx(0)=0
  • d²w/dx²(1)=0
  • d³w/dx³(1)=1

无量纲化的好处不仅仅是数值上更稳定,更重要的是避免了不同物理量纲带来的梯度不平衡问题。在实际计算中,如果直接用真实物理参数,挠度可能只有毫米级别,而坐标x和挠度w之间的量级差会导致损失函数中各项的贡献严重失衡,训练起来非常困难。

2.2 网络结构与激活函数选择

对于这个一维问题,网络结构不需要很复杂。我选用了一个四层全连接网络,每层50个神经元,输入层维度为1(坐标x),输出层维度为1(挠度w)。

激活函数的选择上有讲究。常见的ReLU虽然计算简单,但它的一阶导数是分段常数,二阶导数在分段点处为0,三阶导数几乎到处都是0,根本无法用于四阶方程的损失函数计算。我当时试过ReLU,结果是损失函数根本降不下去,因为网络的导数信息在反向传播中丢失了。

换用tanh激活函数之后,问题得到了解决。tanh是光滑函数,理论上可以无限次求导,而且它在各个阶数的导数都有非零值,能够正常地通过自动微分计算损失函数中的物理残差。

注意:对涉及高阶导数的PINN问题,激活函数必须选择光滑函数。tanh、sigmoid、swish这类函数都可以用,从而避免ReLU族函数导数阶数限制的问题。

网络结构设计的另一个关键点是如何让边界条件更容易被满足。我采用了硬约束的方式,也就是把网络输出做一次变换,使得输出自动满足位移和转角的边界条件:

$$w_{net}(x) = x^3 \cdot NN(x)$$

这里NN(x)是原始网络输出。因为当x=0时x³=0,所以无论如何,w_net在x=0处的值和一阶导数都恒为0。这样就无需在损失函数中再考虑位移和转角的边界约束,只需要约束自由端的弯矩和剪力条件。

这个做法的优点很明显:损失函数中需要约束的项变少了,训练目标更清晰,收敛速度也有明显提升。代价是网络输出和真实解答的表达空间不再完全等价,但在实际测试中,只要网络容量足够,这个变换并不会限制解的精度。

3. 边界条件处理与损失函数构建

3.1 如何把边界条件编码进模型

如前面所说,固定端的位移和转角边界条件通过输出变换被硬编码进网络结构了。这样剩下的只有自由端的两个自然边界条件需要处理。

在损失函数中分别定义自由端的弯矩残差项和剪力残差项:

  • 弯矩条件:网络输出在x=1处的二阶导数应为0
  • 剪力条件:网络输出在x=1处的三阶导数应为1

把这两个残差平方后加入总损失函数中,就实现了对边界条件的软约束。

3.2 损失函数的三项组成

最终的总损失函数由三部分组成:

$$L = L_{PDE} + \lambda_1 L_{moment} + \lambda_2 L_{shear}$$

其中第一项是内部点的控制方程残差,第二项是自由端弯矩残差,第三项是自由端剪力残差。

内部点损失项的计算方法是:在[0,1]区间内随机采样若干点,计算网络输出对这些点的四阶导数,与目标值0做差的平方和。采样点的数量和分布对收敛效果有很大影响,我使用的是均匀随机采样,每个训练轮次重新采样,相当于每次迭代都在使用一组全新的样本点。

权重系数λ的选择也是调参的重点。一开始我将两个边界条件的权重都设为1,训练后发现在边界附近的解答偏差比较大。后来把λ从1逐级调整到10甚至50,才找到了一个合适的平衡点。实际效果是,边界权重越大,边界条件满足得越精确,但过大会导致内部点残差被忽略,整体解的光滑性变差。

4. 代码实现与关键细节

4.1 导入依赖和定义超参数

我用PyTorch作为深度学习框架,因为它的自动微分功能很强大,能方便地计算任意阶数的导数。以下是完整的代码实现。

import torch import torch.nn as nn import numpy as np import matplotlib.pyplot as plt # 设置随机种子,保证可复现性 torch.manual_seed(42) np.random.seed(42) # 超参数定义 n_epochs = 5000 n_samples = 100 # 每次采样的内部点数量 lr = 1e-3 # 学习率 lambda_moment = 10.0 # 弯矩边界权重 lambda_shear = 10.0 # 剪力边界权重

4.2 网络定义

class FCN(nn.Module): def __init__(self, n_input=1, n_output=1, n_hidden=50, n_layer=4): super(FCN, self).__init__() layers = [nn.Linear(n_input, n_hidden)] for _ in range(n_layer - 1): layers.append(nn.Linear(n_hidden, n_hidden)) layers.append(nn.Linear(n_hidden, n_output)) self.activation = nn.Tanh() self.layers = nn.ModuleList(layers) def forward(self, x): for i, layer in enumerate(self.layers): x = layer(x) if i < len(self.layers) - 1: x = self.activation(x) return x

4.3 物理约束与损失计算

def pde_residual(model, x): """计算控制方程残差(四阶导 - 0)""" x = x.requires_grad_(True) w = model(x) # 一阶导 dw_dx = torch.autograd.grad(w, x, grad_outputs=torch.ones_like(w), create_graph=True)[0] # 二阶导 d2w_dx2 = torch.autograd.grad(dw_dx, x, grad_outputs=torch.ones_like(dw_dx), create_graph=True)[0] # 三阶导 d3w_dx3 = torch.autograd.grad(d2w_dx2, x, grad_outputs=torch.ones_like(d2w_dx2), create_graph=True)[0] # 四阶导 d4w_dx4 = torch.autograd.grad(d3w_dx3, x, grad_outputs=torch.ones_like(d3w_dx3), create_graph=True)[0] return d4w_dx4, w, dw_dx, d2w_dx2, d3w_dx3 def compute_loss(model): # 内部点残差 x_interior = torch.rand(n_samples, 1) * 1.0 # [0, 1]区间均匀采样 d4w, _, _, _, _ = pde_residual(model, x_interior) loss_pde = torch.mean(d4w**2) # 自由端x=1处边界条件 x_boundary = torch.ones(1, 1) _, _, _, d2w_b, d3w_b = pde_residual(model, x_boundary) loss_moment = torch.mean(d2w_b**2) # 弯矩应为0 loss_shear = torch.mean((d3w_b - 1.0)**2) # 剪力应为1 # 总损失 loss = loss_pde + lambda_moment * loss_moment + lambda_shear * loss_shear return loss, loss_pde, loss_moment, loss_shear

4.4 训练流程

def train(model): optimizer = torch.optim.Adam(model.parameters(), lr=lr) loss_history = [] for epoch in range(n_epochs): optimizer.zero_grad() loss, loss_pde, loss_moment, loss_shear = compute_loss(model) loss.backward() optimizer.step() loss_history.append(loss.item()) if epoch % 1000 == 0: print(f"Epoch {epoch}: loss={loss.item():.6f}, " f"pde={loss_pde.item():.6f}, " f"moment={loss_moment.item():.6f}, " f"shear={loss_shear.item():.6f}") # 创建模型并训练 model = FCN() train(model)

4.5 解析解对比与可视化

一维悬臂梁在自由端受集中力的挠度解析解表达式为:

$$w_{exact}(x) = \frac{F x^2}{6EI}(3L - x)$$

代入无量纲化参数F=1, EI=1, L=1,得到:

$$w_{exact}(x) = \frac{x^2}{6}(3 - x)$$

用训练好的模型在所有位置推断挠度值,和解析解绘制在同一张图里对比:

# 测试点 x_test = torch.linspace(0, 1, 200).reshape(-1, 1) with torch.no_grad(): w_pred = model(x_test).squeeze().numpy() # 解析解 w_exact = (x_test.squeeze().numpy()**2 * (3 - x_test.squeeze().numpy())) / 6 # 绘图 plt.figure(figsize=(8, 5)) plt.plot(x_test.squeeze().numpy(), w_exact, 'k-', label='Exact Solution') plt.plot(x_test.squeeze().numpy(), w_pred, 'r--', label='PINN Prediction') plt.xlabel('x') plt.ylabel('w(x)') plt.legend() plt.grid(True, alpha=0.3) plt.title('一维悬臂梁挠度:PINN vs 解析解') plt.show() # 计算最大误差 error = np.abs(w_pred - w_exact).max() print(f"最大绝对误差: {error:.6f}")

5. 训练调试与常见问题

5.1 训练过程中的几个关键观察

训练初期,损失函数的下降速度非常快,前几百轮就能看到明显的变化。这是因为硬约束固定端边界后,网络主要需要学习的自由度集中在了如何满足控制方程和剩余边界条件上,寻优空间相对有限。

不过在损失下降到一定程度后,会进入一个平台期。这个阶段损失值在小范围内波动,不再明显下降。这是正常的,PINN的损失函数本质上是多个约束项的加权叠加,各约束之间很难同时达到完美平衡,最终会在某个折中的状态稳定下来。

我在测试中发现,最终模型在内部点的挠度预测最大误差在0.5%左右,对于工程精度要求不高、主要用于方法演示的场景是完全够用的。如果要提高精度,可以增加网络宽度或深度、增加采样点数量、采用自适应采样策略,或者增大边界条件的权重。

5.2 一个容易被忽略的细节:梯度计算

在进行高阶导数计算时,create_graph=True这个参数一定不能忘。如果忘记设置,计算出来的梯度不会参与后续的反向传播,训练时网络参数无法正常更新,损失值会一直卡住不动。我当时第一次跑代码就在这里栽了跟头,找了半天才发现是这个参数的问题。

# 错误写法:无法正常训练 dw_dx = torch.autograd.grad(w, x, grad_outputs=torch.ones_like(w))[0] # 正确写法:保留计算图用于后续反向传播 dw_dx = torch.autograd.grad(w, x, grad_outputs=torch.ones_like(w), create_graph=True)[0]

5.3 四阶导数的数值特性

四阶导数在数值上是比较敏感的。如果网络输出在局部有细微的振荡,四阶导数会被放大得非常明显。我在调试中发现,在训练后期,内部点的四阶导残差虽然在减小,但仍然存在一些小的高频分量。这种情况下可以试试降低学习率,或者对采样点进行平滑处理。

更有效的方法是给损失函数中的PDE残差项换一种归一化方式,比如用相对误差代替绝对误差。我在后续的测试中发现,将PDE残差除以某个参考值后,训练稳定性有比较明显的提升,特别是在网络参数随机初始化的情况下。

5.4 采样点的选择策略是影响收敛速度的关键因素

我最初用的是整个[0,1]区间均匀随机采样。但观察训练过程后发现,靠近自由端的区域对解的准确性影响更大,因为那里有弯矩和剪力的边界条件需要满足。

后来我尝试了一种简单的自适应采样:先训练几百轮,然后统计每个位置的PDE残差分布,在残差较大的区域增加采样点密度。这种方法让同样的训练轮数下最大误差几乎缩小了一半。PINN的采样策略其实挺有讲究,均匀采样、随机采样、自适应采样、带权重的采样各有适用场景,根据问题特性选择合适的方法能事半功倍。

6. 经验总结与扩展方向

6.1 我对PINN求解这个问题的几点体会

用PINN算悬臂梁挠度这个题目,整体跑通不难,但要算得准确、稳定,还是有不少细节需要注意。

第一个体会是边界条件的处理方式直接影响训练效率。硬约束和软约束相比,虽然增加了网络设计的复杂度,但给优化问题排除了一部分不确定性。对于边界条件比较规则的问题,优先考虑能不能把它直接嵌入到网络输出中。

第二个体会是损失函数的权重分配真的很重要。各个物理约束之间天然存在数值尺度差异,如果不做加权调整,训练结果会被数值量级更大的那一项主导。实际操作中建议先用若干个不同权重跑几次短训练,观察各项残差的收敛情况,再根据观察结果调整权重配置。

第三个体会是网络结构和参数量不是越大越好。对于一维简单问题,四层宽度50的网络已经足够。盲目增加网络规模不仅让训练变慢,还可能出现过拟合式的振荡,反而让导数的计算结果变差。

6.2 这个案例还能怎么扩展

这个一维悬臂梁案例虽然简单,但它包含的方程形式、边界条件处理方式、损失函数构建思路,可以作为一系列更复杂问题的起点。

比如可以推广到变截面梁,也就是截面惯性矩I(x)沿梁长变化,控制方程中的系数不再是常数,此时四阶导项会包含位置相关的系数,但在PINN的框架下,修改起来仅仅是损失函数里增加一个系数函数而已。

还可以考虑加上分布荷载、非线性材料本构关系,这类问题在传统数值方法中往往需要反复迭代,但在PINN框架中,只需要调整控制方程的残差表达式即可。

更进一步的,可以用类似框架处理二维弹性力学问题,甚至扩展到基于观测数据做参数辨识。比如用已知的挠度测量值来反推弹性模量E或惯性矩I,这就从正问题走向了反问题,PINN在这个方向上的优势会更明显。

总的来说,PINN为求解物理问题提供了一种全新的思路。传统方法需要网格、离散、组装,而PINN只需要一个神经网络和一套物理约束。它仍然有计算效率、收敛性等方面的不足,但在某些问题上已经展示出独特的价值。这个悬臂梁案例是我入门PINN的第一个完整实现,也让我对这个方向的潜力有了更直观的认识。

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询