規則波浪對水下航行體出水姿態參數的影響

職明洋， 馬貴輝， 任澤宇, 趙建洲

(哈爾濱工程大學船舶工程學院， 哈爾濱 150000)

0 引言

水下航行體出水是一個復雜的流固耦合過程，航行體出水時主要經歷水下航行階段和出水階段。在出水階段，航行體的運動受到波浪的擾動，波浪和航行體的相互作用影響了其出水姿態，以至于影響出水后的二次調整甚至會導致發射失敗，所以研究波浪對航行體出水姿態參數的影響過程及規律對豐富航行體的波浪條件下出水過程有著重要意義。

波浪對水下航行體的擾動屬于流固耦合問題，求解流固耦合問題的數學模型主要包括Navier-Stokes(N-S)模型和理想勢流模型。對于N-S模型的求解，基本思路是針對N-S方程和連續性方程所組成的控制方程組，選擇不同的離散方法和流動模型，其中在出水方面應用較為成功的算法包括流體體積法(VOF)、水平集和流體體積耦合法(CLSVOF)、標記網格法(MAC)、浸沒邊界與流體體積耦合法(IBM-VOF)及光順粒子法(SPH)等。雖然N-S模型可以考慮流體黏性的影響，但其對于計算資源需求大，在工程預報方面，計算效率相對較低。經典的勢流模型在數學上做了簡化處理，大部分計算認為流體是無黏不可壓縮的，Sun等采用基于勢流理論模型的邊界元方法，應用非線性自由表面條件，對船舶在規則波中的響應進行了3D時域分析，當波浪高于船舶甲板高度時，波面最高點和船首甲板高度差的變化周期與入射波的周期相同，的峰值與波高呈正相關；Yu等通過勢流理論和面板單元法，針對張力腿式平臺進行了一系列的系泊分析，以預防平臺在極端海況下失效。Ma等基于勢流理論，建立了完全非線性的數值水池，計算了波浪與垂直圓柱之間的相互作用，通過驗證發現在波浪不發生破碎的情況下，考慮流場黏性與不考慮流場黏性時兩者計算結果差距很小。Ni等采用理想勢流模型對航行體壁面的氣泡行為進行了研究，進一步驗證了理想勢流模型在求解流固耦合問題上的有效性。對于波浪與結構物的耦合作用，大量研究主要針對艦船、水下潛器、海洋平臺等方面，關于波浪擾動在水下航行體出水問題中的影響研究較少。

出水問題是指水下航行體穿過空氣-水界面并釋放到空氣中的過程，涉及非線性變形的復雜問題。Takamure等通過實驗研究了小球出水時自由液面的變形情況以及能量轉化問題，發現當發射深度和小球半徑之比高于一定的臨界值時，動能會以固定比例轉換為勢能，而且此時被小球穿過的自由液面的變形自相似。Wu等設置了多功能出水實驗裝置用于研究小球自由出水和強迫出水問題，結果表明，小球強迫出水后，自由面所形成的水柱高度會隨著傅汝德數的增加而增加，但是，水柱的持續時間會減少，小球自由出水時，受到的流體力主要來自加速度。所以，通過實驗研究發現，結構物穿越自由液面的過程中，涉及自由面變形和水層破裂等難題，學者們一直致力于該問題的解決。Liju等采用約束式運動平臺，研究不同尺度模型在不同運動速度下出水時自由液面的變化，同時通過邊界元方法模擬該出水過程，與實驗對比良好，但其關注的重點在于物體穿越自由液面之前，未計算自由液面的破裂。Rajavaheinthan等采用時間差分迭代求解圓柱出水過程，但是計算結果被迫在自由面破裂之前終止。直到倪寶玉提出一種水層脫落法，當自由面與結構物之間的水層厚度小于某一臨界值時，認為自由面斷裂，結構物出水。目前，結構物出水時自由液面大變形問題的處理方法局限于結構物固定不動或者是運動軌跡已經確定的情況，但是在實際的水下發射過程中，航行體出水姿態傾斜，對于該方面的出水研究較少。

