CN112857805B - 一种基于图相似度特征提取的滚动轴承故障检测方法 - Google Patents

Abstract

本发明公开一种基于图相似度特征提取的滚动轴承故障检测方法，从图相似度特征提取的角度来挖掘非平稳的滚动轴承振动信号中潜藏的特征，并在此基础上实施滚动轴承的故障检测。具体来讲，本发明方法首先将振动信号按照区间转变成能体现信号之间相似度的图矩阵，再利用图特征提取算法对图矩阵进行特征提取。最后，使用提取出的特征来判断滚动轴承是否出现了故障。本发明方法的优势在于：能将振动信号转变成能体现振动信号之间距离相似性的图矩阵，并对图矩阵进行特征提取，从时域的角度实现了滚动轴承振动信号的特征分析与提取，且不涉及复杂的运算，实施操作简单。

Description

一种基于图相似度特征提取的滚动轴承故障检测方法

技术领域

本发明涉及一种轴承故障检测方法，特别涉及一种基于图相似度特征提取的滚动轴承故障检测方法。

背景技术

随着工业技术向着智能化、科技化和集成化发展，现代机械装备的功能要求和使用条件逐步提升，对装备的可靠性也提出了更高的要求，对滚动轴承维修技术的研究逐步获得重视。滚动轴承作为支撑旋转体并降低摩擦系数的关键功能部件而广泛应用于旋转式的机械设备中，其工作性能的好坏直接关系到设备的正常运行。轴承在运转工作过程中不可避免的受到轴向载荷和循环交变载荷作用，易产生疲劳剥落或磨损等失效形式，若不根据使用情况及时更换轴承，严重时将会发生突发事故造成不必要的损失与伤亡。近十年随着人工智能和信息技术的发展，对轴承进行故障检测大多采用数据驱动的方法，但目前数据驱动的故障检测大都是使用振动信号的频域分析特征来进行滚动轴承的故障检测。

对滚动轴承故障诊断进行研究除了具有重要的实际意义之外，还具备深刻的理论意义。在实施滚动轴承故障诊断的过程中，可以推动测控技术、信号处理、模式分类、优化方法等与故障诊断相关的一系列理论方法的发展，不断完善机械故障诊断学科的发展。同时也为下一代新产品的更新换代提供了合理的优化改进依据。目前，滚动轴承的故障诊断技术仍以振动信号分析方法为主。但是由于轴承故障信号是一种非线性非平稳随机信号，传统的信号处理方法难以有效处理，因此，有必要对滚动轴承故障诊断技术进行深入地研究。

早期的轴承振动信号特征提取主要是借助于傅里叶变换进行的，称为经典信号分析方法。基于傅里叶变换的信号分析发展缓慢，主要原因是其计算量大。工程实际中，从机械设备上所测得的振动信号千变万化，其中大量为非平稳、非高斯分布和非线性的随机信号，尤其是当设备出现故障时上述情形显得更为突出。而基于变换的谱分析主要适用于平稳信号时不变信号，因此有必要寻找一种适合于非平稳信号分析的方法，它既能反映信号的频率内容，又能反映该频率内容随时间变化的规律，还能准确地反映出信号能量随时间的分布情况。

发明内容

本发明所要解决的主要技术问题是：如何从图相似度特征提取的角度来挖掘非平稳的滚动轴承振动信号中潜藏的特征，并在此基础上实施滚动轴承的故障检测。具体来讲，本发明方法首先将振动信号按照区间转变成能体现信号之间相似度的图矩阵，再利用图特征提取算法对图矩阵进行特征提取。最后，使用提取出的特征来判断滚动轴承是否出现了故障。

本发明方法解决上述问题所采用的技术方案为：一种基于图相似度特征提取的滚动轴承故障检测方法，包括以下所示步骤：

