田曉潔，劉運祥,劉貴杰，高志明，焦字鵬

(中國海洋大學(xué)工程學(xué)院,山東 青島 266100)

隨著計算機技術(shù)大力發(fā)展、相關(guān)CFD(Computational Fluid Dynamics)技術(shù)的逐漸成熟，出現(xiàn)了使用數(shù)值模擬對水下泳動生物運動進(jìn)行研究分析的方法，研究者使用CFD來仿真魚體周圍的流場和壓力場，進(jìn)而計算出力場，涌現(xiàn)出許多成果[9]。WOLFGANG等[10-11]結(jié)合活魚的實驗數(shù)據(jù)，運用三維面元法分析金槍魚的游動特征，并研究尾跡中的漩渦干擾和控制。蘇玉民和黃勝[12]用三維面元法計算分析了仿生金槍魚剛性及柔性尾鰭的非定常水動力性能。楊亮[13]利用FLUENT軟件的二次開發(fā)，采用動網(wǎng)格技術(shù)詳細(xì)計算了具有金槍魚尾鰭運動規(guī)律的機翼在粘性流場中的水動力性能，計算分析了三維尾鰭的非定常水動力性能。馬堯和趙宏偉對金槍魚在海洋環(huán)境中的游動進(jìn)行仿真，著重分析了機器魚表面壓強，湍流動能以及周圍流體流速等性能指標(biāo)[14]。周澍欣,胡慶松使用fluent將魚體和魚鰭整體作為研究對象在流場中進(jìn)行了受力分析[15]。

金槍魚具有典型流線型外形，其體形特點適合傳感器等的布放。本文以金槍魚魚體特征為對象進(jìn)行分析，將其體形特征進(jìn)行數(shù)學(xué)化描述，將其用于仿生機器魚體的三維建模和網(wǎng)格化，最終在 ANSYS(Fluent模塊)中進(jìn)行來流為定常流時魚體在給定流場下所呈現(xiàn)的力學(xué)特征計算，根據(jù)結(jié)果對比分析流線型魚體在減阻方面的功能，為仿生機器魚的優(yōu)化設(shè)計提供支撐。

1 仿生金槍魚模型

魚的種類繁多，外形更是多種多樣，但一般都可分為頭部、軀干和尾部三部分。在眾多的魚類中，金槍魚魚體的長寬比、軀體厚度、尾部特征等具有較好的比例，適宜進(jìn)行內(nèi)部結(jié)構(gòu)設(shè)計以放置控制系統(tǒng)、驅(qū)動系統(tǒng)和傳感器元件等功能裝置，故以其為仿生對象進(jìn)行分析研究。本文根據(jù)金槍魚的體長、體寬、體高、尾鰭等特征，按文獻(xiàn)[16]的方法建立仿生金槍魚的三維模型，如圖1所示。尾鰭形狀比率為0.14，建立的模型為在自然狀態(tài)下處于直線型的金槍魚。

圖1 仿生機器魚三維模型Fig.1 The 3D model of bionic robot fish

2 數(shù)值計算

2.1 魚體網(wǎng)格的劃分

為便于進(jìn)行數(shù)值分析計算，合理的網(wǎng)格劃分非常重要。通常流體計算域網(wǎng)格劃分一般可將體網(wǎng)格分為結(jié)構(gòu)網(wǎng)格和非結(jié)構(gòu)網(wǎng)格，相對于結(jié)構(gòu)網(wǎng)格，非結(jié)構(gòu)網(wǎng)格雖然生成過程比較復(fù)雜，但有極好的適用性，劃分效率高，方便劃分復(fù)雜形體模型，并且能夠很好地處理復(fù)雜邊界的流場計算問題等優(yōu)點，適用模型中尾鰭等尖端部分的網(wǎng)格劃分，故選取非結(jié)構(gòu)網(wǎng)格主體的劃分方法。

魚體長度l=13 cm、最大寬度b=5 cm、最大高度h=6.8 cm，為了避免外流場邊界引起邊界效應(yīng)，設(shè)定計算域的長度L=50 cm、寬度B=20 cm、高度H=20 cm。并設(shè)定魚體在計算域的中心位置，其中魚頭朝向X軸的正方向，魚尾鰭相對于魚體處在Y軸的正方向。

2.2 湍流模型

