孔祥志，韓宇峰*

天津大學(xué) 高速空氣動力學(xué)研究室，天津 300072

一直以來，氣動加熱在高超聲速飛行器的設(shè)計(jì)中占據(jù)著關(guān)鍵的地位，是工程實(shí)際及理論研究都關(guān)注的問題。飛行器表面的熱防護(hù)材料在經(jīng)過燒蝕后，會發(fā)生一定程度的變形，形成粗糙表面，對流動和氣動加熱產(chǎn)生不可忽略的影響。

在實(shí)際問題中，粗糙壁面形式多樣，開展對粗糙壁面峰值熱流特性的研究，首先要建立描述粗 糙 壁 面 特 征 的 模 型。Colebrook 和White［1］的研究發(fā)現(xiàn)，粗糙密度對于表征粗糙特征是至關(guān)重要的，后來Schlichting［2］量化表征了粗糙密度，用λ=Fr/F 定義粗糙密度，F(xiàn)r表示粗糙壁正面總投影面積，F(xiàn) 表示粗糙與壁面平行的投影面積；Nikuradse［3］基于粗糙密度提出等效砂礫高度的概念，通過試驗(yàn)建立起等效沙礫高度與粗糙壁面特征的關(guān)系；Sigal 和Danberg［4］擴(kuò)展 并改進(jìn)了這一映射關(guān)系，使用Sigal-Danberg 參數(shù)定義了粗糙密度，使等效砂礫高度獲得更為廣泛的應(yīng)用。此外，為了表征粗糙壁面特征，Napoli 等［5］提出有效斜率的概念衡量粗糙壁面的高度及密度，用于對完全粗糙壁面特征的表征，同時研究了有效斜率與粗糙壁面摩阻和熱流的關(guān)系。上述粗糙特征模型由于提出角度的不同，適用于不同壁面粗糙情況［6］。

對于高超聲速粗糙壁的氣動熱規(guī)律，從20 世紀(jì)70 年代以來開展了大量的實(shí)驗(yàn)研究。Bertran等［7］在Jaeck［8］的試驗(yàn)報告基礎(chǔ)上，研究平板波紋型粗糙壁的傳熱特性，重點(diǎn)研究第1 波峰值熱流的演變規(guī)律，發(fā)現(xiàn)了粗糙高度、寬度、間隔等對于熱流的影響規(guī)律；美國PANT 計(jì)劃、文獻(xiàn)［8-14］利用實(shí)驗(yàn)得到的數(shù)據(jù)，擬合得到若干計(jì)算粗糙壁面熱流的經(jīng)驗(yàn)預(yù)測公式，受限于實(shí)驗(yàn)數(shù)據(jù)以及實(shí)驗(yàn)條件，這些公式不具有很大的通用性；Nestler［15］綜述了表面粗糙對傳熱影響的實(shí)驗(yàn)數(shù)據(jù)以及預(yù)測方法，這些表面粗糙的類型包括臺階、空穴、間隙和突起等，通過對該領(lǐng)域關(guān)鍵技術(shù)的評估，指出以往對于粗糙壁氣動熱特性的研究缺乏對傳熱過程的物理認(rèn)知，建立的經(jīng)驗(yàn)相關(guān)性公式是不充分的，需要更完善的實(shí)驗(yàn)數(shù)據(jù)以及新的預(yù)測方法。

除了實(shí)驗(yàn)，學(xué)者們也采用各種數(shù)值方法對高超聲速粗糙壁的熱流進(jìn)行了研究。Polak 等［16］使用隱式有限差分技術(shù)研究不同馬赫數(shù)、雷諾數(shù)下的波紋型粗糙壁面的熱流演變規(guī)律，并對Jaeck等［8］提出的峰值熱流預(yù)測公式進(jìn)行評估，發(fā)現(xiàn)預(yù)測公式對試驗(yàn)結(jié)果存在嚴(yán)重依賴；Egorov 等［17］數(shù)值研究波紋粗糙平板壁面，得出壁面熱流的變化規(guī)律，指出在粗糙壁面峰值附近，熱流會達(dá)到最大；李俊紅等［18］用數(shù)值模擬和理論方法對高超聲速可壓縮湍流中的粗糙壁面熱增量進(jìn)行研究，討論了不同粗糙高度、密度等對熱流的影響規(guī)律，比較了幾種預(yù)測公式［12-14］的結(jié)果與數(shù)值結(jié)果的差異；李玉川和鮑麟［19］研究波紋粗糙壁面上高速氣流的流動和傳熱特性，結(jié)果表明在波狀壁流動中，壁面熱流系數(shù)與總阻力系數(shù)呈現(xiàn)出直接關(guān)聯(lián)。Kumar 和Reddy［20］在 平 板上布 置 三 維 凸起，通過實(shí)驗(yàn)研究凸起附近的熱流，使用不同工況下的實(shí)驗(yàn)數(shù)據(jù)擬合得到用于預(yù)測凸起附近熱流的關(guān)系式；Estruch-Samper［21］在平板上布置短壓縮斜坡粗糙元，研究粗糙偏轉(zhuǎn)角對于峰值熱流等的影響，評估幾種峰值熱流預(yù)測公式的有限適用性；Avallone 等［22］通過實(shí)驗(yàn)研究圓柱形粗糙元的流動以及熱流峰值特征；Thakkar 等［23］使用多組不規(guī)則粗糙表面的模擬數(shù)據(jù)，開發(fā)了具有多個參數(shù)的粗糙壁面特征優(yōu)化模型，但是這一模型的參數(shù)較多，實(shí)際使用起來較為復(fù)雜。

