航空發動機渦輪流動及噪聲數值模擬

嚴 佳, 柴軍生, 任國哲, 吳忱韓, 唐小龍, 楊小權,丁 玨, 翁培奮

(1.上海大學力學與工程科學學院, 上海 200444;2.中國航空發動機集團沈陽航空發動機研究所, 遼寧沈陽 110066;3.上海電力大學能源與機械工程學院, 上海 200090)

渦輪作為航空發動機最核心的熱端部件之一, 其工作環境具有典型的高溫、高壓、高載荷和強腐蝕性等特點.20 世紀90 年代研制成功的推重比10一級發動機的渦輪前溫度已經達到了1 850～1 950 K.在美國開展的“綜合高性能渦輪發動機技術” 和英國開展的“先進核心機軍用發動機”研究計劃中, 發動機渦輪進口溫度更將達到2 200 K 以上.極高的渦輪前溫度對渦輪葉片和渦輪盤的耐熱性是極大的考驗, 常規材料早已無法滿足該耐熱需求, 這對渦輪冷卻技術的發展提出了更高的要求.渦輪旋轉盤腔及旋轉葉片流動換熱數值模擬研究為渦輪冷卻技術的發展提供了強有力的支撐.

對于民用航空, 噪聲水平已經成為推動發動機研制的關鍵性指標, 為此國際民航組織和各國民用航空管理局均制定了嚴苛的民用航空噪聲適航條例.隨著航空發動機噴流、風扇低噪聲設計和控制技術的發展, 渦輪噪聲等曾經不太重要的噪聲分量逐漸變得重要.航空發動機渦輪內部流動結構、發聲機理及噪聲傳播特性十分復雜, 聲學數值模擬逐漸在該領域發揮極大的作用.

為推動渦輪流動、傳熱和噪聲數值模擬方法的發展和工程應用, 本工作對其研究現狀進行了回顧和展望.

1 渦輪流動特性

渦輪作為航空發動機的核心部件, 主要由盤腔和葉片等部分組成, 其可靠性和效率的高低決定了航空發動機的整體性能.渦輪效率與渦輪進口溫度密切相關, 在發動機尺寸不變的情況下, 渦輪進氣口溫度每提高55 K, 發動機推力可提高約10%.然而, 高溫、高壓和強腐蝕性工作環境等對渦輪各部件的耐熱及抗腐蝕性等能力提出了嚴峻的挑戰.渦輪中旋轉盤腔和旋轉葉片等轉動部件在承受高溫高壓燃氣沖擊的情況下, 還受到強大的離心力作用, 其工作環境更為惡劣.為保障渦輪中轉動部件可靠高效地運行, 研究者們開展了大量研究, 結果表明渦輪中流動結構異常復雜, 具體表現為: ①旋轉盤腔和葉片幾何非常復雜;②轉子/靜子相對運動速度高且存在強烈干涉;③冷氣與燃氣存在大密度比和溫度比;④輪轂及機匣附近存在復雜的二次流動;⑤葉片吸力面和葉尖間隙中存在強烈的激波/邊界層干擾;⑥旋轉盤腔和葉片均存在冷熱交替的復雜湍流流動等.早期為解決渦輪復雜流動環境下耐熱極限這一難題, 人們在研制新型耐高溫材料時, 大多通過實驗研究發動機空氣系統對受高溫部件的冷卻作用.隨著計算機技術的發展, 數值模擬在發動機相關研究領域被廣泛應用, 但受限于渦輪復雜流動結構, 準確模擬渦輪流動結構和換熱特性依舊是當前面臨的關鍵性技術難題.渦輪系統中盤腔和葉片是一個整體, 但其流動和傳熱特征又有較大差異, 因此將其分為兩個部分進行討論.

1.1 旋轉盤腔

航空發動機渦輪盤腔根據轉盤的配置不同可以分為轉-靜和轉-轉兩種盤腔系統.轉-靜盤腔系統作為航空發動機渦輪核心部件, 一直備受關注.圖1 為某型號航空發動機渦輪轉-靜盤腔簡化模型.經過大量簡化, 渦輪盤腔的結構依然非常復雜, 而幾何外形的復雜直接導致對腔體內部流動換熱分析難度極大.如何探明盤腔內部流動換熱特性一直是航空發動機渦輪研制中的熱點問題.針對旋轉盤腔流動的研究始于對自由轉盤流動的探索, 基于早期自由轉盤流動的研究結論[1-2], Daily 等[3]通過理論分析和實驗驗證提出了封閉轉-靜盤腔中存在的4 種基本流態, 分別為: 轉子邊界層和靜子邊界層相匯合和不匯合的層流;轉子邊界層和靜子邊界層相匯合和不匯合的湍流.這4 種基本流態為后續旋轉盤腔流動換熱研究奠定了基礎.隨著研究的深入, 理論分析難以滿足復雜腔內流動換熱預測和分析的需求, 而隨著計算機技術的發展, 數值模擬逐漸成為旋轉盤腔內部流動與對流換熱的主要研究手段之一.

