李 明，王奎民，譚 浩，劉 峰

（1.海軍駐錦州地區(qū)軍事代表室，遼寧錦州 121000；2.哈爾濱工程大學(xué) 船舶工程學(xué)院，哈爾濱 150000）

運(yùn)輸型載人潛器的阻力數(shù)值預(yù)報(bào)

李 明1，王奎民1，譚 浩2，劉 峰2

（1.海軍駐錦州地區(qū)軍事代表室，遼寧錦州 121000；2.哈爾濱工程大學(xué) 船舶工程學(xué)院，哈爾濱 150000）

首先對(duì) Myring型回轉(zhuǎn)體的阻力數(shù)值計(jì)算方法進(jìn)行了討論。其中著重對(duì)回轉(zhuǎn)體周?chē)刂朴虻木W(wǎng)格離散方法及湍流模型的選擇進(jìn)行了討論。將回轉(zhuǎn)體數(shù)值計(jì)算得到的結(jié)果與實(shí)驗(yàn)值進(jìn)行對(duì)比，選出適合的湍流模型，同時(shí)驗(yàn)證了計(jì)算方法的可行性。其次應(yīng)用之前提出的數(shù)值計(jì)算方法，對(duì)某型載人潛器的模型進(jìn)行阻力數(shù)值計(jì)算，得到相應(yīng)的阻力系數(shù)。最后通過(guò)換算得到實(shí)艇阻力，并繪制了實(shí)艇的有效功率曲線(xiàn)。研究表明，相比于k-ε模型，k-ω模型能更好地計(jì)算Myring型回轉(zhuǎn)體的阻力，并且在高雷諾數(shù)下有較高的精度，更適宜目標(biāo)潛器的阻力預(yù)報(bào)。從目標(biāo)潛器模型的直航阻力入手，對(duì)實(shí)艇的阻力進(jìn)行預(yù)報(bào)，最后得到其有效功率，為該載人潛器的動(dòng)力、控制等系統(tǒng)的設(shè)計(jì)提供了參考。

載人潛器；回轉(zhuǎn)體艇型；數(shù)值計(jì)算；阻力

0 引言

水動(dòng)力是載人潛器的重要指標(biāo)，對(duì)于其總體性能有著重要的影響。然而，載人潛器的外形復(fù)雜、附體眾多，這些均導(dǎo)致其水動(dòng)力預(yù)報(bào)復(fù)雜且困難。模型實(shí)驗(yàn)法是當(dāng)前解決水動(dòng)力實(shí)際問(wèn)題較為有效的方法，但成本極高。隨著計(jì)算機(jī)技術(shù)的發(fā)展，數(shù)值模擬以其周期短、成本低等優(yōu)點(diǎn)，逐漸成為研究潛器水動(dòng)力性能的有效手段。目前，以“運(yùn)動(dòng)方程為基礎(chǔ)，拘束模型試驗(yàn)為輔”[1]，從而獲得水動(dòng)力系數(shù)的方法已成為潛器水動(dòng)力預(yù)報(bào)的規(guī)范模式，并且得到了大多數(shù)研究人員的認(rèn)可與應(yīng)用。

李佳[2]研究了載人潛器的阻力數(shù)值計(jì)算方法，并與模型艇拖曳實(shí)驗(yàn)進(jìn)行了對(duì)比，證明了數(shù)值模擬的準(zhǔn)確性。龐永杰等[3]提出了新的多面體網(wǎng)格運(yùn)動(dòng)求解器，并模擬了橢球體的平面運(yùn)動(dòng)。李剛[4]對(duì)穿梭艇進(jìn)行了一系列的數(shù)值計(jì)算，最終提出了適用的平方項(xiàng)水動(dòng)力模型。Fuglestad等[5]針對(duì)HUGIN 3000原型的阻力進(jìn)行數(shù)值計(jì)算，通過(guò)計(jì)算值與模型試驗(yàn)數(shù)據(jù)的比對(duì)，驗(yàn)證了數(shù)值計(jì)算的有效性和高雷諾數(shù)阻力預(yù)報(bào)的可能性。涂海文等[6]運(yùn)用雷諾平均 N-S方程，采用RNGk-ε湍流模型，對(duì)四種潛艇模型的阻力進(jìn)行數(shù)值模擬。趙金鑫[7]針對(duì)某型帶吊艙潛器進(jìn)行數(shù)值計(jì)算，通過(guò)比對(duì)試驗(yàn)結(jié)果驗(yàn)證了計(jì)算方法的有效性，并建立了潛器的操縱運(yùn)動(dòng)仿真系統(tǒng)。邱遼原[8]對(duì)潛艇粘性繞流場(chǎng)進(jìn)行了數(shù)值模擬，同時(shí)檢驗(yàn)了數(shù)值計(jì)算結(jié)果，創(chuàng)建了基于該數(shù)學(xué)模型的型線(xiàn)自動(dòng)生成系統(tǒng)。

