施 瑩， 莊 哲， 林建輝

(西南交通大學(xué)牽引動(dòng)力國家重點(diǎn)實(shí)驗(yàn)室 成都, 610031)

引 言

輪對軸承是高速列車走行部的重要部件，它既承擔(dān)了由列車自重及負(fù)載合成的垂向力，還承擔(dān)了輪軌間橫向力。在高速列車長期運(yùn)行過程中，這些外力的聯(lián)合作用使輪對軸承有可能發(fā)生如剝落、點(diǎn)蝕和燒損等故障，這些故障一旦出現(xiàn)，極易在短時(shí)間內(nèi)惡化，最終危及高速列車的運(yùn)行安全。因此，研究高速列車輪對軸承的故障診斷方法對保障列車安全運(yùn)行具有重要意義。

輪對軸承故障檢測的方法主要包括油脂(液)監(jiān)測、聲發(fā)射監(jiān)測、軸溫監(jiān)測及振動(dòng)監(jiān)測[1]。然而，油脂(液)的監(jiān)測方法不適合軸承的實(shí)時(shí)故障監(jiān)測；聲發(fā)射技術(shù)因軸承故障引起的瞬態(tài)彈性波衰減現(xiàn)象使其無法得到廣泛使用；軸溫監(jiān)測是目前較常用的檢測手段：當(dāng)軸箱溫度超過事先設(shè)定的閾值時(shí)，軸溫監(jiān)測裝置發(fā)出警報(bào)，但該方法只對嚴(yán)重故障敏感，無法診斷早期故障[2-3]。軸箱振動(dòng)信號不但易于監(jiān)測，還包含豐富的軸承狀態(tài)信息。因此，振動(dòng)監(jiān)測將成為輪對軸承故障，特別是早期故障的有效監(jiān)測手段。

軸承故障引起的脈沖常與車輪踏面缺陷或輪對不圓時(shí)引起的脈沖混雜在一起，從這些非線性和非平穩(wěn)的具備調(diào)制特性、強(qiáng)噪聲干擾的信號中提取有用故障特征一直是一個(gè)挑戰(zhàn)[4]。

為了有效地應(yīng)對挑戰(zhàn)，眾多學(xué)者提出了許多先進(jìn)的信號處理技術(shù)，包括經(jīng)驗(yàn)?zāi)B(tài)分解[5](empirical mode decomposition,簡稱EMD)，小波變換(wavelet transform,簡稱WT)，并成功地應(yīng)用于旋轉(zhuǎn)機(jī)械的故障診斷中[6-7]。EMD非常適合于分析非線性非平穩(wěn)信號，并能將信號分解為一組本征模函數(shù)(intrinsic mode function,簡稱IMF)和殘余分量，已成功地應(yīng)用于故障檢測、診斷等領(lǐng)域的研究。EMD的缺點(diǎn)(如缺乏理論基礎(chǔ)，噪聲的敏感性，模態(tài)混疊等問題[8])限制了其進(jìn)一步廣泛應(yīng)用。為了解決這些問題，聚合經(jīng)驗(yàn)?zāi)B(tài)分解[9](ensemble empirical mode decomposition,簡稱EEMD)等算法被陸續(xù)提出，但以上問題并沒有從根本上得到解決。WT是另外一種常用的時(shí)頻分析方法。相較于離散小波變換(discrete wavelet transformation,簡稱DWT)的快速算法，連續(xù)小波變換(continuous wavelet transform,簡稱CWT)巨大的計(jì)算成本和冗余系數(shù)給實(shí)際應(yīng)用帶來了巨大困難，小波包變換(wavelet packet transform,簡稱WPT)的提出使信號的高頻部分分解性能得到了較大提升，然而DWT和WPT的分解效果受小波基函數(shù)影響較大。

基于交替方向乘子法[10](alternating direction method of multipliers,簡稱ADMM)的卷積稀疏表示于2015年提出，該方法為軸承故障診斷的特征提取提供了良好前景。然而，CSR提取脈沖的性能受到其懲罰因子的影響較大：當(dāng)懲罰因子過大時(shí)，提取的脈沖小于真正的脈沖數(shù)；當(dāng)懲罰因子過小時(shí)，提取的脈沖大于真正的脈沖數(shù)，這便會(huì)給故障的精確診斷帶來極大困難。