圖1 某型號航空發動機渦輪盤腔簡化模型Fig.1 A simplified model of the rotating cavity in an aero-engine

1985 年, Sambo[4]通過數值求解層流的Navier-Stokes (NS)方程, 分析了旋轉雷諾數、軸向間隙系數、流量系數、盤腔間隙系數及密封冠型式對流動的影響, 預測得到的盤腔切向流動速度與實驗結果吻合較好.Long 等[5]利用有限體積方法和多重網格技術, 對軸向進氣轉-靜盤腔內層流換熱進行了詳細的研究, 討論了不同表面溫度分布對局部努塞爾數的影響, 并給出了相對切向速度和徑向速度的分布曲線, 計算結果與以往實驗測量值相近.旋轉盤腔層流模擬對于簡單構型腔體內流的計算結果與實驗值吻合得較好, 但不適用于存在復雜湍流流動的真實盤腔流動.之后, 湍流模擬成為了旋轉盤腔數值模擬的主流方向.湍流模擬的3 種基本方法為雷諾平均NS 方程(Reynolds-averaged Navier-stokes equations, RANS)方法、大渦模擬(large eddy simulation, LES)方法和直接數值模擬(direct numerical simulation, DNS)方法.RANS方法通過數值求解雷諾平均NS 方程, 依靠理論與經驗的結合, 對平均運動中湍流脈動量引進一系列模型假設, 建立湍流模型, 實現旋轉盤腔流動換熱的數值模擬.由于湍流流動中雷諾應力的主要貢獻來自大尺度脈動, 而大尺度脈動的特性與邊界條件密切相關, 因此雷諾應力的封閉模型難以獲得普適性, 目前尚不存在所有復雜流動都適用的統一封閉模型.

Koosinlin 等[6]通過數值求解普朗特混合長度模型封閉的RANS 方程, 給出了自由旋轉圓盤附近湍流的預測, 結果表明: 在沒有渦流的情況下, 混合長度湍流模型計算結果與實驗數據吻合較好;但在高渦流率下, 由于湍流的各向異性, 計算結果與實驗產生了較大的誤差.Lapworth 等[7]同樣采用混合長度湍流模型, 模擬了施加徑向氣流的非平面盤腔腔內流動和換熱, 并與相應實驗數據進行比較, 得到了較滿意的結果.徐國強等[8]和王美麗等[9]采用SIMPLE 算法和混合長度湍流模型, 模擬了出入流轉-靜盤腔的流場和溫度場.為了抓住腔內氣流的基本規律, 將盤腔內三維非定常可壓縮湍流流動簡化為軸對稱的定常不可壓湍流流動,計算網格采用結構化網格, 并對局部流動變化大的地方進行網格加密, 分析了其盤面平均努塞爾數隨旋轉雷諾數和流量系數的變化, 計算結果與實驗值符合較好.

Chew[10]采用k-ε湍流模型, 給出了同向旋轉盤間徑向出流、靜止盤與旋轉盤間徑向出流以及自由盤流動的數值計算結果, 預測結果顯示傳統的k-ε湍流模型可能不適合腔體內部流動的模擬.Morse[11]通過修正低湍流雷諾數k-ε湍流模型, 對具有徑向出氣口的轉-靜旋轉空腔內的流動進行了預測, 相比于Chew[10]得到了較好的計算結果, 但該方法對復雜流動并不適用.Lacovides 等[12]采用k-ε湍流模型和代數模型(algebraic stress model, ASM)研究了軸對稱不可壓縮等溫流動通過旋轉空腔的數值計算中的湍流建模問題, 結果顯示不同湍流模型對不同流動預測結果差異較大.張靖周等[13]采用柱坐標系結構化網絡, 對具有徑向出流的轉-靜盤腔腔內流動與換熱進行數值模擬, 在盤腔核心區采用標準k-ε模型, 在近壁區采用Launder-Sharma 低雷諾數模型, 討論了旋轉雷諾數和無量綱流動速率對盤腔內流動與換熱的影響, 計算結果與實驗結果相吻合, 而與近壁區采用單方程模型的計算結果相比較, 發現近壁區采用單方程對腔內流動特征模擬的精度要優于兩方程模型.白洛林等[14-15]以羅羅公司未來某型雙軸發動機的高壓渦輪旋轉盤腔為計算模型, 就湍流參數對復雜形狀渦輪盤腔流場的影響進行數值研究, 入口和出口處的面積采用等效面積處理方法, 入口采用質量流量入口, 多個出口間采用分配質量流量的出口條件, 計算采用SA 湍流模型, 為滿足湍流模型對網格的要求,計算網格在壁面處在轉盤的壁面添加20 層邊界層網格, 以保證y+在1 左右.結果表明湍流參數相同的流動, 其流體結構是非常相似的, 并對有預旋進氣的旋轉盤腔腔內湍流流場進行數值模擬, 比較了SA、k-ε、k-ω3 種湍流模型的計算結果, 結果表明SA 湍流模型與實驗結果符合較好.張晶輝等[16]對渦輪盤腔徑向封嚴機理進行數值模擬, 研究了渦輪盤腔徑向封嚴結構內的非定常流動特性.為保證計算的準確性, 計算網格在主流壁面y+<5, 盤腔壁面y+<2,計算采用剪切壓力傳輸(shear stress transport, SST)湍流模型, 主流和盤腔進口設置為理想空氣, 盤腔進口標記為封嚴氣體, 以兩股氣流在封嚴處摻混組分濃度計算封嚴效率, 結果表明:徑向封嚴結構的輪緣搭接邊阻擋了封嚴氣體的出流, 在盤腔內形成了自循環的渦核結構, 封嚴氣流和入侵氣流的相互影響區域局限在封嚴間隙內, 對盤腔內影響較小, 盤腔內渦核結構較為穩定, 速度非定常波動較小.丁水汀等[17]利用自由盤驗證了SSTk-ω湍流模型在模擬旋轉盤腔換熱與流動特性方面的適用性, 結果表明該方法計算結果與分析解、實驗值均吻合較好, 并在此基礎上探究了中心進氣轉靜系盤腔換熱特性對無量綱參數的敏感性, 計算采用結構化網格, 并進行網格無關性驗證, 在入口和出口處為防止產生回流而設置延伸段, 計算結果表明, 在某些工況內, 反映進口耗散效應的?？颂販蕜t與反映轉盤熱邊界的基比切夫準則以及固體的無量綱導熱系數對努塞爾數的影響程度與旋轉雷諾準則量級相同.Liao 等[18]就預旋噴嘴的幾何參數對預旋效率和總壓損失的影響進行了數值研究, 采用SST湍流模型, 選取8種不同壓比和4種轉速, 從預旋效率、總壓損失、預旋噴嘴和接收孔流量系數、噴嘴出口速度和氣流角的不均勻性等研究了預旋系統的流動動力學特性(見圖2).

