李兆飛，任小洪，譚飛，方寧

(1.四川理工學院 自動化與電子信息學院，四川 自貢 643000；2.人工智能四川省重點實驗室，四川 自貢 643000)

隨著對機械設備振動信號非線性機理的認識逐漸加深，許多國內(nèi)外學者用混沌法等多種非線性方法對實際振動信號進行混沌特性的判斷，取得了大量的理論和成果。通常直接計算實測振動信號的非線性特征量來分析設備的非線性特性，多采用關聯(lián)維數(shù)[1-2]、小波包變換頻帶盒維數(shù)[3]、多重分形特征[4-5]、Lyapunov指數(shù)、Kolmogorov熵及相關復雜性測度等方法對振動特性進行研究。然而，一個特征量只能反映振動某方面的特性，而且受信號長度及噪聲影響，使用這些特征量的計算結果會有所偏差，導致判斷錯誤。

研究振動分形特性有助于認識球軸承的內(nèi)在變化規(guī)律，并利用這些規(guī)律對其進行診斷和監(jiān)測。因此，研究了球軸承振動信號的非線性分形特性，通過改進Hurst指數(shù)和多重分形譜的計算方法，對球軸承不同部位及同一部位不同故障狀態(tài)實測振動信號的分形特性進行了研究。

1 Hurst指數(shù)和多重分形譜的計算方法及改進

分形即某些部分以某種方式與其整體相似的集合，分形集就是在仿射變換下某些部分(互不重疊的子集)與該集本身(統(tǒng)計或嚴格地)相通的集合[6]。機械設備振動信號的分形特性主要是振動信號自相似性(統(tǒng)計上)，長期記憶性和狀態(tài)持續(xù)性，以及多重分形性(不均勻的自相似性)。

研究振動信號分形特性主要是計算振動信號的分形特征量。在此，對Hurst指數(shù)和多重分形譜的計算方法加以改進，用以計算振動信號各狀態(tài)下的Hurst指數(shù)和記憶長度及多重分形譜。

1.1 Hurst指數(shù)計算方法的改進

(1)

(2)

取不同的s，可得到不同的(lns,ln(R/S)s)，采用最小二乘法擬合得

ln(R/S)s=Hlns+lnb,

(3)

式中：H為Hurst指數(shù)，0≤H≤1；b為常數(shù)。

1.2 多重分形譜算法的改進

實際工程振動信號應該只在一定的標度范圍內(nèi)才存在分形特性，且不同的局部條件將得到不同局部特征，即表現(xiàn)出多標度分形，也稱多重分形[8]，是對分形結構不規(guī)則和不均勻程度的度量[9]。多重分形奇異譜[10]通過描述定義于分形集上的歸一化測度(分布)的標度性來表現(xiàn)其奇異特性。測度的奇異性由奇異性強度α和相應的測度密度f(α)來表征，分形奇異譜由f(α)及α的取值范圍確定。實際上，一般的集合通常不會得到解析解，故計算多重分形譜較復雜，往往通過數(shù)值算法間接估計。

同樣，原矩法估計多重分形譜在將信號xn(n=1,2,…,N)分劃成長度為ε(ε<1)的一維小盒子，各個盒子的長度為s，N不一定能整除s，信號將有一段剩余,使得計算的結果不夠準確。因此，對算法作如下的改進：

(4)

則第i個子集的概率測度為

(5)

定義配分函數(shù)為

(6)

質(zhì)量指數(shù)為

(7)

則由τ(q)估計得到廣義分維數(shù)為

(8)

該分維數(shù)隨q值不同意義不同。在振動均勻分布時，Dq=1。當q=0時，D0=-τ(0)=D，這表明q=0時，其與物理量的不均勻概率分布p無關，描述對象的空間幾何特征，D0稱作豪斯道夫維。當q=1時，x1(ε)=∑Pi=1，因此有τ(1)=0，這時D1稱為信息維。當q=2時，D2稱關聯(lián)維(其變化表示在數(shù)據(jù)集中點分布的變化，在特定尺度范圍任意選擇兩點的概率)。

最后，計算α和f(α)的公式不變，由統(tǒng)計物理中的勒讓德變換，建立作為描述同一物理對象的3個標度指數(shù)α，f(α)和τ(q)之間的聯(lián)系，即

(9)

通過(9)式測定概率測度，并計算質(zhì)量指數(shù)和配分函數(shù)，就可得具有分形結構的多重分形譜指數(shù)α與f(α)。該算法雖無法得到分維的準確數(shù)值，但可得到奇異強度和對應的分布密度。當N=ms+h且h=0時，即序列所分子集后無剩余時，改進算法與原算法相同。

2 軸承振動信號的采集

試驗數(shù)據(jù)來自Case West Reserve University滾動軸承試驗數(shù)據(jù)庫，包含內(nèi)、外圈和鋼球不同轉(zhuǎn)速、不同故障程度的實測振動數(shù)據(jù)，對這些數(shù)據(jù)進行試驗分析，探索不同故障類型及不同程度故障信號的分形特性。