總之，目前對于粗糙壁面的氣動熱有了深入的認(rèn)知，但多數(shù)研究集中在再入飛行背景下的大尺度（毫米-厘米量級）燒蝕。隨著新型航天航空飛行器的發(fā)展以及其表面精細(xì)化的熱防護(hù)需求，復(fù)合燒蝕材料由于輕量以及高密度等優(yōu)勢，應(yīng)用越來越廣泛，相比于傳統(tǒng)燒蝕材料，復(fù)合燒蝕材料的線燒蝕率更低。例如在10～12 馬赫數(shù)范圍的高超聲速巡航或滑翔飛行器上，機(jī)翼結(jié)構(gòu)采用瓦狀非承重外部隔熱系統(tǒng)，對于低重量的隔熱瓦，會產(chǎn)生百微米量級的波紋型表面粗糙［24］。毫米以下的微粗糙和毫米、厘米量級的大粗糙的流場結(jié)構(gòu)是有顯著差異的。微粗糙尺度遠(yuǎn)小于邊界層厚度的，流場結(jié)構(gòu)較為簡單，壁面熱流增量主要是由于壁面剪切耗散增大而導(dǎo)致［16］；大粗糙尺度與邊界層厚度相當(dāng)，流場結(jié)構(gòu)更為復(fù)雜，粗糙附近壓力波動較大并出現(xiàn)流動分離，熱流會受分離渦影響，熱增量機(jī)理更為復(fù)雜。目前已有的預(yù)測方法，都是針對大粗糙壁的基于實(shí)驗(yàn)擬合得到的半經(jīng)驗(yàn)平均熱流預(yù)測公式。對于微粗糙度的峰值熱流，這些方法是不適用的。此外在這種微燒蝕情況下的粗糙壁面峰值熱流的預(yù)測是很稀缺的［16］。

在實(shí)際中，飛行器表面因燒蝕而產(chǎn)生的粗糙形狀是復(fù)雜的，但是等效砂礫粗糙度概念的提出啟示我們可以利用一種標(biāo)準(zhǔn)粗糙元模型等效地模擬實(shí)際粗糙壁面形狀以獲得壁面熱流變化的規(guī)律。

本文重點(diǎn)關(guān)注微燒蝕條件下出現(xiàn)的百微米量級的波紋型粗糙壁，結(jié)合目前對于粗糙壁熱流影響因素的認(rèn)識以及預(yù)測方法，使用直接數(shù)值模擬計(jì)算不同高度、流動參數(shù)下的流場及熱流。在此基礎(chǔ)上，進(jìn)行理論分析，總結(jié)數(shù)值結(jié)果的規(guī)律，提出定量預(yù)測粗糙壁面峰值熱流的模型。

1 物理模型及數(shù)值方法

1. 1 控制方程

通過直接數(shù)值求解直角坐標(biāo)系下的可壓縮守恒型Navier-Stokes（N-S）方程獲得平板邊界層的流場，其控制方程為

式中：t 為時間變量；x、y 分別為流向和法向的空間坐標(biāo)；U 為守恒型變量；E、F 分別為流向、法向的通量；Ev、Fv為相應(yīng)的黏性項(xiàng)，各項(xiàng)的具體表達(dá)式分別為

其中：u、v 分別為流向、法向的速度分量；p、ρ 和T分別表示壓力、密度和溫度；μ 為黏性系數(shù)；κ 為熱傳導(dǎo)系數(shù)；黏性系數(shù)使用Sutherland 公式計(jì)算：

式中：Tμ為參考溫度，取110.4 K。

熱傳導(dǎo)系數(shù)在普朗特?cái)?shù)以及定壓比熱為常數(shù)的條件下與黏性系數(shù)為線性關(guān)系：

式中：Tκ為參考溫度，取110.4 K。

在計(jì)算過程中需要選擇適當(dāng)?shù)奶卣鞒叨葘Ψ匠踢M(jìn)行無量綱化。上標(biāo)*號表示有量綱物理量，無上標(biāo)*號的量表示無量綱物理量。選取特征長度L*對空間坐標(biāo)進(jìn)行無量綱化：

