壓縮性對渦環(huán)物理特征及其傳播速度的影響規(guī)律

林海燕， 向 陽， 張 斌， 劉 洪

(上海交通大學(xué) 航空航天學(xué)院， 上海 200240)

渦結(jié)構(gòu)廣泛地存在于自然界及工程應(yīng)用中，由于渦環(huán)是3維流場中最基本的渦系結(jié)構(gòu)，眾多學(xué)者開展了對其特征性能的研究.Shariff等[1]指出研究渦環(huán)可以揭示渦動力學(xué)的一些物理特性及基本原理.針對不可壓縮渦環(huán)，研究者從渦環(huán)動力學(xué)的特征[2-3]、極限生長原理[4-5]、多渦環(huán)相互作用[6]等多個方面開展了大量研究.但是，相比于不可壓縮渦環(huán)，可壓縮渦環(huán)的研究十分有限且大多局限于研究可壓縮渦環(huán)的流動結(jié)構(gòu)特征，如針對渦環(huán)夾止和尾跡射流[7]及由激波氣泡相互作用產(chǎn)生的渦環(huán)的生長演化過程[8]的研究，幾乎均未揭示可壓縮性對渦環(huán)物理特征的影響.對于可壓縮流動，一般而言可以分為亞聲速、跨聲速和超聲速3個區(qū)域，其流場存在不同的壓縮性特點.不同的壓縮性會使可壓縮渦環(huán)產(chǎn)生不同的結(jié)構(gòu)特征.

一般采用激波管產(chǎn)生軸對稱的可壓縮渦環(huán)，而通過調(diào)節(jié)激波管驅(qū)動段和被驅(qū)動段之間的壓力比可以改變?nèi)肷浼げǖ鸟R赫數(shù)(Mashock)，即激波到達激波管出口時的馬赫數(shù)，進而得到具有不同可壓縮性的渦環(huán).Elder等[9]在1952年首次通過紋影實驗觀測到了由激波管產(chǎn)生的可壓縮渦環(huán)向下游傳播的過程.Arakeri等[10]利用粒子圖像測速(PIV)技術(shù)研究當Mashock=1.1，1.2，1.3 時的可壓縮渦環(huán)，發(fā)現(xiàn)渦環(huán)完全形成后的傳播速度是激波到達激波管出口時波后速度的0.7倍.Dora等[11]同樣通過PIV技術(shù)研究Mashock=1.27，1.37時的可壓縮渦環(huán)，獲得了渦環(huán)結(jié)構(gòu)參數(shù)的定量數(shù)據(jù)，并驗證了由實驗測得的渦環(huán)傳播速度與理論計算公式求得的結(jié)果相符.此外，部分學(xué)者還關(guān)注了可壓縮渦環(huán)結(jié)構(gòu)中存在的激波.Baird[12]通過差分干涉測量法對可壓縮渦環(huán)進行研究的結(jié)果表明，當Mashock=1.5時，渦環(huán)的徑向區(qū)域存在嵌入后向激波，但該激波在強烈黏性剪切的作用下終止，并不能傳播至渦核中心.Brouillette 等[13]通過電花影圖和紋影拍攝技術(shù)研究不同條件下激波管產(chǎn)生的可壓縮渦環(huán)，發(fā)現(xiàn)能夠產(chǎn)生嵌入激波的Mashock的最小值為1.34.進而，Brouillette等[14]根據(jù)激波的結(jié)構(gòu)特征對可壓縮渦環(huán)進行了更細致的分類：利用紋影拍攝技術(shù)發(fā)現(xiàn)在臨界驅(qū)動長度條件下，當Mashock<1.43時，將會產(chǎn)生沒有嵌入激波的渦環(huán)；當Mashock≥1.43時，則產(chǎn)生帶有嵌入激波的渦環(huán)；當Mashock=1.6時，嵌入激波的逆壓梯度將產(chǎn)生二次渦，這種二次渦被稱為反向旋轉(zhuǎn)渦環(huán).隨后，又有學(xué)者針對反向旋轉(zhuǎn)渦環(huán)開展了更為細致的研究，如：Kontis 等[15]通過高速紋影拍攝技術(shù)觀察到了當Mashock=1.63時，主渦環(huán)下游存在多重反向旋轉(zhuǎn)渦環(huán)；Thangadurai等[16]通過石蠟油煙霧顯示技術(shù)同樣觀測到了反向旋轉(zhuǎn)渦環(huán)并詳細地研究了其形成機理及演變過程.

