基于流固耦合的實體元空心葉片鳥撞數值模擬

2011-04-27 07:45:46劉建明蔣向華馬永峰

航空發(fā)動機 2011年2期

劉建明，蔣向華，武卉，馬永峰，王東

(1.北京航空航天大學能源與動力工程學院，北京 100191；2.中航工業(yè)沈陽發(fā)動機設計研究所，沈陽 110015)

1 引言

鳥撞擊葉片是航空發(fā)動機較為普遍的1種外物損傷事故。對此，發(fā)動機研制技術先進國家開展了許多相關研究，并制定了相應的抗鳥撞的設計規(guī)范和試驗規(guī)范;中國也對飛機抗鳥撞試驗提出了明確要求[1]。

鳥撞研究可分為理論、試驗和數值模擬研究。理論研究包括鳥體切割模型、鳥撞物理過程、鳥體載荷模型等[2]；試驗研究包括國外進行的對葉片模型件進行投鳥試驗[3]，國內采用高壓空氣炮及高速攝影等模擬鳥體對飛機進行撞擊試驗;數值模擬研究主要包括各種算法及各種有限元程序的開發(fā)與應用。目前，國內的研究院所開展了較多鳥撞發(fā)動機葉片的數值模擬研究，如用ANSYS/DYNA軟件，采用接觸沖擊算法進行鳥撞葉片數值模擬[4]；用MSC.DYTRAN軟件，采用流固耦合算法進行鳥撞葉片研究[1]等。在采用流固耦合算法進行葉片鳥撞研究中，普遍采用殼單元劃分葉片，而殼單元在模擬復雜幾何結構方面存在困難。如果采用實體元劃分葉片，由于在失效單元的處理及耦合面的推進等方面存在諸多困難，而且目前流固耦合算法還缺乏完善的理論指導,所以，在現有計算程序中一般只計算殼單元的耦合及失效。這就給具有失效模式的實體元劃分的空心葉片的鳥撞數值模擬帶來了很大挑戰(zhàn)。

本文利用PATRAN軟件建立了鳥撞6面體元劃分的平板葉片的有限元模型，研究了不考慮和考慮平板葉片失效2種情況下的數值模擬方法；建立了鳥撞6面體元劃分的空心葉片轉子級瞬態(tài)動力學有限元模型，研究鳥體密度、長度、半徑、速度的改變對葉片瞬態(tài)響應的影響；最后進行了具有失效模式的空心葉片鳥撞數值模擬。

2 鳥撞平板葉片方法

2.1 不考慮平板失效的方法

為了驗證本文所采用方法的可行性，建立了鳥撞鋁板模型。鋁板尺寸為410 mm×500 mm×10 mm，4邊固支，用分段線性塑性材料模擬，材料參數見表1。

表1 鋁板材料參數

鳥體采用線性流體模擬，其質量M=1.8 kg，密度ρ=928.15 kg/m3，體積模量 K=2200 MPa，用兩端半球中間圓柱實體來模擬，以30°角斜撞向鋁板中央。

鋁板用8節(jié)點6面體元劃分，沿厚度方向劃分2層單元；以鋁板體單元的上、下自由表面的節(jié)點為公共節(jié)點，分別建立了1層厚度為0.001 mm的4節(jié)點4邊形殼元；在板的厚度方向的4個面上建立了1圈啞元，如圖1所示。啞元和殼元形成包圍體單元的封閉的耦合面，如圖2所示。

為了驗證鋪薄殼元對計算結果的影響，計算了體單元表面未鋪和鋪薄殼元時板中心點的位移隨時間的變化情況，如圖3所示。從圖中可見，位移響應曲線幾乎重合，鋪薄殼元對計算結果的影響很小。

2.2 考慮平板失效的方法

建立了1只鳥體撞擊2塊平板葉片的有限元模型，此計算是為模擬鳥體撞擊破壞第1塊板后撞擊第2塊板的情形。參照第2.1節(jié)，2塊板表面由殼單元和啞元分別形成1個包圍平板葉片的耦合面；在上板底端和下板上端之間建立啞元，此啞元和上板下端的殼元以及下板上端的殼元形成1個耦合面。鳥體數據與第2.1節(jié)的相同，平板葉片材料與第2.1節(jié)相比只添加了最大塑性應變破壞準則。定義2塊板的最大塑性應變均為0.01，殼的最大塑性應變?yōu)?.02，計算得到的鳥撞瞬間2塊板應力云圖如圖4所示。從圖中可見，2塊板被撞壞。

