陳宇翔，郜冶，劉乾坤

(哈爾濱工程大學(xué)航天與建筑工程學(xué)院，黑龍江哈爾濱150001)

物體入水現(xiàn)象如船舶航行砰擊、水上飛機(jī)降落、空投魚(yú)雷入水及垮塌橋梁落水廣泛存在于自然界.物體入水問(wèn)題是典型的非線性問(wèn)題，其入水過(guò)程特別是初始瞬間物體會(huì)受到巨大的沖擊載荷，以風(fēng)浪中航行的船舶為例，船體砰擊會(huì)造成船體的震動(dòng)、減緩船速，嚴(yán)重的話還會(huì)導(dǎo)致結(jié)構(gòu)損壞;空投入水的魚(yú)雷則可能因?yàn)槿胨疀_擊造成魚(yú)雷結(jié)構(gòu)損壞、彈道失控等嚴(yán)重問(wèn)題.物體入水問(wèn)題在現(xiàn)代許多軍事及民用領(lǐng)域的廣泛存在及其重要性激發(fā)了眾多科學(xué)家及工程人員的研究熱情［1］.因?yàn)閳A柱、楔等形狀結(jié)構(gòu)簡(jiǎn)單，?？勺鳛樗柩芯课矬w的簡(jiǎn)化外形，所以研究者們常常以它們作為研究對(duì)象.20世紀(jì)90年代早期以前，大部分的研究者［2-4］均基于勢(shì)流理論假設(shè)研究過(guò)楔等形狀物體的入水問(wèn)題.然而，邊界元方法在處理復(fù)雜自由表面問(wèn)題方面不夠強(qiáng)?。?］.自1981年Hirt與 Nichols［6］等提出 VOF方法以來(lái)，帶自由表面的兩相粘性不可壓縮流體運(yùn)動(dòng)的數(shù)值計(jì)算技術(shù)得到了迅速的發(fā)展.近30年來(lái)，VOF方法在數(shù)值計(jì)算物體入水領(lǐng)域得到了廣泛的應(yīng)用.Arai［7］假設(shè)流動(dòng)為無(wú)粘不可壓流動(dòng)，應(yīng)用VOF方法對(duì)二維楔、圓柱、船艏的砰擊入水問(wèn)題進(jìn)行了模擬;Schumann［8］用類似的方法模擬了艏的砰擊入水.Xing-Kaeding［9］基于VOF方法，應(yīng)用Comet軟件對(duì)二維水平圓柱的粘性不可壓出入水運(yùn)動(dòng)進(jìn)行了模擬，取得了與實(shí)驗(yàn)結(jié)果符合較好的結(jié)果.本文應(yīng)用Fluent軟件中的VOF方法結(jié)合動(dòng)網(wǎng)格技術(shù)對(duì)二維圓柱的入水問(wèn)題進(jìn)行了數(shù)值模擬計(jì)算，實(shí)現(xiàn)了圓柱運(yùn)動(dòng)的剛體動(dòng)力學(xué)方程和氣液兩相流動(dòng)的N-S方程的耦合求解.本文數(shù)值模擬計(jì)算結(jié)果與Greenhow和Lin［10］的實(shí)驗(yàn)對(duì)比表明，圓柱入水的重要現(xiàn)象如飛濺現(xiàn)象、空氣墊效應(yīng)［11］等自由表面形狀的變化得到了很好的呈現(xiàn)，圓柱運(yùn)動(dòng)軌跡得到了很好的追蹤.同時(shí)，還對(duì)粘性對(duì)圓柱入水運(yùn)動(dòng)問(wèn)題的影響進(jìn)行了研究.

1 問(wèn)題的描述及數(shù)值水池的建立

