物體自由入水的多參數(shù)影響分析

劉雙， 武雨嫣， 何廣華， 謝濱陽

(1.哈爾濱工業(yè)大學(xué) 機(jī)電工程學(xué)院，黑龍江 哈爾濱 150001; 2.哈爾濱工業(yè)大學(xué)(威海) 船舶與海洋工程學(xué)院，山東 威海 264209)

入水砰擊問題廣泛存在[1],砰擊時(shí)間歷程短、壓力峰值大，易引起結(jié)構(gòu)失效或破壞，對(duì)結(jié)構(gòu)物局部或整體的強(qiáng)度均會(huì)產(chǎn)生較大危害，因此對(duì)結(jié)構(gòu)物入水砰擊問題進(jìn)行研究具有十分重要的意義。關(guān)于結(jié)構(gòu)物入水問題的研究最早始于19世紀(jì)末。Worthington等[2]對(duì)不同表面粗糙度球體的入水過程以及入水噴濺、入水空泡閉合等現(xiàn)象進(jìn)行了系統(tǒng)的研究。Von Karman[3]對(duì)水上飛機(jī)浮箱降落過程中撞擊水面所受的砰擊力進(jìn)行了理論研究。Wagner[4]在Von Karman砰擊理論基礎(chǔ)之上，提出了經(jīng)典的Wagner勢(shì)流砰擊理論。其余還有Dobrovol skaya[5]、Bihnam等[6]、Zhao等[7]、Hu等[8]。陳翔等[9]采用MPS求解器研究了二維楔形體在不同粒子布置方式下的入水問題。張?jiān)狼嗟萚10]采用LS-DYNA軟件對(duì)楔形體及弧形體入水問題進(jìn)行了模擬。其他有張軍等[11]、孫輝等[12]。王永虎等[13]、李云波等[14]采用勢(shì)流方法研究入水問題，張于維等[15]、王易君等[16]使用Fluent軟件研究了結(jié)構(gòu)物入水過程。眾多研究為了問題的簡化，以研究恒速入水現(xiàn)象居多。本研究基于可自動(dòng)生成網(wǎng)格的STAR-CCM+和常用軟件Fluent對(duì)自由入水砰擊問題展開了研究。首先建立了2套數(shù)值模型以實(shí)現(xiàn)自相互驗(yàn)證，通過網(wǎng)格收斂性及速度-壓力耦合算法適用性的研究并結(jié)合試驗(yàn)值進(jìn)行比對(duì)進(jìn)一步驗(yàn)證了數(shù)值方法的準(zhǔn)確性。接著采用2種模型同時(shí)模擬了不同入水速度和底升角等參數(shù)對(duì)楔形體自由入水過程的影響。

1 理論與數(shù)值模型

1.1 二維楔形體砰擊理論

圖1為楔形體入水砰擊模型。楔形體的底升角為β，單位長度楔形體的質(zhì)量為M，入水時(shí)刻速度為V0，入水深度為z時(shí)，對(duì)應(yīng)的下落速度為V，靜水面處半寬為c0，楔形體被液體浸濕部分半寬為c。忽略楔形體重力及浮力的作用，根據(jù)動(dòng)量定理可以得到：

MV0=(M+m)V

(1)

式中m為楔形體單位長度的附加質(zhì)量。

圖1 楔形體入水砰擊模型Fig.1 Water slamming model of wedge body

將式(1)進(jìn)行微分可得某一瞬間作用在楔形體的沖擊力：

(2)

可見，砰擊載荷與入水速度等因素密切相關(guān)。

1.2 基本控制方程

控制方程為：

(3)

式中：ρ代表流體密度；t為時(shí)間，u、v、w分別表示速度沿x、y、z3個(gè)方向上的分量：

(4)

式中：Fbx、Fby、Fbz分別是質(zhì)量力在3個(gè)方向上的分量；pyx為流體內(nèi)應(yīng)力張量的分量。式(3)、(4)分別為連續(xù)性方程和N-S方程的三維形式，本文的入水問題研究主要關(guān)注x-z平面。

2 數(shù)值模擬與結(jié)果分析

2.1 計(jì)算模型的建立

為檢驗(yàn)本數(shù)值模型的準(zhǔn)確性，選用Zhao等[7]模型試驗(yàn)作為參考，圖2是試驗(yàn)中的楔形體模型。在楔形體右側(cè)分布有P1、P2、P3、P4、P55個(gè)傳感器，用來測定測量點(diǎn)的砰擊壓力隨時(shí)間變化情況，同時(shí)通過測定模型下落時(shí)的加速度，計(jì)算出模型整體的砰擊壓力值。表1給出了試驗(yàn)?zāi)Ｐ偷木唧w參數(shù)。圖3是2種CFD計(jì)算模型所采用的網(wǎng)格。