圖2 四種不同轉速下中截面上的總壓和流線分布[18]Fig.2 Total pressure and streamline distributions on the meridian plane at four different rotating speeds[18]

相比于RANS 方法, LES 方法對于可解尺度湍流運動采用數值計算方法直接求解, 而對于小尺度湍流脈動對大尺度湍流運動的作用, 建立亞格子(subgrid scale, SGS)模型, 能夠更好地解析小尺度的脈動對流動的影響, 從而分辨出更多的流動細節.Andersson 等[19]對開放轉-靜盤之間的流動進行了大渦模擬, 研究了間隙比對腔內流動結構的影響, 發現雷諾剪切應力矢量和平均速度梯度矢量在垂直于旋轉軸的平面上趨于一致, 力矩系數在轉子處與經驗相關公式相比有很好的一致性, 在靜子側則相對較低.Gao 等[20]利用大渦模擬方法對轉-靜盤腔封嚴處非定常流動進行計算, 并與RANS 和非定常雷諾平均NS 方程(unsteady Rynolds averaged Navier-Stokes equations, URANS)計算結果對比, 結果顯示只有LES 在非定常計算中表現出與實驗相似的明顯峰值頻率、流動分布和動態特性(見圖3 和4).

圖3 盤腔封嚴處的平均速度分布[20]Fig.3 Mean velocity profiles within the seal gap[20]

圖4 盤腔封嚴處周向渦量和二維速度矢量分布[20]Fig.4 Circumferential vorticity component contours with 2-D velocity vectors on meridian planes[20]

Tuliszka-Sznitko 等[21]利用大渦模擬技術對轉盤加熱的轉-靜盤腔腔內湍流流動進行數值模擬, 局部努賽爾數的分布反映了腔體內部流動結構(見圖5), 結果顯示在湍流程度最高的外沿附近, 傳熱迅速增強.之后Tuliszka-Sznitko 等[22]分析了轉-靜盤腔內雷諾應力張量分量、湍流熱流張量分量、努塞爾數和湍流普朗特數等結構參數的分布, 為建模提供了依據, 計算結果與實驗數據符合較好.Lygren 等[23]利用動態和混合動態亞格子模型對旋轉盤和靜止盤之間的流動進行了大渦模擬, 并將計算結果與直接數值模擬結果進行了比較, 結果表明混合動態模型比動態模型給出了更好的整體預測效果.研究表明, LES 已基本滿足旋轉盤腔流動預測的需求, 為更進一步準確捕捉旋轉盤腔中流動結構, 尤其是流動轉捩等, DNS 方法也逐漸被用于簡單旋轉盤腔問題的研究[24-25].

圖5 無量綱半徑下局部努賽爾數分布[21]Fig.5 Distributions of the local Nusselt numbers[21]

