劉桂珍,于 影,殷寶麟 ,丁海娟,趙炯介,聞邦椿

(1.佳木斯大學(xué) 機械工程學(xué)院,黑龍江 佳木斯 154007; 2.東北大學(xué) 機械工程與自動化學(xué)院,沈陽 110819)

?

激勵頻率影響下的裂紋轉(zhuǎn)子非線性分析

劉桂珍1,于 影1,殷寶麟1,丁海娟1,趙炯介1,聞邦椿2

(1.佳木斯大學(xué) 機械工程學(xué)院,黑龍江 佳木斯 154007; 2.東北大學(xué) 機械工程與自動化學(xué)院,沈陽 110819)

利用拉格朗日方程建立了非穩(wěn)態(tài)油膜力作用下的轉(zhuǎn)子-定子-軸承系統(tǒng)裂紋故障的力學(xué)模型,應(yīng)用數(shù)值分析方法研究了該系統(tǒng)隨激勵頻率變化時的時域波形、軸心軌跡和Poincare截面,揭示了系統(tǒng)分岔特性和進入混沌的途徑.研究結(jié)果表明,當激勵頻率作為唯一控制參數(shù)時,系統(tǒng)呈現(xiàn)出復(fù)雜的非線性動力學(xué)行為.該結(jié)果為大型旋轉(zhuǎn)機械的故障診斷提供了理論依據(jù).

裂紋; 時域波形; 軸心軌跡; Poincare截面

大型旋轉(zhuǎn)機械在高速旋轉(zhuǎn)過程中,由于長期受到交變載荷作用,轉(zhuǎn)軸易出現(xiàn)橫向疲勞裂紋,導(dǎo)致斷軸事故時有發(fā)生,橫向裂紋是大型旋轉(zhuǎn)機械常見故障之一.橫向裂紋對轉(zhuǎn)子系統(tǒng)的安全運行產(chǎn)生巨大威脅,裂紋使轉(zhuǎn)軸剛度減小并發(fā)生周期性變化,導(dǎo)致轉(zhuǎn)子運行不穩(wěn)定.數(shù)十年來,國內(nèi)外學(xué)者在裂紋轉(zhuǎn)子故障診斷方面作了很多工作,其中裂紋轉(zhuǎn)子的運動特性和裂紋征兆已成為各種裂紋轉(zhuǎn)子診斷技術(shù)的基礎(chǔ),被廣泛應(yīng)用于工程實際中[1-3].為了識別早期裂紋的位置和大小,裂紋轉(zhuǎn)子的振動特性研究引起了廣泛的關(guān)注,眾多研究者提出了各不相同的裂紋模型,但至今也沒有一個公認的力學(xué)模型.

筆者以裂紋故障為研究對象,利用拉格朗日方程建立了非穩(wěn)態(tài)油膜力作用下的轉(zhuǎn)子系統(tǒng)裂紋故障的4質(zhì)量8自由度非線性動力學(xué)模型,應(yīng)用數(shù)值分析研究,將激勵頻率作為唯一控制參數(shù)時系統(tǒng)的動態(tài)響應(yīng).通過轉(zhuǎn)子時域波形圖、軸心軌跡和Poincare截面圖等振動信號提供的故障信息,對該系統(tǒng)響應(yīng)的非線性行為和故障機理進行分析,從而為該類轉(zhuǎn)子系統(tǒng)的故障診斷和系統(tǒng)的安全運行提供理論依據(jù).

1 力學(xué)模型與微分方程

1.1 拉格朗日方程描述

設(shè)有n個質(zhì)點組成的質(zhì)點系,受完整的理想約束,具有N個自由度,其位置可由N個廣義坐標方程來確定,則有

(1)

式中:L為拉格朗日函數(shù),L=T-V,T為系統(tǒng)的動能函數(shù),V為系統(tǒng)的勢能函數(shù);R為與系統(tǒng)阻尼相對應(yīng)的耗散函數(shù);Qi為作用在系統(tǒng)上的廣義力;qi為系統(tǒng)獨立的廣義坐標;N為系統(tǒng)的總自由度個數(shù).

1.2 剛度影響系數(shù)的計算

設(shè)轉(zhuǎn)軸半徑為R、長度為L1的無質(zhì)量彈性圓軸,在轉(zhuǎn)軸中央有一深度為a的弓形橫向裂紋,如圖1所示.

圖1 裂紋截面示意圖

