復合材料風扇葉片鳥撞損傷的精細化分析方法

宋長遠，唐 旭，陳大瑋，陳 勇，2

（1.上海交通大學機械與動力工程學院，2.燃氣輪機與民用航空發(fā)動機教育部工程研究中心：上海 200240）

0 引言

近年來，隨著航空發(fā)動機燃油效率的提高，風扇葉片的尺寸也不斷增大，因此將密度低、性能優(yōu)異的碳纖維增強樹脂復合材料投入風扇葉片設計制造已成為當下發(fā)展趨勢。但在鳥撞事故中，復合材料葉片比金屬葉片面臨著更大挑戰(zhàn)。首先，復合材料層間性能較低，延展性與應變能吸收率也低于金屬材料的，對鳥撞等沖擊載荷更加敏感；其次，復合材料風扇葉片受到?jīng)_擊后，損傷主要出現(xiàn)在內部鋪層上，這些肉眼不可見的損傷不但對葉片的結構完整性威脅巨大，也使得復合材料風扇葉片的沖擊損傷試驗研究比金屬葉片的更加困難。因此，針對復合材料風扇葉片鳥撞問題進行數(shù)值仿真，并在鋪層水平上展開精細化分析對復合材料風扇葉片的強度分析尤為重要。

以往的復合材料風扇葉片鳥撞數(shù)值研究在建模時對葉片簡化較多，可提取到的關于葉片的計算結果十分有限。如Chaims等基于鋪層信息和復合材料力學理論，計算單元節(jié)點等效物性，建立了復合材料風扇葉片殼單元模型，但僅分析了葉片的沖擊響應；Siddens將復合材料風扇葉片的鋪層簡化為4層，建立了殼單元模型，分析了葉片損傷及包邊對葉片損傷的影響。近年來，一些學者使用不同方法建立復合材料葉片的精細化有限元模型。Coroneos等、Xiao等分別通過殼單元厚度方向的多個積分點保留了真實葉片的鋪層信息，建立了復合材料風扇葉片的殼單元模型，但使用殼單元計算結果誤差較大，且無法分析層間失效；使用逐層建模的方法，胡殿印等基于Rotor 67風扇葉片建立了包含20層鋪層的復合材料葉片實體模型；Zhang等建立了一種類榫頭結構的實體模型；朱啟晨等提出了針對全尺寸復合材料風扇葉片的精細化有限元建模方法，模型中使用了層合實體單元，準確性得到了試驗驗證；唐旭等、張煜坤等基于文獻[8]中方法分別建立復合材料風扇葉片和榫頭模型，針對鋪層失效進行了預測。復合材料風扇葉片建模技術的發(fā)展為鳥撞損傷的精細化分析提供了條件。

本文針對復合材料風扇葉片鳥撞損傷的分析方法進行了研究。建立了復合材料風扇葉片的精細化有限元模型，提出1套復合材料風扇葉片鳥撞損傷的精細化分析流程，進行中鳥撞擊工況的數(shù)值仿真，分析了不同失效模式發(fā)生位置的規(guī)律。

1 鳥撞損傷精細化分析方法

1.1 鳥撞損傷分析流程

復合材料風扇葉片鳥撞精細化分析流程如圖1所示。復合材料風扇葉片由上百層預浸料鋪層模壓成型，其沖擊損傷可表現(xiàn)為鋪層及鋪層間的失效，因此在針對全尺寸復合材料風扇葉片的鳥撞分析中，葉片的有限元模型需能夠正確體現(xiàn)鋪層形狀、排序、角度等信息。此外，在后處理階段需將鋪層作為獨立單位展開分析。在復合材料結構強度分析中，存在全局、鋪層和材料主軸三種坐標系，在不同坐標系下的應力數(shù)據(jù)不同，復合材料的失效準則基于材料主方向應力建立，因此在進行強度分析時，需要提取材料主軸坐標系下的應力數(shù)據(jù)。在后處理階段的應力數(shù)據(jù)提取、鋪層失效分析和結果可視化環(huán)節(jié)數(shù)據(jù)量較大，需借助腳本程序實現(xiàn)自動化操作。

圖1 復合材料風扇葉片鳥撞精細化分析流程

1.2 復合材料風扇葉片精細化有限元模型

1.2.1 復合材料風扇葉片有限元建模方法

