六自由度三體船操縱運(yùn)動(dòng)自航模擬

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

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

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

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

1 六自由度操縱運(yùn)動(dòng)模擬數(shù)學(xué)模型及方法

1.1 坐標(biāo)系

圖1 操縱運(yùn)動(dòng)坐標(biāo)系Fig.1 The sketch of the co-ordinate system for maneuver

(1)

(2)

其中，

1.2 基本模型

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

(3)

(4)

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

(5)

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

兩相流的混合密度ρ和粘度μ可以用體積分?jǐn)?shù)α表示為：

(6)

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

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

1.3 操縱運(yùn)動(dòng)模型

本文中的三體船由噴水推進(jìn)器提供推進(jìn)力，同時(shí)通過改變噴水推進(jìn)器的噴口方向?qū)崿F(xiàn)三體船的轉(zhuǎn)向，即沒有舵，因此采用六自由度MMG模型建立的操縱運(yùn)動(dòng)模型可表示為：

(7)

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

2 三體船操縱運(yùn)動(dòng)計(jì)算模型及算例

2.1 三體船操縱運(yùn)動(dòng)計(jì)算模型

本文選取三體船型A作為方法驗(yàn)證模型，將模型A的模擬結(jié)果與試驗(yàn)結(jié)果進(jìn)行對(duì)比，從而驗(yàn)證數(shù)值方法的有效性。船型A的三體船噴水推進(jìn)自航模型試驗(yàn)在一個(gè)6 m多深，試驗(yàn)直徑600 m的人工湖上進(jìn)行，試驗(yàn)?zāi)Ｐ退€長6 m。數(shù)值方法得到驗(yàn)證后，選取另一艘與船型A相似的三體船型B，通過在不同傅汝德數(shù)和側(cè)體布局下，對(duì)六自由度的三體船回轉(zhuǎn)和Z形運(yùn)動(dòng)進(jìn)行了數(shù)值模擬，以研究航速和側(cè)體布局對(duì)三體船操縱性能的影響。文中三體船模型的主體和側(cè)體主要參數(shù)如表1所示，表中Ls為側(cè)體水線長，三體船船體布局示意圖如圖3所示。

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

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

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

圖4 計(jì)算域網(wǎng)格Fig.4 Sketch of mesh of domain

2.2 三體船操縱運(yùn)動(dòng)自航模擬方法驗(yàn)證

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

三體船型A的自航模型試驗(yàn)利用差分GPS對(duì)模型運(yùn)動(dòng)軌跡進(jìn)行監(jiān)測，利用光纖陀螺對(duì)垂蕩幅值、橫搖角、俯仰角和航向角進(jìn)行監(jiān)測，所有試驗(yàn)在一個(gè)6 m多深，試驗(yàn)直徑600 m的人工湖上進(jìn)行，三體船型A試驗(yàn)?zāi)Ｐ腿鐖D5所示。

圖5 三體船型A試驗(yàn)?zāi)Ｐ虵ig.5 Sketch of the trimaran A test model

三體船型A六自由度操縱運(yùn)動(dòng)數(shù)值模擬的所有算例都基于一臺(tái)Xeon E5-2650V4集群計(jì)算機(jī)進(jìn)行，采用10個(gè)進(jìn)程并行計(jì)算，每個(gè)算例執(zhí)行50 s。圖6～9所示為數(shù)值預(yù)報(bào)的三體船型A的六自由度回轉(zhuǎn)和Z形運(yùn)動(dòng)的模擬值和試驗(yàn)值對(duì)比?？梢钥闯?種不同側(cè)體布局的回轉(zhuǎn)和Z形操縱運(yùn)動(dòng)的數(shù)值模擬結(jié)果與試驗(yàn)結(jié)果吻合較好，但是回轉(zhuǎn)圈軌跡直徑模擬結(jié)果比試驗(yàn)結(jié)果小6%左右，Z形操縱運(yùn)動(dòng)的運(yùn)動(dòng)周期模擬結(jié)果略大于試驗(yàn)值，這主要是由于數(shù)值模擬中推進(jìn)器的力和力矩是通過經(jīng)驗(yàn)公式獲得，其結(jié)果和試驗(yàn)中的噴水推進(jìn)器產(chǎn)生的力和力矩有所不同。同時(shí)，可以看到在Z形操縱運(yùn)動(dòng)時(shí)三體船發(fā)生了明顯的橫傾。因此，可以證明本文提出的三體船操縱運(yùn)動(dòng)的數(shù)值模擬方法的有效性和正確性，該方法可以很好地模擬三體船靜水中的操縱運(yùn)動(dòng)。

