馬祎蕾, 余 平, 姚世勇

(空間物理重點(diǎn)實(shí)驗(yàn)室, 北京 100076)

0 引 言

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

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

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

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

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

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

1 計(jì)算方法

1.1 控制方程

(1)

式中，

1.2 數(shù)值計(jì)算方法

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

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

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

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

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

(2)

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

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

(3)

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

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

(4)

1.4 算例驗(yàn)證

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

表1 算例工況

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

(a)Streamwise velocity

2 幾何模型與計(jì)算工況

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

圖3 平板鈍三角翼外形及坐標(biāo)系示意圖

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

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

表2 不同壁溫條件下壁溫與來(lái)流總溫比

3 計(jì)算結(jié)果分析

3.1 基本流分析

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

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

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

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

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

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

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

圖6 橫流速度分布

(a)Point1流向速度

3.2 穩(wěn)定性分析

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

(a)Point1

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

(a)Point1

(a)Point1

綜上分析，可知：

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

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

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

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

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

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

表3 Point1點(diǎn)當(dāng)?shù)刈畈环€(wěn)定擾動(dòng)波參數(shù)

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

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

(a)增長(zhǎng)率分布

圖15 N值分布

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

4 結(jié) 論

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

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

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

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