，， ，

(哈爾濱工程大學(xué) 船舶工程學(xué)院，哈爾濱 150001)

現(xiàn)在，常用的三維水動(dòng)力計(jì)算方法主要有自由面格林函數(shù)法和Rankine源法。在應(yīng)用這兩種方法時(shí)都需要對(duì)船體濕表面進(jìn)行網(wǎng)格劃分，而在Rankine源方法中還需要對(duì)自由表面進(jìn)行網(wǎng)格劃分。目前，較為常用的方法是常值面元法，該理論優(yōu)點(diǎn)是便于操作。但缺點(diǎn)是對(duì)于復(fù)雜結(jié)構(gòu)物或在幾何形體存在大曲率變化的地方，需要布置大量的網(wǎng)格，以保證物面的“連續(xù)性”，這勢(shì)必會(huì)增加計(jì)算量。因此，便誕生了一些高階面元法，比如B樣條面元法。B樣條面元法中所用到的插值多項(xiàng)式涉及附近的幾個(gè)單元，而不像一般高階邊界元方法只限制在一個(gè)單元內(nèi)，這樣速度勢(shì)函數(shù)及其導(dǎo)數(shù)在整個(gè)物面上就是連續(xù)的。此外，由于B樣條面元法中相鄰單元的速度勢(shì)控制點(diǎn)是公用的，因此未知數(shù)的個(gè)數(shù)相對(duì)較少，計(jì)算需要的面元數(shù)越多，B樣條面元法計(jì)算優(yōu)勢(shì)越明顯。如何用B樣條方便、快捷地劃分船體水動(dòng)力網(wǎng)格一直是大家關(guān)注的問(wèn)題[1]。為此，提出一種基于B樣條的船體網(wǎng)格劃分及水動(dòng)力計(jì)算方法。

1 三次B樣條曲線和曲面

1.1 插值三次B樣條曲線

為了能描述復(fù)雜形狀和具有局部性質(zhì)，B樣條曲線用一套特殊的基函數(shù)即B樣條基函數(shù)來(lái)表示。這樣，三次B樣條曲線[2]方程可寫(xiě)為

(1)

給定一組數(shù)據(jù)點(diǎn)qi(i=0,1,…,m)，由插值條件給出以n+1個(gè)控制頂點(diǎn)為未知矢量的m+1個(gè)線性方程組成的線性方程組

(2)

將上述線性方程組改寫(xiě)成如下矩陣形式

(3)

解該線性方程組，即可求出全部未知控制頂點(diǎn)，然后利用德布爾遞推插值出三次B樣條曲線，遞推公式為

(4)

1.2 插值雙三次B樣條曲面

B樣條曲面方程比曲線多了一個(gè)參數(shù)v， 雙三次B樣條插值曲面方程寫(xiě)為

(5)

又可改寫(xiě)為

給出類(lèi)似于B樣條曲線方程的表達(dá)式

(6)

(7)

雙三次B樣條曲面可以由兩個(gè)方向的一系列關(guān)聯(lián)的三次B樣條曲線來(lái)表示。

2 水動(dòng)力網(wǎng)格劃分

在應(yīng)用自由面格林函數(shù)法進(jìn)行水動(dòng)力計(jì)算時(shí)，雖然只需在船體濕表面劃分網(wǎng)格[3-4]，但是自由面格林函數(shù)的精確計(jì)算仍然是個(gè)難點(diǎn)。因此，本文運(yùn)用B樣條方法編制相應(yīng)的FORTRAN程序?qū)Υw濕表面網(wǎng)格進(jìn)行劃分，然后采用Rankine源法求解船體水動(dòng)力，船體的水動(dòng)力網(wǎng)格包括船體濕表面和其周?chē)淖杂擅婢W(wǎng)格。

2.1 船體濕表面網(wǎng)格劃分

首先對(duì)船體型值點(diǎn)信息進(jìn)行適當(dāng)處理，成為反算船體表面控制頂點(diǎn)所需要的數(shù)據(jù)點(diǎn)，然后運(yùn)用編制的三次B樣條FORTRAN程序生成計(jì)算船體水動(dòng)力時(shí)所需的網(wǎng)格數(shù)據(jù)信息，包括船體濕表面控制頂點(diǎn)和U、V參數(shù)值。

