扰动引力场中机动突防导弹落点精确预报方法

王磊，周祥，赵卫虎，程先哲，郑伟; WANG Lei; ZHOU Xiang; ZHAO Weihu; CHENG Xianzhe; ZHENG Wei

2025年6月9日 0:41 星期一

网刊加载中。。。

使用Chrome浏览器效果最佳，继续浏览，你可能不会看到最佳的展示效果，

确定继续浏览么?

复制成功，请在其他浏览器进行阅读

扰动引力场中机动突防导弹落点精确预报方法 PDF

- ORCID：
王磊 ¹
✉
- ORCID：
周祥 ²
- ORCID：
赵卫虎 ¹
- ORCID：
程先哲 ¹
- ORCID：
郑伟 ²

1. 国防科技大学信息通信学院, 武汉430000； 2. 国防科技大学空天科学学院, 长沙410073

中图分类号： V412.1

最近更新：2024-12-18

DOI：10.16356/j.1005-2615.2024.06.012

摘要

针对弹道导弹大机动突防后精确制导面临的落点预测需求，提出了一种考虑高阶扰动引力影响的导弹落点解析预报模型。将落点预报问题分解为标准落点预报和落点偏差预报两部分，标准落点预报由二体椭圆轨道理论解析求解，落点偏差预报通过构建的状态空间摄动模型进行求解。基于球谐函数换极法建立高阶扰动引力矢量在轨道柱坐标系中的表达式，并推导得到由F函数和G函数等两类核函数组成的落点偏差预报解析模型以及F函数和G函数的递推公式。该模型无需射前准备工作，相对已有方法具有使用灵活、鲁棒性好等特点。数值仿真结果表明本文提出的落点预测模型残差均值为11.2 m，相对误差小于0.1%，预测耗时小于100 ms，能够为制导算法设计提供支撑，具有一定的工程应用价值。

关键词

弹道导弹; 机动突防; 落点预测; 解析摄动理论; 扰动引力场

机动突防是弹道导弹应对反导拦截系统的有效手段之一。相对传统基于射前规划的突防模式，基于弹载感知系统、智能决策算法和姿轨控动力系统的弹道突防技术能够更加适应高动态对抗态势，是当前和今后导弹突防技术发展的主要趋势。但这种突防模式对弹载制导系统提出了更高的要求，由于突防时机以及机动次数的不确定性，使得弹上制导系统需要具备仅依赖当前状态在线解算制导指令的能力^［

1］。

闭路制导是一种根据目标状态及导弹当前飞行状态，按照控制泛函实时解算制导指令的制导方法，从原理上看，满足机动突防导弹的应用需求。但经典闭路制导算法为降低弹上计算时间，通常引入虚拟目标点的概念以修正地球非球型引力、大气阻力等摄动因素的影响，而且由于弹上计算机性能有限，虚拟目标点一般在射前进行计算并装订到弹上，因此在实际应用中也存在着一定局限性。为克服这一不足，一些学者通过建立考虑地球非球型引力、大气阻力等摄动因素的误差传播解析模型，为实现虚拟目标点的弹上求解提供了可行途径^［

2］。

除了推导解析模型外，还有学者尝试采用神经网络模型^［

3⁃4］建立导弹飞行状态与虚拟目标点的映射关系，在一定应用范围内取得较好的效果。但为了确保足够的拟合精度，同时顾及模型训练的复杂度，神经网络模型在构建时通常会将模型输入的状态空间限定在较小的范围内，一定程度上弱化了模型的泛化能力和对极端情况的适应能力。

从可靠性和可行性的角度考虑，建立弹道落点的解析预测模型仍然是当前解决大范围机动突防导弹高精度制导问题的最有效途径。但当前已有的解析模型仍然存在着如下不足：（1）大多数落点偏差解析预报模型仅考虑至J2项引力的影响，模型精度不足；（2）文献［

5⁃6］分别推导了更高阶引力摄动项的影响，但文献［6］的模型存在计算复杂度高，计算效率低的问题，而文献［5］提出的模型需要依赖一种沿飞行弹道的非球型引力重构模型，使得其对大范围机动导弹的适应能力不足。针对上述问题，本文在文献［5］的基础上，结合球谐函数换极思想，推导出一种无需射前准备工作的导弹落点解析预报模型。

1 基础模型

1.1　导弹落点预测问题描述

