王少純 王哲峰 陳昱辰 劉艷龍 張馨丹

（1.沈陽航空航天大學(xué)航空宇航學(xué)院 沈陽 110000）（2.沈陽大學(xué)環(huán)境工程學(xué)院 沈陽 110044）（3.沈陽航空航天大學(xué)航空發(fā)動機(jī)學(xué)院 沈陽 110000）

1 引言

隨著我國國防力量的發(fā)展，彈道導(dǎo)彈已經(jīng)成為捍衛(wèi)祖國主權(quán)的戰(zhàn)略力量。當(dāng)導(dǎo)彈在飛行時，其質(zhì)心在空間內(nèi)會形成一條運(yùn)動軌跡，這條軌跡就是通常所說的基準(zhǔn)彈道，我們通過對彈道的分析研究可以提高導(dǎo)彈的性能和命中率。

因此在設(shè)計導(dǎo)彈時，要經(jīng)過彈道分析掌握導(dǎo)彈的運(yùn)動規(guī)律，進(jìn)而選擇正確的結(jié)構(gòu)參數(shù)、選擇合適的飛行彈道以及進(jìn)行正確精準(zhǔn)的彈道計算、進(jìn)而評定導(dǎo)彈的基本性能參數(shù)、同時也能夠?yàn)閷?dǎo)彈的飛行試驗(yàn)提供需要的理論彈道參數(shù)數(shù)據(jù)。導(dǎo)彈彈道按照坐標(biāo)系不同劃分為絕對彈道和相對彈道。其中絕對彈道是建立在慣性坐標(biāo)系中的導(dǎo)彈彈道模型。例如地空導(dǎo)彈、岸艦導(dǎo)彈、近程的地地導(dǎo)彈。本文通過對導(dǎo)彈飛行過程中位置參數(shù)和艦船運(yùn)動時的位置參數(shù)進(jìn)行分析，建立模型，確定導(dǎo)彈的軌道模型并解決相關(guān)問題。

2 導(dǎo)彈打擊目標(biāo)的具體分析

導(dǎo)彈的制導(dǎo)過程是一個復(fù)雜的過程，在這一過程中，可能會受到大氣環(huán)境、攻擊參數(shù)和對方反導(dǎo)系統(tǒng)攔截等很多條件的干擾。本文只考慮導(dǎo)彈本身軌跡的改變和攻擊目標(biāo)（航母）的運(yùn)動以及地球的自轉(zhuǎn)及曲率影響。

當(dāng)目標(biāo)為靜止?fàn)顟B(tài)時，通過參數(shù)方程對三維空間內(nèi)的彈道進(jìn)行參數(shù)化分解，將導(dǎo)彈彈道這一復(fù)合的軌跡分解成為三個相互垂直的分運(yùn)動。通過對三維位置點(diǎn)的分解，可以得到三個方向關(guān)于時間的速度參數(shù)方程，并代入地球自轉(zhuǎn)［1]，地主曲率和變軌技術(shù)等產(chǎn)生的影響因素，即可描繪導(dǎo)彈的彈道軌跡。

當(dāng)目標(biāo)具有一定運(yùn)動速度，且運(yùn)動方向與導(dǎo)彈三個分運(yùn)動方向其中一個方向在同一直線上，因此只需要在彈道軌跡的基礎(chǔ)上進(jìn)行相關(guān)方向運(yùn)動參數(shù)的修改即可。在此問題中，我們只需要考慮導(dǎo)彈第二階段的飛行，即發(fā)射中段的導(dǎo)彈飛行模型即可達(dá)到跟蹤的目的。以拋物線的頂點(diǎn)為起點(diǎn)，完善導(dǎo)彈中段的動態(tài)飛行模型與算法。

通過對導(dǎo)彈發(fā)射點(diǎn)和航母起始點(diǎn)在地球表面相對位置的確定，得到如下的相對位置圖。

圖1 導(dǎo)彈發(fā)射點(diǎn)與航母在地球上位置確定

通過二者的坐標(biāo)和在途中的相對位置可以看出，在初始時刻，二者的維度不同。因此，即使航母不運(yùn)動，導(dǎo)彈若命中目標(biāo)，也需要跨越緯度飛行。同時導(dǎo)彈必須要向前飛行和向上方飛行以升高飛行高度。我們可以將導(dǎo)彈的飛行軌跡看成一個三維的曲線［2]。

