弹道助推段误差传播高精度求解方法

王 磊, 郑 伟, 周 祥

(国防科技大学航天科学与工程学院, 湖南 长沙 410073)

0 引 言

高命中精度始终是战略弹道导弹发展的主要目标之一。当前,普遍将影响导弹命中精度的因素划分为制导工具误差和制导方法误差两类[1]。制导工具误差主要包括加速度计和陀螺仪等惯性测量器件的误差,也是影响导弹落点精度最主要的误差源。近年来,随着新型陀螺仪和加速度计的出现,如环形激光陀螺、光纤陀螺、激光加速度计等,惯性测量器件的精度和稳定性逐步提升[2-5](陀螺陀螺漂移量可达10-6°/h,加速度计零偏稳定性可达1×10-5g);同时,针对不同惯性器件的误差补偿模型也愈加精确[6-7],从而使得制导工具误差对弹道导弹落点精度的影响逐渐减小。

相对地,制导方法误差的影响逐渐突显。该类误差主要包括制导算法误差、发动机后效冲量误差、再入段误差、地球引力场模型误差(通常也称为地球扰动引力场)等。其中,制导算法误差和再入段误差本身的影响量级较小,且可采用一定的补偿模型进行修正;后效冲量误差可通过改善发动机性能或增加末修级来减小或者消除;地球扰动引力场对远程弹道导弹落点偏差的影响最大可达1 km量级[8],是制导方法误差中最重要的组成部分,其对弹道影响的误差传播模型及其相应的补偿方法也成为当前亟待解决的课题。

地球扰动引力场对弹道导弹命中精度的影响主要体现在两方面,即导弹发射点垂线偏差引起的惯导系统定向误差与导弹飞行过程中的扰动引力引起的受力误差。将EGM2008模型[9-10](earth gravitational model 2008)等高阶引力场模型直接嵌入到弹载惯导系统中是实现对扰动引力场影响进行补偿最理想的方法,但对弹载计算机的执行效率及容量有极高的要求。为解决这个问题,研究人员提出了一系列地球引力场重构模型来快速逼近真实引力场[11-16],并成功应用于潜艇、舰船、飞机等低速载体的惯导系统中[17-23]。但对弹道导弹等高速飞行器而言,采用地球引力场重构模型对其惯导系统进行直接补偿的实际应用仍未见报道。

基于地球扰动引力场对弹道的影响量,反馈求解弹道特征量(如关机点偏导数)的修正值是实现对地球扰动引力影响进行补偿的另一种常用思路[8]。弹道摄动法是分析系统性扰动对弹道影响最常用的方法。摄动理论最初是由天文学家在分析月球、木星及土星等天体在受到太阳光压及其他引力体摄动后的运动特性的研究中发展起来的。欧拉、拉格朗日、高斯、拉普拉斯等人也对摄动理论的发展也作出了突出的贡献[24-25]。我国学者自20世纪80年代以来也对扰动引力场作用下的弹道摄动问题进行了系统而深入的研究。文献[26]在不考虑视速度偏差的假设下研究了导弹飞行过程中扰动引力对远程弹道导弹惯性制导系统及命中精度的影响特性;文献[27]则对弹道导弹被动段误差传播摄动方程进行了推导,并提出了解析的计算方法;文献[28]推导了垂线偏差、发射点大地经纬度偏差及发射方位角偏差等定位定向误差对弹道助推段关机点状态影响的摄动方程。这些方法在扰动引力场对导弹命中精度影响特性分析的应用中均取得了较好的效果,但由于没有考虑状态偏差与视加速度之间的耦合影响,其精度在一些初始条件下会有不足;另一方面,现有方法中考虑垂线偏差和飞行过程中扰动引力影响的摄动方程各自独立,无法对其影响进行统一分析。

