梁景奇,徐保成,王 瑞,祁曉斌,李瑞杰

(西北機電工程研究所, 陜西咸陽 712099)

0 引言

超空泡航行器大部分表面被空泡包覆,航行阻力顯著減小,航行速度大幅提高。超空泡射彈就是一種利用了超空泡減阻技術的動能武器,具備水下高亞聲速甚至超聲速航行的能力,可以對魚雷、水雷等水下目標實現(xiàn)快速毀傷[1]。

超空泡射彈的重要特征之一是射彈的流體動力受到超空泡與彈體之間相對位置耦合關系的影響。當航速在300～1 200 m/s時,受到初始擾動后,包裹于超空泡內(nèi)的射彈開始周期性地與空泡壁面發(fā)生碰撞,即尾拍運動。Kulkarni[2]等人基于Milwitzky[3]的理論采用加入附加質(zhì)量的方法研究超空泡運動體的“尾拍”力;Euteneuer[4]建立了超空泡運動體的動力學模型,將后體滑行力作為外界高頻擾動看待,忽略了空泡的時間延遲效應;王茂勵[5]基于彈道仿真技術研究了流體動力參數(shù)的攝動對彈道穩(wěn)定性的影響規(guī)律;魏英杰[6]等提出了超空泡射彈流體動力特性的動態(tài)描述方法,建立了射彈尾拍運動的無約束動力學方程;張廣[7]等運用CFX建立了超空泡航行體三維6DOF數(shù)值模型,發(fā)現(xiàn)無控高速超空泡航行體具有運動穩(wěn)定性;時素果[8]等開展了自由飛試驗,獲得了不同預置舵角下超空泡航行體水平運動的彈道特性。水下高速超空泡射彈的彈道穩(wěn)定性與其尾拍運動具有重要關系,彈道失穩(wěn)會造成射擊精度和射程不足,影響射彈的作戰(zhàn)效能。因此,對射彈空泡流場與彈體耦合的尾拍特性的進一步研究極具必要性。

依托CFD軟件FLUENT18.0及其二次開發(fā),結合動網(wǎng)格移動計算域技術,忽略重力因素,在考慮液體壓縮性的基礎上建立了高速射彈剛體運動與空泡流場耦合數(shù)值模型,從理論及實驗方面驗證了數(shù)值模型的準確性;求解超空泡射彈水下縱平面運動,探究水下運動穩(wěn)定規(guī)律,分析了初始速度對射彈尾拍特性的影響。

1 數(shù)值模型

1.1 控制方程

采用VOF多相流模型模擬相界面運動。VOF是一種在固定Euler網(wǎng)格下的界面捕捉法,常用于由兩種及以上不相混液體組成的流體中,適用于多相間有清晰界面的流動。

1.1.1 連續(xù)性方程

(1)

(2)

(3)

式中:ρ是流體混合密度;u是混合速度;n1是相數(shù),由于考慮了不凝氣體,在這里,n1=3;αk、ρk和μk分別為第k相的體積分數(shù)、密度和速度。

1.1.2 動量方程

勻質(zhì)平衡流模型認為各相間不存在速度差,在流體微元中基于平均密度和平均動力粘度的動量守恒方程可描述為:

(4)

式中：ρm、μm分別為流體微元的各相平均密度和動力粘度,且有:

(5)

(6)

1.2 空化模型

采用Schnerr and Sauer空化模型模擬超空泡射彈的水下穩(wěn)定空化繞流,Schnerr and Sauer空化模型將汽相體積分數(shù)和單位體積流體含有的空泡數(shù)量聯(lián)系起來,表達式為:

(7)

(8)

式中:ρl為液體密度;ρv為水蒸氣密度;ρ為混合相密度；pv為飽和為蒸氣壓;p為當前壓力;αnuc為氣核體積分數(shù);RB為氣核空泡直徑;n2為單位體積內(nèi)的空泡數(shù)量。

1.3 湍流模型

Realizabledk-ε湍流模型主要針對充分發(fā)展的湍流,穩(wěn)定性好,近壁面區(qū)域的流動狀態(tài)使用壁面函數(shù)預測,對邊界層網(wǎng)格要求較為寬松;結合尺度化壁面函數(shù),其在不增加計算量和保持模型穩(wěn)定性的前提下增加了模型的適用范圍,模擬效果好。

