sndyuk · January 23, 2025 16:43
diff --git a/qubo_tutorial.ipynb b/qubo_tutorial.ipynb
 {
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "QUBO_Tutorial.ipynb",
      "version": "0.3.2",
      "provenance": [],
      "collapsed_sections": [],
      "include_colab_link": true
    },
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/gist/sndyuk/40cd214b44d14b264ae1798ff0e64964/qubo_tutorial.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "UyJ6ISvF9gwn",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "import numpy as np"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "-avWc61w9gwr",
        "colab_type": "text"
      },
      "source": [
        "# 1. QUBO expression\n",
        "\n",
        "- is written as a 2 dimentional matrix and binary variables.\n",
        "- needs to be upper triangular matrix because it's an one of expression of quadratic form.\n",
        "\n",
        "Sample quadratic form with 3 variables(x, y, z) and the coeficients(10 for (x, x) => x^2, 20 for y^2, 30 for z^2, 2 for (x, y), 1 for (x, z)) .\n",
        "\n",
        "\\begin{equation}\n",
        "A = 10x^2 + 20y^2 + 30z^2 + \\frac{2+2}{2}xy + \\frac{1+1}{2}xz\n",
        "\\end{equation}\n",
        "\n",
        "Let's convert this to matrix.\n",
        "\n",
        "\\begin{equation}\n",
        "  A = \\left(\n",
        "    \\begin{array}{ccc}\n",
        "      10 & 2 & 1 \\\\\n",
        "      2 & 20 & 0 \\\\\n",
        "      1 & 0 & 30\n",
        "    \\end{array}\n",
        "  \\right)\n",
        "  \\left(\n",
        "    \\begin{array}{ccc}\n",
        "      x \\\\\n",
        "      y \\\\\n",
        "      z\n",
        "    \\end{array}\n",
        "  \\right)\n",
        "\\end{equation}\n",
        "\n",
        "Then fold lower off-diagonal to upper side to make it simpler since this is symmetric matrix.\n",
        "\n",
        "\\begin{equation}\n",
        "  A = \\left(\n",
        "    \\begin{array}{ccc}\n",
        "      10 & 4 & 2 \\\\\n",
        "      0 & 20 & 0 \\\\\n",
        "      0 & 0 & 30\n",
        "    \\end{array}\n",
        "  \\right)\n",
        "    \\left(\n",
        "    \\begin{array}{ccc}\n",
        "      x \\\\\n",
        "      y \\\\\n",
        "      z\n",
        "    \\end{array}\n",
        "  \\right)\n",
        "\\end{equation}\n"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Lgvsag8Q9gws",
        "colab_type": "code",
        "outputId": "012459cc-a680-4ed2-b67a-4f1f3479829b",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        }
      },
      "source": [
        "## Quadratic form matrix\n",
        "Q = [\n",
        "    [10, 4, 2],\n",
        "    [0, 20, 0],\n",
        "    [0, 0, 30]\n",
        "]\n",
        "\n",
        "# QUBO consists of the Quadratic form and its variables(= solution) which are all binary.\n",
        "# - The 3x3 QUBO has 3 bits solution space(X).\n",
        "#       4x4 QUBO has 4 bits.\n",
        "#       5x5 QUBO has 5 bits and so on.\n",
        "\n",
        "# The solution is written as vector(or 1(col)xN(row) dimentional matrix).\n",
        "#    x: 0 or 1\n",
        "X = [1, 1, 0]\n",
        "\n",
        "\n",
        "# - The energy of the solution is evaluated by the logic.\n",
        "#   It's a summation of each rows and columns.\n",
        "energy = .0\n",
        "SOLUTION = np.array(X)\n",
        "for i, c in enumerate(SOLUTION):\n",
        "    for j, d in enumerate(SOLUTION[i:]):\n",
        "        if c and d:\n",
        "            energy += Q[i][j + i]\n",
        "print(energy)"
      ],
      "execution_count": 5,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "34.0\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "iO_eZ6uG9gww",
        "colab_type": "code",
        "outputId": "2620a1e0-4529-4540-c378-6e3332e42330",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        }
      },
      "source": [
        "# Also, it can be written as the mathematical formula.\n",
        "# energy = \\inv{X} * Q * X\n",
        "\n",
        "_X = np.matrix(X).T # To 1 col, n row matrix.\n",
        "_Xt = _X.T\n",
        "_Q = np.matrix(Q)\n",
        "\n",
        "energy = _Xt * _Q * _X\n",
        "# OR: (_Xt.dot(_Q)*_X).sum(axis=1)\n",
        "\n",
        "print(float(energy))"
      ],
      "execution_count": 6,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "34.0\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "faKMpNcj9gw0",
        "colab_type": "text"
      },
      "source": [
        "# 2. Graph to QUBO"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "erzBquST9gw0",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "# An element of the QUBO can be a node of the Graph.\n",
        "\n",
        "# In this example, there are 3 nodes:\n",
        "#    node a have weight 10.\n",
        "#    node b have weight 0.\n",
        "#    node c have weight 5.\n",
        "# Every weight of the nodes are written as diagonal matrix:\n",
        "Q = [\n",
        "    [10, 0, 0],\n",
        "    [0, 0, 0],\n",
        "    [0, 0, 5]\n",
        "]\n",
        "\n",
        "# There are some relations between each nodes:\n",
        "#    a-b: 400\n",
        "#    a-c: 3\n",
        "#    b-c: 200\n",
        "# Every relation between each node are written as off-diagonal matrix:\n",
        "Q = [\n",
        "    [10, 400, 3],\n",
        "    [0, 0, 200],\n",
        "    [0, 0, 5]\n",
        "]\n",
        "\n",
        "# Important observation:\n",
        "# If there is the relation below:\n",
        "#    c-a: 500 which means the edge has the specific direction(c->a).\n",
        "# The matrix would be this but QUBO can't evaluate this kind of matrix because QUBO must be upper triangler matrix:\n",
        "Q = [\n",
        "    [10, 400, 3],\n",
        "    [0, 0, 200],\n",
        "    [500, 0, 5]\n",
        "]\n",
        "# It needs to be this:\n",
        "#    a-c: +500 which means the edge doesn't have the specific direction.\n",
        "Q = [\n",
        "    [10, 400, 3 + 500],\n",
        "    [0, 0, 200],\n",
        "    [0, 0, 5]\n",
        "]\n",
        "# we need to expand the solution space to handle the backward direction. E.g,\n",
        "Q2 = [\n",
        "    # Forward direction(a->c)\n",
        "    [10, 400, 3,   0, 0, 0],\n",
        "    [0,  0,   200, 0, 0, 0],\n",
        "    [0,  0,   5,   0, 0, 0],\n",
        "    # Backward direction(c->a)\n",
        "    [0,  0,   0,   0, 0, 500],\n",
        "    [0,  0,   0,   0,  0, 0],\n",
        "    [0,  0,   0,   0,  0, 0],\n",
        "]\n",
        "\n",
        "# The weight and the relation can be called 'cost'. And the summation of the cost can be called 'energy'."
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "em2Bpxjw9gw3",
        "colab_type": "text"
      },
      "source": [
        "# 3.  QUBO expression part.2\n",
        "\n",
        "\\begin{equation}\n",
        "\\sum_{ij} (\\sum_{k=1}^{9} x_{ijk} - 1)^2\n",
        "\\end{equation}\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "KLkSkqUWO0D3",
        "colab_type": "text"
      },
      "source": [
        "### The formula to matrix\n",
        "\n",
        "\\begin{equation}\n",
        "\\sum_{ij} (\\sum_{k=1}^{9} x_{ijk} - 1)^2\\\\\n",
        "= \\sum_{ij} ((\\sum x_{ijk})^2 + (2 \\times (\\sum x_{ijk}) \\times -1) + (-1)^2)\\\\\n",
        "\\end{equation}\n",
        "\n",
        "For example, calculate the first part of the formula$$(\\sum x_{ijk})^2$$ using the sample quadratic form.\n",
        "\n",
        "\n",
        "Sample quadratic form:\n",
        "\\begin{equation}\n",
        "  A = \\left(\n",
        "    \\begin{array}{ccc}\n",
        "      10 & 0 & 0 \\\\\n",
        "      0 & 20 & 0 \\\\\n",
        "      0 & 0 & 30\n",
        "    \\end{array}\n",
        "  \\right)\n",
        "\\end{equation}\n",
        "\n",
        "This can be\n",
        "\n",
        "\\begin{equation}\n",
        "= (10x_{0} + 20x_{1} + 30x_{2})^2\\\\\n",
        "= (10x_{0} + 20x_{1} + 30x_{2})(10x_{0} + 20x_{1} + 30x_{2})\\\\\n",
        "= 100x_{0}^2 + 200x_{0}x_{1} + 300x_{0}x_{2} + ... + 900x_{2}^2\\\\\n",
        "\\end{equation}\n",
        "\n",
        "when x is binary, x^2 and x are equivalent. \n",
        "\n",
        "\\begin{equation}\n",
        "= 100x_{0} + 200x_{0}x_{1} + 300x_{0}x_{2} + ... + 900x_{2}\\\\\n",
        "\\end{equation}\n",
        "\n",
        "To matrix:\n",
        "\n",
        "\\begin{equation}\n",
        "  A = \\left(\n",
        "    \\begin{array}{ccc}\n",
        "      100 & 200 & 300 \\\\\n",
        "      200 & 400 & 600 \\\\\n",
        "      300 & 600 & 900\n",
        "    \\end{array}\n",
        "  \\right)\n",
        "\\end{equation}\n",
        "\n",
        "To upper trianguler matrix:\n",
        "\n",
        "\\begin{equation}\n",
        "  A = \\left(\n",
        "    \\begin{array}{ccc}\n",
        "      100 & 400 & 600 \\\\\n",
        "      0 & 400 & 1200 \\\\\n",
        "      0 & 0 & 900\n",
        "    \\end{array}\n",
        "  \\right)\n",
        "\\end{equation}\n",
        "\n",
        "### To python"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "9mxe3yN-9gw4",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "# This part explain how we build QUBO from the math formula. Using SUDOKU example here: https://www.topcoder.com/challenges/30081256\n",
        "def get_variable_id(N, i, j, k):\n",
        "    return (i * N + j) * N + k\n",
        "\n",
        "N = 9\n",
        "solution_space = N * N * N"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "WA_szQg59gw7",
        "colab_type": "code",
        "outputId": "6af6603b-cf8e-470d-b09b-d9592e109bd6",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 153
        }
      },
      "source": [
        "# Example code 1. It's slow but descriptive.\n",
        "Q = np.zeros((N * N * N, N * N * N))\n",
        "constant = 0\n",
        "\n",
        "# \\sum_{ij}\n",
        "for i in range(N):\n",
        "    for j in range(N):\n",
        "        # \\sum x_{ijk} - 1\n",
        "        Q_tmp = np.zeros((N * N * N, N * N * N))\n",
        "        constant_tmp = 0\n",
        "        for k in range(N):\n",
        "            ijk = get_variable_id(N, i, j, k)\n",
        "            Q_tmp[ijk][ijk] += 1\n",
        "        constant_tmp += -1\n",
        "\n",
        "        # (\\sum x_{ijk} - 1)^2\n",
        "        # = (\\sum x_{ijk})^2 + (2 * x_{ijk} * -1) + -1^2\n",
        "        #   part.1             part.2               part.3\n",
        "\n",
        "        # part.1\n",
        "        diagonal = np.diag(Q_tmp)\n",
        "        A = np.outer(diagonal, diagonal)\n",
        "        A = (np.triu(A, k=1) + np.triu(A.T)) # Fold lower triangular to upper side.\n",
        "        \n",
        "        # part.2\n",
        "        # 2 * x_{ijk} * -1\n",
        "        B = (2 * np.diag(diagonal) * constant_tmp)\n",
        "        Q_tmp = A + B\n",
        "\n",
        "        # part.3\n",
        "        # -1^2\n",
        "        constant_tmp = constant_tmp ** 2\n",
        "        \n",
        "        Q += Q_tmp\n",
        "        constant += constant_tmp\n",
        "\n",
        "print(Q)\n",
        "print('constant:', constant)"
      ],
      "execution_count": 9,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "[[-1.  2.  2. ...  0.  0.  0.]\n",
            " [ 0. -1.  2. ...  0.  0.  0.]\n",
            " [ 0.  0. -1. ...  0.  0.  0.]\n",
            " ...\n",
            " [ 0.  0.  0. ... -1.  2.  2.]\n",
            " [ 0.  0.  0. ...  0. -1.  2.]\n",
            " [ 0.  0.  0. ...  0.  0. -1.]]\n",
            "constant: 81\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "tNc51KU-9gw_",
        "colab_type": "code",
        "outputId": "68c2a0fe-6513-4c3f-e842-f17be90e777c",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 153
        }
      },
      "source": [
        "# Example code 2. It's fast.\n",
        "\n",
        "Q = np.zeros((N * N * N, N * N * N))\n",
        "constant = 0\n",
        "\n",
        "# \\sum_{ij}\n",
        "for i in range(N):\n",
        "    for j in range(N):\n",
        "        # (\\sum x_{ijk} - 1)^2 for each ij\n",
        "        # = (\\sum x_{ijk})^2 + (2 * x_{ijk} * -1) + -1^2\n",
        "        #   part.1             part.2               part.3\n",
        "\n",
        "        for k1 in range(N):\n",
        "            # part.1\n",
        "            ijk1 = get_variable_id(N, i, j, k1)\n",
        "            for k2 in range(N):\n",
        "                ijk2 = get_variable_id(N, i, j, k2)\n",
        "                if ijk1 <= ijk2:\n",
        "                    Q[ijk1][ijk2] += 1\n",
        "                else:\n",
        "                    # Keep it as upper triangular matrix.\n",
        "                    Q[ijk2][ijk1] += 1\n",
        "\n",
        "            # part.2\n",
        "            # 2 * x_{ijk} * -1\n",
        "            Q[ijk1][ijk1] += -2\n",
        "\n",
        "        # part.3\n",
        "        # -1^2\n",
        "        # constant only affects the total energy. It's not a part of QUBO matrix.\n",
        "        constant += 1\n",
        "\n",
        "print(Q)\n",
        "print('constant:', constant)"
      ],
      "execution_count": 10,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "[[-1.  2.  2. ...  0.  0.  0.]\n",
            " [ 0. -1.  2. ...  0.  0.  0.]\n",
            " [ 0.  0. -1. ...  0.  0.  0.]\n",
            " ...\n",
            " [ 0.  0.  0. ... -1.  2.  2.]\n",
            " [ 0.  0.  0. ...  0. -1.  2.]\n",
            " [ 0.  0.  0. ...  0.  0. -1.]]\n",
            "constant: 81\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Kj4-1mli9gxC",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        ""
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "xpEEdqA29gxF",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        ""
      ],
      "execution_count": 0,
      "outputs": []
    }
  ]
 }
	{
	"nbformat": 4,
	"nbformat_minor": 0,
	"metadata": {
	"colab": {
	"name": "QUBO_Tutorial.ipynb",
	"version": "0.3.2",
	"provenance": [],
	"collapsed_sections": [],
	"include_colab_link": true
	},
	"kernelspec": {
	"display_name": "Python 3",
	"language": "python",
	"name": "python3"
	}
	},
	"cells": [
	{
	"cell_type": "markdown",
	"metadata": {
	"id": "view-in-github",
	"colab_type": "text"
	},
	"source": [
	"<a href=\"https://colab.research.google.com/gist/sndyuk/40cd214b44d14b264ae1798ff0e64964/qubo_tutorial.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
	]
	},
	{
	"cell_type": "code",
	"metadata": {
	"id": "UyJ6ISvF9gwn",
	"colab_type": "code",
	"colab": {}
	},
	"source": [
	"import numpy as np"
	],
	"execution_count": 0,
	"outputs": []
	},
	{
	"cell_type": "markdown",
	"metadata": {
	"id": "-avWc61w9gwr",
	"colab_type": "text"
	},
	"source": [
	"# 1. QUBO expression\n",
	"\n",
	"- is written as a 2 dimentional matrix and binary variables.\n",
	"- needs to be upper triangular matrix because it's an one of expression of quadratic form.\n",
	"\n",
	"Sample quadratic form with 3 variables(x, y, z) and the coeficients(10 for (x, x) => x^2, 20 for y^2, 30 for z^2, 2 for (x, y), 1 for (x, z)) .\n",
	"\n",
	"\\begin{equation}\n",
	"A = 10x^2 + 20y^2 + 30z^2 + \\frac{2+2}{2}xy + \\frac{1+1}{2}xz\n",
	"\\end{equation}\n",
	"\n",
	"Let's convert this to matrix.\n",
	"\n",
	"\\begin{equation}\n",
	" A = \\left(\n",
	" \\begin{array}{ccc}\n",
	" 10 & 2 & 1 \\\\\n",
	" 2 & 20 & 0 \\\\\n",
	" 1 & 0 & 30\n",
	" \\end{array}\n",
	" \\right)\n",
	" \\left(\n",
	" \\begin{array}{ccc}\n",
	" x \\\\\n",
	" y \\\\\n",
	" z\n",
	" \\end{array}\n",
	" \\right)\n",
	"\\end{equation}\n",
	"\n",
	"Then fold lower off-diagonal to upper side to make it simpler since this is symmetric matrix.\n",
	"\n",
	"\\begin{equation}\n",
	" A = \\left(\n",
	" \\begin{array}{ccc}\n",
	" 10 & 4 & 2 \\\\\n",
	" 0 & 20 & 0 \\\\\n",
	" 0 & 0 & 30\n",
	" \\end{array}\n",
	" \\right)\n",
	" \\left(\n",
	" \\begin{array}{ccc}\n",
	" x \\\\\n",
	" y \\\\\n",
	" z\n",
	" \\end{array}\n",
	" \\right)\n",
	"\\end{equation}\n"
	]
	},
	{
	"cell_type": "code",
	"metadata": {
	"id": "Lgvsag8Q9gws",
	"colab_type": "code",
	"outputId": "012459cc-a680-4ed2-b67a-4f1f3479829b",
	"colab": {
	"base_uri": "https://localhost:8080/",
	"height": 34
	}
	},
	"source": [
	"## Quadratic form matrix\n",
	"Q = [\n",
	" [10, 4, 2],\n",
	" [0, 20, 0],\n",
	" [0, 0, 30]\n",
	"]\n",
	"\n",
	"# QUBO consists of the Quadratic form and its variables(= solution) which are all binary.\n",
	"# - The 3x3 QUBO has 3 bits solution space(X).\n",
	"# 4x4 QUBO has 4 bits.\n",
	"# 5x5 QUBO has 5 bits and so on.\n",
	"\n",
	"# The solution is written as vector(or 1(col)xN(row) dimentional matrix).\n",
	"# x: 0 or 1\n",
	"X = [1, 1, 0]\n",
	"\n",
	"\n",
	"# - The energy of the solution is evaluated by the logic.\n",
	"# It's a summation of each rows and columns.\n",
	"energy = .0\n",
	"SOLUTION = np.array(X)\n",
	"for i, c in enumerate(SOLUTION):\n",
	" for j, d in enumerate(SOLUTION[i:]):\n",
	" if c and d:\n",
	" energy += Q[i][j + i]\n",
	"print(energy)"
	],
	"execution_count": 5,
	"outputs": [
	{
	"output_type": "stream",
	"text": [
	"34.0\n"
	],
	"name": "stdout"
	}
	]
	},
	{
	"cell_type": "code",
	"metadata": {
	"id": "iO_eZ6uG9gww",
	"colab_type": "code",
	"outputId": "2620a1e0-4529-4540-c378-6e3332e42330",
	"colab": {
	"base_uri": "https://localhost:8080/",
	"height": 34
	}
	},
	"source": [
	"# Also, it can be written as the mathematical formula.\n",
	"# energy = \\inv{X} * Q * X\n",
	"\n",
	"_X = np.matrix(X).T # To 1 col, n row matrix.\n",
	"_Xt = _X.T\n",
	"_Q = np.matrix(Q)\n",
	"\n",
	"energy = _Xt * _Q * _X\n",
	"# OR: (_Xt.dot(_Q)*_X).sum(axis=1)\n",
	"\n",
	"print(float(energy))"
	],
	"execution_count": 6,
	"outputs": [
	{
	"output_type": "stream",
	"text": [
	"34.0\n"
	],
	"name": "stdout"
	}
	]
	},
	{
	"cell_type": "markdown",
	"metadata": {
	"id": "faKMpNcj9gw0",
	"colab_type": "text"
	},
	"source": [
	"# 2. Graph to QUBO"
	]
	},
	{
	"cell_type": "code",
	"metadata": {
	"id": "erzBquST9gw0",
	"colab_type": "code",
	"colab": {}
	},
	"source": [
	"# An element of the QUBO can be a node of the Graph.\n",
	"\n",
	"# In this example, there are 3 nodes:\n",
	"# node a have weight 10.\n",
	"# node b have weight 0.\n",
	"# node c have weight 5.\n",
	"# Every weight of the nodes are written as diagonal matrix:\n",
	"Q = [\n",
	" [10, 0, 0],\n",
	" [0, 0, 0],\n",
	" [0, 0, 5]\n",
	"]\n",
	"\n",
	"# There are some relations between each nodes:\n",
	"# a-b: 400\n",
	"# a-c: 3\n",
	"# b-c: 200\n",
	"# Every relation between each node are written as off-diagonal matrix:\n",
	"Q = [\n",
	" [10, 400, 3],\n",
	" [0, 0, 200],\n",
	" [0, 0, 5]\n",
	"]\n",
	"\n",
	"# Important observation:\n",
	"# If there is the relation below:\n",
	"# c-a: 500 which means the edge has the specific direction(c->a).\n",
	"# The matrix would be this but QUBO can't evaluate this kind of matrix because QUBO must be upper triangler matrix:\n",
	"Q = [\n",
	" [10, 400, 3],\n",
	" [0, 0, 200],\n",
	" [500, 0, 5]\n",
	"]\n",
	"# It needs to be this:\n",
	"# a-c: +500 which means the edge doesn't have the specific direction.\n",
	"Q = [\n",
	" [10, 400, 3 + 500],\n",
	" [0, 0, 200],\n",
	" [0, 0, 5]\n",
	"]\n",
	"# we need to expand the solution space to handle the backward direction. E.g,\n",
	"Q2 = [\n",
	" # Forward direction(a->c)\n",
	" [10, 400, 3, 0, 0, 0],\n",
	" [0, 0, 200, 0, 0, 0],\n",
	" [0, 0, 5, 0, 0, 0],\n",
	" # Backward direction(c->a)\n",
	" [0, 0, 0, 0, 0, 500],\n",
	" [0, 0, 0, 0, 0, 0],\n",
	" [0, 0, 0, 0, 0, 0],\n",
	"]\n",
	"\n",
	"# The weight and the relation can be called 'cost'. And the summation of the cost can be called 'energy'."
	],
	"execution_count": 0,
	"outputs": []
	},
	{
	"cell_type": "markdown",
	"metadata": {
	"id": "em2Bpxjw9gw3",
	"colab_type": "text"
	},
	"source": [
	"# 3. QUBO expression part.2\n",
	"\n",
	"\\begin{equation}\n",
	"\\sum_{ij} (\\sum_{k=1}^{9} x_{ijk} - 1)^2\n",
	"\\end{equation}\n"
	]
	},
	{
	"cell_type": "markdown",
	"metadata": {
	"id": "KLkSkqUWO0D3",
	"colab_type": "text"
	},
	"source": [
	"### The formula to matrix\n",
	"\n",
	"\\begin{equation}\n",
	"\\sum_{ij} (\\sum_{k=1}^{9} x_{ijk} - 1)^2\\\\\n",
	"= \\sum_{ij} ((\\sum x_{ijk})^2 + (2 \\times (\\sum x_{ijk}) \\times -1) + (-1)^2)\\\\\n",
	"\\end{equation}\n",
	"\n",
	"For example, calculate the first part of the formula$$(\\sum x_{ijk})^2$$ using the sample quadratic form.\n",
	"\n",
	"\n",
	"Sample quadratic form:\n",
	"\\begin{equation}\n",
	" A = \\left(\n",
	" \\begin{array}{ccc}\n",
	" 10 & 0 & 0 \\\\\n",
	" 0 & 20 & 0 \\\\\n",
	" 0 & 0 & 30\n",
	" \\end{array}\n",
	" \\right)\n",
	"\\end{equation}\n",
	"\n",
	"This can be\n",
	"\n",
	"\\begin{equation}\n",
	"= (10x_{0} + 20x_{1} + 30x_{2})^2\\\\\n",
	"= (10x_{0} + 20x_{1} + 30x_{2})(10x_{0} + 20x_{1} + 30x_{2})\\\\\n",
	"= 100x_{0}^2 + 200x_{0}x_{1} + 300x_{0}x_{2} + ... + 900x_{2}^2\\\\\n",
	"\\end{equation}\n",
	"\n",
	"when x is binary, x^2 and x are equivalent. \n",
	"\n",
	"\\begin{equation}\n",
	"= 100x_{0} + 200x_{0}x_{1} + 300x_{0}x_{2} + ... + 900x_{2}\\\\\n",
	"\\end{equation}\n",
	"\n",
	"To matrix:\n",
	"\n",
	"\\begin{equation}\n",
	" A = \\left(\n",
	" \\begin{array}{ccc}\n",
	" 100 & 200 & 300 \\\\\n",
	" 200 & 400 & 600 \\\\\n",
	" 300 & 600 & 900\n",
	" \\end{array}\n",
	" \\right)\n",
	"\\end{equation}\n",
	"\n",
	"To upper trianguler matrix:\n",
	"\n",
	"\\begin{equation}\n",
	" A = \\left(\n",
	" \\begin{array}{ccc}\n",
	" 100 & 400 & 600 \\\\\n",
	" 0 & 400 & 1200 \\\\\n",
	" 0 & 0 & 900\n",
	" \\end{array}\n",
	" \\right)\n",
	"\\end{equation}\n",
	"\n",
	"### To python"
	]
	},
	{
	"cell_type": "code",
	"metadata": {
	"id": "9mxe3yN-9gw4",
	"colab_type": "code",
	"colab": {}
	},
	"source": [
	"# This part explain how we build QUBO from the math formula. Using SUDOKU example here: https://www.topcoder.com/challenges/30081256\n",
	"def get_variable_id(N, i, j, k):\n",
	" return (i * N + j) * N + k\n",
	"\n",
	"N = 9\n",
	"solution_space = N * N * N"
	],
	"execution_count": 0,
	"outputs": []
	},
	{
	"cell_type": "code",
	"metadata": {
	"id": "WA_szQg59gw7",
	"colab_type": "code",
	"outputId": "6af6603b-cf8e-470d-b09b-d9592e109bd6",
	"colab": {
	"base_uri": "https://localhost:8080/",
	"height": 153
	}
	},
	"source": [
	"# Example code 1. It's slow but descriptive.\n",
	"Q = np.zeros((N * N * N, N * N * N))\n",
	"constant = 0\n",
	"\n",
	"# \\sum_{ij}\n",
	"for i in range(N):\n",
	" for j in range(N):\n",
	" # \\sum x_{ijk} - 1\n",
	" Q_tmp = np.zeros((N * N * N, N * N * N))\n",
	" constant_tmp = 0\n",
	" for k in range(N):\n",
	" ijk = get_variable_id(N, i, j, k)\n",
	" Q_tmp[ijk][ijk] += 1\n",
	" constant_tmp += -1\n",
	"\n",
	" # (\\sum x_{ijk} - 1)^2\n",
	" # = (\\sum x_{ijk})^2 + (2 * x_{ijk} * -1) + -1^2\n",
	" # part.1 part.2 part.3\n",
	"\n",
	" # part.1\n",
	" diagonal = np.diag(Q_tmp)\n",
	" A = np.outer(diagonal, diagonal)\n",
	" A = (np.triu(A, k=1) + np.triu(A.T)) # Fold lower triangular to upper side.\n",
	" \n",
	" # part.2\n",
	" # 2 * x_{ijk} * -1\n",
	" B = (2 * np.diag(diagonal) * constant_tmp)\n",
	" Q_tmp = A + B\n",
	"\n",
	" # part.3\n",
	" # -1^2\n",
	" constant_tmp = constant_tmp ** 2\n",
	" \n",
	" Q += Q_tmp\n",
	" constant += constant_tmp\n",
	"\n",
	"print(Q)\n",
	"print('constant:', constant)"
	],
	"execution_count": 9,
	"outputs": [
	{
	"output_type": "stream",
	"text": [
	"[[-1. 2. 2. ... 0. 0. 0.]\n",
	" [ 0. -1. 2. ... 0. 0. 0.]\n",
	" [ 0. 0. -1. ... 0. 0. 0.]\n",
	" ...\n",
	" [ 0. 0. 0. ... -1. 2. 2.]\n",
	" [ 0. 0. 0. ... 0. -1. 2.]\n",
	" [ 0. 0. 0. ... 0. 0. -1.]]\n",
	"constant: 81\n"
	],
	"name": "stdout"
	}
	]
	},
	{
	"cell_type": "code",
	"metadata": {
	"id": "tNc51KU-9gw_",
	"colab_type": "code",
	"outputId": "68c2a0fe-6513-4c3f-e842-f17be90e777c",
	"colab": {
	"base_uri": "https://localhost:8080/",
	"height": 153
	}
	},
	"source": [
	"# Example code 2. It's fast.\n",
	"\n",
	"Q = np.zeros((N * N * N, N * N * N))\n",
	"constant = 0\n",
	"\n",
	"# \\sum_{ij}\n",
	"for i in range(N):\n",
	" for j in range(N):\n",
	" # (\\sum x_{ijk} - 1)^2 for each ij\n",
	" # = (\\sum x_{ijk})^2 + (2 * x_{ijk} * -1) + -1^2\n",
	" # part.1 part.2 part.3\n",
	"\n",
	" for k1 in range(N):\n",
	" # part.1\n",
	" ijk1 = get_variable_id(N, i, j, k1)\n",
	" for k2 in range(N):\n",
	" ijk2 = get_variable_id(N, i, j, k2)\n",
	" if ijk1 <= ijk2:\n",
	" Q[ijk1][ijk2] += 1\n",
	" else:\n",
	" # Keep it as upper triangular matrix.\n",
	" Q[ijk2][ijk1] += 1\n",
	"\n",
	" # part.2\n",
	" # 2 * x_{ijk} * -1\n",
	" Q[ijk1][ijk1] += -2\n",
	"\n",
	" # part.3\n",
	" # -1^2\n",
	" # constant only affects the total energy. It's not a part of QUBO matrix.\n",
	" constant += 1\n",
	"\n",
	"print(Q)\n",
	"print('constant:', constant)"
	],
	"execution_count": 10,
	"outputs": [
	{
	"output_type": "stream",
	"text": [
	"[[-1. 2. 2. ... 0. 0. 0.]\n",
	" [ 0. -1. 2. ... 0. 0. 0.]\n",
	" [ 0. 0. -1. ... 0. 0. 0.]\n",
	" ...\n",
	" [ 0. 0. 0. ... -1. 2. 2.]\n",
	" [ 0. 0. 0. ... 0. -1. 2.]\n",
	" [ 0. 0. 0. ... 0. 0. -1.]]\n",
	"constant: 81\n"
	],
	"name": "stdout"
	}
	]
	},
	{
	"cell_type": "code",
	"metadata": {
	"id": "Kj4-1mli9gxC",
	"colab_type": "code",
	"colab": {}
	},
	"source": [
	""
	],
	"execution_count": 0,
	"outputs": []
	},
	{
	"cell_type": "code",
	"metadata": {
	"id": "xpEEdqA29gxF",
	"colab_type": "code",
	"colab": {}
	},
	"source": [
	""
	],
	"execution_count": 0,
	"outputs": []
	}
	]
	}