Qiskit-汉化3.10 Quantum Walk Search Algorithm

量子漫步搜索算法

量子漫步是经典马尔可夫链的量子等价,已成为许多量子算法的关键。在本节中,我们将实现一个量子漫步搜索算法,用于在图中查找标记元素。与经典算法相比,该算法具有二次加速。

1. 经典马尔可夫链

马尔可夫链是一种随机过程,常用于模拟现实生活中的过程。它由状态和相关的转换概率组成,它描述了在每个时间步中状态之间移动的概率。在离散马尔可夫链中,我们在这里研究的,时间步长是离散的。马尔可夫链满足马尔可夫性质,这意味着该过程的下一步只取决于当前,而不是之前的任何步骤。马尔可夫链有一个相关的转移矩阵 P,它描述了在每个状态之间移动的概率。下面我们展示一个马尔可夫链及其相关转移矩阵 P 的例子。

P= \begin{pmatrix} 0.1 & 0.3 & 0.3\\ 0.1 & 0.1 & 0.2 \\ 0.8 & 0.6 & 0.5 \end{pmatrix} \label{eq:matrix_example}

1

给定转移矩阵 P ,我们可以得到 t 时间步长按 P^t 的概率分布。

2. 量子漫步

量子漫步相当于经典的马尔可夫链。由于叠加,在我们测量电路之前,量子漫步将同时采用所有可能的路径。由于量子干涉,某些态会相互抵消。这使得量子漫步算法比经典算法更快,因为我们可以以这样的方式设计它们,错误的答案很快就会抵消。量子漫步有两种常见的模型,硬币态量子漫步和 Szegedy 量子漫步,它们在某些情况下是等价的。硬币态量子漫步是在图的顶点上行走,而 Szegedy 量子漫步发生在图的边缘上。在展示如何实现量子漫步之前,我们将介绍这两个模型。

2.A 硬币态量子漫步

硬币态量子漫步的一个简单例子是在无限整数线上行走。在这种情况下,我们用一个整数 \{\ket{j}: j \ In \mathbb{Z} \} 表示漫步者的位置,因为漫步者可以遍历 \mathbb{Z} 中的所有整数。一枚硬币决定漫步者应该如何移动。在整数线上,漫步者可以向左或向右走。因此,这个硬币的计算基础是 \{\ket{0}, \ket{1}\} ,如果硬币是 \ket{0} ,我们就朝一个方向移动漫步者,如果硬币是 \ket{1} ,我们就朝另一个方向移动漫步者。

硬币态量子是在图中的节点上行走,我们将节点称为状态。漫步者可以在与边相连的状态之间移动。在硬币模型中,我们有两个量子态和两个算符。第一个状态是位置状态,它表示漫步者的位置。对于上面的漫步,这是一个整数,因为漫步可以在整数线上的任何位置。另一种状态是硬币态。硬币状态决定漫步者下一步应该如何移动。我们可以用希尔伯特空间中的向量来表示硬币态和位置态。如果我们可以用 \mathcal{H}_C 表示硬币状态,用 \mathcal{H}_P 表示位置状态,我们就可以用 \mathcal{H} = \mathcal{H}_C \otimes \mathcal{H}_P 表示整个漫步者的量子空间。

正如我们提到的,模型也有两个算符;硬币算符 C 和移位算符 S 。硬币算符在每个时间步骤中作用于 \mathcal{H}_C ,并将漫步者置于叠加状态,因此它同时行走所有可能的路径。对于整数直线上的漫步,这意味着它在每个时间步长中同时向左或向右移动。有不同的硬币算符,但最常见的是Hadamard硬币算符和Grover硬币算符。Hadamard硬币算符是一个哈达玛门,它使漫步者处于一个相等的叠加状态:

H = \frac{1}{\sqrt{2}} \begin{bmatrix} 1 & 1 \\ 1 & -1 \end{bmatrix}

Grover硬币是Grover算法中的Grover扩散算子。我们把它定义为

