基于VOF法的平底結構自由落體入水砰擊載荷模擬

王易君，李明海，張中禮，李上明，史光梅

(中國工程物理研究院總體工程研究所，四川 綿陽 621900)

基于VOF法的平底結構自由落體入水砰擊載荷模擬

王易君，李明海，張中禮，李上明，史光梅

(中國工程物理研究院總體工程研究所，四川 綿陽 621900)

應用動網格技術和VOF方法模擬了平底結構的自由落體入水過程和氣-液兩相流，得到了平底物體自由落體入水過程的運動特性、砰擊載荷特性以及入水空泡演化過程。結果表明：整個入水過程中，平底物體在入水瞬間受到來自液面的砰擊載荷最大，砰擊面的載荷壓力分布中間高兩端低；砰擊載荷隨入水初速的增加而增大；平底物體入水空泡破裂后對上表面造成沖擊壓力。數值計算方法和結果對入水結構安全性和適應性評估具有實際意義。

VOF；入水；砰擊載荷；空泡；數值模擬

入水砰擊載荷和入水空泡問題是入水結構安全性和環境適應性的主要研究內容，其中砰擊載荷大小是決定結構沖擊響應的關鍵參數，而入水空泡的產生、發展和湮滅等演化過程對入水結構的姿態和彈道穩定性具有明顯影響。因此，入水砰擊載荷和入水空泡得到了廣泛的研究[1-3]。近年來，針對入水結構形狀分別開展了圓柱垂直入水空泡[4]、圓柱傾斜入水空泡[5]、含涂層鋼球垂直入水空泡[6]等試驗研究，得到入水過程中空泡的生成、發展和閉合等變化過程。陳偉琪等[3]對物體垂直出入水的非定?？张輸祵W模型進行了研究，得到勻速、空泡壓力不變條件下空泡外形、尺寸隨水深的解析解。何春濤等[7]對球體勻速垂直入水早期空泡形態進行數值計算，得到了垂直勻速入水空泡形態隨時間的變化規律以及空泡表面閉合時間同入水速度的關系。陳宇翔等[8]對水平圓柱垂直入水進行了數值模擬，得到入水初期自由表面的變化，并比較了湍流黏性對入水的影響。

相對楔形體、水平圓柱、彈丸等入水過程，平底結構垂直入水時入水端同初始液面相互平行，在砰擊水面瞬間可視為“面接觸”，因而在入水過程中有明顯的“空氣墊”現象。陳學農等[9]用時間步進法和邊界積分方程對平頭物體垂直入水進行了數值求解，得到垂直入水初期空泡形態。HUERA[10]對平底板以一定角度入水的斜向入水砰擊進行了試驗研究，當入水角度小于5°時，空氣墊現象非常明顯。CHEN等[11]針對平底結構入水砰擊階段的砰擊壓力進行了試驗和數值分析，其中數值分析采用有限元法，在Chuang給出的預估平板入水峰值壓力的關系基礎上，給出了比例系數同入水初速的擬合關系式。楊衡等[12]對彈性結構入水砰擊載荷特性進行了三維數值模擬，基于邊界元法和雙漸近法構建入水砰擊流固耦合的數值模型，但未能精確模擬流體飛濺、自由液面大變形。

在專業的流體分析軟件FLUENT平臺上，基于VOF方法和有限體積法求解氣、液兩相流動的RANS方程，結合動網格技術模擬物體入水運動，研究平底物體自由落體入水過程中的砰擊載荷特性，并對其入水過程的空泡演化過程進行分析。

1 研究對象及計算模型

1.1 研究對象

以1 000 mm×600 mm的2維平底鋼結構為研究對象，其幾何結構示意圖如圖1所示，h為物體下落高度。

圖1 幾何結構示意圖Fig.1 Sketch of the model

1.2 控制方程及邊界條件

VOF模型將控制單元內的流體看作密度可變的單相流體，利用計算網格中流體體積量的變化和網格單元本身體積的比值函數α來確定自由面的位置和形狀。

對于每一個控制單元，其連續性方程為

(1)

(2)

式中：ui為i(i=x,y)方向上速度分量；ρm為單元內流體密度，VOF模型中ρm=α1ρ1+(1-α1)ρ2，其中下標1、2分別表示空氣和水，α1、α2分別為空氣相和水相的體積分數。當α1=0時，該控制單元內全部被水相占據；當α1=1時，該控制單元完全被氣相占據；當0<α1<1時，該單元為包含兩相物質的交界面。

動量方程為

(3)

式中：p為流體壓強，μm為單元動力黏性系數，μm=α1μ1+(1-α1)μ2，μt為湍流黏性系數。

湍流方程采用標準k-ε湍流模型進行求解。