CFD中的湍流模型的引入是解決湍流運動的方程組采用了某種平均而不封閉的問題。數(shù)值模擬的關(guān)鍵是湍流模型的選擇，在眾多的湍流模型中，當(dāng)前最常用的是基于“雷諾平均N-S(RANS)方程組”的k-ε模型。而k-ε模型又可分為標(biāo)準(zhǔn)的k-ε模型和RNG k-ε模型兩種。其中標(biāo)準(zhǔn)的k-ε模型是一種高雷諾數(shù)模型，RNG k-ε模型在標(biāo)準(zhǔn)的k-ε模型的基礎(chǔ)上增加了一個考慮低雷諾數(shù)流動的解析公式并進(jìn)行了相應(yīng)的改進(jìn)，使其在模擬湍流流動中具有更高的精度。本文仿真計算中，由于進(jìn)口的來流速度較小，仿生機器魚特征尺寸較小，因此屬于低雷諾數(shù)范圍，采用RNG k-ε模型能夠?qū)⑽⒎终承阅Ｐ图尤肫渲?，從而得到更接近理論計算值?/p>

RNG k-ε模型的具體形式為：

(1)

(2)

(3)

式中:μt為湍動粘度，ε為湍流耗散率，k為湍流動能，ρ為流體密度，Gk是由于平均速度梯度引起的湍動能k的產(chǎn)生項，Gb是由于浮力引起的湍動能k的產(chǎn)生項，C1ε、C2ε為經(jīng)驗系數(shù)。計算時各經(jīng)驗常數(shù)Cμ=0.084 5，αk=αε=1.39，C1ε=1.42，C2ε=1.68[17]。

庫爾德人分布在土耳其、伊拉克、敘利亞、伊朗、亞美尼亞等國交界處，人口總數(shù)約3000萬?！皫鞝柕滤固埂蔽挥诘鼐壈鍓K的中心，又是中東地區(qū)海拔最高的山地高原，在區(qū)位和地勢上形成對上述地區(qū)的地緣壓制，可以封鎖俄羅斯經(jīng)高加索進(jìn)入中東的戰(zhàn)略通道。如果“庫爾德斯坦”獨立，被大國勢力控制，將對周邊各國尤其是土耳其構(gòu)成致命威脅。近百年來，庫爾德人一直極力追求民族獨立，但遭到土耳其、伊拉克和伊朗等國的強力打壓，尤其是伊拉克薩達(dá)姆政府更甚。

2.3 初始條件

魚類游動的環(huán)境為自然狀態(tài)下的湖泊河流，故選擇模型為不可壓縮流模型，選擇壓力基求解器(Pressure-Based)。穩(wěn)態(tài)計算，考慮重力的影響，重力加速度為9.81 m/s2；求解方程選擇RNG k-ε模型。流體域右側(cè)為速度進(jìn)口(velocity-inlet)邊界條件，左側(cè)為壓力出口(pressure-outlet)邊界條件，其余為壁面(wall)邊界條件。壓力速度耦合方式(Pressure-Velocity)采用SIMPLE算法，因模型是使用3D四面體網(wǎng)格，且研究魚體在流場中的壓力為主要部分，壓力差值采用二階迎風(fēng)格式(Second order upwind)以提高求解后壓力結(jié)果的精度(見圖2)。

在入口處設(shè)置速度為兩個魚身(25 cm/s)及50、100 cm/s的來流。設(shè)置收斂標(biāo)準(zhǔn)為，并進(jìn)行500步的迭代計算。

圖2 流體域網(wǎng)格劃分及邊界條件設(shè)置Fig.2 Set the fluid domain mesh and boundary conditions

3 數(shù)值結(jié)果分析

在本仿真計算中所關(guān)心的問題為仿生機器魚體在來流中所呈現(xiàn)的壓力分布及游動前進(jìn)時的阻力情況。

3.1 魚體在流場中的壓力分布

圖3、4、5分別是流速為25、50、100 cm/s時仿金槍魚體在流場中的壓力分布云圖。從直觀的金槍魚體壓力云圖中可以看出：魚體受壓云圖可大致分為魚頭和魚身與尾鰭三部分，當(dāng)來流速度平行于直線機器魚體長方向時，其頭部所受壓力最大，而體寬最大處則為負(fù)壓，尾鰭中間部分受到的是來自兩側(cè)的正壓，而兩個尖端受到的卻是負(fù)壓。隨著流速的增大，魚體頭部壓力隨之增大，負(fù)壓區(qū)域逐漸減小，但體寬最大處仍為負(fù)壓。綜合分析可以看出仿生機器魚的不同部位在流場中所起到的作用和功能是不同的。在減小某部分運行阻力的設(shè)計時，毫無疑問地必須結(jié)合魚體的流線型設(shè)計與魚鰭的形態(tài)進(jìn)行考慮。

圖3 25 cm/s時魚體壓力分布云圖Fig.3 Pressure distribution nephogram of 25 cm/s

圖4 50 cm/s時魚體壓力分布云圖Fig.4 Pressure distribution nephogram of 50 cm/s

圖5 100 cm/s時魚體壓力分布云圖Fig.5 Pressure distribution nephogram of 100 cm/s

3.2 魚體在流場中的流速分布

