基于CFD的潛艇阻力及流場數值計算

涂海文，孫江龍

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

基于CFD的潛艇阻力及流場數值計算

涂海文，孫江龍

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

運用雷諾平均N－S方程，使用CFD前處理軟件ICEM CFD劃分流場網格，采用RNG k－ε湍流模型，實現了對裸潛體、帶指揮臺圍殼艇體、帶十字尾翼艇體、全附體潛艇4種模型的阻力及粘性流場的數值模擬。通過數值模擬，得到了潛艇表面壓力分布情況和附體附近流場的一些特性，為進一步優化潛艇的艇型和分析潛艇的流噪聲打下了基礎。而阻力的對比在一定程度上驗證了數值模擬的可靠性。

潛艇;CFD;RNG k－ε模型;摩擦阻力

0 引言

潛艇周圍的流場特別是尾流場的流動特性，不僅對潛艇的水動力性能產生直接的影響，而且引發的流動是潛艇水動力噪聲的主要來源之一，對潛艇隱蔽性有重要影響。當航速和排水量一定時，流場的特性與潛艇的主尺度和艇體形狀密切相關。因此，開展潛艇幾何形狀與流場特性的關聯研究，尋求能在已知艇體形狀特征下精確預報流場特性的數值計算方法，為潛艇概念設計階段得出阻力小、噪聲低的最佳潛艇設計方案提供評價依據，一直是船舶水動力學領域中一個重要的前沿研究課題。

長期以來，經典流體理論無法在工程中直接應用，船舶流體力學工作者主要依靠水池實驗來粗略地了解船舶流體性能。但船模實驗耗資巨大，實驗周期長，且水池實驗受到尺度效應、外界干擾等影響，船模實驗的弊端顯而易見。隨著計算機技術的迅猛發展，數值計算方法逐漸成為研究船舶和潛艇流場的另一有效方法［1］。模型試驗和CFD數值仿真模擬各有自己特定的優勢。現在，科學工作者們往往將CFD仿真結果與模型試驗結果相互印證，使這2個方面互相促進和發展。

近年來，潛艇流場的數值模擬得以蓬勃發展。趙峰，周連第［2］運用復雜流場多塊耦合計算方法，對由于指揮臺而造成的潛艇后體不均勻流場作了較成功的數值模擬。張楠等［4］用數值方法預報了研究用美國DARPA潛艇模型SUBOFF與CSSRC潛艇模型SM－x的阻力與尾流場，并與試驗結果進行對比。李艷和姚震球［6］運用數值方法對帶附體潛艇尾流場進行模擬與驗證。操盛文和吳方良［7］采用不同數量的網格分別對不同尺度的SUBOFF模型在高雷諾數條件下的流場進行數值模擬，研究網格數量和艇體主尺度大小對高雷諾數條件下潛艇阻力計算結果的影響。

本文建立了裸艇體、帶指揮臺圍殼艇體、帶十字尾翼艇體、全附體潛艇4個不同的潛艇模型，然后進行CFD的仿真和對比計算。

1 數學模型

1.1 控制方程

不可壓縮流體的基本控制方程由連續性方程和RANS方程組成，其張量形式為:

1.2 湍流模型

上面的方程組是不封閉的，需要采用相應的湍流模型。本文采用RNG κ－ε湍流模型，該湍流模型數學表達可以寫成下面的形式［1］:

湍流動能k方程為

在RNG κ－ε模型中，通過在大尺度運動和修正后的粘度項體現小尺度的影響，而使這些小尺度運動有系統地從控制方程中去除。所得到的k方程和ε方程與標準κ－ε模型非常相似。

與標準κ－ε模型比較發現，RNG κ－ε模型修正了湍動粘度，考慮了平均流動中的旋轉及旋流流動情況，從而可以更好地處理高應變率及流線彎曲程度較大的流動。

需要注意的是，RNG κ－ε模型仍是針對充分發展的湍流有效，即是高Re數的湍流計算模型，而對近壁區內的流動及Re數較低的流動，必須使用壁面函數法或低Re數的κ－ε模型來模擬。

1.3 數值計算方法

采用有限體積法(FVM)對控制方程進行離散，用SIMPLEC方法計算壓力速度關聯方程。時間項采用中心差分格式，對流擴散項采用二階迎風格式。

