clc;
clear all;
tic;
format long;
format compact;

'PSO'

% Choose the problems to be tested. Please note that for test functions F7
% and F25, the global optima are out of the initialization range. For these
% two test functions, we do not need to judge whether the variable violates
% the boundaries during the evolution after the initialization.
% problemSet = [1 :30];
problemSet = [1:30];
for problemIndex = 1:30

    problem = problemSet(problemIndex)
    
    % Record the optimum value.
    optimum = problem*100;
    
    % Define the dimension of the problem
    n = 30;
    lu = [-100 * ones(1, n); 100 * ones(1, n)];

    % Record the best results
    outcome = [];
    outcome_fitness=[];
    
    % set the popsize of the swarm
    if n == 30
      popsize = 40;
    else
      popsize = 60;
    end
    
    % Main body
    time = 1;

    % The total number of runs
    totalTime = 1;
    
    LB = repmat(lu(1,:), popsize, 1);
    UB = repmat(lu(2,:), popsize, 1);


    while time <= totalTime

        rand('seed', sum(100*clock));

        % Initialize the main population
        p = repmat(lu(1, :), popsize, 1) + rand(popsize, n) .* (repmat(lu(2, :) - lu(1, :), popsize, 1));

        % Evaluate the objective function values
        fit =   cec14_func(p', problem);
        fit = fit';
       
        % Set the maximum FES
        FES_max =  n * 10000;    
        
        % Record the history best particle
        p_h_best = p;
        fit_h_best = fit;
        
        % Record the number of function evaluations (FES)
        FES = popsize;
        
        % Initialize the velocity vector
        v = repmat(lu(1, :),  popsize, 1) + rand(popsize, n) .* (repmat(lu(2, :) - lu(1, :), popsize, 1));
        
        while FES < FES_max  
            % PSO-cf
              c1=1.4962;            
              c2=1.4962;  
              w=0.7298;

            % PSO-w
%             c1 = 2.00;
%             c2 = 2.00;
%             w = 0.9  - 0.5 * FES./FES_max;   
%                           
            [ ~, index] = min(fit_h_best);
            
            % Compute the global history best p_best.
            global_best = p_h_best(index,:);
              
            [ ~, indexy] = sort(fit);

            % velocity renew
             for i = 1 : popsize  
                v( i,:) = w .* v(i,:) + c1.*(p_h_best(i,:) - p(i,:)).*rand(1,n) + c2.* (global_best - p(i,:)).*rand(1,n);
             end
            
            % repair the velocity 
            BL  = v <LB; v(BL) =LB(BL);   BU  = v>UB; v(BU) = UB(BU);
             
            p_new = [];
            
            % Renew the particles
             for i = 1 : popsize
               p_new(i ,:) = p(i,:) + v(i,:);
             end
             
            % Repair the new particles
            BL  = p_new <LB;  p_new(BL) = LB(BL); 
            BU  = p_new>UB;  p_new(BU) = UB(BU);
            
            % Evaluate the trial vectors
           fitSet =   cec14_func(p_new', problem);
           fitSet = fitSet';
             
           FES = FES + popsize;
            
           for i = 1 : popsize
                if fitSet(i,:) < fit_h_best(i,:);
                   p_h_best(i,:) = p_new(i,:);
                   fit_h_best(i,:) = fitSet(i,:);
               end
           end   
           p = p_new;
           fit = fitSet;    
          
       end
        outcome = [outcome min(fit_h_best-optimum)];
        time = time + 1;

    end

    sort(outcome)
    mean(outcome)
    std(outcome)
end
toc;