基于URANS 模型的槽道側推器水動力性能數值計算方法

楊星晨，葉金銘，肖昌潤，柯 林，黃 路

(1. 海軍工程大學 艦船與海洋學院，湖北 武漢 430033；2. 國家海洋技術中心漳州基地籌建辦公室，福建 廈門 361007)

0 引 言

船舶在離靠碼頭時有多種方式，其中一種簡便有效的方式是采用側推器，側推器的使用彌補了航速較小時主推進器舵效不足的問題，這使得大型船舶不需要拖船也能自行離靠港。此外，隨著時代的發展出現了功能多樣的載人潛水器、AUV 和ROV，其在海底航行時對操縱性和定位有著更高的要求，因此除了在前后安裝推進器，通常也要安裝側推器，“蛟龍”號深海載人潛水器，為達到深海操縱性要求，在首尾部安裝了6 臺側推器。因此研究側推器的水動力性能具有重要意義。

槽道側推器是出現最早也是目前應用最廣泛的側推器，其主要特點是在船首或船尾橫向開槽，在槽道內部通過安裝T 形軸聯動螺旋槳或者直接采用輪緣驅動式葉片。對于槽道側推器研究主要采用理論和試驗的方法。Taniguch[1]采用敞水模型試驗對定距槳和調距槳的槽道側推器水動力性能進行了較為詳細的研究，采用操縱試驗檢驗側推器的轉向能力。該試驗是目前公布的最為詳細的側推器敞水性能試驗。沈國鑒等[2]根據模型試驗的結果得出了槽道側推器的設計圖譜，據此可初步確定所需功率和轉速，對于側推器的設計有很好的指導性作用。Ridley, Donald E.[3]研究了隧道入口結構對推力產生和船體阻力增強的影響。Beveridge,John L[4]的研究報告中提出船首側推器的水動力性能計算經驗公式，并從設計角度對管道布置、管道形狀和葉輪的選擇做了研究。

20 世紀90 年代，CFD 開始廣泛進入工業設計領域，目前已成為推進器研究設計的方法之一[5]。熊鷹等[6]基于混合網格的RANS 方法準確預報螺旋槳敞水性能。相對于敞水螺旋槳，側推器CFD 方法研究目前較少，近年來國內外學者利用商用軟件CFX，STAR–CCM+，Fluent 等對側推器的不同特征尺寸、工況，船速下的水動力性能和流場特性進行了模擬計算[7–10]。郁程等[11–12]基于MRF 模型和混合網格技術對試驗中[1]部分調距槳做了數值計算，水動力性能誤差在10%以內。Yukun 等[13]采用試驗結合數值仿真的方法對船首側推器進行了研究，主要研究單個首側推器在不同船型，不同來流、不同開口圓角影響下的水動力性能，還研究了多個首側推器互相影響以及開口圓角對葉片壓力分布不均勻的改善作用,利用CFD 方法從流場方面解釋了各參數對側推器水動力影響的原因。

目前采用CFD 方法對側推器仿真計算仍然有很多不確定性，本文對照公開的試驗條件和數據[1]，基于STAR-CCM+中的RANS 方法和SSTk-ω模型，采用高質量結構網格，對調距槳進行仿真計算。通過對計算域、網格數量和時間步長進行比較，驗證CFD 方法的可靠性，并對不同葉數的調距槳側推器流場特性和水動力性能做了計算分析。

1 研究對象

本文以文獻[1]試驗中裝有1308 槳和1308-B 槳的槽道側推器為研究對象，1308 槳是1 個4 葉平板調距槳，1308-B 槳是1 個3 葉槳，其葉片參數與1308 相同。葉片的輪廓如圖1 所示。

圖1 1308 槳葉片輪廓Fig. 1 Blade profile of 1308 propeller

1308 槳和1308-B 槳葉的主要參數如表1 所示。

文獻[1]試驗中對1308 槳和1308-B 測試時使用的均是名義船體，船體槽道壁直徑為203 mm，槽道進出口采用10 mm 的圓弧倒角。試驗的布置如圖2 所示。

圖2 側推器模型試驗布置圖Fig. 2 Arrangement of model test

據此采用Catia 建立試驗的幾何模型，忽略露出水面的試驗體和測量設備，僅對沒入水下的名義船體和螺旋槳建模，如圖3 所示。

圖3 試驗模型的幾何圖Fig. 3 Geometry of the test model

試驗[1]中數據采用無量綱化處理，由于側推器工況屬于系柱工況，進速系數為0，且名義船體推力占比較大，對此采用下列公式處理數據。

