Entanglement entropy calculation

[Faridun]

Dear all,

I have calculated entanglement entropy of a fermionic system with onsite interaction and magnetic and potential fields. I have plotted the results of calculation EA (entanglement entropy in section A) as a function of LA (LA is a size of subsystem).

Why the plotted curve is not smooth? Is this might be a result of non-converged calculations to the ground state?

Here is example of code:

using ITensors

energy_array = Float64[]
energy1_array = Float64[]
energy2_array = Float64[]

energy_difference_array = Float64[]

psi_array = Float64[]
psi_array1 = Float64[]
es = Float64[]
for iteration in 1:500
    N = 96
    Npart = 96
    t1 = 1.0
    b=22
    U = 1
    h = 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
    # Magnetic field term
    for i = 1:N
        hz = 20 * rand() - 10  # [-1, 1]
        println(hz)
        ampo += hz, "Nup", i
        ampo -= hz, "Ndn", i
    end

    # Potential field term
    for i = 1:N
        hE = rand() - 0.5  # [-0.5, 0.5]
        println(hE)
        ampo += hE, "Ntot", i

    end


    H = MPO(ampo, sites)
    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)
    psi1_init = randomMPS(sites, state)
    energy1, psi1 = dmrg(H, [psi0], psi1_init, sweeps)

    # Append energy values to arrays
    push!(energy_array, energy)
    push!(energy_array, energy1)
    orthogonalize!(psi, b)
    U,S,V = svd(psi[b], (linkind(psi, b-1), siteind(psi,b)))
    SvN = 0.0
    for n=1:dim(S, 1)
      p = S[n,n]^2
      SvN -= p * log(p)
    end
    push!(es, SvN)
    orthogonalize!(psi1, b)
    U,S,V = svd(psi1[b], (linkind(psi1, b-1), siteind(psi1,b)))
    SvN = 0.0
    for n=1:dim(S, 1)
      p = S[n,n]^2
      SvN -= p * log(p)
    end
    push!(es, SvN)

end

println("Energy Array: ", energy_array)
println(es)

Thank you in advance.
Faridun

[joshperry3110]

Hi,
I am not sure about why your plotted curves are not smooth, but if you want to check if the DMRG calculation has converged, you may benefit from the Observer system, especially as I believe the entanglement entropy is one of the last properties to converge.
It’s possible to create an Observer that checks for entropy convergence after each sweep of DMRG, up to a specified convergence tolerance, so that you can be more confident that you have indeed reached the correct value.
The documentation for creating custom observers in Julia can be found here: Observer System for DMRG · ITensors.jl
I would suggest using the checkdone! method to reach the desired accuracy.

[miles]

It looks like part of your Hamiltonian has random fields in it (the part where you call rand()). Is that correct? So then is your for loop doing disorder averaging? It may be that you need a very large number of disorder configurations to get a smooth average – try also measuring statistics such as the variance to see how much the value is fluctuating. Before you do that, though, I would try studying various fixed disorder configurations to make sure things are converged in the number of sweeps and DMRG parameters etc. Are you sure 6 sweeps is enough? How did you check?