熊少锋,魏明英,赵明元,熊 华,王卫红,周本春
(1. 北京电子工程总体研究所,北京 100854;2. 中国航天科工集团防御技术研究院,北京 100854;3. 北京航空航天大学自动化科学与电气工程学院,北京 100191)
拦截战术弹道导弹(Tactical ballistic missile,TBM)的制导律设计是一个富有挑战性的课题[1]。TBM为代表的高速机动目标,其再入大气层时的速度可以达到3 km/s,远大于拦截弹的速度,传统的尾追制导方式难以让导弹击中目标。为了提高拦截概率,需采用迎面逆轨拦截的打击方式。理想的逆轨拦截有两个优点:1)导弹速度和目标速度均与弹目视线重合,故无论是导弹速度还是目标速度大小的变化,都无需导弹付出额外的法向过载,只是拦截时间变长或变短而已;2)逆轨拦截时,导弹速度与弹目视线重合,导弹的法向过载与弹目视线垂直,其将全部被用于改变视线转率,效率最高。在制导末端,为使导弹具有较充裕的机动能力,要求在制导末端需用过载尽可能小。
综合而言,逆轨拦截高速机动目标对制导律提出了攻击角度和过载两方面的约束[2-3]。常见的制导律设计方法有比例制导、滑模制导和最优制导,下面分别对其进行介绍。
针对固定目标,文献[4]和[5]分别设计了虚拟目标比例导引律和偏置比例导引律。其中,文献[4]有效地避免了末端控制量发散的问题,而文献[5]通过设计一种盲区控制方案以减少终点处的法向过载。针对高速运动目标,文献[6]将碰撞角约束作为反馈指令,使用正比例系数和负比例系数分别设计了逆轨和顺轨拦截目标的角度约束制导律。
滑模控制,算法简单且对干扰具有鲁棒性,因此被广泛应用于制导律设计。文献[7]和[8]分别采用扩张状态观测器和基于Kriging的无模型有限长单位冲激响应滤波器估计目标机动的加速度,设计了带角度约束的滑模制导律。考虑导弹自动驾驶仪的动态特性,文献[9-10]基于滑模控制设计了角度约束制导律。
最优控制因为对各种约束条件的考虑,在制导律设计方面有着较大的应用潜力[11-12]。针对固定目标,文献[13]和[14]分别应用Schwarz不等式和间接Gauss伪谱法设计了带角度约束的最优制导律。文献[15]考虑中制导末端的视场角约束,提出一种新的虚拟目标设置方法,并基于此设计了满足视场角约束的最优中制导律。文献[16]应用终端投影变换法对高阶制导模型进行降阶处理,根据线性二次型控制和微分对策理论设计了具有角度约束的制导律。在文献[16]的基础上,文献[17]考虑时变的加速度界限,设计了一种打击机动目标的角度约束制导律。考虑导弹自动驾驶仪的二阶动态特性,文献[1]应用高斯伪谱法设计了多约束条件下的三维最优中制导律,且将剩余飞行时间n次幂的倒数作为控制权重系数,使得过载在制导末端趋于一个小值。该制导律需要使用算法软件包SNOPT求解非线性规划问题。
综上所述文献在设计制导律时,一般假设导弹速度恒定不变。但在实际工程中,导弹由于受到发动机推力、空气动力和重力的作用,其速度是时变的,而且变化范围较大,从Ma1以下上升到Ma5~6,故需要研究考虑导弹速度时变的制导律设计。
针对固定目标,文献[18]将剩余飞行距离引入到指标函数,设计了带导引头视场角约束的制导律,使得过载在制导末端趋于小值,在仿真时考虑了导弹速度是时变的,效果良好。针对匀速运动目标,文献[19]应用间接高斯伪谱法设计了考虑导弹速度时变的角度约束制导律。对于低速运动的坦克目标,文献[20]将导弹速度的变化率和目标机动视为总扰动,利用干扰观测器对其进行估计,基于滑模控制设计了角度约束制导律。针对机动目标,文献[21]应用伴随变换法设计了考虑导弹速度时变的最优制导律。文献[22]在文献[16]的基础上,应用线性时变系统理论设计了考虑导弹速度时变的角度约束最优制导律。文献[21-22]的制导律均包含了导弹速度从当前时刻到终端时刻的积分这一项。为解决导弹未来速度未知的问题,文献[21]忽略诱导阻力,给出了预测导弹未来速度的一次性方法;文献[22]则给出了预测导弹未来速度的多步迭代方法。文献[21-22]预测导弹未来速度的方法均需要耗费大量的计算时间。
在传统的最优制导律中[16,22],脱靶量和角度约束的权重系数被设为一个较大的常数,这会使得导弹在初始阶段的过载指令较大。但很多时候,导弹在初始阶段速度较低,还处于加速过程中,难以提供较大的可用过载。为解决这个问题,文献[23]对双曲正切函数进行改造,设计了双曲正切函数的变种,它随时间从零开始增大,且增大的速度越来越快,并将其设为脱靶量和角度约束的权重系数。这么做的思想是,在制导的初始阶段,导弹离目标较远,不需要对脱靶量和攻击角度施加太强的约束,导弹可以“自由”朝目标飞去;随着时间的推移,导弹逐渐接近目标,这时就得逐步强化对脱靶量和攻击角度的约束。双曲正切函数的变种符合这一要求特性,降低导弹在初始阶段的过载需求。
上述文献里的制导律基本都是二维的,但实际中导弹和目标是在三维空间里飞行的,需要研究三维制导律设计。文献[24]基于李亚普诺夫稳定性理论设计了控制撞击角度和时间的三维导引律,其打击对象为固定或者慢速移动目标。
针对导弹逆轨拦截来袭高速机动目标,本文将[23]的制导律从二维扩展到三维。攻击角度定义为导弹速度矢量和目标速度矢量之间的夹角,在三维情况下,难以像文献[25]直接得到期望视线角和攻击角度之间的解析关系式,增加了设计难度。从便于工程应用的角度出发,对三维制导方程进行简化处理,将三维制导模型投影到垂直和横向两个二维平面。在投影过程中,导弹在三维空间里的弹道倾角与其在二维垂直平面内的弹道倾角并不相等,两者之间存在一个等效关系,因此使得垂直平面内制导模型的建立较文献[23]中单纯的二维制导变得复杂,增加了三维制导律的设计难度。截至目前,作者未曾在公开发表的文献中看到类似的逆轨拦截机动目标的三维制导律推导过程,本文中的垂直平面和偏航平面内的制导方程都是作者原创推导的。本文设计的三维制导律具有解析表达式,可实时在线运行,满足工程中导弹拦截TBM目标时的逆轨要求,且导弹在刚开始加速阶段,过载指令不太大。
先定义一下惯性坐标系,惯性坐标系以导弹发射车为原点;x轴正向为发射车的车头方向;y轴在当地铅垂面内,向上为正;z轴在当地水平面内,且与x轴和y轴满足右手定则。
在惯性坐标系中,弹目拦截几何如图1所示:
图1 弹目拦截几何Fig.1 Missile target engagement geometry
在图1中,M和T分别表示导弹和目标,r表示弹目视线,qε和qβ分别表示视线倾角和视线偏角。飞行器(包括导弹和目标)在惯性坐标系中的速度如图2所示。
图2 飞行器速度矢量Fig.2 Velocity vector of aircraft
其中,v,γ和ψ分别表示飞行器的速度、弹道倾角和弹道偏角。导弹速度vm和目标速度vt在三维空间里的投影分别为
(1)
将三维空间分解成铅垂和水平两个二维平面,进而将三维制导律简化为两个二维制导律来设计。在铅垂面内设计制导律时,假设导弹的弹道偏角ψm(t)在t时刻保持瞬时恒定;在水平面内设计制导律时,假设导弹的弹道倾角γm(t)在t时刻保持瞬时恒定。这里需要强调的是,假设导弹的弹道倾角γm和弹道偏角ψm保持瞬时恒定,这样做仅仅只是为了便于将三维制导分解成两个二维制导,从而简化三维制导律的推导过程,但是在执行制导律时,依然还是用γm和ψm随时间变化的实时值。另外,一般来说,弹道倾角和弹道偏角相对于制导指令而言是慢时变量,所以推导过程中假设它们保持瞬时不变。稳定性方面,制导系统的状态变量收敛很快,将三维制导问题解耦成两个二维制导问题不会影响制导系统的稳定性,后面的仿真结果也表明解耦设计的制导律是稳定的。另外,还可以参考文献[26],介绍了三维制导律可以分解成两个二维制导律来设计,而不影响制导系统的稳定性。
在铅垂平面内,弹目相对运动如图3所示。
图3 oxy平面内的弹目拦截几何Fig.3 Missile target engagement geometry in oxy plane
弹目距离沿垂直于初始视线方向的分量y的动力学方程为:
(2)
对式(2)求导可得
(3)
由式(1)可得导弹速度vm和目标速度vt在oxy平面上的投影分别为
(4)
对式(4)求导可得
(5)
根据式(1)可得导弹和目标在oxy平面内的等效弹道倾角分别为
(6)
对式(6)求导可得
(7)
定义如下的两个变量
(8)
则式(7)可以表示为
(9)
将式(4)和式(5)及式(9)代入式(3)可得
(10)
定义如下两个变量为
(11)
将式(11)代入式(10)可得
sin(γmxy-q10)-ηmxyay
(12)
(13)
根据式(12)可建立如下的状态空间方程
(14)
在水平面内,弹目相对运动如图4所示。
图4 oxz平面内的弹目拦截几何Fig.4 Missile target engagement geometry in oxz plane
弹目距离沿垂直于初始视线方向的分量z的动力学方程为:
(15)
对式(15)求导可得
(16)
根据式(1)可得导弹速度vm和目标速度vt在oxz平面上的投影分别为
vmxz=vmcosγm,vtxz=vtcosγt
(17)
对式(17)求导可得
(18)
根据式(1)可得导弹和目标在水平面内的等效弹道偏角分别为
ψmxz=ψm,ψtxz=-ψt±π
(19)
对式(19)求导可得
(20)
将式(18)和(20)代入式(16)可得
(21)
定义一个新的变量ftxz为
(22)
将式(22)代入式(21)可得
cos(ψm-q20)az
(23)
(24)
根据式(23)可建立如下的状态空间方程
(25)
将导弹逆轨拦截机动目标的要求转化为对弹目速度交会角的约束,基于最优控制理论设计制导律。
铅垂面内的最优制导律,设计性能指标函数为
(26)
式中:tf表示终端时刻,tgo表示剩余飞行时间,a为脱靶量的权重系数,b为攻击角度的权重系数。针对式(26),构造Hamiltonian函数为
(27)
根据极小值原理可得最优控制解为
(28)
协态方程为
(29)
对式(29)求解微分方程可得
(30)
(31)
从式(31)可以看到,最优控制解里面包含了终端时刻的状态量x1tf和x4tf,接下来给出这两个终端状态量的求解过程。根据线性系统理论可得
(32)
(33)
求解式(33)可得
(34)
其中,
(35)
(36)
定义一个新的状态向量为
(37)
ZEM(Zero effort miss)和ZEAE(Zero effort angle error)分别为零控脱靶量和零控角误差,它们的数学表达式为
(38)
式中:
(39)
将式(34)代入式(32)右边的第二个部分,并从t到tf积分可得
(40)
式中:
(41)
将式(37)和(40)代入式(32)可得
(42)
求解方程(42)可得
(43)
(44)
将式(43)和式(44)代入式(31)可得最优制导律为
(45)
式中:
(46)
(47)
式(45)中的最优制导律包含了导弹速度从t时刻到tf时刻的积分,但导弹在未来时刻的速度是未知的,这给制导律的执行带来了困难,从工程易于执行的角度出发,简单地假设导弹速度在t时刻到tf时刻是不变的。基于这一假设,可直接得到下面的关系式
(48)
假设导弹和目标的相对位置偏离碰撞三角形较小,也即(q1-q10)<<1,则下列近似关系式成立
y=rxy(q1-q10)
(49)
对式(49)求导可得
(50)
一般来说,目标未来的加速度是难以预知的,因此假设在时间域[t,tf]内目标加速度是不变的,故可得到下列关系式
(51)
将式(48),式(50)~(51)代入式(38)可得零控脱靶量和零控角误差的简化表达式为
(52)
采用同样的方法,可推导设计水平面内的最优制导律,为节省篇幅,此处略去具体的推导过程。
考虑到前面所做的易于工程实施的简化处理,现将y/z向过载指令的计算公式演绎、综合如式(53)和式(54)所示。
(53)
(54)
双曲正切函数的表达式为
(55)
其随时间变化的曲线如图5所示。
图5 双曲正切函数Fig.5 Hyperbolic tangent function
设计双曲正切函数的变种函数为
(56)
令n=2,ζ=2,τ=3,双曲正切函数的变种函数随时间变化的曲线如图6所示。
采用双曲正切函数的变种设计指标函数如下
(57)
通过数字仿真以说明所设计最优制导律的性能。
假设导弹是轴对称的,其动力学方程为式(58)。
(58)
导弹动力学方程(58)阐明了导弹所受到的力,导弹速度的变化关系式,在后续仿真中将由式(58)实时计算导弹的速度,升力、阻力和侧向力。在式(58)中,T为发动机推力,X,Y和Z分别为阻力、升力和侧向力,α和β为攻角和侧滑角,ρ和S为空气密度和特征面积,K,CX0,CYα和CZβ为空气动力系数。
仿真参数:导弹的初始位置为(2000 m, 1000 m, -1000 m),目标的初始位置为(80000 m,75000 m, -10000 m),导弹的初始速度为269.2582 m/s,目标的速度、弹道倾角和弹道偏角的初值分别为3471.3 m/s,-43.7396°和-175.4261°,目标沿x,y和z三轴的加速度分别是-3g,-2g和-2g,导弹的初始质量为734 kg,燃料的燃烧时间和消耗速率分别为10 s和-24.4 kg/s,在燃料燃烧期间发动机的推力T为110 kN,燃料消耗完毕后发动机推力T为0 N,导弹的特征面积S为0.146 m2,重力加速度g=9.8 m/s2,攻角和侧滑角的最大值为35°,仿真步长为0.01 s,空气动力系数K,CX0,CYα和CZβ由表1查表获得。
表1 空气动力系数Table 1 Aerodynamic coefficient
为便于性能对比验证,将脱靶量和角度约束权重系数设为常值的传统最优制导律(Traditional optimal guidance law,TOGL)[16]也一并仿真。本文制导律和文献[16]的制导律,它们的仿真条件是相同的,导弹速度是时变的,其变化关系由式(58)确定。对于本文提出的制导律,控制参数(俯仰通道和偏航通道均采用相同的控制参数)为:N=1,n=2,μ=4,ζ=1。
(59)
TOGL的控制参数为:N=1,a=103,b=107。
导弹的初始弹道倾角和弹道偏角分别为55°和15°,仿真结果如图7~12所示。脱靶量为0.5633 m,导弹的末速度为1274.9 m/s。弹目速度纵向和横向交会角分别为-1.4799°和-179.9103°,与逆轨拦截要求的交会角0°和-180°相差-1.4799°和0.0897°,满足逆轨拦截对弹目交会角的约束。从图10~11看出,初始阶段导弹沿Y轴和Z轴的过载较小,分别为-4.1463gm/s2和2.9191gm/s2,在制导末段,Y向和Z向过载指令均可以趋于零附近的小值,满足制导律对过载指令在初始段和末段不能太大的约束。Y/Z向过载指令在第10 s时存在尖点的现象,这是因为第10 s的时候,导弹燃料消耗完毕,发动机推力由110 kN突降到0 N所致。第18 s时的尖点现象,则是由控制参数τ在第18 s时连续但不可导(不光滑)造成的。在图12中,由于发动机推力的作用,导弹前10 s处于加速过程,第10 s时,导弹速度达到峰值,后面由于发动机推力为零,导弹只受到空气阻力和重力的作用,其速度逐渐减小。
图7 三维空间内的弹目运动轨迹Fig.7 Trajectories of missile and target in threedimensional space
图8 X-Y平面内的弹目运动轨迹Fig.8 Trajectories of missile and target in X-Y plane
图9 X-Z平面内的弹目运动轨迹Fig.9 Trajectories of missile and target in X-Z plane
图10 Y向过载Fig.10 Guidance command along Y axis
图11 Z向过载Fig.11 Guidance command along Z axis
图12 导弹速度Fig.12 Velocity of missile
针对大气层内导弹逆轨拦截高速运动目标的三维制导律设计问题,本文考虑导弹速度时变的情况,运用最优控制理论和双曲正切函数设计了带角度约束的最优制导律。通过合理的简化处理,本文设计的制导律具有解析形式,可实时在线运行。仿真结果表明该最优制导律可以提升导弹的末速度,增强导弹的杀伤动能。