為分析仿金槍魚在流場中的特性，選取魚體所在的中心平面為特征面進(jìn)行分析。

尾跡是一種常見的流動現(xiàn)象，是由各層流體的速度差形成的，當(dāng)中心的流速比四周的流速小，在流場中的物體下游便形成尾跡。在單個流速場中，遠(yuǎn)離魚體的流場沒有受到物體的干擾，流場保持穩(wěn)定；當(dāng)靠近機器魚時，水流向兩側(cè)分離，速度增加，而尾鰭后流速較小，從而形成尾跡。而從圖6、7、8可以看出，由于緊鄰魚體后部的尾跡區(qū)的壓力一般是低于周圍環(huán)境壓力的，所以尾跡在剛開始的區(qū)域，流體有向中心匯聚的趨勢。同時可以看出，水下機器魚前進(jìn)時，對遠(yuǎn)場的影響很小，幾乎為零；對近場的影響隨著速度的增加而逐漸變小。隨著速度的增加，魚體尾部漩渦也逐漸減小，使得分離點靠近模型的尾部，尾跡區(qū)域減小。

3.3 阻力計算結(jié)果分析與對比

仿生機器魚前進(jìn)時的阻力包括壓差阻力、粘性阻力(又稱摩擦阻力)和興波阻力[18]。仿生機器魚由于在水下航行，當(dāng)下潛深度超過1/3體長時，其興波阻力接近于0，而通常情況下仿生機器魚的水下作業(yè)時的深度要遠(yuǎn)遠(yuǎn)超過其體長。故可忽略仿生機器魚產(chǎn)生的興波阻力，分析時僅考慮壓差阻力和粘性阻力。

圖6 25 cm/s時的流速分布云圖Fig.6 The velocity distribution nephogram of 25 cm/s

圖7 50 cm/s時的流速分布云圖Fig.7 The velocity distribution nephogram of 50 cm/s

圖8 100 cm/s時流速分布云圖Fig.8 The velocity distribution nephogram of 100 cm/s

(4)

工程中一般定義粘性阻力系數(shù)Cf和壓差阻力系數(shù)Cp它們的定義式如下：

(5)

則有：

(6)

其中，A為水下機器魚的濕表面積。

不同航行速度下主體的阻力是不同的，為此在不同速度下進(jìn)行了阻力分析，F(xiàn)luent計算結(jié)果如表1所示。

表1 仿生機器魚不同速度時的阻力Table 1 The drags of robot fish at different speeds

此外，本文還建立了與金槍魚模型體長、迎流面積相同的回轉(zhuǎn)體外形水下機器人模型如圖9，對該模型進(jìn)行流場受力數(shù)值計算，與仿生機器魚的外形阻力進(jìn)行對比分析。水下機器人模型的阻力計算結(jié)果如表2所示。

圖9 回轉(zhuǎn)體型水下機器人Fig.9 Rotary type of underwater robot

表2 回轉(zhuǎn)體型水下機器人不同速度時的阻力Table 2 The drags of robot at different speeds

將表1、2中的數(shù)據(jù)綜合分析，可以看出，仿生機器魚可有效的減小其在水中航行的阻力，提高推進(jìn)效率。隨著航行速度的增大，水下機器人受到的摩擦阻力和粘壓阻力都逐漸增大，粘性阻力在總阻力中占有很大份額，總阻力的變化起著決定性影響，壓差阻力只占很小一部分，隨著航速的提高，在總阻力中所占比例也越小。

此外，機器魚在來流平行于直線型魚體的情況下，對受到來流所提供的方向的阻力主要由頭部提供，仿生機器魚魚體的流線型對魚體在該情況下所受的阻力以及升力有著重要的影響，在設(shè)計仿生魚體時，在減小運行阻力方面應(yīng)主要集中在頭部以及魚體的流線型曲面進(jìn)行優(yōu)化設(shè)計。

4 結(jié)語

本文依據(jù)金槍魚體形特征參數(shù)，建立仿生金槍魚的三維仿真模型，重點描述了對流場仿真結(jié)果影響較大的流線型魚體和尾鰭部分。采用四面體結(jié)構(gòu)網(wǎng)格對模型進(jìn)行劃分，設(shè)置合理的流場邊界條件，基于RNG k-ε模型，研究其在特定流場下的力學(xué)特征及所受的阻力情況。計算結(jié)果表明，頭部、尾鰭及魚體流線型部分對仿生機器魚的流場特性影響較大。通過與回轉(zhuǎn)體型水下機器人阻力分析對比可得：仿生機器魚可有效的減小其在水中航行的阻力，提高推進(jìn)效率。因此，仿生機器魚體的流線型外形以及尾鰭結(jié)構(gòu)為水下機器人外形結(jié)構(gòu)的優(yōu)化設(shè)計提供了依據(jù)。