SHARE
TWEET

Untitled

a guest Dec 2nd, 2018 113 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. #include "itensor/all.h"
  2. #include <iostream>
  3. #include <iomanip>
  4.  
  5. using namespace itensor;
  6.  
  7. std::vector<double> test1(int argc, char* argv[], Real lambda){
  8.  
  9.   int N = 50;
  10.   auto sweeps = Sweeps(5);
  11.   sweeps.maxm() = 50,80,100,150,150;
  12.   sweeps.minm() = 10,10,30,50,90;
  13.   sweeps.cutoff() = 1E-6,1E-8,1E-10,1E-12,1E-12;
  14.   sweeps.niter() = 4,3,2,2,2;
  15.   sweeps.noise() = 1E-7,1E-8,1E-10,0,0;
  16.   /*auto sweeps = Sweeps(12);
  17.   sweeps.maxm() = 10,10,10,20,30,50,80,100,150,150,200,250;
  18.   sweeps.minm() = 10,10,10,10,10,10,10,30,50,90,150,200;
  19.   sweeps.cutoff() = 1E-6,1E-6,1E-6,1E-6,1E-6,1E-6,1E-6,1E-8,1E-10,1E-12,1E-12,1E-14,1E-15;
  20.   sweeps.niter() = 4,4,4,4,4,4,3,2,2,2,2,2;
  21.   sweeps.noise() = 1E-7,1E-7,1E-7,1E-7,1E-7,1E-7,1E-8,1E-10,0,0,0,0;*/
  22.  
  23.   bool use_Z3 = true;
  24.   auto sites = SpinHalf(N);
  25.   auto ampo = AutoMPO(sites);
  26.   for(int j = 1; j < N; ++j){
  27.     ampo += -4.0*lambda,"Sz",j,"Sz",j+1;
  28.   }
  29.   for(int j = 1; j <= N; ++j){
  30.     ampo += -2*(1.0-lambda),"Sx",j;
  31.   }
  32.  
  33.   // dmrg
  34.   auto H = MPO(ampo);
  35.   auto psi = MPS(sites);
  36.   auto energy = dmrg(psi,H,sweeps);
  37.  
  38.   // calculate the entanglement entropy
  39.   Real s = 0;
  40.   ITensor S, V;
  41.   int b = N/2;
  42.   psi.position(b);
  43.   ITensor wf = psi.A(b) * psi.A(b+1);
  44.   auto U = psi.A(b);
  45.   auto spectrum = svd(wf,U,S,V);
  46.   for (auto p : spectrum.eigs()){
  47.     if(p > 1E-12) s += -p*log(p);
  48.   }
  49.  
  50.   auto psi1 = MPS(sites);
  51.   auto wfs = std::vector<MPS>(1);
  52.   wfs.at(0) = psi;
  53.   auto energy1 = dmrg(psi1, H, wfs, sweeps, {"Quiet=", true, "Weight=", 20.0});
  54.  
  55.   return std::vector<double>{s, energy, energy1};
  56. }
  57.  
  58. int main(int argc, char* argv[]){
  59.   int N_parameter_sets = 15;
  60.   double lambda_min = 0.01;
  61.   double lambda_max = 0.99;
  62.   double lambda_grad = (lambda_max - lambda_min) / (N_parameter_sets - 1);
  63.  
  64.   std::vector<double> lambdas;
  65.   for (int i = 0; i < N_parameter_sets; i++){
  66.     double lambda = lambda_min + lambda_grad * i;
  67.     lambdas.push_back(lambda);
  68.   }
  69.  
  70.   std::vector<std::vector<double> > results;
  71.   for (auto lambda: lambdas){
  72.     auto result = std::vector<double>();
  73.     result.push_back(lambda);
  74.     for (auto x : test1(argc, argv, lambda)){
  75.       result.push_back(x);
  76.     }
  77.     results.push_back(result);
  78.   }
  79.  
  80.   for (auto r : results){
  81.     std::stringstream ss;
  82.     for (auto y : r){
  83.       //ss << std::fixed << std::setprecision(12) << y << ", ";
  84.       ss << std::setw(12) << std::fixed << std::setprecision(12) << y << ", ";
  85.     }
  86.     /* Output lambda, entropy, E_ground, E_excited */
  87.     std::cout << ss.str() << std::endl;
  88.   }
  89.   return 0;
  90. }
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
 
Top