TDVP for a Mixed Site System

[chemicalfiend]

Hello!

I have a question about running TDVP for a mixed-site system of bosons and fermions. I have a system with nearest-neighbor coupled fermionic sites each coupled locally to bosonic sites (Hamiltonian generated roughly similar to here).

When I run the following block of code :

using ITensors, ITensorMPS

#include("boson.jl")

let
 
Nf = 5         # lattice sites
nboson = 12    # Max boson occupation
N = 2*Nf       # Total number of sites (lattice + boson)
    
s = Vector{Index}(undef, N)
for i in 1:N
    if isodd(i)
        s[i] = siteind("Fermion", i, conserve_qns=true)
    else
        #s[i] = siteind("MyBoson", i,dim=nboson, conserve_qns=true, conserve_number=false,)
        s[i] = Index([QN() => nboson],"Boson,Site,n=$i")
    end
end

mid = findall(i -> isodd(i), 1:N)[ceil(Int, count(isodd, 1:N)/2)]


# initial mps with only middle site populated (other fermionic sites and all bosonic sites are unoccupied)
    
psi0 = MPS(ComplexF64, s, i -> begin
    if isodd(i)
        i == mid ? "1" : "0"
    else
        "0"
    end
end)

# Nearest neighbor electronic site coupling (odd sites) with local bosonic modes (even sites)
# Hamiltonian : H = ω ∑ b†b + g ∑(b† + b) |α><α| - t ∑* |α><β| 

os = OpSum()
t = 1.0
ω = 1.0
g = 1.0

# boson energy (1)
for i in 2:2:N
    os += ω, "N", i
end

# coupling to bosonic modes (2)
for i in 1:2:N-1
    os += g, "N", i, "Adag", i+1
    os += g, "N", i, "A", i+1
end

# inter-site hopping (3)
for i in 1:2:N-2
    os += -t, "Cdag", i, "C", i+2
    os += -t, "Cdag", i+2, "C", i
end


H = MPO(os, s)

psi2 = tdvp(H, -1.0im, psi0; nsweeps=10)

    
println("*****INITIAL POPULATION FERMIONS*****")
println(expect(psi0,"N",sites=1:2:(2*Nf)))

println("*****FINAL POPULATION FERMIONS*****")
println(expect(psi2,"N",sites=1:2:(2*Nf)))

println("*****INITIAL POPULATION BOSONS*****")
println(expect(psi0,"N",sites=2:2:(2*Nf+1)))

println("*****FINAL POPULATION BOSONS*****")
println(expect(psi2,"N",sites=2:2:(2*Nf+1)))
    
end

I get no population change from the initial site whatsoever.

*****INITIAL POPULATION FERMIONS*****
[0.0, 0.0, 1.0, 0.0, 0.0]
*****FINAL POPULATION FERMIONS*****
[0.0, 0.0, 1.0, 0.0, 0.0]
*****INITIAL POPULATION BOSONS*****
[0.0, 0.0, 0.0, 0.0, 0.0]
*****FINAL POPULATION BOSONS*****
[0.0, 0.0, 0.9193953880143915, 0.0, 0.0]

I tried running the same commenting out MPO terms (2) and (3) (coupling to the bosons), and get the same result except for the final bosonic population being 0.0 as well.

However when I change code without mixed sites (just the fermions) :

using ITensors, ITensorMPS

let
N = 5         # Total number of sites 
s = Vector{Index}(undef, N)
for i in 1:N
    s[i] = siteind("Fermion", i, conserve_qns=true)
end

psi0 = MPS(ComplexF64, s, i -> begin  i == div(N,2)+1 ? "1" : "0"  end)


expect(psi0,"N")


os = OpSum()
t = 1.0


# hopping
for i in 1:(N-1)
    os += -t, "Cdag", i, "C", i+1
    os += -t, "Cdag", i+1, "C", i
end

H = MPO(os, s)
psi2 = tdvp(H, -1.0im, psi0; nsweeps=10)

println("*****INITIAL POPULATION FERMIONS*****")
println(expect(psi0, "N", sites=1:N))


println("*****FINAL POPULATION FERMIONS*****")
println(expect(psi2,"N", sites=1:N))

end

I get a change in the populations of the fermions :

*****INITIAL POPULATION FERMIONS*****
[0.0, 0.0, 1.0, 0.0, 0.0]
*****FINAL POPULATION FERMIONS*****
[0.1496491110241341, 0.32475046730718693, 0.05120755314873815, 0.3247385050912503, 0.14965436342869062]

I’m not really sure where I’m going wrong in the implementation of mixed site MPS / MPO. Especially since turning off the coupling to bosons is not making the two codes equivalent. I apologize for the long question, but I was hoping there was a quick fix to the mixed site implementation that would allow this to run.

Thanks in Advance!

For reference : I am running Julia Version 1.11.6 with ITensors v0.9.7

[miles]

My guess is this is most likely because you are doing a quench from an initial state which is a product or zero-entanglement state. TDVP by itself (not just the ITensor implementation but in general) can struggle to grow or adapt the bond dimension when it starts out initially too small. It could be working in the purely fermion-site case because our code defaults to 2-site TDVP which can achieve bond adaptation even in challenging cases, though not always perfectly.

(In contrast, starting from an entangled state with an expressive ‘basis’ like a ground state almost always works very well.)

A solution could be to try our expand function which performs the “global subspace expansion” method by Yang and White. You would need to apply this at some frequency during at least the initial period of time evolution. It’s a bit of an “art” about when to apply it and with what expansion parameters, so it takes some experimentation.

If applicable, another good solution can be to use TEBD for your initial few time steps until the bond dimensions grow to larger values then switch to TDVP for its greater efficiency. If your Hamiltonian is short-ranged enough, you might find TEBD sufficient for your whole problem.

[chemicalfiend]

Thanks for the suggestion Miles! Using expand / TEBD before running TDVP seems to work.

[miles]

Glad to hear it!