{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 概览\n",
    "- 在这个案例中，我们将展示如何通过Paddle Quantum训练量子神经网络来求解量子系统的特征。\n",
    "\n",
    "- 首先，让我们通过下面几行代码引入必要的library和package。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "is_executing": false,
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "import numpy\n",
    "from paddle.complex import matmul, transpose\n",
    "from paddle import fluid\n",
    "from paddle_quantum.circuit import UAnsatz\n",
    "from numpy import array, kron"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 背景\n",
    "- 量子计算中在近期非常有前途的一个量子算法是变分量子特征求解器(VQE, variational quantum eigensolver (VQE)) [1-3].\n",
    "- VQE是量子化学在近期有噪量子设备（NISQ device）上的核心应用之一，其中一个功能比较强大的版本是SSVQE [4]，其核心是去求解一个物理系统的哈密顿量的基态和激发态的性质。数学上，可以理解为求解一个厄米矩阵(Hermitian matrix)的特征值及其对应的特征向量。该哈密顿量的特征值组成的集合我们称其为能谱。\n",
    "- 接下来我们将通过一个简单的例子学习如何通过训练量子神经网络解决这个问题，即求解出给定哈密顿量的能谱。\n",
    "\n",
    "## SSVQE分析物理系统的基态和激发态的能量\n",
    "- 对于具体需要分析的分子，我们需要输入其几何、电荷构型等多项信息。具体的，通过我们内置的量子化学工具包可以利用fermionic-to-qubit映射的技术来输出目标分子的量子比特哈密顿量表示。\n",
    "- 在这里，作为简单的入门案例，我们提供一个简单的2量子位哈密顿量作为例子。 "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "is_executing": false,
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "sigma_I = array([[1, 0], [0, 1]])\n",
    "sigma_X = array([[0, 1], [1, 0]])\n",
    "sigma_Y = array([[0, -1j], [1j, 0]])\n",
    "sigma_Z = array([[1, 0], [0, -1]])\n",
    "H = 0.4 * kron(sigma_Z, sigma_I) + 0.4 * kron(sigma_I, sigma_Z) + 0.2 * kron(sigma_X, sigma_X)\n",
    "hamiltonian = H.astype('complex64')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 搭建量子神经网络（QNN）\n",
    "- 在实现SSVQE的过程中，我们首先需要设计1个量子神经网络QNN（也可以理解为参数化量子电路）。这里，我们提供一个预设的2量子位量子电路。\n",
    "\n",
    "- 我们预设一些该参数化电路的参数，比如宽度为2量子位。\n",
    "\n",
    "- 初始化其中的变量参数，${\\bf{\\theta }}$代表我们量子神经网络中的参数组成的向量，一共有12个参数。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "is_executing": false,
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "N = 2 # 量子神经网络的宽度\n",
    "THETA_SIZE = 12 # 量子神经网络中参数的数量\n",
    "\n",
    "def U_theta(theta, N):\n",
    "    \"\"\"\n",
    "    U_theta\n",
    "    \"\"\"\n",
    "\n",
    "    cir = UAnsatz(N)\n",
    "    # ============== D1=2 ==============\n",
    "    cir.ry(theta[0], 2)\n",
    "    cir.rz(theta[1], 2)\n",
    "    cir.cnot([2, 1])\n",
    "    cir.ry(theta[2], 2)\n",
    "    cir.rz(theta[3], 2)\n",
    "    cir.cnot([2, 1])\n",
    "\n",
    "    # ============== D2=2 ==============\n",
    "    cir.ry(theta[4], 1)\n",
    "    cir.ry(theta[5], 2)\n",
    "    cir.rz(theta[6], 1)\n",
    "    cir.rz(theta[7], 2)\n",
    "    cir.cnot([1, 2])\n",
    "\n",
    "    cir.ry(theta[8], 1)\n",
    "    cir.ry(theta[9], 2)\n",
    "    cir.rz(theta[10], 1)\n",
    "    cir.rz(theta[11], 2)\n",
    "    cir.cnot([1, 2])\n",
    "\n",
    "    return cir.state"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 配置训练模型 - 损失函数\n",
    "- 现在我们已经有了数据和量子神经网络的架构，我们将进一步定义训练参数、模型和损失函数，具体的理论可以参考 [4].\n",
    "- 通过作用量子神经网络$U(\\theta)$在1组正交的初始态上，我们将得到输出态$\\left| {\\psi_k \\left( {\\bf{\\theta }} \\right)} \\right\\rangle $。\n",
    "- 进一步，在SSVQE模型中的损失函数一般由哈密顿量H与量子态$\\left| {\\psi_k \\left( {\\bf{\\theta }} \\right)} \\right\\rangle$的内积的加权求和给出。\n",
    "- 具体的损失函数定义为\n",
    "$$4\\left\\langle {\\psi_1 \\left( {\\bf{\\theta }} \\right)} \\right|H\\left| {\\psi_1 \\left( {\\bf{\\theta }} \\right)} \\right\\rangle + 3\\left\\langle {\\psi_2 \\left( {\\bf{\\theta }} \\right)} \\right|H\\left| {\\psi_2 \\left( {\\bf{\\theta }} \\right)} \\right\\rangle + 2\\left\\langle {\\psi_3 \\left( {\\bf{\\theta }} \\right)} \\right|H\\left| {\\psi_3 \\left( {\\bf{\\theta }} \\right)} \\right\\rangle + \\left\\langle {\\psi_4 \\left( {\\bf{\\theta }} \\right)} \\right|H\\left| {\\psi_4 \\left( {\\bf{\\theta }} \\right)} \\right\\rangle.$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "is_executing": false,
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "SEED = 1\n",
    "\n",
    "class Net(fluid.dygraph.Layer):\n",
    "    \"\"\"\n",
    "    Construct the model net\n",
    "    \"\"\"\n",
    "\n",
    "    def __init__(self, shape, param_attr=fluid.initializer.Uniform(low=0.0, high=2 * numpy.pi, seed=SEED),\n",
    "                 dtype='float32'):\n",
    "        super(Net, self).__init__()\n",
    "\n",
    "        self.theta = self.create_parameter(shape=shape, attr=param_attr, dtype=dtype, is_bias=False)\n",
    "\n",
    "    def forward(self, H, N):\n",
    "        \"\"\"\n",
    "        Args:\n",
    "            input_state: The initial state with default |0..>\n",
    "            H: The target Hamiltonian\n",
    "        Returns:\n",
    "            The loss.\n",
    "        \"\"\"\n",
    "        out_state = U_theta(self.theta, N)\n",
    "\n",
    "        loss_struct = matmul(matmul(\n",
    "            transpose(fluid.framework.ComplexVariable(out_state.real, -out_state.imag), perm=[1, 0]), H),\n",
    "            out_state).real\n",
    "\n",
    "        loss_components = [\n",
    "            loss_struct[0][0],\n",
    "            loss_struct[1][1],\n",
    "            loss_struct[2][2],\n",
    "            loss_struct[3][3]\n",
    "        ]\n",
    "\n",
    "        loss = 4 * loss_components[0] + 3 * loss_components[1] + 2 * loss_components[2] + 1 * loss_components[3]\n",
    "        return loss, loss_components\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 配置训练模型 - 模型参数\n",
    "在进行量子神经网络的训练之前，我们还需要进行一些训练（超）参数的设置，例如学习速率与迭代次数。\n",
    "- 设定学习速率（learning rate）为0.3。\n",
    "- 设定迭代次数为50次。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "is_executing": false,
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "ITR = 50 # 迭代次数\n",
    "\n",
    "LR = 0.3 # 学习速率\n",
    " "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 进行训练\n",
    "\n",
    "- 当训练模型的各项参数都设置完成后，我们将数据转化为Paddle动态图中的变量，进而进行量子神经网络的训练。\n",
    "- 过程中我们用的是Adam Optimizer，也可以调用Paddle中提供的其他优化器。\n",
    "- 我们将训练过程中的每一轮loss可以打印出来。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "is_executing": false,
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "with fluid.dygraph.guard():\n",
    "    # Harmiltonian preparing\n",
    "    hamiltonian = fluid.dygraph.to_variable(hamiltonian)\n",
    "\n",
    "    # net\n",
    "    net = Net(shape=[THETA_SIZE])\n",
    "\n",
    "    # optimizer\n",
    "    opt = fluid.optimizer.AdagradOptimizer(learning_rate=LR, parameter_list=net.parameters())\n",
    "\n",
    "    # gradient descent loop\n",
    "    for itr in range(1, ITR + 1):\n",
    "        loss, loss_components = net(hamiltonian, N)\n",
    "\n",
    "        loss.backward()\n",
    "        opt.minimize(loss)\n",
    "        net.clear_gradients()\n",
    "\n",
    "        print('iter:', itr, 'loss:', '%.4f' % loss.numpy()[0])\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 测试效果\n",
    "我们现在已经完成了量子神经网络的训练，我们将通过与理论值的对比来测试效果。\n",
    "- 理论值由numpy中的工具来求解哈密顿量的各个特征值；\n",
    "- 我们将训练QNN得到的各个能级的能量和理想情况下的理论值进行比对。\n",
    "- 可以看到，SSVQE训练输出的值与理想值高度接近。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "pycharm": {
     "is_executing": false,
     "name": "#%%\n"
    }
   },
   "outputs": [],
   "source": [
    "print('The estimated ground state energy is: ', loss_components[0].numpy())\n",
    "print('The theoretical ground state energy: ', numpy.linalg.eigh(H)[0][0])\n",
    "\n",
    "print('The estimated 1st excited state energy is: ', loss_components[1].numpy())\n",
    "print('The theoretical 1st excited state energy: ', numpy.linalg.eigh(H)[0][1])\n",
    "\n",
    "print('The estimated 2nd excited state energy is: ', loss_components[2].numpy())\n",
    "print('The theoretical 2nd excited state energy: ', numpy.linalg.eigh(H)[0][2])\n",
    "\n",
    "print('The estimated 3rd excited state energy is: ', loss_components[3].numpy())\n",
    "print('The theoretical 3rd excited state energy: ', numpy.linalg.eigh(H)[0][3])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "## 参考文献\n",
    "[1] A. Peruzzo et al., “A variational eigenvalue solver on a photonic quantum processor,” Nat. Commun., vol. 5, no. 1, p. 4213, Dec. 2014.\n",
    "\n",
    "[2] S. McArdle, S. Endo, A. Aspuru-Guzik, S. C. Benjamin, and X. Yuan, “Quantum computational chemistry,” Rev. Mod. Phys., vol. 92, no. 1, p. 015003, Mar. 2020.\n",
    "\n",
    "[3] Y. Cao et al., “Quantum chemistry in the age of quantum computing,” Chem. Rev., vol. 119, no. 19, pp. 10856–10915, 2019.\n",
    "\n",
    "[4] K. M. Nakanishi, K. Mitarai, and K. Fujii, “Subspace-search variational quantum eigensolver for excited states,” Phys. Rev. Res., vol. 1, no. 3, p. 033062, Oct. 2019.\n",
    "\n",
    "\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.7.7"
  },
  "pycharm": {
   "stem_cell": {
    "cell_type": "raw",
    "metadata": {
     "collapsed": false
    },
    "source": []
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}