葉劍平，莊光宇

(1.海軍駐719 所軍事代表室，湖北 武漢430064;2.中國人民解放軍92857 部隊，北京102444)

1 概 述

計算流體力學［1］(computational fluid dynamics，CFD)通過計算機數(shù)值計算和圖像顯示，對包含有流體流動和熱傳導等相關物理現(xiàn)象的系統(tǒng)進行分析。其基本思想可以歸結為:把原來在時間域及空間域上連續(xù)的物理量的場，如速度場和壓力場。用一系列有限個離散點上的變量值的集合來代替，通過一定的原則和方式建立起關于這些離散點上場變量之間關系的代數(shù)方程組，然后求解代數(shù)方程組獲得場變量的近似值。CFD 方法對流體流動進行數(shù)值模擬的過程，通常包括以下幾方面:

1)建立問題中反映各個量之間關系的基本微分控制方程(通常包括質(zhì)量守恒方程、動量守恒方程、能量守恒方程等)及相應的定解條件。

2)尋求離散化控制方程的數(shù)值方法(如有限差分法、有限元法、有限體積法等)，貼著螺旋槳表面建立坐標系，邊界條件的處理等。

3)劃分網(wǎng)格、輸入初始條件和邊界條件、設定控制參數(shù)等。

4)通過圖表等方式顯示計算結果，對質(zhì)量和結果進行檢驗檢查和判斷分析，得出結論。

粘性流體CFD 計算的控制方程是N-S 方程，其主要的湍流數(shù)值模擬可分為直接模擬數(shù)值方法(DNS)、大渦模擬方法 (LES)和雷諾平均方法(RANS)。DNS 方法直接用瞬態(tài)的N-S 方程對湍流進行計算，無需對湍流流動作任何簡化和近似，理論上可以得到比較準確的計算結果;LES 方法采用空間過濾后的瞬態(tài)N-S 方程來模擬湍流中的大尺度渦，而被過濾掉的、比網(wǎng)格尺度小的渦，是通過近似的模型來考慮其對大渦的影響;RANS 方法的核心是不直接求解瞬時的N-S 方程，而是求解時均化的N-S 方程，這樣的近似，既可以避免DNS 方法計算量大的問題，其精度也可滿足工程實際應用的要求。由于螺旋槳周圍流場的復雜性以及目前計算機的性能所限，DNS 方法和LES 方法還很少應用于螺旋槳周圍湍流場的模擬計算，更多的是通過RANS方法來進行數(shù)值模擬。目前，在研究螺旋槳方面，RANS 方法是使用最廣泛、技術最成熟的湍流數(shù)值模擬方法。

1988年，Stern 等［2］應用升力面方法得到螺旋槳載荷分布，代入粘性流體CFD 方法的流場控制方程作為源項，進行了船槳干擾下船舶流場的計算工作。1997年，唐登海等［3］在求解船舶螺旋槳周圍粘性流場的數(shù)值解時，應用直接求解RANS 方程的方法。分析了由于粘性流動引起的螺旋槳周圍及螺旋槳尾流場的一些重要特征。1998年，第22 屆ITTC會議曾組織了一個專題［4］對面元法和RANS 方法進行實例計算比較。通過比較發(fā)現(xiàn)，RANS 方法和面元法都可較準確地預報螺旋槳的敞水性能以及槳葉表面的壓力分布，說明RANS 方法已經(jīng)比較成熟，可成為螺旋槳設計和性能預報的有力工具。2003年，Rhee 等［5］采用混合空泡模型模擬計算了螺旋槳定?？张萘鲃?，計算得到的螺旋槳敞水性能結果同試驗測量結果非常接近，模擬出的空泡初生以及空泡形狀同試驗觀測結果比較相似。2007年，Mikkelsen 等［6］用RANS 方程研究船后螺旋槳性能時，在槳前方用1 個鼓動盤來產(chǎn)生所需的非均勻船體伴流，計算結果顯示這種方法可以用來研究工作于船后尾流場中的螺旋槳。2009年，Vladimir Krasilnikovl 等［7］采用RANS 的方法研究了幾何參數(shù)對螺旋槳誘導流場的影響。計算了不同盤面比、側(cè)斜和葉數(shù)情況下螺旋槳的誘導流場，在此基礎上，進一步分析了這些參數(shù)對螺旋槳尾流場的影響。從計算結果可看出，槳葉盤面比和側(cè)斜的增大都會使螺旋槳的誘導流場變小;葉數(shù)的變化主要影響著盤面內(nèi)螺旋槳誘導流場的周期數(shù)。2010年，王超等［8］結合RNG k-ε 湍流模型，運用滑移網(wǎng)格技術對粘性流場中槳舵相互干擾引起的三維非定常湍流進行了計算，并與試驗測量值作了比較。從對計算結果的分析可知，利用滑移網(wǎng)格技術及RNG k-ε 湍流模型，可以較好地模擬槳舵干擾的水動力性能問題。目前，RANS 方法在組合推進器、空泡性能預報乃至螺旋槳設計等問題上的應用也已相當廣泛。