式中：CT，CF，CTF，分別為螺旋槳，船體，側推器整體的推力系數；CQ為 螺旋槳扭矩系數； η為側推器整體效率；T，F分別為轉子推力和船體推力 kg；Q為轉子扭矩， kg·m 。 ρ 分別表示水的密度， k g·s2/m4；n為轉子轉速， r/s，D為轉子直徑，m。試驗中測試發現當結果以無量綱化方式表示時，轉速對結果基本沒有影響，本文計算時轉速采用20 r /s。

2 計算方法

2.1 計算域劃分和邊界設置

試驗中沒有給出水池的尺度，由于側推器計算工況為系柱工況，為保證計算的收斂性和準確度，需劃分大計算域。計算域初步設置為長度50D，寬度90D，深度40D的長方體，對結果影響最大的主要是長度，在后續計算中改變長度來檢驗計算域的影響。將名義船體置于計算域中央，其上表面與頂面平齊。由于槽道中心距離水面1.25D，螺旋槳浸深高于0.76D時可以不計算自由液面影響[4]，因此外域頂面設置為對稱平面，外域其他面設置為壓力出口，其他表面均使用壁面條件處理，如圖4 所示。

圖4 計算域邊界條件Fig. 4 Boundary conditions of computational domain

圖5 整體外域網格Fig. 5 The mesh of outer domain

圖6 旋轉域網格Fig. 6 The mesh of rotor domain

圖7 槽道口網格分布Fig. 7 The mesh of channel mouth

圖8 葉片網格分布Fig. 8 The mesh of blade

2.2 網格劃分

將計算域分為旋轉域和外域2 部分，采用結構網格對2 個區域進行劃分，區域之間采用交界面進行數據傳遞，在劃分過程中在槳葉表面單獨劃分一層薄層用于邊界層網格數量的控制，由于船體的壓差主要來自于槽道進出口，此處的網格也需要加密。在后續計算中，主要改變這2 處的網格數量來驗證網格數量對于數值計算的影響。

2.3 計算模型

采用STAR-CCM+軟件RANS 方法，在流體不可壓縮的條件下，流場的連續方程和動量方程如下：

式中：ui，uj為各速度分量的時均值 (i,j=1, 2, 3)；P為壓力的時均值； ρ為流體的密度； μ為流體的動力粘性系數；gi為 重力加速度分量；為雷諾應力項。

為封閉控制方程引入SSTk-ω湍流模型，該模型可以對邊界層內的流動和邊界層外的湍流進行更好的模擬，因此能夠對翼型結構和含有分離流的流場進行準確的計算，但這需要較高的網格質量，符合本文計算需求。對于螺旋槳旋轉的模擬目前主要有MRF 和Sliding Mesh 2 種方法，本文先采用MRF 方法對旋轉域的螺旋槳進行定常模擬，待計算穩定后再采用Sliding Mesh 進行非定常計算。在計算中對非定常的時間步設置做了驗證計算。

3 計算結果與分析

3.1 計算域、網格數量、時間步長對側推器仿真計算的影響

對螺距比為0.8 時的1308 槳槽道側推器水動力性能進行計算用來檢驗各因素的影響。首先對計算域的影響進行驗證，構建出3 個長度分別為30D，50D和80D的計算域，采用相同方法和過度率構建旋轉域和外域的結構網格，取相同時間步長，即單位時間旋轉2°，代入計算模型進行計算，結果如表2 所示。

表2 計算域長度對計算結果影響情況Tab. 2 The influence of calculation domain length on calculation results

可看出計算域長度取30D時誤差較大，取80D時誤差較小，但其誤差與50D時的誤差差別不大。因此為減少計算時間，取50D作為計算域的長度，在后續計算中都使用該計算域。

對網格依賴性進行驗證，旋轉域內的葉片表面第1 層網格高度分別取0.0005 mm，0.001 mm，0.005 mm。在能夠反映出葉片形狀的前提下，由低到高改變弦長方向網格數目和過度率，最終生成3 種旋轉域網格；選長度為50D的外域網格，采用相同的方式改變槽道進出口附近的網格數目，并保證區域網格過度平緩，得到3 種外域網格，與旋轉域組合成3 套網格。網格劃分情況如表3 所示。

表3 各組計算域網格劃分Tab. 3 Grid division for computational domain meshes in each group

各組網格的計算結果和試驗進行對比，如表4 所示。

表4 網格劃分對計算結果影響Tab. 4 The influence of the grid on the calculation results

計算結果顯示，當轉子域和外域網格數量增至1 240 萬、1 800 萬時，繼續增加網格數量對計算結果的影響較小，此后產生的誤差與網格數量關系不大。因此采用Mesh2 網格形式構建后續研究對象的網格。

