氣泡高速艇波浪中阻力及運動性能數(shù)值研究

歐勇鵬， 周廣禮， 吳浩

(海軍工程大學(xué) 艦船工程系，湖北 武漢 430033)

目前，船舶氣層減阻技術(shù)已在俄羅斯、烏克蘭、日本、荷蘭等國家獲得了廣泛應(yīng)用，取得了顯著的軍事經(jīng)濟效益[1]。我國自20世紀(jì)80年代開始，在該領(lǐng)域開展了相關(guān)理論及試驗研究，取得了靜水中模型絕對減阻率達(dá)20%～25%的效果[2-4]，但至今尚未實現(xiàn)實艇應(yīng)用，主要原因即是在考慮風(fēng)浪作用下氣泡高速艇的減阻效果及運動性能尚未明確。氣泡高速艇在波浪中是否仍然具備良好的減阻效果？波浪中的艇底氣層是否能保持穩(wěn)定狀態(tài)？氣層作用下高速艇在波浪中的運動性能如何？氣層與搖蕩運動艇體之間的相互作用力學(xué)規(guī)律如何？這一系列問題尚需深入研究和解答。

針對以上問題，國內(nèi)外公開發(fā)表的相關(guān)理論及試驗研究不多。在試驗研究方面，董文才等[5-6]通過拖曳水池試驗探索了波浪中的氣層減阻效果以及氣層對艇體運動的影響規(guī)律，但鑒于試驗條件限制，未能給出艇體運動、氣層形態(tài)、阻力的同步變化過程。在理論研究方面，美國的Jin-Keun Choi等[7]聯(lián)合采用面元法與求解URANS的方法研究了噴氣下瘦長型船的船體興波、氣層形態(tài)及船體運動對氣層分布的影響，但該方法在計算船體興波時采用線性勢流理論，不適用于計算高速艇的噴濺等非線性流動問題。

為此，本文考慮了高速艇流場的強非線性及氣-液兩相流的復(fù)雜性，采用粘性CFD的方法，從波浪中阻力、縱向運動及氣層形態(tài)等角度進(jìn)行研究，試圖通過理論計算揭示波浪中氣層-艇體的相互作用過程，掌握氣泡高速艇在波浪中的航行性能。

1 數(shù)值波浪水池

1.1 研究對象

計算對象為一條底部設(shè)置凹槽的B.H.型氣泡高速艇模型，主尺度參數(shù)如表1所示，圖1給出了艇體及底部凹槽的三維效果圖。其中，Δ為排水量，L為艇體總長，LP為折角線長，BPX為折角線最大寬度，βM為中部艇底斜升角，βT為尾部艇底斜升角；LP/BPX表示折角線長與折角線最大寬度之比，HB/BPX表示艇體凹槽深度與折角線最大寬度之比，SBH/SP表示水平面上的艇底凹槽投影面積與折角線投影面積之比，Ryy/L表示縱向慣性半徑與艇長之比。

表1 計算模型的主尺度參數(shù)

圖1 計算對象三維視圖Fig.1 3-D geometry model of hull

1.2 兩相流模型

氣-液兩相流采用VOF模型(volume of fluid model)，界面追蹤采用高分辨率HRIC方法(high-resolution interface capturing scheme)?？紤]界面壓縮修正的VOF方法的控制方程如式(1)～(5)所示，詳細(xì)可參見文獻(xiàn)[8]。

質(zhì)量守恒方程：

(1)

動量方程：

-αqp+·(μqαquq)+

ρqαqg+FD,q+Fs,q

(2)

式中：ρq、αq和uq分別表示第q相的密度、體積分?jǐn)?shù)和流動速度；FD,q和Fs,q分別為相間曳力與表面張力，本文計算過程中忽略了這兩項的影響；uc為界面壓縮速度，其表達(dá)式為uc=Cα|u|(α/|α|)，式中|u|為流動速度，Cα為銳化因子，取值范圍為0～1。

