非零奇異值數(shù)量的理論分析及其在滑動(dòng)軸承-轉(zhuǎn)子振動(dòng)特征提取應(yīng)用

楊期江， 趙學(xué)智， 湯雅連， 李偉光， 滕憲斌， 郭明軍

(1. 廣州航海學(xué)院 輪機(jī)工程學(xué)院， 廣州 510725; 2． 華南理工大學(xué) 機(jī)械與汽車(chē)工程學(xué)院， 廣州 510640；3. 廣東金融學(xué)院 互聯(lián)網(wǎng)金融與信息工程學(xué)院, 廣州 510521)

奇異值分解(Singular Value Decomposition, SVD)是一種正交變換，其理論本質(zhì)是矩陣的低秩化，使得有用信息集中于非零奇異值序列前部，噪聲信息分布于非零奇異值序列后面，并取非零奇異值序列進(jìn)行矩陣重構(gòu)，濾掉噪聲信息，從而達(dá)到降噪目的。由于奇異值分解具有以上優(yōu)點(diǎn)，被廣泛應(yīng)用于小波變換結(jié)果的后續(xù)處理[1]，信號(hào)消噪[2]、滾動(dòng)軸承故障診斷[3-4]和早期故障微弱信號(hào)分離[5]等方面。

與其他降噪方法相比，SVD具有自身獨(dú)特的優(yōu)勢(shì)，其處理結(jié)果不存在相位漂移，波形失真小[6]。SVD方法不像小波分析那樣依賴(lài)于小波基函數(shù)的選擇，也不像經(jīng)驗(yàn)?zāi)B(tài)分解存在端點(diǎn)效應(yīng)問(wèn)題及模態(tài)混疊現(xiàn)象。但是在應(yīng)用這一方法時(shí)，存在一個(gè)關(guān)鍵問(wèn)題需要解決。那就是非零奇異值的選擇直接決定著消噪效果的好壞[7]。如果非零奇異值數(shù)目選擇過(guò)多，則會(huì)使得處理結(jié)果混入更多的噪聲，如果非零奇異值數(shù)目選擇過(guò)少，卻又會(huì)丟失信號(hào)中的有用頻率成分，造成波形畸變，無(wú)法準(zhǔn)確恢復(fù)我們所需要提取的特定頻率成分的信號(hào)。因此針對(duì)非零奇異值選擇這一關(guān)鍵問(wèn)題，學(xué)者們相繼提出了奇異熵[8]、差分譜理論[9]等，但這些方法并沒(méi)用解釋清楚非零奇異值個(gè)數(shù)背后所蘊(yùn)含著的真正的物理含義。以振動(dòng)信號(hào)為例，當(dāng)含有多個(gè)頻率成分的確定信號(hào)混有噪聲時(shí)，采用SVD方法進(jìn)行消噪，怎么確定非零奇異值的個(gè)數(shù)，既可以完全恢復(fù)原來(lái)含有多個(gè)頻率成分的確定信號(hào)，又能夠準(zhǔn)確地提取某個(gè)或某幾個(gè)頻率成分的信號(hào)。錢(qián)征文等[10]通過(guò)實(shí)驗(yàn)說(shuō)明了利用信號(hào)所構(gòu)造的Hankel矩陣的秩與信號(hào)中的主頻個(gè)數(shù)存在兩倍的關(guān)系，但遺憾的是沒(méi)有直接給出這一結(jié)果的理論推導(dǎo)及證明過(guò)程。趙學(xué)智等[11-12]提出了Hankel 矩陣方式下確定性信號(hào)的非零奇異值和信號(hào)所含頻率數(shù)量之間的關(guān)系，但沒(méi)有從理論上證明非零奇異值和頻率個(gè)數(shù)之間的規(guī)律。

本文通過(guò)反復(fù)仿真實(shí)驗(yàn)的方法分析了Hankel矩陣下非零奇異值數(shù)目與信號(hào)中的頻率個(gè)數(shù)成兩倍的數(shù)量關(guān)系，驗(yàn)證了奇異值成對(duì)出現(xiàn)規(guī)律，當(dāng)構(gòu)造的m×n的Hankel矩陣行數(shù)與列數(shù)充分接近時(shí)，信號(hào)中同一頻率下的兩個(gè)非零奇異值會(huì)是緊密排列在一起的。根據(jù)Hankel矩陣的構(gòu)造方式，從理論上證明了非零奇異值數(shù)目與單個(gè)頻率成2倍的關(guān)系，當(dāng)時(shí)間域信號(hào)是由多個(gè)不同頻率疊加時(shí)，非零奇異值數(shù)目始終是與頻率個(gè)數(shù)成2倍的數(shù)量關(guān)系，且非零奇異值個(gè)數(shù)與幅值和相位無(wú)關(guān)。并將這一規(guī)律應(yīng)用于旋轉(zhuǎn)機(jī)械滑動(dòng)軸承-轉(zhuǎn)子振動(dòng)的特征頻率信號(hào)提取，實(shí)現(xiàn)了對(duì)轉(zhuǎn)子不對(duì)中故障軸心軌跡的準(zhǔn)確提純。

