王東方，鄧小剛，王光學(xué),2，劉化勇

(1.國防科技大學(xué) 航天科學(xué)與工程學(xué)院， 湖南 長沙　410073；2.中山大學(xué) 物理學(xué)院， 廣東 廣州　510006；3.中國空氣動力研究與發(fā)展中心， 四川 綿陽　621000)

?

高超聲速尖雙錐流動高精度數(shù)值模擬*

王東方1，鄧小剛1，王光學(xué)1,2，劉化勇3

(1.國防科技大學(xué) 航天科學(xué)與工程學(xué)院， 湖南 長沙410073；2.中山大學(xué) 物理學(xué)院， 廣東 廣州510006；3.中國空氣動力研究與發(fā)展中心， 四川 綿陽621000)

以25°/55°尖雙錐外形的高超聲速低焓層流流動模擬為例，對高階加權(quán)緊致非線性格式模擬激波/邊界層干擾流動的能力進(jìn)行驗(yàn)證和確認(rèn)。空間離散采用二階MUSCL和三階、五階加權(quán)緊致非線性格式，時間離散采用二階精度雙時間步方法，通量函數(shù)采用混合Roe-Rusanov，AUSMPW+，Van Leer等，對比了不同精度空間離散格式對時間、網(wǎng)格收斂特性和通量函數(shù)耗散特性的影響。數(shù)值模擬結(jié)果表明采用高精度空間離散格式能在較疏的網(wǎng)格上獲得收斂解，并能消除結(jié)果對通量函數(shù)的敏感性，但收斂需要推進(jìn)更久的計算時間。數(shù)值模擬結(jié)果與實(shí)驗(yàn)測量結(jié)果吻合良好，滿足工程精度要求。

驗(yàn)證與確認(rèn)；激波/邊界層干擾；高精度方法；高超聲速流動

激波/邊界層干擾(Shock Wave/Boundary Layer Interaction，SWBLI)在高超聲速飛行器流場中普遍存在并對飛行器的氣動力、熱特性產(chǎn)生重要影響，因而在過去幾十年得到了廣泛研究[1]。干擾激波會在飛行器表面產(chǎn)生很高的局部氣動熱載荷，可能破壞熱防護(hù)系統(tǒng)。激波干擾會誘導(dǎo)邊界層分離，改變飛行器氣動特性，限制其機(jī)動能力，若分離發(fā)生在控制舵前方將會降低控制舵效率。因而在高超聲速飛行器設(shè)計中需要重點(diǎn)考慮激波/邊界層干擾現(xiàn)象。

未來高超聲速飛行器的發(fā)展依賴高保真的數(shù)值模擬工具進(jìn)行設(shè)計和優(yōu)化。目前的計算流體力學(xué)(Computational Fluid Dynamics，CFD)方法已能預(yù)測真實(shí)高超聲速飛行中的大部分流動現(xiàn)象，但準(zhǔn)確模擬SWBLI對CFD來說依然是一個嚴(yán)峻挑戰(zhàn)[1]。因此，十多年來國內(nèi)外一直在持續(xù)不斷地開展CFD預(yù)測SWBLI的確認(rèn)研究活動[2-3]。

25°/55°尖雙錐外形是SWBLI確認(rèn)研究的基準(zhǔn)外形之一，廣泛用于考核數(shù)值方法和物理模型，國際上已對其開展了系統(tǒng)的研究。Holden等[4-6]對該外形進(jìn)行了系列精細(xì)實(shí)驗(yàn)，為CFD確認(rèn)工作提供了高精度的實(shí)驗(yàn)數(shù)據(jù)。通過精心地調(diào)整自由流總焓、雷諾數(shù)和實(shí)驗(yàn)氣體介質(zhì)，突出不同的物理現(xiàn)象。2000年開展的盲比確認(rèn)活動重點(diǎn)關(guān)注沒有化學(xué)反應(yīng)的層流分離流動[4,7-8]。2013年開展的盲比確認(rèn)活動關(guān)注真實(shí)氣體效應(yīng)對層流區(qū)壁面力、熱特性的影響，考核CFD中化學(xué)反應(yīng)和振動松弛等物理模型的準(zhǔn)確性[5]。

