石曉川,李春陽(yáng),孫光曉,仉永超,張祥,王玉帥

1.濰柴動(dòng)力股份有限公司,山東 濰坊 261061;2.濰柴重機(jī)股份有限公司,山東 濰坊,261108

0 引言

雙體客船具有較好的穩(wěn)定性和操作靈活性,在對(duì)抗惡劣海況時(shí)表現(xiàn)優(yōu)異,得到廣泛應(yīng)用。船舶阻力是衡量船舶性能的重要指標(biāo)之一,阻力不僅影響船舶的速度,同時(shí)也與船舶的使用性和經(jīng)濟(jì)性相關(guān)[1]。船舶在定速巡航工況下,主機(jī)的輸出功率主要用于抵抗航行中產(chǎn)生的阻力,因此,船舶阻力的預(yù)測(cè)有利于船主機(jī)的選型匹配,如何更加快速準(zhǔn)確地預(yù)測(cè)船體阻力成為船舶設(shè)計(jì)過(guò)程中的關(guān)鍵。

船舶阻力預(yù)測(cè)的研究方法主要有理論研究、試驗(yàn)研究和數(shù)值模擬。近年來(lái),隨著計(jì)算機(jī)算力的提升,計(jì)算流體動(dòng)力學(xué)(computational fluid dynamics,CFD)成為預(yù)測(cè)阻力的熱門手段[2-3]。趙丙乾等[4]使用CFD軟件計(jì)算了雙艉客船在靜水中的航行阻力,發(fā)現(xiàn)使用軟件對(duì)船體阻力特性的預(yù)測(cè)具有快速性和高精度的優(yōu)點(diǎn),在計(jì)算范圍內(nèi)船舶總阻力的仿真與試驗(yàn)結(jié)果的相對(duì)誤差在3%以內(nèi);Wang等[5]結(jié)合試驗(yàn)研究和數(shù)值模擬的方法分析了高速平面雙體船在不同工況下的阻力,發(fā)現(xiàn)在高速工況下SSTk-ω湍流模型更接近試驗(yàn)結(jié)果;文獻(xiàn)[6-7]對(duì)不同工況下的船舶阻力使用SSTk-ω湍流模型進(jìn)行仿真,所得到的結(jié)果具備參考價(jià)值。仿真模型的不斷優(yōu)化使得計(jì)算結(jié)果的精度得到了進(jìn)一步提高。

RANS法將流體的控制方程時(shí)均化處理,解決了流體在湍流工況下在時(shí)間與空間上的隨機(jī)性。文獻(xiàn)[8-10]通過(guò)RANS法結(jié)合數(shù)值模擬方法對(duì)三維船體阻力進(jìn)行預(yù)測(cè),發(fā)現(xiàn)仿真獲得的自由面興波趨勢(shì)與試驗(yàn)得到的結(jié)果趨勢(shì)相同,具有較好的實(shí)際應(yīng)用價(jià)值。

本文中使用STAR-CCM+軟件,基于RANS法預(yù)測(cè)雙體船阻力,采用二因次換算法將船模阻力換算至實(shí)船阻力,并根據(jù)阻力預(yù)測(cè)結(jié)果對(duì)船主機(jī)進(jìn)行選型匹配。

1 數(shù)學(xué)模型

1.1 控制方程

將海水視為不可壓縮黏性流體,其連續(xù)性方程為:

(1)

式中: ▽為哈密頓算子,ux、uy、uz分別為x、y、z方向上的流體的平均速度,m/s。

流體的動(dòng)量守恒方程為:

(2)

式中:ρ為流體密度,kg/m3;t為時(shí)間,s;p為微元流體壓力,Pa;ν為流體運(yùn)動(dòng)黏度,m2/s;fi為微元流體單位質(zhì)量力,m/s2;下標(biāo)i、j分別代表笛卡爾坐標(biāo)系橫坐標(biāo)與縱坐標(biāo)的方向。

對(duì)于任何形式的流體流動(dòng),均符合Navier-Stokes方程(N-S方程),但隨著雷諾數(shù)增大,流動(dòng)由層流轉(zhuǎn)變?yōu)橥牧?過(guò)度流動(dòng)使得N-S方程的初值與邊值問(wèn)題產(chǎn)生分歧解。為解決該問(wèn)題,目前常采用時(shí)間平均法對(duì)方程進(jìn)行處理,即:

(3)

將式(3)分別帶入到式(1)(2)中可以得到RANS方程[11]

(4)

(5)

1.2 自由表面數(shù)值模型

船體在水面上航行,需要考慮氣-液兩相面對(duì)船舶前進(jìn)阻力的影響[12]。自由面捕捉算法是將液相和氣相同時(shí)計(jì)算,具有較好適用性和靈活性,對(duì)于破碎波等自由液面具有更好的模擬效果。

單一流體的物性參數(shù)由構(gòu)成函數(shù)ci決定:

(6)

式中:V為微元體積,m3;A為微元表面積,m2;ud為微元表面法向矢量上的速度,m/s;n為法向矢量。計(jì)算過(guò)程中,氣相中ci=0,液相中ci=1。

2 物理模型

某雙體客船的長(zhǎng)度為33.4 m,垂線間長(zhǎng)30.9 m。將雙體船的實(shí)船尺寸與模型尺寸按10：1縮小進(jìn)行靜水場(chǎng)航行阻力計(jì)算,實(shí)船與模型的主要參數(shù)如表1所示。

表1 雙體船主要參數(shù)

根據(jù)重力相似準(zhǔn)則,實(shí)船與模型的2個(gè)流動(dòng)系統(tǒng)重力必然相似,其弗勞德數(shù)相等。弗勞德數(shù)[13]

(7)

式中:v為特征速度,m/s;g為自由落體加速度,m/s2;L為船舶特征長(zhǎng)度,m。

根據(jù)式(7)計(jì)算得到不同航速工況下的實(shí)船與模型的速度轉(zhuǎn)換結(jié)果如表2所示。

表2 實(shí)船與模型的速度轉(zhuǎn)換結(jié)果

使用計(jì)算流體力學(xué)軟件STAR-CCM+和SSTk-ω湍流模型進(jìn)行仿真計(jì)算,模擬航速為26～32 kn。雙體船三維模型如圖1所示,計(jì)算時(shí)未考慮附體對(duì)船體的影響。

a)主視圖 b)左視圖

采用RANS法結(jié)合SSTk-ω湍流模型對(duì)虛擬繞流場(chǎng)進(jìn)行求解,計(jì)算域采用長(zhǎng)方體形狀,計(jì)算域的設(shè)置示意如圖2所示。由于船體及流動(dòng)具有對(duì)稱性,所以以半船為計(jì)算對(duì)象,邊界計(jì)算定義如下:船前為船艏前1.5LPP處,速度入口;船后為船艉后5LPP處,壓力出口;船側(cè)為中線面?zhèn)确?LPP處,速度入口;船舯為船體中線面的延展面,對(duì)稱平面;頂部為設(shè)計(jì)水線以上LPP處,速度入口;底部為設(shè)計(jì)水線以下2LPP處,速度入口。

圖2 計(jì)算域與邊界條件的設(shè)置示意圖

計(jì)算域采用切割六面體網(wǎng)格,船體附近六面體網(wǎng)格的密度進(jìn)行加密處理,在靜水面附近采用開爾文波形形狀的密度盒進(jìn)行加密,計(jì)算域內(nèi)網(wǎng)格總數(shù)為220萬(wàn)左右。計(jì)算域網(wǎng)格如圖3所示。

a)計(jì)算域內(nèi)網(wǎng)格 b)局部網(wǎng)格加密

3 模型驗(yàn)證和結(jié)果分析

3.1 模型驗(yàn)證

在進(jìn)行模型數(shù)值模擬之前,需要預(yù)先檢驗(yàn)湍流模型對(duì)于雙體船阻力計(jì)算的準(zhǔn)確性。本文中參考文獻(xiàn)[14]進(jìn)行湍流模型驗(yàn)證,文獻(xiàn)[14]中使用船舶模型為雙體船,湍流模型使用修正Realizablek-ε模型。共選取7種航速工況結(jié)果進(jìn)行對(duì)比,對(duì)比結(jié)果如表3所示。由表3可知:通過(guò)SSTk-ω模型模擬得到的結(jié)果相對(duì)誤差較小,相對(duì)誤差絕對(duì)值最大為5.9%,模擬誤差在可接受范圍內(nèi),SSTk-ω模型以及CFD軟件計(jì)算船舶阻力可靠。

表3 不同湍流模型下總阻力計(jì)算結(jié)果與試驗(yàn)數(shù)據(jù)對(duì)比

3.2 結(jié)果與討論

使用STAR-CCM+對(duì)雙體客船進(jìn)行不同航速下的靜水阻力仿真計(jì)算,航速為26～32 kn,該船在此航速階段航行時(shí),船體姿態(tài)變化明顯。自由面的波形體現(xiàn)出船體的興波特性,數(shù)值計(jì)算收斂后,監(jiān)測(cè)到的自由面波形如圖4、5所示。

a)自由波形圖 b)船側(cè)興波云圖

a)自由波形圖 b)船側(cè)興波云圖

由圖4、5可知:相比于vc=26 kn,當(dāng)vc=32 kn時(shí),隨著弗勞德數(shù)提高,船兩側(cè)波峰提高且范圍增加,興波范圍擴(kuò)大,船艏波峰后移,船艉后興波波峰更加明顯,波峰高度有所提高。

