基于分布式并行伪谱- 神经网络算法的双脉冲导弹多阶段协同轨迹优化

2020-12-08刘超越张成

兵工学报 2020年10期

刘超越，张成

(北京理工大学飞行器动力学与控制教育部重点实验室, 北京 100081)

0 引言

近年来，导弹协同攻击因其在探测机动目标、抵御威胁、穿透防御系统等方面比单个导弹系统具有更好的性能而受到越来越多的关注[1-2]。协同轨迹规划问题可看作是非线性最优控制问题，目前适用于求解最优控制问题的伪谱法已经广泛地应用于无人机的协同轨迹规划问题中。Vera等[3]利用自适应伪谱法生成无人机的协同轨迹，并实现了避障，但其协同飞行时间是预先规定的，局限性较大。白瑞光等[4]利用伪谱法求解了以一种加权和为性能指标的两阶段轨迹规划问题，得到了3个无人机的最优协同航迹，但该研究对无人机的动力学模型进行了较大简化，将其飞行速度设为定值，控制量只改变无人机的飞行方向。

与无人机相比，导弹的弹体动力学具有大范围快速时变的特点：飞行高度可高至几十公里的高空；固体火箭发动机的工作使导弹的飞行速度从0到几马赫之间迅速变化，弹体动压随之变化几百倍；火箭弹在快速穿越大气层时，大气密度、温度等气象参数将剧烈变化；导弹无法像无人机一样作出悬停、后退等机动飞行。这些特性决定了导弹动力学方程的复杂性，使导弹协同轨迹规划问题的收敛域远小于无人机。因此，若将传统的无人机协同轨迹规划方法直接应用到导弹的协同轨迹规划问题中，求解效率很低。王芳等[5]采用伪谱法对3个导弹组成的编队进行了协同轨迹求解，但其以导弹的轴向加速度和法向加速度为控制变量，对动力学模型进行了较大的简化，并未考虑导弹飞行高度变化造成的声速、大气密度的变化，有较大的局限性。

以上协同轨迹规划算法均为基于集中式架构的求解方式，中心控制节点同时联立所有飞行器的动力学方程建立最优控制模型，所有动力学方程基于同一时间变量，因此不需要协调各飞行器的飞行时间。但随着协同飞行飞行器数量的增加，最优控制模型动力学方程的数量急剧增加，算法的求解时间呈指数式上升，因此，集中式方法研究模型中协同飞行的飞行器数量一般不能超过3个。为了提高求解导弹协同轨迹规划问题的效率，需要基于分布式架构进行协同轨迹规划问题的算法研究。

双脉冲导弹具有速度快、突袭性好、射程远等优点。在整个双脉冲导弹的飞行过程中，可以区分为几个飞行阶段。由于双脉冲发动机需要两次点火，发射前必须确定一些关键变量[6-7]，如二级推力的点火时间、初始射角，这些变量会显著影响弹道性能。程仙垒等[8]采用改进的粒子群优化算法，研究了动压、过载以及终端状态等多约束条件下的非连续助推弹道优化问题。但该研究中只优化了滑翔段，上升段都是按照飞行程序俯仰角飞行，并未优化全程弹道。关成启等[9]研究了高斯伪谱法在双脉冲导弹多阶段弹道优化问题中的应用，但二级推力的点火时间、制导启动时间无法随着飞行任务改变，不能很好地发挥双脉冲导弹的优势，对不同的飞行任务适应性较弱。明超等[10]以装备双脉冲发动机的导弹作为研究对象，以二级脉冲推力点火时间和飞行攻角为优化变量，采用高斯伪谱法对双脉冲导弹最大射程进行弹道优化研究，但只对导弹进行竖直平面内的建模，动力学模型简单。

同时，在进行空间发射、导弹试验和武器试验时，需要进行全面的靶场安全分析，以保护人员、资产和公众[11]。靶场安全模板工具包[12-13]结合地理空间信息(如资产位置和人口密度)，可以为制导武器的任务操作规划和安全分析提供伤亡和损害估计，但无法通过主动的轨迹规划将武器碎片(主要是推力段结束后分离的发动机壳体)投送至安全区域。若能提前主动预测分离发动机壳体的落点位置，则可对安全有更好的保证。

