CN103793601A - 基于异常搜索和组合预测的汽轮机组在线故障预警方法 - Google Patents

Abstract

本发明公开了属于电力系统的预警技术领域的一种基于异常搜索和组合预测的汽轮机组在线故障预警方法。该发明包括：负责标准化分割输入参数时间序列并提取序列特征模式的输入初始化处理；通过训练正常状态参数获取异常搜索参考标准的异常特征边界训练；通过搜索超越特征边界而确定异常序列组的异动搜索；利用回归分析识别异常变化趋势，得到异动分布变化规律的异常分析；建立预测模型对异常变化进行趋势预测；根据预测结果，结合异动参数与故障征兆的对应关系预警输出。本发明能够解决传统监测分析仅利用限值理论无法全面识别异常的缺陷，提高了异常预警的精度和深度，为机组故障起因及责任归属问题提供有利证据。

Description

基于异常搜索和组合预测的汽轮机组在线故障预警方法

技术领域

本发明涉及电力系统的预警技术领域，特别涉及一种基于异常搜索和组合预测的汽轮机组在线故障预警方法。具体说是火力发电厂、核能发电厂，需要对汽轮机组或者大型回转机械振动及过程信号进行识别、分析以及预警技术。

背景技术

在激烈的市场竞争中各发电企业谋求更好地发展，力求全面提升企业的综合竞争力，其中最重要的手段就是研发或引进具有国际先进化水平的运行设备。随着中大型企业运行设备结构和功能上的日趋复杂化，企业对设备正常、安全、稳定运行的要求越来越高，在保障设备安全运行方面，若出现故障再进行故障诊断往往为时已晚，因此如何保障大型旋转设备在负荷复杂多变的工况下的安全运行，避免严重的设备运行故障造成巨大经济损失成为企业研究的焦点。

当前，汽轮机组在运行及控制过程中振动及过程信号监测已实现在线采集，并通过简单比较监测参数的测量值与规范阈值来进行机组特性分析，并据此指导机组的运行控制。显然，与现代高精度数据采集技术相比，相对落后的监测分析技术已经严重滞后；传感器监测技术、振动分析诊断技术在相关领域内的应用基本实现了大型旋转机械的状态监测与故障诊断，但故障诊断及排除缺乏预见性，具有一定的功能滞后性，无法实现汽轮机组运行故障的早预警。同时，传统信号采集系统往往专注于对全局信号的采集分析，而缺乏对信号中隐藏的异动数据的二次搜索及挖掘分析。对汽轮机组这类大型高风险设备而言，搜索与异常及故障相关联的异动数据，并根据移动数据建立预测模型进行汽轮机故障预警显然更加具有研究价值。

发明内容

本发明的目的在于提供一种基于异常搜索和组合预测的汽轮机组在线故障预警方法，对机组监测的异常信号进行专业化搜索，并通过异动分析处理对异常进行深入挖掘分析，从而还原异常信号演绎形成的过程，解析其变化的趋势分布，最后通过故障预测模型实现汽轮机组故障预警，其特征在于，包括以下步骤：

1）对输入机组异动数据的参数初始化处理，分割输入参数的时间序列为标准时间子序列模式，分析计算子序列模式中的；

2）所述特征模式包括序列模式高度、模式长度、模式斜率、模式均值和标准差五个特征，并进行标准化；

3）异常特征边界训练，在上述输入参数初始处理基础上，从机组正常运行历史参数数据中，通过训练得到较信任的特征模式值的极大值，作为异常特征边界；利用移动窗口的初始分割，采用一种匹配各个序列点的阀门编码来识别分割子序列的边界点，其中阀门编码采用二进制编码，具有对时间序列极值对位编码进行变异的功能；

4）异动搜索，根据上述异常特征边界训练的边界结果，通过搜索超越边界的序列实现异动序列搜索，形成对应的异动序列数组；

5）异常分析，对上述搜索得到的异动序列数组进行基于确定系数优选的回归分析，识别异常数据的变化趋势，得到异常时间序列演绎下的异动分布变化规律；

6）预警输出，综合上述机组异动数据的发展过程的分析结果，输出异动演绎分布曲线，并结合各类监测参数应用汽轮机组故障预测模型进行故障预测，得到最大可能性故障模式的预警。

所述输入初始化处理过程包括分割输入参数的时间序列为标准时间子序列和计算子序列特征模式并标准化两阶段；其中，前一阶段目的是以期望滑动窗口和设定的阀门编码将预定时间序列划分为不重叠的标准时间子序列；后一阶段目的提供上述划分子序列的模式特征，并进行标准化；故障预警先采用GM(1,1)模型预测和对随机项进行AR(n)预测，最后将两种预测方法进行加权组合，形成灰色加权-AR(n)预测模型。

所述步骤2）各模式特征计算方法如下：

设时间序列X=(<x(t₁),t₁>,<x(t₂),t₂.,…,<x(t_n),t_n>)，第i个子时间序列为X_i=(<x(t_i1),t_i1>,<x(t_i2),t_i2>,…,<x(t_in),t_in>)，则其对应的特征模式定义如下：

模式高度sph    sph=x(t_in)-x(t_i1)

显然，当x(t_in)>x(t_i1)时，sph的值为正；当x(t_in)<x(t_i1)时，sph的值为负。

模式长度spl    spl=in-i1+1

模式斜率spk $\mathrm{spk}=\frac{x\left(t_{\mathrm{in}}\right)-x\left(t_{i1}\right)}{\mathrm{spl}}$

模式均值

$\mathrm{sp}\stackrel{\&OverBar;}{x}=\frac{1}{\mathrm{spl}}\underset{j=t_{i1}}{\overset{t_{\mathrm{in}}}{\mathrm{\&Sigma;}}}x\left(j\right)$

