斜流中艇后螺旋槳水動(dòng)力數(shù)值計(jì)算方法

張文照，肖昌潤(rùn)

(1.海軍工程大學(xué),湖北 武漢 430033； 2.海軍后勤技術(shù)裝備研究所,北京 100072)

0 引 言

斜流條件下，艇后螺旋槳水動(dòng)力性能研究一直是進(jìn)行潛艇操縱性研究的關(guān)鍵，目前主要通過(guò)模型試驗(yàn)的方法來(lái)進(jìn)行研究。但模型試驗(yàn)存在尺度、伴流場(chǎng)的模擬、結(jié)果換算和經(jīng)費(fèi)等眾多問(wèn)題的困擾，受到極大限制。數(shù)值方法可以克服上述試驗(yàn)方法存在的不足，勢(shì)流方法是螺旋槳性能分析中常用的數(shù)值方法。在非均勻流中螺旋槳性能分析也常采用勢(shì)流方法，比較突出的有Kinnas[1]和Hosino[2]。國(guó)內(nèi)有陳家棟[3]、熊鷹[4]、劉小龍[5]等用勢(shì)流方法分析給定不均勻流場(chǎng)中螺旋槳性能。

但采用勢(shì)流方法處理艇后螺旋槳這種多體結(jié)構(gòu)的水動(dòng)力問(wèn)題變得十分困難，并且勢(shì)流方法在處理非均勻進(jìn)流中螺旋槳性能時(shí)，需要處理庫(kù)塔條件和假定尾渦形狀，有一定的人為因素。而粘流CFD方法則更方便處理復(fù)雜的幾何模型，實(shí)現(xiàn)對(duì)流場(chǎng)的真實(shí)模擬。Funeno[6]用粘流CFD方法分析了大側(cè)斜槳非定常性能，WATANABE[7]使用商業(yè)軟件Fluent研究了非均勻伴流中的螺旋槳性能。當(dāng)前粘性CFD的迅速發(fā)展，使得用CFD方法研究艇后螺旋槳的水動(dòng)力性能變得越來(lái)越容易。本文應(yīng)用Fluent軟件計(jì)算了斜流條件下艇后螺旋槳的水動(dòng)力，為開(kāi)展?jié)撏Р倏v性研究奠定了基礎(chǔ)。

1 研究對(duì)象

本文所研究的潛艇模型為加拿大DRDC潛艇水動(dòng)力研究計(jì)劃中的通用艇模[8]，模型長(zhǎng)4.36 m，直徑0.49 m。螺旋槳模型為5葉槳，直徑0.22 m，盤(pán)面比0.5。圖1給出了艇體帶槳的幾何模型。

圖1 計(jì)算所用帶槳潛艇模型Fig.1 Submarine model

2 計(jì)算方法描述

本文采用基于粘流理論的數(shù)值方法模擬艇后螺旋槳周?chē)牧鲃?dòng)，以商用軟件Fluent為求解器求解不可壓三維RANS方程計(jì)算艇后螺旋槳的推力和扭矩。

2.1 控制方程

螺旋槳水動(dòng)力計(jì)算所選取的坐標(biāo)系為固定在螺旋槳上的旋轉(zhuǎn)直角坐標(biāo)系，定義坐標(biāo)軸x為螺旋槳的旋轉(zhuǎn)軸，指向下游為正。控制方程采用Reynolds平均的NS(RANS)方程，其不可壓縮的無(wú)量綱形式為：

(1)

式中：ui和uj為絕對(duì)速度各分量；p為壓力；Re為雷諾數(shù)；

(4)

式中：νt為湍流運(yùn)動(dòng)引起的渦運(yùn)動(dòng)粘性系數(shù)；

采用SSTk-ω湍流模型來(lái)計(jì)算雷諾應(yīng)力項(xiàng)，以實(shí)現(xiàn)控制方程組的封閉。SSTk-ω湍流模型混合了k-ω模型計(jì)算近壁面區(qū)域粘性流動(dòng)的可靠性和k-ε模型計(jì)算遠(yuǎn)場(chǎng)自由流動(dòng)的精確性，并在螺旋槳性能計(jì)算中反復(fù)驗(yàn)證其較高的精度。

