//
// stratégie d'évolution
// variante: adaptation de la matrice de covariance
// variante: stratégie plus ou comma
// choix de fonction: Rastrigin ou Rosenbrock
//
////////////////////////////////////////////////
n=2;   // n=2 si fonction rosen
lambda=10;
mu=2;
Ngen=80;
sigma=1;
alphad=0.2;
alphaq=0.2;
//
/////////////////////////////////////////////////
function y=rast(x)
    y=sum(x.*x-cos(2*%pi*x))+length(x);
endfunction
///////////////////////////////////////////////////
/////////////////////////////////////////////////
function y=rosen(x)
    y=100*(x(2)-x(1)^2)^2+sum(x(1)-1).^2;
endfunction
///////////////////////////////////////////////////

Ap=10*rand(mu,n+1)-5*ones(mu,n+1)// initialisation
for i=1:mu
   Ap(i,n+1)=rosen(Ap(i,1:n))   // evaluation (fonction rosen ou rast) 
end
Ae=zeros(lambda,n+1)
m=1/mu*sum(Ap(:,1:n),'r');
Ye=zeros(lambda,n)
Yp=zeros(mu,n)
d=1/mu*sum(Yp,'r');
BestJ=[]
BestX=[]
Cov=eye(n,n)
//
for ng=1:Ngen
// 
  tau=0;   
  m=m+sigma*d // croisement
  val=rosen(m);    
  for i=1:lambda
     y=grand(1,'mn',zeros(n,1),Cov)
     Ye(i,:)=y';
     Ae(i,1:n)=m+sigma*y' // mutation
     Ae(i,n+1)=rosen(Ae(i,1:n))   // evaluation (fonction rosen ou rast) 
     if (Ae(i,n+1)<val) then
        tau=tau+1/lambda
     end    
  end
//  Atot=[Ap;Ae]; // selection plus (enfants+ parents)
  Atot=Ae        // selection virgule (enfants)
  [a,b]=gsort(-Atot(:,n+1))
  Atot=Atot(b,:)
  Ye=Ye(b,:)
  Ap=Atot(1:mu,:)  
  Yp=Ye(1:mu,:)
  BestJ=[BestJ,Atot(1,n+1)]
  BestX=[BestX;Atot(1,1:n)] 
  sigma=exp(1/3*(tau-1/5)/(1-1/5))*sigma
  d=1/mu*sum(Yp,'r');
  Cd=d'*d
  Cq=zeros(n,n)
  for k=1:mu
    Cq=Cq+1/mu*Yp(k,:)'*Yp(k,:)
  end  
  Cov=(1-alphad-alphaq)*Cov+alphad*Cd+alphaq*Cq
//
H=inv(Cov)
disp([H(1,1)/H(2,2),H(1,1)/H(1,2)]) //comparaison H et inv(Cov)
// valeur cible: 4,01 et 2,005
end 
figure(1)
clf() 
plot(1:Ngen,BestJ)
figure(2)
clf()
for j=1:n
 plot2d(1:Ngen,BestX(:,j),j)
end
