葛明明，曾 明，趙小宇

(國防科技大學(xué) 航天科學(xué)與工程學(xué)院， 湖南 長沙 410073)

?

加權(quán)緊致非線性格式在熱化學(xué)非平衡流數(shù)值模擬中的應(yīng)用*

葛明明，曾 明，趙小宇

(國防科技大學(xué) 航天科學(xué)與工程學(xué)院， 湖南 長沙 410073)

高精度格式因其在較低的計(jì)算資源消耗下仍能保持足夠的精度而在精細(xì)流動(dòng)結(jié)構(gòu)的模擬中具有獨(dú)特優(yōu)勢(shì)。將顯式5階加權(quán)緊致非線性格式(WCNS-E-5)引入二維高溫非平衡流計(jì)算，進(jìn)行自由流速度5.0～5.7 km/s的高焓風(fēng)洞和高空條件下的圓柱非平衡流場(chǎng)數(shù)值模擬，得到正確的流場(chǎng)參數(shù)分布、壁面壓力和壁面熱流等結(jié)果與實(shí)驗(yàn)值吻合較好。計(jì)算結(jié)果表明，WCNS-E-5格式較2階的MUSCL格式具有更好的熱流網(wǎng)格收斂特性。分析并針對(duì)算例條件初步解決了高階格式在剛性化學(xué)反應(yīng)流模擬中的魯棒性問題，為進(jìn)一步采用高階格式開展流場(chǎng)結(jié)構(gòu)更加復(fù)雜的熱化學(xué)非平衡流數(shù)值研究奠定了良好基礎(chǔ)。

高精度格式；熱化學(xué)非平衡流；數(shù)值模擬

高溫條件下的熱化學(xué)非平衡現(xiàn)象是高超聲速流場(chǎng)的重要特征之一，這方面的數(shù)值研究始于20世紀(jì)六七十年代。經(jīng)過半個(gè)多世紀(jì)的發(fā)展，各方面都得到了較為成熟并被普遍接受的理論方法。解算器方面，國內(nèi)外也發(fā)展出了較為完善、功能全面的解算器，如NASA Langley中心發(fā)展的蘭利氣動(dòng)熱力學(xué)迎風(fēng)松弛算法(Langley Aerothermodynamic Upwind Relaxation Algorithm, LAURA)程序[1]，Gaitonde等發(fā)展的通用氣動(dòng)熱力學(xué)模擬程序(General Aerothermodynamic Simulation Program, GASP)[2-3]等。不過目前高超聲速熱化學(xué)非平衡流模擬所采用的數(shù)值算法以2階格式為主，這對(duì)外形較為簡(jiǎn)單的流場(chǎng)模擬和流動(dòng)機(jī)理研究，尚可以達(dá)到精度需求。而21世紀(jì)發(fā)展的新型高超聲速飛行器，外形相對(duì)復(fù)雜，飛行器流場(chǎng)中激波/邊界層、激波/激波相互作用等流動(dòng)現(xiàn)象普遍存在。此時(shí)不僅存在宏觀流動(dòng)和微觀物理化學(xué)過程的緊密耦合，還存在復(fù)雜的流動(dòng)相互作用，多種效應(yīng)并存且彼此制約，互相影響。激波位置、形狀和強(qiáng)度，邊界層厚度和形狀都直接影響到相互作用的性質(zhì)和強(qiáng)度，影響到流動(dòng)分離的特征和分離區(qū)大小，這些又都進(jìn)一步影響到流場(chǎng)各點(diǎn)的熱化學(xué)狀態(tài)。為深入理解存在復(fù)雜波系和相互作用的非平衡流動(dòng)的內(nèi)在機(jī)理，準(zhǔn)確預(yù)測(cè)這類情況下飛行器氣動(dòng)力、氣動(dòng)熱和氣動(dòng)物理特性，要求非平衡流的數(shù)值模擬具有更高的精度和效率。

