超聲速燃?xì)馍淞髁鲌鎏匦缘娜S數(shù)值模擬

張 磊，王 浩，阮文俊，王金龍，張 成

（南京理工大學(xué) 能源與動力工程學(xué)院，南京210094）

超聲速燃?xì)馍淞髟诤娇张c航天、火箭與導(dǎo)彈、渦輪機(jī)、激光切割及化工冶金設(shè)備等許多軍事與民用工程技術(shù)領(lǐng)域中有著廣泛的應(yīng)用［1－3］。固體火箭發(fā)動機(jī)工作時(shí)會產(chǎn)生高溫高壓燃?xì)馍淞?，使單兵火箭發(fā)射過程中存在聲、光、焰和煙等［4］，有時(shí)會對武器系統(tǒng)及人員造成嚴(yán)重的傷害。如何降低這些發(fā)射特征，便于火箭的操作使用，對超音速燃?xì)馍淞髁鲌龅难芯烤哂兄匾饬x［5－7］。Rodionow 等［8］采用時(shí)間推進(jìn)法對發(fā)動機(jī)超聲速射流流場進(jìn)行了數(shù)值模擬，考慮了有限速率化學(xué)反應(yīng)并討論了湍流的作用。于勝春，湯龍生［9］采用流體計(jì)算軟件對某型號導(dǎo)彈發(fā)動機(jī)的噴管和羽流場進(jìn)行了數(shù)值仿真，分析了羽流流場隨飛行高度和速度的變化規(guī)律。李崢等［10］采用NASA熱力學(xué)軟件計(jì)算復(fù)合推進(jìn)劑化學(xué)平衡組分及發(fā)動機(jī)燃燒室溫度，以此為噴管入口邊界條件，模擬了包含二次燃燒及Al2O3顆粒的射流流場，并分析了化學(xué)反應(yīng)和來流速度對射流流場的影響。趙娜、余永剛等［11］采用CFD軟件對小孔高溫燃?xì)馍淞髟诖髿庵械臄U(kuò)展特性進(jìn)行了數(shù)值模擬。大渦模擬能較好地捕捉流場的細(xì)節(jié)問題，對計(jì)算機(jī)內(nèi)存及速度的要求也比較高，但遠(yuǎn)低于DNS對計(jì)算機(jī)資源的要求。本文主要是以單兵火箭發(fā)射為工程背景，采用大渦模擬（LES）對尾噴管超聲速射流進(jìn)行了的三維非穩(wěn)態(tài)數(shù)值計(jì)算，分析不同噴管尺寸對超聲速射流流場的影響規(guī)律。

1 數(shù)學(xué)模型與數(shù)值方法

1.1 大渦模擬理論

大渦模擬（LES）是介于直接數(shù)值模擬（DNS）與Reynolds平均法（RANS）之間的一種湍流數(shù)值模擬方法。其主要思想是：用N－S方程對比網(wǎng)格尺度大的大渦運(yùn)動進(jìn)行直接數(shù)值模擬，而通過建立通用模型來模擬比網(wǎng)格尺度小的小渦運(yùn)動對大尺度渦運(yùn)動的影響，這樣不僅保證了計(jì)算的精度，還減少了大量計(jì)算時(shí)間。隨著這一方法的成熟以及計(jì)算機(jī)能力進(jìn)一步提高，大渦模擬將逐步成為湍流模擬的主要方法。

對N－S方程在波數(shù)空間或物理空間進(jìn)行過濾，得到的LES控制方程為

式中：ρ為流體密度和為濾波后的速度分量，τij為亞格子尺度應(yīng)力（SGS應(yīng)力），式中帶有上劃線的變量為濾波后的場變量。

由于SGS是未知量，為了使式（1）和式（2）可解，采用Smagorinsky－Lilly亞格子模型，假定SGS應(yīng)力為

式中：δij為Kronecker常數(shù)，當(dāng)i≠j時(shí)δij＝0，當(dāng)i＝j(luò)時(shí)δij＝1；μt為亞格子尺度的湍流粘度。

式中：Δ＝ （ΔxΔyΔz）1／3，＝ （）1／2，＝分別為x、y、z軸方向的網(wǎng)格尺寸，Cs為Smagorinsky常數(shù)，取0.1。

1.2 數(shù)值方法

