粒子群算法

  1 function z = Sphere(x)
  2     z = sum(x.^2);
  3 end
  4 
  5 
  6 
  7 function out = PSO(problem, params)
  8     %% Problem Definiton
  9     CostFunction = problem.CostFunction;  % Cost Function
 10     nVar = problem.nVar;        % Number of Unknown (Decision) Variables
 11     VarSize = [1 nVar];         % Matrix Size of Decision Variables
 12     VarMin = problem.VarMin;    % Lower Bound of Decision Variables
 13     VarMax = problem.VarMax;    % Upper Bound of Decision Variables
 14     %% Parameters of PSO
 15     MaxIt = params.MaxIt;   % Maximum Number of Iterations
 16     nPop = params.nPop;     % Population Size (Swarm Size)
 17     w = params.w;           % Intertia Coefficient
 18     wdamp = params.wdamp;   % Damping Ratio of Inertia Coefficient
 19     c1 = params.c1;         % Personal Acceleration Coefficient
 20     c2 = params.c2;         % Social Acceleration Coefficient
 21     % The Flag for Showing Iteration Information
 22     ShowIterInfo = params.ShowIterInfo;    
 23     MaxVelocity = 0.2*(VarMax-VarMin);
 24     MinVelocity = -MaxVelocity;
 25     %% Initialization
 26     % The Particle Template
 27     empty_particle.Position = [];
 28     empty_particle.Velocity = [];
 29     empty_particle.Cost = [];
 30     empty_particle.Best.Position = [];
 31     empty_particle.Best.Cost = [];
 32     % Create Population Array
 33     particle = repmat(empty_particle, nPop, 1);
 34     % Initialize Global Best
 35     GlobalBest.Cost = inf;
 36     % Initialize Population Members
 37     for i=1:nPop
 38         % Generate Random Solution
 39         particle(i).Position = unifrnd(VarMin, VarMax, VarSize);
 40         % Initialize Velocity
 41         particle(i).Velocity = zeros(VarSize);
 42         % Evaluation
 43         particle(i).Cost = CostFunction(particle(i).Position);
 44         % Update the Personal Best
 45         particle(i).Best.Position = particle(i).Position;
 46         particle(i).Best.Cost = particle(i).Cost;
 47         % Update Global Best
 48         if particle(i).Best.Cost < GlobalBest.Cost
 49             GlobalBest = particle(i).Best;
 50         end
 51     end
 52     % Array to Hold Best Cost Value on Each Iteration
 53     BestCosts = zeros(MaxIt, 1);
 54     %% Main Loop of PSO
 55     for it=1:MaxIt
 56         for i=1:nPop
 57             % Update Velocity
 58             particle(i).Velocity = w*particle(i).Velocity ...
 59                 + c1*rand(VarSize).*(particle(i).Best.Position - particle(i).Position) ...
 60                 + c2*rand(VarSize).*(GlobalBest.Position - particle(i).Position);
 61             % Apply Velocity Limits
 62             particle(i).Velocity = max(particle(i).Velocity, MinVelocity);
 63             particle(i).Velocity = min(particle(i).Velocity, MaxVelocity);
 64           
 65             % Update Position
 66             particle(i).Position = particle(i).Position + particle(i).Velocity;
 67             % Apply Lower and Upper Bound Limits
 68             particle(i).Position = max(particle(i).Position, VarMin);
 69             particle(i).Position = min(particle(i).Position, VarMax);
 70             % Evaluation
 71             particle(i).Cost = CostFunction(particle(i).Position);
 72             % Update Personal Best
 73             if particle(i).Cost < particle(i).Best.Cost
 74                 particle(i).Best.Position = particle(i).Position;
 75                 particle(i).Best.Cost = particle(i).Cost;
 76                 % Update Global Best
 77                 if particle(i).Best.Cost < GlobalBest.Cost
 78                     GlobalBest = particle(i).Best;
 79                 end            
 80             end
 81         end
 82         % Store the Best Cost Value
 83         BestCosts(it) = GlobalBest.Cost;
 84         % Display Iteration Information
 85         if ShowIterInfo
 86             disp([‘Iteration ‘ num2str(it) ‘: Best Cost = ‘ num2str(BestCosts(it))]);
 87         end
 88         % Damping Inertia Coefficient
 89         w = w * wdamp;
 90     end
 91     
 92     out.pop = particle;
 93     out.BestSol = GlobalBest;
 94     out.BestCosts = BestCosts;
 95     
 96 end
 97 
 98 
 99 %% Problem Definiton
100 problem.CostFunction = @(x) Sphere(x);  % Cost Function
101 problem.nVar = 5;       % Number of Unknown (Decision) Variables
102 problem.VarMin =  -10;  % Lower Bound of Decision Variables
103 problem.VarMax =  10;   % Upper Bound of Decision Variables
104  
105 %% Parameters of PSO
106 params.MaxIt = 1000;        % Maximum Number of Iterations
107 params.nPop = 50;           % Population Size (Swarm Size)
108 params.w = 1;               % Intertia Coefficient
109 params.wdamp = 0.99;        % Damping Ratio of Inertia Coefficient
110 params.c1 = 2;              % Personal Acceleration Coefficient
111 params.c2 = 2;              % Social Acceleration Coefficient
112 params.ShowIterInfo = true; % Flag for Showing Iteration Informatin
113  
114 %% Calling PSO
115 out = PSO(problem, params);
116  
117 BestSol = out.BestSol;
118 BestCosts = out.BestCosts;
119  
120 %% Results
121 figure;
122 % plot(BestCosts, ‘LineWidth‘, 2);
123 semilogy(BestCosts, ‘LineWidth‘, 2);
124 xlabel(‘Iteration‘);
125 ylabel(‘Best Cost‘);
126 grid on;