基于變厚度殼單元的旋轉葉片固有特性分析

馬 輝，能海強，宋溶澤，張 志，聞邦椿

（東北大學機械工程與自動化學院，沈陽110819）

基于變厚度殼單元的旋轉葉片固有特性分析

馬 輝，能海強，宋溶澤，張 志，聞邦椿

（東北大學機械工程與自動化學院，沈陽110819）

針對1個扭型、變截面旋轉葉片，基于A N SY S有限元軟件，采用變厚度殼單元來模擬真實葉片的方法建模，并與實體單元模型對靜頻和固有振型進行對比，驗證了該建模方法的準確性；在此基礎上分析了科氏力、旋轉軟化、離心剛化以及三者耦合對葉片動頻特性的影響規律；基于殼單元動頻數據，對現在常用的動頻系數經驗公式的適用性進行了簡單評價。研究結果表明：2種模型的動頻規律一致，對系統影響由大到小依次為離心剛化、旋轉軟化和科氏力；葉片A 0振型的動頻系數吻合較好，A 1振型的反之。

旋轉葉片；固有特性；殼單元；A N SY S；動頻

0 引言

旋轉葉片的振動和疲勞損傷故障是航空發動機中較為嚴重的問題。一般說來，研究透平機械葉片的振動問題首先需研究其固有特性，從而據此對葉片進行振動測試和故障診斷，旋轉葉片的固有特性主要指動頻特性和相應的模態振型[1]。葉片高速旋轉時受到離心應力和科氏力的作用，產生旋轉軟化和應力剛化效應，因而在葉片固有特性分析時必須考慮這些因素。國內外許多學者針對葉片固有特性的分析作了大量研究工作。如采用梁單元來模擬葉片進行固有頻率和振動響應分析[2-5]。但由于梁模型忽略了葉片弦向的彎曲變形，順翼展方向的變形并且扭轉也涉及得很少，對于典型的非對稱葉片不能精確計算葉片局部應力和變形，而對于短葉片的分析精度也達不到要求[6]。因此一些學者提出采用板、殼結構來模擬葉片[7-11]。但由于真實葉片一般為空間變截面且帶有預扭角的復雜結構，將葉片考慮為梁或板殼會帶來較大誤差，于是很多學者采用真實葉片的實體有限元模型來計算葉片的動力學特性[12-14]。由于梁忽略了弦向變形，故不適合模擬短葉片，采用實體單元模擬真實葉片需要耗費較大的計算機資源；而采用變厚度的殼單元，不僅可以節省計算機資源，還可以分析帶扭型的變截面葉片。

本文以1個帶扭型的變截面旋轉葉片為研究對象，基于ANSYS有限元軟件，采用變厚度的殼單元建立葉片有限元模型，并與實體單元模型進行對比，通過靜頻、固有振型和動頻的對比結果，驗證建模方法的正確性，并分析不同影響因素對葉片動頻的影響。

1 旋轉葉片振動的有限元動力學方程

旋轉葉片振動的動力學方程為

令C1=C+G，K1=K+S-Kspin，則旋轉葉片的自由振動的動力學方程可化為

當葉片轉速不為0時，由式（2）可知，由于旋轉葉片自由振動方程中存在不對稱的科氏力矩陣，導致方程不能通過經典的振型分解法解耦，利用復模態理論求解系統特征值及對應的特征向量。

令u=Φξ，Φ為與M矩陣正交的特征向量矩陣，ξ為模態坐標向量。將其代入式（2）可得

式中：Λ2為包含前i階特征值ωi的對角陣；ΦTKunsymΦ為不對稱模態剛度矩陣。

通過引入2i維狀態向量的方法，式（3）可寫為

通過QR算法和逆迭代法可計算復模態子空間特征值，而原系統的復特征向量為

2 基于變厚度殼單元葉片有限元模型

