王莎莎，張東東，李新娥，沈大偉，劉奇峰

(中北大學(xué)儀器科學(xué)與動(dòng)態(tài)測(cè)試教育部重點(diǎn)實(shí)驗(yàn)室，山西 太原 030051)

0 引言

外彈道學(xué)是研究彈箭與發(fā)射裝置的作用力脫離之后在空中的運(yùn)動(dòng)規(guī)律和相關(guān)問題的一門學(xué)科[1]。現(xiàn)代戰(zhàn)爭不斷地要求炮彈增加射程，落點(diǎn)精確和毀傷威力，因此，研究外彈道彈丸飛行姿態(tài)和穩(wěn)定性，越來越成為當(dāng)今國際發(fā)展的總趨勢(shì)[2]。

目前，在外彈道的研究上主要使用的方法有風(fēng)洞吹風(fēng)法、工程計(jì)算法、實(shí)驗(yàn)測(cè)量法。比如1 200×1 200 mm2和2 400×2 400 mm2的超音速風(fēng)洞，直徑達(dá)1 000 mm2的高超音速風(fēng)洞為我國新型飛行器的研制提供了充足的技術(shù)保障[3]。但是，風(fēng)洞吹風(fēng)的方法能耗高、周期長、測(cè)試條件復(fù)雜、易受外界干擾、成本高。常用的工程計(jì)算方法主要有二次激波理論、薄翼型理論、二次混合擾動(dòng)、二次激波膨脹、錐形流理論、牛頓理論、小擾動(dòng)線化等多種理論[4]。由于幾何形狀的不同，使用工程計(jì)算的方法得到的數(shù)據(jù)在精度上往往不是很高，對(duì)彈體的某些部分和某些時(shí)段的估算也存在誤差。在彈箭外彈道測(cè)量上，裴東興、岳晗等人提出了使用地磁薄膜線圈磁傳感器測(cè)試彈箭飛行姿態(tài)參數(shù)的方法[5]。沈大偉、邊玉亮等人提出了使用地磁薄膜線圈測(cè)試旋轉(zhuǎn)彈章動(dòng)參數(shù)的方法[6]。其測(cè)試方法是將測(cè)試儀器系統(tǒng)直接放置在實(shí)際運(yùn)動(dòng)著的被測(cè)體內(nèi)或被測(cè)環(huán)境中，使得儀器系統(tǒng)與被測(cè)對(duì)象承受相同的惡劣環(huán)境影響，并實(shí)時(shí)、實(shí)況、準(zhǔn)確地獲取動(dòng)態(tài)參數(shù)，此測(cè)試方法具有體積小、低功耗、穩(wěn)定、抗沖擊、抗高過載等技術(shù)特點(diǎn)。但是此方法一般在測(cè)試時(shí)環(huán)境極其惡劣，并且測(cè)試成本較高。針對(duì)上述問題，提出了基于FLUENT動(dòng)力學(xué)仿真軟件和ADAMS運(yùn)動(dòng)學(xué)仿真軟件相結(jié)合的彈丸外彈道仿真測(cè)試方法。

1 旋轉(zhuǎn)彈外彈道模型

無控旋轉(zhuǎn)彈在空中受到的空氣動(dòng)力不僅影響著彈丸的飛行姿態(tài)，而且影響打擊目標(biāo)的精準(zhǔn)性。因此，建立準(zhǔn)確而又簡單的旋轉(zhuǎn)彈動(dòng)力學(xué)模型，是對(duì)彈丸的外彈道穩(wěn)定性分析和控制系統(tǒng)設(shè)計(jì)的基礎(chǔ)[7]。

1.1 彈丸基本假設(shè)

為了達(dá)到預(yù)期的仿真效果，需要對(duì)彈丸做一定的假設(shè)。這樣的假設(shè)不僅可以將問題變得簡化，同時(shí)也可以將整個(gè)彈丸的運(yùn)動(dòng)規(guī)律展現(xiàn)出來。

1) 彈丸的外形和質(zhì)量相對(duì)于彈軸分布是對(duì)稱的，彈軸為慣性主軸且質(zhì)心位于彈軸上；

