旋轉圓球入水空泡特性與流場結構的大渦模擬研究1)

嚴晨祎 陳 瑛

(上海交通大學船舶海洋與建筑工程學院,水動力學教育部重點實驗室,上海 200240)

引言

物體穿越水面的入水問題與許多軍事應用直接相關,如:艦載發射航行體入水[1]、船舶砰擊[2]、飛機迫降[3]等.物體入水問題的研究已經持續了近一個世紀,這個經典問題的影響因素有很多:物體幾何形狀[4-6]、沖擊速度和角度[5]、物體密度和表面材質[6]等.當引入變化的物體密度[7]和物體入水帶有旋轉時[8],問題變得更加復雜.研究圓球高速旋轉入水過程的空泡特性和精細流場結構,可為將來基于多物體入水[9]的新型航行體入水降載方式提供關鍵性的力學機理認識.

入水問題的關注點主要包括入水沖擊和入水空泡兩個方面.Von Karman[10]提出的附加質量方法為入水沖擊研究提供了先驅性理論指導,其后國際上的理論研究大多圍繞該理論展開.對于圓球入水問題,Watanabe[11]開展的相關試驗為后人理論模型的發展提供了關鍵的驗證數據.入水空泡的研究最早始于19 世紀末Worthington 和Cole[12]的工作,其采用電子閃光攝像捕捉圓球的垂直入水過程,首次系統地闡述了水面飛濺、空泡閉合和射流等現象.Truscott等[8,13-14]對圓球高速旋轉入水作了系統性的試驗探究,包括采用不同入水速度和入水轉速,球面采用不同接觸角的涂層以及采用不同密度的圓球進行試驗.他們發現由于圓球表面的無滑移邊界,使得圓球旋轉會在空泡內部形成橫向楔形射流,并導致圓球彈道會發生明顯的橫向偏轉.

近年來,高性能計算機的發展使得數值模擬成為求解入水問題越來越重要的手段.數值模擬在湍流場細節和渦結構時空演變的精細捕捉方面,相比試驗技術有天然的優勢.而對于入水問題,數值模擬方法有雷諾平均(RANS)類方法,大渦模擬(LES)類方法,也有部分直接數值模擬方法(DNS)[15],還有一些基于拉格朗日體系的無網格方法[16-18]等.

RANS 類方法是最常用的湍流模擬方法,Yu等[19]采用RANS 方法模擬了圓球的入水過程,發現無量綱化之后的沖擊力與入水速度和圓球半徑都無關,只和球的質量有關.周波等[20]使用SST 模型對不同旋轉方向強制運動圓球的傾斜入水過程進行數值模擬,發現球體繞不同軸旋轉傾斜入水時的運動受力特性的差異顯著.

然而傳統RANS 方法采用的渦黏湍流模型對于空泡流這種包含大范圍流動分離和渦旋流動的情形難以獲得較好的計算結果.相對而言大渦模擬類的高精度數值模擬方法(LES,DES,LES/RANS 混合等) 在湍流空泡流領域的應用在近年來取得了長足的進步.Li等[21]建立了六自由度的圓盤撞擊水面空泡演化的數值計算模型,他們發現沖擊生成的空泡在圓盤前進方向存在非對稱與旋轉的現象.Zhao等[22]對低弗勞德數下圓柱斜入水問題進行了數值模擬和試驗探究,他們利用LES 與VOF (volume of fluid)相結合的方法求解湍流場和氣液界面.并且成功捕捉了多種尺度下圓柱周圍的渦結構,同時對多種入水參數進行了數值模擬.

綜上所述,目前入水問題的研究方法眾多,而旋轉圓球入水還有許多復雜的現象值得深入探究,在這一方面采用LES 類方法進行模擬的相關研究較少.本文采用LES 結合動網格與滑移網格技術,基于壁面邊界層八叉樹各向同性網格,對低弗勞德數下旋轉圓球入水的自由運動過程開展數值模擬,并與試驗結果[7-8]進行對比,著重分析空泡形態演化和水中彈道隨入水旋轉速度的變化規律,揭示圓球高速旋轉入水過程的空泡流場精細結構,為多物體入水降載問題提供基礎力學機理的認識.

