中国科技核心期刊      中国指挥与控制学会会刊     军事装备类重点期刊

An algorithm for joint air defense problem based on improved r-NSGA-Ⅱ algorithm

  • JIANG Jia-rui ,
  • XU Guo-liang ,
  • WANG Tuan-tuan
Expand
  • Jiangsu Automation Research Institute, Lianyungang 222061, China

Received date: 2022-10-13

  Revised date: 2022-11-04

  Online published: 2023-02-20

Abstract

Based on the air defense and antimissile operations of traditional kinetic energy weapons and modern high-energy laser weapons, a dynamic weapon target assignment model of heterogeneous defense weapons with multiple tactical targets is established; Considering the commander's choice preference, the optimization framework of weapon target allocation based on improved r-NSGA-Ⅱ algorithm is proposed, and the specific optimization process is given; Based on the specific operational scenario, the design experiment is designed to compare the improved r-NSGA-Ⅱ algorithm with NSGA-Ⅱ, NSGA-Ⅲ, MOEA/D and other classical algorithms in terms of convergence, distribution, time complexity and other performance indicators. The experimental results show that the improved r-NSGA-Ⅱ algorithm can better reflect the commander's operational preferences and intentions while ensuring the computational efficiency.

Cite this article

JIANG Jia-rui , XU Guo-liang , WANG Tuan-tuan . An algorithm for joint air defense problem based on improved r-NSGA-Ⅱ algorithm[J]. Command Control and Simulation, 2023 , 45(1) : 75 -83 . DOI: 10.3969/j.issn.1673-3819.2023.01.013

随着现代信息技术的飞速发展,以导弹为代表的各类武器正在朝着更快、更远、更精确的方向快速进步,如何快速有效地组织我方火力拦截来袭目标这一典型武器目标分配问题(Weapon-Target Assignment, WTA)就显得愈发重要。
WTA问题研究的发展大致遵循从静态到动态、从单目标优化到多目标优化的趋势[1]。Manne[2]于1958年提出了一种以来袭目标存活概率最小化为目标的静态WTA模型。静态WTA问题在建模时未考虑时间因素所产生的影响,优化目标以最小化目标存活概率为主[3,4],模型的区别体现在模型假设上,例如武器同构性与目标同构性等方面。动态WTA模型加入了时间因素,研究着重关注引入时间因素带来的动态性与随机性。Xin等[5]建立了多阶段打击的保护多资源点的动态WTA模型;Murphey[6]使用马尔科夫随机过程描述后续时刻可能新出现的目标数量,并在此基础上建立了动态WTA模型,体现了对引入时间因素后所产生的随机性的研究。
在优化目标数量方面,单目标WTA模型着眼于最重要的战术目标,即最小化来袭目标生存概率,而多目标WTA模型可以同时兼顾防御方打击成本与目标打击时间等因素。由于单目标WTA模型存在资源浪费、不符合实际的缺陷,目前主流WTA模型一般为多目标WTA模型。针对实际作战需求,邵诗佳[7]面向转火时间、火力通道等作战约束,提出了基于最大化毁伤概率以及最小化拦截时间的多阶段动态WTA模型;邱少明等[8]建立了以最小化目标生存概率和最佳费效比为目标函数建立模型,将智能搜索算法引入进化算法以缓解算法进化过程中出现的多样性缺失的问题,经仿真实验证明其存在性能的优化;陈曼[9]建立了以最小化目标生存概率与最小化防御方武器成本为目标函数建立模型,提出粒子群算法改进方案以完成模型的优化。上述研究成果在一定程度上深化了对多目标WTA问题的研究水平,但仍存在两点不足:1)模型并未体现出实际作战中指挥作战人员对不同优化目标的偏好区分以及对作战效果的预期,例如来袭目标的存活概率应该在算法模型的众多优化目标中居于首位,而其他优化目标的重要性应该相对较低;2)目前模型中尚未充分考虑以高能激光武器为代表的新质武器的加入对WTA问题产生的影响。
本文构建的WTA模型为高维多目标优化模型。多目标优化问题所能得到的最优解被称为非支配解(又称Pareto解),对于一个多目标优化问题,其Pareto解往往不只一个。在解空间中,由Pareto解组成的曲面被称为Pareto前沿面。多目标进化算法是求解该类问题的有效方法,尤其是非支配排序进化算法(Non-dominated Sorting Genetic Algorithm,NSGA)及其改进算法。Deb等[10]为解决多目标优化难以区分解的优劣问题,于1995年提出了NSGA算法,随后于2002年对NSGA进行改进,提出NSGA-Ⅱ算法,相较于NSGA,NSGA-Ⅱ算法能够高效地求解多目标优化问题,且生成的解可均匀分布在Pareto近似前沿面上。Said等[11]在2010年基于NSGA-Ⅱ算法提出了基于r-支配的非支配进化算法(Non-r-dominated Sorting Genetic Algorithm Ⅱ, r-NSGA-Ⅱ),该算法将指挥员的抉择偏好体现为参考点与权重向量,并通过这两个参数修改非支配排序进化算法中的支配方式,从而使得算法在一定程度上能够引导种群的进化方向,使种群整体向着权重向量与参考点决定的方向进化,最终收敛在Pareto前沿面与偏好相关部分的附近。作为NSGA-Ⅱ算法的改进算法,r-NSGA-Ⅱ算法在运行效率、收敛性等方面表现也非常优秀。
基于以上论述,本研究将高能激光武器作为防御武器引入WTA模型,在充分考虑打击时间区段、防御方武器资源以及火力冲突约束等信息的基础上,以来袭目标存活概率最小,武器资源消耗最少为优化目标,以0-1决策组建决策向量来构建WTA模型;将指挥员抉择偏好作为影响因子在优化算法中体现,基于改进r-NSGA-Ⅱ算法给出了包含偏好的多目标WTA模型的一种解决方案;最后通过仿真实验对算法的有效性进行了验证。

