噴水推進高速三體滑行艇數值自航研究

仝 博，王永生，楊瓊方，劉承江

（海軍工程大學 動力工程學院，武漢 430033）

0 引 言

相比螺旋槳推進，噴水推進高速船具有高效率、抗空化、低噪聲等優勢。越來越多的水面船舶開始采用噴水推進，尤其在高性能船上得到越來越多的應用[1]。對噴水推進艦船快速性的研究通常采用試驗方法，但由于缺乏統一的標準和規范，其研究發展較為緩慢。近年來，隨著計算流體力學的不斷發展，采用CFD方法研究噴水推進船的流場特性受到越來越多的船舶工作者的青睞。

Keegan等人[2]對某高速噴水推進艇進行“船體+泵”數值計算，計算主要集中在艇體對進流的影響，但未考慮自由液面。Takanor等人[3]實現了“船體+泵”的兩相流數值計算，比較了噴泵裝船后船體阻力的變化，噴泵的作用以體積力代替。Tomohiro[4]在碩士論文中采用高精度的RANS方法對某噴水推進高速艇進行粘性流場的數值模擬，并對球鼻首和流道形狀進行優化，證明了CFD技術進行優化設計的能力。Manivannan等人[5]采用單相Level set方法處理自由面，以體積力模型模擬噴泵作用，對某雙體船模型進行“船+泵”數值自航，重點是對流道內流動進行分析，阻力、流量和推力計算結果與試驗吻合良好。Young等人[6]采用Star-CCM+軟件研究了表面效應船與噴水推進器之間的相互作用以及進流對泵性能的影響。其中，噴水推進器采用真實模型。國內孫存樓[7]采用代替實際泵的邊界條件法對"船體+泵"進行了兩相流的CFD模擬，變非定常計算為定常計算，節約了數值計算時間。并指出中高速噴水推進船推力減額為負值是因為進流破壞了船體邊界層，導致摩擦阻力顯著下降。彭云龍等人[8]運用CFD方法計算了某高速摩托艇“噴泵+流道+格柵+船體”的流場特性，并對流道進行幾何優化改善了渦流，提高了推進效率，計算未考慮自由液面影響。

以上計算大都基于單相流計算，或者以體積力模型代替噴泵的作用，這些簡化在噴水推進船自航流場細節的描述上必然有所欠缺。本文基于CFD方法，采用ITTC推薦的等推力系數法對某高速三體滑行艇進行真實的“船+噴泵”模型的兩相流數值自航。首先，對裸艇阻力進行數值計算，并與試驗比較，驗證CFD方法預報滑行艇阻力性能的精度。然后，對“滑行艇+噴泵”模型進行數值自航，將自航結果換算到實尺度艇，判斷該噴泵能否滿足快速性要求；接著對高速三體滑行艇推力減額為負值作了機理上的分析。最后，對艇底水氣分布進行研究，說明CFD方法模擬滑行艇精細流場的準確度。

1 數值計算方法

滑行艇高速航行時其自由液面非線性程度極高，且艇首與水面接觸的位置飛濺現象嚴重，此時總阻力中飛濺阻力和興波阻力占主要成分。能否很好地捕捉自由液面和飛濺現象直接決定了快速性預報的準確性。

本文的CFD計算對自由液面采用VOF分離求解，即將水和空氣分開求解，可以真實地捕捉到滑行艇附近的復雜流動和飛濺效果[9]。VOF方法通過定義一個流體體積函數F，即水在網格單元中所占的體積比，在對計算域進行離散求解時計算每個網格單元中的水的體積分數，自由液面定義為F=0.5的等值面。其輸運方程為：

根據F的定義，計算區域中的密度和粘性系數可用如下方程統一表示為：

其中：φ（ x，t）代表 ρ（ x，t）、 μ（ x，t）或其它物性參數。

2 研究對象

本文研究對象為某三體滑行艇，以下簡稱M艇。模型的阻力試驗在某研究所進行。船模重166.7 kg，縮尺比1：6。對“艇+噴泵”的數值計算也在模型尺度下進行。專門針對該滑行艇采用三元設計理論和參數化方法設計了噴泵和進水流道，運用已被驗證的CFD方法計算了噴泵和進水流道敞水性能，噴泵敞水效率和進水流道敞水效率分別為91.5%和96.1%。將該噴泵和進水流道裝到艇后，幾何模型如圖1所示。

圖1 “艇體+噴水推進器”模型Fig.1 The model of ‘hull+waterjet’

圖2 計算域Fig.2 Computational region

3 艇體阻力數值計算

3.1 計算域的選取