因此，对于考虑将分离的发动机壳体投送至安全区域的双脉冲导弹协同轨迹规划问题研究，必须保证双脉冲导弹的动力学模型充分反映真实飞行状态，并发挥其双脉冲推力特性的优势，同时提高多导弹协同轨迹规划的求解效率；其次，须对分离发动机壳体的落点位置进行主动预测和落点控制，以提供足够的安全保证。

上述双脉冲导弹的轨迹规划问题是一个终端时刻自由且带有路径约束的非线性最优控制问题。高斯伪谱法是求解最优控制问题数值方法的一种，它采用全局插值多项式在一系列拉格朗日高斯(LG)点上近似状态量和控制量。与一般的配点法相比，高斯伪谱法能够以较少的节点获得较高的求解精度，同时，根据高斯伪谱法共轭变量映射定理[14-15]，采用高斯伪谱法转化得到的非线性规划(NLP)问题的Karush-Kuhn-Tucker(KKT)条件与原最优控制问题1阶最优必要条件的离散形式具有一致性。因此，高斯伪谱方法因其高效、快速收敛和避免显式数值积分等优点，在求解最优控制问题方面得到了广泛的应用[16-17]。导弹分离发动机壳体的落点位置控制问题是一个预测控制问题。基于神经网络的机器学习方法，既可以处理复杂非线性、时变性和不确定性问题，又能够融合启发式经验，计算速度快，已经广泛地应用到了预测控制领域，实现了对制导精度和实时性的提升，从而提高制导系统的智能水平，如自主着陆与返航[18]、自主导航[19]、提高轨道预测精度[20]。

各导弹的飞行轨迹划分为发射段、爬升段、续航段和制导攻击段4个阶段。为了发挥双脉冲导弹可以根据飞行任务确定其二级推力点火时间的推力特性，提高协同飞行的性能，并提高最优协同轨迹的求解效率，提出一种基于分布式并行伪谱- 神经网络(DGPM-ANN)算法，求解考虑将分离的发动机壳体投送至安全区域的双脉冲导弹协同轨迹规划问题。针对伪谱法对初始猜测敏感的问题，在算法的迭代过程中，第1次迭代采用较少的离散点，快速计算最优轨迹；各导弹以得到的最优轨迹作为下一次迭代的参考轨迹进行求解，直至达到指定的精度；各导弹并行地独立求解各自的优化轨迹，并确定各导弹二级推力的点火时间；引入飞行时间下界约束协调飞行时间，以使各导弹的飞行时间随迭代收敛一致，保证协同飞行时间一致性；采用人工神经网络对分离发动机壳体的射程预测函数进行离线拟合，以对发动机落点位置进行快速、精确地预测，形成发动机落点位置约束。

1 导弹协同攻击问题描述

1.1 协同轨迹规划算法

协同轨迹规划算法可以分为集中式、分布式两种算法架构，分别如图1(a)和图1(b)所示。在集中式算法架构中，中心控制节点作为计算中心，需要同时联立所有导弹的动力学方程对协同轨迹进行规划，然后反馈至各导弹进行任务分配，不需要设计协调变量来保证飞行时间的一致性，导弹之间也不需要通信。在分布式算法架构中，各导弹均具备全局信息以及自主轨迹规划能力，相互之间可以通信，不存在主从关系，导弹加入或退出编队也更加灵活。

图1 轨迹规划算法架构Fig.1 Trajectory planning algorithm architecture

本文提出的分布式并行算法结构如图2所示。分布于各导弹的计算单元并行地独立求解各自的轨迹，实现底层的轨迹规划；在上层，导弹之间通过提出的飞行时间一致性协调策略，使各导弹的飞行时间随迭代收敛一致，以保证协同飞行时间一致性，从而实现上层的协调控制。

图2 分布式并行算法结构Fig.2 Decentralized parallel algorithm structure

1.2 轨迹阶段的划分

双脉冲导弹的飞行过程可以根据其发动机的推力特性划分为4个阶段：发射段、爬升段、续航段和制导攻击段，如图3所示。每个阶段的具体表述如下：

1)发射段。始于发射点，结束于一级推力的关机点，为无控飞行状态。一级推力可为导弹在飞行初期提供较高的飞行速度。发射段时长固定，为一级推力的工作时间，设为[t0,t1)。

