張銳戈,李 勁

(1. 三明學(xué)院 機(jī)電工程學(xué)院 福建 三明 365004； 2. 裝備智能福建省高校重點(diǎn)實(shí)驗(yàn)室 福建 三明 365004； 3. 湖南科技大學(xué) 信息與電氣工程學(xué)院 湖南 湘潭 411201)

滾動(dòng)軸承雙譜特征提取及變工況故障診斷

張銳戈1，2,李 勁3

(1. 三明學(xué)院 機(jī)電工程學(xué)院 福建 三明 365004； 2. 裝備智能福建省高校重點(diǎn)實(shí)驗(yàn)室 福建 三明 365004； 3. 湖南科技大學(xué) 信息與電氣工程學(xué)院 湖南 湘潭 411201)

研究工況變化情形下的信號特征提取技術(shù)和軸承狀態(tài)判別方法，以滿足旋轉(zhuǎn)裝備變工況運(yùn)行環(huán)境故障診斷需要。首先分析滾動(dòng)軸承振動(dòng)加速度信號雙譜特征，發(fā)現(xiàn)其在不同故障類型時(shí)差異性顯著，在故障相同但工況不同時(shí)又較為相似，適合提取工況魯棒特征參數(shù)。利用雙譜對稱特性降低數(shù)據(jù)維數(shù)，將其劃分成八個(gè)子區(qū)域后去除冗余數(shù)據(jù)，獲得數(shù)據(jù)量僅為原始雙譜四分之一的優(yōu)選雙譜。再將優(yōu)選雙譜做奇異值分解，依奇異值累積貢獻(xiàn)率選擇相應(yīng)的左奇異向量，計(jì)算元素模值后構(gòu)造特征向量。最后用連續(xù)型隱馬爾可夫模型辨識故障，在輕度、中度和重度故障程度下，使用任意一種工況的數(shù)據(jù)訓(xùn)練診斷模型，都能有效判別四種工況的軸承狀態(tài)，實(shí)現(xiàn)了變載荷、變轉(zhuǎn)速環(huán)境的滾動(dòng)軸承智能故障診斷。

故障診斷;變工況;雙譜;奇異值分解;隱馬爾可夫模型

1 引言

滾動(dòng)軸承是旋轉(zhuǎn)機(jī)械的重要部件，其狀態(tài)監(jiān)測和故障診斷在預(yù)防設(shè)備意外停機(jī)及降低維護(hù)成本等方面具有重要意義[1]。工程實(shí)踐中，旋轉(zhuǎn)機(jī)械運(yùn)行工況常常發(fā)生變化，如數(shù)控機(jī)床為抑制顫振使用變轉(zhuǎn)速工作模式[2]，風(fēng)電設(shè)備轉(zhuǎn)速及載荷隨風(fēng)向和風(fēng)速變化[3]等。變工況情形下，滾動(dòng)軸承振動(dòng)特性比恒定工況時(shí)復(fù)雜，軸承狀態(tài)與信號特征間的對應(yīng)關(guān)系也發(fā)生變化[4]，常用的恒定工況診斷方法難以滿足變工況故障診斷需要。

軸承智能故障診斷是一基于特征參量的模式識別過程。變工況故障診斷，需從工況擾動(dòng)的振動(dòng)信號中提取工況魯棒特征參數(shù)，再用合適的模式識別方法辨識軸承狀態(tài)。旋轉(zhuǎn)機(jī)械變工況運(yùn)行，主要涉及變轉(zhuǎn)速、變載荷、轉(zhuǎn)速和載荷同時(shí)變化三種情形，分別采用不同分析手段或特征提取策略處理信號。

變轉(zhuǎn)速情形下，需消除轉(zhuǎn)速變化對特征參數(shù)影響，有階次跟蹤、階次循環(huán)平穩(wěn)和階次時(shí)頻分析等方法：

(1) 階次跟蹤通過信號重采樣技術(shù)，將時(shí)間域的非平穩(wěn)信號轉(zhuǎn)變?yōu)榻嵌扔虻钠椒€(wěn)信號，進(jìn)而恢復(fù)變轉(zhuǎn)速信號和轉(zhuǎn)速之間的對應(yīng)關(guān)系。西安交通大學(xué)Yi Wang等人先用低通濾波分離出振動(dòng)信號中的轉(zhuǎn)速成分，再用脊線檢測技術(shù)從包絡(luò)信號的短時(shí)譜中提取瞬時(shí)轉(zhuǎn)速信息[5]；北京交通大學(xué)程衛(wèi)東等人，研究了變轉(zhuǎn)速信號等角度重采后的包絡(luò)變形情況，提出修正角度變換值消除包絡(luò)變形的方法[6]。

(2) 階次循環(huán)平穩(wěn)將階次跟蹤和時(shí)域平均技術(shù)相結(jié)合：能消除轉(zhuǎn)速波動(dòng)對振動(dòng)信號影響，還通過時(shí)域平均抑制噪聲和非同步周期成分[7]。印度工學(xué)院的C. Mishra等人先將振動(dòng)信號進(jìn)行短期角度同步平均，再結(jié)合小波去噪和包絡(luò)階次跟蹤提取故障特征信息[8]。