基于此,本文推导了一种可同时考虑发射点垂线偏差和导弹飞行过程扰动引力影响的弹道助推段误差传播求解方法,该方法也考虑了状态偏差引起的视加速度偏差摄动项,求解精度更高。本文正文结构可分为3部分:第一,基于状态空间摄动法推导了弹道导弹在受到扰动引力场作用时的摄动方程,这是推导弹道误差传播模型的理论基础;第二,通过对摄动方程进行分析,深入探索了摄动源的构成及其影响机理;最后,通过数值仿真分析了扰动引力场对弹道助推段的影响量级,并对比验证了本文所提误差传播方法的计算效率和精度。

1 弹道助推段摄动方程

1.1 坐标系的定义

弹道助推段误差传播模型主要在发射惯性坐标系中描述。如图1所示,在考虑发射点垂线偏差的情况下,发射惯性坐标系可以有如下两种定义方式:

(1)以参考椭球法线为定向基准的标准发射惯性系On-xeyeze；

(2)以当地重垂线为定向基准的发射惯性系On-xayaza。

图1 发射惯性系定义示意图Fig.1 Definition of the launch inertial coordinate system

显然,坐标系On-xeyeze与On-xayaza之间的差别由垂线偏差决定,这二者之间的转换关系为

(1)

式中,ρe,ρa,ve,va分别为对应坐标系中的位置和速度矢量;EG为发射系到地心系的方向余弦阵,且

(2)

λ0,B0,A0为大地经纬度和大地方位角;λT,BT,AT为天文经纬度和天文方位角。根据大地测量理论可知它们与垂线偏差子午分量ξ和卯酉分量η之间的关系为

(3)

1.2 摄动方程

在椭球地球假设下,导弹助推段弹道运动方程为

(4)

考虑地球扰动引力场影响时,实际弹道运动方程为

(5)

令DG表示On-xayaza到On-xeyeze的坐标转换矩阵,根据式(1)有

(6)

将式(5)在On-xeyeze中表达,可得

(7)

式中,ge=DGga为在On-xeyeze中表达的真实引力矢量。

求解式(7)与式(4)的等时变分可得

(8)

(9)

式(8)即为扰动引力场影响下弹道助推段运动摄动微分方程。

2 误差传播模型推导

为方便分析,将式(8)改写为矩阵的形式

(10)

其中

(12)

对摄动方程式(10)进行分析可知,扰动引力场对弹道的影响主要包括4部分:

(2)由引力模型误差引起的受力项,即δge;

(3)由状态偏差引起的引力加速度耦合项,即A·[δve,δρe]T;

(4)由状态偏差引起的视加速度耦合项,即B·[δve,δρe]T。

2.1 几何项模型

(13)

数值仿真表明,推力T是产生视加速度的主要来源,忽略气动力及控制力等小量的影响,几何项模型可简化为

(14)

显然,几何项的量级与发射点的垂线偏差及推力的大小直接相关,垂线偏差越大、推力越大,几何项就越大,反之则越小。

2.2 受力项模型

地球的真实引力无法精确获得,通常只能通过模型不断逼近。求解地球扰动引力最常用的方法有点质量法和球谐函数法。前者只能针对固定区域进行求解,后者则可以实现全球范围内任意点的求解。本文仿真中主要采用后者进行计算,且采用EGM2008的模型系数。

扰动引力位在地球外部是一调和函数,满足拉普拉斯方程,基于分离变量法,可将其展开成球谐函数[1]

(15)

扰动引力位对位置的梯度即为扰动引力,即

(16)

由式(16)可知,扰动引力的计算是在地心系内完成,弹道计算时需要将其转换到发射惯性系中,即

(17)

2.3 引力加速度耦合项模型

引力加速度耦合项即是指由状态偏差引起的引力矢量偏差。该项重点在于求解发射惯性系中引力矢量关于速度、位置矢量的雅克比矩阵,根据文献[28]的推导,可知矩阵A的表达式为

(18)

其中

2.4 视加速度耦合项模型

视加速度耦合项即是指由状态偏差引起的视加速度矢量偏差。该项重点在于求解发射惯性系中视加速度矢量关于速度、位置矢量的雅克比矩阵。忽略控制力的影响,视加速度偏差可表示为

(19)

式中,R和T的表达式为

(20)

