圓柱型后體對圓錐垂直入水空泡的影響

黃海龍，王柏秋

(1.哈爾濱工業(yè)大學(xué)土木工程學(xué)院，哈爾濱150090;2.北京宇航系統(tǒng)工程研究所，北京100076;3.哈爾濱工業(yè)大學(xué)建筑學(xué)院，哈爾濱150001)

導(dǎo)彈等物體從空中高速進(jìn)入水中的現(xiàn)象稱為入水問題，關(guān)于入水問題的研究，最早可以追溯到19 世紀(jì)末，A·M·Worthington[1－2]利用高速照相機對液滴和球體入水后所產(chǎn)生的流動現(xiàn)象進(jìn)行了試驗觀察，由此開始，各種針對入水空泡問題的探索不斷的展開.從文獻(xiàn)[3－4]中可以看出，在第二次世界大戰(zhàn)前后，多個國家相繼利用試驗方法和數(shù)值方法開展了針對各種導(dǎo)彈、空投魚雷等武器的入水問題的研究.其中，關(guān)于球體入水流場的試驗研究[5－9]和基于勢流理論研究[10－11]是較多的，其次是關(guān)于圓柱體等平頭物體[12－14]的研究以及具有圓錐頭的柱體模型[15]和射彈模型[16]入水流場的分析.

以上試驗研究和理論研究，分別得到了多種航行體的入水阻力系數(shù)與各種無量綱參數(shù)之間的關(guān)系，以及模型彈道等信息.但是，這些試驗研究并未討論模型后體參數(shù)的變化情況，例如錐柱體入水問題中，柱體段長度對模型入水空泡的影響.并且，由于研究方法和測試手段的局限性，這些研究成果只能反映入水空泡流場的部分參數(shù)信息.因此，隨著計算機技術(shù)的發(fā)展，利用數(shù)值計算，直接求解NS方程的方法逐漸得到了蓬勃發(fā)展[17－18]，這使得所論目標(biāo)更接近真實，能夠獲得的流場信息更豐富.

綜上所述，帶有后體的航行體入水空泡問題研究成果較少，而具有后體的航行體入水產(chǎn)生的空泡流場與不帶后體又具有一定差異.因此，針對后體對入水空泡流場的影響問題，本文利用動網(wǎng)格技術(shù)，采用VOF多相流模型和有限體積法求解RANS方程，對具有不同柱段長度的錐頭圓柱體勻速垂直入水過程進(jìn)行數(shù)值計算，分析模型勻速垂直入水過程中，空泡的生成和發(fā)展過程，以及有、無后體和后體長度對入水空泡的影響.

1 控制方程

1.1 試驗方法連續(xù)性方程和動量方程

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

空氣相連續(xù)性方程:

水蒸氣相連續(xù)性方程:

混合相動量方程:

其中:ui表示混合物速度分量，m/s;下標(biāo)m、l、g和v分別表示混合相、水、空氣和水蒸氣;ρm為混合物密度(kg/m3);μm是混合物動力黏性系數(shù)，m2/s.

其中:ρm和μm分別為:

其式中:αg、αv和 －αg、－αv分別為空氣、水蒸氣和水的體積分?jǐn)?shù).

1.2 空化模型

高速入水物體在入水過程中產(chǎn)生空泡，僅有物體頭部與水接觸，在物體頭部邊緣的分離處會出現(xiàn)較大的壓強梯度變化，當(dāng)壓強降低到飽和蒸汽壓強之下時，可能會出現(xiàn)自然空化現(xiàn)象.本文將這一復(fù)雜過程通過式(7)、(8)中的源項來表示(即Zwart－Gerber－Belamri空化模型).源項表達(dá)式如下:

其中:Re、Rc表示水蒸氣產(chǎn)生和凝結(jié);RB=10－6m為氣泡半徑;αnuc=5×10－4，為不可凝結(jié)氣體體積分?jǐn)?shù);Fvap、Fcond為經(jīng)驗常數(shù)，分別取 50 和 0.001.

1.2 湍流模型

本文計算中采用Realizable k－ε模型模擬湍流，該湍流模型中，湍動能k和湍能耗散率ε的輸運方程分別為:

