CN104075889A - 一种基于故障特征的零空间轴承复合故障诊断方法 - Google Patents

Abstract

一种基于故障特征的零空间轴承复合故障诊断方法，该方法包括根据轴承故障机理建立轴承故障振动模型，依据故障振动模型利用Matlab编程构造基于轴承故障特征的微分算子、对轴承故障非平稳信号进行基于故障特征的零空间微分算子的自适应分解、提取包含故障特征的局部窄带信号、解调得到故障特征。基于轴承故障的零空间微分算子是依据轴承的故障振动模型构造的且不同的故障特征对应算法中的参数不同，进一步对局部窄带信号进行解调便可以提取出不同的轴承故障特征进而实现轴承的复合故障诊断。

Description

一种基于故障特征的零空间轴承复合故障诊断方法

技术领域

本发明涉及一种轴承复合故障诊断方法，属于轴承故障诊断技术领域，尤其涉及一种基于故障特征的零空间轴承复合故障诊断方法。

背景技术

轴承是旋转类机械的重要组成部件，对其运转状态的检测和故障诊断具有很重要的意义。轴承的故障振动信号是一类典型的非平稳信号且含有大量噪声等干扰，较平稳信号而言，非平稳信号的分析是一类十分复杂的非线性问题，尤其在发生复合故障的情况，由于故障间的相互影响、干扰使得故障特征变得复杂进而大大增加了准确检测故障的难度，然而工程实际中所接触的信号往往也是此类信号，所以此类信号的研究对于工程应用具有极其重要的意义。

在传统的信号分析方法中，基本的分析方法包括时域分析和频域分析。时域分析形式简易直观，而频域分析更能反映信号内涵的信息。然而对于复杂的非平稳信号而言，单纯的时域表示或频域表示都不能完整刻画信号富含的特征信息。因此时频分析方法应运而生。典型的时频分析方法有短时傅里叶变化、Wigner-Ville分布、小波变换、Hilbert-Huang变换、EMD分解等，但一般的时频分析方法由于其分解基函数的单一对复杂信号的表达缺乏自适应性。EMD算法虽然具备一定的自适应性，但是在分析故障成分频率接近的复合故障时，会产生模式混叠和虚假分量。

零空间追踪算法是一种基于微分算子的自适应信号分解方法，其中微分算子是基于待分析信号的振动模型构造的。2008年，SL Peng和WL Hwang等人在EMD分解算法的基础之上提出了零空间追踪算法(NSP)。NSP算法主要是通过求解一个最优化问题的正则化方法来实现的，其核心思想是局部窄带信号在奇异局部线性微分算子的作用下“消失”，因此，可以通过奇异局部线性微分算子来抽取待分析信号的局部窄带分量，并将得到的局部窄带信号作为单元信号，用其叠加来逼近原始信号，进而实现信号的自适应分解。2010年，SL Peng对零空间追踪算法进行了进一步优化，使得该算法可以有效地提取调频信号。2011年，胡晰远在此基础上将零空间追踪算法进行了推广，实现了调幅调频信号的提取。2013年肖维维等人提出了自适应三阶微分算子的零空间追踪算法，将算子推广到三阶线性微分算子，它可以提取复指数信号，随后又进一步的提出了四阶线性微分算子。目前零空间追踪算法主要应用在图像处理领域，在机械故障诊断上的应用比较少见。

发明内容

为了解决传统的信号分析方法在轴承故障诊断中的上述技术问题，本发明的目的在于提供了一种基于故障特征的零空间轴承复合故障诊断方法。

为实现上述目的，本发明采用的技术方案为一种基于故障特征的零空间轴承复合故障诊断方法，该方法包括根据轴承故障机理建立轴承故障振动模型，依据故障振动模型利用Matlab编程构造基于轴承故障特征的微分算子、对轴承故障非平稳信号进行基于故障特征的零空间微分算子的自适应分解、提取包含故障特征的局部窄带信号、解调得到故障特征。

所述分解算法包括以下步骤：

(1)根据轴承故障机理建立故障振动模型；

(2)根据故障模型构造基于轴承故障的零空间微分算子；

(3)利用构造的基于轴承故障的零空间微分算子对原始信号进行分解，通过调整算法中的参数得到包含不同故障特征的窄带信号；

(4)对包含不同故障冲击信号成分的窄带信号进行解调分析得到不同故障特征。

轴承振动信号由轴承的旋转运动引起，故障轴承振动信号中还会出现冲击和瞬态振动特征，针对这一特点及动力学可知滚动轴承故障振动形式可以抽象为质量-弹簧-阻尼系统，质量-弹簧-阻尼系统的动力学模型为：