高階格式因其兼具高空間分辨率和計(jì)算效率的特點(diǎn)而在精細(xì)流動(dòng)結(jié)構(gòu)的模擬中具有獨(dú)特優(yōu)勢(shì)[4]，將其應(yīng)用于包含復(fù)雜流場(chǎng)結(jié)構(gòu)的熱化學(xué)非平衡流的數(shù)值模擬很有意義。近期國外應(yīng)用高階格式模擬化學(xué)反應(yīng)流的工作主要是采用Shu[5]的加權(quán)本質(zhì)無振蕩(Weighted Essentially Non Oscillatory, WENO)格式。Wang[6]發(fā)展了一個(gè)基于5階WENO格式的高速化學(xué)非平衡流動(dòng)求解器，通過對(duì)圓柱非平衡流的模擬及與2階總變差減小(Total Variation Diminishing, TVD)格式結(jié)果的對(duì)比驗(yàn)證了WENO求解高超聲速化學(xué)非平衡流動(dòng)的可行性。Lani[7]將5階、7階、9階的WENO格式應(yīng)用于非定常三維可壓縮流求解器ADPDIS 3D，進(jìn)行了高超聲速熱化學(xué)非平衡湍流數(shù)值模擬。Prakash[8]采用一種空間5階、時(shí)間3階精度的激波裝配法進(jìn)行化學(xué)非平衡流計(jì)算，流場(chǎng)捕捉到了不同頻率的擾動(dòng)波疊加效果，結(jié)論指出足夠高的精度對(duì)捕捉擾動(dòng)的必要性。

加權(quán)緊致非線性格式[9-11](Weighted Compact Nonlinear Scheme, WCNS)是鄧小剛等在20世紀(jì)90年代提出的具有激波捕捉能力的高精度有限差分格式，已成功應(yīng)用于低速、超聲速和高超聲速(量熱完全氣體)流場(chǎng)的數(shù)值模擬，在大量復(fù)雜流動(dòng)問題中展示出優(yōu)良性能[12-13]。國外學(xué)者Nonomura[14]指出，采用對(duì)稱守恒網(wǎng)格導(dǎo)數(shù)算法(Symmetrical Conservative Metric Method, SCMM)的WCNS能滿足幾何守恒，因而在復(fù)雜網(wǎng)格問題中更具優(yōu)勢(shì)；Matsukawa[15]評(píng)價(jià)：WCNS既具有高的精度，又具有激波捕捉能力，形式簡(jiǎn)便，計(jì)算耗費(fèi)也較少。

將顯式5階加權(quán)緊致非線性格式(WCNS-E-5)引入高超聲速熱化學(xué)非平衡流的數(shù)值求解，探討高階格式在非平衡流中應(yīng)用的數(shù)值困難并采取相應(yīng)措施，編制采用WCNS-E-5的二維熱化學(xué)非平衡流計(jì)算程序。通過圓柱非平衡流場(chǎng)的模擬初步考察程序的有效性，并研究高階和2階格式在網(wǎng)格收斂性上表現(xiàn)的差異。為進(jìn)一步采用高階格式開展具有復(fù)雜流場(chǎng)結(jié)構(gòu)的熱化學(xué)非平衡流動(dòng)數(shù)值研究奠定基礎(chǔ)。

1 熱化學(xué)模型與控制方程

對(duì)高溫空氣考慮7組元(O2, N2, N, O, NO, NO+, e-)?；瘜W(xué)反應(yīng)模型采用Gupta的7組元6反應(yīng)模型[16]，采用Park的雙溫度模型[17]實(shí)現(xiàn)化學(xué)與熱力學(xué)非平衡的耦合?；瘜W(xué)反應(yīng)方程、振動(dòng)松弛模型及化學(xué)反應(yīng)源項(xiàng)和振動(dòng)松弛源項(xiàng)的計(jì)算方法詳見文獻(xiàn)[18-19]。

采用時(shí)間相關(guān)的二維熱化學(xué)非平衡納維斯托克斯(Navier-Stokes, NS)方程，求解流場(chǎng)的定常解。計(jì)算坐標(biāo)系下的無量綱控制方程為：

(1)

(2)

