一種基于低速、小振幅周期運(yùn)動(dòng)CFD數(shù)值模擬的水下滑翔機(jī)附加質(zhì)量求解方法

張連洪，張宇航，楊紹瓊，馬?偉，牛文棟，侯鵬榮

一種基于低速、小振幅周期運(yùn)動(dòng)CFD數(shù)值模擬的水下滑翔機(jī)附加質(zhì)量求解方法

張連洪1, 2, 3，張宇航1，楊紹瓊1, 2, 3，馬?偉1, 2, 3，牛文棟1, 2, 3，侯鵬榮4

(1. 天津大學(xué)機(jī)構(gòu)理論與裝備設(shè)計(jì)教育部重點(diǎn)實(shí)驗(yàn)室，天津 300350；2. 天津大學(xué)青島海洋技術(shù)研究院，青島 266237；3. 青島海洋科學(xué)與技術(shù)試點(diǎn)國(guó)家實(shí)驗(yàn)室海洋觀測(cè)與探測(cè)聯(lián)合實(shí)驗(yàn)室，青島 266237；4. 大連理工大學(xué)機(jī)械工程學(xué)院，大連 116024)

附加質(zhì)量；水下滑翔機(jī)；計(jì)算流體動(dòng)力學(xué)；水池試驗(yàn)

水下滑翔機(jī)作為一種新型水下航行器，可通過浮力驅(qū)動(dòng)實(shí)現(xiàn)多剖面連續(xù)滑翔運(yùn)動(dòng)，由于其強(qiáng)大的續(xù)航能力和較高的隱蔽性，被廣泛應(yīng)用于海洋調(diào)查和軍事探測(cè)等領(lǐng)域[1]．隨著海洋科學(xué)的不斷發(fā)展，人們對(duì)水下滑翔機(jī)的滑翔性能以及操縱性能要求越來越高．而附加質(zhì)量是水下滑翔機(jī)水動(dòng)力重要影響因素之一，因此，人們對(duì)附加質(zhì)量計(jì)算精度提出了更高的要求．

Lee等[2]分別采用垂直平面運(yùn)動(dòng)機(jī)構(gòu)(vertical plan motion mechanism，VPMM)測(cè)試與計(jì)算流體動(dòng)力學(xué)(computational fluid dynamics，CFD)方法對(duì)水下無人航行器的附加質(zhì)量進(jìn)行求解，驗(yàn)證了VPMM設(shè)備在水下航行器水池試驗(yàn)中的適用性；周景軍等[3]基于相對(duì)運(yùn)動(dòng)法，利用CFD分別模擬水下航行體進(jìn)行非定常運(yùn)動(dòng)與定常運(yùn)動(dòng)，得到流體動(dòng)力的差值來求解附加質(zhì)量；楊磊等[4]采用基于勢(shì)流理論的H-S面元法和基于PMM試驗(yàn)的CFD技術(shù)求得水下滑翔機(jī)附加質(zhì)量，并將結(jié)果進(jìn)行對(duì)比，相對(duì)小于11%；莫慧黠?等[5]基于動(dòng)參考系思想，避免了相對(duì)運(yùn)動(dòng)法添加動(dòng)量源項(xiàng)以及CFD模擬時(shí)動(dòng)網(wǎng)格等問題，附加質(zhì)量求解誤差在10%以內(nèi)；宋麗華[6]基于有限體積法，求解基本速度勢(shì)的拉普拉斯方程計(jì)算附加值量；劉智麗等[7]將物體在靜止來流中的運(yùn)動(dòng)問題轉(zhuǎn)換為物體靜止而來流時(shí)變的計(jì)算問題，避免了在CFD中采用動(dòng)網(wǎng)格技術(shù)，同時(shí)計(jì)算了在黏性流體中橢球體的附加質(zhì)量；此外，還可以采用細(xì)長(zhǎng)體假設(shè)或切片理論求解附加質(zhì)量，但只適用于簡(jiǎn)單三維物體，如細(xì)長(zhǎng)圓柱等[8-9]；對(duì)于復(fù)雜外形水下航行器可使用拖曳水池試驗(yàn)的方法求解附加質(zhì)量，但該方法存在費(fèi)用高、實(shí)驗(yàn)周期長(zhǎng)等問題．目前，最常用的求解形狀復(fù)雜的水下航行器附加質(zhì)量的方法是基于勢(shì)流理論的H-S面元法和數(shù)值模擬平面運(yùn)動(dòng)機(jī)構(gòu)PMM試驗(yàn)的CFD方法．其中，PMM試驗(yàn)法主要包括物體的純升沉運(yùn)動(dòng)、純橫蕩運(yùn)動(dòng)、純俯仰運(yùn)動(dòng)、純搖艏運(yùn)動(dòng)，其難點(diǎn)在于模擬過程涉及到流體流動(dòng)分離、航行器本體和附體之間的干擾作用等．