雷諾數(shù)、馬赫數(shù)分別為

無量綱形式的黏性系數(shù)和熱傳導(dǎo)系數(shù)分別為

1. 2 數(shù)值格式和程序驗(yàn)證

壁面采用無滑移邊界條件，上邊界給定來流條件，入口采用自由來流條件，出口采用線性外推。采用基于有限差分的直接數(shù)值模擬方法開展流場計(jì)算，計(jì)算對流項(xiàng)采用三階WENO（Weighted Essentially Non-Oscillatory）格式，黏性項(xiàng)采用四階中心格式，通量分裂選擇Steger-Warming 分裂，采用四階龍格-庫塔方法在時間方向推進(jìn)求解。

為了驗(yàn)證數(shù)值程序的正確性，與Zhao 等［25］關(guān)于單個正弦粗糙元的高超聲速流動研究中的工況進(jìn)行比對，圖1 展示了流場的壁面壓強(qiáng)分布和粗糙元附近的速度剖面的對比結(jié)果，可以看到，在粗糙元附近，計(jì)算結(jié)果基本吻合，驗(yàn)證了所使用程序在計(jì)算粗糙元時的正確性。

圖1 程序驗(yàn)證Fig.1 Program verification

1. 3 物理模型、網(wǎng)格劃分策略及計(jì)算參數(shù)

數(shù)值模擬的物理模型為全場均布正弦型波紋的平板，如圖2（a）所示，壁面連續(xù)分布有100 個周期的正弦凸起。壁面函數(shù)表達(dá)式為

圖2 正弦型壁面粗糙分布及熱流驗(yàn)證Fig.2 Sinusoidal wall roughness distribution and heat flux verification

式中：h0、kn和L 分別為粗糙元的高度、個數(shù)和寬度。近壁面使用貼體網(wǎng)格。同時為了使粗糙元以外的部分不受到網(wǎng)格彎曲帶來的影響，在距離粗糙元一定高度以上，仍采用平行網(wǎng)格，二者光滑過渡；為了準(zhǔn)確計(jì)算熱流，法向網(wǎng)格尺度使用基于分子平均自由程的網(wǎng)格準(zhǔn)則［26-27］得到，法向第1 層網(wǎng)格尺度在10-6m。在每個粗糙單元尺度上流向分布20 多個點(diǎn)，法向分布100 多個點(diǎn)，計(jì)算域總網(wǎng)格量為100 萬左右，具體網(wǎng)格分布如圖2（a）所示，圖中右上角為x、y 坐標(biāo)比例為1 的示意圖。不同法向第1 層長度（10-7～10-5m）下熱流的計(jì)算結(jié)果見圖2（b），圖中y1表示法向第1 層網(wǎng)格長度，可以看出，10-5m 的結(jié)果偏差很大，5×10-6m 的結(jié)果趨于收斂，但仍有細(xì)微的偏差。而在長度為10-6、5×10-7、10-7m 時，熱流結(jié)果完全一致，說明在法向第1 層網(wǎng)格尺度10-6m下可保證熱流的計(jì)算是準(zhǔn)確的。后文計(jì)算中每個算例都經(jīng)過了網(wǎng)格無關(guān)性檢驗(yàn)，以保證計(jì)算結(jié)果的準(zhǔn)確性。

使用無量綱的熱流系數(shù)ch 來表征壁面熱流：

式中：chsm表示對應(yīng)光滑壁的熱流系數(shù)；計(jì)算的來流參數(shù)選自30 km 的大氣，特征長度選取1 m，具體參數(shù)見表1。

表1 計(jì)算參數(shù)Table 1 Calculation parameters

2 數(shù)值結(jié)果及分析

粗糙壁的幾何特性會影響壁面熱流的演變規(guī)律，Sigal［4］、Bertram［7］、Polak［16］、李俊紅［18］等對粗糙壁高度、寬度、形狀等幾何特征對壁面熱流的影響進(jìn)行了研究，結(jié)果表明，高度和寬度是較為關(guān)鍵的影響因素，兩者通常以比值形式的關(guān)系［4］體現(xiàn)這種影響規(guī)律。

為了研究不同幾何特征以及來流參數(shù)對于粗糙壁峰值熱流的影響規(guī)律，數(shù)值計(jì)算了以下4種工況的流場和熱流，分別是：

1） 來流參數(shù)不變，粗糙壁寬度不變，不同高度下的粗糙壁。

2） 來流參數(shù)不變，粗糙壁高度不變，不同寬度下的粗糙壁。

3） 來流參數(shù)不變，高度寬度均變化，固定寬高比的粗糙壁。

4） 粗糙壁高度和寬度不變，不同馬赫數(shù)以及壁溫下的粗糙壁。

