covering_and_packing.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "covering_and_packing.ipynb",
      "provenance": [],
      "collapsed_sections": [],
      "toc_visible": true,
      "authorship_tag": "ABX9TyMf9u1tj8QuNBuM0ZncvVYl",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/SatyaKuppam/QUBOs-and-Ising-models/blob/master/covering_and_packing.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "C8PbFoPTl334",
        "colab_type": "text"
      },
      "source": [
        "## Covering and packing problems:\n",
        "\n",
        "* Exact cover\n",
        "\n",
        "* Set Packing\n",
        "\n",
        "* Vertex Cover\n",
        "\n",
        "* Minimal Maximal Matching"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "rsFs55oAlrct",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "!pip install dwave-ocean-sdk"
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "cUbhQpiNmMLO",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "import networkx as nx\n",
        "import dwave_networkx as dnx\n",
        "import matplotlib\n",
        "import matplotlib.pyplot as plt\n",
        "import dimod\n",
        "%matplotlib inline\n",
        "\n",
        "from dwave.system.samplers import DWaveSampler\n",
        "from dwave.system.composites import EmbeddingComposite\n",
        "from dwave.cloud import Client\n",
        "from dwave.system.samplers import DWaveSampler\n",
        "from dwave.system.composites import EmbeddingComposite\n",
        "from dwave_qbsolv import QBSolv\n",
        "from collections import defaultdict\n",
        "\n",
        "import networkx as nx\n",
        "import dwave_networkx as dnx\n",
        "import neal"
      ],
      "execution_count": 5,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "ifVPxFGhmciL",
        "colab_type": "text"
      },
      "source": [
        "# 1. Exact Cover\n",
        "\n",
        "*Definition*: An instance $(U, V)$ of a set covering problem consists of a finite set $U$ and a family $H$ of subsets of $U$, such that every element of $U$ belongs to at least one subset in $V$. (From [1]).\n",
        "\n",
        "* The subet $C_i \\in V$ *covers* the elements of $U$. $$V = \\{C_1...C_N\\}$$ \n",
        "$$U = \\bigcup\\limits_{i=1}^{N} C_i$$  \n",
        "\n",
        "* The problem is to find the minimum size subset $C \\subseteq H$ whose members cover all elements of $U$. $$U = \\bigcup\\limits_{C_i \\in C} C_i$$ This is $NP$-hard.\n",
        "\n",
        "* The decision version of this problem is $NP$-complete. In the decision cersion we are asked if a covering exists with size atmost $k$."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "aVdrf0WZJiFi",
        "colab_type": "text"
      },
      "source": [
        "## 1.1 Decision version\n",
        "\n",
        "The energy functions for the decision version and the $NP$-hard version are derived from *Lucas*[2].\n",
        "\n",
        "First the decision version: $H = H_A  + H_B$. $H_A$ ensures that we have an exact cover and $H_B$ ensure that the cover is atmost $k$ which the user provides.\n",
        "\n",
        "$$H_A = \\sum_{\\alpha \\in U} A(1 - \\sum_{\\alpha \\in V_i} x_i)^2=\\sum_{\\alpha \\in U} A(1 + \\sum_{\\alpha \\in V_i}x_i \\sum_{\\alpha \\in V_j}x_j - 2\\sum_{\\alpha \\in V_i}x_i)$$\n",
        "\n",
        "$$H_B = B (\\sum_{i=1}^{|V|}x_i-k)$$\n",
        "\n",
        "So we atleast need $|V|$ spins to encode this problem."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "UED-hdA3mO-N",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "class exact_cover:\n",
        "  def __init__(self, U, V):\n",
        "    self.U = U\n",
        "    self.V = V\n",
        "    self.vartype = dimod.BINARY\n",
        "    self.linear = defaultdict(int)\n",
        "    self.quadratic = defaultdict(int)\n",
        "    self.offset = 0.0\n",
        "    \n",
        "    self.n = len(self.U)\n",
        "    self.B = 1.0\n",
        "    self.A = self.n * self.B + 1.0\n",
        "    print(f\"Setting A:{self.A} and B:{self.B}\")\n",
        "    self.N = len(self.V)\n",
        "\n",
        "    if self.n == 0 or self.N == 0:\n",
        "      raise ValueError(\"U or V cannot be empty.\")\n",
        "\n",
        "    self.solver = None\n",
        "\n",
        "  def _energy_function_h_a(self):\n",
        "    self.offset += self.A\n",
        "\n",
        "    for u in self.U:\n",
        "      for idx_i, v_i in enumerate(self.V):\n",
        "        if u in v_i:\n",
        "          for idx_j, v_j in enumerate(self.V):\n",
        "            if u in v_j:\n",
        "              if idx_i == idx_j:\n",
        "                self.linear[idx_i] += self.A\n",
        "              else:\n",
        "                self.quadratic[(idx_i, idx_j)] += self.A\n",
        "    \n",
        "    for u in self.U:\n",
        "      for idx, v_i in enumerate(self.V):\n",
        "        if u in v_i:\n",
        "          self.linear[idx] -= 2 * self.A\n",
        "\n",
        "  def _energy_function_h_b(self, k=0.0):\n",
        "    self.offset -= self.B * k\n",
        "\n",
        "    for u in U:\n",
        "      for idx, v_i in enumerate(V):\n",
        "        self.linear[idx] += self.B\n",
        "\n",
        "  def decide(self, k, solver=None, **kwargs):\n",
        "    self._energy_function_h_a()\n",
        "    self._energy_function_h_b(k)\n",
        "    return self._solve(solver=solver, **kwargs)\n",
        "\n",
        "  def solve(self, solver=None, **kwargs):\n",
        "    self._energy_function_h_a()\n",
        "    self._energy_function_h_b()\n",
        "    return self._solve(solver=solver, **kwargs)\n",
        "\n",
        "  def _solve(self, solver=None, **kwargs):\n",
        "    if solver == None:\n",
        "      print(\"Setting the solver as QBsolv\")\n",
        "      solver = QBSolv()\n",
        "\n",
        "    print(self.linear)\n",
        "    print(self.quadratic)\n",
        "    print(self.offset)\n",
        "    bqm = dimod.BinaryQuadraticModel(self.linear, self.quadratic, self.offset,\n",
        "                                     self.vartype)\n",
        "    return solver.sample(bqm, **kwargs)\n",
        "  "
      ],
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "9_drd8ftDcU3",
        "colab_type": "text"
      },
      "source": [
        "## 1.2 Testing $H_A$\n",
        "All covers should be equally likely."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "tvKqqBYxV7ZS",
        "colab_type": "code",
        "cellView": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 280
        },
        "outputId": "c92f0da7-6246-4685-cd7f-6b5fccbd25a0"
      },
      "source": [
        "U = [1,2,3,4,5,6]\n",
        "V = [[1], [2], [3],[4],[5],[6],[1,2],[3,4],[5,6]]\n",
        "\n",
        "exact_cover_decision_h_a = exact_cover(U, V)\n",
        "exact_cover_decision_h_a._energy_function_h_a()\n",
        "exact_cover_decision_h_a._solve(num_repeats=1000)"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Setting A:7.0 and B:1.0\n",
            "Setting the solver as QBsolv\n",
            "defaultdict(<class 'int'>, {0: -7.0, 6: -14.0, 1: -7.0, 2: -7.0, 7: -14.0, 3: -7.0, 4: -7.0, 8: -14.0, 5: -7.0})\n",
            "defaultdict(<class 'int'>, {(0, 6): 7.0, (6, 0): 7.0, (1, 6): 7.0, (6, 1): 7.0, (2, 7): 7.0, (7, 2): 7.0, (3, 7): 7.0, (7, 3): 7.0, (4, 8): 7.0, (8, 4): 7.0, (5, 8): 7.0, (8, 5): 7.0})\n",
            "7.0\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "SampleSet(rec.array([([0, 0, 1, 1, 1, 1, 1, 0, 0], -35., 132),\n",
              "           ([0, 0, 0, 0, 0, 0, 1, 1, 1], -35., 119),\n",
              "           ([0, 0, 1, 1, 0, 0, 1, 0, 1], -35., 124),\n",
              "           ([0, 0, 0, 0, 1, 1, 1, 1, 0], -35., 125),\n",
              "           ([1, 1, 0, 0, 0, 0, 0, 1, 1], -35., 129),\n",
              "           ([1, 1, 1, 1, 0, 0, 0, 0, 1], -35., 128),\n",
              "           ([1, 1, 0, 0, 1, 1, 0, 1, 0], -35., 110),\n",
              "           ([1, 1, 1, 1, 1, 1, 0, 0, 0], -35., 134)],\n",
              "          dtype=[('sample', 'i1', (9,)), ('energy', '<f8'), ('num_occurrences', '<i8')]), [0, 1, 2, 3, 4, 5, 6, 7, 8], {}, 'BINARY')"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 6
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "zWvnf2gEDhRU",
        "colab_type": "text"
      },
      "source": [
        "## 1.3 Testing $H_B$\n",
        "No subset should be selected because the minimum energy possible is $-Bk$."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "3s4oKD_T8wSc",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 159
        },
        "outputId": "f01fc35f-1ca0-44ee-c29b-6a48d81b8d42"
      },
      "source": [
        "exact_cover_decision_h_b = exact_cover(U, V)\n",
        "exact_cover_decision_h_b._energy_function_h_b(k=3)\n",
        "exact_cover_decision_h_b._solve(num_repeats=1000)"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Setting A:7.0 and B:1.0\n",
            "Setting the solver as QBsolv\n",
            "defaultdict(<class 'int'>, {0: 6.0, 1: 6.0, 2: 6.0, 3: 6.0, 4: 6.0, 5: 6.0, 6: 6.0, 7: 6.0, 8: 6.0})\n",
            "defaultdict(<class 'int'>, {})\n",
            "-3.0\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "SampleSet(rec.array([([0, 0, 0, 0, 0, 0, 0, 0, 0], -3., 1001)],\n",
              "          dtype=[('sample', 'i1', (9,)), ('energy', '<f8'), ('num_occurrences', '<i8')]), [0, 1, 2, 3, 4, 5, 6, 7, 8], {}, 'BINARY')"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 91
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "q0VCxhBdEhHR",
        "colab_type": "text"
      },
      "source": [
        "## 1.4 Experimentation:\n",
        " $C=[[1,2],[3,4],[5,6]]$ should be selected with high probability."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "xunJvpq7x9Qi",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 211
        },
        "outputId": "7c055620-3e25-4b48-9ad0-d76d7fe17a60"
      },
      "source": [
        "exact_cover_decision = exact_cover(U, V)\n",
        "exact_cover_decision.decide(k=3, num_repeats=1000)"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Setting A:7.0 and B:1.0\n",
            "Setting the solver as QBsolv\n",
            "defaultdict(<class 'int'>, {0: -1.0, 6: -8.0, 1: -1.0, 2: -1.0, 7: -8.0, 3: -1.0, 4: -1.0, 8: -8.0, 5: -1.0})\n",
            "defaultdict(<class 'int'>, {(0, 6): 7.0, (6, 0): 7.0, (1, 6): 7.0, (6, 1): 7.0, (2, 7): 7.0, (7, 2): 7.0, (3, 7): 7.0, (7, 3): 7.0, (4, 8): 7.0, (8, 4): 7.0, (5, 8): 7.0, (8, 5): 7.0})\n",
            "4.0\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "SampleSet(rec.array([([0, 0, 0, 0, 0, 0, 1, 1, 1], -20., 605),\n",
              "           ([0, 0, 1, 1, 0, 0, 1, 0, 1], -14., 160),\n",
              "           ([0, 0, 0, 0, 1, 1, 1, 1, 0], -14., 103),\n",
              "           ([1, 1, 0, 0, 0, 0, 0, 1, 1], -14., 136)],\n",
              "          dtype=[('sample', 'i1', (9,)), ('energy', '<f8'), ('num_occurrences', '<i8')]), [0, 1, 2, 3, 4, 5, 6, 7, 8], {}, 'BINARY')"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 7
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "LQaQ2iZdErpI",
        "colab_type": "text"
      },
      "source": [
        "No solution exists with $k=2$."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "QicNYGeiEKuy",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 211
        },
        "outputId": "94184ac6-d727-49e8-c073-5ad2f3efdd69"
      },
      "source": [
        "exact_cover_decision = exact_cover(U, V)\n",
        "exact_cover_decision.decide(k=2, num_repeats=1000)"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Setting A:7.0 and B:1.0\n",
            "Setting the solver as QBsolv\n",
            "defaultdict(<class 'int'>, {0: -1.0, 6: -8.0, 1: -1.0, 2: -1.0, 7: -8.0, 3: -1.0, 4: -1.0, 8: -8.0, 5: -1.0})\n",
            "defaultdict(<class 'int'>, {(0, 6): 7.0, (6, 0): 7.0, (1, 6): 7.0, (6, 1): 7.0, (2, 7): 7.0, (7, 2): 7.0, (3, 7): 7.0, (7, 3): 7.0, (4, 8): 7.0, (8, 4): 7.0, (5, 8): 7.0, (8, 5): 7.0})\n",
            "5.0\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "SampleSet(rec.array([([0, 0, 0, 0, 0, 0, 1, 1, 1], -19., 604),\n",
              "           ([1, 1, 0, 0, 0, 0, 0, 1, 1], -13., 157),\n",
              "           ([0, 0, 0, 0, 1, 1, 1, 1, 0], -13., 109),\n",
              "           ([0, 0, 1, 1, 0, 0, 1, 0, 1], -13., 134)],\n",
              "          dtype=[('sample', 'i1', (9,)), ('energy', '<f8'), ('num_occurrences', '<i8')]), [0, 1, 2, 3, 4, 5, 6, 7, 8], {}, 'BINARY')"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 8
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "R5Sn2AQwEyxs",
        "colab_type": "text"
      },
      "source": [
        "For $k=4$, we observe that its not the ground state, although solution is possible."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "HsDyp2hKE0sq",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 211
        },
        "outputId": "dff37a50-8363-4265-918b-097dafe001cc"
      },
      "source": [
        "exact_cover_decision = exact_cover(U, V)\n",
        "exact_cover_decision.decide(k=4, num_repeats=1000)"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Setting A:7.0 and B:1.0\n",
            "Setting the solver as QBsolv\n",
            "defaultdict(<class 'int'>, {0: -1.0, 6: -8.0, 1: -1.0, 2: -1.0, 7: -8.0, 3: -1.0, 4: -1.0, 8: -8.0, 5: -1.0})\n",
            "defaultdict(<class 'int'>, {(0, 6): 7.0, (6, 0): 7.0, (1, 6): 7.0, (6, 1): 7.0, (2, 7): 7.0, (7, 2): 7.0, (3, 7): 7.0, (7, 3): 7.0, (4, 8): 7.0, (8, 4): 7.0, (5, 8): 7.0, (8, 5): 7.0})\n",
            "3.0\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "SampleSet(rec.array([([0, 0, 0, 0, 0, 0, 1, 1, 1], -21., 591),\n",
              "           ([0, 0, 0, 0, 1, 1, 1, 1, 0], -15., 100),\n",
              "           ([0, 0, 1, 1, 0, 0, 1, 0, 1], -15., 171),\n",
              "           ([1, 1, 0, 0, 0, 0, 0, 1, 1], -15., 142)],\n",
              "          dtype=[('sample', 'i1', (9,)), ('energy', '<f8'), ('num_occurrences', '<i8')]), [0, 1, 2, 3, 4, 5, 6, 7, 8], {}, 'BINARY')"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 9
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Gjbsum3MFKsD",
        "colab_type": "text"
      },
      "source": [
        "We will find the smallest exact cover. The solution should choose $C=[[1,2],[3,4],[5,6]]$ with high probability."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "i7b9nw0T-z8L",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 211
        },
        "outputId": "099f1fa3-ecdc-4d3c-e276-be8ddafc5062"
      },
      "source": [
        "exact_cover_find = exact_cover(U, V)\n",
        "exact_cover_find.solve(num_repeats=1000)"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Setting A:7.0 and B:1.0\n",
            "Setting the solver as QBsolv\n",
            "defaultdict(<class 'int'>, {0: -1.0, 6: -8.0, 1: -1.0, 2: -1.0, 7: -8.0, 3: -1.0, 4: -1.0, 8: -8.0, 5: -1.0})\n",
            "defaultdict(<class 'int'>, {(0, 6): 7.0, (6, 0): 7.0, (1, 6): 7.0, (6, 1): 7.0, (2, 7): 7.0, (7, 2): 7.0, (3, 7): 7.0, (7, 3): 7.0, (4, 8): 7.0, (8, 4): 7.0, (5, 8): 7.0, (8, 5): 7.0})\n",
            "7.0\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "SampleSet(rec.array([([0, 0, 0, 0, 0, 0, 1, 1, 1], -17., 582),\n",
              "           ([1, 1, 0, 0, 0, 0, 0, 1, 1], -11., 128),\n",
              "           ([0, 0, 0, 0, 1, 1, 1, 1, 0], -11.,  96),\n",
              "           ([0, 0, 1, 1, 0, 0, 1, 0, 1], -11., 198)],\n",
              "          dtype=[('sample', 'i1', (9,)), ('energy', '<f8'), ('num_occurrences', '<i8')]), [0, 1, 2, 3, 4, 5, 6, 7, 8], {}, 'BINARY')"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 10
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "eezdSFa0Nbpd",
        "colab_type": "text"
      },
      "source": [
        "# 2. Set Packing"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "jc3ErNpQNLyH",
        "colab_type": "text"
      },
      "source": [
        "## 2.1 Problem Definition:\n",
        "\n",
        "We have a set of elements $U = \\{1..n\\}$ and a set $H=\\{V_1, V_2,..V_N\\}$ such that $V_i = \\{x|x\\in U\\}$. We need to find $C \\subseteq H$ which is the largest set of subsets $V_i$ such that they are all disjoint. \n",
        "\n",
        "NOTE: $C$ does not have to form a cover.\n",
        "\n",
        "\n",
        "\n",
        "### Equivalence to a graph problem\n",
        "This is equivalent to the maximal inpendent set (MIS) problem, if we have a graph $G=(V,E)$ such that, for every $V_i$ in $H$ we have a vertex in $G$ and an edge exists between the vertices in the graph iff $V_i\\cap V_j \\neq \\emptyset$. Then what is the maximal number of vertices that can be coloured such that no two coloured vertices are connected by an edge."
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "j9-d9LBFl5kl",
        "colab_type": "text"
      },
      "source": [
        "\n",
        "### 2.1.1 Energy Function:\n",
        "\n",
        "* $H = H_A + H_B$\n",
        "* $H_A$ penalizes $V_i$'s that are not disjoint.\n",
        "* $H_B$ ensures that $|C|$ is maximum.\n",
        "* $H_A = A\\sum\\limits_{i,j: V_i \\cap V_j \\neq \\emptyset} x_ix_j$\n",
        "* $H_B = -B\\sum\\limits_{i}x_i$\n",
        "* A > B\n",
        "\n",
        "where $x_i \\in \\{0,1\\}$"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "gxGzneHyAgH3",
        "colab_type": "text"
      },
      "source": [
        "## 2.2 Experimentation\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "9fSTsge6AvRS",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "from itertools import combinations\n",
        "\n",
        "class SetPacking:\n",
        "  def __init__(self, U, V):\n",
        "    self.U = U\n",
        "    self.V = V\n",
        "    self.vartype = dimod.BINARY\n",
        "    self.linear = defaultdict(int)\n",
        "    self.quadratic = defaultdict(int)\n",
        "    self.offset = 0.0\n",
        "    \n",
        "    self.n = len(self.U)\n",
        "    self.B = 1.0\n",
        "    self.A = 10 * self.B\n",
        "    print(f\"Setting A:{self.A} and B:{self.B}\")\n",
        "    self.N = len(self.V)\n",
        "\n",
        "    if self.n == 0 or self.N == 0:\n",
        "      raise ValueError(\"U or V cannot be empty.\")\n",
        "\n",
        "    self.solver = None\n",
        "\n",
        "  def _energy_function_h_a(self):\n",
        "    for idx_i, idx_j in list(combinations(range(len(self.V)), 2)):\n",
        "      if len(list(set(self.V[idx_i]) & set(self.V[idx_j]))) > 0:\n",
        "        self.quadratic[(idx_i, idx_j)] += self.A\n",
        "\n",
        "  def _energy_function_h_b(self):\n",
        "    for idx in range(len(self.V)):\n",
        "      self.linear[idx] -= self.B\n",
        "\n",
        "  def solve(self, solver=None, **kwargs):\n",
        "    self._energy_function_h_a()\n",
        "    self._energy_function_h_b()\n",
        "    return self._solve(solver=solver, **kwargs)\n",
        "\n",
        "  def _solve(self, solver=None, **kwargs):\n",
        "    if solver == None:\n",
        "      print(\"Setting the solver as QBsolv\")\n",
        "      solver = QBSolv()\n",
        "\n",
        "    print(self.linear)\n",
        "    print(self.quadratic)\n",
        "    print(self.offset)\n",
        "    bqm = dimod.BinaryQuadraticModel(self.linear, self.quadratic, self.offset,\n",
        "                                     self.vartype)\n",
        "    return solver.sample(bqm, **kwargs)\n",
        "  "
      ],
      "execution_count": 16,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "08JtNJ9QG38Q",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 139
        },
        "outputId": "369a2726-cd8d-4c85-e641-aa1eb2b278cb"
      },
      "source": [
        "U = [1,2,3,4,5,6]\n",
        "V = [[1,2],[3,4],[5,6], [1,2,3],[3,4,5], [4,5,6]]\n",
        "\n",
        "# test h_a energy function\n",
        "set_packing_h_a = SetPacking(U, V)\n",
        "set_packing_h_a._energy_function_h_a()\n",
        "set_packing_h_a._energy_function_h_b()\n",
        "set_packing_h_a._solve(num_repeats=1000)"
      ],
      "execution_count": 18,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Setting A:10.0 and B:1.0\n",
            "Setting the solver as QBsolv\n",
            "defaultdict(<class 'int'>, {0: -1.0, 1: -1.0, 2: -1.0, 3: -1.0, 4: -1.0, 5: -1.0})\n",
            "defaultdict(<class 'int'>, {(0, 3): 10.0, (1, 3): 10.0, (1, 4): 10.0, (1, 5): 10.0, (2, 4): 10.0, (2, 5): 10.0, (3, 4): 10.0, (4, 5): 10.0})\n",
            "0.0\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "SampleSet(rec.array([([1, 1, 1, 0, 0, 0], -3., 965), ([1, 0, 0, 0, 1, 0], -2.,  37)],\n",
              "          dtype=[('sample', 'i1', (6,)), ('energy', '<f8'), ('num_occurrences', '<i8')]), [0, 1, 2, 3, 4, 5], {}, 'BINARY')"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 18
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "eCqz7J29I8tQ",
        "colab_type": "text"
      },
      "source": [
        "# 3. Vertex Cover"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "b4MFUOvkI-56",
        "colab_type": "text"
      },
      "source": [
        "## 3.1 Problem Definition\n",
        "\n",
        "Given a directed graph $G = (V, E)$ what is the smallest number of vertices that can be coloured such that every edge is incident to a coloured vertex.\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "_eMfP21ZLd2w",
        "colab_type": "text"
      },
      "source": [
        "### 3.1.1 Energy Function.\n",
        "\n",
        "Let $x_v$ be a binary variable for each vertex, that is $0$ when the vertex is not cloured, but $1$ when its is.\n",
        "\n",
        "* $H = H_A + H_B$\n",
        "* $H_A = A \\sum\\limits_{x,v\\in E}(1-x_u)(1-x_v) = A(|E| + \\sum\\limits_{x,v \\in E}x_ux_v-x_u-x_v)$\n",
        "* $H_B = B \\sum\\limits_{v} x_v$"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "J72uLP4YL6Qz",
        "colab_type": "text"
      },
      "source": [
        "## 3.2 Experimentations"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Y_TddxEyMpeD",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "from itertools import combinations\n",
        "\n",
        "class VertexCover:\n",
        "  def __init__(self, G):\n",
        "    self.vartype = dimod.BINARY\n",
        "    self.linear = defaultdict(int)\n",
        "    self.quadratic = defaultdict(int)\n",
        "    self.offset = 0.0\n",
        "\n",
        "    self.B = 1.0\n",
        "    self.A = 10.0 * self.B\n",
        "    print(f\"Setting A:{self.A} and B:{self.B}\")\n",
        "    self.G = G\n",
        "    self.solver = None\n",
        "\n",
        "  def _energy_function_h_a(self):\n",
        "    for u, v in self.G.edges():\n",
        "      self.offset += self.A\n",
        "      self.linear[u] -= self.A\n",
        "      self.linear[v] -= self.A\n",
        "      self.quadratic[(u, v)] += self.A\n",
        "\n",
        "  def _energy_function_h_b(self):\n",
        "    for idx in range(len(self.G.nodes())):\n",
        "      self.linear[idx] += self.B\n",
        "\n",
        "  def solve(self, solver=None, **kwargs):\n",
        "    self._energy_function_h_a()\n",
        "    self._energy_function_h_b()\n",
        "    return self._solve(solver=solver, **kwargs)\n",
        "\n",
        "  def _solve(self, solver=None, **kwargs):\n",
        "    if solver == None:\n",
        "      print(\"Setting the solver as QBsolv\")\n",
        "      solver = QBSolv()\n",
        "\n",
        "    print(self.linear)\n",
        "    print(self.quadratic)\n",
        "    print(self.offset)\n",
        "    bqm = dimod.BinaryQuadraticModel(self.linear, self.quadratic, self.offset,\n",
        "                                     self.vartype)\n",
        "    return solver.sample(bqm, **kwargs)\n",
        "  "
      ],
      "execution_count": 58,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "PES3pbquXNQb",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 335
        },
        "outputId": "b660ef3f-c429-4f00-a09c-79a89b4b6d4b"
      },
      "source": [
        "## Define the graph\n",
        "G_three = nx.Graph()\n",
        "G_three.add_edges_from([(0,1),(0,2),(0,3),(0,4),(1,2),(3,4)])\n",
        "pos = nx.spring_layout(G_three)\n",
        "\n",
        "nx.draw_networkx_nodes(G_three, nodelist=G_three.nodes(), pos=pos, node_color='r')\n",
        "nx.draw_networkx_edges(G_three, pos)\n",
        "nx.draw_networkx_labels(G_three, pos, font_color='w', font_size=12)"
      ],
      "execution_count": 62,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "{0: Text(0.0005584830638797252, 0.001620121988841208, '0'),\n",
              " 1: Text(0.04787180948336928, -1.0, '1'),\n",
              " 2: Text(-0.6539035968808192, -0.7580909027598135, '2'),\n",
              " 3: Text(0.6512507923069597, 0.7580990584101648, '3'),\n",
              " 4: Text(-0.04577748797338965, 0.9983717223608077, '4')}"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 62
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAV0AAADnCAYAAAC9roUQAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deVhU9f4H8PewCAMCEpvisMgWmaIWXs0WKS1RoyJLs9W691Z0rZ9P0oKK2gBibl23eExLblRWptGiXrtZWXnN9aeSGjMMECIq4AIog2zn98f5oYwMq8Occ4b363l8gDlfznwwffvte77nc1SCIICIiKzDTuoCiIh6EoYuEZEVMXSJiKyIoUtEZEUMXSIiK3Jo66C3t7cQHBxspVKIiGzDgQMHygVB8DF3rM3QDQ4Oxv79+7unKqImpaVAZiZw5AhQUQF4eABRUcCzzwI+Zv/cEsmaSqX6s7VjbYYuUbfatw9ITwe2bRO/rqm5emzzZmDePGD8eCApCRg+XJoaiSyMa7okjYwMICYGyM4Ww7Z54AKA0Si+lp0tjsvIkKJKIovjTJesLyMDSEwEqqvbHysI4rjERPHrhITurY2om3GmS9a1b1/HA7e5puDlNQZSOIYuWVd6urh00JqwMPF4VlbLY0aj+P1ECsbQJespLRUvmrXVZGn1anE2bI4gAFu3AmVl3VMfkRUwdMl6MjPbPj5lCnDhArBjR+tjVKr2z0MkYwxdsp4jR1ruUmji5gZotcCrr7Z9DqMRyMmxfG1EVsLQJeupqGj9WEoK8P77wMmT7Z/n/HnL1URkZdwyRtbj4WH+9SFDgLFjgWHDOnSaz7//HivvvBPh4eGIiIi48jEsLAxqtdqCBRNZHkOXrCcqCti0qeUSQ0wMEBwMFBWJX/fuDdjbAwMHArfeajJUUKsxbuZM+I4ZA71eD51Oh99++w06nQ4FBQXw9vY2CeKmjwMGDECvXr2s8mMStUXV1uN6oqOjBfZeIIspLQWCglqGrloNuLtf/ToxUQzhhASgvNx0rLOzGM5mejI0NDTgxIkT0Ol0VwK56WNxcTH69+9vNpADAwNhb29v+Z+XeiyVSnVAEIRoc8c40yXr8fUVeylkZ5tuGzMaTffuXrwoBvO1gatSARMmtNoEx97eHsHBwQgODsZ9991ncqy2thaFhYVXgvjo0aPIzs6GTqdDWVkZgoODzQayv78/VCqVpX4HiDjTJSvbt09cTujsHWkA4OIC7NwJRJudQHSZ0WiEwWAwO0OuqqpCWFiY2UD29vZmIJNZbc10GbpkfZ3pvdDExQVYssTqvRcqKyuRl5dnNpAbGxtbBHF4eDjCw8PRp08fq9ZJ8sLQJflpCl6jse071FQqcc1XgsBtz9mzZ1sEsV6vh16vh1qtNhvIYWFhcHV1lbp06mYMXZKn/fvFXgpbt4rh2nxdV60Ww3jCBLGfroWXFLqTIAg4ffq02UDOz8/HDTfcYDaQQ0JC4OTkJHX5ZAEMXZK3sjLx1t6cHJTn5eGAwYBxiYnAtGk29+SIxsZGnDhxwmwgFxUVwd/f32wgBwUFwcGB172VgqFLiqHX6zF+/Hjk5eVJXYrV1dXVobCw0Gwgnz59GsHBwWYDuX///rCz482lcsLQJcWorq6Gl5cXqquruTOgmZqaGhgMBrOBXFFRgdDQULOB7Ovry99HCTB0SVG8vLyQm5sLb29vqUtRhKqqKuTl5ZkN5Lq6uis7Kq4NZE9PT6lLt1m8OYIURaPRoLi4mKHbQW5ubhg2bBiGmeldce7cuSs7KnQ6HbZs2YJ33nkHer0eTk5OZgM5LCwMvXv3luAn6RkYuiQ7TaE7dOhQqUtRvBtuuAEjRozAiBEjTF4XBAFnzpwxCeTPPvsMOp0OBoMBnp6eZgM5JCQEzs7OEv00toGhS7LTFLrUfVQqFfr27Yu+ffvizjvvNDnW2NiI4uJik0D++eefodPp8Oeff6Jfv35mAzk4OFj5OyxKS8WdNEeOiK1IPTzERk3PPmuxnTQK/x0iW8TQlZadnR0CAwMRGBiIMWPGmByrr6+/ssOiKZC3bt0KnU6HU6dOISgoyGwgazQaee+w2LdP3DO+bZv4dfOmTJs3A/PmiX1DkpKA4cOv6614IY1kZ/369di5cycy+VgeRampqUF+fr5JIDd9PH/+/JUdFtcGsp+fn7Q7LLrh7kheSCNF4UxXmZydnTFw4EAMHDiwxbGLFy9e2WGh1+vxyy+/4IMPPoBOp8Ply5evhPG1gXzDDTd0b9Gd6QMiCOK4xETx6y7els6ZLsnO8ePHER8fjz/++EPqUsgKzp8/fyWMr50hOzg4mDQSav65m5vb9b1xax3vsrKAMWMAV1fg9Glg0SLxUVLNtdPxjvt0SVGqqqrQr18/VFVVcWN/DyYIAsrKykwaCTV9npeXB3d3d7OBHBoa2rHHNj38cMvezoD4xJK8PKC2FrjxRuCnn4CJE4GDB6+OUamA+HjxSShmcHmBFMXNzQ329vaoqKhgi8QeTKVSwdfXF76+vrjjjjtMjjU2NqKkpMQkkHft2gW9Xo+CggL4+fmZDeQBAwbA0dFR3KWwbZv5Ndxjx65+Lgjir9BQ09AVBLFRU1lZp3c1MHRJlprWdRm6ZI6dnR00Gg00Gg3uuecek2P19fUoKioyCeTt27dDr9fj5MmTCAgIQJK9PZ6sq0OrT81bvVpsuOTiIobt1q0tx6hU4vay117rVO0MXZIljUaDEydOYNCgQVKXQgrj4OCAkJAQhISEIDY21uTY5cuXUVBQANcXX0Sv3NzWT/KPfwAvvwzcdpu47nv5cssxRiOQk9Pp+mS8cY56Mu5goO7g5OSEyMhIBHTkIlxjI7BrF6DRtL5T4fz5TtfA0CVZCggIYOhS9/Hw6PhYBwdxTdecLjQNYuiSLHGmS90qKgqCuR4SPj7AlCnidjE7O+C++4CpU4EdO1qOVauBwYM7/dYMXZIlhi51l6qqKqysqsLl5rf6NhEEcSmhuFhcOliyBJgxA/jmG/Njp03r9PszdEmWGLpkaZWVlViwYAFCQ0Px37w81I4ZI+5AaK68XLxw5ul5tdnNunUtT6ZSic/v60ITHIYuyRJDlyzlwoUL0Gq1CA0NxfHjx/Hzzz9jw4YNcE9PF5cIukKtFpvfdAFDl2TJw8MDDQ0NqKyslLoUUqhz585h3rx5CAsLQ35+Pnbt2oWsrCxERkaKA4YPF5cPXFw6d2IXF/H7uviEaoYuyZJKpYJGo8HJkyelLoUUpry8HLNnz0Z4eDhOnjyJPXv2IDMzExERES0HJyRcDd72bjlXqa4Gbheb3QAMXZIxLjFQZ5SWluKNN95AREQEysvLsX//fqxbtw6hrW33apKQIDaviY8HnJ1bLjmo1eLr8fHiuOsIXIB3pJGMMXSpI06fPo0lS5bggw8+wNSpU3Ho0CEEBgZ27iTR0WLzmrIy8dbenBxx94Knp7gtbNo0PjmCbB9Dl9pSUlKCRYsW4cMPP8STTz6JnJwc9O/f//pO6uPT6V4KncXlBZIthi6ZU1xcjJdffhmDBg2CnZ0djh49ihUrVlx/4FoJQ5dki6FLzRUVFSEhIQFRUVFwdnbG8ePHsWzZMvTr10/q0jqFoUuyxdAlACgoKMDzzz+PYcOGoU+fPsjNzcXixYvh5+cndWldwtAl2WLo9mx5eXl47rnnEB0dDV9fX+h0OqSnp8PHQhe0pMLQJdny8vJCdXU1qjvy0ECyGTqdDs888wxGjhyJwMBA5OXlITU1FV5eXlKXZhEMXZItlUqF/v378waJHuL48eN44okncPvttyMsLAx5eXmYP38+PLvQPlHOGLoka1xisH2///47HnvsMYwePRqDBg2CwWBAcnKyzT6qiaFLstb02B6yPYcPH8YjjzyCsWPH4pZbbkF+fj6SkpLg7u4udWndiqFLssaZru05ePAg4uPjERsbi1GjRsFgMOD1119H7969pS7NKhi6JGt8bI/t2LdvH+Li4hAXF4eYmBjk5+fj1Vdfhaurq9SlWRVDl2SNM13l2717NyZMmICHH34YsbGxMBgM+J//+R+ou9rLVuHYe4FkjaGrXL/++iu0Wi1yc3ORlJSEL7/8Ek5OTlKXJTmGLskaQ1d5du7cibfeeguFhYWYNWsWnn76afTq1UvqsmSDoUuy5uPjg4qKCtTU1MDZ3NNbSRYEQcAPP/wArVaLkydPYvbs2XjyySfh6OgodWmyw9AlWbOzs4O/vz9KSkoQEhIidTl0DUEQ8J///AdarRZlZWWYM2cOpk6dCgcHRktr+DtDste0xMDQlQ9BELBt2zZotVpUVlYiOTkZkydPhr29vdSlyR5Dl2SP67ryIQgCvv32W2i1WtTU1CA5ORmTJk1i2HYCQ5dkj6ErvcbGRnz11VdISUlBQ0MD5s6di/j4eNjZcddpZzF0SfY0Gg0MBoPUZfRIjY2N2Lx5M1JSUuDg4IB58+YhLi6OYXsdGLokexqNBjt37pS6jB6loaEBGzduREpKClxdXZGWloaJEydC1d5jyqldDF2SPS4vWE99fT0+++wzpKamok+fPli6dCnGjRvHsLUghi7JHkO3+9XX1+Pjjz9GWloa/Pz8sHLlSowZM4Zh2w0YuiR7ffv2RXl5Oerq6rjZ3sLq6uqQlZWFtLQ0BAYGYs2aNYiJiWHYdiOGLsmevb09/Pz8cOrUKQQGBkpdjk2ora1FZmYm0tPTERoaivXr1+Ouu+6SuqwegaFLitC0xMDQvT6XL1/GBx98gIULF+Kmm27CRx99hNtvv13qsnoUhi4pAtd1r4/RaMS6devw9ttvY8iQIfj8888xYsQIqcvqkRi6pAh8bE/XVFdXY82aNVi8eDGGDx+O7OxsREdHS11Wj8bQJUXgTLdzLl26hIyMDCxduhSjRo3Cli1bMGzYMKnLIvDJEaQQfGxPx1RVVWHhwoUICQnB3r17sX37dmzatImBKyOc6ZIicKbbtoqKCqxatQrLly/H2LFj8cMPP+Dmm2+Wuiwyg6FLisDQNe/ChQtYsWIFVq5cidjYWPz888+IjIyUuixqA5cXSBH69euHM2fOoL6+XupSZOHcuXOYO3cuwsLCUFBQgP/+97/Iyspi4CoAQ5cUwdHREd7e3jhz5ozUpUiqvLwcs2fPRnh4OEpKSrBnzx6sX78e4eHhUpdGHcTQJcXoyUsMpaWleOONNxAREYHy8nIcOHAA69atQ2hoqNSlUScxdEkxemLonj59GjNnzkRkZCQuXryIQ4cOYc2aNQgODpa6NOoihi4pRk8K3ZKSEsyYMQMDBw5EfX09cnJysHr1at4GbQMYuqQYPSF0i4uLMX36dAwaNAh2dnY4evQoli9fjv79+0tdGlkIQ5cUw5ZD988//0RCQgKioqLg4uKC48ePY9myZejXr5/UpZGFMXRJMWwxdAsKCvD3v/8dt9xyC/r06YPc3FwsWrQIfn5+UpdG3YShS4phS6Gbl5eH5557DtHR0ejbty90Oh3S09Ph4+MjdWnUzRi6pBj+/v4oKSlBY2Oj1KV0WW5uLp5++mmMHDkSgYGByMvLQ0pKCry8vKQujayEoUuK4ezsDA8PD5SVlUldSqcdO3YMTzzxBO644w5ERETAYDBg/vz58PT0lLo0sjKGLimK0pYYcnJyMGXKFMTExGDQoEEwGAyYM2cOPDw8pC6NJMLQJUVRSugePnwYjzzyCO69915ER0cjPz8fSUlJcHd3l7o0khhDlxRF7k+QOHjwIB566CHExsZi1KhRMBgMeO2119C7d2+pSyOZYOiSosh1prt3717ExcUhLi4O99xzD/Lz8/Hqq6/C1dVV6tJIZhi6pChyC93du3dj/PjxmDRpEmJjY2EwGPDKK69ArVZLXRrJFJuYk6LI5bE9v/76K9566y3odDrMmjUL2dnZcHJykrosUgCGLimK1DPdn376CVqtFoWFhZg1axaefvpp9OrVS7J6SHkYuqQo/fv3x8mTJyEIAlQqlVXeUxAE/PDDD9BqtSgpKcHs2bPxxBNPwNHR0SrvT7aFoUuK4uLiAhcXF5w9exbe3t7d+l6CIOC7776DVqvF2bNnMWfOHDz22GNwcOBfG+o6/ukhxWlaYuiu0BUEAdu2bYNWq0VlZSWSk5MxefJk2Nvbd8v7Uc/C0CXFaQrdoUOHWvS8giDgm2++gVarxeXLlzF37lxMmjQJdnbc5EOWw9AlxbH0xbTGxkZ89dVX0Gq1EAQBc+fOxUMPPcSwpW7B0CXFsVToNjY2YtOmTUhJSYGjoyPmz5+PBx54wGoX6KhnYuiSspSWIjYnB3UHDgCHDwMeHkBUFPDss0AHe9E2NDRg48aNSElJgaurK9LT0zFhwgSGLVmFShCEVg9GR0cL+/fvt2I5RK3Ytw9ITwe2bUNDYyPsa2uvHlOrAUEAxo8HkpKA4cPNnqK+vh6ffvopUlNT4enpiXnz5mHcuHEMW7I4lUp1QBCEaHPHONMl+cvIABITAaMREAS02ENgNIofs7OB7duBJUuAhIQrh+vr6/Hxxx8jLS0Nfn5+WLVqFcaMGcOwJUkwdEnemgK3urr9sYIgjktMBADU/e1v+PDDD7FgwQIEBgbivffew+jRoxm2JCleniX52rfPfOB6egKbNwMXLwKFhcDUqabHq6tRN2MGJgUF4dNPP8X69evx448/IiYmhoFLkuNMl+QrPf3q0kFzq1cDtbWAnx8wdCiwZYt4Ue3YsStD7GtrsT4iAl7/+Y8VCyZqH2e6JE+lpcC2beKSQXMuLsCkSUByMnDpErBrF/D118BTT5kMswPgtWcPoMDnqZFtY+iSPGVmmn89IgKorwf0+quvHT4M3Hxzy7EqVevnIZIIQ5fk6cgRoKam5eu9ewOVlaavVVQAbm4txxqNQE5O99RH1EUMXZKnigrzr1+8CFz7cEd3d6Cqyvz48+ctWxfRdWLokjy19ohynQ5wcADCwq6+NmQIcPSo+fGenpavjeg6MHRJnqKiIDg7t3y9ulrcLqbVihfVRo0CHnwQyMpqOVatBgYP7v5aiTqBoUuyc/bsWaSfOoXL5tZ0AeCll8RALS0FNmwQ7z5rtl3sCkEApk3r1lqJOouhS7JRVlaGN998ExERESi4dAkN990n7kC41vnzQHy8eFEtKEgM3mupVMCECR1ugkNkLQxdktzp06eRmJiIG2+8EZWVlTh48CDee+89uKamijParlCrxeY3RDLD0CXJlJSUYMaMGRg4cCBqa2tx5MgRvPvuuwgKChIHDB8uNq9xcenciV1cxO+LNtvkiUhSDF2yuhMnTmD69OkYNGgQ7OzscPToUaxYsQIajabl4ISEq8HbXt8Elepq4DbrMkYkJwxdsprCwkK8+OKLGDJkCFxcXHD8+HEsW7YM/fr1a/sbExKAnTvFdVxn55ZLDmq1+Hp8vDiOgUsyxoY31O3y8/OxYMECfPnll3jhhReQm5sLn85e4IqOBjZtEnspZGYCOTn48/BhnLx0CaNeeEHcpcCLZqQADF3qNnq9HgsWLMA333yDl156CTqdDl5eXtd3Uh8f4LXXAAB/bN+OpUuX4rv//5pICRi6ZHF//PEH0tLSsG3bNrz88svQ6/Xw7IY7wzQaDU6cOGHx8xJ1J67pksUcPXoUU6dOxV133YXIyEgYDAbMmzevWwIXAAICAnDixAm09Zw/Irlh6NJ1O3LkCB599FHcc889GDp0KAwGA2bPng2P1vonWIi7uzvs7OxQ0VpzHCIZYuhSlx08eBDx8fEYN24cRo4cifz8fLzxxhtwM9dmsZsEBASguLjYau9HdL0YutRpe/fuRVxcHOLi4hATEwODwYCZM2fC1dXV6rVwXZeUhhfSqMN2794NrVaL33//HW+++SY2btwIZ3OdwKyIM11SGoYuteuXX36BVquFXq9HUlISsrOz4eTkJHVZADjTJeXh8gKZJQgCfvzxR9x999145pln8Nhjj0Gn0+GFF16QTeACnOmS8nCmSyYEQcCOHTug1Wpx6tQpzJkzB48//jgcHR2lLs0sznRJaRi6BEAM2+3bt0Or1eLcuXNITk7GlClT4OAg7z8inOmS0sj7bxR1O0EQsGXLFmi1Wly6dAnJycl49NFHYW9vL3VpHdI00xUEAar2upARyQBDt4dqbGzE119/Da1Wi4aGBiQnJ+Phhx+GnZ2ylvnd3d1hb2+PiooK9OnTR+pyiNrF0O1hGhsbsXnzZqSkpMDe3h5z587FAw88oLiwba5ptsvQJSVg6PYQDQ0N2LhxI1JTU+Hi4oK0tDRMnDjRJv6XvGlddzCf/EsKwNC1cfX19fj000+RmpoKT09PLFmyBOPGjbOJsG3CHQykJAxdG1VXV4ePP/4YaWlp6Nu3L1atWoUxY8bYVNg2aeo2RqQEDF0bU1tbi6ysLCxYsABBQUFYu3YtRo8ebZNh20Sj0eDXX3+VugyiDmHo2ojLly8jMzMT6enpiIiIQGZmJu68806py7IKznRJSRi6CldTU4P3338fCxcuxODBg7FhwwbcdtttUpdlVRqNhjdIkGIwdBWquroaa9euxaJFi3Drrbdi8+bNGD58uNRlSaL5EyRseRmFbINyN2f2UJcuXcLSpUsRGhqKnTt34ptvvsHXX3/dYwMXANzc3ODg4IALFy5IXQpRuxi6ClFVVYW3334bISEh+O233/Dvf/8bmzdvxi233CJ1abLAdV1SCoauzFVUVCAtLQ2hoaE4dOgQduzYgY0bN2LIkCFSlyYrXNclpeCarkxduHABy5cvx8qVKzF+/Hj8/PPPiIyMlLos2eJMl5SCM12ZOXfuHObOnYuwsDAUFhZi9+7dyMrKYuC2gzNdUgqGrkyUl5dj1qxZCA8Px6lTp7B3716sX78e4eHhUpemCJzpklIwdCV25swZvP7664iIiMD58+dx8OBBrF27FiEhIVKXpiic6ZJSMHQlcurUKbz66qu46aabYDQacfjwYWRkZCAoKEjq0hSJM11SCoaulRUXF+OVV17BzTffDEEQ8Pvvv2PlypUICAiQujRFa5rpCoIgdSlEbWLoWklRURFeeuklREVFwcnJCceOHcM777wDf39/qUuzCW5ubnB0dMT58+elLoWoTQzdblZQUIDnn38ew4YNg4eHB3Jzc7F48WL07dtX6tJsDtd1SQkYut0kLy8Pzz33HKKjo+Hn5wedTof09HT4+PhIXZrN4rouKQFvjrCw3NxcpKWlYevWrZg+fTry8vLg6ekpdVk9Ame6pASc6VrIsWPH8Pjjj+OOO+5AREQEDAYD5s+fz8C1Is50SQkYutfpyJEjmDx5Mu6++25ERUUhPz8fc+bMgYeHh9Sl9Tic6ZISMHS76NChQ3j44Ydx33334S9/+QsMBgPefPNNuLm5SV1aj8WZLikBQ7eT9u/fjwcffBATJkzAXXfdhfz8fCQmJqJ3795Sl9bjcaZLSsDQ7aDffvsNEydOxEMPPYR7770XBoMBM2bMgIuLi9Sl0f9rehQ7b5AgOWPotmPXrl0YN24cpkyZgri4OBgMBkyfPh1qtVrq0ugabm5u6NWrF2+QIFlj6LZi586dGDNmDJ566ik8+uij0Ov1ePHFF+Hk5CR1adQGruuS3DF0mxEEATt27MDo0aPx17/+FU8++SRyc3Pxt7/9Db169ZK6POoAruuS3PHmCIhh+91330Gr1aK8vBxz5szB1KlT4eDA3x6l4UyX5K5Hp4ogCNi6dSu0Wi2qqqqQnJyMyZMnw97eXurSqIs40yW565GhKwgCvv76a2i1WtTV1SE5ORmTJk2CnR1XW5QuICAAP/30k9RlELWqR4VuY2MjvvzyS6SkpEClUmHu3Ll48MEHGbY2hDNdkrseEboNDQ344osvkJKSArVajZSUFNx///1QqVRSl0YWxjVdkjubDt36+np89tlnSE1NhYeHBxYvXozY2FiGrQ1r/gQJ/ncmObLJ0K2vr8cnn3yC1NRU+Pr6YsWKFRg7diz/EvYAvXv3hpOTE86dOwcvLy+pyyFqwaZCt66uDllZWViwYAE0Gg3WrFmDmJgYhm0P0zTbZeiSHNnEFaTa2lq89957iIiIwCeffIL3338fP/30E+6++24Gbg/EdV2SM0XPdGtqavDBBx9g4cKFGDhwID766CPcfvvtUpdFEuMOBpIzRYau0WjE2rVrsWjRIgwdOhQbN27EiBEjpC6LZIIzXZIzRS0vXLp0CcuWLUNoaCh++OEHfPXVV/j2228ZuGSiqcUjkRwpYqZ78eJFvPvuu1i2bBnuuOMObN26FUOHDpW6LJKpgIAALi+QbMk6dCsrK7Fq1Sr885//xD333IPvv/8egwYNkroskjnOdEnOLBu6paVAZiZw5AhQUQF4eABRUcCzzwI+Ph0+zYULF7BixQqsXLkS48aNw86dO3HTTTdZtFSyXbxBguTMMqG7bx+Qng5s2yZ+XVNz9djmzcC8ecD48UBSEjB8eKunOXfuHJYvX47Vq1fj/vvvx65duxAREWGREqnn6N27N5ydnXH27Fl4e3tLXQ6Rieu/kJaRAcTEANnZYtg2D1wAMBrF17KzxXEZGS1OUV5ejtmzZyM8PBwnT57Enj17kJmZycClLuO6LsnV9YVuRgaQmAhUVwPtPQxQEMRxiYlXgre0tBRvvPEGIiIiUF5ejgMHDmDdunUIDQ29rrKIuK5LctX10N2372rgNunVC1i3DigsBCorgf/9XyA21vT7qqvROHMm3nniCURGRuLixYs4dOgQ1qxZg+Dg4C6XQ9QcZ7okV11f001PF5cOTM7mAJw4AYweDRQVARMmAJ9/DgweDPz555VhgtGIu/fsweScHPTv37/LJRC1hjNdkquuzXRLS8WLZtcuKVRXA2+9JQasIABbtgAFBcCtt5oMswcw9ORJ9OfDHqmbcKZLctW10M3M7Ng4X18gIgI4erTlMZWq4+ch6iTOdEmuuha6R4603KVwLQcH4OOPgX/9C8jNbXncaARycrr09kTt4UyX5KproVtR0fZxlQrIygJqa4Hp01sddqGwEFVVVV0qgagtzW+QIJKTrl1I8/Bo+/j77wN+fuKFtPr6Vof9kpODx/r2hbu7OxwOQcEAAAbQSURBVCIiIhAeHm7yMTQ0FM7Ozl0qkXo2V1dXqNVq3iBBstO10I2KAjZtMr/EkJEB3HQTMHZs20sQajXiZs1C1cyZKCkpgU6ng16vh06nw6+//gq9Xo+CggL4+fmZDeTg4GA4Ojp2qXzqGZrWdRm6JCeqtv73Kzo6Wti/f3/LA6WlQFBQy1ANDBR3LtTUmM5wX3gB+OQT07HOzuK2sjZ6MtTX16OoqMgkkJs+lpSUICAgwGwgBwQE8LHqhIkTJ+LFF19EXFyc1KVQD6NSqQ4IghBt7ljXZrq+vmIvhexs021jRUXiem77FYlLD+00wXFwcEBISAhCQkIQe81NFpcvX0ZBQcGVID58+DC++OIL6HQ6nD17FiEhIWYDuW/fvmyC0kNwBwPJUddvjkhKArZvN70jraPUavH7r4OTkxMiIyMRGRnZ4tilS5dgMBiuBPKuXbuQmZkJnU4Ho9GI8PDwK7+aBzIfZGhbuIOB5KjroTt8OLBkSctbgdvj4iJ+X7TZmbdFuLq6IioqClFRUS2OVVRUmCxVbN++HatWrYJOp4Odnd2VEG4eyOHh4XB3d++2eql7aDQa7NixQ+oyiExcX2vHhATxY2KiuO+2re05KpU4w12y5Or3ScDDwwPR0dGIvib0BUFAeXm5SSBv2rQJer0eer0ebm5uLWbG4eHhCAsLg1qtluinobZwpktydP39dBMSxFlvejqwdasYrs17MqjVYhhPmCAuKXTjDPd6qFQq+Pj4wMfHB6NGjTI5JgiCyQ4LvV6Pf/3rX9Dr9cjPz4evr6/ZGfKAAQPQi7c6S4ZruiRHXdu90JqyMvHW3pwc4Px5wNNTbHYzbVqnnhyhJA0NDSY7LJrPlIuLi6HRaMwGcmBgIOzt7aUu36ZdunQJXl5eMBqNvHhKVtXW7gXLhi6ZqK2tNdlh0TyQy8rKMGDAALOB7O/vz5CwkBtuuAG5ubnwsdF/9EmeLL9ljDqkV69euPHGG3HjjTe2OFZdXW2yw2LPnj3IysqCXq/HxYsXERYWZjaQvb29Gcid0LSuy9AluWDoSsTFxQWDBw/G4MGDWxyrrKw0mRV///33yMjIgE6nAwCzF/TCw8Ph0d7t2T1Q07rusGHDpC6FCABDV5bc3d1x66234tZr+hALgoCzZ8+aBPKXX355ZenC1dW11R0WLi4uEv000uIOBpIbhq6CqFQqeHt7w9vbG7fddpvJMUEQcOrUKZNAblquyM/Ph7e3t9lADgkJsekdFtzBQHLD0LURKpUK/v7+8Pf3x+jRo02ONTQ04MSJEyaBvGPHDuj1ehQVFUGj0ZgN5KCgIGXvsCgtRWxODi7v3Qv8/rvYHS8qCnj2WZvdTUPyx90LPVxdXR0KCgpaNBTS6/U4c+YMBgwYYDaQ/f395dtUaN8+cd/4tm1oaGyEfW3t1WNN+8bHjxf3jQ8fLl2dZLO4ZYy6xGg0wmAwmA3kyspKhIWFmQ1kHx8f6XZYZGQo6g5Jsk0MXbK4yspK5OXlmQ3khoYGsw2FwsPD0adPn+4rqilwu9ILhMFLFsTQJatq2mFhLpDVarXZQA4LC4Orq2vX33TfPiAmpmXg/uMf4h2RgwcDGzaI67nXcnEBdu6U7S3qpDy8OYKsysvLC15eXhg5cqTJ64Ig4PTp0yaBvGHDBuh0OhgMBnh5ebXYexwREYGQkBA4OTm1/abp6aY9P5qUlACpqcC4ceJygjlGo/j9mzZ18Scm6jjOdEkWGhsbr+ywuHaGXFRUBH9/f7OBHBQUBIdz58w/yaS5lBRAozE/0wU69CQToo7iTJdkz87ODkFBQQgKCsLYsWNNjtXV1aGwsPBKIOfm5uLbb7+FTqfD6dOnkeLmhldqa9HOXLhtKpXYrOm1167nLETtYuiS7Dk6Ol6Z3V6rpqYGNY8+Cqdvv72+NzEaxe54RN2MoUuK5uzsDGdLnez8eUudiahVMt3dTtQJlmr04+lpmfMQtYGhS8oXFSVeCDPH3h5wchI/Nv/8Wmq1uK2MqJsxdEn5pk1r/dicOeKuhqQk4KmnxM/nzGk5ThDaPg+RhTB0Sfl8fcVeCuZuPX7rLfH15r/eest0jEolPsOP28XIChi6ZBuSklq/+aE9arX4/URWwNAl2zB8uNhDobPN2pt6L/AWYLISbhkj29HUtIZdxkjGONMl25KQIDaviY8XdzRcu+SgVouvx8eL4xi4ZGWc6ZLtiY4Wm9eUlYm39ubkiDc+eHqK28KmTeNFM5JMmw1vVCpVGYA/rVcOEZFNCBIEwey/7G2GLhERWRbXdImIrIihS0RkRQxdIiIrYugSEVkRQ5eIyIr+Dy7yr+n/eUR1AAAAAElFTkSuQmCC\n",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": []
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "QSmBdaxcX_PP",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 156
        },
        "outputId": "a6c0cfd8-a43a-4dda-f06b-7703356e8306"
      },
      "source": [
        "vertex_cover = VertexCover(G_three)\n",
        "vertex_cover._energy_function_h_a()\n",
        "vertex_cover._energy_function_h_b()\n",
        "vertex_cover._solve(num_repeats=1000)"
      ],
      "execution_count": 63,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Setting A:10.0 and B:1.0\n",
            "Setting the solver as QBsolv\n",
            "defaultdict(<class 'int'>, {0: -39.0, 1: -19.0, 2: -19.0, 3: -19.0, 4: -19.0})\n",
            "defaultdict(<class 'int'>, {(0, 1): 10.0, (0, 2): 10.0, (0, 3): 10.0, (0, 4): 10.0, (1, 2): 10.0, (3, 4): 10.0})\n",
            "60.0\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "SampleSet(rec.array([([1, 1, 0, 0, 1], 3., 236), ([1, 0, 1, 1, 0], 3., 242),\n",
              "           ([1, 0, 1, 0, 1], 3., 250), ([1, 1, 0, 1, 0], 3., 273)],\n",
              "          dtype=[('sample', 'i1', (5,)), ('energy', '<f8'), ('num_occurrences', '<i8')]), [0, 1, 2, 3, 4], {}, 'BINARY')"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 63
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "KVkEHuLBHnXt",
        "colab_type": "text"
      },
      "source": [
        "# References:\n",
        "\n",
        "[1] Cormen, Thomas H., Charles E. Leiserson, Ronald L. Rivest, and Clifford Stein. *Introduction to algorithms.* MIT press, 2009.\n",
        "\n",
        "[2] Lucas, Andrew. \"Ising formulations of many NP problems.\" *Frontiers in Physics 2* (2014): 5.\n"
      ]
    }
  ]
}