基于BPNN的球艏降阻優(yōu)化模型構(gòu)建研究

張 維 英， 周 俊 秋， 于 博 文， 于 洋

(大連海洋大學(xué) 航海與船舶工程學(xué)院，遼寧 大連 116023)

0 引 言

船舶快速性是船舶設(shè)計(jì)中要考慮的主要問題之一，減少船舶航行過程中的阻力，是提高船舶快速性的主要途徑.近年來，基于CFD的船舶阻力預(yù)報和船型快速性優(yōu)化取得了顯著的成果.葉萌等利用若干設(shè)計(jì)變量樣本點(diǎn)的CFD計(jì)算結(jié)果，通過建立BP神經(jīng)網(wǎng)絡(luò)(BPNN)近似模型成功對載重46 000 t油船進(jìn)行了型線優(yōu)化[1].董震鵬使用CFD對某拖網(wǎng)漁船成功進(jìn)行了球艏形狀降阻優(yōu)化.該文獻(xiàn)指出由于對拖網(wǎng)漁船模型采用的是全船模擬法，導(dǎo)致軟件計(jì)算時間過長，對計(jì)算機(jī)硬件要求較高，此方法不利于在中小型工作站的普遍推廣[2].

對于船型優(yōu)化中產(chǎn)生的大量船型樣本，可根據(jù)分階段采用變精度模型(variable fidelity model，VFM)概念，先后使用低與高兩種精度模型對船體進(jìn)行阻力預(yù)報，經(jīng)過高精度(high fidelity，HF)模型對低精度(low fidelity，LF)模型計(jì)算結(jié)果的補(bǔ)償后，得到逼近實(shí)際值的計(jì)算結(jié)果.采用VFM概念進(jìn)行阻力預(yù)報相比于僅用高精度模型進(jìn)行仿真計(jì)算而言，將大幅減少計(jì)算量，節(jié)約優(yōu)化成本.Robinson等利用VFM在機(jī)翼設(shè)計(jì)和蝙蝠飛行問題的研究上大幅度降低了耗時[3].Toal也驗(yàn)證了使用VFM比只使用高精度模型更能節(jié)約計(jì)算成本[4].

求解N-S方程的CFD軟件具有良好的精度和可信度，可作為高精度模型的計(jì)算工具.高精度模型的計(jì)算值在優(yōu)化過程中不必在每次迭代都進(jìn)行更新，但低精度模型在每次迭代時都需要進(jìn)行更新，故尋求一種使用方便且結(jié)果精確的低精度模型是分階段采用VFM概念進(jìn)行球艏優(yōu)化的關(guān)鍵.林煒杰使用Holtrop法作為低精度模型對船體進(jìn)行阻力計(jì)算，并成功將VFM應(yīng)用于船型初步設(shè)計(jì)[5]，該文獻(xiàn)的Holtrop法阻力預(yù)報是基于高質(zhì)量的三維模型進(jìn)行計(jì)算，每一次優(yōu)化迭代使用Holtrop法進(jìn)行阻力預(yù)報前均需要依據(jù)設(shè)計(jì)點(diǎn)參數(shù)值對船模的尺寸進(jìn)行修改，導(dǎo)致計(jì)算效率低、計(jì)算成本高，不利于VFM應(yīng)用于工程實(shí)踐.本文采用一種將近似技術(shù)與Holtrop法相融合求解低精度模型的思路，得到一種利用數(shù)學(xué)近似函數(shù)表達(dá)式來快速預(yù)報大量設(shè)計(jì)點(diǎn)的低精度模型.近似技術(shù)的出現(xiàn)提供了一種新的優(yōu)化問題求解方法，其基本思想是依據(jù)有限數(shù)量的試驗(yàn)數(shù)據(jù)，構(gòu)建一個既滿足精度要求，計(jì)算成本又相對小的數(shù)學(xué)近似函數(shù)表達(dá)式用于模擬設(shè)計(jì)優(yōu)化問題的輸入輸出關(guān)系，并用于代替實(shí)際仿真程序[6].張舒應(yīng)用BPNN原理，結(jié)合試驗(yàn)資料與圖譜構(gòu)造了神經(jīng)網(wǎng)絡(luò)阻力預(yù)報模型，證明了人工神經(jīng)網(wǎng)絡(luò)在阻力預(yù)報方面的可行性[7].Volpi等使用隨機(jī)徑向基函數(shù)(radial basis functions，RBF)建立動態(tài)徑向基近似模型，對一艘雙體船的阻力和耐波性進(jìn)行優(yōu)化，該方法展現(xiàn)出很高的優(yōu)化精度和效率[8].

