陳時(shí)通,于 勇

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

子母彈是現(xiàn)代武器中很重要的一種大規(guī)模殺傷性武器,它在飛行的過程中從母彈中彈射出多發(fā)子彈,造成大規(guī)模、大面積的毀傷[1]。子母彈的作戰(zhàn)效能由多個(gè)子彈的散布效果決定。子母彈研制過程中面臨的重要問題之一是使多個(gè)子彈正常分離,達(dá)到預(yù)期的散布效果。而要解決這一問題,首先要保證子彈分離過程中姿態(tài)穩(wěn)定。

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

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

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

1 數(shù)值計(jì)算方法

1.1 流場(chǎng)求解

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

(1)

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

1.2 子彈運(yùn)動(dòng)方程

子母彈運(yùn)動(dòng)6DOF方程為

(2)

(3)

(4)

(5)

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

1.3 嵌套網(wǎng)格技術(shù)

隨著工業(yè)技術(shù)的發(fā)展,工程項(xiàng)目模型更加復(fù)雜,網(wǎng)格質(zhì)量的好壞成為了工程模擬能否進(jìn)行的關(guān)鍵。工程項(xiàng)目中,網(wǎng)格的生成與修改往往要占據(jù)80%以上的時(shí)間,項(xiàng)目模型越復(fù)雜,網(wǎng)格生成越困難,需要的時(shí)間也越多。而且,在模擬的過程中,經(jīng)常因?yàn)榫W(wǎng)格質(zhì)量差導(dǎo)致模擬失敗,最終工程項(xiàng)目只能暫停。因此,計(jì)算網(wǎng)格生成技術(shù)的發(fā)展成為了解決問題的關(guān)鍵[16]。

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

1.4 嵌套網(wǎng)格技術(shù)驗(yàn)證

機(jī)彈分離問題是多體分離中的經(jīng)典問題,本節(jié)將對(duì)機(jī)彈分離問題進(jìn)行模擬,并與風(fēng)洞試驗(yàn)進(jìn)行對(duì)比,以此驗(yàn)證計(jì)算方法的正確性。算例來(lái)流靜壓p=20 657 Pa,靜溫T=216.65 K,Ma=1.2,攻角α=0°,仿真時(shí)間步長(zhǎng)Δt=0.001 s。表1為機(jī)彈分離問題參數(shù)表。

表1 機(jī)彈分離問題參數(shù)

表1中,ms,Ls,Ds分別為彈體的質(zhì)量、長(zhǎng)度和直徑;k為彈翼數(shù);ds為質(zhì)心位置;Ix,Iy,Iz為分別x,y,z方向轉(zhuǎn)動(dòng)慣量;FFE,FAE為彈射力。詳細(xì)模型尺寸和實(shí)驗(yàn)結(jié)果可參考文獻(xiàn)[14]。

圖1為仿真模型圖。

圖1 機(jī)彈分離模型

圖2、圖3分別為導(dǎo)彈質(zhì)心位移與姿態(tài)角隨時(shí)間的變化規(guī)律,圖中,D為質(zhì)心位移,λ為姿態(tài)角(姿態(tài)角λ代表了滾轉(zhuǎn)角γ、偏航角ψ和俯仰角θ)。圖中導(dǎo)彈滾轉(zhuǎn)角與實(shí)驗(yàn)值偏差較大,因?yàn)閷?dǎo)彈軸向轉(zhuǎn)動(dòng)慣量較小,模擬較困難,需增加導(dǎo)彈周向的網(wǎng)格節(jié)點(diǎn),但會(huì)極大地降低計(jì)算效率?？傮w來(lái)說,仿真數(shù)據(jù)與實(shí)驗(yàn)數(shù)據(jù)符合較好。機(jī)彈分離的經(jīng)典算例模擬結(jié)果證明了計(jì)算方法的有效性。

圖2 導(dǎo)彈質(zhì)心位移隨時(shí)間的變化

圖3 導(dǎo)彈的姿態(tài)角隨時(shí)間的變化

2 子母彈分離過程的數(shù)值模擬

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

2.1 子母彈模型

子母彈結(jié)構(gòu)設(shè)計(jì)分為2種方案:彈艙靠前設(shè)計(jì)(工況A),彈艙靠后設(shè)計(jì)(工況B)。

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

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

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

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

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

