欠膨脹超聲速沖擊射流的大渦模擬

鄭楓弋， 賴煥新

（ 華東理工大學(xué)機(jī)械與動力工程學(xué)院，上海 200237）

超聲速沖擊射流往往伴隨著強(qiáng)烈的傳熱傳質(zhì)現(xiàn)象，可以提供較大的沖擊力，廣泛應(yīng)用于航空航天等諸多領(lǐng)域[1]。在垂直起降飛行器（STOVL）上，可以利用機(jī)身底部的噴管噴射的高速氣體對地面的沖擊來獲得飛行器所需要的向上推力，從而實(shí)現(xiàn)垂直起降和懸停等動作。超聲速沖擊射流的流場結(jié)構(gòu)復(fù)雜，激波、剪切層與邊界層存在強(qiáng)烈的干涉，其湍流結(jié)構(gòu)是相關(guān)應(yīng)用中流體動力學(xué)的關(guān)鍵。

根據(jù)Mitchell 等[2]的可視化研究，可以簡述欠膨脹超聲速沖擊射流的流場結(jié)構(gòu)，射流噴出后在中心線附近形成一道激波，通常被稱為馬赫盤，馬赫盤周圍有斜激波。在馬赫盤后的亞音速流與斜激波后的超音速流之間會形成剪切層。在馬赫盤后損失的總壓大于通過斜激波的流動損失的總壓，這導(dǎo)致了在平板中心的壓力低于外側(cè)，從而形成一個(gè)再循環(huán)區(qū)。在再循環(huán)區(qū)和馬赫盤后的氣流之間形成一個(gè)接觸面。通過馬赫盤的流體進(jìn)入剪切層，遠(yuǎn)離射流的軸線。斜激波與射流邊界的交點(diǎn)處出現(xiàn)扇形膨脹波，它經(jīng)剪切層反射，形成環(huán)形的壓縮波。這是欠膨脹超聲速沖擊射流可視化得到的平均流場主要特征。

關(guān)于噪聲方面，超聲速沖擊射流會產(chǎn)生強(qiáng)烈的離散頻率噪聲，這種噪聲被稱為沖擊單音，其幅值通常高于沖擊射流中其他噪聲，而成為沖擊射流噪聲的主體。Marsh[3]于1961年首先報(bào)道了射流沖擊平板產(chǎn)生的離散頻率噪聲。Wanger[4]發(fā)展了沖擊射流的不穩(wěn)定性模型，Neuwerth[5]考慮到噴嘴處產(chǎn)生下行壓力波動，對不穩(wěn)定性模型進(jìn)行了改進(jìn)。Tam等[6]基于中性穩(wěn)定波在噴流內(nèi)部主要向上游運(yùn)動的特點(diǎn)，發(fā)展了描述射流沖擊平板的理論模型。Ho等[7]猜測離散頻率的沖擊單音與流場中的反饋環(huán)有關(guān)，是由平板附近的大尺度擬序結(jié)構(gòu)導(dǎo)致。影響沖擊射流流場結(jié)構(gòu)的主要參數(shù)有噴管總壓與環(huán)境壓力之比（Nozzle Pressure Ratio，NPR）、噴口到?jīng)_擊平板的距離以及沖擊平板的大小。Alvi等[8]在實(shí)驗(yàn)中通過粒子圖像測速技術(shù)（Particle Image Velocimetry ，PIV）得到了流場中的速度分布，并測試了沖擊平板上的壓力分布，發(fā)現(xiàn)壓力脈動的主頻與沖擊單音的頻率相同。Powell[9]指出，沖擊單音是受噴嘴與壁面之間的反饋環(huán)控制。然而，由于流場中存在震蕩，并且流-聲關(guān)系十分復(fù)雜，這種反饋控制噪聲的機(jī)理仍然不清晰。