控制方程采用有限體積法離散，有限體積法具有良好的守恒性，能夠采用各種形式的網(wǎng)格以適應(yīng)復(fù)雜的幾何輪廓?？刂品匠讨械膶?duì)流項(xiàng)采用二階迎風(fēng)格式離散，速度與壓力的耦合迭代采用SIMPLE 算法。離散后形成的代數(shù)方程用逐點(diǎn)Gauss-Seidel迭代法求解，并且采用代數(shù)多重網(wǎng)格方法加快求解的收斂速度。

2.2 滑移網(wǎng)格法

在進(jìn)行艇后螺旋槳水動(dòng)力計(jì)算時(shí)，螺旋槳旋轉(zhuǎn)，艇體(包含附體)靜止。對(duì)于這種問(wèn)題，F(xiàn)luent軟件有混合面法、多參考系法、滑移網(wǎng)格法和動(dòng)網(wǎng)格法4種處理方法，前2種方法都假設(shè)流動(dòng)為定常，后2種方法都假定流動(dòng)為非定常。動(dòng)網(wǎng)格法需要配合諸如層生成吞并等技術(shù)，設(shè)置復(fù)雜，而滑移網(wǎng)格法在處理定轉(zhuǎn)子問(wèn)題上更易于實(shí)現(xiàn)。本文對(duì)于艇后螺旋槳周?chē)鲌?chǎng)的計(jì)算采用滑移網(wǎng)格法，包含艇體和附體的部分設(shè)置為靜止區(qū)域，螺旋槳部分為滑動(dòng)區(qū)域，計(jì)算網(wǎng)格隨螺旋槳轉(zhuǎn)動(dòng)。在2個(gè)區(qū)域的交界面上物理量的傳遞采用插值方法處理。

2.3 計(jì)算域與網(wǎng)格劃分

計(jì)算域的進(jìn)口邊界設(shè)置為艇首上游0.8L處(L為艇模長(zhǎng)度)的半球面，出口邊界設(shè)置在螺旋槳下游2L處，柱面外邊界也為進(jìn)口，柱面半徑取1L。整個(gè)計(jì)算域分割成兩大部分，如圖2所示，包含螺旋槳的那一部分區(qū)域?yàn)樾D(zhuǎn)區(qū)域，其他區(qū)域?yàn)殪o止區(qū)域。

圖2 計(jì)算區(qū)域示意圖Fig.2 Computing field

包含螺旋槳的旋轉(zhuǎn)部分由于螺旋槳復(fù)雜的幾何外形，劃分為六面體網(wǎng)格較難實(shí)現(xiàn)，所以對(duì)旋轉(zhuǎn)區(qū)域仍采用四面體網(wǎng)格劃分(見(jiàn)圖3)。

靜止和旋轉(zhuǎn)區(qū)域交界面上的網(wǎng)格劃分見(jiàn)圖4。

圖3 旋轉(zhuǎn)區(qū)域的網(wǎng)格劃分Fig.3 Mesh grid of rotation field

圖4 靜止和旋轉(zhuǎn)區(qū)域交界面處的網(wǎng)格劃分Fig.4 Mesh grid of the interface of static field with rotation field

包含吊艇體、圍殼和十字舵的靜止區(qū)域在進(jìn)行合理的分塊后較易實(shí)現(xiàn)六面體網(wǎng)格劃分。靜止區(qū)域的分塊如圖5所示，分塊會(huì)給網(wǎng)格劃分帶來(lái)一個(gè)問(wèn)題，即指揮臺(tái)圍殼和尾翼表面一層的網(wǎng)格尺寸會(huì)延續(xù)到整個(gè)計(jì)算域的外邊界，然而在計(jì)算域的外邊界，其他方向上的網(wǎng)格劃分已經(jīng)很疏，這樣就會(huì)使得網(wǎng)格長(zhǎng)寬比過(guò)大。于是，在網(wǎng)格劃分用如圖6所示的處理方法，即采用“斜向擴(kuò)展”的方式來(lái)防止網(wǎng)格在計(jì)算域外邊界上長(zhǎng)寬比過(guò)大。