G = \begin{bmatrix} \frac{2}{n} -1 & \frac{2}{n} & \ldots & \frac{2}{n}\\ \frac{2}{n} & \frac{2}{n} - 1 & \ldots & \frac{2}{n} \\ \vdots & \vdots & \ddots & \vdots \\ \frac{2}{n} & \frac{2}{n} & \ldots & \frac{2}{n} -1 \end{bmatrix}

和Hadamard硬币一样,Grover硬币让漫步者处于叠加状态。然而,它的行为有点不同。如果我们将Grover硬币应用于位置为 \ket{000} 的漫步者,我们将获得下图所示的状态向量概率。正如我们所看到的,这个硬币并没有像Hadamard硬币那样把漫步者放在一个相等的叠加位置上。相反, \ket{000} 的概率比其他状态大得多。

2

模型中的另一个算符,即移位算符,作用于 \mathcal{H}_P 并将步行器移动到下一个位置。在整数线上行走时,如果硬币为 \ket{0} ,则移位算符将步行者向左移动;如果硬币为 \ket{1} ,则向右移动。

S \ket{0}\ket{j} = \ket{0}\ket{j+1}
S \ket{1}\ket{j} = \ket{1}\ket{j-1}

通过上面定义的移位算符,我们可以将硬币态量子的一步酉算符 U 表示为

U = SC,

其中C是硬币算符。对于整数线上的量子漫步,我们使用Hadamard硬币,但C可以是Hadamard硬币、Grover硬币或任何其他硬币算符。

我们还可以向前看几步。我们可以将量子态 \ket{\psi} 表示在 t 时间步长之后为

\ket{\psi (t)} = U^t \ket{\psi(0)},

其中 \ket{\psi(0)} 是初始状态,U是漫步[1]中第一步的算子。

硬币态量子漫步最适合于规则图,即所有节点都有相同数量的邻居[2]的图。对于非规则图来说,另一种更方便的量子漫步模型是我们接下来要看的Szegedy模型。

2.B Szegedy量子漫步

硬币态量子漫步是在图的节点上行走,而Szegedy量子漫步是在原始图的二分双覆盖的边缘上行走。双覆盖图是一个顶点数是原始图的两倍的图。当且仅当二分双覆盖图中的两个节点在原图中也连通时,这两个节点与一条边相连。为了创建这个模型,我们从经典行走的转移概率矩阵P开始。如第1节所述,我们用一个转换矩阵 P 表示一个经典的离散时间随机游走。对于任何 N -顶点图, N \ * N 转换矩阵 P ,我们可以将相应的离散时间量子行走定义为希尔伯特空间 \mathcal{H}^N \otimes \mathcal{H}^N 上的酉运算。让 P_{jk} 定义从状态 j 转换到状态 k 的概率。在定义行走之前,我们先定义标准化的状态

\ket{\psi_j} := \sum_{k=1}^N \sqrt{P_{kj}} \ket{j,k}, \; j=1,...,N

以及在 {\ket{\psi_j}}上的投影:j=1,…,N

\Pi := \sum_{j=1}^N \ket{\psi_j} \bra{\psi_j} \label{eq:sz_pi}

我们还引入了移位算子S:

S := \sum_{j,k=1}^N \ket{j,k} \bra{k,j} \label{eq:sz_s}

有了上面定义的 S \Pi ,我们可以引入离散量子漫步的一个步:

U := S(2 \Pi - 1), \label{eq:sz_op}

其中 (2 \Pi - 1) 是反射算符。我们还将行走的 t 步数定义为 U^t [2]。

2.C 硬币态和Szegedy量子漫步的等价性

众所周知,用Grover硬币创造的漫步相当于Szegedy量子漫步。要了解更多细节,请参阅Thomas G. Wong的文章[3],其中他还展示了两个模型中算符之间的等价性。

3.示例:在超立方体上实现量子漫步

