泡狀流中空化與激波相互作用的數(shù)值模擬

董洪輝，王寶壽，陳瑋琪，張珂

（中國船舶科學(xué)研究中心，江蘇無錫214082）

泡狀流中空化與激波相互作用的數(shù)值模擬

董洪輝，王寶壽，陳瑋琪，張珂

（中國船舶科學(xué)研究中心，江蘇無錫214082）

文章針對定?？蓧嚎s泡狀流中回轉(zhuǎn)體的空化繞流進行了數(shù)值模擬，在給定來流速度和環(huán)境壓力的條件下，對不同含氣率下的空化與激波的相互作用進行了研究。首先，對含氣率為0的情形進行計算并同試驗數(shù)據(jù)進行對比，證明所采用的計算模型是可信的。其次，改變流場含氣率從0至0.5進行計算，結(jié)果表明，隨著含氣率的增大，流場的可壓縮性隨之增強，體現(xiàn)為數(shù)值模擬得到了回轉(zhuǎn)體繞流流場的3種波系，包括回轉(zhuǎn)體頭部處的脫體激波、分離面處的膨脹波和空泡尾端的斜激波。而空化與激波的相互作用則體現(xiàn)為：頭激波削弱了空化效應(yīng)，使得空化區(qū)域減??；由于空泡外形的影響，使得空泡尾端出現(xiàn)了斜激波；含氣率超過一定值，空泡已經(jīng)無法閉合在物面上，而是閉合在空泡尾端的激波面上，體現(xiàn)為空泡尾端壁面逆壓梯度趨于平緩。計算結(jié)果揭示了泡狀流中空化與激波相互作用的新的物理現(xiàn)象。

泡狀流；激波；空化

0 引言

均勻混合的氣液兩相流具有可壓縮性，Karplus等人[1]的試驗研究表明，當(dāng)偏離單相極點時，氣液混合物的聲速迅速降低，特別是當(dāng)氣相體積分數(shù)為0.5時，氣液兩相流的聲速接近最小，只有幾十米每秒，如圖1為不同壓力下聲速與含氣率的關(guān)系。因此，在等溫條件下，泡狀流中就可以產(chǎn)生激波，此類激波稱為動力激波。

Eddington[3]對氣液兩相流中的超音速現(xiàn)象進行了研究，得到并證實了當(dāng)氣相體積分數(shù)為0.5時泡狀流中的正激波基本關(guān)系。Jhon等人[4]在等溫條件下，忽略相間速度滑移，提出了一種勢流線性分析方法，對在泡狀流中高速移動的擾動點所產(chǎn)生激波的波型進行了研究。Saito等人[2]使用單相流模型通過求解可壓縮N-S方程組，對NACA0015水翼的三維非定常云空化繞流進行了數(shù)值模擬。陳紅兵[5]使用速度入口條件計算了水翼在氣相體積分數(shù)為0.5的來流中的跨聲速流動，在超聲速來流條件下得到了水翼前端的弓形脫體激波，雖然計算中使用了空化模型但是沒有觀察到空化現(xiàn)象。

圖1 聲速與含氣率的關(guān)系[2]Fig.1 Relation between the speed of sound and the void fraction[2]

可見，對泡狀流動的可壓縮性的理論和數(shù)值研究，多集中于氣相體積分數(shù)為0.5的情形，尚沒有對不同含氣率情況下泡狀流中空化現(xiàn)象開展研究，因而也沒有對空化與激波的相互影響進行研究。本文基于CFD并在給定來流速度的條件下，對回轉(zhuǎn)體在不同氣相體積分數(shù)下的泡狀繞流場中空化與激波的相互作用進行了數(shù)值模擬。

1 多相流基本控制方程

本文采用商用軟件Fluent 12中Mixture多相流模型進行數(shù)值模擬，忽略相間速度滑移[4]，其主要控制方程為[6]：

混合物連續(xù)性方程

混合物動量方程

其中：μt為湍動粘度，混合物粘度

混合物能量方程

氣相體積分數(shù)方程

理想氣體狀態(tài)方程

標(biāo)準(zhǔn)k-ε湍流模型

其中：常數(shù)取值為C1ε=1.44，C2ε=1.92，Cμ=0.09，σk=1.0，σε=1.3。

Zwart-Gerber-Belamri空化模型的汽相輸運方程

其中：Re和Rc分別表示蒸汽泡生成和潰滅的質(zhì)量源項：

