鄒 勁 王瑞宇 孫寒冰 蔣 一
(哈爾濱工程大學(xué) 船舶工程學(xué)院 哈爾濱150001)
高速三體滑行艇是常規(guī)滑行艇、高速多體船和氣膜減阻船的組合船型,綜合以上幾種船型的優(yōu)點(diǎn),在高速滑行時(shí),片體底部與水面接觸,艇體產(chǎn)生的興波與噴濺迅速被吸入槽道內(nèi),槽道內(nèi)會(huì)形成氣相區(qū)域、液相區(qū)域、氣液二相混合區(qū)以及噴濺水層,在一定程度上減小了滑行阻力而且還提供了減震、緩沖的作用[1],但也使其運(yùn)動(dòng)機(jī)理變得較為復(fù)雜,尤其是在高速航態(tài)下的“海豚運(yùn)動(dòng)”。目前,對(duì)于常規(guī)滑行艇縱向運(yùn)動(dòng)穩(wěn)定性的理論判定方法是通過(guò)研究運(yùn)動(dòng)微分方程解的穩(wěn)定性來(lái)判斷系統(tǒng)的穩(wěn)定性[2],但是三體滑行艇底部復(fù)雜的結(jié)構(gòu)使求解其在時(shí)域內(nèi)任意時(shí)刻艇底的水流體動(dòng)升力及力矩變得極為困難。而在CFD 領(lǐng)域,由于主要解決兩相流問(wèn)題,要求在自由液面附近的網(wǎng)格具有很高的分辨率,這一特殊性使商用軟件在船舶六自由度運(yùn)動(dòng)具有較大航態(tài)變化的高性能船舶,如快艇、高速?gòu)?fù)合船方面的模擬具有局限性。然而重疊網(wǎng)格的出現(xiàn)使結(jié)構(gòu)網(wǎng)格的自動(dòng)化成為可能,也使得以上難題獲得解決[3],尤其是在非結(jié)構(gòu)重疊網(wǎng)格技術(shù)推廣到包含任意單元類型子網(wǎng)格重疊的情況之后,非結(jié)構(gòu)重疊網(wǎng)格方法可以既具有非結(jié)構(gòu)網(wǎng)格的靈活性,又具有結(jié)構(gòu)重疊網(wǎng)格的計(jì)算精度[4]。
借助商業(yè)CFD軟件STAR-CCM+,運(yùn)用結(jié)構(gòu)化背景網(wǎng)格與非結(jié)構(gòu)子域網(wǎng)格相結(jié)合的重疊網(wǎng)格方法對(duì)某三體滑行艇在高航速下的“海豚運(yùn)動(dòng)”現(xiàn)象進(jìn)行模擬,隨后比照船模試驗(yàn)結(jié)果制定不同重心位置的計(jì)算工況,借助二分法對(duì)縱向運(yùn)動(dòng)穩(wěn)定性界限曲線進(jìn)行逼近,計(jì)算結(jié)果與試驗(yàn)結(jié)果相符,為類似船型的運(yùn)動(dòng)預(yù)報(bào)以及設(shè)計(jì)優(yōu)化提供了一定參考。
計(jì)算船型為某超高速三體滑行艇,主要參數(shù)如表1所示。圖1為模型剖線示意圖,縮尺比為1∶5.625,運(yùn)用CATIA軟件進(jìn)行建模。
表1 模型參數(shù)
借助商業(yè)CFD軟件STAR-CCM+進(jìn)行仿真計(jì)算,使用改進(jìn)型的Realizablek-ε湍流模型來(lái)封閉RANS方程,該方法可較好地模擬存在流動(dòng)分離和逆壓梯度的復(fù)雜流動(dòng)問(wèn)題。自由液面處用VOF法進(jìn)行追蹤。
圖1 模型橫剖線示意圖
圖2 CATIA模型
運(yùn)動(dòng)模擬中采用自由模方法,釋放船舶的升沉和縱搖兩自由度,完全模擬船模在靜水中的直航狀態(tài)。船模在拖曳前進(jìn)時(shí),會(huì)對(duì)其周圍流場(chǎng)產(chǎn)生影響,使壓力場(chǎng)和剪應(yīng)力場(chǎng)發(fā)生變化,也使船模所受的力和力矩發(fā)生變化:
船模剛體的六自由度控制方程為:
采用重疊網(wǎng)格技術(shù),將控制域劃分為主域與子域兩部分,在計(jì)算中每迭代一次都需重新記錄交叉在一起的不同兩種網(wǎng)格的信息并讓數(shù)值在其中傳遞,雖耗時(shí)較長(zhǎng),卻保證了計(jì)算對(duì)象網(wǎng)格形態(tài)的穩(wěn)定,有利于實(shí)現(xiàn)姿態(tài)大幅變化的運(yùn)動(dòng)模擬[5-6]。
控制域尺寸為:寬為3L,入口距離船尾為3L,出口距離船尾為4L,氣相、液相區(qū)各為1.5L。背景網(wǎng)格采用正交結(jié)構(gòu)化網(wǎng)格,子域網(wǎng)格選擇切割體網(wǎng)格(Trimmer),在自由液面處、船體運(yùn)動(dòng)范圍內(nèi)進(jìn)行網(wǎng)格加密,且加密處網(wǎng)格尺寸與子域網(wǎng)格大致相當(dāng)。背景網(wǎng)格基準(zhǔn)尺寸為0.12倍船長(zhǎng),子域以及加密區(qū)域網(wǎng)格基準(zhǔn)尺寸為0.012倍船長(zhǎng)。
由于三體滑行艇航行速度極快,當(dāng)越過(guò)阻力峰之后雷諾數(shù)Re達(dá)到107,因此壁面第一層網(wǎng)格厚度會(huì)在一定程度上影響計(jì)算機(jī)結(jié)果,邊界層網(wǎng)格的最小尺度可由以下公式[7]求出:
式中:y+代表壁面無(wú)量綱距離;Δy表示壁面第一層網(wǎng)格厚度;L為船長(zhǎng);Re表示雷諾數(shù)。
以=4,Re=2.79×107為例,邊界層層數(shù)設(shè)為三層,增長(zhǎng)率取1.1,選取不同的y+并且生成四套網(wǎng)格,將平衡后得到縱傾角相對(duì)于試驗(yàn)值的誤差作為依據(jù),以此來(lái)確定網(wǎng)格方案。在表2中可以看到方案3誤差較小,因此選擇第三套方案。
表2 網(wǎng)格方案
圖3 控制域
常規(guī)滑行艇在靜水中航行時(shí),隨著速度不斷增加,會(huì)出現(xiàn)一種周期性、有界的垂直平面內(nèi)的運(yùn)動(dòng),這是一種垂蕩和縱搖的耦合運(yùn)動(dòng),稱為“海豚運(yùn)動(dòng)”[8]。其產(chǎn)生原因主要是隨著速度的增大,滑行艇尾部單位長(zhǎng)度上的負(fù)荷逐漸加大。盡管垂蕩及縱搖的小幅度變化只引起浸濕長(zhǎng)度的微小改變,但是由于速度較大,所以整個(gè)水動(dòng)壓力的變化很大,使得滑行艇難以維持原有平衡。
下頁(yè)圖4為當(dāng)xg= 0.561 5 m,速度v= 10.5 m/s時(shí),不同時(shí)刻的艇底壓力分布圖。由圖4(a)可以看到,t= 1.35 s時(shí),水動(dòng)升力集中在艇尾的一小部分區(qū)域,使浸濕長(zhǎng)度減小,整體動(dòng)升力下降,艇體開(kāi)始下沉并伴隨著埋首。當(dāng)動(dòng)浮力點(diǎn)不斷前移并達(dá)到平衡時(shí),由于慣性作用,運(yùn)動(dòng)并不停止,直到動(dòng)浮力集中在船艏附近并形成小范圍高壓區(qū),如圖4(b)所示,此時(shí)艇體又開(kāi)始抬首。這一過(guò)程循環(huán)往復(fù),就形成了周期性的垂蕩與縱搖組合——“海豚運(yùn)動(dòng)”現(xiàn)象。
圖4 艇底壓力分布圖(xg = 0.561 5 m,v = 10.5 m/s)
三體滑行艇的這種縱向不穩(wěn)定性與 和重心位置xg有關(guān)。王慶旭等[9]通過(guò)修改重心縱向位置xg、改變航速v對(duì)該三體滑行艇的船模進(jìn)行靜水直航試驗(yàn),并將得到的縱向失穩(wěn)點(diǎn)進(jìn)行曲線擬合,得到關(guān)于的曲線函數(shù),即縱向穩(wěn)定性界限曲線,其表達(dá)式為:
式中:CB=Δ/(ρV2B2/2),表示動(dòng)負(fù)荷系數(shù)。
由式(6)可見(jiàn),若單純?cè)黾又匦木嚯x尾板的距離xg,失穩(wěn)的臨界會(huì)變大,也就是失穩(wěn)速度會(huì)加大。為了驗(yàn)證數(shù)值計(jì)算結(jié)果,將數(shù)值計(jì)算工況與船模試驗(yàn)對(duì)應(yīng)(文中所有實(shí)驗(yàn)數(shù)據(jù)均參照文獻(xiàn)[9]),運(yùn)用“二分法”對(duì)四種不同重心縱向位置處的失穩(wěn)點(diǎn)進(jìn)行逼近:速度相隔1 m/s,在試驗(yàn)值附近選取三個(gè)工況分別計(jì)算,若在所選工況中不存在未失穩(wěn)點(diǎn)則將工況向下調(diào)整0.5 m/s;若包含,則對(duì)未失穩(wěn)點(diǎn)和下一個(gè)與其相鄰的工況進(jìn)行二分繼續(xù)計(jì)算,精度為0.5 m/s。初始工況如表3所示。
表3 數(shù)值計(jì)算工況
圖5為垂蕩時(shí)歷曲線,下頁(yè)圖6為縱搖時(shí)歷曲線。
圖5 垂蕩時(shí)歷曲線
在圖5(a)、6(a)中可以看到,xg= 0.561 5 m、速度在8.5 m/s時(shí),縱向運(yùn)動(dòng)雖有規(guī)律,但幅度微小,升沉幅度不足0.01 m,縱搖幅度不到0.5°;當(dāng)速度從8.5 m/s提升到9.5 m/s時(shí),升沉以及縱搖幅度明顯加劇,升沉幅度達(dá)到0.03 m,縱搖2.5°左右。不難看出,xg= 0.561 5 m時(shí),縱向運(yùn)動(dòng)失穩(wěn)點(diǎn)在8.5 m/s之前且就在附近。事實(shí)上,將速度降到v= 8.0 m/s時(shí),“海豚運(yùn)動(dòng)”消失。同樣,xg= 0.611 5 m時(shí),v= 10 m/s依然失穩(wěn);當(dāng)速度降到9.5 m/s時(shí),海豚運(yùn)動(dòng)消失。因此,xg= 0.561 5 m、xg= 0.611 5 m時(shí)的臨界失穩(wěn)速度分別在8.0~8.5 m/s、9.5~10 m/s之間。而另外兩個(gè)重心位置,在預(yù)先選取的最小速度點(diǎn)處均未出現(xiàn)“海豚運(yùn)動(dòng)”,對(duì)最小速度點(diǎn)與相鄰速度點(diǎn)進(jìn)行二分,繼續(xù)計(jì)算,分別得到未失穩(wěn)點(diǎn)v= 9.5 m/s、v= 12.5 m/s,將數(shù)值結(jié)果與實(shí)驗(yàn)結(jié)果比較,得到表5。
圖6 縱搖時(shí)歷曲線
表5 模擬結(jié)果與實(shí)驗(yàn)結(jié)果
在表5中,v1表示穩(wěn)定速度上限、v2表示失穩(wěn)速度下限。從計(jì)算結(jié)果中看到,隨著重心位置的前移,縱向失穩(wěn)速度變大,趨勢(shì)與實(shí)驗(yàn)結(jié)果相吻合。與實(shí)驗(yàn)結(jié)果相比,xg= 0.561 5 m時(shí)的v2誤差最大,為19%;xg= 0.731 5 m時(shí)的v2誤差最小,為3.7%;在其他結(jié)果中,除xg= 0.561 5 m時(shí)的v1之外,誤差均小于10%。圖7為縱向穩(wěn)定性界限。
表4 數(shù)值計(jì)算結(jié)果
圖7 縱向穩(wěn)定性界限
如圖7所示,實(shí)曲線以上區(qū)域?yàn)樵囼?yàn)得到的縱向失穩(wěn)區(qū)域(“海豚區(qū)”),以下為穩(wěn)定區(qū),而通過(guò)模擬得到的穩(wěn)定性上限曲線與失穩(wěn)下限曲線的中間地帶即為縱向運(yùn)動(dòng)穩(wěn)定性界限曲線所在位置。不難看出與實(shí)驗(yàn)值相比較,模擬值的穩(wěn)定區(qū)域較小,失穩(wěn)區(qū)域較大,隨著體積傅汝德數(shù)的增加誤差減小,模擬值曲線更加貼近實(shí)驗(yàn)值。
本文運(yùn)用重疊網(wǎng)格技術(shù)對(duì)三體滑行艇在高速航行時(shí)產(chǎn)生的“海豚運(yùn)動(dòng)”進(jìn)行CFD仿真計(jì)算,并結(jié)合二分法對(duì)縱向穩(wěn)定性界限進(jìn)行逼近計(jì)算。結(jié)果顯示:
(1)隨著重心縱向位置的前移,失穩(wěn)臨界速度隨之變大,趨勢(shì)與試驗(yàn)結(jié)果相吻合。
(2)二分法逼近得到的縱向穩(wěn)定性界限位于實(shí)驗(yàn)值下方,即失穩(wěn)區(qū)大于實(shí)驗(yàn)結(jié)果。
(3)在高傅汝德數(shù)條件下,計(jì)算結(jié)果更加貼近實(shí)驗(yàn)值,誤差極小。
綜上所述,結(jié)果顯示文中運(yùn)用的數(shù)值計(jì)算方法可以較好實(shí)現(xiàn)對(duì)三體滑行艇縱向失穩(wěn)運(yùn)動(dòng)的模擬,且失穩(wěn)速度較為準(zhǔn)確,相比于船模試驗(yàn)具有高效能、低成本等優(yōu)勢(shì),可為類似船型的設(shè)計(jì)與優(yōu)化等工作提供可信的參考數(shù)據(jù)。
[1] 孫華偉.三體滑行艇船型與阻力性能研究[D].哈爾濱:哈爾濱工程大學(xué),2010:1-51.
[2] Azcueta R,Caponnetto M,Soding H.Planing boats in waves[C]// Proceedings of IWSH’ 2003.Third International Workshop on Ship Hydrodyamics,Wuhan :Wuhan University of Technology Press,2005,2.
[3] 趙發(fā)明,高成君,夏瓊.重疊網(wǎng)格在船舶CFD中的應(yīng)用研究[J].船舶力學(xué),2011(4):332-341.
[4] 田書玲.基于非結(jié)構(gòu)網(wǎng)格方法的重疊網(wǎng)格算法研究[D].南京:南京航空航天大學(xué),2008:56-80.
[5] Thomas P D,Middlecoff J.F.,Direct Control of the Grid Point Distribution in Meshes Generated by Ellipitic Equations[J].AIAA Journal.1980,18 :652-656.
[6] Royk,Bharats,Rajkeshar S A.Comprehensive Generalized Mesh System for CFD Application[J].Mathe and Computers in Simulation,2008,78:605-617.
[7] 鄧銳.阻流板對(duì)雙體船水動(dòng)力性能影響的數(shù)值研究[D].哈爾濱:哈爾濱工程大學(xué),2010:31-38.
[8] 趙連恩.高性能船舶水動(dòng)力原理與設(shè)計(jì)[M].哈爾濱:哈爾濱工程大學(xué)出版社,2007:195-199.
[9] 王慶旭,鄒勁,史圣哲,等.高速型三體滑行艇簡(jiǎn)介及縱向穩(wěn)定性初步研究[A].第十五屆中國(guó)海洋(岸)工程學(xué)術(shù)討論會(huì)論文集,2011:226-230.