基于分離渦模擬的對轉舵槳水動力性能數值分析

王志勇 范佘明 孫 群 張晨亮

（1.中國船舶及海洋工程設計研究院 上海200011；2.上海市船舶工程重點實驗室 上海200011）

引 言

從幾何外形上，對轉舵槳主要由立柱、下齒輪箱體（從水動力角度舵槳與吊艙槳一致，故下齒輪箱體也可以稱為吊艙），以及立柱前后各一個且旋向相反的兩個螺旋槳組成。對轉舵槳兼具全回轉舵槳和對轉槳的特點，能夠繞立柱軸360°旋轉，具有良好的操縱性。較常規舵槳而言，對轉舵槳的后槳能夠吸收前槳周向旋轉能量，具有較高的推進效率。由于整個推進器的載荷由前后槳分擔，在設計時可以減小兩個槳的直徑，因此對轉舵槳具有良好的空泡和噪聲性能，現已廣泛應用于輪渡等船型上。對海洋平臺而言，風浪條件常常劇烈變化，這就要求海洋平臺作業時必須具有良好的定位能力。對轉舵槳結合動力定位技術能夠很好地滿足工程實踐需求，因此，對轉舵槳具有廣闊的應用前景。

對轉舵槳的幾何結構十分復雜，因此各部件之間的相互作用也非常復雜。但由于專利保護的關系，目前關于對轉舵槳的研究資料非常少。舵槳和吊艙推進器在水動力外形上相似，兩者之間的主要差別僅在于電機的布置方式不同，因此，關于吊艙推進器和常規對轉槳的研究結果都具有參考價值。

近年來，隨著對船舶快速性和節能減排的要求越來越高，研發設計高效的推進器一直是人們研究的重點。在吊艙推進器方面，周劍［1］進行了數值模擬和模型試驗，結果發現吊艙和立柱導致了7.9%的推進效率損失。祝志超［2］對雙槳式吊艙推進器的水動力性能進行了數值預報，給出了在不同偏轉角度下的敞水曲線。在對轉槳方面，前后槳之間的推力和扭矩表現出強烈的非定常性和周期性特征。由于后槳工作在前槳的尾流場中，后槳周圍的流場表現出強烈的復雜性和非定常的特征，因此通常后槳的推力和扭矩難以準確預報。為提高數值模擬的精度，人們對影響數值計算精度的因素進行了大量研究。張濤［3］等采用MRF （moving reference frame）方法對對轉槳的敞水性能進行數值模擬，結果表明：非定常效應不能忽略；當研究對象是對轉槳時，采用準定常方法簡化是不合適的。王展智［4］等研究了時間步長的影響，發現時間步長的選擇取決于前后槳的葉數；當前后槳具有相同的葉數時，時間步長應該相對更小。

RANS模型是螺旋槳敞水性能預報中最常用的湍流模型。大量的研究結果表明，在預報單槳的推力和扭矩時，RANS模型已具有足夠的精度。盡管如此，由于耗散的存在，無法得到槳葉后面流場的準確信息。相較于RANS模型，LES和混合LES/URANS模型精度更高，但需花費更多的計算資源。目前，采用計算資源需求適中的混合LES/URANS模型來取得相對精確的流場信息是一個可行的選擇。

對轉舵槳的水動力性能是工程實踐中關注的焦點。此外，研究對轉舵槳各部件之間相互作用的內在機理對優化設計能夠起到幫助。本文采用IDDES結合滑移網格方法，對一組對轉槳以及對轉舵槳的水動力性能進行數值模擬，對兩者的水動力性能進行對比；進一步對兩者的推力系數進行頻譜分析，分析對轉舵槳振動特性以及各部件之間的相互作用。

1 IDDES方法簡介

IDDES模型是在RANS和LES的基礎上，經過一步步優化發展而來。湍流由各種大小和渦量的不同渦旋疊加而成，最大尺度的渦旋攜帶最多的雷諾應力，必須直接計算，而小尺度的渦旋對雷諾應力貢獻不大，并且接近各向同性，可以采用模型簡化［5］。理論上，可以采用DNS對全尺度范圍的湍流進行直接求解，但受計算資源的限制，目前工程中還不可能采用DNS進行計算。而RANS在處理湍流動能上耗散很大，在預測非定常湍流方面存在不足。LES就誕生在這樣的背景之下，其主要思想是通過空間濾波，濾掉小尺度的渦旋結構，并對其進行模型簡化，對大尺度的渦旋結構進行直接模擬。因此，在解析尺度方面，LES介于RANS和DNS之間。然而，對于壁面附近的流動，黏性作用占主導，不存在慣性子區，因此LES在壁面附近要求的網格量理論上與DNS一致。如果不使用壁面函數，為使邊界層內流場有效求解，物體周圍無因次的第一層網格厚度y+通常要在1以下，這嚴重制約了其在工程上的應用。

