本周汇报目前将分为以下三个部分

量子算法回顾

在周二下午的汇报之后发现对量子傅里叶变换、量子相位估计算法、Shor以及HHL算法理解的还不够透彻，下面列出了一些遗留的疑难问题的讨论和补充。汇报PPT点击此处下载

Grover算法——D变换

首先是Grover算法里的D变换矩阵，满足$D = 2P-I$，$P_{ij} = \frac{1}{N}$，故写成矩阵的形式为

$D = \begin{bmatrix}\frac{2}{N}-1 & \frac{2}{N} & \cdots & \frac{2}{N}\\ \frac{2}{N} & \frac{2}{N}-1 & \cdots & \frac{2}{N} \\ \vdots & \vdots & \vdots & \vdots \\ \frac{2}{N} & \frac{2}{N} & \cdots & \frac{2}{N}-1 \end{bmatrix}$

令$P\bar{v} = T$，其中$T_i= A = \frac{1}{N}\sum\bar{v_i}$，则有$D\bar{v}=(-I+2P)\bar{v}$。故$D\bar{v}$的第$i$个元素为$-\bar{v_i}+2A =A+(A-\bar{v_i})$，因此Grover把这种变换叫做the inversion about average。

设非目标矢量的振幅为$C$，显然，开始时$C$大约为$\frac{1}{\sqrt{N}}$，初始时A(各个矢量的均值)大约为$C$。而目标矢量的振幅是负数，大小也约为$C$，如下图(a)所示 [1]

对量子态进行D变换之后(量子态乘D算符)，非目标矢量完成了$A+(A-C)$的操作，由于初始 $A\approx C$，故D变换之后非目标矢量不发生变换。而目标矢量由于振幅为负，因此完成了$A+(A+C)$的操作，不仅振幅变为正数，而且增长了近$2C$，也就是约为$\frac{2}{\sqrt{N}}$，这就完成了一次$D$变换。目标矢量的振幅大小增加了，非目标矢量的振幅就相应地有所减少。在下一次循环中，A 变小。故目标矢量的增幅不断减少，直至某个极小值。如果一直循环，目标矢量的振幅非但没有增加，反而减小，故循环次数的确定对该算法意义重大。

通过上述过程，对Grover算法中一次G操作的第二次翻转有了更深入的认识，根据公式也能解释为什么需要第一次翻转，是为了制造出目标矢量的负振幅，为第二次对称的D变换服务，至于Grover算法的循环次数计算相信已经在PPT里解释的比较清楚了。

量子相位估计算法(Quantum Phase Estimation)

尽管这部分给出了算法实现的电路，但是相信结合代码能够理解的更清楚，尤其是能将复杂的量子傅里叶变换拆解，参考熊珏婵之前的工作，这部分将主要补充一下二进制相位分析[2]和代码说明[3]

在之前的考虑和讨论之中，忽略了$e^{2\pi} = 1$和$e^{i\pi}=-1$这两个事实(根据欧拉公式可以得出)，因此对公式的推导过程说明的极其模糊

图中$H$表示的是Hardmard门，$R_k$门表示的幺正矩阵形式为$R_{k} \equiv\left[\begin{array}{cc} 1 & 0 \\ 0 & e^{2 \pi i / 2^{k}} \end{array}\right]$。观察发现，对于第一个qubit我们用了一个$H$门和$3-1 = 2$个$R_k$门，对于第二个qubit我们用了一个$H$门和$2-1=1$个$R_k$门，因此总共用了$3+3-1+2-1=6$个门。

结合这两张图所展示的例子，能够把量子傅里叶变换的小数点运算过程解释清楚。这里假设输入的初态是$|110\rang$，对第一个比特位一次进行$H$门、$R_2$和$R_3$受控旋转。接下来便是关键性的一步“去一”，将$e^{2\pi} = 1$和$e^{i\pi}=-1$代入，能够把负号移动到指数上去，就得到了指数求和的形式$\frac{1}{2}+\frac{1}{4}+0$，这样的求和形式恰好能够由小数点的各个二进制位表示。

因此，将上述规律扩展到$N$的情况，就有了以下转换(PPT中有)：