由于真實葉片為變截面扭形葉片，基于ANSYS有限元軟件，通過抽取葉片中面數據，根據葉片型線數據，將葉片的厚度δ擬合為隨空間位置變化的變量，進而將葉片簡化為變厚度的殼單元。采用Shell181殼單元來模擬變厚度扭形葉片，此單元具有4個結點，分別有6個自由度，即沿χ、y、z方向的移動和轉動自由度，適合分析薄及中等厚度的板殼結構零件，并支持線性、大扭轉和大應變、變厚度、非線性等分析。將δ擬合為

式中系數見表1。

表1 葉片厚度表達式系數

3 實體單元和殼單元的靜動頻特性對比分析

3.1 靜頻特性對比

為驗證變厚度殼模型的準確性，本文通過實體單元和殼單元進行靜、動頻特性對比。采用9個截面型線數據來繪制變截面扭形實體葉片，如圖1（a）所示，根據實際葉片數據抽取葉片沿厚度方向的中截面數據，選取7個中截面型線數據，如圖1（b）所示。對于實體葉片的建模采用20節點3維實體單元Soild186，殼葉片的建模采用4節點殼單元Shell181，劃分的有限元網格如圖2所示。約束形式采用葉身根部全固支。本文所選用的葉片數據如下：葉片彈性模量 E=1.25×1011Pa，密度ρ=4370 kg/m3，泊松比υ=0.3，圓盤半徑R=216.52 mm， 葉片工作轉速Ω=15000 r/min。

圖1 葉片幾何模型

不考慮轉速的影響，分別計算實體單元和殼單元模型的固有頻率，見表2。從表2中可見，采用變厚度殼單元和實體單元計算的葉片靜頻具有很好的一致性，前6階固有頻率的誤差均小于1%，且振型完全一致。限于篇幅本文只給出了前3階的固有振型，如圖3、4所示。通過對比可知2種建模所描述的振型完全一致，這也證明了在計算葉片靜頻時，采用變厚度殼單元完全可行。下文考慮轉速的影響，再對比二者的動頻特性。

表2 實體單元和殼單元靜頻對比

圖3 葉片靜止工況下振型(Solid186)

圖4 葉片靜止工況下振型(Solid181)

3.2 動頻特性對比

對模型施加轉速ω=0～20000 r/min(每間隔1000 r/min為1個計算點)，考慮離心剛化效應、旋轉軟化效應、科氏力影響以及三者耦合影響下，分析實體單元和殼單元的前6階動頻特性（f1～f6），如圖5所示。圖中：A表示僅考慮離心剛化的影響，B表示考慮離心剛化、旋轉軟化和科氏力耦合影響，C表示僅考慮科氏力影響，D表示僅考慮旋轉軟化影響；下標1、2分別表示實體單元和殼單元。從圖中可見，對于實體單元和殼單元而言，第1～4和6階所顯示的動頻特性趨勢完全相同，科氏力對動頻幾乎沒有影響；旋轉軟化均使動頻減小，且主要影響第1階動頻；離心剛化對系統動頻影響最大，且主要影響彎曲固有頻率；考慮耦合影響的葉片動頻略小于離心剛化影響下的動頻。第5階實體單元的規律基本與其余5階的規律一致，但殼單元卻出現了不一致：殼單元模型在離心剛化作用下，出現了頻率降低現象。雖然隨轉速變化動頻變化規律不同，但二者實際差值并不大，之所以出現動頻規律差異，主要因為此階振型為彎扭耦合振型，而在計算精度上殼單元要略差于實體單元所致。

圖5 考慮離心剛化、旋轉軟化、科氏力以及三者耦合影響下實體單元和殼單元動頻對比

從以上分析可知，除較為復雜的彎曲耦合振型外，采用變厚度殼單元可以很好地模擬真實葉片的動頻特性，靜頻的前6階固有頻率與實體單元計算誤差都小于1%，證明了變厚度殼單元的精度很高。而且由于殼單元模型的單元數和節點數要比實體單元少很多，所以殼單元模型能節省很多計算時間，因此在高效性上，殼單元要優于實體單元。當然在考慮某些因素對動頻的影響時，實體單元和殼單元還存在一些微小差異，如：考慮旋轉軟化效應時，對比第1階固有頻率，對實體單元的影響更為明顯。

4 常用動頻系數經驗公式適用性分析

