几种粒子群算法的matlab实现

%%%%%%%%%%%%小世界粒子群算法 clear all;close all;clc; Nm=1e3;%循环次数 num=20;%粒子总数 numd=10;%粒子维数 MAXR=100;MINR=-100;%定义域最大值的绝对值 Vmax=2*MAXR/5;%(MAXR-MINR)/2; omiga_max=0.87;omiga_min=0.525;%omiga_max=0.9;omiga_min=0.8; c1=1.1;c2=1.1;%c1=0.45;c2=0.45;c1=1.49618;c2=1.49618 K=5;p=0.2; fid=fopen('小世界PSO.txt','at+'); rand('state',sum(100*clock)); for NN=1:30 TimeStart=cputime;%********************初始化 Z=MAXR-2*MAXR*rand(num,numd);V=Vmax-2*Vmax*rand(num,numd);%初始化粒子位置和速度 %%%%初始化邻域 %---------这段代码可以用在聚类问题上---------------% % for n0=1:num % for n1=n0+1:num % dis(n0,n1)=sum((Z(n0,:)-Z(n1,:)).^2);%dis表示粒子的欧氏距离矩阵 % end % end % for n0=1:num % dis(n0,n0)=0; % end % for n0=2:num % for n1=1:n0-1 % dis(n0,n1)=dis(n1,n0); % end % end % for n0=1:num % [b,ind]=sort(dis(n0,:)); %对欧氏距离矩阵排序并求出相应下标 % neighborhood(n0,:)=ind(1:K);%neighborhood记录粒子n0的邻域粒子的下标 % end %%%%初始化邻域,下标相邻的K个粒子为一个邻域，首尾相接 for n0=1:num for n1=1:num if abs(n1-n0)<=(K-1)/2 || abs(n1-n0-num)<=(K-1)/2 || abs(n1-n0+num)<=(K-1)/2 neighborhood(n0,n1)=1; else neighborhood(n0,n1)=0; end end end for n0=1:num D=fun(Z(n0,:));%计算适应值 Zp(n0,:)=Z(n0,:);%初始化粒子自身最佳位置 Dp(n0)=D; %初始化粒子自身最佳适应值 end for n0=1:num temp=find(neighborhood(n0,:)==1); [Dng(n0),index]=min(Dp(temp));%初始化粒子邻域最佳适应值 ZNg(n0,:)=Zp(temp(index),:);%初始化粒子邻域最佳位置 end for n0=1:num%%%%更新速度和位置 r1=rand(1,numd);r2=rand(1,numd);r3=rand(1,numd); V(n0,:)=omiga_max*V(n0,:)+c1*r1.*(Zp(n0,:)-Z(n0,:))+c2*r2.*(ZNg(n0,:)-Z(n0,:)); V(find(V<-Vmax))=-Vmax; V(find(V>Vmax))=Vmax; Z(n0,:)=Z(n0,:)+V(n0,:); Z(find(Z<MINR))=MINR; Z(find(Z>MAXR))=MAXR; end for n0=1:num%%%%更新粒子的邻域（以概率p断开原有连接并建立新连接） temp=find(neighborhood(n0,:)==1,sum(neighborhood(n0,:)),'first'); for n1=1:sum(neighborhood(n0,:)) if(n0~=temp(n1) && rand<=p) neighborhood(n0,temp(n1))=0; neighborhood(temp(n1),n0)=0; t=find(neighborhood(n0,:)==0,20-sum(neighborhood(n0,:)),'first'); tt=t(ceil(length(t)*rand)); neighborhood(n0,tt)=1;neighborhood(tt,n0)=1; end end end Dg=min(Dng);%Zg=Zng(find(Dng==Dg,1,'first'),:);%全局最佳适应值和最佳位置 BestFitness(1,1)=Dg; for nn=2:Nm%******************************************************循环 for n0=1:num D=fun(Z(n0,:));%计算适应值 if Dp(n0)>=D Dp(n0)=D; %保存粒子自身最佳适应值 Zp(n0,:)=Z(n0,:);%更新粒子自身最佳位置 end end for n0=1:num temp=find(neighborhood(n0,:)==1); [t,index]=min(Dp(temp)); if Dng(n0)>=t Dng(n0)=t;%更新粒子邻域最佳适应值 ZNg(n0,:)=Zp(temp(index),:);%更新粒子邻域最佳位置 end end for n0=1:num%%%%更新速度和位置 r1=rand(1,numd);r2=rand(1,numd);r3=rand(1,numd); V(n0,:)=(omiga_max-(omiga_max-omiga_min)*nn/Nm)*V(n0,:)+c1*r1.*(Zp(n0,:)-Z(n0,:))+c2*r2.*(ZNg(n0,:)-Z(n0,:)); V(find(V<-Vmax))=-Vmax; V(find(V>Vmax))=Vmax; Z(n0,:)=Z(n0,:)+V(n0,:); Z(find(Z<MINR))=MINR; Z(find(Z>MAXR))=MAXR; end for n0=1:num %%%%更新粒子的邻域（以概率p断开原有连接并建立新连接） temp=find(neighborhood(n0,:)==1,sum(neighborhood(n0,:)),'first'); for n1=1:sum(neighborhood(n0,:)) if(n0~=temp(n1) && rand<=p) neighborhood(n0,temp(n1))=0; neighborhood(temp(n1),n0)=0; t=find(neighborhood(n0,:)==0,20-sum(neighborhood(n0,:)),'first'); tt=t(ceil(length(t)*rand)); neighborhood(n0,tt)=1;neighborhood(tt,n0)=1; end end end Dg=min(Dng);%Zg=Zng(find(Dng==Dg,1,'first'),:);%全局最佳适应值和最佳位置 Kc(nn)=sum(sum(neighborhood))/num-1; Cp(nn)=(3*(Kc(nn)-2)*(1-p)^3)/(4*(Kc(nn)-1));%粒子群的聚类系数 Lp(nn)=2*num*(atanh(sqrt((num*Kc(nn)*p/2)/(num*Kc(nn)*p/2+2)))/2*sqrt((num*Kc(nn)*p/2)^2+2*(num*Kc(nn)*p/2)))/Kc(nn);%平均路径长度 BestFitness(1,nn)=Dg; end e=cputime-TimeStart; % for k=2:Nm % temp(1,k-1)=BestFitness(1,k-1)-BestFitness(1,k); % end % subplot(2,1,1),plot(1:Nm,BestFitness); % subplot(2,1,2),plot(temp); [f name]=fun(zeros(1,1)); if NN==1 fprintf(fid,'Function Name:%s K:%d p:%2.2f c1:%2.2f c2:%2.2f 惯性权重:%2.2f-%2.2f 维数:%d 粒子数:%d 搜索代数:%d 搜索范围:[%d,%d]\n%e %f\n',name,K,p,c1,c2,omiga_min,omiga_max,numd,num,Nm,MINR,MAXR,Dg,e); else fprintf(fid,'%e %f\n',Dg,e); end Dglist(NN)=Dg;TimeList(NN)=e; end aveDg=sum(Dglist)/NN stdDg=std(Dglist); aveTime=sum(TimeList)/NN fprintf(fid,'试验次数:%d 平均最优值:%e+-%e 平均耗时:%f\n**************************\n',NN,aveDg,stdDg,aveTime); fclose(fid);