臨近空間高超聲速稀薄流中梯形縫隙繞流分析

李舟波，方蜀州

(北京理工大學(xué)宇航學(xué)院，北京 100081)

0 引 言

航天器的表面實(shí)際上存在很多不連續(xù)的情況，不同幾何形狀的縫隙經(jīng)常出現(xiàn)在各種氣動(dòng)表面上，如航天器殼體、燃?xì)廨啓C(jī)通道、熱交換器和微電子芯片中的肋狀表面[1]。近年來(lái)，縫隙結(jié)構(gòu)更被廣泛應(yīng)用于超燃沖壓發(fā)動(dòng)機(jī)，以延長(zhǎng)超音速空氣和噴油在燃燒室中的停留時(shí)間[2]。但是，縫隙結(jié)構(gòu)會(huì)嚴(yán)重影響氣流的流動(dòng)狀態(tài)與傳熱特性。這是因?yàn)樵诳p隙入口處由于邊界層分離與再附將產(chǎn)生局部高熱流區(qū)；在某些工況下縫隙入口處會(huì)產(chǎn)生激波；縫隙內(nèi)會(huì)卷起旋渦結(jié)構(gòu)等。同時(shí)，旋渦結(jié)構(gòu)會(huì)向縫隙中卷入大量熱量，使得傳熱特性由熱傳導(dǎo)變?yōu)闊醾鲗?dǎo)加對(duì)流傳熱，縫隙壁面熱流隨之大大增加。因此，高超聲速飛行器表面縫隙的流場(chǎng)結(jié)構(gòu)和熱流分布受到了廣泛的關(guān)注。近年來(lái)，國(guó)內(nèi)外許多研究者對(duì)縫隙進(jìn)行了各種數(shù)值和實(shí)驗(yàn)研究[3-15]。

Sudarshan等[3]研究了圓形和拋物形縫隙的傳熱變化和內(nèi)部激波特性，在馬赫數(shù)分別為 10.9和 10.2 的高超聲速風(fēng)洞中進(jìn)行了實(shí)驗(yàn)。結(jié)果發(fā)現(xiàn)，縫隙結(jié)構(gòu)使機(jī)頭表面的熱通量大幅減少，這對(duì)緩解氣動(dòng)加熱是非常有效的。Dhankarghare等[4]將支板組合成一個(gè)支板縫隙，并在不同馬赫數(shù)和長(zhǎng)徑比的情況下，將其與壁面縫隙進(jìn)行比較。結(jié)果表明，兩種縫隙在馬赫數(shù)和寬高比較低時(shí)都表現(xiàn)出內(nèi)部壓力振蕩，而在高馬赫數(shù)和寬高比時(shí)振蕩減小。Krishna等[5]通過(guò)實(shí)驗(yàn)研究了在超聲速條件下，矩形縫隙后壁面幾何形變的影響。研究表明，縫隙壁面坡道結(jié)構(gòu)使得縫內(nèi)整體聲壓和脈動(dòng)壓力分別降低了10 dB和50%。

國(guó)內(nèi)方面，王世芬等[6]、唐貴明[7]分別在激波風(fēng)洞上開(kāi)展了縫隙熱流測(cè)量實(shí)驗(yàn)，研究了平板上的矩形縫隙和單個(gè)橫縫內(nèi)的熱流分布。張昊元等[8]針對(duì)高超聲速飛行器前緣防熱瓦結(jié)構(gòu)，通過(guò)求解三維可壓縮Navier-Stokes方程進(jìn)行了數(shù)值模擬，研究了縫隙誘導(dǎo)形成的三維旋渦的空間分布特征和旋渦運(yùn)動(dòng)對(duì)物面氣動(dòng)加熱的影響規(guī)律。肖志祥等[9]通過(guò)數(shù)值求解雷諾平均的Navier-Stokes方程組，數(shù)值模擬導(dǎo)彈跨、超聲速流動(dòng)，研究導(dǎo)彈環(huán)向凹槽等典型分離流動(dòng)的被動(dòng)減阻措施對(duì)阻力縮減的影響。通過(guò)分析計(jì)算結(jié)果發(fā)現(xiàn)環(huán)向凹槽有利于降低摩擦阻力，但底阻減小不明顯。陸海波等[10]對(duì)高超聲速飛行器鼻錐使用迎風(fēng)凹腔結(jié)構(gòu)作為熱防護(hù)系統(tǒng)時(shí)，凹腔結(jié)構(gòu)的防熱效能進(jìn)行了數(shù)值研究。結(jié)果表明迎風(fēng)凹腔結(jié)構(gòu)能夠有效對(duì)高超聲速飛行器的鼻錐尤其是駐點(diǎn)區(qū)域進(jìn)行冷卻，凹腔越深，其冷卻效果越好。鼻錐氣動(dòng)加熱的最大熱流并不在尖銳唇緣的頂點(diǎn)，而是位于凹腔內(nèi)的側(cè)壁面上，凹腔的深度(L)變化對(duì)最大熱流的出現(xiàn)位置影響很小。除非凹腔很淺(L/D<0.5)，凹腔底面的熱流值都非常小，基本可以忽略。Wang等[11]通過(guò)風(fēng)洞實(shí)驗(yàn)和數(shù)值計(jì)算，揭示了圓形縫隙可壓縮流動(dòng)的換熱特性。在縫隙后緣附近熱流比無(wú)縫隙模型顯著增加，同時(shí)隨著縫隙直徑(D) 與深度(H) 之比的減小，縫隙內(nèi)的升溫幅度減小。