图3 轨迹分段Fig.3 Division of trajectory

2)爬升段。始于一级推力的关机点，结束于二级推力的点火点。一级推力工作结束后，导弹继续爬升，进入爬升段。爬升段连接发射段与续航段，目的是将导弹续航段的开始时刻调整至最佳时机。因此，爬升段的工作时间不固定，设为[t1,t2)。

3)续航段。始于二级推力的点火点，结束于二级推力的关机点。在续航段，二级推力将导弹的姿态和速度调整至最佳状态，为适应战场环境、攻击目标做好充分的准备。续航段时长固定，为二级推力的工作时间，设为[t2,t3)。续航段结束时，发动机结束工作，与弹体分离。

4)制导攻击段。始于二级推力的关机点，结束于轨迹终点。在制导攻击段，导弹不断修正飞行路径，在满足约束的情况下，对地面目标进行攻击。制导攻击段的工作时间不固定，设为[t3,tf]。

1.3 导弹及分离发动机壳体的动力学方程

协同飞行的各导弹在地面坐标系下的质点动力学方程为

(1)

导弹在飞行过程中受到的空气阻力、升力和侧向力的计算公式为

(2)

式中：0.5ρv2为动压；ρ为大气密度；S为导弹参考横截面积；CD、CL和CC分别为阻力系数、升力系数和侧向力系数，均为攻角α、侧滑角β和马赫数Ma的函数。

对于本文中的导弹，阻力系数可看作是关于α和β的二次函数，升力系数可看作关于α的线性函数，侧向力系数可看作关于β的线性函数，表示为

(3)

式中：a1=-0.382;a2=2.415；a3=0.149 8；b1=-0.464;b2=6.466；c1=-0.464;c2=6.466.

标准大气密度ρ、声速a表示为

(4)

式中：T为大气温度，T=292.6-0.010 39y+4.497×10-7y2-6.639×10-12y3.

当续航段完成后，发动机与弹体分离，并以无控状态抛落至地面。分离发动机壳体的动力学方程为

(5)

式中：vb为发动机速度；mb为发动机质量；Fbx为发动机阻力，

(6)

ρb=1.225exp (-yb/7 201)，

(7)

Sb为发动机的参考面积，CbD为发动机的阻力系数；θb为弹道倾角；ψb为弹道偏角角；(xb,yb,zb)为分离发动机壳体的位置坐标。分离发动机壳体的状态向量可表示为xb=[vb,θb,ψb,xb,yb,zb]T.

1.4 约束条件

1)导弹在飞行的每个阶段都要满足相应的阶段约束，在整个轨迹的不同阶段，阶段约束有所不同，导弹的质量随着一级推力以及二级推力的工作发生变化。

发射段：

(8)

爬升段：

(9)

续航段：

(10)

制导攻击段：

(11)

式中：tp为发动机的工作时间，tp1、tp2分别为一级推力发动机和二级推力发动机的工作时间；P1和P2分别为一级推力发动机和二级推力发动机的推力；ξ1和ξ2分别为一级推力发动机和二级推力发动机的质量秒流量；umax和umin分别为控制量的可行域上、下限。

2)变量的边界约束。一般化表示为

Φmin≤Φ(x(t0),t0,x(tf),tf)≤Φmax，

(12)

式中：Φ、Φmax和Φmin分别为变量的一般化表示及其边界约束的上、下限；tf为飞行时间。

3)过载约束。为了确保导弹飞行过程中的结构安全，需要考虑过载约束：

(13)

式中：nL以及nLmax和nLmin分别为过载以及过载约束的上、下限。

4)分离发动机壳体的落点约束。为了保证人员及财产安全，需使分离发动机壳体掉落至安全位置。本文中安全区为位于地面的圆形区域，则对于分离发动机壳体有如下约束：

(14)

式中：(xbf,ybf,zbf)为分离发动机壳体的落点位置坐标；(xe,0,ze)为发动机的落地安全区中心位置坐标；Rs为安全区半径。

5)相邻阶段连接约束。可表示为

(15)

1.5 最优控制问题框架

导弹的多阶段协同轨迹规划问题可转化为最优控制问题P1：

(16)

