馬祎蕾, 余 平, 姚世勇

(空間物理重點實驗室, 北京 100076)

0 引 言

高超聲速飛行器邊界層流場發(fā)生轉(zhuǎn)捩時，飛行器表面熱環(huán)境及飛行器氣動特性發(fā)生顯著變化。正確預(yù)測轉(zhuǎn)捩位置，對于優(yōu)化飛行器熱防護設(shè)計、提升飛行器性能具有十分重要的工程意義[1-2]。

影響高超聲速邊界層轉(zhuǎn)捩的因素眾多，包括轉(zhuǎn)捩雷諾數(shù)、邊界層邊緣馬赫數(shù)、壁面溫度、攻角、頭部鈍度等[1-5]，且作用機理復(fù)雜，許多問題尚未解決。因此，高超聲速邊界層轉(zhuǎn)捩是當前流體力學(xué)領(lǐng)域一個前沿和難點問題，受到廣泛關(guān)注。

壁面溫度是影響高超聲速轉(zhuǎn)捩的一個重要因素。通常情況下，地面風(fēng)洞試驗中壁溫與來流總溫比高，而天上真實飛行條件下的壁溫與來流總溫比低，地面試驗無法模擬天上真實飛行的壁溫條件，這是造成天地差異的一個重要因素。因此，研究壁面溫度條件對邊界層流動穩(wěn)定性及轉(zhuǎn)捩的影響十分重要。

Gary[6]研究表明，轉(zhuǎn)捩雷諾數(shù)隨著壁溫比的降低而增大。Stetson[7]則通過熱線試驗發(fā)現(xiàn)，冷壁的轉(zhuǎn)捩雷諾數(shù)(約3.2×106)低于正常壁溫下的轉(zhuǎn)捩雷諾數(shù)(約4.8×106)。Malik[8]采用eN法分析了5°尖錐的邊界層流場穩(wěn)定性特征，發(fā)現(xiàn)對于絕熱壁邊界層，在來流馬赫數(shù)低于7時為第一模態(tài)導(dǎo)致轉(zhuǎn)捩，而對于冷卻壁，在較低馬赫數(shù)時已變?yōu)榈诙B(tài)導(dǎo)致轉(zhuǎn)捩。Gasperas[9]研究分析了平板和尖錐在不同壁溫下的穩(wěn)定性特征，認為表面溫度沿流向線性降低時，擾動增長率和N值明顯小于絕熱壁。尚慶等[10]研究表明，壁溫升高會使壁面熱流下降，速度邊界層變厚，剖面不飽滿，更易發(fā)生流動分離。蘇彩虹、周恒[11]研究了0°攻角高超聲速鈍錐邊界層，表明絕熱壁邊界層的轉(zhuǎn)捩位置比等溫壁邊界層靠后。曹偉[12]研究發(fā)現(xiàn)，壁面溫度條件對是第一模態(tài)還是第二模態(tài)波決定轉(zhuǎn)捩位置有很大關(guān)系。劉志勇、楊武兵等[13]研究了一定來流工況下的平板邊界層，發(fā)現(xiàn)壁溫升高，初始推遲轉(zhuǎn)捩，隨后變化趨勢發(fā)生反轉(zhuǎn)。Liang等[14]研究了壁溫對馬赫數(shù)7.99的鈍錐邊界層流動穩(wěn)定性影響，研究表明不同壁溫比條件顯著影響了邊界層廣義拐點和壓力梯度分布，導(dǎo)致不同穩(wěn)定性特征，同時在絕熱壁條件下，擾動波增長不再只依賴于第二模態(tài)擾動波。王振清，唐小軍等[15]研究表明，冷壁總體上加速了鈍楔邊界層高頻不穩(wěn)定模態(tài)的發(fā)展。Eppink等[16]通過對比試驗和線性PSE結(jié)果表明，壁溫對橫流行波增長率影響很大，將顯著改變N值分布。