1 多目标WTA模型

假设m个武器资源参与作战拦截n批空中来袭目标,其中,动能武器资源集合W1={wi,i=1,…,m1},高能激光武器资源集合W2={wi,i=m1+1,…,m1+m2},m1+m2=m;武器弹药成本集合C={ci,1,…,m1+m2};目标集合T={ti,i=1,…,n};目标价值集合V={vj,j=1,…,n}。模型假设如下:
1)所有来袭目标的位置、速度、数量等信息已知;
2)高能激光武器与动能武器的毁伤概率不会相互影响,即不考虑激光武器照射目标对动能武器打击效果的影响以及动能武器对目标完成打击后对激光武器打击效果的影响;
3)目标飞过航路捷径垂足点后不再进行打击。航路捷径如图1所示,其中,O点为我方武器位置,P点为目标位置,向量Vm为目标飞行的方向,航路捷径为s,垂足点为T。可以看到,目标飞临垂足点T之前,动能武器弹丸与目标相对飞行,而目标飞过垂足点后弹丸将追逐目标,其毁伤能力将大幅降低,此外,当目标飞离我方时,其尾焰将极大影响激光武器的打击,由此做出该假设。
图1 我舰与目标态势关系水平示意图

Fig.1 Schematic diagram of our ship and target situation relationship

1.1 可行打击时间窗

本文将发现目标的时刻作为WTA决策方案时间轴的原点,将所有来袭目标均飞离防御武器射击范围的时间作为时间轴的终点,以来袭目标相对防御武器的位置与速度、防御武器弹药飞行速度以及防御武器火力范围为依据,可得每型武器wi对每批来袭目标tj的可行打击时间窗区间,有如下三种情况:
Case1:
TimeStartij=TimeEndij=0; if ShortCutij≥DisUpperij
Case2:
T i m e S t a r t i j = m a x 0 , D i s i j 2 - S h o r t C u t i j 2   r a t e j - D i s U p p e r i j 2 - S h o r t C u t i j 2   r a t e j - D i s U p p e r i j v b u l l e t i T i m e E n d i j = m a x 0 , D i s i j 2 - S h o r t C u t i j 2   r a t e j - S h o r t C u t i j v b u l l e t i i f   D i s U p p e r i j S h o r t C u t i j D i s L o w e r i j
Case3:
T i m e S t a r t i j = m a x 0 , D i s i j 2 - S h o r t C u t i j 2 r a t e j - D i s U p p e r i j 2 - S h o r t C u t i j 2   r a t e j - D i s U p p e r i j v b u l l e t i T i m e E n d i j = m a x 0 , D i s i j 2 - S h o r t C u t i j 2   r a t e j - D i s L o w e r i j 2 - S h o r t C u t i j 2   r a t e j - D i s L o w e r i j v b u l l e t i i f   D i s L o w e r i j S h o r t C u t i j
式中:
TimeStartij:武器i打击目标j的最早时间;
TimeEndij:武器i打击目标j的最晚时间;
ShortCutij:目标的航路捷径;
DisLowerij:武器i对目标j的射击近界;
DisUpperij:武器i对目标j的射击远界;
Disij:当前时刻武器i与目标j的直线距离;
ratej:目标j的飞行速度;
vbulleti:武器弹丸飞行速度(高能激光武器对应为无穷)。
在此基础上,确定关于每批来袭目标j对于我方每个武器单元的可行打击时间{[TimeStartij, TimeEndij]}i=1,2,..m j=1,2,…,n,得到所有来袭目标的最早开火时间 m i n i , jTimeStartij,最晚开火时间 m a x i , jTimeEndij

1.2 决策向量

将问题涉及的时间轴[ m i n i , jTimeStartij, m a x i , jTimeEndij]平均划分为T个时隙,每一个时隙可以表示为:
$slo{{t}_{t}}=[\underset{i,j}{\mathop{\text{min}}}\,TimeStar{{t}_{ij}}+(t-1)len,\underset{i,j}{\mathop{\text{min}}}\,TimeStar{{t}_{ij}}+t\cdot len],t=1,\ldots,T$
len=( m a x i , jTimeEndij- m i n i , jTimeStartij)/T
在此基础上,模型假设武器在同一个时隙内至多只能做一次打击决策。如果时隙slott在武器i对于来袭目标j的可行打击时间窗[TimeStartij,TimeEndij]内,则可以做出使用武器i攻击目标j的决策。由于每一个武器-目标对均存在各自的可行打击时间窗,对每个武器-目标对可以做出攻击决策的时隙编号可以记录如下:
Iij={t|t∈{1,…,T},slott⊂[TimeStartij,TimeEndij]}
由于在每一个时隙,决策是否攻击是一个二值变量,模型的决策向量为
X={ x i j t, i=1,…,m j=1,…,n t∈Iij x i j t ∈{0,1}}

1.3 武器毁伤特性分析

1.3.1 武器特性分析