根據(jù)以上4 種情況，在圖3 展示了具體計(jì)算工況的高度和寬度分布，圖3（a）列出了前3 種情況的粗糙壁高度和寬度分布，圖3（b）列出了第4 種情況下的工況分布。共計(jì)算了28 個不同工況下的粗糙壁流場和熱流，使用的網(wǎng)格流向與法向均已無關(guān)。同時，由于平板頭部存在黏性干擾，主要分析0.2～1 m 的流場。在圖3（c）與圖3（d）中給出了高度0.421 2 mm、馬赫數(shù)為6 時的計(jì)算結(jié)果的整體以及局部區(qū)域的速度云圖。

圖3 計(jì)算工況以及典型計(jì)算結(jié)果Fig.3 Calculation conditions and typical calculation results

2. 1 不同高度下的壁面熱流

首先分析第1 類工況下的壁面熱流，即來流參數(shù)不變，粗糙壁寬度不變，計(jì)算不同高度下粗糙壁的基本流場。具體的粗糙元配置：寬度均為10.53 mm，高 度 分 別 為0、0.105 3、0.126 36、0.210 6、0.294 84、0.379 08、0.421 2、0.526 5 mm。

圖4 給出了不同高度下，全場均布粗糙的壁面熱流和光滑平板的對比，同時將壁面形狀給出，方便與壁面熱流做比對。從圖4（a）可以發(fā)現(xiàn)，粗糙壁面的熱流隨壁面形式呈周期性變化，在正弦壁峰附近，熱流達(dá)到峰值；同一流向位置的壁面峰值熱流，隨著高度的增加而增大，這一熱增量不可忽略。由于粗糙元高度相比邊界層厚度很小，因此熱流的變化規(guī)律和壁面幾何形狀類似。圖4（b）給出了局部放大的結(jié)果。圖4（c）給出了較大高度（0.421 2 mm）下平均熱流與光滑壁熱流的比較，粗糙壁平均熱流是對波紋粗糙元一個周期內(nèi)的熱流進(jìn)行平均求得的，而光滑壁面的平均熱流指相同流向區(qū)域的熱流平均值?？梢娫谖⒊叨鹊拇植谙?，壁面平均熱流增量是很小的，即使在較大高度下，微粗糙下的平均熱流增量與光滑壁平均熱流的比值小于3%，因此本文重點(diǎn)關(guān)注峰值熱流。

圖4 不同高度下的壁面熱流Fig.4 Surface heat flux at different heights

為了分析熱流隨著高度的變化規(guī)律，將不同流向位置的壁面峰值熱流提取出來進(jìn)行分析比較，對于每一個工況可以得到100 個峰值點(diǎn)，不同高度下的壁面峰值熱流系數(shù)及相對熱流增量Qratio如圖5 所示。圖5（a）給出了不同高度下的壁面峰值熱流系數(shù)，可以發(fā)現(xiàn)隨著流向位置x 的增大，壁面峰值熱流減小；圖5（b）則給出了不同高度下的相對峰值熱流增量，很明顯，沿流向的相對熱流增量近似保持不變。同時，在同一流向位置，隨著高度增加，壁面峰值熱流和相對峰值熱流增量均增大。壁面峰值熱流沿流向位置的趨勢與反比例函數(shù)的形式近似，隨高度增大的趨勢與線性函數(shù)近似，相對峰值熱流增量與流向位置無關(guān)。為了確定壁面峰值熱流與流向位置的具體關(guān)系，圖6 給出了不同高度下壁面峰值熱流與x-0.5的關(guān)系。從圖中可以看出ch與x-0.5成線性相關(guān)關(guān)系，這一關(guān)系與光滑平板壁面熱流和x-0.5的線性規(guī)律一致。

圖5 不同高度下壁面峰值熱流及相對峰值熱流增量Fig.5 Surface peak heat flux and relative surface peak heat flux increment at different heights

圖6 不同高度下壁面峰值熱流與x-0.5 的關(guān)系Fig.6 Relationship between surface peak heat flux and x-0.5 at different heights

圖7給出了3 個不同流向位置下，壁面峰值熱流隨壁面高度變化的關(guān)系，很明顯，ch與h 具有較強(qiáng)的線性相關(guān)性。因此，對于粗糙壁面峰值熱流的計(jì)算，有理由做出以下假設(shè)：

圖7 不同流向位置壁面峰值熱流與高度的關(guān)系Fig.7 Relationship between surface peak heat flux and height at different flow direction positions

式中：第1 項(xiàng)表示由于粗糙引起的熱流增量，第2 項(xiàng)表示光滑平板的熱流。由式（6）得到相對峰值熱流增量為

式（8）與流向位置x 是無關(guān)的，符合圖5（b）的規(guī)律。

2. 2 不同寬度下的壁面熱流

對于第2 類工況，保持壁面高度不變，計(jì)算不同粗糙寬度的流場。使用的具體粗糙單元配置為：高度為0.421 2 mm，寬度分別為1.755、3.51、5.265、7.02、10.53 mm。對于小寬度，以0.4 m 為中心分布了100 個粗糙，在圖中對應(yīng)的曲線為這100 個粗糙的所在。