本文首先以Myring型回轉(zhuǎn)體為基礎(chǔ)，針對(duì)其數(shù)值計(jì)算方法進(jìn)行了研究與討論，并將計(jì)算得到的阻力系數(shù)與實(shí)驗(yàn)值進(jìn)行對(duì)比，以檢驗(yàn)數(shù)值計(jì)算的精度和驗(yàn)證計(jì)算方法的可行性。之后根據(jù)回轉(zhuǎn)體的數(shù)值計(jì)算方法對(duì)一型載人潛器的阻力進(jìn)行數(shù)值模擬，得到潛器的阻力系數(shù)。最后，繪制實(shí)艇有效功率曲線(xiàn)。

1 計(jì)算方法的驗(yàn)證

1.1 Myring型回轉(zhuǎn)體

Myring型曲線(xiàn)是常用的AUV艏艉形狀，具體形狀見(jiàn)圖1。

本文對(duì)此艇型進(jìn)行數(shù)值模擬，并將計(jì)算值與實(shí)驗(yàn)值對(duì)比，以證明計(jì)算方法的可行性。

Myring型曲線(xiàn)的艏部形狀方程為：

艉部形狀方程為：

式中，a為艏部長(zhǎng)度，b為平行中體長(zhǎng)度，c為艉部長(zhǎng)度，d為中體直徑，x是長(zhǎng)軸上點(diǎn)到艏部頂點(diǎn)的距離，r為x點(diǎn)處的半徑，n和θ分別是控制艏艉曲線(xiàn)飽和程度的參數(shù)。

本文Myring型回轉(zhuǎn)體模型參數(shù)為，平行中體直徑d=2 800 mm，長(zhǎng)度b=724 mm。艏部參數(shù)分別為，a=280 mm，n=2；兩組艉部參數(shù)分別為，c=784 mm，θ=27°。

1.2 Myring型回轉(zhuǎn)體模型的數(shù)值計(jì)算

1.2.1 計(jì)算域的建立及網(wǎng)格劃分

由于計(jì)算模型左右對(duì)稱(chēng)，為減少網(wǎng)格數(shù)量節(jié)約計(jì)算成本，本次數(shù)值計(jì)算僅針對(duì)半個(gè)回轉(zhuǎn)體，采用笛卡爾直角坐標(biāo)系，坐標(biāo)系原點(diǎn)位于回轉(zhuǎn)體頂部。計(jì)算域?yàn)榘雸A柱形，尺寸為：半徑為1.5L，來(lái)流域、去流域分別為1.5L，5L（L為回轉(zhuǎn)體長(zhǎng)度）。

回轉(zhuǎn)體表面網(wǎng)格的離散方式采用混合網(wǎng)格劃分法。在首尾曲率變化較大的區(qū)域采用非結(jié)構(gòu)化網(wǎng)格，平行中體采用結(jié)構(gòu)化網(wǎng)格。為了在保證計(jì)算精度的同時(shí)盡量降低計(jì)算成本，將流體域分為靠近回轉(zhuǎn)體的內(nèi)域和遠(yuǎn)離回轉(zhuǎn)體的外域。內(nèi)域采用細(xì)密的非結(jié)構(gòu)化網(wǎng)格對(duì)回轉(zhuǎn)體周?chē)牧黧w進(jìn)行捕捉，外域采用粗糙的結(jié)構(gòu)化網(wǎng)格，在不影響精度的同時(shí)減少網(wǎng)格數(shù)量。由于選擇了合適的網(wǎng)格離散方法和網(wǎng)格尺度，最終取得了較高的計(jì)算精度并使用了較少的計(jì)算時(shí)間。

