ta=[ 1  1  1  2  2  3  3  4  4  4  5  5  5  5  5  6  6  6  6  6  7  7  7  10 10  15 15 8 11 11 12 12 12 12 13 13 14 14 14 16 16 16 17 17 18 19 21 22 23 23 24];
he=[ 10 19 20  10 8 8 13 9 16 25 12 13 14 17 18 19 21 22  7 20 20 22 24  15 12 11 19 11 12 13 13 14 17 19 14 17 9 16 17 17 24 25 18 23 21 20 22 23 24 25 25];
dl=[40 10 27 22 14 11  4 16 10 28 11 33 24 27 11 20 11 13 25 40 35  8  8 18 42 10 40  5 11 10 21 30 32 36 18 30 10 13 20 10 42 38 12 20  8 38 10 10 11 33 27];

g=make_graph('graf1',1,25,[ta he],[he ta]);

g=add_edge_data(g,'length',[dl dl]);

k=4;   // 4 zlecenia
baza=1
zlecenia=[2 5; 6 8; 3 4; 9 7];  // skad - dokad

// szukamy permutacji wierzchołków od 1 do 2*k pamiętajac że odbiór powinien poprzedzić dostarczenie
// jeżeli dostraczenie wyprzedzi odbiór to interpretujemy jako odbiór

ga_params = init_param();
ga_params = add_param(ga_params ,'minbound',0.00001*ones(1,k*2)); // dlatego 2*k
ga_params = add_param(ga_params ,'maxbound',ones(1,k*2));
ga_params = add_param(ga_params ,'dimension' ,k*2) ;
ga_params = add_param(ga_params , 'beta' ,2) ;
ga_params = add_param(ga_params , 'delta' ,0.1);
ga_params = add_param(ga_params , 'init_func' , init_ga_default);
ga_params = add_param(ga_params , 'crossover_func' , crossover_ga_default);
ga_params = add_param(ga_params , 'mutation_func' , mutation_ga_default) ;
ga_params = add_param(ga_params , 'codage_func' , coding_ga_identity);
ga_params = add_param(ga_params , 'selection_func' ,selection_ga_elitist) ;
ga_params = add_param(ga_params , 'nb_couples' , 110) ;
ga_params = add_param(ga_params , 'pressure' ,0.05 ) ;

function Z=get_wierzcholki(kod)
   stan_zlecen=zeros(k,1);
   wierz=[baza];
   for i=1:(2*k)
     if (kod(1,i)<=k) then
       if (stan_zlecen(kod(1,i),1)==0) then
         wierz=[wierz zlecenia(kod(1,i),1)]; // dodaj wierzcholek skad zlecenia
         stan_zlecen(kod(1,i),1)=1;    // oznacz jako odebrana
       else
         wierz=[wierz zlecenia(kod(1,i),2)]; // dodaj wierzcholek dokad zlecenia  
         stan_zlecen(kod(1,i),2)=1;    // oznacz jako odebrana  
       end    
     else 
       if (stan_zlecen(kod(1,i)-k,1)==1) then  // odebrana   
         wierz=[wierz zlecenia(kod(1,i)-k,2)]; // dodaj wierzcholek dokad zlecenia
         stan_zlecen(kod(1,i)-k,2)=1;    // oznacz jako doreczona           
       else
         wierz=[wierz zlecenia(kod(1,i)-k,1)]; // dodaj wierzcholek skad zlecenia
         stan_zlecen(kod(1,i)-k,1)=1;    // oznacz jako odebrana                      
       end
     end    
   end
   Z=[wierz baza];
endfunction

function [D,T]=oblicz_dlugosc_trasa(lista,g)
    [sx,sy]=size(lista);
    D=0;
    T=[];
    for i=1:(sy-1)
      wp=lista(1,i);
      wk=lista(1,i+1);
      [p,lp]=shortest_path(wp,wk,g,'length');
      ww=path_2_nodes(p,g);
      T=[T ww];
      D=D+lp; 
    end
endfunction


function D=oblicz_dlugosc(kod,g)
    lista=get_wierzcholki(kod);
    [sx,sy]=size(lista);
    D=0;
    for i=1:(sy-1)
      wp=lista(1,i);
      wk=lista(1,i+1);
      [p,lp]=shortest_path(wp,wk,g,'length');
      D=D+lp; 
    end
endfunction

// Fitness Function
function y=f_ocen(KK)
  os=ceil(abs(KK)*k*2);        // zamieniamy ułamki na liczby całkowite z przedziału od 1 do 2*k  
  p=get_permut(os);       // ponieważ nie było gwarancji unikalności liczb to stosujemy tkz. reprezentację porządkową Michalewicz Z. (Algorytmy genetyczne + struktury danych = programy ewolucyjne, rozdział o Zadaniu komiw                           ojażera)  by uzyskać permutację unikalnych liczb
  y=oblicz_dlugosc(p,g);  // oceniamy osobnika opisanego kolejnością p.
endfunction

PopSize = 100;
Proba_cross = 0.7;
Proba_mut = 0.5;  // było 0.1
NbGen = 10;
Log = %T;
nb_disp = 10;

[pop_opt, fobj_pop_opt, pop_init, fobj_pop_init] = optim_ga(f_ocen, PopSize, NbGen, Proba_mut, Proba_cross, Log, ga_params);

Wierz = get_wierzcholki(get_permut(ceil(abs(pop_opt(1))*2*k)));
[Dl,Tr]=oblicz_dlugosc_trasa(Wierz,g);
Wierz
Dl