式中，ρi和ρ分別代表組元i和混合氣體的密度，E為單位體積總能，ev為單位質(zhì)量振動(dòng)能，J為坐標(biāo)變換雅可比行列式。

(3)

(4)

(5)

(6)

bx=uτxx+vτxy+qx

(7)

(8)

(9)

式中，ns為組元個(gè)數(shù)，nm為分子組元個(gè)數(shù)。by，qy，qyv的表達(dá)式類似。

2 數(shù)值方法

2.1 時(shí)間格式

控制方程的時(shí)間離散格式采用Shu[20]提出的三階精度具有TVD 性質(zhì)的Runge-Kutta方法。

(10)

式中，U代表原始量，RHS(U)為控制方程右端項(xiàng)，則：

(11)

2.2 對(duì)流通量

根據(jù)WCNS-E-5對(duì)無粘通量項(xiàng)進(jìn)行空間離散，通量導(dǎo)數(shù)的計(jì)算采用六階精度的中心差分格式。

(12)

(13)

根據(jù)左右原始量得到通量：

(14)

本文嘗試了多種通量格式。AUSMPW+通量格式因其較高的計(jì)算效率與在間斷和邊界層內(nèi)良好的分辨率被廣泛運(yùn)用，對(duì)于非平衡流場(chǎng)熱流的模擬表現(xiàn)較好[21]。但在計(jì)算實(shí)踐中發(fā)現(xiàn)，在含有強(qiáng)激波的算例中，WCNS-E-5與AUSMPW+通量格式的組合可能導(dǎo)致激波附近流場(chǎng)出現(xiàn)振蕩，甚至?xí)褂?jì)算發(fā)散。這可能是因?yàn)楦叻直媛矢袷讲蹲郊げ〞r(shí)耗散不夠，引起了數(shù)值計(jì)算不穩(wěn)定。而Shu研究[22]指出，隨著格式精度的提高，通量格式差別對(duì)于計(jì)算結(jié)果的影響逐漸減小。故本文最后選用了耗散較大的Steger-Warming格式。

左右原始量利用5階加權(quán)非線性插值得到。這里以左值為代表給出具體插值方法，右值類似。

(15)

其中

(16)

(17)

其中

(18)

式中，CLk為最優(yōu)權(quán)重系數(shù)，ISk代表模板k的光滑度量函數(shù)，ε為小量，一般取10-6。

一般情況下，這種非線性加權(quán)方法能夠很好地捕捉激波，并在光滑區(qū)域恢復(fù)精度。但在熱化學(xué)非平衡流動(dòng)中，激波附近微量組元的質(zhì)量分?jǐn)?shù)量級(jí)可能在10-10以下，如果ε仍然取10-6，則光滑度量函數(shù)失效，可能導(dǎo)致微量組元密度在激波附近劇烈振蕩。因此需要減小ε的取值。但有學(xué)者[23]針對(duì)WENO研究指出，ε取值過小也會(huì)降低光滑區(qū)域求解精度。因此，本文采用

ε=max[10-26,10-6×min(IS1,IS2,IS3)]

(19)

這樣針對(duì)間斷處不同物理量的量級(jí)設(shè)定ε，避免在間斷區(qū)間影響光滑度量函數(shù)的效果，又不至于在光滑區(qū)域過小而影響求解精度。

2.3 黏性通量

黏性項(xiàng)的差分格式與對(duì)流項(xiàng)的處理一致，采用6階中心差分。半節(jié)點(diǎn)的通量可以用6階中心插值得到的原始量計(jì)算。

(20)

式中

(21)

