二維平板水漂運動數(shù)值模擬

付曉琴，李陽輝，盧昱錦，肖天航，童明波

南京航空航天大學 航空宇航學院，南京 210016

長期以來，對于水漂運動的機理研究及應用已深入到諸多領域。在炮彈領域，為了能夠避開防魚雷網而輕松摧毀目標，參照二戰(zhàn)時期航空工程師Barnes Neville Wallis設計的“跳彈”，Johnson[1]提出，通過分析炮彈撞擊角度和彈射速度范圍，能夠實現(xiàn)不斷在水上彈跳進而躲避防御的方案。航天方面，宇宙飛船的返回艙再入大氣層時，通過控制其與大氣層的接觸角度，產生“彈跳”的運動過程，能達到降溫的效果[2]。同樣，在水陸兩棲飛機領域，為了保證飛機著水安全性與穩(wěn)定性[3-4]，該問題顯得尤為重要。較差的觸水姿態(tài)可能會導致飛機在著水時發(fā)生“海豚運動”或者跳躍現(xiàn)象[5-6]。因此，對水漂運動的研究具有一定的工程應用價值。

面對水漂這一物理現(xiàn)象，Rosellini[7]、Lyderic[8]以及Clanet[9]等開展了大量的試驗研究，借助高速攝像機記錄圓盤水漂運動過程不同時刻的姿態(tài)與水面變化，探究其中的流體機理與力學性能，結果表明當圓盤與水面的夾角大約為20°時，圓盤打水漂所需的初始投擲速度最小，并且由于圓盤的高度自旋使得姿態(tài)角變化不明顯。此外，Rosellini等[7]還通過簡單的單體理論預測石頭的運動，分析石頭在打水漂過程中受到的重力以及沖擊過程中石頭表面的壓力，證明在高速旋轉的情況下，石頭的反彈是由水動力響應引起，在運動過程中的速度衰減和相對能量損失與自身姿態(tài)角有關。理論分析結果雖然與試驗結果相差不大，但在分析過程中忽略了沖擊時水體的運動，在壓力計算中面積多是采用與流體接觸的石頭面積近似值。Lyderic等[10]則提出了石頭與水碰撞過程的簡化物理模型，依然不考慮由于石頭運動導致的流體運動，結果表明石頭完成水漂的最大次數(shù)與自身速度衰減和穩(wěn)定性有關。通過增大石頭的初始速度可以實現(xiàn)多次的反彈。另外，保持石頭的高度自旋有利于自身的穩(wěn)定性。Do等[11]利用一種基于物理建模的水漂運動仿真方法分別計算水對石頭的作用力和水的變形，得到石頭水漂的運動軌跡。球體能在水面反彈存在一個最大的入射角，如果撞擊速度足夠大，速度與入射角之間的關系可以用一個簡單的常微分方程模型(Ordinary Differential Equation, ODE)[12-13]來表達。Nagahiro等[14]將ODE模型應用到圓盤水漂沖擊中，求得圓盤能彈跳的最小速度和最大投擲角度的解析解。戴巖偉[15]同樣借助水漂運動的物理模型，計算觸水過程中石頭受到的浮力和阻力，并從能量消耗的角度分析石頭能彈跳的最小速度和最大投擲角度。理論分析通常采用簡化模型，得到的結果與實際情況有一定差異。近年來，更多的探索工作集中在水漂運動的數(shù)值模擬研究。Nagahiro等[14]使用光滑粒子流體動力學(Smoothed Particle Hydrodynamics, SPH)方法對Clanet等提出的最佳攻角、最小投擲速度和最大投擲角度進行數(shù)值仿真，計算結果與試驗值相吻合。Yan和Monaghan[16]采用相同的方法對二維平板水漂運動過程進行數(shù)值模擬，論證了三維運動與二維運動的相似性，討論不同姿態(tài)角、投擲角和底面形狀的平板對水漂運動的影響。鄔明[17]和陳詩偉[18]利用任意拉格朗日-歐拉(Arbitrary Lagrange-Euler, ALE)法對圓盤擊水彈跳現(xiàn)象展開研究，分析入水角、俯仰角、初始速度和旋轉角速度的相關性。面對水漂運動可能發(fā)生的多次彈跳，SPH和ALE方法需要構建相當大的計算域，使得計算成本過高，效率低。

