大氣層外彈道導(dǎo)彈中段攔截彈道規(guī)劃

趙蒙， 端軍紅， 王明宇， 殷雙斌

(1.空軍工程大學(xué) 防空反導(dǎo)學(xué)院， 陜西 西安 710051； 2.93792部隊(duì)， 河北 廊坊 065000)

0 引言

彈道導(dǎo)彈具有射程遠(yuǎn)、速度快、精度高、威力大、突防能力強(qiáng)等特點(diǎn)，如何對來襲彈道導(dǎo)彈實(shí)施有效攔截一直是各軍事強(qiáng)國防空能力建設(shè)的重點(diǎn)。2018年隨著美俄相繼退出《中導(dǎo)條約》，反彈道導(dǎo)彈(簡稱反導(dǎo))再次成為世界各國關(guān)注的熱點(diǎn)問題之一。彈道導(dǎo)彈中段攔截又稱中段反導(dǎo)，一般是指對處于飛行中段的彈道導(dǎo)彈進(jìn)行攔截摧毀，其攔截過程如圖1所示。中段反導(dǎo)作戰(zhàn)涉及的關(guān)鍵技術(shù)非常多，但攔截彈道規(guī)劃是中段反導(dǎo)作戰(zhàn)制導(dǎo)、控制攔截彈飛行的前提和依據(jù)，國內(nèi)外有關(guān)于空間軌道目標(biāo)攔截問題的文獻(xiàn)資料相對較多，對于非機(jī)動(dòng)變軌的彈道導(dǎo)彈目標(biāo)，相對于空基武器而言，使用地基攔截彈對目標(biāo)實(shí)施攔截是一種快捷、有效而且技術(shù)復(fù)雜度相對較低的方法，但是目前有關(guān)于地基攔截彈道規(guī)劃的文獻(xiàn)資料相對較少。

本文利用橢圓彈道理論對目標(biāo)飛行軌跡進(jìn)行預(yù)測，在獲得遭遇點(diǎn)的位置后，依據(jù)攔截彈發(fā)射點(diǎn)位置和攔截彈飛行時(shí)間，基于Lambert方程建立攔截彈飛行彈道規(guī)劃模型，并對攔截彈規(guī)劃模型的有效性進(jìn)行仿真驗(yàn)證。該模型為中段反導(dǎo)攔截彈道規(guī)劃設(shè)計(jì)提供一種思路和方法，若攔截彈符合所規(guī)劃的攔截彈道在發(fā)射能量、飛行時(shí)間、攔截高度、飛行程序設(shè)計(jì)等約束條件的要求，理論上可滿足中段反導(dǎo)攔截作戰(zhàn)要求，能夠?qū)崿F(xiàn)對目標(biāo)的攔截。

圖1 反導(dǎo)攔截過程示意圖Fig.1 Schematic diagram of anti-missile interception

1 遭遇點(diǎn)預(yù)測模型

1.1 坐標(biāo)系定義

在建立目標(biāo)與攔截彈的遭遇點(diǎn)以及攔截彈彈道規(guī)劃模型中，主要涉及的坐標(biāo)系有地心慣性坐標(biāo)系、地心固定坐標(biāo)系、站心坐標(biāo)系以及大地坐標(biāo)系，相關(guān)坐標(biāo)系的定義以及各個(gè)坐標(biāo)系之間的轉(zhuǎn)換在文獻(xiàn)[2]中都有詳細(xì)描述和轉(zhuǎn)換公式。本文后續(xù)模型建立過程中如不特殊說明，所有計(jì)算都在地心慣性坐標(biāo)系下進(jìn)行。

在大地坐標(biāo)系下，假設(shè)預(yù)警探測系統(tǒng)探測到彈道導(dǎo)彈發(fā)動(dòng)機(jī)關(guān)機(jī)后，時(shí)刻某一點(diǎn)的位置坐標(biāo)為(,,)、速度矢量為(1,2,3)。通過坐標(biāo)變換可得地心慣性坐標(biāo)系下彈道導(dǎo)彈的位置矢量=(，，)和速度矢量=(，，)。