2) 彈丸的運(yùn)動(dòng)可以看成全部質(zhì)量集中于質(zhì)心的一個(gè)質(zhì)點(diǎn)的運(yùn)動(dòng)；

3) 攻角很小或者不存在，使得線性關(guān)系成立；

4) 地表為平面，重力加速度不變，方向鉛直向下；

5) 氣象條件為標(biāo)準(zhǔn)條件。

1.2 彈丸受力分析

彈丸在出膛后，不會(huì)再受到火藥的壓力，而是受到空氣的阻力。圖1為高速旋轉(zhuǎn)彈在受到空氣阻力和力矩時(shí)的情況。在圖中可以看出總空氣動(dòng)力FN可分為與速度方向相反的水平分力Rx和豎直速度方向的豎直分力Ry(或由沿彈軸的軸向力RA和垂直于彈軸的法向力Rn合成的).由于俯仰力矩Mz的作用，使得彈的軸向方向遠(yuǎn)離速度矢量方向，彈丸產(chǎn)生一個(gè)攻角δ的變化。

圖1 旋轉(zhuǎn)彈受力示意圖Fig.1 Schematic diagram of the force of the rotating bomb

1.2.1不考慮攻角時(shí)的受力分析

彈丸在空中飛行時(shí)，由于空氣具有黏性，彈表會(huì)附帶一層空氣，和彈丸共同運(yùn)動(dòng)，這種阻力稱為摩阻。氣流達(dá)到一定值后，彈丸尾部的流線會(huì)與彈丸分離，形成旋渦，這種阻力稱為渦阻。在超聲速條件下，彈尾將會(huì)出現(xiàn)大量的旋渦，其彈頭和彈尾由于劇烈壓縮空氣，產(chǎn)生接近于錐形的空氣層，此阻力稱為波阻。由此可知，彈丸在飛行時(shí)，同時(shí)受到摩阻、渦阻和波阻的影響。如式(1)、式(2)所示：

(1)

Cx0(Ma)=cxf+cxb+cxw

(2)

此時(shí)所產(chǎn)生的加速度由式(3)所示，ax體現(xiàn)的是空氣阻力對(duì)彈丸質(zhì)心速度大小和方向的改變[8]。

(3)

式(1)—式(3)中，ν為彈丸相對(duì)空氣的速度，ρ為空氣密度，S為特征面積，cx0(Ma)為阻力系數(shù)，cxf為摩阻系數(shù)，cxb為渦阻系數(shù)，cxw為波阻系數(shù)。

1.2.2考慮攻角時(shí)的受力分析

如圖1所示，當(dāng)攻角存在時(shí)，氣流在彈丸的迎風(fēng)側(cè)和背風(fēng)側(cè)的變化關(guān)于彈軸是不對(duì)稱的。迎風(fēng)側(cè)空氣受到彈體擠壓而風(fēng)壓大，背風(fēng)側(cè)空氣受到拉力而風(fēng)壓小，這種壓差隨著速度和攻角的增大導(dǎo)致總空氣動(dòng)力顯著增大，并且它的作用方向既不沿著彈軸或速度的反方向也不通過彈丸瞬時(shí)質(zhì)心，而是在攻角平面內(nèi)偏在彈丸背風(fēng)面與彈軸速度線的一側(cè)[9]。由圖1可以得到Rx、Ry、RA、Rn的關(guān)系為式(4)所示：

(4)

將式(4)換成氣動(dòng)力系數(shù)即為式(5)所示：

(5)

式(5)中，cx為阻力系數(shù)，cy為升力系數(shù)，cA為軸向力系數(shù)，cn為法向力系數(shù)。

總空氣動(dòng)力FN對(duì)彈丸質(zhì)心作用的力矩為式(6)所示：

(6)

將總空氣動(dòng)力FN以作用于壓心的二分力Rx、Ry代替，可得出俯仰力矩Mz和升力Ry、阻力Rx和壓心P到質(zhì)心O之間的距離h之間的關(guān)系。如式(7)、式(8)所示：

Mz=Rycosδ+Rxsinδ

(7)

mz=(cycosδ+cxsinδ)h/l

(8)

從式子中可以看出，其它條件不變的話，俯仰力矩是Ma和攻角δ的函數(shù)，l為參考長度，mz為靜力矩系數(shù)。

