潛艇斜航工況下操舵水動(dòng)力及繞流場數(shù)值研究

郭海鵬，祁江濤，袁文鑫，李廣年

(1.寧波大學(xué) 海運(yùn)學(xué)院，浙江 寧波 315211；2.中國船舶科學(xué)研究中心，江蘇 無錫 214082)

0 引 言

操縱性是衡量潛艇綜合性能優(yōu)良的重要指標(biāo)之一，良好的操縱性能是潛艇安全航行的重要保證[1-2]。潛艇操縱工況下繞流特性極為復(fù)雜，不但有主附體結(jié)合部和流水孔處各種形式渦的產(chǎn)生與演化，而且還涉及艇體、舵以及螺旋槳之間的復(fù)雜相互干擾，導(dǎo)致針對潛艇操縱性的研究具有很高的難度。早期研究主要依靠模型試驗(yàn)方式進(jìn)行，對人力物力要求高而且試驗(yàn)設(shè)施復(fù)雜，難以在潛艇設(shè)計(jì)階段進(jìn)行廣泛應(yīng)用。

隨著計(jì)算機(jī)技術(shù)的發(fā)展，計(jì)算流體動(dòng)力學(xué)（computational fluid dynamics, CFD）方法在水下航行體操縱性預(yù)報(bào)中發(fā)揮越來越重要的作用。張楠等[3-5]提出了適用于水下航行器流場預(yù)報(bào)的數(shù)值方法，并開展了帶自由液面、全附體模型的數(shù)值模擬。孫銘澤等[6]實(shí)現(xiàn)了全附體潛艇平面運(yùn)動(dòng)機(jī)構(gòu)試驗(yàn)數(shù)值模擬，并將計(jì)算得到的水動(dòng)力導(dǎo)數(shù)與試驗(yàn)值以及直線拖曳試驗(yàn)或回轉(zhuǎn)試驗(yàn)的數(shù)值計(jì)算結(jié)果進(jìn)行了對比驗(yàn)證。龐永杰[7]以橢球體為研究對象計(jì)算了其水動(dòng)力導(dǎo)數(shù)，并通過與理論值進(jìn)行對比驗(yàn)證了所采用數(shù)值方法的有效性。焦玉超等[8]對不同X 形尾舵的潛艇繞流場進(jìn)行了數(shù)值模擬，分析了X 舵角度對潛艇性能的影響。翟朔等[9]對某潛艇模型的尾部水平操縱面分別進(jìn)行了共翼型設(shè)計(jì)和非共翼型設(shè)計(jì)，并計(jì)算了2 種操縱面產(chǎn)生的艇體水動(dòng)力和尾流特征。

從目前已有研究可以看出，應(yīng)用CFD 方法開展?jié)撏Р倏v性能預(yù)報(bào)及評估已經(jīng)成為領(lǐng)域內(nèi)的趨勢和研究熱點(diǎn)。本文基于CFD 平臺STAR-CCM+對潛艇斜航工況下的操舵水動(dòng)力及繞流場進(jìn)行數(shù)值研究，以國際上廣泛用于對比驗(yàn)證研究的SUBOFF 潛艇模型為研究對象，分別應(yīng)用雷諾平均（reynolds-averaged navier-stokes,RANS）與分離渦模擬（detached eddy simulation, DES）2 種方法對潛艇不同舵角下的斜拖試驗(yàn)進(jìn)行數(shù)值模擬，據(jù)此分析潛艇縱向力、橫向力、轉(zhuǎn)首力矩以及繞流場隨舵角的變化規(guī)律，揭示潛艇操縱工況下水動(dòng)力特性與繞流場間的內(nèi)在聯(lián)系，為水下航行體大舵角狀態(tài)操縱性能預(yù)報(bào)和評估提供技術(shù)支持。

1 數(shù)值計(jì)算方法

采用STAR-CCM+中的RANS 與DES 求解器開展相關(guān)數(shù)值計(jì)算工作，其中RANS 方法和DES 方法均采用剪切應(yīng)力傳輸模型（shear stress transfer, SST）k-ω湍流模型進(jìn)行構(gòu)造。SSTk-ω兩方程湍流模型為：

式中：ρ為流體密度；uj為速度矢量；xj為位置矢量；k和 ω分別為湍流動(dòng)能和湍流比耗散率；μ和 μt分別為層流和湍流黏性系數(shù)。源項(xiàng)P和函數(shù)F1的具體形式以及系數(shù) β*，γ，σk，σω和σω2可參考SSTk-ω模型的相關(guān)文獻(xiàn)[10]。