界面追蹤采用高分辨率HRIC格式，其基礎(chǔ)是NVD方法(normalized variable diagram)，原理示意圖如圖2所示。圖中，αD、αA、αU分別表示施主單元、受主單元及迎風(fēng)單元的體積分?jǐn)?shù)，αf表示控制體邊界的體積分?jǐn)?shù)，其計算表達(dá)式為：

圖2 NVD方法原理示意圖Fig.2 The schematic diagram of NVD method

圖3 角度θf的定義Fig.3 Definition of angle θf

HRIC方法在求解過程中考慮了界面與網(wǎng)格夾角、當(dāng)?shù)貛炖蕯?shù)等因素的影響，其求解的主要方程和步驟如下[9-10]：

1)采用NVD方法計算歸一化體積分?jǐn)?shù)

(3)

2)考慮角度θf的影響，通過權(quán)重因子γf對界面?zhèn)巫冃芜M(jìn)行修正。

(4)

式中：θf定義為施主單元D指向受主單元A的單元向量d與界面法向n之間的夾角(如圖3所示)，γf=(cosθ)Cθ，C0稱為角度因子。

3)考慮當(dāng)?shù)貛炖蕯?shù)的影響，對計算的不穩(wěn)定性進(jìn)行修正。

(5)

式中：CUL為庫朗數(shù)下界值，CUU為庫朗數(shù)上界值，Cf為當(dāng)?shù)貛炖蕯?shù)。

1.3 湍流模型

湍流模型采用在船舶流體性能計算中被廣泛應(yīng)用和驗證的Realizablek-ε模型[11-12]，其基本數(shù)學(xué)表達(dá)形式如下。

湍動能方程(k方程)：

(6)

耗散率方程(ε方程)：

(7)

(8)

(9)

1.4 數(shù)值造波

數(shù)值波浪采用線性波，無限水深下的波浪參數(shù)滿足：

(10)

式中：ω為圓頻率，k為波數(shù)，a為波幅，U為邊界入口處的流體速度，x和z分別為流體的縱向和垂向坐標(biāo)位置。此外，為保持壓力出口邊界條件的穩(wěn)定性，在距離壓力出口1倍艇長區(qū)域內(nèi)添加阻尼項，動量源的函數(shù)表達(dá)形式為

(11)

1.5 計算域及網(wǎng)格

考慮到流動關(guān)于艇體中縱剖面對稱，僅對流場的一半進(jìn)行網(wǎng)格離散。計算流域分成隨艇體一起運動的重疊網(wǎng)格區(qū)域和固定不動的背景區(qū)域，如圖4所示。這兩個區(qū)域均采用剪切型網(wǎng)格進(jìn)行離散，區(qū)域之間的數(shù)據(jù)傳遞采用重疊網(wǎng)格技術(shù)。

重疊網(wǎng)格區(qū)域前方距艇艏0.25L，后方距艇艉0.3L，上表面距甲板為0.3L，下表面距艇底0.25L，側(cè)面距舷側(cè)為0.3L，均設(shè)置為重疊網(wǎng)格邊界。中縱剖面設(shè)置為對稱面，艇體表面設(shè)置為無滑移物面，噴氣口位于艇底凹槽的艏部，設(shè)置為速度入口。

背景區(qū)域在計算過程中保持靜止，該區(qū)域總長5.5L、總寬2L、總高3.1L，上表面距離艇體甲板為1.3L，下表面距離艇底部為1.7L，側(cè)面距離艇體舷側(cè)為1.8L，均設(shè)置為對稱面；入口距離艇艏為1.5L，設(shè)置為速度入口；出口距離艇艉為3.0L，設(shè)置為壓力出口；中縱剖面設(shè)置為對稱面。