對艇體周圍建立長方體計算域。計算域前端距艇首1倍的艇長，后端距艇尾3倍艇長，頂部距離龍骨線1倍艇長，底部距離龍骨線1.7倍艇長，側邊距離對稱面1倍艇長。坐標原點位于尾部最低點，模型在計算域中的位置及計算域形狀如圖2所示。

圖4 艇體周圍和自由液面加密區網格Fig.4 Refined mesh around hull and the free surface

3.2 網格劃分

本文采用了以六面體為主的切割體網格以及在近壁面設置棱柱層網格。切割體網格通過設置面網格和體網格尺寸以及局部體網格控制來完成，這個過程相對簡單省時。艇體周圍設置棱柱層網格來模擬邊界層流動，棱柱層總厚度設為20 mm，層數8層，延展比1.4。在艇體周圍和自由液面附近進行雙層加密來控制網格過渡。邊界層網格如圖3所示。艇體周圍和自由液面網格如圖4所示。

3.3 邊界條件設置

本次計算采用瞬態分離隱式求解，時間步設為0.000 1 s，兩相流選用VOF方法處理自由液面。艇體表面設為不可滑移壁面，進流面設為速度進口，出流面為分散流出口，計算域頂部和底部為滑移壁面，側邊和對稱面設為對稱邊界。湍流模型選用可實現的k-ε模型。

3.4 阻力計算結果和分析

對M艇模型進行拖曳試驗，拖點位于重心處，拖曳過程中，隨著航速變化艇吃水和縱傾角也會隨著變化，這些變化量都可通過儀器測量到。對該艇的阻力計算在模型尺度下進行，采用固定航行姿態模擬高速航行時的流場。計算工況選取傅如德數在1＜Fr＜1.8 范圍內的五個速度：6 m/s、7 m/s、8 m/s、9.07 m/s和10 m/s，9.07 m/s為設計工況對應模型尺度下的速度。將計算的阻力值按以下公式進行無量綱處理，得到阻力系數。

其中：ρ為水的密度，B為艇寬。阻力和試驗對比結果如圖5所示。

從計算結果看，隨著速度增加，阻力變化規律與試

驗趨勢一致，低航速時，誤差最大，為8.3%，高航速時誤差最小，為0.5%。即使高速航行時存在飛濺、破波現象，阻力數值計算誤差仍能夠達到較高的精度，并且設計點誤差控制在5%以內，說明本文的CFD方法預報高速滑行艇靜水航行的阻力性能是可行的，且準確性較高。

圖5 阻力計算結果Fig.5 Computed resistance

4 “滑行艇+噴泵”數值自航

針對該M艇只進行了阻力試驗，未進行船模自航試驗，因此本文采用CFD方法預報M艇和噴水推進器共同作用下的粘性流場，求取自航因子，預報該噴泵能否達到設計航速要求。

4.1 計算網格

該三體滑行艇采用雙泵推進，其流場具有對稱性，為了節省計算時間，只取半船和單個噴泵進行數值計算。

圖6 艇尾流道網格和泵網格Fig.6 The mesh of duct and pump

圖7 縱剖面網格Fig.7 The mesh of longitudinal section

“艇+噴泵”模型網格劃分為三個部分：“艇+進水流道+噴口”部分網格，網格數量234萬；導葉部分網格，網格數量100萬；葉輪部分網格，網格數量110萬。葉輪和導葉棱柱層網格厚度設為2 mm，層數5層，延展比1.3。三部分網格最后組裝在一起，并設置交界面。艇尾流道網格和泵網格如圖6所示。考慮到前駐點位置艇體和水流碰撞劇烈，噴口射流速度高，落入水中影響艇尾自由面流動，因此專門在這兩處進行網格加密，如圖7所示。

4.2 自航結果

對“艇體+噴泵”數值自航時，采用穩態參考坐標系（MRF）方法，通過設置旋轉坐標系來實現葉輪的旋轉，此時葉輪相對參考坐標系靜止。

由于計算在模型尺度下進行，因此必需滿足一定的相似準則。由于噴水推進器的流道進口與船底齊平，處于船體邊界層內，進口速度分布受船體邊界層影響較大，不能忽略尺度效應，因此不能按照“船+槳”自航所要求的進速系數相等。根據24屆ITTC會議推薦的噴水推進器船模自航試驗相似條件之一：等推力系數法，即保證模型和實尺度泵的推力系數相等。推力系數CT定義為

其中：VS為航速，A為噴口面積。

認為推力減額無尺度效應，則摩擦阻力修正值為

其中：Rm為模型阻力，Rs為實艇阻力，λ為縮尺比。這與螺旋槳船模自航試驗摩擦阻力修正值相同。該M艇在設計航速時實船阻力為67.2 kN，模型阻力317.4 N，計算得到摩擦阻力修正值為6.27 N。