1 Hankel矩陣下奇異值分解信號(hào)消噪理論

當(dāng)應(yīng)用SVD進(jìn)行信號(hào)處理時(shí)，由于信號(hào)是一維的，而SVD針對(duì)的是矩陣，因此需要將一維信號(hào)向量轉(zhuǎn)化成矩陣形式才能利用SVD將信號(hào)分解到不用的向量空間。因此對(duì)一個(gè)確定性信號(hào)y(t),為了對(duì)其進(jìn)行奇異值分解，首先需要利用y(t)構(gòu)造一個(gè)矩陣，這個(gè)矩陣一般為Hankel矩陣，如下式所示：

(1)

式中：y[i]為第i個(gè)采樣點(diǎn)，i=0,1,…，N-1,N是y(t)長(zhǎng)度，H∈Rm×n，矩陣行數(shù)m，列數(shù)n和信號(hào)長(zhǎng)度N之間的約束關(guān)系：m=N-n+1。對(duì)H進(jìn)行奇異值分解

H=USVT

(2)

式中:U、V是單位正交矩陣U∈Rm×m及V∈Rn×n,S=[diag(σ1,σ2,…,σq),0]或其轉(zhuǎn)置，取決于n>m或n

為了說(shuō)明奇異值降噪的原理，我們將U，V用列向量來(lái)表示，U=(u1,u2,…,um)，V=(v1,v2,…,vn)，則式(2)可以改寫(xiě)為列向量ui和vi的表現(xiàn)形式。

(3)

式中:ui∈Rm×1，vi∈Rn×1，i=1,2,…q,q=min(m,n)，式中共有q個(gè)奇異值，但這些奇異值并不都是有效的，對(duì)于確定性的無(wú)噪聲信號(hào)x(t)，它的非零奇異值才是有效奇異值，但是非零奇異值數(shù)量不一定就等于q，因?yàn)槿绻麤](méi)有噪聲，這q個(gè)奇異值中將會(huì)有零奇異值存在；而對(duì)于含噪聲的信號(hào)，在其q個(gè)奇異值中，有很多非零奇異值是由噪聲產(chǎn)生的，因此其有效奇異值的數(shù)目也不等于q。對(duì)于一個(gè)具有確定頻率個(gè)數(shù)的信號(hào)，不管它含噪與否，它的非零奇異值數(shù)目是與確定的頻率個(gè)數(shù)相關(guān)的。近來(lái)有學(xué)者針對(duì)這一問(wèn)題進(jìn)行了實(shí)驗(yàn)上的探索，趙學(xué)智等首先發(fā)現(xiàn)了非零奇異值數(shù)量與信號(hào)中的頻率個(gè)數(shù)成兩倍這一關(guān)系，但沒(méi)有從理論上證明這一結(jié)果。作者經(jīng)過(guò)研究，發(fā)現(xiàn)可以從理論上對(duì)奇異值的數(shù)量規(guī)律進(jìn)行證明。下面先從實(shí)驗(yàn)角度分析奇異值的數(shù)量規(guī)律，然后給出理論證明過(guò)程。

2 Hankel矩陣下非零奇異值數(shù)量規(guī)律的實(shí)驗(yàn)分析

先通過(guò)實(shí)驗(yàn)分析非零奇異值數(shù)量與頻率個(gè)數(shù)之間的關(guān)系，采用反復(fù)實(shí)驗(yàn)的方法。具體為：從含有一個(gè)頻率的信號(hào)出發(fā)，逐步增加信號(hào)中頻率的個(gè)數(shù)，利用信號(hào)構(gòu)造相同維數(shù)的 Hankel 矩陣進(jìn)行 SVD 處理，觀察并分析非零奇異值數(shù)目的變化。采用 8個(gè)信號(hào)進(jìn)行試驗(yàn)，下一個(gè)信號(hào)在前一個(gè)信號(hào)的基礎(chǔ)上增加一個(gè)頻率，且采樣頻率相同，隨著頻率值逐漸增大，其幅值也逐漸加大，而相位隨機(jī)，實(shí)驗(yàn)的8個(gè)信號(hào)如式(4)所示：