my″+cy′+ky＝0

其中m为质量-弹簧-阻尼系统中质量元件的质量，c为质量-弹簧-阻尼系统中阻尼器的阻尼系数，k为质量-弹簧-阻尼系统中弹簧元件的弹性系数，y为系统中质量元件的位移关于时间t函数，y′为y的一阶导数，y″为y的二阶导数，为后续计算方便令n＝c/2m，ξ＝n/ω，ω为有阻尼振动角频率，ξ为阻尼比，n为一常数。

则上式可改写为如下形式：

y″+2ξωy'+ω²y＝0

上式的解为：

y(t)＝y₀e^-ξωtcos(ωt)

其中，y₀为质量元件位移函数当时间t＝0时的初始值。

即轴承故障振动模型为周期出现的指数衰减函数，可表示为：

$x\left(t\right)=\underset{-\∞}{\overset{+\∞}{\mathrm{\Σ}}}{y}_0{e}^{-\mathrm{nt}}\cos \mathrm{\ωt\δ}(t-k{T}_0)$

其中δ(t)为冲击函数T₀为冲击的周期，模型的有效成分为指数衰减信号即：

y(t)＝y₀e^-ntcos(ωt)

综上可知轴承故障模型的有效成分(指数衰减函数)是二阶微分方程y″+2ξωy'+ω²y＝0的解也就是轴承的故障振动模型处在二阶微分算子y″+2ξωy'+ω²y的零空间内。对于轴承故障信号模型：y(t)＝y₀e^-ξωtcos(ωt)容易证明其处在微分算子d²/dt²-2a(t)'d/a(t)dt+ω(t)²+2(a(t)'/a(t))²的零空间内。其中a(t)＝y₀e^-ξωt。这样便可以利用上述零空间微分算子对信号进行自适应分解，由于不同类型的故障在算法中对应的参数不同，因此可以通过改变算法中的参数实现轴承的复合故障诊断。

与现有技术相比，本发明具有如下有益效果。

利用基于轴承故障的零空间微分算子自适应地将传感器采集的含有大量噪声及干扰的轴承复合故障非平稳信号分解为一系列的包含不同故障特征的局部窄带信号的叠加。基于轴承故障的零空间微分算子是依据轴承的故障振动模型构造的且不同的故障特征对应算法中的参数不同，这样通过改变算法中的参数可以得到包含不同故障特征的局部窄带信号，进一步对局部窄带信号进行解调便可以提取出不同的轴承故障特征进而实现轴承的复合故障诊断。

附图说明

图1是基于故障特征的零空间轴承复合故障诊断方法整体流程图。

图2是轴承内外圈复合故障仿真信号的染噪后的时域波形图及频谱图。

图3是轴承内外圈复合故障仿真信号的染噪后解调得到的解调谱图。

图4是本发明中对内外圈复合故障处理后的时域图及频谱图。

图5是本发明中对内外圈复合故障处理后进行解调处理得到的解调谱。

具体实施方式

以下结合附图和实施例对本发明作进一步详细说明。

如图1-5所示，一种基于故障特征的零空间轴承复合故障诊断方法，(1)根据轴承的故障机理建立轴承故障振动模型。针对信号的结构特点及动力学可知滚动轴承故障振动形式可以近似看作质量-弹簧-阻尼系统，质量-弹簧-阻尼系统的动力学模型为：

my″+cy′+ky＝0

其中m为质量-弹簧-阻尼系统中质量元件的质量，c为质量-弹簧-阻尼系统中阻尼器的阻尼系数，k为质量-弹簧-阻尼系统中弹簧元件的弹性系数，y为系统中质量元件的位移关于时间t函数，y′为y的一阶导数，y″为y的二阶导数，为后续计算方便令n＝c/2m，ξ＝n/ω，ω为有阻尼振动角频率，ξ为阻尼比，n为一常数。

则上式可改写为如下形式：

y″+2ξωy'+ω²y＝0

上式的解为：

y(t)＝y₀e^-ξωtcos(ωt)

即轴承故障振动模型为周期出现的指数衰减函数，可表示为：

$x\left(t\right)=\underset{-\∞}{\overset{+\∞}{\mathrm{\Σ}}}{y}_0{e}^{-\mathrm{nt}}\cos \mathrm{\ωt\δ}(t-k{T}_0)$

其中δ(t)为冲击函数T₀为冲击的周期，模型的有效成分为指数衰减信号即：