(3) 階次時(shí)頻分析在時(shí)頻域使用階次跟蹤技術(shù)，能提取非平穩(wěn)信號中的瞬態(tài)時(shí)頻信息。北京交通大學(xué)王天楊博士使用包絡(luò)信號時(shí)頻譜提取等效轉(zhuǎn)頻趨勢線[9]；湖南大學(xué)徐亞軍等人先用線性調(diào)頻小波獲得轉(zhuǎn)速特征頻率[10]，再對信號包絡(luò)重采樣獲得角度域平穩(wěn)信號；西安交通大學(xué)Ming Zhao等人將經(jīng)驗(yàn)?zāi)B(tài)分解和包絡(luò)階次分析相結(jié)合，實(shí)現(xiàn)了復(fù)雜工況下的滾動(dòng)軸承多損傷檢測[11]。

此外，也有學(xué)者開展無需預(yù)濾波和信號重采樣的變轉(zhuǎn)速振動(dòng)信號處理研究。加拿大渥太華大學(xué)的Juanjuan Shi等人提出一種加窗分形維數(shù)變換和廣義解調(diào)的方法，克服了階次方法的重采樣誤差、轉(zhuǎn)頻和瞬時(shí)頻率耦合等問題[12]；國防科技大學(xué)范彬博士先用相空間法構(gòu)建跟蹤矩陣，再用平滑正交分解方法分離軸承損傷和工況變化形成的演化趨勢[13]。

變載荷故障診斷通常選擇載荷變化不敏感的特征值，再用合適的預(yù)處理方法提取載荷魯棒特征參數(shù)。西安電子科技大學(xué)張銳戈博士，發(fā)現(xiàn)振動(dòng)信號短時(shí)過零率對載荷變化不敏感，結(jié)合小波降噪實(shí)現(xiàn)了載荷魯棒的滾動(dòng)軸承早期故障診斷[14]；周勇等人將小波包節(jié)點(diǎn)能量作為特征參數(shù)，經(jīng)獨(dú)立特征選擇后用相關(guān)向量機(jī)實(shí)現(xiàn)變載荷故障診斷[15]。

轉(zhuǎn)速、載荷同時(shí)變化情形，需綜合考慮二者變化對特征參數(shù)的擾動(dòng)。華北電力大學(xué)王珊等人先用階比重采樣技術(shù)，將時(shí)間域非平穩(wěn)信號轉(zhuǎn)變?yōu)榻嵌扔蚱椒€(wěn)信號，再用經(jīng)驗(yàn)?zāi)B(tài)分解方法抽取本征模態(tài)函數(shù)重構(gòu)信號，實(shí)現(xiàn)了風(fēng)電軸承轉(zhuǎn)速、載荷變化情形下的早期故障診斷。華東交通大學(xué)宋京偉等人用小波包分解從振動(dòng)信號時(shí)頻域提取特征參量，結(jié)合連續(xù)型隱馬爾可夫模型實(shí)現(xiàn)變載荷及未知轉(zhuǎn)速下的故障診斷[17]。

電機(jī)是旋轉(zhuǎn)機(jī)械的重要?jiǎng)恿碓?，滾動(dòng)軸承又為電機(jī)關(guān)鍵零部件，工作狀態(tài)會(huì)影響設(shè)備整體性能。針對電機(jī)載荷變化、轉(zhuǎn)速波動(dòng)運(yùn)行環(huán)境，基于雙譜分析提出滾動(dòng)軸承變工況智能故障診斷新方法。與常用的階次方法相比，新方法能在轉(zhuǎn)速、載荷同時(shí)發(fā)生變化情形下實(shí)現(xiàn)變工況故障診斷，且無需進(jìn)行數(shù)據(jù)重采樣和轉(zhuǎn)速信息提取，克服了角度域重采樣誤差和信號成分抽取方法的瞬時(shí)頻率干擾等問題。

利用雙譜分析適合處理非平穩(wěn)、 非高斯隨機(jī)系列，能抑制高斯噪聲、分析周期重復(fù)和調(diào)制特性信號優(yōu)勢[18]提取信號特征。根據(jù)滾動(dòng)軸承振動(dòng)信號雙譜在故障類型不同時(shí)差異性顯著，在軸承狀態(tài)相同、但工況不同時(shí)又具有一定相似性的特點(diǎn)，從雙譜中提取工況魯棒的特征參數(shù)。

基于雙譜對稱特性降低數(shù)據(jù)維數(shù)。雙譜幅值分布滿足y=x和y=-x對稱，過中心原點(diǎn)將雙譜劃分成八個(gè)子區(qū)域后，使用兩個(gè)子區(qū)域數(shù)據(jù)就能表征雙譜信息。將選擇的子區(qū)域進(jìn)行旋轉(zhuǎn)及數(shù)據(jù)合并，得到數(shù)據(jù)量僅為原始雙譜四分之一的優(yōu)選雙譜，且降維過程運(yùn)算量小于主成分分析[19]和灰度共生矩陣[20]等方法。

