04_QAOA.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "ef6fb5e4",
   "metadata": {},
   "source": [
    "# Example: Quantum Approximate Optimization Algorithm (QAOA)\n",
    "\n",
    "### 1. The QAOA algorithm is a quantum algorithm used to solve  Quadratic Unconstrained Binary Optimization (QUBO) problems. \n",
    "\n",
    "### 2. The objective of QUBO problems is to optimize the  quadratic function $C$ subject to a set of constraints. The QUBO problems are typically used to model combinatorial optimization problems, such as the traveling salesman problem, the maximum cut problem (Max-Cut), and the graph coloring problem. \n",
    "\n",
    "\n",
    "### 3. The QAOA was first introduced by E. Farhi, J. Goldstone, and S. Gutmann in [A Quantum Approximate Optimization Algorithm](https://arxiv.org/abs/1411.4028). The algorithm has been shown to be efficient for solving certain optimization problems, such as the Max-Cut problem, which is NP-hard.\n",
    "\n",
    "### In the following, we will implement QAOA for the Max-Cut problem from scratch using Pytorch and its automatic differentiation libraries."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4c707529",
   "metadata": {},
   "source": [
    "# Mathematical formulation of QUBO problems\n",
    "\n",
    "### Mathematically QUBO problem can be defined as follow: Let us denote a set of binary vectors \n",
    "### \\begin{equation}\n",
    "  \\mathbb{S} = \\{ \\vec{s} = (s_1,\\dots,s_L) | L>0,  s_i \\in (0,1) \\}.\n",
    "\\end{equation}\n",
    "\n",
    "### We are given a real-valued upper-triangular matrix $W \\in \\mathbb{R}^{L\\times L}$ whose elements $W_{ij}$ define a weight for  each pair of indices $i,j \\in \\{1,\\dots,L\\}$ within  $\\vec{s}$. We define a QUBO $\\textit{cost}$ function $C(\\vec{s})$ as:\n",
    "### \\begin{equation}\n",
    "C(\\vec{s}) = \\sum_{i=1}^L\\sum_{j=i}^{L} W_{ij}s_is_j.\n",
    "\\end{equation}\n",
    "### The solution of a QUBO problem is a vector $\\vec{s}_*$ that minimizes the cost function, i.e.:\n",
    "### \\begin{equation}\n",
    "\\vec{s}_* = \\arg \\min_{\\vec{s}\\in\\mathbb{S}} C(\\vec{s}).\n",
    "\\end{equation}\n",
    "### As we can see, the search space for the classical optimization algorithm, i.e., size of the vector space $\\mathbb{S}$, grows exponentially with   $L$, $\\dim \\mathbb{S} = 2^L$. Next, we can notice that solution to the considered problem has to belong to a Fock space of the chain of $L$ spins-$1/2$. \n",
    "\n",
    "### The main idea behind the QAOA is to cast a QUBO problem onto the spin-$1/2$ system, where the classical cost function $C$ is represented as a  \\textit{cost} Hamiltonian $\\hat{H}_C$ of the spin-$1/2$ system of size $L$. The optimal solution $s_*$ of the original problem can be inferred from a ground state $|\\psi_*\\rangle$ of the $\\hat{H}_C$. In practice, the solution $\\vec{s}_*$ is given by the spin configuration of a  Fock state vector $\\vec{s}_k$ having the dominant contribution to the $\\hat{H}_C$ ground state $|\\psi_*\\rangle = \\sum_{k} \\alpha_k |\\vec{s}_k\\rangle$:\n",
    "### \\begin{equation}\n",
    " \\begin{split}\n",
    "   k_* &  = \\arg \\max_{k} |\\langle\\vec{s}_k|\\psi_*\\rangle|^2,\\\\\n",
    "   \\vec{s}_* &\\equiv \\vec{s}_{k_*}.\n",
    "   \\end{split}\n",
    "\\end{equation}\n",
    "### As such, the objective is to find a many-body wave function with the lowest possible energy for the cost Hamiltonian $\\hat{H}_C$. In the QAOA we consider an ansatz for the many-body wave function in the form a parametrized quantum circuit\n",
    "### \\begin{equation}\n",
    "|\\psi(\\vec{\\beta},\\vec{\\gamma})\\rangle= \\prod_{j=1}^p \\hat{U}_B(\\beta_j)\\hat{U}_C(\\gamma_j)\\prod_{i=1}^L H_i |0\\rangle,\n",
    "\\end{equation}\n",
    "### where $H_i$ is a Hadamard matrix acting on a $i$-th qubit ($i = 1,\\dots, L)$, while two parametrized unitary operators $\\hat{U}_B$ and $\\hat{U}_C$ are defined as\n",
    "### \\begin{equation}\n",
    "\\begin{split}\n",
    "\\hat{U}_C(\\gamma_j) &= e^{-i\\gamma_j\\hat{H}_C} \\\\\n",
    "\\hat{U}_B(\\beta_j) &= e^{-i\\beta_j\\hat{H}_B},\n",
    "\\end{split}\n",
    "\\end{equation}\n",
    "### where $\\hat{H}_C$ is a considered cost Hamiltonian, and $\\hat{H}_B = \\sum_i \\hat{\\sigma}^x_i$ is a so-called \\textit{mixer} Hamiltonian.\n",
    "### The objective of the algorithm is to find a set of parameters $\\vec{\\beta} = (\\beta_1, \\dots, \\beta_p)$, $\\vec{\\gamma} = (\\gamma_1,\\dots,\\gamma_p)$ ($p\\ge1$ is a hyperparameter of the algorithm) which minimizes the cost function\n",
    "### \\begin{equation}\n",
    "F(\\vec{\\beta},\\vec{\\gamma}) = \\langle\\psi(\\vec{\\beta},\\vec{\\gamma})|\\hat{H}_C|\\psi(\\vec{\\beta},\\vec{\\gamma})\\rangle.\n",
    "\\end{equation}\n",
    "### Finding $\\vec{\\beta}_*$ and $\\vec{\\gamma}_*$ can be done with the classical optimization algorithms in a similar way as in VQE problem. In practice, we run QAOA many times  with random initialization of $\\vec{\\beta}$, and $\\vec{\\gamma}$ and we plot the histogram of the index $k_*$ of the optimal configuration given by the corresponding Fock vector $|\\vec{s}_{k_*}\\rangle$. The solution is given by $|\\vec{s}_{k_*}\\rangle$ corresponding to a peak in the histogram of all $k_*$."
   ]
  },
  {
   "attachments": {
    "Max-cut.svg": {
     "image/svg+xml": [
      "PD94bWwgdmVyc2lvbj0iMS4wIiBlbmNvZGluZz0iVVRGLTgiIHN0YW5kYWxvbmU9Im5vIj8+CjwhLS0gQ3JlYXRlZCB3aXRoIElua3NjYXBlIChodHRwOi8vd3d3Lmlua3NjYXBlLm9yZy8pIC0tPgo8c3ZnCiAgIHhtbG5zOmRjPSJodHRwOi8vcHVybC5vcmcvZGMvZWxlbWVudHMvMS4xLyIKICAgeG1sbnM6Y2M9Imh0dHA6Ly9jcmVhdGl2ZWNvbW1vbnMub3JnL25zIyIKICAgeG1sbnM6cmRmPSJodHRwOi8vd3d3LnczLm9yZy8xOTk5LzAyLzIyLXJkZi1zeW50YXgtbnMjIgogICB4bWxuczpzdmc9Imh0dHA6Ly93d3cudzMub3JnLzIwMDAvc3ZnIgogICB4bWxucz0iaHR0cDovL3d3dy53My5vcmcvMjAwMC9zdmciCiAgIHhtbG5zOnNvZGlwb2RpPSJodHRwOi8vc29kaXBvZGkuc291cmNlZm9yZ2UubmV0L0RURC9zb2RpcG9kaS0wLmR0ZCIKICAgeG1sbnM6aW5rc2NhcGU9Imh0dHA6Ly93d3cuaW5rc2NhcGUub3JnL25hbWVzcGFjZXMvaW5rc2NhcGUiCiAgIHdpZHRoPSIyMDAiCiAgIGhlaWdodD0iMTYwIgogICBpZD0ic3ZnMiIKICAgc29kaXBvZGk6dmVyc2lvbj0iMC4zMiIKICAgaW5rc2NhcGU6dmVyc2lvbj0iMC40NiIKICAgc29kaXBvZGk6ZG9jbmFtZT0ibWF4LWN1dC5zdmciCiAgIGlua3NjYXBlOm91dHB1dF9leHRlbnNpb249Im9yZy5pbmtzY2FwZS5vdXRwdXQuc3ZnLmlua3NjYXBlIgogICB2ZXJzaW9uPSIxLjAiPgogIDxkZWZzCiAgICAgaWQ9ImRlZnM0IiAvPgogIDxzb2RpcG9kaTpuYW1lZHZpZXcKICAgICBpZD0iYmFzZSIKICAgICBwYWdlY29sb3I9IiNmZmZmZmYiCiAgICAgYm9yZGVyY29sb3I9IiM2NjY2NjYiCiAgICAgYm9yZGVyb3BhY2l0eT0iMS4wIgogICAgIGdyaWR0b2xlcmFuY2U9IjEwMDAwIgogICAgIGd1aWRldG9sZXJhbmNlPSIxMCIKICAgICBvYmplY3R0b2xlcmFuY2U9IjEwIgogICAgIGlua3NjYXBlOnBhZ2VvcGFjaXR5PSIwLjAiCiAgICAgaW5rc2NhcGU6cGFnZXNoYWRvdz0iMiIKICAgICBpbmtzY2FwZTp6b29tPSIyLjgyODQyNzEiCiAgICAgaW5rc2NhcGU6Y3g9IjIxMC42NzIxNiIKICAgICBpbmtzY2FwZTpjeT0iMjguODMyNDM3IgogICAgIGlua3NjYXBlOmRvY3VtZW50LXVuaXRzPSJweCIKICAgICBpbmtzY2FwZTpjdXJyZW50LWxheWVyPSJsYXllcjIiCiAgICAgc2hvd2dyaWQ9InRydWUiCiAgICAgaW5rc2NhcGU6d2luZG93LXdpZHRoPSIxMjgwIgogICAgIGlua3NjYXBlOndpbmRvdy1oZWlnaHQ9Ijk2NiIKICAgICBpbmtzY2FwZTp3aW5kb3cteD0iMCIKICAgICBpbmtzY2FwZTp3aW5kb3cteT0iMCI+CiAgICA8aW5rc2NhcGU6Z3JpZAogICAgICAgdHlwZT0ieHlncmlkIgogICAgICAgaWQ9ImdyaWQyMzgzIgogICAgICAgdmlzaWJsZT0idHJ1ZSIKICAgICAgIGVuYWJsZWQ9InRydWUiIC8+CiAgPC9zb2RpcG9kaTpuYW1lZHZpZXc+CiAgPG1ldGFkYXRhCiAgICAgaWQ9Im1ldGFkYXRhNyI+CiAgICA8cmRmOlJERj4KICAgICAgPGNjOldvcmsKICAgICAgICAgcmRmOmFib3V0PSIiPgogICAgICAgIDxkYzpmb3JtYXQ+aW1hZ2Uvc3ZnK3htbDwvZGM6Zm9ybWF0PgogICAgICAgIDxkYzp0eXBlCiAgICAgICAgICAgcmRmOnJlc291cmNlPSJodHRwOi8vcHVybC5vcmcvZGMvZGNtaXR5cGUvU3RpbGxJbWFnZSIgLz4KICAgICAgPC9jYzpXb3JrPgogICAgPC9yZGY6UkRGPgogIDwvbWV0YWRhdGE+CiAgPGcKICAgICBpbmtzY2FwZTpncm91cG1vZGU9ImxheWVyIgogICAgIGlkPSJsYXllcjIiCiAgICAgaW5rc2NhcGU6bGFiZWw9IkVkZ2VzIgogICAgIHN0eWxlPSJkaXNwbGF5OmlubGluZSIKICAgICB0cmFuc2Zvcm09InRyYW5zbGF0ZSgtMTI5LC00MDEuNzAwMDIpIj4KICAgIDxwYXRoCiAgICAgICBzdHlsZT0iZmlsbDpub25lO2ZpbGwtcnVsZTpldmVub2RkO3N0cm9rZTojZDQwMDAwO3N0cm9rZS13aWR0aDo0O3N0cm9rZS1saW5lY2FwOmJ1dHQ7c3Ryb2tlLWxpbmVqb2luOm1pdGVyO3N0cm9rZS1taXRlcmxpbWl0OjQ7c3Ryb2tlLWRhc2hhcnJheTpub25lO3N0cm9rZS1vcGFjaXR5OjEiCiAgICAgICBkPSJNIDE1MCw0ODIuMzYyMTggTCAyMTAsNDQyLjM2MjE4IgogICAgICAgaWQ9InBhdGgyMzk4IiAvPgogICAgPHBhdGgKICAgICAgIHN0eWxlPSJmaWxsOm5vbmU7ZmlsbC1ydWxlOmV2ZW5vZGQ7c3Ryb2tlOiMwMDAwMDA7c3Ryb2tlLXdpZHRoOjI7c3Ryb2tlLWxpbmVjYXA6YnV0dDtzdHJva2UtbGluZWpvaW46bWl0ZXI7c3Ryb2tlLW1pdGVybGltaXQ6NDtzdHJva2UtZGFzaGFycmF5Om5vbmU7c3Ryb2tlLW9wYWNpdHk6MSIKICAgICAgIGQ9Ik0gMTUwLDQ4Mi4zNjIxOCBMIDIxMCw1MjIuMzYyMTgiCiAgICAgICBpZD0icGF0aDI0MDAiIC8+CiAgICA8cGF0aAogICAgICAgc3R5bGU9ImZpbGw6bm9uZTtmaWxsLXJ1bGU6ZXZlbm9kZDtzdHJva2U6I2Q0MDAwMDtzdHJva2Utd2lkdGg6MztzdHJva2UtbGluZWNhcDpidXR0O3N0cm9rZS1saW5lam9pbjptaXRlcjtzdHJva2UtbWl0ZXJsaW1pdDo0O3N0cm9rZS1kYXNoYXJyYXk6bm9uZTtzdHJva2Utb3BhY2l0eToxIgogICAgICAgZD0iTSAyMTAsNDQyLjM2MjE4IEwgMjkwLDQ0Mi4zNjIxOCIKICAgICAgIGlkPSJwYXRoMjQwMiIgLz4KICAgIDxwYXRoCiAgICAgICBzdHlsZT0iZmlsbDpub25lO2ZpbGwtcnVsZTpldmVub2RkO3N0cm9rZTojZDQwMDAwO3N0cm9rZS13aWR0aDo0O3N0cm9rZS1saW5lY2FwOmJ1dHQ7c3Ryb2tlLWxpbmVqb2luOm1pdGVyO3N0cm9rZS1taXRlcmxpbWl0OjQ7c3Ryb2tlLWRhc2hhcnJheTpub25lO3N0cm9rZS1vcGFjaXR5OjEiCiAgICAgICBkPSJNIDIxMCw1MjIuMzYyMTggTCAyOTAsNTIyLjM2MjE4IgogICAgICAgaWQ9InBhdGgyNDA0IiAvPgogICAgPHBhdGgKICAgICAgIHN0eWxlPSJmaWxsOm5vbmU7ZmlsbC1ydWxlOmV2ZW5vZGQ7c3Ryb2tlOiNkNDAwMDA7c3Ryb2tlLXdpZHRoOjQ7c3Ryb2tlLWxpbmVjYXA6YnV0dDtzdHJva2UtbGluZWpvaW46bWl0ZXI7c3Ryb2tlLW1pdGVybGltaXQ6NDtzdHJva2UtZGFzaGFycmF5Om5vbmU7c3Ryb2tlLW9wYWNpdHk6MSIKICAgICAgIGQ9Ik0gMjkwLDQ0Mi4zNjIxOCBMIDI5MCw1MjIuMzYyMTgiCiAgICAgICBpZD0icGF0aDI0MDYiIC8+CiAgICA8cGF0aAogICAgICAgaWQ9InBhdGgyNDA5IgogICAgICAgZD0iTSAyMTAsNDQyLjM2MjE4IEwgMjEwLDUyMi4zNjIxOCIKICAgICAgIHN0eWxlPSJmaWxsOm5vbmU7ZmlsbC1ydWxlOmV2ZW5vZGQ7c3Ryb2tlOiNkNDAwMDA7c3Ryb2tlLXdpZHRoOjQ7c3Ryb2tlLWxpbmVjYXA6YnV0dDtzdHJva2UtbGluZWpvaW46bWl0ZXI7c3Ryb2tlLW1pdGVybGltaXQ6NDtzdHJva2UtZGFzaGFycmF5Om5vbmU7c3Ryb2tlLW9wYWNpdHk6MSIgLz4KICAgIDxwYXRoCiAgICAgICBpZD0icGF0aDMxOTIiCiAgICAgICBkPSJNIDIxMCw0NDIuMzYyMTggTCAyOTAsNDQyLjM2MjE4IgogICAgICAgc3R5bGU9ImZpbGw6bm9uZTtmaWxsLXJ1bGU6ZXZlbm9kZDtzdHJva2U6I2Q0MDAwMDtzdHJva2Utd2lkdGg6NDtzdHJva2UtbGluZWNhcDpidXR0O3N0cm9rZS1saW5lam9pbjptaXRlcjtzdHJva2UtbWl0ZXJsaW1pdDo0O3N0cm9rZS1kYXNoYXJyYXk6bm9uZTtzdHJva2Utb3BhY2l0eToxIiAvPgogIDwvZz4KICA8ZwogICAgIGlua3NjYXBlOmxhYmVsPSJOb2RlcyIKICAgICBpbmtzY2FwZTpncm91cG1vZGU9ImxheWVyIgogICAgIGlkPSJsYXllcjEiCiAgICAgc3R5bGU9ImRpc3BsYXk6aW5saW5lIgogICAgIHRyYW5zZm9ybT0idHJhbnNsYXRlKC0xMjksLTQwMS43MDAwMikiPgogICAgPHBhdGgKICAgICAgIHNvZGlwb2RpOnR5cGU9ImFyYyIKICAgICAgIHN0eWxlPSJvcGFjaXR5OjE7ZmlsbDojZmZmZmZmO2ZpbGwtb3BhY2l0eToxO2ZpbGwtcnVsZTpub256ZXJvO3N0cm9rZTojMDAwMDAwO3N0cm9rZS13aWR0aDoyO3N0cm9rZS1saW5lY2FwOmJ1dHQ7c3Ryb2tlLWxpbmVqb2luOm1pdGVyO21hcmtlcjpub25lO21hcmtlci1zdGFydDpub25lO21hcmtlci1taWQ6bm9uZTttYXJrZXItZW5kOm5vbmU7c3Ryb2tlLW1pdGVybGltaXQ6NDtzdHJva2UtZGFzaGFycmF5Om5vbmU7c3Ryb2tlLWRhc2hvZmZzZXQ6MDtzdHJva2Utb3BhY2l0eToxO3Zpc2liaWxpdHk6dmlzaWJsZTtkaXNwbGF5OmlubGluZTtvdmVyZmxvdzp2aXNpYmxlO2VuYWJsZS1iYWNrZ3JvdW5kOmFjY3VtdWxhdGUiCiAgICAgICBpZD0icGF0aDIzODciCiAgICAgICBzb2RpcG9kaTpjeD0iMjAwIgogICAgICAgc29kaXBvZGk6Y3k9IjQzMi4zNjIxOCIKICAgICAgIHNvZGlwb2RpOnJ4PSIxMCIKICAgICAgIHNvZGlwb2RpOnJ5PSIxMCIKICAgICAgIGQ9Ik0gMjEwLDQzMi4zNjIxOCBBIDEwLDEwIDAgMSAxIDE5MCw0MzIuMzYyMTggQSAxMCwxMCAwIDEgMSAyMTAsNDMyLjM2MjE4IHoiCiAgICAgICB0cmFuc2Zvcm09InRyYW5zbGF0ZSgxMCwxMCkiIC8+CiAgICA8cGF0aAogICAgICAgZD0iTSAyMTAsNDMyLjM2MjE4IEEgMTAsMTAgMCAxIDEgMTkwLDQzMi4zNjIxOCBBIDEwLDEwIDAgMSAxIDIxMCw0MzIuMzYyMTggeiIKICAgICAgIHNvZGlwb2RpOnJ5PSIxMCIKICAgICAgIHNvZGlwb2RpOnJ4PSIxMCIKICAgICAgIHNvZGlwb2RpOmN5PSI0MzIuMzYyMTgiCiAgICAgICBzb2RpcG9kaTpjeD0iMjAwIgogICAgICAgaWQ9InBhdGgyMzg5IgogICAgICAgc3R5bGU9Im9wYWNpdHk6MTtmaWxsOiMwMDAwMDA7ZmlsbC1vcGFjaXR5OjE7ZmlsbC1ydWxlOm5vbnplcm87c3Ryb2tlOiMwMDAwMDA7c3Ryb2tlLXdpZHRoOjI7c3Ryb2tlLWxpbmVjYXA6YnV0dDtzdHJva2UtbGluZWpvaW46bWl0ZXI7bWFya2VyOm5vbmU7bWFya2VyLXN0YXJ0Om5vbmU7bWFya2VyLW1pZDpub25lO21hcmtlci1lbmQ6bm9uZTtzdHJva2UtbWl0ZXJsaW1pdDo0O3N0cm9rZS1kYXNoYXJyYXk6bm9uZTtzdHJva2UtZGFzaG9mZnNldDowO3N0cm9rZS1vcGFjaXR5OjE7dmlzaWJpbGl0eTp2aXNpYmxlO2Rpc3BsYXk6aW5saW5lO292ZXJmbG93OnZpc2libGU7ZW5hYmxlLWJhY2tncm91bmQ6YWNjdW11bGF0ZSIKICAgICAgIHNvZGlwb2RpOnR5cGU9ImFyYyIKICAgICAgIHRyYW5zZm9ybT0idHJhbnNsYXRlKDEwLDkwKSIgLz4KICAgIDxwYXRoCiAgICAgICB0cmFuc2Zvcm09InRyYW5zbGF0ZSg5MCwxMCkiCiAgICAgICBzb2RpcG9kaTp0eXBlPSJhcmMiCiAgICAgICBzdHlsZT0ib3BhY2l0eToxO2ZpbGw6IzAwMDAwMDtmaWxsLW9wYWNpdHk6MTtmaWxsLXJ1bGU6bm9uemVybztzdHJva2U6IzAwMDAwMDtzdHJva2Utd2lkdGg6MjtzdHJva2UtbGluZWNhcDpidXR0O3N0cm9rZS1saW5lam9pbjptaXRlcjttYXJrZXI6bm9uZTttYXJrZXItc3RhcnQ6bm9uZTttYXJrZXItbWlkOm5vbmU7bWFya2VyLWVuZDpub25lO3N0cm9rZS1taXRlcmxpbWl0OjQ7c3Ryb2tlLWRhc2hhcnJheTpub25lO3N0cm9rZS1kYXNob2Zmc2V0OjA7c3Ryb2tlLW9wYWNpdHk6MTt2aXNpYmlsaXR5OnZpc2libGU7ZGlzcGxheTppbmxpbmU7b3ZlcmZsb3c6dmlzaWJsZTtlbmFibGUtYmFja2dyb3VuZDphY2N1bXVsYXRlIgogICAgICAgaWQ9InBhdGgyMzkxIgogICAgICAgc29kaXBvZGk6Y3g9IjIwMCIKICAgICAgIHNvZGlwb2RpOmN5PSI0MzIuMzYyMTgiCiAgICAgICBzb2RpcG9kaTpyeD0iMTAiCiAgICAgICBzb2RpcG9kaTpyeT0iMTAiCiAgICAgICBkPSJNIDIxMCw0MzIuMzYyMTggQSAxMCwxMCAwIDEgMSAxOTAsNDMyLjM2MjE4IEEgMTAsMTAgMCAxIDEgMjEwLDQzMi4zNjIxOCB6IiAvPgogICAgPHBhdGgKICAgICAgIGQ9Ik0gMjEwLDQzMi4zNjIxOCBBIDEwLDEwIDAgMSAxIDE5MCw0MzIuMzYyMTggQSAxMCwxMCAwIDEgMSAyMTAsNDMyLjM2MjE4IHoiCiAgICAgICBzb2RpcG9kaTpyeT0iMTAiCiAgICAgICBzb2RpcG9kaTpyeD0iMTAiCiAgICAgICBzb2RpcG9kaTpjeT0iNDMyLjM2MjE4IgogICAgICAgc29kaXBvZGk6Y3g9IjIwMCIKICAgICAgIGlkPSJwYXRoMjM5MyIKICAgICAgIHN0eWxlPSJvcGFjaXR5OjE7ZmlsbDojZmZmZmZmO2ZpbGwtb3BhY2l0eToxO2ZpbGwtcnVsZTpub256ZXJvO3N0cm9rZTojMDAwMDAwO3N0cm9rZS13aWR0aDoyO3N0cm9rZS1saW5lY2FwOmJ1dHQ7c3Ryb2tlLWxpbmVqb2luOm1pdGVyO21hcmtlcjpub25lO21hcmtlci1zdGFydDpub25lO21hcmtlci1taWQ6bm9uZTttYXJrZXItZW5kOm5vbmU7c3Ryb2tlLW1pdGVybGltaXQ6NDtzdHJva2UtZGFzaGFycmF5Om5vbmU7c3Ryb2tlLWRhc2hvZmZzZXQ6MDtzdHJva2Utb3BhY2l0eToxO3Zpc2liaWxpdHk6dmlzaWJsZTtkaXNwbGF5OmlubGluZTtvdmVyZmxvdzp2aXNpYmxlO2VuYWJsZS1iYWNrZ3JvdW5kOmFjY3VtdWxhdGUiCiAgICAgICBzb2RpcG9kaTp0eXBlPSJhcmMiCiAgICAgICB0cmFuc2Zvcm09InRyYW5zbGF0ZSg5MCw5MCkiIC8+CiAgICA8cGF0aAogICAgICAgdHJhbnNmb3JtPSJ0cmFuc2xhdGUoLTUwLDUwKSIKICAgICAgIHNvZGlwb2RpOnR5cGU9ImFyYyIKICAgICAgIHN0eWxlPSJvcGFjaXR5OjE7ZmlsbDojMDAwMDAwO2ZpbGwtb3BhY2l0eToxO2ZpbGwtcnVsZTpub256ZXJvO3N0cm9rZTojMDAwMDAwO3N0cm9rZS13aWR0aDoyO3N0cm9rZS1saW5lY2FwOmJ1dHQ7c3Ryb2tlLWxpbmVqb2luOm1pdGVyO21hcmtlcjpub25lO21hcmtlci1zdGFydDpub25lO21hcmtlci1taWQ6bm9uZTttYXJrZXItZW5kOm5vbmU7c3Ryb2tlLW1pdGVybGltaXQ6NDtzdHJva2UtZGFzaGFycmF5Om5vbmU7c3Ryb2tlLWRhc2hvZmZzZXQ6MDtzdHJva2Utb3BhY2l0eToxO3Zpc2liaWxpdHk6dmlzaWJsZTtkaXNwbGF5OmlubGluZTtvdmVyZmxvdzp2aXNpYmxlO2VuYWJsZS1iYWNrZ3JvdW5kOmFjY3VtdWxhdGUiCiAgICAgICBpZD0icGF0aDIzOTUiCiAgICAgICBzb2RpcG9kaTpjeD0iMjAwIgogICAgICAgc29kaXBvZGk6Y3k9IjQzMi4zNjIxOCIKICAgICAgIHNvZGlwb2RpOnJ4PSIxMCIKICAgICAgIHNvZGlwb2RpOnJ5PSIxMCIKICAgICAgIGQ9Ik0gMjEwLDQzMi4zNjIxOCBBIDEwLDEwIDAgMSAxIDE5MCw0MzIuMzYyMTggQSAxMCwxMCAwIDEgMSAyMTAsNDMyLjM2MjE4IHoiIC8+CiAgICA8cGF0aAogICAgICAgZD0iTSAyMTAsNDMyLjM2MjE4IEEgMTAsMTAgMCAxIDEgMTkwLDQzMi4zNjIxOCBBIDEwLDEwIDAgMSAxIDIxMCw0MzIuMzYyMTggeiIKICAgICAgIHNvZGlwb2RpOnJ5PSIxMCIKICAgICAgIHNvZGlwb2RpOnJ4PSIxMCIKICAgICAgIHNvZGlwb2RpOmN5PSI0MzIuMzYyMTgiCiAgICAgICBzb2RpcG9kaTpjeD0iMjAwIgogICAgICAgaWQ9InBhdGgzMjAwIgogICAgICAgc3R5bGU9Im9wYWNpdHk6MTtmaWxsOiNmZmZmZmY7ZmlsbC1vcGFjaXR5OjE7ZmlsbC1ydWxlOm5vbnplcm87c3Ryb2tlOiMwMDAwMDA7c3Ryb2tlLXdpZHRoOjI7c3Ryb2tlLWxpbmVjYXA6YnV0dDtzdHJva2UtbGluZWpvaW46bWl0ZXI7bWFya2VyOm5vbmU7bWFya2VyLXN0YXJ0Om5vbmU7bWFya2VyLW1pZDpub25lO21hcmtlci1lbmQ6bm9uZTtzdHJva2UtbWl0ZXJsaW1pdDo0O3N0cm9rZS1kYXNoYXJyYXk6bm9uZTtzdHJva2UtZGFzaG9mZnNldDowO3N0cm9rZS1vcGFjaXR5OjE7dmlzaWJpbGl0eTp2aXNpYmxlO2Rpc3BsYXk6aW5saW5lO292ZXJmbG93OnZpc2libGU7ZW5hYmxlLWJhY2tncm91bmQ6YWNjdW11bGF0ZSIKICAgICAgIHNvZGlwb2RpOnR5cGU9ImFyYyIKICAgICAgIHRyYW5zZm9ybT0idHJhbnNsYXRlKDEwLDEwKSIgLz4KICA8L2c+CiAgPGcKICAgICBpbmtzY2FwZTpncm91cG1vZGU9ImxheWVyIgogICAgIGlkPSJsYXllcjMiCiAgICAgaW5rc2NhcGU6bGFiZWw9IkN1dCIKICAgICBzdHlsZT0iZGlzcGxheTppbmxpbmUiCiAgICAgdHJhbnNmb3JtPSJ0cmFuc2xhdGUoLTEyOSwtNDAxLjcwMDAyKSI+CiAgICA8cGF0aAogICAgICAgc3R5bGU9ImZpbGw6bm9uZTtmaWxsLXJ1bGU6ZXZlbm9kZDtzdHJva2U6IzAwMDAwMDtzdHJva2Utd2lkdGg6MztzdHJva2UtbGluZWNhcDpidXR0O3N0cm9rZS1saW5lam9pbjptaXRlcjtzdHJva2UtbWl0ZXJsaW1pdDo0O3N0cm9rZS1kYXNoYXJyYXk6OSwgMztzdHJva2UtZGFzaG9mZnNldDowO3N0cm9rZS1vcGFjaXR5OjEiCiAgICAgICBkPSJNIDMyMSw0NzguMzYyMTggQyAzMTAuNTYzNTIsNDc3LjY3NDg3IDI5OS42ODAyMyw0NzkuOTYxOTkgMjg5LjU3NzY4LDQ3Ni42NjI1OSBDIDI3Ni40Nzk4OCw0NzEuNzE4MSAyNjUuMTAxNjIsNDYyLjYwMDY4IDI1Ny43ODI1MSw0NTAuNjI3MDEgQyAyNDkuOTgxNzgsNDM5LjMzMTA2IDI0NC4wNDQwOSw0MjUuMTkwNzMgMjMxLjEzMDk0LDQxOC44MTUzNyBDIDIxNy43NDU4NSw0MTQuMzUxMSAyMDIuNzcxMjYsNDEwLjc1NDAyIDE4OC45NjA2Nyw0MTUuMzA2MDUgQyAxODAuMzI5NjgsNDE4Ljk4NTUxIDE3Mi43MDE2NSw0MjcuMzI1NDEgMTczLjc4NjY2LDQzNy4yNTczMiBDIDE3My44NzYyNSw0NDkuNDE5NjcgMTc4LjIyODQyLDQ2MS44MzA5NiAxODcuNDMwMyw0NzAuMDk1MDMgQyAxOTMuNjgxNTksNDc3LjM1NTI1IDIwMC41OTIyMyw0ODMuOTczNyAyMDkuMjY2MDIsNDg4LjI0Mzg2IEMgMjE5LjY1Njc5LDQ5NC43MTAzIDIzMC45NTYsNTAwLjgxMTExIDIzNy44MTA2LDUxMS4zNDU5IEMgMjQzLjU3MTEyLDUxOC4xNzI1OCAyNTAuMDE1NjMsNTI0LjcwMTQ4IDI1My42Njg3Miw1MzIuOTg3MzkgQyAyNTguNTE4ODUsNTM3Ljg4NDcxIDI1OS4wMjU0LDU0NS4xMjE2NCAyNjEsNTUxLjM2MjE4IgogICAgICAgaWQ9InBhdGgyNDEzIiAvPgogIDwvZz4KPC9zdmc+Cg=="
     ]
    }
   },
   "cell_type": "markdown",
   "id": "1c00f3ac",
   "metadata": {},
   "source": [
    "## Max-Cut problem\n",
    "\n",
    "### Let us consider a concrete optimization problem: The Max-Cut problem on a graph. We consider a planar graph $G = (V,E)$, having $L$ vertices (or nodes) enumerated by index $i \\in V$ ($i=1,\\dots,L$), and edges being connected pairs of vertices, i.e. $E = \\{ (i,j) | i,j\\in V, i\\ne j\\}$. Each edge has a weigth $W_{ij}>0$. We consider an undirected graph, thus $W_{ij} = W_{ji}$. \n",
    "\n",
    "![Max-cut.svg](attachment:Max-cut.svg)\n",
    "\n",
    "### We define a $\\textit{cut}$ of a graph $G$, as a division of graph nodes into two parts $S$, $T$ such that $V = S\\cup T$. To each node belonging to subgraph $S$, we assign the value $z_i = +1$, while to each node in subgraph $T$, we assign the value of $z_i = -1$. To check if two nodes belong to the same subgraph, we can define a function\n",
    "### \\begin{equation}\n",
    " c_{ij} = \\frac{1-z_iz_j}{2},\n",
    "\\end{equation}\n",
    "### which takes value $0$ for $z_i = z_j$, i.e. if two nodes belong to the same subset, and $1$ otherwise. Next, we can define a $\\textit{value of a cut}$ as\n",
    "### \\begin{equation}\n",
    "  P = \\sum_{i,j} c_{ij}W_{ij}.\n",
    "\\end{equation}\n",
    "### In a Max-Cut problem,  we aim to  $\\textit{cut}$ maximizing $P$. Rewriting $P$ as\n",
    "### \\begin{equation}\n",
    " P = \\sum_{i<j}W_{ij} - \\sum_{i<j}W_{ij}z_iz_j = W - C(\\vec{z}),\n",
    "\\end{equation} \n",
    "### where $W \\equiv \\sum_{i<j}W_{ij}$ is a positive constant number, the Max-Cut problem reduces to minimize the cost function $C(\\vec{z}) = \\sum_{i<j}W_{ij}z_iz_j$, $\\vec{z}=(z_1,\\dots,z_L)$."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "46531675",
   "metadata": {},
   "source": [
    "## Max-Cut problem expressed as a spin-1/2 problem\n",
    "\n",
    "### The Max-Cut problem can be cast to solving a quantum mechanical problem of spin-1/2 system. To the $i$-th node with the value $z_i$ we assign the state of the spin-1/2, i.e.:\n",
    "\n",
    "### \\begin{equation}\n",
    " \\begin{split}\n",
    "    z_i = +1 \\to & |0\\rangle \\equiv |\\uparrow\\rangle \\\\ \n",
    "    z_i = -1 \\to & |1\\rangle \\equiv |\\downarrow\\rangle.\n",
    " \\end{split}\n",
    "\\end{equation}\n",
    "### As such, minimizing the cost function $C(\\vec{z})$ can be considered as fidning the spin configuration which minimizes the $\\textit{cost Hamiltonian}$, defined as:\n",
    "\n",
    "### \\begin{equation}\n",
    " \\hat{H}_C = \\sum_{ik} w_{ik} \\hat{\\sigma}^z_i\\hat{\\sigma}^z_k,\n",
    "\\end{equation}\n",
    "### where $\\hat{\\sigma}^z_i$ is a Pauli-Z matrix acting on $i$-th qubit.\n",
    "\n",
    "### While we consider a planar graph and enumerate spins with index $i$, the problem is equivalent to considering a one-dimensional chain of spins connected adequately via edges.\n",
    "### For simplicity, let's consider a simple two dimensional grid of size $L_x \\times L_y$, where each node of the grid we $(i_x, i_y)$ have one spin-$1/2$, and $i_x = 1,\\dots,L_x$, $i_y = 1,\\dots,L_y$. We assume $w_{ik} = 1$ where two spins are nearest neighbors and $w_{ik}=0$ otherwise.\n",
    "\n",
    "### Next, to each spin at node $(i_x,i_y)$ we can assign a label $i \\equiv (i_y - 1)L_x + i_x$, , i.e. $(i_x, i_y) \\to i$. With such a mapping, our original 2-dimensional problem can be cast to a one-dimensional problem with $L = L_x L_y$ spins-$1/2$, where total Hilbert space has size ${\\cal D} = 2^{L}$.\n",
    "\n",
    "### The Quantum Approximate Optimization Algorithm is based on a specific form of an ansatz for optimal many-body wave function $|\\psi_*\\rangle$ minimizing the cost Hamiltonian $\\hat{H}_C$. The optimal parameters of an ansatz are obtained in an iterative procedure."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "501af93c",
   "metadata": {},
   "source": [
    "## QAOA ansatz\n",
    "\n",
    "### The QAOA ansatz for the many-body wave-function providing optimal spin configuration minimizing the cost Hamiltonian $\\hat{H}_C$ has the following form:\n",
    "\n",
    "### \\begin{equation}\n",
    " |\\psi(\\vec{\\beta},\\vec{\\gamma}\\rangle = U_B(\\beta_1)U_C(\\gamma_1)\\dots U_B(\\beta_p)U_C(\\gamma_p)|s\\rangle,\n",
    "\\end{equation}\n",
    "\n",
    "### where an initial state $|s\\rangle$ is coherent superposition on spin up/down configuration for each qubit. The state $|s\\rangle$ is generated via acting with the Hadamard operator on $i$-th qubit:\n",
    "### \\begin{equation}\n",
    "    |s\\rangle = \\sum_{i=1}^{L} H_i|0\\rangle \\equiv  \\sum_{i=1}^{L} H_i|\\uparrow\\rangle.\n",
    "\\end{equation}\n",
    "\n",
    "### The variational vectors $\\vec{\\beta} = (\\beta_1,\\dots,\\beta_p)$, and $\\vec{\\gamma} = (\\gamma_1,\\dots,\\gamma_p)$ parametrize two unitaries:\n",
    "### \\begin{equation}\n",
    " \\begin{split}\n",
    "   U_B(\\beta_l) & = e^{-i\\beta_l \\hat{H}_B},\\\\\n",
    "   U_C(\\gamma_l) & = e^{-i\\gamma_l \\hat{H}_C},\n",
    " \\end{split}\n",
    "\\end{equation}\n",
    "### where $\\hat{H}_C$ is our cost Hamiltonian, and $\\hat{H}_B = \\sum_{i} \\hat{\\sigma}^x_i$ is called a $\\textit{mixer Hamiltonian}$. The number of variational coeficients $p$ is a hyperparameter for the QAOA.\n",
    "\n",
    "### In the QOAO  we iteratively update parameters $\\vec{\\beta}$, and $\\vec{\\gamma}$ in order to minimize the cost function, i.e. expectation value of the cost Hamiltonian:\n",
    "### \\begin{equation}\n",
    "F(\\vec{\\beta},\\vec{\\gamma}) = \\langle\\psi(\\vec{\\beta},\\vec{\\gamma})|\\hat{H}_C|\\psi(\\vec{\\beta},\\vec{\\gamma})\\rangle.\n",
    "\\end{equation}\n",
    "### The optimization procedure can be done with the help of automatic differentiation.\n",
    "\n",
    "### After finding the optimal set of parameters $\\vec{\\beta}_*$, $\\vec{\\gamma}_*$ defining the many-body wave function\n",
    "### \\begin{equation}\n",
    " |\\psi_*\\rangle\\equiv |\\psi(\\vec{\\beta}_*,\\vec{\\gamma}_*)\\rangle,\n",
    "\\end{equation}\n",
    "### the optimal spin configuration for our problem corresponds to index $k$ of the Fock vector from our basis, for which the quantity\n",
    "### \\begin{equation}\n",
    " |\\langle v_k|\\psi_*\\rangle|^2\n",
    "\\end{equation}\n",
    "### is maximal, where $|v_k\\rangle$ is a Fock vector of the $L$-qubit Hilbert space:\n",
    "### \\begin{equation}\n",
    " \\begin{split}\n",
    "   |v_1\\rangle & = |\\uparrow \\uparrow \\dots \\uparrow \\rangle \\\\\n",
    "   |v_2\\rangle & = |\\uparrow \\uparrow \\dots \\downarrow \\rangle \\\\\n",
    "   & \\vdots \\\\\n",
    "   |v_D\\rangle & = |\\downarrow \\downarrow \\dots \\downarrow \\rangle \\\\\n",
    " \\end{split}\n",
    "\\end{equation}\n",
    "### In the following we denote: $\\uparrow \\equiv 1$, $\\downarrow \\equiv -1$.\n",
    "\n",
    "### Finally, we repeat the optimization procedure $N$ times, and we plot histogram of the optimal spin configuration indexed by $k$. In the limit of $N\\to\\infty$ the most probable spin configuration is the solution to our classical optimization problem."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c95d835f",
   "metadata": {},
   "source": [
    "## QAOA implementation in Pytorch"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "30d6b69f",
   "metadata": {},
   "source": [
    "### Pauli spin chain operators\n",
    "\n",
    "### The starting point is to implement the Pauli spin chain operators. On a spin-$1/2$ chain with $L$ spins, the Pauli spin operators $X_i,Y_i,Z_i$ acting on $i$-th spin are defined as:\n",
    "### \\begin{equation}\n",
    " \\begin{split}\n",
    " X_i  &= \\mathbb{1}_1\\otimes\\dots\\mathbb{1}_{i-1}\\otimes\\hat{\\sigma}^x\\otimes\\mathbb{1}_{i+1}\\dots\\mathbb{1}_{L},\\\\\n",
    " Y_i  &= \\mathbb{1}_1\\otimes\\dots\\mathbb{1}_{i-1}\\otimes\\hat{\\sigma}^y\\otimes\\mathbb{1}_{i+1}\\dots\\mathbb{1}_{L},\\\\\n",
    " Z_i  &= \\mathbb{1}_1\\otimes\\dots\\mathbb{1}_{i-1}\\otimes\\hat{\\sigma}^z\\otimes\\mathbb{1}_{i+1}\\dots\\mathbb{1}_{L},\n",
    " \\end{split}\n",
    "\\end{equation}\n",
    "### where $\\hat{\\sigma}^{x,y,z}$ are $2\\times2$ Pauli operators.\n",
    "\n",
    "### In a similar manner we can define Hadamard matrix $H_i$ acting on $i$-th spin, i.e.\n",
    "### \\begin{equation}\n",
    " H_i = \\mathbb{1}_1\\otimes\\dots\\mathbb{1}_{i-1}\\otimes H\\otimes\\mathbb{1}_{i+1}\\dots\\mathbb{1}_{L},\n",
    "\\end{equation}\n",
    "### where $H = \\frac{\\hat{\\sigma}^x + \\hat{\\sigma}^z}{\\sqrt{2}}$.\n",
    "\n",
    "### Let's import necessary librarires and define Pauli operators, and Hadamard operator, on a spin chain of lenght $L=L_x\\times L_y$. We consider $L_x = 3$, $L_y = 3$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5624c823",
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import torch as pt\n",
    "from torch import matrix_exp as expm\n",
    "from torch.linalg import eigh as eigh\n",
    "import numpy as np\n",
    "\n",
    "def get_linear_index(i_x, i_y):\n",
    "    return  (i_y-1)*Lx + i_x\n",
    "\n",
    "def get_Identity(k):  # returns k-tensor product of the identity operator, ie. Id^k\n",
    "    Id = id_local\n",
    "    for i in range(0, k-1):\n",
    "        Id = pt.kron(Id, id_local)\n",
    "    return Id\n",
    "         \n",
    "def get_string_operator(A, L, i):\n",
    "    Op = A\n",
    "    if(i == 1):\n",
    "        Op = pt.kron(A,get_Identity(L-1))\n",
    "        return Op\n",
    "    if(i == L):\n",
    "        Op = pt.kron(get_Identity(L-1),A)\n",
    "        return Op\n",
    "    if(i>0 and i<L):\n",
    "        Op = pt.kron(get_Identity(i-1), pt.kron(Op, get_Identity(L-i)))\n",
    "        return Op\n",
    "    \n",
    "id_local = pt.tensor([[1.,0],[0,1.]])\n",
    "sigma_x = pt.tensor([[0,1.],[1.,0]])\n",
    "sigma_y = 1j*pt.tensor([[0,-1.],[1.,0]])\n",
    "sigma_z = pt.tensor([[1.,0],[0,-1.]])\n",
    "hadamard = 1.0/pt.sqrt(pt.tensor(2))*pt.tensor([[1,1],[1,-1]])+1j*0    \n",
    "\n",
    "\n",
    "Lx = 3\n",
    "Ly = 3\n",
    "L = Lx*Ly\n",
    "D = 2**L\n",
    "Id = get_string_operator(id_local, L, 1)\n",
    "X = {}\n",
    "Y = {}\n",
    "Z = {}\n",
    "Hadamard = {}\n",
    "\n",
    "\n",
    "for qubit_i_y in range(1, Ly+1):    # Loop over indices on a 2-dimensional grid (i_x,i_y)\n",
    "    for qubit_i_x in range(1,Lx+1): #\n",
    "        qubit_i = get_linear_index(qubit_i_x, qubit_i_y)\n",
    "        X[qubit_i] = get_string_operator(sigma_x, L, qubit_i)        # Define operator X_i acting on spin (i_x,i_y)\n",
    "        Y[qubit_i] = get_string_operator(sigma_y, L, qubit_i)        # Define operator Y_i acting on spin (i_x,i_y)\n",
    "        Z[qubit_i] = get_string_operator(sigma_z, L, qubit_i)        # Define operator Z_i acting on spin (i_x,i_y)\n",
    "        Hadamard[qubit_i] = get_string_operator(hadamard, L, qubit_i)# Define operator H_i acting on spin (i_x,i_y)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5fa6f464",
   "metadata": {},
   "source": [
    "### Operators are expressed in the diagonal basis of the $Z = \\sum_{i=1}^{L} Z_i$ operator, i.e.\n",
    "\n",
    "### \\begin{equation}\n",
    "\\begin{split}\n",
    "   X & = \\sum_{k,l} <v_k|X|v_l>|v_k\\rangle\\langle v_l| \\\\\n",
    "   Y & = \\sum_{k,l} <v_k|Y|v_l>|v_k\\rangle\\langle v_l| \\\\\n",
    "   Z & = \\sum_{k,l} <v_k|Z|v_l>|v_k\\rangle\\langle v_l| \\delta_{k,l},\n",
    "\\end{split}\n",
    "\\end{equation}\n",
    "and\n",
    "\\begin{equation}\n",
    " \\begin{split}\n",
    "   |v_1\\rangle & = |\\uparrow \\uparrow \\dots \\uparrow \\rangle \\\\\n",
    "   |v_2\\rangle & = |\\uparrow \\uparrow \\dots \\downarrow \\rangle \\\\\n",
    "   & \\vdots \\\\\n",
    "   |v_D\\rangle & = |\\downarrow \\downarrow \\dots \\downarrow \\rangle \\\\\n",
    " \\end{split}\n",
    "\\end{equation}\n",
    "\n",
    "### In the following we denote: $\\uparrow \\equiv 1$, $\\downarrow \\equiv -1$.\n",
    "\n",
    "## Let's us construct our Hilbert space:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1174af8d",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "|v_000> = | 1  1  1  1  1  1  1  1  1 >\n",
      "|v_001> = | 1  1  1  1  1  1  1  1 -1 >\n",
      "|v_002> = | 1  1  1  1  1  1  1 -1  1 >\n",
      "|v_003> = | 1  1  1  1  1  1  1 -1 -1 >\n",
      "|v_004> = | 1  1  1  1  1  1 -1  1  1 >\n",
      "|v_005> = | 1  1  1  1  1  1 -1  1 -1 >\n",
      "|v_006> = | 1  1  1  1  1  1 -1 -1  1 >\n",
      "|v_007> = | 1  1  1  1  1  1 -1 -1 -1 >\n",
      "|v_008> = | 1  1  1  1  1 -1  1  1  1 >\n",
      "|v_009> = | 1  1  1  1  1 -1  1  1 -1 >\n",
      "|v_010> = | 1  1  1  1  1 -1  1 -1  1 >\n",
      "|v_011> = | 1  1  1  1  1 -1  1 -1 -1 >\n",
      "|v_012> = | 1  1  1  1  1 -1 -1  1  1 >\n",
      "|v_013> = | 1  1  1  1  1 -1 -1  1 -1 >\n",
      "|v_014> = | 1  1  1  1  1 -1 -1 -1  1 >\n",
      "|v_015> = | 1  1  1  1  1 -1 -1 -1 -1 >\n",
      "|v_016> = | 1  1  1  1 -1  1  1  1  1 >\n",
      "|v_017> = | 1  1  1  1 -1  1  1  1 -1 >\n",
      "|v_018> = | 1  1  1  1 -1  1  1 -1  1 >\n",
      "|v_019> = | 1  1  1  1 -1  1  1 -1 -1 >\n",
      "|v_020> = | 1  1  1  1 -1  1 -1  1  1 >\n",
      "|v_021> = | 1  1  1  1 -1  1 -1  1 -1 >\n",
      "|v_022> = | 1  1  1  1 -1  1 -1 -1  1 >\n",
      "|v_023> = | 1  1  1  1 -1  1 -1 -1 -1 >\n",
      "|v_024> = | 1  1  1  1 -1 -1  1  1  1 >\n",
      "|v_025> = | 1  1  1  1 -1 -1  1  1 -1 >\n",
      "|v_026> = | 1  1  1  1 -1 -1  1 -1  1 >\n",
      "|v_027> = | 1  1  1  1 -1 -1  1 -1 -1 >\n",
      "|v_028> = | 1  1  1  1 -1 -1 -1  1  1 >\n",
      "|v_029> = | 1  1  1  1 -1 -1 -1  1 -1 >\n",
      "|v_030> = | 1  1  1  1 -1 -1 -1 -1  1 >\n",
      "|v_031> = | 1  1  1  1 -1 -1 -1 -1 -1 >\n",
      "|v_032> = | 1  1  1 -1  1  1  1  1  1 >\n",
      "|v_033> = | 1  1  1 -1  1  1  1  1 -1 >\n",
      "|v_034> = | 1  1  1 -1  1  1  1 -1  1 >\n",
      "|v_035> = | 1  1  1 -1  1  1  1 -1 -1 >\n",
      "|v_036> = | 1  1  1 -1  1  1 -1  1  1 >\n",
      "|v_037> = | 1  1  1 -1  1  1 -1  1 -1 >\n",
      "|v_038> = | 1  1  1 -1  1  1 -1 -1  1 >\n",
      "|v_039> = | 1  1  1 -1  1  1 -1 -1 -1 >\n",
      "|v_040> = | 1  1  1 -1  1 -1  1  1  1 >\n",
      "|v_041> = | 1  1  1 -1  1 -1  1  1 -1 >\n",
      "|v_042> = | 1  1  1 -1  1 -1  1 -1  1 >\n",
      "|v_043> = | 1  1  1 -1  1 -1  1 -1 -1 >\n",
      "|v_044> = | 1  1  1 -1  1 -1 -1  1  1 >\n",
      "|v_045> = | 1  1  1 -1  1 -1 -1  1 -1 >\n",
      "|v_046> = | 1  1  1 -1  1 -1 -1 -1  1 >\n",
      "|v_047> = | 1  1  1 -1  1 -1 -1 -1 -1 >\n",
      "|v_048> = | 1  1  1 -1 -1  1  1  1  1 >\n",
      "|v_049> = | 1  1  1 -1 -1  1  1  1 -1 >\n",
      "|v_050> = | 1  1  1 -1 -1  1  1 -1  1 >\n",
      "|v_051> = | 1  1  1 -1 -1  1  1 -1 -1 >\n",
      "|v_052> = | 1  1  1 -1 -1  1 -1  1  1 >\n",
      "|v_053> = | 1  1  1 -1 -1  1 -1  1 -1 >\n",
      "|v_054> = | 1  1  1 -1 -1  1 -1 -1  1 >\n",
      "|v_055> = | 1  1  1 -1 -1  1 -1 -1 -1 >\n",
      "|v_056> = | 1  1  1 -1 -1 -1  1  1  1 >\n",
      "|v_057> = | 1  1  1 -1 -1 -1  1  1 -1 >\n",
      "|v_058> = | 1  1  1 -1 -1 -1  1 -1  1 >\n",
      "|v_059> = | 1  1  1 -1 -1 -1  1 -1 -1 >\n",
      "|v_060> = | 1  1  1 -1 -1 -1 -1  1  1 >\n",
      "|v_061> = | 1  1  1 -1 -1 -1 -1  1 -1 >\n",
      "|v_062> = | 1  1  1 -1 -1 -1 -1 -1  1 >\n",
      "|v_063> = | 1  1  1 -1 -1 -1 -1 -1 -1 >\n",
      "|v_064> = | 1  1 -1  1  1  1  1  1  1 >\n",
      "|v_065> = | 1  1 -1  1  1  1  1  1 -1 >\n",
      "|v_066> = | 1  1 -1  1  1  1  1 -1  1 >\n",
      "|v_067> = | 1  1 -1  1  1  1  1 -1 -1 >\n",
      "|v_068> = | 1  1 -1  1  1  1 -1  1  1 >\n",
      "|v_069> = | 1  1 -1  1  1  1 -1  1 -1 >\n",
      "|v_070> = | 1  1 -1  1  1  1 -1 -1  1 >\n",
      "|v_071> = | 1  1 -1  1  1  1 -1 -1 -1 >\n",
      "|v_072> = | 1  1 -1  1  1 -1  1  1  1 >\n",
      "|v_073> = | 1  1 -1  1  1 -1  1  1 -1 >\n",
      "|v_074> = | 1  1 -1  1  1 -1  1 -1  1 >\n",
      "|v_075> = | 1  1 -1  1  1 -1  1 -1 -1 >\n",
      "|v_076> = | 1  1 -1  1  1 -1 -1  1  1 >\n",
      "|v_077> = | 1  1 -1  1  1 -1 -1  1 -1 >\n",
      "|v_078> = | 1  1 -1  1  1 -1 -1 -1  1 >\n",
      "|v_079> = | 1  1 -1  1  1 -1 -1 -1 -1 >\n",
      "|v_080> = | 1  1 -1  1 -1  1  1  1  1 >\n",
      "|v_081> = | 1  1 -1  1 -1  1  1  1 -1 >\n",
      "|v_082> = | 1  1 -1  1 -1  1  1 -1  1 >\n",
      "|v_083> = | 1  1 -1  1 -1  1  1 -1 -1 >\n",
      "|v_084> = | 1  1 -1  1 -1  1 -1  1  1 >\n",
      "|v_085> = | 1  1 -1  1 -1  1 -1  1 -1 >\n",
      "|v_086> = | 1  1 -1  1 -1  1 -1 -1  1 >\n",
      "|v_087> = | 1  1 -1  1 -1  1 -1 -1 -1 >\n",
      "|v_088> = | 1  1 -1  1 -1 -1  1  1  1 >\n",
      "|v_089> = | 1  1 -1  1 -1 -1  1  1 -1 >\n",
      "|v_090> = | 1  1 -1  1 -1 -1  1 -1  1 >\n",
      "|v_091> = | 1  1 -1  1 -1 -1  1 -1 -1 >\n",
      "|v_092> = | 1  1 -1  1 -1 -1 -1  1  1 >\n",
      "|v_093> = | 1  1 -1  1 -1 -1 -1  1 -1 >\n",
      "|v_094> = | 1  1 -1  1 -1 -1 -1 -1  1 >\n",
      "|v_095> = | 1  1 -1  1 -1 -1 -1 -1 -1 >\n",
      "|v_096> = | 1  1 -1 -1  1  1  1  1  1 >\n",
      "|v_097> = | 1  1 -1 -1  1  1  1  1 -1 >\n",
      "|v_098> = | 1  1 -1 -1  1  1  1 -1  1 >\n",
      "|v_099> = | 1  1 -1 -1  1  1  1 -1 -1 >\n",
      "|v_100> = | 1  1 -1 -1  1  1 -1  1  1 >\n",
      "|v_101> = | 1  1 -1 -1  1  1 -1  1 -1 >\n",
      "|v_102> = | 1  1 -1 -1  1  1 -1 -1  1 >\n",
      "|v_103> = | 1  1 -1 -1  1  1 -1 -1 -1 >\n",
      "|v_104> = | 1  1 -1 -1  1 -1  1  1  1 >\n",
      "|v_105> = | 1  1 -1 -1  1 -1  1  1 -1 >\n",
      "|v_106> = | 1  1 -1 -1  1 -1  1 -1  1 >\n",
      "|v_107> = | 1  1 -1 -1  1 -1  1 -1 -1 >\n",
      "|v_108> = | 1  1 -1 -1  1 -1 -1  1  1 >\n",
      "|v_109> = | 1  1 -1 -1  1 -1 -1  1 -1 >\n",
      "|v_110> = | 1  1 -1 -1  1 -1 -1 -1  1 >\n",
      "|v_111> = | 1  1 -1 -1  1 -1 -1 -1 -1 >\n",
      "|v_112> = | 1  1 -1 -1 -1  1  1  1  1 >\n",
      "|v_113> = | 1  1 -1 -1 -1  1  1  1 -1 >\n",
      "|v_114> = | 1  1 -1 -1 -1  1  1 -1  1 >\n",
      "|v_115> = | 1  1 -1 -1 -1  1  1 -1 -1 >\n",
      "|v_116> = | 1  1 -1 -1 -1  1 -1  1  1 >\n",
      "|v_117> = | 1  1 -1 -1 -1  1 -1  1 -1 >\n",
      "|v_118> = | 1  1 -1 -1 -1  1 -1 -1  1 >\n",
      "|v_119> = | 1  1 -1 -1 -1  1 -1 -1 -1 >\n",
      "|v_120> = | 1  1 -1 -1 -1 -1  1  1  1 >\n",
      "|v_121> = | 1  1 -1 -1 -1 -1  1  1 -1 >\n",
      "|v_122> = | 1  1 -1 -1 -1 -1  1 -1  1 >\n",
      "|v_123> = | 1  1 -1 -1 -1 -1  1 -1 -1 >\n",
      "|v_124> = | 1  1 -1 -1 -1 -1 -1  1  1 >\n",
      "|v_125> = | 1  1 -1 -1 -1 -1 -1  1 -1 >\n",
      "|v_126> = | 1  1 -1 -1 -1 -1 -1 -1  1 >\n",
      "|v_127> = | 1  1 -1 -1 -1 -1 -1 -1 -1 >\n",
      "|v_128> = | 1 -1  1  1  1  1  1  1  1 >\n",
      "|v_129> = | 1 -1  1  1  1  1  1  1 -1 >\n",
      "|v_130> = | 1 -1  1  1  1  1  1 -1  1 >\n",
      "|v_131> = | 1 -1  1  1  1  1  1 -1 -1 >\n",
      "|v_132> = | 1 -1  1  1  1  1 -1  1  1 >\n",
      "|v_133> = | 1 -1  1  1  1  1 -1  1 -1 >\n",
      "|v_134> = | 1 -1  1  1  1  1 -1 -1  1 >\n",
      "|v_135> = | 1 -1  1  1  1  1 -1 -1 -1 >\n",
      "|v_136> = | 1 -1  1  1  1 -1  1  1  1 >\n",
      "|v_137> = | 1 -1  1  1  1 -1  1  1 -1 >\n",
      "|v_138> = | 1 -1  1  1  1 -1  1 -1  1 >\n",
      "|v_139> = | 1 -1  1  1  1 -1  1 -1 -1 >\n",
      "|v_140> = | 1 -1  1  1  1 -1 -1  1  1 >\n",
      "|v_141> = | 1 -1  1  1  1 -1 -1  1 -1 >\n",
      "|v_142> = | 1 -1  1  1  1 -1 -1 -1  1 >\n",
      "|v_143> = | 1 -1  1  1  1 -1 -1 -1 -1 >\n",
      "|v_144> = | 1 -1  1  1 -1  1  1  1  1 >\n",
      "|v_145> = | 1 -1  1  1 -1  1  1  1 -1 >\n",
      "|v_146> = | 1 -1  1  1 -1  1  1 -1  1 >\n",
      "|v_147> = | 1 -1  1  1 -1  1  1 -1 -1 >\n",
      "|v_148> = | 1 -1  1  1 -1  1 -1  1  1 >\n",
      "|v_149> = | 1 -1  1  1 -1  1 -1  1 -1 >\n",
      "|v_150> = | 1 -1  1  1 -1  1 -1 -1  1 >\n",
      "|v_151> = | 1 -1  1  1 -1  1 -1 -1 -1 >\n",
      "|v_152> = | 1 -1  1  1 -1 -1  1  1  1 >\n",
      "|v_153> = | 1 -1  1  1 -1 -1  1  1 -1 >\n",
      "|v_154> = | 1 -1  1  1 -1 -1  1 -1  1 >\n",
      "|v_155> = | 1 -1  1  1 -1 -1  1 -1 -1 >\n",
      "|v_156> = | 1 -1  1  1 -1 -1 -1  1  1 >\n",
      "|v_157> = | 1 -1  1  1 -1 -1 -1  1 -1 >\n",
      "|v_158> = | 1 -1  1  1 -1 -1 -1 -1  1 >\n",
      "|v_159> = | 1 -1  1  1 -1 -1 -1 -1 -1 >\n",
      "|v_160> = | 1 -1  1 -1  1  1  1  1  1 >\n",
      "|v_161> = | 1 -1  1 -1  1  1  1  1 -1 >\n",
      "|v_162> = | 1 -1  1 -1  1  1  1 -1  1 >\n",
      "|v_163> = | 1 -1  1 -1  1  1  1 -1 -1 >\n",
      "|v_164> = | 1 -1  1 -1  1  1 -1  1  1 >\n",
      "|v_165> = | 1 -1  1 -1  1  1 -1  1 -1 >\n",
      "|v_166> = | 1 -1  1 -1  1  1 -1 -1  1 >\n",
      "|v_167> = | 1 -1  1 -1  1  1 -1 -1 -1 >\n",
      "|v_168> = | 1 -1  1 -1  1 -1  1  1  1 >\n",
      "|v_169> = | 1 -1  1 -1  1 -1  1  1 -1 >\n",
      "|v_170> = | 1 -1  1 -1  1 -1  1 -1  1 >\n",
      "|v_171> = | 1 -1  1 -1  1 -1  1 -1 -1 >\n",
      "|v_172> = | 1 -1  1 -1  1 -1 -1  1  1 >\n",
      "|v_173> = | 1 -1  1 -1  1 -1 -1  1 -1 >\n",
      "|v_174> = | 1 -1  1 -1  1 -1 -1 -1  1 >\n",
      "|v_175> = | 1 -1  1 -1  1 -1 -1 -1 -1 >\n",
      "|v_176> = | 1 -1  1 -1 -1  1  1  1  1 >\n",
      "|v_177> = | 1 -1  1 -1 -1  1  1  1 -1 >\n",
      "|v_178> = | 1 -1  1 -1 -1  1  1 -1  1 >\n",
      "|v_179> = | 1 -1  1 -1 -1  1  1 -1 -1 >\n",
      "|v_180> = | 1 -1  1 -1 -1  1 -1  1  1 >\n",
      "|v_181> = | 1 -1  1 -1 -1  1 -1  1 -1 >\n",
      "|v_182> = | 1 -1  1 -1 -1  1 -1 -1  1 >\n",
      "|v_183> = | 1 -1  1 -1 -1  1 -1 -1 -1 >\n",
      "|v_184> = | 1 -1  1 -1 -1 -1  1  1  1 >\n",
      "|v_185> = | 1 -1  1 -1 -1 -1  1  1 -1 >\n",
      "|v_186> = | 1 -1  1 -1 -1 -1  1 -1  1 >\n",
      "|v_187> = | 1 -1  1 -1 -1 -1  1 -1 -1 >\n",
      "|v_188> = | 1 -1  1 -1 -1 -1 -1  1  1 >\n",
      "|v_189> = | 1 -1  1 -1 -1 -1 -1  1 -1 >\n",
      "|v_190> = | 1 -1  1 -1 -1 -1 -1 -1  1 >\n",
      "|v_191> = | 1 -1  1 -1 -1 -1 -1 -1 -1 >\n",
      "|v_192> = | 1 -1 -1  1  1  1  1  1  1 >\n",
      "|v_193> = | 1 -1 -1  1  1  1  1  1 -1 >\n",
      "|v_194> = | 1 -1 -1  1  1  1  1 -1  1 >\n",
      "|v_195> = | 1 -1 -1  1  1  1  1 -1 -1 >\n",
      "|v_196> = | 1 -1 -1  1  1  1 -1  1  1 >\n",
      "|v_197> = | 1 -1 -1  1  1  1 -1  1 -1 >\n",
      "|v_198> = | 1 -1 -1  1  1  1 -1 -1  1 >\n",
      "|v_199> = | 1 -1 -1  1  1  1 -1 -1 -1 >\n",
      "|v_200> = | 1 -1 -1  1  1 -1  1  1  1 >\n",
      "|v_201> = | 1 -1 -1  1  1 -1  1  1 -1 >\n",
      "|v_202> = | 1 -1 -1  1  1 -1  1 -1  1 >\n",
      "|v_203> = | 1 -1 -1  1  1 -1  1 -1 -1 >\n",
      "|v_204> = | 1 -1 -1  1  1 -1 -1  1  1 >\n",
      "|v_205> = | 1 -1 -1  1  1 -1 -1  1 -1 >\n",
      "|v_206> = | 1 -1 -1  1  1 -1 -1 -1  1 >\n",
      "|v_207> = | 1 -1 -1  1  1 -1 -1 -1 -1 >\n",
      "|v_208> = | 1 -1 -1  1 -1  1  1  1  1 >\n",
      "|v_209> = | 1 -1 -1  1 -1  1  1  1 -1 >\n",
      "|v_210> = | 1 -1 -1  1 -1  1  1 -1  1 >\n",
      "|v_211> = | 1 -1 -1  1 -1  1  1 -1 -1 >\n",
      "|v_212> = | 1 -1 -1  1 -1  1 -1  1  1 >\n",
      "|v_213> = | 1 -1 -1  1 -1  1 -1  1 -1 >\n",
      "|v_214> = | 1 -1 -1  1 -1  1 -1 -1  1 >\n",
      "|v_215> = | 1 -1 -1  1 -1  1 -1 -1 -1 >\n",
      "|v_216> = | 1 -1 -1  1 -1 -1  1  1  1 >\n",
      "|v_217> = | 1 -1 -1  1 -1 -1  1  1 -1 >\n",
      "|v_218> = | 1 -1 -1  1 -1 -1  1 -1  1 >\n",
      "|v_219> = | 1 -1 -1  1 -1 -1  1 -1 -1 >\n",
      "|v_220> = | 1 -1 -1  1 -1 -1 -1  1  1 >\n",
      "|v_221> = | 1 -1 -1  1 -1 -1 -1  1 -1 >\n",
      "|v_222> = | 1 -1 -1  1 -1 -1 -1 -1  1 >\n",
      "|v_223> = | 1 -1 -1  1 -1 -1 -1 -1 -1 >\n",
      "|v_224> = | 1 -1 -1 -1  1  1  1  1  1 >\n",
      "|v_225> = | 1 -1 -1 -1  1  1  1  1 -1 >\n",
      "|v_226> = | 1 -1 -1 -1  1  1  1 -1  1 >\n",
      "|v_227> = | 1 -1 -1 -1  1  1  1 -1 -1 >\n",
      "|v_228> = | 1 -1 -1 -1  1  1 -1  1  1 >\n",
      "|v_229> = | 1 -1 -1 -1  1  1 -1  1 -1 >\n",
      "|v_230> = | 1 -1 -1 -1  1  1 -1 -1  1 >\n",
      "|v_231> = | 1 -1 -1 -1  1  1 -1 -1 -1 >\n",
      "|v_232> = | 1 -1 -1 -1  1 -1  1  1  1 >\n",
      "|v_233> = | 1 -1 -1 -1  1 -1  1  1 -1 >\n",
      "|v_234> = | 1 -1 -1 -1  1 -1  1 -1  1 >\n",
      "|v_235> = | 1 -1 -1 -1  1 -1  1 -1 -1 >\n",
      "|v_236> = | 1 -1 -1 -1  1 -1 -1  1  1 >\n",
      "|v_237> = | 1 -1 -1 -1  1 -1 -1  1 -1 >\n",
      "|v_238> = | 1 -1 -1 -1  1 -1 -1 -1  1 >\n",
      "|v_239> = | 1 -1 -1 -1  1 -1 -1 -1 -1 >\n",
      "|v_240> = | 1 -1 -1 -1 -1  1  1  1  1 >\n",
      "|v_241> = | 1 -1 -1 -1 -1  1  1  1 -1 >\n",
      "|v_242> = | 1 -1 -1 -1 -1  1  1 -1  1 >\n",
      "|v_243> = | 1 -1 -1 -1 -1  1  1 -1 -1 >\n",
      "|v_244> = | 1 -1 -1 -1 -1  1 -1  1  1 >\n",
      "|v_245> = | 1 -1 -1 -1 -1  1 -1  1 -1 >\n",
      "|v_246> = | 1 -1 -1 -1 -1  1 -1 -1  1 >\n",
      "|v_247> = | 1 -1 -1 -1 -1  1 -1 -1 -1 >\n",
      "|v_248> = | 1 -1 -1 -1 -1 -1  1  1  1 >\n",
      "|v_249> = | 1 -1 -1 -1 -1 -1  1  1 -1 >\n",
      "|v_250> = | 1 -1 -1 -1 -1 -1  1 -1  1 >\n",
      "|v_251> = | 1 -1 -1 -1 -1 -1  1 -1 -1 >\n",
      "|v_252> = | 1 -1 -1 -1 -1 -1 -1  1  1 >\n",
      "|v_253> = | 1 -1 -1 -1 -1 -1 -1  1 -1 >\n",
      "|v_254> = | 1 -1 -1 -1 -1 -1 -1 -1  1 >\n",
      "|v_255> = | 1 -1 -1 -1 -1 -1 -1 -1 -1 >\n",
      "|v_256> = |-1  1  1  1  1  1  1  1  1 >\n",
      "|v_257> = |-1  1  1  1  1  1  1  1 -1 >\n",
      "|v_258> = |-1  1  1  1  1  1  1 -1  1 >\n",
      "|v_259> = |-1  1  1  1  1  1  1 -1 -1 >\n",
      "|v_260> = |-1  1  1  1  1  1 -1  1  1 >\n",
      "|v_261> = |-1  1  1  1  1  1 -1  1 -1 >\n",
      "|v_262> = |-1  1  1  1  1  1 -1 -1  1 >\n",
      "|v_263> = |-1  1  1  1  1  1 -1 -1 -1 >\n",
      "|v_264> = |-1  1  1  1  1 -1  1  1  1 >\n",
      "|v_265> = |-1  1  1  1  1 -1  1  1 -1 >\n",
      "|v_266> = |-1  1  1  1  1 -1  1 -1  1 >\n",
      "|v_267> = |-1  1  1  1  1 -1  1 -1 -1 >\n",
      "|v_268> = |-1  1  1  1  1 -1 -1  1  1 >\n",
      "|v_269> = |-1  1  1  1  1 -1 -1  1 -1 >\n",
      "|v_270> = |-1  1  1  1  1 -1 -1 -1  1 >\n",
      "|v_271> = |-1  1  1  1  1 -1 -1 -1 -1 >\n",
      "|v_272> = |-1  1  1  1 -1  1  1  1  1 >\n",
      "|v_273> = |-1  1  1  1 -1  1  1  1 -1 >\n",
      "|v_274> = |-1  1  1  1 -1  1  1 -1  1 >\n",
      "|v_275> = |-1  1  1  1 -1  1  1 -1 -1 >\n",
      "|v_276> = |-1  1  1  1 -1  1 -1  1  1 >\n",
      "|v_277> = |-1  1  1  1 -1  1 -1  1 -1 >\n",
      "|v_278> = |-1  1  1  1 -1  1 -1 -1  1 >\n",
      "|v_279> = |-1  1  1  1 -1  1 -1 -1 -1 >\n",
      "|v_280> = |-1  1  1  1 -1 -1  1  1  1 >\n",
      "|v_281> = |-1  1  1  1 -1 -1  1  1 -1 >\n",
      "|v_282> = |-1  1  1  1 -1 -1  1 -1  1 >\n",
      "|v_283> = |-1  1  1  1 -1 -1  1 -1 -1 >\n",
      "|v_284> = |-1  1  1  1 -1 -1 -1  1  1 >\n",
      "|v_285> = |-1  1  1  1 -1 -1 -1  1 -1 >\n",
      "|v_286> = |-1  1  1  1 -1 -1 -1 -1  1 >\n",
      "|v_287> = |-1  1  1  1 -1 -1 -1 -1 -1 >\n",
      "|v_288> = |-1  1  1 -1  1  1  1  1  1 >\n",
      "|v_289> = |-1  1  1 -1  1  1  1  1 -1 >\n",
      "|v_290> = |-1  1  1 -1  1  1  1 -1  1 >\n",
      "|v_291> = |-1  1  1 -1  1  1  1 -1 -1 >\n",
      "|v_292> = |-1  1  1 -1  1  1 -1  1  1 >\n",
      "|v_293> = |-1  1  1 -1  1  1 -1  1 -1 >\n",
      "|v_294> = |-1  1  1 -1  1  1 -1 -1  1 >\n",
      "|v_295> = |-1  1  1 -1  1  1 -1 -1 -1 >\n",
      "|v_296> = |-1  1  1 -1  1 -1  1  1  1 >\n",
      "|v_297> = |-1  1  1 -1  1 -1  1  1 -1 >\n",
      "|v_298> = |-1  1  1 -1  1 -1  1 -1  1 >\n",
      "|v_299> = |-1  1  1 -1  1 -1  1 -1 -1 >\n",
      "|v_300> = |-1  1  1 -1  1 -1 -1  1  1 >\n",
      "|v_301> = |-1  1  1 -1  1 -1 -1  1 -1 >\n",
      "|v_302> = |-1  1  1 -1  1 -1 -1 -1  1 >\n",
      "|v_303> = |-1  1  1 -1  1 -1 -1 -1 -1 >\n",
      "|v_304> = |-1  1  1 -1 -1  1  1  1  1 >\n",
      "|v_305> = |-1  1  1 -1 -1  1  1  1 -1 >\n",
      "|v_306> = |-1  1  1 -1 -1  1  1 -1  1 >\n",
      "|v_307> = |-1  1  1 -1 -1  1  1 -1 -1 >\n",
      "|v_308> = |-1  1  1 -1 -1  1 -1  1  1 >\n",
      "|v_309> = |-1  1  1 -1 -1  1 -1  1 -1 >\n",
      "|v_310> = |-1  1  1 -1 -1  1 -1 -1  1 >\n",
      "|v_311> = |-1  1  1 -1 -1  1 -1 -1 -1 >\n",
      "|v_312> = |-1  1  1 -1 -1 -1  1  1  1 >\n",
      "|v_313> = |-1  1  1 -1 -1 -1  1  1 -1 >\n",
      "|v_314> = |-1  1  1 -1 -1 -1  1 -1  1 >\n",
      "|v_315> = |-1  1  1 -1 -1 -1  1 -1 -1 >\n",
      "|v_316> = |-1  1  1 -1 -1 -1 -1  1  1 >\n",
      "|v_317> = |-1  1  1 -1 -1 -1 -1  1 -1 >\n",
      "|v_318> = |-1  1  1 -1 -1 -1 -1 -1  1 >\n",
      "|v_319> = |-1  1  1 -1 -1 -1 -1 -1 -1 >\n",
      "|v_320> = |-1  1 -1  1  1  1  1  1  1 >\n",
      "|v_321> = |-1  1 -1  1  1  1  1  1 -1 >\n",
      "|v_322> = |-1  1 -1  1  1  1  1 -1  1 >\n",
      "|v_323> = |-1  1 -1  1  1  1  1 -1 -1 >\n",
      "|v_324> = |-1  1 -1  1  1  1 -1  1  1 >\n",
      "|v_325> = |-1  1 -1  1  1  1 -1  1 -1 >\n",
      "|v_326> = |-1  1 -1  1  1  1 -1 -1  1 >\n",
      "|v_327> = |-1  1 -1  1  1  1 -1 -1 -1 >\n",
      "|v_328> = |-1  1 -1  1  1 -1  1  1  1 >\n",
      "|v_329> = |-1  1 -1  1  1 -1  1  1 -1 >\n",
      "|v_330> = |-1  1 -1  1  1 -1  1 -1  1 >\n",
      "|v_331> = |-1  1 -1  1  1 -1  1 -1 -1 >\n",
      "|v_332> = |-1  1 -1  1  1 -1 -1  1  1 >\n",
      "|v_333> = |-1  1 -1  1  1 -1 -1  1 -1 >\n",
      "|v_334> = |-1  1 -1  1  1 -1 -1 -1  1 >\n",
      "|v_335> = |-1  1 -1  1  1 -1 -1 -1 -1 >\n",
      "|v_336> = |-1  1 -1  1 -1  1  1  1  1 >\n",
      "|v_337> = |-1  1 -1  1 -1  1  1  1 -1 >\n",
      "|v_338> = |-1  1 -1  1 -1  1  1 -1  1 >\n",
      "|v_339> = |-1  1 -1  1 -1  1  1 -1 -1 >\n",
      "|v_340> = |-1  1 -1  1 -1  1 -1  1  1 >\n",
      "|v_341> = |-1  1 -1  1 -1  1 -1  1 -1 >\n",
      "|v_342> = |-1  1 -1  1 -1  1 -1 -1  1 >\n",
      "|v_343> = |-1  1 -1  1 -1  1 -1 -1 -1 >\n",
      "|v_344> = |-1  1 -1  1 -1 -1  1  1  1 >\n",
      "|v_345> = |-1  1 -1  1 -1 -1  1  1 -1 >\n",
      "|v_346> = |-1  1 -1  1 -1 -1  1 -1  1 >\n",
      "|v_347> = |-1  1 -1  1 -1 -1  1 -1 -1 >\n",
      "|v_348> = |-1  1 -1  1 -1 -1 -1  1  1 >\n",
      "|v_349> = |-1  1 -1  1 -1 -1 -1  1 -1 >\n",
      "|v_350> = |-1  1 -1  1 -1 -1 -1 -1  1 >\n",
      "|v_351> = |-1  1 -1  1 -1 -1 -1 -1 -1 >\n",
      "|v_352> = |-1  1 -1 -1  1  1  1  1  1 >\n",
      "|v_353> = |-1  1 -1 -1  1  1  1  1 -1 >\n",
      "|v_354> = |-1  1 -1 -1  1  1  1 -1  1 >\n",
      "|v_355> = |-1  1 -1 -1  1  1  1 -1 -1 >\n",
      "|v_356> = |-1  1 -1 -1  1  1 -1  1  1 >\n",
      "|v_357> = |-1  1 -1 -1  1  1 -1  1 -1 >\n",
      "|v_358> = |-1  1 -1 -1  1  1 -1 -1  1 >\n",
      "|v_359> = |-1  1 -1 -1  1  1 -1 -1 -1 >\n",
      "|v_360> = |-1  1 -1 -1  1 -1  1  1  1 >\n",
      "|v_361> = |-1  1 -1 -1  1 -1  1  1 -1 >\n",
      "|v_362> = |-1  1 -1 -1  1 -1  1 -1  1 >\n",
      "|v_363> = |-1  1 -1 -1  1 -1  1 -1 -1 >\n",
      "|v_364> = |-1  1 -1 -1  1 -1 -1  1  1 >\n",
      "|v_365> = |-1  1 -1 -1  1 -1 -1  1 -1 >\n",
      "|v_366> = |-1  1 -1 -1  1 -1 -1 -1  1 >\n",
      "|v_367> = |-1  1 -1 -1  1 -1 -1 -1 -1 >\n",
      "|v_368> = |-1  1 -1 -1 -1  1  1  1  1 >\n",
      "|v_369> = |-1  1 -1 -1 -1  1  1  1 -1 >\n",
      "|v_370> = |-1  1 -1 -1 -1  1  1 -1  1 >\n",
      "|v_371> = |-1  1 -1 -1 -1  1  1 -1 -1 >\n",
      "|v_372> = |-1  1 -1 -1 -1  1 -1  1  1 >\n",
      "|v_373> = |-1  1 -1 -1 -1  1 -1  1 -1 >\n",
      "|v_374> = |-1  1 -1 -1 -1  1 -1 -1  1 >\n",
      "|v_375> = |-1  1 -1 -1 -1  1 -1 -1 -1 >\n",
      "|v_376> = |-1  1 -1 -1 -1 -1  1  1  1 >\n",
      "|v_377> = |-1  1 -1 -1 -1 -1  1  1 -1 >\n",
      "|v_378> = |-1  1 -1 -1 -1 -1  1 -1  1 >\n",
      "|v_379> = |-1  1 -1 -1 -1 -1  1 -1 -1 >\n",
      "|v_380> = |-1  1 -1 -1 -1 -1 -1  1  1 >\n",
      "|v_381> = |-1  1 -1 -1 -1 -1 -1  1 -1 >\n",
      "|v_382> = |-1  1 -1 -1 -1 -1 -1 -1  1 >\n",
      "|v_383> = |-1  1 -1 -1 -1 -1 -1 -1 -1 >\n",
      "|v_384> = |-1 -1  1  1  1  1  1  1  1 >\n",
      "|v_385> = |-1 -1  1  1  1  1  1  1 -1 >\n",
      "|v_386> = |-1 -1  1  1  1  1  1 -1  1 >\n",
      "|v_387> = |-1 -1  1  1  1  1  1 -1 -1 >\n",
      "|v_388> = |-1 -1  1  1  1  1 -1  1  1 >\n",
      "|v_389> = |-1 -1  1  1  1  1 -1  1 -1 >\n",
      "|v_390> = |-1 -1  1  1  1  1 -1 -1  1 >\n",
      "|v_391> = |-1 -1  1  1  1  1 -1 -1 -1 >\n",
      "|v_392> = |-1 -1  1  1  1 -1  1  1  1 >\n",
      "|v_393> = |-1 -1  1  1  1 -1  1  1 -1 >\n",
      "|v_394> = |-1 -1  1  1  1 -1  1 -1  1 >\n",
      "|v_395> = |-1 -1  1  1  1 -1  1 -1 -1 >\n",
      "|v_396> = |-1 -1  1  1  1 -1 -1  1  1 >\n",
      "|v_397> = |-1 -1  1  1  1 -1 -1  1 -1 >\n",
      "|v_398> = |-1 -1  1  1  1 -1 -1 -1  1 >\n",
      "|v_399> = |-1 -1  1  1  1 -1 -1 -1 -1 >\n",
      "|v_400> = |-1 -1  1  1 -1  1  1  1  1 >\n",
      "|v_401> = |-1 -1  1  1 -1  1  1  1 -1 >\n",
      "|v_402> = |-1 -1  1  1 -1  1  1 -1  1 >\n",
      "|v_403> = |-1 -1  1  1 -1  1  1 -1 -1 >\n",
      "|v_404> = |-1 -1  1  1 -1  1 -1  1  1 >\n",
      "|v_405> = |-1 -1  1  1 -1  1 -1  1 -1 >\n",
      "|v_406> = |-1 -1  1  1 -1  1 -1 -1  1 >\n",
      "|v_407> = |-1 -1  1  1 -1  1 -1 -1 -1 >\n",
      "|v_408> = |-1 -1  1  1 -1 -1  1  1  1 >\n",
      "|v_409> = |-1 -1  1  1 -1 -1  1  1 -1 >\n",
      "|v_410> = |-1 -1  1  1 -1 -1  1 -1  1 >\n",
      "|v_411> = |-1 -1  1  1 -1 -1  1 -1 -1 >\n",
      "|v_412> = |-1 -1  1  1 -1 -1 -1  1  1 >\n",
      "|v_413> = |-1 -1  1  1 -1 -1 -1  1 -1 >\n",
      "|v_414> = |-1 -1  1  1 -1 -1 -1 -1  1 >\n",
      "|v_415> = |-1 -1  1  1 -1 -1 -1 -1 -1 >\n",
      "|v_416> = |-1 -1  1 -1  1  1  1  1  1 >\n",
      "|v_417> = |-1 -1  1 -1  1  1  1  1 -1 >\n",
      "|v_418> = |-1 -1  1 -1  1  1  1 -1  1 >\n",
      "|v_419> = |-1 -1  1 -1  1  1  1 -1 -1 >\n",
      "|v_420> = |-1 -1  1 -1  1  1 -1  1  1 >\n",
      "|v_421> = |-1 -1  1 -1  1  1 -1  1 -1 >\n",
      "|v_422> = |-1 -1  1 -1  1  1 -1 -1  1 >\n",
      "|v_423> = |-1 -1  1 -1  1  1 -1 -1 -1 >\n",
      "|v_424> = |-1 -1  1 -1  1 -1  1  1  1 >\n",
      "|v_425> = |-1 -1  1 -1  1 -1  1  1 -1 >\n",
      "|v_426> = |-1 -1  1 -1  1 -1  1 -1  1 >\n",
      "|v_427> = |-1 -1  1 -1  1 -1  1 -1 -1 >\n",
      "|v_428> = |-1 -1  1 -1  1 -1 -1  1  1 >\n",
      "|v_429> = |-1 -1  1 -1  1 -1 -1  1 -1 >\n",
      "|v_430> = |-1 -1  1 -1  1 -1 -1 -1  1 >\n",
      "|v_431> = |-1 -1  1 -1  1 -1 -1 -1 -1 >\n",
      "|v_432> = |-1 -1  1 -1 -1  1  1  1  1 >\n",
      "|v_433> = |-1 -1  1 -1 -1  1  1  1 -1 >\n",
      "|v_434> = |-1 -1  1 -1 -1  1  1 -1  1 >\n",
      "|v_435> = |-1 -1  1 -1 -1  1  1 -1 -1 >\n",
      "|v_436> = |-1 -1  1 -1 -1  1 -1  1  1 >\n",
      "|v_437> = |-1 -1  1 -1 -1  1 -1  1 -1 >\n",
      "|v_438> = |-1 -1  1 -1 -1  1 -1 -1  1 >\n",
      "|v_439> = |-1 -1  1 -1 -1  1 -1 -1 -1 >\n",
      "|v_440> = |-1 -1  1 -1 -1 -1  1  1  1 >\n",
      "|v_441> = |-1 -1  1 -1 -1 -1  1  1 -1 >\n",
      "|v_442> = |-1 -1  1 -1 -1 -1  1 -1  1 >\n",
      "|v_443> = |-1 -1  1 -1 -1 -1  1 -1 -1 >\n",
      "|v_444> = |-1 -1  1 -1 -1 -1 -1  1  1 >\n",
      "|v_445> = |-1 -1  1 -1 -1 -1 -1  1 -1 >\n",
      "|v_446> = |-1 -1  1 -1 -1 -1 -1 -1  1 >\n",
      "|v_447> = |-1 -1  1 -1 -1 -1 -1 -1 -1 >\n",
      "|v_448> = |-1 -1 -1  1  1  1  1  1  1 >\n",
      "|v_449> = |-1 -1 -1  1  1  1  1  1 -1 >\n",
      "|v_450> = |-1 -1 -1  1  1  1  1 -1  1 >\n",
      "|v_451> = |-1 -1 -1  1  1  1  1 -1 -1 >\n",
      "|v_452> = |-1 -1 -1  1  1  1 -1  1  1 >\n",
      "|v_453> = |-1 -1 -1  1  1  1 -1  1 -1 >\n",
      "|v_454> = |-1 -1 -1  1  1  1 -1 -1  1 >\n",
      "|v_455> = |-1 -1 -1  1  1  1 -1 -1 -1 >\n",
      "|v_456> = |-1 -1 -1  1  1 -1  1  1  1 >\n",
      "|v_457> = |-1 -1 -1  1  1 -1  1  1 -1 >\n",
      "|v_458> = |-1 -1 -1  1  1 -1  1 -1  1 >\n",
      "|v_459> = |-1 -1 -1  1  1 -1  1 -1 -1 >\n",
      "|v_460> = |-1 -1 -1  1  1 -1 -1  1  1 >\n",
      "|v_461> = |-1 -1 -1  1  1 -1 -1  1 -1 >\n",
      "|v_462> = |-1 -1 -1  1  1 -1 -1 -1  1 >\n",
      "|v_463> = |-1 -1 -1  1  1 -1 -1 -1 -1 >\n",
      "|v_464> = |-1 -1 -1  1 -1  1  1  1  1 >\n",
      "|v_465> = |-1 -1 -1  1 -1  1  1  1 -1 >\n",
      "|v_466> = |-1 -1 -1  1 -1  1  1 -1  1 >\n",
      "|v_467> = |-1 -1 -1  1 -1  1  1 -1 -1 >\n",
      "|v_468> = |-1 -1 -1  1 -1  1 -1  1  1 >\n",
      "|v_469> = |-1 -1 -1  1 -1  1 -1  1 -1 >\n",
      "|v_470> = |-1 -1 -1  1 -1  1 -1 -1  1 >\n",
      "|v_471> = |-1 -1 -1  1 -1  1 -1 -1 -1 >\n",
      "|v_472> = |-1 -1 -1  1 -1 -1  1  1  1 >\n",
      "|v_473> = |-1 -1 -1  1 -1 -1  1  1 -1 >\n",
      "|v_474> = |-1 -1 -1  1 -1 -1  1 -1  1 >\n",
      "|v_475> = |-1 -1 -1  1 -1 -1  1 -1 -1 >\n",
      "|v_476> = |-1 -1 -1  1 -1 -1 -1  1  1 >\n",
      "|v_477> = |-1 -1 -1  1 -1 -1 -1  1 -1 >\n",
      "|v_478> = |-1 -1 -1  1 -1 -1 -1 -1  1 >\n",
      "|v_479> = |-1 -1 -1  1 -1 -1 -1 -1 -1 >\n",
      "|v_480> = |-1 -1 -1 -1  1  1  1  1  1 >\n",
      "|v_481> = |-1 -1 -1 -1  1  1  1  1 -1 >\n",
      "|v_482> = |-1 -1 -1 -1  1  1  1 -1  1 >\n",
      "|v_483> = |-1 -1 -1 -1  1  1  1 -1 -1 >\n",
      "|v_484> = |-1 -1 -1 -1  1  1 -1  1  1 >\n",
      "|v_485> = |-1 -1 -1 -1  1  1 -1  1 -1 >\n",
      "|v_486> = |-1 -1 -1 -1  1  1 -1 -1  1 >\n",
      "|v_487> = |-1 -1 -1 -1  1  1 -1 -1 -1 >\n",
      "|v_488> = |-1 -1 -1 -1  1 -1  1  1  1 >\n",
      "|v_489> = |-1 -1 -1 -1  1 -1  1  1 -1 >\n",
      "|v_490> = |-1 -1 -1 -1  1 -1  1 -1  1 >\n",
      "|v_491> = |-1 -1 -1 -1  1 -1  1 -1 -1 >\n",
      "|v_492> = |-1 -1 -1 -1  1 -1 -1  1  1 >\n",
      "|v_493> = |-1 -1 -1 -1  1 -1 -1  1 -1 >\n",
      "|v_494> = |-1 -1 -1 -1  1 -1 -1 -1  1 >\n",
      "|v_495> = |-1 -1 -1 -1  1 -1 -1 -1 -1 >\n",
      "|v_496> = |-1 -1 -1 -1 -1  1  1  1  1 >\n",
      "|v_497> = |-1 -1 -1 -1 -1  1  1  1 -1 >\n",
      "|v_498> = |-1 -1 -1 -1 -1  1  1 -1  1 >\n",
      "|v_499> = |-1 -1 -1 -1 -1  1  1 -1 -1 >\n",
      "|v_500> = |-1 -1 -1 -1 -1  1 -1  1  1 >\n",
      "|v_501> = |-1 -1 -1 -1 -1  1 -1  1 -1 >\n",
      "|v_502> = |-1 -1 -1 -1 -1  1 -1 -1  1 >\n",
      "|v_503> = |-1 -1 -1 -1 -1  1 -1 -1 -1 >\n",
      "|v_504> = |-1 -1 -1 -1 -1 -1  1  1  1 >\n",
      "|v_505> = |-1 -1 -1 -1 -1 -1  1  1 -1 >\n",
      "|v_506> = |-1 -1 -1 -1 -1 -1  1 -1  1 >\n",
      "|v_507> = |-1 -1 -1 -1 -1 -1  1 -1 -1 >\n",
      "|v_508> = |-1 -1 -1 -1 -1 -1 -1  1  1 >\n",
      "|v_509> = |-1 -1 -1 -1 -1 -1 -1  1 -1 >\n",
      "|v_510> = |-1 -1 -1 -1 -1 -1 -1 -1  1 >\n",
      "|v_511> = |-1 -1 -1 -1 -1 -1 -1 -1 -1 >\n"
     ]
    }
   ],
   "source": [
    "D = 2**L\n",
    "basis = pt.zeros((D,L)) \n",
    "for v_i in range(0,D):\n",
    "    fock_state = pt.zeros(D) \n",
    "    fock_state[v_i] = 1\n",
    "    for i in range(1,L+1):\n",
    "        basis[v_i,i-1] = pt.vdot(fock_state, Z[i]@fock_state)\n",
    "\n",
    "for v_i in range(0,D):  \n",
    "    string_fock_vector = \"|v_\" + \"{:03d}\".format(v_i) + \"> = |\"\n",
    "    for i in range(1,L+1):\n",
    "        tmp = int(basis[v_i,i-1].item())\n",
    "        if(tmp==1):\n",
    "            string_plus_minus = \" {:1d}\".format(tmp)\n",
    "        if(tmp==-1):\n",
    "            string_plus_minus = \"{:1d}\".format(tmp)\n",
    "            \n",
    "         \n",
    "        string_fock_vector = string_fock_vector + string_plus_minus + \" \"\n",
    "    string_fock_vector = string_fock_vector + \">\"\n",
    "    print(string_fock_vector)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e7173404",
   "metadata": {},
   "source": [
    "### The many-body wave-function $|\\psi\\rangle$ is naturally decomposed in the Fock basis:\n",
    "\n",
    "### \\begin{equation}\n",
    "|\\psi\\rangle = \\sum_{k} \\psi_k |v_l\\rangle,\n",
    "\\end{equation}\n",
    "### with $\\sum_k |\\psi_k|^2 = 1$.\n",
    "\n",
    "### As such, the state of the system is encoded in a vector:\n",
    "\n",
    "### \\begin{equation}\n",
    "   |\\psi\\rangle    \\Leftrightarrow \\begin{pmatrix} \\psi_1 \\\\ \n",
    "                               \\psi_2 \\\\\n",
    "                               \\vdots \\\\ \n",
    "                               \\psi_D \n",
    "                \\end{pmatrix}\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8e711830",
   "metadata": {},
   "source": [
    "## Cost and mixer Hamiltonians\n",
    "\n",
    "### Our cost Hamiltonian for Max-Cut problem has the following form:\n",
    "\n",
    "### \\begin{equation}\n",
    " \\hat{H}_C = \\frac{1}{2}\\sum_{\\langle i,j\\rangle} Z_i Z_j,\n",
    "\\end{equation}\n",
    "\n",
    "### where $\\langle i,j\\rangle$ means sum over nearest neighbors. Note, that each node $i$ has two coordinates, i.e. $ i \\equiv (i_x, i_y)$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6526a153",
   "metadata": {},
   "outputs": [],
   "source": [
    "H_C = pt.zeros((D,D)) + 0*1j\n",
    "for qubit_i_x in range(1,Lx+1):\n",
    "    for qubit_i_y in range(1,Ly+1):\n",
    "        qubit_i = get_linear_index(qubit_i_x,qubit_i_y)\n",
    "        for qubit_j_x in range(1,Lx+1):\n",
    "            for qubit_j_y in range(1,Ly+1):\n",
    "                qubit_j = get_linear_index(qubit_j_x,qubit_j_y)\n",
    "                # Sum over nearest neighbors only\n",
    "                if(np.abs(qubit_i_x - qubit_j_x)==1 and np.abs(qubit_i_y - qubit_j_y)==0):\n",
    "                    H_C = H_C + Z[qubit_i]@Z[qubit_j]\n",
    "                if(np.abs(qubit_i_x - qubit_j_x)==0 and np.abs(qubit_i_y - qubit_j_y)==1):\n",
    "                    H_C = H_C + Z[qubit_i]@Z[qubit_j]\n",
    "H_C = 0.5*H_C                    "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c9f1eaf9",
   "metadata": {},
   "source": [
    "### The mixer Hamiltonian has the following form:\n",
    "\n",
    "### \\begin{equation}\n",
    "\\hat{H}_B = \\sum_i X_i\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "55f87351",
   "metadata": {},
   "outputs": [],
   "source": [
    "H_B = pt.zeros((D,D)) + 0*1j\n",
    "for qubit_i_x in range(1,Lx+1):\n",
    "    for qubit_i_y in range(1,Ly+1):\n",
    "        qubit_i = get_linear_index(qubit_i_x,qubit_i_y)\n",
    "        H_B = H_B + X[qubit_i]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fcc56be3",
   "metadata": {},
   "source": [
    "### We define two unitary evolution operators for cost, and mixer Hamiltonians, parametrized by $\\gamma$ and $\\beta$\n",
    "\n",
    "### \\begin{equation}\n",
    "\\begin{split}\n",
    " U_B(\\beta) & = e^{-i\\beta \\hat{H}_B} \\\\\n",
    " U_C(\\gamma) & = e^{-i\\gamma \\hat{H}_C}.\n",
    "\\end{split}\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d753bfc3",
   "metadata": {},
   "outputs": [],
   "source": [
    "def U_B(H_B,beta):\n",
    "    return expm(-1j*beta*H_B)\n",
    "\n",
    "def U_C(H_C,gamma):\n",
    "    return expm(-1j*gamma*H_C) "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "93df42b5",
   "metadata": {},
   "source": [
    "### QAOA ansatz\n",
    "\n",
    "### Finally, we can define the QAOA ansatz, parametrized by two set of parameters, $\\vec{\\beta}$, and $\\vec{\\gamma}$, having $p$ elements:\n",
    "### \\begin{equation}\n",
    " |\\psi(\\vec{\\beta},\\vec{\\gamma})\\rangle = \\prod_{l = 1}^p U_B(\\beta_l)U_C(\\gamma_l)|s\\rangle,\n",
    "\\end{equation}\n",
    "### where state $|s\\rangle$ is coherent superposition of states $|0\\rangle$ adn $|1\\rangle$ of each qubit, i.e.\n",
    "### \\begin{equation}\n",
    " |s\\rangle = \\prod_{i = 1}^{L} H_i|0\\rangle,\n",
    "\\end{equation}\n",
    "### where $H_i$ is Hadamard matrix acting on $i$-th qubit initialized in $|0\\rangle$ state.\n",
    "\n",
    "### Length $p$ of vectors $\\vec{\\beta}$, and $\\vec{\\gamma}$ is one of the hyperparameters in $QAOA$ anstaz:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6a09d31a",
   "metadata": {},
   "outputs": [],
   "source": [
    "def psi(parameters):\n",
    "    beta_vec  = parameters[0,:]\n",
    "    gamma_vec = parameters[1,:]\n",
    "    psi_tmp = pt.zeros(D) + 0*1j\n",
    "    psi_tmp[0] = 1 # in basis of total Z operator Z = \\sum_i Z_i this elements corresponds  \n",
    "                   # to configurations with all spins pointing up\n",
    "   \n",
    "    # prepare Hadamard operation on each qubit\n",
    "    for qubit_i_x in range(1,Lx+1):\n",
    "        for qubit_i_y in range(1,Ly+1):\n",
    "            qubit_i = get_linear_index(qubit_i_x,qubit_i_y)\n",
    "            psi_tmp = Hadamard[qubit_i]@psi_tmp\n",
    "    \n",
    "    for i in range(0,p):\n",
    "        psi_tmp = U_C(H_C, gamma_vec[i])@psi_tmp\n",
    "        psi_tmp = U_B(H_B, beta_vec[i])@psi_tmp\n",
    "    \n",
    "    return psi_tmp"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4d321fd9",
   "metadata": {},
   "source": [
    "### The cost function is as a expectation value of the cost Hamilonian $F(\\vec{\\beta},\\vec{\\gamma}) = \\langle \\psi(\\vec{\\beta},\\vec{\\gamma})|\\hat{H}_C|\\psi(\\vec{\\beta},\\vec{\\gamma})\\rangle$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "08a9f51a",
   "metadata": {},
   "outputs": [],
   "source": [
    "def F(parameters):\n",
    "    return pt.vdot(psi(parameters), H_C@psi(parameters)).real"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "03e5989f",
   "metadata": {},
   "source": [
    "## Optimization loop\n",
    "\n",
    "### In the optimization loop, we choose a number of layers $p$, and initialize the $\\vec{\\beta}$, and $\\vec{\\gamma}$ as a random number. Next, with the help of automatic differentiation, we find a set of parameters that minimizes the expectation value of the cost Hamiltonian $\\hat{H}_C$. When we find the optimal set of parameters, our solution to the problem (i.e. spin configuration) corresponds to the Fock state with the maximal probability of the final wave-function, i.e. $v_{\\rm optimal} = argmax(\\rho)$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "47504e47",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Optimal spins configuration after single run corresponds to basis vector with index k = 170\n"
     ]
    }
   ],
   "source": [
    "p = 3\n",
    "parameters = pt.rand((2,p),requires_grad=True)   # first row corresponds to (beta_1,...,beta_p)\n",
    "                                                 # second row corresponds to (gamma_1,...,gamma_p)\n",
    "optimizer = pt.optim.Adam([parameters],lr = 1e-2)\n",
    "loss_vs_epoch = []\n",
    "for i in range(0,100):\n",
    "    loss = F(parameters)\n",
    "    optimizer.zero_grad()\n",
    "    loss.backward()\n",
    "    optimizer.step()\n",
    "    loss_vs_epoch.append(loss.item())\n",
    "\n",
    "psi_optimal = psi(parameters)\n",
    "rho = np.abs(psi_optimal.detach().numpy())**2\n",
    "i_optimal = np.argmax(rho)\n",
    "print('Optimal spins configuration after single run corresponds to basis vector with index k = ' + str(i_optimal))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "42c543ee",
   "metadata": {},
   "source": [
    "### We have to repeat an optimization procedure and make histogram of optimal configurations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "efcda2af",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0.5, 0, 'Fock state index')"
      ]
     },
     "execution_count": null,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAioAAAGwCAYAAACHJU4LAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAA9hAAAPYQGoP6dpAAArnElEQVR4nO3deXTU9b3/8dckJJOFJBAQk5AAEQQKQWRVCC4Ui+xFbUstnlLb01t2EKqAeg24ELytFqwCR22TWn4IvQWstiWVsorIZQmU9YQiAYIkpigmQZqQ5fP7w8OUIQnOhJn5fmd4Ps6ZczLfZeadj6N5+dnGYYwxAgAAsKEwqwsAAABoDEEFAADYFkEFAADYFkEFAADYFkEFAADYFkEFAADYFkEFAADYVjOrC7gedXV1Onv2rOLi4uRwOKwuBwAAeMAYo4qKCqWkpCgs7Np9JkEdVM6ePau0tDSrywAAAE1QVFSk1NTUa14T1EElLi5O0le/aHx8vMXVAAAAT5SXlystLc31d/xagjqoXB7uiY+PJ6gAABBkPJm2wWRaAABgWwQVAABgWwQVAABgWwQVAABgWwQVAABgWwQVAABgWwQVAABgWwQVAABgWwQVAABgW0G9My0QTGrrjHYVfq7Sikq1iYtS//REhYfxZZoAcC0EFSAA8g4Va8F7R1RcVuk6lpwQpazR3TQsI9nCygDA3hj6Afws71CxJq3IdwspklRSVqlJK/KVd6jYosoAwP4IKoAf1dYZLXjviEwD5y4fW/DeEdXWNXQFAICgAvjRrsLP6/WkXMlIKi6r1K7CzwNXFAAEEYIK4EelFY2HlKZcBwA3GoIK4Edt4qJ8eh0A3GgIKoAf9U9PVHJClBpbhOzQV6t/+qcnBrIsAAgaBBXAj8LDHMoa3a3Bc5fDS9bobuynAgCNIKgAfjYsI1nLHuktZzP3f92SEqK07JHe7KMCANfAhm9AAAzLSNa9nT/R3458Kkl6+6d3sjMtAHiAoAIEiMPxn1AyoGMrCysBgODB0A8AALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtggoAALAtS4NKTU2Nnn76aaWnpys6Olq33HKLnn32WdXV1VlZFgAAsIlmVr75iy++qOXLl+t3v/udunfvrj179ujRRx9VQkKCZsyYYWVpAADABiwNKh999JG+/e1va+TIkZKkDh066O2339aePXsavL6qqkpVVVWu5+Xl5QGpEwAAWMPSoZ9BgwZp48aNOnbsmCTpH//4h7Zv364RI0Y0eH12drYSEhJcj7S0tECWCwAAAszSHpU5c+aorKxMXbt2VXh4uGpra/XCCy/o4YcfbvD6efPmadasWa7n5eXlhBUAAEKYpUFl9erVWrFihVauXKnu3btr//79mjlzplJSUjRhwoR61zudTjmdTgsqBQAAVrA0qDz++OOaO3euvv/970uSevTooVOnTik7O7vBoAIAAG4sls5RuXjxosLC3EsIDw9neTIAAJBkcY/K6NGj9cILL6hdu3bq3r279u3bp5dfflk//vGPrSwLAADYhKVB5de//rX++7//W5MnT1ZpaalSUlL0s5/9TM8884yVZQEAAJuwNKjExcVp8eLFWrx4sZVlAAAAm+K7fgAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAAAgG0RVAALGGOsLgEAggJBBQAA2BZBBbAAHSoA4BmCCmABcgoAeIagAgAAbIugAliAybQA4BmCCmABYgoAeIagAgAAbIugAliAkR8A8AxBBbCAYfAHADxCUAEAALZFUAEswNAPAHjG8qDyySef6JFHHlGrVq0UExOj22+/XXv37rW6LAAAYAPNrHzz8+fPKzMzU4MHD9b69evVpk0bffzxx2rRooWVZQEAAJuwNKi8+OKLSktLU05OjutYhw4drCsICBCGfgDAM5YO/bz77rvq27evvvvd76pNmzbq1auX3njjjUavr6qqUnl5udsDCEas+gEAz1gaVE6cOKFly5bp1ltv1d/+9jdNnDhR06dP11tvvdXg9dnZ2UpISHA90tLSAlwxAAAIJIex8EtHIiMj1bdvX+3YscN1bPr06dq9e7c++uijetdXVVWpqqrK9by8vFxpaWkqKytTfHx8QGoGmmri7/cq73CJJOnwgvsV67R05BUALFNeXq6EhASP/n5b2qOSnJysbt26uR37xje+odOnTzd4vdPpVHx8vNsDCEYM/ACAZywNKpmZmSooKHA7duzYMbVv396iigAAgJ1YGlQee+wx7dy5UwsXLtTx48e1cuVKvf7665oyZYqVZQF+Z+GIKwAEFUuDSr9+/bRu3Tq9/fbbysjI0HPPPafFixdr/PjxVpYF+B0xBQA8Y/lsvlGjRmnUqFFWlwEAAGzI8i30gRsRIz8A4BmCCmAFggoAeISgAliAnWkBwDMEFQAAYFsEFcACzFEBAM8QVAALkFMAwDMEFQAAYFsEFcAC7EwLAJ4hqAAWIKYAgGcIKgAAwLYIKkCAXLl3CiM/AOAZggoQIFeGEzZ8AwDPEFQAAIBtEVSAADGNPgEANIagAgSI+9APAMATBBUAAGBbBBUgYFj1AwDeIqgAAcKqHwDwHkEFAADYltdBJS8vT9u3b3c9f+2113T77bfrBz/4gc6fP+/T4oBQcmUfCkM/AOAZr4PK448/rvLycknSwYMHNXv2bI0YMUInTpzQrFmzfF4gECqu/CJCcgoAeKaZtzcUFhaqW7dukqQ1a9Zo1KhRWrhwofLz8zVixAifFwgAAG5cXveoREZG6uLFi5Kkv//97xo6dKgkKTEx0dXTAqA+96Ef+lQAwBNe96gMGjRIs2bNUmZmpnbt2qXVq1dLko4dO6bU1FSfFwiECrdVP+QUAPCI1z0qr776qpo1a6Y//vGPWrZsmdq2bStJWr9+vYYNG+bzAgEAwI3L6x6Vdu3a6c9//nO947/61a98UhAQquhEAQDved2jEh4ertLS0nrHP/vsM4WHh/ukKCAUua36IbUAgEe8DiqNTQKsqqpSZGTkdRcEAABwmcdDP6+88ookyeFw6M0331Tz5s1d52pra7Vt2zZ17drV9xUCIYgt9AHAMx4HlctzUIwxWr58udswT2RkpDp06KDly5f7vkIgRLDqBwC853FQKSwslCQNHjxYa9euVcuWLf1WFAAAgNSEVT+bN2/2Rx1AyLtyuIcOFQDwjNdBpba2Vrm5udq4caNKS0tVV1fndn7Tpk0+Kw4IJe5DP0QVAPCE10FlxowZys3N1ciRI5WRkSGHw+GPugAAALwPKqtWrdIf/vAHvoAQ8JJbj4p1ZQBAUGnSlxJ26tTJH7UAIc1tjgpJBQA84nVQmT17tpYsWcIYOwAA8Duvh362b9+uzZs3a/369erevbsiIiLczq9du9ZnxQGhxD3bE/QBwBNeB5UWLVrogQce8EctQEi7MprQIQkAnvE6qOTk5PijDgAAgHq8nqMCoIlY9QMAXvO6RyU9Pf2ae6ecOHHiugoCQhWrfgDAe14HlZkzZ7o9r66u1r59+5SXl6fHH3/cV3UBAAA0bWfahrz22mvas2fPdRcEhCr3Dd/oUgEAT/hsjsrw4cO1Zs0aX70cEHJY9QMA3vNZUPnjH/+oxMREX70cAACA90M/vXr1cptMa4xRSUmJ/vWvf2np0qU+LQ4IJVfu5kyPCgB4xuugMnbsWLfnYWFhuummm3Tvvfeqa9euvqoLCDluQz/MUQEAj3gdVLKysvxRBwAAQD1eBxVJqq2t1TvvvKOjR4/K4XCoW7duGjNmjMLDw31dHxAy3Fb90KECAB7xOqgcP35cI0aM0CeffKIuXbrIGKNjx44pLS1Nf/nLX9SxY0d/1AkEPbIJAHjP61U/06dPV8eOHVVUVKT8/Hzt27dPp0+fVnp6uqZPn+6PGgEAwA3K6x6VrVu3aufOnW5LkVu1aqVFixYpMzPTp8UBIYVVPwDgNa97VJxOpyoqKuodv3DhgiIjI31SFBCKWPUDAN7zOqiMGjVK//Vf/6X/+7//kzFGxhjt3LlTEydO1JgxY/xRIwAAuEF5HVReeeUVdezYUQMGDFBUVJSioqKUmZmpTp06acmSJf6oEQgJrPoBAO95PUelRYsW+tOf/qTjx4/r6NGjMsaoW7du6tSpkz/qA0LGlcM95BQA8IxXQaW8vFzNmzdXWFiYOnXq5AondXV1Ki8vV3x8vF+KBAAANyaPh37WrVunvn37qrKyst65yspK9evXT++9955PiwNCifvQD30qAOAJj4PKsmXL9MQTTygmJqbeuZiYGM2ZM0evvvqqT4sDQolbULGuDAAIKh4HlUOHDunee+9t9Pzdd9+tgwcP+qImAAAASV4ElfPnz6umpqbR89XV1Tp//rxPigJCkds+KnSpAIBHPA4qHTp00J49exo9v2fPHrVv394nRQGhyH1eCkkFADzhcVB58MEH9dRTT+nTTz+td66kpERPP/20HnroIZ8WBwAAbmweB5W5c+cqLi5Ot956qyZPnqwlS5bolVde0aRJk9S5c2c1b95cc+fObXIh2dnZcjgcmjlzZpNfAwgWDP0AgGc83kclLi5OH374oebNm6fVq1e75qO0bNlSjzzyiBYuXKi4uLgmFbF79269/vrruu2225p0PxAMWPUDAN7zagv9hIQELV26VOfOndOnn36qkpISnTt3TkuXLlWLFi2aVMCFCxc0fvx4vfHGG2rZsmWTXgMAAIQmr7/rR5IcDoduuukmtWnTRg6H47oKmDJlikaOHKn77rvva6+tqqpSeXm52wMIFm5b6NOlAgAe8fq7fnxp1apVys/P1+7duz26Pjs7WwsWLPBzVYB/sDMtAHivST0qvlBUVKQZM2ZoxYoVioqK8uieefPmqayszPUoKiryc5WA77A4GQC8Z1mPyt69e1VaWqo+ffq4jtXW1mrbtm169dVXVVVVpfDwcLd7nE6nnE5noEsFAAAWsSyoDBkypN6W+48++qi6du2qOXPm1AspQLC7criHkR8A8EyTgsrUqVP17LPPKjExsclvHBcXp4yMDLdjsbGxatWqVb3jQChwH/ohqQCAJzyeo3LmzBnXzytXrtSFCxckST169GCuCAAA8AuPe1S6du2qVq1aKTMzU5WVlSoqKlK7du108uRJVVdX+6SYLVu2+OR1AFtiNi0AeM3jHpWysjL97//+r/r06aO6ujqNGDFCnTt3VlVVlf72t7+ppKTEn3UCQY+cAgDe8zioVFdXq3///po9e7aio6O1b98+5eTkKDw8XL/97W/VsWNHdenSxZ+1AgCAG4zHQz/x8fHq1auXMjMzdenSJV28eFGZmZlq1qyZVq9erdTUVO3atcuftQJBjVU/AOA9j3tUzp49q6efflpOp1M1NTXq27ev7rrrLl26dEn5+flyOBwaNGiQP2sFghqrfgDAex4HldatW2v06NHKzs5WTEyMdu/erWnTpsnhcOjnP/+54uPjdc899/izVgAAcINp8hb6CQkJ+t73vqeIiAht2rRJhYWFmjx5si9rA0KK+3f9WFcHAASTJm34duDAAbVt21aS1L59e0VERCgpKUnjxo3zaXFAKHH79mQL6wCAYNKkoJKWlub6+dChQz4rBgAA4EqWfXsycKNxH/qhTwUAPEFQAQLELahYVwYABBWCCgAAsC2CCmAFulQAwCMEFSBA3HamJakAgEcIKgAAwLYIKkCAuG2hT4cKAHiEoAIECDvTAoD3CCoAAMC2CCpAgLCFPgB4j6ACBAg70wKA9wgqAADAtggqQICYRn4GADSOoAIECKt+AMB7BBUAAGBbBBUgYBj8AQBvEVSAAGHoBwC8R1ABAAC2RVABAoSBHwDwHkEFCJArN3lj6AcAPENQAQAAtkVQAQLEfeiHLhUA8ARBBQgQVv0AgPcIKgAAwLYIKkCAuE2mtbAOAAgmBBUgQNzmqDD2AwAeIagAAADbIqgAgUInCgB4jaACBIj70I9lZQBAUCGoAAAA2yKoAAHivuqHLhUA8ARBBQgQhn4AwHsEFQAAYFsEFSBA2EIfALxHUAEC5Mp5KeQUAPAMQQUAANgWQQUIEPehH/pUAMATBBUgQEwjPwMAGkdQAQAAtkVQAQKFLhUA8BpBBQgQ91U/JBUA8ARBBQAA2BZBBQgQNnwDAO8RVIAAYYoKAHiPoAIAAGyLoAIEyJWbvDH0AwCeIagAAeI+9ENSAQBPEFQAAIBtEVSAAGHVDwB4j6ACWICcAgCeIagAAVDv25LpUgEAjxBUAACAbRFUgACo16FiTRkAEHQIKkAAXB1MGPkBAM8QVAAAgG0RVIAAuHoybb3JtQCABlkaVLKzs9WvXz/FxcWpTZs2Gjt2rAoKCqwsCfCLekM/llQBAMHH0qCydetWTZkyRTt37tSGDRtUU1OjoUOH6ssvv7SyLAAAYBPNrHzzvLw8t+c5OTlq06aN9u7dq7vvvrve9VVVVaqqqnI9Ly8v93uNgC+wjQoANI2t5qiUlZVJkhITExs8n52drYSEBNcjLS0tkOUBTXb1lxCSUwDAM7YJKsYYzZo1S4MGDVJGRkaD18ybN09lZWWuR1FRUYCrBAAAgWTp0M+Vpk6dqgMHDmj79u2NXuN0OuV0OgNYFeAb9Yd+6FMBAE/YIqhMmzZN7777rrZt26bU1FSrywEAADZhaVAxxmjatGlat26dtmzZovT0dCvLAQAANmNpUJkyZYpWrlypP/3pT4qLi1NJSYkkKSEhQdHR0VaWBvgUq34AoGksnUy7bNkylZWV6d5771VycrLrsXr1aivLAnyu/qofkgoAeMLyoR8AAIDG2GZ5MhDKGPoBgKYhqAABwHf9AEDTEFQAAIBtEVSAALh6PhZDPwDgGYIKEAD1h35IKgDgCYIKAACwLYIKEACs+gGApiGoAIFAMAGAJiGoAAAA2yKoAAFQbwt9xn4AwCMEFSAAmKMCAE1DUAEAALZFUAECgC30AaBpCCpAALAzLQA0DUEFAADYFkEFCAC20AeApiGoAAHAqh8AaBqCCgAAsC2CChAA9TZ8s6gOAAg2BBUgEOpNUiGqAIAnCCoAAMC2CCpAALDhGwA0DUEFCABW/QBA0xBUAACAbRFUgACov+qHLhUA8EQzqwsAbgQM/cCfauuMdhV+rtKKSrWJi1L/9ESFhzmsLgvwCYIKAASxvEPFWvDeERWXVbqOJSdEKWt0Nw3LSLawMsA3GPoBAoBVP/CHvEPFmrQi3y2kSFJJWaUmrchX3qFiiyoDfIegAgSAuWqsh6EfXK/aOqMF7x1pMPRePrbgvSOqrePDhuBGUAGAILSr8PN6PSlXMpKKyyq1q/DzwBUF+AFBBQiAepNpGfzBdSqtaDykNOU6wK4IKoAVyCm4Tm3ionx6HWBXBBUACEL90xOVnBClxhYhO/TV6p/+6YmBLAvwOYIKEAD1h36A6xMe5lDW6G4NnrscXrJGd2M/FQQ9ggoQAPV2pmXZD3xgWEaylj3SW1dnkaSEKC17pDf7qCAksOEbAASxYRnJahPnVEl5lSTp7Z/eyc60CCkEFSAA2EIf/vWfUDKgYysL6wB8j6EfIADYmRYAmoagAgAAbIugAgQAW+gDQNMQVIAAqD/0Q1IBAE8QVAAAgG0RVIAAYNUPADQNQQUICJIJ/MfBlikIYQQVAAgh7HqMUENQAQKg/tAPf0zgH7V1fLYQWggqQACw4RsCpbqWTxdCC0EFAELIpdo6q0sAfIqgAgQAq37gT3VXfKCqCSoIMQQVIACu3uCNDd/gSzW1BBWELoIKEAD0qMCfrhzuqa7hw4XQQlABgCB3ZS8Kc1QQaggqQADU61GxpgyEqGqGfhDCCCpAANSbo0JSgY/U1hm3vVMu1RBUEFoIKgAQxK7uQaFHBaGGoAIEQP0eFLpU4BtXz0lhjgpCDUEFsABDP/CV6pqre1T4cCG0EFQAIIhdHUyuDi5AsCOoAAHAPirwF+aoINQRVIAAYGda+AtzVBDqCCoAEMTq96gQghFaCCpAADD0A3+5est8hn4QappZXYAd1dYZ7fjnOf0xv0hnzv9bzmZhat3cKYfjP9cYY3TuwiVV1tQqqlm42/lrnePeG/Pe8xcvuX3Gth/7l8a//pFa2bhm7g2Oez/70v2z9f92ntTOj8+5HbNbzdwbHPdKksPhUNuW0RrYsbXuvKWVwsOuOBkgDmOs/X+7pUuX6he/+IWKi4vVvXt3LV68WHfddZdH95aXlyshIUFlZWWKj4/3ST15h4o16w//0MVLtT55PQAAQkGLmAgterCHhmUkX/drefP329Khn9WrV2vmzJl66qmntG/fPt11110aPny4Tp8+bUk9eYeKNXFFPiEFAICrfHGxWhNX5CvvUHFA39fSHpU77rhDvXv31rJly1zHvvGNb2js2LHKzs7+2vt92aNSW2c0MPvv+rTi0tdfDADADSo5IUrb53zzuoaBgqJH5dKlS9q7d6+GDh3qdnzo0KHasWNHg/dUVVWpvLzc7eEruwo/J6QAAPA1issqtavw84C9n2VB5dy5c6qtrdXNN9/sdvzmm29WSUlJg/dkZ2crISHB9UhLS/NZPaUVlT57LQAAQlkg/2ZavjzZ4XDvOjLG1Dt22bx581RWVuZ6FBUV+ayONnFRPnstAABCWSD/Zlq2PLl169YKDw+v13tSWlpar5flMqfTKafT6Zd6+qcn6ua4SIZ/AAC4huSEKPVPTwzY+1nWoxIZGak+ffpow4YNbsc3bNiggQMHBrye8DCHFnw7I+DvCwBAMMka3S2g+6lYOvQza9Ysvfnmm/rtb3+ro0eP6rHHHtPp06c1ceJES+oZlpGs5Y/0VkxkuCXvDwCAXbWMidDyR3r7ZB8Vb1i6M+24ceP02Wef6dlnn1VxcbEyMjL017/+Ve3bt7espmEZyfpWtyR2puVen9/rcDh0c4JTFf+uUWl5lf59qYadabnXZzuIXosda+Ze+98rsTPtdfPHzrQAAMC/gmIfFQAAgK9DUAEAALZFUAEAALZFUAEAALZFUAEAALZFUAEAALZFUAEAALZFUAEAALZFUAEAALZl6Rb61+vyprrl5eUWVwIAADx1+e+2J5vjB3VQqaiokCSlpaVZXAkAAPBWRUWFEhISrnlNUH/XT11dnc6ePau4uDg5PP12Lg+Vl5crLS1NRUVFfI+QH9C+/kX7+hft63+0sX9Z3b7GGFVUVCglJUVhYdeehRLUPSphYWFKTU3163vEx8fzL4kf0b7+Rfv6F+3rf7Sxf1nZvl/Xk3IZk2kBAIBtEVQAAIBtEVQa4XQ6lZWVJafTaXUpIYn29S/a179oX/+jjf0rmNo3qCfTAgCA0EaPCgAAsC2CCgAAsC2CCgAAsC2CCgAAsC2CSgOWLl2q9PR0RUVFqU+fPvrggw+sLikobNu2TaNHj1ZKSoocDofeeecdt/PGGM2fP18pKSmKjo7Wvffeq8OHD7tdU1VVpWnTpql169aKjY3VmDFjdObMmQD+FvaVnZ2tfv36KS4uTm3atNHYsWNVUFDgdg1t3HTLli3Tbbfd5toAa8CAAVq/fr3rPG3rW9nZ2XI4HJo5c6brGG18febPny+Hw+H2SEpKcp0P2vY1cLNq1SoTERFh3njjDXPkyBEzY8YMExsba06dOmV1abb317/+1Tz11FNmzZo1RpJZt26d2/lFixaZuLg4s2bNGnPw4EEzbtw4k5ycbMrLy13XTJw40bRt29Zs2LDB5Ofnm8GDB5uePXuampqaAP829nP//febnJwcc+jQIbN//34zcuRI065dO3PhwgXXNbRx07377rvmL3/5iykoKDAFBQXmySefNBEREebQoUPGGNrWl3bt2mU6dOhgbrvtNjNjxgzXcdr4+mRlZZnu3bub4uJi16O0tNR1Pljbl6Bylf79+5uJEye6HevatauZO3euRRUFp6uDSl1dnUlKSjKLFi1yHausrDQJCQlm+fLlxhhjvvjiCxMREWFWrVrluuaTTz4xYWFhJi8vL2C1B4vS0lIjyWzdutUYQxv7Q8uWLc2bb75J2/pQRUWFufXWW82GDRvMPffc4woqtPH1y8rKMj179mzwXDC3L0M/V7h06ZL27t2roUOHuh0fOnSoduzYYVFVoaGwsFAlJSVubet0OnXPPfe42nbv3r2qrq52uyYlJUUZGRm0fwPKysokSYmJiZJoY1+qra3VqlWr9OWXX2rAgAG0rQ9NmTJFI0eO1H333ed2nDb2jX/+859KSUlRenq6vv/97+vEiROSgrt9g/pLCX3t3Llzqq2t1c033+x2/Oabb1ZJSYlFVYWGy+3XUNueOnXKdU1kZKRatmxZ7xra350xRrNmzdKgQYOUkZEhiTb2hYMHD2rAgAGqrKxU8+bNtW7dOnXr1s31H2na9vqsWrVK+fn52r17d71zfH6v3x133KG33npLnTt31qeffqrnn39eAwcO1OHDh4O6fQkqDXA4HG7PjTH1jqFpmtK2tH99U6dO1YEDB7R9+/Z652jjpuvSpYv279+vL774QmvWrNGECRO0detW13natumKioo0Y8YMvf/++4qKimr0Otq46YYPH+76uUePHhowYIA6duyo3/3ud7rzzjslBWf7MvRzhdatWys8PLxeciwtLa2XQuGdyzPPr9W2SUlJunTpks6fP9/oNZCmTZumd999V5s3b1ZqaqrrOG18/SIjI9WpUyf17dtX2dnZ6tmzp5YsWULb+sDevXtVWlqqPn36qFmzZmrWrJm2bt2qV155Rc2aNXO1EW3sO7GxserRo4f++c9/BvVnmKByhcjISPXp00cbNmxwO75hwwYNHDjQoqpCQ3p6upKSktza9tKlS9q6daurbfv06aOIiAi3a4qLi3Xo0CHaX1/9X83UqVO1du1abdq0Senp6W7naWPfM8aoqqqKtvWBIUOG6ODBg9q/f7/r0bdvX40fP1779+/XLbfcQhv7WFVVlY4ePark5OTg/gxbMYPXzi4vT/7Nb35jjhw5YmbOnGliY2PNyZMnrS7N9ioqKsy+ffvMvn37jCTz8ssvm3379rmWdi9atMgkJCSYtWvXmoMHD5qHH364waVxqamp5u9//7vJz8833/zmNy1fGmcXkyZNMgkJCWbLli1uyw8vXrzouoY2brp58+aZbdu2mcLCQnPgwAHz5JNPmrCwMPP+++8bY2hbf7hy1Y8xtPH1mj17ttmyZYs5ceKE2blzpxk1apSJi4tz/f0K1vYlqDTgtddeM+3btzeRkZGmd+/eruWfuLbNmzcbSfUeEyZMMMZ8tTwuKyvLJCUlGafTae6++25z8OBBt9f497//baZOnWoSExNNdHS0GTVqlDl9+rQFv439NNS2kkxOTo7rGtq46X784x+7/r2/6aabzJAhQ1whxRja1h+uDiq08fW5vC9KRESESUlJMQ8++KA5fPiw63ywtq/DGGOs6csBAAC4NuaoAAAA2yKoAAAA2yKoAAAA2yKoAAAA2yKoAAAA2yKoAAAA2yKoAAAA2yKoAAAA2yKoAGjQyZMn5XA4tH//fqtL8RmHw6F33nnnul4jNzdXLVq08Ek9AL4eQQUIYj/60Y/kcDjqPY4fP251afXMnz9ft99+u9f3+TIYFBcXa/jw4T55LQCB0czqAgBcn2HDhiknJ8ft2E033WRRNfZ2+avuAQQPelSAIOd0OpWUlOT2CA8PlyRt3bpV/fv3l9PpVHJysubOnauamhrXvXV1dXrxxRfVqVMnOZ1OtWvXTi+88EKD71NXV6ef/vSn6ty5s06dOtXgNVu2bFH//v0VGxurFi1aKDMzU6dOnVJubq4WLFigf/zjH65en9zcXEnSyy+/rB49eig2NlZpaWmaPHmyLly44Hq9Rx99VGVlZa775s+fL+mrr6h/4okn1LZtW8XGxuqOO+7Qli1brtlWVw79XB7aWrt2rQYPHqyYmBj17NlTH330kds9ubm5ateunWJiYvTAAw/os88+q/e67733nvr06aOoqCjdcsstWrBggaudn332WaWkpLjdN2bMGN19992qq6u7Zr0AJL49GQhiEyZMMN/+9rcbPHfmzBkTExNjJk+ebI4ePWrWrVtnWrdubbKyslzXPPHEE6Zly5YmNzfXHD9+3HzwwQfmjTfeMMYYU1hYaCSZffv2maqqKvPQQw+Z22+/3Xz66acNvl91dbVJSEgwP//5z83x48fNkSNHTG5urjl16pS5ePGimT17tunevbspLi42xcXF5uLFi8YYY371q1+ZTZs2mRMnTpiNGzeaLl26mEmTJhljjKmqqjKLFy828fHxrvsqKiqMMcb84Ac/MAMHDjTbtm0zx48fN7/4xS+M0+k0x44da7S9JJl169a5/X5du3Y1f/7zn01BQYH5zne+Y9q3b2+qq6uNMcbs3LnTOBwOk52dbQoKCsySJUtMixYtTEJCgus18/LyTHx8vMnNzTUff/yxef/9902HDh3M/PnzjTHG1NTUmAEDBpixY8caY4xZtmyZSUhIMCdPnrzGP1kAlxFUgCA2YcIEEx4ebmJjY12P73znO8YYY5588knTpUsXU1dX57r+tddeM82bNze1tbWmvLzcOJ1OVzC52uU/5B988IG57777TGZmpvniiy8areWzzz4zksyWLVsaPJ+VlWV69uz5tb/TH/7wB9OqVSvX85ycHLdgYIwxx48fNw6Hw3zyySdux4cMGWLmzZvX6Gs3FFTefPNN1/nDhw8bSebo0aPGGGMefvhhM2zYMLfXGDdunFs9d911l1m4cKHbNb///e9NcnKy6/nHH39s4uLizJw5c0xMTIxZsWJF4w0AwA1zVIAgN3jwYC1btsz1PDY2VpJ09OhRDRgwQA6Hw3UuMzNTFy5c0JkzZ1RSUqKqqioNGTLkmq//8MMPKzU1VRs3blRMTEyj1yUmJupHP/qR7r//fn3rW9/Sfffdp+9973tKTk6+5utv3rxZCxcu1JEjR1ReXq6amhpVVlbqyy+/dP0uV8vPz5cxRp07d3Y7XlVVpVatWl3z/a522223uX6+XGtpaam6du2qo0eP6oEHHnC7fsCAAcrLy3M937t3r3bv3u02ZFZbW6vKykpdvHhRMTExuuWWW/TLX/5SP/vZzzRu3DiNHz/eqxqBGxlBBQhysbGx6tSpU73jxhi3kHL5mPTVXI3o6GiPXn/EiBFasWKFdu7cqW9+85vXvDYnJ0fTp09XXl6eVq9eraefflobNmzQnXfe2eD1p06d0ogRIzRx4kQ999xzSkxM1Pbt2/WTn/xE1dXVjb5PXV2dwsPDtXfvXtd8nMuaN2/u0e91WUREhOvny+11ee7I5fa6lrq6Oi1YsEAPPvhgvXNRUVGun7dt26bw8HCdPHlSNTU1ataM//wCnmAyLRCiunXrph07drj9sd2xY4fi4uLUtm1b3XrrrYqOjtbGjRuv+TqTJk3SokWLNGbMGG3duvVr37dXr16aN2+eduzYoYyMDK1cuVKSFBkZqdraWrdr9+zZo5qaGr300ku688471blzZ509e9btmobu69Wrl2pra1VaWqpOnTq5PXy5sqdbt27auXOn27Grn/fu3VsFBQX16ujUqZPCwr76T+zq1au1du1abdmyRUVFRXruued8ViMQ6oj0QIiaPHmyFi9erGnTpmnq1KkqKChQVlaWZs2apbCwMEVFRWnOnDl64oknFBkZqczMTP3rX//S4cOH9ZOf/MTttaZNm6ba2lqNGjVK69ev16BBg+q9X2FhoV5//XWNGTNGKSkpKigo0LFjx/TDH/5QktShQwcVFhZq//79Sk1NVVxcnDp27Kiamhr9+te/1ujRo/Xhhx9q+fLlbq/boUMHXbhwQRs3blTPnj0VExOjzp07a/z48frhD3+ol156Sb169dK5c+e0adMm9ejRQyNGjPBJG06fPl0DBw7U//zP/2js2LF6//333YZ9JOmZZ57RqFGjlJaWpu9+97sKCwvTgQMHdPDgQT3//PM6c+aMJk2apBdffFGDBg1Sbm6uRo4cqeHDhzfa0wTgCpbOkAFwXa616scYY7Zs2WL69etnIiMjTVJSkpkzZ45rRYsxxtTW1prnn3/etG/f3kRERJh27dq5JoZeuernspdeesnExcWZDz/8sN57lZSUmLFjx5rk5GQTGRlp2rdvb5555hlTW1trjDGmsrLSPPTQQ6ZFixZGksnJyTHGGPPyyy+b5ORkEx0dbe6//37z1ltvGUnm/PnzrteeOHGiadWqlZHkWrV06dIl88wzz5gOHTqYiIgIk5SUZB544AFz4MCBRttDDUymvfL3O3/+vJFkNm/e7Dr2m9/8xqSmppro6GgzevRo88tf/rLe5N68vDwzcOBAEx0dbeLj403//v3N66+/burq6syQIUPM/fff7zap+bHHHjMdO3Z0rWAC0DiHMR4MwgIAAFiAOSoAAMC2CCoAAMC2CCoAAMC2CCoAAMC2CCoAAMC2CCoAAMC2CCoAAMC2CCoAAMC2CCoAAMC2CCoAAMC2CCoAAMC2/j/AzdtyibKegwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "N_runs = 10\n",
    "histogram = np.zeros(D)\n",
    "for run in range(0,N_runs):\n",
    "    parameters = pt.rand((2,p),requires_grad=True) \n",
    "    optimizer = pt.optim.Adam([parameters],lr = 1e-2)\n",
    "    loss_vs_epoch = []\n",
    "    for i in range(0,200):\n",
    "        loss = F(parameters)\n",
    "        optimizer.zero_grad()\n",
    "        loss.backward()\n",
    "        optimizer.step()\n",
    "        loss_vs_epoch.append(loss.item())\n",
    "    psi_optimal = psi(parameters)\n",
    "    rho = np.abs(psi_optimal.detach().numpy())**2\n",
    "    histogram[np.argmax(rho)] += 1\n",
    "    \n",
    "plt.plot(histogram, '-o')\n",
    "plt.ylabel(\"# Counts\")\n",
    "plt.xlabel(\"Fock state index\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c9e0a474",
   "metadata": {},
   "source": [
    "### The optimal configuration corresponds to:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "afcaa397",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "|v_511> = | 1 -1  1 -1  1 -1  1 -1  1 >\n"
     ]
    }
   ],
   "source": [
    "k_optimal = np.argmax(histogram)\n",
    "string_fock_vector = \"|v_\" + \"{:03d}\".format(v_i) + \"> = |\"\n",
    "for i in range(1,L+1):\n",
    "    tmp = int(basis[k_optimal,i-1].item())\n",
    "    if(tmp==1):\n",
    "        string_plus_minus = \" {:1d}\".format(tmp)\n",
    "    if(tmp==-1):\n",
    "        string_plus_minus = \"{:1d}\".format(tmp)\n",
    "    string_fock_vector = string_fock_vector + string_plus_minus + \" \"\n",
    "string_fock_vector = string_fock_vector + \">\"\n",
    "print(string_fock_vector)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}