基于振动信号的水电机组状态劣化在线评估方法研究

刘 东，赖 旭，胡 晓，肖志怀

（1.武汉大学 水资源与水电工程科学国家重点实验室，湖北 武汉 430072；2.武汉大学 水力机械过渡过程教育部重点实验室，湖北 武汉 430072）

1 研究背景

水电机组作为水电能源转换的关键设备，其安全性和稳定性一直是电力行业关注的热点。随着机组不断向大型化、复杂化发展，集成化程度越来越高，结构也日趋复杂。因此，为确保水电机组安全稳定运行，提高设备利用率，避免重大经济损失和人员伤亡，必须对水电机组运行状态进行合理的监测与评估。

目前对水电机组进行状态监测和健康评估的主要方法是针对机组振动信号进行分析［1-2］，有效的特征提取方法是信号分析的常用手段。近年来，学者们研究了多种信号处理方法用于水电机组振动信号的特征提取，如经验模态分解［3］，小波奇异值分解［4］等，通过对提取出的特征向量进行分类，可以准确判断振动信号所反映的机组运行状态是正常还是故障。而水电机组的故障发生为一个渐变过程，特征提取与分类相结合的故障诊断方法侧重于对机组故障发生时的判断，难以体现机组从正常转变为故障的演化过程，无法为机组提供主动有效的故障预警保护。随着传感器技术、计算机技术的不断进步，水电行业开展了大量状态监测方面的工作，建立了比较完善的机组实时状态监测系统，为实现机组健康评估和故障预警提供了有效的技术支撑。现有的故障预警方法大多通过振动监测值与振动阈值的比较来实现，而振动阈值通常根据相关国家标准［5］、行业标准［6］，以及厂家建议值等来确定，这种单一限值比较法未能充分反映机组自身特点，忽略了机组工况对振动的影响，且没有对状态监测系统产生的海量历史数据进行有效利用，其预警可靠度有待提高。

针对目前水电机组状态监测与健康评估存在的不足，部分学者开展了考虑机组特性和运行工况的故障诊断方法研究［7-8］，取得了一定成效。本文提出一种融合振动信号时域和频域特征的水电机组劣化趋势在线评估方法。首先，以机组健康状态下工况X（水头、开度等）和振动信号时域特征Y为健康样本，利用最小二乘支持向量机建立机组状态健康模型Y=f(X)。将机组状态监测系统中的实时工况数据代入健康模型中，预测该工况下振动信号时域特征健康值，计算信号时域特征实际值与健康值之间的相对误差，作为代表时域方向的劣化评估指标。在信号时域特征选取方面，使用数理统计中的检测指数选取对故障最为敏感的时域特征，避免特征选取的盲目性。其次，结合小波变换与奇异值分解理论，对机组健康状态振动信号进行特征提取，得到包含信号不同频率成分的奇异值特征向量，利用数据挖掘中的K 均值聚类算法得到信号特征向量的健康聚类中心，计算实测信号特征向量与健康聚类中心之间的相对欧式距离，作为当前频域方向的劣化指标。进而将时域劣化指标与频域劣化指标结合，构建综合劣化指标，实时评估机组劣化程度。文中以实际电站机组运行案例对该方法反映机组劣化趋势的有效性进行了检验。

2 理论基础

2.1 检测指数检测指数（Detection Index，DI）是数理统计中的概念［9］，本文将其重新定义，用于水电机组振动信号时域特征参数的选取。假设x1、x2分别为状态1 和状态2 下所采集信号的同一种时域特征值，并且分别符合正态分布N(μ1，σ1)、N(μ2，σ2)，μ和σ分别为特征参数的均值和方差，结合概率统计理论可以认为：|x2-x1|越大，那么该特征参数对于不同状态的区分度就越强，即具有更高的敏感性。大量研究表明，水电机组振动信号的统计特征有着明显的正态分布规律［8-10］，满足检测指数的适用条件。

令z=x2-x1，且μ2＞μ1，则z同样服从正态分布N(μ2-μ1，σ1+σ2)，其概率密度函数为：

