基于极光优化算法（Polar Lights Optimization, PLO）的多无人机协同三维路径规划（提供MATLAB代码）

一、极光优化算法介绍

极光优化算法（Polar Lights Optimization, PLO）是2024年提出的一种新型的元启发式优化算法，它从极光这一自然现象中汲取灵感。极光是由太阳风中的带电粒子在地球磁场的作用下，与地球大气层中的气体分子碰撞而产生的光显示。PLO算法通过模拟这些带电粒子的运动轨迹和动力学过程，提出了一种有效的优化策略，旨在解决复杂的优化问题。

1.PLO算法的核心思想

PLO算法的核心思想是模拟带电粒子在地球磁场中的运动，包括它们的旋转运动（gyration motion）、在极光椭圆区域内的行走（aurora oval walk），以及粒子间的碰撞（particle collision）。这些运动策略不仅模拟了自然界中粒子的行为，而且被设计用来平衡算法的全局探索能力和局部开发能力。

2.PLO算法的主要组成部分

1. 旋转运动（Gyration Motion）：

   • 模拟带电粒子在地球磁场作用下沿磁力线的旋转运动。这种运动有助于算法在局部区域内进行细致的搜索，提高解的精度。
     

2. 极光椭圆行走（Aurora Oval Walk）：

   • 模拟带电粒子在极光椭圆区域内的自由移动，这种策略有助于算法进行全局探索，寻找更广泛的解空间。
     

3. 粒子碰撞策略（Particle Collision Strategy）：

   • 模拟带电粒子间的相互碰撞，这种策略有助于算法跳出局部最优解，增强解的多样性。
     

3.PLO算法的实现步骤

1. 初始化：

   • 生成一个初始粒子群，每个粒子代表一个潜在的解。

2. 适应度评估：

   • 对每个粒子计算其适应度值，即评估它们在优化问题中的表现。

3. 更新最优解：

   • 根据适应度值，确定当前找到的最优解。

4. 迭代过程：

   • 重复以下步骤直到满足停止条件（如达到最大迭代次数或达到预定的适应度阈值）：
     a. 更新粒子速度和位置：根据旋转运动和极光椭圆行走策略更新粒子的速度和位置。
     b. 粒子碰撞：以一定的概率执行粒子碰撞策略，进一步更新粒子的位置。
     c. 适应度评估：对更新后的粒子群再次进行适应度评估。
     d. 更新最优解：如果找到更好的解，则更新全局最优解。

5. 输出最优解：

   • 迭代结束后，输出找到的最优解或最优解集合。

参考文献：
[1]Yuan, Chong, Dong Zhao, Ali Asghar Heidari, Lei Liu, Yi Chen and Huiling Chen. “Polar lights optimizer: Algorithm and applications in image segmentation and feature selection.” Neurocomputing 607 (2024): 128427.

二、无人机（UAV）三维路径规划

单个无人机三维路径规划数学模型参考如下文献：

Phung M D , Ha Q P . Safety-enhanced UAV Path Planning with Spherical Vector-based Particle Swarm Optimization[J]. arXiv e-prints, 2021.

每个无人机的目标函数由路径长度成本，安全性与可行性成本、飞行高度成本和路径平滑成本共同组成：

2.1路径长度成本

路径长度成本由相邻两个节点之间的欧氏距离和构成，其计算公式如下：

2.2路径安全性与可行性成本

路径安全性与可行性成本通过下式计算：

2.3路径飞行高度成本

飞行高度成本通过如下公式计算所得：

2.4路径平滑成本

投影向量通过如下公式计算：

转弯角度的计算公式为：

爬坡角度的计算公式为：

平滑成本的计算公式为：

2.5总成本（目标函数）

总成本由最优路径成本，安全性与可行性成本、飞行高度成本和路径平滑成本的线性加权所得。其中，b为加权系数。

三、实验结果

在三维无人机路径规划中，无人机的路径由起点，终点以及起始点间的点共同连接而成。因此，自变量为无人机起始点间的各点坐标，每个无人机的目标函数为总成本（公式9）。本文研究3个无人机协同路径规划，总的目标函数为3个无人机的总成本之和。

Xmin=[Xmin0,Xmin1,Xmin2];
Xmax=[Xmax0,Xmax1,Xmax2];
dim=dim0+dim1+dim2;
fobj=@(x)GetFun(x,fobj0,fobj1,fobj2);%总的目标函数
pop=50;
maxgen=1500;

[fMin ,bestX,Convergence_curve]=eco(pop,maxgen,Xmin,Xmax,dim,fobj);%Trajectories,fitness_history, population_history
% save bestX bestX
BestPosition1 = SphericalToCart(bestX(1:dim/3),model);% 第一个无人机得到的路径坐标位置
BestPosition2 = SphericalToCart(bestX(1+dim/3:2*dim/3),model1);% 第二个无人机得到的路径坐标位置
BestPosition3 = SphericalToCart(bestX(1+2*dim/3:end),model2);% 第三个无人机得到的路径坐标位置

gca1=figure(1);
gca2=figure(2);
gca3=figure(3);
PlotSolution(BestPosition1,model,gca1,gca2,gca3);% 画第一个无人机
PlotSolution1(BestPosition2,model1,gca1,gca2,gca3);% 画第二个无人机
PlotSolution2(BestPosition3,model2,gca1,gca2,gca3);% 画第三个无人机

figure
plot(Convergence_curve,'LineWidth',2)
xlabel('Iteration');
ylabel('Best Cost');
grid on;