本文應用體積分數(shù) (VOF) 模型結合整體動網格方法對二維平板水漂運動進行數(shù)值模擬，加入了有限體積法的優(yōu)勢，規(guī)避了多次彈跳過程導致的數(shù)值水池區(qū)域過大的問題，在實現(xiàn)精確捕捉水面的同時提高了計算效率與精度。確定能完成水漂運動的最小投擲速度和最大投擲角度，并從能量損失的角度分析不同姿態(tài)角、投擲角、投擲速度對水漂運動的影響，從動量變化的角度探究發(fā)生“彈跳”的原因。進而，引入速度入口造波法研究波浪水面水漂運動過程，掌握平板運動姿態(tài)和能量變化在同一波形不同位置和不同波高工況下的運動規(guī)律，為水陸兩棲飛機著水或類水漂運動的研究工作提供技術支持和數(shù)據(jù)支撐。

1 數(shù)值計算方法

對于平板水漂運動的數(shù)值模擬計算方法主要包括以下4個部分：流體控制方程及求解，自由液面捕捉方法，速度入口造波方法，整體動網格方法。流體控制方程主要為數(shù)值模擬提供理論支持，自由液面采用VOF方法進行精確捕捉，波浪數(shù)值水池的建立依賴速度入口造波法，整體動網格方法最終實現(xiàn)平板運動過程的模擬。

1.1 流體控制方程求解

流體控制方程為非定常雷諾平均Navier-Stokes(URANS)方程，湍流模型選用k-εRNG。SIMPLE算法用于求解壓力-速度耦合計算；有限體積法實現(xiàn)控制方程對空間的離散；在時間離散中，采用格林-高斯格點方法求解梯度，選擇PRESTO！離散壓力項；利用幾何重構方法描述體積分數(shù)；其余均選用二階迎風格式離散。

1.2 自由液面捕捉方法

在研究剛體與水面撞擊的過程中，對于水面變化和液體飛濺的精確捕捉有利于還原實際的物理過程，同時提高仿真結果的準確性。VOF方法可以實時追蹤自由液面的形狀和位置，利用網格單元中流體體積和網格單元自身體積的比值函數(shù)來確定自由液面的狀態(tài)[19]。比值函數(shù)范圍為0～1。數(shù)值為0則表示此網格不存在該流體，數(shù)值為1說明此網格被該流體充滿，數(shù)值介于0～1為自由液面。

1.3 速度入口造波方法

速度入口造波方法是將波面方程(式(1))和水質點速度(式(2))施加在速度入口邊界。隨著計算時間的推進，波浪從入口向計算域遠場傳播開來。為了避免波浪傳播到計算域出口產生反射波，對平板黏性流場和運動的計算產生影響，需要將壓力出口邊界附近的一段區(qū)域設置為阻尼消波區(qū)。文獻[20-21]中采用速度入口造波法數(shù)值模擬了線性規(guī)則波與不規(guī)則波，與試驗結果吻合。該方法的造波精度論證工作在前期已經完成[22-23]，這里不做贅述。本文采用線性規(guī)則波實現(xiàn)平板在波浪水面水漂運動的數(shù)值模擬。波浪參數(shù)及波形位置說明如圖1所示，波面方程表達式為

(1)

式中：η為波面高度；H為波高；x為水質點的x坐標；k為波數(shù)；ω為圓頻率。

波浪中水質點(x,y)的水平速度分量u與垂直速度分量v的表達式為

(2)

式中：φ為速度勢函數(shù)；T為周期；d為水深；θx為波浪傳播方向與x軸正向的夾角。

阻尼消波采用動量源項消波的方法，對于距離壓力出口位置兩個波長的長度區(qū)域進行消波，消波方式為在動量方程中增加動量衰減項[22]，其表達式為

(3)

式中：x0為在消波區(qū)中的計算域x坐標。消波區(qū)示意如圖1所示，可以看到，添加動量衰減項之后，波幅明顯變小。

圖1 波浪參數(shù)與波浪不同位置

1.4 整體動網格方法