該泵設計轉速為ns=1 315 rpm。若按進速系數相等，換算到模型泵轉速為53.7 rps，以此為初始轉速對“艇+泵”模型進行數值自航。得到計算結果后對轉速進行調整，多次反復計算最終確定T+FD-R=0時轉速為nm=50.3 rps。自航結果見表1。

表1 自航計算結果Tab.1 Calculated results of self-propulsion

由于模型和實尺推力系數相等，通過換算可求得實尺推力為59.7 kN，假定模型和實尺推力減額相等，則實尺推力可克服的阻力為67.96 kN，該值大于實艇阻力67.2 kN，因此該噴水推進器可以滿足快速性要求。圖8反應的是設計點時艇尾部的波形和噴泵射流。

4.3 負推力減額機理分析

從表2中可看出，自航模型的摩擦阻力比裸艇的摩擦阻力小，但減小不多，而壓阻力卻減小了21%，這是導致總阻力下降12%的主要原因。摩擦阻力變化不大原因主要有兩個：一是因為本次計算未考慮滑行艇裝上噴水推進器后的姿態變化，濕表面積和吃水基本沒改變，所以粘性引起的摩擦阻力也變化不大；二是雖然進流破壞了船底邊界層，但其影響僅集中在流道進口附近，對粘性阻力影響不大。

圖8 設計點時尾部波形和噴泵射流Fig.8 The wave profiles and jet flow in the wake of craft at design point

表2 M艇阻力成分Tab.2 The resistance of M boat

為了探尋壓阻力明顯下降的原因，取y=0.26L縱剖面處的艇底壓力分布，如圖9所示。艇中部到艇首范圍內裸艇和自航模型的壓力分布是非常接近的，而在距尾板約0.9 m的區間內（流道進口所在區域），自航模型的艇底壓力出現了波峰和波谷兩個極值，與裸艇壓力分布相差較大。從圖10自航模型的艇底壓力分布圖可以看出，出現波谷是因為流道進口前方處水流由于泵的抽吸作用得到加速，此處出現低壓區；出現波峰是因為進水流道進口兩側由于泵的抽吸，來流縱向速度銳減為零，并改變了速度方向，被吸入流道內，此處形成了高壓區。該高壓區的壓力絕對值比低壓區的壓力絕對值大，這樣導致首尾壓差阻力減小，在摩擦阻力變化不大的情況下，總阻力下降，推力減額為負值。

圖9 y=0.26L處縱剖面艇底壓力分布Fig.9 The pressure distribution on the bottom at the longitudinal section of y=0.26L

圖10 自航模型艇底壓力分布Fig.10 The pressure distribution on the bottom of self-propulsion model

4.4 艇底水氣分布

船舶航行時水線面以下的船體表面不與空氣接觸，全由水覆蓋，但在裸艇阻力計算和自航模擬時均發現一個不真實的現象：艇底表面存在不合理的水氣分布，這與Brizzolara[10-11]的模擬結果相似。圖11反映的是速度為6 m/s和9 m/s時裸艇體表面水氣分布。紅色表示水，藍色表示空氣，介于兩者之間的顏色表示氣水混合物。從圖中可以看出，兩個速度下艇底表面都出現了不正常的水氣分布，尤其在龍骨線附近，氣水混合層由前駐點（即艇首產生破波的位置）開始，沿著龍骨一直延伸到艇尾板。速度低時氣水層較窄，速度高時氣水層較寬。

圖11 不同速度時艇底水氣分布Fig.11 The VOF distribution on the bottom in different speeds

圖12 速度為9 m/s時加密前后艇底水氣分布Fig.12 The VOF distribution on the bottom at speed of 9 m/s before and after refinement

推斷氣水混合層的出現可能是由于高速航行時前駐點產生破波，流動劇烈，自由液面嚴重扭曲，現有的網格數量不足以捕捉這種非線性程度較高的自由液面流動。因此在前駐點和龍骨附近進行網格加密，將整體網格數量由100萬增加至287萬。圖12反應的是速度為9 m/s工況下，網格加密前后艇底水氣分布情況。從圖中可以看出，前駐點和龍骨附近網格加密后艇底氣水層變窄了，說明網格局部加密對改善氣水層是有效的，但隨著進一步加密，氣水層并沒有完全消失。

在對前駐點和龍骨附近網格加密后，總阻力誤差從7%降為4%，這是因為艇底氣水層的存在一定程度上影響了邊界層內的流動，造成邊界層內的速度場和壓力場的計算不準確，進而影響摩擦阻力的計算精度，這可能是滑行艇阻力計算精度難以控制的原因之一。