如同文獻(xiàn)［10］中Greenhow和Lin設(shè)計(jì)的實(shí)驗(yàn)，一個(gè)直徑11 cm的零浮力(neutral buoyancy)水平圓柱由靜水上方某一確定位置自由落體進(jìn)入水池中.而且當(dāng)圓柱在接觸水面的瞬間，其速度為3.13 m/s，隨后圓柱繼續(xù)向下運(yùn)動(dòng)直至到達(dá)水深為30 cm處.

由以上描述的問(wèn)題，本文建立了一個(gè)寬60 cm，高60 cm(水深30 cm)的二維數(shù)值水池.為準(zhǔn)確地捕捉圓柱入水沖擊引起的自由表面變化等重要現(xiàn)象，整個(gè)計(jì)算區(qū)域采用正交四邊形網(wǎng)格，為此在圓柱邊界采用階梯近似［12］處理，并對(duì)圓柱附近的網(wǎng)格進(jìn)行了加密處理，見(jiàn)圖1.圓柱邊界應(yīng)用無(wú)滑移壁面邊界條件，水池上邊界應(yīng)用壓力出口邊界.為了節(jié)省計(jì)算網(wǎng)格并盡可能模擬真實(shí)流動(dòng)區(qū)域，采用文獻(xiàn)［9］中的方法，水池的左右邊界(包括靜水和空氣區(qū)域)應(yīng)用UDF定義邊界條件使壓力保持為自由表面未變化前的靜壓;水池下邊界應(yīng)用速度為零的速度入口邊界條件.

圖1 圓柱入水問(wèn)題計(jì)算網(wǎng)格Fig.1 Numerical grid for the circular cylinder waterentry case

2 理論基礎(chǔ)和計(jì)算方法

圓柱入水問(wèn)題是復(fù)雜的兩相流體流動(dòng)和固體運(yùn)動(dòng)的耦合過(guò)程.圓柱的入水沖擊會(huì)突然改變兩相流場(chǎng)中水的壓力場(chǎng)和速度場(chǎng)，而流場(chǎng)的改變同時(shí)又會(huì)反作用來(lái)影響圓柱的運(yùn)動(dòng).本文中，有限體積法被用于計(jì)算存在自由表面的粘性不可壓兩相流體流動(dòng);二維圓柱視為剛體，其在落水過(guò)程中為一自由度(1DOF)運(yùn)動(dòng)，流場(chǎng)形狀隨時(shí)間變化，需要應(yīng)用動(dòng)網(wǎng)格.

2.1 流體流動(dòng)與剛體運(yùn)動(dòng)方程

式中:V是流體運(yùn)動(dòng)速度矢量，Vb是有限體積即網(wǎng)格單元表面的運(yùn)動(dòng)速度.動(dòng)量方程為

對(duì)于每一個(gè)有限體積單元，連續(xù)性方程為

式中:T為粘性應(yīng)力張量，I為單位張量，b為體積力矢量.當(dāng)網(wǎng)格單元移動(dòng)或改變形狀時(shí)，有限體積單元滿足空間守恒律:

由于自由表面的形狀和位置是未知的，需要由求解過(guò)程給出，因此數(shù)值求解帶自由表面的兩相流動(dòng)問(wèn)題變得十分困難.為了模擬有自由表面的兩相流流動(dòng)，本文應(yīng)用VOF(volume of fluid)方法.VOF是一種捕捉自由表面的方法，它同時(shí)對(duì)水和空氣建立控制域，通過(guò)求解附加的方程來(lái)捕捉自由面，這種方法在處理復(fù)雜的自由表面流動(dòng)問(wèn)題具有很大優(yōu)勢(shì).引入的附加體積分?jǐn)?shù)輸運(yùn)方程:

式中:c=Vl/V為液相體積分?jǐn)?shù).VOF方法將0＜c＜1的兩相混合流體處理為單相流，認(rèn)為其有同樣的壓力和速度值.其物性參數(shù)如密度處理為

由牛頓第二定律，流體中運(yùn)動(dòng)的圓柱剛體動(dòng)力學(xué)方程為

式中:VC為運(yùn)動(dòng)時(shí)的圓柱重心的速度矢量.