整體動網格方法是一種在物體運動過程中將整個計算域(包括網格單元和邊界)與物體相固連的方法，即整個計算域與物體有相同的運動規(guī)律[24]。在兩相流計算中，由于網格區(qū)域跟隨物體一起運動，為了保持水面位置不發(fā)生變化。需要對邊界條件進行特殊處理。在靜水面的計算中，水面保持固定高度，需給定合理的體積分數(shù)邊界條件，將水面之上的部分賦予空氣的屬性，水面之下的部分賦予水的屬性。在波浪水面的計算中，保證波浪平衡位置為一定值，速度入口邊界的波面高度隨時間做余弦函數(shù)變化，以維持穩(wěn)定的波形。在體積分數(shù)邊界條件的設置中，水的體積分數(shù)α1∈[0，1]。通過判定當前網格與水面位置的關系，給出合理的體積分數(shù)。α1=0表示網格中水體積分數(shù)為0，此網格在水面上方；α1=1表示網格中水體積分數(shù)為1，此網格在水面下方；α1=0.5表示網格中水體積分數(shù)為0.5，此網格為水氣交界面。通過對網格的分類，實現(xiàn)不同介質的壓力初始化。位于水面以上為空氣域，壓力即為大氣壓101 325 Pa，位于水面以下為水域，壓力與深度成正比，P=ρghw，其中，g=9.8 kg/m2，ρ=998.2kg/m3，hw為當前網格與水面的垂直距離。此方法在課題組之前的工作中已經取得良好的成果[22-23]。

2 數(shù)值模型與網格設置

2.1 平板水漂模型

參照Rosellini等[7]的水漂試驗圓板尺寸構建數(shù)值模擬的二維平板模型。采用的二維鋁制板如圖2所示，其中，α為平板姿態(tài)角，β為平板投擲角，R為平板半長度，h為平板厚度,U為投擲速度。該平板長度L為2R=50 mm，厚度為h=2.75 mm。

圖2 二維平板模型

2.2 網格劃分

為保證遠場不受影響，計算域長為16L，高為10L。采用結構網格劃分流場計算域，如圖3(a)所示。為了更精確地捕捉水氣兩相流的情況，同時考慮到平板在運動過程中垂向位移較大，需要對平板與水面附近區(qū)域的網格進行局部加密，垂向加密區(qū)寬度為2L(-0.07～0.03 m)。根據(jù)k-ε湍流模型適用的y+值的范圍為30～300，得到最小網格尺寸為0.000 2 m。取網格增長率為1.2，實現(xiàn)加密區(qū)的網格向外圍網格的光滑過渡，得到最大網格尺寸為0.006 m，網格量約為13萬。靜水面工況計算域邊界條件的設置如圖3(b)所示，平板表面為無滑移壁面邊界。在波浪水面的數(shù)值模擬中，由于波浪傳播到計算域出口會產生回流，為了避免對流場造成影響，計算域長為30L，高為10L，網格劃分與上文一致，邊界條件如圖3(c)所示。

圖3 二維網格與邊界條件

3 算例驗證

為了驗證數(shù)值方法及整體動網格技術在平板水漂數(shù)值模擬過程中的可行性和準確性，采用文獻[7]中的平板水漂試驗進行算例驗證。計算域x、y向尺寸與二維網格一致，z向長度為10L；采用全結構網格劃分，最小網格尺寸為0.000 2 m，網格增長率為1.2，得到的總網格量約為470萬，如圖4(a)所示。邊界條件設置與坐標軸的定義如圖4(b)所示，平板初始最低位置與水面重合，通過數(shù)值模擬得到三維圓板在不同姿態(tài)角下能完成水漂運動所需的最小速度與最大投擲角。在最小速度的計算過程中，平板投擲角為20°；在數(shù)值模擬最大投擲角時，平板投擲速度為3.5 m/s。

圖4 三維網格與邊界條件

圖5(a)給出平板完成水漂運動所需最小投擲速度Umin的數(shù)值模擬結果與試驗值[6]和理論解[14]的對比。可以看出，數(shù)值模擬值與試驗值相比偏小，但一直都在理論值A和B之間，三者具有相同的變化規(guī)律。具體來看，當a<5°時，數(shù)值模擬得到的最小速度比理論值A、B和試驗值低，隨著α的增大，數(shù)值模擬得到的最小速度與試驗值逐漸接近，直到α>20°時，數(shù)值模擬值與理論值A和試驗值吻合度越來越高。同時，通過數(shù)值模擬可以看出，當平板姿態(tài)角為20°左右時，其能進行水漂運動所需的投擲速度最小，這與試驗得到的結論相吻合。

圖5(b)給出平板能進行水漂運動的最大投擲角度βmax的數(shù)值模擬結果與試驗值[7]和理論解[14]的對比?？梢钥闯?，數(shù)值模擬值與試驗值和理論值具有相同的變化規(guī)律，雖然相比試驗值稍大，但一直在理論值A、B形成的區(qū)域內。具體表現(xiàn)為，當姿態(tài)角較小時，數(shù)值模擬得到的值與理論值B吻合度較高，與試驗值相比偏大；隨著姿態(tài)角的增大，開始與理論值A和試驗值相吻合。

