基于質(zhì)量源法的全墊升氣墊船興波模擬研究

陳熙， 朱仁傳， 顧孟瀟， 高嵩

(1.上海交通大學(xué) 海洋工程國(guó)家重點(diǎn)實(shí)驗(yàn)室，上海 200240；2 上海船舶及海洋工程設(shè)計(jì)研究院，上海 200011)

氣墊船是以靜態(tài)空氣壓力支撐的一類高性能船[1]，其中依靠全部靜態(tài)空氣壓力支撐且速度較高的類型為全墊升式氣墊船ACV。該類型氣墊船在航行中的阻力絕大部分來(lái)自于興波阻力，探究氣墊面興波的實(shí)質(zhì)及與一般排水型船舶的異同，厘清其與航速、主尺度、氣墊壓力、噴氣流量等參數(shù)之間的關(guān)系，建立求解興波阻力的方法，對(duì)于氣墊船快速性優(yōu)化、及后續(xù)圍裙動(dòng)力學(xué)和氣墊空氣動(dòng)力學(xué)[2]的研究具有重要意義。

研究氣墊船興波問(wèn)題的方法主要有3大類：模型試驗(yàn)、理論方法及數(shù)值方法[3]。Newman等[4]提出了氣墊做勻速直線運(yùn)動(dòng)的興波阻力計(jì)算方法，將氣墊簡(jiǎn)化為一個(gè)作用在水面上的壓力面，面上壓力分布均勻，并給出了相應(yīng)工況下的興波阻力系數(shù)圖譜，后續(xù)學(xué)者研究多基于此來(lái)參照。Barratt[5]針對(duì)氣墊底部不同壓力分布時(shí)的情況，提出了淺水中的興波阻力計(jì)算方法，結(jié)果表明興波阻力的峰值隨著流域的加深而提前發(fā)生，但是因其假設(shè)興起的波浪是微幅波，其結(jié)果均較大于試驗(yàn)結(jié)果。Newman等[4]提供的圖譜在中高速時(shí)與試驗(yàn)結(jié)果吻合良好，然而針對(duì)低速情況下的圖譜結(jié)果出現(xiàn)了稠密的峰谷波動(dòng)，與試驗(yàn)不符。Doctors[6]針對(duì)該問(wèn)題提出了新的氣墊底部壓力分布方式，其假定底部為雙曲正切分布，而非之前給定的均布?jí)毫?，?jì)算表明在阻力峰出現(xiàn)之前沒(méi)有出現(xiàn)稠密的峰谷波動(dòng)，彌補(bǔ)了Newman圖譜的不足。同時(shí)，Doctors采用該方法研究了氣墊船在淺水中直航的阻力性能，結(jié)果表明淺水會(huì)使興波阻力增大。Lazauskas[7]采用線性理論方法預(yù)報(bào)了氣墊船非勻速航行時(shí)的阻力性能。Mileewski[8]運(yùn)用流體仿真軟件計(jì)算分析了氣墊船航行時(shí)船身四周的興波波形。Nikseresht等[9]模擬了全墊升式氣墊船氣室內(nèi)部氣墊的空氣流動(dòng)，指出氣墊的壓力分布、航速等因素極大的影響了其興波阻力與興波波形。Bhushan[10]等運(yùn)用線性理論和非線性理論2種方法對(duì)深淺水下不同壓力面形狀的氣墊興波阻力進(jìn)行了計(jì)算分析，結(jié)果表明非線性理論更有利于計(jì)算淺水中的阻力性能。Kevin等[11]對(duì)線性方法和有限體積法這2種方法進(jìn)行了對(duì)比，指出2種方法對(duì)計(jì)算較高速段的興波阻力的準(zhǔn)確性很高，而因低速時(shí)非線性影響的增強(qiáng)，線性方法不再穩(wěn)定。劉寧[12]對(duì)全墊升式氣墊船在波浪上航行的運(yùn)動(dòng)響應(yīng)進(jìn)行了研究，考慮若干非線性因素的影響，分析了圍裙觸水的影響。

前人的研究多是給定氣墊面壓力分布模擬興波情況，而不同航速下適宜的壓力分布情況較難把握，且實(shí)際情況中氣墊船無(wú)法直接控制底部壓力，而是通過(guò)控制風(fēng)機(jī)流量來(lái)調(diào)整船體。本文介紹了求解全墊升式氣墊船興波阻力的基本假設(shè)及簡(jiǎn)化方法，建立了氣室模型，并引入質(zhì)量源模擬風(fēng)機(jī)供氣，運(yùn)用STAR-CCM+商用軟件，對(duì)興波阻力進(jìn)行了計(jì)算，并與多個(gè)文獻(xiàn)計(jì)算結(jié)果進(jìn)行對(duì)比，驗(yàn)證了該方法的有效性。此外，采用該方法研究了淺水中氣墊船的興波阻力變化情況。

