基于剛性網(wǎng)格法的潛艇水動(dòng)力數(shù)值計(jì)算

劉瑞杰，劉　洋，劉一函，賈　地

（1. 中國(guó)人民解放軍92857部隊(duì)，北京 100161；2. 中國(guó)船舶重工集團(tuán)公司第七一四研究所，北京 100101；3. 中國(guó)石油天然氣股份有限公司北京油氣調(diào)控中心，北京 100012）

基于剛性網(wǎng)格法的潛艇水動(dòng)力數(shù)值計(jì)算

劉瑞杰1，劉洋2，劉一函3，賈地1

（1. 中國(guó)人民解放軍92857部隊(duì)，北京 100161；2. 中國(guó)船舶重工集團(tuán)公司第七一四研究所，北京 100101；3. 中國(guó)石油天然氣股份有限公司北京油氣調(diào)控中心，北京 100012）

通過(guò)對(duì)潛艇六自由度空間運(yùn)動(dòng)方程進(jìn)行分析，得到了直航阻力運(yùn)動(dòng)、垂直面變漂角運(yùn)動(dòng)、水平面變漂角運(yùn)動(dòng)的空間運(yùn)動(dòng)方程。利用 Suboff 模型，基于剛性網(wǎng)格法，通過(guò)設(shè)定計(jì)算域的運(yùn)動(dòng)形式，實(shí)現(xiàn)了潛艇直航阻力運(yùn)動(dòng)、垂直面變漂角運(yùn)動(dòng)、水平面變漂角運(yùn)動(dòng)的數(shù)值模擬，通過(guò)對(duì)計(jì)算結(jié)果進(jìn)行擬合得到了線性水動(dòng)力系數(shù)和非線性水動(dòng)力系數(shù)，通過(guò)與相關(guān)試驗(yàn)結(jié)果進(jìn)行對(duì)比可知，計(jì)算誤差在 7% 左右。

ALE 方程；剛性網(wǎng)格法；變漂角運(yùn)動(dòng)；空間運(yùn)動(dòng)方程

0　引　言

潛艇作為海軍重要的作戰(zhàn)武器，對(duì)于維護(hù)國(guó)家海權(quán)、協(xié)同水面艦船作戰(zhàn)、進(jìn)行戰(zhàn)略威懾具有重要意義［1］。獲取潛艇水動(dòng)力系數(shù)是評(píng)價(jià)其操縱性的重要途徑。計(jì)算流體力學(xué)（CFD）作為重要的途徑，可以節(jié)省經(jīng)費(fèi)、縮短研究周期，且由于獲得的水動(dòng)力系數(shù)越來(lái)越準(zhǔn)確，故越來(lái)越受到認(rèn)可。

文獻(xiàn)［2-5］以 Suboff 模型為研究對(duì)象對(duì)潛艇定常運(yùn)動(dòng)進(jìn)行了計(jì)算，分析了網(wǎng)格和湍流模型對(duì)計(jì)算精度的影響，并獲得了部分線型和非線性水動(dòng)力系數(shù)。文獻(xiàn)［6］通過(guò)MRF方法將回轉(zhuǎn)運(yùn)動(dòng)計(jì)算轉(zhuǎn)化為定常計(jì)算，并進(jìn)行了數(shù)值模擬，分析了回轉(zhuǎn)運(yùn)動(dòng)的原理。文獻(xiàn)［7-9］通過(guò)自編程序?qū)剞D(zhuǎn)運(yùn)動(dòng)進(jìn)行了研究，分析了湍流模型、網(wǎng)格結(jié)構(gòu)、邊界層等對(duì)計(jì)算精度的影響。

以上文獻(xiàn)涉及到的數(shù)值計(jì)算均是基于固聯(lián)于艇體的運(yùn)動(dòng)坐標(biāo)系，并將非定常問(wèn)題通過(guò) UDF 功能轉(zhuǎn)化為定常計(jì)算。本文以潛艇 Suboff 模型為研究對(duì)象，基于大地坐標(biāo)系，借助剛性網(wǎng)格方法和Mesh Motion 功能，通過(guò)流體計(jì)算軟件 Fluent，對(duì)潛艇直航阻力運(yùn)動(dòng)、水平面和垂直面變漂角運(yùn)動(dòng)、懸臂回轉(zhuǎn)運(yùn)動(dòng)進(jìn)行數(shù)值計(jì)算，并將計(jì)算得到的水動(dòng)力系數(shù)與試驗(yàn)值進(jìn)行對(duì)比，驗(yàn)證方法的有效性和準(zhǔn)確性。