使用奇異值分解提取特征參數(shù)。將優(yōu)選雙譜做奇異值分解，依奇異值累積貢獻(xiàn)率選擇奇異向量，計(jì)算元素模值后獲得故障診斷特征向量。與雙譜各種切片譜[21-23]、雙譜倒譜[24]、雙譜能量等參數(shù)相比，奇異值分解提取特征參數(shù)具有性能穩(wěn)健和抑制干擾等優(yōu)點(diǎn)。最后用連續(xù)型隱馬爾可夫模型辨識故障，憑借非平穩(wěn)、高冗余及重復(fù)性低時(shí)間序列建模優(yōu)勢，利用觀察符號可見和狀態(tài)隱藏特點(diǎn)，從工況擾動(dòng)的特征參量中判別軸承狀態(tài)。

2 滾動(dòng)軸承振動(dòng)信號雙譜

基于滾動(dòng)軸承振動(dòng)信號模型，推導(dǎo)雙譜理論表征，分析軸承狀態(tài)、載荷和轉(zhuǎn)速等要素如何影響振動(dòng)信號及其雙譜。

2.1 滾動(dòng)軸承振動(dòng)信號模型

依據(jù)滾動(dòng)軸承振動(dòng)模型，當(dāng)滾動(dòng)體周期性通過故障部位時(shí)，可抽象成一周期性脈沖作用于軸承系統(tǒng)，產(chǎn)生的振動(dòng)加速度信號可以表征為[26]：

(1)

式(1)中，T0和τi分別為沖擊激勵(lì)周期及其擾動(dòng)，s(t)為周期性激勵(lì)產(chǎn)生的衰減振蕩，Ai為幅度調(diào)制成分，n(t)為加性噪聲。s(t)和Ai分別表示為：

(2)

(3)

式(2)中fn和β分別為軸承系統(tǒng)的固有頻率和衰減系數(shù)。式(3)中AM和AT取值由徑向固有載荷和變化載荷決定，fr為轉(zhuǎn)頻。

若軸承狀態(tài)相同、但運(yùn)行工況發(fā)生變化，載荷和轉(zhuǎn)速對振動(dòng)信號產(chǎn)生下述影響：

(1) 隨著軸向載荷增加，系統(tǒng)共振頻率fn會(huì)逐漸變大，但上升幅度會(huì)逐漸減小[27]。

(2) 轉(zhuǎn)速除了影響衰減振蕩成分s(t)的激勵(lì)周期T0，對應(yīng)的轉(zhuǎn)頻fr還調(diào)制幅度分量Ai。

因此，滾動(dòng)軸承振動(dòng)信號的幅度特性受轉(zhuǎn)速調(diào)制，衰減振蕩成分則受載荷和轉(zhuǎn)速雙重因素影響。

2.2 振動(dòng)信號雙譜表征

將信號x(t)作零均值預(yù)處理，其三階累積量可以表示為：

(4)

將式(1)代入式(4)，有：

(5)

計(jì)算(5)式數(shù)學(xué)期望：

C3(τ1,τ2)=

(6)

(7)

雙譜定義為三階累積量的二維傅立葉變換：

τ2)e-j2πf1τ1e-j2πf2τ2dτ1dτ2dt.

(8)

將式(6)代入式(8)，求得振動(dòng)信號x(t)的雙譜：

(9)

為驗(yàn)證上述理論推導(dǎo)，使用凱斯西儲(chǔ)大學(xué)軸承研究中心數(shù)據(jù)，計(jì)算了不同故障類型、不同工況相同故障類型的振動(dòng)信號雙譜。其中0.355 6mm故障點(diǎn)直徑的滾動(dòng)體、內(nèi)圈和外圈故障雙譜如圖1所示。0.533 4mm直徑外圈故障，四種不同工況下的振動(dòng)信號雙譜如圖2所示。

圖1 不同故障狀態(tài)雙譜

圖2 不同工況外圈故障雙譜

圖1和圖2結(jié)果表明，不同故障狀態(tài)雙譜分布特性和幅值大小皆不相同，而不同工況、同種故障類型雙譜的幅值大小和分布特性較為相似。亦即不同軸承狀態(tài)雙譜差異性顯著，不同工況、相同故障狀態(tài)雙譜相似性顯著，驗(yàn)證了前述理論推導(dǎo)結(jié)論。

3 雙譜特征參數(shù)提取

3.1 基于對稱特性降低雙譜維數(shù)

振動(dòng)信號x(t)為實(shí)序列，其傅立葉變換滿足共軛對稱性S(f)=S*(-f)，式(9)具有下述特征：

B(f1,f2)=B(f2,f1)=B*(-f1,-f2)=