對于簡單構型旋轉盤腔內流動結構, 研究人員進行了十分全面的研究和分析, 數值模擬中對盤腔結構和邊界條件進行了大量的簡化, 且計算多采用結構化網格.而實際工程中旋轉盤腔結構非常復雜, 數值模擬中結構化網格的使用受到較大的限制, 同時對盤腔構型的邊界條件的簡化使得數值模擬結果只能從定性上分析腔體內部流動結構和換熱特性, 對實際盤腔涉及跨尺度旋流及流熱耦合問題的精細化研究還不夠.目前, 針對真實構型的旋轉盤腔流動換熱相關數值模擬研究相對較少.在數值模擬方法方面, RANS 和LES 是旋轉盤腔流動與換熱的主要研究手段.就RANS 而言, 如何得到一種普適性好、精度高的湍流模型是其難點之一.現有的湍流模型中, 普朗特混合長度模型僅適用于局部平衡狀態的湍流, 無法模擬對流和擴散作用的影響;標準k-ε模型中假設湍流黏度各項同性, 但旋轉盤腔內部強旋流中湍流黏度表現出明顯的各項異性;改進的k-ε模型對旋流的預測精度有所提高, 但對跨雷諾數量級的湍流模擬難度依舊較大, 而旋轉盤腔中切向和主流方向的雷諾數差別較大.由此可見, 目前湍流模型在旋轉盤腔數值模擬中還存在一定的局限性.而相較于RANS, LES 通過模化小尺度渦結構的形式, 能夠得到更加精確的計算結果和更加豐富的流場信息, 這也就要求更小的網格尺寸和更高的網格質量, 從而增加了數值模擬計算成本.如何較好地平衡計算精度和計算成本之間的關系, 是真實盤腔流動與換熱面臨的一大難點, 而隨著計算機技術的發展, 目前LES 方法已經表現出了較強的工程適用性, 在未來復雜流動換熱問題數值模擬方向具有較好的應用前景, 而DNS 應用于旋轉盤腔的數值模擬仍受限于其巨大的計算成本.

1.2 旋轉葉片

旋轉葉片作為渦輪主要熱端部件之一, 一直是發動機相關方向的研究熱點.隨著數值計算能力的發展, 渦輪葉片研究對象已由早期簡化的平面葉柵向全尺寸全環旋轉葉片轉變, 研究方向從原先的主要考慮葉片氣動性能發展至同時考慮氣動、傳熱和噪聲等多種性能的氣-熱-聲多場耦合分析.這促進了旋轉葉片冷卻技術的迅速發展, 1960 年以前不采用冷卻技術, 渦輪入口最高溫度不超過1 300 K, 發展至今冷卻量達到500～600 K, 渦輪前燃氣溫度達到1 950 K,其中旋轉葉片流動換熱數值模擬方法的發展和由此推動的流動換熱機理研究起了關鍵性作用.

Boyle[26]通過求解準靜態三維NS 方程, 湍流模型選取修正的BL 模型, 模擬7 種不同外形葉片的表面流動換熱特性, 并與實驗結果相比較, 同時研究了近壁面對計算結果的影響.Boyle 等[27]采用多步龍格庫塔方法和隱式殘值光順, 湍流模型為混合長度模型, 對葉片表面流動與換熱特性進行定常求解, 研究了旋轉、葉頂間隙和進口邊界層厚度變化對葉片和端壁傳熱的影響.Boyle 等[28]研究網格正交性對葉片表面流動與換熱的影響, 計算采用5種不同的C型網格, 計算結果顯示, 網格的不同會影響葉片尾跡區壓力分布, 而對葉片表面換熱特性影響不大.此外, Ameri 等[29]通過求解三維RANS 方程, 并在BL 湍流模型中加入Mayle 提出的過渡模型, 結果顯示, 加入過渡模型后的計算結果與實驗結果符合更好.Tutar 等[30]使用基于Smagorinsky-Lilly 亞格子模型的大渦模擬方法對葉片表面流動與換熱進行計算, 結果表明相比于求解RANS 方程, LES 能夠在預測渦輪葉柵前氣熱流動特性方面提供更準確、可靠的計算結果(見圖6).

圖6 計算模型和不同出口雷諾數對葉片表面換熱系數的影響[30]Fig.6 Model and the influence of different outlet Reynolds number on blade surface heat transfer coefficient[30]