流形學(xué)習(xí)是一種非線性降維方法，目前廣泛應(yīng)用的流形學(xué)習(xí)算法一般分為兩種:基于全局的流形學(xué)習(xí)分析方法以及基于局部的流形學(xué)習(xí)分析方法?；谌值姆治龇椒ㄖ饕械染嘤成鋄11]和最大方差展開算法[12](maximum variance unfolding,簡稱MVU)等?；诰植康姆治龇椒ㄖ饕芯植壳锌臻g分析[13](local tangent space alignment,簡稱LTSA)、局部線性嵌入[14](locally linear embedding,簡稱LLE)和拉普拉斯特征映射算法[15](Laplacian eigenmaps,簡稱LE)等。

筆者提出了一種基于卷積稀疏表示、希爾伯特變換和流形學(xué)習(xí)的新的故障檢測方法。首先，通過在不同懲罰因子下的卷積稀疏表示提取不同稀疏特征的脈沖；然后，針對提取的一系列脈沖進(jìn)行希爾伯特變換，構(gòu)造脈沖包絡(luò)空間；最后，利用等距映射流形學(xué)習(xí)算法對脈沖包絡(luò)空間求解低維本征包絡(luò)，以實(shí)現(xiàn)故障診斷。

1 基本原理

1.1 卷積稀疏表示脈沖提取

傳統(tǒng)的信號稀疏表示理論基于正交線性變換，通過將信號展開為一系列同形態(tài)基函數(shù)線性加權(quán)和的形式。然而，在工程實(shí)際中大部分信號都不是單一形態(tài)的成分，因而很難用單一的形態(tài)基函數(shù)稀疏表示。因此，就需要利用兩種或兩種以上的不同形態(tài)基函數(shù)對信號加以表示，其基本思想是用超完備的冗余基函數(shù)原子庫取代正交的基函數(shù)庫，原子庫可以由任意基函數(shù)構(gòu)成，從而最優(yōu)地匹配原始信號的結(jié)構(gòu)。從原子庫中尋求信號的最佳線性原子組合來表示信號，被稱為信號的稀疏表示。

卷積稀疏表示理論將信號表示為一系列原子與稀疏系數(shù)的卷積和，其定義[10]為

(1)

其中：sk∈Rn為第k項(xiàng)待分析的振動(dòng)信號；n為信號長度；dm∈Rp為一系列原子；p為原子長度；xk,m∈Rn-p+1為與給定信號sk和原子dm相關(guān)的稀疏系數(shù)；λ∈R+為懲罰因子；K為給定信號sk的個(gè)數(shù)；M為原子dm的個(gè)數(shù)。

當(dāng)給定一組信號和原子，式(1)中的稀疏系數(shù)可由平移不變稀疏編碼[10](shift invariant sparse coding,簡稱SISC)得到，該方法將信號表示為多個(gè)原子與系數(shù)的卷積形式，自適應(yīng)地學(xué)習(xí)原子與系數(shù)達(dá)到對原信號的稀疏表達(dá)，其表達(dá)式為

(2)

式(2)中，不同稀疏系數(shù)xk,m針于不同信號sk進(jìn)行解耦，稀疏系數(shù)xm表示為

(3)

當(dāng)信號sk及稀疏系數(shù)xk,m已知，式(1)中的原子dm可由平移不變性字典學(xué)習(xí)[10](shift invariant dictionary learning,簡稱SIDL)得到，該算法為變量交替更新算法，即輪流更新x和d，表達(dá)式為

(4)

所以，dm=1?m。

(5)

1.2 希爾伯特脈沖包絡(luò)提取

關(guān)于軸承故障的特征信息包含在提取的脈沖的包絡(luò)中，提取的脈沖的希爾伯特變換定義[16]為

(6)

(7)

(8)

(9)

1.3 脈沖包絡(luò)空間流形降維

流形學(xué)習(xí)將式(9)所示的高維脈沖包絡(luò)空間IESm映射至一個(gè)低維本征脈沖包絡(luò)中。文中采用等距映射流形算法進(jìn)行降維計(jì)算。

