魏斯行，馬 寧，顧解忡，王 楠，陸明鋒

(1.上海交通大學(xué)海洋工程國(guó)家重點(diǎn)實(shí)驗(yàn)室，上海200240；2.上海交通大學(xué)船舶海洋與建筑工程學(xué)院，上海200240；3.高新船舶與深海開(kāi)發(fā)裝備協(xié)同創(chuàng)新中心，上海200240；4.南通中遠(yuǎn)海運(yùn)川崎船舶工程有限公司，江蘇南通226000)

0 引 言

多用途船以其分艙靈活的特點(diǎn)，正逐步成為除集裝箱船、油船和散貨船三大主流船型以外的第4種主流船型。因此針對(duì)多用途船進(jìn)行型線優(yōu)化設(shè)計(jì)，以改善多用途船的阻力性能，對(duì)節(jié)能降耗和增加營(yíng)運(yùn)效益來(lái)說(shuō)都具有重要意義。傳統(tǒng)的型線設(shè)計(jì)模式主要依賴(lài)設(shè)計(jì)經(jīng)驗(yàn)以及船舶型線數(shù)據(jù)庫(kù)，耗時(shí)耗力，且方案非最優(yōu)。而近些年提出的SBD(Simulation Based Design)技術(shù)，則是結(jié)合了最優(yōu)化技術(shù)以及計(jì)算流體動(dòng)力學(xué)，使得整個(gè)型線設(shè)計(jì)流程更加高效。

近十幾年來(lái)，國(guó)內(nèi)外船型優(yōu)化尤其在以阻力最優(yōu)為目標(biāo)的船體型線優(yōu)化方面得到了迅猛的發(fā)展，主要工作集中在對(duì)中高速船的船首或船的中前部型線進(jìn)行優(yōu)化，并獲得了較好的減阻效果[1–3]。而對(duì)于低速肥大型船舶而言，相較于阻力性能，其尾部流場(chǎng)的品質(zhì)也越來(lái)越受到設(shè)計(jì)研究者的重視[4–5]。

本文以某低速肥大型的多用途船為研究目標(biāo)，以其阻力以及伴流場(chǎng)不均勻度為優(yōu)化目標(biāo)，通過(guò)RBF插值方法實(shí)現(xiàn)船體曲面變形，進(jìn)而采用Shipflow軟件計(jì)算船舶的總阻力和伴流場(chǎng)不均勻度，結(jié)合非劣分類(lèi)遺傳算法以及Kriging代理模型快速求解得到性能最優(yōu)的型線，并通過(guò)STAR CCM+對(duì)優(yōu)化結(jié)果進(jìn)行驗(yàn)證，根據(jù)優(yōu)化結(jié)果為該多用途船提出了首尾型線趨勢(shì)性建議。

1 船體型線優(yōu)化

1.1 優(yōu)化目標(biāo)

以某多用途船作為研究對(duì)象，分別考慮船首和船尾形狀改變對(duì)船舶水動(dòng)力性能的影響。優(yōu)化工況選取設(shè)計(jì)航速13.77 kn，吃水11m。船首形狀改變主要影響興波阻力，因此以興波阻力為優(yōu)化目標(biāo)；船尾型線優(yōu)化時(shí)，則以總阻力和伴流場(chǎng)不均勻度為優(yōu)化目標(biāo)。該船的主尺度信息以及優(yōu)化時(shí)的約束條件見(jiàn)表1。

表1 多用途船參數(shù)及約束條件Tab.1 Parameters and restrictions of multipurpose ship

1.2 船體曲面變形方法

基于徑向基函數(shù)[6]（RBF）的曲面變形方法屬于船體曲面局部變形方法，該方法主要通過(guò)局部控制點(diǎn)的位移實(shí)現(xiàn)船體曲面的變形，并保證變形后船體曲面的連續(xù)性與光滑性。

船體表面節(jié)點(diǎn)的位移可以通過(guò)插值函數(shù)S(x)描述，其中x=(xj,yj,zj)為船體表面網(wǎng)格節(jié)點(diǎn)。插值函數(shù)S(x)有如下表述:

式中：xj=(xj,yj,zj) 為徑向基函數(shù)中心；p(x)為多項(xiàng)式；N為控制點(diǎn)的個(gè)數(shù)； ‖x‖表 示歐式距離；?為對(duì)應(yīng)于歐氏距離的徑向基函數(shù)。本文選用W end land's方程，其具有2階連續(xù)性：

p(x)

