simulateannealing.asv

在matlab环境下用人工智能算法--模拟退火算法求函数极植的matlab程序

% 使用模拟退火法求函数 f(x,y)=5sin(xy)+ x^2 + y^2 的最小值, 

% 根据题意,我们设计冷却表进度表为:
% 初温度为100
% 衰减参数为0.95
% 马可夫链长度为10000
% 使用Metropolis接受准则进行模拟,Metropolis的步长为0.02
% 结束条件为根据上一个最优解与最新的一个最优解的之差小于某个容差.
function r=SimulateAnnealing()
  format long
  tic
  %搜索的最大区间
  XMAX=4;
  YMAX=4;
  
  %冷却表参数
  MarkovLength=10000;   %马可夫链长度
  DecayScale=0.95;      %衰减参数
  StepFactor=0.02;      %步长因子
  Temperature=100;      %初始温度
  Tolerance=1e-8;       %容差
  
  %PreX,NextX  分别为prior and next value of x
  %PreY,NextY  分别为prior and next value of y
  %PreBextX,PreBestY 是上一个最优解
  %BestX,BestY   是最终解
  AcceptPoints=0;     %Metropolis过程中总接受点
  
  %随机选点
  PreX=-XMAX*rand;
  PreY=-YMAX*rand;
  PreBestX=PreX;
  PreBestY=PreY;
  BestX=PreX;
  BestY=PreY;
  
  %每迭代一次退火一次(降温),直到满足迭代条件为止
 %Decay=Temperature*DecayScale;
  while(1)
     Temperature=Temperature*DecayScale;
     AcceptPoints=0;
     
     %当前温度下迭代loop(即MARKO链长度)次
     for i=1:MarkovLength
         
        % 1) 在此点附近随机选下一点
        while(1)
            NextX=PreX+StepFactor*XMAX*(rand-0.5);
            NextY=PreY+StepFactor*YMAX*(rand-0.5);
            if (NextX>=-XMAX & NextX<=XMAX & NextY>=-YMAX & NextY<=YMAX)
                break;
            end
        end    

        % 2) 是否全局最优解
        if ObjectFunction(BestX,BestY)>ObjectFunction(NextX,NextY)
            %保留上一个最优解
            PreBestX=BestX;
            PreBestY=BestY;
            %此为新的最优解
            BestX=NextX;
            BestY=NextY;
        end              
        
        % 3) Metropolis过程              
        if ObjectFunction(PreX,PreY)-ObjectFunction(NextX,NextY)>0
            % 接受,此处lastPoint即下一个迭代的点以新接受的点开始
            PreX=NextX;
            PreY=NextY;
            AcceptPoints=AcceptPoints+1;
        else
            change=-1*(ObjectFunction(NextX,NextY)-ObjectFunction(PreX,PreY))/Temperature;
            if exp(change)>rand
                PreX=NextX;
                PreY=NextY;
                AcceptPoints=AcceptPoints+1;
            end %不接受,保存原解
        end   
     end
     %符合终止条件退出
     if (abs(ObjectFunction(PreBestX,PreBestY)-ObjectFunction(BestX,BestY))<Tolerance)
            BestX
            BestY
            BestF=ObjectFunction(BestX,BestY)
            break;
     end
  end    
  toc
  
  %graph
  grid on
  [X,Y]=meshgrid(-4:4,-4:4);
  Z=5*sin(X*Y)+X^2+Y^2;
  plot3(X,Y,Z)
  plot3(BestX,BestY,BestF,'*r')
  
  figure:clf
  mesh(X,Y,Z)
  grid on
  
  figure;clf;
  grid on
  surf(X,Y,Z)
  hold on
  plot3(BestX,BestY,BestF,'*r')
  hold off
  
function f=ObjectFunction(x,y)
f=5*sin(x*y)+x^2+y^2;