基于γ-Reθ轉(zhuǎn)捩模型的高超聲速 復(fù)雜構(gòu)型轉(zhuǎn)捩模擬

易淼榮, 趙慧勇, 樂嘉陵

(中國空氣動力研究與發(fā)展中心 高超聲速沖壓發(fā)動機(jī)技術(shù)重點實驗室, 四川 綿陽 621000)

0 引 言

邊界層轉(zhuǎn)捩現(xiàn)象對航空領(lǐng)域中的摩阻和熱流具有重要影響。對于高超聲速吸氣式飛行器來說，若流入進(jìn)氣道的邊界層為層流，則有可能因逆壓梯度和激波邊界層干擾而引起流動分離，但湍流邊界層在抵抗分離方面具有優(yōu)勢。因此，為了提高飛行穩(wěn)定性，同時增加燃燒室內(nèi)燃料和氧氣的混合，經(jīng)常會在進(jìn)氣道前體上安裝強(qiáng)制轉(zhuǎn)捩裝置，以促進(jìn)邊界層在進(jìn)入唇口之前轉(zhuǎn)捩為湍流。然而，盡管對邊界層轉(zhuǎn)捩現(xiàn)象的研究已經(jīng)長達(dá)一個多世紀(jì)，人們?nèi)匀徊荒艽_切地給出什么情況下邊界層會由層流轉(zhuǎn)捩成為湍流，特別是在高超聲速復(fù)雜構(gòu)型下，由于受到太多因素(湍流度、雷諾數(shù)、馬赫數(shù)、粗糙度、鈍度、壁溫、幾何構(gòu)型等)及各種因素之間互相耦合作用的影響，使得轉(zhuǎn)捩預(yù)測更加困難。

現(xiàn)階段，有多種數(shù)值方法可以用來進(jìn)行邊界層轉(zhuǎn)捩研究。DNS和LES被用來研究典型狀態(tài)下的基本轉(zhuǎn)捩機(jī)制和詳細(xì)的轉(zhuǎn)捩過程[1-2]；LST和PSE 則被用來研究非穩(wěn)定波的形成和發(fā)展歷程；轉(zhuǎn)捩準(zhǔn)則則在總結(jié)大量實驗數(shù)據(jù)的基礎(chǔ)上，用簡單的形式給出轉(zhuǎn)捩位置[3-4]，如NASP的尖前緣平板轉(zhuǎn)捩準(zhǔn)則：Reθ/Mae=305[5]。目前來看，半經(jīng)驗的eN方法在工程上的應(yīng)用最為廣泛，但是eN方法也有幾個致命缺陷：由于其基于線性穩(wěn)定性理論，故難以處理較為復(fù)雜的構(gòu)型，更無法模擬粗糙顆粒誘導(dǎo)的轉(zhuǎn)捩；同時N因子的確定也需要大量的試驗數(shù)據(jù)，而且其在不同馬赫數(shù)、不同構(gòu)型下的具體數(shù)值也不一樣；而且，來流湍流度本是影響轉(zhuǎn)捩的一個重要因素，但在eN方法中卻沒有任何體現(xiàn)。相反，基于RANS模型發(fā)展而來的工程轉(zhuǎn)捩模型，則因其具有較好的魯棒性且易于與現(xiàn)有的湍流模型框架融合而在最近越來越得到人們的重視，其中，Langtry和Menter[6-8]基于SST湍流模型，提出了一種基于經(jīng)驗耦合的間歇因子轉(zhuǎn)捩模型——γ-Reθ模型。在該模型中，所有的計算都是基于局部變量，沒有任何沿流線積分的操作，因此可以輕易地對復(fù)雜構(gòu)型進(jìn)行模擬，Langtry將該模型應(yīng)用于大量低速工程算例。由于γ-Reθ轉(zhuǎn)捩模型只是為把基于經(jīng)驗關(guān)系式的模型結(jié)合進(jìn)一般CFD方法中提供了一個框架，而具體的物理機(jī)制均包含在經(jīng)驗關(guān)系式中，因此很多研究者在其框架下針對各種不穩(wěn)定機(jī)制對模型進(jìn)行修正，從而適應(yīng)不同情況下的轉(zhuǎn)捩模擬。如Watanabe[9]通過類比G?rtler不穩(wěn)定性和橫流不穩(wěn)定性，引入類似于G?rtler數(shù)的Kohama參數(shù)，Krumbein[10-11]針對三維邊界層，通過求解Falkner-Skan-Cooke近似解引入橫流轉(zhuǎn)捩動量厚度雷諾數(shù)，各自提出了基于該框架下能夠模擬橫流模態(tài)的轉(zhuǎn)捩模型。Langel[12]則通過引入“粗糙幅度” 的輸運方程發(fā)展了能夠反映壁面粗糙度影響的修正版本。在高超聲速領(lǐng)域，Krause[13]和Zhang[14]分別提出了可以應(yīng)用于高超聲速情況下轉(zhuǎn)捩預(yù)測的經(jīng)驗關(guān)系式及相應(yīng)的修正，并都將其應(yīng)用于Ma8.3的雙楔平板模型，得到了與試驗吻合較好的壁面壓力分布。Cheng[15]則添加了一個針對壓力梯度的修正因子，對Ma7.93的7°尖錐進(jìn)行了模擬，研究了來流湍流度和網(wǎng)格對計算結(jié)果的影響。Bensassi[16]模擬了Ma8的7°尖錐模型，得到了與試驗一致的壁面熱流分布。但是，目前并沒有見到該模型在更復(fù)雜構(gòu)型上的應(yīng)用，如光滑壁面或安裝有轉(zhuǎn)捩帶的高超聲速進(jìn)氣道構(gòu)型。

