EPR pairs (Bell states) preparation circuit in SymPy

Question

I am trying to recretate this circuit from the book: https://www.elsevier.com/books/quantum-information-processing-quantum-computing-and-quantum-error-correction/djordjevic/978-0-12-821982-9

The author does the folowing analysis:

\begin{equation} \begin{aligned} \operatorname{CNOT}(H \otimes I)\left|\psi_{1}\right\rangle\left|\psi_{2}\right\rangle &=\left[\begin{array}{cccc} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 1 \\ 0 & 0 & 1 & 0 \end{array}\right] \otimes\left(\begin{array}{c} 1 \\ \sqrt{2} \end{array}\left[\begin{array}{cc} 1 & 1 \\ 1 & -1 \end{array}\right] \otimes\left[\begin{array}{cc} 1 & 0 \\ 0 & 1 \end{array}\right]\right)\left(\left[\begin{array}{l} a_{1} \\ b_{1} \end{array}\right] \otimes\left[\begin{array}{l} a_{2} \\ b_{2} \end{array}\right]\right) \\ &=\frac{1}{\sqrt{2}}\left[\begin{array}{llll} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 1 \\ 0 & 0 & 1 & 0 \end{array}\right]\left[\begin{array}{cccc} 1 & 0 & 1 & 0 \\ 0 & 1 & 0 & 1 \\ 1 & 0 & -1 & 0 \\ 0 & 1 & 0 & -1 \end{array}\right]\left[\begin{array}{l} a_{1} a_{2} \\ a_{1} b_{2} \\ b_{1} a_{2} \\ b_{1} b_{2} \end{array}\right] \\ &=\frac{1}{\sqrt{2}}\left[\begin{array}{l} a_{1} a_{2}+b_{1} b_{2} \\ a_{1} b_{2}+b_{1} b_{2} \\ b_{1} a_{2}-b_{1} b_{2} \\ a_{1} a_{2}-b_{1} a_{2} \end{array}\right] \end{aligned} \end{equation}

Now I am using the folowing code in SymPy Quantum:

from sympy import *
from sympy import init_printing

init_printing(use_latex=True)
from sympy.physics.quantum import TensorProduct

Pauli_Z = Matrix(2, 2, [1, 0, 0, -1])
U_I_2 = Matrix(2, 2, [1, 0, 0, 1])
U_H = 1 / sqrt(2) * (Matrix(2, 2, [1, 1, 1, -1]))
Spin_up = Matrix([1, 0])
X_Spin_up = 1 / sqrt(2) * (Matrix([1, 0]) + Matrix([0, 1]))
U_I = Matrix(2, 2, [1, 0, 0, 1])

Cnot = Matrix(4, 4, [1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0])

xx = TensorProduct(X_Spin_up, Spin_up)
HU = TensorProduct(U_H, U_I_2)

Cnot, HU, xx

TensorProduct(Cnot, TensorProduct(U_H, U_I_2) * TensorProduct(X_Spin_up, Spin_up))

But the SIZE of the resulting tensor product vector is wrong:

My input is: $\psi_{1} = \frac{1}{\sqrt{2}}(|0\rangle+|1\rangle)$ and: $\psi_{2} = |0\rangle$

and if I print $Cnot, HU, xx$, I get:

\begin{equation} \left(\left[\begin{array}{cccc} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 1 \\ 0 & 0 & 1 & 0 \end{array}\right],\left[\begin{array}{cccc} \frac{\sqrt{2}}{2} & 0 & \frac{\sqrt{2}}{2} & 0 \\ 0 & \frac{\sqrt{2}}{2} & 0 & \frac{\sqrt{2}}{2} \\ \frac{\sqrt{2}}{2} & 0 & -\frac{\sqrt{2}}{2} & 0 \\ 0 & \frac{\sqrt{2}}{2} & 0 & -\frac{\sqrt{2}}{2} \end{array}\right],\left[\begin{array}{c} \frac{\sqrt{2}}{2} \\ 0 \\ \frac{\sqrt{2}}{2} \\ 0 \end{array}\right]\right) \end{equation}

What am I doing wrong?

Thanks!

Answer 1

The problem that i can see is that you are using TensorProduct between $CNOT$ and $H \otimes I$ matrices. But they should undergo normal matrix multiplication, not TensorProduct. The reason for that is that whilst $H$ and $I$ are in parallel, $\textit{CNOT}$ and $H \otimes I$ are not in parallel but in series. That's why, there should be a tensor product between the gates $H$ and $I$ but a normal matrix product between $CNOT$ and $H\otimes I$.

A SymPy code could be:

from sympy import Matrix, symbols, sqrt, init_printing
from sympy.physics.quantum import TensorProduct
from IPython.display import display_pretty

init_printing(use_latex=True)

U_I = Matrix([[1,0],
              [0,1]])
U_H = 1/sqrt(2)*Matrix([[1, 1],
                        [1,-1]])

Cnot=Matrix([[1,0,0,0],
             [0,1,0,0],
             [0,0,0,1],
             [0,0,1,0]])

U_HI = TensorProduct(U_H,U_I)

a1 = symbols('a_1')
b1 = symbols('b_1')
a2 = symbols('a_2')
b2 = symbols('b_2')

psi1 = Matrix(2,1,[a1,b1])
psi2 = Matrix(2,1,[a2,b2])
psi12 = TensorProduct(psi1,psi2)

Cnot*U_HI*psi12

And the output would be:

$\displaystyle \left[\begin{matrix}\frac{\sqrt{2} a_{1} a_{2}}{2} + \frac{\sqrt{2} a_{2} b_{1}}{2}\\\frac{\sqrt{2} a_{1} b_{2}}{2} + \frac{\sqrt{2} b_{1} b_{2}}{2}\\\frac{\sqrt{2} a_{1} b_{2}}{2} - \frac{\sqrt{2} b_{1} b_{2}}{2}\\\frac{\sqrt{2} a_{1} a_{2}}{2} - \frac{\sqrt{2} a_{2} b_{1}}{2}\end{matrix}\right]$