基于流固耦合方法的跨聲速葉片氣動和強度性能研究

汪松柏， 李紹斌,2 ， 宋西鎮

(1.北京航空航天大學 能源與動力工程學院 航空發動機氣動熱力國家級重點實驗室，北京 100191；2．先進航空發動機協同創新中心，北京 100191)

?

基于流固耦合方法的跨聲速葉片氣動和強度性能研究

汪松柏1， 李紹斌1,2， 宋西鎮1

(1.北京航空航天大學 能源與動力工程學院 航空發動機氣動熱力國家級重點實驗室，北京 100191；2．先進航空發動機協同創新中心，北京 100191)

為研究流固耦合對跨聲速壓氣機葉片氣動和強度性能的影響，采用商用軟件ANSYS-CFX/Multiphysics研究了跨聲速葉片在氣動力和離心力共同作用下的氣動性能變化規律，通過單向和雙向兩種耦合方法對比分析了葉片的強度性能。結果表明：葉片在氣動力和離心力共同作用下，葉尖前緣變形量最大，氣動性能與冷態葉片相比有明顯變化，葉尖激波位置前移，堵塞流量增大。雙向流固耦合與單向耦合相比，葉片的總體變形量增大約1.1%，最大等效應力增大0.8%左右。研究結果表明工程設計應采用雙向流固耦合方法進行葉片氣動性能設計和強度校核，以提高壓氣機的安全性和可靠性。

跨聲速壓氣機；流固耦合；氣動性能；強度性能

高性能航空發動機的發展趨勢對壓縮系統的穩定性和可靠性提出越來越高的要求。對于未來推比15一級的渦扇發動機，壓氣機單級增壓比高達2.2左右[1]，葉片級負荷越來越高，葉尖切線速度越來越大，使得高速旋轉葉片在氣動力和離心力共同作用下容易產生較大的幾何變形，葉片幾何形狀細微的變化將使得壓氣機偏離理想的工作狀態，甚至引起葉片的結構強度和疲勞壽命問題[2]。因此，研究跨聲速壓氣機葉片在氣動力和離心力共同作用下氣動性能變化規律和結構強度性能具有重要的意義，為壓氣機的氣動設計和結構優化提供參考依據。

葉片工作過程中承受氣動力、離心力和非均勻分布的溫度載荷，NASA Lewis Research Center[3]用多學科耦合分析方法研究了GE公司E3大涵道比風扇葉片在氣動力、離心力和非均勻溫度載荷作用下的幾何變形，葉片最大變形出現在葉尖前緣，氣動和離心載荷引起的風扇幾何變形量遠遠比溫度載荷大。WILSON等[4]用非線性氣動彈性模型研究了跨聲速風扇在不同工作狀態下的幾何變形對性能的影響及葉尖相對安裝角的變化。KALLESOEE等[5-6]對不同轉速下風力渦輪機葉片大變形的氣動彈性穩定性問題進行了研究，指出葉片變形后會降低其顫振極限。HOU等[7]研究表明不考慮氣動力和慣性力對葉片變形的影響會使得固有頻率偏差達-1%～5%。鄭赟等[8-9]通過流固耦合的方法對跨聲速壓氣機風扇的靜氣動彈性問題進行了研究，分析不同模態下風扇變形規律及對氣動性能影響，指出氣動力變形對氣動性能的影響在非設計轉速下大于設計轉速。李彬等[10]采用雙向流固耦合方法對透平長葉片的氣動和強度性能進行了分析，發現葉片的變形作用會改變表面氣動載荷分布和出口氣流角。姜偉等[11]利用流固耦合方法預測了汽輪機末級葉片的顫振，發現高葉高區處于不穩定區域容易發生顫振，葉片受力與位移之間的相位關系對顫振的發作有重要影響。譚禎等[12]用不同湍流模型研究了旋轉葉片氣耦合條件下的動力學特性，發現不同湍流模型對葉片結構響應影響不大。