步骤(1)：利用加速度振动传感器采集滚动轴承在健康运行状态下的振动信号数据，记为振动信号向量z∈R^N×1，并同时记录加速度振动传感器采集信号的频率f和滚动轴承的转动频率r；其中，R^N×1表示N×1维的实数向量，N为振动信号数据的个数，频率的单位为赫兹(Hz)，R表示实数集。

步骤(2)：先设置移动窗口长度L等于f除以r的商，再从振动信号向量z∈R^N×1的第一个元素开始，利用移动窗口依次将z划分成n个L×1维的窗口向量x₁，x₂，…，x_n，具体的实施过程如下所示：

步骤(2.1)：初始化i＝1。

步骤(2.2)：将振动信号向量z中的第(i-1)×L+1个元素至第i×L个元素组建成第i个窗口向量x_i∈R^L×1；其中，R^L×1表示L×1维的实数向量。

步骤(2.3)：判断是否满足条件：(i+1)×L＞N；若否，则设置i＝i+1后，返回步骤(2.2)；若是，则得到n个L×1维的窗口向量x₁，x₂，…，x_n。

步骤(3)：根据如下所示公式计算图矩阵Z_i∈R^L×L中的第a行第b列元素Z_i(a，b)，重复本步骤直至得到图矩阵Z₁，Z₂，…，Z_n；其中，a∈{1，2，…，L}，b∈{1，2，…，L}，i∈{1，2，…，n}。

上式中，e为自然常数，x_i(a)和x_i(b)分别表示窗口向量x_i中的第a个元素和第b个元素，δ表示图矩阵参数，其取值范围是5L≤δ≤500L。由此可见，图矩阵Z_i是一个对称矩阵。

步骤(4)：根据如下所示步骤(4.1)至步骤(4.6)，依次逐个提取特征向量s₁，s₂，…，s_K，并保留左载荷向量u₁，u₂，…，u_K和右载荷向量v₁，v₂，…，v_K。

步骤(4.1)：初始化k＝1。

步骤(4.2)：分别初始化左载荷向量u_k和右载荷向量v_k为任意的L×1维的非零实数向量。

步骤(4.3)：若k＜2，则根据

计算矩阵G_u；若k≥2，则根据公式

计算矩阵G_u；再根据公式u_k＝G_uu_k/||G_uu_k||更新u_k；其中，

表示计算G_uu_k的长度，上标号T表示矩阵或向量的转置符号。

步骤(4.4)：根据

计算矩阵G_v后，再根据公式v_k＝G_vv_k/||G_vv_k||更新v_k。

步骤(4.5)：判断v_k是否收敛；若否，则返回步骤(4.3)；若是，则根据

计算特征向量s_k∈R^n×1中的第i个元素s_k(i)；其中，R^n×1表示n×1维的实数向量。

步骤(4.6)：判断特征向量s_k中所有元素的方差是否小于阈值

若否，则设置k＝k+1后，返回步骤(4.2)；若是，则设置K＝k-1后，得到K个特征向量s₁，s₂，…，s_K，并保留左载荷向量u₁，u₂，…，u_K和右载荷向量v₁，v₂，…，v_K。

值得注意的是，上述步骤(4.4)中的实施过程为图矩阵特征提取算法的数值求解方法。为了提取图矩阵Z₁，Z₂，…，Z_n中的特征向量，需通过左载荷向量u∈R^L×1和右载荷向量v∈R^L×1对各个图矩阵进行特征提取，从而使如下所示的目标函数最大化。

s.t.u^Tu＝1，v^Tv＝1

上式②的求解可使用经典的朗格朗日乘子法，即通过拉格朗日乘子λ和η构造如下所示的朗格朗日函数J。

其中，

然后，分别计算J相对于u和v的偏微分，并设置偏微分等于0。

由此可得到：

这两个公式是标准的特征值与特征向量的求解问题。若是在

的等号两边同时左乘u^T，则可以得到：

由于上式②是求解最大化问题，因此需要计算最大特征所对应的特征向量。

