基于三維輻射能量法的船舶波浪增阻計算分析

洪 亮，朱仁傳，繆國平，董國祥

（1.上海交通大學(xué) 船舶海洋與建筑工程學(xué)院，海洋工程國家重點實驗室，高新船舶與深海開發(fā)裝備協(xié)同創(chuàng)新中心，上海 200240；2.上海船舶運輸科學(xué)研究所，上海 200135）

0 引 言

船舶在波浪中航行時，船體所受到的阻力和船舶所需的動力相對于靜水中的情況會有所增加，這無疑會影響到船舶在海上保持航速的能力。所以船舶在波浪中的阻力增加現(xiàn)象是在船舶設(shè)計階段就必須考慮的重要因素之一。

對于船舶在波浪中阻力增加的研究始于20世紀中葉。早期的研究方法多是基于二維切片理論開展的。Maruo[1]在1957年提出了勢流求解方法，并指出波浪中船舶阻力增加主要是由縱搖和升沉運動引起的；Joosen[2]根據(jù)對波浪漂移力的分析引申出了波浪中船舶阻力增加的方法；Gerritsma和Beukelman[3]通過計算輻射阻尼波能量的大小來確定船舶在波浪中的阻力增加；Salvesen[4]給出了計入細長體的擾動力和力矩二階項時波浪中阻力平均增值的計算方法；Faltinsen[5]通過對船體濕表面上的直接壓力積分計算出阻力增加；Kihara[6]基于高速船切片理論提出了一種非線性時域數(shù)值方法。

隨著三維勢流理論在船舶耐波性中的應(yīng)用，三維頻域和時域的邊界元法越來越多地用于計算船舶水動力系數(shù)和運動預(yù)報的計算中。Papanikolaou等[7]采用了頻域三維脈動源格林函數(shù)建立邊界積分方程并求解，獲得了船舶運動和波浪載荷，并將其應(yīng)用于有航速的船舶計算中。Zakaria[8]則采用三維移動脈動源格林函數(shù)求解有航速船舶在波浪中的流場的速度勢分布，并通過對船體濕表面的直接壓力積分得到了有航速船舶在波浪中的阻力增加值。此外也有不少學(xué)者，如Kim[9]采用三維時域方法來計算船舶在波浪中的阻力增加。

總的來說，現(xiàn)有的阻力增加的計算方法可以分為兩類，即近場法和遠場法。遠場法考慮了船體的輻射波和繞射波在無窮遠處的能量和動量通量，通過對總動量的變化率的分析得到穩(wěn)定的阻力增加值。目前已有的遠場計算結(jié)果大多是基于二維的切片理論，切片理論通常是基于細長體假定且適用于高頻低航速的工況，這影響了切片理論在計算某些肥大型船舶或者是高航速工況時的計算準確性。近場法則是通過對船體濕表面上二階壓力的直接積分來計算阻力增加，而且這種方法可以通過對一階速度勢及其偏導(dǎo)數(shù)的精確計算來實現(xiàn)，這對三維船體表面速度勢的計算精度要求很高。理論上來說，兩種方法應(yīng)該可以獲得一致的結(jié)果，但實際計算中，由于兩種方法要求在近場和遠場求解流體運動速度勢，精度難以保證。船舶波阻增加是二階力范疇，相對于一階力，量級不大，精度要求高，即使采用完全相同的勢流理論及相關(guān)的數(shù)值方法來求解線性耐波性問題，上述兩種方法計算出的阻力增加結(jié)果依然有所差別。

本文基于三維頻域勢流理論采用移動脈動源格林函數(shù)法計算了流體運動速度勢及船體水動力系數(shù)和運動，根據(jù)輻射能量原理通過頻域水動力系數(shù)獲得了航行船舶一個遭遇周期內(nèi)的輻射能量，并通過能量與阻力的關(guān)系計算出船體波阻增加。分別就不同航行工況下的Wigley III、S175船以及某型油船的波阻增加進行了計算與分析，理論結(jié)果與試驗比較吻合良好。同時將本文方法與工程上常用的基于二維切片理論計算結(jié)果進行了比較與分析。研究表明了不同船型不同航速下本文方法的準確性和良好的適用性，為船舶波阻增加提供了快速有效的手段。