2.2 網(wǎng)格運(yùn)動(dòng)方案

為了適應(yīng)由于圓柱的剛體運(yùn)動(dòng)而不斷改變的流場(chǎng)形狀，需要應(yīng)用動(dòng)網(wǎng)格模型.為了獲得很好的數(shù)值解，特別是在對(duì)于圓柱入水沖擊引起的的強(qiáng)非線性自由表面變化模擬的情況下，需要很好的網(wǎng)格質(zhì)量，因此本文采用整體網(wǎng)格移動(dòng)而不變形的方案［13］.本文應(yīng)用Fluent的動(dòng)態(tài)層模型，由定義整個(gè)流動(dòng)區(qū)域和圓柱為剛體區(qū)域，而將數(shù)值水池上下邊界定義為靜止區(qū)域來(lái)實(shí)現(xiàn)網(wǎng)格整體運(yùn)動(dòng)而不變形.每一當(dāng)前時(shí)間步的網(wǎng)格運(yùn)動(dòng)速度為上一時(shí)間步網(wǎng)格速度與當(dāng)前時(shí)間步圓柱所受的兩相流體壓力、粘性力以及自身重力所引起的速度增量之和.

2.3 計(jì)算方法

本文中，PISO算法用來(lái)耦合速度場(chǎng)和壓力場(chǎng);VOF兩相流模型的自由面重構(gòu)方法用來(lái)捕捉自由表面的變化情況;對(duì)于流場(chǎng)中的湍流的描述則采用標(biāo)準(zhǔn)k-ε兩方程模型和壁面函數(shù)法;壓強(qiáng)差值采用PRESTO!差分格式.

3 計(jì)算結(jié)果及分析

本文應(yīng)用之前描述的Fluent動(dòng)網(wǎng)格技術(shù)結(jié)合VOF方法對(duì)二維零浮力水平圓柱入水問(wèn)題進(jìn)行了計(jì)算，對(duì)比分析了數(shù)值計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果，并且對(duì)于粘性對(duì)圓柱入水運(yùn)動(dòng)的影響進(jìn)行了分析.

3.1 自由表面的變化

圓柱入水的瞬間，改變了其周?chē)鷫毫Φ攘鲌?chǎng)參數(shù)，圓柱的速度突然下降，而同時(shí)與其接觸的水相應(yīng)地被突然加速，在圓柱兩側(cè)激起兩股射流，射流向斜上方上升直至其變得不穩(wěn)定.由于重力的影響，隨著時(shí)間的發(fā)展，射流形狀由直線逐漸彎曲向下.

圖2 零浮力圓柱入水自由表面變化:實(shí)驗(yàn)(左)和數(shù)值(右)結(jié)果對(duì)比Fig.2 Free surface deformation of a neutrally buoyant cylinder water-entry:simulation(right)and experiment(left)

在圓柱入水的過(guò)程中，一些空氣也隨同圓柱一起被壓入水中，會(huì)形成所謂的＂空氣墊＂.由于圓柱底部?jī)A角較大，導(dǎo)致底部的空氣逃逸較多，可壓縮空氣墊主要存在于圓柱兩側(cè)，將圓柱兩側(cè)表面與水隔開(kāi)，隨著圓柱入水深度的不斷增加，空氣不斷排出，圓柱的沾濕面積不斷增加.

為了方便比較觀察實(shí)驗(yàn)照片和數(shù)值計(jì)算結(jié)果，圖2給出了5個(gè)不同時(shí)刻實(shí)驗(yàn)和數(shù)值模擬的水池自由表面變化情況對(duì)比圖.由圖2可知，本文計(jì)算的零浮力圓柱入水過(guò)程自由表面變化情況和實(shí)驗(yàn)對(duì)比符合很好，成功地模擬顯示了圓柱入水沖擊引起的水射流形成、行進(jìn)、彎曲以及入水圓柱周?chē)目諝鈮|等自由表面的變化現(xiàn)象.