其中：蒸汽泡半徑RB=10-6m，氣核體積分數(shù)αnuc=5×10-4，汽化系數(shù)Fvap=50，凝結(jié)系數(shù)Fvap=0.001。

2 計算網(wǎng)格、求解器設(shè)置及驗證

2.1 計算網(wǎng)格與求解器設(shè)置

研究對象為回轉(zhuǎn)體，其頭部接近圓錐，后體為半徑R=100 mm、長度L=1 000 mm的圓柱，不計及尾部影響，將圓柱尾端延伸至流場出口邊界。由于流動具有對稱性，因此計算域為二維軸對稱網(wǎng)格，其中入口邊界距回轉(zhuǎn)體頭部50R，圓柱長度L=100R，出口邊界距回轉(zhuǎn)體頭部100R，如圖2。使用ICEM CFD 12.0進行結(jié)構(gòu)化網(wǎng)格劃分，共劃分71 631個四邊形單元，回轉(zhuǎn)體附近網(wǎng)格采用C型網(wǎng)格，如圖2所示。

圖2 計算域與回轉(zhuǎn)體附近網(wǎng)格Fig.2 Computational domain and mesh around the body

陳紅兵[5]的數(shù)值試驗表明，當(dāng)泡狀流中氣相本身為低亞聲速流動或氣相馬赫數(shù)小于0.1時，速度入口條件可用于泡狀流的無窮遠或入口邊界條件。因此，入口和外部邊界均采用速度入口條件，給定來流速度為V∞＝30 m/s；回轉(zhuǎn)體表面采用壁面無滑移無穿透條件；出口采用壓力出口條件；計算域的中軸線設(shè)置為軸對稱邊界條件；流場初始溫度T=300K，工作壓力P=138 540 Pa，空化數(shù)σ=0.3。

采用基于壓力—速度耦合的SIMPLE算法進行求解，計算過程是非定常的，采用Mixture多相流模型，忽略各相間的速度滑移，第一相為水，第二相分別為空氣和水蒸氣，其中空氣為理想氣體。以圓柱直徑為特征長度的雷諾數(shù)Re=3×106，采用標(biāo)準(zhǔn)k-ε湍流模型，近壁處理采用標(biāo)準(zhǔn)壁面函數(shù)。邊界湍流參數(shù)設(shè)置同褚學(xué)森等人[7]所采用的一致，為Turbulent Intensity=0.5%，Turbulent Viscosity Ratio=5。

空化模型采用Zwart-Gerber-Belamri模型，飽和蒸氣壓為3 540 Pa。

2.2 計算模型驗證

來流速度V∞＝30 m/s、工作壓力P=138 540 Pa（即空化數(shù)σ=0.3），對含氣率αg=0的情形進行了計算，并同權(quán)曉波等人[8]的試驗結(jié)果進行了對比，如圖3所示。

由圖3的回轉(zhuǎn)體繞流流場的密度云圖可以清晰地看見空泡外形以及尾部的回射流結(jié)構(gòu)，而壁面壓力系數(shù)曲線與試驗結(jié)果吻合得較好，這表明所采用的計算模型是可信的。

圖3 回轉(zhuǎn)體密度云圖和壁面壓力系數(shù)曲線Fig.3 Axis-symmetrical body’s density color contour and pressure coefficient curve on wall

3 計算結(jié)果及討論

基于以上經(jīng)過驗證的計算模型，本文保持來流速度V∞＝30 m/s、工作壓力P=138 540 Pa，計算了氣液均勻混合介質(zhì)含氣率0～0.5的流動狀況，以期得到含氣率導(dǎo)致的可壓縮性對回轉(zhuǎn)體繞流流場的影響。

3.1 含氣率對繞流流場的影響

圖4至圖6分別給出了各含氣率下的水相體積分數(shù)云圖、絕對壓力云圖和各含氣率下的來流馬赫數(shù)和無量綱空泡長度。由圖4可見，不同含氣率下，水相在流場中的體積分布差異很大。而造成這種差異的原因正是含氣率導(dǎo)致的流場可壓縮性的不同：從圖5可以看到回轉(zhuǎn)體繞流流場的壓力分布，由于流場含氣，氣相在壓力高的區(qū)域受壓縮，在壓力低的區(qū)域膨脹，從而導(dǎo)致回轉(zhuǎn)體頭部駐點附近的高壓區(qū)水相體積分數(shù)最高，而分離面后的低壓區(qū)水相體積分數(shù)最低；含氣率的升高也導(dǎo)致了流場可壓縮性的增強，由圖6可見，來流速度不變的情況下，隨著含氣率的增大，來流馬赫數(shù)逐漸增大，當(dāng)含氣率增至0.3附近，來流馬赫數(shù)開始大于1，故在圖4（e）、（f）和圖5（e）、（f）中可以清晰地看到回轉(zhuǎn)體頭部前端存在脫體激波。

