Heisenberg model with an external field

I am a student who has just started researching DMRG, and I don’t know much about the subject yet. I would appreciate the guidance and help from more experienced researchers. Currently, I am reading a paper related to DMRG with the DOI: 10.1103/PhysRevB.110.014413. It seems that this paper also uses ITensor for implementation, and I want to try to replicate its results to better understand how to use ITensor.

However, I encountered difficulties when trying to reproduce the DMRG-related figures, particularly Figure 3 from the paper. The Hamiltonian in the paper for (\theta = 0) is given by:

H_{\theta=0} = \sum_{i} \left(-\sigma_i^x \sigma_{i+1}^x - \sigma_i^y \sigma_{i+1}^y + \sigma_i^z \sigma_{i+1}^z\right) - (h_z + 2) \sum_{i} \sigma_i^z

Based on this Hamiltonian, I calculated

m_z = \frac{1}{N} \sum_i \langle S_i^z \rangle,

but the final result remains around 0.5, without large change. However, the results in the paper show that the power law continuously rises and ultimately stabilizes at 0.5.

using ITensors, ITensorMPS
N = 257
sites = siteinds("S=1/2",N)
os = OpSum()

for hz=0:0.1:2.5
    for j=1:N-1
      os -= 0.5,"S+",j,"S-",j+1
      os -= 0.5,"S-",j,"S+",j+1
      os += "Sz",j,"Sz",j+1
    end
    for j=1:N
        os -=(hz+2),"Sz",j
    end
    H = MPO(os,sites)
    nsweeps = 5 # number of sweeps is 5
    maxdim = [10,20,100,100,200,380] # gradually increase states kept
    cutoff = [1E-10] # desired truncation error
    psi0 = random_mps(sites;linkdims=2)
    energy,psi = dmrg(H,psi0;nsweeps,maxdim,cutoff)
    magz = expect(psi,"Sz")
    mz = sum(magz) / N
    println("$hz $mz")
end

Could anyone please help me identify any errors in my approach or calculations? Thank you!

Glad you are trying ITensor. It’s a little hard to know, even for experts, why the results you are obtaining are not matching what you saw in the paper. Usually this will be a combination of things like that your DMRG calculation may not be fully converged, or there are some finite-size or open-boundary effects, or other practical considerations like how observables are obtained and processed, which could all be playing a role.

I think the three main things you should investigate here are:

  1. do you get different results if you increase nsweeps ?
  2. do you get different results for different system sizes N?
  3. what does \langle S^z_j \rangle look like, visually? You should definitely make plots of this quantity and look at them to reason about if you are e.g. measuring boundary or bulk effects, and just get a general intuition about the behavior of your system.

We also provide a detailed FAQ page about checking if DMRG is converged and how to ensure convergence.

Thank you for your response. I will further explore these questions and conduct additional tests. I will also share the relevant results afterwards.