對導(dǎo)彈的運(yùn)動建立三維模型，將其運(yùn)動軌跡分解為豎直（z）方向，東西（x）方向，南北（y）方向。

若從上空觀測導(dǎo)彈的飛行軌跡，會發(fā)現(xiàn)如圖2（a）中鎖死的彈道軌跡，從某一緯度所在平面觀測，會發(fā)現(xiàn)如圖2（b）中所示的彈道軌跡。

圖2 不同角度的導(dǎo)彈彈道軌跡

導(dǎo)彈在發(fā)射后的運(yùn)動軌跡可分為三部分［3]，即發(fā)射段、中段和末段。如圖3所示，導(dǎo)彈在發(fā)射段受重力和推進(jìn)器產(chǎn)生的推力作用（空氣阻力忽略不計）做拋物線式運(yùn)動；在發(fā)射中段，導(dǎo)彈的運(yùn)動可再分為兩部分，一是只受重力和空氣阻力作用的自由彈道，二是在衛(wèi)星調(diào)控的基礎(chǔ)上，受到重力和空氣阻力共同作用的再入彈道；當(dāng)航母進(jìn)入導(dǎo)彈的打擊范圍，導(dǎo)彈進(jìn)入末段軌跡，自主打擊目標(biāo)。

圖3 導(dǎo)彈發(fā)射全程示意圖

將導(dǎo)彈在x，y和z方向上的運(yùn)動均看成為多次變加速運(yùn)動的合運(yùn)動，得到如下的運(yùn)動方程。

將導(dǎo)彈的運(yùn)動分為若干個函數(shù)，式中的an，bn，cn等分別為導(dǎo)彈在第n段中x、y、z方向上的加速度，tn為導(dǎo)彈在第n段所經(jīng)過的時間［4]。

以導(dǎo)彈發(fā)射車為坐標(biāo)原點(diǎn)，xi為原點(diǎn)與地球的切線，方向指向正東方向，xi軸和zi軸相互垂直，方向豎直向上。導(dǎo)彈在發(fā)射段時，由于只受燃料推力，空氣阻力和重力的作用，其軌跡可以簡化為拋物線（y等于0）。在o-xizi坐標(biāo)系內(nèi)可以建立其彈道方程：

式中v0，h0分別為導(dǎo)彈發(fā)動機(jī)關(guān)機(jī)時的速度和高度，xt0為導(dǎo)彈在發(fā)動機(jī)關(guān)機(jī)時，在平面坐標(biāo)系o-xizi上的橫坐標(biāo)，vx0，vz0分別為導(dǎo)彈發(fā)動機(jī)關(guān)機(jī)處，速度在坐標(biāo)系上沿oxi軸和ozi軸上的分量，vx1，vz1和xt1，zt1分別為導(dǎo)彈于發(fā)射段時，在坐標(biāo)系上的速度和位移分量，θ0為導(dǎo)彈發(fā)射段開始時的傾角（極?。瑃為導(dǎo)彈在發(fā)射段的運(yùn)行的時間，a為導(dǎo)彈在上升時所受到的加速度矢量和，H為導(dǎo)彈飛行高度［5]。

在導(dǎo)彈的中段，導(dǎo)彈可分為自由彈道和再入彈道兩個階段，自由彈道是近似拋物線的運(yùn)動軌跡，再入彈道為導(dǎo)彈從脫離拋物線到再一次進(jìn)入到大氣層后的運(yùn)動軌跡。

導(dǎo)彈的自由彈道可以看成是導(dǎo)彈質(zhì)心運(yùn)動和繞心運(yùn)動的結(jié)合。在中段，導(dǎo)彈的運(yùn)動微分方程應(yīng)由運(yùn)動學(xué)方程，控制方程，質(zhì)心運(yùn)動方程和繞質(zhì)心運(yùn)動的動力學(xué)方程共同組成。我們同樣以導(dǎo)彈發(fā)射車為坐標(biāo)原點(diǎn)，研究導(dǎo)彈姿態(tài)運(yùn)動方程［6]。

導(dǎo)彈質(zhì)心運(yùn)動方程可由動力學(xué)方程和質(zhì)心運(yùn)動方程積分求得：

導(dǎo)彈姿態(tài)方程可以利用導(dǎo)彈轉(zhuǎn)動角速度ω1在坐標(biāo)系的投影方程：

結(jié)合導(dǎo)彈發(fā)射的實(shí)體情況來看，ψ和γ都是很小的量，可以近似認(rèn)為 sinψ≈ψ，sinγ≈γ，cosψ≈cosγ≈1。