將RANS 湍流模型中的長度尺度使用DES 的長度尺度代替，就獲得了基于該湍流模型的混合RANS/LES方法即DES 方法。DES 的長度尺度由RANS 和LES 的長度尺度以如下方式混合得到，即

LES 長度尺度由網(wǎng)格間距得到，即

式中，Δ為局部網(wǎng)格間距，一般取網(wǎng)格單元與相鄰網(wǎng)格單元中心連線長度的最大值。

RANS 長度尺度為最近壁面距離，即

由DES 模型的構(gòu)造可知，當(dāng)RANS 長度尺度激活時(shí)，湍流模型為傳統(tǒng)的RANS 湍流模型，而當(dāng)LES 長度尺度激活時(shí)，湍流模型轉(zhuǎn)換為經(jīng)典的Smagorisnky 亞格子模型，從而實(shí)現(xiàn)RANS 和LES 方法的轉(zhuǎn)換。

2 數(shù)值計(jì)算設(shè)置

2.1 研究對象

以美國DARPA 潛艇模型SUBOFF-8 為研究對象，該潛艇由裸艇、指揮臺圍殼及4 個(gè)完全相同的尾翼（左右為升降舵，上下為方向舵）組成，主艇體總長Loa=4.356 m，垂線間長L=4.261 m，其中前體長1.016 m，平行中體長2.229 m，后體長1.111 m；艇體最大直徑D=0.508 m；指揮臺圍殼長0.368 m，高0.205 m，上部有2∶1 的橢圓形橫截面的頂蓋；4 個(gè)相同的尾翼截面為NACA0020，呈十字布置在尾部，頂端弦長為0.152 m，詳細(xì)幾何特征由文獻(xiàn)[11]提供，三維幾何形狀如圖1所示。

圖1 SUBOFF 幾何模型Fig.1 Geometric model of SUBOFF

2.2 計(jì)算工況

計(jì)算工況具體參數(shù)如表1 所示?？紤]潛艇在漂角β=-4°、舵角0°～40°范圍內(nèi)的水動(dòng)力性能，其中δ在0°～15°范圍內(nèi)數(shù)值計(jì)算結(jié)果可與Roddy[12]所做的公開實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對比，從而確定數(shù)值方法的可靠性。在STAR-CCM+中，漂角的處理方法是將SUBOFF 以艇體質(zhì)心所在的Z軸為旋轉(zhuǎn)軸，進(jìn)行逆時(shí)針旋轉(zhuǎn)4°，其他保持不變。

表1 工況設(shè)定Tab.1 Setup of study condition

2.3 計(jì)算域設(shè)置

計(jì)算域總長為7 倍艇長，首部距離進(jìn)口1.5 倍，尾部距離出口4.5 倍，左右及上下邊界距艇體均為2.0 倍艇長。對于邊界條件的設(shè)定，除計(jì)算域出流面采用壓力出口條件以外，其他邊界均采用速度入口條件，艇體所有表面采用無滑移壁面邊界條件。計(jì)算域尺度及邊界條件如圖2 所示。

圖2 計(jì)算域尺度及邊界條件Fig.2 Dimensions and boundary conditions of the computational domain

計(jì)算域采用STAR-CCM+中的切割體網(wǎng)格進(jìn)行離散。切割體網(wǎng)格屬于非結(jié)構(gòu)網(wǎng)格，適用于具有復(fù)雜外形模型的網(wǎng)格劃分。切割體網(wǎng)格在生成時(shí)可以通過控制體進(jìn)行任意區(qū)域的網(wǎng)格加密，同時(shí)還可以對幾何模型中的控制線、控制面進(jìn)行線加密及面加密。為了保證潛艇附近的網(wǎng)格質(zhì)量，對其表面進(jìn)行了加密；通過使用切割體網(wǎng)格單元，在艇體附近及其尾流進(jìn)行網(wǎng)格密化，從而更好地捕捉流動(dòng)的細(xì)節(jié)，獲得更準(zhǔn)確的流場信息。同時(shí)，在近壁面生成棱柱層網(wǎng)格以便于邊界層流動(dòng)的求解。通過調(diào)整棱柱層網(wǎng)格的總厚度、層數(shù)及增長比率來控制近壁面第一層網(wǎng)格高度，使近壁面Y+值在30～60 范圍內(nèi)。網(wǎng)格數(shù)量在500 萬左右，網(wǎng)格劃分情況如圖3 所示。

