Randomize sweep order in sfdp_layout()

This fixes some strange bugs where deterministic position updates lead
to strange artifacts.

src/graph/layout/graph_sfdp.cc

@@ -33,7 +33,7 @@ void sfdp_layout(GraphInterface& g, boost::any pos, boost::any vweight,
                  boost::any eweight, boost::any pin, python::object spring_parms,
                  double theta, double init_step, double step_schedule,
                  size_t max_level, double epsilon, size_t max_iter,
-                 bool adaptive, bool verbose)
+                 bool adaptive, bool verbose, rng_t& rng)
 {
     typedef ConstantPropertyMap<int32_t,GraphInterface::vertex_t> vweight_map_t;
     typedef ConstantPropertyMap<int32_t,GraphInterface::edge_t> eweight_map_t;
@@ -68,12 +68,13 @@ void sfdp_layout(GraphInterface& g, boost::any pos, boost::any vweight,
     run_action<graph_tool::detail::never_directed>()
         (g,
          std::bind(get_sfdp_layout(C, K, p, theta, gamma, mu, mu_p, init_step,
-                              step_schedule, max_level, epsilon,
+                                   step_schedule, max_level, epsilon,
                                    max_iter, adaptive),
                    placeholders::_1, g.GetVertexIndex(), placeholders::_2,
                    placeholders::_3, placeholders::_4,
-              pin_map.get_unchecked(num_vertices(g.GetGraph())),
-              groups.get_unchecked(num_vertices(g.GetGraph())), verbose),
+                   pin_map.get_unchecked(num_vertices(g.GetGraph())),
+                   groups.get_unchecked(num_vertices(g.GetGraph())), verbose,
+                   std::ref(rng)),
          vertex_floating_vector_properties(), vertex_props_t(), edge_props_t())
         (pos, vweight, eweight);
 }

src/graph/layout/graph_sfdp.hh

@@ -200,10 +200,10 @@ struct get_sfdp_layout
     bool simple;
 
     template <class Graph, class VertexIndex, class PosMap, class VertexWeightMap,
-              class EdgeWeightMap, class PinMap, class GroupMap>
+              class EdgeWeightMap, class PinMap, class GroupMap, class RNG>
     void operator()(Graph& g, VertexIndex vertex_index, PosMap pos,
                     VertexWeightMap vweight, EdgeWeightMap eweight, PinMap pin,
-                    GroupMap group, bool verbose) const
+                    GroupMap group, bool verbose, RNG& rng) const
     {
         typedef typename property_traits<PosMap>::value_type pos_t;
         typedef typename property_traits<PosMap>::value_type::value_type val_t;
@@ -215,6 +215,7 @@ struct get_sfdp_layout
 
         vector<pos_t> group_cm;
         vector<vweight_t> group_size;
+        vector<size_t> vertices;
 
         int i, N = num_vertices(g), HN=0;
         for (i = 0; i < N; ++i)
@@ -223,6 +224,8 @@ struct get_sfdp_layout
                 vertex(i, g);
             if (v == graph_traits<Graph>::null_vertex())
                 continue;
+            if (!pin[v])
+                vertices.push_back(v);
             pos[v].resize(2, 0);
             size_t s = group[v];
 
@@ -280,21 +283,19 @@ struct get_sfdp_layout
                 qt.put_pos(pos[v], vweight[v]);
             }
 
+            std::shuffle(vertices.begin(), vertices.end(), rng);
+
             pos_t diff(2, 0), pos_u(2, 0), ftot(2, 0), cm(2, 0);
 
             size_t nmoves = 0;
+            N = vertices.size();
             #pragma omp parallel for default(shared) private(i)  \
                 firstprivate(Q, diff, pos_u, ftot, cm) \
                 reduction(+:E, delta, nmoves) schedule(static) if (N > 100)
             for (i = 0; i < N; ++i)
             {
                 typename graph_traits<Graph>::vertex_descriptor v =
-                    vertex(i, g);
-                if (v == graph_traits<Graph>::null_vertex())
-                    continue;
-
-                if (pin[v])
-                    continue;
+                    vertex(vertices[i], g);
 
                 ftot[0] = ftot[1] = 0;
 
@@ -357,8 +358,8 @@ struct get_sfdp_layout
                         target(*e, g);
                     if (u == v)
                         continue;
+                    #pragma omp critical
                     {
-                        #pragma omp critical
                         for (size_t l = 0; l < 2; ++l)
                             pos_u[l] = pos[u][l];
                     }
@@ -424,8 +425,8 @@ struct get_sfdp_layout
 
                 E += power(norm(ftot), 2);
 
+                #pragma omp critical
                 {
-                    #pragma omp critical
                     for (size_t l = 0; l < 2; ++l)
                     {
                         group_cm[group[v]][l] *= group_size[group[v]];

src/graph_tool/draw/__init__.py

@@ -507,7 +507,7 @@ def coarse_graph_stack(g, c, coarse_stack, eweight=None, vweight=None,
 
 def sfdp_layout(g, vweight=None, eweight=None, pin=None, groups=None, C=0.2,
                 K=None, p=2., theta=0.6, max_level=11, gamma=1., mu=0., mu_p=1.,
-                init_step=None, cooling_step=0.9, adaptive_cooling=True,
+                init_step=None, cooling_step=0.95, adaptive_cooling=True,
                 epsilon=1e-2, max_iter=0, pos=None, multilevel=None,
                 coarse_method="hybrid", mivs_thres=0.9, ec_thres=0.75,
                 coarse_stack=None, weighted_coarse=False, verbose=False):
@@ -549,7 +549,7 @@ def sfdp_layout(g, vweight=None, eweight=None, pin=None, groups=None, C=0.2,
         connected component, or group assignment.
     init_step : float (optional, default: ``None``)
         Initial update step. If not provided, it will be chosen automatically.
-    cooling_step : float (optional, default: ``0.9``)
+    cooling_step : float (optional, default: ``0.95``)
         Cooling update step.
     adaptive_cooling : bool (optional, default: ``True``)
         Use an adaptive cooling scheme.
@@ -634,7 +634,7 @@ def sfdp_layout(g, vweight=None, eweight=None, pin=None, groups=None, C=0.2,
         K = _avg_edge_distance(g, pos)
 
     if init_step is None:
-        init_step = 10 * max(_avg_edge_distance(g, pos), K)
+        init_step = 2 * max(_avg_edge_distance(g, pos), K)
 
     if multilevel is None:
         multilevel = g.num_vertices() > 1000
@@ -671,8 +671,8 @@ def sfdp_layout(g, vweight=None, eweight=None, pin=None, groups=None, C=0.2,
                               max_iter=max_iter,
                               cooling_step=cooling_step,
                               adaptive_cooling=False,
-                              init_step=max(2 * K,
-                                            _avg_edge_distance(u, pos) / 10),
+                              # init_step=max(2 * K,
+                              #               _avg_edge_distance(u, pos)),
                               multilevel=False,
                               verbose=False)
             #graph_draw(u, pos)
@@ -699,7 +699,7 @@ def sfdp_layout(g, vweight=None, eweight=None, pin=None, groups=None, C=0.2,
                                      (C, K, p, gamma, mu, mu_p, _prop("v", g, groups)),
                                      theta, init_step, cooling_step, max_level,
                                      epsilon, max_iter, not adaptive_cooling,
-                                     verbose)
+                                     verbose, _get_rng())
     return pos
 
 def radial_tree_layout(g, root, weighted=False, r=1.):