近自由液面氣泡與沖擊波的相互作用

葉曦，初文華，陳林，張阿漫

1中國船舶及海洋工程設計研究院，上海200011

2上海海洋大學海洋科學學院，上海201306

3上海海洋大學國家遠洋漁業工程技術研究中心，上海201306

4上海海洋大學大洋漁業資源可持續開發教育部重點實驗室，上海201306

5哈爾濱工程大學船舶工程學院，黑龍江哈爾濱150001

近自由液面氣泡與沖擊波的相互作用

葉曦1，初文華2，3，4，陳林1，張阿漫5

1中國船舶及海洋工程設計研究院，上海200011

2上海海洋大學海洋科學學院，上海201306

3上海海洋大學國家遠洋漁業工程技術研究中心，上海201306

4上海海洋大學大洋漁業資源可持續開發教育部重點實驗室，上海201306

5哈爾濱工程大學船舶工程學院，黑龍江哈爾濱150001

［目的］為了研究自由液面、氣泡與沖擊波三者之間的相互作用，［方法］基于間斷迦遼金法，結合Level Set與Real Ghost Fluid方法，分析復雜流場中沖擊波傳播特性及氣泡運動特性，描述流場內各種波系的生成與發展過程。［結果］結果表明：在相互作用過程中，流場生成的復雜波系中包含多個稀疏波和沖擊波。自由液面減緩了氣泡的潰滅速度，而入射沖擊波則加快了氣泡的潰滅速度，并使自由液面的上拱運動增大。［結論］所得結果可為水下爆炸對艦船結構的毀傷機理提供參考。

氣泡；自由液面；沖擊波；間斷迦遼金法

Abstract：［Objectives］This paper presents research into the interaction between a free surface，shock wave and bubble.［Methods］Based on the Discontinuous Galerkin method and combined with the Level Set method and Real Ghost Fluid method，the characteristics of shock wave propagation and bubble motion are analyzed，the generation and evolution of waves in fluid fields are described in detail.［Results］The results show the presence of complex waves in the fluid field after the interaction，including multiple rarefaction waves and shock waves.The collapse speed of the bubble is slowed down with the existence of the free surface，while the shock wave of the incident accelerates the collapse speed of the bubble and increases the up warp motion of the free surface.［Conclusions］The conclusions drawn from this paper can be used as reference points for further research into the damage mechanisms of ship structures subjected to underwater explosions.

Key words：bubble；free surface；shock wave；Discontinuous Galerkin method

0 引 言

在船舶與海洋工程領域，經常會出現沖擊波、氣泡與自由液面間的相互作用，例如，水下爆炸或結構反射的沖擊波與爆炸產生的氣泡，以及船舶航行過程中在自由液面附近產生的空泡（船舶尾流中存在大量空泡）之間的相互作用。這種相互作用會改變沖擊波的傳播特性并影響氣泡的運動特性，由此導致艦船結構載荷特性發生變化。因此，分析沖擊波、自由液面以及氣泡之間的相互作用對于深入研究水下爆炸工況的艦船結構毀傷機理具有重要意義。

由水中氣泡界面兩側介質差異造成的流場間斷將給數值模擬帶來一定的困難。Johnsen等［1-2］采用體積分數法結合有限體積法進行加權本質無振蕩（Weighted Essentially Non-Oscillatory，WENO）格式重構，計算了水中氣泡在強沖擊波和弱沖擊波作用下的演化歷程，并考慮了氣泡下游存在壁面的情況。王革和關奔［3］采用水平集（Level Set）方法結合虛擬流體（Ghost Fluid）方法，并采用有限體積法求解Euler方程組，分析了水中沖擊波與氣泡的相互作用。此外，自由拉格朗日法［4］、邊界元法［5］、間斷迦遼金法［6-10］等也適用于該物理現象的求解。然而，目前尚無自由液面、氣泡與沖擊波三者相互作用的研究成果。

本文將基于間斷迦遼金法，結合LevelSet方法［11］、Real Ghost Fluid方法［12］以及自適應網格加密（Adaptive Mesh Refinement，AMR）技術［13］研究水下氣泡、沖擊波、自由液面之間的相互作用，分析自由液面對氣泡運動及流場特性的影響，詳細闡述新波系的產生與發展過程。

1 基本理論與方法