本文以一艘參數(shù)化試驗(yàn)油船模型球艏的長度與高度快速性優(yōu)化為例，通過基于拉丁超立方設(shè)計(jì)法選取的樣本點(diǎn)及Holtrop法阻力預(yù)報響應(yīng)值構(gòu)建基于BPNN的球艏降阻優(yōu)化低精度模型，并使用該模型對油船進(jìn)行大量設(shè)計(jì)點(diǎn)的快速阻力預(yù)報，通過分析設(shè)計(jì)變量對總阻力的影響特性，得到逼近最優(yōu)解的設(shè)計(jì)點(diǎn).

1 標(biāo)準(zhǔn)模型阻力預(yù)報研究

本文將一艘75 000 DWT散貨船作為標(biāo)準(zhǔn)模型，采用Rhino軟件對其進(jìn)行三維建模，并采用Grasshopper(簡稱GH)軟件的船舶設(shè)計(jì)插件Nemo 1.1 Beta 的Holtrop法阻力預(yù)報功能對75 000 DWT散貨船標(biāo)準(zhǔn)模型進(jìn)行一系列航速下的阻力預(yù)報.75 000 DWT散貨船主要參數(shù)見表1，船模的縮尺比為1∶1.圖1展示了該75 000 DWT散貨船的三維模型.

圖1 75 000 DWT散貨船三維模型

表1 75 000 DWT散貨船主要參數(shù)

Holtrop依據(jù)334個模型的試驗(yàn)數(shù)據(jù)，加上系列64的高速船數(shù)據(jù)進(jìn)行回歸，使得Holtrop法運(yùn)用范圍擴(kuò)展至弗勞德數(shù)為0.55以上.艾亞法適用范圍較廣，常應(yīng)用于中、低速商船，但其統(tǒng)計(jì)資料僅代表20世紀(jì)40年代以前的船型，故用于近代船型與大型豐滿船型時通常誤差較大.蘭潑凱勒法將船的剩余阻力和摩擦阻力分開計(jì)算，剩余阻力由曲線圖譜查得，適用于單螺旋槳大型船舶.圖2為75 000 DWT散貨船的Holtrop法阻力預(yù)報值與艾亞法阻力預(yù)報值、蘭潑凱勒法阻力預(yù)報值及試驗(yàn)值的對比[9].

圖2 多方法的阻力值比較

比較圖2中各方法的阻力預(yù)報值，參照各方法的性質(zhì)，可知Holtrop法的誤差與其他兩種傳統(tǒng)經(jīng)驗(yàn)公式法相差不大，且船型適用范圍更廣，Holtrop法可作為低精度模型的求解方法.

2 基于Rhino與GH的船體建模

通過Rhino及其可視化編程插件GH建立本次試驗(yàn)的優(yōu)化對象——一艘油船三維模型，其基準(zhǔn)設(shè)計(jì)模型的主尺度參數(shù)如表2所示.

表2 試驗(yàn)油船模型主尺度參數(shù)

試驗(yàn)油船基準(zhǔn)設(shè)計(jì)模型的建立是依據(jù)表3所展示的油船船型參數(shù)取值范圍[10].試驗(yàn)油船基準(zhǔn)設(shè)計(jì)模型如圖3所示.

表3 油船船型參數(shù)取值范圍

圖3 油船三維模型

3 神經(jīng)網(wǎng)絡(luò)的選用及其應(yīng)用機(jī)理

3.1 神經(jīng)網(wǎng)絡(luò)的選用

本文在低精度模型中引入近似技術(shù)，提出一種將BPNN作為球艏降阻優(yōu)化低精度模型的方法.

工程問題中較為常用的近似模型形式有響應(yīng)面模型(RSM)、Kriging模型和神經(jīng)網(wǎng)絡(luò)模型.文獻(xiàn)[11]對以上3種近似模型形式進(jìn)行了研究與對比.由于神經(jīng)網(wǎng)絡(luò)模型適用范圍廣，預(yù)測空間相對平滑，擬合效率較高、耗時短，且阻力計(jì)算采用神經(jīng)網(wǎng)絡(luò)模型精度最高，其性能指標(biāo)的平均絕對百分比誤差(MAPE)較其余近似模型低，故本文選用在船型優(yōu)化設(shè)計(jì)領(lǐng)域應(yīng)用較為成熟的BPNN作為本次優(yōu)化問題的近似模型形式.

3.2 神經(jīng)網(wǎng)絡(luò)的應(yīng)用機(jī)理

只具有一個神經(jīng)元的單層感知機(jī)模型如圖4所示，單層感知機(jī)功能十分有限.其原理是通過均勻試驗(yàn)設(shè)計(jì)法選取多組樣本點(diǎn)，通過訓(xùn)練神經(jīng)網(wǎng)絡(luò)，計(jì)算得出滿足誤差要求的權(quán)重，則可以生成數(shù)學(xué)近似函數(shù)表達(dá)式用于模擬問題的輸入輸出關(guān)系.

圖4 單層感知機(jī)模型

設(shè)輸入為X=(x1x2…xn)T∈Rn，相應(yīng)權(quán)值為W=(w1w2…wn)T∈Rn，則單層感知機(jī)的輸出y可表示為[12]

