{ "cells": [ { "cell_type": "markdown", "metadata": { "tags": [ "remove_cell" ] }, "source": [ "# Estimating pi ($\\pi$) using Quantum Phase Estimation Algorithm" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 1. Quick overview of the [Quantum Phase Estimation Algorithm](https://qiskit.org/textbook/ch-algorithms/quantum-phase-estimation.html)\n", "\n", "Quantum Phase Estimation (QPE) is a quantum algorithm that forms the building block of many more complex quantum algorithms. At its core, QPE solves a fairly straightforward problem: given an operator $U$ and a quantum state $\\vert\\psi\\rangle$ that is an eigenvalue of $U$ with $U\\vert\\psi\\rangle = \\exp\\left(2 \\pi i \\theta\\right)\\vert\\psi\\rangle$, can we obtain an estimate of $\\theta$?\n", "\n", "The answer is yes. The QPE algorithm gives us $2^n\\theta$, where $n$ is the number of qubits we use to estimate the phase $\\theta$.\n", "\n", "## 2. Estimating $\\pi$\n", "\n", "In this demo, we choose\n", "$$U = p(\\theta), \\vert\\psi\\rangle = \\vert1\\rangle$$\n", "where \n", "$$\n", "p(\\theta) = \\begin{bmatrix}\n", "1 & 0\\\\ 0 & \\exp(i\\theta)\n", "\\end{bmatrix}\n", "$$\n", "is one of the quantum gates available in Qiskit, and\n", "$$p(\\theta)\\vert1\\rangle = \\exp(i\\theta)\\vert1\\rangle.$$ \n", "\n", " By choosing the phase for our gate to be $\\theta = 1$, we can solve for $\\pi$ using the following two relations:\n", "\n", "1. From the output of the QPE algorithm, we measure an estimate for $2^n\\theta$. Then, $\\theta = \\text{measured} / 2^n$ \n", "2. From the definition of the $p(\\theta)$ gate above, we know that $2\\pi\\theta = 1 \\Rightarrow \\pi = 1 / 2\\theta$\n", "\n", "Combining these two relations, $\\pi = 1 / \\left(2 \\times (\\text{(measured)}/2^n)\\right)$. \n", "\n", "For detailed understanding of the QPE algorithm, please refer to the chapter dedicated to it in the Qiskit Textbook located at [qiskit.org/textbook](https://qiskit.org/textbook/ch-algorithms/quantum-phase-estimation.html)." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 3. Time to write code" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We begin by importing the necessary libraries." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "## import the necessary tools for our work\n", "from IPython.display import clear_output\n", "from qiskit import *\n", "from qiskit.visualization import plot_histogram\n", "import numpy as np\n", "import matplotlib.pyplot as plotter\n", "from qiskit.tools.monitor import job_monitor\n", "# Visualisation settings\n", "import seaborn as sns, operator\n", "sns.set_style(\"dark\")\n", "\n", "pi = np.pi" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The function qft_dagger computes the inverse Quantum Fourier Transform. For a detailed understanding of this algorithm, see the dedicated chapter for it in the [Qiskit Textbook](https://qiskit.org/textbook/ch-algorithms/quantum-fourier-transform.html)." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "## Code for inverse Quantum Fourier Transform\n", "## adapted from Qiskit Textbook at\n", "## qiskit.org/textbook\n", "\n", "def qft_dagger(circ_, n_qubits):\n", " \"\"\"n-qubit QFTdagger the first n qubits in circ\"\"\"\n", " for qubit in range(int(n_qubits/2)):\n", " circ_.swap(qubit, n_qubits-qubit-1)\n", " for j in range(0,n_qubits):\n", " for m in range(j):\n", " circ_.cp(-np.pi/float(2**(j-m)), m, j)\n", " circ_.h(j)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The next function, qpe_pre, prepares the initial state for the estimation. Note that the starting state is created by applying a Hadamard gate on the all but the last qubit, and setting the last qubit to $\\vert1\\rangle$." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "## Code for initial state of Quantum Phase Estimation\n", "## adapted from Qiskit Textbook at qiskit.org/textbook\n", "## Note that the starting state is created by applying \n", "## H on the first n_qubits, and setting the last qubit to |psi> = |1>\n", "\n", "def qpe_pre(circ_, n_qubits):\n", " circ_.h(range(n_qubits))\n", " circ_.x(n_qubits)\n", "\n", " for x in reversed(range(n_qubits)):\n", " for _ in range(2**(n_qubits-1-x)):\n", " circ_.cp(1, n_qubits-1-x, n_qubits)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next, we write a quick function, run_job, to run a quantum circuit and return the results." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "## Run a Qiskit job on either hardware or simulators\n", "\n", "def run_job(circ, backend, shots=1000, optimization_level=0):\n", " t_circ = transpile(circ, backend, optimization_level=optimization_level)\n", " qobj = assemble(t_circ, shots=shots)\n", " job = backend.run(qobj)\n", " job_monitor(job)\n", " return job.result().get_counts()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Then, load your account to use the cloud simulator or real devices." ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "tags": [ "uses-hardware" ] }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "ibmqfactory.load_account:WARNING:2021-02-18 11:32:59,136: Credentials are already in use. The existing account in the session will be replaced.\n" ] } ], "source": [ "## Load your IBMQ account if \n", "## you'd like to use the cloud simulator or real quantum devices\n", "my_provider = IBMQ.load_account()\n", "simulator_cloud = my_provider.get_backend('ibmq_qasm_simulator')\n", "device = my_provider.get_backend('ibmq_16_melbourne')" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [], "source": [ "simulator = Aer.get_backend('qasm_simulator')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Finally, we bring everything together in a function called get_pi_estimate that uses n_qubits to get an estimate for $\\pi$." ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [], "source": [ "## Function to estimate pi\n", "## Summary: using the notation in the Qiskit textbook (qiskit.org/textbook),\n", "## do quantum phase estimation with the 'phase' operator U = p(theta) and |psi> = |1>\n", "## such that p(theta)|1> = exp(2 x pi x i x theta)|1>\n", "## By setting theta = 1 radian, we can solve for pi\n", "## using 2^n x 1 radian = most frequently measured count = 2 x pi\n", "\n", "def get_pi_estimate(n_qubits):\n", "\n", " # create the circuit\n", " circ = QuantumCircuit(n_qubits + 1, n_qubits)\n", " # create the input state\n", " qpe_pre(circ, n_qubits)\n", " # apply a barrier\n", " circ.barrier()\n", " # apply the inverse fourier transform\n", " qft_dagger(circ, n_qubits)\n", " # apply a barrier\n", " circ.barrier()\n", " # measure all but the last qubits\n", " circ.measure(range(n_qubits), range(n_qubits))\n", "\n", " # run the job and get the results\n", " counts = run_job(circ, backend=simulator, shots=10000, optimization_level=0)\n", " # print(counts) \n", "\n", " # get the count that occurred most frequently\n", " max_counts_result = max(counts, key=counts.get)\n", " max_counts_result = int(max_counts_result, 2)\n", " \n", " # solve for pi from the measured counts\n", " theta = max_counts_result/2**n_qubits\n", " return (1./(2*theta))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now, run the get_pi_estimate function with different numbers of qubits and print the estimates." ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "scrolled": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Job Status: job has successfully run\n", "2 qubits, pi ≈ 2.0\n", "Job Status: job has successfully run\n", "3 qubits, pi ≈ 4.0\n", "Job Status: job has successfully run\n", "4 qubits, pi ≈ 2.6666666666666665\n", "Job Status: job has successfully run\n", "5 qubits, pi ≈ 3.2\n", "Job Status: job has successfully run\n", "6 qubits, pi ≈ 3.2\n", "Job Status: job has successfully run\n", "7 qubits, pi ≈ 3.2\n", "Job Status: job has successfully run\n", "8 qubits, pi ≈ 3.1219512195121952\n", "Job Status: job has successfully run\n", "9 qubits, pi ≈ 3.1604938271604937\n", "Job Status: job has successfully run\n", "10 qubits, pi ≈ 3.1411042944785277\n", "Job Status: job has successfully run\n", "11 qubits, pi ≈ 3.1411042944785277\n", "Job Status: job has successfully run\n", "12 qubits, pi ≈ 3.1411042944785277\n" ] } ], "source": [ "# estimate pi using different numbers of qubits\n", "nqs = list(range(2,12+1))\n", "pi_estimates = []\n", "for nq in nqs:\n", " thisnq_pi_estimate = get_pi_estimate(nq)\n", " pi_estimates.append(thisnq_pi_estimate)\n", " print(f\"{nq} qubits, pi ≈ {thisnq_pi_estimate}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "And plot all the results." ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", "\n" ], "text/plain": [ "