本文研究對象為一種鋪層復合材料風扇葉片，其葉尖直徑為1828.8 mm，葉片高度為727.7 mm，葉尖弦長為337.9 mm，采用直燕尾榫頭；葉片中線上的最大厚度為47.5 mm，位于榫頭上；最小厚度為8.8 mm，氣動葉身前、尾緣的厚度分別為3.4、2.1 mm。葉片鋪層材料為碳纖維增強中溫高韌性環(huán)氧樹脂預浸料，鋪層厚度為0.125 mm。葉片鋪層方案為[0°/45°/0°/-45°]，含242層結構層、206層插入層，鋪層總數(shù)為448層，葉片及部分鋪層形狀如圖2所示。

圖2 葉片及部分鋪層形狀

根據(jù)文獻[8]，首先在Fibersim軟件中以葉片中面為鋪貼面，利用體積填充功能分別向吸力面和壓力面填充鋪層形狀的鋪層，以滿足葉片的幾何外形。然后根據(jù)文獻[11]中的復合材料風扇葉片鋪層設計準則，調整鋪層角度及排序，并將鋪層信息儲存到HDF5文件中，鋪層命名原則見表1。

在ANSYS Workbench軟件中，使用ACP（Pre）模塊預先定義復合材料預浸料單向帶的材料屬性，劃分鋪貼面網(wǎng)格，最小網(wǎng)格尺寸為5 mm，復合材料環(huán)氧樹脂工程彈性常數(shù)及單向帶工程彈性常數(shù)見表2、3。

表1 鋪層命名原則

表2 環(huán)氧樹脂工程彈性常數(shù)

表3 單向帶工程彈性常數(shù)

在ACP主程序中，以鋪貼面網(wǎng)格為基準，通過法向拉伸生成葉片的有限元模型。風扇葉片鋪層使用層合實體Solid185單元建立，定義單元最大厚度為1.25 mm，即單元在厚度方向上至多包含10層鋪層信息，該單元類型沿厚度方向包含多個積分點，每層鋪層通過3個積分點表示。鋪層邊緣、葉片厚度漸變區(qū)域存在Solid185退化而成的四面體單元和楔形單元。葉片單元劃分、退化單元以及層合單元中鋪層信息定義如圖3所示。

1.2.2 復合材料風扇葉片有限元模型驗證

使用試驗與數(shù)值仿真對比的方法驗證1.2.1節(jié)中復合材料風扇葉片有限元模型的準確性。葉片裝夾及測量儀器設置如圖4所示。采用錘擊法對葉片進行模態(tài)測試，按照發(fā)動機實際運行時的葉片與輪盤的相互作用方式，使用專用振動夾具對葉片施加約束；根據(jù)葉片運行時產(chǎn)生的離心力，對夾具底部的5個頂緊螺栓各施加65 N·m的預緊力矩，并將試驗裝置固定到隔振底座上以降低外界干擾。

圖3 葉片單元劃分、退化單元以及層合單元中鋪層信息定義

使用PCB模態(tài)力錘敲擊葉片，選用金屬錘頭，敲擊點位置如圖4所示。在葉片前6階模態(tài)振型幅值較大位置選取6個測點，使用Polytec PSV-500型激光多普勒測振儀測量各測點的時域響應信號，進行FFT變換，獲得葉片前6階固有頻率。

在本次試驗測量中，各測點得到的固有頻率一致性較高，此處選取1號測點的測量結果與仿真分析結果進行對比。根據(jù)夾具與葉片之間的相互作用方式，在葉片有限元模型相應位置施加節(jié)點位移約束，進行模態(tài)分析，獲取葉片的前6階固有頻率，榫頭節(jié)點約束如圖5所示。模態(tài)分析得到的前6階固有頻率與測試結果一致性較高，測量數(shù)據(jù)與有限元數(shù)值計算結果對比，見表4。從表中可見，前3階的誤差僅為1%左右，最大誤差為7.71%。驗證了復合材料風扇葉片精細化有限元模型的準確性，模型滿足后續(xù)計算要求。

圖4 敲擊點位置

圖5 榫頭節(jié)點約束

表4 測點1處測量數(shù)據(jù)與有限元數(shù)值結果對比

1.3 復合材料風扇葉片精細化后處理方法