式中：多項(xiàng)式 在本文中采用線性多項(xiàng)式以實(shí)現(xiàn)旋轉(zhuǎn)平移，即

基于上述方法，為了實(shí)現(xiàn)船首和船尾的變形，首先需要對(duì)可變控制點(diǎn)進(jìn)行布置。布置的原則是可變控制點(diǎn)可以使得首尾形狀UV度在一個(gè)較大的范圍內(nèi)變化，對(duì)于船首型線，可變控制點(diǎn)應(yīng)該布置在船首的進(jìn)流段，對(duì)于船尾型線，應(yīng)該布置在船尾的去流段，所有可變控制點(diǎn)的變化方向均為Y方向。為了保證船中型線不受影響，還應(yīng)在船中部分添加固定控制點(diǎn)。具體布置情況見(jiàn)圖1，可變控制點(diǎn)的詳細(xì)布置信息見(jiàn)表2。

圖1 可變控制點(diǎn)的位置Fig.1 Location of variable control points

表2 可變控制點(diǎn)的參數(shù)Tab.2 Parametersof variable control points

1.3 阻力和伴流場(chǎng)計(jì)算方法

圖2 Shipflow計(jì)算區(qū)域劃分Fig.2 Computational domain in Shipflow

借助流體力學(xué)分析軟件Shipflow進(jìn)行水動(dòng)力性能快速計(jì)算。該軟件對(duì)不同的阻力成分，在計(jì)算和分析時(shí)會(huì)采用不同的理論模型：對(duì)于興波阻力，計(jì)算時(shí)采用Rankine源非線性勢(shì)流理論方法（船體表面壓力積分和船尾能量法）；對(duì)于摩擦阻力，采用邊界層動(dòng)量積分方法進(jìn)行求解；對(duì)于低速肥大型船中占總阻力成分較小的粘壓阻力，計(jì)算分析時(shí)采用k?ε湍流模型來(lái)求解RANS方程，計(jì)算域從船體中部開(kāi)始。不同方法的計(jì)算域劃分如圖2所示。Shipflow通過(guò)對(duì)不同阻力成分進(jìn)行分區(qū)計(jì)算，雖然可以大大縮短計(jì)算周期，但和完全基于RANS方程求解的CFD方法相比，計(jì)算精度卻下降了。表3為該多用途船的阻力性能實(shí)驗(yàn)結(jié)果與Shipflow計(jì)算結(jié)果和完全基于RANS方法的流體力學(xué)計(jì)算軟件STAR CCM+的計(jì)算結(jié)果對(duì)比。Shipflow對(duì)于總阻力的預(yù)報(bào)結(jié)果誤差在5%左右，STAR CCM+的精度則更高，誤差僅有0.51%。因此，對(duì)于船型優(yōu)化來(lái)說(shuō)，可以采用預(yù)報(bào)時(shí)間更短的Shipflow計(jì)算船舶水動(dòng)力性能，采用預(yù)報(bào)精度更高的STAR CCM+來(lái)驗(yàn)證優(yōu)化結(jié)果的可靠性。

表3 多用途船阻力數(shù)值計(jì)算值與試驗(yàn)結(jié)果對(duì)比Tab.3 Comparison between numerical simulation and test results of resistanceof multipurpose ship

對(duì)于多用途船的尾部伴流品質(zhì)的評(píng)估，本文采用的是荷蘭MARIN水池提出的伴流目標(biāo)函數(shù)WOF(wake object function)[7]。WOF越小說(shuō)明伴流越均勻，反之則均勻性越差。WOF的計(jì)算公式如下：

其中：R為螺旋槳半徑；φ值為0°～350°，間隔10°；Vaxial為水流軸向速度。

1.4 優(yōu)化算法

采用帶精英策略的非支配排序遺傳算法NSGA-II[8]在設(shè)計(jì)空間搜索滿足排水量約束條件的最優(yōu)船型。考慮到時(shí)間成本與計(jì)算精度，船首優(yōu)化時(shí)選取的初始種群個(gè)體數(shù)為30，迭代次數(shù)設(shè)置為30；船尾優(yōu)化時(shí)選取的初始種群個(gè)體數(shù)為50，迭代次數(shù)設(shè)置為50。

2 Kriging近似模型

通過(guò)引入代理模型，可以解決船型優(yōu)化中傳統(tǒng)CFD計(jì)算耗時(shí)的問(wèn)題。首先通過(guò)拉丁超立方采樣方法[9]，生成100個(gè)樣本集作為輸入值；然后使用Ship flow對(duì)這些樣本點(diǎn)所對(duì)應(yīng)的船型進(jìn)行水動(dòng)力性能計(jì)算，計(jì)算結(jié)果作為輸出值；最后利用這些輸入輸出值構(gòu)造Kriging近似模型，并驗(yàn)證近似模型的可信度。