目前，壓氣機葉片設計過程中氣動設計和強度校核是分別進行的，沒有考慮葉片在各種載荷同時作用下的變形問題。流固耦合技術的發展可以將氣動設計和強度校核結合起來，使得設計分析的結果更加接近葉片實際的工作狀態。本文采用時域的雙向流固耦合方法對NASA Rotor67跨聲速壓氣機轉子葉片進行三維數值模擬，通過固體域和流體域連續區域的離散化和構造插值函數求解各自的動力學方程，在流固耦合的交界面上實現物理信息的傳遞，分析葉片在氣動力和離心力共同作用下氣動性能的變化規律。并采用有限元程序ANSYS對葉片進行單向流固耦合研究，探討了兩種耦合方法在分析壓氣機葉片強度性能方面的差異，為預防事故發生和提高航空發動機的可靠性提供參考和分析方法。

1 數值方法

1.1 流體動力學數值求解

流固耦合計算中，由流體動力學求解提供葉片表面的非定常氣動載荷。采用流體動力學計算軟件CFX求解了雷諾平均N-S方程，其方程組形式如下：

(1)

(2)

(3)

1.2 葉片運動結構動力學模型

流固耦合計算必須在每一個物理時間步內對葉片運動方程進行求解，以得到更新的葉片外型。葉片的瞬態響應求解結構動力學方程為：

(4)

式中:[M]為質量矩陣，[C]為阻尼矩陣，[K]為剛度矩陣；X表示葉片響應的位移向量；Fa為作用在葉片表面的載荷向量。

1.3 時域雙向流固耦合數值方法

圖1 時域雙向流固耦合分析框圖Fig.1 Schematics of time domain two-way fluid-structure interaction analysis

葉片在氣動載荷和離心載荷共同作用下的幾何變形是一個典型的雙向流固耦合問題，流場的氣動力加載到葉片表面引起葉片的幾何變形，葉片的變形作為邊界條件會改變周圍流場，同時流場變化的壓力分布又影響結構的幾何。本文采用雙向瞬態的流固耦合分析方法，分別用有限體積法計算流場和有限元法計算結構場。流體域和固體域采用內外循環迭代求解的方法，其耦合求解過程如圖1所示。外循環控制流固耦合計算的時間，內循環控制流體域和固體域交界面上的插值和數據傳遞。當內循環流體域和結構域的載荷和位移都達到收斂狀態時迭代結束，進入下一個時間步的外循環，整個時間推進過程中求解動網格上的非定常N-S方程。

流固耦合交界面上，利用ANSYS多場求解器MFX(Multi-Field Solver)定義流固耦合邊界，采用守恒插值法將氣動載荷和結構邊界的位移等信息在交界面上實現傳遞。網格變形技術是實現壓氣機葉片流固耦合計算的關鍵，本文采用多塊網格變形技術，葉片的變形通過指數衰減的方式傳遞到各塊邊界上，網格較密和變形較大的區域增加網格剛度，以保證葉片變形后的網格質量較高。對于機匣處的網格，網格隨葉片變形在機匣表面上滑移，實現葉尖間隙內網格的運動。

2 數值模擬實例及分析

2.1 計算模型及網格離散

高校思政實踐課實踐主角是學生，過程性管理者(主要由指導教師和朋輩導師擔任)是配角，只有充分調動主角的實踐積極性，才能確保實踐實效性的不斷提升。

本文以NASA Lewis Research Center的跨聲速壓氣機Rotor67轉子葉片為研究對象。轉子的設計轉速為16 043 r/min，總壓比為1.63，葉片數為22，展弦比為1.58。葉片采用鈦合金制造，其材料物理性能參數為：密度ρ=4 500 kg/m3，彈性模量E=140 GPa，泊松比μ=0.34。

圖2 流體通道網格及局部放大圖Fig.2 Fluid passage mesh and detail views

圖3 固體計算網格Fig.3 Structure computational grid

圖2和圖3分別給出流固耦合計算的流體通道網格和固體網格。其中，流場計算網格采用了多塊結構化網格，環繞葉片采用正交性較好的O型網格，其他區域采用H型網格，葉尖間隙內采用蝶型網格。葉尖間隙為1.016 mm，間隙內網格層數為9，整個流體域網格總節點數約26萬，網格近壁面區域的Y+≈10，滿足湍流模型的要求。葉片結構相對簡單，有限元模型采用Solid185實體單元，單元通過8個節點來定義，每個節點有3個沿著xyz方向平移的自由度，單元具有超彈性、應力鋼化、蠕變、大變形和大應變能力，在葉片受力變形的模擬中可靠性高。整個固體網格采用六面體結構化單元，葉片的葉高、弦長和葉厚方向的單元數分別為28、36、2，整個固體網格總節點數約3 000個。