與導(dǎo)彈進(jìn)入自由彈道相比，導(dǎo)彈在進(jìn)入再入彈道后所受的力主要為大氣產(chǎn)生的阻力和導(dǎo)彈自身所受的重力。所以在以導(dǎo)彈發(fā)射車為坐標(biāo)原點(diǎn)的三維坐標(biāo)系中，可以建立再入彈道的微分方程：

式中：m為導(dǎo)彈去掉推進(jìn)器后的彈頭質(zhì)量；vt為導(dǎo)彈質(zhì)心在三維坐標(biāo)系上速度的模；XF為空氣阻力，并且：

式中ρ為空氣密度；s為導(dǎo)彈橫截面積Rb為導(dǎo)彈底部半徑；d為導(dǎo)彈最大直徑；vl為音速，v為導(dǎo)彈速度；i為彈性系數(shù)，該系數(shù)與導(dǎo)彈外形有關(guān)；為馬赫數(shù)對應(yīng)的標(biāo)準(zhǔn)阻力系數(shù)；F為阻力函數(shù)，另：

H為無因次空氣比重，其跟隨導(dǎo)彈高度變化的經(jīng)驗(yàn)公式為

當(dāng)高度小于9300m時：

當(dāng)高度大于等于9300m時：

式中：k=0.0001。

在導(dǎo)彈中段的自由彈道結(jié)束后，導(dǎo)彈馬上進(jìn)入再入彈道，結(jié)合自由彈道和再入彈道，可以建立方程為

式中：h1為導(dǎo)彈結(jié)束自由彈道時的高度；θz為再入角；tz為導(dǎo)彈結(jié)束自由彈道時的時間，θ1為導(dǎo)彈進(jìn)入再入彈道時的傾角。

考慮到在制導(dǎo)導(dǎo)彈的打擊過程中，對方的反導(dǎo)系統(tǒng)會進(jìn)行攔截，故找出有效的行進(jìn)軌跡來對導(dǎo)彈進(jìn)行變軌處理，來最大化降低攔截的可能性十分必要［7]。機(jī)動變軌技術(shù)是導(dǎo)彈在飛行中可隨時改變其彈道，以躲避敵方反導(dǎo)防御系統(tǒng)攔截的一種突防技術(shù)，通常分為全彈道變軌和彈道中末段變軌兩種。因?yàn)閷?dǎo)彈發(fā)射初期是在己方國土內(nèi)，所以敵人很難如此迅速的檢測到發(fā)射信號，并且導(dǎo)彈在彈道末端是以二十秒內(nèi)保持1000m/s的超音速進(jìn)行俯沖，在這種高速短時狀態(tài)下敵人很難精確對導(dǎo)彈進(jìn)行攔截。因此，為節(jié)約發(fā)動機(jī)燃料和實(shí)現(xiàn)最高效的精確打擊，我們將變軌技術(shù)只應(yīng)用在中段?，F(xiàn)階段常用的變軌軌跡為蛇形機(jī)動和s形機(jī)動，其變軌原理均是有自帶的小火箭在短時間內(nèi)產(chǎn)生推力而完成。

變軌過程加速度是由推力T大小，方向單位向量uT和導(dǎo)彈的質(zhì)量m決定的，即

式中導(dǎo)彈的質(zhì)量隨燃料燃燒而減小，T可以表示為

由上可得，推力方向也就是變軌方向是決定推力大小等參數(shù)的決定因素。所以我們研究出了基于大風(fēng)級下的海浪波動模型。此模型通過模擬在大風(fēng)級自然風(fēng)下海浪的波形，來形成導(dǎo)彈的變軌軌跡。此軌跡具有隨機(jī)性強(qiáng)，風(fēng)級固定下穩(wěn)定的特性，所以對最終打擊精度不會有太大誤差。我們將海浪高度在坐標(biāo)系內(nèi)的變化用一個三元函數(shù)來表示［8]：

其中，ai，ki，wi，εi分別為第i次諧波的波幅，波數(shù)，角頻率，初相位。

我們知道導(dǎo)彈變軌是由導(dǎo)彈攜帶的多個矢量火箭發(fā)動機(jī)，產(chǎn)生的脈沖推力來完成，在變軌過程中由于減少了前進(jìn)方向的速度分量，會增加飛行時間。為了降低導(dǎo)彈負(fù)擔(dān)，和避免貽誤戰(zhàn)機(jī)，并我們將海浪模型簡化成二維問題。

