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： 各代最佳路线的长度
%=========================================================================
%%第一步：变量初始化
m=31;Alpha=1;Beta=5;Rho=0.1;NC_max=200;Q=100; %总共31只蚂蚁，200次迭代，那么每次迭代时每个城市随机分配之多一个蚂蚁
load('C:\Users\sobigdog\Documents\MATLAB\C.mat') % 这个路径是C矩阵的路径，事先建好的，表示n个城市的坐标，n×2矩阵。
n=size(C,1);% n表示问题的规模（城市个数），size(C,1)表示输出的是矩阵C的行数，即输出为n即城市的个数
D=zeros(n,n); % D表示完全图的赋权邻接矩阵，D为n行n列的空矩阵，那么从下面的定义可知该矩阵储存的是城市之间的距离
for i=1:n
for j=1:n
if i~=j % ~=在Matlab中表示“不等于”
D(i,j)=((C(i,1)-C(j,1))^2+(C(i,2)-C(j,2))^2)^0.5; %城市i和城市j之间的距离,C为n×2的矩阵
else
D(i,j)=eps; %eps应该是matlab中的一个极小的数，在matlab中输入eps，其值为2.2204e-16;
end
D(j,i)=D(i,j); %矩阵D为对称矩阵，该矩阵中保存两两城市之间的距离，注意不只是相邻城市之间的距离
end
end %计算完两两城市之间的距离
Eta=1./D; %Eta为启发因子，这里设为距离的倒数，注意的是1./D是矩阵的除法，所以Eta也是一个n×n的矩阵
Tau=ones(n,n); %Tau为信息素矩阵
Tabu=zeros(m,n); %存储并记录路径的生成，m行n列的零矩阵，表示m只蚂蚁的各自的解
NC=1; %迭代计数器，第一次迭代
R_best=zeros(NC_max,n); %各代最佳路线，之所以是n列，因为这n列表示出一个解即n个城市的排列顺序，从这两行来看，在每一次迭代后都会保存最优路线及其长度，并不会再下次迭代时将其覆盖掉
L_best=inf.*ones(NC_max,1); %各代最佳路线的长度，在Matlab中inf表示的是infinite即无穷大的意思，之所以这样设置的原因是什么？？？？。
L_ave=zeros(NC_max,1); %各代路线的平均长度
while NC<=NC_max %停止条件之一：达到最大迭代次数
%%第二步：将m只蚂蚁放到n个城市上
Randpos=[]; %随机的位置为空
for i=1:(ceil(m/n)) %ceil函数的作用是朝正无穷方向取整，比如m/n=3.12，则ceil(m/n)=4
Randpos=[Randpos,randperm(n)]; %在Matlab中randperm的最常用的用法是返回一个从1~n的包含n个数的随机排列（每个数字只出现一次），因此每次循环都会产生n个随机数，相当于随机n个城市序列，比如randperm(5)———ans= 3 5 1 2 4
%再运行一次randperm(5)———ans= 1 2 3 4 5
end
Tabu(:,1)=(Randpos(1,1:m))'; %将m个蚂蚁随机，每个蚂蚁放到前面产生的城市序列中，每个蚂蚁一个城市，需要m个，Tabu(:,1)表示的就是所有蚂蚁的初始城市，所以括号内为1
%第三步：在第二步中m只蚂蚁都随机选择了初始城市，接下来开始游历，即m只蚂蚁按概率函数选择下一座城市，完成各自的周游
for j=2:n %思考为什么j从2开始到n？？？因为我们不知道m个蚂蚁的随机初始城市的编号是什么，可能是1，也可能是10，或者37，
for i=1:m %每只蚂蚁都要循环一次
visited=Tabu(i,1:(j-1)); %已访问的城市，每只蚂蚁都是从城市1开始判断，注意的是Tabu(i,1:(j-1))的输出并不是从1,2,3,···这样排列，因为城市的编号是随机排列的
J=zeros(1,(n-j+1)); %待访问的城市的编号（也就是待访问城市的个数）
P=J; %待访问城市的选择概率分布，
Jc=1; %定义一个初值，该值是用在表示待访问的城市的个数
for k=1:n % k表示城市
if length(find(visited==k))==0 % find函数用于返回所需要元素的所在位置，因此find（）语句找到visited中等于k的元素在数组visited中的位置，length可以理解为统计个数，如果说a=[1 2 3 4 5 2]，find(a==2)→ans =2 6，而此时length(find(a==2))→ans=2,也就是表示有2个元素满足条件，因此如果：
% if length(find(visited==k))==0，那么表示没有访问过城市k。
J(Jc)=k; %注意上面已经指出J是待访问的城市，所以，如果判断上式满足的话，就将城市k列入到J数组中作为未访问的城市
Jc=Jc+1; %J中元素是未访问过的城市的编码，J（Jc）中的Jc是1~n-j+1，表示有n-j+1个未访问过的城市，而且J中元素是从小到大的顺序排列的
end
end %上上步是随机分别m只蚂蚁到n个城市，而这两个for循环完成的是计算每只蚂蚁还没有访问的城市，将这些城市序号保存在数组J中
%下面计算待选城市的概率分布
for k=1:length(J) %此处k就表示待访问的城市的编号，比如说有10个待访问的城市，即length(J)=10，那么k的取值就是1~10，该for循环计算的是最后一个访问过的城市（也就是蚂蚁现在的城市）与所有待访问的城市间的访问概率
P(k)=(Tau(visited(end),J(k))^Alpha)*(Eta(visited(end),J(k))^Beta);
% 此处注意该式中的结果还不是蚂蚁从城市visited(end)到城市J(k)的概率，因为我们知道：
end
P=P/(sum(P)); %如果除以了sum(P)，才是我们想要得到的概率
Pcum=cumsum(P); %cumsum虽然也是求和，但是通常用于计算一个数组各行的累加值。因此该式就是对概率累加求和，P是1行，n-j+1的数组，但是cumsum(P)并不是我们想象的是一个数字，经过查询可知，当a是向量时，返回向量中各元素分别是该元素在a中位置之前所有元素之和，当a是矩阵时，默认b中每一列都是对应列前面各行元素之和 ，例如a=[1 2 3 4 5 2]，则cumsum(a)→ans=[1 3 6 10 15 17]，因此，Pcum应该也就是1行，n-j+1的，即列数等于未过城市的数目。
Select=find(Pcum>=rand);
to_visit=J(Select(1)); %选中的城市幅值给to_visit，之所以取J(Select(1))，重点是Select(1)，是因为概率肯定越加越大，根据赌轮盘选择原理，肯定要选择首先满足条件的那个区段，也就是说选择出下一个访问概率较大的城市，注意的是只是选中的概率大但是并不是一定选中
Tabu(i,j)=to_visit; %该行的作用就是将to_visit的值作为第i只蚂蚁在该次迭代中第j个城市，即确定了下一个去访问的城市。这行命令很重要，之前我们不了解为什么第46行命令中j要从2开始到n，j从2到n就代表 Tabu(i,j)中的j
end
end %到这里就对该次迭代中的所有蚂蚁的路线确定了
if NC>=2 %我们没有发现迭代次数NC会逐渐增加？？？？？？？，注意在第87中才出现了 NC=NC+1;
Tabu(1,:)=R_best(NC-1,:); %更新一下？？？？为什么将上一次迭代结果的最优路线作为第1只蚂蚁的解？？？？？？？？？？
end
%%第四步：记录本次迭代最佳路线
L=zeros(m,1); % L记录的是m只蚂蚁的路线的长度
for i=1:m
R=Tabu(i,:); % R应该表示的是第i只蚂蚁的城市
for j=1:(n-1)
L(i)=L(i)+D(R(j),R(j+1)); %在前面已经定义过了函数D（i,j）的作用是计算城市i和城市j之间的距离，那么可知L定义的也是距离，那么上式L(i)输出的就是第i只蚂蚁的解的总共的长度（没有加上返回路线长度）
end
L(i)=L(i)+D(R(1),R(n)); %之所以还要加上第1和城市和最后一个城市间的距离，是因为蚂蚁还要从最后一个城市回到出发点即第一个城市
end
L_best(NC)=min(L); %该min也是Matlab中的一个函数，很显然求的是最小值，该行的作用就是求出最优的一条路线并作为该代中最佳路线的长度
pos=find(L==L_best(NC)); %可能是考虑到有多条路线的长度同时最短，find函数的输出对应的应该是蚂蚁的编号，该编号就是路线最短的蚂蚁的编号
R_best(NC,:)=Tabu(pos(1),:); %选出满足最小长度的蚂蚁中的第一只蚂蚁的路线作为最优路线
L_ave(NC)=mean(L); %mean函数求解所有蚂蚁的路线长度的平均值
NC=NC+1; %下一步迭代，该行命令仍然在第37行的while NC<=NC_max的循环中
%%第五步：更新信息素，注意到这是在完成一次迭代循环后再更新信息素
Delta_Tau=zeros(n,n); %开始时信息素为n行n列的零矩阵，思考问题：如果每次都设置为零，那么是否每次迭代后信息素不会叠加吗？？？？？？？，对于该问题的回答就是：该行只是信息素的释放，而我们知道信息素的更新包括信息素的挥发和信息素的释放
for i=1:m
for j=1:(n-1)
Delta_Tau(Tabu(i,j),Tabu(i,j+1))=Delta_Tau(Tabu(i,j),Tabu(i,j+1))+Q/L(i); %上式是两个相邻城市间路线的信息量增量。但注意到L(i)是整条线路的长度，有些类似蚁周系统，但是蚁周系统中又没有这个系数Q，该如何理解？？？？？？
end
Delta_Tau(Tabu(i,n),Tabu(i,1))=Delta_Tau(Tabu(i,n),Tabu(i,1))+Q/L(i); %这是终点城市和起点城市间直接连线上的信息量的增量
end %用两个for循环完成了m个蚂蚁的信息量的更新
Tau=(1-Rho).*Tau+Delta_Tau;
%%第六步：禁忌表清零，个人理解Tabu中储存的是m只蚂蚁的路线，将Tabu清零是为了下次重新迭代做准备
Tabu=zeros(m,n); %随机分布m个蚂蚁到n个城市时，有的说法是每个城市至多分布1个蚂蚁，该如何理解？？
%清空禁忌表是为了让m个蚂蚁重新游历�
