1 files changed, 48 insertions, 58 deletions

diff --git a/biroli-mezard.cpp b/biroli-mezard.cpp
index 6f390e3..58c1e6a 100644
--- a/biroli-mezard.cpp
+++ b/biroli-mezard.cpp
@@ -151,6 +151,7 @@ public:
   unsigned occupiedNeighbors;
   unsigned maximumNeighbors;
   bool marked;
+  bool visited;
 
   Vertex() {
     occupiedNeighbors = 0;
@@ -277,6 +278,16 @@ public:
     return o;
   }
 
+  std::list<std::reference_wrapper<Vertex<D>>> overlaps(const std::list<std::reference_wrapper<Vertex<D>>>& vs) {
+    std::list<std::reference_wrapper<Vertex<D>>> o;
+
+    for (Vertex<D>& v : vs) {
+      o.splice(o.begin(), overlaps(v));
+    }
+
+    return o;
+  }
+
   bool canReplaceWith(const Vertex<D>& v, unsigned t) const {
     if (t == Empty) {
       return true;
@@ -446,76 +457,55 @@ public:
     }
   }
 
-  unsigned flipCluster(const Transformation<D>& R, Vertex<D>& v0) {
+  unsigned flipClusterHelper(const Transformation<D>& R, Vertex<D>& v, bool initial = false) {
     unsigned n = 0;
 
-    Vertex<D>& v0New = vertices[vectorToIndex(R.apply(v0.position))];
-
-    if (&v0New != &v0) {
-      std::list<std::reference_wrapper<Vertex<D>>> q;
-
-      v0.marked = true;
-      v0New.marked = true;
-
-      swap(v0, v0New);
-
-      for (Vertex<D>& vn : overlaps(v0)) {
-        if (!vn.marked) {
-          q.push_front(vn);
-        }
-      }
-      for (Vertex<D>& vn : overlaps(v0New)) {
-        if (!vn.marked) {
-          q.push_front(vn);
-        }
-      }
-
-      while (!q.empty()) {
-        Vertex<D>& v = q.front();
-        q.pop_front();
-
-        if (!v.marked && !overlaps(v).empty()) {
-          v.marked = true;
-          Vertex<D>& vNew = vertices[vectorToIndex(R.apply(v.position))];
+    if (!v.marked && (!overlaps(v).empty() || initial)) {
+      v.marked = true;
+      Vertex<D>& vNew = vertices[vectorToIndex(R.apply(v.position))];
 
-          if (&vNew != &v) {
-            vNew.marked = true;
+      if (&vNew != &v) {
+        vNew.marked = true;
 
-            swap(v, vNew);
-
-            for (Vertex<D>& vn : overlaps(v)) {
-              if (!vn.marked) {
-                q.push_front(vn);
-              }
-            }
-
-            for (Vertex<D>& vn : overlaps(vNew)) {
-              if (!vn.marked) {
-                q.push_front(vn);
-              }
-            }
-
-            n += 2;
-          } else {
-            n += 1;
+        swap(v, vNew);
+        /*
+        while (!overlaps({v, vNew}).empty()) {
+          for (Vertex<D>& vn : overlaps({v, vNew})) {
+            flipClusterHelper(R, vn);
           }
         }
-
-        // For some reason, the process above misses frustrated states. This is
-        // an inelegant stopgap.
-        if (q.empty()) {
-          for (Vertex<D>& vv : vertices) {
-            if (!vv.marked && !overlaps(vv).empty()) {
-              std::cerr << "Found bad vertex." << std::endl;
-              q.push_front(vv);
-            }
+        */
+        for (Vertex<D>& vn : overlaps({v, vNew})) {
+          if (!vn.marked) {
+            n += flipClusterHelper(R, vn);
           }
         }
+
+        n += 2;
+      } else {
+        n += 1;
       }
     } else {
-      n = 1;
     }
 
+    return n;
+  }
+
+  unsigned flipCluster(const Transformation<D>& R, Vertex<D>& v0) {
+    unsigned n = flipClusterHelper(R, v0, true);
+
+    bool badVertex = true;
+    while (badVertex) {
+      badVertex = false;
+      for (Vertex<D>& vv : vertices) {
+        if (!vv.marked && !overlaps(vv).empty()) {
+          std::cerr << "Found bad vertex." << std::endl;
+          n += flipClusterHelper(R, vv);
+          badVertex = true;
+        }
+      }
+    } 
+
     if (n > occupancy() / 2) {
       orientation = R.apply(orientation);
     }