模拟退火二元法

%matlab程序实现模拟退火算法程序函数求极值(引用后修改，感谢ARMYLAU)%使用模拟退火法求函数f(x,y)=3*cos(xy)+x+y2的最小值%解：根据题意，我们设计冷却表进度表为：%即初始温度为30%衰减参数为0.95%马可夫链长度为10000%Metropolis的步长为0.02%结束条件为根据上一个最优解与最新的一个最优解的之差小于某个容差。%使用METROPOLIS接受准则进行模拟，程序如下%*日期:2012—11-29%*作者:steven%*EMAIL:hxs2004@126.com%*结束条件为两次最优解之差小于某小量function[BestX,BestY]=SimulateAnnealing1clear;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));whilemm>ToleranceTemperature=DecayScale*Temperature;AcceptPoints=0.0;%在当前温度T下迭代loop(即MARKOV链长度)次fori=0:MarkovLength:1%1)在此点附近随机选下一点p=0;whilep==0NextX=PreX+StepFactor*XMAX*(rand-0.5);NextY=PreY+StepFactor*YMAX*(rand-0.5);ifp==(〜(NextX>=-XMAX&&NextX<=XMAX&&NextY>=-YMAX&&NextY<=YMAX))p=1；endend%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;elsechanger=-1*(ObjectFunction(NextX,NextY)-ObjectFunction(PreX,PreY))/Temperature;rnd=rand;p1=exp(changer);double(p1);ifp1>rand%//不接受，保存原解PreX=NextX;PreY=NextY;AcceptPoints=AcceptPoints+1;endendendmm=abs(ObjectFunction(BestX,BestY)-ObjectFunction(PreBestX,PreBestY));enddisp('最小值在点:');BestXBestYdisp('最小值为:{0}');ObjectFunction(BestX,BestY)end****************************************************子函数，目标函数值计算functionvalue=ObjectFunction(x,y)value=3*cos(x*y)+x+y*y;end%使用模拟退火法求函数f(x,y)=sin(x*y)+xA2+yA2的最小值formatlongXMAX=4;%搜索的最大区间YMAX=4;%搜索的最大区间MarkovLength=10000;%马可夫链长度DecayScale=0.95;%衰减参数0.95StepFactor=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链长度)次fori=0:1:MarkovLengthwhile(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%退出循环endend%判断新产生点与原来最优点哪个更优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;elsechange=-1*(minfunction(NextX,NextY)-minfunction(PreX,PreY))/Temperature;%不接受，保存原解if(exp(change)>rand)PreX=NextX;PreY=NextY;AcceptPoints=AcceptPoints+1;endendend%结束条件为根据上一个最优解与最新的一个最优解的之差小于某个容差if(~(abs(minfunction(BestX,BestY)-minfunction(PreBestX,PreBestY))