本研究所假定的防御武器主要包括动能武器及高能激光武器,动能武器又包含制导弹药武器以及非制导类武器。依据武器特性可知,非制导类动能武器对目标的毁伤能力主要受武器系统误差以及弹丸散布误差等因素影响,其毁伤概率可依据现有国军标进行计算;制导类动能武器主要受其制导能力影响,其对目标的毁伤能力可近似认为仅与目标类型有关,而与打击时间、射程范围以内的距离无关;高能激光武器主要通过将激光束持续、稳定地投射在目标表面某具体位置上,在能量累积达到目标表面材质毁伤阈值后,目标结构发生损坏,从而达到毁伤目的,其毁伤概率计算特点如下:
1)高能激光光束不与来袭目标靶面接触时,基本认为无法毁伤;
2)命中后,毁伤情况与激光束照射时间长短直接相关;
3)高精度跟瞄系统和光传播的瞬时性,使得高能激光武器对目标同一部位实现多次打击成为可能,命中后毁伤概率具有累计增加的可能性。

1.3.2 基于马尔科夫决策过程的高能激光武器毁伤概率计算

本研究中构建了动态WTA模型,对于作战过程时间采用离散化方式处理,由于本文在构建多目标武器目标分配时,对于模型的动态性采用时间离散化的方式,高能激光武器的打击决策过程可通过离散马尔科夫决策过程进行描述。马尔科夫决策过程假设系统下一时刻的状态仅取决于当前系统状态以及当前所选择的行为,与系统以往过程无关。马尔科夫决策过程包含5个模型要素:状态、动作、策略、奖励和回报[12],针对高能激光武器对目标的打击过程,状态集可以描述为S={si,i=1,…6},动作集可以描述为U={u1,u2},元素含义表示如下:
s1:打击过程开始;
s2:打击过程结束,目标存活;
s3:能量累积完成,目标摧毁;
s4:打击未命中或能量累积未完成;
s5:累积能量仍有效;
s6:累积能量无效;
u1:开始打击或继续打击目标;
u2:停止打击目标。
高能激光武器是通过向目标投射足量的能量完成对目标打击的武器。由于打击持续时长较短,模型假设激光武器的打击过程可近似为绝热过程,能量投射到目标上后不会逸散,之后的照射打击在目标原本的照射点附近时,能量累积的进程有可能被继承,基于绝热过程的假设,能量累积的概率与时间无关,激光武器对目标的连续打击与间断打击是等价的。因此,在描述激光武器打击来袭目标的马尔科夫决策过程中,策略仅设定为m次连续的打击指令,而本节目标是计算m次连续打击后目标的毁伤概率,不另外设定奖励函数与回报函数。高能激光武器打击目标马尔科夫决策过程的图模型如图2
图2 高能激光武器打击过程图模型

Fig.2 High energy laser weapon strike process

其中:
p1:高能激光首次打击命中并毁伤目标的概率;
p2:后续打击到达时,累积能量仍有效的概率;
p3:累积能量有效的情况下,高能激光摧毁目标的概率;
p4:累积能量无效的情况下,高能激光摧毁目标的概率;
由上述马尔科夫决策过程可得到高能激光武器连续m次打击目标的毁伤概率为
prij(m)=   p 1 ,                       i f m = 1 p 1 + ( 1 - p 1 ) β 1 - α m - 1 1 - α , i f m 2
α=p2(1-p3)+(1-p2)(1-p4)β=p2p3+(1-p2)p4
在上述推导中,累积能量有效的概率p2与能量留存情况下毁伤目标的概率p3对于总体毁伤概率存在重要影响,其值与前次照射后留存的能量有较大关系。考虑实际作战情况下,随着打击次数m的提升,p2p3的值均会增大,可以将p2p3的取值设定为线性递增,也可以根据实际量测情况求取均值。

1.4 目标函数

本研究构建的多WTA模型为多目标优化模型,其优化目标:来袭目标突防概率最小以及武器资源消耗最少。
以来袭目标突防概率最小为目标建立的函数如下:
f1= j = 1 nvj i = 1 m 1 t = 1 T(1-pij ) x i j t× i = m 1 + 1 mprij t = 1 T x i j t  / j = 1 nvj
以武器资源消耗最少为目标建立函数如下:
f2= j = 1 n i = 1 m t = 1 Tci x i j t

1.5 约束条件

