李 源，王 盟，劉金勇，江興隆

（中國船舶集團(tuán)有限公司第七一〇研究所，湖北 宜昌 443003）

0 引言

建設(shè)海洋強(qiáng)國離不開海洋信息的獲取，海洋信息的獲取與傳輸離不開水下無線通信。為了完成某水下電磁波通信實(shí)驗(yàn)，采用水下實(shí)驗(yàn)平臺搭載磁性天線作為接收天線完成通訊。該水下實(shí)驗(yàn)平臺利用氣囊式浮力調(diào)節(jié)系統(tǒng)達(dá)到升沉運(yùn)動(dòng)的目的，其升沉方向流體阻力影響著其無動(dòng)力下潛與上升的運(yùn)動(dòng)過程，因此對水下實(shí)驗(yàn)平臺升沉方向流體阻力分析有著重要的意義。

當(dāng)前，國內(nèi)外使用CFD方法進(jìn)行水動(dòng)力性能研究的案例非常多。張洪彬等[1]利用AUV模型，使用 ANSYS CFX軟件在不同攻角下計(jì)算阻力系數(shù)，證明了計(jì)算結(jié)果的可靠性和CFD方法的有效性。孫鵬飛[2]以某全海深載人潛水器為研究對象，進(jìn)行直航工況下的數(shù)值模擬，發(fā)現(xiàn)研究潛器前進(jìn)方向進(jìn)行數(shù)值模擬使用Realizablek–ε湍流模型更加合適，后退方向使用Standardk–ε湍流模型。劉峰等[3]針對橢球體進(jìn)行數(shù)值模擬，求解相關(guān)的水動(dòng)力系數(shù)，發(fā)現(xiàn)CFD方法可以用于水動(dòng)力性能研究。張赫[4]針對某長航程潛水器，先后采用勢流方法、數(shù)值模擬和模型實(shí)驗(yàn)計(jì)算水動(dòng)力系數(shù)，并進(jìn)行了深入的分析和比較。孫麗[5]采用不同的網(wǎng)格模型對AUV進(jìn)行三維水動(dòng)力仿真，得出最佳網(wǎng)格模型與最佳湍流計(jì)算模型。Amiri等[6]則利用STAR-CCM+對 SUBOFF潛艇在不同水深、不同航速下直航運(yùn)動(dòng)時(shí)自由面對潛艇阻力、升力及力矩的影響進(jìn)行了分析。Toxopeus[7]從網(wǎng)格劃分、不確定度分析、敏感度分析等方面對 Suboff潛艇模型直航及偏航時(shí)的流場計(jì)算進(jìn)行了詳細(xì)討論。

本文通過CFD與經(jīng)驗(yàn)公式的計(jì)算結(jié)果對比，對不同尺寸的網(wǎng)格進(jìn)行計(jì)算分析，選出一組與經(jīng)驗(yàn)公式計(jì)算結(jié)果吻合度較高的網(wǎng)格作為后續(xù)分析的基礎(chǔ)，最終擬合出阻力系數(shù)與速度的函數(shù)，從而對實(shí)驗(yàn)平臺無動(dòng)力下潛進(jìn)行仿真，通過與海試數(shù)據(jù)對比，證明了CFD計(jì)算結(jié)果較為可靠。

1 實(shí)驗(yàn)平臺建模及網(wǎng)格劃分

1.1 求解原理

CFD數(shù)值模擬的核心思想為原來在時(shí)間域、空間域的連續(xù)場（如壓力場、速度場、密度場等）離散化，通過一定的數(shù)學(xué)方法建立關(guān)于這些離散點(diǎn)場變量關(guān)系的方程組，求解并獲得場變量的近似值[8]。一切流動(dòng)傳熱問題都必須遵從物理守恒定律，即三大守恒方程（質(zhì)量守恒、動(dòng)量守恒、能量守恒），其中質(zhì)量守恒方程又稱連續(xù)性方程，動(dòng)量守恒方程又稱N–S方程。FLUENT所求解的便有這三大方程，利用計(jì)算機(jī)強(qiáng)大的計(jì)算能力揭示了流體的運(yùn)動(dòng)特性[9]。

連續(xù)性方程：

式中：ρ是密度；t是時(shí)間；u、v、w是速度矢量u在x、y、z三方向上的分量。

N–S方程：在x、y、z三方向的方程為