數值模擬時計算域選取x方向10倍研究對象特征尺寸，y方向10倍研究對象特征尺寸的進行計算。計算域空氣上方為壓力出口邊界，水域底部和計算域兩側為固壁邊界條件。初始時刻，物體距離未擾動水面高度為h。

1.3 數值計算

采用GAMBIT建立二維仿物理數值入水模型，結構化網格進行劃分，并在平底結構附近進行網格加密。由于運動速度較低，空氣和水均視為不可壓縮流體。編寫用戶自定義函數(UDF)描述物體自由落體運動及入水后運動速度變化，并結合動態分層算法動網格技術實現物體入水過程。物體自由落體以6DOF模型實現，通過設置物體質量、轉動慣量從而實現物體自由下落。求解基于VOF多相流模型的RANS方程，其空間、時間離散采用一階迎風格式，壓力項采用PRESTO！離散格式，動量方程和湍流方程均采用一階迎風格式，自由液面重構采用Geo-Reconstruct格式。

2 數值計算結果

2.1 數值方法驗證

同文獻[13]中平底體模型入水試驗數據進行比較分析，下落高度為0.3 m。建立了最小網格尺寸分別為40 mm、20 mm、10 mm的二維平底結構入水計算模型，其網格數量分別為9 588、39 890、155 520，計算得到入水砰擊過程中入水面中間點的壓力時程曲線，如圖2所示，其中pMax為峰值壓力。同試驗結果相比，計算所得的砰擊峰值壓力相對誤差分別為11.4%、17.9%、11.0%，所采用的網格尺寸對計算結果的敏感性較低。

圖2 入水過程中入水端面壓力時程曲線Fig.2 The curve of pressure-time on slamming side

圖3為最小網格尺寸為20 mm，不同時間步長和計算砰擊壓力的曲線。時間步長分別為0.001 s、0.000 5 s、0.000 2 s、0.000 1 s、0.000 05 s，通過對不同時間步長和砰擊壓力曲線進行擬合，得到當時間步長趨于0時，砰擊壓力為pΩ。從圖3可知，隨著時間步長的減小，砰擊壓力大小逐漸增加，最終收斂為一常值。

圖3 砰擊相對壓力與時間步長Fig.3 The curve of relative pressure and timesteps

平底結構在砰擊瞬間入水面中心軸線上壓力分布如圖4所示，其中pMid為中心位置處所受壓力。從圖4可知，計算得到的壓力分布規律同試驗結果一致。

采用提出的計算方法得到的平底結構入水時，結構的砰擊壓力曲線同試驗結果一致，砰擊面壓力分布規律同試驗結果一致，但砰擊壓力值大于試驗結果。說明計算方法有效，但計算結果同試驗相比存在一定的誤差。在網格尺寸為結構長度的1%～4%時，計算結果受網格尺寸影響較小。計算結果隨時間步長減小而收斂為一常值，由于時間步長的減小導致計算成本的增加，因此后文以最小網格為20 mm、時間步長為0.001 s，對平底結構入水砰擊壓力大小、砰擊加速度以及入水空泡等規律進行分析。

圖4 入水端面壓力分布Fig.4 Sketch of pressure on slamming side

2.2 砰擊面砰擊壓力

圖5為入水前后入水端面壓力分布，其中τ為砰擊時刻。從圖5中可以看出，入水砰擊瞬間對結構作用壓力最大，此后會迅速下降。砰擊發生時，砰擊面的壓力分布為中間高、兩端低，隨著結構逐漸浸入水域，砰擊面壓力不斷減小，壓力分布逐漸均勻。

砰擊過程中，砰擊面中間位置受到的砰擊壓力以及壓力變化量都大于結構邊緣位置。砰擊面與流體之間的空氣只能經結構和自由液面間隙向空氣域流動，使得中間位置受到壓力更高，當空氣被完全排出后，砰擊面壓力分布一致。

圖5 入水前后砰擊面壓力分布Fig.5 Sketch of pressure on slamming side during water-entry

2.3 入水過程加速度

對于同一結構，力和加速度成正比，因此以砰擊加速度作為參數來分析結構入水過程受到的砰擊作用力。數值計算獲得了結構在不同高度入水的加速度曲線，如圖6所示。在撞擊液面之前，空氣墊的存在產生向上的加速度，此時加速度曲線略有上升；發生撞擊液面的瞬間，加速度極速增加。由此判斷結構在入水砰擊時刻受到來自液面的砰擊作用力最大，在撞擊水面后迅速下降。在砰擊加速度發生后，可能還存在結構和流體的二次作用，其原因可能是砰擊發生后，原來的自由液面受到擾動，從而產生變化，使得液面和砰擊面出現空氣間隙，隨著結構向下運動，會再次受到流體的砰擊。結合圖2分析，二次砰擊受網格劃分影響，局部網格越小，二次砰擊作用越明顯，但不同的網格使得二次砰擊具有一定的隨機性。

