function [bestX, bestFitness, bestFitnessEvolution]= ICEFO_v1(options)
%--------------------------------------------------------------------------
% Improved Chaotic Electromagnetic Field Optimization (ICEFO)
% Dr Houssem BOUCHEKARA
% 26/07/2019
%--------------------------------------------------------------------------
% For the initial version of the EFO:
% Abedinpourshotorban, H., Shamsuddin, S. M., Beheshti, Z., & Jawawi, D. N.
% (2015). Electromagnetic field optimization: A physics-inspired
% metaheuristic optimization algorithm. Swarm and Evolutionary
% Computation.
%
% For the Improved Chaotic Electromagnetic Field Optimization (ICEFO)
% algorithm
% Bouchekara, H. Neural Comput & Applic (2019).
% https://doi.org/10.1007/s00521-019-04298-3
%--------------------------------------------------------------------------
bestFitnessEvolution=[];
%--------------------------------------------------------------------------
N_var = options.N_var ;
N_emp = options.N_emp ;
Max_gen = options.Max_gen ;
minval = options.minval ;
maxval = options.maxval ;
R_rate = options.R_rate ;
Ps_rate = options.Ps_rate ;
P_field = options.P_field ;
N_field = options.N_field ;
ObjFunction = options.ObjFunction ;
phi = (1 + sqrt(5))/2; %golden ratio
%initializatin
em_pop=repmat(minval,N_emp,1)+rand(N_emp,N_var).*repmat((maxval - minval),N_emp,1); % random population
fit= feval(ObjFunction,em_pop(1:N_emp,1:N_var));
em_pop = [em_pop(:,1:N_var) fit'];
em_pop = sortpop(em_pop,N_var+1);
em_pop1=em_pop;
em_pop2=em_pop;
em_pop3=em_pop;
em_pop4=em_pop;
%random vectors (this is to increase the calculation speed instead of
%determining the random values in each iteration we allocate them in the
%beginning before algorithm start
r_index1 = randi([1 (round(N_emp.* P_field))],[N_var Max_gen]); %random particles from positive field
r_index2 = randi([(round(N_emp.*(1-N_field))) N_emp],[N_var Max_gen]); %random particles from negative field
r_index3 = randi([(round(N_emp.* P_field)+1) (round(N_emp.*(1-N_field))-1)],[N_var Max_gen]);%random particles from neutral field
ps= rand(N_var,Max_gen);%= probability of selecting electromagnets of generated particle from the positive field
r_force = rand(1,Max_gen);%random force in each generation
%--------------------------------------------------------------------------
% Chaos
% 1: 'Chebyshev map'
% 2: 'Circle map'
% 3: 'Gauss/mouse map'
% 4: 'Iterative map'
% 5: 'Logistic map'
% 6: 'Piecewise map'
% 7: 'Sine map'
% 8: 'Singer map'
% 9: 'Sinusoidal map'
% 10: 'Tent map'
Index_c=2;
r_force_c1=chaos(1,0.7,Max_gen);
r_force_c2=chaos(2,0.7,Max_gen);
r_force_c3=chaos(4,0.7,Max_gen);
r_force_c4=chaos(5,0.7,Max_gen);
r_force1=r_force_c1;
r_force2=r_force_c2;
r_force3=r_force_c3;
r_force4=r_force_c4;
%--------------------------------------------------------------------------
rp = rand(1,Max_gen);%some random numbers for checking randomness probability in each generation
randomization = rand(1,Max_gen);%coefficient of randomization when generated electro magnet is out of boundary
RI=1;%index of the electromagnet (variable) which is going to be initialized by random number
generation=N_emp;
new_emp = zeros(1,N_var+1); %temporary array to store generated particle
ITER=1;
counter=0;
while (generation <= Max_gen)
r1 = r_force1(1,generation);
r2 = r_force2(1,generation);
r3 = r_force3(1,generation);
r4 = r_force4(1,generation);
% 1------------------------------------------------------------------------------------------------
for i=1:N_var
if (ps(i,generation) > Ps_rate)
new_emp1(i) = em_pop1(r_index3(i,generation), i) +...
phi * r1 * (em_pop1(r_index1(i,generation), i) - em_pop1(r_index3(i,generation), i)) +...
r1 * (em_pop1(r_index3(i,generation), i) - em_pop1(r_index2(i,generation), i));
else
new_emp1(i) = em_pop1 (r_index1(i,generation), i);
end
%checking whether the generated number is inside boundary or not
if ( new_emp1(i) >= maxval(i) || new_emp1(i) <= minval(i) )
new_emp1(i) = minval(i) + (maxval(i) - minval(i)) .* randomization(1, generation);
end
end
%replacement of one electromagnet of generated particle with a random number (only for
%some generated particles) to bring diversity to the population
if ( rp(1,generation) < R_rate)
new_emp1(RI) = minval(RI) + (maxval(RI) - minval(RI)) .* randomization(1, generation);
RI=RI+1;
if (RI > N_var)
RI=1;
end
end
new_emp1(N_var+1) = feval(ObjFunction,new_emp1(1:N_var));
%updating the population if the fitness of the generated particle is better than worst fitness in
%the population (because the population is sorted by fitness, the last particle is the worst)
if ( new_emp1(N_var+1) < em_pop1(N_emp , N_var+1) )
position=find(em_pop1(:,N_var+1) > new_emp1(N_var+1));
em_pop1=insert_in_pop(em_pop1,new_emp1,position(1));
end
%-------------------------------------------------------------------------------------------------
% 2-----------------------------------------------------------------------------------------------
for i=1:N_var
if (ps(i,generation) > Ps_rate)
new_emp2(i) = em_pop2(r_index3(i,generation), i) +...
phi * r2 * (em_pop2(r_index1(i,generation), i) - em_pop2(r_index3(i,generation), i)) +...
r2 * (em_pop2(r_index3(i,generation), i) - em_pop2(r_index2(i,generation), i));
else
new_emp2(i) = em_pop2 (r_index1(i,generation), i);
end
%checking whether the generated number is inside boundary or not
if ( new_emp2(i) >= maxval(i) || new_emp2(i) <= minval(i) )
new_emp2(i) = minval(i) + (maxval(i) - minval(i)) .* randomization(1, generation);
end
end
%replacement of one electromagnet of generated particle with a random number (only for
%some generated particles) to bring diversity to the population
if ( rp(1,generation) < R_rate)
new_emp2(RI) = minval(RI) + (maxval(RI) - minval(RI)) .* randomization(1, generation);
RI=RI+1;
if (RI > N_var)
RI=1;
end
end
new_emp2(N_var+1) = feval(ObjFunction,new_emp2(1:N_var));
%updating the population if the fitness of the generated particle is better than worst fitness in
%the population (because the population is sorted by fitness, the last particle is the worst)
if ( new_emp2(N_var+1) < em_pop2(N_emp , N_var+1) )
position=find(em_pop2(:,N_var+1) > new_emp2(N_var+1));
em_pop2=insert_in_pop(em_pop2,new_emp2,position(1));
end
%-------------------------------------------------------------------------------------------------
% 3-----------------------------------------------------------------------------------------------
for i=1:N_var
if (ps(i,generation) > Ps_rate)
new_emp3(i) = em_pop3(r_index3(i,generation), i) +...
phi * r3 * (em_pop3(r_index1(i,generation), i) - em_pop3(r_index3(i,generation), i)) +...
r3 * (em_pop3(r_index3(i,generation), i) - em_pop3(r_index2(i,generation), i));
else
new_emp3(i) = em_pop3 (r_index1(i,generation), i);
end
%checking whether the generated number is inside boundary or not
if ( new_emp3(i) >= maxval(i) || new_emp3(i) <= minval(i) )
new_emp3(i) = minval(i) + (maxval(i) - minval(i)) .* randomization(1, generation);
end
end
%replacement of one electromagnet of generated particle with a random number (only for
%some generated particles) to bring diversity to the population
if ( rp(1,generation) < R_rate)