則彈道目標(biāo)在點(diǎn)的極徑、速度以及速度矢量與當(dāng)?shù)厮矫骈g的夾角的大小分別為

(1)

(2)

(3)

1.2 橢圓彈道參數(shù)

彈道導(dǎo)彈的飛行彈道按照發(fā)動(dòng)機(jī)的工作狀態(tài)，可分主動(dòng)段和被動(dòng)段兩大部分。根據(jù)彈道導(dǎo)彈飛行過程中受力的不同，被動(dòng)段又可分為在大氣層以外的自由飛行段和進(jìn)入大氣層后的再入飛行段。其中在自由飛行段由于彈道導(dǎo)彈始終處于大氣層以外飛行，可認(rèn)為導(dǎo)彈只受到地心引力的影響，則彈道導(dǎo)彈與地球就構(gòu)成一個(gè)二體系統(tǒng)，其飛行軌跡近似為以地心為焦點(diǎn)的橢圓軌道。

彈道導(dǎo)彈目標(biāo)自由飛行段的彈道方程為

(4)

式中：為地心矢徑的模；為橢圓彈道的半通徑；為橢圓彈道的偏心率；為橢圓彈道的真近點(diǎn)角。

由目標(biāo)狀態(tài)矢量、可得橢圓軌道的半通徑、偏心率及橢圓半長軸具體計(jì)算公式如下：

(5)

(6)

(7)

式中：=，為萬有引力常數(shù)，為地球質(zhì)量。

則時(shí)刻目標(biāo)在點(diǎn)的真近點(diǎn)角和偏近地點(diǎn)角，可由(8)式、(9)式唯一確定：

(8)

(9)

式中：2與2總是處于同一象限內(nèi)。

1.3 目標(biāo)軌跡預(yù)測

根據(jù)二體理論，時(shí)刻后直至彈道導(dǎo)彈再入大氣層這段時(shí)間內(nèi)，彈道導(dǎo)彈目標(biāo)任意時(shí)刻的位置矢量可由橢圓軌道方程解算。設(shè)經(jīng)過Δ=-后，時(shí)刻的偏近點(diǎn)角為，則可利用牛頓迭代法求解下式獲得：

(10)

則時(shí)刻目標(biāo)在點(diǎn)的位置矢量=(1,1,1)、速度矢量=(1,1,1)預(yù)測模型為

(11)

(12)

式中：Δ=-。

2 攔截彈道規(guī)劃模型

2.1 問題描述

中段攔截武器系統(tǒng)的攔截彈飛行過程與彈道導(dǎo)彈基本類似。攔截彈的彈道按照發(fā)動(dòng)機(jī)工作狀態(tài)可分為主動(dòng)段和被動(dòng)段。主動(dòng)段主要依靠發(fā)動(dòng)機(jī)提供的推力按照程序飛行，在攔截彈發(fā)動(dòng)機(jī)關(guān)機(jī)后就進(jìn)入被動(dòng)飛行段，若不考慮中段制導(dǎo)控制情況下，其后續(xù)彈道飛行軌跡完全由主動(dòng)段結(jié)束時(shí)刻攔截彈的位置矢量和速度矢量所決定。攔截彈在主動(dòng)段的受力情況非常復(fù)雜，直接建模非常困難，根據(jù)反導(dǎo)攔截彈主動(dòng)段的飛行程序，應(yīng)用彈道方程可比較容易獲取攔截彈在主動(dòng)段結(jié)束時(shí)刻的運(yùn)動(dòng)狀態(tài)參數(shù)，為了簡化建模的復(fù)雜性，在被動(dòng)段近似認(rèn)為攔截彈的軌跡為以地心為焦點(diǎn)的橢圓彈道，從而將復(fù)雜的彈道規(guī)劃問題轉(zhuǎn)化為經(jīng)典的Lambert問題進(jìn)行處理。

