pso1.m

一个很好的pso程序 粒子群优化法的一个好例子

%Initialization of PSO parameters
wmax=0.9;
wmin=0.4;
itmax=200; %Maximum iteration number
c1=1.4;
c2=1.4;

for iter=1:itmax
    W(iter)=wmax-((wmax-wmin)/itmax)*iter;  % W[100]啊吝摹 痢痢 累酒咙
end

%**********************************************************

%Initialization of positions of agents
% agents are initialized between -5,+5 randomly
a=0;  
b=5;
N=20;
D=2;

x=a+(b-a)*rand(N,D,1);  % x[20,2,2] 檬扁拳

%Initialization of velocities of agents
%Between -5 , +5, (which can also be started from zero)

m=0;
n=1;
V=m+(n-m)*rand(N,D,1);  % V[20,2,2] 檬扁拳

%**********************************************************
%Function to be minimized. Scheaffer F6
%x(i,2,1)=x2    % y axis
%x(i,1,1)=x1    % x axis
%f=100*(x2-x1^2)^2+(x1-1)^2
for i=1:N;
    F(i,1,1)=100*(x(i,2,1)-x(i,1,1)^2)^2+(x(i,1,1)-1)^2;
end
%**********************************************************
[C,I]=min(abs(F(:,1,1)));

B(1,1,1)=C;
XX(1,1,1)=I;
gbest(1,1,1)=x(I,1,1);  % x axis
gbest(1,2,1)=x(I,2,1);  % y axis

%********************************************************
%Matrix composed of gbest vector 
for p=1:N
    for r=1:D  
        G(p,r,1)=gbest(1,r,1);
    end
end
 
Fbest(1,1,1)=100*(G(1,2,1)-G(1,1,1)^2)^2+(G(1,1,1)-1)^2;

for i=1:N;
    pbest(i,:,1)=x(i,:,1);
end
    

V(:,:,2)=W(1)*V(:,:,1)+c1*rand*(pbest(:,:,1)-x(:,:,1))+c2*rand*(G(:,:,1)-x(:,:,1));

x(:,:,2)=x(:,:,1)+V(:,:,2);

 %       plot(x(:,1,1),x(:,2,1),'o') % Particle狼 檬扁 困摹 Display
%clear
        %getframe;              % 弊覆 capture

Fb(1,1,1)=100*(gbest(1,2,1)-gbest(1,1,1)^2)^2+(gbest(1,1,1)-1)^2;

%******************************************************
for j=2:itmax-1
% Calculation of new positions
    for i=1:N;
        F(i,1,j)=100*(x(i,2,j)-x(i,1,j)^2)^2+(x(i,1,j)-1)^2;
    end
    [C,I]=min(abs(F(:,:,j)));   %C绰 弥家蔼, I绰 割锅掳啊 弥家蔼牢瘤 唱鸥辰促
    B(1,1,j)=C;
    gbest(1,1,j)=x(I,1,j);
    gbest(1,2,j)=x(I,2,j);
        
    Fb(1,1,j)=100*(gbest(1,2,j)-gbest(1,1,j)^2)^2+(gbest(1,1,j)-1)^2;

    [C,I]=min(Fb(1,1,:));

    if Fb(1,1,j)<=C
        gbest(1,1,j)=gbest(1,1,j);
        gbest(1,2,j)=gbest(1,2,j);
    else
        gbest(1,1,j)=gbest(1,1,I);
        gbest(1,2,j)=gbest(1,2,I);
    end

    %Matrix composed of gbest vector 
    for p=1:N
        for r=1:D  
            G(p,r,j)=gbest(1,r,j);
        end
    end

    Fbest(1,1,j)=100*(G(1,2,j)-G(1,1,j)^2)^2+(G(1,1,j)-1)^2;

 
    for i=1:N;
        [C,I]=min(F(i,1,:));
        if F(i,1,j)<=C
            pbest(i,:,j)=x(i,:,j);
        else
            pbest(i,:,j)=x(i,:,I);
        end
    end

    V(:,:,j+1)=W(j)*V(:,:,j)+c1*rand*(pbest(:,:,j)-x(:,:,j))+c2*rand*(G(:,:,j)-x(:,:,j));

    x(:,:,j+1)=x(:,:,j)+V(:,:,j+1);
    
end

%plot(x(:,1,200),x(:,2,200),'o',x(:,1,190),x(:,2,190),'x')
 
for i=1:200;
    plot(x(:,1,i),x(:,2,i),'o')
    axis([0 5 0 5])
    getframe;              % 弊覆 capture
    naver_delay(0.01);
end
clear