周生通， 李鴻光， 張 龍， 周新建

(1.華東交通大學(xué) 機(jī)電與車輛工程學(xué)院，南昌 330013；2.上海交通大學(xué) 機(jī)械系統(tǒng)與振動國家重點(diǎn)實驗室，上海 200240)

?

基于嵌入式譜隨機(jī)有限元法的轉(zhuǎn)子系統(tǒng)隨機(jī)不平衡響應(yīng)特性分析

周生通1， 李鴻光2， 張 龍1， 周新建1

(1.華東交通大學(xué) 機(jī)電與車輛工程學(xué)院，南昌 330013；2.上海交通大學(xué) 機(jī)械系統(tǒng)與振動國家重點(diǎn)實驗室，上海 200240)

原始不平衡量在轉(zhuǎn)軸中沿軸線呈任意空間分布，因而實際轉(zhuǎn)子的不平衡量需用分布不平衡量和點(diǎn)不平衡量兩種方式共同表達(dá)。考慮兩種不平衡量的隨機(jī)性，以隨機(jī)場和隨機(jī)變量的方式模擬轉(zhuǎn)子中存在的兩類不平衡量，并采用嵌入式譜隨機(jī)有限元方法建立隨機(jī)不平衡量下的轉(zhuǎn)子系統(tǒng)不平衡響應(yīng)隨機(jī)分析模型。結(jié)果顯示，提供的隨機(jī)表征方法能較好地模擬兩類轉(zhuǎn)子不平衡量的隨機(jī)性；不平衡量的隨機(jī)性影響著轉(zhuǎn)子不平衡響應(yīng)的變異性，可以發(fā)現(xiàn)：① 共振頻率附近的位移響應(yīng)幅值標(biāo)準(zhǔn)差也較大；② 位移響應(yīng)幅值的變異系數(shù)受低頻段內(nèi)共振頻率的影響較小，但隨著共振頻率階次的增高影響也逐漸增大；③ 位移響應(yīng)相位角的標(biāo)準(zhǔn)差變化情況與位移幅值的變異系數(shù)變化情況是基本類似的。

質(zhì)量不平衡；譜隨機(jī)有限元方法；隨機(jī)不平衡響應(yīng)

一個轉(zhuǎn)子在設(shè)計上一般都使它相對于旋轉(zhuǎn)軸線是軸對稱的。但是由于工藝上的一系列因素，最后裝配完畢的轉(zhuǎn)子總是不能做到動力上的完全軸對稱，也就是存在一定的不平衡量，這種不平衡量通常稱之為原始不平衡量[1]。造成原始不平衡量的因素主要有：轉(zhuǎn)子材質(zhì)的不均勻性，鍵槽不對稱引起的不平衡，轉(zhuǎn)子加工中總是產(chǎn)生一些圓度偏差和偏心等。對于汽輪機(jī)轉(zhuǎn)子尚有各個葉片之間的差別，葉片鎖口及末葉片的不對稱等不平衡因素。對于壓縮機(jī)尚有葉輪的不平衡量。所有這些因素造成的不平衡量都屬于隨機(jī)性質(zhì)的[1]。