(1)

當(dāng)采集的是高度非線性的樣本點(diǎn)時，則需要多層感知機(jī)BPNN模型來訓(xùn)練樣本點(diǎn)，根據(jù)Kolmogorov 定理，只具有一個隱含層的三層BPNN就可以實(shí)現(xiàn)對任意精度的逼近.BPNN由輸入層、隱含層和輸出層組成，通過將輸出層的校正誤差沿原信號傳遞通路反傳，并逐層調(diào)整輸出層到隱含層、隱含層到輸入層的權(quán)值.其模型圖如圖5所示.

圖5 三層BPNN模型圖

4 基于BPNN的船舶阻力預(yù)報方法

4.1 試驗(yàn)設(shè)計(jì)方法的選取及其應(yīng)用

BPNN的學(xué)習(xí)訓(xùn)練需要依據(jù)樣本點(diǎn)數(shù)據(jù)，而均勻的樣本點(diǎn)需要利用試驗(yàn)設(shè)計(jì)方法選取.本文采用具有高靈活性和廣泛適用性的拉丁超立方設(shè)計(jì)法選取構(gòu)建BPNN的樣本點(diǎn)[13].

拉丁超立方抽樣是一種從多元參數(shù)分布中近似隨機(jī)抽樣的方法，屬于分層抽樣法.該方法根據(jù)樣本點(diǎn)數(shù)量k的要求，將n維設(shè)計(jì)空間中每一維設(shè)計(jì)變量平均劃分為k個子空間，所有設(shè)計(jì)變量的子空間組合形成kn個子空間.再根據(jù)每一個設(shè)計(jì)變量的每一個子空間只出現(xiàn)一次的原則，隨機(jī)選出k個子空間，從而保證取點(diǎn)的均勻性，最后在選出的每一個子空間中隨機(jī)取一個樣本點(diǎn)，形成k個樣本點(diǎn).

4.2 Holtrop法阻力預(yù)報原理

本文選用GH的船舶設(shè)計(jì)插件Nemo 1.1 Beta軟件計(jì)算樣本點(diǎn)響應(yīng)值.Nemo軟件的阻力預(yù)報模塊是根據(jù)Holtrop在1978～1984年發(fā)表的資料為計(jì)算依據(jù)的.Holtrop法總阻力表示為[14-16]

Rt=Rf(1+k1)+Rapp+Rw+Rb+Rtr+Ra

(2)

式中：Rt是總阻力；Rf是摩擦阻力；1+k1是船型黏性阻力因子；Rapp是附體阻力；Rw是興波阻力；Rb是球艏附加壓阻力；Rtr是艉封板浸水產(chǎn)生的附加壓阻力；Ra是模型與實(shí)船相關(guān)修正因子.

Holtrop已在文獻(xiàn)[14-16]中對上式中各阻力的回歸計(jì)算公式進(jìn)行了詳細(xì)說明.因?yàn)榇驮O(shè)計(jì)初期只對裸船體進(jìn)行研究，故附體阻力為0.通過對其余各阻力回歸計(jì)算公式的參數(shù)歸納，得到總阻力的一般函數(shù)表達(dá)式：

Rt=f(v，Shull，Lwl，B，T，，Cp，Lcb，Cm，

Cb，Cw，Abt，Atr，Cstern，Hb，Tf，ie)

(3)

式中：v是航速(kn)；Shull是船體濕表面積；Lwl是設(shè)計(jì)水線長；B是設(shè)計(jì)水線寬；T是平均吃水；是型排水體積；Cp是棱形系數(shù)；Lcb是浮心縱向位置；Cm是舯橫剖面系數(shù)；Cb是方形系數(shù)；Cw是水線面系數(shù)；Abt是艏垂線處球艏橫剖面面積；Atr是零航速下的艉封板浸水面積；Cstern是船艉形狀系數(shù)；Hb是Abt的形心到基線的高度；Tf是艏吃水；ie是半進(jìn)流角.

孫源等使用基于Holtrop法的Maxsurf軟件分析DTMB5415標(biāo)準(zhǔn)模型得出不同航速下球艏縱向位置的合理變化可使總阻力降低1.26%～10.31%[17]，說明改變球艏線性尺寸的方法在降阻方式中有著重要地位.

式(3)中受球艏前端點(diǎn)距艏垂線距離lb、球艏前端點(diǎn)距設(shè)計(jì)水線距離hb影響的變量有以下4個：船體濕表面積Shull、型排水體積、艏垂線處球艏橫剖面面積Abt及其形心到基線的高度Hb，故通過改變設(shè)計(jì)變量lb、hb的值可改變上述4個變量的值，進(jìn)而改變由Holtrop法計(jì)算的總阻力Rt的值.在優(yōu)化球艏過程中保持油船總長、型深及型寬等總尺度不變，艏垂線位置隨球艏線性尺寸變化而做較小的改變.本文的設(shè)計(jì)變量lb、hb示意見圖6.