可以看到，在臨近空間(60～100 km)中的高超聲速縫隙結(jié)構(gòu)的研究仍然很少，這可能是高超聲速條件下可靠的實(shí)驗(yàn)數(shù)據(jù)測(cè)量困難和臨近空間稀薄氣體效應(yīng)的原因。但隨著計(jì)算機(jī)技術(shù)的發(fā)展，尤其是直接模擬蒙特卡洛(Direct simulation Monte Carlo,DSMC) 方法的計(jì)算效率增強(qiáng)，對(duì)稀薄縫隙繞流的研究也越來(lái)越多。Jin等[12]采用DSMC方法對(duì)一組具有圓形分離/再附角的二維矩形縫隙模型進(jìn)行了高超聲速稀薄流動(dòng)數(shù)值模擬。結(jié)果表明，對(duì)于圓形分離角，外部氣體分子更容易進(jìn)入縫隙，從而促進(jìn)了與縫隙內(nèi)氣體分子的動(dòng)量和能量交換。而圓形再附角則會(huì)使得氣體分子從縫隙內(nèi)逸出，削弱動(dòng)量和能量的交換。Guo等[13]采用直接模擬蒙特卡洛法,研究了帶后掠角縫隙的流動(dòng)特性，模擬條件為海拔高度20～60 km,馬赫數(shù) 6～20。可以觀察到，增加高度，縫隙內(nèi)二次回流區(qū)的旋轉(zhuǎn)方向會(huì)發(fā)生逆轉(zhuǎn)。增加馬赫數(shù)，縫隙內(nèi)的速度、壓力和溫度均增加。Palharini等[14]采用DSMC方法研究了80 km高度稀薄過(guò)渡流區(qū)高超聲速縫隙繞流問(wèn)題，探索縫隙幾何尺寸(長(zhǎng)深比)對(duì)縫隙內(nèi)部流場(chǎng)結(jié)構(gòu)、壁面熱流和壓力的影響規(guī)律。結(jié)果表明，在80 km高度稀薄過(guò)渡流區(qū)，幾何尺寸對(duì)高超聲速縫隙繞流的影響與連續(xù)流區(qū)存在不同，主要表現(xiàn)為，當(dāng)長(zhǎng)深比等于4時(shí)縫隙內(nèi)部就出現(xiàn)了雙渦結(jié)構(gòu)和剪切層在縫隙底面再附，而連續(xù)流區(qū)產(chǎn)生類似現(xiàn)象對(duì)應(yīng)的長(zhǎng)深比為14。靳旭紅等[15]采用DSMC方法研究了70 km,75 km,80 km和90 km等4個(gè)飛行高度下稀薄流區(qū)高超聲速縫隙流動(dòng)問(wèn)題。結(jié)果表明：上述飛行高度下,外部流動(dòng)的分離和再附在縫隙內(nèi)部形成一個(gè)充滿腔體的單渦結(jié)構(gòu)；稀薄氣體效應(yīng)對(duì)縫隙內(nèi)部流動(dòng)結(jié)構(gòu)和壁面熱流影響明顯，隨著高度的增加，主渦渦心上移,其形狀逐漸變得“扁長(zhǎng)”，右上角逐漸變尖，熱流越來(lái)越集中分布于縫隙下游側(cè)面的頂部區(qū)域。

從目前的研究成果來(lái)看，上述文獻(xiàn)的研究模型主要為不發(fā)生形變的矩形結(jié)構(gòu)，并且研究?jī)?nèi)容主要集中在流場(chǎng)特性方面，而對(duì)表面氣動(dòng)參數(shù)的研究較少。但是在實(shí)際飛行中，高超聲速流場(chǎng)中縫隙內(nèi)表面會(huì)受到嚴(yán)重的氣動(dòng)加熱，縫隙形狀結(jié)構(gòu)尤其是兩側(cè)壁面會(huì)被強(qiáng)烈的表面熱流燒蝕而發(fā)生變化。因此，對(duì)形變縫隙進(jìn)行流場(chǎng)結(jié)構(gòu)和氣動(dòng)加熱分析顯得十分重要。

本文采用直接蒙特卡洛方法進(jìn)行仿真計(jì)算，將縫隙模型簡(jiǎn)化為壁面可變傾斜角的梯形模型。研究了在稀薄高超聲速環(huán)境中壁面傾角改變對(duì)流場(chǎng)結(jié)構(gòu)和氣動(dòng)熱特性分布的影響規(guī)律，對(duì)飛行器表面優(yōu)化設(shè)計(jì)提供進(jìn)一步的理論幫助。

1 數(shù)值方法

