葉劍平，莊光宇

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

1 概 述

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

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

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

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

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

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

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

2 計(jì)算模型的建立

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

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

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

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

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

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

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

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

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

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

3 計(jì)算結(jié)果與分析

3.1 敞水性能計(jì)算

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

進(jìn)速系數(shù)

槳推力系數(shù)

槳扭矩系數(shù)

推進(jìn)效率

3.2 敞水性能計(jì)算結(jié)果及與試驗(yàn)結(jié)果對比

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

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

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

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

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

以DTRC4119 槳為例，將壓力分布的CFD 計(jì)算結(jié)果與試驗(yàn)測量結(jié)果［10］以及Hoshino［11］采用面元法的計(jì)算結(jié)果進(jìn)行比較，進(jìn)速系數(shù)J 取0.833，比較結(jié)果見圖4 ～圖6。在0.3R 剖面處，CFD 計(jì)算結(jié)果和實(shí)驗(yàn)結(jié)果相比較整體偏高，但是二者的形狀基本相似，整體趨勢基本相當(dāng)，與面元法計(jì)算結(jié)果的精度大致相同;在0.7R 和0.9R 剖面處，除了在導(dǎo)邊和隨邊的附近CFD 計(jì)算結(jié)果和實(shí)驗(yàn)結(jié)果稍有差別外，其他地方的壓力分布吻合很好，雖然導(dǎo)邊附近和隨邊附近的壓力分布計(jì)算結(jié)果和試驗(yàn)結(jié)果有差別，但整個葉切面壓力分布的整體趨勢基本相當(dāng)。說明CFD 的計(jì)算方法能較好地計(jì)算螺旋槳定常狀態(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 計(jì)算的螺旋槳DTRC4119 在x/R=0.295，r/R=0.7，J=0.806 處的速度分布與實(shí)驗(yàn)結(jié)果以及面元法的比較［12］?？梢钥闯觯趶较蛩俣确植忌?，CFD 計(jì)算值與實(shí)驗(yàn)值有良好的一致性，比面元法準(zhǔn)確度略高一些;在切向速度分布上，CFD 計(jì)算值仍然與實(shí)驗(yàn)值保持良好的一致性，特別是在速度峰值處的模擬值要優(yōu)于面元法;在軸向速度分布上，CFD 計(jì)算值雖然與實(shí)驗(yàn)值有相似的變化規(guī)律，但量值上存在一定的差距，基本滿足工程上的應(yīng)用。本文的CFD 方法雖然數(shù)值模擬了粘性的影響，但由于采用的是時均化的RANS 方法，此方法是人為采取的平均化方法，并且引入的湍流方程也是一種近似化的假設(shè)。另外，流場中各參量在計(jì)算中由于網(wǎng)格間的數(shù)據(jù)傳遞而存在耗散現(xiàn)象，因此還無法完全體現(xiàn)出螺旋槳尾流場中尖銳的誘導(dǎo)速度峰值。

圖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 清晰地顯示出螺旋槳稍渦的存在。應(yīng)用CFD 數(shù)值模擬，可以直觀模擬螺旋槳尾流場，便于對螺旋槳的尾流場進(jìn)行分析。

圖10 槳后軸向誘導(dǎo)速度Fig.10 Induced axle velocity behind propeller

圖11 槳后橫向誘導(dǎo)速度分布Fig.11 Induced transverse velocity behind propeller

4 結(jié) 語

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

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

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

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

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

［1］王福軍.計(jì)算流體動力學(xué)分析——CFD 軟件原理與應(yīng)用［M］.北京:清華大學(xué)出版社，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ù)值預(yù)報(bào)與流場分析［J］.水動力學(xué)研究與進(jìn)展(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］.船舶力學(xué)，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］.哈爾濱:哈爾濱工程大學(xué)出版社，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.