式中：Ψ(·)为性能指标中的终端函数；L(·)为性能指标中的被积函数；Nm为导弹的个数。

2 多阶段最优控制问题离散化

高斯伪谱法将最优控制问题的状态量和控制量在一系列LG点上离散，并以离散点为节点构造拉格朗日插值多项式来逼近状态量和控制量。通过对全局插值多项式求导来近似状态量对时间的导数，从而将微分方程约束转换为一组代数约束，性能指标函数中的积分项由高斯求积公式近似，终端状态由初始状态加右端函数在整个过程的积分获得。经上述变换，可将最优控制问题转化为具有一系列代数约束的NLP问题：

1)时域变换。在每个飞行阶段中，t∈[t0,tf]被划分成K个时间区间[tk-1,tk)，k=1,2,…,K. 采用高斯伪谱法计算最优控制问题时，需要通过(17)式将每个时间区间转换到[-1,1]：

(17)

2)状态变量与控制变量的插值近似。高斯伪谱法需要在一系列离散点上对状态量和控制量进行全局插值多项式逼近。在[tk-1,tk)中，采用拉格朗日插值多项式作为基函数来近似状态量和控制量：

(18)

(19)

3)微分方程的约束转换为代数约束。对全局插值多项式(18)式微分后代入运动学微分方程中，并在LG点上离散，可得

(20)

4)约束条件离散化。边界约束(12)式表示为

(21)

在LG点对路径约束(13)式进行离散化处理，得

(22)

式中：C、Cmax和Cmin分别为路径约束的一般化表示及其约束的上、下限。

最后，最优控制问题P1可用高斯求积公式进行离散化表示为P2：

(23)

至此，多阶段连续最优控制问题被转化为离散参数的NLP问题。为了解决NLP问题，本文采用了开源的NLP求解器IPOPT. IPOPT采用内点算法，可以处理具有大量等式和不等式约束的大规模稀疏非凸问题。

3 射程预测函数和人工神经网络拟合

本节基于分离发动机壳体的动力学方程(5)式，构造了分离发动机壳体的射程预测函数，实现了发动机由分离位置到着陆位置的距离预测，然后建立数据集，利用人工神经网络对该函数进行了离线拟合。

3.1 分离发动机壳体的射程预测函数

以t3时刻导弹的部分状态量组成分离发动机壳体的积分初值[vb0,θb0,ψb0,0,yb0,0]，采用定步长数值积分算法对(5)式进行数值积分。取积分终止条件为yb≤0 km，得到积分终止时刻tbf的分离发动机壳体x轴、z轴坐标sbfx和sbfz：

(24)

令射程预测函数的输入向量为sb0=[vb0,θb0,ψb0,yb0]，输出向量为sbf1=[sbfx,sbfz]，基于(5)式构建分离发动机壳体射程预测函数，表示为

sbf1=fb(sb0).

(25)

3.2 人工神经网络训练模型

人工神经网络可以通过学习样本集自动提取规则，实现从输入到输出的映射，适用于拟合复杂非线性函数，因此，它可以根据训练样本预测输出结果。本文利用人工神经网络的拟合能力，对分离发动机壳体射程预测函数(25)式进行函数逼近，以提高预测速度，实现分离发动机壳体落点的实时预测，为导弹的在线轨迹规划提供条件。因此，需要训练得到一个神经网络拟合函数sbf2=Fb(sb0)，使得对于任意的sb0:

‖sbf1-sbf2‖≤δn，

(26)

式中:δn为拟合误差向量的上限，δn=0.5 km. 则分离发动机壳体的落点位置约束(14)式可表示为

‖(x(t3)+sbf2(1)-xe,z(t3)+sbf2(2)-ze)‖2≤Rs.

(27)

导弹的多阶段协同轨迹规划问题可以表示为P3：

(28)

本文采用前馈神经网络进行函数拟合，包括1个输入层、1个隐含层和1个输出层。图4展示了神经网络的结构，含有4个输入，隐含层有35个神经元，2个输出。

图4 人工神经网络结构Fig.4 Artificial neural network structure

