基于CFD數值的模擬船/槳相互干擾流動

陳 科，賴明雁

（1.中遠海運發展股份有限公司，上海 200135；2.上海船舶研究設計院，上海 201203）

0 引 言

隨著計算機性能的提升和數值方法的發展，計算流體力學 (Computational Fluid Dynamics, CFD) 在船舶海洋工程領域的應用越來越廣，目前對許多流動問題的預報精度已達到較高的水平。從Tokyo 2015 CFD Workshop的計算結果來看，多家機構的阻力、自航預報精度已在5%的誤差范圍內[1]，在某種程度上能滿足工程設計要求。同時，CFD在耐波性和操縱性數值模擬方面也顯示出較高的可靠性，并可為一些非線性流動現象的數值計算提供參考[2-6]。近幾年，隨著綠色船舶和節能減排概念的深入發展，船舶設計者越來越重視船舶推進效率的提升和新型節能裝置的開發，迫切需要得到能較為可靠地預報船/槳相互干擾問題的數值處理方法。CFD因具有較高的靈活性和巨大的潛力而成為人們研究的熱點。

然而，船/槳相互干擾是一種比較復雜的非線性流動現象。從船體周圍的流動特點來看，艉部附近流線曲率變化較快，常伴隨流動分離現象，導致槳盤面處伴流較為復雜，使螺旋槳的推進性能有別于敞水。同時，螺旋槳的旋轉會產生抽吸現象，影響船體附近的流動，造成阻力增額。船/槳干擾問題的復雜性給CFD計算帶來很大挑戰。為研究船/槳干擾，相關研究人員開發出不同的數值處理方法，這些方法各有優缺點，具體區別體現在對螺旋槳的處理上，主要分為體積力法和船/槳整體建模法2類。體積力法不考慮螺旋槳的幾何形狀，而是利用虛擬圓盤代替螺旋槳，在圓盤所在區域的控制方程中施加體積力源項。體積力法建模簡單、計算效率高，因此得到廣泛應用[7-8]。然而，體積力法的計算精度有限，不能很好地體現出螺旋槳幾何形狀改變對流場的影響，且對槳盤面所在區域的流場預報不準。隨著研究的不斷深入，相關研究人員開始嘗試對船/槳整體建立模型來求解[9]。這類方法中常用的是Moving Reference Frame（MRF）法和滑移網格法。在利用整體建模法計算時，須建立一個包含螺旋槳的小域。在利用 MRF法計算時，螺旋槳所在的區域處在一個相對旋轉的坐標系內，但網格是不動的；在利用滑移網格法計算時，螺旋槳按照相應的轉速旋轉，流場處于瞬態變化中。整體來說，MRF法和滑移網格法各有優缺點，其中MRF法的計算效率較高且穩定，而滑移網格法則能更真實地模擬船/槳之間的非定常流動現象。目前，已有一些學者對這2種方法進行探討，如ZHANG[10]和楊春蕾等[11]從船體阻力、螺旋槳受力及伴流等方面分析MRF法與滑移網格法的區別。這些研究豐富了對船/槳相互干擾現象的認識，但對自航因子預報及螺旋槳梢渦等流動特征方面的分析還比較少。另外，船/槳相互干擾計算的可靠性還需更多的驗證和分析。

本文針對國際標準船模韓國船舶與海洋工程研究所的集裝箱船（KRISO Container Ship，KCS），采用MRF法和滑移網格法對船/槳相互干擾的流動問題進行計算，將計算結果與試驗結果和其他數值結果相比較，從伴流、船體表面壓力、螺旋槳梢渦、船體和螺旋槳受力及自航因子等方面分析MRF法和滑移網格法的預報精度及船/槳之間相互干擾的流動特點，為工程實際應用提供參考。

1 控制方程及數值方法

數值模擬采用商業CFD軟件STAR-CCM+，計算中控制方程采用Reynolds平均Navier-Stokes方程（即RANS方程），湍流模型采用SST k-ω兩方程湍流模型，自由面基于VOF（Volume Of Fluid）方法求解。對流項離散采用二階精度的迎風格式。

控制方程包括連續性方程和動量守恒方程，對于不可壓縮流體，有

式(1)和式(2)中：i和j分別為變量在坐標系中第i方向及第j方向的分量；P為壓力；u為速度矢量；ρ為流體密度；f為質量力。

SST k-ω模型是一種結合k-ε模型和原始k-ω模型的混合方法，湍動能k和湍動耗散率ω的求解方程為

VOF控制方程為

2 數值計算結果及分析

2.1 船型參數