模式标准差spσ $\mathrm{sp\&sigma;}=\sqrt{\frac{1}{\mathrm{spl}}\underset{j=t_{\mathrm{il}}}{\overset{t_{\mathrm{in}}}{\mathrm{\&Sigma;}}}{(x\left(j\right)-\mathrm{sp}\stackrel{\&OverBar;}{x})}^2}$

其中，根据监测参数的变化规律，一般模式高度sph和模式长度spl模式单独使用仅用于具有周期变化的参数。

所述步骤2）将特征模式标准化，尽管这些特征值本身的值域不同，但它们在衡量时间序列的异常情况时具有同等重要的地位，因此有必要将其规范化，设c=c₁,c₂,…,c_n为其中一组特征值，则规范化定义为：

$\mathrm{norm}\left(c_i\right)=\frac{c_i-c_{\min }}{c_{\max }-c_i}$

其中，c_min为c=c₁,c₂,…,c_n中的最小值；c_max为c=c₁,c₂,…,c_n中的最大值，这样c中的每一个特征值就被规范到0-1之间。

所述步骤4）搜索的异动参数数据，按照异常演绎的可能规律进行回归分析拟合，有七种异常可能演绎趋势：持续恒定型、(陡升/陡降)跳跃型、线性(增长/下降)型、半岭形上升型、指数(增长/下降)型、浴盆曲线型、(凸型/凹型)抛物线型，利用最小二乘法估计出表达式的参数。

所述步骤4）中七种异常可能演绎趋势，首先假设

表示回归函数对应估计异常参数的拟合值，x(t)表示提取异常时间序列A_λ中的数据，则异常演绎的七类回归模型描述如下：

（1）持续恒定型

机组的异动搜索数据序列持续高于正常水平，且比较恒定，其时间序列模型为： $\hat{y}=x\left(t\right)=a_0$

其中a₀为常数且a₀≥0；

（2）跳跃型

机组的异动搜索数据突然从较低（高）水平跃至较高（低）水平，变化急速，称为陡升（陡降）分布，其时间序列模型为：