試驗軸承為深溝球軸承6205-2RS JEM SKF，電動機空載，采樣頻率fs=12 kHz，轉(zhuǎn)速為1 797 r/min，軸承故障由電火花加工機在球軸承內(nèi)、外圈及鋼球上模擬，故障直徑分別為0.18，0.36，0.54和0.72 mm(圖中對應表示為故障1～4)，故障深度為0.28 mm。正常狀態(tài)和各故障狀態(tài)振動信號的時域波形如圖1和圖2所示。

圖1 正常狀態(tài)振動信號時域波形圖

圖2 故障狀態(tài)下內(nèi)圈振動信號的時域波形

3 試驗及分析

3.1 Hurst指數(shù)分析

正常和故障狀態(tài)下振動信號的Hurst指數(shù)如圖3和圖4所示，圖中擬合直線的斜率即為Hurst指數(shù)(正常時H=0.680 55，故障1時H=0.309 79，故障2時H=0.452 47，故障3時H=0.391 01，故障4時H=0.321 58)，振動信號的Hurst指數(shù)在0～1之間。

圖3 正常狀態(tài)Hurst指數(shù)的計算

圖4 故障狀態(tài)Hurst指數(shù)的計算

如果振動信號具有Brown運動特性，給定一個從-t到0的過去增量Z0-Z-t，與未來增量Zt-Z0的相關系數(shù)為

(10)

因E[(Z0-Z-t)(Zt-Z0)]2=E[(Zt-Z-t)]2，代入上式整理有

(11)

則可使用H的大小表征振動信號的統(tǒng)計特性，即：

l)當0

2)當H=0.5時，b(t)=0。表明信號不相關，信號為標準的隨機信號，即過去與將來不存在相關性，振動信號為完全獨立過程。

3)當0.5

Vs統(tǒng)計量最初用來檢驗R/S分析的穩(wěn)定性[7]，經(jīng)改進后被用來分析序列的記憶長度，其計算式為

(12)

為保證Hurst指數(shù)的魯棒性，在計算大樣本數(shù)據(jù)時，還需去除數(shù)據(jù)的馬爾可夫短期記憶性。結果在數(shù)據(jù)長度大于10時，馬爾可夫短期記憶被成功去除，也可通過Vs統(tǒng)計出非周期循環(huán)的精確度量。

Vs由上升轉(zhuǎn)為顯著常數(shù)或下降的分界點，即為振動長記憶的消失點，相應s就是信號去相關性的時刻。各狀態(tài)振動信號的Hurst指數(shù)和記憶長度的計算結果見表1。

表1 不同狀態(tài)振動信號的Hurst指數(shù)和記憶長度

由表可知，正常狀態(tài)振動信號的Hurst指數(shù)大于0.5，而各故障狀態(tài)振動信號的Hurst指數(shù)均小于0.5，說明正常狀態(tài)振動信號存在持續(xù)的長程相關性，而故障狀態(tài)振動信號存在反持續(xù)性；正常狀態(tài)的記憶長度比故障狀態(tài)的大很多，說明正常狀態(tài)振動的平均影響范圍比故障狀態(tài)大，這是因為正常狀態(tài)的振動平穩(wěn)有序，不會產(chǎn)生較大的沖擊振動，會影響到較大的范圍。Hurst及記憶長度估計了不同狀態(tài)振動波動的影響范圍。

3.2 多重分形性分析

由改進算法計算正常及故障狀態(tài)下振動信號的多重分形譜，根據(jù)前面的分析，尺寸ε取值為2，4，…,10時，取q=-50～50，步長取1，分析繪制軸承不同狀態(tài)的相關圖像，結果如圖5和圖6所示。

圖5 正常狀態(tài)信號ln xq(ε)-ln ε擬合直線

圖6 故障狀態(tài)信號ln xq(ε)-ln ε擬合直線

從圖中可以看出，在q>0或q<0的取值區(qū)間上，球軸承不同振動狀態(tài)配分函數(shù)與測量單元尺寸的自然對數(shù)均有較好的線性關系對應，也就是符合標度不變性。說明振動的標度分布為多重分形。而利用單一分形或正態(tài)分布對其進行刻畫是不合理的。

τ(q)與q的關系如圖7所示，從圖中可以看出，正常狀態(tài)的線性關系很好，而故障狀態(tài)時，τ(q)為上凸函數(shù)，表明τ(q)與q間為非線性關系，說明軸承故障狀態(tài)有更強的多重分形特征。

圖7 正常狀態(tài)及各故障狀態(tài)信號的τ(q)-q關系圖