其中:σk和σε分別為k和ε的Prandtl數(shù)，取值分別為 σk=1.0，σε=1.2;為由平均速度梯度所引起的湍動能生成項;經(jīng)驗常數(shù)C1和C2取值分別為1.44 和1.9.

2 數(shù)值求解

2.1 計算域和邊界條件

本文選取130°圓錐頭型為研究對象，模型后體圓柱直徑 D=10 mm.根據(jù)A·M·Worthington[1]的研究，軸對稱物體垂直入水所產(chǎn)生的空泡流場同樣具有軸對稱性質(zhì)，因此，本文數(shù)值計算均采用二維軸對稱計算模型，計算域設(shè)置和模型表面網(wǎng)格分布示意圖如圖1所示.圖1中，空氣域高度為60 D，水域深度為400 D，計算域徑向為80 D，使用四邊形網(wǎng)格劃分計算域.不失一般性，在保證計算精度和正確性的基礎(chǔ)上，選取網(wǎng)格密度和計算量適中的網(wǎng)格分布作為本文數(shù)值模型.

圖1 圓錐及計算域分布示意圖

2.2 數(shù)值方法

本文采用有限體積法對控制方程進(jìn)行空間離散，應(yīng)用PISO方法對動量方程方程和連續(xù)性方程聯(lián)立求解.計算過程中，對流項的離散采用二階迎風(fēng)格式;水蒸氣相的離散采用QUICK格式;湍流輸運方程的離散采用二階迎風(fēng)格式;對時間的離散采用一階顯示格式.計算時間步長的選取根據(jù)網(wǎng)格尺度、運動速度等因素綜合確定.通過動網(wǎng)格技術(shù)實現(xiàn)模型的運動過程與計算域內(nèi)的網(wǎng)格更新.

由于本文計算采用結(jié)構(gòu)化網(wǎng)格，網(wǎng)格主要是被拉伸或被壓縮變形，因此，采用動態(tài)層分裂法進(jìn)行網(wǎng)格更新.動態(tài)層分裂法的基本思想是當(dāng)網(wǎng)格被拉伸到大于設(shè)定的理想最大尺度((1+αs)hideal)時，該網(wǎng)格分裂成兩層新的網(wǎng)格，當(dāng)網(wǎng)格被壓縮到小于設(shè)定的理想最小尺度(αchideal)時，網(wǎng)格與其相鄰的網(wǎng)格合并成新的網(wǎng)格，即:

引入動態(tài)網(wǎng)格后，相當(dāng)于網(wǎng)格按一定速度作對流運動，在實際計算中應(yīng)減掉網(wǎng)格的對流效應(yīng).因而引入動網(wǎng)格后控制體V對變量φ的統(tǒng)一形式的守恒方程可寫為:

其中:ρm為混合物密度為速度矢量為動網(wǎng)格的運動速度;Γ為擴散系數(shù);Sφ為標(biāo)量φ的源項;?V表示控制體積V的邊界.

3 數(shù)值方法對比與驗證

為驗證方法的正確性和可信性，首先對直徑為25.4 mm球體勻速垂直入水過程進(jìn)行數(shù)值計算，并將計算結(jié)果的空泡輪廓與Albert May[6]等人提出的勻速垂直入水理想空泡模型進(jìn)行對比如下.

在二維軸對稱坐標(biāo)系下，垂直入水早期空泡形態(tài)的理想空泡模型表達(dá)式如式(14)所示，該入水空泡模型主要用于描述入水過程中，空泡閉合之前發(fā)展較為充分時期的外部輪廓.

其中:CD為入水物體的阻力系數(shù);d為入水物體的特征尺寸.

模型從空氣域勻速垂直入水，入水速度30 m/s，空氣域環(huán)境壓強為101 325 Pa，空氣、水、和水蒸氣的密度分別為 1.225、998.2、0.554 kg/m3.選取入水空泡充分發(fā)展時的空泡形態(tài)，與相同時刻下，由式(14)的理想空泡模型擬合結(jié)果進(jìn)行對比，結(jié)果如圖2所示，其中，阻力系數(shù)取CD=0.20.

圖2 無后體圓錐入水空泡輪廓計算結(jié)果與經(jīng)驗公式對比