在非平衡流動(dòng)中，跨過激波后的化學(xué)反應(yīng)使部分組元質(zhì)量分?jǐn)?shù)發(fā)生量級(jí)變化，此時(shí)高階線性插值會(huì)在激波附近產(chǎn)生非物理振蕩，甚至插值出負(fù)值，導(dǎo)致計(jì)算失敗。文獻(xiàn)[13]提出了一種應(yīng)用于完全氣體條件下的非線性插值格式，它利用3個(gè)線性模板，通過構(gòu)造合理的權(quán)重，來消除間斷導(dǎo)致的振蕩，并且在光滑區(qū)域理論上可達(dá)到6階精度。本文則采用了類似2.2節(jié)中給出的對(duì)流項(xiàng)插值的方法，首先得到半節(jié)點(diǎn)的左右值，對(duì)二者取平均得到半節(jié)點(diǎn)的原始量，進(jìn)而計(jì)算出黏性通量。算例計(jì)算實(shí)踐中發(fā)現(xiàn)，對(duì)于非平衡流動(dòng)，本文做法的魯棒性更強(qiáng)。

3 算例與結(jié)果分析

根據(jù)上述數(shù)值方法編制了應(yīng)用WCNS-E-5的二維熱化學(xué)非平衡流計(jì)算程序，開展了高焓條件下的圓柱非平衡流場(chǎng)數(shù)值模擬，通過與實(shí)驗(yàn)和文獻(xiàn)計(jì)算結(jié)果的比較，初步考察程序的有效性。

這里給出三個(gè)代表性算例。算例1為德國航空宇航中心(Deutsches zentrum für Luftund Raumfahrt, DLR)進(jìn)行的高焓風(fēng)洞(High Enthalpy shock tunnel Gōttingen, HEG)圓柱繞流實(shí)驗(yàn)[24]，重點(diǎn)對(duì)比分析實(shí)驗(yàn)測(cè)得的和本文計(jì)算得到的壁面壓力、熱流分布結(jié)果。算例2為中科院力學(xué)所高焓風(fēng)洞自由流條件[18]下的圓柱繞流，主要比較2階(采用MUSCL格式，記作MUSCL2)和高階格式的計(jì)算結(jié)果及二者關(guān)于駐點(diǎn)熱流的網(wǎng)格收斂性。算例3采用與算例2自由流總焓相同的高空自由流條件，與算例1、算例2(均為高焓風(fēng)洞非平衡自由流條件)不同的是，算例3的自由流中僅含N2和O2組元，跨過激波后的化學(xué)反應(yīng)引起原子和離子組元的質(zhì)量分?jǐn)?shù)發(fā)生很大的量級(jí)變化。這對(duì)應(yīng)用高階格式計(jì)算程序的魯棒性提出了更高要求。

3.1 算例條件與計(jì)算網(wǎng)格

三個(gè)算例的來流參數(shù)見表1。

對(duì)于算例1，文獻(xiàn)[25]采用了5組元(不含NO+和e-)的高溫空氣模型，本文將自由流的NO+質(zhì)量分?jǐn)?shù)設(shè)定為10-20。算例3為高空自由流條件，只含有N2和O2組元，本文將其他組元的質(zhì)量分?jǐn)?shù)人為設(shè)定為10-10。

算例1的圓柱半徑為45 mm，算例2和算例3的均為35 mm。計(jì)算域的外邊界在駐點(diǎn)區(qū)離開壁面0.5個(gè)半徑長度。計(jì)算網(wǎng)格量為121×121或121×61，在壁面附近采用指數(shù)拉伸的方法加密，采用了多個(gè)不同的網(wǎng)格拉伸因子生成多套網(wǎng)格(對(duì)應(yīng)壁面處網(wǎng)格雷諾數(shù)在1～1000范圍內(nèi))，以驗(yàn)證計(jì)算結(jié)果的網(wǎng)格收斂性。壁面網(wǎng)格雷諾數(shù)ReΔη采用式(22)計(jì)算：

(22)

式中：Δη1為駐點(diǎn)無量綱法向第一層網(wǎng)格高度；Rn為圓柱半徑；ρw，Vs,w，μw分別為壁面處混合氣體密度、聲速以及黏性系數(shù)。下面對(duì)三個(gè)算例分別給出ReΔη=8，10，8的121×121網(wǎng)格對(duì)應(yīng)的計(jì)算結(jié)果。對(duì)于算例2，還以駐點(diǎn)熱流為代表，詳細(xì)考察了計(jì)算結(jié)果的網(wǎng)格收斂性。

