# AutoMPO for fermions and DMRG

Hi,

I try creating a Hubbard model Hamiltonian in the two different ways using AutoMPO:

In the first case, I used Jordan-Wigner strings:

    auto ampoGS = AutoMPO(sites);
ampoGS += -1.0,"Adagup",b,"F",b,"Aup",b+1; //Hopping for spin up


In the second case, I used fermion operators assuming AutoMPO handles this fine.

    auto ampo = AutoMPO(sites);
ampo += -1.0,"Cdagup",b,"Cup",b+1; //Hopping for spin up
ampo += -1.0,"Cdagdn",b,"Cdn",b+1; //Hopping for spin down
ampo += -1.0,"Cdn",b,"Cdagdn",b+1;


Next, I found wildly different results performing the DMRG operation on the Hamiltonian constructed using HGS = IQMPO(ampoGS) vs. H = IQMPO(ampo). For free fermions at half-filling, I checked that the DMRG with HGS gives the correct results, while H does not. For the interacting case with Hubbard U, I found that dmrg method on H gives results that vary from run to run?

Ultimately, I want to do time evolution, and for this I need toExpH which flags an error when working with HGS because it contains more than "2 operator" terms. (And I can't trust "H" but maybe I am making a simple error working with AutoMPO?)

Any help would be really appreciated! Thanks!

my complete code is:

auto sites = Hubbard(N);
auto state = InitState(sites);
for(int i = 1; i <= N; ++i)
{
if(i%2 == 1)
state.set(i,"Up");
else
state.set(i,"Dn");
}
auto psioriginal = IQMPS(state);

// Create the Target Hamiltonian and find the Ground State Correlations

auto ampo = AutoMPO(sites);
auto ampoGS = AutoMPO(sites);

auto ampoLocalEnergy = AutoMPO(sites);

//AutoMPO version with C operators doesn't work correctly. Tested for ground state energy at half-filling expecting -1.27 and got something else.

for (int b = 1; b < N; b++)
{

ampoGS += -1.0,"Adagup",b,"F",b,"Aup",b+1; //Hopping for spin up

ampo += -1.0,"Cdagup",b,"Cup",b+1; //Hopping for spin up
ampo += -1.0,"Cdagdn",b,"Cdn",b+1; //Hopping for spin down
ampo += -1.0,"Cdn",b,"Cdagdn",b+1;

ampoGS += U,"Nupdn",b; //Interactions between upspin/downspin electrons on the same site
ampo += U,"Nupdn",b; //Interactions between upspin/downspin electrons on the same site

}
ampoGS += U,"Nupdn",N;
ampo += U,"Nupdn",N;

auto HGS = IQMPO(ampoGS);

auto sweeps = Sweeps(5); //number of sweeps is 5
sweeps.maxm() = 10,20,100,100,200; //gradually increase states kept
sweeps.cutoff() = 1E-10; //desired truncation error
double Finalgroundstateenergydensity = dmrg(psioriginal,HGS,sweeps,{"Quiet=",true})/N; // psioriginal is now the ground state of the System.


/* auto H = IQMPO(ampo);

// auto sweeps = Sweeps(5); //number of sweeps is 5
// sweeps.maxm() = 10,20,100,100,200; //gradually increase states kept
// sweeps.cutoff() = 1E-10; //desired truncation error
// double Finalgroundstateenergydensity = dmrg(psioriginal,H,sweeps,{"Quiet=",true})/N; // psioriginal is now the ground state of the System.

std::cout << Finalgroundstateenergydensity << std::endl;


selected

Hi, it looks like your Hamiltonian is not Hermitian for ampo. That's probably why its different every time you run it. Assuming you directly copied and pasted:

ampo += -1.0,"Cup",b,"Adagup",b+1;


is the offending line. The Adagup should be Cdagup

Once you fix that, they still don't agree and that's because you didn't write your C operators in normal ordering, so they are off by a commutator I think? When I do

ampo += -1.0,"Cdagup",b,"Cup",b+1; //Hopping for spin up
ampo += -1.0,"Cdagup",b+1,"Cup",b; // <-
ampo += -1.0,"Cdagdn",b,"Cdn",b+1; //Hopping for spin down
ampo += -1.0,"Cdagdn",b+1,"Cdn",b;// <-


I get the same value for both cases using this change

commented by (160 points)
Whoops, should have definitely caught these myself. :O Thanks for your help!