湍流強度k方程:

(9)

湍流耗散率ε的方程:

(10)

式中:μt為湍動粘度;μ是流體的時均速度;σk、σε分別為k、ε方程的湍流能量普朗特數(shù);C1、C2為經(jīng)驗常數(shù);E為時均應變率;ν為運動粘度;xi、xj為各方向距離。

1.4 可壓縮液體模型

Tait方程是通過采用非線性回歸的方法,對能夠反應p-v-T三者關系的試驗數(shù)據(jù)進行擬合,而得到的液體狀態(tài)方程,廣泛應用于描述可壓縮液體的物性。沒有溫度修正的簡化Tait液體狀態(tài)方程可描述為:

(11)

K=K0+n3ΔP

(12)

Δp=p-p0

(13)

(14)

式中:p0為參考壓力;ρ0為參考壓力下的液體密度;K0為參考壓力下的液體體積彈性模型;n3為密度指數(shù);p為當前壓力;ρ為當前壓力下的液體密度;K為當前壓力下的液體體積彈性模量;c為水中聲速。

1.5 動網(wǎng)格耦合運動

耦合運動即射彈運動與空泡流場的計算相互耦合。已知前一時刻剛體重心位置和偏轉(zhuǎn)角,通過對物體表面壓力和剪切應力積分得到流體動力和力矩,再根據(jù)剛體運動方程計算物體運動的平移速度和角速度,然后重新計算重心位置和偏轉(zhuǎn)角。

(15)

(16)

式中:G為變換矩陣,剛體的位置與方向根據(jù)每時間步線速度和角速度的變化而變化,c.g.表示質(zhì)心。

采用移動計算域技術模擬彈體運動,計算過程中僅涉及計算域的移動,不存在網(wǎng)格的變形與重構,計算效率高,結果一致性好。

2 計算模型及方法驗證

2.1 模型、網(wǎng)格劃分及邊界條件

文中研究的射彈模型如圖1所示,射彈采用圓盤空化器,前部為兩段錐段,中間為圓柱段,柱段尾部安裝有6片尾翼。前端空化器直徑3.2 mm,圓柱段最大直徑15 mm,質(zhì)量為0.23 kg。

圖1 射彈模型圖

采用圓柱形計算域,直徑取50倍尾截面空泡直徑,計算域軸向長度為11倍彈長,前端邊界距離空化器4倍彈長,后端邊界距離彈尾6倍彈長,該計算域徑向尺度可以忽略空泡阻塞效應[9]。

針對所建立的三維計算域,采用ICEM軟件的O-Block技術劃分全結構化網(wǎng)格,如圖2所示在彈體周圍3 mm范圍內(nèi)的流域劃分外O-block用于設置邊界層網(wǎng)格,近壁面添加邊界層網(wǎng)格,并根據(jù)y+值對網(wǎng)格進行優(yōu)化。劃分網(wǎng)格時特別注意在空泡兩相交界面位置進行網(wǎng)格加密,最終劃分的網(wǎng)格總數(shù)約80萬,網(wǎng)格質(zhì)量均在0.6以上。

圖2 剖面網(wǎng)格劃分圖

計算域四周邊界均設置為壓力入口條件(pressure-inlet),設置靜壓值且靜壓值隨深度變化,射彈表面的邊界條件設置為壁面(wall),并且壁面與臨界網(wǎng)格相對靜止。計算域和邊界條件設置如圖3所示。

圖3 計算域及邊界條件圖

2.2 數(shù)值模型驗證

在超空泡射彈運動中,運動特性取決于流體動力特性,流體動力取決于彈體與空泡的相對位置關系,因此空泡形態(tài)計算準確性是運動模擬準確性的前提。將文獻[10]經(jīng)驗公式和文獻[11]試驗結果與數(shù)值模型的空泡外形作對比,驗證數(shù)值模型準確性。

2.2.1 經(jīng)驗公式驗證

Savechenko[10]等在烏克蘭國家科學院水動力實驗室進行了空泡截面獨立擴張原理有效性驗證試驗,給出了空泡外形描述經(jīng)驗公式:

(17)

(18)

(19)