1.2.2 數(shù)值計(jì)算參數(shù)的設(shè)定與選取

求解器采用分離式求解器，基于有限體積法來(lái)離散控制方程，求解穩(wěn)態(tài)三維RANS方程，單相流模型。邊界條件：入口條件為速度進(jìn)口，出口條件為outflow，物體壁面條件選為無(wú)滑移條件，入口邊界條件的湍流參數(shù)項(xiàng)為湍流強(qiáng)度好湍流粘度比。

湍流強(qiáng)度和湍動(dòng)粘度比的設(shè)置參考相關(guān)文獻(xiàn)[1]設(shè)為I=0.5%，μr/μ。壓力速度耦合方式選擇SIMPLEC，壓力插值選擇PRESTO。對(duì)流項(xiàng)及擴(kuò)散項(xiàng)的選擇根據(jù)文獻(xiàn)[2]“應(yīng)用二階迎風(fēng)格式其結(jié)果比一階迎風(fēng)格式更接近試驗(yàn)值，并且這種結(jié)論不隨網(wǎng)格數(shù)目的增減而變化”，采用二階迎風(fēng)格式離散。

湍流模型的選擇需通過(guò)對(duì)比試驗(yàn)進(jìn)行定量分析。在其余設(shè)置相同的情況下，分別采用Realizablek-ε湍流模型和SSTk-ω湍流模型對(duì)回轉(zhuǎn)體進(jìn)行數(shù)值模擬。之后比對(duì)實(shí)驗(yàn)值選擇適宜的湍流模型。

1.2.3 湍流模型

1）Realizablek-ε模型

標(biāo)準(zhǔn)k-ε模型對(duì)時(shí)均變率特別大的情形，有可能導(dǎo)致負(fù)的正壓力[2]。為使流動(dòng)符合湍流的物理定律，需要對(duì)正應(yīng)力進(jìn)行某種數(shù)學(xué)約束。湍動(dòng)粘度計(jì)算式中的系數(shù)Cμ不應(yīng)該是常數(shù)，而應(yīng)與應(yīng)變率聯(lián)系起來(lái)。從而提出了Realizablek-ε模型。在Realizablek-ε模型中，關(guān)于k和ε輸運(yùn)方程如下：

其中，

式中，k為湍動(dòng)能；ε為湍流耗散率；ρ為密度；t為時(shí)間；ui(i=1,2,3)對(duì)應(yīng)x、y、z三個(gè)方向；μ為動(dòng)力粘度；Gk為湍動(dòng)能k的產(chǎn)生項(xiàng)；σk和σε為Prandtl數(shù)，分別與湍動(dòng)能k和耗散率ε對(duì)應(yīng)；μt為各向異性的張量；C1和C2為經(jīng)驗(yàn)值。

2）SSTk-ω方程模型

SSTk-ω模型由Menter發(fā)展而來(lái)，以便其在廣泛的領(lǐng)域中可以?xún)?yōu)于k-ε模型，并且在近壁面自由流中k-ω模型有廣泛的應(yīng)用范圍和精度。為了達(dá)到此目的，k-ε模型被修改為k-ω模型。SSTk-ω模型和標(biāo)準(zhǔn)k-ω模型相似，但SSTk-ω模型比標(biāo)準(zhǔn)k-ε模型在廣泛的流動(dòng)領(lǐng)域有更高的精度和可信度。

湍動(dòng)能k的輸運(yùn)方程為：

湍流比耗散率ω的輸運(yùn)方程為：

其中，湍流粘度為：

Gk表示由于平均速度梯度引起的湍動(dòng)能的產(chǎn)生項(xiàng)：

Gω表示ω的產(chǎn)生項(xiàng)：

SSTk-ω湍流模型考慮了湍流剪應(yīng)力的輸運(yùn)特性，能夠準(zhǔn)確地預(yù)報(bào)由于逆壓梯度導(dǎo)致的流動(dòng)分離點(diǎn)和分離區(qū)域。

1.2.4 計(jì)算結(jié)果與實(shí)驗(yàn)值對(duì)比分析

圖5為實(shí)驗(yàn)值與不同湍流模型所測(cè)得的數(shù)值對(duì)比圖。