圖8 顯示了不同寬度下，粗糙壁面的熱流結(jié)果與光滑壁面（x=0.4～0.41 m）的比較?？梢钥闯觯?dāng)粗糙帶高度固定時，寬度越大，壁面峰值熱流越小。

圖8 不同寬度下的壁面熱流Fig.8 Surface heat flux at different widths

圖9 給出了粗糙壁峰值熱流及相對峰值熱流增量沿流向的演變規(guī)律?？梢钥吹剑S著粗糙寬度的增加，壁面峰值熱流逐漸減小，但是整體并不成線性關(guān)系；同時，相對峰值熱流增量也隨寬度增加而減小。

圖9 不同寬度下的壁面峰值熱流及相對峰值熱流增量Fig.9 Surface peak heat flux and relative surface peak heat flux increment at different widths

2. 3 固定寬高比下的壁面熱流

在得到上述關(guān)于熱流與寬度成負(fù)相關(guān)的規(guī)律后，需要進(jìn)一步確定壁面熱流與寬度的定量關(guān)系。文獻(xiàn)中常用Sigal-Danberg 參數(shù)來表征粗糙帶特征，這一參數(shù)以高度與寬度的比值來定量刻畫粗糙度特征。本節(jié)使用第3 類工況，目的是研究高度和寬度對粗糙壁面峰值熱流的影響強(qiáng)度以及這兩個因素之間比值規(guī)律，具體的粗糙壁參數(shù)見表2。

表2 固定寬高比的計(jì)算參數(shù)Table 2 Fixed aspect ratio calculation parameters

所取工況的粗糙寬高比均固定為25，假如寬度和高度對壁面熱流有相同的影響強(qiáng)度，那么在固定寬高比的情況下，幾種不同配置的粗糙壁面峰值熱流在同一流向位置應(yīng)是相同的。圖10（a）給出了5 個工況下的壁面熱流。可以看到，在固定寬高比的情況下，粗糙壁面峰值熱流并不相同。圖10（b）給出了每個工況的壁面峰值熱流隨流向x 的變化，可以發(fā)現(xiàn)高度越大，壁面峰值熱流越大。這說明熱流受高度的影響更大。同時，圖10（c）進(jìn)一步顯示了相對峰值熱流增量與流向位置的無關(guān)性。通過上述分析可以認(rèn)識到，壁面峰值熱流受高度的影響大于受寬度的影響，熱流與寬度的冪次關(guān)系應(yīng)該在-1～0 這一范圍。

圖10 相同寬高比下的壁面熱流Fig.10 Surface heat flux with the same aspect ratio

2. 4 不同馬赫數(shù)下的壁面熱流

將粗糙壁面的寬度取為10.53 mm，高度取為0.379 08 mm，分析了不同馬赫數(shù)（4、5、6、7、8、10）下的壁面峰值熱流。圖11（a）給出了粗糙壁面峰值熱流結(jié)果的比較?？梢园l(fā)現(xiàn)，同一流向位置，在不同馬赫數(shù)下，壁面峰值熱流隨馬赫數(shù)的變化并非單調(diào)。這是由于光滑壁的熱流本身隨馬赫數(shù)的變化就是非單調(diào)的。圖11（b）給出了不同馬赫數(shù)下，相對峰值熱流增量的演變規(guī)律，可以發(fā)現(xiàn)同一流向位置，在不同馬赫數(shù)下，相對峰值熱流增量隨馬赫數(shù)增大而單調(diào)減小，但與流向位置x 無關(guān)，這一規(guī)律與前文相同。

圖11 不同馬赫數(shù)下壁面峰值熱流及其相對峰值熱流增量Fig.11 Surface peak heat flux and relative surface peak heat flux increment at different Mach numbers

2. 5 不同壁溫下的壁面熱流

同樣，研究壁溫對于粗糙壁面熱流的影響時，固定粗糙壁面的為寬度10.53 mm、高度0.379 08 mm，計(jì)算不同壁溫（350、400、500、600、650、700、800 K）下的粗糙壁熱流。

圖12 顯示了不同壁溫下，全場均布粗糙的壁面峰值熱流和相對峰值熱流增量，可以發(fā)現(xiàn)，在同一流向位置，壁溫越高，壁面峰值熱流和相對峰值熱流越小；相對峰值熱流增量近似與流向位置無關(guān)。

圖12 不同壁溫下的壁面峰值熱流及相對峰值熱流增量Fig.12 Surface peak heat flux and relative surface peak heat flux increment at different wall temperatures

3 熱流預(yù)測模型

第2 節(jié)中不同高度、不同寬度、相同寬高比、不同馬赫數(shù)以及不同壁溫下，粗糙壁面峰值熱流的定性變化規(guī)律可以概括為

