繆宇躍 孫江龍

華中科技大學(xué) 船舶與海洋工程學(xué)院，湖北 武漢 430074

1 引 言

由于數(shù)值模擬相對于實驗研究具有獨特的優(yōu)點，如成本低、周期短，能獲得完整的數(shù)據(jù)，能模擬出實際運行過程中各種測量數(shù)據(jù)的狀態(tài)，因而目前計算流體動力學(xué)（CFD）技術(shù)被廣泛應(yīng)用于工程領(lǐng)域。

在運用CFD動網(wǎng)格方法研究推進(jìn)方面，Liefvendahl和 Troeng［1］以潛艇螺旋槳為例運用并行計算模擬了大范圍邊界移動變形問題；文獻(xiàn)［2］對比了定?；旌厦鏍顟B(tài)與非定常動網(wǎng)格狀態(tài)下螺旋槳的推力和扭矩，并研究了非均勻伴流的影響；喻欣［3］運用滑移網(wǎng)格和動網(wǎng)格技術(shù)研究了螺旋槳在升沉運動中的水動力性能；李輝［4］對全方位推進(jìn)器運用動網(wǎng)格模擬進(jìn)行了分析并與MRF技術(shù)模擬結(jié)果予以了對比。他們的研究方法效果良好，值得借鑒。

本文采用CAD及Gambit軟件，根據(jù)螺旋槳的投影原理及其型值參數(shù)，建立螺旋槳的三維模型，分區(qū)劃分混合網(wǎng)格，使用CFD軟件（Fluent），分別采用動網(wǎng)格技術(shù)和MRF技術(shù)，對敞水螺旋槳的水動力性能進(jìn)行了數(shù)值模擬，得出各種不同進(jìn)速系數(shù)下槳葉的敞水性能曲線，計算結(jié)果與試驗測量值吻合良好。

2 數(shù)學(xué)模型

2.1 控制方程

RANS方程是粘性流體運動學(xué)和動力學(xué)的普適性控制方程，本文以之為求解螺旋槳水動力性能計算的基本方程。其形式如下：

式中，ρ為流體密度；P為靜壓；fi為單位質(zhì)量的質(zhì)量力；ui，uj為速度分量。

湍流脈動動能方程（k方程）為：

2.2 建立模型

DTMB P4119槳是一種無側(cè)斜無后傾分布的三葉螺旋槳，被ITTC選為考證數(shù)值方法預(yù)報精度的標(biāo)準(zhǔn)，本文將采用該槳進(jìn)行研究，其尺寸如表1所示［5］。

根據(jù)螺旋槳的投影原理（圖 1）［6］，首先將坐標(biāo)系O1X1Y1Z1中的型值點經(jīng)過一次旋轉(zhuǎn)得到坐標(biāo)系O1UVW 中坐標(biāo)，然后由坐標(biāo)系 O′X′Y′Z′與O1UVW的關(guān)系以及圖1b中平面坐標(biāo)與圖1a中柱坐標(biāo)的關(guān)系，可推導(dǎo)出螺旋槳葉切面處局部坐標(biāo)系O1X1Y1Z1至全局坐標(biāo)系OXYZ的坐標(biāo)變換公式為：

表1 DTMB P4119螺旋槳的幾何參數(shù)Tab.1 Geometric features of DTMB P4119 propeller

式中，θ為縱傾角，φ為螺距角；全局坐標(biāo)系OXYZ的 OYZ 平面與槳轂端面平行，坐標(biāo)系 O′X′Y′Z′與OXYZ平行，坐標(biāo)系O1X1Y1Z1的O1X1軸經(jīng)過葉切面的最厚處，坐標(biāo)系O1UVW由O1X1Y1Z1經(jīng)過一次旋轉(zhuǎn)得到；L為O′O1的長度；r為葉切面所在圓柱面的半徑。