1 數學模型

1.1 基本控制方程

本文采用均質多相流的VOF 模型,氣液混合介質由各相體積分數比例混合而成,并將這種多相混合物看作一種單一的變密度流體.各相介質共享相同的速度、壓力、湍流量等各種流動物理量,混合物的物性參數則通過各相加權平均獲得.

混合介質的流場控制方程如下

1.2 大渦模擬濾波的控制方程

在大渦模擬(LES)中,物理量的求解是通過濾波核函數G(x,x′) 在空間上對上述方程進行濾波來得到的.濾波過程能很好的拆分尺度大于濾波長度的渦與尺度小于濾波長度的渦.濾波之后的變量如下定義

濾波后小尺度變量為 ?′=?-?ˉ,而在有限體積法的數值求解過程中,網格尺度V是隱式的濾波尺度,即G(x,x′)=1/V通過對式(1)和式(2)進行濾波運算之后,能得到以下濾波之后的控制方程

其中τij為亞格子應力,定義如下

1.3 亞格子尺度應力模型(SGS)

對N-S 方程濾波之后的亞格子應力項是未知的,所以需要對其進行建模來求解.采用Boussinesq假定[23],亞格子應力計算方法如下

其中μt為亞格子渦黏系數.采用常規濾波方法,所以各向同性部分τkk可以與濾波之后的壓力項合并[24],為濾波之后的應變率張量定義如下

而對于亞格子渦黏系數μt,采用Nicoud 與Ducros[25]提出的WALE (wall-adapting local eddyviscosity)模型.該模型優勢在于剪切層流的μt=0,相比傳統的Smagorinsky-Lilly 模型[26]能較好地返回流體區域中層流的特性,并且在壁面擁有更好的對數律y3性質.該模型具體形式如下

其中馮卡門常數κ=0.418 7,d為單元中心到壁面距離,ΔV為單元體積.WALE 常數采用Cw=0.325,這個值對于大多數流動擁有更好的適應性.

1.4 表面張力模型

當氣液固三相接觸時,在三相交界處產生接觸線與接觸角.當溫度和壓強給定,在平衡狀態下這個接觸角被稱為靜態接觸角.然而通常在試驗中能觀察到移動的接觸線,這個接觸線對應的接觸角為動態接觸角(如圖1 所示).而在入水過程中動態接觸角對空泡的產生與形態有很重要的影響[27].

圖1 入水參數與浸潤性示意圖Fig.1 The schematic description for the impact and wettability parameters

針對這個問題,本文采用 Brackbill等[28]提出的CSS (continuum surface stress)表面張力模型描述球體表面的潤濕性,將表面張力作為VOF 動量方程的源項,具體模型如下

CSS 模型避免顯式計算界面曲率,將表面張力影響作為各向異性的毛細應力來計算,從而能對界面上的尖角得到較好的模擬結果.而對于氣液界面與球體壁面的接觸線,也引入了壁面依附約束模型(wall adhesion) 來處理.將流體與壁面的接觸角用于計算距離壁面的1 個流體單元界面法向的條件,如下所示

其中θw為壁面依附角,為壁面法向單位矢量,為壁面切向單位矢量.

2 數值算法與驗證

本文采用有限體積法對流動控制方程進行時空離散,同時采用基于壓力的分離式求解器(SIMPLE)對控制方程中耦合的非線性項進行解耦迭代求解.LES 需要解析流場中大部分高波數的湍流脈動,從而需要高精度的數值格式以及高精度的網格,而偶階次的數值格式天生不具有耗散性,很好的規避了對SGS 模型在耗散上的影響.所以采用中心差分格式能很好的適應LES 方法,但是中心差分格式會對流場引入數值震蕩,所以最終我們選擇有界中心差分格式[29]作為動量的離散格式.而對于方程中變量的梯度項,采用最小二乘格式來進行計算.由于入水過程為低速所以忽略空化效應的影響.本文采用了VOF 多相流模型對氣液界面進行捕捉,而迎風格式對耗散會過估計,中心差分則有數值震蕩都無法很好的對界面進行捕捉.為了克服這些問題,采用一種高分辨率界面捕捉算法(HRIC)進行空泡面捕捉[30],相比QUICK 與其他二階格式擁有更高的計算精度,相比傳統的幾何重構格式計算速度能更快.圓球的三自由度運動采用動網格方法中的動框架來實現計算域平移運動控制,圓球帶空泡在水中平移和旋轉運動采用相對于背景計算域的滑移網格來控制.