Bhaskaran 等[31]采用高階數紙法對跨音速葉片葉柵進行了大渦模擬, 研究了入流湍流作用下葉片表面的傳熱和邊界層的發展, 并與Arts 等的實驗進行定量比較, 二者表現出極好的一致性.王鵬等[32]通過對影響流熱耦合計算精度的若干因素的敏感性分析, 得出一套合理的參數設置, 并對三維MarkⅡ葉型進行數值模擬研究, 探討了雷諾數對渦輪葉片表面換熱的影響.李虹楊等[33]提出對γ ?Reθ轉捩模型的粗糙度修正方法, 并對MarkⅡ渦輪導葉進行數值模擬.為驗證數值算法的準確性, 首先對粗糙平板和某高壓渦輪導葉(high-pressure turbine vane, HPTV)進行數值模擬, 結果顯示模擬預測的努賽爾數分布與實驗值吻合較好, 同時轉捩位置的預測也與實驗值接近.MarkⅡ導葉數值模擬網格壁面處y+保證在0.4～0.9, 邊界條件為標準實驗5441 工況, 計算分析導葉表面粗糙度對渦輪葉片流動轉捩以及傳熱特性的影響,預測結果顯示修正后的γ-Reθ模型預測的導葉表面溫度、熱流量、轉捩位置等參數的變化趨勢與真實流動情況相符.江立軍等[34]采用尺度自適應模擬方法對燃燒室及渦輪進行了數值模擬, 計算網格經網格無關性驗證后, 燃燒室采用六面體和四面體混合網格, 網格數約為1 570萬, 燃氣渦輪均采用六面體網格, 其中燃氣渦輪導葉網格數為188 萬, 動葉網格數為633 萬.燃燒室和渦輪的交界面采用凍結轉子方法進行處理, 保證計算域間聯通及數據傳遞, 燃氣渦輪導葉與動葉采用型轉換方法進行數據傳遞.計算結果表明, 尺度自適應模擬方法具有較高的計算精, 能夠精準捕捉燃燒室及渦輪流道中復雜的流動渦系(見圖7).此外, 吳忱韓等通過大渦模擬研究了內冷渦輪葉片的流動換熱機理(見圖8), 并提出了波紋內冷通道等高性能內冷換熱策略.張國花等[35]在給定通道雷諾數的條件下, 通過實驗和數值計算相結合的方式, 研究了矩形內冷通道中截斷肋片在不同排布方式下的換熱特性, 并基于流動特征深入分析了其中的對流換熱機理.

圖7 燃燒室渦輪交互作用尺度自適應模擬[34]Fig.7 Vortex structure in the combustor and turbine[34]

圖8 帶氣冷的渦輪葉片流動換熱數值模擬Fig.8 Simulation of flow and heat transfer in turbine with cooling pipes

以上數值研究為旋轉葉片冷卻結構設計提供了方向.隨著旋轉葉片冷卻技術的發展, 葉片表面流動換熱特性發生了巨大的變化, 也提高了相關研究的難度.氣膜冷卻是渦輪葉片重要的冷卻技術之一, 最早出現在機翼的防凍問題上, 隨后被應用于葉片的冷卻保護, 氣膜冷卻的數值研究也得到了相應的發展.氣膜冷卻基本原理是由葉片表面氣膜孔向主流注入冷氣, 從而在葉片表面形成溫度較低的冷氣膜, 將葉片表面與高溫主流燃氣隔離, 同時隨流動帶走部分熱量, 從而對葉片起到良好的保護作用.通常影響氣膜冷卻效果的因素有包括氣膜孔幾何參數、氣膜孔氣動參數、氣膜孔前邊界層發展情況、主流流動特性等.

McGovern 等[36]較早地對帶有氣膜孔的葉片開展了計算研究, 得到了單列冷卻孔平板氣膜冷卻的流動細節, 數值模擬沒有考慮到固體溫度邊界對流場的影響.Garg[37]通過數值計算發現, 旋轉葉片氣膜孔冷氣出口方向對氣膜的冷卻效率有很大的影響, 且在氣膜孔出口處具有很強的三維效應.Voigt 等[38]研究數值模擬方法的在對氣膜冷卻進行氣熱耦合計算時的精度影響, 同時在計算中比較不同湍流模型的優劣.計算結果表明, 通過求解URANS 方程的非穩態計算方法計算精度要比求解RANS 方程的穩態計算方法高, 而穩態計算中k-ω湍流模型計算結果比SST, SGS 模型更符合實驗結果.Takahahshi 等[39]則采用基于脫體渦模擬方法,對葉片前緣氣膜冷卻效率進行計算, 計算結果明顯優于基于S-A 和SST 湍流模型的定常計算, 能夠更好地捕捉氣膜冷卻的湍流結構和溫度波動(見圖9).姚玉等[40]和楊蓓潔等[41]利用SIMPLE 算法和二階迎風離散格式, 采用RNGk-ε湍流模型和非平衡的壁面函數, 計算網格采用非結構混合網格, 并在氣膜孔和壁面處進行網格加密, 計算了在不同葉柵通道主流雷諾數和二次流吹風比情況下, 在吸力面3 個典型弦向位置處開設收斂縫形孔對葉柵通道損失系數和氣膜冷卻效果的影響等, 此外還對肋條結構對氣膜冷卻凹槽葉尖流動換熱的影響進行研究,結果表明在葉尖增加肋條結構能夠起到減小葉尖高表面傳熱系數區, 提高葉尖平均氣膜冷卻效率的作用.孟通等[42]采用四面體網格, 對平板上4 種氣膜冷卻結構孔內渦流場進行數值模擬, 設置主流入口為速度入口, 二次流入口為質量入口, 空間離散選擇二階迎風格式, 通過對計算結果的分析得到了氣膜孔內渦結構對氣膜冷卻效率的作用規律, 在此基礎上通過對孔內渦結構的優化, 提出兩種氣膜孔改型結構能夠有效地提高氣膜冷卻效率.此外, 張晶輝等[16]和Wang 等[43]采用數值模擬方法研究了旋轉工況下吸力側葉尖結構對氣膜冷卻性能的影響, 結果表明吸力側葉尖區域的氣膜性能主要受葉尖泄漏渦和通道渦的影響: 葉尖結構對流場的影響很大;劉言明等[44]基于逆向工程方法對某內冷旋轉渦輪葉片進行了詳細的三維氣熱耦合計算, 分析了葉片溫度分布規律, 準確預測了渦輪最高溫度出現在前緣葉頂處, 為葉片冷卻結構設計優化提出了參考.

