不確定性激勵下碰摩轉(zhuǎn)子的振動響應(yīng)識別

王 威， 甘春標(biāo)

(浙江大學(xué) 機(jī)械工程學(xué)院，杭州 310027)

工程實(shí)際中，由于載荷、摩擦及其他因素的影響，多數(shù)旋轉(zhuǎn)機(jī)械系統(tǒng)呈現(xiàn)出非線性特性。特別是當(dāng)系統(tǒng)處于故障狀況下，這類非線性系統(tǒng)的振動特性變得更加復(fù)雜，呈現(xiàn)出復(fù)雜的動力學(xué)特性甚至混沌特性；此外，轉(zhuǎn)子系統(tǒng)在實(shí)際運(yùn)行過程中，不確定性擾動的影響是客觀存在的[1]。因此，對此類旋轉(zhuǎn)機(jī)械轉(zhuǎn)子系統(tǒng)的雜亂響應(yīng)信號進(jìn)行分析及識別，有效的判斷非線性特性及混沌特性，對于故障診斷有著重要意義。

現(xiàn)有的信號處理方法如小波分析[2]、時頻分析[3]等被廣泛應(yīng)用于故障診斷領(lǐng)域，用以對系統(tǒng)非線性信號進(jìn)行分析，但其對系統(tǒng)混沌特性缺乏足夠的反映。針對此問題，非線性動力學(xué)方法可以有效的克服其缺點(diǎn)，現(xiàn)有的判別混沌特性的方法主要有兩類：一類是針對信號的直接判別方法，如Poincaré截面[4]、最大Lyapunov指數(shù)[5]、相關(guān)維數(shù)等，但其對數(shù)據(jù)要求較高，且結(jié)果受噪聲的影響較大[6]；另一類是針對信號的間接判別方法，主要通過系統(tǒng)的非線性檢驗(yàn)對系統(tǒng)混沌特性進(jìn)行判別。替代數(shù)據(jù)方法作為間接判別方法的一種，由Theiler等[7]于1992年首次提出，由于有效地避免直接判別法的局限性，其后被廣泛應(yīng)用于各領(lǐng)域，用于對混沌時間序列進(jìn)行分析[8-11]。Small等[12]針對具有明顯周期特征的時間序列動力特性的特點(diǎn)，發(fā)展了替代數(shù)據(jù)方法，提出了一種偽周期替代數(shù)據(jù)算法(Pseudo-Periodic Surrogate, PPS),通過對R?ssler系統(tǒng)的響應(yīng)信號和人體心電圖的實(shí)驗(yàn)數(shù)據(jù)進(jìn)行分析，證明了該算法可以有效的區(qū)分混沌與噪聲周期(noisy periodic)信號。該算法已被廣泛應(yīng)用于動力學(xué)系統(tǒng)的非線性特性分析[13-15],但較少地應(yīng)用于故障診斷領(lǐng)域[16]。