將以型值表形式表達(dá)的螺旋槳葉切面的局部坐標(biāo)轉(zhuǎn)換為全局坐標(biāo)后，在CAD及Gambit軟件中進(jìn)行三維建模，坐標(biāo)軸的定義為：X軸與螺旋槳旋轉(zhuǎn)軸一致，指向下游；Y軸與槳葉參考線一致；Z軸服從右手系，計算采用全尺度模型。

3 網(wǎng)格劃分及計算策略

網(wǎng)格劃分是CFD模擬過程中最為耗時的環(huán)節(jié)，也是直接影響模擬精度和效率的關(guān)鍵因素之一。網(wǎng)格過疏或過密都會極大地影響計算結(jié)果。網(wǎng)格若過疏，往往會得到不精確甚至完全錯誤的結(jié)果；網(wǎng)格若過密，又會使計算量增大，使計算難以收斂。本文在劃分網(wǎng)格時使用了局部加密的方法，對于出口段的網(wǎng)格，將其密度給予了適當(dāng)降低，以便于控制總網(wǎng)格數(shù)。這樣，在網(wǎng)格模型總節(jié)點數(shù)一定的情況下，不僅可以提高計算精度，還可以避免流場變化平緩區(qū)域的計算資源浪費［7］。

螺旋槳表面網(wǎng)格劃分如圖2所示，是采用適應(yīng)性好的三角網(wǎng)格對槳葉和槳轂進(jìn)行加密，尺寸大小控制在直徑的1%左右，以滿足計算要求和精度。整個計算流域應(yīng)取一個足夠大的圓柱體，直徑和長度應(yīng)取螺旋槳直徑的數(shù)倍，以模擬無限流場。借鑒計算經(jīng)驗，遠(yuǎn)流場對螺旋槳的影響十分微小，而動網(wǎng)格技術(shù)對模型、網(wǎng)格和步長等參數(shù)及計算機配置的要求則較苛刻，且極為耗時。為減少網(wǎng)格，節(jié)約時間，本文流域取3倍直徑的圓柱，螺旋槳置于其中，距入口1倍直徑，距出口4倍直徑。

同時，在緊鄰螺旋槳的區(qū)域建立一個圓柱包裹螺旋槳（圖3中深黑色小圓柱），為網(wǎng)格加密區(qū)，加密區(qū)采用非結(jié)構(gòu)網(wǎng)格來適應(yīng)復(fù)雜的螺旋槳形狀。該區(qū)域的網(wǎng)格尺寸應(yīng)較小，以接近槳葉面網(wǎng)格尺寸，并向外平緩增大。為此，加密小圓柱外表面網(wǎng)格尺寸，取為槳葉面網(wǎng)格尺寸的3倍。該區(qū)域之外的流域網(wǎng)格采用稀疏的結(jié)構(gòu)網(wǎng)格，網(wǎng)格較大，向外增長也較大，但要平穩(wěn)過渡。為此，大圓柱外表面的網(wǎng)格尺寸取為槳葉面網(wǎng)格尺寸的6倍。尾流段按一定比例，采取由面生成體及網(wǎng)格方法拖出。這樣，既可保證計算精度，又可減少網(wǎng)格總數(shù)。流域如圖3所示，網(wǎng)格總數(shù)約100×104。

Fluent軟件有3種動邊界控制實現(xiàn)方法：彈簧近似光滑模型、動態(tài)分層模型和局部重劃模型?？紤]到螺旋槳復(fù)雜的槳葉幾何形狀和大的網(wǎng)格變化范圍，本文采用彈簧近似光滑模型與局部重劃模型相結(jié)合的方法控制槳葉的旋轉(zhuǎn)［8］。

計算策略為：入口設(shè)為速度入口，出口設(shè)為水流出口，流域壁面設(shè)為速度入口以模擬無限流場。動量、湍流動能和湍流耗散率均采用二階迎風(fēng)格式進(jìn)行離散，其余保持默認(rèn)。固定轉(zhuǎn)速為20 r／s，通過改變來流速度來實現(xiàn)不同的進(jìn)速系數(shù)，注意需設(shè)置合適的步長。

4 計算結(jié)果及分析

4.1 流場分析