為解決RANS在預測非定常湍流方面的不足，同時解決LES在近壁面附近網格數量要求非常高的問題，將兩者結合起來是一個較好的解決方案。假定流體不可壓縮，經過空間濾波后的質量和動量方程如下：

式中：ρ代表密度，kg/m3；上劃線代表空間濾波。而采用時間平均之后N-S方程與濾波之后的N-S方程具有相同的數學形式，這是LES和RANS能夠得以結合起來的數學基礎。

SPALART等學者利用LES和RANS的控制方程在數學形式上相同的特點，發揮兩者的優點，同時回避兩者的缺點，巧妙地將LES和RANS結合起來，并于1997年提出了DES（detached eddy simulation）模型，通常稱為DES97。DES97的主要思想是在近壁面和亞格子區域采用Spalart-Allmaras（S-A）模型進行計算，在遠場采用LES進行計算，LES/RANS之間的切換通過引入一個長度尺度來實現。DES97對長度尺度 定義如下：

式中：d為離壁面最近的距離，分別是三個方向上的局部網格尺寸。CDES為常數，通常取0.65。

但在實踐中，人們發現RANS與LES之間的切換，完全取決于網格尺度大小，在某些情況下，RANS會提前切換到LES，出現模化應力耗散（modeled stress depletion， MSD）的現象，從而導致網格誘導分離（grid-induced separation， GID）的問題。此外，DES97還存在速度型偏離對數率（log layer mismatch， LLM）的問題。針對GID，SPALART在 2006年又進一步提出 DDES［6］（delayed detached eddy simulation）模型。DDES對長度尺度進行修改：

式中：fd為延遲過渡函數。

fd能避免LES在邊界層內求解，從而解決GID問題。

針對LLM問題，TARVIN［7］在2006年又提出IDDES（improved delayed detached eddy simulation）模型，并在2008年，由SHUR［8］進一步發展出一種融合DDES和WMLES（wall modeled LES）的方法。

IDDES不僅可以與S-A湍流模型結合，也可以與其他的兩方程模型（如SSTk-ω模型）結合。本文采用的是基于S-A湍流模型的IDDES模型。

2 幾何模型和邊界條件

2.1 幾何模型

數值模擬的對象是MILLER［9］設計的一組4X4 的對轉槳，記為 CRP6。CRP6 由 DTMB3686和 DTMB3687A 組成。前后槳之間的軸間距為0.043 2 m。另外，本文設計了一組吊艙和立柱，與CRP6共同組成了對轉舵槳（CRRP）系統。立柱是剖面為NACA66mod并且無拱度的機翼。CRP6和CRRP的外形見圖1。

圖1 計算模型

CRP6的主要幾何參數如表1所示。

表1 CRP6主要幾何參數

2.2 計算網格和邊界條件

網格利用ANSYS ICEM 16.0 生成。如圖2所示，整個計算域分為三部分：繞前后槳周圍的兩個旋轉區域，包含立柱和吊艙以及遠場的靜止區域。靜止區域是一個正方體，其長、寬、高分別為10DF（前槳直徑）、6DF，、6DF。吊艙中心距上游入口為4DF，距下游出口為 6DF。CRP6的計算域與CRRP的計算域基本一致。時考慮到吊艙、立柱和兩個旋轉區域。這兩個因素使分塊的布置較困難。為獲得網格質量高的六面體網格，分塊必須盡可能布置為正六面體。常規的Y形切分由于會在近旋轉域附近產生負網格，因此不能滿足網格質量要求。一種可行的方法是對吊艙進行兩次O形切分，對立柱進行兩次Y形切分，分成6個六面體。圖3提供了一種合適的分塊布置，CRP6和CRRP的計算網格如圖4 所示。

圖2 計算域

計算域的各部分均采用結構化六面體網格。由于前后槳所在的旋轉域在生成靜止區域網格的時候必須刪除，在布置吊艙和立柱周圍的分塊時必須同

圖3 立柱和吊艙周圍的分塊布置

圖4 計算網格

在物面附近，y+取40，邊界層網格為10層。每片槳葉所在區域網格在100萬左右，立柱和吊艙所在的靜止區域網格數目約為200萬，總網格數量為1 000萬左右。

邊界條件為：入口速度設為來流速度，壓力梯度為0；出口速度設為0梯度，壓力設為0；槳葉、吊艙以及立柱表面設為不可穿透物面條件；除入口和出口外，遠場邊界設為對稱性邊界。

3 計算結果與分析

3.1 CRP6

首先進行CRP6的敞水性能計算，并與MILLER在1976年試驗報告提供的試驗數據［9］作比較，以驗證數值方法的可行性。計算時，前后槳的轉速均取為720 r/min，與試驗時相同。時間步長設為0.000 115 74 s，對應于每轉720個時間步。進速系數為0.5～1.1，其通過改變來流速度而改變。進速系數、推力和扭矩系數等參數定義如下：

