Some simplifying refactoring.

1 files changed, 74 insertions, 62 deletions

diff --git a/rbmp.cpp b/rbmp.cpp
index 1d0c47e..3a689b9 100644
--- a/rbmp.cpp
+++ b/rbmp.cpp
@@ -14,28 +14,32 @@ using Rng = randutils::random_generator<pcg32>;
 class HalfEdge;
 
 class Vertex {
-private:
-  unsigned index;
-
 public:
-  std::vector<HalfEdge*> neighbors;
+  unsigned index;
+  std::vector<std::reference_wrapper<HalfEdge>> neighbors;
+  std::array<unsigned, 2> coordinate;
 
-  void setIndex(unsigned i) {index = i;}
-  unsigned getIndex() const {return index;}
   void addEdge(HalfEdge& e) {
-    neighbors.push_back(&e);
+    neighbors.push_back(e);
   }
 };
 
+class Edge;
+
 class HalfEdge {
 private:
   Vertex* tail;
   Vertex* head;
 
 public:
+  const Edge& edge;
   double oldX;
   double X;
-  double weight;
+
+  HalfEdge(const Edge& e) : edge(e) {
+    X = 0;
+    oldX = X;
+  }
   void setVertices(Vertex& vt, Vertex& vh) {
     tail = &vt;
     head = &vh;
@@ -54,78 +58,86 @@ public:
   }
 };
 