圖6 球艏幾何參數(shù)

5 基于BPNN的球艏降阻優(yōu)化

5.1 球艏降阻優(yōu)化BPNN的構(gòu)建流程

本文所構(gòu)建的球艏降阻優(yōu)化低精度模型是一種基于一定數(shù)量樣本點(diǎn)及其響應(yīng)值所訓(xùn)練的BPNN.圖7給出了構(gòu)建球艏降阻優(yōu)化BPNN的流程.通過將均勻試驗(yàn)設(shè)計(jì)法選取的樣本點(diǎn)及其響應(yīng)值送入BPNN進(jìn)行訓(xùn)練，生成可代替船舶性能計(jì)算軟件進(jìn)行快速阻力預(yù)報的低精度模型.再通過大量設(shè)計(jì)點(diǎn)的BPNN總阻力計(jì)算結(jié)果，預(yù)報總阻力隨設(shè)計(jì)變量的變化趨勢，分析最優(yōu)設(shè)計(jì)點(diǎn)的大致位置，確定最優(yōu)設(shè)計(jì)點(diǎn)的近似解，驗(yàn)證球艏降阻優(yōu)化BPNN作為低精度模型的可行性.

圖7 BPNN構(gòu)建流程圖

5.2 第一批樣本點(diǎn)的選取及球艏變形

依據(jù)文獻(xiàn)[10]，只有當(dāng)航速大于某一數(shù)值時，球艏才會顯示其降阻作用，該航速值稱為球艏的界限速度，并將該界限速度作為是否選用球艏的衡量標(biāo)準(zhǔn)[10]，該試驗(yàn)油船界限速度vk對應(yīng)的弗勞德數(shù)計(jì)算公式如下：

Frb=0.644-0.641Cb

(4)

對3 500～7 000 DWT的油船，其服務(wù)航速通常為11.5～13 kn.通過將試驗(yàn)油船的Cb代入式(4)中，得出該油船球艏的臨界弗勞德數(shù)Frb，進(jìn)而算出該臨界弗勞德數(shù)對應(yīng)的界限速度vk=12.26 kn.經(jīng)查閱船型資料，本文以油船在v=12.5 kn時的降阻為優(yōu)化目標(biāo).

本文分別以GH中試驗(yàn)油船基準(zhǔn)設(shè)計(jì)模型球艏的初始長度lb=4 693 mm和初始高度hb=5 015 mm為兩設(shè)計(jì)變量的基準(zhǔn)點(diǎn).lb的取值范圍取以基準(zhǔn)點(diǎn)為中心，以1 000 mm為半徑的閉區(qū)間；hb的取值范圍取以基準(zhǔn)點(diǎn)為中心，以500 mm 為半徑的閉區(qū)間.設(shè)計(jì)空間如下：

3 693 mm≤lb≤5 693 mm

4 515 mm≤hb≤5 515 mm

(5)

基于拉丁超立方設(shè)計(jì)法在設(shè)計(jì)空間內(nèi)選取20個樣本點(diǎn)作為BPNN的第一批學(xué)習(xí)樣本點(diǎn)，樣本點(diǎn)的二維分布如圖8所示，樣本點(diǎn)數(shù)值見表4.

表4 基準(zhǔn)點(diǎn)與部分第一批樣本點(diǎn)

圖8 第一批樣本點(diǎn)的分布

依據(jù)20個樣本點(diǎn)，在GH中對試驗(yàn)油船模型進(jìn)行球艏變形，生成20個不同球艏幾何參數(shù)的三維模型，樣本1～5的球艏模型示意如圖9所示.

圖9 部分球艏變形示意圖

5.3 第一批樣本點(diǎn)的總阻力計(jì)算及相關(guān)分析

將依據(jù)基準(zhǔn)點(diǎn)及第一批樣本點(diǎn)所生成的不同尺寸球艏的試驗(yàn)油船依次拾取到GH的船舶設(shè)計(jì)插件Nemo 1.1 Beta中進(jìn)行阻力預(yù)報.圖10展示了Nemo的船體輸入模塊，在Float區(qū)域中設(shè)置設(shè)計(jì)水線的高度，Brep電池用于拾取Rhino中的船體曲面.

圖10 Nemo的船體輸入模塊

輸入船體曲面后，Nemo即對船體進(jìn)行流體靜力學(xué)計(jì)算，圖11展示了Nemo的流體靜力學(xué)計(jì)算模塊.船體數(shù)據(jù)經(jīng)電池左側(cè)Hull_Brep接口輸入，通過電池運(yùn)算后輸出于其右側(cè)的HH_Json接口.

圖11 Nemo的流體靜力學(xué)計(jì)算模塊

