三體船興波問題的數值計算

周利蘭 ， 高 高 ， 尹 巍

（武漢理工大學 a.交通學院；b.高速船舶工程教育部重點實驗室，武漢430063）

1 引 言

近年來三體船的研究漸趨活躍。三體船的水下部分由主體和兩個小側體組成，通過連接橋將側體與主體連接成一體，這種船型構造使三體船的興波阻力小，兩個側體又能提供足夠的穩性,連接橋具有提高總縱強度的功能，利于形成寬闊甲板面，該船型還具有優良的耐波性。三體船用作軍船還有防護能力、破損穩性以及隱蔽性好等優點。盡管高速三體船與常規單體船相比也有制造工藝復雜、制造成本高的缺點，但性能方面的諸多優勢使其應用前景廣闊。

三體船除單個片體自身首尾波系之間的興波干擾外，其主體和兩個側體之間的波系干擾導致復雜的興波波形，因而其阻力特性比單體船的要復雜得多。Zafer Elcin[1]采三維Rankine源面元法系統地研究了側體位置對三體船興波阻力的影響。Yang[2]等分別采用包括新細長體理論在內的三種線性理論方法進行三體船側體位置的優化分析。陳京普等[3]采用改進的Dowson方法，開發了可用于雙體船、三體船興波阻力的數值預報程序。近年發展的船舶興波問題的基于NURBS的廣義高階面元法算法，對排水型單體船和雙體船的計算都獲得了成功[4]，本文討論將其應用于三體船興波問題計算。

2 基本方程與邊界條件

在勢流假定下，擾動速度勢φ滿足Laplace方程：

物面邊界條件：

自由表面運動學邊界條件：

自由表面動力學邊界條件：

無窮遠條件：

再加上輻射條件即構成定解問題。實際計算中自由表面邊界條件采用（3）式與（4）式的耦合條件：

3 數值方法

以Rankine源為Green函數的勢流解法由于其函數形式簡單且因此易于推廣應用于非線性計算問題，所以從它問世到現在的二十多年間一直都是船舶興波勢流計算的主流方法。船舶流體力學問題數值計算中基于NURBS高階面元法的引入，起到了計算流體力學CFD與船舶工程CAD之間的橋梁作用。

與通常邊界元法不同的是，在以上述速度勢表達式通過邊界條件構造的求解方程組中，直接的變量不是源強密度而是源強分布曲面的控制頂點，在邊界面上選取的配置點數目應與源強控制頂點數目相同。

輻射條件的滿足是以Rankine源為奇點的面元方法的主要困難之一。本文采用自由面上置源，從而可方便地以配置點移動法數值上滿足輻射條件。高高曾分析過有關數值誤差，相應參數的取值可參見有關文獻[5]。縱向采用均勻網格，配置點前移可由自由面上的源分布曲面后移來實現。最后，自由面上的源分布曲面表為：

式中的δx和zs分別代表自由面上的源分布面沿縱向及垂向的偏移量。此處δx取縱向網格間距，zs取縱向網格間距的1.5倍。

約定以上標F表示相應物理量為自由面上的，上標B表示物面上的，記源點的參數坐標為（u0,v0），則源分布對計算域內任一點（u,v）產生的誘導速度為：

將以上速度勢的表達式代入物面邊界條件、自由表面邊界條件后，即構成了以曲面控制頂點及為未知數的代數方程組，未知數的總數為由于NURBS曲面積分難以解析進行，文中采用分片高斯積分。解得及后，由（9）式決定速度勢，進而求得波高及興波阻力。

4 數值求解與結果

本文數值計算采用兩種船型，一種是文獻[6]所使用的數學船型；另一種采用文獻[7]給出的方尾船型及布局，側體的尺度比主體的小很多，主體、側體主尺度要素見下表1。計算域取主體船前0.5水線長，主體船后1.5水線長，主體船側1.0水線長，主體和側體均取21×6個配置點，自由面取61×20個配置點。對于三體船，側體流場的不對稱性會對興波計算帶來不同程度的影響，對于船型一，注意到比較的參考文獻[6]未計及此影響，本文亦在未計及該影響的條件下進行計算比較；而對于方尾船型二，由于本文在方尾處滿足以下方尾邊界條件[7]：