這些研究顯示雙錐流動模擬對數(shù)值方法提出的挑戰(zhàn)主要有3個方面：數(shù)值耗散、時間積分和物理模型。數(shù)值耗散主要體現(xiàn)為數(shù)值模擬結(jié)果對空間格式精度、通量函數(shù)和網(wǎng)格分辨率的敏感程度，采用高分辨率低耗散的數(shù)值方法更利于獲得準(zhǔn)確的流場結(jié)構(gòu)。高效的時間積分方法用于解決激波/激波干擾和分離區(qū)的耦合振蕩導(dǎo)致的收斂困難[9]。在高焓條件下，流場中存在熱化學(xué)非平衡，現(xiàn)已開展了很多研究來確認(rèn)和改進(jìn)熱化學(xué)模型[10-12]。總的來說，來流焓值越高，CFD模擬結(jié)果與實(shí)驗(yàn)結(jié)果差異越大，尤其是實(shí)驗(yàn)介質(zhì)包含氧氣的時候。此外壁面催化也對熱流結(jié)果產(chǎn)生重要影響。

Druguet等[13]基于二階精度算法系統(tǒng)研究了數(shù)值方法對雙錐流場模擬的影響，發(fā)現(xiàn)模擬結(jié)果對采用的通量函數(shù)和斜率限制器很敏感。粗網(wǎng)格結(jié)果顯示分離區(qū)大小強(qiáng)烈依賴于數(shù)值方法的耗散，低耗散的通量函數(shù)和斜率限制器模擬的分離區(qū)更大。網(wǎng)格收斂性研究發(fā)現(xiàn)分離區(qū)大小對網(wǎng)格加密程度很敏感。因此，雙錐流動模擬應(yīng)采用高分辨率低耗散的數(shù)值方法，同時要進(jìn)行仔細(xì)的網(wǎng)格收斂性研究以確保網(wǎng)格分辨率足夠。

Gaitonde等[14]的時間收斂性研究表明雙錐流動模擬至少要推進(jìn)100個特征時間長度才能達(dá)到收斂。為確保計算收斂，除了監(jiān)測殘差，還要監(jiān)測分離區(qū)大小和氣動力系數(shù)。隱式時間方法常用來加速收斂，應(yīng)注意局部時間步長可能會導(dǎo)致非物理解，此外較大的時間步長可能使計算發(fā)生振蕩而不能收斂[9,13]。

本文重點(diǎn)考察數(shù)值方法因素對雙錐流動模擬的影響。數(shù)值模擬結(jié)果對限制器、通量函數(shù)以及網(wǎng)格敏感，其主要與二階精度方法耗散較大而不能準(zhǔn)確模擬分離區(qū)大小有關(guān)。高精度方法分辨率高耗散低，對分離渦等多尺度結(jié)構(gòu)具有良好的捕捉能力。Tissera[15]采用高階加權(quán)本質(zhì)無振蕩(Weighted Essentially NonOscillatory,WENO)格式模擬尖雙錐流動的研究表明采用高精度方法比加密網(wǎng)格能更高效地改善模擬結(jié)果。因而有必要對高精度方法模擬雙錐流動的特性進(jìn)行系統(tǒng)研究。

近年來，高精度方法研究已受到廣泛關(guān)注，其在直接數(shù)值模擬(Direct Numerical Simulations，DNS)、大渦模擬(Large Eddy Simulations，LES)、計算聲學(xué)(Computational AeroAcoustics，CAA)等領(lǐng)域獲得了大量應(yīng)用[16]。但應(yīng)用于高超聲速流動模擬時，高精度方法必須具備健壯的激波捕捉能力。常通過引入加權(quán)技術(shù)構(gòu)造非線性的高精度方法，如WENO格式[17]和加權(quán)緊致非線性格式(Weighted Compact Nonlinear Schemes, WCNS)[18]。WCNS由Deng等[18]基于緊致非線性格式(Compact Nonlinear Schemes，CNS)[19]構(gòu)造，目前已成功應(yīng)用于各種類型亞跨聲速復(fù)雜外形流動模擬[20-22]，其不但能滿足幾何守恒律(Geometric Conservation Law，GCL)[23-24]，還具有曲線網(wǎng)格下的自由流和渦保持特性[25]。典型的高超聲速流動模擬表明WCNS還具有健壯的激波捕捉能力，本文重點(diǎn)考核其在高超聲速激波/邊界層干擾流動模擬中的特性。

1　數(shù)值方法

1.1控制方程

在笛卡爾坐標(biāo)系下，不考慮體積力和外部熱源，無量綱形式的三維可壓縮Navier-Stokes方程為：

(1)

式中：Q=[ρ,ρu,ρv,ρw,ρe]T；E,F,G和Eν,Fν,Gν為無黏和黏性通量。方程式(1)轉(zhuǎn)換至曲線坐標(biāo)系下變?yōu)椋?/p>