在航天动力学领域，落点预测问题可归结为经典初值问题，即给定飞行器当前状态 $S_{0}$ ，预测其未来某一时刻状态 $S_{T}$ ，预测的关键是找到 $S_{0}$ 到 $S_{T}$ 的映射关系。

弹道导弹落点预报可形式化的描述为

f : S_{0} \to S_{T}

（1）

式中： $S_{0} \in R^{6 \times 1}$ 表示导弹当前时刻运动状态的集合， $S_{T} \in R^{3 \times 1}$ 表示落点位置参数的集合。

1.2　轨道柱坐标系中的落点精确预测模型

轨道柱坐标系有助于描述空间无动力飞行器的运动状态，尤其便于分析摄动因素对飞行器的影响。本文后续推导的解析模型均以此为基础，其定义如图1所示。其中，导弹状态 $S_{0}$ 参数确定的二体轨道为标准椭圆弹道，标准椭圆弹道所在的平面为标准弹道平面，坐标原点 $O$ 位于地心。 $r$ 轴位于标准弹道平面内且沿地心矢径方向； $β$ 轴同样位于标准弹道平面内，垂直于 $r$ 轴且朝向导弹运动的方向； $z$ 轴垂直于标准弹道平面，且与 $r$ 轴和 $β$ 轴组成右手坐标系。

图1 轨道柱坐标系定义示意图

Fig.1 Definition of orbit cylindrical coordinate frame

在轨道柱坐标系中，选取 $v = [v_{r}, v_{β}, v_{z}]^{T}$ ， $x = {[r, z]}^{T}$ 以及地心角 $ξ$ 组成状态变量，建立导弹运动方程为