式中:Rn為空化器半徑;Rc,max、Lc,max分別為超空泡最大截面半徑和空泡全長;x為空泡截面到空化器的距離;σ為空化數(shù),σ=0.012～1.057。

對比同工況下經(jīng)驗公式和數(shù)值模型計算的空泡輪廓,驗證數(shù)值模型的準確性。取水深1 m,航速100 m/s,空化器半徑1.6 mm,前段空泡輪廓對比如圖4所示。

圖4 空泡輪廓對比圖

圖4對比結果顯示,同工況下經(jīng)驗公式得到的空泡輪廓在徑向尺度略大于CFD計算的空泡輪廓,計算誤差小于12%,可以證明數(shù)值模型的準確性。

2.2.2 試驗結果驗證

Hrubes[11]采用試驗方法,在水深4 m,速度970 m/s條件下水平發(fā)射平頭錐形射彈,利用高速攝像機捕捉射彈入水后的空泡形態(tài)。試驗結果如圖5所示,射彈完全被空泡包裹。運用建立的數(shù)值模型,采用與文獻[11]相同的彈體模型及初始條件,模擬射彈水下直航運動,空泡輪廓對比如圖6所示。

圖5 970 m/s初速試驗與數(shù)值模擬空泡形態(tài)圖

圖6 數(shù)值模擬與試驗空泡輪廓對比圖

圖6給出了Hrubes試驗數(shù)據(jù)點和數(shù)值模擬的空泡外形曲線。從中可以看出,數(shù)值模擬結果與試驗數(shù)據(jù)點基本吻合,均能有效包裹射彈,最大誤差不超過10%,表明文中數(shù)值模型合理可行。

3 計算結果及分析

采用變參數(shù)法探究不同初速下超空泡射彈尾拍的流體動力特性及運動特性,水深恒定1 m,忽略重力、海流和波浪影響,考慮液體壓縮性,模擬定深繞質(zhì)心旋轉(zhuǎn)的縱平面尾拍運動[12]。先定常計算初始空泡流場,再基于三維動網(wǎng)格模型并結合二次開發(fā)技術模擬自由減速運動,探究運動穩(wěn)定性及流體動力、彈道參數(shù)的變化規(guī)律,初始條件如表1所示。

表1 初始條件表

3.1 流體動力特性

探究1 m水深、30 rad/s角速度擾動條件下不同初速射彈的流體動力變化特性,以計算時長0.015 s為計算終止條件,時間步長根據(jù)計算殘差在1e-08到1e-05之間逐漸增加,計算結果如圖7、圖8所示。

圖7 尾拍升力對比圖

圖8 尾拍阻力對比圖

圖7、圖8分別為3種速度工況尾拍升、阻力對比曲線。初始擾動使彈體軸線轉(zhuǎn)動,并與空泡壁面發(fā)生周期性碰撞,彈體的尾拍升、阻力均呈現(xiàn)周期性振蕩特性。根據(jù)曲線特征,可將受擾動后超空泡射彈的振蕩過程分為發(fā)散階段和收斂階段。在第一階段,由于彈體穿刺空泡深度不斷增加,尾拍升、阻力曲線均呈發(fā)散特性,振蕩幅值迅速增大并達到最大,振蕩周期不斷減小;在第二階段,由于速度衰減,升、阻力曲線呈收斂特性,振蕩幅值逐漸減小,振蕩周期緩慢增大,兩個階段升、阻力曲線的幅值和周期均同步變化,同一時刻升力幅值約為阻力幅值的3.5倍。

尾拍力與速度的平方正相關,初速越大,振蕩周期越小,尾拍升、阻力越快達到最大幅值,最大幅值越大,均值越大;隨著速度的衰減,3種工況下阻力和升力均值的差距逐漸減小;0.01 s后時間步長增至1e-05,振蕩曲線呈階梯形,但依然能客觀反映尾拍力振蕩特性。

3.2 空泡形態(tài)

以1 200 m/s工況為例,分別取振蕩發(fā)散階段和振蕩收斂階段的一個尾拍周期T,分析尾拍過程中空泡形態(tài)變化及彈體與空泡的位置關系,結果如圖9、圖10所示。

圖9 振蕩發(fā)散階段尾拍動態(tài)圖(第一周期)