1） 壁面峰值熱流的變化與粗糙壁的高度成正比，與流向位置x 的0.5 次成反比；相對峰值熱流增量與流向位置近似無關(guān)。

2） 壁面峰值熱流的變化與粗糙壁的寬度成反比。

3） 壁面峰值熱流受高度的影響大于受寬度的影響。

4） 相對峰值熱流增量隨馬赫數(shù)增高而減小。

5） 相對峰值熱流增量隨壁溫增高而減小。

在以上定性規(guī)律認(rèn)知的基礎(chǔ)上，為了進(jìn)一步分析得到定量預(yù)測粗糙壁面峰值熱流的模型，首先從光滑平板壁面熱流的理論解出發(fā)給出熱流預(yù)測公式的形式。

對于可壓縮平板邊界層相似性解，存在理論解用來計(jì)算熱流，其具體推導(dǎo)過程見文獻(xiàn)［28］。對于熱流系數(shù)ch，理論解為

式（10）即光滑平板壁面熱流的理論解。在之前的數(shù)值分析中，對粗糙壁面峰值熱流做了的假設(shè)，當(dāng)式（7）中的h=0 時，即表示光滑壁面的熱流，應(yīng)與式（10）相同，那么由此可得：

對于式（7）中第1 項(xiàng)粗糙壁引起的峰值熱流增量，從第2 節(jié)的數(shù)值結(jié)果來看，顯然，它的系數(shù)c1應(yīng)與來流參數(shù)（馬赫數(shù)、壁溫比）和粗糙帶的特征有關(guān)。結(jié)合第2 節(jié)中對粗糙帶高度及寬度的數(shù)值分析結(jié)果，作出以下假設(shè)：

相應(yīng)地，相對峰值熱流增量為

式（13）中需要確定的參數(shù)有a、b 和C。使用2.2 節(jié)、2.4 節(jié)與2.5 節(jié)中得到的全場波紋壁上的峰值熱流數(shù)值實(shí)驗(yàn)數(shù)據(jù)，計(jì)算其相對峰值熱流增量（此量僅于馬赫數(shù)、壁溫比、粗糙壁特征有關(guān)，與流向位置無關(guān)）。對流向取平均，使用冪次關(guān)系擬合這個值與寬度、馬赫數(shù)、壁溫比的關(guān)系（分別對應(yīng)2.2 節(jié)、2.4 節(jié)與2.5 節(jié)中的工況）。擬合的結(jié)果如圖13 所示，圖13（a）～圖13（c）分別給出了對寬度、馬赫數(shù)以及壁溫的擬合結(jié)果。由擬合公式可以最終確定，C 為0.5，a 為2.67，b 為1.38。

圖13 熱流預(yù)測公式冪次關(guān)系的修正Fig.13 Modification of exponential relationship of heating prediction formula

圖14 給出了參數(shù)確定后，利用熱流式（13）反求得到的不同寬度、馬赫數(shù)、壁溫情況下的Δ。可以發(fā)現(xiàn)這一常數(shù)整體集中在300～400 之間，根據(jù)本文的假設(shè)，即Δ 為表示粗糙類型的常數(shù)，通過將所有工況（均為正弦粗糙）先沿流向位置平均，計(jì)算得到的Δ 再取平均，最終確定正弦壁面的這一常數(shù)為322。

圖14 熱流預(yù)測公式中Δ的計(jì)算Fig.14 Calculation of Δ in heating prediction formula

這樣可以得到預(yù)測粗糙帶壁面峰值熱流的公式：

對于正弦波紋壁Δ可取322。那么，在已知粗糙壁的形狀、位置、來流參數(shù)的前提下，可以使用式（15）預(yù)測由粗糙壁導(dǎo)致的壁面峰值熱流。

4 預(yù)測模型驗(yàn)證

本文提出的壁面峰值熱流預(yù)測公式需要進(jìn)一步進(jìn)行以下驗(yàn)證：一是基于相似性解得到的熱流公式，與實(shí)際直接數(shù)值模擬得到的熱流是否一致；二是式（15）在實(shí)際運(yùn)用到粗糙壁的熱流預(yù)測時，其誤差及適用范圍有多大。

4. 1 光滑壁面的熱流驗(yàn)證

圖15 給出了光滑壁面，在不同馬赫數(shù)下基于相似性解通過變換得到的熱流和熱流預(yù)測公式（15）的比較，結(jié)果顯示，兩者得到的熱流結(jié)果是一致的，這驗(yàn)證了理論推導(dǎo)過程的正確性。

圖15 相似性解熱流結(jié)果與預(yù)測公式結(jié)果Fig.15 Similarity surface heat flux results and prediction formula results