集群目标火力分配模型的约束条件主要体现在火力资源的约束与火力冲突的约束。其中,火力资源约束包括武器弹药消耗方面的约束与火力资源在时间上的约束;火力冲突约束包括动能武器之间的火力冲突约束与动能武器与高能激光武器之间的火力冲突约束。综合考虑模型的性质并将各约束条件离散化后,可表示如下:
1)在武器弹药消耗方面,动能武器存在弹药总量的火力资源约束,假设每型武器弹余量为ri;高能激光武器的消耗主要体现在电能上,若有电能使用限制,可以从电能的限制与高能激光武器每次打击的平均能量消耗计算得到该激光武器的剩余打击次数ri。弹药消耗方面的火力资源约束可以表示如下:
j = 1 n t = 1 T x i j t≤ri,∀i
2)在时间方面,基于对时隙的假设,武器i在同一时隙slott内至多只能进行一次打击决策。
j = 1 n x i j t≤1,∀i,t
3)武器存在连续射击时间限制。另外,各武器也存在转火时间限制。假设连续时间约束为tlim1i,转火时间约束为tlim2i,其表示如下:
x i j t 1 x i j t 2 ≠1, if |t 1-t2|≤tlim1i,∀t1≠t2
x i j 1 t 1 x i j 2 t 2 ≠1, if |t 1-t2|≤tlim2i,∀t1≠t2,j1≠j2
4)动能武器之间存在火力冲突现象。如果两型武器在较短的时间间隔内发射制导弹药攻击同一个目标,可能会导致制导弹药在打击过程中相互干扰甚至追踪错误目标;如果两型武器在较短的时间间隔内分别发射非制导弹药和制导弹药,后者可能会被前者摧毁。因此建立约束,使得不同武器不能在时间约束tlim3内攻击同一个目标。
x i 1 j t 1 x i 2 j t 2 ≠1, if |t 1-t2|≤tlim3,∀i1≠i2
5)激光武器与动能武器存在时间与空间耦合的火力冲突现象。由于动能武器弹丸飞行过程中产生的尾焰会对高能激光武器命中产生极大影响,与该动能武器弹丸飞行与爆炸影响的区域存在时空重叠的目标不能被高能激光武器选为打击的目标;此外,由于动能武器在发射时产生的粉尘会形成烟幕干扰高能激光武器的瞄准,同时粉尘会折射高能激光,使得平台光辐射浓度超标,伤害平台上的战斗人员和设备,当武器平台上有动能武器发射时,高能激光武器的打击将受到限制。假设动能武器wi飞行至打击远界的时间为tUpperi,上述约束条件表示如下:
x i j t=0,i∈[m1+1,m], if j = 1 n i = 1 m 1 x i j t≥1
x i 1 j t 1=0,i1∈[m1+1,m], if x i 2 j t 2=1,i2∈[1,m1], 0≤t1-t2<tUpperi

1.6 数学模型

综上所述,本研究构建的多目标WTA模型为
min(f1,f2), s.t.   j = 1 n t = 1 T x i j t r i j = 1 n x i j t 1 , i , t x i j t 1   x i j t 2 1 , i f | t 1 - t 2 | t l i m 1 i , t 1 t 2 x i j 1 t 1   x i j 2 t 2 1 , i f | t 1 - t 2 | t l i m 2 i , t 1 t 2 , j 1 j 2 x i 1 j t 1   x i 2 j t 2 1 , i f | t 1 - t 2 | t l i m 3 , i 1 i 2 x i j t = 0 , i [ m 1 + 1 , m ] ,   i f j = 1 n i = 1 m 1 x i j t 1 x i 1 j t 1 = 0 , i 1 [ m 1 + 1 , m ] , i f x i 2 j t 2 = 1 , i 2 [ 1 , m 1 ] , 0 t 1 - t 2 < t U p p e r i
该模型为多目标优化的离散非线性模型,解模型的关键在于Pareto前沿面的求解。指挥员抉择偏好将使得Pareto前沿面上不同位置的解处于不同的地位,最终导致最优解区域集中在前沿面的某些部位。指挥员抉择偏好可以体现为权重向量与参考点,前者用于表示不同优化目标的重要性,后者用于表示指挥员对于不同优化目标的解的预期,采用的寻优算法需要能体现不同的偏好。

2 算法简介

基于r-支配的非支配排序进化算法(r-NSGA-Ⅱ)可以很好地解决带偏好的多目标优化问题。算法将DM的偏好体现为参考点与权重向量并以此改进了NSGA-Ⅱ算法中的支配关系,使算法能够在一定程度上引导种群进化的方向。另外,针对多目标WTA问题,本文对该r-NSGA-Ⅱ算法进行了改进,设计了一种改进r-NSGA-Ⅱ算法。

2.1 NSGA-Ⅱ算法简介

NSGA-Ⅱ算法的主要步骤包括种群初始化,子代种群生成,父子种群合并,非支配排序,拥挤距离计算,新一代种群选择这六个步骤。算法的运行流程如图3所示。
图3 NSGA-Ⅱ算法流程

Fig.3 Algorithm flow of NSGA-Ⅱ

2.1.1 种群初始化

种群初始化分为编码和生成初始化种群两个部分,生成的初始种群记为A0:
1)种群初始化采用0-1编码方式,种群个体为武器资源打击目标的打击方案。
2)初始化种群需要兼顾多样性,均匀性和可行性,生成大量可行解并均匀选取进入初代种群。

2.1.2 非支配排序

构建新的种群St+1,将种群Rt+1划分为若干非支配层(F1,F2,…),可以认为,非支配层数i较小的个体的优先度更高。Deb给出了一种快速进行非支配排序的算法,该算法的复杂度为O(nm)。

2.1.3 拥挤距离计算

拥挤距离计算用于在同一个非支配层Fl中比较个体的优先度。出于基因的多样性考虑,种群的进化过程中倾向于选择具有独特基因的个体,因此拥挤距离较大的个体具备更高的优先级。个体i的拥挤度:
$\begin{matrix} & {{I}_{crowd}}=\overset{m}{\mathop{\underset{k=1}{\mathop \sum }\,}}\,(I[i+1].k-I[i-1].k)/(f_{k}^{\text{max}}-f_{k}^{\text{min}}) \\ & I[i+1].k=\underset{f{{\left[ j \right]}_{k}}\ge f{{\left[ i \right]}_{k}},i\ne j,i,j\in {{F}_{l}}}{\mathop{\text{min}}}\,f{{[j]}_{k}} \\ & I[i-1].k=\underset{f{{\left[ j \right]}_{k}}\le f{{\left[ i \right]}_{k}},i\ne j,i,j\in {{F}_{l}}}{\mathop{\text{max}}}\,f{{[j]}_{k}} \\ \end{matrix}$
其中,m为优化目标数量,f[j]k,jFl表示非支配层Fl中的个体j在第k个目标函数上的取值。