水下滑翔機(jī)依靠機(jī)翼提供的升力實(shí)現(xiàn)鋸齒形剖面運(yùn)動(dòng)，依靠尾舵實(shí)現(xiàn)其滑翔的滑翔穩(wěn)定性[10-11]．在進(jìn)行PMM試驗(yàn)時(shí)，機(jī)翼會(huì)產(chǎn)生較大的流動(dòng)分離現(xiàn)象．因此本文提出了一種采用CFD技術(shù)模擬水下滑翔機(jī)低速、小振幅周期運(yùn)動(dòng)的試驗(yàn)方法，以減小機(jī)翼與流體的相互作用．同時(shí)，本文采用基于勢(shì)流理論的H-S面元法編寫了附加質(zhì)量計(jì)算程序，用以驗(yàn)證CFD方法的準(zhǔn)確性．此外，通過對(duì)水下滑翔機(jī)等比例模型進(jìn)行水池試驗(yàn)，將得到的結(jié)果與CFD數(shù)值模擬的結(jié)果進(jìn)行對(duì)比，以驗(yàn)證這種方法的可行性．

1?附加質(zhì)量求解原理

1.1?研究對(duì)象

本文研究對(duì)象為天津大學(xué)自主研發(fā)的海燕-L長(zhǎng)航程水下滑翔機(jī)，如圖1所示．海燕-L主要由耐壓艙、機(jī)翼、前后導(dǎo)流罩和天線構(gòu)成．耐壓艙內(nèi)部包括浮力調(diào)節(jié)系統(tǒng)、姿態(tài)調(diào)節(jié)系統(tǒng)、控制系統(tǒng)、通信和導(dǎo)航系統(tǒng)等[13]．表1為海燕-L的主要特征參數(shù)．

圖1?海燕-L水下滑翔機(jī)

表1?海燕-L主要特征參數(shù)

Tab.1?Main characteristic parameters of Petrel-L

圖2?海燕-L體坐標(biāo)系

1.2?求解原理

采用空間固定坐標(biāo)系，平面為豎直面，軸根據(jù)右手定則垂直于豎直面方向．設(shè)定海燕-L在流體域中沿、、軸同時(shí)進(jìn)行低速度、小振幅的周期性升沉和橫蕩運(yùn)動(dòng)．圖3為水下滑翔機(jī)運(yùn)動(dòng)軌跡放大后的示意圖，實(shí)際模擬過程中振幅很小．

圖3?海燕-L運(yùn)動(dòng)軌跡示意

海燕-L運(yùn)動(dòng)方程為

則加速度為

將式(1)、式(2)代入式(3)，可得

在相同流體域中，分別使水下滑翔機(jī)繞體坐標(biāo)系的、、軸做純搖艏運(yùn)動(dòng)．設(shè)定水下滑翔機(jī)運(yùn)動(dòng)方程為

將式(6)代入式(7)，可得

2?利用CFD計(jì)算水下滑翔機(jī)附加質(zhì)量

2.1?CFD模型建立

本文采用CFD商用軟件STAR-CCM＋對(duì)水下滑翔機(jī)低速、小振幅周期運(yùn)動(dòng)進(jìn)行模擬．STAR-CCM+具有強(qiáng)大的網(wǎng)格支持能力，支持動(dòng)網(wǎng)格、滑移網(wǎng)格以及重疊網(wǎng)格等．其中，重疊網(wǎng)格綜合了滑移網(wǎng)格與動(dòng)網(wǎng)格的優(yōu)點(diǎn)，既能夠很好地描述部件運(yùn)動(dòng)，又能夠保證物體運(yùn)動(dòng)過程中的網(wǎng)格質(zhì)量，可以對(duì)任意不規(guī)則運(yùn)動(dòng)的物體進(jìn)行求解．因此本文采用重疊網(wǎng)格技術(shù)對(duì)水下滑翔機(jī)運(yùn)動(dòng)進(jìn)行模擬．