在光滑平板上，圖16 給出了相似性解的熱流結(jié)果與實(shí)際DNS 數(shù)值結(jié)果的比較，可以看出，兩者的熱流計(jì)算結(jié)果差距小于3%，吻合的較好，這驗(yàn)證了使用相似性解來計(jì)算熱流的可行性。

圖16 相似性解熱流結(jié)果與DNS 數(shù)值模擬結(jié)果Fig.16 Similarity solution surface heat flux results and DNS numerical simulation results

4. 2 全場均布波紋粗糙壁面的熱流驗(yàn)證

完成對于式（15）理論依據(jù)的驗(yàn)證后，對其在預(yù)測粗糙壁面峰值熱流時的誤差進(jìn)行分析。首先使用控制變量法，對高度、寬度、壁溫、馬赫數(shù)分別進(jìn)行驗(yàn)證。使用相對誤差以及絕對誤差來衡量其大小。理論公式預(yù)測的壁面峰值熱流與實(shí)際DNS 數(shù)值結(jié)果的相對誤差：

式中：chDNS為直接數(shù)值模擬得到的熱流系數(shù)。絕對誤差：

式中：Qw 為有量綱的熱流；Qwpredict為預(yù)測公式得到的有量綱的熱流。

圖17～圖21 給出了在不同高度、不同寬度、相同寬高比、不同馬赫數(shù)以及不同壁溫下的壁面峰值熱流的絕對誤差和同一流向位置實(shí)際與預(yù)測結(jié)果的比較。同時將其相對誤差在圖中標(biāo)出?？梢钥闯?，在大多數(shù)工況，除了在平板頭部由于黏性干擾導(dǎo)致誤差較大，下游的相對誤差都在5%以下。值得注意的是，在不同寬度下的相對誤差較大（＜10%），其原因可能為模型使用的用于衡量粗糙特征的冪次表達(dá)有一定的局限性，但絕對誤差也可保持在一個較小的范圍內(nèi)。相同寬高比、不同馬赫數(shù)以及不同壁溫下的預(yù)測結(jié)果誤差都在3%之內(nèi)。表明在已有工況下，式（15）預(yù)測效果不錯。

圖17 不同高度下全場波紋壁預(yù)測公式與DNS 結(jié)果壁面峰值熱流誤差Fig.17 Surface peak heat flux error between full-field corrugated wall prediction formula and DNS results at different heights

圖18 不同寬度下全場波紋壁預(yù)測公式與DNS 結(jié)果壁面峰值熱流誤差Fig.18 Surface peak heat flux error between full-field corrugated wall prediction formula and DNS results at different widths

圖19 相同寬高比下的全場波紋壁預(yù)測公式與DNS 結(jié)果壁面峰值熱流誤差Fig.19 Surface peak heat flux error between full-field corrugated wall prediction formula and DNS results with same aspect ratio

圖20 不同馬赫數(shù)下全場波紋壁預(yù)測公式與DNS 結(jié)果壁面峰值熱流誤差Fig.20 Surface peak heat flux error between full-field corrugated wall prediction formula and DNS results at different Mach numbers

圖21 不同壁溫下全場波紋壁預(yù)測公式與DNS 結(jié)果壁面峰值熱流誤差Fig.21 Surface peak heat flux error between full-field corrugated wall prediction formula and DNS results at different wall temperatures

4. 3 測試工況的驗(yàn)證

為了進(jìn)一步驗(yàn)證式（15）對于更為復(fù)雜的工況，即多變量均變化情況下的適用性，對更多的測試工況進(jìn)行了驗(yàn)證。進(jìn)行了高度、壁溫、馬赫數(shù)均變化的工況的驗(yàn)證，計(jì)算工況見表3。

表3 測試工況的計(jì)算參數(shù)Table 3 Calculation parameters of test conditions

圖22 給出3 個驗(yàn)證工況下的壁面峰值熱流以及相對誤差結(jié)果，可以看出，在多變量變化的情況下，所提出的預(yù)測模型的誤差可以維持在較低的范圍內(nèi)（不超過3%），這也證明了式（15）對于多變量變化的復(fù)雜情況的適用性。

圖22 測試工況驗(yàn)證結(jié)果Fig.22 Test case verification results

4. 4 極限工況的驗(yàn)證

對建立預(yù)測模型工況的壁溫、馬赫數(shù)的最大和最小值，即極限的幾個工況，模型的適用性進(jìn)行進(jìn)一步驗(yàn)證的。使用的工況見表4。

表4 極限工況的計(jì)算參數(shù)Table 4 Calculation parameters of limited conditions

圖23 給出壁面峰值熱流比較以及模型預(yù)測的相對誤差?？梢园l(fā)現(xiàn)，實(shí)際數(shù)值模擬的壁面峰值熱流與預(yù)測公式的結(jié)果基本吻合；其中，馬赫數(shù)為4、壁溫為800 K 時的相對誤差較大，原因在于此時壁溫接近絕熱壁溫（841.51 K），熱流較小，相對偏差會更大一些。