圖5 靜止區(qū)域網(wǎng)格劃分Fig.5 Mesh grid of static field

圖6 指揮臺(tái)圍殼網(wǎng)格斜向擴(kuò)展法Fig.6 Oblique extension method

圖7 圍殼和尾翼端面網(wǎng)格劃分Fig.7 Mesh grid of shell and end face

此外，對(duì)于圍殼和尾翼端面網(wǎng)格的劃分，將上端面先分成如圖7所示的2個(gè)端部面和1個(gè)中間面，當(dāng)支柱上端面邊的節(jié)點(diǎn)固定后，這樣的分塊能有效控制網(wǎng)格層數(shù)。

2.4 邊界條件和初始條件的設(shè)置

計(jì)算域的進(jìn)口為圖2中右側(cè)的圓面和徑向的柱面，進(jìn)口設(shè)置為速度進(jìn)口，給定來(lái)流速度的大小和方向。計(jì)算域的出口為左側(cè)的圓面，出口設(shè)置為壓力出口，給定壓力大小。螺旋槳、艇體、指揮臺(tái)圍殼和尾附體表面設(shè)為物面無(wú)滑移條件，靜止區(qū)域與旋轉(zhuǎn)區(qū)域相重疊的面設(shè)置為交界面，完成2個(gè)區(qū)域的數(shù)據(jù)傳遞。用進(jìn)口速度初始化整個(gè)流場(chǎng)。

3 計(jì)算結(jié)果的比較與分析3.1 0°攻角不同進(jìn)速非定常計(jì)算結(jié)果

首先分別計(jì)算了艇體0°攻角下，螺旋槳艇后進(jìn)速Jb為0.8，1.0，1.2，1.4，1.6工況下，螺旋槳轉(zhuǎn)動(dòng)1周過(guò)程中推力和扭矩的變化值。選取設(shè)計(jì)點(diǎn)Jb=1.0，畫(huà)出該工況下螺旋槳轉(zhuǎn)動(dòng)1周所產(chǎn)生的推力和扭矩的變化曲線(xiàn)(見(jiàn)圖8和圖9)。通過(guò)曲線(xiàn)可以看出螺旋槳轉(zhuǎn)動(dòng)1周，推力和扭矩呈現(xiàn)出規(guī)律的周期性，周期數(shù)等于螺旋槳的頁(yè)數(shù)。螺旋槳轉(zhuǎn)速通常很高，因此推力和扭矩的這種周期性的脈動(dòng)相對(duì)于潛艇運(yùn)動(dòng)而言十分快速，因此在研究潛艇操縱運(yùn)動(dòng)時(shí)，給定艇體運(yùn)動(dòng)狀態(tài)的情況下(如潛艇航速、縱傾角等)，可以忽略螺旋槳推力和扭矩的這種周期脈動(dòng)，而用一個(gè)周期的平均值作為艇體該運(yùn)動(dòng)狀態(tài)下螺旋槳產(chǎn)生的推力和力矩計(jì)算值。

表1和表2分別列出了螺旋槳艇后進(jìn)速Jb為0.8，1.0，1.2，1.4，1.6工況下計(jì)算得到平均推力系數(shù)Kt和平均扭矩系數(shù)10Kq，并和試驗(yàn)值做比較，計(jì)算誤差在3%以?xún)?nèi)。

圖8 Jb=1.0推力系數(shù)1周內(nèi)的變化Fig.8 Trust coefficient in a period when Jb=1.0

圖9 Jb=1.0扭矩系數(shù)1周內(nèi)的變化Fig.9 Torque coefficient in a period when Jb=1.0

進(jìn)速JbKt計(jì)算值Kt試驗(yàn)值相對(duì)誤差/%080275028019100242024408120211020816140175017123160136013325

表2 不同進(jìn)速螺旋槳平均扭矩系數(shù)

3.2 艇體攻角對(duì)螺旋槳水動(dòng)力性能影響計(jì)算

