軸對稱短艙自然層流優(yōu)化及轉捩敏感性分析

曹凡，張美芳，胡驍，唐智禮

南京航空航天大學 航空學院，非定?？諝鈩恿W與流動控制工業(yè)和信息化部重點實驗室，南京 210016

發(fā)動機是飛機的核心部件之一，為了保護發(fā)動機，降低阻力和噪聲，短艙應運而生。短艙主要由進氣道、整流罩和內、外涵道噴管組成，是推進系統(tǒng)的重要組成部分[1]。隨著渦扇發(fā)動機的發(fā)展，涵道比不斷增大，短艙直徑也越來越大，勢必帶來阻力的上升，影響飛行效率。面積越大，短艙的設計難度就越高，大涵道比發(fā)動機短艙是綜合了氣動設計、結構設計、材料和工藝等多種學科的復雜系統(tǒng)。為了減小短艙面積增大帶來的摩擦阻力，通過氣動設計的手段對大涵道比發(fā)動機短艙進行外形層流化設計是近年來興起的有效手段[2]。

相同來流條件下，層流的摩阻遠小于湍流摩阻，如果將層流技術應用于機翼、尾翼和短艙，能夠減少整機阻力15%以上，因此，層流減阻技術在新型民機的氣動設計中應用前景廣闊。自然層流技術是一種層流主動控制方法，該技術通過對短艙型面的設計，使得短艙外表面保持較大的順壓梯度區(qū)域，從而抑制流向T-S(Tollmien-Schlichting)波的增長，推遲轉捩的發(fā)生保持大范圍層流區(qū)域[3]。目前，機翼和機身的優(yōu)化減阻技術比較成熟，但對短艙的研究較少，在現(xiàn)有的經(jīng)驗和技術基礎上對短艙進行優(yōu)化，實現(xiàn)有效的減阻仍然是主要攻關問題。

近年來，國內學者對NLF短艙氣動外形設計也開展了一系列研究，但相關的風洞和飛行試驗卻很少。西北工業(yè)大學的何小龍等[12]研究了帶動力的自然層流短艙氣動優(yōu)化設計；孟曉軒等[13]基于雙eN方法研究了CRM(Common Research Model)短艙流動轉捩的影響因素，認為馬赫數(shù)和攻角對轉捩位置影響較大,而雷諾數(shù)和湍流度相對較小。2019年,復旦大學Wang等[14]應用差分進化算法優(yōu)化NLF短艙，將層流范圍提升了16.64%，獲得良好減阻效果;杜璽等[15]針對設計的自然層流短艙在德國ETW(European Transonic Wind)跨聲速風洞進行了對比驗證，結果表明低雷諾下轉捩預測結果吻合較好，而高雷諾數(shù)下轉捩提前發(fā)生與CFD(Computational Fluid Dynamics)結果差別明顯。目前來看，針對NLF短艙的氣動外形設計已經(jīng)取得較好的進展，但實際設計過程中，識別對層流影響顯著的變量更加重要，因此需進行參數(shù)靈敏度分析。高雷諾數(shù)下，保持大范圍的層流流動極其困難，在短艙層流化設計中，外形幾何控制參數(shù)對轉捩位置的影響和氣動參數(shù)與層流范圍的關系，尚有待研究。

以CRM短艙為研究對象，在高雷諾數(shù)跨聲速狀態(tài)滿足幾何約束條件下利用發(fā)展的進化/確定性混合算法進行氣動外形優(yōu)化，分別以總阻力最小化和層流最大化設計軸對稱NLF短艙，驗證層流短艙減阻的可能性。為了分析幾何控制參數(shù)、流量系數(shù)和馬赫數(shù)對層流范圍的影響，基于層流最大化的短艙構型開展轉捩敏感性研究，分析層流與阻力系數(shù)的關系，以期為后續(xù)自然層流短艙的設計和魯棒優(yōu)化提供依據(jù)。

1 層流轉捩預測

(1)

Langtry和Menter[17]采用的關聯(lián)函數(shù)為

(2)

(3)

(4)