(2)

(3)

1.2空間格式

空間離散包含三部分：加權(quán)插值、數(shù)值通量計算、通量差分。由于顯式差分格式計算高效、便于矢量化和并行化，Deng等[27]建議在通量差分采用顯式差分格式。以ξ方向離散為例說明本文的空間離散方法。

四階顯式差分格式WCNS-E4為：

(4)

六階顯式差分格式WCNS-E6為：

(5)

(6)

(7)

(8)

式(7)可表達(dá)成式(8)的線性組合：

(9)

式中，

(10)

稱為最優(yōu)權(quán)或線性權(quán)。

非線性插值是將三個子模板的線性插值進(jìn)行非線性加權(quán)組合：

(11)

在光滑區(qū)，非線性權(quán)應(yīng)趨于線性權(quán)以達(dá)到最優(yōu)精度；在間斷附近，包含間斷的子模板的非線性權(quán)應(yīng)趨于零以避免跨間斷插值。 依此原則定義的非線性權(quán)為：

(12)

式中，ε=10-6是避免分母為零的小量，βk為各子模板的光滑因子：

(13)

數(shù)值模擬發(fā)現(xiàn)式(12)所示的權(quán)函數(shù)耗散較大，并且在光滑區(qū)極值點(diǎn)處不能恢復(fù)為最優(yōu)權(quán)。 為此，Liu等[28]在WCNS中引入了Borges等[29]提出的最新的權(quán)函數(shù)：

(14)

(15)

則三階非線性插值為：

(16)

非線性權(quán)定義為：

(17)

式(4)和式(16)標(biāo)記為WCNS-E4E3，式(5)和式(11)標(biāo)記為WCNS-E6E5。為了對比，還采用了二階迎風(fēng)格式MUSCL(monotone upstream-centred scheme for conservation laws)和Minmod限制器，標(biāo)記為MUSCL-MM。三種空間格式的近似色散分析(Approximate Dispersion Relation，ADR)[31-32]結(jié)果如圖1所示(其中φ為波數(shù)；ΦR,ΦI分別為修正波數(shù)的實(shí)部和虛部，分別反映格式的色散和耗散特性)，可以看出隨著格式精度提高，格式耗散特性明顯改善。

圖1　不同格式的近似色散分析結(jié)果Fig.1　Approximate dispersion relation for various schemes

1.3時間格式和邊界條件

根據(jù)文獻(xiàn)[1]中Candler的建議，采用二階精度雙時間步方法[33]，子迭代采用隱式LU-SGS方法[34]。

本文將軸對稱的尖雙錐流動按三維問題求解。首先生成二維結(jié)構(gòu)網(wǎng)格，再沿對稱軸旋轉(zhuǎn)生成三維計算域，周向取5個計算單元，因此兩側(cè)取對稱邊界條件。上邊界取超聲速入口條件，右側(cè)出口按零梯度外推處理。壁面取無滑移等溫壁條件。

2　計算條件

2.1雙錐幾何參數(shù)

雙錐模型幾何參數(shù)如圖2所示。前面是兩個半錐角分別為25°和55°的錐體，后面是圓柱體。

圖2　雙錐模型幾何參數(shù)Fig.2　Double cone configuration

2.2自由流條件

選取Holden等[4-5]的兩個低焓實(shí)驗(yàn)狀態(tài)作為自由流條件，以最小化真實(shí)氣體效應(yīng)的影響。具體參數(shù)如表1所示，其中u∞,ρ∞,T∞分別為實(shí)驗(yàn)段來流的速度、密度和溫度，Twall為實(shí)驗(yàn)?zāi)Ｐ捅砻鏈囟取?/p>

表1　雙錐實(shí)驗(yàn)自由流條件

2.3網(wǎng)格

雙錐流動對網(wǎng)格很敏感，因而生成了四套逐漸加密的網(wǎng)格進(jìn)行網(wǎng)格收斂性研究(流向×法向)：128×128，256×128，512×256，1024×256。網(wǎng)格分布如圖3所示。邊界層內(nèi)加密以更好刻畫邊界層和壁面熱流特性，流向上在錐前緣、分離和再附點(diǎn)處都適當(dāng)加密。對512×256的網(wǎng)格，近壁第一層網(wǎng)格高度為1×10-5m，計算結(jié)果顯示所有位置的y+都小于1，滿足氣動熱計算要求。