復合材料風扇葉片鋪層數(shù)量多，為解決現(xiàn)有有限元軟件中鋪層應力數(shù)據(jù)提取操作及鋪層失效分析結果顯示過程的繁復性，參考文獻[12]中后處理程序的實現(xiàn)邏輯與架構，編寫ANSYS Workbench軟件的二次開發(fā)ACT腳本及MATLAB程序，實現(xiàn)精細化后處理。在瞬態(tài)計算過程的時間點上，使用Workbench腳本自動化提取鋪層6個材料主軸應力（圖1），對3個積分點上的計算結果取平均值，得到鋪層上的結果?；诓牧现鬏S應力和相關失效準則，使用MATLAB程序對鋪層進行失效分析，并根據(jù)鋪層的單元節(jié)點信息，重構出鋪層有限元網(wǎng)格，生成Tecplot FEM格式數(shù)據(jù)文件，實現(xiàn)鋪層失效的可視化。

在失效準則選用方面，復合材料損傷起始和失效機理復雜，沒有一種廣泛適用的理論判據(jù)適用于復合材料所有失效工況，針對失效準則的適用性的分析不在本文研究范圍內。選用Hou等改進后的Chang-Chang失效準則，區(qū)分其中纖維失效模式為纖維拉伸失效和纖維壓縮失效，將其中一些需要大量材料試驗才能獲得的面外橫向基體開裂剪切強度、分層剪切強度等強度極限做簡化，應用于后續(xù)分析。在材料主軸坐標系下，1、2、3方向分別對應鋪層纖維方向、面內橫向和面外法向，材料主軸應力及對應強度極限如圖6所示，可得到如下失效準則：

（1）纖維失效。

當≥0時，拉伸失效是否發(fā)生的判斷準則為

當＜0時，壓縮失效是否發(fā)生的判斷準則為

（2）基體失效。

當≥0時，基體開裂是否發(fā)生的判斷準則為

當0時，基體壓潰是否發(fā)生的判斷準則為

（3）分層。

僅當≥0時，分層是否發(fā)生的判斷準則為

圖6 鋪層材料主軸應力及對應強度極限

2 中鳥撞擊工況下算例分析

2.1 載荷工況

本文研究的風扇葉片對應型號發(fā)動機的進氣道喉道面積為3.02 m，風扇轉子葉片總數(shù)=18，轉速=3600 r/min。根據(jù)適航條款CCAR33.76中的相關要求，需以不低于起飛決斷速度吸入1只1.15 kg和6只0.7 kg的鳥來驗證發(fā)動機的適航性。

圖7 風扇葉片切鳥及最大鳥切片形狀[17]

選取發(fā)動機在風扇葉片75%葉高處吸入1只1.15 kg的中鳥作為本次研究的載荷工況。使用Storace等提出的最大鳥切片載荷模型及參考文獻[18]的載荷-時間關系對葉片橢圓形區(qū)域施加均布的沖擊載荷。風扇葉片切鳥及最大鳥切片形狀如圖7所示。參考文獻[19]中的參數(shù)設置，取鳥的密度=938 kg/m，建立長徑比為2∶1的平頂圓柱模型，長度=184.14 mm，直徑=92.07 mm。沖擊載荷參數(shù)見表5。本次分析中、、方向分別為葉片的周向、展向和軸向。

表5 沖擊載荷參數(shù)

中鳥受相鄰風扇葉片切割后，撞擊在葉片壓力面上，形成橢圓形沖擊加載區(qū)域。加載區(qū)域及節(jié)點力施加方式如圖8（a）所示。根據(jù)總沖擊加載時長定義無量綱時間，沖擊力-無量綱時間的函數(shù)曲線如圖8（b）所示。在榫頭上施加了與圖5中相同的節(jié)點位移約束。

圖8 沖擊載荷施加區(qū)域及時間-載荷

2.2 復合材料風扇葉片動態(tài)響應分析

以1 ms為間隔，選取10 ms內葉片吸力面50%弦長位置曲線在-平面上的投影（如圖9所示），分析葉片的周向彎曲變形特性。圖中數(shù)據(jù)進行了無量綱化處理。從圖9（a）中可見，從初始時刻到最大加載時刻，葉片以載荷作用區(qū)域附近的局部變形為主；在后續(xù)的響應中，葉片整體彎曲變形占據(jù)主導。從圖9（b）中可見，在約10 ms內，葉片完成1次周向彎曲-回彈運動，振動頻率與葉片的1階相固有頻率相近，在彎曲過程中，氣動葉身的彎掠造型有展成平直的趨勢。

圖9 吸力面50%弦長曲線在x-y平面投影