圖13為網格加密前后噴水推進器縱剖面水氣分布。從圖中可看出，網格加密前，一些空氣被吸入流道，聚集在背部，經過葉輪旋轉和導葉整流后從噴口噴出，使得噴口射流中夾雜著空氣。這勢必會導致流量不足，泵的做功能力下降，揚程、功率、推力和航速均會達不到設計值。在前駐點進行網格加密后，流道內部則充滿了水，此時噴泵的性能可以完全發揮出來。因此得出，艇底不合理的水氣分布會影響噴泵的宏觀特性，通過局部網格加密改善艇底水氣分布后，可以避免這種不利影響。

圖13 加密前后噴水推進器縱剖面水氣分布Fig.13 The VOF distribution on the longitudinal section of waterjet

5 結 論

基于CFD方法，數值計算了某噴水推進高速三體滑行艇阻力及推進性能，對自航結果和流場細節的模擬進行分析，有以下結論：

（1） 在高傅汝德數（1＜ Fr＜1.8）工況下，隨著航速增加，三體滑行艇阻力變化趨勢與試驗一致，最大誤差8.3%，最小誤差0.5%，說明本文的數值方法預報滑行艇的阻力性能可以達到較好的精度。

（2）采用等推力系數法，在模型尺度下進行“滑行艇+噴泵”的數值自航，并將結果推算到實尺度艇，結果表明該噴泵可以克服實艇阻力，滿足快速性要求。說明基于CFD技術的模型自航方法可快速準確地檢驗噴水推進船能否達到設計航速。

（3）噴水推進滑行艇高速航行時，進水流道進口處存在高壓區和低壓區，而且高壓區壓力絕對值大于低壓區壓力絕對值，尾部總體壓力較裸艇航行時大，這使得首尾壓差阻力減小，這是推力減額為負值的主要原因。

（4）考慮自由液面情況下對滑行艇進行自航模擬，艇底會出現不正常的水氣分布，這會影響到滑行艇的阻力性能和噴泵的推進性能。通過局部網格加密可以避免空氣進入進水流道，從而避免對噴泵造成不利影響，但艇底氣水層難以完全消除，這可能是CFD方法預報滑行艇阻力精度難以控制的原因之一。因此作者認為，采用較優的網格策略可以對三體滑行艇進行較為準確的快速性預報，但在局部流場細節的模擬精度有待進一步研究和試驗驗證。

參 考 文 獻：

[1]丁江明.船舶噴水推進技術國內外研究與應用現狀[C]//船舶水動力學學術會議論文集.中國 西安,2013.

[2]Keegan D,Martin D.Use of RANS for waterjet analysis of a high-speed sealift concept vessel[C]//Kourosh K.First International Symposium on Marine Propulsors.Norway,2009.

[3]Takanori H,Kunihide O.Numerical simulation of flow around a waterjet propeller ship[C]//First International Symposium on Marine Propulsors.Norway,2009.

[4]Tomohiro T.Simulation based design for high speed sea lift with waterjet by high fidelity RANS approach[D].Iowa：The University of Iowa,2010.

[5]Manivannan K,Svetlozar G.Numerical and experimental evaluation of waterjet propelled delft catamarans[C]//11th International Conference on Fast Sea Transportation,September 2011.Honolulu,USA,2011.

[6]Young Y L,Savander B R.Numerical investigation of the impact of SES-Waterjet interactions and flow non-uniformity on pump performance[C]//11th International Conference on Fast Sea Transportation,September 2011.Honolulu,USA,2011.

[7]孫存樓.噴水推進器與螺旋槳混合推進艦船推進性能研究[D].武漢：海軍工程大學,2011.Sun Cunlou.Research on propulsion performance of hybrid propulsion ship of waterjet and propeller[D].Wuhan：Naval University of Engineering,2011.

[8]彭云龍,王永生.高速摩托艇噴水推進器性能分析及其改進[J].大連海事大學學報,2013,39(2)：13-17.Peng Yunlong,Wang Yongsheng.Performance analysis and its implementation of water jet propeller equipped on motorboat[J].Journal of Dalian Maritime University,2013,39(2)：13-17.

[9]王 碩,蘇玉民.滑行艇靜水直航及波浪中運動的數值模擬[J].華南理工大學學報(自然科學版),2013,41(4)：119-126.Wang Shuo,Su Yumin.Numerical simulation of planning crafts sailing in calm water and in waves[J].Journal of South China University of Technology(Natural Science Edition),2013,41(4)：119-126.

[10]Brizzolara S,Serra F.CFD simulation of planing hulls[C]//Seventh International Conference on High-Performance Marine Vehicles,13-15 October 2010.Melbourne,Florida,USA,2010.

[11]Villa D,Brizzolara S.Systematic CFD analysis of flaps/interceptors hydrodynamic performance[C]//10th International Conference on Fast Sea Transportation.Athens,2009.