DSMC方法的基本思想是：用有限個(gè)仿真分子代替真實(shí)氣體分子，并在計(jì)算機(jī)中存儲(chǔ)仿真分子的位置坐標(biāo)、速度分量以及內(nèi)能等信息，并在模擬分子運(yùn)動(dòng)、碰撞以及與邊界相互作用的過(guò)程中不斷改變這些信息的具體數(shù)值，最后通過(guò)統(tǒng)計(jì)網(wǎng)格內(nèi)仿真分子的運(yùn)動(dòng)狀態(tài)實(shí)現(xiàn)對(duì)真實(shí)氣體流動(dòng)問(wèn)題的模擬，流場(chǎng)的宏觀參數(shù)可以由統(tǒng)計(jì)平均求得。該方法的關(guān)鍵是在時(shí)間步長(zhǎng)內(nèi)將分子的運(yùn)動(dòng)與碰撞解耦。在時(shí)間步長(zhǎng)內(nèi)讓所有分子運(yùn)動(dòng)一定的距離并考慮在邊界的反射，然后計(jì)算此時(shí)間步長(zhǎng)內(nèi)具有代表性的分子間的碰撞。

本文數(shù)值計(jì)算工具選擇的是美國(guó)桑迪亞國(guó)家實(shí)驗(yàn)室的開(kāi)源SPARTA程序，并且使用北京并行云科學(xué)技術(shù)計(jì)算中心的工作站來(lái)進(jìn)行數(shù)值模擬計(jì)算。模擬分子采用可變徑軟球模型，并考慮到了空氣分子的轉(zhuǎn)動(dòng)能和平動(dòng)能之間的能量傳遞過(guò)程。分子內(nèi)能和平動(dòng)能之間的能量傳遞采用標(biāo)準(zhǔn)Larsen-Borgnakke模型。另外，本文選擇的時(shí)間步長(zhǎng)為當(dāng)?shù)刈杂煞肿悠骄鲎矔r(shí)間的三分之一。模擬分子數(shù)在每個(gè)網(wǎng)格內(nèi)最少為10個(gè)。

同時(shí)，由于采用DSMC方法進(jìn)行數(shù)值計(jì)算時(shí)，模擬分子數(shù)和網(wǎng)格劃分引起的計(jì)算量非常大，需要開(kāi)展并行計(jì)算。一般并行算法采用區(qū)域分解，根據(jù)計(jì)算的處理器數(shù)(或CPU核心數(shù))的多少將計(jì)算區(qū)域劃分為等量的子區(qū)域，每個(gè)處理器在其分配的子區(qū)域內(nèi)部獨(dú)立地計(jì)算模擬分子的碰撞和遷移，離開(kāi)子區(qū)域的模擬分子把攜帶的信息傳遞給對(duì)應(yīng)子區(qū)域的處理器。但DSMC方法的一個(gè)顯著特征是模擬分子在流場(chǎng)中的分布是不斷變化的，計(jì)算時(shí)間依賴于子區(qū)域內(nèi)的分子個(gè)數(shù)、分子的碰撞次數(shù)以及分子與物面的碰撞次數(shù)。因?yàn)槟M前不可能預(yù)測(cè)到不同子區(qū)域內(nèi)的分子個(gè)數(shù)，因此在區(qū)域分解時(shí)將各子區(qū)域的計(jì)算負(fù)載劃分得比較均衡就特別困難。如果簡(jiǎn)單地將流場(chǎng)均勻劃分成若干同等大小的子區(qū)域，就有可能造成處理器之間的負(fù)載分布不均勻，負(fù)載最少的處理器就要花費(fèi)較長(zhǎng)的時(shí)間等待負(fù)載最大的處理器，造成同步等待的時(shí)間顯著增加。這種方式并不適合DSMC并行模擬。

本文發(fā)展了一套靜態(tài)隨機(jī)負(fù)載平衡方法用于解決不同處理器之間的計(jì)算時(shí)間同步問(wèn)題，該方法基于概率近似原理，將計(jì)算區(qū)域的全部網(wǎng)格平均分配給指定的所有處理器。由于采用相同的概率隨機(jī)選取流場(chǎng)網(wǎng)格，按照均分后的數(shù)量分配給每個(gè)處理器，因此當(dāng)計(jì)算區(qū)域的網(wǎng)格數(shù)量較多時(shí)，每個(gè)處理器包含近似相等的物體邊界、高密度流動(dòng)區(qū)域和稀薄氣體區(qū)域的網(wǎng)格數(shù)，這樣每個(gè)處理器的計(jì)算負(fù)載也非常接近。

2 仿真計(jì)算

為了驗(yàn)證SPARTA程序的可行性，首先對(duì)二維的平板模型進(jìn)行了模擬，平板的長(zhǎng)度為100 mm，厚度為5 mm，并將模擬計(jì)算得到的結(jié)果與風(fēng)洞實(shí)驗(yàn)[16]測(cè)量得到的結(jié)果進(jìn)行比較。模擬仿真的來(lái)流條件設(shè)置為風(fēng)洞的實(shí)驗(yàn)條件：來(lái)流氣體為氮?dú)?，自由?lái)流的速度為1504 m/s，平均自由程為1.6 mm，溫度為13.32 K，氣體分子數(shù)密度3.716×1020m-3。

