How to calculate the Schmidt decomposition of a state without SVD

Question

I have this state of two qubits here: $$ |\psi_{AB}\rangle = \frac{1}{2}(|0\rangle_A |0\rangle_B + |1\rangle_A |1\rangle_B + |1\rangle_A |0\rangle_B - |0\rangle_A |1\rangle_B) $$ Which means that the density matrix (with order $|00\rangle$, $|01\rangle$, $|10\rangle$, $|11\rangle$) is: $$ \rho_{AB} = \frac{1}{4}\begin{pmatrix}1&-1&1&1\\-1&1&-1&-1\\1&-1&1&1\\1&-1&1&1\end{pmatrix} $$ I have to calculate the Schmidt decomposition WITHOUT USING THE SVD. What I did was calculating the reduced density matrix which at the end are: $$ \rho_A = \text{Tr}_B(\rho_{AB}) = \frac{1}{2}I \qquad\qquad \rho_B = \text{Tr}_A(\rho_{AB}) = \frac{1}{2}I $$ With this I get the Schmidt coefficients which are $\left\{\frac{1}{2}, \frac{1}{2}\right\}$ with eigenstates $\left\{|0\rangle_A, |1\rangle_A\right\}$ for $\rho_A$ and $\left\{|0\rangle_B, |1\rangle_B\right\}$ for $\rho_B$.

So my decomposition is (in the computational basis) $$ |\psi\rangle_{AB} = \frac{1}{2} (|0\rangle_A |0\rangle_B + |1\rangle_A |1\rangle_B) $$ But this is different than the first state! I would expect the correct decomposition to be: $$ |\psi\rangle_{AB} = \frac{1}{2} (|+\rangle_A |0\rangle_B - |-\rangle_A |1\rangle_B) $$ Is it wrong to assume that the eigenstates of the reduced density matrix are part of the Schmidt basis? How can I calculate the Schmidt basis in a correct way?

Answer 1

Why do you need to calculate the Schmidt decomposition without SVD? Doing it in other ways will just lead to the same result.

That aside, the problem with the calculation here is that knowing the reduced states doesn't uniquely identify the full state. The full state must be some purification of the reduced ones, but you can't know which one of those. In other words, studying the eigenvalues of the reduced states correctly gives you the Schmidt coefficients, but does not necessarily tell you the associated states in the decomposition. More precisely, knowing that a reduced state is $\rho_A$ only tells you that the original bipartite state has the form $\operatorname{vec}(\sqrt{\rho_A}V^\dagger)$ for some isometry $V$.

The correct approach to get the Schmidt vectors would be to compute the eigendecomposition of $\Psi\Psi^\dagger$ and $\Psi^\dagger\Psi$, with $\Psi$ the original state thought of as a matrix (that is, the linear operator such that $|\Psi\rangle=\operatorname{vec}(\Psi)$). But of course, this is just the same thing as doing the SVD of $|\Psi\rangle$.

It might be interesting to note that the approach of computing the eigendecomposition of reduced states does work, sometimes. More precisely, you get the correct Schmidt vectors iff the Schmidt coefficients are non-degenerate (equivalently, iff the SVD of the original state are non-degenerate). To see it, write the SVD of the original state as $$|\Psi\rangle = \sum_k \sqrt{p_k} |u_k\rangle\otimes|v_k\rangle$$ for some pair of orthonormal vectors $\{|u_k\rangle\}_k$ and $\{|v_k\rangle\}_k$. Then the reduced states are $$\rho_A = \sum_k p_k |u_k\rangle\!\langle u_k|.$$ This relation might lead you to believe that the eigendecomposition of $\rho_A$ gives you back the vectors $|u_k\rangle$. After all, the above equation is an eigendecomposition for $\rho_A$. However, when the singular values $(p_k)$ are not all distinct, it is not the only such decomposition. For exactly the same reason why the eigenvectors of a linear operator are not uniquely defined when there are degenerate eigenvalues. Whatever method you use to compute the eigenvectors will necessarily make some choice, and won't generally give you back $|u_k\rangle$, but rather some (unitary) linear combination of them.

You can also see it directly yourself with your example: try adding an arbitrarily small perturbation and apply the same procedure to the state $$|+,0\rangle+(1+\epsilon)|-,1\rangle$$ and you should get the correct results.

Answer 2

When a question says "without using the SVD" it may be implying that there are easier ways, at least in this specific case. It's always worth just having a stare at the state for a minute (as I suspect you've done given that you have spotted the answer). I generally start just grouping things based on the state of the first qubit: $$ \frac12(|0\rangle(|0\rangle+|1\rangle)+|1\rangle(|0\rangle-|1\rangle)) $$ which is the same as $$ \frac{1}{\sqrt{2}}(|0+\rangle+|1-\rangle) $$ (make sure you're careful with normalisation). Emphasise that this is the Schmidt decomposition because of the orthonormal basis on each subsystem, and so you can just read off what the Schmidt coefficients are.

(I realise this doesn't answer your question about why your method was failing, but @glS did that already)