1 三維頻域勢流理論

三維頻域勢流理論已經(jīng)發(fā)展得比較成熟，相關(guān)的文獻資料有很多，這里只對主要公式進行簡要的介紹。

1.1 坐標系

航行船舶的運動計算所采用的坐標系o-xyz如圖1所示，坐標原點在靜止自由面上，oxy平面與靜水面重合，z軸豎直向上，船舶以定速U0航行，在波浪激勵下作六自由度運動。

圖1 船體運動坐標系Fig.1 The coordinate system for ship motion

1.2 速度勢與運動方程

假設(shè)入射波為微幅規(guī)則波，則流場內(nèi)非定常速度勢Φ=φe-iωt，定常速度勢可表達為

其中：ω 為波浪頻率，ωe為遭遇頻率，φI為入射波速度勢，φD為繞射勢，φj與Xˉj分別為第 j個運動模態(tài)的規(guī)范化速度勢與運動幅值。三維頻域勢流理論所求解的邊值問題為：

本文采取源偶混合分布法求解速度勢，邊界積分方程為：

其中：P為場點，Q為源點，α為場點所在位置的空間立體角，SB為船體濕表面，G（ P，Q ）為格林函數(shù)，WL為水線。本文采用移動脈動源格林函數(shù)來求解邊界積分方程，其表達如下：

其中：U0為航速，L1，L2分別為不同區(qū)間上的積分路徑，（x，y，z）為場點坐標，（ζ，η，ξ）為源點坐標，r=為被積函數(shù)分母的零點。Havelock型格林函數(shù)及其偏導(dǎo)數(shù)的計算可參考文獻[10]。需要說明的是，由于移動脈動源格林函數(shù)在自由面上不收斂，所以在求解邊界積分方程的時候，實際上忽略了水線積分項的影響。

對邊界積分方程離散求解后獲得各個運動模態(tài)的附加質(zhì)量、阻尼系數(shù)以及波浪力，并以此建立船舶6自由度的運動方程：

其中：Mkj為廣義質(zhì)量矩陣，Akj（ω ）為附加質(zhì)量，Bkj（ω ）為阻尼系數(shù)，Ckj為回復(fù)力矩陣，F(xiàn)wk（ω）為k方向的波浪力。

1.3 船舶在波浪中阻力增加

從能量的角度出發(fā)，在船速為U0，波速為C的情況下，于船后某一位置處一個遭遇周期內(nèi)波浪的輻射能量與船舶阻力平均增值所做的功相等，即：

其中：E為在一個波浪遭遇周期內(nèi)整個船體所受阻尼力輻射的能量，Te為船舶波浪的遭遇周期，RAW為船舶在波浪中的阻力增加值，λ為入射波波長。

一個波浪遭遇周期內(nèi)船舶輻射能量為各模態(tài)輻射力所做的功之和，則能量E的可表達如下：

其中：Fj為 j模態(tài)的輻射力別為船體j模態(tài)的位移、速度和加速度，將其代入（7）式則可得到船舶在波浪中阻力增加的最終表達式：

其中：XkR，XkI分別為第k模態(tài)的運動幅值的實部與虛部。（8）式可以看出，船舶在波浪中輻射運動引起的阻力增加是由兩部分組成的，一是與附加質(zhì)量相關(guān)的增阻RAM，一是與阻尼系數(shù)相關(guān)的增阻RDP。在無航速時，由于水動力系數(shù)存在Ajk=Akj的關(guān)系，可以證明RAM=0，這樣（8）式則可表達如下：

（9）式與Joosen法[2]一致。而Joosen在推廣到有航速的情況時僅僅是用遭遇頻率ωe代替ω0而忽略了附加質(zhì)量對波浪中阻力增加的影響。

需要說明的是，本文采用了Salvensen法公式中與繞射相關(guān)項對阻力增加進行了修正，詳細做法可參考文獻[4]。后面給出了對輻射增阻和繞射增阻的比較。

2 數(shù)值計算結(jié)果的比較與分析