轉(zhuǎn)子系統(tǒng)的碰摩故障是一種強(qiáng)非線性現(xiàn)象，系統(tǒng)運(yùn)動會出現(xiàn)跳躍、分岔以及混沌等非線性動力學(xué)行為。國內(nèi)外學(xué)者對該類系統(tǒng)進(jìn)行了深入研究，很好地解釋了碰摩轉(zhuǎn)子系統(tǒng)中一些復(fù)雜的非線性動力學(xué)現(xiàn)象[17-20]。Adams等研究了碰摩轉(zhuǎn)子系統(tǒng)的混沌運(yùn)動，通過改變間隙，獲得了豐富的次諧波、擬周期和混沌運(yùn)動。羅躍綱等利用分岔圖對含有裂紋-碰摩耦合故障的雙跨轉(zhuǎn)子系統(tǒng)的復(fù)雜非線性運(yùn)動進(jìn)行了分析，并討論了裂紋、碰摩及其耦合故障對系統(tǒng)非線性動力學(xué)響應(yīng)的影響。在工程實(shí)際中，不確定性激勵對轉(zhuǎn)子系統(tǒng)的影響已越來越被國內(nèi)外學(xué)者所重視[21]。Leng等[22]對隨機(jī)激勵作用下裂紋故障轉(zhuǎn)子系統(tǒng)的非線性特性進(jìn)行了研究，并發(fā)現(xiàn)隨機(jī)激勵在擬周期解和混沌解的分岔區(qū)間對轉(zhuǎn)子系統(tǒng)的振動響應(yīng)影響比較顯著，且隨機(jī)激勵越大，其影響也越顯著。高喆等[23]以Jeffcort轉(zhuǎn)子為模型，研究了碰摩轉(zhuǎn)子-軸承系統(tǒng)的隨機(jī)分岔與混沌特性，分析了白噪聲對轉(zhuǎn)子系統(tǒng)響應(yīng)的影響。游震洲等[24]對含碰摩故障的轉(zhuǎn)子-基礎(chǔ)系統(tǒng)受外部隨機(jī)激勵作用下動力學(xué)特性進(jìn)行了研究。分析了隨機(jī)外激勵強(qiáng)度和譜寬的變化以及周期外激勵頻率的改變對系統(tǒng)振動模式的影響?，F(xiàn)有的研究主要集中于應(yīng)用直接判別方法對系統(tǒng)混沌特性進(jìn)行分析，實(shí)際工程應(yīng)用中的信號數(shù)據(jù)因受到不確定激勵的影響，通常難以滿足其要求，而替代數(shù)據(jù)方法可以有效的克服此不足。

本文應(yīng)用偽周期替代數(shù)據(jù)算法，對受不確定性激勵下碰摩故障轉(zhuǎn)子系統(tǒng)的振動響應(yīng)進(jìn)行識別和分類。首先介紹由Small等提出的偽周期替代數(shù)據(jù)方法；其次，建立碰摩故障Jeffcott轉(zhuǎn)子系統(tǒng)動力學(xué)模型以及碰摩力與油膜力模型，并引入有界噪聲表征不確定性激勵形式；然后利用數(shù)值方法進(jìn)行分析，選擇相關(guān)維數(shù)作為判別統(tǒng)計(jì)量對系統(tǒng)輸出信號的混沌特性進(jìn)行判別，同時應(yīng)用原始信號的Poincaré截面圖及最大Lyapunov指數(shù)對判別結(jié)果進(jìn)行驗(yàn)證；最后對本文工作進(jìn)行總結(jié)。

1 偽周期替代數(shù)據(jù)算法

偽周期替代數(shù)據(jù)算法(PPS)能夠檢測具有明顯周期特征的時間序列的動力特征。PPS算法是建立在局部建模技術(shù)的基礎(chǔ)上，通過引入噪聲來生成替代時間序列。生成的替代時間序列與原時間序列的偏差性取決于噪聲引入的級別。

PPS算法的基本算法步驟如下：

式中： 參數(shù)ρ為噪聲半徑。

步驟4令si+1=sj+1， 并且i=i+1。

步驟5重復(fù)步驟3和步驟4直到滿足i>N。

通過以上算法，雖然引入噪聲干擾，但得到的替代數(shù)據(jù)序列保持了原時間序列隱藏的動力學(xué)特征。其中參數(shù)ρ的選取對替代數(shù)據(jù)能否保留原序列隱藏的動力學(xué)特征有著重要影響。當(dāng)參數(shù)ρ的取值過小，則得到的替代序列與原序列幾乎一樣；當(dāng)參數(shù)ρ的取值過大，則產(chǎn)生的替代數(shù)據(jù)近似于獨(dú)立同分布的噪聲數(shù)據(jù)。因此本文依據(jù)Small等關(guān)于參數(shù)ρ的取值法則，即令原始時間序列與替代數(shù)據(jù)的長度2的短節(jié)段的期望值最大。

2 轉(zhuǎn)子動力學(xué)模型