圖4 各含氣率下水相體積分數(shù)云圖Fig.4 Water volume color contours of each gas fraction

圖5 各含氣率下絕對壓力云圖Fig.5 Absolute pressure color contours of each gas fraction

圖4給出的是水相體積分數(shù)云圖，分離面后水相體積分數(shù)低的區(qū)域即是空泡，若將空泡界面定義為ρ=0.9ρm，則圖6給出了不同含氣率下無量綱空泡長度。由圖6可見，隨著含氣率的增大，空泡長度先增大后減小，其原因可以由空泡的氣體成分來解釋：空泡內(nèi)氣體由水蒸汽和空氣兩部分組成，隨著含氣率的增大，由圖5可見回轉(zhuǎn)體分離面后的低壓區(qū)逐漸減小，則空化效應(yīng)被減弱，因而空泡內(nèi)的水蒸汽部分逐漸變少；但是流場中均勻混合的空氣對空泡可起到通氣作用，隨著含氣率的增大其作用愈發(fā)明顯；含氣率較低時，通氣作用對空泡長度的增益占優(yōu)，含氣率較高時，空化效應(yīng)減弱對空泡長度的影響占優(yōu)，兩種因素綜合影響下，導(dǎo)致了空泡長度先增大后減小。

圖6 不同含氣率下來流馬赫數(shù)和空泡無量綱長度Fig.6 Incoming flow Mach number and dimensionless cavitylength of each gas fraction

3.2 空化與激波的相互作用

上一節(jié)中，我們已經(jīng)知道隨著含氣率的升高，流場可壓縮性增強，流場中將出現(xiàn)激波。現(xiàn)在，讓我們來研究當(dāng)流場可壓縮性增強至出現(xiàn)激波時，激波與空化的相互作用。

圖7 各含氣率下Ma數(shù)云圖Fig.7 Mach number color contours of each gas fraction

圖7給出了各含氣率下的Ma數(shù)云圖，從中可以對比流場的可壓縮性（其中為了便于進行對比，云圖標(biāo)尺是相同的）。圖7（a）表明，含氣率為0時，流場基本是不可壓縮的；圖7（a）至（f）表明，隨著含氣率的增大，來流馬赫數(shù)逐漸增大，當(dāng)αg=0.5時，可以清楚地看到回轉(zhuǎn)體頭部出現(xiàn)脫體激波，這說明隨著含氣率的升高，來流流場的可壓縮性逐漸增強；由圖7（b）-（f），可以看到頭部分離面處都出現(xiàn)了膨脹波，而空泡尾部出現(xiàn)了激波，這也從另一個側(cè)面說明氣液摻混的混合介質(zhì)具有了可壓縮性。對比圖7（b）和（f），發(fā)現(xiàn)含氣率為0.1時的Mamax比含氣率為0.5時要大，這是由于空化效應(yīng)和含氣率的相互作用導(dǎo)致的：如圖1所示，氣液混合物的聲速是壓力和含氣率的函數(shù)，壓力越低、含氣率越接近于0.5，聲速越低；當(dāng)含氣率為0.1時，回轉(zhuǎn)體頭部還未出現(xiàn)頭激波，分離面后的空化效應(yīng)比較明顯，該處壓力僅為水的飽和蒸汽壓量級并且局部含氣率接近0.5；當(dāng)含氣率為0.5時，回轉(zhuǎn)體頭部存在頭激波，波后壓力較高；因此，含氣率為0.1時流場中的Mamax比含氣率為0.5時要大。另外，由圖7（b）也驗證了氣液混合物的聲速與含氣率的關(guān)系—在回轉(zhuǎn)體頭部駐點處水相占支配地位，而在回轉(zhuǎn)體分離面后空化區(qū)域蒸氣相占支配地位，因而兩處的Ma數(shù)均較?。欢诳张萁缑娓浇?，各相摻混導(dǎo)致局部含氣率接近0.5，故Ma數(shù)最大。

