案例4：使用残差分析进行离心泵故障诊断

本案例延续了案例3中基于稳态实验的故障诊断方法的情况，所不同的是本案例中的水泵运行于各种不同工况中，即转速不再限于2900 RPM

第一步：多速泵系统建模——基于残差分析的诊断

当水泵运行于转速变化快且变化范围大时，稳态泵扬程和扭矩方程不再能得到精确结果。此时，摩擦和其他损失的影响变得更加明显，而模型参数又依赖水泵转速。可行的方案便是建立一个黑盒模型，模型参数不必具有物理意义。该模型作为水泵已知行为的模拟器，模型的输出为从对应的测量信号中提取得到的残差。残差的特性，如均值、方差和功率被用来区分正常和故障操作工况

由如下系统简图可知，使用静态泵扬程方程，动态泵-管道系统方程，可以计算出图中标红的4种残差

系统图中四个灰色模型的表示如下：

  • 静态泵模型：Δp̂(t) = θ1ω2(t) + θ2ω(t)
  • 动态管道模型：$\hat Q(t)=\theta_3+\theta_4\sqrt{\Delta p(t)}+\theta_5\hat Q(t-1)$
  • 动态泵-管道模型：$\hat{\hat Q}(t)=\theta_3+\theta_4\sqrt{\Delta\hat p(t)}+\theta_5\hat{\hat Q}(t-1)$
  • 动态逆泵模型：M̂motor(t) = θ6 + θ7ω(t) + θ8ω(t − 1) + θ9ω2(t) + θ10M̂motor(t − 1)

模型参数θ1, ..., θ10与水泵转速相关。

本案例中，模型参数将通过分段线性估计方法得到，转速区间分为以下3个档：

  1. ω ≤ 900RPM
  2. 900 < ω ≤ 1500RPM
  3. ω > 1500RPM

正常情况下，水泵运行于由闭环控制器控制的0-3000 RPM参考转速区间内，该参考转速输入属于修正伪随机位序列信号。电机扭矩、泵扭矩、泵转速和压强数据在10 Hz的采样频率下采集

以下代码用于加载测量得到的这些信号并且画出参考转速和实际转速的图像：

以下代码定义泵转速划分的三种区间范围：

第二步：模型识别

静态泵模型识别

使用泵转速ω(t)和压差Δp(t)作为静态泵模型方程的输入输出得到模型参数θ1和θ2的估计值，此处函数staticPumpEst用来执行估计过程：

动态管道模型识别

类似于静态泵模型参数估计中的做法，给定方程$\hat{Q}(t)=\theta_3+\theta_4\sqrt{\Delta p(t)}+\theta_5\hat Q(t-1)$，通过使用流出率Q(t)和压差Δp(t)的测量值作为方程的输入输出，可以估计得到参数θ3，θ4，θ5的估计值。此处函数dynamicPipeEst用来执行方程估计过程：

不同于静态泵模型，动态管道模型具有对液体流出率动态依赖的关系。为了模拟该模型在不同转速区间下的表现，从控制系统工具箱中的LPV System模块开发得到一个分段线性模型。该模型被命名为LPV_pump_pipe，函数simulatePumpPipeModel用来运行该模拟过程：

动态泵-管道模型识别

动态泵-管道模型使用同以上部分相同的参数（θ3,θ4,θ5），不同的是模型模拟使用静态泵模型中估计得到的压差，而不是测量得到的压差，其他步骤如同动态管道模型识别过程中的操作

该结果与使用压差测量值得到的结果几乎一致（比较上图与上上图）

动态逆泵模型识别

参数θ6,…,θ10的估计符合同样的准则，回归分析测量得到的扭矩值，前一时刻的扭矩值以及速度值得到参数。然而，此时分速度区间分段线性模型并不能很好的拟合数据，因此使用了一种不同的黑盒方法。该方法使用有理回归器识别得到一个经数据拟合后的非线性的自回归各态历经模型identifyNonlinearARXModel，过程如下：

得到的非线性的自回归各态历经模型如下：

使用该模型模拟估计得到的电机扭矩值和测量得到的电机扭矩值绘图如下：

第三步：残差生成

由上述四种方法得到了四个模型的估计值，它们与对应测量值的差便可作为残差用作后续分析的输入

可视化得到的残差：

第四步：残差特征提取

从第三步中得到的残差，可以用来做特征提取（如常规的信号处理方法，最大峰值、信号方差等），得到的特征被用于故障类型分离

