郭 鵬， 劉志遠(yuǎn)， 張桂昌， Reza Hedayati， 張俊紅

(1. 天津理工大學(xué) 天津市先進(jìn)機(jī)電系統(tǒng)設(shè)計與智能控制重點(diǎn)實驗室，天津 300384 2.天津大學(xué) 內(nèi)燃機(jī)燃燒學(xué)國家重點(diǎn)實驗室，天津 300350；3. 中國民航大學(xué) 航空工程學(xué)院，天津 300300；4. 代爾夫特理工大學(xué) 航空航天工程學(xué)院，荷蘭 代爾夫特 2600AA)

鳥撞對飛機(jī)運(yùn)行安全有重要影響。根據(jù)美國聯(lián)邦航空管理局(federal aviation administration, FAA)統(tǒng)計[1]，從1988年到2018年間，野生動物與飛機(jī)相撞已造成超過263架飛機(jī)損毀，287人喪生。圖1(a)為從1990-2017年間野生動物撞擊飛機(jī)的事故統(tǒng)計，2017年較前年增加1 069起，其中的鳥撞所占的比例達(dá)到95%。圖1(b)為飛機(jī)各部件發(fā)生撞鳥的比例[1]，其中發(fā)動機(jī)比例最高達(dá)到37.21%。鳥撞時，風(fēng)扇最先受到撞擊，其抗鳥撞能力成為影響飛行安全的重要因素。研究鳥撞風(fēng)扇過程對于發(fā)展風(fēng)扇安全強(qiáng)度設(shè)計理論具有重要的科學(xué)意義，對于風(fēng)扇葉片適航認(rèn)證具有主要的工程意義。

圖1 野生動物撞擊飛機(jī)次數(shù)統(tǒng)計及各部件發(fā)生撞鳥比例

國內(nèi)外在鳥撞模擬中常使用圓柱體、半球端圓柱體、橢圓體和球體來替代真實的鳥體，該方法可以反映鳥軀干的主要質(zhì)量和形狀對碰撞過程的影響[2]，但忽略了其他部位對碰撞過程的影響。近年來，國內(nèi)外開始使用更加復(fù)雜的鳥模型進(jìn)行研究。Hedayati等[3]建立了真實綠頭鴨的SPH模型，研究分析鳥的幾何形狀和撞擊方向?qū)B撞平板過程的影響，該模型相比于傳統(tǒng)的替代鳥模型更接近于試驗結(jié)果。此后Hedayati等[4]使用多個簡單幾何體對白枕鵲鴨的形狀進(jìn)行簡化，建立了一種在幾何形狀類似于白枕鵲鴨的新的替代鳥模型，并且發(fā)現(xiàn)半球端圓柱體模型更接近真實鳥模型以尾部方向撞擊平板的計算結(jié)果，橢圓體模型可以作為腹部撞擊方向的替代模型。McCallum等[5-6]利用SPH方法，基于生物識別和已發(fā)布的CT掃描數(shù)據(jù)建立了一種多材料鳥模型，該模型的數(shù)值模擬結(jié)果與經(jīng)典流體動力學(xué)理論非常吻合。Zhang等[7-8]建立帶有更多幾何特征的帆背潛鴨模型，該模型沖擊剛性靶板的計算結(jié)果與試驗吻合較好。并使用帆背潛鴨模型撞擊靜止風(fēng)扇，其研究表明鳥的幾何以及撞擊方向?qū)B撞風(fēng)扇有著重要且不可忽略的影響。寇建峰等[9]針對雷達(dá)罩以及機(jī)翼前緣開展不同姿態(tài)鳥體的鳥撞研究，發(fā)現(xiàn)鳥體姿態(tài)對結(jié)構(gòu)撞擊響應(yīng)有很大的影響，并且提出結(jié)構(gòu)的抗鳥撞安全評估必須考慮不同鳥體姿態(tài)的影響。楊杰[10]模擬鳥撞整級旋轉(zhuǎn)風(fēng)扇的不同位置，結(jié)果表明撞擊位置越高造成的損傷越大，損傷葉片損失的能量越大。