Wang等基于動態網格技術建立了水下航行體三維出水彈道模型，該模型的建立為出水軌跡的研究提供了參考。朱坤等采用Mixture模型研究了波浪相位對水下航行體出水過程中水動力學特性的影響，發現波浪相位會導致航行體肩部空泡空間幾何的不對稱，進而引起航行體表面水動力分布的不對稱性。權曉波等將二階Stokes波的波面運動速度耦合到入射邊界條件中，采用有限體積法研究了五級海況下水下航行體出水時的位置偏移，研究發現相對于靜水出水，五級海況下水下航行體完全出水后彈體尾部偏移1 m。陳超倩等基于獨立膨脹原理對水下航行體的水下彈道進行了仿真，所得結果和實驗數據對比良好，驗證了計算的有效性，發現發射潛艇的速度使彈道軌跡出現彎曲，航行體俯仰角偏差持續變大，不利于航行體姿態的穩定，且艇速越大，出水姿態角越大，在計算過程中雖然考慮了空化，但是空化模型對水下彈道的影響不大。張重先運用Morison公式，對波浪影響下水下航行體的出水過程進行了仿真，波浪相位對波浪擾動的大小和方向均有影響，且影響效果最為顯著，水下航行體的小尺寸性和快速性有利于降低波浪對水下航行體的擾動作用。雖然前人對波浪影響下水下航行體的出水姿態做了一定的研究，但是波浪參數對水下航行體出水姿態的影響缺少系統、全面的論述。

本文基于邊界元方法，采用規則波浪模型計及波浪入射勢，建立了波浪作用下的航行體出水模型，采用截平面的方法將航行體出水模型做了簡化處理，并與實驗數據進行了對比，驗證了模型的有效性。模擬了航行體在不同波浪條件下的出水姿態，研究了浪向、波高、初相位以及周期對水下航行體的出水姿態的影響，對水下航行體出水姿態參數進行了時域分析。

1 理論與方法

本章應用勢流理論建立水下航行體出水的數值模型，內容包括：1)基本的控制方程及邊界條件；2)波浪模型及入射勢；3)網格劃分及單元映射過程；4)航行體出水過程模型簡化；5)有效性驗證。圖1為本文物理模型的示意圖。

圖1 物理模型示意圖Fig.1 Schematic diagram of physical model

航行體的出水過程是一個復雜的非線性流固耦合過程，本文基于邊界元方法建立了波浪作用下的航行體出水模型，具體計算流程如下：

1)在初始條件下將航行體表面單元、節點的相關信息輸入到三維計算程序中，計算得到速度勢；

2)將速度勢輸入到式(8)中，可以得到航行體表面壓力；

3)通過航行體運動學方程計算受力后航行體整體的運動響應；

4)得到更新后的時間步時刻的航行體表面單元、節點的相關信息；

5)重復迭代計算以上步驟，直到航行體出水過程結束。

1.1 控制方程及邊界條件

假定流體為理想流體，無旋，不考慮黏性作用，計算過程中考慮波浪的影響，將總速度勢由波浪入射勢和結構表面擾動勢組成，即

=+

(1)

滿足拉普拉斯方程

(2)

根據勢流基本理論，方程滿足以下邊界條件:

結構表面滿足不可穿透條件

(3)

式中，為航行體表面法向，為邊界運動速度矢量。

流場無窮遠處滿足邊界條件

(4)

即認為航行體擾動勢在無窮遠處對流場無影響。

自由液面(在=(,,)處)邊界條件

dd=0

(5)

代入伯努利方程并且攝動展開后保留一階項得到

(6)

式(6)滿足平面進行波的運動方程，是線性化的自由液面條件。

對于擾動勢，根據格林公式可以用其邊界上的速度勢及其法向導數來確定，所以需要采用邊界積分方法求解流場擾動勢

()=

(7)

式中，為航行體表面邊界面，和分別為場點和源點。自由空間的格林函數為(,)=1|-|，和分別為場點和源點的位置矢量。

將求解得出的速度勢帶入非定常的伯努利方程中，對水下航行體在運動過程中所受的壓力進行求解，求解過程中需要計算流場靜壓和大氣壓力，求解公式如式(8)所示