在進(jìn)行阻力預(yù)報前，要先定義流體的參數(shù)和航速范圍.本文設(shè)定流體的溫度為15 ℃，類別為海水.航速上下限分別為15 kn和12 kn，梯度為0.5 kn.Nemo的流體參數(shù)及航速設(shè)置模塊如圖12所示.

圖12 Nemo的流體參數(shù)及航速設(shè)置模塊

圖13是Nemo的Holtrop法阻力預(yù)報模塊.

圖13 Nemo的Holtrop法阻力預(yù)報模塊

圖14是Nemo軟件阻力預(yù)報模塊全局布置圖.含有不同功能的電池模塊通過正確的邏輯連線耦合，使數(shù)據(jù)在各模塊間傳遞，用于計(jì)算響應(yīng)值.

圖14 Nemo的阻力預(yù)報模塊全局布置圖

圖15展示了油船在球艏尺寸基準(zhǔn)點(diǎn)處航速為12.5 kn時GH顯示面板的阻力計(jì)算結(jié)果等數(shù)據(jù).

圖15 Nemo的阻力數(shù)值顯示面板

使用船舶設(shè)計(jì)軟件Nemo在試驗(yàn)油船球艏尺寸基準(zhǔn)點(diǎn)及20個樣本點(diǎn)處測得的模型總阻力如表4所示.

總阻力在設(shè)計(jì)空間的分布如圖16所示.通過顏色標(biāo)尺觀察圖16中的散點(diǎn)，可以發(fā)現(xiàn)總阻力朝著球艏高度最小值和球艏長度最小值方向減小.

圖16 第一批樣本點(diǎn)的總阻力空間分布

為了量化變量間的相關(guān)程度，本文對第一批樣本點(diǎn)的兩個自變量與總阻力進(jìn)行Pearson相關(guān)分析.相關(guān)系數(shù)可以量化變量間的相關(guān)程度，反映變量之間的變化趨勢.變量x和y的Pearson簡單相關(guān)系數(shù)r定義如下：

(6)

式中：k為樣本容量，xj、yj為進(jìn)行相關(guān)分析的兩變量所對應(yīng)的樣本值.

數(shù)據(jù)正態(tài)分布是Pearson相關(guān)分析的前提之一.通過對3組變量進(jìn)行正態(tài)性檢驗(yàn)，得出3組變量正態(tài)性檢驗(yàn)顯著性均大于0.05，可直接進(jìn)行Pearson相關(guān)分析.hb與Rt的相關(guān)系數(shù)為0.929，Sig值小于0.001，故hb與Rt存在極顯著的強(qiáng)正相關(guān)關(guān)系.lb與Rt的相關(guān)系數(shù)為0.052，表明lb與Rt存在弱正相關(guān)關(guān)系，Sig值為0.827(Sig值大于0.05表示變量間相關(guān)性不顯著)，說明lb與Rt之間的相關(guān)性不顯著.但根據(jù)Pearson簡單相關(guān)系數(shù)的檢驗(yàn)t統(tǒng)計(jì)量公式：

(7)

可知對任意值的r，樣本容量k足夠多就能使統(tǒng)計(jì)量足夠大，通過顯著性檢驗(yàn)的可能性就高.故本文將通過對BPNN計(jì)算的大量設(shè)計(jì)點(diǎn)變量值進(jìn)行相關(guān)分析，來說明整個設(shè)計(jì)空間內(nèi)變量間的相關(guān)關(guān)系.

觀察圖16的散點(diǎn)在ZX平面和YZ平面的投影，可發(fā)現(xiàn)hb與Rt的散點(diǎn)較lb與Rt更接近一條直線，說明hb與Rt的線性相關(guān)性強(qiáng)于lb與Rt.但hb數(shù)值變化時為了保證總長等主尺度值不變，從而使艏垂線位置發(fā)生較小的改變，進(jìn)而影響了lb的數(shù)值.為進(jìn)一步分析總阻力的變化規(guī)律，并基于變化規(guī)律擴(kuò)充用于訓(xùn)練BPNN的樣本點(diǎn)數(shù)量，本文將在第一批樣本點(diǎn)的阻力最小值點(diǎn)附近開展探究球艏長度與總阻力關(guān)系的單因素試驗(yàn)，對球艏長度與總阻力進(jìn)行Pearson相關(guān)分析.

5.4 球艏長度與總阻力的相關(guān)分析

在表4樣本點(diǎn)的總阻力最小值處，將球艏高度值(hb=4 620 mm)恒定不變，取[4 019，4 419]為球艏長度的研究范圍，采用均勻設(shè)計(jì)試驗(yàn)法生成10組球艏長度試驗(yàn)樣本，在GH里對試驗(yàn)油船模型基于10組樣本生成10個船體三維模型，并使用Nemo軟件對10個模型進(jìn)行阻力預(yù)報.數(shù)據(jù)見表5.

表5數(shù)據(jù)所對應(yīng)的散點(diǎn)圖及擬合情況如圖17所示，分析圖17得出：球艏長度單因素試驗(yàn)數(shù)據(jù)基本符合線性相關(guān)關(guān)系.