該條件等價于以下流動在方尾下緣切向光滑脫體的條件[8]（從某種意義上類似于無界域問題中的Kutta條件）：

式中：ζ為波高；h為方尾邊沿高度；xT,y為方尾處縱、橫向坐標；x+T的上標“＋”，表示從x的正向（自由面上）趨近于xT（在自由面上）。

為表達三體船側體布置，記：主體船中至側體船中的縱向間距為Xs（側體在主體中前為正），主體船中至側體船中的橫向間距為Ys，如圖1所示。

圖1 三體船側體與主體位置關系Fig.1 Configuration of trimaran

文獻[6]主體和側體均采用如下的數學船型：

文中主體船長、船寬及吃水均三倍于側體。計算取Xs為-1/3L，Ys為0.2L，圖2給出了興波阻力系數隨Fn數的變化曲線,圖3給出了Xs=-1/3L，Ys=0.2L，Fn=0.5下的三維波形圖。計算所得興波阻力系數與試驗結果符合較好，優于文獻[6]的計算結果。

表1 三體船主體和側體主要要素Tab.1 Main dimensions of trimaran

圖 2 計算 Cw與文獻 Cw之比較（Xs=-1/3L，Ys=0.2L）Fig.2 Calculated wave-making resistance coefficient Cwand referenced Cw(Xs=-1/3L,Ys=0.2L)

圖 3 三維波形圖（Xs=-1/3L，Ys=0.2L，Fn=0.5）Fig.3 Three dimensional wave pattern of trimaran(Xs=-1/3L,Ys=0.2L,Fn=0.5)

對于船型二，選取不同的側體布置參數Xs、Ys，由Xs、Ys組合成不同的側體布置方案，本文數值計算中分別取 Xs為-0.3L，-0.2L、0L 和-0.1L，Ys為 0.086 5L、0.092 3L、0.105 8L 和 0.111 7L。 圖 4－9分別反映了Xs為-0.3L，-0.2L和-0.1L，Ys為0.092 3L和0.111 7L這6種情況下的興波阻力系數隨Fn數的變化，并與模型試驗結果得到的剩余阻力系數進行了比較（此剩余阻力系數含有主體、側體的粘性干擾影響）[9]，圖10、圖11給出了不同橫向偏距和縱向偏距共8種情況下興波阻力預報曲線，圖12給出了Xs=-0.3L、Ys=0.073 1的三體船阻力系數圖，圖13、圖14分別Xs=-0.2L，Ys=9 m（0.086 5L），Fn=0.6的船舶興波的三維波形圖和波形等高線圖。注意到模型試驗結果得到的剩余阻力系數除含有粘壓阻力、曲率對摩擦阻力的影響，主體、側體的興波干擾影響外，還包含主體、側體的粘性干擾影響及本文興波計算采用的是線性自由面邊界條件等因素，與模型試驗得到的剩余阻力系數相比，計算興波阻力系數隨Fn數的變化合理。計算結果表明：縱向偏距Xs對三體船Cw有顯著影響，高速時橫向偏距也對Cw有較大影響，這與其他學者對高速三體船側體的橫向偏距和縱向偏距對高速三體船興波阻力影響規律的研究結果基本一致[10-11]；圖12中與文獻[12]的興波阻力系數計算結果及試驗所得剩余阻力系數的比較，結果表明本文計算結果較優于文獻[12]的計算結果。

圖4 計算的興波阻力系數Cwp與實驗Cr之比較（Xs=-0.1L，Ys=0.092 3L）Fig.4 Calculated wave-making resistance coefficient Cwp and experimental residual resistance coefficient Cr(Xs=-0.1L,Ys=0.092 3L)

圖5 計算的興波阻力系數Cwp與實驗Cr之比較（Xs=-0.2L，Ys=0.092 3L）Fig.5 Calculated wave-making resistance coefficient Cwp and experimental residual resistance coefficient Cr(Xs=-0.2L,Ys=0.092 3L)

圖6 計算的興波阻力系數Cwp與實驗Cr之比較（Xs=-0.3L，Ys=0.092 3L）Fig.6 Calculated wave-making resistance coefficient Cwp and experimental residual resistance coefficient Cr(Xs=-0.3L,Ys=0.092 3L)

