量子近似优化算法
概述
量子近似优化算法(Quantum Approximate Optimization Algorithm,QAOA)是利用量子计算机来近似解决组合优化问题的量子算法,最早由Farhi等人于2014年提出。在本教程里,我们将利用QAOA算法来解决最大割问题(Max-Cut),来熟悉MindQuantum中量子线路的搭建和训练。
本文档适用于CPU环境。
环境准备
如果没有安装MindQuantum,可通过如下命令安装MindQuantum。
[ ]:
!pip install https://ms-release.obs.cn-north-4.myhuaweicloud.com/1.4.0/MindQuantum/any/mindquantum-0.3.0-py3-none-any.whl -i https://pypi.tuna.tsinghua.edu.cn/simple
本教程所需要的额外库:
networkx
networkx
是创建、操作和研究复杂网络的结构、动态和功能库。可通过pip3 install networkx
来进行安装。
Max-Cut问题描述
Max-Cut问题是图论中的一个NP-complete问题,它需要将一个图中的顶点分成两部分,并使得两部分被切割的边最多。如下图(a),一个图由五个顶点构成,相互连接的边为(0, 1), (0, 2), (1, 2), (2, 3), (3, 4), (0, 4)
。为了使得被切割的边最多,我们尝试通过(b)图的分割,将1、2、4分为一组,0、3分成另一组,因此可得到被切割的边有5条。当图中顶点增多时,我们很难找到有效的经典算法来解决Max-Cut问题。下面,我们介绍怎么将Max-Cut问题转化为一个哈密顿量的基态能力求解问题。
Max-Cut问题量子化
这里我们将图中的每个顶点赋予一个量子比特,当顶点被分到左边时,我们将该顶点上的量子比特设置为\(\left|0\right>\)态,同理,右边为\(\left|1\right>\)态,当两个顶点被分到不同的集合中时,这两个顶点上的比特将处于不同的量子态。例如对于第0个顶点和第1个顶点,当其连线被切割是,两个顶点上的比特对应的量子态可以为\(\left|\psi\right>=\left|0_11_0\right>\)或\(\left|\psi\right>=\left|1_10_0\right>\),其中下角标表示顶点的序号。此时,我们选择哈密顿量\(H=(Z_1Z_0-1)/2\),这里\(Z\)为泡利\(Z\)算符。不难发现:
而当顶点被分到同一集合中是,不难验证此时:
因此,我们只用按照上面的规则,写出图对应的哈密顿量\(H\),利用量子计算机求得\(H\)的基态能量与基态,我们就可以得到该图的Max-Cut切割方案与最大切割边数。我们记所有边的集合为\(C\),所有边个数为\(c\),则哈密顿量可写为:
导入相关依赖
[1]:
from mindquantum.circuit import Circuit, UN, StateEvolution
from mindquantum.gate import Hamiltonian, H, ZZ, RX
from mindquantum.nn import MindQuantumAnsatzOnlyLayer
from mindquantum.ops import QubitOperator
import networkx as nx
import mindspore.nn as nn
import numpy as np
import matplotlib.pyplot as plt
搭建所需求解的图
通过add_path
可在图中添加边。最后画出图的结构。
[2]:
g = nx.Graph()
nx.add_path(g, [0,1])
nx.add_path(g, [1,2])
nx.add_path(g, [2,3])
nx.add_path(g, [3,4])
nx.add_path(g, [0,4])
nx.add_path(g, [0,2])
nx.draw(g,with_labels=True, font_weight='bold')
如上如,我们得到一个由5个节点和6条边构成的图结构。
搭建QAOA量子线路
线路搭建
这里我们采用量子绝热近似算法,经过演化将量子态从\(X^{\otimes n}\)的本征态演化到图对应哈密的量的基态。
搭建图对应哈密顿量的含时演化线路:
[3]:
def build_hc(g,para):
hc = Circuit()
for i in g.edges:
hc += ZZ(para).on(i)
return hc
搭建\(X^{\otimes n}\)含时演化的量子线路:
[4]:
def build_hb(g, para):
hc = Circuit()
for i in g.nodes:
hc += RX(para).on(i)
return hc
为了使得最后优化的结果足够准确,我们需要将量子线路重复多次,因此我们通过如下函数搭建多层的训练网络:
[5]:
def build_ansatz(g, p):
c = Circuit()
for i in range(p):
c += build_hc(g,f'g{i}')
c += build_hb(g,f'b{i}')
return c
构建图对应的哈密顿量:
[6]:
def build_ham(g):
hc = QubitOperator()
for i in g.edges:
hc += QubitOperator(f'Z{i[0]} Z{i[1]}')
return hc
生成完整的量子线路和图所对应的哈密顿量
这里我们选择p = 4
,表示选用4层的QAOA量子线路,ansatz
是求解该问题的量子线路,init_state_circ
是将量子态制备到均匀叠加态上的量子线路。
[7]:
p = 4
ham = Hamiltonian(build_ham(g))
ansatz = build_ansatz(g, p)
init_state_circ = UN(H, g.nodes)
搭建待训练量子神经网络
由于该问题不需要编码层量子线路,我们这里使用MindQuantumAnsatzOnlyLayer
作为待训练的量子神经网络,并采用Adam
优化器。
[8]:
net = MindQuantumAnsatzOnlyLayer(ansatz.para_name, init_state_circ+ansatz, ham)
opti = nn.Adam(net.trainable_params(), learning_rate=0.05)
train_net = nn.TrainOneStepCell(net, opti)
训练并展示结果
[9]:
for i in range(600):
if i%10 == 0:
print("train step:", i, ", cut:", (len(g.edges)-train_net())/2)
train step: 0 , cut: [[3.0059216]]
train step: 10 , cut: [[3.3262742]]
train step: 20 , cut: [[3.7228582]]
train step: 30 , cut: [[3.983411]]
train step: 40 , cut: [[4.135832]]
train step: 50 , cut: [[4.216693]]
train step: 60 , cut: [[4.2141833]]
train step: 70 , cut: [[4.2036085]]
train step: 80 , cut: [[4.260594]]
train step: 90 , cut: [[4.373112]]
train step: 100 , cut: [[4.4853263]]
train step: 110 , cut: [[4.5553446]]
train step: 120 , cut: [[4.587566]]
train step: 130 , cut: [[4.611128]]
train step: 140 , cut: [[4.637698]]
train step: 150 , cut: [[4.6584387]]
train step: 160 , cut: [[4.66508]]
train step: 170 , cut: [[4.663408]]
train step: 180 , cut: [[4.6678705]]
train step: 190 , cut: [[4.6875486]]
train step: 200 , cut: [[4.7206187]]
train step: 210 , cut: [[4.7580614]]
train step: 220 , cut: [[4.7893686]]
train step: 230 , cut: [[4.8074245]]
train step: 240 , cut: [[4.8116426]]
train step: 250 , cut: [[4.8077316]]
train step: 260 , cut: [[4.803544]]
train step: 270 , cut: [[4.8039436]]
train step: 280 , cut: [[4.8088512]]
train step: 290 , cut: [[4.8154163]]
train step: 300 , cut: [[4.821649]]
train step: 310 , cut: [[4.8281393]]
train step: 320 , cut: [[4.8366113]]
train step: 330 , cut: [[4.847317]]
train step: 340 , cut: [[4.858108]]
train step: 350 , cut: [[4.865946]]
train step: 360 , cut: [[4.8693476]]
train step: 370 , cut: [[4.869488]]
train step: 380 , cut: [[4.868954]]
train step: 390 , cut: [[4.8695197]]
train step: 400 , cut: [[4.8711824]]
train step: 410 , cut: [[4.8730283]]
train step: 420 , cut: [[4.874686]]
train step: 430 , cut: [[4.8768916]]
train step: 440 , cut: [[4.880748]]
train step: 450 , cut: [[4.8865013]]
train step: 460 , cut: [[4.8930907]]
train step: 470 , cut: [[4.898922]]
train step: 480 , cut: [[4.9031305]]
train step: 490 , cut: [[4.906122]]
train step: 500 , cut: [[4.9088955]]
train step: 510 , cut: [[4.9119415]]
train step: 520 , cut: [[4.9149566]]
train step: 530 , cut: [[4.9175825]]
train step: 540 , cut: [[4.920064]]
train step: 550 , cut: [[4.9228735]]
train step: 560 , cut: [[4.925872]]
train step: 570 , cut: [[4.9282985]]
train step: 580 , cut: [[4.929679]]
train step: 590 , cut: [[4.930426]]
根据上面的训练结果我们发现,该问题哈密顿量的基态能量对应的边切割数趋近与5。
量子态展示
前面我们通过训练得到了量子线路中参数的最优值,下面,我们通过StateEvolution
类的final_state
来输出量子线路在最优参数时的量子态,其中ket
参数表示是否将最终量子态表示为右矢形式。
[10]:
pr = dict(zip(ansatz.para_name, net.weight.asnumpy()))
print(StateEvolution(init_state_circ+ansatz).final_state(pr, ket=True))
(0.017737679183483124-0.03180303797125816j)¦00000⟩
(-0.02683155983686447+0.0012889178469777107j)¦00001⟩
(0.011993971653282642+0.006973826792091131j)¦00010⟩
(-0.014608755707740784-0.003942559473216534j)¦00011⟩
(-0.02683155983686447+0.0012889178469777107j)¦00100⟩
(0.00725862430408597+0.10942266136407852j)¦00101⟩
(-0.014608755707740784-0.003942559473216534j)¦00110⟩
(0.008969870395958424-0.004171415697783232j)¦00111⟩
(0.00950924027711153-0.00026544960564933717j)¦01000⟩
(-0.37196943163871765-0.3156493902206421j)¦01001⟩
(-0.040885526686906815+0.037214867770671844j)¦01010⟩
(-0.37196943163871765-0.3156493902206421j)¦01011⟩
(-0.03160367161035538+0.009305878542363644j)¦01100⟩
(-0.040885526686906815+0.037214867770671844j)¦01101⟩
(-0.03160367161035538+0.009305878542363644j)¦01110⟩
(0.00950924027711153-0.00026544960564933717j)¦01111⟩
(0.00950924027711153-0.00026544960564933717j)¦10000⟩
(-0.03160367161035538+0.009305878542363644j)¦10001⟩
(-0.040885526686906815+0.037214867770671844j)¦10010⟩
(-0.03160367161035538+0.009305878542363644j)¦10011⟩
(-0.37196943163871765-0.3156493902206421j)¦10100⟩
(-0.040885526686906815+0.037214867770671844j)¦10101⟩
(-0.37196943163871765-0.3156493902206421j)¦10110⟩
(0.00950924027711153-0.00026544960564933717j)¦10111⟩
(0.008969870395958424-0.004171415697783232j)¦11000⟩
(-0.014608755707740784-0.003942559473216534j)¦11001⟩
(0.00725862430408597+0.10942266136407852j)¦11010⟩
(-0.02683155983686447+0.0012889178469777107j)¦11011⟩
(-0.014608755707740784-0.003942559473216534j)¦11100⟩
(0.011993971653282642+0.006973826792091131j)¦11101⟩
(-0.02683155983686447+0.0012889178469777107j)¦11110⟩
(0.017737679183483124-0.03180303797125816j)¦11111⟩
概率图
我们画出最终量子态在计算基矢下的概率分布
[11]:
def show_amp(state):
amp = np.abs(state)**2
n_qubits = int(np.log2(len(amp)))
labels = [bin(i)[2:].zfill(n_qubits) for i in range(len(amp))]
plt.bar(labels, amp)
plt.xticks(rotation=45)
plt.show()
state = StateEvolution(init_state_circ+ansatz).final_state(pr)
show_amp(state)
根据概率分布图我们发现,该Max-Cut问题具有四个简并解,每个解对应的概率大概为25%。
总结
这里我们通过量子近似优化算法来解决了Max-Cut问题,并得到了案例中的图对应的最大切割方案。
参考文献
[1] Edward Farhi, Jeffrey Goldstone, and Sam Gutmann. A Quantum Approximate Optimization Algorithm