(8)

式中,為大氣壓。

式(8)即航行體受力迭代式，可以認為其和邊界條件(3)構成了本文數值方法的基礎。

1.2 波浪模型及入射勢

考慮到本文所有工況計算中波高相對于水深為小量，因此選取的波浪模型為滿足重力波周期性解的規則微幅波模型，進而采用微幅波理論對水下航行體在規則波浪中的運動進行分析。微幅波模型如圖2所示，圖中，為波長，為波高，為波幅。

圖2 微幅波波面Fig.2 Small amplitude wave surface

微幅波面方程如式(9)所示。

(,)=cos(-+)

(9)

式中，代表波浪頻率，代表波數，代表初相位，即初始時刻相位。

微幅波理論中，水質點運動緩慢，同時注意到微幅波中勢函數位移偏導較時間偏導為小量，故自由表面非線性邊界條件可以簡化為線性邊界條件，將微幅波問題簡化為線性問題。因此，波浪向前傳播時，波浪的入射勢為

(10)

假定波浪傳播方向與固定坐標系存在夾角，則微幅波經過坐標轉換，可得到新的數學模型

(11)

在本文中航行體水平速度方向和波浪傳播方向相同時為順浪，此時=180°；相反時為逆浪，此時=0°；兩個方向相互垂直時為橫浪，此時=90°，范圍為0°～180°。

根據入射勢，流體質點運動速度為

(12)

1.3 網格劃分及單元映射過程

首先對航行體結構表面離散化處理，離散過程如圖3(a)所示，由于結構表面三角形單元形態各異，需要把所有單元映射到平面上的直角三角形內，以便將各單元以同一標準參數化。

(a)結構表面離散

(b)單元映射過程圖3 網格處理Fig.3 Grid processing

映射如圖3(b)所示，分別把和映射到坐標軸和上，三角形的3個頂點分別落在(1，0)，(0，1)和(0，0)上，線性插值分別用于表示幾何位置函數，速度勢函數以及速度的法向分量=??，三角形的3個頂點分別用1,2,3來表示，則插值公式如下所示

(13)

三角形在公式(7)的第一項積分為

=++

(14)

式中，是三角形從整體坐標系到局部坐標系(,)上的雅克比轉換。

另外，積分公式右端在所選單元上的積分可表示為

=′+′+′

(15)

式中,=·(-)，為三角形單元的法向量。

令=′+，其中,=1,2,…，，可將積分方程(7)簡化為矩陣形式

(16)

1.4 水下航行體出水過程模型簡化

本文對于航行體穿越自由液面時涉及的相變以及網格重塑問題進行了簡化處理，如圖4所示，對于節點的重塑，超過自由液面的節點將會被祛除，并且在截取出來的面上設置一個節點作為中心點，節點信息的確定采用以下思路：對于單元的祛除，三角形單元的3個節點只要有1個節點超過自由液面那么該單元就會被刪除，但是有一種情況特殊，當三角形單元的3個節點有2個節點在自由液面以下，另外1個節點稱之為特殊節點()(=0,1,2…,,為所有該類型節點的總數)在自由液面以上時，node()的節點信息將會被用于節點信息的更新，即通過截面上方所有特殊節點()的運動參數表示截面新節點的運動信息

圖4 平面截取后的水下航行體Fig.4 Underwater vehicle after plane interception

(17)

式中，和分別為點和特殊點的坐標，和分別為點和特殊點的節點速度，自由液面以下的航行體部分繼續計算，直到航行體尾部出水。

1.5 有效性驗證

(18)

為驗證算法的有效性，將在開放水池測得的試驗結果與數值模擬結果對比。試驗為本文工況之一：航行體初始水平速度為，垂直速度1666，波面狀態為靜水。對比結果如圖5所示。

(a) 無量綱俯仰角數值模擬和實驗對比

(b) 無量綱俯仰角速度數值模擬和實驗對比圖5 水下航行體數值模擬和實驗的俯仰姿態對比Fig.5 Comparison of pitching attitude between numerical simulation and experiment

