模擬退火算法matlab源程序

1、function MinD,B estPath=M ainAneal(CityPosition,pn)function MinD,Best Path=MainAneal2(Cit yPosition,pn)% 此題以中國 31 省會城市的最短旅行路徑為例，給出 TSP 問題的模擬退火程序%Cit yPosition_31=1304 2312;3639 1315;4177 2244;3712 1399;3488 1535;3326 1556;. %3238 1229;4196 1044;4312?790;4386?570;3007 1970;2562 1756;.%2788 1491;2381

2、1676;1332?695;3715 1678;3918 2179;4061 2370;.%3780 2212;3676 2578;4029 2838;4263 2931;3429 1908;3507 2376;.%3394 2643;3439 3201;2935 3240;3140 3550;2545 2357;2778 2826;2370 2975;%T0=clockglobal path p2 D; m,n=size(CityPosition);% 生成初始解空間，這樣可以比逐步分配空間運行快一些 TracePath=zeros(1e3,m);Dist ance=inf*zeros(1,

3、1e3);D = sqrt(Cit yPosition( :,ones(1,m) - CityPosition( :,ones(1,m)').A2 +(Cit yPosition( : ,2*ones(1,m) - CityPosition( :,2*ones(1,m)').A2 );% 將城市的坐標矩陣轉(zhuǎn)換為鄰接矩陣(城市間距離矩陣)for i=1:pnpath(i,:)=randperm(m);% 構(gòu)造一個初始可行解endt=zeros(1,pn);p2=zeros(1,m);iter_max=100;%input(' 請輸入固定溫度下最大迭代次數(shù) iter_max

4、=' );m_max=5;%input(' 請輸入固定溫度下目標函數(shù)值允許的最大連續(xù)未改進次數(shù) m_nax=' ) ; % 如果考慮到降溫初期新解被吸收概率較大，容易陷入局部最優(yōu)% 而隨著降溫的進行新解被吸收的概率逐漸減少，又難以跳出局限% 人為的使初期 iter_max,m_max 較小，然后使之隨溫度降低而逐步增大 ,可能% 會收到到比較好的效果T=1e5;N=1;tau=1e-5;%input(' 請輸入最低溫度 tau=' ); %nn=ceil(log10(t au/T)/log10(0.9);while?T>=tau%&m_nu

5、m<m_max? ? ? ?iter_num=1;% 某固定溫度下迭代計數(shù)器 m_num=1;% 某固定溫度下目標函數(shù)值連續(xù)未改進次數(shù)計算器%it er_max=100; %m_max=10;%ceil(10+0.5*nn-0.3*N);while m_num<m_max&it er_num<iter_max%MRRTT(Metropolis, Rosenbluth, Rosenbluth, Teller, Teller) 過程 :% 用任意啟發(fā)式算法在 path 的領(lǐng)域 N(path) 中找出新的更優(yōu)解for i=1:pnLen1(i)=sum(D(path(i,1

6、:m-1)+m*(path(i,2:m)-1) D(path(i,m)+m*(path(i,1)-1);% 計算一次行遍所有城市的總路程path2(i,: )=ChangePath2(path(i,: ),m);%更新路線Len2(i)=sum(D(path2(i,1:m-1)+m*(path2(i,2:m)-1) D(path2(i,m)+m*(path2(i,1)-1);end%Len1%Len2%if Len2-Len1<0|exp(Len1-Len2)/(T)>randR=rand(1,pn);%Len2-Len1<t|exp(Len1-Len2)/(T)>Ri

7、f find(Len2-Len1<t|exp(Len1-Len2)/(T)>R)=0)path(find(Len2-Len1<t|exp(Len1-Len2)/(T)>R)=0), : )=path2(find(Len2-Len1<t|exp(Len1-Len 2)/(T)>R)=0), : );Len1(find(Len2-Len1<t|exp(Len1-Len2)/(T)>R)=0)=Len2(find(Len2-Len1<t|exp(Len1-Len2)/( T)>R)=0);TempMinD,T empIndex=min(Le

8、n1);%T empMinDTracePath(N,: )=path(T empIndex,: );Dist ance(N,: )=TempMinD;N=N+1;T=T*0.9 ；m_num=0;elsem_num=m_num+1;enditer_num=it er_num+1;endT=T*0.9%m_num,iter_num,NendMinD,Index=min(Distance);BestPath=TracePath(Index,: );disp(MinD)%T1=clock% 更新路線子程序function p2=ChangePath2(p1,CityN um)global p2;while(1)R=unidrnd(CityNum,1,2);if abs(R(1)-R(2)>1break;endendR=unidrnd(CityNum,1,2);I=R(1);J=R(2);%len1=D(p(I),p(J)+D(p(I+1),p(J+1);%len2=D(p(I),p(I+1)+D(p(J),p(J+1);if I<Jp2(1:I)=p1(1:I);p2