1　基本理論與控制方程

1.1ALE 方程

在很多運(yùn)動(dòng)物體模擬中，幾何體的運(yùn)動(dòng)可以分解為整個(gè)計(jì)算域的剛體運(yùn)動(dòng)和局部位置的相對(duì)運(yùn)動(dòng)的疊加。本文用 ALE 方程來(lái)描述這種計(jì)算域的剛性運(yùn)動(dòng)。其在直角坐標(biāo)系下的積分形式為：

式中的解向量為：

式中：Ω 為控制體，?Ω 為控制體單元的邊界；dV 為體積微元；為控制體邊界外法向單元向量；dS 為面積微元；ρ，u，v，w，e 分別為流體的密度，x，y，z軸的流體速度和單位體積流體的總內(nèi)能。

1.2剛性網(wǎng)格法

對(duì)一些特殊的運(yùn)動(dòng)形式可以將物體作為剛性物體處理，此時(shí)，計(jì)算域網(wǎng)格作為剛性網(wǎng)格與物體固連，所以不用在每個(gè)時(shí)間步重新生成新的網(wǎng)格，基于 ALE方程便可描述計(jì)算域的運(yùn)動(dòng)形式。剛性運(yùn)動(dòng)網(wǎng)格計(jì)算方法與固定計(jì)算域相比只是在通量計(jì)算時(shí)考慮了網(wǎng)格的運(yùn)動(dòng)速度。網(wǎng)格運(yùn)動(dòng)速度和網(wǎng)格在每個(gè)時(shí)刻點(diǎn)的位置由如下方法確定。

設(shè)物體的轉(zhuǎn)動(dòng)軸的方程為：

2　水動(dòng)力系數(shù)數(shù)值計(jì)算實(shí)現(xiàn)

水動(dòng)力系數(shù)數(shù)值計(jì)算通過(guò)以下數(shù)值試驗(yàn)實(shí)現(xiàn)：直航阻力數(shù)值試驗(yàn)、垂直面變漂角數(shù)值試驗(yàn)、水平面變漂角數(shù)值試驗(yàn)。通過(guò) Suboff 模型，利用 Fluent 軟件實(shí)現(xiàn)各種運(yùn)動(dòng)。

2.1數(shù)值計(jì)算的準(zhǔn)備

通過(guò)專業(yè)網(wǎng)格劃分軟件 ICEM CFD，采用全局映射式六面體網(wǎng)格離散方法建立網(wǎng)格。計(jì)算區(qū)域的大小是：距離潛艇前部有1 L 距離，尾部有2 L 距離，上下左右為1 L 距離。根據(jù) RNG k-ε 湍流模型對(duì)于近壁面網(wǎng)格處理方法的要求，Y+ 主要集中在（0，10），網(wǎng)格數(shù)目為360 萬(wàn)。

圖1　計(jì)算域網(wǎng)格截面圖Fig. 1　Mesh graph of computational domain

2.2直航阻力數(shù)值試驗(yàn)

根據(jù)縱向直航阻力試驗(yàn)的特點(diǎn)設(shè)定計(jì)算域運(yùn)動(dòng)速度 u0=0.5 m/s，1 m/s，2 m/s，2.5 m/s，3.045 m/s，5.144 m/s，6.091 m/s，7.160 m/s，8.230 m/s，9.254m/s，Turbulent Intensity和Turbulent Viscosity Ratio 均設(shè)定為1%。根據(jù)直航阻力運(yùn)動(dòng)特點(diǎn)，u=u0，v=0，w=0，p=0，q=0，r=0，則六自由度空間運(yùn)動(dòng)方程［10］可簡(jiǎn)化為：

其中 u，v，w，p，q，r 分別為計(jì)算域橫向速度、縱向速度、垂向速度、橫傾角速度、縱傾角速度及偏航角速度。通過(guò) Matlab 對(duì)計(jì)算結(jié)果進(jìn)行二次樣條線擬合如圖2所示。圖3是流場(chǎng)的速度剖面圖。

圖2　直航阻力數(shù)值試驗(yàn)擬合曲線Fig. 2　Fitted curves of direct resistance numerical experiments

圖3　流場(chǎng)縱剖面速度分布圖Fig. 3　The longitudinal velocity distribution of flow field

2.3垂直面變漂角數(shù)值試驗(yàn)