式中：U為當?shù)厮俣?；θ為動量厚度；s為流線弧長；Ui為速度張量分量；xi為笛卡爾坐標張量分量；μ為層流黏性系數(shù)；μt為湍流黏性系數(shù)；S為流線弧長；k為湍動能；Tu為湍流度；λθ為壓力梯度；ρ為密度；Pγ、Eγ為生成項、耗散項；σr=1。

(5)

(6)

(7)

式中：Fwake為確保Fθ t在翼型下游尾流區(qū)為0的參數(shù)；y為靠近壁面最近距離；δ為邊界層厚度；cθ t為源項系數(shù)；σθ t為擴散系數(shù)，取2.0；ce2為耗散項系數(shù);t為時間。

σθ t的大小決定了臨界動量厚度雷諾數(shù)對流動歷史和當?shù)亓鲃訁?shù)的依賴程度，σθ t越大,轉捩對流動歷史越不敏感，更多地取決于湍流度和壓力梯度等當?shù)亓鲃訁?shù)。擴散系數(shù)和源項系數(shù)的變化雖然會帶來轉捩位置的改變，但與氣動參數(shù)和幾何參數(shù)的變化相比，其影響程度相對較小[18-19]。

混合函數(shù)Fθ t在邊界層外為0，在邊界層內為1，起到關閉Pθ t的作用。為了耦合k-ωSST湍流方程，引入有效間歇因子γeff與湍動能方程進行交互，表達式為

γeff=max(γ,γsep)

(8)

Fθ t

(9)

式中：γsep用于捕獲并模擬分離誘導轉捩；Freattach為一種開關函數(shù)，表示當黏度比大到足以引起層流再附著時，禁止其修改；Rev為渦量雷諾數(shù)；Reθc為臨界動量厚度雷諾數(shù)。

(10)

(11)

圖1 NLF(1)-0416摩擦阻力系數(shù)曲線

表1 NLF(1)-0416翼型轉捩結果對比[20-21]

DLR-F5機翼前緣后掠角20°，采用超臨界對稱翼型。計算工況:Ma=0.82，α=2°，Re=1.5×106，參考弦長0.15 m，湍流強度0.5%，湍流黏性比10，物面表面網(wǎng)格y+<1，第1層網(wǎng)格高度L=1.5×10-6m，網(wǎng)格量約為247萬，如圖2所示。

圖2 DLR-F5計算網(wǎng)格

圖3 DLR-F5試驗結果與數(shù)值模擬結果對比

2 高雷諾數(shù)下軸對稱自然層流短艙優(yōu)化設計

大型民用運輸機短艙的流動雷諾數(shù)一般高于2×107，為了實現(xiàn)高雷諾數(shù)下的層流流動，需要保持更加強烈的順壓梯度，這遠比在低雷諾數(shù)下的優(yōu)化更加困難。通過精確的參數(shù)化方法、高精度流場求解、可靠的轉捩預測結合高效的優(yōu)化算法進行高雷諾數(shù)下的自然層流短艙優(yōu)化設計。

2.1 軸對稱短艙幾何設計

為了設計滿足幾何約束的短艙，同時對比優(yōu)化前后的層流范圍，首先對CRM標模短艙[22]進行自然通氣狀態(tài)下仿真。設計工況：巡航高度11 000 m，基于短艙參考弦長雷諾數(shù)3.5×107，自然來流馬赫數(shù)Ma=0.85，攻角α=0°，湍流強度Tu=0.05%，湍流黏性比為10，流場計算采用有限體積方法，空間離散為二階精度，通量計算采用 Roe 格式，時間推進為 LU-SGS (Lower Upper Symmetric Gauss Seidel) 隱式方法。

計算網(wǎng)格如圖4所示，劃分六面體結構網(wǎng)格共850萬，其中物面流向和周向各201個網(wǎng)格節(jié)點，第1層網(wǎng)格高度4×10-6m。

圖4 CRM短艙物面網(wǎng)格

計算結果如圖5所示。CRM短艙在自然通氣狀態(tài)下，其層流范圍很小，且流動在靠近唇口的位置便開始轉捩發(fā)生流動分離。

