1 files changed, 0 insertions, 142 deletions
diff --git a/langevin.cpp b/langevin.cpp
deleted file mode 100644
index ea28d65..0000000
--- a/langevin.cpp
+++ /dev/null
@@ -1,142 +0,0 @@
-#include <getopt.h>
-#include <limits>
-#include <unordered_map>
-#include <list>
-
-#include "Eigen/Dense"
-#include "Eigen/src/Eigenvalues/ComplexEigenSolver.h"
-#include "Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h"
-#include "complex_normal.hpp"
-#include "p-spin.hpp"
-#include "dynamics.hpp"
-#include "stokes.hpp"
-
-#include "collectStokesData.hpp"
-
-#include "pcg-cpp/include/pcg_random.hpp"
-#include "randutils/randutils.hpp"
-#include "unsupported/Eigen/CXX11/src/Tensor/TensorFFT.h"
-
-#define PSPIN_P 3
-const int p = PSPIN_P; // polynomial degree of Hamiltonian
-using Complex = std::complex<Real>;
-using RealVector = Vector<Real>;
-using RealMatrix = Matrix<Real>;
-using RealTensor = Tensor<Real, p>;
-using ComplexVector = Vector<Complex>;
-using ComplexMatrix = Matrix<Complex>;
-using ComplexTensor = Tensor<Complex, p>;
-
-using Rng = randutils::random_generator<pcg32>;
-
-int main(int argc, char* argv[]) {
-  // model parameters
-  unsigned N = 10; // number of spins
-  Real T = 1;    // temperature
-  Real Rκ = 0;   // real part of distribution parameter
-  Real Iκ = 0;   // imaginary part of distribution parameter
-  // simulation parameters
-  Real ε = 1e-15;
-  Real εJ = 1e-5;
-  Real δ = 1;   // threshold for determining saddle
-  Real Δ = 1e-3;
-  Real γ = 1e-1;   // step size
-  unsigned t = 1000; // number of Langevin steps
-  unsigned M = 100;
-  unsigned n = 100;
-
-  int opt;
-
-  while ((opt = getopt(argc, argv, "N:M:n:T:e:r:i:g:t:E:")) != -1) {
-    switch (opt) {
-    case 'N':
-      N = (unsigned)atof(optarg);
-      break;
-    case 't':
-      t = (unsigned)atof(optarg);
-      break;
-    case 'T':
-      T = atof(optarg);
-      break;
-    case 'e':
-      δ = atof(optarg);
-      break;
-    case 'E':
-      ε = atof(optarg);
-      break;
-    case 'g':
-      γ = atof(optarg);
-    case 'r':
-      Rκ = atof(optarg);
-      break;
-    case 'i':
-      Iκ = atof(optarg);
-      break;
-    case 'n':
-      n = (unsigned)atof(optarg);
-      break;
-    case 'M':
-      M = (unsigned)atof(optarg);
-      break;
-    default:
-      exit(1);
-    }
-  }
-
-  Rng r;
-
-  collectStokesData<2, 4>(N, r.engine(), 1e-15, 1.0, 0.01);
-
-  pSpinModel<Real, 2, 4>pSpin(N, r.engine(), 1, 0.01);
-
-  std::normal_distribution<Real> Red(0, 1);
-
-  RealVector zMin = randomMinimum(pSpin, Red, r.engine(), ε);
-  Real Hr;
-  RealVector dHr;
-  RealMatrix ddHr;
-  std::tie(Hr, dHr, ddHr, std::ignore) = pSpin.hamGradHess(zMin);
-  RealMatrix M1 = ddHr - (dHr * zMin.transpose() + zMin.dot(dHr) * Matrix<Real>::Identity(zMin.size(), zMin.size()) + (ddHr * zMin) * zMin.transpose()) / (Real)zMin.size();
-  RealMatrix M2 = ddHr - zMin.dot(dHr) * Matrix<Real>::Identity(zMin.size(), zMin.size()) / Real(N);
-  Eigen::EigenSolver<Matrix<Real>> eigenS(M2);
-  for (unsigned i = 0; i < N; i++) {
-    RealVector zNew = normalize(zMin + 1e-5 * eigenS.eigenvectors().col(i).real());
-    std::cout << eigenS.eigenvalues()(i) << " " << 2.0 * (pSpin.getHamiltonian(zNew) - Hr) / 1e-10 << " " << real(eigenS.eigenvectors().col(i).dot(zMin)) << " " << eigenS.eigenvectors().col(0).dot(eigenS.eigenvectors().col(i)) << std::endl;
-  }
-  std::cout << std::endl;
-  getchar();
-
-  complex_normal_distribution<Real> d(0, 1, 0);
-
-  pSpinModel complexPSpin = pSpin.cast<Complex>();;
-  ComplexVector zSaddle = zMin.cast<Complex>();
-
-  ComplexVector zSaddleNext;
-  bool foundSaddle = false;
-  while (!foundSaddle) {
-    ComplexVector z0 = normalize(zSaddle + δ * randomVector<Complex>(N, d, r.engine()));
-    try {
-      zSaddleNext = findSaddle(complexPSpin, z0, ε);
-      Real saddleDistance = (zSaddleNext - zSaddle).norm();
-      if (saddleDistance / N > 1e-2) {
-        std::cout << saddleDistance << std::endl;
-        foundSaddle = true;
-      }
-    } catch (std::exception& e) {}
-  }
-
-  Complex H1 = complexPSpin.getHamiltonian(zSaddle);
-  Complex H2 = complexPSpin.getHamiltonian(zSaddleNext);
-
-  Real φ = atan2( H2.imag() - H1.imag(), H1.real() - H2.real());
-  std::cerr << (zSaddle - zSaddleNext).norm() / (Real)N << " " << φ << " " << H1 * exp(Complex(0, φ)) << " " << H2 * exp(Complex(0, φ)) << std::endl;
-
-  complexPSpin *= exp(Complex(0, φ));
-
-  Cord test(complexPSpin, zSaddle, zSaddleNext, 3);
-  test.relaxNewton(10, 1, 1e4);
-
-  std::cout << test.totalCost(10) << " " << test.totalCost(20) << std::endl;
-
-  return 0;
-}