基于罰函數(shù)的滑行艇阻力性能優(yōu)化方法研究*

許蘊(yùn)蕾

(海軍駐上海地區(qū)艦艇設(shè)計(jì)研究軍事代表室 上海 200011)

基于罰函數(shù)的滑行艇阻力性能優(yōu)化方法研究*

許蘊(yùn)蕾

(海軍駐上海地區(qū)艦艇設(shè)計(jì)研究軍事代表室 上海 200011)

滑行艇;阻力;罰函數(shù);優(yōu)化方法

滑行艇的阻力性能優(yōu)化是設(shè)計(jì)滑行艇艇體的重要內(nèi)容之一,如何減小艇體阻力是設(shè)計(jì)師優(yōu)先考慮的目標(biāo)。文章從傳統(tǒng)的SIT阻力估算方法著手,分析滑行艇受到的力和力矩,并引入Savitsky對(duì)艇體阻力的修正,即考慮因噴濺而產(chǎn)生的摩擦阻力。由于滑行艇的重心縱向位置對(duì)阻力性能有很大影響,故通過(guò)構(gòu)造罰函數(shù)法,把有約束問(wèn)題化為無(wú)約束問(wèn)題,利用MATLAB優(yōu)化工具箱中的fm incon函數(shù)計(jì)算出最佳的重心縱向位置。最后,通過(guò)某滑行艇模型的水池拖曳試驗(yàn)驗(yàn)證了該優(yōu)化方法的可行性,具有一定的工程參考價(jià)值。

0 引 言

滑行艇水動(dòng)力性能的研究歷來(lái)是研究的難點(diǎn),其中重要的一項(xiàng)內(nèi)容是對(duì)阻力性能的研究。阻力是船舶最重要的性能之一,如果能從理論上設(shè)計(jì)出最小阻力船型,這將為整個(gè)船舶設(shè)計(jì)周期節(jié)省大量的勞力,時(shí)間和費(fèi)用。隨著“時(shí)間”價(jià)值的不斷提高,人們對(duì)各種運(yùn)輸工具速度的要求也日益提高。對(duì)民用船舶來(lái)說(shuō),航速的提高意味著貨物周轉(zhuǎn)周期的縮短以及經(jīng)濟(jì)效益的提高;對(duì)軍用艦艇來(lái)說(shuō)則意味著戰(zhàn)斗力的增強(qiáng)[1]。

滑行艇阻力性能優(yōu)化主要包括艇體型線優(yōu)化、艇體參數(shù)優(yōu)化和添加附體等。艇體型線優(yōu)化通常可以應(yīng)用數(shù)學(xué)方法對(duì)型線進(jìn)行光順,但是必須以艇體的布置、水動(dòng)力與結(jié)構(gòu)性能的要求為目標(biāo)函數(shù)。在以艇體阻力性能為型線優(yōu)化對(duì)象時(shí),只能在某些約束條件下完成[2]。海軍工程大學(xué)船舶與海洋工程系與水動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室,通過(guò)對(duì)三種艇型及不同噴氣方式的模型試驗(yàn),研究了斷階滑行艇模型氣層減阻的實(shí)施途徑及減阻效果,取得了總阻力減少25%以上的結(jié)果,提出了一種適合于采用氣層減阻技術(shù)且阻力性能優(yōu)良的艇型[3]。崔焰、盧曉平對(duì)方尾水面艦船,采用Hohrop阻力計(jì)算公式,以總阻力為目標(biāo)函數(shù),據(jù)遺傳算法,對(duì)主要的船型參數(shù)進(jìn)行了優(yōu)化計(jì)算[4]?；型У闹匦目v向位置對(duì)阻力性能有很大影響。從減小阻力考慮,重心后移是有利的,它對(duì)應(yīng)較大的有利寬度,使滑行效率提高;但對(duì)避免海豚運(yùn)動(dòng)發(fā)生和波浪中的運(yùn)動(dòng)響應(yīng)和沖擊力都會(huì)帶來(lái)不利影響。因此,如何確定艇體的重心縱向位置是設(shè)計(jì)滑行艇阻力性能的重要內(nèi)容之一。

1 SIT阻力估算方法