圖2 Zhao等[7]的楔形體試驗(yàn)?zāi)Ｐ虵ig.2 Wedge test model of Zhao et al[7]

表1 楔形體的主要參數(shù)Table 1 Principle parameters of wedge

為減小邊界對(duì)計(jì)算結(jié)果的影響，本文計(jì)算域?qū)挾热⌒ㄐ误w寬度的10倍，為5 000 mm；高度為5 200 mm；其中空氣域高度2 500 mm；液體高度2 700 mm；楔形體的底升角為30°；寬度500 mm；頂端在自由液面上方100 mm處。計(jì)算域上方邊界為壓力出口，其余邊界為標(biāo)準(zhǔn)壁面，楔形體在重力加速作用下以一定初速度由空氣落入水中。在STAR-CCM+中沒有模擬純二維情況，楔形體取單位厚度。為與試驗(yàn)保持一致，楔形體尖端觸水速度均為6.15 m/s。

圖3 楔形體流場網(wǎng)格劃分Fig.3 Computational mesh for water entry of wedge

本Fluent模型基于網(wǎng)格重構(gòu)技術(shù)，采用VOF及6-DOF方式追蹤自由液面及楔形體運(yùn)動(dòng)，湍流模型采用Realizablek-ω模型；而在STAR-CCM+里，基于重疊網(wǎng)格技術(shù)和DFBI方法模擬大尺度運(yùn)動(dòng)問題；通過壁面函數(shù)法對(duì)近壁區(qū)Realizablek-ε模型進(jìn)行處理；采用歐拉多相流模型對(duì)空氣和水進(jìn)行建模；采用VOF法捕捉自由液面。

2.1.1 收斂性驗(yàn)證

首先，對(duì)2種CFD計(jì)算模型的網(wǎng)格收斂性展開了驗(yàn)證。表2、3分別給出了2種計(jì)算模型中采用的3種網(wǎng)格。Fluent是通過改變最小網(wǎng)格尺寸來控制網(wǎng)格數(shù)量；STAR-CCM+中是通過改變網(wǎng)格基礎(chǔ)尺寸來控制網(wǎng)格數(shù)量。

表2 Fluent中采用的網(wǎng)格Table 2 Number of grid used in Fluent

表3 STAR-CCM+中采用的網(wǎng)格Table 3 Number of grid used in STAR-CCM+

圖4分別給出了Fluent及STAR-CCM+中網(wǎng)格密度對(duì)楔形體入水過程速度衰減的影響，同時(shí)與試驗(yàn)結(jié)果進(jìn)行比對(duì)。

圖4 速度變化曲線圖Fig.4 Speed curve

由速度變化曲線可以看出，F(xiàn)luent及STAR-CCM+中預(yù)報(bào)的速度變化趨勢(shì)和試驗(yàn)結(jié)果整體上均趨勢(shì)一致，證明本文所采用的2種CFD計(jì)算模型具有較好的可靠性。為了保障計(jì)算的精度，在接下來的模擬中，2種方法均選用較為精細(xì)的Mesh C網(wǎng)格。整體來看，STAR-CCM+的結(jié)果與試驗(yàn)更為接近。原因是在Fluent中的二維模型無法考慮到三維效應(yīng)，而STAR-CCM+中考慮了三維效應(yīng)，更接近試驗(yàn)的實(shí)際情況。但同時(shí)STAR-CCM+的網(wǎng)格數(shù)量也會(huì)遠(yuǎn)大于Fluent，致使計(jì)算量偏大。

為了更好地捕捉楔形體入水時(shí)的流場信息，2種計(jì)算模型均選取了較小時(shí)間步長0.000 1 s，均滿足收斂性要求，在此不再贅述。

2.1.2 速度-壓力耦合求解算法的選擇

對(duì)速度和壓力的耦合關(guān)系進(jìn)行求解時(shí)，常用算法有SIMPLE、SIMPLEC、PISO和Coupled。在STAR-CCM+中采用的是默認(rèn)的SIMPLE算法；在Fluent中為了對(duì)比計(jì)算精度，選擇更為穩(wěn)定的算法，分別采用不同算法進(jìn)行模擬得到了砰擊力曲線，并與STAR-CCM+及試驗(yàn)的結(jié)果進(jìn)行對(duì)比，見圖5。

圖5 各方法砰擊壓力結(jié)果比較Fig.5 Comparison of slamming pressure of each method