流體域計算的邊界條件如下：壓氣機進口氣動條件為標況，總溫為288.15 K，總壓為101 325 Pa，軸向進氣，計算域出口采用簡單徑向平衡方程，設置平均靜壓，轉子域給定旋轉速度。流固耦合計算過程中，轉子葉片旋轉通過單個通道的非定常物理步數設為20步，時間步長為8.5×10-6s，每個子步內流體和結構耦合迭代10次。對每一個時間步的葉尖位置位移進行監測，在一個物理時間步內監測點的最大位移隨內時間步迭代的變化小于1×10-5m時，即認為達到了變形后的平衡位置。固體域計算時葉根固定約束，壓力面、吸力面和葉尖給定流固耦合的交界面。

2.2 葉片變形對氣動性能的影響

計算過程中同時監測總參數和殘差水平，隨迭代步數的推進，總參數保持水平并且殘差降低至小于1×10-5，表明結果收斂。設計轉速下，通過改變出口背壓獲得轉子氣動特性如圖4所示，其中“Steady Simulation”表示不考慮葉片變形(冷態)的計算結果，“FSI Simulation”表示葉片在氣動和離心慣性力載荷沖擊作用變形后(熱態)的雙向流固耦合計算結果，“Experiment”為NASA公開發表的實驗結果數據[13]。從圖中可以看出，數值結果和實驗結果曲線趨勢一致，葉片變形后的特性線相對于冷態結果整體向右偏移。近失速點工況冷態葉型的總壓比為1.67，葉片變形后的總壓比提升約1.8%，與實驗總壓比結果吻合很好。最高效率點處，不考慮變形的最大效率為90.24%，葉片幾何變形后的最大效率為90.25%，兩個數值結果都較實驗最大效率偏低2.9%左右，這是由于數值湍流黏性偏大引起的，對效率計算值的偏差量在工程應用中是可接受的。實驗獲得的堵塞流量為34.96 kg/s，冷態葉型計算的堵塞流量比實驗值偏低0.92%，ZHENG等[14]計算Rotor67冷態葉型的堵塞流量比實驗值偏低1.2%。葉片變形后的堵塞流量為35.3 kg/s，比實驗結果增大0.97%，葉片變形主要改變了上半葉高前緣的進口葉型角，使通道的流通能力顯著增強。

圖4 葉片的氣動特性Fig.4 Aerodynamic characteristics of blade

壓氣機一般工作在高效率點附近，因此本文針對最大效率工況進行分析。圖5給出最大效率工況90%和70%葉高截面變形前后的靜壓分布對比，其中“Steady computation” 表示冷態的數值結果，“Static deflection”表示熱態的雙向流固耦合數值結果。從圖中可以看出，葉片變形明顯改變了氣動載荷分布，葉尖前緣氣動載荷的變化最大。葉片變形后使得前緣壓力面壓力升高，吸力面的壓力降低，葉尖前緣的進氣攻角增大。從90%葉高截面壓力分布可知，激波位于葉片吸力面約65%弦長位置，葉片變形后明顯引起激波位置前移，激波強度變化不大。70%葉高處葉片的變形量相對減小，葉片幾何變形引起的壓力變化主要集中前半弦長，對激波位置和激波強度的影響較小。從70%和90%葉高截面變形前后氣動載荷的對比可以看出，葉片的幾何變形對氣動性能的影響隨葉高的增加而增大，氣動性能的改變主要集中在高葉高區域。

圖5 最大效率點不同葉高處的靜壓對比Fig.5 Comparison of the pressure distribution at various span locations at peak efficiency point

2.3 葉片變形對激波的影響

對于跨聲速壓氣機，葉尖相對馬赫數較高，葉片通道經過激波的強壓縮產生較高的增壓比。同時，激波與葉片附面層的干擾作用產生較大的流動損失。因此，跨聲速壓氣機中激波結構直接影響壓氣機的氣動性能。