SIT法又稱(chēng)Savitsky法[5],它是根據(jù)美國(guó)史蒂文斯實(shí)驗(yàn)室水池試驗(yàn)結(jié)果提出的。自1964年Savitsky提出該方法以來(lái),已有幾種不同的版本,即使Savitsky版也是Krovkin-Kroukovsky和Murray等人研究工作的延續(xù)和修改。以多種形式出現(xiàn)的SIT法可能是在美國(guó)應(yīng)用得最多的滑行艇阻力預(yù)報(bào)方法。不同的版本常給出不同的結(jié)果,不過(guò),它們的核心卻是相同的。該方法認(rèn)為動(dòng)浮力與攻角不應(yīng)該是線性關(guān)系,而是經(jīng)滑行參數(shù)表示的升力、扭矩方程(以壓力中心函數(shù)的形式)。

1.1 Savitsky的力和力矩方程

滑行艇的阻力包括壓差阻力Rp(或稱(chēng)剩余阻力)和摩擦阻力Rf兩大部分。壓差阻力又包括飛濺阻力和興波阻力,它們分別等于動(dòng)水力壓力合力P和靜水壓力P0在運(yùn)動(dòng)方向的分量;而P和P0在鉛錘方向的分量之和等于艇重。因此可知壓差阻力與艇重的關(guān)系為:

滑行艇的阻力:

一般地,傾角較小(cosτ≈1),艇底來(lái)流速度近似地取等于航速,這樣R為傾角τ和濕長(zhǎng)度l的函數(shù),欲確定τ和l,可以通過(guò)滑行艇運(yùn)動(dòng)時(shí)的受力及力矩平衡方程來(lái)求解。

不計(jì)推力與摩擦力的鉛垂分量及這些力對(duì)尾緣點(diǎn)的力矩,并考慮到小攻角的條件,得出鉛垂方向力的平衡方程和對(duì)尾緣點(diǎn)力矩的平衡方程:

式中:lG——重心至尾緣距離(沿龍骨方向);

Mp——?jiǎng)訅毫狭?duì)尾緣的力矩;

Mp0——靜壓力合力對(duì)尾緣的力矩。

SIT法的無(wú)因次力平衡方程

式中:第一項(xiàng)為無(wú)因次動(dòng)浮力項(xiàng),第二項(xiàng)為無(wú)因次靜浮力項(xiàng);傾角的單位為“°”。當(dāng)有橫傾角時(shí),需要進(jìn)行以下修正:

SIT法的無(wú)因次力矩平衡方程,以壓力中心的方式表達(dá):

式中:

為求解λ和τ,可采用如下計(jì)算步驟:

1)已知滑行艇排水量Δ、航速V、艇寬b和斜升角β,計(jì)算寬度傅汝德數(shù)Fb及艇的實(shí)際負(fù)載系數(shù)

由(5)式可以求得相當(dāng)平板的動(dòng)負(fù)荷系數(shù)Cb0。

2)令lp=lG,則(6)式可轉(zhuǎn)化為:

3)聯(lián)立式(4)、(8),采用擬牛頓法編程求解這一非線性方程組。

1.2 噴濺阻力

1964年的Savitsky版本中未考慮噴濺引起的濕面積對(duì)摩擦阻力的影響。2006年,Savitsky重新修改了1964年版本,使其能夠計(jì)及噴濺引起的濕面積所產(chǎn)生的摩擦阻力。濕面積由兩部分組成:濕底壓力區(qū)和噴濺區(qū)?；型?dòng)升力完全由濕底壓力區(qū)提供,噴濺區(qū)只產(chǎn)生摩擦阻力[7]。

當(dāng)滑行艇高速滑行時(shí),由于受到流體動(dòng)力作用,出現(xiàn)有明顯的飛濺現(xiàn)象,艇底流動(dòng)模型如圖1所示。在兩舷側(cè)形成一種連續(xù)噴射的圓錐面水膜,稱(chēng)為膜狀飛濺或主飛濺;此外,在駐點(diǎn)線前還有向兩舷外側(cè)噴射的,由很細(xì)小的水滴組成的水流束,稱(chēng)為須狀飛濺。這兩種飛濺損耗的能量就產(chǎn)生飛濺阻力,而須狀飛濺還將增加艇底浸濕面積,以致對(duì)摩擦阻力有所影響。由于噴濺速度難以測(cè)量,噴濺層還不同于常規(guī)邊界層,常規(guī)邊界層不能用于準(zhǔn)確求解噴濺對(duì)摩擦阻力的影響,這造成噴濺對(duì)總阻力的影響難以準(zhǔn)確確定。

