模拟退火二元法.docx

%matlab 程序实现 模拟退火算法程序 函数求极值（引用后修改，感谢 ARMYLAU）
%使用模拟退火法求函数f(x,y) = 3*cos(xy) + x + y2的最小值
%解：根据题意，我们设计冷却表进度表为：
%即初始温度为30
%衰减参数为0.95
%马可夫链长度为10000
%Metropolis的步长为0.02
%结束条件为根据上一个最优解与最新的一个最优解的之差小于某个容差。
%使用METROPOLIS接受准则进行模拟, 程序如下
%* 日期:2012－11-29
%* 作者:steven
%*　EMAIL:hxs2004@
%* 结束条件为两次最优解之差小于某小量
function [BestX,BestY]=SimulateAnnealing1
clear;
clc;
%// 要求最优值的目标函数,搜索的最大区间
  XMAX= 4;
  YMAX = 4;
%冷却表参数
 MarkovLength = 10000; %// 马可夫链长度
 DecayScale = 0.95; %// 衰减参数
 StepFactor = 0.02; %// 步长因子
 Temperature=30; %// 初始温度
 Tolerance = 1e-8; %// 容差
AcceptPoints = 0.0; %// Metropolis过程中总接受点
rnd =rand;
% 随机选点 初值设定
PreX = -XMAX * rand ;
PreY = -YMAX * rand;
PreBestX  = PreX;
PreBestY = PreY;
PreX = -XMAX * rand ;
PreY = -YMAX * rand;
BestX = PreX;
 BestY = PreY;
% 每迭代一次退火一次(降温), 直到满足迭代条件为止
mm=abs( ObjectFunction( BestX,BestY)-ObjectFunction (PreBestX, PreBestY));
while mm > Tolerance
Temperature=DecayScale*Temperature;
AcceptPoints = 0.0;
% 在当前温度T下迭代loop(即MARKOV链长度)次
for i=0:MarkovLength:1
% 1) 在此点附近随机选下一点
p=0;
while p==0
     NextX = PreX + StepFactor*XMAX*(rand-0.5);
    NextY = PreY + StepFactor*YMAX*(rand-0.5);
    if p== (~(NextX >= -XMAX && NextX <= XMAX && NextY >= -YMAX && NextY <= YMAX))
        p=1;
    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
 changer = -1 * ( ObjectFunction(NextX,NextY) - ObjectFunction(PreX,PreY) ) / Temperature ;
rnd=rand;
p1=exp(changer);
double (p1);
if p1 > rand             %// 不接受, 保存原解
    PreX=NextX;
    PreY=NextY;
    AcceptPoints=AcceptPoints+1;
end
end
end
mm=abs( ObjectFunction( BestX,BestY)-ObjectFunction (PreBestX, PreBestY));
end
disp('最小值在点:');
BestX
BestY
disp( '最小值为:{0}');
ObjectFunction(BestX, BestY)
end

****************************************************
子函数，目标函数值计算

function  value=ObjectFunction(x,y)
          value=3*cos(x*y)+x+y*y;
end

%使用模拟退火法求函数f(x,y)=sin(x*y)+x^2+y^2的最小值
format long
XMAX=4; %搜索的最大区间
YMAX=4; %搜索的最大区间
MarkovLength=10000; %马可夫链长度
DecayScale=0.95; %衰减参数0.95
StepFactor=0.02; %步长因子
Temperature=100; %初始温度
Tolerance=1e-8; %容差
AcceptPoints=0.0; %Metropolis过程中总接受点
PreX=-XMAX*rand; %初始的搜索值
PreY=-YMAX*rand; %初始的搜索值
PreBestX=PreX; %上一个最优解
PreBestY=PreY; %上一个最优解
BestX=PreX; %最终解
BestY=PreY; %最终解
while(1)
Temperature=Temperature*DecayScale; %每迭代一次退火一次（降温），直到满足迭代条件为止
AcceptPoints=0.0;
%在当前温度下迭代（即MARKOV链长度）次
for i=0:1:MarkovLength
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
%判断新产生点与原来最优点哪个更优
if(minfunction(BestX,BestY)>minfunction(NextX,NextY))
PreBestX=BestX; %保留上一个最优解
PreBestY=BestY;
BestX=NextX; %新的最优解
BestY=NextY;
end
%接受新产生的点为下一迭代的开始点
if(minfunction(PreX,PreY)-minfunction(NextX,NextY)>0)
PreX=NextX;
PreY=NextY;
AcceptPoints=AcceptPoints+1;
else
change=-1*(minfunction(NextX,NextY)-minfunction(PreX,PreY))/Temperature;
%不接受，保存原解
if(exp(change)>rand)
PreX=NextX;
PreY=NextY;
AcceptPoints=AcceptPoints+1;
end
end
end
%结束条件为根据上一个最优解与最新的一个最优解的之差小于某个容差
if(~(abs(minfunction(BestX,BestY)-minfunction(PreBestX,PreBestY))>Tolerance))
break

end
end
a=BestX
b=BestY
c=minfunction(BestX,BestY)

%%%%%%%%%%%%%子程序%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function minf=minfunction(x,y)
minf=sin(x*y)+x*x+y*y;%求目标函数的值