From 22f0e9aee93c724209ee0c8b0d079b0d191e767c Mon Sep 17 00:00:00 2001
From: Jaron Kent-Dobias <jaron@kent-dobias.com>
Date: Thu, 12 Dec 2019 09:02:51 -0500
Subject: added a stage for burnout steps to address superheating

---
 hadamard_mcmc.hpp | 22 +++++++++++++---------
 1 file changed, 13 insertions(+), 9 deletions(-)

(limited to 'hadamard_mcmc.hpp')

diff --git a/hadamard_mcmc.hpp b/hadamard_mcmc.hpp
index 5ef062f..9f0d8c0 100644
--- a/hadamard_mcmc.hpp
+++ b/hadamard_mcmc.hpp
@@ -156,7 +156,7 @@ public:
     E = M.energy();
   }
 
-  bool step(Givens& g) {
+  bool step(Givens& g, bool dry = false) {
     double ΔE = g.tryRotation();
 
     bool accepted = ΔE < 0 || exp(-β * ΔE) > rng.uniform((double)0.0, 1.0);
@@ -166,11 +166,13 @@ public:
       g.acceptRotation();
     }
 
-    A.after_step(accepted, g, θ0, E, ΔE, M);
+    if (!dry) {
+      A.after_step(accepted, g, θ0, E, ΔE, M);
+    }
     return accepted;
   }
 
-  double sweep() {
+  double sweep(bool dry = false) {
     unsigned total = 0;
     unsigned accepted = 0;
 
@@ -179,9 +181,9 @@ public:
         Givens g1(M, false, i, j, θ0, rng);
         Givens g2(M, true, i, j, θ0, rng);
 
-        if (this->step(g1))
+        if (this->step(g1, dry))
           accepted++;
-        if (this->step(g2))
+        if (this->step(g2, dry))
           accepted++;
 
         total += 2;
@@ -193,7 +195,7 @@ public:
 
   void tune(unsigned N, double ε) {
     for (unsigned i = 0; i < N; i++) {
-      double ratio_accepted = this->sweep();
+      double ratio_accepted = this->sweep(true);
       if (ratio_accepted > 0.5) {
         θ0 *= 1 + ε;
       } else {
@@ -202,10 +204,12 @@ public:
     }
   }
 
-  void run(unsigned N) {
+  void run(unsigned N, bool dry = false) {
     for (unsigned i = 0; i < N; i++) {
-      this->sweep();
-      A.after_sweep(θ0, E, M);
+      this->sweep(dry);
+      if (!dry) {
+        A.after_sweep(θ0, E, M);
+      }
     }
   }
 };
-- 
cgit v1.2.3-70-g09d2