y1(t)=2.8sin(t),

y2(t)=2.8sin(t)+4.6sin(4t+0.56),

y3(t)=2.8sin(t)+4.6sin(4t+0.56)+

6.6sin(12t+1.35)

y4(t)=2.8sin(t)+4.6sin(4t+0.56)+

6.6sin(12t+1.35)+8.4sin(20t+2.58)

y5(t)=2.8sin(t)+4.6sin(4t+0.56)+

6.6sin(12t+1.35)+8.4sin(20t+2.58)+

10.6sin(30t+3.12)

y6(t)=2.8sin(t)+4.6sin(4t+0.56)+

6.6sin(12t+1.35)+8.4sin(20t+2.58)+

10.6sin(30t+3.12)+14sin(60t+2.79)

y7(t)=2.8sin(t)+4.6sin(4t+0.56)+

6.6sin(12t+1.35)+8.4sin(20t+2.58)+

10.6sin(30t+3.12)+14sin(60t+2.79)+

17sin(86t+1.29)

y8(t)=2.8sin(t)+4.6sin(4t+0.56)+

6.6sin(12t+1.35)+8.4sin(20t+2.58)+

10.6sin(30t+3.12)+14sin(60t+2.79)+

17sin(86t+1.29)+21.4sin(104t+3.12)

(4)

設(shè)每個(gè)信號(hào)在[0, 2π]內(nèi)均采集1 024點(diǎn)，利用采樣點(diǎn)構(gòu)造512×513的大小的Hankel矩陣，并計(jì)算其奇異值，得到各個(gè)信號(hào)的奇異值結(jié)果具體如圖1所示，圖中從左至右依次為y1(t)～y8(t)的奇異值曲線(xiàn)。從圖 1可見(jiàn)，當(dāng)信號(hào)每增加一個(gè)頻率成分，不管這個(gè)幅值、相位如何，非零奇異值的數(shù)目總是只增加兩個(gè)。在下一個(gè)信號(hào)的奇異值曲線(xiàn)中，除了新增加的兩個(gè)奇異值較大外，其余的奇異值均和上一個(gè)信號(hào)的奇異值大小幾乎相等，而且新增的兩個(gè)奇異值之間幾乎相等。采用反復(fù)實(shí)驗(yàn)方法，在仿真信號(hào)(4)中，如果采用其它的頻率和幅值，也會(huì)得到同樣結(jié)果：即不管頻率大小及其幅值如何，一個(gè)頻率總是只產(chǎn)生兩個(gè)非零奇異值。

為了分析新增的兩個(gè)奇異值之間幾乎兩兩相等的原因。我們構(gòu)造一個(gè)含有4個(gè)頻率的信號(hào)，具體如式(5)所示：

y(t)=sin(5t)+2.6sin(9t+0.35)+

3.7sin(19t+1.26)+5.4sin(32t+0.98)

(5)

對(duì)此信號(hào)在[0, 2π]內(nèi)均采集1 024點(diǎn)，分別構(gòu)造不同維數(shù)大小的Hankel 矩陣，并計(jì)算其奇異值。在不同維數(shù)的Hankel矩陣下各個(gè)信號(hào)的奇異值結(jié)果具體如圖2所示。從圖 2可見(jiàn)，隨著構(gòu)造的m×n的Hankel矩陣維數(shù)增加，奇異值大小都有不同程度的增大。當(dāng)行數(shù)與列數(shù)充分接近時(shí)，會(huì)出現(xiàn)信號(hào)中同一頻率下的兩個(gè)奇異值幾乎相等。

圖2 非零奇異值與矩陣維數(shù)之間的關(guān)系

