超聲速流場中子母彈分離與子彈姿態變化規律的數值模擬

陳時通,于 勇

(北京理工大學 宇航學院,北京 100081)

子母彈是現代武器中很重要的一種大規模殺傷性武器,它在飛行的過程中從母彈中彈射出多發子彈,造成大規模、大面積的毀傷[1]。子母彈的作戰效能由多個子彈的散布效果決定。子母彈研制過程中面臨的重要問題之一是使多個子彈正常分離,達到預期的散布效果。而要解決這一問題,首先要保證子彈分離過程中姿態穩定。

國外關于子母彈的研究資料較少。Edge等[2]自編CFD程序,采用嵌套網格,研究了子母彈拋撒過程中的激波相互作用。Panneerselvam等[3]通過風洞實驗研究了子母彈在彈艙內和彈艙外的法向力系數、力矩系數和壓心位置的變化。Corder等[4]使用自編程序研究了帶彈翼和不帶彈翼的子彈在不同攻角下的法向力系數和俯仰力矩系數的變化。Dietz[5]使用軟件Overflow模擬了24枚子彈拋撒前的馬赫數、壓力分布以及拋撒后各子彈的位移和姿態變化。 Deep等[6]通過自編程序研究了有、無彈艙的情況下,子彈在不同馬赫數及不同攻角情況下的法向力系數及俯仰力矩系數并與實驗對比,吻合良好。

國內對子母彈的研究起步較晚。鄒德坤等[7]使用CFX模擬了子彈在不同攻角、不同馬赫數下的升力系數和阻力系數。王金龍等[8]將氣囊實驗的結果運用于FLUENT數值仿真,研究了在燃氣干擾下不同來流馬赫數對子母彈的流場結構和子彈的運動特性的影響。也有學者研究母彈的殼片在分離過程的運動特性。蔣增輝等[9]使用風洞投放模型試驗技術研究了母彈處于大迎角(25°)狀態下殼片的分離運動規律。王巍[10]利用彈簧近似和網格重構相結合的非結構動網格技術,耦合求解N-S方程及彈道方程,模擬了不同的殼片分離速度和角速度的拋殼過程。陶如意等利用AUSM+格式求解采用k-ωSST湍流模型的雷諾平均N-S方程,研究了子母彈時序拋撒過程中的子彈與子彈和子彈與母彈間的激波相互作用[11];通過風洞實驗研究了子彈在母彈流場中不同位置、不同攻角下的阻力系數、升力系數和力矩系數的變化[12]。

實驗研究可以給出子母彈拋撒過程中的氣動參數及流場狀況,但實驗成本大,周期長。與實驗相比,數值模擬成本低,周期短。已有大量的學者通過數值模擬研究子母彈拋撒過程,研究多采用動網格技術,使用非結構網格進行網格重構。但是,所帶來的問題是網格數量多,網格重構計算量大,而且重構的網格質量會逐漸降低,嚴重時甚至使計算中斷。因此,本文在計算方法上試圖使用嵌套動網格技術,采用結構化網格。在網格邊界運動時,不需要重新生成網格,從而提高網格質量,減少計算時間。同時耦合動網格技術和6DOF方程計算分離過程中的子彈軌跡與姿態,利用這種非定常的模擬方法研究了母彈上不同彈倉位置對子彈拋撒后的4枚子彈姿態的影響。

1 數值計算方法

1.1 流場求解

雷諾平均N-S方程在笛卡爾坐標系下的描述為

(1)

式中:Ω為控制體體積;S為控制體表面積;n為表面S的外法向單位矢量;dV為體積分微元;dS為面積分微元;Q為守恒變量,F(Q)為對流項,H(Q)為黏性通量。

1.2 子彈運動方程

子母彈運動6DOF方程為

(2)

(3)

(4)

(5)

本文使用的是Fluent18.0,采取k-ω湍流模型,耦合6DOF方程(弱耦合),使用的求解方法是隱式的Roe-FDS方法,離散方程是二階迎風格式,選取密度求解器。

1.3 嵌套網格技術

隨著工業技術的發展,工程項目模型更加復雜,網格質量的好壞成為了工程模擬能否進行的關鍵。工程項目中,網格的生成與修改往往要占據80%以上的時間,項目模型越復雜,網格生成越困難,需要的時間也越多。而且,在模擬的過程中,經常因為網格質量差導致模擬失敗,最終工程項目只能暫停。因此,計算網格生成技術的發展成為了解決問題的關鍵[16]。

嵌套網格技術(Fluent中稱為Overset)能夠很好地解決模型復雜、網格質量不高的問題,同時還具有部件網格生成容易、網格數量少,網格質量高且高可移植性的特點。嵌套網格技術將模型拆解成一個背景網格和多個簡單的組件網格。每套網格需單獨劃分結構化網格,并在Fluent中進行插值,網格重疊的邊界處采用三線性差值的方式傳遞數據。組件網格運動的過程中,每運動一次都需要重新生成新的嵌套網格,保證生成網格具有高質量。同時,背景網格和組件網格單獨劃分,可根據工程需要隨時更改組件形狀,提高網格生成效率[13]。