圖5給出了重疊網(wǎng)格區(qū)域及背景區(qū)域的網(wǎng)格劃分情況。對艇底凹槽、艏部噴濺及艉部興波等區(qū)域進(jìn)行了網(wǎng)格加密，艇體附近的網(wǎng)格尺度為5 mm，噴濺處的網(wǎng)格尺度為3 mm，邊界層網(wǎng)格厚度為3 mm，設(shè)置4層網(wǎng)格，沿法向的增長率為1.25。對自由液面附近的網(wǎng)格進(jìn)行了加密，垂直自由液面方向的網(wǎng)格尺度為4 mm，同時為保證艇體姿態(tài)變化過程中重疊區(qū)的網(wǎng)格尺度基本相同，自由液面加密區(qū)的高度在重疊網(wǎng)格區(qū)域中為0.06 m，在背景區(qū)域中為0.1 m，且在兩者相互重疊的位置處進(jìn)行適當(dāng)增加(如圖5(b)所示)。網(wǎng)格總數(shù)為120萬，其中重疊網(wǎng)格區(qū)域的網(wǎng)格數(shù)為90萬，背景區(qū)域的網(wǎng)格數(shù)為30萬。

圖4 區(qū)域劃分與邊界條件Fig.4 Calculation region and boundary condition

圖5 計算區(qū)域的網(wǎng)格劃分Fig.5 The meshes of calculation region

2 數(shù)值方法驗證

艇體縱搖和垂蕩的計算采用6-DOF剛體運動模型，僅計及沿垂直方向的平動及繞橫軸方向的轉(zhuǎn)動，限制了其余自由度的運動。非穩(wěn)態(tài)時間離散采用二階格式，時間步長Δt=5.0×10-4s。空氣的密度為1.184 kg/m3，動力粘性系數(shù)為1.855×10-5Pa·s；水的密度為997.561 kg/m3，動力粘性系數(shù)為8.887×10-4Pa·s。

在計算過程中為避免氣-液相在艇體表面發(fā)生偽擴散現(xiàn)象，考慮模型尺度、航速、時間步長等因素的影響，HRIC-VOF方法的主要參數(shù)取值為Cа=0.7、Cθ=0.5、CUL=0.5、CUU=6.0，其中模型試驗在中國特種飛行器的高速拖曳水池中進(jìn)行，試驗結(jié)果可參見文獻(xiàn)[14]。

2.1 數(shù)值波浪

圖6給出了波長λ=4.0 m、波高h(yuǎn)w=0.06 m時充分發(fā)展后的波形圖像以及距離區(qū)域入口1.5L處(艇艏所在位置)數(shù)值計算波升與理論結(jié)果的對比。從圖6可以看出：采用本文的數(shù)值造波方法可以清晰表達(dá)所需的數(shù)值波浪，計算所得波形與理論值吻合較好，驗證了本文數(shù)值造波方法的有效性。

2.2 靜水阻力

圖7給出了飽和噴氣下的阻力、縱傾計算值與試驗值的對比，圖中Rt/Δ為單位排水量阻力。從圖7可以看出：采用該數(shù)值方法計算所得阻力、縱傾值與試驗結(jié)果吻合較好，阻力最小相對偏差可達(dá)1.48%，最大相對偏差為4.59%；縱傾值的最小相對偏差可達(dá)3.10%，絕對偏差值為0.14°?？梢?，該數(shù)值計算方法計算艇體阻力及航行姿態(tài)且具備較好的精度。圖8給出了FV=2.80時，艇體興波及氣層形態(tài)圖像的計算結(jié)果。通過試驗過程中的觀察發(fā)現(xiàn)，計算所得興波及艇底氣層在形態(tài)上與試驗結(jié)果基本一致。

2.3 頂浪中的縱向運動

圖9給出了噴氣流量Q=12 m3/h、容積傅汝德數(shù)FV=2.24、波長船長比λ/Lp=2.171、波高船長比hw/Lp=0.021 7時艇體垂蕩與縱搖的時歷過程。圖中k為入射波的波數(shù)，ζa為波幅，ZG為垂蕩值，θ為縱搖角。從圖中可以看出，當(dāng)時間t>15 s后，艇體運動處于穩(wěn)定波動狀態(tài)。