1 氣墊船興波模型及數(shù)值計(jì)算

氣墊興波數(shù)值模擬在一個(gè)三維數(shù)值水池中進(jìn)行的，在出口處和遠(yuǎn)離船體的邊界處設(shè)有人工阻尼消波區(qū)。以兩相流假定為基礎(chǔ)，自由面沒(méi)有擾動(dòng)時(shí)，其上部為空氣，下部為水。

為模擬全墊升式氣墊船在正常航行中的興波情況，將其放置在一個(gè)以兩相流假定為基礎(chǔ)的三維數(shù)值水池中進(jìn)行。采用固定在船上并隨船以定常速度航行的參考坐標(biāo)系來(lái)模擬計(jì)算興波情況。因全墊升氣墊船是一個(gè)多自由度非線性復(fù)雜系統(tǒng)，其中耦合了氣墊、圍裙、水面等多方面的影響，本文僅對(duì)興波情況進(jìn)行模擬，作如下簡(jiǎn)化假定：

1)氣墊船高速航行中，其水動(dòng)阻力主要由氣墊下表面與圍裙觸水的興波阻力及圍裙?jié)癖砻娴哪Σ磷枇?部分構(gòu)成，本文假設(shè)船體正常航行，圍裙不觸水，即水動(dòng)阻力只有氣墊的興波阻力。同時(shí)將圍裙視作剛性，簡(jiǎn)化氣墊船上層建筑及船體結(jié)構(gòu)。

2)所研究的工況為氣墊船穩(wěn)定航行，不考慮船體姿態(tài)變化，且氣墊底部壓力分布達(dá)到穩(wěn)態(tài)，在氣室風(fēng)機(jī)口加入質(zhì)量源來(lái)模擬供氣，從而達(dá)到墊升船體、作用水面的效果。

3)因上述假定，氣墊船穩(wěn)定航行，船體與水沒(méi)有接觸，不存在水動(dòng)摩擦阻力，因此不考慮粘性，視為理想流體，Bhushan等[9]的研究也表明在研究氣墊船純興波問(wèn)題時(shí)，不考慮粘性的結(jié)果與RANS方法差異相當(dāng)小。

1.1 控制方程

本計(jì)算模型中的流場(chǎng)是三維不可壓縮理想流體，流場(chǎng)應(yīng)滿足連續(xù)性方程和歐拉方程：

(1)

(2)

采用流體體積函數(shù)(volume of fluid，VOF)法追蹤自由面?？刂企w內(nèi)水的體積占控制體總體積之比定義為水的體積分?jǐn)?shù)α，空氣的體積分?jǐn)?shù)即為(1-α)。該模型假定水和空氣共享速度和壓力場(chǎng)，可將其看作一種等效流體，故只要處理和單相流相同形式的質(zhì)量、動(dòng)量守恒方程，該等效流體的密度為：

ρ=αρw+(1-α)ρa(bǔ)

(3)

式中下標(biāo)w、a分別指代水和空氣。體積分?jǐn)?shù)滿足連續(xù)性方程：

(4)

方程(1)、(2)、(4)組成了氣墊興波模擬數(shù)值水池的控制方程。該方程組與一般的N-S方程相比，忽略了粘性項(xiàng)的影響，同時(shí)引入了附加質(zhì)量源項(xiàng)q。

1.2 質(zhì)量源模擬風(fēng)機(jī)系統(tǒng)

全墊升氣墊船正常航行過(guò)程中，墊升風(fēng)機(jī)向船體氣道供氣，氣體經(jīng)過(guò)圍裙進(jìn)入氣室，形成氣墊將船體全部托離水面，因墊升后船體距離水面有一段距離，氣體會(huì)通過(guò)囊指下側(cè)縫隙溢出，因此供氣需持續(xù)穩(wěn)定。圖1為氣墊船橫剖面的示意圖。

如前述假設(shè)2，為模擬出氣墊施加于水面的效果，簡(jiǎn)化氣墊船結(jié)構(gòu)，僅對(duì)氣室區(qū)域進(jìn)行建模，如圖所示為氣室橫剖面圖。采用三維質(zhì)量源方法[13]模擬風(fēng)機(jī)不斷給氣室供氣。

圖2 簡(jiǎn)化氣室橫剖面Fig.2 Cross section diagram of simplified chamber

