clear all

close all

clc

%%%%%%%%%%初始化参数%%%%%%%%%

NP=20;

FoodNumber=NP/2;

Limit=100;

maxCycle=2500;

%%%%%定义全局变量%%%%%%%%%

global NC;  %总的城市数

global d;   %各城市之间的距离

NC=14;

x(1)=16.47;x(2)=16.47;x(3)=20.09;x(4)=22.39;x(5)=25.23;x(6)=22.00;x(7)=20.47;

x(8)=17.20;x(9)=16.30;x(10)=14.05;x(11)=16.53;x(12)=21.52;x(13)=19.41;x(14)=20.09;

y(1)=96.10;y(2)=94.44;y(3)=92.54;y(4)=93.37;y(5)=97.24;y(6)=96.05;y(7)=97.02;y(8)=96.29;

y(9)=97.38;y(10)=98.12;y(11)=97.38;y(12)=95.59;y(13)=97.13;y(14)=94.55;

%计算各城市之间的距离

for m=1:14

for n=1:14

d(m,n)=sqrt((x(m)-x(n))^2+(y(m)-y(n))^2);

end

end

%画初始路径图

figure(1);

for m=1:NC

scatter(x(m),y(m),50,'r')

hold on

end

plot([x(1),x(2),x(3),x(4),x(5),x(6),x(7),x(8),x(9),x(10),x(11),x(12),x(13),x(14),x(1)],...

[y(1),y(2),y(3),y(4),y(5),y(6),y(7),y(8),y(9),y(10),y(11),y(12),y(13),y(14),y(1)]);

%title('初始路径')；

xlable('城市x坐标');

ylable('城市y坐标');

runtime=10;  %通过修改runtime的值，改变程序的运行次数，用以算法的健壮性

GlobalMins=zeros(1,runtime);

for r=1:runtime

%初始化蜂群

for i=1:FoodNumber

Foods(i,:)=initial();

end

%计算适应度函数值

for i=1:FoodNumber

route=Foods(i,:);

Fitness(i)=calculateFitness(route);

end

%初始化搜索次数，用于Limit比较

trial=zeros(1,FoodNumber);

%找出适应度函数值得最小值

BestInd=find(Fitness==min(Fitness));

BestInd=BestInd(end);  %避免有两个相同位置，只取其一

GlobalMin=Fitness(BestInd);

GlobalParams=Foods(BestInd,:);

%迭代开始

iter=1;

j=1;  %用以初始化结果

while ((iter<=maxCycle))

%%%采蜜蜂模式%%

%随机交换3个城市的顺序，改变此时的路线

route_next=Foods(i,:);

l=round(2+12*rand());

m=round(2+12*rand());

n=round(2+12*rand());

temp=route_next(1);

route_next(l)=route_next(m);

route_next(m)=route_next(n);

route_next(n)=temp;

%计算新蜜源的适应度函数值

FitnessSol=calculateFitness(route_next);

%使用贪婪原则，寻找最优蜜源

if (FitnessSol

Foods(i,:)=route_next;

Fitness(i)=FitnessSol;

trial(i)=0;

else

trial(i)=trial(i)+1;  %超过设定的limit次，则该蜂成为侦查蜂

end

end

%%%%根据适应度值计算采蜜蜂被跟随的概率%%%%%%

pro=(0.9*Fitness./max(Fitness))+0.1;

%%%观察蜂%%%%%%%

i=1;          %要跟随的采蜜蜂

t=0;           %标记观察蜂

while (t

if(rand

t=t+1;

%随机交换3个城市的顺序，改变此时的路线

route_next=Foods(i,:);

l=round(2+12*rand());

m=round(2+12*rand());

n=round(2+12*rand());

temp=route_next(1);

route_next(l)=route_next(m);

route_next(m)=route_next(n);

route_next(n)=temp;

%计算新蜜源的适应度函数值

FitnessSol=calculateFitness(route_next);

%使用贪婪原则，寻找最优蜜源

if (FitnessSol

Foods(i,:)=route_next;

Fitness(i)=FitnessSol;

trial(i)=0;

else

trial(i)=trial(i)+1;  %超过设定的limit次，则该蜂成为侦查蜂

end

end

i=i+1;  %要跟随的下一个采蜜蜂

if (i==(FoodNumber)+1)

i=1;

end

end

%记录此时更好的解

ind=find(Fitness==min(Fitness));

ind=ind(end);

if (Fitness(ind)

GlobalMin=Fitness(ind);

GlobalParams=Foods(ind,:);

end

%%%侦查蜂模式%%%

ind=find(trial==max(trial));

ind=ind(end);

if (trial(ind)>Limit)  %若搜索次数超过极限值，则进行随机搜索产生新的解

Bas(ind)=0;

for m=1:1:NC

route_new(m)=m;

end

route_new(NC+1)=1;

FitnessSol=calculateFitness(route_new);

Foods(ind,:)=route_new;

Fitness(ind)=FitnessSol;

end

%%%%%%建立一个次数和最优解的矩阵，以便于画图%%%%%

Cishu(j)=iter;  %迭代次数行向量

Zuiyou(j)=GlobalMin;  %每次迭代得到的最优解

j=j+1;

iter=iter+1;

end

while(iter<=maxCycle)

GlobalMins(r)=GlobalMin;  %程序运行完一次，记录这次的最优路径的长度

disp(['第',num2str(r),'次运行得到的最优路径是:',num2str(GlobalParams),'此条路径的长度是:',num2str(GlobalMins(r))]);

end

%%图曲线图%%

figure(2);

plot(Cishu,Zuiyou,'b');

%title('优化曲线')；

xlable('迭代次数');

ylable('路径长度');

figure(3);

for m=1:NC

scatter(x(m),y(m),50,'r')

hold on

end

plot([x(GlobalParams(1)),x(GlobalParams(2)),x(GlobalParams(3)),x(GlobalParams(4)),x(GlobalParams(5)),x(GlobalParams(6)),

x(GlobalParams(7)),x(GlobalParams(8)),x(GlobalParams(9)),x(GlobalParams(10)),x(GlobalParams(11)),x(GlobalParams(12)),

x(GlobalParams(13)),x(GlobalParams(14)),x(GlobalParams(1)),y(GlobalParams(1)),y(GlobalParams(2)),y(GlobalParams(3)),

y(GlobalParams(4)),y(GlobalParams(5)),y(GlobalParams(6)),y(GlobalParams(7)),y(GlobalParams(8)),y(GlobalParams(9)),

y(GlobalParams(10)),y(GlobalParams(11)),y(GlobalParams(12)),y(GlobalParams(13)),y(GlobalParams(14)),y(GlobalParams(1))]);

%titlt('优化路径')；

xlable('城市x坐标');

ylable('城市y坐标');