由對(duì)比圖可看出，當(dāng)流速較低時(shí)k-ε模型與實(shí)驗(yàn)值更為接近，最小誤差甚至達(dá)到了0.7%。雖然低速時(shí)k-ω模型比k-ε模型精度低，但是最大誤差為6.6%、最小誤差為4.8%，這種誤差工程上可以接受。

但隨著速度的增加，k-ε模型所得到的結(jié)果與實(shí)驗(yàn)值的誤差越來(lái)越大，最大可達(dá)20.03%。而k-ω模型在高速時(shí)仍能較好的模擬流場(chǎng)，最小誤差達(dá)到2.8%。綜上，k-ω模型更適于模擬這種流動(dòng)，并且在低、高速均具有較為滿(mǎn)意的精度，所以之后載人潛器的數(shù)值計(jì)算采用k-ω湍流模型。

2 目標(biāo)載人潛器阻力預(yù)報(bào)

2.1 潛器模型的網(wǎng)格劃分及設(shè)置

針對(duì)目標(biāo)艇的網(wǎng)格劃分，艇體面網(wǎng)格劃分方法采用非結(jié)構(gòu)化網(wǎng)格；體網(wǎng)格內(nèi)域采用 Tet/Hybrid，外域采用Hex/Map網(wǎng)格類(lèi)型。

數(shù)值計(jì)算的初始設(shè)置與 Myring型回轉(zhuǎn)體的設(shè)置相同。在湍流模型的選擇方面，采用之前得到的結(jié)論，使用SSTk-ω湍流模型。

2.2 潛器模型阻力預(yù)報(bào)

載人潛器模型直航阻力數(shù)值計(jì)算結(jié)果見(jiàn)表1。

表1 載人潛器模型直航阻力計(jì)算

2.3 潛器實(shí)艇有效功率預(yù)報(bào)

由于本次研究的潛器作業(yè)地點(diǎn)遠(yuǎn)離水面，所以總阻力只包含粘性阻力。因此模型與實(shí)體之間的阻力轉(zhuǎn)換采用雷諾相似。因?yàn)槟Σ磷枇ο禂?shù)是Re的函數(shù)，故相同雷諾數(shù)下潛器模型與實(shí)體的阻力系數(shù)相同。

再根據(jù)相同外形潛器的粘壓阻力系數(shù)恒定的假定可知，相同雷諾數(shù)下潛器模型與實(shí)體的阻力系數(shù)相同。根據(jù)雷諾數(shù)相等條件（Rem=Rep）可知：

λ=為縮尺比。（m為模型，p為實(shí)艇）當(dāng)ν相同時(shí)，由Re數(shù)相似準(zhǔn)則，得出

由此推斷出，若要進(jìn)行拖曳驗(yàn)證實(shí)驗(yàn)，拖車(chē)速度需要達(dá)到實(shí)艇速度的6倍。這種要求明顯不能實(shí)現(xiàn)?？紤]到該潛器主體部分為半球體與圓柱體的結(jié)合，在水下直航運(yùn)動(dòng)時(shí)。當(dāng)Re數(shù)超過(guò)臨界Re數(shù)時(shí)，其阻力系數(shù)曲線(xiàn)將趨于水平。因此，可通過(guò)此法計(jì)得到實(shí)艇阻力。

根據(jù)圖7的數(shù)據(jù)，得出本潛器模型臨界雷諾數(shù)下的阻力系數(shù)，Cd=0.020 4。因此目標(biāo)實(shí)艇在臨界雷諾數(shù)下的阻力系數(shù)

3 結(jié)論

本文通過(guò)對(duì)某運(yùn)輸型載人潛器的阻力性能進(jìn)行研究，得到的結(jié)論如下：

1）針對(duì)Myring型回轉(zhuǎn)體的數(shù)值計(jì)算，采用SSTk-ω湍流模型能更好的模擬流場(chǎng)情況，計(jì)算值與實(shí)驗(yàn)值更為接近。

2）將Myring型回轉(zhuǎn)體的計(jì)算值與實(shí)驗(yàn)值進(jìn)行對(duì)比。由于計(jì)算誤差在工程允許范圍內(nèi)，本次計(jì)算方法的可行性得以驗(yàn)證。

