海流對出水空泡演化過程影響機(jī)理數(shù)值研究

劉元清，崔 軍

(北京宇航系統(tǒng)工程研究所，北京 100076)

0 引言

海洋中的海水沿一定路徑大規(guī)模流動形成了海流[1]，海流流向、路徑和速度均隨時間和空間呈現(xiàn)出一定的變化。航行體水中及出水過程空泡形態(tài)受到海流作用，會造成空泡不對稱性發(fā)生變化，從而影響航行體水中運(yùn)動軌跡和出水姿態(tài)[1-3]，是水下航行體研制過程中應(yīng)重點(diǎn)關(guān)注的問題。

根據(jù)流動形式不同，海流分為近海岸潮流與遠(yuǎn)海岸環(huán)流。近海岸潮流指與潮汐運(yùn)動關(guān)聯(lián)的海水運(yùn)動，有往復(fù)流和旋轉(zhuǎn)流；流速的大小和方向隨時間存在周期性變化，變化周期一般小于12 h；受近海岸地形影響，不同深度的流速存在切變。遠(yuǎn)海岸環(huán)流指在海洋中相對穩(wěn)定的流動；流速和流切變一般較小，流動變化周期較長。海流參數(shù)一般包括：流動速度、流動方向、隨時間的周期性及隨空間的切變。

由于物理問題的復(fù)雜性，考慮浪流與非定常、非線性空泡流動耦合影響的航行體水下發(fā)射研究工作開展相對較少，水下航行體對海洋環(huán)境的適應(yīng)性研究，主要集中在波浪對航行體流體動力特征及運(yùn)動姿態(tài)的影響方面，多在流場邊界采用Stokes非線性造波法開展波浪浪級、浪向及相位影響研究[4-5]。對于海流條件影響的分析較少，一般采用理論計算方法，通過施加隨深度線性變化的水流速度，開展水下航行體出水姿態(tài)參數(shù)研究[6-8]，近年來，有學(xué)者采用機(jī)理模型實驗及三維數(shù)值仿真的方法，研究橫流對水下發(fā)射航行體流體動力的影響[9-10]。

本文采用數(shù)值仿真方法，重點(diǎn)針對海流對空泡演化過程的影響開展研究，形成考慮海流影響的空泡多相流理論與數(shù)值計算方法，分析了海流流向與流速的影響，獲得了海流對空泡發(fā)展演化的作用機(jī)制及對航行體流體動力特征的影響規(guī)律。

1 考慮海流影響的航行體出水空泡多相流數(shù)學(xué)模型

1.1 物理模型及網(wǎng)格劃分

建立計算模型時，需將航行體對流的速度V12分解為航行體對地的速度V10與地對流的速度V02的矢量和

V12=V10+V02

(1)

式中，海流運(yùn)動即地對流的速度V02通過設(shè)置計算域邊界條件給定海流的流速、流向及沿深度方向的流切變，航行體平動即航行體對地的速度V10則通過設(shè)置動框架沿流向、垂直方向機(jī)橫向移動的速度給定，而航行體轉(zhuǎn)動及俯仰平面內(nèi)運(yùn)動角速度則通過設(shè)置球面網(wǎng)格轉(zhuǎn)動速度進(jìn)行給定。

計算域網(wǎng)格劃分如圖1所示，整個流體計算域被劃分為兩部分，一部分為包裹航行體的球形區(qū)域，按照指定彈道進(jìn)行法向運(yùn)動；另一部分為遠(yuǎn)場流體域，計算過程中保持靜止。兩部分區(qū)域在內(nèi)部邊界相接，每迭代計算一步即在內(nèi)部邊界通過差值進(jìn)行數(shù)據(jù)更新，算例網(wǎng)格總數(shù)約1.05×106。

圖1 計算域網(wǎng)格劃分Fig.1 Computational domain and local refined grid

1.2 多相流模型

為精確捕捉流動現(xiàn)象，采用VOF多相流模型，體積分?jǐn)?shù)控制方程采用具有截面重構(gòu)性質(zhì)的PLIC格式。針對水、蒸汽和空氣等多相流場，認(rèn)為不同組分的流體介質(zhì)具有相同的運(yùn)動速度，用同一速度、壓力變量描述流場。航行體出水過程的流場滿足質(zhì)量守恒、動量守恒、能量守恒規(guī)律，可構(gòu)筑如下控制方程。