圖10為艇體運動穩(wěn)定狀態(tài)下無因次垂蕩、縱搖與試驗結(jié)果的對比。從圖中可以看出，垂蕩與縱搖的計算結(jié)果與試驗結(jié)果吻合較好，垂蕩幅值與試驗結(jié)果的偏差為6.4%，縱搖幅值與試驗結(jié)果的偏差為2.8%。說明該數(shù)值計算方法求解艇體在波浪中的搖蕩運動具備較好的可靠性。

圖6 數(shù)值造波及結(jié)果分析Fig.6 Numerical wave patterns and its Comparison with theoretical value

圖7 數(shù)值計算結(jié)果與試驗值的對比Fig.7 Comparison of numerical results with experimental results

圖8 艇體興波及艇底氣層數(shù)值計算與試驗觀測的對比Fig.8 Comparison of spray and air cavity pattern of numerical results with experimental results

圖9 垂蕩與縱搖的時歷計算過程Fig.9 the time history of heave and pitch

圖10 垂蕩與縱搖曲線與試驗結(jié)果的對比Fig.10 Comparison of heave and pitch numerical results with experimental results

3 氣層對波浪中阻力及縱向運動性能的影響

3.1 氣層對波浪中阻力的影響

圖11給出了FV=2.24時波浪中噴氣前后的總阻力變化，此時波長λ/LP=2.171、波高h(yuǎn)w/LP=0.021 7、氣流量Q=12 m3/h。從圖中可以看出：在波浪作用下，艇體阻力隨時間呈周期性波動，氣層可使波浪中的阻力大幅度降低。

圖11 噴氣對波浪中阻力的影響Fig.11 The effect of air cavity on hull resistance in wave

表2給出了不同波長下噴氣前后艇體平均阻力及減阻率的變化，從表中可以看出：在不同波長下，氣層均可使阻力大幅度降低，減阻率可達(dá)27.53%～30.41%；從表2中還可以看出，隨著波長增長，減阻率呈增大趨勢。

表2不同波長下噴氣對阻力及減阻率的影響

Table2Theeffectofairinjectiononhullresistanceatdifferentwavelength

λ/LPRt/ΔQ=12m3/hQ=0m3/h減阻率/%靜水0.16920.117430.621.450.17430.126327.531.810.17620.128227.242.170.17550.125028.772.890.17540.122730.023.620.17410.121130.41

圖12給出了不噴氣及噴氣下波浪中阻力相對靜水阻力的增加量隨波長的變化，圖中ΔR/R表示波浪中的阻力增加與靜水阻力之比。從圖中可以看出：噴氣及不噴氣條件下艇體阻力增值隨波長的變化規(guī)律相類似；但氣層作用下艇體在波浪中的阻力增值有所增大；隨著波長增加，兩者阻力增值的百分比趨于相同。

圖12 相對靜水的阻力增加隨波長的變化Fig.12 The add resistance in different wave length

3.2 氣層對縱向運動的影響

圖13給出了噴氣前后艇體垂蕩、縱搖時歷曲線的對比，此時航速FV=2.237，波長λ/LP=2.17，波高h(yuǎn)w/LP=0.021 7，噴氣時的氣流量Q=12 m3/h。圖14給出的是艇體運動處于穩(wěn)定波動狀態(tài)下時中心化時歷曲線(去均值)的對比情況。

從圖13看出：噴氣使得垂蕩和縱搖運動的波動中心發(fā)生變化，氣層使得垂蕩運動的平均值增大，而縱搖運動的平均值略有降低。去除波動運動的中心值后發(fā)現(xiàn)：氣層對垂蕩運動幅值的影響不大，使縱搖運動幅值略有降低(如圖14所示)。

圖15給出了不同波長下，噴氣前后艇體垂蕩與縱搖響應(yīng)的變化。從圖15中可以看出：頂浪規(guī)則波中，艇體氣層對垂蕩運動響應(yīng)的影響不大，對縱搖運動響應(yīng)有所影響；噴氣使得縱搖運動響應(yīng)有所減小，且隨著波長增長，氣層導(dǎo)致的縱搖響應(yīng)減小量呈增大趨勢。

4 艇體運動對氣層形態(tài)及其穩(wěn)定性的影響