總體來說，三種格式基于1024×256網(wǎng)格模擬的流場特性均達(dá)到網(wǎng)格收斂性要求。但對Run28狀態(tài)，二階格式?jīng)]有明顯的網(wǎng)格收斂現(xiàn)象，仍需要進(jìn)一步加密網(wǎng)格。

圖3　256×128網(wǎng)格分布Fig.3　Structure of grid 256×128

3　結(jié)果分析

本節(jié)先討論雙錐流場結(jié)構(gòu)，再分析不同空間格式對時間和網(wǎng)格收斂性的影響，隨后分析計算效率和通量函數(shù)的影響，最后與實(shí)驗(yàn)數(shù)據(jù)對比評估高精度方法的準(zhǔn)確性。

為便于和實(shí)驗(yàn)對比將數(shù)據(jù)進(jìn)行無量綱化。長度按前錐軸向長度L=92.08 mm無量綱化，壁面壓力p和熱流率Qw特性按式(18)無量綱化：

(18)

3.1流場結(jié)構(gòu)分析

這里結(jié)合高精度數(shù)值模擬結(jié)果對雙錐流場結(jié)構(gòu)做系統(tǒng)總結(jié)。以Run28狀態(tài)下五階WCNS-E6E5格式和AUSMPW+通量函數(shù)在1024×256網(wǎng)格上計算的結(jié)果為例分析雙錐流場結(jié)構(gòu)，數(shù)值紋影結(jié)果如圖4所示。

圖4　Run28流場結(jié)構(gòu)Fig.4　Flowfield structure for Run28

前錐產(chǎn)生的附著激波與后錐產(chǎn)生的脫體激波相互干擾并形成局部高壓區(qū)，逆壓梯度在兩錐交界附近形成分離區(qū)，分離區(qū)將形成分離激波。前錐激波與分離激波干擾將強(qiáng)化分離激波并形成接觸間斷。分離激波與后錐脫體激波干擾屬于V 類激波干擾，在三叉點(diǎn)處形成入射激波和接觸間斷。該接觸間斷將后錐脫體激波后的亞聲速區(qū)和近壁附近的超聲速區(qū)分開。分離區(qū)上方的區(qū)分流線有輕微下凹，其形成的壓縮波在激波層中間匯聚，與三叉點(diǎn)處發(fā)出的入射激波干擾形成馬赫盤，其后方是兩道入射激波包圍著的高壓亞聲速區(qū)。上傳的入射激波與三叉點(diǎn)處發(fā)出的接觸間斷干擾反射為膨脹波，下傳的入射激波在分離區(qū)尾部終止。分離區(qū)尾部的亞聲速區(qū)快速膨脹，促使分離區(qū)上方流動加速沖向壁面，壓縮邊界層導(dǎo)致局部壓力和熱流劇烈增加，隨后形成壁面再壓縮激波，其在后錐上方超聲速區(qū)域內(nèi)不斷反射：壓縮波從壁面反射為壓縮波，從接觸間斷反射為膨脹波；膨脹波從壁面反射為膨脹波，從接觸間斷反射為壓縮波?？偟牧鲃咏Y(jié)構(gòu)與Gnoffo等[9,11,35]分析的結(jié)果一致，但熱流和壓力峰值產(chǎn)生的原因是分離區(qū)再附點(diǎn)后的流動沖擊壁面導(dǎo)致的，而非三叉點(diǎn)處的入射激波導(dǎo)致。

3.2時間收斂性

本小節(jié)討論空間格式和網(wǎng)格因素對時間收斂特性的影響。根據(jù)文獻(xiàn)[1]中Candler的建議采用二階時間精度的雙時間步方法計算以避免非物理振蕩。Druguet等[13-14]將以自由流速度經(jīng)過特征長度的時間定義為特征時間(tchar=l/u∞)，以非定常計算推進(jìn)時間與特征時間的比值TC來衡量時間收斂性，這里特征長度取0.2 m。通過檢測殘差和阻力系數(shù)判斷收斂。

