文章目录

• 一、理论基础
• 二、案例背景
• • 1、问题描述
  • 2、算法流程
  • 3、算法实现
• 三、MATLAB程序实现
• 四、延伸阅读
• 五、参考文献
• 六、完整程序

一、理论基础

标准粒子群算法通过追随个体极值和群体极值完成极值寻优，虽然操作简单，且能够快速收敛，但是随着迭代次数的不断增加，在种群收敛集中的同时，各粒子也越来越相似，可能在局部最优解周围无法跳出。混合粒子群算法摒弃了传统粒子群算法中的通过跟踪极值来更新粒子位置的方法，而是引进了遗传算法中的交叉和变异操作，通过粒子同个体极值和群体极值的交叉以及粒子自身变异的方式来搜寻全局最优解。

二、案例背景

1、问题描述

旅行商问题(traveling salesman problem,TSP)又称为推销员问题、货郎担问题，该问题是最基本的路线问题。该问题寻求单一旅行者由起点出发，通过所有给定的需求点之后，最后再回到起点的最小路径成本。最早的旅行商问题的数学模型是由Dantzig(1959)等人提出的。旅行商问题是车辆路线问题(VRP)的特例，已证明旅行商问题是NP难题。

2、算法流程

基于混合粒子群算法的TSP算法流程如图1所示。

图1 混合粒子群算法流程

其中，种群初始化模块初始化粒子群种群；适应度值计算模块计算粒子群个体的适应度值；更新粒子模块则根据粒子适应度值更新个体最优粒子和群体最优粒子；个体最优交叉即把个体和个体最优粒子进行交叉得到新粒子；群体最优交叉则把个体和群体最优粒子进行交叉得到新粒子；粒子变异是指粒子自身变异得到新粒子。

3、算法实现

• 个体编码
  粒子个体编码采用整数编码的方式，每个粒子表示历经的所有城市，比如当历经的城市数为10，个体编码为[9 4 2 1 3 7 6 10 8 5]，表示城市遍历从9开始，经过4、2、1、3、⋯ \dotsm⋯最终返回起点城市9，从而完成TSP遍历。

• 适应度值
  粒子适应度值表示为遍历路径的长度，计算公式为f i t n e s s ( i ) = ∑ i , j = 1 n p a t h i , j (1) fitness(i)=\sum_{i,j=1}^{n}path_{i,j}\tag{1}fitness(i)=i,j=1∑n​pathi,j​(1)其中，n nn为城市数量；p a t h i , j path_{i,j}pathi,j​为城市i , j i,ji,j间路径长度。

• 交叉操作
  个体通过个体极值和群体极值交叉来更新，交叉方法采用整数交叉法。首先选择两个交叉位置，然后把个体和个体极值或个体与群体极值进行交叉，假定随机选取的交叉位置为3和5，操作方法如下：
  个体[9 4 2 1 3 7 6 10 8 5] 极值[9 2 1 6 3 7 4 10 8 5]
  交叉为新个体[9 4 1 6 3 7 6 10 8 5]
  产生的新个体如果存在重复位置则进行调整，调整方法为用个体中未包括的城市代替重复包括的城市，如下所示：
  [9 4 1 6 3 7 6 10 8 5]调整为[9 4 2 1 3 7 6 10 8 5]
  对得到的新个体采用了保留优秀个体策略，只有当新粒子适应度值好于旧粒子时才更新粒子。

• 变异操作
  变异方法采用个体内部两位互换方法，首先随机选择变异位置p o s 1 pos1pos1和p o s 2 pos2pos2，然后把两个变异位置互换，假设选择的变异位置为2和4，变异操作如下所示：
  [9 4 2 1 3 7 6 10 8 5] 变异为[9 1 2 4 3 7 6 10 8 5]
  对得到的新个体采用了保留优秀个体策略，只有当新粒子适应度值好于旧粒子时才更新粒子。

三、MATLAB程序实现

根据混合粒子群算法原理，在MATLAB中编程实现基于混合粒子群的TSP搜索算法。

• 适应度函数
  适应度函数计算个体适应度值，个体适应度值为路线总长度，代码如下：

function indiFit = fitness(x, cityCoor, cityDist)
%% 该函数用于计算个体适应度值
% x           input      个体
% cityCoor    input      城市坐标
% cityDist    input      城市距离
% indiFit     output     个体适应度值