同样我们采用“去一”的操作可以对上述图片中的（6）式进行进一步的化简(引入欧拉方程，将相位用二进制小数表示)，就能得到如下结果

$=\frac{1}{\sqrt{N}}\left(|0\rangle+e^{\left.2 \pi i \mid 0 . x_{n}\right]}|1\rangle\right) \otimes\left(|0\rangle+e^{2 \pi i\left[0 . x_{n-1} x_{n}\right]}|1\rangle\right) \otimes \cdots \otimes\left(|0\rangle+e^{2 \pi i\left[0 . x_{1} x_{2} \ldots x_{n}\right]}|1\rangle\right) (1)$

并且该式子还可以由如下矩阵形式(Hilbert空间中的算符)进行表示

$\frac{1}{\sqrt{N}}\left[\begin{array}{cccccc} 1 & 1 & 1 & 1 & \ldots & 1 \\ 1 & \omega^{1} & \omega^{1 \cdot 2} & \omega^{1 \cdot 3} & \ldots & \omega^{1 \cdot(N-1)} \\ 1 & \omega^{2 \cdot 1} & \omega^{2 \cdot 2} & \omega^{2 \cdot 3} & \ldots & \omega^{2 \cdot(N-1)} \\ 1 & \omega^{3 \cdot 1} & \omega^{3 \cdot 2} & \omega^{3 \cdot 3} & \ldots & \omega^{3 \cdot(N-1)} \\ \ldots & \ldots & \ldots & \ldots & \omega^{j k} & \ldots \\ 1 & \omega^{(N-1) \cdot 1} & \omega^{(N-1) \cdot 2} & \omega^{(N-1) \cdot 3} & \ldots & \omega^{(N-1) \cdot(N-1)} \end{array}\right]$

其中$w_n:=e^{\frac{2\pi i}{2^n}}$，$N := 2^n$。此外，类比3qubit的电路进行分析，可以发现电路中采用了$N$个$H$门，第一个qubit采用了$N-1$个$R_k$旋转门，第二个qubit采用了$N-2$个$R_k$旋转门，以此类推，一共采用了$N+N-1+N-2+\cdots+1=\frac{N(N+1)}{2}$个门，其电路图如下图所示

当然图中省略了一些SWAP运算(用来在计算后交换qubits的位置，应该第一行与最后一行交换，第二行与倒数第二行交换…全部换完后，和式（1）中的最终算式的顺序保持一致，即靠前的qubit存相位信息少的量子态)，最多有n/2个SWAP运算，每个SWAP可以用3个CNOT门实现，因此QFT所用的门的数量是$O(N^2)$。

经过对量子门的分析，那么就能更容易去理解代码里的逻辑，基于IBM的开源Python库Qiskit代码实现QFT

from qiskit import QuantumCircuit, ClassicalRegister, QuantumRegister
from qiskit import execute
from qiskit import BasicAer
import math
n = 4  # 第一组寄存器x
m = 1  # 第二组寄存器
# 创建n qubit的量子寄存器
c = ClassicalRegister(n+m, 'c')
# 创建量子电路
measure = QuantumCircuit(q, c)
# 构建量子傅里叶变换
for j in range(n):
    for k in range(j):
        # 受控门Rk，第k个qubit控制第j个qubit
        register.cu1(math.pi/float(2**(j-k)), q[j], q[k])
    # 第j个qubit的h门运算
    register.h(q[j])

# 定义qubit交换规则，即三个C-not门
def swap(qc, q, i, j):
    qc.cx(q[i], q[j])
    qc.cx(q[j], q[i])
    qc.cx(q[i], q[j])

# 控制n/2个swap门的操作，第一个与最后一个对应，以此类推
if n%2==1:
    for i in range(int((n-1)/2)):
        swap(register, q, i, n-i-1)
else:
    for i in range(int(n/2)):
        swap(register, q, i, n-i-1)
register.draw(output="mpl")

通过以上代码能够实现QFT电路，关键之处都给出了清晰的注释，而所谓逆傅里叶变换只不过改变受控旋转门$R_k$的方向，绘制出的电路图如下所示

