楊釗,李杰,牛笑天
西北工業大學 航空學院,西安 710072
當前,在氣動力相關性分析研究中一個主流趨勢是將CFD、風洞試驗和飛行試驗這3種用于確定飛機氣動特性的技術手段緊密結合起來為飛行器的設計提供有力的支撐。圖1顯示了3種方法各自的特點和它們之間的區別。CFD和風洞試驗一般用于在設計研制階段初步確定飛行器的基本氣動特性,而飛行試驗則通常用來在最后階段對飛行器真實飛行性能進行綜合評估。通過建立3種分析手段之間的關聯模型,便可以利用CFD手段輔助數據修正,發展從地面風洞試驗數據向天上真實飛行數據的外推方法,從而盡可能在飛行試驗開展之前完善飛行器的總體設計[1],同時也能夠為飛行試驗計劃提供參考。對于CFD計算分析而言,可以針對全尺寸飛機幾何外形開展數值模擬研究,狀態可以涵蓋完整的飛行工況,并且不存在雷諾數效應問題,能夠作為縮比風洞試驗和飛行試驗之間的橋梁,在風洞試驗數據的修正和飛行結果的對比校驗中都發揮著相當重要的作用[2]。

圖1 確定飛機氣動特性的技術手段和方法
相關性分析中最重要的一個方面就是對飛行器雷諾數效應的分析和修正。對于一般的附著流動來說,雷諾數的大小影響模型表面上附面層的性質,從而改變附面層的厚度、附面層轉捩位置、表面摩擦阻力以及與氣體黏性有關的氣流分離情況。這樣,在試驗雷諾數與飛行雷諾數相差較大的情況下必然導致風洞試驗數據與飛行情況不一致[3-4]。換言之,兩種雷諾數之間的差異使風洞試驗數據的準確度受到了影響。雷諾數效應又可分為直接雷諾數效應和間接雷諾數效應[1]。直接雷諾數效應是與恒定壓力分布有關的效應,而間接雷諾數效應則與雷諾數變化時壓力分布的變化有關。通常依賴于間接雷諾數效應的氣動特性是:升力和俯仰力矩、波阻、阻力發散和抖振邊界;而依賴于直接雷諾數效應的特性主要有:黏性阻力、邊界層分離特性和抖振邊界。眾所周知,邊界層的流動轉捩位置也會隨雷諾數變化而發生改變。由于層流邊界層更容易分離,而湍流邊界層在一定雷諾數下具有更高的壁面摩阻,因此飛機的阻力和升力都會隨著轉捩位置的改變而發生不同程度的變化。因此,雷諾數的變化不僅會直接影響到全機的主要氣動力特性,還能夠通過改變表面流動轉捩位置來間接對升、阻力及力矩特性等產生影響[1]。
在風洞試驗中,模型尺寸遠小于真實飛機,因此試驗雷諾數遠小于飛行雷諾數。這種由雷諾數差異所引起的“雷諾數效應”對附面層的結構和飛機的氣動特性往往有著不可忽略的影響。一般來說,低雷諾數風洞試驗結果主要會存在以下問題[4]:
1) 大迎角時氣流分離比真實飛行提前,影響飛機的最大升力系數以及其他與分離有關的氣動特性。
2) 飛機表面摩擦系數與真實飛行相比偏大,雷諾數的差異會顯著影響飛機的阻力系數(最小阻力系數、零升阻力系數、升致阻力因子等均是雷諾數的敏感項)。
3) 風洞試驗中飛機表面附面層狀態與真實情況下有很大不同,從而使表面各種小的突出物對全機阻力的影響也不同。
4) 雷諾數會對物面邊界層的流動轉捩特性產生顯著影響,進而影響到全機的升阻特性和壓力分布特征。通常來說,風洞試驗中所得到的層流流動轉捩位置與真實情況相比會有所推遲。
近年來,國內外利用高雷諾數或變雷諾數試驗[5-10]針對各種類型的飛機開展了大量的關于雷諾數效應的試驗研究,取得了豐碩的成果。但是大量的風洞試驗所帶來的高昂費用是設計人員們不得不面臨的一個嚴峻的問題。因此,通過數值模擬方法來研究飛行器的雷諾數效應逐漸成為研究熱點。Pettersson等開展了運輸機構型雷諾數影響的數值研究[11];Said等通過數值模擬方法研究了雷諾數對三角翼流場和氣動特性的具體影響[12];Rudnik等利用CFD方法研究了雷諾數的差異對某增升構型低速氣動特性的影響規律[13];Curtin等則深入探究了跨聲速飛行狀態下雷諾數效應對雙發運輸機構型氣動性能的影響[14]。近年來,國內同樣在這方面開展了很多的研究工作。周林[15]、馬明生[16]等通過數值計算研究了雷諾數效應對于運輸機構型氣動性能的具體影響規律;張培紅等[17]采用CFD方法,同時結合風洞試驗數據開展了雷諾數對某飛機氣動性能的影響研究;張彥軍等[18]采用試驗和數值計算兩種手段分別研究了雷諾數差異對于超臨界翼型氣動性能的影響;李中武等[19]開展了低速增升裝置變雷諾數的計算研究。通過這些研究進一步加深了研究人員對于雷諾數效應的認識,也更加清晰地得到了雷諾數對于飛機氣動性能的具體影響規律。
圖2中的階段1部分給出了本文的主要研究路線和內容,同時也表明本文研究工作在整個層流機翼設計、驗證和評估項目中的重要性。文中利用數值計算方法針對特殊布局形式的層流機翼飛行驗證平臺低速起降和高速巡航構型在試驗和飛行雷諾數下進行數值模擬分析,以得到全機氣動力系數和流場結果。通過計算結果與風洞試驗數據的對比分析,進行已有數據的相互校驗。同時,深入剖析雷諾數效應對飛行驗證平臺高、低速氣動性能、低速失速分離特性和高速層流轉捩特性等的具體影響規律,并據此對低速和高速試驗數據進行修正,為后期階段2中飛行試驗的設計提供數據支撐,協助完成國內首次基于高亞聲速無人驗證平臺開展的高空層流流動驗證工作,實現階段3中層流機翼設計效果的綜合評估。