2.2 自由面網(wǎng)格劃分

由于在計(jì)算船體水動(dòng)力時(shí)要滿足自由面和遠(yuǎn)方輻射條件，所以應(yīng)使船體周?chē)淖杂擅娼苓h(yuǎn)疏，以便較好地模擬船體周?chē)牧鲌?chǎng)。本文采用“拉伸變換”[5]來(lái)輔助船體周?chē)杂擅婢W(wǎng)格的生成，拉伸變化公式為

(8)

圖1、2為Wigley船型、某拖船和它們周?chē)杂擅婢W(wǎng)格劃分實(shí)例。

3 實(shí)船水動(dòng)力計(jì)算

對(duì)一艘50萬(wàn)t，船體艏艉曲率變化較大的大型礦砂船采用編制的FORTRAN程序，進(jìn)行網(wǎng)格劃分及水動(dòng)力計(jì)算，圖3為50萬(wàn)t礦砂船周?chē)杂擅婢W(wǎng)格。其中包括船體廣義附加質(zhì)量、廣義附加阻尼、入射波波浪激勵(lì)力、繞射波波浪激勵(lì)力、剖面載荷和運(yùn)動(dòng)響應(yīng)(見(jiàn)圖4～9)的求解，并將結(jié)果與文獻(xiàn)[6]進(jìn)行了對(duì)比。

圖1 Wigley船型周?chē)杂擅婢W(wǎng)格(網(wǎng)格總數(shù)2 320)

圖2 某拖船周?chē)杂擅婢W(wǎng)格(網(wǎng)格總數(shù)1 320)

圖3 50萬(wàn)t礦砂船周?chē)杂擅婢W(wǎng)格(網(wǎng)格總數(shù)1 722)

圖4 無(wú)因次附加質(zhì)量隨ω變化

圖5 無(wú)因次阻尼力隨ω變化

圖6 入射波波浪力F1隨ω變化

圖7 繞射波波浪力Fd隨ω變化

圖8 船舯剖面載荷隨ω變化

從求解的水動(dòng)力結(jié)果可以看出，本文的結(jié)果與文獻(xiàn)[6]的結(jié)果非常吻合，由此可見(jiàn)，用B樣條方法進(jìn)行船體水動(dòng)力網(wǎng)格的劃分，并用于實(shí)船水動(dòng)力計(jì)算是可行的。

圖9 運(yùn)動(dòng)響應(yīng)η1隨ω變化

4 結(jié)論

B樣條面元法是一種高階面元法，它不像一般高階邊界元方法那樣在物體表面以及物體表面的速度勢(shì)都不連續(xù)，在單元之間有跳躍以及基本源的高階導(dǎo)數(shù)存在奇異性。本文采用的雙三次B樣條面元法保證了物體表面的光滑性，速度勢(shì)函數(shù)及其一階、二階導(dǎo)數(shù)在整個(gè)物面上也是連續(xù)的，從而提高了水動(dòng)力計(jì)算結(jié)果的準(zhǔn)確性和該方法的適用性。

[1] PIEGL L, TILLER W. The NURBS Book[M]. New York:Springer，1996.

[2] 施法中.計(jì)算機(jī)輔助幾何設(shè)計(jì)與非均勻有理B樣條[M].北京:高等教育出版社,2001.

[3] 張海彬，任慧龍，宋競(jìng)正，等.一種船體三維濕表面網(wǎng)格自動(dòng)生成方法[J].中國(guó)造船，2001，42(4)：61-65.

[4] 陳賓康.用B樣條擬合與繪制船體型線[J].武漢水運(yùn)工程學(xué)院學(xué)報(bào),1989,13(4):1-10.

[5] 戴愚志，余建星.一種船體及周?chē)杂擅娴木W(wǎng)格自動(dòng)生成方法[J].船舶工程,2006,28(5):1-4.

[6] 劉日明.基于B樣條面元法的浮體二階水動(dòng)力計(jì)算[D].哈爾濱:哈爾濱工程大學(xué),2009.