本文基于Langtry和Menter的γ-Reθ模型，進(jìn)行一些必要的修正，建立一套適用于高超聲速情況下的轉(zhuǎn)捩預(yù)測方法。通過較為簡單的三維構(gòu)型自然轉(zhuǎn)捩和粗糙顆粒誘導(dǎo)強(qiáng)制轉(zhuǎn)捩的算例對模型進(jìn)行了驗證。將修正的模型應(yīng)用于X-51A前體邊界層自然轉(zhuǎn)捩和強(qiáng)制轉(zhuǎn)捩的研究中，通過與試驗結(jié)果進(jìn)行比較，表明本文建立的方法能夠較好地進(jìn)行高超聲速復(fù)雜構(gòu)型中的自然轉(zhuǎn)捩和強(qiáng)制轉(zhuǎn)捩預(yù)測。

1 數(shù)值計算方法

1.1 流場求解器

本文的γ-Reθ模型是基于一套有限體積大規(guī)模并行結(jié)構(gòu)網(wǎng)格Navier-Stokes求解器AHL3D建立的。該求解器由中國空氣動力研究與發(fā)展中心吸氣式高超聲速技術(shù)研究中心開發(fā)，能夠?qū)硐霘怏w和化學(xué)非平衡狀態(tài)下的混合氣體進(jìn)行求解[17-18]。時間積分采用 LU-SGS，無粘通量采用重構(gòu)-推進(jìn)方法，重構(gòu)采用 3 階MUSCL插值，推進(jìn)采用Steger-Warming和AUSMPW+格式，粘性通量采用改進(jìn)的Gauss定理計算。湍流模型為Menter的SST雙方程模型。

1.2 γ-Reθ轉(zhuǎn)捩模型

Langtry和Menter[6]提出的γ-Reθ轉(zhuǎn)捩模型主要基于間歇因子γ和轉(zhuǎn)捩起始動量厚度雷諾數(shù)Reθt的輸運方程。γ方程用來控制轉(zhuǎn)捩程序，其源項由臨界動量厚度雷諾數(shù)Reθc和最大渦雷諾數(shù)Rev之間的關(guān)系來控制，而Reθc又與Reθt存在對應(yīng)關(guān)系。Reθt則被當(dāng)做一個輸運標(biāo)量，主要思想就是使用經(jīng)驗耦合關(guān)系式得到邊界層外的Reθt，并且允許Reθt通過輸運方程擴(kuò)散到邊界層內(nèi)去。具體形式如下：

(1)

(2)

Pγ=Flengthca1ρS[γFonset]0.5(1-ce1γ)

(3)

Eγ=ca2ρΩγFturb(ce2γ-1)

(4)

(5)

(6)

Reθ t=F(Tu)F(λθ)

(7)

(8)

(9)

(10)

(11)

(12)

1.3 針對模型的修正