如圖1所示，計算模型采用zor軸對稱坐標系。在初始狀態，氣泡與自由液面相對靜止，氣泡半徑為R，氣泡中心與自由液面的初始距離為df。強度為M的沖擊波波前與自由液面平行，與氣泡中心的距離為ds，自下而上沖擊氣泡。本文所有參數均進行無量綱化處理，參考量分別為氣泡初始半徑、標準大氣壓和空氣密度。由于本文主要研究流場中的沖擊波特性，故假設流場無粘且初始無旋，并采用軸對稱Euler方程描述流場：

其中

式（1）～式（2）中：t為時間；ρ為流場密度；p為流場壓力；u，v分別為z，r軸方向的速度；，為總能量（內能與動能之和），其中ein為單位質量的內能。

圖1 計算模型Fig.1 Numerical model

氣泡內部采用理想氣體狀態方程［14］描述：

式中：γg=1.4，為氣體比熱比常數。

氣泡外流場采用Tait狀態方程［6］進行描述：

式中：N和Pw為介質常數，N=7.15，Pw= 3.31×108Pa。

采用間斷迦遼金法求解Euler方程，半離散格式為

式中：pm為多項式的最高次數；Ω為單元面積；Φl(l=0，1，…，pm)和Ψk(k=0，1，…，pm)為試函數；Ek為第k次多項式中Ψk的系數；和為單元邊界的法向量；A為單元邊界。

選取Legendre正交多項式作為試函數，則流場中的變量E為

式（5）采用三階Runge-Kutta顯式格式進行時間步進，采用HLLC格式［15］計算數值通量F和G，通過斜率限制器［16］抑制高精度離散格式帶來的間斷處非物理震蕩。為了提高計算精度與效率，采用AMR技術構建多層網格（圖2），在流場中壓力、密度等物理量出現間斷的位置進行局部加密。

圖2 AMR多層網格示意圖Fig.2 Diagram for AMR multiple levels grid

2 數值方法驗證

圖3所示為水中沖擊波與柱狀氣泡相互作用的數值計算結果對比。介質界面兩側擁有較大的密度與壓力差，沖擊波強度較大，對間斷捕捉算法的計算穩定性及精度要求較高。計算區域取為24 mm×24 mm（此處仍采用帶量綱的物理量，以便與文獻結果［17-19］進行對比），底層網格數為150× 150，采用3層AMR網格。柱狀氣泡半徑為3 mm，氣泡中心位于（12 mm，12 mm）處，沖擊波馬赫數為1.72。氣泡密度為1 kg/m3，壓力為1 bar，氣體常數γg=1.4；氣泡外部水的密度為1 000 kg/m3，壓力為1 bar，介質常數N=4.4，Pw=6×108Pa。沖擊波初始位置與氣泡中心相距54 mm。計算域的左邊界為入口邊界，右邊界和上、下邊界均為無反射邊界。本文將計算結果與文獻［17-19］的結果進行了對比，雖然本文采用的底層網格數和AMR優化層數均小于文獻［17］的結果，但仍能較好地模擬水中氣泡與沖擊波之間的相互作用。

圖3 水中沖擊波與氣泡相互作用的數值計算結果對比Fig.3 Results comparison for numerical calculation about shock-bubble interaction in water

3 數值計算分析

為了簡化說明，設SWi為入射波；SWb為水下爆炸產生的沖擊波；SWir為SWi在自由液面的反射波；SWir*為SWi在氣泡上游表面的反射波；SWirr為SWir在氣泡下游表面的反射波；SWj為生成環狀氣泡時產生的沖擊波；SWjr為SWj在自由液面的反射波；SWt為氣泡內部的透射沖擊波。

首先，考慮自由液面附近沖擊波與初始內外壓平衡的氣泡之間的相互作用。氣泡中心與自由液面的初始距離為df=1.5，沖擊波波前與氣泡中心的距離為ds=2。計算域大小為3×6，采用3層 AMR網格，最底層網格數量為120×240，下邊界為入口邊界條件，上邊界為出口邊界條件，左、右邊界為壁面邊界條件。各介質內的物理量記為［ρ，p，u，v，γ］，則氣泡內部參數為［1，1，0，0，1.4］，外部流場中水域參數為［1 000，1，0，0，7.15］，空氣域參數與氣泡內部相同，沖擊波強度M=1.1，沖擊波與近自由面氣泡的相互作用如圖4所示。由于流場中未采用空化模型，因此壓力存在負值。