式中,Cx、Cy、Cz分别为阻力系数、升力系数和侧力系数;ρ为大气密度;v为飞行器相对大气的速度;Sm为弹体最大横截面积;m为弹体质量;ue为排气速度;Se为喷口截面积;pe为排气端面压力;pH为当地大气压力。

记Mv为R对速度矢量的偏导数,Mr为R对位置矢量的偏导数,具体表达式为

(21)

其中

M1=-Cxρvx;M2=-Cxρvy;M3=-Cxρvz

且有

(22)

记Nr为T对位置矢量的偏导数,具体为

(23)

可得视加速度相对导弹状态矢量的雅克比矩阵为

(24)

即

(25)

3 扰动引力场对弹道影响特性分析

3.1 扰动引力场对助推段弹道的影响特性分析

以民兵-3导弹总体参数、气动参数和我国标准大气参数为基础构建了导弹误差传播仿真系统。基于该仿真系统首先分析了摄动模型不同摄动项的影响量级,进而采用摄动法与弹道求差法对弹道误差传播特性进行仿真分析,并对比验证本文所提摄动法的计算效率及精度。

发射点参数及垂线偏差信息如表1所示,导弹飞行过程中的扰动引力采用360阶的球谐函数进行计算。

表1 导弹发射点参数

图2为4类摄动项的量级及其随时间的变化曲线。可以看出,在表1所示的初始条件下,几何项在所有摄动项中所占比例最大,最高可达500 mgal (1 mgal=10-5m/s2), 且其变化曲线与发动机推力曲线基本一致,符合前面的分析结论;受力项的量级仅次于几何项,且随着时间的增加而减小;视加速度耦合项随时间先增大后减小,当导弹飞出大气层外则减小为零;在整个助推段时间周期内,引力加速度耦合项的量级均在2 mgal以下,因而本文后续仿真中均忽略了该项。

图2 不同摄动项随时间变化曲线图Fig.2 Curves of different perturbation accelerations over time

下面对扰动引力场影响下的弹道误差传播特性进行仿真分析。弹道求差法的仿真流程可描述为:在由λ0,B0,A0确定的发射惯性系内和由λT,BT,AT确定的发射惯性系内采用同一时间步长分别进行弹道积分,且前者在积分时引力模型只考虑到J2项,后者则考虑到360阶扰动引力项。将两条弹道状态量(速度和位置)在同一坐标系内等时求差即可获得弹道真实误差传播结果。该方法的误差源仅为计算机本身累积误差,可以作为本文所提摄动模型精度评估的参照。为保证积分精度,弹道积分步长应足够小,通常使其与制导步长一致,一般取0.02 s。

与弹道求差法不同,摄动法的积分变量为导弹状态的偏差,其量级小且变化缓慢,因而可以采用较大的时间步长进行积分,其步长大小基于如下流程确定:

步骤1设置其步长初值Δt0=0.02 s,并按照当前步长采用摄动法求解助推段关机点状态偏差ΔX0;

步骤2按照0.01 s的增量,逐步增加步长值,同时求解相应步长下的关机点状态偏差ΔXk;

步骤3计算ΔXk相对ΔX0的误差百分比ε,并判断ε是否大于设定的相对误差容许值,若成立则终止循环,最终步长取为循环中上一步的步长,若不成立,则返回到步骤2。

本文仿真中,将相对误差容许值设置为0.001,并依据上述流程求得摄动法可采用的积分步长为1 s。同时,根据是否忽略视加速度耦合项将摄动法分为两类,即摄动法-1和摄动法-2。

图3为弹道求差法、摄动法-1和摄动法-2这3种方法求解弹道助推段误差传播的结果对比曲线,其中图3(a)～图3(c)依次为位置矢量偏差随时间的变化曲线,图3(d)～图3(f)为速度矢量偏差随时间的变化曲线。可以看出,本文提出的误差传播模型能较好反应地球扰动引力场影响下的弹道助推段误差传播特性;同时可知,若忽略摄动方程中的视加速度耦合项,会使误差传播模型的精度有明显的下降。

图3 沿助推段弹道的导弹状态偏差计算结果Fig.3 Calculation results of the state error in the boost phase