攔截彈彈道規(guī)劃的Lambert問題描述是：已知攔截彈發(fā)射點(diǎn)的位置矢量=(2,2,2)、遭遇點(diǎn)的位置矢量=(3,3,3)以及攔截彈的飛行時(shí)間，規(guī)劃一條橢圓半長軸為、偏心率為、半通徑為的橢圓軌道，實(shí)現(xiàn)攔截彈在時(shí)刻從點(diǎn)發(fā)射，時(shí)刻與目標(biāo)在點(diǎn)遭遇。完成彈道規(guī)劃計(jì)算后，則可解得攔截彈發(fā)射后的位置矢量和速度矢量，根據(jù)攔截彈發(fā)動(dòng)機(jī)關(guān)機(jī)點(diǎn)高度，利用(13)式可求得攔截彈發(fā)動(dòng)機(jī)關(guān)機(jī)點(diǎn)的真近點(diǎn)角為

(13)

式中：為地球赤道半徑大小。

將(13)式代入(8)式、(9)式，可得偏近點(diǎn)角以及攔截彈從發(fā)射后至發(fā)動(dòng)機(jī)關(guān)機(jī)點(diǎn)的飛行時(shí)間Δ，最終可求解得到攔截彈在發(fā)動(dòng)機(jī)關(guān)機(jī)點(diǎn)處的位置矢量和速度矢量。求解出、后，則可通過工程方法設(shè)計(jì)出攔截彈的主動(dòng)段飛行程序。

2.2 彈道規(guī)劃模型

根據(jù)上述描述，在已知、情況下，實(shí)質(zhì)上只要求出攔截彈發(fā)射點(diǎn)或者遭遇點(diǎn)處的速度矢量(2,2,2)、(3,3,3)問題即可解決。根據(jù)Lambert定理，引入普適變量，從而建立飛行時(shí)間與的函數(shù)關(guān)系，根據(jù)文獻(xiàn)[18]，攔截彈的徑向速度滿足

(14)

式中：為徑向速度；為周向速度。

定義滿足

(15)

(16)

式中：為積分常數(shù)。

令=，則飛行時(shí)間的普適公式可表示為

(17)

式中：

(18)

(19)

(20)

利用、()和()表示的拉格朗日系數(shù)和分別為

(21)

(22)

(23)

式中：Δ表示兩點(diǎn)真近點(diǎn)角之差，Δ=-；·=23+23+23。

為求出的值，建立函數(shù)

(24)

將代入(24)式中，通過牛頓迭代法可最終求出的值，將求得的值代入拉格朗日系數(shù)中，由此可得出與為

(25)

(26)

再利用(4)式～(12)式即可得到攔截彈發(fā)射后的位置矢量和速度矢量為

(27)

(28)

即攔截彈的飛行狀態(tài)信息。

2.3 相關(guān)約束條件

時(shí)間約束，由橢圓方程可知對于給定的、兩點(diǎn)，通過、兩點(diǎn)的橢圓彈道有無窮多條。對于反導(dǎo)攔截，需要考慮的是飛行時(shí)間的約束，主要是指目標(biāo)和攔截彈在同一時(shí)刻同時(shí)運(yùn)動(dòng)到預(yù)測遭遇點(diǎn)位置。因此給定了攔截彈的飛行時(shí)間之后，則可唯一地確定出通過攔截點(diǎn)和遭遇點(diǎn)的彈道。

初始段飛行程序約束，由于攔截彈在大氣層內(nèi)飛行過程受力非常復(fù)雜，而且初始段攔截彈的飛行彈道一般是按預(yù)先設(shè)定的控制程序飛行，初始段的彈道難以按照理論規(guī)劃彈道飛行，實(shí)際的飛行彈道與攔截彈初始段的飛行程序有關(guān)。

3 仿真計(jì)算

3.1 模型驗(yàn)證

在大地慣性坐標(biāo)系下，假設(shè)時(shí)刻獲取彈道導(dǎo)彈發(fā)動(dòng)機(jī)關(guān)機(jī)后點(diǎn)的狀態(tài)位置矢量為(2.945 36×10m，5.998 132×10m，2.967 413×10m)，速度矢量為(-2.659 2×10m/s，2.172 93×10m/s，3.326 04×10m/s)。假設(shè)攔截彈的發(fā)射點(diǎn)坐標(biāo)為(-2.882 326×10m，4.541 819×10m，3.416 752×10m)。以時(shí)刻為起始時(shí)間，則目標(biāo)落地時(shí)刻為1 011 s，下面對攔截彈道規(guī)劃模型進(jìn)行仿真驗(yàn)證。其中圖2所示為飛行時(shí)間分別為500 s、600 s、700 s，不同遭遇點(diǎn)位置下攔截彈道的規(guī)劃結(jié)果。圖3所示為飛行時(shí)間分別為500 s、600 s、700 s，相同遭遇點(diǎn)位置下攔截彈道的規(guī)劃結(jié)果。圖4為飛行時(shí)間為500 s，不同遭遇點(diǎn)位置下攔截彈道的規(guī)劃結(jié)果。