2.1 計算域及邊界條件

為了校驗當前數值模擬方案的可靠性,與Aristoff等[3]的垂直入水試驗進行了對比.球直徑為D=0.025 4 m,半徑為R,分別對比了鋼和聚丙烯兩種不同密度比m*=ρl/ρs=7.86,0.86 的情況.如圖2 所示,整個計算域采用球直徑R進行標準化,垂直長度為100R,水平長度為40R,寬度為20R,球心距離頂部壓力邊界設置為80R,距離側面速度邊界為20R,在寬度方向選取了對稱面,同時在球面上設置為無滑移邊界.水面位置距離球為3R.坐標軸原點同水面保持一致如圖3 所示,其中Vl代表圓球入水時左側壁面速度大小,Vr代表圓球入水時右側壁面速度大小,圓球各項初始化參數如表1 所示.

表1 Aristoff等[7]垂直入水參數Table 1 Aristoff et al.[7] experimental condition

圖2 計算區域示意圖Fig.2 Computational domain of sphere

圖3 對稱面坐標軸示意圖Fig.3 Coordinate system in cut view

具體計算網格采用結構化與非結構化相結合的方式,盡可能采用結構化網格以減小數值誤差(如圖4 所示).為了實現精度較高的流場大渦模擬,近壁面邊界層附近采用盡量各向同性和均勻的網格,以減小大渦模擬低通濾波算法中濾波函數的交換誤差.計算網格采用盡量合理的網格疏密分布.在物體表面布置致密網格,網格的精細度達到能夠刻畫湍流邊界層結構和典型湍流擬序結構的程度.除了滿足數值模擬對網格質量的一般性要求,即保證無量綱壁面垂向距離y+＜1 之外,還要盡量使壁面流向和展向的無量綱網格單元大小滿足Δx+≤100 和Δz+≤30.對于當前模擬工況Δx+,Δz+列在表2 中,選取的時間為圓球剛剛沉入水面以下的時刻(t=10 ms 左右).

圖4 近壁面網格示意圖Fig.4 Illumination of near wall mesh

表2 流向與展向平均無量綱網格尺度表(t=10 ms)Table 2 Averaged nondimensional grid sizes of the finest mesh in the streamwise and spanwise (t=10 ms)

而另一種普遍采用的衡量LES 解析度的方式是Benard等[31]建議整體網格量能計算出流動區域內至少80%的湍動能.將未解析湍動能比率定義為

其中參數C選取為100.

雖然湍流脈動速度和湍動能ER在LES 中不是顯式求解獲得的物理量,但是對于那些宏觀時均流動幾乎定常的流場,可以通過數據采樣統計的方法獲得和ER.雖然圓球入水空泡流動并不是宏觀定常流場,但是可以采用同樣的網格對圓球在水中重力場中無空泡垂直降落過程這樣類似的時均定常流場進行LES 模擬,進行數據采樣統計從而獲得脈動速度并計算湍動能,以檢驗所采用的網格分辨率的可靠性.所以模擬了圓球初始速度Vo=5.45 m/s 不帶空泡在水中下沉,其中t∈[0.07,0.18] (1100 時間步數據采樣),開啟數據采樣,得到均方根(RMS)的來計算ER與Esgs.

他逢到小小的失意也就哭，沒有一點隱藏?；ㄉ追涞降厣狭耍约航劳瓷囝^了，小貓不肯吃糕了，都是他哭的原因。有一天，他自己失手，把外婆送給他的泥人跌碎了，他哭得異常悲傷，和旁人死了親人沒有兩樣。