圖17 總阻力的二維分布

表5 球艏長度樣本及響應(yīng)值

通過對兩變量進(jìn)行正態(tài)性檢驗(yàn)，得出球艏長度和總阻力正態(tài)性檢驗(yàn)顯著性均大于0.05，可直接使用Pearson簡單相關(guān)系數(shù)量化球艏長度與總阻力間的相關(guān)程度，經(jīng)計(jì)算得到其相關(guān)分析結(jié)果：lb與Rt的相關(guān)系數(shù)為0.631，Sig值為0.05，故在lb的[4 019，4 419]范圍內(nèi)，lb與Rt之間存在顯著的強(qiáng)正相關(guān)關(guān)系.

5.5 第二批樣本點(diǎn)的選取及總阻力計(jì)算

為了增強(qiáng)BPNN在尋優(yōu)區(qū)間的預(yù)報精度，本文在第一批樣本點(diǎn)的基礎(chǔ)上增添40個樣本點(diǎn).

依據(jù)圖16中散點(diǎn)的分布趨勢和單因素試驗(yàn)得出的變量相關(guān)分析結(jié)果，判斷出最優(yōu)設(shè)計(jì)點(diǎn)位于整個設(shè)計(jì)空間的球艏長度最小值和球艏高度最小值附近.故采用拉丁超立方設(shè)計(jì)法在該范圍內(nèi)選取40個樣本點(diǎn).第二批樣本點(diǎn)的選取范圍如下：

3 693 mm≤lb≤4 293 mm

4 515 mm≤hb≤4 815 mm

(8)

使用Nemo軟件對選取的第二批樣本點(diǎn)進(jìn)行阻力預(yù)報，部分第二批樣本點(diǎn)及總阻力如表6所示.

表6 部分第二批樣本點(diǎn)

第二批樣本點(diǎn)的總阻力在設(shè)計(jì)空間的分布見圖18.

圖18 第二批樣本點(diǎn)的總阻力空間分布

5.6 球艏降阻優(yōu)化BPNN的訓(xùn)練

將基準(zhǔn)點(diǎn)及第一批和第二批樣本點(diǎn)的集合作為輸入變量，將其響應(yīng)值即總阻力作為輸出變量，在邁實(shí)軟件里對球艏降阻優(yōu)化BPNN進(jìn)行訓(xùn)練.

圖19展示了該BPNN訓(xùn)練過程的訓(xùn)練曲線，從圖中可以看出，其均方差σ很快收斂于某極小值處，說明該神經(jīng)網(wǎng)絡(luò)的訓(xùn)練達(dá)到預(yù)期效果.本次試驗(yàn)將該極小值設(shè)置為1×10-10，最大迭代數(shù)設(shè)置為500 000，其意義是當(dāng)全局誤差達(dá)到1×10-10時或達(dá)到最大迭代數(shù)即500 000時，該神經(jīng)網(wǎng)絡(luò)的訓(xùn)練終止.

圖19 BPNN的訓(xùn)練曲線

通過合理的參數(shù)設(shè)置，可以得到精度高、誤差小的神經(jīng)網(wǎng)絡(luò)，進(jìn)而應(yīng)用于工程優(yōu)化.當(dāng)使用一個隱含層就可以滿足求解精度的時候，則建立只含一個隱含層的神經(jīng)網(wǎng)絡(luò)，這樣有利于神經(jīng)網(wǎng)絡(luò)的穩(wěn)定性，故本神經(jīng)網(wǎng)絡(luò)只建立一個隱含層.通過調(diào)節(jié)該BPNN的隱含層節(jié)點(diǎn)個數(shù)，確定了該BPNN的隱含層節(jié)點(diǎn)為6個，此時訓(xùn)練效果較優(yōu).

5.7 球艏降阻優(yōu)化BPNN的偏差分析

球艏降阻優(yōu)化BPNN原值與計(jì)算值的偏差曲線見圖20.

圖20 BPNN的偏差曲線

本文采用相對均方根誤差(RMSE)和平均絕對百分比誤差(MAPE)作為衡量近似模型精度的評價指標(biāo)，依據(jù)第一批與第二批樣本點(diǎn)的Nemo軟件原值與BPNN計(jì)算值而得出的RMSE和MAPE分別為171.6和0.092.

由偏差曲線可以看出，該BPNN計(jì)算值能夠比較完好地接近于原值，特別是在10～21號樣本點(diǎn)之間，計(jì)算值與原值幾乎重合，在其余樣本點(diǎn)上該BPNN計(jì)算值也能較好地追蹤原值，并反映出總阻力的變化趨勢.BPNN的兩個評價指標(biāo)均較小，說明該BPNN對船體總阻力的預(yù)報精度較高.該球艏降阻優(yōu)化BPNN可以在滿足精度要求的前提下接替Nemo軟件對船體進(jìn)行總阻力預(yù)報，滿足成為下一步基于變精度模型的球艏降阻優(yōu)化研究中的低精度模型的要求.