圖6為最高效率工況不同葉高截面的相對馬赫數云圖對比，從熱態計算結果、冷態計算結果和NASA的實驗結果中通道相對馬赫數分布的對比來看，變形前后相同葉高截面的波系結構基本一致，變形前的通道激波傾斜更大,變形后的激波位置明顯前移，與圖5一致。葉片變形后，激波前的相對馬赫數增大，通道相對馬赫數分布更接近實驗結果。90%葉高處通道進口段的相對馬赫數實驗結果大于1.35，葉片變形后的吸力面一側的高馬赫數區域明顯增大，與實驗結果吻合更好。與90%葉高的相對馬赫數對比，70%葉高截面變形前后通道的相對馬赫數變化較小。在葉片尾跡區域，數值計算得到了較細致的尾跡虧損分布，而實驗的穩態測量未給出詳細分布。

3 強度和模態計算結果及分析

3.1 各載荷作用下葉片變形規律

高速旋轉的葉片承受氣動力和離心力，單向流固耦合僅考慮載荷對葉片結構的影響，忽略結構幾何變化對流場的作用。雙向流固耦合方法考慮了流體和固體相互耦合作用過程，各載荷作用下葉片的變形如圖7所示。從圖中可以看出，氣動載荷和離心載荷使得葉片產生彎曲變形，葉片變形量隨葉高增加而增大，最大變形量出現在葉尖前緣，離心載荷造成的最大變形幅值是氣動載荷的2.75倍。葉片在氣動載荷和離心載荷共同作用下，雙向流固耦合計算的總體變形幅值比單向耦合大1.1%，這種差異反映了葉片受載荷作用變形后反過來對流場的影響，變形作用使得流場作用于葉片的氣動載荷增強，體現了流固之間相互影響的耦合作用。對于壓氣機中廣泛存在的非線性行為，需要同時考慮結構和流場之間的相互作用才能較準確的反映這種非線性行為。

圖6 最高效率工況不同葉高處相對馬赫數分布Fig.6 Relative Mach number distribution at various span locations at peak efficiency point

圖7 各載荷作用下的葉片變形量(m)Fig.7 Contour of blade deflection under various forces(m)

氣動力和離心力共同作用下葉片的總體變形量及三個方向的變形分量如圖8所示，表1為各載荷作用下總體變形量和三個分量最大值對比的結果， 其值均為絕對值。壓氣機葉片從葉根到葉尖彎扭積疊使得各截面重心不重合，在離心載荷的作用下，葉片并不是單一的徑向變形。從表1可以看出，周向方向的變形量最大，軸向其次，徑向方向變形量最小，氣動力引起的變形占總變形量的26%～30%。葉片由于徑向變形，使得工作狀態的葉尖間隙和設計值有所不同，表1中兩種流固耦合方法計算的徑向變形量相差0.6%。壓氣機工作過程中徑向變形使得葉尖可能與機匣發生不均勻碰磨，導致葉片疲勞出現葉尖掉角等故障[15]。實際工程設計階段較準確的預估葉片變形對壓氣機的安全性和可靠性具有重要的意義。葉片在氣動力和離心力共同作用下的周向變形量達2.473 mm，熱態葉片出現扭轉恢復現象，這種扭轉恢復隨著葉高的增加而增大，使得葉尖區域的扭轉角偏離設計值，進而導致壓氣機設計點的氣動性能發生變化。

表1 各載荷作用下葉片變形量最大值

圖8 葉片總體變形量及三個變形分量(m)Fig.8 Total and three component deformation of blade(m)

3.2 強度計算結果及分析

葉片在各種載荷作用下的等效應力與應變最大值見表2，圖9為各載荷作用下葉片的等效應力云圖。從圖9可以看出，氣動力作用葉片后產生的最大等效應力出現在吸力面中部下半葉高區域，氣動力和離心力共同作用下的等效應力最大值產生于壓力面靠近前緣葉根區。從表2最大等效應力和應變值看出，葉片的氣動力加劇了葉片的變形，氣動力的作用使得葉片最大等效應力增大約5.3%，可見氣動力對結構強度的影響是不能忽略的。氣動載荷和離心載荷共同作用下的最大應力約為431.3 MPa，比單向耦合的最大應力大0.8%，其最大應力值遠低于鈦合金的屈服極限。

表2 不同載荷下葉片等效應力與應變最大值

圖9 各載荷作用下的葉片等效應力(Pa)Fig.9 Von-Mises stress contour of blade under various forces(Pa)

