王 中,盧曉平,王 瑋

(1.海軍工程大學(xué) 船舶與動力學(xué)院,湖北 武漢 430033;2.海軍裝備研究院 艦船所,北京 100073)

采用三維Rankine源方法計算船舶興波阻力時,自由面網(wǎng)格的生成是一個必要步驟,也是一個關(guān)鍵步驟[1].近年來三體船成為國內(nèi)外研究的熱點[2-4].由于側(cè)體的存在,自由面網(wǎng)格的劃分比單體船更加復(fù)雜、繁瑣.Dawson型方法自由面采用流線網(wǎng)格,三體船側(cè)體附近的流線需要反復(fù)迭代計算、調(diào)整才能確定,費時費力;另外,更為合理的自由面條件不能完整地用流線的導(dǎo)數(shù)表示,如改進的Dawson法自由面條件,包括Rankine源方法各種階次非線性自由面條件以及全非線性自由面條件.在這一類興波阻力計算方法中,舍棄流線網(wǎng)格,采用自由面上的“貼水線”網(wǎng)格勢在必行[5-6].但對三體船來說,快速生成自由面上的“貼水線”網(wǎng)格也并非易事.側(cè)體以及方尾的存在使得網(wǎng)格生成變得相當(dāng)復(fù)雜,采用商用軟件(如Gambit、ICEM等)來劃分此類網(wǎng)格固然可以解決問題,但過程繁瑣,此類軟件生成的網(wǎng)格文件必須經(jīng)過進一步的解析轉(zhuǎn)換才能供自己編制的興波阻力計算程序使用.因此,針對三體船構(gòu)型特點開發(fā)適用于三體船的自由面網(wǎng)格快速生成系統(tǒng),成為提高三體船興波問題計算效率的重要技術(shù)環(huán)節(jié).

1 網(wǎng)格生成方法概述

實用的自由面網(wǎng)格快速生成系統(tǒng)應(yīng)符合輸入簡單、修改靈活、存儲方便、過程可視等要求,自由面網(wǎng)格快速生成系統(tǒng)界面如圖 1所示,具備上述特點.通過該系統(tǒng),可將整個輸入?yún)?shù)及網(wǎng)格劃分結(jié)果保存于數(shù)據(jù)文件中,需要修改時直接打開該文件,改變一些參數(shù)即可快速重新生成需要的網(wǎng)格;同時,可將生成的網(wǎng)格數(shù)據(jù)導(dǎo)出成文本文件,直接供其他程序使用,導(dǎo)出的網(wǎng)格文件分塊存儲了主區(qū)域、主體方尾后側(cè)區(qū)域、側(cè)體方尾后側(cè)區(qū)域 3部分網(wǎng)格頂點坐標(biāo)信息.

在網(wǎng)格劃分前首先需要確定必要的信息.按如下方式定義系統(tǒng)的輸入?yún)?shù),如圖 1、2所示:主體靜水面水線型值點、側(cè)體靜水面水線型值點(側(cè)體局部坐標(biāo)系)及側(cè)體橫向偏距Dy和縱向偏距Dx、上游長度P1P12、下游長度P13P2、橫向?qū)挾?P1P11、傾角∠P11P1P8、縱向單位船長網(wǎng)格數(shù)、橫向總網(wǎng)格數(shù)、主體與側(cè)體間橫向網(wǎng)格數(shù)、橫向網(wǎng)格密度控制因子、主體方尾橫向網(wǎng)格數(shù)、側(cè)體方尾橫向網(wǎng)格數(shù),其中上游長度、下游長度、橫向長度以主體長度 L為單位.

圖1 軟件系統(tǒng)主界面Fig.1 Window of the software

圖2 自由面網(wǎng)格參數(shù)定義Fig.2 Free-surfacemesh parameter definition

如圖 2所示,由輸入的網(wǎng)格控制參數(shù)可以確定邊界點:P1、P2、P3、P4、P5、P6、P7、P8、P9、P10.這些邊界點與輸入的船體水線型值點可生成供插值的NURBS縱向邊界線:l1、l2、l3、l4、l5、l6、l7、l8、l9,參見圖2.這些邊界線是由一條或多條NURBS曲線組成的,不論是 2點組成的線段還是多點組成的曲線,均統(tǒng)一寫成NURBS曲線形式,邊界線可以將自由面劃分成 S1～S5共 5個子區(qū)域,然后根據(jù)定義的縱向和橫向網(wǎng)格數(shù)在各個子區(qū)域分別劃分網(wǎng)格.NURBS曲線表示三體船自由面網(wǎng)格(子)區(qū)域的邊界線是本文貼體坐標(biāo)代數(shù)法生成自由面網(wǎng)格的特點之一,這種方法便于自由面網(wǎng)格邊界線分割,從而較好地控制網(wǎng)格線密度.

