I calculate two-site reduced density matrix from a MPS , just like

//index linking i to i+1:

auto ir = commonIndex(psi.A(i),psi.A(i+1));

auto rho = psi.A(i)*dag(prime(psi.A(i),Site,ir));

for(int k = i+1; k < j; ++k)

{

rho *= psi.A(k);

rho *= dag(prime(psi.A(k),Link));

}

rho *= psi.A(j);

rho *= dag(prime(psi.A(j)));

for(int k = j+1; k <= psi.N(); ++k)

{

rho *= psi.A(k);

rho *= dag(prime(psi.A(k),Link));

}

How can I obtain eigenvalue of \rho