圖5 展示了對稱面上Qc與ER大小,可以發現在ER較高的圓球尾流區Qc遠遠小于0.2.所以可以認為在流場中絕大多數區域Qc＜0.2 的要求都能被滿足.

圖5 對稱面上(a)大渦模擬解析湍動能ER,(b)湍動能未解析率Qc,采用數據采樣方法,采樣點1100 (Vo=5.45 m/s,ωo=0 rad/s)Fig.5 (a) Resolved turbulence kinetic energy ER and (b) its resolution Qc at the symmetry,obtained by the data sampling method using 1100 sampling points (Vo=5.45 m/s,ωo=0 rad/s)

根據經典湍流理論,湍流場能量在慣性子區滿足-5/3 次方律(Kolmogorov 的K41 理論)[32]

而湍流壓力脈動的能譜滿足文獻[33]的-7/3 次方律

這也是LES 研究論文中經常采用的驗證方法.為了進一步驗證所采用的網格分辨率的可靠性,實時監測獲取穩定之后流場中的壓力測點(在隨球運動的動坐標系內,測點坐標x=0,y=3R,z=0)的壓力信息,采用快速傅里葉變換得到歸一化的功率譜密度估計PSD (見圖6).其中采樣頻率Fs=104,NFFT=1024.可以直接測量如下的雙對數坐標PSD (power spectral density)圖中的峰值包絡線斜率,發現與上述經典湍流理論中的-7/3 次方律非常接近.

圖6 歸一化的測點壓力功率譜密度及其與經典湍流理論-7/3 次方律的對比(Vo=5.45 m/s,ωo=0 rad/s)Fig.6 Normalized PSD of the pressure in comparison with the-7/3 power spectrum of Kolmogorov’s K41 theory

2.2 垂直入水與試驗對比驗證結果

圖7 給出了圓球垂直入水數值模擬結果與文獻[7]的試驗結果的對比.從圖中可以看出,當前模擬方案在總體上對圓球深度有些低估,但總體上模擬與試驗結果是非常接近的.表3 選擇了幾個特定時間,對比了圓球模擬深度與試驗深度,可以發現相對誤差不大于10%,總體來說誤差水平是可以接受的.

表3 垂直入水圓球深度試驗[7]與模擬對比表(Vo=2.17 m/s)Table 3 Comparison between the experimental[7] and simulation vertical depths for the case at Vo=2.17 m/s

圖7 垂直入水圓球深度試驗[7]與模擬對比(Vo=2.17 m/s)Fig.7 Comparison between the experimental[7] and simulation vertical depths for the case at Vo=2.17 m/s

圖8 展示了數值模擬和試驗泡型的對比.此處初始時刻定義同文獻[7]維持一致,即將圓球球心通過水面位置作為初始時刻.雖然在圖8(a) 68.9 ms 與圖8 (b) 62.2 ms 接近空泡斷裂時上半部分空泡略有偏小.但總體對于定性的泡形,當前的模擬方案也能和試驗吻合較好.特別是在空泡斷裂之后(圖8 (a)75.9 ms 與圖8 (b) 68.7 ms),還能捕捉到空泡面的收縮與向上的“Worthington”射流.

圖8 圓球垂直入水試驗結果(上)與數值模擬結果(下)的對比Fig.8 Experimental (up) and simulation (down) results of the evolution of cavity shapes from the side view

圖8 圓球垂直入水試驗結果(上)與數值模擬結果(下)的對比(續)Fig.8 Experimental (up) and simulation (down) results of the evolution of cavity shapes from the side view (continued)

3 不同轉速對圓球旋轉入水影響

本文主要考察相同入水速度Vo=5.45 m/s 下不同入水轉速ωo=50,100,150,199 rad/s 對于入水空泡演化、圓球彈道、流場各個物性參數的影響.所有初始時刻t=0 ms 選取在圓球剛剛接觸自由液面的時刻.

3.1 入水角速度對圓球運動學性質影響

