基于 MindQuantum 计算受控量子门的矩阵

基于 MindQuantum 计算受控量子门的矩阵

MindQuantum 当前版本 (0.9.0) 中采用了一种间接计算的方法：先将受控量子门构造成量子线路，然后在后端计算量子线路的矩阵，从而得到与之等价的受控量子门的矩阵。 本文介绍一种方法，来直接计算受控量子门的矩阵。原理上和 MindQuantum 类似，虽然不如后者优雅，但在计算效率上会有一定优势。

先来看看 MindQuantum 的实现方式。如果直接使用 API: gate.matrix，会只输出被控量子门的矩阵。比如对于一个控制为 0， 受控位为 1 的 CNOT 门：

from mindquantum import *
import numpy as np

X.on(1,0).matrix()

array([[0, 1],
       [1, 0]])

可见，只返回了单个量子比特下的 X 门的矩阵。下面将其封装在 circuit 中来计算。

circ = Circuit()
circ += X.on(1,0)
res0 = circ.matrix()
print(res0.real)

[[1. 0. 0. 0.]
 [0. 0. 0. 1.]
 [0. 0. 1. 0.]
 [0. 1. 0. 0.]]

其计算利用了矩阵列的定义：矩阵的第 $i$ 列为该矩阵作用在第 $i$ 个基矢上得到的新矢量，即 $|v_i⟩=A|i⟩$。所以，整个矩阵即为

$$A=[|v_0⟩,|v_1⟩,...,|v_{n-1}⟩]$$

以下为具体的计算流程：

n_qubits = 2
sim = Simulator('mqvector', n_qubits)
res1 = np.zeros((1<<n_qubits, 1<<n_qubits), dtype=complex)

for i in range(1<<n_qubits):
    basis = np.zeros((1<<n_qubits), dtype=complex) # 第 i 个基矢
    basis[i] = 1
    sim.set_qs(basis) 
    sim.apply_circuit(circ)
    col = sim.get_qs() # 线路作用在第 i 个基矢后得到的量子态
    res1[:,i] = col # 该量子态作为矩阵的第 i 列
    
print(res0==res1)

[[ True  True  True  True]
 [ True  True  True  True]
 [ True  True  True  True]
 [ True  True  True  True]]

关于这个 CNOT 的矩阵为什么和教科书上的不一样，是因为 MindQuantum 中采用了 little-endian 的表示方法，具体介绍可以参考 我之前的一篇文章。

当控制位为 1，而受控位为 0 时，矩阵为：

circ = Circuit()
circ += X.on(1,0)
circ.matrix().real

array([[1., 0., 0., 0.],
       [0., 0., 0., 1.],
       [0., 0., 1., 0.],
       [0., 1., 0., 0.]])

接下来，介绍我们的方法：基于子空间计算量子门的矩阵。关于子空间的概念及如何在子空间下对量子态的施加操作以节省计算量可以参考 我之前的一篇文章。

我们考虑一个受控 U 门

$$U=\left(\begin{array}{cc}{a}_{00}& a_{01}\\ a_{10}& a_{11}\end{array}\right)$$

我们考虑双量子比特的四个基矢，$|q_1{q}_0⟩$ 并且假设 $q_1$ 为控制量子比特，而 $q_0$ 为受控比特。

1. $|00⟩$，状态不变。$|0⟩\otimes |0⟩\to |0⟩\otimes |0⟩$，受控量子门的矩阵第一列由 $\left(\begin{array}{c}1\\ 0\end{array}\right)\otimes \left(\begin{array}{c}1\\ 0\end{array}\right)\to \left(\begin{array}{c}1\\ 0\end{array}\right)\otimes \left(\begin{array}{c}1\\ 0\end{array}\right)$，即 $[1,0,0,0{]}^T\to [1,0,0,0{]}^T$；

2. $|01⟩$，状态不变。$|0⟩\otimes |1⟩\to |0⟩\otimes |1⟩$，受控量子门的矩阵第二列由 $\left(\begin{array}{c}1\\ 0\end{array}\right)\otimes \left(\begin{array}{c}0\\ 1\end{array}\right)\to \left(\begin{array}{c}1\\ 0\end{array}\right)\otimes \left(\begin{array}{c}0\\ 1\end{array}\right)$，即 $[0,1,0,0{]}^T\to [0,1,0,0{]}^T$；