如果只考慮裂紋處彎矩的作用,根據(jù)斷裂力學(xué)理論,由于裂紋的存在將在裂紋局部產(chǎn)生附加角位移.設(shè)在η,ξ方向彎矩作用下的無量綱局部柔度系數(shù)分別為[4]

(2)

對式(2)進行數(shù)值積分,得到無量綱的局部柔度系數(shù)與無量綱裂紋深度的關(guān)系曲線,如圖2所示.ε,δ為僅與裂紋深度a有關(guān)的相關(guān)剛度參數(shù),其表達式如下:

(4)

式中:L1為兩端支承間軸的長度;R為軸頸;ν為泊松比.

圖2 柔度與無量綱裂紋深度關(guān)系

1.3 非穩(wěn)態(tài)油膜力模型[5]

非線性油膜力模型采用短軸承假設(shè)下的Ca-pone非線性油膜力模型,該模型有較好的精度和收斂性.在短軸承油膜力假設(shè)條件下的無量綱雷諾方程為

由式(5)可得無量綱油膜壓力為

(6)

式中:D為軸承直徑.

(7)

1.4 運動微分方程

轉(zhuǎn)子-定子-軸承系統(tǒng)中間軸上具有一橫向裂紋的力學(xué)模型,如圖3所示,兩端由滑動軸承支撐,滑動軸承直徑為D,長度為L.兩軸承之間為一無質(zhì)量彈性軸,其半徑為R,長度為L1,轉(zhuǎn)軸中央有一對稱布置的圓盤,深度為a的弓形橫向裂紋.O1為轉(zhuǎn)子的幾何中心;Oc為轉(zhuǎn)子質(zhì)心;mi(i=1,2,3,4)分別為轉(zhuǎn)子、軸承、定子和轉(zhuǎn)軸在軸承處的半集中質(zhì)量(kg);ki(i=1,2,3)分別為轉(zhuǎn)軸、軸承支撐處和定子基礎(chǔ)的剛度系數(shù)(N·m-1);ci(i=1,2,3)分別為轉(zhuǎn)軸、軸承支撐處和基礎(chǔ)對定子的阻尼系數(shù)(N·sm-1).

圖3 非穩(wěn)態(tài)油膜力的轉(zhuǎn)子-軸承系統(tǒng)裂紋故障的力學(xué)模型

設(shè)轉(zhuǎn)盤在水平方向、垂直方向移動的位移為x1,y1,左端軸承在水平方向、垂直方向移動的位移為x2,y2,定子在水平方向、垂直方向移動的位移為x3,y3,右端軸承座處水平方向、垂直方向移動的位移為x4,y4.

采用耗散能量形式的拉格朗日方程推導(dǎo)出對應(yīng)裂紋故障的轉(zhuǎn)子-定子-軸承系統(tǒng)的質(zhì)點運動微分方程:

(y1-y4)sin2τ]=m1eω2cosωt

(y1-y4)sin2τ]=m1eω2sinωt-m1g

k3y3+m3g=0

(y4-y1)sin2τ]=

(8)

式中:xi,yi(i=1,2,3,4)分別為應(yīng)各質(zhì)量的位移坐標;b為軸頸的偏心量(mm);Fx,Fy分別為水平、鉛垂方向的油膜力(N).

為研究問題的方便性,將上述方程無量綱化.設(shè):

則上述方程轉(zhuǎn)化為

ρcosωt

(y1-y4)sin2τ]+ρsinωt-G

η3y3-G

(9)

2 數(shù)值模擬

運用四階Runge-Kutta法對數(shù)值進行求解.在計算中為了能夠較快地得到穩(wěn)定解,應(yīng)將步長選得盡量小且周期足夠多.為了消除瞬態(tài)響應(yīng)的影響,舍棄前40個周期,計算軌跡圖時取后10～20個周期.選取系統(tǒng)參數(shù)如下:m1=4.0kg,m2=32.1kg,m3=50.0kg,m4=20.0kg,c1=1.05kN·sm-1,c2=2.1kN·sm-1,c3=2.1kN·sm-1,k1=250kN·m-1,k2=250kN·m-1,k3=25 000kN·m-1,kr=1 000kN·m-1,R=0.025m,L=0.570m,δ2=0.2mm;η=0.018MPa,R1=0.015m,f=0.2.通過計算得出轉(zhuǎn)子系統(tǒng)的3級固有頻率分別為f1=13.974 7Hz,f2=43.586 4Hz,f3=113.109 6Hz.利用激勵頻率作為控制參數(shù),分別考察轉(zhuǎn)子系統(tǒng)在臨界轉(zhuǎn)速以下、臨界轉(zhuǎn)速附近以及臨界轉(zhuǎn)速以上3種情況下,裂紋轉(zhuǎn)子系統(tǒng)非線性動力學(xué)響應(yīng).如圖4～11分別表示對應(yīng)轉(zhuǎn)子系統(tǒng)三階臨界轉(zhuǎn)速時的轉(zhuǎn)子系統(tǒng)動態(tài)響應(yīng).