圖6 回轉(zhuǎn)圈軌跡對(duì)比Fig.6 Comparison of turning circle trajectories

2.3 計(jì)算結(jié)果及分析

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

圖7 Z形操縱運(yùn)動(dòng)軌跡對(duì)比Fig.7 Comparison of zigzag maneuvering trajectories

圖8 Z形操縱運(yùn)動(dòng)艏向角/舵角時(shí)歷曲線對(duì)比Fig.8 Comparison of zigzag maneuvering heading/rudder angle time curves

圖9 Z形操縱運(yùn)動(dòng)橫傾角時(shí)歷曲線對(duì)比Fig.9 Comparison of zigzag maneuvering roll angle time curves

圖10 回轉(zhuǎn)圈軌跡對(duì)比Fig.10 Comparison of turning circle trajectories

以上結(jié)果說明，三體船回轉(zhuǎn)運(yùn)動(dòng)的回轉(zhuǎn)直徑受到航速和側(cè)體布局變化的影響明顯，但并不是單調(diào)的增加或減小，而是存在明顯的波動(dòng)變化。側(cè)體橫向間距增大，回轉(zhuǎn)直徑增大，回轉(zhuǎn)性能降低。在低速段Fn低于0.35時(shí)，回轉(zhuǎn)直徑隨航速增大而增大，隨側(cè)體縱向位置前移而減小，在高速段則由于主側(cè)體間的興波干擾和波形位置的改變導(dǎo)致回轉(zhuǎn)性能發(fā)生改變，回轉(zhuǎn)直徑在側(cè)體靠近船中時(shí)反而增大，此時(shí)相對(duì)低速時(shí)的回轉(zhuǎn)直徑也不再是單調(diào)增加。

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

圖11 Z形操縱運(yùn)動(dòng)艏向角/舵角時(shí)歷曲線對(duì)比Fig.11 Comparison of zigzag maneuvering heading/rudder angle time curves

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

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

圖12 回轉(zhuǎn)運(yùn)動(dòng)仿真波形Fig.12 Waveform of turning motion simulation

圖13 z形操縱運(yùn)動(dòng)仿真波形圖Fig.13 Waveform of zigzag maneuvering motion simulation

3 結(jié)論

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

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

3)對(duì)于航向改變性能，側(cè)體橫向間距增大，初轉(zhuǎn)期、超越時(shí)間和周期顯著增大，航向改變性能降低。側(cè)體縱向位置前移，超越角明顯增大，航向改變性能降低。

由于當(dāng)前數(shù)值模擬中推進(jìn)器的力和力矩是通過經(jīng)驗(yàn)公式獲得，其結(jié)果和試驗(yàn)中的噴水推進(jìn)器產(chǎn)生的力和力矩有所不同，這也導(dǎo)致了目前的數(shù)值預(yù)報(bào)存在一定的誤差，將來的工作將基于更為精確的模擬方法進(jìn)行三體船操縱性數(shù)值模擬研究，獲得更高精度的數(shù)值預(yù)報(bào)模擬。同時(shí)，對(duì)三體船操縱性能受到的航速和側(cè)體布局變化的復(fù)雜影響，從流場和優(yōu)化方法角度進(jìn)行更系統(tǒng)的研究。