2 計算模型

2.1 計算對象

建立如圖1所示的帶全附體的潛艇計算模型［3］。潛艇的指揮臺機翼采用的是3:2橢圓首部和NACA0020機翼尾部的復合機翼，則相應的機翼厚度為330 mm，潛艇的十字尾翼采用了NACA0020尾翼。潛艇總長度為995 mm，最大半徑為45.6 mm，指揮臺機翼距潛艇首部長度為220 mm，機翼后面的平行中體長度為240 mm，指揮臺機翼的弦長為140 mm，指揮臺高60 mm，十字尾翼最高點距軸中心線55 mm，弦長分布在軸向910～960 mm之間。取潛艇中心線與頭部的交點為原點，沿潛艇長度方向定義為x軸正向，機翼附體的高度方向為y軸正向，建立笛卡爾直角坐標系［1］。

圖1 全附體潛艇計算模型Fig.1 The computing model of submarine with all appendages

計算域尺寸選取參照圖2所示。為了使數值計算結果更接近無限空間內的三維繞流，計算域應取得足夠大。為此，上游邊界距原點0.3倍艇長，下游邊界距原點1.7倍艇長，上、下側邊界距軸線400 mm。

圖2 計算域xy平面圖Fig.2 Calculation domain of xy plane

光艇體模型、帶指揮臺圍殼艇體模型、帶十字尾翼艇體模型均可由全附體模型得到，它們的計算域與全附體計算域相同。

2.2 計算網格

劃分網格時主要遵循的原則是:靠近潛艇指揮臺、十字尾翼和潛艇表面及交接部等區域采取比較細密的網格，由于壁面及附體附近流場比較復雜，含有分離、旋渦及其脫落等現象，因此在這一區域需要布置較為密集的網格節點，以有利于研究潛艇指揮臺交接部和十字尾翼處的流場分析。而在遠離壁面及附體的區域，則可布置得較稀疏一些，以節省整個計算的工作量。

網格工具使用的是ICEMCFD網格劃分軟件。采用非結構性的四面體網格和三棱柱網格，四面體網格適合對結構復雜的幾何模型進行快速高效的網格劃分，三棱柱網格可以更好地模擬邊界層效應。在靠近潛艇指揮臺和十字尾翼的區域必須通過網格加密形成更為細密的網格，以便于更好地對這些區域的流場進行模擬。網格圖如圖3～圖5所示。

2.3 邊界條件

對于邊界條件的定義，為了更好地模擬無限流場的環境，在整個流域中，除了出口，其余2個邊界均設置為速度入口。數字計算域的前端面I和側表面M為來流入口，后端面O為來流出口，見圖6。

1)層流邊界條件設置

潛艇、指揮臺和十字尾翼的壁面設置為固壁條件，滿足無滑移條件。

采取層流計算模型時，在非定常計算過程中，每個時步取0.01 s，計算時長為10 s。

2)湍流邊界條件設置

計算域來流邊界設置為速度入口velocity-inlet，u=U0，v=w=0。其中初始來流速度U0=6.6 m/s，同樣以潛艇長作為特征長度，則對應的雷諾數Re=6.6×106。計算域出口處邊界條件設置為自由出流outflow。

潛艇、指揮臺和十字尾翼設置為壁面邊界條件wall，滿足無滑移條件，近壁區修正采用加強型壁面函數 Enhanced Wall Treatment。

在非定常計算過程中，每個時步取0.003 s，計算時長為3 s。

3 計算結果分析

3.1 摩擦阻力分析

潛艇運動中所受的總阻力Rt是所有流體作用力沿運動方向的合力亦即潛艇表面上所有微面積ds上切應力τ和壓力p在運動方向的合力，即

式中:前1項積分表示由作用在潛艇表面上切應力所造成的阻力，即摩擦阻力Rf;第2項積分表示由作用在潛艇表面上壓力所造成的阻力，即粘壓阻力Rpv。(τ，x)表示表面剪切應力與 x方向的夾角，(p，x)表示表面壓力與x方向的夾角。工程中一般定義摩擦阻力系數及粘壓阻力系數如下:

在工程實踐中，對于船舶摩擦阻力的預報通常使用“相當平板假定”計算，即利用與船體表面積相當的平板的摩擦阻力來代替船舶的摩擦阻力。

