1 files changed, 32 insertions, 2 deletions
diff --git a/langevin.cpp b/langevin.cpp
index 870879b..5ca2d72 100644
--- a/langevin.cpp
+++ b/langevin.cpp
@@ -83,11 +83,41 @@ int main(int argc, char* argv[]) {
   ComplexVector z0 = normalize(randomVector<Complex>(N, d, r.engine()));
 
   ComplexVector zSaddle = findSaddle(J, z0, ε);
-  ComplexVector zSaddlePrev = ComplexVector::Zero(N);
   ComplexVector z = zSaddle;
 
+  std::function<double(const ComplexTensor&, const ComplexVector&)> energyNormGrad = []
+    (const ComplexTensor& J, const ComplexVector& z) {
+      double W;
+      std::tie(W, std::ignore) = WdW(J, z);
+      return W;
+    };
+
+  double aGoal = 1e3;
+
+  std::function<double(const ComplexTensor&, const ComplexVector&)> energyInvA = [aGoal]
+    (const ComplexTensor& J, const ComplexVector& z) {
+      double a = z.squaredNorm();
+      if (a > aGoal) {
+        return -aGoal;
+      } else {
+        return -z.squaredNorm();
+      }
+    };
+
+  while (zSaddle.squaredNorm() < aGoal) {
+    std::tie(std::ignore, z) = metropolis(J, z, energyInvA, T, γ, 100, d, r.engine());
+    try {
+      std::cerr << "Starting descent from " << z.squaredNorm() << "." << std::endl;
+      zSaddle = findSaddle(J, z, ε);
+    } catch (std::exception& e) {
+    }
+    std::cerr << "Current saddle is of size " << zSaddle.squaredNorm() << "." << std::endl;
+  }
+
+  ComplexVector zSaddlePrev = ComplexVector::Zero(N);
+
   while (δ < (zSaddle - zSaddlePrev).norm()) { // Until we find two saddles sufficiently close...
-    std::tie(std::ignore, z) = langevin(J, z, T, γ, M, d, r.engine());
+    std::tie(std::ignore, z) = metropolis(J, z, energyNormGrad, T, γ, M, d, r.engine());
     try {
       ComplexVector zSaddleNext = findSaddle(J, z, ε);
       if (Δ < (zSaddleNext - zSaddle).norm()) { // Ensure we are finding distinct saddles.