圖2 研究路線圖
本文的研究對象是航空工業第一飛機設計研究院研制的層流機翼飛行驗證平臺,是國內首個為開展層流飛行試驗設計的高亞聲速無人飛行平臺。該平臺氣動布局采用雙機身、π型尾翼設計,外翼下吊裝4臺小型渦噴發動機。機翼分為內、外翼兩段,內翼段為層流驗證段,上表面光滑,且無運動翼面。外翼段后緣內側為單縫簡單偏轉式襟翼,供起飛、著陸時增升使用,外側為副翼控制面。尾翼為π字型結構,雙垂尾底部分別和雙機身尾部相連,頂部則與平尾相連。垂尾和平尾后緣分別設有方向舵和升降舵。該特殊布局層流機翼飛行驗證平臺的主要幾何參數如表1中所示。

表1 飛行驗證平臺幾何參數
針對該層流機翼驗證平臺的低速風洞試驗是在航空工業空氣動力研究院的FL-8低速風洞中進行的,風洞試驗中采用的是1∶3.25的縮比模型,試驗雷諾數大概在150萬量級。風洞試驗模型示意圖如圖3所示。該低速試驗采用腹撐式天平進行測力,并通過圖 4所示的對稱天平測量試驗來定量扣除天平裝置氣動外形對飛機流場和氣動力的影響。

圖3 低速風洞試驗模型示意圖

