cleanup, and change to initalization with walsh matrices

1 files changed, 72 insertions, 66 deletions

diff --git a/hadamard.cpp b/hadamard.cpp
index 6b9a64a..15f73b1 100644
--- a/hadamard.cpp
+++ b/hadamard.cpp
@@ -1,57 +1,54 @@
 
 #include "hadamard_pt.hpp"
 
-#include <iostream>
 #include <fstream>
+#include <iostream>
 
 class MeasureEnergy : public Measurement {
-  public:
-    unsigned N;
-    double totalE;
-    MeasureEnergy() {
-      N = 0;
-      totalE = 0;
-    }
+public:
+  unsigned N;
+  double totalE;
 
-    void after_sweep(double, double E, const Orthogonal&) override {
-      totalE += E;
-      N++;
-    }
+  MeasureEnergy() {
+    N = 0;
+    totalE = 0;
+  }
 
-    double energy() const {
-      return totalE / N;
-    }
+  void after_sweep(double, double E, const Orthogonal&) override {
+    totalE += E;
+    N++;
+  }
+
+  double energy() const { return totalE / N; }
 };
 
 class MeasureTransitionRates : public ParallelMeasurement {
-  public:
-    std::vector<std::vector<unsigned>> nAccepted;
-    unsigned total_steps;
-
-    MeasureTransitionRates(unsigned n) : nAccepted(n - 1) {
-      total_steps = 0;
-      for (unsigned i = 0; i < n - 1; i++) {
-        nAccepted[i].resize(i + 1);
-      }
+public:
+  std::vector<std::vector<unsigned>> nAccepted;
+  unsigned total_steps;
+
+  MeasureTransitionRates(unsigned n) : nAccepted(n - 1) {
+    total_steps = 0;
+    for (unsigned i = 0; i < n - 1; i++) {
+      nAccepted[i].resize(i + 1);
     }
+  }
 
-    void after_step(bool accepted, unsigned i, unsigned j, double, double, const MCMC&, const MCMC&) override {
-      if (accepted) {
-        nAccepted[j - 1][i]++;
-      }
-    }
+  void after_step(bool accepted, unsigned i, unsigned j, double, double, const MCMC&,
+                  const MCMC&) override {
+    if (accepted)
+      nAccepted[j - 1][i]++;
+  }
 
-    void after_sweep(const std::vector<MCMC>&) override {
-      total_steps++;
-    }
+  void after_sweep(const std::vector<MCMC>&) override { total_steps++; }
 };
 
 int main(int argc, char* argv[]) {
-  unsigned trials = 1e6;
+  unsigned n_tuning = 1e2;
   double β0 = 1.00;
   double βf = 40.00;
   unsigned num = 40;
-  unsigned n = 20;
+  unsigned k = 2;
   double ε = 0.01;
 
   unsigned M = 10;
@@ -59,37 +56,39 @@ int main(int argc, char* argv[]) {
 
   int opt;
 
-  while ((opt = getopt(argc, argv, "d:b:c:n:t:N:M:e:")) != -1) {
+  while ((opt = getopt(argc, argv, "k:b:c:n:t:N:M:e:")) != -1) {
     switch (opt) {
-      case 'd':
-        n = atoi(optarg);
-        break;
-      case 'b':
-        β0 = atof(optarg);
-        break;
-      case 'c':
-        βf = atof(optarg);
-        break;
-      case 'e':
-        ε = atof(optarg);
-        break;
-      case 'n':
-        num = atof(optarg);
-        break;
-      case 't':
-        trials = (unsigned)atof(optarg);
-        break;
-      case 'N':
-        N = (unsigned)atof(optarg);
-        break;
-      case 'M':
-        M = (unsigned)atof(optarg);
-        break;
-      default:
-        exit(1);
+    case 'k':
+      k = atoi(optarg);
+      break;
+    case 'b':
+      β0 = atof(optarg);
+      break;
+    case 'c':
+      βf = atof(optarg);
+      break;
+    case 'e':
+      ε = atof(optarg);
+      break;
+    case 'n':
+      num = atof(optarg);
+      break;
+    case 't':
+      n_tuning = (unsigned)atof(optarg);
+      break;
+    case 'N':
+      N = (unsigned)atof(optarg);
+      break;
+    case 'M':
+      M = (unsigned)atof(optarg);
+      break;
+    default:
+      exit(1);
     }
   }
 
+  unsigned n = pow(2, k);
+
   std::vector<Measurement*> As(num);
   for (Measurement*& A : As) {
     A = new MeasureEnergy();
@@ -98,13 +97,19 @@ int main(int argc, char* argv[]) {
 
   PT p(β0, βf, num, n, B, As);
 
-  std::cout << "Beginning " << trials << " tuning steps for " << n << ".\n";
-  p.tune(trials, ε);
-  std::cout << "Finished tuning, beginning " << N << " tempering updates of " << M << " sweeps each.\n";
+  for (MCMC& sim : p.Ms) {
+    sim.M = walsh(k);
+  }
+
+  std::cout << "Beginning " << n_tuning << " tuning steps for " << n << ".\n";
+  p.tune(n_tuning, ε);
+  std::cout << "Finished tuning, beginning " << N << " tempering updates of " << M
+            << " sweeps each.\n";
   p.run(N, M);
   std::cout << "Finished " << n << ".\n";
 
-  std::string filename = "probs_" + std::to_string(n) + "_" + std::to_string(β0) + "_" + std::to_string(βf) + "_" + std::to_string(num) +  ".dat";
+  std::string filename = "probs_" + std::to_string(n) + "_" + std::to_string(β0) + "_" +
+                         std::to_string(βf) + "_" + std::to_string(num) + ".dat";
   std::ifstream file(filename);
 
   unsigned N_old = 0;
@@ -141,7 +146,8 @@ int main(int argc, char* argv[]) {
 
   file_out.close();
 
-  std::string efilename = "energies_" + std::to_string(n) + "_" + std::to_string(β0) + "_" + std::to_string(βf) + "_" + std::to_string(num) +  ".dat";
+  std::string efilename = "energies_" + std::to_string(n) + "_" + std::to_string(β0) + "_" +
+                          std::to_string(βf) + "_" + std::to_string(num) + ".dat";
   std::ifstream efile(efilename);
 
   unsigned Ne_old = 0;