本文建立的偏心-碰摩耦合故障Jeffcott轉(zhuǎn)子系統(tǒng)如圖1所示，在模型中，軸承采用滑動軸承，轉(zhuǎn)子考慮為集中參數(shù)模型，軸承力以集中力的方式施加于轉(zhuǎn)子兩端，轉(zhuǎn)盤與轉(zhuǎn)子之間也通過集中約束力來實(shí)現(xiàn)連接。O1,O2,O3為軸承幾何中心、轉(zhuǎn)子幾何中心、轉(zhuǎn)子質(zhì)心；δ為轉(zhuǎn)子圓盤和定子的間隙；k,c為轉(zhuǎn)軸與圓盤之間的連接剛度和阻尼；m2為轉(zhuǎn)子圓盤質(zhì)量；m1為左右滾動軸承外圈質(zhì)量；kc為碰摩接觸剛度；e為偏心量；Fx,Fy為左端軸承的支承反力；Px,Py分別為碰摩力對圓盤的x和y向作用力；frx,fry分別為作用于圓盤x和y向的不確定性激勵。

圖1 Jeffcott轉(zhuǎn)子模型Fig.1 The Jeffcott rotor model

本文動力學(xué)模型考慮滑動軸承油膜力、全周碰摩力與不確定性激勵作用，并考慮不平衡偏心的影響。由此，通過Lagrange方程可得轉(zhuǎn)子運(yùn)動方程，其模型為

(1)

對式(1)進(jìn)行無量綱化，令

則原方程可以化為

(2)

以下分別對軸承油膜力、全周碰摩力以及隨機(jī)激勵模型進(jìn)行分析。

2.1 油膜力

本文選用短軸承油膜力模型，其無量綱表達(dá)式為[26]

(3)

式中：x,y為軸承無量綱位移；fx,fy分別為無量綱油膜力在x,y方向的分量，重疊點(diǎn)表示該變量的時間導(dǎo)數(shù)，其中

2.2 全周碰摩力

設(shè)靜止時轉(zhuǎn)子和定子之間的間隙為δ,kc為定子徑向剛度，并且轉(zhuǎn)子和定子間的摩擦符合庫倫摩擦定律，摩擦因數(shù)為f, 則碰摩時徑向碰撞力PN及切向摩擦力PT分別可以表達(dá)為

PN=(d-δ)kc,PT=f·PN
(d≥δ)

(4)

(5)

式中：Px,Py分別為碰摩力在x,y方向的分量。

圖2 全周碰摩模型Fig.2 The schematic diagram of the full annular rub

將式(4)代入式(5)得

(6)

則碰摩力在x-y坐標(biāo)系上相對于無量綱位移X2,Y2的表達(dá)式

(7)

2.3 不確定性激勵

在實(shí)際工況中，旋轉(zhuǎn)機(jī)械系統(tǒng)通常會受到不確定性激勵的作用，如來自地面的傳動振動激勵、氣流引起的非平穩(wěn)激勵、由工件表面粗糙度引起的不確定性切削力等。本文引入的不確定性激勵直接作用于轉(zhuǎn)子圓盤以模擬不確定性切削力，并假設(shè)x和y向受到的激勵相同，即：frx=fry=fran。

以下應(yīng)用有界噪聲表征系統(tǒng)所受到的不確定性激勵，并假設(shè)有界噪聲信號可表示為

ξ(t)=Aisin(Ωit+σiBi(t)+γi)

(8)

式中：Ai為噪聲幅值；Bi(t)為單位Wiener過程；γi為[0,2π)上均勻分布的隨機(jī)變量；Ωi為中心頻率；σi為常數(shù)。有界噪聲的均值為零，其相關(guān)函數(shù)為式

(9)

其雙邊功率譜密度為

(10)

由文獻(xiàn)[27]可知， 有界噪聲的任一物理可識化樣本η(t)可以由式(11)逼近

(11)

圖3為通過諧波疊加法模擬出來的有界噪聲。根據(jù)Shinozuka等的研究，當(dāng)模擬點(diǎn)數(shù)大于一定值(N>20 000)時，模擬出來的有界噪聲的功率譜與理論值基本相同。

因此，可假設(shè)不確定性激勵fran為

fran(t)=ε0f0ξ(t)

(12)

式中：f0為常數(shù)； 0<ε0?1為小量。

圖3 有界噪聲樣本的時間歷程Fig.3 A sample time series of the bounded noise

3 振動響應(yīng)識別