B*(-f2,-f1).

(10)

因而雙譜分布滿足下述對稱特性：

(11)

過雙譜中心原點(diǎn)，依直線x=0、y=0、y=x、和y=-x，將雙譜劃分為如圖3所示的8個(gè)子區(qū)域。

圖3 雙譜分布特性

式(11)和圖3表明雙譜幅度分布具有較大信息冗余，利用對稱特性降低數(shù)據(jù)維數(shù)，步驟如下：

(3) 將區(qū)域(III)數(shù)據(jù)以原點(diǎn)為軸心，向右旋轉(zhuǎn)90度后到達(dá)區(qū)域(I)位置，使降維后的數(shù)據(jù)具有正方性質(zhì)。數(shù)據(jù)旋轉(zhuǎn)通過下式完成：

(12)

式(12)中x和y為區(qū)域(III)數(shù)據(jù)的原始坐標(biāo)位置，x'和y'為旋轉(zhuǎn)后的坐標(biāo)位置。

圖4 優(yōu)選雙譜

3.2 基于奇異值分解提取特征參數(shù)

將優(yōu)選雙譜Φ作奇異值分解：

Φ=UΣVT=

(13)

式(13)中，矩陣U和V的列分別稱為左、右奇異向量，σi(i=1,…,M)為與Ui、Vi相對應(yīng)的奇異值。奇異值大小能表征奇異向量重要程度，其分布特性如圖5所示：

圖5 奇異值分布特性

圖5表明優(yōu)選雙譜奇異值快速減小，因而只要選擇K個(gè)大奇異值對應(yīng)的奇異向量，就能保留優(yōu)選雙譜的大部分信息并有效降低數(shù)據(jù)維數(shù)。K個(gè)奇異向量具有的信息量用奇異值累積貢獻(xiàn)率表征[29]：

(14)

基于奇異值分解提取特征參數(shù)過程如下：

(1) 根據(jù)奇異值累積貢獻(xiàn)率ξ，計(jì)算所需的奇異向量數(shù)目K；

(3) 計(jì)算矩陣ψ各元素ψij的模值，作為特征參數(shù)構(gòu)造M×K維特征向量O：

(15)

依上述過程得到的觀察向量O，歷經(jīng)雙譜優(yōu)選和奇異值分解，繼承大部分雙譜信息前提下大幅度降低數(shù)據(jù)維數(shù)，實(shí)現(xiàn)振動(dòng)信號雙譜特征提取。

4 基于隱馬爾可夫模型智能故障診斷

4.1 連續(xù)型隱馬爾可夫模型概述[30]

將特征向量O的組成元素Oi,j作為觀察符號，依聚類原理分成N類，每一個(gè)類別屬于隱馬爾可夫模型的一個(gè)狀態(tài)，得到模型狀態(tài)集S：

S={S1,S2,…,SN}.

(16)

用狀態(tài)初始概率分布矩陣描述模型初始狀態(tài)：

π=[π1,π2,…,πN].

(17)

元素πi表示狀態(tài)的初始分布概率：

(18)

式中q0表示模型初始時(shí)刻的狀態(tài)。

時(shí)間序列建模過程中，模型狀態(tài)隨觀察向量的變化而改變，變化過程用狀態(tài)轉(zhuǎn)移概率矩陣描述：

(19)

aij表示模型從狀態(tài)i轉(zhuǎn)移到狀態(tài)j的概率：

(20)

隱馬爾可夫模型的轉(zhuǎn)移概率aij無法直接觀察，但能用觀察符號概率間接表征，通過高斯混合獲得：

(21)

式中bj(Ot)表示t時(shí)刻觀察符號Ot屬于狀態(tài)j的概率。η為一高斯概率密度函數(shù)，參數(shù)μ和U為混合均值矩陣和協(xié)方差矩陣，M是混合數(shù)目；cjm為混合權(quán)系數(shù)，N個(gè)狀態(tài)的權(quán)系數(shù)構(gòu)成混合權(quán)矩陣C。

連續(xù)型隱馬爾可夫模型表示為：

λ=(A,π,C,μ,U).

(22)

實(shí)際使用時(shí)，基于可見的觀察向量，使用最大期望算法確定模型參數(shù)值。

4.2 故障診斷過程描述

變工況智能故障診斷為一基于特征參數(shù)的模式識別過程，由模型訓(xùn)練和故障辨識兩個(gè)環(huán)節(jié)構(gòu)成。主要步驟如下：

(1) 提取信號特征。模型訓(xùn)練和故障辨識兩個(gè)環(huán)節(jié)，皆需估計(jì)振動(dòng)信號雙譜和優(yōu)選雙譜，并依奇異值分解獲得特征向量。