m = size(x, 1);
n = size(cityCoor, 1);
indiFit = zeros(m, 1);
for i = 1:m
    for j = 1:n-1
        indiFit(i) = indiFit(i)+cityDist(x(i, j), x(i, j+1));
    end
    indiFit(i) = indiFit(i)+cityDist(x(i, 1), x(i, n));
end

• 粒子初始化
  粒子初始化用于初始化粒子，计算粒子适应度值，并根据适应度值确定个体最优粒子和群体最优粒子。程序代码如下：

nMax = 200;                      % 进化次数
indiNumber = 1000;               % 个体数目
individual = zeros(indiNumber,n);
% 初始化粒子位置
for i = 1:indiNumber
    individual(i, :) = randperm(n);    
end
%% 计算种群适应度
indiFit = fitness(individual, cityCoor, cityDist);
[value, index] = min(indiFit);
tourPbest = individual;                          % 当前个体最优
tourGbest = individual(index, :) ;               % 当前全局最优
recordPbest = inf*ones(1, indiNumber);           % 个体最优记录
recordGbest = indiFit(index); % value            % 群体最优记录
xnew1 = individual;

• 交叉操作
  交叉操作把粒子同个体极值和群体极值进行交叉，从而得到较好的个体，交叉操作代码如下：

%% 交叉操作
function a = Recombin(a, b)
%% 交叉
% 输入：
% a为待交叉的个体，b为交叉的个体极值
% 输出：
% a为交叉后得到的新个体
n = length(a);
c1 = unidrnd(n-1);         % 产生交叉位
c2 = unidrnd(n-1);         % 产生交叉位
while c1 == c2
    c1 = round(rand*(n-2))+1;
    c2 = round(rand*(n-2))+1;
end
chb1 = min(c1, c2);
chb2 = max(c1, c2);
cros = b(chb1:chb2);        % 交叉区域元素
ncros = length(cros);      % 交叉元素个数
% 删除与交叉区域相同元素
for j=1:ncros
    for k=1:n
        if a(k) == cros(j)
            a(k) = 0;
            for t = 1:n-k
                temp = a(k+t-1);
                a(k+t-1) = a(k+t);
                a(k+t)=temp;
            end
        end
    end
end
%插入交叉区域
a(n-ncros+1:n) = cros;        
% a0 = a; b0 = b;
% for i = chb1:chb2
%     a1 = a; b1 = b;
%     a(i) = b0(i);
%     b(i) = a0(i);
%     x = find(a == a(i));
%     y = find(b == b(i));
%     i1 = x(x ~= i);
%     i2 = y(y ~= i);
%     if ~isempty(i1)
%         a(i1) = a1(i);
%     end
%     if ~isempty(i2)
%         b(i2) = b1(i);
%     end
% end

• 变异操作
  变异操作对自身进行变异，从而得到更好的个体。变异操作代码如下：

function x = Mutate(x)
% 输入：x        个体
% 输出：x        变异后的个体
n = length(x);
while true
    c1 = round(rand*(n-1))+1;       % 产生变异位
    c2 = round(rand*(n-1))+1;       % 产生变异位
    if c1 ~= c2
        break;
    end
end
temp = x(c1);
x(c1) = x(c2);
x(c2) = temp;

• 仿真结果
  采用混合粒子群算法规划TSP路径，各城市的初始位置如图2所示。
  

图2 城市初始位置

混合粒子群算法的进化次数为100，种群规模为100，算法进化过程中最优粒子适应度值变化和规划出的最优路径如图3和图4所示。

图3 适应度值变化

图4 规划出的最优路径

四、延伸阅读

使用混合粒子群算法规划其他城市模型的最短路径。
城市分布图如图5所示，混合粒子群算法的进化次数为200，种群规模为1000。采用混合粒子群算法规划的路径如图6所示。

图5 城市分布图

图6 规划路径

从图2~图6可以看出，基于混合粒子群算法的TSP算法可以解决规模较小的旅行商问题，对于规模较大的旅行商问题，混合粒子群算法也可以得到较优路径。

五、参考文献

[1] Kennedy J . Particle swarm optimization[J]. Proc. of 1995 IEEE Int. Conf. Neural Networks, (Perth, Australia), Nov. 27-Dec. 2011, 4(8):1942-1948 vol.4.
[2] 王伟. 混合粒子群算法及其优化效率评价[J]. 中国水运:学术版, 2007, 7(006):102-103.
[3] 郁磊等. MATLAB智能算法30个案例分析(第2版)[M].北京航空航天大学出版社.2015年.

代码添加QQ3381151092