www.pudn.com > spso.rar > spso.m, change:2005-04-07,size:6788b


function main() 
global POPSIZE DIMENSION W C1 C2 VMAX XMIN XMAX P PBEST individual 
POPSIZE = 20; 
W = 0.73; 
C1 = 2.05; 
C2 = 2.05; 
VMAX = 5; 
% clusters spherical coordinates bounds 
LJmin = [0,-1.0000,-3.0000,-6.0000,-9.1038,-12.7120,-16.5053,-19.8214,-24.1133,-28.4225,-32.7659,-37.9676,-44.3268]; 
% N = 3; radius = (3*N/(4*pi*sqrt(2)))^(1/3); 
% unit = [0 2*pi;0 pi;0 radius];bounds=unit; 
% for i=1:N-1 bounds = cat(1,bounds,unit); end 
% clear radius unit 
bounds = ones(2,1) * [-5.12 5.12]; 
XMIN = bounds(:,1)'; 
XMAX = bounds(:,2)'; 
DIMENSION = length(XMIN); 
% DIMENSION = 2; 
% XMIN = -5.12; 
% XMAX = 5.12; 
P=[]; % 1:DIMENSION 
PBEST=[]; 
individual = struct([]);%struct('number','best','bestfitness','fitness','speed'); 
 
% fprintf(1,'Minimization of LJ clusters: N = %d\n',N); 
total = 0; 
sum = 0; 
tic 
for j=1:1 
    initiate; 
    for i=1:2000 
        opt = PBEST; 
        flag = -1; 
        if(round(rand)) 
            direc = 1; 
        else  
            direc = 1; 
        end 
        for k=1:POPSIZE 
            if((individual(k).fitness == individual(k).bestfitness) && (individual(k).fitness == PBEST)) 
                flag = k; 
                break; 
            end 
        end 
        if(flag ~= -1) 
            generate(flag); 
        end        
        for k=1:POPSIZE 
            if(k ~= flag) 
                for t=1:DIMENSION 
                    individual(k).speed(t) = W * ( individual(k).speed(t) + direc * ( C1 * rand * (individual(k).best(t) - individual(k).number(t)) + ... 
                                                                                      C2 * rand * (P(t) - individual(k).number(t))) ); 
                    if (individual(k).speed(t) > VMAX)	individual(k).speed(t) = VMAX; end 
                    individual(k).number(t) = individual(k).number(t) + individual(k).speed(t); 
%                     if(individual(k).number(t) < XMIN(t))  
%                         individual(k).number(t) = XMIN(t); 
%                     elseif(individual(k).number(t) > XMAX(t))  
%                         individual(k).number(t) = XMAX(t); 
%                     end 
                end 
                calculation(k); 
                localbest(k); 
            end 
%                     direc = calc_div(direc,i); 
        end 
        globalbest(1); 
        if(abs(PBEST)<1e-3) 
%         if(abs(PBEST-LJmin(N))<1e-3)  
            total = total + 1; 
            sum = sum + i; 
            break 
        elseif((opt-PBEST)>1e-3) 
            fprintf(1,'\n%d %.4f\n',i,PBEST); 
        else 
            fprintf(1,'%d ',i);            
        end 
    end  
    fprintf(1,'\n%d\t%.4f\n',i,PBEST);     
end 
fprintf(1,'Total number is %d\nSum is %d\nAverage is %.1f\n',total,sum,sum/total); 
toc 
 
function calculation(num) 
global POPSIZE DIMENSION W C1 C2 VMAX XMIN XMAX P PBEST individual 
x = individual(num).number; 
% x1 = individual(num).number(2);%x0 = -1.3071;x1 = -1.4249; 
% s = 0; 
% h = 0; 
% for i=1:5 
%     s = s + i * cos((i-1)*x0+i); 
%     h = h + i * cos((i+1)*x1+i); 
% end 
% individual(num).fitness = s * h + (x0+1.42513)^2.0 + (x1+0.80032)^2.0; 
% s = (x0+x1+1.0)^2.0 * (19.0-14.0*x0+3.0*x0^2.0-14.0*x1+6.0*x0*x1+3.0*x1^2.0) + 1.0; 
% h = 30.0 + (2.0*x0-3.0*x1)^2.0 * (18.0-32.0*x0+12.0*x0^2.0+48.0*x1-36.0*x0*x1+27.0*x1^2.0); 
% individual(num).fitness = s * h; 
% Rastrigin's function (-5.12 5.12) 
len=length(x);A=8; 
fval=0; 
for i=1:len 
    fval=fval+x(i)^2-A*cos(2*pi*x(i)); 