由2塊平板葉片的鳥撞破壞模擬可知，要進行具有破壞模式、實體元劃分的空心葉片的鳥撞數值模擬，可以采用在空心葉片表面鋪殼元的方法。

3 不考慮失效的葉片鳥撞模擬

3.1 有限元模型

8節(jié)點6面體空心葉片模型采用PATRAN建模。所建空心葉片模型每個葉片沿葉高有14層節(jié)點，共13層單元。為了節(jié)省計算時間，模型只建立了在發(fā)動機中與鳥體有初始碰撞行為的第1級葉片，由13片葉片組成，每片葉片有節(jié)點2044個、單元1122個；為了進行流固耦合計算，以葉片體單元自由表面上的節(jié)點為公共節(jié)點，建立覆蓋在體表面的4節(jié)點的殼元，其厚度為0.001 mm。所建有限元模型如圖5所示。

3.2 相關數據

葉片材料的靜態(tài)應變強化效應采用雙線性應力-應變模型模擬，應力循環(huán)效應采用等向強化模型模擬，動態(tài)應變率效應采用Cowper-Symonds模型模擬。其硬化方程為

式中：σd為動態(tài)應力；σy為靜態(tài)屈服應力；ε˙為等效應變率；D、P為表征動態(tài)硬化的材料系數，計算時取D=100、P=10。

葉片的材料密度ρ=4430 kg/m3，泊松比μ=0.34，彈性模量E=112500MPa，靜態(tài)屈服應力σy=900MPa，靜態(tài)硬化模量EH=14286MPa。葉片的轉速為12000 r/min。

鳥體的幾何尺寸（包括半徑和長度）和初始速度見表2。考慮到鳥體在碰撞過程中的流變假設，即當鳥體相對速度大于50 m/s時，會發(fā)生流變；鳥撞一般發(fā)生在飛機起降階段，飛機與鳥的相對速度不大。因此，選取了2種速度（70、100 m/s）進行對比計算。

表2 鳥體的幾何尺寸及速度

鳥體采用線性流體模型，選取3種材料計算，見表3。

表3 鳥體材料

3.3 葉片的預應力分析

鳥撞事故發(fā)生時，發(fā)動機一般處于穩(wěn)定旋轉狀態(tài)，葉片由于旋轉引起的慣性力已經具有一定應力，為了得到更準確的分析結果，必須考慮預應力的影響，在瞬態(tài)分析開始前給葉片模型加上一定的應力場。這部分由NASTRAN計算得到。

3.4 鳥撞數值模擬

選取第2.2節(jié)中給定材料的葉片、鳥體密度1、鳥體半徑1、鳥體長度1、鳥體速度1進行計算，并以該計算工況作為1個標準工況。計算時間為0.002 s。

3.4.1 鳥體密度對葉片應力的影響

選取鳥體密度2、3進行對比計算。在3種密度鳥體撞擊下葉片的應力（極大值應力處，下同）及其峰值變化曲線分別如圖6、7所示。從圖6中可見，在遭受鳥體撞擊后，葉片會在較短時間內產生應力值較大的沖擊波，其峰值隨著鳥體密度的增大而增大，這是因為鳥體高速撞擊使葉片發(fā)生應變硬化及應變率硬化，硬化程度隨鳥體密度增大而增大，葉片的應力隨著其硬化程度的增大而增大；從圖7中可見，鳥體密度對葉片應力峰值的影響不是線性的，因為當鳥體密度增大到一定程度后，由于葉片的可變形能力和硬化能力已充分發(fā)揮，葉片的應力峰值的增加速度就會變慢。

圖7 不同密度鳥體撞擊下的葉片應力峰值曲線

3.4.2 鳥體長度對葉片應力的影響

選取鳥體長度2、3進行對比計算。3種長度鳥體撞擊下葉片的應力及其峰值曲線分別如圖8、9所示。從圖中可見，鳥體長度的影響與鳥體密度的影響相類似。鳥體長度對葉片應力峰值的影響不大，這是因為葉片應力峰值的產生主要來自鳥體前端與葉片發(fā)生碰撞的部位，由于應力波的作用，鳥體后端的速度明顯減弱，對葉片碰撞的影響也就大大減弱。