張 軍,蔡曉偉,宣建明,王英霖
(1.南京航空航天大學 航空學院 非定常空氣動力學與流動控制工業和信息化部重點實驗室,江蘇 南京 210016;2.中國船舶科學研究中心 水動力學重點實驗室,江蘇 無錫 214082;3.南京理工大學 能源與動力工程學院,江蘇 南京 210094)
潛射導彈作為一種可以進行水下發射的武器系統,具有射程遠、速度快、反應時間短、命中精度高、隱蔽性強等諸多優點,具有重要的軍事與戰略意義[1-2]。在高緯度以及部分中緯度海域,如我國黃海、渤海部分海域,由于氣候的原因,存在著一定時間的冰期[3]。海面具有大面積的冰體覆蓋,可以進一步提高潛射導彈發射的隱蔽性與突然性。
國內外均已開展了水下發射導彈出水過程的數值與試驗研究[4-7],但對于細長導彈穿越冰水混合物出水過程的研究較少,其流場的變化比單純的出水過程更加劇烈,成為導彈冰下發射成功與否的關鍵因素之一。ANSYS/LS-DYNA軟件基于有限元算法,對流固耦合問題具有良好的適應性[8-9]。本文利用ANSYS/LS-DYNA有限元分析軟件,對細長彈體穿越冰水混合區域出水過程進行數值模擬研究,獲得這一過程中流場的變化規律以及彈體質心位移變化等結果,研究成果可以為導彈冰下發射研究提供技術支持。
采用ALE算法[10],先執行一個或者幾個Lagrange時間步計算;當單元網格產生變形之后,再執行ALE時間步計算,對內部單元進行網格重新劃分,將變形網格中的單元變量和節點矢量值插值到重新劃分后的網格中。基于ALE方法,在LS-DYNA中可以方便地將Euler網格與Lagrange網格進行耦合,以處理流體與結構之間的相互作用問題。本文采用此方法進行彈體穿越冰水氣混合物流場的數值模擬研究。
采用流固耦合算法時,通常需要對Lagrange結構進行約束,將結構的相關參量傳遞給流體單元。本文采用罰函數約束[11]的方法進行耦合計算,罰函數耦合系數考慮拉格朗日節點(彈體,即從物質)和歐拉流體(水與空氣,即主物質)之間的相對位移。檢查每一個節點對主物質表面的貫穿,如果從節點不出現貫穿,就不進行任何操作;如果發生從節點對主節點的貫穿,界面力F就分布到歐拉流體的節點上。
流體區域包含水與空氣2種流體介質,水與空氣采用NULL材料模型,水介質的狀態方程采用Gruneisen狀態方程:
(1)

空氣介質采用線性多項式狀態方程:
p=C0+C1θ+C2θ2+C3θ3+(C4+C5θ+C6θ2)E
(2)

對于彈體,在出水的過程中,彈體的變形量很小,為了加快運算的速度,減少計算時間,將彈體作為剛體處理,彈體的材料模型選擇RIGID模型。表1給出了冰體材料參數,水、空氣、彈體、冰體對應的材料模型與狀態方程[12]如表2所示。

表1 冰體材料參數

表2 材料模型與狀態方程
計算模型包括細長彈體(見圖1)、冰體,以及由空氣域和水域組成的流體區域,彈體的底部直徑為118 mm,彈體長1 247 mm(質心距離頭部頂點676.12 mm),彈體材質為鋁合金,殼體厚度為5 mm。由于本文主要關注的是流場的變化以及彈體的運動姿態,因此忽略彈體的變形,將彈體作為剛體進行模擬計算,彈體采用Lagrange單元。彈體的初始速度為20 m/s,運動方向為豎直向上,攻角為0°;初始時刻,彈的運動方向和z軸正方向重合。所建立的三維有限元網格模型如圖2所示。

圖1 彈體外形側視圖

圖2 彈體出水網格模型
冰體的分布方式以及彈體相對位置對出水過程存在一定的影響,本文研究了3種不同的分布方式(規則分布、菱形分布、隨機分布,見圖3)及2種不同的相對位置(冰體正下方和冰體間隙)對質心運動軌跡的影響,共有6種計算工況,具體工況定義見表3。

圖3 冰體分布方式
對于以上6種不同工況,細長彈體穿越冰水混合區域出水時,彈體質心在x、y2個方向上的位移如圖4所示。表4給出了t=1 s時,彈體質心位移。從計算結果可以發現,彈體的初始相對位置對彈體出水后x方向質心位移影響最大,彈體從冰體間隙之間出水相較于彈體從冰體正下方出水,在x方向上的偏移要大數倍至數百倍。在同樣的冰體分布方式下,彈體質心在y方向上的位移同樣受到相對位置的影響,彈體從冰體間隙之間出水時,在y方向上的位移要明顯小于從冰體正下方出水時的位移。相對位置相同、冰體分布方式不同時,彈體質心在x方向位移差距不大,但是彈體質心在y方向位移存在著較大的偏差。

表3 計算工況條件定義

圖4 6種工況下彈體質心位移

表4 t=1 s時彈體質心的位移
對工況1進行深入的流場結果分析,圖5~圖10給出了具體的流場云圖和速度隨時間變化曲線。
彈體在離開水域時(t=1 s),流場的密度等值線圖如圖5所示,從圖中可以看出,彈體在穿越冰水混合區域后,水域流場中密度的等值線在各個方向上均有分布。彈體在出水時,水域與空氣域界面上,圍繞彈體周圍的密度等值線分布較為密集,說明此位置流體的密度梯度較大,變化較為劇烈。隨著到水-氣界面距離的增加,水域中的密度等值線變得稀少,同時水域流場的四周出現了因為空泡而形成的密集等值線圓環。浮在水面的冰體受水與彈體的作用,位置發生改變,位于彈體正上方的冰體被推離出水面,同時出現較大程度的破碎。

