网上已经有很多matlab实现的蚁群算法的源代码,也有详细的注释,但是有些注释有误,在这里更正,如有错误,希望各位批评指正。
以下是解放军信息工程大学一个老师编的matlab程序,请尊重原作者劳动,引用时请注明出处。注释有一部分来源于网络,对于明显错误的地方我已经更正
function [R_best,L_best,L_ave,Shortest_Route,Shortest_Length]=ACATSP(C,NC_max,m,Alpha,Beta,Rho,Q)
%% 主要符号说明
%% C n个城市的坐标,n×2的矩阵
%% NC_max 最大迭代次数
%% m 蚂蚁个数
%% Alpha 表征信息素重要程度的参数,信息重要程度的参数如果不是固定的,而是动态变化的会怎么样
%% Beta 表征启发式因子重要程度的参数,同理启发因子也可能会动态变化
%% Rho 信息素蒸发系数
%% Q 信息素增加强度系数
%% R_best 各代最佳路线
%% L_best 各代最佳路线的长度
%%第一步:变量初始化
n=size(C,1);%n表示问题的规模(城市个数)
D=zeros(n,n);%D表示完全图的赋权邻接矩阵,D(i,j)表示城市i与城市j之间的距离
%生成n个城市的邻接矩阵
for i=1:n
for j=1:n
if i~=j
D(i,j)=((C(i,1)-C(j,1))^2+(C(i,2)-C(j,2))^2)^0.5;%欧氏距离
else
D(i,j)=eps; %i=j时不计算,应该为0,但后面的启发因子要取倒数,用eps(浮点相对精度)表示
end
D(j,i)=D(i,j); %对称矩阵
end
end
Eta=1./D; %Eta为启发因子,这里设为距离的倒数
Tau=ones(n,n); %Tau为信息素矩阵
Tabu=zeros(m,n); %记录每一次迭代中每一只蚂蚁的旅游路径
NC=1; %迭代计数器,记录迭代次数
R_best=zeros(NC_max,n); %各代最佳路线,记录每一次迭代的最佳路径
L_best=inf.*ones(NC_max,1); %各代最佳路线的长度
L_ave=zeros(NC_max,1); %各代路线的平均长度
%迭代
%%
%学习本算法最弄明白以下问题,思考以下问题,把这几个问题弄明白了,蚁群算法你就可以灵活运用
%1:各城市的邻接矩阵,其值为个城市之间的欧氏距离,对称矩阵(现实应用中也可以是不对称的)
%2:启发因子矩阵,其值为各城市距离邻接矩阵的倒数,在整个算法流程中,启发因子矩阵是如何发挥作用的,迭代过程中可不可以不断发展变化?,还是一层不变??%3:信息素矩阵,信息素为每一个蚂蚁在履行过程中留下的信息,所以一开始初始化时所有的值相同,信息素在迭代过程中如何起作用,如何发展变化
%4:在某些特殊且复杂的任务中,可能存在一只蚂蚁需要反复访问某一个节点的情况(例如某一个城市为众多城市的中转站,某些城市之间不能直达)
%5:有的任务不需要返回出发点,只需要访问所有的城市
%6:有些任务是不需要遍历所有的城市,只需要遍历给定的城市,(从一个城市到达另一个城市,中途必须经过哪些城市,有的城市可选(有的城市必须经过,有的城市不是必须经过的))
%%
while NC<=NC_max %停止条件之一:达到最大迭代次数,停止
%%第二步:将m只蚂蚁放到n个城市上,为每一只蚂蚁随机指定一个开始的城市
Randpos=[]; %随即存取
for i=1:(ceil(m/n))%m只蚂蚁除以n个城市,向上取整
Randpos=[Randpos,randperm(n)];%一个列向量,每一次生成一个1:n的乱序整数序列,重复ceil(m/n)次
end
Tabu(:,1)=(Randpos(1,1:m))';%随机指定每一只蚂蚁一开始所在的城市
%%第三步:m只蚂蚁按概率函数选择下一座城市,完成各自的周游
%概率函数如何给定???
for j=2:n %遍历还没有旅行过的城市
for i=1:m %遍历所有的蚂蚁
visited=Tabu(i,1:(j-1)); %记录已访问的城市,避免重复访问
J=zeros(1,(n-j+1)); %初始化待访问的城市列表
P=J; %待访问城市的选择概率分布
Jc=1;%记录每一只蚂蚁还没有访问过的城市个数
for k=1:n%此处遍历将每一只蚂蚁还没有访问过的城市查找出来,编程中可以使用队列优化此处
if length(find(visited==k))==0 %满足条件表示第k个城市没有被第i只蚂蚁访问过
J(Jc)=k; %将第k个城市添加到待访问城市列表中
Jc=Jc+1; %待访问的城市个数自加1
end
end
%下面计算待选城市的概率分布
for k=1:length(J)%遍历每一只蚂蚁还没有访问过的城市
%Alpha信息素指数,启发因子指数
%访问概率为信息素的Alpha次方与启发因子的Beta次方的乘积(启发式因子及遇到障碍物改变方向的可能性)
P(k)=(Tau(visited(end),J(k))^Alpha)*(Eta(visited(end),J(k))^Beta);%visited(end)表示上一次访问过的城市,J(k)表示还没有访问过的城市,Tau(visited(end),J(k))表示两个城市之间的信息素
end
P=P/(sum(P));
%按概率原则选取下一个城市
Pcum=cumsum(P); %cumsum,按列增量求和,即从第一个元素到当前元素的累加和
Select=find(Pcum>=rand); %若计算的概率大于原来的就选择这条路线,rand每次返回一个0到1之间的随机数
to_visit=J(Select(1));%随机选择下一个要访问的城市
Tabu(i,j)=to_visit;%第i只蚂蚁,第j次访问的城市为to_visit,为什么这样选出来的下一个城市是合理的???????????????
end
end
if NC>=2%如果不是第一次迭代,
Tabu(1,:)=R_best(NC-1,:);
end
%%第四步:记录本次迭代最佳路线
L=zeros(m,1); %开始距离为0,m*1的列向量,L存储每一只蚂蚁本次迭代周游所有城市的路径长度
for i=1:m
R=Tabu(i,:);%R为第i只蚂蚁的旅游路线
for j=1:(n-1)
L(i)=L(i)+D(R(j),R(j+1)); %原距离加上第j个城市到第j+1个城市的距离
end
L(i)=L(i)+D(R(1),R(n)); %一轮下来后走过的距离,从最后一个点回到原点的距离
end
L_best(NC)=min(L); %最佳距离取最小,第NC次迭代的最短路径长度
pos=find(L==L_best(NC));
R_best(NC,:)=Tabu(pos(1),:); %此轮迭代后的最佳路线
L_ave(NC)=mean(L); %此轮迭代后的平均距离
NC=NC+1 %迭代继续
%%第五步:更新信息素
Delta_Tau=zeros(n,n);%开始时信息素为n*n的0矩阵,本次遍历对信息素矩阵的增量
for i=1:m%遍历每一只蚂蚁
for j=1:(n-1)
%Tabu第i行为第i只蚂蚁在本次便利所有城市时访问城市的顺序列表,L(i)为第i只蚂蚁本次旅行路程的总长度,Q为一个正数,路程越短,对信息素的增量越大(Q/L(i))
Delta_Tau(Tabu(i,j),Tabu(i,j+1))=Delta_Tau(Tabu(i,j),Tabu(i,j+1))+Q/L(i);
%此次循环在路径(i,j)上的信息素增量
end
Delta_Tau(Tabu(i,n),Tabu(i,1))=Delta_Tau(Tabu(i,n),Tabu(i,1))+Q/L(i);%从最后一个城市返回第一个城市(有些任务是不需要返回原点,只需要遍历所有城市,要求路径最短)
%此次循环在整个路径上的信息素增量
end
%Rho为信息素蒸发系数
Tau=(1-Rho).*Tau+Delta_Tau; %考虑信息素挥发,更新后的信息素
%%第六步:禁忌表清零
Tabu=zeros(m,n); %%直到最大迭代次数
end
%%第七步:输出结果
Pos=find(L_best==min(L_best)); %找到最佳路径(非0为真)
Shortest_Route=R_best(Pos(1),:) %最大迭代次数后最佳路径
Shortest_Length=L_best(Pos(1)) %最大迭代次数后最短距离
subplot(1,2,1) %绘制第一个子图形
DrawRoute(C,Shortest_Route) %画路线图的子函数
subplot(1,2,2) %绘制第二个子图形
plot(L_best)
hold on %保持图形
plot(L_ave,'r')
title('平均距离和最短距离') %标题
function DrawRoute(C,R)
%%=========================================================================
%%
%% 画路线图的子函数
%%-------------------------------------------------------------------------
%% C Coordinate 节点坐标,由一个N×2的矩阵存储
%% R Route 路线
%%=========================================================================
N=length(R);
scatter(C(:,1),C(:,2));
hold on
plot([C(R(1),1),C(R(N),1)],[C(R(1),2),C(R(N),2)],'g')
hold on
for ii=2:N
plot([C(R(ii-1),1),C(R(ii),1)],[C(R(ii-1),2),C(R(ii),2)],'g')
hold on
end
title('旅行商问题优化结果 ')