本文參考Henderson等[10]的實(shí)驗(yàn)工況與劉海等[11-12]的數(shù)值模擬，對欠膨脹沖擊射流進(jìn)行數(shù)值模擬分析，著重研究噴嘴到壁面的反饋機(jī)制，考察湍流場中大尺度擬序結(jié)構(gòu)的產(chǎn)生和發(fā)展演化過程以及沖擊射流中的周期性現(xiàn)象。

1 計(jì)算對象與模型

欠膨脹沖擊射流的大渦模擬（LES）使用經(jīng)過Favre濾波后的Navier-Stokes方程作為控制方程。為使方程中的亞網(wǎng)格應(yīng)力項(xiàng)封閉，本文使用了WALE亞網(wǎng)格模型[13]，其亞網(wǎng)格黏度被?；癁椋?/p>

式中：Ls=min(κd,CwV1/3) ；κ為von Kármán常數(shù)；d為網(wǎng)格點(diǎn)到壁面的距離；Cw=0.325，為WALE常數(shù)；V為控制體的體積。

本文計(jì)算參考Henderson的試驗(yàn)工況，其中收縮噴嘴的出口直徑Dj=25.4 mm，噴口距壁面的高度與噴口直徑之比為h/Dj=2.08，來流總壓與環(huán)境壓強(qiáng)的比值ptotal/p0=4.03 ，其完全膨脹馬赫數(shù)為1.55，噴管外廓錐角為30°，唇厚度為1.27 mm。該工況對應(yīng)的完全膨脹馬赫數(shù)為1.55。計(jì)算域的設(shè)置如圖1所示。計(jì)算網(wǎng)格數(shù)為228萬，在沖擊平板上壁面量綱為一距離y+＜5。對于噴管入口給定總壓邊界條件，對于遠(yuǎn)場出口給定無反射邊界條件，沖擊壁面和噴嘴壁面使用絕熱無滑移邊界條件。

圖1 計(jì)算域設(shè)置Fig. 1 Computational domain

使用Fluent UDF在噴口前Dj處添加周向多模態(tài)渦環(huán)擾動，以觸發(fā)流動轉(zhuǎn)捩[14-15]。

計(jì)算中使用AUSM（Advection Upstream Splitting Method）通量分裂，空間上使用3階MUSCL（Monotone Upstream-Centered Schemes for Conservation Laws）格式，時(shí)間上使用2階精度隱式格式。初場由RANS Realizablek-ε湍流模型計(jì)算得出。大渦模擬的時(shí)間步長?t1=1×10?6s ，可保證全場庫朗數(shù)CFL=u＜1，Δx為網(wǎng)格間距。經(jīng)過時(shí)長為400Dj/Uj的發(fā)展，認(rèn)為已經(jīng)排除了初始流場的影響，開始進(jìn)行變量的采樣統(tǒng)計(jì)。采樣間隔為Δt2=5×10?6s，采樣時(shí)間持續(xù)5×10?3s。

圖2 不同渦環(huán)擾動計(jì)算結(jié)果對比Fig. 2 Comparison of calculation results by different vortex ring perturbations

2 結(jié)果與分析

2.1 時(shí)均流場

為驗(yàn)證網(wǎng)格無關(guān)性，本文分別對網(wǎng)格數(shù)為124萬（Coarse），228萬（Mid），342萬（Fine）的3組網(wǎng)格進(jìn)行雷諾時(shí)均模擬。圖3示出了各組網(wǎng)格中心線上壓力分布，從計(jì)算結(jié)果中可以看出，在網(wǎng)格數(shù)228萬與342萬的算例中，計(jì)算結(jié)果在峰值處的偏差3.8%，可滿足網(wǎng)格無關(guān)性的要求。

圖3 各組網(wǎng)格中心線上壓力分布Fig. 3 Pressure distribution on centerline of each mesh

圖4 中心線上物理量分布與實(shí)驗(yàn)對比Fig. 4 Comparison of the value on centerline with experiment

圖5 沖擊射流激波的比較Fig. 5 Comparison of the shock between the calculation and experiment in impinging jet

