# TEBD for 2-leg Hubbard Model

**URL:** <https://itensor.discourse.group/t/tebd-for-2-leg-hubbard-model/139>\
**Category:** DMRG and Numerical Methods\
**Tags:** julia, excited\_states\
**Created:** [May 10, 2022, 9:21pm UTC](https://itensor.discourse.group/t/tebd-for-2-leg-hubbard-model/139 "2022-05-10T21:21:14Z")\
**Posts on this page:** 3\
**Page:** 1

<div class="post-metadata">

**Author:** ![nticea](https://avatars.discourse-cdn.com/v4/letter/n/eb8c5e/32.png) [@nticea](https://itensor.discourse.group/u/nticea)\
**Post date:** [May 10, 2022, 9:21pm UTC](https://itensor.discourse.group/t/tebd-for-2-leg-hubbard-model/139/1 "2022-05-10T21:21:15Z")

</div>

Hi there! I am trying to implement time evolving block decimation for the Hubbard model. I believe my DMRG works just fine (the state convergences, my equilibrium correlations look good, and the ground state energy matches up with benchmark values), but there is a bug in how I am constructing the Trotter gates for my TEBD evolution.

Specifically: when I calculate the inner product of the time-evolved ground state with its original self, the overlap decreases with time. I encountered this previously in the case of the 1D TEBD, and ultimately realized that it was due to a mistake in how I was constructing my Trotter gates. The only difference between the 1D case and the 2-leg ladder, as far as I am aware, is the insertion of fSWAP gates before and after the application of a two-site operator between non-adjacent sites. Here is my code for constructing the Trotter gates:

```julia
  gates = ITensor[]
  for b in lattice
      s1_idx = min(b.s1, b.s2)
      s2_idx = max(b.s1, b.s2)

      # Check if we need to apply any swaps (sites are not adjacent) 
      if abs(s2_idx-s1_idx) > 1
          # Apply SWAP to get to new configuration
          SWAP_forward = []
          s_last = s1_idx
          for s_next in s1_idx+1:(s2_idx-1)
              push!(SWAP_forward, op("fSWAP", sites[s_last], sites[s_next]))
              s_last = s_next
          end
          append!(gates, SWAP_forward)

          # Apply gate
          hj_twosite = make_twosite(sites[s_last], sites[s2_idx])
          push!(gates,hj_twosite)

          # Apply SWAP to get back 
          SWAP_backward = []
          for s_next in s_last-1:s1_idx
              push!(SWAP_backward, op("fSWAP", sites[s_last], sites[s_next]))
              s_last = s_next
          end
          append!(gates, SWAP_backward)
          #append!(gates, reverse(SWAP_forward))

      else # no need to apply swaps 
          hj_twosite = make_twosite(sites[s1_idx], sites[s2_idx])
          push!(gates,hj_twosite)
      end

  end

  # One-site terms (no swap needed)
  for j=1:N-1
      hj_onesite = make_onesite(sites[j])     
      push!(gates,hj_onesite)
  end

  # End site (no swap needed)
  hn = make_endsite(sites[N])
  push!(gates,hn)

  # Add all the gates in reverse 
  append!(gates,reverse(gates))

```

I am fairly certain that my `make_twosite(...)`, `make_onesite(...)`, and `make_endsite(...)` gates are correct, because when I run this code with only one leg in the ladder, I am getting good overlap between the initial ground state and the time-evolved ground state. I’ve also attached a screenshot of my fSWAP operator in (what I believe to be) the correct basis.

 ![Screen Shot 2022-05-10 at 2.19.49 PM](https://global.discourse-cdn.com/free1/uploads/itensor/original/1X/787cd9712a1b5da04a79a9f4db25d6feb0bbaef1.png)

Thank you so much for all the continued support! Hope you have a great week.

---

<div class="post-metadata">

**Author:** ![mtfishman](https://yyz2.discourse-cdn.com/free1/user_avatar/itensor.discourse.group/mtfishman/32/11_2.png) [@mtfishman](https://itensor.discourse.group/u/mtfishman)\
**Post date:** [May 10, 2022, 9:24pm UTC](https://itensor.discourse.group/t/tebd-for-2-leg-hubbard-model/139/2 "2022-05-10T21:24:39Z")

</div>

You should not need to implement your own (fermionic) swap gates, the `apply` function inserts them automatically. You can just input nonlocal gates into `apply` and it should just work.

---

<div class="post-metadata">

**Author:** ![nticea](https://avatars.discourse-cdn.com/v4/letter/n/eb8c5e/32.png) [@nticea](https://itensor.discourse.group/u/nticea)\
**Post date:** [May 11, 2022, 5:20am UTC](https://itensor.discourse.group/t/tebd-for-2-leg-hubbard-model/139/3 "2022-05-11T05:20:28Z")

</div>

Hi Matt,

Thanks so much for the help! Your solution ended up fixing my problem. When I tried adding phonons (I have my own Hubbard Holstein site type) I got an out-of-bounds indexing error that I traced back to the `_fermionic_swap(...)` function. I think this is due to the fact that there is a bug in that function at the part where we fill in the diagonals of the blocks of the fSWAP operator. Here is what I did to fix the problem:

```julia
function _fermionic_swap(s1::Index, s2::Index)
  T = ITensor(QN(), s1', s2', dag(s1), dag(s2))
  for b in nzblocks(T) 
    dval = 1.0 
    # Must be a diagonal block
    ((b[1] ≠ b[3]) || (b[2] ≠ b[4])) && continue
    n1, n2 = b[1], b[2]
    if _isodd_fermionic_parity(s1, n1) && _isodd_fermionic_parity(s2, n2)
      dval = -1.0 
    end
    Tb = ITensors.blockview(ITensors.tensor(T), b) 
    mat_dim = prod(dims(Tb)[1:2])  
    Tbr = reshape(Tb, mat_dim, mat_dim)
    for i in diagind(Tbr)
      # NDTensors.setdiagindex!(Tbr, dval, i) <-- OLD CODE
      setindex!(Tbr, dval, i) # <-- MY MODIFICATION  
    end
  end
  return T
end 

```

Not sure if this is actually a bug, but changing that line avoids the out-of-bounds error and (seemingly) results in the correct time evolution, as measured by the fact that the overlap between |\psi\_{GS}\rangle and |\psi\_{GS}(t)\rangle remains close to 1.

Thanks again for all the help!