不同時刻葉片頂部在平面投影如圖10所示。從圖中可見，葉片頂部同樣經(jīng)過約10 ms時間回到初始位置，從葉頂?shù)奈恢米兓芍?，葉片在整體上還發(fā)生軸向彎曲運動和扭轉運動；同時，葉頂?shù)男螤畎l(fā)生了變化，葉片橫截面呈現(xiàn)出“中間厚、兩邊薄”的特點，在鳥撞響應中，葉片前緣、尾緣相對于葉片弦向中部區(qū)域更加活躍。

圖10不同時刻葉片頂部在x-z平面投影

圖11 前緣、尾緣附近相對于葉片50%弦長位置的振動

在有限元模型中選取測點，分析在葉片不同高度上前緣、尾緣和葉片弦向中點（50%弦長）位置的振動，如圖11（a）、（b）所示。從圖中可見，氣動葉身頂部、中部的前、尾緣測點的位移曲線相對于相同高度的弦向中部測點的曲線有明顯波動。若以弦向中點的位移描述葉片整體位移，則在葉身上部區(qū)域的前緣、尾緣附近發(fā)生了較為顯著的局部振動，其中以測點LE_MID附近的局部振動最為明顯；在伸根段頂部SHANK測點組中，前緣、尾緣局部振動不明顯，如圖11（c）所示。由于伸根段及以下區(qū)域剛度較高，前緣尾緣與葉片中部相對厚度差較小，前緣、尾緣的活躍程度會進一步降低。提取高階模態(tài)振型，葉片第14、17、21階模態(tài)振型如圖12所示。從圖中可見，一些高階模態(tài)在葉片前緣、尾緣處活躍，而在伸根段及榫頭處振動幅值小，體現(xiàn)出與葉片響應相似的特點。因此對于復合材料風扇葉片的外物沖擊問題，葉片高階模態(tài)對響應和損傷的影響不容忽視。

圖12 葉片高階模態(tài)振型

2.3 復合材料風扇葉片鋪層失效分析

2.3.1 不同失效模式發(fā)生位置

使用第1.3節(jié)中方法提取鋪層的材料主軸應力并將鋪層失效模式可視化。由于鋪層復合材料結構層間性能較低，首先針對分層失效模式展開分析。

發(fā)生分層的位置主要受結構特點和葉片周向彎曲運動影響，集中發(fā)生在榫頭、伸根段內部的50%弦長位置附近，其他位置鋪層失效較少。葉片受到?jīng)_擊后，在結構慣性影響下，葉片的變形和運動從受載區(qū)域向下傳導，最終對榫頭及伸根段形成了拉-彎-扭組合加載模式。變厚度區(qū)域在拉伸載荷下的分層如圖13所示。榫頭、伸根段厚度變化快，其中存在較多的遞減鋪層，在受載過程中，不同高度鋪層之間的位移一致性較低，導致此處層間應力水平較高，且在榫頭、伸根段內部存在鋪層遞減，遞減鋪層的邊緣樹脂淤積形成了力學性能較差的樹脂囊，易成為分層的起點，致使伸根段、榫頭部位分層嚴重。

圖13 變厚度區(qū)域在拉伸載荷下的分層

在加載結束時刻，伸根段彎曲幅值增大，分層失效及發(fā)生位置如圖14所示。從圖中可見，內部鋪層在50%弦長位置附近出現(xiàn)肉眼不可見的分層。在=4~8 ms內，葉身彎曲幅值較大，但由于葉片的氣動葉身厚度較薄，遞減鋪層少，層間應力水平較低，因此分層未向氣動葉身發(fā)展，在=8~10 ms內，葉身中部的彎掠造型恢復，期間層間剪切應力方向突變，在=9 ms時，葉盆吸力面?zhèn)蠕亴覵S_P043_45和SS_P041_0與SS_P023_45p之間出現(xiàn)分層。

圖14 分層失效及發(fā)生位置

發(fā)生纖維失效的鋪層較少，在葉片完成1次回彈期間，受葉片的軸向彎曲運動及扭轉運動的影響，±45°鋪層在伸根段前緣位置發(fā)生纖維拉伸失效。壓力面?zhèn)蠕亴釉?0%葉高以下部位的拉伸失效指標較高；吸力面?zhèn)瓤拷~片表面的0°、-45°結構層壓縮失效指標較高，而45°鋪層的失效指標仍然指向拉伸失效模式。