本節將詳細討論圓球高速旋轉入水的水動力學特性,主要關注在圓球的彈道(圖9)、速度(圖10)、加速度(圖11) 以及圓球受力情況(圖12) 隨時間t的變化.

圖9 不同入水轉速圓球入水彈道對比(ωo=199,150,100,50 rad/s,Vo=5.45 m/s)Fig.9 Trajectories of all spinning cases (ωo=199,150,100,50 rad/s,Vo=5.45 m/s)

圖10 圓球質心速度與角速度隨時間t 變化Vo=5.45 m/sFig.10 The velocity of sphere as a function of time with Vo=5.45 m/s

圖11 圓球加速度與角加速度隨時間t 變化 Vo=5.45 m/sFig.11 The acceleration of sphere as a function of time with Vo=5.45 m/s

圖12 圓球所受升阻力系數隨時間變化Vo=5.45 m/sFig.12 Time-history of the simulation coefficients with Vo=5.45 m/s

圖9 展示了當圓球入水速度為Vo=5.45 m/s 時,不同入水旋轉角速度ωo=50,100,150,199 rad/s 對于圓球彈道的影響.其中為了和文獻[8]的試驗結果進行對比,使用密度比m*=1.74,半徑R=28.6 mm 的圓球進行模擬.可以發現ω0=199 rad/s 的模擬結果與文獻[8]的試驗結果吻合較好.圓球入水后在馬格努斯力的作用下彈道發生偏轉,在入水的初始階段不同轉速ωo對彈道的影響并不是很大,但是隨著圓球進入水面以下(t=10 ms 之后,見圖13),入水轉速越大的球,相應的彈道偏轉也越大.當空泡發生深閉合時,空泡斷裂時圓球所在深度都非常近似.

為了進一步探究彈道偏轉與入水轉速的關系,在圖10 中展示了圓球質心速度與角速度隨時間變化關系圖,其中采用入水速度Vo與角速度ωo對速度和角速度進行歸一化.通過比較歸一化之后的結果可以發現,入水轉速的改變對水平方向速度vx/Vo的影響很大(見圖10 (a)).相應的轉速越大,水平方向速度大小也越大,從而使得彈道偏轉也越大.對于ωo=50 rad/s,vx/Vo近似于線性下降,而對于ωo=100,150,199 rad/s 都隨時間t,先快速下降(t=10～20 ms),然后趨于常數(t＞ 60 ms).通過對比x方向加速度ax隨時間變化(見圖11(a))可以發現,入水砰擊階段(t＜ 5 ms)ωo=50 rad/s 的加速度峰值不到ωo=100,150,199 rad/s 的30%,并且當砰擊階段結束,ax也趨于平穩,這就導致了vx/Vo是近似的線性下降.而對于垂直方向速度vy/Vo,并不受到ωo改變的影響(見圖10 (b)),所有工況下都呈現在入水砰擊階段(t＜ 5 ms),下降較快,之后下降趨勢逐漸放緩的情況.而圓球轉速ω/ωo隨時間t變化總體不大(見圖10 (c),min (ω/ωo)=0.94).

在圖11 中展示了圓球加速度與角加速度隨時間t的變化,在入水砰擊階段 (見圖11(a)和圖11(b)t=0～4 ms),y方向加速度ay到達峰值之后快速下降,然后趨于不變并且不同ωo對于ay影響不大.但是對于ωo=199,150,100 rad/s 時,ax繼續變大直到圓球大部分進入自由液面以下(t=6 ms)而到達峰值.峰值的大小隨ωo變大而變大,而ωo=199,150 rad/s 兩個情況峰值大小非常接近,說明入水砰擊階段ax峰值大小受到Vo大小的限制.在t=20 ms 之后ax趨于穩定直到空泡斷裂(t＞ 100 ms).圓球角加速度β(見圖11 (c)),在砰擊階段震蕩過后都基本維持不變,β的大小也是隨著ωo的變大而變大.結合圓球表面未沾濕面積之比Aair/Asphere隨時間t變化(見圖14),可以發現ωo越大,未沾濕比例越小,水依附在球面的面積越大,從而導致球面上受到黏性力也而越大,所以β隨著ωo的變大而變大.過了入水砰擊階段(t＞ 5 ms)Aair/Asphere基本維持不變,使得β維持不變.

