基于歐拉法的自尋的炮射導(dǎo)彈彈道建模與仿真

白 卓，趙河明，楊晉偉

（1.中北大學(xué)機(jī)電工程學(xué)院，太原030051；2.淮海工業(yè)集團(tuán)，山西長(zhǎng)治046012）

常規(guī)武器發(fā)射平臺(tái)在使用和維護(hù)不變的情況下，遇到要求精確打擊時(shí)，可直接使用制導(dǎo)彈藥。這種擴(kuò)展使常規(guī)武器兼?zhèn)淞顺Ｒ?guī)武器的傳統(tǒng)功能及過(guò)去只有制導(dǎo)武器才具有的精確打擊功能。從20 世紀(jì)的幾次局部戰(zhàn)爭(zhēng)來(lái)看，精確制導(dǎo)武器在現(xiàn)代戰(zhàn)爭(zhēng)中起著決定性的作用。據(jù)有關(guān)資料統(tǒng)計(jì)，過(guò)去平均使用250 發(fā)155 mm 的常規(guī)彈藥，只可能擊毀一輛坦克，現(xiàn)在使用精確制導(dǎo)技術(shù)的炮彈僅需1～2 發(fā)即可，其效能提高了125～259倍。現(xiàn)在世界各國(guó)都在爭(zhēng)相研制的炮射導(dǎo)彈就是一種精確制導(dǎo)武器[1-7]。裝有精確制導(dǎo)系統(tǒng)的炮射導(dǎo)彈的出現(xiàn)使現(xiàn)代火炮進(jìn)入了一個(gè)新的時(shí)代，實(shí)現(xiàn)了彈炮集合，由無(wú)控向有控的轉(zhuǎn)變[8]。相對(duì)于機(jī)載導(dǎo)彈，炮射導(dǎo)彈目前還有很大的發(fā)展?jié)摿Γ貏e是炮射自尋的防空導(dǎo)彈技術(shù)，使坦克擁有抗衡武裝直升機(jī)，甚至攻擊機(jī)的“撒手锏”[9]。

大多數(shù)炮射導(dǎo)彈采用直瞄攻擊模式，不僅耗時(shí)長(zhǎng)，而且精確度與毀傷無(wú)法達(dá)到預(yù)想效果。目前，西方國(guó)家制導(dǎo)彈藥多采用制導(dǎo)體制，自尋的模式代表著制導(dǎo)彈藥的發(fā)展方向。

1 炮射導(dǎo)彈系統(tǒng)數(shù)學(xué)模型

由理論力學(xué)可知，導(dǎo)彈運(yùn)動(dòng)在每一瞬時(shí)，把導(dǎo)彈當(dāng)作一個(gè)質(zhì)量不變的，在氣動(dòng)力、推力、操縱力等作用下運(yùn)動(dòng)的剛體來(lái)處理[9]。

通過(guò)這些假定，來(lái)研究導(dǎo)彈質(zhì)心的移動(dòng)。并且，在利用質(zhì)心運(yùn)動(dòng)的動(dòng)量定理時(shí)，應(yīng)將彈道坐標(biāo)系作為動(dòng)坐標(biāo)系，這樣會(huì)使方程式簡(jiǎn)單明了，更便于編寫程序，求解方程組。

1.1 導(dǎo)彈的受力分析

導(dǎo)彈在地球周圍大氣里飛行，會(huì)受到重力、空氣動(dòng)力和發(fā)動(dòng)機(jī)推力的共同作用。在導(dǎo)彈上的重力為地心引力和離心慣性力的矢量和。

式（1）、（2）中：cx、cy、cz分別為阻力系數(shù)、升力系數(shù)和側(cè)向力系數(shù)；mx、my、mz分別稱為滾動(dòng)力矩系數(shù)、偏航力矩系數(shù)和俯仰力矩系數(shù)；L為特征長(zhǎng)度；S為特征面積。

發(fā)動(dòng)機(jī)推力是由發(fā)動(dòng)機(jī)內(nèi)的燃?xì)饬饕愿咚賴娚涑龆a(chǎn)生的反作用力。

1.2 導(dǎo)彈數(shù)學(xué)模型建立

本文的彈道設(shè)計(jì)為方案飛行彈道，末端采用比例導(dǎo)引法制導(dǎo)，實(shí)現(xiàn)炮射導(dǎo)彈的“自尋的”。導(dǎo)彈的飛行路線設(shè)計(jì)：炮口發(fā)射→無(wú)控階段→程序爬升階段→平飛階段→比例導(dǎo)引段毀傷目標(biāo)。圖1為導(dǎo)彈和目標(biāo)的運(yùn)動(dòng)學(xué)幾何關(guān)系。