圖4 天平影響測量試驗
高速風洞試驗是在航空工業空氣動力研究院的FL-2風洞中進行的,風洞試驗中采用的是1∶7的全金屬通氣縮比模型,試驗雷諾數大概在300萬量級。風洞試驗模型示意圖如圖5所示。全模測力主支撐形式為中央翼段下部畸變的尾支撐形式,畸變帶來的支撐干擾通過雙支桿支撐方式進行扣除。試驗模型在風洞試驗中的具體支撐形式如圖6所示。

圖5 高速風洞試驗模型示意圖

圖6 試驗模型在風洞中的支撐示意圖
數值模擬計算主要采用NASA的CFL3D求解器。其作為一種多塊結構網格求解器,計算效率高,計算精度基本滿足工程應用要求,適合作為飛機設計研制過程中的計算分析和設計輔助手段。
CFL3D求解器采用積分形式的雷諾平均Navier-Stokes方程作為流場控制方程:
(1)
式中:V為控制體體積;S為控制體表面積;Q為守恒變量向量;f為通過控制體表面的通量向量,包含3個方向的無黏和黏性通量;n為控制體表面的外法向單位向量。
采用格心格式的有限體積法對控制方程進行離散,其中黏性項的空間離散采用二階中心差分格式,無黏項則采用二階Roe迎風通量差分格式離散。時間項采用雙時間隱式近似因子分解法(A-F)進行推進,計算中應用了多重網格技術加速流場收斂。

(2)
γ間歇因子輸運方程的作用是進行轉捩的判定和對轉捩區的發展進行預測,具體形式如下:
(3)
為了模擬分離誘導轉捩,對輸運方程做了特殊處理,使得γ在流動分離處快速增長,以模擬層流分離/湍流再附現象,主要的修正如下:
(4)
(5)
γeff=max(γ,γsep)
(6)
上公式中涉及到相關變量的具體定義詳見文獻[20-22]。
CFD計算采用如圖7中所示的真實尺寸全機模型,模型中設置有襟翼、副翼、升降舵和方向舵等部件,舵面與其他部件之間留有縫隙。計算網格采用多塊結構化網格劃分,計算區域的遠場邊界流向前后取30倍平均氣動弦長,展向和豎直方向各取為20倍。附面層第1層網格高度為5.0×10-6m,法向增長率為1.15,物面網格的y+值均小于1[23]。y+是與壁面第1層網格間距和具體流動參數直接相關的無量綱壁面距離,表征的是第1層網格點在附面層分區結構中的位置。對于中央翼段進行流向加密,以滿足高速狀態層流計算的要求,最終總網格量為3 000萬左右,計算中所采用的驗證平臺表面網格如圖8所示。

圖7 計算模型三視圖

圖8 計算網格
針對層流機翼驗證機低速起降構型(襟翼下偏30°)和高速干凈構型(襟翼下偏0°)不同馬赫數飛行狀態,分別開展試驗和飛行雷諾數下的數值計算分析,所有的CFD計算狀態如表2所示。通過氣動力系數計算結果與試驗數據的對比驗證了方法的計算能力和可靠性。同時利用不同雷諾數下的計算結果分析總結雷諾數的差異在低速和高速狀態下對全機氣動性能的影響程度和具體規律。

表2 CFD計算狀態
針對低速起降構型飛行雷諾數和試驗雷諾數下的數值計算結果(Cal)與試驗數據(Exp)的對比如圖9所示。可以看出,-4°~8°迎角之間計算所得升力系數CL與試驗數據吻合較好;計算與試驗零升力矩相差不大,但試驗結果中的縱向穩定度略大于數值預測結果。低速狀態下雷諾數的差異對線性段小迎角狀態的升力和力矩系數Cm影響不大,主要會對失速迎角附近的全機升力特性產生較為明顯的影響。高雷諾數狀態全機升力線性段由低雷諾數狀態下的8°迎角增加至10°以上,最大升力系數也由低雷諾數時的1.28增加至1.35左右。可見,在低速時雷諾數的增加會明顯提升大迎角狀態下機翼表面流動的穩定性,進而改善層流機翼飛行驗證平臺的失速特性。