圖1 艇底流動(dòng)模型

在滑行狀態(tài)下,噴濺基線同駐點(diǎn)線的縱向距離非常小,如圖1所示,在噴濺基線處的水流速度沿駐點(diǎn)線的方向。目前,估算船體總粘性阻力時(shí),以噴濺基線后的面積作為總的濕面積。

圖1中,被駐點(diǎn)線OC、浸濕龍骨(長(zhǎng)度為L(zhǎng)k)、兩條浸濕舭緣線(長(zhǎng)度為L(zhǎng)C),以及尾端圍成的兩塊區(qū)域?yàn)榈撞繅毫^(qū)域Ap、As為噴濺區(qū)域面積,噴濺區(qū)域的大小同斜升角β和縱傾角τ有關(guān)。對(duì)于棱柱型船,Savitsky給出了浸濕長(zhǎng)度同斜升角以及縱傾角之間的關(guān)系:

Bowles和Denny給出了預(yù)測(cè)滑行艇艏部附近的水表面擾動(dòng)模型公式:

因水動(dòng)壓力產(chǎn)生的噴濺,同船體構(gòu)型(如舭部升高、縱傾以及舭部濕/干的區(qū)域)和工作條件(如滑行艇的速度和所受的自由表面波浪)和成比例關(guān)系。因此,采用經(jīng)驗(yàn)公式來(lái)確定噴濺阻力:

其中K1由滑行艇的速度決定,它有以下關(guān)系:

噴濺處的壓力由下面的公式進(jìn)行計(jì)算:

其中P為噴濺基線處的壓力,系數(shù)K2可近似確定為:

其中b(y)是噴濺處距縱中剖面的距離,如圖2所示:

圖2 噴濺基線處壓力系數(shù)的確定

因而噴濺所產(chǎn)生的噴濺阻力和升力可用以下公式進(jìn)行計(jì)算:

1.3 總阻力

考慮到噴濺引起的濕面積對(duì)摩擦阻力的影響,滑行艇實(shí)際受到的阻力應(yīng)再加上噴濺引起的摩擦阻力Rsp,則滑行艇總阻力表示為:

(15)式中的排水量通常已知,而攻角τ未知。在摩擦阻力Rf中,壓力區(qū)浸濕面積未知,它取決于浸濕長(zhǎng)度;噴濺引起的摩擦阻力Rsp也與浸濕長(zhǎng)度有關(guān),而浸濕長(zhǎng)度是隨艇體排水量、航速、艇寬、重心縱向位置和斜升角的變化而變化的。

一般來(lái)說(shuō),根據(jù)任務(wù)書(shū)來(lái)設(shè)計(jì)滑行艇的阻力性能,其艇體排水量、艇寬、最高航速、巡航速度通常都已經(jīng)確定,可以改變的是重心縱向位置和斜升角。而如何確定最佳的重心縱向位置和斜升角,以達(dá)到滑行艇在最高航速和巡航速度具有較小的功率,則是本文的研究目標(biāo)。

2 罰函數(shù)的艇體阻力優(yōu)化方法

2.1 基本原理

罰函數(shù)法即序列無(wú)約束極小化方法,它的基本原理是將有約束問(wèn)題化為無(wú)約束問(wèn)題,亦即將原來(lái)的目標(biāo)函數(shù)和約束函數(shù),按一定方式構(gòu)成一個(gè)新的函數(shù),當(dāng)這個(gè)新的函數(shù)向函數(shù)逼近,它的最優(yōu)解也就是原問(wèn)題的最優(yōu)解。罰函數(shù)法又分為:內(nèi)點(diǎn)罰函數(shù)法,外點(diǎn)罰函數(shù)法和混合罰函數(shù)法[6]。

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

