* fichier : ordo_2.dgibi

************************************************************************
*  TEST DE L OPERATEUR ORDO 'COUT'
*  POUR LE CALCUL DE LA PERMUTATION OPTIMISANT UN COUT
*  BP, 2016-06-24
*  mot-cles : mathematiques, permutation, arrangement
************************************************************************



************************************************************************
* petite procedure utile pour le calcul a la main
************************************************************************
debp coutperm lili*'LISTENTI' jperm*'LISTENTI';
  Scout = 0;
  repe bi (dime jperm);
    j = extr jperm &bi;
    k = ((&bi - 1) * n) + j;
    Scout = Scout + (ABS (extr lili k));
  fin bi;
finp Scout;


************************************************************************
* donnees
************************************************************************

* ml1 = lect 
*    6   3   7   1   8
*    5   2   4   4   7
*    2   5   3   9   4
*    6   3   4   5   4
*    2   2   6   5   8 ;
ml1 = brui 'BLAN' 'POIS' 1000 (10**2);
list ml1;



************************************************************************
* on verifie a la main sur toutes les permutations
************************************************************************

n2 = dime ml1;
n  = enti proche (n2**0.5);

nfacto = factorie n;
mess 'n!=' nfacto;


si (nfacto < 1000) ;


Tperm = table;
Tcout = lect; 
* combinaison initiale
perm1 = lect 1 pas 1 n;
cout1 =  coutperm ml1 perm1;
Tperm . 1 = (perm1 + 0);
Tcout = Tcout et cout1;
mess '----- Combinaison 1 de cout = '   cout1   '-----';
list perm1;
mess '--------------------------------------------------' ;

* boucle sur les combinaisons possibles
repe bcomb (nfacto - 1);
      icomb = &bcomb + 1;

*     i=n-1
      i = n - 1;
      
*  10 if(a(i).lt.a(i+1)) go to 20
*     i=i-1
*     if(i.eq.0) go to 20
*     go to 10
      repe b10;
        si ((extr perm1 i) < (extr perm1 (i+1))); quit b10; finsi;
        i = i - 1;
        si (i ega 0); quit b10; finsi;
      fin  b10;
      
*  20 j=i+1
*     k=n
      j = i + 1;
      k = n;
      
*    30 t=a(j)
*       a(j)=a(k)
*       a(k)=t
*       j=j+1
*       k=k-1
*       if(j.lt.k) go to 30
*       j=i
*       if(j.ne.0) go to 40
*       nextp=.false.
*       return
      repe b30;
*       swap
        t = extr perm1 j;
        REMPLACER perm1 j (extr perm1 k);
        REMPLACER perm1 k t;
        j = j + 1;
        k = k - 1;
        si (j < k); iter b30; finsi;
        j = i;
        si (j neg 0); quit b30; finsi;        
        si (j ega 0); quit bcomb; finsi;
      fin b30;
      
*    40 j=j+1
*       if(a(j).lt.a(i)) go to 40
*       t=a(i)
*       a(i)=a(j)
*       a(j)=t
*       nextp=.true.
*       end
      repe b40;
        j = j + 1;
        si ((extr perm1 j) < (extr perm1 i)); iter b40; finsi;
*       swap
        t = extr perm1 i;
        REMPLACER perm1 i (extr perm1 j);
        REMPLACER perm1 j t;
        quit b40;
      fin  b40;
      
      cout1 =  coutperm ml1 perm1;
      Tperm . icomb = (perm1 + 0);
      Tcout = Tcout et cout1;
      mess '----- Combinaison 'icomb' de cout = '   cout1   '-----';
      list perm1;
      mess '--------------------------------------------------' ;
      
        
fin  bcomb;

jlist = lect 1 pas 1 nfacto;
Tcout2 jlist2 = ordo Tcout jlist;
jmin = extr jlist2 1;
cout2 = extr Tcout2 1;
mess '----- Combinaison ' jmin ' de cout mini = ' cout2   '-----';
list Tperm . jmin;

finsi;
   
************************************************************************
* calcul via ORDO 'COUT' LISTENTI
************************************************************************

* calcul de la permutation et du cout associe

*------> methode Hongroise :
jperm = lect 1 pas 1 n;
temp zero;
c_ordo p_ordo = ORDO 'COUT' 'HONG' ml1 jperm;
temp 'SGAC' 'IMPR';
temp ;
mess '>>>>> Cout mini =' c_ordo '<<<<<';
list p_ordo;

*------> methode ou lon calcule tout :
temp zero;
c_ordo2 p_ordo2 = ORDO 'COUT'  ml1 jperm;
temp 'SGAC' 'IMPR';
temp ;
mess '>>>>> Cout mini =' c_ordo2 '<<<<<';
list p_ordo2;

************************************************************************
* calcul via ORDO 'COUT' LISTREEL
************************************************************************

*------> methode Hongroise :
Xperm = lect 1 pas 1 n; list Xperm;
XML1 = FLOT ml1; list XML1;
temp zero;
* opti surv 2121211;
c_ordX p_ordX = ORDO 'COUT' 'HONG' Xml1 Xperm;
temp 'SGAC' 'IMPR';
temp ;
mess '>>>>> Cout mini =' c_ordX '<<<<<';
list c_ordX;
toto = p_ordX + 1;
list toto;
list p_ordX;

*------> methode ou lon calcule tout :
temp zero;
c_ordX2 p_ordX2 = ORDO 'COUT'  Xml1 Xperm;
temp 'SGAC' 'IMPR';
temp ;
mess '>>>>> Cout mini =' c_ordX2 '<<<<<';
list p_ordX2;



************************************************************************
* test de non regression
************************************************************************

si (nfacto < 1000) ;
  errcout = abs (cout2 - c_ordo);
  errperm = maxi (Tperm . jmin - p_ordo) 'ABS';
sinon;
  errcout = abs (c_ordo2 - c_ordo);
  errperm = maxi (p_ordo2 - p_ordo) 'ABS';
finsi;

SI ((errcout > 0) OU (errperm > 0)) ;
   ERRE 5 ;
SINO ;
   ERRE 0 ;
FINSI ;

FIN ;