轉子單位時間步長內旋轉角度分別設置為1°，2°，5°時，內部迭代采用10 步，采用二階時間離散，計算結果如表5 所示。

表5 時間步長對計算結果影響Tab. 5 The influence of time step on calculation results

結果顯示時間步長采用5°時的誤差最小，時間步長1°和2°的計算結果接近。分析其原因，5°誤差小主要是因為時間步長過大對于非定常流場的計算粗糙，為了節省計算時間，在后續研究中采用2°時間步長。

3.2 不同螺距比的側推器水動力計算結果

試驗中沒有明確調距方式，在數值計算中采用工程實踐中的調距方式，從設計螺距比 γ=0出發，調整到虛擬螺距角 γ′， 虛擬螺距變化為 Δγ=γ′-γ，令槳葉繞設計轉動軸轉動 Δα=Δγ ，并認為轉動 Δα后，新調距槳是具有螺距角 γ′及對應螺距比(P/D)的槳[14]。據此對1308 槳和1308-B 槳葉的0.8，1.0，1.2 螺距進行建模和計算，計算域、網格和計算設置采用已驗證的結果進行構建，計算結果如圖9 所示。

圖9 1308 和1308-B 槽道側推器不同螺距的水動力系數Fig. 9 Hydrodynamic coefficients for different pitches of 1308 and 1308-B side thrusters

對比CFD 計算結果和試驗數據，發現CFD 計算的水動力性能與試驗基本相符，隨著螺距的增大，水動力系數呈現上升的趨勢。不同螺距的螺旋槳推力系數CT誤差在-2%～3%左右，船體的推力系數CF誤差在4%～6%，1 0CQ的誤差在-2%～4%左右，可見誤差主要來自于船體推力的計算。

根據試驗報告中采用的效率計算公式，對計算結果進行換算，CFD 計算效率和試驗測試效率結果如圖10 所示。

圖10 1308 和1308-B 不同螺距下的效率Fig. 10 Efficiency of different pitches for 1308 and 1308-B

可以看出CFD 的效率計算結果與試驗基本相同，1308-B 和1308 槳隨著螺距的增加效率逐漸增加。不同葉片數的敞水螺旋槳效率在系柱情況下基本相同[15]，與敞水螺旋槳不同的是試驗中葉片數多的槽道側推器效率較高，這主要是因為葉片數量的增減影響了槽道內部的流速。雖然葉型相同，而1308 槳的槽道側推器流速較快，因此具有更高的效率點。

3.3 不同葉數的側推器流場特性分析

取槽道中部橫截面的軸向速度進行分析，如圖11所示。可以看出，1308 槽道側推器不同螺距下的流速均稍高于1308-B 槽道側推器的流速，這使得1308 槽道側推器的推進效率略高于1308-B 的效率。另外可以看出1308 槽道側推器的流場不均勻性更大，其產生稍隙泄流的流速更大。

取1308 槽道側推器和1308-B 槽道側推器不同螺距下葉片的吸力面壓力分布進行對比分析，如圖12 所示。可以看出，葉片的壓力分布呈“波浪式”分布，低壓集中于高半徑的導邊附近.這主要是因為槽道側推器采用的葉片是無拱度的葉片，且該槳葉剖面的螺距角均相同。這使得葉片高半徑處的負載較重，因此在導邊高半徑部分的低壓明顯。對比可知1308 槽道側推器不同螺距下的吸力面最低壓力都低于1308-B，特別是螺距P/D=1.2 時，最低壓力達到-270 kPa 左右，容易產生明顯的面空化問題。

4 結 語

本文驗證CFD 方法計算槽道側推器水動力性能的可靠性，研究螺距比和葉片數對側推器水動力性能的影響，并對流場特性進行了分析，得到以下結論：

1）槽道側推器的水動力性能計算對計算域的要求較高，需選擇較長的計算域，長度應不小于30D，采用結構網格時在對葉片、葉間隙、槽道口附近網格進行加密后，選擇2°以內時間步長的非定常計算可以較為準確計算出槽道側推器的水動力性能。本文計算誤差主要來自于船體推力，誤差在6%以內。

2）1308 和1308-B 槽道側推器的水動力系數隨螺距變化趨勢相同，但1308 的效率較高，主要是因為相同螺距下1308 槳的槽道內部流速較快，使得葉片處于更高的效率點。

3）采用等螺距和無拱度葉型的槽道側推器葉片壓力分布呈“波浪式”分布特點，其低壓集中于導邊的高半徑部分，容易產生空化問題。葉梢部分的流動不均性較大，葉頂間隙有明顯的回流速度。

總體而言，CFD 計算中URANS 模型能夠準確預報槽道側推器水動力性能和流場特性，可以作為設計研究側推器的輔助工具。