圖4 轉(zhuǎn)動角速度ω=80 rad/s時的轉(zhuǎn)、定子時域圖、相圖、Poincare截面和頻譜圖

圖4～6表示轉(zhuǎn)子系統(tǒng)位于一階臨界轉(zhuǎn)速附近的3組裂紋故障轉(zhuǎn)子系統(tǒng)的響應(yīng).觀察此時 Poincare截面圖、相圖、時間歷程、幅值譜圖等,表明一個穩(wěn)態(tài)周期響應(yīng)隨系統(tǒng)參數(shù)變化的過程.在很短的參數(shù)區(qū)間內(nèi),系統(tǒng)穩(wěn)態(tài)響應(yīng)快速地完成了從擬周期運動→周期運動→倍周期運動的整個過程.從對應(yīng)的Poincare截面圖可以看到從周期吸引子(有限個不動點)→孤點→P2周期點的演化過程;從對應(yīng)的幅值譜圖可以看到,系統(tǒng)響應(yīng)的次諧波響應(yīng)能不斷分岔產(chǎn)生新的諧波分量,這些諧波分量持續(xù)分岔又產(chǎn)生了大量的低頻頻率成分,因而使定子的時域波形圖表現(xiàn)出比轉(zhuǎn)子時域波形圖變化劇烈.

圖5 ω=90 rad/s時的轉(zhuǎn)、定子時域圖、相圖、Poincare截面和頻譜圖

圖6 ω=120 rad/s時的轉(zhuǎn)、定子時域圖、相圖、Poincare截面和頻譜圖

圖7 ω=250 rad/s時的轉(zhuǎn)、定子時域圖、相圖、Poincare截面和頻譜圖

圖8 ω=273 rad/s時的轉(zhuǎn)、定子時域圖、相圖、Poincare截面和頻譜圖

圖7～9表示位于轉(zhuǎn)子系統(tǒng)二階臨界轉(zhuǎn)速附近的3組裂紋故障轉(zhuǎn)子的響應(yīng).Poincare截面圖、相圖、時域圖和幅值譜圖等多種圖示,表明系統(tǒng)在很短的參數(shù)區(qū)間內(nèi),快速地完成了從周期運動→混沌運動→倍周期運動的整個過程.從對應(yīng)的Poincare截面圖可以看到從一個吸引子→擬周期→二周期點的快速演化過程;從對應(yīng)的幅值譜圖可以看到,系統(tǒng)響應(yīng)的次諧波響應(yīng)能不斷分岔產(chǎn)生新的諧波分量,這些諧波分量持續(xù)分岔又產(chǎn)生了大量的低頻頻率成分,因而使定子上的時間歷程呈現(xiàn)出劇烈變化的特性.

圖10～11是轉(zhuǎn)子系統(tǒng)位于三階臨界轉(zhuǎn)速附近時裂紋故障轉(zhuǎn)子的響應(yīng).按50 rad/s間隔取點過程發(fā)現(xiàn),當ω位于(500,700) rad/s和ω位于(750,1 000) rad/s區(qū)間時,系統(tǒng)呈發(fā)散狀態(tài),因而將450 rad/s作為三階臨界轉(zhuǎn)速的臨界值進行討論.

圖9 ω=300 rad/s時的轉(zhuǎn)、定子時域圖、相圖、Poincare截面和頻譜圖

圖10 轉(zhuǎn)動角速度ω=450 rad/s時的轉(zhuǎn)、定子時域圖、相圖、Poincare截面和頻譜圖

當ω=450 rad/s時,Poincare截面圖、相圖、時間歷程組和幅值譜圖等多種圖示,表明系統(tǒng)處于擬周期運動,隨后系統(tǒng)在(500,700)rad/s區(qū)間內(nèi)呈現(xiàn)發(fā)散狀態(tài),當ω=710 rad/s時,系統(tǒng)快速恢復(fù)周期運動,在ω在(750,1 000)rad/s區(qū)間時,系統(tǒng)重新呈現(xiàn)發(fā)散狀態(tài).