式(1)、(2)中引入的質(zhì)量源項(xiàng)q表達(dá)式為：

(5)

式中：Ωs為質(zhì)量源區(qū)域，即圖中虛線框住的部分；Δz為質(zhì)量源區(qū)域的厚度；vair為模擬風(fēng)機(jī)吹氣口的空氣速度；q的單位為s-1。通過(guò)在該區(qū)域添加與空氣進(jìn)氣速度有關(guān)的質(zhì)量源項(xiàng)，產(chǎn)生供氣效果。

在本文的驗(yàn)證工作中，前人研究的前提多是給定壓力分布情況，因此需要確定質(zhì)量源項(xiàng)q的大小與氣墊底部壓長(zhǎng)比的關(guān)系。

Bhushan等[14]給出了由伯努利方程得到的流量公式：

(6)

式中：pc為氣墊面上實(shí)際壓強(qiáng)減去標(biāo)準(zhǔn)大氣壓的壓強(qiáng)由氣墊壓長(zhǎng)比可得氣墊相對(duì)壓強(qiáng)；Q為模擬風(fēng)機(jī)的出口流量，m3/s；Cd為泄漏系數(shù)，F(xiàn)altine建議取為0.6～1.0；Ai為簡(jiǎn)化氣室的四周邊界最低點(diǎn)與水面的間隙形成的總面積；ρa(bǔ)ir為空氣密度。質(zhì)量源項(xiàng)可通過(guò)式q=Q/VΩ求出，其中VΩ為質(zhì)量源區(qū)域的體積。

1.3 興波阻力系數(shù)計(jì)算方法

全墊升氣墊船在正常航行中，船體并不與水面接觸，難以采用排水型船體的表面壓力積分法求得阻力。本文將簡(jiǎn)化氣室和氣墊視為一個(gè)整體，通過(guò)壓力施加到水面，對(duì)氣室正下方與氣墊直接接觸水面的形狀及壓力進(jìn)行處理，求得內(nèi)水面所受力，由牛頓第三定律即可得氣墊船受到的興波阻力。

取氣室正下方的水面上一矩形微元面，其長(zhǎng)為dL，寬為dB，如圖3所示。

該微元面與水平面夾角為αw，取直角坐標(biāo)系ξOη,其原點(diǎn)置于前端點(diǎn)O，Oξ軸沿航行方向向后為正，Oη軸垂直于Oξ軸，向下為正。因壓力作用與作用面垂直，則微元面上的壓力為：

Fη=pc(x,y)dBdL

(7)

該力沿船舶方向航行的分力為：

(8)

式中：dξw為微元面的波高；dx為微元面水平面投影的x方向距離。整個(gè)氣室正下方水面上沿航行方向受到的總分力Rx可積分得到：

(9)

水面受到的力為氣墊施加于其上的力，由牛頓第三定律，氣室航行時(shí)的興波阻力Rw為Rx的反作用力為：

(10)

具體求解采用將積分式離散求和的方法，表達(dá)式為：

(11)

在得到興波阻力后，可根據(jù)紐曼-波爾的氣墊船興波阻力系數(shù)計(jì)算出興波阻力系數(shù)：

(12)

2 氣墊船的興波波形及阻力計(jì)算驗(yàn)證

影響氣墊船阻力最主要的參數(shù)為壓長(zhǎng)比，其決定興波阻力特別是越峰阻力的大小。許多學(xué)者在研究氣墊船興波問(wèn)題時(shí)，多是指定氣墊底面的壓力分布，針對(duì)不同平面形狀的氣墊展開(kāi)了廣泛研究，其中較為廣泛對(duì)照的是紐曼-波爾圖譜，該圖譜給出了指定壓長(zhǎng)比的氣墊面在不同傅汝德數(shù)下的興波阻力曲線。本文取一種壓長(zhǎng)比的氣墊船模型在靜水中航行的興波情況進(jìn)行研究，以驗(yàn)證該方法的有效性。

2.1 計(jì)算域及邊界條件

選取氣墊船氣室模型長(zhǎng)LC=3 m，寬BC=2 m，氣墊壓長(zhǎng)比PC=0.012 7ρgLC，如圖4(a)，并置于圖4(b)所示的數(shù)值水池中。計(jì)算域原點(diǎn)位于氣墊體的垂直中心線和靜水面的交點(diǎn)處，x軸指向航行方向，沿船舶航行方向?yàn)檎?，y軸指向船寬方向，以船舶右舷為正，z軸向上為正。計(jì)算域范圍為：-7LC

圖4 氣室及計(jì)算域區(qū)域Fig.4 Air chamber and calculation domain