$\hat{y}x\left(t\right)=\left\{\begin{array}{c}{a}_0\&PlusMinus;e^t,t<t_0\\ a_t,t\&GreaterEqual;t_0\end{array}\right.,a_t>>a_0+e^t\mathrm{or}{a}_t<<a_0-e^t$

其中a₀、e、a_t为常数且a₀≥0，a₀±e^t>0，a_t>0；

（3）线性

机组的异动数据呈线性增长是机组振动长期运行常见的参数变化规律之一，振动异动呈线性增加（降低）变化是机组故障渐变过程的表现，其时间序列变化模型为：

$\hat{y}=x\left(t\right)=a_0\&PlusMinus;\mathrm{kt}$

其中a₀、k为常数且a₀≥0，k>0；

（4）半岭形

机组的异动数据随着时间持续上升到某一定值a₀后不再变化，称为半岭型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=\left\{\begin{array}{cc}0& t\&le;a_1\\ a_0\sin \frac{\mathrm{\&pi;}(t-a_1)}{2(a_2-a_1)}& a_1\&le;t\&le;a_2\\ a_0& t\&GreaterEqual;a_2\end{array}\right.$

其中a₀、a₁、a₂为常数a₀≥0，a₂≥a₁≥0；

（5）指数型

机组异动数据的非线性上升（下降）多呈指数型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=a\&PlusMinus;{\mathrm{be}}^t$

其中a、b为常数且a≥0，b≥0；

（6）浴盆型

机组异动数据的发展过程呈现由高到低，继而稳定一段时间，再由低到高，称为浴盆型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=\left\{\begin{array}{cc}{a}_0+a_{1}t+a_2{t}^2,& t\&le;t_1\mathrm{ort}\&GreaterEqual;t_2\\ a_3& t_1\&le;t\&le;t_2\end{array}\right.$

其中a₀、a₁、a₂、a₃为常数；

（7）抛物型

机组异动数据发展呈现由高（低）到低（高），再由低（高）到高（低）的过程，称为凹（凸）形抛物线型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=a_0+a_{1}t+a_2{t}^2$

其中a₀、a₁、a₂为常数且a₂>0；

各回归分析结果模型的确定系数R，其表达式为：

$R=\frac{\mathrm{\&Sigma;}{({\hat{x}}_i-\stackrel{\&OverBar;}{x})}^2}{\mathrm{\&Sigma;}{(x_i-\stackrel{\&OverBar;}{x})}^2},i=\mathrm{0,1,2}\&CenterDot;\&CenterDot;\&CenterDot;n$

其中x_i为当前异常序列数据，n为序列长度，

所述异常序列平均值，

为回归函数对应估计异常参数，选取确定系数最大的模型作为异动参数演绎分布形式。

所述步骤6）预警输出包括：在综合上述机组异动数据的发展过程的分析结果中，采用异动参数演绎分布形式，提供显示确定的异常分布趋势线图，触发给出异动等级作为预警输出技术方案；进行振动信号划分并分析振动参数的主要特征频率成分通过上述组合预测模型，预测得到特征频率各成分的变化趋势，由此预测特征频率进行故障预警。

所述预测特征频率进行故障预警是将汽轮机组常见的几种故障分为基频故障、广谱频故障、高频故障和低频故障，其中基频故障包括原始质量不平衡、转子部件脱落、转子热弯曲和支撑松动；高频故障包括不对中故障和轴裂纹故障；低频故障包括油膜失稳和蒸汽激振；动静碰磨属于广谱频故障；通过预测的特征频率与故障类型的对应确定故障模式，实现汽轮机组故障预警。

本发明的有益效果是本发明借鉴数据挖掘中的搜索技术思想，建立更加专业化的异常信号搜索模式，采用回归选优分析算法针对异常数据集进行演绎拟合及趋势分析，并配合异常演绎趋势结果及曲线显示进行详细预警显示，给出预警性指导；并采用灰色-AR(n)组和预测模型进行汽轮机组故障预警。本发明不但可以使得振动及过程监测数据得到更深化利用，并且大幅度提高振动及过程参数预警精度以及内容的深度，并为机组故障起因及责任归属问题提供更有利的判据线索。

附图说明

图1为本发明汽轮机组故障预警过程各步骤示意图。

图2为本发明汽轮机组振动及过程参数异常搜索分析方法中输入初始化处理过程的工作原理及流程图。

图3为本发明汽轮机组振动及过程参数异常搜索分析方法中异常边界训练过程的工作原理及流程图。

图4为本发明汽轮机组振动及过程参数异常搜索分析方法中异常分析过程的工作原理及流程图。

图5为本发明汽轮机组预警机制规则示意图。

图6为本发明汽轮发电机组典型故障振动特征频率属性识别模型流程图。

图7为本发明汽轮发电机组检测参数异常属性识别功能结构图。

具体实施方式

本发明提供一种基于异常搜索和组合预测的汽轮机组在线故障预警方法，下面结合附图予以说明。

图1所示为本发明汽轮机组故障预警过程各步骤示意图。

图中，基于汽轮机组振动及过程参数异常搜索分析的汽轮机组故障预警步骤，具体包括汽轮机组振动检测、常见故障故障模式积累与征兆因素异常分析的关联性分析、振动参数时间序列的初始化处理、输入初始处理、异常特征边界训练、异动搜索、异常分析、导入并训练预测模型、故障属性识别以及预警输出等步骤。

在上述汽轮机组振动及过程参数异常搜索分析方法中，输入初始处理技术方案包括标准输入时间子序列模式分割和提取子模式序列特征两阶段，前一阶段是以期望滑动窗口和规则将待检参数时间序列划分为不重叠的标准时间子序列模式；后一阶段包含序列模式高度、模式长度、模式斜率、模式均值和标准差五个特征计算程序。具体原理及流程参见图2。所述过程实施方法包括如下步骤：

步骤101：设置时间序列X=(<x(t₁),t₁>,<x(t₂),t₂>,…,<x(t_n),t_n>)的输入个数n，确定用于分割时间序列的检测窗口m。其中，对于时间序列输入个数n充分考虑算法以及参数数据采集周期，根据监测或预警期望感应时间进行设置，推荐选取1-5分钟的采集数为宜。对于检测窗口m的确定，通过先验数据的仿真实验表明：对于序列数据间变化比较平缓的数据集，m的值在9-11之间变化时都能够取得比较理想的效果；而对于序列数据间变化比较陡峭的数据集，则m的值在5-7之间变化，可根据现场实际情况据此设定。

步骤102：通过检测窗口内时间序列中的极值点来划分时间子序列。其中窗口序列极值点的统计通过每个序列点设置的阀门编码来识别，阀门编码定义如下：

阀门编码为一组对应时间序列点的二进制编码，初始值为零，即M_1×n=(0,0,…,0)，则规定若对应位上序列点为检测窗口内的单属性极值，则相应阀门值进行变异操作，即逻辑非操作，上述过程表达式如下：

$M\left(i\right)=\left\{\begin{array}{c}\stackrel{\&OverBar;}{M\left(i\right)}\\ M\left(i\right)\end{array}\right.<x\left(i\right),t_i>\&Element;<x_{\max },t>+<x\left(i\right),t_i>\&Element;<x_{\min },t>$

式中i--第i检测窗口；

〈x_max,t〉——在所有检测窗口中检测得到的最大序列值序列；

〈xmin,t〉——在所有检测窗口中检测得到的最小序列值序列。

在完成全部检测窗口阀门编码操作后，将阀门编码为1的序列点设为子序列边界点，时间序列的首尾两个序列点也被定义为边界点，因此时间序列的子序列就表现为从相邻的极小(大)点到极大(小)点之间的序列段。

步骤103：分别计算每一段时间子序列的五个特征模式，提供异常搜索引擎元素：模式高度sph、模式长度spl、模式斜率sps、模式均值

和模式标准差spσ，并将其标准化。其中各模式特征计算方法如下：

设时间序列X=(<x(t₁),t₁>,<x(t₂),t₂>,…,〈x(t_n),t_n〉)，第i个子时间序列为X_i=(<x(t_i1),t_i1>,<x(t_i2),t_i2>,…,<x(t_in),t_in>)，则其对应的特征模式定义如下：

模式高度sph    sph=x(t_in)-x(t_i1)

显然，当x(t_in)>x(t_i1)时，sph的值为正；当x(t_in)<x(t_i1)时，sph的值为负。

模式长度spl    spl=in-i1+1

模式斜率spk

模式均值

模式标准差spσ $\mathrm{sp\&sigma;}=\sqrt{\frac{1}{\mathrm{spl}}\underset{j=t_{\mathrm{il}}}{\overset{t_{\mathrm{in}}}{\mathrm{\&Sigma;}}}{(x\left(j\right)-\mathrm{sp}\stackrel{\&OverBar;}{x})}^2}$

其中，根据监测参数的变化规律，一般模式高度sph和模式长度spl模式单独使用仅用于具有周期变化的参数。

步骤104：将特征模式标准化。尽管这些特征值本身的值域不同，但它们在衡量时间序列的异常情况时具有同等重要的地位，因此有必要将其规范化。设c=<c₁,c₂,…,c_n>为其中一组特征值，则规范化定义为：

$\mathrm{norm}\left(c_i\right)=\frac{c_i-c_{\min }}{c_{\max }-c_i}$

其中，c_min为c=c₁,c₂,…,c_n中的最小值；c_max为c=c₁,c₂,…,c_n中的最大值，这样c中的每一个特征值就被规范到0-1之间。

在上述汽轮机组振动及过程参数异常搜索分析方法中，异常特征边界训练技术方案在上述输入初始处理基础上，以机组正常状态的历史参数时间序列为训练对象，将子序列集中的各特征模式的极大值通过训练由大到小排序，以期望最低出现率为指标选取期望数量极值的均值作为其特征边界。其具体原理及计算流程见图3所示。

步骤105：确定训练空间的大小，给出历史参数时间序列组数N和序列维数k。对于拥有丰富历史数据的现场而言，训练的组数越多，得到的异常特征边界越接近实际，一般选取N≥200为佳；序列维数k可仿照上述输入初始处理步骤101的方法确定。

步骤106：按照上述输入初始处理步骤102的方法，对N组历史参数时间序列经初始化处理，形成N组包含5类特征模式值的向量集：

步骤107：依次从N组特征模式向量集中抽取各类特征模式的极大值，即

并按从大到小排列，形成特征模式极值排序序列 $Z_i=(z_i^1,z_i^2,\&CenterDot;\&CenterDot;\&CenterDot;,z_i^j,\&CenterDot;\&CenterDot;\&CenterDot;z_i^N),z_i^j>z_i^{j+1},$ 其中i代表特征模式类别。

步骤108：设定异常边界的期望最低出现率ε，在Z_i中统计满足ε约束的最小极大值数目p，即求解如下不等式：

其中N为Z_i中的极大值序列长度。

将Z_i中的

作为第i特征模式的异常特征边界，也可根据需要将前p位极大值的平均值

作为结果。ε一般根据训练空间大小N进行取值，对于N=200的训练空间而言，ε取5%～10%为佳。

在上述汽轮机组振动及过程参数异常搜索分析方法中，异动搜索技术方案是在上述输入初始处理基础上，采用上述异常特征边界训练的边界结论，通过搜索超越边界的序列实现异动序列搜索。一般规定子序列中，任意特征模式超过异常边界即视为异常序列，并存储到异常序列组A中，以便异常分析使用，上述搜索过程即可表达为：

$A=\left\{A_i\right|\&ForAll;{\mathrm{sp}}_j>z_j,j\&Element;\left(\mathrm{1,2,3,4,5}\right),0\&le;i\&le;n\}$

其中A_i为组成异常序列组的异常子序列；

sp_j为子序列的特征模式值；

z_j为特征模式sp_j对应的异常特征边界；

n为输入时间序列X在输入初始化处理过程中分割的子序列数。

在上述汽轮机组振动及过程参数异常搜索分析方法中，异常分析技术方案采用基于确定系数分析的异常发展趋势回归选优算法，分析拟合出异常序列组中异常信号的演绎分布规律，提供分布趋势图，具体原理及流程参见图4。所述过程实施方法包括如下步骤：

步骤109：将搜索的异动参数数据，按照异常演绎的可能规律进行回归分析拟合，并利用最小二乘法估计出表达式的参数。本发明提供七种异常可能演绎趋势：持续恒定型、(陡升/陡降)跳跃型、线性(增长/下降)型、半岭形上升型、指数(增长/下降)型、浴盆曲线型、(凸型/凹型)抛物线型。

假设

表示回归函数对应估计异常参数的拟合值，x(t)表示提取异常时间序列Aλ中的数据，则上述方法中提供的异常演绎的七类回归模型可以描述为：

（8）持续恒定型

机组的异动搜索数据序列持续高于正常水平，且比较恒定，其时间序

列模型为： $\hat{y}=x\left(t\right)=a_0$

其中a₀为常数且a₀≥0；

（9）跳跃型

机组的异动搜索数据突然从较低（高）水平跃至较高（低）水平，变化急速，称为陡升（陡降）分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=\left\{\begin{array}{c}{a}_0\&PlusMinus;e^t,t<t_0\\ a_t,t\&GreaterEqual;t_0\end{array}\right.,a_t>>a_0+e^t\mathrm{or}{a}_t<<a_0-e^t$

其中a₀、e、a_t为常数且a₀≥0，a₀±e^t>0，a_t>0；

（10）  线性

机组的异动数据呈线性增长是机组振动长期运行常见的参数变化规律之一。振动异动呈线性增加或降低的变化是机组故障渐变过程的表现，其时间序列变化模型为：

$\hat{y}=x\left(t\right)=a_0\&PlusMinus;\mathrm{kt}$

其中a₀、k为常数且a₀≥0，k>0；

（11）  半岭形

机组的异动数据随着时间持续上升到某一定值a₀后不再变化，称为半岭型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=\left\{\begin{array}{cc}0& t\&le;a_1\\ a_0\sin \frac{\mathrm{\&pi;}(t-a_1)}{2(a_2-a_1)}& a_1\&le;t\&le;a_2\\ a_0& t\&GreaterEqual;a_2\end{array}\right.$

其中a₀、a₁、a₂为常数a₀≥0，a₂≥a₁≥0；

（12）  指数型

机组异动数据的非线性上升或下降，多呈指数型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=a\&PlusMinus;{\mathrm{be}}^t$

其中a、b为常数且a≥0，b≥0；

（13）浴盆型

机组异动数据的发展过程呈现由高到低，继而稳定一段时间，再由低到高，称为浴盆型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=\left\{\begin{array}{cc}{a}_0+a_{1}t+a_2{t}^2,& t\&le;t_1\mathrm{ort}\&GreaterEqual;t_2\\ a_3& t_1\&le;t\&le;t_2\end{array}\right.,$

其中a₀、a₁、a₂、a₃为常数；

（14）抛物型

机组异动数据发展呈现由高（低）到低（高），再由低（高）到高（低）的过程，称为凹（凸）形抛物线型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=a_0+a_{1}t+a_2{t}^2$

其中a₀、a₁、a₂为常数且a₂>0。

步骤110：计算109步骤中各回归分析结果模型的确定系数R，其表达式为：

$R=\frac{\mathrm{\&Sigma;}{({\hat{x}}_i-\stackrel{\&OverBar;}{x})}^2}{\mathrm{\&Sigma;}{(x_i-\stackrel{-}{x})}^2}i=\mathrm{0,1,2}\&CenterDot;\&CenterDot;\&CenterDot;n$

其中x_i为当前异常序列数据，n为序列长度，

所述异常序列平均值，

为回归函数对应估计异常参数。

步骤111：选取确定系数最大的模型作为异动参数演绎分布形式。

在上述汽轮机组振动及过程参数异常搜索分析方法中，预警输出技术方案采用与上述异常分析结果对应的预警机制，其中所述预警机制规则参见图5。提供显示确定的异常分布趋势线图，触发给出异动等级。

步骤112：进行振动信号划分并分析振动参数的主要特征频率成分。以连续轴振动采集信号的时间子序列的整数倍N作为分析单元，，取N=1min。对各分析单元进行快速傅里叶变换（FFT），求取此区间振动信号的各频率成分，一般情况下，用成分比例来代替各频率成分分布。

本过程通过快速傅里叶变换（FFT）对输入异常分析单元的时间序列进行频域处理，并根据分辨率求取工频(X)、2倍频(2X)、高倍频（≥3倍频，nX)、1/2倍频(1/2X)、低倍频（<1/2倍频，1/nX）对应的最大振幅Ai，i=1，2，3，4，5。设定分析信号序列为X={x(1),x(2)),...,x(N)}，且满足N=2^m点，对其进行FFT变换，即Y=FFT(X)。具体流程参见图6。

对于偶数序号输出采用基2算法，表示为

$x\left(2r\right)=\underset{n=0}{\overset{\frac{N}{2}-1}{\mathrm{\&Sigma;}}}[x\left(n\right)+x(n+\frac{N}{2})]W_{N/2}^{\mathrm{nr}}$

式中r=0,1,…,(N/2-1)。

对于奇数序号输出采用基4算法，表示为

$x(4r+1)=\underset{n=0}{\overset{\frac{N}{4}-1}{\mathrm{\&Sigma;}}}\left\{\right[x\left(n\right)-x(n+\frac{N}{2})]-j[x(n+\frac{N}{4})-x(n+\frac{3N}{4})\left]\right\}{W}_N^n{W}_{N/2}^{\mathrm{nr}}$

$x(4r+3)=\underset{n=0}{\overset{\frac{N}{4}-1}{\mathrm{\&Sigma;}}}\left\{\right[x\left(n\right)-x(n+\frac{N}{2})]-j[x(n+\frac{N}{4})-x(n+\frac{3N}{4})\left]\right\}{W}_N^{3n}{W}_{N/2}^{\mathrm{nr}}$

式中r=0,1,…,(N/4-1)。

在上述变换基础上，设分析单元内振动采集信号分辨率为Δf，信号频率f，则对应的1X、2X、nX、1/2X、1/nX对应的时域采样点编号为：

k(1X)=f/Δf；k(2X)=2f/Δf；{k(nX)}=(3:5)f/Δf；

k(1/2X)=f/2Δf；{k(1/nX)}=(1/3:1/5)f/2Δf。

对应1X、2X、nX、1/2X、1/nX的振幅为：

A1=2*|Y(k(1X))|；A2=2*|Y(k(2X))|；A3=Σ2*|Y({k(nX)})|；

A4=2*|Y(k(1/2X))|；A5=Σ2*|Y({k(1/nX)})|。

则频率成分可用成分比率表示：

$\mathrm{zj}=\left(A_j\right)/\underset{i=1}{\overset{5}{\mathrm{\&Sigma;}}}{A}_i,j=\mathrm{1,2,3,4,5}$

对K个分析单元分别重复上述操作，则得到各频率成分比率集合：

Z=(Z₁,Z₂,Z₃,Z₄,Z₅)=({z₁(1),z₁(2),…,z₁(K)},…,{z₅(1),z₅(2),…,z₅(K)})

步骤113：进行级比检验，使振动频率成分比例序列满足预测模型要求：原始时间序列

计算级比σ(i)=x(i-1)/x(i)，i=2,…,N；求区间

判断σ(i)是否落于区间

内，若大多数点落入该区域内，可建模。否则，对原始序列利用数据平滑方法进行数据预处理，具体方法如下：

x⁽⁰⁾(i)=[x(i-1)+2x(i)+x(i+1)]4

两端点数据处理方法如下：

对新得到的数列再次进行级比检验，判断是否可以建模。

步骤114：进行振动频率成分比例趋势项的预测。进行GM(1,1)模型预测，其中包括数据序列的生成、建立灰色模型、精度检验、趋势项预测四部分内容：

(1)数据序列的生成。设原始序列为经过AGO-1后为 $X^{\left(1\right)}={\left\{x^{\left(1\right)}\left(i\right)\right\}}_{i=1}^N,$ 其中 $x^{\left(1\right)}\left(i\right)=\underset{j=1}{\overset{i}{\mathrm{\&Sigma;}}}{x}^{\left(0\right)}\left(j\right),i=2,\&CenterDot;\&CenterDot;\&CenterDot;,N\&CenterDot;$

(2)建立灰色模型GM(1,1)并求解参数。GM(1,1)模型的灰微分方程为：x⁽⁰⁾(i)+aZ⁽¹⁾(i)=b，其对应的白化方程为：

其中参数

由下式确定：

$\hat{a}=\left[\begin{array}{c}a\\ b\end{array}\right]={\left(B^{T}B\right)}^{-1}{B}^T\mathrm{\&omega;}{Y}_N$

其中B、Y_N和ω分别为：

$B=\left[\begin{array}{cc}-\frac{1}{2}[x^{\left(1\right)}\left(1\right)+x^{\left(1\right)}\left(2\right)]& 1\\ -\frac{1}{2}[x^{\left(1\right)}\left(2\right)+x^{\left(1\right)}\left(3\right)]& 1\\ \&CenterDot;& \&CenterDot;\\ \&CenterDot;& \&CenterDot;\\ \&CenterDot;& \&CenterDot;\\ -\frac{1}{2}[x^{\left(1\right)}(N-1)+x^{\left(1\right)}\left(N\right)]& 1\end{array}\right]$

$Y_N=\left[\begin{array}{c}{X}^{\left(0\right)}\left(2\right)\\ X^{\left(0\right)}\left(3\right)\\ \&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\\ X^{\left(0\right)}\left(N\right)\end{array}\right]$

其中ω_i=R^i-1，i=1,2,…,n，R=1.5

得到GM(1,1)模型的白化响应式为：

${\hat{x}}^{\left(1\right)}(i+1)=(x^{\left(0\right)}\left(1\right)-\frac{b}{a})e^{-\mathrm{ai}}+\frac{b}{a}$

还原得原序列值：

${\hat{x}}^{\left(0\right)}(i+1)={\hat{x}}^{\left(1\right)}(i+1)-{\hat{x}}^{\left(i\right)}=(1-e^a)(x^{\left(0\right)}\left(1\right)-\frac{b}{a})e^{-\mathrm{ai}},i=\mathrm{1,2},\&CenterDot;\&CenterDot;\&CenterDot;,N$

(3)精度检验。残差 $\left(i\right)=x^{\left(0\right)}\left(i\right)-{\hat{x}}^{\left(0\right)}\left(i\right),i=\mathrm{1,2},\&CenterDot;\&CenterDot;\&CenterDot;,N,$ 残差向量e=(e(1),e(2),…e(N))，记原始序列X⁽⁰⁾以及残差序列e的方差分别为

和

则有：

${\stackrel{\&OverBar;}{X}}^{\left(0\right)}=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{x}^{\left(0\right)}\left(i\right),{\mathrm{\&sigma;}}_1^2=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}(x^{\left(0\right)}\left(i\right)-{\stackrel{\&OverBar;}{x}}^{\left(0\right)}{)}^2$

$\stackrel{\&OverBar;}{e}=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}e\left(i\right),{\mathrm{\&sigma;}}_2^2=\frac{1}{N}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{(e\left(i\right)-\stackrel{\&OverBar;}{e})}^2$

计算后验比值C=σ₂/σ₁和小误差概率

若C，p均符合要求，则可按照操作(2)进行趋势项预测，否则进行残差修正，直到符合要求为止。

(4)趋势项预测。利用步骤114中的公式

${\hat{x}}^{\left(0\right)}(i+1)={\hat{x}}^{\left(1\right)}(i+1)-{\hat{x}}^{\left(1\right)}\left(i\right)=(1-e^a)(x^{\left(0\right)}\left(1\right)-\frac{b}{a})e^{-\mathrm{ai}}$

计算

模型初始序列

变为：

$X^{\left(0\right)}=\{{\left\{x^{\left(0\right)}\left(i\right)\right\}}_{i=2}^N,{\hat{x}}^{\left(0\right)}(N+1)\}$

转至(1)数列的生成，对新序列求AGO并求参数a、b，计算

同理依次求预测值 ${\hat{X}}^{\left(0\right)}(N+3),={\hat{x}}^{\left(0\right)}(N+4),...,+{\hat{x}}^{\left(0\right)}(N+l).$

步骤115：对随机项进行AR(n)预测。其中包括平稳性及零均值性检验、计算自相关系数、计算各阶AICC数值、AR模型预测四部分内容。

(1)平稳性、零均值性检验。对GM(1,1)模型提取趋势项所产生的残差e，首先进行平稳性、零均值性检验。

1)计算残差序列e的逆总序数。对于e(i)(i=1,2,…,N-1)来说，若后面有一个值大于它，也就是说e_j>e_i(j>i)，则称e_j为一个逆序，并记逆序数为A_i，逆序总数为：

2)统计检验

$A=\underset{i=1}{\overset{N-1}{\mathrm{\&Sigma;}}}{A}_i$

$E\left(A\right)=\frac{1}{4}N(N-1)$

$D\left(A\right)=\frac{N({2N}^2+3N-5)}{72}$

统计量 $Z=\frac{[A+\frac{1}{2}-E\left(A\right)]}{\sqrt{D\left(A\right)}}$

在显著性水平α=0.05情况下，若|Z|<1.96，则认为残差序列是平稳的，否则认为序列是非平稳的。对非平稳序列一般用差分法进行平稳化。

ΔX_t=X_t-X_t-1

进行d次差分形式

Δ^dX_t=Δ(Δ^d-1X_t-Δ^d-1X_t-1)=…

3)序列平均值化。将满足统计检验的序列数据减去其均值，然后用得到的零均值序列建模。