表4給出了數(shù)值計(jì)算和試驗(yàn)得到的砰擊力峰值。結(jié)合圖5與表4發(fā)現(xiàn)，在Fluent中，SIMPLEC算法和Coupled算法結(jié)果比較相近，但Coupled算法耗時(shí)較長，PISO算法雖然測得砰擊力與試驗(yàn)誤差最小，但計(jì)算結(jié)果較不穩(wěn)定，因此在后續(xù)Fluent計(jì)算中采用SIMPLEC算法。

表4 數(shù)值結(jié)果與試驗(yàn)結(jié)果對(duì)比Table 4 Comparison of numerical and experimental results

整體來說，F(xiàn)luent與STAR-CCM+中的計(jì)算結(jié)果均與試驗(yàn)值吻合較好，均具有較高精度；具體來看，STAR-CCM+計(jì)算結(jié)果要略小于Fluent，對(duì)于砰擊壓力極值的預(yù)報(bào)，STAR-CCM+誤差略大.

2.1.3 測量點(diǎn)砰擊壓力曲線對(duì)比

按照Zhao等[7]的試驗(yàn)，設(shè)定楔形體的質(zhì)量為241 kg，入水時(shí)刻速度為6.15 m/s，選擇入水時(shí)刻為時(shí)間零點(diǎn)。根據(jù)上面2節(jié)的分析，F(xiàn)luent中采用1.25 mm網(wǎng)格及SIMPLEC計(jì)算模型進(jìn)行計(jì)算，在STAR-CCM+中采用基礎(chǔ)尺寸60 mm網(wǎng)格進(jìn)行模擬。對(duì)楔形體P1～P55個(gè)測量點(diǎn)處的砰擊壓強(qiáng)進(jìn)行監(jiān)測，圖6是STAR-CCM+、Fluent計(jì)算結(jié)果以及試驗(yàn)結(jié)果的時(shí)歷曲線圖，僅展示P1、P3、P53個(gè)點(diǎn)。

從圖6可見，2種計(jì)算模型所得各點(diǎn)數(shù)值結(jié)果均可實(shí)現(xiàn)自相互驗(yàn)證，且與試驗(yàn)結(jié)果吻合較好，展現(xiàn)出了相同的規(guī)律：測量點(diǎn)入水后所受砰擊壓力由0迅速增加到最大值，之后緩慢降低；隨著各點(diǎn)入水深度的增加，砰擊壓力極值先增加后減少。相比來看，STAR-CCM+中各點(diǎn)砰擊壓力的極值結(jié)果略小于Fluent，但與試驗(yàn)值更為接近。圖7為楔形體入水過程中2種計(jì)算模型的自由液面云圖比較。

圖6 各點(diǎn)砰擊壓力比較Fig.6 Comparison of every point slamming pressure

圖7 入水自由液面云圖對(duì)比Fig.7 Comparison of free surface cloud image

從圖7可見，兩者均能清楚地觀察到液面上升以及射流時(shí)的液面分離現(xiàn)象。故為了簡便，后續(xù)對(duì)于流場云圖的分析將主要展示Fluent監(jiān)測的結(jié)果。

綜上所述，F(xiàn)luent與STAR-CCM+對(duì)于楔形體入水問題的模擬結(jié)果較為接近，并與試驗(yàn)值吻合較好。由于STAR-CCM+中楔形體相比Fluent中純二維情況多了一層厚度，更為接近真實(shí)情況。故STAR-CCM+結(jié)果與試驗(yàn)更為吻合，但同時(shí)其網(wǎng)格數(shù)較多，計(jì)算耗時(shí)長，需要根據(jù)實(shí)際情況進(jìn)行選擇。

2.2 速度對(duì)入水砰擊的影響

2.2.1 入水速度對(duì)砰擊壓力的影響

設(shè)定楔形體的底升角為20°，質(zhì)量為250 kg，分別以4、6、8、10 m/s速度入水。取楔形體尖端觸水時(shí)刻為時(shí)間零點(diǎn)，監(jiān)測楔形體從時(shí)間零點(diǎn)到完全入水過程中的各個(gè)參數(shù)變化情況。

圖8給出了楔形體以4種不同速度入水時(shí)的砰擊壓力時(shí)歷曲線。可得：STAR-CCM+與Fluent計(jì)算結(jié)果整體吻合較好，僅在極值點(diǎn)附近略有差別；STAR-CCM+得到的壓力結(jié)果略低于Fluent結(jié)果。不同速度入水的楔形體砰擊壓力曲線呈現(xiàn)相同的趨勢(shì)：在入水初期，隨著入水深度的增加，砰擊壓力逐漸增加，在發(fā)生液面分離時(shí)砰擊力達(dá)到最大值，之后緩慢降低。

