基于SPH方法的M型快艇入水載荷仿真分析

楊一凡，陳三桂，周維星，張濤*

1 華中科技大學(xué) 船舶與海洋工程學(xué)院，湖北 武漢 430074

2 中國艦船研究設(shè)計(jì)中心，湖北 武漢 430064

0 引 言

隨著船舶向高速、輕量化方向發(fā)展，具有吃水淺、速度快和靈活性等多項(xiàng)優(yōu)點(diǎn)的M型快艇越來越受到國際造船界的重視。M型快艇由主船體、輔助船體和連接橋3部分組成。其中，M型快艇的中間區(qū)域?yàn)橹鞔w，主要作用是通過排水為船體提供浮力；兩側(cè)的圍壁為輔助船體，主要起到密封作用；主船體與輔助船體之間則通過連接橋連接。M型快艇航行時(shí)，主船體、輔助船體、連接橋與水面圍成的區(qū)域中的空氣受到壓縮而形成空氣氣墊，其氣墊效應(yīng)使船體抬升并快速且平穩(wěn)地航行。同時(shí)，在M型快艇高速航行中，主船體與輔助船體之間的空間內(nèi)會(huì)反復(fù)進(jìn)出水，船底也會(huì)受到波浪反復(fù)的砰擊，嚴(yán)重影響M型快艇的安全性。因此，研究M型快艇入水過程中船體受到的砰擊載荷具有重要意義，其空氣墊效應(yīng)也是不可忽視的。

Zhang等[1]通過拖船試驗(yàn)研究了M型快艇在不同波浪作用下的砰擊載荷，但總體上國內(nèi)外學(xué)者對(duì)M型快艇入水問題的研究仍較少。相反，對(duì)于與M型快艇十分相似的三體船結(jié)構(gòu)入水問題的研究比較多。例如，Davis等[2]最早開始采用實(shí)驗(yàn)和附加質(zhì)量理論研究三體船入水問題，還開展了雙體船不同截面結(jié)構(gòu)的入水實(shí)驗(yàn)，運(yùn)用數(shù)值方法分析入水過程，并針對(duì)結(jié)構(gòu)入水時(shí)雙體船縫隙中殘留有空氣的問題，采用軟連接方式將水與結(jié)構(gòu)物的附加質(zhì)量聯(lián)系起來。此外，曹正林等[3]對(duì)高速三體船入水時(shí)影響連接橋砰擊壓力峰值的因素也進(jìn)行了相關(guān)研究。

對(duì)于結(jié)構(gòu)物入水問題的研究，最初基本上采用的都是試驗(yàn)方法。例如，Lin等[4]通過平板入水試驗(yàn)，發(fā)現(xiàn)平板在砰擊水面時(shí)存在空氣層，并比較了不同高度下平板自由落體入水產(chǎn)生的加速度和壓力。在數(shù)值計(jì)算方面，最具代表性的是源于馮·卡曼相關(guān)研究的楔形物體入水問題[5]。然而，除了研究著陸過程中的水對(duì)水上飛機(jī)的影響外，并未考慮空氣對(duì)其的影響。Zhao等[6]基于非線性邊界元理論，研究了具有任意形狀的二維結(jié)構(gòu)的入水問題，在計(jì)算模型中檢查了質(zhì)量、動(dòng)量和能量的守恒性。之后，有學(xué)者在此基礎(chǔ)上也進(jìn)行了廣泛和深入的研究[7-9]。不僅如此，針對(duì)結(jié)構(gòu)物入水問題的研究，還發(fā)展出了不同的理論。然而，有學(xué)者發(fā)現(xiàn)無論是使用邊界元法（BEM）還是有限元法（FEM），在模擬結(jié)構(gòu)物砰擊水面時(shí)，評(píng)估自由表面的非線性演化過程仍然具有挑戰(zhàn)性。

