简介

及早发现和隔离机器运行中的异常现象，有助于减少事故，减少停机时间，从而节省操作成本。该方法涉及处理来自系统操作的实时测量，标记任何可能指向新出现的故障的意外行为。

本例探讨了以下故障诊断方面的内容:

识别系统行为的动态模型

在基于模型的检测方法中，首先使用测量的输入和输出数据建立相关系统的动态模型。一个好的模型能够准确地预测系统在未来一定时间范围内的响应。当预测不好时，残差可能很大，可能包含相关性。利用这些方面来检测故障的发生。

考虑一个受冲击和振动影响的建筑。根据系统的不同，振动的来源可以是不同类型的刺激，如阵风、与运行中的发动机和涡轮机的接触，或地面振动。这些冲击是对系统进行脉冲碰撞试验的结果，这些试验是为了充分激发系统而添加的。万博1manbetx仿真软件模型pdmMechanicalSystem.slx是这种结构的一个简单示例。这种激励来自于周期性的颠簸以及经过过滤的白噪声模拟的地面振动。系统的输出由一个受测量噪声影响的传感器收集。该模型能够模拟结构处于健康或损坏状态的各种情况。

sysA =“pdmMechanicalSystem”;open_system (sysA)%将模型设置为正常运行模式set_param ([sysA,/机械系统的]，“LabelModeActiveChoice”，“正常”）模拟系统并记录响应数据sim(sysA) ynormal = logsout.getElement(“y”) . values;

输入信号未测量;我们所记录的只是反应ynormal．因此，我们利用“盲识别”技术建立了系统的动态模型。特别地，我们建立了记录信号的ARMA模型作为系统的表示。当输入信号被假设为(滤波后)白噪声时，这种方法有效。由于数据受到周期性波动的影响，我们将数据分成几个部分，每个部分从波动的发生处开始。这样，每个数据段包含对一个碰撞和随机激励的响应——这种情况可以使用时间序列模型捕获，其中碰撞的影响归因于适当的初始条件。

Ts = 1/256;%数据采样时间Nr = 10;信号中碰撞的百分比N = 512;碰撞之间的数据长度%Znormal = cell(nr,1);为1:nr ysegment = ynormal.Data((ct-1)*N+(1:500));z = iddata(ysegment，[]，Ts);Znormal {ct} = z;每个节段只有一个凸起结束情节(znormal {}):绘制记录片段的采样图。标题(“测量反应段”）

将数据分成估计和验证两部分。

Ze = merge(znormal{1:5});Zv = merge(znormal{6:10});

用状态空间形式估计7阶时间序列模型党卫军()命令。通过交叉验证(检查与验证数据的拟合)和残差分析(检查残差不相关)来选择模型顺序。

Nx = 7;模型= sest(ze, nx，“形式”，“规范”，“t”、Ts);礼物(模型)查看具有参数不确定性的模型方程

离散时间识别状态空间模型:x (t + Ts) = x (t) + K e (t) y (t) = C x (t) + e (t) = (x1, x2) x3 x1 0 1 0 x2 0 0 1 x3 0 0 0 x4 0 0 0 x5 0 0 0 x6 0 0 0 x7 0.5548 + / - 0.04606 - -2.713 + / - 0.2198 - 5.885 + / - 0.4495 x4 x5 x6 x1 0 0 0 x2 0 0 0 x3 1 0 0 x4 0 1 0 0 0 1 x5 x6 0 0 0 x7 -8.27 + / - 0.5121 - 9.234 + / - 0.3513 - -7.956 + / - 0.1408 x7 x1 0 x2 0 x3 0 x4 0 x5 0 x6 1 x7 4.263 + / - 0.02599 C = (x1, x2) x3 x4 x5 x6 x7 y1 1 0 0 0 0 0 0 K = 1.025 + / - 0.01401日元x1 x2 1.444 + / - 0.0131 x3 1.907 + / - 0.01271 x4 2.385 + / - 0.01203 x5 2.857 + / - 0.01456 x63.26 +/- 0.0222 x7 3.552 +/- 0.0336采样时间:0.0039062秒参数化:CANONICAL形式，索引:7。参数及其不确定性用“idssdata”、“getpvec”、“getcov”表示。状态:终止条件:接近(本地)最小值，(范数(g) < tol)..迭代次数:7，函数计算次数:15在时域数据“ze”上使用SSEST进行估计。拟合估计数据:[99.07 99.04 99.15 99.05 99.04]%(预测焦点)FPE: 0.6242, MSE:[0.5974 0.6531 0.5991 0.5871 0.6496]模型“报告”属性中的更多信息。

模型显示参数估计的不确定性相对较小。我们可以通过计算测量信号估计频谱上的1-sd(99.73%)置信限来确定可靠性。

H =谱图(模型);showConfidence (h, 3)

尽管在较低频率的响应中有大约30%的不确定度，但置信区域很小。验证的下一步是看模型对验证数据集中的响应预测得有多好zv．我们使用超前25步的预测范围。

比较(zv，模型，25)%针对一个数据集的验证

图显示，该模型能够预测未来验证数据集25个时间步(= 0.1秒)的第一次实验中的响应，准确率为> 85%。若要查看与数据集中其他实验的适合度，请使用图轴的右键上下文菜单。

验证模型的最后一步是分析由它产生的残差。对于一个好的模型，这些残差应该是白色的，即对非零滞后显示统计上不显著的相关性:

zv渣油(模型)

残差在非零滞后时大多不相关。在导出了正常行为的模型之后，我们继续研究如何使用该模型来检测故障。

基于健康状态模型的残差分析故障检测

故障检测是对系统观察结果中不需要的或意外的变化进行标记。故障是由于系统的逐渐磨损或传感器故障或部件损坏引起的突然变化而引起的系统动力学变化。当故障出现时，在正常工作条件下得到的模型无法预测观测到的响应。这导致测量响应和预测响应之间的差值(残差)增加。这种偏差通常用较大的残差平方和或相关性的存在来标记。

将Simulin万博1manbetxk模型放在受损系统变体中进行仿真。我们使用单个凸点作为输入，因为残余测试需要白色输入，由于初始条件可能是瞬态的。

set_param ([sysA,/机械系统的]，“LabelModeActiveChoice”，“DamagedSystem”）;set_param ([sysA,/脉冲的]，“时间”，“5120”）%来强制只进行一次碰撞sim(sysA) y = logsout.getElement(“y”) . values;

resid(model, y.Data) set_param([sysA，/脉冲的]，“时间”，“512”）%还原原始

残差现在更大，并在非零滞后时显示相关性。这是故障检测背后的基本思想——创建一个剩余度量，并观察它如何随着每一组新的测量而变化。这里使用的是基于1步预测误差的简单残差。实际上，更高级的残差是根据应用程序的需要而定制的。

使用正常状态和恶化状态模型进行故障检测

一种更详细的故障检测方法是识别系统故障(损坏)状态的模型。然后我们可以分析哪个模型更有可能解释来自系统的实时测量。这种安排可以推广到各种类型故障的模型，因此不仅用于检测故障，而且用于确定是哪一个故障(“隔离”)。在本例中，我们采用以下方法:

我们已经模拟了系统的正常运行模式。现在我们模拟模型pdmMechanicalSystem在"生命尽头"模式下。在这种情况下，系统已经恶化到允许操作的最终状态。

set_param ([sysA,/机械系统的]，“LabelModeActiveChoice”，“DamagedSystem”）;sim(sysA) y = logsout.getElement(“y”) . values;Zfault = cell(nr,1);为ct = 1: nr z = iddata (y.Data ((ct-1) * N + (1:50 0)), [], Ts);Zfault {ct} = z;结束

现在我们创建了一组模型，每个模型对应一个数据段。如前所述，我们建立了状态空间形式的七阶时间序列模型。关闭协方差计算以提高速度。

mNormal = cell(nr,1);mFault = cell(nr, 1);Nx =阶(型号);opt = ssestOptions(“EstimateCovariance”, 0);为1:nr mNormal{ct} = ssest(znormal{ct}， nx，“形式”，“规范”，“t”， Ts, opt);mFault{ct} = sest(zfault{ct}， nx，“形式”，“规范”，“t”， Ts, opt);结束

验证模型mFault都很好地代表了错误的操作模式:

compare(merge(zfault{:})， mFault{:}， 25)

正常和错误的估计光谱如图所示。

Color1 =“k”;Color2 =“r”;ModelSet1 = cat(2,mNormal,repmat({Color1}，[nr, 1]))';ModelSet2 = cat(2,mFault,repmat({Color2}，[nr, 1]))';频谱(ModelSet1{:}，ModelSet2{:})轴([1 1000 -45 40])title(输出光谱(黑色:正常，红色:故障)）

谱图显示了两者的区别:受损模式的主共振被放大，但光谱在其他方面重叠。接下来，我们创建了一种定量区分正常状态和故障状态的方法。我们可以使用数据聚类和分类方法，例如:

在本例中，我们使用支持向量机分类技术。万博1manbetx来自两类模型的信息聚类(mNormal而且mFault)可以基于这些模型可以提供的不同类型的信息，如它们的极点和零点的位置，它们的峰值共振的位置或它们的参数列表。在这里，我们根据与两个共振相对应的极点位置对模式进行分类。为了聚类，我们将健康状态模型的极点标记为“good”，将故障状态模型的极点标记为“faulty”。

ModelTags = cell(nr*2,1);% nr为数据段数ModelTags(1:nr) = {‘好’};ModelTags(nr+1:end) = {“错误”};ParData = 0 (nr*2,4);plist = @ (p)[真正的(p(1)),图像放大(p(1)),真正的(p(3)),图像放大(p (3)));占主导共振的%极点为ct = 1:nr ParData(ct，:) = plist(esort(极点(mNormal{ct})));ParData (nr + ct:) = plist (esort(杆(mFault {ct})));结束cl = fitcsvm(ParData,ModelTags，“KernelFunction”，“rbf”，.．.“BoxConstraint”正,“类名”, {‘好’，“错误”});cl.ConvergenceInfo.Converged

符合逻辑的1

cl是分离训练数据的SVM分类器ParData分为好区和坏区。使用预测该方法的分类器1可以分配一个输入nx × 1的向量到两个区域之一。

现在，我们可以测试分类器的预测(正常vs损坏)，从一个参数以连续方式从健康(mode = ' normal ')到完全损坏(mode = 'DamagedSystem')的方式变化的系统中收集数据批。为了模拟这个场景，我们将模型置于“恶化系统”模式中。

set_param ([sysA,/机械系统的]，“LabelModeActiveChoice”，“DeterioratingSystem”）;sim(sysA) ytv = logsout.getElement(“y”) . values;ytv =挤压(ytv. data);PredictedMode = cell(nr,1);为ct = 1: nr zSegment = iddata (ytv ((ct-1) * 512 + (1:50 0)), [], Ts);mSegment = sest(zSegment, nx，“形式”，“规范”，“t”、Ts);PredictedMode(ct) = predict(cl, plist(esort(pole(mSegment))));结束I = strcmp(PredictedMode，‘好’）;Tags = ones(nr,1);Tags(~I) = -1;t = (0:5120)'*Ts;模拟时间%时间= t(1:512:end-1);情节(时间(I)、标签(我),“g *”、时间(~我),标签(~我),的r *，“MarkerSize”12)网格在轴([0 20 -2])“绿色:正常，红色:故障”)包含(“数据评估时间”) ylabel (“预测”）

图中显示，分类器预测直到中点的行为是正常的，此后处于故障状态。

基于模型参数在线自适应的故障检测

前面的分析使用了在系统运行期间不同时间收集的批量数据。监视系统运行状况的另一种通常更方便的方法是创建系统行为的自适应模型。新的测量值被连续处理，并被用于以递归的方式更新模型的参数。磨损或故障的影响通过模型参数值的变化来表示。

再次考虑损耗场景。随着系统的老化，会有更大的“咔啦咔啦”，这表现为几个共振模态的激发以及系统峰值响应的上升。这个场景在模型中描述pdmDeterioratingSystemEstimation的“变质系统”模式相同pdmMechanicalSystem除了为离线识别而添加的脉冲凸点不存在。系统的响应被传递给一个“递归多项式模型估计器”块，该块被配置为估计ARMA模型结构的参数。实际系统在正常状态下启动，但在200秒的时间跨度内会恶化到生命结束状态。

Initial_model = translatecov(@(x)idpoly(x)，模型);sysB =“pdmDeterioratingSystemEstimation”;open_system (sysB);

在转换为多项式(ARMA)格式后，使用前一节中导出的正常行为估计模型的参数和协方差数据对“ARMA模型”块进行了初始化。的translatecov ()函数的使用，使参数协方差数据也转换。该块使用“遗忘因子”算法，遗忘因子设置为略小于1，以在每个采样瞬间更新参数。遗忘因素的选择会影响系统更新的速度。一个小的值意味着更新将有很高的方差，而一个大的值将使估计器难以适应快速的变化。

利用模型参数估计来更新输出谱及其3-sd置信域。当频谱的置信区域在相关频率处不与健康系统的置信区域重叠时，系统将明显发生了变化。故障检测阈值在图中使用黑线表示在某些频率下允许的最大增益。随着系统中变化的积累，光谱会越过这条线。这可以作为故障的可视化指示器，在现实生活中的系统中可以使用它来要求修复。

运行模拟并观察频谱图的更新。

sim (sysB)

模型参数的运行估计也用于计算系统极点位置，然后将其输入支持向量机分类器，以预测系统是否处于“良好”或“故障”状态。这一决定也体现在剧情上。当预测归一化得分小于.3时，认为该决策是试探性的(接近区分边界)。见脚本pdmARMASpectrumPlot.m有关如何计算频谱和分类器预测的运行估计的详细信息。

可以在Simulink外部实现自适应估计和绘制过程万博1manbetxrecursiveARMA ()函数。“递归多项式模型估计器”块以及recursiveARMA ()函数支持用于部署目万博1manbetx的的代码生成。

该分类方案可以推广到有几种已知失效模式的情况。为此，我们需要多组分类器，其中模式指的是某种类型的故障。本例中没有讨论这些方面。

结论

这个示例展示了结合数据聚类和分类方法的系统识别方案如何帮助检测和隔离故障。讨论了顺序批处理分析和在线适应方案。确定了测量输出信号的ARMA结构模型。在既可以访问输入信号又可以访问输出信号的情况下，可以采用类似的方法，并希望采用其他类型的模型结构，如状态空间或Box-Jenkins多项式模型。

在这个例子中，我们发现: