sanneal.m

本程序是基本的MATLB退火模拟算法

function s=sanneal()

%------基本模拟退火优化算法-----------

%------名称：基本模拟退火优化算法（SANNEAL）

%------作用：求解函数优化问题

%------说明：全局性，并行性，高效的智能算法

%------作者：张弓长

%------单位：湖北大学

%------时间：2008年8月17日 

%------------------------------------------------------------------

%------初始格式化--------------------------------------------------


%clear all;

clc;

format long;
tic
%------给定初始化条件----------------------------------------------
T=1e-1;

eps=10^(-6);           %设置精度(在已知最小值时候用)

a=-100;b=100;D=2;
L=10;dt=0.915;scale=0.1;k=0.5;K=195;
%f='tripod';D=2;%Y but undefinite
f='Foxhole';D=2;%unkonn

j=0;i=0;

x1=rand(1,D).*(b-a)+a;u1=feval(f, x1);x=x1;
while j<K
    
    for l=1:L
        r=ceil(D*rand);
        scale=0.1*exp(-j/(K*20));
        x(l,r)=x1(r)+(-1+2*rand)*scale*(b-a);
        if x(l,r)<a
            x(l,r)=b-(a-x(l,r));
        end
        if x(l,r)>b
            x(l,r)=a+(x(l,r)-b);
        end
        u2=feval(f, x(l,:));
        df=abs(u2-u1);
        if u2<u1 
            u1=u2;x1=x(l,:);
        elseif exp(-df/(T))>=rand
            u1=u2;x1=x(l,:);
            
        end
        
        i=i+1;
        
    end
    T=T*dt;j=j+1;
end 

while j<0%K
    for l=1:L
        r=ceil(D*rand);
        scale=0.1*exp(-j/(K*20));
        x(r)=x1(r)+(-1+2*rand)*scale*(b-a);
        if x1(r)<a
            x(r)=b-(a-x(r));
        end
        if x(r)>b
            x(r)=a+(x(r)-b);
        end
        u2=feval(f, x);
        df=abs(u2-u1);
        if df<0
            u1=u2;x1=x; 
        elseif df>=0
            if rand<exp(-df/T)
                u1=u2;x1=x; 
            else
                u1=u1;x1=x1;
            end
        end
    end
    j=j+1;T=dt*T;
end
u1
s=[i,u1,x1];
format;
t1=toc
%%%------算法结束----------------------------------------------





function [out]=tripod(in)

x1=in(:,1);
x2=in(:,2);

px1=((x1) >= 0);
px2=((x2) >= 0);

out= ( px2.*(1+px1) ...
    + abs(x1 + 50*px2.*(1-2*px1))...
    + abs(x2 + 50*(1-2.*px2)) );


function [out]=Foxhole(in)
x=in(:,1);
y=in(:,2);

% location of foxholes, you can change this, these are what DeJong used
% if you change 'em, a{1} and a{2} must each have same # of elements
a{1} = [... 
        -32 -16   0  16  32 ;...
        -32 -16   0  16  32 ;...
        -32 -16   0  16  32 ;...
        -32 -16   0  16  32 ;...
        -32 -16   0  16  32 ;...
    ];

a{2} = [...
        -32 -32 -32 -32 -32 ;...
        -16 -16 -16 -16 -16 ;...
        0   0   0   0   0 ;...
        16  16  16  16  16 ;...
        32  32  32  32  32 ;...
    ];

term_sum=0;   
for j=1:prod(size((a{1})))
    ax=a{1} (j);
    ay=a{2} (j);
    term_num = (x - ax).^6 + (y - ay).^6;
    term_sum=term_sum+ 1./(j+term_num);
end
out = 1./(.002 + term_sum);