高超聲速激波/邊界層干擾流動(dòng)數(shù)值模擬研究

龔安龍，劉 周，楊云軍，周偉江

(中國航天空氣動(dòng)力技術(shù)研究院，北京 100074)

0 引言

近年來，全世界范圍內(nèi)掀起了高超聲速飛行器研制的熱潮，其中高超聲速流動(dòng)中激波/邊界層干擾現(xiàn)象的準(zhǔn)確預(yù)測(cè)成為亟待解決的空氣動(dòng)力學(xué)難題。高超聲速雙錐分離流動(dòng)是一類典型的激波/邊界層干擾流動(dòng)。其流動(dòng)現(xiàn)象廣泛存在于高超聲速飛行器氣動(dòng)控制面、局部凸起物等與飛行器本體形成的氣流壓縮區(qū)域，可能帶來極其嚴(yán)重的氣動(dòng)力、熱問題。高超聲速雙錐流動(dòng)已經(jīng)得到了廣泛的實(shí)驗(yàn)和數(shù)值模擬研究，并成為CFD代碼確認(rèn)的重要算例。Holden等人［1-4］開展了雙錐高超聲速流動(dòng)的系列實(shí)驗(yàn)研究，獲得了多個(gè)狀態(tài)下準(zhǔn)確的流動(dòng)特性和壓力、熱流分布數(shù)據(jù)。在數(shù)值模擬方面，Candler、Roy、Gnoffo、Nompelis 等人［5-24］做出了很多貢獻(xiàn)，他們研究了高溫真實(shí)氣體效應(yīng)以及稀薄氣體效應(yīng)對(duì)風(fēng)洞來流條件和雙錐壁面流動(dòng)特性的影響，也研究了網(wǎng)格收斂性、數(shù)值格式等方面對(duì)分離流動(dòng)模擬可靠性的影響情況。本文在基于Roe格式的全 Navier-Stokes方程計(jì)算流體力學(xué)(CFD)代碼基礎(chǔ)上引入了一種局部熵修正方法，用于進(jìn)行雙錐高超聲速分離流動(dòng)的模擬驗(yàn)證，并分析了網(wǎng)格和熵修正方法耗散性對(duì)流動(dòng)模擬準(zhǔn)確性的影響特性。

1 雙錐模型及來流條件

雙錐模型來源于Holden等人的高超聲速實(shí)驗(yàn)，這里選用了尖頭的前錐25°和后錐55°模型，如圖1所示。Holden等人通過多個(gè)狀態(tài)的實(shí)驗(yàn)，詳細(xì)研究了高超聲速雙錐流動(dòng)的機(jī)理，建立了其流動(dòng)結(jié)構(gòu)。如圖2所示，雙錐高超聲速流動(dòng)包含了前錐上激波與粘性分離流動(dòng)的相互干擾、兩個(gè)錐面上激波的相互干擾以及由此形成的輸運(yùn)激波侵入第二錐底層流動(dòng)等復(fù)雜的流動(dòng)結(jié)構(gòu)。兩錐之間形成的分離區(qū)大小、輸運(yùn)激波侵入位置及當(dāng)?shù)氐膲毫εc熱流分布等跟來流的條件非常敏感。本文選用了兩個(gè)典型來流狀態(tài)——Run 28和 Run 35——進(jìn)行研究，來源于 holden等人［3-4］的實(shí)驗(yàn)，見表1。

表1 雙錐流動(dòng)的風(fēng)洞實(shí)驗(yàn)來流條件Table 1 Flow conditions of double-cone experiment in wind tunnel

2 網(wǎng)格生成

復(fù)雜的高超聲速雙錐分離流動(dòng)對(duì)于CFD數(shù)值模擬來說極具挑戰(zhàn)性，其中一個(gè)重要原因就是計(jì)算結(jié)果的網(wǎng)格敏感性非常強(qiáng)［7-9］，即網(wǎng)格達(dá)到極密的情況下才能得到滿足網(wǎng)格收斂性條件的流場(chǎng)解，甚至對(duì)于某些數(shù)值方法無法達(dá)到網(wǎng)格收斂性。本文采用了三套由稀到密的網(wǎng)格(二維軸對(duì)稱網(wǎng)格)用以考察網(wǎng)格的收斂性:(1)coarse——稀網(wǎng)格300×100(流向×法向，下同);(2)medium——中等網(wǎng)格 600×200;(3)fine——密網(wǎng)格1200×400。圖3顯示了中等網(wǎng)格規(guī)模的網(wǎng)格分布情況，在前錐和兩錐的壓縮拐角處網(wǎng)格進(jìn)行了加密處理以保證分離區(qū)流動(dòng)的精細(xì)捕捉。