對于超音速粘性流動的數(shù)值計(jì)算是采用有限體積法對控制方程進(jìn)行了離散，為了保證計(jì)算的準(zhǔn)確性和更好地捕捉膨脹壓縮波等流場細(xì)節(jié)，使用具有間斷分辨率高、穩(wěn)定性好的AUSM＋格式求解對流項(xiàng)，利用三階精度的MUSCL格式對無粘通量進(jìn)行求解，而粘性通量采用標(biāo)準(zhǔn)的二階中心差分格式。時(shí)間步長采用LU－SGS隱式推進(jìn)法提高計(jì)算效率。

2 計(jì)算模型及網(wǎng)格劃分

固體火箭發(fā)動機(jī)在工作過程中會產(chǎn)生高溫高壓燃?xì)?，通過噴管加速流動，形成超聲速燃?xì)馍淞鲊娤蛲獠看髿猸h(huán)境。計(jì)算模型采用4種不同尺寸的噴管進(jìn)行模擬，研究了不同噴管尺寸和噴管出口馬赫數(shù)對射流流場的影響。噴管的幾何參數(shù)如表1所示。表中Dt為喉部直徑，De為出口直徑，Ma為出口馬赫數(shù)。

表1 噴管尺寸

2.1 計(jì)算域設(shè)置

圖1為噴管及外流場的計(jì)算區(qū)域，計(jì)算區(qū)域長為噴管出口直徑的50倍，寬為出口直徑的30倍。在本文計(jì)算中，不考慮火箭發(fā)動機(jī)燃燒室內(nèi)的流動，直接在噴管入口給出燃?xì)馍淞鳁l件。燃?xì)饨茷榭蓧嚎s理想氣體，粘性系數(shù)由Sutherland公式確定。

①壓強(qiáng)入口條件。

噴管入口采用壓強(qiáng)入口邊界條件，設(shè)置總壓p0為1MPa，總溫T0為3 000K；

②壓強(qiáng)出口條件。

大氣環(huán)境邊界采用壓強(qiáng)出口邊界，設(shè)置壓強(qiáng)p1為0.1MPa，溫度T1為300K；

③壁面條件。

噴管壁面采用絕熱、無滑移壁面條件，選用標(biāo)準(zhǔn)壁面函數(shù)處理邊界湍流。

圖1 計(jì)算區(qū)域及邊界條件

2.2 網(wǎng)格劃分

整個(gè)流場計(jì)算區(qū)域均為結(jié)構(gòu)化網(wǎng)格，如圖2所示，為了使結(jié)構(gòu)網(wǎng)格保持較好的正交性，將整個(gè)計(jì)算區(qū)域分割成多個(gè)子區(qū)域。采用六面體網(wǎng)格，對噴管喉部及出口附近的網(wǎng)格進(jìn)行加密處理。通過對網(wǎng)格的無關(guān)性檢測發(fā)現(xiàn)，對噴管1來說，當(dāng)網(wǎng)格數(shù)目達(dá)到約120萬時(shí)，噴管射流流場的變化極小。因此，本文對4種計(jì)算區(qū)域選取的網(wǎng)格數(shù)為120～180萬。

圖2 計(jì)算區(qū)域三維網(wǎng)格劃分

3 數(shù)值模擬結(jié)果及討論

圖3、圖4分別給出了4種不同尺寸噴管燃?xì)馍淞髦行妮S線上的壓強(qiáng)和速度vx的分布曲線。由圖3和圖4可以看出，由于噴管尺寸的不同，燃?xì)馍淞鹘?jīng)過了不同強(qiáng)度的膨脹壓縮過程。在噴管出口附近膨脹壓縮波的強(qiáng)度大，因此壓強(qiáng)和速度的波動比較明顯。隨著典型膨脹—壓縮—膨脹過程的循環(huán)，由于粘性效應(yīng)及環(huán)境壓力的作用，膨脹壓縮強(qiáng)度逐漸減弱，在射流軸線上壓力和速度呈振蕩衰減變化。