2.1.4 子代生成

种群采用锦标赛选择机制产生基因池,采用单点交叉和单点随机变异机制生成子代种群,第t次进化生成的子代种群记为Bt
锦标赛选择机制反复从基因池中随机选取两个个体P1P2,通过依次比较两个个体的非支配层级与拥挤距离,选择优先度更高的个体加入基因池,反复选择直到基因池中个体数量为n
单点交叉机制从基因池中随机选取两个个体P1P2,按交叉概率p1进行交叉操作,随机选取一个单点,将两个个体在断点出断开,重新组合成为一个新的个体。
单点随机变异机制对完成交叉操作后的个体进行操作,按照变异概率p2进行变异操作,随机选取一个单点,并对对应位置的值取反。

2.1.5 父代种群与子代种群合并

将父代种群At与生成的子代种群Bt+1合并为种群Rt+1

2.2 r-NSGA-Ⅱ算法简介

r-NSGA-Ⅱ算法在结构上与NSGA-Ⅱ算法保持一致,因此能继承NSGA-Ⅱ算法优秀的效率与准确性。r-NSGA-Ⅱ算法将NSGA-Ⅱ算法中用于衡量两个个体之间的支配关系的Pareto支配重构为r-支配,其余部分保持不变。在下文中将个体x在Pareto意义上支配个体y记为xy,将个体x r-支配个体y记为xryxry指下面两种情况之一发生:
1)xy
2)xy在Pareto意义上等价且D(x,y,g)<-δ
其中:
Dist(x,g)= i = 1 m w i f i ( x ) - f i ( g ) f i m a x - f i m i n   2
D(x,y,g)= D i s t ( x , g ) - D i s t ( y , g ) D i s t m a x - D i s t m i n
wi∈(0,1), i = 1 mwi=1,x∈R
Distmax= m a x z RDist(z,g)
Distmin= m i n z RDist(z,g)
R表示种群整体,m表示优化目标的数量,fi(x)表示个体x在等i个目标函数上的取值, f i m a x 表示种群R中所有个体在第i个目标函数上的最大取值, f i m i n 表示种群R中所有个体在第i个目标函数上的最小取值。g表示参考点,wi表示不同目标函数的权重,δ∈[01],表示r-非支配限界。
通过重构个体间的支配关系,算法可以通过设置参考点g和权重w来引导种群进化的方向,通过设置r-非支配限界δ来改变种群收敛速度。可以证明,当δ取值为1时,r-支配退化为Pareto支配。

2.3 针对多目标WTA问题的改进