圖3 為尾部控制點(diǎn)變量1（Y1）和變量2（Y2）在取值范圍內(nèi)的分布，可以觀察到變量是隨機(jī)且均勻分布在取值范圍內(nèi)的。

圖3 部分變量在二維空間投影Fig.3 Projection of part of variables

K riging[10]代理模型對(duì)于待測(cè)點(diǎn)的響應(yīng)預(yù)測(cè)值表達(dá)如下：

其中：Y為由n個(gè)樣本點(diǎn)所組成的n維列向量；f為常值列向量，其元素全部為1；R為n階對(duì)稱(chēng)方陣，其元素為表示任意2個(gè)樣本點(diǎn) (xi,xj)之間空間相關(guān)性的相關(guān)函數(shù)R(xi,xj)；β?為回歸模型參數(shù)；r為待測(cè)點(diǎn)和樣本點(diǎn)(x,xi) 之間空間相關(guān)性相關(guān)函數(shù)R(xi,xj)所組成的n維向量。

圖4 顯示了船首優(yōu)化時(shí)，K riging模型響應(yīng)預(yù)測(cè)值與實(shí)際CFD計(jì)算結(jié)果的對(duì)比，同時(shí)計(jì)算兩者的均方根誤差可以得到RMS=0.0176?？梢钥闯觯瑑烧呓Y(jié)果吻合良好。該代理模型可用于船型優(yōu)化。

3 優(yōu)化結(jié)果分析

3.1 首部型線優(yōu)化結(jié)果

興波阻力對(duì)船首形狀變化敏感，因此將興波阻力設(shè)為船首優(yōu)化的目標(biāo)函數(shù)。船首優(yōu)化過(guò)程中興波阻力優(yōu)化迭代歷程如圖5所示?？梢钥闯觯芭d波阻力在迭代的過(guò)程中逐步地減小并逐步收斂至一穩(wěn)定值。圖6為多用途船船首優(yōu)化前后船側(cè)的波切圖對(duì)比，可以觀察到優(yōu)化船在靠近船首的第一個(gè)波谷值的峰值有所下降，這是興波阻力得到改善的原因。

圖4 Kriging近似模型預(yù)測(cè)值與CFD計(jì)算值結(jié)果對(duì)比Fig. 4 Comparison between Kriging approximation model predictive valueand CFD results

圖5 興波阻力優(yōu)化迭代歷程Fig.5 Iteration process of wave resistance coefficient

圖6 船首優(yōu)化前后波切圖對(duì)比Fig. 6 Comparison of wave cut profile before and after bow optim ization

根據(jù)最終的Pareto解，得到4個(gè)可變控制點(diǎn)的變化值，從而得到了優(yōu)化后的型線。通過(guò)對(duì)新的型線重新建模，并在高精度的CFD求解器STAR CCM+中進(jìn)行驗(yàn)證，驗(yàn)證后的阻力成分對(duì)比見(jiàn)表4。結(jié)果顯示，阻力優(yōu)化趨勢(shì)與Shipflow中的計(jì)算結(jié)果基本一致，其中剩余阻力減少2.47%，摩擦阻力減少0.37%，總阻力減少0.80%。

表4 船首優(yōu)化結(jié)果的STAR CCM+驗(yàn)證Tab.4 STAR CCM +verification of bow optimization results

船舶型線的變化如圖7所示。從型線變化的趨勢(shì)來(lái)看，船首橫剖線的形狀由U型過(guò)渡為V型，V型船首較于U型船首濕面積小，因此可減小摩擦阻力。同時(shí)V型艏水下部分較瘦，易于使水流沿著縱剖線方向流動(dòng)，從而減少舭部產(chǎn)生的旋渦，所以對(duì)剩余阻力性能也是有利的。

圖7 船首優(yōu)化前后橫剖線對(duì)比Fig.7 Cross-section lines beforeand after bow optim ization

因此對(duì)該多用途船來(lái)說(shuō)，在保證排水量前提下，為了獲得更好的阻力性能，首部橫剖線形狀應(yīng)該更趨向于V型。

3.2 尾部型線優(yōu)化結(jié)果