像葉片和葉輪位置處的不平衡量是轉(zhuǎn)子系統(tǒng)不平衡量的主要貢獻(xiàn)源。在理論分析時，一般將葉片和葉輪簡化為剛性圓盤，相應(yīng)的不平衡量則以力的形式施加在圓盤的幾何中心處，用不平衡質(zhì)量、偏心距和相位角三個量確定不平衡的大小和方向。不過，如果將連續(xù)轉(zhuǎn)軸分割為微小厚度的圓片，那么正如單個圓盤一樣，每一個圓片上都存在一定的不平衡量，且大小和方向各不相同。也就是說，轉(zhuǎn)軸中的不平衡量分布函數(shù)實為一條任意的、隨機(jī)的空間分布曲線[1-2]。雖然不平衡量的空間任意分布事實很早就被人們認(rèn)識到，但直到1993年才被LEE等[3]在傳遞矩陣方法中予以考慮。后來，SHIH等[4]又提出了識別柔性轉(zhuǎn)子不平衡量分布情況的方法。該方法將分段連續(xù)的質(zhì)量偏心曲線用傅里葉級數(shù)展開，但這種方法不太適用于轉(zhuǎn)子動力學(xué)的有限元分析程序，為此YANG等[5]提出采用多項式曲線模擬轉(zhuǎn)子不平衡質(zhì)量偏心曲線，以便借助有限元程序更容易地識別轉(zhuǎn)子中存在的不平衡量分布。DEEPTHIKUMAR等[6]在此基礎(chǔ)上進(jìn)一步發(fā)展了同時存在分布不平衡和軸彎曲故障的柔性轉(zhuǎn)子模態(tài)平衡方法。不過，上述這些研究都是從確定性的角度處理轉(zhuǎn)子不平衡量的。

現(xiàn)有對轉(zhuǎn)子不平衡量隨機(jī)性的考慮多是針對圓盤不平衡量的[7-8]，即將圓盤的偏心距、質(zhì)量或相位角作為隨機(jī)變量處理。不過，STOCKI等[9]在轉(zhuǎn)子振動離散性分析中則考慮了分布不平衡量的不確定性，其采用的思路是將分布不平衡量表示為轉(zhuǎn)子前四階振型的加權(quán)疊加，權(quán)值設(shè)為隨機(jī)變量，并使用非嵌入式的譜方法開展了相應(yīng)的隨機(jī)分析工作。本文同樣從隨機(jī)角度考慮轉(zhuǎn)子系統(tǒng)中的不平衡量，并以隨機(jī)場和隨機(jī)變量的形式將它們表達(dá)為兩種形式，即集中不平衡量(Lumped Unbalance)和分布不平衡量(Distributed Unbalance)。最后，對某一50 MW汽輪機(jī)轉(zhuǎn)子同時計及這兩種不平衡量，并基于嵌入式譜隨機(jī)有限元方法建立該轉(zhuǎn)子的隨機(jī)不確定性分析模型，實現(xiàn)其不平衡響應(yīng)的隨機(jī)分析。

1 轉(zhuǎn)子不平衡量的隨機(jī)表征

1.1 分布不平衡量

圖1 轉(zhuǎn)子中的不平衡量Fig.1 Rotor unbalances

將分布不平衡量在單個轉(zhuǎn)子單元內(nèi)產(chǎn)生的動能代入拉格朗日方程，即可得到對應(yīng)的不平衡力單元向量：

(1)

(2)

為了處理積分中的隨機(jī)量，將一維隨機(jī)場u(x,ω)和φ(x,ω)分別作Karhunen-Loeve展開并分別保留至L1階和L2階，有：

(3)

(4)

(5)

式(5)由于包含了連乘的指數(shù)函數(shù)，進(jìn)一步處理會比較繁瑣。這里將相位角隨機(jī)場φ(x,ω)退化為正態(tài)隨機(jī)變量，即?。?/p>

(6)

這樣，式(5)就可簡化為：

(7)

進(jìn)一步，取指數(shù)函數(shù)的級數(shù)展開近似表達(dá)：

(8)

則由式(7)表示的轉(zhuǎn)子空間隨機(jī)不平衡量產(chǎn)生的單元力向量就可整理為如下形式：

(9)

式中:Qkj是維數(shù)為8×1的列矩陣，表達(dá)式為：

(10)

更特別地，若不平衡量分布密度u(x,ω)也同樣退化為正態(tài)隨機(jī)變量，即：

(11)

則有：

(12)

式中:

(13)

1.2 集中不平衡量

