0 votes
asked by (180 points)
edited by

I am comparing the two methods for time-evolution.

First, I write the Hamiltonian in MPO form using autoMPO, and then time evolve using toExpH and applyMPO.

Secondly, I implement a TEBD type algorithm by creating a sequence of two-site gates using BondGate and time evolving with gateTEvol.

However, both give very different results for the dynamics for a Fermi-Hubbard model. I have tested this for a spin-1/2 chain and this works fine, I am wondering if the trotter gate method does not properly take into account the fermionic statistics?

Here is sample code:

#include "itensor/all.h"

using namespace itensor;

int main(int argc, char* argv[])
{
    int BD = 1024;

    int t1 = 1;

    int N = 20;
    auto sites = Electron(N);
    auto state = InitState(sites);
    int p=1;
    for(int i = N; i >= 1; --i)
    {
        if(p == 1)
        {
            state.set(i,"Up");
            p = 2;
        }
        else
        {
            state.set(i,"Dn");
            p = 1;
        }
    }
    auto psi = MPS(state);

    //
    // Create the Hamiltonian using AutoMPO
    //
    auto ampo = AutoMPO(sites);
    for(int b = 1; b < N; ++b)
    {
        ampo += -t1,"Cdagup",b,"Cup",b+1;
        ampo += -t1,"Cup",b,"Cdagup",b+1;
        ampo += -t1,"Cdagdn",b,"Cdn",b+1;
        ampo += -t1,"Cdn",b,"Cdagdn",b+1;
    }
    auto H = toMPO(ampo);
    auto tau = 0.01;
    auto expH = toExpH(ampo,tau*Cplx_i);


    auto ttotal = 1.0;
    auto nt = int(ttotal/tau+(1e-9*(ttotal/tau)));

    std::ofstream ofGup("DATA/MPO_TE_Test.txt");

    for(int n = 1; n <= nt; ++n)
    {
        cpu_time sw_time;

        psi = applyMPO(expH,psi,{"Method=","DensityMatrix","MaxDim=",BD,"Cutoff=",1E-15});
        psi.noPrime().normalize();

        double updo;

        for(int j = 1; j <= N; ++j)
        {
            psi.position(j);
            updo = (eltC(dag(prime(psi(j),"Site"))*op(sites,"Nup",j)*psi(j))).real();
            printfln(ofGup,"%d ",updo);
        }

        printfln("MPO: Time %d",n*tau);
        printfln("Maximum MPS bond dimension is %d",maxLinkDim(psi));

        auto sm = sw_time.sincemark();
        printfln("CPU time = %s ",showtime(sm.time));
    }



    auto sites2 = Electron(N);
    auto state2 = InitState(sites2);
    p=1;
    for(int i = N; i >= 1; --i)
    {
        if(p == 1)
        {
            state2.set(i,"Up");
            p = 2;
        }
        else
        {
            state2.set(i,"Dn");
            p = 1;
        }
    }
    auto psi_TEBD = MPS(state2);

    Print(totalQN(psi_TEBD));

    auto gates = vector<BondGate>();

    //Create the gates exp(-i*tstep/2*hterm)
    //and add them to gates
    for(int b = 1; b <= N-1; ++b)
    {
        auto hterm = -t1*op(sites2,"Cdagup",b)*op(sites2,"Cup",b+1);
        hterm += -t1*op(sites2,"Cup",b)*op(sites2,"Cdagup",b+1);
        hterm += -t1*op(sites2,"Cdagdn",b)*op(sites2,"Cdn",b+1);
        hterm += -t1*op(sites2,"Cdn",b)*op(sites2,"Cdagdn",b+1);

        auto g = BondGate(sites2,b,b+1,BondGate::tReal,tau/2.,hterm);
        gates.push_back(g);
    }
    //Create the gates exp(-i*tstep/2*hterm) in reverse
    //order (to get a second order Trotter breakup which
    //does a time step of "tstep") and add them to gates
    for(int b = N-1; b >= 1; --b)
    {
        auto hterm = -t1*op(sites2,"Cdagup",b)*op(sites2,"Cup",b+1);
        hterm += -t1*op(sites2,"Cup",b)*op(sites2,"Cdagup",b+1);
        hterm += -t1*op(sites2,"Cdagdn",b)*op(sites2,"Cdn",b+1);
        hterm += -t1*op(sites2,"Cdn",b)*op(sites2,"Cdagdn",b+1);

        auto g = BondGate(sites2,b,b+1,BondGate::tReal,tau/2.,hterm);
//        auto g = BondGate(sites,b,b+1);
        gates.push_back(g);
    }


    Real cutoff = 1E-15; //truncation error cutoff when restoring MPS form
    std::ofstream ofGup2("DATA/TEBD_TE_Test.txt");

    for(int n = 1; n <= nt; ++n)
    {
        cpu_time sw_time;

        //Time evolve, overwriting psi when done
        gateTEvol(gates,tau,tau,psi_TEBD,{"Cutoff=",cutoff,"Verbose=",true,"MaxDim=",BD});
        psi_TEBD.noPrime().normalize();

        double updo;

        for(int j = 1; j <= N; ++j)
        {
            psi_TEBD.position(j);
            updo = (eltC(dag(prime(psi_TEBD(j),"Site"))*op(sites2,"Nup",j)*psi_TEBD(j))).real();
            printfln(ofGup2,"%d ",updo);
        }

        printfln("Time %d",n*tau);
        printfln("Maximum MPS bond dimension is %d",maxLinkDim(psi_TEBD));

        auto sm = sw_time.sincemark();
        printfln("CPU time = %s ",showtime(sm.time));
    }





    return 0;
}

What have I done wrong??
Note that this is not an issue with the time-step or bond dimension.

Thanks,

Stuart

3 Answers

+1 vote
answered by (180 points)
edited by

This is still an open issue.....

0 votes
answered by (6.1k points)

Hi Stuart,

I believe the issue is that you should define your Hamiltonian as follows:

for(int b = 1; b < N; ++b)
  {
  ampo += -t1,"Cdagup",b,"Cup",b+1;
  ampo += -t1,"Cdagup",b+1,"Cup",b;
  ampo += -t1,"Cdagdn",b,"Cdn",b+1;
  ampo += -t1,"Cdagdn",b+1,"Cdn",b;
  }

to correctly account for the Fermionic sign (take a look at: https://itensor.org/docs.cgi?page=tutorials/fermions for more information).

Let me know if that helps.

-Matt

0 votes
answered by (36.6k points)

In addition to Matt's answer, which I also agree looks very likely to be the issue, I'd like to note that we have been finding that the toExpH approach to performing time evolution has not always been the most reliable unfortunately (though it depends on details of the Hamiltonian and time step size etc.). So in general you might not be able to get as good of results with that method as, say, the Trotter method for cases where both can be used.

In the near future, we are planning to offer alternative time evolution methods which will supersede the toExpH approach.

Miles

Welcome to ITensor Support Q&A, where you can ask questions and receive answers from other members of the community.

Formatting Tips:
  • To format code, indent by four spaces
  • To format inline LaTeX, surround it by @@ on both sides
  • To format LaTeX on its own line, surround it by $$ above and below
  • For LaTeX, it may be necessary to backslash-escape underscore characters to obtain proper formatting. So for example writing \sum\_i to represent a sum over i.
If you cannot register due to firewall issues (e.g. you cannot see the capcha box) please email Miles Stoudenmire to ask for an account.

To report ITensor bugs, please use the issue tracker.

Categories

...