原始的轉(zhuǎn)捩模型是針對低速情況提出的，因此不能直接應(yīng)用于高超聲速情況下的轉(zhuǎn)捩預(yù)測，這個模型一經(jīng)提出，很多研究者均在其框架下針對各種情況對模型的Reθc和Flength進(jìn)行修正，前者用來修正轉(zhuǎn)捩起始位置，后者用來修正轉(zhuǎn)捩區(qū)域長度。如Krause和張曉東均提出了該框架下能夠應(yīng)用于高超聲速情況下的經(jīng)驗關(guān)系式。Krause[13]的修正如下：

(13)

(14)

將來流粘性比RT∞設(shè)置成0.001，因此Tu在到達(dá)模型前緣時均衰減到了一個很小的值，通過將Tuinlet加入到經(jīng)驗關(guān)系式中來體現(xiàn)不同來流湍流度對轉(zhuǎn)捩的影響規(guī)律，經(jīng)過系列低速平板的標(biāo)定后，應(yīng)用于Ma8.3下雙楔模型試驗的轉(zhuǎn)捩預(yù)測中，得到的壁面壓力系數(shù)與試驗吻合較好，顯示出高超聲速情況下預(yù)測轉(zhuǎn)捩的潛力。尚慶[19]也采用該修正關(guān)系式對雙楔模型進(jìn)行了更進(jìn)一步的研究，其結(jié)果顯示基于該修正關(guān)系式的γ-Reθ轉(zhuǎn)捩模型能夠在壓力系數(shù)、熱流Stanton數(shù)和分離區(qū)大小上與試驗吻合良好。

張曉東[14]的修正則如下所示：

張曉東還基于Reshotko[20]對高超聲速試驗中Reθ t與邊界層外緣馬赫數(shù)Mae關(guān)系的總結(jié)，利用當(dāng)?shù)伛R赫數(shù)對模型進(jìn)行了馬赫數(shù)修正，以適應(yīng)不同Ma數(shù)下的轉(zhuǎn)捩預(yù)測情況：

Ma∈[0,12]

(17)

張曉東也對雙楔模型進(jìn)行了與試驗吻合較好的模擬。

本文利用文獻(xiàn)[21]中的進(jìn)氣道構(gòu)型，采用二維計算對轉(zhuǎn)捩模型進(jìn)行考核及修正。進(jìn)氣道外壓縮形式為四波系平面頂壓，流道全長為550mm，本文模擬其唇口之前壓縮面的情況，模型構(gòu)型及網(wǎng)格如圖1所示。來流條件如表1所示，表中case3模擬的是在距前緣88mm的第一道壓縮面上安裝了高度為1mm的鉆石型轉(zhuǎn)捩裝置的試驗結(jié)果，轉(zhuǎn)捩帶的具體參數(shù)可參考文獻(xiàn)[22]，由于是二維模擬，轉(zhuǎn)捩帶構(gòu)型模化為矩形。計算時為了更好地與試驗對比熱流值，壁面溫度根據(jù)實驗測量結(jié)果分段給定。對比用到的試驗數(shù)據(jù)來自于文獻(xiàn)[21-22]。

圖1 高超聲速進(jìn)氣道二維模型及網(wǎng)格示意圖

CaseMa總壓/MPa總溫/K迎角/(°)Re/(107·m-1)轉(zhuǎn)捩類型15.962.00473.21.01.82自然24.971.01357.7-0.62.21自然35.962.00473.21.01.82強(qiáng)制

試驗在FL-31風(fēng)洞中進(jìn)行。FL-31屬于常規(guī)噪聲風(fēng)洞，根據(jù)一般的認(rèn)識，常規(guī)噪聲風(fēng)洞的來流湍流度約為1%的量級。在case1情況下，Menter、Krause及張曉東三種版本的經(jīng)驗關(guān)系式在不同來流湍流度下的計算結(jié)果與試驗結(jié)果比較如圖2所示。從圖中可以看出，Menter版本和Krause版本在來流Tu=0.03時(圖中Tu0.03表示來流湍流度為0.03%，RT1表示來流湍流粘性比μt/μ為1)，轉(zhuǎn)捩位置仍早于試驗值，明顯過早地預(yù)測了轉(zhuǎn)捩的出現(xiàn)。因此，在本程序框架下，Menter版本和Krause版本將無法正確模擬FL-31風(fēng)洞湍流度水平下的轉(zhuǎn)捩，更無法模擬靜音風(fēng)洞或飛行試驗條件等低湍流度情況下的轉(zhuǎn)捩，而張曉東版本則能在本程序框架下正確反映FL-31風(fēng)洞的來流湍流度的量級，因此選取張曉東提出的經(jīng)驗關(guān)系式作為本文采用的版本。