進速系數：

前槳推力系數：

前槳扭矩系數：

后槳推力系數：

后槳扭矩系數：

推進器推力系數：

推進器扭矩系數：

推進器敞水效率：

上述式中：VA為進速，m/s；n為轉速，r/min；T為推力，N；Q為扭矩，N·m；下標F和A分別為前槳和后槳。為了方便作圖，KTM和KQM分別定義為總推力和總扭矩的一半。

圖5 CRP6 敞水特性曲線

圖5是CRP6敞水性能的計算結果與試驗結果的對比。前槳的計算結果與試驗結果吻合良好，推力和扭矩的誤差均在2%左右。后槳推力較試驗結果低6%左右，扭矩高4%左右。整個推進器的推力誤差在3%左右，扭矩誤差在2%左右，效率的誤差在5%左右，與試驗結果吻合較好，證實本文采用的數值方法基本可行。

3.2 CRRP

MILLER設計CRP6主要是為研究對轉槳前后槳之間的相互干擾，前后槳的槳葉數目相同時，更容易產生共振。為研究CRRP各部件之間的相互干擾，設計一套吊艙和立柱，加上CRP6，這四個部分，共同組成整個CRRP系統。

CRRP的前后槳轉速與CRP6相同，均為720 r/min。時間步長設為0.000 115 74 s，對應于每轉720個時間步。進速系數為0.5～1.1，通過改變來流速度來調節進速系數。推力系數、扭矩系數等參數定義基本與CRP6相同，但由于吊艙和立柱的存在，其中一些系數存在差別，其定義如下：

吊艙單元軸向力：

吊艙單元橫向力：

推進器推力系數：

式中：下標pod代表吊艙和立柱。

圖6是CRRP的敞水性能曲線。由圖6（a）可以看出：隨著進速系數增大，前槳和后槳的推力系數和扭矩系數都減小，并且后槳的推力和扭矩系數減小得更快。吊艙和立柱單元的橫向力和軸向力系數如圖6（b）所示。在0°舵角工況下，吊艙和立柱單元的橫向力是一個接近于0的小量。隨著進速系數增加，吊艙和立柱所受的軸向力基本不變。由于和CRRP推力的方向相反，這部分力造成了CRRP系統的效率損失。

圖6 CRRP的敞水特性曲線

吊艙和立柱表面壓力分布如下頁圖7所示。由圖7（a）可以看出，立柱前端同時存在一個高壓中心和低壓中心，這個現象是由于前槳尾流場的周向旋轉運動造成的。由于后槳的抽吸作用，立柱的后端產生一個低壓中心，如圖7（b）所示。立柱和吊艙的壓力差同時作用，即產生壓差阻力。壓差阻力是吊艙和立柱軸向力的主要部分，粘性阻力只占其中很小一部分。

圖7 J= 0.9時吊艙和立柱表面的壓力分布

渦旋結構采用Q準則［10］定義，Q=500等值面如圖8所示。從圖中圓圈標出的部分可以看出，在立柱的前端存在明顯的回流，這就解釋了為什么在立柱的前端同時存在高壓中心和低壓中心。這種流動結構不僅會造成CRRP系統的推進效率損失，嚴重情況下，在低壓中心附近還可能產生空泡，在工程實踐中，應該盡量避免這種情況發生。隨著渦結構靠近后槳盤面，前槳的梢渦逐漸收縮，在經過后槳之后，前槳的梢渦被打碎，形成網狀渦。后槳的梢渦受前槳尾流的影響，存在一定程度的變形，但整體形狀并未被破壞。

圖8 Q = 500等值面

3.3 CRP6 和CRRP的敞水性能比較

CRP6和CRRP 的敞水性對比如表2所示。CRRP前槳的推力和扭矩系數明顯高于CRP6，后槳的推力系數兩者基本一致，而CRRP的扭矩系數略大于CRP6。整個推進器方面，盡管CRRP前槳推力系數大于后槳，但由于吊艙和立柱阻力的存在，CRRP整個推進器的單元推力小于CRP，扭矩系數則大于CRP。CRRP推進器推進效率比CRP小6%左右。

表2 CRP6和CRRP的敞水性能比較

3.4 CRP6 和 CRRP 各部件之間的相互干擾

圖9 CRP6和CRRP推力扭矩系數變化曲線