通過(guò)上述方法的反復(fù)實(shí)驗(yàn)，我們發(fā)現(xiàn)了非零奇異值數(shù)目與頻率個(gè)數(shù)之間的數(shù)量關(guān)系規(guī)律以及奇異值成對(duì)出現(xiàn)規(guī)律：① 對(duì)于一個(gè)含有固定頻率數(shù)量的確定性信號(hào)，利用其構(gòu)造m×n的 Hankel 矩陣，當(dāng)m和n中最小值大于信號(hào)中頻率數(shù)量的兩倍之后，非零奇異值的數(shù)目始終為信號(hào)中頻率數(shù)量的兩倍，稱(chēng)之為奇異值數(shù)量規(guī)律。② 當(dāng)構(gòu)造的m×n的Hankel矩陣行數(shù)與列數(shù)充分接近時(shí)，會(huì)出現(xiàn)信號(hào)中同一頻率下的兩個(gè)非零奇異值是緊密排列在一起的，即非零奇異值成對(duì)出現(xiàn)規(guī)律。Hankel矩陣下非零奇異值與信號(hào)中頻率成分個(gè)數(shù)之間呈2倍的奇妙的數(shù)量關(guān)系規(guī)律。通過(guò)研究，我們發(fā)現(xiàn)這一現(xiàn)象背后存在理論上的必然性，下面從理論上來(lái)證明這一規(guī)律。

3 非零奇異值數(shù)量規(guī)律的理論證明

設(shè)有時(shí)間域連續(xù)的確定性信號(hào)y(t)包含I個(gè)頻率成分，由下式表示：

(6)

以采樣間隔為T(mén)s，采樣點(diǎn)數(shù)為j，那么時(shí)間離散為t=Tsj，因此將原始連續(xù)信號(hào)離散為y[j]。根據(jù)y[j]構(gòu)建一個(gè)m×n的Hankel矩陣H，顯然H可以表示為單個(gè)頻率分量構(gòu)造的Hankel矩陣Hi之和[1]，即：

(7)

式中Hi是頻率分量yi[j]=aisin(ωiTs·j+φi)構(gòu)造的m×n的Hankel矩陣，如下式所示[2]：

(8)

式中:m=N-n+1，而N是xi(t)離散后的長(zhǎng)度。

現(xiàn)在分三種情況分別對(duì)Hankel矩陣H進(jìn)行討論，第一種情況為確定性信號(hào)y(t)僅含有1個(gè)頻率成分，第二種情況為確定性信號(hào)y(t)含有I個(gè)頻率成分，第三種情況為確定性信號(hào)y(t)含有噪聲時(shí)。

(1) 當(dāng)確定性信號(hào)y(t)僅含有1個(gè)頻率成分時(shí)

此時(shí)原始信號(hào)y(t)僅表示為：

y(t)=A1sin(ω1t+φ1)

(9)

那么其離散之后構(gòu)建的m×n的Hankel矩陣為：

(10)

根據(jù)三角函數(shù)關(guān)系式，可得：

sin(ω1t+φ1)=a1sin(ω1t)+b1cos(ω1t)

(11)

將式(11)代入式(10)，可得：

(12)

由歐拉公式ejω1Ts=cosω1Ts+jsinω1Ts，得：

(13)

將式(13)代入式(12)中得：

(14)

通過(guò)化簡(jiǎn)為：

(15)

可將H1改寫(xiě)為兩個(gè)Hankel矩陣之和的形式：

(16)

其中：

(2) 當(dāng)確定性信號(hào)y(t)含有I個(gè)頻率成分時(shí)

此時(shí)原始信號(hào)y(t)表示含單個(gè)頻率時(shí)域信號(hào)之和，即為：

(17)

那么其離散之后構(gòu)建的m×n的Hankel矩陣H為單個(gè)頻率分量構(gòu)造的Hankel矩陣Hi之和：

(18)

同理，將Hi改寫(xiě)為兩個(gè)Hankel矩陣之和的形式：

(19)

式中：

將式(19)代入到式(18)中，可得含有I個(gè)頻率成分信號(hào)總Hankel矩陣可寫(xiě)成如下形式：

(20)

由式(20)可知，矩陣H由I個(gè)Hankel矩陣Hi求和得到，其中第i個(gè)頻率成分信號(hào)的Hankel矩陣Hi的秩為2，僅與頻率ωi有關(guān)，與幅值A(chǔ)i，相位角φi都無(wú)關(guān)。由于ω1≠ω2≠ω3≠…≠ωI，可知每一個(gè)頻率信號(hào)所構(gòu)成的Hankel矩陣Hi之間是互相線(xiàn)性無(wú)關(guān)的，因此可知當(dāng)原始信號(hào)含有I個(gè)頻率成分時(shí)，總矩陣H的秩為2I，當(dāng)且僅當(dāng)m行×n列的總矩陣H滿(mǎn)足min(m,n)≥2I時(shí)，總矩陣H秩才為2I。由非零奇異值數(shù)量等于矩陣的秩，因此I個(gè)頻率成分信號(hào)構(gòu)造的m行×n列Hankel矩陣滿(mǎn)足min(m,n)≥2I時(shí)，通過(guò)SVD分解之后的非零奇異值個(gè)數(shù)也應(yīng)該為2I個(gè)，且非零奇異值個(gè)數(shù)與幅值和相位無(wú)關(guān)。