圖5 CRM短艙計算結果

經(jīng)測量CRM短艙，其最大弦長6.06 m，半徑1.75 m，為了滿足幾何約束(見表2)，對CRM短艙剖面進行CST(Class Shape Transformation)參數(shù)化[23]。

表2 短艙設計約束條件

CST方法由一個類別函數(shù)和一個形狀函數(shù)構成，具有精度高和設計參數(shù)少等特點，其設計參數(shù)具有明確的物理意義，能夠直接反映翼型的變化過程。其數(shù)學表達式為

(12)

式中：c為翼型弦長；Rle為翼型前緣半徑；Δzte/c為后緣點縱坐標；β為后緣傾角；bi為形狀函數(shù)控制參數(shù)(i=1,2,…,n-1)。采用基于Bernstein多項式的CST方法時，為了提高參數(shù)化精度,且防止過高的階數(shù)使得參數(shù)化趨于病態(tài)，取n=6。此時，參數(shù)化短艙外剖面共需要8個設計變量(Rle、Δzte、β、b1、b2、b3、b4、b5)。

對于流量系數(shù)的約束，采用了閆海津和杜璽[24]提出的在短艙出口處添加堵錐的方式代替動力短艙實現(xiàn)流量系數(shù)的連續(xù)調節(jié)，初始流量系數(shù)(Mass Flow Ratio,MFR)MFR=0.7，最終形成的滿足約束的軸對稱短艙半模如圖6所示。

圖6 帶堵錐的軸對稱短艙

2.2 自然層流短艙優(yōu)化設計方法

通常更大范圍的層流要求短艙表面維持較大范圍的順壓梯度，這會導致短艙后半段的壓力快速恢復產生強激波，增大壓阻和波阻從而對總阻力的減小產生負面效應。本文著眼于自然層流短艙的優(yōu)化設計，期望在短艙的表面實現(xiàn)更大范圍的層流，但對于具體的可優(yōu)化程度沒有經(jīng)驗參考，伴隨摩擦阻力減小對激波阻力與總阻力的影響也未知。因此，為了更好地分析自然層流短艙的局部敏感性，分別以轉捩位置最大化和總阻力最小化為優(yōu)化目標進行優(yōu)化，分析不同參數(shù)影響下的層流范圍，驗證層流減阻的可能性。希望獲得一款能夠保持較大層流范圍的短艙進行敏感性分析，有利于研究層流范圍變化規(guī)律。

2.2.1 進化/確定性混合優(yōu)化算法

解決氣動優(yōu)化問題可以用3種不同的方法：① 無導數(shù)搜索算法(Derivative-Free Search Algorithns,DFSAs)；② 基于梯度的算法(Gradient-Based Algorithms,GBAs)；③ 進化算法(Evolntionan Algorithms,EAs)。每一種算法都有各自的優(yōu)缺點，可以結合不同的算法[25]，把其優(yōu)點同時發(fā)揮出來，比如:EAs的魯棒性和全局搜索,以及GBAs的快速收斂能力。通常使用基于梯度的搜索方法，計算目標函數(shù)的導數(shù)以確定最佳的下降方向，梯度法已與遺傳算法(Genetic Algorithm,GA)結合用來設計翼型和機翼[26-27]。

以GA為主算法，用擬牛頓方法(BFGS算法進行)輔助加速收斂，發(fā)展了一種新的多級混合優(yōu)化算法。文獻[28]對該混合算法的原理進行了詳細介紹。進化和確定性算法在這種混合優(yōu)化方法中交替進行，并成功應用于高超聲速飛行器氣動外形優(yōu)化，提升了優(yōu)化效率。GA/BFGS算法流程如下：

Step1初始化種群，并計算適應值。定義初始種群數(shù)量Np進化代數(shù)Ke，局部搜索選擇頻率為f。

Step2執(zhí)行GA算法，產生較優(yōu)種群。產生一個隨機整數(shù)Ns(Ns≤Np)，并從當前種群中選中并標記Ns個個體。