定义识别率（Discrimination Rate，DR）为同一种时域特征值对水电机组不同运行状态的区分度，计算方式为：

令：，利用换元积分法，结合式（1）、式（2）可得：

式中，检测指数DI表征了机组时域特征参数对机组不同运行状态的敏感性，计算方式为：

或者

式中，、为x1、x2的均值。很明显，由式（3）可知，DI值越大，则识别率DR越大，那么特征参数对于不同状态的敏感性就越强，同理当μ1＞μ2时，也会得到同样的结果。本文正式选用的DI值为：

2.2 最小二乘支持向量机理论最小二乘支持向量机（Least Squares Support Vector Machine，LS-SVM）是支持向量机的变体与改进，它结合最小二乘法理论，将支持向量机分类与回归中的不等式约束条件转换为等式约束，并用误差平方和损失函数代替其不敏感损失函数，因此，支持向量机的最优解问题也由二次规划问题转变为求解线性方程组，降低了计算的复杂度，具有良好的泛化能力［11-12］。

给定包含N组样本的训练集，其中xk、yk分别为第k组输入、输出样本，且xk∈R n，yk∈R。则支持向量机解决回归分类问题的方法是构造如下形式的数学模型：

式中：ψ·(，·)为映射函数，将低维输入向量映射到高维空间；αk、b均为常数。

最小二乘支持向量机通过引入最小二乘法，其优化问题为：

式中：ω 为权向量矩阵；e为误差，称为松弛因子；γ为正则化参数，又称惩罚因子。该优化问题的约束条件为：

式中，φ(xk)为非线性映射函数。

为了将带约束条件的优化问题转化为无约束条件的优化问题，引入拉格朗日方程，并根据Ka⁃rush-Kuhn-Tucker 条件得到最小二乘支持向量机的线性回归方程［13］：

因此，式（5）中支持向量机的回归与分类问题转变为求解线性方程，参数α与b通过最小二乘法求出，映射函数ψ·(，·)选用径向基核函数，表达式为：

式中，σ为核参数，与式（6）中的惩罚因子γ共同影响着最小二乘支持向量机的回归精度。

2.3 小波奇异值分解理论

2.3.1 小波理论 小波变换具有良好的时频局部化特性，可有效地将信号分解为包含不同频率成分的信号分量，目前广泛应用在信号降噪［14］、信号特征提取［15］等领域。

对信号f(t)进行连续小波变换，其表达式为：

将式（10）小波基函数ψ(∙)的尺度因子a和时间平移因子τ在同一尺度下离散化处理，按幂级数展开得a=a0j，τ=ka0j(j，k∈z)，取a0=2［16］，即得到信号的离散小波变换形式：

式中：j为离散小波变换中的尺度参数；k为时间平移参数；小波基函数ψ(t)通过j与k的伸缩平移得到小波信号ψj，k(t)，计算ψj，k(t)与信号f(t)的内积，最终得到信号展开系数DWT f(j，k)，实现对信号的离散小波变换。

使用多分辨率分析理论［17］将信号f(t)分解到不同的小波空间和尺度空间，为了便于使用计算机实现多分辨率分析，需要对分解系数进行离散化处理。对于能量有限的信号，采用有限精度分解方法，如下式所示：

式中：cj，k为尺度系数；ϕj，k(t)为尺度空间；dj，k为小波系数；ψj，k(t)为小波空间；变量j0 与j含义与式（11）相同，为尺度参数，表示小波空间和尺度空间对原信号的分辨率，在实际应用中j0 可理解为将信号f(t)累计应分解的层数，j为当前分解层数；k为时间平移参数，意义与式（11）相同。为了区别尺度参数j和时间平移参数k，展开系数cj，k和dj，k一般写为cj[k]、dj[k]，cj[k]与dj[k]可以通过Mallat塔式分解算法得到，计算式为：

式中，l[n]与h[n]分别为低通数字滤波器和高通数字滤波器。

该算法的原理如图1所示。