由此可见，u₁和v₁分别是

和

的最大特征值所对应的特征向量，并且满足单位长度要求。由于u₁和v₁的求解相互耦合，因此需要通过往复迭代的策略进行求解，也就是步骤(4.5)中有返回步骤(4.3)的实施过程。

当u₁和v₁求解完成后，还需求解更对的载荷向量u₂，u₃，…，u_K和v₂，v₃，…，v_K。为了保证载荷向量的多样性，需要在上式②增加约束条件，即：

上式⑥中，k≥2。同样可使用拉格朗日乘子法对上式⑥进行求解，以如下所示的朗格朗日函数相对于u_k的偏微分为例。

上式⑦可等价推理出：

也就是说，u_k为矩阵

的最大特征值所对应的特征向量。

值得说明的是，对于特征值求解问题：Aα＝λα，其最大特征值所对应的特征向量的数值计算方法通常是：随机初始化特征向量α为任意非零实数向量后，不断的通过α＝Aα后再实施归一化处理α＝α/||α||的方式直至α收敛，即可得到最大特征值所对应的特征向量。这就是为什么步骤(4.3)和步骤(4.4)中分别使用u_k＝G_uu_k/||G_uu_k||和v_k＝G_vv_k/||G_vv_k||更新u_k和v_k，而未曾出现求解特征值的实施过程。

步骤(5)：根据公式Λ＝S^TS/(n-1)计算协方差矩阵Λ∈R^n×n后，再计算检测指标向量

并将

中最大值元素记录为φ_max；其中，特征矩阵S＝[s₁，s₂，…，s_K]，diag{}表示将大括号内的矩阵对角线元素转变成列向量的操作。

步骤(6)：利用加速度振动传感器在线实时采集滚动轴承运行时的振动信号数据，当振动信号数据的采集个数达到L时，将这L个振动信号数据组建成窗口向量z_t∈R^L×1，并根据如下所示公式计算图矩阵X_t中的第a行第b列元素X_t(a，b)；其中，t表示窗口向量z_t对应的采样时刻。

上式中，z_t(a)和z_t(b)分别表示z_t中的第a个元素和第b个元素。

步骤(7)：根据公式y_t＝diag{U^TX_tV}计算特征向量y_t∈R^K×1后，再根据φ_t＝y_t ^TΛ^{- 1}y_t计算检测指标φ_t；其中，左载荷矩阵U＝[u₁，u₂，…，u_K]，右载荷矩阵V＝[v₁，v₂，…，v_K]。

步骤(8)：判断是否满足条件：φ_t≤φ_max；若是，则滚动轴承运行状态健康，并返回步骤(6)继续利用新采样时刻的窗口向量实施滚动轴承故障检测；若否，则执行步骤(9)。

步骤(9)：返回步骤(6)继续利用最新采集的L个振动信号数据实施故障检测，若连续6个采样时刻的检测指标皆大于φ_max，则触发滚动轴承的故障警报；反之，则滚动轴承运行状态健康。

通过以上所述实施步骤，本发明方法的优势介绍如下。

首先，本发明方法将振动信号转变成能体现振动信号之间距离相似性的图矩阵，并对图矩阵进行特征提取，从时域的角度实现了滚动轴承振动信号的特征分析与提取。其次，相比于利用频域分析实施滚动轴承故障检测的技术，本发明方法不涉及复杂的运算，实施操作简单。

附图说明

图1为本发明方法的实施流程示意图。

图2为本发明方法的故障检测结果图。

具体实施方式

下面结合附图和具体实施方式对本发明进行详细说明。

如图1所示，本发明公开了一种基于图相似度特征提取的滚动轴承故障检测方法，下面结合一个具体应用实例来说明本发明方法的具体实施方式。