對于集中不平衡量，只需用隨機(jī)變量就可完全表征不平衡量大小和相位角的隨機(jī)性。在文獻(xiàn)[8]中已經(jīng)給出了當(dāng)兩者都是正態(tài)隨機(jī)變量但相位角均值為零時的不平衡量表征公式。這里則直接將其擴(kuò)展為不平衡量大小和相位分別是具有任意分布參數(shù)的正態(tài)隨機(jī)變量的情況，對應(yīng)的集中不平衡量表征公式為：

(14)

式中:Qkj是維數(shù)為4×1的列矩陣，表達(dá)式為：

(15)

2 不平衡響應(yīng)的嵌入式譜隨機(jī)有限元方程

譜隨機(jī)有限元方程是依據(jù)多項式混沌基的正交性質(zhì)建立的。從轉(zhuǎn)子隨機(jī)不平衡量的表征公式(9)、(12)和(14)中可以看到，作為右端項的不平衡力Fu是以多項式的形式表達(dá)的，無法直接利用多項式混沌基的正交性。需要將Fu的表達(dá)式改寫成以多項式混沌展開的形式，即：

(16)

將原來Fu的表達(dá)式改寫為式(16)的形式后，即可建立起轉(zhuǎn)子系統(tǒng)不平衡響應(yīng)的譜隨機(jī)有限元方程：

(17)

式中:Z(ω)是轉(zhuǎn)子系統(tǒng)的動剛度矩陣。在不考慮其他隨機(jī)因素的情況下，式(17)最終可轉(zhuǎn)化為如下P個線性方程組：

Z(ω)δi=Fi,i=0,1,…,P-1

(18)

解方程組即可求得不平衡響應(yīng)的多項式混沌展開中的廣義坐標(biāo)值δi,i=0,…,P-1，從而表達(dá)出隨機(jī)不平衡響應(yīng)Δ，即為：

(19)

(20)

3 數(shù)值算例

求解某一50 MW汽輪機(jī)轉(zhuǎn)子動平衡系統(tǒng)(如圖2所示)的隨機(jī)不平衡響應(yīng)。設(shè)汽輪機(jī)轉(zhuǎn)軸中的不平衡量分布密度的大小u(x,ω)服從高斯隨機(jī)場且具有指數(shù)協(xié)方差函數(shù)，相應(yīng)的相位角φ服從正態(tài)隨機(jī)量；同時，取作用在圓盤10和圓盤19上的不平衡量U10和U19的大小和相位角亦均服從正態(tài)分布。具體的隨機(jī)量設(shè)定值如表1所示。

圖2 某50 MW汽輪機(jī)轉(zhuǎn)子動平衡系統(tǒng)Fig.2 A 50 MW steam turbine rotor dynamic balancing system

元件屬性名稱分布類型均值標(biāo)準(zhǔn)差備注圓盤10不平衡量U10正態(tài)分布0.190.038相位角φ10正態(tài)分布0π/15圓盤19不平衡量U19正態(tài)分布0.30.06相位角φ19正態(tài)分布π/3π/15轉(zhuǎn)軸相位角φ正態(tài)分布π/6π/60分布不平衡量大小u(x)一維高斯隨機(jī)場0.010.001設(shè)轉(zhuǎn)軸分布不平衡量隨機(jī)場具有指數(shù)協(xié)方差函數(shù),相關(guān)長度6.36m,區(qū)間為[0,6.36]

該模型中，單個集中不平衡量采用2維5階精度描述，分布不平衡量采用5維5階精度描述，故整個譜隨機(jī)有限元模型將形成一個具有9維5階的多項式混沌展開模型，共計P=2 002項多項式混沌基函數(shù)。進(jìn)而，依據(jù)式(18)～(20)即可得到轉(zhuǎn)子上任意節(jié)點(diǎn)位移的隨機(jī)響應(yīng)情況。限于篇幅，這里僅給出圓盤10位置處的位移響應(yīng)概率統(tǒng)計信息，如圖3～7所示。圖3(a)和(b)所示分別為圓盤10處的水平和豎直位移頻響函數(shù)幅值的均值、上下包絡(luò)線和標(biāo)準(zhǔn)差曲線；圖4給出了對應(yīng)水平和豎直位移頻響幅值的變異系數(shù)隨激勵頻率的變化情況；圖5和圖6則類似地給出了水平和豎直位移相位角的統(tǒng)計信息；圖7(a)和(b)給出了在特定轉(zhuǎn)速頻率(25.5 Hz和40.0 Hz)下水平位移的直方圖和基于核密度估計的概率密度信息，其中頻率25.5 Hz靠近轉(zhuǎn)子系統(tǒng)一階臨界轉(zhuǎn)速，頻率40.0 Hz則遠(yuǎn)離共振頻率。

