Fix 3-site SU dt and reduce artificial C4v breaking

[Yue-Zhengyuan]

This PR aims to solve issues of the 3-site simple update discovered in #215 (credit to @ogauthe).

• Now each su3site_iter actually performs evolution over a period of dt. In particular, the exponentiation of the Hamiltonian is corrected. For each ⅃-shaped cluster

     r-1         M3 
                 | 
                 ↓ 
     r   M1 -←- M2 
         c      c+1 

the involved Hamiltonian terms are $H = (H_{12} + H_{23}) / 2 + H_{13}$ (the 1/2 factor is due to repeated counting from L-shaped clusters). Previously I mistakably exponentiated each H-term first and add them up:

$$\exp({-\epsilon H_{12}})/2 + \exp({-\epsilon H_{23}})/2 + \exp({-\epsilon H_{13}}) \quad \text{(wrong)}$$

But this is approximately

$$\begin{aligned} & \frac{1}{2}(1 - \epsilon H_{12}) + \frac{1}{2}(1 - \epsilon H_{23}) + (1 - \epsilon H_{13}) + O(\epsilon^2) \\\ &= 2 - \epsilon \left[ \frac{1}{2}(H_{12} + H_{23}) + H_{13} \right] + O(\epsilon^2) = 2 \, \exp({-\epsilon H / 2}) \end{aligned}$$

with a useless extra factor of 2 in the front, and it reduced the time step to $\epsilon / 2$. Now I directly calculate $\exp({-\epsilon H})$ for each cluster.

• The C4v symmetry is better preserved by more thorough rotation of the iPEPS during evolution. Previously, I only did rotl90 -> rotr90 in each iteration. Now I do rotl90 four times. As a result now each NN bond is being updated 4 times in each iteration, so the exponentiation of the Hamiltonian is done with alg.dt / 4.

[codecov]

Codecov Report

All modified and coverable lines are covered by tests ✅

Files with missing lines	Coverage Δ
src/algorithms/time_evolution/evoltools.jl	95.58% <100.00%> (+0.20%)	⬆️
src/algorithms/time_evolution/simpleupdate.jl	100.00% <100.00%> (ø)
src/algorithms/time_evolution/simpleupdate3site.jl	100.00% <100.00%> (ø)
src/operators/localoperator.jl	87.61% <100.00%> (+9.03%)	⬆️

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[lkdvos]

I realize I left a bunch of slightly unrelated comments here, but it might be good to fix this while we're at it. Otherwise looks good to go for me!

[lkdvos]

Update: I realized none of these comments are really relevant to the specific PR, so I'm okay with merging this as is and leaving the rest for follow-up PRs

[Yue-Zhengyuan]

I'll make some final clean up today, and we can merge afterwards.

[lkdvos]

Additional comment: it seems @ogauthe still cannot get the results to agree so maybe we should hold off for a second

[Yue-Zhengyuan]

I think I get dt wrong again... I should divide dt by 2 instead of 4. Hopefully except the size of dt there are no other serious issues.

[Yue-Zhengyuan]

I took a closer look of the J1-J2 SU + AD example. It turns out that the SU result there (using master branch) also have C4v symmetry breaking (see here). I thought I solved it in #171 before but didn't. We then need to regenerate the J1-J2 example before merging this PR.

[lkdvos]

As an idea for a follow-up PR:

I think part of the confusion here is/was facilitated by the use of LocalOperator, which represents a sum of local terms. In principle, before going to simple update, it would be nice to have an explicit (possibly configurable) exponentiate or trotterize call that takes in a LocalOperator, and converts it to a new datastructure which holds the exponentiated gates. This is not really for performance, more as a generic organizational thing to really capture the difference, which would make it slightly more obvious that you cannot sum things the same way.

Otherwise, is there any way we can get a testcase that would be sensitive to these kinds of factors, to try and better detect if we have it right?

[Yue-Zhengyuan]

Indeed we have discussed before that storing the Trotter gates as a LocalOperator is a little bit weird. But the only difference is that exp(-dt * H) is the product of each gate (up to 1st order Trotter decomposition), instead of their sum. So I think the data structure that stores the Trotterized Hamiltonian will be essentially the same as the current LocalOperator. I just mentally shift to a different interpretation of each term when doing time evolution.

To check each evolution step actually corresponds dt, we can use either real time evolution (I use this in the full update PR) or finite temperature calculation (like Olivier's test).

[lkdvos]

I agree that they would be similar, except the product has an order: we are now explicitly imposing this order in the simple update code by taking first horizontal and then vertical edges, but in principle you could have different schemes, for example horizontal * vertical * diagonal * ..., which could be expressed by using a dictionary for the LocalOperator and a Vector or tuple for the gates. It's definitely not urgent, but it might be quite elegant if we ever want different trotterization schemes.

[lkdvos]

Good to go for me if tests pass

[Yue-Zhengyuan]

@ogauthe Are you satisfied with fixes in this PR?

[ogauthe]

So I just checked, the factor 2 is now correct and I get the same results for J2=0.0001 and for J2=0. If I find some other issue I will open another issue, I think this can be merged. Thank you for the fix @Yue-Zhengyuan !