通過(guò)計(jì)算不同攻角條件下潛艇的受力和力矩，可推算相關(guān)水動(dòng)力系數(shù)。假定潛艇的速度為U=3.045 m/s，攻角 α 為± 2°，± 4°，± 6°，± 8°，± 10°，± 12°，± 15°， ± 20°，± 25°，則計(jì)算域速度計(jì)算結(jié)果如圖3所示。

圖4　垂直面變漂角斜航試驗(yàn)計(jì)算結(jié)果Fig. 4　Test results of vertical drift angle experiment

由于 Suboff 模型對(duì)稱性，故 Y0′=0，K0′=0，N0′=0。而 Mww′ 數(shù)值較小可約等于 0。垂向力Z和縱傾力矩 M 關(guān)于垂向速度 w 的函數(shù)為：

Z、M 關(guān)于不同速度 w 可得如下方程組：

解方程組并進(jìn)行無(wú)因次化，可得水動(dòng)力系數(shù)。

2.4水平面變漂角數(shù)值試驗(yàn)

假定潛艇的速度為U=3.045 m/s，漂角 β 為± 2°，± 4°，± 6°，± 8°，± 10°，± 12°，± 14°，± 16°，由垂直面變漂角 β 試驗(yàn)特點(diǎn)可知，u ≠ 0，v ≠ 0，w=0，p=0， q=0，r=0。計(jì)算域運(yùn)動(dòng)速度

由于 Suboff 模型對(duì)稱性，故 Y0′=0，N0′=0。橫向力 Y和偏航力矩 N 關(guān)于橫向速度 v 的函數(shù)為：

Y，N 關(guān)于不同速度 v 可得如下方程組：

解方程組并進(jìn)行無(wú)因次化，可得水動(dòng)力系數(shù)。

2.5計(jì)算結(jié)果校驗(yàn)

通過(guò)直航阻力數(shù)值試驗(yàn)及水平面變漂角數(shù)值試驗(yàn)獲得的水動(dòng)力導(dǎo)數(shù)計(jì)算結(jié)果誤差為3.7%，滿足計(jì)算工程需求。通過(guò)垂直面變漂角數(shù)值試驗(yàn)獲得的水動(dòng)力導(dǎo)數(shù)計(jì)算結(jié)果誤差為7%，這是由于當(dāng)攻角較大時(shí)，通過(guò)湍流模型計(jì)算不能準(zhǔn)確捕捉到流動(dòng)分離；而且全附體潛艇大攻角運(yùn)動(dòng)時(shí)，潛艇指揮臺(tái)圍殼導(dǎo)致背流區(qū)分離，而且有主附體干擾的存在，同時(shí)指揮臺(tái)圍殼對(duì)流場(chǎng)的干擾也會(huì)對(duì)翼產(chǎn)生影響，影響舵效及伴流場(chǎng)。

3　結(jié)　語(yǔ)

本文根據(jù)每種運(yùn)動(dòng)的不同，對(duì)潛艇六自由度空間運(yùn)動(dòng)方程進(jìn)行化簡(jiǎn)，獲得了直航阻力運(yùn)動(dòng)、垂直面變漂角運(yùn)動(dòng)、水平面變漂角運(yùn)動(dòng)的空間運(yùn)動(dòng)方程。利用Suboff 模型，基于剛性網(wǎng)格法，通過(guò)設(shè)定計(jì)算域的運(yùn)動(dòng)形式，實(shí)現(xiàn)了潛艇直航阻力運(yùn)動(dòng)、垂直面變漂角運(yùn)動(dòng)、水平面變漂角運(yùn)動(dòng)的數(shù)值模擬，通過(guò)對(duì)計(jì)算結(jié)果進(jìn)行擬合得到了線性水動(dòng)力系數(shù)和非線性水動(dòng)力系數(shù)，通過(guò)與相關(guān)試驗(yàn)結(jié)果進(jìn)行對(duì)比可知，計(jì)算誤差在7% 左右。

［1］肖昌潤(rùn)， 劉瑞杰， 許可， 等. 潛艇旋臂回轉(zhuǎn)試驗(yàn)數(shù)值模擬［J］. 江蘇科技大學(xué)（自然科學(xué)版）， 2014， 28（4）： 313-316. XIAO Chang-run， LIU Rui-jie， XU Ke， et al. Simulation for submarine rotating-arm tests［J］. Journal of Jiangsu University of Science and Technology （Natural Science Edition）， 2014， 28（4）：313-316.

