# Time evolution of local green function using two applyMPO for spinless-fermions.

+1 vote
edited

Hi,
I am trying to calculate the local green function Gii(t)= <\psi0|Cdagi(t)Ci|\psi0> at site i=N, using two applyMPO for one-dimensional Kitaev model [equation 1 of PRB 88, 161103(R) (2013)]. For this I am doing following steps:
(i) I am calculating |\phi>=C
i |\psi0> by applying local operator Ci to the ground state |psi0>.
(ii) Time evolution of |\phi> using applyMPO: |phi(t)>=exp(-iHt)|\phi>.
(iii) Application of local operator Cdag
i to |\phi(t): Cdagi|phi(t)>.
(iv)Time evolution of |psi
0> using applyMPO: |psi(t)>=exp(-iHt)|\psi_0>.
(v) Overlap of <psi(t)|phi(t)> using innerC.

After running the code, I am getting non-zero values for the odd time steps and zero for even number of time steps:
0.001 -0.000605309 -0.671164
0.002 0 -0
0.003 0.00300302 0.671158
0.004 0 -0
0.005 -0.00540069 -0.671143

Here is the code:

//-------------Creating |phi> = C_i1|Psig>--------

 auto i1=N;  // for site i1=N.
psig.position(i1);
auto newpsi = noPrime(psig(i1)*op(sites,"A",i1));
psig.set(N, newpsi);


//--------------Jordan-Wigner string----------------------

 for(int k = i1-1; k >=1; k--)
{
psig.position(k);
auto newpsi1 = noPrime(psig(k)*op(sites,"F",k));
psig.set(k, newpsi1);
}
psig.noPrime().normalize();


auto tau=0.001;
auto ii = Complexi;
auto args = Args("Method=","DensityMatrix","Cutoff=",1E-14,"MaxDim=",7000);
auto expH = toExpH(ampo,tstep*Cplx
i);
//-------------Time Evolution--------------------------

auto nt = int(ttotal/tau+(1e-9*(ttotal/tau)));
for(int n = 1; n <= nt;++n)
{
psig = applyMPO(expH,psig,args);
psig.noPrime().normalize();


//--Jordan-Wigner string

 for(int k = 1; k <i1; ++k)
{
psig.position(k);
auto newpsi2 = noPrime(psig(k)*op(sites,"F",k));
psig.set(k, newpsi2);
}


// ----------------- Cdag|phi(t)>--------------------

    psig.position(i1);
psig.set(i1, newpsi3);
psig.noPrime().normalize();


//------------------psi(t)>=exp(-iHt)|psig>-----------------------

      psi = applyMPO(expH,psi,args);
psi.noPrime().normalize();
auto result2 = -ii*innerC(psi,psig);
file2<<float(n*tstep)<<' '<< result2.real()<<' ' <<result2.imag() <<std::endl;
`

+1 vote

Hi, thanks for the good question and your patience with a slow reply.

While I'm not sure why you getting the value of zero on every other time step, I think there are a few major improvements to your approach which you should consider anyway, and after you make them this problem may disappear in the process of changing your code.

Here are the improvements or fixes:
1. since one of the exp(-i Ht) operators acts on the ground state of H, you can rigorously replace it by a scalar exp(-i E0 t) where E0 is the ground state energy <psi0|H|psi0> and so skip that step

1. a much better choice for doing the time evolution would be the TDVP algorithm, for which a high-quality implementation in ITensor is now available here:
https://github.com/ITensor/TDVP