(3) 當(dāng)確定性信號(hào)y(t)中含有噪聲時(shí)

噪聲的頻譜是分布在整個(gè)頻率軸上的，因此噪聲將產(chǎn)生很多非零奇異值，但是噪聲產(chǎn)生的非零奇異值只是影響確定性信號(hào)的非零奇異值大小，而不會(huì)影響確定性信號(hào)的非零奇異值個(gè)數(shù)。而且非零奇異值的大小是由頻率的幅值決定的，某頻率的幅值越大，其對(duì)應(yīng)的非零奇異值也越大，因?yàn)槠娈愔敌蛄惺墙敌蚺帕械?，如果確定性信號(hào)的幅值較大，經(jīng)過(guò)SVD處理后，噪聲產(chǎn)生的非零奇異值就會(huì)排列在確定性信號(hào)的非零奇異值后面，因此，根據(jù)確定性信號(hào)中的頻率數(shù)量，選擇前面合適的奇異值進(jìn)行SVD重構(gòu)，就可以消除噪聲的影響，得到干凈的信號(hào)。

綜上理論分析可得到如下重要結(jié)論：采用Hankel矩陣構(gòu)造一維信號(hào)，如果一維信號(hào)中所含互不相等的頻率成分個(gè)數(shù)為I，當(dāng)矩陣當(dāng)中的行數(shù)與列數(shù)的最小值大于頻率成分個(gè)數(shù)的兩倍時(shí)，通過(guò)SVD之后的非零奇異值個(gè)數(shù)應(yīng)該為2I個(gè)，可表述為每?jī)蓚€(gè)非奇異值對(duì)應(yīng)信號(hào)中的一個(gè)頻率成分；且非零奇異值個(gè)數(shù)是與幅值和相位無(wú)關(guān)的，即非零奇異值數(shù)量規(guī)律。

對(duì)一維信號(hào)所構(gòu)造的Hankle矩陣進(jìn)行奇異值分解，得到從大到小依次排列的非零奇異值序列譜。根據(jù)以上的實(shí)驗(yàn)及理論分析可知，振動(dòng)信號(hào)的特征頻率成分就隱含在這一系列的非零奇異值序列譜當(dāng)中。下面通過(guò)非零奇異值和頻率的關(guān)系來(lái)提取滑動(dòng)軸承-轉(zhuǎn)子振動(dòng)信號(hào)特征頻率。

4 對(duì)滑動(dòng)軸承-轉(zhuǎn)子振動(dòng)信號(hào)特征的提取

滑動(dòng)軸承是大重型旋轉(zhuǎn)機(jī)械的核心支撐部件。軸承間隙中所形成的楔形油膜支撐轉(zhuǎn)子高速旋轉(zhuǎn)。相對(duì)于滾動(dòng)軸承，滑動(dòng)軸承與轉(zhuǎn)子之間存在便于形成油膜的間隙，大型機(jī)組在進(jìn)行安裝調(diào)試時(shí)，聯(lián)軸器與轉(zhuǎn)子之間的對(duì)中精度往往是制約轉(zhuǎn)子高速運(yùn)行的關(guān)鍵因素。因此機(jī)組在進(jìn)行現(xiàn)場(chǎng)安裝前，需要利用實(shí)驗(yàn)臺(tái)調(diào)試模擬實(shí)際安裝環(huán)境，通過(guò)對(duì)滑動(dòng)軸承-轉(zhuǎn)子進(jìn)行振動(dòng)測(cè)試提取軸心軌跡等特征信號(hào)以判斷滑動(dòng)軸承-轉(zhuǎn)子的運(yùn)行狀況。實(shí)際測(cè)試過(guò)程中，機(jī)組試驗(yàn)臺(tái)的振動(dòng)受到各個(gè)部件耦合振動(dòng)的影響，所測(cè)試得到的實(shí)驗(yàn)數(shù)據(jù)噪聲干擾較大，導(dǎo)致軸心軌跡等特征信號(hào)模糊不清，難以辨別。

(a) 試驗(yàn)臺(tái)架實(shí)物圖

(b) 電渦流傳感器安裝圖