(1)連續(xù)性方程

(2)

式中，ρ為網(wǎng)格單元流體的平均密度，ui為速度分量，其中坐標(biāo)變量x1=x,x2=y,x3=z。

(2)動量方程

(3)

式中，gi為重力分量，p為壓力，τ為黏性應(yīng)力張量，對于牛頓流體來說，黏性應(yīng)力張量為

(4)

式中，μ為動力黏度系數(shù)。

(3)能量方程

(5)

式中，E為流體微團(tuán)能量,T為溫度,J為介質(zhì)的擴(kuò)散通量。

(6)

在以上各方程中，物性參數(shù)密度ρ和黏性系數(shù)μ為各相體積分?jǐn)?shù)α的加權(quán)平均值

ρ=ρw(1-αg-αv)+ρgαg+ρvαv

(7)

μ=μw(1-αg-αv)+μgαg+μvαv

(8)

式中，下標(biāo)w，g，v分別代表水、空氣和蒸汽相。

1.3 空化模型

本文空化模型采用兩相Yuan模型。該模型的輸運(yùn)源項基于微氣泡的泡動力學(xué)模型，因此模型對壓力的變化較為敏感，在非定?？张萘饔嬎阒芯哂辛己帽憩F(xiàn)，適用于出水中環(huán)境壓力發(fā)生迅速變化空泡流的計算。同時，對于有清晰交界面的出水問題，VOF算法在捕捉交界面時相對其他算法具有明顯的優(yōu)勢，而Yuan模型空泡模型對VOF算法具有較好的適用性。

假定出水過程中，空泡與大氣沒有發(fā)生聯(lián)通，因此認(rèn)為在同一計算單元不會同時存在蒸汽和大氣?？栈话l(fā)生在水汽填充的區(qū)域，在αg不為零的單元，認(rèn)為不會發(fā)生空化。采用的空化模型還是基于水汽兩相模型，蒸汽相連續(xù)性方程的右端源相Sv的表達(dá)式為

(9)

其中

(10)

式中n0為單位體積水內(nèi)含有的氣核數(shù)，一般為給定常數(shù)，本文中為1.5×1010；Rn為氣核半徑。將式(10)代入(9)，得到

(11)

氣核半徑隨時間的變化率dRn/dt可以由Rayleigh-Plesset方程求解，但更方便的處理是采用下式直接給出

(12)

2 海流影響數(shù)值方法校驗

為了驗證文中所采用的計算方法與計算模型，采用上述的計算方法對實驗相同工況進(jìn)行數(shù)值模擬研究，獲取航行體表面壓力系數(shù)分布、空泡壓力與推進(jìn)長度，并與實驗數(shù)據(jù)進(jìn)行對比分析，如圖2～4所示。研究結(jié)果表明，本文采用的計算方法與數(shù)學(xué)模型不僅能夠比較準(zhǔn)確地預(yù)測空泡形態(tài)變化規(guī)律，而且對于空泡泡壓的計算也可以取得令人滿意的精度，因此，該數(shù)值方法與數(shù)學(xué)模型可以滿足文中研究需要。

(a)本文計算結(jié)果

(b)實驗測量結(jié)果圖2 航行體表面壓力對比Fig.2 Comparison of pressure distribution

圖3 空泡內(nèi)壓力對比Fig.3 Comparison of predicted and measured ventilated cavity pressure

圖4 空泡推進(jìn)速度Fig.4 Comparison of predicted and measured ventilated cavity length

3 海流對出水空泡發(fā)展及流體動力特征影響分析

3.1 考慮航行體俯仰運(yùn)動下航向海流對出水空泡的影響研究

在航行過程中，航行體姿態(tài)是不斷變化的，航行體的運(yùn)動速度方向也在不斷變化，因此造成航行體質(zhì)心攻角的變化。而海流的空間變化及航行體的轉(zhuǎn)動使不同位置的局部攻角也不相同。局部攻角可以由下式求得

(13)

式中,Vs為海流速度，Vx為航行體運(yùn)動垂向速度，如果Vx按照30 m/s計算,由海流引起的局部攻角變化如表1所示，而流體與航行體運(yùn)動的相對攻角是影響空泡在出水過程中發(fā)展規(guī)律的主要因素之一。

表1 海流速度引起的相對攻角Tab.1 Attack angle induced by different sea route velocities