圖9 Ma=0.2不同雷諾數下的全機氣動力系數結果
圖10給出了Ma=0.2,不同迎角下兩種雷諾數狀態計算所得全機表面壓力分布云圖和流線的對比。從圖中可以明顯看出雷諾數對于大迎角翼面流動狀態的具體影響規律:在8°~12°迎角范圍內高、低雷諾數狀態下中央翼段上表面流動狀態總體差異不大,而外翼段表面的流線形態差別則非常明顯。低雷諾數狀態外翼段外側在8°迎角時即存在顯著流動分離現象,并隨著迎角的增大分離區域進一步擴大;而高雷諾數狀態下直到12°迎角時翼面上大部分區域依舊保持著較好的附著流動,只在翼根部分出現了局部的流動分離現象。這充分表明雷諾數的增加能夠顯著提升大迎角狀態機翼表面流動的穩定性。

圖10 Ma=0.2不同迎角下的表面壓力分布云圖及流線
風洞試驗雷諾數計算狀態,8°迎角下機翼表面出現局部流動分離,升力系數進入非線性區域,針對8°迎角下的表面壓力分布和流線結果開展進一步分析。圖11給出了機翼展向壓力分布截面所對應的站位信息。圖 12則顯示了8°迎角下自由轉捩計算所得不同雷諾數下機翼展向不同站位截面壓力分布曲線的對比。可以看出,試驗雷諾數下中央翼段和外翼上表面前部均存在明顯的鋸齒狀,表明在此處預測到前緣層流分離泡現象;而對于飛行雷諾數計算結果,其前緣壓力分布曲線相對較為光滑,未能從中觀察到明顯的層流分離泡現象。圖13兩種雷諾數下中央翼段前緣表面壓力分布云圖及流線可以更為直觀的看出不同計算結果的差異。試驗雷諾數下前緣分離泡的尺寸和分離程度要明顯大于飛行雷諾數下的計算結果。

圖11 機翼展向壓力分布截面站位示意圖

圖12 α=8°機翼展向不同站位壓力分布曲線對比

圖13 Ma=0.2不同雷諾數下的中央翼段前緣表面壓力分布云圖及流線
圖13中所表現出的現象與圖14中央翼段前緣附近表面摩擦力系數云圖的結果相一致,試驗雷諾數下上表面前緣附近存在明顯的分離區。從圖15中展向15%站位和65%站位表面摩擦力系數分布曲線的對比中能夠更加直觀的看出兩種雷諾數下機翼前緣流動分離特性的差異。在試驗雷諾數下兩個站位處均有明顯的層流分離泡存在,尤其是在外翼段,其前緣分離區域長度將近局部弦長的25%;而對于飛行雷諾數,15%站位處前緣分離極小,而在65%站位處表面分離區長度也不及局部弦長的10%。

圖14 Ma=0.2不同雷諾數下的中央翼段前緣表面摩擦力系數分布云圖

圖15 α=8°機翼展向不同站位表面摩擦力系數分布曲線對比
針對干凈構型高速不同馬赫數狀態的計算氣動力系數與試驗結果的對比如圖16中所示。可以看出,雷諾數對不同馬赫數下的全機升力系數影響相對較小。從300~1 000萬雷諾數以上線性段升力系數均與試驗數據吻合良好。雷諾數的差異對于全機的阻力和縱向力矩系數影響較為明顯。由于采用自由轉捩方式進行計算,雷諾數對阻力的影響規律與全湍流計算狀態下有所不同,試驗雷諾數300萬狀態下小迎角阻力系數最小,中等和高雷諾數狀態阻力系數差異不大,均略高于低雷諾數狀態。這主要是因為雷諾數增加會導致翼面流動轉捩特性發生變化,進而改變全機表面摩擦阻力量值。由于驗證機在試驗和計算中所采用短艙模型略有差異,試驗和計算所得阻力量值整體差距偏大,但是數值計算結果中依舊能夠反映出雷諾數對整體阻力特性的影響規律。從全機縱向力矩特性曲線來看,不同馬赫數狀態試驗雷諾數下的計算結果在絕對量值和縱向穩定度方面均與試驗數據吻合得更好。對于馬赫數0.6和0.8狀態試驗雷諾數下的力矩結果與中、高雷諾數差異較大;而馬赫數0.7下3種雷諾數的計算結果之間則差異相對較小。