(2)计算自相关系数。

${\widehat{\mathrm{\&rho;}}}_i=\frac{{\widehat{\mathrm{\&gamma;}}}_i}{{\widehat{\mathrm{\&gamma;}}}_0}=\frac{\frac{1}{N}\underset{t=i+1}{\overset{N}{\mathrm{\&Sigma;}}}{X}_t{X}_{t-i}}{\frac{1}{N}\underset{t=1}{\overset{N}{\mathrm{\&Sigma;}}}{X}_t^2}=\frac{\underset{t=i+1}{\overset{N}{\mathrm{\&Sigma;}}}{X}_t{X}_{t-i}}{\underset{t=1}{\overset{N}{\mathrm{\&Sigma;}}}{X}_t^2}$

(3)计算各阶AICC数值，利用AICC准则判定AR模型的阶数。计算偏自相关系数

(4)AR模型预测。判断AR模型是否稳定，若模型稳定，用该模型对残差序列进行一步预测，加上平均化处理的均值，如需要再进行差分反运算，得到预测值，可以应用等维灰数递补模型继续进行预测。

步骤116：实现组合预测。将残差序列的时序AR预测结果加到灰色趋势项预测值上，得到最终的灰色加权-AR(n)预测结果。

灰色-时序组合预测：将残差序列的时序AR预测结果加到灰色趋势项预测值上，得到最终的灰色加权-AR(n)预测模型式子为：

