近年来，多无人机(UAV)编队控制已经发展成为无人机控制领域的重要研究方向^[1-6]，无论是在军事应用还是民事应用上，其都吸引了大量的关注。与单机执行任务相比，多无人机编队执行任务在时间消耗、方案选择、任务完成率等方面具有显著优势。但是，当战场态势环境或者编队任务改变时，需要变换多无人机编队队形或者调整编队的松散程度，即编队重构，因此研究一种快速、安全和实用的多无人机编队重构控制方法非常重要^[7-8]。

当前，国内外对于多无人机编队重构进行了广泛的研究。针对战场威胁环境条件与机载不同载荷的无人机编队重构问题，文献[9]设计了一种基于纳什议价的分布式模型预测控制(DMPC)方法，能够使编队无人机协同躲避各种威胁以完成战术战略任务，但是整个飞行过程不考虑无人机编队的通信距离。在考虑各种约束条件(终端状态、安全防撞距离、通信距离等)的前提下，文献[10]采用控制作用参数化与时间离散化方法将寻找连续控制输入和终端时刻的无人机编队重构时间最优问题近似为寻找最优常数参数几何和分段时间问题进行求解，但是并没有涵盖障碍环境。通过使用编队中无人机与障碍的视线角以及无人机之间的相对速度，文献[11]提出了一种基于微分几何的控制律以实现无人机编队的避障与重构，但是这种方法很可能会产生过多的输入控制量从而使得结果不是很理想。文献[12]提出了一种用于四旋翼无人机编队重构的控制律，以使无人机编队稳定快速实现各种编队构型的重构切换，但是没有涉及避障问题。文献[13]将监控容错机制引入无人机编队控制律以实现重构，但是存在与文献[12]同样的问题。

在考虑各种约束条件以及队形切换的情况下，针对无人机数学模型、编队模型及碰撞环境威胁约束，本文首先通过引入邻居集这一概念并结合分布式协同策略构建多无人机编队的重构代价函数，进而提出了改进量子粒子群优化(RQPSO)算法对其进行求解，并与采用粒子群优化(PSO)算法求解的结果进行了对比。仿真结果表明，本文算法能安全、有效地完成无人机编队自主重构，实现隐蔽突防。

1 预备知识 1.1 模型预测控制

模型预测控制(MPC)主要包括3个方面，即预测模型、滚动优化与反馈校正。预测模型是指通过被控对象的过去信息和未来输入状态预测未来输出状态；滚动优化是指通过在线改变具有固定形式的优化性能指标的时间区域以实现优化；反馈校正用于修正模型失配。由于MPC采用的是一种有限时域优化，而不是全局优化策略，因此其抗干扰性能与鲁棒性能都非常强^[13-15]。

依据图 1所示，典型的MPC算法可以简单描述如下：

1) k时刻系统状态x(k)，使用系统过去的输入信息预测其未来的输出状态。

2) 选取某个性能指标进行在线优化计算，得到预测时域N内系统开环最优控制输入序列u^*(k)=[u^*(k∣k), u^*(k+1∣k), …, u^*(k+N-1∣k)]。

3) 取u^*(k)的第1项u^*(k∣k)作为k时刻的输入控制。

4) k+1时刻在新的时域上重复上述步骤。

1.2 量子粒子群优化算法

PSO算法自诞生以来就一直被用于求解各种复杂的优化问题^[16-17]，然而PSO算法本身具有搜索速度有限这一缺陷使得其不可能以概率1搜索全局最优解。科研人员经过研究发现，粒子在量子空间中具有聚集性的性质，通过将这一性质应用到PSO算法中就可以克服其存在的固有缺陷，这就是量子粒子群优化(QPSO)算法^[18]。在QPSO算法中，可以通过波函数ψ(X, t)来对粒子的运动状态进行描述^[19-20]：

式中：Q为概率密度函数，其满足如下条件：

在量子空间中，粒子的运动能够满足如下方程：

式中：ħ为普朗克常数；Ĥ为哈密顿算子，可以描述为如下形式：

