Update.
authorFrancois Fleuret <francois@fleuret.org>
Tue, 21 Aug 2012 03:49:15 +0000 (20:49 -0700)
committerFrancois Fleuret <francois@fleuret.org>
Tue, 21 Aug 2012 03:49:15 +0000 (20:49 -0700)
miniksp.cc

index f562423..4ec9423 100644 (file)
@@ -20,7 +20,7 @@
 // END_IP_HEADER                                                         //
 ///////////////////////////////////////////////////////////////////////////
 
-#define VERBOSE
+// #define VERBOSE
 
 #include <iostream>
 #include <fstream>
@@ -42,6 +42,8 @@ typedef float scalar_t;
 #define ASSERT(x)
 #endif
 
+// Adds to all the edge length the length of the shortest (which can
+// be negative)
 void raise_es(int nb_edges, scalar_t *el) {
   scalar_t min_es = el[0];
   for(int e = 1; e < nb_edges; e++) {
@@ -52,30 +54,33 @@ void raise_es(int nb_edges, scalar_t *el) {
   }
 }
 
+// Add to every edge length the differential of the psi function on it
 void add_dpsi_es(int nb_edges, scalar_t *el, int *ea, int *eb, scalar_t *psi) {
   for(int e = 0; e < nb_edges; e++) {
     el[e] += psi[ea[e]] - psi[eb[e]];
   }
 }
 
+// Find the shortest path in the graph and return in return_edge_back,
+// for each vertex the edge to follow back from it to reach the source
+// with the shortest path, and in return_dist the distance to the
+// source. The edge lengths have to be positive.
 void find_shortest(int nb_vertices,
                    int nb_edges, scalar_t *el, int *ea, int *eb,
                    int source, int sink,
-                   int *edge_back, scalar_t *dist) {
+                   int *return_edge_back, scalar_t *return_dist) {
   for(int v = 0; v < nb_vertices; v++) {
-    // dist[v] = FLT_MAX;
-    dist[v] = 999;
+    dist[v] = FLT_MAX;
+    edge_back[v] = -1;
   }
 
   dist[source] = 0;
 
-  for(int v = 0; v < nb_vertices; v++) {
-    edge_back[v] = -1;
-  }
-
+#ifdef DEBUG
   for(int e = 0; e < nb_edges; e++) {
-    ASSERT(el[e] >= 0);
+    if(el[e] < 0) abort();
   }
+#endif
 
   int nb_changes;
   scalar_t d;
@@ -98,43 +103,62 @@ void find_best_paths(int nb_vertices,
                      int *edge_occupation) {
   scalar_t *dist = new scalar_t[nb_vertices];
   int *edge_back = new int[nb_vertices];
-  scalar_t *tmp_el = new scalar_t[nb_edges];
+  scalar_t *positive_el = new scalar_t[nb_edges];
+  scalar_t s;
+  int v;
 
   for(int e = 0; e < nb_edges; e++) {
-    tmp_el[e] = el[e];
+    positive_el[e] = el[e];
     edge_occupation[e] = 0;
   }
 
-  raise_es(nb_edges, tmp_el);
+  raise_es(nb_edges, positive_el);
 
-  scalar_t s;
   do {
-    find_shortest(nb_vertices, nb_edges, tmp_el, ea, eb, source, sink, edge_back, dist);
-    add_dpsi_es(nb_edges, tmp_el, ea, eb, dist);
+    find_shortest(nb_vertices, nb_edges, positive_el, ea, eb, source, sink, edge_back, dist);
+    add_dpsi_es(nb_edges, positive_el, ea, eb, dist);
     s = 0.0;
+
+    // If the new path reaches the sink, we will backtrack on it to
+    // compute its score and invert edges
+
     if(edge_back[sink] >= 0) {
-      int v = sink;
+
+      v = sink;
       while(v != source) {
         int e = edge_back[v];
         ASSERT(eb[e] == v);
         v = ea[e];
-#ifdef VERBOSE
-        cout << "INVERTING " << ea[e] << " -> " << eb[e] << " [" << el[e] << "]" << endl;
-#endif
         s += el[e];
-        int k = eb[e];
-        eb[e] = ea[e];
-        ea[e] = k;
-        tmp_el[e] = - tmp_el[e];
-        edge_occupation[e] = 1 - edge_occupation[e];
       }
+
+      // We found a good path (score < 0), we revert the edges along
+      // the path and note that they are occupied
+
+      if(s < 0) {
+        v = sink;
+#ifdef VERBOSE
+        cout << "FOUND A PATH OF LENGTH " << s << endl;
+#endif
+        while(v != source) {
+          int e = edge_back[v];
+          ASSERT(eb[e] == v);
+          v = ea[e];
 #ifdef VERBOSE
-      cout << "FOUND A PATH OF LENGTH " << s << endl;
+          cout << "INVERTING " << ea[e] << " -> " << eb[e] << " [" << el[e] << "]" << endl;
 #endif
+          int k = eb[e];
+          eb[e] = ea[e];
+          ea[e] = k;
+          positive_el[e] = - positive_el[e];
+          el[e] = - el[e];
+          edge_occupation[e] = 1 - edge_occupation[e];
+        }
+      }
     }
   } while(s < 0);
 
-  delete[] tmp_el;
+  delete[] positive_el;
   delete[] dist;
   delete[] edge_back;
 }
@@ -146,19 +170,20 @@ int main(int argc, char **argv) {
     exit(EXIT_FAILURE);
   }
 
-  ifstream file(argv[1]);
+  ifstream *file = new ifstream(argv[1]);
 
   int nb_edges, nb_vertices;
   int source, sink;
 
-  if(file.good()) {
-    file >> nb_vertices >> nb_edges;
-    file >> source >> sink;
+  if(file->good()) {
+
+    (*file) >> nb_vertices >> nb_edges;
+    (*file) >> source >> sink;
 
-    cout << "nb_edges = " << nb_edges << endl;
-    cout << "nb_vertices = " << nb_vertices << endl;
-    cout << "source = " << source << endl;
-    cout << "sink = " << sink << endl;
+    cout << "INPUT nb_edges " << nb_edges << endl;
+    cout << "INPUT nb_vertices " << nb_vertices << endl;
+    cout << "INPUT source " << source << endl;
+    cout << "INPUT sink " << sink << endl;
 
     scalar_t *el = new scalar_t[nb_edges];
     int *ea = new int[nb_edges];
@@ -166,10 +191,8 @@ int main(int argc, char **argv) {
     int *edge_occupation = new int[nb_edges];
 
     for(int e = 0; e < nb_edges; e++) {
-      file >> ea[e] >> eb[e] >> el[e];
-#ifdef VERBOSE
-      cout << ea[e] << " -> " << eb[e] << " [" << el[e] << "]" << endl;
-#endif
+      (*file) >> ea[e] >> eb[e] >> el[e];
+      cout << "INPUT_EDGE " << ea[e] << " " << eb[e] << " " << el[e] << endl;
     }
 
     find_best_paths(nb_vertices, nb_edges, el, ea, eb, source, sink,
@@ -177,7 +200,7 @@ int main(int argc, char **argv) {
 
     for(int e = 0; e < nb_edges; e++) {
       if(edge_occupation[e]) {
-        cout << ea[e] << " " << eb[e] << endl;
+        cout << "RESULT_OCCUPIED_EDGE " << ea[e] << " " << eb[e] << endl;
       }
     }
 
@@ -189,9 +212,12 @@ int main(int argc, char **argv) {
   } else {
 
     cerr << "Can not open " << argv[1] << endl;
+
+    delete file;
     exit(EXIT_FAILURE);
 
   }
 
+  delete file;
   exit(EXIT_SUCCESS);
 }