表1 來流參數(shù)

計(jì)算中除5階的WCNS-E-5外，還采用了2階的MUSCL格式，對(duì)比分析二者結(jié)果。通量格式為Steger-Warming格式，熵修正系數(shù)取0.05，MUSCL2格式中采用minmod限制器。三個(gè)算例均設(shè)為等壁溫(300 K)、全催化壁條件。

3.2 算例1

將本文計(jì)算結(jié)果分別與文獻(xiàn)[24]實(shí)驗(yàn)測(cè)量結(jié)果及文獻(xiàn)[25]計(jì)算進(jìn)行對(duì)比，文獻(xiàn)[25]是應(yīng)用GASP程序[2]，采用2階的MUSCL格式，通量格式為AUSMPW+格式，使用minmod限制器，網(wǎng)格量為161×151。從表面壓力分布看，本文結(jié)果、文獻(xiàn)[25]的數(shù)值結(jié)果均與實(shí)驗(yàn)[24]測(cè)量值吻合很好。從表面熱流分布(如圖1所示)看，本文的WCNS-E-5與2階MUSCL格式結(jié)果基本一致，較文獻(xiàn)[25]的數(shù)值結(jié)果更加接近實(shí)驗(yàn)[24]值。另外，本文的5階和2階格式熱流結(jié)果在駐點(diǎn)區(qū)存在一定差異，2階格式得到的熱流分布曲線的光滑性不如5階格式。

圖1 算例1壁面熱流分布Fig.1 Surface heat flux distribution for case 1

(a) p (b) T (c) Tv

(d) CO (e) CNO+圖2 算例2流場(chǎng)參數(shù)分布Fig. 2 Flow field parameters distribution for case 2

3.3 算例2

圖2分別給出了本文采用2階和5階格式得到的流場(chǎng)壓力、溫度和O，NO+組元質(zhì)量分?jǐn)?shù)分布，圖中上半部分為2階格式結(jié)果，下半部分為5階格式結(jié)果。在此算例條件下，兩種格式的結(jié)果基本一致。注意到5階格式得到的流場(chǎng)壓力等值線在圖中標(biāo)出的位置存在抖動(dòng)，這一現(xiàn)象在文獻(xiàn)[26]的量熱完全氣體流場(chǎng)結(jié)果中也存在，是由于激波與網(wǎng)格走向不一致造成的，文獻(xiàn)[26]通過優(yōu)化網(wǎng)格來消除。

2階和5階格式得到的駐點(diǎn)線上組元質(zhì)量分?jǐn)?shù)分布總體吻合，只有CNO+在激波附近和過激波后的峰值存在微小差別，如圖3所示。

圖3 算例2駐點(diǎn)線上組元NO+質(zhì)量分?jǐn)?shù)分布Fig.3 NO+mass fraction distribution along the stagnation line for case 2

下面對(duì)比分析5階WCNS-E-5和2階MUSCL格式在熱流計(jì)算的網(wǎng)格收斂性上的表現(xiàn)。圖4給出了二者得到的駐點(diǎn)熱流隨ReΔη的變化曲線。在5%的誤差范圍內(nèi)，5階格式在ReΔη達(dá)到300左右基本收斂，而2階格式則要求ReΔη低至60。

圖4 算例2駐點(diǎn)熱流網(wǎng)格收斂性曲線Fig.4 Grid convergence curve of the heat flux at stagnation point for case 2

可見同樣精度要求下高階格式對(duì)網(wǎng)格密度要求明顯降低，因此可提高計(jì)算效率。圖5給出了5階格式采用121×61(網(wǎng)格雷諾數(shù)ReΔη=90)，121×61(ReΔη=60)和2階格式采用121×121(ReΔη=30)網(wǎng)格的殘差收斂曲線。5階格式可比2階格式減少網(wǎng)格量，從而減少每個(gè)時(shí)間步的計(jì)算時(shí)間；壁面附近網(wǎng)格也可更疏，這樣能顯著加大每個(gè)時(shí)間步推進(jìn)的Δt，因此得到定常結(jié)果需要的時(shí)間推進(jìn)步數(shù)大大減少。雖然網(wǎng)格量相同時(shí)5階格式的單步計(jì)算時(shí)間比2階格式長(該算例條件下高22%)，但殘差達(dá)到10-12量級(jí)時(shí)5階格式采用ReΔη=90網(wǎng)格的總計(jì)算時(shí)間僅為2階格式采用ReΔη=30網(wǎng)格的22.5%。