圖14 圓球表面未沾濕面積Aair/Asphere 比例隨時間t 變化Vo=5.45 m/sFig.14 Time-history of non-wetted area ratio Aair/Asphere with Vo=5.45 m/s

基于上述運動學的討論,圓球質心的受力可以做如下分析

其中ms為小球質量,Fhx,Fhy為球面上所受黏性力與壓力的合力,球面上還應有表面張力Fst對球質心的作用,根據式(15)表面張力大小的最大值不會超過Fst=πDσ,其大小不到重力1%:6σ/(D2ρsg) ≈ 10-3,ρs為小球密度.所以為了簡化分析,忽略了表面張力對球受力的影響[15,34].

圓球所受升力Fl,阻力Fd的方向定義如圖15所示,故根據式(21)和式(22)可以將升阻力化簡為如下形式

圖15 圓球所受升阻力方向簡圖Fig.15 Schematic illumination of direction of drag and lift forces

其中彈道切線角度γ=arctan(-vx/vy),從而得到相應的升阻力系數

其中A=0.5πR2為圓球截面面積.

根據上述定義在圖12 中展示了圓球所受升力系數Cl,阻力系數Cd隨時間t變化.可以發現Cl與Cd的變化趨勢分別與ax,ay比較接近.Cd在度過入水砰擊階段之后(t＞5 ms)從Cd=0.2 開始隨時間線性下降直到空泡斷裂,并且ωo對其影響并不是很明顯.但是對于升力系數Cl,ωo越大入水初始階段(t=0～10 ms)的峰值越大.

3.2 入水角速度對空泡形態的影響

本節展示了ωo=199,150,100,50 rad/s 的空泡泡形(圖13),對稱面上壓力、速度云圖,以及Q準則等值面,其中Q準則為速度梯度張量的第二不變量定義如下[35].

圖13 圓球入水空泡演化對比圖(Vo=5.45 m/s)Fig.13 Time evolution of cavity during water entry of the spinning spheres from the longitudinal view (Vo=5.45 m/s)

圖13 展示了空泡隨時間的演化序列,其中球面上為水氣混合密度的云圖,即不同顏色的接觸線代表了圓球的水氣接觸線,在圖13(a) 中定性比較了文獻[8]的試驗與當前的模擬結果.可以發現模擬泡形相比試驗結果略小一些,這在空泡斷裂之后較為明顯(圖13 (a)實驗111 ms 與模擬110.8 ms).這可能是由于隨著圓球的下沉網格逐漸稀疏導致的.同時在試驗中飛濺冠在液面以上先行閉合,而在當前模擬中無法很好的模擬水面以上的飛濺冠使得這一現象無法觀察到.在圖13 和圖16 可以看到有橫向楔形劈尖在入水的初始階段(t=10 ms 附近)開始生成侵入空泡內部,結合圖11中ax隨t變化可以發現,橫向劈尖在球面上的攀升是ax變大的主要原因.當t＞20 ms 時,橫向劈尖在球面上不再繼續延展 (見圖13(a)～圖13(c))ax也趨于穩定.隨著圓球的下沉,橫向劈尖逐漸發展直到擊穿泡面(ωo=199,150 rad/s).而相同時刻,ωo越大橫向劈尖也越大(見圖16).而對于ωo=100 rad/s,上半部分楔形劈尖還未到達另一邊空泡面空泡已經發生斷裂.

圖16 圓球入水空泡演化俯視圖(Vo=5.45 m/s)Fig.16 Time evolution of cavity during water entry of the spinning spheres from the top view ( Vo=5.45 m/s)

對于ωo=50 rad/s 的情況(圖13,圖16 (d)),能觀察到空泡在自由液面上發生表面閉合,也沒有明顯的橫向劈尖產生(圖16 (d)),這是由于ωo=50 rad/s 對應的無量綱旋轉參數S=0.262＜1,對應Vr＜0,也就是接觸線大部分都在圓球赤道以下,所以沒有橫向劈尖的產生,與文獻[8]的結論保持一致.