6 球艏設(shè)計(jì)變量的局部優(yōu)化

6.1 第一批設(shè)計(jì)點(diǎn)的選取及尋優(yōu)

使用訓(xùn)練完成的球艏降阻優(yōu)化BPNN在整個設(shè)計(jì)空間預(yù)報第一批由拉丁超立方設(shè)計(jì)法選取的500個設(shè)計(jì)點(diǎn)總阻力.BPNN的總阻力預(yù)報值空間分布如圖21所示.

圖21 第一批設(shè)計(jì)點(diǎn)的總阻力空間分布

觀察圖21的散點(diǎn)在ZX平面和YZ平面的投影，可以判斷第一批設(shè)計(jì)點(diǎn)的設(shè)計(jì)變量與總阻力數(shù)據(jù)基本符合線性相關(guān)關(guān)系.

通過對3組變量進(jìn)行正態(tài)性檢驗(yàn)，得出3組變量正態(tài)性檢驗(yàn)顯著性均小于0.05，故不可直接采用Pearson相關(guān)分析，可使用Spearman相關(guān)系數(shù)進(jìn)行研究.Spearman相關(guān)系數(shù)通常被認(rèn)為是排列后變量之間的Pearson線性相關(guān)系數(shù)，其用在數(shù)據(jù)的分布使得Pearson線性相關(guān)系數(shù)不能用來描述或是用來描述將導(dǎo)致錯誤的結(jié)論時，評價兩個統(tǒng)計(jì)變量的相關(guān)性.

經(jīng)計(jì)算得到其相關(guān)分析結(jié)果：lb與Rt的相關(guān)系數(shù)為0.214，Sig值小于0.001，故在整個設(shè)計(jì)空間內(nèi)lb與Rt存在極顯著的弱正相關(guān)關(guān)系；hb與Rt的相關(guān)系數(shù)為0.921，Sig值小于0.001，故在整個設(shè)計(jì)空間內(nèi)hb與Rt存在極顯著的強(qiáng)正相關(guān)關(guān)系.

6.2 第二批設(shè)計(jì)點(diǎn)的選取及尋優(yōu)

依據(jù)整個設(shè)計(jì)空間內(nèi)lb、hb與Rt的相關(guān)分析結(jié)論，通過縮小在設(shè)計(jì)空間的尋優(yōu)范圍，在存在最優(yōu)解的局部設(shè)計(jì)空間：

3 693 mm≤lb≤4 500 mm

4 515 mm≤hb≤5 000 mm

(9)

使用拉丁超立方設(shè)計(jì)法選取第二批500個設(shè)計(jì)點(diǎn)，并使用球艏降阻優(yōu)化BPNN進(jìn)行阻力預(yù)報.第二批設(shè)計(jì)點(diǎn)總阻力預(yù)報值空間分布如圖22所示.

觀察圖22的散點(diǎn)在ZX平面和YZ平面的投影，可以判斷第二批設(shè)計(jì)點(diǎn)的設(shè)計(jì)變量與總阻力數(shù)據(jù)基本符合線性相關(guān)關(guān)系.通過對3組變量進(jìn)行正態(tài)性檢驗(yàn)，得出3組變量正態(tài)性檢驗(yàn)顯著性均小于0.05，故不可直接采用Pearson相關(guān)分析，可使用Spearman相關(guān)系數(shù)對第二批500個設(shè)計(jì)點(diǎn)變量進(jìn)行相關(guān)分析.

圖22 第二批設(shè)計(jì)點(diǎn)的總阻力空間分布

經(jīng)計(jì)算得到其相關(guān)分析結(jié)果：lb與Rt的相關(guān)系數(shù)為0.318，Sig值小于0.001，故在整個設(shè)計(jì)空間內(nèi)lb與Rt存在極顯著的弱正相關(guān)關(guān)系；hb與Rt的相關(guān)系數(shù)為0.784，Sig值小于0.001，故在整個設(shè)計(jì)空間內(nèi)hb與Rt存在極顯著的強(qiáng)正相關(guān)關(guān)系.該結(jié)論與整個設(shè)計(jì)空間的相關(guān)分析結(jié)論一致.

比較各設(shè)計(jì)點(diǎn)及樣本點(diǎn)的總阻力，得出總阻力最小值點(diǎn)(lb=3 770 mm，hb=4 638 mm)，該點(diǎn)的總阻力為144 257.6 N，該點(diǎn)是局部優(yōu)化中逼近最優(yōu)解的設(shè)計(jì)點(diǎn).優(yōu)化前后的球艏幾何形狀對比如圖23所示.

圖23 優(yōu)化前后球艏變形示意圖

7 球艏降阻優(yōu)化結(jié)果分析