图1 Mallat 塔式分解算法

式（13）与图1表明，通过Mallat 算法，第j层的小波展开系数cj[k]与dj[k]都可以通过上一层展开系数cj-1[k]经过数字滤波器组l[n]和h[n]处理后得到，同理，可以从cj[k]按此步骤进一步得到展开系数cj+1[k]与dj+1[k]，如此重复，最终得到一系列不同尺度下的小波展开系数。

在水电机组实际状态监测中，振动信号的采集一般满足奈奎斯特采样定理，因此采集的离散振动信号f(n)可近似等于初始系数c0[k]，即：

2.3.2 奇异值分解理论 奇异值分解是一种重要的矩阵分解，其定义为：设定矩阵A(A ∈Crm×n)，r=rank(A)，存在酉矩阵U (U ∈Cm×m)与酉矩阵V (V ∈Cn×n)，使得：

式中：Σ=diag(σ1，σ2，…，σr)为对角阵;σ1、σ2、…、σr为矩阵A的奇异值，且σ1≥σ2≥…≥σr ＞0。

令：

式中，向量S 称为矩阵A 的奇异值特征向量。

3 水电机组状态劣化评估方法

3.1 基于检测指数与支持向量机的机组时域劣化评估

3.1.1 时域特征敏感性分析 不同时域特征参数对于机组故障的敏感性不同，本文选取了信号分析中具有代表性的6 个时域特征量，它们是：均值、标准差、峭度、偏度、峰峰值、均方差，具体计算方式见表1。

表1中，xi为信号x第i个采样点值，N为信号的长度。计算机组不同状态下振动信号时域特征参数的均值与方差，代入式（4）得到各时域特征参数的检测指数。关于时域特征个数的选取，若仅选一个最大值，容易忽略其他对区分机组状态有价值的时域特征；若选多个值时，又容易包含对机组状态变化不敏感的时域特征。鉴于检测指数最大的前两个时域特征其检测指数值显著大于其他时域特征的检测指数，且所选检测指数最大的前两个时域特征（标准差、峰峰值）已能够明显区分出机组三种不同运行状态，为此选取检测指数最大的两个时域特征参数代表振动信号的时域特征。

表1 常用时域特征参数

3.1.2 时域劣化评估指标计算 以健康状态下水电机组历史工况数据Xt（水头，开度等）作为输入样本，所选振动信号时域特征Y1t、Y2t作为输出样本，利用LS-SVM 建立机组状态健康模型：

式中，t为时间变量。

将机组状态监测系统采集到的实时工况数据代入健康模型即式（17）中，得到当前工况下振动信号时域特征健康值。计算信号时域特征实际值与健康值之间的相对误差，进行加权计算，作为当前时刻机组时域方向的劣化评估指标，计算方式为：

式中：S1t、S2t为振动信号时域特征健康值；Y1t、Y2t为振动信号时域特征实际值，本文选取的时域特征为标准差和峰峰值，由专家经验判定两者权重占比相同，均为0.5。

3.2 基于小波奇异值与聚类分析的机组频域劣化评估虽然振动信号时域特征在一定程度上能反映机组的状态劣化趋势，但实际中往往会存在振动信号时域特征没有明显改变，但在频域结构上信号特征发生明显变化的情况。因此，有必要构建机组频域方向的劣化评估指标。

利用小波变换对机组健康状态下的历史振动信号进行分解。根据式（13）和图1中的Mallat 分解算法得到包含信号不同频率成分的小波分解系数。将各分支小波系数进行差值重构，使其长度与原信号保持一致。重构后的小波系数作为奇异值输入矩阵A ，代入式（15），得到机组健康状态下反映振动信号频域特征的奇异值特征向量，即式（16）。结合数据挖掘技术中的K 均值算法［17］得到奇异值特征向量的健康聚类中心C，该聚类中心表征了机组在健康状态下，振动信号应有的标准奇异值特征向量，因此可作为频域方向评估机组运行状态的基准值。计算实时监测信号的奇异值特征向量与聚类中心C 之间的相对欧氏距离，作为频域方向的劣化评估指标，计算方式为：

