clear;
clc;
tic;
%数据导入
D=xlsread('E:\LXHAI.xlsx','A2:I10');
S=xlsread('E:\LXHAI.xlsx','A15:I15');% 卸货时间
q=xlsread('E:\LXHAI.xlsx','A13:I13');%吨
a=xlsread('E:\LXHAI.xlsx','A17:I17');
b=xlsread('E:\LXHAI.xlsx','A19:I19');
%初始化参数
B1=zeros(200,9);
B2=zeros(200,9);
B3=zeros(200,9);
tabu=[];
v=50;
n=9;
m=3;
t=zeros(m,n); %蚂蚁游历的时间
Time=zeros(1,n+1);
alpha=1;
beta=5;
rho=0.1;
Q=1;
Eta=1./D;
weight=zeros(m,n); %蚂蚁的载重量
X=zeros(m,n);
Tau=ones(n,n);
Table=zeros(m,n);
iter=1;
iter_max=200;
Route_best=zeros(iter_max,n);
Length_best=zeros(iter_max,1);
%迭代寻找最佳路径
while iter<=iter_max
Table(:,1)=[1;1;1];
citys_index=1:n;
%逐个蚂蚁路径选择
for i=1:m
%逐个城市路径选择
for j=2:n
PP=zeros(n,2);
if(i==1)
tabu=Table(i,1:(j-1)); %已访问城市记录表（禁忌表）
else
if((Table(i,(j-1)))~=1)
tabu=[tabu,Table(i,(j-1))];
end
end
allow_index=~ismember(citys_index,tabu);%产生待访问城市位置对应的逻辑数值
allow=citys_index(allow_index);%待访问城市集合
%%%%%%%%%%%%%%%%
%计算 城市间转移概率
Time(1,1)=t(i,j-1);
u=1;
Time1=zeros(m,n);
weight1=zeros(m,n);
P=zeros(1,n);
for k=1:length(allow)
Time(1,k+1)=Time(1,1)+D(Table(i,j-1),allow(k))/v;
if((Time(1,k+1)>=a(allow(k)))&&(Time(1,k+1)<=b(allow(k)))) %%%t(i,j)=t(i,j-1)+Time(1,k+1);
X(i,j)=1;
Time(1,k+1)= Time(1,k+1)+S(allow(k));
else
X(i,j)=0;
end
P(k)=Tau(tabu(end),allow(k))^alpha*(Eta(tabu(end),allow(k))*X(i,j))^beta;
if(P(k)~=0)
PP(u,1)=P(k);
PP(u,2)=allow(k);%%存储城市名
Time1(i,u)=Time(1,k+1);
weight1(i,u)=weight(i,j-1)+q(allow(k)); %%%%%
u=u+1;
end
end
P=zeros(1,n);
g=0;
l=1;
L=zeros(n,2);%%1表示城市信息，2表示
for d1=1:length(PP(:,1))
if(weight1(i,d1)<=8)
P(l)=PP(d1,1);
L(l,1)=PP(d1,2);%%存储城市名
L(l,2)=d1; %调用时间、载重时有用
l=l+1;
else
g=g+1;
end
end
if(L(1,1)==0)
break;
end
if(g==length(PP(:,1)))%%%length
Table(i,j)=1;
break;
else
P=P/sum(P);
%用轮盘赌法选择下一个访问位置
Pc=cumsum(P);
target_index=find(Pc>=rand);
target=L(target_index(1));
Table(i,j)=target;
weight(i,j)= weight(i,j-1)+q(target);
t(i,j)= t(i,j-1)+D(Table(i,j-1),target)/v+S(target);
end
end
end
if(~isempty(allow))
continue;
end
B1(iter,:)=Table(1,:);
B2(iter,:)=Table(2,:);
B3(iter,:)=Table(3,:);
%%计算各个蚂蚁的路径距离
Length=zeros(m,1);
for i=1:m
Route=Table(i,:);
for j=1:(sum(Table(i,:)>0)-1)
Length(i)=Length(i)+D(Route(j),Route(j+1));
end
Length(i)=Length(i)+D(Route(sum(Table(i,:)>0)),Route(1));
end
%%%%%
Length1=sum(Length);
if iter==1
Length_best(iter)=Length1;
Route_best=Table;
else
if(Length_best(iter-1)<Length1)
Route_best=Table;
end
end
%更新信息素
Delta_Tau=zeros(n,n);
%逐个蚂蚁计算
for i=1:m
%逐个城市计算
for j=1:(sum(Table(i,:)>0)-1)
Delta_Tau(Table(i,j),Table(i,j+1))=Delta_Tau(Table(i,j),Table(i,j+1))+Q/Length(i);
end
Delta_Tau(Table(i,sum(Table(i,:)>0)),Table(i,1))= Delta_Tau(Table(i,sum(Table(i,:)>0)),Table(i,1))+Q/Length(i);
end
Tau=(1-rho)*Tau+Delta_Tau;
%迭代次数加1，清空路径记录表
iter=iter+1;
Table=zeros(m,n);
end
disp(['最短路程:' num2str(Length_best(1))]);
disp(['最短路径:' num2str([Route_best(1,:),Route_best(2,:),Route_best(3,:)])]);
toc;