圖2 不同計算方法得到的熱流分布(case1)

圖3 自然轉(zhuǎn)捩熱流分布(case2)

圖4 強(qiáng)制轉(zhuǎn)捩熱流分布(case3)

圖5 自然轉(zhuǎn)捩熱流分布(case1，修正后的模型)

Fig.5Heatfluxdistributionofnaturaltransition(case1,modifiedmodel)

圖6 自然轉(zhuǎn)捩熱流分布(case2，修正后的模型)

Fig.6Heatfluxdistributionofnaturaltransition(case2,modifiedmodel)

圖7 強(qiáng)制轉(zhuǎn)捩熱流分布(case3，修正后的模型)

Fig.7Heatfluxdistributionofforcedtransition(case3,modifiedmodel)

修改后的k-w方程為：

(18)

其中：

α1=1.0,α2=0.4,α3=0.2

兩個經(jīng)驗關(guān)系式為：

(19)

(20)

此外，對Reθ t的修正改為：

Ma∈[0,12]

(21)

1.4 網(wǎng)格無關(guān)性驗證

將修正后的轉(zhuǎn)捩模式在1.3節(jié)中的二維構(gòu)型中進(jìn)行網(wǎng)格無關(guān)性驗證。選取1.3節(jié)中的case2，并以1.3節(jié)中的網(wǎng)格為base網(wǎng)格，網(wǎng)格在第一至第四道壓縮面的流向上分別布置了121、61、41、65個網(wǎng)格點，在法向上布置了141個網(wǎng)格點，第一層網(wǎng)格距離壁面0.001mm?？偟木W(wǎng)格量約為6萬。再在此基礎(chǔ)上在法向加密1.5倍得到v1.5網(wǎng)格，在法向和流向各加密1.5倍得到v1.5-s1.5網(wǎng)格。3套網(wǎng)格的計算結(jié)果如圖2所示。從中可以看出，3套網(wǎng)格計算得到的熱流幾乎重合，可認(rèn)為轉(zhuǎn)捩模型此時已經(jīng)達(dá)到網(wǎng)格無關(guān)性。

圖8 多套網(wǎng)格計算結(jié)果(case2，修正后的模型)

2 計算結(jié)果與分析

由于上一節(jié)對模型進(jìn)行的修正都是以二維計算為基礎(chǔ)，為了驗證其是否能夠適應(yīng)更寬的馬赫數(shù)、雷諾數(shù)范圍及三維復(fù)雜構(gòu)型的轉(zhuǎn)捩預(yù)測，選取了Ames的全尺寸升力體構(gòu)型在有迎角情況下的轉(zhuǎn)捩，驗證模型對高超聲速三維構(gòu)型自然轉(zhuǎn)捩的預(yù)測能力；通過Ma6平板上單個粗糙顆粒誘導(dǎo)轉(zhuǎn)捩的算例，驗證模型對粗糙顆粒誘導(dǎo)強(qiáng)制轉(zhuǎn)捩的預(yù)測能力。最終，將模型應(yīng)用于X-51A進(jìn)氣道構(gòu)型的轉(zhuǎn)捩模擬當(dāng)中。

2.1 Ma7.4 下 Ames 的 all-body 試驗

Lockman[24]對Ames研究中心的all-body三維升力體提供了一個完整的高超聲速三維實驗，模型是一個用于高超聲速巡航的通用升力體外形。實驗是在Ames的3.5英寸高超聲速風(fēng)洞開展的。模型示意圖見圖9。模型全長0.9144m，具有75°后掠角的尖三角翼形狀。全體是一個主軸比4的橢圓錐，后體的橫截面是具有直邊的橢圓形。實驗條件如表2所示。網(wǎng)格拓?fù)淙鐖D10所示，保持壁面第一層y+<1，網(wǎng)格量約為150萬。來流湍流度設(shè)為 0.5%。

表2 全尺寸升力體模型風(fēng)洞試驗來流條件Table 2 Test conditions for the all-body model

圖9 全尺寸升力體模型示意圖

圖10 全尺寸升力體模型網(wǎng)格示意圖

試驗測量了迎風(fēng)面和背風(fēng)面表征熱流的Stanton數(shù)，即St，其定義為:

ρ∞和U∞分別表示來流密度和速度，Hw為壁面上流體的焓值，Ht表示來流總焓。