式中， Qt为当前t时段实测振动信号的奇异值特征向量。

3.3 综合劣化评估将指标Tt与指标Ft加权组合，构建综合劣化评估指标：

式中，a和b为对应的权重系数，这里认为水电机组振动信号的时域成分与频域成分同等重要，故在本文中，取a=b=0.5。

水电机组振动状态综合劣化评估模型结构如图2所示。

图2 水电机组振动状态综合劣化评估模型

4 应用实例

以国内某电站3 号机组为研究对象，进行实例分析。已知该机组为轴流转桨式机组，2015年8月底机组上机架、转轴等处出现明显的振动及异音，经停机检修发现故障原因为机组转轮室里衬脱落，机组经维修后恢复正常运行。将所提出的水电机组状态劣化评估方法用于该机组故障发生前后的状态劣化评估，验证其有效性。

4.1 时域特征筛选对3 号机组轴向振动信号进行采样分析，根据轴向振动信号频谱特征，将机组故障前后的状态分为正常、故障预警和故障三种状态，每种状态各40 组信号数据。信号采样频率为458 Hz，每组信号长度为4096 个点。各状态的信号时域波形图如图3所示。

图3 三种状态轴向振动信号波形

计算每组数据时域特征参数，利用检测指数DI值进行特征敏感性分析，每种时域特征参数对应的DI值如表2所示。

表2 不同时域特征参数对应的DI 值

表2中，1、2、3 分别指机组正常、故障预警和故障3 种状态，从表2可以看出，标准差和峰峰值对应的DI值总体大于其他时域特征参数所对应的DI值，故选择标准差与峰峰值作为反映信号时域特征的特征量。在信号处理中，标准差与峰峰值也多用来作为信号分析的关键指标［19-20］，这也体现了DI值所选取特征的代表性。3 种状态下轴向振动波形的各时域特征参数值如图4所示，图中，信号序列1 ～40 为正常状态采集的40 组轴向振动信号编号，41 ～80 为故障预警状态采集的40 组轴向振动信号编号，81 ～120 为故障状态采集的40 组轴向振动信号编号。从图4可以看出，相较于其他时域特征参数，3 种不同状态之间的信号波形标准差与峰峰值差异明显，进一步印证了利用DI值作为敏感性指标选取时域特征的可靠性。

图4 不同运行状态下机组振动信号时域特征值

4.2 机组状态劣化评估

4.2.1 劣化评估健康样本构建 由于该机组为轴流转桨式机组，选择导叶开度、桨叶开度和工作水头作为机组工况数据。在发生故障前后时间段内，3 号机组的主要运行工况集中在导叶开度63%，桨叶开度20%，水头50 ～55 m，故选取该工况附近的数据点进行建模分析。

选取机组未发生故障时的正常数据（包含导叶开度X1t、桨叶开度X2t、水头H、轴向振动信号标准差σ与峰峰值f）共计190 组，利用LS-SVM 建立如式（20）所示的机组状态健康模型：

随机选取150 组数据作为训练集，20 组数据作为验证集，20 组数据作为测试集，对健康模型进行训练，所选数据覆盖导叶开度61% ～65%，桨叶开度16% ～20%，水头51 ～55 m 这一工况范围。考虑到在水电生产领域，监测系统采样得到的信号为实时值，对数据进行归一化不利于对水电机组进行实时状态评估，且LS-SVM 对于不同数据量级的数据依然具有良好的数学拟合特性，所以未对实验数据进行归一化处理。

参数优化方面，使用遗传算法［21］对LS-SVM 的核参数σ与惩罚因子γ进行寻优，遗传算法优化LS-SVM的算法流程如图所示。

图5中，LS-SVM 模型参数寻优范围为［1，500］，［1，500］，选取平均相对误差（MAPE）作为遗传算法的适应度函数。LS-SVM 预测效果评价方面，使用相对误差（Relative Error，RE）、平均绝对百分比误差（Mean Absolute Percentage Error，MAPE）、均方根误差（Root Mean Square Error，RMSE）作为评价指标。评价指标表达式为：