y(t)＝y₀e^-ntcos(ωt)

(2)根据轴承故障模型构造基于轴承故障的零空间微分算子。由(1)易知轴承故障模型的有效成分(指数衰减函数)是二阶微分方程y″+2ξωy'+ω²y＝0的解也就是轴承的故障振动模型处在二阶微分算子y″+2ξωy'+ω²y的零空间内。对于轴承故障信号模型：y(t)＝y₀e^-ξωtcos(ωt)容易证明其处在微分算子d²/dt²-2a(t)'d/a(t)dt+ω(t)²+2(a(t)'/a(t))²的零空间内。其中a(t)＝y₀e^-ξωt。

(3)输入原始信号S,终止阈值ε，λ₂ ⁰，γ⁰和λ₁ ⁰的初始值；

令j＝0，U_j＝0，λ^j ₁＝λ⁰ ₁，γ^j＝γ⁰；

(4)按下式计算

$\widehat{\mathrm{\Φ}}=-{(A^{T}A+{\mathrm{\λ}}_2{M_2}^T{M}_2)}^{-1}{A}^T{D}_2(S-U);$

(5)按下式计算参数

${\widehat{\mathrm{\λ}}}^{j+1}=\frac{1}{1+{\mathrm{\γ}}^j}\frac{S^{T}M{({{\mathrm{\λ}}_1}^j,{\mathrm{\γ}}^j,T)}^{T}S}{S^{T}M{({{\mathrm{\λ}}_1}^j,{\mathrm{\γ}}^j,T)}^{T}M({{\mathrm{\λ}}_1}^j,{\mathrm{\γ}}^j,T)S};$

(6)按下式计算

${\hat{U}}_{j+1}={(T^{T}T+(1+{\mathrm{\γ}}^j){{\mathrm{\λ}}_1}^{j+1}E)}^{-1}(T^T\mathrm{TS}+{\mathrm{\γ}}^j{{\mathrm{\λ}}_1}^{j+1}S);$

(7)按下式计算γ^j+1：

${\mathrm{\γ}}^{j+1}=\frac{(S-{\hat{U}}_{j+1})S}{{\left|\right|S-{\hat{U}}_{j+1}\left|\right|}^2}-1;$

(8)判断是否满足 $\left|\right|{\hat{U}}_{j+1}-{\hat{U}}_j\left|\right|<\mathrm{\&epsiv;}\left|\right|S\left|\right|,$ 若满足，则 $\hat{U}={\hat{U}}_{j+1},{\widehat{\mathrm{\&lambda;}}}_1={{\widehat{\mathrm{\&lambda;}}}_1}^{j+1},\widehat{\mathrm{\&gamma;}}={\mathrm{\&gamma;}}^{j+1},\widehat{\mathrm{\&alpha;}}={\widehat{\mathrm{\&alpha;}}}_j$ 输出及U＝S-R；否则令；否则令j＝j+1回到步骤(5)；

其中， $A=[A_{D_1(S-U)},A_{S-U}],$ M₁＝[D₁ ^T，E^T]^T， $M_2=\left[\begin{array}{cc}{D}_2& \\ & E\end{array}\right],$ p和q可以从 $\widehat{\mathrm{\&Phi;}}{[p^T,q^T]}^T$ 中获得， $T=D_2+B_{\widehat{\mathrm{\&Phi;}}}{M}_1,M({{\mathrm{\&lambda;}}_1}^j,{\mathrm{\&gamma;}}^j,T)={(T^{T}T+(1+\widehat{\mathrm{\&gamma;}}){\mathrm{\&lambda;}}_{1}E)}^{-1}.$

矩阵A_x是对角元素为x向量的对角矩阵。E是单位矩阵，D₁和D₂代表一阶和二阶微分矩阵，是拉格朗日参数，是保留在算子T的零空间中确定S-R的信息量的参数，代表满足条件的窄带信号，代表信号的瞬时频率。算法中的参数的初始值是根据经验设置的。由于不同的故障特征对应着算法中不同的参数，因此可以通过改变算法中的参数来实现轴承复合故障的诊断。

图2是本发明中轴承内外圈复合故障信号的染噪后的时域波形图及频谱图。采样频率为15369Hz，采样点数为4096个点，内圈特征频率为123Hz，外圈特征频率为78Hz。图中可以看出在噪声影响下冲击并不明显。

图3是本发明中轴承内外圈复合故障仿真信号的染噪后解调得到的解调谱图。可以看到既包含内圈特征频率又包含外圈特征频率。