為檢驗上述方法的正確性，本文選取了具有試驗結(jié)果的三類不同類型船舶的阻力增加進行計算，三類船分別為數(shù)學(xué)船型Wigley III、細長型集裝箱船S175和某肥大型油船。

2.1 Wigley III

Wigley III型船模的型線如圖2所示，其主尺度見表1。對Wigley III船模在Fn=0.2、0.3和0.4三個航速下進行了頻域數(shù)值計算，浪向角為180°迎浪。

圖3-5分別給出了傅汝德數(shù)Fn=0.2、0.3和0.4時Wigley III的垂蕩、縱搖和阻力增加的計算結(jié)果，圖中DUT試驗數(shù)據(jù)以及切片理論的數(shù)值計算結(jié)果均來自文獻[16]。與試驗結(jié)果相比，采用三維移動脈動源格林函數(shù)所計算出的Wigley III的垂蕩運動結(jié)果略好于切片理論，而縱搖運動的結(jié)果則比切片理論的有不小的改善。從圖4中可以看出，隨著航速的提高，切片理論所計算的縱搖運動RAO的峰值頻率與試驗結(jié)果向高頻段的偏移程度在增加。圖6給出了Fn=0.3時Wigley III水動力系數(shù)的計算結(jié)果，可以看到，兩種方法所計算出的垂蕩模態(tài)的水動力系數(shù)均與試驗吻合良好，而三維移動脈動源格林函數(shù)所計算出的縱搖模態(tài)水動力系數(shù)比切片理論更好。這也解釋了二者計算出的縱搖RAO有所差別的原因。通過（8）式可知，船舶在波浪中的阻力增加是與其運動幅值及水動力系數(shù)有關(guān)的。圖5的結(jié)果可以看出三維移動脈動源格林函數(shù)計算出的波浪中阻力增加的結(jié)果在3個航速下均與試驗?zāi)茌^好地吻合，僅在Fn=0.2的峰值上偏??；而切片理論由于其計算出的水動力系數(shù)和RAO的偏差使其阻力增加的計算結(jié)果也產(chǎn)生了偏差，結(jié)果中的共振頻率較之試驗向高頻偏移。

圖2 Wigley III船型線圖Fig.2 The body plan of Wigley III

圖3 Wigley III垂蕩運動RAOFig.3 RAO of heave motion of Wigley III

圖4 Wigley III縱搖運動RAOFig.4 RAO of pitch motion of Wigley III

圖5 Wigley III阻力增加Fig.5 Added resistance of Wigley III

表1 Wigley III主尺度Tab.1 Principle dimension of Wigley III

續(xù)表1

圖6 Wigley III水動力系數(shù)Fig.6 Hydrodynamic coefficients of Wigley III

圖7 S175船型線圖Fig.7 The body plan of S175

2.2 S175集裝箱船

S175集裝箱船是船舶耐波性試驗的標準模型之一。其實船主尺度見表2。圖7所示為S175船型線圖。計算航速分別為Fn=0.15、0.2、0.25 和 0.3，浪向角為 180°迎浪。

圖8給出S175船在迎浪工況下以4種不同航速航行時的阻力增加數(shù)值結(jié)果。采用了Fuji[17]、Takahashi[18]和ITTC報告中Nakamura[19]的試驗結(jié)果進行了比較，同時給出了基于切片理論的Joosen法(以下簡稱J法)和Gerritsma和Beukelman法（以下簡稱G&B法）的計算結(jié)果?？梢钥闯?，在所給的4種航速下，除Fn=0.2時J法的計算效果較差外，所給的三種數(shù)值方法可較為準確地計算出S175船舶在波浪中航行時阻力增加的峰值位置和大小。在低頻長波區(qū)域即圖8中大于峰值位置的區(qū)域，G&B法的計算結(jié)果比試驗結(jié)果要偏大，本文方法和J法的計算結(jié)果均和試驗吻合較好且本文方法略優(yōu)于J法。在高頻短波區(qū)域即圖8中小于峰值位置的區(qū)域，J法和G&B法的結(jié)果均小于試驗結(jié)果，G&B法在高頻區(qū)域出現(xiàn)了振蕩現(xiàn)象，J法的高頻結(jié)果基本趨近于0，而本文方法在高頻區(qū)域的計算結(jié)果依然可以和試驗吻合較好。綜合來看，本文方法在全部的頻率區(qū)間上的計算結(jié)果均優(yōu)于基于二維切片理論的G&B法和J法。圖9給出了S175船在Fn=0.25時的波浪增阻的成分。圖中可以看出輻射增阻是船舶在波浪中阻力增加的最主要成分，繞射增阻隨入射波波長的增加而緩慢減小。在中低頻率的入射波中，繞射增阻在總增阻中所占的比重很小，此時輻射增阻是主要的；在高頻率的入射波中，繞射增阻在總增阻中的比重變大，而輻射增阻的比重變小，二者均為總增阻的主要成分。