由前面的分析知，回轉(zhuǎn)體繞流流場可能存在頭部脫體激波、分離面處的膨脹波和空泡尾部的壓縮波，圖8給出了典型情況下回轉(zhuǎn)體空泡繞流流場的波系示意圖，其中空泡界面定義為αl=0.1。下面對回轉(zhuǎn)體泡狀繞流流場的波系產(chǎn)生的原因進行分析。

圖8 回轉(zhuǎn)體空泡繞流流場的波系Fig.8 Wave system of flow field around axis-symmetrical body

圖9為含氣率αg分別為0～0.5時的壓力等值線圖。對比可以發(fā)現(xiàn)，當(dāng)含氣率增大至一定程度時（由圖6，αg≈0.3時，來流Ma≥1），在回轉(zhuǎn)體頭部會出現(xiàn)脫體激波，而產(chǎn)生脫體激波后，波后壓力升高，這也是圖6中當(dāng)αg＞0.1后空化長度減小的原因之一；對比有無含氣的情形可以發(fā)現(xiàn)，由于含氣率導(dǎo)致的流場大范圍的可壓縮性，使得空泡尾端無法閉合在物面上，而是閉合在激波面上，如圖9（b）-（f）所示。

圖9 各含氣率下壓力等值線圖Fig.9 Absolute pressure contour line of each gas fraction

圖10為αg=0.1時蒸氣相體積分數(shù)云圖和Ma數(shù)云圖，其中圖10（b）為便于顯示Ma≥1的區(qū)域，云圖標(biāo)尺下限定為Ma=1。可以發(fā)現(xiàn)，雖然含氣率αg=0.1時，來流Ma數(shù)小于1，但是由于在蒸汽泡界面附近，空隙率接近0.5，由圖1可知，則空泡界面處馬赫數(shù)大于1為超聲速流，故如圖10（b）所示在流經(jīng)凸鈍角時產(chǎn)生了扇形膨脹波；空泡界面附近的流場膨脹加速后依舊為超聲速流，而由于空泡外形的影響，空泡尾端相當(dāng)于形成了凹鈍角，超聲速流動流經(jīng)該處形成了斜激波。此外，由于繞流流場的壓力分布和含氣率分布不均勻，導(dǎo)致各地聲速的分布也不均勻，最終導(dǎo)致的結(jié)果就是空泡尾端的激波面為曲面，如圖10（b）所示。

圖10 αg=0.1時蒸氣相體積分數(shù)云圖和Ma數(shù)云圖Fig.10 Vapor volume color contour and Mach number color contour of gas fraction 0.1

3.3 不同含氣率對空化水動力的影響

為研究不同含氣率對回轉(zhuǎn)體繞流流場空化水動力的影響，對含氣率αg=0～0.5工況的壁面壓力系數(shù)進行了對比，如圖11所示，其中壓力系數(shù)的定義為

圖11 各含氣率下壁面壓力系數(shù)曲線Fig.11 Pressure coefficient curve on wall of each gas fraction

圖11表明，隨著含氣率的增大，壓力系數(shù)最小的區(qū)域（即可能產(chǎn)生空化的區(qū)域）逐漸變短，其原因為隨著含氣率增大，混合介質(zhì)的可壓縮性增強（特別是頭部脫體激波的產(chǎn)生），導(dǎo)致分離面后的低壓也隨之升高，從而導(dǎo)致低壓力的空化區(qū)域減小。因此，隨著含氣率的增大，流場的空化效應(yīng)被削弱。另外，由圖11也可以看出，隨著流場含氣率的增大，空泡尾端的逆壓梯度逐漸減??；當(dāng)含氣率αg≥0.2時，逆壓梯度已變得相當(dāng)平緩；這印證了前面的結(jié)論，由于流場含氣帶來的大范圍可壓縮性的影響，空泡已經(jīng)無法閉合在物面上，而是閉合在空泡尾端的激波面上，從而體現(xiàn)為空泡尾端壁面逆壓梯度減小。

4 結(jié)論