對Run28和Run1狀態(tài)，基于1024×256網(wǎng)格分別采用不同精度的空間格式進(jìn)行計算，通量函數(shù)統(tǒng)一采用混合Roe-Rusanov。阻力系數(shù)的收斂歷史如圖5所示。對Run28狀態(tài)，二階MUSCL格式收斂需要推進(jìn)150個特征時間，而三階和五階WCNS需要推進(jìn)250個以上特征時間。對Run1狀態(tài)，三種格式推進(jìn)100個特征時間都能收斂，彼此差異不大。這說明計算需要的收斂時間強(qiáng)烈依賴流動狀態(tài)，采用高階格式時推進(jìn)需更多的時間。對Run1狀態(tài)還對比了三階WCNS在不同網(wǎng)格上的收斂特性，如圖6所示。由圖可見密網(wǎng)格收斂需要更多的計算時間。流動收斂的過程也是分離區(qū)演化至穩(wěn)態(tài)的過程，提高格式精度和加密網(wǎng)格都會減小數(shù)值耗散，不利于消除流場中的擾動，因而需要更久的時間才能收斂。

(a) Run28

(b) Run1圖5　阻力系數(shù)收斂歷史Fig.5　Drag coefficient convergence histories

圖6　不同網(wǎng)格的阻力收斂特性Fig.6　Drag coefficient convergence characteristics on different grids

3.3網(wǎng)格收斂性

本小節(jié)討論空間格式對網(wǎng)格收斂性的影響。通過對比壁面特性隨網(wǎng)格的變化來判斷網(wǎng)格收斂性，這里只給出壁面熱流特性。Run28和Run1狀態(tài)不同空間格式的網(wǎng)格收斂性分別如圖7、圖8所示?？梢婋S著網(wǎng)格加密，分離區(qū)逐漸增大，干擾點(diǎn)產(chǎn)生的峰值熱流位置逐漸后移。

綜合壁面熱流特性和分離區(qū)大小進(jìn)行分析可見，格式精度越高越容易在較粗的網(wǎng)格上獲得收斂。對Run1狀態(tài)，二階格式在1024×256網(wǎng)格上收斂，三階格式在512×256網(wǎng)格上收斂，五階格式在256×128網(wǎng)格上收斂。Run28狀態(tài)較難收斂，就五階格式在512×256網(wǎng)格上結(jié)果而言，盡管分離區(qū)大小還在變化，但壁面熱流已看不出明顯差異，可認(rèn)為達(dá)到網(wǎng)格收斂。二階格式在1024×256網(wǎng)格上計算的分離區(qū)大小與高階收斂結(jié)果相比明顯偏小，因而認(rèn)為未達(dá)到網(wǎng)格收斂，還需進(jìn)一步加密網(wǎng)格。

3.4計算效率

對推進(jìn)10 000步的計算時間統(tǒng)計表明，WCNS-E4E3和WCNS-E6E5要比MUSCL-MM格式分別多消耗5%和25%的CPU時間。若以二階格式基于128×128網(wǎng)格推進(jìn)一個特征時間所需的CPU時間為標(biāo)準(zhǔn)，將各階格式獲得收斂結(jié)果所需的CPU時間歸一化，結(jié)果如表2所示，可見提高格式精度能更高效得到收斂解。

(a) 二階格式(a) MUSCL-MM

(b) 三階格式(b) WCNS-E4E3

(c) 五階格式(c) WCNS-E6E5圖7　Run28壁面熱流網(wǎng)格收斂性Fig.7　Grid convergence of heat transfer rate for Run28

(a) 二階格式(a) MUSCL-MM

(b) 三階格式(b) WCNS-E4E3

(c) 五階格式(c) WCNS-E6E5圖8　Run1壁面熱流網(wǎng)格收斂性Fig.8　Grid convergence of heat transfer rate for Run1

3.5通量函數(shù)影響

通量函數(shù)的耗散是數(shù)值耗散的重要部分。本小節(jié)評估采用不同精度空間格式對通量函數(shù)的影響。在512×256網(wǎng)格上對Run28狀態(tài)采用不同的空間格式和通量函數(shù)組合進(jìn)行模擬，如Van Leer通量[36]、Mod.SW通量[13]、AUSMPW+通量[37]、混合Roe-Rusanov通量[38]。分離區(qū)大小如表3所示?？梢娡亢瘮?shù)耗散越小，相應(yīng)的分離區(qū)越大，與Druguet等[13]的結(jié)論一致。但是不同通量間分離區(qū)大小的差異隨著空間格式精度的提高而越來越小。

五階WCNS和不同通量函數(shù)計算的壁面熱流率特性如圖10所示。可見不同通量函數(shù)計算的峰值熱流幾乎一致。這說明空間格式精度的提高使計算結(jié)果對通量函數(shù)的敏感性大大降低。

3.6結(jié)果確認(rèn)