3.2 圓柱的剛體運(yùn)動(dòng)

進(jìn)入水池之前，圓柱在空氣中自由下落，在氣流場(chǎng)中受到的阻力較小.圓柱接觸水面的瞬間，兩相流場(chǎng)中的水由于圓柱的沖擊被大大加速，壓力場(chǎng)、速度場(chǎng)特別是圓柱附近的壓力分布突然改變.相應(yīng)地，圓柱受到流場(chǎng)很大的阻力(包括壓差阻力、摩擦阻力)，速度大大降低.圖3給出了圓柱入水過(guò)程的無(wú)量綱阻力(流場(chǎng)數(shù)據(jù)積分得到的入水阻力與重力的比值)與時(shí)間的關(guān)系，零時(shí)刻圓柱受到很大的一個(gè)沖擊力，而隨著時(shí)間的進(jìn)行，進(jìn)入水中的圓柱周?chē)膲毫Ψ植贾饾u均勻，圓柱速度越來(lái)越慢，阻力的變化逐漸平緩.

圖3 圓柱入水過(guò)程的無(wú)量綱阻力時(shí)間曲線Fig.3 Time histories of non-dimensional drag force during cylinder water-entry

圖4 圓柱豎直運(yùn)動(dòng)時(shí)間曲線Fig.4 Time histories of vertical motion of the cylinder

圖4給出了圓柱入水過(guò)程沒(méi)入水中的圓柱位置(y軸坐標(biāo)為入水后圓柱底部最低點(diǎn)的位置)和入水時(shí)間關(guān)系曲線，隨著圓柱在水中的下降，圓柱受到的阻力不斷減小，速度也逐漸降低，因此圓柱位置的變化曲線也逐漸平緩.由圖4中的曲線與實(shí)驗(yàn)記錄的數(shù)據(jù)點(diǎn)對(duì)比可知，數(shù)值解同樣和Greenhow和Lin的實(shí)驗(yàn)相當(dāng)吻合(圖中偏差較大、標(biāo)記問(wèn)號(hào)的一個(gè)點(diǎn)同時(shí)也是實(shí)驗(yàn)中認(rèn)為存在問(wèn)題的數(shù)據(jù)點(diǎn)).

3.3 粘性對(duì)入水圓柱運(yùn)動(dòng)的影響

由于水的粘性，入水過(guò)程中圓柱除了受到由于兩相流場(chǎng)的壓力場(chǎng)分布不均帶來(lái)的壓差阻力，還受到粘性引起的摩擦阻力.

圖5 圓柱入水豎直運(yùn)動(dòng)過(guò)程時(shí)間曲線(湍流和無(wú)粘流)Fig.5 Time histories of vertical motion of the cylinder during cylinder water-entry(turbulence and without turbulence model)

圖6 圓柱入水過(guò)程的無(wú)量綱阻力(湍流和無(wú)粘流對(duì)比)Fig.6 Time histories of non-dimensional drag force during cylinder water-entry(turbulence and without turbulence model)

為了考察粘性對(duì)圓柱入水問(wèn)題的影響，本文計(jì)算了無(wú)粘流模型下的圓柱入水過(guò)程，并與之前基于湍流模型的數(shù)值計(jì)算結(jié)果進(jìn)行比較.圖5、6分別給出了基于湍流模型和無(wú)粘流模型計(jì)算的圓柱入水過(guò)程圓柱運(yùn)動(dòng)軌跡和圓柱所受阻力時(shí)間曲線.可以看出，2種模型的數(shù)值模擬結(jié)果基本一致，所以入水過(guò)程中圓柱主要受氣液兩相流場(chǎng)壓力分布不均引起的壓差阻力，而其所受的摩擦阻力對(duì)圓柱入水過(guò)程受力和圓柱運(yùn)動(dòng)軌跡影響較小.

4 結(jié)束語(yǔ)