圖8 計算的興波阻力系數Cwp與實驗Cr之比較（Xs=-0.2L，Ys=0.111 7L）Fig.8 Calculated wave-making resistance coefficient Cwp and experimental residual resistance coefficient Cr(Xs=-0.2L,Ys=0.111 7L)

圖7 計算的興波阻力系數Cwp與實驗Cr之比較（Xs=-0.1L，Ys=0.111 7L）Fig.7 Calculated wave-making resistance coefficient Cwp and experimental residual resistance coefficient Cr(Xs=-0.1L,Ys=0.111 7L)

圖9 計算的興波阻力系數Cwp與實驗Cr之比較（Xs=-0.3L，Ys=0.111 7L）Fig.9 Calculated wave-making resistance coefficient Cwp and experimental residual resistance coefficient Cr(Xs=-0.3L,Ys=0.111 7L)

圖10 Ys=9 m興波阻力預報曲線Fig.10 Calculated wave-making resistance coefficient curves for Ys=9 m

圖11 Ys=11 m興波阻力預報曲線Fig.11 Calculated wave-making resistance coefficient curves for Ys=11 m

圖12 三體船阻力系數比較（Xs=-0.3L，Ys=0.073 1L）Fig.12 Resistance coefficients of trimaran(Xs=-0.3L,Ys=0.073 1L)

圖13 三維波形圖（Xs=-0.2L，Ys=9 m，Fn=0.6）Fig.13 Three dimensional wave pattern of triamaran(Xs=-0.2L,Ys=9 m,Fn=0.6)

5 結 語

本文所采用基于NURBS的廣義高階面元法對于三體船興波波形與興波阻力的數值計算可以反映三體船的興波特性，得到較為合理且優于薄船理論的計算結果。對于高速方尾三體船側體的橫向偏距和縱向偏距對興波阻力影響的初步計算結果與其他學者的研究結果基本一致。今后將對有關的非線性問題作進一步的研究探討。

圖14 三維波形等高線圖（Xs=-0.2L，Ys=9 m，Fn=0.6）Fig.14 Contour map of three dimensional wave pattern for triamaran(Xs=-0.2L,Ys=9 m,Fn=0.6)

[1]Elcin Z.Wave Making Resistance Characteristics of Trimaran Hulls[D].Monterey,CA:Master’s Thesis,Naval Postgraduate School,2003.

[2]Yang C,Noblesse F,Lohner R.Practical hydrodynamic optimization of a trimaran[C]//SNAME Annual Meeting.Orlando,FL,2001:185-196.

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

Chen Jingpu,Zhu Dexiang,He Shulong.Research on numerical prediction method for wavemaking resistance of catamaran/trimaran[J].Journal of Ship Mechanics,2006,10(2):23-29.

[4]高 高,馬 嶺,鄒早建.排水型高速雙體船興波問題的數值計算[J].武漢理工大學學報,2003,25(6):586-596.

[5]高 高.船波問題的輻射條件、自由面條件的研究及面元法誤差分析[D].武漢:武漢交通科技大學,1997.

[6]Suzuki K,Kai H.Numerical estimation of wave making characteristics of high-speed ships taking sinkage and trim effects into account[C]//International Conference on Fast Sea Transportation,FAST’2005,June 2005.St.Petersburg,Russia,2005.

[7]Gao Gao.Numerical implementation of transom conditions for high-speed displacement ships[J].Journal of Ship Mechanics,2006,10(3):1－9.

[8]Nakos D E,Sclavounos P D.Kelvin wake and wave resistance of cruiser-and transom-stern ships[J].Journal of Ship Research,1994,38(1):9-29.

[9]李培勇,裘泳銘,顧敏童,王文富.三體船阻力模型試驗[J].中國造船,2002,43(4):6-12.

[10]Tuck E O,Lazauskas L.Optimum hull spacing of a family of multihulls[J].Applied Mathematics Publications,1998,45:180-195.

[11]酈 云,盧曉平.高速三體船阻力性能研究[J].船舶力學,2007,11(2):191-198.

Li Yun,Lu Xiaoping.An investigation on the resistance of high speed trimarans[J].Journal of Ship Mechanics,2007,11(2):191-198.