超立方体是 n 维的3维立方体的对应物。所有节点的度都是 n ,超立方体总共有 n =2^n 个节点。我们可以用 n -元组的二进制数来表示超立方体图中的节点。一个节点的邻居的二进制表示将仅相差一个二进制数。例如,在4维超立方体中, 0000 的邻居是 0001 0010 0100 1000 。因此,一个节点连接到汉明距离为1的所有节点。边缘也被标记。第a:th位上不同的两个相邻节点由标记为 a 的边连接。

表示超立方体上硬币态量子漫步的希尔伯特空间是 \mathcal{H} = \mathcal{H}^n \otimes \mathcal{H}^{2^n} ,其中 \mathcal{H}^n 表示硬币空间, \mathcal{H}^{2^n} 表示漫步者的位置。计算基础是

\big\{ \ket{a,\vec{v}}, 0 \leq a \leq n-1, \vec{v} \in \{(00...00), (00...01), ....., (11...11 )\} \big\}.

硬币计算基础 a 的价值,与边 a 相关联,决定漫步者应该移动到哪里。如果 a=0 , 漫步者将转到第一个二进制值与当前节点不同的节点。如果 a=1 , 漫步者将移动到第二个值与当前值不同的节点,以此类推。设 \vec{e}_a 为n元组,其中除索引为 a 的值外,所有二进制值均为 0 。然后移位算符 S 将漫步者从状态 \ket{\vec{v} \oplus \vec{e}_a} 移动:

S \ket{a} \ket{\vec{v}} = \ket{a} \ket{\vec{v} \oplus \vec{e}_a}.

我们用Grover硬币, G ,来进行这次行走。因此,演化算符为

U = SG.

现在我们将展示如何在四维超立方体上实现量子漫步。我们需要实现硬币算符和移位算符。首先从Qiskit导入所有必要的库。

# Importing standard Qiskit libraries
from qiskit import QuantumCircuit, execute, Aer, IBMQ, QuantumRegister, ClassicalRegister
from qiskit.compiler import transpile, assemble
from qiskit.tools.jupyter import *
from qiskit.visualization import *
from qiskit.circuit.library import QFT
from numpy import pi
from qiskit.quantum_info import Statevector
from matplotlib import pyplot as plt
import numpy as np
# Loading your IBM Q account(s)
provider = IBMQ.load_account()

电路将有6个量子比特,4个代表位置,2个代表硬币。如前所述,这枚硬币是Grover硬币,它是Grover算法中的扩散器。我们从实现这个开始。

one_step_circuit = QuantumCircuit(6, name=' ONE STEP') 
# Coin operator
one_step_circuit.h([4,5])
one_step_circuit.z([4,5])
one_step_circuit.cz(4,5)
one_step_circuit.h([4,5])
one_step_circuit.draw() 

3

现在,让我们实现移位算符。我们知道漫步者只能移动到邻近节点,并且所有邻近节点仅相差1位。我们想要根据硬币移动漫步者,我们通过对其中一个节点量子比特应用NOT门来移动漫步者。如果硬币的状态为 \ket{11} ,我们将漫步者移动到第一个节点量子比特不同的状态。如果硬币是 \ket{10} \ket{01} ,步行者将分别移动到第二个和第三个量子比特不相同的状态。最后,如果Grover硬币为 \ket{00} ,我们翻转第四个量子比特。我们在Grover币之后用CCNOT和NOT门实现了这一点。它们合在一起,就是量子漫步在四维超立方体上的一步。

# Shift operator function for 4d-hypercube
def shift_operator(circuit):
    for i in range(0,4):
        circuit.x(4)
        if i%2==0:
            circuit.x(5)
        circuit.ccx(4,5,i)

shift_operator(one_step_circuit)

one_step_gate = one_step_circuit.to_instruction() 
one_step_circuit.draw() 

4

4. 量子漫步搜索算法

