other/BasicQuantumAlgorithms.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Quantum Fourier Transform and Quantum Phase Esitmation with Tequila"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this tutorial, we introduce the quantum fourier transform (QFT) and the quamtum phase estimation (QPE) algorithm, derive the circuit, and implement it using Tequila.\n",
    "\n",
    "## Contents\n",
    "1. [Introduction](#introduction)\n",
    "2. [Quantum Fourier Transform](#qft)     \n",
    "    2.1 [Mathematical Description of the Quantum Fourier Transform](#math)         \n",
    "    2.2 [Expansion of the Quantum Fourier Transform](#expansion)  \n",
    "    2.3 [The Derivation of the Circuit Implementation](#deriv)  \n",
    "    2.4 [Tequila Implementation of the Quantum Fourier Transform](#impl)  \n",
    "    2.5 [Simulating the QFT Circuit](#sim)\n",
    "3. [Quantum Phase Estimation](#qpe)     \n",
    "    3.1 [Mathematical Description of the Quantum Phase Estimation](#math2)  \n",
    "    3.2 [Quantum Phase Estimation Circuit](#sim2)\n",
    "4. [References](#references)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 1. Introduction <a id='introduction'></a>\n",
    "\n",
    "The Fourier transform occurs in many different places throughout classical computing, ranging from signal processing applications to data compression to complexity theory. The quantum Fourier transform (QFT) is the quantum counterpart of the discrete Fourier transform over the amplitudes of a wavefunction. It is part of many quantum algorithms, such as Shor's factoring algorithm and quantum phase estimation. As an example of where quantum Fourier transform is used, we will also take a look at the quantum phase estimation algorithm, and how it utilizes the QFT. Quantum phase estimation is also one of the most important subroutines in quantum computation. It serves as a central building block for many quantum algorithms, such as Shor's algorithm and Quantum Counting."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 2. Quantum Fourier Transform <a id='qft'></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2.1 Mathematical Description of the Quantum Fourier Transform <a id='math'></a>\n",
    "The Quantum Fourier Transform acts on a quantum state \n",
    "$$ \n",
    "\\lvert x \\rangle = \\sum_{i = 0}^{N-1} x_i \\lvert i \\rangle \n",
    "$$\n",
    "to give the quantum state \n",
    "$$\n",
    "\\lvert y \\rangle = \\sum_{i = 0}^{N-1} y_i \\lvert i \\rangle\n",
    "$$\n",
    "such that \n",
    "$$\n",
    "y_i = \\frac{1}{\\sqrt{N}} \\sum_{j = 0}^{N-1}x_je^{2\\pi i \\frac{ji}{N}} = \\frac{1}{\\sqrt{N}} \\sum_{j = 0}^{N-1}x_jw^{ji}_{N}\n",
    "$$\n",
    "We note that this mapping is equivalent to applying to each basis state $\\lvert i \\rangle$ the following mapping:\n",
    "$$\n",
    "\\lvert i \\rangle \\mapsto \\frac{1}{\\sqrt{N}}\\sum_{y = 0}^{N-1}w^{iy}_N\\lvert y \\rangle\n",
    "$$\n",
    "To see this equivalance, expand and combine liked terms, then it should be clear. Moreover, this mapping correponds to the following Unitary Matrix $U_{QFT}$, where:\n",
    "$$\n",
    "U_{QFT} = \\frac{1}{\\sqrt{N}}\\sum_{i=0}^{N-1}\\sum_{y=0}^{N-1} w^{iy}_N \\lvert y \\rangle \\langle i \\rvert \n",
    "$$\n",
    "To see this that this matrix represents the transformation as defined above and is Unitary, we note that:\n",
    "\\begin{align*}\n",
    "U_{QFT}\\lvert j \\rangle &= \\Big(\\frac{1}{\\sqrt{N}}\\sum_{i=0}^{N-1}\\sum_{y=0}^{N-1} w^{iy}_N \\lvert y \\rangle \\langle i \\rvert\\Big) \\lvert j \\rangle\\\\\n",
    "&= \\frac{1}{\\sqrt{N}}\\sum_{y=0}^{N-1} \\lvert y \\rangle \\sum_{i=0}^{N-1} w^{iy}_N  \\langle i \\mid j \\rangle\\\\\n",
    "&= \\frac{1}{\\sqrt{N}}\\sum_{y=0}^{N-1} \\lvert y \\rangle \\sum_{i=0}^{N-1} w^{iy}_N  \\delta_{ij} && \\text{Since $(\\lvert 0 \\rangle, \\dots, \\lvert {N-1} \\rangle)$ is an orthonormal basis}\\\\\n",
    "&= \\frac{1}{\\sqrt{N}}\\sum_{y=0}^{N-1} w^{jy}_N \\lvert y \\rangle \\\\\n",
    "\\end{align*}\n",
    "\n",
    "as wanted. To see that $U_{QFT}$ is Unitary, we note that:\n",
    "$$\n",
    "U_{QFT}^{\\dagger} = \\frac{1}{\\sqrt{N}}\\sum_{i'=0}^{N-1}\\sum_{y'=0}^{N-1} w^{(-1)i'y'}_N \\lvert i' \\rangle \\langle y' \\rvert \n",
    "$$\n",
    "then we get that:\n",
    "\\begin{align*}\n",
    "U_{QFT}^{\\dagger}U_{QFT} &= \\frac{1}{N} \\sum_{i, y, i', y'} e^{2\\pi i (i'y' - iy)/N}\\rvert y \\rangle \\langle y' \\rvert \\delta_{ii'}\\\\\n",
    "&= \\frac{1}{N} \\sum_{y, i, y'} e^{2\\pi i (y' - y)i/N}\\rvert y \\rangle \\langle y' \\lvert \\\\\n",
    "&= \\frac{1}{N} \\sum_{y, i, y'} \\rvert y \\rangle \\langle y' \\lvert \\delta_{yy'}\\\\\n",
    "&= \\sum_y\\lvert y \\rangle \\langle y \\rvert = I\n",
    "\\end{align*}\n",
    "as desired. Hence, $U_{QFT}$ is unitary. \n",
    "\n",
    "Now we attempt to develop a Quantum Circuit in Python that will perform the $U_{QFT}$ transformation."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2.2 Expansion of the Quantum Fourier Transform$^1$ <a id='expansion'></a>\n",
    "To see why this circuit implements the Quantum Fourier Transform, we first consider the expansion of $U_{QFT}$ as it operates on the state $\\lvert x \\rangle = \\lvert x_1\\dots x_n \\rangle$, where each $x_i$ represents a qubit, thus $N = 2^n$:\n",
    "$$\n",
    "\\begin{aligned}\n",
    "QFT_N\\vert x \\rangle & = \\frac{1}{\\sqrt{N}} \\sum_{y=0}^{N-1}\\omega_N^{xy} \\vert y \\rangle \n",
    "\\\\\n",
    "& = \\frac{1}{\\sqrt{N}} \\sum_{y=0}^{N-1} e^{2 \\pi i xy / 2^n} \\vert y \\rangle ~\\text{since}\\: \\omega_N^{xy} = e^{2\\pi i \\frac{xy}{N}} \\:\\text{and}\\: N = 2^n \n",
    "\\\\\n",
    "& = \\frac{1}{\\sqrt{N}} \\sum_{y=0}^{N-1} e^{2 \\pi i \\left(\\sum_{k=1}^n y_k/2^k\\right) x} \\vert y_1 \\ldots y_n \\rangle \\:\\text{rewriting in fractional binary notation}\\: y = y_1\\ldots y_n, y/2^n = \\sum_{k=1}^n y_k/2^k \n",
    "\\\\\n",
    "& = \\frac{1}{\\sqrt{N}} \\sum_{y=0}^{N-1} \\prod_{k=1}^n e^{2 \\pi i x y_k/2^k } \\vert y_1 \\ldots y_n \\rangle \\:\\text{after expanding the exponential of a sum to a product of exponentials} \n",
    "\\\\\n",
    "& = \\frac{1}{\\sqrt{N}} \\bigotimes_{k=1}^n  \\left(\\vert0\\rangle + e^{2 \\pi i x /2^k } \\vert1\\rangle \\right) \\:\\text{after rearranging the sum and products, and expanding} \n",
    "\\sum_{y=0}^{N-1} = \\sum_{y_1=0}^{1}\\sum_{y_2=0}^{1}\\ldots\\sum_{y_n=0}^{1} \n",
    "\\\\\n",
    "& = \\frac{1}{\\sqrt{N}}\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2}x} \\vert1\\rangle\\right) \n",
    "\\otimes\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2^2}x} \\vert1\\rangle\\right) \n",
    "\\otimes  \n",
    "\\ldots\n",
    "\\otimes\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2^{n-1}}x} \\vert1\\rangle\\right) \n",
    "\\otimes\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2^n}x} \\vert1\\rangle\\right) \n",
    "\\end{aligned}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2.3 The Derivation of the Circuit Implementation <a id='deriv'></a>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.3.1 Hadamard and Controlled Rotational Gates\n",
    "In this section, we will look at how the circuit implemented above equates the the $QFT_N$ matrix mathematically. Looking at the code, we see that there are only two gates that are used: The single-qubit Hadamard gate, $H$, and the two-qubit controlled rotation $CROT_k$ gate. The Hadamard gate is the following:\n",
    "$$\n",
    "H = \\frac{1}{\\sqrt{2}}\\begin{bmatrix} 1 & 1 \\\\ 1 & -1 \\end{bmatrix}\n",
    "$$\n",
    "and controlled rotation $CROT_k$ given in block-diagonal form as: \n",
    "\n",
    "$$CROT_k =\n",
    "\\begin{bmatrix} 1 & 0 & 0 & 0 \\\\ 0 & 1 & 0 & 0 \\\\ 0 & 0 & 1 & 0 \\\\ 0 & 0 & 0 & \\exp\\left(\\frac{2\\pi i}{2^k}\\right)\\end{bmatrix}\n",
    "$$\n",
    "Now we note that applying the $H$ operator to a single-qubit state $\\vert\\psi\\rangle = \\alpha \\vert 0 \\rangle + \\beta \\vert 1 \\rangle$, we get a new state:\n",
    "$$\\frac{1}{\\sqrt{2}}(\\alpha + \\beta) \\vert 0 \\rangle + \\frac{1}{\\sqrt{2}}(\\alpha - \\beta)  \\vert 1 \\rangle $$\n",
    "letting $\\vert\\psi\\rangle = \\vert0\\rangle$, we get $\\alpha = 1, \\beta = 0$, therefore:\n",
    "$$\n",
    "H\\vert\\psi\\rangle = H\\vert0\\rangle = \\frac{1}{\\sqrt{2}}(\\vert0\\rangle + \\vert1\\rangle)\n",
    "$$\n",
    "and similarly when  $\\vert\\psi\\rangle = \\vert1\\rangle$, we get:\n",
    "$$\n",
    "H\\vert\\psi\\rangle = H\\vert1\\rangle = \\frac{1}{\\sqrt{2}}(\\vert0\\rangle - \\vert1\\rangle)\n",
    "$$\n",
    "Hence, in general, we can say:\n",
    "$$\n",
    "H\\vert \\psi \\rangle = \\frac{1}{\\sqrt{2}}\\left(\\vert0\\rangle + \\exp\\left(\\frac{2\\pi i}{2}\\psi\\right)\\vert1\\rangle\\right)\n",
    "$$\n",
    "Now we look at how the $CROT_k$ gate acts on the $2-$qubit system $\\lvert\\psi\\rangle = \\alpha\\lvert00\\rangle + \\beta\\lvert01\\rangle + \\gamma\\lvert10\\rangle + \\delta\\lvert11\\rangle:$\n",
    "$$\n",
    "CROT_k\\lvert\\psi\\rangle = \\lvert\\psi\\rangle = \\alpha\\lvert00\\rangle + \\beta\\lvert01\\rangle + \\gamma\\lvert10\\rangle + \\exp\\left(\\frac{2\\pi i}{2^k}\\right)\\delta\\lvert11\\rangle\n",
    "$$\n",
    "Then it is easy to see that:\n",
    "$$CROT_k\\vert 0\\psi_j\\rangle = \\vert 0\\psi_j\\rangle$$\n",
    "\n",
    "\n",
    "and\n",
    "\n",
    "\n",
    "$$CROT_k\\vert 1\\psi_j\\rangle = \\exp\\left( \\frac{2\\pi i}{2^k}\\psi_j \\right)\\vert 1\\psi_j\\rangle$$\n",
    "\n",
    "Let us apply these two gates to get an equivalent tranform as the Unitary Matrix $U_{QFT}$ as defined above.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.3.2 Equivalance of Circuit and $U_{QFT}$ Transformation Proof$^1$\n",
    "The circuit operates as follows. We start with an n-qubit input state $\\vert x_1x_2\\ldots x_n\\rangle$.\n",
    "\n",
    "\n",
    "After the first Hadamard gate on qubit 1, the state is transformed from the input state to \n",
    "\n",
    "$$\n",
    "H_1\\vert x_1x_2\\ldots x_n\\rangle = \n",
    "\\frac{1}{\\sqrt{2}}\n",
    "\\left[\\vert0\\rangle + \\exp\\left(\\frac{2\\pi i}{2}x_1\\right)\\vert1\\rangle\\right]\n",
    "\\otimes\n",
    "\\vert x_2x_3\\ldots x_n\\rangle\n",
    "$$\n",
    "\n",
    "<li> After the $UROT_2$ gate on qubit 1 controlled by qubit 2, the state is transformed to\n",
    "\n",
    "$$\n",
    "\\frac{1}{\\sqrt{2}}\n",
    "\\left[\\vert0\\rangle + \\exp\\left(\\frac{2\\pi i}{2^2}x_2 + \\frac{2\\pi i}{2}x_1\\right)\\vert1\\rangle\\right]\n",
    "\\otimes\n",
    "\\vert x_2x_3\\ldots x_n\\rangle\n",
    "$$\n",
    "\n",
    "<li> After the application of the last $UROT_n$ gate on qubit 1 controlled by qubit $n$, the state becomes\n",
    "\n",
    "$$\n",
    "\\frac{1}{\\sqrt{2}}\n",
    "\\left[\\vert0\\rangle + \n",
    "\\exp\\left(\n",
    "\\frac{2\\pi i}{2^n}x_n + \n",
    "\\frac{2\\pi i}{2^{n-1}}x_{n-1} + \n",
    "\\ldots + \n",
    "\\frac{2\\pi i}{2^2}x_2 + \n",
    "\\frac{2\\pi i}{2}x_1\n",
    "\\right)\n",
    "\\vert1\\rangle\\right]\n",
    "\\otimes\n",
    "\\vert x_2x_3\\ldots x_n\\rangle\n",
    "$$\n",
    "\n",
    "Letting\n",
    "\n",
    "$$\n",
    "x = 2^{n-1}x_1 + 2^{n-2}x_2 + \\ldots + 2^1x_{n-1} + 2^0x_n\n",
    "$$\n",
    "\n",
    "we can write the above state as \n",
    "\n",
    "$$\n",
    "\\frac{1}{\\sqrt{2}}\n",
    "\\left[\\vert0\\rangle + \n",
    "\\exp\\left(\n",
    "\\frac{2\\pi i}{2^n}x \n",
    "\\right)\n",
    "\\vert1\\rangle\\right]\n",
    "\\otimes\n",
    "\\vert x_2x_3\\ldots x_n\\rangle\n",
    "$$\n",
    "\n",
    "<li> After the application of a similar sequence of gates for qubits $2\\ldots n$, we find the final state to be:\n",
    "\n",
    "$$\n",
    "\\frac{1}{\\sqrt{2}}\n",
    "\\left[\\vert0\\rangle + \n",
    "\\exp\\left(\n",
    "\\frac{2\\pi i}{2^n}x \n",
    "\\right)\n",
    "\\vert1\\rangle\\right]\n",
    "\\otimes\n",
    "\\frac{1}{\\sqrt{2}}\n",
    "\\left[\\vert0\\rangle + \n",
    "\\exp\\left(\n",
    "\\frac{2\\pi i}{2^{n-1}}x \n",
    "\\right)\n",
    "\\vert1\\rangle\\right]\n",
    "\\otimes\n",
    "\\ldots\n",
    "\\otimes\n",
    "\\frac{1}{\\sqrt{2}}\n",
    "\\left[\\vert0\\rangle + \n",
    "\\exp\\left(\n",
    "\\frac{2\\pi i}{2^{2}}x \n",
    "\\right)\n",
    "\\vert1\\rangle\\right]\n",
    "\\otimes\n",
    "\\frac{1}{\\sqrt{2}}\n",
    "\\left[\\vert0\\rangle + \n",
    "\\exp\\left(\n",
    "\\frac{2\\pi i}{2^{1}}x \n",
    "\\right)\n",
    "\\vert1\\rangle\\right]\n",
    "$$\n",
    "\n",
    "which is exactly the QFT of the input state as derived above but in reversed order, which is why the swaps in the circuit is applied.\n",
    "    \n",
    "## 2.4 Tequila Implementation of the Quantum Fourier Transform <a id='impl'></a>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "     ┌─────────┐                                 ┌─────────┐                 »\n",
      "q_0: ┤ RZ(π/8) ├─────────────■────────────────■──┤ RZ(π/4) ├─────────────────»\n",
      "     ├─────────┤             │                │  └─────────┘                 »\n",
      "q_1: ┤ RZ(π/4) ├─────────────┼────────────────┼───────────────■──────────────»\n",
      "     └──┬───┬──┘┌─────────┐┌─┴─┐┌──────────┐┌─┴─┐┌─────────┐┌─┴─┐┌──────────┐»\n",
      "q_2: ───┤ H ├───┤ RZ(π/8) ├┤ X ├┤ RZ(-π/8) ├┤ X ├┤ RZ(π/4) ├┤ X ├┤ RZ(-π/4) ├»\n",
      "        └───┘   └─────────┘└───┘└──────────┘└───┘└─────────┘└───┘└──────────┘»\n",
      "c: 3/════════════════════════════════════════════════════════════════════════»\n",
      "                                                                             »\n",
      "«                                                ┌───┐   \n",
      "«q_0: ───────────────────────■────────────────■──┤ H ├─X─\n",
      "«          ┌───┐┌─────────┐┌─┴─┐┌──────────┐┌─┴─┐└───┘ │ \n",
      "«q_1: ──■──┤ H ├┤ RZ(π/4) ├┤ X ├┤ RZ(-π/4) ├┤ X ├──────┼─\n",
      "«     ┌─┴─┐└───┘└─────────┘└───┘└──────────┘└───┘      │ \n",
      "«q_2: ┤ X ├────────────────────────────────────────────X─\n",
      "«     └───┘                                              \n",
      "«c: 3/═══════════════════════════════════════════════════\n",
      "«                                                        \n"
     ]
    }
   ],
   "source": [
    "import tequila as tq\n",
    "from numpy import pi\n",
    "\n",
    "def qft_rotations(n):\n",
    "    \"\"\"Given the number of qubit registers <n>, create the appropriate Controlled Rotation gates\"\"\"\n",
    "    circuit = tq.gates.H(target=0)\n",
    "    for i in range(1, n + 1):\n",
    "        for qubit in reversed(range(i)):\n",
    "            circuit = tq.gates.Phase(phi = pi/(2 ** (i - qubit)), target = i, control = qubit) + circuit\n",
    "        circuit = tq.gates.H(target=i) + circuit\n",
    "    return circuit\n",
    "\n",
    "def qft_swap(circuit, n):\n",
    "    \"\"\"Given the number of qubit registers <n> and the circuit <circuit>, create the appropriate Controlled Rotation gates,\n",
    "    perform appropriate register swaps\"\"\"\n",
    "    for qubit in range(n//2):\n",
    "        circuit += tq.gates.SWAP(qubit, n-qubit - 1)\n",
    "    return circuit\n",
    "\n",
    "def qft(n):\n",
    "    \"\"\"Given the number of qubit registers <n>, return the QFT circuit.\"\"\"\n",
    "    return qft_swap(qft_rotations(n - 1), n)\n",
    "\n",
    "circuit = qft(3)\n",
    "tq.draw(circuit, backend='qiskit')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2.5 Simulating the QFT Circuit <a id='sim'></a>\n",
    "\n",
    "In order to check the correctness of our circuit, we will take the Quantum Fourier state of the bits representing the number $4$ and then take the inverse of our circuit and simulate the circuit on real hardware. If we see that the highest probability is placed on $3$ (ie. $011$), then we can be sure our circuit is correct. Let us first define our qft_inverse function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "           ┌───┐    ┌──────────┐                     ┌──────────┐            »\n",
      "q_0: ─X────┤ H ├────┤ RZ(-π/4) ├──■───────────────■──┤ RZ(-π/8) ├────────────»\n",
      "      │ ┌──┴───┴───┐└──────────┘┌─┴─┐┌─────────┐┌─┴─┐└──┬───┬───┘┌──────────┐»\n",
      "q_1: ─┼─┤ RZ(-π/4) ├────────────┤ X ├┤ RZ(π/4) ├┤ X ├───┤ H ├────┤ RZ(-π/4) ├»\n",
      "      │ ├──────────┤            └───┘└─────────┘└───┘   └───┘    └──────────┘»\n",
      "q_2: ─X─┤ RZ(-π/4) ├─────────────────────────────────────────────────────────»\n",
      "        └──────────┘                                                         »\n",
      "c: 3/════════════════════════════════════════════════════════════════════════»\n",
      "                                                                             »\n",
      "«                                                                \n",
      "«q_0: ───────────────────────────────────■───────────────■───────\n",
      "«                                        │               │       \n",
      "«q_1: ──■───────────────■────────────────┼───────────────┼───────\n",
      "«     ┌─┴─┐┌─────────┐┌─┴─┐┌──────────┐┌─┴─┐┌─────────┐┌─┴─┐┌───┐\n",
      "«q_2: ┤ X ├┤ RZ(π/4) ├┤ X ├┤ RZ(-π/8) ├┤ X ├┤ RZ(π/8) ├┤ X ├┤ H ├\n",
      "«     └───┘└─────────┘└───┘└──────────┘└───┘└─────────┘└───┘└───┘\n",
      "«c: 3/═══════════════════════════════════════════════════════════\n",
      "«                                                                \n"
     ]
    }
   ],
   "source": [
    "tq.draw(qft(3).dagger())\n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we get the Quantum Fourier state for $011$, which rotates the most significant qubit by $\\frac{3}{4\\pi}$, the middle qubit by $\\frac{3}{2\\pi}$ and the least significant qubit by $3\\pi$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "     ┌───┐┌──────────┐\n",
      "q_0: ┤ H ├┤ RZ(3π/4) ├\n",
      "     ├───┤├──────────┤\n",
      "q_1: ┤ H ├┤ RZ(3π/2) ├\n",
      "     ├───┤└┬────────┬┘\n",
      "q_2: ┤ H ├─┤ RZ(3π) ├─\n",
      "     └───┘ └────────┘ \n",
      "c: 3/═════════════════\n",
      "                      \n"
     ]
    }
   ],
   "source": [
    "n = 3\n",
    "qt_state = tq.gates.H(target=0) + tq.gates.H(target=1) + tq.gates.H(target=2) + tq.gates.Phase(target=0, phi=n* pi/4) \\\n",
    "+ tq.gates.Phase(target=1, phi=n * pi/2) + tq.gates.Phase(target=2, phi=n * pi)\n",
    "tq.draw(qt_state)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Simulating the $QFT^\\dagger$ circuit should get us back to the original state of $\\lvert 011 \\rangle$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "+100.0000|011> \n"
     ]
    }
   ],
   "source": [
    "print(tq.simulate(qt_state + qft(3).dagger(), samples=100, read_out_qubits=[2,1,0]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "After a $100$ sample simulations, we got the state  $\\lvert 011 \\rangle$ every time, as desired."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 3. Quantum Phase Estimation <a id='qpe'></a>\n",
    "\n",
    "## 3. 1 Mathematical Description of the Quantum Phase Estimation <a id='math2'></a>\n",
    "\n",
    "This circuit is aimed at estimating the phase of a unitary operator $U$. It estimates $\\theta$ in $U\\vert\\psi \\rangle =e^{\\boldsymbol{2\\pi i} \\theta }|\\psi \\rangle$, where $|\\psi\\rangle$ is an eigenvector and $e^{\\boldsymbol{2\\pi i}\\theta}$ is the corresponding eigenvalue. The circuit operates in the following steps:\n",
    "\n",
    "\n",
    "\n",
    " \n",
    "\n",
    "<li> First we apply a $n$-bit Hadamard gate operation $H^{\\otimes n}$ on the $n$-counting register: \n",
    "\n",
    "\n",
    "\n",
    "$$ \\psi_1 = {\\frac {1}{2^{\\frac {n}{2}}}}\\left(|0\\rangle +|1\\rangle \\right)^{\\otimes n} \\lvert \\psi \\rangle$$\n",
    "\n",
    "\n",
    "\n",
    "<li> $\\textbf{Controlled Unitary Operations}$: Like the controlled phase rotation gate, the controlled unitary $C-U$ applies the unitary operator $U$ on the target register if and only if its corresponding control bit is $|1\\rangle$. Since $U$ is a unitary operatory with eigenvector $|\\psi\\rangle$ such that $U|\\psi \\rangle =e^{\\boldsymbol{2\\pi i} \\theta }|\\psi \\rangle$, we get the following:\n",
    "\n",
    "\n",
    "$$U^{2^{j}}|\\psi \\rangle =U^{2^{j}-1}U|\\psi \\rangle =U^{2^{j}-1}e^{2\\pi i\\theta }|\\psi \\rangle =\\cdots =e^{2\\pi i2^{j}\\theta }|\\psi \\rangle$$\n",
    "\n",
    "\n",
    "\n",
    "Applying all the $n$ controlled operations $C − U^{2^j}$ with $0\\leq j\\leq n-1$, and using the relation $|0\\rangle \\otimes |\\psi \\rangle +|1\\rangle \\otimes e^{2\\pi i\\theta }|\\psi \\rangle =\\left(|0\\rangle +e^{2\\pi i\\theta }|1\\rangle \\right)\\otimes |\\psi \\rangle$, carefully noting the relationship between the powers of $e^{\\boldsymbol{2\\pi i} \\theta}$ and the binary represention of each qubit, we get the following:\n",
    "\n",
    "\\begin{aligned}\n",
    "\\psi_{2} =\\frac {1}{2^{\\frac {n}{2}}} \\left(|0\\rangle+{e^{\\boldsymbol{2\\pi i} \\theta 2^{n-1}}}|1\\rangle \\right) \\otimes \\cdots \\otimes \\left(|0\\rangle+{e^{\\boldsymbol{2\\pi i} \\theta 2^{1}}}\\vert1\\rangle \\right) \\otimes \\left(|0\\rangle+{e^{\\boldsymbol{2\\pi i} \\theta 2^{0}}}\\vert1\\rangle \\right) \\otimes |\\psi\\rangle  = \\frac{1}{2^{\\frac {n}{2}}}\\sum _{k=0}^{2^{n}-1}e^{\\boldsymbol{2\\pi i} \\theta k}|k\\rangle \\otimes \\vert\\psi\\rangle\n",
    "\\end{aligned}\n",
    "where $k$ denotes the integer representation of n-bit binary numbers. \n",
    "\n",
    "<li> $\\textbf{Inverse Fourier Transform}$: Recall from above that QFT maps an n-qubit input state $\\vert x\\rangle$ into an output as\n",
    "\n",
    "$$\n",
    "QFT\\vert x \\rangle = \\frac{1}{2^\\frac{n}{2}}\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2}x} \\vert1\\rangle\\right) \n",
    "\\otimes\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2^2}x} \\vert1\\rangle\\right) \n",
    "\\otimes  \n",
    "\\ldots\n",
    "\\otimes\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2^{n-1}}x} \\vert1\\rangle\\right) \n",
    "\\otimes\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2^n}x} \\vert1\\rangle\\right) \n",
    "$$\n",
    "\n",
    "Replacing $x$ by $2^n\\theta$ in the above expression gives exactly the expression derived in step 2 above. Therefore, to recover the state $\\vert2^n\\theta\\rangle$, apply an inverse Fourier transform on the n-bit counting registers. Doing so get us:\n",
    "\n",
    "$$\n",
    "QFT^\\dagger QFT\\vert 2^n \\theta \\rangle =  \\vert 2^n \\theta \\rangle = QFT^\\dagger(\\frac{1}{2^\\frac{n}{2}}\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2}x} \\vert1\\rangle\\right) \n",
    "\\otimes\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2^2}x} \\vert1\\rangle\\right) \n",
    "\\otimes  \n",
    "\\ldots\n",
    "\\otimes\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2^{n-1}}x} \\vert1\\rangle\\right) \n",
    "\\otimes\n",
    "\\left(\\vert0\\rangle + e^{\\frac{2\\pi i}{2^n}x} \\vert1\\rangle\\right)) \n",
    "$$\n",
    "\n",
    "<li> $\\textbf{Measurement}$: Hence, cleary the above set of actions when measured should peak at $x = 2^n\\theta$. For the case when $2^n\\theta$ is an integer, measuring in the computational basis gives the phase in the n-bit counting registers with high probability. For the case when $2^n\\theta$ is not an integer, it can be shown that the above expression still peaks near $x = 2^n\\theta$ with probability better than $4/\\pi^2 \\approx 40\\%$.\n",
    "\n",
    "## 3.2 Quantum Phase Estimation Circuit <a id='sim2'></a>\n",
    "\n",
    "Let us consider an example where we apply the $T$-gate:\n",
    "$$ T|1\\rangle = \n",
    "\\begin{bmatrix}\n",
    "1 & 0\\\\\n",
    "0 & e^\\frac{i\\pi}{4}\\\\ \n",
    "\\end{bmatrix}\n",
    "\\begin{bmatrix}\n",
    "0\\\\\n",
    "1\\\\ \n",
    "\\end{bmatrix}\n",
    "= e^\\frac{i\\pi}{4}|1\\rangle $$\n",
    "\n",
    "thus the expected results of the QPE should give us $\\theta = \\frac{1}{8}$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let us implement and apply QPE:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "     ┌───┐┌─────────┐┌───┐┌──────────┐┌───┐                                 »\n",
      "q_0: ┤ H ├┤ RZ(π/8) ├┤ X ├┤ RZ(-π/8) ├┤ X ├─────────────────────────────────»\n",
      "     ├───┤├─────────┤└─┬─┘└──────────┘└─┬─┘           ┌───┐┌──────────┐┌───┐»\n",
      "q_1: ┤ H ├┤ RZ(π/8) ├──┼────────────────┼─────────────┤ X ├┤ RZ(-π/8) ├┤ X ├»\n",
      "     ├───┤├─────────┤  │                │             └─┬─┘└──────────┘└─┬─┘»\n",
      "q_2: ┤ H ├┤ RZ(π/8) ├──┼────────────────┼───────────────┼────────────────┼──»\n",
      "     ├───┤├─────────┤  │                │  ┌─────────┐  │                │  »\n",
      "q_3: ┤ X ├┤ RZ(π/8) ├──■────────────────■──┤ RZ(π/8) ├──■────────────────■──»\n",
      "     └───┘└─────────┘                      └─────────┘                      »\n",
      "c: 4/═══════════════════════════════════════════════════════════════════════»\n",
      "                                                                            »\n",
      "«                                                                        »\n",
      "«q_0: ───────────────────────────────────────────────────────────────────»\n",
      "«     ┌─────────┐┌───┐┌──────────┐┌───┐┌──────────┐                      »\n",
      "«q_1: ┤ RZ(π/8) ├┤ X ├┤ RZ(-π/8) ├┤ X ├┤ RZ(-π/4) ├──────────────────────»\n",
      "«     └─────────┘└─┬─┘└──────────┘└─┬─┘└──────────┘┌───┐┌──────────┐┌───┐»\n",
      "«q_2: ─────────────┼────────────────┼──────────────┤ X ├┤ RZ(-π/8) ├┤ X ├»\n",
      "«     ┌─────────┐  │                │  ┌─────────┐ └─┬─┘└──────────┘└─┬─┘»\n",
      "«q_3: ┤ RZ(π/8) ├──■────────────────■──┤ RZ(π/8) ├───■────────────────■──»\n",
      "«     └─────────┘                      └─────────┘                       »\n",
      "«c: 4/═══════════════════════════════════════════════════════════════════»\n",
      "«                                                                        »\n",
      "«                                                                       »\n",
      "«q_0: ──────────────────────────────────────────────────────────────────»\n",
      "«                                                                       »\n",
      "«q_1: ──────────────────────────────────────────────────────────────────»\n",
      "«     ┌─────────┐┌───┐┌──────────┐┌───┐┌─────────┐┌───┐┌──────────┐┌───┐»\n",
      "«q_2: ┤ RZ(π/8) ├┤ X ├┤ RZ(-π/8) ├┤ X ├┤ RZ(π/8) ├┤ X ├┤ RZ(-π/8) ├┤ X ├»\n",
      "«     ├─────────┤└─┬─┘└──────────┘└─┬─┘├─────────┤└─┬─┘└──────────┘└─┬─┘»\n",
      "«q_3: ┤ RZ(π/8) ├──■────────────────■──┤ RZ(π/8) ├──■────────────────■──»\n",
      "«     └─────────┘                      └─────────┘                      »\n",
      "«c: 4/══════════════════════════════════════════════════════════════════»\n",
      "«                                                                       »\n",
      "«                                            ┌───┐    ┌──────────┐     »\n",
      "«q_0: ──────────────────────────────────X────┤ H ├────┤ RZ(-π/4) ├──■──»\n",
      "«                                       │    └───┘    └──────────┘┌─┴─┐»\n",
      "«q_1: ──────────────────────────────────┼─────────────────────────┤ X ├»\n",
      "«     ┌─────────┐┌───┐┌──────────┐┌───┐ │ ┌──────────┐            └───┘»\n",
      "«q_2: ┤ RZ(π/8) ├┤ X ├┤ RZ(-π/8) ├┤ X ├─X─┤ RZ(-π/4) ├─────────────────»\n",
      "«     ├─────────┤└─┬─┘└──────────┘└─┬─┘   └──────────┘                 »\n",
      "«q_3: ┤ RZ(π/8) ├──■────────────────■──────────────────────────────────»\n",
      "«     └─────────┘                                                      »\n",
      "«c: 4/═════════════════════════════════════════════════════════════════»\n",
      "«                                                                      »\n",
      "«                     ┌──────────┐                                             »\n",
      "«q_0: ─────────────■──┤ RZ(-π/8) ├─────────────────────────────────────────────»\n",
      "«     ┌─────────┐┌─┴─┐└──┬───┬───┘┌──────────┐                                 »\n",
      "«q_1: ┤ RZ(π/4) ├┤ X ├───┤ H ├────┤ RZ(-π/4) ├──■───────────────■──────────────»\n",
      "«     └─────────┘└───┘   └───┘    └──────────┘┌─┴─┐┌─────────┐┌─┴─┐┌──────────┐»\n",
      "«q_2: ────────────────────────────────────────┤ X ├┤ RZ(π/4) ├┤ X ├┤ RZ(-π/8) ├»\n",
      "«                                             └───┘└─────────┘└───┘└──────────┘»\n",
      "«q_3: ─────────────────────────────────────────────────────────────────────────»\n",
      "«                                                                              »\n",
      "«c: 4/═════════════════════════════════════════════════════════════════════════»\n",
      "«                                                                              »\n",
      "«                               \n",
      "«q_0: ──■───────────────■───────\n",
      "«       │               │       \n",
      "«q_1: ──┼───────────────┼───────\n",
      "«     ┌─┴─┐┌─────────┐┌─┴─┐┌───┐\n",
      "«q_2: ┤ X ├┤ RZ(π/8) ├┤ X ├┤ H ├\n",
      "«     └───┘└─────────┘└───┘└───┘\n",
      "«q_3: ──────────────────────────\n",
      "«                               \n",
      "«c: 4/══════════════════════════\n",
      "«                               \n"
     ]
    }
   ],
   "source": [
    "from typing import Any\n",
    "\n",
    "def qft_dagger(n):\n",
    "    \"\"\"Given the number of qubit registers <n>, return the QFT^dagger circuit.\"\"\"\n",
    "    return qft(n).dagger()\n",
    "\n",
    "def qpe(psi : 'QCircuit', n : int, angle : float):\n",
    "    \"\"\"\n",
    "    Given an eigenvector <psi> and the number of counting registers <n>,\n",
    "    perform the Quantum Phase Algorithm and phase rotation angle <angle>\n",
    "    \"\"\"\n",
    "    circuit = psi\n",
    "    for i in range(n):\n",
    "        circuit += tq.gates.H(target=i)\n",
    "    repetitions = 1\n",
    "    for counting_qubit in range(n):\n",
    "        for i in range(repetitions):\n",
    "            circuit += tq.gates.Phase(angle, counting_qubit, n)\n",
    "        repetitions *= 2\n",
    "    circuit += qft_dagger(3)\n",
    "    return circuit\n",
    "\n",
    "psi = tq.gates.X(target=3)\n",
    "angle = pi/4\n",
    "circuit = qpe(psi, 3, angle)\n",
    "\n",
    "tq.draw(circuit)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "+100.0000|001> \n"
     ]
    }
   ],
   "source": [
    "print(tq.simulate(circuit, samples=100, read_out_qubits=[2,1,0]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So the results of the QPE algorithm gave us, with hundred percent accuracy, the state $\\lvert 001 \\rangle$. Thus, we get that $2^n\\theta = 2^3\\theta = b(001) = 1$, which gives us:\n",
    "$$ \\theta = \\frac{1}{8}$$\n",
    "as expected."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 4. References <a id='references'></a>\n",
    "1. Quantum Fourier Transform. (2020). Retrieved 21 December 2020, from https://qiskit.org/textbook/ch-algorithms/quantum-fourier-transform.html\n",
    "2. Quantum Fourier Tranform. (2020). Retrieved 21 December 2020, from https://viterbi-web.usc.edu/~tbrun/Course/lecture13.pdf\n",
    "3. Quantum Phase Estimation. (2020). Retrieved 31 January 2021, from https://qiskit.org/textbook/ch-algorithms/quantum-fourier-transform.html"
   ]
  }
 ],
 "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.8.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}