tutorials/quantum_simulation/ClassicalShadow_Application_CN.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "4f560a55",
   "metadata": {},
   "source": [
    "# 基于经典影子的量子态性质估计"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a894cbba",
   "metadata": {},
   "source": [
    "<em> Copyright (c) 2021 Institute for Quantum Computing, Baidu Inc. All Rights Reserved. </em>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e4cdaa03",
   "metadata": {},
   "source": [
    "## 概览"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3b3ea48a",
   "metadata": {},
   "source": [
    "在[未知量子态的经典影子](./ClassicalShadow_Intro_CN.ipynb)教程中，我们介绍了经典影子的一些理论知识并展示了如何构建一个量子态 $\\rho$ 的经典影子。根据经典影子的理论推导可以看出，它十分适合用来估计量子态线性的性质，目前其基本应用场景有如下四个：量子态的保真度估计、纠缠验证、局部可观测量期望的估计、全局可观测量期望的估计 [1]。其中可观测量期望值的估计广泛地出现在当前的量子算法中，例如用于估计复杂分子哈密顿量（Hamiltonian）基态能量的[变分量子本征求解器](./VQE_CN.ipynb)（variational quantum eigensolver, VQE）。接下来我们将重点讨论基于经典影子的可观测量期望值估计算法，并展示如何使用量桨中的 shadow 功能来进行可观测量期望值的估计。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "db8b954f",
   "metadata": {},
   "source": [
    "## 可观测量期望值估计"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e2869ea5",
   "metadata": {},
   "source": [
    "### 问题描述"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a6adea50",
   "metadata": {},
   "source": [
    "在量子化学领域，核心任务之一是求解一个量子尺度上封闭物理系统的哈密顿量 $\\hat{H}$ 的基态能量及其对应的基态，主要的实现方法是通过在量子设备上准备一个参数化的试探波函数 $|\\Psi(\\theta)\\rangle$，然后结合经典机器学习中的优化算法（例如梯度下降法）去不断地调整、优化参数 $\\theta$ 使得期望值 $\\langle\\Psi(\\theta)|\\hat{H}| \\Psi(\\theta)\\rangle$ 最小化。这套方案的基本原理是基于 Rayleigh-Ritz 变分原理。\n",
    "\n",
    "$$\n",
    "E_{0}=\\min _{\\theta}\\langle\\Psi(\\theta)|\\hat{H}| \\Psi(\\theta)\\rangle, \\tag{1}\n",
    "$$\n",
    "\n",
    "其中 $E_{0}$ 表示该系统的基态能量。从数值分析的角度来看，该问题可以被理解为求解一个离散化哈密顿量 $\\hat{H}$ （埃尔米特矩阵）的最小本征值 $\\lambda_{\\min }$ 和其对应的本征向量 $\\left|\\Psi_{0}\\right\\rangle$ 。经典影子发挥作用的场景就是每一次优化中计算 $\\langle\\Psi(\\theta)|\\hat{H}| \\Psi(\\theta)\\rangle = \\operatorname{tr}(\\hat{H}\\rho )$ 的部分（其中 $\\rho = | \\Psi(\\theta)\\rangle\\langle\\Psi(\\theta)|$）。\n",
    "\n",
    "于是问题转化为：对于一个 $n$ 个量子比特的量子态 $\\rho$ 和一个可写成一组泡利算子$\\{I,X,Y,Z\\}^{\\otimes n}$线性组合的可观测量（哈密顿量）$\\hat{H}$，\n",
    "\n",
    "$$\n",
    "\\hat{H}=\\sum_{Q \\in\\{I, X, Y, Z\\} ^{\\otimes n}} \\alpha_{Q} Q, \\quad \\text{where} \\quad \\alpha_{Q} \\in  \\mathbb{R}, \\tag{2}\n",
    "$$\n",
    "\n",
    "如何用经典影子来估计可观测量期望值 $\\operatorname{tr}(\\hat{H}\\rho )$ ？\n",
    "\n",
    "最直观的方法是将哈密顿量的每一项分别作为测量基，对量子态 $\\rho$ 进行相应的泡利测量，并对每一项的测量进行一定次数的重复，再统计测量结果得到估计值。在这里，我们称此方法为逐项测量的方法。\n",
    "\n",
    "读者可以看到当哈密顿量 $\\hat{H}$ 项数较少，$n$ 也较小时，我们可以通过逐项测量的方法来得到 $\\operatorname{tr}(\\hat{H}\\rho )$，但当 $\\hat{H}$ 项数增多，且 $n$ 较大时，逐项测量的方法所需要的代价将大大增加。而将要介绍的基于经典影子的方法，可以用更少的代价得到同等精度的 $\\operatorname{tr}(\\hat{H}\\rho )$ 的估计。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ff8af318",
   "metadata": {},
   "source": [
    "### 基于经典影子的改进算法"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a3db3328",
   "metadata": {},
   "source": [
    "在经典影子的构建中，一个关键的步骤是从固定集合中均匀随机采样酉变换，[未知量子态的经典影子](./ClassicalShadow_Intro_CN.ipynb)教程展示了所选集合为 Clifford 群时的情况。当所选的集合是作用在单量子比特上的 Clifford 群时，构建时的采样与测量步骤就相当于对量子态做泡利测量。量桨中提供了使用随机泡利测量的经典影子算法（Classical shadows using random Pauli measurements，CS）。简单来说，在 CS 算法中，我们重复地为每个量子比特均匀随机选择一个泡利基来测量量子态 $\\rho$，根据测量结果估计可观测量期望值，具体的原理，读者可参考  [1-2] 学习。进一步地，同样是选取泡利测量基，当选取的方式不再是均匀随机时，基于经典影子的改进算法被先后提出 [2-3]。量桨中也提供了相关的算法功能：局部偏置的经典影子算法（Locally-biased classical shadows，LBCS）[2]，自适应泡利影子算法（Adaptive Pauli shadows，APS）[3]。感兴趣的读者可以参考 [1-3] 来详细学习这些算法。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1b772de7",
   "metadata": {},
   "source": [
    "## Paddle Quantum 代码实现"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e094e752",
   "metadata": {},
   "source": [
    "在量桨中，我们基于经典影子提供了 shadow 功能，主要包含两个函数，支持用户使用上述基于经典影子的三种算法来估计可观测量的期望值，以及获取未知量子态的经典影子数据。下面我们将展示如何基于量桨中的 shadow 功能来实现氢分子（H$_{2}$）和氢化锂（LiH）基态能量估计。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "efc769b4",
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# 导入需要的包\n",
    "import numpy as np\n",
    "import paddle\n",
    "import paddle_quantum as pq\n",
    "from paddle_quantum import Hamiltonian\n",
    "from paddle_quantum.VQE.chemistrysub import H2_generator\n",
    "from paddle_quantum.loss import ExpecVal\n",
    "from paddle_quantum.state import zero_state, to_state\n",
    "from paddle_quantum.qinfo import shadow_trace"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a58e8871",
   "metadata": {},
   "source": [
    "### 估计氢分子（H$_{2}$）基态能量"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ebda5aeb",
   "metadata": {},
   "source": [
    "导入拥有 4 个量子比特的氢分子（H$_{2}$）的哈密顿量（用户具体可以参考[变分量子本征求解器](./VQE_CN.ipynb) 教程，来获得氢分子（H$_{2}$）哈密顿量）。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "70b8fdef",
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "H2 hamiltonian =  -0.04207897647782277 I\n",
      "0.17771287465139946 Z0\n",
      "0.1777128746513994 Z1\n",
      "-0.2427428051314046 Z2\n",
      "-0.24274280513140462 Z3\n",
      "0.17059738328801055 Z0, Z1\n",
      "0.04475014401535163 Y0, X1, X2, Y3\n",
      "-0.04475014401535163 Y0, Y1, X2, X3\n",
      "-0.04475014401535163 X0, X1, Y2, Y3\n",
      "0.04475014401535163 X0, Y1, Y2, X3\n",
      "0.12293305056183797 Z0, Z2\n",
      "0.1676831945771896 Z0, Z3\n",
      "0.1676831945771896 Z1, Z2\n",
      "0.12293305056183797 Z1, Z3\n",
      "0.1762764080431959 Z2, Z3\n"
     ]
    }
   ],
   "source": [
    "# 导入量桨中预处理好的氢分子哈密顿量\n",
    "H2_pauli_str, H2_qubit = H2_generator()\n",
    "# 根据 H2_pauli_str 创建哈密顿量类\n",
    "H2_hamiltonian = Hamiltonian(H2_pauli_str)\n",
    "print('H2 hamiltonian = ', H2_hamiltonian)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0d6a6532",
   "metadata": {},
   "source": [
    "为了展示如何用基于经典影子的算法来估计基态能量，我们首先通过量桨中的 VQE 来估计氢分子（H$_{2}$）的基态，并得到基态所对应的量子电路。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "40036ae4",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "iter: 20 loss: -1.0562\n",
      "iter: 20 estimated ground state energy: -1.0562 Ha\n",
      "iter: 40 loss: -1.1200\n",
      "iter: 40 estimated ground state energy: -1.1200 Ha\n",
      "iter: 60 loss: -1.1344\n",
      "iter: 60 estimated ground state energy: -1.1344 Ha\n",
      "iter: 80 loss: -1.1360\n",
      "iter: 80 estimated ground state energy: -1.1360 Ha\n",
      "\n",
      "训练后的电路：\n",
      "--Ry(6.286)----*--------------x----Ry(-0.00)----*--------------x----Ry(3.159)--\n",
      "               |              |                 |              |               \n",
      "--Ry(6.277)----x----*---------|----Ry(-0.21)----x----*---------|----Ry(6.269)--\n",
      "                    |         |                      |         |               \n",
      "--Ry(3.137)---------x----*----|----Ry(6.286)---------x----*----|----Ry(6.272)--\n",
      "                         |    |                           |    |               \n",
      "--Ry(6.277)--------------x----*----Ry(6.290)--------------x----*----Ry(3.145)--\n",
      "                                                                               \n"
     ]
    }
   ],
   "source": [
    "ITR = 80  # 设置训练的总迭代次数\n",
    "LR = 0.4   # 设置学习速率\n",
    "D = 2      # 设置量子神经网络中重复计算模块的深度 Depth\n",
    "N = H2_hamiltonian.n_qubits # 获取量子比特数\n",
    "\n",
    "paddle.seed(42)\n",
    "# 初始化待训练的量子电路\n",
    "cir_H_2 = pq.ansatz.Circuit(N)\n",
    "cir_H_2.real_entangled_layer(depth=D)\n",
    "cir_H_2.ry()\n",
    "\n",
    "# 定义初始态\n",
    "init_state = zero_state(N)\n",
    "\n",
    "# 定义损失函数\n",
    "loss_func = ExpecVal(H2_hamiltonian)\n",
    "\n",
    "# 一般来说，我们利用Adam优化器来获得相对好的收敛，\n",
    "opt = paddle.optimizer.Adam(learning_rate=LR, parameters=cir_H_2.parameters())\n",
    "\n",
    "# 记录优化结果\n",
    "summary_iter, summary_loss = [], []\n",
    "\n",
    "# 优化循环\n",
    "for itr in range(1, ITR + 1):\n",
    "\n",
    "    # 前向传播计算损失函数\n",
    "    state = cir_H_2(init_state)\n",
    "    loss = loss_func(state)\n",
    "\n",
    "    # 反向传播极小化损失函数\n",
    "    loss.backward()\n",
    "    opt.minimize(loss)\n",
    "    opt.clear_grad()\n",
    "\n",
    "    # 更新优化结果\n",
    "    summary_loss.append(loss.numpy())\n",
    "    summary_iter.append(itr)\n",
    "\n",
    "    # 打印结果\n",
    "    if itr % 20 == 0:\n",
    "        print(\"iter:\", itr, \"loss:\", \"%.4f\" % loss.numpy())\n",
    "        print(\"iter:\", itr, \"estimated ground state energy:\", \"%.4f Ha\" \n",
    "                                            % loss.numpy())\n",
    "    if itr == ITR:\n",
    "        print(\"\\n训练后的电路：\") \n",
    "        print(cir_H_2)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a6c9d23f",
   "metadata": {},
   "source": [
    "#### shadow 功能介绍"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c9891733",
   "metadata": {},
   "source": [
    "此时，我们获得了生成氢分子（H$_{2}$）的基态所对应的量子电路。在此电路上，我们可以直接运行 `shadow_trace` 函数，来获得使用经典影子算法估计得到的基态能量。\n",
    "在 `shadow_trace` 函数中，我们的输入为需要估计的哈密顿量、采样次数、以及所选择使用的采样算法。用户可以通过指定 `method` 参数来选择想要使用的采样算法。其中 CS 的适用范围更广，速度最快，但是其估计精度可能稍差；LBCS 的精度更高，但在哈密顿量的项数较高时运行偏慢；APS 的精度也会更高，但是量子比特数目较大时运行偏慢。 "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "3b31e57d",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "H2 ground state energy =  -1.1359991\n",
      "H2 ground state energy CS=  -1.1532653284115573\n",
      "H2 ground state energy LBCS=  -1.1219509441210997\n",
      "H2 ground state energy APS=  -1.18267941594884\n"
     ]
    }
   ],
   "source": [
    "# 估计出的基态对应的实际能量值\n",
    "H2_energy = loss_func(cir_H_2()).numpy()[0]\n",
    "# 采样次数\n",
    "sample = 1500 \n",
    "# 分别采用三种算法估计可观测量期望值，即基态能量\n",
    "H2_energy_CS = shadow_trace(cir_H_2(), H2_hamiltonian, sample, method=\"CS\")\n",
    "H2_energy_LBCS = shadow_trace(cir_H_2(), H2_hamiltonian, sample, method=\"LBCS\")\n",
    "H2_energy_APS = shadow_trace(cir_H_2(), H2_hamiltonian, sample, method=\"APS\")\n",
    "\n",
    "print('H2 ground state energy = ', H2_energy)\n",
    "print('H2 ground state energy CS= ', H2_energy_CS)\n",
    "print('H2 ground state energy LBCS= ', H2_energy_LBCS)\n",
    "print('H2 ground state energy APS= ', H2_energy_APS)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1cb856da",
   "metadata": {},
   "source": [
    "现在让我们采用对哈密顿量逐项测量的传统方法来估计氢分子（H$_{2}$）的基态能量，"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "28af249f",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "3ae0242c",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "H2 ground state energy traditional =  -1.13771533304257\n"
     ]
    }
   ],
   "source": [
    "# 使用哈密顿量逐项测量方法来估计基态能量\n",
    "H2_energy_traditional = cir_H_2().expec_val(H2_hamiltonian, shots=100)\n",
    "print('H2 ground state energy traditional = ', H2_energy_traditional)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b2e1ad02",
   "metadata": {},
   "source": [
    "我们可以看到在 1500 次采样下，三种算法估计的基态能量与 VQE 算法估计出的基态的实际能量已经十分接近，而逐项测量的方法是针对哈密顿量中的每一项各进行 100 次测量，氢分子（H_$2$）哈密顿量有 15 项，相当于测量了 1500 次，得到的结果与 VQE 算法估计出的基态的实际能量的差距也较小。在这种小规模的情况，基于经典影子的算法与逐项测量的方法比并没有体现出明显的优势，但在大规模量子系统场景下，该类算法需要的采样次数仅仅是关于哈密顿量项数的常数级别的增长，逐项测量或已有的一些方法则需要关于哈密顿量项数呈多项式甚至指数级别增长的采样次数，来达到同样精度 [1]。事实上，[2] 中指出对于 CS 算法以及 LBCS 算法，我们得到的估计的平均误差 $\\epsilon$，方差 $\\operatorname{var}(\\nu)$ 以及采样次数 $S$ 有如下关系关系：\n",
    "\n",
    "$$\n",
    "S = O(\\epsilon^{-2} \\operatorname{var}(\\nu) ),  \\tag{3}\n",
    "$$\n",
    "\n",
    "其中方差 $\\operatorname{var}(\\nu)$ 与采样次数相互独立，且方差与哈密顿量的项数有关。于是根据我们想要达到的精度（平均误差），可以计算出实际所需要的采样次数。同样，根据采样次数，也可以定义我们的平均误差：\n",
    "\n",
    "$$\n",
    "\\epsilon = \\sqrt{\\frac{\\operatorname{var}}{S}}.  \\tag{4}\n",
    "$$\n",
    "\n",
    "可以看到我们上述的实验中，CS 算法和 LBCS 算法得到的误差，均在文献 [3] 理论估计的精度之内。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5e0e174d",
   "metadata": {},
   "source": [
    "同时，量桨提供了与经典影子有关的采样函数 `shadow_sample` ，支持对未知量子态的事先采样，方便读者探索经典影子的其他应用，具体使用方式如下："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "fc6b465c",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[('zyzx', '1001'), ('xyzx', '0100'), ('yxyx', '1100'), ('yxzx', '0001'), ('zxxx', '1010'), ('xxyx', '0000'), ('zyxy', '1011'), ('xxzz', '0100'), ('xxzz', '1100'), ('yxzx', '0100')]\n"
     ]
    }
   ],
   "source": [
    "from paddle_quantum.shadow import shadow_sample\n",
    "\n",
    "# 在态矢量模式下运行电路，获得当前电路输出态\n",
    "H2_rho = to_state(cir_H_2().data)\n",
    "# 获得经典影子的数据，输出成 list 的形式\n",
    "H2_sample_data_CS = shadow_sample(H2_rho, H2_qubit, sample_shots=10, mode='density_matrix', \n",
    "                                  hamiltonian=H2_hamiltonian, method='CS')\n",
    "print(H2_sample_data_CS)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e74ffd18",
   "metadata": {},
   "source": [
    "### 估计氢化锂（LiH）基态能量"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e4480cb6",
   "metadata": {},
   "source": [
    "接下来考虑氢化锂（LiH）的基态能量，这里我们直接读取预先计算好的文件，以此来生成具有 12 个量子比特的氢化锂（LiH）的泡利形式的分子哈密顿量。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "6916c8a4",
   "metadata": {},
   "outputs": [],
   "source": [
    "with open('./LiH_hamiltonian.txt', 'r') as lih_file:\n",
    "    unprocessed_pauli_str = lih_file.read()\n",
    "    LiH_pauli_str = [term.split(maxsplit=1) for term in unprocessed_pauli_str.split('\\n')]\n",
    "    LiH_pauli_str = [[float(term[0]), term[1]] for term in LiH_pauli_str]\n",
    "    LiH_hamiltonian = Hamiltonian(LiH_pauli_str)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "61d8fd59",
   "metadata": {},
   "source": [
    "接下来，我们同样可以通过运行 VQE 电路估计其分子哈密顿量的基态。这里，由于该分子哈密顿量较大，VQE 所需要的训练时间较长，我们直接提供了已经训练后的 VQE 电路的参数，用户可以通过它直接获得估计的 LiH 的基态并在其上测试基于经典影子的方法。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "fe326409",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "预训练 VQE 得到的基态能量为：-7.7720 \n"
     ]
    }
   ],
   "source": [
    "# 读取事先训练好的参数\n",
    "pretrained_parameters = paddle.load('LiH_VQE_parameters.pdtensor')\n",
    "N = LiH_hamiltonian.n_qubits\n",
    "\n",
    "# 按照量子比特数量/网络宽度初始化量子神经网络\n",
    "cir_LiH = pq.ansatz.Circuit(N)\n",
    "\n",
    "# 内置的 {R_y + CNOT} 电路模板\n",
    "cir_LiH.real_entangled_layer(depth=D)\n",
    "cir_LiH.ry()\n",
    "cir_LiH.update_param(pretrained_parameters)\n",
    "\n",
    "# 计算哈密顿量为 LiH 时的期望值\n",
    "expec_val_LiH = cir_LiH().expec_val(LiH_hamiltonian, shots=0)\n",
    "\n",
    "print('预训练 VQE 得到的基态能量为：%.4f ' % expec_val_LiH)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8f8ff501",
   "metadata": {},
   "source": [
    "得到了估计的氢化锂分子（LiH）基态所对应的电路后，我们直接使用 `shadow_trace` 函数进行随机测量即可。同时，由于该分子哈密顿量有 631 项，为了保证两类方法测量次数一致，我们规定函数 `shadow_trace` 的 `sample = 1262`，函数 `expecval` 的 `shots = 2`。\n",
    "\n",
    "又因为 LiH 基态的量子比特数为 12，所以对 LiH 的基态做不同的泡利测量时，共有 $3^{12}$ 种可能的测量组合，那么仅仅进行 1262 次采样从而得到估值，具有随机性。于是，我们分别运行 20 次上述四种方法，取这 20 个样本数据的均值作为做为各个算法的估计值，并计算样本方差，对算法进行简单的比较。（运行下述代码块需要的时间可能较长，通常需要 20 分钟左右）"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "8f3684ce",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "LiH ground state energy =  -7.77198\n",
      "ave LiH ground state energy CS =  -7.850611450035954\n",
      "ave LiH ground state energy LBCS =  -7.793319461201525\n",
      "ave LiH ground state energy APS =  -7.774484195247901\n",
      "ave LiH ground state energy traditional =  -7.622710303123105\n",
      "time =  1003.3228180408478\n"
     ]
    }
   ],
   "source": [
    "import time\n",
    "\n",
    "begin = time.time()\n",
    "estimator_CS = []\n",
    "estimator_LBCS = []\n",
    "estimator_APS = []\n",
    "estimator_traditional = []\n",
    "\n",
    "# 估计出的基态对应的实际能量值\n",
    "LiH_energy = cir_LiH().expec_val(LiH_hamiltonian, shots=0)\n",
    "\n",
    "# 运行算法次数\n",
    "n = 20 \n",
    "\n",
    "for i in range(n):\n",
    "    \n",
    "    LiH_energy_CS = shadow_trace(cir_LiH(), LiH_hamiltonian, 1262, method=\"CS\")\n",
    "    LiH_energy_LBCS = shadow_trace(cir_LiH(), LiH_hamiltonian, 1262, method=\"LBCS\")\n",
    "    LiH_energy_APS = shadow_trace(cir_LiH(), LiH_hamiltonian, 1262, method=\"APS\")\n",
    "    LiH_energy_traditional = cir_LiH().expec_val(LiH_hamiltonian, shots=2)\n",
    "\n",
    "    estimator_CS.append(LiH_energy_CS) \n",
    "    estimator_LBCS.append(LiH_energy_LBCS) \n",
    "    estimator_APS.append(LiH_energy_APS) \n",
    "    estimator_traditional.append(LiH_energy_traditional) \n",
    "\n",
    "ave_LiH_energy_CS = np.mean(estimator_CS)\n",
    "ave_LiH_energy_LBCS = np.mean(estimator_LBCS)\n",
    "ave_LiH_energy_APS = np.mean(estimator_APS)\n",
    "ave_LiH_energy_traditional = np.mean(estimator_traditional)\n",
    "end = time.time() \n",
    "\n",
    "print(\"LiH ground state energy = \", LiH_energy)\n",
    "print(\"ave LiH ground state energy CS = \", ave_LiH_energy_CS)\n",
    "print(\"ave LiH ground state energy LBCS = \", ave_LiH_energy_LBCS)\n",
    "print(\"ave LiH ground state energy APS = \", ave_LiH_energy_APS)\n",
    "print('ave LiH ground state energy traditional = ', ave_LiH_energy_traditional)\n",
    "print('time = ', end-begin)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "89469273",
   "metadata": {},
   "source": [
    "从结果来看，基于经典影子算法得到的均值比逐项测量的更接近 VQE 算法估计出的 LiH 基态的实际能量，且算法的误差均在文献 [3] 理论估计的精度之内。那么各算法的样本方差又是怎样的呢？"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "742f4f15",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "LiH variance CS =  0.20268237699326175\n",
      "LiH variance LBCS =  0.016956271213630907\n",
      "LiH variance APS =  0.0028097294191431854\n",
      "LiH variance traditional =  0.11731067295981745\n"
     ]
    }
   ],
   "source": [
    "# 计算样本方差\n",
    "variance_CS = []\n",
    "variance_LBCS = []\n",
    "variance_APS = []\n",
    "variance_traditional = []\n",
    "\n",
    "for i in range(n):\n",
    "    variance_CS.append((estimator_CS[i] - ave_LiH_energy_CS) ** 2)\n",
    "    variance_LBCS.append((estimator_LBCS[i] - ave_LiH_energy_LBCS) ** 2)\n",
    "    variance_APS.append((estimator_APS[i] - ave_LiH_energy_APS) ** 2)\n",
    "    variance_traditional.append((estimator_traditional[i] - ave_LiH_energy_traditional) ** 2)\n",
    "\n",
    "var_CS = sum(variance_CS)/(n-1)\n",
    "var_LBCS = sum(variance_LBCS)/(n-1)\n",
    "var_APS = sum(variance_APS)/(n-1)\n",
    "var_traditional = sum(variance_traditional)/(n-1)\n",
    "\n",
    "print('LiH variance CS = ', var_CS)\n",
    "print('LiH variance LBCS = ', var_LBCS)\n",
    "print('LiH variance APS = ', var_APS)\n",
    "print('LiH variance traditional = ', var_traditional)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d40ef806",
   "metadata": {},
   "source": [
    "可以看到，APS 算法的样本方差是最低的，其次是 LBCS 算法，接着是逐项测量的方法，最后是 CS 算法。据此，我们可以发现哈密顿量的项数规模增大后，基于经典影子的算法与逐项测量的方法相比，在同等代价下精度更高，且更加稳定。其中 APS 算法是最稳定的。\n",
    "\n",
    "值得一提的是，对于经典影子算法来说 12 个量子比特的场景仍不能较好地展现出其与现有一些算法相比的巨大优势。在具有更多量子比特的大规模系统中，其在算法代价上的优势才能更好地被展现 [1]。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fa9e6a26",
   "metadata": {},
   "source": [
    "## 总结"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "91a60770",
   "metadata": {},
   "source": [
    "本教程讨论了如何用基于经典影子的改进算法来得到可观测量期望值的估计，并展示了如何使用量桨中的 shadow 功能。可以看到，基于经典影子的改进算法可以对可观测量期望得到很好的估计。相比逐项测量的方法，在采样次数一致的情况下，它的估计值更精确，且算法更稳定。在大规模量子系统场景下，经典影子方法在一些问题中所需要的采样次数与系统大小无关，故对于系统大小来说仅仅是常数级别的增长 [1]。所以它在 NISQ（noisy intermediate-scale quantum）时代所能发挥的作用将继续地被不断挖掘。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5fcdcad3",
   "metadata": {},
   "source": [
    "_______\n",
    "\n",
    "## 参考文献"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "947b52db",
   "metadata": {},
   "source": [
    "[1] Huang, Hsin-yuan, R. Kueng and J. Preskill. “Predicting many properties of a quantum system from very few measurements.” [Nature Physics (2020): 1-8.](https://www.nature.com/articles/s41567-020-0932-7?proof=t)\n",
    "\n",
    "[2] Hadfield, Charles, et al. \"Measurements of quantum hamiltonians with locally-biased classical shadows.\" [arXiv preprint arXiv:2006.15788 (2020).](https://arxiv.org/abs/2006.15788)\n",
    "\n",
    "[3] Hadfield, Charles. \"Adaptive Pauli Shadows for Energy Estimation.\" [arXiv preprint arXiv:2105.12207 (2021).](https://arxiv.org/abs/2105.12207)"
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "4261e4eef114648d37e4a637967bd8d2966507f48b194e5e336ba3366b740269"
  },
  "kernelspec": {
   "display_name": "Python 3.8.13 ('pq-ns-icode')",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.16"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}