在神经网络的训练阶段，为了能够训练得到足够精确的拟合函数，应该建立充分的训练数据集，多次考虑不同状态量初值来求解射程预测函数。各状态量初值是以在飞行期间分离发动机壳体的每个状态量可以假定的最大值和最小值来定义的。设分离发动机壳体质量mb=50 kg，阻力系数CbD=0.149，参考面积Sb=0.292 m2. 各状态量初值的上、下界以及步长由表1所示。

表1 状态量初值边界及步长Tab.1 Boundaries and step sizes of initial state values

以表1所示的状态量初值按相应步长进行网络划分，得到4维初值网格空间。以网格空间内的所有节点对应的值作为求解射程预测函数的初值，求解得到4 290 783个对应的预测射程向量，即产生4 290 783组输入输出样本。在完成预测射程向量数据采集步骤后，利用人工神经网络对采集到的数据集进行离线训练，生成学习模型，并将模型的权值和阈值表示的参数存储在数据库中，从而得到射程预测函数的人工神经网络拟合函数Fb(sb0)。图5显示了数据训练模型。

图5 数据集训练模型Fig.5 Dataset training model

定义拟合结果的误差百分比ex和ez分别为

(29)

为了分析人工神经网络对以不同积分步长求解射程预测函数得到的样本数据集进行拟合的精度以及拟合函数的求解速度，在积分步长分别为0.6 s、0.8 s和1.0 s的条件下，分别求解射程预测函数得到相应的数据集，然后利用人工神经网络分别进行离线拟合，拟合结果参数如表2所示。从表2中可以看出：随着积分步长的增大，射程预测函数的积分节点减少，平均计算时间也相应减小，但预测射程的精确度也势必有所减小；同时，3种情况下的神经网络拟合函数的平均计算时间以及拟合函数预测射程x轴方向分量和z轴方向分量的平均误差百分比几乎一致，且拟合函数的预测射程非常接近射程预测函数的预测射程。积分步长分别为0.6 s、0.8 s和1.0 s的情况下，采用神经网络拟合函数预测射程，预测速度分别提高了42倍、30倍和24倍。基于以上分析，本文以0.6 s作为积分步长生成样本数据集，可以得到较高的射程预测函数求解精度，并可以进一步得到较高的神经网络拟合函数预测精度和求解速度。

表2 不同积分步长对应的数据集人工神经网络拟合结果Tab.2 ANN fitting results of data sets corresponding todifferent integral step sizes

样本数据经人工神经网络拟合后，拟合结果x轴方向分量和z轴方向分量的误差百分比分别如图6(a)和图6(b)所示。经验证，拟合结果满足误差约束(26)式。由此可见，采用人工神经网络拟合预测模型，可以提供较高的精度和更快的求解速度。

图6 训练误差Fig.6 Training error

4 分布式并行伪谱- 神经网络算法

求解轨迹规划问题的流程图如图7所示。流程图左上方的3个方框表示了离线状态下分离发动机壳体的射程预测函数人工网络拟合过程。右上方两个方框表示将多阶段轨迹规划的最优控制问题转化为非线性规划问题的过程。底部的方框表示用结合序列高斯伪谱法和人工神经网络的DGPM-ANN算法对轨迹规划问题进行求解，从而得到原最优控制问题的解。

图7 轨迹规划求解流程图Fig.7 Flowchart of trajectory planning

本文提出的DGPM-ANN算法是一种基于迭代策略解耦形式的协同轨迹求解方法。在第1次迭代中，采用较少的离散点，快速计算最优轨迹；以得到的最优轨迹作为下一次迭代的参考轨迹进行求解，直至达到指定的精度；每一轮次的迭代中各导弹并行地独立求解轨迹规划问题，并利用神经网络拟合函数对分离发动机壳体的落点位置进行计算，以对发动机落点位置进行约束；在每次迭代中引入飞行时间下界约束来协调各导弹的飞行时间，以使所有导弹的飞行时间随迭代趋于一致；所有导弹的飞行时间达到一致后，DGPM-ANN算法完成。

由于在分布式算法的每一次迭代中，各导弹均独立求解最优轨迹，所以不同导弹的飞行时间可能不同。为了保证不同导弹飞行时间的一致性，迭代收敛时必须满足以下条件：

(30)

(31)

图8 DGPM-ANN算法过程

5 数值仿真