其中，位移均值(實線)、標(biāo)準(zhǔn)差(實線)、上包絡(luò)(點(diǎn)線)、下包絡(luò)(虛線)和1 000個Monte Carlo樣本(灰線)圖3 圓盤10的位移幅值均值和標(biāo)準(zhǔn)差信息Fig.3 Mean and standard deviation of displacementamplitudes at disk 10

圖4 圓盤10處水平和豎直位移幅值的變異系數(shù)曲線Fig.4 The coefficient of variation curves of horizontal and vertical displacement amplitudes at disk 10

從圖3(a)和(b)中可以看到，受到轉(zhuǎn)子隨機(jī)不平衡量的影響，圓盤10處的水平和豎直位移頻響幅值曲線表現(xiàn)出了離散特性，并圍繞著均值曲線上下隨機(jī)波動，而偏離均值的大小與激勵頻率有關(guān)，當(dāng)激勵頻率與系統(tǒng)固有頻率一致，即發(fā)生共振時，偏離量最大，即在共振處出現(xiàn)有較大的標(biāo)準(zhǔn)差峰值。但此處的標(biāo)準(zhǔn)差大小并不能準(zhǔn)確反映位移頻響幅值隨機(jī)離散的程度，圖4則采用變異系數(shù)表達(dá)離散程度以消除均值的影響。從圖4中可以看到，水平和豎直位移頻響函數(shù)幅值的變異系數(shù)在小于40 Hz的頻段內(nèi)基本沒有變化，而這一頻段內(nèi)就包含有轉(zhuǎn)子系統(tǒng)的第一階臨界轉(zhuǎn)速頻率，且此臨界頻率處的響應(yīng)幅值標(biāo)準(zhǔn)差就比較大。隨著激勵頻率的增加，當(dāng)達(dá)到其他高階臨界轉(zhuǎn)速時，變異系數(shù)就開始有較大的浮動，可以看到階數(shù)越高變異系數(shù)峰值也就越大。這說明轉(zhuǎn)子的不平衡隨機(jī)因素對系統(tǒng)高階共振頻率處位移幅值的離散程度影響較大，而對低階影響則較小。如在第一階臨界轉(zhuǎn)速頻率處的影響幾乎可以忽略，因為它與附近非共振頻率處的離散程度幾乎一致。

其中，位移相位角均值(實線)、上包絡(luò)(點(diǎn)線)、下包絡(luò)(虛線)和1000個Monte Carlo樣本(灰線)圖5 圓盤10的位移相位角信息Fig.5 Phase angles of displacement at disk 10

圖5(a)和(b)則分別給出了隨機(jī)不平衡量影響下的水平和豎直位移相位角的離散情況，圖6(a)和(b)中顯示的是對應(yīng)相位角的標(biāo)準(zhǔn)差和變異系數(shù)情況。可以看到，與位移幅值離散特性表達(dá)不同，位移相位角的離散程度用標(biāo)準(zhǔn)差表征更好，且此時的水平和豎直位移相位角標(biāo)準(zhǔn)差的變化情況與對應(yīng)的水平和豎直位移幅值變異系數(shù)的變化情況類似，如圖4所示。其中的原因主要是由于位移幅值均值及其標(biāo)準(zhǔn)差的數(shù)值大小在不同的頻率范圍內(nèi)存在有量級上的較大差別，當(dāng)表示離散特性時不同量級上的標(biāo)準(zhǔn)差之間不易反映整體離散程度的大小，此時就需要剔除均值的影響用變異系數(shù)來表示離散程度，但對于位移相位角來說，其標(biāo)準(zhǔn)差不存在量級差別問題，而均值在一定的范圍內(nèi)(不超過360°)變化，此時用標(biāo)準(zhǔn)差表示離散程度會更好，若采用變異系數(shù)則效果上不佳(如圖6(a)和(b)虛線所示)。