在艇型優(yōu)化之前,需要明確滑行艇阻力性能設(shè)計(jì)的目標(biāo),即在巡航速度和最高航速時(shí),阻力較小。從發(fā)動(dòng)機(jī)有效功率來(lái)看,優(yōu)化目標(biāo)可以定義為:滑行艇在最高速度航行時(shí)的有效功率與巡航速度時(shí)的有效功率之比最小,并且以最高速度航行時(shí)的有效功率盡可能小。

2.3 設(shè)計(jì)變量

如前所述,影響滑行艇阻力的幾個(gè)主要要素為:排水量Δ、艇寬B、航速V、重心縱向位置lG和斜升角β。由于艇寬和排水量受到設(shè)計(jì)指標(biāo)的限定,這里均設(shè)為定值。

2.4 目標(biāo)函數(shù)

由修正后的SIT預(yù)報(bào)模型(15)式可知:

假設(shè)滑行艇巡航速度和最高航速分別為V1、V2,則優(yōu)化的目標(biāo)函數(shù)為:

引入中間變量τ1,τ2,λ1,λ2,上式可寫(xiě)作:

2.5 約束函數(shù)

(1)由SIT法的無(wú)因次力平衡方程得

(2)對(duì)于相當(dāng)平板的動(dòng)負(fù)荷系數(shù)Cb0有

(3)由SIT法的無(wú)因次力矩平衡方程得

(4)當(dāng)攻角＞12°時(shí),已非滑行狀態(tài)了,故有

(5)斜升角≤30°

(6)重心位置應(yīng)處在艇中之后,以艇中靠后3%～4%艇長(zhǎng)L(折角線長(zhǎng))為最佳[6]。

在(1)、(2)、(3)、(5)之條件下,和約束條件完全相同,故一共有15組。

2.6 優(yōu)化計(jì)算

本文滑行艇的阻力性能優(yōu)化模型屬于非線性約束的最優(yōu)化問(wèn)題,采用SQP算法。因此,這里可以使用MATLAB優(yōu)化工具箱中的fmincon函數(shù)[8]。MATLAB中SQP法的實(shí)現(xiàn)分3步,即:

(1)拉格朗日函數(shù)Hess矩陣的更新;

(2)二次規(guī)劃問(wèn)題求解;

(3)一維搜索和目標(biāo)函數(shù)的計(jì)算。

3 基于模型艇的計(jì)算實(shí)例

根據(jù)上述罰函數(shù)的滑行艇阻力優(yōu)化方法,對(duì)某滑行艇模型進(jìn)行優(yōu)化計(jì)算,其模型參數(shù)如表1所示,通過(guò)母型船換算得到的巡航速度和最高航速分別為6.0 m/s、10.0m/s。首先,在MATLAB程序中是建立目標(biāo)函數(shù)的m文件,即編寫(xiě)目標(biāo)函數(shù)程序并保存在objfun.m文件中;然后建立非線性約束條件的m文件并保存在confun.m文件中;最后編寫(xiě)優(yōu)化程序。程序執(zhí)行過(guò)程是通過(guò)優(yōu)化程序中的fmincon函數(shù)調(diào)用目標(biāo)函數(shù)和約束條件函數(shù),并用表1中的各種參數(shù)作為優(yōu)化初始值實(shí)現(xiàn)優(yōu)化運(yùn)算。通過(guò)MATLAB程序優(yōu)化模型的計(jì)算,可得到最佳的重心縱向位置lG=0.71m,即重心距船舯后0.39m。

表1 某滑行艇船模參數(shù)

3.1 模型的水池拖曳試驗(yàn)驗(yàn)證

為了確定艇體各種參數(shù)對(duì)滑行艇阻力性能的影響,在某所高速水動(dòng)力試驗(yàn)室進(jìn)行了船模的水池拖曳試驗(yàn),通過(guò)測(cè)試該模型不同的重量、重心位置、附體等情況下的阻力、升沉、縱傾角等參數(shù),從而為滑行艇的艇型改進(jìn)和實(shí)艇的設(shè)計(jì)提供技術(shù)依據(jù)。該模型為木質(zhì),表面打磨光滑、噴漆,經(jīng)檢驗(yàn)符合“滑行艇船模阻力試驗(yàn)方法(CB/Z 244-88)”中的允差標(biāo)準(zhǔn)。模型主要參數(shù)如表1所示。

