{
  "cells": [
    {
      "cell_type": "markdown",
      "id": "492d7cec",
      "metadata": {},
      "source": [
        "---\n",
        "title: CHSH inequality\n",
        "description: Run an experiment on a quantum computer to demonstrate the violation of the CHSH inequality with the Estimator primitive.\n",
        "---\n",
        "\n",
        "{/* cspell:ignore zorder Tsirelson */}\n",
        "\n",
        "# CHSH inequality\n",
        "\n",
        "*Usage estimate: Two minutes on a Heron r3 processor (NOTE: This is an estimate only. Your runtime might vary.)*\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "a57d9139",
      "metadata": {},
      "source": [
        "## Learning outcomes\n",
        "\n",
        "After completing this tutorial, you can expect to understand the following information:\n",
        "\n",
        "* How to construct a parameterized Bell-state CHSH circuit and measure the four expectation values that make up the CHSH witnesses.\n",
        "* How to compute expectation values of multiple observables on a parameter sweep in a single call to the [`EstimatorV2`](/docs/api/qiskit-ibm-runtime/estimator-v2) primitive.\n",
        "* How to validate a quantum workflow on a noisy local simulator with `AerSimulator.from_backend` before submitting to hardware.\n",
        "* How to scale a CHSH experiment into a device-wide entanglement benchmark by running many independent Bell pairs in parallel on IBM Quantum® hardware.\n",
        "\n",
        "## Prerequisites\n",
        "\n",
        "It is recommended that you familiarize yourself with these topics:\n",
        "\n",
        "* [Entanglement in action](/learning/courses/basics-of-quantum-information/entanglement-in-action/chsh-game), a course lesson on Bell states and the CHSH game.\n",
        "* [`SparsePauliOp`](/docs/api/qiskit/qiskit.quantum_info.SparsePauliOp) and the Qiskit primitives [introduction](/docs/guides/primitives).\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "f956eda1",
      "metadata": {},
      "source": [
        "## Background\n",
        "\n",
        "In this tutorial, you will run an experiment on a quantum computer to demonstrate the violation of the CHSH inequality with the Estimator primitive.\n",
        "\n",
        "The CHSH inequality, named after Clauser, Horne, Shimony, and Holt, is used to experimentally test Bell's theorem (1969). The theorem asserts that local hidden-variable theories cannot account for some consequences of entanglement in quantum mechanics. Demonstrating a violation of the CHSH inequality shows that quantum mechanics is incompatible with local hidden-variable theories, an experiment that is foundational to our understanding of quantum mechanics.\n",
        "\n",
        "The 2022 Nobel Prize for Physics was awarded to Alain Aspect, John Clauser, and Anton Zeilinger in part for their pioneering work in quantum information science, and in particular, for their experiments with entangled photons demonstrating violation of Bell's inequalities.\n",
        "\n",
        "For this experiment, we will create an entangled pair on which we measure each qubit in two different bases. We will label the bases for the first qubit $A$ and $a$ and the bases for the second qubit $B$ and $b$. This allows us to compute the CHSH quantity $S_1$:\n",
        "\n",
        "$$\n",
        "S_1 = A(B-b) + a(B+b).\n",
        "$$\n",
        "\n",
        "Each observable is either $+1$ or $-1$. Clearly, one of the terms $B\\pm b$ must be $0$, and the other must be $\\pm 2$. Therefore, $S_1 = \\pm 2$. The average value of $S_1$ must satisfy the inequality:\n",
        "\n",
        "$$\n",
        "|\\langle S_1 \\rangle|\\leq 2.\n",
        "$$\n",
        "\n",
        "Expanding $S_1$ in terms of $A$, $a$, $B$, and $b$ gives:\n",
        "\n",
        "$$\n",
        "|\\langle S_1 \\rangle| = |\\langle AB \\rangle - \\langle Ab \\rangle + \\langle aB \\rangle + \\langle ab \\rangle| \\leq 2.\n",
        "$$\n",
        "\n",
        "You can define another CHSH quantity $S_2$:\n",
        "\n",
        "$$\n",
        "S_2 = A(B+b) - a(B-b),\n",
        "$$\n",
        "\n",
        "which leads to another inequality:\n",
        "\n",
        "$$\n",
        "|\\langle S_2 \\rangle| = |\\langle AB \\rangle + \\langle Ab \\rangle - \\langle aB \\rangle + \\langle ab \\rangle| \\leq 2.\n",
        "$$\n",
        "\n",
        "If quantum mechanics could be described by local hidden-variable theories, these inequalities would always hold. As demonstrated in this tutorial, they can be violated on a quantum computer, so quantum mechanics is not compatible with local hidden-variable theories.\n",
        "\n",
        "We create the entangled pair by preparing the Bell state $|\\Phi^+\\rangle = \\frac{|00\\rangle + |11\\rangle}{\\sqrt{2}}$. Using the Estimator primitive, we obtain the expectation values $\\langle AB \\rangle, \\langle Ab \\rangle, \\langle aB \\rangle$, and $\\langle ab \\rangle$ directly, without reconstructing them from raw counts. We measure the second qubit in the $Z$ and $X$ bases. The first qubit is measured in orthogonal bases as well, but with a rotation angle $\\theta$ that we sweep between $0$ and $2\\pi$. The Estimator primitive evaluates this parameter sweep in a single [primitive unified bloc (PUB)](/docs/guides/primitive-input-output).\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c7346a67",
      "metadata": {},
      "source": [
        "## Requirements\n",
        "\n",
        "Before starting this tutorial, be sure you have the following installed:\n",
        "\n",
        "* Qiskit SDK v2.0 or later, with [visualization](/docs/api/qiskit/visualization) support\n",
        "* Qiskit Runtime v0.40 or later (`pip install qiskit-ibm-runtime`)\n",
        "* Qiskit Aer v0.17 or later (`pip install qiskit-aer`)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "d5ffeb42",
      "metadata": {},
      "source": [
        "## Setup\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "id": "323b27db",
      "metadata": {
        "execution": {
          "iopub.execute_input": "2026-05-19T22:38:49.177326Z",
          "iopub.status.busy": "2026-05-19T22:38:49.177203Z",
          "iopub.status.idle": "2026-05-19T22:38:49.804694Z",
          "shell.execute_reply": "2026-05-19T22:38:49.804154Z"
        }
      },
      "outputs": [],
      "source": [
        "# General\n",
        "import numpy as np\n",
        "\n",
        "# Qiskit imports\n",
        "from qiskit import QuantumCircuit\n",
        "from qiskit.circuit import Parameter\n",
        "from qiskit.quantum_info import SparsePauliOp\n",
        "from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager\n",
        "\n",
        "# Qiskit Runtime imports\n",
        "from qiskit_ibm_runtime import QiskitRuntimeService\n",
        "from qiskit_ibm_runtime import EstimatorV2 as Estimator\n",
        "\n",
        "# Qiskit Aer for local noisy simulation\n",
        "from qiskit_aer import AerSimulator\n",
        "\n",
        "# Plotting routines\n",
        "import matplotlib.pyplot as plt\n",
        "import matplotlib.ticker as tck"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "fd284a0f",
      "metadata": {
        "execution": {
          "iopub.execute_input": "2026-05-19T22:38:49.806301Z",
          "iopub.status.busy": "2026-05-19T22:38:49.806113Z",
          "iopub.status.idle": "2026-05-19T22:38:53.208796Z",
          "shell.execute_reply": "2026-05-19T22:38:53.208025Z"
        }
      },
      "outputs": [
        {
          "data": {
            "text/plain": [
              "'ibm_pittsburgh'"
            ]
          },
          "execution_count": 2,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Select an IBM Quantum backend.\n",
        "service = QiskitRuntimeService()\n",
        "backend = service.least_busy(\n",
        "    min_num_qubits=127, operational=True, simulator=False\n",
        ")\n",
        "backend.name"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "ca035271",
      "metadata": {},
      "source": [
        "## Small-scale simulator example\n",
        "\n",
        "Before submitting a hardware job, we validate the entire workflow on a local noisy simulator. We use `AerSimulator.from_backend(backend)` to build a simulator that inherits the noise model and coupling map of the backend you selected, so the simulator response is qualitatively similar to what we expect from hardware.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "0539a796",
      "metadata": {},
      "source": [
        "### Step 1: Map classical inputs to a quantum problem\n",
        "\n",
        "We write the CHSH circuit with a single parameter $\\theta$, which sweeps the measurement basis of the first qubit. The [`Estimator`](/docs/api/qiskit-ibm-runtime/estimator-v2) primitive simplifies the analysis: it returns expectation values of observables directly, and it can evaluate a parameterized circuit at many parameter values in a single call.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "id": "c3f57d25",
      "metadata": {
        "execution": {
          "iopub.execute_input": "2026-05-19T22:38:53.210931Z",
          "iopub.status.busy": "2026-05-19T22:38:53.210786Z",
          "iopub.status.idle": "2026-05-19T22:38:53.349189Z",
          "shell.execute_reply": "2026-05-19T22:38:53.348782Z"
        }
      },
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/chsh-inequality/extracted-outputs/c3f57d25-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 3,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "theta = Parameter(r\"$\\theta$\")\n",
        "\n",
        "chsh_circuit = QuantumCircuit(2)\n",
        "chsh_circuit.h(0)\n",
        "chsh_circuit.cx(0, 1)\n",
        "chsh_circuit.ry(theta, 0)\n",
        "chsh_circuit.draw(output=\"mpl\", idle_wires=False, style=\"iqp\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "af999dde",
      "metadata": {},
      "source": [
        "Next, we create a list of 21 phase values from $0$ to $2\\pi$ at which to evaluate the parameterized circuit ($0$, $0.1\\pi$, $0.2\\pi$, ..., $1.9\\pi$, $2\\pi$).\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 4,
      "id": "9d587e2f",
      "metadata": {
        "execution": {
          "iopub.execute_input": "2026-05-19T22:38:53.350400Z",
          "iopub.status.busy": "2026-05-19T22:38:53.350337Z",
          "iopub.status.idle": "2026-05-19T22:38:53.352045Z",
          "shell.execute_reply": "2026-05-19T22:38:53.351729Z"
        }
      },
      "outputs": [],
      "source": [
        "number_of_phases = 21\n",
        "phases = np.linspace(0, 2 * np.pi, number_of_phases)\n",
        "# Phases need to be expressed as a list of lists for the Estimator PUB\n",
        "individual_phases = [[ph] for ph in phases]"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "af917afb",
      "metadata": {},
      "source": [
        "Finally we define the observables. The first qubit is measured along axes rotated by $\\theta$; the second qubit is measured in $Z$ and $X$. With those choices, the four CHSH correlators map to the Pauli operators $ZZ$, $ZX$, $XZ$, and $XX$:\n",
        "\n",
        "$$\n",
        "\\langle S_1 \\rangle = \\langle ZZ \\rangle - \\langle ZX \\rangle + \\langle XZ \\rangle + \\langle XX \\rangle,\n",
        "$$\n",
        "\n",
        "$$\n",
        "\\langle S_2 \\rangle = \\langle ZZ \\rangle + \\langle ZX \\rangle - \\langle XZ \\rangle + \\langle XX \\rangle.\n",
        "$$\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 5,
      "id": "21ae6238",
      "metadata": {
        "execution": {
          "iopub.execute_input": "2026-05-19T22:38:53.353007Z",
          "iopub.status.busy": "2026-05-19T22:38:53.352948Z",
          "iopub.status.idle": "2026-05-19T22:38:53.355217Z",
          "shell.execute_reply": "2026-05-19T22:38:53.354857Z"
        }
      },
      "outputs": [],
      "source": [
        "# <S_1> = <ZZ> - <ZX> + <XZ> + <XX>\n",
        "observable1 = SparsePauliOp.from_list(\n",
        "    [(\"ZZ\", 1), (\"ZX\", -1), (\"XZ\", 1), (\"XX\", 1)]\n",
        ")\n",
        "\n",
        "# <S_2> = <ZZ> + <ZX> - <XZ> + <XX>\n",
        "observable2 = SparsePauliOp.from_list(\n",
        "    [(\"ZZ\", 1), (\"ZX\", 1), (\"XZ\", -1), (\"XX\", 1)]\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "2d571aa2",
      "metadata": {},
      "source": [
        "### Step 2: Optimize problem for quantum hardware execution\n",
        "\n",
        "V2 primitives only accept circuits and observables that conform to the instructions and connectivity supported by the target system (instruction set architecture, or ISA, circuits and observables). We build the `AerSimulator` from the backend and transpile against the simulator's target so the same pass manager is exercised end-to-end.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 6,
      "id": "3e139a89",
      "metadata": {
        "execution": {
          "iopub.execute_input": "2026-05-19T22:38:53.356197Z",
          "iopub.status.busy": "2026-05-19T22:38:53.356140Z",
          "iopub.status.idle": "2026-05-19T22:38:56.713188Z",
          "shell.execute_reply": "2026-05-19T22:38:56.712683Z"
        }
      },
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/chsh-inequality/extracted-outputs/3e139a89-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 6,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Build a noisy simulator from the ibm_pittsburgh backend\n",
        "aer_sim = AerSimulator.from_backend(backend)\n",
        "\n",
        "pm = generate_preset_pass_manager(target=aer_sim.target, optimization_level=3)\n",
        "chsh_isa_circuit = pm.run(chsh_circuit)\n",
        "chsh_isa_circuit.draw(output=\"mpl\", idle_wires=False, style=\"iqp\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "2817b069",
      "metadata": {},
      "source": [
        "We also transform the observables to match the transpiled circuit's qubit layout using `SparsePauliOp.apply_layout`.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 7,
      "id": "d1faa7ab",
      "metadata": {
        "execution": {
          "iopub.execute_input": "2026-05-19T22:38:56.714405Z",
          "iopub.status.busy": "2026-05-19T22:38:56.714288Z",
          "iopub.status.idle": "2026-05-19T22:38:56.716319Z",
          "shell.execute_reply": "2026-05-19T22:38:56.715992Z"
        }
      },
      "outputs": [],
      "source": [
        "isa_observable1 = observable1.apply_layout(layout=chsh_isa_circuit.layout)\n",
        "isa_observable2 = observable2.apply_layout(layout=chsh_isa_circuit.layout)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "e734572e",
      "metadata": {},
      "source": [
        "### Step 3: Execute using Qiskit primitives\n",
        "\n",
        "Run the parameter sweep with `EstimatorV2` in `aer_sim` mode. The Estimator `run()` method takes an iterable of PUBs. Each PUB has the format `(circuit, observables, parameter_values, precision)`. We pass both observables together so they share the same parameter sweep.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 8,
      "id": "7f954aeb",
      "metadata": {
        "execution": {
          "iopub.execute_input": "2026-05-19T22:38:56.717549Z",
          "iopub.status.busy": "2026-05-19T22:38:56.717473Z",
          "iopub.status.idle": "2026-05-19T22:38:58.217539Z",
          "shell.execute_reply": "2026-05-19T22:38:58.217000Z"
        }
      },
      "outputs": [],
      "source": [
        "# Use the AerSimulator-backed Estimator to validate the workflow locally\n",
        "estimator_sim = Estimator(mode=aer_sim)\n",
        "\n",
        "pub = (\n",
        "    chsh_isa_circuit,  # ISA circuit\n",
        "    [[isa_observable1], [isa_observable2]],  # ISA observables\n",
        "    individual_phases,  # Parameter values\n",
        ")\n",
        "\n",
        "sim_result = estimator_sim.run(pubs=[pub]).result()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "9c5230fc",
      "metadata": {},
      "source": [
        "### Step 4: Post-process and return result in desired classical format\n",
        "\n",
        "The Estimator returns expectation values for both observables. We plot them against $\\theta$ together with the classical bound ($\\pm 2$) and the Tsirelson bound ($\\pm 2\\sqrt{2}$). The shaded grey regions mark the gap between the two. Points that lie inside these bands violate the CHSH inequality.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 9,
      "id": "c8fd5140",
      "metadata": {
        "execution": {
          "iopub.execute_input": "2026-05-19T22:38:58.219407Z",
          "iopub.status.busy": "2026-05-19T22:38:58.219306Z",
          "iopub.status.idle": "2026-05-19T22:38:58.285202Z",
          "shell.execute_reply": "2026-05-19T22:38:58.284609Z"
        }
      },
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/chsh-inequality/extracted-outputs/c8fd5140-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "chsh1_sim = sim_result[0].data.evs[0]\n",
        "chsh2_sim = sim_result[0].data.evs[1]\n",
        "\n",
        "\n",
        "def plot_chsh(phases, chsh1, chsh2, title):\n",
        "    fig, ax = plt.subplots(figsize=(10, 6))\n",
        "\n",
        "    ax.plot(\n",
        "        phases / np.pi, chsh1, \"o-\", label=r\"$\\langle S_1 \\rangle$\", zorder=3\n",
        "    )\n",
        "    ax.plot(\n",
        "        phases / np.pi, chsh2, \"o-\", label=r\"$\\langle S_2 \\rangle$\", zorder=3\n",
        "    )\n",
        "\n",
        "    # classical bound +-2\n",
        "    ax.axhline(y=2, color=\"0.9\", linestyle=\"--\")\n",
        "    ax.axhline(y=-2, color=\"0.9\", linestyle=\"--\")\n",
        "\n",
        "    # quantum bound, +-2*sqrt(2)\n",
        "    ax.axhline(y=np.sqrt(2) * 2, color=\"0.9\", linestyle=\"-.\")\n",
        "    ax.axhline(y=-np.sqrt(2) * 2, color=\"0.9\", linestyle=\"-.\")\n",
        "    ax.fill_between(phases / np.pi, 2, 2 * np.sqrt(2), color=\"0.6\", alpha=0.7)\n",
        "    ax.fill_between(\n",
        "        phases / np.pi, -2, -2 * np.sqrt(2), color=\"0.6\", alpha=0.7\n",
        "    )\n",
        "\n",
        "    ax.xaxis.set_major_formatter(tck.FormatStrFormatter(\"%g $\\\\pi$\"))\n",
        "    ax.xaxis.set_major_locator(tck.MultipleLocator(base=0.5))\n",
        "\n",
        "    ax.set_xlabel(r\"$\\theta$\")\n",
        "    ax.set_ylabel(\"CHSH witness\")\n",
        "    ax.set_title(title)\n",
        "    ax.legend()\n",
        "    plt.show()\n",
        "\n",
        "\n",
        "plot_chsh(\n",
        "    phases,\n",
        "    chsh1_sim,\n",
        "    chsh2_sim,\n",
        "    \"CHSH witnesses from AerSimulator (ibm_pittsburgh noise model)\",\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "e9744b84",
      "metadata": {},
      "source": [
        "The simulator's CHSH witnesses already exceed the classical bound of $\\pm 2$ at several values of $\\theta$, even with the backend's noise model. The peaks fall just short of the Tsirelson bound $\\pm 2\\sqrt{2}$ because of simulated device noise. With the workflow validated, we move on to actual hardware.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "747e948d",
      "metadata": {},
      "source": [
        "## Large-scale hardware example\n",
        "\n",
        "A CHSH test is intrinsically a *two-qubit* experiment, so it does not scale by making one circuit bigger. Instead, it scales by running **many tests in parallel**. Here we tile the backend with as many disjoint Bell pairs as its connectivity allows (a *matching* of the coupling map) and run an independent CHSH sub-circuit on every pair, all in a single job.\n",
        "\n",
        "This turns CHSH into a **device-wide benchmark of entanglement quality**: rather than a single hand-picked pair, we test entanglement across a large fraction of the chip at once, under realistic conditions where every pair contends with its neighbors' crosstalk and parallel-gate errors. Violating the inequality on every pair *simultaneously* certifies that genuine entanglement is available everywhere on the device.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 10,
      "id": "3376bc73",
      "metadata": {
        "execution": {
          "iopub.execute_input": "2026-05-19T22:38:58.286726Z",
          "iopub.status.busy": "2026-05-19T22:38:58.286656Z",
          "iopub.status.idle": "2026-05-19T22:38:59.198086Z",
          "shell.execute_reply": "2026-05-19T22:38:59.197690Z"
        }
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Tiling ibm_pittsburgh with 64 parallel Bell pairs (128 of 156 qubits)\n",
            "Job ID: d86efd5g7okc73el0rp0\n",
            "63/64 Bell pairs violated the CHSH inequality (mean peak witness 2.75, classical bound 2)\n"
          ]
        },
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/chsh-inequality/extracted-outputs/3376bc73-1.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "# -------------------------Step 1: Map classical inputs to a quantum problem-------------------------\n",
        "# A CHSH test is bipartite, so we scale up by running one independent CHSH\n",
        "# experiment on every disjoint Bell pair the device can host. A greedy\n",
        "# matching of the coupling map gives a set of edges that share no qubits.\n",
        "num_qubits = backend.num_qubits\n",
        "used = set()\n",
        "pairs = []\n",
        "for qa, qb in backend.coupling_map.get_edges():\n",
        "    if qa not in used and qb not in used:\n",
        "        pairs.append((qa, qb))\n",
        "        used.update((qa, qb))\n",
        "num_pairs = len(pairs)\n",
        "print(\n",
        "    f\"Tiling {backend.name} with {num_pairs} parallel Bell pairs \"\n",
        "    f\"({2 * num_pairs} of {num_qubits} qubits)\"\n",
        ")\n",
        "\n",
        "# One parameterized CHSH sub-circuit per pair, all sharing the angle theta\n",
        "theta = Parameter(r\"$\\theta$\")\n",
        "chsh_circuit = QuantumCircuit(num_qubits)\n",
        "for qa, qb in pairs:\n",
        "    chsh_circuit.h(qa)\n",
        "    chsh_circuit.cx(qa, qb)\n",
        "    chsh_circuit.ry(theta, qa)\n",
        "\n",
        "# Embed the two CHSH observables onto each pair's qubits (identity elsewhere)\n",
        "obs1 = SparsePauliOp.from_list([(\"ZZ\", 1), (\"ZX\", -1), (\"XZ\", 1), (\"XX\", 1)])\n",
        "obs2 = SparsePauliOp.from_list([(\"ZZ\", 1), (\"ZX\", 1), (\"XZ\", -1), (\"XX\", 1)])\n",
        "observables = []\n",
        "for qa, qb in pairs:\n",
        "    observables.append([obs1.apply_layout([qa, qb], num_qubits)])\n",
        "    observables.append([obs2.apply_layout([qa, qb], num_qubits)])\n",
        "\n",
        "number_of_phases = 21\n",
        "phases = np.linspace(0, 2 * np.pi, number_of_phases)\n",
        "individual_phases = [[ph] for ph in phases]\n",
        "\n",
        "# -------------------------Step 2: Optimize problem for quantum hardware execution-------------------------\n",
        "pm = generate_preset_pass_manager(target=backend.target, optimization_level=3)\n",
        "chsh_isa_circuit = pm.run(chsh_circuit)\n",
        "isa_observables = [\n",
        "    [o[0].apply_layout(chsh_isa_circuit.layout)] for o in observables\n",
        "]\n",
        "\n",
        "# -------------------------Step 3: Execute using Qiskit primitives-------------------------\n",
        "estimator_hw = Estimator(mode=backend)\n",
        "estimator_hw.options.environment.job_tags = [\"TUT_CI\"]\n",
        "\n",
        "pub = (chsh_isa_circuit, isa_observables, individual_phases)\n",
        "job = estimator_hw.run(pubs=[pub])\n",
        "print(f\"Job ID: {job.job_id()}\")\n",
        "hw_result = job.result()\n",
        "\n",
        "# -------------------------Step 4: Post-process and return result in desired classical format-------------------------\n",
        "# evs has shape (2 * num_pairs, number_of_phases); rows alternate S1, S2\n",
        "evs = np.asarray(hw_result[0].data.evs)\n",
        "chsh1_all = evs[0::2]\n",
        "chsh2_all = evs[1::2]\n",
        "\n",
        "# A pair \"violates\" CHSH if its strongest witness exceeds the classical bound\n",
        "peak = np.maximum(\n",
        "    np.abs(chsh1_all).max(axis=1), np.abs(chsh2_all).max(axis=1)\n",
        ")\n",
        "n_violate = int(np.sum(peak > 2))\n",
        "print(\n",
        "    f\"{n_violate}/{num_pairs} Bell pairs violated the CHSH inequality \"\n",
        "    f\"(mean peak witness {peak.mean():.2f}, classical bound 2)\"\n",
        ")\n",
        "\n",
        "fig, ax = plt.subplots(figsize=(10, 6))\n",
        "\n",
        "# Faint individual per-pair curves\n",
        "for row in chsh1_all:\n",
        "    ax.plot(phases / np.pi, row, color=\"#1f77b4\", alpha=0.2, lw=1)\n",
        "for row in chsh2_all:\n",
        "    ax.plot(phases / np.pi, row, color=\"#ff7f0e\", alpha=0.2, lw=1)\n",
        "\n",
        "# Bold mean curves across all pairs\n",
        "ax.plot(\n",
        "    phases / np.pi,\n",
        "    chsh1_all.mean(axis=0),\n",
        "    color=\"#1f77b4\",\n",
        "    lw=2.5,\n",
        "    label=r\"$\\langle S_1 \\rangle$ (mean)\",\n",
        ")\n",
        "ax.plot(\n",
        "    phases / np.pi,\n",
        "    chsh2_all.mean(axis=0),\n",
        "    color=\"#ff7f0e\",\n",
        "    lw=2.5,\n",
        "    label=r\"$\\langle S_2 \\rangle$ (mean)\",\n",
        ")\n",
        "\n",
        "# classical bound +-2 and Tsirelson bound +-2*sqrt(2)\n",
        "ax.axhline(y=2, color=\"0.9\", linestyle=\"--\")\n",
        "ax.axhline(y=-2, color=\"0.9\", linestyle=\"--\")\n",
        "ax.axhline(y=np.sqrt(2) * 2, color=\"0.9\", linestyle=\"-.\")\n",
        "ax.axhline(y=-np.sqrt(2) * 2, color=\"0.9\", linestyle=\"-.\")\n",
        "ax.fill_between(phases / np.pi, 2, 2 * np.sqrt(2), color=\"0.6\", alpha=0.7)\n",
        "ax.fill_between(phases / np.pi, -2, -2 * np.sqrt(2), color=\"0.6\", alpha=0.7)\n",
        "\n",
        "ax.xaxis.set_major_formatter(tck.FormatStrFormatter(\"%g $\\\\pi$\"))\n",
        "ax.xaxis.set_major_locator(tck.MultipleLocator(base=0.5))\n",
        "ax.set_xlabel(r\"$\\theta$\")\n",
        "ax.set_ylabel(\"CHSH witness\")\n",
        "ax.set_title(\n",
        "    f\"CHSH witnesses for {num_pairs} parallel Bell pairs on {backend.name}\"\n",
        ")\n",
        "ax.legend()\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "3c32f01c",
      "metadata": {},
      "source": [
        "The faint curves are the individual Bell pairs and the bold curves are their mean across the device. Every pair traces the same sinusoid predicted by quantum mechanics, and the spread between the faint curves reflects the variation in noise from pair to pair. Wherever a curve enters the grey bands, it has crossed the classical bound of $\\pm 2$, and the printed summary confirms that essentially every pair violates the CHSH inequality at the same time.\n",
        "\n",
        "The peaks fall short of the Tsirelson bound $\\pm 2\\sqrt{2}$ because of device noise, but the conclusion is unambiguous: the backend sustains genuine entanglement across the whole chip simultaneously, not just on a single hand-picked pair. This is the sense in which the CHSH experiment \"scales\": not as one larger circuit, but as a parallel benchmark that certifies entanglement everywhere at once.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "86e635a3",
      "metadata": {},
      "source": [
        "## Next steps\n",
        "\n",
        "<Admonition type=\"tip\" title=\"Recommendations\">\n",
        "  If you found this work interesting, you might be interested in the following material:\n",
        "\n",
        "  * [Entanglement in action](/learning/courses/basics-of-quantum-information/entanglement-in-action/chsh-game): a course lesson by John Watrous on Bell states and the CHSH game.\n",
        "  * [Get started with the Estimator primitive](/docs/guides/primitives): a guide on PUBs and parameter sweeps.\n",
        "  * [Real-time benchmarking for qubit selection](/docs/tutorials/real-time-benchmarking-for-qubit-selection): another way to characterize qubit and entanglement quality across a device.\n",
        "  * [`SparsePauliOp` API reference](/docs/api/qiskit/qiskit.quantum_info.SparsePauliOp).\n",
        "</Admonition>\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
}