{
  "cells": [
    {
      "cell_type": "markdown",
      "id": "2114652b-8f1f-4ba1-817b-48e98e3c6053",
      "metadata": {},
      "source": [
        "---\n",
        "title: Pauli correlation encoding to reduce max-cut requirements\n",
        "description: Use Pauli correlation encoding to encode optimization problems into qubits with greater efficiency for quantum computation.\n",
        "---\n",
        "\n",
        "{/* cspell:ignore lbrack setminus coloneqq rbrack binom rhobeg nfev PCEFQ */}\n",
        "\n",
        "# Pauli correlation encoding to reduce max-cut requirements\n",
        "\n",
        "*Usage estimate: 35 minutes on an Eagle r3 processor (NOTE: This is an estimate only. Your runtime might vary.)*\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "7239458e-d833-490e-8462-eaf2b5a115d4",
      "metadata": {},
      "source": [
        "## Learning outcomes\n",
        "\n",
        "After going through this tutorial, users should expect the following outcomes:\n",
        "\n",
        "* Understand the theoretical principles behind Pauli Correlation Encoding (PCE), including how multi‑body Pauli strings enable polynomial compression of classical optimization problems.\n",
        "* Implement PCE in practice to encode and solve large‑scale optimization tasks on near‑term quantum hardware.\n",
        "\n",
        "## Prerequisites\n",
        "\n",
        "We recommend familiarity with the following topics before going through this tutorial:\n",
        "\n",
        "* [Variational quantum algorithms](/learning/courses/variational-algorithm-design)\n",
        "* [QAOA and max-cut](/docs/tutorials/quantum-approximate-optimization-algorithm)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "8a343d18-ab46-431f-899d-0664c2f99cc0",
      "metadata": {},
      "source": [
        "## Background\n",
        "\n",
        "This tutorial presents *Pauli Correlation Encoding* (PCE) [\\[1\\]](#references), an approach designed to encode optimization problems into qubits with greater efficiency for quantum computation. PCE maps classical variables in optimization problems to multi-body Pauli-matrix correlations, resulting in a polynomial compression of the problem's space requirements. By employing PCE, the number of qubits needed for encoding is reduced, making it particularly advantageous for near-term quantum devices with limited qubit resources. Furthermore, it is analytically demonstrated that PCE inherently mitigates barren plateaus, offering super-polynomial resilience against this phenomenon. This built-in feature enables unprecedented performance in quantum optimization solvers.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c133831c-1eff-4913-aff6-8c0d82df9d61",
      "metadata": {},
      "source": [
        "### Overview\n",
        "\n",
        "The PCE approach consists of three main steps, as illustrated in Figure 1 from [\\[1\\]](#references) in below:\n",
        "\n",
        "1. Encoding the optimization problem into a Pauli correlation space.\n",
        "2. Solving the problem using a quantum-classical optimization solver.\n",
        "3. Decoding the solution back to the original optimization space.\n",
        "   The PCE approach is adaptable to any quantum optimization solver capable of processing Pauli correlation matrices.\n",
        "\n"
      ]
    },
    {
      "attachments": {},
      "cell_type": "markdown",
      "id": "40a636a3-c0c0-42ba-a68a-2abcb8a7187b",
      "metadata": {},
      "source": [
        "![Overview of PCE.](https://quantum.cloud.ibm.com/docs/images/tutorials/solving-maxcut-with-reduced-qubit-requirements-using-pauli-correlation-encoding/af2cb835-88db-4a3d-9c86-51424b1a4bd3.avif)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "1f28c0b3-b6dd-4627-90e6-db70b9114cd0",
      "metadata": {},
      "source": [
        "In Figure 1 from [\\[1\\]](#references), the [max-cut](/docs/tutorials/quantum-approximate-optimization-algorithm) problem is used as an example to illustrate the PCE approach. The max-cut problem with $m=9$ nodes is encoded into a Pauli correlation space, representing the optimization problem as a correlation matrix — specifically, two-body Pauli-matrix correlations across $n=3$ qubits $(Q_1, Q_2, Q_3)$. Node colors indicate the Pauli string used for each encoded node.\n",
        "For example, node 1, which corresponds to binary variable $x_1$, is encoded by the expectation value of $Z_1 \\otimes Z_2 \\otimes I_3$, while $x_8$ is encoded by $I_1 \\otimes Y_2 \\otimes Y_3$.\n",
        "This corresponds to compressing the problem's $m$ variables into $ n = O(m^{1/2})$ qubits. More broadly, $k $-body correlations enable polynomial compressions of order $k$, with $k>1$. The chosen Pauli set comprises three subsets of mutually-commuting Pauli strings, allowing all $m$ correlations to be experimentally estimated with only three measurement settings.\n",
        "\n",
        "A loss function $\\mathcal{L}$ of Pauli expectation values that imitates the original max-cut objective function is constructed. The loss function is then optimized using a quantum-classical optimization solver, such as the  [Variational Quantum Eigensolver (VQE)](/learning/courses/quantum-diagonalization-algorithms/vqe).\n",
        "\n",
        "Once the optimization is complete, the solution is decoded back to the original optimization space, yielding the optimal max-cut solution.\n",
        "\n",
        "## Requirements\n",
        "\n",
        "Before starting this tutorial, be sure you have the following installed:\n",
        "\n",
        "* Qiskit SDK v1.0 or later, with [visualization](/docs/api/qiskit/visualization) support\n",
        "* Qiskit Runtime v0.22 or later (`pip install qiskit-ibm-runtime`)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "adfa1c9b-0dd3-42d0-afd9-bce648cf668e",
      "metadata": {},
      "source": [
        "## Setup\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "5abb33b8-080d-4375-ac16-7788f2f1516a",
      "metadata": {},
      "outputs": [],
      "source": [
        "from itertools import combinations\n",
        "\n",
        "import numpy as np\n",
        "import rustworkx as rx\n",
        "import networkx as nx\n",
        "\n",
        "from scipy.optimize import minimize, OptimizeResult\n",
        "\n",
        "from qiskit.circuit.library import efficient_su2\n",
        "from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager\n",
        "from qiskit.quantum_info import SparsePauliOp\n",
        "from qiskit_ibm_runtime import EstimatorV2 as Estimator\n",
        "from qiskit_ibm_runtime import QiskitRuntimeService\n",
        "from qiskit_ibm_runtime import Session\n",
        "from rustworkx.visualization import mpl_draw\n",
        "from qiskit_aer import AerSimulator"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "id": "89c2999b-d309-4cf4-820e-9ef8a5cd5807",
      "metadata": {},
      "outputs": [],
      "source": [
        "def calc_cut_size(graph, partition0, partition1):\n",
        "    \"\"\"Calculate the cut size of the given partitions of the graph.\"\"\"\n",
        "\n",
        "    cut_size = 0\n",
        "    for edge0, edge1 in graph.edge_list():\n",
        "        if edge0 in partition0 and edge1 in partition1:\n",
        "            cut_size += 1\n",
        "        elif edge0 in partition1 and edge1 in partition0:\n",
        "            cut_size += 1\n",
        "    return cut_size"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "19cc2e86-204f-4ebf-b1a3-942025cc7015",
      "metadata": {},
      "source": [
        "## Small-scale simulator example\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "id": "af0c1db0-3b69-459c-8013-5149ede620c2",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "We are using the aer_simulator_from(ibm_pittsburgh)\n"
          ]
        }
      ],
      "source": [
        "service = QiskitRuntimeService()\n",
        "real_backend = service.least_busy(\n",
        "    operational=True, simulator=False, min_num_qubits=156\n",
        ")\n",
        "backend = AerSimulator.from_backend(real_backend)\n",
        "print(f\"We are using the {backend.name}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "c8084430-5386-4788-97c3-c6e4fb7cb191",
      "metadata": {},
      "source": [
        "### Step 1: Map classical inputs to a quantum problem\n",
        "\n",
        "#### The max-cut problem\n",
        "\n",
        "The max-cut problem is a combinatorial optimization problem that is defined on a graph $G = (V, E)$, where $V$ is the set of vertices and $E$ is the set of edges. The goal is to partition the vertices into two sets, $S$ and $V \\setminus S$, such that the number of edges between the two sets is maximized.\n",
        "For the detailed description of the max-cut problem, please refer to the [Quantum approximate optimization algorithm](/docs/tutorials/quantum-approximate-optimization-algorithm) tutorial.\n",
        "The max-cut problem is also used as an example in the [Advanced techniques for QAOA](/docs/tutorials/advanced-techniques-for-qaoa) tutorial.\n",
        "In those tutorials, the QAOA algorithm is used to solve the max-cut problem.\n",
        "\n",
        "#### Graph -> Hamiltonian\n",
        "\n",
        "Let us first consider a random graph with 100 nodes.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 4,
      "id": "37edb718-2bab-49d7-ad66-5f2f67d2aeff",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/pauli-correlation-encoding-for-qaoa/extracted-outputs/37edb718-2bab-49d7-ad66-5f2f67d2aeff-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "num_nodes = 100  # Number of nodes in graph\n",
        "seed = 42\n",
        "graph = rx.undirected_gnp_random_graph(num_nodes, 0.1, seed=seed)\n",
        "mpl_draw(graph)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "eb7a80dc-74ea-472b-a13f-11cb8d4c0ca9",
      "metadata": {},
      "outputs": [],
      "source": [
        "nx_graph = nx.Graph()\n",
        "nx_graph.add_nodes_from(range(num_nodes))"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 6,
      "id": "63451877-908a-4e80-9e11-23eb0d288bfc",
      "metadata": {},
      "outputs": [],
      "source": [
        "for edge in graph.edge_list():\n",
        "    nx_graph.add_edge(edge[0], edge[1])"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 7,
      "id": "515e7220-586f-4e2a-82b6-3885e3e38566",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Initial cut size: 345\n"
          ]
        }
      ],
      "source": [
        "curr_cut_size, partition = nx.approximation.one_exchange(nx_graph, seed=1)\n",
        "print(f\"Initial cut size: {curr_cut_size}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "57e3804e-61ae-45c5-9d94-3665cf01784b",
      "metadata": {},
      "source": [
        "We encode the graph with 100 nodes into two-body Pauli-matrix correlations across nine qubits (see the explanation below). The graph is represented as a correlation matrix, where each node is encoded by a Pauli string. The sign of the expectation value of the Pauli string indicates the partition of the node. For example, node 0 is encoded by a Pauli string, $\\prod_0 = I_{8} \\otimes ... I_2 \\otimes X_1 \\otimes X_0$. The sign of the expectation value of this Pauli string indicates the partition of node 0. We define a *Pauli-correlation encoding* (PCE) relative to $\\prod$ as\n",
        "\n",
        "$x_i \\coloneqq \\textit{sgn}(\\langle\\prod_i \\rangle),$\n",
        "\n",
        "where $x_i$ is the partition of node $i$ and $\\langle \\prod_i \\rangle \\coloneqq  \\langle \\psi |\\prod_i| \\psi \\rangle $ is the expectation value of the Pauli string encoding node $i$ over a quantum state $|\\psi \\rangle$.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "678e00fb-5bb9-450e-8229-0eab5676057a",
      "metadata": {},
      "source": [
        "Now, let's encode the graph into a Hamiltonian using PCE.\n",
        "We divide the nodes into three sets: $S_1$, $S_2$, and $S_3$.\n",
        "Then, we encode the nodes in each set using the Pauli strings with $X$, $Y$, and $Z$, respectively.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "4d5d48e5-ebe8-4b36-9496-288f40174a7b",
      "metadata": {},
      "source": [
        "We need to extract a relationship between the number of nodes and qubits that we will need to encode all the nodes. Using all possible permutations for the encoding yields to:\n",
        "\n",
        "$$\n",
        "m=3\\binom{n}{k}.\n",
        "$$\n",
        "\n",
        "In this example we consider $k=2$, hence,\n",
        "\n",
        "$$\n",
        "m  = \\frac{3}{2} n(n-1).\n",
        "$$\n",
        "\n",
        "Therefore, the number of qubits $n$ needed to express a certain number of nodes $m$ read as:\n",
        "\n",
        "$$\n",
        "n = \\left\\lceil \\frac{1 + \\sqrt{1 + \\tfrac{8}{3}m}}{2} \\right\\rceil.\n",
        "$$\n",
        "\n",
        "*Note that the $\\lceil \\cdot \\rceil$ symbol represents the ceiling function, which rounds any real number up to the next integer. This ensures that the number of qubits is an integer.*\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 13,
      "id": "8ea1f545-3e9f-4620-bde8-755178ad3ec9",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Number of qubits: 9\n",
            "List 1: [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32]\n",
            "List 2: [33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65]\n",
            "List 3: [66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99]\n"
          ]
        }
      ],
      "source": [
        "num_qubits = int(np.ceil((1 + np.sqrt(1 + (8 / 3) * num_nodes)) / 2))\n",
        "\n",
        "list_size = num_nodes // 3\n",
        "node_x = [i for i in range(list_size)]\n",
        "node_y = [i for i in range(list_size, 2 * list_size)]\n",
        "node_z = [i for i in range(2 * list_size, num_nodes)]\n",
        "\n",
        "print(f\"Number of qubits: {num_qubits}\")\n",
        "print(\"List 1:\", node_x)\n",
        "print(\"List 2:\", node_y)\n",
        "print(\"List 3:\", node_z)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 14,
      "id": "d2649acb-7857-4cf4-88dc-4ef381a8552f",
      "metadata": {},
      "outputs": [],
      "source": [
        "def build_pauli_correlation_encoding(pauli, node_list, n, k=2):\n",
        "    pauli_correlation_encoding = []\n",
        "    for idx, c in enumerate(combinations(range(n), k)):\n",
        "        if idx >= len(node_list):\n",
        "            break\n",
        "        paulis = [\"I\"] * n\n",
        "        paulis[c[0]], paulis[c[1]] = pauli, pauli\n",
        "        pauli_correlation_encoding.append((\"\".join(paulis)[::-1], 1))\n",
        "\n",
        "    hamiltonian = []\n",
        "    for pauli, weight in pauli_correlation_encoding:\n",
        "        hamiltonian.append(SparsePauliOp.from_list([(pauli, weight)]))\n",
        "\n",
        "    return hamiltonian\n",
        "\n",
        "\n",
        "pauli_correlation_encoding_x = build_pauli_correlation_encoding(\n",
        "    \"X\", node_x, num_qubits\n",
        ")\n",
        "pauli_correlation_encoding_y = build_pauli_correlation_encoding(\n",
        "    \"Y\", node_y, num_qubits\n",
        ")\n",
        "pauli_correlation_encoding_z = build_pauli_correlation_encoding(\n",
        "    \"Z\", node_z, num_qubits\n",
        ")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "4ce84c90-8b64-4959-a315-2380482801ad",
      "metadata": {},
      "source": [
        "### Step 2: Optimize problem for quantum hardware execution\n",
        "\n",
        "#### Quantum circuit\n",
        "\n",
        "Here, the state $|\\psi \\rangle$ is parameterized with $\\mathbf{\\theta}$, and we optimize these parameters $\\mathbf{\\theta}$ using a variational approach.\n",
        "This tutorial employs the `efficient_su2` ansatz for our variational algorithm due to its expressive capabilities and ease of implementation.\n",
        "We also use the relaxed loss function, which will be introduced later in this tutorial.\n",
        "As a result, we can address large-scale problems with fewer qubits and shallower circuit depths.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 15,
      "id": "035f6b4a-4de0-452a-b60f-7260f9e3103a",
      "metadata": {},
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<Image src=\"/docs/images/tutorials/pauli-correlation-encoding-for-qaoa/extracted-outputs/035f6b4a-4de0-452a-b60f-7260f9e3103a-0.avif\" alt=\"Output of the previous code cell\" />"
            ]
          },
          "execution_count": 15,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# Build the quantum circuit\n",
        "qc = efficient_su2(num_qubits, su2_gates=[\"ry\", \"rz\"], reps=2)\n",
        "qc.draw(\"mpl\")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 16,
      "id": "162f1384-98f5-406e-b5ae-0e12d0ad4b59",
      "metadata": {},
      "outputs": [],
      "source": [
        "# Optimize the circuit\n",
        "\n",
        "pm = generate_preset_pass_manager(optimization_level=3, backend=backend)\n",
        "qc = pm.run(qc)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "644c9317-7988-427e-979a-975c0a616f52",
      "metadata": {},
      "source": [
        "#### Loss function\n",
        "\n",
        "For the loss function $\\mathcal{L}$, we use a relaxation of the max-cut objective function as described in [\\[1\\]](#references), which is defined as $\\mathcal{V}(\\mathbf{x}) \\coloneqq \\sum_{(i, j) \\in E} W_{i, j}(1-x_i x_j)$. Here, $W_{i, j}$ denotes the weight of the edge $(i, j)$, and $x_i$ represents the partition of node $i$.\n",
        "The loss function  $\\mathcal{L}$ is given by:\n",
        "\n",
        "$\\mathcal{L}\\coloneqq \\sum_{(i, j) \\in E} W_{i, j} \\text{tanh} (\\alpha \\langle\\prod_i \\rangle) \\text{tanh} (\\alpha \\langle\\prod_j \\rangle) + \\mathcal{L}^{(\\text{reg})},$\n",
        "\n",
        "where the max-cut objective function is replaced by the smooth hyperbolic tangents of the expectation values of the Pauli strings encoding the nodes. The regularization term $\\mathcal{L}^{(\\text{reg})}$ and the rescaling factor $\\alpha$, proportional to the number of qubits, are introduced to improve the solver's performance.\n",
        "\n",
        "The regularization term is defined as:\n",
        "\n",
        "$\\mathcal{L}^{(\\text{reg})}$ is defined as $\\mathcal{L}^{(\\text{reg})} \\coloneqq \\beta \\nu \\lbrack \\frac{1}{m} \\sum_{i \\in V} \\text{tanh} (\\alpha \\langle\\prod_i \\rangle)^2 \\rbrack ^2$\n",
        "\n",
        "where $\\beta=1/2$, $\\nu = |E|/2 + (m -1) /4$, $|E|$ is the number of edges, and $m$ is the number of nodes in the graph.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "fe608e6a-08ce-493d-9b0a-eb9d6e0028ff",
      "metadata": {},
      "outputs": [],
      "source": [
        "def loss_func_estimator(x, ansatz, hamiltonian, estimator, graph):\n",
        "    \"\"\"\n",
        "    Calculates the specified loss function for the given ansatz, Hamiltonian,\n",
        "    and graph.\n",
        "\n",
        "    The expectation values of each Pauli string in the Hamiltonian are first\n",
        "    obtained by running the ansatz on the quantum backend. These\n",
        "    expectation values are then passed through the nonlinear function\n",
        "    tanh(alpha * prod_i). The loss function is\n",
        "    subsequently computed from these transformed values.\n",
        "    \"\"\"\n",
        "    job = estimator.run(\n",
        "        [\n",
        "            (ansatz, hamiltonian[0], x),\n",
        "            (ansatz, hamiltonian[1], x),\n",
        "            (ansatz, hamiltonian[2], x),\n",
        "        ]\n",
        "    )\n",
        "    result = job.result()\n",
        "\n",
        "    # calculate the loss function\n",
        "    node_exp_map = {}\n",
        "    idx = 0\n",
        "    for r in result:\n",
        "        for ev in r.data.evs:\n",
        "            node_exp_map[idx] = ev\n",
        "            idx += 1\n",
        "\n",
        "    loss = 0\n",
        "    alpha = num_qubits\n",
        "    for edge0, edge1 in graph.edge_list():\n",
        "        loss += np.tanh(alpha * node_exp_map[edge0]) * np.tanh(\n",
        "            alpha * node_exp_map[edge1]\n",
        "        )\n",
        "\n",
        "    regulation_term = 0\n",
        "    for i in range(len(graph.nodes())):\n",
        "        regulation_term += np.tanh(alpha * node_exp_map[i]) ** 2\n",
        "    regulation_term = regulation_term / len(graph.nodes())\n",
        "    regulation_term = regulation_term**2\n",
        "    beta = 1 / 2\n",
        "    v = len(graph.edges()) / 2 + (len(graph.nodes()) - 1) / 4\n",
        "    regulation_term = beta * v * regulation_term\n",
        "\n",
        "    loss = loss + regulation_term\n",
        "\n",
        "    global experiment_result\n",
        "    print(f\"Iter {len(experiment_result)}: {loss}\")\n",
        "    experiment_result.append({\"loss\": loss, \"exp_map\": node_exp_map})\n",
        "    return loss"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "49943596-7f90-4226-a900-d5890deefbd9",
      "metadata": {},
      "source": [
        "### Step 3: Execute using Qiskit primitives\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "583ba786-e2e9-4593-bf75-00723b589d78",
      "metadata": {},
      "source": [
        "In this tutorial, we set `max_iter=50` in the optimization loop for demonstration purposes. If we increase the number of iterations, we can expect better results.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 18,
      "id": "8d203dd5-8b72-4b78-a36e-c10fdef3ebc3",
      "metadata": {},
      "outputs": [],
      "source": [
        "pce = []\n",
        "pce.append(\n",
        "    [op.apply_layout(qc.layout) for op in pauli_correlation_encoding_x]\n",
        ")\n",
        "pce.append(\n",
        "    [op.apply_layout(qc.layout) for op in pauli_correlation_encoding_y]\n",
        ")\n",
        "pce.append(\n",
        "    [op.apply_layout(qc.layout) for op in pauli_correlation_encoding_z]\n",
        ")"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "8d9d6313-9bcd-4ffb-b40c-361d18c68afe",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Iter 0: 159.88755362682548\n",
            "Iter 1: 113.46202580636677\n",
            "Iter 2: 56.76494226400048\n",
            "Iter 3: 32.63357946896002\n",
            "Iter 4: 21.517837239610117\n",
            "Iter 5: 30.96034960483569\n",
            "Iter 6: 20.780475923938027\n",
            "Iter 7: 24.54251816279811\n",
            "Iter 8: 27.834486461763042\n",
            "Iter 9: 16.705460776812693\n",
            "Iter 10: 18.020587887236864\n",
            "Iter 11: 12.252379762741352\n",
            "Iter 12: 5.253885750886939\n",
            "Iter 13: 6.985984759592262\n",
            "Iter 14: 6.908717244584757\n",
            "Iter 15: 12.915466016863858\n",
            "Iter 16: 4.105776920457279\n",
            "Iter 17: 11.707504530740305\n",
            "Iter 18: 7.154360511076546\n",
            "Iter 19: 10.3890865704735\n",
            "Iter 20: 10.376147647857252\n",
            "Iter 21: 2.533430195296697\n",
            "Iter 22: 3.8612421907795462\n",
            "Iter 23: 6.103735057461906\n",
            "Iter 24: -1.1190368234312347\n",
            "Iter 25: 6.125915279494738\n",
            "Iter 26: 11.086280445482455\n",
            "Iter 27: 10.102569882302827\n",
            "Iter 28: -0.02664415648133822\n",
            "Iter 29: 7.621887727398785\n",
            "Iter 30: 5.967346615554497\n",
            "Iter 31: 3.85345716014828\n",
            "Iter 32: 4.5494846149011\n",
            "Iter 33: 10.006668112637232\n",
            "Iter 34: -3.1927138938527877\n",
            "Iter 35: 2.8829882366285116\n",
            "Iter 36: 3.3130087521654144\n",
            "Iter 37: -4.907566569808272\n",
            "Iter 38: -4.980134722109894\n",
            "Iter 39: -2.990457463896541\n",
            "Iter 40: -5.938401817344579\n",
            "Iter 41: -2.1807712386469724\n",
            "Iter 42: -1.0945774380342126\n",
            "Iter 43: -4.7548102593556685\n",
            "Iter 44: -3.8762362299208144\n",
            "Iter 45: -4.9348321021624\n",
            "Iter 46: -6.487722842864011\n",
            "Iter 47: 0.7064210113389331\n",
            "Iter 48: -2.3428323031772216\n",
            "Iter 49: -2.626032270380895\n",
            " message: Return from COBYLA because the objective function has been evaluated 50 times.\n",
            " success: False\n",
            "  status: 3\n",
            "     fun: -2.626032270380895\n",
            "       x: [ 1.375e+00  1.951e+00 ...  9.395e-01  8.948e-01]\n",
            "    nfev: 50\n"
          ]
        }
      ],
      "source": [
        "max_iter = 50\n",
        "counter = {\"i\": 0}\n",
        "last_x = {\"value\": None}\n",
        "last_fun = {\"value\": None}\n",
        "\n",
        "with Session(backend=backend) as session:\n",
        "    estimator = Estimator(mode=session)\n",
        "\n",
        "    experiment_result = []\n",
        "\n",
        "    def loss_func(x):\n",
        "        last_x[\"value\"] = x.copy()\n",
        "        if counter[\"i\"] + 1 > max_iter:\n",
        "            return last_fun[\"value\"]\n",
        "        counter[\"i\"] += 1\n",
        "        val = loss_func_estimator(\n",
        "            x, qc, [pce[0], pce[1], pce[2]], estimator, graph\n",
        "        )\n",
        "        last_fun[\"value\"] = val\n",
        "        return val\n",
        "\n",
        "    np.random.seed(seed)\n",
        "    initial_params = np.random.rand(qc.num_parameters)\n",
        "\n",
        "    result = minimize(\n",
        "        loss_func, initial_params, method=\"COBYLA\", options={\"rhobeg\": 1.0}\n",
        "    )\n",
        "\n",
        "    if counter[\"i\"] >= max_iter:\n",
        "        result = OptimizeResult(\n",
        "            message=f\"Return from COBYLA because the objective function \"\n",
        "            f\"has been evaluated {max_iter} times.\",\n",
        "            success=False,\n",
        "            status=3,\n",
        "            fun=last_fun[\"value\"],\n",
        "            x=last_x[\"value\"],\n",
        "            nfev=counter[\"i\"],\n",
        "        )\n",
        "\n",
        "print(result)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "5c9dd9f8-ed15-4008-8841-e44caf11cf99",
      "metadata": {},
      "source": [
        "### Step 4: Post-process and return result in desired classical format\n",
        "\n",
        "The partitions of the nodes are determined by evaluating the sign of the expectation values of the Pauli strings that encode the nodes.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 20,
      "id": "fa2db108-754b-4036-af98-87f5390a9c11",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "{0, 2, 3, 8, 9, 11, 12, 13, 17, 18, 20, 22, 23, 24, 25, 26, 27, 30, 35, 37, 38, 40, 43, 46, 48, 49, 50, 51, 53, 57, 61, 62, 63, 66, 67, 68, 70, 71, 74, 77, 81, 82, 83, 84, 87, 88, 94, 96, 99} {1, 4, 5, 6, 7, 10, 14, 15, 16, 19, 21, 28, 29, 31, 32, 33, 34, 36, 39, 41, 42, 44, 45, 47, 52, 54, 55, 56, 58, 59, 60, 64, 65, 69, 72, 73, 75, 76, 78, 79, 80, 85, 86, 89, 90, 91, 92, 93, 95, 97, 98}\n"
          ]
        }
      ],
      "source": [
        "# Calculate the partitions based on the final expectation values\n",
        "# If the expectation value is positive, the node belongs to partition 0 (par0)\n",
        "# Otherwise, the node belongs to partition 1 (par1)\n",
        "def get_partitions(experiment_result):\n",
        "    par0, par1 = set(), set()\n",
        "    best_index = min(\n",
        "        range(len(experiment_result)),\n",
        "        key=lambda i: experiment_result[i][\"loss\"],\n",
        "    )\n",
        "    for i in experiment_result[best_index][\"exp_map\"]:\n",
        "        if experiment_result[best_index][\"exp_map\"][i] >= 0:\n",
        "            par0.add(i)\n",
        "        else:\n",
        "            par1.add(i)\n",
        "    return par0, par1, best_index\n",
        "\n",
        "\n",
        "par0, par1, best_index = get_partitions(experiment_result)\n",
        "print(par0, par1)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "6f82e256-36f1-4ac1-badb-abdc98a26b23",
      "metadata": {},
      "source": [
        "We can calculate the cut size of the max-cut problem using the partitions of the node.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 21,
      "id": "bd85ceae-ef8b-4e21-b447-e92ca92e06eb",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Cut size: 268\n"
          ]
        }
      ],
      "source": [
        "cut_size = calc_cut_size(graph, par0, par1)\n",
        "print(f\"Cut size: {cut_size}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "ef94b209-ea19-4292-aedc-f77edf77a2eb",
      "metadata": {},
      "source": [
        "Once the training is complete, we perform one round of single-bit swap search to improve the solution as a classical post-processing step.\n",
        "In this process, we swap the partitions of two nodes and evaluate the cut size. If the cut size is improved, we keep the swap. We repeat this process for all possible pairs of nodes connected by an edge.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 22,
      "id": "b0df38ef-98d0-4d8c-bbdf-75d14e4680f7",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "[1, 0, 1, 1, 0, 0, 0, 0, 1, 1, 0, 1, 1, 1, 0, 0, 0, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 0, 0, 1, 0, 0, 0, 0, 1, 0, 1, 1, 0, 1, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1, 1, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 0, 1, 1, 0, 0, 1, 0, 0, 1, 0, 0, 0, 1, 1, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 1]\n"
          ]
        }
      ],
      "source": [
        "cur_bits = []\n",
        "\n",
        "for i in experiment_result[best_index][\"exp_map\"]:\n",
        "    if experiment_result[best_index][\"exp_map\"][i] >= 0:\n",
        "        cur_bits.append(1)\n",
        "    else:\n",
        "        cur_bits.append(0)\n",
        "print(cur_bits)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 23,
      "id": "232658b0-a86b-4cb5-aff8-c8f2423491b6",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "279 [1, 0, 1, 1, 0, 0, 0, 0, 1, 0, 0, 1, 1, 1, 0, 0, 0, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 0, 0, 1, 0, 0, 0, 0, 1, 0, 1, 1, 0, 1, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 1, 0, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 0, 1, 1, 0, 0, 1, 0, 0, 1, 0, 0, 0, 1, 1, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 1]\n"
          ]
        }
      ],
      "source": [
        "# Swap the partitions and calculate the cut size\n",
        "\n",
        "\n",
        "def swap_partitions(graph, cur_bits):\n",
        "    best_cut = 0\n",
        "    best_bits = []\n",
        "    for edge0, edge1 in graph.edge_list():\n",
        "        swapped_bits = cur_bits.copy()\n",
        "        swapped_bits[edge0], swapped_bits[edge1] = (\n",
        "            swapped_bits[edge1],\n",
        "            swapped_bits[edge0],\n",
        "        )\n",
        "\n",
        "        cur_partition = [set(), set()]\n",
        "        for i, bit in enumerate(swapped_bits):\n",
        "            if bit > 0:\n",
        "                cur_partition[0].add(i)\n",
        "            else:\n",
        "                cur_partition[1].add(i)\n",
        "        cut_size = calc_cut_size(graph, cur_partition[0], cur_partition[1])\n",
        "        if best_cut < cut_size:\n",
        "            best_cut = cut_size\n",
        "            best_bits = swapped_bits\n",
        "    return best_cut, best_bits\n",
        "\n",
        "\n",
        "best_cut, best_bits = swap_partitions(graph, cur_bits)\n",
        "print(best_cut, best_bits)"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "5142e8d0-ef6f-4abb-974e-ed25e5bfe692",
      "metadata": {},
      "source": [
        "# Large-scale hardware example\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "id": "d91c4ac2-3c5b-4afb-aa2c-42d1b50fa674",
      "metadata": {},
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "We are using 33 qubits\n",
            "We are using the ibm_pittsburgh\n",
            "Iter 0: 57399.57543902076\n",
            "Iter 1: 56458.787143794\n",
            "Iter 2: 40778.45608998947\n",
            "Iter 3: 35571.58511146131\n",
            "Iter 4: 33861.6835761173\n",
            "Iter 5: 39697.22637736274\n",
            "Iter 6: 34984.77893767163\n",
            "Iter 7: 32051.882157096858\n",
            "Iter 8: 26134.153216063707\n",
            "Iter 9: 24914.322627065787\n",
            "Iter 10: 24030.21227315425\n",
            "Iter 11: 23047.463945514\n",
            "Iter 12: 22629.42866110748\n",
            "Iter 13: 17374.859132614685\n",
            "Iter 14: 18020.11637762458\n",
            "Iter 15: 17924.7066364044\n",
            "Iter 16: 15825.1992250984\n",
            "Iter 17: 16553.346711978447\n",
            "Iter 18: 12393.565736512377\n",
            "Iter 19: 11994.021456089155\n",
            "Iter 20: 11199.994322735669\n",
            "Iter 21: 9624.895532927634\n",
            "Iter 22: 9073.811130188606\n",
            "Iter 23: 9836.721241931278\n",
            "Iter 24: 10555.925186133794\n",
            "Iter 25: 9179.1179493286\n",
            "Iter 26: 8495.394826965305\n",
            "Iter 27: 8913.688189840399\n",
            "Iter 28: 7830.448471810181\n",
            "Iter 29: 7757.430542422075\n",
            "Iter 30: 6796.187594518731\n",
            "Iter 31: 7307.985913766867\n",
            "Iter 32: 7340.225833330675\n",
            "Iter 33: 7064.731899380469\n",
            "Iter 34: 7632.270657372515\n",
            "Iter 35: 7049.154710767935\n",
            "Iter 36: 7486.118442084411\n",
            "Iter 37: 6302.12602219333\n",
            "Iter 38: 6244.934230209166\n",
            "Iter 39: 7154.9748739261395\n",
            "Iter 40: 6482.109600054041\n",
            "Iter 41: 5718.475169152395\n",
            "Iter 42: 5693.008457857462\n",
            "Iter 43: 4869.782667921923\n",
            "Iter 44: 4957.625304450959\n",
            "Iter 45: 5582.240637063214\n",
            "Iter 46: 4983.90082772116\n",
            "Iter 47: 5416.268575648202\n",
            "Iter 48: 4809.98398457807\n",
            "Iter 49: 5092.527306646118\n",
            " message: Return from COBYLA because the objective function has been evaluated 50 times.\n",
            " success: False\n",
            "  status: 3\n",
            "     fun: 5092.527306646118\n",
            "       x: [ 1.375e+00  1.951e+00 ...  7.259e-01  8.971e-01]\n",
            "    nfev: 50\n",
            "Cut size: 56152\n",
            "The best max-cut value achieved for a graph with 1500 nodes on 33 qubits is 56219\n",
            "and the specific partition we obtained is [1, 0, 0, 0, 1, 1, 0, 1, 1, 0, 0, 0, 0, 0, 1, 1, 0, 1, 0, 0, 1, 0, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 1, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 1, 0, 0, 0, 1, 0, 1, 1, 1, 1, 0, 1, 0, 0, 1, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 1, 0, 0, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 1, 1, 1, 1, 1, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 0, 0, 1, 0, 0, 0, 1, 1, 0, 1, 1, 0, 1, 0, 0, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 1, 1, 1, 0, 1, 0, 1, 0, 0, 1, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 1, 1, 0, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 0, 1, 1, 1, 0, 1, 0, 1, 1, 1, 0, 1, 0, 0, 0, 1, 1, 1, 1, 1, 1, 0, 1, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 0, 1, 0, 0, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 0, 1, 0, 1, 0, 1, 0, 1, 1, 1, 0, 0, 1, 1, 1, 0, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 1, 1, 0, 1, 0, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 0, 1, 1, 1, 0, 1, 1, 0, 0, 0, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 1, 1, 0, 0, 0, 1, 0, 0, 1, 0, 0, 1, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 1, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 0, 0, 0, 1, 1, 1, 0, 1, 1, 1, 1, 0, 0, 1, 0, 1, 0, 1, 1, 0, 0, 0, 1, 0, 1, 0, 1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 0, 1, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, 0, 1, 1, 0, 1, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 1, 0, 1, 1, 1, 1, 0, 1, 1, 0, 1, 0, 0, 0, 0, 1, 0, 0, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, 1, 0, 0, 0, 1, 1, 1, 0, 0, 1, 0, 0, 1, 0, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 1, 0, 1, 1, 0, 1, 1, 1, 1, 0, 0, 1, 1, 0, 1, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 1, 1, 1, 0, 1, 1, 0, 0, 0, 1, 0, 1, 0, 1, 1, 1, 0, 1, 1, 1, 0, 1, 1, 1, 1, 0, 0, 1, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 0, 1, 1, 1, 1, 1, 0, 0, 1, 0, 1, 0, 1, 0, 1, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 1, 1, 0, 1, 0, 0, 1, 0, 1, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 0, 0, 1, 1, 1, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0, 1, 0, 1, 0, 1, 0, 1, 1, 0, 0, 1, 1, 0, 1, 1, 0, 0, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 0, 1, 0, 0, 1, 1, 0, 1, 1, 1, 0, 1, 1, 0, 1, 1, 1, 0, 1, 0, 0, 0, 1, 0, 0, 1, 1, 1, 0, 0, 1, 0, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 0, 0, 0, 1, 0, 1, 0, 1, 1, 0, 0, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 1, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 0, 0, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, 1, 1, 0, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 0, 0, 1, 0, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 0, 1, 0, 1, 0, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 1, 1, 1, 1, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 0, 1, 0, 0, 1, 1, 0, 0, 0, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 1, 0, 1, 1, 0, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 0, 1, 1, 1, 0, 0, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]\n"
          ]
        }
      ],
      "source": [
        "# -------------------------Step 1-------------------------\n",
        "\n",
        "num_nodes = 1500  # Number of nodes in graph\n",
        "graph = rx.undirected_gnp_random_graph(num_nodes, 0.1, seed=seed)\n",
        "nx_graph = nx.Graph()\n",
        "nx_graph.add_nodes_from(range(num_nodes))\n",
        "for edge in graph.edge_list():\n",
        "    nx_graph.add_edge(edge[0], edge[1])\n",
        "\n",
        "num_qubits = int(np.ceil((1 + np.sqrt(1 + (8 / 3) * num_nodes)) / 2))\n",
        "\n",
        "list_size = num_nodes // 3\n",
        "node_x = [i for i in range(list_size)]\n",
        "node_y = [i for i in range(list_size, 2 * list_size)]\n",
        "node_z = [i for i in range(2 * list_size, num_nodes)]\n",
        "\n",
        "pauli_correlation_encoding_x = build_pauli_correlation_encoding(\n",
        "    \"X\", node_x, num_qubits\n",
        ")\n",
        "pauli_correlation_encoding_y = build_pauli_correlation_encoding(\n",
        "    \"Y\", node_y, num_qubits\n",
        ")\n",
        "pauli_correlation_encoding_z = build_pauli_correlation_encoding(\n",
        "    \"Z\", node_z, num_qubits\n",
        ")\n",
        "print(f\"We are using {num_qubits} qubits\")\n",
        "\n",
        "# -------------------------Step 2-------------------------\n",
        "backend = real_backend\n",
        "print(f\"We are using the {backend.name}\")\n",
        "qc = efficient_su2(num_qubits, [\"ry\", \"rz\"], reps=2)\n",
        "pm = generate_preset_pass_manager(optimization_level=3, backend=backend)\n",
        "qc = pm.run(qc)\n",
        "# -------------------------Step 3-------------------------\n",
        "pce = []\n",
        "pce.append(\n",
        "    [op.apply_layout(qc.layout) for op in pauli_correlation_encoding_x]\n",
        ")\n",
        "pce.append(\n",
        "    [op.apply_layout(qc.layout) for op in pauli_correlation_encoding_y]\n",
        ")\n",
        "pce.append(\n",
        "    [op.apply_layout(qc.layout) for op in pauli_correlation_encoding_z]\n",
        ")\n",
        "\n",
        "# Run the optimization using a session.\n",
        "max_iter = 50\n",
        "counter = {\"i\": 0}\n",
        "with Session(backend=backend) as session:\n",
        "    estimator = Estimator(mode=session)\n",
        "    estimator.options.environment.job_tags = [\"TUT_PCEFQ\"]\n",
        "    experiment_result = []\n",
        "\n",
        "    def loss_func(x):\n",
        "        last_x[\"value\"] = x.copy()\n",
        "        if counter[\"i\"] + 1 > max_iter:\n",
        "            return last_fun[\"value\"]\n",
        "        counter[\"i\"] += 1\n",
        "        val = loss_func_estimator(\n",
        "            x, qc, [pce[0], pce[1], pce[2]], estimator, graph\n",
        "        )\n",
        "        last_fun[\"value\"] = val\n",
        "        return val\n",
        "\n",
        "    np.random.seed(seed)\n",
        "    initial_params = np.random.rand(qc.num_parameters)\n",
        "    result = minimize(\n",
        "        loss_func, initial_params, method=\"COBYLA\", options={\"rhobeg\": 1.0}\n",
        "    )\n",
        "    if counter[\"i\"] >= max_iter:\n",
        "        result = OptimizeResult(\n",
        "            message=\"Return from COBYLA because the objective function \"\n",
        "            \"has been evaluated {max_iter} times.\",\n",
        "            success=False,\n",
        "            status=3,\n",
        "            fun=last_fun[\"value\"],\n",
        "            x=last_x[\"value\"],\n",
        "            nfev=counter[\"i\"],\n",
        "        )\n",
        "print(result)\n",
        "\n",
        "# -------------------------Step 4-------------------------\n",
        "\n",
        "par0, par1, best_index = get_partitions(experiment_result)\n",
        "cut_size = calc_cut_size(graph, par0, par1)\n",
        "print(f\"Cut size: {cut_size}\")\n",
        "\n",
        "best_bits = []\n",
        "cur_bits = []\n",
        "for i in experiment_result[best_index][\"exp_map\"]:\n",
        "    if experiment_result[best_index][\"exp_map\"][i] >= 0:\n",
        "        cur_bits.append(1)\n",
        "    else:\n",
        "        cur_bits.append(0)\n",
        "best_cut, best_bits = swap_partitions(graph, cur_bits)\n",
        "# Print final solution\n",
        "\n",
        "print(\n",
        "    f\"The best max-cut value achieved for a graph with {num_nodes} nodes \"\n",
        "    f\"on {num_qubits} qubits is {best_cut}\"\n",
        ")\n",
        "print(f\"and the specific partition we obtained is {best_bits}\")"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "9d68b120-bf43-4070-b897-013652c824d7",
      "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",
        "  * [Advanced techniques for QAOA](/docs/tutorials/advanced-techniques-for-qaoa)\n",
        "  * [Combine error mitigation options with the Estimator primitive](/docs/tutorials/combine-error-mitigation-techniques)\n",
        "</Admonition>\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "99f8259d-fb92-4fae-ab68-96b1e50929a0",
      "metadata": {},
      "source": [
        "## References\n",
        "\n",
        "\\[1] Sciorilli, M., Borges, L., Patti, T. L., García-Martín, D., Camilo, G., Anandkumar, A., & Aolita, L. (2024). Towards large-scale quantum optimization solvers with few qubits. arXiv preprint [arXiv:2401.09421](https://arxiv.org/abs/2401.09421).\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "a4d97730-3c7f-4ce9-8bb9-e4e1c3801c10",
      "metadata": {},
      "source": [
        "## Tutorial survey\n",
        "\n",
        "Please take this short survey to provide feedback on this tutorial. Your insights will help us improve our content offerings and user experience.\n",
        "\n",
        "[Link to survey](https://your.feedback.ibm.com/jfe/form/SV_8ANZAlsKSFf6DA2)\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "id": "0c934e9b-2864-4292-94c9-7fa9b5bce007",
      "metadata": {},
      "source": [
        "© IBM Corp. 2024-2026\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"
    },
    "hours": 1.5,
    "qpuSeconds": 2100
  },
  "nbformat": 4,
  "nbformat_minor": 4
}