本文對回轉(zhuǎn)體的可壓縮泡狀繞流進行了數(shù)值計算，其中含氣率αg=0的計算結(jié)果與試驗結(jié)果進行了對比，驗證了本文計算模型的合理性。各含氣率下的計算結(jié)果表明，隨著含氣率的增大，流場的可壓縮性也相應(yīng)增強，體現(xiàn)為數(shù)值模擬得到了回轉(zhuǎn)體繞流流場的3種波系，包括回轉(zhuǎn)體頭部處的脫體激波、分離面處的膨脹波和空泡尾端的斜激波。而空化與激波的相互作用則體現(xiàn)為：頭激波削弱了空化效應(yīng)，使得空化區(qū)域減??；由于空泡外形的影響，使得空泡尾端出現(xiàn)了斜激波；當(dāng)含氣率超過一定值時，空泡已經(jīng)無法閉合在物面上，而是閉合在空泡尾端的激波面上，體現(xiàn)為空泡尾端壁面逆壓梯度趨于平緩。本文的計算結(jié)果揭示了可壓縮泡狀流動的新的物理現(xiàn)象，有助于實驗結(jié)果的驗證。

[1]Karplus H B.The velocity of sound in a liquid containing gas bubbles[J].Armour Research Foundation of Illinois Institute of Technology,1958:C00-248.

[2]Saito Yoshinori,Takami Rieko,Nakamori Ichiro,Ikohagi Toshiaki.Numerical analysis of unsteady behavior of cloud cavitation around a NACA0015 foil[J].Computational Mechanics,2007,40(1):85-96.

[3]Eddington.Investigation of supersonic shock phenomena in a two-phase(liquid-gas)tunnel[J].AIAA,1970,08(01):65-74.

[4]Grant J R,Kirschner I N.High-speed motion in bubbly flows[C]//Fifth International Symposium on Cavitation(CAV2003). Osaka,Japan,2003:1-10.

[5]陳紅兵.定?？蓧嚎s泡狀流動激波的數(shù)值仿真[J].船舶力學(xué),2012,12(10):1107-1114. Chen Hongbing.Numerical study on the shocks in thesteady compressible bubbly flows[J].Journal of Ship Mechanics, 2012,12(10):1107-1114.

[6]ANSYS Inc.ANSYS FLUENT 12.0 Theory Guide[K].2009.

[7]褚學(xué)森,王志,顏開.自然空化流動數(shù)值模擬中參數(shù)取值影響的研究[J].船舶力學(xué),2007,11(1):32-39. Chu Xuesen,Wang Zhi,Yan Kai.Parametric study on numerical simulation of natural cavitation flow[J].Journal of Ship Mechanics,2007,11(1):32-39.

[8]權(quán)曉波,李巖,魏海鵬,王寶壽,孔德才.大攻角下軸對稱航行體空化流動特性試驗研究[J].水動力學(xué)研究與進展A輯,2008,23(6):662-667. Quan Xiaobo,Li Yan,Wei Haipeng,Wang Baoshou,Kong Decai.An experiment study on cavitation of underwater vehicle’s surface at large angles of attack[J].Journal of Hydrodynamics,Ser.A,2008,23(6):662-667.

Numerical simulation on interaction between cavitation and shock wave in compressible bubbly flows

DONG Hong-hui,WANG Bao-shou,CHEN Wei-qi,ZHANG Ke
(China Ship Scientific Research Center,Wuxi 214082,China)

The numerical simulation on the cavitating flow around an axis-symmetrical body in the steady compressible bubbly flow was simulated,and the interaction between cavitation and shock wave under different air volume fractions with the constant flow velocity and operating pressure was investigated.Firstly,to verify the used calculation model,the case which air volume fraction takes 0 in the incoming flow was calculated and shows a good agreement with the experimental result.Secondly,more cases with air volume fractions varying from 0 to 0.5 were calculated.The results show that with the increase of air volume fractions, flow field becomes more compressible.This result was confirmed by the appearance of three kinds of wave systems,including detached shock wave before the head,expansion wave at separation surface and oblique shock wave after cavity tail.The interaction between cavitation and shock wave is reflected as follows:bow shock before caviting region weakens cavitation effects,making caviting region decrease;cavity generated by cavitation results in occurrence of oblique shock;cavity closure point cannot locate on the body surface but on the shock surface,which can be confirmed by the adverse pressure gradient at the end of cavity.These results revealed a new physical phenomenon which existed in the interaction process between the cavity and the shock wave in compressible bubbly flows.

bubbly flow;shock wave;cavitation

U661.3

A

10.3969/j.issn.1007-7294.2015.11.001

1007-7294（2015）11-1295-09

2015-08-15

董洪輝（1989-），男，碩士研究生，E-mail：donghonghui@163.com；

王寶壽（1963-），男，研究員。