文獻[15-16]給出動頻fd和靜頻fs之間的關系式

式中：B為葉片動頻系數；fr為轉動頻率。

對于葉片的A0型(1階彎曲)振動的動頻系數常用公式有[15-16]式中：Dm為葉片的平均直徑；l為葉片長度其中 β0、β1分別為葉根、葉頂振動平面與葉輪平面的夾角。

文獻[17]根據試驗測試結果，給出了葉根剛性固定的變截面葉片A0型動頻率系數

對于高頻彎曲振動，由于振動頻率較高，振幅較小，所以離心力對振動頻率影響不大。高階彎曲振動的動頻系數計算公式較少，僅文獻[18]給出了葉片的A1型振動（2階彎曲）的動頻系數公式

本文葉片的幾何參數為：Dm=521.64mm，l=88.6 mm，β=36°，R0=216.52 mm，A1=77.72 mm2，A0=263.33 mm2。將葉片參數代入式（9）～（11）可獲得葉片動頻率系數分別為4.417、5.148和31.319。

通過對殼單元僅考慮離心剛化條件下計算的第1階（A0型振動：1階彎曲振動）和第3階（A1型振動：2階彎曲振動）仿真數據，利用公式計算動頻系數，如圖6所示。從圖中可見，計算結果為2條直線，即第1、2階彎曲振動的動頻系數為一定值，直線斜率即為動頻系數，從而確定前2階動頻系數分 別 為 B1=5.755 和B2=17.491。

圖6 動靜頻和轉頻關系曲線

對比仿真結果和經驗公式可知，對于A0型動頻系數，二者吻合較好；而對于A1型動頻系數，二者差別很大。因而還需要仔細確認其適用范圍和適用條件。

5 結論

本文主要基于ANSYS有限元軟件，提出了采用變厚度殼單元來模擬真實變截面、帶扭型的葉片建模方法，得到以下主要結論。

（1）本文所建模型與實體單元模型通過對靜頻和固有振型進行對比分析，驗證了本建模方法的正確性。

（2）分析了科氏力、旋轉軟化、離心剛化以及三者耦合對實體建模和殼建模葉片動頻特性的影響規律，結果表明，2種模型的動頻規律一致，離心剛化對系統的影響最大，其次為旋轉軟化效應，而科氏力對動頻的影響最小，在計算時可忽略。

（3）基于殼單元離心剛化的動頻數據，與常用的動頻系數經驗公式進行了簡單對比，A0型的動頻系數吻合較好，而A1型的動頻系數吻合較差。對其適用范圍和適用條件還需進一步探討。

[1]王建軍，卿立偉，李其漢.旋轉葉片頻率轉向與振型轉換特性[J].航空動力學學報，2007，22（1）：8-11. WANG Jianjun,QING Liwei,LIQihan.Frequency veering and mode shape interaction for rotating blades [J].Journal of Aerospace Power,2007，22（1）：8-11.(in Chinese)

[2]Hashemi S M,Richard M J.Natural frequencies of rotating uniform beams with coriolis effects[J].Journal of Vibration and Acoustics,2001,123（4）：444-455.

[3]Huang B W,Kuang J H.Mode localization in a rotating mistuned turbo disk with Coriolis effect[J].International Journal of Mechanical Sciences, 2001, 43 （7）：1643-1660.

[4]Lee C L,Alsalems M F,Woehrle T G.Natural frequency measurements for rotating spanwise uniform cantilever beams[J].Journal of Sound and Vibration,2001,240 （5）：957-961.

[5]鄧峰巖，和興鎖，張娟，等.旋轉Timoshenko梁的動力學分析[J].航空學報，2006，27（6）：1092-1096. DENG Fengyan,HE Xingsuo,ZHANG Juan,et al. Dynamic analysis of rotating timoshenko beams[J].Acta Aeronautica etAstronautica Sinica,2006,27（6）：1092-1096.(in Chinese)

