alpha=0.75; // parametre de croisement
B = 10; // largeur du domaine
T=50; // nombre d'iterations
N=100; // nombre de solutions
pm=0.1; // probabilité de mutation
pc=0.7; // proba de croisement

// Affichage de la fonction
Domaine=[0:0.001:B]; length(Domaine);
Image=[]; // vecteur
for i=1:length(Domaine)
    //im = sin*sin
    im=10*sin(0.3*Domaine(i))*sin(1.3*Domaine(i)^2+0.00001*Domaine(i)^4+0.2*Domaine(i)+80)
    Image=[Image,im]
end
xset("window",1); xset("thickness",2);
plot2d(Domaine,Image,2);xtitle("La fonctin f ")

// ***Debut de l'algorithme génétique*** //
P0=grand(1,N,'unf',0,B); // tirage de N solutions de la loi uniforme entre 0 et B
f0=[];
for i=1:N // fs = les images des abcusses stoqués dans P0
    f0=[f0, 10*sin(0.3*P0(i))*sin(1.3*P0(i)^2+0.00001*P0(i)^4+0.2*P0(i)+80) ];
end

[f0,indice]=gsort(f0); // ordonner  en décroissant les valeurs de f0 et repèrer l'indice
P0=P0(indice); // réordonner selon l'indice
MeilleurI=[]; MeilleurV=[];
MeilleurI=[MeilleurI,P0(N)]; // abscisse de l'individu avec la solution minimale
MeilleurV=[MeilleurV,f0(N)]; // image de l'individu avec la solution minimale

// Croisement
for t =1:T
    Bernoulli=grand(1,N/2,'bin',1,pc) // N/2 valeurs suivant bernoulli [1,pc]
    P1=ones (1,N) // initiation d'une population de N individus
    for i=0:(N/2)-1 // par paire pour croisement!
    if Bernoulli(i+1)==1   then
        // croisement par pair qui ont la meme qualité (similarité par barycentre)
        P1(2*i+1)= alpha* P0(2*i+1)+(1-alpha)*P0(2*i+2);
        P1(2*i+2)= alpha* P0(2*i+2)+(1-alpha)*P0(2*i+1);
    else // conservation des parents
        P1(2*i+1)=P0(2*i+1)
        P1(2*i+2)=P0(2*i+2)
    end,
end

// Mutation
P2=P0;
Bernoulli=grand(1,N,'bin',1,pm); // N/2 valeurs suivant bernoulli [1,pm]
for i=1:N
    // U[-B/4;B/4] c'est le rayon de mutation
    // Faire la mutation ssi bernoulli=1 sinon P2(i)=P(i)
    P2(i)=  P0(i)+ grand(1,'unf',-B/4,B/4) *Bernoulli(i) ;
end

P3=[];
P3=[P3,P1,P2]; // P1 : croisé, P2:muté
f=[]; // contient les images de P3 pour les évaluer
for i=1:(2*N)
    f=[f, 10*sin(0.3*P3(i))*sin(1.3*P3(i)^2+0.00001*P3(i)^4+0.2*P3(i)+80) ];
end

//transformer f et P3 selon une permutaton aléatoire
permut=grand(1,'prm',[1:(2*N)]'); //2N car P3 contient P1 et P2 de tailles N
f=f(permut);
P3=P3(permut);

P4=[] //on stocke les abcisses des meilleures solutions dans P4
for i=0:(N-1)
    if  f(2*i+1)>=f(2*i+2) & grand(1,'unf',0,1)>=0.5 then // mauvaise solution
        P4=[P4,P3(2*i+1)]; // on l'accepte avec une probaba 0.5
    else
        P4=[P4,P3(2*i+2)];
    end
end

//calculer les images de la population P4
f=[]
for i=1:(N)
    f=[f,10*sin(0.3*P4(i))*sin(1.3*P4(i)^2+0.00001*P4(i)^4+0.2*P4(i)+80)];
end

[f,V]=gsort(f); // tri des valeurs de f; V : vecteur d'indices
P4=P4(V); // tri des abcisses selon le Tri de f pour prendre les meilleurs individus

// l'individu qui a fait l'image minimale : (minimiser sinus = maximiser -sinus)
// fi = - Fi + cste
MeilleurI=[MeilleurI,P4(N)]; //abcisse
MeilleurV=[MeilleurV,f(N)]; //individu (dernière valeur que je veu minimiser)
P0=P4;
end

// Affichage des meilleurs individus ayant le minimum des valeurs
disp(MeilleurV);xset("window",3); plot2d([1:1:T+1],MeilleurV,2);
xtitle("Fluctuation du minimum : ","axe du temps","")