由于本文建立的模型是多目标约束优化模型,优化算法需要处理模型的约束问题。Deb[13]提出的NSGA-Ⅱ算法使用的是可行性法则这一最为经典的约束处理机制。作为NSGA-Ⅱ算法的变种算法,r-NSGA-Ⅱ算法也可以选用可行性法则这一机制来处理模型中的约束条件。可行性法则的基本思想如下:
当解x与解y均为可行解或均为不可行解但违约度相同时,按照两个解的目标函数比较占优关系;
当解xy均为不可行解且违约度不同时,按照两个解的违约度比较占优关系;
可行解优于不可行解。
可行性准则原理简单,方便执行,但完全忽略了不可行解的信息,当存在多个不连续的可行域时,算法难以穿过不可行域。由于多目标WTA模型的约束条件较多,可行域的间断是不可避免的。
为了充分利用不可行解的信息,引导进化中的种群穿过不可行域,实现全局搜索,本文将ε约束法[14]引入r-NSGA-Ⅱ算法。将解xε约束法意义下优于解y记为xεyxεy指下面三种情况之一发生:
1)xryϕ(x),ϕ(y)<ε;
2)xryϕ(x)=ϕ(y);
3)ϕ(x)<ϕ(y)且ϕ(x),ϕ(y)<ε不成立。
ε约束法的关键问题是ε值的设定。当ε=0时,ε约束法退化为可行性约束,当ε=∞时,ε约束法退化为r-支配,将完全无视约束条件。按照进化算法的常规流程,种群会在进化的过程中逐渐优化,不可行解的比例逐渐减少,不可行解中的有效信息自然也会逐渐变少。因此ε值应该随着进化过程推荐逐渐变小。传统的ε设置方式如下[15]:
ε(0)=φ(xθ)
ε(G)= ε ( 0 ) ( 1 - G / T c ) c p . 0 < G < T c 0 . G T c
其中,G是进化代数,xθ是初始种群个体按照违约度ϕ(x)降序排列后的第θ个个体,cp是控制ε下降速度的参数。经过试验证明,θ的值通常为0.05*N,N为种群个体数量,cp∈[2,10],Tc通常设置在区间[0.1*Gmax,0.8*Gmax],Gmax为最大进化代数[16]
对于多目标WTA模型,个体x的违约度ϕ(x)设置为个体x约束函数归一化后的加权求和:
φ(x)= k = 1 7akg'k (x)
g'k (x)= 0 , i f   g k ( x ) = 0 , x R g k ( x ) / m a x x R g k ( x )
其中,R表示当前种群,A=(ak),k=1,…,7为权重向量,代表各约束条件的重要程度,gk(x)为个体x的约束函数。按照多目标WTA模型的约束条件,约束函数定义如下:
$\begin{matrix} & {{g}_{1}}(x)=\overset{m}{\mathop{\underset{i=1}{\mathop \sum }\,}}\,max\left( 0,\overset{n}{\mathop{\underset{j=1}{\mathop \sum }\,}}\,\overset{T}{\mathop{\underset{t=1}{\mathop \sum }\,}}\,x_{ij}^{t}-{{r}_{i}} \right) \\ & {{g}_{2}}(x)=\overset{m}{\mathop{\underset{i=1}{\mathop \sum }\,}}\,\overset{T}{\mathop{\underset{t=1}{\mathop \sum }\,}}\,max\left( 0,\overset{n}{\mathop{\underset{j=1}{\mathop \sum }\,}}\,x_{ij}^{t}-1 \right) \\ & {{g}_{3}}(x)=\overset{m}{\mathop{\underset{i=1}{\mathop \sum }\,}}\,\overset{n}{\mathop{\underset{j=1}{\mathop \sum }\,}}\,\underset{{{t}_{1}}\ne {{t}_{2}},x_{ij}^{{{t}_{1}}}\text{ }\!\!~\!\!\text{ }x_{ij}^{{{t}_{2}}}=1}{\mathop \sum }\,max(0,tlim{{1}_{i}}-\left| {{t}_{1}}-{{t}_{2}} \right|) \\ & {{g}_{4}}(x)=\overset{m}{\mathop{\underset{i=1}{\mathop \sum }\,}}\,\underset{\forall {{t}_{1}},{{t}_{2}},{{j}_{1}}\ne {{j}_{2}},x_{ij}^{{{t}_{1}}}\text{ }\!\!~\!\!\text{ }x_{ij}^{{{t}_{2}}}=1}{\mathop \sum }\,max(0,tlim{{2}_{i}}-\left| {{t}_{1}}-{{t}_{2}} \right|) \\ & {{g}_{5}}(x)=\overset{m}{\mathop{\underset{i=1}{\mathop \sum }\,}}\,\underset{\forall {{t}_{1}},{{t}_{2}},{{j}_{1}}\ne {{j}_{2}},x_{ij}^{{{t}_{1}}}\text{ }\!\!~\!\!\text{ }x_{ij}^{{{t}_{2}}}=1}{\mathop \sum }\,max(0,tlim{{3}_{i}}-\left| {{t}_{1}}-{{t}_{2}} \right|) \\ & {{g}_{6}}(x)=\overset{m}{\mathop{\underset{i={{m}_{1}}+1}{\mathop \sum }\,}}\,\overset{n}{\mathop{\underset{j=1}{\mathop \sum }\,}}\,\overset{T}{\mathop{\underset{t=1}{\mathop \sum }\,}}\,x_{ij}^{t}\times \overset{n}{\mathop{\underset{j=1}{\mathop \sum }\,}}\,\overset{{{m}_{1}}}{\mathop{\underset{i=1}{\mathop \sum }\,}}\,x_{ij}^{t}-1 \\ & {{g}_{7}}(x)=\overset{n}{\mathop{\underset{j=1}{\mathop \sum }\,}}\,\left[ \overset{m}{\mathop{\underset{{{i}_{1}}={{m}_{1}}+1}{\mathop \sum }\,}}\,\overset{T}{\mathop{\underset{{{t}_{1}}=1}{\mathop \sum }\,}}\,x_{{{i}_{1}}j}^{{{t}_{1}}}\times \overset{{{m}_{1}}}{\mathop{\underset{{{i}_{2}}=1}{\mathop \sum }\,}}\,\underset{{{t}_{2}},0<{{t}_{1}}-{{t}_{2}}<tUpper}{\mathop \sum }\,x_{{{i}_{2}}j}^{{{t}_{2}}}\text{ }\!\!~\!\!\text{ } \right] \\ \end{matrix}$
通过将ε约束法引入r-NSGA-Ⅱ算法,改进的r-NSGA-Ⅱ算法放宽了约束条件,能更好地利用不可行解,更加适用于多目标WTA问题。

3 仿真计算与结果分析

下面将模拟平台典型对空作战场景下的火力协同决策过程,采用NSGA-Ⅱ、NSGA-Ⅲ、MOEA/D、改进r-NSGA-Ⅱ算法这四类多目标进化算法进行仿真求解,并对仿真结果进行对比分析。

3.1 性能评价指标

本文采用Mohammadi[17]提出的基于复合前沿的反转世代距离评价指标(Inverted Generational Distance based on Composite Front, IGD-CF)用于衡量算法在体现指挥员决策偏好方面的能力。IGD-CF是一个综合性能的评价指标,可以综合评价算法在偏好区域内的收敛性和分布性,IGD-CF值越小,说明算法在偏好区域内的收敛性与分布性越好。
该指标评价过程如下:
1)生成复合前沿。将要比较的算法的解集进行合并并进行非支配排序,将得到的非支配解集作为复合前沿。
2)为参考点生成偏好区域。计算复合前沿中所有点与参考点之间的距离,将距离参考点最近的解作为中心点,在决策者指定半径r后,以中心点为圆心,r为半径的区域被视为偏好区域。
3)计算IGD值。
IGD-CF(P*,Q)= v P * d ( v , Q ) | P * |
其中,P*为偏好区域内的复合前沿,|P*|为偏好区域内复合前沿的解集个数,Q为算法获得的最优Pareto解集,d(v,Q)为偏好区域内复合前沿上的个体到种群Q的最小欧氏距离。