圖3 計(jì)算域網(wǎng)格劃分情況Fig.3 Gird generation in the computational domain

3 結(jié)果分析與討論

3.1 水動(dòng)力及力矩

圖4 給出了RANS 與DES 方法獲得的作用在潛艇上的縱向力X、橫向力Y及轉(zhuǎn)首力矩N。為了便于分析，均采用無因次形式表述：

圖4 數(shù)值結(jié)果與模型試驗(yàn)數(shù)據(jù)對比Fig.4 Comparison between the numerical results and the experimental data

式中：ρ為流體密度；V為航行速度；L為潛艇總長。

相對X′而言，無論RANS 還是DES 對阻力的預(yù)報(bào)精度都在±3%以內(nèi)；δ的變化會對阻力產(chǎn)生較大的影響，δ在0°～30°范圍內(nèi)，隨著舵角的增大，艇體所受的阻力以拋物線的趨勢增大；隨著δ的進(jìn)一步增大，RANS 計(jì)算結(jié)果在30°～40°范圍內(nèi)基本呈線性增大，DES 計(jì)算結(jié)果先有略微減小然后再增大。相對Y′而言，RANS 計(jì)算結(jié)果與試驗(yàn)值的吻合程度略好于DES 方法，RANS 計(jì)算的誤差均在10%以內(nèi)，DES 計(jì)算除了δ=15°時(shí)誤差偏大，其余工況均在12% 以內(nèi)；δ在0°～30°范圍內(nèi)，DES 與RANS 計(jì)算得到的橫向力隨著舵角的增大而近似呈現(xiàn)線性增大；δ在30°～35°范圍內(nèi)，RANS 及DES 計(jì)算的橫向力顯著下降，發(fā)生了明顯的舵失速現(xiàn)象；進(jìn)一步增大舵角，艇體所受的橫向力有微幅的增加，這表明舵效并未隨舵角進(jìn)一步增大而發(fā)生明顯改善。相對于N′而言，DES 計(jì)算在各工況下的計(jì)算精度均在±12% 以內(nèi)，而RANS 方法在δ=10°和δ=15°計(jì)算誤差較大，均超過20%；N′與Y′的變化規(guī)律保持基本一致。通過對水動(dòng)力分析可知，就對比現(xiàn)有的公開試驗(yàn)數(shù)據(jù)而言，RANS 方法在X′與Y′相關(guān)的計(jì)算精度比DES 方法略好，但是DES 方法對N′的計(jì)算精度明顯優(yōu)于RANS 方法。

3.2 泄出渦系結(jié)構(gòu)

圖5 給出了由RANS 及DES 方法獲得的不同舵角下艇尾區(qū)域的渦系結(jié)構(gòu)，其中渦系結(jié)構(gòu)采用Q準(zhǔn)則等值面表示，并以無量綱縱向速度Ux/U0著色。可以看出，隨著方向舵舵角的變化，艇尾區(qū)域渦系發(fā)展呈現(xiàn)出顯著的差異。在小舵角（δ=10°）下，舵背風(fēng)面主要以附著流為主，舵梢部有輕微渦泄出。隨著舵角增大（δ=20°），方向舵附近產(chǎn)生的渦系結(jié)構(gòu)明顯增多。方向舵背風(fēng)面靠近尾部區(qū)域開始產(chǎn)生泄出渦系。同時(shí)，由于偏轉(zhuǎn)方向舵的阻滯效應(yīng)，在艇身與方向舵的交匯處形成了明顯的馬蹄渦結(jié)構(gòu)。隨著舵角的進(jìn)一步增大（δ=30°），舵背風(fēng)面完全由泄出渦系覆蓋，而且馬蹄渦結(jié)構(gòu)也更為顯著。當(dāng)舵角達(dá)到δ=40°時(shí)，舵背風(fēng)面的泄出渦完全脫離舵表面發(fā)展，在下游形成了復(fù)雜的渦系結(jié)構(gòu)，與舵失速狀態(tài)相對應(yīng)。通過對比可知，在小舵角工況下，RANS 方法與DES 方法得到的艇尾區(qū)域渦系結(jié)構(gòu)差異不明顯，隨著舵角增大，二者的差異越發(fā)明顯。特別是發(fā)生大規(guī)模流動(dòng)分離的工況下，DES方法捕捉到的渦系結(jié)構(gòu)更加精細(xì)。