图5 遗传算法优化LS-SVM 算法流程

式中：yi为实际输出值；y′i为期望值。

LS-SVM 的参数寻优与训练结果如表3和图6、图7所示。

表3 LS-SVM 参数寻优与回归预测结果

从表3和图6、图7可以看出，训练后的LS-SVM 预测效果较好，误差（MAPE）分别为1.04%、2.86%，对应精度为98.96%、97.14%，可以准确的反映机组正常运行状态下所处工况与振动信号时域特征之间的映射关系，为进行水电机组时域方向的状态劣化评估提供了精度保证。

图6 预测参数为标准差的LS-SVM 训练结果

图7 预测参数为峰峰值的LS-SVM 训练结果

将同样190 组机组振动信号通过小波奇异值进行特征提取。根据3.2 节算法，利用离散小波变换将每组信号分解为包含不同频率成分的小波分解系数，分解小波基选用“DB8”小波，分解层数为4层。将得到的小波分解系数进行差值重构，使其长度与原信号保持一致，重构后的各分支小波分解系数如图8所示。图8中，5 组小波分解系数由上到下依次代表了信号由低频到高频的频率成分。将小波分解系数矩阵作为奇异值输入矩阵A ，此时A ∈，代入式（15）中进行奇异值分解，得到190 组奇异值特征向量，利用K 均值聚类得到奇异值向量的健康聚类中心，如表4所示。

图8 正常状态振动信号小波重构分解系数

4.2.2 机组状态实时跟踪评估 对3 号机组故障发生前到发生时这段时间的运行状态进行劣化评估分析，按时间顺序选取该时段内同工况下（导叶开度61% ～65%，桨叶开度16% ～20%，水头51 ～55 m）270 组数据。一方面，将每组工况数据输入已训练好机组状态健康模型，即式（21），计算当前振动信号标准差与峰峰值的健康值，求取与实际值之间的相对误差，代入式（18）得到机组时域劣化指标Tt；另一方面，运用小波奇异值理论对每组信号进行分解，提取奇异值特征向量，计算与表4健康聚类中心之间的相对欧式距离，根据式（19）得到机组频域劣化度指标Rt。将Tt与Rt代入式（20），最终得到反映3 号机组当前健康状态的综合劣化评估指标Dt。

表4 正常运行状态下奇异值特征向量及健康聚类中心

图9为3 号机组综合劣化评估指标随时间的变化趋势图，可以看出，随着时间的推移，机组状态劣化趋势明显，呈曲折上升态势，反映了机组由正常到故障的演变过程。通过与电站运行人员交流并查看机组在该时间段内的运行检修记录，发现在2015年8月15 至8月20 号运行期间，机组轴系各部位的振动、摆度、温度、气隙等各种运行参数幅值均处于正常水平，振摆信号的频率成分均无异常，且机组在运行时无异常响动，处于正常运行状态，从图9也可以看出，在该时段内机组综合劣化指标一直维持在较低水平且波动平稳。点C 处的综合劣化度首次突破0.2，可认为是机组故障开始的标志。根据该电站的事后分析报告，3 号机组在2015年8月28日开机过程中，机组上机架、水车室、蜗壳以及尾水管等处均有明显异常声音，通过检修处理发现故障。而根据图9的评估分析可以看出，该机组早在8月23日前后即出现发生故障的趋势，因此通过综合劣化评估指标可更好地实现机组故障预警，有效避免事故发生。

图9 机组状态劣化评估结果

为进一步印证图9所反映的机组振动劣化过程，从图9中选取6 组有代表性的数据A ～F 进行频谱分析，查看其1 ～5倍频、6倍频（转轮叶片数倍频）、12倍频（导叶片数倍频）幅值变化情况，如表5所示。