表2 S175主尺度Tab.2 Principle dimension of S175

圖8 S175船阻力增加Fig.8 Added resistance of S175 containership

2.3 油船

油船的試驗是在上海船舶運輸科學(xué)研究所進行的，表3給出了某型油船的主尺度。由于不便給出此船的型線圖，圖10給出了油船數(shù)值計算網(wǎng)格的示意圖，由圖中可以看出此船屬于肥大型船舶，其方形系數(shù)Cb=0.83。

圖9 S175船阻力增加成分Fig.9 Components of added resistance

圖10 油船計算網(wǎng)格示意圖Fig.10 Calculation mesh of tanker

圖11 油船阻力增加Fig.11 Added resistance of tanker

表3 油船主尺度Tab.3 Principle dimension of tanker

圖11給出了四種航速下油船的波浪中阻力增加的結(jié)果，浪向角為180°迎浪。圖中可以看出本文方法能夠計算得到較為準確的波浪增阻。四種不同航行工況下，理論計算的波阻增加的幅值、趨勢以及峰值與共振頻率都與試驗結(jié)果吻合良好。隨著航速的增加，波浪增阻的共振頻率基本沒有變化，但是峰值略有增加。結(jié)果說明本文方法在計算油船這種肥大型船舶時也有著令人滿意的計算準確度。

2.4 結(jié)果分析

從（8）式給出的阻力增加計算表達來看，本文方法基于船舶三維水動力系數(shù)和運動可以簡單快捷地計算出船舶在波浪中的阻力增加值，避免了由于物面速度勢梯度的計算誤差帶來的影響。上述給出的3種不同船型的數(shù)值結(jié)果均采用的是均值面元法對速度勢進行求解，計算網(wǎng)格少且效率高，計算結(jié)果令人滿意。

船舶在波浪中阻力增加的計算準確度與其水動力系數(shù)和運動直接相關(guān)。本文方法的計算結(jié)果比G&B法和Joosen法更優(yōu)。一方面是由于本文方法是基于三維有航速船舶水動力計算的結(jié)果，其水動力系數(shù)和運動結(jié)果均優(yōu)于二維切片理論；另一方面，本文方法考慮了有航速時船舶附加質(zhì)量對阻力增加的影響，這是G&B法和Joosen法的表達中所忽略的。

本文的計算對象是三種不同船型，計算涵蓋了Fn從0.15至0.4的工況，從計算結(jié)果來看，可以較為準確地計算不同航速下船舶在波浪中的阻力增加值及共振頻率和峰值。研究表明，無論是細長體船型S175還是肥大型油船 ，都可以獲得和試驗值相吻合的結(jié)果，說明本文方法有著廣泛的適用性。

3 結(jié) 論

本文基于三維頻域勢流理論采用移動脈動源格林函數(shù)法對有航速船舶水動力進行計算，并采用輻射能量原理計算船舶在波浪中的阻力增加。對不同航速下Wigley III，S175船以及某型油船的運動及阻力增加進行了頻域數(shù)值計算與分析，有如下結(jié)論：

（1）本文方法的數(shù)值計算結(jié)果與試驗吻合良好，可以較好地對航行船舶在頂浪中的阻力增加進行數(shù)值計算。

（2）相較于工程上常用基于切片理論的方法，本文方法直接通過整船的三維水動力系數(shù)及運動結(jié)果計算獲得阻力增加值，以三維完整的船體為計算對象，在保證較好的計算效果的同時有著很高的計算效率。