由圖3、圖4還可發(fā)現(xiàn)，噴管出口馬赫數(shù)越大，射流軸線上的靜壓和速度變化越劇烈，且經(jīng)歷的膨脹壓縮循環(huán)次數(shù)增加，對發(fā)動機(jī)后面的大氣環(huán)境影響區(qū)域變大。對于出口馬赫數(shù)Ma為2的噴管，燃?xì)馍淞鬏S線上壓強(qiáng)和速度經(jīng)歷4～5個(gè)膨脹壓縮波后趨于穩(wěn)定變化，而出口馬赫數(shù)Ma為2.5時(shí)要經(jīng)歷9～10個(gè)膨脹壓縮波后，射流軸線上離噴口較遠(yuǎn)位置的靜壓才逐漸和環(huán)境靜壓保持一致。當(dāng)噴管出口馬赫數(shù)相同時(shí)，噴管尺寸的變化對射流軸線上一定范圍內(nèi)的各參數(shù)影響較小。

圖3 中心軸線上的靜壓分布

圖4 中心軸線上的速度分布

圖5為距離噴管出口200mm處，不同徑向點(diǎn)速度vy的分布曲線圖。

圖5 徑向上的速度分布

通過圖5可以觀察到，射流軸線上的速度最大，距離軸線較近的核心區(qū)內(nèi)速度變化幅度很小。核心區(qū)就是射流的起始段，并隨著馬赫數(shù)的增大而增大。到達(dá)射流剪切層邊界時(shí)，徑向速度的衰減梯度很大，在y＝15mm處后，速度變化比較平緩。離噴管出口距離相同時(shí)，隨噴管尺寸增大，射流核心區(qū)在徑向上增加。這些分析結(jié)果與射流理論相符合。

圖6為4種不同尺寸噴管燃?xì)馇放蛎洺曀偕淞髁鲌龅乃俣确植荚茍D，由于燃?xì)馍淞髁鞒鰢娍跁r(shí)的靜壓高于噴口外周圍環(huán)境的靜壓，射流一出噴口首先進(jìn)行膨脹，形成Prandtl－Meyer流。燃?xì)饫^續(xù)膨脹加速，使中心區(qū)域壓強(qiáng)降低，當(dāng)射流外部區(qū)域壓強(qiáng)低于環(huán)境壓強(qiáng)時(shí)，由于受到外部環(huán)境的壓縮作用，從而形成了上、下相交射流激波，相交于軸線上，射流在剪切層邊界也形成反射的激波。這樣反復(fù)循環(huán)，在射流流場形成了系列的膨脹壓縮波結(jié)構(gòu)，且膨脹壓縮波結(jié)構(gòu)逐漸崩解，這種循環(huán)過程隨著噴管尺寸的增大而增多，與圖3、圖4中軸線上的參數(shù)相對應(yīng)。在超聲速初始條件相同時(shí)，僅改變噴管尺寸，對射流流場的總體結(jié)構(gòu)基本無影響。隨著膨脹壓縮波崩解后，由于Kelvin－Helmholtz不穩(wěn)定效應(yīng)，導(dǎo)致超聲速欠膨脹射流中出現(xiàn)激烈的紊亂現(xiàn)象，在噴管尺寸最大時(shí)表現(xiàn)最為明顯。

圖6 射流流場速度云圖

4 結(jié)論

從數(shù)值模擬結(jié)果可以得到以下結(jié)論：

①超聲速燃?xì)馍淞髁鲌鲋徐o壓和速度在射流軸線上具有大的波動，噴管出口馬赫數(shù)越大，波動越劇烈。出口馬赫數(shù)相同時(shí)，噴管尺寸的變化對波動幅度影響較小。隨著離噴口距離的增大，波動幅度逐漸減小。

②在距噴管出口200mm處，速度沿徑向在射流核心區(qū)基本不變，隨徑向距離的增加而迅速衰減。隨著噴管尺寸的增大，射流核心區(qū)在徑向上增加。

③超聲速欠膨脹燃?xì)馍淞髋c周圍大氣劇烈摻混，形成了典型的膨脹—壓縮—膨脹循環(huán)過程。隨著噴管尺寸的增大，射流經(jīng)歷的膨脹壓縮循環(huán)次數(shù)就越多，但射流流場總體結(jié)構(gòu)基本不變。由于粘性效應(yīng)及環(huán)境壓力的作用，膨脹壓縮波結(jié)構(gòu)沿軸向逐漸崩解。