將海流參數(shù)(包括流速、流向)抽象為5種典型海流特征開展數(shù)值仿真計算，見表2和圖5。

表2 海流速度的數(shù)學(xué)表達(dá)Tab.2 Calculation cases with different sea route mode

圖5 給定海流速度分布Fig.5 Sea route velocity distribution of different cases

圖6 各工況迎流面空泡發(fā)展比較Fig.6 Comparison of predicted ventilated cavity length in different cases

通過提取并對比出水時刻各工況下迎流面空泡發(fā)展過程(見圖6)，從空泡推進(jìn)規(guī)律看，海流形式對空泡發(fā)展的影響不顯著。海流流速增大，迎水面空泡泡長略短，背水面泡長略有增長。主要是因為各海流的差別較小，以水深0 m時的海流速度相比，前4種工況的差別不超過0.2 m/s，其引起的攻角變化不超過0.5°。海流作用機(jī)制主要為海流疊加引起的局部攻角變化使得空泡不對稱性加劇。

通過對航行體表面各壓力監(jiān)測點(diǎn)空泡潰滅峰值出現(xiàn)時間進(jìn)行對比分析，5種工況迎流面沿航行體表面各測點(diǎn)的潰滅時間見圖7。從空泡潰滅特征看，計算結(jié)果表現(xiàn)出一定的階梯特征。從潰滅峰值時間來看，工況5中的來流速度最大，迎流面的空泡也相對更薄一些，因此各截面潰滅出現(xiàn)的時間也相對更早。

圖7 各工況迎流面測點(diǎn)潰滅時間Fig.7 Comparison of cavity collapse time

出水時，不同海流工況下空泡的形態(tài)見圖8。比較工況4與工況5，隨海流速度增大，迎流面的空泡厚度略薄。工況4中，泡內(nèi)回射水體更多地出現(xiàn)在背流面。

圖8 航向海流對空泡形態(tài)影響Fig.8 Cavity shape impacted by different sea route mode in flow direction

3.2 不同橫向海流對出水空泡的影響研究

為研究不同橫向來流速度對空泡形態(tài)的影響，橫流速度分別取0，0.25，0.5，0.75 m/s，分別給出了不同橫流條件下空泡形態(tài)變化及空泡泡長發(fā)展演化對比分析。

不同橫向來流作用下的空泡形態(tài)差別不大，隨橫向海流速度增加，空泡對稱性發(fā)生變化。受橫向來流影響，空泡沿橫流方向由迎流面向背流面發(fā)生偏斜，由圖9，10可以看出，背流面空泡泡長較長、空泡較厚。通過對比空泡厚度最大橫截面迎、背流面空泡形態(tài)可知，來流對空泡的橫向擠壓作用是造成空泡橫向偏轉(zhuǎn)的主要原因。

(a) 橫向來流速度(0 m/s)

(b) 橫向來流速度(0.75 m/s)圖9 某時刻不同橫流速度下橫向空泡形態(tài)對比(視角為沿來流方向)Fig.9 Cavity shape impacted by different sea route velocity in cross direction

(a)橫向來流0 m/s

(b)橫向來流0.25 m/s

(c)橫向來流0.5 m/s

(d)橫向來流0.75 m/s圖10 不同橫向來流空泡最大橫截面形狀Fig.10 Cavity shape in cross-section impacted by different sea route velocity in cross direction

4 結(jié)論

通過建立考慮海流影響的航行體垂直及航行體出水過程空泡多相流數(shù)值計算模型，得到海流對出水空泡演化過程的作用特點(diǎn)，對其流體動力特征進(jìn)行深入分析，并研究流速、流向等重要海流參數(shù)對出水空泡演化過程的影響機(jī)理，并得到以下結(jié)論：

1) 在航行體帶空泡出水過程中，海流的不同會影響航行體運(yùn)動狀態(tài)和姿態(tài)，因此，其直接和間接地影響肩部的局部攻角變化，從而影響空泡的發(fā)展過程。

2) 航向海流對潰滅階段的空泡形態(tài)變化影響較為明顯，不同的海流狀態(tài)下，泡內(nèi)水、氣分布有所不同，高壓產(chǎn)生次序及潰滅過程的流動結(jié)構(gòu)變化也因此不同。

3) 不同橫向來流作用下的空泡形態(tài)差別不大，受到合成來流的方向發(fā)生偏轉(zhuǎn)的影響，迎流面向橫向來流方向偏斜。