圖6示出了沖擊區(qū)域時(shí)均速度場計(jì)算與Henderson實(shí)驗(yàn)的比較，圖6(a)是本文計(jì)算的速度云圖，可以看到在噴口與激波之間的膨脹區(qū)域中，流速不斷增大，并且在激波與壁面間存在低速的再循環(huán)區(qū)。這種趨勢與圖6(b)中Henderson實(shí)驗(yàn)得到的時(shí)均速度場基本一致。說明本文的模擬結(jié)果可以有效地模擬超聲速沖擊射流的流場。

圖6 沖擊射流時(shí)均速度計(jì)算與實(shí)驗(yàn)的比較Fig. 6 Comparison of average velocity between the calculation and experiment in impinging jet

2.2 瞬態(tài)結(jié)果與反饋機(jī)制

圖7示出了軸截面瞬時(shí)合速度分布，各幅圖是等時(shí)間間隔的，展示了1個(gè)完整的流場震蕩周期(T=2×10-4s)內(nèi)流場的變化，可以看到流場中激波的位置與強(qiáng)弱存在周期性。圖7(a)和7(b)反映了波系沿軸向向下游運(yùn)動的過程。在馬赫盤與斜激波向下游運(yùn)動的過程中，環(huán)形激波隨之生成并沿軸向向下游發(fā)展，進(jìn)入壁面射流區(qū)。在激波下行過程中，噴口噴出的欠膨脹氣體有了更充分的空間進(jìn)行膨脹，可以達(dá)到更高的速度，這導(dǎo)致了穿過馬赫盤與斜激波的流體速度差增大，使得剪切作用增強(qiáng)，再循環(huán)區(qū)的影響范圍增大。當(dāng)再循環(huán)區(qū)向下游運(yùn)動時(shí)，剪切層和斜激波后的超聲速流動也發(fā)生畸變，斜激波后的環(huán)形激波隨再循環(huán)區(qū)向下游運(yùn)動而消失，其環(huán)形激波生成和消失的周期與馬赫盤的震蕩周期相同。圖7(c)和圖7(d)示出了馬赫盤向上游運(yùn)動的過程中，斜激波隨馬赫盤向上游擺動，在斜激波后方的流體速度差減小，再循環(huán)區(qū)也變小。在外剪切層與壁面射流的交界處（近壁面處y=0.8Dj）的位置產(chǎn)生渦結(jié)構(gòu)，并有波向上游傳播。

圖7 一個(gè)流動震蕩周期內(nèi)的瞬態(tài)速度分布（Δt/T=0.25）Fig. 7 Transient velocity distribution during a flow oscillation period (Δt/T=0.25)

Q判據(jù)被廣泛用于渦的辨識，其定義式為

圖8 瞬態(tài)流場的Q判據(jù)等值面Fig. 8 Iso-surface of Q-criterion in unsteady flow field

圖9 一個(gè)流動震蕩周期內(nèi)的瞬態(tài)速度脈動分布Fig. 9 Transient velocity fluctuation distribution during a flow oscillation period

速度散度 ?·U為連續(xù)性方程中的一項(xiàng)，可以反映流場中的膨脹與壓縮現(xiàn)象，也可以描述聲波的傳播，圖10示出了一個(gè)流動震蕩周期內(nèi)的瞬態(tài)速度散度分布?？梢钥闯觯淞饕唤?jīng)流出噴口，就發(fā)生高強(qiáng)度的膨脹，當(dāng)膨脹波從噴管唇部傳播到軸線附近時(shí)與對稱側(cè)的膨脹波交匯會帶來更加劇烈的膨脹，在軸線上膨脹程度最高處的靜壓值低于環(huán)境壓力。為了與射流外界的環(huán)境壓力平衡，會形成環(huán)形的斜激波，氣體在通過斜激波時(shí)會發(fā)生壓縮。但由于本工況欠膨脹程度過高，斜激波的角度超過了激波正常反射條件，會出現(xiàn)類似于馬赫反射的情況，在軸線附近形成了正激波馬赫盤，攔截激波與正激波在交點(diǎn)處產(chǎn)生反射激波。在初始時(shí)刻，靠近壁面y=0.8Dj的位置處向上游輻射出一道較強(qiáng)的波，波會在噴管唇部發(fā)生反射，反射波在剪切層附近向下游傳播，與新生成的波發(fā)生干涉并繼續(xù)傳播至沖擊平板，從而形成反饋環(huán)。