鳥撞葉片過程是一個高度非線性的瞬態(tài)沖擊過程，上述研究表明撞擊姿態(tài)以及撞擊位置對葉片損傷有著不可忽略的影響，同時葉片的扭轉(zhuǎn)結(jié)構(gòu)使得鳥撞風(fēng)扇葉片過程更加復(fù)雜。現(xiàn)有研究中，對真實鳥與真實航空發(fā)動機(jī)風(fēng)扇葉片碰撞過程的研究較少，且簡化鳥體撞擊葉片過程中僅考慮簡化鳥體的吸入速度，未考慮風(fēng)扇旋轉(zhuǎn)速度造成鳥體-葉片相對運(yùn)動速度的影響。本文基于真實綠頭鴨SHP模型，模擬鳥撞擊固定轉(zhuǎn)速下的某航空渦扇發(fā)動機(jī)整級風(fēng)扇過程，分析不同撞擊位置以及不同撞擊姿態(tài)對風(fēng)扇的動力學(xué)響應(yīng)及損傷變形情況的影響。

1 數(shù)值模型

1.1 網(wǎng)格模型建立

FAA統(tǒng)計數(shù)據(jù)表明[1]，在1990—2017年綠頭鴨撞擊飛機(jī)的總次數(shù)達(dá)到1 003次，為與民航發(fā)動機(jī)發(fā)生碰撞次數(shù)最頻繁的鳥類之一。因而本文將綠頭鴨作為鳥撞模擬的原型，鴨重820 g，體長40 cm，翼展45 cm。在醫(yī)學(xué)成像中心對綠頭鴨進(jìn)行CT掃描，如圖2(a)，獲取1 566張醫(yī)學(xué)數(shù)字成像(DICOM)圖像，如圖2(b)。利用DICOM圖像及SPH方法對綠頭鴨進(jìn)行建模，建模過程中考慮綠頭鴨的內(nèi)部空腔(空腔部位沒有SPH單元)，最終模型由41 685個SPH單元組成，每個單元的質(zhì)量為0.019 1 g，如圖2(c)。

圖2 綠頭鴨建模過程

國內(nèi)外鳥撞數(shù)值模擬中替代鳥體模型目前普遍采用半球端圓柱體模型，但是半球端圓柱體模型只是對鳥體軀干的簡化，而沒有將鳥體其他部位考慮進(jìn)來，不是對鳥體真實外形的建模。為了驗證真實鳥模型相比于傳統(tǒng)替代鳥體模型具有更準(zhǔn)確的數(shù)值模擬結(jié)果，本文選取傳統(tǒng)替代鳥體模型(半球端圓柱體模型)作為真實綠頭鴨模型的對比。半球端圓柱體的長徑比選取為2[11]，根據(jù)鳥的質(zhì)量0.8 kg以及鳥的密度938 kg/m3可以得到半球端圓柱體的直徑D=0.086 9 m和長度L=0.173 8 m。利用SPH法對半球端圓柱體離散化，粒子距離為3.2 mm，模型由32 038個SPH粒子組成每個粒子質(zhì)量為0.024 97 g,如圖3。

圖3 傳統(tǒng)替代鳥體SPH模型

本文研究的發(fā)動機(jī)風(fēng)扇由24塊窄弦葉片組成，葉身高度為603.2 mm,初始扭轉(zhuǎn)角為61.3°。用3 mm六面體單元對整級風(fēng)扇進(jìn)行有限元網(wǎng)格劃分，獲得1 520 208網(wǎng)格單元，如圖4所示。

圖4 渦扇發(fā)動機(jī)風(fēng)扇有限元網(wǎng)格模型

1.2 材料模型

鳥材料采用Null材料模型結(jié)合Gruneisen狀態(tài)方程來描述；在空材料參數(shù)中提供模型的本構(gòu)關(guān)系計算黏性應(yīng)力，使用狀態(tài)方程來計算壓力，其中材料密度設(shè)為為938 kg/m3。Gruneisen狀態(tài)方程按下式定義壓縮材料的壓力：

(1)

對于膨脹材料，Gruneisen狀態(tài)方程按下式定義壓力：