通過對比發現數值模擬和實驗的俯仰姿態結果對比良好，但是在航行體出水階段，由于對航行體出水條件的簡化，采用了以截平面節點代替出水航行體的方法，忽略了航行體在空氣中的運動和結構變化，導致該階段的偏差較水下運動階段的誤差大，航行體尾部觸水時俯仰姿態相對誤差如表1所示，實際誤差較小，可以認為計算方法有效。

表1 航行體尾部觸水時數值模擬和實驗數據的俯仰姿態參數對比Tab.1 Comparison of pitching attitude between numerical simulation and experiment when the tail of an underwater vehicle touches the water

2 結果與討論

2.1 靜水條件下出水姿態參數分析

在分析波浪參數對水下航行體出水姿態參數影響之前，需要對本文研究的出水姿態參數變化有一個較為完整的認識，因此忽略波浪的作用，分析水下航行體在靜水條件下姿態參數變化。靜水條件下水下航行體的無量綱俯仰角和無量綱俯仰角速度變化曲線如圖6所示。

圖6 靜水條件下無量綱俯仰角和無量綱俯仰角速度變化Fig.6 Comparison of dimensionless pitch angle and dimensionless pitch angular velocity under static water conditions

頭部觸水(HT)之前為水下運動階段，運動中的航行體受到的力主要是豎直方向的重力與浮力以及水平方向的流體作用力，由于其流體作用力的等效作用點位于質心前方，在轉矩的作用下導致航行體做順時針旋轉運動。航行體上升階段做俯仰運動，其受力不平衡性進一步加劇，從而導致航行體的無量綱俯仰角速度不斷地增大以及無量綱俯仰角不斷地減小。頭部觸水和尾部觸水(TT)之間為出水階段，航行體濕表面積不斷減小，整體所受浮力不斷下降，流體作用力的作用點不斷下移并逐漸位于質心下方，同時航行體受到的流體作用力也在不斷減小，其產生的作用可以使得航行體有逆時針轉動的趨勢，所以可以看到航行體無量綱俯仰角速度在尾部觸水之前不斷減小。

表2 靜水條件下頭部接觸水時刻和尾部接觸水時刻的無量綱俯仰角和無量綱俯仰角速度Tab.2 Dimensionless pitch angle and dimensionless pitch angular velocity at the moment when the head touches the water surface and the moment when the tail touches the water surface under static water conditions

2.2 波高對出水姿態參數的影響

(a) 無量綱俯仰角曲線

(b) 無量綱俯仰角速度曲線圖7 逆浪條件下波高對出水姿態參數的影響Fig.7 Influence of wave height on water-exiting attitude parameters under the condition of head seas

在水下航行體實際發射過程中，往往會更加注重頭部觸水時刻和尾部觸水時刻的出水姿態參數，為了更好探究兩個時刻的出水姿態參數同靜水條件下的出水姿態參數的關系，將逆浪條件下，頭部和尾部觸水的無量綱俯仰角和無量綱俯仰角速度與靜水條件下的做比較，得到俯仰角偏差值和俯仰角速度偏差值

(19)

得到的不同波高影響下的出水姿態參數偏差如圖8所示，可以發現俯仰角偏差呈負數，俯仰角速度偏差呈正數。根據式(19)可以發現，在逆浪條件下，相比于靜水，波高會增加航行體的俯仰程度，不利于發射成功，對水下發射呈負影響；隨著波高的增加，俯仰角和俯仰角速度偏差值近似呈線性增加，說明隨著波高的增加，波浪力也逐漸增加并對航行體尾部觸水無量綱俯仰角的影響越來越大。

(a) 俯仰角偏差

(b) 俯仰角速度偏差圖8 逆浪條件下波高對俯仰角和俯仰角速度偏差值的影響Fig.8 Influence of wave height on pitch angle and pitch angular velocity deviation under the condition of head seas

如圖9所示，順浪條件下波高對出水姿態參數的作用規律和逆浪條件下的相反，雖然出水姿態參數的變化規律仍然不變，但是該工況下的無量綱俯仰角大于靜水的無量綱俯仰角，無量綱俯仰角速度小于靜水的無量綱俯仰角速度。即相比于靜水，波高能夠抑制航行體的俯仰變化，對水下發射呈正影響，而且隨著波高的增加，對航行體俯仰變化的抑制效果更加顯著。

(a) 無量綱俯仰角曲線

(b) 無量綱俯仰角速度曲線圖9 順浪條件下波高對出水姿態參數的影響Fig.9 Influence of wave height on water-exiting attitude parameters under the condition of following seas

將順浪條件下，頭部和尾部觸水的無量綱俯仰角和無量綱俯仰角速度與靜水條件下的做比較，如圖10所示,俯仰角偏差>0，俯仰角速度偏差<0，順浪條件下不同波高對航行體的俯仰程度呈正影響，俯仰角和俯仰角速度偏差仍然呈線性變化。

(a) 俯仰角偏差

(b) 俯仰角速度偏差圖10 逆浪條件下波高對俯仰角和俯仰角速度偏差值的影響Fig.10 Influence of wave height on pitch angle and pitch angular velocity deviation under the condition of head seas

2.3 周期對出水姿態參數的影響

周期為兩個相鄰波峰經過海面上同一點的時間間隔，周期越長的波浪完成該過程所需的時間就越長，由于頻率=2π，且由式(10)可得出水相位

(20)

在本節計算中，浪向為逆浪，波高取1.88 m，初相位0°，周期取1～13 s，間隔1 s，得到頭部和尾部觸水的無量綱俯仰角和無量綱俯仰角速度與靜水條件下的偏差，所得結果如圖11所示。相同浪向、波高和初相位的情況下，周期較小時，由式(20)可知，相位改變量較大，出水相位由相位改變量和初相位共同決定，所以圖11中俯仰角偏差和俯仰角速度偏差在小周期階段有較大起伏；而周期較大時(≥10 s)，頻率對時間的敏感度降低，加之水下航行體在水中停留時間極短，相位改變量較小，出水相位和初相位近似相等，各周期工況下航行體出水相位大致接近，此時周期對出水相位的影響較小，所以圖11中俯仰角偏差和俯仰角速度偏差在大周期階段的曲線近乎水平。

(a) 俯仰角偏差

(b) 俯仰角速度偏差圖11 周期對出水姿態偏差的影響Fig.11 Influence of period on water-exiting attitude deviation

因此，周期通過影響航行體水下航行階段這一過程中波浪的相位變化，進而間接影響了航行體出水時刻的相位。周期的變化對相位的影響呈負相關，周期較小工況對航行體出水姿態影響大于周期較大工況對其姿態的影響。同時注意到實際工程中波浪周期較長，且航行體出水速度較快，認為在出水過程中相位近似一致，故周期對航行體出水姿態的影響較小。

2.4 出水相位對出水姿態參數的影響

在海浪運動中，浪和流是存在相對運動的，一方面，流體質點的運動速度和波浪的相速度是不同的；另一方面，波浪水質點運動方向和波浪傳播方向在不同位置也存在差異。這會導致航行體在不同相位出水時，航行體出水姿態參數將會受到不同的影響。

本節計算浪向取逆浪，波高取1.88 m，為了降低頻率對時間的敏感度，周期取為100 s，從而可以近似認為初相位和航行體的出水相位相等；同時航行體出水速度快，認為頭部和尾部出水時刻相位一致。即通過改變初相位來代替頭尾出水相位的變化。工況選擇中出水相位變化為0°～315°，間隔45°。出水姿態參數變化曲線如圖12所示，無論是無量綱俯仰角還是無量綱俯仰角速度，其變化趨勢與靜水曲線是相同的，差別在變化幅度不同，觀察圖12兩圖的尾部放大曲線并結合一個周期內余弦曲線的形狀，可以將尾部出水時刻的姿態參數大致分為3類：1)波峰范圍(0°，45°，315°)，該位置出水對航行體的俯仰程度呈負影響，不利于發射成功；2)波谷范圍(135°，180°，225°)，該位置出水對航行體的俯仰程度呈正影響，有利于發射成功；3)波節范圍(90°，270°)，該位置出水時刻，波浪對俯仰參數的影響介于波谷和波峰位置。