\{\begin{array}{l} {\dot{v}}_{r} = \frac{v_{β}^{2}}{r} - \frac{μ}{r^{2}} c o s^{3} η + δ a_{r} & \dot{r} = v_{r} \\ {\dot{v}}_{β} = - \frac{v_{r} v_{β}}{r} + δ a_{β} & \dot{ξ} = \frac{v_{β}}{r} \\ {\dot{v}}_{z} = - \frac{μ}{r^{2}} s i n η c o s^{2} η + δ a_{z} & \dot{z} = v_{z} \end{array}

（2）

式中： $μ$ 表示地球引力常数， $η = a r c t a n (z / r)$ ， $δ a_{r}$ 、 $δ a_{β}$ 和 $δ a_{z}$ 表示在轨道柱坐标系3个方向上除地球中心引力外的摄动加速度。

若摄动力已知，采用式（2）可精确获得导弹落点位置。另外，由式（2）的结构可知，导弹的实际落点位置 $S_{T}$ 可以表示为两部分之和，即标准落点 $S_{T}^{*}$ 和由摄动因素引起的落点偏差 $Δ S_{T}^{}$ 。用形式化语言可表示为 $f = f_{1} + f_{2}$ ，且

f_{1} : S_{0}^{} \to S_{T}^{*}

（3）

f_{2} : S_{0}^{} \to Δ S_{T}^{}

（4）

式中函数映射 $f_{1}$ 可通过二体轨道理论解析求解，不再赘述。本文重点聚焦函数 $f_{2}$ 的解析求解。

1.3　导弹落点偏差预测的状态空间摄动模型

状态空间摄动模型最早由我国科学家任萱教授于20世纪80年代初提出^［

7］，是一种求解摄动模型解析解的通用框架和理论，近几年国防科技大学课题组在这方面进行了更加深入系统的研究，提出一系列模型算法^{［参考文献 8⁃13}8⁃13］。国内一些其他学者也探讨了状态空间摄动理论在弹道预报方面的应用问题^{［参考文献 14⁃15}14⁃15］。运用状态空间摄动理论求解导弹落点偏差的基本思想如图2所示。首先基于变分原理，在轨道柱坐标系中建立以轨道偏差

δ X

为状态量，以真近点角

ξ

为自变量的摄动微分方程，保留一阶量；其次推导该一阶微分方程的状态转移矩阵

Φ (t, τ)

，得到积分形式的摄动方程半解析解；而后对摄动力进行函数变换，将其转换为自变量的函数，并代入摄动方程半解析解进行函数积分，得到完整解析解。

图2 状态空间摄动理论框架

Fig.2 Theoretical framework for state space perturbation

定义导弹状态偏差为

Y = [Δ v_{r}, Δ r, Δ v_{β}, Δ t, Δ v_{z} {, Δ z]}^{T}

（5）

基于状态空间摄动理论，可建立关于 $Y$ 的积分方程

\{\begin{array}{l} Δ v_{r} = \frac{1}{\sqrt[]{μ p_{ξ}}} \int_{0}^{β} r_{ξ}^{2} [λ_{11} δ g_{r} + λ_{13} δ g_{β}] d ξ \\ Δ r = \frac{1}{\sqrt[]{μ p_{ξ}}} \int_{0}^{β} r_{ξ}^{2} [λ_{21} δ g_{r} + λ_{23} δ g_{β}] d ξ \\ Δ v_{β} = \frac{1}{\sqrt[]{μ p_{ξ}}} \int_{0}^{β} r_{ξ}^{2} [λ_{31} δ g_{r} + λ_{33} δ g_{β}] d ξ \\ Δ t = \frac{1}{\sqrt[]{μ p_{ξ}}} \int_{0}^{β} r_{ξ}^{2} [λ_{41} δ g_{r} + λ_{43} δ g_{β}] d ξ \\ Δ v_{z} = \frac{1}{\sqrt[]{μ p_{ξ}}} \int_{0}^{β} r_{ξ}^{2} λ_{55} δ g_{z} d ξ \\ Δ z = \frac{1}{\sqrt[]{μ p_{ξ}}} \int_{0}^{β} r_{ξ}^{2} λ_{65} δ g_{z} d ξ \end{array}

（6）

式中： $p_{ξ}$ 表示标准二体轨道半通径，为常量； $r_{ξ}$ 表示标准二体轨道地心距； $λ_{i j}$ 为状态转移矩阵对应的元素，其中 $i, j = 1,2, \dots, 6$ ，其具体表达式可详见文献［

5］；

δ g_{r}

，

δ g_{β}

和

δ g_{z}

表示轨道柱坐标系中高阶扰动引力矢量三分量。文献［5］通过引入沿飞行弹道的扰动引力重构模型，将高阶扰动引力拟合为分段3阶多项式，并基于式（6）推导得到了考虑高阶扰动引力摄动的落点偏差解析预测模型。

2 改进的导弹落点偏差解析预测模型

尽管已有模型^［

5］在求解高阶扰动引力场中导弹落点预报方面表现出较好的性能，但该方法需要在导弹发射前基于标准弹道构建沿飞行弹道的扰动引力重构模型，显然无法应用于导弹大范围机动突防的情况。本节介绍一种摄动模型的直接推导方法。

2.1　改进的半解析模型

根据式（6）可知，摄动模型完整解析解的推导首先需要将摄动量表示为自变量的函数。在地球物理与大地测量领域，高阶地球引力模型一般采用球谐函数表示，自变量为地心距、地心经度和地心纬度，显然无法直接代入式（6）进行解析解的推导。文献［

8］提出的球谐函数换极模型设计了一种可将球谐函数自变量变为以任意大圆为“新赤道平面”的坐标系中，为解决本文的问题提供了思路。将换极后的球谐函数表达式代入式（6），且保留与落点偏差计算相关的项，可得

\begin{array}{l} Δ r = \sqrt[]{\frac{μ}{p_{ξ}}} \int_{0}^{β} [- λ_{21} \sum_{n = 2}^{s} (n + 1) {(\frac{a_{e}}{r_{ξ}})}^{n} \sum_{k = 0}^{n} ({\bar{A}}_{n k} c o s k ξ + \\ {\bar{B}}_{n k} s i n k ξ) \cdot {\bar{P}}_{n k} (0) - λ_{23} \sum_{n = 2}^{s} {(\frac{a_{e}}{r_{ξ}})}^{n} \cdot \sum_{k = 0}^{n} k \cdot ({\bar{A}}_{n k} s i n k ξ - {\bar{B}}_{n k} c o s k ξ) {\bar{P}}_{n k} (0)] d ξ \end{array}

（7）

\begin{array}{l} Δ t = \sqrt[]{\frac{μ}{p_{ξ}}} \int_{0}^{β} [- λ_{41} \sum_{n = 2}^{s} (n + 1) {(\frac{a_{e}}{r_{ξ}})}^{n} \sum_{k = 0}^{n} ({\bar{A}}_{n k} c o s k ξ + \\ {\bar{B}}_{n k} s i n k ξ) \cdot {\bar{P}}_{n k} (0) - λ_{43} \sum_{n = 2}^{s} {(\frac{a_{e}}{r_{ξ}})}^{n} \cdot \sum_{k = 0}^{n} k \cdot ({\bar{A}}_{n k} s i n k ξ - {\bar{B}}_{n k} c o s k ξ) {\bar{P}}_{n k} (0)] d ξ \end{array}

（8）

\begin{array}{l} Δ z = \sqrt[]{\frac{μ}{p_{ξ}}} \int_{0}^{β} λ_{65} \sum_{n = 2}^{s} {(\frac{a_{e}}{r_{ξ}})}^{n} \sum_{k = 0}^{n} ({\bar{A}}_{n k} c o s k ξ + \\ {\bar{B}}_{n k} s i n k ξ) \cdot {\bar{P}}_{n k}^{'} (0) d ξ \end{array}

（9）

式中： $a_{e}$ 表示地球赤道平均半径； ${\bar{A}}_{n k}$ 和 ${\bar{B}}_{n k}$ 为换极后的正则化球谐系数； ${\bar{P}}_{n k} (0)$ 、 ${\bar{P}}_{n k}^{'} (0)$ 分别为缔合勒让德函数及其一阶导数； $n 、 k$ 分别为球谐的阶数和级数； $s$ 为截断阶数。由于 ${\bar{P}}_{n k} (0)$ 和 ${\bar{P}}_{n k}^{'} (0)$ 仅包 $ξ$ 一个自变量，因此式（7~9）原则上可积分得到完整解析解。

2.2　两类核函数的积分求解

通过对式（7~9）进行整理合并可知，求解以上问题的关键是处理以下2个核函数

\{\begin{matrix} F (n, m, p, q) = \int_{0}^{β} \frac{c o s (m ξ) c o s^{p} ξ s i n^{q} ξ}{r^{n}} d ξ \\ G (n, m, p, q) = \int_{0}^{β} \frac{s i n (m ξ) c o s^{p} ξ s i n^{q} ξ}{r^{n}} d ξ \end{matrix}

（10）

式中 $m$ 、 $n$ 、 $p$ 、 $q$ 均表示非负整数。

由于F函数和G函数中的各项均为 $ξ$ 的初等函数形式，因此以上积分理论上均可积分。若逐项推导，特别是当 $n$ 、 $m$ 、 $p$ 、 $q$ 较大时，积分表达式十分复杂。因此，寻求以上积分的递推或其他简化表达形式是解决问题的关键。

通过对式（10）进行深入剖析和推导，得到了一组递推公式如下：

（1）当 $n = p = q = 0$ 时

\begin{matrix} F (0, m, 0,0) = \int_{0}^{β} c o s m ξ d ξ = {[\frac{1}{m} s i n m ξ]}_{0}^{β} \\ G (0, m, 0,0) = \int_{0}^{β} s i n m ξ d ξ = {[- \frac{1}{m} c o s m ξ]}_{0}^{β} \end{matrix}

（11）

（2）当 $p = 1, n = q = 0$ 时

\begin{array}{l} F (0, m, 1,0) = \int_{0}^{β} c o s m ξ c o s ξ d ξ = \\ {[\frac{1}{2 m + 2} s i n (m + 1) ξ + \frac{1}{2 m - 2} s i n (m - 1) ξ]}_{0}^{β} \\ G (0, m, 1,0) = \int_{0}^{β} s i n m ξ c o s ξ d ξ = \\ {[- \frac{1}{2 m + 2} c o s (m + 1) ξ - \frac{1}{2 m - 2} c o s (m - 1) ξ]}_{0}^{β} \end{array}

（12）

（3）当 $p \geq 2, n = q = 0$ ，且 $m \neq p$ 时

\begin{matrix} F (0, m, p, 0) = η_{0,0} + η_{0,1} F (0, m, p - 2,0) \\ G (0, m, p, 0) = μ_{0,0} + μ_{0,1} G (0, m, p - 2,0) \end{matrix}

（13）

式中： $η_{0,0} = \frac{p}{p^{2} - m^{2}} c o s m ξ c o s^{p - 1} ξ s i n ξ - \frac{m}{p^{2} - m^{2}} \cdot$ $s i n m ξ c o s^{p} ξ$ ， $μ_{0,0} = \frac{p}{p^{2} - m^{2}} s i n m ξ c o s^{p - 1} ξ s i n ξ + \frac{m}{p^{2} - m^{2}} c o s m ξ c o s^{p} ξ$ ， $η_{0,1} = μ_{0,1} = \frac{p (p - 1)}{p^{2} - m^{2}}$ 。

（4）当 $p \geq 2, n = q = 0$ ，且 $m = p$ 时

\begin{array}{l} F (0, m, p, 0) = η_{0,0}^{'} + η_{0,1}^{'} F (0, m - 1, p + 1,0) \\ G (0, m, p, 0) = μ_{0,0}^{'} + μ_{0,1}^{'} G (0, m - 1, p + 1,0) \end{array}

（14）

式中： $η_{0,0}^{'} = \frac{1}{p + 1} s i n ((m - 1) ξ) c o s^{p + 1} ξ$ ， $μ_{0,0}^{'} = - \frac{1}{p + 1} c o s ((m - 1) ξ) c o s^{p + 1} ξ$ ， $η_{0,1}^{'} = μ_{0,1}^{'} = \frac{p - m + 2}{p + 1}$ 。

（5）当 $p \in N, n = 0, q = 1$ 时

\begin{matrix} F (0, m, p, 1) = η_{1,0} + η_{1,1} G (0, m, p + 1,0) \\ G (0, m, p, 1) = μ_{1,0} + μ_{1,1} F (0, m, p + 1,0) \end{matrix}

（15）

式中： $η_{1,0} = - \frac{1}{p + 1} c o s (m ξ) c o s^{p + 1} ξ$ ， $μ_{1,0} = - \frac{1}{p + 1} s i n (m ξ) c o s^{p + 1} ξ$ ， $η_{1,1} = - \frac{m}{p + 1}$ ， $μ_{1,1} = \frac{m}{p + 1}$ 。

（6）当 $p \in N, n = 0, q \geq 2$ 时

\begin{array}{l} F (0, m, p, q) = η_{2,0} + η_{2,1} F (0, m, p + 2, q - 2) + \\ η_{2,2} G (0, m, p + 1, q - 1) \\ G (0, m, p, q) = μ_{2,0} + μ_{2,1} G (0, m, p + 2, q - 2) + \\ μ_{2,2} F (0, m, p + 1, q - 1) \end{array}

（16）

式中： $η_{2,0} = - \frac{1}{p + 1} c o s m ξ c o s^{p + 1} ξ s i n^{q - 1} ξ$ ， $μ_{2,0} = - \frac{1}{p + 1} s i n m ξ c o s^{p + 1} ξ s i n^{q - 1} ξ$ ， $η_{2,1} = μ_{2,1} = \frac{q - 1}{p + 1}$ ， $η_{2,2} = μ_{2,2} = - \frac{m}{p + 1}$ 。

（7） $m$ 方向递推

\{\begin{array}{l} F (n, m, p, q) = F (n, m - 1, p + 1, q) - \\ G (n, m - 1, p, q + 1) \\ G (n, m, p, q) = G (n, m - 1, p + 1, q) + \\ F (n, m - 1, p, q + 1) \end{array}

（17）

（8） $n$ 方向递推

\{\begin{array}{l} F (n, m, p, q) = β_{0} F (n - 1, m, p, q) + β_{1} F (n - 1, \\ m, p + 1, q) + β_{2} F (n - 1, m, p, q + 1) \\ G (n, m, p, q) = β_{0} G (n - 1, m, p, q) + β_{1} G (n - 1, \\ m, p + 1, q) + β_{2} G (n - 1, m, p, q + 1) \end{array}

（18）

式中： $β_{0} = \frac{1}{p_{ξ}}$ ， $β_{1} = \frac{p_{ξ} - ‖r_{0}‖}{p_{ξ} ‖r_{0}‖}$ ， $β_{2} = - \frac{r_{0} \cdot v_{0}}{\sqrt[]{μ p_{ξ}} ‖r_{0}‖}$ ，且 $‖r_{0}‖$ 为预测起始点地心距， $r_{0} 、 v_{0}$ 分别为预测起始点位置矢量和速度矢量。

2.3　完整解析模型求解

经整理，落点偏差解析预测模型（式（7~9））的完整表达式为

\begin{array}{l} Δ r = \frac{1}{h} [r λ_{23} U_{n}]_{0}^{β} + K_{2} \sum_{k = 0}^{n} a_{e}^{n} {\bar{P}}_{n k} (0) \cdot [α_{20} {\bar{A}}_{n k} F (n - \\ 2, k, 0,0) + α_{20} {\bar{B}}_{n k} G (n - 2, k, 0,0) + \\ α_{21} {\bar{A}}_{n k} F (n - 2, k, 1,0) + α_{22} {\bar{A}}_{n k} F (n - 2, k, \\ 0,1) + α_{21} {\bar{B}}_{n k} G (n - 2, k, 1,0) + \\ α_{22} {\bar{B}}_{n k} G (n - 2, k, 0,1)] \end{array}

（19）

\begin{array}{l} Δ t = \frac{1}{h} [r λ_{43} U_{n}]_{0}^{β} + K_{4} \sum_{k = 0}^{n} a_{e}^{n} {\bar{P}}_{n k} (0) \cdot [α_{40} {\bar{A}}_{n k} F (n - \\ 2, k, 0,0) + α_{40} {\bar{B}}_{n k} G (n - 2, k, 0,0) + \\ α_{41} {\bar{A}}_{n k} F (n - 2, k, 1,0) + α_{42} {\bar{A}}_{n k} F (n - 2, k, \\ 0,1) + α_{41} {\bar{B}}_{n k} G (n - 2, k, 1,0) + \\ α_{42} {\bar{B}}_{n k} G (n - 2, k, 0,1)] \end{array}

（20）

\begin{array}{l} Δ z = \frac{r_{1} μ}{h^{2}} \sum_{k = 0}^{n} a_{e}^{n} {\bar{P}}_{n k}^{'} (0) \cdot [s i n β {\bar{A}}_{n k} F (n - 1, k, 1,0) - \\ c o s β {\bar{A}}_{n k} F (n - 1, k, 0,1) + s i n β {\bar{B}}_{n k} G (n - 1, \\ k, 1,0) - c o s β {\bar{B}}_{n k} G (n - 1, k, 0,1)] \end{array}

（21）

式中 $h$ 为标准二体轨道的动量矩，为常量。

K_{2} = \frac{r_{1}^{2}}{h^{2} p_{ξ}^{2}}

，

K_{4} = \frac{r_{1}^{2}}{h p_{ξ}^{2}}

\begin{array}{l} {\bar{P}}_{n k} (0) = \\ \{\begin{array}{l} \begin{matrix} 0 & m o d (n - k, 2) = 1 \end{matrix} \\ \begin{array}{l} {(- 1)}^{\frac{n - k}{2}} \sqrt[]{\frac{υ (2 n + 1) (n - k)! (n + k)!}{2^{2 n}}} l_{n k} \\ m o d (n - k, 2) = 0 \end{array} \end{array} \\ υ = \{\begin{matrix} 1 & k = 0 \\ 2 & k \neq 0 \end{matrix}, l_{n k} = \frac{1}{(\frac{n - k}{2})! (\frac{n + k}{2})!} \end{array}

\begin{array}{l} {\bar{P}}_{n k}^{'} (0) = \\ \{\begin{array}{l} \sqrt[]{\frac{n (n + 1)}{2}} {\bar{P}}_{n 1}^{'} (0) & k = 0 \\ \sqrt[]{(n - k) (n + k - 1)} {\bar{P}}_{n, k + 1}^{'} (0) & 1 \leq k \leq s \end{array} \end{array}

α_{20} = (1 - 2 n) e s i n f_{1}

\begin{array}{l} α_{21} = (1 - n) s i n β - n e^{2} s i n β + 2 (n - 1) \\ e s i n f_{0} - e^{2} s i n f_{1} s i n f_{0} \end{array}

\begin{array}{l} α_{22} = (n - 1) c o s β + n e^{2} c o s β + 2 (n - 1) \\ e c o s f_{0} - e^{2} c o s f_{1} c o s f_{0} \end{array}

α_{40} = (1 - 2 n)

α_{41} = (n - 1) e c o s f_{1} c o s β + 2 (n - 1) c o s β - n e c o s f_{0}

α_{42} = (n - 1) e c o s f_{1} s i n β + 2 (n - 1) s i n β + n e s i n f_{0}

式中 $e$ 为标准椭圆轨道的偏心率。

式（19~21）即为本文推导的落点偏差预测模型，可以看出，落点偏差的预测仅依赖预测点起始处的状态参数。同时由于采用了递推形式的核函数积分公式，本文推导的预测模型解析解与文献［

5］中模型相比，在形式上更加简洁。

3 仿真实验

本节通过数值仿真来验证上述导弹落点解析预报模型的精度。分别采用数值积分法、文献［

5］的方法以及本文推导的模型预报导弹落点，并将数值积分法的计算结果作为评价解析模型精度的基准。仿真中扰动引力采用72×72阶EGM2008模型进行计算。导弹突防时刻的位置参数如表1所示。

表1 预报点状态参数

Table 1 Parameters of the initial condition

变量	取值
高度 $H_{0}$ /km	700
速度 $v_{0}$ /（m·s^-1）	7 000
速度倾角 $θ_{0}$ /（°）	30
方位角 $α_{0}$ /（°）	0~180
经度 $λ_{0}$ /（°）	0
纬度 $φ_{0}$ /（°）	0

仿真结果如图3~5所示。图（3（a），4（a），5（a））分别表示由高精度数值模型计算出的72×72阶地球不规则引力引起的导弹落点偏差。可以看出，对射程超过10 000 km的导弹，地球不规则引力可引起接近20 km的落点偏差。图（3（b），4（b），5（b））分别表示文献［

5］所提模型以及本文所提模型在用于求解导弹落点偏差时相对高精度数值解的残差。可以看出，文献［5］模型和本文所提模型径向残差绝对均值分别为6.3和5.9 m，侧向残差绝对均值分别为6.1和7.3 m，时间项残差绝对均值分别为0.005 8 s和0.005 3 s，精度相当。但文献［5］的方法需要在导弹发射前基于标准弹道构建地球扰动引力场重构模型，而本文所提方法则不需要在导弹发射前做任何计算，这正是本文方法的核心优势。对于大范围机动导弹，由于无法在射前明确其飞行轨迹，因而也无法构建覆盖其飞行轨迹的扰动引力重构模型，从而使得文献［5］所提方法在此类情形下不可用。可以看出，相对文献［5］的模型，本文方法能够适应导弹大范围机动变轨等情况，灵活性和鲁棒性更好。

图3 纵向位置预报误差

Fig.3 Prediction error of longitudinal position

图4 时间项预报误差

Fig.4 Prediction error of time term

图5 侧向位置预报误差

Fig.5 Prediction error of lateral position

表2所示为3种计算模型的计算效率对比。仿真时发射方位角固定为40°，其他参数与表1保持一致。仿真结果表明，本文所提方法计算效率与所考虑的扰动引力阶次呈正相关，阶次越高，耗时越长。当考虑72阶引力时，解析模型的计算耗时为0.096 8 s。与本文方法不一样，文献［

5］方法的效率与引力阶次关系不大，基本维持在7 ms左右，这是由于文献［5］的模型通过有限元网格来拟合扰动引力矢量，重构模型阶次是固定的，不随所考虑的扰动引力阶次而改变。

表2 不同模型在不同引力阶次条件下计算效率的对比结果

Table 2 Comparative results of computing efficiency of different models under different gravity order conditions

方法	EGM2×2		EGM8×8		EGM24×24		EGM48×48		EGM72×72
方法	精度/m	耗时/s	精度/m	耗时/s	精度/m	耗时/s	精度/m	耗时/s	精度/m	耗时/s
数值法	9.5E-4	0.092	9.5E-4	0.663	9.5E-4	4.21	9.5E-4	15.84	9.5E-4	37.65
文献［5］方法	15.210	0.007 3	15.644 2	0.007 2	15.983 6	0.007 2	15.223 3	0.007 9	15.272 2	0.007 8
本文方法	14.670	1.56E-4	14.681 2	5.91E-4	14.679 6	0.006 2	14.661 3	0.031 2	14.672 2	0.096 8

如前文所述，本文所提方法的核心优势是不依赖标准轨迹，无需在射前做任何计算，但这种优势也是有代价的，正如表2中所示，该方法在计算效率方面存在一定程度的降低，但随着弹载计算机性能的提升，通过适当提高模型复杂度来进一步提升方法的应用灵活性和鲁棒性有较明显的工程应用价值。

4 结论

本文针对弹道导弹大机动突防后精确制导面临的落点预测需求，基于状态空间摄动理论提出了一种考虑高阶扰动引力影响的导弹落点解析预报模型。仿真结果显示了该模型在计算精度和效率方面的优良性能。特别地，相比已有方法，该方法仅依靠导弹当前状态即可实现考虑高阶扰动引力影响的导弹落点偏差精确预报，无需射前准备工作。这一特性使得该方法能够较好地适应导弹大范围机动突防制导的应用需求。需要说明的是，当前的仿真分析均是在实验室仿真环境下进行的，下一步仍须验证并优化模型在弹载计算环境下的性能，从而为未来先进制导系统设计与开发提供支撑。

参考文献

唐学海，王军，肖志明. 弹道式再入航天器落点预报技术［J］.南京航空航天大学学报， 2019， 51（S）： 145-148. [百度学术]

TANG Xuehai， WANG Jun， XIAO Zhijun. Landing point prediction ballistic reentry spacecraft［J］. Journal of Nanjing University of Aeronautics & Astronautics， 2019， 51（S）： 145-148. [百度学术]

刘伯阳，唐学海，孙中兴. 跳跃式返回航天器落点预报技术研究［J］. 南京航空航天大学学报， 2021， 53（S）： 45-50. [百度学术]

LIU Boyang， TANG Xuehai， SUN Zhongxing. Research on prediction of skip entry spacecraft’s landing site［J］. Journal of Nanjing University of Aeronautics & Astronautics， 2021， 53（S）： 45-50. [百度学术]

任乐亮，鲜勇，李少朋，等．基于改进二阶优化器并行学习的弹道导弹神经网络落点预测方法［J］. 航空学报， 2023， 44（14）： 327964. [百度学术]

REN Leliang， XIAN Yong， LI Shaopeng， et al. A neural network model for impact point prediction of ballistic missile based on improved second-order optimizer with parallel learning［J］. Acta Aeronautica et Astronautica Sinica， 2023， 44（14）： 327964. [百度学术]

PENG H， BAI X. Improving orbit prediction accuracy through supervised machine learning［J］. Advances in Space Research， 2018， 61（10）： 2628-2646. [百度学术]

王磊. 基于状态空间摄动法的战略导弹弹道快速预报与制导方法研究［D］. 长沙：国防科技大学， 2018. [百度学术]

WANG Lei. Fast trajectory prediction and guidance algorithm for strategic missiles based on state space perturbation method［D］. Changsha： National University of Defense Technology， 2018. [百度学术]

VINTI J P. Theory of an accurate intermediary orbit for satellite astronomy［J］. Journal of Research of the National Bureau of Standards， 1961， 65（3）： 169-201. [百度学术]

任萱. 自由飞行时摄动方程的状态转移矩阵的解析解［J］. 中国空间科学技术， 1983， 1： 1-16. [百度学术]

REN Xuan. Analytical solution of state transition matrix of perturbation equation in free flight phase［J］. Chinese Space Science and Technology， 1983， 1： 1-16. [百度学术]

WANG Lei， LIN Zhiqiang， PENG Yongjun. Progress on error propagation and correction of long-range rockets in disturbing gravity field［J］. Journal of Aerospace Engineering， 2022， 236（5）： 1693-1704. [百度学术]

WANG Lei， ZHENG Wei. Approximate analytic solution of motion for spacecraft around an oblate planet［J］.Journal of Aerospace Engineering， 2019， 233（4）： 1423-1434. [百度学术]

ZHANG H， LI B. Velocity-to-be-gained deorbit guidance law using state space perturbation method［J］. Journal of Aerospace Engineering， 2018， 31（2）： 04017099. [百度学术]

WANG L， ZHENG W. A novel feedback correction strategy for dosed⁃loop guidance of long⁃range rockets［c］//Proceedings of International Conference on Guidance， Navigation and Control. Tianjin：［s.n.］2020. [百度学术]

ZHOU Huan， AN Xueying， ZHENG Wei， et al. Fast local representation of gravity anomaly along flight trajectory［J］. Journal of Aerospace Engineering， 2016， 230（7）： 1201-1215. [百度学术]

WANG Lei， ZHENG Wei. Fast trajectory deviation prediction algorithm with arbitrary degree/order gravitational potential［J］. Journal of Guidance， Control and Dynamics， 2019， 42（7）： 1595-1600. [百度学术]

DENG Zexiao， LIU Luhua. State-space perturbation analytical solution for the dynamics of launch-vehicle boost phase［J］. Journal of Aerospace Engineering， 2024， 37（5）： 04024049. [百度学术]

WU Yansheng， WANG Zongqiang， ZHANG Bing. Comprehensive compensation method for the influence of disturbing gravity field on long-range rocket guidance computing［J］. Chinese Journal of Aeronautics， 2022， 35（3）： 408-418. [百度学术]

请使用 Firefox、Chrome、IE10、IE11、360极速模式、搜狗极速模式、QQ极速模式等浏览器，其他浏览器不建议使用!

扰动引力场中机动突防导弹落点精确预报方法 PDF

摘要

关键词

1 基础模型

1.1 导弹落点预测问题描述

1.2 轨道柱坐标系中的落点精确预测模型

1.3 导弹落点偏差预测的状态空间摄动模型