非線性水動力導數的數值計算與研究*

趙小仨 徐海祥

(高性能船舶技術教育部重點實驗室1) 武漢 430063) (武漢理工大學交通學院2) 武漢 430063)

非線性水動力導數的數值計算與研究*

趙小仨1,2)徐海祥1,2)

(高性能船舶技術教育部重點實驗室1)武漢 430063) (武漢理工大學交通學院2)武漢 430063)

針對船舶的非線性運動難以界定和非線性運動難以預報的問題，以供應船為研究對象，采用CFD商用軟件FLUENT，結合動網格技術對大振幅平面運動機構試驗進行數值模擬，通過對比不同工況的流場壓力云圖，分析得出供應船水動力達到非線性的振幅范圍.設計供應船非線性運動的試驗方案，分別模擬不同頻率時的大漂角斜航運動及大振幅的純縱蕩、純橫蕩、純首搖、組合運動，擬合得到接近零頻率的非線性水動力導數.

非線性水動力；大振幅PMM試驗；數值計算；供應船

0 引 言

船舶操縱性與船舶航行安全緊密相關，是重要的水動力性能之一.近些年，國際海事組織(international maritime organization, IMO)前后頒布了A.751(18)和MSC.137(76)號決議，針對船舶操縱性的問題提出了明確的要求，并建議各國政府機構按要求執行.SIMMAN 2008和SIMMAN 2014的研討會，評估了CFD(computational fluid dynamics, CFD)方法預報船舶操縱性的能力.第25屆ITTC(international towing tank conference, ITTC)操縱會議對現有的船舶操縱性預報方法做了總結.總之，船舶操縱性能越來越受到造船界的重視[1-3].

水動力導數對船舶操縱性的預報至關重要.目前，通過平面運動機構試驗(planar motion mechanism test, PMM)確定船舶水動力導數是最可靠的方法之一.從SIMMAN 2008發布了針對三個標準船模進行的PMM試驗的實驗數據以來，國內外學者開始對CFD模擬PMM試驗進行探究.Turnock等[4-12]用CFD軟件模擬小振幅PMM試驗，求取線性水動力導數；Toxopeus等[13-16]建立非線性水動力模型，模擬大振幅PMM試驗，求取非線性水動力導數.

雖然許多學者對數值模擬PMM試驗做了大量研究工作，但是迄今沒有一個定量的標準來判斷船模的運動是否達到非線性范疇，相關文章也較少.評判船舶的運動是否達到非線性，不僅取決于運動參數，還與船型等因素有關.文中將以供應船為研究對象[17]，通過數值模擬船模不同運動幅值的PMM試驗，分析出供應船水動力達到非線性的運動幅值范圍.在此基礎上，設計試驗工況，計算零頻率附近的非線性水動力導數.

1 數學模型

研究船舶在大振幅下的操縱運動，用線性水動力模型很難準確的表達船舶所受到的水動力，為了更準確的描述船舶的運動，須考慮運動狀態的非線性項[18].根據經驗，在粘性類流體動力和力矩的泰勒級數展開式中保留至三階項，對描述船舶在常速域中的運動已足夠精確.

1) 流體慣性力(矩)XI，YI，MI

(1)

式中：u，v，r分別為船舶的縱向、橫向速度與轉首角速度.

2) 粘性類流體動力(矩)Xh，Yh，Mh

(2)

3) 水動力導數的無因次化 采用MMG模型系統建議的以Ld為參考面積進行無因次化.

2 平面運動機構試驗

PMM試驗是約束模試驗的一種.是通過測量船模所受到的水動力和力矩，求得計算船舶操縱運動所需的各種加速度導數、速度導數和耦合導數.PMM試驗有小振幅和大振幅之分，前者只能測定線性水動力導數，后者可以測定非線性水動力導數.本文主要模擬PMM試驗的以下幾種運動形式：定漂角斜航；純橫蕩運動；純縱蕩運動；純首搖運動；組合運動.

3 數值計算方法

3.1 計算模型

計算模型為一艘75 m供應船，縮尺比為1∶20.船模幾何參數見表1.三維模型見圖1.

表1 供應船模參數

圖1 供應船模型

3.2 坐標系

坐標系見圖2，O點位于船舯；X軸指向船首為正；Y軸指向左舷為正；Z軸正向依據右手定則確定.