不同工況下的船底壓力分布如圖6所示。由圖6可知:海水壓力主要作用于2個(gè)分離船體上且高壓力點(diǎn)集中在船體的中后方;隨著航速提高,船體受到的壓力逐漸增大,航速為26 kn時(shí),最大船舶壓力為1 456.76 Pa,航速為32 kn時(shí),最大船舶壓力為2 123.11 Pa,相比于航速為26 kn時(shí)增大約為45.7%。隨著海水壓力增大,作用于船體的阻力隨之增大。

圖6 雙體船在不同工況下的壓力分布對(duì)比

計(jì)算出模型阻力后,使用二因次換算法將船模阻力換算至實(shí)船,規(guī)定艉傾的縱傾角為正,艏傾的縱傾角為負(fù);升沉用上升或下降高度與吃水的比值表示,下沉為負(fù),上升為正?？紤]到平板摩擦阻力因數(shù)區(qū)別于真正的摩擦阻力因數(shù),二因次法忽略了附體及風(fēng)阻影響,工程上通常采用的實(shí)船總阻力補(bǔ)償因數(shù)為10%[15]。實(shí)船有效功率

Pc=wcR′tc,

(8)

計(jì)算得到的實(shí)船總阻力及有效功率如表4所示,不同航速下實(shí)船總阻力及總阻力因數(shù)變化曲線如圖7所示。

表4 阻力計(jì)算結(jié)果

a)總阻力 b)總阻力系數(shù)

由表4及圖7可知:基于RANS法得到的實(shí)船總阻力隨實(shí)船航速的增大而增大;航行速度為32 kn時(shí),補(bǔ)償后的實(shí)船總阻力達(dá)到122.77 kN,相較于航速為26 kn時(shí)的阻力提升了約22.4%。隨著航速提高,興波阻力對(duì)船體的影響顯著,在總阻力中影響占比增大;實(shí)船總阻力因數(shù)明顯隨航速增大而減小。

船舶航行工況復(fù)雜,船-機(jī)-槳達(dá)到最佳匹配狀態(tài)可以使整個(gè)系統(tǒng)的推進(jìn)效果最優(yōu)[16]。船舶的有效功率

Pc=nη0δPs,

(9)

式中:Ps為主機(jī)額定功率,kW;n為發(fā)動(dòng)機(jī)數(shù);η0為船機(jī)槳總效率;δ為發(fā)動(dòng)機(jī)負(fù)荷率,是主機(jī)實(shí)際功率與額定功率之比。

船舶航速計(jì)算采用插值法,計(jì)算式為:

(10)

式中:uc為實(shí)際航速,kn,下標(biāo)1、2表示不同航速對(duì)應(yīng)的工況點(diǎn)。

此雙體客船采用雙機(jī)、雙槳的推進(jìn)系統(tǒng)配置,主機(jī)的額定功率為1 545 kW,船-機(jī)-槳總效率約為64%。

根據(jù)數(shù)值模擬計(jì)算得到船舶總阻力結(jié)果,匹配發(fā)動(dòng)機(jī)的負(fù)載特性,可以計(jì)算得到不同發(fā)動(dòng)機(jī)功率下的船舶航速。按照發(fā)動(dòng)機(jī)70%～100%負(fù)荷率(1 082～1 545 kW),共計(jì)算7個(gè)發(fā)動(dòng)機(jī)功率在上述條件下對(duì)應(yīng)的船舶航速,如表5所示。由表5可知:在發(fā)動(dòng)機(jī)額定功率時(shí),基本滿足客戶最大設(shè)計(jì)航速為32 kn的要求;該船主機(jī)與實(shí)船匹配效果較好,發(fā)動(dòng)機(jī)負(fù)載可以覆蓋航行常用速度,可以達(dá)到節(jié)能增效的目的。

表5 不同發(fā)動(dòng)機(jī)功率下所對(duì)應(yīng)的船舶航速

4 結(jié)論

通過(guò)RANS法對(duì)某雙體客船的船舶阻力進(jìn)行預(yù)測(cè),并應(yīng)用二因次轉(zhuǎn)換法獲得了航速為26～32 kn時(shí)的實(shí)船總阻力。

1)阻力隨著航速的提升而增大,當(dāng)航行速度為32 kn時(shí),補(bǔ)償后的實(shí)船總阻力達(dá)到122.77 kN;總阻力因數(shù)隨著航速的提升而減小,當(dāng)航行速度為32 kn時(shí),總阻力因數(shù)降低至2.91×10-3;在高航速工況下,雙體船具有較好的快速性能。

2)使用STAR-CCM+計(jì)算流體力學(xué)軟件,針對(duì)雙體船的阻力特性預(yù)測(cè)具有效率高、精度高和可靠性強(qiáng)等特點(diǎn),模擬計(jì)算的阻力結(jié)果可以更好地進(jìn)行船舶早期設(shè)計(jì)的船機(jī)匹配,確保整船效率得到提升。