4.1 氣層面積

圖16給出了穩(wěn)定運動狀態(tài)下，艇底氣層面積覆蓋率隨時間的變化，此時的波高、噴氣量、波浪等參數(shù)與上一節(jié)計算工況完全相同。其中，氣層面積覆蓋率是指氣層面積SA與艇底凹槽面積SBH之比，表征了氣層在凹槽中所占有的面積比率。

圖13 噴氣及不噴氣下艇體縱向運動時歷曲線對比Fig.13 The comparison of longitudinal motions between with and without air cavity

圖14 中心化后的縱向運動波動曲線對比Fig.14 The comparison of longitudinal motions after subtraction of mean value

圖15 不同波長下氣層對垂蕩與縱搖響應(yīng)的影響Fig.15 The effect of air cavity on heave and pitch at different wave length

從圖16可以看出：在波浪作用下，氣層面積覆蓋率呈周期性波動，這表征了艇體運動對氣層形態(tài)的干擾作用；波長不同，氣層面積覆蓋率的波動情況略有不同，但差別不大。

表3給出了不同波長下，氣層平均覆蓋率、波動幅值的統(tǒng)計結(jié)果，從表中可以看出：在不同波長下，氣層在凹槽中的平均覆蓋率可達(dá)78.67%～83.76%，且波動幅值僅為6.92%～9.54%，說明氣層在縱向規(guī)則波中未發(fā)生嚴(yán)重破碎現(xiàn)象，具備良好的穩(wěn)定性，這正是B.H.型氣泡高速艇在頂浪中仍然具備良好減阻效果的主要原因。

4.2 氣層形態(tài)

圖17給出了波長λ/LP=2.171下艇底氣層在一個運動周期中的變化過程，航速、波長、波高、噴氣量等參數(shù)與上一節(jié)的計算相一致。圖中，t1時刻表示縱搖角位于波谷位置，t2時刻表示縱搖角位于從下往上穿過平衡位置，t3表示縱搖位于波峰位置，t4表示縱搖位于從上往下穿過平衡位置。從圖中可以看出：在一個縱搖運動周期過程中，氣層基本上均能填滿整個底部凹槽，其形態(tài)與艇底凹槽形狀接近；除了氣層尾部發(fā)生小量變化外，其余變化不大。這說明采用設(shè)置凹槽的艇底構(gòu)型，可較好地保持氣層在頂浪航行狀態(tài)下的穩(wěn)定性。

表3波浪中氣層覆蓋率的平均值及波動情況

Table3Theaveragevalueofaircavitycoverratioonbottomhollow

λ/LP平均覆蓋率/%波動幅值/%1.4578.678.551.8179.888.442.1781.249.542.8982.448.533.6283.766.92

圖16 氣層面積覆蓋率隨時間的變化曲線Fig.16 Time history of air cavity cover ratio on bottom hollow

圖17 艇底氣層的形態(tài)變化Fig.17 The change history of air cavity pattern within a period of wave

5 結(jié)論

1)采用HRIC-VOF方法、Overset網(wǎng)格技術(shù)及剛體6-DOF模型相結(jié)合，可有效計算氣層作用下高速艇波浪中的阻力、垂蕩、縱搖及氣層形態(tài)。

2)B.H.型氣泡高速艇在頂浪中仍可取得良好的氣層減阻效果，波長對減阻率有所影響，但影響不大。

3)氣層對頂浪條件下B.H.型氣泡高速艇垂蕩的影響甚微；氣層對縱搖有所影響，噴氣可使縱搖響應(yīng)減小。

4)頂浪條件下，氣層在凹槽中的覆蓋率可達(dá)78.67%～83.76%，氣層面積變化不大，這是該型艇在波浪中仍具備良好噴氣減阻效果的主要原因。

[1] SVERCHKOV A V．Application of air cavities on high-speed ships in Russia[C]// Proceeding of the International Conference on Ship Drag Reduction SMOOTH-SHIPS. Istanbul, Turkey, 2010: 20-24．