${\hat{X}}^{\left(0\right)}(N+l)={\hat{x}}^{\left(0\right)}(N+l)+{\hat{e}}^{\left(0\right)}(N+l)$

通过上述组合预测模型，可以预测得到特征频率各成分的变化趋势。

步骤117：根据对特征频率的预测进行故障预警。将汽轮机组常见的几种故障分为基频故障、广谱频故障、高频故障和低频故障，其中基频故障包括原始质量不平衡、转子部件脱落、转子热弯曲和支撑松动；高频故障包括不对中故障和轴裂纹故障；低频故障包括油膜失稳和蒸汽激振；动静碰磨属于广谱频故障。通过预测的特征频率与故障类型的对应确定故障模式，实现汽轮机组故障预警。

本发明描述的汽轮机组振动异常搜索分析和故障预警方法可针对机组振动及任意过程参数信号的分析，利用专业搜索方法对其进行异动搜索，并通过异常分析对异常信号进行回归分析，最终通过故障预测模型实现汽轮机组故障预警。本发明提高了电厂对汽轮机组异常信号的利用和挖掘，丰富了预警内容，为保证电厂汽轮机组的安全稳定运行提供了有力的依据。

对于本发明各个步骤中所阐述的方法，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围内。

Claims (8)

1.一种基于异常搜索和组合预测的汽轮机组在线故障预警方法，对机组监测的异常信号进行专业化搜索，并通过异动分析处理对异常进行深入挖掘分析，从而还原异常信号演绎形成的过程，解析其变化的趋势分布，最后通过故障预测模型实现汽轮机组故障预警，其特征在于，包括以下步骤：