觀察對比圖4和圖5（縱坐標(biāo)為壓力值），可以看到在進(jìn)速系數(shù)相同的情況下，葉面的表面壓力要大于葉背表面壓力，這與旋轉(zhuǎn)中的螺旋槳產(chǎn)生往前的推力相符。實際上，螺旋槳產(chǎn)生的推力就是螺旋槳所受到的壓差力與阻力的代數(shù)和。從葉根到葉梢，壓力不斷增加，在0.7r左右達(dá)到最大，然后逐漸減小，在葉梢處達(dá)到最小值。而且，可以明顯地看出壓力由隨邊向?qū)н呏饾u增大，在導(dǎo)邊處達(dá)到最大。

再結(jié)合圖6和圖7（縱坐標(biāo)為壓力值），可以看到水流集中在螺旋槳葉間，在螺旋槳盤面區(qū)域形成高壓區(qū)和高速區(qū)。從中可看出，螺旋槳對上游槳盤外的水流有明顯的抽吸現(xiàn)象，使流體由槳盤外流向槳盤內(nèi)。在葉梢處，形成了低壓區(qū)，產(chǎn)生了梢渦，這是因為該處的水流較快，且葉面和葉背存在著較大的壓差。

對比圖8和圖9可看出，從槳轂開始，速度和壓力都在下降，雖然經(jīng)過一段流域后會有所回復(fù)，但仍然比周圍流域低，可以看到尾渦的影響；從槳葉往后的水流壓力和速度都高于外流域，這是螺旋槳向后推水，水反作用于螺旋槳，由此產(chǎn)生推力的結(jié)果。

同時，從這些圖中還可看到，較大的變化都高度集中于螺旋槳盤面區(qū)域及其尾流段，盤面區(qū)以外的變化相對較小，隨著遠(yuǎn)離葉梢，外流場的變化和影響迅速減小，可見計算模型是可靠的。

圖10 分別為進(jìn)速系數(shù) J＝0.5，0.7 和 0.9 時螺旋槳的尾流情況。

通過比較可以看出，尾渦外半徑處的螺距是隨J的增大而增大，其直徑小于螺旋槳直徑。文獻(xiàn)［9］采用PIV測量四葉側(cè)斜槳尾渦，也揭示了相同的情況。圖中也顯示出，內(nèi)半徑處尾流在螺旋槳的作用下向內(nèi)聚攏，其直徑逐漸減小，反映了螺旋槳對流體的抽吸。

4.2 敞水計算結(jié)果與試驗的對比

本文同時使用MRF技術(shù)進(jìn)行了螺旋槳敞水模擬，以求對比于動網(wǎng)格技術(shù)的準(zhǔn)確性，如表2～4所示。設(shè)置計算域內(nèi)的流體按MRF模型，繞軸以角速度20 r／s旋轉(zhuǎn)，其他設(shè)置默認(rèn)。通過Fluent模擬計算，得出不同進(jìn)速系數(shù)J情況下的螺旋槳推力與扭矩，進(jìn)而求出槳的推力系數(shù)Kt、轉(zhuǎn)矩系數(shù)Kq以及敞水效率η。

表2 動網(wǎng)格與實驗的對比Tab.2 Comparisons of the dynamic mesh and experiments

表3 MRF原模型與實驗的對比Tab.3 Comparisons of the original model and experiments

表4 MRF新模型與實驗的對比Tab.4 Comparisons of the new model and experiments

由于基于原模型使用MRF技術(shù)計算的結(jié)果相對于試驗值［10］其誤差并沒有降低，沒有表現(xiàn)出相對動網(wǎng)格技術(shù)的精確優(yōu)勢，因此本文對原模型進(jìn)行了加密處理，主要是增加了包裹螺旋槳的小區(qū)域的網(wǎng)格數(shù)量，最終網(wǎng)格總數(shù)約為180×104，然后，對新模型進(jìn)行MRF模擬，得到結(jié)果并進(jìn)行對比，由表2與表4的數(shù)據(jù)與實驗數(shù)據(jù)得到圖11。