圖1 導(dǎo)彈和目標(biāo)的運(yùn)動(dòng)學(xué)幾何關(guān)系Fig.1 The kinematic geometry relationshipbetween missile and target

1）無(wú)控升段。在無(wú)控彈道段為自由飛行，忽略其他干擾，只考慮重力和氣動(dòng)力的作用。

2）程序爬升階段。對(duì)于程序爬升段，目的是提升彈道高度，增加射程和為實(shí)現(xiàn)大落腳準(zhǔn)備。對(duì)于爬升段的方案飛行，考慮到技術(shù)上的可行性，實(shí)現(xiàn)的難易程度以及精確性，本文采用給定俯仰角變化規(guī)律??(t)的飛行方案。這段時(shí)間內(nèi)，發(fā)動(dòng)機(jī)開始工作，相比于無(wú)控段的彈道方程組，爬升段要引入一個(gè)新的變量mc（質(zhì)量秒流量）。

3）平飛階段。導(dǎo)彈在平飛段不會(huì)進(jìn)行大的機(jī)動(dòng)飛行（速度、高度恒定），并且在飛行期間彈道傾角不會(huì)改變。對(duì)于平飛階段，常采用超低空飛行或大高度飛行，以避開敵人雷達(dá)偵查和防空火力的射擊。

式（3）中：v為導(dǎo)彈的速度；P為發(fā)動(dòng)機(jī)推力；X為軸向力；g為重力加速度；θ為彈道傾角；ωz為俯仰角角速度；Mz為靜俯仰力矩；Jz為赤道轉(zhuǎn)動(dòng)慣量；-mc為質(zhì)量流量；?為俯仰角；α為平衡攻角。

導(dǎo)彈—目標(biāo)相對(duì)位置：

導(dǎo)彈—目標(biāo)相對(duì)速度：

導(dǎo)彈—目標(biāo)相對(duì)距離：

導(dǎo)彈—目標(biāo)視線角速度：

導(dǎo)彈—目標(biāo)相對(duì)距離變化率：

2 彈道方程數(shù)值解法

彈道方程組是一階變系數(shù)聯(lián)立方程組，只能用數(shù)值方法求得數(shù)值解，僅在一些特定條件下經(jīng)過(guò)適當(dāng)?shù)暮?jiǎn)化才能求得近似解析解。

本節(jié)通過(guò)對(duì)比龍格—庫(kù)塔法（Runge—Kutta）、阿當(dāng)姆斯法以及歐拉法，選出合適的方法求解彈道微分方程[13-16]。

2.1 龍格—庫(kù)塔法

龍格—庫(kù)塔法實(shí)質(zhì)上是以函數(shù)y(x)的臺(tái)勞級(jí)數(shù)為基礎(chǔ)的一種改進(jìn)方法。最常用的是4階龍格—庫(kù)塔法，若已知在點(diǎn)n處的值(tn,y1n,y2n,…,ymn) ，則求點(diǎn)n+1處的函數(shù)值的公式可由龍格—庫(kù)塔公式推導(dǎo)出：

式（9）中：

2.2 阿當(dāng)姆斯預(yù)報(bào)—校正法

2.3 歐拉法

3 自尋的炮射導(dǎo)彈彈道仿真

3.1 仿真方法

本文采用歐拉法來(lái)求解彈道方程組，在初始條件給定的情況下，程序可以得到以時(shí)間t為自變量的彈道上的任意諸元。對(duì)于剛體彈道方程，時(shí)間的步長(zhǎng)h必須小于0.005 s。否則，計(jì)算發(fā)散。故在使用歐拉法求解彈道方程組時(shí)，選取0.001作為積分步長(zhǎng)，既可便捷的求解彈道方程，又可滿足誤差要求。比例導(dǎo)引方法的比例系數(shù)K，應(yīng)選擇在1～∞的范圍內(nèi)，通常在2～6 的范圍內(nèi)選擇。本文選擇K=4。

采用Matlab 編寫程序來(lái)實(shí)現(xiàn)自尋的炮射導(dǎo)彈的各段彈道，程序以時(shí)間t作為全局變量，將計(jì)算出的數(shù)值存儲(chǔ)在數(shù)組中，便于以后作圖[17-20]。

3.2 仿真結(jié)果及分析

仿真初始條件為：炮射導(dǎo)彈以初速度300 m/s 發(fā)射，初始攻角α=0°，初始俯仰角?=18°，初始彈道的傾角θ=18°，得到仿真結(jié)果如圖2所示。

圖2 彈道曲線Fig.2 Ballistic curve