其中：m为粒子质量；V(X)为粒子所在的势能；为拉普拉斯算子。该模型称为δ势阱模型。QPSO算法就是基于δ势阱模型被提出来的：假设在N维的目标搜索空间中存在M个粒子，种群中粒子组成的群体XL(t)=[XL₁(t), XL₂(t), …, XL_M(t)]为问题解, 在t时刻，XL_i(t)=[XL_{i, 1}(t), XL_{i, 2}(t), …, XL_{i, N}(t)](i=1, 2, …, M)表示第i个粒子的位置。整个种群中个体最好位置可以表示为：P_i(t)=[P_{i, 1}^l(t), P_{i, 2}^l(t), …, P_{2, N}^l(t)], 全局最好位置可以表示为: P^g=(P₁^g, P₂^g, …, P_N^g)。在QPSO算法中，每个粒子收敛于一点，即吸引子，p_i=(p_{i, 1}, p_{i, 2}, …, p_{i, N})。

式中：r_{1, j}和r_{2, j}为[0,1]之间随机数；c₁和c₂分别为自身因子和全局因子。

QPSO算法中构建了平均最好位置这一概念，其定义为种群中所有粒子个体最好位置的平均值，即

式中：M为种群规模。

QPSO算法的进化方程为

式中：u为[0,1]之间随机数；b为收缩-膨胀因子，一般，b取值为从1.0线性减小到0.5时，可以达到比较好的效果，即

其中：T为最大迭代次数；t_d为当前迭代次数。与PSO算法相比，在QPSO算法中，对于聚集在P^g附近的粒子，由于吸引子p_i点是在P_{i, j}^l和P_j^g之间，所以粒子以很大的概率在P^g附近出现^[20-25]。

2 无人机编队重构模型 2.1 无人机运动模型

假设无人机编队中共有N_v架无人机，则无人机编队集合可表示为V_c={uv_i|i=1, 2, …, N_v}，无人机自身能够实现速度与滚转角的相互解耦，能够始终保持飞行高度不变，则无人机近似离散化后的质心运动模型可以描述为^[26-27]

式中：x_i、y_i和v_i分别为无人机的坐标和速度；χ_i和γ_i分别为无人机的航迹方位角和航迹倾斜角；v_i^c和γ_i^c分别为速度指令和航迹倾斜角指令；α_v和α_γ分别为速度时间常数和航迹倾斜角时间常数；g为重力加速度；τ为采样周期。

取编队第i架无人机uv_i在k时刻状态变量为: X_i(k)=[x_i(k)，y_i(k)，v_i(k)，χ_i(k)，γ_i(k)]^T，第i架无人机uv_i在k时刻的控制变量为:U_i(k)=[v_i^c(k)，γ_i^c(k)]^T，这样可将无人机速度v_i运动方程简化为式(10)，并可得到约束条件如式(11) 所示。

通常存在一条已知参考轨迹以使无人机编队跟踪飞行，假设参考轨迹满足如下方程^[28]：

式中:x_r、y_r和χ_r分别为参考轨迹在地面坐标系的坐标和方位角；v_r和ω_r分别为参考轨迹的速度和角速度，v_r与ω_r满足分段连续且一直有界，满足如下约束方程：

2.2 编队模型

多无人机编队控制策略主要包括以下3种：领航-跟随法、虚拟领航法以及行为控制法^[29-30]。由于领航-跟随法具有误差传导等问题，而行为控制法具有数学模型不易建立等缺点，所以本文选用虚拟领航法。假设存在一无人机跟随轨迹的参考点，编队中的无人机能够提前或者通过无线通信实时获取参考点的轨迹，并且时刻保持与参考点的相对距离和角度，无人机之间不存在参考误差。轨迹参考点坐标系X_rO_rY_r固连于参考点O_r(见图 2)，可以得到编队中所有无人机的期望位置：

式中：x_i^dr和y_i^dr(k)分别为k时刻无人机uv_i期望位置与虚拟轨迹参考点位置之间的相对距离。

3 编队重构问题描述 3.1 威胁描述