上述研究結(jié)果均表明，對于可壓縮渦環(huán)，可壓縮性的強弱對渦環(huán)流動結(jié)構(gòu)的影響有很大的不同.然而，有關(guān)渦環(huán)自身的可壓縮性描述及可壓縮性對渦環(huán)自身結(jié)構(gòu)參數(shù)影響的較為系統(tǒng)的研究尚未見報道.激波結(jié)構(gòu)的存在與否是一個可對渦環(huán)結(jié)構(gòu)的穩(wěn)定性產(chǎn)生重要影響的因素，而渦環(huán)自身的結(jié)構(gòu)參數(shù)(如渦核結(jié)構(gòu)，渦環(huán)半徑，傳播速度)對于研究由渦環(huán)結(jié)構(gòu)引起的夾帶、混合、燃燒等現(xiàn)象均十分重要.本文通過有限體積法對激波管產(chǎn)生的可壓縮渦環(huán)進行數(shù)值計算，以研究渦環(huán)自身的可壓縮性及其在不同可壓縮性條件下對渦環(huán)半徑、渦核半徑、渦環(huán)傳播速度以及渦核內(nèi)渦量分布的影響.

1 研究方法

1.1 計算方法和模型

采用課題組自行研發(fā)的軟件進行數(shù)值計算；采用有限體積法對3維可壓縮Navier-Stokes方程進行空間離散；采用5階WENO(Weighted Essentially Non-Oscillatory)格式重構(gòu)對流項；時間推進采用3階Runge-Kutta法；針對激波管產(chǎn)生渦環(huán)這個物理現(xiàn)象，為了方便驗證，采用與文獻[17]相同的物理模型.激波管計算模型示意圖及各組成部分的結(jié)構(gòu)尺寸如圖1所示.激波管由驅(qū)動段A和被驅(qū)動段B組成，內(nèi)徑為64 mm；C和D為環(huán)境區(qū)域，即渦環(huán)的形成區(qū)域及運動區(qū)域；被驅(qū)動段和環(huán)境區(qū)域的初始壓力為101.325 kPa；整個計算區(qū)域的初始溫度為300 K；改變驅(qū)動段的壓力即可得到不同驅(qū)動段與被驅(qū)動段的壓力比(Rp)，進而得到不同的Mashock.定義入射激波到達激波管出口時時間t=0；激波管出口處x=0.

通過改變渦環(huán)生長演化最主要區(qū)域(0 mm

表1 網(wǎng)格信息

為了驗證計算方法的正確性，在Rp=10及驅(qū)動段長度L=230 mm的條件下進行計算TS1-1，并將計算結(jié)果與文獻[18]中相同條件下算例的計算結(jié)果(TS2)進行比較，渦環(huán)生長過程中不同時刻渦環(huán)截面密度梯度云圖對比結(jié)果如圖3所示.由圖3可見：TS1-1與TS2的渦核中心位置以及嵌入激波位置均較為一致(見圖3(a))；TS1-1與TS2在對稱軸附近產(chǎn)生的滑移線位置也較為一致(圖3(b))；TS1同樣得到了反向旋轉(zhuǎn)渦環(huán)結(jié)構(gòu)(圖3(c)).進而在Rp=7，L=165 mm的條件下進行計算(TS1-2)，并將獲得的計算結(jié)果與文獻[17]的計算結(jié)果(TS3)進行對比.當t=560 μs時，流場中當?shù)亓黧w沿x軸橫向穿過渦環(huán)中心線上的徑向速度(v)的分布規(guī)律如圖4所示，其中TS3-m為TS3的實測值.由圖4可見，TS1-2與TS3的實驗及計算結(jié)果無論是趨勢還是拐點數(shù)據(jù)均較為吻合，從而驗證了所提計算方法的準確性.

