闞瑾瑜 韓 峰 劉佳侖 秦江濤*
(武漢理工大學(xué)船海與能源動(dòng)力工程學(xué)院1) 武漢 430063)海裝沈陽局駐大連地區(qū)第一軍事代表室2) 大連 116000)(武漢理工大學(xué)智能交通系統(tǒng)研究中心3) 武漢 430063)(武漢理工大學(xué)國(guó)家水運(yùn)安全工程技術(shù)研究中心4) 武漢 430063)
隨著計(jì)算機(jī)科學(xué)的不斷發(fā)展,數(shù)值模擬在船舶阻力性能的研究方面得到了越來越廣泛的應(yīng)用[1-4].近年來許多學(xué)者結(jié)合數(shù)值模擬方法對(duì)船舶實(shí)船阻力的預(yù)報(bào)方法進(jìn)行了研究.倪崇本等[5]提出了一種基于計(jì)算流體力學(xué) (computational fluid dynamics, CFD) 的實(shí)船總阻力工程預(yù)報(bào)方法.該方法通過對(duì)疊模船體繞流場(chǎng)的數(shù)值計(jì)算得到形狀因子,再基于歐拉方程的求解獲得船體興波阻力,并在此基礎(chǔ)上對(duì)實(shí)船阻力進(jìn)行預(yù)報(bào);經(jīng)過驗(yàn)證,該方法預(yù)報(bào)精度滿足工程需要,具有較強(qiáng)的工程實(shí)用性.沈興榮等[6]基于數(shù)值模擬和模型試驗(yàn)提出了一種實(shí)船阻力三因次預(yù)報(bào)方法,采用勢(shì)流方法計(jì)算興波阻力和對(duì)應(yīng)航速下的航態(tài),根據(jù)船舶姿態(tài)采用疊模黏性流場(chǎng)求解的方法計(jì)算黏壓阻力,依據(jù)1957-ITTC公式計(jì)算摩擦阻力獲得形狀因子,從而預(yù)報(bào)實(shí)船阻力.研究中模型試驗(yàn)數(shù)據(jù)和數(shù)值模擬結(jié)果對(duì)比良好,但由于缺乏實(shí)船試驗(yàn)數(shù)據(jù),無法評(píng)估預(yù)報(bào)的精度.張恒等[7]使用CFD方法,通過調(diào)整流體介質(zhì)的黏性系數(shù)進(jìn)行不同尺度船舶的力學(xué)全相似仿真,計(jì)算結(jié)果表明該方法求解的不同尺度船模運(yùn)動(dòng)、受力特性及流場(chǎng)分布較為相似,從而在滿足全相似情況下避免了尺度效應(yīng).蔡博奧等[8]以NPL高速圓舭排水船模型為研究對(duì)象,使用CFD方法探究了其阻力成分隨著尺度的變化規(guī)律,提出了一種新的三體船阻力預(yù)報(bào)方法并對(duì)其合理性進(jìn)行驗(yàn)證.這些研究進(jìn)一步證明了CFD方法對(duì)于實(shí)船仿真的準(zhǔn)確性與實(shí)用性,但寬扁肥大駁船類似船型與傳統(tǒng)的運(yùn)輸船舶阻力特性較為不同,且關(guān)于此船型阻力性能尺度效應(yīng)的相關(guān)研究相對(duì)較少[9-10].
文中以一寬扁肥大海上自升式平臺(tái)工作船為研究對(duì)象,使用數(shù)值模擬方法考察其不同尺度在相當(dāng)速度下的阻力性能,并結(jié)合當(dāng)前常規(guī)阻力換算方法進(jìn)一步探討寬扁駁船船型的實(shí)船阻力預(yù)報(bào)方法.
研究對(duì)象海上自升式平臺(tái)工作船船體幾何模型見圖1.該船長(zhǎng)寬比較小(L/B=1.62),方形系數(shù)大(Cb=0.895),艏艉船型無明顯收縮過渡,與常規(guī)駁船船型較為相似.布置有樁靴、吊艙等非流線型附體以及三個(gè)不對(duì)稱分布的月池.
圖1 海上自升式平臺(tái)工作船船體幾何模型
通過求解RANS方程對(duì)黏性不可壓縮流體進(jìn)行數(shù)值模擬,該方程的張量形式為
(i、j=1,2,3)
(1)
與N-S方程相比,RANS方程多出的“雷諾應(yīng)力項(xiàng)”導(dǎo)致了方程組的不封閉,因此需要建立湍流模型模擬雷諾應(yīng)力.k-ε湍流模型為船舶行業(yè)普遍采用的兩方程模型,采用的Realizablek-ε湍流模型為該模型的改型,并結(jié)合壁面函數(shù)方法進(jìn)行近壁面流動(dòng)的處理.采用壁面函數(shù)方法時(shí),y+一般要求在30~300,經(jīng)過嘗試,取y+為240.
由于研究對(duì)象布置有不對(duì)稱的月池,因此在進(jìn)行數(shù)值模擬時(shí)選擇整船及其繞流場(chǎng)為對(duì)象.計(jì)算域的上邊界、下邊界、兩個(gè)側(cè)邊界及入口邊界設(shè)定為速度入口,出口邊界設(shè)定為壓力出口邊界,船體壁面設(shè)置為不可滑移壁面邊界,即壁面處法向速度為0且切向無滑移.
戰(zhàn)斗也驚動(dòng)了云浮族的族人們,他們奔到天葬場(chǎng),望著滿地的狼藉,無不驚訝得變了顏色。有人眼尖,看到了尸體脖頸上插著的竹葉鏢,交頭接耳地議論:“看到了嗎,是天使殺死的蜘蛛精……”并不時(shí)偷眼瞟著唐玉煙。
在進(jìn)行數(shù)值模擬時(shí)需要確定合適大小的計(jì)算域,以保證數(shù)值結(jié)果精度的前提下節(jié)約計(jì)算資源.計(jì)算域大小的選取經(jīng)過嘗試最終確定為船前約1.5倍船長(zhǎng),船側(cè)及船底約2.5倍船長(zhǎng),船后約3.5倍船長(zhǎng),船上方約1倍船長(zhǎng),船底約2.5倍船長(zhǎng).
選擇商業(yè)軟件STAR-CCM+自帶的網(wǎng)格處理模塊作為網(wǎng)格生成器,該軟件可以自動(dòng)化地按照設(shè)置的網(wǎng)格劃分原則生成高質(zhì)量的非結(jié)構(gòu)化網(wǎng)格.在劃分網(wǎng)格時(shí),需要對(duì)一些曲率變化較大的地方或者需要準(zhǔn)確捕捉流動(dòng)的地方進(jìn)行適當(dāng)?shù)募用?在本文的研究中,研究對(duì)象船的船體較為規(guī)則,但布置有不同種類的附體,部分附體形狀較為復(fù)雜且尺度較小,需要進(jìn)行一定的區(qū)域加密以捕捉小尺度幾何的形狀,部分附體加密的網(wǎng)格布置見圖2.
圖2 部分附體加密情況
依據(jù)ITTC規(guī)程估算數(shù)值模擬的不確定度USN,進(jìn)行數(shù)值模擬方法的有效性驗(yàn)證.參照規(guī)程,以網(wǎng)格數(shù)量為參數(shù)進(jìn)行了網(wǎng)格收斂性研究,網(wǎng)格尺度比rG取為2.8,針對(duì)粗、中、細(xì)三套網(wǎng)格采用相同的CFD方法在設(shè)計(jì)航速Fr=0.119 8的工況下進(jìn)行數(shù)值模擬.三重系列網(wǎng)格的參數(shù)及數(shù)值結(jié)果見表1.
表1 不同密度系列網(wǎng)格參數(shù)及阻力數(shù)值結(jié)果
總阻力系數(shù)的變化趨勢(shì)為單調(diào)收斂情況,具體驗(yàn)證過程見表2.
表2 總阻力系數(shù)驗(yàn)證
使用驗(yàn)證后的CFD方法對(duì)平臺(tái)工作船從試驗(yàn)?zāi)P统叨鹊綄?shí)船尺度共五個(gè)尺度工況在設(shè)計(jì)航速Fr=0.119 8下進(jìn)行了數(shù)值模擬,各個(gè)尺度工況的具體信息見表3.
表3 平臺(tái)工作船不同尺度數(shù)值模擬工況
各個(gè)阻力系數(shù)數(shù)值結(jié)果及摩擦阻力經(jīng)驗(yàn)公式結(jié)果隨尺度變化趨勢(shì)見圖3.摩擦阻力經(jīng)驗(yàn)公式見式(2).
圖3 平臺(tái)工作船不同尺度下阻力系數(shù)變化規(guī)律
(2)
由圖3可知,摩擦阻力系數(shù)的數(shù)值模擬結(jié)果與ITTC經(jīng)驗(yàn)公式結(jié)果較為吻合,兩條曲線近乎重疊.隨著尺度的不斷增大,即雷諾數(shù)的不斷增大,總阻力系數(shù)和摩擦阻力系數(shù)都呈現(xiàn)下降的趨勢(shì),并且二者的變化趨勢(shì)曲線近似平行.剩余阻力系數(shù)的變化相對(duì)較小近似為常數(shù),處于小范圍內(nèi)波動(dòng)的狀態(tài).綜上可以認(rèn)為尺度效應(yīng)主要由摩擦阻力導(dǎo)致.
不同尺度工況下數(shù)值模擬的航態(tài)結(jié)果見圖4.其中縱傾角定義尾傾為正,首傾為負(fù),升沉量依據(jù)船長(zhǎng)進(jìn)行量綱一的量化,定義上浮為正,下沉為負(fù).
圖4 平臺(tái)工作船不同尺度下航態(tài)變化規(guī)律
由圖4可知,各個(gè)尺度下升沉量的量綱一的量值是較為相近的,縱傾角的差距和吃水相比為一小量,因此可以近似認(rèn)為各個(gè)尺度下船舶的航態(tài)是保持一致的.
圖5為船艉處以LAMD2準(zhǔn)則辨識(shí)的漩渦結(jié)構(gòu).
圖5 船艉漩渦結(jié)構(gòu)及壓力分布
結(jié)合船艉的壓力系數(shù)分布可知,該漩渦結(jié)構(gòu)出現(xiàn)在船艉的低壓區(qū).正是由于該漩渦結(jié)構(gòu)影響了船艉的壓力分布,才導(dǎo)致船艉出現(xiàn)低壓區(qū).
圖6為不同尺度工況下船側(cè)剪切應(yīng)力系數(shù)分布情況,通過式(3)進(jìn)行了量綱一的量化.
圖6 平臺(tái)工作船不同尺度下船側(cè)剪切應(yīng)力系數(shù)圖(側(cè)視圖)
(3)
式中:τ為船體所受剪切力;V為航速.
由圖6可知,在尺度較小的情況下,船體的剪切應(yīng)力系數(shù)相對(duì)較高.隨著尺度的增加,船體的剪切應(yīng)力系數(shù)在不斷降低.
圖7為不同尺度時(shí)中縱剖面與船體波形圖,首柱處橫坐標(biāo)為1,尾柱處為0,波高通過式(4)進(jìn)行無量綱化.
圖7 平臺(tái)工作船不同尺度下中縱剖面與船體波形圖
(4)
式中:h為波面高度.
由圖7可知,各個(gè)尺度下的波形總體相近,縮尺比在4.33和2時(shí),尾柱處的波形與其余三個(gè)尺度有所不同,但二者之間較為接近.
基于水池試驗(yàn)相同模型尺度的數(shù)值模擬結(jié)果,使用傳統(tǒng)的二因次換算法和三因次換算法對(duì)實(shí)船阻力進(jìn)行預(yù)報(bào),并將兩者的結(jié)果與實(shí)尺度數(shù)值模擬的阻力結(jié)果進(jìn)行比較.模型尺度時(shí)不同航速的阻力數(shù)值結(jié)果見表4和圖8,其中剩余阻力部分借助疊模的數(shù)值模擬分解為黏壓阻力和興波阻力.
表4 平臺(tái)工作船不同航速下阻力數(shù)值模擬結(jié)果
圖8 平臺(tái)工作船阻力數(shù)值結(jié)果隨航速變化規(guī)律
由圖8可知,隨著航速的增加,各個(gè)阻力成分都呈現(xiàn)下降的趨勢(shì),摩擦阻力系數(shù)和粘壓阻力系數(shù)的變化趨勢(shì)較為接近,兩者的變化趨勢(shì)近乎平行;興波阻力系數(shù)隨著航速先降低后緩慢增加.
二因次換算法中剩余阻力系數(shù)采用模型尺度設(shè)計(jì)航速的數(shù)值模擬結(jié)果,摩擦阻力系數(shù)依據(jù)ITTC經(jīng)驗(yàn)公式進(jìn)行估算.三因次換算法中形狀因子取表4中的平均值2.517.三種實(shí)船阻力預(yù)報(bào)方法對(duì)不同尺度平臺(tái)工作船的阻力預(yù)報(bào)比較見表5和圖9.
表5 三種實(shí)船阻力預(yù)報(bào)方法結(jié)果比較
圖9 不同換算方法預(yù)報(bào)不同尺度總阻力系數(shù)
由圖9可知,三種阻力預(yù)報(bào)方法所得到的結(jié)果趨勢(shì)是相似的,總阻力系數(shù)均隨著尺度的增大而降低,并且逐漸趨于平穩(wěn).但定量分析來說,數(shù)值模擬方法的結(jié)果更接近于二因次換算法的預(yù)報(bào)結(jié)果,兩者的差距小與5%,而和三因次換算法的結(jié)果相差較大.
1) 在同航速時(shí),隨著尺度的不斷增大,總阻力系數(shù)與摩擦阻力系數(shù)都呈現(xiàn)逐漸減小的趨勢(shì),剩余阻力系數(shù)較為穩(wěn)定且基本為常數(shù).摩擦阻力系數(shù)與ITTC經(jīng)驗(yàn)公式計(jì)算值吻合情況良好.尺度的變化對(duì)于航態(tài)和流場(chǎng)的影響較小,不同尺度下的航態(tài)較為接近,中縱剖面與船體波形基本重合.
2) 通過比較數(shù)值模擬及兩種傳統(tǒng)的實(shí)船阻力預(yù)報(bào)方法對(duì)于不同尺度下航行阻力的預(yù)報(bào)結(jié)果,發(fā)現(xiàn)三者的變化趨勢(shì)相似,數(shù)值模擬結(jié)果與二因次換算法的結(jié)果較為接近,因此認(rèn)為對(duì)于本文所研究的平臺(tái)工作船及類似寬扁肥大船型,二因次換算法對(duì)于實(shí)船阻力的預(yù)報(bào)相對(duì)更為準(zhǔn)確.