一种基于ICEEMDAN-SVD的滚动轴承故障诊断方法

一种基于iceemdan-svd的滚动轴承故障诊断方法
技术领域
1.本发明涉及旋转机械故障类型识别和诊断技术领域，具体涉及一种基于iceemdan-svd的滚动轴承故障诊断方法。


背景技术：

2.随着近年全球风力发电行业发展，我国风电装机量进一步增加。风力机发电机组绝大多数运行在工况恶劣多变的露天环境，运行工况复杂多变，各种相关因素导致风电机组容易发生故障。齿轮箱作为风力发电机传动链中关键部件之一，齿轮箱组成部件多，容易发生故障并且齿轮箱故障造成的停机时间和维护成本最高。齿轮箱由多个齿轮、轴承和轴构成。高速轴位于齿轮箱末端与发电机相连。风力发电机在发电期间高速轴的转速在每分钟1500到1800转。在高速轴和发电机之间，发生任何的微小振动都可能导致轴承损坏。提高在噪声环境下的滚动轴承故障诊断精度非常重要，可以有效降低设备维修和维护成本，提高机械设备运行的经济性，保证安全生产。


技术实现要素：

3.本发明的目的在于提供一种基于iceemdan-svd的滚动轴承故障诊断方法，用以解决emd、eemd、ceemdan等方法存在的伪模态和残留噪声的问题，同时结合相关系数－峭度筛选最优分量，通过svd分解，将信号重构，最后进行包络谱分析。
4.为解决上述技术问题，本发明所采用的技术方案是：
5.一种基于iceemdan-svd的滚动轴承故障诊断方法，包括以下步骤：
6.步骤1：获取风力机高速轴承故障数据；
7.步骤2：计算高速轴承峭度指标，判断轴承是否发生故障；
8.步骤3：对异常信号进行iceemdan分解，得到多个模态分量；
9.步骤4：对模态分量计算person相关系数和峭度，根据相关系数－峭度准则选择最优分量；
10.步骤5：对最优分量进行svd奇异值分解，得到重构信号；
11.步骤6：对svd处理信号进行包络解调分析，提取故障特征对比滚动轴承理论计算值，诊断得到滚动轴承故障类型。
12.本发明技术方案的进一步明细在于：所述步骤2包括：
13.峭度系数计算公式：
[0014][0015]
式中，q为分析的轴承振动信号，μ为信号q的均值；σ为信号q的标准差。
[0016]
本发明技术方案的进一步明细在于：所述步骤3包括以下步骤：
[0017]
步骤3.1：定义x(t)为带分解信号，ek(
·
)表示由局部均值分解产生的k阶模态分量，n(
·
)表示产生信号的局部均值，ω(i)代表高斯噪声，imfi表示第i个模态分量；
[0018]
步骤3.2：向待处理信号x(t)中分别加入多次白噪声ω(i)，构造序列，得到第一组残差：
[0019]
r1＝(n(x(i)))
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(2)
[0020]
步骤3.3：计算第一模态分量：
[0021]
imf1＝x-r1ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(3)
[0022]
步骤3.4：继续添加白噪声，利用局部均值分解计算第二组残差：
[0023]
r1+β1e1(ω(i))
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(4)
[0024]
步骤3.5：定义第二个模态分量imf2：
[0025]
imf2＝r
1-r2＝imf
1-(n(imf1+β1e1(ω(i))))
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(5)
[0026]
步骤3.6：计算第k个残差：
[0027]rk
＝(n(r
k-1
+β
k-1ek-1
(ω(i))))
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(6)
[0028]
第k个模态分量：
[0029]
imfk＝r
k-1-rkꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(7)
[0030]
步骤3.7：直至计算分解结束，得到所有模态分量与残差数。
[0031]
本发明技术方案的进一步明细在于：所述步骤4包括以下步骤：
[0032]
步骤4.1：person相关系数计算公式如下：
[0033][0034]
式中，z和y分别为两个相关变量，e()为数学期望；
[0035]
步骤4.2：峭度系数计算公式如下：
[0036][0037]
式中，q为分析的轴承振动信号，μ为信号q的均值；σ为信号q的标准差；
[0038]
步骤4.3：根据person相关系数筛选出相关系数大于0.5的模态分量且选择峭度系数最大的分量为最优分量。
[0039]
本发明技术方案的进一步明细在于：所述步骤5包括以下步骤：
[0040]
步骤5.1：对信号构造hankel矩阵，设p＝(p(1),p(2),
…
,p(l))是长度为l的振动信号，构造hankel矩阵：
[0041][0042]
式中，特征矩阵a∈rm×n，1《n《l；l＝n+m-1；
[0043]
步骤5.2：对矩阵a进行svd分解；
[0044]
a＝usv
t
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(11)
[0045]
式中，u和v是正交矩阵；s为非负对角矩阵，s的对角元素成为奇异值λi(i＝1,2,
…
k),且奇异值存在如下关系：λ1≥λ2≥
…
λk≥0；
[0046]
步骤5.3：确定svd重构分量个数，为进行合适降噪，采用奇异值差分谱方法，其定义如下：bi＝λ
i-λ
i+1
(i＝1,2,
…
,n-1)，则序列b＝(b1,b2,
…
,b
p-1
)即为奇异值的差分谱；整个差分谱中存在一个最大峰值bk，bk即为最大突变点；最大突变点代表着理想信号和噪声的分界。
[0047]
本发明技术方案的进一步明细在于：所述步骤6具体包括：滚动轴承内圈特征频率计算公式如下：
[0048][0049]
式中，bpfi表示内圈特征频率的理论值，z代表滚动体个数，d表示滚动体直径，d表示节圆直径，α表示轴承接触角，fr表示转频。
[0050]
由于采用了上述技术方案，本发明相对现有技术来说，取得的技术进步是：
[0051]
本发明提供一种基于iceemdan-svd的滚动轴承故障诊断方法，与现有技术相比，可以对实际的振动数据进行有效的降低噪声，提取有效的故障特征信息，在轴承故障早期能够发现故障特征，相比于其他故障诊断方法，iceemdan解决了经验模态分解、集合经验模态分解、完备经验模态分解存在的残留噪声和伪模态问题，svd分解可以去除数据中的噪声和冗余信息，通过选择最重要的奇异值进行重构信号，基于iceemdan-svd的滚动轴承故障判断方法不依靠个人经验选取参数，符合工业简易、实用和可移植性强的准则。
附图说明
[0052]
图1为本发明提供的基于iceemdan-svd的滚动轴承故障诊断方法系统图；
[0053]
图2为本发明的风力发电机高速轴承时域数据图；
[0054]
图3为本发明的风力发电机高速轴承峭度系数图；
[0055]
图4为本发明的奇异值差分谱图；
[0056]
图5为本发明的奇异值重构信号图；
[0057]
图6为本发明的信号包络谱图。
具体实施方式
[0058]
下面结合实施例对本发明做进一步详细说明：
[0059]
实施例1
[0060]
如图1所示，本实施例提供一种基于iceemdan-svd的滚动轴承故障诊断方法，包括以下步骤：
[0061]
步骤1：获取风力机高速轴承故障数据。
[0062]
步骤2：峭度系数计算公式：
[0063][0064]
式中，q为分析的轴承振动信号，μ为信号q的均值；σ为信号q的标准差。
[0065]
步骤3.1：定义x(t)为带分解信号，ek(
·
)表示由emd分解产生的k阶模态分量，n(
·
)表示产生信号的局部均值，ω(i)代表高斯噪声。
[0066]
步骤3.2：向待处理信号x(t)中分别加入多次白噪声ω(i)，构造序列，得到第一组残差：
[0067]
r1＝(n(x(i)))(2)
[0068]
步骤3.3：计算第一模态分量：
[0069]
imf1＝x-r1(3)
[0070]
步骤3.4：继续添加白噪声，利用局部均值分解计算第二组残差：
[0071]
r1+β1e1(ω(i))
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(4)
[0072]
步骤3.5：定义第二个模态分量imf2：
[0073]
imf2＝r
1-r2＝imf
1-(n(imf1+β1e1(ω(i))))
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(5)
[0074]
步骤3.6：计算第k个残差：
[0075]rk
＝(n(r
k-1
+β
k-1ek-1
(ω(i))))
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(6)
[0076]
第k个模态分量：
[0077]
imfk＝r
k-1-rkꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(7)
[0078]
步骤3.7：直至计算分解结束，得到所有模态与残差数。
[0079]
步骤4包括：
[0080]
步骤4.1：person相关系数计算公式如下：
[0081][0082]
式中，z和y分别为两个相关变量，e()为数学期望；
[0083]
步骤4.2：峭度系数计算公式如下：
[0084][0085]
式中，q为分析的轴承振动信号，μ为信号q的均值；σ为信号q的标准差。
[0086]
步骤4.3：根据person相关系数筛选出相关系数大于0.5的本征模态分量且选择峭度系数最大的分量为最优分量。
[0087]
步骤5.1：对信号构造hankel矩阵，设p＝(p(1),p(2),
…
,p(l))是长度为l的振动信号，构造hankel矩阵：
[0088][0089]
式中，特征矩阵a∈rm×n，1《n《l；l＝n+m-1；
[0090]
步骤5.2：对矩阵a进行svd分解。
[0091]
a＝usv
t
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(11)
[0092]
式中，u和v是正交矩阵；s为非负对角矩阵，s的对角元素成为奇异值λi(i＝1,2,
…
k)，且奇异值存在如下关系：λ1≥λ2≥
…
λk≥0。
[0093]
步骤5.3：确定svd重构分量个数，为进行合适降噪，采用奇异值差分谱方法，其定
义如下：bi＝λ
i-λ
i+1
(i＝1,2,
…
,n-1)。则序列b＝(b1,b2,
…
,b
p-1
)即为奇异值的差分谱。整个差分谱中存在一个最大峰值bk，即为最大突变点。最大突变点代表着理想信号和噪声的分界。
[0094]
步骤6.1：滚动轴承内圈特征频率计算公式如下：
[0095][0096]
式中，bpfi表示内圈特征频率的理论值，z代表滚动体个数，d表示滚动体直径，d表示节圆直径，α表示轴承接触角，fr表示转频。
[0097]
如图2所示2mw风力发电机高速轴承振动数据集，该数据集每天采样6秒，采样频率为97656hz，共采集了50天的数据。最后发现风力发电机高速轴有内圈故障。该轴承类型为skf32222j2型圆锥滚子轴承，具体参数如表1。风力发电机组高速轴转动频率大约为30hz。通过表1的轴承结构参数可以计算出轴承内圈理论特征频率bpfi＝284.1hz。
[0098]
表1 skf轴承参数
[0099][0100]
如图3所示，该图为数据集的峭度系数，同时根据时域峭度图可以将时域数据分为三个阶段：正常、退化和失效三个阶段。根据峭度图可以发现该轴承在0到25天是正常状态，在26天到48天之间轴承是处于退化状态，在49天以后轴承处于故障状态。
[0101]
选取退化阶段数据进行iceemdan处理，得到若干本征模态分量，根据相关系数－峭度准则重构分量，如表2所示相关系数前6个本征模态分量。根据准则选取imf1为最优分量进行处理。将最优分量构造hankel矩阵，进行svd分解，如图4所示，奇异值差分谱，选择最大奇异点进行重构信号，选择第6点重构信号，其余奇异值全部置零。如图5所示为重构信号。将重构信号进行包络谱处理。如图6所示，内圈故障特征频率计算为bpfi＝283hz，故障特征频率明显，存在3个明显谱峰，谱图干净，无任何干扰成分，故障特征频率及其倍频成分被准确提取出来。
[0102]
表2本征模态分量相关系数－排列熵值
[0103][0104]
本发明的有益效果是，与现有技术相比，可以对实际的振动数据进行有效的降低噪声，提取有效的故障特征信息，在轴承故障早期能够发现故障特征。相比于其他故障诊断方法，iceemdan解决了经验模态分解、集合经验模态分解、完备经验模态分解存在的残留噪声和伪模态问题。svd分解可以去除数据中的噪声和冗余信息。通过选择最重要的奇异值进行重构信号。基于iceemdan-svd的滚动轴承故障判断方法不依靠个人经验选取参数，符合工业简易、实用和可移植性强的准则。
[0105]
上文一般性的对本发明做了详尽的描述，但在本发明基础上，可以对之做一些修改或改进，这对于技术领域的一般技术人员是显而易见的。因此，在不脱离本发明思想精神的修改或改进，均在本发明的保护范围之内。

