胡聲超，鮑福廷，王 中，蔡 強(qiáng)

(1.西北工業(yè)大學(xué)燃燒、熱結(jié)構(gòu)與內(nèi)流場(chǎng)重點(diǎn)實(shí)驗(yàn)室，西安 710072;2.西安近代化學(xué)研究所，西安 710065)

0 引言

大推力火箭和導(dǎo)彈在發(fā)射過(guò)程中，由噴管燃?xì)馍淞髟斐傻脑胍艨蛇_(dá)160 dB以上，從而對(duì)飛行器的載荷、結(jié)構(gòu)及地面設(shè)施造成巨大危害[1]。近20～30年，國(guó)外在這方面做了較多的基礎(chǔ)性研究工作，Seiner與Ponton對(duì)不同溫度空氣的超音速射流進(jìn)行了試驗(yàn)研究[2]。Freund利用自己開(kāi)發(fā)的DNS代碼，分別對(duì)馬赫數(shù)0.9及1.92的高溫空氣射流流場(chǎng)及聲場(chǎng)進(jìn)行了數(shù)值模擬[3－4]。Lupoglazoff N 與 Biancherin A 利用法國(guó)航空航天實(shí)驗(yàn)室(ONERA)的MSD代碼對(duì)亞音速及超音速的高溫空氣射流噪音進(jìn)行了詳細(xì)的數(shù)值模擬[5－6]。由于試驗(yàn)成本昂貴，數(shù)值模擬可大大減少該項(xiàng)費(fèi)用。

固體火箭發(fā)動(dòng)機(jī)燃?xì)馍淞魇且粋€(gè)極其復(fù)雜的流動(dòng)過(guò)程，要建立一個(gè)包括所有因素的數(shù)學(xué)模型非常困難。為了簡(jiǎn)化模型，本文忽略氣-固兩相、多組分和化學(xué)反應(yīng)的影響，對(duì)噴管超音速高溫空氣射流流場(chǎng)及聲場(chǎng)進(jìn)行數(shù)值模擬，建立一套有效的計(jì)算方法，為如何降低噴管射流噪音的工作提供基礎(chǔ)。

本文分兩部分對(duì)射流聲場(chǎng)進(jìn)行計(jì)算，首先針對(duì)噴管超音速高溫空氣射流流場(chǎng)進(jìn)行三維非穩(wěn)態(tài)數(shù)值模擬，將非穩(wěn)態(tài)流場(chǎng)數(shù)據(jù)作為輸入?yún)?shù)，進(jìn)行氣動(dòng)聲學(xué)計(jì)算，得到其近場(chǎng)及遠(yuǎn)場(chǎng)聲學(xué)特性，并將數(shù)值計(jì)算結(jié)果與相關(guān)文獻(xiàn)中的試驗(yàn)結(jié)果對(duì)比，驗(yàn)證該方法的可行性。

1 三維非穩(wěn)態(tài)流場(chǎng)數(shù)值模擬

1.1 數(shù)學(xué)模型

由于忽略氣-固兩相、多組分和化學(xué)反應(yīng)的影響，因此僅對(duì)三維純氣相的控制方程進(jìn)行求解:

式中 ρ為密度;→u為速度矢量;Γ為關(guān)于φ的交換系數(shù);Sφ為源項(xiàng)。

由于聲學(xué)計(jì)算對(duì)精度要求較高，雷諾時(shí)均N-S方程無(wú)法準(zhǔn)確描述聲音的產(chǎn)生及傳播，DNS雖然不采用任何簡(jiǎn)化模型對(duì)流動(dòng)進(jìn)行計(jì)算，精度高，但其對(duì)網(wǎng)格質(zhì)量、數(shù)量(N≈Re9/4)要求較高，特別是在高雷諾數(shù)情況下，給運(yùn)算帶來(lái)很大不便，綜合考慮精度及計(jì)算成本，本文采用大渦模擬(LES)模型，式(1)即是過(guò)濾后的N-S方程:

τij為亞格子應(yīng)力:

LES通過(guò)方程(1)濾掉小尺寸渦，然后利用亞格子模型將亞格子應(yīng)力模型化，本文中選取的亞格子模型為Smagorinskv模型，使整個(gè)方程封閉，最后對(duì)其進(jìn)行數(shù)值求解。

1.2 物理模型及網(wǎng)格劃分

為了方便與試驗(yàn)相對(duì)比，本文采用Seiner＆Ponton噴管模型[2]，由于該噴管詳細(xì)的幾何模型沒(méi)有給出，因此在文中定義了一個(gè)簡(jiǎn)單的錐形噴管，滿足出口尺寸De=91.44 mm及出口馬赫數(shù)Ma=2的設(shè)計(jì)條件。外場(chǎng)軸向長(zhǎng)度取70De，徑向取45De(如圖1所示)。

圖1 計(jì)算區(qū)域及邊界條件Fig.1 Com putational domain and boundary conditions

整個(gè)計(jì)算區(qū)域均采用結(jié)構(gòu)化網(wǎng)格，在噴管出口以及聲源面內(nèi)部對(duì)網(wǎng)格加密，網(wǎng)格總數(shù)在400萬(wàn)左右，網(wǎng)格劃分如圖2所示。

1.3 邊界條件

如圖1所示，入口處采用壓強(qiáng)入口邊界條件，給定入口總壓、總溫，由于外場(chǎng)直徑取值較大，可近似認(rèn)為射流對(duì)其無(wú)影響，選用壓強(qiáng)遠(yuǎn)場(chǎng)邊界，外場(chǎng)出口采用壓強(qiáng)出口邊界?？紤]邊界上壓強(qiáng)波動(dòng)的反射對(duì)計(jì)算結(jié)果的影響，遠(yuǎn)場(chǎng)及外場(chǎng)出口均附加上無(wú)反射邊界條件，內(nèi)部聲源面采用內(nèi)部邊界條件。各邊界條件具體參數(shù)值如表1所示。

圖2 軸向截面網(wǎng)格劃分Fig.2 Axial section of the grid

表1 邊界條件Table 1 Boundary conditions

1.4 流場(chǎng)計(jì)算過(guò)程及結(jié)果

整個(gè)計(jì)算過(guò)程分成2部分，首先計(jì)算15 000個(gè)時(shí)間步(時(shí)間步長(zhǎng)Δt=2.5×10－5)，使計(jì)算過(guò)程經(jīng)歷約2個(gè)聲學(xué)階段(聲波完全穿過(guò)整個(gè)計(jì)算區(qū)域所用的時(shí)間);然后，開(kāi)啟參數(shù)統(tǒng)計(jì)功能，對(duì)整個(gè)流場(chǎng)的平均數(shù)據(jù)進(jìn)行統(tǒng)計(jì)、記錄，與此同時(shí)，每2個(gè)時(shí)間步，記錄聲源面內(nèi)部非定常流動(dòng)參數(shù)，為遠(yuǎn)場(chǎng)聲場(chǎng)的計(jì)算提供必要的數(shù)據(jù)支持。

圖3和圖4是計(jì)算完成后流場(chǎng)的結(jié)構(gòu)云圖。由圖3、圖4(a)、(b)可見(jiàn)，噴管的射流流場(chǎng)內(nèi)部是由一個(gè)復(fù)雜的膨脹壓縮波系組成。由于噴管出口壓強(qiáng)略大于環(huán)境壓強(qiáng)，燃?xì)庠趪姽艹隹谑紫冗M(jìn)行膨脹，噴管出口處出現(xiàn)膨脹波束，膨脹使中心區(qū)域壓強(qiáng)降低，外部區(qū)域壓強(qiáng)接近于環(huán)境壓強(qiáng)，同時(shí)燃?xì)馑俣认蜉S外偏移;當(dāng)外部區(qū)域壓強(qiáng)低于環(huán)境壓強(qiáng)時(shí)，由于外部環(huán)境壓強(qiáng)作用，燃?xì)馑俣确较蛴指淖優(yōu)橄蜉S內(nèi)偏移，燃?xì)膺M(jìn)行壓縮，射流流場(chǎng)中出現(xiàn)壓縮波束。因此，噴管出口流場(chǎng)結(jié)構(gòu)為先出現(xiàn)膨脹波，再出現(xiàn)壓縮波。壓縮結(jié)果使中心區(qū)域壓強(qiáng)升高，升高到一定程度時(shí)，燃?xì)庑柙俅芜M(jìn)行膨脹，形成膨脹波，反復(fù)循環(huán)，從而出現(xiàn)膨脹波和壓縮波的交替過(guò)程。因此，射流流場(chǎng)形成了一個(gè)膨脹壓縮波系結(jié)構(gòu)。此外還可看出，當(dāng)經(jīng)歷7～8個(gè)膨脹壓縮過(guò)程后，氣流變得紊亂，出現(xiàn)激烈的湍流現(xiàn)象，考慮到氣動(dòng)噪音的主要聲源來(lái)自湍流，因此選取的聲源面必須包含該區(qū)域。