本算例需要建立兩個(gè)流體域，如圖4所示，分別對(duì)這兩個(gè)域劃分網(wǎng)格，即組件網(wǎng)格和背景網(wǎng)格．其中，組件網(wǎng)格為任意長(zhǎng)方體與其內(nèi)部水下滑翔機(jī)外表面之間區(qū)域的網(wǎng)格(該網(wǎng)格可由長(zhǎng)方體域與內(nèi)部水下滑翔機(jī)進(jìn)行布爾減運(yùn)算得到，長(zhǎng)方體尺寸略大于水下滑翔機(jī)整體尺寸)，背景網(wǎng)格為流體域網(wǎng)格，長(zhǎng)度為3(2.6m)，寬、高分別為1.5．本文將流體域劃分為切割體網(wǎng)格單元．

圖4?邊界條件

根據(jù)第1.2節(jié)所述方法，模擬水下滑翔機(jī)進(jìn)行低速小振幅運(yùn)動(dòng)，監(jiān)控水下滑翔機(jī)所受水動(dòng)力和力矩，即可得到水下滑翔機(jī)的8個(gè)附加質(zhì)量．

2.2?網(wǎng)格數(shù)量與計(jì)算步長(zhǎng)設(shè)置

圖5?網(wǎng)格數(shù)量設(shè)置

可以看出，不同網(wǎng)格數(shù)量對(duì)計(jì)算結(jié)果的影響不大，附加質(zhì)量上下波動(dòng)可能與計(jì)算步長(zhǎng)的選取有關(guān)．為使附加質(zhì)量結(jié)果達(dá)到一定精度，并節(jié)省計(jì)算資源，本文進(jìn)行仿真時(shí)采用85×104的網(wǎng)格數(shù)量．

圖6?計(jì)算步長(zhǎng)設(shè)置

可以看出當(dāng)計(jì)算步長(zhǎng)減小為0.02s，即一個(gè)周期計(jì)算200次時(shí)，計(jì)算結(jié)果已經(jīng)不再隨著計(jì)算步長(zhǎng)的減小而發(fā)生明顯的變化．為節(jié)省計(jì)算資源，提高附加質(zhì)量計(jì)算效率，本文采用計(jì)算步長(zhǎng)0.02s計(jì)算水下滑翔機(jī)附加質(zhì)量，此時(shí)能夠達(dá)到很高的精度．

2.3?速度、振幅影響分析

圖7?速度幅值與附加質(zhì)量的關(guān)系

圖8?角速度幅值與附加質(zhì)量的關(guān)系

2.4?計(jì)算結(jié)果

本文最終使用網(wǎng)格數(shù)為85×104的模型，計(jì)算步長(zhǎng)為0.02s，對(duì)水下滑翔機(jī)低速、小振幅周期運(yùn)動(dòng)進(jìn)行數(shù)值模擬，速度幅值設(shè)置為0.01m/s，角速度幅值設(shè)置為0.005rad/s，得到水下滑翔機(jī)8個(gè)附加質(zhì)量值如表2所示．

表2?基于CFD的附加質(zhì)量計(jì)算值

Tab.2?Value of additionalmassbased on CFD

3?結(jié)果驗(yàn)證及分析

3.1?基于H-S面元法的水下滑翔機(jī)附加質(zhì)量計(jì)算

H-S面元法是基于勢(shì)流理論，采用分布奇點(diǎn)法來求解物體附加質(zhì)量的，所針對(duì)的是物體在沒有黏性、不可壓縮的無旋流體中的運(yùn)動(dòng)問題[13]．

假設(shè)物體在靜止的理想流體中以速度作非定常直線運(yùn)動(dòng)，則速度勢(shì)滿足控制方程為

邊界條件為

由勢(shì)流理論可知，流場(chǎng)內(nèi)任意點(diǎn)的速度勢(shì)為

將式(13)兩端取法相導(dǎo)數(shù)，帶入式(12)可以得到

則物體的附加質(zhì)量為