在本实施案例中，用于测试故障检测性能的轴承故障具体包括：轴承内圈缺陷和轴承外圈缺陷。轴承在正常运行状态时，加速度振动传感器的采集信号频率是f＝12000Hz，轴承的转动频率是r＝60Hz。利用这些数据来验证本发明方法用于轴承故障检测可行性，具体包括如下所示步骤。

步骤(1)：利用加速度振动传感器采集滚动轴承在健康运行状态下的振动信号数据，记为振动信号向量z∈R^N×1，并同时记录加速度振动传感器采集信号的频率f和滚动轴承的转动频率r。

步骤(2)：先设置移动窗口长度L等于f除以r的商，再从振动信号向量z∈R^N×1的第一个元素开始，利用移动窗口依次将z划分成n个L×1维的窗口向量x₁，x₂，…，x_n，具体的实施过程如前述步骤(2.1)至步骤(2.3)所示。

步骤(3)：根据上述公式①计算图矩阵Z_i∈R^L×L中的第a行第b列元素Z_i(a，b)，重复本步骤直至得到图矩阵Z₁，Z₂，…，Z_n。

步骤(4)：根据前述步骤(4.1)至步骤(4.6)，依次逐个提取特征向量s₁，s₂，…，s_K，并保留左载荷向量u₁，u₂，…，u_K和右载荷向量v₁，v₂，…，v_K。

步骤(5)：根据公式Λ＝S^TS/(n-1)计算协方差矩阵Λ∈R^n×n后，再计算检测指标向量

并将

中最大值元素记录为φ_max。

上述步骤(1)至步骤(5)为本发明方法的离线建模阶段，完成后，即可利用实时的移动窗口数据来检测轴承是否出现了故障。

步骤(6)：利用加速度振动传感器在线实时采集滚动轴承运行时的振动信号数据，当振动信号数据的采集个数达到L时，将这L个振动信号数据组建成窗口向量z_t∈R^L×1，并根据如下所示公式计算图矩阵X_t中的第a行第b列元素X_t(a，b)。

步骤(7)：根据公式y_t＝diag{U^TX_tV}计算特征向量y_t∈R^K×1后，再根据φ_t＝y_t ^TΛ^{- 1}y_t计算检测指标φ_t；其中，左载荷矩阵U＝[u₁，u₂，…，u_K]，右载荷矩阵V＝[v₁，v₂，…，v_K]。

步骤(8)：判断是否满足条件：φ_t≤φ_max；若是，则滚动轴承运行状态健康，并返回步骤(6)继续利用新采样时刻的窗口向量实施滚动轴承故障检测；若否，则执行步骤(9)。

步骤(9)：返回步骤(6)继续利用最新采集的L个振动信号数据实施故障检测，若连续6个采样时刻的检测指标皆大于φ_max，则触发滚动轴承的故障警报；反之，则滚动轴承运行状态健康。

在图2中，本发明方法可以成功的检测出轴承内圈故障和外圈故障。由此可见，本发明方法在时域特征提取的基础上，即可实现轴承的故障检测。

Claims (1)

1.一种基于图相似度特征提取的滚动轴承故障检测方法，其特征在于，包括以下步骤：

步骤(1)：利用加速度振动传感器采集滚动轴承在健康运行状态下的振动信号数据，记为振动信号向量z∈R^N×1，并同时记录加速度振动传感器采集信号的频率f和滚动轴承的转动频率r；其中，R^N×1表示N×1维的实数向量，R表示实数集，N为振动信号数据的个数，频率的单位为赫兹；

步骤(2)：先设置移动窗口长度L等于f除以r的商，再从振动信号向量z∈R^N×1的第一个元素开始，利用移动窗口依次将z划分成n个L×1维的窗口向量x₁，x₂，…，x_n，具体的实施过程如步骤(2.1)至步骤(2.3)所示：

步骤(2.1)：初始化i＝1；

步骤(2.2)：将振动信号向量z中的第(i-1)×L+1个元素至第i×L个元素组建成第i个窗口向量x_i∈R^L×1；其中，R^L×1表示L×1维的实数向量；