圖5 數(shù)值模擬與試驗值和理論值對比

出現(xiàn)上述現(xiàn)象的原因是：當平板以較小的姿態(tài)角觸水時，與水面的夾角較小，平板下方的空氣容易被壓入水中，導致明顯的空氣墊[25]現(xiàn)象?？諝馄扑樾纬蓺馀萑菀自斐伤娴淖冃?，它的狀態(tài)與流動形式不穩(wěn)定，其中水氣比例不斷變化，很難采用合理的理論模型與數(shù)值模擬方法進行準確的描述，因此導致數(shù)值模擬得到的結果在小姿態(tài)角的條件下與試驗值和理論值具有一定差距。當平板姿態(tài)角增大時，在平板下方的空氣更容易逃逸，空氣墊現(xiàn)象不明顯，此時數(shù)值模擬得到的值與試驗值和理論值A誤差更小?？紤]到網格精度有限，數(shù)值模擬得到的值與試驗值和理論解仍然有所差別。但總體來看，數(shù)值模擬的結果與試驗值和理論解有相同的變化規(guī)律，且誤差在可接受范圍內。

4 二維與三維平板的運動特性相似性

本文主要對二維平板進行水漂運動數(shù)值模擬，在外形上與實際三維圓板有所區(qū)別，故需要對二維與三維平板的運動特性相似性進行說明。

平板與水面接觸碰撞期間，流體作用力是決定平板能否反彈的重要因素，根據(jù)平板垂向動力學方程：

(4)

(5)

式中：Fz為垂向流體作用力；Swet為平板觸水面積；CL為升力系數(shù)；CD為阻力系數(shù)，令Cn=CLcosα-CDsinα，其只與模型形狀和α相關，由于水漂運動選用的α較小，故Cn為常數(shù)。

分3種情況對Swet進行討論：三維圓板、三維方板、二維方板；其中，圓板直徑、三維方板底部邊長以及二維方板邊長均為2R，三者高度均為h。

Swet=

(6)

由式(6)可知，當模型為三維和二維方板時，Swet與z呈線性關系；當模型為三維圓板時，Swet與z呈非線性關系。如圖6所示，為R∈[15 mm,30 mm]，α=20°的三維圓板Swet與z的關系曲線簇。

圖6 α=20°時Swet與z的關系曲線

利用多項式線性擬合可以得到Swet與z的近似線性關系式：

Swet≈-(1.779 7R)kz-0.174 6R2-

0.916 4R+7.107

(7)

對原始曲線和擬合曲線在R∈[15 mm,30 mm]，α∈[10°，30°]范圍內進行擬合優(yōu)度計算，得到的決定系數(shù)r2>0.995，決定系數(shù)r2[26-28]越接近1代表兩條曲線吻合程度越高；故式(7)能很好地代替式(6)中的三維圓板關系曲線。

以R=25 mm，α=20°為例，Swet的理論曲線與通過式(6)得到的擬合曲線的決定系數(shù)r2=0.995 7， 曲線如圖7所示，在初始觸水與觸水到最大深度時，擬合曲線與理論曲線有一定誤差；在其他位置，該線性曲線可以很好地描述圓板實際觸水的運動特征。

圖7 Swet與z的理論與擬合曲線

當模型為圓板時，Swet與z的關系式為

Swet≈-0.889 85(L/sinα)z+φ

(8)

式中：L=2R；φ=-0.174 6R2-0.916 4R+7.107。

綜上，

(9)

故式(4)可以改寫為

(10)

式中：q為常數(shù)項，

其中：Cn1、Cn2、Cn3為常數(shù)，由物體形狀決定。

因此，式(10)為受定值外力作用的簡諧運動方程式，可以改寫為

(11)

(12)

其中：

代入式(12)，得到

(13)

圖8 觸水時間τ隨的變化曲線

5 計算結果及分析

5.1 二維平板靜水面水漂數(shù)值模擬

參照文獻[16]，由于平板在水漂過程中高速自旋，姿態(tài)角α在整個碰撞過程中保持不變，因此在數(shù)值模擬過程中只放開x方向與y方向平移自由度。本節(jié)通過改變平板初始姿態(tài)角α、投擲角β、投擲速度U對平板水漂運動過程進行數(shù)值模擬。計算直到平板沉入水中。模型如圖2所示。通過分析平板水平速度、垂直速度、垂直位移、能量變化，確定上述參數(shù)對平板水漂運動過程的影響。能量僅考慮平板機械能，其定義為

