Hi miles

I find one example.The hamiltonian is H=He+Hph+He-ph

$$

H_{e}=-t\sum_{i,\sigma}(c^{\dag}_{i\sigma}+h.c.)+U\sum_{i}n_{i\uparrow}n_{i\downarrow} \\

Hph=\Omega b^{dag}b \\

He-ph=-\frac{g}{\sqrt(N)}(b^{\dag}+b)\sum_{i,\sigma}(-1)^{i}n_{i\sigma}

$$

And I write a code. When I execute the code, the energy of every iteration is not always descending. It will be fluctuating, as -10,-11,-10,-9.... I try different parameters. However the result is same.

#include "itensor/all.h"

using namespace itensor;

using Holstein = MixedSiteSet<ElectronSite,BosonSite>;

int

main()

{

auto N = 20;

auto U = 0.5;

auto J = 1.0;

auto g = 1.0;

auto sigma = 1.0;

auto sites = Holstein(N,{"ConserveNf=",true,

"ConserveNb=",false,

"MaxOcc=",50});

auto ampo = AutoMPO(sites);

for(int j=1;j <= N-2; j+=2)

{

ampo += -J,"Cdagup",j+2,"Cup",j;

ampo += -J,"Cdagup",j,"Cup",j+2;

ampo += -J,"Cdagdn",j+2,"Cdn",j;

ampo += -J,"Cdagdn",j,"Cdn",j+2;

}

ampo += -J,"Cdagup",1,"Cup",19;

ampo += -J,"Cdagup",19,"Cup",1;

ampo += -J,"Cdagdn",1,"Cdn",19;

ampo += -J,"Cdagdn",19,"Cdn",1;

for(int j=1;j < N; j += 2)

{

ampo += U,"Nup",j,"Ndn",j;

}

for(int j=1;j<N;j +=4)

{

ampo += g,"Adag",2,"Nup",j;

ampo += g,"Adag",2,"Ndn",j;

}

for(int j=3;j<N;j +=4)

{

ampo += -g,"Adag",2,"Nup",j;

ampo += -g,"Adag",2,"Ndn",j;

}

ampo += sigma,"Adag",2,"A",2;

auto H = toMPO(ampo);

auto sweeps = Sweeps(20);

//very important to use noise for this model

sweeps.noise() = 1E-6,1E-6,1E-8,1E-12;

sweeps.maxdim() = 10,20,100,100,200;

sweeps.cutoff() = 1E-10;

auto state = InitState(sites);

for(auto j : range1(N))

{

if(j%4==1) state.set(j,"Up");

else if(j%4==3) state.set(j,"Dn");

else state.set(j,"0");

}

auto psi0 = MPS(state);

auto [energy,psi] = dmrg(H,psi0,sweeps,{"Quiet=",true});

printfln("Ground State Energy = %.12f",energy);