圖2 隨船坐標系

3.3 劃分計算域

計算域尺寸船首上游取1.5倍船長，船尾下游取3倍船長，船兩側取2倍船長，水深方向取8.3倍吃水，見圖3.

圖3 供應船的計算域

3.4 劃分網格以及驗證收斂性

為了保證網格質量，采用分塊全結構化網格，并在首尾部以及呆木處進行適當加密.

以表1中的船模為研究對象，采用3種不同數量的網格模擬相同工況下的純首搖運動，進行網格收斂性的驗證.網格數分別為100萬，200萬和300萬.圖4給出了不同網格數計算得到的縱向力X、橫向力Y和轉首力矩M在一個周期內的曲線.由圖4可知，網格數量從100萬增至200萬，計算結果變化明顯；當網格數量從200萬增至300萬時，計算結果幾乎不變.因此選擇250萬左右的網格，既能保證收斂性，又能節省計算資源和計算時間.

圖4 不同網格數的計算結果對比

3.5 數值計算方法

1) 邊界條件 ①入流與出流邊界條件，入流面設為速度入口；出流面設為自由流出口，權重為1；②船體表面，在船體表面施加無滑移壁面條件；③自由面，考慮到供應船舶在定位過程中時低速航行，忽略自由面興波的影響，將自由面設為對稱面.

2) 定義動網格(dynamic mesh) 編寫純橫蕩與純首搖運動的用戶自定義函數(UDF)，并使之與FLUENT相關聯；網格更新方法(mesh methods)選擇網格光順方法(smoothing)和動態層方法(layering).

3) 離散格式和求解算法 非穩態流動；壓強插值格式選用標準格式；空間離散采用二階迎風差分格式；時間積分方案采用一階隱式；求解算法采用基于速度-壓力耦合的SIMPLEC算法.

4 試驗設計及計算結果分析

4.1 供應船非線性模型的振幅范圍

計算工況見表2～3.其中，v′=va/Uc和r′=ra·L/Uc分別為純橫蕩和純首搖無因次化的速度幅值.

表2 純橫蕩運動計算工況

表3 純首搖運動計算工況

圖5～6給出了不同運動幅值時的純橫蕩與純首搖運動所受橫向力Y在一個周期的變化曲線.圖7～8給出了純橫蕩與純首搖運動橫向力Y高階量的實際值與三階擬合值的對比.從圖中可以得出以下結論.

1) 隨著橫蕩幅值和首搖角幅值的增大，橫向力幅值明顯增大，但橫向力的相位幾乎沒有發生變化.

2) 當振蕩幅度超過某一幅值，橫向力已經不再是一階正(余)弦的形式，即運動幅值越大，非線性表現的越強.

3) 由圖6可知，當初始首向角為8°時，純首搖運動表現為一階正余弦的形式；當初始首向角為18°時，純首搖運動已經不再是一階正(余)弦的形式.

4) 由圖7～8可知，當橫蕩運動的橫向幅值y超過0.9 m(運動過程中的最大漂角約為31°)，當首搖運動的初始首向角超過30°時，橫向水動力高階量的實際值與三階擬合值間的誤差越來越大.

圖5 75 m供應船純橫蕩運動橫向力Y隨橫蕩幅值的變化

圖6 75 m供應船純首搖運動橫向力Y隨初始首向角的變化

圖7 純橫蕩運動Y方向高階量的擬合值與實際值的對比

圖8 純首搖運動Y方向高階量的擬合值與實際值的對比

因此，對于75 m供應船模，其非線性水動力模型的振幅范圍約為漂角β大于10°小于30°.

4.2 流場分析

圖9～10分別給出了75 m供應船首搖角幅值為8°和25°的純首搖運動一個周期內的船體表面壓力云圖.壓力變化的總體趨勢是：縱向來看，從船首向船中壓力逐漸遞減，從船中向船尾壓力逐漸遞增；橫向來看，隨著首向角轉動位置與方向的不同，最大壓力的分布左右交替出現.通過對比圖9～10可知，首搖幅值為25°時的壓差比首搖幅值為8°的壓差要大得多.

圖9 75 m供應船首搖幅值為8°的船體表面壓力分布圖