1.2.3壓力中心

由圖2可以看出壓力中心的位置處于彈丸質(zhì)心和頭部之間的連線上。實(shí)驗(yàn)證明，彈丸的壓心位置不僅隨馬赫數(shù)變化而變化，也隨攻角δ大小而改變，根據(jù)旋轉(zhuǎn)彈丸在風(fēng)洞中吹風(fēng)的實(shí)驗(yàn)結(jié)果，壓心隨著馬赫數(shù)的增大向彈底移動(dòng)，使得實(shí)際的壓心、質(zhì)心距h減小[10]。

圖2 彈丸各部尺寸圖Fig.2 Dimensions of the various parts of the projectile

對(duì)于旋轉(zhuǎn)彈，壓力距離質(zhì)心的距離可用高巴爾公式進(jìn)行估算，如式(9)所示：

h=h0+0.37ht-0.16d

(9)

式(9)中，h0表示頭部底端至質(zhì)心的距離，ht表示彈丸頭部的長度。

除上節(jié)所介紹的彈丸受力外，彈丸在外彈道飛行過程中，還受到重力的影響。標(biāo)準(zhǔn)重力加速度為：北緯45°，零海拔：gNo=9.806 65 m/s2，g=gNof(y,Λ)。y表示高度，Λ表示緯度。

2 外彈道仿真

2.1 彈丸外流場(chǎng)仿真

FLUENT是一個(gè)用于模擬和分析復(fù)雜幾何區(qū)域內(nèi)的流體流動(dòng)和傳熱現(xiàn)象專用軟件，具有廣泛的物理模型，能夠快速、準(zhǔn)確的獲得CFD(流體動(dòng)力學(xué))分析結(jié)果。在仿真的過程中，必要時(shí)還要對(duì)網(wǎng)格進(jìn)行細(xì)化，改變數(shù)值和物理模型，使得仿真條件更加接近真實(shí)工況，得到更為精確的結(jié)果。

首先對(duì)無控旋轉(zhuǎn)彈建立一個(gè)長方體的計(jì)算域，如圖3所示。在GAMBIT中進(jìn)行彈丸模型的網(wǎng)格劃分，如圖4所示。計(jì)算出不同工況下(不同馬赫數(shù)和不同攻角)的外部繞流場(chǎng)。在獲得流場(chǎng)的壓力分布之后，F(xiàn)LUENT會(huì)自動(dòng)對(duì)彈丸的表面進(jìn)行積分，從而求得彈體和各個(gè)位置表面受到的空氣動(dòng)力。為了簡化分析，只觀察空氣阻力和俯仰力矩的數(shù)值。除此之外，彈丸運(yùn)動(dòng)速度通常都超過了聲速，甚至達(dá)到幾倍音速，因此彈丸周圍的流場(chǎng)為超音速氣流，密度的變化不可忽略。因此，將仿真定義為三維彈丸的外部超聲速可壓縮繞流流場(chǎng)的數(shù)值分析。

圖3 計(jì)算域的建立Fig.3 Establishment of the calculation domain

圖4 彈丸的網(wǎng)格劃分Fig.4 Meshing of the projectile

由于繞流場(chǎng)的計(jì)算要求在彈丸附近布置的網(wǎng)格足夠密，才能準(zhǔn)確的獲得所需要的參數(shù)，同時(shí)計(jì)算域又要足夠的大，以適應(yīng)遠(yuǎn)場(chǎng)邊界條件，避免流場(chǎng)中壓力波等邊界反射引起數(shù)值計(jì)算的不穩(wěn)定[11]。在計(jì)算域的中間位置再設(shè)置一個(gè)小圓柱體，將里邊的網(wǎng)格進(jìn)行加密處理，外邊的網(wǎng)格劃分上較為稀疏，在保證了計(jì)算精度的同時(shí)，又加快了仿真進(jìn)度。

將劃分好網(wǎng)格的模型文件導(dǎo)入到FLUENT中，本次仿真選用密度基耦合顯示求解器，S-A湍流模型，在流體物理屬性定義上，選擇了空氣，滿足薩蘭德定律計(jì)算粘性，將定壓比熱和熱傳導(dǎo)率都設(shè)置為常數(shù)。在密度基的求解過程中，由于會(huì)涉及接下來的壓力輸入和數(shù)據(jù)對(duì)比，為了與邊界條件的設(shè)定值表述習(xí)慣相一致，將參考?jí)毫χ昧?。設(shè)置的邊界條件參數(shù)如表1所示。