計(jì)算域?qū)挾仍O(shè)置為180 mm，水平方向設(shè)置360個(gè)網(wǎng)格；高度設(shè)置為205 mm，垂直方向設(shè)置410個(gè)網(wǎng)格。每個(gè)網(wǎng)格的長(zhǎng)寬都是0.5 mm。計(jì)算域與平板示意圖如圖1所示。平板壁面溫度為290 K。模擬粒子的仿真比例FNUM設(shè)置為4.56×1011，即一個(gè)模擬粒子代表4.56×1011個(gè)真實(shí)分子。時(shí)間步長(zhǎng)設(shè)置為3.102×10-7s。同時(shí)，采用自適應(yīng)結(jié)構(gòu)網(wǎng)格，分子碰撞采用VHS模型，NTC方法進(jìn)行取樣，能量交換采用Larsen-Borgnakke模型，壁面恒溫，采用漫反射模型。不考慮氣體化學(xué)反應(yīng)。

圖1 平板測(cè)試算例示意圖Fig.1 Calculation example of plate tests

平板上平面的壓力仿真結(jié)果與文獻(xiàn)[17]中DAC方法的仿真結(jié)果以及風(fēng)洞實(shí)驗(yàn)結(jié)果對(duì)比，如圖2所示。SPARTA方法得到的結(jié)果與DAC的結(jié)果比較吻合，但是都比實(shí)驗(yàn)結(jié)果略大，但壓力分布的趨勢(shì)與實(shí)驗(yàn)結(jié)果大致相同。計(jì)算結(jié)果表明本文選擇的數(shù)值方法和計(jì)算程序?qū)τ谶@類高超聲速流場(chǎng)的計(jì)算模擬具有較好的可信度。

圖2 仿真結(jié)果比較Fig.2 Comparison of simulation results

3 模型和流場(chǎng)參數(shù)

3.1 幾何定義

臨近空間航天器上的表面不連續(xù)結(jié)構(gòu)或縫隙遠(yuǎn)遠(yuǎn)小于航天器本身的尺寸，因此將航天器表面上的縫隙結(jié)構(gòu)簡(jiǎn)化為無(wú)限長(zhǎng)的平面上的矩形縫隙。本文僅對(duì)二維情況的縫隙模型進(jìn)行分析，模型的幾何示意圖如圖3所示。

圖3 縫隙模型外形示意圖Fig.3 Schematic diagram of the gap model

傾角θ表示前壁面和底面的夾角，傾角θ分別取30°,45°,60°,75°,90°(即矩形縫隙)。定義縫隙的長(zhǎng)度為L(zhǎng)=9 mm，深度為D=3 mm?？p隙上游的平板長(zhǎng)度為3D，下游平板的長(zhǎng)度為D。

3.2 來(lái)流條件

本文選取的高度為70 km，來(lái)流的速度、壓強(qiáng)、溫度、密度、粒子數(shù)密度以及平均自由程均采用美國(guó)標(biāo)準(zhǔn)大氣(1976) 確定，計(jì)算工況如表1所示。氣體成分為78%的N2和22%的O2。根據(jù)70 km自由來(lái)流的平均自由程和縫隙模型的特征尺寸，得到全局克努森數(shù)Kn大約為0.3，該縫隙模型的繞流問(wèn)題屬于稀薄流范疇。本文設(shè)定縫隙模型的壁面溫度Tw為220 K。

表1 自由流條件Table 1 Free stream conditions

4 流場(chǎng)和壁面參數(shù)計(jì)算公式

4.1 流場(chǎng)特性分析

在本文的研究中主要關(guān)注的流場(chǎng)特性參數(shù)有流場(chǎng)中的馬赫數(shù)、壓力和溫度。DSMC方法中流場(chǎng)的宏觀參數(shù)是由計(jì)算域內(nèi)每個(gè)網(wǎng)格單元中的模擬粒子的參數(shù)計(jì)算平均得到的。

每個(gè)網(wǎng)格單元的質(zhì)心速度矢量的計(jì)算方法為：

(1)

式中：N,mj和cj分別為每個(gè)網(wǎng)格單元中粒子的數(shù)量、質(zhì)量和速度矢量。粒子的熱運(yùn)動(dòng)速度或擾動(dòng)速度c′定義為

c′≡c-c0

(2)

網(wǎng)格中平動(dòng)溫度TT的計(jì)算公式為：

(3)

轉(zhuǎn)動(dòng)溫度TR與網(wǎng)格單元的轉(zhuǎn)動(dòng)能量εR和轉(zhuǎn)動(dòng)自由度ζR有關(guān)，計(jì)算公式為：

(4)

每個(gè)網(wǎng)格單元的壓力的計(jì)算方法為：

(5)

每個(gè)網(wǎng)格單元的密度的計(jì)算方法為：

(6)

4.2 壁面參數(shù)

