Dear ITensor Team,
Thanks for your reply. Recently, I want to repeat entanglement entropy calculation as Fig.1 in Phys. Rev. Lett. 110, 115701 by using ITensor v3, which is just simulate 1/4 doped hubbard chain. I followed the tutorial http://itensor.org/docs.cgi?vers=cppv3&page=formulas/entanglement_mps
But i find my results(entropy and energy) are not converge ( we get different results every time )when i increase bond dimension. Is there any bug in my naive code? Please help me. Many thanks!
There is my naive code:
int main()
{
int Nx = 24;
auto N = Nx;
auto sites = Electron(N);
auto t = 1.0;
auto ampo = AutoMPO(sites);
for(int j=1; j<N; ++j)
{
ampo += -t, "Cdagup",j,"Cup",j+1;
ampo += -t, "Cdagdn",j,"Cdn",j+1;
}
for(int j=1; j<=N;++j)
{
ampo +=UU, "Nupdn",j;
}
auto H = toMPO(ampo);
auto state = InitState(sites,"Emp");
// for this part, naively, i want to doped my hubbard chain, is that correct?
for(int i=1; i<=Nx/4; ++i)
{
state.set(i,"Up");
}
for(int i=Nx/4+1; i<=Nx/2; ++i)
{
state.set(i,"Dn");
}
auto sweeps = Sweeps(10);
sweeps.maxdim() = 10, 20, 50,500;
sweeps.noise() = 1E-7, 1E-8, 1E-10, 0;
sweeps.cutoff() = 1E-10;
PrintData(sweeps);
auto psi0 = MPS(state);
auto [energy,psi] = dmrg(H,psi0,sweeps,{"Quiet=",true});
PrintData(Nx);
PrintData(UU);
PrintData(t);
PrintData(totalQN(psi));
PrintData(maxLinkDim(psi));
PrintData(energy);
//Given an MPS called "psi",
//and some particular bond "b" (1 <= b < length(psi))
//across which we want to compute the von Neumann entanglement
auto b = Nx/2;
//"Gauge" the MPS to site b
psi.position(b);
//SVD this wavefunction to get the spectrum
//of density-matrix eigenvalues
auto l = leftLinkIndex(psi,b);
auto s = siteIndex(psi,b);
auto [U,S,V] = svd(psi(b),{l,s});
auto u = commonIndex(U,S);
//Apply von Neumann formula
//to the squares of the singular values
Real SvN = 0.;
for(auto n : range1(dim(u)))
{
auto Sn = elt(S,n,n);
auto p = sqr(Sn);
if(p > 1E-12) SvN += -plog(p);
}
printfln("Across bond b=%d, SvN = %.10f",b,SvN-2log(Nx)/3.0);
printfln("Across bond b=%d, SvN = %.10f",b,SvN);
return 0;
}