為了研究艇體攻角對(duì)艇后螺旋槳推力和扭矩的影響，計(jì)算設(shè)計(jì)進(jìn)速Jb=1.0， 艇速u(mài)=0.7 m/s，艇體處于0°，5°，10°，15°，30°攻角工況下艇后螺旋槳平均推力、扭矩系數(shù)(見(jiàn)圖10)。由于整個(gè)計(jì)算對(duì)象幾何對(duì)稱(chēng)特性，負(fù)攻角與正攻角下的槳模水動(dòng)力對(duì)稱(chēng)。

與試驗(yàn)值相對(duì)比，計(jì)算值的相對(duì)誤差隨著攻角的增大而增大，總體控制在7%以?xún)?nèi)，該誤差是由于在大攻角時(shí)，艇體流體分離比較嚴(yán)重，而計(jì)算網(wǎng)格不能有效地捕捉到該現(xiàn)象導(dǎo)致的。

圖10 Jb=1.0時(shí)艇后螺旋槳推力、扭矩系數(shù)隨艇體攻角的變化Fig.10 Behind-boat trust of the propeller with different pitch angle when Jb=1.0

從圖10可看出，螺旋槳推力和扭矩是隨著艇體攻角呈現(xiàn)出先減小后增大的趨勢(shì)。究其原因，是因?yàn)楫?dāng)艇體產(chǎn)生攻角時(shí)，螺旋槳進(jìn)流速度會(huì)減小，伴流分?jǐn)?shù)會(huì)增大，當(dāng)攻角比較小的時(shí)候伴流分?jǐn)?shù)的增大導(dǎo)致的螺旋槳實(shí)效進(jìn)速的提高占主要作用，因此在小攻角時(shí)螺旋槳推力和扭矩均隨著攻角的增加而降低。但是，隨著攻角的不斷變大，攻角導(dǎo)致的螺旋槳進(jìn)流速度的降低越來(lái)越嚴(yán)重，進(jìn)而超過(guò)伴流分?jǐn)?shù)減小的影響，此后，螺旋槳的推力和扭矩隨著攻角的增大而有所增加。

4 結(jié) 語(yǔ)

本文采用CFD方法對(duì)艇后螺旋槳的水動(dòng)力性能進(jìn)行計(jì)算，并對(duì)計(jì)算結(jié)果進(jìn)行比較和分析后得出結(jié)論如下：

1)本文采用的數(shù)值計(jì)算方法對(duì)帶攻角的艇后螺旋槳推力和扭矩是有效的；

2)隨著攻角的增大，螺旋槳推力和扭矩的計(jì)算誤差有所增加，總體控制在7%以?xún)?nèi)；

3)在同一進(jìn)速下，艇后螺旋槳推力和扭矩隨著艇體攻角的增大呈現(xiàn)出先減小后增大的趨勢(shì)。

[1] KINNAS S A,HSIN C J, KEENAN D.A potential based panel method for the unsteady flow around open and ducted propellers[C].Proceedings of the 8thSymposium on Naval Hydrodynamics,1990.

[2] HOSHINO T.Hydrodynamic analysis of propeller in unsteady flow using a surface panel method[C].Journal of Society of Naval Architectures of Japan,1993.

[3] 陳家棟,董世湯.非定常螺旋槳表面壓力面元法預(yù)報(bào)計(jì)算[J].中國(guó)造船，1998(1)：21-26.

[4] 熊鷹.非均勻流中螺旋槳空泡和脈動(dòng)壓力地?cái)?shù)值和試驗(yàn)研究[D]，武漢：武漢理工大學(xué)，2002.

[5] 劉小龍,王國(guó)強(qiáng).螺旋槳非定常性能預(yù)估的面元法[J].船舶力學(xué),2006,10(2).

[6] FUNENO I. Analysis of unsteady viscous flows around a highly skewed propeller(in Japanese)[J]. Kansai Society of Naval Architects,2002,237：39-45.

[7] WATANABE T.Simulation of steady and unsteady cavitation on a marine propeller using a RANS code，F(xiàn)ifth International[J].Symposium on Cavitation (CAV2003)，Osaka, Japan,2003：1-4.

[8] FOURNIER E.Wind tunnel investigation of the DREA po-wered MkI submarine model[R].CANADA:DREA,1993.