3 結(jié)論

通過對帶有裂紋故障的轉(zhuǎn)子-定子-軸承系統(tǒng)在不同轉(zhuǎn)速比下的運動特性的研究,得出以下結(jié)論:

圖11 轉(zhuǎn)動角速度ω=710 rad/s時的轉(zhuǎn)、定子時域圖、相圖、Poincare截面和頻譜圖

(1) 當轉(zhuǎn)子出現(xiàn)裂紋時,在臨界轉(zhuǎn)速為kωc/n(k=1,2,…,;n=1,2,…,h)附近振動增大,定子在系統(tǒng)響應(yīng)中出現(xiàn)高次諧波分量,且當轉(zhuǎn)軸轉(zhuǎn)速為臨界轉(zhuǎn)速的倍數(shù)時,(n/k)×分量達到最大值,即存在分數(shù)次頻共振現(xiàn)象.

(2) 反映在轉(zhuǎn)子、定子的時域波形上看,低速時(ω≤300 rad/s)時定子振動幅頻值變化更為劇烈;高速時(ω≥450 rad/s),轉(zhuǎn)子、定子時域圖形出現(xiàn)同步振動.

(3) 因定子的幅頻特性曲線具有豐富的分數(shù)次頻和倍頻成分,因而定子的相圖比轉(zhuǎn)子的相圖更能體現(xiàn)出組合頻率的特征.

[1] 丁千,陳予恕.非線性轉(zhuǎn)子-密封系統(tǒng)的亞諧共振失穩(wěn)機理研究[J].振動工程學(xué)報,1997(4):404-412.

DING Qian,CHEN Yushu.Study on the mechanism of subharmonic instability of nonlinear rotor/seal system[J].Journal of Vibration Engineering,1997(4):404-412.

[2] 曹樹謙,丁千,陳予恕.具有滑動軸承的穩(wěn)態(tài)轉(zhuǎn)子系統(tǒng)有限元建模分析[J].汽輪機技術(shù),1999,41(6):347-350.

CAO Shuqian,DING Qian,CHEN Yushu.Analysis on modeling steady rotor system with sliding bearings by using FEM[J].Turbine Technology,1999,41(6):347-350.

[3] 鄭惠萍,陳予恕.滑動軸承轉(zhuǎn)子系統(tǒng)抗擾動穩(wěn)定裕度的研究[J].汽輪機技術(shù),2001,43(1):35-37.

ZHENG Huiping,CHEN Yushu.Research on stability margin of rotor system with sliding bearing[J].Turbine Technology,2001,43(1):35-37.

[4] GASCH R.A survey of the dynamic behavior of a simple rotating shaft with a transverse crack[J].Journal of Sound and Vibration,1993,160(2):313-332.

[5] 劉桂珍,于影,李憲芝,等.轉(zhuǎn)子系統(tǒng)非穩(wěn)態(tài)油膜力的非線性分析[J].機床與液壓,2013,41(11):31-34.

LIU Guizhen,YU Ying,LI Xianzhi,et al.The nonlinear analysis of unsteady oil film force in rotor system[J].Machine Tool & Hydraulics,2013,41(11):31-34.

Nonlinear analysis oncrack rotor via excitation frequency

LIU Guizhen1,YU Ying1,YIN Baolin1,DING Haijuan1,ZHAO Jiongjie1,WEN Bangchun2

(1.College of Mechanical Engineering,Jiamusi University,Jiamusi 154007,Heilongjiang,China;2.School of Mechanical Engineering & Automation,Northeastern University,Shenyang 110819,China)

Established the unsteady oil-film force with crack fault of rotor-stator-bearing system mechanical models by using Lagrange equation.The numerical analysis method is used to study system with excitation frequency of the time-domain waveform and orbit trace and Poincare section,revealing the system into chaos and bifurcation characteristic ways.The results show that when the excitation frequency as the only control parameters,system showing a complex nonlinear dynamics behavior,the results for the follow-up study of large rotating machinery fault diagnosis provides a theoretical basis.

crack; time-domain waveform; axis orbit; Poincare section

黑龍江省自然科學(xué)基金資助項目(E2016066);佳木斯大學(xué)精密制造技術(shù)研究團隊資助項目(cxtdpy-2016-04)

劉桂珍(1965—),女，副教授,博士.E-mail:jmslgz@163.com

TH 133

A

1672-5581(2017)02-0170-08