圖1 尖雙錐實(shí)驗(yàn)?zāi)Ｐ统叽鐖DFig.1 Size of the double-cone experimental model

圖2 高超聲速尖雙錐分離流動(dòng)結(jié)構(gòu)示意圖Fig.2 Diagrammatic sketch of the hypersonic sharp double-cone separated flow

圖3 CFD數(shù)值模擬采用的中等規(guī)模網(wǎng)格(600×200)Fig.3 Medium grid adopted in CFD numerical simulations(600×200)

3 數(shù)值計(jì)算方法簡(jiǎn)介

本文采用的CFD數(shù)值計(jì)算方法［25］以可壓縮全Navier-Stokes方程為控制方程;空間離散采用有限體積方法，無粘通量采用Roe格式;時(shí)間推進(jìn)采用隱式SGS方法;定常解獲得的定義為在典型位置(分離點(diǎn)和再附點(diǎn)附近等)的物面監(jiān)控點(diǎn)和分離區(qū)流場(chǎng)參數(shù)不再發(fā)生改變?yōu)橹?。由于所研究的風(fēng)洞實(shí)驗(yàn)條件為層流狀態(tài)，因此CFD數(shù)值模擬中粘性流動(dòng)采用層流處理。

對(duì)于Roe格式，在高超聲速流動(dòng)計(jì)算中需要引入熵修正來抑制紅玉(Carbuncle)現(xiàn)象和保證計(jì)算的穩(wěn)定性。但是通常采用的熵修正方法在空間是各項(xiàng)同性的，并且在全部計(jì)算單元都是同一個(gè)值，這樣就會(huì)在弱梯度方向和粘性主控流動(dòng)區(qū)域產(chǎn)生多余甚至錯(cuò)誤的耗散。多余的耗散將改變分離區(qū)流動(dòng)特性，從而使得粘性摩擦力和熱流計(jì)算結(jié)果出現(xiàn)異常。為此，本文引入了一種局部熵修正方法［26］，用于克服常規(guī)全局熵修正的不足。該局部熵修正方法采用基于壓力的輔助限制器。對(duì)于通量Jacobian矩陣的每一個(gè)特征值λ，需要由最小特征值λmin來限定，即:

其中λmin是當(dāng)?shù)亓鲌?chǎng)參數(shù)的函數(shù)，表示為:

其中V和c分別是當(dāng)?shù)貑卧乃俣群吐曀?，Φ是?dāng)?shù)赜?jì)算單元所有面上的壓力比值的函數(shù)［26］。

該局部熵修正方法在捕捉強(qiáng)激波間斷(單元面上的壓力值存在較大的跳躍)時(shí)能夠有效的抑制紅玉現(xiàn)象發(fā)生，而在其它弱壓力梯度的流動(dòng)區(qū)域(如粘性邊界層流動(dòng)區(qū)域等)不引入熵修正，即不產(chǎn)生多余耗散，很好的克服了全局熵修正方法的弊端。

4 模擬結(jié)果與討論

4.1 網(wǎng)格敏感性分析

對(duì)于Run 28實(shí)驗(yàn)狀態(tài)，采用三套不同密度的網(wǎng)格進(jìn)行了計(jì)算。圖4給出了不同密度網(wǎng)格下計(jì)算獲得的激波邊界線的分布情況比較，可以看到:網(wǎng)格過稀會(huì)導(dǎo)致計(jì)算的分離區(qū)減小;中等網(wǎng)格與密網(wǎng)格計(jì)算結(jié)果已經(jīng)非常吻合。圖5進(jìn)一步比較了壓力系數(shù)和反映熱流的Stanton數(shù)沿流向的分布情況，不難發(fā)現(xiàn):稀網(wǎng)格計(jì)算的流動(dòng)分離區(qū)小，而壓力和熱流峰值高;中等網(wǎng)格已經(jīng)非常接近密網(wǎng)格計(jì)算結(jié)果，即已經(jīng)滿足網(wǎng)格收斂性。