［2］劉志華， 熊鷹， 韓寶玉. 潛艇流場(chǎng)數(shù)值計(jì)算網(wǎng)格與湍流模型選?。跩］. 華中科技大學(xué)學(xué)報(bào)（自然科學(xué)版）， 2009， 39（7）： 98-101. LIU Zhi-hua， XIONG Ying， HAN Bao-yu. Computational grid and turbulent model for calculating submarine viscous flow field［J］. Journal of Huazhong University of Science and Technology （Natural Science Edition）， 2009， 37（7）： 98-101.

［3］孫銘澤， 王永生， 張志宏， 等. 基于網(wǎng)格變形技術(shù)的全附體潛艇操縱性計(jì)算［J］. 武漢理工大學(xué)學(xué)報(bào)（交通科學(xué)與工程版），2013， 37（2）： 420-424. SUN Ming-ze， WANG Yong-sheng， ZHANG Zhi-hong， et al. Numerical simulation of submarine maneuverability based on mesh deformation technology［J］. Journal of Wuhan University of Technology （Transportation Science & Engineering）， 2013，37（2）： 420-424.

［4］柏鐵朝， 梁中剛， 周軼美， 等. 潛艇操縱性水動(dòng)力數(shù)值計(jì)算中湍流模式的比較與運(yùn)用［J］. 中國(guó)艦船研究， 2010， 5（2）： 22-28. BAI Tie-chao， LIANG Zhong-gang， ZHOU Yi-mei， et al. Comparison and application of turbulence modes in submarine maneuvering hydrodynamic forces computation［J］. Chinese Journal of Ship Research， 2010， 5（2）： 22-28.

［5］吳方良， 吳曉光， 許建， 等. 全附體潛艇三維粘性流場(chǎng)數(shù)值計(jì)算方法研究［J］. 船舶工程， 2009， 31（4）： 8-12. WU Fang-liang， WU Xiao-guang， XU Jian， et al. Study on numerical calculation method of the 3D viscous flow field over a submarine with full appendages［J］. Ship Engineering， 2009，31（4）： 8-12.

［6］劉帥， 葛彤， 趙敏. 基于源項(xiàng)法的潛艇旋臂試驗(yàn)?zāi)M［J］. 大連海事大學(xué)學(xué)報(bào)， 2011， 37（2）： 1-4. LIU Shuai， GE Tong， ZHAO Min. Simulation for submarine rotating-arm test based on added momentum source method［J］. Journal of Dalian Maritime University， 2011， 37（2）： 1-4.

［7］GREGORY P A， JOUBERT P N， CHONG M S. Flow over a body of revolution in a steady turn［R］. Rockingham： Defence Science and Technology Organisation Victoria Platform Science Lab， 2004.

［8］TOXOPEUS S， ATSAVAPRANEE P， WOLF E， et al. Collaborative CFD exercise for a submarine in a steady turn［C］//Proceedings of the OMAE ASME 31st International Conference on Ocean Offshore and Arctic Engineering. Rio de Janeiro：ASME， 2012： 761-772.

［9］MARSHALLSAY P G， ERIKSSON A M. Use of computational fluid dynamics as a tool to assess the hydrodynamic performance of a submarine［C］//Australasian Fluid Mechanics Conference. Launceston， 2012.

［10］施生達(dá). 潛艇操縱性［M］. 北京： 國(guó)防工業(yè)出版社， 1995.

［11］PAN Y C， ZHANG H X， ZHOU Q D. Numerical prediction of submarine hydrodynamic coefficients using CFD simulation［J］. Journal of Hydrodynamics， Series B， 2012，24（6）： 840-847.

Research on the electromechanical actuation system for submarine outboard

LIU Rui-jie1， LIU Yang2， LIU Yi-han3， JIA Di1
（1. No. 92857 Unit of PLA， Beijing 100161， China；2. The 714 Research Institute of CSIC， Beijing 100101， China；3. Beijing Oil and Gas Control Center， Beijing 100012， China）

Through analysis of submarine six degrees of freedom space motion equation， the direct resistance movement floated， vertical motion， horizontal angle of the space motion equation of angular motion. SUBOFF model， based on the rigid grid method， through setting calculation domain form of exercise， has realized the submarine direct resistance movement， vertical plane change drift angle motion， the horizontal drift Angle motion of the numerical simulation， based on the calculation results for fitting the hydrodynamic coefficient of linear and nonlinear hydrodynamic coefficients， by comparing with relevant experiment results， the calculation error is about 7%.

ALE equation；rigid grid；drift angle motion；the space motion equation

U675.6

A

1672-7619（2016）09-0065-04

10.3404/j.issn.1672-7619.2016.09.012

2016-06-12

劉瑞杰（1987-），男，研究方向?yàn)闈撏Р倏v性能。