first working version

diff --git a/lib/neatogen/CMakeLists.txt b/lib/neatogen/CMakeLists.txt
index 6ccb9455a7c8981fbe4eb40f9d85db01b3ab7025..eeac6ac82f3f3f738bb8a198aadd932cc47d8907 100644 (file)
--- a/lib/neatogen/CMakeLists.txt
+++ b/lib/neatogen/CMakeLists.txt
@@ -45,6 +45,7 @@ set(sources
     sparsegraph.h
     stress.h
     voronoi.h
+       sgd.h
 
     # Source files
     adjust.c
@@ -85,6 +86,7 @@ set(sources
     stuff.c
     stress.c
     voronoi.c
+       sgd.c
 )
 
 if (with_ipsepcola)

diff --git a/lib/neatogen/Makefile.am b/lib/neatogen/Makefile.am
index 039e7f712d538962f69ad3739337d6d637e80165..feea48a6ab357d6bdbc42d37586a0e087b4e654c 100644 (file)
--- a/lib/neatogen/Makefile.am
+++ b/lib/neatogen/Makefile.am
@@ -21,7 +21,8 @@ noinst_HEADERS = adjust.h edges.h geometry.h heap.h hedges.h info.h mem.h \
        bfs.h closest.h conjgrad.h defs.h dijkstra.h embed_graph.h kkutils.h \
        matrix_ops.h pca.h stress.h quad_prog_solver.h digcola.h \
     overlap.h call_tri.h \
-       quad_prog_vpsc.h delaunay.h sparsegraph.h multispline.h fPQ.h
+       quad_prog_vpsc.h delaunay.h sparsegraph.h multispline.h fPQ.h \
+       sgd.h
 
 IPSEPCOLA_SOURCES = constrained_majorization_ipsep.c \
        mosek_quad_solve.c mosek_quad_solve.h quad_prog_vpsc.c 
@@ -37,6 +38,7 @@ libneatogen_C_la_SOURCES = adjust.c circuit.c edges.c geometry.c \
        conjgrad.c pca.c closest.c bfs.c constraint.c quad_prog_solve.c \
        smart_ini_x.c constrained_majorization.c opt_arrangement.c \
     overlap.c call_tri.c \
-       compute_hierarchy.c delaunay.c multispline.c $(WITH_IPSEPCOLA_SOURCES)
+       compute_hierarchy.c delaunay.c multispline.c $(WITH_IPSEPCOLA_SOURCES) \
+       sgd.c
 
 EXTRA_DIST = $(IPSEPCOLA_SOURCES) gvneatogen.vcxproj*

diff --git a/lib/neatogen/neatoinit.c b/lib/neatogen/neatoinit.c
index 4da2a477f8a79b22d2132dc2d44642676eb42f63..c99c3e4c5bb2ae40d6d8b4ab7d3f4b032f830bd4 100644 (file)
--- a/lib/neatogen/neatoinit.c
+++ b/lib/neatogen/neatoinit.c
@@ -28,6 +28,7 @@
 #endif
 #include "kkutils.h"
 #include "pointset.h"
+#include "sgd.h"
 
 #ifndef HAVE_SRAND48
 #define srand48 srand
@@ -1335,41 +1336,6 @@ static void kkNeato(Agraph_t * g, int nG, int model)
     }
     solve_model(g, nG);
 }