建立2個(gè)坐標(biāo)系:母彈坐標(biāo)系為定坐標(biāo)系,以母彈質(zhì)心為原點(diǎn),以母彈彈軸為X1軸,彈尾指向彈頭方向?yàn)閄1軸正方向。Y1軸在母彈模型縱向?qū)ΨQ面且與X1軸垂直,指向上與重力方向相反為Y1軸正方向。Z1軸由右手法則確定。子彈坐標(biāo)系為動(dòng)坐標(biāo)系,以子彈質(zhì)心為原點(diǎn),子彈彈軸為X2軸,彈尾指向彈頭方向?yàn)閄2軸正方向。Y2軸始終垂直于X2軸且指向上為正方向,初始時(shí)刻Y2軸方向與重力方向相反。Z2軸由右手法則確定。4枚子彈分別位于母彈的不同周向位置,圖8為子母彈裝配圖。

圖8 子母彈裝配圖

2.2 邊界條件

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

圖9 計(jì)算域形狀

2.3 計(jì)算網(wǎng)格

圖10是定常狀態(tài)和非定常狀態(tài)下某一時(shí)刻的Overset網(wǎng)格。

圖10 子母彈Overset網(wǎng)格

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

3 結(jié)果分析

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

3.1 子母彈壓力場(chǎng)

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

圖11 子母彈壓力場(chǎng)

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

3.2 彈艙位置對(duì)子彈姿態(tài)的影響

3.2.1 子彈1、子彈4的姿態(tài)變化

圖12為子彈1的姿態(tài)角隨時(shí)間的變化曲線。從圖12可知,工況A時(shí)子彈1的姿態(tài)角變化較工況B更大。因?yàn)閺椗撐恢每壳?接近母彈激波,母彈激波附近壓力較大,對(duì)子彈1的姿態(tài)影響較大。子彈1的近激波面壓力較大,使子彈1的姿態(tài)發(fā)生偏轉(zhuǎn)。在母彈激波作用下,子彈1的偏航角正向增大到0.7°,俯仰角負(fù)向增大到-1.7°。工況B時(shí),由于子彈1位置遠(yuǎn)離母彈激波,在子彈2激波和空氣動(dòng)力的作用下,子彈1的偏航角和俯仰角負(fù)向增大到-0.5°左右,小于工況A。

圖12 子彈1的姿態(tài)變化

圖13為子彈4的姿態(tài)角隨時(shí)間變化曲線。從圖13中能夠得出與子彈1類似的結(jié)論。

3.2.2 子彈2、子彈3的姿態(tài)變化

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

圖13 子彈4的姿態(tài)變化

圖14 子彈2的姿態(tài)變化

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

圖15 子彈3的姿態(tài)變化

圖16(a)、16(b)為工況A子母彈在t=0.01 s和t=0.02 s時(shí)的壓力場(chǎng),截面經(jīng)過子彈質(zhì)心,方向?yàn)閄1軸方向。從圖中可以看出,子母彈壓力場(chǎng)對(duì)稱,因此,子彈1、子彈4和子彈2、子彈3的姿態(tài)運(yùn)動(dòng)規(guī)律對(duì)稱。

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

圖16 工況A、工況B不同時(shí)刻子母彈壓力場(chǎng)

4 結(jié)論

采用嵌套網(wǎng)格的方法對(duì)子母彈的分離進(jìn)行仿真與分析,由本文分析可知,采用嵌套網(wǎng)格技術(shù)和k-ω湍流模型,同時(shí)耦合6DOF方程,能夠很好地模擬子母彈分離的姿態(tài)變化過程及其干擾流場(chǎng)。

子母彈分離過程中,若彈艙位置在前部,子彈分離過程始終出現(xiàn)子彈激波與母彈激波的相交,子彈受到母彈激波影響姿態(tài)變化比彈艙位置在后部更大。若彈艙位置在后部,子母彈分離初期出現(xiàn)子彈激波的相交,中間位置對(duì)稱分布的子彈在拋撒過程中受到子彈激波和不對(duì)稱壓力場(chǎng)的作用,姿態(tài)會(huì)發(fā)生不對(duì)稱的變化。無(wú)論彈艙位于前部還是后部,上下位置對(duì)稱分布的子彈的壓力場(chǎng)是對(duì)稱的,運(yùn)動(dòng)姿態(tài)也對(duì)稱。