1）对输入机组异动数据的参数初始化处理，分割输入参数的时间序列为标准时间子序列模式，分析计算子序列模式中的；

2）所述特征模式包括序列模式高度、模式长度、模式斜率、模式均值和标准差五个特征，并进行标准化；

3）异常特征边界训练，在上述输入参数初始处理基础上，从机组正常运行历史参数数据中，通过训练得到较信任的特征模式值的极大值，作为异常特征边界；利用移动窗口的初始分割，采用一种匹配各个序列点的阀门编码来识别分割子序列的边界点，其中阀门编码采用二进制编码，具有对时间序列极值对位编码进行变异的功能；

4）异动搜索，根据上述异常特征边界训练的边界结果，通过搜索超越边界的序列实现异动序列搜索，形成对应的异动序列数组；

5）异常分析，对上述搜索得到的异动序列数组进行基于确定系数优选的回归分析，识别异常数据的变化趋势，得到异常时间序列演绎下的异动分布变化规律；

6）预警输出，综合上述机组异动数据的发展过程的分析结果，输出异动演绎分布曲线，并结合各类监测参数应用汽轮机组故障预测模型进行故障预测，得到最大可能性故障模式的预警。

2.根据权利要求1所述一种基于异常搜索和组合预测的汽轮机组在线故障预警方法，其特征在于，所述输入初始化处理过程包括分割输入参数的时间序列为标准时间子序列和计算子序列特征模式并标准化两阶段；其中，前一阶段目的是以期望滑动窗口和设定的阀门编码将预定时间序列划分为不重叠的标准时间子序列；后一阶段目的提供上述划分子序列的模式特征，并进行标准化；故障预警先采用GM(1,1)模型预测和对随机项进行AR(n)预测，最后将两种预测方法进行加权组合，形成灰色加权-AR(n)预测模型。