圖10 75 m供應船首搖幅值為25°的船體表面壓力分布圖

由圖9可知，前半周期船體表面的壓力分布情況與后半周期是關于船體中縱剖面對稱的.圖9b)表示首搖運動3/8周期時的壓力分布圖，它的左舷壓力分布恰好與圖9c)中右舷壓力分布相同.出現該種現象的主要原因是：純首搖運動的前半周期與后半周期的轉首角速度大小相同方向相反，而且當首搖幅值為8°時，屬于小振幅振動，符合線性假定，因此由運動輸出的力在前半周期與后半周期是左右舷對稱的.

由圖10可知，前半周期船體表面的壓力分布情況與后半周期已不再是關于船體中縱剖面對稱.其主要原因是：當首搖幅值為25°時，屬于大振幅振動，運動的輸入與輸出已不是線性關系.

4.3 非線性水動力導數的計算結果

以75米供應船為研究對象，根據4.1節所得結論，確定非線性水動力的計算工況.進而通過模擬直航、大漂角斜航,以及頻率為0.2～0.4 rad/s的大振幅純縱蕩、純橫蕩、純首搖和組合運動，計算各運動模態下的水動力和水動力矩，然后求取不同頻率時非線性模型的水動力導數，最終得到75 m供應船在零頻率附近的各水動力導數[19].

1) 數據處理方案 以純橫蕩運動為例，介紹本文的數據處理方案.其他運動形式的數據處理方法類似.根據純橫蕩的運動規律，非線性水動力數學模型可變形為：

(3)

分別計算不同振動頻率時船體所受的水動力X，Y，N，再用最小二乘法把水動力按上式擬合，分離出系數Ys1，Yc1，Yc3，Ms1，Mc1，Mc3，進而得到相應的水動力導數.

2) 直航運動 通過模擬0.4m/s～0.8m/s的直航運動，可求取縱向速度導數Xu，Xuu，Xuuu.圖11給出了不同航速時縱向力的變化曲線.用最小二乘法對圖11中的曲線進行擬合，再根據縱向水動力的三階非線性表達式，可以直接得出3個縱向水動力導數，結果見表4.

圖11 直航運動縱向力的擬合曲線

無因次量X'uX'uuX'uuu數值×10-3-0.42-2.370.21

3) 斜航運動 通過模擬10°～30°不同漂角下的斜航運動，可求取速度導數Xvv，Yv，Yvvv，Nv，Nvvv.圖12給出了不同漂角時縱向力、橫向力和力矩的變化曲線.與直航運動求取水動力導數的方法相同，得出斜航運動對應的五個水動力導數，結果見表5.

圖12 斜航運動中縱向力、橫向力及首搖力矩的擬合曲線

無因次量X'vvY'vM'vY'vvvM'vvv數值×10-20.50-2.95-1.63-7.820.27

表6 純縱蕩計算結果

表7 純橫蕩計算結果

7) 組合運動 通過模擬不同頻率下首搖幅值為25°漂角為12°的組合運動，可求取水動力導數Xvr，Xrr，Yvvr，Yvrr，Mvvr，Mvrr.計算結果見表10～12.

表8 純首搖橫向力計算結果

表9 純首搖轉首力矩計算結果

表10 組合運動縱向力計算結果

表11 組合運動橫向力計算結果

表12 組合運動轉首力矩計算結果

5 結 束 語

文中以供應船為研究對象，在FLUENT中數值模擬供應船模在不同運動幅值的PMM試驗，進而對比各流場壓力云圖，分析得出供應船水動力達到非線性的振幅范圍約為漂角β大于10°小于30°.最后，設計供應船非線性運動的試驗方案，分別模擬不同頻率時的大漂角斜航運動以及大振幅的純縱蕩、純橫蕩、純首搖、組合運動，擬合得到接近零頻率的非線性水動力導數.該結論和方法可為供應船舶操縱性能的預報提供一定的參考.

[1]吳寶山.第七屆船舶力學學術委員會操縱性學組工作總結報告[C].第七屆船舶力學學術委員會全體會議論文集，2010:51-61.

[2]鄒早建.船舶操縱性研究進展[C].第六屆船舶力學學術委員會全體會議專集，2006:54-67.

[3]ITTC. The Manoeuvring Committee[C]. Proceedings of 25th ITTC-Volume I,2008:143-203.