步骤(2.3)：判断是否满足条件：(i+1)×L＞N；若否，则设置i＝i+1后，返回步骤(2.2)；若是，则得到n个L×1维的窗口向量x₁，x₂，…，x_n；

步骤(3)：根据如下所示公式计算图矩阵Z_i∈R^L×L中的第a行第b列元素Z_i(a，b)，重复本步骤直至得到图矩阵Z₁，Z₂，…，Z_n；其中，a∈{1，2，…，L}，b∈{1，2，…，L}，i∈{1，2，…，n}；

上式中，e为自然常数，x_i(a)和x_i(b)分别表示窗口向量x_i中的第a个元素和第b个元素，δ表示图矩阵参数，其取值范围是5L≤δ≤500L；

步骤(4)：根据如下所示步骤(4.1)至步骤(4.6)，依次逐个提取特征向量s₁，s₂，…，s_K，并保留左载荷向量u₁，u₂，…，u_K和右载荷向量v₁，v₂，…，v_K；

步骤(4.1)：初始化k＝1；

步骤(4.2)：分别初始化左载荷向量u_k和右载荷向量v_k为任意的L×1维的非零实数向量；

步骤(4.3)：若k＜2，则根据

计算矩阵G_u；若k≥2，则根据公式

计算矩阵G_u；再根据公式u_k＝G_uu_k/||G_uu_k||更新u_k；其中，

表示计算G_uu_k的长度，上标号T表示矩阵或向量的转置符号；

步骤(4.4)：根据

计算矩阵G_v后，再根据公式v_k＝G_vv_k/||G_vv_k||更新v_k；

步骤(4.5)：判断v_k是否收敛；若否，则返回步骤(4.3)；若是，则根据

计算特征向量s_k∈R^n×1中的第i个元素s_k(i)；其中，R^n×1表示n×1维的实数向量；

步骤(4.6)：判断特征向量s_k中所有元素的方差是否小于阈值

若否，则设置k＝k+1后，返回步骤(4.2)；若是，则设置K＝k-1后，得到K个特征向量s₁，s₂，…，s_K，并保留左载荷向量u₁，u₂，…，u_K和右载荷向量v₁，v₂，…，v_K；

步骤(5)：根据公式Λ＝S^TS/(n-1)计算协方差矩阵Λ∈R^n×n后，再计算检测指标向量

并将

中的最大值记录为φ_max；其中，特征矩阵S＝[s₁，s₂，…，s_K]，diag{}表示将大括号内的矩阵对角线元素转变成列向量的操作；

步骤(6)：利用加速度振动传感器在线实时采集滚动轴承运行时的振动信号数据，当振动信号数据的采集个数达到L时，将这L个振动信号数据组建成窗口向量z_t∈R^L×1，并根据如下所示公式计算图矩阵X_t中的第a行第b列元素X_t(a，b)；其中，t表示窗口向量z_t对应的采样时刻；

上式中，z_t(a)和z_t(b)分别表示z_t中的第a个元素和第b个元素；

步骤(7)：根据公式y_t＝diag{U^TX_tV}计算特征向量y_t∈R^K×1后，再根据φ_t＝y_t ^TΛ^-1y_t计算检测指标φ_t；其中，左载荷矩阵U＝[u₁，u₂，…，u_K]，右载荷矩阵V＝[v₁，v₂，…，v_K]；

步骤(8)：判断是否满足条件：φ_t≤φ_max；若是，则滚动轴承运行状态健康，并返回步骤(6)继续利用新采样时刻的窗口向量实施滚动轴承故障检测；若否，则执行步骤(9)；

步骤(9)：返回步骤(6)继续利用最新采集的L个振动信号数据实施故障检测，若连续6个采样时刻的检测指标皆大于φ_max，则触发滚动轴承的故障警报；反之，则滚动轴承运行状态健康。

Images