[6]韓廣才，吳艷紅，王寅超，等.旋轉葉片剛柔耦合系統動力學分析[J].哈爾濱工程大學學報，2011，32（6）：736-741. HAN Guangcai,WU Yanhong,WANG Yinchao,et al. Dynamic analysis of a rotating blade with a rigid-flexible coupling,pre-twisted angle,and variable cross-section [J].Journal of Ha rbin Engineering University,2011,32 （6）：736-741.(in Chinese)

[7]Young T H,Liou G T.Coriolis effect on the vibration of a cantilever plate with time-varying rotating speed[J]. Journal of Vibration and Acoustics,1992,114（2）：232-241.

[8]Yoo H H,Park J H.Vibration analysis rotating pre-twisted blades[J].Computers and Structures,2001, 79：1811-1819.

[9]Sunil S K,Turner K E.Natural frequencies of a pre-twisted blade in a centrifugal force field[J].Journal of Sound and Vibration,2011,330：2655-2681.

[10]李紅影，謝里陽，郭星輝.殼板葉片的固有特性分析[J].東北大學學報（自然科學版），2011，32（11）：1627-1630. LIHongying,XIE Liyang,GUO Xinghui.Inherent characteristics of the cylindrical panel blade[J].Journal of Northeastern University (Natural Science),2011,32 （11）：1627-1630.(in Chinese)

[11]Subrahmanyam K B,Kaza K R V.Vibration and buckling ofrotating,pretwisted,preconed beams including Corioliseffects [J].Journal of Vibration, Acoustics,Stress and Reliability in Design,1986,108：140-149.

[12]Bhat M M,Ramamurti V,Sujatha C.Studies on the determination ofnaturalfrequencies ofindustrial turbine blades [J].Journal of Sound and Vibration, 1996,196（5）：681-703.

[13]Tsai G C.Rotating vibration behavior of the turbine blades with different groups of blades[J].Journal of Sound and Vibration,2004,271（3）：547-575.

[14]Nikolic M,Petrov E P,Ewins D J.Coriolis forces in forced response analysis of mistuned bladed disks[J]. Journal of Turbomachinery,2007,129（4）：730-739.

[15]劉尚明，孟繁娟，何光新.對葉片動頻系數經驗公式的評價[J].汽輪機技術，2001，43（1）：14-16. LIU Shangming,MENG Fanjuan,HE Guangxin.The evaluation on experimental formula of blade rotational frequency coefficient[J].Turbine Technology,2001,43 （1）：14-16.(in Chinese)

[16]華東電力試驗研究所，江蘇電力試驗研究所，南京工學院.汽輪機葉片的振動特性和調整 [M].北京：電力工業出版社，1981：115-118. East China Electric Power Test and Research Institute, Jiangsu Electric Power Test,Research Institute and Nanjing Institute of Technology. Vibration characteristics and adjustment of the steam turbine

Natural Characteristic Analysis of Rotating Blade Based on Variable Thickness Shell Element

MA Hui,NAI Hai-qiang,SONG Rong-ze,ZHANG Zhi,WEN Bang-chun
（School of Mechanical Engineering&Automation,Northeastern University,Shenyang 110819,China）

Taking a rotating variable cross-section twisted blade as the research object,the finite element model of the blade was established using variable thickness shell element based on ANSYS software.The comparison between static frequencies and their model shapes was conducted to validate the accuracy of modeling method.Then the influences law of Coriolis force,spin softening,centrifugal stiffening and coupling of them on dynamic frequencies of the blade were analyzed based on the model.Finally,the common empirical formulas of dynamic frequency coefficients were evaluated based on the dynamic frequency data of the shell element model.The study shows that the law of dynamic frequencies

from the shell element model is in good agreement with that from a solid element model.The effects are listed from large to small as follows:centrifugal stiffening,spin softening and Coriolis force.In addition,dynamic frequency coefficients of the A0 vibration mode are in good agreement with the results from other references,whereas the result of A1 vibration mode is just the reverse of that of A0.

rotating blade;natural characteristics;shell element;ANSYS;dynamic frequency

馬輝（1978），男，副教授，主要研究方向為故障診斷、非線性振動和轉子動力學。

中央高?；究蒲袠I務費專項基金（N100403008）、教育部新世紀人才支撐計劃（NCET-11-0078）資助

2012-12-28