三階和五階WCNS都計算得到了網(wǎng)格收斂的結(jié)果，壁面壓力和熱流率特性與實(shí)驗(yàn)數(shù)據(jù)[4,6]對比如圖11、圖12所示。

總的來說，前錐和分離區(qū)壓力特性與實(shí)驗(yàn)符合很好。前錐的壁面熱流率特性與實(shí)驗(yàn)有差異。Run28狀態(tài)計算的干擾點(diǎn)比實(shí)驗(yàn)靠后，峰值大小與實(shí)驗(yàn)接近。Run1狀態(tài)計算的干擾點(diǎn)位置與實(shí)驗(yàn)差異不大，但峰值偏大。

(a) Run28

(b) Run1圖9　不同網(wǎng)格分離區(qū)大小對比Fig.9　Size of separation zone on different grids

MUSCL-MMWCNS-E4E3WCNS-E6E5Run281×32×150=48001.05×16×300=50401.25×8×300=3000Run11×16×100=16001.05×8×150=12601.25×2×150=375

表3　不同通量函數(shù)計算的分離區(qū)大小

圖10　不同通量壁面熱流率特性對比Fig.10　Wall heat transfer rate computed with different flux functions

(a) 壁面壓力特性(a) Surface pressure

(b) 壁面熱流率特性(b) Heat transfer rate圖11　Run28狀態(tài)計算與實(shí)驗(yàn)結(jié)果對比Fig.11　Comparison between experimental measurements and numerical predictions for Run28

(a) 壁面壓力特性(a) Surface pressure

(b) 壁面熱流率特性(b) Heat transfer rate圖12　Run1狀態(tài)計算與實(shí)驗(yàn)結(jié)果對比Fig.12　Comparison between experimental measurements and numerical predictions for Run1

4　結(jié)論

采用二階MUSCL、三階和五階WCNS格式模擬了高超聲速25°/55°尖雙錐流場，研究了空間格式對時間收斂性、網(wǎng)格收斂性和通量耗散特性的影響。隨著空間格式精度提高，收斂需要推進(jìn)的計算特征時間會增加，需要的網(wǎng)格量會減少，對通量函數(shù)的敏感性會降低。不同來流狀態(tài)收斂需要的計算時間和網(wǎng)格不一致，這可能與流場中的分離區(qū)大小有關(guān)。分離區(qū)越大，計算越難收斂，網(wǎng)格要求越高。因而在推廣到相關(guān)問題模擬時仍有必要開展時間和網(wǎng)格收斂性研究。

高階格式盡管需要推進(jìn)更多的特征時間才能獲得收斂解，但需要的網(wǎng)格量少會使計算量減小，綜合來看高階格式獲得收斂解的計算效率還是優(yōu)于二階格式的。

數(shù)值模擬結(jié)果與實(shí)驗(yàn)對比整體符合較好，滿足工程精度要求。但干擾區(qū)峰值特性預(yù)測還與實(shí)驗(yàn)存在較大差異，需要進(jìn)一步研究。

References)

[1]Babinsky H, Harver J K. Shock wave-boundary-layer interactions[M].USA:Cambridge University Press, 2011.

[2]Kinight D. RTO WG10: test cases for CFD validation of hypersonic flight[C]//Proceedings of 40th AIAA Aerospace Sciences Meeting & Exhibit, 2006.

[3]Gaitonde D V. Progress in shock wave/boundary interactions[J]. Progress in Aerospace Sciences, 2015, 72: 80-99.

[4]Harvey J K, Holden M S, Wadhams T P. Code validation study of laminar shock/boundary layer and shock/shock interactions in hypersonic flow. Part a: experimental measurements[C]//Proceedings of 39th Aerospace Sciences Meeting and Exhibit, AIAA 2001-1031, 2001.

[5]Holden M S, Wadhams T P, MacLean M G, et al. Measurements of real gas effects on regions of laminar shock wave/boundary layer interaction in hypersonic flows for blind code validation studies [C]//Proceedings of 21st AIAA Computational Fluid Dynamics Conference, AIAA 2013-2837, 2013.

[6]MacLean M G, Holden M S, Dufrene A. Measurements of real gas effects on regions of laminar shock wave/boundary layer interaction in hypervelocity flows[C]. AIAA Aviation 2014, Atlanta, GA, 2014.

[7]Candler G V, Nompelis I, Holden M S. Computational analysis of hypersonic laminar viscous-inviscid interactions[C]//Proceedings of 38th Aerospace Sciences Meeting and Exhibit, AIAA 2000-0532, 2000.