圖16 不同馬赫數狀態全機氣動力系數曲線
高速飛行狀態下雷諾數的差異不僅會影響全機基本氣動力系數,還會對全機表面,尤其是中央層流驗證段上表面的層流轉捩特性產生明顯的影響。風洞試驗中來流湍流度大致為0.3%~0.4% 之間,計算采用的自由來流湍流度為0.3%,黏性比為10,計算中維持流動湍流度沿流向不發生明顯衰減。圖17中展示了Ma=0.7、迎角為2°時不同雷諾數下的全機表面摩擦阻力系數云圖。可以看出,隨著雷諾數的增加層流轉捩位置呈現出逐漸提前的趨勢,當雷諾數達到11.90×106時中央翼段上表面僅剩22%倍弦長的層流區長度。圖18給出了不同雷諾數下預測所得中央驗證段上表面轉捩位置與試驗結果的對比。在試驗雷諾數的范圍內試驗結果與數值預測結果有一定的差異,但整體絕對量值差異不大,并且計算中所得轉捩位置隨雷諾數的變化趨勢與試驗結果基本一致。

圖17 Ma=0.7、α=2°時不同雷諾數下的全機表面摩擦阻力系數云圖

圖18 Ma=0.7、α=2°時預測所得中央翼段上表面轉捩位置與試驗結果的對比
圖19展示了馬赫數0.7、迎角2°狀態不同雷諾數下中央翼段上表面監測點處湍動能沿法向的分布。監測點設置在對稱面上,沿流向距離機翼前緣0.2倍弦長處。橫坐標k代表湍動能的絕對量值,縱坐標代表沿著壁面法向的各個點與壁面的絕對距離h。可以看出,在該監測點處隨著流動雷諾數的增大邊界層的厚度逐漸變薄,但其內部湍動能的峰值卻會明顯增加。因此,雷諾數的增加會使得邊界層內部流動動能和不穩定性均大幅提升,導致流動轉捩更容易發生。

圖19 Ma=0.7、α=2°時不同雷諾數下預測所得中央翼段上表面轉捩位置與試驗結果的對比
圖20給出了全機表面壓力分布截面站位示意,在內、外翼段及平尾上各取多個展向截面,并對所得不同雷諾數下的壓力分布曲線進行對比,分別如圖21和圖22中所示。由圖可以看出,在2°迎角狀態3個不同雷諾數下,內、外翼段表面各個站位截面壓力分布曲線差異相對較小,主要集中在上表面前部和中部區域;并且,雷諾數的差異對于機翼外段壓力分布的影響更為明顯一些。雷諾數的改變對平尾截面壓力分布形態也有一定程度的影響,但從壓力分布結果的對比中未觀察到顯著的規律性。而雷諾數的差異對平尾截面壓力分布形態的影響則稍大一些。

圖20 壓力分布截面位置示意

圖21 Ma=0.7、α=2°時不同雷諾數下的內外翼段展向不同站位截面壓力分布曲線對比

