0 引言
具有快速弹道优化制导能力是现在弹道优化的重要发展方向之一。弹道优化方法可分为解析法和数值解法。解析法是针对实际物理模型,基于最优控制原理推导出解析解。文献[1]中利用庞特里亚金极小值原理,设计了滑翔制导炮弹以最大射程为指标的最优控制参数。解析法可以进行直接求解,对于简单系统较为有效,但难以求解复杂的多约束非线性优化问题。数值解法则是通过某种离散方法,对时间、状态和控制量进行离散,将非线性最优控制问题转化为有限参数规划问题后,采用合适的参数规划算法求解。直接打靶法是数值解法中常用的一种[2]。直接打靶法只离散控制量,其优化参数少,便于理解,但是需要使用高阶积分算法求解状态变量,计算量大,效率不高且对初始猜测要求很高[3]。伪谱法是另一种常用的数值解法。伪谱离散与序列二次规划法结合的伪谱法在求解非线性最优化问题时十分有效[4-5]。但是在目标函数及约束较复杂的情况下其求解速度会较慢,其主要原因是序列二次规划方法在每一次迭代中不仅要更新Jacobian矩阵,还需要更新Hessian矩阵。
凸优化方法由于收敛速度快,同时能考虑到各类约束,在轨迹优化、预测控制等最优控制问题中凸显优势,近年来逐渐受到重视。对于一个凸问题,局部最优解即为全局最优解,算法的最优性得到了保证[6]。但是一般的非线性问题具有高度的非凸性[7]。针对这类非线性规划问题,一般基于直接法求解,如序列二次规划等方法[8],也可以利用高斯伪谱法求解基于线性协方差的弹道优化问题[9]。然而,这类方法在使用中常存在初值敏感、收敛较慢的问题。因此,需要将原非线性规划问题转化为一系列近似的线性规划子问题,通过迭代逐步逼近原问题的解,这一过程称为序列凸优化[10]。基于序列凸优化,能够在多约束条件下以较快的速度收敛且对初值不敏感[11-12]。对于在线弹道优化,着重考虑其实时性与收敛性,采用多段序列凸优化,将航点约束转化为各段轨迹相连的平滑轨迹凸优化问题进行求解[13]。
一般制导律设计是最速接近为性能指标,在制导中会导致制导指令频繁调整,造成更多能量损耗。在终端条件约束变化不大时,中制导段合理优化弹道有利于提高飞行品质,从而有效的命中目标[14]。考虑到实际飞行存在扰动,控制律参数难以整定[15],设计了在飞行过程中利用可变的控制率参数进行循环计算得到当前路径的方法。基于某制导火箭弹滑翔弹道模型,以控制能量最优轨迹规划模型为目标函数,利用Radau伪谱法离散连续变量,线性凸化非线性动态方程,建立标准凸优化模型,并进行模拟仿真验算。
1 建模
考虑到滑翔段属于无动力飞行,为了提高远程打击能力以及末端存速,增加毁伤效果,需要制导火箭弹在滑翔的过程中尽量节省能量以实现一定的打击效能。本文中以控制能量最优为性能指标,优化制导火箭弹打击固定目标的滑翔弹道,并要求落角不大于给定上界值,末速度应不小于给定下界值。
1.1 制导火箭弹滑翔段弹道模型
针对方案弹道规划问题,在初步设计阶段,为了抓住问题的本质,便于研究滑翔弹的飞行弹道及飞行特性,可将滑翔弹视为一个可操纵的质点来研究,简化后的运动模型为
(1)
式中: V为速度,θ为弹道倾角。
FD、FL计算公式如下所示:
(2)
(3)
式(2)、式(3)中:s为特征面积,Cx为零升阻力系数,kc为诱导阻力系数,
为升力系数导数,α为攻角,β为侧滑角。
1.2 模型约束
假设制导火箭弹在滑翔起控点速度为V0,起控时的弹道倾角为θ0,起控点坐标为(x0,y0),目标点坐标为(xt,yt)。为了实现精确打击,需要制导火箭弹的末端位置与目标重合。
建立初始状态约束为
(4)
同时,为了保证最终的攻击效果,需要对末端存速进行约束。末速度约束下界为VfL,落角约束上界为θfU。
建立末端状态约束为
(5)
式(5)中: xT0为理想弹道的终点,xT为计算得到的弹道的终点,δ为弹着点的允许偏差率(δ≤100 m)。
在弹道模型中,控制变量为攻角α。考虑到制导火箭弹的控制性能及弹体的稳定性,同时避免弹体失速,需要对控制变量取值进行约束。
建立控制变量约束为
|α|≤αmax
(6)
1.3 优化模型
对于不同打击任务需求,可选择不同的优化指标。考虑到滑翔段无动力飞行,为了提高远程打击能力以及末端存速增加毁伤效果,需要制导火箭弹在滑翔的过程中尽量节省能量[5]。因此,对应的方案弹道的性能指标函数为
J=
α2dt
(7)
式(7)中:初始时间t0与终止时间tf是可变的,也是该优化问题的一组决策变量。
综上,制导火箭弹的滑翔段轨迹优化问题模型为
(8)
2 伪谱离散与凸化
2.1 模型一般化
为了方便后续公式推导,本文中先建立制导火箭弹滑翔控制段能量最优的一般化模型。一般化模型中状态变量、控制变量以及状态函数都以向量形式表示。
状态向量x∈Rn×1,控制变量u∈Rm×1,本文中,取n=4,m=1。
状态向量为:x=(v θ x y)T
控制向量为: u=α
于是动态方程可以写为
(9)
式(9)中: f∈Rn×1为状态函数向量。
目标函数可以写为
J=
(L2(u))2dt
(10)
式(10)中: L2(u)为控制向量u的Euclidean范数。状态向量 、控制向量和时间变量的上下界约束可以写为
xL≤x≤xU
(11)
uL≤u≤uU
(12)
tL≤t≤tU
(13)
式(11)—式(13)中: L表示约束下界,U表示约束上界。
综上所述,制导火箭弹的滑翔段轨迹优化问题模型为
(14)
2.2 Radau伪谱离散
伪谱法采用全局多项式插值,相对于传统的等距梯形离散方法,伪谱离散可以获得更高的精度[18]。Radau 伪谱法利用N阶Radau勒让德多项式的根为节点,其节点是非对称、非等间距的,采用它的根作为插值节点可以避免等距高阶多项式插值经常出现的龙格现象[19]。相比Gauss伪谱离散,Radau伪谱离散没有附加的积分约束,更为简洁[20]。
Radau多项式是N阶勒让德多项式lN和N-1阶勒让德多项式lN-1的差:
RN(τ)=lN(τ)-lN-1(τ), τ∈[-1,1]
N阶Radau多项式有N个根,这N个根包含τ0=-1, τN<1。可以看出,配置点包含初始点,但不包含终止点。
定义τ∈[-1,1]为伪谱时间, t=[t0,tf]为物理时间。伪谱时间与物理时间存在以下的映射关系:
(15)
Radau伪谱法的节点包含N阶Radau多项式的N个根和一个附加终止点,共计N+1个节点。
定义Li为第i个节点对应的拉格朗日多项式函数的基函数,为
于是,状态向量x的第j个分量xj(τ)可以用经过这N+1个节点的拉格朗日多项式函数来近似。
(16)
对两边求导,可得:
于是第j个状态变量在第ik个配置点τik处微分的值为
令
则
(17)
式(17)中: X∈R(N+1)×n为状矩阵,Xij为第j个状态变量在第i个离散节点的状态。
定义控制矩阵为U∈RN×m,则Uij为第j个控制变量在第i个离散点的控制量。微分矩阵D∈RN×(N+1),其每一行对应相应的配置点,每一列对应相应的拉格朗日基函数。
对伪谱映射关系求导,并代入动态方程可得:
定义伪谱域函数:
(18)
于是:
离散形式为
(19)
式(19)中: F∈RN×n
Fij= fτj(X,U,t0,tf)
(21)
式(20)中: 下标i对应第i个离散点。
可得:
DX=F(X, U, t0, tf)
(21)
可以看出,该最优化问题的决策变量包含:状态变量、控制变量、起始时间和结束时间。为了便于数值求解,基于标准数值规划模型,将所有决策变量转化为列向量。
定义决策变量列向量为
xdv=((vec(X))T (vec(U))T t0 tf)T
(22)
式(22)中,vec表示矩阵列变换, xdv∈Rndv×1, ndv为决策变量总数,表示为
ndv=(N+1)n+Nm+2
(23)
式(23)中,(N+1)n对应状态变量,Nm对应控制变量,2对应起始时间和结束时间。
Radau伪谱法的离散求积公式为
(24)
式(24)中: φi为被积函数在第i个配置点的值;ωi为对应配置点的积分权系数。
于是目标函数的离散形式为
(25)
式(25)中: ju表示第ju个控制变量。
引入辅助变量ηi建立二阶锥约束:
(26)
式(26)中:ηi为二次函数
的上镜图。可以看出,该二次函数为凸函数。
可知,凸函数上镜图的最优解即为该函数的最优解。在上式中表现为,当收敛到最优解时,小于等于号可以近似取等号,达到了对目标函数的无损转换:
(27)
2.3 线性凸化
为了应用凸优化方法得到优化弹道,需要保证目标函数与约束条件均是凸函数。因此需要对问题进行凸化,使得凸函数的性质得到展现,即凸函数最优解也是原问题的全局最优解。文中采用了序列凸优化来逐步逼近原问题的最优解。
设优化的初始状态为
考虑到末端时间tf未知,设虚拟时间τ∈[0,1],有t=t0+(tf-t0)τ
令(tf-t0)=ut,则有
可设优化的控制量序列为
u(k)=[α(k) β(k) ut(k)]T
对方程进行多元函数离散化处理得到:
A(k)(z(k)u(k))+B1(k)(z(k)α(k)u(k))·α+
B2(k)(z(k)α(k)u(k))·ut(k)+C(k)
式(28)中,
C(k)=f-A(k)z(k)-B1(k)α(k)-B2(k)ut(k)
3 弹道优化
在火箭弹飞行的时间T内,在每次间隔一定的时间段Tk后调用凸优化程序进行弹道优化。其中,间隔时间Tk是远小于飞行时间T的,同时为了保证能得到收敛的最优弹道,间隔时间Tk必须大于程序计算时间。
本节介绍本文中提出的弹道优化方法的实现方法,数值实现的步骤总结如下:
步骤1：选择以积分方法得到的弹道参数装定为优化的初始值。
步骤2：在滑翔段时对弹道进行模拟计算,求解凸优化的子问题,得到当前的优化弹道。
步骤3：自解算后飞行一段时间Tk后,再次求解凸问题,得到新的优化弹道。
步骤4：制导火箭弹落地击中目标,则计算结束,否则重复步骤3。
图1 弹道解算流程
Fig.1 Trajectory solution process
4 数值仿真与结果分析
本节以制导火箭弹滑翔段弹道规划问题为例,将本文中提出的方法进行仿真计算并与传统的伪谱法进行对比。仿真计算中使用非线性问题最优化通用工具GPOPS2进行计算。
其中,制导火箭弹的参数如表1所示,表2为决策变量的上下边界。
表1 制导火箭弹参数
Table 1 Parameters of guided rockets
m/kgkcCx0CyS/m2ρ/(kg·m-3)70350.510.00.01170.8723
表2 决策变量的上下边界
Table 2 Upper and lower boundaries of decision variables
v/(m·s-1)θψαβt/s上界/π/2π/266300下界200-π/2-π/2-6-60
在弹丸飞行过程中,风是弹道产生偏差的重要影响之一,文中以考虑风对弹道产生偏差,然后对弹道进行修正。随机风由于其风速大小不确定,不能写出其具体的模型公式,此处的随机风利用Matlab语言随机生成[19-20],如图2所示。风速w为大气相对地面的速度,通常将水平风速分解为纵风与横风:
wx=-wcos(αw-αN)
wz=-wsin(αw-αN)
图2 随机风模型曲线
Fig.2 Random wind model curve
其中,wx为纵风,wz为横风,aN为发射方向与正北方的夹角,aw为风向与正北方的夹角。
在火箭弹滑翔控制段内,假设在开始第1次计算的位置,存在一定误差,使此时的速度略大于预定速度而高度未达到预定值,通过本文提到的方法得到的优化弹道与原弹道的对比如图3所示。可以发现,通过文中的制导方法可以在不回到原弹道的情况下依旧到达了原定的目标,表示了在具有误差影响的情况下,该方法能够满足制导火箭弹的制导律要求,在制导火箭弹中的具有可行性的。
图3 弹道优化轨迹
Fig.3 Trajectory optimization
如上图所示,图4为攻角随时间的变化曲线,图5为弹道倾角随时间的变化曲线,图6为速度与时间的变化曲线。在其中,攻角、倾角和速度选取开始进行优化的位置为时间起点。分别在20、40、60 s时取当前参数重新优化弹道。从图中可以发现由于风对弹道的影响整体偏弱,得到攻角曲线等并不会产生突变。与原弹道的优化方案相比较可以发现,在满足毁伤任务的前提下,本文中的方法可以做到在弹丸飞行过程中进行弹道优化,修正弹道偏差,且其在理论上可以做到满足约束的能量最优。
图4 攻角随时间变化
Fig.4 Trajectory change over time
图5 弹道倾角随时间的变化
Fig.5 Trajectory inclination changes with time
图6 速度随时间的变化
Fig.6 Velocity changes with time
5 结论
本文中提出基于凸优化的制导方法,用于在线计算制导火箭弹滑翔段的优化弹道,计算结果表明:
1) 所提出的在线弹道优化方法能够做到在满足约束的情况下收敛得到结果,即在存在误差的情况下能回到原来的位置,证明了该制导方法是可行的。
2) 该方法能够以较快的速度收敛得到结果,并能够在飞行的过程中以当前的状态为基础进行计算,说明该方法可以用来进行在线优化计算。
[1] 史金光,王中原,孙洪辉,等.制导炮弹滑翔弹道优化设计方法研究[J].南京理工大学学报,2011,35(5):610-613,620.SHI Jinguang,WANG Zhongyuan,SUN Honghui,et al.Optimization design method for glide trajectory of guided projectiles[J].Journal of Nanjing University of Science and Technology,2011,35(5):610-613,620.
[2] 陈思远.助推-滑翔导弹弹道优化与制导方法研究[D].北京:北京理工大学,2016.LI Yu.Study of trajectory optimization and guidance method for boost-glide missile[D].Beijing:Beijing Institute of Technology,2016.
[3] FABIEN B C.Some tools for the direct solution of optimal control problems[J].Advances in Engineering Software,1998,29(1):45-61.
[4] 黄诘,张友安,刘永新,等．基于Radau伪谱方法的 轨迹优化[J].海军航空工程学院学报,2015,30(4):315-320.HUANG Jie,ZHANG Youan,LIU Yongxin,etal．Trajectory optimization based on Radau pseudospectral method[J].Journal of Naval Aeronautical and Astronautical University,2015,30(4):315-320.
[5] 陈琦.滑翔增程制导炮弹弹道优化及制导控制方法研究[D].南京:南京理工大学,2017.CHEN Qi.Investigation on trajectory optimization and guidance &control scheme for a type of gliding guided projectiles[D].Nanjing:Nanjing University of Science and Technology,2017.
[6] MAO Y Q, SZMUK M, ACIKMESE B.Successive convexification of non-convex optimal control problems and its convergence properties[C]//2016 IEEE 55th Conference on Decision and Control(CDC).Las Vegas,NV.IEEE,2016:3636-3641.
[7] LIU X F,LU P.Solving nonconvex optimal control problems by convex optimization[J].Journal of Guidance,Control,and Dynamics,2014,37(3):750-765.
[8] ALLEN L,BEIJERSBERGEN M W,SPREEUW R J C,et al.Orbital angular momentum of light and the transformation of Laguerre-Gaussian laser modes[J].Physical Review A,1992,45(11):8185-8189.
[9] MATSUMOTO N,ANDO T,INOUE T,et al.Generation of high-quality higher-order Laguerre-Gaussian beams using liquid-crystal-on-silicon spatial light modulators[J].Journal of the Optical Society of America A,2008,25(7):1642.
[10] LIU X F,SHEN Z J,LU P.Solving the maximum-crossrange problem via successive second-order cone programming with a line search[J].Aerospace Science and Technology,2015,47:10-20.
[11] 李臻,许冰青,李庆波,等.基于序列凸优化算法的飞行器轨迹规划[J].空天防御,2021(4):50-56,79.LI Zhen,XU Bingqing,LI Qingbo,etal.Gliding vehicle trajectory planning based on sequential convex optimization algorithm[J].Space and Space Defense,2021(4):50-56,79.
[12] 姜欢,赵建博,熊芬芬,等.多约束下基于凸优化的协同制导研究[C]//中国自动化学会控制理论专业委员会.第36届中国控制会议论文集.北京理工大学宇航学院;飞行器动力学与控制教育部重点实验室,2017:1302-1307.JIANG Huan,ZHAO Jianbo,XIONG Fenfen,et al.Research on cooperative guidance based on convex optimization under multiple constraints[C]//Committee on Control Theory of the Chinese Society of Automation.Proceedings of the 36th China Control Conference (D).College of Astronautics,Beijing Institute of Technology;Key Laboratory of Vehicle Dynamics and Control,Ministry of Education;,2017:1302-1307.
[13] 禹春梅.飞行器在线轨迹规划与制导控制方法研究[D].长沙:国防科学技术大学,2018.YU Chunmei.Research on online trajectory planning and guidance control for aircraft[D].Changsha:National Defense University of Science and Technology,2018.
[14] 高俊.无人飞行器控制参数离线优化与在线调整研究[D].长沙:国防科学技术大学,2015.GAO Jun.Research on control parameter offline optimize and online adjust of unmanned aerial vehicles[D].Changsha:National Defense University of Science and Technology,2015.
[15] 路钊.高超声速飞行器再入末段轨迹在线优化[D].哈尔滨:哈尔滨工业大学,2014.LU Zhao.Online trajectory optimization for the terminal stage of reentry hypersonic vehicles[D].Harbin:Harbin Institute of Technology,2014.
[16] 钱杏芳,林瑞雄,赵亚男．导弹飞行力学[M]．北京:北京理工大学出版社,2000.QIAN Xingfang,LIN Ruixiong,ZHAO Yanan． Flight mechanics of guided missile[M].Beijing:Beijing Institute of Technology Press, 2000.
[17] SAGLIANO M.Pseudospectral convex optimization for powered descent and landing[J].Journal of Guidance,Control,and Dynamics,2018,41(2):320-334.
[18] 李卫丽,严洪森,张维琴.风干扰下某型导弹的弹道仿真[J].计算机技术与发展,2011,21(1):246-249.
LI Weili,YAN Hongsen,ZHANG Weiqin.Missile’s trajectory simulation under the impact of wind[J].Computer Technology and Development,2011,21(1):246-249.
[19] 阚菲菲,白晨曦,张学锋.随机风对某型火箭布雷弹落点散布的影响[J].常州工学院学报,2016,29(3):27-32.KAN Feifei,BAI Chenxi,ZHANG Xuefeng.Shell distribution of a mine-laying bomb under the influence of stochastic wind[J].Journal of Changzhou Institute of Technology,2016,29(3):27-32.
[20] 李庆扬,王能超,易大义．数值分析[M].北京:清华大学出版社,2008.LI Qingyang,WANG Nengchao,YI Dayi． Numerical analysis[M].Beijing:TsingHua University Press,2008．
[21] GARG D,PATTERSON M A,FRANCOLIN C,et al．Direct trajectory optimization and costate estimation of finite-horizon and infinite-horizon optimal control problems using a Radau pseudospectral method[J].Computational Optimization and Applications,2011,49(2):335-358.
[22] 李宁波,雷虎民,周觐,等.反临拦截弹中制导弹道在线优化设计[J].控制与决策,2017,32(12):2195-2200.LI Ningbo,LEI Humin,ZHOU Jie,et al.Design of midcourse trajectory online optimization for near space interceptor[J].Control and Decision-making,2017,32 (12):2195-2200.
[23] 王俊波,曲鑫,任章.基于在线轨迹规划的混合再入制导方法[J].宇航学报,2012,33(9):1217-1224.WANG Junbo,QU Xin,REN Zhang.Hybrid reentry guidance based on the online trajectory planning[J].NASA,2012,33(9):1217-1224.
[24] 姜欢.凸优化方法在导弹轨迹优化与制导控制中的应用[D].北京:北京理工大学,2017.JIANG Huan.Application convex optimization theory in trajectory optimization,guidance and control of missile[D].Beijing:Beijing Institute of Technology,2017.