p=ρ0C2μ+(γ0+aμ)E

(2)

在上述Gruneisen狀態(tài)方程中C是vs(vp)曲線的截距；S1，S2，和S3是vs(vp)曲線的斜率系數(shù)[12]；γ0是Gruneisen常數(shù)，a是的一階體積修正量。

(3)

式中：C=1 480，S1=1.92，S2=0，S3=0。

風(fēng)扇葉片采用的材料是Ti-6Al-4V鈦合金材料，采用塑性動力學(xué)模型，材料參數(shù)見表1。

表1 葉片材料參數(shù)

1.3 邊界條件

沿葉根向上選取1/6、2/6、3/6、4/6、5/6葉高五個撞擊高度，研究撞擊位置對鳥撞風(fēng)扇葉片過程的影響，如圖5所示。

圖5 鳥體撞擊位置

分別選取15種鳥體吸入姿態(tài)研究鳥體撞擊姿態(tài)對撞擊過程的影響。如圖6所示，自然飛行頭向風(fēng)扇為0°姿態(tài)，鳥體繞Y軸依次旋轉(zhuǎn)45°得到Y(jié)-45°、Y-90°、Y-135°、Y-180°、Y-225°、Y-270°、Y-315°姿態(tài)；鳥體繞Z軸依次旋轉(zhuǎn)45°得到Z-45°、Z-90°、Z-135°、Z-180°、Z-225°、Z-270°、Z-315°姿態(tài)。

圖6 撞擊姿態(tài)示意圖

模擬中，鳥的飛行速度為116 m/s，風(fēng)扇轉(zhuǎn)速為3 800 r/min。以往鳥撞風(fēng)扇研究大多未考慮風(fēng)扇的轉(zhuǎn)動，鳥沿風(fēng)扇軸線撞向葉片，鳥撞過程中只與一或兩個葉片接觸，發(fā)生正向碰撞。實際上，鳥撞向旋轉(zhuǎn)的風(fēng)扇，會與多個葉片發(fā)生接觸，并且風(fēng)扇葉片和鳥的相互作用包括撞擊和對鳥體的切割兩部分。為更加真實的模擬鳥撞風(fēng)扇的情況，鳥的速度設(shè)為該位置處相對于旋轉(zhuǎn)風(fēng)扇的相對速度[13]，如圖7。相對速度的大小計算公式為：

圖7 鳥相對葉片速度

(4)

v葉片=rω

(5)

式中：ω為風(fēng)扇的轉(zhuǎn)速；r為撞擊位置旋轉(zhuǎn)半徑；v相對為風(fēng)扇葉片撞擊位置處的切向速度。

2 模型驗證

本文通過與Wilbeck試驗數(shù)據(jù)[14]進(jìn)行對比來進(jìn)行模型準(zhǔn)確性驗證。Wilbeck采用壓縮空氣炮裝置將質(zhì)量為1 kg的鳥體分別以116 m/s、225 m/s和253 m/s的速度射出，與剛性目標(biāo)板撞擊，在剛性平板中心位置安裝壓力傳感器以獲取撞擊中心位置的壓力時間曲線，實驗簡圖，見圖8。

1.氣罐及空氣炮； 2.鳥彈及彈底板； 3.炮管； 4.壓力傳感器； 5.剛性目標(biāo)板； 6.彈底板止動裝置。

進(jìn)行半球端圓柱體鳥模型與綠頭鴨真實鳥模型撞擊平板的模擬，分別獲取116 m/s、225 m/s、253 m/s的不同速度撞向固定平板時平板中心的壓力曲線，并與對應(yīng)速度下Wilbeck的試驗數(shù)據(jù)進(jìn)行比較。

本文建立的平板模型為直徑0.6 m厚度0.06 m的圓盤，采用六面體實體網(wǎng)格對圓盤進(jìn)行離散化，獲得18 000網(wǎng)格單元。該數(shù)值模擬中，綠頭鴨以尾部垂直撞向平板。通過固連接觸將殼單元節(jié)點(diǎn)約束、限定在圓盤外表面節(jié)點(diǎn)上，以獲取撞擊中心的壓力曲線。平板模型和殼單元的布置位置及其相對于鳥模型的位置，如圖8。