由圖4（a）和圖4（b）可知，SWi傳播至氣泡表面后分為3部分：一部分被反射形成稀疏波；另一部分透射至氣泡內部并引起氣泡內部氣體的高速運動；剩下的一部分則繞過氣泡繼續向前傳播。

圖4 沖擊波（M=1.1）與近自由面氣泡的相互作用（分圖中左側為壓力分布，右側為馬赫場分布）Fig.4 Shock-bubble interaction near free surface for M=1.1（In figures：left are pressure distribution，and right are Mach number field distribution）

由圖4（c）和圖4（d）可知，當SWi到達自由液面后，一方面將從水中入射空氣，由于水的特征阻抗高于空氣，因此SWi以稀疏波的形式被反射回來，形成SWir；另一方面，反射回來的SWir作用于氣泡表面，也從水中入射空氣，因此稀疏波SWir將以沖擊波的形式被反射回來，形成SWirr。以此類推，在氣泡與自由液面之間，沖擊波與稀疏波將不斷地往復傳播。此時，兩側自由液面在沖擊波的作用下將出現上拱，而氣泡上游表面將內凹形成射流，且在射流附近出現局部高壓區域。

由圖4（e）和圖4（f）可知，此時氣泡在沖擊波誘導下形成的射流已突破下游表面形成環狀氣泡，形成新的沖擊波SWj，并向四周傳播。由于SWj的強度遠高于SWi，因此SWj在自由液面處反射形成的稀疏波SWjr所造成的負壓也遠高于SWir，此時，自由液面中心在SWj的作用下出現上拱。而在環狀氣泡的上游則出現了若干局部高壓區域，這是由氣泡形成射流時在上游產生的局部高壓與向上游傳播的SWj相交所產生。此外，當SWi和 SWj作用于自由液面時，有低強度沖擊波透射入氣相，將引起液面另一側的氣體運動。

圖5所示為不同時刻r=0處的壓力分布。圖5（a）中單向箭頭與各曲線相交的時間點從左至右分別為t=0.402，0.702，0.772，0.816，0.852，0.916，1.030，1.145；圖5（b）中t1=0.702，t2=0.744，t3=0.775，t4=0.809，t5=0.863，t6=0.954，t7=0.982，t8= 1.016。SWi到達自由面后將以稀疏波SWir的形式反射回來，并形成局部負壓。由圖5（b）可知，當稀疏波SWir作用在氣泡下游表面后，將反射回沖擊波SWirr，使附近流場逐漸改變為正壓狀態。而SWirr將再次到達自由面，并反射回稀疏波作用于氣泡下游表面，最終再以沖擊波的形式被反射。此時，氣泡內部的沖擊波SWt已接近氣泡下游表面。圖5（b）中，t8時刻曲線即為氣泡射流沖擊下游表面后在沖擊處附近產生的局部高壓，并形成沖擊波SWj。當SWj作用于自由液面時也會形成相應的稀疏波SWjr，并形成局部負壓。由于SWj的強度遠高于SWi，則SWjr所造成的負壓也遠高于SWir，如圖5（a）所示（t=1.145）。

圖5 r=0處的壓力分布Fig.5 Pressure distribution at r=0

圖6所示為有、無自由液面時的氣泡動態特性。自由液面反射的稀疏波作用于氣泡表面時會使氣泡下游界面發生朝向下游的運動，導致相應時刻的氣泡體積增大，則形成最小體積的時間更長。由于反射波的作用時間較短，因此其產生的影響也較小。由于射流是由背向自由液面的氣泡上游表面所形成，故自由液面的反射波對其影響較小，因此有、無自由液面時，氣泡的射流頂點速度基本一致。

圖6 有、無自由液面時的氣泡動態特性Fig.6 Dynamics of bubble with and without free surface

水下多點爆炸比較常見，在產生氣泡的同時向四周傳播沖擊波，該沖擊波與其他方向傳播而來的沖擊波相互耦合，進而對氣泡運動及流場特性造成影響。沖擊波（M=1.1）與初始內部高壓的氣泡在自由液面附近的相互作用如圖7所示，計算域大小為8×12，采用4層AMR網格，最底層網格數為120×180。氣泡內部介質參數為［1 630，1 500，0，0，1.25］，其他參數與前文相同。

圖7 沖擊波（M=1.1）與初始內部高壓的氣泡在自由液面附近的相互作用（分圖中左側為壓力分布，右側為馬赫場分布）Fig.7 Shock-bubble interaction with initial high inner pressure near free surface for M=1.1（In figures：left are pressure distribution，and right are Mach number field distribution）