本文基于上述方法使用Matlab編寫了水下滑翔機(jī)附加質(zhì)量計(jì)算程序．與CFD數(shù)值模擬法不同的是，CFD方法是對(duì)流體域進(jìn)行體網(wǎng)格劃分，而H-S面元法是對(duì)水下滑翔機(jī)進(jìn)行表面網(wǎng)格劃分．考慮到劃分網(wǎng)格時(shí)邊長(zhǎng)相同的四邊形面元較三角形面元會(huì)節(jié)省大量計(jì)算資源，而三角形面元相比四邊形面元更能貼合物體表面，因此本文采用四邊形占優(yōu)來劃分網(wǎng)格．本文編寫的程序能兼容三角形和四邊形兩種形狀面元，從而達(dá)到既貼合物體表面又能節(jié)省計(jì)算資源的目的．

水下滑翔機(jī)三維模型的坐標(biāo)原點(diǎn)位置應(yīng)與浮心重合．使用workbench的mesh模塊對(duì)水下滑翔機(jī)表面進(jìn)行網(wǎng)格劃分，如圖9所示，最小網(wǎng)格尺寸為15mm，表面網(wǎng)格數(shù)量為17293．

圖9?Petrel-L表面網(wǎng)格

將軟件輸出的網(wǎng)格文件導(dǎo)入程序，即可得到海燕-L水下滑翔機(jī)的附加質(zhì)量．

3.2?水池試驗(yàn)

本文針對(duì)海燕-L水下滑翔機(jī)進(jìn)行了水池試驗(yàn)，以驗(yàn)證上述CFD數(shù)值模擬方法的可行性．

試驗(yàn)地點(diǎn)為哈爾濱工程大學(xué)，試驗(yàn)設(shè)備采用水平型循環(huán)水槽，它與普通水池有所不同，被測(cè)試的模型沒有前進(jìn)速度，而周圍流場(chǎng)是持續(xù)運(yùn)動(dòng)的，因此試驗(yàn)時(shí)間不受場(chǎng)地限制，可實(shí)現(xiàn)長(zhǎng)時(shí)間的自動(dòng)化測(cè)量，極大縮減了試驗(yàn)周期和成本．

長(zhǎng)航程水下滑翔機(jī)的水池試驗(yàn)主要通過力天平測(cè)力傳感器，在循環(huán)水槽的均勻流場(chǎng)中，使用平面運(yùn)動(dòng)機(jī)構(gòu)(PMM)驅(qū)動(dòng)水下滑翔機(jī)做指定方式的運(yùn)動(dòng). 試驗(yàn)過程包括水平面定長(zhǎng)變漂角六分力測(cè)量、水平面純搖艏簡(jiǎn)諧運(yùn)動(dòng)六分力測(cè)量、水平面純橫蕩簡(jiǎn)諧運(yùn)動(dòng)六分力測(cè)量、垂直面定常變攻角六分力測(cè)量、重直面純俯仰簡(jiǎn)諧運(yùn)動(dòng)六分力測(cè)量以及垂直面純升沉簡(jiǎn)諧運(yùn)動(dòng)六分力測(cè)量．同時(shí)，采集水下滑翔機(jī)運(yùn)動(dòng)過程中所受的力和力矩，繪制并擬合受力曲線，求取水下滑翔機(jī)附加質(zhì)量．圖10為水下滑翔機(jī)模型的加工與裝配圖．圖11為試驗(yàn)過程．

圖10?海燕-L模型加工與裝配

圖11?試驗(yàn)過程

3.3?結(jié)果與分析

表3?試驗(yàn)結(jié)果

Tab.3?Experimental results

(1) 如圖11(a)所示，天平支撐桿與滑翔機(jī)殼體連接處存在較大縫隙，在進(jìn)行水池試驗(yàn)時(shí)會(huì)產(chǎn)生較大的內(nèi)流場(chǎng)；

(2) 受到航道尺度的限制，水下滑翔機(jī)純升沉簡(jiǎn)諧運(yùn)動(dòng)受到邊界的干擾較為嚴(yán)重；

(3) 水下滑翔機(jī)機(jī)翼材質(zhì)為碳纖維，在進(jìn)行純升沉簡(jiǎn)諧運(yùn)動(dòng)時(shí)機(jī)翼會(huì)發(fā)生較大的抖動(dòng)．

