%%=========================================================================
%% ACATSP.m
%% Ant Colony Algorithm for Traveling Salesman Problem
%% ChengAihua,PLA Information Engineering University,ZhengZhou,China
%% Email:aihuacheng@gmail.com
%% All rights reserved
%%-------------------------------------------------------------------------
%% 主要符号说明
%% C n个城市的坐标,n×2的矩阵
%% NC_max 最大迭代次数
%% m 蚂蚁个数
%% Alpha 表征信息素重要程度的参数
%% Beta 表征启发式因子重要程度的参数
%% Rho 信息素蒸发系数
%% Q 信息素增加强度系数
%% R_best 各代最佳路线
%% L_best 各代最佳路线的长度
%% L_ave 各代路线的平均长度
%%=========================================================================
%%第一步:变量初始化
point=[-67.15 18.42
-67.13 18.21
-67.19 17.95
-67.02 18.08
-66.88 18.30
-66.78 18.48
-66.56 18.44
-66.74 18.18
-66.93 17.97
-66.61 17.99
-66.32 18.99];
n=11;
m=6;
Alpha=1;
Beta=2;
vol=0.2;
Q=10;
for i=1:n
for j=1:n
D(i,j)=sqrt(sum((point(i,:)-point(j,:)).^2));
end
end
for i=1:n
D(i,i)=eps;
end
Heu_F=1./D;%Eta为启发因子,这里设为距离的倒数
Tau=ones(n,n);%Tau为信息素矩阵
Table=zeros(m,n);%存储并记录路径的生成
iter=1;%迭代计数器
iter_max=100;
R_best=zeros(iter_max,n);%各代最佳路线
L_best=inf.*ones(iter_max,1);%各代最佳路线的长度
L_ave=zeros(iter_max,1);%各代路线的平均长度
while iter<=iter_max%停止条件之一:达到最大迭代次数
%%第二步:将m只蚂蚁放到n个城市上
start=zeros(m,1);
for i=1:m
temp=randperm(n);
start(i)=temp(1);
end
Table(:,1)=start;
%
citys_index=1:n;
%%第三步:m只蚂蚁按概率函数选择下一座城市,完成各自的周游
for i=1:m
for j=2:n
tabu=Table(i,1:(j-1));
allow_index=~ismember(citys_index,tabu);%已访问的城市
allow=citys_index(allow_index);%待访问的城市
P=allow;%待访问城市的选择概率分布
%
for k=1:length(allow)
P(k)=(Tau(tabu(end),allow(k))^Alpha)*(Heu_F(tabu(end),allow(k))^Beta);%(信息素^信息素系数)*(启发因子^启发因子系数)
end
P=P/sum(P);
%按概率原则选取下一个城市
Pc=cumsum(P);
target_index=find(Pc>=rand);
target=allow(target_index(1));
Table(i,j)=target;
end
end
if iter>=2
Table(1,:)=R_best(iter-1,:);
end
%%第四步:记录本次迭代最佳路线
L=zeros(m,1);
for i=1:m
R=Table(i,:);
for j=1:(n-1)
L(i)=L(i)+D(R(j),R(j+1));
end
L(i)=L(i)+D(R(1),R(n));
end
L_best(iter)=min(L);
pos=find(L==L_best(iter));
R_best(iter,:)=Table(pos(1),:);
L_ave(iter)=mean(L);
iter=iter+1
%%第五步:更新信息素
Delta_Tau=zeros(n,n);
for i=1:m
for j=1:(n-1)
Delta_Tau(Table(i,j),Table(i,j+1))=Delta_Tau(Table(i,j),Table(i,j+1))+Q/L(i);
end
Delta_Tau(Table(i,n),Table(i,1))=Delta_Tau(Table(i,n),Table(i,1))+Q/L(i);
end
Tau=(1-vol).*Tau+Delta_Tau;
%%第六步:禁忌表清零
Table=zeros(m,n);
end
%%第七步:输出结果
Pos=find(L_best==min(L_best));
Shortest_Route=R_best(Pos(1),:)
Shortest_Length=L_best(Pos(1))