//
// algorithme génétique réel
//
////////////////////////////////////////////////
//
// selection par le rang
//
function Asel=rang(A)
[Npop,n1]=size(A);
[s,p]=gsort(A(:,n1));
A=A(p,:);
//
Asel=[];
p=(1:Npop)/sum(1:Npop);ps=cumsum(p);
for i=1:Npop
  u=rand();isel=1;
  while (u>ps(isel)) 
     isel=isel+1;
  end
Asel=[Asel;A(isel,:)];
end  
endfunction
///
/////////////////////////////////////////////////
// mutation
//
function Amut=mutation(A)
        [Npop,n1]=size(A);
        Amut=A;
       for i=1:Npop
         if rand()<pm then
           Amut(i,:)=A(i,:)+alpha*grand(1,n1,'nor',0,1)
         end 
        end       
endfunction
/////////////////////////////////////////////////
// croisement
//
function Acrois=croisement(A)
       [Npop,n1]=size(A);
       Acrois=A;
       for i=1:Npop/2
        u=1+floor(rand()*Npop);
        v=1+floor(rand()*Npop);
        r=rand()
        if (rand()<pc)
         Acrois(2*i-1,:)=r*A(u,:)+(1-r)*A(v,:)
         Acrois(2*i,:)=(1-r)*A(u,:)+r*A(v,:) 
        end       
      end
endfunction
/////////////////////////////////////////////////
//
function y=rastrigin(x)  // the function to optimize
n1=length(x);
y=n1+sum(x.^2-cos(2*%pi*x));
endfunction
////////////////////////////////////////////////////
Npop=60;Ngen=50;pm=0.9;pc=0.6;alpha=0.05;n=2;
////////////////////////////////////////////////////
// initialisation
A=10*rand(Npop,n+1)-5*ones(Npop,n+1)
//////////////////////////////////////////////////////
Jmin=[];minJ=1E10
for i=1:Ngen
    //eval
    for j=1:Npop
       A(j,n+1)=rastrigin(A(j,1:n)) 
    end   
    c=gsort(A(:,n+1));
    //
    Jmin=[Jmin,c(Npop)]
    if c(Npop)<minJ then bestX=A(Npop,1:n),minJ=c(Npop),end
    // selection
    Asel=rang(A)    
    // croisement
    Acrois=croisement(Asel)
    // mutation
    Amut=mutation(Acrois)
A=Amut;    
end
figure(1)
clf()
plot(1:Ngen,Jmin)
xtitle('meilleur valeur obtenue pour la fonction J')
disp('meilleur élément:');disp(A(Npop,1:n))

