基于等距映射流形方法是一種非迭代、非線性的數(shù)據(jù)降維算法，其基本原理是：將任意兩個(gè)數(shù)據(jù)點(diǎn)之間的測地線距離(也就是兩點(diǎn)間的最短距離)當(dāng)成流形的幾何表述，通過多維縮放算法(multiple dimensional scaling,簡稱MDS)算法保持這個(gè)點(diǎn)與點(diǎn)之間的最短距離，形成局部測地距離和低維非線性特征空間的映射關(guān)系，進(jìn)而獲取最佳的低維流形結(jié)果，還可以滿足數(shù)據(jù)測地距離在低維流形空間維持不變的特點(diǎn)。

等距映射流形算法的優(yōu)化函數(shù)[11]如下

(10)

τ(DG)-τ(Dγ)L2

(11)

那么Isomap的最優(yōu)化解Y=(y1,y2,…,ym)可以通過τ(DG)的最大特征值表示的向量得到。

1.4 估測原子類型數(shù)量

在式(1)中，原子的數(shù)目是不確定的，而是基于軸承故障振動(dòng)特性來確定[17]。脈沖可以用共振頻率和阻尼參數(shù)兩個(gè)參數(shù)來完全描述，共振頻率可以定義不同的原子，因此利用原子的主頻率來估計(jì)原子數(shù)，考慮到軸承的故障特性，在輪對軸承故障檢測中，一般預(yù)先設(shè)定為4，然后通過CSR的迭代得到M的準(zhǔn)確值。每個(gè)原子的主頻率通過傅里葉變換提取fm(m=1,2,…,M)。如果任何兩個(gè)主頻率的差大于頻率分辨率fsp-1(fs表示采樣頻率)，則表明有原子；否則M=M-1，循環(huán)執(zhí)行CSR并且計(jì)算主頻率直到任何兩個(gè)主頻率的差大于頻率分辨率。輸出M的終值作為原子數(shù)目估測值。

1.5 基于卷積稀疏表示及流形學(xué)習(xí)故障診斷總體框架

基于卷積稀疏表示及流形學(xué)習(xí)的軸承故障診斷流程如圖1所示。

圖1 基于卷積稀疏表示及流形學(xué)習(xí)的故障診斷流程Fig.1 Flowchart of the fault detection based on CSR and the manifold

具體步驟包括：

1) 準(zhǔn)備分析信號；

2) 估測使用主頻分析的原子類型數(shù)目；

3) 根據(jù)式(1)～(5)提取不同稀疏特性的脈沖；

4) 根據(jù)式(9)構(gòu)建脈沖包絡(luò)空間；

5) 針對脈沖包絡(luò)空間展開流形降維，求解低維本征脈沖包絡(luò)；

6) 計(jì)算本征10包絡(luò)的包絡(luò)譜判斷軸承故障。

2 仿真驗(yàn)證

為了證明所提方法的有效性，設(shè)置如下仿真信號

(12)

其中：L為模擬脈沖數(shù)目，為39；Al為l的故障脈沖振幅，為5×10-9m；Tp為故障特征頻率fp的倒數(shù)；fp設(shè)為49.1 Hz；β為結(jié)構(gòu)阻尼系數(shù)，為1 200 Ns/m；ω為固有頻率，為2 000 Hz。

圖2表示仿真信號原始時(shí)域波形及信噪比為-7 dB的含白噪聲信號波形。

圖2 模擬信號Fig.2 Simulation signals

圖3(a)表示基于CSR-HT-Isomap方法提取的低維本征脈沖包絡(luò)，其中，懲罰因子λ=6∶1∶15，Isomap的近鄰點(diǎn)個(gè)數(shù)[11]取10，可以發(fā)現(xiàn)，該方法成功地提取了全部39條脈沖包絡(luò)。圖3(b)表示該低維本征脈沖包絡(luò)的包絡(luò)譜，可以發(fā)現(xiàn)該方法將故障特征頻率fp及其20次諧波均成功提取，證明了筆者所提方法的有效性。

圖3 CSR-HT-Isomap故障特征提取Fig.3 Fault feature extraction by the CSR-HT-Isomap method