圖5 算例2密度殘差收斂曲線Fig.5 Convergence histories of density residuals for case 2

3.4 算例3

(a) CO，CN，CNO

(b) CNO+圖6 算例3駐點(diǎn)線組元質(zhì)量分?jǐn)?shù)分布Fig.6 Species mass fraction distribution along the stagnation line for case 3

與算例2類似，5階和2階格式得到的壓力和溫度分布沒有差異。圖6給出了兩種格式得到的駐點(diǎn)線組元質(zhì)量分?jǐn)?shù)分布。與算例2不同的是，兩種格式得到的微量組元分布在激波附近有細(xì)小差別，峰值也略有不同。這是因?yàn)楦呖諄砹鳁l件下微量組元跨過激波變化更大，而高階格式在同樣網(wǎng)格下對(duì)梯度的分辨精度更高，能夠更好地捕捉梯度。以圖6(b)給出的CNO+分布為例，5階格式曲線上升的起始位置略滯后于2階格式，而峰值則比2階格式稍高。其實(shí)兩種格式結(jié)果差異的這個(gè)特點(diǎn)在其他參數(shù)分布上也有體現(xiàn)，但由于相對(duì)差值很小，表現(xiàn)不明顯。

5階和2階格式在熱流計(jì)算的網(wǎng)格收斂性上表現(xiàn)與算例2類似。采用ReΔη等于320，160，80，40，20，8的121×121網(wǎng)格的計(jì)算結(jié)果表明，5階格式在ReΔη低于160后駐點(diǎn)熱流結(jié)果的變化就在2%以內(nèi)，而2階格式需要ReΔη低至40。

4 結(jié)論

將顯式五階加權(quán)緊致非線性格式WCNS-E-5引入非平衡流的數(shù)值求解，針對(duì)高階格式在化學(xué)反應(yīng)流應(yīng)用中的魯棒性問題采取措施，編制了采用WCNS-E-5格式的二維熱化學(xué)非平衡流計(jì)算程序，開展了高焓條件下的圓柱流場(chǎng)模擬。得到以下結(jié)論：

1)關(guān)于非線性加權(quán)插值方法的改進(jìn)(即針對(duì)間斷處不同物理量的量級(jí)設(shè)定權(quán)系數(shù)計(jì)算中的小量)，確保了微量化學(xué)組元在間斷區(qū)的光滑度量函數(shù)效果；黏性項(xiàng)插值計(jì)算中借鑒對(duì)流項(xiàng)非線性插值方法，解決了激波附近組元質(zhì)量分?jǐn)?shù)量級(jí)變化引起的插值振蕩問題。

2)在高焓風(fēng)洞條件下的圓柱繞流計(jì)算中，5階WCNS-E-5格式結(jié)果優(yōu)于2階的MUSCL格式。HEG風(fēng)洞算例中本文結(jié)果比文獻(xiàn)計(jì)算結(jié)果更接近實(shí)驗(yàn)值，5階格式得到的熱流分布曲線比2階格式光滑。根據(jù)駐點(diǎn)熱流考察的網(wǎng)格收斂性，5階格式也明顯優(yōu)于2階格式。

3)高空自由流條件下，跨過激波后的化學(xué)反應(yīng)使來流中微量組元質(zhì)量分?jǐn)?shù)發(fā)生量級(jí)變化。算例結(jié)果表明，5階格式對(duì)組元梯度的捕捉比2階格式精細(xì)，在熱流計(jì)算的網(wǎng)格收斂性上也有明顯優(yōu)勢(shì)。