考虑对使用伪随机位序列作为输入，实现的一组50次实验的数据采集。每组实验分别按照以下10种模式进行：

  1. 健康水泵
  2. 故障1：间隙距磨损
  3. 故障2：叶轮出口小沉积
  4. 故障3：叶轮进口沉积
  5. 故障4：叶轮出口物理磨损
  6. 故障5：叶片折断
  7. 故障6：气穴现象
  8. 故障7：转速传感器偏差
  9. 故障8：流量计偏差
  10. 故障9：压力传感器偏差

下面加载实验数据（50 × 1 cell型数据）：

接下来生成针对每种实验工况下每组数据得到的残差值，可使用以下helperComputeEnsembleResidues函数

基于此函数的基础上，可以出10组实验得到的残差值，如健康水泵模型的残差计算如下所示：

因针对所有数据的计算较为耗时，这里直接导入计算后的结果，如下：

每组实验得到数据形式为50 × 4 cell型的残差数据，4表示共有4种残差值，对于每个cell中的元素，它们是形式为550 × 1时序数据

得到残差后，能有效区分故障类型的特征无法直观得到，因此使用函数helperGenerateFeatureTable生成多种常用的备选特征，诸如均值、最大峰值、方差、谷值、1-范数等，当然在此之前所有残差数据都应使用MinMaxScaler归一化到健康模式数据的范围内

每个残差特征表形式为50 × 20 table数据，20表示所有4种残差得到的特征（每种残差提取5种特征），50表示50次实验

下一步，将所有特征数据打上标签并拼接到一起，操作如下：

第五步：分类器设计

使用散点绘图可视化模式分离

有了残差特征训练数据后，可以通过散点绘图目视的方法观察哪些特征适合用来分离故障。例如对于故障1：间隙距磨损，可以生成一幅标记上“Healthy”和“ClearanceGapWear”标签的特征散点图：

依稀可以看出，第1列和第17列的特征可以有效地分离健康和故障1工况

接下来仔细分析这两个特征：

继续查看这两个特征构成的散点图：

绘图显示，特征1：“Mean1”和特征17：“OneNorm1”可以被用来区分健康和间隙距故障工况，同样的分析方法可以应用于其他故障模式分析。对于所有的故障模式，都可以找到一组所有特征的子集特征来区分这些故障模式。然而，由于统一特征可能会被多种故障模式影响，所以故障分离会更加困难。例如特征1和特征17就并不是对所有故障类型都能区分，相反对于传感器偏差类故障来说，故障分离就简单一些，因为针对它们的诊断，有很多的有效特征可选

对于定义的三种传感器偏差类故障，可以从散点图的手动查看中选出合适的特征

画出这五种特征在三种传感器偏差故障和健康时的散点图：

图中可以看出，这四种工况模式可以被单个或多个特征区分

接下来可以用选定的特征子集，训练一个含有20棵子树的装袋树分类器，来区分这几种传感器偏差故障

可以看出，最终错分率低于3%，证明可以选择一部分特征子集来区分故障子类（如归属于传感器偏差类的故障）

使用Classification Learner APP进行多分类

上述部分主要通过人工检测的方式选择特征划分故障类型，这种方法比较原始并且也不一定能覆盖所有的故障类型。一种更时髦的方法是使用Matlab统计和机器学习工具箱中的Classification Learner APP，以便快速建模和比较结果，操作过程如下：

导入数据并定义验证方式（此处选择20%样本留出）

训练并分析结果

从结果可看出，

  • 三次支持向量机具有最好的分类效果，准确率达到了86%
  • 间隙距磨损故障易和叶轮出口沉积故障的识别混淆
  • 物理磨损易和叶片折断和健康状态两种模式混淆，推测可能是跟物理磨损的程度有关，轻微时更像健康状态，严重时更像叶片折断

要想导出生成的结果模型，可直接使用“导出”功能

第六步：总结

本案例介绍了基于比较测量数据和模型估计得到参数之间残差的方式建模分析离心泵故障的过程。当建模任务教复杂，且模型参数依赖于操作工况时，该方法是对参数估计和基于轨迹方法的很好补充

基于残差的故障诊断方法遵循以下步骤：

  1. 使用物理知识和黑盒系统识别技术对测量数据和系统输出的动态过程建模
  2. 计算测量信号和模型输出信号的残差，该残差在方差较大时可能需要经过平滑化、滤波处理
  3. 对每种残差提取残差特征，如峰值、功率、谷值等信号特征
  4. 选择可用特征，使用异常检测和分类技术进行故障检测
  5. 并不是所有残差和由它导出的特征都对每种故障的识别有效。对特征绘制散点图和柱状图有利于找到哪些特征对哪些特定的故障类型分离有用。对可被用来分离故障的特征挑选和评估过程属于可迭代过程（即可不断优化）

参考：https://ww2.mathworks.cn/help/predmaint/ug/fault-diagnosis-of-centrifugal-pumps-using-residual-analysis.html