圖2 不同攔截點(diǎn)不同飛行時(shí)間仿真結(jié)果Fig.2 Simulation results of different interception points and different flight durations

圖3 相同攔截點(diǎn)不同飛行時(shí)間仿真結(jié)果Fig.3 Simulation results of different flight durations at the same interception point

圖4 不同攔截點(diǎn)相同飛行時(shí)間仿真結(jié)果Fig.4 Simulation results of the same flight duration at different intercept points

假設(shè)攔截彈主動(dòng)段結(jié)束時(shí)刻即發(fā)動(dòng)機(jī)關(guān)機(jī)點(diǎn)距離地面高度為100 km，遭遇點(diǎn)為目標(biāo)從起始時(shí)刻飛行700 s的位置(-1 487 113 m，5 535 938 m，4 280 217 m)，設(shè)攔截彈飛行時(shí)間為300～700 s，仿真計(jì)算攔截彈主動(dòng)段發(fā)動(dòng)機(jī)關(guān)機(jī)點(diǎn)以及與目標(biāo)遭遇點(diǎn)的速度隨飛行時(shí)間的變化曲線如圖5所示。由圖5可知，在同一遭遇點(diǎn)，攔截彈發(fā)射時(shí)間越早，所需的關(guān)機(jī)點(diǎn)速度越小。假設(shè)攔截彈的速度受到最大發(fā)射能量約束，可進(jìn)一步約束飛行時(shí)間的取值范圍。

圖5 攔截彈發(fā)動(dòng)機(jī)關(guān)機(jī)點(diǎn)、遭遇點(diǎn)速度隨tf的變化曲線Fig.5 Curves of variation of velocity at the interceptor’s shutdown point and impact point with tf

利用圖2中的仿真參數(shù)，根據(jù)模型可計(jì)算出關(guān)機(jī)點(diǎn)處的位置矢量和速度矢量，具體參數(shù)如表1所示，以表1中攔截彈關(guān)機(jī)點(diǎn)的參數(shù)可作為主動(dòng)段飛行程序的設(shè)計(jì)指標(biāo)。同樣也可解得遭遇點(diǎn)處目標(biāo)與攔截彈的運(yùn)動(dòng)參數(shù)，具體的參數(shù)如表2所示。對比分析表2中的遭遇點(diǎn)參數(shù)可知，在遭遇點(diǎn)1、遭遇點(diǎn)2、遭遇點(diǎn)3處攔截彈與目標(biāo)的距離偏差分別為|Δ|=925 m，|Δ|=956 m，|Δ|=1 092 m。一般情況攔截彈的末制導(dǎo)機(jī)動(dòng)范圍可達(dá)幾千米，在攔截彈末制導(dǎo)控制下攔截彈有能力修正飛行彈道的偏差，實(shí)現(xiàn)攔截彈與目標(biāo)的遭遇，所以規(guī)劃的攔截彈道在理論上能夠滿足中段反導(dǎo)攔截作戰(zhàn)要求。

表1 攔截彈關(guān)機(jī)點(diǎn)位置和速度

表2 遭遇點(diǎn)參數(shù)

3.2 攔截解算