3. $|10⟩$，状态变为 $|0⟩\otimes (a_{00}|0⟩+a_{10}|1⟩)$。受控量子门的矩阵第三列由 $\left(\begin{array}{c}0\\ 1\end{array}\right)\otimes \left(\begin{array}{c}1\\ 0\end{array}\right)\to \left(\begin{array}{c}0\\ 1\end{array}\right)\otimes \left(\begin{array}{c}{a}_{00}\\ a_{10}\end{array}\right)$，即 $[0,0,1,0{]}^T\to [0,0,a_{00},a_{10}{]}^T$；

4. $|11⟩$，状态变为 $|0⟩\otimes (a_{01}|0⟩+a_{11}|1⟩)$。受控量子门的矩阵第四列由 $\left(\begin{array}{c}0\\ 1\end{array}\right)\otimes \left(\begin{array}{c}0\\ 1\end{array}\right)\to \left(\begin{array}{c}0\\ 1\end{array}\right)\otimes \left(\begin{array}{c}{a}_{01}\\ a_{11}\end{array}\right)$，即 $[0,0,0,1{]}^T\to [0,0,a_{01},a_{11}{]}^T$。

于是，量子门的矩阵即为

$$\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& a_{00}& a_{01}\\ 0& 0& a_{10}& a_{11}\end{array}\right)$$
比如对于 X 门由于 $X=\left(\begin{array}{cc}{a}_{00}& a_{01}\\ a_{10}& a_{11}\end{array}\right)=\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right)$

所以，X.on(0,1) 的矩阵即为
$$\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& 0& 1\\ 0& 0& 1& 0\end{array}\right)$$

同理，可以构造 $q_0$ 为控制量子比特，$q_1$ 为受控量子比特的受控 $U$ 门的矩阵为
$$\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& a_{00}& 0& a_{01}\\ 0& 0& 1& 0\\ 0& a_{10}& 0& a_{11}\end{array}\right)$$
而 X.on(1,0) 的矩阵即为
$$\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 0& 0& 1\\ 0& 0& 1& 0\\ 0& 1& 0& 0\end{array}\right)$$

下面我们实现一下，以 Y 门为例

from mindquantum import *
circ = Circuit()
circ += Y.on(1,0)
a = circ.matrix()
print(a)

[[1.+0.j 0.+0.j 0.+0.j 0.+0.j]
 [0.+0.j 0.+0.j 0.+0.j 0.-1.j]
 [0.+0.j 0.+0.j 1.+0.j 0.+0.j]
 [0.+0.j 0.+1.j 0.+0.j 0.+0.j]]

import numpy as np

matrix = np.array([[0,-1j], [1j,0]]) # Y 门的矩阵

one = np.array([0,1],dtype=complex) # |1> 态

n_qubits = 2
ctrl_qubit = 0
obj_qubit = 1

blank = np.identity(1<<n_qubits,dtype=complex)

for i in range(1<<n_qubits): # 对所有列都轮流进行操作
    index_str = bin(i)[2:].zfill(n_qubits)
    reversed_index_str = index_str[::-1] # 需要反过来，因为我们这里采用 little-endian 的表示方法
        
    if reversed_index_str[ctrl_qubit] == '1': # 当控制量子比特的状态为 1 的时候，才施加 U 门操作
        col = 1
        if reversed_index_str[obj_qubit] == '0': # 如果此时受控位为 0 则其状态变为 U 的第一列
            for j in range(n_qubits-1, -1, -1):
                if j == obj_qubit:
                    col = np.kron(col, matrix[:,0]) # 直积操作
                else:
                    col = np.kron(col, one)
        else:                                   # 如果此时受控位为 1 则其状态变为 U 的第二列
            for j in range(n_qubits-1, -1, -1):
                if j == obj_qubit:
                    col = np.kron(col, matrix[:,1]) # 直积操作
                else:
                    col = np.kron(col, one)
            
        blank[:,i] = col
                        
print(blank)

[[1.+0.j 0.+0.j 0.+0.j 0.+0.j]
 [0.+0.j 0.+0.j 0.+0.j 0.-1.j]
 [0.+0.j 0.+0.j 1.+0.j 0.+0.j]
 [0.+0.j 0.+1.j 0.+0.j 0.+0.j]]

检验一下：

a == blank

array([[ True,  True,  True,  True],
       [ True,  True,  True,  True],
       [ True,  True,  True,  True],
       [ True,  True,  True,  True]])

可见，本文方法得到的结果是正确的。