模擬退火算法matlab源程序

1、function MinD,BestPath=MainAneal(CityPosition,pn) function MinD,BestPath=MainAneal2(CityPosition,pn)%此題以中國31省會城市的最短旅行路徑為例，給出TSP問題的模擬退火程序%CityPosition_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 1676;1332?

2、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);Distance=inf*zeros(1,1e3);D = sqrt(

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

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

5、定溫度下迭代計數(shù)器 m_num=1;%某固定溫度下目標(biāo)函數(shù)值連續(xù)未改進次數(shù)計算器 %iter_max=100; %m_max=10;%ceil(10+0.5*nn-0.3*N); while m_num<m_max&iter_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:m-1)+m*(path(i,2:m)-1) D(path(i,m)+m*(

6、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)>R if find(Len2-Len1<t|exp(Len1-Le

7、n2)/(T)>R)=0)path(find(Len2-Len1<t|exp(Len1-Len2)/(T)>R)=0), : )=path2(find(Len2-Len1<t|exp(Len1-Len2)/(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,TempIndex=min(Len1);%TempMinDTracePath(N,: )=path(Te

8、mpIndex,: );Distance(N,: )=TempMinD;N=N+1;T=T*0.9；m_num=0; elsem_num=m_num+1; end iter_num=iter_num+1;endT=T*0.9%m_num,iter_num,Nend MinD,Index=min(Distance);BestPath=TracePath(Index,: );disp(MinD)%T1=clock%更新路線子程序function p2=ChangePath2(p1,CityNum)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)