依據(jù)本文構(gòu)建的攔截彈道模型，在相關(guān)約束條件下對反導(dǎo)攔截作戰(zhàn)的整個(gè)作戰(zhàn)過程進(jìn)行仿真分析。攔截彈的發(fā)射點(diǎn)位置為(-2.882 326×10m，4.541 819×10m，3.416 752×10m)，以美軍地基中段防御系統(tǒng)的相關(guān)攔截能力為參照，假設(shè)攔截彈的攔截低界為200 km、攔截高界為2 000 km、攔截斜距為4 000 km，最大速度增量為5 200 m/s，下面對攔截彈的發(fā)射窗口、可攔截弧段進(jìn)行仿真計(jì)算。仿真結(jié)果如圖6、圖7所示。其中，圖6為攔截彈在發(fā)射時(shí)間和飛行時(shí)間兩個(gè)維度上的發(fā)射窗口，藍(lán)色區(qū)域?yàn)閿r截窗口，表示攔截彈在相應(yīng)的飛行時(shí)間和發(fā)射時(shí)間約束下發(fā)射可實(shí)現(xiàn)對目標(biāo)的攔截；圖7為在攔截彈在攔截能力約束下對彈道導(dǎo)彈目標(biāo)的可攔截弧段解算結(jié)果。根據(jù)攔截彈攔截能力約束下的發(fā)射窗口和可攔截弧段參數(shù)，可進(jìn)一步指導(dǎo)修正攔截彈的設(shè)計(jì)。

圖6 攔截彈二維發(fā)射窗口Fig.6 Two-dimensional launch window of the interceptor

圖7 可攔截弧段解算結(jié)果Fig.7 Interceptable arc calculation results

4 結(jié)論

1)建立了反彈道導(dǎo)彈中段攔截彈道規(guī)劃模型，并通過仿真計(jì)算驗(yàn)證了彈道規(guī)劃模型的有效性。

2)所建立的攔截彈彈道模型能夠在發(fā)射時(shí)間和飛行轉(zhuǎn)移時(shí)間維度上對彈道導(dǎo)彈攔截窗口進(jìn)行解算，并且可以依據(jù)攔截彈的攔截能力對彈道導(dǎo)彈的可攔截弧段進(jìn)行計(jì)算。

3)該模型為攔截彈道規(guī)劃設(shè)計(jì)提供了一種思路和方法，對中段攔截彈道規(guī)劃設(shè)計(jì)具有一定的參考和借鑒意義。

[1] 高昂，胡延蘇. 基于概率模型的彈道導(dǎo)彈防御系統(tǒng)攔截方案優(yōu)化與評估[J].兵工學(xué)報(bào)，2016，37(2): 379-384.

GAO A，HU Y S. The optimization and evaluation of ballistic missile defense system interception scheme based on probability model[J]. Acta Armamentarii, 2016, 37(2): 379-384. (in Chinese)

[2] 張毅，楊輝耀，李俊莉.彈道導(dǎo)彈彈道學(xué)[M].長沙:國防科技大學(xué)出版社，1999.

ZHANG Y，YANG H Y，LI J L. Ballistic missile ballistics[M].Changsha: National University of Defense Technology Press，1999. (in Chinese)

[3] 劉瑤，張占月，王磊. 中段反導(dǎo)火力單元優(yōu)化部署研究[J].火力與指揮控制，2018，43(5):157-161.

LIU Y，ZHANG Z Y，WANG L. Research on the optimal deployment of the mid-section antimissile firepower unit[J]. Firepower and Command Control，2018，43(5): 157-161. (in Chinese)

[4] 楊光宇，陳國光，王捷，等. 基于PSO-SQP算法的空地導(dǎo)彈彈道優(yōu)化[J].兵器裝備工程學(xué)報(bào), 2021, 42(5): 153-157.

YANG G Y, CHEN G G, WANG J, et al. Trajectory optimization of air to ground missile based on PSO-SQP algorithm[J]. Journal of Ordnance Equipment Engineering, 2021, 42(5): 153-157. (in Chinese)

[5] 張浩，張奕群，張鵬飛. 攔截主動(dòng)防御目標(biāo)的微分對策制導(dǎo)律[J].系統(tǒng)工程與電子技術(shù), 2021, 43(5): 1335-1345.

ZHANG H, ZHANG Y Q, ZHANG P F. Differential game guidance law for intercepting the active defense target[J].Systems Engineering and Electronic Technology, 2021, 43(5): 1335-1345. (in Chinese)

[6] 柴華.基于可達(dá)集的大氣層外導(dǎo)彈防御系統(tǒng)不確定性交戰(zhàn)分析[D].長沙: 國防科技大學(xué)，2015.