3）在潛器直航阻力計(jì)算方面，通過(guò)對(duì)模型艇進(jìn)行數(shù)值計(jì)算，得到了相應(yīng)的阻力系數(shù)。同時(shí)采用“臨界雷諾數(shù)法”，得到實(shí)艇的阻力曲線(xiàn)及有效功率曲線(xiàn)。

[1]王小波.輕載HOV水動(dòng)力系數(shù)數(shù)值預(yù)報(bào)[D].哈爾濱：哈爾濱工程大學(xué), 2013.

[2]WANG Xiaobo Numerical Prediction for Hydrodynamic Coefficient of Light Load HOV[D].Harbin： Harbin Engineering University, 2013.

[3]李佳.載人潛器阻力性能的數(shù)值和試驗(yàn)預(yù)報(bào)及外形優(yōu)化研究[D].哈爾濱： 哈爾濱工程大學(xué), 2009.

[4]LI Jia.Prediction with Numerical Calculation and Model Test for Drag Performance and Form Optimization for Drag Reduction of a Manned Submersible[D].Harbin Engineering University, 2009.

[5]龐永杰, 楊路春, 李宏偉, 等.潛體水動(dòng)力導(dǎo)數(shù)的CFD計(jì)算方法研究[J].哈爾濱工程大學(xué)學(xué)報(bào),2009,30(8)： 903-908.

[6]李剛.穿梭潛器水動(dòng)力特性的數(shù)值模擬和試驗(yàn)研究[D].哈爾濱： 哈爾濱工程大學(xué), 2011.

[7]Fuglestad A L., Grahl Madsen M.Computational Fluid Dynamics Applied on an Autonomous Underwater Vehicle[C]// Proceedings of the 23rd International Conference on Offshore Mechanics and Arctic Engineering.Vancouver BC Cana-da, 2004, USA：American Society of Mechanical Engineers, 2004.

[8]涂海文, 孫江龍.基于 CFD的潛艇阻力及流場(chǎng)數(shù)值計(jì)算[J].艦船科學(xué)技術(shù), 2012,34(3)： 19-25.

[9]趙金鑫.某潛器水動(dòng)力性能計(jì)算及運(yùn)動(dòng)仿真[D].哈爾濱： 哈爾濱工程大學(xué), 2011.

[10]邱遼原.潛艇粘性流場(chǎng)的數(shù)值模擬及其阻力預(yù)報(bào)的方法研究[D].武漢： 華中科技大學(xué), 2006.

[11]楊卓懿, 宋磊.回轉(zhuǎn)體艇型阻力近似計(jì)算方法與試驗(yàn)研究[J].船舶工程, 2015(6)： 18-21.

[12]劉峰, 韓端鋒, 王小波, 等.載人潛器水動(dòng)力性能研究[J].科技導(dǎo)報(bào), 2014,32(6)： 29-33.

Numerical Prediction on Resistance of Transportation Manned Submersible

LI Ming1, WANG Kuimin1, TAN Hao2, LIU Feng2
(1.Navy Military Representative Office, Liaoning Jinzhou 121000, China; 2.College of Shipbuilding Engineering,Harbin Engineering University, Harbin 150000, China)

Firstly, the numerical method of Myring type rotary body’s resistance is discussed.The method of grid discrete on control domain surrounding the axisymmetric body and the selection of turbulence model are emphatically analyzed.The calculation results are obtained and compared with the experimental data.Thus, an available turbulence model would be selected.The calculation results are compared with the experimental data.The feasibility of the calculation method is verified.Hereafter, a series of corresponding resistance coefficients can be obtained by calculating the model of a specific type manned submarine’s resistance.At last, the resistance of real submarine is got by conversion.And the effective power curve of the real submarine is drew.The results show that, comparing k-ε model, the drag of Myring type rotary body can be calculated better by k-ω model, and more accurate data can be received under the high Reynolds number.So k-ω model is more plausible for the predication of the aim submersible vehicle’s drag.Starting with the direct navigation resistance of the aim submersible vehicle model, the drag of real marine is estimated.Finally the effective power is obtained and the paper provides references for the design of the power system and control system of manned submarines.

manned submersible; rotatory submersible; numerical calculation; resistance

U661.31

A

10.14141/j.31-1981.2017.06.015

李明（1971—），男，高級(jí)工程師，研究方向：艦艇總體設(shè)計(jì)與管理。