這種“二元不規(guī)則波”可表示為

3 反導(dǎo)導(dǎo)彈殺傷區(qū)模型

從航空母艦發(fā)射的反導(dǎo)導(dǎo)彈的殺傷區(qū)是其不低于某一給定概率殺傷目標(biāo)的空間范圍，它主要是有反導(dǎo)系統(tǒng)的火控系統(tǒng)的參數(shù)決定，并且需要同時滿足以下四個不等式。這表示出的區(qū)域就是反導(dǎo)系統(tǒng)的有效殺傷區(qū)。其中我們不考慮時間對相遇的影響［9]。

式中：Qm為OX與正北方向夾角。

我們可以把地球看成是一個繞地軸自西向東旋轉(zhuǎn)的球體，由于它的旋轉(zhuǎn)運(yùn)動，地球表面上從兩極到赤道上，各點(diǎn)的運(yùn)動速度是隨當(dāng)量半徑的增長而增大的。越靠近地球兩端，表面的運(yùn)動線速度越大，反之越小。而且，在同一維度上，海拔高度的不同也會導(dǎo)致線速度的偏差，比如說摩天大樓樓頂?shù)木€速度要大于底座的線速度。所以計算導(dǎo)彈軌跡時，如果不計算自轉(zhuǎn)產(chǎn)生的影響，那么在導(dǎo)彈的長距離跨經(jīng)緯度的飛行過程中，就會因?yàn)檫@兩方面因素的影響，最終會導(dǎo)致無法精確打擊目標(biāo)。

這其中一個重要參數(shù)為偏近角E

平進(jìn)點(diǎn)角：

式中：M0=E0-esinE0相當(dāng)于T=0時的平近點(diǎn)角。

從開普勒方程可知，如果已知偏近點(diǎn)角E求M非常容易，但是在已知M的條件下求E，這是一個超越方程，用微分迭代比較容易求解。

即

收斂條件：

之后為精確計算制導(dǎo)導(dǎo)彈的軌道和打擊點(diǎn)，我們在兩個慣性參考系中，選擇對發(fā)射段結(jié)束的分離點(diǎn)進(jìn)行分析。A坐標(biāo)系原點(diǎn)為地球質(zhì)心，X，Y，Z軸滿足右手法則且固定在地球上；B坐標(biāo)系的X，Y，Z三軸的指向和A相同，但是它們隨地球自轉(zhuǎn)而旋轉(zhuǎn)，即三軸的指向在不斷變化。由上述兩個坐標(biāo)系的定義可得：1）在發(fā)射段結(jié)束時刻，A坐標(biāo)系中導(dǎo)彈的位置矢量與B坐標(biāo)系的導(dǎo)彈位置矢量相同［10]，即：rKA=rKB，rKA和rKB的兩個位置分量分別表示關(guān)機(jī)點(diǎn)經(jīng)度和關(guān)機(jī)點(diǎn)緯度，K表示關(guān)機(jī)點(diǎn)。2）導(dǎo)彈在坐標(biāo)系A(chǔ)，B中的速度不一樣。速度矢量差為地球自轉(zhuǎn)速度在該高度引起的牽連速度，用VKe來表示。導(dǎo)彈在A坐標(biāo)系中的速度矢量VKA（Vx，Vy，Vz）為已知參數(shù)，則在B坐標(biāo)系中的速度矢量VKB為

地球自轉(zhuǎn)轉(zhuǎn)向?yàn)樽晕飨驏|，其角速度為w，所以，其引起的牽連速度VKe的方向沿K點(diǎn)自西向東，大小為ωrKcosΘk，這里Θk為K點(diǎn)緯度，rK為K點(diǎn)

賦初值：

迭代公式：到地心的距離。

因此，考慮地球自轉(zhuǎn)的導(dǎo)彈速度為

XK，YK，ZK是rKA的直角坐標(biāo)形式。因此，當(dāng)不考慮地球自轉(zhuǎn)時，關(guān)機(jī)點(diǎn)的速度參數(shù)為rKA，rKA，而上述驗(yàn)證的考慮地球自轉(zhuǎn)時的關(guān)機(jī)點(diǎn)參數(shù)為rKA，rKB。

如果把地球視為一個半徑為6300km的球體，那么它的曲率不能忽略。

圖4 地球曲率對導(dǎo)彈打擊點(diǎn)的影響