（3）本文方法在計算具有細長體船型的S175船和肥大型油船時都有著令人滿意的計算效果，說明本文方法有著廣泛的適用性，為波浪中船舶阻力增加的預(yù)報提供了快速有效的計算手段。

參 考 文 獻：

[1]Maruo H.The excess resistance of a ship in a rough seas[J].International Shipbuilding Progress,1957,4(85)：337-345.

[2]Jossen W P A.Added resistance in waves[C]//6th Symposium on Naval Hydro-dynamics.Washington,USA,1966.

[3]Gerritsma J,Beukelman W.Analysis of the resistance increase in waves of a fast cargo-ship[J].International Shipbuilding Progress,1972,18(217)：285-293.

[4]Salvensen N.Added resistance of ships in waves[J].Journal of Hydronautics,1978,12(1)：24-34.

[5]Faltinsen O M,Minsaas K J,Liapis N,Skj?rdal S.Prediction of resistance and propulsion of a ship in a seaway[C]//Proceedings of the 13th Symposium on Naval Hydrodynamics.Tokyo,Japan,1980.

[6]Kihara H,Naito S,Sueyoshi M.On the practical prediction of added resistance using nonlinear slender body theory[C]//Proceeding of 4th Osaka Colloquium on Sea-keeping Performance of Ships.Osaka,Japan,2000.

[7]Papanikolaou A,Zaraphonitis G.On an improved method for the evaluation of second-order motions and loads on 3D floating bodies in waves[J].Journal of Schiff-stechnik,1987,34：170-211.

[8]Zakaria N M G.Prediction of added resistance of ships by 3-D Green function method and its comparison with semi-3-D and 2-D approaches[J].IE(I)Journal M R,2008,88：16-25.

[9]Kim K H,Kim Y.Numerical study on added resistance of ships by using a time-domain rankine panel method[J].Ocean Engineering,2011,38：1357-1367.

[10]洪 亮,朱仁傳,繆國平,范 菊.三維頻域有航速格林函數(shù)的數(shù)值計算與分析[J].水動力學(xué)研究與進展A輯,2013,4：423-430.Hong Liang,Zhu Renchuan,Miao Guoping,Fan Ju.Numerical calculation and analysis of 3-D Green’s function with forward speed in frequency domain[J].Journal of Hydrodynamics,2013,4：423-430.

[11]Tsujimoto M,Kuroda M.On a calculation of decrease of ship speed in actual seas[J].Journal of the Japan Society of Naval Architects and Ocean Engineers,2009,9：79-85.

[12]Tsujimoto M,Takagi K.A calculation method of 10 model index for ships[J].Journal of the Japan Society of Naval Architects and Ocean Engineers,2009,10：97-104.

[13]Kuroda M,Tsujimoto M,Takagi K.Influence of bow flare shapes on added resistance in waves[J].Journal of the Japan Society of Naval Architects and Ocean Engineers,2010,11：101-107.

[14]段文洋,李傳慶,劉 亮.關(guān)于船舶繞射增阻計算方法的評注[J].中國造船,2011,S1：38-46.Duan Wenyang,Li Chuanqing,Liu Liang.A comment on calculation methods of ship added resistance due to wave diffraction[J].Shipbuilding of China,2011,S1：38-46.

[15]Alexandersson M.A Study of methods to predict added resistance in waves[D].KTH Centre for Naval Architecture,2009.

[16]Journee J M J.Experiments and calculations on 4 Wigley hull forms in head waves[DB/OL].(1992-05-01)[2003-11-26].http：//www.shipmotions.nl/DUT/PapersReports/0909-DUT-92.pdf.

[17]Fuji H,Takahashi T.Experimental study on the resistance increase of a ship in regular oblique waves[C].Proc.14th ITTC,1975,4：351.

[18]Takahashi T.A practical prediction method for added resistance of a ship in waves and the direction of its application to hull form design[J].Transactions on the West-Japan Society of Naval Architects,1988,75：75-95.

[19]TTC Seakeeping Committee.Prediction of added resistance[C].In：Proceedings of the 17th ITTC,1984：500-503.