六自由度三體船操縱運動自航模擬

姜帆， 李云波， 龔家燁 (.哈爾濱工程大學 船舶工程學院，黑龍江 哈爾濱 5000；.上海海事大學 海洋科學與工程學院， 上海 0306)

三體船是公認的高性能船舶之一，因其2個側體與主船體共享1個甲板及上層建筑，具有較大的甲板面積和良好的橫搖穩定性。目前三體船的水動力性能研究主要集中在阻力[1-3]和耐波性能[4-6]，而操縱性能對于高速三體船在高速航行中的危險規避和航行安全又極為重要，又由于三體船比傳統船型較寬，隨著技術的成熟和實船建造數量的增多，它的機動性、回轉性能和航向的糾偏能力研究和預報越來越緊迫。同時，國際海事組織[7-9]對船舶的操縱能力提出了更嚴格的要求，有必要對三體船操縱性做深入研究和特性預報分析，以指導后續的研究和實船設計。

目前，船舶操縱性能預報方法主要有數值方法和試驗方法。然而，試驗方法由于場地和設備成本高昂、試驗過程復雜，以及操縱運動中模型周圍流場無法觀測等原因，同時隨著數值方法和計算機能力的不斷發展，人們開始越來越多的采用數值方法代替試驗方法進行船舶操縱性能研究。自20世紀初船舶操縱運動方程[10]和船舶運動數學模型[11-15]的提出，數值方法開始在船舶操縱性研究中得到發展。平野等[16]提出了一種考慮橫傾運動的大橫傾船舶回轉運動的數值模擬方法。文獻[17-19]采用勢流方法對三體船型操縱運動進行了仿真計算。張豪[20]采用STAR-CCM+軟件對三體船斜拖運動和旋臂運動進行了數值模擬，對不同側體布局的三體船回轉運動進行了預報及研究。郭海鵬等[21-22]使用STAR-CCM+對ONR Tumblehome船模型的回轉、靜態漂移和橫傾試驗進行了數值模擬，并基于四自由度 MMG模型預報了該船型在靜水中的回轉和Z形操縱運動。文獻[23-24]利用STAR-CCM+軟件分別對KCS集裝箱船和KVLCC2油輪的約束模型試驗進行了數值模擬，并基于操縱運動數學模型預報了船舶的回轉和Z形操縱運動。

綜上所述，當前三體船操縱運動數值預報方法還是以勢流方法為主，同時船舶操縱運動的數值預報研究通常采用三自由度或四自由度的運動數學模型，但這對于三體船的操縱性預報而言卻忽略了航行姿態對預報結果的影響，而這可能導致預報結果和實際情況存在較大偏差，因此有必要開展六自由度的三體船操縱運動的精確預報。本文將考慮航態的影響，研究不同航速及片體布局的三體船在靜水中水平面內的操縱運動。基于開源CFD求解器OpenFOAM，利用六自由度MMG模型，對三體船在靜水中的回轉和Z形的操縱運動進行了數值模擬。為了驗證數值方法的有效性，以某三體船型A為例進行了數值模擬并與試驗結果進行對比。利用驗證后的理論方法對另一三體船型B在不同傅汝德數和片體布局下的操縱運動進行了預報。

1 六自由度操縱運動模擬數學模型及方法

1.1 坐標系

圖1 操縱運動坐標系Fig.1 The sketch of the co-ordinate system for maneuver

(1)

(2)

其中，

1.2 基本模型

本文計算域流場求解的控制方程用三維不可壓縮非定常雷諾平均N-S(URANS)方程表示：

(3)

(4)

湍流模型采用SSTk-ω模型[25]。自由面求解采用帶有人工可壓縮項的流體體積函數(volume of fluid，vof)方法[26]，僅在空氣與水的交界面處保持連續和非零，體積輸運用平流函數表示為：

(5)

式中：Ur為用于壓縮界面的速度場，Ur=min{cα|U|,max(|Uf|)}；cα為壓縮系數，取默認值1；α為歐拉兩相流的體積分數，代表液體部分所占體積百分比，在空氣中，α=0;水中α=1；在水平面0<α<1。

兩相流的混合密度ρ和粘度μ可以用體積分數α表示為：

(6)

本文的耦合離散方程采用PIMPLE算法求解，計算域尺寸和邊界條件如圖2所示，計算域總長度為4.5Lm，寬度為3Lm，水深為1.5Lm，自由面以上高度為0.5Lm，Lm為主體水線長。邊界條件：船體表面、計算域底部和頂部設置為固壁邊界條件，計算域后部設置為壓力出口，其他計算域外面的3個位置設置為速度入口。在2個側面、出口和入口處，各設置了0.5Lm的人工數值阻尼消波區，來降低數值震蕩，由于模擬的是靜水條件下的回轉和Z形運動，因此在入口處的速度設置為0，出口處的壓力梯度設置為0，以確保不會發生回流。由于本文中的模擬沒有風，所以對于α<0的部分，初始速度設置為0。