战争中敌方的重要目标一般处于各种保护之中，为了实施战略战术纵深打击，多无人机编队需要隐蔽突防敌方防御区域，到达指定区域才能完成这一任务。在突防的过程中，多无人机编队会面临多种威胁，本文考虑预警雷达、防空雷达、禁飞区(包括固定障碍物、固定威胁以及移动威胁)3种威胁，假设编队中的无人机都能够探测到敌方预警雷达和防空雷达，且无人机携带干扰设备，能够对预警雷达与防空雷达实施一定程度的干扰。

1) 预警雷达威胁

如图 3所示，假设预警雷达位置为(x_ewr, y_ewr)，干扰设备位置为(x_je, y_je)，干扰设备的有效干扰距离为R_ewr，预警雷达最大探测半径为R_max，当无人机对预警雷达进行干扰时，其探测半径R可以描述为^{[9, 31]}

式中：

2) 防空雷达威胁

如图 4所示，假设防空雷达位置为(x_adr, y_adr)，干扰设备的位置为(x_je, y_je)，无人机uv_i位置为(x_i, y_i)，干扰设备的有效干扰距离为R_adr，防空雷达最大探测半径和最大探测角度分别为R_max和β，则无人机处于安全区域的条件为^{[9, 31]}

式中：;

3) 禁飞区威胁

如图 5所示，多无人机编队一旦遇到禁飞区或者障碍物只能够绕行避开，为了简化问题，本文统一将上述不规则图形或者矩形的最小外接圆作为禁飞区或者障碍物来讨论。假设禁飞区最小外接圆圆心位置和半径分别为(x_irc, y_irc)和R_irc，则无人机编队处于安全区域的条件为

3.2 编队无人机约束条件构造

一方面，由于无人机编队面临各种威胁(预警雷达威胁、防空雷达威胁、禁飞区或障碍威胁)，只有躲避这些威胁，无人机编队实现隐蔽突防飞抵指定区域才能实施对敌方目标的有效打击；另一方面，在无人机编队执行战略战术任务的过程中，需要避免编队内部的碰撞以及保持正常的实时通信。针对上述情况，可以构造如下约束条件：

1) 预警雷达威胁约束

假设k时刻编队内无人机探测到N_ewr个敌预警雷达ewr_j(j=1, 2, …, N_ewr)，预警雷达的位置为(x_ewrj, y_ewrj)，被干扰后的探测半径为R_j，根据图 3可知，编队内的无人机只有处于被干扰后预警雷达探测距离之外时才能保证其不被预警雷达发现，由此可以构造第i架无人机uv_i面临第j个预警雷达的威胁约束条件为

2) 防空雷达威胁约束

假设在k时刻编队内无人机探测到N_adr个敌防空雷达adr_j(j=1, 2, …, N_adr)，防空雷达的位置为(x_adrj, y_adrj)，被干扰后的安全距离和安全角余弦分别为d_sj和cs_sj，根据图 4可知，编队内的无人机只有处于安全区域内时才能保证其不被防空雷达发现，由此可以构造第i架无人机uv_i面临第j个防空雷达的威胁约束条件为

3) 禁飞区威胁约束

假设在k时刻编队内无人机探测到N_nfz个敌防空雷达nfz_j(j=1, 2, …, N_nfz)，禁飞区最小外接圆圆心的位置为(x_ircj, y_ircj)，根据图 5可知，编队内的无人机只有不飞抵禁飞区时才能保证其安全飞行，由此可以构造第i架无人机uv_i飞临第j个禁飞区的威胁约束条件为

4) 无人机防碰撞约束

为了防止编队内无人机发生碰撞，无人机之间的距离一定要大于防撞安全距离D_safe，由此可以构造k时刻uv_i与uv_j之间的防碰撞距离约束条件为

式中：。

5) 无人机通信距离约束

多无人机编队执行任务期间需要接收各种指令及信息，如上层发送的改变任务指令或者编队无人机各自的状态信息等，所以无人机之间必须能够进行正常实时的通信，由于通信设备的通信距离有限，所以需要保证编队内无人机之间的距离小于通信保障距离D_com，由此可以构造k时刻uv_i与uv_j之间的通信距离约束条件为

3.3 编队重构过程