3.2 仿真结果

典型的对空作战场景一般会面临超声速反舰导弹、亚声速导弹、无人机混合来袭的战场态势,因此实验假定中存在三类来袭目标,共9批,来袭目标位置速度等信息已知。假定防空平台共搭载了三型动能武器以及一型高能激光武器。
实验设定的来袭目标与武器资源参数如下:表1 为威胁目标参数取值情况,表2 为武器资源参数取值情况,表3为武器对抗来袭目标的特征参数。
表1 威胁目标参数取值

Tab.1 Threat target parameter values

编号 类型 目标价值 位置/km 速度/Ma
1 A 3 (23,22,7) (2.1,2.2,0.9)
2 A 3 (27,20,8) (1.5,1.9,0.9)
3 A 3 (15,16,7) (3.2,2.2,1)
4 B 2 (17,16,6) (1.4,0.9,0.3)
5 B 2 (17,13,5) (0.6,0.6,0.3)
6 B 2 (12,14,4) (0.8,0.8,0.3)
7 C 1.5 (9,9,6) (0.4,1.1,0.2)
8 C 1.5 (7.7,11,6) (0.3,0.2,0.2)
9 C 1.5 (9.8,9.3,6.5) (0.25,0.62,0.2)
表2 武器资源参数

Tab.2 Weapon resource parameters

武器
编号
弹药
库存
武器
成本
弹药飞行
速度
连射
时间
转火
时间
1 60 50 2.5 1.5 0
2 60 20 2 1.5 0
3 100 5 1.5 1 0
4 100 1 Inf 0 3
表3 武器对抗来袭目标的毁伤概率与远近界

Tab.3 Destruction probability of weapons against incoming targets with near and far boundaries

目标类型 武器1 武器2 武器3 高能激光武器
A (0.9,[40,120]) (0.65,[30,90]) (0.55,[0.3,15]) p1=0.2,
B (0.0,[0,0]) (0.82,[25,90]) (0.62,[0.3,20]) p2=0.7
C (0.0,[0,0]) (0.77,[15,40]) (0.75,[0.3,20]) p3=0.4;[3,10]
全部被测算法中,种群数量同一设定为100,最大迭代次数设定为5 000,交叉概率设为1,变异期望设为1。所有算法与仿真实验均运行在一台使用Intel(R)Core(TM) i7-7700HQ CPU(4核,2.8 GHz),8 GB内存的电脑上。
实验对比了多种不同的多目标优化算法在典型对空作战场景模型中的优化性能表现,参与比较的指标包括:算法运行时间(Runtime),超体积(Hyper Volume, HV)[18]以及基于复合前沿的反转世代距离评价指标(IGD-CF)。其中HV指标综合反映了算法的收敛性与分布性,其值越大,性能表现越好;IGD-CF指标反映了算法在反映指挥员偏好方面的能力,其值越小,性能表现越好。
为消除随机误差,实验对典型防空反导场景进行了15次仿真。各算法具体表现如表4表5所示。
表4 最优/均值表现分析

Tab.4 Analysis of optimal/mean performance

目标函数 改进r-NSGA-Ⅱ NSGA-Ⅱ NSGA-Ⅲ MOEA/D
f1 最优 0.18 0.19 0.19 0.21
平均 0.28 0.3 0.3 0.32
f2 最优 0.15 0.15 0.19 0.28
平均 0.28 0.27 0.31 0.32
表5 试验参数分析

Tab.5 Analysis of test parameters

评价指标 改进r-NSGA-Ⅱ NSGA-Ⅱ NSGA-Ⅲ MOEA/D
最好 0 0 0 0
IGD-CF 最差 0.077 0.094 0.108 0.126
平均 0.028 0.033 0.051 0.062
最好 0.733 0.766 0.728 0.712
HV 最差 0.618 0.593 0.561 0.466
平均 0.682 0.672 0.647 0.636
最好 14.46 14.76 15.36 15.02
runtime 最差 16.5 15.72 14.34 12.96
平均 15.24 15.12 17.64 15.95
表5展现了各多目标优化算法在运行时间、HV指标以及IGD-CF指标上的性能表现。其中,对于超体积HV指标,由于本问题的实际Pareto前沿未知,实验采用基于蒙特卡洛方法的超体积分析方法[19],超体积分析选用种群中所有个体各个优化目标的最优值作为下界,选用{1,1}作为上界,每次仿真产生的采样点数量为1 000 000,采样点分布采用随机分布。
表4中,f1表示目标突防概率,f2表示武器打击成本。从表4可知,改进r-NSGA-Ⅱ算法能够找到具备最优突防概率与打击成本的参数组合。
表5可知,在全部的测试算法中,改进r-NSGA-Ⅱ算法在各种情况下均具备最好的IGD-CF指标性能表现,这表明改进r-NSGA-Ⅱ算法能够充分地体现指挥员的决策偏好,其种群在偏好点附近分布较为密集,能够提供更多符合指挥员偏好的决策方案;在全部的测试算法中,改进r-NSGA-Ⅱ算法在平均与最坏情况下具备最好的超体积指标性能表现,这表示在平均与最坏情况下,改进r-NSGA-Ⅱ算法能够覆盖最大范围的目标空间,具备良好的收敛性与分布性;改进r-NSGA-Ⅱ算法在运行时间方面与NSGA-Ⅱ,NSGA-Ⅲ算法基本持平,慢于MOEA/D算法,这表示改进r-NAGA-Ⅱ算法具备不弱于部分传统算法的运行速度,并在这方面依然具备改进空间。
基于表4表5的测试结果可得,改进r-NSGA-Ⅱ算法在测试的全部多目标优化算法中具备最好的优化效能表现。