4?結(jié)?語

本文以海燕-L水下滑翔機(jī)為研究對(duì)象，提出了一種采用CFD數(shù)值模擬法求解水下滑翔機(jī)附加質(zhì)量的新方法，該方法只需要建模一次，在同一個(gè)流體域中模擬水下滑翔機(jī)的4種運(yùn)動(dòng)方式，監(jiān)控其所受的水動(dòng)力和力矩，即可得到海燕-L的8個(gè)附加質(zhì)量，相比文獻(xiàn)[4]中所述的基于PMM試驗(yàn)的CFD方法，極大提高了附加質(zhì)量求解效率．同時(shí)，本文采用基于勢(shì)流理論的H-S面元法，對(duì)海燕-L水下滑翔機(jī)的附加質(zhì)量進(jìn)行求解，得到兩種方法計(jì)算結(jié)果較為一致，相對(duì)誤差基本可以控制在5%以內(nèi)．最后，本文通過水池試驗(yàn)對(duì)等比例水下滑翔機(jī)模型進(jìn)行試驗(yàn)驗(yàn)證，相對(duì)誤差基本可以控制在6%以內(nèi)，可以看出本文提出的CFD方法具有較高的可行性和準(zhǔn)確性．由于水池試驗(yàn)對(duì)試驗(yàn)設(shè)備的要求較高，且費(fèi)用高、周期長(zhǎng)，因此，該CFD數(shù)值模擬方法可以在一定程度上代替水池試驗(yàn)對(duì)水下滑翔機(jī)附加質(zhì)量進(jìn)行求解，且該方法能適用于任意尺寸、形狀的水下滑翔機(jī)，對(duì)不同型號(hào)水下滑翔機(jī)水動(dòng)力外形設(shè)計(jì)、動(dòng)力學(xué)建模都具有較大的參考價(jià)值．

［1］ 何柏巖，杜金輝，楊紹瓊，等. 基于VMD-LSSVM的水下滑翔機(jī)深平均流預(yù)測(cè)[J]. 天津大學(xué)學(xué)報(bào)(自然科學(xué)與工程技術(shù)版)，2021，54(4)：388-396.

He Baiyan，Du Jinhui，Yang Shaoqiong，et al. Prediction of underwater glider depth-averaged current velocity based on VMD-LSSVM[J]. Journal of Tianjin University(Science and Technology)，2021，54(4)：388-396(in Chinese).

［2］ Lee S K，Joung T H，Cheon S J，et al. Evaluation of the added mass for a spheroid-type unmanned underwater vehicle by vertical planar motion mechanism test[J]. International Journal of Naval Architecture and Ocean Engineering，2011，3(3)：174-180.

［3］ 周景軍，李育英，項(xiàng)慶睿. 一種水下航行體附加質(zhì)量數(shù)值計(jì)算方法[J]. 魚雷技術(shù)，2013，21(4)：246-249.

Zhou Jingjun，Li Yuying，Xiang Qingrui. A numerical computation method of additional mass for underwater vehicle[J]. Torpedo Technology，2013，21(4)：246-249(in Chinese).

［4］ 楊?磊，曹軍軍，姚寶恒，等. 水下滑翔機(jī)附加質(zhì)量數(shù)值計(jì)算[J]. 艦船科學(xué)技術(shù)，2016，38(23)：116-120.

Yang Lei，Cao Junjun，Yao Baoheng，et al. Numerical calculation of additional mass of underwater glider[J]. Ship Science and Technology，2016，38(23)：116-120(in Chinese).

［5］ 莫慧黠，黨建軍，羅?凱，等. 水下航行器附加質(zhì)量數(shù)值計(jì)算方法[J]. 水下無人系統(tǒng)學(xué)報(bào)，2017，25(4)：250-255.

Mo Huixia，Dang Jianjun，Luo Kai，et al. Numerical calculation method of additional mass of underwater vehicle[J]. Journal of Underwater Unmanned Systems，2017，25(4)：250-255(in Chinese).

［6］ 宋麗華. 水下航行體附加質(zhì)量的計(jì)算方法[J]. 水雷戰(zhàn)與艦船防護(hù)，2017，25(4)：58-61，54.