圖9 無量綱溫度染色下的瞬時渦結構和溫度等值面[39]Fig.9 Instantaneous vortex structures, colored by local non-dimensional temperature and isosurfaces of non-dimensional temperature[39]

綜上研究可以發現, 通過改變葉片構型實現改變葉片表面換熱特性的研究為早期葉片設計提供了很大的幫助, 但該方式起到的葉片冷卻效果有限.相比之下, 利用氣膜冷卻技術可以大幅提高渦輪葉片的承溫能力, 其在將來很長時間內都將是葉片冷卻設計的關鍵技術之一.但氣膜冷卻也給數值模擬方法的發展帶來了巨大的挑戰, 這要求數值模擬方法在如下情況下依舊保持魯棒高效: ①帶復雜冷卻通道的葉片;②轉子/靜子干涉;③冷氣與燃氣大密度比和溫度比;④輪轂及機匣附近復雜二次流動;⑤葉片吸力面和葉尖間隙附近強烈的激波/邊界層干擾;⑥冷熱交替的復雜湍流流動等.

綜合渦輪旋轉盤腔和旋轉葉片流動與對流換熱數值模擬研究進展可以發現, 目前對盤腔和葉片換熱特性的數值模擬大都會對模型進行適當的簡化, 同時對于發動機內部復雜的流動邊界進行近似處理, 使流場能夠收斂到一個相對合適的解, 從而實現對渦輪內部流場的數值預測和分析.而對于存在復雜流動與換熱問題的真實構型, 同時包含盤腔和葉片的數值模擬相對較少, 其中原因包括如下3個方面.

(1) 巨大的計算網格量和計算成本。真實渦輪存在眾多的細小結構, 為準確刻畫這些細小結構對流場的影響, 需使用更小尺寸的網格, 從而大幅增加計算成本, 同時會在網格布局和保證網格質量上花費大量時間.

(2) 無法獲取準確的邊界條件.渦輪位于發動機燃燒室之后, 其主流入口為燃燒室噴出的高溫高壓氣體, 其流動相當復雜, 對于數值模擬而言, 如何給定合適的邊界條件是準確模擬渦輪流動與對流換熱的一大難點.

(3) 缺少高精度的動靜交界面信息傳遞數值方法.對于存在相對運動的物理問題, 動靜交界面信息傳遞一直是不可避免的一大問題, 現常用的處理方法包括凍結轉子法和滑移網格方法, 二者均需要處理旋轉域和靜止域網格交界面上的數據傳遞, 而目前的數值方法還不能夠很好地保證交界面上的計算精度.

2 渦輪氣動噪聲研究現狀

渦輪作為典型的旋轉機械, 其氣動噪聲可分為純音和寬頻噪聲兩個部分.渦輪噪聲產生的原因主要包括: ①渦輪轉子靜子干涉;②葉片與來流湍流、不均勻燃氣相互作用;③葉片表面湍流邊界層引起的壓力脈動及與尾緣相互作用;④葉片表面渦脫落引起的隨機脈動, 及所導致葉片局部升阻力變化;⑤葉尖間隙流動及與于葉片尾緣、下游葉片相互作用等.航空發動機渦輪噪聲分量在20 世紀90 年代前后達到一個低點[45].但在新型航空發動機中, 為降低發動機重量和維護成本, 提升可靠性, 渦輪葉片稠度顯著降低, 同時緊湊的發動機設計使得轉靜葉盤間距縮小, 渦輪噪聲相對于上一代渦輪不僅沒有顯著下降, 部分新設計渦輪的噪聲反而有較大幅度提升[45-46], 這導致渦輪噪聲的重要性明顯增強.

