[发明专利]一种基于温度数据的风机齿轮箱子空间故障预测方法

权利要求书

1.一种基于温度数据的风机齿轮箱子空间故障预测方法，其特征在于，包括以下步骤：

A1温度数据的预处理；利用回归分析方法对温度数据进行单步预测，得到实际值和预测值之间的差值，称为残差，将残差作为随机状态空间模型的观测量Y;

A11多元线性回归模型

多元线性回归模型的一般形式如下：

y＝β₀+β₁x₁+…+β_px_p+ε   (1)

式(1)中，β₀,β₁,…,β_p是未知参数，β₀为回归常数，β₁,…,β_p为回归系数；y为因变量；x₁,x₂,…,x_p为自变量，这里是与因变量相关的监测量；ε为随机误差。若已知参数的估计值则可以实现温度的预测：

y^=β^0+β^1x1+...+β^pxp---(2)]]>

式(2)为经验回归方程，称为y的预测值;

随机误差ε，常假定服从正态分布：E(ε)＝0，Var(ε)＝σ²;

A12参数估计

已知n组监测数据(x_i1,x_i2,…,x_ip;y_i)，样本预测的误差：

yi-y^i=yi-β0-β1xi1-β2xi2...-βpxip---(3)]]>

采用最小二乘法估计多元线性回归模型的参数，即使式(4)取最小值时的解;

Q=Σi=1n(yi-y^i)2=Σi=1n(yi-β0-β1xi1-β2xi2-...-βpxip)2---(4)]]>

分别对β₀,β₁,…,β_p求偏导数，并令其等于零，然后联立求解即可求得回归参数的估计值

A13温度预测及残差求取

要预测k时刻的齿轮箱温度值需要与齿轮箱温度相关的监测量作为自变量：k时刻之前两个时刻的环境温度（T_e(k-1)、T_e(k-2)）、k时刻之前两个时刻的齿轮箱油温度（T_o(k-1)、T_o(k-2)）、k时刻之前两个时刻的齿轮轴承温度（T_b(k-1)、T_b(k-2)）。则建立模型的自变量为：X_(k)＝(T_e(k-1),T_e(k-2),T_o(k-1),T_o(k-2),T_b(k-1),T_b(k-2));

齿轮箱温度T_k，采用齿轮箱油温度和齿轮轴承温度的加权和来表示，即T_k＝0.5T_o(k)+0.5T_b(k);

上述参数估计用到的(x_i1,x_i2,…,x_ip;y_i)已经知道（X_(k)、T_k），由最小二乘法估计出模型的参数那么齿轮箱温度的预测值为温度的残差为ek=Tk-T^k;]]>

A2随机状态空间模型的识别

随机子空间的线性状态空间模型描述如下：

Xk+1=AXk+wkYk=CXk+vk---(5)]]>

方程(5)是随机状态空间模型，其中，X是状态量，一般不可测量，只是为了便于描述这个系统动态的数学对象；Y是观测量，在这里就是温度的残差；w是系统噪声，由于建模不精确和一些干扰造成的；v是测量噪声。这些量均为相应维数的列向量；A是系统矩阵；C是输出矩阵；

A21正交投影

定义由温度残差组成的分块Hankel矩阵Y：

其中，矩阵Y_p的维数q×N，矩阵Y_f的维数(p+1)×N，N值一般很大，除以的意义是对测量值进行标准化处理;

对Y阵进行重新分块，如式(7)所示：

其中，为Y_f的第一行移到了Y_p的末行之后的矩阵，为Y_f去掉了第一行后的矩阵;

首先，将Y_f正交投影到Y_p的空间上，由正交投影的定义可以计算出O_i的值;

由于O_i一般是非常庞大的矩阵，在实际计算时，我们为了保证数值的稳定性，通常先对Y阵进行分块和LQ分解，得到一个稀疏的下三角阵，也对其进行分块，如式(9)所示，经过推导，投影矩阵O_i亦可以由式(10)求得；又根据随机子空间识别理论，投影矩阵O_i可分解为可观矩阵Γ_i和Kalman滤波状态序列的乘积。

Oi=Yf/Yp=L21L31Q1=ΓiX^i---(10)]]>

Γi=CCACA2...CAi-1---(11)]]>

同样的，对Y阵进行二次分块和LQ分解，并将正交投影到的空间上，得到了O_i-1的表达式(12);其中，Γ_i-1为Γ_i去掉最后一行CA^i-1后的矩阵;

Oi-1=Yf-/Yp+L31L32Q1Q2=Γi-1X^i+1---(13)]]>

A22奇异值分解

对式(10)中的L阵进行奇异值分解：

L21L31=U1U2S100S2V1TV2T=U1S1V1T---(14)]]>

其中U、V为酉矩阵；S₁为对角阵，对角阵元素由大到小排列，系统的阶数为非零奇异值的个数，但在实际计算中由于噪声的影响，比较高次的奇异值不会完全等于0，而是和0比较接近的数，在这种情况下，通常奇异值会有一个跳跃，即σ_i+1<<σ_i，选择这个相对很小的奇异值σ_i+1，令这个奇异值近似等于0，并且从这个奇异值开始之后的值都近似处理为0即（σ_i+1≈σ_i+2≈…≈0）;

根据奇异值分解(14)，得到可观矩阵Γ_i和状态量的表达式：

Γi=U1S11/2---(15)]]>

X^i=S11/2V1TQ1---(16)]]>

再由式(13)，得到的表达式，其中为求Moore-Penrose逆;

A23估计系数矩阵A、C

已经估计出状态量根据随机状态空间模型，得到下面的方程;

X^i+1Yi|i=AC(X^i)+ρwρv---(18)]]>

式中的Y_i|i就是之前出现过的Y_q+1|q+1，为了表示方便，这里对下标进行了改变;ρ_w、ρ_v是残差向量，与不相关;利用最小二乘思想，即令残差项最小时可以解出这个方程，得到A、C阵的估计值;

A3齿轮箱的故障预测

利用大量历史温度残差数据按照A2中所述的识别方法，得到很多组阵（k＝1,2…），并求取阵的特征值，在正常情况下，每一次计算所得特征值相差不大，也就是说系统的特征是基本一致的，均为正常状态;在这里定义参考特征值它是历次阵特征值的平均;

有了参考特征值，利用它与实时数据所求的特征值进行比较来判断齿轮箱的状态好坏;在这里定义离散度R作为评价指标，即由实时数据所求每个特征值与对应参考特征值之间距离的平均值，公式如下：

R=1nΣi=1n(xi-xoi)2+(yi-yoi)2---(20)]]>

其中，n为系统阶数；x_i、y_i为由实时温度残差数据求得的特征值的实部和虚部；x_oi、y_oi为参考特征值的实部和虚部;

每次计算都得到一个R值，取一段正常运行时间所求得的离散度值，计算其均值μ和根方差σ，定义预警门槛值R_k为：

R_k＝μ+3σ   (21)

有了门槛值R_k就可以通过监测R的变化情况，来实时评估齿轮箱的状态好坏，如果R越过所设定的门槛值，则说明齿轮箱有发生故障的风险。