为了验证所提DGPM-ANN算法的高效性和有效性，通过两个仿真算例，利用一种集中式伪谱- 人工神经网络(CGPM-ANN)算法与本文提出的DGPM-ANN算法分别求解双脉冲导弹的多阶段协同轨迹规划问题，然后将两种算法的仿真结果进行比较。CGPM-ANN算法同时对所有导弹进行协同轨迹规划问题求解，并采用与CGPM-ANN算法中相同的神经网络拟合函数，以续航段末端状态为分离发动机壳体的初始状态，预测发动机落点。采用通用伪谱最优控制软件GPOPS实现算法，GPOPS软件使用默认的NLP求解程序IPOPT，仿真在一台使用英特尔i7-4790@3.6GHz处理器的计算机上运行。仿真参数和约束参数分别如表3和表4所示，表中：αmin和αmax分别为攻角的可行域上、下限；βmin和βmax分别为侧滑角的可行域上、下限；θ0min和θ0max分别为初始弹道倾角的可行域上、下限，各导弹根据任务分别确定其最优初始弹道倾角θ0. 建立仿真场景如下：4枚导弹同时从表5所示的不同地点发射，对坐标位置为(120 km,0 km,60 km)的地面目标进行协同攻击。DGPM-ANN算法与CGPM-ANN算法中各导弹的每个阶段离散点个数均为40，GPOPS软件采用默认参数。

表3 仿真参数表Tab.3 Simulation parameters

表4 约束参数Tab.4 Constraint parameters

表5 发射点位置Tab.5 Launch positions km

图9 算例1最优轨迹Fig.9 Optimal trajectory in Case 1

算例1以飞行时间最短为仿真目标，DGPM-ANN算法和CGPM-ANN算法得到的最优轨迹分别如图9(a)～图9(d)所示，相邻阶段的连接节点用黑点标记。导弹轨迹用实线表示；以求解得到的最优轨迹导弹续航段末端状态为分离发动机壳体初始状态，根据(5)式进行数值积分，得到分离的发动机壳体轨迹，由虚线表示；发动机落点安全区由黑色圆圈表示。由图9可知，DGPM-ANN算法与CGPM-ANN算法均实现了导弹对指定地面目标的协同攻击，且各分离发动机壳体均落到了安全区域。

DGPM-ANN算法得到的算例1最优飞行参数如图10所示。由图10可知，为了使飞行时间达到最短，发动机分离后，各导弹攻角迅速减小，使导弹获得负升力，同时在侧滑角的控制下，各导弹的飞行方向立即指向目标位置，制导攻击段导弹无明显侧向机动，以快速攻击地面目标。整个飞行过程中的过载满足约束。

算例2以控制量最省为仿真目标，DGPM-ANN算法和CGPM-ANN算法得到的最优轨迹分别如图11(a)～图11(d)所示，相邻阶段的连接节点用黑点标记。由图11可知，两种算法均实现了导弹对指定地面目标的协同攻击，且各分离发动机壳体均落到了安全区域。

算例2 DGPM-ANN算法得到的最优飞行参数如图12所示。由图12可知：为了使各导弹控制量最省，在飞行阶段的前期，导弹只在攻角和侧滑角作用下进行较小机动，以保证分离发动机壳体被抛落至安全区域；制导攻击段的后期，导弹进入空气密度较大的高度后，在攻角和侧滑角作用下对目标进行协同攻击。整个飞行过程中的过载满足约束。

两个算例中各导弹的优化参数，分别如表6和表7所示。由表6和表7可知，两个算例中DGPM-ANN算法的最优性能指标分别为CGPM-ANN算法的99.5%和92.4%，获得了更好的性能指标。这是因为采用CGPM-ANN算法求解该多阶段的最优控制问题时，所有导弹的动力学方程被联立为大规模方程组，所有方程共用同一时间变量，各导弹无法独立地确定其最适合协同飞行的最优二级推力点火时间，因此所有导弹的各阶段时间节点均相同。作为对比，DGPM-ANN算法求解时，各导弹可根据飞行任务和协同飞行任务的性能指标独立地优化其最优轨迹，因此各导弹拥有独立的时间变量，可以确定各自更适合协同飞行任务的最优二级推力点火时间，从而更好地实现各导弹发动机参数与弹道参数的共同优化，使协同飞行获得更好的性能指标。