圖11為中心線St數(shù)分布，對于迎風(fēng)面，測量到的St數(shù)顯示轉(zhuǎn)捩大概在距前緣0.25m的位置發(fā)生，但是原始的γ-Reθ模型給出的轉(zhuǎn)捩位置卻在距前緣0.06m左右，這是由于高馬赫數(shù)和可壓縮性導(dǎo)致了轉(zhuǎn)捩起始動量厚度雷諾數(shù)的增加，而這一增加效果在原始的轉(zhuǎn)捩模型中被忽略掉了，因此轉(zhuǎn)捩位置被大大前移。所以原始的轉(zhuǎn)捩模型和全湍流下的 SST 模型均在x=0.25m的上游區(qū)域給出了遠(yuǎn)高于試驗值的熱流。而本文修改的模型，不僅準(zhǔn)確地捕捉到了轉(zhuǎn)捩的位置，壁面熱流也與試驗吻合較好。轉(zhuǎn)捩結(jié)束位置的熱流大約是層流區(qū)的3倍，并且要稍高于全湍流計算的結(jié)果。對于背風(fēng)面，由于壁面溫度與恢復(fù)溫度比較接近，熱流值都比較小，因此不同模型得到的熱流相差不大。圖12所示為迎風(fēng)面的St數(shù)分布和摩擦力線。原始的γ-Reθ和全湍流 SST 模型給出了相似的St數(shù)分布。本文修改的模型計算結(jié)果則顯示：在上游區(qū)域，熱流情況基本與層流結(jié)果一致, 0.23m

圖11 中心線St數(shù)分布

Fig.12Stantonnumberdistributionandfrictionforcelinesofthewindward

2.2 Ma6平板上單個粗糙顆粒誘導(dǎo)轉(zhuǎn)捩

Tirtey和Chazot[25-26]在VKI H3高超聲速風(fēng)洞中對單個粗糙顆粒誘導(dǎo)平板轉(zhuǎn)捩情況進(jìn)行了一系列試驗研究。粗糙顆粒構(gòu)型包括圓柱形、鉆石形、斜坡形和半球形，本文選取圓柱形和鉆石形粗糙顆粒在高雷諾數(shù)下的情況進(jìn)行模擬。模型如圖13所示，一塊長 290mm、寬100mm，前緣鈍度 0.5mm的平板，在距前緣60mm處安裝了一個粗糙顆粒，粗糙顆粒細(xì)節(jié)如圖13(b)所示，其中a=4mm，對于圓柱形粗糙顆粒，k=1.02mm，對于鉆石形粗糙顆粒，k=0.83mm。來流馬赫數(shù)為 6，來流單位雷諾數(shù)為2.7×107/m，來流總壓3.1MPa，總溫為500K，壁溫在試驗開始時均為(292±2)K，試驗結(jié)束后層流區(qū)最大溫升為4K，湍流區(qū)最大溫升為12K，由于這里壁溫變化都不大，因此選取300K為壁面溫度。粗糙顆粒附近的網(wǎng)格如圖14所示，第一層壁面y+<1，網(wǎng)格量為500萬左右。來流湍流度設(shè)為 1.5%。

(a) (b)

圖14 粗糙顆粒附近網(wǎng)格示意圖

此算例中St數(shù)定義為：

