#
# Programme developpe par Guy Tremblay, Professeur au dept. d'informatique
# dans le cadre du cours INF7440.
# (Automne 2002, 2003, 2004)
# Modification par Cedric Chauve (Automne 2004)


resource TSP_profondeur()
  const int INFINI = high(int);

  #############################################################
  # Matrices des couts pour le graphe.
  #############################################################

  # Lecture et verification du nombre d'elements a generer et a traiter
  int n; getarg(1, n); 

  if (n <= 0) {
    write("*** Erreur: n <= 0" );
    stop(1);
  }

  int germe; getarg(2, germe); 
  if (germe <= 0) {
    write("*** Erreur: germe <= 0" );
    stop(1);
  }

  int W[n][n];

  # Generation (aleatoire) des elements du tableau a.
  procedure generer( ref int a[*,*], int n ) 
  {
    seed(1.0 * germe); 
    for [i = 1 to n, j = 1 to n st i != j] {
      a[i][j] = int(random(1, 10*n));
    }
    for [i = 1 to n] {
      W[i][i] = 0;
    }
  }

  # Impression des elements d'un tableau.
  procedure imprimerCouts( int a[*][*], int n ) {
    for [i = 1 to n] {
      for [j = 1 to n-1] {
        writes( a[i][j], ", " );
      }
      write( a[i][n] );
    }
  }

  generer( W, n );
  write( "W = " );
  imprimerCouts( W, n );
  write();

  #############################################################
  # Types et procedures pour les sequences.
  #############################################################

  const int N = 100;
  type Sequence    = ptr SequenceRec;
  type SequenceRec = rec( int elems[N]; int ln; );

  procedure creer() returns Sequence s
  # POSTCONDITION
  #   s = []
  { s = new(SequenceRec); s^.ln = 0; }

  procedure cloner( Sequence s1 ) returns Sequence s2
  # POSTCONDITION
  #   s2 est un nouvel objet qui est une copie (profonde) de s1.
  { s2 = new(SequenceRec); s2^ = s1^; }

  procedure longueur( Sequence s ) returns int l
  # POSTCONDITION
  #   l = length(s)
  { l = s^.ln; }
 
  procedure ajouterEnQueue( Sequence s, int elem )
  # POSTCONDITION
  #   s = s' || [elem]
  { s^.ln += 1; s^.elems[s^.ln] = elem; }
  
  procedure supprimerQueue( Sequence s )
  # PRECONDITION
  #   s != []
  # POSTCONDITION
  #   s = s'[1..length(s')-1]
  { s^.ln -= 1; }
  
  procedure tete( Sequence s ) returns int t
  # PRECONDITION
  #   s != []
  # POSTCONDITION
  #   t = s[1]
  { t = s^.elems[1]; }

  procedure queue( Sequence s ) returns int q
  # PRECONDITION
  #   s != []
  # POSTCONDITION
  #   q = s[length(s)]
  { q = s^.elems[s^.ln]; }

  procedure element( Sequence s, int i ) returns int e
  # PRECONDITION
  #   i IN domain(s)
  # POSTCONDITION
  #   e = s[i]
  { e = s^.elems[i]; }

  procedure estElement( Sequence s, int e ) returns bool r 
  # POSTCONDITION
  #   r <=> e IN s
  { 
    r = false;
    for [i = 1 to s^.ln] {
      if( s^.elems[i] == e ) { r = true; return; }
    }
  }

  procedure imprimer( Sequence s ) 
  # POSTCONDITION
  #   Chacun des elements de s a ete imprime sur stdout.
  {
    writes("[");
    for [i = 1 to s^.ln-1] {
        writes( s^.elems[i], ", " );
    }
    writes( s^.elems[s^.ln], "]" );
  }

  #############################################################
  # Types et procedures pour les files de priorite.
  #############################################################

  const int MAXELEMS   = 5000;
  type ElemFP          = rec( int priorite; Sequence elem );
  type FilePriorite    = ptr FilePrioriteRec;
  type FilePrioriteRec = rec( ElemFP elements[MAXELEMS]; int nbElements );

  procedure creerFP() returns FilePriorite fp
  # POSTCONDITION
  #   fp = {}
  { fp = new(FilePrioriteRec); fp^.nbElements = 0; }

  procedure ajouter( FilePriorite fp, int priorite, Sequence elem )
  # POSTCONDITION 
  #  fp = fp' U {(priorite, elem)}
  { 
    fp^.nbElements += 1;
    fp^.elements[fp^.nbElements].priorite = priorite;
    fp^.elements[fp^.nbElements].elem  = elem;
  }

  procedure estVide( FilePriorite fp ) returns bool r
  # POSTCONDITION 
  #   r <=> fp = {}
  { r = fp^.nbElements == 0; }

  procedure obtenirMin( FilePriorite fp, ref int prio, ref Sequence elem )
  # PRECONDITION
  #   ~estVide(fp)
  # POSTCONDITION
  #   (prio, elem) IN fp' & prio = MINIMUM( p :: (p, e) IN fp' )
  #   fp = fp' - {(prio, elem)}
  {
    prio = high(int);
    int posMin;
    for [i = 1 to fp^.nbElements] {
      if (fp^.elements[i].priorite < prio) {
        prio = fp^.elements[i].priorite;
        posMin = i;
      }
    }
    elem = fp^.elements[posMin].elem;
    fp^.elements[posMin] = fp^.elements[fp^.nbElements];
    fp^.nbElements -= 1;
  }

  #############################################################
  # Procedures pour le parcours du graphe. 
  #############################################################

  int nbCircuitsComplets = 0;
  int nbNoeudsParcourus  = 0;

  procedure areteExiste( int i, int j ) returns bool r 
  # ROLE
  #   Determine si une arete existe du sommet i vers le sommet j dans le graphe.
  # PRECONDITION
  #   i et j denote des numeros de sommets valides pour la matrice W.
  # POSTCONDITION
  #   r <=> W[i][j] != INFINI
  { r = W[i][j] != INFINI; }


  procedure dejaVisite( Sequence chemin, int sommet ) returns bool r
  # ROLE
  #   Determine si le sommet a deja ete visite, c'est-a-dire, 
  #   appartient deja au chemin indique.  
  # POSTCONDITION 
  #   r <=> estElement(chemin, sommmet)
  { r = estElement(chemin, sommet); nbNoeudsParcourus++; }


  procedure coutChemin( Sequence chemin ) returns int cout
  # POSTCONDITION 
  #   cout = sommme des couts des aretes reliant les sommets du chemin indique.
  {
    cout = 0;
    for [i = 1 to longueur(chemin)-1] {
      cout += W[element(chemin, i)][element(chemin, i+1)];
    }  
  }

  procedure coutMin( int couts[*], Sequence chemin, int aExclure ) returns int cout
  # POSTCONDITION
  #   cout = minimum des couts indiques, mais en ignorant le sommet aExclure
  #          de meme que les sommets deja visites de chemin (mais a l'exception 
  #          du sommet de depart, pour creer le cycle).
  {
    cout = high(int);
    for [i = 1 to ub(couts) st couts[i] > 0  
                            & i != aExclure 
                            & (~dejaVisite(chemin, i) | element(chemin, 1) == i) ] {
      cout = min(couts[i], cout);
    }
  }

  procedure coutEstime( Sequence chemin, int sommet, int n ) returns int cout
  # ROLE
  #   Determiner un estime (une borne inferieure) pour le cout d'un circuit
  #   obtenu en ajoutant le sommet indique au chemin existant.
  # POSTCONDITION
  #   cout = somme des couts du chemin + somme de l'arete reliant le dernier sommet
  #          du chemin a sommet + somme des couts minimum pour 
  #          aller vers les sommets pas visites puis vers le sommet de depart.
  #
  {
    cout = coutChemin(chemin);
    cout += W[queue(chemin)][sommet];
    cout += coutMin( W[sommet][*], chemin, queue(chemin) );
    for [s = 1 to n st s != sommet & ~dejaVisite(chemin, s)] {
      cout += coutMin( W[s][*], chemin, sommet );
    }
  }

  #########################################################################
  # Variables (globales) et procedures pour mise a jour du chemin minimum.
  #########################################################################

  Sequence cheminMin;
  int coutCheminMin = INFINI;

  procedure majCheminMin( Sequence chemin, int cout )
  {
    if (cout < coutCheminMin) {
      cheminMin = chemin;
      coutCheminMin = cout;
    } 
  }

  procedure explorerCircuits( Sequence chemin, int n )
  {
    if (longueur(chemin) == n) {

      if ( areteExiste(queue(chemin), tete(chemin)) ) {
        int cout = coutChemin(chemin) + W[queue(chemin)][tete(chemin)];
        if (cout < coutCheminMin) {
          # On vient de trouver un circuit ayant un cout inferieur.
          nbCircuitsComplets += 1;
          majCheminMin(chemin, cout);
        }
      } else {  
        # On ne fait rien (pas de cycle vers le sommet de depart).
      }

    } else {
      # Certains sommets n'ont pas encore ete visites.

      FilePriorite aExplorer = creerFP();
      # Pour la position suivante, on explore, parmi tous les sommets n'ayant
      # pas encore ete visites, uniquement ceux pour lesquels une arete existe.
      for [s = 1 to n st ~dejaVisite(chemin, s)] {
        int estime = coutEstime(chemin, s, n);
        if ( areteExiste(queue(chemin), s) & estime < coutCheminMin ) {
          # Une arete existe: on l'ajoute dans la file de priorite des noeuds a examiner.
          Sequence ch = cloner(chemin); ajouterEnQueue(ch, s);
          ajouter(aExplorer, estime, ch);
        } else {
          # Aucune arete => cul-de-sac.
        }
      }
      # On explore les sommets en ordre croissant de couts estimes.
      while (~estVide(aExplorer)) {
        int estime;
        Sequence chemin;
        obtenirMin(aExplorer, estime, chemin);
        if ( estime < coutCheminMin ) {
          # On explore uniquement si l'estime est encore interessant
          explorerCircuits(chemin, n);
        }
      }
    }
  }

  procedure trouverCircuitMin( int sommetDepart, int n )
  {
    Sequence chemin = creer();

    # On explore l'ensemble des circuits.
    ajouterEnQueue(chemin, sommetDepart);
    explorerCircuits(chemin, n);

    # On imprime le resultat.
    imprimer( cheminMin ); 
    write( " =>", coutCheminMin );
 
  }

  trouverCircuitMin( 1, n );
  write( nbCircuitsComplets, " circuits complets ont ete examines et ", nbNoeudsParcourus, "noeuds ont ete parcourus" ); 
end