(14)

式中：E為物體機械能；m為物體質量；U為運動速度大??；g=9.8 kg/m2；H為高度差。

5.1.1 姿態(tài)角α的影響

平板投擲角為20°，投擲速度為4 m/s，初始時刻平板最低點與水面重合。數(shù)值模擬α=10°，20°，30°時平板的水漂運動。圖9為不同姿態(tài)角的平板第1次水漂運動過程中不同時刻的水面變化，可以看出，隨著姿態(tài)角的增大，平板觸水時形成的水面凹陷越大，其前端形成的射流越大，完成一次水漂的耗時更長。

圖9 不同姿態(tài)角條件不同時刻水面變化

不同姿態(tài)角的平板在水漂運動過程中速度變化曲線以及平板重心位置的垂直位移和能量變化曲線見圖10。如圖10(a)所示，平板在觸水的過程中，水平速度持續(xù)衰減。隨著姿態(tài)角的增大，水平速度衰減幅度增大。此外，圖10(b)垂直速度變化曲線與圖10(c)垂直位移變化曲線顯示：在第1次水漂運動的過程中，姿態(tài)角越大，平板垂直速度變化范圍越大，反彈的高度也就越高。而后續(xù)過程，垂直速度變化范圍逐漸縮小，圖10(d)給出平板水漂運動過程的能量變化曲線，表1反映了不同姿態(tài)角的平板在前3次觸水過程中的相對能量損失:

表1 不同姿態(tài)角平板前3次水漂運動相對能量損失

圖10 不同姿態(tài)角條件對比曲線

顯然，姿態(tài)角越大的平板，一次水漂運動過程的相對能量損失越大。當姿態(tài)角從30°減小到10°，在第1次水漂運動的過程中，相對能量損失從55.56%降低到25.56%。同時，姿態(tài)角為10°、20°的平板在后續(xù)的運動過程中，相對能量損失呈現(xiàn)遞減趨勢，而姿態(tài)角為30°的平板相對能量損失基本不變，導致30°姿態(tài)角的平板水漂次數(shù)要小于其余兩個姿態(tài)角的平板。分析原因：姿態(tài)角小的平板在觸水后，與水接觸面積迅速增大，導致升力增大，平板在入水過程中減速的更快，平板浸水深度小，在水面向前運動所需推動的流體量小，所受到的阻力小，能量的損失小。平板在水漂運動過程中能否從水面彈起取決于在外載荷作用下的y方向動量改變量，即y方向外載荷F對平板的沖量。根據(jù)沖量的定義：

(15)

不同姿態(tài)角下，平板質量m相同，入水時y向速度Uy1相同，故出水時y向速度Uy2只與沖量相關。圖11給出平板在水漂運動過程中受到的沖量變化，可以看到，在第1次水漂運動時，10°姿態(tài)角的平板所受到的沖量最小，僅為30°姿態(tài)角的78%左右，故平板出水速度小，從水面反彈的高度也相應變小。而姿態(tài)角越大，平板完成一次觸水過程的所受沖量大，因此出水速度大，平板能彈起更高的高度；但平板初始觸水時，與水面小的接觸面積導致的小升力使得平板更容易浸入水中，產生更大的附加質量，能量損失有所增大。當姿態(tài)角過大時，能量損失過大引起水平速度快速衰減，垂直速度變化范圍迅速變小，平板無法離開水面，很快沉沒。綜上所述，對比3個不同的姿態(tài)角，可以看出以相對能量損失最小的10°姿態(tài)角進行水漂運動可以得到最多的水漂次數(shù)。

圖11 y方向沖量變化曲線

5.1.2 投擲角β的影響

平板姿態(tài)角固定為20°，投擲速度為4 m/s，初始時刻平板最低處與水面重合，對比分析投擲角β=10°，20°，30°時的仿真結果，總結出投擲角對平板水漂的影響規(guī)律。圖12為不同投擲角的平板第1次水漂運動過程中不同時刻的水面變化情況?？梢钥闯?，采用相同的平板，投擲角較小時，平板觸水深度小，在運動過程中激起的水花也小，平板很快從水面反彈。而當投擲角較大時則相反，更容易在平板前部形成射流，激起的水花更多，完成一次水漂的耗時也相對較長。

圖12 不同投擲角條件不同時刻水面變化