在参照虚拟参考轨迹的情况下，编队中所有无人机在各自预定的相对位置飞行。当探测并定位出敌方预警雷达或者防空雷达威胁后，编队可以选择靠近威胁的某个无人机采用干扰，其余无人机进行队形调整，位于进行干扰无人机所产生的安全区域内，从而实现隐蔽突防。自动形成期望队形，根据威胁进行自动变换，以提高编队的安全性以及协同突防的有效性。

假设无人机编队可以对进行干扰任务的无人机自动分配，本文不进行研究。无人机编队能够根据环境威胁以及飞行过程中需要执行的任务进行构型自动选择。由于预警雷达和防空雷达系统2类威胁的作战能力和行为不同，因此针对这2类威胁应该选择不同的构型。

4 无人机编队重构控制问题的求解 4.1 邻居集定义

定义1 无人机uv_i的邻居集定义为与其具有耦合约束关系的所有无人机的集合。如用κ(i)表示无人机uv_i的邻居集，则无人机uv_i与无人机uv_j∈κ(i)存在耦合约束关系。另外，ς(i)表示与无人机uv_i存在紧约束的无人机集合。显然不是所有的耦合约束都是紧约束。

结合编队重构控制的特点以及MPC的预测特性，将邻居集内的无人机再进行分类，将在MPC预测时域内会违背最小间距约束的无人机视为紧约束邻居，在进行在线滚动优化时将紧约束邻居无人机的代价一起进行优化。这样做的好处是：不仅与邻居无人机进行了协同，而且避免了集中式优化问题中考虑邻居无人机内所有无人机的代价，降低了决策空间的维数，求解时间将缩短。如图 6所示，无人机uv_i的邻居集内拥有j₁、j₂、j₃、j₄和j₅ 5架无人机，但只有j₁和j₂属于紧约束邻居(图 6中标红色)，无人机uv_i在进行轨迹优化时，将考虑这2架无人机的代价，可以看出这一集合是时变的，这里假设每个时刻所有无人机能得到一致的紧约束邻居集合。

4.2 基于分布式模型预测控制的协同策略

在传统的集中式编队队形重构控制方法中，重构控制器需要在线优化所有无人机的控制输入，计算量非常大，不具有扩展性^[32]。

本文提出的分布式编队重构控制方法原理框图如图 7所示。无人机uv_i和无人机uv_j之间相互传递协同信息，并以某种性能指标最小为控制目标。常见的一种协同策略是重构完成时刻的协调，即相互协同确定重构完成时刻T_r，每个无人机的主要任务就是规划一条轨迹使得T_r时刻编队重构任务是可以保证完成的，并且还要保证自身的安全性。T_r是协调变量，每个无人机必须将其作为航迹生成的约束^[7-9]。

当耦合约束不是紧约束时，每个无人机“冻结”其他无人机的决策行为，然后求解自身的局部优化问题。当存在紧耦合约束时，每个无人机将考虑存在紧约束邻居无人机的局部代价以及编队总体代价，以在必要时需要牺牲自身性能来满足整体编队的性能指标。

当无人机uv_i在无人机uv_j之后求解优化问题时，对耦合约束进行修改，如果是紧约束，有

显然，当参数η>0时松弛了耦合约束，可能会使得局部性能变差，但对整个编队的代价有利。这样，无人机uv_i在求解自身优化问题时也考虑到了其他无人机的性能指标，即对整个编队的性能指标也有利。

利用无人机之间耦合约束的稀疏特性，决策空间将明显减小，计算量不会太大，这种协同优化策略的优点在于：① 存在紧耦合约束时不“冻结”其他无人机的行为；② 缩小了决策空间，局部优化问题求解的复杂性大大降低。

图 8所示为3种不同在线优化形成的航迹。图 8(a)所示为非协同优化航迹，无人机uv_i和无人机uv_j均优化自身的局部代价，则很容易违背最小间距约束；图 8(b)所示为基于优先级的协调航迹，无人机uv_j将无人机uv_i视为移动的障碍物进行规避，但由于进行了非协同的规避机动，无人机uv_j到达自身期望位置的时间较长，从而编队总体代价较大；图 8(c)所示为进行协同优化的航迹，无人机uv_i和无人机uv_j以编队总体代价为优化目标，无人机uv_i牺牲自身局部代价获得了编队总体代价最优。