式中：p是作用在微元體上的壓力；τxx、τyx、?τzx是分子粘性應(yīng)力的分量；Fx、Fy、Fz是作用在微元體上的體力。

對于不可壓縮流，如果熱交換很小可忽略時(shí)，可不考慮能量方程，故在本文中，可不考慮能量方程的影響。本文選用Standardk–ε湍流模型[10-11]進(jìn)行分析。

1.2 實(shí)驗(yàn)平臺幾何建模

為提高網(wǎng)格劃分質(zhì)量與計(jì)算收斂，本文對水下實(shí)驗(yàn)平臺進(jìn)行適當(dāng)簡化，如圖1，主要參數(shù)如表1所示。圖1中第3部分為圓柱形的氣囊保護(hù)罩，平行于氣囊保護(hù)罩母線方向?yàn)閷?shí)驗(yàn)平臺的升沉方向。后文中所研究阻力方向均為實(shí)驗(yàn)平臺的升沉方向。

圖1 水下實(shí)驗(yàn)平臺示意圖Fig.1 Schematic diagram of underwater experimental platform

表1 水下實(shí)驗(yàn)平臺主要參數(shù)Table 1 Main parameters of underwater experimental platform

1.3 網(wǎng)格尺寸對升沉阻力的分析

本文使用的前處理軟件是ICEM CFD。網(wǎng)格劃分是影響水下實(shí)驗(yàn)平臺模擬精度和計(jì)算效率的因素之一。在實(shí)際計(jì)算中，需要對網(wǎng)格進(jìn)行加密以提高精度，在對水下設(shè)備進(jìn)行CFD分析時(shí)，網(wǎng)格適當(dāng)程度的加密會(huì)提高計(jì)算精度，但網(wǎng)格并不是越細(xì)越好，由于CFD軟件分析時(shí)存在截?cái)嗾`差，任何一種網(wǎng)格劃分，過粗或過細(xì)都會(huì)容易導(dǎo)致計(jì)算偏差甚至失敗，所以對網(wǎng)格情況的研究是很有意義的，因此選取合適的網(wǎng)格是重要的環(huán)節(jié)。

考慮到本文中的水下實(shí)驗(yàn)平臺形狀不規(guī)則，本文對于外流場的劃分采用非結(jié)構(gòu)化四面體網(wǎng)格，四面體網(wǎng)格對于復(fù)雜不規(guī)則的物體有較好的適應(yīng)性，網(wǎng)格劃分如圖2所示，計(jì)算域?yàn)?.6 m*9.6 m*13.9 m的矩形體。為了比較網(wǎng)格尺寸對計(jì)算結(jié)果的影響，現(xiàn)采用5種不同的網(wǎng)格尺寸進(jìn)行計(jì)算，壁面網(wǎng)格尺寸與實(shí)驗(yàn)平臺的表面網(wǎng)格尺寸如表2所示。

圖2 計(jì)算域網(wǎng)格劃分Fig.2 Grid generation of computing domain

表2 網(wǎng)格尺寸Table 2 Grid size

以下潛工況為例，采用Standardk–ε湍流模型，來流速度v取為0.1 m/s，用以上5種不同的網(wǎng)格尺寸通過FLUENT仿真軟件對阻力性能進(jìn)行分析，并用造船界常用的摩擦阻力公式桑海公式和普朗特–許立汀公式對摩擦阻力系數(shù)進(jìn)行計(jì)算，將結(jié)果進(jìn)行對比，如表3所示。其中：

表3 Standard k–ε模型下水下實(shí)驗(yàn)平臺阻力計(jì)算表Table 3 Resistance calculation table of underwater experimental platform based on standard k–ε model

桑海公式為

普朗特–許立汀為

從表3中可以看出，在Standardk–ε模型下，采用網(wǎng)格1的劃分方式計(jì)算出的摩擦阻力系數(shù)與經(jīng)驗(yàn)公式的計(jì)算值相比偏小，網(wǎng)格2計(jì)算出的水下實(shí)驗(yàn)平臺摩擦阻力系數(shù)在經(jīng)驗(yàn)公式計(jì)算范圍之內(nèi)，而網(wǎng)格3計(jì)算出的水下實(shí)驗(yàn)平臺摩擦阻力系數(shù)略大于經(jīng)驗(yàn)公式計(jì)算值，網(wǎng)格4、5的劃分方式計(jì)算出的摩擦阻力系數(shù)與經(jīng)驗(yàn)公式的計(jì)算值相比偏小?？梢姴煌木W(wǎng)格劃分情況對實(shí)驗(yàn)平臺阻力的計(jì)算是有影響的，任何一種網(wǎng)格劃分，過粗或過細(xì)都會(huì)容易導(dǎo)致計(jì)算偏差甚至失敗，所以對網(wǎng)格情況的研究是很有意義的。Standardk–ε模型的計(jì)算中，通過比較可以得出，網(wǎng)格2是比較合適的網(wǎng)格類型，所以在本文的后續(xù)仿真中采用網(wǎng)格2進(jìn)行分析。

2 計(jì)算結(jié)果與分析

雷諾數(shù)的物理意義為流體慣性力與粘性力的無因次比，其定義為Re=ρVL/μ，其中ρ為流體密度，μ為動(dòng)力粘性系數(shù)，仿真中取20℃海水的物理性質(zhì)，ρ為 1 024.7 kg/m3，μ為 1.08×10–3Pa·s，V為流體運(yùn)動(dòng)速度，L為特征長度，取電池艙與天線艙兩翼板距離的平均值為0.34 m。

對來流速度0.01～0.19 m/s的情況利用FLUENT進(jìn)行了大量的計(jì)算，邊界條件設(shè)置如下。

流體：根據(jù)實(shí)驗(yàn)平臺的工作環(huán)境，設(shè)置流體材料的物理屬性為20℃海水的物理性質(zhì)；

入口：采用速度入口邊界條件，給定速度u=V，對于入口條件湍流定義方法的設(shè)置，選取湍流強(qiáng)度I（I=0.16Re–0.125）和湍流粘性比率（取默認(rèn)值 10）；

出口：采用 outflow自由出流邊界條件；

壁面邊界：采用無滑移邊界條件。

最終計(jì)算結(jié)果如表4所示。

表4 不同速度下水下實(shí)驗(yàn)平臺的阻力及阻力系數(shù)表Table 4 Resistance and resistance coefficient table of underwater experimental platform at different speeds

將表4中的數(shù)據(jù)整理分析，得到壓差阻力系數(shù)、摩擦阻力系數(shù)及總阻力系數(shù)隨雷諾數(shù)變化的曲線圖，如圖3所示。

圖3 不同雷諾數(shù)阻力系數(shù)變化圖Fig.3 Variation of drag coefficient at different Reynolds numbers

圖3分別給出了不同雷諾數(shù)下水下實(shí)驗(yàn)平臺的總阻力系數(shù)、壓差阻力系數(shù)和摩擦阻力系數(shù)的變化曲線圖，從圖中可以看出阻力系數(shù)為壓差阻力系數(shù)與摩擦阻力系數(shù)之和，壓差阻力系數(shù)項(xiàng)占總阻力系數(shù)的大部分，摩擦阻力系數(shù)影響很小，隨著雷諾數(shù)的變大，阻力系數(shù)都在變小，變小的速度隨雷諾數(shù)的增大而變慢，阻力隨雷諾數(shù)增大快速變大。

圖4為不同速度下阻力系數(shù)變化圖，將已計(jì)算的點(diǎn)進(jìn)行擬合，擬合成的函數(shù)為

圖4 不同速度下阻力系數(shù)變化圖Fig.4 Variation of resistance coefficient at different speeds

圖4中，R2為0.996，擬合程度較好，為后文仿真水下實(shí)驗(yàn)平臺下潛時(shí)間做準(zhǔn)備。

由于本文中水下實(shí)驗(yàn)平臺只在升沉方向上運(yùn)動(dòng)，對升沉方向上列微分方程有：

式中：f為浮力，N；Rt為實(shí)驗(yàn)平臺下潛所受阻力，N；m為總質(zhì)量，kg；g為重力加速度，m/s2；y為下潛深度，m；t為時(shí)間，s；Ct為阻力系數(shù)；A為實(shí)驗(yàn)平臺在垂直來流方向上投影面積，m2；ρ為海水的密度，kg/m3。

為了預(yù)測實(shí)驗(yàn)平臺無動(dòng)力下潛的時(shí)間，在Simulink中搭建出實(shí)驗(yàn)平臺下潛時(shí)的模型，如圖5所示。其中 MATLAB function為上文擬合出的速度阻力系數(shù)的四次函數(shù)，將下潛的深度與速度用示波器輸出。為方便表述，如圖1所示氣囊保護(hù)罩在艙體之下時(shí)沿紙面向下，為正向下潛，反之如圖8所示，氣囊保護(hù)罩在艙體之上沿紙面向下為反向下潛。同以上分析步驟，對反向下潛進(jìn)行仿真，將速度時(shí)間圖與正向下潛進(jìn)行對比，如圖6所示，正向下潛在13 s左右達(dá)到勻速狀態(tài)，速度為0.059 m/s，反向下潛在 15 s左右達(dá)到勻速狀態(tài)，速度為0.063 m/s。由于本實(shí)驗(yàn)平臺用于磁性天線通信實(shí)驗(yàn)，一方面為防止實(shí)驗(yàn)平臺下潛速度過快對磁性天線通信實(shí)驗(yàn)有所影響。另一方面，由于配重塊安裝在氣囊保護(hù)罩的底部，正向下潛有利于水下實(shí)驗(yàn)平臺重心的降低，從而提高實(shí)驗(yàn)平臺的穩(wěn)性，所以在海試中選用正向下潛的方式。

圖5 水下實(shí)驗(yàn)平臺下潛Simulink模型Fig.5 Simulink model of experimental platform during submergence

圖6 速度時(shí)間正反向下潛對比圖Fig.6 Comparison chart of speed and time during forward and reverse diving

同擬合正向下潛阻力系數(shù)與速度的函數(shù)一樣，反向下潛的阻力系數(shù)與速度的擬合函數(shù)為

與水下實(shí)驗(yàn)平臺正向下潛 Simulink搭建的模型類似，同樣的可以對水下實(shí)驗(yàn)平臺反向下潛進(jìn)行仿真。

將水下實(shí)驗(yàn)平臺正向下潛到10 m過程中的仿真與海試數(shù)據(jù)對比，如圖7所示，Simulink模型計(jì)算得出下潛到10 m水深需要171 s，海試數(shù)據(jù)可以看出下潛到10 m水深為173 s，仿真中達(dá)到勻速時(shí)為0.059 m/s，海試中為0.062 m/s，二者吻合度較高，從而也間接證明了不同速度下阻力系數(shù)計(jì)算的準(zhǔn)確度較高。圖8為實(shí)驗(yàn)平臺的實(shí)物圖，圖9為實(shí)驗(yàn)平臺的海試照片。

圖7 深度時(shí)間海試仿真對比圖Fig.7 Comparison of depth and time in sea trial simulation

圖8 實(shí)驗(yàn)平臺實(shí)物Fig.8 Physical object of experimental platform

圖9 實(shí)驗(yàn)平臺海試照片F(xiàn)ig.9 Sea trial photo of experimental platform

3 結(jié)束語

本文以水下實(shí)驗(yàn)平臺為研究對象，通過 CFD方法對實(shí)驗(yàn)平臺下潛工況進(jìn)行數(shù)值模擬，得出主要結(jié)論如下。

1）對本文中的水下實(shí)驗(yàn)平臺進(jìn)行五組網(wǎng)格的計(jì)算分析，通過與經(jīng)驗(yàn)公式計(jì)算的摩擦阻力系數(shù)對比，選出第2種網(wǎng)格為較為合適的網(wǎng)格，作為后續(xù)分析的條件。

2）通過FLUENT計(jì)算，得到該水下實(shí)驗(yàn)平臺速度從0.01～0.19 m/s的阻力與阻力系數(shù)，繪制出阻力與阻力系數(shù)在對應(yīng)雷諾數(shù)范圍內(nèi)的變化圖。

3）通過不同速度下阻力系數(shù)變化圖擬合出阻力系數(shù)與速度的函數(shù)，根據(jù)擬合出的函數(shù)與水下實(shí)驗(yàn)平臺升沉方向的微分方程搭建出Simulink模型，得到了正向下潛與反向下潛的速度時(shí)間圖，為水下實(shí)驗(yàn)平臺海試的下潛方式提供了依據(jù)。最終將正向下潛的深度時(shí)間圖與海試結(jié)果進(jìn)行了對比，水下實(shí)驗(yàn)平臺無動(dòng)力下潛到10 m的仿真與海試結(jié)果的水深分別為171 s和173 s，仿真與實(shí)驗(yàn)結(jié)果吻合度較高，從而也證明了 CFD計(jì)算結(jié)果與擬合出的速度阻力系數(shù)函數(shù)較為可靠。