(2) 訓(xùn)練模型集。在模型訓(xùn)練階段，先將隱馬爾可夫模型參數(shù)在取值范圍內(nèi)隨機(jī)賦值。再依當(dāng)前特征向量用Baum-Welch算法[30]更新模型參數(shù)，獲得當(dāng)前軸承狀態(tài)隱馬爾可夫模型的初始分布概率、狀態(tài)轉(zhuǎn)移概率、混合權(quán)系數(shù)、均值矩陣和協(xié)方差矩陣。每一種工況，都計(jì)算正常狀態(tài)、滾動(dòng)體故障、內(nèi)圈和外圈故障的隱馬爾可夫模型，將其標(biāo)識為λ1、λ2、λ3、λ4，形成診斷模型集。

(3) 軸承狀態(tài)判別。故障辨識階段，將任意工況、未知軸承狀態(tài)獲取的特征向量輸入模型集，并計(jì)算不同軸承狀態(tài)模型下的似然概率P(O|λi)，取值最大的模型狀態(tài)被判別為待診斷的軸承狀態(tài)。

為保證故障診斷有效性，模型訓(xùn)練和故障辨識兩個(gè)環(huán)節(jié)數(shù)據(jù)沒有重疊。具體流程如圖6所示。

圖6 基于隱馬爾可夫模型的變工況智能故障診斷

5 實(shí)驗(yàn)研究

5.1 數(shù)據(jù)描述

研究對象為安裝于電機(jī)驅(qū)動(dòng)端的6205-2RS軸承，通過調(diào)節(jié)測功機(jī)扭矩和電機(jī)轉(zhuǎn)速模擬運(yùn)行工況，詳情如表1所示。

表1 工況參數(shù)

工況載荷/kW轉(zhuǎn)速/(r/min)工況101797工況20.7351772工況31.4711750工況42.2061730

用電火花加工點(diǎn)蝕故障模擬不同的故障類型。每一類故障，都加工深度為0.279 4mm，直徑分別為0.177 8mm，0.355 6mm和0.533 4mm的故障點(diǎn)，來表征不同的故障程度。數(shù)據(jù)由凱斯西儲(chǔ)大學(xué)軸承研究中心采集，加速度傳感器用磁力吸座安裝于軸承座正上方，信號采樣頻率為12kHz。每一類故障程度，皆有4種工況、4種軸承狀態(tài)振動(dòng)加速度信號用于變工況故障診斷研究。

5.2 參數(shù)設(shè)置

(1) 雙譜計(jì)算。用直接法估計(jì)，序列長度取1 024數(shù)據(jù)點(diǎn)，數(shù)據(jù)片段取256點(diǎn)、75%重復(fù)度，頻域用寬度為5的Rao-Gabr窗平滑。

(2) 奇異向量數(shù)目選擇。奇異值累積貢獻(xiàn)率ξ在0.177 8mm故障時(shí)設(shè)為50%，0.355 6mm故障時(shí)為80%，0.533 4故障時(shí)為85%，對應(yīng)的奇異向量個(gè)數(shù)分別取4，16和12。

(3) 模型訓(xùn)練。隱馬爾可夫模型的狀態(tài)和混合個(gè)數(shù)設(shè)為2，參數(shù)估計(jì)的收斂門限值及最大迭代次數(shù)分別設(shè)為0.001和50。

5.3 實(shí)驗(yàn)結(jié)果

實(shí)驗(yàn)依故障程度不同分3組進(jìn)行。每組實(shí)驗(yàn)中，先訓(xùn)練不同工況的隱馬爾可夫模型集，再用其辨識任意工況軸承狀態(tài)。模型訓(xùn)練和故障辨識使用不同數(shù)據(jù)以保證診斷有效性，每種狀態(tài)用1 000組數(shù)據(jù)驗(yàn)證診斷效果，正確識別的數(shù)據(jù)組數(shù)與輸入數(shù)據(jù)組數(shù)間的比值稱為辨識精度。三種故障程度診斷結(jié)果見表2、表3和表4，表中最左列表示待辨識數(shù)據(jù)，后面4列表示不同工況診斷模型的辨識結(jié)果，詳情如表2—表4所示。

表2 0.177 8 mm故障程度辨識精度 /%

表3 0.355 6mm故障程度辨識精度 /%

表4 0.533 4mm故障程度辨識精度 /%

5.4 實(shí)驗(yàn)分析

實(shí)驗(yàn)結(jié)果表明，在三種不同的故障程度下，使用任意一種工況數(shù)據(jù)訓(xùn)練診斷模型，都能正確辨識四種工況的軸承狀態(tài)。即診斷方法克服了轉(zhuǎn)速和載荷變化對特征參量的擾動(dòng)，是一種適應(yīng)工況變化的魯棒診斷方法。與當(dāng)前常用的階次方法相比，文中方法無需數(shù)據(jù)重采樣，也不要提取設(shè)備轉(zhuǎn)速信息，是一種穩(wěn)健、高效的智能故障診斷方法。

表2、表3和表4的診斷結(jié)果有下述特點(diǎn)：