本文在研究模型的壁面特性時(shí)主要關(guān)注表面壓力系數(shù)和表面熱流系數(shù)。

模型表面的壓力pw由每個(gè)時(shí)間步內(nèi)入射和反射的粒子在物面上的法向動(dòng)量之和得出，表達(dá)式如下:

(7)

式中：A為物面面積；N為單位時(shí)間內(nèi)在單位面積上與物面碰撞的粒子的個(gè)數(shù)；v是j粒子在物面法線法向的速度分量；i,r分別表示入射和反射粒子。

壓力系數(shù)Cp的定義如下：

(8)

熱流qw是通過(guò)計(jì)算碰撞壁面的模擬粒子的凈能量通量來(lái)計(jì)算的，與入射、反射粒子的平動(dòng)能、轉(zhuǎn)動(dòng)能和振動(dòng)能有關(guān)，表達(dá)式如下：

(9)

(10)

式中：eRj和eVj分別代表粒子的轉(zhuǎn)動(dòng)能和振動(dòng)能。

傳熱系數(shù)Ch的定義如下：

(11)

5 流場(chǎng)和壁面參數(shù)計(jì)算公式

5.1 流場(chǎng)特性分析

5.1.1流線分析

為了清晰地描述縫隙的流動(dòng)特征，縫隙內(nèi)回流區(qū)的變化是研究的重點(diǎn)。流線圖能較好反應(yīng)出其變化。圖4顯示了不同側(cè)壁傾角下的流線分布。

圖4 縫隙內(nèi)流線分布圖Fig.4 Distribution of streamline traces inside the gaps

從圖中可以看出，當(dāng)θ=90°時(shí)，縫隙內(nèi)部存在連在一起的兩個(gè)旋渦，此時(shí)仍然認(rèn)為其是一個(gè)回流區(qū)且填滿整個(gè)縫隙，上方流線直接跨過(guò)縫隙。當(dāng)θ=75°時(shí)，隨著傾角的減小，左、右旋渦強(qiáng)度都開(kāi)始減小，在縫隙上方的外部流體逐漸進(jìn)入到縫隙內(nèi)部，縫隙中間的旋渦受到剪切層的擠壓作用而往縫隙底部凹陷。在這兩種工況下，旋渦仍然是連在一起，是上一種情況的延續(xù)，在縫隙內(nèi)表現(xiàn)為單個(gè)旋渦的特性。但兩個(gè)渦心的距離開(kāi)始增大。當(dāng)θ=60°時(shí)，在縫隙前、后側(cè)出現(xiàn)了明顯分開(kāi)的兩個(gè)旋渦，此時(shí)剪切層也深入縫隙底部?？p隙右側(cè)旋渦尺寸變化不大，但縫隙的左側(cè)旋渦明顯變小。θ=45°時(shí)，左側(cè)角落出現(xiàn)非常小的旋渦。對(duì)于右側(cè)旋渦區(qū)域，它在每種工況下，其姿態(tài)和流向長(zhǎng)度變化不大。當(dāng)傾角繼續(xù)減小，θ=30°時(shí)，縫隙的左側(cè)角落幾乎被剪切層充滿，左旋渦徹底消失，在縫隙左側(cè)壁面上邊界層沒(méi)有流動(dòng)分離。

綜合比較圖4，當(dāng)傾角從90°變化到30°時(shí)，縫隙內(nèi)部的回流區(qū)發(fā)生了明顯的變化，并以60°和 30°作為分界點(diǎn)。當(dāng)θ≤60°時(shí)，縫隙內(nèi)部出現(xiàn)了兩個(gè)回流區(qū)，且上方的流體能夠進(jìn)入到縫隙的底面；θ≤30°，左側(cè)的旋渦區(qū)域完全消失。同時(shí)，可以看出，旋渦區(qū)域的旋轉(zhuǎn)方向均為順時(shí)針?？p隙右側(cè)的角部渦流與前壁的傾角無(wú)關(guān)，都保持大致相同的尺寸?？梢酝茢?，改變?chǔ)鹊拇笮?duì)縫隙后壁面區(qū)域的影響較小，這會(huì)在下文表面參數(shù)分析中進(jìn)行驗(yàn)證。

最后，通過(guò)傾斜縫隙后壁面也可以消除右側(cè)的旋渦區(qū)域。為此，以30°的角度傾斜縫隙前后壁面，然后使用DSMC方法計(jì)算。圖5顯示了兩側(cè)壁面傾角都為30°時(shí)的流線分布圖。通過(guò)與圖4中的結(jié)果進(jìn)行比較，發(fā)現(xiàn)當(dāng)后壁面傾角為30°時(shí)，縫隙右側(cè)旋渦區(qū)域也完全消失。因此，一個(gè)重要的結(jié)論是，通過(guò)以適當(dāng)?shù)慕嵌葍A斜縫隙前后壁面，可以完全消除縫隙內(nèi)的旋渦區(qū)域。

圖5 流線分布圖Fig.5 Distribution of streamline traces

5.1.2馬赫數(shù)分析