Song Lihua. Calculation method of additional mass of underwater vehicle[J]. Mine Warfare and Ship Protection，2017，25(4)：58-61，54(in Chinese).

［7］ 劉智麗，林獻(xiàn)武. 有黏流場(chǎng)中旋成橢球體的軸向附加質(zhì)量系數(shù)計(jì)算[J]. 廈門大學(xué)學(xué)報(bào)(自然科學(xué)版)，2019，58(3)：408-414.

Liu Zhili，Lin Xianwu. Calculation of axial additional mass coefficient of rotating ellipsoid in viscous flow field[J]. Journal of Xiamen University(Natural Science Edition)，2019，58(3)：408-414(in Chinese).

［8］ Wang K，Zhang X，Zhang Z Q，et al. Numercal analysis of added mass and damping of floating production，storage and offloading system[J]. Acta Mechanica Sinica，2012，28(3)：870-876.

［9］ Jensen R L，Watts P B，Lawrence J E，et al. Vertical hand force and forearm emg during a high-step rock-on climbing move with and without added mass[C]// Proceedings of 23rd International Symposium on Biomechanics in Sports. Beijing，China，2005：466-469.

［10］沈新蕊，王延輝，楊紹瓊，等. 水下滑翔機(jī)技術(shù)發(fā)展現(xiàn)狀與展望[J]. 水下無人系統(tǒng)學(xué)報(bào)，2018，26(2)：89-106.

Shen Xinrui，Wang Yanhui，Yang Shaoqiong，et al. Development of underwater gliders：An overview and prospect [J]. Journal of Unmanned Undersea Systems，2018，26(2)：89-106(in Chinese).

［11］張?亮，李云波. 流體力學(xué)[M]. 哈爾濱：哈爾濱工程大學(xué)出版社，2001.

Zhang Liang，Li Yunbo. Fluid Mechanics[M]. Harbin：Harbin Engineering University Press，2001(in Chinese).

［12］ Yang Ming，Wang Yanhui，Wang Shuxin，et al. Motion parameter optimization for gliding strategy analysis of underwater gliders[J]. Ocean Engineering，2019，191：106502-1-106502-12.

［13］戴遺山，段文洋. 船舶在波浪中運(yùn)動(dòng)的勢(shì)流理論[M]. 北京：國(guó)防工業(yè)出版社，2008.

Dai Yishan，Duan Wenyang. The Potential Flow Theory of Ship Moving in Waves[M]. Beijing：National Defense Industry Press，2008(in Chinese).

Method for Solving Added Mass of Underwater Glider Based on CFD of Low-Speed and Small-Amplitude Periodic Motion

Zhang Lianhong1, 2, 3，Zhang Yuhang1，Yang Shaoqiong1, 2, 3，Ma Wei1, 2, 3，Niu Wendong1, 2, 3，Hou Pengrong4

(1. Key Laboratory of Mechanism Theory and Equipment Design of Ministry of Education，Tianjin University，Tianjin 300350，China；2. Qingdao Institute for Ocean Engineering of Tianjin University，Qingdao 266237，China；3. The Joint Laboratory of Ocean Observing and Detection，Pilot National Laboratory for Marine Science and Technology(Qingdao)，Qingdao 266237，China；4. School of Mechanical Engineering，Dalian University of Technology，Dalian 116024，China)

added mass；underwater glider；computational fluid dynamics(CFD)；pool test

U674.941

A

0493-2137(2022)03-0239-08

10.11784/tdxbz202010039

2020-10-21；

2020-11-07.

張連洪（1963—??），男，博士，教授，zhanglh@tju.edu.cn.

馬?偉，wei.ma@tju.edu.cn.

國(guó)家重點(diǎn)研發(fā)計(jì)劃資助項(xiàng)目(2019YFC0311803，2019YFC0311701)；國(guó)家自然科學(xué)基金資助項(xiàng)目(11902219)；天津市自然科學(xué)基金資助項(xiàng)目(18JCJQJC46400).

Supported by the National Key Research and Development Program of China(No.2019YFC0311803，No.2019YFC0311701)，the National Natural Science Foundation of China(No.11902219)，the Natural Science Foundation of Tianjin，China (No.18JCJQJC46400).

(責(zé)任編輯：王曉燕)