Step3對標記的Ns個個體執(zhí)行確定性優(yōu)化算法BFGS，BFGS將這Ns個個體每個迭代Kb步，產生Ns個新個體,并將其加入種群替代原有個體。

Step4最優(yōu)選擇(即從GA種群中選擇前Ns個最優(yōu)個體)，新的搜索空間簡化為收斂的較小空間。隨著搜索空間的縮小，種群規(guī)模和搜索頻率也可以減少。

Step5繼續(xù)對簡化的搜索空間和局部搜索頻率進行優(yōu)化，直至收斂。

多級方法在這個交替的過程中耦合。為了保持基因的多樣性，可以在進化過程開始時使用大種群規(guī)模，以便選擇更好的基因并傳給下一代。該算法交替過程如圖7所示。

圖7 GA/BFGS混合優(yōu)化算法示意圖

為驗證算法的有效性，對(GA/BFGS算法)與純梯度(BFGS)算法[29]、差分進化算法(Differential Evdution,DE)[30]、GA算法[31]在Ackley多峰函數(shù)上進行了對比驗證，如圖8所示。種群個體均為100、代數(shù)5 000、交叉概率0.9、變異概率0.1，差分算法中縮放因子取0.5。

Ackley多峰函數(shù)為

(13)

Ackley函數(shù)取值范圍為x∈[-32, 32]，維度d=30，具有很多局部極小值，雖然梯度算法進化效率很高，但容易陷入局部最優(yōu)解，而GA、DE這種進化算法能夠利用交叉、變異算子跳出局部最優(yōu)解，進化效率相對較慢。本文所使用的混合算法在進化前期效率不如梯度算法，如圖8所示，但避免了陷入局部最優(yōu)解，很好地平衡了進化效率和尋優(yōu)能力。

圖8 Ackley 函數(shù)收斂曲線

2.2.2 軸對稱短艙自然層流設計

自然層流短艙設計過程中，目標是在激波前保持順壓梯度，抑制T-S波轉捩，從而獲得大范圍層流區(qū)域。采用GA/BFGS混合優(yōu)化算法結合六階CST(形狀類別變化函數(shù))參數(shù)化方法，以總阻力最小化和轉捩位置最大化為目標分別進行單點優(yōu)化。

短艙計算工況：巡航高度11 000 m，基于短艙最大長度雷諾數(shù)Re=2.5×107，流量系數(shù)MFR=0.7，湍流度Tu=0.05%，黏性比為10。劃分結構網(wǎng)格如圖9所示，第1層網(wǎng)格高度5×10-3m，網(wǎng)格前緣用C型網(wǎng)格保證正交性，并進行加密處理，軸向分布200個網(wǎng)格點，網(wǎng)格數(shù)量約11萬。

圖9 軸對稱短艙剖面網(wǎng)格

多層混合交替初始階段種群規(guī)模Np=40，局部搜索選擇頻率為f=2，進化代數(shù)100代，保證優(yōu)化過程中的基因多樣性，隨著搜索空間的縮小，種群規(guī)模也可以減少，這里減少到20個，局部搜索頻率減少到1，繼續(xù)進化10代，確保最后幾代個體能夠收斂到全局最優(yōu)解。

計算結果如表3所示，Baseline為短艙設計工況，minCD是以總阻力系數(shù)最小為優(yōu)化目標的工況，maxXtr是以轉捩位置最大為優(yōu)化目標的工況。ΔCD為短艙內外剖面的總阻力系數(shù)相對于基準外形的變化量，由于此次設計不改變內型線，其值變化可近似代表外剖面的阻力變化；Xtr為短艙轉捩位置；Cf-up為短艙外剖面的摩擦阻力系數(shù)，參考面積為383.7 m2；Dmax為短艙外罩最大直徑；Rle為短艙前緣半徑。

表3 軸對稱短艙優(yōu)化前后計算結果對比