4 结束语

为实现现代战场下防空反导武器-目标分配的准确建模求解,本文基于动态时间建立了多武器多目标决策分配模型,在充分考虑打击时间区段、打击模式和开火约束等信息的基础上,将高能激光作为防御武器引入模型,提高了模型的有效性与适用性,较以往模型有较大改进。本文提出改进r-NSGA-Ⅱ算法,将指挥员偏好体现为参考点与权重引入算法,在进化算法中引导种群进化方向,经过仿真实验,改进r-NSGA-Ⅱ算法具有很好的效果,算法具有不弱于传统算法的运算性能,并能为指挥人员提供更符合偏好的选择方案。
[1]
Alexander Kline, Darryl Ahner, Raymond Hill. The weapon-target assignment problem[J]. Computers & Operations Research, 2019(105): 226-236.

[2]
Manne, Alan S. A target-assignment problem[J]. Operations Research, 1958, 6(3): 346-351.

DOI

[3]
Kline, A, Ahner, D, Pachter, M. A greedy hungarian algorithm for the weapon-target assignment problem[R]. Air Force Institute of Technology Center for Operational Analysis, 2017.

[4]
Karasakal, O. Air defense missile-target allocation models for a naval task group[J]. Comput. Oper. Res, 2008, 35 (6):1759-1770.

DOI

[5]
XIN B, CHEN J, PENG Z, et al. An efficient rule-based constructive heuristic to solve dynamic weapon-target assignment problem[J]. IEEE Trans. Syst. ManCybern, 2011, 41 (3):598-606.

[6]
Robert A. Murphey. An approximate algorithm for a weapon target assignment stochastic program[J]. Approximation and Complexity in Numerical Optimization, 2000(42):406-421.

[7]
邵诗佳. 基于智能算法的武器目标分配问题研究[D]. 哈尔滨: 哈尔滨工程大学, 2019.

SHAO S J. Research on weapon target assignment based on intelligent algorithm[D]. Harbin: Harbin Engineering University, 2019.

[8]
邱少明, 冯江惠, 杜秀丽, 等. 基于改进多目标HQPSOGA求解武器目标分配问题[J]. 计算机应用与软件, 2021, 38(11):255-262.

QIU S M, FENG J H, DU X L, et al. Weapon target assignment based on improved multi-objective HQPSOGA[J]. Computer Applications and Software, 2021, 38(11):255-262.

[9]
陈曼. 粒子群算法在舰载武器-目标分配中的研究与应用[D]. 武汉: 武汉科技大学, 2018.

CHEN M. Research and application of particle swarm optimization in ship-borne weapon target assignment[D]. Wuhan: Wuhan University of Science and Technology, 2018.

[10]
K. Deb, A. Pratap, S. Agarwal, T. Meyarivan. A fast and elitist multiobjective genetic algorithm: NSGA-Ⅱ[J]. IEEE Transactions on Evolutionary Computation, 2002, 6(2):182-197.

DOI

[11]
L. Ben Said, S. Bechikh, K. Ghedira. The r-Dominance: a new dominance relation for interactive evolutionary multicriteria decision making[J]. IEEE Transactions on Evolutionary Computation, 2010, 14(5):801-818.

DOI

[12]
张波, 商豪. 应用随机过程[M].2版. 北京: 中国人民大学出版社, 2009.

ZHANG B, SHANG H. Applied stochastic processes[M]. Second Edition. Beijing: Press of Chinese People University, 2009.

[13]
K.Deb. An efficient constraint handling method for genetic algorithms[J]. Computer Methods in Applied Mechanics and Engineering, 2000, 186(2): 311-338.

DOI

[14]
Takahama, Tetsuyuki, Setsuko Sakai. Efficient constrained optimization by the ε constrained rank-based differential evolution[C]. 2012 IEEE Congress on Evolutionary Computation, 2012: 1-8.

[15]
YukiTakahama, Tetsu, Setsuko Sakai. Efficient constrained optimization by the ε constrained differential evolution with rough approximation using kernel regression[C]. 2013 IEEE Congress on Evolutionary Computation, 2013: 1334-1341.

[16]
李笠, 等. 约束进化算法及其应用研究综述[J]. 计算机科学, 2021, 48(4):1-13.

LI L. Survey of constrained evolutionary algorithms and their applications[J]. Computer Science, 2021, 48(4):1-13.

DOI

[17]
A. Mohammadi, M. N. Omidvar, X. Li. A new performance metric for user-preference based multi-objective evolutionary algorithms[C]. 2013 IEEE Congress on Evolutionary Computation, 2013:2825-2832.

[18]
Eckart Zitzler, Lothar Thiele. Multiobjective evolutionary algorithms: a comparative case study and the strength pareto approach[J]. IEEE Transactions on Evolutionary Computation, 1999, 3(4):257-271.

DOI

[19]
Bringmann K, Friedrich T. Approximating the least hypervolume contributor: NP-hard in general, but fast in practice[C]. Proceedings of the International Conference on Evolutionary Multi-Criterion Optimization. Berlin,Germany, 2009: 6-20.

Outlines

/