［1］張福祥.火箭燃?xì)馍淞鲃恿W(xué)［M］.哈爾濱：哈爾濱工程大學(xué)出版社，2004.ZHANG Fu－xiang.Dynamics of the rocket jet flow［M］.Harbin：Harbin Engineering University Press，2004.（in Chinese）

［2］趙承慶，姜毅.氣體射流動力學(xué)［M］.北京：北京理工大學(xué)出版社，1998.ZHAO Cheng－qing，JIANG Yi.Gas jet dynamics［M］.Beijing：Beijing Institute of Technology Press，1998.（in Chinese）

［3］鄒寧.超聲速噴管設(shè)計(jì)及其數(shù)值模擬和實(shí)驗(yàn)研究［D］.南京：南京航空航天大學(xué)，2008.ZOU Ning.Supersonic nozzle design and its investigation with numerical and experimental methods［M］.Nanjing：Nanjing University of Aeronautics and Astronautics，2008.（in Chinese）

［4］武瑞清，阮文俊，李昕，等.單兵火箭燃?xì)庀嬖O(shè)計(jì)的試驗(yàn)研究與數(shù)值模擬［J］.彈箭與制導(dǎo)學(xué)報(bào)，2009，29（3）：153－156.WU Rui－qing，RUAN Wen－jun，LI Xin，et al.The experimental research and numerical simulation of flame damper design of shoulder launched rocket［J］.Journal of Projectiles，Rockets，Missiles and Guidance，2009，29（3）：153－156.（in Chinese）

［5］常見虎，周長省，李軍.高低空環(huán)境下火箭發(fā)動機(jī)射流流場的數(shù)值分析［J］.系統(tǒng)仿真學(xué)報(bào)，2007，19（16）：3 672－3 675.CHANG Jian－h(huán)u，ZHOU Chang－sheng，LI Jun.Numerical analysis around flow field of rocket in high－low altitude［J］.Journal of System Simulation，2007，19（16）：3 672－ 3 675.（in Chinese）

［6］陳二云，馬大為，趙改平，等.燃?xì)庾杂缮淞鞯母呔乳g斷有限元數(shù)值模擬［J］.彈道學(xué)報(bào)，2009，21（1）：79－82.CHEN Er－yun，MA Da－wei，ZHAO Gai－ping，et al.Numerical simulation of gas free－jet using high－precision discontinuous finite element method［J］.Journal of Ballistics，2009，21（1）：79－82.（in Chinese）

［7］孫振華，徐東來，何國強(qiáng).飛行參數(shù)對導(dǎo)彈發(fā)動機(jī)羽流的影響［J］.固體火箭技術(shù)，2005，28（3）：188－191.SUN Zhen－h(huán)ua，XU Dong－lai，HE Guo－qiang.Effect of flight status parameters on missile motor plume［J］.Journal of Solid Rocket Technology，2005，28（3）：188－191.（in Chinese）

［8］RODIONOW A V.New space－marching technique for exhaust plume simulation，AIAA 2000－3390［R］.2000.

［9］于勝春，湯龍生.固體火箭發(fā)動機(jī)噴管及羽流流場的數(shù)值分析［J］.固體火箭技術(shù)，2004，27（2）：95－97.YU Sheng－chun，TANG Long－sheng.Numerical analysis of nozzle and plume flow field of a solid rocket motor［J］.Journal of Solid Rocket Technology，2004，27（2）：95－97.（in Chinese）

［10］李崢，向紅軍，張小英.復(fù)合推進(jìn)劑固體火箭發(fā)動機(jī)噴流流場數(shù)值模擬［J］.固體火箭技術(shù)，2014，37（1）：37－42.LI Zheng，XIANG Hong－jun，ZHANG Xiao－ying.Numerical simulation of composite solid propellant rocket motor exhaust plume［J］.Journal of Solid Rocket Technology，2014，37（1）：37－42.（in Chinese）

［11］趙娜，余永剛，劉東堯，等.小孔流量發(fā)生器噴口流場特性的數(shù)值模擬［J］.彈道學(xué)報(bào)，2010，22（2）：81－85.ZHAO Na，YU Yong－gang，LIU Dong－yao，et al.Numerical simulation of nozzle flowfield characteristic of orifice flow generator［J］.Journal of Ballistics，2010，22（2）：81－85.（in Chinese）