{
  "cells": [
    {
      "cell_type": "markdown",
      "id": "44ef87b3",
      "metadata": {},
      "source": [
        "---\n",
        "title: Simulation of kicked Ising Hamiltonian with dynamic circuits\n",
        "description: Tutorial demonstrating utility-scale dynamic circuits using a hexagonal kicked Ising model simulation\n",
        "---\n",
        "\n",
        "{/* cspell:ignore hcords ycords xcords fontsize ncol Krsulich Lishman */}\n",
        "\n",
        "# Simulation of kicked Ising Hamiltonian with dynamic circuits\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "2ae52bd4",
      "metadata": {},
      "source": [
        "*Usage estimate: 7.5 minutes on a Heron r3 processor. (NOTE: This is an estimate only. Your runtime may vary.)*\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c76f2627",
      "metadata": {},
      "source": [
        "Dynamic circuits are circuits with classical feedforward - in other words, they are mid-circuit measurements followed by classical logic operations that determine quantum operations conditioned on the classical output. In this tutorial, we simulate the kicked Ising model on a hexagonal lattice of spins and use dynamic circuits to realize interactions beyond the physical connectivity of the hardware.\n",
        "\n",
        "The Ising model has been studied extensively across areas of physics. It models spins that undergo Ising interactions between lattice sites, as well as kicks from the local magnetic field on each site. The Trotterized time evolution of the spins considered in this tutorial, taken from [\\[1\\]](#references), is given by the following unitary:\n",
        "\n",
        "$$\n",
        "U(\\theta)=\\left(\\prod_{\\langle j, k\\rangle} \\exp \\left(i \\frac{\\pi}{8} Z_j Z_k\\right)\\right)\\left(\\prod_j \\exp \\left(-i \\frac{\\theta}{2} X_j\\right)\\right)\n",
        "$$\n",
        "\n",
        "To probe the spin dynamics, we study the average magnetization of the spins at each site as a function of Trotter steps. Hence, we construct the following observable:\n",
        "\n",
        "$$\n",
        "\\langle O\\rangle =  \\frac{1}{N} \\sum_i \\langle Z_i \\rangle\n",
        "$$\n",
        "\n",
        "To realize the ZZ interaction between lattice sites, we present a solution using the dynamic circuit feature, leading to a significantly shorter two-qubit depth compared to the standard routing method with SWAP gates. On the other hand, the classical feedforward operations in dynamic circuits typically have longer execution times than quantum gates; hence, dynamic circuits have limitations and tradeoffs. We also present a way to add a dynamical decoupling sequence on idling qubits during the classical feedforward operation using the [stretch](/docs/guides/stretch) duration.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "a6c0af3b",
      "metadata": {},
      "source": [
        "## Requirements\n",
        "\n",
        "Before starting this tutorial, ensure that you have the following installed:\n",
        "\n",
        "* Qiskit SDK v2.0 or later with [visualization](/docs/api/qiskit/visualization) support\n",
        "* Qiskit Runtime v0.37 or later with visualization support (`pip install 'qiskit-ibm-runtime[visualization]'`)\n",
        "* Rustworkx graph library (`pip install rustworkx`)\n",
        "* Qiskit Aer (`pip install qiskit-aer`)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "03860584",
      "metadata": {},
      "source": [
        "## Set up\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "id": "88a21408",
      "metadata": {},
      "outputs": [],
      "source": [
        "import numpy as np\n",
        "from typing import List\n",
        "import rustworkx as rx\n",
        "import matplotlib.pyplot as plt\n",
        "from rustworkx.visualization import mpl_draw\n",
        "from qiskit.circuit import (\n",
        "    Parameter,\n",
        "    QuantumCircuit,\n",
        "    QuantumRegister,\n",
        "    ClassicalRegister,\n",
        ")\n",
        "from qiskit.transpiler import CouplingMap\n",
        "from qiskit.quantum_info import SparsePauliOp\n",
        "from qiskit.circuit.classical import expr\n",
        "from qiskit.transpiler.preset_passmanagers import (\n",
        "    generate_preset_pass_manager,\n",
        ")\n",
        "from qiskit.transpiler import PassManager\n",
        "from qiskit.circuit.library import RZGate, XGate\n",
        "from qiskit.transpiler.passes import (\n",
        "    ALAPScheduleAnalysis,\n",
        "    PadDynamicalDecoupling,\n",
        ")\n",
        "\n",
        "from qiskit.transpiler.basepasses import TransformationPass\n",
        "from qiskit.circuit.measure import Measure\n",
        "from qiskit.transpiler.passes.utils.remove_final_measurements import (\n",
        "    calc_final_ops,\n",
        ")\n",
        "from qiskit.circuit import Instruction\n",
        "\n",
        "from qiskit.visualization import plot_circuit_layout\n",
        "from qiskit.circuit.tools import pi_check\n",
        "\n",
        "from qiskit_aer import AerSimulator\n",
        "from qiskit_aer.primitives import SamplerV2 as Aer_Sampler\n",
        "\n",
        "from qiskit_ibm_runtime import (\n",
        "    QiskitRuntimeService,\n",
        "    Batch,\n",
        "    SamplerV2 as Sampler,\n",
        ")\n",
        "from qiskit.providers.exceptions import QiskitBackendNotFoundError\n",
        "from qiskit_ibm_runtime.visualization import (\n",
        "    draw_circuit_schedule_timing,\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "df65631b",
      "metadata": {},
      "source": [
        "## Step 1: Map classical inputs to a quantum circuit\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "188f83ae",
      "metadata": {},
      "source": [
        "We start by defining the lattice to simulate. We choose to work with the honeycomb (also called hexagonal) lattice, which is a planar graph with nodes of degree 3. Here, we specify the size of the lattice, the relevant circuit parameters of interest in the Trotterized dynamics. We simulate the Trotterized time evolution under the Ising model under three different $\\theta$ values of the local magnetic field.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "id": "e8bc54ef",
      "metadata": {},
      "outputs": [],
      "source": [
        "hex_rows = 3  # specify lattice size\n",
        "hex_cols = 5\n",
        "depths = range(9)  # specify Trotter steps\n",
        "zz_angle = np.pi / 8  # parameter for ZZ interaction\n",
        "max_angle = np.pi / 2  # max theta angle\n",
        "points = 3  # number of theta parameters\n",
        "\n",
        "θ = Parameter(\"θ\")\n",
        "params = np.linspace(0, max_angle, points)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "id": "0b12f364",
      "metadata": {},
      "outputs": [],
      "source": [
        "def make_hex_lattice(hex_rows=1, hex_cols=1):\n",
        "    \"\"\"Define hexagon lattice.\"\"\"\n",
        "    hex_cmap = CouplingMap.from_hexagonal_lattice(\n",
        "        hex_rows, hex_cols, bidirectional=False\n",
        "    )\n",
        "    data = list(hex_cmap.physical_qubits)\n",
        "    graph = hex_cmap.graph.to_undirected(multigraph=False)\n",
        "    edge_colors = rx.graph_misra_gries_edge_color(graph)\n",
        "    layer_edges = {color: [] for color in edge_colors.values()}\n",
        "    for edge_index, color in edge_colors.items():\n",
        "        layer_edges[color].append(graph.edge_list()[edge_index])\n",
        "    return data, layer_edges, hex_cmap, graph"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "b2286ebb",
      "metadata": {},
      "source": [
        "Let's start with a small test example:\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 4,
      "id": "c011bc1a",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/c011bc1a-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "hex_rows_test = 1\n",
        "hex_cols_test = 2\n",
        "\n",
        "data_test, layer_edges_test, hex_cmap_test, graph_test = make_hex_lattice(\n",
        "    hex_rows=hex_rows_test, hex_cols=hex_cols_test\n",
        ")\n",
        "\n",
        "# display a small example for illustration\n",
        "node_colors_test = [\"lightblue\"] * len(graph_test.node_indices())\n",
        "pos = rx.graph_spring_layout(\n",
        "    graph_test,\n",
        "    k=5 / np.sqrt(len(graph_test.nodes())),\n",
        "    repulsive_exponent=1,\n",
        "    num_iter=150,\n",
        ")\n",
        "mpl_draw(graph_test, node_color=node_colors_test, pos=pos)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "8e5f862a",
      "metadata": {},
      "source": [
        "We will use the small example for illustration and simulation. Below we also construct a large example to show the workflow can be extended to large sizes.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 5,
      "id": "ba481bd4",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "num_qubits = 46\n"
          ]
        },
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/ba481bd4-1.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "data, layer_edges, hex_cmap, graph = make_hex_lattice(\n",
        "    hex_rows=hex_rows, hex_cols=hex_cols\n",
        ")\n",
        "num_qubits = len(data)\n",
        "print(f\"num_qubits = {num_qubits}\")\n",
        "\n",
        "# display the honeycomb lattice to simulate\n",
        "node_colors = [\"lightblue\"] * len(graph.node_indices())\n",
        "pos = rx.graph_spring_layout(\n",
        "    graph,\n",
        "    k=5 / np.sqrt(num_qubits),\n",
        "    repulsive_exponent=1,\n",
        "    num_iter=150,\n",
        ")\n",
        "mpl_draw(graph, node_color=node_colors, pos=pos)\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "6eb03e83",
      "metadata": {},
      "source": [
        "### Build unitary circuits\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "f4b421fa",
      "metadata": {},
      "source": [
        "With the problem size and parameters specified, we are now ready to build the parametrized circuit that simulates the Trotterized time evolution of $U(\\theta)$ with different Trotter steps, specified by the `depth` argument. The circuit we construct has alternating layers of `Rx`($\\theta$) gates and `Rzz` gates. The `Rzz` gates realize the ZZ interactions between coupled spins, which will be placed between each lattice site specified by the `layer_edges` argument.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 6,
      "id": "9e2dc428",
      "metadata": {},
      "outputs": [],
      "source": [
        "def gen_hex_unitary(\n",
        "    num_qubits=6,\n",
        "    zz_angle=np.pi / 8,\n",
        "    layer_edges=[\n",
        "        [(0, 1), (2, 3), (4, 5)],\n",
        "        [(1, 2), (3, 4), (5, 0)],\n",
        "    ],\n",
        "    θ=Parameter(\"θ\"),\n",
        "    depth=1,\n",
        "    measure=False,\n",
        "    final_rot=True,\n",
        "):\n",
        "    \"\"\"Build unitary circuit.\"\"\"\n",
        "    circuit = QuantumCircuit(num_qubits)\n",
        "    # Build trotter layers\n",
        "    for _ in range(depth):\n",
        "        for i in range(num_qubits):\n",
        "            circuit.rx(θ, i)\n",
        "        circuit.barrier()\n",
        "        for coloring in layer_edges.keys():\n",
        "            for e in layer_edges[coloring]:\n",
        "                circuit.rzz(zz_angle, e[0], e[1])\n",
        "        circuit.barrier()\n",
        "    # Optional final rotation, set True to be consistent with Ref. [1]\n",
        "    if final_rot:\n",
        "        for i in range(num_qubits):\n",
        "            circuit.rx(θ, i)\n",
        "    if measure:\n",
        "        circuit.measure_all()\n",
        "\n",
        "    return circuit"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "24235b4a",
      "metadata": {},
      "source": [
        "Visualize the small test circuit:\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 7,
      "id": "268e6999",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/268e6999-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 7,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "circ_unitary_test = gen_hex_unitary(\n",
        "    num_qubits=len(data_test),\n",
        "    layer_edges=layer_edges_test,\n",
        "    θ=Parameter(\"θ\"),\n",
        "    depth=1,\n",
        "    measure=True,\n",
        ")\n",
        "circ_unitary_test.draw(output=\"mpl\", fold=-1)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "0a8abb0d",
      "metadata": {},
      "source": [
        "Similarly, construct the unitary circuits of the large example at different Trotter steps and the observable to estimate the expectation value.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 8,
      "id": "6c9e388a",
      "metadata": {},
      "outputs": [],
      "source": [
        "circuits_unitary = []\n",
        "for depth in depths:\n",
        "    circ = gen_hex_unitary(\n",
        "        num_qubits=num_qubits,\n",
        "        layer_edges=layer_edges,\n",
        "        θ=Parameter(\"θ\"),\n",
        "        depth=depth,\n",
        "        measure=True,\n",
        "    )\n",
        "    circuits_unitary.append(circ)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 9,
      "id": "695e2bad",
      "metadata": {},
      "outputs": [],
      "source": [
        "observables_unitary = SparsePauliOp.from_sparse_list(\n",
        "    [(\"Z\", [i], 1 / num_qubits) for i in range(num_qubits)],\n",
        "    num_qubits=num_qubits,\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "84d2cd91",
      "metadata": {},
      "source": [
        "### Build dynamic circuit implementation\n",
        "\n",
        "This section demonstrates the main dynamic circuit implementation to simulate the same Trotterized time evolution. Note that the honeycomb lattice we want to simulate does not match the heavy lattice of the hardware qubits. One straightforward way to map the circuit to the hardware is to introduce a series of SWAP operations to bring interacting qubits next to each other, to realize the ZZ interaction. Here we highlight an alternative approach using dynamic circuits as a solution, which illustrates that we can use the combination of quantum and real-time classical computation within a circuit in Qiskit to realize interactions beyond nearest-neighbor.\n",
        "\n",
        "In the dynamic circuit implementation, the ZZ interaction is effectively implemented by using ancilla qubits, mid-circuit measurement, and feedforward. To understand this, note that the ZZ rotations apply a phase factor $e^{i\\theta}$ to the state based on its parity. For two qubits, the computational basis states are $|00\\rangle$, $|01\\rangle$, $|10\\rangle$, and $|11\\rangle$. The ZZ rotation gate applies a phase factor to states $|01\\rangle$ and $|10\\rangle$ whose parity (the number of ones in the state) is odd and leaves even-parity states unchanged. The following describes how we can effectively implement ZZ interactions on two qubits using dynamic circuits.\n",
        "\n",
        "1. Compute parity into an ancilla qubit: instead of directly applying ZZ to two qubits, we introduce a third qubit, the ancilla qubit, to store the parity information of the two data qubits. We entangle the ancilla with each data qubit using CX gates from the data qubit to the ancilla qubit.\n",
        "\n",
        "2. Apply a single-qubit Z rotation to the ancilla qubit: this is because the ancilla has the parity information of the two data qubits, which effectively implements the ZZ rotation on the data qubits.\n",
        "\n",
        "3. Measure the ancilla qubit in the X basis: this is the key step that collapses the state of the ancilla qubit, and the measurement outcome tells us what has happened:\n",
        "\n",
        "   * Measure 0: when a 0 outcome is observed, we have in fact correctly applied a $ZZ(\\theta)$ rotation to our data qubits.\n",
        "\n",
        "   * Measure 1: when a 1 outcome is observed, We've applied $ZZ(\\theta + \\pi)$ instead.\n",
        "\n",
        "4. Apply correction gate when measuring 1: If we measured 1, we apply Z gates to the data qubits to \"fix\" the extra $\\pi$ phase.\n",
        "\n",
        "The resulting circuit is the following:\n",
        "\n",
        "![dynamic implementation](https://quantum.cloud.ibm.com/docs/images/tutorials/dc-hex-ising/circuit-1.avif)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "9872fed2",
      "metadata": {},
      "source": [
        "When we adopt this approach to simulate a honeycomb lattice, the resulting circuit embeds perfectly into the hardware with a heavy-hex lattice: all data qubits reside on the degree-3 sites of the lattice, which forms a hexagonal lattice. Every pair of data qubits shares an ancilla qubit residing on a degree-2 site. Below, we construct the qubit lattice for the dynamic circuit implementation, introducing ancilla qubits (shown in the darker purple circles).\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 10,
      "id": "c0e7afd2",
      "metadata": {},
      "outputs": [],
      "source": [
        "def make_lattice(hex_rows=1, hex_cols=1):\n",
        "    \"\"\"Define heavy-hex lattice and corresponding lists of data and ancilla nodes.\"\"\"\n",
        "    hex_cmap = CouplingMap.from_hexagonal_lattice(\n",
        "        hex_rows, hex_cols, bidirectional=False\n",
        "    )\n",
        "    data = list(hex_cmap.physical_qubits)\n",
        "\n",
        "    heavyhex_cmap = CouplingMap()\n",
        "    for d in data:\n",
        "        heavyhex_cmap.add_physical_qubit(d)\n",
        "\n",
        "    # make coupling map\n",
        "    a = len(data)\n",
        "    for edge in hex_cmap.get_edges():\n",
        "        heavyhex_cmap.add_physical_qubit(a)\n",
        "        heavyhex_cmap.add_edge(edge[0], a)\n",
        "        heavyhex_cmap.add_edge(edge[1], a)\n",
        "        a += 1\n",
        "    ancilla = list(range(len(data), a))\n",
        "    qubits = data + ancilla\n",
        "\n",
        "    # color edges\n",
        "    graph = heavyhex_cmap.graph.to_undirected(multigraph=False)\n",
        "    edge_colors = rx.graph_misra_gries_edge_color(graph)\n",
        "    layer_edges = {color: [] for color in edge_colors.values()}\n",
        "    for edge_index, color in edge_colors.items():\n",
        "        layer_edges[color].append(graph.edge_list()[edge_index])\n",
        "\n",
        "    # construct observable\n",
        "    obs_hex = SparsePauliOp.from_sparse_list(\n",
        "        [(\"Z\", [i], 1 / len(data)) for i in data],\n",
        "        num_qubits=len(qubits),\n",
        "    )\n",
        "\n",
        "    return (data, qubits, ancilla, layer_edges, heavyhex_cmap, graph, obs_hex)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "5c39eeab",
      "metadata": {},
      "source": [
        "Visualize the heavy-hex lattice for data qubits and ancilla qubits at a small scale:\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 11,
      "id": "2d7224ef",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "number of data qubits = 46\n",
            "number of ancilla qubits = 60\n"
          ]
        },
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/2d7224ef-1.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "(data, qubits, ancilla, layer_edges, heavyhex_cmap, graph, obs_hex) = (\n",
        "    make_lattice(hex_rows=hex_rows, hex_cols=hex_cols)\n",
        ")\n",
        "\n",
        "print(f\"number of data qubits = {len(data)}\")\n",
        "print(f\"number of ancilla qubits = {len(ancilla)}\")\n",
        "\n",
        "node_colors = []\n",
        "for node in graph.node_indices():\n",
        "    if node in ancilla:\n",
        "        node_colors.append(\"purple\")\n",
        "    else:\n",
        "        node_colors.append(\"lightblue\")\n",
        "\n",
        "pos = rx.graph_spring_layout(\n",
        "    graph,\n",
        "    k=1 / np.sqrt(len(qubits)),\n",
        "    repulsive_exponent=2,\n",
        "    num_iter=200,\n",
        ")\n",
        "\n",
        "# Visualize the graph, blue circles are data qubits and purple circles are ancillas\n",
        "mpl_draw(graph, node_color=node_colors, pos=pos)\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "bf9177d2",
      "metadata": {},
      "source": [
        "Below, we construct the dynamic circuit for the Trotterized time evolution. The `RZZ` gates are replaced with the dynamic circuit implementation using the steps described above.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 12,
      "id": "4c98664f",
      "metadata": {},
      "outputs": [],
      "source": [
        "def gen_hex_dynamic(\n",
        "    depth=1,\n",
        "    zz_angle=np.pi / 8,\n",
        "    θ=Parameter(\"θ\"),\n",
        "    hex_rows=1,\n",
        "    hex_cols=1,\n",
        "    measure=False,\n",
        "    add_dd=True,\n",
        "):\n",
        "    \"\"\"Build dynamic circuits.\"\"\"\n",
        "    (data, qubits, ancilla, layer_edges, heavyhex_cmap, graph, obs_hex) = (\n",
        "        make_lattice(hex_rows=hex_rows, hex_cols=hex_cols)\n",
        "    )\n",
        "    # Initialize circuit\n",
        "    qr = QuantumRegister(len(qubits), \"qr\")\n",
        "    cr = ClassicalRegister(len(ancilla), \"cr\")\n",
        "    circuit = QuantumCircuit(qr, cr)\n",
        "\n",
        "    for k in range(depth):\n",
        "        # Single-qubit Rx layer\n",
        "        for d in data:\n",
        "            circuit.rx(θ, d)\n",
        "        circuit.barrier()\n",
        "\n",
        "        # CX gates from data qubits to ancilla qubits\n",
        "        for same_color_edges in layer_edges.values():\n",
        "            for e in same_color_edges:\n",
        "                circuit.cx(e[0], e[1])\n",
        "        circuit.barrier()\n",
        "\n",
        "        # Apply Rz rotation on ancilla qubits and rotate into X basis\n",
        "        for a in ancilla:\n",
        "            circuit.rz(zz_angle, a)\n",
        "            circuit.h(a)\n",
        "        # Add barrier to align terminal measurement\n",
        "        circuit.barrier()\n",
        "\n",
        "        # Measure ancilla qubits\n",
        "        for i, a in enumerate(ancilla):\n",
        "            circuit.measure(a, i)\n",
        "        d2ros = {}\n",
        "        a2ro = {}\n",
        "        # Retrieve ancilla measurement outcomes\n",
        "        for a in ancilla:\n",
        "            a2ro[a] = cr[ancilla.index(a)]\n",
        "\n",
        "        # For each data qubit, retrieve measurement outcomes of neighboring\n",
        "        # ancilla qubits\n",
        "        for d in data:\n",
        "            ros = [a2ro[a] for a in heavyhex_cmap.neighbors(d)]\n",
        "            d2ros[d] = ros\n",
        "\n",
        "        # Build classical feedforward operations (optionally add DD on idling\n",
        "        # data qubits)\n",
        "        for d in data:\n",
        "            if add_dd:\n",
        "                circuit = add_stretch_dd(circuit, d, f\"data_{d}_depth_{k}\")\n",
        "\n",
        "            # # XOR the neighboring readouts of the data qubit;\n",
        "            # if True, apply Z to it\n",
        "            ros = d2ros[d]\n",
        "            parity = ros[0]\n",
        "            for ro in ros[1:]:\n",
        "                parity = expr.bit_xor(parity, ro)\n",
        "            with circuit.if_test(expr.equal(parity, True)):\n",
        "                circuit.z(d)\n",
        "\n",
        "        # Reset the ancilla if its readout is 1\n",
        "        for a in ancilla:\n",
        "            with circuit.if_test(expr.equal(a2ro[a], True)):\n",
        "                circuit.x(a)\n",
        "        circuit.barrier()\n",
        "\n",
        "    # Final single-qubit Rx layer to match the unitary circuits\n",
        "    for d in data:\n",
        "        circuit.rx(θ, d)\n",
        "\n",
        "    if measure:\n",
        "        circuit.measure_all()\n",
        "    return circuit, obs_hex\n",
        "\n",
        "\n",
        "def add_stretch_dd(qc, q, name):\n",
        "    \"\"\"Add XpXm DD sequence.\"\"\"\n",
        "    s = qc.add_stretch(name)\n",
        "    qc.delay(s, q)\n",
        "    qc.x(q)\n",
        "    qc.delay(s, q)\n",
        "    qc.delay(s, q)\n",
        "    qc.rz(np.pi, q)\n",
        "    qc.x(q)\n",
        "    qc.rz(-np.pi, q)\n",
        "    qc.delay(s, q)\n",
        "    return qc"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c8b3b632",
      "metadata": {},
      "source": [
        "#### Dynamical decoupling (DD) and support for `stretch` duration\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "1f08bda1",
      "metadata": {},
      "source": [
        "One caveat of using the dynamic circuit implementation to realize the ZZ interaction is that the mid-circuit measurement and the classical feedforward operations typically take a longer time to execute than quantum gates. To suppress qubit decoherence during the idling time for the classical operations to happen, we added a [dynamical decoupling](/docs/guides/error-mitigation-and-suppression-techniques#dynamical-decoupling) (DD) sequence after the measurement operation on the ancilla qubits, and before the conditional Z operation on the data qubit, before the `if_test` statement.\n",
        "\n",
        "The DD sequence is added by the function `add_stretch_dd()`, which uses the [`stretch` durations](/docs/guides/stretch) to determine the time intervals between the DD gates. A `stretch` duration is a way to specify a stretchable time duration for the `delay` operation such that the delay duration can grow to fill up the qubit idling time. The duration variables specified by `stretch` are resolved at compile time into desired durations that satisfy a certain constraint. This is very useful when the timing of DD sequences is essential to achieve good error suppression performance. For more details on the `stretch` type, see the [OpenQASM](https://openqasm.com/language/delays.html#duration-and-stretch-types) documentation. Currently, support for the `stretch` type in Qiskit Runtime is experimental. For details on its usage constraints, please refer to [the limitations section](/docs/guides/stretch#qiskit-runtime-limitations) of the `stretch` documentation.\n",
        "\n",
        "Using the functions defined above, we build the Trotterized time evolution circuits, with and without DD, and the corresponding observables.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "997419c8",
      "metadata": {},
      "source": [
        "We start by visualizing the dynamic circuit of a small example:\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 13,
      "id": "b6e2e76c",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/b6e2e76c-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "hex_rows_test = 1\n",
        "hex_cols_test = 1\n",
        "\n",
        "(\n",
        "    data_test,\n",
        "    qubits_test,\n",
        "    ancilla_test,\n",
        "    layer_edges_test,\n",
        "    heavyhex_cmap_test,\n",
        "    graph_test,\n",
        "    obs_hex_test,\n",
        ") = make_lattice(hex_rows=hex_rows_test, hex_cols=hex_cols_test)\n",
        "\n",
        "node_colors = []\n",
        "for node in graph_test.node_indices():\n",
        "    if node in ancilla_test:\n",
        "        node_colors.append(\"purple\")\n",
        "    else:\n",
        "        node_colors.append(\"lightblue\")\n",
        "pos = rx.graph_spring_layout(\n",
        "    graph_test,\n",
        "    k=5 / np.sqrt(len(qubits_test)),\n",
        "    repulsive_exponent=2,\n",
        "    num_iter=150,\n",
        ")\n",
        "\n",
        "# display a small example for illustration\n",
        "node_colors_test = [\"lightblue\"] * len(graph_test.node_indices())\n",
        "mpl_draw(graph_test, node_color=node_colors, pos=pos)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 14,
      "id": "735e590a",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/735e590a-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 14,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "circuit_dynamic_test, obs_dynamic_test = gen_hex_dynamic(\n",
        "    depth=1,\n",
        "    θ=Parameter(\"θ\"),\n",
        "    hex_rows=hex_rows_test,\n",
        "    hex_cols=hex_cols_test,\n",
        "    measure=False,\n",
        "    add_dd=False,\n",
        ")\n",
        "circuit_dynamic_test.draw(\"mpl\", fold=-1)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 15,
      "id": "5de9381a",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/5de9381a-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 15,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "circuit_dynamic_dd_test, _ = gen_hex_dynamic(\n",
        "    depth=1,\n",
        "    θ=Parameter(\"θ\"),\n",
        "    hex_rows=hex_rows_test,\n",
        "    hex_cols=hex_cols_test,\n",
        "    measure=False,\n",
        "    add_dd=True,\n",
        ")\n",
        "circuit_dynamic_dd_test.draw(\"mpl\", fold=-1)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "956dd43e",
      "metadata": {},
      "source": [
        "Similarly, construct the dynamic circuits for the large example:\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 16,
      "id": "bd7b2be0",
      "metadata": {},
      "outputs": [],
      "source": [
        "circuits_dynamic = []\n",
        "circuits_dynamic_dd = []\n",
        "observables_dynamic = []\n",
        "for depth in depths:\n",
        "    circuit, obs = gen_hex_dynamic(\n",
        "        depth=depth,\n",
        "        θ=Parameter(\"θ\"),\n",
        "        hex_rows=hex_rows,\n",
        "        hex_cols=hex_cols,\n",
        "        measure=True,\n",
        "        add_dd=False,\n",
        "    )\n",
        "    circuits_dynamic.append(circuit)\n",
        "\n",
        "    circuit_dd, _ = gen_hex_dynamic(\n",
        "        depth=depth,\n",
        "        θ=Parameter(\"θ\"),\n",
        "        hex_rows=hex_rows,\n",
        "        hex_cols=hex_cols,\n",
        "        measure=True,\n",
        "        add_dd=True,\n",
        "    )\n",
        "    circuits_dynamic_dd.append(circuit_dd)\n",
        "    observables_dynamic.append(obs)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "0f0fdf70",
      "metadata": {},
      "source": [
        "## Step 2: Optimize problem for hardware execution\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "86642d2b",
      "metadata": {},
      "source": [
        "We are now ready to transpile the circuit to the hardware. We will transpile both the unitary standard implementation and the dynamic circuit implementation to the hardware.\n",
        "\n",
        "To transpile to hardware, we first instantiate the backend. If available, we will choose a backend where the [`MidCircuitMeasure`](/docs/guides/measure-qubits) (`measure_2`) instruction is supported.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 17,
      "id": "ec904b18",
      "metadata": {},
      "outputs": [],
      "source": [
        "service = QiskitRuntimeService()\n",
        "try:\n",
        "    backend = service.least_busy(\n",
        "        operational=True,\n",
        "        simulator=False,\n",
        "        use_fractional_gates=True,\n",
        "        filters=lambda b: \"measure_2\" in b.supported_instructions,\n",
        "    )\n",
        "except QiskitBackendNotFoundError:\n",
        "    backend = service.least_busy(\n",
        "        operational=True,\n",
        "        simulator=False,\n",
        "        use_fractional_gates=True,\n",
        "    )"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "610622b4",
      "metadata": {},
      "source": [
        "### Transpilation for dynamic circuits\n",
        "\n",
        "First, we transpile the dynamic circuits, with and without adding the DD sequence. To ensure we use the same set of physical qubits in all circuits for more consistent results, we first transpile the circuit once, and then use its layout for all subsequent circuits, specified by [`initial_layout`](/docs/api/qiskit/qiskit.transpiler.TranspileLayout#initial_layout) in the pass manager. We then construct the [primitive unified blocs](/docs/guides/primitive-input-output) (PUBs) as the Sampler primitive input.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 18,
      "id": "aa653907",
      "metadata": {},
      "outputs": [],
      "source": [
        "pm_temp = generate_preset_pass_manager(\n",
        "    optimization_level=3,\n",
        "    backend=backend,\n",
        ")\n",
        "isa_temp = pm_temp.run(circuits_dynamic[-1])\n",
        "dynamic_layout = isa_temp.layout.initial_index_layout(filter_ancillas=True)\n",
        "\n",
        "pm = generate_preset_pass_manager(\n",
        "    optimization_level=3, backend=backend, initial_layout=dynamic_layout\n",
        ")\n",
        "\n",
        "dynamic_isa_circuits = [pm.run(circ) for circ in circuits_dynamic]\n",
        "dynamic_pubs = [(circ, params) for circ in dynamic_isa_circuits]\n",
        "\n",
        "dynamic_isa_circuits_dd = [pm.run(circ) for circ in circuits_dynamic_dd]\n",
        "dynamic_pubs_dd = [(circ, params) for circ in dynamic_isa_circuits_dd]"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "6979ad68",
      "metadata": {},
      "source": [
        "We can visualize the qubit layout of the transpiled circuit below. The black circles show the data qubits and the ancilla qubits used in the dynamic circuit implementation.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 19,
      "id": "65b2df70",
      "metadata": {},
      "outputs": [],
      "source": [
        "def _heron_coords_r2():\n",
        "    cord_map = np.array(\n",
        "        [\n",
        "            [\n",
        "                0,\n",
        "                1,\n",
        "                2,\n",
        "                3,\n",
        "                4,\n",
        "                5,\n",
        "                6,\n",
        "                7,\n",
        "                8,\n",
        "                9,\n",
        "                10,\n",
        "                11,\n",
        "                12,\n",
        "                13,\n",
        "                14,\n",
        "                15,\n",
        "                3,\n",
        "                7,\n",
        "                11,\n",
        "                15,\n",
        "                0,\n",
        "                1,\n",
        "                2,\n",
        "                3,\n",
        "                4,\n",
        "                5,\n",
        "                6,\n",
        "                7,\n",
        "                8,\n",
        "                9,\n",
        "                10,\n",
        "                11,\n",
        "                12,\n",
        "                13,\n",
        "                14,\n",
        "                15,\n",
        "                1,\n",
        "                5,\n",
        "                9,\n",
        "                13,\n",
        "                0,\n",
        "                1,\n",
        "                2,\n",
        "                3,\n",
        "                4,\n",
        "                5,\n",
        "                6,\n",
        "                7,\n",
        "                8,\n",
        "                9,\n",
        "                10,\n",
        "                11,\n",
        "                12,\n",
        "                13,\n",
        "                14,\n",
        "                15,\n",
        "                3,\n",
        "                7,\n",
        "                11,\n",
        "                15,\n",
        "                0,\n",
        "                1,\n",
        "                2,\n",
        "                3,\n",
        "                4,\n",
        "                5,\n",
        "                6,\n",
        "                7,\n",
        "                8,\n",
        "                9,\n",
        "                10,\n",
        "                11,\n",
        "                12,\n",
        "                13,\n",
        "                14,\n",
        "                15,\n",
        "                1,\n",
        "                5,\n",
        "                9,\n",
        "                13,\n",
        "                0,\n",
        "                1,\n",
        "                2,\n",
        "                3,\n",
        "                4,\n",
        "                5,\n",
        "                6,\n",
        "                7,\n",
        "                8,\n",
        "                9,\n",
        "                10,\n",
        "                11,\n",
        "                12,\n",
        "                13,\n",
        "                14,\n",
        "                15,\n",
        "                3,\n",
        "                7,\n",
        "                11,\n",
        "                15,\n",
        "                0,\n",
        "                1,\n",
        "                2,\n",
        "                3,\n",
        "                4,\n",
        "                5,\n",
        "                6,\n",
        "                7,\n",
        "                8,\n",
        "                9,\n",
        "                10,\n",
        "                11,\n",
        "                12,\n",
        "                13,\n",
        "                14,\n",
        "                15,\n",
        "                1,\n",
        "                5,\n",
        "                9,\n",
        "                13,\n",
        "                0,\n",
        "                1,\n",
        "                2,\n",
        "                3,\n",
        "                4,\n",
        "                5,\n",
        "                6,\n",
        "                7,\n",
        "                8,\n",
        "                9,\n",
        "                10,\n",
        "                11,\n",
        "                12,\n",
        "                13,\n",
        "                14,\n",
        "                15,\n",
        "                3,\n",
        "                7,\n",
        "                11,\n",
        "                15,\n",
        "                0,\n",
        "                1,\n",
        "                2,\n",
        "                3,\n",
        "                4,\n",
        "                5,\n",
        "                6,\n",
        "                7,\n",
        "                8,\n",
        "                9,\n",
        "                10,\n",
        "                11,\n",
        "                12,\n",
        "                13,\n",
        "                14,\n",
        "                15,\n",
        "            ],\n",
        "            -1\n",
        "            * np.array([j for i in range(15) for j in [i] * [16, 4][i % 2]]),\n",
        "        ],\n",
        "        dtype=int,\n",
        "    )\n",
        "\n",
        "    hcords = []\n",
        "    ycords = cord_map[0]\n",
        "    xcords = cord_map[1]\n",
        "    for i in range(156):\n",
        "        hcords.append([xcords[i] + 1, np.abs(ycords[i]) + 1])\n",
        "\n",
        "    return hcords"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 20,
      "id": "98d402e0",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/98d402e0-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 20,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "plot_circuit_layout(\n",
        "    dynamic_isa_circuits_dd[8],\n",
        "    backend,\n",
        "    qubit_coordinates=_heron_coords_r2(),\n",
        "    view=\"virtual\",\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "1bff4d50",
      "metadata": {},
      "source": [
        "<Admonition type=\"note\">\n",
        "  If you get errors about `neato` not found from `plot_circuit_layout()`, make sure you have the `graphviz` package installed and available in your PATH. If it installs into a non-default location (for example, using `homebrew` on MacOS), you may need to update your `PATH` environment variable. This can be done inside this notebook using the following:\n",
        "\n",
        "  ```python\n",
        "  import os\n",
        "  os.environ['PATH'] = f\"path/to/neato{os.pathsep}{os.environ['PATH']}\"\n",
        "  ```\n",
        "</Admonition>\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 21,
      "id": "82fb6fa8",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/82fb6fa8-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 21,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "dynamic_isa_circuits[1].draw(fold=-1, output=\"mpl\", idle_wires=False)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 22,
      "id": "99ad295c",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/99ad295c-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 22,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "dynamic_isa_circuits_dd[1].draw(fold=-1, output=\"mpl\", idle_wires=False)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "47d43b9e",
      "metadata": {},
      "source": [
        "#### Transpile using `MidCircuitMeasure`\n",
        "\n",
        "`MidCircuitMeasure` is an addition to the available measurement operations, calibrated specifically to perform [mid-circuit measurements](/docs/guides/execute-dynamic-circuits#midcircuit). The `MidCircuitMeasure` instruction maps to the `measure_2` instruction supported by the backends. Note that `measure_2` is not supported on all backends. You can use `service.backends(filters=lambda b: \"measure_2\" in b.supported_instructions)` to find backends that support it. Here, we show how to transpile the circuit so that the mid-circuit measurements defined in the circuit are executed using the `MidCircuitMeasure` operation, if the backend supports it.\n",
        "\n",
        "Below, we print the duration for the `measure_2` instruction and the standard `measure` instruction.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 23,
      "id": "de870864",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Mid-circuit measurement `measure_2` duration: 1.3800000000000003 μs\n",
            "Terminal measurement `measure` duration: 2.1800000000000006 μs\n"
          ]
        }
      ],
      "source": [
        "print(\n",
        "    f'Mid-circuit measurement `measure_2` duration: '\n",
        "    f'{backend.instruction_durations.get('measure_2',0) * backend.dt * 1e9/1e3} μs'\n",
        ")\n",
        "print(\n",
        "    f'Terminal measurement `measure` duration: '\n",
        "    f'{backend.instruction_durations.get('measure',0) * backend.dt *1e9/1e3} μs'\n",
        ")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 24,
      "id": "a1bdc9e7",
      "metadata": {},
      "outputs": [],
      "source": [
        "\"\"\"Pass that replaces terminal measures in the middle of the circuit with\n",
        "MidCircuitMeasure instructions.\"\"\"\n",
        "\n",
        "\n",
        "class ConvertToMidCircuitMeasure(TransformationPass):\n",
        "    \"\"\"This pass replaces terminal measures in the middle of the circuit with\n",
        "    MidCircuitMeasure instructions.\n",
        "    \"\"\"\n",
        "\n",
        "    def __init__(self, target):\n",
        "        super().__init__()\n",
        "        self.target = target\n",
        "\n",
        "    def run(self, dag):\n",
        "        \"\"\"Run the pass on a dag.\"\"\"\n",
        "        mid_circ_measure = None\n",
        "        for inst in self.target.instructions:\n",
        "            if isinstance(inst[0], Instruction) and inst[0].name.startswith(\n",
        "                \"measure_\"\n",
        "            ):\n",
        "                mid_circ_measure = inst[0]\n",
        "                break\n",
        "        if not mid_circ_measure:\n",
        "            return dag\n",
        "\n",
        "        final_measure_nodes = calc_final_ops(dag, {\"measure\"})\n",
        "        for node in dag.op_nodes(Measure):\n",
        "            if node not in final_measure_nodes:\n",
        "                dag.substitute_node(node, mid_circ_measure, inplace=True)\n",
        "\n",
        "        return dag\n",
        "\n",
        "\n",
        "pm = PassManager(ConvertToMidCircuitMeasure(backend.target))\n",
        "\n",
        "dynamic_isa_circuits_meas2 = [pm.run(circ) for circ in dynamic_isa_circuits]\n",
        "dynamic_pubs_meas2 = [(circ, params) for circ in dynamic_isa_circuits_meas2]\n",
        "\n",
        "dynamic_isa_circuits_dd_meas2 = [\n",
        "    pm.run(circ) for circ in dynamic_isa_circuits_dd\n",
        "]\n",
        "dynamic_pubs_dd_meas2 = [\n",
        "    (circ, params) for circ in dynamic_isa_circuits_dd_meas2\n",
        "]"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "a1e982ee",
      "metadata": {},
      "source": [
        "### Transpilation for unitary circuits\n",
        "\n",
        "To establish a fair comparison between the dynamic circuits and their unitary counterpart, we use the same set of physical qubits used in the dynamic circuits for the data qubits as the layout for transpiling the unitary circuits.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 25,
      "id": "87962c09",
      "metadata": {},
      "outputs": [],
      "source": [
        "init_layout = [\n",
        "    dynamic_layout[ind] for ind in range(circuits_unitary[0].num_qubits)\n",
        "]\n",
        "\n",
        "\n",
        "pm = generate_preset_pass_manager(\n",
        "    target=backend.target,\n",
        "    initial_layout=init_layout,\n",
        "    optimization_level=3,\n",
        ")\n",
        "\n",
        "\n",
        "def transpile_minimize(circ: QuantumCircuit, pm: PassManager, iterations=10):\n",
        "    \"\"\"Transpile circuits for specified number of iterations and return the one\n",
        "    with smallest two-qubit gate depth\"\"\"\n",
        "    circs = [pm.run(circ) for i in range(iterations)]\n",
        "    circs_sorted = sorted(\n",
        "        circs,\n",
        "        key=lambda x: x.depth(lambda x: x.operation.num_qubits == 2),\n",
        "    )\n",
        "    return circs_sorted[0]\n",
        "\n",
        "\n",
        "unitary_isa_circuits = []\n",
        "for circ in circuits_unitary:\n",
        "    circ_t = transpile_minimize(circ, pm, iterations=100)\n",
        "    unitary_isa_circuits.append(circ_t)\n",
        "\n",
        "unitary_pubs = [(circ, params) for circ in unitary_isa_circuits]"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "e7da0161",
      "metadata": {},
      "source": [
        "We visualize the qubit layout of the transpiled unitary circuits. The black circles indicate the physical qubits used to transpile the unitary circuits and their indices correspond to the virtual qubit indices. By comparing this with the layout plotted for the dynamic circuits, we can confirm that the unitary circuits use the same set of physical qubits as the data qubits in the dynamic circuits.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 26,
      "id": "8c3c633f",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/8c3c633f-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 26,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "plot_circuit_layout(\n",
        "    unitary_isa_circuits[-1],\n",
        "    backend,\n",
        "    qubit_coordinates=_heron_coords_r2(),\n",
        "    view=\"virtual\",\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "2fd70904",
      "metadata": {},
      "source": [
        "We now add the DD sequence to the transpiled circuits and construct the corresponding PUBs for job submission.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 27,
      "id": "383ba663",
      "metadata": {},
      "outputs": [],
      "source": [
        "pm_dd = PassManager(\n",
        "    [\n",
        "        ALAPScheduleAnalysis(target=backend.target),\n",
        "        PadDynamicalDecoupling(\n",
        "            dd_sequence=[\n",
        "                XGate(),\n",
        "                RZGate(np.pi),\n",
        "                XGate(),\n",
        "                RZGate(-np.pi),\n",
        "            ],\n",
        "            spacing=[1 / 4, 1 / 2, 0, 0, 1 / 4],\n",
        "            target=backend.target,\n",
        "        ),\n",
        "    ]\n",
        ")\n",
        "\n",
        "unitary_isa_circuits_dd = pm_dd.run(unitary_isa_circuits)\n",
        "unitary_pubs_dd = [(circ, params) for circ in unitary_isa_circuits_dd]"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c4980277",
      "metadata": {},
      "source": [
        "### Compare two-qubit gate depth of unitary and dynamic circuits\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 28,
      "id": "36f1d72d",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<matplotlib.legend.Legend at 0x12628b0e0>"
            ]
          },
          "execution_count": 28,
          "metadata": {},
          "output_type": "execute_result"
        },
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/36f1d72d-1.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "# compare circuit depth of unitary and dynamic circuit implementations\n",
        "unitary_depth = [\n",
        "    unitary_isa_circuits[i].depth(lambda x: x.operation.num_qubits == 2)\n",
        "    for i in range(len(unitary_isa_circuits))\n",
        "]\n",
        "\n",
        "dynamic_depth = [\n",
        "    dynamic_isa_circuits[i].depth(lambda x: x.operation.num_qubits == 2)\n",
        "    for i in range(len(dynamic_isa_circuits))\n",
        "]\n",
        "\n",
        "plt.plot(\n",
        "    list(range(len(unitary_depth))),\n",
        "    unitary_depth,\n",
        "    label=\"unitary circuits\",\n",
        "    color=\"#be95ff\",\n",
        ")\n",
        "plt.plot(\n",
        "    list(range(len(dynamic_depth))),\n",
        "    dynamic_depth,\n",
        "    label=\"dynamic circuits\",\n",
        "    color=\"#ff7eb6\",\n",
        ")\n",
        "plt.xlabel(\"Trotter steps\")\n",
        "plt.ylabel(\"Two-qubit depth\")\n",
        "plt.legend()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c566243b",
      "metadata": {},
      "source": [
        "The main benefit of the measurement-based circuit is that, when implementing multiple ZZ interactions, the CX layers can be parallelized, and measurements can occur simultaneously. This is because all ZZ interactions commute, so the computation can be performed with measurement depth 1. After transpiling the circuits, we observe that the dynamic circuit approach yields a significantly shorter two-qubit depth than the standard unitary approach, with the caveat that the additional mid-circuit measurement and the classical feedforward themselves take time and introduce their own sources of errors.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "206a7278",
      "metadata": {},
      "source": [
        "## Step 3: Execute using Qiskit primitives\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "53a3e876",
      "metadata": {},
      "source": [
        "#### Local testing mode\n",
        "\n",
        "Before submitting the jobs to the hardware, we can run a small test simulation of the dynamic circuit using the [local testing mode](/docs/guides/local-testing-mode).\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 29,
      "id": "cdfd9576",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Simulated average magnetization at trotter step = 1 at three theta values\n"
          ]
        },
        {
          "data": {
            "text/plain": [
              "array([ 0.16666667,  0.01529948, -0.14290365])"
            ]
          },
          "execution_count": 29,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "aer_sim = AerSimulator()\n",
        "pm = generate_preset_pass_manager(backend=aer_sim, optimization_level=1)\n",
        "circuit_dynamic_test.measure_all()\n",
        "isa_qc = pm.run(circuit_dynamic_test)\n",
        "with Batch(backend=aer_sim) as batch:\n",
        "    sampler = Sampler(mode=batch)\n",
        "    result = sampler.run([(isa_qc, params)]).result()\n",
        "\n",
        "print(\n",
        "    \"Simulated average magnetization at trotter step = 1 at three theta values\"\n",
        ")\n",
        "result[0].data[\"meas\"].expectation_values(obs_dynamic_test[0])"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "59a3a8a0",
      "metadata": {},
      "source": [
        "#### MPS simulation\n",
        "\n",
        "For large circuits, we can use the `matrix_product_state` (MPS) simulator, which provides an approximate result to the expectation value according to the chosen bond dimension. We later use the MPS simulation results as the baseline to compare the results from the hardware.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 30,
      "id": "6fab4462",
      "metadata": {},
      "outputs": [],
      "source": [
        "# The MPS simulation below took approximately 7 minutes to run on a\n",
        "# laptop with Apple M1 chip\n",
        "\n",
        "mps_backend = AerSimulator(\n",
        "    method=\"matrix_product_state\",\n",
        "    matrix_product_state_truncation_threshold=1e-5,\n",
        "    matrix_product_state_max_bond_dimension=100,\n",
        ")\n",
        "mps_sampler = Aer_Sampler.from_backend(mps_backend)\n",
        "\n",
        "shots = 4096\n",
        "\n",
        "data_sim = []\n",
        "for j in range(points):\n",
        "    circ_list = [\n",
        "        circ.assign_parameters([params[j]]) for circ in circuits_unitary\n",
        "    ]\n",
        "\n",
        "    mps_job = mps_sampler.run(circ_list, shots=shots)\n",
        "    result = mps_job.result()\n",
        "\n",
        "    point_data = [\n",
        "        result[d].data[\"meas\"].expectation_values(observables_unitary)\n",
        "        for d in depths\n",
        "    ]\n",
        "\n",
        "    data_sim.append(point_data)  # data at one theta value\n",
        "\n",
        "data_sim = np.array(data_sim)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "ebd28d97",
      "metadata": {},
      "source": [
        "With the circuits and observables prepared, we now execute them on hardware using the Sampler primitive.\n",
        "\n",
        "Here we submit three jobs for `unitary_pubs`, `dynamic_pubs`, and `dynamic_pubs_dd`. Each is a list of parametrized circuits corresponding to nine different Trotter steps with three different $\\theta$ parameters.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 31,
      "id": "76b5e07e",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "unitary: d96s4b52su3c739hakrg\n",
            "unitary_dd: d96s4bt2su3c739haksg\n",
            "dynamic: d96s4c0tcv6s73dk55mg\n",
            "dynamic_dd: d96s4ckqp3as739qvid0\n",
            "dynamic_meas2: d96s4csqp3as739qvie0\n",
            "dynamic_dd_meas2: d96s4daf47jc73a5v8s0\n"
          ]
        }
      ],
      "source": [
        "shots = 10000\n",
        "\n",
        "with Batch(backend=backend) as batch:\n",
        "    sampler = Sampler(mode=batch)\n",
        "\n",
        "    sampler.options.experimental = {\n",
        "        \"execution\": {\n",
        "            \"scheduler_timing\": True\n",
        "        },  # set to True to retrieve circuit timing info\n",
        "    }\n",
        "\n",
        "    job_unitary = sampler.run(unitary_pubs, shots=shots)\n",
        "    print(f\"unitary: {job_unitary.job_id()}\")\n",
        "\n",
        "    job_unitary_dd = sampler.run(unitary_pubs_dd, shots=shots)\n",
        "    print(f\"unitary_dd: {job_unitary_dd.job_id()}\")\n",
        "\n",
        "    job_dynamic = sampler.run(dynamic_pubs, shots=shots)\n",
        "    print(f\"dynamic: {job_dynamic.job_id()}\")\n",
        "\n",
        "    job_dynamic_dd = sampler.run(dynamic_pubs_dd, shots=shots)\n",
        "    print(f\"dynamic_dd: {job_dynamic_dd.job_id()}\")\n",
        "\n",
        "    job_dynamic_meas2 = sampler.run(dynamic_pubs_meas2, shots=shots)\n",
        "    print(f\"dynamic_meas2: {job_dynamic_meas2.job_id()}\")\n",
        "\n",
        "    job_dynamic_dd_meas2 = sampler.run(dynamic_pubs_dd_meas2, shots=shots)\n",
        "    print(f\"dynamic_dd_meas2: {job_dynamic_dd_meas2.job_id()}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "273fa3a0",
      "metadata": {},
      "source": [
        "## Step 4: Post-process and return results in desired classical format\n",
        "\n",
        "After the jobs have completed, we can retrieve the circuit duration from the job results metadata and visualize the circuit schedule information. To read more about visualizing a circuit's scheduling information, refer to [this page](/docs/guides/qiskit-runtime-circuit-timing).\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 32,
      "id": "bc16418f",
      "metadata": {},
      "outputs": [],
      "source": [
        "# Circuit durations is reported in the unit of `dt`\n",
        "# which can be retrieved from `Backend` object\n",
        "unitary_durations = [\n",
        "    job_unitary.result()[i].metadata[\"compilation\"][\"scheduler_timing\"][\n",
        "        \"circuit_duration\"\n",
        "    ]\n",
        "    for i in depths\n",
        "]\n",
        "\n",
        "dynamic_durations = [\n",
        "    job_dynamic.result()[i].metadata[\"compilation\"][\"scheduler_timing\"][\n",
        "        \"circuit_duration\"\n",
        "    ]\n",
        "    for i in depths\n",
        "]\n",
        "\n",
        "dynamic_durations_meas2 = [\n",
        "    job_dynamic_meas2.result()[i].metadata[\"compilation\"][\"scheduler_timing\"][\n",
        "        \"circuit_duration\"\n",
        "    ]\n",
        "    for i in depths\n",
        "]\n",
        "\n",
        "result_dd = job_dynamic_dd.result()[1]\n",
        "circuit_schedule_dd = result_dd.metadata[\"compilation\"][\"scheduler_timing\"][\n",
        "    \"timing\"\n",
        "]\n",
        "\n",
        "# to visualize the circuit schedule, one can show the figure below\n",
        "fig_dd = draw_circuit_schedule_timing(\n",
        "    circuit_schedule=circuit_schedule_dd,\n",
        "    included_channels=None,\n",
        "    filter_readout_channels=False,\n",
        "    filter_barriers=False,\n",
        "    width=1000,\n",
        ")\n",
        "\n",
        "# Save to a file since the figure is large\n",
        "fig_dd.write_html(\"scheduler_timing_dd.html\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "bee08b77",
      "metadata": {},
      "source": [
        "We plot the circuit durations for unitary circuits and the dynamic circuits. From the plot below, we can see that, despite the time needed for mid-circuit measurements and classical operations, dynamic circuit implementation with `measure_2` results in comparable circuit durations as the unitary implementation.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 33,
      "id": "639221e6",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<matplotlib.legend.Legend at 0x12bfde270>"
            ]
          },
          "execution_count": 33,
          "metadata": {},
          "output_type": "execute_result"
        },
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/639221e6-1.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "# visualize circuit durations\n",
        "\n",
        "\n",
        "def convert_dt_to_microseconds(circ_duration: List, backend_dt: float):\n",
        "    dt = backend_dt * 1e6  # dt in microseconds\n",
        "    return list(map(lambda x: x * dt, circ_duration))\n",
        "\n",
        "\n",
        "dt = backend.target.dt\n",
        "plt.plot(\n",
        "    depths,\n",
        "    convert_dt_to_microseconds(unitary_durations, dt),\n",
        "    color=\"#be95ff\",\n",
        "    linestyle=\":\",\n",
        "    label=\"unitary\",\n",
        ")\n",
        "plt.plot(\n",
        "    depths,\n",
        "    convert_dt_to_microseconds(dynamic_durations, dt),\n",
        "    color=\"#ff7eb6\",\n",
        "    linestyle=\"-.\",\n",
        "    label=\"dynamic\",\n",
        ")\n",
        "plt.plot(\n",
        "    depths,\n",
        "    convert_dt_to_microseconds(dynamic_durations_meas2, dt),\n",
        "    color=\"#ff7eb6\",\n",
        "    linestyle=\"-.\",\n",
        "    marker=\"s\",\n",
        "    mfc=\"none\",\n",
        "    label=\"dynamic w/ meas2\",\n",
        ")\n",
        "\n",
        "plt.xlabel(\"Trotter steps\")\n",
        "plt.ylabel(r\"Circuit durations in $\\mu$s\")\n",
        "plt.legend()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "3c2bd06f",
      "metadata": {},
      "source": [
        "After the jobs have completed, we retrieve the data below and compute the average magnetization estimated by the observables `observables_unitary` or `observables_dynamic` we constructed earlier.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 34,
      "id": "f4947049",
      "metadata": {},
      "outputs": [],
      "source": [
        "runs = {\n",
        "    \"unitary\": (\n",
        "        job_unitary,\n",
        "        [observables_unitary] * len(circuits_unitary),\n",
        "    ),\n",
        "    \"unitary_dd\": (\n",
        "        job_unitary_dd,\n",
        "        [observables_unitary] * len(circuits_unitary),\n",
        "    ),\n",
        "    # Omitting Dyn w/o DD and Dynamic w/ DD plots for better readability\n",
        "    # \"dynamic\": (job_dynamic, observables_dynamic),\n",
        "    # \"dynamic_dd\": (job_dynamic_dd, observables_dynamic),\n",
        "    \"dynamic_meas2\": (job_dynamic_meas2, observables_dynamic),\n",
        "    \"dynamic_dd_meas2\": (\n",
        "        job_dynamic_dd_meas2,\n",
        "        observables_dynamic,\n",
        "    ),\n",
        "}"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 35,
      "id": "583d417e",
      "metadata": {},
      "outputs": [],
      "source": [
        "data_dict = {}\n",
        "for key, (job, obs) in runs.items():\n",
        "    data = []\n",
        "    for i in range(points):\n",
        "        data.append(\n",
        "            [\n",
        "                job.result()[ind].data[\"meas\"].expectation_values(obs[ind])[i]\n",
        "                for ind in depths\n",
        "            ]\n",
        "        )\n",
        "    data_dict[key] = data"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "507073d9",
      "metadata": {},
      "source": [
        "Below we plot the spin magnetization as a function of the Trotter steps at different $\\theta$ values, corresponding to different strengths of the local magnetic field. We plot both the pre-computed MPS simulation results for the unitary ideal circuits, together with the experimental results from the following:\n",
        "\n",
        "1. running the unitary circuits with DD\n",
        "2. running the dynamic circuits with DD and `MidCircuitMeasure`\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 36,
      "id": "662239cf",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/dc-hex-ising/extracted-outputs/662239cf-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "plt.figure(figsize=(10, 6))\n",
        "\n",
        "colors = [\"#0f62fe\", \"#be95ff\", \"#ff7eb6\"]\n",
        "for i in range(points):\n",
        "    plt.plot(\n",
        "        depths,\n",
        "        data_sim[i],\n",
        "        color=colors[i],\n",
        "        linestyle=\"solid\",\n",
        "        label=f\"θ={pi_check(i*max_angle/(points-1))} (MPS)\",\n",
        "    )\n",
        "    # plt.plot(\n",
        "    #     depths,\n",
        "    #     data_dict[\"unitary\"][i],\n",
        "    #     color=colors[i],\n",
        "    #     linestyle=\":\",\n",
        "    #     label=f\"θ={pi_check(i*max_angle/(points-1))} (Unitary)\",\n",
        "    # )\n",
        "\n",
        "    plt.plot(\n",
        "        depths,\n",
        "        data_dict[\"unitary_dd\"][i],\n",
        "        color=colors[i],\n",
        "        marker=\"o\",\n",
        "        mfc=\"none\",\n",
        "        linestyle=\":\",\n",
        "        label=f\"θ={pi_check(i*max_angle/(points-1))} (Unitary w/DD)\",\n",
        "    )\n",
        "\n",
        "    # Omitting Dyn w/o DD and Dynamic w/ DD plots for better readability\n",
        "    # plt.plot(\n",
        "    #     depths,\n",
        "    #     data_dict[\"dynamic\"][i],\n",
        "    #     color=colors[i],\n",
        "    #     linestyle=\"-.\",\n",
        "    #     label=f\"θ={pi_check(i*max_angle/(points-1))} (Dyn w/o DD)\",\n",
        "    # )\n",
        "    # plt.plot(\n",
        "    #     depths,\n",
        "    #     data_dict[\"dynamic_dd\"][i],\n",
        "    #     marker=\"D\",\n",
        "    #     mfc=\"none\",\n",
        "    #     color=colors[i],\n",
        "    #     linestyle=\"-.\",\n",
        "    #     label=f\"θ={pi_check(i*max_angle/(points-1))} (Dynamic w/ DD)\",\n",
        "    # )\n",
        "\n",
        "    # plt.plot(\n",
        "    #     depths,\n",
        "    #     data_dict[\"dynamic_meas2\"][i],\n",
        "    #     color=colors[i],\n",
        "    #     marker=\"s\",\n",
        "    #     mfc=\"none\",\n",
        "    #     linestyle=':',\n",
        "    #     label=f\"θ={pi_check(i*max_angle/(points-1))} (Dynamic w/ MidCircuitMeas)\",\n",
        "    # )\n",
        "\n",
        "    plt.plot(\n",
        "        depths,\n",
        "        data_dict[\"dynamic_dd_meas2\"][i],\n",
        "        color=colors[i],\n",
        "        marker=\"*\",\n",
        "        markersize=8,\n",
        "        linestyle=\":\",\n",
        "        label=f\"θ={pi_check(i*max_angle/(points-1))} \"\n",
        "        f\"(Dynamic w/ DD & MidCircuitMeas)\",\n",
        "    )\n",
        "\n",
        "\n",
        "plt.xlabel(\"Trotter steps\", fontsize=16)\n",
        "plt.ylabel(\"Average magnetization\", fontsize=16)\n",
        "plt.xticks(rotation=45)\n",
        "handles, labels = plt.gca().get_legend_handles_labels()\n",
        "plt.legend(\n",
        "    handles,\n",
        "    labels,\n",
        "    loc=\"upper right\",\n",
        "    bbox_to_anchor=(1.46, 1.0),\n",
        "    shadow=True,\n",
        "    ncol=1,\n",
        ")\n",
        "plt.title(\n",
        "    f\"{hex_rows}x{hex_cols} hex ring, {num_qubits} data qubits, \"\n",
        "    f\"{len(ancilla)} ancilla qubits \\n{backend.name}: Sampler\"\n",
        ")\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "aead1034",
      "metadata": {},
      "source": [
        "When we compare the experimental results with the simulation, we see that the dynamic circuit implementation (dotted line with stars) overall has better performance than the standard unitary implementation (dotted line with circles). In summary, we present dynamic circuits as a solution for simulating Ising spin models on a honeycomb lattice, a topology that is not native to the hardware. The dynamic circuit solution allows ZZ interactions between qubits that are not nearest neighbors, with a shorter two-qubit gate depth than using SWAP gates, at the cost of introducing extra ancilla qubits and classical feedforward operations.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "52fb43e1",
      "metadata": {},
      "source": [
        "## References\n",
        "\n",
        "\\[1] Quantum computing with Qiskit, by Javadi-Abhari, A., Treinish, M., Krsulich, K., Wood, C.J., Lishman, J., Gacon, J., Martiel, S., Nation, P.D., Bishop, L.S., Cross, A.W. and Johnson, B.R., 2024. arXiv preprint [arXiv:2405.08810 (2024)](https://arxiv.org/abs/2405.08810)\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"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 5
}