圖5 不同舵角工況下的Q=100 s-2 等值面圖Fig.5 Iso-surface of Q=100 s-2 under different rudder angles

3.3 流線分布

圖6 給出了DES 與RANS 方法獲得的舵葉中間位置處的舵附近流線分布。可知，在舵角δ=10°工況下，舵剖面附近流動(dòng)主要以附著流為主，而且2 種方法得到的流場沒有明顯的差異。在舵角δ=20°工況下，舵剖面背流面發(fā)生了明顯的流動(dòng)分離現(xiàn)象，其中RANS方法呈現(xiàn)出的流動(dòng)分離更為明顯，導(dǎo)致該區(qū)域的速度明顯下降。在舵角δ=30°工況下，舵剖面的整個(gè)背流面均發(fā)生了流動(dòng)分離，2 種數(shù)值方法獲得的流線分布差異也更加明顯。RANS 方法中背流面處的流動(dòng)分離更加明顯，而DES 方法捕捉到的流動(dòng)分離并未隨舵角明顯增強(qiáng)。另外，2 種方法捕捉得到流場均在舵下游呈現(xiàn)出明顯的低速區(qū)。在舵角δ=40°工況下，舵剖面背流面處的流動(dòng)分離顯著增強(qiáng)，對附近流場產(chǎn)生的擾動(dòng)也十分明顯，在下游形成了更為顯著的低速區(qū)，這一流場特征與舵失速工況下的阻力增加及升力下降相對應(yīng)。綜上，DES 方法在流場捕捉方面具有明顯優(yōu)勢，能清晰地捕捉艇尾區(qū)域的復(fù)雜流場信息。

圖6 不同舵角工況下舵附近流線分布Fig.6 Streamlines around the rudder under different rudder angles

3.4 壓力分布

圖7 給出了RANS 方法與DES 方法獲得的不同舵角工況下的舵剖面壓力分布情況。壓力采用無因次系數(shù)Cp表示，Cp=p/(0.5ρV2)。橫坐標(biāo)由無因次量x/C表示，坐標(biāo)原點(diǎn)位于該舵前緣處，無因次量“1”表示舵尾緣。在舵角δ=10°工況下，舵剖面壓力呈現(xiàn)出明顯的吸力面和壓力面，正壓和負(fù)壓的峰值均位于前緣附近，而且2 種方法獲得的壓力分布較為一致。在舵角δ=20°工況下，DES 方法得到的吸力面負(fù)壓區(qū)較RANS方法更為顯著，但整體分布差異很小。在舵角δ=30°工況下，吸力面負(fù)壓區(qū)更為顯著，同時(shí)發(fā)生流動(dòng)分離的區(qū)域壓力分布較為平緩。DES 方法得到的吸力面負(fù)壓區(qū)明顯強(qiáng)于RANS 方法，這也與該舵角工況下DES計(jì)算得到的舵橫向力高于RANS 計(jì)算結(jié)果相對應(yīng)。在舵角δ=40°工況下，吸力面負(fù)壓區(qū)峰值明顯下降，而且發(fā)生流動(dòng)分離的區(qū)域壓力分布呈平直狀，導(dǎo)致吸力面與壓力面之間的壓力差顯著下降，這與舵發(fā)生失速導(dǎo)致的舵效下降相對應(yīng)?？傮w而言，2 種數(shù)值方法得到的壓力分布差異相對較小，與潛艇操舵水動(dòng)力性能的一般規(guī)律較為吻合。

圖7 不同舵角工況下舵剖面壓力分布圖Fig.7 Pressure distribution on the rudder section under different rudder angles

4 結(jié) 語

本文針對潛艇斜航工況下的操舵水動(dòng)力性能及繞流場進(jìn)行數(shù)值模擬研究，主要得到以下結(jié)論：

1）潛艇斜航工況下的縱向力、橫向力及轉(zhuǎn)首力矩隨著舵角的變化呈現(xiàn)復(fù)雜的變化規(guī)律，在潛艇操舵水動(dòng)力預(yù)報(bào)精度方面DES 方法整體上優(yōu)于RANS 方法；

2）DES 方法對于潛艇操舵工況下艇尾區(qū)域的繞流場捕捉更為精細(xì)，在潛艇操縱工況下的水動(dòng)力學(xué)機(jī)理探討方面有明顯的優(yōu)勢，在潛艇復(fù)雜操縱運(yùn)動(dòng)分析和性能評估方面具有很好的工程應(yīng)用潛力。