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

(1.四川理工學(xué)院 自動(dòng)化與電子信息學(xué)院，四川 自貢 643000；2.人工智能四川省重點(diǎn)實(shí)驗(yàn)室，四川 自貢 643000)

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

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

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

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

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

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

(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 多重分形譜算法的改進(jìn)

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

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

(4)

則第i個(gè)子集的概率測(cè)度為

(5)

定義配分函數(shù)為

(6)

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

(7)

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

(8)

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

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

(9)

通過(9)式測(cè)定概率測(cè)度，并計(jì)算質(zhì)量指數(shù)和配分函數(shù)，就可得具有分形結(jié)構(gòu)的多重分形譜指數(shù)α與f(α)。該算法雖無(wú)法得到分維的準(zhǔn)確數(shù)值，但可得到奇異強(qiáng)度和對(duì)應(yīng)的分布密度。當(dāng)N=ms+h且h=0時(shí)，即序列所分子集后無(wú)剩余時(shí)，改進(jìn)算法與原算法相同。

2 軸承振動(dòng)信號(hào)的采集

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

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

圖1 正常狀態(tài)振動(dòng)信號(hào)時(shí)域波形圖

圖2 故障狀態(tài)下內(nèi)圈振動(dòng)信號(hào)的時(shí)域波形

3 試驗(yàn)及分析

3.1 Hurst指數(shù)分析

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

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

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

如果振動(dòng)信號(hào)具有Brown運(yùn)動(dòng)特性，給定一個(gè)從-t到0的過去增量Z0-Z-t，與未來增量Zt-Z0的相關(guān)系數(shù)為

(10)

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

(11)

則可使用H的大小表征振動(dòng)信號(hào)的統(tǒng)計(jì)特性，即：

l)當(dāng)0

2)當(dāng)H=0.5時(shí)，b(t)=0。表明信號(hào)不相關(guān)，信號(hào)為標(biāo)準(zhǔn)的隨機(jī)信號(hào)，即過去與將來不存在相關(guān)性，振動(dòng)信號(hào)為完全獨(dú)立過程。

3)當(dāng)0.5

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

(12)

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

Vs由上升轉(zhuǎn)為顯著常數(shù)或下降的分界點(diǎn)，即為振動(dòng)長(zhǎng)記憶的消失點(diǎn)，相應(yīng)s就是信號(hào)去相關(guān)性的時(shí)刻。各狀態(tài)振動(dòng)信號(hào)的Hurst指數(shù)和記憶長(zhǎng)度的計(jì)算結(jié)果見表1。

表1 不同狀態(tài)振動(dòng)信號(hào)的Hurst指數(shù)和記憶長(zhǎng)度

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

3.2 多重分形性分析

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

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

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

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

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

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

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

圖8 q-α曲線

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

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

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

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

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

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

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

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

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

4 結(jié)論

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

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

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

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