造船界常用的摩擦阻力公式為(ITTC57):

從表1可看出，數值模擬出來的摩擦阻力系數與用經驗公式算出來的摩擦阻力系數相比，誤差較小，它們的差別在于數值模擬中考慮了潛艇三維船體表面彎曲度的影響，而經驗公式未考慮此影響，因此本文的數值模擬結果是合理及可以接受的。從另一方面也驗證了數值模擬的可靠性。大，在指揮臺圍殼和十字尾翼處，壓力產生了突變。而對于帶指揮臺圍殼艇體模型，在湍流和層流2種情況下的差別在于湍流的最大壓力系數比層流的要大一些，再就是指揮臺尾部的“高”壓力區域比層流的要小，也就是說湍流情況下壓力恢復至潛艇光體相同的壓力比層流的要快。從圖中還可以看出，4種模型模擬的情況也比較吻合。

1)裸艇體模型靜壓力系數分布圖

?

3.2 壁面靜壓力系數分析

圖7～圖14依此描述了裸艇體模型、帶指揮臺圍殼艇體模型、帶十字尾翼艇體模型、全附體模型4種模型分別在湍流和層流狀態下y=0剖面上下壁面沿x向靜壓力系數分布情況。從圖中可以看出，湍流和層流的計算結果比較相似，均是潛艇頭部壓力比較

2)帶指揮臺圍殼艇體模型靜壓力系數分布圖

圖9 湍流下y=0剖面沿X向靜壓力系數分布Fig.9 The distribution of static pressure coefficient on y=0 profile along the X direction under the turbulence

圖10 層流下y=0剖面沿X向靜壓力系數分布Fig.10 The distribution of static pressure coefficient on y=0 profile along the X direction under the laminar

3)帶十字尾翼艇體靜壓力系數分布圖

4)全附體潛艇模型靜壓力系數分布圖

圖13 湍流下y=0剖面沿X向靜壓力系數分布Fig.13 The distribution of static pressure coefficient on y=0 profile along the X direction under the turbulence

圖14 層流下y=0剖面沿X向靜壓力系數分布Fig.14 The distribution of static pressure coefficient on y=0 profile along the X direction under the laminar

圖15和圖16將湍流和層流2種情況下全附體潛艇模型靜壓力系數的分布云圖形式清晰，變化的過程很直觀。

3.3 流場分析

圖17～圖19展示了指揮臺圍殼處的流線分布。從圖中可以看出，在圍殼頭部附近流體出現分流現象，速度方向發生改變，然后在圍殼尾部合流，這就是指揮臺圍殼頭部壓力很大的原因。

圖20展示了十字尾翼在y=0剖面處速度分布云圖。圖21和圖22是十字尾翼在y=0剖面和z=0剖面處的流線圖。在y=0剖面u速度分布圖上觀察到了負速度區，而在流線圖21和圖22上，表現為明顯的渦。

4 結語

1)運用CFD數值模擬技術很好地預報了潛艇的摩擦阻力系數，而通過對潛艇壁面壓力的數值預報，可以發現在潛艇的頭部，指揮臺圍殼的頭部以及十字尾翼的前端面壓力均比較大，從而掌握了潛艇上一些顯著地方的壁面壓力分布情況，這些都為潛艇快速性和潛艇艇型優化提供了非常好的參考。

2)從層流和湍流2種狀態下潛艇的壓力分布對比當中可以看出，在2種不同狀態下潛艇的壓力分布曲線的走勢具有一致性，只是壓力的大小和一些局部地方的壓力分布有所不同。

3)從流場分析中可以看出，在附體的附近流線比較密集，在潛艇的尾翼處有明顯的渦現象產生，這就為后面潛艇流噪聲的分析作了很好的鋪墊。

［1］陳小鄒.翼型附體與平板交接部流場分析［J］.艦船科學技術，2010，32(10):8 －11.

CHEN Xiao-zou.Junction flows analysis around flat plate mounted a foilappendage［J］.Ship Science and Technology，2010，32(10):8 －11.

［2］趙峰，周連第.潛艇含指揮臺附體區域周圍粘性流場的多塊耦合計算［J］.水動力學研究與進展(A輯)，1996，11(4):448－458.

