Excited states of hubbard model

Dear all,
I am trying to calculate excited states of spinful fermi hubbard model as follows

using ITensors


let
  N = 20
  Npart = 10
  t1 = 1.0
  U  = 1.0

  sites = siteinds("Electron", N;
                   conserve_qns = true)

  ampo = AutoMPO()
  for f=1:N-1
    ampo += -t1, "Cdagup", f,    "Cup", f+1
    ampo += -t1, "Cdagup", f+1,  "Cup", f
    ampo += -t1, "Cdagdn", f,    "Cdn", f+1
    ampo += -t1, "Cdagdn", f+1,  "Cdn", f
  end

  for i=1:N
    ampo += U, "Nupdn", i
  end
  H = MPO(ampo,sites)

  #sweeps = Sweeps(6)
  #nsweeps = 30
  #maxdim!(sweeps,50,100,200,400,800,800)
  #cutoff!(sweeps,1E-12)
  nsweeps = 30
  maxdim = [10,10,10,20,20,40,80,100,200,200]
  cutoff = [1E-8]
  noise = [1E-6]
  #@show sweeps

  #state = ["Emp" for n=1:N]
  state = InitState(sites);
  p = Npart
  for i=N:-1:1
    if p > i
      println("Doubly occupying site $i")
      state[i] = "UpDn"
      p -= 2
    elseif p > 0
      println("Singly occupying site $i")
      state[i] = (isodd(i) ? "Up" : "Dn")
      p -= 1
    end
  end
  # Initialize wavefunction to be bond 
  # dimension 10 random MPS with number
  # of particles the same as `state`
  psi0 = MPS(state)

  # Check total number of particles:
  #@show flux(psi0)

  # Start DMRG calculation:
  energy, psi = dmrg(H, psi0, nsweeps, maxdim, cutoff, noise)
  #psi1_init = randomMPS(sites, state, 10)
  #energy1,psi1 = dmrg(H,[psi0],psi1_init; nsweeps)

  println("\nGround State Energy = $energy")
  #println("DMRG energy gap = ",energy1-energy);
end

What I should to use instead InitState in Julia? According to Exited states of Hubbard model - ITensor Support Q&A randomMPS does not work in a quantum number (QN) conserving calculation. Is any ideas how to manage this problem?

You can use this Tutorial about setting up a quantum-number conserving DMRG calculation in Julia ITensor:
https://itensor.github.io/ITensors.jl/dev/tutorials/QN_DMRG.html

There the initial state is chosen to be a product state, which you make by calling psi = MPS(sites,state).

Also, fortunately randomMPS does in fact work for QN conserving calculations, you are just required to provide an array like state telling the randomMPS function what product state to “base” the randomMPS from (i.e. to use as an initializer which is then scrambled while staying in the same QN sector as the product state provided).
Please see this documentation on randomMPS.

As I understood, I should use randomMPS:

####previous code
  state = [isodd(n) ? "Up" : "Dn" for n=1:N]
  psi0 = randomMPS(sites,state)
  energy, psi = dmrg(H, psi0, nsweeps, maxdim, cutoff, noise)

because using only MPS(sites,states) it is initial state and I got MethodError: no method matching MPS.
But in case of randomMPS(sites,state) I got another error MethodError: no method matching dmrg(::MPO, ::MPS, ::Int64, ::Vector{Int64}, ::Vector{Float64}, ::Vector{Float64})

Yes, please review the code example above closely. The reason for your error is not related to randomMPS but it is that nsweeps, maxdim, and onward are named keyword arguments so you will need to use a semicolon to separate them from the first (positional) arguments. This is a standard Julia syntax and not specific to ITensor.

Thank you for your reply.
But separating them by using semicolon changed nothing. It works only in case:

  state = [isodd(n) ? "Up" : "Dn" for n=1:N]
  psi0 = randomMPS(sites,state)
  @show flux(psi0)
  sweeps = Sweeps(6)
  maxdim!(sweeps,50,100,200,400,800,800)
  cutoff!(sweeps,1E-12)
  @show sweeps
  energy, psi = dmrg(H, psi0, sweeps)

But I need to calculate excited states also:

  psi1_init = randomMPS(sites,state)
  energy1,psi1 = dmrg(H,[psi0],psi1_init, sweeps)
  println(energy, energy1)

In this case I am getting the same values of energy and energy1.
Is this correct? How can I calculate the correct first excited state in this case?

The comment about adding the semicolon was to remove the “MethodError: no method matching dmrg(::MPO, ::MPS, ::Int64, ::Vector{Int64}, ::Vector{Float64}, ::Vector{Float64})”. Did that error go away?

Regarding excited states, yes the excited-state DMRG is not guaranteed to work without adjusting and experimenting with some parameters. The most important one is the “weight” parameter which adjusts the amount of penalty for the new state to be orthogonal to the previous ones provided. For more information about the weight parameter, please see this code example section of the documentation:
https://itensor.github.io/ITensors.jl/dev/examples/DMRG.html#Compute-excited-states-with-DMRG
Also, after your DMRG run, I would suggest putting a line like

@show inner(psi1, psi0)

to check the overlap of the two states (is it close to zero or to one?).