圖10 一個(gè)流動震蕩周期內(nèi)的瞬態(tài)速度散度分布Fig. 10 Transient velocity divergence distribution during a flow oscillation period

為精確探究流場震蕩與反饋回路的周期，使用帶有Hanning窗的快速傅里葉變換對流場中的脈動壓力進(jìn)行處理，樣本數(shù)N=1024，采樣間隔與數(shù)值模擬采樣間隔Δt相同。圖11和圖12分別為快速傅里葉變換得到的功率譜密度（PSD）在中心線和唇線上的分布，圖中示出了流場中的脈動壓力在fs=5212 Hz和倍頻2fs=10425 Hz時(shí)存在強(qiáng)烈的離散頻率。這與Henderson實(shí)驗(yàn)中測得的沖擊單音頻率（10253 Hz）相符合。在圖11示出的中心線脈動壓力頻率分布和圖12示出的唇線脈動壓力頻率分布（圖中使用對數(shù)標(biāo)尺）中可以看到，在fs以及其各倍頻處都有明顯的強(qiáng)烈離散頻率，該頻率對應(yīng)的周期為T=1.92×10?4s，與前面觀測得到的周期T=2×10?4s基本相同。由于中心線上存在超聲速流動，在超聲速流動中脈動量較小，因此中心線上的數(shù)據(jù)存在一段低谷，而唇線上各處都有較高的脈動水平。

圖11 中心線脈動壓力頻率分布Fig. 11 Frequency distribution of pressure fluctuation on centerline

圖12 唇線脈動壓力頻率分布Fig. 12 Frequency distribution of pressure fluctuation on lip line

圖13、圖14分別示出了對應(yīng)頻率為fs和2fs的脈動壓力分布，圖中示出了在正激波馬赫盤與反射激波的后方有頻率為fs脈動壓力分布，其位于軸線上x=1.6Dj處。并且，在壁面射流區(qū)中的渦環(huán)結(jié)構(gòu)也以fs為主頻。在正激波馬赫盤與反射激波的前后均有頻率為2fs壓力脈動分布，集中于軸線上x=1.26Dj與x=1.6Dj處，并且存在反饋環(huán)的區(qū)域中，也分布有2fs的壓力脈動信號，這可以直接說明，反饋環(huán)的頻率與沖擊單音的頻率是相關(guān)的。

圖13 5 212 Hz脈動壓力幅值Fig. 13 Amplitude of pressure fluctuation at 5 212 Hz

圖14 10 405 Hz脈動壓力幅值Fig. 14 Amplitude of pressure fluctuation at 10 405 Hz

2.3 本征正交分解分析

本征正交分解（Proper Orthogonal Decomposition，POD），也被稱為Karhunen-Loeve展開或主成分分析方法（Principal Components Analysis），是一種基于矩陣論的數(shù)據(jù)統(tǒng)計(jì)分析方法，廣泛應(yīng)用于數(shù)據(jù)降維，流場分析等。POD的主要目的是為流場提供分析手段且保持其物理意義，同時(shí)著眼于使用動態(tài)系統(tǒng)的概念來預(yù)測流場，可以實(shí)現(xiàn)對復(fù)雜的非線性系統(tǒng)的線性降維處理?？紤]脈動分量u′(x,ti) ，將其分解為時(shí)間系數(shù)an(t)與POD的空間模態(tài) φni(x) ，即

本文使用Snapshot方法進(jìn)行POD分解[17-18]，對仿真計(jì)算得到由M行空間點(diǎn)、N列采樣時(shí)間構(gòu)成的原始數(shù)據(jù)矩陣AM×N，定義相關(guān)矩陣RN×N：

求解相關(guān)矩陣RN×N的特征值問題

