粒子群优化算法求解旅行商问题

1粒子群优化算法求解旅行商问题深圳大学信息工程学院黄彩玲2005年6月16日SZU－TIDSPLab2粒子群优化算法求解旅行商问题参照：粒子群优化算法求解旅行商问题黄岚等吉林大学学报（理学版）2003年10月SZU－TIDSPLab3五个定义1设n个节点的TSP问题的解序列为s=(ai),I=1…n.定义交换子SO(i1,i2)为交换解S中的点ai1和ai2,则S’=S+SO(i1,i2)为解S经算子SO(i1,i2)操作后的新解。这里的＋的含义是执行交换操作。2一个或多个交换子的有序队列就是交换序，记作SS，SS=(SO1,SO2,…SON),SO1,SO2等是交换子，之间的顺序是有意义的。作用于一个TSP问题是意味着所有的交换子依次作用于该解上。3不同的交换序作用于同一解上可能产生相同的新解，所有有相同效果的交换序的集合称为交换序的等价集。4若干个交换序可以合并成一个新的交换序，定义＋为两个交换序的合并算子。5在交换序等价集中，拥有最少交换子的交换序称为该等价集的基本交换序。SZU－TIDSPLab4算式Vid’=Vid+alpha(Pid-Xid)+beta(Pgd-Xid)（式1）Xid’=Xid+Vid（式2）－alpha和beta为（0，1）之间的随机数。－(Pid-Xid)和(Pgd-Xid)是基本交换序，表示Xid经过交换得到Pid和Pgd。－alpha(Pid-Xid)表示基本交换序(Pid-Xid)中的所有交换子以概率alpha保留。beta(Pgd-Xid)同理。SZU－TIDSPLab5算法步骤1初始化粒子群，给每个粒子一个初始解和随机的交换序。2如果满足结束条件，转步骤5。3根据粒子当前位置X计算下一新解X’:1）计算（Pid-Xid）交换序。2）计算(Pgd-Xid)。3）计算式1，并将Vid’转换为基本交换序。4）计算式2，更新Xid。5）如果找到一个更好得解，更新Pid。4如果整个群体找到一个更好的解，更新Pgd，转步骤2。5显示结果。SZU－TIDSPLab6程序分析主要数据结构：种群大小(PopSize)空间维数(NDim)最大迭代次数(MaxIter)城市之间的距离(nodeDistance)各粒子当前适应度值(fvalue)更新前各粒子适应度值(fpbest)各粒子位置(population)各粒子速度(velocity)各粒子的最佳位置(pbest)全局最佳粒子位置(gbest)全局最佳粒子序号(index)得到相近适应度值的迭代次数(samecounter)计算samecounter的临时变量(oldbestval)SZU－TIDSPLab7主要函数BeginWith1：使每个路径都从节点1开始排起。（这个命名不好）wholeDistance：计算路径即一个循环后的距离。（改为caculateWholeDistance）ExchangeIndex：计算（Pid-Xid）等交换序。（大小写要统一,CaculateWholeDistance,ChangeIndex）HoldByAlpha：计算以概率保留交换序。changeIndex：进行交换操作。SZU－TIDSPLab8初始化各主要数据（设求14点的TSP）flag=0;%停止程序标志oldbestval=0;%记录旧的适应度值samecounter=0;%记录得到相同适应度值的迭代次数iteration=0;%迭代次数MaxIter=2000;%最大迭代次数PopSize=200;%种群大小alpha=0.8;%概率beta=0.4;w=0.4;NDim=14;population=ones(NDim,PopSize);%初始化各粒子，即产生路径种群fori=1:PopSizepopulation(:,i)=randperm(NDim)';endpopulation=BeginWith1(population);%以1为起点重排每个路径SZU－TIDSPLab9初始化各主要数据node=[16.4796.10;16.4794.44;20.0992.54;22.3993.37;25.2397.24;...22.0096.05;20.4797.02;17.2096.29;16.3097.38;14.0598.12;...16.5397.38;21.5295.59;19.4197.13;20.0994.55];％节点坐标nodeDistance=zeros(NDim,NDim);%计算每个城市之间的距离fori=1:NDimforj=1:NDimnodeDistance(i,j)=sqrt((node(i,1)-node(j,1))^2+(node(i,2)-node(j,2))^2);endEndfori=1:PopSize%计算各路径的距离fvalue(i)=wholeDistance(population(:,i)',nodeDistance);Endvelocity=zeros(NDim,PopSize);%产生交换序fori=1:PopSizevelocity(:,i)=round(rand(1,NDim)'*NDim);endSZU－TIDSPLab10初始化各主要数据pbest=population;%记录各粒子的个体极值点位置,即个体找到的最短路径fpbest=fvalue;%记录最佳适应度值,即个体找到的最短路径的长度[fbestval,index]=min(fvalue);%找出全局极值和相应的序号SZU－TIDSPLab11主程序while(flag==0)&(iterationMaxIter)%寻找最优主程序iteration=iteration+1;%迭代次数递增fori=1:PopSize%更新全局极值点位置,这里指路径gbest(:,i)=population(:,index);endpid_xid=ExchangeIndex(population,pbest);%求pid-xid,pgd-xid交换序pid_xid=HoldByAlpha(pid_xid,alpha);%以概率alpha保留交换序pgd_xid=ExchangeIndex(population,gbest);pgd_xid=HoldByAlpha(pgd_xid,beta);velocity=HoldByAlpha(velocity,w);%以概率w保留交换序population=changeIndex(population,velocity);%根据交换序进行路径交换population=changeIndex(population,pid_xid);population=changeIndex(population,pgd_xid);SZU－TIDSPLab12主程序fori=1:PopSize%更新各路径距离fvalue(i)=wholeDistance(population(:,i)',nodeDistance);endchangeColumns=fvaluefpbest;%更新后的距离优于更新前的，记录序号pbest(:,find(changeColumns))=population(:,find(changeColumns));%更新个体最佳路径fpbest=fpbest.*(~changeColumns)+fvalue.*changeColumns;%更新个体最佳路径距离[fbestval,index]=min(fvalue);%更新全局最佳路径,记录相应的序号iffbestval==oldbestval%比较更新前和更新后的适应度值;samecounter=samecounter+1;%相等时记录加一;elseoldbestval=fbestval;%不相等时更新适应度值，并记录清零;samecounter=0;endSZU－TIDSPLab13主程序ifsamecounter=20%多次迭代的适应度值相近时程序停止flag=1;endBest(iteration)=fbestval;%输出及描出找到的全局极值endSZU－TIDSPLab14BeginWith1()函数functionpopulation=BeginWith1(population)[xy]=size(population);[NO1,index]=min(population',[],2);%找最小值1fori=1:ypop=population(:,i);temp1=pop([1:index(i)-1]);temp2=pop([index(i):x]);population(:,i)=[temp2'temp1']';endSZU－TIDSPLab15changeIndex（）函数functionpopulation=changeIndex(population,velocity)[xy]=size(population);fori=1:ya=velocity(:,i);%取出其中一组交换序pop=population(:,i);%取出对应的粒子forj=1:x%取出其中一个交换算子作交换ifa(j)~=0pop1=pop(j);%temp';%temp(a(j));pop(j)=pop(a(j));pop(a(j))=pop1;endendpopulation(:,i)=pop;endSZU－TIDSPLab16ExchangeIndex（）函数functionpgid_xid=ExchangeIndex(population,pgbest);[xy]=size(population);pgid_xid=zeros(x,y);fori=1:ypop=pgbest(:,i);%从pgbest取出一个顺序pop1=population(:,i);%从粒子群中取出对应的顺序forj=1:x%从pgbest的顺序中取出一个序号NoFrompgbest=pop(j);fork=1:x%从对应的粒子顺序中取出一个序号NoFromPopulation=pop1(k);if(NoFrompgbest==NoFromPopulation)&&(j~=k)%两序号同且不在同一位置pgid_xid(j,i)=k;%交换算子pop1(k)=pop1(j);pop1(j)=NoFromPopulation;endendendendSZU－TIDSPLab17HoldByAlpha（）函数functionvelocity=HoldByAlpha(velocity,w)[x,y]=size(velocity);fori=1:xforj=1:yifrandwvelocity(i,j)=0;endendendSZU－TIDSPLab18wholeDistance（）函数functiondist=wholeDistance(path,nodeDistance)L=length(path);%path为一个循环的节点顺序dist=0;fori=1:L-1dist=dist+nodeDistance(path(i),path(i+1));enddist=dist+nodeDistance(path(1),path(L));%加上首尾节点的距离SZU－TIDSPLab19运行结果节点数alphabetaw粒子数理想达优值opt程序运行次数结果不变最大迭代次数限制达优次数平均值ave平均值平均运行时间err140.80.40.420030.87855020531.868531.60943.135s2.37%Err=(ave-opt)/optSZU－TIDSPLab20运行结果运行五十次，每次得的最短距离如图所示0510152025303540455030.53131.53232.53333.53434.535SZU－TIDSPLab21改进微粒群优化算法求解旅行商问题参照：改进微粒群优化算法求解旅行商问题肖健梅等计算机工程与应用2004.35该论文在上述论文基础上作出改