計算模型為國際標模KCS，該模型的原型是韓國KRISO（Koera Research Institute of Ships and Ocean Engineering）設計的集裝箱船，船體幾何形狀見圖1。KCS被國際學術界推薦用來進行CFD驗證和分析，有較為充分的試驗數據作為參考，包括靜水阻力、自航、耐波性及操縱性等。船體和槳的主要參數見表1。KCS實船長230m，設計航速為24kn，模型縮尺比為31.6。本文計算采用的航速為設計航速，對應Fr=0.26。

圖1 KCS船體幾何形狀

表1 KCS船體和槳參數

2.2 邊界條件及網格

在 CFD計算中，計算域的選取不僅影響計算效率，還會影響計算精度。一般來說，計算域和邊界條件的指定依賴于計算問題本身。對于自航計算，由于船體左右流動不對稱，計算域寬度比常規阻力計算大。本文的計算域取為長方體，大小為：船前1倍船長，兩側和底部1.5倍船長。船前兩側和上下面均為來流邊界條件，取x軸負向為壓力出口邊界。計算網格為非結構六面體網格，為更好地捕捉船體附近的流動特征，特別是艉部船/槳相互干擾區域的流動細節，對艉部和螺旋槳附近區域進行加密；為更好地模擬螺旋槳泄渦，對槳葉的導邊和隨邊均進行加密。計算域網格劃分見圖2，網格總數為271萬個，其中船體部分網格為220萬個，螺旋槳所在區域網格為51萬個。體網格類型主要以六面體為主，為更好地捕捉物面附近的流動，在船體和螺旋槳表面生成邊界層網格。

圖2 計算域網格劃分

2.3 計算結果及分析

為分析船/槳干擾問題及自航因子計算，首先對KCS在裸船條件下的靜水阻力進行計算，計算網格配置參數與“2.2”節中的敘述相似。表2給出阻力的計算結果和試驗結果，二者能較好地吻合。圖3為KCS靜水狀態航行時波形，對比了船體附近的自由面興波，可看到通過數值計算得到的自由面波形能與試驗結果較好地吻合，數值計算得到的波峰和波谷位置及波形都能與試驗結果相一致，表明本文采用的數值方法具有一定的可靠性。

表2 KCS靜水阻力計算結果

KCS自航試驗是在實船自航點下進行的，因此在船模尺度計算中需添加一項恒定的強制力，或稱摩擦力修正力（Skin Friction Correction，SFC，以Fc表示），用于修正尺度效應引起的摩擦力誤差。Fc定義為

圖3 KCS靜水狀態航行時波形

式(6)中：分別為船模和實船的摩擦阻力系數；k為形狀因子，對于KCS而言，k+1=1.1；U0為航速；SW為濕表面積。

計算中螺旋槳的轉速與模型試驗一致，即n=9.50r/s。表3給出分別用MRF法和滑移網格法得到的推進及自航因子，并與試驗結果和其他數值結果相對比，其中CASTRO等[12]的計算采用IOWA-CFDSHIP求解器，并在計算中使用DES（Detached Eddy Simulation）模擬湍流流動，螺旋槳旋轉采用over-set網格方式來實現。需說明的是，本文的自航因子計算采用的是強制自航法，即：首先計算 2個轉速（n=9.50r/s和n=9.75r/s）下的自航，記錄相應的強制力、螺旋槳推力系數KT和扭矩系數KQ；然后通過強制力插值得到自航時的螺旋槳轉速、推力和扭矩等信息，插值原則為保證強制力等于試驗的強制力（即30.25N）。自航因子采用等推力法計算得到，需說明的是，螺旋槳敞水數據采用試驗測量結果。

從轉速的計算結果來看，滑移網格法精度較高，MRF法給出的轉速偏高。這主要是因為利用MRF法計算得到的船體阻力較大，因此在進行強制力插值時得到的轉速較大。利用MRF法和滑移網格法得到的KT均與試驗結果相近，KQ相差略大，即通過數值模擬得到的螺旋槳推進效率較低，CASTRO等[12]的計算結果存在同樣的現象。這可能是因為艉流場流動比較復雜，存在邊界層的分離，數值模擬未能很好地捕捉這些流動特征。由表3中的自航因子可知，利用滑移網格法得到的自航因子與試驗結果能較好地吻合，相對而言，利用MRF法計算得到的推力減額偏大，推進效率偏低。對于伴流分數和相對旋轉效率等因子，利用MRF法和滑移網格法得到的結果都能與試驗結果相一致。

表3 利用MRF法和滑移網格法得到的自航計算結果