表1 邊界條件參數(shù)設(shè)定值

在表1中，由于在模型的建立過程中，坐標(biāo)系X軸正向的選擇是彈尾朝向彈頭，與來流氣體相對(duì)彈丸的運(yùn)動(dòng)方向剛好相反，因此橫坐標(biāo)添加負(fù)號(hào)。

完成上述的所有步驟后，對(duì)控制參數(shù)進(jìn)行相關(guān)的設(shè)定，完成整個(gè)仿真。得到1.5Ma和4°攻角下的外流場(chǎng)仿真，進(jìn)而得出速度云圖，如圖5所示。

從云圖中可以看到，高速飛行的彈丸，其彈表和彈尾受到摩阻、渦阻的影響(彈尾波顯示)。由于彈丸周圍的流場(chǎng)為超音速氣流，彈頭和彈尾有相似于錐形的空氣層。這與之前的空氣阻力分析相一致。并且由于攻角的存在，使得彈丸上下兩表面產(chǎn)生了明顯的不對(duì)稱，這與分析和實(shí)際測(cè)試相一致。

圖5 速度云圖Fig.5 Speed cloud

在FLUENT軟件建立的模型基礎(chǔ)上，通過修改邊界條件，分別計(jì)算出不同攻角和馬赫數(shù)下彈丸的空氣阻力和俯仰力矩。如表2、表3所示。

從表2、表3中可以看出，空氣阻力的大小受到攻角的影響不大，主要受馬赫數(shù)值的影響；而俯仰力矩同時(shí)受到攻角和馬赫數(shù)的影響。

2.2 彈丸運(yùn)動(dòng)學(xué)仿真

在FLUENT仿真的基礎(chǔ)上，利用ADAMS運(yùn)動(dòng)學(xué)仿真軟件的仿真過程可視化特性，可以對(duì)彈的整個(gè)外彈道過程進(jìn)行跟蹤觀察和測(cè)量，并得出相應(yīng)的運(yùn)動(dòng)曲線圖。首先建立彈體的三維模型，將模型導(dǎo)入ADAMS中得到圖6。

表2 不同攻角和馬赫數(shù)下的空氣阻力值

表3 不同攻角和馬赫數(shù)下的俯仰力矩

圖6 ADAMS中的模型圖Fig.6 Model diagram in ADAMS

由于本次實(shí)驗(yàn)所用的測(cè)試裝置為彈載測(cè)試儀，因此在仿真時(shí)也需加入測(cè)試裝置。為了在不影響整個(gè)彈丸的質(zhì)量和質(zhì)心的基礎(chǔ)上，需要將測(cè)試儀做到小型化。因此在仿真時(shí)用一個(gè)小圓柱體模型代替，如圖6彈頭小方塊顯示。

在ADAMS中對(duì)彈丸進(jìn)行參數(shù)設(shè)定，總質(zhì)量為5.86 kg，質(zhì)心位置為118.5 mm，質(zhì)心坐標(biāo)為(0,0,0),初速為980 m/s，轉(zhuǎn)速為2 513.27 rad/s,射角為21°。重力加速度g恒定為9.806 65 m/s2，方向鉛直向下，得到設(shè)置完參數(shù)的模型圖7。

圖7 設(shè)置完參數(shù)的模型Fig.7 Model for setting parameters

2.3 彈丸仿真曲線分析

在ADAMS將設(shè)置完參數(shù)的模型施加之前分析的載荷后進(jìn)行仿真，就可以得到彈丸在外彈道飛行過程中的速度和加速度仿真曲線，如圖8和圖9所示。

圖8 加速度仿真曲線Fig.8 Acceleration simulation curve

圖9 速度仿真曲線Fig.9 Speed simulation curve