圖6展示了不同前壁傾角θ下縫隙模型的流場(chǎng)馬赫數(shù)等值線圖。

從圖中可以看出，當(dāng)θ=90°時(shí)，在縫隙內(nèi)的大部分區(qū)域，其馬赫數(shù)小于 1，說(shuō)明在內(nèi)部回流區(qū)的氣流速度很小，且由于壁面的黏性作用，在縫隙中間的流體速度大于兩側(cè)的流體，縫隙底面中間位置受到的侵蝕最大。當(dāng)θ減小時(shí)，縫隙外剪切層從縫隙中部緩慢地深入，導(dǎo)致最大流速區(qū)往前移動(dòng)，亞聲速區(qū)被縫隙外剪切層往下壓縮到縫隙底部，流場(chǎng)內(nèi)的整體速度也開(kāi)始增大。這是由于縫隙前壁面往上游傾斜時(shí)，一方面使縫隙具有更寬的開(kāi)口，一方面讓縫隙入口處流體的再附點(diǎn)從縫隙底面移動(dòng)到縫隙前壁面，外部流動(dòng)可以更深地作用縫隙內(nèi)部流動(dòng)，因此可以看到縫隙內(nèi)流體速度開(kāi)始增加。當(dāng)θ=30°時(shí)，亞聲速區(qū)只存在于縫隙右半?yún)^(qū)域，縫隙外超聲速流場(chǎng)進(jìn)入縫隙左側(cè)區(qū)域。

5.1.3密度云圖

圖7為不同前壁面傾角時(shí)縫隙內(nèi)部的無(wú)量綱密度(當(dāng)?shù)孛芏扰c來(lái)流密度之比)等值線圖。

從圖中可以看出，當(dāng)θ=90°時(shí)，縫隙左側(cè)完全被低密度區(qū)占據(jù)，其值小于來(lái)流密度。又由于氣流在后壁面被壓縮，縫隙右側(cè)流體密度明顯比左側(cè)高很多，密度最高處在縫隙右側(cè)角落。隨著傾角θ的減小，流場(chǎng)內(nèi)的整體密度變化不大，只是由于左側(cè)壁面向前傾斜，縫隙內(nèi)的低密度區(qū)范圍也隨之往上游延伸。而在縫隙右側(cè)，隨著傾角減小，密度只在縫隙右側(cè)角落處略微增大。可以看出前壁面傾角的改變對(duì)流場(chǎng)密度的影響較小。

5.1.4壓力云圖

圖8為不同前壁面傾角下縫隙內(nèi)部的無(wú)量綱壓力(當(dāng)?shù)貕毫εc來(lái)流壓力之比)等值線圖。

從圖中可以看出，當(dāng)θ=90°時(shí)，在縫隙左側(cè)的壓強(qiáng)較低。這是因?yàn)楫?dāng)來(lái)流經(jīng)過(guò)縫隙的前壁面時(shí)，相當(dāng)于流經(jīng)一個(gè)后向臺(tái)階，流場(chǎng)的面積突然增加，自由來(lái)流在這里產(chǎn)生流動(dòng)分離，因此在前壁面頂點(diǎn)處附近出現(xiàn)了局部壓強(qiáng)較低的區(qū)域。隨著前壁傾角的減小，縫隙左側(cè)的低壓區(qū)域被顯著壓縮。隨著縫隙前壁面的傾斜，縫外高壓氣流加速流入，壓力呈上升趨勢(shì)。

在縫隙右側(cè)，還可以觀察到在后壁面頂點(diǎn)的左上方流場(chǎng)處出現(xiàn)了一個(gè)局部高壓區(qū)域。從流線圖中可以看出，這是高速流體撞擊在后壁面頂點(diǎn)區(qū)域所致。隨著θ的減小，該高壓區(qū)域開(kāi)始緩慢增大。這也是因?yàn)榍氨诿鎯A斜，進(jìn)入縫隙內(nèi)的流體速度開(kāi)始增大。其最終會(huì)在右側(cè)壁面開(kāi)始重新被壓縮，并且壓縮效應(yīng)隨速度的增大而增大，所以，前壁傾角的減小會(huì)增大縫隙右側(cè)的壓力。又因?yàn)橄”⌒?yīng)，減弱了縫隙內(nèi)的氣流速度對(duì)氣體壓縮的影響，所以，其改變幅度也并不顯著。

圖6 縫隙內(nèi)馬赫數(shù)分布Fig.6 Mach number distribution in the gaps

最后，還可以明顯發(fā)現(xiàn)前壁傾角的改變對(duì)流場(chǎng)壓強(qiáng)的影響比對(duì)流場(chǎng)密度的影響更大。

5.1.5溫度云圖

圖9展示了不同前壁傾角下縫隙內(nèi)部的無(wú)量綱溫度(當(dāng)?shù)販囟扰c來(lái)流溫度之比)等值線圖。