圖9為振蕩發(fā)散階段(第一個尾拍周期)內(nèi)的空泡動態(tài)變化過程。彈體向下偏轉(zhuǎn),尾翼局部穿刺空泡下壁面后產(chǎn)生二次空泡,升、阻力達到最大,彈體回彈并觸碰上壁面,最終回到起始位置。在該階段,尾拍角幅值隨穿刺深度的增加呈增大趨勢,尾拍角為0°時,彈體除空化器外不沾濕,空泡軸線與彈體軸線基本重合。

圖10為振蕩收斂階段(第6個尾拍周期)的空泡動態(tài)變化過程,在該階段,彈體穿刺空泡深度更大,空泡形態(tài)發(fā)生顯著改變。到達最大尾拍角時,除尾翼外,彈體第二錐段和靠近尾翼的圓柱段也發(fā)生沾濕,且兩次沾濕區(qū)域具有不對稱性,而尾拍角為0°時彈體依然存在沾濕情況,這是由于空泡變化與彈體振蕩周期存在相位差,空泡軸線變化滯后于射彈軸線變化,即空泡時間延遲效應。

3.3 運動特性

以地面系為參考坐標系,分析3種初速射彈受擾動后的尾拍運動特性,俯仰角、俯仰角速度、軸向速度、質(zhì)心位移對比曲線分別如圖11～圖14所示。

圖11、圖12為俯仰角和俯仰角速度隨時間變化曲線。不同速度下的俯仰角、俯仰角速度曲線均呈周期性振蕩特性,同樣可分為振蕩發(fā)散和振蕩收斂階段,變化規(guī)律與流體動力曲線相似。由俯仰角曲線可知,初速越大,振蕩周期越小,振蕩發(fā)散階段俯仰角幅值越小;在振蕩收斂階段,3種工況俯仰角幅值均接近4°,基本保持恒定,這表明,振蕩收斂階段彈體俯仰角幅值受初速影響極小。

圖11 俯仰角隨時間變化曲線

圖12 俯仰角速度隨時間變化曲線

由俯仰角速度曲線可知,初速越大,俯仰角速度曲線振蕩周期越小,越快達到最大幅值,最大幅值越大,1 200 m/s工況最大俯仰角約為270 rad/s,增幅高達9倍;達到最大幅值后,俯仰角速度振蕩幅值均隨時間逐步衰減。通過與圖7、圖8對比可知,尾拍力曲線俯仰角曲線存在相位差,亦是由空泡延遲效應所致。

圖13為X向速度隨時間變化曲線,3條曲線均呈單調(diào)遞減趨勢,初速度越大,阻力越大,因此速度衰減越快,0.014 s內(nèi)1 200 m/s、1 000 m/s、800 m/s工況的速度衰減幅度分別為260 m/s、190 m/s和130 m/s;隨著時間的增加3種工況航速差距逐漸縮小,由此可知通過增加初速增加射程效果有限。

圖13 X向速度隨時間變化曲線

圖14為0.015 s內(nèi)3種工況射彈質(zhì)心位移變化曲線,質(zhì)心隨著彈體的振蕩均產(chǎn)生了側向位移,且位移量不斷增加。以1 200 m/s為例,軸向位移15 m時,側向位移為0.02 m,與軸向位移量相比,側向位移量可忽略不計,彈體基本保持了直線航行,這表明超空泡射彈的水下彈道具有動態(tài)穩(wěn)定性。

圖14 質(zhì)心位移變化曲線

4 結論

文中基于建立的流場運動耦合數(shù)值模型,數(shù)值仿真了超空泡射彈的尾拍運動,分析了初速對尾拍力的特性及運動特性的影響,結論如下:

1)受擾動后超空泡射彈的水下運動具有動態(tài)穩(wěn)定特性,尾拍力、俯仰角、俯仰角速度等參數(shù)均呈現(xiàn)準周期振蕩特性,可將該過程分為振蕩發(fā)散階段和振蕩收斂階段；

2)尾拍運動使空泡軸線產(chǎn)生偏移,空泡軸線變化滯后于射彈軸線變化,同一尾拍周期的兩次沾濕區(qū)域具有不對稱性；

3)初速度越大,尾拍周期越小,彈體尾拍升、阻力及俯仰角速度振蕩幅值越大,俯仰角振蕩幅值受初速影響極小。