圖13給出不同投擲角的平板在水漂運動過程中速度變化的曲線以及垂直位移和能量變化曲線。如圖13(a)所示，投擲角越大，平板在水漂運動過程中水平速度衰減幅度越大。圖13(b)垂直速度變化曲線與圖13(c)垂直位移變化曲線表明，更大投擲角度使平板擁有更大的垂直速度，導致平板浸沒深度更大。這意味著在水漂過程中平板會產生更大的附加質量和能量損失。在完成第一次水漂后，20°投擲角的平板反而是彈起高度最高的。反映在能量損失方面，如圖13(d)和表2所示，投擲角越大的平板在第1次水漂過程入水深度也越大，由于附加質量的原因，平板消耗在水中的能量也越大。當平板投擲角β從10°增加到30°，在完成第1次水漂后，相對能量損失從24.59%上升到57.75%。正因如此，投擲角β=10°的平板雖然垂直方向速度最小，但得益于其每次相對能量損失為3種工況中最小，因此實現(xiàn)的水漂次數(shù)最多。

圖13 不同投擲角條件對比曲線

表2 不同投擲角平板前3次水漂運動相對能量損失

5.1.3 投擲速度U的影響

平板姿態(tài)角為20°，投擲角為20°，初始時刻平板最低點與水面重合，分析對比U=3，4，6 m/s時的仿真結果，得到投擲速度對平板水漂的影響。圖14分別不同投擲速度的平板第1次水漂運動過程中不同時刻的水面變化。投擲速度小的平板在水漂運動過程中水面形成的凹陷小，在運動過程中所激起的水花與平板前端的射流不明顯；當投擲速度變大時，水面形成的凹陷變大，在運動過程中所激起的水花與平板前端的射流變得明顯。

圖14 不同投擲速度條件不同時刻水面變化

不同投擲速度條件下平板在水漂運動過程中速度變化曲線以及垂直位移和能量變化曲線如圖15所示。圖15(a)水平速度變化曲線顯示平板投擲速度越大，初始水平速度越大，衰減的程度趨于一致，因此投擲速度小的平板更快停下來。圖15(b)垂直速度變化曲線與圖15(c)垂直位移變化曲線表明，平板投擲速度越大，初始垂直速度越大，導致平板在水面反彈的高度高，完成一次水漂運動的時間長。從圖15(d)平板水漂運動過程的能量變化曲線以及表3中平板前3次水漂運動的相對能量損失數(shù)值可以看出，在第1次水漂運動的過程中，投擲速度大的平板的相對能量損失略小于投擲速度小的平板。具體而言，平板速度大小從速度為3 m/s增加到6 m/s時，相對能量損失從41.27% 下降到38.75%。在后續(xù)水漂運動的過程中，相對能量損失的比例也都維持在30%左右。分析原因：由于平板的初始速度有差異，投擲速度大的平板初始水平速度與垂向速度大，因此在水漂運動過程中，平板的垂向速度與垂向位移變化范圍大。而根據(jù)伯努利定理，水平速度大的平板在觸水時帶動附近的水運動，與下層靜水面的壓差大，造成平板在觸水時受到的升力大，平板在垂向的加速度大，減速更快。但正是由于此時初始垂直速度也大，因此投擲速度對觸水深度影響不大，相對能量損失隨其變化規(guī)律也不明顯。且對比前兩個工況，平板初始能量相等，改變投擲速度時，平板的初始能量存在差異。投擲速度越大初始能量較多，而每次水漂運動的相對能量損失比例相當，因此投擲速度小的平板在水面水漂運動的次數(shù)更少。綜上所述，投擲速度大的平板能完成更多次數(shù)的水漂運動。

圖15 不同投擲速度條件對比曲線

表3 不同投擲速度平板前3次水漂運動相對能量損失

5.2 二維平板波浪水面水漂數(shù)值模擬

平板在波浪水面水漂運動時，從水面彈起后無法預測下一次觸水點的位置，因此本節(jié)重點研究平板在波浪水面的第1次水漂過程，討論不同相位位置和不同波高對平板水漂運動的影響。

5.2.1 同一波形不同相位觸水對平板水漂運動的影響