图4是本发明中对内外圈复合故障处理后的时域图及频谱图。可以发现冲击相较图2更加明显。

图5是本发明中对内外圈复合故障处理后进行解调处理得到的解调谱。从图中可以清楚的发现内外圈故障特征频率及倍频被分别的提取出来了。

Claims (2)

1.一种基于故障特征的零空间轴承复合故障诊断方法，其特征在于：该方法包括根据轴承故障机理建立轴承故障振动模型，依据故障振动模型利用Matlab编程构造基于轴承故障特征的微分算子、对轴承故障非平稳信号进行基于故障特征的零空间微分算子的自适应分解、提取包含故障特征的局部窄带信号、解调得到故障特征；

所述分解算法包括以下步骤：

(1)根据轴承故障机理建立故障振动模型；

(2)根据故障模型构造基于轴承故障的零空间微分算子；

(3)利用构造的基于轴承故障的零空间微分算子对原始信号进行分解，通过调整算法中的参数得到包含不同故障特征的窄带信号；

(4)对包含不同故障冲击信号成分的窄带信号进行解调分析得到不同故障特征；

轴承振动信号由轴承的旋转运动引起，故障轴承振动信号中还会出现冲击和瞬态振动特征，针对这一特点及动力学可知滚动轴承故障振动形式可以抽象为质量-弹簧-阻尼系统，质量-弹簧-阻尼系统的动力学模型为：

my″+cy′+ky＝0

其中m为质量-弹簧-阻尼系统中质量元件的质量，c为质量-弹簧-阻尼系统中阻尼器的阻尼系数，k为质量-弹簧-阻尼系统中弹簧元件的弹性系数，y为系统中质量元件的位移关于时间t函数，y′为y的一阶导数，y″为y的二阶导数，为后续计算方便令n＝c/2m，ξ＝n/ω，ω为有阻尼振动角频率，ξ为阻尼比，n为一常数；

则上式可改写为如下形式：

y″+2ξωy'+ω²y＝0

上式的解为：

y(t)＝y₀e^-ξωtcos(ωt)

其中，y₀为质量元件位移函数当时间t＝0时的初始值；

即轴承故障振动模型为周期出现的指数衰减函数，可表示为：

$x\left(t\right)=\underset{-\&infin;}{\overset{+\&infin;}{\mathrm{\&Sigma;}}}{y}_0{e}^{-\mathrm{nt}}\cos \mathrm{\&omega;t\&delta;}(t-k{T}_0)$

其中δ(t)为冲击函数T₀为冲击的周期，模型的有效成分为指数衰减信号即：

y(t)＝y₀e^-ntcos(ωt)

综上可知轴承故障模型的有效成分(指数衰减函数)是二阶微分方程y″+2ξωy'+ω²y＝0的解也就是轴承的故障振动模型处在二阶微分算子y″+2ξωy'+ω²y的零空间内；对于轴承故障信号模型：y(t)＝y₀e^-ξωtcos(ωt)容易证明其处在微分算子d²/dt²-2a(t)'d/a(t)dt+ω(t)²+2(a(t)'/a(t))²的零空间内；其中a(t)＝y₀e^-ξωt；这样便可以利用上述零空间微分算子对信号进行自适应分解，由于不同类型的故障在算法中对应的参数不同，因此可以通过改变算法中的参数实现轴承的复合故障诊断。

2.根据权利要求1所述的一种基于故障特征的零空间轴承复合故障诊断方法，其特征在于：(1)根据轴承的故障机理建立轴承故障振动模型；针对信号的结构特点及动力学可知滚动轴承故障振动形式可以近似看作质量-弹簧-阻尼系统，质量-弹簧-阻尼系统的动力学模型为：

my″+cy′+ky＝0

其中m为质量-弹簧-阻尼系统中质量元件的质量，c为质量-弹簧-阻尼系统中阻尼器的阻尼系数，k为质量-弹簧-阻尼系统中弹簧元件的弹性系数，y为系统中质量元件的位移关于时间t函数，y′为y的一阶导数，y″为y的二阶导数，为后续计算方便令n＝c/2m，ξ＝n/ω，ω为有阻尼振动角频率，ξ为阻尼比，n为一常数；

则上式可改写为如下形式：

y″+2ξωy'+ω²y＝0

上式的解为：

y(t)＝y₀e^-ξωtcos(ωt)

即轴承故障振动模型为周期出现的指数衰减函数，可表示为：