圖4 不同網(wǎng)格密度下計(jì)算獲得的激波位置(Run 28)Fig.4 Results of shock position from different grid levels(Run 28)

圖5 不同網(wǎng)格密度下壓力系數(shù)和Stanton數(shù)計(jì)算結(jié)果(Run 28)Fig.5 Results of pressure coefficient and Stanton number from different grid levels(Run 28)

4.2 熵修正方法耗散性影響分析

為考察本文引入的局部熵修正方法是否有效，這里比較分析了全局和局部熵修正方法在不同修正系數(shù)(反映修正強(qiáng)弱程度)下計(jì)算得到的壓力和熱流參數(shù)，如圖6所示。采用全局熵修正方法，選用較大的修正系數(shù)將帶來很大的耗散，從而導(dǎo)致計(jì)算的分離區(qū)偏小;而如果選用較小的修正系數(shù)，則很容易出現(xiàn)紅玉現(xiàn)象和計(jì)算的不穩(wěn)定。而采用局部熵修正方法，選用不同的修正系數(shù)計(jì)算獲得的分離區(qū)大小幾乎沒有差異，但較大的修正系數(shù)能夠很好的抑制紅玉現(xiàn)象。

4.3 與實(shí)驗(yàn)結(jié)果的比較與分析

為了進(jìn)一步驗(yàn)證本文計(jì)算方法和結(jié)果的可靠性，將局部熵修正方法計(jì)算的結(jié)果與Holden的實(shí)驗(yàn)結(jié)果以及Candler［5-8］采用全局熵修正方法計(jì)算得到的結(jié)果進(jìn)行了比較。兩種實(shí)驗(yàn)條件下的計(jì)算和測(cè)量獲得的壓力系數(shù)和Stanton數(shù)沿流向的分布曲線見圖7和圖8。通過比較發(fā)現(xiàn)，兩種計(jì)算方法獲得的流動(dòng)分離區(qū)相差不大，均與實(shí)驗(yàn)結(jié)果比較吻合;而采用局部熵修正方法獲得的壓力系數(shù)和Stanton數(shù)的分布及其峰值與實(shí)驗(yàn)結(jié)果更加一致。

圖6 不同熵修正方法的壓力系數(shù)和Stanton數(shù)計(jì)算結(jié)果(Run 28)Fig.6 Results of pressure coefficient and Stanton number from different entropy fix methods(Run 28)

圖7 計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果的比較(Run 28)Fig.7 Results comparison between computation and experiment(Run 28)

圖8 計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果的比較(Run 35)Fig.8 Results comparison between computation and experiment(Run 35)

5 結(jié)論

本文采用基于全Navier-Stokes方程的CFD數(shù)值模擬技術(shù)研究了一類典型的高超聲速激波/邊界層干擾流動(dòng)——雙錐分離流動(dòng)。通過研究數(shù)值預(yù)測(cè)結(jié)果的網(wǎng)格敏感性、熵修方法耗散性影響等，獲得如下一些結(jié)論:

(1)激波干擾下的分離流動(dòng)特性及流場(chǎng)物理參數(shù)對(duì)計(jì)算網(wǎng)格比較敏感，過稀的網(wǎng)格將導(dǎo)致分離區(qū)減小，壓力和熱流峰值增大，本文采用的中等網(wǎng)格已經(jīng)達(dá)到網(wǎng)格收斂性要求;

(2)在Roe格式中引入局部熵修正方法，既能夠有效的抑制紅玉現(xiàn)象，也保證在弱壓力梯度的粘性層不引入多余耗散，是一種能夠更加準(zhǔn)確模擬復(fù)雜高超聲速激波/邊界層干擾流動(dòng)的方法。

［1］HOLDEN M S.Shock interaction phenomena in hypersonic flows［R］.AIAA Paper 98-2751，1998.

［2］HOLDEN M S.Experimental studies of laminar separated flows induced by shock wave/boundary layer and shock/shock interaction in hypersonic flows for CFD validationp［R］.AIAA Paper 2000-0930，2000.

［3］HOLDEN M S，WADHAMS T P.Code validation study of laminar shock/boundary layer and shock/shock interactions in hypersonic flows.part A:experimental measurements［R］.AIAA Paper2001-1031，2001.

［4］HOLDEN M S.Experimental studies in LENS shock and expansion tunnel to examine real-gas effects in hypervelocity flows［R］.AIAA Paper 2004-0916，2004.