由圖2可以看出，數(shù)值計算結(jié)果與理想空泡模型擬合結(jié)果具有較好的一致性.在空泡輪廓兩端具有一定的誤差，產(chǎn)生該誤差的主要原因是由于理想空泡模型本身忽略了自由表面處和分離點處的空泡輪廓的不規(guī)則性.通過圖2的對比驗證可見，本文所采用的數(shù)值計算方法是正確可信的，在此基礎(chǔ)之上，下面針對帶有不同尺度圓柱后體的圓錐模型，展開勻速垂直入水空泡的數(shù)值計算.

4 計算結(jié)果分析

計算模型外形如圖3所示，頂部錐角 a=130°，直徑 D=10mm，后體長度 L分別為:0(圓錐)，0.5D，1.0D，2.0D 和 4.0D.為考察后體對入水空泡的影響，將各模型入水計算結(jié)果的某一時刻流場分布顯示于圖4中.

圖3 不同長度后體模型

圖4 入水空泡內(nèi)部組分分布

有后體后，可以發(fā)現(xiàn)入水空泡的最大直徑隨著后體尺度有一定變化.圖5表示垂直入水空泡最大直徑隨后體長度的變化過程，入水速度為200 m/s.從圖5可以看出，模型后體長度較小時(例如圖5中的無后體模型和后體長度為L=0.5D模型)，入水空泡最大直徑略高于其余帶后體模型，隨著模型后體長度的增加，入水空泡最大直徑逐漸降低，當(dāng)后體長度增加到一定程度后，入水空泡最大直徑趨于穩(wěn)定，不再隨后體長度變化而變化.

圖5 空泡最大直徑隨模型后體的變化

圖4中的入水空泡尾部的閉合情況清晰表明入水空泡閉合時間隨后體的出現(xiàn)而有所變化.圖6表示模型勻速垂直入水過程中，入水空泡尾部無量綱閉合時間隨后體長度的變化情況，模型入水速度為200 m/s.

圖6 入水空泡閉合時間隨后體的變化

從圖6可以看出，模型后體長度較小時，空泡閉合時間較長，尤其當(dāng)后體長度L=0時，閉合時間最長;當(dāng)模型后體長度較大時，空泡閉合時間趨于一致，例如當(dāng)L＞2.0時，空泡閉合時間不再隨后體長度變化而變化.由圖5、6可以看出，帶后體模型產(chǎn)生的入水過程，當(dāng)模型后體長度大于某一數(shù)值時，流場是相近的.為考察穩(wěn)定的入水過程，下面以L=4.0D的模型為例，進(jìn)行分析.

圖7為長后體模型(L=4.0D)勻速垂直入水過程中，水蒸氣和空氣在空泡內(nèi)部的空間分布的變化情況，其中，模型入水速度為200 m/s.從圖7可以看出，在入水初期，模型肩部的流動分離點附近和空泡入口處的噴濺內(nèi)部邊緣有少量的水蒸氣.隨著入水過程的進(jìn)行，空泡入口的噴濺液體迅速向內(nèi)彎曲，空泡開始趨于閉合.此時，在空泡入口處的彎曲液壁內(nèi)側(cè)，由于液面迅速彎曲，附近流體流速增高，靜壓力急劇下降，導(dǎo)致液體迅速空化，從而使得該處水蒸氣含量隨著空泡閉合而迅速增加，隨著時間的推移，空化逐漸沿著空泡壁面向空泡內(nèi)部傳播.在模型附近，既空泡頭部附近，空泡輪廓線穩(wěn)定，不隨時間變化，因此，該處流場內(nèi)部壓力和流速均穩(wěn)定，發(fā)生空化的概率相對較低，該處只存在少量入水沖擊產(chǎn)生的水蒸氣或者不存在水蒸氣.

另外，由圖7的入水過程可知，入水空泡內(nèi)部流體壓力低于外界大氣壓和周圍流體內(nèi)部靜壓，當(dāng)空泡壁運動能量不足時，將導(dǎo)致空泡壁向軸線坍縮，從而擠壓內(nèi)部空氣，最終可能導(dǎo)致一股沿空泡軸線向上的液體噴濺的產(chǎn)生.

圖7 不同時刻入水空泡組分分布

5 結(jié)論