[8]Harvey J K, Holden M S, Wadhams T P. Code validation study of laminar shock/boundary layer and shock/shock interaction in hypersonic flow. Part b: comparison with Navier-Stokes and DSMC solutions [C]//Proceedings of 39th Aerospace Sciences Meeting and Exhibit, AIAA 2001-1031, 2001.

[9]Nompelis I. Computational study of hypersonic double-cone experiments for code validation[D]. USA: University of Minnesota, 2004 .

[10]Nompelis I, Candler G V, MacLean M, et al. Numerical investigation of double-cone flow experiments with high-enthalpy effects [C]//Proceedings of 48th AIAA Aerospace Sciences Meeting Including the New Horizons Forum and Aerospace Exposition, AIAA 2010-1283, 2010.

[11]Knight D, Longo J, Drikakis D, et al. Assessment of CFD capability for prediction of hypersonic shock interactions[J].Progress in Aerospace Science, 2012, 48: 8-26.

[12]Nompelis I, Candler G V. US3d predictions of double-cone and hollow cylinder-flare flows at high enthalpy[C]//Proceedings of 44th AIAA Fluid Dynamics Conference, AIAA 2014-3366, 2014.

[13]Druguet M C, Candler G V, Nompelis I. Effects of numerics on Navier-Stokes computations of hypersonic double-cone flows[J]. AIAA Journal, 2005, 43(3): 616-623.

[14]Gaitonde D V, Canupp P W, Holden M S. Heat transfer predictions in a laminar hypersonic viscous/inviscid interaction [J]. Journal of Thermophysics and Heat Transfer, 2002, 16(4): 481-489.

[15]Tissera S. Assessment of high-resolution methods in hypersonic real-gas flows[D].USA: Cranfield University, 2010.

[16]Wang Z J, Fidkowski K, Abgrall R, et al. High-order CFD methods: current status and perspective[J]. International Journal for Numerical Methods in Fluids, 2013, 72(8): 811-845.

[17]Jiang G S, Shu C W. Efficient implementation of weighted ENO schemes[J]. Journal of Computational Physics, 1996, 126: 202-228.

[18]Deng X G, Zhang H X. Developing high-order weighted compact nonlinear schemes [J]. Journal of Computational Physics, 2000, 165(1): 22-44.

[19]Deng X G, Maekawa H. Compact high-order accurate nonlinear schemes[J].Journal of Computational Physics, 1997, 130(1): 77-91.

[20]Liu X, Deng X G, Mao M L. High-order behaviors of weighted compact fifth-order nonlinear schemes [J]. AIAA Journal, 2007, 45(8): 2093-2097.

[21]Deng X G, Mao M L, Tu G H, et al. Extending weighted compact nonlinear schemes to complex grids with characteristic-based interface conditions[J].AIAA Journal, 2010, 48(12): 2840-2851.

[22]Deng X G, Mao M L, Tu G H, et al. High-order and high accurate CFD methods and their applications for complex grid problems [J].Communications in Computational Physics, 2012, 11(4): 1081-1102.

[23]Deng X G, Mao M L, Tu G H, et al. Geometric conservation law and applications to high-order finite difference schemes with stationary grids[J]. Journal of Computational Physics, 2011, 230(4): 1100-1115.

[24]Deng X G, Min Y B, Mao M L, et al. Further studies on geometric conservation law and applications to high-order finite difference schemes with stationary grids[J]. Journal of Computational Physics, 2013, 239: 90-111.

[25]Nonomura T, Iizuka N, Fujii K. Freestream and vortex preservation properties of high-order WENO and WCNS on curvilinear grids[J]. Computer & Fluids, 2010, 39(2): 197-214.

[26]Deng X G, Zhu H J, Min Y B, et al. Symmetric conservation metric method: a link between high order finite-difference and finite-volume schemes for flow computations around complex geometries[C]//Proceedings of Eighth International Conference on Computational Fluid Dynamics (ICCFD8), ICCFD8-2014-0005, Chengdu, China, 2014.

[27]Deng X G, Liu X, Mao M L, et al. Investigation on weighted compact fifth-order nonlinear scheme and applications to complex flows [C]//Proceedings of 17th AIAA Computational Fluid Dynamics Conference, AIAA 2005-5246, 2005.

[28]Liu H Y, Ma Y K, Yan Z G, et al. A shock-capturing methodology based on high order compact interpolation[C]//Proceedings of Eighth International Conference on Computational Fluid Dynamics (ICCFD8), ICCFD8-2014-0082, Chengdu, China, 2014.