4)在同樣精度要求下，采用5階格式可比2階格式減少網(wǎng)格量，并且顯著放寬對(duì)壁面附近網(wǎng)格密度的要求，從而可加大每個(gè)時(shí)間步推進(jìn)的Δt。雖然同樣網(wǎng)格量時(shí)5階格式單步計(jì)算時(shí)間更長，但得到定常結(jié)果需要的總計(jì)算時(shí)間減少。在本文算例2中可降至2階格式的1/4。

5)未來還需在多方面努力，包括較寬范圍速度高度條件下流場(chǎng)激波附近網(wǎng)格的調(diào)試，組元連續(xù)方程和振動(dòng)能方程的特殊處理，進(jìn)一步將高階格式應(yīng)用于存在復(fù)雜干擾現(xiàn)象的熱化學(xué)非平衡流場(chǎng)數(shù)值模擬，發(fā)揮高階格式在捕捉精細(xì)流場(chǎng)結(jié)構(gòu)方面的優(yōu)勢(shì)。

References)

[1] Cheatwood F M, Gnoffo P A. User′s manual for the langley aerothermodynamic upwind relaxation algorithm (LAURA)[R]. NASA TM-4674, 1996.

[2] Knight D, Longo L, Drikakis D, et al. Assessment of CFD capability for prediction of hypersonic shock in interactions[J]. Progress in Aerospace Sciences, 2012, S48/49(2): 8-26.

[3] Olynick D, Tam T. Trajectory based validation of the shuttle heating environment [C]//Proceedings of 31st Thermophysics Conference, AIAA 96-1891, 1996.

[4] 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.

[5] Jiang G S, Shu C W. Efficient implementation of weighted ENO schemes [J]. Journal of Computational Physics, 1995, 126(1): 202-228.

[6] Wang X W, Zhong X L. Nonequilibrium and reactive high-speed flow simulations with a fifth-order WENO scheme[C]//Proceedings of 39th AIAA Fluid Dynamics Conference, AIAA 2009-4041, 2009.

[7] Lani A, Sjogreen B, Yee H C, et al. High-order simulation of hypersonic nonequilibrium flows on overset grids [M]//Center for Turbulence Research Annual Research Briefs, 2010.

[8] Prakash A, Parsons N, Wang X, et al. High-order shock-fitting methods for direct numerical simulation of hypersonic flow with chemical and thermal nonequilibrium [J]. Journal of Computational Physics, 2011, 230(23): 8474-8507.

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

[10] Deng X G, Mao M L. Weighted compact high-order nonlinear schemes for the Euler equations[C]//Proceedings of 13th Computational Fluid Dynamics Conference, AIAA 97-1941, 1997.

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

[12] Deng X G, Wang G X, Tu G H, et al. Applications of high-order weighted compact nonlinear scheme for complex transonic flows[C]//Proceedings of 49th AIAA Aerospace Sciences Meeting including the New Horizons Forum and Aerospace, AIAA 2011-364, 2011.

[13] Liu H Y, Deng X G, Mao M L, et al. High order nonlinear schemes for viscous terms and the application to complex conguration problems[C]//Proceedings of 49th AIAA Aerospace Sciences Meeting including the New Horizons Forum and Aerospace, AIAA 2011-369, 2011.

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

[15] Matsukawa Y. Implicit large eddy simulation of a supersonic flat-plate boundary layer flow by weighted compact nonlinear scheme [J]. International Journal of Computational Fluid Dynamics, 2011, 25(2): 47-57.

[16] Gupta R N, Yos J M, Thompson R A. Review of reaction rates and thermodynamic and transport properties for an 11-species air model for chemical and thermal nonequilibrium calculations to 30000 K[R]. NASA RP-1232, 1990.

[17] Park C. Problems of rate chemistry in the flight regimes of aeroassisted orbital transfer vehicles [C]//Proceedings of 19th Thermophysics Conference, AIAA 1984-1730, 1984.