對優(yōu)化結(jié)果按船舶阻力理論結(jié)合Holtrop法原理進(jìn)行分析.在Holtrop法各阻力成分中，試驗(yàn)油船優(yōu)化前后變化最大的分別是興波阻力Rw、摩擦阻力Rf及球艏附加壓阻力Rb.優(yōu)化后該油船的興波阻力發(fā)生了較明顯的減小，摩擦阻力及球艏附加壓阻力發(fā)生了增大，由于Rw的減少量大于其他阻力成分增量的總和，故總阻力Rt減小，達(dá)到了提高船舶快速性的目的.由Holtrop法計(jì)算的優(yōu)化前后兩種球艏的總阻力及各阻力成分?jǐn)?shù)值對照如表7所示.

表7 優(yōu)化前后球艏阻力對比

該試驗(yàn)油船摩擦阻力在優(yōu)化后發(fā)生了增大.優(yōu)化后試驗(yàn)油船排水體積為6 880.658 m3，較試驗(yàn)油船基準(zhǔn)設(shè)計(jì)模型的排水體積增加了 84.587 m3，根據(jù)船體濕表面積Shull的估算公式[10]

(10)

(11)

得出：由于B/T和舯橫剖面系數(shù)Cm在優(yōu)化前后未發(fā)生改變，故Cs視為常數(shù)，由于該試驗(yàn)油船在優(yōu)化前后設(shè)計(jì)水線長的改變較小，故可視Lwl未改變.且排水體積較優(yōu)化前增大，故Shull增大.根據(jù)

Rf=0.037 5ρv2Shull/(lgRe-2)2

(12)

(13)

得出：由于水密度ρ、航速v及設(shè)計(jì)水線長Lwl在球艏優(yōu)化前后未發(fā)生改變，水的運(yùn)動黏性系數(shù)ν≈1.188×10-6，故摩擦阻力Rf與Shull呈正比例關(guān)系[18].優(yōu)化后的Rf較優(yōu)化前增大0.349%.

優(yōu)化后球艏幾何形狀的改變使得接近水表面處的球艏附加壓阻力Rb發(fā)生較大幅度增加，根據(jù)球艏附加壓阻力Rb的Holtrop法計(jì)算公式[15]

(14)

(15)

(16)

得出：在航速v和艏吃水Tf為定值時，Rb僅與艏垂線處球艏橫剖面面積Abt及其形心到基線的高度Hb有關(guān).優(yōu)化后Abt與Hb均較優(yōu)化前增大，經(jīng)過船舶設(shè)計(jì)軟件Nemo計(jì)算，得出優(yōu)化后Rb較優(yōu)化前增大1 363.5%.圖24分別是艏垂線處球艏橫剖面在優(yōu)化前及優(yōu)化后的Abt及Hb的變化情況示意圖.

船艏在航行時破波會產(chǎn)生破波阻力，破波阻力本質(zhì)上是一種興波阻力，若采用尖瘦的水線和尖劈型的橫剖面組成的球艏可減緩波浪破碎，從而減少破波能量的虧損.觀察圖24中橫剖面的形狀，優(yōu)化后的艏垂線處球艏橫剖面較優(yōu)化前更為瘦削，故優(yōu)化后的球艏可在一定程度上減少破波阻力[18].

圖24 優(yōu)化前后艏垂線處球艏橫剖面參數(shù)變化示意圖

該試驗(yàn)油船的興波阻力在優(yōu)化后發(fā)生了較明顯的減小.由于設(shè)計(jì)水線長Lwl在優(yōu)化前后未改變，故弗勞德數(shù)Fr=0.221在優(yōu)化前后未改變，依據(jù)Holtrop法在0

Rw=C1C2C5ρgexp[m1Frd+m4cos(λFr-2)]

(17)

改變設(shè)計(jì)變量lb、hb的數(shù)值，進(jìn)而引起Abt、Hb及的變化是該油船總阻力Rt改變的主要原因.且優(yōu)化后的船舶排水量較優(yōu)化前增加了1.245%，故該試驗(yàn)油船的船舶載貨能力也將增大.

8 結(jié) 語

本文以油船球艏降阻優(yōu)化為目標(biāo)構(gòu)建了基于BPNN的球艏降阻優(yōu)化模型，該近似模型可以在保證精度的前提下快速預(yù)報大量設(shè)計(jì)點(diǎn)的總阻力.模型的構(gòu)建為下一步基于變精度模型的研究中神經(jīng)網(wǎng)絡(luò)的應(yīng)用以及進(jìn)行阻力、排水量等多目標(biāo)優(yōu)化提供了經(jīng)驗(yàn)與支持.

在構(gòu)建過程中得出兩個結(jié)論：lb與Rt存在極顯著的弱正相關(guān)關(guān)系；hb與Rt存在極顯著的強(qiáng)正相關(guān)關(guān)系.可為同類型船舶球艏形狀優(yōu)化提供新思路與參考.