(1) 不同故障程度的辨識精度存在差異。如表2中辨識精度皆大于99%，而表3和表4的辨識精度在某些軸承狀態(tài)下稍低一些。其原因是不同故障程度時(shí)產(chǎn)生的周期性激勵(lì)不一致，而某些情形的激勵(lì)使雙譜產(chǎn)生較大擾動(dòng)影響模型辨識精度。

(2) 訓(xùn)練數(shù)據(jù)與辨識數(shù)據(jù)間的工況差異程度影響辨識精度。如表3中工況1模型對工況4的內(nèi)圈故障診斷，表3中工況4模型對工況1和工況2的滾動(dòng)體故障診斷，以及表4中工況4模型對工況1的滾動(dòng)體故障診斷，辨識精度要比工況接近時(shí)低一些。其原因是雙譜特征受工況差異程度影響，工況接近時(shí)雙譜較為相似，反之則差異較大。圖7(a)和圖7(b)為表3中工況1和工況2的優(yōu)選雙譜，其幅度和分布特性比較接近。圖7(c)為表3中工況4優(yōu)選雙譜，其幅度和分布特性與前兩種工況差別較大。雙譜差異導(dǎo)致特征參數(shù)擾動(dòng)，干擾隱馬爾可夫模型似然概率估值并影響軸承狀態(tài)辨識。

圖7 三種工況內(nèi)圈故障優(yōu)選雙譜，0.355 6 mm故障點(diǎn)直徑

5.5 診斷方法可適性驗(yàn)證

本節(jié)對另一型號、不同使用環(huán)境滾動(dòng)軸承開展變工況故障診斷，以驗(yàn)證診斷方法的普適性。軸承型號為6203-2RS，安裝于電機(jī)風(fēng)扇端，用電火花加工點(diǎn)蝕故障模擬不同故障類型，用不同直徑故障點(diǎn)模擬不同故障程度。安裝時(shí)將外圈故障點(diǎn)置于6點(diǎn)鐘位置，電機(jī)的轉(zhuǎn)速和載荷數(shù)據(jù)與表1相同，故障點(diǎn)直徑與4.1節(jié)的6205-2RS軸承相同。數(shù)據(jù)由凱斯西儲(chǔ)大學(xué)軸承研究中心提供，0.177 8mm和0.355 6mm故障點(diǎn)直徑時(shí)包含四種運(yùn)行工況數(shù)據(jù)，0.533 4mm故障點(diǎn)直徑時(shí)包含工況2、工況3和工況4的數(shù)據(jù)。

三種故障程度的奇異向量數(shù)目分別為12、24和12，雙譜估計(jì)、隱馬爾可夫模型訓(xùn)練參數(shù)與4.2節(jié)一致，辨識精度用1 000組數(shù)據(jù)計(jì)算。6203-2RS軸承的辨識精度表明診斷方法具有較好的適用性，詳見表5、表6和表7。

表5 6203-2RS軸承0.177 8mm故障程度辨識精度/%

故障類型使用的診斷模型工況1工況2工況3工況4工況1,正常100100100100工況1,滾動(dòng)體98.996.298.199.8工況1,內(nèi)圈10010010085工況1,外圈99.399.990.897.1工況2,正常100100100100工況2,滾動(dòng)體99.899.897.799.6工況2,內(nèi)圈98.110010099.8工況2,外圈92.299.810098.2工況3,正常100100100100工況3,滾動(dòng)體10010099.9100工況3,內(nèi)圈80.510010099工況3,外圈79.710097.799.7工況4,正常100100100100工況4,滾動(dòng)體98.798.798.798.9工況4,內(nèi)圈76.8100100100工況4,外圈91.510097.699.8

表6 6203-2RS軸承0.355 6mm故障程度辨識精度/%

故障類型使用的診斷模型工況1工況2工況3工況4工況1,正常100100100100工況1,滾動(dòng)體96.684.392.698.8工況1,內(nèi)圈100100100100工況1,外圈10010098.2100工況2,正常100100100100工況2,滾動(dòng)體96.883.49899.9工況2,內(nèi)圈100100100100工況2,外圈100100100100工況3,正常100100100100工況3,滾動(dòng)體75.578.3100100工況3,內(nèi)圈100100100100工況3,外圈100100100100工況4,正常100100100100工況4,滾動(dòng)體8375.6100100工況4,內(nèi)圈98.9100100100工況4,外圈100100100100

表7 6203-2RS軸承0.533 4mm故障程度辨識精度/%

故障類型使用的診斷模型工況2工況3工況4工況2,正常100100100工況2,滾動(dòng)體98.499.796.8工況2,內(nèi)圈99.799.9100工況2,外圈100100100工況3,正常100100100工況3,滾動(dòng)體98.599.795.9工況3,內(nèi)圈100100100工況3,外圈100100100工況4,正常100100100工況4,滾動(dòng)體8397.397.7工況4,內(nèi)圈99100100工況4,外圈99.399100

6 結(jié)論