以船舶航行方向視角來(lái)定義計(jì)算域邊界條件，前邊界面及上下面設(shè)定為速度入口，左右面設(shè)定為對(duì)稱邊界，后邊界面為壓力出口，氣室四壁及頂面為壁面。同時(shí)在左右面及后邊界面設(shè)阻尼消波區(qū)，設(shè)定阻尼波長(zhǎng)度大致為定常興波波長(zhǎng)的2倍，以消除波反射影響結(jié)果。

計(jì)算網(wǎng)格的劃分如圖5所示，為在保證計(jì)算精度的同時(shí)提高計(jì)算效率，加密部分網(wǎng)格。在凱爾文波系范圍內(nèi)進(jìn)行網(wǎng)格加密，精準(zhǔn)捕捉波形，防止數(shù)值耗散，自由面垂向方向要保證足夠網(wǎng)格，氣室內(nèi)部特別是質(zhì)量源區(qū)域網(wǎng)格要加密，提高氣體流場(chǎng)信息的準(zhǔn)確度。

圖5 計(jì)算域網(wǎng)格劃分Fig.5 Grid division of calculation domain

2.2 興波阻力系數(shù)計(jì)算與對(duì)比

對(duì)氣墊模型進(jìn)行了12個(gè)傅汝德數(shù)下的靜水航行興波模擬計(jì)算，在確定質(zhì)量源q大小時(shí)，先根據(jù)式(6)估算出大致范圍，因氣室與水面之間圍合區(qū)域的面積隨著自由面形狀不同而異，在計(jì)算每種傅汝德數(shù)的情況時(shí)，需微調(diào)q的值，使氣墊底面平均壓強(qiáng)大致滿足壓長(zhǎng)比的要求。

表1給出了不同傅汝德數(shù)下氣墊船靜水航行中的興波阻力及其系數(shù)，以及模擬所采用的質(zhì)量源項(xiàng)的大小。將興波阻力系數(shù)隨傅汝德數(shù)變化的數(shù)據(jù)繪制成圖6。

表1 氣墊船興波阻力及其系數(shù)計(jì)算結(jié)果

圖6 氣墊船興波阻力系數(shù)隨傅汝德數(shù)變化的曲線圖Fig.6 Curve of wave making resistance coefficient of ACV changing with Froude number

由圖6可知，曲線有2個(gè)較為明顯的阻力峰，第1個(gè)阻力峰發(fā)生在Fr為0.35時(shí)，第2個(gè)阻力峰發(fā)生在Fr為0.6附近時(shí)。氣墊船底部與水面接觸部分的壓力減去大氣壓力，僅為大氣壓力的2%～6%，壓力沖量相比一般排水型船舶產(chǎn)生的沖量很小。此外，隨著航速提高，水面受到的壓力也不會(huì)劇增，因而興波阻力不會(huì)隨之大量增加，反而可能會(huì)減弱，這是全墊升氣墊船具有高速性能的原因之一。

圖6中給出了其他學(xué)者在相同工況下的模擬結(jié)果。本文的結(jié)果與Bhushan[7]的CFD方法較為吻合；與Everst等[15]的試驗(yàn)結(jié)果相比，曲線變化趨勢(shì)尤其是2個(gè)阻力峰發(fā)生的傅汝德數(shù)基本一致?？傮w來(lái)看，本文計(jì)算氣墊船興波阻力的方法具有一定有效性和可靠性。

2.3 興波波形分析

圖7給出了5個(gè)不同傅汝德數(shù)時(shí)氣墊船靜水航行時(shí)的自由表面興波波形?？梢钥闯觯退贂r(shí)，凱爾文波系中的橫波較為明顯，隨航速的增加，橫波不再明顯，凱爾文角逐漸減小，氣墊壓力中心逐漸后移。

圖7 氣墊船自由面興波Fig.7 Waveform inducing by ACV

圖8給出了4種工況下計(jì)算域中縱剖面上的自由面升高圖。圖中的橫坐標(biāo)x/LC為縱向長(zhǎng)度和氣墊船長(zhǎng)的比值，0為氣墊船縱向中點(diǎn)的位置。結(jié)合圖8，在航速較低(Fr=0.3)時(shí)，在氣墊船區(qū)域有2個(gè)較為明顯的波谷，隨著航速增大，氣墊區(qū)域第1個(gè)波谷的位置不斷后移，當(dāng)Fr為0.7時(shí)已處于船體后部。圖中也可看出氣墊深度隨航速的變化，在Fr為0.7附近時(shí)有極大值，當(dāng)傅汝德數(shù)偏離0.7時(shí)，氣墊深度有減小的趨勢(shì)。