圖6 圓盤10的位移相位角標(biāo)準(zhǔn)差和變異系數(shù)頻譜曲線Fig.6 The standard deviation and coefficient of variation curves of displacement phase angle at disk 10

圖7 圓盤10的水平位移在不同轉(zhuǎn)速頻率下的均值、標(biāo)準(zhǔn)差、概率密度函數(shù)直方圖和核密度估計曲線Fig.7 Means, standard deviations, histograms of probability density functions, kernel density estimation curves of the horizontal displacement at disk 10 under different rotating frequencies

進(jìn)一步，從圖7(a)和(b)中可以看到某一特定激勵頻率下水平位移頻響幅值的概率分布比較接近正態(tài)分布，這與我們在例子中所采用的隨機(jī)量類型均為正態(tài)型的有較大關(guān)系。其中，在第一階臨界轉(zhuǎn)速頻率附近(25.5 Hz處)時，其對應(yīng)位移幅值的均值為1.214×10-3m、標(biāo)準(zhǔn)差為2.006×10-4m；在非共振區(qū)域(40 Hz處)時，均值為4.728×10-5m、標(biāo)準(zhǔn)差為7.960×10-6m。從量級上對比，可以看到兩個頻率處的標(biāo)準(zhǔn)差在數(shù)值上相差有兩個數(shù)量級，但這也只能表明系統(tǒng)在共振頻率處的響應(yīng)離散幅值更大；而從幅值離散程度上對比，可以計算出兩頻率處的幅值變異系數(shù)分別為0.165 2和0.168 4，比較接近，說明實際中這兩個頻率處的隨機(jī)離散程度差別并不大，與圖4結(jié)果一致。

4 結(jié) 論

本文采用分布不平衡量和集中不平衡量兩種概念描述轉(zhuǎn)子中存在的不平衡量。在考慮不平衡量隨機(jī)性的基礎(chǔ)上，推導(dǎo)了分布不平衡量和集中不平衡量的力表達(dá)式，并利用嵌入式譜隨機(jī)有限元方法建立了轉(zhuǎn)子系統(tǒng)不平衡響應(yīng)的隨機(jī)分析模型，求解了某一50 MW汽輪機(jī)轉(zhuǎn)子在動平衡擺架支撐系統(tǒng)上的隨機(jī)不平衡響應(yīng)概率信息。總結(jié)算例結(jié)果，可以看到：

(1) 受到轉(zhuǎn)子隨機(jī)不平衡量的影響，系統(tǒng)不平衡響應(yīng)(位移幅值頻響函數(shù))表現(xiàn)出了較大的隨機(jī)離散特性，尤其是在共振頻率附近，其對應(yīng)的響應(yīng)幅值標(biāo)準(zhǔn)差要明顯高于非共振頻率區(qū)；

(2) 但從位移頻響幅值的變異系數(shù)上分析，在低頻段內(nèi)各頻率上的變異系數(shù)基本保持不變(即使是在第一階臨界轉(zhuǎn)速頻率處)，而隨著激振頻率的增大，變異系數(shù)的浮動也開始增大，尤其是在高階次的共振頻率處，變異系數(shù)出現(xiàn)了較大峰值，這說明了轉(zhuǎn)子不平衡隨機(jī)因素對更高階共振頻率處的位移幅值離散程度影響更為明顯；

