[ITensors] [BUG] QN error when applying MPO to MPS

[emstoudenmire]

Description of bug

When applying certain MPOs to MPS, there is an error coming from eigen from within the apply(::MPO,::MPS) function. It seems to only happen when:

1. the MPO has non-trivial flux
2. the MPS has a bond dimension of 1

Looking into it some more, it could actually be a problem with randomMPS (and a similar problem occurring with productMPS), which is that the link indices seem to have the wrong quantum numbers in the case of making a product state. See below for the product state I made.

Minimal code demonstrating the bug or unexpected behavior

using ITensors

let
  N = 4
  s = siteinds("Boson",N; conserve_qns=true)

  terms = OpSum()
  terms += "Adag",3
  W = MPO(terms,s)

  psi = randomMPS(s,[isodd(n) ? "1" : "0" for n=1:N]; linkdims=1)
  @show psi

  Wpsi = apply(W,psi)
  @show norm(Wpsi)

  return
end

Actual output or behavior

Output of minimal runnable code

ERROR: LoadError: AssertionError: flux(A) == QN()
Stacktrace:
  [1] macro expansion
    @ ~/.julia/dev/ITensors/src/decomp.jl:229 [inlined]
  [2] macro expansion
    @ ~/.julia/dev/ITensors/src/global_variables.jl:177 [inlined]
  [3] eigen(A::ITensor, Linds::Vector{Index{Vector{Pair{QN, Int64}}}}, Rinds::Vector{Index{Vector{Pair{QN, Int64}}}}; kwargs::Base.Pairs{Symbol, Any, Tuple{Symbol, Symbol}, NamedTuple{(:ishermitian, :tags), Tuple{Bool, TagSet}}})
    @ ITensors ~/.julia/dev/ITensors/src/decomp.jl:227
  [4] contract(::ITensors.Algorithm{:densitymatrix}, A::MPO, ψ::MPS; kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ ITensors ~/.julia/dev/ITensors/src/mps/mpo.jl:631
  [5] contract
    @ ~/.julia/dev/ITensors/src/mps/mpo.jl:587 [inlined]
  [6] contract(A::MPO, ψ::MPS; alg::String, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ ITensors ~/.julia/dev/ITensors/src/mps/mpo.jl:531
  [7] contract
    @ ~/.julia/dev/ITensors/src/mps/mpo.jl:516 [inlined]
  [8] #apply#922
    @ ~/.julia/dev/ITensors/src/mps/mpo.jl:507 [inlined]
  [9] product(A::MPO, ψ::MPS)
    @ ITensors ~/.julia/dev/ITensors/src/mps/mpo.jl:507
 [10] top-level scope
    @ ~/software/troubleshooting/miles/ampo_redesign/apply_bug.jl:14
 [11] include(fname::String)
    @ Base.MainInclude ./client.jl:451
 [12] top-level scope
    @ REPL[4]:1
in expression starting at /Users/mstoudenmire/software/troubleshooting/miles/ampo_redesign/apply_bug.jl:3

Version information

• Output from versioninfo():

julia> versioninfo()
Julia Version 1.7.2
Commit bf53498635 (2022-02-06 15:21 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin19.5.0)
  CPU: Apple M1 Max
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-12.0.1 (ORCJIT, westmere)
Environment:
  JULIA_DIR = /Applications/Julia-1.7.app/Contents/Resources/julia

Using latest development branch of ITensors.jl

[mtfishman]

Interesting, looks like somehow the density matrix ends up with a nonzero flux. We should probably actually print flux(A) in that error message.

[emstoudenmire]

Yes and I think the bug could really be originating back with one of the MPS constructors

[mtfishman]

@emstoudenmire I'm trying to clear out some issues. Is this still an issue?

[emstoudenmire]

Appears to still be an issue, though the error message has changed. I checked that with the current version of everything that the inputs to apply seem reasonable (i.e. it's not failing before that). Both the W MPO in the example and the psi MPS look to be valid and well-formed, unless there is something subtly wrong with one of them.

It's something I could look into as well. Anyway yes let's keep this open.