CHAI H. Analysis on uncertain combat of exoatmospheric missile defense system based on reachable set[D]. Changsha: University of Defense Science and Technology，2015. (in Chinese)

[7] 王莎，劉進(jìn)忙，李超，等.中段高精度彈道預(yù)測方法研究[J].科學(xué)技術(shù)與工程，2013，13(23):6683-6687.

WANG S，LIU J M，LI C，et al. Research on high precision trajectory prediction method in midcourse[J].Science，Technology and Engineering，2013，13(23):6683-6687.(in Chinese)

[8] 楊少春，劉永蘭，劉進(jìn)忙.采用分段思想的中段彈道預(yù)測算法[J].空軍工程大學(xué)學(xué)報(bào)(自然科學(xué)版)，2014，15(3):33-36.

YANG S H，LIU Y L，LIU J M. The middle trajectory prediction algorithm using segmented thinking[J]. Journal of Air Force Engineering University (Natural Science Edition), 2014, 15(3): 33-36. (in Chinese)

[9] 周啟航，張劉，霍明英.彈道導(dǎo)彈中段突防彈道設(shè)計(jì)與驗(yàn)證[J].光學(xué)精密工程，2015，23(9): 2645-2655.

ZHOU Q H，ZHANG L，HUO M Y. Design and verification of ballistic missile midcourse penetration trajectory[J]. Optical Precision Engineering, 2015, 23(9): 2645-2655. (in Chinese)

[10] 方群，李新國，朱戰(zhàn)霞，等.航天飛行動(dòng)力學(xué)[M].西安:西北工業(yè)大學(xué)出版社，2015.

FANG Q, LI X G, ZHU Z X, et al. Space flight dynamics [M]. Xi’an: Northwestern Polytechnic University Press, 2015. (in Chinese)

[11] BATTIN R H. Lambert’s problem revisited[J]. AIAA Journal，1977，15(5):705-713．

[12] 柴華，仲明，梁彥剛.采用狀態(tài)轉(zhuǎn)移矩陣的攔截中段制導(dǎo)方法[J].國防科技大學(xué)學(xué)報(bào)，2015，37(4):137-142.

CHAI H，ZHONG M，LIANG Y G. Intercepting midcourse guidance method using state transition matrix[J].Journal of National University of Defense Science and Technology，2015，37(4): 137-142. (in Chinese)

[13] 吳啟星，張為華.彈道導(dǎo)彈中段防御的攔截窗口分析[J].飛行力學(xué)，2005，23(2):56-58.

WU Q X，ZHANG W H. Analysis of intercepting window of ballistic missile midcourse defense [J]. Flight Mechanics，2005，23(2):56-58. (in Chinese)

[14] AVANZINI G，PALMAS A，VELLUTINI E．Solution of low-thrust Lambert problem with perturbative expansions of equinoctial elements[J]. Journal of Guidance, Control, and Dynamics, 2015, 38(8): 1585-1601.

[15] SCHUMACHERJR P W, SABOL C, HIGGINSON C, et al. Uncertain Lambert problem[J]. Journal of Guidance，Control, and Dynamics, 2015, 38(7): 1573-1584．

[16] 徐利民，張濤，陶佳偉. 能量最優(yōu)與燃料最優(yōu)Lambert交會(huì)問題[J].北京航空航天大學(xué)學(xué)報(bào)，2018，44(9):1888-1893．

XU L M，ZHANG T，TAO J W. Energy optimal and fuel optimal Lambert rendezvous problem[J]. Journal of Beijing University of Aeronautics and Astronautics，2018，44(9):1888-1893. (in Chinese)

[17] WEN C，ZHAO Y，SHI P．Derivative analysis and algorithm modification of transverse-eccentricity-based Lambert problem[J].Journal of Guidance, Control, and Dynamics, 2014, 37(4): 1195-1201.

[18] 茅永興，何晶.航天器軌道確定的單位矢量法[M].北京:國防工業(yè)出版社，2009.

MAO Y X, HE J. Unit vector method for spacecraft orbit determination[M]. Beijing: National Defense Industry Press，2009.(in Chinese)