而量子相位估计算法其实就是在上面的QFT基础上进行的。给定一个算符$U$，该算法能够有效估计$U|\psi\rang=e^{2\pi i \theta}|\psi\rang$中的相位$\theta$，其中$e^{2\pi i \theta}$是特征向量$\psi$对应的特征值。控制$U$门能够使qubit成比例的旋转$e^{2\pi i \theta}$，通过连续的控制$U$门重复此旋转适当的次数，直到将相位$\theta$编码到$0$和$2^t$之间，那么我们就能用基底进行表示，其电路图为

显然量子相位估计算法QPE就是在逆傅里叶变换的基础上加了$H$门和控制$U$门，那么只需要在上述傅里叶变换代码之前加上如下的结构即可

for j in range(n-1, -1, -1):
    register.h(q[j])
    # 当受控位q[j]为1时，对量子位使用幺正算子U
    register.ul(math.pi/float(2**(j)), q[j])

同时之前代码里的register.cu1(math.pi/float(2**(j-k)), q[j], q[k])要换成register.cu1(-math.pi/float(2**(j-k)), q[j], q[k])，逆傅里叶变换就这么搞定了，于是QPE算法的电路图就显现出来了。

在图像中进行了简单标注，以便与QPE算法电路图进行比对理解，一目了然。

HHL算法

再来首先回顾一下HHL算法，我们要求解的问题是$A|x\rang=|b\rang$，由于$A$是一个厄米算符，因此它能谱分解为如下形式

$A=\sum_{j=0}^{N-1} \lambda_{j}\left|u_{j}\right\rangle\left\langle u_{j}\right|, \quad \lambda_{j} \in \mathbb{R}$

其中$|u_j\rang$是第$j^{th}$个特征值对应的特征向量，那么有

$A^{-1}=\sum_{j=0}^{N-1} \lambda_{j}^{-1}\left|u_{j}\right\rangle\left\langle u_{j}\right|$

而待求解等式的右侧也可以这样的形式

$|b\rangle=\sum_{j=0}^{N-1} b_{j}\left|u_{j}\right\rangle, \quad b_{j} \in \mathbb{C}$

HHL算法从寄存器中读出目标状态就结束

$|x\rangle=A^{-1}|b\rangle=\sum_{j=0}^{N-1} \lambda_{j}^{-1} b_{j}\left|u_{j}\right\rangle$

HHL的算法电路图如下图所示

在HHL算法中，令QPE算法中的$U=e^{iAt}$，其中$A$是与要求解的系统相关联的矩阵。那么对于特征值为$e^{i\lambda_{j}t}$特征向量$|u_j\rang_{nb}$，QPT将输出 $\left|\tilde{\lambda}_{j}\right\rangle_{n l}\left|u_{j}\right\rangle_{n b}$。其中$\tilde{\lambda}_{j}$表示$n_l$-bit二进制近似等于$2^{n l} \frac{\lambda j t}{2 \pi}$。因此，每个$\lambda_j$也能用$n_l$个bit精确表示。

$\mathrm{QPE}\left(e^{i A t}, \sum_{j=0}^{N-1} b_{j}|0\rangle_{n l}\left|u_{j}\right\rangle_{n_{b}}\right)=\sum_{j=0}^{N-1} b_{j}\left|\lambda_{j}\right\rangle_{n l}\left|u_{j}\right\rangle_{n b}$

换句话说，通过第一阶段的相位估计，中间的寄存器q[1]会存储一系列A的特征值$j$（存储在基态$|j\rang$中），而底部的的寄存器存储的输入$|b\rang$会在A的特征空间中进行分解，表示为$|b\rang=\beta_j|u_j\rang$

而第二阶段则是通过附加的的量子比特q[0]将基态值的倒数$\frac{1}{\lambda}$按比例提取到了对应基态的概率幅上（控制R门），用公式表示为

$\sum_{j=0}^{N-1} b_{j}\left|\lambda_{j}\right\rangle_{n l}\left|u_{j}\right\rangle_{n b}\left(\sqrt{1-\frac{C^{2}}{\lambda_{j}^{2}}}|0\rangle+\frac{C}{\lambda_{j}}|1\rangle\right)$