圖2 不同網(wǎng)格條件下Γ隨t的變化規(guī)律

圖3 TS1與TS2的渦環(huán)截面密度梯度云圖對比結(jié)果

為了研究不同可壓縮性條件下產(chǎn)生的渦環(huán)，設(shè)計5個具有不同Rp值的算例，每個算例對應(yīng)的Mashock(由激波到達出口時根據(jù)正激波關(guān)系計算獲得)取值如表2所示.

圖4 徑向速度分布規(guī)律對比

表2 5個算例的激波馬赫數(shù)

1.2 流函數(shù)確定渦邊界法

因此，在渦環(huán)的隨體坐標系下，流函數(shù)ψ0為

(1)

式中：u0為渦環(huán)傳播速度，由各時刻渦量最大值點的所在位置確定.為了消除取渦量最大值所在位置時的測量誤差，對時間及位置曲線進行高斯擬合[21]，得到u0的平均值.

ψ0的計算步驟如下：

(1) 定義在對稱軸(r=0)上，ψ0=0；以ψ0=0的點為起點，根據(jù)速度場(u,v)進行積分，從而求得整個計算平面的ψ0值.

(2) 獲得整個區(qū)域ψ0值的曲線分布，可知ψ0=0(不包括對稱軸r=0)形成一個封閉曲線，該曲線即為渦環(huán)邊界.根據(jù)算例5計算得到的流函數(shù)分布如圖5所示.從圖中可以清晰地看到渦環(huán)邊界ψ0=0.鑒于渦環(huán)的對稱性，在ψ0=0的上半部分及r=0所包圍的區(qū)域計算環(huán)量.

圖5 流函數(shù)分布示意圖

1.3 渦核內(nèi)渦量分布的無量綱化

渦核內(nèi)的渦量分布是一個對理論分析(數(shù)學(xué)角度)及渦環(huán)的特征描述(物理角度)非常重要的參數(shù).目前，已有多個用于描述渦核內(nèi)渦量分布的模型被提出.Saffman[22]的研究表明，在真實流體中不可壓縮渦環(huán)的渦核內(nèi)渦量分布滿足高斯分布.

對渦量分布和位置的無量綱計算公式如下：

ωn=ω/ω0

(2)

yn=(y-y0)/rc

(3)

式中：ω為渦量；ω0為峰值渦量；y為ω對應(yīng)的徑向位置；y0為ω0對應(yīng)的徑向位置；rc為渦核半徑.

1.4 可壓縮渦環(huán)傳播速度計算

可壓縮渦環(huán)由于自誘導(dǎo)速度而向下游傳播.Moore[23]提出了該傳播速度的理論計算公式：

(4)

式中：?？赏ㄟ^對由流函數(shù)方法確定的邊界區(qū)域的渦量進行積分得到；rr為渦環(huán)半徑；Mavortex為渦馬赫數(shù)，可用于表征渦環(huán)的可壓縮性，

(5)

c∞為無窮遠處的聲速.

2 結(jié)果與討論

2.1 可壓縮渦環(huán)的形成及其演化過程

激波離開激波管出口后向下游運動的同時發(fā)生衍射，導(dǎo)致激波管出口壓力減小，產(chǎn)生向激波管上游傳播的膨脹波，使得激波管內(nèi)流體速度增大；被加速的流體在激波管出口形成剪切層后卷起形成渦環(huán)，渦環(huán)不斷生長并向下游傳播；渦環(huán)完成生長后，從激波管流出的流體會形成尾跡射流.尾跡射流的長度及持續(xù)時間取決于激波管驅(qū)動段的長度和壓力.

圖6 當Rp=3時 的演化過程云圖

圖7 當Rp=6時 的演化過程云圖

2.2 渦環(huán)可壓縮性的定量表征