[18] 曾明. 高焓風(fēng)洞流場(chǎng)測(cè)量的數(shù)值重建和非平衡效應(yīng)的數(shù)值分析[D].北京:中國科學(xué)院力學(xué)研究, 2007.

ZENG Ming. Numerical rebuilding of free-stream measurement and analysis of nonequilibrium effects in high-enthalpy tunnel [D]. Beijing: Institute of Mechanics, Chinese Academy of Sciences, 2007. (in Chinese)

[19] Zeng M, Xu D, Liu J. Novel method to calculate vibrational thermal conduction in hypersonic nonequilibrium flow [J]. Journal of Thermophysics and Heat Transfer, 2016, 30(1): 12-24.

[20] Shu C W, Osher S. Efficient implementation of essentially non-oscillatory shock capturing schemes[J]. Journal of Computational Physics, 1988, 77(2): 493-471.

[21] 柳軍. 熱化學(xué)非平衡流及其輻射現(xiàn)象的實(shí)驗(yàn)和數(shù)值計(jì)算研究[D].長沙:國防科學(xué)技術(shù)大學(xué), 2004.LIU Jun. Experimental and numerical research on thermo-chemical nonequilibrium flow with radiation phenomenon[D]. Changsha: National University of Defense Technology, 2004. (in Chinese)

[22] Shu C W. Essentially non-oscillatory and weighted essentially non-oscillatory schemes for hyperbolic conservation laws [R]. NASA/CR-97-206253, 1997: 19-26

[23] 侯中喜. 超聲速復(fù)雜流場(chǎng)并行數(shù)值分析及高階格式研究[D].長沙: 國防科學(xué)技術(shù)大學(xué), 2000.

HOU Zhongxi. Parallel numerical analysis of hypersonic complex flow and research of high-oder schemes[D]. Changsha: National University of Defense Technology, 2000. (in Chinese)

[24] Karl S, Martinez-Schramm J, Hannemann K. High enthalpy cylinder flow in HEG: a basis for CFD validation [C]//Proceedings of 33rd AIAA Fluid Dynamics Conference and Exhibit,AIAA 2003-4252, 2003.

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

[26] 董義道. WCNS高階格式在典型流動(dòng)狀態(tài)下的應(yīng)用研究[D]. 長沙:國防科學(xué)技術(shù)大學(xué), 2014.

DONG Yidao. Application research of weighted compact high-oder nonlinear schemes in some typical flow conditions[D]. Changsha: National University of Defense Technology, 2014.(in Chinese)

Applications of high-order weighted compact nonlinear scheme for thermo-chemical nonequilibrium flow

GE Mingming， ZENG Ming， ZHAO Xiaoyu

(College of Aerospace Science and Engineering, National University of Defense Technology, Changsha 410073, China)

High-order scheme has advantages in the simulations of complex flow for its good behavior of preserving high spatial resolution with lower cost than the second order methods. WCNS-E-5(weighted compact nonlinear scheme) was applied in the simulations of thermo-chemical nonequilibrium flow around a cylinder with the freestream velocity from 5 km/s to 5.7 km/s. The distribution of correct flow field parameter was acquired, the pressure and heat flux on the surface agreed well with the experimental data. WCNS-E-5 also exhibited better performance on the grid convergence of heat flux than the second order MUSCL scheme. Results show that the present attempt to combine high order scheme with the calculations of stiff chemically reacting flow is basically successful, and it has established a good foundation for the further study.

high-order scheme; thermo-chemical nonequilibrium flow; numerical simulation

10.11887/j.cn.201605026

http://journal.nudt.edu.cn

2016-03-21

國家自然科學(xué)基金資助項(xiàng)目(11572348)；國防科技大學(xué)重大應(yīng)用基礎(chǔ)研究資助項(xiàng)目(ZDYYJCYJ20140101)

葛明明(1991—)，男，江蘇南通人，博士研究生，E-mail：owen2024@sina.cn；曾明(通信作者)，女，副教授，博士，碩士生導(dǎo)師，E-mail：ming_z@163.com

V211.3

A

1001-2486(2016)05-163-07