ZHAO Feng，ZHOU Lian-di.Multiblock coupled computation of the viscous flow field around a submarine sail zone［J］.Journal of Hydrodynamics，1996，11(4):448 －458.

［3］LEE Y T，HUANG T T，DAVID S C H.Flow prediction about an appended body of revolution［J］.Taylor Research Center Bethesda，MD 20084 －5000.

［4］張楠，沈泓萃，姚惠之.潛艇阻力與流場的數值模型與驗證及艇型的數值優化研究［J］.船舶力學，2005，9(1):1－13.

ZHANG Nan，SHEN Hong-cui，YAO Hui-zhi.Validation of numerical simulation on resistance and flow field of submarine and numerical optimization of submarine hull form［J］.Journal of Ship Mechanics，2005，9(1):1 －13.

［5］張楠，沈泓萃，姚惠之，等.帶流水孔潛體流場數值模擬［J］.船舶力學，2004，8(1):1 －11.

ZHANG Nan，SHEN Hong-cui，YAO Hui-zhi，et al.Numerical simulation of flow around a submerged body with flood hole［J］.Journal of Ship Mechanics，2004，8(1):1 －11.

［6］李艷，姚震球.帶附體潛艇尾流場的數值模擬與驗證［J］.江蘇科技大學學報(自然科學版)，2006，20(2):7－12.

LI Yan，YAO Zhen-qiu.Numerical simulation and validation of wake field of submarine with stern appendage［J］.Journal of Jiangsu University of Science and Technology(Natural Science Edition)，2006，20(2):7 －12.

［7］操盛文，吳方良.尺度效應對全附體潛艇阻力數值計算結果的影響［J］.中國艦船研究，2009，4(1):33 －37.

CAO Sheng-wen，WU Fang-liang.Investigation of scaling effects on numerical computation of submarine resistance［J］.Chinese Journal of Ship Research，2009，4(1):33 －37.

［8］YANG Chi，LOHNER R.Prediction of flows over an axisymmetric body with appendages［A］.The 8 th International Conferenceon NumericalShip Hydrodynamics，Busan，Korea，September，2003.

［9］李新汶，陳源，王文琦.CFD在潛艇艇附體與艇體連接形式研究中的應用［J］.船舶力學，2003，7(5):28－32.

LI Xin-wen，CHEN Yuan，WANG Wen-qi.Applicationof CFD to junction form of strern appendages.［J］.Journal of Ship Mechanics，2003，7(5):28 －32.

［10］王志博，姚惠之，張楠.指揮臺圍殼對潛艇尾流影響的計算研究［J］.船舶力學，2009，13(2):196 －202.

WANG Zhi-bo，YAO Hui-zhi，ZHANG Nan.Calculation and analysis of the effects of sail form on submarine wake［J］.Journal of Ship Mechanics，2009，13(2):196 －202.

［11］盧云濤，張懷新，潘徐杰.全附體潛艇的流場和流噪聲的數值模擬［J］.振動與沖擊，2008，27(9):142 －189.

LU Yun-tao，ZHANG Huai-xin，PAN Xu-jie.Numerical simulation of flow-field and flow-noise of a fully appendage submarine［J］.Journal of Vibration and Shock，2008，27(9):142－189.

Numerical analysis of resistance and flow field of submarine based on CFD

TU Hai-wen，SUN Jiang-long
(College of Naval Architecture and Ocean Engineering，Huazhong University of Science and Technology，Wuhan 430074，China)

The article numerically simulated the resistance and flow field of bare submarine body with different appendages and submarine body full-appended with RNG k － ε model.And solving solutions with RANS equations.The girds were divided by ICEM CFD of CFD pretreatment software.Through numerical simulation，it obtained the submarine surface pressure distribution and some characteristics of the flow field near appendages.This is further to lay the foundation for optimizing the boat type and analyzing noise of submarine that induced by flow.And resistance contrast verified the reliability of numerical simulation to some degree.

submarine;RNG k－ε model;ICEM CFD;frictional resistance

U674.76;TB53

A

1672－7649(2012)03－0019－07

10.3404/j.issn.1672－7649.2012.03.004

2011－05－10;

2011－06－28

涂海文(1986－)，男，碩士研究生，主要從事船舶與海洋工程研究。