超聲速欠膨脹沖擊射流數(shù)值模擬

李恩義，樂貴高，馬大為，張英琦，何 強

（南京理工大學機械工程學院，南京 210094）

0 引言

沖擊射流是一種既有工程應用背景，又有理論研究價值的獨特的流動現(xiàn)象［1］。其中，欠膨脹沖擊射流在航空航天領域中有著重要影響，如推進系統(tǒng)、表面溫度、熱應力、制造、特別是飛行器系統(tǒng)安全。在超聲速欠膨脹沖擊射流流場中，有著復雜的流場結構，會出現(xiàn)間斷激波、反射激波、馬赫盤、滯止泡等。此外，由于射流的不穩(wěn)定性，沖擊區(qū)域的流場內可能產生卷吸作用、質量擴散、化學反應、傳熱和聲學特性［2］?？紤]到實驗的局限性，運用CFD方法分析高速沖擊射流不僅提供了流場結構特征，而且有助于更好地理解流場結構和進行定量分析。

近幾十年來，國內外對超聲速沖擊射流展開了廣泛深入的研究。在國外，S.Mackie和R.Taghavi［3］采用Wind計算程序分析了二維、三維矩形超聲速沖擊射流，對比分析了湍流模型、噴射距離、壓比和噴管總溫對壁面射流速度、卷吸特征以及沖擊壁面的壓力和表面摩擦力的影響。H.Lou和C.Shih［4］采用PIV實驗的方法分析了軸對稱沖擊射流的流場結構和聲學特性，并研究了微射流對流場聲學特性的影響，得出微射流有助于減弱噪聲反饋環(huán)和減小沖擊射流的不穩(wěn)定性。M.D.Limaye［5］采用薄金屬膜技術對單圓形可壓縮沖擊射流展開實驗研究，分析了馬赫數(shù)和噴射距離對沖擊壁面?zhèn)鳠崴俾实挠绊?。在國內，何楓、謝峻石［6］等利用有限體積法對軸對稱可壓縮N-S程進行離散，對超聲速垂直沖擊射流進行了數(shù)值模擬，得到詳細的波系結構，激波形狀和位置與實驗相吻合。劉海等［7］采用大渦模擬方法對高度欠膨脹超聲速沖擊射流流場進行數(shù)值分析，結果顯示內剪切層的大尺度渦結構的形成與馬赫盤的振蕩相關，在內外剪切層的作用下形成了壁面射流區(qū)域內外交錯的渦結構。姚朝輝等［8］采用PIV技術對高速沖擊射流的流場結構和渦結構進行深入研究，探討流場與聲場的相關性，發(fā)現(xiàn)沖擊單音的存在與否及強弱與渦結構的存在與否及強弱大小相對應，且沖擊單音隨壓比、沖擊距離、噴嘴唇厚等參數(shù)變化規(guī)律也與渦結構與這些因素的變化規(guī)律相一致。劉小軍等［9］等基于RNG k-ε湍流模型，建立了超音速燃氣射流垂直沖擊平板和沖擊浸沒平板的計算模型，分析了不同沖擊條件下努賽爾數(shù)分布規(guī)律和溫度分布、論述了超音速射流傳熱的特性及影響傳熱特性的因素。盡管對沖擊射流傳熱問題已有大量的研究成果，但是多是針對亞聲速、不可壓縮射流或是采用SA、k-ε、k-ω等湍流模型，其與實驗值有著一定的偏差［10-11］。

1 計算方法

1.1 控制方程

對于可壓縮粘性流動的控制方程可寫為：

式中，ρ為密度，ui為xi方向上的速度分量，p為熱力學壓強，ij為剪切應力張量，qi為熱通量向量，f為標量變量，ET為總能量密度：

其中，e為單位質量所含內能，黏性系數(shù)μ根據(jù)Sutherland公式確定

其中，Re為雷諾數(shù)，Sc為施密特數(shù)，具體符號意義可參見文獻［12］。

1.2 湍流模型

k-l湍流模型的優(yōu)勢體現(xiàn)在：1）避免了ε中非齊次壁面邊界條件；2）克服了k-ε封閉模型中出現(xiàn)的反向壓力梯度。k-l湍流模型的輸運方程為：

式中，μt為渦流粘度：

動能產生項基于線性（Boussinesq）模型：

σk=1.0，σl=0.67，α=0.011 8，C2=0.42，

Cμ=0.09，k=0.42

