OutOfMemoryError() - Githubissues

Some code that worked perfectly fine 24 hours ago is now giving me this weird OutOfMemoryError. As far as I can tell, the only thing that has changed is that I downloaded the updated version of the code yesterday morning on 04/23/20 (following pull request #319), whereas the code I was using before was the version from 04/16/20 (following pull request #287).

Here is the relevant part of my code. It is attempting to run DMRG on the spin-1/2 lattice formulation of the Schwinger model. Up until a day ago, my code worked fine even up to N = 900. Yes, it ran very slowly, because this Hamiltonian has non-local interactions leading to a quadratic number of terms in the MPO. However, it still always worked, with no hint of anything resembling an OutOfMemoryError.

N = 200
x = 100
sites = siteinds("S=1/2",N, conserve_qns=true) # need to stay in total z-spin = 0 sector

# Input operator terms which define a Hamiltonian matrix, and convert these terms to an MPO tensor network
# Refer to paper to see exactly how these terms correspond to Schwinger model Hamiltonian
ampo = AutoMPO()
trace = N*μ/2 + N^2/8 # constant multiple of identity subtracted from Hamiltonian to make it traceless
add!(ampo,2*2 * trace,"Sz*Sz",1) # cheat designed to add trace to Hamiltonian
for j = 1:N
  if j != N
    add!(ampo,x,"S+",j,"S-",j+1) # S+ and S- are the same as σ+ and σ-
    add!(ampo,x,"S-",j,"S+",j+1)
    # S operators are 1/2 of Pauli matrices, so I have put factors of 2 to compensate
    add!(ampo,2 * trunc((N-j+1)/2)/2,"Sz",j)
    for i=1:j-1
      add!(ampo,2*2 * (N-j)/2,"Sz",i,"Sz",j)
    end
  end
  add!(ampo,2 * μ*(-1)^(j-1)/2,"Sz",j)
end
H = MPO(ampo,sites)

sweeps = Sweeps(10)
# make max number dimensions larger with larger chain length to maintain accuracy
maxdim!(sweeps, 10,20,100,100,200)
cutoff!(sweeps, 1E-10)
init_state = [isodd(n) ? "Up" : "Dn" for n=1:N]
psi_init = productMPS(sites, init_state)

# Run the DMRG algorithm, returning energy
# (dominant eigenvalue) and optimized MPS
ground_energy, psi_ground = dmrg(H, psi_init, sweeps; outputlevel=0)
ground_energy_density = ground_energy / (2*x*N)
println("Ground state energy = $ground_energy")
println("Ground state energy density = $ground_energy_density")

However, now, when I attempt to run the above code, I get the following error:

ERROR: LoadError: OutOfMemoryError()
Stacktrace:
 [1] Array at .\boot.jl:406 [inlined]
 [2] Array at .\boot.jl:414 [inlined]
 [3] similar at .\array.jl:334 [inlined]
 [4] similar at D:\buildbot\worker\package_win64\build\usr\share\julia\stdlib\v1.3\LinearAlgebra\src\adjtrans.jl:197 [inlined]
 [5] similar at .\abstractarray.jl:627 [inlined]
 [6] _unsafe_getindex(::IndexCartesian, ::LinearAlgebra.Adjoint{Float64,Array{Float64,2}}, ::UnitRange{Int64}, ::UnitRange{Int64}) at .\multidimensional.jl:682
 [7] _getindex at .\multidimensional.jl:670 [inlined]
 [8] getindex(::LinearAlgebra.Adjoint{Float64,Array{Float64,2}}, ::UnitRange{Int64}, ::UnitRange{Int64}) at .\abstractarray.jl:981
 [9] #qn_svdMPO#333(::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(ITensors.qn_svdMPO), ::AutoMPO, ::Array{Index{Array{Pair{QN,Int64},1}},1}) at C:\Users\user\.julia\packages\ITensors\JA9pG\src\physics\autompo.jl:614
 [10] qn_svdMPO at C:\Users\user\.julia\packages\ITensors\JA9pG\src\physics\autompo.jl:531 [inlined]
 [11] #MPO#348 at C:\Users\user\.julia\packages\ITensors\JA9pG\src\physics\autompo.jl:794 [inlined]
 [12] MPO(::AutoMPO, ::Array{Index{Array{Pair{QN,Int64},1}},1}) at C:\Users\user\.julia\packages\ITensors\JA9pG\src\physics\autompo.jl:789
 [13] top-level scope at D:\Julia Files\Schwinger Model\Schwinger_DMRG_paper_replication.jl:69        
in expression starting at D:\Julia Files\Schwinger Model\Schwinger_DMRG_paper_replication.jl:16

The line numbers in [12] and [13] won't match, but the important part is that the error occurs at the line "H = MPO(ampo,sites)." Because of this, it can't be due to DMRG. There could be some issue on my desktop, but since this worked fine a day ago on my desktop, I am inclined to believe that it is an issue with some update between 04/16/20 and 04/23/20. Any help with this would be greatly appreciated.

Hi Sujay, Writing you back here since I think you said you’re not getting notifications about ITensor Support Forum posts. I took a bit more time to think about your specific Hamiltonian and I’m fairly convinced about two things:

this is a tough case for AutoMPO, the way it currently works, so I’m not surprised it eats up a lot of time and memory
but the good news is: I’m fairly sure your Hamiltonian has a very compact MPO representation

As far as I can see from a brief reading of the code, the only term which is troublesome and non-local is the one where the site values are named “i” and “j” and where i < j but can otherwise take any value between 1 and j-1.

For a term like that, there ought to be a simple way to construct an MPO, by using the “finite state automaton” picture of how MPOs work. (If you don’t know this picture, I could send you some references.) In your case, the automaton would have the following internal states (this is just for making that one type of term, not the others):

we haven’t placed any operators yet
we have placed the “i” operator now, but not the “j” one yet
we have now placed “j”, and are done

For the off-diagonal elements in the MPO that do the transitions, the one to take the most care with would be the elements that do the transitions from state 2 to state 3, as the coefficients go like (N-j). But since the MPO has “knowledge” of what site its on (in the sense that each MPO tensor can be different for a finite system), you can just set the coefficient to be proportional to (N-j) by hand.

So I’m pretty convinced you could make a nice, custom MPO for your case if you research it a bit. You could also still use AutoMPO to make an MPO for the rest of your local terms, and use ITensor’s feature which lets you pass multiple MPOs to the dmrg function and which “lazily” treats them as summed together in an efficient way (without actually ever summing the MPOs together).

Let me know if you have any questions -

Miles

Hello Miles,

Thank you for letting me know! It would be great if you could send me some of those resources. I would definitely be interested in trying to figure out how to find a compact MPO representation for the Hamiltonian that I am working on.

By the way, a few days ago, I used the old version of the repository (following pull request #287 ) for my code, and I can confirm that it worked nicely. So for now I will stick with the old version of the code (until I can figure out how to implement your recommendation), but that may help localize what change made my code in particular stop working. I'll keep looking into that if I have time and let you know if I find anything more specific. I'll also try some of your simpler recommendations (like passing multiple MPOs to the dmrg function) one of these days.

Oh yeah, and now that I know I'm not getting notifications, I'm planning to at least do a quick check on the ITensor Support Q&A page occasionally so that I don't miss any replies for an extended period of time. But thank you for taking the time to put your reply both here and there!

Thanks, Sujay

Hi Sujay, Thanks for letting us know about the behavior of the older version. Not much has changed in the AutoMPO code since we added it to the Julia version, so if its behavior changed a lot in terms of memory usage, it could be something subtle such as how something internal to the ITensor type itself works.

Regarding resources for designing MPOs manually, the most detailed discussion of the "automaton" picture of MPS and MPO's is, I believe, in this article: https://arxiv.org/abs/0708.1221

Another nice article that briefly discusses some examples of MPOs and how they work is this one: https://arxiv.org/abs/0804.2509

Finally, I have some talk slides that go over such constructions very briefly, but the visual aids may help you: http://tensornetwork.org/reviews_resources.html (See the "Tensor Networks and Applications" talk slides by me which are linked there.) Specific slide deck discussing MPOs: https://itensor.org/miles/BrazilLectures/TNAndApplications03.pdf

ITensor / ITensors.jl

OutOfMemoryError() #322