for a MPS state. But I don’t know how to get the eigenvalues of a general MPO. My guess is something like svd(A::MPO, s). If A is my MPO density operator and sis the corresponding sites with total length N.
Hi good question. Getting the von Neumann entropy of a mixed state \rho which is in the form of an MPO is not an easy task, and how best to do it is somewhat of an open research question. But obtaining the Renyi entropy, especially the second Renyi entropy S_2 is something straightforwardly computable with MPOs.
For others reading this discussion, recall that S_2 = -\log(\text{Tr}[ \rho_A^2]) where \rho_A = \text{Tr}_B[\rho] is the partial trace of \rho taken over region B which is the complement of region A.
Here is then a sketch of \rho, \rho_A and \text{Tr}[ \rho_A^2] for the case where A is the first half of the sites of the MPO:
Then you just take the negative log of the result.
One thing I like about this tensor network viewpoint of the second Renyi entropy is that it makes it clear what people are talking about when they say it corresponds to a “double sheeted path integral” or some similar language to that. Really it just means the diagram above where the two density matrices are traced individually on region B, but then together on region A.
Hi Miles. Thanks for the answer. For other order Renyi entropy or von Neumann entropy. Is it possible to directly calculate it by the spectrum of a MPO (suppose it represents a density matrix). I guess this will be much slower but it may be acceptable if my system sizes are not very large. I tried to get the spectrum by svd(rho,s) with rho as my density operator and s as my siteinds. But it seems that this only works for an ITensor type object. Can I transfer a MPO to an ITensor type object?
The answer is that yes, you could get the von Neumann entropy of an MPO but you’d have to pay an exponential cost. If you want to turn an MPO into a single ITensor, all you need to do is just contract all of the MPO tensors together, like this:
# Say you have an MPO called M
T = ITensor(1.)
for j=1:length(M)
T *= M[j]
end
A shorter way is to use the prod function like this: T = prod(M).
Then you can just SVD apart the unprimed site indices from the primed site indices to get the spectrum of the operator, thought of as a map from the unprimed space to the primed space.
Say that the array of site indices of M is s, then you can do:
U,S,V = svd(T,s)
Whether there is a reliable way to get the same spectrum information in a way that doesn’t scale exponentially is more of an open question!
Hi Miles. A small question is that if I want to calculate the square of a MPO (my density operator) \rho^2, I can use apply(rho,rho) right? Also, is the quantity \text{tr}(\rho^2) equal to inner(rho',rho)? Thanks.
It’s a good question: actually inner(A,B) for two MPOs A and B will transpose (Hermitian conjugate) the first MPO A, thought of as a large matrix, before contracting it with B. So in general it is not the same as multiplying A times B, which is what apply(A,B) does, and then tracing the result AB.
But these would be the same for the case that A is Hermitian.
I would recommend the following for you: first use tr(apply(rho,rho)) because you can be sure that computes \text{Tr}[\rho^2]. Then verify that inner(rho,rho) gives the same result using the assumption that rho is Hermitian. Then you might want to use inner because it is faster.
Hello, if i want to calculate a renyi entropy of a MPS tensor, please help me check my code,
function entanglement(psi, b)
orthogonalize!(psi, b)
U,S,V = svd(psi[b], (linkind(psi, b-1), siteind(psi, b)) )
SVN = 0.0
for n=1: dim(S, 1)
p = S[n, n]^2
if p > 0.0
SVN = (SVN - log(p^2) )
end
end
SVN= SVN
return SVN
end
The second-order renyi entropy is S=-\text{Tr}[\log \rho^2].
Also, please be specific about the question you are asking (versus asking us just to read your code). What is the part of the code you are unsure about?
What have you done to test the code and if it didn’t work as you expect, what happened?