圖2 計算域和邊界條件Fig.2 Diagram of domain size and boundaries

1.3 操縱運動模型

本文中的三體船由噴水推進器提供推進力，同時通過改變噴水推進器的噴口方向實現三體船的轉向，即沒有舵，因此采用六自由度MMG模型建立的操縱運動模型可表示為：

(7)

式中：m為船舶質量；(Ix′,Iy′,Iz′)為船在局部坐標系下的慣性矩；(FHx′,FHy′,FHz′)和(MHx′,MHy′,MHz′)分別為船在局部坐標系下產生的水動力和力矩；(FPx′,FPy′,FPz′)和(MPx′,MPy′,MPz′)則分別為噴水推進器在局部坐標系下產生的力和力矩。操縱運動的力和力矩是通過噴水推進器噴口轉動實現的，所以沒有舵的力和力矩項。本文沒有模擬實際的噴水推進器，力和力矩是通過經驗公式計算[27-29]，忽略速度損失對噴水推進器效率的影響。

2 三體船操縱運動計算模型及算例

2.1 三體船操縱運動計算模型

本文選取三體船型A作為方法驗證模型，將模型A的模擬結果與試驗結果進行對比，從而驗證數值方法的有效性。船型A的三體船噴水推進自航模型試驗在一個6 m多深，試驗直徑600 m的人工湖上進行，試驗模型水線長6 m。數值方法得到驗證后，選取另一艘與船型A相似的三體船型B，通過在不同傅汝德數和側體布局下，對六自由度的三體船回轉和Z形運動進行了數值模擬，以研究航速和側體布局對三體船操縱性能的影響。文中三體船模型的主體和側體主要參數如表1所示，表中Ls為側體水線長，三體船船體布局示意圖如圖3所示。

表1 模型主尺度Table 1 Principal particulars of models

圖3 三體船船體及布局Fig.3 Schematic diagram of hull and layout of trimaran

數值計算中，所有域網格均通過snappyHexMesh生成，以背景網格為基礎，對船體表面周圍區域進行二次網格細化，以確保船舶運動的興波的傳播足夠平滑。船體表面采用棱柱層網格，從而精確求解流場，無量綱化的壁面到第一層網格的距離y+取25。網格示意圖如圖4所示，計算網格總數量為143萬。

圖4 計算域網格Fig.4 Sketch of mesh of domain

2.2 三體船操縱運動自航模擬方法驗證

為了驗證數值方法的有效性，本文分別對三體船型A的回轉和Z形運動進行了由噴水推進器推進的三體船自航模型試驗和數值模擬，試驗和模擬工況選取2種不同的側體位置，主側體橫向間距p=d1/Lm=0.10，主側體縱向間距a=d2/Lm=0，0.23，操舵角δ=25°，船舶初始航速Fn=0.35，合計4個工況，進而將試驗結果和數值模擬結果進行比較。

三體船型A的自航模型試驗利用差分GPS對模型運動軌跡進行監測，利用光纖陀螺對垂蕩幅值、橫搖角、俯仰角和航向角進行監測，所有試驗在一個6 m多深，試驗直徑600 m的人工湖上進行，三體船型A試驗模型如圖5所示。

圖5 三體船型A試驗模型Fig.5 Sketch of the trimaran A test model

三體船型A六自由度操縱運動數值模擬的所有算例都基于一臺Xeon E5-2650V4集群計算機進行，采用10個進程并行計算，每個算例執行50 s。圖6～9所示為數值預報的三體船型A的六自由度回轉和Z形運動的模擬值和試驗值對比。可以看出2種不同側體布局的回轉和Z形操縱運動的數值模擬結果與試驗結果吻合較好，但是回轉圈軌跡直徑模擬結果比試驗結果小6%左右，Z形操縱運動的運動周期模擬結果略大于試驗值，這主要是由于數值模擬中推進器的力和力矩是通過經驗公式獲得，其結果和試驗中的噴水推進器產生的力和力矩有所不同。同時，可以看到在Z形操縱運動時三體船發生了明顯的橫傾。因此，可以證明本文提出的三體船操縱運動的數值模擬方法的有效性和正確性，該方法可以很好地模擬三體船靜水中的操縱運動。

圖6 回轉圈軌跡對比Fig.6 Comparison of turning circle trajectories

2.3 計算結果及分析

利用上文經過驗證的數值模擬方法，對三體船型B，在3種側體布局(p=0.10,a=0；p=0.15,a=0；p=0.10,a=0)和3種航速(Fn=0.24, 0.35, 0.47)下操舵角為25°時的六自由度回轉和Z形操縱運動進行了數值模擬，探討了不同航速和側體布局對三體船操縱性能的影響。圖10所示為三體船型B數值預報的靜水中回轉運動在不同側體布局和航速下的運動軌跡對比圖。從圖中可以看出，相較于側體位置p=0.10,a=0，側體橫向間距增大，回轉直徑增大25%～29%；側體縱向位置前移，回轉直徑減小35%左右，但在高航速Fn為0.47時，回轉直徑反而增大10%。相同側體布局不同航速間比較時，航速Fn由0.24到0.35，各側體布局的回轉直徑增加8%～14%，當在高航速Fn為0.47時，相比于Fn為0.24時的回轉直徑并不是簡單的增加，當側體位于艉部，即a=0時，回轉直徑減小7%左右，當側體前移，即a=0.30時，回轉直徑則增加60%。