end 
fval=fval+len*A; 
% spherical coordinates 
% x = individual(num).number; 
% N = DIMENSION / 3; 
% fval=0; 
% for i=1:N-1 
%     for j=i+1:N 
%         r=(x(3*j)*cos(x(3*j-2))*cos(x(3*j-1))-x(3*i)*cos(x(3*i-2))*cos(x(3*i-1)))^2 + ... 
%           (x(3*j)*cos(x(3*j-2))*sin(x(3*j-1))-x(3*i)*cos(x(3*i-2))*sin(x(3*i-1)))^2 + ... 
%           (x(3*j)*sin(x(3*j-2))-x(3*i)*sin(x(3*i-2)))^2; 
%         fval = fval + r^(-6) - 2*r^(-3); 
%     end 
% end 
individual(num).fitness = fval; 
return 
 
function generate(num) 
global POPSIZE DIMENSION W C1 C2 VMAX XMIN XMAX P PBEST individual 
individual(num).number = rand(1,DIMENSION).*(XMAX-XMIN) + XMIN; 
calculation(num); 
individual(num).bestfitness = individual(num).fitness; 
individual(num).best = individual(num).number; 
return 
 
function globalbest(number) 
global POPSIZE DIMENSION W C1 C2 VMAX XMIN XMAX P PBEST individual 
if(number == 0) 
    [PBEST,flag] = max([individual.fitness]); 
    P = individual(flag).number; 
else 
    for i=1:POPSIZE 
        if(individual(i).bestfitness < PBEST) 
            P = individual(i).best; 
            PBEST = individual(i).bestfitness; 
        end 
    end 
end 
return 
 
function initiate() 
global POPSIZE DIMENSION W C1 C2 VMAX XMIN XMAX P PBEST individual 
for i=1:POPSIZE 
    individual(i).number = rand(1,DIMENSION).*(XMAX-XMIN) + XMIN; 
    individual(i).speed = rand(1,DIMENSION).* VMAX; 
    individual(i).best = individual(i).number; 
    calculation(i); 
    individual(i).bestfitness = individual(i).fitness; 
end 
globalbest(0); 
return 
 
function localbest(num) 
global POPSIZE DIMENSION W C1 C2 VMAX XMIN XMAX P PBEST individual 
if(individual(num).bestfitness > individual(num).fitness) 
    individual(num).best = individual(num).number; 
end 
individual(num).bestfitness = individual(num).fitness; 
return 
 
function Stretch_Tech 
global POPSIZE DIMENSION W C1 C2 VMAX XMIN XMAX P PBEST individual 
gma1=1e+4/2;gma2=1/2;meu=1e-10; 
for i=1:POPSIZE 
    if(abs(individual(i).fitness-PBEST)>1e-5) 
        x = individual(i).number; 
        fval = individual(i).fitness; 
        G = fval + gma1 * sqrt(sum((x-P).^2)) * (sign(fval-PBEST)+1); 
        H = G + gma2 * (sign(fval-PBEST)+1) / tanh(meu*(G-PBEST)); 
        individual(i).fitness = H; 
    end 
end 
return 
 
function direc = calc_div(direc,iter) 
global POPSIZE DIMENSION W C1 C2 VMAX XMIN XMAX P PBEST individual 
D_low = 5e-6 * exp(1-iter/2000); 
D_high = 0.25 * exp(1-iter/2000); 
div = 0; 
x = [individual.number]; 
x = reshape(x,DIMENSION,POPSIZE); 
x = x'; 
for i=1:POPSIZE 
    for j=1:DIMENSION 
        div = div + (individual(i).number(j) - mean(x(:,j)))^2; 
    end 
end 
div = div / ( POPSIZE * max(XMAX-XMIN) ); 
if(direc>0 && div<D_low) 
    direc = -1; 
elseif(direc<0 && div>D_high) 
    direc = 1; 
end 
return 
 
function direc = calc_dir2(t) 
global POPSIZE DIMENSION W C1 C2 VMAX XMIN XMAX P PBEST individual 
D_low = 5e-6; 
D_high = 0.25; 
div = 0; 
x = [individual.number]; 
x = reshape(x,DIMENSION,POPSIZE); 
x = x'; 
for i=1:POPSIZE 
    div = div + (individual(i).number(t) - mean(x(:,t)))^2; 
end 
div = div / (POPSIZE * (XMAX(t)-XMIN(t))); 
if(div<D_low) 
    direc = -1; 
elseif(div>D_high) 
    direc = 1; 
elseif round(rand)  
    direc = -1; 
else  
    direc = 1; 
end 
return