Gauge improvements

examples/vumps/vumps_hubbard_extended.jl

@@ -56,6 +56,7 @@ subspace_expansion_kwargs = (cutoff=cutoff, maxdim=maxdim)
 
 println("\nRun VUMPS on initial product state, unit cell size $N")
 ψ = vumps_subspace_expansion(H, ψ; outer_iters, subspace_expansion_kwargs, vumps_kwargs)
+ψ = orthogonalize(ψ.AL, :; tol=1e-14) # ensure translation invariance
 
 # Check translational invariance
 println("\nCheck translational invariance of optimized infinite MPS")

src/infinitecanonicalmps.jl

@@ -249,21 +249,37 @@ function InfMPS(s::CelledVector, f::Function)
   return ψ = InfiniteCanonicalMPS(ψL, ψC, ψR)
 end
 
-function finite_mps(ψ::InfiniteCanonicalMPS, range::AbstractRange)
+function finite_mps(
+  ψ::InfiniteCanonicalMPS,
+  range::AbstractRange;
+  ortho_lims::UnitRange=last(range):last(range),
+)
   @assert isone(step(range))
+  @assert first(ortho_lims) == last(ortho_lims) # TODO: variable ortho_lims
+  @assert first(ortho_lims) ∈ range
+
   N = length(range)
-  ψ_finite = ψ.AL[range]
-  ψ_finite[N] *= ψ.C[last(range)]
+  ψ_finite = Vector{ITensor}(undef, N)
+  for i in first(range):first(ortho_lims)
+    ψ_finite[i] = ψ.AL[i]
+  end
+  ψ_finite[first(ortho_lims)] *= ψ.C[first(ortho_lims)]
+
+  for i in (last(ortho_lims) + 1):last(range)
+    ψ_finite[i] = ψ.AR[i]
+  end
+
   l0 = linkind(ψ.AL, first(range) - 1 => first(range))
-  l̃0 = sim(l0)
   lN = linkind(ψ.AR, last(range) => last(range) + 1)
+  l̃0 = sim(l0)
   l̃N = sim(lN)
   δl0 = δ(dag(l̃0), l0)
   δlN = δ(dag(l̃N), lN)
   ψ_finite[1] *= δl0
   ψ_finite[N] *= dag(δlN)
-  ψ_finite = MPS([dag(δl0); [ψ_finiteᵢ for ψ_finiteᵢ in ψ_finite]; δlN])
-  set_ortho_lims!(ψ_finite, (N + 1):(N + 1))
+  ψ_finite = MPS(
+    [dag(δl0); ψ_finite; δlN]; ortho_lims=(first(ortho_lims) + 1):(last(ortho_lims) + 1)
+  )
   return ψ_finite
 end
 function ITensorMPS.expect(ψ::InfiniteCanonicalMPS, o::String, n::Int)

src/orthogonalize.jl

@@ -1,7 +1,13 @@
+using KrylovKit: schursolve, Arnoldi
 # TODO: call as `orthogonalize(ψ, -∞)`
 # TODO: could use commontags(ψ) as a default for left_tags
 function right_orthogonalize(
-  ψ::InfiniteMPS; left_tags=ts"Left", right_tags=ts"Right", tol::Real=1e-12
+  ψ::InfiniteMPS;
+  left_tags=ts"Left",
+  right_tags=ts"Right",
+  tol::Real=1e-12,
+  eager=true,
+  ishermitian_rtol::Real=tol * 100,
 )
   # A transfer matrix made from the 1st unit cell of the infinite MPS
   T = TransferMatrix(ψ)
@@ -12,9 +18,20 @@ function right_orthogonalize(
   # Start by getting the right eivenvector/eigenvalue of T
   # TODO: make a function `right_environments(::InfiniteMPS)` that computes
   # all of the right environments using `eigsolve` and shifting unit cells
-  λ⃗₁ᴿᴺ, v⃗₁ᴿᴺ, eigsolve_info = eigsolve(T, v₁ᴿᴺ, 1, :LM; tol=tol)
+
+  # original eigsolve function, switch to schur which enforces real
+  #λ⃗₁ᴿᴺ, v⃗₁ᴿᴺ, eigsolve_info = eigsolve(T, v₁ᴿᴺ, 1, :LM; tol, eager)
+  TT, v⃗₁ᴿᴺ, λ⃗₁ᴿᴺ, info = schursolve(T, v₁ᴿᴺ, 1, :LM, Arnoldi(; tol, eager))
   λ₁ᴿᴺ, v₁ᴿᴺ = λ⃗₁ᴿᴺ[1], v⃗₁ᴿᴺ[1]
 
+  if info.converged == 0
+    @warn "orthogonalize not converged after $(info.numiter) iterations"
+  end
+
+  if size(TT, 2) > 1 && TT[2, 1] != 0
+    @warn("Non-unique largest eigenvector of transfer matrix found")
+  end
+
   if imag(λ₁ᴿᴺ) / norm(λ₁ᴿᴺ) > 1e-15
     @show λ₁ᴿᴺ
     error(
@@ -24,11 +41,11 @@ function right_orthogonalize(
 
   # Fix the phase of the diagonal to make Hermitian
   v₁ᴿᴺ .*= conj(sign(v₁ᴿᴺ[1, 1]))
-  if !ishermitian(v₁ᴿᴺ; rtol=tol)
-    @show λ₁ᴿᴺ
-    @show v₁ᴿᴺ
+  if !ishermitian(v₁ᴿᴺ; rtol=ishermitian_rtol)
+    #@show λ₁ᴿᴺ
+    #@show v₁ᴿᴺ
     @show norm(v₁ᴿᴺ - swapinds(dag(v₁ᴿᴺ), reverse(Pair(inds(v₁ᴿᴺ)...))))
-    error("v₁ᴿᴺ not hermitian")
+    @warn("v₁ᴿᴺ is not hermitian within rtol=$ishermitian_rtol")
   end
   if norm(imag(v₁ᴿᴺ)) / norm(v₁ᴿᴺ) > 1e-13
     println(
@@ -105,11 +122,16 @@ end
 function mixed_canonical(
   ψ::InfiniteMPS; left_tags=ts"Left", right_tags=ts"Right", tol::Real=1e-12
 )
-  _, ψᴿ, _ = right_orthogonalize(ψ; left_tags=ts"", right_tags=ts"Right")
-  ψᴸ, C, λ = left_orthogonalize(ψᴿ; left_tags=ts"Left", right_tags=ts"Right")
+  _, ψᴿ, _ = right_orthogonalize(ψ; left_tags=ts"", right_tags)
+  ψᴸ, C, λ = left_orthogonalize(ψᴿ; left_tags, right_tags)
   if λ ≉ one(λ)
     error("λ should be approximately 1 after orthogonalization, instead it is $λ")
   end
+  # tags have been added, remove them now
+  # This is broken until we fix the tags for C
+  #removetags!(ψᴸ, left_tags)
+  #removetags!(ψᴿ, right_tags)
+  #removetags!(C,addtags(right_tags,left_tags))
   return InfiniteCanonicalMPS(ψᴸ, C, ψᴿ)
 end