Guide-to-cuda-q-backends.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# SPDX-License-Identifier: Apache-2.0 AND CC-BY-NC-4.0\n",
    "#\n",
    "# Licensed under the Apache License, Version 2.0 (the \"License\");\n",
    "# you may not use this file except in compliance with the License.\n",
    "# You may obtain a copy of the License at\n",
    "#\n",
    "# http://www.apache.org/licenses/LICENSE-2.0\n",
    "#\n",
    "# Unless required by applicable law or agreed to in writing, software\n",
    "# distributed under the License is distributed on an \"AS IS\" BASIS,\n",
    "# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.\n",
    "# See the License for the specific language governing permissions and\n",
    "# limitations under the License."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Accelerating Quantum Computing: A Step-by-Step Guide to Expanding Simulation Capabilities and Enabling Interoperability of Quantum Hardware\n",
    "  \n",
    "## Overview of methods of accelerating quantum simulation with GPUs\n",
    "\n",
    "This notebook includes the following: \n",
    "* Introduction to CUDA-Q through two Hello World examples using `sample` and `observe` calls.  \n",
    "\n",
    "* Guide to different backends for executing quantum circuits, emphasizing a variety of patterns of parallelization: \n",
    "    * Statevector memory over multiple processors for simulation\n",
    "    * Circuit sampling over multiple processors\n",
    "    * Hamiltonian batching\n",
    "    * Circuit cutting\n",
    "    * Quantum hardware\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Hello World Examples "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Uncomment and execute this cell if necessary\n",
    "#!pip install cuda-quantum==0.8.0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "import cudaq\n",
    "from cudaq import spin\n",
    "from typing import List\n",
    "import numpy as np\n",
    "import sys"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the first example below, we demonstrata how to define and sample a quantum kernel that encodes a quantum circuit. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Example 1 - defining, drawing, and sampling a quantum kernel\n",
    "\n",
    "##############################################################\n",
    "#  1. Select a backend for kernel execution\n",
    "cudaq.set_target(\"qpp-cpu\")\n",
    "##############################################################\n",
    "\n",
    "##############################################################\n",
    "# 2. Define a kernel function \n",
    "@cudaq.kernel\n",
    "def kernel(qubit_count: int):\n",
    "    # Allocate our `qubit_count` to the kernel.\n",
    "    qvector = cudaq.qvector(qubit_count)\n",
    "\n",
    "    # Apply gates to the qubit indexed by 0.\n",
    "    # CUDA-Q has several built in gates beyond the few examples below\n",
    "    # For a full list see https://nvidia.github.io/cuda-quantum/latest/api/default_ops.html\n",
    "    z(qvector[0])\n",
    "    z(qvector[0])\n",
    "    s(qvector[0])\n",
    "    t(qvector[0])\n",
    "    s(qvector[0])\n",
    "    h(qvector[0])\n",
    "        \n",
    "    # Apply gates to all qubits\n",
    "    x(qvector)\n",
    "    \n",
    "    # Apply a Controlled-X gate between qubit 0 (acting as the control)\n",
    "    # and each of the remaining qubits.  \n",
    "    for i in range(1, qubit_count):\n",
    "        x.ctrl(qvector[0], qvector[i])\n",
    "\n",
    "    # Measure the qubits\n",
    "    # If we don't specify measurements, all qubits are measured in\n",
    "    # the Z-basis by default.\n",
    "    mz(qvector)\n",
    "\n",
    "##############################################################\n",
    "# 3. Call the kernel function with the variable qubit_count set to 2 and sample the outcomes\n",
    "qubit_count = 2\n",
    "result = cudaq.sample(kernel, qubit_count, shots_count=1000)\n",
    "print(result)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now it's your turn to try it out. \n",
    "\n",
    "**Exercise 1:** Edit the code below to create a kernel that produces a circuit for the GHZ state with 4 qubits"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Exercise 1 - Edit the code below to create a kernel that produces a\n",
    "# circuit for the GHZ state with 4 qubits\n",
    "\n",
    "##############################################################\n",
    "#  1. Select a backend for kernel execution\n",
    "cudaq.set_target(\"qpp-cpu\")\n",
    "##############################################################\n",
    "\n",
    "##############################################################\n",
    "# 2. Define a kernel function \n",
    "@cudaq.kernel\n",
    "def kernel(qubit_count: int):\n",
    "    # Allocate our `qubit_count` to the kernel.\n",
    "    qvector = cudaq.qvector(qubit_count)\n",
    "\n",
    "    ##############################################################\n",
    "    # Edit code below\n",
    "    \n",
    "    # Apply a Hadamard gate to the qubit indexed by 0.\n",
    "    \n",
    "    \n",
    "    # Apply a Controlled-X gate between qubit 0 (acting as the control)\n",
    "    # and each of the remaining qubits.  \n",
    "\n",
    "    # Measure the qubits\n",
    "\n",
    "    # Edit code above\n",
    "    ##############################################################\n",
    "    \n",
    "##############################################################\n",
    "# 3. Call the kernel function with the variable qubit_count set to 2 and sample the outcomes\n",
    "qubit_count = 4\n",
    "result = cudaq.sample(kernel, qubit_count, shots_count=1000)\n",
    "\n",
    "print(result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The next example illustrates a few things:\n",
    "* Kernels can be used to define subcircuits.  \n",
    "* `cudaq.draw` can produce ascii or LaTeX circuit diagrams\n",
    "* We can define Hamiltonians with `spin` operators and compute expecation values with `observe`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Example 2 - Expectation value calculations\n",
    "\n",
    "# Define a quantum kernel function to apply a CNOT gate between a control qubit and a each qubit in a list of target qubits\n",
    "\n",
    "@cudaq.kernel\n",
    "def cnot_kernel(control: cudaq.qubit, targets: cudaq.qview):\n",
    "    # Apply a Controlled-X gate between qubit 0 (acting as the control)\n",
    "    # and each of the remaining qubits.  \n",
    "    for i in range(len(targets)):\n",
    "        x.ctrl(control, targets[i])\n",
    "\n",
    "# Define a quantum kernel function to generate a GHZ state on multpile qubits\n",
    "@cudaq.kernel\n",
    "def kernel(qubit_count: int):\n",
    "    # Allocate our `qubit_count` to the kernel.\n",
    "    control_qubit = cudaq.qubit()\n",
    "    target_qubits = cudaq.qvector(qubit_count-1)\n",
    "\n",
    "    # Apply a Hadamard gate to the qubit indexed by 0.\n",
    "    h(control_qubit)\n",
    "    # Apply a Controlled-X gate between qubit 0 (acting as the control)\n",
    "    # and each of the remaining qubits.  \n",
    "    cnot_kernel(control_qubit, target_qubits)\n",
    "\n",
    "# Define a Hamiltonian in terms of Pauli Spin operators.\n",
    "hamiltonian = spin.z(0) + 2*spin.y(1) - spin.x(0) * spin.z(1) - spin.i(2)\n",
    "\n",
    "# Compute the expectation value given the state prepared by the kernel.\n",
    "qubit_count = 3\n",
    "result = cudaq.observe(kernel, hamiltonian, qubit_count, shots_count = 1000).expectation()\n",
    "\n",
    "print(cudaq.draw(kernel, qubit_count)) \n",
    "#print(cudaq.draw('latex', kernel, qubit_count)) # Use the 'latex' option to generate LaTeX code for\n",
    "print(hamiltonian)\n",
    "print('<psi|H|psi> =', result)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Guide to Different Simulation Targets\n",
    "\n",
    "\n",
    "The figure below illustrates a few options for accelerating statevector simulations of single quantum processor kernel executions on one CPU, one GPU, or a multi-node, multi-GPU system. \n",
    "\n",
    "![](https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/single-processor-backends.jpg)\n",
    "\n",
    "In the Hello World examples in the previous section, we saw statevector simulations of a QPU on a CPU.  When GPU resources are available, we can use a single-GPU or multi-node, multi-GPU systems for fast statevector simulations. The `nvidia` target accelerates statevector simulations through `cuStateVec` library. This target offers a variety of configuration options:\n",
    "\n",
    "* **Single-precision GPU simulation** (default): The default of setting the target to `nvidia` through the command `cudaq.set_target('nvidia')` provides single (`fp32`) precision statevector simulation on one GPU.\n",
    "\n",
    "* **Double fp64 precision on a single-GPU**: The option `cudaq.set_target('nvidia', option='fp64')` increases the precision of the statevector simulation on one GPU.\n",
    "\n",
    "* **Multi-node, multi-GPU simulation**: To run the `cuStateVec` simulator on multiple GPUs, set the target to `nvidia` with the `mgpu` option (`cudaq.set_target('nvidia', option='mgpu,fp64')`) and then run the python file containing your quantum kernels within a `MPI` context: `mpiexec -np 2 python3 program.py`. Adjust the `-np` tag according to the number of GPUs you have available.\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, we'll cover a few of the ways you can organize the distribution of quantum simulations over multiple GPU processors, whether you are simulating a single quantum processing unit (QPU) or multiple QPUs."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Single-QPU Statevector Simulations\n",
    "\n",
    " \n",
    "\n",
    "\n",
    "In some cases, the memory required to hold the entire statevector for a simulation exceeds the memory of a single GPU. In these cases, we can distribute the statevector across multiple GPUs as the diagram in the image below suggests.  \n",
    "\n",
    "\n",
    "![](https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/statevector-distribution.png)\n",
    "\n",
    "This is handled automatically within the `mgpu` option when the number of qubits in the statevector exceeds 25.  By changing the environmental variable `CUDAQ_MGPU_NQUBITS_THRESH` prior to setting the target, you can change the threshold at which the statevector distribution is invoked.\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Simulating Parallel QPU computaiton\n",
    "\n",
    "Future quantum computers will accelerate performance by linking multiple QPUs for parallel processing. Today, you can simulate and test programs for these systems using GPUs, and with minimal changes to the target platform, the same code can be executed on multi-QPU setups once they are developed.\n",
    "\n",
    "![](https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/multi-qpus.png)\n",
    "\n",
    "We'll examine a few multi-QPU parallelization patterns here:\n",
    "\n",
    "* Circuit sampling distributed over multiple processors\n",
    "* Hamiltonian batching\n",
    "* Circuit cutting\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### Circuit Sampling\n",
    "\n",
    "One method of parallelization is to sample a circuit over several processors as illustrated in the diagram below.\n",
    "\n",
    "![](https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/circuit-sampling.png)\n",
    "\n",
    "Check out the [documentation](https://nvidia.github.io/cuda-quantum/latest/using/backends/platform.html#remote-mqpu-platform) for code that demonstrates how to launch asynchronous sampling tasks using `sample_async` on multiple virtual QPUs, each simulated by a tensornet simulator backend using the `remote-mqpu` target."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### Hamiltonian Batching\n",
    "Another method for distributing the computational load in a simulation is Hamiltonian batching. In this approach, the expectation values of the Hamiltonian's terms are calculated in parallel across several virtual QPUs, as illustrated in the image below.\n",
    "\n",
    "![](https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/Hamiltonian-batching.png)\n",
    "\n",
    "The `nvidia-mqpu`option of the `nvidia` target along with the `execution=cudaq.parallel.thread` option in the `observe` call handles the distribution of the expectation value computations of a multi-term Hamiltonian across multiple virtual QPUs for you.  Refer to the example below to see how this is carried out:  \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cudaq.set_target(\"nvidia\", option=\"mqpu\")\n",
    "target = cudaq.get_target()\n",
    "num_qpus = target.num_qpus()\n",
    "print(\"Number of QPUs:\", num_qpus)\n",
    "\n",
    "\n",
    "# Define spin ansatz.\n",
    "@cudaq.kernel\n",
    "def kernel(angle: float):\n",
    "    qvector = cudaq.qvector(2)\n",
    "    x(qvector[0])\n",
    "    ry(angle, qvector[1])\n",
    "    x.ctrl(qvector[1], qvector[0])\n",
    "\n",
    "\n",
    "# Define spin Hamiltonian.\n",
    "hamiltonian = 5.907 - 2.1433 * spin.x(0) * spin.x(1) - 2.1433 * spin.y(\n",
    "    0) * spin.y(1) + .21829 * spin.z(0) - 6.125 * spin.z(1)\n",
    "\n",
    "exp_val = cudaq.observe(kernel,\n",
    "                        hamiltonian,\n",
    "                        0.59,\n",
    "                        execution=cudaq.parallel.thread).expectation()\n",
    "print(\"Expectation value: \", exp_val)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the above code snippet, since the Hamiltonian contains four non-identity terms, there are four quantum circuits that need to be executed. When the `nvidia-mqpu` platform is selected, these circuits will be distributed across all available QPUs. The final expectation value result is computed from all QPU execution results.  \n",
    "\n",
    "An alternative method for orchestrating Hamiltonian batching is to use the MPI context and multiple GPUs.  You can read more about this [here](https://nvidia.github.io/cuda-quantum/latest/using/backends/platform.html#nvidia-mqpu-platform)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### Circuit cutting\n",
    "\n",
    "Circuit cutting is a widely used technique for parallelization. One way to conceptualize circuit cutting is through the Max Cut problem. In this scenario, we aim to approximate the Max Cut of a graph using a divide-and-conquer strategy, also known as QAOA-in-QAOA or QAOA². This approach breaks the graph into smaller subgraphs and solves the Max Cut for each subgraph in parallel using QAOA (see references such as [arXiv:2205.11762v1](https://arxiv.org/abs/2205.11762), [arxiv.2101.07813v1](https://arxiv.org/abs/2101.07813), [arxiv:2304.03037v1](https://arxiv.org/abs/2304.03037), [arxiv:2009.06726](https://arxiv.org/abs/2009.06726), and [arxiv:2406:17383](https://arxiv.org/abs/2406.17383)). By doing so, we effectively decompose the QAOA circuit for the larger graph into smaller QAOA circuits for the subgraphs.\n",
    "\n",
    "![](https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/qaoa-cut.png)\n",
    "\n",
    "To complete the circuit cutting, we'll need to merge the results of QAOA on the subgraphs into a result for the entire graph.  This requires solving another smaller optimization problem, which can also be tackled with QAOA.  You can read about that in more detail in a series of [interactive labs](https://github.com/NVIDIA/cuda-q-academic/tree/main/qaoa-for-max-cut).\n",
    "\n",
    "![](https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/circuit-cutting.png)\n",
    "\n",
    "This example illustrates how to use the `MPI` context to orchestrate running `@cudaq.kernel` decorated functions in parallel. Additionally, a few exercises are built into this longer example to provide some practice with the CUDA-Q commands introduced earlier in this notebook. Solutions to these exercises appear in the [solutions-sc24.ipynb](solutions-sc24.ipynb) file, but we encourage you to first attempt the exercises out yourself."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First we need to define a graph and subgraphs.  Let's start with the graph drawn below.\n",
    "\n",
    "<img src=\"https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/example-graph.png\" width=\"500\">\n",
    "\n",
    "For this demonstration, we'll divide our example graph into the five subgraphs depicted below:\n",
    "\n",
    "<img src=\"https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/subgraphs.png\" width=\"500\">\n",
    "\n",
    "Execute the cell below to generate subgraphs for the divide-and-conquer QAOA. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Identify subgraphs, separating out the edges as source and target nodes\n",
    "num_subgraphs = 5 # Number of subgraphs\n",
    "nodeCountList = [8,7,6,5,4] # Number of nodes in each subgraph\n",
    "nodeList : List[int] = [] # List of nodes in each of the subgraphs \n",
    "edgeListSources : List[int] = [] # List of edge sources in each subgraph\n",
    "edgeListTargets : List[int] = [] # List of edge targets in each subgraph\n",
    "\n",
    "# subgraph0 data\n",
    "nodeList.append([3, 6, 9, 10, 13, 14, 21, 22])\n",
    "edgeListSources.append([3,3,3,3,6,6,9,14])\n",
    "edgeListTargets.append([14,9,10,13,22,13,21,22])\n",
    "\n",
    "# subgraph1 data\n",
    "nodeList.append([8, 11, 12, 15, 16, 25, 26])\n",
    "edgeListSources.append([8, 8, 11, 11, 11, 11, 12, 15, 16, 16, 25])\n",
    "edgeListTargets.append([25, 12, 26, 25, 15, 12, 15, 16, 25, 26, 26])\n",
    "\n",
    "# subgraph2 data\n",
    "nodeList.append([4, 5, 7, 18, 20, 24])\n",
    "edgeListSources.append([4, 4, 5, 7, 18, 20])\n",
    "edgeListTargets.append([5, 24, 7, 24, 20, 24])\n",
    "\n",
    "# subgraph3 data\n",
    "nodeList.append([0, 19, 27, 28, 29])\n",
    "edgeListSources.append([0, 0, 19, 19, 27, 27])\n",
    "edgeListTargets.append([19, 28, 27, 29, 29, 28])\n",
    "\n",
    "# subgraph4 data\n",
    "nodeList.append([1, 2, 17, 23])\n",
    "edgeListSources.append([1, 1, 2, 17])\n",
    "edgeListTargets.append([23, 2, 17, 23])\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, we need a helper function that will be used to map graph nodes to qubits."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# We'll create this function to rename the nodes to be sequential integers\n",
    "# beginning with 0 as a way to map the graph nodes to qubits. We take this\n",
    "# approach because we can't use the `.index` option for lists\n",
    "# within a cudaq.kernel.\n",
    "def rename_nodes(edge_src, edge_tgt, nodes):\n",
    "    \"\"\" \n",
    "    Parameters\n",
    "    ----------\n",
    "    edges_src: List[int]\n",
    "        List of the first (source) node listed in each edge of the graph, when the edges of the graph are listed as pairs of nodes\n",
    "    edges_tgt: List[int]\n",
    "        List of the second (target) node listed in each edge of the graph, when the edges of the graph are listed as pairs of nodes\n",
    "    nodes: List[int]\n",
    "        List of nodes of the graph\n",
    "    \n",
    "    Returns\n",
    "    -------\n",
    "    new_edge_src : List[int]\n",
    "        List of the first (source) node listed in each edge of the graph after renaming nodes to be sequential integers beginning with 0, \n",
    "        when the edges of the graph are listed as pairs of nodes\n",
    "    new_edge_tgt : List[int]\n",
    "        List of the second (target) node listed in each edge of the graph after renaming nodes to be sequential integers beginning with 0, \n",
    "        when the edges of the graph are listed as pairs of nodes\n",
    "    \"\"\"\n",
    "    new_edge_src = []\n",
    "    new_edge_tgt = []\n",
    "    for i in range(len(edge_src)):\n",
    "        new_edge_src.append(nodes.index(edge_src[i]))\n",
    "        new_edge_tgt.append(nodes.index(edge_tgt[i]))\n",
    "    return new_edge_src, new_edge_tgt    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next let's create kernels to combine into the QAOA circuit:\n",
    "\n",
    "* `qaoaProblem` kernel adds the gate sequence depicted below for each edge in the graph\n",
    "\n",
    "<img src=\"https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/problem-kernel.png\" width=\"500\">\n",
    "\n",
    "* `qaoaMixer` applies a parameterized `rx` gate to all the qubits, highlighted in green in the diagram below\n",
    "\n",
    "* `kernel_qaoa` builds the QAOA circuit drawn below using the `qaoaProblem` and `qaoaMixer`\n",
    "\n",
    "<img src=\"https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/qaoa-circuit.png\" width=\"500\">\n",
    "\n",
    "**Exercise 2:**  The `kernel_qaoa` kernel has been defined for you.  Your task is edit the two `###FIX_ME###`s in the code below to complete the `qaoaProblem` and `qaoaMixer` kernels.  \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Exercise 2 Edit the two ###FIX_ME###s in the code below. One is in `qaoaProblem` and the other is in `qaoaMixer`\n",
    "\n",
    "# Problem Kernel\n",
    "@cudaq.kernel\n",
    "def qaoaProblem(qubit_0 : cudaq.qubit, qubit_1 : cudaq.qubit, alpha : float):\n",
    "    \"\"\"Build the QAOA gate sequence between two qubits that represent an edge of the graph\n",
    "    Parameters\n",
    "    ----------\n",
    "    qubit_0: cudaq.qubit\n",
    "        Qubit representing the first vertex of an edge\n",
    "    qubit_1: cudaq.qubit\n",
    "        Qubit representing the second vertex of an edge\n",
    "    alpha: float\n",
    "        Free variable\n",
    "\n",
    "    \"\"\"\n",
    "    x.ctrl(qubit_0, qubit_1)\n",
    "    ###FIX_ME###\n",
    "    x.ctrl(qubit_0, qubit_1)\n",
    "\n",
    "# Mixer Kernel\n",
    "@cudaq.kernel\n",
    "def qaoaMixer(###FIX_ME###):\n",
    "    \"\"\"Build the QAOA gate sequence that is applied to each qubit in the mixer portion of the circuit\n",
    "    Parameters\n",
    "    ----------\n",
    "    qubit_0: cudaq.qubit\n",
    "        Qubit\n",
    "    beta: float\n",
    "        Free variable\n",
    "\n",
    "    \"\"\"\n",
    "    rx(2.0*beta, qubits)\n",
    "\n",
    "\n",
    "# We now define the kernel_qaoa function which will build the QAOA circuit for our graph\n",
    "@cudaq.kernel\n",
    "def kernel_qaoa(qubit_count :int, layer_count: int, qubits_src: List[int], qubits_tgt: List[int], thetas : List[float]):\n",
    "    \"\"\"Build the QAOA circuit for max cut of the graph with given edges and nodes\n",
    "    Parameters\n",
    "    ----------\n",
    "    qubit_count: int\n",
    "        Number of qubits in the circuit, which is the same as the number of nodes in our graph\n",
    "    layer_count : int\n",
    "        Number of layers in the QAOA kernel\n",
    "    edges_src: List[int]\n",
    "        List of the first (source) node listed in each edge of the graph, when the edges of the graph are listed as pairs of nodes\n",
    "    edges_tgt: List[int]\n",
    "        List of the second (target) node listed in each edge of the graph, when the edges of the graph are listed as pairs of nodes\n",
    "    thetas: List[float]\n",
    "        Free variables to be optimized\n",
    "    \"\"\"\n",
    "    # Let's allocate the qubits\n",
    "    qreg = cudaq.qvector(qubit_count)\n",
    "\n",
    "    # And then place the qubits in superposition\n",
    "    h(qreg)\n",
    "    \n",
    "    # Each layer has two components: the problem kernel and the mixer\n",
    "    for i in range(layer_count):\n",
    "        # Add the problem kernel to each layer\n",
    "        for edge in range(len(qubits_src)):\n",
    "            qubitu = qubits_src[edge]\n",
    "            qubitv = qubits_tgt[edge]\n",
    "            qaoaProblem(qreg[qubitu], qreg[qubitv], thetas[i])\n",
    "        # Add the mixer kernel to each layer\n",
    "        qaoaMixer(qreg,thetas[i+layer_count])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We'll need a Hamiltonian to encode the cost function:   $$H= \\frac{1}{2}\\sum_{(u,v)\\in E} (Z_uZ_v-II),$$ where $E$ is the set of edges of the graph."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define a function to generate the Hamiltonian for a max cut problem using the graph G\n",
    "\n",
    "def hamiltonian_max_cut(sources : List[int], targets : List[int]):\n",
    "    \"\"\"Hamiltonian for finding the max cut for the graph with edges defined by the pairs generated by sources and targets\n",
    "\n",
    "    Parameters\n",
    "    ----------\n",
    "    sources: List[int]\n",
    "        list of the source vertices for edges in the graph\n",
    "    targets: List[int]\n",
    "        list of the target vertices for the edges in the graph\n",
    "\n",
    "    Returns\n",
    "    -------\n",
    "    cudaq.SpinOperator\n",
    "        Hamiltonian for finding the max cut of the graph defined by the given edges\n",
    "    \"\"\"\n",
    "    hamiltonian = 0\n",
    "   \n",
    "    # Since our vertices may not be a list from 0 to n, or may not even be integers,\n",
    "    # we need to map the vertices to the list of integers 0 to qubit_count -1\n",
    "    \n",
    "    for i in range(len(sources)):\n",
    "        # Add a term to the Hamiltonian for the edge (u,v)\n",
    "        qubitu = sources[i]\n",
    "        qubitv = targets[i]\n",
    "        hamiltonian += 0.5*(spin.z(qubitu)*spin.z(qubitv)-spin.i(qubitu)*spin.i(qubitv))\n",
    "\n",
    "    return hamiltonian"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's put this all together in a function that finds the the optimal parameters for QAOA of a given subgraph."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "metadata": {},
   "outputs": [],
   "source": [
    "def find_optimal_parameters(qubit_src : List[int], qubit_tgt : List[int], qubit_count : int, layer_count: int, seed :int):\n",
    "    \"\"\"Function for finding the optimal parameters of QAOA for the max cut of a graph\n",
    "    Parameters\n",
    "    ----------\n",
    "    qubit_src: List[int]\n",
    "    qubit_tgt: List[int]\n",
    "        Sources and targets defining the edges of the graph\n",
    "    nodes: List[int]\n",
    "        Integer labels of the nodes of the graph\n",
    "    qubit_count: int\n",
    "        qubit_count is the number of nodes in the graph\n",
    "    layer_count : int\n",
    "        Number of layers in the QAOA circuit\n",
    "    seed : int\n",
    "        Random seed for reproducibility of results\n",
    "\n",
    "    Returns\n",
    "    -------\n",
    "    list[float]\n",
    "        Optimal parameters for the QAOA applied to the given graph\n",
    "    \"\"\"\n",
    "    \n",
    "    # Each layer of the QAOA kernel contains 2 parameters\n",
    "    parameter_count : int = 2*layer_count\n",
    "\n",
    "    # Specify the optimizer and its initial parameters.\n",
    "    optimizer = cudaq.optimizers.COBYLA()\n",
    "    np.random.seed(seed)\n",
    "    optimizer.initial_parameters = np.random.uniform(-np.pi, np.pi,\n",
    "                                                     parameter_count)\n",
    "    \n",
    "    # Pass the kernel, spin operator, and optimizer to `cudaq.vqe`.\n",
    "    optimal_expectation, optimal_parameters = cudaq.vqe(\n",
    "        kernel=kernel_qaoa,\n",
    "        spin_operator=hamiltonian_max_cut(qubit_src, qubit_tgt),\n",
    "        argument_mapper=lambda parameter_vector: (qubit_count, layer_count, qubit_src, qubit_tgt, parameter_vector),\n",
    "        optimizer=optimizer,\n",
    "        parameter_count=parameter_count)\n",
    "\n",
    "    return optimal_parameters"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Before running this function in parallel, let's execute it sequentially. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[2.760265752934724, -1.1781480121902836], [2.863390666464872, -1.2262588372574177], [2.7568258048030376, -1.1779328954753276], [2.8169010481390853, -1.2171055596613192], [2.7488320354701976, -1.178054032241366]]\n"
     ]
    }
   ],
   "source": [
    "# Testing the find_optimal_parameters function\n",
    "layer_count = 1\n",
    "seed = 123\n",
    "\n",
    "new_src = []\n",
    "new_tgt = []\n",
    "optimal_parameters = []\n",
    "\n",
    "for i in range(num_subgraphs):\n",
    "    src, tgt= rename_nodes(edgeListSources[i], edgeListTargets[i], nodeList[i])\n",
    "    new_src.append(src)\n",
    "    new_tgt.append(tgt)\n",
    "    optimal_parameters.append(find_optimal_parameters(src, tgt, nodeCountList[i], layer_count, seed))\n",
    "\n",
    "print(optimal_parameters)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, we'll need to sample the `kernel_qaoa` circuit with the optimal parameters to find approximate max cut solutions to each of the subgraphs.\n",
    "\n",
    "**Exercise 3:** Edit the `FIX_ME` in the code block below to sample the QAOA circuits for each of the subgraphs using the optimal parameter found above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Exercise 3\n",
    "# Sampling the QAOA circuits with the optimal parameters to identify an appoximate max cut of the subgraphs\n",
    "\n",
    "shots = 10000\n",
    "\n",
    "for i in range(num_subgraphs):\n",
    "    counts = FIX_ME(kernel_qaoa, nodeCountList[i], layer_count, new_src[i], new_tgt[i], optimal_parameters[i], shots_count=shots)\n",
    "    print('subgraph ',i,' has most_probable outcome = ',counts.most_probable())\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "That completes the \"conquer\" stage of the divide-and-conquer algorithm. To learn more about how the results of the subgraph solutions are merged together to get a max cut approximation of the original graph, check out the 2nd notebook of this [series of interactive tutorials](https://github.com/NVIDIA/cuda-q-academic/tree/main/qaoa-for-max-cut). For the purposes of today's tutorial we'll set that step aside and examine how we could parallelize the divide-and-conquer QAOA algorithm using CUDA-Q and `MPI`. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The diagram below illustrates a strategy for implementing the divide-and-conquer QAOA across 4 processors (which could be distinct GPUs or separate processes on a single GPU). The approach involves storing subgraph data in a dictionary, where the keys represent subgraph names. These dictionary keys are distributed among the 4 processors, with each processor responsible for solving the QAOA problem for the subgraphs corresponding to its assigned keys.\n",
    "\n",
    "![](https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/parallel-workflow.png)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We've coded this strategy up in the one-step-qaoa.py file for you.  The command line below executes the [qaoa-divide-and-conquer.py](qaoa-divide-and-conquer.py) file on 4 processors in parallel.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Uncomment and execute this cell to install mpi4py if necessary\n",
    "#%pip install mpi4py"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "# MPI call\n",
    "# Uncomment if you have OpenMPI installed with a GPU\n",
    "#print(sys.executable)\n",
    "#python_path = sys.executable\n",
    "#!mpiexec -np 4 --oversubscribe --allow-run-as-root {python_path} divide-and-conquer-qaoa.py\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The animation below captures a small instance of a recursive divide-and-conquer QAOA running on a CPU versus a GPU in parallel. The lineplots on the top depict the error between the calculated max cut solution and the true max cut of the graph over time. The graphs on the bottom represent the max cut solutions as various subgraph problems are solved and merged together. The green graphs on the right show the parallelization of solving subgraph problems simultaneously.\n",
    "\n",
    "![](https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/maxcut_ani.gif)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Beyond Statevector Simulations\n",
    "\n",
    "#### Other simulators\n",
    "\n",
    "When using CUDA-Q on NVIDIA GPU with available CUDA runtime libraries, the default target is set to `nvidia`. This will utilize the cuQuantum single-GPU statevector simulator. On CPU-only systems, the default target is set to `qpp-cpu` which uses the OpenMP CPU-only statevector simulator.\n",
    "\n",
    "For many applications, it's not necessary to simluate and access the entire statevector. The default simulator can be overridden by the environment variable CUDAQ_DEFAULT_SIMULATOR where tensor network, matrix product state simulators can be selected. Please refer to the table below for a list of backend simulator names along with its multi-GPU capability.\n",
    "\n",
    "![](https://raw.githubusercontent.com/NVIDIA/cuda-q-academic/refs/heads/main/images/accelerating/backends.png)\n",
    "\n",
    "For more information about all the simulator backends available on [this documentation page](https://nvidia.github.io/cuda-quantum/latest/using/backends/simulators.html#tensor-network-simulators).\n",
    "\n",
    "#### Quantum processing units\n",
    "In addition to executing simulations, CUDA-Q is equipped to run quantum kernels on quantum processing units.  For more information on how to execute CUDA-Q code on quantum processing units, check out the [documentation](https://nvidia.github.io/cuda-quantum/latest/using/backends/hardware.html)."
   ]
  }
 ],
 "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.10.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}