一是抓牢組建工作，夯實黨組織基礎。組建黨的組織是搞好非公有制經濟組織黨建工作的基礎，發(fā)展黨員工作又是組建黨的組織的關鍵點和基礎。首先，要始終抓住發(fā)展黨員工作不動搖，把一線生產能手培養(yǎng)成中層骨干，把中層骨干培養(yǎng)成入黨積極分子，把入黨積極分子培養(yǎng)成合格黨員。要健全黨員能進能出機制，優(yōu)化黨員隊伍結構，重視從青年工人、農民工和知識分子中發(fā)展黨員。要創(chuàng)新黨組織設置模式，全面推進非公經濟組織黨組織建設的“全覆蓋”。積極落實十八大報告精神，以服務群眾、做群眾工作為主要任務，加強基層服務型黨組織建設，以黨的基層組織建設帶動其他各類基層組織建設。

詳細說明可參見文獻［13］。

1.3 邊界條件和求解方法

本文采用3種不同的邊界條件：沖擊壁面采用粘性無滑移壁面；收斂噴管（出口Ma=1）采用給定的總溫總壓；外部來流和側邊出口采用無反射邊界條件。且為更好地捕捉流場特征和減少工作量，對能夠體現(xiàn)流場主要特征的結構化多塊網格進行細化加密。

基于有限體積法離散控制方程，對空間項采用具有二階精度的TVD離散格式和HLLC黎曼求解器以及Min-mod限制器，時間項采用隱式特性的點隱式算法。同時采用變Courant數(shù)和多重網格法［14］來達到加速收斂的目的。

2 數(shù)值計算

2.1 亞聲速速傳熱驗證

傳熱系數(shù)定義為沖擊平面溫度與絕熱壁面溫度之差的形式：

式中：q為通過壁面的熱流密度；Tw和Taw分別是沖擊平面溫度與絕熱壁面溫度；h稱為表面?zhèn)鳠嵯涤址Q對流換熱系數(shù)。

為獲得具有一般意義的傳熱特性，射流沖擊的傳熱特性采用無量綱量—努塞爾數(shù)表示：

式中：h為對流換熱系數(shù)，k為導熱系數(shù)，Lref為特征長度。

為驗證k-l湍流模型的有效性，對亞音速氣體射流垂直沖擊平板的模型進行數(shù)值計算并與Baugh［15］和 Yan［16］實驗結果對比。亞音速氣流經過長為50D的長直管充分發(fā)展，從管口噴出，沖擊到等溫壁面。其中，噴口與平板間距為2D，p=105Pa，T=293 K，Re=23 000，Tw=314.9 K，D=0.040 3 m。

圖1 沖擊平板努賽爾數(shù)對比圖

圖1給出了RNG k-ε湍流模擬、k-l湍流模型和實驗測得的沖擊平板徑向上的努賽爾對比圖。從圖中可以得到，文獻［9］中采用的RNG k-ε湍流模擬中最大誤差約為10.7%，而如圖1所示，采用k-l湍流模型的計算結果與實驗值吻合良好，更為準確地描述了努塞爾數(shù)二次峰值現(xiàn)象，驗證了該模型用于傳熱計算分析的有效性。針對努塞爾數(shù)二次峰值現(xiàn)象，Behnia等［17］指出是由于該區(qū)域內湍動能增大造成的。當壁面射流形成一個急劇增加的速度峰值時，該區(qū)域靜溫減小促使臨近壁面?zhèn)鳠?，進而使其湍動能增加。

2.2 超聲速欠膨脹沖擊射流流場結構

本文所采用的計算工況是Alvi［18］經典實驗，實驗條件：采用不含升力板的聲速噴管，NPR即射流總壓與環(huán)境壓力之比（NPR=P0/Pa）為5，P0=506 625 Pa，Pa=101 325 Pa，噴管出口直徑D=25.4 mm，射流總溫等于環(huán)境溫度T0/Ta=293 K，沖擊距離是噴管直徑的3倍即h/D=3。

沖擊射流的流場結構按其流動特性一般可以分為3個區(qū)域：1）主射流區(qū)，2）沖擊區(qū)，3）壁面射流區(qū)。下頁圖2描述流場結構示意圖，圖3是實驗陰影圖。主射流與自由射流相似，包含壓縮波和膨脹波，經過平板激波后速度降低，進入沖擊區(qū)。沖擊區(qū)內的回流稱為滯止泡，由于滯止點附近的梯度變化較大使該區(qū)域流場結構變化較大。在沖擊區(qū)射流改變方向，沿徑向流動稱之為壁面射流。

圖4是計算所得密度云圖，對比圖3可以得出采用k-l湍流模型能夠很好地捕捉到流場的結構，馬赫盤、內部剪切層、外部射流邊界、滯止泡等。圖5為計算所得馬赫數(shù)云圖，出噴管后射流壓強大于外界壓強，產生膨脹波和壓縮波交替出現(xiàn)的流場，且最高馬赫數(shù)可達2.8。圖6為計算與實驗測量值對比圖。其中，Cp為平面壓力系數(shù)無量綱量，Cp=（Ps-P∞）/（P0-P∞）且Ps為沖擊平面靜壓。由圖中可得，計算值與測量值吻合良好，可見該計算方法的合理性。從整體上看，平面壓力較小，對于環(huán)形峰值的出現(xiàn)與回流區(qū)滯止泡的出現(xiàn)相關。