我们现在将实现一个量子漫步搜索算法,在图中找到一个标记的顶点。首先,我们描述了算法,然后我们讨论了它的实现。量子漫步搜索算法解决了通过量子行走在图中寻找标记顶点的问题。也就是说,我们标记一些顶点集合 |M| ,从图中的任意节点开始,根据行走轨迹移动,直到找到标记的节点。量子漫步搜索算法中的基态有两个寄存器,一个对应于当前节点,另一个对应于前一个节点。也就是说,基态对应于图中的边。用 \mathcal{H} 上的酉运算 W(p) 表示基于经典马尔可夫链的具有转移矩阵 p 的量子漫步。我们还定义了 \ket{p_x} = \sum_y \sqrt{P_{xy}}\ket{y} 作为节点 x 的邻居的一致叠加。设 \ket{x}\ket{y} 为基态。如果 x 是一个标记的节点,我们将基状态 \ket{x}\ket{y} 定义为“好”。否则,我们称之为“坏”。现在我们介绍“好”和“坏”状态:

\ket{G} = \frac{1}{\sqrt{|M|}} \sum_{x \in M} \ket{x} \ket{p_x}, \; \ket{B} = \frac{1}{\sqrt{N-|M|}} \sum_{x \notin M} \ket{x} \ket{p_x},

也就是好基态和坏基态的叠加。接下来,让我们定义 \epsilon = |M|/N \theta = \arcsin(\sqrt{\epsilon})
简而言之,该算法由三步组成:

  1. 设置初始状态 \ket{U} = \frac{1}{\sqrt{N}} \sum_{x} \ket{x} \ket{p_x} = \sin{\theta} \ket{G} + \cos{\theta} \ket{B} ,所有边的均匀叠加

  2. 重复执行 O(1/\sqrt{\epsilon}) 次:

    (a) 通过 \ket{B} 进行反射

    (b) 通过 \ket{U} 进行反射

  3. 在计算的基础上做测量

我们可以很容易地用哈达玛门实现步骤 1 ,通过 \ket{B} 反射,使用相位 oracle ,如果 x 在第一个寄存器中,则移位 x 的相位,否则保持电路不变。

步骤2(b)等价于找到一个酉 R(P) ,执行如下映射:

\begin{align} \label{eq:mapping_1} \ket{U} &\mapsto \ket{U}, \: \text{和} \\ \ket{\psi} &\mapsto -\ket{\psi}, \: \forall \ket{\psi} \text{的特征向量张成的空间 $ W(P) $ 正交于 $ \ket{U} $ } \label{eq:mapping_2} \end{align}

为了找到这个算子,我们对 W(P) 应用相位估计。上面我们定义了 W(P) 作为随机游走的演化算子。
正如我们在 2.A 节中看到的,这是一个酉算子酉算符。由此得出 W(P) 的特征值有范数 1 。因此,我们可以将 W(P) 的特征值写成 e^{\pm 2i\theta_j} 的形式。酉 W(P) 有一个特征向量和对应的特征值
1 ,即 \ket{U} 。这由 \theta_1=0 给出, R(P) 将通过添加一个带有辅助量子比特的寄存器来找到这个向量 \ket{U} ,并以 O(1/\sqrt{\delta}) 的精度执行相位估计,其中 δ P 的谱隙。为此,我们需要应用 W(P) O(1/\sqrt{\delta}) 次。设 \ket{w} w (P) 的特征向量,特征值 e^{\pm 2i\theta_j} 。假设 \tilde{\theta_j} 是由相位估计给出的 \tilde{\theta_j} 的最佳近似。操作 R(P)
为步骤2(b)中的 \ket{w} 执行映射,由[4]给出

\ket{w} \ket{0} \mapsto \ket{w} \ket{\tilde{\theta_j}} \mapsto (-1)^{|\tilde{\theta_j} \neq 0|} \ket{w} \ket{\tilde{\theta_j}} \mapsto (-1)^{|\tilde{\theta_j} \neq 0|} \ket{w} \ket{0}

5. 示例:4维超立方体上的量子漫步搜索

量子漫步搜索算法可以在 O(1/\sqrt{\epsilon}) 步中找到一组标记的节点, \epsilon = |M|/N ,其中 M 是标记的节点数, N 是节点总数。该算法最初用于Szegedy量子漫步,其中我们使用两个节点寄存器来表示量子态。然而,使用Grover硬币的硬币态漫步相当于Szegedy量子漫步,因为硬币态漫步的实现通常不太复杂,所以我们选择使用硬币态漫步来实现算法。我们将使用在第3节中实现的4维超立方体。