圖4為靜水阻力計算和自航計算得到的船體表面動壓力分布對比。由圖4可知，在自航狀態下，受螺旋槳抽吸影響，船體表面的動壓力比靜水阻力低，導致自航時船體阻力較大，在功率換算中等價于推力減額。利用MRF法和滑移網格法得到的計算結果都能較好地反映該船/槳干擾的流動現象。

圖5為利用試驗、MRF法和滑移網格法得到的螺旋槳下游x/LPP=0.9941處的伴流對比。從整體上看，通過數值模擬得到的伴流場與通過試驗測量得到的結果相一致；但從細節上看，利用滑移網格法得到的結果與試驗結果比較吻合。這里需指出的是，在利用MRF法計算時螺旋槳不動，因此螺旋槳槳葉的相位會對結果有一定影響。圖6為艉部流場中Q等值面。對比采用MRF法和滑移網格法得到的計算結果可看到，由于MRF法僅在螺旋槳旋轉域內部施加體積力，因此得到的螺旋槳梢渦在小范圍內呈旋轉狀態，而在距離螺旋槳稍遠處隨來流向下游水平泄出，在螺旋槳計算子域和外域交界處存在渦的非物理變形。相對而言，利用滑移網格法得到的渦結構與物理直觀更為接近。

圖4 船體表面動壓分布

圖5 螺旋槳下游的軸向速度分布（x/LPP=0.9941)

圖6 艉部流場中Q等值面

3 結 語

本文基于MRF法和滑移網格法對國際標模KCS的船/槳干擾問題進行數值研究。從計算結果來看，螺旋槳的存在會改變船體表面的壓力分布，造成船體阻力增大，即出現推力減額現象。利用MRF法和滑移網格法都能較好地反映該流動現象。從自航因子結果來看，滑移網格法給出的自航因子能與試驗結果和其他數值結果相吻合，而采用MRF法所得結果的誤差較大，這主要是因為利用MRF法得到的船體阻力偏大。從伴流結果來看，滑移網格法與試驗吻合較好。此外，在螺旋槳附近，利用MRF法和利用滑移網格法得到的槳葉梢渦形狀基本一致；在遠離螺旋槳處，利用滑移網格法得到的計算結果更符合物理直觀。

【 參 考 文 獻 】

[1] LARSSON L, STERN F, VISONNEAU M, et al.Proceedings of CFD workshop Tokyo 2015[C].National Maritime Research Institute, 2015.

[2] CARRICA P M, FU H, STERN F.Computations of self-propulsion free to sink and trim and of motions in head waves of the KRISO Container Ship (KCS) model[J].Applied Ocean Research, 2011, 33: 309-320.

[3] SHEN Z R, WAN D C.RANS computations of added resistance and motions of a ship in head waves[J].International Journal of Offshore and Polar Engineering, 2013, 23 (4): 263-271.

[4] SHEN Z R, YE H X, WAN D C.URANS simulations of ship motion responses in long-crest irregular waves[J].Journal of Hydrodynamics, 2014, 26 (3): 436-446.

[5] WANG J H, ZHAO W W, WAN D C.Free maneuvering simulation of ONR Tumblehome using overset grid method in naoe-FOAM-SJTU solver[C]//31stSymposium on Naval Hydrodynamics, September 11-16, 2016, California.

[6] HU C H, SUEYOSHI M, KASHIWAGI M.Numerical simulation of strongly nonlinear wave-ship interaction by CIP/Cartesian grid method[J].International Journal of offshore and Polar Engineering, 2010: 81-87.

[7] 吳召華，陳作鋼，代燚.基于體積力法的船體自航性能數值預報[J].上海交通大學學報，2013, 47 (6): 943-949.

[8] SIMONSEN C D, STERN F.RANS maneuvering simulation of esso osaka with rudder and a body-force propeller[J].Journal of Ship Research, 2005, 49(2):98-120.

[9] CARRICA P M, SADAT-HOSSEINI H, STERN F.CFD analysis of broaching for a model surface combatant with explicit simulation of moving rudders and rotating propellers[J].Computers & Fluids, 2012, 53: 117-132.

[10] ZHANG Z R.Verification and validation for RANS simulation of KCS container ship without/with propeller[J].Journal of Hydrodynamics, 2010, 22 (5): 932-939.

[11] 楊春蕾，朱仁傳，繆國平，等.基于CFD方法的船/槳/舵干擾數值模擬[J].水動力學研究與進展，2011, 26 (6): 667-673.

[12] CASTRO A M, CARRICA P M, STERN F.Full scale self-propulsion computations using discretized propeller for the KRISO Container Ship KCS[J].Computers & Fluids, 2011, 51: 35-47.