對比分析發(fā)現(xiàn)在總阻力系數(shù)最小時，轉捩位置推遲了3.6%，摩擦阻力也略有降低，但短艙外罩最大直徑大大降低，短艙剖面面積減小約28.5%，前緣半徑變小，如圖10所示。當層流面積最大化時，短艙外剖面的摩擦阻力減小24%，表明延遲轉捩，增大層流范圍極大地降低摩擦阻力，但也使得總阻力明顯增加。這表明合理的平衡層流面積最大化和總阻力最小化2個目標，才能實現(xiàn)真正的短艙減阻目的。

圖10 優(yōu)化前后短艙外形對比

歸一化后Cp、Cf曲線如圖11、圖12所示，原始外形由于唇口吸力峰的存在過早轉捩，層流范圍最大化的外形在前半段維持了良好的順壓梯度，可達到總長的27%左右，但是在短艙尾部出現(xiàn)由激波導致的強逆壓梯度，壓力快速恢復，出現(xiàn)相對強烈的激波，這是導致阻力略有增加的原因。阻力最小化外形后半段外剖面迅速向內剖面靠近，壓縮了短艙面積。

圖11 優(yōu)化前后壓力系數(shù)曲線對比

圖12 優(yōu)化前后摩阻系數(shù)曲線對比

3 軸對稱自然層流短艙局部敏感性分析

基于單點優(yōu)化設計，獲得了2種不同的短艙構型，研究表明如果以總阻力為優(yōu)化目標最后的優(yōu)化結果的確會使總阻力減小，但是層流轉捩位置推遲較小，層流范圍增大不明顯。單點優(yōu)化的設計結果在偏離設計點時，往往不能保證飛行性能。同時，為了提取層流短艙的設計規(guī)律，分析了幾何參數(shù)和飛行參數(shù)的敏感性。為了更好地分析單變量影響下軸對稱自然層流短艙的轉捩變化規(guī)律，選取層流范圍最大的短艙構型為初始構型，較大范圍的層流更容易對比分析不同參數(shù)帶來的影響。

3.1 局部單變量幾何敏感性

由分析可知，短艙剖面采用六階CST參數(shù)化方法進行參數(shù)化，考慮幾何約束去掉后緣厚度Δzte后還有7個參數(shù)，分別為前緣半徑Rle、后緣傾角β、短艙剖面形狀函數(shù)的控制參數(shù)(b1、b2、b3、b4、b5),控制范圍如圖13所示。在保持飛行狀態(tài)不變的情況下，對CST參數(shù)化設計變量初始值的(-15%, 15%)區(qū)間如表4所示，取6組數(shù)據(jù)進行單一變量的幾何敏感性分析，并給出Cp以及歸一化后轉捩位置,如表5、圖14所示。

表4 設計變量取值范圍

通過對比表5的數(shù)據(jù)，按對轉捩位置影響大小進行排序，轉捩位置隨b1、b2、β的變化最為明顯，其次為b5、Rle，影響最小的參數(shù)為b3、b4。其中b1對層流范圍的影響最為顯著，當其降低10%時，層流范圍降至10%以下。由圖13可知,b1主要控制短艙截面的頭部區(qū)域，在短艙的設計加工中需要重點考慮。b2、β變化下，層流范圍在10%、20%內波動較大，是次要考慮的設計點。而b3、b4的變化對轉捩位置的影響不敏感，層流范圍都可以保持在20%以上。

圖13 形狀控制參數(shù)控制范圍

表5 不同幾何參數(shù)變化下的轉捩位置

通過分析層流短艙的初始構型壓力分布如圖14 中所示，在短艙截面前半部分保持較長范圍的順壓梯度，緊接著出現(xiàn)由激波導致的強逆壓梯度，壓力迅速恢復?？缏曀贍顟B(tài)下，流動從翼型前緣開始，穿過一系列膨脹波加速到超聲速狀態(tài)，而后保持超聲速流動直到后緣附近。由于壓力需要恢復到后駐點處的壓力值，同時由超聲速狀態(tài)變?yōu)閬喡曀贍顟B(tài)，因此在層流短艙截面后緣附近產生一道激波，通過逆壓梯度的大小可以評估激波的強度。