圖8 氣墊各航速下中縱剖面自由面升高曲線Fig.8 Free-form surface elevation curve of mid-longitudinal section at each speed of air cushion

進(jìn)一步分析波形，計(jì)算出各個(gè)航速下氣墊正下方興波水表面的平均縱傾角θw，如圖9所示。

圖9 氣墊下方水表面平均縱傾角隨航速的變化曲線Fig.9 Curve of average pitch angle of water surface under air cushion with speed

對(duì)比興波阻力系數(shù)變化圖(圖6)和平均縱傾角變化圖(圖9)可見(jiàn)，兩者變化趨勢(shì)十分相似，平均縱傾角在Fr=0.35及Fr=0.6附近各有1個(gè)峰值。

因此，在氣墊船底部壓力一定的情況下，興波阻力的主要影響因素之一是氣墊下方水表面的平均縱傾角。氣墊壓力的強(qiáng)度與航速關(guān)系不密切，航速更多地影響底部波形分布，因此在后續(xù)對(duì)氣墊船阻力性能研究中應(yīng)著重關(guān)注設(shè)計(jì)航速下的氣墊底部波形分布。

3 淺水中的氣墊船興波阻力變化

氣墊船因具有兩棲性能，對(duì)淺水中氣墊船的阻力特性進(jìn)行研究具有研究意義。前文針對(duì)深水情況的氣墊船興波情況進(jìn)行了計(jì)算研究，驗(yàn)證了該方法的有效性與可靠性，可運(yùn)用該方法對(duì)淺水中氣墊船的興波阻力變化情況進(jìn)行研究。

針對(duì)近似無(wú)限水深(d/Lc=4)和極淺水(d/Lc=0.25)2種情況，分別計(jì)算了10個(gè)傅汝德數(shù)下靜水航行的興波阻力系數(shù)。氣墊及流域參數(shù)如表2所示。

表2 氣墊船及流域參數(shù)Table 2 Parameter of ACV and flow field

與第2節(jié)中控制底部壓力分布情況不同，本節(jié)控制質(zhì)量源項(xiàng)，即各航速下的進(jìn)氣流量保持一致，如此設(shè)置更接近真實(shí)情況，試驗(yàn)中控制風(fēng)機(jī)流量較控制底部壓力也更為便捷，方便對(duì)比結(jié)果。模擬得到的結(jié)果如圖10所示。

圖10 2種不同水深下的氣墊船興波阻力系數(shù)曲線Fig.10 Curve of wave resistance coefficient of ACV under two different water depths

由圖10可知，水域的深淺對(duì)興波阻力影響很大，具有較為顯著的淺水效應(yīng)。淺水(d/L=0.25)與深水(d/L=4.0)下氣墊船均有2個(gè)阻力峰，第1個(gè)阻力峰的位置兩者基本相同，淺水的第1個(gè)阻力峰值較低，而淺水下的第2個(gè)阻力峰發(fā)生在Fr=0.5附近，深水發(fā)生在Fr=0.7左右，淺水的第2個(gè)阻力峰值顯著大于深水的值。在越過(guò)阻力峰后，淺水下的阻力值較深水下要低。

因此氣墊船的淺水效應(yīng)大致與常規(guī)船型的規(guī)律一致，在臨界速度附近阻力會(huì)出現(xiàn)極大值，一般為氣墊船興波的第2個(gè)阻力峰，在相關(guān)設(shè)計(jì)時(shí)需考慮該處影響。

4 結(jié)論

1)文中計(jì)算了不同航速下給定平均壓力值的氣墊船興波阻力，結(jié)果準(zhǔn)確給出了氣墊船興波阻力隨速度變化的2個(gè)明顯的阻力峰，與相關(guān)實(shí)驗(yàn)結(jié)果對(duì)比吻合良好。

2)運(yùn)用同一方法對(duì)淺水中的氣墊船興波阻力變化情況進(jìn)行了研究，表明氣墊船同樣有淺水效應(yīng)，在臨界速度時(shí)興波阻力會(huì)出現(xiàn)極大值，也是第2個(gè)阻力峰出現(xiàn)的位置。

綜上可知，論文提出的質(zhì)量源方法結(jié)合CFD平臺(tái)可以模擬移動(dòng)氣墊引起的水表面變形，方法簡(jiǎn)單、計(jì)算穩(wěn)定、高效，能夠獲得準(zhǔn)確的興波阻力，較為客觀反映實(shí)際氣墊船移動(dòng)興波特性。