在建立耦合系統(tǒng)模型和各故障模型的基礎(chǔ)上，本文利用定步長Runge-Kutta數(shù)值積分方法獲取了系統(tǒng)的非線性動力響應(yīng)，計(jì)算步長取2π/1 000，為使系統(tǒng)得到穩(wěn)定解，略去前300個周期，取后面的300個周期進(jìn)行分析。

在未特別說明時，系統(tǒng)模型參數(shù)取值為

m1=4.0 kg,m2=32.1 kg,R=25 mm,L=12 mm,c=0.11 mm,μ=0.018 Pa·s,c1=1 050 N·s/m,c2=2 100 N·s/m,k=2.5×107N/m,e=0.05 mm,kc=3.5×106N/m,δ=0.11 mm,f=0.1,ε0f0=0.1,σ=0.01,A0=0.005,Ω0=0.2,σ0=0.8,ωl=0.2,ωr=2 000,N=20 000。

取系統(tǒng)轉(zhuǎn)速為分岔參數(shù)，不確定性激勵強(qiáng)度系數(shù)分別取值σ=0，σ=0.01， 其余參數(shù)取給定參數(shù)畫x方向振動量的分岔圖，由圖4(b)可以看出，不確定性激勵作用下系統(tǒng)的動力學(xué)特性非常復(fù)雜，經(jīng)歷了周期占優(yōu)、倍周期占優(yōu)和混沌占優(yōu)等多個狀態(tài)，系統(tǒng)響應(yīng)非線性特性明顯。轉(zhuǎn)子轉(zhuǎn)速較低時，不確定性激勵對分岔區(qū)域的影響更為明顯。

圖4 系統(tǒng)振動響應(yīng)隨轉(zhuǎn)速變化的分岔圖Fig.4 The bifurcation diagram of the vibration response of the rotor under the change of rotating speed

應(yīng)用偽周期替代數(shù)據(jù)算法對系統(tǒng)的隨機(jī)動力學(xué)特性進(jìn)行分析。圖5～圖7中，分別對式(2)求解獲得的系統(tǒng)樣本時間序列及其30組替代數(shù)據(jù)樣本進(jìn)行關(guān)聯(lián)維數(shù)的計(jì)算，其中關(guān)聯(lián)維數(shù)的計(jì)算參考Judd[28]提出的算法，根據(jù)文獻(xiàn)中提出的算法，將關(guān)聯(lián)維數(shù)dc看作是觀察尺度系數(shù)ε0的函數(shù)，這樣求解的關(guān)聯(lián)維數(shù)為一條曲線，而不只是一個數(shù)值。同時應(yīng)用Rosenstein等[29]提出的最大Lyapunov指數(shù)算法對原始樣本進(jìn)行分析，用以驗(yàn)證所分析結(jié)果的有效性。

圖5為轉(zhuǎn)速ω取800 rad/s時隨機(jī)激勵下系統(tǒng)動力學(xué)響應(yīng)，由圖5(a)和圖5(b)可以看出，系統(tǒng)呈現(xiàn)倍周期占優(yōu)；對其進(jìn)行替代數(shù)據(jù)檢驗(yàn)，取噪聲半徑ρ=0.001，圖5(d)為原樣本與替代數(shù)據(jù)的關(guān)聯(lián)維數(shù)曲線的比較，由圖5(d)可以看出，原樣本與替代數(shù)據(jù)樣本關(guān)聯(lián)維數(shù)曲線不能被區(qū)分，含噪周期序列的零假設(shè)檢驗(yàn)不能被否定。

圖5 周期占優(yōu)振動響應(yīng)(ω=800 rad/s)Fig.5 The periodic-dominant vibration response (ω=800 rad/s)

