%**************************************************************************************************
%Reference:  1) J. Brest, S. Greiner, B. Boskovic, M. Mernik, and V. Zumer, ˇ°Self-adapting
%                     control parameters in differential evolution: A comparative study on numerical
%                     benchmark problems,ˇ± IEEE Trans. Evolut. Comput., vol. 10, no. 6,
%                     pp. 646¨C657, Dec. 2006.
%                     2) J. Zhang and A. C. Sanderson, ˇ°JADE: adaptive differential evolution
%                     with optional external archive,ˇ± IEEE Trans. Evolut. Comput., vol. 13,
%                     no. 5, pp. 945-958, 2009.
%
% Note: We obtained the MATLAB source code from Dr.Jane Jing Liang, and did some
%           minor revisions in order to solve the 30 cec2014 benchmark test functions,
%           however, the main body was not changed.
%**************************************************************************************************
warning off;
clear all
clc;
tic;

'cec2014'
'jDE'

% Set the dimension of the problem
n = 30;
% Set the population size
popsize = 100;
% Set two parameters in jDE
tau1 = 0.1;
tau2 = 0.1;

% Define the search space
lu = [-100 * ones(1, n); 100 * ones(1, n)];

problemSet = [1 : 30];

for problemIndex= 1 : 30
    
    problem = problemSet(problemIndex)
    
    % Set the global optimum
    optimum = problem * 100.0;
    
    % record the best results
    outcome = [];
    
    %Main body which was provided by Dr. J. Zhang
    time = 1;
    
    % The total number of runs
    totalTime = 5;
    
    while time <= totalTime
        
        rand('seed', sum(100 * clock));
        
        % Initialize the main population
        popold = repmat(lu(1, :), popsize, 1) + rand(popsize, n) .* (repmat(lu(2, :)-lu(1, :), popsize, 1));
        
        % Evaluate the initial population
        valParents = cec14_func(popold', problem);
        valParents = valParents';
        
        % The settings of F and CR
        F = 0.5 * ones(popsize, 1);
        CR = 0.9 * ones(popsize, 1);
        
        FES = popsize;
        
        while FES < n * 10000
            
            pop = popold;      % the old population becomes the current population
            
            Fold = F;
            CRold = CR;
            
            % Generate F and CR
            IF  = rand(popsize, 1) < tau1;
            ICR = rand(popsize, 1) < tau2;
            
            F(IF)  = 0.1 + 0.9 * rand(sum(IF), 1);
            CR(ICR) = 0.0 + 1.0 * rand(sum(ICR), 1);
            
            r0 = [1:popsize];
            [r1, r2, r3] = gnR1R2R3(popsize, r0);
            
            % == == == == == == == == ==Mutation == == == == == == == == ==
            vi  = pop(r1, :) + F(:, ones(1, n)) .* (pop(r2, :) - pop(r3, :));
            
            vi = boundConstraint(vi, lu);
            
            %== == == == == == == == == Crossover == == == == == == == == =
            
            mask = rand(popsize, n) > CR(:, ones(1, n));                   % mask is used to indicate which elements of ui comes from the parent
            rows = (1:popsize)'; cols = floor(rand(popsize, 1) * n) + 1;   % choose one position where the element of ui doesn't come from the parent
            jrand = sub2ind([popsize n], rows, cols); mask(jrand) = false;
            ui = vi;  ui(mask) = pop(mask);
            
            % Evaluate the trial vectors
            valOffspring = cec14_func(ui', problem);
            valOffspring = valOffspring';
            
            FES = FES + popsize;
            
            %  == == == == == == == == == Selection == == == == == == == ==
            % I == 1: the parent is better; I == 2: the offspring is better
            [valParents, I] = min([valParents, valOffspring], [], 2);
            popold = pop;
            
            popold(I == 2, :) = ui(I == 2, :);
            
            F(I == 1) = Fold(I == 1);
            CR(I == 1) = CRold(I == 1);
            
        end
        
        outcome = [outcome min(valParents-optimum)];
        
        time = time + 1;
        
    end
    
    sort(outcome)
    mean(outcome)
    std(outcome)
    
end
toc;