圖2為自尋的炮射導(dǎo)彈的彈道曲線。從圖2可以得出，導(dǎo)彈以初速300 m/s 發(fā)射，在發(fā)動(dòng)機(jī)不工作的情況下爬升一段距離。此時(shí)，彈道傾角逐漸減小；隨后，發(fā)動(dòng)機(jī)開始工作，導(dǎo)彈在發(fā)動(dòng)機(jī)的推動(dòng)下爬升，爬升過(guò)程中隨著彈道傾角的增加，爬升率下降；當(dāng)t=22 s時(shí)，導(dǎo)彈爬升到最高點(diǎn)；然后，經(jīng)歷一段平飛，導(dǎo)彈迅速調(diào)整姿態(tài)，彈道傾角在平飛段為0，平飛段的目的為使導(dǎo)彈飛到目標(biāo)上空，增加導(dǎo)彈射程；當(dāng)彈目距離達(dá)到預(yù)定值之后，導(dǎo)彈加速下落，比例導(dǎo)引法擊毀目標(biāo)。

平衡攻角和俯仰角與彈道傾角之間存在著簡(jiǎn)單的幾何關(guān)系：α=?-θ。

如圖3 所示，在導(dǎo)彈發(fā)射時(shí)，平衡攻角α=0°，平衡攻角α一直保持為0°，導(dǎo)彈保持平衡飛行，直至彈道轉(zhuǎn)為比例導(dǎo)引段，導(dǎo)彈姿態(tài)調(diào)整，進(jìn)入末端制導(dǎo)，實(shí)際攻角開始減小，其收斂于平衡攻角，導(dǎo)彈飛行穩(wěn)定。顯然，圖示攻角一直保持為0，僅為理想狀態(tài)下攻角曲線，實(shí)際的攻角曲線應(yīng)如圖4所示。

圖4 所示為某無(wú)控彈道的攻角的仿真曲線，舵偏角δz=0°，與此對(duì)應(yīng)的平衡攻角α也為0。導(dǎo)彈在飛行過(guò)程中攻角的變化會(huì)產(chǎn)生嚴(yán)重的振蕩，因?yàn)閺楏w阻尼力矩的存在，攻角的變化會(huì)收斂于平衡攻角。

圖3 平衡攻角-時(shí)間曲線Fig.3 Angle of attack-time curve

圖4 平衡攻角-時(shí)間曲線Fig.4 Angle of attack-time curve

由圖5 速度時(shí)間變化的曲線可以看出，速度曲線近似為3段折線。導(dǎo)彈經(jīng)發(fā)射之后，在0～30 s，導(dǎo)彈在較高的切向加速度加速飛行；在30 s 之后，導(dǎo)彈到達(dá)預(yù)定飛行高度，速度保持在625 m/s 左右，開始以勻速飛行；150 s 之后，導(dǎo)彈進(jìn)入末端制導(dǎo)段，速度開始減小，飛向并擊毀目標(biāo)。

依據(jù)圖6 的彈道傾角隨時(shí)間的變化曲線可以看出，在無(wú)控階段，彈道傾角逐漸減小，直到發(fā)動(dòng)機(jī)開始工作，彈道的傾角開始增大；導(dǎo)彈轉(zhuǎn)至平飛段后，彈道傾角急速減小到0，并在平飛段一直保持在0；在比例導(dǎo)引階段，導(dǎo)彈加速向目標(biāo)飛去，并且在此過(guò)程中，彈道傾角不斷減小。

圖5 速度-時(shí)間曲線Fig.5 Velocity-time curve

圖6 彈道傾角-時(shí)間曲線Fig.6 Ballistic inclination-time curve

圖7 俯仰角-時(shí)間曲線Fig.7 Pitch angle-time curve

4 結(jié)論

根據(jù)以上一系列圖示的仿真結(jié)果，炮射導(dǎo)彈可以實(shí)現(xiàn)自尋的模式。本文在前人提出的“自尋的”炮射導(dǎo)彈的基礎(chǔ)之上，設(shè)計(jì)了上文所示的可以實(shí)現(xiàn)頂部攻擊的、對(duì)靜止目標(biāo)或慢速移動(dòng)目標(biāo)實(shí)現(xiàn)高毀傷打擊的炮射導(dǎo)彈彈道模型。

從仿真的結(jié)果來(lái)看，當(dāng)初速等參數(shù)一定時(shí)，就可以確定一條彈道曲線。所以，只要改變模型中的相應(yīng)參數(shù)，導(dǎo)彈飛行時(shí)間、射程、彈道高度、落點(diǎn)速度都會(huì)做出相應(yīng)的改變。本文所提到的研究方法、思路和分析方法可以為未來(lái)實(shí)現(xiàn)炮射導(dǎo)彈智能化提供借鑒，具有一定的實(shí)際意義。