简而言之,我们将实现如下算法。通过对节点量子比特和硬币量子比特应用哈达玛门,我们实现了第1步,在所有边缘上均匀叠加。对于步骤2(a),我们实现了一个阶段oracle。步骤2(b)通过在超立方体上量子行走的一个步骤上进行相位估计来实现,然后标记 \theta \neq 0 的所有量子态。我们通过旋转一个辅助量子比特来做到这一点。在这一步的最后一部分,我们反转相位估计。θ量子比特的数量取决于 \theta 的精度。

下面,我们在 4 维超立方体上实现量子漫步搜索算法。

对于这个算法,我们需要使用之前实现的一步门的逆运算。我们通过使用内置电路函数 inverse() 得到它。

one_step_circuit.inverse().draw() 

5

反向一步门将用于后面的相位估计。我们需要根据我们在第 3 节中实现的一步门及其逆过程来制作受控门。我们稍后将根据控制量子比特的值使用它们。

# Make controlled gates
inv_cont_one_step = one_step_circuit.inverse().control()
inv_cont_one_step_gate = inv_cont_one_step.to_instruction()
cont_one_step = one_step_circuit.control()
cont_one_step_gate = cont_one_step.to_instruction()

受控一步门和受控逆一步门都将用于相位估计。我们将在相位估计中使用的另一件事是量子傅立叶变换。Qiskit 有一个函数 QFT,它实现了量子傅立叶变换。相位估计使用了量子傅里叶逆变换,但我们还需要使用普通的 QFT 来反转相位估计。

inv_qft_gate = QFT(4, inverse=True).to_instruction()  
qft_gate = QFT(4, inverse=False).to_instruction()

QFT(4, inverse=True).decompose().draw("mpl")

6

在我们实现相位估计之前,我们实现了一个标记状态 1011 和 1111 的相位 oracle。然后,我们将它做成一个电路。这是算法的步骤 2(a)。

phase_circuit =  QuantumCircuit(6, name=' phase oracle ')
# Mark 1011
phase_circuit.x(2)
phase_circuit.h(3)
phase_circuit.mct([0,1,2], 3)
phase_circuit.h(3)
phase_circuit.x(2)
# Mark 1111
phase_circuit.h(3)
phase_circuit.mct([0,1,2],3)
phase_circuit.h(3)
phase_oracle_gate = phase_circuit.to_instruction()
# Phase oracle circuit
phase_oracle_circuit =  QuantumCircuit(11, name=' PHASE ORACLE CIRCUIT ')
phase_oracle_circuit.append(phase_oracle_gate, [4,5,6,7,8,9])
phase_circuit.draw() 

7

我们现在将实现一个门,如果其他量子比特非零,则旋转辅助量子比特。我们将在相位估计中使用此门,如果 \theta \neq 0 ,它将旋转辅助量子比特。

# Mark q_4 if the other qubits are non-zero 
mark_auxiliary_circuit = QuantumCircuit(5, name=' mark auxiliary ')
mark_auxiliary_circuit.x([0,1,2,3,4])
mark_auxiliary_circuit.mct([0,1,2,3], 4)
mark_auxiliary_circuit.z(4)
mark_auxiliary_circuit.mct([0,1,2,3], 4)
mark_auxiliary_circuit.x([0,1,2,3,4])

mark_auxiliary_gate = mark_auxiliary_circuit.to_instruction()
mark_auxiliary_circuit.draw()

8

现在,我们将实现算法的第2(b)步。这一步包括相位估计,量子行走的一步,接着是一个辅助量子比特,如果 \theta \neq 0 ,我们就旋转它。为此,我们使用我们刚刚创建的mark_auxiliary_gate。然后,我们反转相位估计。

