Universal Quantum Computing II#
Show code cell content
import mmf_setup; mmf_setup.nbinit()
This cell adds /home/docs/checkouts/readthedocs.org/user_builds/physics-555-quantum-technologies/checkouts/latest/src to your path, and contains some definitions for equations and some CSS for styling the notebook. If things look a bit strange, please try the following:
- Choose "Trust Notebook" from the "File" menu.
- Re-execute this cell.
- Reload the notebook.
Two-Qubit Operations#
In the previous assignment, you learned how to implement arbitrary single-qubit operations by controlling magnetic fields. A universal set of gates can be constructed out of these and a single two-qubit operation such as the ᴄɴᴏᴛ gate:
\[\begin{gather*}
\mathsf{\small CNOT} =
\begin{pmatrix}
1 & 0 & 0 & 0\\
0 & 1 & 0 & 0\\
0 & 0 & 0 & 1\\
0 & 0 & 1 & 0
\end{pmatrix}.
\end{gather*}\]
Implementing this gate is significantly more challenging. Here we will explore some possible ways of doing this.
Continuing with out example spin-½ qubits, consider the dipole-dipole interaction:
\[\begin{gather*}
\op{V} = \mu_{D} \vect{\op{S}}_{A} \cdot \vect{\op{S}}_{B}
= \mu_{D}\begin{pmatrix}
1\\
& -1 & 2\\
& 2 & -1\\
&&& 1
\end{pmatrix}.
\end{gather*}\]
This interaction describes the behavior of one qubit in the presence of the magnetic field generated by the spin of the other qubit and vice-versa. To simplify our results, we shall choose units so that \(\hbar = \mu_D = 1\).
Show that, in these simplified units,
\[\begin{gather*} \op{V} = \frac{1}{4} \begin{pmatrix} 1\\ & -1 & 2\\ & 2 & -1\\ &&& 1 \end{pmatrix} \end{gather*}\]by explicitly computing this from the Pauli matrices:
\[\begin{gather*} \vec{\op{S}_{A}} = \frac{\vec{\mat{σ}}}{2}\otimes \mat{1}, \qquad \vec{\op{S}_{B}} = \mat{1}\otimes \frac{\vec{\mat{σ}}}{2}. \end{gather*}\]I.e. show that
\[\begin{gather*} 4\op{V} = \sum_{i=\{x, y, z\}} (\mat{σ}_i\otimes \mat{1})(\mat{1}\otimes \mat{σ}_i)\\ = \begin{pmatrix} 0 & 1\\ 1 & 0 \\ && 0 & 1\\ && 1 & 0\\ \end{pmatrix} \begin{pmatrix} && 1 & 0\\ && 0 & 1 \\ 1 & 0\\ 0 & 1\\ \end{pmatrix} + \cdots\\ = \begin{pmatrix} 1\\ & -1 & 2\\ & 2& -1\\ &&& 1 \end{pmatrix} \end{gather*}\]Diagonalize this matrix so that you can compute the matrix exponential and find the form of the corresponding 2-qubit unitary operations. (Hint: you only need to diagonalize a 2×2 block, and this block is diagonalized by the matrix corresponding to the Hadamard gate.)
Answer
The answer can be expressed in terms of the Hadamard gate \(H\):
\[\begin{gather*} \begin{pmatrix} -1 & 2\\ 2& -1 \end{pmatrix} = H \begin{pmatrix} 1\\ & -3 \end{pmatrix} H \end{gather*}\]I.e. show that
\[\begin{gather*} \mat{U} = \exp\left(\mat{V}t/\I\hbar\right) = \begin{pmatrix} e^{\I\theta}\\ & H\begin{pmatrix} e^{\I\theta}\\ & e^{-3\I\theta} \end{pmatrix} H\\ && e^{\I\theta} \end{pmatrix}, \end{gather*}\]where \(\theta = \hbar \mu_D t/4\).
Can you combine this with single qubit gates to form ᴄɴᴏᴛ?
# Brute-force optimization to see if we can do this.
# Apparently not.
import mmf_setup; mmf_setup.nbinit()
import numpy as np
import scipy.linalg
from scipy.linalg import expm, logm, eigh
from scipy.optimize import minimize
from numpy import kron
from phys_555_2022.utils import sigmas
H = np.array([[1, 1], [1, -1]]) / np.sqrt(2)
X, Y, Z = sigmas
CNOT = np.eye(4)
CNOT[2:, 2:] = X.real
CZ = np.eye(4)
CZ[2:, 2:] = Z.real
V = np.diag([1, -1, -1, 1])
V[1,2] = V[2, 1] = 2
I = np.eye(2)
sigma_mus = np.asarray([I] + sigmas.tolist())
def unpack(y):
theta = y[-1]
bs = np.reshape(y[:-1], (4, 4))
Hs = [np.einsum('a,a...', _b, sigma_mus) for _b in bs]
Us = [expm(_H/2j) for _H in Hs]
return Us, theta
def get_residue(y):
(Ui1, Ui2, Uo1, Uo2), theta = unpack(y)
M = np.zeros((4, 4), dtype=complex)
M[1:3, 1:3] = H @ np.diag([np.exp(1j*theta), np.exp(-3j*theta)]) @ H
M[0, 0] = M[-1, -1] = np.exp(1j*theta)
U = kron(Ui1, Ui2) @ M @ kron(Uo1, Uo2)
#U *= np.exp(-1j*np.angle(U[0,0]))
assert np.allclose(U.conj().T @ U, np.eye(4))
return (abs(U - CNOT)**2).sum()
def get_res(rng=np.random.default_rng(seed=1)):
y0 = rng.random(17) * 2*np.pi
res = minimize(get_residue, x0=y0)
return res.fun
def find_match(U_goal, H, rng=np.random.default_rng(seed=1)):
"""Try to find a set of 4 single-qubit gates transforming U to U_goal.
Arguments
---------
U_goal : 4x4 array
Target 2-qubit gate.
H : 4x4 array
Starting 2-qubit Hamiltonian. Single qubit operations will be applied to
the inputs and outputs to try to match U_goal with exponentiation of H.
"""
assert np.allclose(H, H.T.conj())
D, V = eigh(H)
assert np.allclose(V @ np.diag(D) @ V.T.conj(), H)
def unpack(y):
theta = y[-1]
bs = np.reshape(y[:-1], (4, 4))
Hs = [np.einsum('a,a...', _b, sigma_mus) for _b in bs]
Us = [expm(_H/2j) for _H in Hs]
return Us, theta
def get_U(y):
(Ui1, Ui2, Uo1, Uo2), theta = unpack(y)
M = V @ np.diag(np.exp(1j*theta*D)) @ V.T.conj()
U = kron(Ui1, Ui2) @ M @ kron(Uo1, Uo2)
assert np.allclose(U.conj().T @ U, np.eye(4))
return U
def get_residue(y):
return (abs(get_U(y) - U_goal)**2).sum()
y0 = rng.random(17) * 2*np.pi
res = minimize(get_residue, x0=y0)
U = get_U(res.x)
H = logm(U)/1j
return res.fun, U, H
import numpy as np
import scipy.linalg
from scipy.linalg import expm, logm, eigh
from scipy.optimize import minimize
from numpy import kron
from phys_555_2022.utils import sigmas
H = np.array([[1, 1], [1, -1]]) / np.sqrt(2)
X, Y, Z = sigmas
CNOT = np.eye(4)
CNOT[2:, 2:] = X.real
CZ = np.eye(4)
CZ[2:, 2:] = Z.real
I = np.eye(2)
sigma_mus = np.asarray([I] + sigmas.tolist())
def find_match(U_goal, H, y0=None, rng=np.random.default_rng(seed=1)):
"""Try to find a set of 4 single-qubit gates transforming U to U_goal.
Arguments
---------
U_goal : 4x4 array
Target 2-qubit gate.
H : 4x4 array
Starting 2-qubit Hamiltonian. Single qubit operations will be applied to
the inputs and outputs to try to match U_goal with exponentiation of H.
"""
assert np.allclose(H, H.T.conj())
D, V = eigh(H)
assert np.allclose(V @ np.diag(D) @ V.T.conj(), H)
def unpack(y):
theta = y[-1]
bs = np.reshape(y[:-1], (4, 4))
Hs = [np.einsum('a,a...', _b, sigma_mus) for _b in bs]
Us = [expm(_H/2j) for _H in Hs]
return Us, theta
def get_U(y):
(Ui1, Ui2, Uo1, Uo2), theta = unpack(y)
M = V @ np.diag(np.exp(1j*theta*D)) @ V.T.conj()
U = kron(Ui1, Ui2) @ M @ kron(Uo1, Uo2)
assert np.allclose(U.conj().T @ U, np.eye(4))
return U
def get_residue(y):
return (abs(get_U(y) - U_goal)**2).sum()
if y0 is None:
y0 = rng.random(17) * 2*np.pi
res = minimize(get_residue, x0=y0)
U = get_U(res.x)
H = logm(U)/1j
return res, U, H
---------------------------------------------------------------------------
ModuleNotFoundError Traceback (most recent call last)
Cell In[2], line 6
4 from scipy.optimize import minimize
5 from numpy import kron
----> 6 from phys_555_2022.utils import sigmas
8 H = np.array([[1, 1], [1, -1]]) / np.sqrt(2)
9 X, Y, Z = sigmas
ModuleNotFoundError: No module named 'phys_555_2022'
rng = np.random.default_rng(seed=3)
V = rng.normal(size=(4, 4*2)).view(complex)
V += V.T.conj()
V = V.reshape((2, 2, 2, 2))
V = (V + np.einsum('abcd->badc', V)).reshape((4, 4))
V = np.diag([1, -1, -1, 1])
V[1,2] = V[2, 1] = 2
res, U, UH = find_match(CNOT, V, rng=rng)
for n in range(10):
res, U, UH = find_match(CNOT, UH, rng=rng)
print(res.fun)
#res, U, UH = find_match(CNOT, UH, y0=res.x, rng=rng)
#print(res.fun)
#for n in range(10):
# res, U, UH = find_match(CNOT, UH, y0=0*res.x, rng=rng)
# print(res.fun)
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Cell In[3], line 10
7 V = np.diag([1, -1, -1, 1])
8 V[1,2] = V[2, 1] = 2
---> 10 res, U, UH = find_match(CNOT, V, rng=rng)
12 for n in range(10):
13 res, U, UH = find_match(CNOT, UH, rng=rng)
NameError: name 'find_match' is not defined
(abs(U - CNOT)**2).sum()
---------------------------------------------------------------------------
NameError Traceback (most recent call last)
Cell In[4], line 1
----> 1 (abs(U - CNOT)**2).sum()
NameError: name 'U' is not defined