[29]Borges R, Carmona M, Costa B, et al. An improved weighted essentially non-oscillatory scheme for hypersonic conservation laws[J]. Journal of Computational Physics, 2008, 227(6): 3191-3211.

[30]燕振國， 劉化勇， 毛枚良， 等. 三階HWCNS格式的構(gòu)造及在高超流動中的應(yīng)用[J].航空學(xué)報， 2015, 36(5): 1460-1470.

YAN Zhenguo, LIU Huayong, MAO Meiliang, et al. Development of 3rd order HWCNS schemes and its application in hypersonic flow simulation[J]. Acta Aeronautica et Astronautica Sinica, 2015, 36(5): 1460-1470. (in Chinese)

[31]Pirozzoli S. On the spectral properties of shock-capturing schemes[J].Journal of Computational Physics, 2006, 219(2): 489-497.

[32]毛枚良， 燕振國， 劉化勇， 等. 高階加權(quán)非線性格式的擬線性頻譜分析方法研究[J]. 空氣動力學(xué)學(xué)報， 2015, 33(1): 1-9.

MAO Meiliang, YAN Zhenguo, LIU Huayong, et al. Study of quasi-linear spectral analysis method of high-order weighted nonlinear schemes[J].Acta Aerodynamica Sinica, 2015, 33(1): 1-9.(in Chinese)

[33]Pulliam T H. Time accuracy and the use of implicit methods[C]//Proceedings of 11th Computational Fluid Dynamics Conference, AIAA 1993-3360, 1993.

[34]Yoon S, Jameson A. Lower-upper symmetric-Gauss-Seidel method for the Euler and Navier-Stokes equations[J]. AIAA Journal, 1988, 26(9): 1025-1026.

[35]Gnoffo P A. CFD validation studies for hypersonic flow prediction [C] //Proceedings of 39th Aerospace Sciences Meeting and Exhibit, AIAA 2001-1025, 2001.

[36]Van Leer B. Flux-vector splitting for the Euler equations[C]//Proceedings of Eighth International Conference on Numerical Methods in Fluid Dynamics, Springer, 1982: 507-512 .

[37]Kim K H, Kim C, Rho O H. Methods for the accurate computations of hypersonic flows: I. AUSMPW+ scheme[J]. Journal of Computational Physics, 2001, 174(1): 38-80.

[38]Wang D F, Deng X G, Wang G X, et al. Developing a hybrid flux function suitable for hypersonic flow simulation with high-order methods[J].International Journal for Numerical Methods in Fluids, 2016, 81(5): 309-327.

Numerical simulation of hypersonic double cone flows with high-order methods

WANG Dongfang1, DENG Xiaogang1, WANG Guangxue1,2, LIU Huayong3

(1. College of Aerospace Science and Engineering, National University of Defense Technology, Changsha 410073, China；2. School of Physics, Sun Yat-sen University, Guangzhou 510006, China; 3. China Aerodynamics Research and Development Center, Mianyang 621000, China)

Hypersonic low enthalpy laminar flows of double cone with 25°/55° geometry were simulated by using high-order WCNS (weighted compact nonlinear schemes), and their capabilities to accurately predict laminar shock wave/boundary layer interaction were examined. The simulations were performed through adopting the second order MUSCL, the third-order and the fifth-order WCNS as spatial discretization schemes, employing the second-order dual time-stepping approach for time integration and using different flux functions, such as hybrid Roe-Rusnov, AUSMPW+ and Van Leer, for comparison. The effects of high-order methods on time and grid convergence, as well as the dissipation characteristics of flux functions, were analyzed. The numerical simulation results indicate that the high-order methods can obtain well-resolved results on coarse grid and eliminate the sensitivity of flux functions. However, the high-order methods need longer computational time to reach convergence. The computed results show good agreement with the experimental data, and the computational accuracy may be characterized as reasonable for most engineering purposes.

verification and validation; shock wave/boundary layer interaction; high-order methods; hypersonic flow

10.11887/j.cn.201604009http://journal.nudt.edu.cn

2016-04-06

國防科學(xué)技術(shù)大學(xué)科研計劃資助項目(ZDYYJCYJ20140101)

王東方(1987—)，男，河南周口人，博士研究生，E-mail：dongfanghnnj313@163.com；鄧小剛(通信作者)，男，教授，博士，博士生導(dǎo)師，E-mail：xgdeng2000@vip.sina.com

V211.3

A

1001-2486(2016)04-054-10