圖7 Z形操縱運動軌跡對比Fig.7 Comparison of zigzag maneuvering trajectories

圖8 Z形操縱運動艏向角/舵角時歷曲線對比Fig.8 Comparison of zigzag maneuvering heading/rudder angle time curves

圖9 Z形操縱運動橫傾角時歷曲線對比Fig.9 Comparison of zigzag maneuvering roll angle time curves

圖10 回轉圈軌跡對比Fig.10 Comparison of turning circle trajectories

以上結果說明，三體船回轉運動的回轉直徑受到航速和側體布局變化的影響明顯，但并不是單調的增加或減小，而是存在明顯的波動變化。側體橫向間距增大，回轉直徑增大，回轉性能降低。在低速段Fn低于0.35時，回轉直徑隨航速增大而增大，隨側體縱向位置前移而減小，在高速段則由于主側體間的興波干擾和波形位置的改變導致回轉性能發生改變，回轉直徑在側體靠近船中時反而增大，此時相對低速時的回轉直徑也不再是單調增加。

圖11所示為三體船型B數值預報的靜水中Z形運動在不同側體布局和航速下的艏向角/舵角時歷曲線對比圖。從圖中可以看出，相較于側體位置p=0.10,a=0，側體橫向間距增大，初轉期、超越時間和周期都隨之增大，增大幅度達29%～66%，超越角在Fn為0.24時變化不大，在0.35時增大32%，在0.47時則減小9%；側體縱向位置前移，超越角則明顯增大，且航速越低，增大幅度越明顯，航速Fn由0.24到0.47，超越角分別增大了78%、57%和20%。相同側體布局不同航速間比較時，隨著航速增大，初轉期、超越時間和周期承減小趨勢，超越角則是除了側體位置p=0.15,a=0，其他位置也均呈減小趨勢。

圖11 Z形操縱運動艏向角/舵角時歷曲線對比Fig.11 Comparison of zigzag maneuvering heading/rudder angle time curves

以上結果說明，三體船Z形運動的主要參數受航速和側體布局的影響明顯，但和回轉運動一樣并非是單一變化。側體橫向間距增大，初轉期、超越時間和周期顯著增大，超越角則在不同航速下有不同的變化，航向改變性能降低。側體縱向位置前移，超越角明顯增大。

圖12和圖13所示為三體船型B在靜水中回轉和Z形操縱運動過程中，自由面興波仿真結果圖，分別對應側體位置p=0.10,a=0和位置p=0.10,a=0.30，航速Fn=0.35的仿真結果。從圖中可以看出，在回轉和Z形操縱運動時，三體船船艏興波會延伸到側體位置，這種變化在航速發生變化時也同時存在，航速和側體位置的變化會導致側體遭遇的艏部興波波形和干擾情況發生變化，這也是三體船回轉和Z形操縱運動受到航速和側體布局變化的影響明顯，但并不是單調的增加或減小，而是存在明顯的波動變化的主要原因之一。同時三體船兩側波面存在的高度差導致的船體兩側壓力分布不均，也使得三體船操縱性能發生變化。

圖12 回轉運動仿真波形Fig.12 Waveform of turning motion simulation

圖13 z形操縱運動仿真波形圖Fig.13 Waveform of zigzag maneuvering motion simulation

3 結論

1)三體船操縱性能受航速和側體布局變化的影響明顯，但這種影響并不是單調變化，而是存在明顯波動的變化。這種變化是由于三體船航行過程中側體遭遇的艏部興波波形和干擾會受到航速和側體布局變化的影響，同時三體船兩側波面高度差導致的壓力分布不均的原因。

2)對于回轉性能，側體橫向間距增大，回轉直徑增大，回轉性能降低。在低速段Fn低于0.35時，隨側體縱向位置前移而減小，回轉性能增加，在高速段時回轉直徑在側體靠近船中時增大，回轉性能降低。

3)對于航向改變性能，側體橫向間距增大，初轉期、超越時間和周期顯著增大，航向改變性能降低。側體縱向位置前移，超越角明顯增大，航向改變性能降低。

由于當前數值模擬中推進器的力和力矩是通過經驗公式獲得，其結果和試驗中的噴水推進器產生的力和力矩有所不同，這也導致了目前的數值預報存在一定的誤差，將來的工作將基于更為精確的模擬方法進行三體船操縱性數值模擬研究，獲得更高精度的數值預報模擬。同時，對三體船操縱性能受到的航速和側體布局變化的復雜影響，從流場和優化方法角度進行更系統的研究。