[4]TURNOCK S R, PHILLIPS A B, Furlong M. Urans simulations of static drift and dynamic manouveres of the kvlcc2 tanker[C]. Proceedings of Workshop on Verification and Validation of Ship Maneuvering Simulation Methods (SIMMAN2008),2008.

[5]PAN Y C, ZHANG H X, ZHOU Q D. Numerical prediction of submarine hydrodynamic coefficients using CFD simulation[J]. Science Direct Journal of Hydrodynamics,2012,24(6):840-847.

[6]張赫，龐永杰，李曄.基于FLUENT軟件模擬平面運動機構試驗[J].系統仿真學報，2010，22(3):566-569.

[7]石愛國，聞虎，李理,等.船舶淺水水動力導數的數值計算[J].中國航海，2011，34(3):69-82.

[8]鄒早建，楊勇，張晨曦.深淺水中KVLCC船體橫蕩運動水動力數值計算[J].水動力學研究與進展，2011，26(1):85-92.

[9]楊勇，鄒早建，樓鵬宇.小振幅平面運動機構試驗粘性流場數值模擬[C].第二十三屆全國水動力學研討會第十屆全國水動力學學術會議文集，2011:190-195.

[10]楊勇.非定常操縱運動船體水動力數值計算[D].上海：上海交通大學，2011.

[11]邢磊.三體船水動力導數及操縱性能預報研究[D].哈爾濱：哈爾濱工程大學，2012.

[12]程捷，張志國，蔣奉兼,等.平面運動機構試驗的數值模擬[J].水動力學研究與進展，2013，28(4):460-464.

[13]TOXOPEUS S L. Using CFD calculations to improve predictions of ship manoeuvres[J]. RINA Developments in Marine CFD,2011(3),13-25.

[14]NOBUAKI S, PABLO M, CARRICA, et al. URANS simulations of static and dynamic maneuvering for surface combatant: verification and validation for forces, moment, and hydrodynamic derivatives[J]. Journal of Marine Science and Technology,2012,17(4):422-445.

[15]李冬荔.粘性流場中船舶操縱水動力導數計算[J].哈爾濱工程大學學報，2010，31(4):421-427.

[16]劉山.基于CFD技術數值模擬平面運動機構試驗[D].武漢：武漢理工大學，2012.

[17]林洪波.供應船的穩性特點[J].天津航海,1987(3):27-29.

[18]賈欣樂，楊鹽生.船舶運動數學模型-機理建模與辨識建模[M].大連：大連海事大學出版社，1999.

[19]趙小仨，羅薇，柯梟冰，等.呆木對供應船線性水動力導數影響[J].船舶工程，2014，36(3):33-36.

[20]趙小仨，羅薇，柯梟冰，等.供應船縱向加速度導數的估算方法[J].中國航海，2014，37(4):88-91.

Numerical Simulation and Research of Nonlinear Hydrodynamic Derivatives

ZHAO Xiaosa1,2)XU Haixiang1,2)

(KeyLaboratoryofHighPerformanceShipTechnologyofMinistryofEducation,WuhanUniversityofTechnology,Wuhan430063,China)1)(Collegeoftransportation,WuhanUniversityofTechnology,Wuhan430063,China)2)

For the problem that the nonlinear dynamics of ships are difficult to define and forecast, this article uses CFD commercial software FLUENT in combination with dynamic mesh to simulate the large amplitude planar motion mechanism test for a supply vessel. The range of nonlinear motion amplitude of the supply vessel is obtained by comparing the pressure nephograms of different scenarios. Besides, the test scheme of the supply vessel’s nonlinear motion is designed. The large drift angle of oblique motion and large amplitude of pure surge, pure sway, pure yaw, combination movement with several frequencies are simulated, respectively. Finally, the nonlinear hydrodynamic derivatives close to zero frequency are calculated.

nonlinear hydrodynamic derivatives; large amplitude pmm test; numerical simulation; supply vessel

2016-12-14

*國家自然科學基金項目 (61301279, 51479158)、中央高校基本科研業務費專項資金項目(163102006)資助

U661.33

10.3963/j.issn.2095-3844.2017.01.014

趙小仨(1989—)：女，工學碩士，實驗員，主要研究領域為船舶水動力研究