由仿真曲線可以看出，飛行階段的加速度曲線呈現(xiàn)負(fù)值，這是由于施加載荷時(shí)，將指向彈丸頭部方向定為正向，因此加速度方向與速度方向相反。外彈道加速度值由最開始的-117.10 m/s2快速下降到最高點(diǎn)的-8.18 m/s2后又平緩的下降至29.6 s的-3.263 m/s2,后期增加到落地前的-3.498 m/s2，彈丸外彈道的飛行速度由最開始的980 m/s下降到26.4 s處248.2 m/s之后速度小幅度的增加到267.2 m/s，這是由于彈在飛行的后期，飛行處于平穩(wěn)階段，速度降低，所受的飛行阻力降低，重力彈在飛行爬升階段的阻力變成了下降階段的動(dòng)力，使速度有所增加。

3 彈丸實(shí)測(cè)曲線分析

靶場(chǎng)進(jìn)行發(fā)射實(shí)驗(yàn)時(shí)，將彈載測(cè)試儀置于實(shí)際彈丸中，彈載測(cè)試儀底部中心安裝加速度傳感器，用于測(cè)量彈丸在外彈道飛行過程中所受的空氣阻力引起的加速度變化[12]。通過加速度計(jì)測(cè)得加速度實(shí)測(cè)曲線，如圖10所示。對(duì)所測(cè)的加速度實(shí)測(cè)曲線進(jìn)行濾波后得到加速度濾波曲線，如圖11所示。通過對(duì)加速度曲線積分，得到速度曲線，如圖12所示。

圖10 加速度實(shí)測(cè)曲線Fig.10 Acceleration measurement curve

圖11 加速度濾波曲線Fig.11 Acceleration filter curve

由圖10可以看到，實(shí)測(cè)曲線中存在很多的雜波，對(duì)其進(jìn)行濾波處理后可以看到加速度值由最開始的-105.7 m/s2快速下降到最高點(diǎn)的-10.2 m/s2后又平緩的下降至29.6 s的-7.436 m/s2,最后變?yōu)?5.23 m/s2。

將加速度實(shí)測(cè)曲線進(jìn)行積分后，得到速度曲線。由圖12可以看出，彈丸速度由最開始的981.3 m/s下降到26.4 s處262.3 m/s，最后變?yōu)?39.1 m/s。

圖12 速度曲線Fig.12 Speed curve

由實(shí)測(cè)曲線和仿真曲線進(jìn)行對(duì)比，發(fā)現(xiàn)彈在飛行爬升階段持續(xù)時(shí)間低于下降的階段，整個(gè)外彈道的持續(xù)階段為37.616 s。彈丸在16.2 s飛到彈道最高點(diǎn)，由于加速度計(jì)本身的原因，使得在彈丸飛行過程中的爬升段在0～16.2 s內(nèi)，加速度實(shí)測(cè)值小于實(shí)際仿真值，在彈丸飛行過程中的下降段16.2～37.616 s內(nèi)，加速度實(shí)測(cè)值大于實(shí)際仿真值，因此導(dǎo)致加速度的測(cè)試曲線與實(shí)際仿真曲線不完全相同，進(jìn)而導(dǎo)致速度的測(cè)試曲線與實(shí)際仿真曲線不完全重合。

由實(shí)測(cè)加速度、速度曲線和仿真得到的加速度、速度曲線可以看出，彈丸在實(shí)際飛行中的加速度、速度的變化與仿真得到的變化規(guī)律基本保持一致，都是前期曲線比較陡，后期比較平緩。

4 結(jié)論

本文提出了無控旋轉(zhuǎn)彈外彈道仿真方法，該方法在建立彈丸受到的空氣力和力矩的基礎(chǔ)上，運(yùn)用FLUENT動(dòng)力學(xué)仿真軟件對(duì)彈丸飛行過程進(jìn)行仿真。得到了不同馬赫數(shù)和不同攻角下彈丸受到的空氣阻力和俯仰力矩，將這些參數(shù)帶入到ADAMS運(yùn)動(dòng)學(xué)仿真軟件中，在一定的初始條件下，完成對(duì)無控旋轉(zhuǎn)彈的外彈道過程的仿真，得到速度和加速度仿真曲線。并與實(shí)測(cè)曲線進(jìn)行對(duì)比，發(fā)現(xiàn)仿真曲線和實(shí)測(cè)曲線基本一致。