圖6 不同高度入水過程Y方向加速度曲線Fig.6 The curve of Y acceleration with different height in water-entry

砰擊時刻加速度記為a，以峰值壓力的80%、90%作為砰擊持續載荷，其對應的砰擊持續時間分別記為T80、T90。表1給出了不同高度自由下落的結構在入水過程中的參數。從表1中可以看出，下落初速越大，結構受到的砰擊作用力越強，砰擊面受到的壓力載荷也越大。砰擊作用的峰值壓力范圍越小，其持續時間也越短，但整體認為，不同高度下自由入水，砰擊持續的時間約為1 ms。

表1 不同高度入水過程參數

2.4 入水空泡

結構在撞擊水面后，會形成同大氣相通的開空泡，入水空泡的產生、發展等演化過程對入水結構的姿態或者彈道穩定性有影響，尤其是空泡閉合以及破裂時刻。

圖7為下落高度為1.4 m的入水過程液面變化情況，其中初始時刻為物體到達未擾動液面位置時刻。從圖7得出，結構撞擊液面改變了初始的流場參數，使得原來的液面發生分離。水自結構撞擊入水端面處向外排擠，并在物體兩側形成射流，射流初始方向朝上，隨后向外傾斜。在結構浸入水域的過程中，射流由初始的連續流體逐漸變得不穩定，形成分離飛濺的水花，并在重力的作用下下落匯入水中。圖7同文獻[13]中平底體模型砰擊時產生的飛濺水花相似。

圖7 入水過程中氣液界面變化(h=1.4 m)Fig.7 The evolution of gas-liquid interface when flatted-bottom body entering water(h=1.4 m)

在入水過程中，空氣同物體一起被壓入水中，會形成“空氣墊”。由于平底結構的擊水面同初始液面間不存在傾角，因此初始階段結構擊水面的空氣逃逸較少，隨著浸入水深的增加，空氣才不斷由入水端面排出，經結構兩側進入大氣。

從圖7中t=0.05 s到t=0.73 s時刻可以看出在平底結構入水過程中，以入水端面擊水面為界，最初會形成與大氣相通的開空泡，隨著物體浸入水域深度的增加，在流體壓力的作用下，被物體排向兩側的流體向中間流動并發生翻卷，隨后阻斷了開空泡在物體上方形成閉合?？张蓍]合后，物體上方的流體相向流動，發生撞擊，撞擊處的流體一部分向上方空氣域產生射流，一部分流體向結構上表面流動并對其產生沖擊。

圖8為空泡閉合前后的速度矢量圖，從圖8可以看出，隨著水深的增加，流體壓力迫使物體入水端面、兩側以及上方的流場進一步發生變化。其中物體入水端面流體受到擠壓并進一步被排開，結合圖7可以推斷出物體兩側的空泡內部氣體被擠向物體上表面，導致側向空泡變窄，側向的沾濕面積增加。結構上方的流體相向運動，撞擊點位于上方中心位置處，撞擊后流體在沿著結構上表面流動、翻卷，并將空泡分割開。

圖8 閉合空泡沖擊上表面的速度矢量圖Fig.8 The vector of speed of the cavity when water impact on the top of the flatted-bottom body

3 結 論

利用FLUENT平臺中基于VOF方法和動網格技術，對平底鋼結構自由落體入水過程的砰擊載荷特性和氣液兩相流動進行了求解，數值計算得到入水過程運動規律、入水過程的載荷特性、壓力分布規律以及入水空泡演變過程，其具體結論如下：

(1)平底結構自由落體入水過程中，由于液面對結構的砰擊作用，使得撞擊液面時刻產生明顯的砰擊加速度，入水初速越大，物體受到液面的砰擊載荷越大；

(2)平底結構入水過程中，入水端面受到液面的砰擊壓力呈現中間高兩端低的分布特點，隨著入水時間的增加，入水端面空氣被排出，其壓力分布逐漸趨于平緩；

(3)平底結構自由落體入水會形成入水空泡，該空泡經歷開空泡、閉合、破裂等過程?？张蓍]合后，物體上方的流體在壓力作用下會對物體上表面產生沖擊，并在上表面形成翻卷的流體，從而導致空泡被割裂。

[1] WAUGH J G,STUBSTAD G W. 水彈道學模擬[M].陳九錫,張開榮,譯. 北京: 國防工業出版社,1979.

[2] 王永虎,石秀華. 入水沖擊問題研究的現狀與進展[J]. 爆炸與沖擊, 2008,28 (3): 276-282. WANG Yonghu, SHI Xiuhua. Review on research and development of water-entry inpact problem[J]. Explosion and Shock Wave,2008,28(3):276-282.