一旦按上述方法確定了自由面網(wǎng)格的邊界線,即可以分區(qū)采用貼體坐標(biāo)代數(shù)法生成自由面網(wǎng)格.自由面橫向網(wǎng)格和縱向網(wǎng)格密度的控制對興波阻力Rankine源方法計算結(jié)果的影響不可忽視,在貼體坐標(biāo)代數(shù)法生成自由面網(wǎng)格過程中可采用各種方法控制網(wǎng)格密度,本文縱向網(wǎng)格(即橫向網(wǎng)格線)取為等間距的,重點考慮了橫向網(wǎng)格(即縱向網(wǎng)格線)密度的控制.經(jīng)分析和試算,確定在 l4邊界以內(nèi)采用較為密集的等間距橫向網(wǎng)格,而在 l4邊界以外的 S3區(qū)域采用等比數(shù)列方法控制橫向網(wǎng)格密度,假設(shè)對一段NURBS曲線,要按照等比因子k分成N份,即將參數(shù)域u∈[0,1]按k分成N份,然后根據(jù)得到的u0、u1、u2、…、uN計算NURBS曲線上的點P(ui),根據(jù)得到的這些點即可得到網(wǎng)格的邊界線段,如圖 3所示.

圖3 NURBS曲線分段示意圖Fig.3 NURBS curve subsection

設(shè)第i個分段參數(shù)變化量:

則有

由式(1)可以計算得到

以上橫向網(wǎng)格(即縱向網(wǎng)格線)密度的控制是本文貼體坐標(biāo)代數(shù)法生成自由面網(wǎng)格的特點之二.歸納以上自由面網(wǎng)格生成方法,可以給出本系統(tǒng)劃分網(wǎng)格的流程如圖4所示.

圖 4 自由面網(wǎng)格快速生成流程圖Fig.4 Flow chart of free-surface mesh generation

以上探討的是方尾三體船型自由面網(wǎng)格的劃分,而其他船型,如不含方尾的三體、雙體、單體船型,自由面網(wǎng)格劃分均屬于它的子集,更容易生成.

2 網(wǎng)格生成示例

圖 5所示為采用該系統(tǒng)生成的自由面網(wǎng)格,其中圖5(a)是某單體方尾船型自由面網(wǎng)格劃分軟件屏幕截圖,自由面計算區(qū)域參數(shù)設(shè)置與圖 1相同,縱向單位船長網(wǎng)格數(shù)為 20,橫向總網(wǎng)格數(shù)為 20,橫向網(wǎng)格密度控制因子為 1.07,方尾橫向劃分 2個網(wǎng)格,總網(wǎng)格數(shù)為1 260.

圖5(b)給出的是某5 000 t級三體船自由面網(wǎng)格劃分截圖,自由面計算區(qū)域同圖1參數(shù),縱向單位船長網(wǎng)格數(shù)為 30,橫向總網(wǎng)格數(shù)為 30,主體與側(cè)體間橫向網(wǎng)格數(shù)為 8,傾角為 45°,橫向網(wǎng)格密度控制因子為 1.1,半個主體方尾橫向網(wǎng)格數(shù)為 3,側(cè)體方尾橫向網(wǎng)格數(shù)為 2,總網(wǎng)格數(shù)為 2 935.由圖 5可以看出,所開發(fā)的自由面網(wǎng)格快速生成系統(tǒng)能夠快速生成較復(fù)雜的自由面網(wǎng)格,簡化了網(wǎng)格劃分的過程;而且通過合理控制橫向網(wǎng)格密度,節(jié)省了自由面單元數(shù),提高了自由面網(wǎng)格生成效率.

圖 5 自由面網(wǎng)格劃分示例截圖Fig.5 Free-surfacemesh examp les

3 數(shù)值實驗

采用以上自由面網(wǎng)格生成方法基于Wigley數(shù)學(xué)三體船[8]實施興波阻力數(shù)值計算,并對計算結(jié)果進行分析.Wigley數(shù)學(xué)三體船主體和側(cè)體船型函數(shù)表達式在各自坐標(biāo)系下均為

式中:L/B=10,B/T=1.6,側(cè)體長度為中體長度的一半.根據(jù)作者們多年研究發(fā)現(xiàn),當(dāng)側(cè)體位于主體中部時,不論是線性薄船理論還是Dawson型方法,都能得到與試驗結(jié)果吻合較好的興波阻力曲線,而當(dāng)側(cè)體位于主體后部時,線性薄船理論得到的結(jié)果與試驗結(jié)果差別較大,因此選擇這種側(cè)體布局作為計算實例:該三體船側(cè)體縱向偏距為0.25 L,橫向偏距為0.2 L.

