Handle diagonal C - ensure TensorKit v0.14 compatibility (#223)

* ensure C is not diagonal

* improve entanglement spectrum/entropy

Author: lkdvos

src/algorithms/changebonds/svdcut.jl

@@ -91,7 +91,8 @@ function changebonds(ψ::InfiniteMPS, alg::SvdCut)
     ψ = if space(ncr, 1) != space(copied[1], 1)
         InfiniteMPS(copied)
     else
-        InfiniteMPS(copied, complex(ncr))
+        C₀ = TensorMap(complex(ncr))
+        InfiniteMPS(copied, C₀)
     end
     return normalize!(ψ)
 end

src/algorithms/toolbox.jl

@@ -5,10 +5,20 @@ Calculate the Von Neumann entanglement entropy of a given MPS. If an integer `si
 given, the entropy is across the entanglement cut to the right of site `site`. Otherwise, a
 vector of entropies is returned, one for each site.
 """
-entropy(state::InfiniteMPS) = map(c -> -tr(safe_xlogx(c * c')), state.C);
+entropy(state::InfiniteMPS) = map(Base.Fix1(entropy, state), 1:length(state))
 function entropy(state::Union{FiniteMPS,WindowMPS,InfiniteMPS}, loc::Int)
-    return -tr(safe_xlogx(state.C[loc] * state.C[loc]'))
-end;
+    S = zero(real(scalartype(state)))
+    tol = eps(typeof(S))
+    for (c, b) in entanglement_spectrum(state, loc)
+        s = zero(S)
+        for x in b
+            x < tol && break
+            s += (x * log(x))^2
+        end
+        S += oftype(S, dim(c) * s)
+    end
+    return S
+end
 
 # function infinite_temperature(H::MPOHamiltonian)
 #     return [permute(isomorphism(storagetype(H[1, 1, 1]), oneunit(sp) * sp,
@@ -84,9 +94,8 @@ matrix to the right of a given site. This is a dictionary mapping the charge to
 values.
 """
 function entanglement_spectrum(st::Union{InfiniteMPS,FiniteMPS,WindowMPS}, site::Int=0)
-    @assert site <= length(st)
-    _, S, = tsvd(st.C[site])
-    return TensorKit.SectorDict(c => real(diag(b)) for (c, b) in blocks(S))
+    checkbounds(st, site)
+    return LinearAlgebra.svdvals(st.C[site])
 end
 
 """

src/utility/utility.jl

@@ -148,11 +148,6 @@ function randomize!(a::AbstractBlockTensorMap)
     return a
 end
 
-function safe_xlogx(t::AbstractTensorMap, eps=eps(real(scalartype(t))))
-    (U, S, V) = tsvd(t; alg=SVD(), trunc=truncbelow(eps))
-    return U * S * log(S) * V
-end
-
 """
     tensorexpr(name::Symbol, ind_out, [ind_in])