[2] 董文才．滑行艇及平板氣層減阻的研究[D]．武漢：海軍工程大學(xué)，2003．

DONG Wencai. Study on the effect of air layer on planning craft and flat plate[D]. Wuhan: Naval University of Engineering, 2003.

[3] 董文才，郭日修，陳小玲，等．滑行艇氣層減阻試驗[J]．中國造船， 2002， 43(4)： 13-18．

DONG Wencai，GUO Rixiu，CHEN Xiaoling，et al．Experimental study on resistance reduction of planning craft by air[J]．Shipbuilding of China， 2002， 43(4)： 13-18．

[4] 唐桂林，倪其軍，王麗艷，等．高速氣泡艇阻力數(shù)值模擬及氣泡減阻效果分析[J]．船舶力學(xué)， 2014， 8(18)： 882-888．

TANG Guilin, NI Qijung, WANG Liyan, et al. Numerical evaluation of the resistance reduction effect for three-dimensional high-speed air cavity craft[J]. Journal of ship mechanics, 2014， 8(18)： 882-888．

[5] DONG Wencai，OU Yongpeng．Experimental study on resistance and longitudinal motion of high-speed air cavity craft[J]．Journal of ship mechanics， 2011， 9(15): 950-959．

[6] 歐勇鵬，董文才，許勇．B.H型氣泡高速艇規(guī)則波中縱向運動試驗研究[J]. 中國造船, 2010, 51(3)： 1-10．

OU Yongpeng, DONG Wencai, XU Yong. Experimental study on the longitudinal motion performances of bottom hollowed high speed air craft in regular wave[J]．Shipbuilding of China， 2010， 51(3)： 1-10．

[7] JIN-KEUN C，CHAO-TSUNG H，GEORGES L．Numerical studies on the hydrodynamic performance and the startup stability of high speed ship hulls with air plenums and air tunnels[C]//Ninth International Conference on Fast Sea Transportation FAST2007.Shanghai，China, 2007： 476-484．

[8] MUZAFERIJA S, PERIC M, SAMES P, et al. A two-fluid Navier-Stokes solver to simulate water entry[C]//Proc Twenty-Second Symposium on Naval Hydrodynamics. 1998.

[9] SHONIBARE O Y,WARDLE K E．Numerical investigation of vertical plunging jet using a hybrid multifluid-VOF multiphase CFD solver[C]//International Journal of Chemical Engineering.[S.l.], 2015, 925639.

[10] WACLAWCZYK T, KORONOWICZ T. Comparison of CICSAM and HRIC high-resolution schemes for interface capturing[J]. Journal of theoretical and applied mechanics, 2008, 46(2): 325-345.

[11] SHIH T H, LIOU W W, SHABBIR A, et al. A new k-ε eddy-viscosity model for high Reynolds number turbulent flows-model development and validation[J]. Computers & fluids, 1995， 24(3)： 227-238.

[12] YONG Xua, ZHANG Guoqing. Numerical calculation for the flow in the air-thrust bearings[C]//Procedia Engineering. [S.l.], 2011: 922-927.

[13] FENTON J D. A fifth-order stokes theory for steady waves[J]. J waterway, port, coastal and ocean engineering, 1985,111(2): 216-234.

[14] 歐勇鵬．氣泡高速艇水動力性能研究[D]．武漢：海軍工程大學(xué)， 2009．

OU Yongpeng. Study on hydrodynamic performance of high speed air cavity craft[D]. Wuhan: Naval University of Engineering, 2009.

本文引用格式：

歐勇鵬， 周廣禮， 吳浩. 氣泡高速艇波浪中阻力及運動性能數(shù)值研究[J]. 哈爾濱工程大學(xué)學(xué)報， 2017， 38(12): 1849-1857.

OU Yongpeng, ZHOU Guangli, WU Hao. Numerical investigation of the effects of air cavity on the resistance and longitudinal motion of a high-speed air cavity craft in regular waves[J]. Journal of Harbin Engineering University， 2017， 38(12): 1849-1857.