(1) 雙譜能有效表征變工況情形的滾動(dòng)軸承狀態(tài)，故障類型不同時(shí)具有顯著差異性，故障相同但工況不同時(shí)又具有較高相似性；滾動(dòng)軸承振動(dòng)信號雙譜具有較大信息冗余，基于對稱特性剔除冗余數(shù)據(jù)能大幅度降低數(shù)據(jù)維數(shù)。

(2) 奇異值分解是一種有效的特征提取方法，經(jīng)奇異值累積貢獻(xiàn)率選擇奇異向量，能獲得可靠的特征信息并降低數(shù)據(jù)維數(shù)。

(3) 隱馬爾可夫模型是一種穩(wěn)健的時(shí)間序列建模工具，能適應(yīng)軸承運(yùn)行工況變化引起的觀察向量擾動(dòng)，并依據(jù)可見的特征參數(shù)辨識隱藏的故障狀態(tài)；使用期望最大算法計(jì)算隱馬爾可夫模型參數(shù)，訓(xùn)練過程存在收斂于局部最優(yōu)值情形，最優(yōu)的隱馬爾可夫模型初始參數(shù)選擇還有待于進(jìn)一步研究。

[1]El-ThaljiI,JantunenE.Asummaryoffaultmodellingandpredictivehealthmonitoringofrollingelementbearings[J].MechanicalSystemsandSignalProcessing, 2015, 60-61: 252-272.

[2] 李勤良，汪博，趙斌，等. 考慮非線性遲滯力的機(jī)床顫振系統(tǒng)穩(wěn)定性研究[J]. 機(jī)械工程學(xué)報(bào), 2013(11): 43-49.

[3]VillaLF,ReonesA,PeránJR,etal.Angularresamplingforvibrationanalysisinwindturbinesundernon-linearspeedfluctuation[J].MechanicalSystemsandSignalProcessing, 2011, 25(6): 2157-2168.

[4] 林京，趙明. 變轉(zhuǎn)速下機(jī)械設(shè)備動(dòng)態(tài)信號分析方法的回顧與展望[J]. 中國科學(xué):技術(shù)科學(xué), 2015(7): 669-686.

[5]WangY,XuG,ZhangQ,etal.Rotatingspeedisolationanditsapplicationtorollingelementbearingfaultdiagnosisunderlargespeedvariationconditions[J].JournalofSoundandVibration, 2015, 348: 381-396.

[6] 程衛(wèi)東，王天楊，WangJin-jiang，等. 包絡(luò)變形對階比分析結(jié)果的影響及消除方法[J]. 振動(dòng)工程學(xué)報(bào), 2015(3): 470-477.

[7]WangT,LiangM,LiJ,etal.Bearingfaultdiagnosisunderunknownvariablespeedviagearnoisecancellationandrotationalordersidebandidentification[J].MechanicalSystemsandSignalProcessing, 2015, 62-63: 30-53.

[8]MishraC,SamantarayAK,ChakrabortyG.Rollingelementbearingdefectdiagnosisundervariablespeedoperationthroughanglesynchronousaveragingofwaveletde-noisedestimate[J].MechanicalSystemsandSignalProcessing, 2016, 72-73: 206-222.

[9] 王天楊，李建勇，程衛(wèi)東. 基于低次故障特征階比系數(shù)的變轉(zhuǎn)速滾動(dòng)軸承等效轉(zhuǎn)頻估計(jì)算法[J]. 機(jī)械工程學(xué)報(bào), 2015(3): 121-128.

[10] 徐亞軍，于德介，孫云嵩，等. 滾動(dòng)軸承故障診斷的階比多尺度形態(tài)學(xué)解調(diào)方法[J]. 振動(dòng)工程學(xué)報(bào), 2013(2): 252-259.

[11]MingZ,JingL,XiaoqiangX.Multi-FaultDetectionofRollingElementBearingsunderHarshWorkingConditionUsingIMF-BasedAdaptiveEnvelopeOrderAnalysis[J].Sensors, 2014, 14(11): 20320-20346.

[12]ShiJ,LiangM,GuanY.Bearingfaultdiagnosisundervariablerotationalspeedviathejointapplicationofwindowedfractaldimensiontransformandgeneralizeddemodulation:Amethodfreefromprefilteringandresampling[J].MechanicalSystemsandSignalProcessing, 2016, 68-69: 15-33.

[13] 范彬，胡雷，胡蔦慶. 變工況下旋轉(zhuǎn)機(jī)械故障跟蹤的相空間曲變方法[J]. 物理學(xué)報(bào), 2013(16): 63-70.

[14]WanX,WangD,TsePW,etal.Acriticalstudyofdifferentdimensionalityreductionmethodsforgearcrackdegradationassessmentunderdifferentoperatingconditions[J].Measurement, 2016, 78: 138-150.

[15] 周勇，何創(chuàng)新. 基于獨(dú)立特征選擇與相關(guān)向量機(jī)的變載荷軸承故障診斷[J]. 振動(dòng)與沖擊,2012(3): 157-161.