［5］OLEJNICZAK J，CANDLER G V，HORMUNG H G.Computation of double-cone experiments in high enthalpy nitrogen［R］.AIAA Paper 97-2549，1997.

［6］WRIGHT M J，SINHA K，OLEJNICZAK J，et al.Numerical and experimental investigation of double-cone shock interactions［J］.AIAA Journal，2000，38(12):2268-2276.

［7］CANDLER G V，NOMPELIS I，HOLDEN M S.Computational analysis of hypersonic laminar viscous-inviscid interactions［R］.AIAA Paper 2000-0532，2000.

［8］CANDLER G V，NOMPELIS I，DRUGUET M C.Navier-Stokes predictions of hypersonic double-cone and cylinder-flare flow fields［R］.AIAA Paper 2000-1024，2000.

［9］ROY C J，GALLIS M A，BARTEL T J，et al.Navier-Stokes and direct simulation monte carlo predictions for laminar hypersonic separation［R］.AIAA Paper 2001-1024，2001.

［10］GNOFFO P A.CFD Validation studies for hypersonic flow prediction［R］.AIAA Paper 2001-1025，2001.

［11］CANDLER G V，NOMPELIS I.CFD validation for hypersonic flight:real gas effects［R］.AIAA Paper 2002-0434，2002.

［12］CANDLER G V，NOMPELIS I，DRUGUET M C，et al.CFD validation for hypersonic flight:hypersonic double-cone flow simulations［R］.AIAA Paper 2002-0581，2002.

［12］ROY C J，GALLIS M A，BARTEL T J，et al.Navier-Stokes and DSMC simulations for hypersonic laminar shock-shock interaction flows［R］.AIAA Paper 2002-0737，2002.

［13］DRUGUET M C，CANDLER G V，NOMPELIS I.Simulations of viscous hypersonic double-cone flows:influence of numerics［R］.AIAA Paper 2003-3548，2003.

［14］NOMPELIS I，CANDLER G V，HOLDEN M S，et al.Computational investigation of hypersonic viscous/inviscid interactions in high enthalpy flows［R］.AIAA Paper 2003-3642，2003.

［15］WADHAMS T P，HOLDEN M S.Experimental studies and numerical predictions for hypervelocity vehicle design and code validation［R］.AIAA Paper 2003-6905，2003.

［16］NOMPELIS I，CANDLER G V，MACLEAN M，et al.Numerical investigation of high enthalpy chemistry on hypersonic double-cone experiments［R］.AIAA Paper 2005-0584，2005.

［17］DRUGUET M C，CANDLER G V，NOMPELIS I.Effect of numerics on Navier-Stokes computations of hypersonic double-cone flows［J］.AIAA Journal，2005，43(3):616-623.

［18］NOMPELIS I，CANDLER G V，MACLEAN M，etal.Numerical investigation of double-cone flow experiments with high-enthalpy effects［R］.AIAA Paper 2010-1283，2010.

［19］NOMPELIS I，CANDLER G V，HOLDEN M S，et al.Real gas effects on hypersonic shock wave laminar boundary layer interactions［R］.AIAA Paper 2003-0443，2003.

［20］NOMPELIS I，CANDLER G V，HOLDEN M S.Effect of vibrational nonequilibrium on hypersonic double-cone experiments［J］.AIAA Journal，2003，41(11):2162-2169.

［22］NOMPELIS I，CANDLER G V，WADHAMS T P，et al.Numerical simulation of high-enthalpy experiments in the LENS-X expansion tube facility［R］.AIAA Paper 2004-1000，2004.

［23］MACLEAN M，CANDLER G V，HOLDEN M S.Numerical evaluation of flow conditions in the LENS reflected shock-tunnel facilities［R］.AIAA Paper 2005-0903，2005.

［24］MACLEAN M，HOLDEN M S，WADHAMS T P，etal.A computational analysis of thermochemical studies in the LENS facilities［R］.AIAA Paper 2007-0121，2007.

［25］LIU Z.Delayed detached eddy simulation for static and forced oscillating airfoil at high angle of attack［R］.IAC 2013，2013.

［26］GNOFFO P A.Updates to multi-dimensional flux reconstruction for hypersonic simulations on tetrahedral grids［R］.AIAA Paper2010-1271，2010.