编队重构分布式协同优化的步骤如下：

1) k时刻，无人机uv_i优化的目标函数中增加紧约束邻居无人机的局部代价以及编队总体代价。

式中：α>0；和可以体现出整个编队重构时间代价，即要求每个无人机不要偏离最终期望位置太远，尽快完成重构控制任务^[31]。

2) 根据预测时域对无人机最小间距约束进行单调性松弛(松弛因子为η)。

3) 优化求解得到U_i^*和U_j∈ς(i)^*，实施U_i^*的第1项U_i^*(k∣k)，并将预测状态向量X_i^*发送给无人机uv_j。

4) k+1时刻，再根据新的状态和新的邻居集j∈ς(i)重复上面过程，得到k+1时刻的控制输入。

综上，可得无人机编队的整体重构优化模型为

4.3 RQPSO-DMPC算法求解编队重构问题

2002年，Riget等^[32-33]将多样性引导策略引入PSO算法中，提出吸引-排斥PSO算法。多样性的大小能够改变搜索模式(即是进行全局搜索还是进行局部搜索)。其中吸引对应的是收敛阶段，而排斥则对应的是发散阶段。

本文选用粒子与中心点之间的平均欧几里得距离表示多样性，分别将粒子当前位置组成的种群和个体最好位置组成的种群表示为

这样可以得到QPSO算法的多样性度量为

式中：∣A∣为种群搜索空间中粒子最大值与最小值之间的一种关系表示；N为求解问题的空间维数。

时间对于战场来说至关重要，然而在QPSO算法中，种群多样性将随着搜索时间而下降，因此当进入搜索后期时，粒子多样性降低，全局搜索能力也将随之大大减弱，为了有效避免早熟现象，在文献[21]的基础上，设置多样性的下限d_low，本文将多样性引导的QPSO算法应用于多无人机编队重构控制。

多样性控制的量子粒子群优化(DCQPSO)算法中，diversity(S_XL)被用于引导算法搜索方向(即收敛方向还是发散方向)。一旦预先设定的下限d_low大于diversity(S_XL)，整个粒子群将由收敛状态进入到发散状态，此时多样性将会出现暂时性增加的情况，直到多样性diversity(S_XL)再次大于d_low。

当diversity(S_XL)或diversity(S_P)降低到d_low以下时，进行如下变异操作：

式中：P_j^g为处于最好位置的粒子；ε为服从标准正态分布的随机数；为变异参数。

通过施加变异操作，可以使得P_i(t)与M_best之间距离增加，从而增大diversity(S_P)。与此同时，由于P^g(t)的变动导致M_best改变，这样粒子当前位置与M_best之间距离也会增大，势必引起粒子产生一定程度的发散，从而增大diversity(S_XL)。因此可以通过增加2种多样性(diversity(S_P)或者diversity(S_XL))而进行变异，前一种方法称为DCQPSO_P，后一种方法称为DCQPSO_X。

除了对粒子群进行变异操作以增加种群多样性，本文还通过在群体全局最好位置P^g上施加选择操作以丰富种群的多样性，从而避免种群过早陷入局部最优解。即p_i不再通过P^g(t)和P_i^l决定，而是由P_i^l和某个随机选择出的其他粒子的个体最好位置(其目标函数值要比P_i^l好)决定。在从粒子群中随机选择某个其他粒子P_k^l(t)后，可以进行如下操作：

当选择出来的P_k^l(t)的目标函数值好于P_i^l(t)，则p_i由P_k^l(t)和P_i^l决定，否则P_i(t)仍由P^g(t)和P_i^l(t)决定。因此，通过上述选择操作选出P^g后，P_i(t)的坐标为

引入多样性与选择策略可以得到改进的QPSO算法，即RQPSO算法，其流程描述如下^[34-35]：

步骤1 对粒子进行编码(见图 9)，设置算法基本参数，随机初始化粒子群。