3.根据权利要求1所述一种基于异常搜索和组合预测的汽轮机组在线故障预警方法，其特征在于，所述步骤2）各模式特征计算方法如下：

设时间序列X=(<x(t₁),t₁>,<x(t₂),t₂>,…,<x(t_n),t_n>)，第i个子时间序列为X_i=(<x(t_i1),t_i1>,<x(t_i2),t_i2>,…,<x(t_in),t_in>)，则其对应的特征模式定义如下：

模式高度sph：    sph=x(t_in)-x(t_i1)，

显然，当x(t_in)>x(t_i1)时，sph的值为正；当x(t_in)<x(t_i1)时，sph的值为负，

模式长度spl：    spl=in-i1+1，

模式斜率spk： $\mathrm{spk}=\frac{x\left(t_{\mathrm{in}}\right)-x\left(t_{i1}\right)}{\mathrm{spl}},$

模式均值

$\mathrm{sp}\stackrel{\&OverBar;}{x}=\frac{1}{\mathrm{spl}}\underset{j=t_{i1}}{\overset{t_{\mathrm{in}}}{\mathrm{\&Sigma;}}}x\left(j\right),$

模式标准差spσ： $\mathrm{sp\&sigma;}=\sqrt{\frac{1}{\mathrm{spl}}\underset{j=t_{\mathrm{il}}}{\overset{t_{\mathrm{in}}}{\mathrm{\&Sigma;}}}{(x\left(j\right)-\mathrm{sp}\stackrel{\&OverBar;}{x})}^2},$

其中，根据监测参数的变化规律，一般模式高度sph和模式长度spl模式单独使用仅用于具有周期变化的参数。