圖2 沖擊射流流場結構圖

圖3 實驗陰影圖

圖4 計算密度云圖

圖5 計算馬赫數(shù)云圖

圖6 計算與測量壓力系數(shù)對比圖

2.3 超聲速欠膨脹沖擊射流傳熱

對于傳熱特性的研究，所采用計算模型參考文獻［19］中數(shù)據(jù)，噴管出口直徑25.4 mm，設計馬赫數(shù)為1。噴管出口射流總壓與外部環(huán)境壓力之比（P0/Pa）從2.0到4.4。射流滯止溫度分別為493 K、580 K、591 K。沖擊距離（h/D）是噴管直徑的3、6、10倍。表1給出詳細說明。文中比熱比γ取1.4，壓比滿足下式：

表1 計算模型參數(shù)

圖7 文獻與計算溫度對比圖

圖7給出了文獻［19］和本文對于h/d=3下3種不同壓比和射流滯止溫度下的對比，其中，橫坐標為徑向距離與噴管出口直徑之比，縱坐標為沖擊平面溫度與射流滯止溫度之比。從圖中可以看出在沖擊平面徑向距離r/D小于4時有一定誤差，在沖擊平面徑向距離r/D大于4時吻合很好，同時也驗證了所采用的數(shù)值方法的有效性。結合流場結構可以得出，射流在沖擊區(qū)內，速度被滯止，溫度和壓強升高。在壁面射流區(qū)，氣體得到加速，隨之溫度和壓強降低。此外，在射流沖擊區(qū)內沖擊平面溫度最大值略大于射流滯止溫度，這是由于數(shù)值耗散造成的。與Case1和Case2不同，Case3出現(xiàn)了環(huán)形峰值，類似于圖6，這是由于Case1和Case2屬于中等程度欠膨脹射流，而Case3屬于高度欠膨脹射流，其流場出現(xiàn)了滯止泡，影響了溫度和壓強的分布。

圖8～圖10分別給出了射流滯止溫度為493 K、591 K和580 K下，沖擊距離h/d等于3、6、10時，沖擊平面溫度與射流滯止溫度之比的對比圖。從圖中可以得出，隨著沖擊距離的增加，傳熱效果減弱，沖擊平面溫度減小。在r/D小于1區(qū)域內，隨著壓比的增加圖中曲線越陡峭，說明溫度變化越大，且該區(qū)域是流場傳熱特征的主要影響區(qū)。值得注意的是，在z/d=6的工況下，在r/D大于2的區(qū)域內溫度下降幅值小于z/d=3和z/d=10的工況，特別是隨著壓比的增加，這一特征更加明顯。這表明沖擊平面溫度并不是一定隨沖擊距離的增大而得到相應比例的減小。因此，合理選擇沖擊距離對沖擊射流工程應用具有重要的實踐指導意義。

圖8 493 K下，不同沖擊距離的溫度比對比圖

圖9 591 K下，不同沖擊距離的溫度比對比圖

圖10 580 K下，不同沖擊距離的溫度比對比圖

3 結論

利用數(shù)值計算方法，采用k-l湍流模型，對超聲速欠膨脹沖擊射流進行研究分析，得出如下結論：

1）通過k-ε湍流模型和k-l湍流模型與實驗值對比，可得k-l湍流模型在亞聲速射流數(shù)值傳熱問題上的求解具有更好的效果。

2）k-l湍流模型能夠很好地捕捉到超聲速欠膨脹沖擊射流的流場結構，如馬赫盤、內部剪切層、外部射流邊界、滯止泡等，且獲得的平面壓力系數(shù)與實驗值吻合良好，驗證了k-l湍流模型在超聲速欠膨脹沖擊射流流場結構數(shù)值模擬中具有合理性和有效性。

3）通過對計算的平面溫度與射流滯止溫度之比與文獻值對比，進一步驗證了k-l湍流模型在超聲速欠膨脹沖擊射流流場傳熱數(shù)值模擬中具有合理性和有效性。

4）當壓比增加形成高度欠膨脹沖擊射流時，沖擊區(qū)內滯止泡的出現(xiàn)會影響沖擊平面上的壓力和溫度分布。

5）當噴管出口射流總壓與外部環(huán)境壓力之比和射流滯止溫度一定時，沖擊平面部分區(qū)域溫度并不會隨著沖擊距離的增加而得到相應程度的減小。因此，合理選擇沖擊距離對沖擊射流工程應用具有重要的實踐指導意義。