根据第2.4节的分析可知,视加速度耦合项反应了导弹状态偏差引起的气动力及发动机推力偏差,且主要引起气动力偏差。图2表明该摄动项主要在10～70 s的飞行时间内起作用,这是因为导弹在该飞行时间段内主要处于稠密大气层内,气动力对导弹状态的变化比较敏感;当导弹飞出大气层后,该摄动项则几乎不再起作用。图3表明,为了获得较高精度的误差传播模型,视加速度耦合项不应该被忽略。

3.2 扰动引力场对导弹射程的影响特性分析

为充分验证本文所提弹道助推段误差传播模型的效率、精度及其对不同初始条件的适应性,首先采用遍历法仿真分析不同导弹发射方位角条件下两种摄动模型相对弹道求差法的精度,其次采用蒙特卡罗打靶法分析不同垂线偏差条件下两种摄动模型相对弹道求差法的精度及其误差统计特性。

结果分析时将助推段关机点状态偏差导致的导弹射程偏差ΔL设为待评价指标,且ΔL的计算公式为

(26)

式中,等式右边行向量为导弹射程关于导弹关机点状态的偏导数,通常称为关机点射程偏导数,列向量为导弹关机点状态偏差。本文仿真中关机点射程偏导数取文献[27]中给定的值,即有

3.2.1 发射方位角遍历

仿真中,方位角由-180°～175°按照5°等间隔遍历,共计算72条弹道,其余发射点参数与表1所示一致。仿真结果如图4和图5所示。图4显示了弹道求差法和摄动法-2计算耗时对比曲线。显然,完成一次给定初始条件下弹道助推段误差传播仿真计算,求差法平均耗时1 025.3 s,而摄动法-2只需要23.4 s。其耗时之比为50.2,与这两种方法积分步长之比的倒数基本一致。实际上,采用球谐函数计算导弹飞行过程中的扰动引力是整个积分计算中最耗时的部分,几乎占据了总耗时的99%,且计算时间随着球谐阶数级数的增加呈幂级数增加。这也成为制约弹道求差法计算效率的最主要原因。

图4 计算时间对比Fig.4 Comparison of time-consuming between the two methods

图5显示了扰动引力场引起的射程偏差及摄动法-2的计算精度随着方位角的变化情况。

图5 不同方位角条件下的射程偏差求解精度对比Fig.5 Accuracy comparison of range error calculation under different azimuth conditions

3.2.2 垂线偏差随机打靶

设置发射点处垂线偏差的子午分量及卯酉分量都为满足正态分布的随机数,且其均值和方差均取5″。仿真中,在发射方位角为-160°、20°和140° 3种条件下分别进行500次的蒙特卡罗打靶,并统计摄动法-2计算所得射程偏差的相对误差百分比。

仿真结果如图6所示,图6(a),图6(c)和图6(e)分别为3种发射方位角条件下的摄动法-2相对误差百分比的散点图,图6(b),图6(d)和图6(f)则分别为其对应的直方图统计图。结果表明,3种方位角条件下摄动法-2的相对误差百分比都小于10%,且均值分别为1.834 0%,6.663 1%和2.525 4%。

图6 垂线偏差随机条件下的射程偏差求解精度对比Fig.6 Accuracy comparison of range error calculation under different deflection of the vertical conditions

4 结 论

地球扰动引力场对远程弹道导弹落点偏差的影响在百米到1 km量级。导弹发射前快速分析扰动引力场对拟飞行弹道的影响是实现扰动引力场影响精确补偿的有效手段之一。本文主要针对这一需求,基于状态空间摄动法提出一种可同时考虑发射点垂线偏差和导弹飞行过程扰动引力影响的弹道助推段误差传播高精度求解方法。仿真结果表明,考虑视加速度耦合项的误差传播模型在计算精度方面较不考虑视加速度耦合项的情况有明显提升;同时,相对弹道求差法而言,考虑视加速度耦合项的误差传播模型计算效率提升了约50倍,计算相对误差均小于10%。可以为实现扰动引力场影响快速补偿提供理论依据和方法支撑。