船尾形狀改變主要影響船舶粘性阻力以及尾部伴流場(chǎng)，因此船尾優(yōu)化時(shí)以總阻力以及槳盤(pán)處伴流不均勻度作為優(yōu)化目標(biāo)。船尾優(yōu)化過(guò)程中船舶總阻力與槳盤(pán)面伴流不均勻度迭代歷程如圖8所示。從優(yōu)化結(jié)果上來(lái)看，船尾型線變化導(dǎo)致總阻力和艉部伴流場(chǎng)變化，且總體上呈現(xiàn)負(fù)相關(guān)的趨勢(shì)，即船舶總阻力越小，伴流不均勻度越大。

圖8 船尾優(yōu)化結(jié)果Pareto解集Fig.8 Pareto solution of Stern optim ization result

從圖8 Pareto解集中挑選出2組具有代表性的解，對(duì)K riging近似模型得到的優(yōu)化結(jié)果進(jìn)行驗(yàn)證，通過(guò)Shipflow計(jì)算結(jié)果見(jiàn)表5?？梢杂^察到，當(dāng)尾流場(chǎng)不均勻度大大減少了近29%時(shí)（opt1），船舶總阻力卻增加了4.27%，這顯然不是希望得到的優(yōu)化船型。而方案opt2中的WOF值減小23.53%，同時(shí)總阻力也減小了2.79%。因此選擇op t2中的優(yōu)化船舶進(jìn)行模型重建，并在STAR CCM+中進(jìn)行驗(yàn)證。

STAR CCM+驗(yàn)證后的阻力成分對(duì)比如表6所示。結(jié)果顯示，船舶總阻力增加了0.89%，其中剩余阻力增加了5%左右，摩擦阻力減少了0.16%。這說(shuō)明船舶阻力的改變主要是因?yàn)槭Ｓ嘧枇Φ淖兓鸬摹?/p>

船舶型線的變化如圖9所示。從優(yōu)化結(jié)果來(lái)看，船尾型線在接近推進(jìn)器時(shí)，形狀逐漸由V型朝著U型變化。圖10為尾部型線優(yōu)化前后槳盤(pán)處速度云圖。原始船V型尾的存在使得船底下部的流體順著船體舭部流向螺旋槳盤(pán)面，由于邊界層的影響，在螺旋槳盤(pán)面的中心垂直平面易形成高伴流區(qū)；優(yōu)化船所對(duì)應(yīng)的U型尾，在靠近尾軸處上方型線呈現(xiàn)“內(nèi)凹”狀，船尾舭部半徑變小，因此水流在舭部更容易發(fā)生邊界層分離，形成舭渦，從而削弱了槳盤(pán)面上部產(chǎn)生的伴流峰值，從而使得伴流場(chǎng)均勻性得到改善。但由于舭渦的存在會(huì)造成能量損失，因此優(yōu)化船的粘壓阻力較之前會(huì)有所增大。這正是優(yōu)化船槳盤(pán)面處不均勻度得到改善，但總阻力有所增大的原因。

表5 優(yōu)化方案及總阻力與伴流不均勻度變化Tab.5 Optimization scheme and variation of total resistance and wake non-uniformity

表6 船尾優(yōu)化結(jié)果的STAR CCM+驗(yàn)證Tab.6 STAR CCM + verification of stern optimization results

圖9 船尾優(yōu)化前后橫剖線對(duì)比Fig.9 Cross-section lines beforeand after stern optimization

圖10 船尾優(yōu)化前后槳盤(pán)面處速度分布對(duì)比Fig.10 Velocity distribution at propeller disk before and after stern optim ization

因此對(duì)該多用途船來(lái)說(shuō)，為了獲得更均勻伴流場(chǎng)，同時(shí)阻力不至于惡化太多，可以使船尾型線在接近推進(jìn)器時(shí)，形狀逐漸由V型朝著U型變化。

4 結(jié)語(yǔ)

本文將多用途船作為研究對(duì)象，采用徑向基函數(shù)插值方法實(shí)現(xiàn)船體曲面變形，利用Kriging近似模型替代耗時(shí)的CFD計(jì)算以提高優(yōu)化效率，分別對(duì)多用途船船首和船尾形狀進(jìn)行優(yōu)化。優(yōu)化結(jié)果表明：對(duì)于該多用途船船首形狀，采用V型船首可以獲得更好的阻力性能；對(duì)于船尾形狀，由于尾部伴流場(chǎng)不均勻度與阻力性能呈現(xiàn)矛盾的關(guān)系，因此為了獲得更均勻的伴流場(chǎng)同時(shí)阻力不至于惡化太多，可以使船尾橫剖線在接近推進(jìn)器時(shí)，形狀逐漸由V型朝著U型變化。