本文試驗(yàn)研究采用的軸承-轉(zhuǎn)子試驗(yàn)臺(tái)為團(tuán)隊(duì)自主研發(fā)搭建的大型滑動(dòng)軸承轉(zhuǎn)子試驗(yàn)臺(tái)系統(tǒng)。最高轉(zhuǎn)速5 000 r/min；轉(zhuǎn)子重量約5T，具體如圖3(a)所示。圖3(b)是實(shí)際機(jī)組試驗(yàn)臺(tái)電渦流位移傳感器安裝圖。本文轉(zhuǎn)子試驗(yàn)轉(zhuǎn)速為2 600 r/min，圖4是滑動(dòng)軸承-轉(zhuǎn)子近電機(jī)端實(shí)際運(yùn)行過(guò)程中采集到的軸心軌跡，通過(guò)該軸心軌跡難以判斷機(jī)組試驗(yàn)臺(tái)運(yùn)行狀態(tài)。

圖4 滑動(dòng)軸承-轉(zhuǎn)子振動(dòng)軸心軌跡

為了更進(jìn)一步地分析影響軸心軌跡的各個(gè)頻率成分，首先分別對(duì)近電機(jī)端1號(hào)測(cè)點(diǎn)(x方向)與2號(hào)測(cè)點(diǎn)(y方向)上采集得到位移振動(dòng)原始信號(hào)進(jìn)行快速傅里葉變化(FFT)，可得到其頻率幅值譜如圖5所示。從圖5(b)，(d)頻域圖可知，1倍頻與2倍頻幅值最為突出，且2倍頻與1倍頻幅值接近，同時(shí)存在高倍轉(zhuǎn)頻，根據(jù)轉(zhuǎn)子振動(dòng)理論，二倍轉(zhuǎn)頻的存在并伴有高倍轉(zhuǎn)頻這種現(xiàn)象主要是由于轉(zhuǎn)子對(duì)中不良所致，其中還有明顯的電力系統(tǒng)的工頻干擾，分別是50 Hz與150 Hz。如果只通過(guò)單個(gè)測(cè)點(diǎn)信號(hào)的頻譜分析難以判別轉(zhuǎn)子不對(duì)中的程度，需要引入軸心軌跡以進(jìn)一步識(shí)別，而由于高倍轉(zhuǎn)頻以及電力工頻等頻率成分疊加導(dǎo)致原始信號(hào)所繪制的軸心軌跡模糊，難以辨識(shí)。

根據(jù)Hankel矩陣下非零奇異值數(shù)目與頻率個(gè)數(shù)之間的數(shù)量關(guān)系，以及奇異值成對(duì)出現(xiàn)規(guī)律，對(duì)原始信號(hào)1倍頻與2倍頻成分進(jìn)行特征提取。原信號(hào)的數(shù)據(jù)長(zhǎng)度為4 096，故利用圖5(a)、(c)的時(shí)域信號(hào)構(gòu)造兩個(gè)2 048×2 049的Hankel矩陣并對(duì)其進(jìn)行奇異值分解，得到兩個(gè)原信號(hào)的奇異值如圖6所示(圖中只列出前40個(gè)奇異值)。

原信號(hào)中主要有1倍頻與2倍頻這兩個(gè)幅值較大的頻率成分，它們產(chǎn)生4個(gè)奇異值，而由于頻率的幅值對(duì)奇異值的影響最大[11]，因此這2個(gè)幅值較大的頻率對(duì)應(yīng)的很可能就是前4個(gè)奇異值，現(xiàn)在利用第1對(duì)和第2對(duì)奇異值進(jìn)行 SVD 重構(gòu)，即有重構(gòu)矩陣為[12]：

(21)

(a) 近電機(jī)端1號(hào)測(cè)點(diǎn)時(shí)域圖

(b) 近電機(jī)端1號(hào)測(cè)點(diǎn)頻域圖

(c) 近電機(jī)端2號(hào)測(cè)點(diǎn)時(shí)域圖

(d) 近電機(jī)端2號(hào)測(cè)點(diǎn)頻域圖

圖5 滑動(dòng)軸承轉(zhuǎn)子振動(dòng)信號(hào)及其幅值譜

Fig.5 Vibration signal of sliding bearing-rotor and its spectrum

(a) 近電機(jī)端1號(hào)測(cè)點(diǎn)

(b) 近電機(jī)端2號(hào)測(cè)點(diǎn)