2 計算模型的建立

以DTRC4119 槳為研究對象，建立模型及進行網(wǎng)格劃分、邊界條件設定，給出湍流模型的選擇和求解控制參數(shù)的設置，梳理螺旋槳的水動力性能預報、螺旋槳的尾流場模擬與分析的整個過程。DTRC4119 槳的主參數(shù)［9］見表1。

表1 DTRC4119 槳的主參數(shù)Tab.1 The primary parameter of DTRC4119 propeller

2.1 計算網(wǎng)格劃分與數(shù)值模擬參數(shù)設定

采用CFD 數(shù)值模擬螺旋槳的敞水試驗。鑒于來流均勻和螺旋槳幾何上的周期性，為減少計算時間，提高計算效率，可只取單個槳葉所在的單通道(周向展開120°)作為計算域。計算模型見圖1。在建模過程中使用的是直角坐標系O-XYZ。計算域的內(nèi)邊界取在槳轂和葉片表面上，其中槳轂模擬敞水試驗;外邊界面取在直徑約為螺旋槳直徑5 倍的圓柱體表面上。為計算的需要，又把整個域分成1 個小域和1 個大域，這樣便于在劃分網(wǎng)格時進行局部加密，提高計算結果的準確度。網(wǎng)格劃分是CFD 模擬過程中的重要環(huán)節(jié)?？拷鼧~和槳轂表面的這層非結構網(wǎng)格尺度為直徑的1/200 左右，以保證網(wǎng)格y+值的取值范圍(30 ＜y+＜300)。同時，在劃分網(wǎng)格時使用了局部加密的方法，對于槳葉與槳轂的連接處以及葉梢部分等進行局部加密，以便捕捉到重要的流場信息。螺旋槳表面及周圍流場網(wǎng)格劃分如圖2所示。

圖1 單通道數(shù)值模擬模型Fig.1 Single chanel numerical simulate model

圖2 網(wǎng)格劃分示意圖Fig.2 Sketch map of reseau partition

計算區(qū)域包括的邊界條件有:入口邊界、出口邊界、槳葉和槳轂表面、扇形柱狀側(cè)面。本文選用Fluent 軟件提供的運動參考坐標系模型(即MRF 模型)，進口邊界和扇形柱體圓弧面設定為速度進口條件，出口邊界設定為壓力出口，槳葉和槳轂表面為無滑移物面邊界條件，扇形柱狀大域、小域兩側(cè)面設定為旋轉(zhuǎn)周期性邊界條件，轉(zhuǎn)動區(qū)域(小域)與靜止區(qū)域(大域)的交界面設為交換面。計算域內(nèi)的流體則按MRF 模型設置為繞x 軸以角速度N/rpm 旋轉(zhuǎn)。

2.2 湍流模型的選擇與求解控制參數(shù)的設置

計算時，湍流模型的選擇與求解控制參數(shù)的設置見表2。一般迭代1 000 步左右即可收斂。

表2 各類數(shù)值格式和選項的設置Tab.2 The numeric format and the option settings

3 計算結果與分析

3.1 敞水性能計算

對于艦船螺旋槳需要計算:槳推力T，槳扭矩Q，螺旋槳轉(zhuǎn)速n，水流進速VA。為了便于比較分析，通常均以槳直徑無因次化。

進速系數(shù)

槳推力系數(shù)

槳扭矩系數(shù)

推進效率

3.2 敞水性能計算結果及與試驗結果對比

進速系數(shù)J 分別取為0.5，0.7，0.833，0.9 和1.0，螺旋槳轉(zhuǎn)速為一定值n=900 r/min。表3 分別列出了在不同進速系數(shù)下，得到的螺旋槳推力、扭矩和推進效率計算值與試驗值的比較結果。圖3 繪出敞水特征曲線結果比較。

表3 DTRC4119 敞水性能比較Tab.3 The hydrodynamic open water performance of DTRC4119 propeller

圖3 螺旋槳敞水性能比較Fig.3 The hydrodynamic open water performance of propeller