隨著入水速度增加，達(dá)到完全入水所需時(shí)間變短，砰擊壓力極值逐漸變大，出現(xiàn)砰擊壓力極值的時(shí)間也逐漸提前。入水速度為10 m/s的工況，其砰擊壓力極值遠(yuǎn)大于其他3種入水速度，這是由于入水初期，速度大的楔形體在表面形成的射流區(qū)域較大，沾濕面積也較大，產(chǎn)生的砰擊壓力極值也隨之增大。

在楔形體右側(cè)均勻布置P1～P44個(gè)壓力測量點(diǎn)；P1、P2、P3、P4分別距離楔形體尖端53.2 mm、106.4 mm、159.6 mm、212.8 mm。監(jiān)測楔形體在不同入水速度時(shí)砰擊壓力的變化情況，圖9為STAR-CCM+與Fluent中的計(jì)算結(jié)果。在同一測點(diǎn)由于各速度下出現(xiàn)砰擊力峰值的時(shí)間接近；為了易觀察，采用曲線整體平移的方式：將4 m/s曲線整體向右平移0.015 s，6 m/s曲線向右平移0.01 s，8 m/s曲線向右平移0.005 s，10 m/s曲線保持不動(dòng)。

圖8 砰擊壓力計(jì)算結(jié)果Fig.8 Calculation results of sniping pressure

圖9 各點(diǎn)砰擊壓力比較Fig.9 Comparison of sniper pressure at each point

由圖9可以看出，兩者展現(xiàn)出了相同的變化趨勢(shì)，STAR-CCM+模擬得到的砰擊壓力極值略小于Fluent數(shù)值結(jié)果。隨著入水速度增加，各監(jiān)測點(diǎn)的砰擊壓力極值也隨之增加；在同一入水速度下，P2、P3測量點(diǎn)處的砰擊壓力極值總是大于P1和P4測量點(diǎn)。說明楔形體入水時(shí)的最大砰擊壓力點(diǎn)出現(xiàn)在兩側(cè)中部區(qū)域，故在進(jìn)行實(shí)際操作時(shí)需要對(duì)楔形體兩側(cè)中間區(qū)域進(jìn)行材料加強(qiáng)處理。

2.2.2 入水速度對(duì)楔形體運(yùn)動(dòng)特性的影響

以楔形體尖端觸水時(shí)刻為時(shí)間零點(diǎn)，楔形體完全入水時(shí)結(jié)束，圖10給出了楔形體速度時(shí)歷曲線。

由圖10可見，2種數(shù)值模型計(jì)算結(jié)果吻合較好，展現(xiàn)了共性的變化規(guī)律：入水初期沾濕面積小，故阻力小，速度緩慢衰減；隨著入水深度增加，阻力作用面積增大，砰擊壓力逐漸增長導(dǎo)致速度衰減加快；入水后期砰擊壓力減小又導(dǎo)致阻力減小，速度的衰減再次放緩，曲線出現(xiàn)拐點(diǎn)。STAR-CCM+的速度減少相對(duì)Fluent更加緩慢。楔形體入水速度越大，衰減速度越快，衰減幅度越大。結(jié)合圖8可知，在入水初期，高速入水的楔形體會(huì)在較短的時(shí)間內(nèi)產(chǎn)生巨大的砰擊壓力，且砰擊壓力極值遠(yuǎn)大于其他低速情況。因此入水速度越大，速度衰減幅度越大、速率越快，且由于楔形體入水速度較高時(shí)，砰擊壓力在峰值過后往往率先開始衰減，故導(dǎo)致阻力率先減小，相應(yīng)的，曲線出現(xiàn)拐點(diǎn)的時(shí)間會(huì)提前。

圖11為底升角20°楔形體以6 m/s入水時(shí)的相圖。可以發(fā)現(xiàn)隨著入水深度的增加，楔形體射流區(qū)域越來越明顯，沾濕面積越來越大；到楔形體完全入水時(shí)液面飛濺形成水花。

圖10 速度變化趨勢(shì)比較Fig.10 Comparison of speed trends

圖11 底升角20°楔形體以6 m/s入水相圖Fig.11 Phase diagram of 6 m/s with β=20°

2.3 楔形體底升角對(duì)入水砰擊的影響

在進(jìn)行船型設(shè)計(jì)研究時(shí)，需要確定船體剖面的載荷分布，故研究船體底升角變化對(duì)砰擊壓力的影響具有非常重要的意義。研究對(duì)象為底升角20°、30°和40°的楔形體，設(shè)定楔形體的質(zhì)量為250 kg，保證入水時(shí)刻速度為10 m/s，監(jiān)測楔形體由尖端觸水到完全入水過程中的各參數(shù)變化情況。