圖22 Ma=0.7、α=2°時不同雷諾數下平尾展向不同站位截面壓力分布曲線對比
鑒于文中變雷諾數CFD計算結果所體現出的雷諾數影響規律基本符合相關物理機理,氣動力及流動狀態差異可以得到合理解釋,同時計算結果與高低速狀態風洞試驗結果之間也具備很好的數據相關性,采用如圖23所示的雷諾數效應修正思路,將CFD計算中獲取的雷諾數差異引起的氣動力系數差量疊加在對應的基礎風洞試驗數據上,從而獲得真實飛行狀態氣動力數據預測值。對于高速狀態試驗迎角范圍比計算迎角范圍小的問題,首先利用試驗迎角范圍之外試驗雷諾數下的CFD計算結果,以試驗迎角區間邊界點上的各個氣動力系數的吹風結果為基準,進行數據平移以將風洞試驗結果向外拓展;而后,在此基礎上通過差量疊加的方式進行雷諾數效應的修正。圖24 和圖25分別給出了低速馬赫數0.2和高速馬赫數0.7狀態下風洞試驗數據修正實例。真實飛行狀態氣動力數據預測值反映出了本文通過數值計算所得到的雷諾數對驗證平臺低速和高速氣動力系數的基本影響規律,具備一定的合理性。

圖23 雷諾數效應修正思路

圖24 Ma=0.2基于自由轉捩CFD計算結果的風洞試驗氣動力系數修正

圖25 Ma=0.7風洞試驗氣動力系數雷諾數修正
采用基于當地流動變量的轉捩預測方法,針對層流機翼飛行驗證平臺低速起降和高速巡航構型分別在試驗和飛行雷諾數下進行數值模擬分析。通過計算結果與風洞試驗數據的對比分析,驗證了結果的可靠性和準確性,歸納總結了雷諾數效應對驗證平臺高低速全機氣動特性、低速大迎角流動分離狀態以及高速層流轉捩特性的具體影響規律和機理如下:
1) 低速飛行狀態,試驗和飛行雷諾數下計算所得升力系數在線性段范圍內均與試驗值吻合良好,在非線性段則較試驗值整體偏高,且變化形態更加和緩;試驗所得最大升力系數與試驗雷諾數下預測值相差不大,但明顯低于飛行雷諾數下的計算值;對于縱向力矩特性,兩種雷諾數下計算所得縱向靜穩定裕度相近,但較試驗結果明顯偏低。
2) 高速飛行狀態,不同雷諾數下升力系數計算結果均與試驗值吻合良好;線性段阻力系數計算結果呈現出隨雷諾數增加而逐漸收斂的趨勢,但其絕對量值均明顯低于試驗值;對于縱向力矩系數,其試驗值整體上與試驗雷諾數下計算結果的吻合程度更好。計算和試驗所得中央翼段表面轉捩位置絕對量值上略有差異,但其隨雷諾數的變化規律基本一致,均隨雷諾數增加逐步前移。
3) 以上結果表明,低速和高速風洞試驗所采用的試驗雷諾數尚未完全進入“自準區[24]”范圍。考慮到該驗證平臺尺寸較小,且試驗與計算氣動力結果在大部分情況下吻合良好,因而在未進行變雷諾數風洞試驗的前提下,利用變雷諾數CFD計算結果輔助進行風洞試驗數據外推修正有助于獲得較為合理的修正結果。
4) 雷諾數的改變會導致流動中慣性和黏性影響的占比以及邊界層特性發生變化。雷諾數對于低速升力和失速特性的影響主要在于其改變了邊界層的厚度和發展速度;而雷諾數對于高速翼面轉捩位置的影響則在于流動中黏性力占比的相對變化改變了邊界層流動對于擾動的抑制能力。
文中基于CFD結果開展了驗證機高低速狀態試驗氣動力數據的修正,預測的飛行狀態氣動力結果反映出了雷諾數效應的具體影響,具備一定的合理性。文中的研究工作也進一步表明CFD計算所得到的規律能夠幫助飛機設計人員在飛行驗證之前更加清晰的了解和認識驗證平臺真實的氣動性能,從而為飛行試驗的設計提供數據支撐。