在葉片彎曲過程中，伸根段、榫頭因拉伸-扭轉載荷發(fā)生基體失效。其中，壓力面?zhèn)纫恍?°鋪層在榫頭約束處、伸根段75%弦長位置，以及±45°鋪層在伸根段前緣發(fā)生了基體開裂，少量0°鋪層在伸根段前緣出現(xiàn)基體壓潰；吸力面?zhèn)仁т亴虞^少，少量0°鋪層在榫頭約束部位、少量±45°鋪層在伸根段前緣發(fā)生了基體開裂，靠近葉片表面的0°鋪層在伸根段25%葉高位置發(fā)生了基體壓潰。此外，葉片彎曲導致葉盆處的葉型曲度降低。在葉片回彈過程中，在葉盆處葉型曲度快速回復，導致吸力面?zhèn)蠕亴酉蚁蚴芾?，一些高度超過50%葉高的±45°鋪層在葉盆位置發(fā)生基體開裂，0°鋪層無明顯變化?？傮w而言，葉片基體失效的分布較為廣泛，且基體開裂較基體壓潰更加嚴重。

2.3.2 失效位置規(guī)律性及影響因素分析

對于纖維失效和基體失效，鋪層失效模式及位置具有一定的規(guī)律性，相同角度鋪層失效函數(shù)分布的相似性如圖15所示。從圖中可見，在某一時刻，中面一側大部分相同角度鋪層失效函數(shù)分布相似，僅在裁切鋪層的邊緣存在一定差異，并且由于在彎曲變形過程中靠近葉片表面的鋪層普遍承受更大的拉力，從葉片中面到葉片表面，鋪層上失效函數(shù)值整體水平逐漸升高。

不同葉片損傷模式及位置如圖16所示。從圖中可見，葉片的遠場失效更加嚴重，并且內部失效比外部失效更加嚴重，其中伸根段前緣位置的±45°鋪層發(fā)生了纖維斷裂失效，嚴重威脅結構葉片完整性，為本次分析中葉片最為薄弱的環(huán)節(jié)。

圖15 相同角度鋪層失效函數(shù)分布的相似性

圖16 葉片損傷模式及位置

在精細化后處理過程中發(fā)現(xiàn)，雖然大部分鋪層的失效函數(shù)的分布存在規(guī)律性，失效位置易于預測，在復合材料結構的沖擊損傷分析中，由于響應中高階模態(tài)的作用以及鋪層形狀、角度的復雜性，存在一些特殊情況。例如，在=9 ms時刻，鋪層PS_P017_0、PS_P002_0在葉頂中部出現(xiàn)明顯的基體開裂，與同角度鋪層的失效特性存在明顯差異。值得注意的是，在Siddens的分析中，不帶包邊復合材料風扇葉片經(jīng)歷鳥撞后，同樣于葉頂中部發(fā)生了基體開裂，但受限于葉片建模及后處理的精細化程度，未能指出發(fā)生基體開裂的具體鋪層。

針對鋪層PS_P017_0展開分析，結合圖10（b）可知此時葉片在失效位置產(chǎn)生彎折，根據(jù)式（3）包含的主軸應力，分析、、對失效的影響。PS_P017_0基體開裂失效相關應力如圖17所示。從圖中可見，對于此處失效的作用最大，其高應力區(qū)與失效位置重合，最大值達126.22 MPa，遠超鋪層的橫向拉伸強度；次之，在失效區(qū)域應力水平較高，極值的絕對值為41.48 MPa，略低于強度極限；的作用較小，僅在葉頂應力較高，最大為10.60 MPa，不足強度極限的1/3。高階模態(tài)造成的葉身上部區(qū)域的彎曲運動使PS_P017_0這一0°鋪層在失效位置產(chǎn)生了較大的橫向拉應力，進而引發(fā)了基體開裂。

圖17 PS_P017_0基體開裂失效相關應力

3 結論

（1）葉片的動態(tài)響應包括周向彎曲、軸向彎曲、扭轉，在葉身上部區(qū)域的響應中，高階模態(tài)的貢獻顯著；在失效分析中發(fā)現(xiàn)，分層受結構特點的影響大，主要發(fā)生在伸根段、榫頭內部、吸力面葉盆內部一些裁切鋪層的邊緣與相鄰鋪層之間。發(fā)生纖維失效的鋪層較少，±45°鋪層在伸根段前緣位置纖維拉伸失效顯著，為本次分析中葉片最薄弱的環(huán)節(jié)。

（2）對于纖維失效和基體失效，相同角度鋪層的失效函數(shù)值分布的相似性高，鋪層失效位置具有一定規(guī)律性；存在一些鋪層在葉片高階模態(tài)影響下于特殊位置失效，必須通過精細化分析才能得到有效預測，如在部分鋪層葉頂中部發(fā)生的基體開裂。