1.4 嵌套網格技術驗證

機彈分離問題是多體分離中的經典問題,本節將對機彈分離問題進行模擬,并與風洞試驗進行對比,以此驗證計算方法的正確性。算例來流靜壓p=20 657 Pa,靜溫T=216.65 K,Ma=1.2,攻角α=0°,仿真時間步長Δt=0.001 s。表1為機彈分離問題參數表。

表1 機彈分離問題參數

表1中,ms,Ls,Ds分別為彈體的質量、長度和直徑;k為彈翼數;ds為質心位置;Ix,Iy,Iz為分別x,y,z方向轉動慣量;FFE,FAE為彈射力。詳細模型尺寸和實驗結果可參考文獻[14]。

圖1為仿真模型圖。

圖1 機彈分離模型

圖2、圖3分別為導彈質心位移與姿態角隨時間的變化規律,圖中,D為質心位移,λ為姿態角(姿態角λ代表了滾轉角γ、偏航角ψ和俯仰角θ)。圖中導彈滾轉角與實驗值偏差較大,因為導彈軸向轉動慣量較小,模擬較困難,需增加導彈周向的網格節點,但會極大地降低計算效率。總體來說,仿真數據與實驗數據符合較好。機彈分離的經典算例模擬結果證明了計算方法的有效性。

圖2 導彈質心位移隨時間的變化

圖3 導彈的姿態角隨時間的變化

2 子母彈分離過程的數值模擬

本文計算模型為1枚母彈和8枚子彈裝配,由于子彈對稱分布,因此模型只考慮右邊部分。本文采用嵌套網格對子母彈分離過程進行數值模擬,子母彈工作高度為8 km,母彈攻角為0°,且無旋轉,飛行馬赫數為2.43,靜壓為35 652 Pa,靜溫為236.2 K,子彈彈射速度為50 m/s(彈射沿徑向方向均勻作用在子彈質心上),仿真時間為0.03 s。初始狀態下彈艙已打開,子彈暴露在氣流中。非定常計算時,給定子彈初始彈射速度,子彈做無控運動。

2.1 子母彈模型

子母彈結構設計分為2種方案:彈艙靠前設計(工況A),彈艙靠后設計(工況B)。

彈艙靠前設計及母彈尺寸如圖4所示,子彈布局如圖5所示,子彈尺寸如圖6所示。彈艙靠后設計及母彈尺寸如圖7所示。

圖4 彈艙靠前設計母彈尺寸(單位:mm)

圖5 子彈布局(單位:mm)

圖6 子彈尺寸(單位:mm)

圖7 彈艙靠后設計及母彈尺寸(單位:mm)

建立2個坐標系:母彈坐標系為定坐標系,以母彈質心為原點,以母彈彈軸為X1軸,彈尾指向彈頭方向為X1軸正方向。Y1軸在母彈模型縱向對稱面且與X1軸垂直,指向上與重力方向相反為Y1軸正方向。Z1軸由右手法則確定。子彈坐標系為動坐標系,以子彈質心為原點,子彈彈軸為X2軸,彈尾指向彈頭方向為X2軸正方向。Y2軸始終垂直于X2軸且指向上為正方向,初始時刻Y2軸方向與重力方向相反。Z2軸由右手法則確定。4枚子彈分別位于母彈的不同周向位置,圖8為子母彈裝配圖。

圖8 子母彈裝配圖

2.2 邊界條件

數值計算中,母彈速度為超聲速,外部空氣邊界條件為壓力遠場,母彈設置為壁面,外部空氣邊界條件設置為Overset,子彈設置為壁面。由于母彈的彈頭形狀是小圓弧,彈尾是圓柱形,因此,計算域形狀與母彈形狀類似,計算域頭部為半圓形,尾部為圓柱形,計算域半徑為2 m,長度為11 m。計算域形狀如圖9所示。

圖9 計算域形狀

2.3 計算網格

圖10是定常狀態和非定常狀態下某一時刻的Overset網格。

圖10 子母彈Overset網格

由于計算采用的是嵌套網格,因此母彈和子彈需分別劃分結構化網格。子彈和母彈靠近壁面的第1層網格高度為1 mm,子彈部件網格數量為21萬,母彈網格數量為170萬。母彈周圍的網格為背景網格,子彈周圍的網格為組件網格。在Fluent中選定背景網格和組件網格進行插值,重疊區域插值方式為最小二乘插值。計算網格壁面y+約為200,時間步長為10-4s,殘差收斂達到10-3。在非定常狀態下,每一個時間步組件網格要更新移動,嵌套網格也會進行重構,保證每一次迭代的網格質量。

3 結果分析

本文分別對子彈彈艙在前部(工況A)或后部(工況B)進行數值模擬。

3.1 子母彈壓力場