$x\left(t\right)=\underset{-\&infin;}{\overset{+\&infin;}{\mathrm{\&Sigma;}}}{y}_0{e}^{-\mathrm{nt}}\cos \mathrm{\&omega;t\&delta;}(t-k{T}_0)$

其中δ(t)为冲击函数T₀为冲击的周期，模型的有效成分为指数衰减信号即：

y(t)＝y₀e^-ntcos(ωt)

(2)根据轴承故障模型构造基于轴承故障的零空间微分算子；由(1)易知轴承故障模型的有效成分(指数衰减函数)是二阶微分方程y″+2ξωy'+ω²y＝0的解也就是轴承的故障振动模型处在二阶微分算子y″+2ξωy'+ω²y的零空间内；对于轴承故障信号模型：y(t)＝y₀e^-ξωtcos(ωt)容易证明其处在微分算子d²/dt²-2a(t)'d/a(t)dt+ω(t)²+2(a(t)'/a(t))²的零空间内；其中a(t)＝y₀e^-ξωt；

(3)输入原始信号S,终止阈值ε，λ₂ ⁰，γ⁰和λ₁ ⁰的初始值；

令j＝0，U_j＝0，λ^j ₁＝λ⁰ ₁，γ^j＝γ⁰；

(4)按下式计算

$\widehat{\mathrm{\&Phi;}}=-{(A^{T}A+{\mathrm{\&lambda;}}_2{M_2}^T{M}_2)}^{-1}{A}^T{D}_2(S-U);$

(5)按下式计算参数

${\widehat{\mathrm{\&lambda;}}}^{j+1}=\frac{1}{1+{\mathrm{\&gamma;}}^j}\frac{S^{T}M{({{\mathrm{\&lambda;}}_1}^j,{\mathrm{\&gamma;}}^j,T)}^{T}S}{S^{T}M{({{\mathrm{\&lambda;}}_1}^j,{\mathrm{\&gamma;}}^j,T)}^{T}M({{\mathrm{\&lambda;}}_1}^j,{\mathrm{\&gamma;}}^j,T)S};$

(6)按下式计算

${\hat{U}}_{j+1}={(T^{T}T+(1+{\mathrm{\&gamma;}}^j){{\mathrm{\&lambda;}}_1}^{j+1}E)}^{-1}(T^T\mathrm{TS}+{\mathrm{\&gamma;}}^j{{\mathrm{\&lambda;}}_1}^{j+1}S);$

(7)按下式计算γ^j+1：

${\mathrm{\&gamma;}}^{j+1}=\frac{(S-{\hat{U}}_{j+1})S}{{\left|\right|S-{\hat{U}}_{j+1}\left|\right|}^2}-1;$

(8)判断是否满足 $\left|\right|{\hat{U}}_{j+1}-{\hat{U}}_j\left|\right|<\mathrm{\&epsiv;}\left|\right|S\left|\right|,$ 若满足，则 $\hat{U}={\hat{U}}_{j+1},{\widehat{\mathrm{\&lambda;}}}_1={{\widehat{\mathrm{\&lambda;}}}_1}^{j+1},\widehat{\mathrm{\&gamma;}}={\mathrm{\&gamma;}}^{j+1},\widehat{\mathrm{\&alpha;}}={\widehat{\mathrm{\&alpha;}}}_j$ 输出及U＝S-R；否则令；否则令j＝j+1回到步骤(5)；

其中， $A=[A_{D_1(S-U)},A_{S-U}],$ M₁＝[D₁ ^T，E^T]^T， $M_2=\left[\begin{array}{cc}{D}_2& \\ & E\end{array}\right],$ p和q可以从 $\widehat{\mathrm{\&Phi;}}{[p^T,q^T]}^T$ 中获得， $T=D_2+B_{\widehat{\mathrm{\&Phi;}}}{M}_1,M({{\mathrm{\&lambda;}}_1}^j,{\mathrm{\&gamma;}}^j,T)={(T^{T}T+(1+\widehat{\mathrm{\&gamma;}}){\mathrm{\&lambda;}}_{1}E)}^{-1};$

矩阵A_x是对角元素为x向量的对角矩阵；E是单位矩阵，D₁和D₂代表一阶和二阶微分矩阵，是拉格朗日参数，是保留在算子T的零空间中确定S-R的信息量的参数，代表满足条件的窄带信号，代表信号的瞬时频率；算法中的参数的初始值是根据经验设置的；由于不同的故障特征对应着算法中不同的参数，因此可以通过改变算法中的参数来实现轴承复合故障的诊断。