在波浪條件下，由于平板姿態(tài)角為20°，計算得到平板在垂向投影長度為0.017 m，為了防止波高太高，平板直接沖入水中，無法離開水面，因此選擇線性規(guī)則波浪參數(shù)為：波高H=0.02 m，波長l=0.3 m。利用速度入口造波法得到的波浪初始水面如圖16所示。初始時刻模型最低點與波浪相應位置重合。計算過程中平板姿態(tài)角α=20°，投擲角β=20°，投擲速度大小U=4 m/s，波浪從左向右傳播，平板迎浪投擲。選擇波峰、波谷、平衡位置-上行速度和下行速度4個觸水位置(見圖1) 數(shù)值模擬平板第1次水漂運動過程。

圖16 波浪初始相位

圖17為平板在波浪水面的不同位置第1次水漂運動的水面變化情況。可以看出，平板在波峰觸水與水面撞擊形成的凹陷最小，在波谷處最大。在波谷位置和波浪平衡位置-上行速度位置觸水的平板，需要跟隨波浪運動到波峰才能離開水面。

圖17 平板在波浪不同位置觸水時水面變化

圖18給出平板在波浪不同位置觸水時的速度、垂直位移以及平板能量變化曲線。可以看出，平板在波浪水面水漂運動與在靜水面水漂運動有相似的運動趨勢。在第1次水漂運動的過程中，不同觸水位置水平速度衰減由大到小依次為波谷、波浪平衡位置-上行速度、靜水面、波浪平衡位置-下行速度、波峰位置。其中，最大衰減幅度約為25.5%，最小約為 14%。對于垂直速度的變化，卻呈現(xiàn)出與水平速度近乎相反的規(guī)律，在波谷處觸水的平板最大垂直速度可達到1.2 m/s，在波峰處觸水的平板僅約為0.2 m/s。因此導致不同位置觸水的平板從水面反彈的高度差異巨大。圖18(d)以及表4從能量角度解析平板運動特征的差異。平板在波浪水面水漂時，除波峰位置觸水的相對能量損失較小，其余條件下相對能量損失大致相似，在40%左右，與5.1節(jié)靜水面條件下水漂運動相對能量損失相似。具體來看，在波浪平衡位置-上行速度觸水的平板相對能量損失最大，約占40.95%，波峰處最小，僅為最大值的86.8%左右。分析原因：平板所有的能量損失可以理解為水漂過程中水相獲得的能量。平板與水相的能量交換與平板和撞擊點波浪運動切線的夾角θ以及接觸時間有關，其中，夾角θ的定義如圖19所示，夾角越大，平板觸水的面積越大，受到的阻力越大，導致水平方向速度衰減幅度大，造成更多的能量損失；對于上述4個觸水點位置，在波浪平衡位置-上行速度觸水，此時平板迎浪撞擊波面，與波浪切線的夾角接近180°，因此導致其相對能量損失大；當平板在波峰處觸水時，與波浪接觸的時間最短，同時接觸面積也較小，造成平板的能量損失??；在波浪的波谷位置觸水時，雖然此時夾角不是最大，由于波浪和平板存在相對運動，平板必須貼著波面運動到波峰位置，與水面的接觸時間最長，導致平板的能量損失較大。最終呈現(xiàn)出上述運動姿態(tài)與能量變化規(guī)律。

圖18 波浪不同位置觸水對比曲線

表4 平板在不同位置觸水的相對能量損失

圖19 平板和波浪在撞擊點切線的夾角

5.2.2 不同波高同一觸水位置對平板水漂運動的影響

參照5.2.1節(jié)得到的結論，選擇能量損失最大的波浪平衡位置-上行速度與能量損失最小的波峰位置觸水，通過改變波高為0.01 m、0.02 m、0.03 m，實現(xiàn)在不同波高對于平板水漂運動的影響。

1) 波浪平衡位置-上行速度觸水

圖20給出不同波高條件下在波浪平衡位置-上行速度觸水的水面變化圖。不同波高明顯改變了水與平板的接觸面積，且平板均隨波浪運動到波峰位置才能離開水面。

圖20 平衡位置-上行速度觸水水面變化

平板在不同波高的波浪平衡位置-上行速度觸水時的速度、垂直位移和能量變化曲線見圖21?？梢钥闯?，平板運動姿態(tài)和能量的變化與波高線性相關，但影響不大。在波高為0.01 m的波形觸水，平板水平速度衰減最小，約為23.6%，垂直速度變化量最小，約為2.17 m/s，相對應的垂直位移也最小，此時的相對能量損失約為39.1%。隨著波高的遞增，水平速度衰減幅度增量約為1.99%，垂直速度變化增量在0.07 m/s左右，相對能量損失增量約為1.35%。分析原因：不同波高參數(shù)改變了平板與接觸點波浪運動切線的夾角。當平板在波浪平衡位置-上行速度觸水時，波高越高，意味著平板與波浪切線的夾角越大，此時平板觸水時與水面貼合的更好，接觸面積更大，因此平板水平速度衰減幅度更大，由于此時受到的垂向力更大，平板能反彈更高的高度。同時，平板在此處觸水后需要運動到波峰的位置才能離開水面，因此在波高越高的工況下平板完成一次水漂運動的觸水時間長，造成的相對能量損失也就相應更大。