3.3 頻率計算結果及分析

各載荷作用下葉片的前5階頻率以及無載荷下葉片的固有頻率見表3。由表3可知，由于離心力的剛化作用，葉片的1階固有頻率相對無載荷時增大約30%，而高階頻率的增加比例相對1階頻率減小。氣動力對葉片的固有頻率影響較小，1階固有頻率與無載荷時一樣，高階固有頻率的變化范圍在-0.06%～0.075%之間，葉片固有頻率基本不受氣動載荷影響。

表3 各載荷作用下葉片的頻率

對于旋轉機械，如果葉片的固有頻率與系統在旋轉過程中產生的激振頻率相同或者接近，將會引起葉片共振，產生較大的振幅和振動應力，容易造成高周或者低周疲勞斷裂。本文研究的壓氣機葉片旋轉過程產生的激振頻率為5 882.4，11 764.8 Hz，…。從表3葉片的固有頻率看出，葉片的固有頻率比激振頻率小很多，能有效的避開共振現象的發生。

4 結 論

采用時域的雙向流固耦合方法數值研究了NASA Rotor 67跨聲速轉子葉片在氣動載荷和離心載荷共同作用下的氣動性能變化規律，對兩種流固耦合方法計算獲得的葉片強度性能結果進行了對比分析，其主要結論如下：

(1) 跨聲速葉片在氣動力和離心力共同作用下氣動性能變化明顯，葉片變形后導致激波位置前移，堵塞流量增大，流場參數變化主要集中在上半葉高區，實際工程設計中應該考慮氣動載荷和離心載荷共同作用后對葉片變形及氣動性能的影響。

(2) 基于時域的雙向流固耦合能更加真實反映葉片和流場之間相互作用關系，比僅考慮流場對結構影響更準確。通過該分析方法，可以獲得葉片最大應力、應變位置和變形后的結構幾何，分析葉片總體和各個方向的變形分量大小，為葉片氣動和結構設計提供參考依據和方法。

[1] 王占學，劉增文，蔡元虎，等. 推重比15一級發動機關鍵技術及分析[J]. 航空發動機，2010, 36(1): 58-62.

WANG Zhanxue, LIU Zengwen, CAI Yuanhu, et al. Key technologies and analysis of aeroengine with thrust to weight ratio up to level of 15[J]. Aeroengine, 2010, 36(1): 58-62.

[2] 陳懋章，劉寶杰. 風扇/壓氣機氣動設計技術發展趨勢——用于大型客機的大涵道比渦扇發動機[J]. 航空動力學報，2008, 23(6): 961-975.

CHEN Maozhang, LIU Baojie. Fan/compressor aero design trend and challenge on the development of high bypass ratio turbofan[J]. Journal of Aerospace Power, 2008, 23(6): 961-975.

[3] MAHAJAN A J, STEFKO G L. An iterative multidisciplinary analysis for rotor blade shape determination[R]. AIAA/SAE/AS-ME/ASEE 29th Joint Propulsion Conference and Exhibit 1993.

[4] WILSON M J, IMREGUN M, SAYMA A I. The effect of stagger variability in gas turbine fan assemblies[J]. Journal of Turbomachinery, 2007, 129(2): 404-411.

[5] KALLES?E B S, HANSEN M H. Some effects of large blade deflections on aeroelastic stability[C]// 47th AIAA Aerospace Sciences Meeting Including the New Horizons Forum and Aerospace Exposition. Orlando, FL, 2009: 5-8.

[6] KALLES?E B S. Equations of motion for a rotor blade, including gravity, pitch action and rotor speed variations[J]. Wind Energy, 2007, 10(3): 209-230.

[7] HOU J, WICKS B J. Root flexibility and untwist effects on vibration characteristics of a gas turbine blade[R]. Defence Science and Technology Organization Victoria (Australia) Platform Sciences Lab, 2002.

[8] 鄭赟，田曉，楊慧. 跨聲速風扇葉片變形對氣動性能的影響[J]. 航空動力學報，2011, 26(7): 1621-1627.

ZHENG Yun, TIAN Xiao, YANG Hui. Impact of blade deflection on aerodynamic performance of transonic fan[J]. Journal of Aerospace Power, 2010, 26(7): 1621-1627.