构建无人机编队重构问题与RQPSO算法粒子之间的映射关系，即粒子编码。本文采用的非线性模型中控制输入指令为速度和航向，每个时刻无人机uv_i优化的变量维数并不一样，因为邻居集是时变的，所以粒子群的维数也是时变的。如图 9所示，如果无人机uv_i只优化自身局部代价，粒子群的维数为M×2N，即M个粒子，每个粒子包含2N个变量。考虑紧耦合约束无人机后，无人机uv_i在线优化的变量个数为2(n+1)N(n为ς(i)内元素的个数)，此时粒子群的维数为M×2(n+1)N。

步骤2 当迭代次数小于T时，执行以下步骤。

步骤3 计算粒子群的平均最好位置M_best。

步骤4 计算粒子的当前适应值，并与前一次迭代的适应值比较，如果f(XL_i(t+1)) < f(P_i(t))，则P_i(t+1)=XL_i(t+1)。

步骤5 计算群体当前的全局最优位置P^g。

步骤6 将当前时刻全局最优位置与前一次迭代的全局最优位置进行比较，选择更优值作为群体的全局最优位置。

步骤7 对每一个粒子，选择点P^g。对粒子的每一维，根据式(35) 计算p_{i, j}。

步骤8 根据式(7) 计算粒子的新的位置。

步骤9 计算粒子群的多样性diversity(S_XL)(DCQPSO_X)或diversity(S_P)(DCQPSO_P)。

步骤10 通过式(8) 计算收缩-扩张因子(收敛模式)。

步骤11 判断diversity(S_XL)或diversity(S_P)，如果小于d_low，则进行变异操作(发散模式)。

步骤12 按照QPSO算法进化方程对种群中所有粒子进行位置更新。

步骤13 返回步骤2。

5 仿真与分析

为了验证本文所提算法的有效性，使用MATLAB编译环境进行仿真验证，分别选用PSO、RQPSO 2种算法对无人机编队重构控制进行仿真实验，粒子群规模都设置为20。对于PSO算法，c₁=c₂=2，w=0.5，T=800；对于RQPSO算法，b的值将随着算法的运行从1.0线性减小到0.5，d_low=0.000 4， =0.000 012。

假设编队中的无人机之间不存在通信时延；无人机的速度、航向不受任何其他因素干扰影响。各无人机的初始运动以及参考轨迹参数如表 1和表 2所示，各无人机与参考轨迹相对距离如表 3所示，预警雷达、防空雷达、禁飞区坐标及作用距离如表 4所示。

仿真结果如图 10~图 15所示。

从图 10可以看出，无人机编队在整个飞行过程中基本没有偏离虚拟参考轨迹，当发现敌方搜索雷达，需要进行突防时，无人机编队中的uv₅被选择对敌搜索雷达进行干扰，其余无人机能够以uv₅与虚拟参考轨迹为参考自动进行编队构型重构，以实施安全突防；当发现敌方防空导弹跟踪雷达时，无人机编队中的uv₃被选择对敌跟踪雷达进行干扰，此时编队中的uv₃与uv₅相互协作，同时还要参考虚拟参考轨迹以便选择出最佳干扰位置，其余无人机也自动进行编队构型重构，以实现安全突防；当整个编队突防成功后，需要根据此时的位置和期望的编队构型在此进行编队重构。从图 11和图 12可以看出，uv₁的速度和航向控制输入满足控制约束。从图 14可以看出，无人机之间的相对距离也满足约束条件，从而验证了所设计的编队构型控制算法是有效的。图 15表明，RQPSO算法进行编队重构的效果要明显优于PSO算法。

6 结论

本文提出了一种基于RQPSO算法的多无人机编队重构控制算法，并验证了其有效性与安全性。

1) 相比传统的PSO算法，RQPSO算法的变异操作和选择操作使得整个算法的搜索能力更强大，且不易于陷入局部最优值，而使得到的解不是最优的。

2) 通过引入邻居集降低编队重构代价问题的复杂性，提高了优化求解的效率，更有利于对实时性要求高的多无人机编队自主重构。

3) 与以往的各种多无人机编队重构约束条件相比，本文考虑的约束条件更加全面具体，具有更多的现实意义。