圖21 平衡位置-上行速度觸水變化對比曲線

2) 波峰位置觸水

圖22給出不同波高條件下在波浪波峰位置觸水的水面變化圖。可以看出，不同波高條件下平板水漂運動時水面形成的凹陷大致相同。平板從觸水到離水的時間隨著波高的增加而變小。

圖22 波峰位置觸水水面變化

圖23為平板在波高為0.01 m、0.02 m、0.03 m的波浪波峰觸水時的速度、垂直位移和能量變化曲線?？梢钥闯?，其變化趨勢與波浪平衡位置-上行速度位置觸水相反，且受波高影響更明顯。其中，在波高為0.01 m的波形觸水，平板水平速度衰減最大，約為17.5%，垂直速度變化量最大，約為1.928 m/s，相對應的垂直位移也最大；此時相對能量損失約為37.8%。隨著波高的遞增，水平速度衰減幅度減小約4%，垂直速度變化減小值在0.25 m/s左右，相對能量損失減小值約為2.36%。分析原因：波高參數(shù)越大，波浪越陡，波峰形狀也更尖銳。因此，在波高越高的波峰位置觸水時，平板越容易離開水面，由于接觸的時間短，平板水平速度衰減幅度小，能量損失也就相應的變小。同時由于越尖銳的波峰位置處波浪內流速度變化越小，平板在觸水過程中受到的沖量小，導致垂直速度變化小，垂直位移相應變小。

圖23 波峰位置觸水變化對比曲線

6 結 論

1) 采用整體動網格結合VOF方法可以對平板水漂運動過程進行很好的數(shù)值模擬，計算得到的結果與試驗值和理論值的誤差均在合理范圍內。

2) 當姿態(tài)角為20°左右時，平板能以最小的投擲速度完成水漂運動。當投擲速度一定時，姿態(tài)角為10°的平板的相對能量損失小，隨著投擲角的增大相對能量損失也隨之增大，因此10°姿態(tài)角的平板水漂運動的次數(shù)最多；平板投擲角與投擲速度也影響著水漂運動的次數(shù)，投擲角為30°的平板相對能量損失最大，隨著投擲角的減小，平板的相對能量損失減少，水漂運動的次數(shù)增多；投擲速度則影響平板的初始能量。投擲速度越大，平板初始能量越大，雖然相對能量損失比例相當，投擲速度大的平板可以完成更多次數(shù)的水漂運動。

3) 研究平板在波浪水面不同位置觸水時，在波浪平衡位置-上行速度觸水的平板相對能量損失最大，約為40.95%，其次為在波谷觸水，上述兩種條件均比在相同條件下的靜水面工況的相對能量損失大。在波浪平衡位置-下行速度位置觸水的相對能量損失較小，最小發(fā)生在波峰位置觸水，僅為最大值的86.8%左右。這是由于不同的觸水位置決定著平板與水面的接觸面積和觸水時間，接觸面積越大，觸水時間越長，相對能量損失越大。

4) 探究波高對平板水漂運動的影響，平板運動姿態(tài)與能量損失與波高線性相關。波高對平板在波浪平衡位置-上行速度觸水影響較小，當波高為0.01 m時，平板水平速度衰減最小，約為23.6%，垂直速度變化量最小，約為2.17 m/s，相對應的垂直位移最小，相對能量損失約為39.1%。隨著波高的遞增，水平速度衰減幅度增量約為1.99%。垂直速度變化增量在0.07 m/s左右，相對能量損失增量約為1.35%。波高對平板在波峰位置觸水影響稍大，且與在波浪平衡位置-上行速度觸水工況的規(guī)律相反。當波高為0.01 m時，平板水平速度衰減最大，約為17.5%，垂直速度變化量最大，約為1.928 m/s，相對應的垂直位移也最大，相對能量損失約為37.8%。隨著波高的遞增，水平速度衰減幅度減小約4%。垂直速度變化減小值在0.25 m/s左右，相對能量損失減小值約為2.36%。