-class BipartiteGraph {
+class Edge {
+public:
+  std::array<HalfEdge, 2> halfedges;
+  double weight;
+
+  Edge() : halfedges{*this, *this} {};
+  void setVertices(Vertex& vt, Vertex& vh) {
+    halfedges[0].setVertices(vt, vh);
+    halfedges[1].setVertices(vh, vt);
+    vt.addEdge(halfedges[0]);
+    vh.addEdge(halfedges[1]);
+  }
+  double belief() const {
+    return halfedges[0].X + halfedges[1].X;
+  }
+  bool active() const {
+    return belief() >= 0;
+  }
+};
+
+class Graph {
 public:
-  std::vector<Vertex> redVertices;
-  std::vector<Vertex> blueVertices;
-  std::vector<HalfEdge> redEdges;
-  std::vector<HalfEdge> blueEdges;
-
-  BipartiteGraph(unsigned n, Rng& r) : redVertices(n * (n + 1)), blueVertices(n * (n + 1)), redEdges(pow(2 * n, 2)), blueEdges(pow(2 * n, 2)) {
-    for (unsigned i = 0; i < redVertices.size(); i++) {
-      redVertices[i].setIndex(i);
-      blueVertices[i].setIndex(i);
+  std::vector<Vertex> vertices;
+  std::vector<Edge> edges;
+
+  Graph(unsigned n, Rng& r) : vertices(2 * n * (n + 1)), edges(pow(2 * n, 2)) {
+    for (unsigned i = 0; i < vertices.size() / 2; i++) {
+      vertices[i].index = i;
+      vertices[i].coordinate = {2 * (i % (n + 1)), 2 * (i / (n + 1)) + 1};
+      vertices[vertices.size() / 2 + i].index = i + vertices.size() / 2;
+      vertices[vertices.size() / 2 + i].coordinate = {2 * (i % n) + 1, 2 * (i / n)};
     }
-    for (unsigned i = 0; i < redEdges.size(); i++) {
-      Vertex& blueVertex = blueVertices[(i % (2 * n)) / 2 + n * (((i + 2*n) / 4) / n)];
-      Vertex& redVertex = redVertices[(1+(i % (2 * n))) / 2 + (n + 1) * ((i / 4) / (n))];
-      redEdges[i].setVertices(redVertex, blueVertex);
-      blueEdges[i].setVertices(blueVertex, redVertex);
-      redVertex.addEdge(redEdges[i]);
-      blueVertex.addEdge(blueEdges[i]);
-      redEdges[i].X = 0;
-      blueEdges[i].X = 0;
-      redEdges[i].oldX = 0;
-      blueEdges[i].oldX = 0;
-      double w = r.variate<double, std::exponential_distribution>(1);
-      redEdges[i].weight = w;
-      blueEdges[i].weight = w;
+    for (unsigned i = 0; i < edges.size(); i++) {
+      Vertex& blueVertex = vertices[vertices.size() / 2 + (i % (2 * n)) / 2 + n * (((i + 2*n) / 4) / n)];
+      Vertex& redVertex = vertices[(1+(i % (2 * n))) / 2 + (n + 1) * ((i / 4) / (n))];
+      edges[i].setVertices(redVertex, blueVertex);
+      edges[i].weight = r.variate<double, std::exponential_distribution>(1);
     }
   }
 
-  void propagateBeliefs() {
-    for (HalfEdge& e : redEdges) {
-      double Xt = std::numeric_limits<double>::infinity();
-      for (const HalfEdge* en : e.getHead().neighbors) {
-        if (e.getTail().getIndex() != en->getHead().getIndex()) {
-          Xt = std::min(en->weight - en->oldX, Xt);
+  double propagateBeliefs() {
+    for (Edge& e : edges) {
+      for (HalfEdge& he : e.halfedges) {
+        double Xt = std::numeric_limits<double>::infinity();
+        for (const HalfEdge& en : he.getHead().neighbors) {
+          if (he.getTail().index != en.getHead().index) {
+            Xt = std::min(en.edge.weight - en.oldX, Xt);
+          }
         }
+        he.X = Xt;
       }
-      e.X = Xt;
     }
-    for (HalfEdge& e : blueEdges) {
-      double Xt = std::numeric_limits<double>::infinity();
-      for (const HalfEdge* en : e.getHead().neighbors) {
-        if (e.getTail().getIndex() != en->getHead().getIndex()) {
-          Xt = std::min(en->weight - en->oldX, Xt);
-        }
+
+    double ΔX = std::numeric_limits<double>::infinity();
+
+    for (Edge& e : edges) {
+      ΔX = std::min(ΔX, e.belief());
+
+      for (HalfEdge& he : e.halfedges) {
+        he.oldX = he.X;
       }
-      e.X = Xt;
-    }
-    for (HalfEdge& e : redEdges) {
-      e.oldX = e.X;
-    }
-    for (HalfEdge& e : blueEdges) {
-      e.oldX = e.X;
     }
+
+    return ΔX;
   }
 };
 
 int main() {
   Rng r;
-  BipartiteGraph G(10, r);
+  Graph G(1000, r);
 
-  for (unsigned i = 0; i < 1000; i++) {
-    G.propagateBeliefs();
-//    std::cout << G.redEdges[40].X << " " << G.blueEdges[40].X << std::endl;
-  }
+  double Δx = 0;
 
-  for (unsigned i = 0; i < G.blueEdges.size(); i++) {
-    std::cout << G.blueEdges[i].getTail().getIndex() << " " << G.blueEdges[i].getHead().getIndex() + G.blueVertices.size() << std::endl;
+  while (Δx < 1) {
+    Δx = G.propagateBeliefs();
   }
-  for (unsigned i = 0; i < G.blueEdges.size(); i++) {
-    if (G.blueEdges[i].X + G.redEdges[i].X >= 0) {
-      std::cout << G.blueEdges[i].getTail().getIndex() << " " << G.blueEdges[i].getHead().getIndex() + G.blueVertices.size() << std::endl;
+
+  for (unsigned i = 0; i < G.edges.size(); i++) {
+    if (G.edges[i].active()) {
+      std::cout << G.edges[i].halfedges[0].getTail().coordinate[0] << " " << G.edges[i].halfedges[0].getTail().coordinate[1] << " " << G.edges[i].halfedges[0].getHead().coordinate[0] << " " << G.edges[i].halfedges[0].getHead().coordinate[1] << std::endl;
     }
   }