近年來，一種新的數(shù)值計(jì)算方法——光滑粒子流體動(dòng)力學(xué)（smooth particle dynamics，SPH）算法逐漸受到學(xué)者們的關(guān)注。SPH算法是一種無網(wǎng)格方法，因沒有大變形引起的網(wǎng)格變形，便于跟蹤材料的軌跡，界面處理也靈活，所以SPH算法非常適合處理結(jié)構(gòu)物入水過程中的液面大變形問題。近年來，該方法成為了研究結(jié)構(gòu)物入水問題的重要方法。例如，閆蕊等[10]基于SPH算法研究了空氣對(duì)平板入水中的影響， Oger等[11]提出一種在SPH算法模擬過程中提取近固體邊界局部壓力的新方法，并模擬了楔形物體入水過程，所得結(jié)果與實(shí)驗(yàn)數(shù)據(jù)吻合較好。此外，Yan等[12]使用SPH單相流模型模擬了三體船的部分入水情況。

基于現(xiàn)有研究，本文擬采用SPH液?氣兩相流算法，對(duì)M型快艇典型截面結(jié)構(gòu)入水過程中的砰擊載荷進(jìn)行研究。通過模擬平板和弓形結(jié)構(gòu)的入水過程，將仿真結(jié)果與文獻(xiàn)試驗(yàn)結(jié)果進(jìn)行對(duì)比，以驗(yàn)證基于SPH液?氣兩相流理論計(jì)算結(jié)構(gòu)物入水結(jié)果的準(zhǔn)確性。此外，通過SPH兩相流算法模擬不同斜升角和初始入水速度下M型快艇典型截面結(jié)構(gòu)的入水過程，分析砰擊載荷的變化，討論斜升角對(duì)船舶入水過程的影響。

1 SPH液?氣兩相流理論

SPH方法是一種無網(wǎng)格的拉格朗日算法，其開源的DualSPHysics代碼是一種硬件加速的“光滑粒子流體動(dòng)力學(xué)”代碼，已被證實(shí)可用于模擬當(dāng)前結(jié)構(gòu)，解決自由表面流動(dòng)問題[13]。

1.1 光滑核函數(shù)

采用SPH算法建立的模型其性能很大程度上取決于光滑核函數(shù)的選擇，此函數(shù)在SPH算法計(jì)算中的作用至關(guān)重要。此外，偏微分方程所有解都需轉(zhuǎn)化為核函數(shù)，這在一定程度上也決定了計(jì)算的精確性。

目前，核函數(shù)種類繁多，為了提高計(jì)算效率，新的核函數(shù)不斷涌現(xiàn)。而常用的則是Wendland核函數(shù)[14]，具體形式如下：