4.根据权利要求1所述一种基于异常搜索和组合预测的汽轮机组在线故障预警方法，其特征在于，所述步骤2）将特征模式标准化，尽管这些特征值本身的值域不同，但它们在衡量时间序列的异常情况时具有同等重要的地位，因此有必要将其规范化，设c=c₁,c₂,…,c_n为其中一组特征值，则规范化定义为：

$\mathrm{norm}\left(c_i\right)=\frac{c_i-c_{\min }}{c_{\max }-c_i},$

其中，c_min为c=c₁,c₂,…,c_n中的最小值；c_max为c=c₁,c₂,…,c_n中的最大值，这样c中的每一个特征值就被规范到0-1之间。

5.根据权利要求1所述一种基于异常搜索和组合预测的汽轮机组在线故障预警方法，其特征在于，所述步骤4）搜索的异动参数数据，按照异常演绎的可能规律进行回归分析拟合，有七种异常可能演绎趋势：持续陡升或陡降恒定型、跳跃型、线性增长或线性下降型、半岭形上升型、指数增长或下降型、浴盆曲线型、凸型或凹型抛物线型，利用最小二乘法估计出表达式的参数。

6.根据权利要求1所述一种基于异常搜索和组合预测的汽轮机组在线故障预警方法，其特征在于，所述步骤4）中七种异常可能演绎趋势，首先假设

表示回归函数对应估计异常参数的拟合值，x(t)表示提取异常时间序列Aλ中的数据，则异常演绎的七类回归模型描述如下：

（1）持续恒定型：

机组的异动搜索数据序列持续高于正常水平，且比较恒定，其时间序列模型为： $\hat{y}=x\left(t\right)=a_0,$

其中a₀为常数且a₀≥0；

（2）跳跃型：

机组的异动搜索数据突然从较低或较高水平跃至较高或较低水平，变化急速，称为陡升或陡降分布，其时间序列模型为：

$\hat{y}x\left(t\right)=\left\{\begin{array}{c}{a}_0\&PlusMinus;e^t,t<t_0\\ a_t,t\&GreaterEqual;t_0\end{array}\right.,a_t>>a_0+e^t\mathrm{or}{a}_t<<a_0-e^t$

其中a₀、e、a_t为常数且a₀≥0，a₀±e^t>0，a_t>0；

（3）线性：

机组的异动数据呈线性增长是机组振动长期运行常见的参数变化规律之一，振动异动呈线性增加（降低）变化是机组故障渐变过程的表现，其时间序列变化模型为：

$\hat{y}=x\left(t\right)=a_0\&PlusMinus;\mathrm{kt},$

其中a₀、k为常数且a₀≥0，k>0；

（4）半岭形：

机组的异动数据随着时间持续上升到某一定值a₀后不再变化，称为半岭型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=\left\{\begin{array}{cc}0& t\&le;a_1\\ a_0\sin \frac{\mathrm{\&pi;}(t-a_1)}{2(a_2-a_1)}& a_1\&le;t\&le;a_2\\ a_0& t\&GreaterEqual;a_2\end{array}\right.,$

其中a₀、a₁、a₂为常数a₀≥0，a₂≥a₁≥0；

（5）指数型：

机组异动数据的非线性上升或下降多呈指数型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=a\&PlusMinus;{\mathrm{be}}^t,$

其中a、b为常数且a≥0，b≥0；

（6）浴盆型：

机组异动数据的发展过程呈现由高到低，继而稳定一段时间，再由低到高，称为浴盆型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=\left\{\begin{array}{cc}{a}_0+a_{1}t+a_2{t}^2,& t\&le;t_1\mathrm{ort}\&GreaterEqual;t_2\\ a_3& t_1\&le;t\&le;t_2\end{array}\right.,$

其中a₀、a₁、a₂、a₃为常数；

（7）抛物型：

机组异动数据发展呈现由高低或由低到高，再由低到高或由到高到低的过程，称为凹或凸形抛物线型分布，其时间序列模型为：

$\hat{y}=x\left(t\right)=a_0+a_{1}t+a_2{t}^2,$

其中a₀、a₁、a₂为常数且a₂>0；

各回归分析结果模型的确定系数R，其表达式为：

$R=\frac{\mathrm{\&Sigma;}{({\hat{x}}_i-\stackrel{\&OverBar;}{x})}^2}{\mathrm{\&Sigma;}{(x_i-\stackrel{\&OverBar;}{x})}^2},i=\mathrm{0,1,2}\&CenterDot;\&CenterDot;\&CenterDot;n,$

其中x_i为当前异常序列数据，n为序列长度，

所述异常序列平均值，为回归函数对应估计异常参数，选取确定系数最大的模型作为异动参数演绎分布形式。

7.根据权利要求1所述一种基于异常搜索和组合预测的汽轮机组在线故障预警方法，其特征在于，所述步骤6）预警输出包括：在综合上述机组异动数据的发展过程的分析结果中，采用异动参数演绎分布形式，提供显示确定的异常分布趋势线图，触发给出异动等级作为预警输出技术方案；进行振动信号划分并分析振动参数的主要特征频率成分通过上述组合预测模型，预测得到特征频率各成分的变化趋势，由此预测特征频率进行故障预警。

8.根据权利要求7所述一种基于异常搜索和组合预测的汽轮机组在线故障预警方法，其特征在于，所述预测特征频率进行故障预警是将汽轮机组常见的几种故障分为基频故障、广谱频故障、高频故障和低频故障，其中基频故障包括原始质量不平衡、转子部件脱落、转子热弯曲和支撑松动；高频故障包括不对中故障和轴裂纹故障；低频故障包括油膜失稳和蒸汽激振；动静碰磨属于广谱频故障；通过预测的特征频率与故障类型的对应确定故障模式，实现汽轮机组故障预警。

Images