圖3 軸向截面瞬時(shí)馬赫數(shù)云圖Fig.3 Contour of instantaneousmach number in axial section

圖4 軸向截面射流流場(chǎng)壓強(qiáng)與速度結(jié)構(gòu)云圖Fig.4 Contours of pressure and velocity in axial sectio

圖4(c)、(d)是射流平均場(chǎng)結(jié)構(gòu)云圖。由平均參數(shù)可計(jì)算出不同時(shí)間點(diǎn)聲壓(瞬時(shí)壓強(qiáng)與平均壓強(qiáng)的差值)及速度的瞬時(shí)波動(dòng)值，進(jìn)而得到近場(chǎng)聲場(chǎng)結(jié)構(gòu)。

2 聲場(chǎng)的數(shù)值計(jì)算

2.1 近場(chǎng)聲場(chǎng)計(jì)算

流場(chǎng)計(jì)算完成后，將記錄的統(tǒng)計(jì)數(shù)據(jù)與瞬時(shí)數(shù)據(jù)相比較、運(yùn)算，可得出近場(chǎng)聲場(chǎng)的結(jié)構(gòu)。

2.2 遠(yuǎn)場(chǎng)聲場(chǎng)計(jì)算

計(jì)算遠(yuǎn)場(chǎng)聲場(chǎng)時(shí)，利用之前記錄的非定常流動(dòng)參數(shù)，采用Ffowcs Williams and Hawkings(FW-H)方程，選取遠(yuǎn)場(chǎng)觀測(cè)點(diǎn)進(jìn)行計(jì)算(遠(yuǎn)場(chǎng)觀測(cè)點(diǎn)分布見(jiàn)圖5)。

FW-H方程具體表達(dá)式如下:

式中ui為i方向流體速度分量;vi為i方向面速度分量;un為積分面法向流體速度;vn為積分面法向面速度;δ(f)為狄拉克(Dirac delta)函數(shù);H(f)為亥維賽(Heaviside)函數(shù);T ij為萊特希爾(Lighthill)應(yīng)力張量:

從方程的結(jié)構(gòu)可看出，求解方程主要由求解面積分及體積分構(gòu)成。其中，面積分對(duì)應(yīng)單極子、雙極子以及積分面內(nèi)的四極子聲源;體積分對(duì)應(yīng)積分面外的四極子聲源。考慮到體積分過(guò)于復(fù)雜、計(jì)算量大，因此本文去掉了體積分以簡(jiǎn)化計(jì)算模型。

圖5 遠(yuǎn)場(chǎng)觀測(cè)點(diǎn)分布示意圖Fig.5 Location of observers in the far field

3 計(jì)算結(jié)果與試驗(yàn)結(jié)果對(duì)比分析

圖6是經(jīng)過(guò)計(jì)算后得到的近場(chǎng)聲壓云圖。由圖6可看出，在z=1.5 m左右，即位于7～8個(gè)膨脹壓縮過(guò)程之后的區(qū)域，聲壓絕對(duì)值很大。由此可見(jiàn)，聲源中心位于此區(qū)域，也就是湍流現(xiàn)象最激烈的區(qū)域，聲音由該聲源向四周傳播，其中z軸正方向的聲壓比負(fù)方向大。