# Phase estimation
phase_estimation_circuit = QuantumCircuit(11, name=' phase estimation ')
phase_estimation_circuit.h([0,1,2,3])
for i in range(0,4):
    stop = 2**i
    for j in range(0,stop):
        phase_estimation_circuit.append(cont_one_step, [i,4,5,6,7,8,9])

# Inverse fourier transform
phase_estimation_circuit.append(inv_qft_gate, [0,1,2,3])

# Mark all angles theta that are not 0 with an auxiliary qubit
phase_estimation_circuit.append(mark_auxiliary_gate, [0,1,2,3,10])

# Reverse phase estimation
phase_estimation_circuit.append(qft_gate, [0,1,2,3])   

for i in range(3,-1,-1):
    stop = 2**i
    for j in range(0,stop):
        phase_estimation_circuit.append(inv_cont_one_step, [i,4,5,6,7,8,9])
phase_estimation_circuit.barrier(range(0,10))
phase_estimation_circuit.h([0,1,2,3])

# Make phase estimation gate
phase_estimation_gate = phase_estimation_circuit.to_instruction()
phase_estimation_circuit.draw() 

现在我们使用之前制作的门来实现整个量子漫步搜索算法。我们首先将哈达玛门门应用于节点和硬币量子比特,这是算法的第一步。然后,我们迭代地应用相位oracle门和相位估计门(步骤2(a)和2(b))。我们需要 O(1/\sqrt{\epsilon}) 迭代,如第4节中对算法的描述所述。最后,我们测量节点量子比特。

# Implementation of the full quantum walk search algorithm
theta_q = QuantumRegister(4, 'theta')
node_q = QuantumRegister(4, 'node')
coin_q = QuantumRegister(2, 'coin')
auxiliary_q = QuantumRegister(1, 'auxiliary')
creg_c2 = ClassicalRegister(4, 'c')
circuit = QuantumCircuit(theta_q, node_q, coin_q, auxiliary_q, creg_c2)
# Apply Hadamard gates to the qubits that represent the nodes and the coin
circuit.h([4,5,6,7,8,9])
iterations = 2

for i in range(0,iterations):
    circuit.append(phase_oracle_gate, [4,5,6,7,8,9])
    circuit.append(phase_estimation_gate, [0,1,2,3,4,5,6,7,8,9,10])

circuit.measure(node_q[0], creg_c2[0])
circuit.measure(node_q[1], creg_c2[1])
circuit.measure(node_q[2], creg_c2[2])
circuit.measure(node_q[3], creg_c2[3])
circuit.draw()

最后,我们在qasm模拟器上运行实现。我们看到电路在大多数情况下都坍缩到标记状态。

backend = Aer.get_backend('qasm_simulator') 
job = execute( circuit, backend, shots=1024 ) 
hist = job.result().get_counts() 
plot_histogram( hist )

11

6. 参考文献

  1. Renato Portugal. Quantum Walks and Search Algorithms. New York, NY: Springer New York, 2013
  2. Markus G. Kuhn.Some Introductory Notes on Quantum Computing. Apr. 2000
  3. Thomas G. Wong. “Equivalence of Szegedy’s and coined quantum walks”. In: Quantum InformationProcessing 16.9 (July 2017). ISSN: 1573-1332. DOI:10.1007/s11128-017-1667-y. URL:http://dx.doi.org/10.1007/s11128-017-1667-y.37
  4. Ronald de Wolf. Quantum Computing: Lecture Notes. 2021. arXiv:1907.09415 [quant-ph]
import qiskit.tools.jupyter
%qiskit_version_table

版本信息

Qiskit Software Version
qiskit-terra 0.18.0
qiskit-aer 0.8.2
qiskit-ignis 0.6.0
qiskit-ibmq-provider 0.15.0
qiskit-aqua 0.9.4
qiskit 0.28.0
System information
Python 3.7.7 (default, May 6 2020, 04:59:01) [Clang 4.0.1 (tags/RELEASE_401/final)]
OS Darwin
CPUs 8
Memory (Gb) 32.0

Tue Oct 05 13:51:31 2021 BST