圖9 鳥體撞擊平板

Heimbs研究結(jié)果表明，鳥撞剛性目標(biāo)板撞擊中心壓力響應(yīng)分為四個階段。第一階段為初始撞擊階段，初始撞擊階段的峰值壓力稱為雨貢紐壓力，第二階段為壓力衰減階段，第三階段為恒定流動階段，第四階段為壓力結(jié)束階段。不同撞擊速度下，碰撞中心壓力曲線仿真與試驗對比如圖10所示，可以清晰的看出相比于半球端圓柱體模型，綠頭鴨模型各階段壓力響應(yīng)與Wilbeck實驗結(jié)果吻合更好。尤其是在初始撞擊階段吻合更好，對比116 m/s、225 m/s以及253 m/s撞擊速度下的撞擊中心峰值壓力，綠頭鴨鳥體模型的峰值壓力都要比半球端圓柱體模型更加接近Wilbeck實驗的峰值壓力。造成這個結(jié)果的主要原因是綠頭鴨模型與實際鳥外形相似，并且同時考慮到綠頭鴨的內(nèi)部空腔，對空腔進(jìn)行無粒子化處理，使得模型更加接近真實的鳥，而半球端圓柱體只是對鳥軀干的一個近似的模擬，并且半球端圓柱體模型內(nèi)部并沒有空腔，因而綠頭鴨SPH模型結(jié)果相比于傳統(tǒng)的半球端圓柱體模型更加接近Wlibeck的實驗結(jié)果，真實綠頭鴨SPH模型可更真實反應(yīng)鳥撞動力學(xué)特性，可進(jìn)行下一步研究。

圖10 撞擊中心處的壓力曲線

3 結(jié)果與討論

參考Zhang等分析鳥撞風(fēng)扇過程中風(fēng)扇損傷的判據(jù)，為了分析撞擊位置以及撞擊姿態(tài)對葉片不同位置損傷的影響，本文也將討論分析葉根及前緣的等效應(yīng)力、撞擊力及前緣最大位移。葉根等效應(yīng)力用于判斷葉根處得損傷狀況，前緣等效應(yīng)力以及前緣最大位移用于判斷前緣損傷狀態(tài)，前緣是撞擊過程中出現(xiàn)大變形的位置因而增加前緣最大位移來表示前緣的變形情況。

3.1 撞擊位置對風(fēng)扇葉片損傷的影響

圖11為撞擊過程中不同時刻鳥體狀態(tài)，撞擊過程中鳥體被葉片切割成鳥塊，隨后鳥塊與葉片發(fā)生碰撞，鳥塊在與葉片高速的碰撞過程中發(fā)生流變。

圖11 鳥撞風(fēng)扇過程中不同時刻鳥的狀態(tài)

圖12(a)為鳥與風(fēng)扇的接觸力時間曲線，每條接觸力-時間曲線有多個波峰且每個波峰的高度不相同，接觸力-時間曲線波峰數(shù)量與鳥被切割成鳥塊的數(shù)量相同，且每條接觸力-時間曲線每個波峰的大小與葉片撞擊的鳥塊大小相關(guān)。撞擊位置為4/6、5/6葉高時整體接觸力要小于撞擊位置為1/6、2/6、3/6葉高的接觸力，如表2所示，1/6、2/6、3/6、4/6、5/6位置接觸力峰值為119 kN、155 kN、109 kN、67.2 kN、35.2 kN。從葉根位置至2/6葉高位置葉片的扭轉(zhuǎn)幅度很小，接觸力主要受撞擊速度影響，撞擊高度越大，撞擊位置葉片線速度越大，相對撞擊速度越大，撞擊接觸力越大。從2/6葉高位置至葉尖位置，葉片扭轉(zhuǎn)幅度增大，接觸力主要受葉片扭轉(zhuǎn)角度影響，撞擊高度越大,葉片扭轉(zhuǎn)角越大，葉片對鳥體的切割作用更加明顯，葉片法向的拍擊減弱使得接觸力逐漸減小。