1975 年, Mathews 等[47]綜述了渦輪噪聲預測方法, 主要包括半經驗模型[48]和解析方法[49]等.在此之前, 研究者們提出了管道聲模態理論[50]、陣風翼型干涉理論[51-52]等解析模型, 但實際應用中以半經驗模型為主.該時期的半經驗模型尚無法考慮渦輪中的非定常流動及相關復雜流動結構, 模型變量通常選為質量流率、溫度、壓比和渦輪轉速等宏觀參數, 非常依賴幾何和工況相似性, 因此產生了以各個主要發動機供應商或研究機構命名的預測模型, 如Rolls-Royce 方法、Pratt &Whitney 方法、NASA 方法和General Electric 方法等.在渦輪噪聲傳播和輻射方面, 研究者們引進了風扇噪聲的管道聲模態傳播、透射和截止相關理論, 但渦輪噪聲有其獨特傳播環境特性, 如高溫高壓、高載荷等, 同時噪聲與噴流剪切層相互作用將導致噴流噪聲寬頻化(hay-stacking)等.最關鍵的是噪聲預測模型精度亟需進一步提升.2011 年,Hultgren[45]總結了2003 年AARC 核心機噪聲研討會和2008 年AARC 渦輪噪聲研討會的討論結果, 并對渦輪噪聲進行了評論, 提出了渦輪噪聲存在問題和未來研究方向.Hultgren 認為20 世紀90 年代以前的發動機中大稠度渦輪葉片本身可以有效降低渦輪噪聲, 結合管道聲模態截止和聲學處理, 渦輪噪聲已經大幅降低.但是, 近年來渦輪設計趨勢是降低葉片數和轉子靜子間距, 這導致: ①增加了葉片載荷, 純音噪聲增強;②降低了葉片稠度, 增加透射傳播噪聲;③降低了葉片通過頻率(blade passing frequency, BPF), 增加噪聲響度和計權聲壓級.同時,針對熱端部件設計聲學處理裝置依舊面臨很大挑戰.

在數值模擬方法廣泛發展起來之前, 研究者們在解析方法上開展了大量的工作.1977 年,Cumpsty 等[53]基于波長大于葉片弦長和緊致渦輪盤假設, 從理論上分析了熵波、渦波和聲波與渦輪葉片的相互作用, 并提出了低頻渦輪噪聲的二維預測方法, 結果表明低速噴流情況下,低頻渦輪和間接燃燒噪聲很重要.2012 年, Duran 等[54]考慮了經過渦輪葉片的焓降影響, 進一步發展了Cumpsty 等提出的二維解析方法, 并與Doyle 等[55]的實驗結果進行了對比, 驗證了其在渦輪噪聲和傳播損失預測中的準確性.

與其他氣動噪聲預測方法相同, 渦輪噪聲的數值模擬方法主要包括直接計算方法和混合計算方法兩大類[56].直接計算方法最有代表性的包括直接數值模擬、大渦模擬和脫體渦模擬(detached eddy simulation, DES)混合方法等三類.直接計算方法通過采用高精度、低耗散、低色散的數值方法和準確的邊界條件, 直接求解可壓縮流體控制方程, 同時獲得流場和聲場信息.該方法的好處是不依賴于聲源簡化模型, 同時由于流場信息和聲場信息是同時獲得的, 故可考慮聲傳播過程中流場和聲場的相互作用.混合計算方法是一種分區求解的方法, 可以分成兩個或三個主要物理過程: 聲源區, 主要是復雜流動變化產生的流場脈動;中間聲場的噪聲傳播, 包括: 聲散射、反射、折射等聲波與環境發生相互作用的過程;遠場的噪聲輻射區.近場聲源主要通過計算流體力學方法求解獲得, 其難點是湍流的數值模擬.近場聲源通常通過諸如DNS、LES、DES 或URANS 等方法獲得.噪聲源通常需要采用高精度、低耗散、低色散的數值格式進行求解.中間聲場通?；谇蠼庵T如線性化歐拉方程(linearized Euler equation,LEE)[57]和聲擾動方程(acoustic perturbation equation, APE)[58]等聲傳播方程獲得.LEE 和APE 的基本思想是將流場脈動量與平均場進行線性分解, 利用高精度數值方法進行聲源求解,進而獲得速度和壓力脈動在背景流場下的分布.求解聲學方程開展聲場計算的難點是噪聲源的精確求解, 通常聲學方程的聲源可以基于DNS 或者LES 方法給出, 也可以通過隨機噪聲產生和傳播(stochastic noise generation and radiation, SNGR)方法[59]和隨機粒子網格(random particle mesh, RPM)[60]方法得到.這兩種方法通過求解RANS 獲得湍流平均信息生成噪聲源項, 能夠比較準確預測寬頻噪聲.遠場聲輻射的求解通常是基于聲比擬理論的積分方法.這種方法將近場求解獲得的速度和壓力脈動等流動信息轉化成聲源信號代入到聲輻射方程中去求解遠場噪聲, 其中Kirchhoff方程[61]和FW-H 方程是最具有代表性的求解包含固體壁面問題的聲輻射方程的主要形式.相比于Kirchhoff方程, FW-H 積分方程的顯著優勢在于考慮固體邊界對聲場分布的影響, 從物理角度清晰地解釋聲波輻射的物理本質, 因此成為當前簡單高效的遠場聲輻射計算方法.