技术特征：
1.一种基于iceemdan-svd的滚动轴承故障诊断方法，其特征在于：该基于iceemdan-svd的滚动轴承故障诊断方法，包括以下步骤：步骤1：获取风力机高速轴承故障数据；步骤2：计算高速轴承峭度指标，判断轴承是否发生故障；步骤3：对异常信号进行iceemdan分解，得到多个模态分量；步骤4：对模态分量计算person相关系数和峭度，根据相关系数－峭度准则选择最优分量；步骤5：对最优分量进行svd奇异值分解，得到重构信号；步骤6：对svd处理信号进行包络解调分析，提取故障特征对比滚动轴承理论计算值，诊断得到滚动轴承故障类型。2.根据权利要求1所述的一种基于iceemdan-svd的滚动轴承故障诊断方法，其特征在于：所述步骤2包括：峭度系数计算公式：式中，q为分析的轴承振动信号，μ为信号q的均值；σ为信号q的标准差。3.根据权利要求1所述的一种基于iceemdan-svd的滚动轴承故障诊断方法，其特征在于：所述步骤3包括以下步骤：步骤3.1：定义x(t)为带分解信号，ek(
·
)表示由局部均值分解产生的k阶模态分量，n(
·
)表示产生信号的局部均值，ω(i)代表高斯噪声，imfi表示第i个模态分量；步骤3.2：向待处理信号x(t)中分别加入多次白噪声ω(i)，构造序列，得到第一组残差：r1＝(n(x
(i)
))(2)步骤3.3：计算第一模态分量：imf1＝x-r1ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(3)步骤3.4：继续添加白噪声，利用局部均值分解计算第二组残差：r1+β1e1(ω(i))
ꢀꢀꢀꢀꢀꢀ
(4)步骤3.5：定义第二个模态分量imf2：imf2＝r
1-r2＝imf
1-(n(imf1+β1e1(ω
(i)
)))
ꢀꢀꢀꢀꢀꢀ
(5)步骤3.6：计算第k个残差：r
k
＝(n(r
k-1
+β
k-1
e
k-1
(ω
(i)
)))
ꢀꢀꢀꢀꢀꢀꢀ
(6)第k个模态分量：imf
k
＝r
k-1-r
k
ꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(7)步骤3.7：直至计算分解结束，得到所有模态分量与残差数。4.根据权利要求1所述的一种基于iceemdan-svd的滚动轴承故障诊断方法，其特征在于：所述步骤4包括以下步骤：步骤4.1：person相关系数计算公式如下：
式中，z和y分别为两个相关变量，e()为数学期望；步骤4.2：峭度系数计算公式如下：式中，q为分析的轴承振动信号，μ为信号q的均值；σ为信号q的标准差；步骤4.3：根据person相关系数筛选出相关系数大于0.5的模态分量且选择峭度系数最大的分量为最优分量。5.根据权利要求1所述的一种基于iceemdan-svd的滚动轴承故障诊断方法，其特征在于：所述步骤5包括以下步骤：步骤5.1：对信号构造hankel矩阵，设p＝(p(1),p(2),
…
,p(l))是长度为l的振动信号，构造hankel矩阵：式中，特征矩阵a∈r
m
×
n
，1<n<l；l＝n+m-1；步骤5.2：对矩阵a进行svd分解；a＝usv
t
ꢀꢀꢀꢀꢀꢀꢀ
(11)式中，u和v是正交矩阵；s为非负对角矩阵，s的对角元素成为奇异值λi(i＝1,2,
…
k),且奇异值存在如下关系：λ1≥λ2≥
…
λ
k
≥0；步骤5.3：确定svd重构分量个数，为进行合适降噪，采用奇异值差分谱方法，其定义如下：b
i
＝λ
i-λ
i+1
(i＝1,2,
…
,n-1)，则序列b＝(b1,b2,
…
,b
p-1
)即为奇异值的差分谱；整个差分谱中存在一个最大峰值b
k
，b
k
即为最大突变点；最大突变点代表着理想信号和噪声的分界。6.根据权利要求1所述的一种基于iceemdan-svd的滚动轴承故障诊断方法，其特征在于：所述步骤6具体包括：滚动轴承内圈特征频率计算公式如下：式中，bpfi表示内圈特征频率的理论值，z代表滚动体个数，d表示滚动体直径，d表示节圆直径，α表示轴承接触角，fr表示转频。

技术总结
本发明公开了一种基于ICEEMDAN-SVD的滚动轴承故障诊断方法，涉及旋转机械故障类型识别和诊断技术领域，包括获取风力机高速轴承故障数据；计算高速轴承峭度指标，判断轴承是否发生故障；对异常信号进行改进自适应完备噪声模态分解，得到多个模态分量；对模态分量计算Person相关系数和峭度，根据相关系数－峭度准则选择最优分量；对最优分量进行SVD奇异值分解，根据奇异值差分谱进行重构信号；对重构信号进行包络解调分析，提取故障特征对比滚动轴承理论计算值，诊断得到滚动轴承故障类型；本发明可以解决模态分解所造成的模态混淆和残留噪声的问题，不用人为选定参数。符合工业界简易、实用和可移植性强的准则，具有较高应用价值。价值。价值。


技术研发人员：张文广 陈松 陈文华
受保护的技术使用者：华北电力大学
技术研发日：2023.06.21
技术公布日：2023/9/19