噴管射流流場及噪音的數值模擬①

胡聲超，鮑福廷，王 中，蔡 強

(1.西北工業大學燃燒、熱結構與內流場重點實驗室，西安 710072;2.西安近代化學研究所，西安 710065)

0 引言

大推力火箭和導彈在發射過程中，由噴管燃氣射流造成的噪音可達160 dB以上，從而對飛行器的載荷、結構及地面設施造成巨大危害[1]。近20～30年，國外在這方面做了較多的基礎性研究工作，Seiner與Ponton對不同溫度空氣的超音速射流進行了試驗研究[2]。Freund利用自己開發的DNS代碼，分別對馬赫數0.9及1.92的高溫空氣射流流場及聲場進行了數值模擬[3－4]。Lupoglazoff N 與 Biancherin A 利用法國航空航天實驗室(ONERA)的MSD代碼對亞音速及超音速的高溫空氣射流噪音進行了詳細的數值模擬[5－6]。由于試驗成本昂貴，數值模擬可大大減少該項費用。

固體火箭發動機燃氣射流是一個極其復雜的流動過程，要建立一個包括所有因素的數學模型非常困難。為了簡化模型，本文忽略氣-固兩相、多組分和化學反應的影響，對噴管超音速高溫空氣射流流場及聲場進行數值模擬，建立一套有效的計算方法，為如何降低噴管射流噪音的工作提供基礎。

本文分兩部分對射流聲場進行計算，首先針對噴管超音速高溫空氣射流流場進行三維非穩態數值模擬，將非穩態流場數據作為輸入參數，進行氣動聲學計算，得到其近場及遠場聲學特性，并將數值計算結果與相關文獻中的試驗結果對比，驗證該方法的可行性。

1 三維非穩態流場數值模擬

1.1 數學模型

由于忽略氣-固兩相、多組分和化學反應的影響，因此僅對三維純氣相的控制方程進行求解:

式中 ρ為密度;→u為速度矢量;Γ為關于φ的交換系數;Sφ為源項。

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

τij為亞格子應力:

LES通過方程(1)濾掉小尺寸渦，然后利用亞格子模型將亞格子應力模型化，本文中選取的亞格子模型為Smagorinskv模型，使整個方程封閉，最后對其進行數值求解。

1.2 物理模型及網格劃分

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

圖1 計算區域及邊界條件Fig.1 Com putational domain and boundary conditions

整個計算區域均采用結構化網格，在噴管出口以及聲源面內部對網格加密，網格總數在400萬左右，網格劃分如圖2所示。

1.3 邊界條件

如圖1所示，入口處采用壓強入口邊界條件，給定入口總壓、總溫，由于外場直徑取值較大，可近似認為射流對其無影響，選用壓強遠場邊界，外場出口采用壓強出口邊界。考慮邊界上壓強波動的反射對計算結果的影響，遠場及外場出口均附加上無反射邊界條件，內部聲源面采用內部邊界條件。各邊界條件具體參數值如表1所示。

圖2 軸向截面網格劃分Fig.2 Axial section of the grid

表1 邊界條件Table 1 Boundary conditions

1.4 流場計算過程及結果

整個計算過程分成2部分，首先計算15 000個時間步(時間步長Δt=2.5×10－5)，使計算過程經歷約2個聲學階段(聲波完全穿過整個計算區域所用的時間);然后，開啟參數統計功能，對整個流場的平均數據進行統計、記錄，與此同時，每2個時間步，記錄聲源面內部非定常流動參數，為遠場聲場的計算提供必要的數據支持。

圖3和圖4是計算完成后流場的結構云圖。由圖3、圖4(a)、(b)可見，噴管的射流流場內部是由一個復雜的膨脹壓縮波系組成。由于噴管出口壓強略大于環境壓強，燃氣在噴管出口首先進行膨脹，噴管出口處出現膨脹波束，膨脹使中心區域壓強降低，外部區域壓強接近于環境壓強，同時燃氣速度向軸外偏移;當外部區域壓強低于環境壓強時，由于外部環境壓強作用，燃氣速度方向又改變為向軸內偏移，燃氣進行壓縮，射流流場中出現壓縮波束。因此，噴管出口流場結構為先出現膨脹波，再出現壓縮波。壓縮結果使中心區域壓強升高，升高到一定程度時，燃氣需再次進行膨脹，形成膨脹波，反復循環，從而出現膨脹波和壓縮波的交替過程。因此，射流流場形成了一個膨脹壓縮波系結構。此外還可看出，當經歷7～8個膨脹壓縮過程后，氣流變得紊亂，出現激烈的湍流現象，考慮到氣動噪音的主要聲源來自湍流，因此選取的聲源面必須包含該區域。

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

圖4 軸向截面射流流場壓強與速度結構云圖Fig.4 Contours of pressure and velocity in axial sectio

圖4(c)、(d)是射流平均場結構云圖。由平均參數可計算出不同時間點聲壓(瞬時壓強與平均壓強的差值)及速度的瞬時波動值，進而得到近場聲場結構。

2 聲場的數值計算

2.1 近場聲場計算

流場計算完成后，將記錄的統計數據與瞬時數據相比較、運算，可得出近場聲場的結構。

2.2 遠場聲場計算

計算遠場聲場時，利用之前記錄的非定常流動參數，采用Ffowcs Williams and Hawkings(FW-H)方程，選取遠場觀測點進行計算(遠場觀測點分布見圖5)。

FW-H方程具體表達式如下:

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

從方程的結構可看出，求解方程主要由求解面積分及體積分構成。其中，面積分對應單極子、雙極子以及積分面內的四極子聲源;體積分對應積分面外的四極子聲源。考慮到體積分過于復雜、計算量大，因此本文去掉了體積分以簡化計算模型。

圖5 遠場觀測點分布示意圖Fig.5 Location of observers in the far field

3 計算結果與試驗結果對比分析

圖6是經過計算后得到的近場聲壓云圖。由圖6可看出，在z=1.5 m左右，即位于7～8個膨脹壓縮過程之后的區域，聲壓絕對值很大。由此可見，聲源中心位于此區域，也就是湍流現象最激烈的區域，聲音由該聲源向四周傳播，其中z軸正方向的聲壓比負方向大。

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

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

圖7 不同角度下遠場的聲強分布Fig.7 The SPL in different angles

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

4 結論

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

圖8 不同角度下觀測點處的聲音頻譜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.