由表 2～4及圖 11可看出，除了在 J＝1.1時出現(xiàn)了明顯的誤差外，用動網(wǎng)格計算的kt，10kq與η的平均誤差分別為 3.52%，5.63%，-2.19%；在 J不是很大的情況下，曲線吻合較好，驗證了該數(shù)值方法模擬的可靠性。

對模型加密前，采用MRF技術(shù)計算的結(jié)果其誤差均高于動網(wǎng)格技術(shù)的誤差；而對模型加密后，采用MRF技術(shù)計算的結(jié)果其平均誤差分別為0.6%，4.97%和-4.16%，除 η 外，其它在 5%以內(nèi)，低于動網(wǎng)格技術(shù)的平均誤差，體現(xiàn)了MRF技術(shù)的精確性。由此可見，網(wǎng)格數(shù)量對MRF技術(shù)模擬具有重要影響，可為網(wǎng)格劃分提供參考依據(jù)。

5 結(jié) 論

1）根據(jù)螺旋槳投影原理，利用CAD軟件建立三維模型，實現(xiàn)了數(shù)據(jù)到實體模型的轉(zhuǎn)換。

2）通過CFD動網(wǎng)格技術(shù)數(shù)值模擬，分析了螺旋槳敞水壓力速度分布及尾流情況，為螺旋槳的研究提供了依據(jù)。

3）探討了運用動網(wǎng)格技術(shù)實現(xiàn)三維復(fù)雜模型回轉(zhuǎn)運動的動網(wǎng)格計算。結(jié)合動網(wǎng)格和MRF技術(shù)，將得到的敞水結(jié)果與實驗結(jié)果進(jìn)行對比，驗證了數(shù)值模擬方法的可靠性，滿足工程要求。

4）建模過程中對模型進(jìn)行了簡化處理，對于高速旋轉(zhuǎn)下產(chǎn)生的空泡對槳葉壓力的影響、模擬條件與實驗情況差別以及網(wǎng)格劃分和數(shù)量等，都是影響數(shù)值模擬準(zhǔn)確性的因素，對以后進(jìn)行深入研究具有參考價值。

［1］LIEFVENDAHL M，TROENG C.Parallelization and load balancing of a dynamic mesh method for moving boundary CFD ［J/OL］.Swedish Defense Research Agency，F(xiàn)OI，November 2008，［2011-01-17］//powerlab·fsb·hr/ped/kturbo/…/oscic-08-Liefvendah/Martin.pdf.

［2］GUO C Y，MA N，YANG C J.Numerical simulation of a podded propulor in viscous flow ［J］.Journal of Hydrodynamics，2009，21（1）：71-76.

［3］喻欣.考慮升沉運動的螺旋槳水動力性能計算［D］.哈爾濱：哈爾濱工程大學(xué)，2008.

［4］李輝.基于CFD方法的全方向推進(jìn)器水動力性能研究［D］.哈爾濱：哈爾濱工程大學(xué)，2009.

［5］YOUNG Y L.A BEM for prediction of unsteady midchord face and／or back propeller cavitation［J］.Journal of Fluids Engineering，2001，6（123）：311-319.

［6］王國強，董世湯.船舶螺旋槳理論與應(yīng)用［M］.哈爾濱：哈爾濱工程大學(xué)出版社，2007.

［7］王超，黃勝，解學(xué)參.基于CFD方法的螺旋槳水動力性能預(yù)報［J］.海軍工程大學(xué)學(xué)報，2008，20（4）：107-112.

［8］孫存樓，王永生，李堅波.基于CFD的調(diào)距槳調(diào)距過程槳葉水動力性能研究 ［J］.哈爾濱工程大學(xué)學(xué)報，2010，31（2）：149-153.

［9］DIFELICE F.Experimental investigation of the propeller wake at different loading conditions by particle image velocimetry［J］.Journal of Ship Research，2004，48（2）：168-190.

［10］YAMAK K.Comparative calculations of propeller by surface panel method ［C］//Workshop of 20thITTC Propulsor Committee，Ship Research Institute，Supplement No 15，1993：57-66.