由圖7（a）可知，氣泡釋放的沖擊波SWb與入射波SWi在相匯處疊加，形成高壓區。SWi與氣泡上游表面接觸后發生反射形成稀疏波SWir*，與先前的高壓區疊加形成高壓帶狀區并繼續向四周傳播。在氣泡與自由液面之間的區域，SWb在自由液面處以稀疏波的形式被反射，隨后在氣泡下游表面以沖擊波的形式被反射，如此往復，形成了間隔的沖擊波區域與稀疏波區域。在初始時刻，氣泡除了向四周釋放沖擊波SWb外，還在氣泡內部形成內聚稀疏波，最終在氣泡中心形成高壓區，如圖7（b）所示。

圖8所示為不同初始條件下氣泡（氣泡內部初始壓力1 500）與自由液面形態隨時間的變化。與沒有入射沖擊波時相比，在SWi的作用下，同時刻的氣泡體積減小，氣泡膨脹的同時向自由液面偏移，自由液面兩側發生明顯的上拱運動。當入射沖擊波SWi強度較大時，自由液面的上拱運動更明顯，氣泡下游表面較快地形成明顯的尖峰，并與自由液面發生融合。

圖8 不同初始條件下的氣泡與自由液面形態Fig.8 Shape of bubble and free surface with different initial conditions

4 結 論

本文采用間斷迦遼金法，結合Level Set和Real Ghost Fluid方法求解近自由液面氣泡與沖擊波間的相互作用，詳細分析了氣泡的運動特性、各類沖擊波與稀疏波的生成和發展過程以及流場中壓力、速度等物理量的分布特性。主要結論如下：

1）在氣泡、沖擊波與自由液面三者相互作用的過程中，流場內存在多個復雜的稀疏波系和沖擊波系。

2）當氣泡初始內外壓平衡時，氣泡在入射沖擊波的作用下潰滅形成環狀氣泡，在此過程中產生的沖擊波強度遠高于入射沖擊波。在該沖擊波作用下，自由液面中心形成向上凸起的水冢。當自由液面存在時，氣泡的運動過程與自由場中相似，但是在自由液面反射的稀疏波作用下，氣泡的潰滅速度減緩。

3）當氣泡內部存在一定的初始壓力時，氣泡運動開始時釋放的沖擊波與入射沖擊波發生耦合，形成更為復雜的流場分布。在入射沖擊波作用下，氣泡并未形成環狀氣泡，且與無入射波的情況相比，氣泡向自由液面偏移，體積減小，自由液面的上拱運動增大。

［1］JOHNSEN E，COLONIUS T.Shock-induced collapse of a gas bubble in shockwave lithotripsy［J］.The Journal of the Acoustical Society of America，2008，124（4）：2011-2020.

［2］JOHNSEN E.Numerical simulation of non-spherical bubble collapse：with applications to shockwave lithotripsy［D］.California：California Institute of Technology，2007.

［3］王革，關奔.激波作用下R22氣泡射流現象研究［J］.力學學報，2013，45（5）：707-715.WANG G，GUAN B.A study on jet phenomenon of R22gas cylinder under the impact of shock wave［J］.Chinese Journal of Theoretical and Applied Mechanics，2013，45（5）：707-715（in Chinese）.

［4］JAMALUDDIN A R，BALL G J，TURANGAN C K，et al.The collapse of single bubbles and approximation of the far-field acoustic emissions for cavitation induced by shock wave lithotripsy［J］.Journal of Fluid Mechanics，2011，677：305-341.

［5］KLASEBOER E，FONG S W，TURANGAN C K，et al.Interaction of lithotripter shockwaves with single inertial cavitation bubbles［J］.Journal of Fluid Mechanics，2007，593：33-56.

［6］PARK J.A coupled runge kutta discontinuous galerkin-directghostfluid（RKDG-DGF） method to near-field early-time underwater explosion（UNDEX）simulations［D］.Blacksburg：Virginia Polytechnic Institute and State University，2008.

［7］COCKBURN B，HOU S C，SHU C W.The runge-kutta local projection discontinuous Galerkin finite element method for conservation laws IV：the multidimensional case［J］.Mathematics of Computation，1990，54（190）：545-581.