吊艙和立柱的存在改變了前后槳之間的相互作用。圖9是J= 0.9時CRP6和CRRP的推力和扭矩系數在旋轉一周內的時歷曲線，這些系數取自CFD計算結果。由圖9（a）可以看出，CRP6前槳的推力扭矩系數有8個周期，而CRRP前槳的推力扭矩系數僅有4個周期。CRP6推力和扭矩的變化幅值都比CRRP大。在圖9（b）中，CRP6和CRRP后槳的推力扭矩均具有8個周期，CRP6的系數隨時間變化更快，但變化幅值比CRRP小。圖9（c）是吊艙和立柱所受的軸向力和橫向力系數的時歷曲線，這些系數在槳旋轉一周內僅有4個周期。為進一步分析前后槳之間的相互作用，對CRP6和CRRP的前槳和后槳的推力系數進行了快速傅里葉變換（FFT）。根據STRASBERG［11］的研究結果，CRP推力和扭矩系數的波動頻率由式（16）決定：

式中：下標F和A分別代表前槳和后槳；f0代表軸頻；Z表示槳葉數目；m為整數，m和Z必須滿足方程；mF ZF=mA ZA。

下頁圖10是CRP6和CRRP推力系數的頻譜特性曲線。在圖10（a）中，CRP6前槳推力系數峰值的頻率為8倍軸頻的整數倍，滿足公式（16）。而CRRP前槳推力系數峰值的頻率為4倍軸頻的整數倍。從峰值上來看，CRRP前槳的峰值遠小于CRP6，并且8倍軸頻處峰值基本可以忽略。這意味著由于吊艙和立柱的存在，后槳對前槳的作用被弱化了，并且后槳周向誘導速度分量的影響可以忽略，軸向誘導速度分量是否有影響還需要進一步研究。由圖10（b）可以看出，CRP6后槳推力系數幅值在6～7階時仍較大，其基頻為8倍軸頻。而CRRP后槳推力系數幅值在5階以后變得很小，其基頻為4倍軸頻。2階幅值最大，頻率為8倍（即前后槳葉數之和）軸頻，說明前槳尾流的周向誘導速度分量對后槳有影響，兩者之間產生共振。除前槳尾流導致的共振峰之外，后槳本身的運動也是其推力脈動幅值的主要部分。CRP6前槳和后槳推力系數的峰值分別為0.028 89和0.005 16，而CRRP前槳和后槳推力系數的峰值分別為0.006 86和0.013 44，吊艙和立柱的軸向力峰值為0.004 73。

圖10 CRP6 和CRRP推力系數的頻譜分析

4 結 論

因此，盡管CRRP后槳的推力系數脈動幅值比CRP6大，但就推進器整體而言，CRRP推力產生的振動要比CRP6小得多。

本文設計了一套吊艙和立柱，與Miller設計的對轉槳結合在一起，組成CRRP系統。運用滑移網格方法，結合IDDES模型，對CRP6和CRRP的水動力性能進行了數值模擬。通過比較CRP6的敞水性能計算值與試驗值，驗證數值方法的可行性。另外，為研究推進器各部件之間的相互干擾，對CRP6和CRRP的推力系數進行FFT變換。通過計算和分析，得到以下初步結論：

（1）比較CRP6的計算值與試驗值，整個推進器單元的推力誤差約3%、扭矩誤差約2%、效率的誤差約5%，與試驗結果吻合較好，證明了本文采用的數值方法的可靠性。滑移網格方法結合IDDES模型，能夠用來對CRP這種相互作用十分明顯的復雜推進器的水動力性能進行數值預報。

（2）對CRRP進行水動力性能數值預報，并與CRP6進行對比。結果發現，CRRP前槳的推力和扭矩比CRP6前槳的推力和扭矩大，而兩者后槳的推力扭矩則基本一致。吊艙和立柱阻力的存在使整個CRRP推進器單元的推力比CRP6小，而CRRP的扭矩則比CRP6大，導致CRRP的推進效率比CRP6小6%左右。

（3）吊艙和立柱阻力的主要成分是由其前后端的壓力差產生壓差阻力。在立柱前端存在明顯回流，從而使立柱的前端同時存在一個高壓中心和低壓中心。這不僅會減小CRRP系統的推進效率，還有可能在低壓中心附近產生空泡，工程實踐中應該避免出現這種情況。

（4）吊艙和立柱的存在改變了CRRP前后槳的頻譜特性曲線。前槳推力脈動幅值的頻率不再滿足公式（16），而是4倍軸頻的整數倍。CRRP前槳推力脈動幅值比CRP6前槳明顯減小。后槳周向誘導速度分量對前槳的影響基本可以忽略，但軸向誘導速度是否有影響則還需作進一步研究。就后槳而言，前后槳之間相互作用產生的共振依然存在，其自身運動產生的振動也是后槳脈動幅值的主要成分。但從整個推進器系統來看，CRRP推力的脈動幅值明顯小于CRP6。相對于CRP，采用CRRP能夠減小由軸承力導致的振動。

本文工作得到了本院吳瓊和于海兩位老師的指導，在此一并表示感謝。