圖5 流場密度等值線圖

圖6 流體密度變化云圖
為了更清楚地分析彈體出水過程,選取4個時間節點,獲得截面的流場屬性云圖,見圖6~圖9。
截面上的流體密度云圖如圖6所示。從圖中可以發現:在t=5 ms時,如圖6(a)所示,在彈體的底部,流體的密度開始變小,出現了空化的現象;隨著彈體的運動,當t=50 ms時,如圖6(b)所示,隨著彈體在水中的運動,彈底部空泡體積逐步增大,彈體側面也出現了空化現象;隨著時間的進一步推進,彈體周圍空化現象加劇,使得彈體側面與彈底的空泡連成一體,如圖6(c)所示;當彈體出水時,部分的水被帶入空氣中,如圖6(d)所示。
彈體出水過程中的流場應力云圖如圖7所示。從圖中可以看到:當t=5 ms時,如圖7(a)所示,彈底處有較大的應力值,同時彈體側面也有較大的應力值,但是和彈底處相比,值要小一些;隨著彈體的運動,流場出現較大應力值的范圍逐漸擴大,逐步向遠離彈體的區域擴散,但是應力峰值開始減小,如圖7(b)所示;當t=500 ms時,彈體部分進入空氣中,此時彈頭以及彈體側面與空氣接觸的位置出現較大的應力值,其他位置的應力值較小, 如圖7(c)所示;當t=1 s時,彈體完全進入空氣中,彈體周圍出現了明顯的高應力值,且出現了較大的應力峰值,彈后空泡區域的應力值變得極小,如圖7(d)所示。

圖7 流體應力變化云圖
流體區域的壓力云圖如圖8所示。從圖中可以發現:在彈體運動的初期t=5 ms時,如圖8(a)所示,流場中的壓力主要集中在彈體所在位置且集中于彈體的后部,最大壓力約為1.192 MPa;當t=50 ms時,如圖8(b)所示,流場中壓力集中分布的區域發生變化,最大壓力值減小,約為0.268 7 MPa,當彈體還未出水時,流場中壓力主要集中在彈體所在位置;當t=500 ms時,如圖8(c)所示,彈體部分進入空氣,空氣中的壓力分布較彈體在未出水前變化劇烈,最大壓力位置位于彈體頭部上方。由于空氣阻力較水中阻力小得多,此時流場中最大壓力為25.150 kPa,空氣域出現部分負壓區域;彈體完全出水進入空氣后,流場區域的壓力分布變得較為均勻,僅彈體頭部附近存在最大壓力為1 806 Pa的壓力分布,如圖8(d)所示。

圖8 流體壓力變化云圖
冰水混合區域流場的z方向速度(即豎直向上)變化云圖如圖9所示。
從圖9(a)中可以看出,由于受到彈體運動的影響,彈體周圍的流體開始流動,流速最大的區域位于彈底部,最大速度約為13.27 m/s。當t=50 ms時,如圖9(b)所示,彈體產生的擾動向遠處傳播,導致更大范圍的流體開始流動,彈體周圍的流體速度約為1.94~2.90 m/s。當彈體進入空氣中后,空氣受到彈體的影響開始運動,彈體上方冰體附近的空氣也產生流動現象,如圖9(c)所示。當彈體離開水域后,彈體附近空氣的流速約為1.5 m/s,水域中的流速減小,如圖9(d)所示。

圖9 流體z方向速度變化云圖
在彈體出水的過程中,對彈體的速度進行一定的濾波處理,得到如圖10所示的速度隨時間變化曲線。從圖中可以看出:在運動初始階段,由于受到水的阻力作用,彈體的速度變化劇烈,速度從20 m/s迅速下降;隨后速度出現波動,波動幅度不斷減小,大約在2 ms后速度的波動幅度明顯變緩,逐步接近勻速運動。

圖10 彈體在z方向的速度隨時間的變化曲線
采用LS-DYNA軟件,對冰下發射的導彈出水運動過程進行了數值模擬,獲得如下結果:
①彈體在水下運動過程中,首先在彈底形成空泡,隨著彈體的運動,空泡逐漸擴大,空化現象愈發顯著。
②彈體在運動過程中,水域流場中應力主要集中于彈底以及彈體側面,隨著彈體的運動,流場存在應力的范圍不斷向彈體遠處傳播,且峰值不斷減小。當彈體出水進入空氣中后,彈體周圍的空氣流場區域存在較大的應力,但應力峰值較小。流場應力的存在對彈體的運動產生一定的影響。
③水域流場中,彈體附近流體壓力出現較大的梯度,隨著彈體逐漸離開水域,壓力峰值不斷減小。當彈體進入空氣域時,最大壓力位置位于彈體頭部上方,彈體附近會出現部分負壓區域。彈體完全出水后,水域壓力分布變得較為均勻。
④彈體在出水運動時,彈體周圍流體的流速較大。隨著彈體的運動,流場擾動會逐漸向遠離彈體的地方傳播,導致更大范圍的流體出現流動現象,但是速度峰值逐漸變小。當彈體進入空氣中,空氣受到擾動,靜止的空氣開始運動,彈體與冰體附近的速度值較大。彈體在運動的初期,速度驟降,隨后速度降低程度變緩,最終彈體的速度趨于穩定。
⑤不同的冰體分布方式以及初始彈體相對位置,會造成彈體質心位移不同。