对于粒子群等优化算法,将物理背景抽象为数学模型,我认为主要在两方面:
(1)如何将研究对象抽象为粒子,怎样将粒子的属性和粒子所属空间用数学的方式去表达。
(2)如何将现实约束抽象为目标函数,只有有明确的目标函数,才会有优化的目标,粒子才能按照这个目标去寻优。
卫星任务地面站资源调度问题
阅读文献:
樊慧晶, 章文毅, 田妙苗,等.基于粒子群算法的卫星任务地面站资源调度方法[J]. 中国科学院大学学报, 2022, 39(6) : 801-808.
樊慧晶. 基于卫星数传与测控任务的地面站资源调度研究[D]. 北京: 中国科学院大学, 2021.
卫星地面站资源调度问题是指根据卫星的任务应用需求,在有限的时间范围和特定的约束条件下,为卫星的数传、测控等任务配置相应的地面站资源( 如天线、记录器等) ,得到满足卫星任务需求的资源分配调度方案。
若不考虑信息传输的具体细节,将背景简化一下,就是调度卫星与地面基站的任务分配。这里涉及到两种任务:数传任务和测控任务。数传指,卫星将获取到的数据下传到地面站(记录器);测控指,地面站天线向卫星发射信号,对卫星的有效载荷提出指令。
地面站资源调度问题涉及的变量:
(1)卫星任务集合 J = { j1,j2,j3,…,ji,…,jn } ;
(2)卫星集合 S = { s1,s2,s3,…,sp,…,sz} ,任务 ji 的卫星用 Sji 表示;
(3)地面站集合 G = { g1,g2,g3,…,gq,…,gm } ;
(4)天线集合 A = { a1,a2,a3,…,ak,…,ax } ,任务 ji 所用天线用 Aji 表示。
(5)记录器集合 R = { r1,r2,r3,…,rl,…,ry } ,任务 ji 所用记录器用 Rji 表示。
采用离散粒子群算法,设计离散粒子群的矩阵 C,N行(N个粒子),n列(粒子维数为n,每一维表示每一个卫星任务的完成情况),矩阵元素的值 cji(0或1)表示任务 ji 是否安排了符合约束的设备。
优化目标:
(1)任务执行率 f1 最大
cji 表示任务 ji 是否执行,若任务 ji 安排了符合约束的设备,那么 cji = 1,否则 cji = 0。
任务执行率:
m
a
x
f
1
=
∑
i
=
1
i
=
n
C
j
i
n
maxf_1= \frac {\sum\limits^{i=n}_{i=1}C_{j_{i}}} n
maxf1=ni=1∑i=nCji
(2)调度方案中所有任务的总时长 f2 最大
S
T
j
i
=
m
i
n
(
F
S
S
j
i
,
F
C
S
j
i
)
ST_{j_{i}}=min(FSS_{j_i},FCS_{j_i})
STji=min(FSSji,FCSji)
E T j i = m a x ( F S E j i , F C E j i ) ET_{j_i}=max(FSE_{j_i},FCE_{j_i}) ETji=max(FSEji,FCEji)
m a x f 2 = ∑ i = 1 i = n ( E T j i − S T j i ) max\ f_2=\sum\limits^{i=n}_{i=1}(ET_{j_i}-ST_{j_i}) max f2=i=1∑i=n(ETji−STji)
FSSji 、FCSji 、FSEji 、FCEji 分别表示任务的实际数传开始时间、实际测控开始时间、实际数传结束时间、实际测控结束时间。
(3)尽量为任务安排其优选设备
m
a
x
f
3
=
∑
i
=
1
i
=
n
(
(
E
T
j
i
−
S
T
j
i
)
×
(
1
/
S
A
(
S
j
i
,
A
j
i
)
)
+
(
F
S
E
j
i
−
F
S
S
j
i
)
×
(
1
/
S
R
(
S
j
i
,
R
j
i
)
)
)
max\ f_3 = \sum\limits^{i=n}_{i=1}((ET_{j_i}-ST_{j_i})×(1/SA(S_{j_i},A_{j_i}))+(FSE_{j_i}-FSS_{j_i})×(1/SR(S_{j_i},R_{j_i})))
max f3=i=1∑i=n((ETji−STji)×(1/SA(Sji,Aji))+(FSEji−FSSji)×(1/SR(Sji,Rji)))
f3 表示调度方案每个任务的实际时长与其设备首选程度的乘积,其值越大,代表调度方案的设备越优选。
其中,SA( Sji ,Aji ) 表示任务 ji 的卫星对其天线的首选性值,SR( Sji ,Rji ) 表示任务 ji 的卫星对其记录器的首选性值。由于设备首选性数值与其首选程度呈负相关,因此,对设备首选性数值取倒数表示设备的首选程度。
此外,在算法的过程中还需注意一些现实约束条件,这里列举主要几条:
(1)设备能力约束
天线能力与任务类型相匹配。即任务 ji 的天线 Aji 的能力 ATAji要满足任务类型 Tji,如果任务 ji 的类型 Tji 为单一测控任务,就只需要安排天线资源,不需要对该任务安排记录器。
(2)设备时间占用约束
①某一时刻,一部天线只能安排一颗卫星的任务:即在[STji ,ETji ]时间段内,卫星任务 ji 的天线 Aji 不能被其他任务所用。
②当一部天线完成一个任务,需要一定的切换时间 t 才能供下一任务使用:即在[ETji ,ETji + t]时间段内,卫星任务 ji 的天线 Aji 不能被其他卫星利用。
(还有其他的约束条件,比如记录器的逻辑码速率等等,这里暂不列举。)
这是一个多目标优化问题,作者没有采取帕累托最优方法,而是根据问题自身的特点采取了另外的策略:
由于问题的前两个优化目标(任务执行率最大和任务总时长最大),进化范围比较有限,如果通过算法进行进化或搜索,不仅会增加计算的时间成本,效果也并不明显。所以,作者先随机生成一定数量的初始解,并按照前两个优化目标(任务执行率最大和任务总时长最大)对其进行筛选,得到较优的初始种群,然后对优化目标3(尽量为任务安排其优选设备)进行算法的进化与寻优。(我认为这一点比较巧)
(这两篇文章花了一些时间才看明白,因为文章里没有明确说建立了的什么样的矩阵、算法的具体过程是什么,一开始不知道要将哪个变量抽象为粒子,也搞不明白这么多纷繁的目标和约束是怎么回事。反复看了几遍大概知道了,粒子群矩阵所抽象的是任务分配情况(0或1),其他的数据比如卫星数、地面站数、时长、首选性等等,是提前设定好的参数(应该)。但其实这个就是课题的难点之一所在,属性很多,约束很多,要将它们梳理好并反映在算法中,估计代码里有很大一部分是在处理这些约束条件。)
无人机协同多任务分配问题
文献阅读:王峰, 张衡, 韩孟臣, 等.基于协同进化的混合变量多目标粒子群优化算法求解无人机协同多任务分配问题[J].计算机学报, 2021, 44(10): 1967-1983.
无人机协同多任务分配:多架无人机通过协同合作分别对多个军事目标进行观测、打击和评估任务,最终完成整个军事行动。
一个无人机协同多任务分配模型可表示为一个四元组< U,T,M,C >,其中,U 为无人机集合,T为军事区域内存在的目标集合,M 是所有目标上的任务集合,C 为问题模型中的约束集合。无人机、目标和任务是模型的三个重要对象。
任务包含三种类型:观测任务(Observe)、打击任务(Attack)、打击结果评估任务(Evaluate)。
每一个军事目标包含以上三种任务。
无人机分为两种,战斗无人机和侦察无人机;其中,战斗无人机只能执行打击任务,侦察无人机可以执行观测任务和评估任务。无人机集合包含战斗无人机和侦察无人机。
无人机完成任务都需要消耗不同数量和类型的资源。战斗无人机完成打击任务需要消耗一定数目的导弹,侦察无人机完成观测和评估任务需要消耗一定的侦察时间。
无人机任务分配计划举例:
对于两个目标(T1、T2),分配 3 架侦察无人机和 4 架战斗无人机。
从表中可以看出,第一个目标的观测任务被分配给第 1 架和第 2 架侦察无人机执行,两架无人机执行目标 T1 的观测任务时一共消耗了 30 个单位时间. 对第一个目标的打击任务则由第 1 架和第 3 架战斗无人机协同完成,两架无人机在执行该任务时一共消耗了 4 枚导弹。
粒子编码
采用混合编码的方式,每个粒子的编码包含三个部分,即 P = ( MN,UN,C ). 其中,第一部分 MN 是离散变量,记录任务编号;第二部分 UN 是离散变量,记录执行对应任务的无人机编号;第三部分 C 是混合变量,记录无人机执行任务时的资源消耗。
编码举例:
对于一个包含 2 个目标,1 架战斗无人机和 1 架侦察无人机的任务分配模型,目标 T1上存在三种任务 ( M1,M2,M3 ),目标 T2 上存在三种任务 ( M4,M5,M6 ),无人机 U1 为侦察无人机,无人机 U2 为战斗无人机。
粒子的位置向量的前两行表示了任务分配结果。任务M1,M3,M4 和 M6 被分配给无人机 U1 执行.。此时,对应的任务分配计划向量中的X1k =1 ( k ∈{ 1,3,4,6 } );任务 M2 和 M5 被分配给了无人机 U2 执行,此时,对应的任务分配计划向量中的 X1k = 1 ( k ∈{ 2,5 } )。
粒子的位置向量的第三行则表示了无人机执行每个任务时的资源消耗结果,该向量中既有表示侦察无人机执行观测任务和打击结果评估任务时间消耗的连续变量,也有表示战斗无人机执行打击任务消耗导弹数目的离散变量。
此外,为了表示无人机的路径规划结果,粒子的位置向量中任务编号出现的顺序代表了无人机执行任务的顺序。黑色实心箭头表示无人机 U1 执行任务依次为 M1 M4 M6 M3 ,黑色空心箭头表示无人机 U2 的执行任务依次为 M2 M5。
在初始化时,无人机和无人机执行的任务将随机地从所有任务集合和所有无人机集合中选取。
约束条件
(1)战斗无人机最大携弹数目约束
单个战斗无人机执行打击任务消耗的导弹数要小于其最大携弹数目。
∑
k
=
1
N
m
C
i
k
X
i
k
≤
L
o
a
d
i
\sum\limits^{Nm}_{k=1}C_i^kX_i^k\leq Load_i
k=1∑NmCikXik≤Loadi
(2)无人机类型约束
侦察无人机只能完成观测任务和评估任务,而战斗无人机只能完成打击任务。
(3)任务执行时序约束
军事任务要求打击任务在观测任务完成后才能开始执行,而打击结果评估任务需要在打击任务完成后开始。
(4)多机协同约束
为了保证针对目标的军事任务顺利完成,该目标的每种类型的任务需要分配给至少一架无人机执行,即:
∑
i
=
1
S
X
i
k
⋀
∑
i
=
S
+
1
S
=
F
X
i
k
+
1
⋀
∑
i
=
1
S
X
i
k
+
2
=
1
\sum\limits^S_{i=1}X_i^k \bigwedge \sum\limits^{S=F}_{i = S+1}X_i^{k+1}\bigwedge\sum\limits^S_{i=1}X_i^{k+2}=1
i=1∑SXik⋀i=S+1∑S=FXik+1⋀i=1∑SXik+2=1
(5)任务需求约束
为了保证每个任务被成功完成,执行任务的无人机消耗的资源需要满足任务完成所要求的条件。
对于打击任务,执行该打击任务的战斗无人机消耗的导弹总数需要大于该任务完成所需导弹数;对于观测任务和打击结果评估任务,侦察无人机执行任务的时间消耗需要大于该任务完成所需执行时间。
(6)最大航程约束
无人机执行任务时的总飞行距离应小于无人机的最远飞行距离。
优化目标
同时考虑无人机飞行总航程和任务完成时间两个优化指标。
(1)无人机飞行总航程 f1
要求无人机执行任务时的飞行总航程最短,保证无人机执行任务时消耗的飞行资源最少。
设无人机 Ui 依次执行的任务序列为 Seqi = { M1i ,M2i …,M iNi},Dis 为无人机 Ui 执行第 k 个任务到第 k + 1个任务的飞行距离,则无人机飞行总航程为
f
1
=
∑
i
=
1
N
u
∑
M
k
i
∈
S
e
q
i
D
i
s
M
k
i
M
k
+
1
i
f1=\sum\limits^{Nu}_{i=1}\sum\limits_{M_k^i \in Seq^i}Dis_{M_k^i}^{M_{k+1}^i}
f1=i=1∑NuMki∈Seqi∑DisMkiMk+1i
(2)任务完成时间 f2
要求所有任务的完成时间最短,保证整个任务能被快速执行完毕。
f
2
=
m
a
x
1
≤
k
≤
N
m
t
k
e
f2=\mathop{max}_{1\leq k \leq Nm}t_k^e
f2=max1≤k≤Nmtke
其中,tke 表示任务 Mk 的完成时间。所有目标任务的最大完成时间即代表军事行动的任务完成时间。
对于多目标优化,文章采用基于帕累托支配的方法选择最优个体。
(个人感觉这篇文章条理、逻辑、内容组织非常清晰,其中,粒子编码采用向量的形式,可以将大量的信息囊括在一个粒子中。)
多导弹动态武器目标分配问题
文献阅读:刘攀, 徐胜利, 张迪, 等.基于粒子群优化的多导弹动态武器目标分配算法[J].南京航空航天大学学报, 2023, 55(1): 108 - 115.
武器目标分配(Weapon target assignment, WTA)是根据目标威胁的大小和运动状态,按照期望目标(如最小化对防御资源毁伤概率或最大化对目标拦截概率)进而确定武器 - 目标之间的分配关系,形成对目标的分配策略。
进攻方战机组成战斗编队,对防御方战机造成实质性威胁;防御方战机作为空中移动武器发射平台,其每架战机均携带有数量一定的对空导弹,可以对进攻方战机进行拦截。
在目标分配之前,需要对作战双方进行态势评估,通过建立一套科学合理的评判标准,实时反映双方的威胁和优势信息。然后综合威胁和优势信息进行武器目标分配。
粒子编码
将粒子群中的每个粒子的位置视为一个导弹 - 进攻方战机的目标分配方案。
设置粒子群的粒子个数为 m,粒子的空间维数为 Z。第 i 个粒子的位置向量 Xi、速度向量 Vi 分别表示为
{
X
i
=
(
x
i
1
,
.
.
.
,
x
i
r
,
.
.
.
,
x
i
Z
)
V
i
=
(
v
i
1
,
.
.
.
,
v
i
r
,
.
.
.
,
v
i
Z
)
\begin{cases} X_i=(x_{i1}, ..., x_{ir}, ..., x_{iZ}) \\ V_i\ = (v_{i1}, ..., v_{ir}, ..., v_{iZ}) \end{cases}
{Xi=(xi1,...,xir,...,xiZ)Vi =(vi1,...,vir,...,viZ)
式中,位置向量 Xi 与多导弹目标分配方案对应,位置向量中具体数值 xir 为 1~ J 的正整数,即为第 r 枚导弹分配的进攻方战机编号 ,Z 为导弹总数;i ∈( 1,2,⋯,m ),r ∈( 1,2,⋯,Z )。
优化目标
综合考虑威胁和优势信息,每次以进攻方战机生存概率最小的原则进行防御方战机导弹的分配,则第 k 次分配中目标优化函数为
E
k
(
π
)
=
a
r
g
m
i
n
∑
j
=
1
J
k
∑
i
=
1
I
k
{
T
j
i
k
⋅
[
∏
i
=
1
I
k
(
1
−
S
i
j
k
)
X
i
j
k
]
}
,
π
∈
Ω
k
E^k(\pi)=arg\ min\sum\limits^{J^k}_{j=1}\sum\limits^{I^k}_{i=1}\lbrace{T^k_{ji}·[\prod^{I^k}_{i=1}(1-S^k_{ij})X^k_{ij}]}\rbrace \ ,\ \pi \in \Omega^k
Ek(π)=arg minj=1∑Jki=1∑Ik{Tjik⋅[i=1∏Ik(1−Sijk)Xijk]} , π∈Ωk
式中,Ωk 为当前所有满足约束条件的防御方战机 - 进攻方战机分配方案构成的集合。Tkji 和 Skji 分别为防御方战机 Bi 和进攻方 Rj 之间的威胁和优势指数。Ik 和 Jk 分别为每次分配时防御方和进攻方战机的总数。
文章中给出了威胁指数 Tji 的计算方法:
T j i = w 1 ⋅ T a j i + w 2 ⋅ T v j i + w 3 ⋅ T h j i + w 4 ⋅ T r j i + w 5 ⋅ T i n j i T_ji=w_1·Ta_{ji}+w_2·Tv_{ji}+w_3·Th_{ji}+w_4·Tr_{ji}+w_5·Tin_ji Tji=w1⋅Taji+w2⋅Tvji+w3⋅Thji+w4⋅Trji+w5⋅Tinji
其中,Taji、Tvji、Thji、Trji、Tinji 分别表示角度威胁因子、速度威胁因子、高度威胁因子、距离威胁因子和隐身性能威胁因子。文中有具体计算方法,这里不再展开。wi (i = 1, 2, 3, 4, 5) 为各威胁因子权值,且由层次分析法,满足下式要求:
∑ i = 1 5 w i = 1 w i > 0 \sum\limits^5_{i=1}w_i=1\ \ \ w_i>0 i=1∑5wi=1 wi>0
以及关于优势指数 Sij 的解释:第 k 次分配时防御方战机 Bi 对进攻方战机 Rj 的优势指数 Skij 等同于对进攻方战机 Rj 的毁伤概率,则第 k 次分配到下次分配前,进攻方战机的生存概率为
A j k = ∏ i = 1 I k ( 1 − S i j k ) X i j k A^k_j=\prod^{I^k}_{i=1}(1-S^k_{ij})X^k_{ij} Ajk=i=1∏Ik(1−Sijk)Xijk
Ik 表示第 k 阶段防御方可参与目标分配的战机总数,Xkij 为布尔值。当防御方战机 Bi 攻击进攻方战机 Rj 时,Xkij = 1;当防御方战机 Bi 不攻击进攻方战机 Rj 时,Xkij = 0。每次目标分配完成后进行战场评估,计算进攻方战机的累计生存概率,若大于生存阈值,则判定该战机被击中;否则仍将其作为待拦截目标。
当进攻方战机剩余数量为零或防御方战机剩余导弹数为零时,不再进行武器目标分配决策,空战结束。
约束条件
(1)打击时间窗约束
假设第 j 架进攻方战机在时间 tij 被第 i 架防御方战机拦截,打击时间约束如下
t
i
m
e
i
n
≤
t
i
j
≤
t
i
m
e
o
u
t
time_{in}\leq t_{ij} \leq time_{out}
timein≤tij≤timeout
式中,timein 为第 j 架进攻方战机抵达防御方第 i 架战机射程最远处的时间;timeout 为第 j 架进攻方战机飞离防御方第 i 架战机最大射程的时间。防御方战机对进攻方战机拦截时间应满足这一约束条件。
(2)转火时间窗约束
防御方战机在发射一次导弹后需经过一定的时间以完成武器转火,之后才能参与下一次分配。
无人机群协同航迹规划问题
阅读文献:陈霁月.基于改进粒子群算法和 PRM 算法的无人机群协同航迹规划研究[D].哈尔滨: 哈尔滨工业大学, 2022.
无人机的航迹规划是针对无人机复杂约束条件和任务需求进行多路径规划,规划出的航迹需要满足无人机本身的物理机能约束限制,同时不与环境中的地形障碍或其他威胁发生碰撞,保证无人机沿规划航迹飞行时的安全性和对任务的完成率,即要求无人机能够以最小的代价安全完成指定任务。
粒子编码
设粒子群中共包含 N 个粒子,无人机航迹规划期望得到的航迹点数量为 d,则第 i 个粒子的位置可以表示为
x
i
=
(
x
i
,
1
,
x
i
,
2
,
x
i
,
d
)
T
=
(
(
x
,
y
,
z
)
i
,
1
,
(
x
,
y
,
z
)
i
,
2
,
.
.
.
,
(
x
,
y
,
z
)
i
,
d
)
T
x_i=(x_{i, 1},x_{i, 2},x_{i, d})^T=((x,y,z)_{i, 1},(x, y,z)_{i,2},...,(x,y,z)_{i,d})^T
xi=(xi,1,xi,2,xi,d)T=((x,y,z)i,1,(x,y,z)i,2,...,(x,y,z)i,d)T
其中,(x, y, z)i, j, i =1, 2, …, N, j = 1, 2, …, d 表示第 i 个粒子中第 j 个航迹点的位置坐标。
第 i 个粒子的速度可以表示为
v
i
=
(
v
i
,
1
,
v
i
,
2
,
.
.
.
,
v
i
,
d
)
T
=
(
(
v
x
,
v
y
,
v
z
)
i
,
1
,
(
v
x
,
v
y
,
v
z
)
i
,
2
,
.
.
.
,
(
v
x
,
v
y
,
v
z
)
i
,
d
)
T
v_i=(v_{i, 1},v_{i,2},...,v_{i,d})^T=((v^x,v^y,v^z)_{i,1},(v^x,v^y,v^z)_{i,2},...,(v^x,v^y,v^z)_{i,d})^T
vi=(vi,1,vi,2,...,vi,d)T=((vx,vy,vz)i,1,(vx,vy,vz)i,2,...,(vx,vy,vz)i,d)T
其中,(vx, vy, vz)i, j, i = 1, 2, …, N, j = 1, 2, …, d 表示第 i 个粒子中第 j 个航迹点的在三维空间中的速度。
约束条件
(作者对山区环境进行了随机三维地形建模,在粒子运动过程中,会实时接收到周围环境信息,下述的俯仰角、水平转弯角、转弯半径是由粒子所处地形而计算得到的数据)
(1)飞行高度约束
设满足避障的最低飞行高度为 Hmin ,满足安全性的最高飞行高度为 Hmax ,无人机在各个航迹段的飞行过程中,高度 Hi 应满足
H
m
i
n
≤
H
i
≤
H
m
a
x
,
(
1
=
1
,
2
,
3
,
.
.
.
,
N
)
H_{min} \leq H_i \leq H_{max},\ (1=1, 2,3,...,N)
Hmin≤Hi≤Hmax, (1=1,2,3,...,N)
(2)俯仰角约束
无人机飞行姿态的改变通常发生在航迹点处,设当前航迹点为 i,该航迹点处飞行方向发生改变的俯仰角为 γi ,最大俯仰角为 γmax,则应满足
γ
i
<
γ
m
a
x
,
(
i
=
1
,
2
,
3
,
.
.
.
N
)
\gamma_i<\gamma_{max}, \ (i=1,2,3,... N)
γi<γmax, (i=1,2,3,...N)
(3)水平转弯角约束
设当前航迹点为 i,该航迹点处飞行方向发生改变的水平转弯角为 αi,最大转弯角为 αmax,则应满足
α
i
<
α
m
a
x
,
(
i
=
1
,
2
,
3
,
.
.
.
,
N
)
\alpha_i<\alpha_{max}, \ (i = 1, 2, 3, ..., N)
αi<αmax, (i=1,2,3,...,N)
(4)转弯半径约束
设无人机在航迹点处发生转弯的半径为 ri,受到物理限制无人机自身可以完成转弯的最小半径为 rmin ,则有
r
i
≥
r
m
a
x
,
(
i
=
2
,
3
,
.
.
.
,
N
−
1
)
r_i \geq r_{max}, \ (i=2,3,...,N-1)
ri≥rmax, (i=2,3,...,N−1)
优化目标
无人机的航迹代价是对航迹优劣的评估,代价值越低说明规划的航迹越好。
航迹代价可以概括为无人机自身带来的代价约束和外界带来的代价,其中无人机自身产生的代价包括航迹能耗代价、高度代价和转向代价,外界威胁带来的代价主要包括地形代价、雷达威胁代价、多无人机之间的协同避碰代价和通信代价,综合这些代价可以得到综合航迹代价模型。
(1)航迹能耗代价 fL
f
L
=
∑
i
=
1
N
−
1
a
i
L
i
f_L=\sum\limits^{N-1}_{i=1}a_iL_i
fL=i=1∑N−1aiLi
其中,ai 为能耗系数,Li 表示第 i 段航迹。
(2)高度威胁代价 fHi
f
H
i
=
{
H
i
−
H
m
i
n
H
m
a
x
−
H
m
i
n
H
m
i
n
<
H
i
<
H
m
a
x
I
其他
f_{H_i}=\begin{cases} \frac {H_i-H_{min}} {H_{max}-H_{min}} \ \ \ H_{min}<H_i<H_{max} \\ I \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ 其他 \end{cases}
fHi={Hmax−HminHi−Hmin Hmin<Hi<HmaxI 其他
其中,Hmin为最低航迹高度,Hmax为最高航迹高度,Hi为当前航迹高度;I 为常数,为高度威胁代价的惩罚因子,惩罚因子的大小远大于无人机的航迹代价。若无人机飞行高度不在约束范围内,则会获得大小为 I 的代价值,由于其远大于航迹代价,此航迹将不可行。
(3)转向代价 fαi
f
α
i
=
{
0
θ
≤
ϕ
/
2
k
θ
i
ϕ
/
2
<
θ
≤
ϕ
f_{\alpha_i}=\begin{cases} 0 \ \ \ \ \ \ \ \theta \leq \phi/2 \\ k\theta_i \ \ \ \phi/2 < \theta \leq \phi \end{cases}
fαi={0 θ≤ϕ/2kθi ϕ/2<θ≤ϕ
其中,k 为常数,为转向角度代价的惩罚因子,且 k > 1;θi 表示第 i 个节点处发生的转向角度,包括俯仰角和水平转弯角,ϕ 为最大转向角度。
(4)转弯代价 fri
f
r
i
=
{
0
r
i
>
r
m
i
n
C
其他
f_{r_i}=\begin{cases} 0 \ \ \ \ \ \ \ r_i > r_{min} \\ C \ \ \ \ \ 其他 \end{cases}
fri={0 ri>rminC 其他
其中,惩罚因子 C 为常数,其大小远大于无人机的航迹代价。
外部威胁代价:
(5)雷达威胁代价 fR
f
R
=
∑
i
=
2
N
∑
j
=
1
n
A
i
,
j
f_R=\sum\limits^N_{i=2}\sum\limits^n_{j=1}A_{i, j}
fR=i=2∑Nj=1∑nAi,j
A i , j = { ( δ / B i , j ) 4 B i , j ≤ R k 0 其他 A_{i, j}=\begin{cases} (\delta/B_{i, j})^4 \ \ \ B_{i, j} \leq R_k \\ 0 \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ 其他 \end{cases} Ai,j={(δ/Bi,j)4 Bi,j≤Rk0 其他
B i , j = ( x i , j − x k ) 2 + ( y i , j − y k ) 2 + ( z i , j − z k ) 2 B_{i,j}=\sqrt{(x_{i,j}-x_k)^2+(y_{i,j}-y_k)^2+(z_{i,j}-z_k)^2} Bi,j=(xi,j−xk)2+(yi,j−yk)2+(zi,j−zk)2
其中, δ 表示雷达强度。
(6)多无人机避碰代价 fdis
f
d
i
s
=
∑
I
=
1
N
−
1
∑
j
=
2
N
C
i
,
j
f_{dis}=\sum\limits^{N-1}_{I=1}\sum^N_{j = 2}C_{i,j}
fdis=I=1∑N−1j=2∑NCi,j
C i , j = { 1 d i j ≤ d s a f e 0 其他 C_{i,j}=\begin{cases} 1 \ \ \ d_{ij} \leq d_{safe} \\ 0 \ \ \ 其他 \end{cases} Ci,j={1 dij≤dsafe0 其他
(7)多无人机通信代价 fcom
f
c
o
m
=
∑
I
=
1
N
−
1
∑
j
=
2
N
D
i
,
j
f_{com}=\sum\limits^{N-1}_{I=1}\sum^N_{j = 2}D_{i,j}
fcom=I=1∑N−1j=2∑NDi,j
D i , j = { 0 可通信 1 不可通信 D_{i,j}=\begin{cases} 0 \ \ \ 可通信 \\ 1 \ \ \ 不可通信\end{cases} Di,j={0 可通信1 不可通信
其中,N 表示总航迹点个数,每个航迹点之间被分为 n 段,Di, j表示无人机 i 与 j 每个航迹点之间的通信代价。
则无人机的综合航迹代价为
f
=
k
1
f
L
+
k
2
∑
i
=
1
N
f
H
i
+
k
3
∑
i
=
1
N
f
α
i
+
k
4
∑
i
=
1
N
f
r
i
+
k
5
f
R
+
k
6
f
d
i
s
+
k
7
f
c
o
m
f = k_1f_L+k_2\sum\limits^N_{i=1}f_{H_i}+k_3\sum\limits^N_{i=1}f_{\alpha_i}+k_4\sum\limits^N_{i=1}f_{r_i}+k_5f_R+k_6f_{dis}+k_7f_{com}
f=k1fL+k2i=1∑NfHi+k3i=1∑Nfαi+k4i=1∑Nfri+k5fR+k6fdis+k7fcom
其中,k1、k2、k3、k4、k5、k6 和 k7 分别为各项代价的比例因子,通过调整比例因子可以设置对不同约束条件的需求程度,进而改变航迹的收敛效果。
认知无线电(Cognitive Radio)系统自适应资源分配
阅读文献:Yi Yang, Qinyu Zhang, Ye Wang, et al. Adaptive Resources Allocation Algorithm Based on Modified PSO for Cognitive Radio System[J]. CHINA COMMUNICATIONS, 2019, 16(5): 83-92.
认知无线电是一种智能的频谱感知和分配技术,能够智能地感知无线通信环境,并根据一定的学习决策算法自适应地改变频谱分配参数,动态地检测和有效地利用空闲频谱。频谱分配是认知无线电的主要部分,主要解决了多个用户间不矛盾的选择频率并最大化频谱效率的问题。
粒子编码
有 N 个子载波,将其中的部分同时分配给 K 个次级用户(Secondary User,SU),也就是 K <= N。
采用离散粒子群优化算法,令群体由 m 个粒子组成,每个粒子的维数为 K,每一维表示分配给了哪个频谱(子载波)。
约束条件
(1)每个用户的传输功率不能超过系统所要求的总传输功率限制
∑
k
=
1
K
∑
n
=
1
N
p
k
,
n
≤
P
t
o
t
a
l
\sum\limits^K_{k=1}\sum\limits^N_{n=1}p_{k,n} \leq P_{total}
k=1∑Kn=1∑Npk,n≤Ptotal
(2)每个子载波的传输功率为正
p
k
,
n
≥
0
p_{k, n} \geq 0
pk,n≥0
(3)每个子载波只能分给一个用户
(4)保证用户之间的公平约束
优化目标
(1)最大化信道容量 R
R
k
(
n
)
=
ρ
k
,
n
N
l
o
g
2
(
1
+
p
k
,
n
h
k
,
n
2
Γ
N
0
B
/
N
)
R_k(n)=\frac {\rho_{k, n}} N log_2(1+\frac {p_{k,n}h^2_{k,n}} {\Gamma N_0B/N})
Rk(n)=Nρk,nlog2(1+ΓN0B/Npk,nhk,n2)
m a x ∑ k = 1 K ∑ n = 1 N R k ( n ) max \ \sum\limits^K_{k=1}\sum\limits^N_{n=1}R_k(n) max k=1∑Kn=1∑NRk(n)
其中,ρk,n 为[0, 1]离散二进制值,值为 0 表示未分配,值为 1 表示分配;pk,n 为传输功率;B 是频谱带宽;hk,n 是信道增益;Γ 是信噪比间隙;N0 是高斯白噪声功率谱密度。
(2)用户之间公平性最佳
公平性表达为
Φ
=
(
∑
k
=
1
K
R
k
)
2
/
K
∑
k
=
1
K
(
R
k
)
2
\Phi=(\sum\limits^K_{k=1}R_k)^2/K\sum^K_{k=1}(R_k)^2
Φ=(k=1∑KRk)2/Kk=1∑K(Rk)2
(感觉看外文文献还是比较有难度的,里面有不少不太能看得懂,需要结合中文去理解)
频率指配(Frequency Assignment,FA)问题
阅读一
文献阅读:薛寒. 战场频率指配算法研究[D]. 长沙: 国防科学技术大学, 2015.
常见复杂电磁环境下频率指配的问题可以表述为,在已给定的有限可用频率点(段)的限制条件下,为一定数量的用频设备找到一种干扰代价最小的频率指配方案,以尽可能多地满足现有的频率约束关系。
粒子编码
采用整数编码的方式,将所有可用频点(假设有 M 个)按照频率值由小到大的顺序进行整数编号。
假设指配问题模型中有 N 台用频设备,每个用频设备 i 的可用频率为一个正整数集合 fi (i=1,2,…,N),他们的可用频率范围形成 N 个整数子集。
每一个指配方案 f = (f1, f2, …, fN)表示为一个正整数序列。每一个指配方案 f 视为一个粒子,即每一个粒子都是一个正整数序列。
例如,可用的频点为 M=30 个,编号 1~30,地域内用频设备 N=5 台,每台设备的可用频率均为全频点,即 5 个整数子集都是相同的{1,2,…,30},指配方案 f=(2,12,14,25,7)就表示 1 号~5 号的用频设备分别被指配了编号 2,12,14,25,7 的频率点。
也就是在算法中,构建这样的粒子矩阵,n行(n个粒子),N列(粒子维数为N,每一维表示该设备所分配的频率点是哪一个)。
粒子所属的空间为可用频率点组成的空间,限制了粒子的位置最大值和位置最小值。
目标函数
目标:干扰代价最小。
常见的干扰类型有三类:同频干扰、邻频干扰、互调干扰。
分别建立相应的干扰约束矩阵,用矩阵中的元素值表示其行和列所代表的设备间的干扰约束关系,得到同频干扰、邻频干扰和互调干扰的三个约束矩阵 A、B、C。
假设同频约束矩阵 A 中某元素 aij 为 1,则表示设备 i 和 j 在使用相同频率时将存在同频干扰;反之,若元素 aij 为 0,则表示设备 i 和 j 可以指配相同频率且不会产生干扰。
在某一指配方案 s 下,用 E(s) 表示指配方案 s 的干扰代价( 其中 α、β、γ 为分别为同频干扰、邻频干扰、互调干扰三种约束的惩罚因子,即各种干扰的加权值):
E
(
s
)
=
α
∑
i
=
1
N
∑
j
=
i
+
1
N
A
i
j
+
β
∑
i
=
1
N
∑
j
=
i
+
1
N
B
i
j
+
γ
∑
i
=
j
N
∑
j
=
m
N
∑
m
=
1
N
C
i
j
m
E(s)=\alpha \sum \limits^N_{i = 1} \sum\limits^N_{j=i+1}A_{ij}+\beta\sum\limits^N_{i=1} \sum\limits^N_{j=i+1}B_{ij}+\gamma\sum \limits^N_{i=j}\sum\limits^N_{j=m}\sum\limits^N_{m=1}C_{ijm}
E(s)=αi=1∑Nj=i+1∑NAij+βi=1∑Nj=i+1∑NBij+γi=j∑Nj=m∑Nm=1∑NCijm
(1)同频干扰约束代价函数Aij
同频干扰指落在设备接收机通带内的干扰信号与有用信号落在了同样的频率点上。同频约束即一对用频设备间除非有足够的间距、合格的屏蔽条件或者合适的地形因素,否则不能指配相同的频率。
A
i
j
=
{
a
i
j
,
f
(
i
)
=
f
(
j
)
且
a
i
j
≠
0
0
,
其他
A_{ij}=\begin{cases} a_{ij},\ f(i)=f(j) 且 a_{ij}\neq 0 \\ 0, \ \ \ 其他\end{cases}
Aij={aij, f(i)=f(j)且aij=00, 其他
(2)邻频干扰约束代价函数Bij
邻频干扰指来自其他用频设备使用相邻或相近频道造成的干扰。邻频约束即一对设备除非处在合适的规避因素下,否则不能指配相邻或相近的频率。
B
i
j
=
{
b
i
j
,
∣
f
(
i
)
−
f
(
j
)
∣
<
k
且
b
i
j
≠
0
,
k
为频道间隔
0
,
其他
B_{ij}=\begin{cases} b_{ij}, \ |f(i)-f(j)|<k且b_{ij}\neq 0,\ k为频道间隔 \\ 0,\ \ \ 其他 \end{cases}
Bij={bij, ∣f(i)−f(j)∣<k且bij=0, k为频道间隔0, 其他
(3)互调干扰约束代价函数Cij
主要由设备接收机中的非线性器件引起,若有两个或多个不同频率的干扰信号同时进入接收机,将可能在非线性器件的作用下产生加性或减性的新频率,这一频率落入接用频设备收机通带与有用信号重叠则造成互调干扰。
C
i
j
m
=
{
c
i
j
m
,
2
f
(
i
)
−
f
(
j
)
=
f
(
m
)
且
c
i
j
m
≠
0
0
,
其他
C_{ijm} = \begin{cases} c_{ijm}, \ 2f(i)-f(j)=f(m)且c_{ijm}\neq 0 \\ 0, \ \ \ \ \ 其他 \end{cases}
Cijm={cijm, 2f(i)−f(j)=f(m)且cijm=00, 其他
干扰代价 E(s) 即在指配方案 s 的情况下违反同频、邻频、互调三种约束条件的代价函数的加权和。
(标记一个问题,之后再具体看,频率点是离散的整数,而标准粒子群算法是适用于连续空间的问题的,这个怎么处理?)
阅读二
文献阅读:刘铭, 查淞, 黄纪军, 等. 基于多目标优化的联合作战用频规划方法[J]. 电波科学学报, 2022, 37(3): 434-442.
文章在引言里指出:然而,目前在用频规划问题相关文献中,建模通常只设有频域一维的变量,得到的用频规划方案仅包含各用频设备频点的指配方案,无法实现对用频设备时域、空域、频域、能域等多维度的统筹调度。同时,针对联合作战场景的用频规划问题通常基于单目标或多目标加权建模,难以精确描述决策者对实际战场上用频规划方案执行预期的考量与要求,可能导致解的优劣判定的偏差,造成理想用频规划方案的丢失。
(正好对应了上一篇文章里的不足。)
文章建立了多目标的的联合作战用频规划模型,并采用带精英策略的非支配排序蚁群算法求解问题。
问题编码
战斗集团由 Np 个作战平台组成,各作战平台分别列装一台或多台各型通信、导航、雷达、制导、电子战等用频装备。战斗集团的作战过程可依次序拆解为 Nr 个作战阶段。(将战斗集团的作战过程依次序拆解为 NT 个作战 阶段,目的是使得同一频点可以被指配给不在同一作战阶段使用的多台用频装备,从而实现频谱资源时间维度上的复用。同一台用频装备,原则上从作战开始到终止的全部作战阶段均使用同一频点,尽可能避免换频的操作。)
问题中包含的信息有:
(1)作战平台标识。P = { p1 , p2 , … , p(Np) } ,共 Np 个作战平台。
(2)用频装备标识。W = { w1, w2 , … , w(Nw) },共 Nw 台用频装备。
(3)用频装备所处平台序号。S = { s1, s2, … , s(Nw) } 表示各用频装备所处作战平台标识号的集合。
(4)用频需求优先级。 Ri 表示用频装备 wi 的用频需求优先级。
(5)干扰门限。 Di 表示用频装备 wi 的接收机干扰门限 (单位:dB)。
(6)作战阶段标识。T = { t1, t2, … , t(NT) } 表示 NT 个作战阶段。
(7)可用频谱资源。F = { f1, f2, … , f(NF) }表示全体可用频点,NF为全体可用频点数量。
(8)可选频率信息。Sfi = { sf1, sf2, …, sf(N(sfi)) } 表示用频装备 wi 的可选频点集合,N(sfi) 为用频装备 wi 可选频点数量。
(9)可选发射功率信息。Spi = { sp1, sp2, …, sp(N(spi)) } 表示用频装备 wi 的可选功率值集合,N(spi) 为用频装备 wi 可选功率值数量。
(10)可选地理坐标信息。Sci = { sc1, sc2, …, sc(N(Sci)) } 表示作战平台 pi 的可选横纵地理坐标集合,N(Sci) 为作战平台 pi 可选横纵地理坐标数量。
每只蚂蚁走过的城市用元任务表示。
设置元任务 Uijkv 用以描述单个用频装备的用频规划方案,采用以下的六元组来表示元任务:
U
i
j
k
v
=
{
E
q
u
i
p
m
e
n
t
I
d
,
S
t
a
g
e
I
d
,
F
r
e
q
u
e
n
c
y
,
P
o
w
e
r
,
C
o
r
d
x
,
C
o
r
d
y
}
U_{ijkv}=\lbrace EquipmentId,StageId,Frequency,Power,Cord_x,Cord_y \rbrace
Uijkv={EquipmentId,StageId,Frequency,Power,Cordx,Cordy}
其中,EquipmentId 为用频装备标识;StageId 为用频装备所处平 台标识;Frequency 表示该用频装备的使用频点;Power 表示该用频装备的发射功率;Cordx 表示该用频装备的 x 轴坐标位置,Cordy 表示该用频装备的 y 轴坐标位置。
即元任务表示了,安排用频装备 wi 使用频点 sfi,其发射功率为spk,所在平台的横纵地理坐标为 scv。
设 xijkv 表示元任务 Uijkv 的完成状态,若该元任务得到执行,则 xijkv = 1,否则 xijkv = 0。
则问题的决策变量可以表示为
约束条件
(1)每个用频装备只能且必须执行一次频率、功率、地理坐标的综合调度。
∑
j
=
1
N
(
S
f
i
)
∑
k
=
1
N
(
S
p
i
)
∑
v
=
1
N
(
S
c
i
)
x
i
j
k
v
=
1
,
∀
i
=
1
,
2
,
.
.
.
,
N
w
\sum\limits^{N_{(Sf_i)}}_{j=1}\sum\limits^{N_{(Sp_i)}}_{k=1}\sum\limits^{N_{(Sc_i)}}_{v=1}x_{ijkv}=1,\ \forall i=1, 2,...,N_w
j=1∑N(Sfi)k=1∑N(Spi)v=1∑N(Sci)xijkv=1, ∀i=1,2,...,Nw
(2)用频规划方案中用频装备选 定的频点应在可用频谱资源与该用频装备可选频率信息的交集中。
(3)处在同一作战平台的不同用频装备地理坐标应当保持一致性。
优化目标
(1)干扰冲突最少
用频装备间造成同频干扰的总次数最少
min
D
T
=
∑
i
=
1
N
w
∑
j
=
1
N
w
d
t
i
j
\min DT=\sum\limits^{N_w}_{i=1}\sum\limits^{N_w}_{j=1}dt_{ij}
minDT=i=1∑Nwj=1∑Nwdtij
式中,dtij = 0 表示用频装备 wi 与用频装备 wj 之间不存在同频干扰,dtij = 1 表示存在同频干扰。
(2)需求满足最高
用频需求得到满足的用频装备乘以其优先级系数后的总量最高。用频需求得到满足的定义为其他任意用频装备均不对该用频装备产生同频干扰。
min
S
A
=
∑
i
=
1
N
w
(
s
a
i
×
R
i
)
\min SA = \sum\limits^{N_w}_{i=1}(sa_i×R_i)
minSA=i=1∑Nw(sai×Ri)
式中,sai = 0 表示存在至少一台用频装备对用频装备 wi 产生同频干扰,sai = 1 表示全体用频装备均不对用频装备 wi 产生同频干扰。
(3)邻频风险最低
计算用频规划方案中全体具有电磁干扰可能性的成对用频装备之间邻频干扰风险的总和。邻频干扰风险是两台用频装备间产生邻频干扰的概率及干扰程度的量化表示。
min
A
R
=
∑
i
=
1
N
w
∑
j
=
1
N
w
a
i
j
×
a
r
i
j
\min AR = \sum\limits^{N_w}_{i=1}\sum\limits^{N_w}_{j=1}a_{ij}×ar_{ij}
minAR=i=1∑Nwj=1∑Nwaij×arij
式中,arij = 1/|Frequencyi −Frequencyj|,arij 表示用频装备 wi 与用频装备 wj 间的邻频干扰风险,频点相同时取 arij = 0。
综上,文章所要解决的用频规划问题可概括为,在满足频谱资源约束、同平台位置一致性约束、用频装备唯一性约束等约束条件的前提下,结合联合作战用频需求与用频装备冲突干扰分析,对我方战斗集团中的作战平台与用频装备进行综合调度,为用频装备指配可选频点、明确发射功率以及为作战平台分配地理坐标,制定优化的用频规划方案,以最大化实现干扰冲突最少、需求满足最高、邻频风险最低等多目标的总体最优化。
为了实现“对用频设备时域、空域、频域、能域等多维度的统筹调度”,将用频设备、频点、发射功率、位置坐标信息都包含在信息素矩阵中:建立大小为 Nw × max(NSf) × max(NSp) × max(NSc) 的信息素矩阵 PHM。在初始化阶段,如果用频装备 wi 不使用可选频率信息 Sfij、可选发射功率信息 Spik 与可选地理坐标信息 Sciv,则将 PHM(i, j, k, v) 初始化为 0;否则将 PHM(i, j, k, v) 初始化为 1。
蚂蚁选择下一个城市(指派下一个元任务)的概率:
P
r
(
i
,
j
,
k
,
v
)
=
P
H
M
(
i
,
j
,
k
,
v
)
∑
m
=
1
N
(
S
f
i
)
∑
n
=
1
N
(
S
p
i
)
∑
r
=
1
N
(
S
c
i
)
P
H
M
(
i
,
m
,
n
,
r
)
Pr(i,j,k,v)=\frac {PHM(i,j,k,v)} {\sum\limits^{N_{(Sf_i)}}_{m=1}\sum\limits^{N_{(Sp_i)}}_{n=1}\sum\limits^{N_{(Sc_i)}}_{r=1}PHM(i,m,n,r)}
Pr(i,j,k,v)=m=1∑N(Sfi)n=1∑N(Spi)r=1∑N(Sci)PHM(i,m,n,r)PHM(i,j,k,v)
上式表示用频装备 wi 使用频点 sfj、功率 spk,其所在平台的横纵地理坐标为 scv 的元任务 Uijkv 被成功指配的概率。
信息素的更新方式:
其中,Q 表示信息素矩阵的信息素增量系数;符号 < > 表示四舍五入取整;η 为信息素更新分布系数,其数值随迭代次数 iter 的增长而自适应调整。
为了克服单目标或多目标加权建模不足,采用非支配排序的方法获取 Pareto 最优解,使指挥员可依据自身决策偏好在解集中选择最满意的用频规划方案。
(如果将问题用粒子群算法来实现的话,在粒子编码部分,可以将每个粒子设计为多维矩阵,粒子的每一维不仅包括该用频设备所分配的频点这一一维向量,还包括功率信息、所处地理坐标,形成一个Nw 维的矩阵(假设有 Nw 个用频装备)。)