表2 不同撞擊位置的最大接觸力、葉根等效應(yīng)力和前緣等效應(yīng)力

文獻(xiàn)[15]計算了圓柱鳥模型以102 m/s速度撞擊轉(zhuǎn)速3 500 r/min風(fēng)扇3/6葉高過程中葉片瞬態(tài)沖擊響應(yīng)特性，結(jié)果表明，葉片峰值應(yīng)力1.396 GPa, 在本文研究中鳥體以116 m/s速度撞擊轉(zhuǎn)速3 800 r/min風(fēng)扇3/6葉高過程中葉片峰值應(yīng)力為1.53 GPa。由于葉片模型、鳥體模型及計算工況的差異，本研究葉片峰值應(yīng)力計算結(jié)果與文獻(xiàn)[15]計算結(jié)果略有差異但處于同一量級且數(shù)值接近。因此可以認(rèn)為本研究中等效應(yīng)力計算結(jié)果合理可信，可用于進(jìn)一步討論分析。

圖12(b)為不同撞擊位置下葉根應(yīng)力-時間曲線，1/6、2/6、3/6、4/6、5/6葉高位置葉根最大等效應(yīng)力為1.23 GPa、1.24 GPa、1.10 GPa、0.969 GPa、與0.479 GPa。圖12(c)為不同撞擊位置最大鳥塊撞擊葉片的前緣等效應(yīng)力-時間曲線，1/6、2/6、3/6、4/6、5/6葉高位置前緣最大等效應(yīng)力為1.18 GPa、1.62 GPa、1.53 GPa、1.03 GPa與0.479 GPa。不同撞擊位置下葉根及葉片前緣處最大等效應(yīng)力變化規(guī)律與接觸力變化規(guī)律相似，其產(chǎn)生原因也與接觸力變化原因相同。

圖12 接觸力、葉根應(yīng)力及前緣應(yīng)力

除5/6葉高位置外，其他撞擊位置葉根及葉片前緣處最大等效應(yīng)力都大于材料的屈服極限0.95 GPa，會在葉根及葉片前緣處產(chǎn)生塑性變形造成不可逆損傷，撞擊位置為2/6葉高時，葉片葉根、前緣處損傷達(dá)最大。

3.2 撞擊姿態(tài)對風(fēng)扇葉片損傷的影響

圖13為鳥以Y-315°姿態(tài)撞擊葉片過程中不同時刻的應(yīng)力分布，葉片最大應(yīng)力主要出現(xiàn)在前后緣、前后葉根位置。相應(yīng)區(qū)域在撞擊過程中最可能發(fā)生塑性變形、撕裂等損傷[16]。

圖13 不同時刻葉片應(yīng)力分布圖

圖14(a)，為繞Y軸旋轉(zhuǎn)得到的7個姿態(tài)的鳥撞擊風(fēng)扇過程中的接觸力時間曲線。繞Y-0°姿態(tài)的接觸力峰值為105 kN，在所有姿態(tài)中最小。Y-270°和Y-315°姿態(tài)的接觸力峰值要選大于其他姿態(tài)，接觸力峰值都是229 kN。圖14(b)，為繞Z軸旋轉(zhuǎn)得到的7個姿態(tài)的鳥撞擊風(fēng)扇過程中的接觸力時間曲線，鳥體Z軸姿態(tài)對接觸力峰值影響較小，Z-45°、90°、135°、180°、225°、270°、315°姿態(tài)的接觸力峰值接近，分別為153 kN、138 kN、165 kN、145 kN、144 kN、158 kN、145 kN。對比15個姿態(tài)的接觸力峰值，Y-270°和Y-315°姿態(tài)的接觸力最大。

圖14 不同姿態(tài)接觸力