為了確定艇體重心位置對(duì)滑行艇阻力性能的影響,分別進(jìn)行重心距船舯后350 mm、370 mm、390 mm的水池拖曳試驗(yàn),對(duì)測(cè)量的阻力制成圖表如圖3所示:

由圖3可知,當(dāng)Vm＜5.0m/s時(shí),艇體重心距船舯后Xgm=350 mm的模型阻力較小;當(dāng)Vm＞5.0 m/s以后,Xgm=390 mm的模型阻力較小,阻力減小最大處達(dá)2.78%左右。這說(shuō)明在相同排水量下,滑行艇開(kāi)始滑行后,艇體重心后移將使縱傾角增大、濕表面積減小,阻力減小。

圖3 不同重心縱向位置的阻力變化曲線

4 結(jié) 語(yǔ)

本文根據(jù)Savitsky的阻力估算方法分析了滑行艇受到的力和力矩,并在總阻力中計(jì)入噴濺引起的摩擦阻力;通過(guò)利用罰函數(shù)的優(yōu)化方法,確定了阻力優(yōu)化的目標(biāo),設(shè)計(jì)了優(yōu)化變量、目標(biāo)函數(shù)和約束函數(shù),然后利用MATLAB優(yōu)化工具箱中的fmincon函數(shù)對(duì)艇體阻力進(jìn)行優(yōu)化計(jì)算,得出最佳的艇體重心縱向位置,并通過(guò)某滑行艇的模型水池拖拽試驗(yàn)驗(yàn)證了該方法的可行性。

[1] 董祖舜.快艇動(dòng)力學(xué)[M].武漢:華中理工大學(xué)出版社,1994.

[2] 黃金森.基于最小阻力的船型優(yōu)化研究[D].天津大學(xué)工學(xué)碩士學(xué)位論文,2005.12.

[3] 董文才,郭日修,陳小玲,呂巖松.滑行艇氣層減阻試驗(yàn)[J].中國(guó)造船,2002,43(4).

[4] 崔焰,盧曉平.基于Holtrop公式與遺傳算法的船型參數(shù)優(yōu)化[J].海軍工程大學(xué)學(xué)報(bào),2009,21(6):26-29.[5] 董文才,郭日修.滑行艇阻力研究進(jìn)展[J].船舶力學(xué),2000,4(4):68-81.

[6] 朱凱.滑行艇的型線優(yōu)化設(shè)計(jì)及性能分析[D].哈爾濱工程大學(xué)工學(xué)碩士學(xué)位論文,2007.

[7] 孫華偉.三體滑行艇船型與快速性研究[D].哈爾濱:哈爾濱工程大學(xué)碩士論文,2010.1

[8] 褚洪生,杜增吉,閻金華等編著.MATLAB7.2優(yōu)化設(shè)計(jì)實(shí)例指導(dǎo)教程[M].機(jī)械工業(yè)出版社,2007:91-119.

Optimization on Resistance Capability of Planing Craft Based on Penalty Function

Xu Yunlei

planing craft;resistance;penalty function;optimization

Resistance capability optimization is one of the important issues to design a planing craft hull,and the top-priority aim for designer is how to reduce the hull resistance.Firstly,the SIT resistance estimation method is introduced,and the force and moment on planning craft are analyzed,then the Savitsky modification of hull resistance is proposed,namely the frictional resistance caused by spray.Due to the great effect of the barycenter position to resistance capability,the penalty function is implemented to change the constrained problem to unconstrained problem.The best barycenter position is calculated by the fmincon function in the MATLAB optimized kit.Finally,the feasibility of the optimization method is validated by the tank daggle trial of a planing craft model.

U 661.33

A

1001-9855(2010)05-0009-05

2010-06-25

許蘊(yùn)蕾(1981.06-),女,漢族,吉林人,助理工程師,主要從事船舶設(shè)計(jì)審查工作。