子彈彈射出艙之前,需要計算定常狀態下子彈在母彈艙內的繞流流場,并以此為初始流場計算子母彈的分離過程。圖11(a)為定常狀態下子母彈的壓力場,圖中,p為靜壓。從圖中可以看出,在定常狀態下,母彈彈頭處壓力較高,頭部形成馬赫波,從母彈頭部到尾部,壓力降低,在母彈彈艙后部位置形成馬赫波。子彈尾部位于彈艙馬赫波處,壓力較高。

圖11 子母彈壓力場

圖11(b)為工況A非定常狀態下子彈流場壓力分布,仿真時間t=0.03 s。圖中截面為子彈2、子彈3的縱向截面,下同。從圖中可以看出,子彈彈頭處壓力較高,形成激波且子彈形成的激波與母彈形成的激波相交。圖11(c)為工況B非定常狀態下子彈流場壓力分布,仿真時間t=0.03 s。從圖中可以看出,子彈彈頭處也形成激波,但由于彈艙位置遠離母彈激波,沒有出現與母彈激波相交的狀況。圖11(d)為工況B非定常狀態下子彈流場壓力分布,仿真時間t=0.01 s。從圖中可以看出,由于仿真時間較短,子彈位移較小,因此發生了子彈間的激波相交。子母彈分離的壓力場表明,若彈艙位置在后部,子母彈分離初期子彈間發生激波相交;子母彈分離后期,由于子彈間距較大,沒有發生激波相交。

3.2 彈艙位置對子彈姿態的影響

3.2.1 子彈1、子彈4的姿態變化

圖12為子彈1的姿態角隨時間的變化曲線。從圖12可知,工況A時子彈1的姿態角變化較工況B更大。因為彈艙位置靠前,接近母彈激波,母彈激波附近壓力較大,對子彈1的姿態影響較大。子彈1的近激波面壓力較大,使子彈1的姿態發生偏轉。在母彈激波作用下,子彈1的偏航角正向增大到0.7°,俯仰角負向增大到-1.7°。工況B時,由于子彈1位置遠離母彈激波,在子彈2激波和空氣動力的作用下,子彈1的偏航角和俯仰角負向增大到-0.5°左右,小于工況A。

圖12 子彈1的姿態變化

圖13為子彈4的姿態角隨時間變化曲線。從圖13中能夠得出與子彈1類似的結論。

3.2.2 子彈2、子彈3的姿態變化

圖14為子彈2的姿態角隨時間變化曲線。子彈2的位移變化規律與子彈1相同。從圖14可知,工況A時母彈激波對子彈2的姿態影響較大。因此,在母彈激波的作用下,子彈2的偏航角正向增大到1.7°,俯仰角負向增大到-0.7°。工況B時,子彈2受到子彈1、子彈3的激波和空氣動力的作用,俯仰角正向增大到1.8°。

圖13 子彈4的姿態變化

圖14 子彈2的姿態變化

圖15為子彈3的姿態角隨時間變化曲線。子彈3的位移變化規律與子彈1相同。從圖15中可以看出,在工況A時子彈3變化規律與子彈2相同。在工況B時子彈3偏航角與俯仰角變化都較大,與子彈2不同。工況A時,母彈激波對子彈3的姿態影響較大,使子彈3的偏航角正向增大到1.7°,俯仰角正向增大到0.7°。工況B時,子彈3受到子彈2和子彈4的激波的相交作用和空氣動力作用,導致子彈3的俯仰角正向增大到1.8°,偏航角負向增大到-1.8°。

圖15 子彈3的姿態變化

圖16(a)、16(b)為工況A子母彈在t=0.01 s和t=0.02 s時的壓力場,截面經過子彈質心,方向為X1軸方向。從圖中可以看出,子母彈壓力場對稱,因此,子彈1、子彈4和子彈2、子彈3的姿態運動規律對稱。

圖16(c)、16(d)為工況B子母彈在t=0.01 s和t=0.02 s時的壓力場。從圖中可以看出,子彈1、子彈4的壓力場軸對稱,子彈2、子彈3的壓力場不是軸對稱。子彈2、子彈3的背激波面的低壓區方向不同,造成了子彈2、子彈3的姿態變化不同。

圖16 工況A、工況B不同時刻子母彈壓力場

4 結論

采用嵌套網格的方法對子母彈的分離進行仿真與分析,由本文分析可知,采用嵌套網格技術和k-ω湍流模型,同時耦合6DOF方程,能夠很好地模擬子母彈分離的姿態變化過程及其干擾流場。

子母彈分離過程中,若彈艙位置在前部,子彈分離過程始終出現子彈激波與母彈激波的相交,子彈受到母彈激波影響姿態變化比彈艙位置在后部更大。若彈艙位置在后部,子母彈分離初期出現子彈激波的相交,中間位置對稱分布的子彈在拋撒過程中受到子彈激波和不對稱壓力場的作用,姿態會發生不對稱的變化。無論彈艙位于前部還是后部,上下位置對稱分布的子彈的壓力場是對稱的,運動姿態也對稱。