在縫隙左側(cè)，當(dāng)θ=90°時(shí)，低溫區(qū)占據(jù)縫隙左下方區(qū)域，也正好是左旋渦區(qū)域的范圍。當(dāng)傾角減小時(shí)，縫隙前壁面附近的低溫區(qū)開(kāi)始減小。高溫區(qū)開(kāi)始前移，縫隙內(nèi)的整體溫度也開(kāi)始有所提高。當(dāng)θ=30°時(shí)，隨著左側(cè)旋渦區(qū)域消失，縫隙外的高溫氣流進(jìn)入縫隙的左側(cè)角落，區(qū)域溫度明顯升高，遠(yuǎn)大于縫隙右側(cè)角落。同樣在縫隙右側(cè)上方，由于壓縮作用形成了一個(gè)高溫區(qū)。當(dāng)傾角減小時(shí)，右側(cè)的溫度分布基本沒(méi)有變化，只是縫隙右側(cè)上方形成的高溫區(qū)范圍略微增大。這驗(yàn)證了前壁傾角的改變對(duì)縫隙后壁面附近的影響較小。

最后，為了更深入了解其變化，研究其縫隙內(nèi)部的在稀薄氣體流動(dòng)中的熱不平衡效應(yīng)，對(duì)轉(zhuǎn)動(dòng)溫度也進(jìn)行研究。圖10展示了不同前壁傾角下縫隙內(nèi)部的無(wú)量綱轉(zhuǎn)動(dòng)溫度(當(dāng)?shù)剞D(zhuǎn)動(dòng)溫度與來(lái)流轉(zhuǎn)動(dòng)溫度之比)等值線圖。

綜合比較圖9，可以看出，當(dāng)θ=90°時(shí)，在縫隙左側(cè)，由于整體溫度都較低，平動(dòng)溫度TT與轉(zhuǎn)動(dòng)溫度TR的差距并不明顯，該處的熱不平衡性較小。而在縫隙右側(cè)，由于氣流在后壁面處受到強(qiáng)壓縮，該區(qū)域的平動(dòng)溫度TT遠(yuǎn)大于轉(zhuǎn)動(dòng)溫度TR，存在強(qiáng)烈的熱不平衡效應(yīng)。而當(dāng)傾角減小時(shí)，平動(dòng)溫度和轉(zhuǎn)動(dòng)溫度都開(kāi)始增加，但轉(zhuǎn)動(dòng)溫度增大的幅度更大。當(dāng)θ=30°時(shí)，平動(dòng)溫度和轉(zhuǎn)動(dòng)溫度間的差距小于θ=90°時(shí)兩者間的差值。這說(shuō)明傾角減小能減弱縫隙內(nèi)流場(chǎng)中的熱不平衡效應(yīng)。

圖7 縫隙內(nèi)密度分布Fig.7 Density ratio distribution in the gaps

5.2 表面參數(shù)分析

5.2.1表面熱流系數(shù)分析

圖11分別展示了不同前壁傾角下的縫隙的前壁面、底面、后壁面三個(gè)表面上的熱流系數(shù)分布。

從圖中可以看出，在前壁面，當(dāng)θ=90°時(shí)，熱流系數(shù)從表面頂點(diǎn)往下緩慢增大，在表面中部位置附近達(dá)到最大值。之后，熱流系數(shù)往下再緩慢下降直到表面底端。當(dāng)θ=75°和 60°時(shí)，在表面的上半部分熱流分布與θ=90°區(qū)別不大。在表面下半部分，熱流系數(shù)隨前壁傾角減小而增大。這是因?yàn)?，?dāng)前壁傾角減小時(shí)，其底部附近的流體速度也加快，傳熱速度也隨之增大。當(dāng)θ繼續(xù)減小到 45°和 30°時(shí)，前壁面上整體熱流系數(shù)有了大幅提高。熱流系數(shù)從表面頂點(diǎn)往下單調(diào)增大，并在底部位置達(dá)到最大值。

在縫隙底面上，對(duì)于不同前壁傾角，沿來(lái)流方向，表面熱流系數(shù)都是先增加后減少，因此存在一個(gè)熱流峰值，峰值對(duì)應(yīng)在底面上的x軸坐標(biāo)都是約為 0.8，該位置也是縫隙右側(cè)旋渦的起始位置。同時(shí)，可以看到，隨著前壁傾角的減小，在底面的前半部分，整體熱流系數(shù)也隨之增大。這是因?yàn)殡S著傾角的減小，流體越來(lái)越深入縫隙，造成熱流系數(shù)的整體增加。當(dāng)傾角越小，增加的速度也越大。在底面的后半部分，不同工況下的熱流系數(shù)開(kāi)始趨于定值。

在后壁面上，各工況下的表面熱流系數(shù)均在頂點(diǎn)處達(dá)到最大值，且從頂點(diǎn)沿后壁面向下越來(lái)越小。同時(shí)，隨著傾角的減小，在后壁面的同一位置處，表面熱流系數(shù)雖然逐漸增大，但幅度很小。這說(shuō)明改變前壁面傾角θ，并不會(huì)對(duì)右側(cè)壁面產(chǎn)生較大影響。前壁傾角改變對(duì)縫隙內(nèi)表面的影響主要集中在前壁面以及底面的前半部分。

