# [C++ ITensor 3] MPO time evolution does not work for large interactions

+1 vote

Hi,

It seems the MPO time evolution gives wrong results for large interactions. The sample code is below.

#include "itensor/all.h"
using namespace std;
using namespace itensor;

int
main(int argc, char* argv[])
{
int N = 10;
auto sites = Boson(N,{"MaxOcc=",1,"ConserveQNs=",false});
auto state = InitState(sites);
for(int i = 1; i <= N; ++i)
{
state.set(i,"0");
}
auto psi0 = MPS(state);
auto psi = MPS(state);

auto ampo = AutoMPO(sites);
//Make the Hamiltonian
for(int b = 1; b < N; ++b)
{
ampo += 10.0,"N",b,"N",b+1;
}
for(int b = 1; b <= N; ++b)
{
ampo += 0.5,"A",b;
ampo += -2.0,"N",b;
}
auto H = toMPO(ampo);
auto tau = 0.1;
//auto expH = toExpH(ampo,tau);
auto expH = toExpH(ampo,tau*Cplx_i);

auto args = Args("Method=","DensityMatrix","Cutoff=",1E-9,"MaxDim=",3000);
auto ttotal = 3.0;
auto nt = int(ttotal/tau+(1e-9*(ttotal/tau)));

for(int n = 1; n <= nt; ++n)
{
psi = applyMPO(expH,psi,args);
psi.noPrime().normalize();
auto ovlp = norm(innerC(psi,psi0));
auto bb = n*tau;
printfln("\nOverlap at time %.4f %.20f",bb,ovlp);
}
return 0;
}


When the NN interaction is large, like 10.0, the overlap between the time-evolved psi and the initial psi will quickly become "0.00000000000000000000" and then keeps zero forever. This makes no sense. Is this due to the algorithm itself or some bug in the code? Thanks.

Jin

commented by (550 points)
edited by
commented by (1.2k points)
Thanks JR. That works. But I think I did not bring out my real concern. If the NN interaction is very large, tau needs to be very small and the calculation will be very slow to reach the target final time. I was looking for an efficient solution for large interaction. My naive idea is to delete NN interaction term and instead putting a projector onto zero states on the neighboring sites of Adag and A. This results in a 3-site operator which is not allowed in the built-in MPO time evolution method.
commented by (550 points)
Yes. that’s a problem. It is hard to reach the large time scale. please pay attention to the half-block entanglement entropy of system(sometimes  it increases ~lnt or ~t). If it is large, discard weight would be large, and you need keeping much more states.           The current C++ version ITensors does not support three-site interaction, and  you can try the Julia verison.
cheers!

+1 vote

The NN interaction is large and time step tau=0.1. The real time step tau*10=1, is very large for MPO time evolution method. You can try tau=0.01 or smaller.
Best wishes!

+1 vote

Hi Jin,
I believe JR is correct that a smaller tau could help here, but of course it would be nice if you can use a method that allows a larger tau.

Have you thought about using the TDVP method instead? It has all the advantages of the MPO approach you are using, but over the last handful of years has turned out to be much superior in terms of accuracy and ability to reliably use larger time steps.

There is a high-quality TDVP code available for ITensor through this code repository:
https://github.com/ITensor/TDVP
and it also includes a sample code to help you get started.

There's one caution about TDVP, as you may know, that it historically has trouble 'adapting' the bond dimension automatically, either in the sense that it doesn't do so at all (1-site TDVP) or in the sense that it does it in a sub-optimal way (2-site TDVP). But the code linked above includes a novel 'basis expansion' method which you'll see used in the sample that improves the MPS basis to alleviate problems with the growth of the bond dimension in a proper way.

As with any time evolution code, I'd recommend comparing to exact results either for small systems or in limits where an exact solution may be available.

Best regards,
Miles

commented by (1.2k points)
Hi Miles, thanks for your answer. The TDVP algorithm looks great. The sample code only adds basis in the first two steps. Then the bond-dimension does not change after that. For longer time evolution, the entanglement entropy grows, should we add basis every a few times? How do we keep the bond dimension optimal (the most economy MPS representation) during the time evolution?
commented by (70.1k points)
Hi, so I have not actually used the TDVP method myself a whole lot recently, so please cross-check my answers with other literature and by testing things on small systems against other parameter values and ideally other methods like using Trotter gates. But here are my thoughts:

1.  the sample code might only be doing the basis expansion in early steps because it is doing imaginary time evolution (I think, just from glancing at it). For imaginary time evolution, after a long enough time the bond dimension saturates to the ground-state value. For real-time evolution it may be important to do the basis expansion more often.

2. the basis expansion is only strictly necessary when doing 1-site TDVP which cannot grow the bond dimension.

3. there is a two-site TDVP algorithm which is available by setting "NumCenter" equal to 2. Two-site TDVP does adapt the bond dimension, which you control by setting a truncation error cutoff and maxdim, very similarly to two-site DMRG.

4. the relative benefits and drawbacks of 1-site versus 2-site TDVP are complicated. As mentioned, two-site TDVP is basis adaptive so it would seem better but it is slower than one-site TDVP, and more importantly, might be less accurate in some cases (for a given bond dimension). See the following paper for more information: https://arxiv.org/abs/1809.01400

I hope to eventually put some systematic tests of TDVP and Trotter time evolution on the ITensor website so we can make our own recommendations of when to use which method and for what parameters.

For some recent work I was involved in using METTS for the Hubbard model, starting off with Trotter gates (TEBD) for the initial imaginary time steps, then swtiching to two-site TDVP for a while and finally one-site TDVP worked very well in terms of balancing accuracy and speed.

Best,
Miles