式中：W(r,h）為Wendland核函數(shù)，其中r和h分別表示任意粒子a和b的距離以及光滑長度；在二維空間，αD= 7/(4πh2)，在三維空間，αD=21/(16πh3)；q為r與h的比值。

1.2 動(dòng)量方程

SPH算法中的動(dòng)量方程用于移動(dòng)粒子，而人工黏度方案[15]則是SPH算法中進(jìn)行流體仿真的常用方法，故可以寫為：

式中：Pa和ρa(bǔ)分別為粒子a的壓力和密度；va為粒子a的速度；Pb和ρb分別為粒子b的壓力和密度；mb為粒子b的質(zhì)量；?a為粒子a的哈密爾頓算子；Wab為核函數(shù)；g為重力加速度；Πab為黏度項(xiàng)，并由式（3）給出。

式中：α為引入適當(dāng)?shù)南⒍枰{(diào)整的系數(shù)；為 平均聲速，=0.5(ca+cb)， 其中ca，cb分別為粒子a和b的 聲 速； μab=vab·rab/(rab2+η2)， 其中vab表 示粒子a和b的速度差值，rab表示粒子a和b的距離， 且vab=va?vb，rab=ra?rb， η2=0.01h2。

1.3 狀態(tài)方程

SPH方法中將流體視為弱可壓縮，其壓力由基于粒子密度來確定，如下所示[16-17]：

式中：γ為常數(shù)項(xiàng)，對(duì)于水γ取 7，對(duì)于空氣γ取 1.4；ρ為密度；ρ0和c0為設(shè)置的密度和聲速。

1.4 連續(xù)性方程

在SPH模擬弱可壓縮流體的整個(gè)過程中，通過求解質(zhì)量守恒方程或連續(xù)性方程來計(jì)算密度的變化，具體如下所示：

1.5 時(shí)間步

粒子a的速度、加速度和密度變化率分別表示為：

式中：ra為粒子a的距離；Fa為粒子a的加速度；Da是粒子a的密度變化率?；赩erlet方案[18]，對(duì)這些方程進(jìn)行積分。

1.6 液?氣兩相流

在單相中使用一致形式的壓力梯度及速度散度來更新液相的動(dòng)量方程和連續(xù)性方程。 對(duì)于氣相，因界面處的密度有較大的不連續(xù)性，故不適合采用相同的SPH公式，而使用如下形式[19]：

式中：ui，pi，ρi，?i分別為粒子i的速度、壓力、密度和哈密爾頓算子；mj，pj，ρj分別為粒子j的質(zhì)量、壓力和密度；Πij和Wij分別為黏度項(xiàng)和光滑核函數(shù)；a為沖擊加速度。

1.7 邊界條件

邊界處理的核心思想是由被視為滿足不可滲透條件的單獨(dú)一組流體粒子來描述邊界。

2 算法驗(yàn)證

2.1 平板模型

為驗(yàn)證SPH兩相流算法的精確性，本文選擇文獻(xiàn)[4]試驗(yàn)測試的平板入水結(jié)果作為依據(jù)。該試驗(yàn)采用的平板參數(shù)如下：平板質(zhì)量12.5 kg，尺寸0.45 m×0.2 m，平板從距自由液面0.1 m的高度自由下落，砰擊入水。

圖1所示為本文參照文獻(xiàn)[4]入水試驗(yàn)的平板建立的SPH模型，采用此二維模型模擬平板的入水過程。設(shè)定平板寬為0.45 m，密度為13 888.89 kg/m3，厚度為0.01 m。計(jì)算域的水域和空氣域尺寸分別如圖所示。試驗(yàn)時(shí)，平板初始位置距水面的高度為0.1 m，初始入水速度為0 m/s，重力加速度為9.81 m/s2。水及空氣的人工黏度分別為0.02和0.05，粒子的間距為2 mm，計(jì)算粒子的總數(shù)為1 028 905。圖1中，x軸對(duì)應(yīng)水域的底部邊緣，y軸代表模型的三維方向，z軸則對(duì)應(yīng)模型的中心軸線。

圖1 采用SPH兩相流算法的平板入水模型Fig.1 Simulation model of flat plate during water entry using twophase SPH method

經(jīng)模擬，分別得到平板模型入水的加速度及受到的壓力分布結(jié)果。圖2為仿真結(jié)果與文獻(xiàn)[4]試驗(yàn)結(jié)果的對(duì)比，圖3所示為模擬仿真的壓力云圖。結(jié)合兩圖可以發(fā)現(xiàn)，在130 ms時(shí)，盡管平板未接觸水面，但因平板與水面間存在空氣，空氣在平板和水面的壓縮作用下對(duì)平板產(chǎn)生壓力，導(dǎo)致平板的加速度開始緩慢上升。從圖3還可以明顯看出，當(dāng)平板接近水面時(shí)，平板底部的空氣迅速逸出，但仍有部分空氣殘留在平板底部，形成了“空氣墊”，該空氣墊將在入水過程中起到緩沖作用。

圖2 平板入水加速度數(shù)值仿真與試驗(yàn)結(jié)果對(duì)比Fig.2 Comparison of acceleration between simulations and experimental results of flat plate during water entry

圖3 平板入水過程中的壓力云圖Fig.3 Pressure contours of flat plate during water entry

文獻(xiàn)[4]試驗(yàn)結(jié)果針對(duì)的是初始高度為0.1 m的平板入水過程，得到的最大壓力為19.796 kPa，即對(duì)應(yīng)于文獻(xiàn)[4]給出的入水速度1.40 m/s的結(jié)果。根據(jù)文獻(xiàn)[4]中壓力測點(diǎn)的位置布置方式，本文提取了SPH仿真計(jì)算中平板中心點(diǎn)的壓力結(jié)果，如圖4所示。由圖可見，其最大壓力峰值為17.940 kPa。

圖4 平板試驗(yàn)中心點(diǎn)的測試壓力[4]Fig.4 Pressures at centre point of flat plate in test[4]

從圖2和圖4可以看出，在兩相流下通過SPH模擬獲得的加速度和最大壓力值略小于文獻(xiàn)[4]得到的試驗(yàn)值，且模擬的加速度脈沖寬度略大于試驗(yàn)結(jié)果。其原因在于，在仿真中空氣無法從y軸方向溢出，導(dǎo)致仿真的緩沖效果比試驗(yàn)的結(jié)果更明顯。總體上，兩相流下通過SPH模擬得到的結(jié)果與試驗(yàn)結(jié)果擬合得較好。

2.2 弓形模型

為進(jìn)一步驗(yàn)證SPH兩相流算法的精確性，針對(duì)Aarsnes[20]弓形模型的入水進(jìn)行了模擬。該文獻(xiàn)中所采用的測試設(shè)備的總重261 kg，弓形模型總長1 m，測試段長0.1 m。對(duì)于入水過程，采用壓力傳感器P1～P3測點(diǎn)（圖5（a）中所示位置）記錄壓力值。弓形模型截面的碰撞速度由速度傳感器測得，而相應(yīng)的加速度則通過速度差獲取。對(duì)于給定的沖擊加速度，沖擊力由式（10）給出。

式中：Fv為沖擊力；m為截面質(zhì)量。

圖5（b）所示為本文相應(yīng)建立的SPH模型，模型質(zhì)量為261 kg/m，水和空氣的人工黏度為0.02和0.05。為了保證計(jì)算的精度，SPH模型中的粒子間距為2 mm。在本節(jié)中，由于重點(diǎn)不是研究弓形模型的入水，而是對(duì)SPH算法的精確性進(jìn)行驗(yàn)證， 因此，僅模擬了斜升角為0°且沖擊速度為0.61 m/s時(shí)弓形模型的入水過程。模擬結(jié)果與文獻(xiàn)[20]的試驗(yàn)結(jié)果進(jìn)行了對(duì)比，如圖6和圖7所示。其中，圖6所示為弓形截面所受沖擊力，圖7所示為P1～P3測點(diǎn)的壓力。不僅如此，還將Sun[21]使用BEM方法獲得的計(jì)算結(jié)果，以及Wang[22]使用LS-DYNA有限元方法獲得的模擬結(jié)果一起進(jìn)行了比較。

圖5 弓形結(jié)構(gòu)截面及SPH模型Fig.5 Arched structure's cross-section and SPH model

從圖6可以看出，與BEM和LS-DYNA方法得到的結(jié)果相比，本文采用SPH兩相流算法計(jì)算的結(jié)果更接近于文獻(xiàn)[4]獲得的壓力峰值及其對(duì)應(yīng)的時(shí)間。

圖6 弓形結(jié)構(gòu)截面的砰擊載荷圖Fig.6 Slamming loads on the arched cross-section

從圖7可以看出，本文采用兩相流SPH算法所得結(jié)果與文獻(xiàn)[4]中3個(gè)測點(diǎn)的記錄結(jié)果相吻合，這證明了本文采用的兩相流SPH算法具有較高的準(zhǔn)確度?？傊疚脑趦上嗔飨碌腟PH仿真結(jié)果與文獻(xiàn)[4]所給出的試驗(yàn)結(jié)果擬合較好，驗(yàn)證了兩相流下SPH模擬計(jì)算方法的可靠性。

圖7 P1～P3測點(diǎn)記錄的壓力Fig.7 Pressures recorded at three measuring points

3 M型快艇結(jié)構(gòu)仿真模型

圖8所示為本文采用的M型快艇截面幾何形狀（單位：mm）。該截面長度為846.8 mm，高度為305.6 mm，重量為80 kg/s。為了分析典型位置的壓力，選擇布置了4個(gè)壓力測點(diǎn)。P1～P4測點(diǎn)分別位于點(diǎn)A、AB和CD、DE段中心處。仿真時(shí)，本文選擇了如圖9所示4種主船體斜升角α（13°，18°，23°，28°），以研究該斜升角對(duì)M型快艇船體結(jié)構(gòu)入水過程的影響。研究中，不同斜升角時(shí)的模型下落高度及質(zhì)量與文獻(xiàn)[4]相同。

圖8 M型快艇截面示意圖Fig.8 Geometry of high-speed M-boat cross-section

圖9 不同斜升角下M型快艇截面示意圖Fig.9 Diagram of high-speed M-boat cross-section with various dead-rise angles

圖10所示為相對(duì)應(yīng)的二維M型快艇截面的SPH模型。空氣域和水域尺寸分別為3 000 mm×1 000 mm和3 000 mm×1 500 mm。水及空氣的人工黏度分別為0.02和0.05，快艇截面模型的質(zhì)量為80 kg/s。同時(shí)，考慮到計(jì)算時(shí)間的影響，SPH數(shù)值仿真中粒子間距設(shè)為2 mm，后文將詳細(xì)討論粒子間距的影響。圖10中，x軸對(duì)應(yīng)于水域底部邊緣，y軸代表模型三維方向，z軸則對(duì)應(yīng)于模型中心軸線。除了考慮主船體斜升角的影響外，本文還將討論入結(jié)構(gòu)入水速度的影響，具體工況如表1所示。

表1 模擬工況Table 1 The modelling conditions

圖10 M型快艇截面的SPH模型Fig.10 The SPH model of igh-speed M-boat cross-section

4 結(jié)果分析

4.1 粒子間距的影響

本文選擇粒子間距為5，3.3，2.5，2和1.5 mm。圖11給出了斜升角為23°時(shí)計(jì)算的M型快艇截面加速度， 模擬中所有初始入水速度均取為5 m/s。從圖11可以發(fā)現(xiàn)，隨著粒子間距的減小，加速度曲線將逐漸穩(wěn)定, 且粒子間距分別為2和1.5 mm的結(jié)果吻合較好，兩者的差別較小，說明繼續(xù)縮小粒子間距對(duì)結(jié)果的影響不大。因此，考慮到計(jì)算時(shí)間的合理性，本文選擇了2 mm作為初始粒子的間距。

圖11 不同粒子間距下M型快艇截面加速度（α=23°）Fig.11 Calculated acceleration of high-speed M-boat crosssection with various particle spacings (α =23°)

4.2 加速度結(jié)果

圖12所示為在SPH兩相流模型計(jì)算中，不同斜升角和入水速度下得到的M型快艇加速度的計(jì)算結(jié)果。模擬中，船舶截面斜升角分別為13°，18°，23°和28°，對(duì)于每個(gè)斜升角，分別選擇了5，7.5和10 m/s這3種不同初始入水速度。

從圖12中可以很好地觀察到，船舶加速度曲線有2個(gè)峰值。在船舶開始運(yùn)動(dòng)后，加速度將迅速增加，達(dá)到第1個(gè)峰值后快速減小，然后再較緩慢地增加并達(dá)到第2個(gè)峰值，最后又慢慢降低。根據(jù)加速度結(jié)果，可以推測出船舶截面在入水過程中遭受了2次砰擊。關(guān)于這2次砰擊的具體位置及其影響在下文詳細(xì)說明。

圖12 不同斜升角下的加速度計(jì)算結(jié)果Fig.12 Calculation results of acceleration for various dead-rise angles

保持船舶截面斜升角不變，在初始入水速度快的工況下，加速度峰值都更大，這表明受到的砰擊載荷更大，且峰值出現(xiàn)的時(shí)間都更早。為進(jìn)一步分析斜升角對(duì)船舶入水加速度的影響，提取出不同工況下加速度的2處峰值并予以對(duì)比，結(jié)果如圖13所示。

圖13 不同初始入水速度下加速度計(jì)算結(jié)果與斜升角間的關(guān)系Fig.13 Relationship between calculated acceleration and dead-rise angle at various initial water entry velocities

從圖13可以看出，在相同初始入水速度下，隨著斜升角的增大，加速度第1個(gè)峰值減小，也即第1次砰擊載荷降低，而第2個(gè)峰值基本沒有太大變化。但是，若斜升角增大到一定程度時(shí)，可見第2個(gè)峰值顯著增加，也即第2次砰擊載荷突然增大，具體原因后文將進(jìn)行分析。

4.3 測點(diǎn)壓力結(jié)果

圖14～圖17顯示的是圖8中所描述的壓力測點(diǎn)P1～P4的壓力監(jiān)測結(jié)果。

從圖14～圖17中可以看出，在斜升角確定的情況下，隨著模型初始入水速度的增加，各測點(diǎn)的壓力峰值增大，且除測點(diǎn)P1外，其余測點(diǎn)壓力到達(dá)峰值的時(shí)間明顯提前。同時(shí)，還可以明顯看出，隨著模型初始入水速度的增加，測點(diǎn)P1壓力的脈寬減小。這說明船舶初始入水速度越快，受到的壓力載荷會(huì)越大，入水所需時(shí)間越短。

圖14 斜升角為13°時(shí)測點(diǎn)記錄的壓力變化Fig.14 Pressure histories recorded at measuring points when deadrise angle is 13°

圖15 斜升角為18°時(shí)測點(diǎn)記錄的壓力變化Fig.15 Pressure histories recorded at measuring points when deadrise angle is 18°

圖16 斜升角為23°時(shí)測點(diǎn)記錄的壓力變化Fig.16 Pressure histories recorded at measuring points when deadrise angle is 23°

圖17 斜升角為28°時(shí)測點(diǎn)記錄的壓力變化Fig.17 Pressure histories recorded at measuring points when deadrise angle is 28°

根據(jù)上節(jié)中加速度曲線分析中2次砰擊及對(duì)應(yīng)的時(shí)間可以推斷，P1～P2測點(diǎn)的壓力在第1次砰擊過程中達(dá)到峰值，在第2次砰擊過程中卻波動(dòng)不大。而P3～P4測點(diǎn)的壓力則相反，在第2次砰擊過程中達(dá)到峰值。同時(shí)，結(jié)合測點(diǎn)位置可以得出結(jié)論，即第1次砰擊發(fā)生在船舶主船體的龍骨處，而第2次砰擊則發(fā)生在連接橋和輔助船體之間。

4.4 M型快艇入水過程

為了更詳細(xì)了解M型快艇入水過程及過程中加速度和壓力變化，本文選取初始入水速度為5 m/s、斜升角分別為13°和28°的算例，提取典型時(shí)刻的壓力云圖，如圖18所示。

從圖18（a）并結(jié)合圖19加速度曲線可以看出：

圖18 不同斜升角下自由液面壓力云圖及高度Fig.18 Pressure contours and free surface elevation at various dead-rise angle

在t=0 ms時(shí)，船舶開始向下運(yùn)動(dòng), 約7 ms時(shí)，對(duì)應(yīng)于圖19中的加速度曲線第1個(gè)峰值達(dá)到最大值，此時(shí)，船舶結(jié)構(gòu)上的最大壓力位于截面的龍骨處。這進(jìn)一步表明加速度的第1個(gè)峰值是由水面與龍骨之間的相互作用引起的，即第1次砰擊發(fā)生在船舶的龍骨位置。

在t=21.4 ms的瞬間，對(duì)應(yīng)于圖19中A處，水射流砰擊到連接橋，使得加速度曲線出現(xiàn)了一個(gè)小尖峰。

圖19 v=5 m/s時(shí)不同斜升角的加速度結(jié)果Fig.19 Acceleration results for various dead-rise angle at v=5 m/s

在t=25.2 ms時(shí)，船的輔助船體接觸水面，之后殘留于主船體與輔助船體間的空氣因無法溢出而受到壓縮，進(jìn)而對(duì)船舶施加了一個(gè)作用力，使得加速度曲線開始緩慢上升，并在31 ms時(shí)刻瞬間達(dá)到峰值，形成加速度的第2個(gè)峰值。這表明加速度的第2個(gè)峰值發(fā)生在船的輔助船體接觸水面后，即第2次砰擊發(fā)生在連接橋和輔助船體之間，也印證了壓力分析的結(jié)果。

同時(shí)，從圖18（a）中還可以看出，空氣對(duì)于船舶入水有著非常重要的影響，特別是在M型快艇的側(cè)板接觸水面之后的入水過程。

對(duì)于圖18（b），即主船體體斜升角為28°時(shí)的M型快艇截面壓力云圖，與圖18（a）所示斜升角為13°時(shí)的壓力結(jié)果相似。但是，由于斜升角的增大，主船體底部側(cè)面的長度增加，水面淹沒主船體底部的時(shí)間變長，導(dǎo)致主船體船底射流砰擊連接橋的過程發(fā)生在輔助船體接觸水面之后，時(shí)間對(duì)應(yīng)于圖19中的A處和圖18(b)中的t=24 ms時(shí)刻。這也是主船體斜升角為28°時(shí)加速度曲線的第2個(gè)峰值突然增加的原因所在。

圖20所示為斜升角為28°、初始入水速度v=5 m/s、有/無空氣時(shí)M型快艇入水計(jì)算得到的加速度結(jié)果。從圖20可以發(fā)現(xiàn)，加速度曲線的第1個(gè)峰吻合較好，而對(duì)于第2個(gè)峰值有空氣作用的結(jié)果明顯小于沒有空氣作用的結(jié)果，同時(shí)有空氣作用的加速度曲線脈寬更寬。這說明空氣在第2個(gè)峰值形成過程中起到緩沖作用，在輔助船體接觸水面后，受到壓縮空氣的作用，減緩了模型的下落，使得船體受到的砰擊載荷更均勻，起到類似于“空氣墊”的作用。

圖20 斜升角為28°時(shí)有/無空氣的加速度結(jié)果Fig.20 Acceleration results with or without air influence when dead-rise angle is 28°

在第2次砰擊發(fā)生時(shí)，因部分空氣無法及時(shí)溢出，在主船體周圍形成了氣穴，如圖21所示。這類似于文獻(xiàn) [2]在三體船的入水測試中觀察到的氣穴結(jié)果，如圖22所示。根據(jù)文獻(xiàn)[2]的測試結(jié)果，隨著開始入水，空氣流向形成氣墊的龍骨位置。該結(jié)果表明，SPH兩相流算法可以輕松模擬出結(jié)構(gòu)入水過程空氣中形成的渦旋，而有限元方法則很模擬得到。

圖21 氣穴Fig.21 Air pocket

圖22 三體船的入水測試圖[2]Fig.22 View of water entry test of trimaran [2]

5 結(jié) 論

本文采用SPH兩相流算法研究了M型快艇入水過程，利用該算法模擬了平板和弓形模型的入水問題，并將仿真結(jié)果與文獻(xiàn)中已有試驗(yàn)結(jié)果進(jìn)行了對(duì)比，從而驗(yàn)證了所用算法的精確性，通過仿真模擬還觀察到了“空氣墊”現(xiàn)象。

通過對(duì)M型快艇截面SPH模型的仿真模擬，發(fā)現(xiàn)船舶在入水過程中受到了二次砰擊作用，分別發(fā)生于主船體與水面的接觸過程和側(cè)體接觸水面之后。其中，第2次砰擊發(fā)生的原因主要是，船舶與水面對(duì)空氣的壓縮作用形成壓力，并與主船體底部射流砰擊連接橋共同作用。受到壓縮的空氣在第2個(gè)峰值形成過程中起到了緩沖作用，使得側(cè)體接觸水面后，主船體受到的砰擊載荷更均勻，其作用類似于一個(gè)“空氣墊”。

斜升角對(duì)M型快艇的入水過程有很大的影響。對(duì)于第1次砰擊，因龍骨幾何形狀是楔形的，故對(duì)于較大的斜升角，龍骨處受到的砰擊載荷相對(duì)較小。對(duì)于第2次砰擊，在主船體斜升角大到一定程度時(shí)，結(jié)構(gòu)受到的砰擊載荷會(huì)顯著增加，這是因?yàn)镸型快艇和水面對(duì)空氣的壓縮以及主船體底部射流對(duì)連接橋砰擊的共同作用所致。因此，建議設(shè)計(jì)M型快艇時(shí)選擇合適的斜升角。

此外，第2次砰擊發(fā)生時(shí)，因部分空氣無法及時(shí)溢出會(huì)在主船體附近形成氣穴，而SPH兩相流算法能夠非常好地對(duì)此現(xiàn)象進(jìn)行模擬和分析。