圖14 自然層流短艙幾何參數(shù)變化下的Cp曲線

隨著b1的增大，曲線吸力峰值減小，外形表面的壓力變化更為平緩，轉捩位置變化較小，但b1減小時吸力峰變強使得提前發(fā)生轉捩，層流范圍變??；b2對前緣吸力峰的影響也較大。同時由于轉捩位置變化較大，導致后緣處激波位置和強度變化也較大；b3、b4、b5、β在吸力峰后的壓力變化趨于平緩，對壓力分布的影響較??；隨著前緣半徑Rle增大前緣吸力峰變強，對激波強度的影響較小。

3.2 流量系數(shù)敏感性

MFR是指進入進氣道的實際空氣質量流量與同馬赫數(shù)下進入進氣道的理論最大質量流量之比[24]。本文優(yōu)化后的軸對稱層流短艙MFR=0.699，在(0.67，0.73)范圍內通過移動堵錐實現(xiàn)流量系數(shù)的調節(jié)，轉捩位置變化如圖15所示。

圖15 轉捩位置隨流量系數(shù)變化關系

流量系數(shù)增加的過程中，轉捩位置的大小波動明顯，MFR=0.678時，轉捩位置前移到0.4 m以下?？梢园l(fā)現(xiàn)，巡航狀態(tài)下進氣流量對層流影響很大，此構型下盡量保持MFR>0.69才能獲得較大范圍層流。

3.3 馬赫數(shù)敏感性

短艙能否在不同馬赫數(shù)下產生較為強烈的順壓梯度以保持穩(wěn)定的層流范圍以適應飛行狀況的變化是非常關鍵的。此次馬赫數(shù)的單變量影響分析圍繞巡航點依次在(0.8, 0.9)內取了10個工況，并在亞聲速區(qū)間取4個工況,分析該自然層流短艙的變化規(guī)律，結果如圖16～圖18所示。

圖16 轉捩位置隨馬赫數(shù)變化關系

圖17 Cp隨馬赫數(shù)變化曲線

圖18 Cf隨馬赫數(shù)變化曲線

壓力和阻力系數(shù)曲線都對馬赫數(shù)變化體現(xiàn)出了明顯的敏感性。需要特別指出的是當Ma=0.87 時，轉捩位置出現(xiàn)強烈的波動，經(jīng)重復計算排除了計算意外。除此之外，馬赫數(shù)越大順壓梯度越強，吸力峰小幅降低，激波增強且后移。尤其是在飛行速度相差較大的情況，轉捩位置在亞音速區(qū)間比跨聲速區(qū)間更加敏感?？缏曀贍顟B(tài)下，隨著馬赫數(shù)增大，短艙整體上能夠保持層流范圍大于20%，層流范圍波動較小。

4 結 論

在幾何約束條件下,針對短艙剖面開展軸對稱自然層流短艙優(yōu)化設計研究，并利用CST參數(shù)化后帶有物理意義的幾何參數(shù)，分析了幾何擾動、流量系數(shù)和馬赫數(shù)的轉捩敏感性。主要得到以下結論：

2) 基于層流最大化的軸對稱短艙構型，開展了局部幾何控制參數(shù)敏感性研究。外形控制參數(shù)b1、b2對壓力分布和層流范圍影響最大，而前緣半徑Rle、后緣傾角β對壓力分布影響很小，是層流范圍的次要影響因素。b3、b4、b5的敏感性相對較低，其主要影響短艙外表面湍流段的壓力分布。為保持大范圍層流，要重點關注短艙剖面頭部外形的設計。形成了基于CST控制參數(shù)的短艙剖面外形設計準則。

3) 轉捩位置會隨著流量系數(shù)和馬赫數(shù)的增大呈波動狀態(tài)，敏感性較高。在設計點(MFR=0.7，Ma=0.85)附近，層流范圍能夠保持在20%以上，降低了摩擦阻力，除了Ma=0.87以外，設計點附近均能保持較大范圍的層流，為后續(xù)魯棒性優(yōu)化設計提供指導。