From 122c2b4aa525df69a2e95d47e9b1202e02fc4301 Mon Sep 17 00:00:00 2001
From: Jaron Kent-Dobias <jaron@kent-dobias.com>
Date: Mon, 10 Oct 2022 21:39:20 +0200
Subject: Some refactoring, and example code for comparing low-T probabilities
 to zero-T.

---
 uniform.cpp | 45 +++++++++++++++++++++++++++++++++++++--------
 1 file changed, 37 insertions(+), 8 deletions(-)

(limited to 'uniform.cpp')

diff --git a/uniform.cpp b/uniform.cpp
index c6f7869..7e2d583 100644
--- a/uniform.cpp
+++ b/uniform.cpp
@@ -1,6 +1,4 @@
-#include <functional>
 #include <fstream>
-#include <stack>
 
 #include "rbmp.hpp"
 
@@ -30,16 +28,13 @@ int main(int argc, char* argv[]) {
   std::string filename = "order_" + std::to_string(n) + "_" + std::to_string(T) + ".dat";
 
   Rng r;
-  AztecDiamond a(n, r);
+  AztecDiamond a(n);
 
   std::vector<Real> avgProbabilities(a.edges.size());
 
   for (unsigned i = 0; i < m; i++) {
-    for (AztecDiamond::Edge& e : a.edges) {
-      e.weights.push(exp(- r.variate<Real, std::exponential_distribution>(1) / T));
-      e.probability = 0;
-    }
-    a.computeWeights();
+    a.setWeights(r);
+    a.computeWeights(T);
     a.computeProbabilities();
 
     for (unsigned j = 0; j < a.edges.size(); j++) {
@@ -68,5 +63,39 @@ int main(int argc, char* argv[]) {
 
   output.close();
 
+  std::vector<long int> data_x0(a.vertices.size());
+  std::vector<long int> data_y0(a.vertices.size());
+
+  PerfectMatching pm(a.vertices.size(), a.edges.size());
+
+  for (const AztecDiamond::Edge& e : a.edges) {
+    pm.AddEdge(e.head->index, e.tail->index, e.weight);
+  }
+
+  pm.options.verbose = false;
+  pm.Solve();
+
+  for (unsigned i = 0; i < a.vertices.size() / 2; i++) {
+    unsigned j = pm.GetMatch(i);
+
+    data_x0[i] += a.vertices[i].coordinate[0];
+    data_y0[i] += a.vertices[i].coordinate[1];
+    data_x0[i] -= a.vertices[j].coordinate[0];
+    data_y0[i] -= a.vertices[j].coordinate[1];
+
+    data_x0[j] += a.vertices[i].coordinate[0];
+    data_y0[j] += a.vertices[i].coordinate[1];
+    data_x0[j] -= a.vertices[j].coordinate[0];
+    data_y0[j] -= a.vertices[j].coordinate[1];
+  }
+
+  std::ofstream output2("order_0.dat");
+
+  for (unsigned i = 0; i < a.vertices.size(); i++) {
+    output2 << data_x0[i] << " " << data_y0[i] << std::endl;
+  }
+
+  output.close();
+
   return 0;
 }
-- 
cgit v1.2.3-70-g09d2