Hi,
I would like to compute the operator entanglement of an MPO, defined as S_{op} = - \sum_{\alpha} \lambda_\alpha^2 \log_2{\lambda_\alpha^2} where \lambda_\alpha are the Schmidt values across a cut of the system (see this paper). This is the same definition as the entanglement entropy for MPS states but extended to MPOs. Unlike the MPS case, the operator entanglement is not a genuine measure of entanglement. For example, the operator entanglement across bond b (with N_A on the left and N_B sites on the right) of the maximally mixed state \rho_{mm} = \frac{1}{2^N} \mathbb{I}^{\otimes N} is
since the reduced density matrix on \rho_{N_A} is
where | \alpha_A \rangle is an orthonormal basis, making the schmidt values \lambda_\alpha ^2 = \frac{1}{2^{N_A}} the Schmidt values.
I thought I could compute the operator entanglement the same was as computing the von Neumann entanglement for the MPS, but I’m not getting the expected result for the maximally mixed state. Do you know what I’m doing wrong?
using ITensors, ITensorMPS
function operator_entanglement(state::Union{MPS, MPO}, bond::Int)
state = orthogonalize(state, bond)
U,S,V = svd(state[bond], (linkinds(state, bond-1)..., siteinds(state, bond)...))
SvN = 0.0
for n=1:dim(S, 1)
p = S[n,n]^2
SvN -= p * log2(p)
end
return SvN
end
let
N = 4
s = siteinds("Qubit", N)
# Operator entanglement of the maximally mixed state
max_mixed = MPO(s, "Id")/2^N
Na = 3
Sop_mix = operator_entanglement(max_mixed, Na+1)
@show maxlinkdim(max_mixed)
println("Operator Entanglement = ", Sop_mix)
println("Exact Operator Entanglement = ", Na*log2(2))
println("Does it match the expected value?", Sop_mix == Na * log2(2) )
end
gives the following output:
maxlinkdim(max_mixed) = 1
Operator Entanglement = 0.25
Exact Operator Entanglement = 3.0
Does it match the expected value?false