从公式中可以看出，当对辅助寄存器观测值为1时，原来的寄存器由一系列$|j\rang$的叠加变为$f(j)*|0\rang$的叠加。经过上述操作以后，基态$|j\rang$中的值按照$f(j)$的比例被提取到了基态$|j\rang$的概率幅上(注：$f(j) = \frac{1}{j}$)。

可以看出相位估计使得寄存器和输入量子态纠缠到了一起，而由于控制R门的存在，使得对辅助寄存器q[0]的影响q[1]和q[2]的纠缠态$|\lambda_j\rang|u_j\rang$。而逆相位估计就是将q[1]和q[2]解纠缠，即将所有q[1]中的$\lambda_j$都还原为$|00\cdots\rang$，这样通过控制R门提取的比例最终影响的是q[2]寄存器中的$|u\rang$。 别忘了我们上面化简后的HHL求解表达式

$|x\rangle=A^{-1}|b\rangle=\sum_{j=0}^{N-1} \lambda_{j}^{-1} b_{j}\left|u_{j}\right\rangle$

显然，经过电路的处理，从输入寄存器q[2]就能得到我们上式$|x\rang$的结果，至此HHL算法结束。在整个过程中，**受控旋转操作只涉及了第一寄存器和第二寄存器，但最终却影响了第三寄存器的值，而这种影响就借助了量子计算特性中的纠缠性。**上述算法确实不好理解，在学习过程中我最深的理解还是来源于段小佳的HHL算法笔记和IBM的Qiskit量子计算，感谢这些作者。

MERA代码

本周也学些了TEBD for MPS的相关代码，在学习过程中又接触到一个新的张量网络收缩库ncon，能够用负数表示张量网络的Open boundry，从而完成对整个张量网络的收缩，官方的库地址点击这里，下面给出代码示例

from ncon import ncon
a = np.random.randn(3, 4)
b = np.random.randn(4, 5)
ab_ncon = ncon([a, b], ((-1, 1), (1, -2)))
ab_np = np.dot(a, b)
assert np.allclose(ab_ncon, ab_np)

这一小段代码就完成了张量a和张量b的收缩，可以看到张量a的第一个腿3和张量b的第二个腿5被保留下来，最终乘完是一个3*5的张量，相当于进行了矩阵乘积，同时该库的运算结果和numpy.dot以及einsum都是一致的，而einsum还需要指定索引，显然这个库的张量网络收缩使用起来更方便一些。

TEBD的代码还没有完全看完，但是梳理一下整个处理的思路就是

• Step1: 首先随机初始化MPS
• Step2: 采用Lanzcos方法求解左右环境张量
• Step3: 将MPS转换为正交形式
• Step4: 在A-B之间应用时间演化门，然后在B-A链接之间应用时间演化门。
• Step5: 重复Step2-Step4，在迭代中减少到较小的时间步长（提高准确率）

论证Covid-19生成型张量网络对X射线图像的分类中文完善

存在的问题及解决方案

（1）特征映射前后没有说明x所表示的含义

添加了x表示图片的像素，以及取值范围

（2）经过特征映射怎么表示到Hilbert空间

补充了Hilbert空间的定义，以及张量网络和Hilbert空间的关联性。同时也加了数学公式来描述feature map的作用，更清晰的表示了两个空间数据的转换。

（3）怎么得到的$\psi$

发现这个地方写的有些问题，其实输入的图片样本进行映射之后是一个Hilbert空间的向量$v$，而不是$\psi$，$\psi$是MPS量子态来的。

（4）生成型张量网络模型是什么样子的，怎么去训练的这个模型，梳理一套详细的流程图介绍过程

（5）对实验参数的介绍，参数怎么设置，怎么初始化，然后迭代次数，等等，举个例子来说明

把参数举例说明加到了第四部分结果的第二段，也就是拿实际实验参数作为例子，避免重复

Reference

1. 孙吉贵等 量子搜索算法
2. 熊珏婵IBMQiskit学习
3. Learn Quantum Computation using Qiskit
4. 矩阵低秩近似