[9] 鄭赟，王彪，楊慧. 跨聲速風扇葉片的靜態氣動彈性問題[J].航空動力學報，2013, 28(11): 2475-2482.

ZHENG Yun, WANG Biao, YANG Hui. Static aeroelastic problems of transonic fan blade[J]. Journal of Aerospace Power, 2013, 28(11): 2475-2482.

[10] 李彬，宋立明，李軍. 基于雙向流固耦合的長葉片氣動和強度性能的數值研究[J]. 推進技術, 2014, 35(2): 202-207.

LI Bin, SONG Liming, LI Jun. Numerical investigations on aerodynamic and mechanical performance of long blade based on two-way fluid-structure coupling approach[J]. Journal of Propulsion Technology, 2014,35(2): 202-207.

[11] 姜偉，謝誕梅，陳暢，等. 基于時域分析法的汽輪機末級葉片顫振預測及分析[J]. 振動與沖擊，2015, 34(11): 194-199.

JIANG Wei, XIE Danmei, CHEN Chang, et al. Flutter prediction and analysis for a steam turbine last-stage blade based on time domain analysis method[J]. Journal of Vibration and Shock, 2015, 34(11): 194-199.

[12] 譚禎，李朝峰，太興宇，等. 不同湍流模型旋轉葉片氣固耦合動力學特性研究[J]. 振動與沖擊, 2013, 32(11): 46-50.

TAN Zhen, LI Chaofeng, TAI Xingyu, et al. Dynamic characteristics of rotating blades with gas-structure interaction of various turbulence models[J]. Journal of Vibration and Shock, 2013, 32(11): 46-50.

[13] STRAZISAR A J, WOOD J R, HATHAWAY M D, et al. Laser anemometer measurements in a transonic axial-flow fan rotor[R]. NASA Technical Paper 2879, 1989.

[14] ZHENG R, XIANG J, SUN J. Blade geometry optimization for axial flow compressor[C]//ASME Turbo Expo 2010: Power for Land, Sea, and Air. Glasgow, UK: American Society of Mechanical Engineers, 2010: 633-644.

[15] 范順昌，唐曉輝，張銀東，等. 航空發動機高壓壓氣機三級轉子葉片掉角分析[J]. 失效分析與預防，2014, 9(2): 110-114.

FAN Shunchang, TANG Xiaohui, ZHANG Yindong, et al. Failure analysis of third-stage rotor blade of high-pressure compressor in aero-engine[J]. Failure Analysis and Prevention, 2014, 9(2): 110-114.

Investigations on aerodynamic and mechanical performance of a transonic blade based on a fluid-structure interaction method

WANG Songbai1， LI Shaobin1,2, SONG Xizhen1

(1. National Key Laboratory of Science and Technology on Aero-thermodynamics, School of Energy and Power Engineering, Beijing University of Aeronautics and Astronautics, Beijing 100191, China;2. Collaborative Innovation Center for Advanced Aero-Engine, Beijing 100191, China)

In order to study the influence of the fluid-structure interaction for the aerodynamic and mechanical performance of a transonic compressor blade, a time domain two-way fluid-structure interaction numerical method with ANSYS-CFX/Multiphysics was applied to examine the aerodynamic performance of the transonic blade under aerodynamic and centrifugal forces. The mechanical performance was analyzed by comparing the results between one-way and two-way fluid-structure interaction methods. The results show that the maximum displacement of the blade appears at the tip near leading edge, and that the blade deformation has significant effects on the aerodynamic performance in comparison with the cold blade, which makes passage shock wave position forward and choke mass flow rate increase. Comparing the results between one-way and two-way FSI numerical computations in mechanical performance, the amplitude of the overall blade deformation increases by 1.1% and the maximum Von-Mises stress increases by 0.8%. The overall results indicate that the industrial practice should adopt a two-way fluid-structure interaction method to guide the aerodynamic design and check the strength, and to improve the safety and reliability of the compressor.

transonic compressor; fluid-structure interaction (FSI); aerodynamic performance; mechanical performance

國家自然科學基金項目(51106003)

2015-06-02 修改稿收到日期：2015-09-14

汪松柏 男，碩士，1989年生

李紹斌 男，博士，講師，1980年生

E-mail: lee_shaobin@buaa.edu.cn

V231.3

A

10.13465/j.cnki.jvs.2016.20.017