2.3.1 底升角對(duì)砰擊壓力的影響

圖12展了底升角為20°、30°和40°的楔形體受砰擊力隨時(shí)間變化曲線。從圖12可知：不同底升角的砰擊壓力變化呈現(xiàn)相同的趨勢(shì)。在入水初期，隨著入水深度的增加，砰擊壓力逐漸增加；在發(fā)生液面分離以后，砰擊力迅速降低，在液面分離時(shí)砰擊力達(dá)到最大值。底升角越小，砰擊壓力極值越大，出現(xiàn)時(shí)間越早，這是由于小底升角楔形體在入水初期會(huì)產(chǎn)生較大的射流區(qū)域，使壓力迅速增加。與前文一致的是，STAR-CCM+中的壓力模擬結(jié)果略低于Fluent中模擬結(jié)果。

圖13是不同底升角的楔形體在浸沒時(shí)的相圖。

圖12 不同底升角入水砰擊壓力的比較Fig.12 Comparison of the water pressure among three different bottom angle

圖13 不同底升角的楔形體入水相圖Fig.13 Phase diagram of wedges at different bottom angles

從圖13可以看出，底升角越小的楔形體，沾濕長度越長，產(chǎn)生的射流越明顯。由于底升角小導(dǎo)致液面升高時(shí)液體所受擠壓增加，在大氣壓力和重力等因素的共同作用下，砰擊壓力極值會(huì)增加。

為研究不同底升角時(shí)砰擊壓力極值分布規(guī)律，取3種底升角的楔形體在相同下落時(shí)刻(t=0.004 s)的砰擊壓力分布云圖進(jìn)行對(duì)比。如圖14所示。可以發(fā)現(xiàn)，底升角越小，砰擊壓力極值的產(chǎn)生范圍越小，這是由于射流在楔形體表面形成的角度變小而造成的。圖15給出了底升角為20°的楔形體在入水t=0.004 s時(shí)刻的流場速度矢量圖，可見氣體在楔形體的邊緣會(huì)形成漩渦，最大速度產(chǎn)生在噴射區(qū)，壓力峰值出現(xiàn)在楔形體與自由液面接觸區(qū)域。

圖14 楔形體壓力云圖Fig.14 Wedge pressure cloud diagram

圖15 底升角20°楔形體流場Fig.15 Flow field at 20° bottom angle

表5為t=0.004 s時(shí)刻下不同底升角楔形體的壓力峰值。分析表5，在相同時(shí)刻下，不同底升角的楔形體壓力峰值隨底升角的增加而減小。

表5 在t=0.004 s時(shí)刻不同底升角的楔形體壓力峰值Table 5 Peak pressure with different bottom angles at t=0.004 s

2.3.2 底升角對(duì)速度變化的影響

圖16為速度變化曲線；可見，STAR-CCM+與Fluent趨勢(shì)一致：入水初期，阻力小，速度減少緩慢；隨著入水深度增加，砰擊力加大，加速度增大，速度減少加快；之后砰擊力又變小，速度衰減再次放緩. 比較各底升角的情況：底升角越小，入水時(shí)間越短，速度衰減越快。由于STAR-CCM+模擬結(jié)果中砰擊壓力略小，故速度衰減比Fluent緩慢。

3 結(jié)論

1)STARCCM+與Fluent對(duì)于楔形體自由入水模擬可形成自相互驗(yàn)證，并均與試驗(yàn)結(jié)果吻合較好。STAR-CCM+模擬的楔形體速度降低速率及所受砰擊壓力值均略低于Fluent結(jié)果。整體上STAR-CCM+結(jié)果更加準(zhǔn)確，但計(jì)算資源消耗偏大。

2)楔形體入水速度越大，速度衰減越快，砰擊壓力極值越大，砰擊壓力極值出現(xiàn)時(shí)間越早。楔形體入水時(shí)砰擊力逐漸增加，在發(fā)生液面分離時(shí)達(dá)到最大；極值出現(xiàn)在楔形體兩側(cè)1/2高度附近。

3)底升角越小，入水時(shí)噴射區(qū)域根部液體擠壓越嚴(yán)重，產(chǎn)生的砰擊力極值越大，速度衰減越快。

本研究可為后續(xù)采用STAR-CCM+和Fluent進(jìn)行更為復(fù)雜三維結(jié)構(gòu)入水砰擊載荷預(yù)報(bào)提供參考。