-/* sgdNeato:
- * Solve using stochastic gradient descent a la Zheng, Pawar, Goodman.
- */
-static void sgdNeato(Agraph_t * g, int nG, int model)
-{
-    fprintf("hello, love from sgd");
-
-    if (model == MODEL_SUBSET) {
-        subset_model(g, nG);
-    } else if (model == MODEL_CIRCUIT) {
-        if (!circuit_model(g, nG)) {
-            agerr(AGWARN,
-                 "graph %s is disconnected. Hence, the circuit model\n",
-                 agnameof(g));
-            agerr(AGPREV,
-                 "is undefined. Reverting to the shortest path model.\n");
-            agerr(AGPREV,
-                 "Alternatively, consider running neato using -Gpack=true or decomposing\n");
-            agerr(AGPREV, "the graph into connected components.\n");
-            shortest_path(g, nG);
-        }
-    } else if (model == MODEL_MDS) {
-        shortest_path(g, nG);
-        mds_model(g, nG);
-    } else
-        shortest_path(g, nG);
-    initial_positions(g, nG);
-    diffeq_model(g, nG);
-    if (Verbose) {
-        fprintf(stderr, "Solving model %d iterations %d tol %f\n",
-               model, MaxIter, Epsilon);
-        start_timer();
-    }
-    solve_model(g, nG);
-}
 
 /* neatoLayout:
  * Use stress optimization to layout a single component
@@ -1394,7 +1360,7 @@ neatoLayout(Agraph_t * mg, Agraph_t * g, int layoutMode, int layoutModel,
        if (layoutMode == MODE_KK)
        kkNeato(g, nG, layoutModel);
        else if (layoutMode == MODE_SGD)
-       sgdNeato(g, nG, layoutModel)
+       sgd(g, nG, layoutModel);
        else
        majorization(mg, g, nG, layoutMode, layoutModel, Ndim, MaxIter, am);
 }

diff --git a/lib/neatogen/sgd.c b/lib/neatogen/sgd.c
new file mode 100644 (file)
index 0000000..a8d91ca
--- /dev/null
+++ b/lib/neatogen/sgd.c
@@ -0,0 +1,86 @@
+#include "sgd.h"
+#include "neatoprocs.h"
+#include <math.h>
+
+typedef struct term {
+    node_t *i, *j;
+    float d, w;
+} term;
+
+void fisheryates_shuffle(term *terms, int nC2) {
+    int i;
+    for (i=nC2-1; i>=1; i--) {
+        // unsigned j = rk_interval(i, &rstate);
+        int j = (int)(drand48() * (i+1)); // TODO: better rng
+        fprintf(stderr, "%d %d %f %f\n", i, j, terms[i].d, terms[i].w);
+
+        term temp = terms[i];
+        terms[i] = terms[j];
+        terms[j] = temp;
+    }
+}
+
+void sgd(graph_t *G, /* input graph */
+         int n, /* number of nodes */
+         int model /* distance model */)
+{
+    // calculate shortest paths (from neatoprocs)
+    shortest_path(G, n);
+
+    // initialise array of terms
+    int nC2 = (n*(n-1))/2;
+    term *terms = N_NEW(nC2, term);
+    int i, j;
+    int ij = 0;
+    for (i=0; i<n; i++) {
+        for (j=0; j<i; j++) {
+            terms[ij].i = GD_neato_nlist(G)[i];
+            terms[ij].j = GD_neato_nlist(G)[j];
+            terms[ij].d = GD_dist(G)[ND_id(terms[ij].i)][ND_id(terms[ij].j)];
+            terms[ij].w = 1 / (terms[ij].d*terms[ij].d);
+            ij++;
+        }
+    }
+
+    // initialise annealing schedule
+    float w_min = terms[0].w, w_max = terms[0].w;
+    for (ij=1; ij<nC2; ij++) {
+        if (terms[ij].w < w_min)
+            w_min = terms[ij].w;
+        if (terms[ij].w > w_max)
+            w_max = terms[ij].w;
+    }
+    float eta_max = 1 / w_min;
+    float eta_min = .01 / w_max; // TODO: read in epsilon
+    float lambda = log(eta_max/eta_min) / (30-1); // TODO: read in t_max
+
+    // initialise starting positions (from neatoprocs)
+    initial_positions(G, n);
+
+    // perform optimisation
+    int t;
+    for (t=0; t<30; t++) {
+        fisheryates_shuffle(terms, nC2);
+        float eta = eta_max * exp(-lambda * t);
+        for (ij=0; ij<nC2; ij++) {
+            // cap step size
+            float mu = eta * terms[ij].w;
+            if (mu > 1)
+                mu = 1;
+            
+            float dx = ND_pos(terms[ij].i)[0] - ND_pos(terms[ij].j)[0];
+            float dy = ND_pos(terms[ij].i)[1] - ND_pos(terms[ij].j)[1];
+            float mag = sqrt(dx*dx + dy*dy);
+            
+            float r = (mu * (mag-terms[ij].d)) / (2*mag);
+            float r_x = r * dx;
+            float r_y = r * dy;
+
+            ND_pos(terms[ij].i)[0] -= r_x;
+            ND_pos(terms[ij].i)[1] -= r_y;
+            ND_pos(terms[ij].j)[0] += r_x;
+            ND_pos(terms[ij].j)[1] += r_y;
+        }
+    }
+    free(terms);
+}

diff --git a/lib/neatogen/sgd.h b/lib/neatogen/sgd.h
new file mode 100644 (file)
index 0000000..b99263e
--- /dev/null
+++ b/lib/neatogen/sgd.h
@@ -0,0 +1,8 @@
+#ifndef SGD_H
+#define SGD_H
+
+#include "neato.h"
+
+extern void sgd(graph_t *, int, int);
+
+#endif /* SGD_H */

diff --git a/lib/neatogen/stuff.c b/lib/neatogen/stuff.c
index 9f60854e9e170aadcf0766d9e92576d392e52043..8d031650da7be9d5970f44ab563ae462f68d0ed1 100644 (file)
--- a/lib/neatogen/stuff.c
+++ b/lib/neatogen/stuff.c
@@ -247,7 +247,7 @@ int scan_graph_mode(graph_t * G, int mode)
     nE = agnedges(G);
 
     lenx = agattr(G, AGEDGE, "len", 0);
-    if (mode == MODE_KK) {
+    if (mode == MODE_KK || mode == MODE_SGD) {
        Epsilon = .0001 * nV;
        getdouble(G, "epsilon", &Epsilon);
        if ((str = agget(G->root, "Damping")))
@@ -276,7 +276,7 @@ int scan_graph_mode(graph_t * G, int mode)
     else
        Initial_dist = total_len / (nE > 0 ? nE : 1) * sqrt(nV) + 1;
 
-    if (!Nop && (mode == MODE_KK)) {
+    if (!Nop && (mode == MODE_KK || mode == MODE_SGD)) {
        GD_dist(G) = new_array(nV, nV, Initial_dist);
        GD_spring(G) = new_array(nV, nV, 1.0);
        GD_sum_t(G) = new_array(nV, Ndim, 1.0);