得到的矩陣H′也是一個(gè)2 048×2 049的矩陣，采用平均法[12]從H′中恢復(fù)出信號(hào)，結(jié)果如圖7(b)、(d)所示。從幅值譜可見(jiàn)這兩個(gè)重構(gòu)信號(hào)的頻率成分確實(shí)為轉(zhuǎn)頻和二倍轉(zhuǎn)頻，它們的幅值分別為1.542×10-5m和1.586×10-5m、1.266×10-5m和2.152×10-5m，而在兩個(gè)原信號(hào)的幅值譜中 (圖5(b)、(d)所示)，這兩個(gè)頻率的幅值分別為1.539×10-5m和1.594×10-5m、1.269×10-5m和2.15×10-5m，SVD提取的結(jié)果與它們只存在微小的差異，這種差異是噪聲的影響造成的，是不可避免的，可以說(shuō)，SVD把轉(zhuǎn)頻和二倍轉(zhuǎn)頻基本完整地提取出來(lái)了。從提取出的時(shí)域波形(圖7(a)、(c)所示)可以看出1倍頻與二倍頻所對(duì)應(yīng)的時(shí)域波形都沒(méi)有發(fā)生較大的波動(dòng)，因此轉(zhuǎn)子不對(duì)中狀態(tài)是比較穩(wěn)定的。從提取出的頻域波形(圖7(b)、(d)所示)可以看出該頻域圖是非常干凈的，只有1倍轉(zhuǎn)頻和2倍頻兩個(gè)譜峰存在，尤其令人注意的是，不但這兩個(gè)譜峰之間的電力系統(tǒng)50 Hz的工頻干擾也完全消失了，而且達(dá)到了對(duì)高次諧波干擾完美的濾波效果，這種效果是一般帶通濾波器達(dá)不到的?？梢钥闯鯯VD分離特定頻率成分的效果是非常優(yōu)異的，一般的濾波器很難做到這一點(diǎn)，濾波器總是會(huì)存在一定的帶寬，并且還有過(guò)渡帶的影響，除非兩個(gè)頻率隔得特別遠(yuǎn)，否則濾波器在分離單個(gè)頻率時(shí)總是會(huì)受到臨近頻率的影響。

(a) 近電機(jī)端1號(hào)測(cè)點(diǎn)時(shí)域圖

(b) 近電機(jī)端1號(hào)測(cè)點(diǎn)頻域圖

(c) 近電機(jī)端2號(hào)測(cè)點(diǎn)時(shí)域圖

(d) 近電機(jī)端2號(hào)測(cè)點(diǎn)頻域圖

圖7 SVD提取的故障波形特征及其幅值譜

Fig.7 Waveform feature of fault extracted by SVD and its spectrum

分別以圖7(a)、(c)時(shí)域波形為x軸、y軸得到滑動(dòng)軸承-轉(zhuǎn)子近電機(jī)端振動(dòng)軸心軌跡，具體如圖8(a)所示。同理，采用上述方法可以得到滑動(dòng)軸承-轉(zhuǎn)子遠(yuǎn)電機(jī)端振動(dòng)軸心軌跡如圖8(b)所示。根據(jù)轉(zhuǎn)子不對(duì)中故障的特征分析，當(dāng)不對(duì)中比較輕微時(shí)，軸心軌跡呈橢圓狀；當(dāng)不對(duì)中故障達(dá)到中等程度時(shí)，軸心軌跡呈香蕉形；當(dāng)不對(duì)中故障較嚴(yán)重時(shí)，軸心軌跡呈外“8”字形。從圖8可以看出，滑動(dòng)軸承-轉(zhuǎn)子近電機(jī)端振動(dòng)軸心軌跡呈香蕉狀，不對(duì)中故障為中等程度；滑動(dòng)軸承-轉(zhuǎn)子遠(yuǎn)電機(jī)端振動(dòng)軸心軌跡呈外“8”字形，不對(duì)中故障較嚴(yán)重，所以遠(yuǎn)電機(jī)端的不對(duì)中程度要大于近電機(jī)端。

(a) 近電機(jī)端

(b) 遠(yuǎn)電機(jī)端