圖15為各撞擊姿態(tài)下，與最大鳥塊相撞的葉片前葉根的等效應(yīng)力，各撞擊姿態(tài)下前葉根最大等效應(yīng)力如表3所示，Y-135°、Y-270°、Y-315°、Z-135°和Z-315°撞擊姿態(tài)下，撞擊過程中前葉根最大等效應(yīng)力分別為1.99 GPa、1.98 GPa、2.00 GPa、1.98 GPa和2.04 GPa，高于其他姿態(tài)。圖16為Y軸8種姿態(tài)和Z軸8種姿態(tài)下，與最大鳥塊相撞的葉片后葉根的等效應(yīng)力時間曲線，各撞擊姿態(tài)下后葉根在撞擊過程中最大等效應(yīng)力如表3所示，Z-135°姿態(tài)在撞擊過程中后葉根受到的等效應(yīng)力最大，為1.94 GPa。在撞擊過程中前葉根受到的應(yīng)力普遍大于后葉根，因而在撞擊過程中前葉根受到的損傷大于后葉根。

表3 不同姿態(tài)最大接觸力、前后葉根等效應(yīng)力與前緣接觸處最大位移

圖15 不同姿態(tài)前葉根等效應(yīng)力

圖16 不同姿態(tài)后葉根等效應(yīng)力

選用與最大鳥塊相撞的葉片前緣接觸處節(jié)點(diǎn)最大位移表征鳥撞風(fēng)扇過程中與最大鳥塊撞擊的葉片前緣的變形情況。通過比較發(fā)現(xiàn)，Y-270°姿態(tài)葉片前緣的變形要明顯大于其他所有位置，前緣節(jié)點(diǎn)位移最大值為52.3 mm。

在鳥撞風(fēng)扇葉片過程中，鳥體姿態(tài)主要影響鳥體被切割次數(shù)、鳥體被切割后各鳥塊的質(zhì)量、鳥塊與葉片的接觸位置以及接觸面積。這些因素進(jìn)一步影響葉片的應(yīng)力響應(yīng)及撞擊過程對葉片的損傷大小。Y-135°、Y-270°、Y-315°、Z-135°及Z-315°撞擊姿態(tài)對前葉根損傷最大，Z-135°撞擊姿態(tài)對后葉根損傷最大，Y-270°撞擊姿態(tài)對前緣的損傷最大。0°撞擊姿態(tài)下，接觸力、前葉根應(yīng)力、后葉根應(yīng)力及前緣位移為最小，對葉片前后葉根及前緣的損傷最小。在撞擊過程中前葉根受到的應(yīng)力要普遍大于后葉根，因而在撞擊過程中前葉根相對于后葉根受到的損傷要更大。

4 結(jié) 論

本文采用真實綠頭鴨SPH模型，模擬鳥撞整級旋轉(zhuǎn)風(fēng)扇葉片過程，對鳥撞位置和撞擊姿態(tài)對葉片損傷的影響進(jìn)行研究，得到如下結(jié)論：

(1) 鳥撞擊高速旋轉(zhuǎn)風(fēng)扇的過程中會與多個葉片發(fā)生相互作用，被多個葉片切割為不同大小的鳥塊，鳥體被切割后與葉片發(fā)生碰撞發(fā)生流變。與最大鳥塊發(fā)生碰撞的葉片應(yīng)力最大。鳥撞過程中，葉片最大等效應(yīng)力出現(xiàn)位置為葉片前緣、后緣以及葉根。

(2) 從葉根到2/6葉高區(qū)域，撞擊的位置越高，撞擊力、葉根處應(yīng)力以及葉片前緣應(yīng)力越大。從2/6葉高到葉尖區(qū)域，撞擊的位置越低，葉片受到的撞擊力、葉根處及前緣處的應(yīng)力越大，葉根及前緣受到損傷越大。

(3)Y-135°、Y-270°、Y-315°、Z-135°及Z-315°撞擊姿態(tài)對前葉根損傷最大，Z-135°撞擊姿態(tài)對后葉根損傷最大，Y-270°撞擊姿態(tài)對前緣的損傷最大。0°撞擊姿態(tài)下，接觸力、前葉根應(yīng)力、后葉根應(yīng)力及前緣位移為最小，對葉片前后葉根及前緣的損傷最小。

(4) 在撞擊過程中前葉根應(yīng)力要普遍大于后葉根，因而在撞擊過程中前葉根相對于后葉根受到的損傷更大。