以往關(guān)于高超聲速邊界層轉(zhuǎn)捩壁溫影響的研究多基于二維平板或軸對稱的鈍錐等簡單外形，研究結(jié)果表明了壁面溫度條件對轉(zhuǎn)捩影響的復(fù)雜性。然而，對于采用升力體外形的高超聲速飛行器來說，其三維幾何特征明顯，包括鈍化的端頭和翼前緣、大后掠角的翼面形狀、較扁平的迎風(fēng)面構(gòu)型等。平板鈍三角翼外形具有鈍化的鈍頭、鈍化的翼前緣、大的翼后掠角，有攻角下，其迎風(fēng)面為簡單的平板構(gòu)型。顯然，平板鈍三角翼相比二維平板和軸對稱鈍錐更能表征升力體外形的高超聲速流動，包括橫流效應(yīng)、前緣效應(yīng)和熵層效應(yīng)等，研究平板鈍三角翼外形壁溫條件對邊界層流動穩(wěn)定性及轉(zhuǎn)捩的影響既具有學(xué)術(shù)價值，也具有工程價值。但是，對此類外形的相關(guān)研究工作很少，為此，有必要開展相關(guān)的壁溫影響研究，為高超聲速飛行器的轉(zhuǎn)捩分析與預(yù)示奠定理論基礎(chǔ)。

本文針對平板鈍三角翼外形，圍繞典型高超聲速風(fēng)洞試驗狀態(tài)開展流場數(shù)值計算和線性穩(wěn)定性分析，研究不同壁溫條件下該外形邊界層的流動穩(wěn)定性特征，并采用eN法獲取不同壁溫條件下平板鈍三角翼表面的N值分布，研究壁溫條件對邊界層轉(zhuǎn)捩的影響。

1 計算方法

1.1 控制方程

(1)

式中，

1.2 數(shù)值計算方法

N-S方程數(shù)值求解中，采用有限體積法，對流項采用三階迎風(fēng)Roe-FDS格式離散，黏性項采用中心差分格式離散，時間推進為AF-LU隱式方法。

不考慮對底部流場的模擬，物面采用無滑移邊界條件，出口設(shè)為超聲速外插邊界條件。計算域為半模，網(wǎng)格量為201×281×137，如圖1，法向網(wǎng)格布置點數(shù)較多，原因在于后續(xù)穩(wěn)定性分析需要較為精細的邊界層流場。

圖1 計算域及網(wǎng)格

1.3 穩(wěn)定性計算方法

根據(jù)線性穩(wěn)定性理論，二維局部平行流假設(shè)下，分析O-S方程可將小擾動表示為行進波形式：

(2)

對于空間模式，式(2)中波數(shù)α、β一般為復(fù)數(shù)，其虛部的負值-αi、-βi分別為擾動在x向、z向的空間增長率，ω為擾動波頻率。

擾動向下游傳播時增長率將發(fā)生變化。目前工程上常用的一種轉(zhuǎn)捩預(yù)測方法eN法[17-19]就是基于線性穩(wěn)定性理論，將增長率沿擾動傳播路徑積分,得到擾動幅值放大因子N：

(3)

式中，(x0,z0)為頻率為ω的擾動波開始失穩(wěn)處(αi、βi由正轉(zhuǎn)負的中性點)，其位置在中性曲線上；增長率-αi、-βi取當?shù)仡l率為ω的最不穩(wěn)定擾動波(增長率最大)的增長率，擾動傳播路徑方向即該擾動波的群速度方向。

預(yù)設(shè)轉(zhuǎn)捩判據(jù)為Ntr(通常由經(jīng)驗給出)，即對于不同頻率ω的擾動波，若在(xtr,ztr)處幅值放大因子達到Ntr，預(yù)示轉(zhuǎn)捩發(fā)生。不同頻率擾動波達到Ntr的位置不同，一般取最上游位置作為轉(zhuǎn)捩發(fā)生位置，對應(yīng)N值稱為包絡(luò)N值。轉(zhuǎn)捩位置曲線應(yīng)滿足：

(4)

1.4 算例驗證

文獻[20]針對HIFiRE-1飛行工況進行了數(shù)值模擬，鈍錐半錐角7°，頭部半徑2.5 mm，計算工況如表1，壁面條件取恒溫壁Tw=410 K(Tw/T0,∞=0.168)。

表1 算例工況

為考核三維流場計算的可靠性，特別是邊界層流場和橫流效應(yīng)，本文采用前述計算方法對文獻[20]中算例進行了計算，得到的流向x=850 mm、背風(fēng)面θ=135°位置點的基本流速度剖面與文獻[20]結(jié)果比較如圖2所示，其中橘色速度剖面為計算結(jié)果，與算例基本相符，驗證了邊界層流場計算的正確性。

(a)Streamwise velocity

2 幾何模型與計算工況

幾何模型為一平板鈍三角翼外形[21]，長度取400 mm，鈍度半徑3 mm，后掠角75°。定義坐標系原點O為鈍三角板頂點，x、z軸指向如圖3所示，y軸符合右手坐標系。