α與q的曲線關系如圖8所示，故障狀態(tài)越強，α與q的非線性越強，且α是q的非線性減函數(shù)。因此，當q為極大值時，α為最小值αmin=α(+∞)。由(6)式可知，xq(ε)由Pi(ε)的最大值決定，即αmin反映最大速度。再由ε<1及f(α)與Nα(ε)的正比關系可知，f(αmin)反映最大幅值的出現(xiàn)次數(shù)。相應的，αmax反映最小速度，f(αmax)反映最小幅值的出現(xiàn)次數(shù)。

圖8 q-α曲線

廣義分形維數(shù)的函數(shù)關系如圖9所示。從圖中可以看出，正常狀態(tài)時Dq幾乎為常數(shù)1，表明振動均勻分布；而故障狀態(tài)時Dq不為常數(shù)1，表示振動非均勻分布。當q=0時，所有的廣義分維數(shù)都是1。隨q增大而單調(diào)下降曲線為故障狀態(tài)廣義分維數(shù)，也存在多重分形特性，廣義分維數(shù)在常數(shù)1周圍波動越大則故障越嚴重，表明故障狀態(tài)與重分形性存在對應關系。

圖9 廣義分形維數(shù)Dq-q關系圖

軸承振動信號的多重分形譜如圖10所示，為開口向下的一條拋物線。由重分形譜的分析，多重分形譜寬Δα=αmax-αmin表征了所有分形結構的不規(guī)則分布度，描述了振動幅度大??；Δf=fmax-fmin表征了物理量生成子集中元素個數(shù)在最小、最大處的比例，同理，可反映振動信號大、小峰值占有的比例。fmax有相同概率的單元數(shù)和ε的變化關系則反映了振動大、小峰值的速度變化及所占的比例，可以反映振動的劇烈程度。

圖10 多重分形譜f(α)-α關系圖

由圖10可以看出：不同狀態(tài)時fmax均為1，這是由于軸承的轉(zhuǎn)速相同，具有相同的變化速率，所有的區(qū)域均有分布，即多重分形譜圖fmax相同。而多重分形譜中5條曲線的寬度及形狀各異，則反映了軸承相異狀態(tài)的內(nèi)在動力學特性，不同故障狀態(tài)間的分形特性和復雜動力學行為及分形特征的差異可由多重分形反映。

不同狀態(tài)振動信號多重分形譜的Δα和Δf計算結果見表2。Δα越大，故障程度越大。表示在此標度內(nèi)，較大故障的振動信號有較顯著的多重分形特性，說明故障狀態(tài)下振動分布不均勻性最大，此時振動比較劇烈。正常狀態(tài)球軸承振動的多重分形性最小，變化波動的奇異性也最弱，說明正常狀態(tài)下振動比較均勻。且故障程度較強振動重分形譜圖呈左鉤狀，即Δf>0，奇異譜的頂部相對較圓潤，振動大多在波峰。表明指數(shù)在小幅值的機會比處于大幅值的幾率小。

由表2可知，不同故障狀態(tài)有相應的Δα和Δf，表明振動狀態(tài)與分形譜參數(shù)間有一定的對應關系，不同振動狀態(tài)擁有不同的分形特性，Δα及Δf不光進行了定性的說明還進行了定量的比較。

表2 不同狀態(tài)振動信號的Δα和Δf值

另外，對軸承外圈和鋼球不同故障狀態(tài)進行了相同的分析，結果表明：不同故障程度時振動信號的長程相關性及分布的均勻性具有完全相似的結果，正常狀態(tài)Hurst指數(shù)均大于0.5，故障狀態(tài)Hurst指數(shù)均小于0.5，正常狀態(tài)記憶長度均大于故障狀態(tài)；振動狀態(tài)與振動多重分形譜之間存在對應關系，只是多重分形譜寬及振動的大、小峰值的增長速度和所占的比例不同而已。

4 結論

1)研究了滾動軸承振動信號的非線性分形特征，改進了Hurst指數(shù)和多重分形譜2種分形特征算法，通過2次分割充分利用了剩余數(shù)據(jù)。

2)通過Hurst指數(shù)研究了滾動軸承振動信號的長程相關性，其正常狀態(tài)振動信號具有正的持續(xù)性，故障狀態(tài)振動信號具有反持續(xù)性。通過記憶長度確定振動的影響范圍，故障狀態(tài)振動信號的影響波動范圍要比正常狀態(tài)振動信號的小。

3)通過多重分形譜分析滾動軸承振動信號的均勻性。振動狀態(tài)與其多重分形譜之間存在對應關系，即振動的劇烈度由多重分形譜寬Δα刻畫，振動大、小峰值的增長速度及所占比例由最大、最小概率子集分維差Δf刻畫。正常狀態(tài)的Δα和Δf最小，故障狀態(tài)越大，Δα和Δf越大。

另外，有關系統(tǒng)狀態(tài)對滾動軸承振動多重分形特性內(nèi)在機理的影響還需要更多的實證研究分析，上述研究結果為下一步利用軸承振動信號分形方法進行故障特征提取、診斷及檢測奠定了基礎。