(a) 無量綱俯仰角曲線

(b) 無量綱俯仰角速度曲線圖12 出水相位對出水姿態參數的影響Fig.12 Influence of water-exiting phase on water-exiting attitude

出水相位對出水俯仰角和俯仰角速度偏差的影響如圖13所示。圖13(a)中，分界線以下的部分為對無量綱俯仰角的負影響，不利于發射成功；以上部分為對無量綱俯仰角的正影響，有利于發射成功。圖13(b)中，分界線以上的部分為對無量綱俯仰角的負影響，不利于發射成功;以下部分對無量綱俯仰角的正影響，有利于發射成功?？梢园l現，兩圖中的偏差圖線近似為余弦曲線，與規則波浪的波形近似，這是因為流體(波浪)質點和波浪的傳播存在相對運動，流體質點運動速度的方向隨著波浪相位的改變而發生改變。通過圖13可以發現，在該計算工況下，俯仰角偏差變化隨著出水相位呈余弦變化規律，因為隨著出水相位的變化，流體質點的軌跡速度也發生改變，流體質點和波浪的傳播存在相對運動，受到波浪相位的影響，流體質點運動速度的方向近似呈余弦變化。浪流的相對運動如圖14所示，在波峰位置，波浪傳播方向和流體質點的運動方向是相同的；在波谷位置，波浪傳播方向和流體質點的運動方向相反。波峰位置和波谷位置浪流相對運動的不同，造成兩個相位下航行體尾部觸水無量綱俯仰角相反的情況：波峰位置出水，流體質點始終沿著軸正方向運動，此時彈體所受的附加慣性力偏向于航行體俯仰方向，波浪對尾部觸水無量綱俯仰角的負影響達到最大；波谷位置出水，流體質點會沿軸負方向運動，此時彈體所受波浪的附加慣性力對航行體俯仰運動呈抑制作用，波浪對尾部觸水無量綱俯仰角的正影響達到最大。

(a) 俯仰角偏差

(b) 俯仰角速度偏差圖13 出水相位對出水姿態偏差的影響Fig.13 Influence of water-exiting phase on water-exiting attitude deviation

圖14 浪流相對運動Fig.14 Relative motion of waves and fluid particles

波浪對水下航行體的作用力可分為兩部分：1)由波浪壓力場及流場引起的壓差力； 2)由于波浪質點的軌跡速度引起的附加慣性力。浪級主要影響波浪力幅值的大小，對航行體出水姿態參數有著直接影響。而浪級、初相位和周期主要影響波浪質點的運動方向，改變附加慣性力，進而對出水姿態參數產生影響(或負影響)。因此，出水相位(近似初相位)的改變對出水姿態的影響主要體現在改變了波浪質點與航行體的相對運動，進而改變了附加慣性力的作用，通過與壓差力共同作用，最終影響了航行體的出水姿態參數。

3 結論

本文通過邊界元計算方法實現了對水下航行體水下發射的數值模擬，分析了波浪參數(浪向、波高、初相位和周期)對出水姿態參數的影響，得到以下結論：

1)波高越大，對航行體出水的姿態參數影響越大，但浪向和出水相位影響波高對出水姿態參數的影響效果。

2)周期對航行體出水姿態參數的影響實際上是影響了航行體水下航行階段過程波浪的相位變化，進而改變了航行體出水時刻的相位。隨著周期增大，頻率對時間的敏感度減小，周期對航行體的出水姿態參數的影響逐漸減小。但在實際工程中波浪周期較長，航行體出水為瞬態，認為整個出水過程相位一致，可以不計周期對航行體出水姿態的影響。

3)出水相位主要通過改變波浪水質點與航行體的相對運動影響出水姿態。逆浪時，波浪傳播和水質點運動同向的波峰位置對俯仰程度的負影響最大，不利于發射成功；波浪傳播和水質點運動反向的波谷位置對俯仰程度的正影響最大，有利于發射成功；波節位置對俯仰程度的影響介于波峰和波谷位置之間。