{
  "cells": [
    {
      "cell_type": "markdown",
      "id": "39a1fd9b-a17f-486e-aa7a-38d46b867984",
      "metadata": {},
      "source": [
        "---\n",
        "title: Utility I\n",
        "description: In this lesson, we run a utility-scale circuit which appears in what we informally call, the utility paper, published in Nature Vol 618, 15 June 2023.\n",
        "---\n",
        "\n",
        "# Utility-scale experiment I\n",
        "\n",
        "<Admonition type=\"note\">\n",
        "  Tamiya Onodera (5 July 2024)\n",
        "\n",
        "  [Download the pdf](https://ibm.ent.box.com/public/static/k38xhk5yaqy9b3gjdmu0ynhnd9i4c3e1.zip) of the original lecture. Note that some code snippets might become deprecated since these are static images.\n",
        "\n",
        "  *Approximate QPU time to run this experiment is 45 seconds.*\n",
        "</Admonition>\n",
        "\n"
      ]
    },
    {
      "attachments": {},
      "cell_type": "markdown",
      "id": "f0d73c71-11ba-4dcb-b230-ecea500e8732",
      "metadata": {},
      "source": [
        "## 1. Introduction to the utility paper\n",
        "\n",
        "In this lesson, we run a utility-scale circuit which appears in what we informally call [\"the utility paper\"](https://www.nature.com/articles/s41586-023-06096-3) published in Nature Vol 618, 15 June 2023. The paper deals with the time evolution of the 2D transverse-field Ising model. In particular, they consider time dynamics of the Hamiltonian,\n",
        "\n",
        "$$\n",
        "H = H_{ZZ} + H_X = - J \\sum_{(i,j)} Z_i Z_j + h \\sum_{i} X_i\n",
        "$$\n",
        "\n",
        "in which $J > 0$ is the coupling of nearest-neighbor spins with $i < j$ and $h$ is the global transverse field.\n",
        "\n"
      ]
    },
    {
      "attachments": {},
      "cell_type": "markdown",
      "id": "c92cf6bb-df3e-4020-b980-edf17188e096",
      "metadata": {
        "jp-MarkdownHeadingCollapsed": true
      },
      "source": [
        "They simulate spin dynamics from an initial state by means of first-order Trotter decomposition of the time-evolution operator,\n",
        "\n",
        "$$\n",
        "\\begin{aligned}\n",
        " \\exp(-i H_{ZZ} \\delta t) &= \\prod_{(i,j)} \\exp (i J \\delta t Z_i Z_j)  = \\prod_{(i,j)} \\mathrm{R}_{Z_i Z_j} ( - 2 J \\delta t) \\\\\n",
        " \\exp(-i H_X \\delta t) &= \\prod_{i} \\exp (-i h \\delta t X_i )  = \\prod_{i} \\mathrm{R}_{X_i} ( 2 h \\delta t)\n",
        "\\end{aligned}\n",
        "$$\n",
        "\n",
        "in which the evolution time $T$ is discretized into $T / \\delta t$ Trotter steps and $\\mathrm{R}_{Z_i Z_j}(\\theta_J)$ and $\\mathrm{R}_{X_i}(\\theta_h)$\n",
        "are $ZZ$ and $X$ rotation gates, respectively.\n",
        "\n",
        "They ran experiments on an IBM Quantum® Eagle processor, which is a 127-qubit device with heavy-hex connectivity, applying $X$ interactions to all the qubits and $ZZ$ interactions for all the edges of the coupling map.  Note that all the $ZZ$ interactions cannot be applied simultaneously due to \"data dependance\".  Hence, they **color the coupling map** to group them into layers. Those in a layer are assigned the same color, and can be applied in parallel.\n",
        "\n",
        "In addition, for experimental simplicity, they focused on the case $\\theta_J=-\\pi /2$.\n",
        "\n",
        "The novel contribution of the paper is that they built quantum circuits at a scale beyond statevector simulation, ran them on noisy quantum computers, and succeeded in extracting reliable results. That is, they demonstrated the utility of noisy quantum computers. In doing so, they applied zero-noise extrapolation (ZNE) with probabilistic error amplification (PEA) to mitigate errors from noisy devices.\n",
        "\n",
        "From then on, we called such experiments and circuits \"utility-scale\".\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "fcb359b4-469a-4e37-b260-8bf5044d1f60",
      "metadata": {},
      "source": [
        "### 1.1 Your goal\n",
        "\n",
        "Your goal in this lesson is to build a utility-scale circuit and run it on an Eagle processor.  It is beyond the scope of this notebook to extract reliable results, partly because PEA is an experimental feature of Qiskit at this time of writing and partly because applying ZNE with PEA will take a fair amount of time.\n",
        "\n",
        "Concretely, you are asked to build and run the circuit corresponding to Figure 4b of the paper, and plot the \"unmitigated\" points of your own.  As you see, it is a 127-qubit $\\times$ 60-layer (20 Trotter steps) circuit with the $\\langle Z_{62} \\rangle$ as an observable.\n",
        "\n"
      ]
    },
    {
      "attachments": {},
      "cell_type": "markdown",
      "id": "a4aecf28-fc89-42db-ad8f-67e8afacff87",
      "metadata": {},
      "source": [
        "![image.png](https://quantum.cloud.ibm.com/learning/images/courses/utility-scale-quantum-computing/utility-i/utility1.avif)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "b8c5a5ae-fc75-4fa9-9c70-2de36652f49a",
      "metadata": {},
      "source": [
        "Sound like a big deal?    Don't worry. The last three lessons of this course provide stepping stones. To begin, we will demonstrate a smaller-scale experiment which is to build and run on a fake device a 27-qubit $\\times$ 6-layer (2 Trotter steps) circuit with the $\\langle Z_{13} \\rangle$ as an observable.\n",
        "\n",
        "That's all for the introduction. Let us go on a utility-scale adventure!\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "id": "ad57179f-240e-441d-b2e4-e531cd192b45",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "'2.0.2'"
            ]
          },
          "execution_count": 1,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "import qiskit\n",
        "\n",
        "qiskit.__version__"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "id": "35fff60f-7c73-4492-ace9-ecdc8fd33b0c",
      "metadata": {},
      "outputs": [],
      "source": [
        "#!pip install qiskit_ibm_runtime\n",
        "#!pip install qiskit_aer"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "id": "2f7464cc-4190-400f-b947-4f38100b43d2",
      "metadata": {},
      "outputs": [],
      "source": [
        "import matplotlib.pyplot as plt\n",
        "import numpy as np\n",
        "import rustworkx as rx\n",
        "\n",
        "from qiskit import QuantumCircuit, transpile\n",
        "from qiskit.circuit import Parameter\n",
        "from qiskit.circuit.library import YGate\n",
        "from qiskit.quantum_info import SparsePauliOp\n",
        "from qiskit_ibm_runtime import (\n",
        "    QiskitRuntimeService,\n",
        "    fake_provider,\n",
        "    EstimatorV2 as Estimator,\n",
        ")\n",
        "from qiskit_aer import AerSimulator"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "acd628f1-e010-44b7-ba77-cceb3b4198c6",
      "metadata": {},
      "outputs": [],
      "source": [
        "service = QiskitRuntimeService()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "06556d62-08bc-4c43-a7d0-382b73ec239a",
      "metadata": {},
      "source": [
        "## 2. Preparation\n",
        "\n",
        "### 2.1 Construct RZZ(-$\\pi$ / 2)\n",
        "\n",
        "First, observe that the RZZ gate in general requires two $CX$ gates.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 5,
      "id": "d0b96520-3a4b-43c4-ade8-85cec68e3e09",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/learning/images/courses/utility-scale-quantum-computing/utility-i/extracted-outputs/d0b96520-3a4b-43c4-ade8-85cec68e3e09-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 5,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "from qiskit.circuit.library import RZZGate\n",
        "\n",
        "θ_h = Parameter(\"$\\\\theta_h$\")\n",
        "qc1 = QuantumCircuit(2)\n",
        "qc1.append(RZZGate(θ_h), [0, 1])\n",
        "qc1.decompose(reps=1).draw(\"mpl\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "2770e3d2-d9d1-4bd5-8f00-8d8401f10438",
      "metadata": {},
      "source": [
        "As mentioned above, we focus on the RZZ gate with a specific angle, -$\\pi$ / 2, for this experiment.  As shown in the paper, it can be realized with just one $CX$ gate.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 6,
      "id": "14b0ecb2-a8da-484a-bf49-606450bb074a",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/learning/images/courses/utility-scale-quantum-computing/utility-i/extracted-outputs/14b0ecb2-a8da-484a-bf49-606450bb074a-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 6,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "qc2 = QuantumCircuit(2)\n",
        "\n",
        "qc2.sdg([0, 1])\n",
        "qc2.append(YGate().power(1 / 2), [1])\n",
        "qc2.cx(0, 1)\n",
        "qc2.append(YGate().power(1 / 2).adjoint(), [1])\n",
        "\n",
        "qc2.draw(\"mpl\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "269a5094-b4e6-4249-b9d5-ea25ded1767e",
      "metadata": {},
      "source": [
        "We define a gate in terms of this circuit for future reference.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 7,
      "id": "2b5fbb1a",
      "metadata": {},
      "outputs": [],
      "source": [
        "rzz = qc2.to_gate(label=\"RZZ\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "d89cccac-13a8-49b5-9dab-c1e98f9735c8",
      "metadata": {},
      "source": [
        "Let us make a random use of the newly-defined `rzz`.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 8,
      "id": "6753b9a4-90b1-4981-82b2-8177ac258af2",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/learning/images/courses/utility-scale-quantum-computing/utility-i/extracted-outputs/6753b9a4-90b1-4981-82b2-8177ac258af2-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "qc3 = QuantumCircuit(3)\n",
        "qc3.append(rzz, [0, 1])\n",
        "qc3.append(rzz, [0, 2])\n",
        "display(qc3.draw(\"mpl\"))\n",
        "# display(qc.decompose(reps=1).draw(\"mpl\"))"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "ed9c4589-c996-4199-8dd7-19e8d144a384",
      "metadata": {},
      "source": [
        "Before using this further, let's verify the logical equivalence of `qc1` (the RZZ gate) for -pi/2 and our newly-defined `rzz` or `qc2` gate:\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 9,
      "id": "b84c6f00-d50a-4874-9fa4-0c81742654e8",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "True"
            ]
          },
          "execution_count": 9,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "from qiskit.quantum_info import Operator\n",
        "\n",
        "op1 = Operator(qc1.assign_parameters([-np.pi / 2]))\n",
        "op2 = Operator(qc2)\n",
        "\n",
        "op1.equiv(op2)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "ace1a98d-2ef7-4fe8-b67d-7799e7566f5a",
      "metadata": {},
      "source": [
        "### 2.2 Color the coupling map\n",
        "\n",
        "Let us study how we color the coupling map of a backend.  This is needed to group $ZZ$ interactions into layers.\n",
        "\n",
        "To begin with, let us visualize the coupling map of a backend.  Note the coupling maps are heavy-hexagonal for all the current IBM Quantum devices.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 10,
      "id": "84e11b8e-bc74-4f3a-8538-d4d34d069e0f",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/learning/images/courses/utility-scale-quantum-computing/utility-i/extracted-outputs/84e11b8e-bc74-4f3a-8538-d4d34d069e0f-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 10,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "backend = service.least_busy(operational=True, simulator=False)\n",
        "\n",
        "backend.coupling_map.draw()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "a5fbfce4-c8cc-49fe-831d-31ad59f84c08",
      "metadata": {},
      "source": [
        "To color a coupling map, we use [rustworkx](https://www.rustworkx.org/), which is a Python package for working with graphs and complex networks.  It provides multiple coloring algorithms, which are all heuristic and thus aren't guaranteed to find a minimal coloring.\n",
        "\n",
        "Having said this, since heavy-hex graphs are bipartite, we pick up [graph\\_bipartite\\_edge\\_color](https://www.rustworkx.org/apiref/rustworkx.graph_bipartite_edge_color.html), which should find a minimal coloring for these graphs.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 11,
      "id": "79d3af17-efa0-4e48-89c5-ad9a962cf83b",
      "metadata": {
        "scrolled": true
      },
      "outputs": [],
      "source": [
        "def color_coupling_map(backend):\n",
        "    graph = backend.coupling_map.graph\n",
        "    undirected_graph = graph.to_undirected(multigraph=False)\n",
        "    edge_color_map = rx.graph_bipartite_edge_color(undirected_graph)\n",
        "    if edge_color_map is None:\n",
        "        edge_color_map = rx.graph_greedy_edge_color(undirected_graph)\n",
        "    # build a map from color to a list of edges\n",
        "    edge_index_map = undirected_graph.edge_index_map()\n",
        "    color_edges_map = {color: [] for color in edge_color_map.values()}\n",
        "    for edge_index, color in edge_color_map.items():\n",
        "        color_edges_map[color].append(\n",
        "            (edge_index_map[edge_index][0], edge_index_map[edge_index][1])\n",
        "        )\n",
        "    return edge_color_map, color_edges_map"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "1a409faf-cb52-4c1c-95bd-beeca2cfd379",
      "metadata": {},
      "source": [
        "Heavy-hexagonal graphs should be painted in three colors.  Let us check this for the above coupling map.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 12,
      "id": "5cf7f215-2155-498e-a8e8-acc51ccba363",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "ibm_strasbourg, 127-qubit device, 3 colors assigned.\n"
          ]
        }
      ],
      "source": [
        "edge_color_map, color_edges_map = color_coupling_map(backend)\n",
        "print(\n",
        "    f\"{backend.name}, {backend.num_qubits}-qubit device, {len(color_edges_map.keys())} colors assigned.\"\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "574bb71f-436e-499c-bb23-1adbab9445db",
      "metadata": {},
      "source": [
        "Yes, it is!\n",
        "\n",
        "For fun, let us paint the coupling map to the coloring obtained, using the feature of rustworks visualization.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 13,
      "id": "fb9b6915-5ace-4104-b81d-91666e689f0d",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/learning/images/courses/utility-scale-quantum-computing/utility-i/extracted-outputs/fb9b6915-5ace-4104-b81d-91666e689f0d-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 13,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "color_str_map = {0: \"green\", 1: \"red\", 2: \"blue\"}\n",
        "\n",
        "undirected_graph = backend.coupling_map.graph.to_undirected(multigraph=False)\n",
        "for i in undirected_graph.edge_indices():\n",
        "    undirected_graph.get_edge_data_by_index(i)[\"color\"] = color_str_map[\n",
        "        edge_color_map[i]\n",
        "    ]\n",
        "\n",
        "rx.visualization.graphviz_draw(\n",
        "    undirected_graph, method=\"neato\", edge_attr_fn=lambda edge: {\"color\": edge[\"color\"]}\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "4ae34f8d-6a46-4d81-bd0d-66a0d6fca588",
      "metadata": {},
      "source": [
        "## 3. Solve the Trotterized time evolution of a 2D Ising model.\n",
        "\n",
        "Let us define a routine to construct a circuit of the utility paper for the time evolution of a 2D Ising model.  The routine takes three parameters, a backend, an integer indicating the number of Trotter steps, and a Boolean controlling the barrier insertion.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 14,
      "id": "b9d5dfd9-997f-4d8c-af1f-7927df31bffb",
      "metadata": {},
      "outputs": [],
      "source": [
        "def get_utility_circuit(backend, num_steps: int, barrier: bool = False):\n",
        "    num_qubits = backend.num_qubits\n",
        "    _, color_edges_map = color_coupling_map(backend)\n",
        "    θ_h = Parameter(\"$\\\\theta_h$\")\n",
        "    qc = QuantumCircuit(num_qubits)\n",
        "\n",
        "    for i in range(num_steps):\n",
        "        qc.rx(θ_h, range(num_qubits))\n",
        "\n",
        "        for _, edge_list in color_edges_map.items():\n",
        "            for edge in edge_list:\n",
        "                qc.append(rzz, edge)\n",
        "\n",
        "        if barrier:\n",
        "            qc.barrier()\n",
        "    return qc"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "e6689c89-8321-44e8-8983-82d6b0735751",
      "metadata": {},
      "source": [
        "Please note that we have already manually performed qubit mapping and routing for the circuit constructed.  Thus, when we transpile the circuit later, we do not (**should not**) ask the transpiler to do qubit mapping and routing. As you will soon see, we invoke it with the optimization level as 1 and the layout method as \"trivial\".\n",
        "\n",
        "Next we define a simple routine to get the information about the circuit constructed for a quick check.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 15,
      "id": "5c002120",
      "metadata": {},
      "outputs": [],
      "source": [
        "def get_circuit_info(qc: QuantumCircuit, reps: int = 0):\n",
        "    qc0 = qc.decompose(reps=reps)\n",
        "    return (\n",
        "        f\"{qc0.num_qubits} qubits × {qc0.depth(lambda x: x.operation.num_qubits == 2)} layers ({qc0.depth()}-depth)\"\n",
        "        + \",  \"\n",
        "        + f\"\"\"Gate breakdown: {\", \".join([f\"{k.upper()} {v}\" for k, v in qc0.count_ops().items()])}\"\"\"\n",
        "    )"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "30b69c52-09f0-4d5a-855c-e8f90acd0984",
      "metadata": {},
      "source": [
        "Let us exercise these routines.  You should see a circuit of 27 qubit $\\times$ 15 layers (5 Trotter steps). Since the fake device has 28 edges, there should be 28\\*5 entangling gages.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 16,
      "id": "af7bdb29-3bdb-4f6e-a53c-76ce935dff03",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/learning/images/courses/utility-scale-quantum-computing/utility-i/extracted-outputs/af7bdb29-3bdb-4f6e-a53c-76ce935dff03-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        },
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "27 qubits × 15 layers (20-depth),  Gate breakdown: CIRCUIT-165 140, RX 135, BARRIER 5\n",
            "27 qubits × 15 layers (60-depth),  Gate breakdown: SDG 280, UNITARY 280, CX 140, R 135, BARRIER 5\n"
          ]
        }
      ],
      "source": [
        "backend = fake_provider.FakeTorontoV2()\n",
        "num_steps = 5\n",
        "qc = get_utility_circuit(backend, num_steps, True)\n",
        "\n",
        "display(qc.draw(output=\"mpl\", fold=-1))\n",
        "print(get_circuit_info(qc, reps=0))\n",
        "print(get_circuit_info(qc, reps=1))"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "bf038a4d-4482-4368-8171-8cb4cd828c51",
      "metadata": {},
      "source": [
        "## 4. Solve the 27-qubit version of the problem.\n",
        "\n",
        "We now demonstrate a smaller-scale version of the utility experiment. We build a 27-qubit $\\times$ 6-layer (2 Trotter steps) circuit with $\\langle Z_{13} \\rangle$ as an observable, and run it both on AerSimulator and a fake device.\n",
        "\n",
        "Of course, we follow our four-step workflow, \"Qiskit pattern\", which consists of Map, Optimize, Execute, and Post-Process. More concretely,\n",
        "\n",
        "* Map classical inputs to a quantum computation.\n",
        "* Optimize circuits for quantum computation.\n",
        "* Execute circuits using primitives.\n",
        "* Post-process and return results in classical format.\n",
        "\n",
        "In the following, we have the Map step for creating a circuit for a smaller-scale experiment. We then have one set of Optimize and Execute for AerSimulator and another for a fake device.  Finally, we have the Post-Process step to plot the results.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "ad5f3abd-0f25-47ec-9541-16fe83b7e710",
      "metadata": {},
      "source": [
        "### 4.1 Step 1: Map\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 17,
      "id": "cd9c6774",
      "metadata": {},
      "outputs": [],
      "source": [
        "backend = fake_provider.FakeTorontoV2()  # a 27 qubit fake device.\n",
        "num_steps = 2\n",
        "qc = get_utility_circuit(backend, num_steps)\n",
        "obs = SparsePauliOp.from_sparse_list(\n",
        "    [(\"Z\", [13], 1)], num_qubits=backend.num_qubits\n",
        ")  # Falcon\n",
        "angles = [\n",
        "    0,\n",
        "    0.1,\n",
        "    0.2,\n",
        "    0.3,\n",
        "    0.4,\n",
        "    0.5,\n",
        "    0.6,\n",
        "    0.7,\n",
        "    0.8,\n",
        "    1.0,\n",
        "    np.pi / 2,\n",
        "]  # We try 11 angles for theta_h."
      ]
    },