clear all
%设定城市坐标
Coord=[0.6783 0.4195 0.6 0.2949 0.2727 0.4233 0.4174 0.7752 0.6778 0.7488 0.1621;
0.2236 0.1634 0.9439 0.2463 0.1293 0.561 0.6424 0.4536 0.7219 0.4609 0.4781];
t0=1; %初始温度
iLK=30; %内循环迭代次数
oLK=60; %外循环最大迭代次数
lam=0.9;
istd=0.001; %若内循环函数值方差小于istd则停止
ostd=0.001; %若外循环函数值方差小于ostd则停止
ilen=5; %内循环保存的目标函数值个数
olen=5; %外循环保存的目标函数值个数
%程序主体
m=length(Coord); %城市的个数
fare=citydistance(Coord); %路径费用
path=1:m; %初始路径
pathfar=pathfare(fare,path);%路径费用
ores=zeros(1,olen); %外循环保存的目标函数值
e0=pathfar; %能量初始值
t=t0; %温度
for out=1:oLK %外循环模拟退火过程
ires=zeros(1,ilen); %内循环保存的目标函数
for in=1:iLK %内循环模拟热平衡过程
[newpath,~]=swap(path,1); %产生新状态
e1=pathfare(fare,newpath);%新状态能量
%Metropolis抽样稳定准则
r=min(1,exp(-(e1-e0)/t));
if rand<r
path=newpath; %更新最新状态
e0=e1; %%%%%%%%%%%%% 注意:e1有可能大于e0,这一步是跳出局部极小值的关键
end
ires=[ires(2:end) e0]; %保存新状态能量
%内循环终止准则:连续ilen个状态能量波动小于istd
if std(ires,1)<istd
break;
end
end
%保存新状态能量
ores=[ores(2:end) e0];
%外层循环终止准则:连续olen个状态能量波动小于ostd
if std(ores,1)<ostd
break;
end
t=lam*t;
end
pathfar=e0;
%输入结果
fprintf('近似最优路径是:\n')
disp(path)
fprintf('近似最优路径长度是\t pathfar=')
disp(pathfar);
myplot(path,Coord,pathfar);