可壓縮渦環(huán)根據(jù)結(jié)構(gòu)特征可分為不含嵌入激波的渦環(huán)及含有嵌入激波的渦環(huán).Kontis等[15]提出根據(jù)Mashock是否大于1.43來定量判斷是否存在嵌入激波，但Mashock只能作為判斷渦環(huán)是否會產(chǎn)生的一種參考值，對于表征渦環(huán)可壓縮性的相對強弱并不具有物理意義.Baird[12]指出把含有嵌入激波的渦環(huán)稱為超聲速渦環(huán)，并根據(jù)由渦環(huán)對稱軸上最大軸向速度計算得到的局部馬赫數(shù)(Malocal)定量區(qū)分亞聲速和超聲速渦環(huán).此外，Moore[23]通過理論推導(dǎo)提出Mavortex也可用于表征渦環(huán)的可壓縮性.

傳統(tǒng)的空氣動力學(xué)把Ma<0.8稱為亞聲速，Ma=0.8～1.2稱為跨聲速，Ma>1.2則稱為超聲速.而在可壓縮渦環(huán)的研究中，并沒有系統(tǒng)地研究過渦環(huán)自身定量的可壓縮性.在不同Rp條件下，3種表征可壓縮性的馬赫數(shù)(Mashock，Malocal，Mavortex)的取值變化趨勢如圖9所示.由圖9可見：①Mashock與Rp近似成正比，且都大于1，表明都為可壓縮渦環(huán).② 當Rp≥4時，Malocal均大于1.根據(jù)文獻[12]的定義，Rp≥4產(chǎn)生的皆為超聲速渦環(huán)，但當Rp=4時Malocal只是略大于1，結(jié)合圖8可知，此時產(chǎn)生的渦環(huán)中并不存在嵌入激波，說明產(chǎn)生的不是傳統(tǒng)意義上的超聲速渦環(huán)，而更趨向于跨聲速渦環(huán).③ 當Rp=4時，Mavortex略小于1，表明此時的渦環(huán)處于跨聲速渦環(huán)的范圍；當Rp=5時，結(jié)合圖8可知該條件下的渦環(huán)存在嵌入激波，Malocal也遠大于1，表明此時產(chǎn)生的是超聲速渦環(huán)，但其Mavortex的值反而小于1.導(dǎo)致這種現(xiàn)象的原因可能是渦核受到了激波的影響，軸向被拉長，徑向被壓縮，使得渦核的半徑產(chǎn)生變化，導(dǎo)致Mavortex的數(shù)值受到了影響，但這種影響在Rp=6之后就不再顯著.

通過以上分析，為了使渦環(huán)可壓縮性的定量表征能夠與傳統(tǒng)的流動可壓縮性相對應(yīng)，將渦環(huán)也區(qū)分為亞聲速，跨聲速和超聲速3種特征，通過Malocal和Mavortex可以很明顯地區(qū)分具有弱壓縮性的亞聲速特征渦環(huán)和具有強壓縮性的超聲速特征渦環(huán).在跨聲速區(qū)域，由于從是否存在嵌入激波這一現(xiàn)象上無法區(qū)分可壓縮渦環(huán)的特征，只能通過Malocal和Mavortex的值作近似判斷.由于Malocal的計算較為簡便，故可通過其對渦環(huán)壓縮性的強弱進行預(yù)判；而Mavortex是基于整個渦環(huán)的結(jié)構(gòu)參數(shù)計算獲得的，因此更能夠代表渦環(huán)自身的可壓縮性，同時也更具有實用意義.

圖8 不同Rp下主渦環(huán)ω及 的分布云圖

圖9 3種表征渦環(huán)壓縮性的Ma隨Rp的變化

3 可壓縮性對渦環(huán)物理特征的影響

本節(jié)分析Mavortex與可壓縮渦環(huán)的相關(guān)結(jié)構(gòu)參數(shù).描述渦環(huán)的重要參數(shù)包括ω，rc，rr及u0，下面將對這些參數(shù)進行定量分析.

圖10 不同Mavortex條件下ωn分布及其相應(yīng)的高斯分布對比

3.1 渦核內(nèi)渦量分布

渦核內(nèi)無量綱渦量沿yn的分布(ωn)規(guī)律以及在對應(yīng)條件下的高斯分布情況如圖10所示.其中，yn=0為渦核中心.由圖10可見：隨著Mavortex取值的增加，渦量的聚集程度也有所增加；可壓縮渦環(huán)中渦核內(nèi)渦量的分布不滿足高斯分布，并且隨著Mavortex的增加，渦核內(nèi)渦量更加集中，也愈加偏離高斯分布.

3.2 渦環(huán)半徑和渦核半徑

Norbury[24]根據(jù)渦環(huán)的無量綱半徑對渦環(huán)進行分類.可壓縮渦環(huán)直徑的定義如圖11所示.渦環(huán)直徑取渦量絕對值最大的兩點之間的距離，渦環(huán)半徑則取為該距離的1/2；渦核直徑由穿越渦核中心線(平行于y軸)上的軸向速度分布中最大值和最小值的所在位置確定(見圖12)， 渦核半徑則取為渦核直徑的1/2.

圖11 通過渦環(huán)截面渦量分布定義的渦環(huán)直徑

圖12 通過渦核中心u的分布定義的渦核直徑

渦環(huán)半徑和渦核半徑隨Mavortex的變化趨勢如圖13所示.由圖13可見：隨著Mavortex的增加，渦環(huán)半徑呈逐漸增加的趨勢；渦核半徑在Mavortex較小時呈逐漸增加的趨勢，而在Mavortex>1.2時，即出現(xiàn)相對較強的嵌入激波后，由于激波的拉伸作用，使得渦核徑向半徑略微變??；在Mavortex接近1的區(qū)域，Mavortex對兩個半徑的影響突然加劇，表明可壓縮性開始對渦環(huán)產(chǎn)生更劇烈的影響，也可以認為渦環(huán)進入超聲速特征的區(qū)域.

3.3 渦環(huán)傳播速度

圖14 u0及隨Mavortex的變化情況

4 結(jié)論

本文結(jié)合有限體積法求解3維Navier-Stockes方程，對激波管產(chǎn)生的可壓縮渦環(huán)進行數(shù)值模擬,計算結(jié)果與文獻的相關(guān)實驗結(jié)果及計算結(jié)果相符合；分析渦環(huán)的演化過程及其相關(guān)結(jié)構(gòu)；給出了定量表征渦環(huán)可壓縮性的參數(shù)，并定量分析可壓縮性對渦環(huán)參數(shù)的重要影響；揭示渦環(huán)在可壓縮性作用下自身物理特征的變化情況.研究主要得到以下結(jié)論：

(1) 可壓縮渦環(huán)的可壓縮性可以通過Malocal及Mavortex定量表征.由于Malocal的計算方法較為簡單，所以可通過其對渦環(huán)壓縮性的強弱進行預(yù)判；由于Mavortex更能代表渦環(huán)自身的可壓縮性，所以更具有實用意義.

(2) 激波管產(chǎn)生的可壓縮渦環(huán)能夠表現(xiàn)出亞聲速、跨聲速和超聲速3種特征，根據(jù)Malocal和Mavortex的取值可以很明顯地區(qū)分具有弱壓縮性的亞聲速特征渦環(huán)和具有強壓縮性的超聲速特征渦環(huán).在跨聲速區(qū)域，由于從嵌入激波這一現(xiàn)象上無法區(qū)分可壓縮渦環(huán)的特征，只能通過Malocal和Mavortex的值作近似判斷；而Mashock只能表征產(chǎn)生渦環(huán)的激波的強度，對于可壓縮渦環(huán)的特征判定不具備實際意義.

(3) 隨著可壓縮性的增加，渦環(huán)半徑呈上升趨勢，渦核半徑由于嵌入激波的存在使得其在超聲速特征區(qū)域略微有所減小，而渦環(huán)的傳播速度逐漸增加；在3種可壓縮性特征區(qū)域計算得到的渦環(huán)傳播速度和理論公式計算結(jié)果相符合，表明該理論計算公式可同時適用于3種特征區(qū)域，從而可以更為準確地預(yù)測渦環(huán)的傳播速度；可壓縮渦環(huán)渦核內(nèi)渦量的分布不符合高斯分布，隨著可壓縮性的增加，渦量更加集中.

致謝感謝上海交通大學(xué)高性能計算中心超級計算機Л提供的計算幫助.