本文以Fluent軟件為計(jì)算工具，應(yīng)用VOF結(jié)合動(dòng)網(wǎng)格技術(shù)的數(shù)值方法對(duì)水平圓柱入水過(guò)程氣液兩相流動(dòng)和圓柱剛體運(yùn)動(dòng)問(wèn)題進(jìn)行了模擬，數(shù)值模擬的自由表面變化過(guò)程以及圓柱的剛體運(yùn)動(dòng)軌跡與經(jīng)典文獻(xiàn)的實(shí)驗(yàn)符合得很好.從中可以看出，本文應(yīng)用的數(shù)值方法可以用來(lái)模擬物體入水特別是伴隨著強(qiáng)非線性自由表面變化的問(wèn)題，而之前常應(yīng)用的邊界元方法不能精確模擬此類問(wèn)題;本文的數(shù)值方法可以應(yīng)用來(lái)計(jì)算物體入水所受的沖擊力，進(jìn)而對(duì)物體如船體結(jié)構(gòu)的設(shè)計(jì)提出指導(dǎo)意見(jiàn)和提供參考數(shù)據(jù).

［1］王永虎，石秀華.入水沖擊問(wèn)題研究的現(xiàn)狀與進(jìn)展［J］.爆炸與沖擊，2008，28(3):276-282.

WANG Yonghu，SHI Xiuhua.Review on research and development of water-entry impact problem［J］.Explosion and Shock Waves，2008，28(3):276-282.

［2］Von KARMAN T.The impact on seaplane floats during landing.NACA TN 321［R］.Washington:National Advisory Committee for Aeronautics，1929.

［3］WAGNER H.Uber stoss-und gleitvorgange an der oberflache von flussigkeiten［J］.Zeitschrift Für Angewandte Mathematik und Mechanik，1932，12(4):193-215.

［4］ZHAO R，F(xiàn)ALTINSEN O M.Water-entry of two-dimensional bodies［J］.Journal of Fluid Mechanics，1993，246:593-612.

［5］ZHU Xinying，F(xiàn)ALTINSEN O M，HU Changhong.Water entry and exit of a horizontal circular cylinder［J］.Journal of Offshore Mechanics and Arctic Engineering，2007，129:253-264.

［6］HIRT C W，NICHOLS B D.Volume of fluid(VOF)method for dynamics of free boundaries［J］.Journal of Computational Physics，1981，39:201-225.

［7］ARAI M，CHENG L Y，INOUE Y.A computing method for the analysis of water impact of arbitrary shaped bodies［J］.Journal of the Society of Naval Architects of Japan，1995，176:233-240.

［8］SCHUMANN C.Volume-of-fluid computations of water entry of bow sections［C］//Proceedings of Euromech 374.Poitiers，1998:209-217.

［9］XING-KAEDING Y，JENSEN G，PERIC M.Numerical simulation of water-entry and water-exit of a horizontal circular cylinder［C］//Proceedings of the 6th International Conference on Hydrodynamics.Perth，2004:663-669.

［10］GREENHOW M，LIN W M.Nonlinear free-surface effects:Experiments and theory.MIT report 83-19［R］.Cambridge:Massachusetts Institute of Technology，1983.

［11］陳震，肖熙.空氣墊在平底結(jié)構(gòu)入水砰擊中作用的仿真分析［J］.上海交通大學(xué)學(xué)報(bào)，2005，39(5):670-673.

CHEN Zhen，XIAO Xi.Simulation analysis on the role of air cushion in the slamming of a flat-bottom structure［J］.Journal of Shanghai Jiao Tong University，2005，39(5):670-673.

［12］XING T，STERN F.Grid generation and post-processing for computational fluid dynamics(CFD)［Z］.［s.l.］:IIHR-Hydroscience and Engineering，2006.

［13］HAD?ICA I，HENNIGB J，PERICM，et al.Computation of flow-induced motion of floating bodies［J］.Applied Mathematical Modelling，2005，29(12):1196-1210.