
    {
      "cell_type": "markdown",
      "id": "81b39024-a7f4-4a7a-a543-59185def648b",
      "metadata": {},
      "source": [
        "### 4.2 Steps 2 and 3: Optimize and Execute (Simulator)\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 18,
      "id": "fa307165-fc74-4b6c-88c5-a2ac36341a3c",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "27 qubits × 6 layers (23-depth),  Gate breakdown: SDG 112, UNITARY 112, CX 56, R 54\n",
            "27 qubits × 6 layers (16-depth),  Gate breakdown: U3 80, CX 56, R 54, U1 32, U 28\n"
          ]
        }
      ],
      "source": [
        "backend_sim = AerSimulator()\n",
        "transpiled_qc_sim = transpile(\n",
        "    qc, backend_sim, optimization_level=1, layout_method=\"trivial\"\n",
        ")\n",
        "transpiled_obs_sim = obs.apply_layout(layout=transpiled_qc_sim.layout)\n",
        "\n",
        "print(get_circuit_info(qc, reps=1))\n",
        "print(get_circuit_info(transpiled_qc_sim, reps=1))"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "0a101b40-5555-41d9-96b2-a8c986ea9903",
      "metadata": {},
      "source": [
        "One user ran the next cell using a MacBook Pro with 2.3 GHz quad-core Intel Core i7 Processor equipped with 32GB 3LPDDR4X RAM, which is running macOS 14.5.  It took 161ms in wall time. Each laptop will be slightly different.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 19,
      "id": "178dd8c4-439b-4324-ba1f-687acc57a678",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "CPU times: user 231 ms, sys: 186 ms, total: 417 ms\n",
            "Wall time: 111 ms\n"
          ]
        }
      ],
      "source": [
        "%%time\n",
        "params = [[p] for p in angles]\n",
        "estimator = Estimator(mode=backend_sim)\n",
        "pub = (transpiled_qc_sim, transpiled_obs_sim, params)\n",
        "result_sim = estimator.run([pub]).result()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "48b684f9-061a-49fd-acb1-3ed7b4bab1f2",
      "metadata": {},
      "source": [
        "### 4.3 Steps 2 and 3: Optimize and execute (fake device)\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 20,
      "id": "f944eeae-8dc5-4ad6-97ad-92f95159956b",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "27 qubits × 6 layers (23-depth),  Gate breakdown: SDG 112, UNITARY 112, CX 56, R 54\n",
            "27 qubits × 6 layers (49-depth),  Gate breakdown: SDG 324, U1 274, H 162, CX 56, U3 14\n"
          ]
        }
      ],
      "source": [
        "backend_fake = fake_provider.FakeTorontoV2()\n",
        "transpiled_qc_fake = transpile(\n",
        "    qc, backend_fake, optimization_level=1, layout_method=\"trivial\"\n",
        ")\n",
        "transpiled_obs_fake = obs.apply_layout(layout=transpiled_qc_fake.layout)\n",
        "\n",
        "print(get_circuit_info(qc, reps=1))\n",
        "print(get_circuit_info(transpiled_qc_fake, reps=1))"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "59adfdfb-0954-42f7-9dd2-c819d04264b6",
      "metadata": {},
      "source": [
        "When the same user ran the next cell with the same environment as above, it took 2min 19s in Wall Time.  Executing a circuit on a fake device invokes noisy simulation which takes much longer time than exact simulation.   We recommend you not execute a larger circuit (such as a 27-qubit $\\times$ 9 layers with 3 Trotter steps) on a fake device.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 21,
      "id": "62d0e0bb-9acb-474a-8d19-b132013131a0",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "CPU times: user 4min 42s, sys: 9.35 s, total: 4min 51s\n",
            "Wall time: 38.3 s\n"
          ]
        }
      ],
      "source": [
        "%%time\n",
        "params = [[p] for p in angles]\n",
        "estimator = Estimator(mode=backend_fake)\n",
        "pub = (transpiled_qc_fake, transpiled_obs_fake, params)\n",
        "result_fake = estimator.run([pub]).result()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "181e5f61-2830-485c-bd22-a1d70b443ed9",
      "metadata": {},
      "source": [
        "### 4.4 Step 4: Post-process\n",
        "\n",
        "We plot the results from exact and noisy simulations. You see the severe effects of noise on FakeToronto.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 22,
      "id": "4cf3cad0-e720-469c-a9ed-e397dae0f2f6",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/learning/images/courses/utility-scale-quantum-computing/utility-i/extracted-outputs/4cf3cad0-e720-469c-a9ed-e397dae0f2f6-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "plt.plot(angles, result_fake[0].data.evs, \"o\", label=\"Fake Device\")\n",
        "plt.plot(angles, result_sim[0].data.evs, \"o\", label=\"AerSimulator\")\n",
        "plt.xlabel(\"$\\\\mathrm{R_x}$ angle $\\\\theta_h$\")\n",
        "plt.title(\"$\\\\langle Z_{13} \\\\rangle$\")\n",
        "plt.legend()\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "1b05966c-efe2-4eb4-b63c-43c07cf81f41",
      "metadata": {},
      "source": [
        "## 5. Solve the 127-qubit version of the problem\n",
        "\n",
        "It is your goal to run the utility-scale experiment as mentioned at the beginning. You will create and execute a 127-qubit and 60-layer (20 Trotter steps) circuit with  $\\langle Z_{62} \\rangle$ as an observable. We recommend you try to do this yourself, using the code for the 27-qubit version when appropriate. But the solution is provided here.\n",
        "\n",
        "**Solution:**\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "a18a42dd-5891-4b5f-8c3b-1bc7e1ae4a5e",
      "metadata": {},
      "source": [
        "### 5.1 Step 1: Map\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "8b900d9a-0351-426a-b25b-19e631fd257e",
      "metadata": {},
      "outputs": [],
      "source": [
        "# backend_map = service.backend(\"ibm_brisbane\")\n",
        "backend_map = service.least_busy(operational=True, simulator=False)\n",
        "\n",
        "num_steps = 20\n",
        "qc = get_utility_circuit(backend_map, num_steps)\n",
        "obs = SparsePauliOp.from_sparse_list(\n",
        "    [(\"Z\", [62], 1)], num_qubits=backend_map.num_qubits\n",
        ")  # Eagle\n",
        "angles = [0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 1.0, np.pi / 2]"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "51ea0288-d8e9-4ee1-bbe8-dc24826a5026",
      "metadata": {},
      "source": [
        "### 5.2 Steps 2 and 3: Optimize and execute\n",
        "\n",
        "We note that the coupling map of the Eagle processor has 144 edges.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "f9fadf76-48e3-43a3-b773-8fa3b5fd5f3d",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "156 qubits × 60 layers (221-depth),  Gate breakdown: SDG 7040, UNITARY 7040, CX 3520, R 3120\n",
            "156 qubits × 60 layers (201-depth),  Gate breakdown: RZ 11933, SX 6240, CZ 3520\n"
          ]
        }
      ],
      "source": [
        "# backend = service.backend(\"ibm_brisbane\")\n",
        "backend = backend_map\n",
        "\n",
        "transpiled_qc = transpile(qc, backend, optimization_level=1, layout_method=\"trivial\")\n",
        "transpiled_obs = obs.apply_layout(layout=transpiled_qc.layout)\n",
        "\n",
        "print(get_circuit_info(qc, reps=1))\n",
        "print(get_circuit_info(transpiled_qc))"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 27,
      "id": "2ddff1e1-c868-4300-bdbb-cdd3b6932302",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "job id=d1479n6qf56g0081sxa0\n"
          ]
        }
      ],
      "source": [
        "params = [[p] for p in angles]\n",
        "estimator = Estimator(mode=backend)\n",
        "pub = (transpiled_qc, transpiled_obs, params)\n",
        "job = estimator.run([pub])\n",
        "\n",
        "job_id = job.job_id()\n",
        "print(f\"job id={job_id}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c497fc9d-aa58-4eb2-851e-369c7209da97",
      "metadata": {},
      "source": [
        "### 5.3 Post-process\n",
        "\n",
        "We provide the values for the \"mitigated\" points in Figure 4b of the utility paper.   Plot these together with your results.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 28,
      "id": "1c4ede0c-e865-41fe-9b38-8ecf39afa64f",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/learning/images/courses/utility-scale-quantum-computing/utility-i/extracted-outputs/1c4ede0c-e865-41fe-9b38-8ecf39afa64f-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "result_paper = [\n",
        "    1.0171,\n",
        "    1.0044,\n",
        "    0.9563,\n",
        "    0.9602,\n",
        "    0.8394,\n",
        "    0.8120,\n",
        "    0.5466,\n",
        "    0.4556,\n",
        "    0.1953,\n",
        "    0.0141,\n",
        "    0.0117,\n",
        "]\n",
        "\n",
        "# REPLACE WITH YOUR OWN JOB ID\n",
        "job = service.job(job_id)\n",
        "\n",
        "plt.plot(angles, job.result()[0].data.evs, \"o\", label=f\"{job.backend().name}\")\n",
        "plt.plot(angles, result_paper, \"o\", label=\"Utility Paper\")\n",
        "plt.xlabel(\"$\\\\mathrm{R_x}$ angle $\\\\theta_h$\")\n",
        "plt.title(\"$\\\\langle Z_{62} \\\\rangle$\")\n",
        "plt.legend()\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "7fb3e103-af1d-4a9f-8016-1a270c5706be",
      "metadata": {},
      "source": [
        "Are your results similar to the \"unmitigated\" in Figure 4b?   They could be very different, depending on the device and its condition at the time of experiment.  Do not worry about the results themselves. What we will check is whether you have done the coding correctly.  If you have done so, congratulations, you have reached the starting line of the utility era.\n",
        "\n",
        "As in the Utility paper, scientists across the globe have put forth tremendous ingenuity to extracting meaningful results even in the presence of noise. The end goal of this collective effort is quantum advantage: a state in which quantum computers can solve some problems of use in industry faster, with higher fidelity, or more cheaply than classical computers. This is not likely to be a single event, but rather an era in which classical reproduction of quantum results takes progressively longer, until at some point that quantum lead time becomes critically important. One thing is clear about quantum advantage: we get there only through the utility-scale experiments.  If this course results in your joining the quest, which is full of challenge and fun, we would be more than happy.\n",
        "\n"
      ]
    },
    {
      "attachments": {},
      "cell_type": "markdown",
      "id": "a346d104-299e-4152-a8ea-5e297d05b6a1",
      "metadata": {},
      "source": [
        "## Reference\n",
        "\n",
        "* Kim, Y., Eddins, A., Anand, S. et al. Evidence for the utility of quantum computing before fault tolerance. Nature 618, 500–505 (2023). [https://doi.org/10.1038/s41586-023-06096-3](https://doi.org/10.1038/s41586-023-06096-3)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "id": "a1b8767d",
      "source": "© IBM Corp., 2017-2026"
    }
  ],
  "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"
    },
    "widgets": {
      "application/vnd.jupyter.widget-state+json": {
        "state": {},
        "version_major": 2,
        "version_minor": 0
      }
    }
  },
  "nbformat": 4,
  "nbformat_minor": 5
}