[3] 陳偉琪,王寶壽,顏開,等. 物體垂直入水的非定?？张輸祵W模型[J]. 應用數學和力學,2013,34(11):1130-1140. CHEN Weiqi,WANG Baoshou, YAN Kai,et al. Model of the unsteady vertical water-entry and water-exit cavities[J]. Applied Mathematics and Mechanics,2013,34(11):1130-1140.

[4] 何春濤,王聰,何乾坤,等. 圓柱體低速入水空泡試驗研究[J]. 物理學報,2012,61(13):1-8. HE Chuntao,WANG Cong, HE Qiankun,et al. Lowspeed water-entry of cylinderical projectile[J]. Acta Physica Sinica,2012, 61(13):1-8.

[5] WEI Zhaoyu, SHI Xiuhua. Evolution of three-dimensional cavitation following water entry of an inclined cylinder[J]. Theoretical & Applied Mechanics Letters,2012,2(2):1-3.

[6] ARISTOFF J M,TRUSCOTT T T,TECHET A H,et al. The water-entry cavity formed by low bond number impacts[J]. Physics of Fluid, 2008,20:091111.

[7] 何春濤,王聰,閔景新,等. 回轉體勻速垂直入水早期空泡數值模擬研究[J]. 工程力學,2012,29(4):237-243. HE Chuntao,WANG Cong,MIN Jingxin,et al. Numerical simulation of early air-cavity of cylinder cone with vertical water-entry[J].Engineering Mechanics,2012, 29(4):237-243.

[8] 陳宇翔,郜冶,劉乾坤. 應用VOF方法的水平圓柱入水數值模擬[J]. 哈爾濱工程大學學報, 2011,32(11):1339-1442. CHEN Yuxiang,GAO Ye, LIU Qiankun. Numerical simulation of water-entry in a horizontal circular cylinder using the volume of fluid(VOF) method[J]. Journal of Harbin Engineering University, 2011,32(11):1339-1442.

[9] 陳學農，何友聲.平頭物體三維帶空泡入水的數值模擬[J].力學學報,1990,22(3):129-138. CHEN Xuenong, HE Yousheng. Numerical simulation of 3D water entry of a blunt cylinder with a ventilated cavity[J].Acta Mechanica Sinica , 1990,22(2):129-138.

[10] HUERA-HUARTE F J, JEON D,GHARIB M. Experimental investigation of water slamming loads on panels[J].Ocean Engineering,2011,38(11/12):1347-1355.

[11] CHEN Zhen, XIAO Xi, WANG Deyu. Simulation of flat-bottom structure slamming[J]. China Ocean Engineering,2005,19(3):385-394.

[12] 楊衡,孫龍泉,龔小超,等.彈性結構入水砰擊載荷特性三維數值模擬研究[J].振動與沖擊,2014,33(19):28-34. YANG Heng, SUN Longquan, GONG Xiaochao,et al. 3D numeriacal simulation of slamming load character for water entry of an elastic structure[J].Journal of Vibration and Shock,2014,33(19):28-34.

[13] 黃震球,張文海.減小平底體砰擊的試驗研究[J]. 華中工學院學報,1986，14(5): 725-730. HUANG Zhenqiu,ZHANG Wenhai.Experimental investigation on the reduction of flat-botton body slamming[J].Journal of Huazhong University of Science and Technology，1986,14(5):725-730.

Numerical simulation on the slamming load in the water-entry process of flatted-bottom body based on the method VOF

WANG Yijun, LI Minghai, ZHANG Zhongli, LI Shangming, SHI Guangmei

(Insititue of Engineering system of CAEP, Mianyang 621900, China)

The volume of fluid(VOF) method coupled with the dynamic mesh method was used to simulate the slamming load on a flatted-bottom body and the liquid-gas multiphase flow during the water-entry process. The unsteady water-entry cavity produced by the flatted-bottom body was analyzed. The flatted-bottom body motion, slamming load, pressure distribution, free surface deformation, and air cushion effect were captured. The evolutionary process of cavity such as its movement and closure was presented. The maximal slamming load occurs at the very beginning and the pressure on the middle part of the bottom surface is higher than that on the edge part. The slamming load increases with the initial velocity of water-entry. Moreover, the water rushes on the surface after the cavity collapse. The numerical calculation method and results are of practical significance for the safety and adaptability assessment of flatted-bottom body during water-entry.

VOF; water-entry; slamming load; cavity; numerical simulation

中國工程物理研究院“雙百人才工程”基金項目(ZX04001)；中國工程物理研究院總體工程研究所科技專項基金項目(2013KJZ01)

2015-06-17 修改稿收到日期：2016-01-12

王易君 女，碩士，工程師，1986年生

李明海 男，研究員，博士生導師，1970年生

U663.2

A

10.13465/j.cnki.jvs.2017.02.030