从表5可以看出，A - F 振动信号1 ～5 倍频均无明显差别，未能有效反映机组故障，6 倍频在F处大幅升高，12 倍频从A 到F 逐步递增，说明两者为故障特征频率，且变化趋势与图9所反映的结果相吻合，进一步表明了本文所提出的劣化评估方法能真实有效的反映机组劣化过程。

表5 不同信号各频率成分幅值对比

图10为3 号机组2017年7月某段时间的状态劣化趋势变化情况，此时机组处于正常运行状态，选取该时间段内同工况下10 组监测数据，利用劣化评估方法进行机组劣化状态评估。从图10可以看出，该时段内机组综合劣化度维持在一个很低水平，且始终小于0.2，符合正常运行状态下应有的劣化度特征，因此可将Di=0.2 作为机组综合劣化度上限值，若Di＞0.2 ，则说明机组出现故障征兆，需停机检查。

图10 机组状态劣化评估结果

4.2.3 对比实验 为进一步突出该评估指标的有效性，选取表1中轴向振动信号均值等6 个时域特征量和信号1 倍频、6 倍频、12 倍频幅值作为对比，绘制3 号机组轴向振动信号各特征量在15年8月故障发生前后随时间的变化趋势如图11所示。

从图11（a）（c）（g）可以看出，信号均值、均方差、1 倍频幅值随时间变化趋势不明显，未能反映机组运行状态的变化；图11（b）（h）中，信号峭度指标和6 倍频幅值虽在8月26日出现下降或升高，但在之前一段时间内并无明显变化，不能有效反映机组由正常到故障的劣化渐变过程；图11（d）（e）（f）（i）中，信号标准差、峰峰值、偏度和12 倍频幅值则随着时间的推移不断增大，反映了与图9机组综合劣化指标相同的机组劣化趋势，计算各指标的缓变率即变化程度，衡量指标对机组状态变化的灵敏度。分别选取图11（d）（e）（f）（i）中指标幅值变化比较明显的A、C、E 三个时刻（与图9中的A、C、E 点时刻相同），已知A 时刻点指标为正常值，分别计算C、E 时刻指标幅值相对于A 时刻的变化量，作为指标缓变率，计算方式如下式所示：

式中，yA、yC、yE分别为信号标准差、峰峰值、偏度和12 倍频在时刻A、C、E 的幅值。

各指标的缓变率如表6所示。从表6可以看出，综合劣化评估指标的缓变率远大于其他几个指标，说明该指标对机组运行状态变化的反应更为敏感，比传统的时域和频域指标更能体现出机组由正常到故障的劣化过程。

图11 机组信号时、频量随时间变化趋势

表6 指标缓变率

5 结论

本文提出了一种全新的水电机组状态劣化评估方法，在考虑运行工况的前提下，以机组振动信号为分析对象，利用检测指数理论、最小二乘支持向量机理论和小波奇异值理论，构建同时涵盖机组时频特征的状态综合劣化评估模型，最终通过计算劣化指标的形式实时评估机组的健康状态。结合电站机组实例进行验证分析，得到以下结果：

检测指数可有效反映机组振动信号各时域特征对不同运行状态的敏感性，从而选取最能反映信号特征的时域指标，避免了特征选取的盲目性。

该水电机组状态劣化评估方法能很好地对机组振动劣化状态进行实时跟踪评估，及时发现故障征兆，实现故障预警，并且充分利用机组监测系统历史数据，体现了数据挖掘的思想。

该方法能充分反映机组特色，所确定的综合劣化评估指标健康限值根据机组类型和不同运行工况会有所不同。针对本文所研究的3 号机组以及对应工况，其正常运行区间可确定为[0，0.2]。当劣化指标小于0.2 时，该机组处于正常运行状态；当指标大于0.2 时，则认为机组出现异常，监测系统报警，提醒运行人员关注。

与传统信号监测指标相比，本文所提出的综合劣化评估指标对水电机组故障的发生更加敏感，可更明显地反映出机组由正常到故障的状态变化。

实际电站机组运行案例的分析结果表明，该方法能有效反映机组劣化趋势，可为实现水电机组故障预警提供准确可靠的指导。