圖23 極限工況驗(yàn)證結(jié)果Fig.23 Limited case verification results

4. 5 其他粗糙壁面形式的適用性

對于建立的熱流預(yù)測模型式（15），詳細(xì)驗(yàn)證了其對于正弦型連續(xù)波紋壁面的可用性，確定了其適用范圍，為了進(jìn)一步研究此模型對于其他壁面形式的適用性，數(shù)值計(jì)算分散的正弦型壁面以及橢圓型壁面，壁面形式如圖24 所示。圖24（a）為分散的正弦壁面，在0～1 m 的范圍內(nèi)，布置了50 個高度為0.210 6 mm、寬度為10.53 mm 的正弦型突起，相鄰粗糙元之間的間距為10.53 mm，在粗糙元附近進(jìn)行了局部加密；圖24（b）為分散的橢圓壁面，在0～1 m 的范圍內(nèi)，布置了40 個高度為0.210 6 mm、寬度為7.02 mm 的橢圓型凸起，相鄰粗糙元之間間距為14.04 mm，粗糙元附近進(jìn)行了局部加密。兩套網(wǎng)格流向分布2 501 個點(diǎn)，法向分布401 個點(diǎn)。

圖24 壁面形狀Fig.24 Wall shape

使用上述配置進(jìn)行數(shù)值計(jì)算，在計(jì)算分散的正弦壁面時，Δ取322，與之前的正弦波紋壁是一致的，這是符合Δ表征粗糙帶類型的常數(shù)這一假設(shè)的；而計(jì)算橢圓時，首先需要確定Δ，使用的方法與第4 節(jié)中方法相同，使用數(shù)值得到的相對熱流增量求得Δ，發(fā)現(xiàn)其分布為1 條與流向無關(guān)的直線，符合對于Δ的定義，對所有流向位置得到的Δ取平均，得到Δ為300。圖25 給出了不同壁面形式下粗糙附近的熱流分布，不同壁面形式下，粗糙附近的熱流變化是不同的。圖26 給出了2 種形式下壁面峰值熱流的變化規(guī)律，同時比對了其與預(yù)測值的差異，發(fā)現(xiàn)兩者是基本吻合的。圖27 則給出了壁面峰值熱流預(yù)測模型的相對誤差，可以看出，相對誤差在2%以下，說明預(yù)測公式對于其他2 種壁面形狀也具有一定的適用性。

圖25 不同壁面形式下的壁面熱流Fig.25 Surface heat flux with different wall shapes

圖26 不同壁面形式下的壁面峰值熱流及其預(yù)測Fig.26 Surface peak heat flux and its prediction with different wall shapes

圖27 不同壁面形式下的壁面峰值熱流預(yù)測誤差Fig.27 Prediction error of surface peak heat flux with different wall shapes

5 結(jié) 論

針對高超聲速平板層流邊界層，研究了表面均布高度小于600 μm 正弦型波紋凸起的粗糙導(dǎo)致的壁面熱流的變化規(guī)律，并提出了一種可用于定量預(yù)測由波紋粗糙壁引起的壁面峰值熱流的公式，詳細(xì)驗(yàn)證了其適用性。

使用直接數(shù)值模擬研究粗糙壁在不同馬赫數(shù)、高度、寬度以及壁溫下的壁面峰值熱流，發(fā)現(xiàn)壁面峰值熱流的變化與粗糙壁的高度成正比，而與粗糙壁的寬度成反比；壁面峰值熱流受高度的影響大于受寬度的影響；壁面相對峰值熱流增量隨馬赫數(shù)增大而減??；壁面相對峰值熱流增量隨壁溫增大而減小。

基于數(shù)值結(jié)果得到的定性規(guī)律，從光滑壁面熱流的理論解出發(fā)，假設(shè)熱流增量與粗糙高度、寬度、馬赫數(shù)和壁溫成冪次關(guān)系，得到預(yù)測粗糙壁面峰值熱流的公式，并根據(jù)數(shù)據(jù)擬合得到了各冪次關(guān)系中的系數(shù)。首先通過控制變量法，對于寬度、高度、馬赫數(shù)、壁溫等單一變量進(jìn)行了驗(yàn)證。然后驗(yàn)證了這些參數(shù)均變化時公式的有效性，并在極限工況下測試了其適用范圍。最后在2 種不同的壁面粗糙形式下進(jìn)一步探究了公式的可用性。整體上，預(yù)測模型的誤差維持在8%以下。

本文關(guān)注的粗糙元都是緩變的，對于高度與寬度大小相當(dāng)?shù)拇植谠?，本文的定性結(jié)論應(yīng)該能夠適用，但定量預(yù)測公式是否適用，還有待進(jìn)一步研究。

致 謝

感謝中國空氣動力研究與發(fā)展中心國義軍研究員有啟發(fā)性的幫助。