根據(jù)對稱性,半(右)側(cè)自由面計算區(qū)域縱向取為[-L,2L],橫向取為[0,L],即自由面網(wǎng)格船首方向伸展0.5 L,船尾方向伸展1.5 L,船寬方向伸展L.為對比分析,采用疊模流線網(wǎng)格和本文方法生成的自由面網(wǎng)格這兩類自由面網(wǎng)格進行計算,疊模流線網(wǎng)格即原始的Dawson型方法所采用的自由面網(wǎng)格,根據(jù)疊模計算得到的結(jié)果按流線劃分自由面網(wǎng)格,物理量的導(dǎo)數(shù)采用流線上的差分來表示.

流線網(wǎng)格取為均勻網(wǎng)格,L長度內(nèi)劃分38個網(wǎng)格.由自由面疊模流線與平行于y軸的直線組成,如圖6所示,各流線在起點處y方向上間隔相等.三體船側(cè)體的存在會使得側(cè)體內(nèi)外流線 y方向上間隔不均勻,為此盡可能地保證同一寬度下各縱向位置處在y方向上的網(wǎng)格數(shù)接近相等,主體與側(cè)體之間橫向劃分網(wǎng)格數(shù)為 7.船體表面與自由面保持相同的網(wǎng)格縱向尺度,即主體縱向劃分 38個四邊形網(wǎng)格,側(cè)體縱向劃分 19個四邊形網(wǎng)格,主體和側(cè)體沿吃水方向上分別劃分10個和8個網(wǎng)格.

圖6 流線自由面網(wǎng)格Fig.6 Stream line free-surfacemesh

采用本文自由面網(wǎng)格快速生成系統(tǒng)生成兩種不同自由面上的“貼水線”網(wǎng)格進行計算,兩者區(qū)別為傾角不同,一個傾角為 0°,即垂直網(wǎng)格(橫向網(wǎng)格線垂直于x軸線),類似圖6的流線網(wǎng)格;另一個橫向網(wǎng)格線傾角為 45°,類似圖 1所示網(wǎng)格,其他參數(shù)保持相同,即縱向單位船長網(wǎng)格數(shù)為 38,橫向總網(wǎng)格數(shù)為 38,主體與側(cè)體間橫向網(wǎng)格數(shù)為 7,橫向網(wǎng)格密度控制因子為 1.0.如此參數(shù)設(shè)置是盡可能與流線網(wǎng)格保持一致,以客觀比較.

圖 7給出的是興波阻力系數(shù)計算結(jié)果與試驗結(jié)果[8]的比較.圖中流線網(wǎng)格曲線是指采用流線網(wǎng)格興波阻力程序[9]根據(jù)圖 6所示流線網(wǎng)格計算得到的結(jié)果;0°網(wǎng)格曲線、45°網(wǎng)格曲線是指采用貼體網(wǎng)格興波阻力計算程序按上述貼體網(wǎng)格計算得到的結(jié)果,縱向?qū)?shù)采用迎風(fēng)格式計算,橫向?qū)?shù)采用中心差分格式計算,導(dǎo)數(shù)計算采用物理平面與計算平面的轉(zhuǎn)換關(guān)系在計算平面上進行計算,假設(shè)計算平面坐標(biāo)(ξ,η)與物理平面坐標(biāo)(x,y)之間為一一對應(yīng)關(guān)系:

將物理空間的導(dǎo)數(shù)變換為對計算空間 ξ、η的導(dǎo)數(shù).

其中:

在計算平面上對自由面上物理量的導(dǎo)數(shù)進行計算是本文三體船興波阻力數(shù)值計算方法的第 3個特點,如前所述,前 2個特點是自由面網(wǎng)格生成技術(shù)方面的特點.此外,圖 7還給出了根據(jù)線性薄船理論興波阻力計算結(jié)果[10-11]加以比較.

圖 7 興波阻力系數(shù)計算結(jié)果對比Fig.7 Wave drag comparison between the calculated results and experimental resu lts

從圖 7可以看出,與線性薄船理論計算結(jié)果相比,線性自由面條件的三維Rankine源方法與試驗結(jié)果更加接近,薄船理論結(jié)果與模型試驗結(jié)果偏差較大,由于算例側(cè)體位于主體后部,如前文所述,側(cè)體位于主體后部,伴有較強的噴濺應(yīng)該是計算結(jié)果與模型試驗結(jié)果偏差的原因之一.

