Version dynamique de l'algo TSP

tsp_dyn.c

double tsp_prog_dyn(point *V, int n, int *Q) {
  /*
    Version programmation dynamique du TSP. Le résultat (la tournée
    optimale) doit être écrit dans la permutation Q, tableau qui doit
    être alloué avant l'appel. On renvoie la valeur de la tournée
    optimale ou -1 s'il y a eut une erreur (pression de 'q' pour sortir
    de l'affichage). Une fois que votre programme arrivera à calculer
    une tournée optimale (à comparer avec tsp_brute_force() ou
    tsp_brute_force_opt()), il sera intéressant de dessiner à chaque
    fois que possible le chemin courant avec drawPath(V,n,Q). La
    variable "running" indique si l'affichage est actif (ce qui peut
    ralentir le programme), la pression de 'q' faisant passer running à
    faux.

    Pour résumer, D est un tableau de cellules à deux dimensions indexé
    par un sommet t (int) et un ensemble S (un int aussi):

    - D[t][S].length = longueur minimum d'un chemin allant de V[n-1] à
      V[t] qui visite tous les points d'indice dans S

    - D[t][S].pred = l'indice du point précédant V[t] dans le chemin
      ci-dessus de longueur D[t][S].length

    NB: ne pas lancer tsp_prog_dyn() avec n>31 car sinon les entiers
    (int sur 32 bits) ne seront pas assez grands pour stocker les
    sous-ensembles. De plus, pour n=32, n*2^n / 10^9 = 137 secondes
    est un peu long pour un ordinateur à 1 GHz. On pourrait
    outrepasser cette limitation en utilisant des "long" sur 64 bits
    pour coder les ensembles. Mais, outre le temps de calcul, l'espace
    mémoire (le malloc() pour la table D) risque d'être problématique:
    32 * 2^32 * sizeof(cell) représente déjà 1,536 To de mémoire.
  */

  int t, S, T, x;

  // déclaration de la table D[t][S] qui comportent (n-1)*2^(n-1) cellules
  cell **D = malloc((n - 1) * sizeof(cell *)); // n-1 lignes
  for (t = 0; t < n - 1; t++)
    D[t] = malloc((1 << (n - 1)) * sizeof(cell)); // 2^{n-1} colonnes

  // cas de base, lorsque S={t}, pour t=0..n-2: D[t][S]=d(V[n-1],V[t])
  S = 1; // S=00...001 et pour l'incrémenter utiliser NextSet(S,n-1)
  for (t = 0; t < n-1; t++) {
    D[t][S].length = dist(V[n-1], V[t]);
    D[t][S].pred = n - 1;
    S = NextSet(S, n - 1);
  }

  // cas |S|>1 (cardinalité de S > 1): D[t][S] = min_x { D[x][S\{t}] +
  // d(V[x],V[t]) } avec t dans S et x dans S\{t}. NB: pour faire T =
  // S\{t}, utilisez T=DeleteSet(S,t); et pour savoir si x appartient
  // à l'ensemble S testez tout simplement si (S != DeleteSet(S,x)).
  
  // S = 2;
  do {
    for (t = 0; t < n - 1; t++) {
      T = DeleteSet(S, t);
      if (T != S) {
        double dmin = DBL_MAX;
        int vmin;
        for (x = 0; x < n - 1; x++) {
          if (T != DeleteSet(T, x)) {
            double d = D[x][T].length + dist(V[x], V[t]);
            if (d < dmin) {
              dmin = d;
              vmin = x;
            }
          }
        }
        D[t][S].length = dmin;
        D[t][S].pred = vmin;

      }
    }

    // ici D[t][S].length et D[t][S].pred viennent d'être calculés.
    // Il reste a extraire le chemin de V[t] à V[n-1] et le mettre
    // dans Q, c'est le rôle de la fonction ExtractPath().

    int k = ExtractPath(D, t, S, n, Q);
    drawPath(V, n, Q, k);
    // }
    S = NextSet(S, n - 1);
  } while (S && running);

  double w = 0; // valeur de la tournée optimale, 0 par défaut
  S = (1 << (n - 1)) - 1;
  if (running) {
    // on a terminé correctement le calcul complet de la table D. Il
    // faut alors calculer la longueur w de la tournée optimale et
    // éventuellement construire la tournée pour la mettre dans Q. NB:
    // si le calcul a été interrompu (pression de 'q' pendant
    // l'affichage) alors ne rien faire et renvoyer 0
    double min = DBL_MAX;
    for(int t=0;t<n-1;t++){
        double tmp = D[t][S].length+dist(V[t],V[n-1]);
        if(tmp<min) {
            min=tmp;
            int k = ExtractPath(D, t, S, n, Q);
            drawPath(V, n, Q, k);
        }
    }
    w= min;
  }

  for (t = 0; t < n - 1; t++)
    free(D[t]);
  free(D);
  return w;
}