圖3 平板鈍三角翼外形及坐標系示意圖

基本的計算工況選取典型的地面高超聲速風(fēng)洞試驗條件：單位來流雷諾數(shù)Re∞=7.12×106/m，馬赫數(shù)M∞=6，總溫T0,∞=420 K，攻角、側(cè)滑角均為0°。選取壁溫條件時，考慮了以下四種典型情況：

1)/真實飛行狀態(tài)下，來流總溫高，壁溫比小，約為0.3；2)在飛行器再入段，由于前期的氣動加熱和隨后的飛行速度降低，有些情況下會出現(xiàn)真實壁溫高于絕熱壁溫的情況，即壁面向流場放熱，此時壁溫比很高，超過絕熱壁溫比；3)高超聲速風(fēng)洞試驗中，來流總溫較低，壁溫比通常為0.7左右；4)在溫敏漆等試驗中，由于采用非金屬材料，壁面條件可看作是絕熱壁。于是，本文選取了四個典型的壁溫條件，包括一個絕熱壁條件和三個等溫壁條件(150 K，300 K，400 K)。不同壁溫條件下，壁溫與來流總溫比如表2所示。

表2 不同壁溫條件下壁溫與來流總溫比

3 計算結(jié)果分析

3.1 基本流分析

為便于比較，本文結(jié)果分析采用如下的無量綱化形式：

V*=V/V，T*=T/T

式中，Lx、Lz分別為鈍三角翼外形的長度和底部半寬度，δ為當?shù)剡吔鐚雍穸龋琕∞為來流速度(860.838 m/s)，T0,∞為來流總溫(420 K)。

圖4為絕熱壁條件下平板三角翼不同流向位置處的馬赫數(shù)分布云圖，由于工況為0°攻角狀態(tài)，迎/背風(fēng)流場對稱，同一流向截面上從側(cè)前緣至中心子午線發(fā)展，邊界層逐漸增厚，并在中心區(qū)域形成流向渦結(jié)構(gòu)。圖5為三角翼表面壓力系數(shù)及壁面極限流線分布，圖6為無量綱最大橫流速度分布，圖中Z/Lz=0對應(yīng)對稱中心線位置，Z/Lz=1對應(yīng)鈍三角翼底部最大展向位置?？芍行淖游缇€流向渦區(qū)域和側(cè)前緣區(qū)域外，翼身上壁面極限流線方向基本平行，無復(fù)雜流場結(jié)構(gòu)。

圖4 不同流向位置處的馬赫數(shù)云分布

圖5 表面壓力系數(shù)分布及流線

由于中心對稱線區(qū)域流動復(fù)雜，存在大尺度流向渦結(jié)構(gòu)，流場變化劇烈，因此，后續(xù)采用線性穩(wěn)定性理論分析時應(yīng)避開該區(qū)域，且平板三角翼外形需考慮三維效應(yīng)，故基于絕熱壁狀態(tài)流場，特別是中心流向渦尺度特征和邊界層橫流流動特征(如圖6)，選擇當?shù)亓骶€近乎平行、橫流效應(yīng)較弱位置Point1(x*=0.88，z*=0.2)，和靠近前緣、橫流效應(yīng)較強點Point2(x*=0.88，z*=0.8)兩點進行研究，分析不同壁面溫度條件對不同模態(tài)流動穩(wěn)定性特征影響。

圖6 橫流速度分布

(a)Point1流向速度

3.2 穩(wěn)定性分析

在獲取各壁溫條件基本流流場的基礎(chǔ)上，針對平行性較好的流動區(qū)域(中心流向渦區(qū)域和側(cè)前緣區(qū)域除外)開展了線性穩(wěn)定性分析計算，獲得了邊界層流動的穩(wěn)定性特征。

(a)Point1

圖9、圖10分別為不同壁面溫度條件下Point1和Point2兩點當?shù)孛總€頻率的最大增長率分布和中性曲線。圖中，低頻和高頻兩個分支分別對應(yīng)第一模態(tài)和第二模態(tài)的擾動。ωr為擾動波的無量綱頻率，βr為擾動波的展向波數(shù)。

(a)Point1

(a)Point1

綜上分析，可知：