接著通過(guò)FW-H方程對(duì)遠(yuǎn)場(chǎng)觀測(cè)點(diǎn)的聲壓進(jìn)行計(jì)算，得出不同角度下觀測(cè)點(diǎn)處的總聲壓級(jí)(OLSPL)，并與文獻(xiàn)[2]中的試驗(yàn)結(jié)果進(jìn)行對(duì)比(如圖7所示)。從圖7可看出，兩者的整體趨勢(shì)幾乎一致，在40°～50°左右，達(dá)到最高聲壓級(jí)，但數(shù)值計(jì)算得到的聲壓級(jí)最大值比實(shí)際值大。經(jīng)分析，主要有兩方面的原因:首先，由于本文計(jì)算工況的雷諾數(shù)較高(Re=106)，捕捉的頻域范圍為1～10 kHz，因此對(duì)網(wǎng)格的質(zhì)量要求很高，要想減少該因素的影響，可考慮降低雷諾數(shù)或改善網(wǎng)格的質(zhì)量;其次，實(shí)際試驗(yàn)過(guò)程中，各種外界因素及噴管詳細(xì)結(jié)構(gòu)的不同，也會(huì)對(duì)計(jì)算結(jié)果造成影響?？傮w來(lái)說(shuō)，該方法較好地預(yù)示了噴管射流遠(yuǎn)場(chǎng)聲場(chǎng)的分布。

圖6 軸向截面聲壓云圖Fig.6 Contour of sound pressure in axial section

圖7 不同角度下遠(yuǎn)場(chǎng)的聲強(qiáng)分布Fig.7 The SPL in different angles

最后，分別對(duì)20°、50°處觀測(cè)點(diǎn)處的聲壓進(jìn)行計(jì)算對(duì)比(見(jiàn)圖8)可得到，在低頻區(qū)域(0～4 000 Hz)，計(jì)算結(jié)果與試驗(yàn)結(jié)果較為吻合，但在高頻區(qū)域(4 000～10 000 Hz)相差較大。經(jīng)過(guò)分析，主要是由于網(wǎng)格的質(zhì)量造成的。要準(zhǔn)確地捕獲不同頻率下的聲音，必須在該頻率對(duì)應(yīng)的1個(gè)波長(zhǎng)內(nèi)包含一定數(shù)量的網(wǎng)格節(jié)點(diǎn)，頻率越高，波長(zhǎng)越小，相應(yīng)的網(wǎng)格也需要更密。

4 結(jié)論

利用大渦模擬(LES)對(duì)噴管超音速高溫射流聲場(chǎng)進(jìn)行了數(shù)值模擬。結(jié)果表明，該方法得到的計(jì)算結(jié)果與試驗(yàn)結(jié)果較吻合，但由于網(wǎng)格限制，遠(yuǎn)場(chǎng)的高頻聲譜存在一定誤差。今后工作重點(diǎn)是根據(jù)尾噴管噴射噪音分布特點(diǎn)，提出相應(yīng)降噪方案，并用該數(shù)值模擬方法對(duì)其進(jìn)行驗(yàn)證，降低試驗(yàn)成本，為工程實(shí)踐提供參考。

圖8 不同角度下觀測(cè)點(diǎn)處的聲音頻譜Fig.8 Sound spectrum for different angles

[1] Sutherland L C.Progress and problems in rocket noise prediction for ground facilities[R].AIAA 93-4383.

[2] Seiner JM，Ponton M K.The effects of temperature on supersonic jet noise emission[R].DGLR/AIAA 92-02-046.

[3] Freund JB，Lele SK，Moin P.Direct simulation of a Mach 1.92 jet and its sound field[R].AIAA 98-2291.

[4] Freund JB.Noise sources in a low Reynolds number turbulent jet atMach 0.9[J].J.Fluid Mech.，2001，438:277-305.

[5] Lupoglazoff N，Biancherin A.Comprehensive 3D unsteady simulations of subsonic and supersonic hot jet flow field.part I:aerodynamic analysis[R].AIAA-2002-2599.

[6] Lupoglazoff N，Biancherin A.Comprehensive 3D unsteady simulations of subsonic and supersonic hot jet flow field.part II:acoustics analysis[R].AIAA-2002-2600.