question about asymmetric hubbard model

[luckyzy]

Hello everyone, I have been working on my two-dimensional asymmetric Hubbard model and have noticed that the computation speed is very slow, and it consistently fails to converge. I would like to ask if this issue might be caused by my program. Below are my parameters:

let
  Ny = 6
  Nx = 6
  U=2
  tupup=1
  tdndn=0.2
  N = Nx * Ny

  sites = siteinds("Electron",N; conserve_sz=true)
  lattice = square_lattice(Nx, Ny; yperiodic=false)

  os = OpSum()
  for b in lattice
    os += tupup, "Cdagup", b.s1, "Cup", b.s2
    os += tupup, "Cup", b.s1, "Cdagup", b.s2
    os += tdndn, "Cdagdn", b.s1, "Cdn", b.s2
    os += tdndn, "Cdn", b.s1, "Cdagdn", b.s2
    os += U, "Nup", b.s1, "Ndn", b.s1
  end

  H = MPO(os, sites)
  state = [isodd(n) ? "Up" : "Dn" for n in 1:N]
  psi0 = random_mps(sites, state; linkdims=30)
  
  nsweeps = 20
  maxdim = [ 400, 800, 1600, 1800 ,2000 ,2200]
  cutoff = [1E-8]
  noise = [1E-5 ,1E-6, 1E-7, 1E-8, 0.0]

  energy, psi = dmrg(H, psi0; nsweeps, maxdim, cutoff,  noise)
  H2 = inner(H,psi,H,psi)
  E = inner(psi',H,psi)
  var = H2-E^2
  @show var
  return nothing

H = t_{\uparrow} \sum_{\langle i,j \rangle} \left( C_{i\uparrow}^{\dagger} C_{j\uparrow} + C_{j\uparrow}^{\dagger} C_{i\uparrow} \right)+t_{\downarrow} \sum_{\langle i,j \rangle} \left( C_{i\downarrow}^{\dagger} C_{j\downarrow} + C_{j\downarrow}^{\dagger} C_{i\downarrow} \right) + U \sum_{i} n_{i\uparrow} n_{i\downarrow}

Could everyone please help me verify if my Hamiltonian for the two-dimensional DMRG is correct? Feel free to adjust any part of it!

[miles]

My feedback would mainly be to say that Ny = 6 is rather ambitious for the Hubbard model with DMRG. State-of-the-art DMRG calculations for Ny = 4 cylinders needed to keep up to 20,000 states (maxdim=20,000) to converge and even then there were some questions about the accuracy. See this paper:
https://www.science.org/doi/full/10.1126/science.aal5304

(While that paper is about about a t-t'-U model, I would think that the general point would apply to your case too.)

In general, for the Hubbard model one has to approach wider cylinders very gradually and carefully with DMRG. This includes studying smaller cylinders first, trying various initial states and plotting the results to visualize the behavior e.g. are holes or particles ‘sticking’ to the boundaries or phase separating, and also considering if there are different / modified Hamiltonians which might show similar qualitative behavior but have less entanglement.

[luckyzy]

Thanks, Miles! I think the way you mentioned works really well.

[corbett5]

I’ve been using ITensors to study the Ferm-Hubbard Hamiltonian as well and can confirm that 6x6 is ambitious. The one other remark I have is that I think your interaction term is wrong.

You have (ignoring the other terms)

for b in lattice
  os += U, "Nup", b.s1, "Ndn", b.s1
end

but here for b in lattice is iterating over all the lattice bonds not sites. Even worse the number of times b.s1 will appear is not constant from site to site so you’ll wind up with different interaction strengths.

julia> lattice = square_lattice(2, 2; yperiodic=false)
4-element Vector{LatticeBond}:
 LatticeBond(1, 3, 1.0, 1.0, 2.0, 1.0, "")
 LatticeBond(1, 2, 1.0, 1.0, 1.0, 2.0, "")
 LatticeBond(2, 4, 1.0, 2.0, 2.0, 2.0, "")
 LatticeBond(3, 4, 2.0, 1.0, 2.0, 2.0, "")

julia> for b in lattice
       @show b.s1, b.s2
       end
(b.s1, b.s2) = (1, 3)
(b.s1, b.s2) = (1, 2)
(b.s1, b.s2) = (2, 4)
(b.s1, b.s2) = (3, 4)

I would imagine there’s a way to use the lattice system to do this but I would just do

for i in 1:N
  os += U, "Nup", i, "Ndn", i
end

Also, about the computational speed. Currently you’re only conserving the total spin, do you also want to conserve the total number of electrons (I think this is the most common case). Have you tried enabling the block sparse multithreading? I would try increasing the number of sweeps at lower bond dimensions substantially as a way to speed convergence at the larger bond dimensions. I would also recommend saving the MPS to disk so that you can resume the calculation and not have to redo it every time (here is an example of how to do this).