圖6為轉(zhuǎn)速ω取1 200 rad/s時隨機(jī)激勵下系統(tǒng)動力學(xué)響應(yīng)，由圖6(a)和圖6 (b)可以看出，系統(tǒng)進(jìn)入混沌狀態(tài)，系統(tǒng)Poincaré截面呈現(xiàn)兩個混沌小島；由圖6(c)可以看出系統(tǒng)響應(yīng)最大Lyapunov指數(shù)為正值。對其進(jìn)行替代數(shù)據(jù)檢驗(yàn)，取噪聲半徑ρ=0.005， 圖6(d)為原樣本與替代數(shù)據(jù)的關(guān)聯(lián)維數(shù)曲線的比較，由圖6(d)可以看出，原樣本與替代數(shù)據(jù)樣本關(guān)聯(lián)維數(shù)曲線被顯著區(qū)分，其零假設(shè)檢驗(yàn)應(yīng)被否定。PPS算法可以有效的對系統(tǒng)混沌特性進(jìn)行識別。

圖7為轉(zhuǎn)速ω取500 rad/s時隨機(jī)激勵下系統(tǒng)動力學(xué)響應(yīng)，這里為獲得隨機(jī)響應(yīng)，取偏心量e=0，系統(tǒng)僅受隨機(jī)外激勵作用；由圖7(a)～圖7(c)可以看出，系統(tǒng)呈現(xiàn)隨機(jī)占優(yōu)特性，其最大Lyapunov指數(shù)在初始階段產(chǎn)生跳變，其最大值趨向于無窮大。對其進(jìn)行替代數(shù)據(jù)檢驗(yàn)，取噪聲半徑ρ=0.005，圖7(d)為原樣本與替代數(shù)據(jù)的關(guān)聯(lián)維數(shù)曲線的比較，由圖7(d)可以看出，原樣本與替代數(shù)據(jù)樣本關(guān)聯(lián)維數(shù)曲線不能被區(qū)分，隨機(jī)時間序列的零假設(shè)檢驗(yàn)不能被否定，隨機(jī)時間序列不能被識別。

綜上分析，應(yīng)用偽周期替代數(shù)據(jù)算法可以有效的判定系統(tǒng)振動響應(yīng)的周期占優(yōu)與混沌占優(yōu)特性，對不確定性激勵下碰摩故障轉(zhuǎn)子系統(tǒng)的振動響應(yīng)的特征識別與分類有較好的適用性，可以很好的應(yīng)用于對工程實(shí)際中含噪信號的混沌判別。

圖6 混沌占優(yōu)振動響應(yīng)(ω=1 200 rad/s)Fig.6 The chaotic-dominant vibration response (ω=1 200 rad/s)

圖7 隨機(jī)占優(yōu)振動響應(yīng)(ω=500 rad/s)Fig.7 The random-dominant vibration response (ω=500 rad/s)

4 結(jié) 論

工程實(shí)際中，故障轉(zhuǎn)子系統(tǒng)呈現(xiàn)出復(fù)雜的非線性動力學(xué)特性，并且隨機(jī)擾動的影響是客觀存在的，對不確定性激勵下故障轉(zhuǎn)子系統(tǒng)的非線性特性進(jìn)行分析，有效的判斷系統(tǒng)動力學(xué)特性，對于故障診斷有著重要意義。針對此問題，本文建立了一類受不確定性激勵的碰摩故障轉(zhuǎn)子系統(tǒng)的動力學(xué)模型，并應(yīng)用偽周期替代數(shù)據(jù)法對其非線性特性進(jìn)行分析，對故障信號進(jìn)行了識別與分類。

由仿真結(jié)果分析可知，不確定性激勵下碰摩轉(zhuǎn)子系統(tǒng)動力學(xué)特性復(fù)雜，經(jīng)歷了周期占優(yōu)、混沌占優(yōu)等多個狀態(tài)，系統(tǒng)響應(yīng)非線性特性明顯。轉(zhuǎn)子轉(zhuǎn)速較低時，隨機(jī)擾動的影響明顯，并且對分岔點(diǎn)區(qū)域的影響更為顯著；轉(zhuǎn)子轉(zhuǎn)速較高時，應(yīng)用偽周期替代數(shù)據(jù)法能夠有效的對含噪混沌時間序列進(jìn)行識別。今后將有效結(jié)合時頻分析方法，對不確定性激勵下故障轉(zhuǎn)子系統(tǒng)的振動信號進(jìn)行比較分析、識別及分類，為振動系統(tǒng)狀態(tài)預(yù)測與故障預(yù)警提供新的方法和思路。