入水初始階段(t=10 ms 附近)與空泡斷裂后(t=110 ms 附近)不同ωo的壓力云圖展示在圖17 中,其中實線代表混合密度ρm=400 kg/m3自由液面,可以看出由于馬格努斯效應圓球一側形成低壓區并且隨著ωo變大,低壓區也在變大,這也使得彈道的橫向位移也越大.同時ωo的變化對于空泡斷裂后的壓力場有顯著的影響.對于ωo=100 rad/s (見圖17(c)),在空泡的斷裂位置有明顯的高壓區,但是對于ωo=199 rad/s (見圖17(a)),并不能觀察到明顯的高壓區.這可以如下解釋:橫向劈尖侵入泡內從而影響空泡深閉合的斷裂,而ωo越大橫向楔形越大,從而使得水進入泡內更快.最終結果就是空泡因靜水壓夾斷時,產生的高壓區更小乃至幾乎沒有.

圖17 圓球對稱面壓力云圖(Vo=5.45 m/s)Fig.17 Pressure contours of the spinning spheres from symmetric view(Vo=5.45 m/s)

圖18 展示了ωo=199,150,100 rad/s 在100 ms 附近時的對稱面速度云圖與矢量場,對比圖13可以發現雖然此時空泡還未斷裂,但在對稱面上由于橫向楔形劈尖的存在,已經接近斷裂.觀察速度云圖可以發現ωo=100 rad/s 在自由液面上沒有出現明顯的橫向氣流射出與ωo=199,150 rad/s 有所不同,這可以被解釋為ωo較小,橫向劈尖速度小,使得自由液面變形速度較慢從而抑制了橫向氣流的生成.

圖18 圓球對稱面速度大小與矢量場(Vo=5.45 m/s)Fig.18 Velocity and vector field of the spinning spheres from the symmetric view ( Vo=5.45 m/s)

為了探究渦結構對空泡泡形的影響,圖19 展示了Q=100 000 s-2的等值面,可以發現渦結構生成的時間點有兩個:入水初始階段(t=11 ms 附近)與空泡斷裂階段(t=90～111 ms).這一點與垂直入水的情況保持一致,但橫向劈尖的存在使得空泡上半部分始終能觀察到相應的渦結構.ωo的增加使入水初始階段的渦結構更加趨向于非對稱,但是到了空泡斷裂階段,ωo越大渦結構越少.這可以被解釋為小球橫向位移的增加,使得泡內空氣與外界交換受到阻礙,抑制了渦結構的生成.

圖19 旋轉圓球入水Q=100 000 s-2 等值面 (Vo=5.45 m/s)Fig.19 Vortex Q iso-surface of the spinning spheres from longitudinalview (Vo=5.45 m/s)

4 結論

本文采用大渦模擬對圓球高速旋轉入水建立了數值模型.研究了圓球在相同入水速度Vo=5.45 m/s 下不同入水轉速ωo=50,100,150,199 rad/s的空泡演化與水動力學特性,主要結論如下.

(1)數值模擬結果通過與垂直入水,高速旋轉入水試驗對比,在圓球彈道,空泡泡形上與試驗結果有較好的一致性,驗證了數值模型的正確性與可靠性.

(2)入水轉速的改變對圓球的水動力學特性有顯著的影響.隨著入水轉速的增加,圓球彈道水平位移更大,圓球水平方向獲得的速度更大,在入水砰擊階段圓球所獲得的升力也更大,但升力的峰值受到入水速度的限制.而入水轉速對于圓球垂直方向特性影響不大,包括垂直方向速度,加速度,空泡斷裂時圓球所在深度.

(3)入水轉速對于空泡形態與流場特性得到詳細的探究.在模擬中觀察到ωo=50 rad/s 的表面閉合ωo=100,150,199 rad/s 的深閉合.對于深閉合,入水轉速的增加帶來了更強的橫向楔形射流并且抑制了空泡斷裂帶來的高壓以及相應渦結構的生成.