Tt表示總溫，Tw表示壁溫。圓柱形粗糙顆粒尾跡區(qū)的紅外熱試驗和計算的St數(shù)分布對比如圖15、16所示，試驗和本文修正的γ-Reθ轉(zhuǎn)捩模型均顯示粗糙顆粒誘導(dǎo)的尾跡區(qū)的熱流明顯超出了非尾跡區(qū)，而對于原始的γ-Reθ模型，尾跡區(qū)和非尾跡區(qū)均出現(xiàn)了轉(zhuǎn)捩現(xiàn)象，當(dāng)來流湍流度被修改為 0.5% 時，非尾跡區(qū)在下游階段仍出現(xiàn)了轉(zhuǎn)捩現(xiàn)象，而尾跡區(qū)的轉(zhuǎn)捩現(xiàn)象卻比試驗更靠后。這表明原始的γ-Reθ轉(zhuǎn)捩模型不僅對自然轉(zhuǎn)捩的預(yù)測比實際情況早，而且對粗糙顆粒對轉(zhuǎn)捩的促進(jìn)作用模擬不夠。而修正后的模型克服了以上兩個缺點，在尾跡區(qū)和非尾跡區(qū)均給出了與試驗吻合較好的結(jié)果。將修正后的模型計算得到的尾跡區(qū)和非尾跡區(qū)的St數(shù)取出，與試驗值以及非尾跡區(qū)的全湍流值進(jìn)行對比，如圖17所示。當(dāng)流體流過粗糙顆粒后，尾跡區(qū)的壁面熱流出現(xiàn)了一個快速增長區(qū)，之后由于分離渦的作用，其值降低到一個較低值，而后再慢慢增長到全湍流值。St數(shù)在尾跡區(qū)的值是非尾跡區(qū)的2.5～3倍。由于全湍流計算結(jié)果中粗糙顆粒后緣的尾跡區(qū)也有一個熱流值較低的區(qū)域，因此可認(rèn)為轉(zhuǎn)捩在粗糙顆粒后緣就已經(jīng)完成。試驗與計算的對比顯示：修正的模型可以模擬粗糙顆粒誘導(dǎo)轉(zhuǎn)捩的整個過程，也能較好地預(yù)測壁面熱流。但是，在試驗中，隨著尾跡區(qū)的流體向下游流動，非穩(wěn)定結(jié)構(gòu)在展向發(fā)展明顯，因此轉(zhuǎn)捩區(qū)域越來越寬，而在本文的計算中，轉(zhuǎn)捩區(qū)域在展向的發(fā)展卻比較緩慢。圖16(a)、(b)和圖17(b)顯示的是修改后的模型模擬鉆石形粗糙顆粒的情況，結(jié)果與圓柱形粗糙顆粒比較類似。

(a) 紅外熱試驗結(jié)果

(b) 計算結(jié)果

(c) 計算結(jié)果

(d) 計算結(jié)果

圖15 紅外熱試驗和計算得到的St數(shù)分布(圓柱形粗糙顆粒誘導(dǎo)轉(zhuǎn)捩)

Fig.15Infra-redvisualizationandcalculatedStantonnumberdistributionforcylinderroughnesselements

(a) 紅外熱試驗結(jié)果

(b) 計算結(jié)果

圖16 紅外熱試驗和計算得到的St數(shù)分布(鉆石形粗糙顆粒誘導(dǎo)轉(zhuǎn)捩)

Fig.16Infra-redvisualizationandcalculatedStantonnumberdistributionfordiamondroughnesselements

(a) 圓柱形粗糙顆粒

(b) 鉆石形粗糙顆粒

Fig.17ComparisonofStnumberinandoutofwakeregionforcylinderanddiamondroughnesselements

圖18則為粗糙顆粒下游流向渦的發(fā)展情況，其中展向切片上顯示的是流向渦量的大小|ωx|(只顯示了|ωx|>2000s-1的情況)。從圖中可知，對于兩種形式的粗糙顆粒，下游均發(fā)展形成了由粗糙顆粒后緣發(fā)展起來的主渦和粗糙顆粒兩側(cè)形成的次渦，并在繼續(xù)向下游的發(fā)展過程中慢慢耗散掉，流動則轉(zhuǎn)捩為湍流。但本文的計算結(jié)果顯示其流向渦破碎之后并沒有往展向有明顯的拓展，因此轉(zhuǎn)捩區(qū)域在展向也并未出現(xiàn)明顯的展向拓展，考慮到尾跡導(dǎo)致的轉(zhuǎn)捩位置區(qū)域很靠前(圖15(a))和很靠后(圖15(c))時，尾跡往展向的發(fā)展速度均很緩慢，因此這可能是由于計算方法中耗散過大致使未能捕捉到足夠的尾跡流動細(xì)節(jié)造成的。而過大的耗散則是因雷諾平均引入的湍流粘性和MUSCL插值造成較大的數(shù)值耗散共同導(dǎo)致的。

(a) 圓柱形粗糙顆粒

(b) 鉆石形粗糙顆粒

Fig.18Developmentofstreamwisevortexdownstreamoftheroughnesselements

2.3 20% 縮比X-51A前體模型在BAM6QT中的轉(zhuǎn)捩