得到相關(guān)矩陣RN×N的特征值λi與特征向量 φi。若原始數(shù)據(jù)矩陣AM×N為脈動速度，則λn可以表示第n階模態(tài)對流場中湍動能的貢獻(xiàn)。

定義第n階模態(tài)的能量貢獻(xiàn)率為

圖15 POD能量貢獻(xiàn)率Fig. 15 Energy contribution rate of POD modes

圖16 POD累計(jì)能量貢獻(xiàn)率Fig. 16 Cumulated energy contribution rate of POD modes

圖17示出了脈動速度的1、2、3、4、8階速度模態(tài)（各模態(tài)均以該模態(tài)最大幅值進(jìn)行歸一化處理并記為 φˉ ）。1階模態(tài)是流場中能量占比最高的模態(tài)，可以反映噴流流場中主要的大尺度相干結(jié)構(gòu)，具有較強(qiáng)的對稱性。在1階模態(tài)中，射流邊界、馬赫盤與斜激波后、再循環(huán)區(qū)以及壁面射流均具有較強(qiáng)的相關(guān)性。再循環(huán)區(qū)在斜激波后超聲速流的剪切作用下表現(xiàn)為對稱的環(huán)形渦結(jié)構(gòu)，再循環(huán)區(qū)前有明顯的分界面。超聲速流與壁面射流交界處產(chǎn)生的渦結(jié)構(gòu)沿徑向發(fā)展，并將外部環(huán)境流體卷吸到壁面射流中。在1、2階模態(tài)的壁面射流區(qū)域中捕捉到了沿徑向的模態(tài)對。3階模態(tài)主要反映了在射流邊界，再循環(huán)區(qū)與壁面射流的交界存在較強(qiáng)相關(guān)性。其主要現(xiàn)象為沿徑向射出的壁面射流導(dǎo)致在靠近壁面y=0.8Dj～1.6Dj處產(chǎn)生渦結(jié)構(gòu)，卷吸流體回到再循環(huán)區(qū)與壁面射流的交界處。而4階模態(tài)表現(xiàn)為斜激波后的流動直接匯入壁面射流，這與前文提到的在再循環(huán)區(qū)與壁面射流的交界處的渦結(jié)構(gòu)具有周期性相符合。8階模態(tài)捕捉到了在反饋環(huán)控制區(qū)域存在的相干結(jié)構(gòu)，與圖10示出的反饋環(huán)相對應(yīng)。

圖17 1、2、3、4、8階脈動速度模態(tài)Fig. 17 Velocity fluctuation plots of 1st, 2nd, 3rd, 4th, 8th POD mode

3 結(jié) 論

(1)馬赫盤在軸向震蕩，激波下行過程中，噴出的欠膨脹氣體有更充分的膨脹空間，因此可以達(dá)到更高的速度，從而導(dǎo)致穿過馬赫盤與斜激波的流體速度差增大，使得剪切作用增強(qiáng)，再循環(huán)區(qū)的影響范圍增大。

(2)在靠近壁面y=0.8Dj的位置輻射出一道較強(qiáng)的波動并向上游傳播，并經(jīng)噴管唇部發(fā)生反射后，在剪切層附近向下游傳播，與新生成的波發(fā)生干涉并繼續(xù)傳播至沖擊平板，從而形成反饋環(huán)，它主導(dǎo)了沖擊單音的產(chǎn)生，本文通過對脈動壓力進(jìn)行快速傅里葉變換，證實(shí)反饋環(huán)與沖擊單音具有相同的頻率。

(3)在壁面射流區(qū)中，超聲速流與壁面射流交界處受反饋?zhàn)饔卯a(chǎn)生大尺度環(huán)形渦結(jié)構(gòu)，在其沿徑向發(fā)展的過程中，外部環(huán)境流體不斷被卷吸到壁面射流中，并且大尺度環(huán)形渦結(jié)構(gòu)隨沖擊發(fā)生破碎，生成更小的渦旋結(jié)構(gòu)。