ARTICLE · INTELLIGENCE

战地情报 · 详情页

来自尧图项目组的一线实战观察与深度解析

物理约束神经网络PINN故障诊断:MATLAB完整实现与损失函数设计

物理约束神经网络PINN故障诊断:MATLAB完整实现与损失函数设计 简介面向工业故障诊断与智能运维领域研究人员与工程师这份MATLAB项目围绕物理约束神经网络PINN实现故障诊断分类预测将物理机理与数据驱动深度融合重点解决小样本、噪声干扰与模型可解释性不足等问题。资源包含1个docx文档约92KB文档从物理建模、数据生成、网络搭建、复合损失函数设计到训练优化与GUI开发逐层展开并配有完整可运行代码及代码详解。已有275人学习下载。阅读后可掌握通过自动微分将物理方程嵌入损失函数的实现技巧、多目标优化平衡策略以及如何构建集数据处理、模型训练、预测评估与可视化于一体的诊断系统适合具备一定MATLAB和机器学习基础、希望将AI与物理模型结合落地的开发者参考。 做故障诊断的人这两年应该没少被“物理约束神经网络”这个词刷屏。我第一次看到PINNPhysics-Informed Neural Network是在流体力学方向的论文里当时第一反应是“这玩意跟故障诊断有什么关系”。直到后来做旋转机械的小样本故障分类纯数据驱动的模型频频在变工况下翻车我才意识到把物理规律硬塞进神经网络可能才是工业场景里真正能落地的解法。这篇就用一个完整的MATLAB项目实例把PINN怎么跟故障诊断分类结合、损失函数怎么写、GUI怎么做全部拆开讲透。项目已跑通代码结构完整既有全连接网络作为骨干又把系统动力学方程的残差作为物理约束项加入训练适合有MATLAB基础、想在故障诊断方向引入PINN的工程师和研究生参考。1. 为什么故障诊断需要物理约束纯数据模型的三个死穴先聊聊我踩过的坑。以前做轴承故障分类用的是CNNLSTM那套经典组合数据集是西储大学的公开数据训练集测试集随机划分准确率能到99%以上看着挺美。但一换到实测数据准确率直接掉到85%以下。问题出在三个地方第一随机划分数据意味着同一工况下的样本被同时分进了训练集和测试集模型其实“记住”了工况特征而不是故障特征。第二纯数据驱动模型完全不理解轴承的动力学特性它学到的决策边界在特征空间中非常脆弱稍微有点噪声就翻车。第三小样本场景下数据量根本撑不起复杂网络的学习需求过拟合是必然的。PINN解决的就是第三个问题同时捎带缓解前两个。它的核心思路是别让网络自己瞎学把系统的物理规律——比如转子动力学方程、故障特征频率公式——写成残差项直接加到损失函数里。网络输出的分类结果不仅要符合训练数据的标签还要满足物理方程的约束。这就相当于给网络加了一个“物理正则化器”让它在数据稀疏的地方也能往正确的方向收敛。用大白话说普通神经网络是个“考试全靠刷题”的学生PINN是个“既刷题又背物理定律”的学生。题目没见过的时候后者至少能靠物理规律蒙对方向。2. PINN的核心原理与损失函数构造2.1 网络输出不再是裸的分类概率在标准分类任务里网络的最后一层接softmax输出一个概率分布然后直接用交叉熵算损失反向传播完事。PINN不是这么干的。它的损失函数长这样L_total L_data λ * L_physics其中L_data就是普通的交叉熵分类损失L_physics是物理约束损失λ是权重系数。但问题来了物理方程通常是描述连续时间动力学的比如二阶微分方程而分类任务的输出是离散的故障类别。这两者怎么对齐这里有个设计巧思。我的做法是网络主干输出一个特征向量接两个分支。分支一是softmax分类头输出故障类别概率分支二是物理头的输入把中间特征经过一个解码层还原成系统的物理量——比如振动位移、速度、加速度。然后把这个还原出来的位移序列代入已知的动力学方程计算残差残差的均方误差就是L_physics。换句话说网络中间层学到的特征不仅要能区分故障类型还要能重构出符合物理规律的信号。这两件事共享同一个特征空间互相约束特征就比纯分类任务学到的扎实得多。2.2 物理残差的具体形式以最常见的单自由度转子系统为例振动微分方程是m * x(t) c * x(t) k * x(t) F(t)其中m是质量c是阻尼k是刚度F(t)是激振力。当系统发生不同故障时这些参数会变化。比如裂纹会导致刚度k下降不对中会产生二倍频激振力松动会引入非线性项。物理残差定义为residual m * d2x_dt2 c * dx_dt k * x - F;如果网络重构出的x(t)符合动力学规律residual就接近零。如果网络为了分类而强行扭曲了特征residual就会变大体现在L_physics上就是梯度惩罚。2.3 MATLAB里怎么算二阶导数这可能是很多人在MATLAB里实现PINN时卡住的第一关。用符号微分太慢。用数值差分精度不够。正确做法是自动微分对应MATLAB的dlgradient函数。% x是dlarray类型t是时间 % 网络输出重构位移 x_pred model(parameters, t); % 一阶导数 dx_dt dlgradient(sum(x_pred), t); % 二阶导数 d2x_dt2 dlgradient(sum(dx_dt), t);注意dlgradient的书写姿势第一个参数必须是标量或对每个元素求梯度后reduction成一个标量不能直接对向量求导。我一开始就是没做sum直接dlgradient(x_pred, t)报错报了半天才看明白文档。3. MATLAB完整实现训练循环与物理约束注入3.1 数据集准备与预处理我用的是仿真加实验混合数据。先按三类故障生成仿真信号正常状态、轴承外圈故障、轴承内圈故障。每类500个样本每个样本是1024个采样点的振动信号。采样频率12kHz故障特征频率按理论公式设置。% 生成正常状态信号 fs 12000; % 采样频率 t (0:1023) / fs; x_normal sin(2*pi*50*t) 0.1 * randn(size(t)); % 外圈故障在特征频率处加冲击 f_outer 107.3; x_outer sin(2*pi*50*t) 0.8 * exp(-20*mod(t, 1/f_outer)) .* sin(2*pi*2000*t) 0.1 * randn(size(t)); % 内圈故障冲击幅值受转频调制 f_inner 162.2; x_inner sin(2*pi*50*t) 0.8 * (1 0.5*sin(2*pi*30*t)) .* exp(-20*mod(t, 1/f_inner)) .* sin(2*pi*2000*t) 0.1 * randn(size(t));特征提取这里要多说两句。直接拿1024点原始信号喂网络也可以但训练效率低而且物理约束不好设计。我把每段信号转成频域特征和时域统计特征拼接峰值、均方根、峭度、频域重心、2倍频能量占比。一组特征总共12维既保留故障特征又方便后续物理约束设计。3.2 网络结构定义用dlnetwork定义网络。结构不复杂输入层12维中间两层全连接各64个神经元ReLU激活然后是特征层32维接两个分支。% 主网络层 layers [ featureInputLayer(12, Normalization, zscore) fullyConnectedLayer(64, Name, fc1) reluLayer(Name, relu1) fullyConnectedLayer(64, Name, fc2) reluLayer(Name, relu2) fullyConnectedLayer(32, Name, feature) reluLayer(Name, relu3) ]; % 分类分支 classify_layers [ fullyConnectedLayer(3, Name, class_out) softmaxLayer(Name, class_softmax) ]; % 物理分支从特征层重构动力学参数和质量/阻尼/刚度 physics_layers [ fullyConnectedLayer(16, Name, phy_fc) reluLayer(Name, phy_relu) fullyConnectedLayer(3, Name, phy_out) % 输出 [m, c, k] ]; lgraph layerGraph(layers); lgraph addLayers(lgraph, classify_layers); lgraph addLayers(lgraph, physics_layers); % 把分类分支接到feature层 lgraph connectLayers(lgraph, relu3, class_out); % 把物理分支也接到feature层 lgraph connectLayers(lgraph, relu3, phy_fc); dlnet dlnetwork(lgraph);这里的物理分支输出不是重构位移而是输出动力学方程的系数m、c、k。好处是不用额外重构整个信号只要判别网络学出来的“故障特征”能对应到合理的物理参数范围就行。比如正常状态下k应该接近设计刚度值故障状态下k应该明显变小。这比硬构造位移再求导好实现得多训练也更稳定。3.3 自定义损失函数训练循环里最关键的是损失函数。损失分三块function [loss, grad] modelLoss(dlnet, dlX, dlT, dlY, lambda, m_true, c_true, k_true) % 前向传播 dlY_pred forward(dlnet, dlX); % 分类损失 class_probs dlY_pred{1}; loss_class crossentropy(class_probs, dlY); % 物理参数输出 phy_params dlY_pred{2}; m_pred phy_params(:, 1); c_pred phy_params(:, 2); k_pred phy_params(:, 3); % 物理约束损失预测参数与理论参数之间的偏差 % 再加上参数之间的物理约束关系 loss_phy mean((m_pred - m_true).^2) ... mean((c_pred - c_true).^2) ... mean((k_pred - k_true).^2); % 总损失 loss loss_class lambda * loss_phy; % 自动梯度 grad dlgradient(loss, dlnet.Learnables); end物理分支的监督信号从哪里来这里要说明一下。我做的是半监督式设计在训练集里除了故障标签外还额外附带了每个样本对应的近似物理参数。这些参数可以从振动信号的包络谱解算出来。这样物理约束有明确的目标值网络收敛更快。如果实际场景里标注不出物理参数也可以用另一种模式——只约束参数之间的相对关系比如“故障状态下k比正常值低至少10%”用约束不等式变体来写但训练难度会明显增加。3.4 完整训练循环训练循环的标准写法如下注意几个细节dlarray数据格式、mini-batch的构造、dlfeval的使用。% 数据准备 X feature_matrix; % [12, N]的double矩阵 T categorical(labels); % 故障标签 params_true param_matrix; % 每个样本对应的[m, c, k]真实近似值 % 转为dlarray dlX dlarray(X, CB); dlT dlarray(onehotencode(T, 2), CB); % 超参数 numEpochs 200; miniBatchSize 32; lambda 0.1; % 物理约束权重 learnRate 0.001; gradDecay 0.9; sqGradDecay 0.99; avgGrad []; avgSqGrad []; % 训练批处理 for epoch 1:numEpochs % 每个epoch打乱数据 idx randperm(size(dlX, 2)); dlX dlX(:, idx); dlT dlT(:, idx); params_true params_true(:, idx); for i 1:miniBatchSize:size(dlX, 2) idxBatch i:min(iminiBatchSize-1, size(dlX, 2)); dlXBatch dlX(:, idxBatch); dlTBatch dlT(:, idxBatch); paramsBatch params_true(:, idxBatch); [loss, grad] dlfeval(modelLoss, dlnet, dlXBatch, dlTBatch, ... paramsBatch, lambda, m_true, c_true, k_true); % Adam更新 [dlnet, avgGrad, avgSqGrad] adamupdate(dlnet, grad, ... avgGrad, avgSqGrad, epoch, learnRate, gradDecay, sqGradDecay); end % 每个epoch打印损失 if mod(epoch, 20) 0 fprintf(Epoch %d, Loss: %.4f\n, epoch, extractdata(loss)); end end这里m_true、c_true、k_true是常值参数如果是对应具体系统的标定值那没问题。如果是故障状态下的参数其实每个batch的样本对应的参数目标值不同直接用样本自带的paramsBatch就行。4. 实验对比加不加物理约束的差异4.1 小样本场景下的对比我用三组实验做的对比全部训练数据下不加物理约束、小样本下不加物理约束、小样本下加物理约束。前面两个是参照组第三个是PINN方案。小样本设置每类只取30个样本总共90个训练样本。测试集固定每类200个样本。实验结果如下表方案训练样本数/类测试准确率训练损失收敛情况全量数据 纯分类网络50098.7%正常收敛小样本 纯分类网络3078.3%过拟合明显小样本 PINNλ0.013084.7%轻微过拟合小样本 PINNλ0.13091.2%正常收敛小样本 PINNλ1.03086.5%欠拟合可以看到λ0.1时效果最好准确率比纯数据驱动提高了将近13个百分点。但λ也不是越大越好到1.0时物理约束太强网络反而没法拟合分类任务了。物理约束和分类损失之间的平衡在实际调参时非常重要。4.2 变工况泛化性测试再做了一个更有说服力的实验训练时只用0.15mm和0.21mm两档负载的数据测试时加上0.30mm负载。纯分类网络的准确率从98%跌到76%PINN模型只跌到89%。这个结果说明物理约束学到的特征不像纯数据驱动那样“贴着”训练集分布而是更接近系统本身的固有规律。5. GUI设计与交互逻辑5.1 界面布局用MATLAB App Designer设计了一个故障诊断交互界面整体分四个区域左侧数据加载区中间模型配置区右侧训练状态区底部结果显示区。界面的控件配置列表如下区域控件类型作用数据加载区按钮 编辑框选择数据文件夹、显示路径模型配置区下拉框 编辑框选择网络结构、设置λ、学习率、迭代轮数训练状态区坐标轴 文本区实时显示损失曲线和准确率曲线结果显示区表格 坐标轴显示混淆矩阵、分类结果、物理参数预测值5.2 核心回调函数训练按钮的回调是关键。这里要注意一点App Designer的代码模型里训练循环中要更新UI控件必须在回调函数里调用drawnow或直接更新组件属性否则界面会卡死不动。function TrainButtonPushed(app, ~) % 读取参数 lambda_val app.LambdaEdit.Value; lr app.LrEdit.Value; epochs app.EpochsEdit.Value; % 禁用训练按钮防止重复点击 app.TrainButton.Enable off; % 训练循环中实时更新损失曲线 for epoch 1:epochs % 训练一个epoch [loss, acc] trainOneEpoch(app.dlnet, ...); % 更新坐标轴 hold(app.LossAxes, on); plot(app.LossAxes, epoch, loss, b-, LineWidth, 1.5); drawnow; % 更新准确率显示 app.AccuracyLabel.Text sprintf(当前准确率: %.2f%%, acc*100); % 检查用户是否点击了停止按钮 if app.StopFlag break; end end % 重新启用训练按钮 app.TrainButton.Enable on; end这里有一个开发过程中的小坑如果你的训练循环里有耗时操作而界面没做异步处理App会显示“未响应”。解决方案有两种一是像上面的代码一样用drawnow强制刷新二是在后台使用parfeval做异步调用。训练规模小的话drawnow就够了异步调用处理数据回调时反而更绕。5.3 模型导入导出GUI里加了一个模型保存和加载功能用save和load函数序列化网络结构。% 保存模型 function SaveModelButtonPushed(app, ~) [filename, pathname] uiputfile(*.mat, 保存模型文件); if filename ~ 0 fullpath fullfile(pathname, filename); save(fullpath, dlnet, lambda, fs, feature_names); uialert(app.UIFigure, 模型保存成功, 提示); end end6. 常见问题与调试经验先说三行总体经验物理约束不是银弹λ的调节是最主要的功夫自动微分用不对最常见小样本本身有噪声PINN提升的幅度和问题的物理规律强度直接相关。下面列我踩过的坑6.1 dlgradient报错或结果为NaN最常见的原因是你没用dlfeval包裹损失函数调用或者输入不是dlarray类型。dlgradient必须配合dlfeval使用这是MATLAB自动微分的语法要求。% 错误写法 [loss, grad] modelLoss(dlnet, x, y); % 正确写法 [loss, grad] dlfeval(modelLoss, dlnet, x, y);还要检查你的dlarray维度标注格式批量维用B通道/特征维用C。维度标错的话很多计算方法都会出奇怪的结果。6.2 物理约束权重λ到底应该取多少我的建议是先做一个基线实验λ0纯分类网络记录准确率和损失曲线。然后从λ0.001开始以10倍步长递增分别训练观察。选择标准是在验证集准确率最高的λ此时训练集的分类损失应该和基线在同一个量级不能差太多。如果分类损失明显高于基线说明物理约束太强需要调低或者缩减物理约束项的范围。6.3 小样本下物理参数输出不准物理分支直接输出参数值时如果参数值的量级差异大比如k是几万c是几百网络收敛就会慢。建议对物理分支的输出做归一化处理比如输出相对标定标准值的比例。% 输出相对残差而不是绝对值 k_pred k_std * (1 tanh(phy_params(:, 3)));这样输出范围被限制在k_std的[0, 2]倍区间训练稳定性大幅提升。6.4 训练过程中损失曲线震荡PINN的损失函数是两个目标的加权和两个目标的梯度方向不一致时优化过程很容易震荡。解决办法一个是对两个损失分别做梯度裁剪二是在训练初期让物理约束权重从0逐渐增大到目标值课程学习策略。% 课程学习策略先学分类再逐步加入物理约束 lambda_current lambda_max * min(1, epoch/warmupEpochs);我在实际项目中用了这个策略epoch达到50时λ才升到最大值震荡明显减少。6.5 混淆矩阵里某两类总是互相分错这个不一定是网络的问题先检查数据。轴承外圈故障和内圈故障的频域特征在峰值位置上有重叠如果特征提取时没有加入相位信息仅仅靠幅值谱网络很难区分。建议在特征里增加包络谱的相位特征或者小波包分解后的子频带能量分布。7. 实际项目中的物理约束设计经验最后分享一个更贴近工程的观点很多人一上来就想把完整的系统动力学方程二阶、非线性、多自由度写进损失函数这通常不是好主意。完整的物理约束表达能力虽强但方程本身对参数误差非常敏感你的模型参数估计稍微偏一点物理残差就大得离谱然后分类任务直接被压垮。我推荐的做法是分层设计。第一层刚性的已知物理规律。比如故障特征频率与转频、轴承几何尺寸的关系这类公式是教科书级的、绝对可信可以直接写死。第二层带有参数估计的弱约束。比如“正常状态下振动能量主要集中在转频及其整数倍处”这个规律本身没问题但具体比例受工况影响较大需要设置容差带。第三层完全数据驱动。不要试图用物理规律约束一切把所有角落都填满是不可行的。PINN的失败案例绝大多数都是物理约束写得太满、太绝对没有给数据驱动留余地。我做过的成功调优往往是“物理约束把大方向框住数据驱动在细节上自由发挥”这跟带新人是一个道理——定方向给框架具体路径让他自己走。这个项目的所有代码和数据生成脚本已经整理成型核心训练循环不到100行GUI部分约500行。有基础的读者照着第四节的内容半天就能复现出基础版本。物理约束神经网络在故障诊断这条路上我觉得还是在早期但方向值得押注。本文还有配套的精品资源点击获取
RELATED READING

延伸阅读

更多一线实战笔记与深度复盘,助您持续精进