目前已有的轉(zhuǎn)捩預(yù)測方法雖然很多，但是研究的多為平板、圓錐等簡單構(gòu)型，而針對高超聲速進(jìn)氣道構(gòu)型進(jìn)行完整模擬的還比較少。王亮[27-28]用其提出的三方程轉(zhuǎn)捩模型對X-51A 20%縮比模型在普渡大學(xué)Ma6靜音風(fēng)洞中的自然轉(zhuǎn)捩試驗結(jié)果進(jìn)行了模擬，結(jié)果顯示該轉(zhuǎn)捩模型抓住了靜音狀態(tài)下流動一直保持為層流的特征，也較準(zhǔn)確地刻畫了噪聲狀態(tài)下的轉(zhuǎn)捩起始位置，但是在轉(zhuǎn)捩結(jié)束位置上與試驗結(jié)果還存在較大差別。Orlik[29]則采用了風(fēng)洞試驗與數(shù)值模擬相結(jié)合的方法，研究了一個高超聲速前體的自然轉(zhuǎn)捩情況，其計算采用的是基于線性穩(wěn)定性分析的eN方法，由于對擾動增長率的積分并不是從中性線上開始的，而是從距前緣1.33mm的固定位置開始，因此積分起始位置的N值是未知的，造成了對N值的刻畫出現(xiàn)較大的不確定性。本文希望能夠針對進(jìn)氣道模型同時對自然轉(zhuǎn)捩和強(qiáng)制轉(zhuǎn)捩進(jìn)行模擬，從而考核本文建立的轉(zhuǎn)捩模型對高超聲速復(fù)雜構(gòu)型的轉(zhuǎn)捩預(yù)測能力。

Brog[30]在 Boeing/AFOSR Mach-6靜音風(fēng)洞中采用 20% 的縮比模型對X-51A進(jìn)氣道轉(zhuǎn)捩情況進(jìn)行了一系列試驗研究。主要研究了噪聲對光滑壁面的自然轉(zhuǎn)捩和粗糙顆粒誘導(dǎo)的強(qiáng)制轉(zhuǎn)捩的影響。模型如圖19所示。模型總長34.42cm。為評估本文的模型對復(fù)雜構(gòu)型轉(zhuǎn)捩的模擬能力，計算了3個算例，并將計算結(jié)果與試驗數(shù)據(jù)進(jìn)行了對比。計算設(shè)置如表3所示。由于沒有偏航角，為節(jié)省計算量，只計算了一半的模型，在中心面上采用對稱邊界條件。光滑壁面模型的網(wǎng)格為400萬左右，安裝了鉆石形轉(zhuǎn)捩帶的模型網(wǎng)格量為600萬左右。壁面第一層網(wǎng)格y+小于 1。試驗結(jié)果顯示壁面溫度在整個試驗過程中變化很小，因此將壁溫設(shè)為300K。湍流度在噪聲條件下為3.0%，靜音條件下為0.05%。由于3個算例的馬赫數(shù)和雷諾數(shù)都非常接近，其壁面摩阻和熱流差別也很小，因此在計算全層流和全湍流結(jié)果時，只計算了case5，并用其來判斷3個算例的轉(zhuǎn)捩位置。

圖19 20%縮比X-51A前體構(gòu)型示意圖

CaseMaRe/m-1總壓/kPa總溫/K來流湍流度轉(zhuǎn)捩類型46.006.59×1065864180.05%自然55.807.40×1066214243.00%自然65.787.35×1066144243.00%強(qiáng)制

圖20為噪聲情況下模型迎風(fēng)面測量得到的壁面溫度與計算的摩阻系數(shù)。試驗與計算結(jié)果均在壓縮拐角之后開始上升，但測量的壁面溫度顯示轉(zhuǎn)捩位置在中間最靠后，在兩側(cè)區(qū)域最靠前；而計算結(jié)果顯示轉(zhuǎn)捩區(qū)域在展向分布較為均勻。從圖21可知，安裝了轉(zhuǎn)捩帶后，測量的壁面溫度在壓縮拐角之前就明顯大于光滑壁面情況，計算得到的摩阻系數(shù)則在轉(zhuǎn)捩帶之后立即上升到較高值，表明轉(zhuǎn)捩位置位于壓縮拐角和轉(zhuǎn)捩帶之間。

圖20 迎風(fēng)面測量得到的壁面溫度和計算得到的摩阻系數(shù)(case5)

Fig.20Measuredsurfacetemperatureandcalculatedskinfrictionforwindwardsurface(case5)