[16] 王珊，蘇璐瑋，顧煜炯，等. 變工況特性下的風(fēng)電軸承早期故障診斷方法[J]. 電力科學(xué)與工程,2014(03): 39-44.

[17] 宋京偉，梅秀莊，郭云. 變轉(zhuǎn)速變載荷下滾動(dòng)軸承故障的小波-HMM診斷[J]. 軸承，2005(12): 28-32.

[18]DongG,ChenJ,ZhaoF.Afrequency-shiftedbispectrumforrollingelementbearingdiagnosis[J].JournalofSoundandVibration. 2015, 339: 396-418.

[19] 張銳戈，譚永紅. 雙譜主成分分析的滾動(dòng)軸承智能故障診斷[J]. 振動(dòng)工程學(xué)報(bào)，2014(05): 763-769.

[20] 林勇，胡夏夏，朱根興，等. 基于振動(dòng)譜圖像識別的智能故障診斷[J]. 振動(dòng)、測試與診斷，2010, 30(2): 175-180.

[21]ZhouY,ChenJ,DongGM,etal.Applicationofthehorizontalsliceofcyclicbispectruminrollingelementbearingsdiagnosis[J].MechanicalSystemsandSignalProcessing. 2012, 26(0): 229-243.

[22] 張琳，黃敏. 基于EMD與切片雙譜的軸承故障診斷方法[J]. 北京航空航天大學(xué)學(xué)報(bào)，2010(3): 287-290.

[23] 周宇，陳進(jìn)，董廣明，等. 基于循環(huán)雙譜的滾動(dòng)軸承故障診斷[J]. 振動(dòng)與沖擊. 2012, 31(09): 78-81.

[24] 李輝，鄭海起，唐力偉. 基于倒雙譜分析的軸承故障診斷研究[J]. 振動(dòng)、測試與診斷，2010, 30(4): 353-356.

[25] 黃晉英，潘宏俠，畢世華，等. 基于高階累量譜的軸承故障診斷[J]. 火炮發(fā)射與控制學(xué)報(bào)，2007(2): 56-59.

[26] Cong F, Chen J, Dong G, et al. Vibration model of rolling element bearings in a rotor-bearing system for fault diagnosis[J]. Journal of Sound and Vibration. 2013, 332(8): 2081-2097.

[27] 趙聯(lián)春，馬家駒，馬純青. 載荷對球軸承振動(dòng)特性的影響[J]. 軸承，2003(05): 38-41.

[28] Yang D M, Stronach A F, Macconnell P, et al. Third-order spectral techniques for the diagnosis of motor bearing condition using artificial neural netwroks[J]. Mechanical Systems and Signal Processing. 2002, 16(2-3): 391-411.

[29] He Q, Yan R, Kong F, et al. Machine condition monitoring using principal component representations[J]. Mechanical Systems and Signal Processing. 2009, 23(2): 446-466.

[30] Rabiner L R. A tutorial on hidden Markov models and selected applications in speech recognition[J]. Readings in speech recognition. 1990, 53(3): 267-296.

Fault diagnosis of variety operation conditions for rolling element bearings based on bispectrum feature extraction

ZHANGRui-ge1,2,LIJing3

(1. College of Mechanical and Electronic Engineering, Sanming University, Sanming 365004,China; 2. Key Laboratory of Intelligent Equipment of Fujian Province Universities, Sanming 365004,China; 3. School of Information and Electrical Engineering, Hunan University of Science and Technology, Xiangtan 411201,China)

This paper presents a novel method for the feature extraction and bearing status identification under the varying loads and rotation speed, which can be applied in the fault diagnosis of the rotating equipment of variety operation conditions. Firstly, we explored the bispectrum of the vibration acceleration signals, and found the bispectrum distribution was of significant differences under different bearing status, but very similar to the same fault type under different operation conditions. Thus, it is appropriate to extract the robust features from the bispectrum for the fault diagnosis of variety operation conditions. Then, the bispectrum was divided into eight sub-zones according to the symmetric characteristics. After the redundant data were removed, the retained data were called the preferred bispectrum, which only had a quarter data of the original bispectrum. The feature extraction applied the singular value decomposition, the number of the singular vectors was also determined by the accumulated contribution rate by the singular value, and the modulus of the selected singular vectors was as the features. Finally, we conducted a number of experiments under four operation conditions of three fault severities to validate the efficiency of the proposed method. The results show that the proposed method effectively identify various fault types and is robust to the varying operating conditions.

fault diagnosis;variety operation conditions;bispectrum;singular value decomposition;Hidden Markov model

2017-05-02.

張銳戈(1977-)，男，博士，副教授，E-mail:ruig_zhang@126.com.

國家自然科學(xué)基金資助項(xiàng)目(61474042)；福建省科技廳重點(diǎn)項(xiàng)目(2013N0032).

2095-7386(2017)02-0057-10

10.3969/j.issn.2095-7386.2017.02.012

TN 911.6;TP 206.3

A