(3) 對于位移相位角來說，其也表現(xiàn)出了較大的隨機(jī)離散特性，但其離散程度更適合用相位角標(biāo)準(zhǔn)差來表征。

[1] 周仁睦. 轉(zhuǎn)子動平衡——原理、方法和標(biāo)準(zhǔn)[M]. 北京: 化學(xué)工業(yè)出版社, 1992.

[2] 王悅武，田杜平，徐錫林. 高速平衡技術(shù)裝備[M]. 上海: 上海交通大學(xué)出版社, 2013.[3] LEE A, SHIH Y, KANG Y. The analysis of linear rotor-bearing systems: a general transfer matrix method[J]. Journal of Vibration and Acoustics,1993, 115(4): 490-497.

[4] SHIH Y, LEE A. Identification of the unbalance distribution in flexible rotors[J]. International Journal of Mechanical Sciences,1997, 39(7): 841-857.

[5] YANG T, LIN C. Estimation of distributed unbalance of rotors[J]. Journal of Engineering for Gas Turbines and Power,2002,124(4): 976-983.

[6] DEEPTHIKUMAR M B, SEKHAR A S, SRIKANTHAN M R. Modal balancing of flexible rotors with bow and distributed unbalance[J]. Journal of Sound and Vibration,2013,332(24):6216-6233.

[7] ZHANG Y, WEN B, LEUNG A Y T. Reliability analysis for rotor rubbing[J]. Journal of Vibration and Acoustics,2001, 124(1): 58-62.

[8] DIDIER J, FAVERJON B, SINOU J. Analysing the dynamic response of a rotor system under uncertain parameters by polynomial chaos expansion[J]. Journal of Vibration and Control,2012, 18(5): 712-732.

[9] STOCKI R, LASOTA R, TAUZOWSKI P, et al. Scatter assessment of rotating system vibrations due to uncertain residual unbalances and bearing properties[J]. Computer Assisted Methods in Engineering and Science,2012, 19(2): 95-120.

Stochastic unbalance response characteristics of rotor systems based on intrusive spectral stochastic finite element method

ZHOU Shengtong1, LI Hongguang2, ZHANG Long1, ZHOU XinJian1

(1. School of Mechatronics & Vehicle Engineering, East China Jiaotong University, Nanchang 330013, China;2. State Key Lab of Mechanical Systems & Vibration, Shanghai Jiao Tong University, Shanghai 200240, China)

Mass unbalance in a flexible shaft is distributed arbitrarily along its axis, so distributed mass unbalance and lumped mass unbalance should be considered and presented, respectively for mass unbalance simulation of real rotor systems. Taking the random nature of these two unbalances in a flexible shaft into account, they are represented as random fields and random variables. Here, stochastic analysis of unbalance response of rotor systems with random mass unbalances were performed using the intrusive spectral stochastic finite element method. Numerical results illustrated that the presented stochastic representation methods of distributed and lumped mass unblances are suitable to simulate the random nature of mass unbalances. Results showed that the variations of unblance responses are largely affected by the randomness of unbalances, for example, ① the standard devitations of displacement amplitudes near the resonance frequencies become larger; ② the effects of resonance frequences with in a lower frequency band on the variation coefficient of displacement amplitudes are less, but the effects increase with increase in resonance frequency order; ③ the change curve of the standard deviation of displacement phase angle is similar to that of the variation coefficient of displacement amplitudes.

mass unbalance; spectral stochastic finite element method; stochastic unbalance response

國家自然科學(xué)基金(51505146)；江西省自然科學(xué)基金(20161BAB216135;20122BAB206027)

2015-07-02 修改稿收到日期：2015-09-27

周生通 男，博士，講師，1984年生

李鴻光 男，博導(dǎo)，教授，1972年生

E-mail: hgli@sjtu.edu.cn;zhoust@ecjtu.edu.cn

TH113

A

10.13465/j.cnki.jvs.2016.19.008