圖8 縫隙內(nèi)壓強(qiáng)分布Fig.8 Pressure distribution in the gaps

綜合分析三個(gè)壁面上的熱流系數(shù)，熱流系數(shù)最大的區(qū)域是在縫隙的后壁面，這是因?yàn)榱黧w在縫隙后壁面產(chǎn)生壓縮。同時(shí)，隨著θ的減小，前壁面和底面的熱流系數(shù)大幅提高，這會(huì)導(dǎo)致前壁面和底面更容易燒蝕，并進(jìn)一步減小θ值。

5.2.2表面壓力系數(shù)分析

圖12分別展示了不同前壁傾角下的縫隙的前壁面、底面、后壁面三個(gè)表面上的壓力系數(shù)分布。

當(dāng)傾角θ從 90°減小到 45°時(shí)，在前壁面底部的壓力系數(shù)變化不大，但在底部上方的大部分區(qū)域內(nèi)，壓力系數(shù)值隨著前壁傾角的減小而減小。當(dāng)θ=30°時(shí)，在表面上方部分的壓力系數(shù)分布曲線和θ=45°時(shí)基本相同，但在表面底部，壓力系數(shù)遠(yuǎn)大于θ=45°時(shí)的分布。這是因?yàn)楫?dāng)θ=30°時(shí)，縫隙左側(cè)角落旋渦開(kāi)始徹底消失，流體與表面完全貼合，消除了縫隙角落處旋渦的影響，使得壓力系數(shù)在底部突然增大。

在縫隙底面上，各工況的表面壓力系數(shù)都有一個(gè)沿來(lái)流方向持續(xù)增大的趨勢(shì)，其最大值所對(duì)應(yīng)的位置均在底面最右端取得。同時(shí)，底面上整體壓力系數(shù)隨著前壁傾角的減小而增大，不過(guò)增大的幅度很小。

在縫隙后壁面上，壓力系數(shù)的變化規(guī)律與熱流系數(shù)相同。對(duì)于不同的工況，表面壓力系數(shù)在后壁面的頂點(diǎn)處取得最大值，從頂點(diǎn)沿后壁面向下，表面壓力系數(shù)迅速減小。到表面中部區(qū)域時(shí)，減速放緩。

圖9 縫隙內(nèi)溫度分布Fig.9 Temperature distribution in the gaps

圖10 縫隙內(nèi)轉(zhuǎn)動(dòng)溫度分布Fig.10 Rotational temperature distribution in the gaps

圖12 表面壓力系數(shù)分布Fig.12 Pressure coefficient distribution on the surface

并且在表面的底部區(qū)域，表面壓力系數(shù)基本穩(wěn)定。同時(shí),在后壁面的同一位置處，傾角越小，對(duì)應(yīng)的表面壓力系數(shù)越大。這是因?yàn)椋S著θ減小，縫隙內(nèi)氣流速度增大，對(duì)后壁面的壓縮效應(yīng)也隨之增強(qiáng)。

6 結(jié) 論

本文聚焦高超聲速稀薄流航天器表面的縫隙，使用直接模擬蒙特卡洛方法進(jìn)行數(shù)值模擬，對(duì)可變壁面傾角的梯形縫隙的流場(chǎng)情況及壁面參數(shù)進(jìn)行分析，結(jié)果表明：

1) 側(cè)壁傾斜式變形對(duì)縫隙內(nèi)部流場(chǎng)形成了高溫引流和加速導(dǎo)流兩種效應(yīng)。前者可以誘導(dǎo)更多高溫流體更深入縫隙內(nèi)部，流場(chǎng)參數(shù)整體有所提高。后者能夠加快縫內(nèi)氣體的切向流速，使其更快溢流。兩種效應(yīng)都將促使縫隙底部熱流密度升高，而這里是最接近飛行器內(nèi)部零件的位置，因此壁面傾斜形變將給飛行器縫隙內(nèi)的熱防護(hù)帶來(lái)極為嚴(yán)峻的挑戰(zhàn)。

2) 當(dāng)傾角改變時(shí)，縫隙內(nèi)部的回流區(qū)發(fā)生了明顯的變化，并以60°和30°作為分界點(diǎn)。當(dāng)θ≤60°時(shí)，縫隙內(nèi)部從一個(gè)大旋渦徹底分裂為兩個(gè)旋渦；當(dāng)θ≤30°，左側(cè)的旋渦區(qū)域完全消失，在前壁面消除了流動(dòng)分離。同時(shí)，通過(guò)改變縫隙兩側(cè)壁面的角度大小，能完全消除縫隙內(nèi)的旋渦區(qū)域，改善其氣動(dòng)特性。

3) 對(duì)于縫隙內(nèi)的表面而言，前壁面傾角的改變主要影響前壁面和底面前半部分。隨著傾角θ的減小，在前壁面，傳熱系數(shù)顯著增大，而壓強(qiáng)系數(shù)減小。而在底面上，傳熱系數(shù)隨θ減小而增大，壓力系數(shù)也小幅增大。在后壁面上，隨著θ減小，熱流系數(shù)和壓力系數(shù)都有增大，但幅度很小。