流線網(wǎng)格與 0°自由面網(wǎng)格計算結(jié)果差別甚微,這相互驗證了流線網(wǎng)格、貼體網(wǎng)格興波阻力計算結(jié)果的正確性,同時也說明采用自由面上的“貼水線”網(wǎng)格計算方法可以得到與流向網(wǎng)格方法一致的結(jié)果.流線網(wǎng)格生成費時費力,尤其是三體船興波阻力計算過程中自由面網(wǎng)格生成更加復(fù)雜,由于側(cè)體的存在,側(cè)體兩側(cè)緊貼側(cè)體表面的流線需要反復(fù)計算和調(diào)整才能得到,采用貼體網(wǎng)格預(yù)報三體船興波阻力計算過程則較為簡潔方便.

45°自由面網(wǎng)格計算結(jié)果與流線網(wǎng)格、0°自由面網(wǎng)格計算結(jié)果有一定差別,與試驗結(jié)果比較來看似乎 45°自由面網(wǎng)格計算結(jié)果與試驗結(jié)果更加接近,分析可能原因是 45°自由面網(wǎng)格更符合興波的物理現(xiàn)象,但是否說明斜網(wǎng)格一定優(yōu)于直網(wǎng)格尚待進一步研究.

4 結(jié)束語

作者成功開發(fā)了適用于三體船興波阻力計算的自由面網(wǎng)格快速生成系統(tǒng),該系統(tǒng)僅需簡單輸入必要的船體型值和幾個控制參數(shù)即可快速生成三體船興波阻力計算的自由面貼體網(wǎng)格,并實現(xiàn)網(wǎng)格密度變化控制,同時可應(yīng)用于方尾三體船,避免了商用網(wǎng)格生成軟件生成三體船興波阻力計算自由面網(wǎng)格的繁瑣過程.對一艘Wigley數(shù)學(xué)船型興波阻力的計算說明,采用系統(tǒng)生成的貼體網(wǎng)格計算三體船興波阻力,在計算過程簡潔方便的同時還可以得出良好的計算結(jié)果.

自由面網(wǎng)格傾角變化對興波阻力計算結(jié)果有一定影響,需要進一步研究.另外,方尾三體船興波阻力的計算與非方尾船略有不同,除需要該文已實現(xiàn)的生成方尾后流動區(qū)域網(wǎng)格分塊外,還需在阻力計算程序中對方尾做相應(yīng)的特殊處理,這將在后續(xù)研究中考慮.

[1]高 斌,繆國平.關(guān)于Dawson方法中網(wǎng)格劃分的研究[J].海洋工程,2000,18(1):20-28.

GAO Bin,MIAO Guoping.Investigation of grid generation in Dawson method[J].Ocean Engineering,2000,18(1): 20-28.

[2]SAUNDERS JR P.An investigation of the resistance properties of a modern trimaran combatant ship base on Taylor standard series and series64[D].CA:naval Postgraduate School,1995:1-9.

[3]盧曉平,潘雨村.高速三體船興波阻力與片體布局優(yōu)化研究[J].水動力學(xué)研究與進展(A輯),2004,19(3): 347-359.

LU Xiaoping,PAN Yucun.An investigation of wave resistance on the high speed trimarans and their piece hull position optim ization[J].Journal of Hydrodynamics(Ser A),2004,19(3):347-359.

[4]WANG Hu,ZOU Zaojian.Numerical research on wavemaking resistance of trimaran[J].J Shanghai Jiao Tong University,2008,13(3):348-351.

[5]劉應(yīng)中.船舶興波阻力理論[M].北京:國防工業(yè)出版社,2003:89-90.

[6]陳京普,朱德祥,何術(shù)龍.雙體船/三體船興波阻力數(shù)值預(yù)報方法研究[J].船舶力學(xué),2006,10(2):23-29.

CHEN Jingpu,ZHU Dexiang,HEShulong.Research on numerical prediction method for wavemaking resistanceof catamaran/trimaran[J].Journal of Ship Mechanics,2006,10 (2):23-29.

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

[8]BATTISTIN D,DANIELLIA,ZOTTI I.Numericaland experimental investigations on wave resistanceof trimaran congurations[C]//Proceedings of the 9th International Maritime Association ofMediterranean.Ischia,Italy,2000.

[9]DAWSON CW.A practical computer method for solving ship wave problems[C]//Proc of 2nd Intl Con f on Num Ship Hydrodyn.Berkeley,Calif,1977.

[10]王 中,盧曉平,鐘士崗.基于單元柯欽函數(shù)精確積分的多體船興波阻力計算[J].哈爾濱工程大學(xué)學(xué)報, 2009,30(6):602-606.

WANGZhong,LU Xiaoping,ZHONG Shigang.Multi-hull ship wave making resistance calculation based on Kochin function integral analysis expression on the surface panels [J].Journal of Harbin Engineering University,2009,30 (6):602-606.

[11]PENG Hongxuan.Numerical computation of multi-hull ship resistance and motion[D].Halifax Dalhousie University,2001:9-38.