20 世紀90 年代, 德國航空發動機公司的Kahl 等[62]針對葉片顫振等預測開發了時間線化歐拉方程求解器(Lin3D);Kennepohl 等[63]對該求解器進行了改進并用于渦輪轉子靜子干涉噪聲預測, 其研究結果表明黏性湍流尾跡干涉噪聲是渦輪噪聲中的主要聲源, 相比無黏預測聲壓級高10 dB 左右.Broszat 等[64]和Korte 等[65]使用Lin3D 求解器捕捉到了低壓渦輪噪聲的模態截止, 過適當的邊界條件調整后Lin3D 可以較好預測渦輪噪聲, 但對模態截止預測誤差較大, 他們指出多級渦輪中已經截止的模態可能再次截通(cut-on), 因此需要多級模擬.除了Lin3D, NASA 開發的TURBO 求解器也被用于渦輪噪聲的研究中, Zante 等[66]使用TURBO求解器研究了噪聲在單級渦輪中的透射傳播特性.

2013 年, 加拿大Sherbrooke 大學的Wang 等[67]和Papadogiannis 等[68]基于所開發的TurboAVBP 求解器, 通過大渦模擬研究了MT1 高壓渦輪流場, 測試了不同亞格子模型和網格分辨率對渦輪流動預測的影響, 結果表明: ①滿足RANS 分辨率情況下, 平均流場分布基本不受亞格子模型和網格分辨率影響;②激波強度和湍流結構顯著受到亞格子模型和網格分辨率影響.2015 年, Papadogiannis 等[69]進一步發展了TurboAVBP 求解器, 加入了熵波邊界條件, 并將其用于渦輪噪聲、熵波-渦輪相互作用引起的間接燃燒噪聲等.其研究表明渦輪噪聲在管道內以純音為主, 三維LES 模擬結果表明間接燃燒噪聲透射傳播效率低于二維模擬和緊致聲源理論模型.2016 年, Wang 等[70]使用TurboAVBP 求解器, 通過大渦模擬和動態模擬分解(dynamic mode decomposition, DMD)等探索了高壓渦輪噪聲中的主要發聲機理, 主要具體包括: 轉子-靜子干涉噪聲、鈍尾緣噪聲和激波-渦干涉噪聲.Papadogiannis 等的研究表明渦輪噪聲預測中仍未解決的關鍵問題是: ①來流湍流的引入;②機匣壁面分辨率需要進一步提升以提高間隙流動預測精度等.

3 結論與展望

為推動航空發動機渦輪流-熱-聲數值模擬方法的發展, 本文面向工程需求, 對航空發動機渦輪流動換熱與噪聲數值模擬研究的現狀進行了回顧.在此基礎上, 基于航空發動機旋轉盤腔/旋轉葉片熱流耦合求解器的開發經歷, 對渦輪流-熱-聲數值方法和軟件開發的發展做出了展望.

渦輪數值模擬研究在旋轉盤腔和旋轉葉片流動、換熱、噪聲預測及評估中已經取得了豐碩的成果.為應對新一代航空發動機渦輪冷卻設計、低噪聲設計中面臨的復雜流-熱-聲耦合問題帶來的挑戰, 數值模擬方法的發展及工程應用仍需在如下方面更進一步.

(1) 發展精度更高、色散耗散更小的數值格式, 發展魯棒性更好、效率更高、適用于極端流動換熱環境的計算方法.渦輪流動換熱數值模擬中需要考慮的因素遠多于風扇/壓氣機數值模擬, 除了轉子/靜子交界面處理、復雜二次流動和激波/邊界層干擾模擬等因素外, 還包括數值方法對帶復雜冷卻通道葉片幾何的適應性, 大密度比和溫度比冷熱交替湍流模擬等.這些額外的極端流動換熱環境對數值方法的計算精度、色散耗散、魯棒性和計算效率提出了更高的要求.

(2) 面向工程應用的求解器開發.目前國內無論高校還是工業部門, 對航空發動機渦輪數值模擬的研究大多數依然依靠國外成熟的商用軟件, 缺少自主可控的國產渦輪流動與換熱數值模擬軟件是目前國內航空領域面臨的一大難題, 發展形成面向工程應用的自主可控專用軟件是推動數值模擬深度參與渦輪設計的關鍵, 有助于進一步推動國家航空發動機渦輪設計技術的發展.

(3) 進一步提升渦輪噪聲預測能力.渦輪低噪聲設計的一個關鍵性難題是對渦輪噪聲的準確預測, 提升渦輪噪聲預測能力主要包括如下幾個方面: ①發展更接近真實環境的渦輪噪聲理論模型;②建立經過驗證可以考慮噪聲寬頻化的渦輪噪聲半經驗預測模型;③發展更準確的聲模態預測和分解方法;④進一步提升非設計點工況噪聲預測能力;⑤建立可靠的寬頻噪聲高精度數值方法;⑥建立符合物理機制的間接燃燒噪聲傳輸預測模型等.