［8］COCKBURN B，SHU C W.The runge-kutta discontinuous Galerkin method for conservation laws V：multidimensional systems［J］.Journal of Computational Physics，1998，141（2）：199-224.

［9］QIU J X，LIU T G，KHOO B C.Simulations of compressible two-medium flow by runge-kutta discontinuous Galerkin methods with ghost fluid method［J］.Communications in Computational Physics，2008，3（2）：479-504.

［10］葉曦.水下氣泡運動的噪聲及近場載荷特性研究［D］.哈爾濱：哈爾濱工程大學，2014.YE X.Study on noise and near-field load characteristics induced by bubble motion［D］.Harbin：Harbin Engineering University，2014（in Chinese）.

［11］OSHER S，FEDKIW R.Level Set methods and dynamic implicit surfaces［M］.New York：Springer，2002.

［12］WANG C W，LIU T G，KHOO B C.A real ghost fluid method for the simulation of multimedium compressible flow［J］.SIAM Journal on Scientific Computing，2006，28（1）：278-302.

［13］BERGER M J，COLELLA P.Local adaptive mesh refinement for shock hydrodynamics［J］.Journal of Computational Physics，1989，82（1）：64-84.

［14］LANEY C B.Computational gasdynamics［M］.Oxford：Cambridge University Press，1998.

［15］HARTEN A，LAX P D，LEER B V.On upstream differencing and Godunov-type schemes for hyperbolic conservation laws［J］.SIAM Review，1983，25（1）：35-61.

［16］COCKBURN B，LIN S Y，SHU C W.TVB runge-kutta local projection discontinuous Galerkin finite element method for conservation laws III：one-dimensional systems［J］.Journal of Computational Physics，1989，84（1）：90-113.

［17］NOURGALIEV R R，DINH T N，THEOFANOUS T G.Adaptive characteristics-based matching for compressible multifluid dynamics［J］.Journal of Computational Physics，2006，213（2）：500-529.

［18］BALL G J，HOWELL B P，LEIGHTON T G，et.al.Shock-induced collapse of a cylindrical air cavity in water：a free-lagrange simulation［J］.Shock Wave，2000，10（4）：265-276.

［19］BOURNE N K，FIELD J E.Shock-induced collapse of single cavities in liquids［J］.Journal of Fluid Mechanics，1992，244（1）：225-240.

Interaction between bubble near free surface and shock wave

YE Xi1，CHU Wenhua2，3，4，CHEN Lin1，ZHANG Aman5
1 Marine Design and Research Institute of China，Shanghai 200011，China
2 College of Marine Sciences，Shanghai Ocean University，Shanghai 201306，China
3 National Engineering Research Center for Oceanic Fisheries，Shanghai Ocean University，Shanghai 201306，China
4 Key Laboratory of Sustainable Exploitation of Oceanic Fisheries Resources，Ministry of Education，Shanghai Ocean University，Shanghai 201306，China
5 School of Shipbuilding Engineering，Harbin Engineering University，Harbin 150001，China

U661.144

A

10.3969/j.issn.1673-3185.2017.05.011

2016-10-31< class="emphasis_bold">網絡出版時間：

時間：2017-9-26 10:43

國家自然科學基金資助項目（51479041，51279038）；國家自然科學基金青年科學基金資助項目（11402143）；上海市高校青年教師培養資助計劃項目（A1-2035-14-0010-18）；大洋漁業資源可持續開發教育部重點實驗室開放基金資助項目（A1-0203-16-2007-5）

葉曦，男，1987年生，博士，工程師。研究方向：艦船振動噪聲與沖擊。E-mail：yexi0527@163.com

初文華（通信作者），女，1986年生，博士，講師。研究方向：流體動力學。E-mail：whchu@shou.edu.cn

陳林，男，1987年生，碩士，工程師。研究方向：艦船振動噪聲與沖擊。E-mail：1021763445@qq.com

張阿漫，男，1981年生，博士，教授。研究方向：艦船抗爆抗沖擊。E-mail：zhaman@163.com

http://kns.cnki.net/kcms/detail/42.1755.TJ.20170926.1043.016.html期刊網址：www.ship-research.com

葉曦，初文華，陳林，等.近自由液面氣泡與沖擊波的相互作用［J］.中國艦船研究，2017，12（5）：90-96.

YE X，CHU W H，CHEN L，et al.Interaction between bubble near free surface and shock wave［J］.Chinese Journal of Ship Research，2017，12（5）：90-96.