不考慮地球曲率時，導(dǎo)彈平在飛段的彈道傾角θ≈0°，而實(shí)際上導(dǎo)彈存在彈道傾角θ'：

式中，d為導(dǎo)彈的飛行距離。

當(dāng)導(dǎo)彈飛行340km時，θ'-θ=-3.09°，所以地球曲率對導(dǎo)彈打擊點(diǎn)的影響不可忽略。

圖5 俯仰角示意圖

因?yàn)榻虒?dǎo)彈的主軸不跟蹤地球曲面［11]，導(dǎo)彈發(fā)射時決定其運(yùn)動狀態(tài)，而且導(dǎo)彈的各方向角是以各向主軸作為基準(zhǔn)測出的，所以導(dǎo)彈的各方向角的定義取決于以上所定義的坐標(biāo)系。因?yàn)閷?dǎo)彈發(fā)射時的水平面與主軸始終成一定角度，所以導(dǎo)彈的傾角始終是導(dǎo)彈發(fā)射時的水平面與彈軸的夾角(?)，而不是和所處位置水平面的夾角(?')，如圖3所示。所以只需假設(shè)原始坐標(biāo)系是慣性參考系，根據(jù)導(dǎo)彈距離地表的高度(h)及其變化率(h')來模擬無線電高度表測得的導(dǎo)彈高度及其變化率，將其加入到導(dǎo)彈的高度控制回路中，并且把重力加速度改為g'，建立考慮地球曲率的導(dǎo)彈彈道模型：

綜上我們可以得出導(dǎo)彈打擊模型為

通過對導(dǎo)彈運(yùn)動的三維建模進(jìn)行分析求解，我們已經(jīng)得出了導(dǎo)彈打擊靜止點(diǎn)的運(yùn)動方程。當(dāng)航母向正南方向（y方向）以速度Vk勻速運(yùn)動時，在上一個模型基礎(chǔ)上，給導(dǎo)彈運(yùn)動方程添加一個y方向的速度Vk即可得到問題二的模型［12]。

導(dǎo)彈質(zhì)心運(yùn)動方程可由動力學(xué)質(zhì)心運(yùn)動方程積分得：

導(dǎo)彈姿態(tài)方程可以利用導(dǎo)彈轉(zhuǎn)動角速度ω1在坐標(biāo)系的投影方程：

與問題一相同，結(jié)合導(dǎo)彈發(fā)射的實(shí)體情況，ψ和γ都是很小的量，可以近似認(rèn)為sinψ≈ψ，sinγ≈γ，cosψ≈cosγ≈1。

在目標(biāo)移動時，導(dǎo)彈的中段自由彈道運(yùn)動方程為

導(dǎo)彈進(jìn)入再入段后，除y方向上加的勻速運(yùn)動外，其運(yùn)動方程和問題一相同，即：

式中：m為導(dǎo)彈去掉推進(jìn)器后的彈頭質(zhì)量；vt為導(dǎo)彈質(zhì)心在三維坐標(biāo)系上速度的模；XF為空氣阻力，并且：式中ρ為空氣密度；s為導(dǎo)彈橫截面積Rb為導(dǎo)彈底部半徑；d為導(dǎo)彈最大直徑；vl為音速，v為導(dǎo)彈速度；i為彈性系數(shù)，該系數(shù)與導(dǎo)彈外形有關(guān)為馬赫數(shù)對應(yīng)的標(biāo)準(zhǔn)阻力系數(shù)；F為阻力函數(shù)，另：

H為無因次空氣比重，其大小跟隨導(dǎo)彈高度變化的公式為

結(jié)合自由彈道和再入彈道，建立方程：

式中：h1為導(dǎo)彈結(jié)束自由彈道時的高度；θz為再入角；tz為導(dǎo)彈結(jié)束自由彈道時的時間，θ1為導(dǎo)彈進(jìn)入再入彈道時的傾角

所以在目標(biāo)航母以Vk的速度向南方向行駛時導(dǎo)彈的中段運(yùn)動方程為

4 結(jié)語

該模型針對打擊靜止目標(biāo)時的軌跡和打擊靜止目標(biāo)時的軌跡均進(jìn)行了數(shù)字化分析，并結(jié)合于公式，建立了全面的導(dǎo)彈打擊和預(yù)防攔截模型，增加了精準(zhǔn)率并減小了偏差。同時，考慮到到溫度，天氣，空氣質(zhì)量等對本模型的干擾，均將其設(shè)為干擾因子，計算比例后帶入模型。