2)由于鈍三角翼存在復(fù)雜的三維邊界層流動，不同展向位置的邊界層穩(wěn)定性特征不同，壁溫的影響也出現(xiàn)差異，這與二維平板和0°攻角軸對稱鈍錐不同。Point2(靠近前緣區(qū)域)的橫流模態(tài)擾動增長率與第二模態(tài)增長率幅值相當，且當?shù)財_動增長率明顯高于Point1(靠近對稱中心線)，并先于Point1出現(xiàn)模態(tài)相交現(xiàn)象。

3.3 包絡(luò)N值及轉(zhuǎn)捩分析

四種壁溫條件下整體包絡(luò)N值在鈍三角翼表面的分布情況見圖12。

圖13 轉(zhuǎn)捩判據(jù)為Ntr=4.5時轉(zhuǎn)捩線位置

為探討轉(zhuǎn)捩反轉(zhuǎn)現(xiàn)象的內(nèi)在機理，下面分析過Point1點最不穩(wěn)定擾動波沿勢流方向的變化情況。

首先提取不同壁溫條件下Point1點處的包絡(luò)N值及對應(yīng)該包絡(luò)N值的擾動波頻率(本文稱之為Point1點的最不穩(wěn)定波頻率)和擾動增長率，如表3所示。

表3 Point1點當?shù)刈畈环€(wěn)定擾動波參數(shù)

過Point1點對應(yīng)最不穩(wěn)定波頻率的擾動波其增長率沿勢流方向的分布如圖14(a)所示，橫軸ls*為無量綱勢流方向距離。圖14(b)、圖14(c)分別為流向ls*≈0.0065位置(靠近擾動失穩(wěn)起始點)和Point1點擾動增長率的頻域分布，四個圓點分別為表3中提取的Point1點的最不穩(wěn)定波頻率對應(yīng)增長率。

由圖14(a)可知，對應(yīng)Point1點最不穩(wěn)定波頻率的擾動增長率沿勢流方向明顯具有兩個峰值，ls*=0.0065點和Point1點則分別位于前后這兩個峰值區(qū)域內(nèi)。

(a)增長率分布

圖15 N值分布

由以上分析可知，轉(zhuǎn)捩反轉(zhuǎn)的內(nèi)在機理在于大后掠角平板鈍三角翼邊界層流動獨特的橫流特征：在靠近前緣區(qū)域，橫流效應(yīng)強，擾動以橫流模態(tài)和第一模態(tài)為主，形成一個擾動增長率峰值區(qū)域；而在靠近對稱中心線的區(qū)域，橫流效應(yīng)弱，擾動以第二模態(tài)為主，形成另一個擾動增長率峰值區(qū)域。壁溫比低于絕熱壁值時，壁溫比的增加主要影響第二個峰值區(qū)域，導(dǎo)致N值降低，轉(zhuǎn)捩延遲；當壁溫比超過絕熱壁值后，壁溫比對第一個峰值區(qū)域的影響超過了對第二個峰值區(qū)域的影響，導(dǎo)致N值增大，轉(zhuǎn)捩前移。

4 結(jié) 論

本文研究了在地面風(fēng)洞試驗條件(M∞=6，Re∞=7.12×106/m，α=β=0°)下，壁面溫度變化對大后掠角平板鈍三角翼外形高超聲速邊界層流動穩(wěn)定性及轉(zhuǎn)捩的影響。所得結(jié)論如下：

2)壁溫條件的變化對不同模態(tài)的擾動產(chǎn)生顯著影響。壁溫比升高促進橫流模態(tài)和第一模態(tài)擾動增長，同時第二模態(tài)的擾動受到抑制。

3)隨著壁溫比的增加，受第一、第二模態(tài)聯(lián)合作用影響，鈍三角翼表面N值分布呈現(xiàn)先減小后增大的特點，預(yù)示存在轉(zhuǎn)捩反轉(zhuǎn)現(xiàn)象，反轉(zhuǎn)點壁溫比大約在絕熱壁壁溫比(約0.8)附近。

高超聲速靜風(fēng)洞轉(zhuǎn)捩試驗中，來流總溫低，壁溫比的量值高(約0.7)，不少試驗采用非金屬模型，壁面可視為絕熱壁。而在實際的高超聲速飛行試驗中，來流總溫高，壁溫比的量值低(量級大致在0.3)。本文的研究結(jié)果表明，這兩種狀態(tài)的穩(wěn)定性特征差異大，這預(yù)示著可能存在顯著的天地差異。顯然，直接基于地面的靜風(fēng)洞試驗結(jié)果來預(yù)示飛行條件下的轉(zhuǎn)捩可能是不合適的，需要開展相關(guān)的天地差異分析和天地相關(guān)性研究。