Solution (source code)

= Solution

Diagonalize the <density operator> in the <Hadamard gate> basis $|\pm\rangle=(|0\rangle\pm|1\rangle)/\sqrt2$. Its eigenvalues are $2/3,1/3$, so
$$
\rho_A=\frac23|+\rangle\langle+|+\frac13|-\rangle\langle-|.
$$
Choose orthonormal records on a second qubit to construct the first <purification of a density operator>:
$$
\boxed{|\Omega_1\rangle=\sqrt{\frac23}|+\rangle_A|0\rangle_B+\sqrt{\frac13}|-\rangle_A|1\rangle_B.}
$$
It is normalized. In its outer product the cross terms have partial trace zero because $\langle0|1\rangle_B=0$, while the diagonal terms give precisely $\rho_A$.

A different purification follows from <unitary freedom of purification>. Apply the <Hadamard gate> to $B$:
$$
\boxed{|\Omega_2\rangle=(I\otimes H)|\Omega_1\rangle
=\sqrt{\frac23}|+\rangle_A|+\rangle_B+\sqrt{\frac13}|-\rangle_A|-\rangle_B.}
$$
The same partial-trace calculation uses $\langle+|-\rangle_B=0$ and gives $\rho_A$. The two global states are distinct, not related by a global phase, although their reduced density operators agree. Explicitly the second is
$$
|\Omega_2\rangle=a(|00\rangle+|11\rangle)+b(|01\rangle+|10\rangle),\qquad
a=\frac{\sqrt2+1}{2\sqrt3},\quad b=\frac{\sqrt2-1}{2\sqrt3}.
$$
Its reduction has diagonal entries $a^2+b^2=1/2$ and off-diagonal entries $2ab=1/6$, checking all the required matrix entries.