對比表3 和圖3 中螺旋槳水動力性能計算值與試驗值可知，計算所得的Kt 曲線與10Kq 曲線都比試驗值略高，這應該是由于空泡的存在所導致。在試驗中，吸力面上低于飽和壓力的區(qū)域出現(xiàn)了空泡，空泡區(qū)的壓力保持為飽和蒸汽壓力;而數(shù)值模擬時，由于未使用特定的空化模型，所以導致模擬得到的Kt 與10Kq 值存在一定的誤差?？偟膩碚f，在設計工況J=0.833 處吻合得最理想，在考察的進速系數(shù)范圍內(nèi)(0.5 ～1.1)，推力系數(shù)Kt、扭矩系數(shù)Kq 以及敞水效率η 的計算結果與試驗結果吻合較好，推力系數(shù)Kt 誤差在2.4%以內(nèi)，扭矩系數(shù)Kq 誤差在1%以內(nèi)。在高進速系數(shù)下，由于考察值的絕對量已經(jīng)很小，所以計算值微小的變化都會造成誤差相應增大。說明在有關螺旋槳敞水性能的計算中，采用CFD 軟件進行的數(shù)值模擬結果可以滿足工程應用的要求。

3.3 螺旋槳的剖面壓力分布比較

以DTRC4119 槳為例，將壓力分布的CFD 計算結果與試驗測量結果［10］以及Hoshino［11］采用面元法的計算結果進行比較，進速系數(shù)J 取0.833，比較結果見圖4 ～圖6。在0.3R 剖面處，CFD 計算結果和實驗結果相比較整體偏高，但是二者的形狀基本相似，整體趨勢基本相當，與面元法計算結果的精度大致相同;在0.7R 和0.9R 剖面處，除了在導邊和隨邊的附近CFD 計算結果和實驗結果稍有差別外，其他地方的壓力分布吻合很好，雖然導邊附近和隨邊附近的壓力分布計算結果和試驗結果有差別，但整個葉切面壓力分布的整體趨勢基本相當。說明CFD 的計算方法能較好地計算螺旋槳定常狀態(tài)下的壓力分布。

圖4 P4119 r/R=0.3 處剖面壓力分布比較Fig.4 P4119 r/R=0.3 section pressure distributions

圖5 P4119 r/R=0.7 處剖面壓力分布比較Fig.5 P4119 r/R=0.7 section pressure distributions

圖6 P4119 r/R=0.9 處剖面壓力分布比較Fig.6 P4119 r/R=0.9 section pressure distributions

3.4 螺旋槳尾流場的數(shù)值模擬

圖7 ～圖9 為CFD 計算的螺旋槳DTRC4119 在x/R=0.295，r/R=0.7，J=0.806 處的速度分布與實驗結果以及面元法的比較［12］?？梢钥闯?，在徑向速度分布上，CFD 計算值與實驗值有良好的一致性，比面元法準確度略高一些;在切向速度分布上，CFD 計算值仍然與實驗值保持良好的一致性，特別是在速度峰值處的模擬值要優(yōu)于面元法;在軸向速度分布上，CFD 計算值雖然與實驗值有相似的變化規(guī)律，但量值上存在一定的差距，基本滿足工程上的應用。本文的CFD 方法雖然數(shù)值模擬了粘性的影響，但由于采用的是時均化的RANS 方法，此方法是人為采取的平均化方法，并且引入的湍流方程也是一種近似化的假設。另外，流場中各參量在計算中由于網(wǎng)格間的數(shù)據(jù)傳遞而存在耗散現(xiàn)象，因此還無法完全體現(xiàn)出螺旋槳尾流場中尖銳的誘導速度峰值。

圖7 徑向速度分布Fig.7 Radial velocity distributions

圖8 切向速度分布Fig.8 Favorite velocity distributions

圖9 軸向速度分布Fig.8 Axle velocity distributions

圖10 和圖11 分別為x/R=0.295 處截面的軸向速度和橫剖面內(nèi)的速度(徑向速度與切向速度的合速度)的分布情況。從中可以清楚地看出尾渦面的位置及周圍流場的速度分布，圖11 清晰地顯示出螺旋槳稍渦的存在。應用CFD 數(shù)值模擬，可以直觀模擬螺旋槳尾流場，便于對螺旋槳的尾流場進行分析。

圖10 槳后軸向誘導速度Fig.10 Induced axle velocity behind propeller

圖11 槳后橫向誘導速度分布Fig.11 Induced transverse velocity behind propeller

4 結 語

本文在介紹CFD 概念的基礎上，利用CFD 計算方法預報了螺旋槳水動力性能并數(shù)值模擬了螺旋槳槳后粘性尾流場。通過對獲得的數(shù)值結果和試驗結果進行比較，得出以下結論:

1)網(wǎng)格質(zhì)量對CFD 計算精度和計算效率有重要的影響。網(wǎng)格主要分為結構網(wǎng)格和非結構網(wǎng)格。結構網(wǎng)格使用規(guī)則單元便于計算，但適用性較差;非結構網(wǎng)格具有良好的靈活性，尤其對復雜問題比較有效，但在計算上相對耗時。因此，一般使用在較簡單的幾何模型上生成單塊結構化網(wǎng)格和多塊結構化網(wǎng)格，在復雜幾何模型上生成非結構化網(wǎng)格，在結構化網(wǎng)格和非結構化網(wǎng)格之間使用楔形網(wǎng)格過渡的混合網(wǎng)格。并且利用局部加密技術，獲取特定位置的精確信息。

2)湍流模型的不同對螺旋槳水動力性能的數(shù)值預報存在明顯差異。目前幾種常用的湍流模型，如標準k-ε 湍流模型、RNG k-ε 湍流模型和SST k-ω模型等都無法完全模擬不同狀況下的真實流場，只是在不同應用中有各自的優(yōu)點。因此，選擇合適的湍流模型才能夠?qū)β菪龢恼承粤鲌鲞M行更準確的數(shù)值模擬。

3)近壁面區(qū)域處理方式的選取影響數(shù)值模擬的結果。對于近壁面區(qū)域的處理，F(xiàn)luent 選用的默認方法實際上是一組半經(jīng)驗公式，用于將壁面上的物理量與湍流核心區(qū)域內(nèi)待求的未知量直接聯(lián)系起來。這里應該結合網(wǎng)格的壁面精細劃分，保證y+值在有效范圍內(nèi)，使流動處于對數(shù)律層。

4)發(fā)現(xiàn)CFD 對螺旋槳尾流場以及梢渦的形成與結構的數(shù)值模擬有較高的準確性，不但能從定性上而且能從定量上預報螺旋槳粘流場。CFD 計算結果與試驗結果相比，總的來看吻合較好，采用CFD軟件進行的數(shù)值模擬結果可以滿足工程應用的要求。但也應看到，在細節(jié)部分尚存在一定的偏差，在某些位置仍有較大的誤差，這都需要在以后的工作中進一步改進與完善。

［1］王福軍.計算流體動力學分析——CFD 軟件原理與應用［M］.北京:清華大學出版社，2004.

［2］STERN F，KIM H，PATEL V C，et al.A viscous-flow approach to the computation of propeller-hull interaction［J］.Journal of Ship Research，1988，32(4):246-184.

［3］唐登海，董世湯.船舶螺旋槳周圍粘性流場數(shù)值預報與流場分析［J］.水動力學研究與進展(A)，1997，12(4):426-436.

［4］Proceedings of 22ndITTC Propulsion Committee Propeller RANS/PANEL Method Workshop.Grenoble，F(xiàn)rance，April，1998，5-6.

［5］RHEE S，JOSHI S.CFD validation for a marine propeller using an unstructured mesh based RANS method［C］.Proceedings of FEDSM’03，Honolili，USA，2003.

［6］MIKKELSEN R，ANDERSEN P，SORENSEN J N.Modeling of behind condition wake flow in RANS computation on a conventional and high skew propeller［C］.Proceedings of 10thNumerical Towing Tank Symposium，Hamburg，Germany，September，2007.

［7］KRASILNIKOV V，SUN J，HALSE K H.CFD investigation in scale effect on propellers with different magnitude of skew in turbulent flow［C］.First International Symposium on Marine Propulsors SMP′ 09，Trondheim，Norway，June，2009.

［8］王超，黃勝，常欣，等.基于滑移網(wǎng)格與RNG k-ε 湍流模型的槳舵干擾性能研究［J］.船舶力學，2011，15(7):715-721.

WANG Chao，HUANG Sheng，CHANG Xin，et al.Research on the hydrodynamics performance of propeller-rudder interaction based on sliding mesh and RNG k-ε model［J］.Journal of Ship Mechanics，2011，15(7):715-721.

［9］蘇玉民，黃勝.船舶螺旋槳理論［M］.哈爾濱:哈爾濱工程大學出版社，2003.163-169.

［10］JESSUP S D.An experimental investigation of viscous aspects of propeller blade flow［D］.Washington:School of Engineering and Architecture，The Catholic University of America，1989.

［11］HOSHINO.Data for the 22ndITTC 98 workshop on panel method and RANS solution［C］.Proceedings of the 22ndITTC，Grenoble，F(xiàn)rance，1998.