經(jīng)驗(yàn)?zāi)B(tài)分解(Empirical Mode Decomposition，EMD)是一種廣泛應(yīng)用的信號(hào)分解方法，也可以分離出信號(hào)中的各個(gè)頻率。近年來(lái)發(fā)展了EEMD(總體平均經(jīng)驗(yàn)?zāi)B(tài)分解)、CEEMD(補(bǔ)充的 EEMD)等不同的改進(jìn)方法[13]。對(duì)于本文的滑動(dòng)軸承-轉(zhuǎn)子振動(dòng)信號(hào)，利用CEEMD方法自適應(yīng)分解得到一系列從高頻到低頻的具有不同物理意義內(nèi)稟模態(tài)分量(Intrinsic Mode Function, IMF)，具體如圖9(a)、(c)所示。通過(guò)快速傅里葉變換(FFT)得到各IMF的幅值譜，具體如圖9(b)、(d)所示。從圖9(b)可知近電機(jī)端1號(hào)測(cè)點(diǎn)信號(hào)經(jīng)CEEMD方法分解之后得到的IMF3的頻譜圖中主要是以2倍轉(zhuǎn)頻為主，但是沒(méi)有能夠分離電力系統(tǒng)的工頻干擾成分，如1倍工頻50 Hz，3倍工頻150 Hz；IMF4主要是以1倍轉(zhuǎn)頻為主，但是仍然存在模態(tài)混疊現(xiàn)象，如2倍頻的存在。從圖9(d)可知近電機(jī)端2號(hào)測(cè)點(diǎn)信號(hào)CEEMD分解之后得到的IMF3的頻譜圖中主要是以2倍轉(zhuǎn)頻為主，但同樣還伴有電力系統(tǒng)的工頻干擾，如1倍工頻50 Hz，3倍工頻150 Hz，同時(shí)仍然存在模態(tài)混疊現(xiàn)象，如1倍頻的存在；IMF4主要是以1倍轉(zhuǎn)頻為主。因此利用內(nèi)稟模態(tài)分量IMF3與IMF4進(jìn)行CEEMD重構(gòu)，即有重構(gòu)信號(hào)為：

(22)

式中：ci為第i個(gè)IMF。

可得到近電機(jī)端1號(hào)、2號(hào)測(cè)點(diǎn)提取的故障波形，分別以這兩個(gè)重構(gòu)之后的時(shí)域波形為x軸、y軸得到滑動(dòng)軸承-轉(zhuǎn)子近電機(jī)端振動(dòng)軸心軌跡，具體如圖10(a)所示。同理，采用上述方法可以得到滑動(dòng)軸承-轉(zhuǎn)子遠(yuǎn)電機(jī)端振動(dòng)軸心軌跡如圖10(b)所示。與圖8對(duì)比可見(jiàn)，SVD方法提取的軸心軌跡比EMD方法提取的軸心軌跡清楚很多。

5 結(jié) 論

(1) 從實(shí)驗(yàn)中發(fā)現(xiàn)了Hankel矩陣下非零奇異值成對(duì)出現(xiàn)規(guī)律：當(dāng)構(gòu)造的m×n的Hankel矩陣行數(shù)與列數(shù)充分接近時(shí)，會(huì)出現(xiàn)信號(hào)中同一頻率下的兩個(gè)非零奇異值是緊密排列在一起的。

(2) 從理論上證明了Hankel矩陣下非零奇異值數(shù)量規(guī)律：對(duì)于一個(gè)含有固定頻率數(shù)量的信號(hào)，利用其構(gòu)造m×n的 Hankel 矩陣，當(dāng)矩陣維數(shù)大于信號(hào)中頻率數(shù)量的兩倍之后，非零奇異值的數(shù)目與信號(hào)中頻率個(gè)數(shù)始終成兩倍的數(shù)量關(guān)系。且非零奇異值個(gè)數(shù)是與幅值和相位無(wú)關(guān)的。

(3) SVD分離特定頻率成分的效果是非常明顯的，利用Hankel矩陣下非零奇異值數(shù)量規(guī)律，可實(shí)現(xiàn)滑動(dòng)軸承-轉(zhuǎn)子振動(dòng)不對(duì)中故障1倍轉(zhuǎn)頻與2倍頻的特征提取，實(shí)現(xiàn)了對(duì)轉(zhuǎn)子不對(duì)中故障軸心軌跡的準(zhǔn)確提純。

(a) 近電機(jī)端1號(hào)測(cè)點(diǎn)CEEMD固有模態(tài)分量

(b) 近電機(jī)端1號(hào)測(cè)點(diǎn)CEEMD固有模態(tài)分量幅值譜

(c) 近電機(jī)端2號(hào)測(cè)點(diǎn)CEEMD固有模態(tài)分量

(d) 近電機(jī)端2號(hào)測(cè)點(diǎn)CEEMD固有模態(tài)分量幅值譜

圖9 CEEMD提取的故障波形特征及其幅值譜

Fig.9 Waveform feature of fault extracted by CEEMD and its spectrum

(a) 近電機(jī)端

(b) 遠(yuǎn)電機(jī)端

圖10 CEEMD提純后的軸心軌跡

Fig.10 Sliding bearing-rotor vibration orbit extracted by CEEMD