現(xiàn)比較所提方法與傳統(tǒng)的基于EEMD和WPT包絡(luò)空間的Isomap特征提取效果，其中，EEMD中添加的白噪聲的幅值選擇為0.2σ，分解次數(shù)為400。WPT中母小波為Db8小波。以IMF和多尺度包絡(luò)空間為基礎(chǔ)的本征脈沖包絡(luò)如圖4所示。IMF包絡(luò)空間發(fā)現(xiàn)10次諧波的故障特征頻率，多尺度包絡(luò)空間發(fā)現(xiàn)5次諧波，而且通過基于IMF和多尺度包絡(luò)空間的包絡(luò)幅值和包絡(luò)譜幅值均小于圖3。這些表明，筆者所提方法要優(yōu)于IMF包絡(luò)和多尺度包絡(luò)空間。

圖4 本征包絡(luò)流形提取Fig.4 The intrinsic envelope extracted by the manifold with Isomap and their envelope spectra

3 臺(tái)架實(shí)驗(yàn)驗(yàn)證

為了獲得高速列車輪對軸承振動(dòng)信號，采用高速列車輪對跑合實(shí)驗(yàn)臺(tái)進(jìn)行臺(tái)架實(shí)驗(yàn)，如圖5所示。跑合實(shí)驗(yàn)臺(tái)由電機(jī)、驅(qū)動(dòng)輪、加載裝置、輪對及軸箱組成，實(shí)驗(yàn)中，對輪對軸承做了人工傷處理，設(shè)置了外圈故障與滾柱故障，外圈故障的深度、寬度皆為1 mm，長度為5 mm，滾柱故障的深度、寬度皆為1 mm。如圖6所示，加速度傳感器安置于軸箱上，軸承參數(shù)見表1。

圖5 輪對跑合實(shí)驗(yàn)臺(tái)及加速度傳感器Fig.5 The test bench and measurement sensor

圖6 高速列車輪對軸承人工傷Fig.6 Man-made fault photos of the wheelset of high speed train

滾子數(shù)量滾子直徑/mm節(jié)圓直徑/mm接觸角/rad轉(zhuǎn)速/(r·s-1)1926.91800.157 120.57

3.1 外圈故障檢測

當(dāng)輪對軸承外圈存在圖6(a)所示故障時(shí)，軸箱振動(dòng)信號如圖7所示。

圖7 外圈故障軸箱振動(dòng)信號Fig.7 Vibration signals of bearing outer race fault

外圈的故障特征頻率fBPFO表示為

(13)

其中:Nb為滾子數(shù)量;Bd為滾子直徑;Pd為節(jié)圓直徑;φ為接觸角;fw為轉(zhuǎn)速；經(jīng)計(jì)算fBPFO=166.6 Hz。

利用CSR-HT-Isomap方法提取的本征脈沖包絡(luò)如圖8所示，其中圖8(b)成功地發(fā)現(xiàn)了外圈的故障特征頻率及其7階諧波。

圖8 CSR-HT-Isomap故障特征提取Fig.8 Fault feature extraction by the CSR-HT-Isomap method

現(xiàn)比較筆者所提的基于CSR-HT-Isomap方法與基于傳統(tǒng)的EEMD和WPT包絡(luò)空間的Isomap特征提取效果，其中，EEMD中添加的白噪聲的幅值選擇為信號標(biāo)準(zhǔn)差，分解次數(shù)為400。WPT中母小波為Db8小波。以IMF和多尺度包絡(luò)空間為基礎(chǔ)的本征脈沖包絡(luò)流形提取如圖9所示。IMF包絡(luò)空間發(fā)現(xiàn)了故障信息中7次諧波的故障特征頻率，如圖9(b)所示，多尺度包絡(luò)空間僅僅發(fā)現(xiàn)2次諧波的故障特征頻率，如圖9(d)所示，而且通過基于IMF和多尺度包絡(luò)空間得到的包絡(luò)譜幅值均小于基于CSR-HT-Isomap方法提取的本征脈沖包絡(luò)空間包絡(luò)譜幅值，如圖8(b)、圖9(b)及圖9(d)所示。這些表明在提取由軸承外圈故障特征時(shí)，脈沖包絡(luò)空間的表現(xiàn)要優(yōu)于IMFs和多尺度包絡(luò)空間，無論是在放大故障特征頻率的諧波數(shù)方面，還是在強(qiáng)化本征包絡(luò)譜方面，前者均優(yōu)于后兩者。