圖22為迎風(fēng)面 case4～case6的壁面熱流及 case5的全層流和全湍流情況。靜音條件下，熱流分布與層流情況基本一致，表示轉(zhuǎn)捩并沒有發(fā)生；噪聲情況下，拐角之后出現(xiàn)了明顯的熱流抬升，即轉(zhuǎn)捩現(xiàn)象。而安裝了鉆石形轉(zhuǎn)捩帶之后，轉(zhuǎn)捩位置被推到了轉(zhuǎn)捩帶和壓縮拐角之間的位置。熱流結(jié)果與摩阻系數(shù)得到的結(jié)論吻合。摩擦力線則顯示：在層流和靜音模式下，由于壓力梯度的作用，進(jìn)氣道下表面有明顯的橫流溢流，會導(dǎo)致進(jìn)入進(jìn)氣道的流量降低；而對于全湍流或者轉(zhuǎn)捩的狀況，這樣的橫向流動得到了較為明顯的抑制，有利于提升進(jìn)氣道捕獲流量。在圖23中，噪聲條件下，計算的中心線熱流顯示轉(zhuǎn)捩一經(jīng)過壓縮拐角就開始，而在x=0.28m的位置完成。測量的壁面溫度顯示轉(zhuǎn)捩起始于壓縮拐角，結(jié)束于x=0.27m的位置。而靜音條件下計算與試驗結(jié)果都顯示沒有出現(xiàn)轉(zhuǎn)捩現(xiàn)象。在靜音和噪聲情況下，試驗與計算得到的轉(zhuǎn)捩起始和結(jié)束位置均吻合較好。圖24中，在中心線位置，計算的摩阻系數(shù)和測量的壁面溫度在壓縮拐角的位置均已經(jīng)超過了光滑壁面的情況，進(jìn)一步表明邊界層在壓縮拐角之前就已經(jīng)轉(zhuǎn)捩為湍流。從計算結(jié)果可知，轉(zhuǎn)捩帶將轉(zhuǎn)捩完成位置前推到了距前緣0.08m的位置。

圖21 迎風(fēng)面測量得到的壁面溫度和計算得到的摩阻系數(shù)(case6)

Fig.21Measuredsurfacetemperatureandcalculatedskinfrictionforwindwardsurface(case6)

圖22 計算得到的迎風(fēng)面熱流對比

(a) 測量得到的壁面溫度

(b) 計算得到的摩阻系數(shù)

圖23 中心線上測量得到的壁面溫度和計算得到的摩阻系數(shù)(case4, case5)

Fig.23Measuredsurfacetemperatureandcalculatedskinfrictionofthecenterline(case4,case5)

(a) 測量得到的壁面溫度

(b) 計算得到的摩阻系數(shù)

圖24 中心線上測量得到的壁面溫度和計算得到的摩阻系數(shù)(case6)

Fig.24Measuredsurfacetemperatureandcalculatedskinfrictionofthecenterline(case6)

3 結(jié) 論

在 RANS 求解器框架下，通過添加可壓縮性修正，修改分離誘導(dǎo)轉(zhuǎn)捩參數(shù)，并對經(jīng)驗耦合關(guān)系式進(jìn)行修改，搭建了針對高超聲速情況下的γ-Reθ轉(zhuǎn)捩模型。用原始的γ-Reθ模型和修正后的模型對Ma7.4的 Ames 全尺寸模型和Ma6的單個粗糙顆粒誘導(dǎo)平板轉(zhuǎn)捩進(jìn)行了計算。 結(jié)果顯示：原始模型在計算自然轉(zhuǎn)捩時，轉(zhuǎn)捩位置過于靠前，而且不能捕捉粗糙顆粒對轉(zhuǎn)捩的促進(jìn)作用；而修正模型克服了原始模型的這兩個缺點，得到了與試驗比較一致的結(jié)果。最后，將修正后的模型應(yīng)用于20%縮比X-51A前體進(jìn)氣道模型的自然轉(zhuǎn)捩和強(qiáng)制轉(zhuǎn)捩計算，結(jié)果顯示模型能夠反映高湍流度和粗糙顆粒對轉(zhuǎn)捩的促進(jìn)作用，表明該模型有望應(yīng)用于高超聲速復(fù)雜構(gòu)型自然轉(zhuǎn)捩和強(qiáng)制轉(zhuǎn)捩的預(yù)測研究。

當(dāng)然，由于該模型目前僅將來流湍流度作為來流擾動條件，而高超聲速來流擾動中的湍流度、擾動頻率、擾動模態(tài)均可影響轉(zhuǎn)捩位置，因此，如何使模型更有效反映來流擾動，特別是來流壓力擾動發(fā)展成為第二模態(tài)后對轉(zhuǎn)捩的影響，有待進(jìn)一步研究。