本文對具有不同長度圓柱型后體的圓錐勻速垂直入水過程進(jìn)行了數(shù)值模擬研究，分析了有、無后體和不同長度后體模型入水后，空泡的發(fā)展過程及空泡內(nèi)部流場的分布，得到了以下結(jié)論:

1)模型在勻速垂直入水后形成的入水空泡內(nèi)，無后體模型流場內(nèi)的空化作用弱于有后體模型流場，入水空泡內(nèi)的水蒸氣含量明顯偏低.

2)帶圓柱后體模型，后體的存在有利于模型肩部流動分離點附近自然空化的發(fā)生，并有少量的水蒸氣分布在分離點附近的空泡壁面.

3)無后體模型入水空泡的閉合明顯晚于帶圓柱后體模型入水空泡的閉合，當(dāng)圓柱后體長度小于1倍圓柱體直徑時，空泡最大直徑隨著后體長度的增加而減小，且空泡最大直徑略大于具有更長帶后體的模型;當(dāng)后體長度大于2倍以上圓柱體直徑時，入水空泡最大直徑與后體長度無關(guān).

4)無后體時，入水空泡的閉合時間較長;有后體時，后體長度大于兩倍圓柱直徑時，入水空泡的閉合時間與后體長度無關(guān).

通過本文研究可以看出，模型后體的存在對入水空泡形態(tài)及相關(guān)參數(shù)有重要影響，為了使入水問題研究結(jié)果在魚雷和導(dǎo)彈等武器中得到更好的應(yīng)用，應(yīng)加強對模型后體影響的研究.

[1]WORTHINGTON A M.A study of splashes[M].New York:Longmans Green and Company，1908.

[2]WORTHINGTON A M，COLE R S.Impact with a liquid surface[J].Phil.Trans.Roy.Soc.，1897，189A:137 －148.

[3]WHALLEY I A.Project upkeep－a review of the WWⅡdambubuster weapon[C]//Indianapolis:38th AIAA Joint propulsion Conf，2002.

[4]ALBERT M.Water entry and the cavity－running behavior of missiles[R].AD A020429，1975.

[5]MALLOCK A.Sounds produced by drops falling on water[J].Pro.R.Soc.Lond.，1918，A95:138 －143.

[6]BELL G E.On the impact of a solid sphere with a fluid surface

[J].Phil.Mag.j.Sci，1924，48:753 － 765.

[7]ASHLEY S.Warp drive underwater[J].Sci.Amer，2001，284(5):70－79.

[8]MAY A，WOODHULL J C.Drag coefficients of steel spheres entering water vertically[J].Journal of Applied Physics，1948，19(12):1109－1121.

[9]TRUSCOTT T T，TECHET A H.Water entry of spinning spheres[J].J.Fluid Mech，2009，625:135 －165.

[10]Yan H，Liu Y，KOMINIARCZUK J，et al.Cavity dynamics in water entry at low Froude numbers[J].J.Fluid Mech，2009，641:441－461

[11]ARISTOFF J M，BUSH J W M.Water entry of small hydrophobic spheres[J].J.Fluid Mech，2009，619:45 －78.

[12]SHI H，ITOH M，TAKAMI T.Optical observation of the supercavitation induced by high－speed water entry[J].Journal of Fluids Engineering，2000，122:806－810.

[13]何春濤，王 聰，何乾坤，等.圓柱體低速入水空泡試驗研究[J].物理學(xué)報，2012，61(13):1－8.

[14]陳學(xué)農(nóng)，何友聲.平頭物體三維帶空泡入水的數(shù)值模擬[J].力學(xué)學(xué)報，1990，22(2):129－138.

[15]王 聰，何春濤，權(quán)曉波，等.空氣壓強對垂直入水空泡影響的數(shù)值研究[J].哈爾濱工業(yè)大學(xué)學(xué)報，2012，44(5):14－19.

[16]MAY A.Effect of surface condition of a sphere on its waterentry cavity[J].Journal of Applied Physics，1951，22:1219－1222.

[17]HE C，WANG C，WEI Y，et al.Numerical simulation of pressure distribution in vertical water－ entry cavity[J].Journal of Ship Mechanics，2011，15(9):960 －968.

[18]何春濤，王 聰，閔景新，等.回轉(zhuǎn)體勻速垂直入水早期空泡數(shù)值模擬研究[J].工程力學(xué)，2012，29(4):237－243.