圖9 本征包絡(luò)流形提取Fig.9 The intrinsic envelope extracted by the manifold with Isomap and their envelope spectra

3.2 滾子故障檢測

當(dāng)輪對軸承外圈存在圖6(b)所示故障時(shí)，軸箱振動(dòng)信號如圖10所示。

圖10 滾子故障軸箱振動(dòng)信號Fig.10 Vibration signals of a bearing roller fault

滾子的故障特征頻率fBSF表示為

(14)

計(jì)算后，fBSF=67.9 Hz。由CSR-HT-Isomap算法提取的軸承滾子故障振動(dòng)本征脈沖包絡(luò)如圖11所示。

圖11 CSR-HT-Isomap故障特征提取Fig.11 Fault feature extraction by the CSR-HT-Isomap method

根據(jù)圖12所示的滾子和滾道的運(yùn)動(dòng)關(guān)系。在輪對軸承滾子負(fù)載區(qū)(load zone，簡稱LZ)中有從脈沖B1經(jīng)E1到H1有7個(gè)脈沖。在圖11(a)中有7個(gè)脈沖在滾子負(fù)載區(qū)，在輪對軸承滾子非負(fù)載區(qū)域(none load zone，簡稱NLZ)，由于在滾子非負(fù)載區(qū)域缺乏接觸力導(dǎo)致作用在滾子故障處的力太小而無法激發(fā)脈沖。在脈沖B1到脈沖H1這7個(gè)脈沖中，脈沖B1和脈沖H1的位置在滾子負(fù)載區(qū)的邊界，當(dāng)有缺陷的滾子處于脈沖E1的位置上時(shí)，在滾道和有缺陷的滾子之間出現(xiàn)了一個(gè)缺口，相當(dāng)于脈沖F1和D1的兩個(gè)正常滾子形成的對稱結(jié)構(gòu)造成缺口的產(chǎn)生。所以脈沖B1，H1和E1的包絡(luò)邊緣小于脈沖C1，D1，F(xiàn)1和G1的包絡(luò)邊緣。通過CSR-HT-Isomap方法獲得的本征包絡(luò)成功地發(fā)現(xiàn)了滾子故障時(shí)軸承的振動(dòng)特點(diǎn)，同時(shí)本征包絡(luò)的包絡(luò)譜成功地顯示了故障特征頻率fBSF和它的2次諧波。仔細(xì)觀察還可發(fā)現(xiàn)，包絡(luò)譜的頻率呈現(xiàn)周期性波動(dòng)的特點(diǎn)，除第一部分頻率峰值圍繞在轉(zhuǎn)頻左右，其他的峰值皆圍繞在滾子故障特征頻率及其倍頻頻率附近，這為故障的識(shí)別提供了有力的證據(jù)。

圖12 滾子與滾道間的運(yùn)動(dòng)關(guān)系Fig.12 The relative motion between the rollers and races

比較基于IMFs包絡(luò)空間和多尺度包絡(luò)空間的本征脈沖包絡(luò)如圖13所示?？梢园l(fā)現(xiàn)，其包絡(luò)空間無法體現(xiàn)滾子故障時(shí)滾子處于不同位置時(shí)的軸承振動(dòng)特點(diǎn)，包絡(luò)譜的幅值較小，且故障特征頻率淹沒在噪聲頻率中，無法反應(yīng)軸承故障狀態(tài)。

圖13 本征包絡(luò)流形提取Fig.13 The intrinsic envelope extracted by the manifold with Isomap and their envelope spectra:

4 結(jié)束語

筆者提出了一種基于CSR-HT-Isomap相結(jié)合的故障檢測方法，通過仿真數(shù)據(jù)以及臺(tái)架實(shí)驗(yàn)數(shù)據(jù)表明：該輪對軸承故障診斷方法可以很好地提取仿真信號、軸承內(nèi)圈及滾動(dòng)體故障特征，通過與基于聚合經(jīng)驗(yàn)?zāi)B(tài)分解和小波包變換的包絡(luò)空間算法進(jìn)行比較，證明本方法在提取本征包絡(luò)、強(qiáng)化本征包絡(luò)譜以及放大故障特征頻率的諧波數(shù)方面具備較大優(yōu)勢，具有一定的工程應(yīng)用前景。