Skip to content

Instantly share code, notes, and snippets.

@mohzeki222

mohzeki222/chapter2.ipynb

Created May 19, 2021 09:48
Show Gist options
  • Save mohzeki222/45b897e59eb84a6ebb71ddbeb7237722 to your computer and use it in GitHub Desktop.
Save mohzeki222/45b897e59eb84a6ebb71ddbeb7237722 to your computer and use it in GitHub Desktop.
Download ZIP
chapter2.ipynb
Display the source blob
Display the rendered blob
Raw
chapter2.ipynb
{
"nbformat": 4,
"nbformat_minor": 0,
"metadata": {
"colab": {
"name": "chapter2.ipynb",
"provenance": [],
"collapsed_sections": [],
"include_colab_link": true
},
"kernelspec": {
"name": "python3",
"display_name": "Python 3"
},
"language_info": {
"name": "python"
}
},
"cells": [
{
"cell_type": "markdown",
"metadata": {
"id": "view-in-github",
"colab_type": "text"
},
"source": [
"<a href=\"https://colab.research.google.com/gist/mohzeki222/45b897e59eb84a6ebb71ddbeb7237722/chapter2.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "G1hCcMPMBmei"
},
"source": [
"## **はじめに**\n",
"\n",
"---\n",
"\n",
"まずはいつも通りGoogle colaboratory上にD-Wave Systemsが提供しているOcean SDKをインストールするところから始めましょう。\n"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "FQb5Otb7BviP",
"outputId": "0791f546-7ea0-424b-91e6-91e70c1ca59b"
},
"source": [
"pip install dwave-ocean-sdk"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"Requirement already satisfied: dwave-ocean-sdk in /usr/local/lib/python3.7/dist-packages (3.3.0)\n",
"Requirement already satisfied: dimod==0.9.13 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.9.13)\n",
"Requirement already satisfied: dwave-inspector==0.2.5 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.2.5)\n",
"Requirement already satisfied: dwave-system==1.4.0 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (1.4.0)\n",
"Requirement already satisfied: dwave-hybrid==0.6.1 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.6.1)\n",
"Requirement already satisfied: penaltymodel-mip==0.2.4; platform_machine == \"x86_64\" or platform_machine == \"amd64\" or platform_machine == \"AMD64\" in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.2.4)\n",
"Requirement already satisfied: dwave-qbsolv==0.3.2 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.3.2)\n",
"Requirement already satisfied: dwave-neal==0.5.7 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.5.7)\n",
"Requirement already satisfied: dwave-tabu==0.3.1 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.3.1)\n",
"Requirement already satisfied: penaltymodel-cache==0.4.3 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.4.3)\n",
"Requirement already satisfied: penaltymodel==0.16.4 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.16.4)\n",
"Requirement already satisfied: dwave-cloud-client==0.8.4 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.8.4)\n",
"Requirement already satisfied: dwave-greedy==0.1.2 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.1.2)\n",
"Requirement already satisfied: dwave-networkx==0.8.8 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.8.8)\n",
"Requirement already satisfied: minorminer==0.2.5 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.2.5)\n",
"Requirement already satisfied: penaltymodel-lp==0.1.4 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.1.4)\n",
"Requirement already satisfied: pyqubo==1.0.10 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (1.0.10)\n",
"Requirement already satisfied: dwavebinarycsp==0.1.2 in /usr/local/lib/python3.7/dist-packages (from dwave-ocean-sdk) (0.1.2)\n",
"Requirement already satisfied: numpy<2.0.0,>=1.17.3 in /usr/local/lib/python3.7/dist-packages (from dimod==0.9.13->dwave-ocean-sdk) (1.19.5)\n",
"Requirement already satisfied: importlib-resources>=3.2.0; python_version < \"3.9\" in /usr/local/lib/python3.7/dist-packages (from dwave-inspector==0.2.5->dwave-ocean-sdk) (5.1.2)\n",
"Requirement already satisfied: Flask>=1.1.1 in /usr/local/lib/python3.7/dist-packages (from dwave-inspector==0.2.5->dwave-ocean-sdk) (1.1.2)\n",
"Requirement already satisfied: homebase<2.0.0,>=1.0.0 in /usr/local/lib/python3.7/dist-packages (from dwave-system==1.4.0->dwave-ocean-sdk) (1.0.1)\n",
"Requirement already satisfied: networkx<3.0,>=2.0 in /usr/local/lib/python3.7/dist-packages (from dwave-system==1.4.0->dwave-ocean-sdk) (2.5.1)\n",
"Requirement already satisfied: click>5 in /usr/local/lib/python3.7/dist-packages (from dwave-hybrid==0.6.1->dwave-ocean-sdk) (8.0.0)\n",
"Requirement already satisfied: plucky>=0.4.3 in /usr/local/lib/python3.7/dist-packages (from dwave-hybrid==0.6.1->dwave-ocean-sdk) (0.4.3)\n",
"Requirement already satisfied: ortools<9.0.0,>=6.6.4659 in /usr/local/lib/python3.7/dist-packages (from penaltymodel-mip==0.2.4; platform_machine == \"x86_64\" or platform_machine == \"amd64\" or platform_machine == \"AMD64\"->dwave-ocean-sdk) (8.2.8710)\n",
"Requirement already satisfied: six<2.0.0,>=1.11.0 in /usr/local/lib/python3.7/dist-packages (from penaltymodel-cache==0.4.3->dwave-ocean-sdk) (1.15.0)\n",
"Requirement already satisfied: requests[socks]>=2.18 in /usr/local/lib/python3.7/dist-packages (from dwave-cloud-client==0.8.4->dwave-ocean-sdk) (2.23.0)\n",
"Requirement already satisfied: python-dateutil>=2.7 in /usr/local/lib/python3.7/dist-packages (from dwave-cloud-client==0.8.4->dwave-ocean-sdk) (2.8.1)\n",
"Requirement already satisfied: decorator<5.0.0,>=4.1.0 in /usr/local/lib/python3.7/dist-packages (from dwave-networkx==0.8.8->dwave-ocean-sdk) (4.4.2)\n",
"Requirement already satisfied: fasteners in /usr/local/lib/python3.7/dist-packages (from minorminer==0.2.5->dwave-ocean-sdk) (0.16)\n",
"Requirement already satisfied: scipy in /usr/local/lib/python3.7/dist-packages (from minorminer==0.2.5->dwave-ocean-sdk) (1.4.1)\n",
"Requirement already satisfied: Deprecated>=1.2.10 in /usr/local/lib/python3.7/dist-packages (from pyqubo==1.0.10->dwave-ocean-sdk) (1.2.12)\n",
"Requirement already satisfied: zipp>=0.4; python_version < \"3.8\" in /usr/local/lib/python3.7/dist-packages (from importlib-resources>=3.2.0; python_version < \"3.9\"->dwave-inspector==0.2.5->dwave-ocean-sdk) (3.4.1)\n",
"Requirement already satisfied: Jinja2>=2.10.1 in /usr/local/lib/python3.7/dist-packages (from Flask>=1.1.1->dwave-inspector==0.2.5->dwave-ocean-sdk) (2.11.3)\n",
"Requirement already satisfied: itsdangerous>=0.24 in /usr/local/lib/python3.7/dist-packages (from Flask>=1.1.1->dwave-inspector==0.2.5->dwave-ocean-sdk) (2.0.0)\n",
"Requirement already satisfied: Werkzeug>=0.15 in /usr/local/lib/python3.7/dist-packages (from Flask>=1.1.1->dwave-inspector==0.2.5->dwave-ocean-sdk) (2.0.0)\n",
"Requirement already satisfied: protobuf>=3.14.0 in /usr/local/lib/python3.7/dist-packages (from ortools<9.0.0,>=6.6.4659->penaltymodel-mip==0.2.4; platform_machine == \"x86_64\" or platform_machine == \"amd64\" or platform_machine == \"AMD64\"->dwave-ocean-sdk) (3.17.0)\n",
"Requirement already satisfied: absl-py>=0.11 in /usr/local/lib/python3.7/dist-packages (from ortools<9.0.0,>=6.6.4659->penaltymodel-mip==0.2.4; platform_machine == \"x86_64\" or platform_machine == \"amd64\" or platform_machine == \"AMD64\"->dwave-ocean-sdk) (0.12.0)\n",
"Requirement already satisfied: idna<3,>=2.5 in /usr/local/lib/python3.7/dist-packages (from requests[socks]>=2.18->dwave-cloud-client==0.8.4->dwave-ocean-sdk) (2.10)\n",
"Requirement already satisfied: urllib3!=1.25.0,!=1.25.1,<1.26,>=1.21.1 in /usr/local/lib/python3.7/dist-packages (from requests[socks]>=2.18->dwave-cloud-client==0.8.4->dwave-ocean-sdk) (1.24.3)\n",
"Requirement already satisfied: certifi>=2017.4.17 in /usr/local/lib/python3.7/dist-packages (from requests[socks]>=2.18->dwave-cloud-client==0.8.4->dwave-ocean-sdk) (2020.12.5)\n",
"Requirement already satisfied: chardet<4,>=3.0.2 in /usr/local/lib/python3.7/dist-packages (from requests[socks]>=2.18->dwave-cloud-client==0.8.4->dwave-ocean-sdk) (3.0.4)\n",
"Requirement already satisfied: PySocks!=1.5.7,>=1.5.6; extra == \"socks\" in /usr/local/lib/python3.7/dist-packages (from requests[socks]>=2.18->dwave-cloud-client==0.8.4->dwave-ocean-sdk) (1.7.1)\n",
"Requirement already satisfied: wrapt<2,>=1.10 in /usr/local/lib/python3.7/dist-packages (from Deprecated>=1.2.10->pyqubo==1.0.10->dwave-ocean-sdk) (1.12.1)\n",
"Requirement already satisfied: MarkupSafe>=0.23 in /usr/local/lib/python3.7/dist-packages (from Jinja2>=2.10.1->Flask>=1.1.1->dwave-inspector==0.2.5->dwave-ocean-sdk) (2.0.0)\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "0rgQk54CJvxX"
},
"source": [
"またシミュレーターとしてOpenJijを利用する場合にはこちらも用意しておきましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "ZNDyMuUqJu31",
"outputId": "1bca201f-6cae-41d8-89fc-dcd1fbc9323f"
},
"source": [
"pip install openjij"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"Requirement already satisfied: openjij in /usr/local/lib/python3.7/dist-packages (0.3.3)\n",
"Requirement already satisfied: scipy in /usr/local/lib/python3.7/dist-packages (from openjij) (1.4.1)\n",
"Requirement already satisfied: requests in /usr/local/lib/python3.7/dist-packages (from openjij) (2.23.0)\n",
"Requirement already satisfied: dimod>=0.9.1 in /usr/local/lib/python3.7/dist-packages (from openjij) (0.9.13)\n",
"Requirement already satisfied: jij-cimod>=1.2.0 in /usr/local/lib/python3.7/dist-packages (from openjij) (1.2.1)\n",
"Requirement already satisfied: numpy>=1.18.4 in /usr/local/lib/python3.7/dist-packages (from openjij) (1.19.5)\n",
"Requirement already satisfied: certifi>=2017.4.17 in /usr/local/lib/python3.7/dist-packages (from requests->openjij) (2020.12.5)\n",
"Requirement already satisfied: urllib3!=1.25.0,!=1.25.1,<1.26,>=1.21.1 in /usr/local/lib/python3.7/dist-packages (from requests->openjij) (1.24.3)\n",
"Requirement already satisfied: chardet<4,>=3.0.2 in /usr/local/lib/python3.7/dist-packages (from requests->openjij) (3.0.4)\n",
"Requirement already satisfied: idna<3,>=2.5 in /usr/local/lib/python3.7/dist-packages (from requests->openjij) (2.10)\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "-TN6A5YvBmjX"
},
"source": [
"そしてご自身がお持ちのAPI tokenと接続先のendpointを設定します。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "iH93qZmAB9tE"
},
"source": [
"token = '**'\n",
"endpoint = 'https://cloud.dwavesys.com/sapi/'"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "3FeUTX8BBmlm"
},
"source": [
"さて今回は様々な組合せ最適化問題に対応するQUBO行列を作ってみましょう！\n",
"正直に言って量子アニーリングマシンはQUBO行列さえ用意すれば様々な問題を解くことができますよ！と言っても、どうやってQUBO行列にすれば良いのか、とかどういうふうに考えたら良いのか、とか。\n",
"組合せ最適化問題について勉強してからかな？とか、かなりハードル高いです。\n",
"\n",
"そこで今回は量子アニーリングマシンの定式化でよくやる基本作法を学び、それらを組み合わせて色々な問題に対応させていくということを紹介します。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "z93p0hyYDUhg"
},
"source": [
"### 量子アニーリングマシンで解ける組合せ最適化問題"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "dp7bpZhwBmoB"
},
"source": [
"量子アニーリングマシンで考える組合せ最適化問題は、マシンの性質によってQUBO行列を利用した二次関数のみ解くことができます。\n",
"\\begin{equation}\n",
"E({\\bf x}) = \\sum_{i=1}^N \\sum_{j=1}^N Q_{ij} x_i x_j\n",
"\\end{equation}"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "sqOjRswAD88U"
},
"source": [
"量子アニーリングの解説やWebの情報等を調べると必ずそうした記述があります。そのせいでQUBO行列に意識が行きがちです。それよりも**大事なことは、二次関数である**ということです。\n",
"二次関数と聞けば、僕にでも私にでもできるのではないか？と考えることができるのではないでしょうか。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "oNbhr8kmDhEj"
},
"source": [
"そのため最初は**できるだけ二乗を作る意識**で、解くべき組合せ最適化問題を考えましょう。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "sXnDEfHmFah0"
},
"source": [
"二次関数で教わったことといえば最大値や最小値は頂点のところにあるよ、ということを思い出すことができるのではないでしょうか。\n",
"最小化問題を考えるという場合であれば、下に凸の次のような二次関数であれば簡単そうです。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "vDQMC2uKFLxn"
},
"source": [
"import numpy as np\n",
"a = 0.5\n",
"x = np.linspace(0,10,11)\n",
"y = a*x**2"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "XAvXEjB2HKBO"
},
"source": [
""
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "vT8qWKYfFuta"
},
"source": [
"ここでnp.linspace()を用いて、$0$から始まり、$10$でおわる$11$個の数列を用意しました。\n",
"つまり$x=(0,1,2,3,,,,10)$です。\n",
"その数列をxに格納しています。\n",
"それに対してyはaをかけて、xを二乗していますね。aは$0.5$としていますが、他の値にしたりして形を変えても結構ですよ。\n",
"このxとyによるグラフを見てみましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/",
"height": 265
},
"id": "uffGaFOwGRlq",
"outputId": "a5f6bbb3-62ee-4eac-ac15-6bdfac5bc696"
},
"source": [
"import matplotlib.pyplot as plt\n",
"plt.plot(x,y)\n",
"plt.show()"
],
"execution_count": null,
"outputs": [
{
"output_type": "display_data",
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAD4CAYAAAD1jb0+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAeLElEQVR4nO3deXxU9b3/8dcnIRC2hJ1AWAKCYAiyhUWo1r1YcasVlwqoINSlrvdWbXur1y5Xe1uX1v6uslU2WepS0OuG20V/F4GENQEFBAIJCQlbAoGs8/39kZEfUpaQzMyZybyfj4ePzJxMct6jyduT8/2e8zXnHCIiEnlivA4gIiJ1owIXEYlQKnARkQilAhcRiVAqcBGRCNUolDtr166dS0lJCeUuRUQiXmZm5l7nXPsTt4e0wFNSUsjIyAjlLkVEIp6Z5Zxsu06hiIhEKBW4iEiEUoGLiEQoFbiISIRSgYuIRKhazUIxsx3AIaAaqHLOpZtZG2AhkALsAMY65w4EJ6aIiJzobI7AL3HODXTOpfufPw587JzrDXzsfy4iIiFSn1Mo1wGz/I9nAdfXP46ISMOy73A5T7+9kaMV1QH/3rUtcAd8aGaZZjbZv62jcy7f/7gA6HiyLzSzyWaWYWYZRUVF9YwrIhI5Kqt93DNvNfNW5LB9b2nAv39tr8T8nnMuz8w6AEvN7KvjP+mcc2Z20pUhnHNTgakA6enpWj1CRKLG029vZOX2/bxw80BSOycE/PvX6gjcOZfn/1gIvAUMA/aYWScA/8fCgKcTEYlQC1buZM6XOUy+qCfXD0oOyj7OWOBm1tzMWn77GLgSyAKWABP8L5sALA5KQhGRCJOZs59/W5zFhb3b8djovkHbT21OoXQE3jKzb1//mnPufTNbBSwys4lADjA2aClFRCJEfvFRpsxZTXKrprx062BiYyxo+zpjgTvntgEDTrJ9H3BZMEKJiESisspqpszJ5GhFFfPvHk5is7ig7i+kt5MVEWmonHP84s0NrM8tZuq4IfTu2DLo+9Sl9CIiATDji+28uSaPR644lyv7JYVknypwEZF6+nxLEb9/dxOj+yVx/yW9QrZfFbiISD3k7Cvl/tfW0LtDS/40dgAxQRy0PJEKXESkjg6XV3H37AzMYNr4dJo3Ce2wogpcRKQOfD7Ho4vWsrXwMC/dOphubZuFPIMKXESkDv7yyVY+yN7DL69O5Xu923mSQQUuInKWPsgu4PmPNnPj4C7cNSrFsxwqcBGRs7B5zyEeWbiWAV0S+d0NafivUveEClxEpJYOHqng7tkZNGvSiFfGpRMfF+tpHhW4iEgtVFX7+Nn8NeQfLOPl24eQlBjvdSRdSi8iUhvPvv8Vn2/Zy7M39mdI99ZexwF0BC4ickZvrcll2ufbmXBBd24e2s3rOMeowEVETmN97kEee2MDI3q24VdjUr2O8x0qcBGRUyg6VM6UOZm0b9GEv942mLjY8KpMnQMXETmJiiof98zN5MCRCt64ZyRtWzTxOtI/UYGLiJzEk0uyycg5wEu3DaJf50Sv45xUeP09ICISBuZ+mcP8lTu59+JzGHN+Z6/jnJIKXETkOCu27eOpJdlc2rcDj17Zx+s4p6UCFxHxyzt4lHvnraZb22a8cMvAoC5IHAgqcBER4GhFNZNnZ1BR5WPa+HQS4oO7IHEgaBBTRKKec47H3ljPxvwSZkxI55z2LbyOVCs6AheRqPfKsm0sWbebf7myD5f27eh1nFpTgYtIVPvs60Keff8rxpzfiXsvPsfrOGdFBS4iUWtb0WF+Nn8NfZMS+MOPz/f03t51oQIXkah0qKySyXMyiYuNYeq4ITRrHHlDgpGXWESknnw+x8ML17J9bylzJw6na5vQL0gcCDoCF5Go8/xHm/loUyFPXpPKBee09TpOnanARSSqvLshn798spWb07sybkR3r+PUiwpcRKLGpvwSHl20jsHdWvH09f0ibtDyRCpwEYkK+0trFiROaNqIl28fQpNG3i5IHAi1LnAzizWzNWb2jv95DzNbYWZbzWyhmTUOXkwRkborq6zmnrmZFB4q55Vx6XRI8H5B4kA4myPwB4FNxz1/FnjeOdcLOABMDGQwEZFAqPY5Hlm0lhXb9/OfPz6fgV1beR0pYGpV4GbWBbgamO5/bsClwOv+l8wCrg9GQBGRunLO8dSSbN7dUMCvrj6P6wYmex0poGp7BP4C8HPA53/eFjjonKvyP88FTvpvxswmm1mGmWUUFRXVK6yIyNn4yydbmfNlDlO+35NJF/b0Ok7AnbHAzWwMUOicy6zLDpxzU51z6c659Pbt29flW4iInLXXVuzkuaWbuXFwFx4f3dfrOEFRmysxRwHXmtkPgXggAXgRaGVmjfxH4V2AvODFFBGpvQ+yC/jVPzZwSZ/2PHNj/4ifLngqZzwCd8494Zzr4pxLAW4BPnHO/QT4FPix/2UTgMVBSykiUksrtu3jZ/PXMKBrK/76k8HExTbc2dL1eWePAY+Y2VZqzonPCEwkEZG6+aqghEmzM+jauikzJwyNyBtUnY2zenfOuc+Az/yPtwHDAh9JROTs7dp/hPEzVtK8cSNmTxxO6+YN/9KUhvu3hYhEjf2lFUyYuZKyympm3TWM5FZNvY4UEg377wsRafBKy6u489VV5B08ytxJw+mT1NLrSCGjI3ARiViV1T7umbeaDbkHeem2wQxNaeN1pJDSEbiIRCSfz/Hz19ezbHMRz97YnytSI2cx4kDREbiIRKRn3v+Kt9bk8a8/6MPNQ7t5HccTKnARiThTl33D1GXbuGNkSsStJB9IKnARiShvrs7l9+9+xdXnd+LXY1Ib7FWWtaECF5GI8enXhfz89fWMPKctz40dQExM9JY3qMBFJEKs2XmAe+eupk9SS14Z1zBW1KkvFbiIhL2thYe569VVdEhowqt3DqNlfJzXkcKCClxEwlpBcRkTZq4kNsaYfdcw2rds4nWksKF54CIStoqPVDJh5kqKj1ayYPIIurdt7nWksKIjcBEJS2WV1dw9O4Pte0uZOm4IacmJXkcKOzoCF5GwU1Xt44H5a1iVs5+/3DqIkb3aeR0pLOkIXETCinOOf1ucxYcb9/DkmFTGnN/Z60hhSwUuImHl+aWbmb9yF/df0os7RvXwOk5YU4GLSNiYvXwHf/5kKzend+XRK8/1Ok7YU4GLSFj47/X5PLkkm8vP68jvbkiL6kvka0sFLiKe+9+te3l44VrSu7fmpdsG0agBL0QcSPq3JCKeysorZvKcTFLaNWP6+KHEx+kS+dpSgYuIZ3buO8Idf1tFQnwjZt01jMRmukT+bKjARcQTRYfKGTdzBVU+H7MnDqNTYnQsRBxIKnARCbnD5VXc+epK9pSUMfOOofTqED0LEQeSrsQUkZAqr6pmypwMNuUfYvr4dAZ3a+11pIilI3ARCRmfz/HoonX83637+MON53NJ3w5eR4poKnARCYlqn+NfXl/HO+vzefyqvtw4pIvXkSKeTqGISNBVVvt4ZNE63l63m0evOJeffj96FyIOJBW4iARVRVXNnQXfzy7giav6MkXlHTAqcBEJmvKqau6bt5qPNhXy6zGp3PU93ZwqkFTgIhIUZZXVTJmTyf9sLuI316cxbkR3ryM1OGccxDSzeDNbaWbrzCzbzP7dv72Hma0ws61mttDMGgc/rohEgiMVVUyctYplW4p49sb+Ku8gqc0slHLgUufcAGAgMNrMRgDPAs8753oBB4CJwYspIpHicHkVd/xtFcu/2cefbhrAzUO7eR2pwTpjgbsah/1P4/z/OOBS4HX/9lnA9UFJKCIRo6SskvEzVpCZc4AXbhnEjwZrqmAw1WoeuJnFmtlaoBBYCnwDHHTOVflfkgskn+JrJ5tZhpllFBUVBSKziISh4iOVjJu+gvW5xbx06yCuHaCl0IKtVgXunKt2zg0EugDDgL613YFzbqpzLt05l96+ffs6xhSRcLa/tILbpn/JpvxDvHz7EK7q38nrSFHhrGahOOcOmtmnwAVAKzNr5D8K7wLkBSOgiIS3vYfLuX36CrbvLWXq+CFc3EeXx4dKbWahtDezVv7HTYErgE3Ap8CP/S+bACwOVkgRCU+FJWXcMvVLduwrZeYdQ1XeIVabI/BOwCwzi6Wm8Bc5594xs43AAjP7LbAGmBHEnCISZvKLj3LbtBXsKSnj1TuHMaJnW68jRZ0zFrhzbj0w6CTbt1FzPlxEosyu/Ue4bfqXHCytZM7EYQzp3sbrSFFJV2KKyFnJ2VfKbdNWcKiskrmThjOgayuvI0UtFbiI1No3RYf5ybQVlFdV89rdI0hLTvQ6UlRTgYtIrWzZc4hbp63AOcf8ySPom5TgdaSopwIXkTPalF/C7dNXEBNjzL97BL07ag3LcKAVeUTktLLyirl12pfExcawcLLKO5zoCFxETmntroOMn7GClvFxzL97BN3aNvM6khxHBS4iJ5WxYz93/G0VbZo35rW7h9Oltco73OgUioj8k+Xf7GP8zJV0aNmEhVNGqLzDlI7AReQ7vtiyl0mzV9G1dTPmTRpOh4R4ryPJKajAReSYT78uZMqcTHq2a87cScNp16KJ15HkNFTgIgLA0o17uG/eanp3bMHcicNp3VyrJIY7FbiI8O6GfB6Yv4Z+yYnMvnMYic3ivI4ktaACF4lyi9fm8ciidQzs2opX7xxKy3iVd6TQLBSRKPZ6Zi4PL1xLevfWzL5rmMo7wugIXCRKLVi5kyfe2sCoc9oxbXw6TRvHeh1JzpIKXCQKzfxiO0+/s5GL+7Tn5duHEB+n8o5EKnCRKFJV7ePpdzYye3kOP+jXkT/fOogmjVTekUoFLhIlSsoquW/eaj7fspcpF/Xk56P7EhtjXseSelCBi0SBnfuOcNesVezYW8qzN/bn5qHdvI4kAaACF2ngVu3Yz5Q5mVT7HLMnDmPkOe28jiQBogIXacDeXJ3L429sILl1U2ZMSKdn+xZeR5IAUoGLNEA+n+O5pZt56dOtXNCzLf91+2BaNdOl8Q2NClykgTlaUc2jf1/LuxsKuDm9K7+5Po3GjXTNXkOkAhdpQApLypg0O4MNecX88ofnMenCHphppklDpQIXaSCydxczaVYGB49U8srtQ7iyX5LXkSTIVOAiDcDSjXt4cMEaEpvG8fefXkBacqLXkSQEVOAiEcw5x/TPt/P79zbRPzmR6ePTtYJOFFGBi0Soiiofv16cxYJVu/hh/yT+dNNA3ZAqyqjARSLQwSMV3DN3Ncu37eP+S3rxyBXnEqPL4qOOClwkwmzfW8rEV1eRe+Aoz40dwI8Gd/E6knhEBS4SQZZ/s4+fzs0kNsaYd/dwhqa08TqSeOiMs/vNrKuZfWpmG80s28we9G9vY2ZLzWyL/2Pr4McViV4LV+1k3IwVtG/ZhH/cO0rlLbVaUq0KeNQ5lwqMAO4zs1TgceBj51xv4GP/cxEJsGqf4z/e3cRjb2zggnPa8sY9I+nWtpnXsSQMnPEUinMuH8j3Pz5kZpuAZOA64GL/y2YBnwGPBSWlSJQqLa/ioYVrWbpxD7eP6MZT1/SjUawui5caZ3UO3MxSgEHACqCjv9wBCoCOp/iaycBkgG7ddA9ikdrKLz7KxFcz+KqghKeuSWXCyBRdFi/fUesCN7MWwBvAQ865kuN/kJxzzszcyb7OOTcVmAqQnp5+0teIyHetzz3IpFkZHKmoZsYdQ7mkTwevI0kYqtXfYmYWR015z3POvenfvMfMOvk/3wkoDE5Ekejy7oZ8xr6ynLjYGN64Z6TKW06pNrNQDJgBbHLOPXfcp5YAE/yPJwCLAx9PJHo45/jrp1u5d95qzuuUwD/uG0WfpJZex5IwVptTKKOAccAGM1vr3/YL4BlgkZlNBHKAscGJKNLwlVdV88SbG3hzdR7XDujMH358PvFxuixeTq82s1C+AE41cnJZYOOIRJ/9pRVMmZPBqh0HePjyc3ngsl4arJRa0ZWYIh7aWniIu17NoKCkjD/fOohrB3T2OpJEEBW4iEeWrNvNL9/cQJO4GBZMHsHgbrqYWc6OClwkxA6VVfLkkmzeXJ3HoG6t+Mutg+jSWldWytlTgYuE0OqdB3howVpyDxzhgct688ClvXRlpdSZClwkBKp9NVMEX/x4C0kJ8SyccoFuRiX1pgIXCbLcA0d4eOFaVu04wLUDOvOb69NIbBrndSxpAFTgIkG0ZN1ufvnWBpyD528ewPUDkzVFUAJGBS4SBIfLq/j14qxjA5Uv3jxIt4CVgFOBiwSYBiolVFTgIgGigUoJNRW4SABooFK8oAIXqafjByqfGzuAGwZpoFJCQwUuUkcaqBSvqcBF6mDNzgM8qIFK8ZgKXOQsVPsc/+fTrbyggUoJAypwkVrSQKWEGxW4SC1ooFLCkQpc5DQ0UCnhTAUucgrfGai8tBc/u6w3cRqolDCiAhc5gQYqJVKowEWOk3fwKA8vWMvKHfu5ZkBnfquBSgljKnARv7fX7eYXGqiUCKICl6h3qKySp5Zs5I3VuRqolIiiApeo5fM53lqTxzPvf8W+w+UaqJSIowKXqLRu10GeXJLN2l0HGdi1FdPGpzOwayuvY4mcFRW4RJWiQ+X85wdfsSgjl3YtmvDHmwbwo0HJxMToXLdEHhW4RIWKKh+zl+/gxY+2UFZVzZSLenL/pb1oGa8ZJhK5VODS4P3P5iKefjubb4pKubhPe349JpWe7Vt4HUuk3lTg0mDl7CvlN+9s4qNNe0hp24yZd6Rzad+OXscSCRgVuDQ4peVV/PXTrUz/fDtxscbjV/XlzlEpNGkU63U0kYA6Y4Gb2UxgDFDonEvzb2sDLARSgB3AWOfcgeDFFDkz5xyL1+7mP97bxJ6Scn40KJnHrupLx4R4r6OJBEVtJry+Cow+YdvjwMfOud7Ax/7nIp7JyivmppeX89DCtXRoGc8b94zkuZsHqrylQTvjEbhzbpmZpZyw+TrgYv/jWcBnwGMBzCVSK/sOl/PHD79mwapdtGnWmGdv7M9NQ7pqWqBEhbqeA+/onMv3Py4ATjkyZGaTgckA3bp1q+PuRL6rstrHnOU5PP/RZo5WVHPXqB48cFlv3XhKokq9BzGdc87M3Gk+PxWYCpCenn7K14nU1hdb9vLvb2ezpfAwF/Zux6/HpNK7Y0uvY4mEXF0LfI+ZdXLO5ZtZJ6AwkKFETmbX/iP89r838kH2Hrq2acrUcUO4IrWj7hgoUauuBb4EmAA84/+4OGCJRE5wpKKKlz/7hpeXbSPWjH/9QR8mfq8H8XGaFijRrTbTCOdTM2DZzsxygSepKe5FZjYRyAHGBjOkRCfnHO+sz+f3724iv7iMawd05okf9qVTYlOvo4mEhdrMQrn1FJ+6LMBZRI7ZuLuEp97OZuX2/ZzXKYEXbxnEsB5a1kzkeLoSU8LKgdIK/rT0a15bsZPEpnH87oY0bhnajVhNCxT5JypwCQtHK6pZsGonL3y0hcPlVYy/IIWHLu9Nq2aNvY4mErZU4OKpfYfLmb08hzlf5rC/tIILerblyWtT6ZuU4HU0kbCnAhdPbN9byvTPt/F6Zi7lVT4uP68Dd1/Yk2E92mhaoEgtqcAlpDJz9jN12TY+3LiHuJgYfjQ4mUkX9qBXB12II3K2VOASdNU+x9KNe5i67BtW7zxIYtM47ru4F+NHdqdDS91sSqSuVOASNGWV1byemcuML7azfW8pXVo35alrUhk7tCvNGutHT6S+9FskAXfiwOT5XRJ56bZBjO6XRKPY2tzBWERqQwUuAXPiwORlfTsw+SINTIoEiwpc6k0DkyLeUIFLnXw7MDnt821k5hzQwKSIB1TgclY0MCkSPvQbJ7Wy73A5c77MYfZyDUyKhAsVuJzW9r2lzPhiG3/P0MCkSLhRgctJZeYcYNqybXywsYC4mBhuGJTM3RdpYFIknKjA5Zj9pRV8mF3A3zNzNTApEgFU4FGusKSMD7ILeC+rgC+37cPnIKVtM566JpWb0rvSvIl+RETClX47o1DugSO8n1XA+1kFZO48gHNwTvvm3HtxL0anJdGvc4LOb4tEABV4lNi+t5T3svJ5P6uA9bnFAJzXKYGHLz+Xq9KS6N1R57ZFIo0KvIFyzrGl8DDvbSjgvax8vio4BMCALok8NrovV6UlkdKuuccpRaQ+VOANiHOO7N0lvJeVz3tZBWwrKsUM0ru35t/GpDI6LYnkVlrRXaShUIFHOJ/PsTb3IO9tyOf97AJ27T9KbIwxvEcb7hyZwg/6JdEhQTNIRBoiFXgEqvY5Vu3Yf2wgsqCkjLhYY1Svdtx/SS+uSE2iTXMtBizS0KnAI0RltY/l3+zjvawClm4sYO/hCpo0iuH757bnsf59uLRvRxKbxnkdU0RCSAUexsoqq/liy17eyyrgo017KD5aSbPGsVzStwNXpSVxSZ8OmqctEsX02x9GCg+VkZ1XQvbuYtbnFvO/3+zjcHkVLeMbccV5HRmdlsRF57YnPi7W66giEgZU4B5wzrG7uIysvGKy84rJ2l1CVl4xhYfKj72mR7vmXN2/E1f1T2LkOe1o3Eh3/BOR71KBB5nP58jZf4SsvGKydhez0V/WB45UAhBj0KtDC77Xqx39khNJ65xAaucEWsbrfLaInJ4KPICqqn1s21taU9Z5JccK+3B5FQBxsUafpJb8oF/SsbLum5RA08Y6JSIiZ08FXkflVdVs2XP42JF19u4SNuWXUFbpAyA+LobUTgncMCiZtOQE+nVO5NyOLXUqREQCRgVeC0crqtlUUFJzvtp/ZL15zyEqqx0ALZs0IrVzAj8Z3p205ATSOifSo11zrVQjIkFVrwI3s9HAi0AsMN0590xAUgWJc46ySh8lZZUUH62k5GjNx///uOrY577dtvdwOdv3luKr6WpaN4sjLTmRSRf2JK1zIv06J9CtTTNiYnT3PhEJrToXuJnFAn8FrgBygVVmtsQ5tzFQ4U7G53McKq86Vr7fKeFjxVx10m0lRyupqPad9vs3bxxLYtM4Evz/9OrQgqvP70xa5wTSkhPplBivW62KSFiozxH4MGCrc24bgJktAK4DAl7gv3hrA8s2F1FytJJD5VU4d+rXxhjHCjixaRwJ8XF0TmzqL+RGJB63/djjY9sa6bSHiESM+hR4MrDruOe5wPATX2Rmk4HJAN26davbjlo1Jb176+8UbsKJRdyspoBbNGmkI2QRiQpBH8R0zk0FpgKkp6ef5tj51O67pFdAM4mINAT1OV+QB3Q97nkX/zYREQmB+hT4KqC3mfUws8bALcCSwMQSEZEzqfMpFOdclZndD3xAzTTCmc657IAlExGR06rXOXDn3LvAuwHKIiIiZ0Fz5kREIpQKXEQkQqnARUQilApcRCRCmTvddemB3plZEZBTxy9vB+wNYJxIoPccHfSeG776vt/uzrn2J24MaYHXh5llOOfSvc4RSnrP0UHvueEL1vvVKRQRkQilAhcRiVCRVOBTvQ7gAb3n6KD33PAF5f1GzDlwERH5rkg6AhcRkeOowEVEIlREFLiZjTazr81sq5k97nWeYDKzrmb2qZltNLNsM3vQ60yhYmaxZrbGzN7xOksomFkrM3vdzL4ys01mdoHXmYLNzB72/1xnmdl8M4v3OlOgmdlMMys0s6zjtrUxs6VmtsX/sXUg9hX2BX7c4slXAanArWaW6m2qoKoCHnXOpQIjgPsa+Ps93oPAJq9DhNCLwPvOub7AABr4ezezZOABIN05l0bNbahv8TZVULwKjD5h2+PAx8653sDH/uf1FvYFznGLJzvnKoBvF09ukJxz+c651f7Hh6j5pU72NlXwmVkX4GpgutdZQsHMEoGLgBkAzrkK59xBb1OFRCOgqZk1ApoBuz3OE3DOuWXA/hM2XwfM8j+eBVwfiH1FQoGfbPHkBl9oAGaWAgwCVnibJCReAH4O+LwOEiI9gCLgb/7TRtPNrLnXoYLJOZcH/BHYCeQDxc65D71NFTIdnXP5/scFQMdAfNNIKPCoZGYtgDeAh5xzJV7nCSYzGwMUOucyvc4SQo2AwcB/OecGAaUE6M/qcOU/73sdNf/z6gw0N7PbvU0Veq5m7nZA5m9HQoFH3eLJZhZHTXnPc8696XWeEBgFXGtmO6g5RXapmc31NlLQ5QK5zrlv/7p6nZpCb8guB7Y754qcc5XAm8BIjzOFyh4z6wTg/1gYiG8aCQUeVYsnm5lRc150k3PuOa/zhIJz7gnnXBfnXAo1/30/cc416CMz51wBsMvM+vg3XQZs9DBSKOwERphZM//P+WU08IHb4ywBJvgfTwAWB+Kb1mtNzFCIwsWTRwHjgA1mtta/7Rf+9UelYfkZMM9/YLINuNPjPEHlnFthZq8Dq6mZbbWGBnhJvZnNBy4G2plZLvAk8AywyMwmUnNL7bEB2ZcupRcRiUyRcApFREROQgUuIhKhVOAiIhFKBS4iEqFU4CIiEUoFLiISoVTgIiIR6v8BVaINbGH9sxwAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"tags": [],
"needs_background": "light"
}
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "rj1Q4ba7GXuE"
},
"source": [
"おなじみの放物線ができました。片側だけなので寂しいですか？\n",
"それでは頂点を動かしてみましょう。\n"
]
},
{
"cell_type": "code",
"metadata": {
"id": "rpLrCG3aGhRI"
},
"source": [
"b = 3\n",
"y = a*(x-b)**2"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "NkHDingzGgmN"
},
"source": [
"ここでbの値によって二次関数が平行移動するようにします。b=3ですから右に3ズレるはずです。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/",
"height": 265
},
"id": "rRXMtj3iGrLW",
"outputId": "a133f5a7-e3cc-40da-f38a-ffe7d322d979"
},
"source": [
"plt.plot(x,y)\n",
"plt.show()"
],
"execution_count": null,
"outputs": [
{
"output_type": "display_data",
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAD4CAYAAAD1jb0+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3dd3xUZb7H8c+PAAmhJyS0JHTpPYCK67Kii2LBXQXBxlpW73q9u3t18Vquu+p6r671WlEUF2y4gtgrgq4FFRKkh16S0BIIhEBInef+kegGBQnJzJwp3/frlVcmZ8r5DuWbZ8488xxzziEiIuGngdcBRESkblTgIiJhSgUuIhKmVOAiImFKBS4iEqYaBnNnbdq0cZ07dw7mLkVEwl5mZuZu51zSD7cHtcA7d+5MRkZGMHcpIhL2zGzrkbbrEIqISJg6ZoGbWaqZfWJmq81slZn9oXr7HWa2zcyWVn+NDXxcERH5Tm0OoVQANzrnlphZcyDTzOZVX/ewc+6BwMUTEZGjOWaBO+d2ADuqLxeZWRbQMdDBRETkpx3XMXAz6wwMBr6p3nS9mS03s+fMrPVR7nONmWWYWUZ+fn69woqIyL/UusDNrBnwGvBH59x+YCrQDRhE1Qj9wSPdzzk3zTmX7pxLT0r60SwYERGpo1oVuJk1oqq8X3LOzQVwzu1yzlU653zAM8DwwMUUEZEfqs0sFAOmA1nOuYdqbG9f42a/Alb6P56ISHgrKa/kjrdWsaPwkN8fuzazUEYClwErzGxp9bZbgUlmNghwwBbgWr+nExEJc9M+28SMhVv4Zd+2tG/ZxK+PXZtZKF8AdoSr3vNrEhGRCLN93yGe/HQDZ/Vrx8nd2vj98fVJTBGRALn3/TX4HNw6tndAHl8FLiISAIu3FPDWsu1ce2pXUhPiA7IPFbiIiJ9V+hx3vLWKdi3i+N2obgHbjwpcRMTPZmfksGr7fm4Z24v4xoFb9FUFLiLiR4WHyrn/w7Wkd2rNeQM7BHRfQV0PXEQk0j02fz0FxWXMOHc4VR+jCRyNwEVE/GRD3gFmLNzChKGp9E9pGfD9qcBFRPzk7ndX06RRDH8a0zMo+1OBi4j4wYI1u/h0bT6/H92DpOaxQdmnClxEpJ7KKnz89Z0surZpyuSTOwdtvypwEZF6mrFwM5t3H+T2c/rQuGHwalUFLiJSD3lFJTw6fwO/6JnEL3olB3XfKnARkXp44MO1lJRXcvs5fYK+bxW4iEgdLc/dx+zMXK4Y2ZmuSc2Cvn8VuIhIHThXtd5JYtPG/MfoHp5kUIGLiNTBm0u3syR7H1PG9KRFXCNPMqjARUSO08HSCu55P4v+HVsyfmiqZzm0FoqIyHGa+ulGdu0v5clLhtCgQWDXO/kpGoGLiByH7D3FTPt8E+MGdWBopwRPs6jARUSOw/++l0WMGTef1cvrKCpwEZHaWrhhNx+s2sl1o7r5/QzzdaECFxGphYpKH3e+vZqU1k347aldvY4DqMBFRGrl5UXZrN1VxG1jexPXKMbrOIAKXETkmPYeLOPBj9ZxUtdEzuzXzus431OBi4gcw8Mfr6OopJy/nNcn4KdJOx4qcBGRn7Bm535e/Horl4zoRK92LbyOcxgVuIjIUTjnuOvt1TSPa8QNZ5zgdZwfUYGLiBzFh6t2snDjHm444wRaN23sdZwfUYGLiBxBSXkld7+bxQltm3HJiDSv4xyR1kIRETmCZz/fRO7eQ7x09QgaxoTmWDc0U4mIeGhnYQlPfLKRMX3bMrJ7G6/jHNUxC9zMUs3sEzNbbWarzOwP1dsTzGyema2v/t468HFFRALv3vezqHSO28YG/zRpx6M2I/AK4EbnXB/gRODfzawPcDMw3znXA5hf/bOISFjL3FrAG0u389ufdSEtMd7rOD/pmAXunNvhnFtSfbkIyAI6AuOAmdU3mwmcH6iQIiLB4PM57nx7NW1bxHLdqO5exzmm4zoGbmadgcHAN0Bb59yO6qt2Am39mkxEJMjmZOayPLeQm8/qRdPY0J/jUesCN7NmwGvAH51z+2te55xzgDvK/a4xswwzy8jPz69XWBGRQCkqKee+D9cwJK0V5w/q6HWcWqlVgZtZI6rK+yXn3NzqzbvMrH319e2BvCPd1zk3zTmX7pxLT0pK8kdmERG/e2zBBnYfKOMv5/YNqfVOfkptZqEYMB3Ics49VOOqt4DJ1ZcnA2/6P56ISOBtyj/A37/czPihKQxMbeV1nFqrzUGekcBlwAozW1q97VbgXuBVM7sK2ApMCExEEZHAuvvdLGIbxjDlzJ5eRzkuxyxw59wXwNFeT4z2bxwRkeD6ZG0eC9bkcctZvUhuHud1nOOiT2KKSNQqq/Dx13dW06VNU64Y2cXrOMdNBS4iUev5r7awKf8g/312bxo3DL86DL/EIiJ+sPtAKY98vJ6fn5DEab2SvY5TJypwEYlKD3y4lkPlldx+TmidJu14qMBFJOqs3FbIPzJymHxyZ7onN/M6Tp2pwEUkqjjnuOOtVSTEN+b3o3t4HadeVOAiElXeXr6DjK17+dOYnrRs0sjrOPWiAheRqFFcVsE972XRt0MLJqSneh2n3lTgIhI1nvp0IzsKS/jLuX2JaRCeb1zWpAIXkaiQU1DM059t4tyBHRjeJcHrOH6hAheRiOdc1YkazOCWs3p5HcdvVOAiEvHmZObycdYubjyjJx1aNfE6jt+owEUkouXuLebOt1czvEsCV54Sfuud/BQVuIhELJ/PMWX2cpxzPDh+YES8cVmTClxEItaMhVv4atMebj+nD6kJoX2G+bpQgYtIRNqQd4C/fbCG03olc9Gw8J/zfSQqcBGJOBWVPm58dSlNGsdw76/7h+1iVcdSm1OqiYiElSc/3ciy3EKeuHgIyS3C6yw7x0MjcBGJKCtyC3l0/nrOG9iBswe09zpOQKnARSRilJRXcsOrS0ls1pi7xvX1Ok7A6RCKiESMBz9ay/q8A8y8cjit4ht7HSfgNAIXkYjw9aY9PPvFZi4ZkcbPT0jyOk5QqMBFJOwdKK3gT7OXkZYQz61je3sdJ2h0CEVEwt7d76xm275DzL72JJrGRk+taQQuImFtwZpdvLI4h2tP7UZ658hYJra2VOAiErYKDpZx05wV9GrXnP88I7zPb1kX0fNaQ0QiinOO299YSeGhMp6/cjixDWO8jhR0GoGLSFh6a9l23l2xgz+efgJ9OrTwOo4nVOAiEnZ2FpZw+xsrGZzWimtP7ep1HM+owEUkrDjnuOm15ZRXOh6aMIiGMdFbY9H7zEUkLL30TTafrcvn1rG96NKmqddxPKUCF5GwsWX3Qf7n3Sx+1qMNl57Yyes4njtmgZvZc2aWZ2Yra2y7w8y2mdnS6q+xgY0pItGu0ue4cfYyGsYY9104IGLX+D4etRmBzwDOPML2h51zg6q/3vNvLBGRw037bBOZW/dy17i+tG8ZOWeWr49jFrhz7jOgIAhZRESOKGvHfh6et46z+rXj/EEdvY4TMupzDPx6M1tefYil9dFuZGbXmFmGmWXk5+fXY3ciEo3KKnzc8OoyWjRpxN3n99OhkxrqWuBTgW7AIGAH8ODRbuicm+acS3fOpSclRccSjyLiP4/MX0fWjv3c++v+JDaL9TpOSKlTgTvndjnnKp1zPuAZYLh/Y4mIQObWvUz9dCMT0lM4vU9br+OEnDoVuJnVPNHcr4CVR7utiEhdFJdVrfHdvmUTbj+nj9dxQtIxF7Mys1nAKKCNmeUCfwFGmdkgwAFbgGsDmFFEotC9769h8+6DzPrtiTSPa+R1nJB0zAJ3zk06wubpAcgiIgLA5+vzef6rrVw5sgsndUv0Ok7I0icxRSSkFBaXM2X2crolNeWmM3t6HSekqcBFJKTc8fYq8g+U8vBFg4hrFH1rfB8PFbiIhIz3V+zg9W+3cf0vujMgpZXXcUKeClxEQkJ+USm3vr6C/h1bcv1p3b2OExZU4CLiOecct8xdzsGySh6aMJBGUbzG9/HQn5KIeG52Zi4fZ+Vx05ie9Gjb3Os4YUMFLiKeyiko5q63VzOiSwJXjuzidZywogIXEc/4fI4pc5YB8MD4gTRooIWqjocKXEQ88/eFW/h6UwF/PqcPqQnxXscJOypwEfHEhrwi7vtgDaN7JTM+PcXrOGFJBS4iQVdeWbXGd3zjGO65oL/W+K6jY66FIiLib09+spHluYVMvWQIyc3jvI4TtjQCF5GgWp67j8cWrOf8QR04q3/7Y99BjkoFLiJBU1xWwQ2vLqNNs1juPK+f13HCng6hiEhQ+HyOP76ylE35B3j+yhG0jNca3/WlEbiIBMX9H63lo9W7uP2cPpzSo43XcSKCClxEAu61zFymfrqRi0ek8ZuTO3sdJ2KowEUkoBZvKeCWuSs4uVsid57XV1MG/UgFLiIBk1NQzLUvZNKxdROevGSIVhn0M/1pikhAFJWUc9XMxVT6HNMnp9MqvrHXkSKOClxE/K7S5/j9rG/ZlH+QqZcMoWtSM68jRSRNIxQRv/vf97L4ZG0+//OrfpzcXTNOAkUjcBHxq1mLspn+xWauGNmZS0Z08jpORFOBi4jfLNy4m9vfWMnPT0jitrG9vY4T8VTgIuIXm3cf5HcvLqFLm6Y8dvFgGmrGScDpT1hE6q2wuJyrZiwmpoExffIwWsTpY/LBoAIXkXopr/Rx3cuZ5Owt5qlLh5KWqDPrBItmoYhInTnnuPPtVXy5YQ/3XziA4V0SvI4UVTQCF5E6m7lwCy9+nc21P+/K+PRUr+NEHRW4iNTJp2vzuOud1ZzRpy3/NaaX13GikgpcRI7b+l1F/MfL39KzXQv+76JBNGigBaq8cMwCN7PnzCzPzFbW2JZgZvPMbH3199aBjSkioaLgYBlXzcwgtlEMz05Op2ms3krzSm1G4DOAM3+w7WZgvnOuBzC/+mcRiXBlFT7+7cVMdu4v4ZnLh9KxVROvI0W1Yxa4c+4zoOAHm8cBM6svzwTO93MuEQkxzjlue30FizYXcP+FAxicphfeXqvrMfC2zrkd1Zd3Am2PdkMzu8bMMswsIz8/v467ExGvPfP5JmZn5vL70T0YN6ij13EEP7yJ6ZxzgPuJ66c559Kdc+lJSUn13Z2IeGDe6l3c8/4azh7Qnj+O7uF1HKlW1wLfZWbtAaq/5/kvkoiEkqwd+/nDK9/Sv2NLHrhwoGachJC6FvhbwOTqy5OBN/0TR0RCSX5RKVfPzKBFXCOeuTydJo1jvI4kNdRmGuEs4Cugp5nlmtlVwL3AGWa2Hji9+mcRiSAl5ZVc80IGBQfLeHZyOm1bxHkdSX7gmBM4nXOTjnLVaD9nEZEQ4Zzjv15bzrfZ+3jq0iH069jS60hyBPokpoj8yOMLNvDm0u1MGdOTM/u19zqOHIUKXEQO896KHTw4bx2/HtyR60Z18zqO/AQVuIh8b3nuPm54dSlDO7Xmngv6Y6YZJ6FMBS4iAOwsLOG3z2eQ2DSWpy8bSmxDzTgJdSpwEeFQWSVXP7+YAyUVTP9NOm2axXodSWpBy4iJRDmfz3HDq0tZtX0/0yen06tdC68jSS1pBC4S5R6at473V+7ktrG9Oa3XUZc1khCkAheJYm98u43HP9nAxGGpXHVKF6/jyHFSgYtEqcyte7npteWc2DWBu8b104yTMKQCF4lCOQXFXPtCBh1axjH1kqE0bqgqCEd6E1MkymzKP8Bl0xdRWuHjlWuG0bppY68jSR2pwEWiyMpthUx+bhEAs357It2Tm3mcSOpDBS4SJRZtLuCqGYtpHteQF64eQbcklXe4U4GLRIEFa3bxuxeXkNK6CS9cNYIOOhlxRFCBi0S4N77dxp9mL6N3+xbMuGIYifqUZcRQgYtEsJkLt/CXt1ZxYtcEnrk8neZxjbyOJH6kAheJQM45Hp2/gYc/Xsfpvdvy+MWDiWukxakijQpcJML4fI673lnNjIVbuGBICn+7oD8NYzTPOxKpwEUiSHmlj/+as5y5327jypFd+O+ze+ss8hFMBS4SIUrKK7n+5SV8nJXHjWecwPWnddfH4yOcClwkAhSVlHP1zAy+2VzAXeP6cvlJnb2OJEGgAhcJc3sOlDL574tYs6OIRyYOYtygjl5HkiBRgYuEse37DnHp9G/YtvcQ0y4fqvW8o4wKXCRMbcw/wGXPfkNRSQUvXDWC4V0SvI4kQaYCFwlDK7cVcvlzi2hgMOuaE+nXsaXXkcQDKnCRMPP1pj1cPTODlk0a8eLVI+jSpqnXkcQjKnCRMDJv9S7+/eUlpCXE88JVw2nfUotSRTMVuEiYmLsklylzltOvQwv+fsVwEnQihqinAhcJA899sZm73lnNyd0SmXZ5Os1i9V9XVOAiIc05x8Mfr+fR+esZ07ctj0zUolTyLypwkRDl8znufHsVM7/ayvihKdzzay1KJYerV4Gb2RagCKgEKpxz6f4IJRLtyit9TJm9jDeWbue3P+vCrWN7a10T+RF/jMB/4Zzb7YfHERGqFqW67qUlLFiTx5QxPbluVDeVtxyRDqGIhJD9JeVcPSODxVsLuPv8flx6YievI0kIq+8BNQd8ZGaZZnaNPwKJRKvdB0qZNO1rlmTv5dGJg1Xeckz1HYGf4pzbZmbJwDwzW+Oc+6zmDaqL/RqAtLS0eu5OJDLl7i3m8umL2F54iGcnpzOqZ7LXkSQM1GsE7pzbVv09D3gdGH6E20xzzqU759KTkpLqszuRiLQhr4jxT31F/oFSXrxqhMpbaq3OBW5mTc2s+XeXgV8CK/0VrKY9B0rJ3VsciIcW8dTy3H2Mf+oryisd/7jmJNI7a0VBqb36jMDbAl+Y2TJgEfCuc+4D/8Q63P++t4YxD3/GC19vxedzgdiFSFBV+hzPfr6JCU9/RdPYhsz5t5Po06GF17EkzNT5GLhzbhMw0I9Zjuo/z+jBrv0l3P7GSt5bvoP7LhxAakJ8MHYt4ncb8w9w05zlZG7dy+heydxzQX+Sm8d5HUvCkDkXvBFtenq6y8jIqNN9nXP8Y3EOd7+bhc85bj6rF5eO6KQzbkvYqPQ5pn+xiQc/WkdcoxjuOK8P5w/qqDneckxmlnmkD0qGzTxwM2Pi8DR+dkISt8xdwZ/fXMW71aPxTolaD1lC24a8IqbMWc632fv4ZZ+23H1+P5JbaNQt9RM2I/CanHPMzszlr++sprzSx01jevGbkztrNC4hp6LSxzOfb+bhj9fRtHEMd47rx7kD2mvULccl7EfgNZkZE9JT+VmPNtw6dwV3vbOa91fu4L4LB+rsJBIy1u0qYsrsZSzLLeTMvu346/n9SGoe63UsiSBhOQKvyTnH3CXbuPPtVZRW+JgypidXjOxCjEbj4pGKSh9Pf7aJRz5eT7O4htw1ri9n99eoW+ouokbgNZkZFwxN4ZQebbjt9RXc/W4W763Ywf3jB9ItqZnX8STKrNm5nymzl7NiWyFnD2jPXef1JbGZRt0SGBGzuHDbFnE8c3k6/3fRIDbmH+SsRz7n6X9upFLzxiUIyit9PDp/Pec+9gU7Cg8x9ZIhPHHxEJW3BFTYj8BrMjPOH9yRk7sn8t+vr+Se99fw/sqdPDB+AN2Tm3sdTyLUqu2FTJm9nNU79nPewA7ccV5fna9SgiJiRuA1JTeP4+nLhvLopMFs3XOQsY9+wdRPN1JR6fM6mkSQsgofD89bx7jHvySvqJSnLq36N6fylmCJqBF4TWbGeQM7cFLXRP785kr+9sEaPqieqdKznUbjUj8rtxXyp9nLWLOziF8N7sifz+lDaxW3BFlEjsBrSmoey9RLh/LExUPI2XuIcx/7gscXrKdco3Gpg9KKSh78aC3jnviSgoNlPHt5Og9fNEjlLZ6I2BH4D509oD0ndk3gz2+t4oGP1vHBqp3cf+FAerfXAkJSO8tz9zFl9nLW7irigiEp/PmcPrSMb+R1LIliET8CrymxWSxPXDyEqZcMYWdhCec9/gWPztdoXH5aaUUl932whl89uZB9h8p47jfpPDhhoMpbPBc1I/CazurfnhFdE7nz7VU8NG8dH6zcyf3jB9C3Q0uvo0mIWZqzjymzl7E+7wAT0lO47ew+tGyi4pbQEFUj8JoSmjbmkYmDefqyoeQVlTLu8S95aN46yio0GpeqM8Pf834Wv37ySw6UVjDjimHcd+FAlbeElKgcgdc0pm87RnRJ4M63V/Po/PV8tGonD4wfSL+OGo1Hq8yte7lpzjI25h9k0vBUbhnbmxZxKm4JPVE7Aq+pVXxjHr5oEM9enk7BwTLGPfElD360ltKKSq+jSRCVlFfyP++u5sKnFlJS7uP5K4dzz68HqLwlZEX9CLym0/u0ZVjnBP767moeW7CBD6tH4wNSWnkdTQIsY0sBN81ZzqbdB7lkRBo3n9WL5ipuCXFhvxphoHyyJo9b5q4gr6iEMX3bMWl4Gqd0b6M1xyNIpc/xz3V5vPxNDvPX7KJjqyb87YIBjOzexutoIoc52mqEKvCfsL+knCcWbGB2Zi4FB8tITWjCxGFpjB+aorOphLHt+w7xakYOry7OYXthCW2axTJxWCq/G9WNprF6USqhRwVeD6UVlXy0ahezFmWzcOMeYhoYo3slM2lEGqf2SNLa42GgotLHJ2vzmbUom0/X5uGAn/VI4uLhqYzu3ZZGMXo7SEKXCtxPtuw+yCuLc5iTmcPuA2V0bNWEi4alMiE9lXYtNSoPNTkFxVWj7Ywcdu0vJbl5LBPSU7loWCqpCfFexxOpFRW4n5VV+Pg4q2pU/vn63TQwOK1XMpOGp/HzE5JoqBGdZ8orfczPymPWomw+W58PwKgTkpg0PI3TeiXr70bCjgo8gLL3FPOPjGxezcglv6iU9i3jmJCeyoRhqXRs1cTreFEje08xryzOZnZm1d9DuxZxVa+O9PcgYU4FHgTfjfxeWZzNP9flY8DPNfILKL0SkmigAg+y3L3FvLo4h39UH3tt26Lq2OuEdB179Yctuw8ya3E2r2XmHvZexPj0FNq31GhbIosK3CPfzX54ZVE2n2j2Q71oNpBEKxV4CPhu/vE/Fuewo3r+8fj0FCYOS6VTYlOv44WsjfkHeGVRNq8t2UbBwTJSWjdh0nDNx5fooQIPIZU+x2fr8nl5UTYL1uRR6XOc0r0Nk4ancUaftjRuqFF5SXklH67aycvfZPPN5gIaNjDO6NNWn4iVqKQCD1E7C0uYnZHDK4tz2LbvEIlNG3NhegoTh6XRpU30jcrX7ypi1qIc5n6by77icjolxjNxWBoXDk0hqXms1/FEPKECD3GVPsfn6/N5ZVEO87J2Uelz9GrXnE6J8aS2jictMZ7UhHjSEuJJad2E2IYxXkeus5LySnL3HiKnoJjsGl9bdh9kfd4BGsUYv+zbjouHp3FS10SNtiXqHa3AtfBDiIhpYIzqmcyonsnk7S9hzpJclmzdy+bdB/nnunxKyv91ogkzaNci7vtC/+7ru5/bNGuMmXel55wjv6j0sHLOLij+vrB37S897PZNGsWQlhBPp8SmjE9P4YIhKSQ202hb5FjqVeBmdibwCBADPOucu9cvqaJccos4rhvV/fufnXPkHyj914h1z6HvC/GL9bvZub/ksPt/V4ipCU2OWPJxjeo/ei8uqyCn4NBhxfz9973FP/qF0776F86pPZKqsiTGk9I6NH7hiISrOhe4mcUATwBnALnAYjN7yzm32l/hpIqZkdw8juTmcQztlPCj6492SCKnoJiFG/dQXHb4iSmSm8f+aNSeVn2oJrl5LA0aGD6fY1dRCdl7Di/nqq9D7D5w+Ci6WWxDUhPi6ZrUlFE9kw577I5hfshHJFTVZwQ+HNjgnNsEYGavAOMAFXiQxTWKoXtyM7onN/vRdc45Cg6W/egwRnZBMd9sLuD1pduo+TZI44YNSGoWS35RKWWV/xpFNzDo0KoJaQnxnN47mdSa5Z8QT+v4RhpFiwRZfQq8I5BT4+dcYMQPb2Rm1wDXAKSlpdVjd1IXZkZis1gSm8UyOK31j64vrahk+76Swwo+b38JbVvGHXbopUOrJvrQkUiICfibmM65acA0qJqFEuj9yfGJbRhDlzZNo3LKoki4q8+QahuQWuPnlOptIiISBPUp8MVADzPrYmaNgYnAW/6JJSIix1LnQyjOuQozux74kKpphM8551b5LZmIiPykeh0Dd869B7znpywiInIcNK1ARCRMqcBFRMKUClxEJEypwEVEwlRQl5M1s3xgax3v3gbY7cc44UDPOTroOUeH+jznTs65pB9uDGqB14eZZRxpPdxIpuccHfSco0MgnrMOoYiIhCkVuIhImAqnAp/mdQAP6DlHBz3n6OD35xw2x8BFRORw4TQCFxGRGlTgIiJhKiwK3MzONLO1ZrbBzG72Ok+gmVmqmX1iZqvNbJWZ/cHrTMFgZjFm9q2ZveN1lmAws1ZmNsfM1phZlpmd5HWmQDOz/6z+N73SzGaZWZzXmfzNzJ4zszwzW1ljW4KZzTOz9dXff3x6rDoI+QKvcfLks4A+wCQz6+NtqoCrAG50zvUBTgT+PQqeM8AfgCyvQwTRI8AHzrlewEAi/LmbWUfg90C6c64fVctQT/Q2VUDMAM78wbabgfnOuR7A/Oqf6y3kC5waJ092zpUB3508OWI553Y455ZUXy6i6j92R29TBZaZpQBnA896nSUYzKwlcCowHcA5V+ac2+dtqqBoCDQxs4ZAPLDd4zx+55z7DCj4weZxwMzqyzOB8/2xr3Ao8COdPDmiy6wmM+sMDAa+8TZJwP0fcBPg8zpIkHQB8oG/Vx82etbMIvrEpM65bcADQDawAyh0zn3kbaqgaeuc21F9eSfQ1h8PGg4FHrXMrBnwGvBH59x+r/MEipmdA+Q55zK9zhJEDYEhwFTn3GDgIH56WR2qqo/7jqPql1cHoKmZXeptquBzVXO3/TJ/OxwKPCpPnmxmjagq75ecc3O9zhNgI4HzzGwLVYfITjOzF72NFHC5QK5z7rtXVnOoKvRIdjqw2TmX75wrB+YCJ3ucKVh2mVl7gOrvef540HAo8Kg7ebKZGVXHRrOccw95nSfQnHO3OOdSnHOdqfr7XeCci+iRmXNuJ5BjZj2rN40GVnsYKRiygRPNLL763yKnzA0AAACYSURBVPhoIvyN2xreAiZXX54MvOmPB63XOTGDIUpPnjwSuAxYYWZLq7fdWn0OUokc/wG8VD0w2QRc4XGegHLOfWNmc4AlVM20+pYI/Ei9mc0CRgFtzCwX+AtwL/CqmV1F1ZLaE/yyL32UXkQkPIXDIRQRETkCFbiISJhSgYuIhCkVuIhImFKBi4iEKRW4iEiYUoGLiISp/wfzBDZsDZbaOQAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"tags": [],
"needs_background": "light"
}
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "FtZ9cvw8GsrK"
},
"source": [
"確かに右にずれて、下に凹んだ頂点が見えるようになりました。左右に伸びた放物線も見ると都ができました。\n",
"もっと鋭いものを描きたかったらaの値を変えてみましょう。\n",
"こうしたグラフを描くのもプログラミングの威力ですぐに描くことができます。\n",
"この場合は最小値が頂点のところの$3$にくるはずというのはすぐにわかります。\n",
"\n",
"量子アニーリングマシンでもそうした結果を出すことはできるでしょうか。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "cA7tJ5wFHy4t"
},
"source": [
"### 量子アニーリングマシンで整数値を扱う"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "6z1WQY_RH5qC"
},
"source": [
"ただ問題は量子アニーリングマシンでは$x_i=0$または$1$をとる二値の変数のみが利用できました。それでは答えが$3$のような整数になるということはありません。\n",
"よくよく考えてみると$0$と$1$しか扱えないというのはなかなか不自由です。\n",
"整数値くらい使えるようにならないものでしょうか。\n",
"そのために幾つかの変数を足し合わせて、その合計値で整数を示すようにしましょう。\n",
"例えば\n",
"\\begin{equation}\n",
"x = \\sum_{i=1}^N x_i\n",
"\\end{equation}\n",
"というように二値をとる$x_i$を足し算していけば、$x$は最小で$0$（全ての$x_i$が$0$をとる場合）となり、$x$は最大で$N$（全ての$x_i$が$1$をとる場合）となります。\n",
"これで最小$0$、最大$N$の整数を量子アニーリングマシンの中で動作する量子ビットで表すことができました。\n",
"この結果をそのまま、上記の二次関数に入れてみましょう。\n",
"\\begin{equation}\n",
"y = a\\left(\\sum_{i=1}^N x_i - b \\right)^2\n",
"\\end{equation}\n",
"ここで二乗というのは同じものを二回かけるという意味であるから次のような意味を持ちます。\n",
"\\begin{equation}\n",
"y = a\\left(\\sum_{i=1}^N x_i - b \\right)\\left(\\sum_{j=1}^N x_j - b \\right)\n",
"\\end{equation}\n",
"ここで和をとる記号のシグマ記号が$i$が$1$から$N$までその値を変えて足し算されるということを思い出すと、別に文字は$i$でなくてもなんでも良い。どうせ$1$から$N$までの値に変わるので。そこで$i$の代わりに2回目の登場時には$j$を利用することにしました。\n",
"これをそのまま展開すると、\n",
"\\begin{equation}\n",
"y = a \\sum_{i=1}^N \\sum_{j=1}^N x_ix_j - 2ab \\sum_{i=1}^N x_i + ab^2\n",
"\\end{equation}\n",
"という結果を得ます。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "fUnqiSkRKdHm"
},
"source": [
"この結果を見てみると、QUBO行列が現れていることに気づきます。\n",
"第一項はそのままですね。係数$a$がQUBO行列に代わります。\n",
"さらに第二項から$-2ab$が$x_i^2=x_i$であることに注意すると、$x_ix_i$の係数として追加されます。第三項はあくまで定数項で最小化には関係しないもの、$x_i$には依存しないもので放っておきましょう。\n",
"\n",
"さてそれではQUBO行列を作るコードを書いてみましょう。ここでは$N=10$とします。\n",
"つまり最小$0$、最大$10$となるような整数の範囲で、二次関数の最小化を考えます。\n",
"答えはもちろん$b=3$になるはずです。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "Eqxy6VC1LEho"
},
"source": [
"N = 10\n",
"QUBO = np.zeros(N**2).reshape(N,N)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "JIQSKNKsLWH7"
},
"source": [
"np.zeros()で全成分が$0$の数列を100個作ります。それを$10 \\times 10$の行列に変形して並べます。\n",
"さてこれに第一項の結果を入れてみましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "pLIiOEEcLiX-"
},
"source": [
"for i in range(N):\n",
" for j in range(N):\n",
" QUBO[i][j] = a"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "gJxFlec7Lsuf"
},
"source": [
"次は**第二項を追加**しましょう。\n",
"$x_ix_i$の係数だけですから、QUBO[i][i]と指定して、対角成分だけに追加されるようにしましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "MWWPTtqGLryh"
},
"source": [
"for i in range(N):\n",
" QUBO[i][i] = QUBO[i][i] - 2*a*b\n"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "bWiNO09BL_1i"
},
"source": [
"これで準備完了です。二次関数の展開からQUBO行列が決まり、量子アニーリングマシンに投じる形に処理することができました。\n",
"\n",
"それでは前回の復習です。量子アニーリングマシンに問題を投じましょう。\n",
"まずはDWaveSampler()で準備をします。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "usBLnckrMN5H"
},
"source": [
"from dwave.system import DWaveSampler, EmbeddingComposite\n",
"\n",
"dw_sampler = DWaveSampler(solver='Advantage_system1.1', token=token)\n",
"sampler = EmbeddingComposite(dw_sampler)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "zUkwY7mAMmbR"
},
"source": [
"ここでは埋め込みが自動的に決まり、そのまま量子アニーリングマシンに問題が投じられるEmbeddingComposite()を利用しました。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "1Xi-nZ7ONdYn",
"outputId": "59242d16-82fc-4a02-9bbe-2724c6ff3478"
},
"source": [
"sampleset = sampler.sample_qubo(QUBO,num_reads=10)\n",
"print(sampleset.record)"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"[([1, 0, 0, 1, 0, 1, 0, 0, 0, 0], -4.5, 1, 0.)\n",
" ([0, 0, 1, 0, 0, 0, 0, 1, 0, 1], -4.5, 1, 0.)\n",
" ([1, 1, 0, 0, 0, 0, 0, 0, 1, 0], -4.5, 1, 0.)\n",
" ([0, 1, 0, 0, 1, 0, 0, 1, 0, 0], -4.5, 1, 0.)\n",
" ([0, 0, 1, 0, 0, 1, 0, 0, 1, 0], -4.5, 1, 0.)\n",
" ([0, 0, 0, 1, 1, 0, 1, 0, 0, 0], -4.5, 2, 0.)\n",
" ([0, 0, 0, 1, 0, 0, 0, 1, 1, 0], -4.5, 1, 0.)\n",
" ([0, 0, 0, 0, 1, 0, 1, 1, 1, 0], -4. , 1, 0.)\n",
" ([0, 1, 0, 0, 0, 1, 1, 1, 0, 0], -4. , 1, 0.)]\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "4Ki7pqTqOW-J"
},
"source": [
"結果を見ると、$x_i=1$となっているのは合計して$3$つになっているはずです。\n",
"たまにしくじっているものもあるかもしれません。\n",
"ちゃんと二次関数の頂点のところに、最小値を持つところに結果を出してくれています。\n",
"-4.5というエネルギー値を出しています。それは$9a$と足算をするとぴったり$0$になります。\n",
"つまり二次関数の最小値である$0$です。ちゃんと答えになっていますね。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "BqKycg8TPy9u"
},
"source": [
"### 数式からQUBO行列にするライブラリpyQUBO\n"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "wMz07Q2jP4Vt"
},
"source": [
"二次関数からなる数式からQUBO行列を得ることができるということがわかりました。\n",
"展開して丁寧にプログラムに起こすことでQUBO行列に相当する部分を取り出すことができます。\n",
"しかし計算間違いや、プログラムの間違いがあり得ますので、それを排除するために便利なライブラリがあります。\n",
"リクルートコミュニケーションズ作成のpyQUBOです。\n",
"D-Wave Ocean SDKにも搭載されているため、そのまま利用することができます。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "BriNdvAwPc7u"
},
"source": [
"pyQUBOの使い方は、まずどんな文字で数式を表すのか、ということを宣言します。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "ifqiNggYP1OJ",
"colab": {
"base_uri": "https://localhost:8080/"
},
"outputId": "c0338ae7-d124-4c13-cd89-2fc9fcb19906"
},
"source": [
"from pyqubo import Array\n",
"\n",
"# バイナリ変数\n",
"x = Array.create(name='x', shape=(N), vartype='BINARY')\n",
"print(x)"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"Array([Binary(x[0]), Binary(x[1]), Binary(x[2]), Binary(x[3]), Binary(x[4]), Binary(x[5]), Binary(x[6]), Binary(x[7]), Binary(x[8]), Binary(x[9])])\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "W2g4Rl1uP91g"
},
"source": [
"まずはpyQUBOからArrayというクラスを呼び出します。\n",
"Array.create()で数式で使われる文字式として、xを用意します。\n",
"そのxはshape=NとありますようにN個の文字式を用意します。その文字は'BINARY'ですから、2値の$0$または$1$をとるよ、と教えています。ここで'SPIN'という選択をすると$x_i=-1$または$x_i=+1$を取るタイプにも変えられます。このタイプはイジング変数と呼びます。\n",
"狙い通りN=10個のBinary(x[i])が用意できたかと思います。これが$x_i$のことだと考えてください。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "tOAlQduaRikS"
},
"source": [
"これを用いて数式を書くかのようにプログラムすることができます。\n",
"例えば先程の二次関数を書いてみましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "KC4tllBeP8Zb"
},
"source": [
"y = a*(sum(x)-b)**2"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "OAia1VowSiP8"
},
"source": [
"数式のままですよね。sum(x）は$\\sum_{i=1}^N x_i$のことです。他はPythonの数式の書き方ですね。$a$をかけたり、$b$をずらしたり、二乗したり。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "qCtjsG1uSxk4"
},
"source": [
"このコスト関数からQUBO行列への変換は、以下の通りです。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "cgnfP-iJS3_m"
},
"source": [
"model = y.compile()\n",
"qubo, offset = model.to_qubo()"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "GRQUobN6TEaK"
},
"source": [
"最初のy.compile()でQUBO行列の変換を行う準備をして、model.to_qubo()により、QUBO行列に変換をすることができます。\n",
"実際に見てみましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "f0c0t3H7TQKQ",
"outputId": "64038436-1b1a-40f3-9184-e631a2aa29a3"
},
"source": [
"print(qubo)"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"{('x[8]', 'x[8]'): -2.5, ('x[7]', 'x[7]'): -2.5, ('x[1]', 'x[1]'): -2.5, ('x[0]', 'x[0]'): -2.5, ('x[4]', 'x[4]'): -2.5, ('x[2]', 'x[2]'): -2.5, ('x[8]', 'x[2]'): 1.0, ('x[6]', 'x[2]'): 1.0, ('x[8]', 'x[1]'): 1.0, ('x[8]', 'x[0]'): 1.0, ('x[2]', 'x[0]'): 1.0, ('x[8]', 'x[5]'): 1.0, ('x[7]', 'x[5]'): 1.0, ('x[7]', 'x[1]'): 1.0, ('x[0]', 'x[1]'): 1.0, ('x[9]', 'x[7]'): 1.0, ('x[9]', 'x[6]'): 1.0, ('x[7]', 'x[4]'): 1.0, ('x[3]', 'x[2]'): 1.0, ('x[7]', 'x[2]'): 1.0, ('x[8]', 'x[6]'): 1.0, ('x[4]', 'x[2]'): 1.0, ('x[5]', 'x[4]'): 1.0, ('x[3]', 'x[0]'): 1.0, ('x[3]', 'x[3]'): -2.5, ('x[4]', 'x[3]'): 1.0, ('x[9]', 'x[0]'): 1.0, ('x[6]', 'x[1]'): 1.0, ('x[7]', 'x[3]'): 1.0, ('x[6]', 'x[0]'): 1.0, ('x[8]', 'x[3]'): 1.0, ('x[6]', 'x[4]'): 1.0, ('x[3]', 'x[1]'): 1.0, ('x[9]', 'x[2]'): 1.0, ('x[5]', 'x[0]'): 1.0, ('x[2]', 'x[1]'): 1.0, ('x[5]', 'x[5]'): -2.5, ('x[9]', 'x[8]'): 1.0, ('x[9]', 'x[3]'): 1.0, ('x[5]', 'x[2]'): 1.0, ('x[6]', 'x[5]'): 1.0, ('x[6]', 'x[6]'): -2.5, ('x[4]', 'x[0]'): 1.0, ('x[9]', 'x[4]'): 1.0, ('x[4]', 'x[1]'): 1.0, ('x[9]', 'x[5]'): 1.0, ('x[5]', 'x[1]'): 1.0, ('x[9]', 'x[1]'): 1.0, ('x[8]', 'x[7]'): 1.0, ('x[9]', 'x[9]'): -2.5, ('x[5]', 'x[3]'): 1.0, ('x[7]', 'x[0]'): 1.0, ('x[6]', 'x[3]'): 1.0, ('x[8]', 'x[4]'): 1.0, ('x[7]', 'x[6]'): 1.0}\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "MlFj7ojsTTXo"
},
"source": [
"これを見ると、どの文字とどの文字の間に係数がかかっているのかというdict形式による記述がされています。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "3N7QfAS4TZAK"
},
"source": [
"これをそのままD-Waveマシンに投じることができます。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "8Hgu6CnwJYMY"
},
"source": [
"dw_sampler = DWaveSampler(solver='Advantage_system1.1', token=token)\n",
"sampler = EmbeddingComposite(dw_sampler)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "3FJfb_ZoJXlE"
},
"source": [
"どのマシンやソルバーを利用するか、埋め込み方法が必要であれば指定したのちに計算を実行してください。\n",
"もちろんマシンタイムの節約のためにシミュレータを使うのも良いでしょう。その場合はこちらを使いましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "TN0eukJCJkZ3"
},
"source": [
"from openjij import SQASampler\n",
"sampler = SQASampler()"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "sa3FOh0eKDw0"
},
"source": [
"pyQUBOではdict形式で問題が記述されているので、そのままOpenJijによるシミュレーターに問題を投入することができます。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "exB_lG9STdlu",
"outputId": "3e44eebc-90ab-43c9-d0c2-a8613cfc77f3"
},
"source": [
"sampleset = sampler.sample_qubo(qubo,num_reads=10)\n",
"print(sampleset.record)"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"[([1, 0, 1, 0, 0, 1, 0, 0, 0, 0], -4.5, 1, 0.)\n",
" ([0, 0, 0, 0, 0, 0, 1, 1, 0, 1], -4.5, 1, 0.)\n",
" ([0, 0, 0, 0, 0, 1, 0, 1, 0, 1], -4.5, 3, 0.)\n",
" ([0, 0, 0, 0, 0, 0, 1, 1, 1, 0], -4.5, 1, 0.)\n",
" ([1, 0, 0, 0, 0, 1, 0, 1, 0, 1], -4. , 1, 0.)\n",
" ([0, 1, 0, 0, 0, 0, 1, 1, 1, 0], -4. , 1, 0.)\n",
" ([0, 0, 1, 0, 1, 1, 0, 0, 1, 0], -4. , 1, 0.)\n",
" ([0, 0, 1, 0, 0, 1, 0, 1, 1, 0], -4. , 1, 0.)]\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "EyldySpFT0hi"
},
"source": [
"同じような結果が出てきましたね。頂点である$b=3$付近のところで答えが出ていると思います。\n",
"ちなみにoffsetというのは定数部分として$x$に関係しない定数の値が得られています。みてみると$ab^2=4.5$となるはずです。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "rOpN9t5xUD2P",
"outputId": "b38a8c83-c713-4436-a699-7d2e49d09022"
},
"source": [
"print(offset)"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"4.5\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "BA63mR7FRAxc"
},
"source": [
"はいバッチリですね！"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "Y2ZPqJoADhGy"
},
"source": [
"### 罰金法"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "k65FLIVtOc18"
},
"source": [
"こうした二次関数からQUBO行列を作ることができるのだ、と知っていれば$x_i=0$または$1$をとる二値変数に対して、制限を加えることができます。先程の例では、みんなが知っているような二次関数を利用すれば$x_i$の合計した値が$b$になるようにできるということがわかりました。\n",
"逆説的に、合計が$b$となるようにしたい場合には、次のような二次関数を与えれば良いということがわかります。\n",
"\\begin{equation}\n",
"E(x) = \\frac{a}{2}\\left( \\sum_{i=1}^N x_i - b \\right)^2\n",
"\\end{equation}"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "DIkpIAQhfdYP"
},
"source": [
"つまりこうした二次関数は、役割を持っていて、$\\sum_{i=1}^N x_i$を$b$にするという条件を満たすためのものであるということになります。\n",
"最適化問題では、**何かと何かを等しくしておいてくださいという等式制約条件**を課すことが往々にしてあります。\n",
"こうした二次関数の頂点を利用して**等式制約条件を実装する方法を罰金法**と呼びます。\n",
"\n",
"ちなみに足し算ですから、どこが$1$になって所望の数値を与えても良いことに注目すると次のような問題に利用することができます。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "GRwQAXxmiHBc"
},
"source": [
"### 例題"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "JkAjeKOViKuL"
},
"source": [
"荷物が複数あり、それらには重さが適当に決められているとする。\n",
"その荷物を決められた個数を取り出すとします。\n",
"そのときに荷物の総量が最大となる組み合わせを考えよう。\n",
"\n",
"まず荷物の個数を$N=10$としましょう。\n",
"そして重さを適当に乱数値で決めておきましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "XXhpPE0wFH4m",
"outputId": "73973570-8e70-4c1b-c9d1-5f2a8ca6fa3b"
},
"source": [
"N = 10\n",
"w = np.random.rand(N)\n",
"print(w)"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"[0.28804661 0.92110721 0.72182371 0.13426707 0.72502264 0.93964647\n",
" 0.08257674 0.29841768 0.14978464 0.8095039 ]\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "0YAZUYpzi_JG"
},
"source": [
"この荷物の中から重いものを探してくれば答えは自ずからわかりますが、まずは例題として、これを量子アニーリングマシンで解くことを考えてみましょう。\n",
"まずは荷物の選択を表す$x_i$という変数を用意しましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "f1xNZqcijpLP"
},
"source": [
"x = Array.create(name='x', shape=N,vartype='BINARY')"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "65sGcTBmjxCn"
},
"source": [
"次に決められた個数というのを$K=3$としましょう。\n",
"3つの荷物を取り出すというわけです。そこで罰金法を利用しましょう。\n",
"\\begin{equation}\n",
"a\\left(\\sum_{i=1}^N x_i - K \\right)^2\n",
"\\end{equation}\n",
"これをpyQUBOを利用して、次のようにコードを書いておきます。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "p1u5xuBkF2Uv"
},
"source": [
"K = 3\n",
"constr = (sum(x)-K)**2"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "50khq8q8j-wk"
},
"source": [
"このconstrに罰金法で利用する二次関数を準備しました。$x_i$の合計がKになるようにしています。適当な係数である$a$は後でかけることにします。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "WjBiKNFMkOKR"
},
"source": [
"次に重さが最大となるように、という点を考慮しましょう。\n",
"基本的に最小化問題を解くのが量子アニーリングマシンですが、最大化をしたい場合は、符号をマイナスにします。\n",
"罰金法を含め、次のようなコスト関数はいかがでしょうか。\n",
"\\begin{equation}\n",
"E(x) = - \\sum_{i=1}^N w_i x_i + a \\left(\\sum_{i=1}^N x_i - K \\right)^2\n",
"\\end{equation}\n",
"第1項では$x_i=1$になるものだけが計算にかかり、その合計が大きくなればなるほど、コスト関数が下がるようになります。つまり重いほどコスト関数的には嬉しいわけです。\n",
"第2項は罰金法ですね。\n",
"それではコスト関数の第1項をコードで書いてみましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "XmKM3XkWjY5N"
},
"source": [
"cost = - np.dot(w,x)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "2Be6QSkYlj4d"
},
"source": [
"と書いても構いませんし、"
]
},
{
"cell_type": "code",
"metadata": {
"id": "KC3G1PgYleNb"
},
"source": [
"cost = 0\n",
"for i in range(N):\n",
" cost = cost - w[i]*x[i]"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "Bhvtz9p5lyB6"
},
"source": [
"と書いても構いません。数式で見たまんまの書き方をすることができるのがpyQUBOの利点です。\n",
"さてこれらを集めて、QUBO行列を作りましょう。\n",
"ここで少し便利な機能を利用しましょう。pyQUBOからConstrainおよびPlaceholderという機能を利用します。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "-qmfYsZclsC8"
},
"source": [
"from pyqubo import Constraint, Placeholder\n",
"cost_func = cost + Placeholder('a')*Constraint(constr, label='one-hot')"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "jbcYn2bemPpl"
},
"source": [
"ここでPlaceholder('a')というのは、罰金法の二次関数にかかる係数です。\n",
"こうしておくことで後で係数を変えるときに、計算時間がかからずにスムーズに繰り返し実験が行えるようになります。\n",
"またConstraint()で挟むことにより、挟まれた罰金法による制約条件を満たしているかどうかを簡単に識別することができます。量子アニーリングマシンでは、そのマシンの性能により、制約条件を必ずしも満たすとは限りません。\n",
"そのため制約条件を違反しているかどうかは重要な評価となります。\n",
"\n",
"準備は整いましたね。\n",
"このコスト関数からQUBO行列を計算する準備をまずは整えます。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "wN1rVylAlswe"
},
"source": [
"model = cost_func.compile()"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "QCuUPCCKm6m2"
},
"source": [
"この操作の後QUBO行列を出してもらうのですが、ここでPlaceholder('a')としてまだ定めていない係数値を決めます。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "EJEEZhjxm5x8"
},
"source": [
"feed_dict = {'a': 2.0}\n",
"qubo, offset = model.to_qubo(feed_dict=feed_dict)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "aGIhyTQfnHqm"
},
"source": [
"そのためにfeed_dict = {'a':2.0}という書き方をしていますが、これは係数$a$を$2.0$にしてくださいね、という命令をしています。\n",
"もしも制約条件を違反する答えが出てしまう場合には、ここの値を変えることになります。\n",
"その手前までの数式の形は変わらないので、いじるのはここだけになります。\n",
"そのための機能がPlaceholder('a')です。\n",
"QUBO行列を作る際に、そうした係数の指定を反映させるためにmodel.to_qubo(feed_dict=**feed_dict**)と書いてあります。\n",
"これで完成。あとは量子アニーリングマシンに投入するだけです。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "5i6zsZgxKtig"
},
"source": [
"まずはどのマシン・ソルバーを使うか、そして埋め込みを実行して準備をします。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "yjZXRfhLKsoH"
},
"source": [
"dw_sampler = DWaveSampler(solver='Advantage_system1.1', token=token)\n",
"sampler = EmbeddingComposite(dw_sampler)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "A69tpRGcKPhb"
},
"source": [
"シミュレータを利用する場合にはコチラ。\n"
]
},
{
"cell_type": "code",
"metadata": {
"id": "-RJoXEtaKQF2"
},
"source": [
"sampler = SQASampler()"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "9gB_isJzKmrh"
},
"source": [
"準備が整いましたらいつも通り答えを出してもらいましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "kcZqmnAgnGOy",
"outputId": "b841064a-c34a-4844-d1d7-285729efe9fc"
},
"source": [
"sampleset = sampler.sample_qubo(qubo,num_reads=10)\n",
"print(sampleset.record)"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"[([1, 1, 0, 0, 0, 0, 0, 0, 0, 1], -20.01865773, 1, 0.)\n",
" ([0, 1, 0, 0, 1, 0, 0, 1, 0, 0], -19.94454753, 1, 0.)\n",
" ([0, 0, 1, 0, 0, 0, 0, 1, 0, 1], -19.82974529, 1, 0.)\n",
" ([0, 1, 0, 0, 1, 0, 0, 0, 1, 0], -19.79591449, 2, 0.)\n",
" ([0, 0, 0, 0, 0, 1, 1, 1, 0, 0], -19.32064089, 1, 0.)\n",
" ([0, 0, 0, 1, 0, 0, 0, 0, 1, 1], -19.09355561, 1, 0.)\n",
" ([0, 0, 0, 1, 0, 0, 1, 0, 0, 1], -19.0263477 , 1, 0.)\n",
" ([1, 0, 0, 0, 0, 1, 0, 1, 1, 0], -17.6758954 , 1, 0.)\n",
" ([1, 0, 0, 0, 0, 0, 0, 0, 1, 0], -16.43783125, 1, 0.)]\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "bXoAW8YxokCc"
},
"source": [
"答えがいくつか出てきたかもしれませんが、ちゃんと$K=3$つ$1$が立っている事を確認してください。その中で一番低いコスト関数値を出した結果で何番目に$1$が立っているかをみましょう。\n",
"それで荷物の重さが重い順に選ばれたかどうかをみてみましょう。大丈夫ですかね。\n",
"もしもうまくいっていなかったら、もう一度やってみたり、feed_dict={'a': 2.0}における値をもう少し大きくしたりしてみてください。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "N5qU0AnuGpYK"
},
"source": [
"制約条件を満たしているかどうかのチェックをするためには次の手順を踏んでいきます。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "T_D-2tLeGnjf"
},
"source": [
"decoded_samples = model.decode_sampleset(sampleset=sampleset, feed_dict=feed_dict)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "jUwPfcqeGugo"
},
"source": [
"得られたsamplesetに対して、decode_sampleset（）を適用します。\n",
"その結果をdecoded_samplesとします。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "VrxiCqPBHbHN"
},
"source": [
"これには答えの情報が格納されており一つ一つ取り出して、.constraints(only_broken=True)とすれば、制約条件を違反しているかどうかを見ることができます。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "9ciM2YNjG1t-",
"outputId": "4ca8ed7a-0774-42c0-cfa9-552b4e285424"
},
"source": [
"for sample in decoded_samples:\n",
" print(sample.constraints(only_broken=True))"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"{}\n",
"{}\n",
"{}\n",
"{}\n",
"{}\n",
"{}\n",
"{}\n",
"{'one-hot': (False, 1.0)}\n",
"{'one-hot': (False, 1.0)}\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "v-b2SUGuHaN9"
},
"source": [
"Falseとあるのは制約条件に違反しており、いくつ違反しているのかという個数も見ることができます。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "LC_1tqc9yvpA"
},
"source": [
"### 量子アニーリングマシンが解くべき問題"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "Z-C6iv4Myy0f"
},
"source": [
"さて上記の例題では、重い順から選べば、わざわざ量子アニーリングマシンを利用しないでも解けるような問題でした。\n",
"しかし荷物の形状によっては、一緒に持つのが難しいものだったり、持ちやすくなったりするものがあるはずです。そうした情報をもしもどこかしらで知ることができて、コスト関数に組み込むことができれば、より**現実的な問題となる**でしょう。\n",
"ある荷物を持ち($x_i=1$)、それと同時に異なる荷物を持っているとき$x_j=1$、コスト関数に登場するような寄与は次のような形です。\n",
"\\begin{equation}\n",
"\\sum_{i \\neq j} W_{ij} x_ix_j\n",
"\\end{equation}"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "0-cRb7iy3Dnr"
},
"source": [
"この$W_{ij}$が同時に持ちやすいか持ちにくいかを示す量ということになります。\n",
"ここでは乱数でそのような量を用意してみましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "acWd-WYPoS9h"
},
"source": [
"W = 0.1*np.random.randn(N**2).reshape(N,N)\n",
"W = W - np.diag(np.diag(W))"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "Po_SQcaZ3blt"
},
"source": [
"ここで行列Wを正負の値の両方をとりうる乱数値にしました。\n",
"対角成分を取り出すnp.diag(W)から対角行列を作るnp.diag()により、Wの対角成分を$0$にするという操作を行っています。\n",
"\n"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "g0hbeL5x391I",
"outputId": "fe1d3c41-b4f1-4d1c-c8ea-63877e919795"
},
"source": [
"print(W)"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"[[ 0. 0.07858091 0.18983013 -0.04886906 -0.12782428 0.09670302\n",
" 0.22284804 -0.12923688 0.07252795 0.00617773]\n",
" [ 0.05648087 0. 0.02300602 -0.07359657 0.00454248 -0.08407988\n",
" 0.21148559 0.1152086 0.05466151 0.0874889 ]\n",
" [ 0.04203252 0.12042396 0. -0.03557292 0.04142517 -0.02474132\n",
" -0.07501519 0.11324322 0.11675406 -0.0258398 ]\n",
" [-0.01135026 -0.04086623 -0.24498328 0. 0.07660923 -0.10099807\n",
" 0.01619057 -0.02053899 -0.08437451 -0.17968568]\n",
" [ 0.15042138 -0.08266102 -0.08705548 -0.06078506 0. 0.19604097\n",
" -0.03289579 0.05366866 -0.02201612 -0.02816382]\n",
" [ 0.01003466 -0.00794726 -0.06399696 -0.19841438 -0.09426893 0.\n",
" -0.18342792 0.10586101 -0.10893612 0.07632177]\n",
" [-0.0657326 0.05548889 -0.12918447 -0.11942815 0.14032939 -0.16042424\n",
" 0. 0.07490986 0.1117644 0.1112174 ]\n",
" [ 0.24678508 0.00340811 0.05875284 0.00688242 -0.05727437 0.09665674\n",
" -0.06822112 0. -0.09048486 -0.13662666]\n",
" [-0.0251032 -0.13767392 0.05599694 0.02218825 -0.1225431 -0.00962018\n",
" 0.06307695 -0.0281562 0. 0.05594002]\n",
" [-0.07856222 -0.04762884 0.12460813 -0.09173484 -0.00320353 -0.13512872\n",
" 0.07859983 -0.09599088 -0.15241805 0. ]]\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "sBJwcsiV3-ae"
},
"source": [
"こうした行列を可視化するためには画像で表示することも直感的でみやすい時があります。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/",
"height": 265
},
"id": "eFD3NFbQ3a3d",
"outputId": "aa3359a3-73f2-4bec-8740-9686ba4de6e2"
},
"source": [
"import matplotlib.pyplot as plt\n",
"plt.imshow(W)\n",
"plt.colorbar()\n",
"plt.show()"
],
"execution_count": null,
"outputs": [
{
"output_type": "display_data",
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAATAAAAD4CAYAAABrN7qeAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAT+ElEQVR4nO3df5BdZX3H8fcnu4SQBfKDMJH8UNIaYVJtxdkiiKKVOIbqkP6hLTjQoGj6QxSsMxqlox3+6GB1rM6UOq5RhwpKNeKQ0UjEqLU/BAmEEZMYEwMmGxbyS35FmmR3v/3j3tB1Sfae5T7n3vvc83nNnMm59558z3c3m+8+z3Oe5xxFBGZmOZrS7gTMzF4oFzAzy5YLmJllywXMzLLlAmZm2eotI2jPqX3RO2t28ri904aTxwSYN+2J5DEH989JHhOg99mSrhrPLed7u+jkA6XEfXjX3OQxR84cTR4TYORwT/KYwwcPMvLMITUT481/0hcHDo4UOvb+nx1eHxHLmjlfGUopYL2zZjPvg9cnjzv7ZQeTxwT4h3PXJo/54S+8K3lMgDM2l1Noeq9/rJS4X37ZV0uJe9XfXpc85lPveSp5TIAnH5mZPObQJz/TdIwDB0f46foXFzq256zt5fxGblIpBczMOl8Ao5TT6mwVFzCzigqCo1GsC9mpXMDMKswtMDPLUhCMZL6U0AXMrMJGcQEzswwFMJJ5ASs0kVXSMknbJO2QtKrspMysNUaJQlunatgCk9QD3Ay8CRgE7pO0NiK2lJ2cmZUngKOZj4EVaYGdD+yIiJ0RcQS4HVheblpmVrYgGCm4daoiY2Dzgd1jXg8Crx5/kKSVwEqAnlmzkiRnZiUKGOnc2lRIssXcETEQEf0R0d/T15cqrJmVpDYTv9jWqYq0wPYAC8e8XlB/z8yyJkZoaj142xUpYPcBiyUtola4LgfeUWpWZla62iB+lxewiBiWdC2wHugBvhQRm0vPzMxKVZsH1uUFDCAi1gHrSs7FzFpstNtbYGbWnbqhBeZbSptVVCBGmFJoa6TRah1Jfydpi6SfSdog6SUpvgYXMLMKGw0V2iYyZrXOpcAS4ApJS8Ydtgnoj4g/BNYA/5Qif3chzSoqEEciyf36n1utAyDp2Gqd55YbRsQPxxx/D3BlihO7gJlVVG0ia+FO2BxJG8e8HoiIgfp+odU6Y1wDfLfoiSdSSgHT1FF65/82edyDv0z/pCOAa/dclTzmzKeThwRgcGk5vf5Zt88vJe4Tf1/O78i5q36VPObQYDnfg3nn7E0ec/+0o0niTGIQf39E9Dd7PklXAv3A65uNBW6BmVVWhBiJJL8QC63WkbQUuAF4fUQcTnFiFzCzChtNM42i4WodSecBnweWRUSyJqkLmFlF1Qbxmy8BJ1qtI+lGYGNErAU+CZwKfEMSwK6IuKzZc7uAmVXUJAfxJ451nNU6EfGxMftLk5xoHBcwswob8VIiM8vRsZn4OXMBM6uw0TRXIdvGBcysomqLuV3AzCxDgTiaZilR27iAmVVUBKkmsraNC5hZZSnVRNa2cQEzq6jALTAzy5gH8c0sS0HjmxV2Ohcws4qqPVYt7xKQd/Zm1oRqPNjWzLpQ4Jn4ZpYxt8DMLEsRcgvMzPJUG8T3UiIzy1Kye+K3TSkFLAKGj6Sv7Bde8IvkMQF+cs+5yWMePiN5yJozkjwL4XmeefG0UuJe/YkPlBL3996xPXnMVy0cTB4T4MEN5ySPOXzopKZj1AbxPQZmZpnyTHwzy5Jn4ptZ1lI91KNdXMDMKioCjo66gJlZhmpdSBcwM8uUZ+KbWZa6YRpFw/ajpIWSfihpi6TNkq5rRWJmVrZaF7LI1qmKtMCGgQ9GxAOSTgPul3R3RGwpOTczK1nX3xM/IoaAofr+05K2AvMBFzCzjNWuQlZoLaSks4HzgHuP89lKYCVAzxkzEqRmZmXqhomshTu3kk4FvglcHxFPjf88IgYioj8i+ntO70uZo5mVZLT+aLVGW6cq1AKTdBK14nVbRNxRbkpm1grdcBWyYQGTJOCLwNaI+HT5KZlZq3TyFcYiirTALgKuAh6S9GD9vY9GxLry0jKzskWI4W4vYBHxX9DBnWAze8G6vgtpZt2pG8bA8m4/mllTRkOFtkYkLZO0TdIOSauO8/nFkh6QNCzpbanydwEzq6hj88CaLWCSeoCbgUuBJcAVkpaMO2wXcDXw1ZRfg7uQZhWWaI7X+cCOiNgJIOl2YDljVutExCP1z0ZTnPCYcgrYkSlM2ZP+IRGbto4v6mm89B//J3nMHZ++IHlMgJkzD5US96RtJ5cS9/GLkv68Pueh/1icPOZpjyQPCcCie36TPObjvxlpOkYEDBe/oeEcSRvHvB6IiIH6/nxg95jPBoFXN51gAW6BmVXYJAbx90dEf5m5vBAuYGYVlXAt5B5g4ZjXC+rvlc6D+GYVFqFCWwP3AYslLZI0FbgcWFt68riAmVVaisXcETEMXAusB7YCX4+IzZJulHQZgKQ/ljQIvB34vKTNKfJ3F9KsoiLSTWStLy1cN+69j43Zv49a1zIpFzCzyhIjfqyameWqwPhWR3MBM6uoblgL6QJmVlVRGwfLmQuYWYV18u2ii3ABM6uo8CC+meXMXUgzy5avQppZliJcwMwsY55GYWbZ8hiYmWUpEKO+Cmlmucq8AeYCZlZZHsQ3s6xl3gRzATOrMLfAjuOMGU/zjkt/nDzud/7l4uQxAQ5cc2HymKfuKucH4wlmlxJ34MYvlBL32lv+qpS4I0ueSR7zafqSxwQ4tCD9v9nhz/c0HSOA0VEXMDPLUQBugZlZrjwPzMzy5QJmZnkq9Mi0juYCZlZlboGZWZYCwlchzSxfeRewwis5JfVI2iTp22UmZGYtFAW3DjWZpejXUXtsuJl1iyoUMEkLgLcAq8tNx8xa5thE1iJbhyo6BvYZ4EPAaSc6QNJKYCXA6Wed0nxmZla63CeyNmyBSXorsDci7p/ouIgYiIj+iOjvmzU1WYJmVqJRFds6VJEW2EXAZZL+FJgGnC7p1oi4stzUzKxs6vYWWER8JCIWRMTZwOXAD1y8zLpA0QH8Di5yngdmVlmdPUBfxKQKWET8CPhRKZmYWet1cOuqCLfAzKpstN0JNMcFzKyquuCGhnk/FM7MmqIotjWMIy2TtE3SDkmrjvP5yZL+vf75vZLOTpG/C5hZlSW4CimpB7gZuBRYAlwhacm4w64BfhMRLwX+GfhEivRdwMysWecDOyJiZ0QcAW4Hlo87ZjlwS31/DXCJpKb7r6WMgR14+jS+8uPXJo/bszh5SACm70k/DrDgzj3JYwLsu3heKXFX8p5S4p7zzQOlxJ3yxqeSxzz4renJYwIcvjL992DK9OEkcSYxkXWOpI1jXg9ExEB9fz6we8xng8Crx/39546JiGFJTwJnAPsnm/NYHsQ3q6pgMsuE9kdEf4nZvCDuQppVWZqZ+HuAhWNeL6i/d9xjJPUCM4Cmm6YuYGYVlugq5H3AYkmLJE2ltuRw7bhj1gIr6vtvo7YkselptO5CmlVZgpn49TGta4H1QA/wpYjYLOlGYGNErAW+CHxF0g7gILUi1zQXMLMqS7SUKCLWAevGvfexMfv/C7w9zdn+nwuYWUUVnaTayVzAzKqsg29WWIQLmFmFuQVmZvlyATOzLHkMzMyy5gJmZrlS5jc09Ex8M8uWW2BmVeYupJllyYP4ZpY1FzAzy5YLmJnlSOR/FdIFzKyqPAZmZllzATOzbLmAPV/PYTh9e0/yuM9e9EzymAAzfzwtecynzntR8pgAr3v/vaXEfeh9rygl7rZrZpUSd0b6hxJx9nt/lT5oSR7uGUkSx11IM8uXC5iZZSl8FdLMcuYWmJnlymNgZpYvFzAzy1Kxp253NBcws4oS+XchC93QUNJMSWsk/ULSVkkXlp2YmZXv2LMhG22dqmgL7LPAXRHxNklTgekl5mRmrdLBxamIhgVM0gzgYuBqgIg4AhwpNy0za4nMC1iRLuQiYB/wZUmbJK2W1Df+IEkrJW2UtHH42UPJEzWzxAp2Hzu5C1mkgPUCrwI+FxHnAYeAVeMPioiBiOiPiP7eU55X38ysE0XBrUMVKWCDwGBEHFtFvIZaQTOzzGm02NapGhawiHgM2C3pnPpblwBbSs3KzFoi9y5k0auQ7wNuq1+B3Am8s7yUzKwlOrx7WEShAhYRDwL9JediZq1WhQJmZt2nG2biu4CZVZhG865gLmBmVdUFY2CF1kKaWXdqxVVISbMl3S1pe/3P4z4oQdJdkp6Q9O2isV3AzKqsNRNZVwEbImIxsIHjTISv+yRw1WQCl9KFXPKiffx01b8mj/v7X//r5DEB9q5Iv/Tp8FA5693vvu2CUuKeeeNgKXFPu7OcVRlPkv5pR7sfmp08JsCpQ8PJY44MTU0Sp0WD+MuBN9T3bwF+BHx4/EERsUHSG8a/PxGPgZlVWfECNkfSxjGvByJioODfnRsRQ/X9x4C5hc/agAuYWVVN7qlE+yPihHNBJX0fON7DUG/4nVNGhJSu3ecCZlZRKeeBRcTSE55HelzSWRExJOksYG+as3oQ36zaIoptzVkLrKjvrwDubDbgMS5gZhXWosXcNwFvkrQdWFp/jaR+Saufy0X6T+AbwCWSBiW9uVFgdyHNqqpFE1kj4gC1u9iMf38j8O4xr1832dguYGYV1sn3+irCBcyswlzAzCxPQYoB+rZyATOrMN9Ox8zy5QJmZjnyDQ3NLF8RvqGhmWUs7/rlAmZWZe5CmlmeAnAX0syylXf9cgEzqzJ3Ic0sW74KaWZ56oLHqpVSwH6+/0zOXf03yeO++L+PJo8J8Oi7lDxmnFzOKtmZO8v5idt1//xS4k59w5PlxH1oRvKYB19Rzr/ZU8sOJ495dFvzPwe1iax5VzC3wMyqzHejMLNcuQVmZnnyGJiZ5ctrIc0sZ+5CmlmWJvdg247kAmZWZZm3wAo9F1LSByRtlvRzSV+TNK3sxMysBaLg1qEaFjBJ84H3A/0R8XKgB7i87MTMrHwaHS20daqiXche4BRJR4HpwKPlpWRmLRFkP5G1YQssIvYAnwJ2AUPAkxHxvfHHSVopaaOkjSOHDqXP1MySEoGi2NapinQhZwHLgUXAPKBP0pXjj4uIgYjoj4j+nr6+9JmaWXoRxbYOVWQQfynwcETsi4ijwB3Aa8pNy8xaIvMCVmQMbBdwgaTpwLPAJcDGUrMys/J1wRhYwwIWEfdKWgM8AAwDm4CBshMzs/J18hXGIgpdhYyIjwMfLzkXM2upzu4eFuGZ+GZVFbiAmVnG8u5BuoCZVVknz/EqwgXMrMoyL2CFFnObWReKgJHRYlsTJM2WdLek7fU/Zx3nmFdK+kn9phE/k/QXRWKX0gKbchT69qSPu++VU9MHBV7zkm3JY27+7suTxwR49LXl/MZcsGG4lLiH/qCcJ0lNKeHpTKf/Ov3TgwB2Lp+ePGYcSdT2aE0LbBWwISJukrSq/vrD4475LfCXEbFd0jzgfknrI+KJiQK7BWZWZa2Zib8cuKW+fwvwZ89PI34ZEdvr+48Ce4EzGwX2GJhZVQVQ/J74cySNXYEzEBFFJ7TPjYih+v5jwNyJDpZ0PjAV+FWjwC5gZpUVEIXHt/ZHRP+JPpT0feBFx/noht85Y0RIOmHVlHQW8BVgRUTj5FzAzKoqaHqA/rlQEUtP9JmkxyWdFRFD9QK19wTHnQ58B7ghIu4pcl6PgZlVWWvGwNYCK+r7K4A7xx8gaSrwLeDfImJN0cAuYGZV1poCdhPwJknbqd2e6yYASf2SVteP+XPgYuBqSQ/Wt1c2CuwupFlltWYxd0QcoHYbrvHvbwTeXd+/Fbh1srFdwMyqKoAq3E7HzLpU5kuJXMDMKiuSXYVsFxcws6oKKDDVqqO5gJlVWfGZ+B3JBcysyjwGZmZZivBVSDPLmFtgZpanIEZG2p1EU1zAzKpqcrfT6UguYGZV5mkUZpajAMItMDPLUkzqhoYdyQXMrMJyH8RXlHAZVdI+4NcFDp0D7E+eQHlyyjenXCGvfDsh15dERMOHXkxE0l3UvpYi9kfEsmbOV4ZSCljhk0sbJ7rPdqfJKd+ccoW88s0p127nO7KaWbZcwMwsW+0uYEWfK9cpcso3p1whr3xzyrWrtXUMzMysGe1ugZmZvWAuYGaWrbYVMEnLJG2TtEPSqnbl0YikhZJ+KGmLpM2Srmt3TkVI6pG0SdK3253LRCTNlLRG0i8kbZV0YbtzmoikD9R/Dn4u6WuSprU7pyprSwGT1APcDFwKLAGukLSkHbkUMAx8MCKWABcA7+3gXMe6Dtja7iQK+CxwV0ScC/wRHZyzpPnA+4H+iHg50ANc3t6sqq1dLbDzgR0RsTMijgC3A8vblMuEImIoIh6o7z9N7T/Y/PZmNTFJC4C3AKsbHdtOkmZQexrzFwEi4khEPNHerBrqBU6R1AtMBx5tcz6V1q4CNh/YPeb1IB1eFAAknQ2cB9zb3kwa+gzwIaDTV+ouAvYBX653d1dL6mt3UicSEXuATwG7gCHgyYj4XnuzqjYP4hck6VTgm8D1EfFUu/M5EUlvBfZGxP3tzqWAXuBVwOci4jzgENDJ46GzqPUUFgHzgD5JV7Y3q2prVwHbAywc83pB/b2OJOkkasXrtoi4o935NHARcJmkR6h1zd8o6db2pnRCg8BgRBxr0a6hVtA61VLg4YjYFxFHgTuA17Q5p0prVwG7D1gsaZGkqdQGQte2KZcJSRK1MZqtEfHpdufTSER8JCIWRMTZ1L6vP4iIjmwlRMRjwG5J59TfugTY0saUGtkFXCBpev3n4hI6+KJDFbTlfmARMSzpWmA9tSs5X4qIze3IpYCLgKuAhyQ9WH/voxGxro05dZP3AbfVf5HtBN7Z5nxOKCLulbQGeIDa1elNeFlRW3kpkZlly4P4ZpYtFzAzy5YLmJllywXMzLLlAmZm2XIBM7NsuYCZWbb+DwEX+ceanPHrAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 2 Axes>"
]
},
"metadata": {
"tags": [],
"needs_background": "light"
}
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "zTbVU6Bk4Ih1"
},
"source": [
"それではこの効果を取り入れた新しいコスト関数を用意しましょう。\n",
"pyQUBOの形式にするために次のようなコードを書きます。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "AZDV5bUT32-5"
},
"source": [
"cost2 = 0\n",
"for i in range(N):\n",
" for j in range(N):\n",
" cost2 = cost2 + W[i][j]*x[i]*x[j]"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "HKDCmGfR5I7I"
},
"source": [
"先程のコスト関数から、今用意した項も追加することにしましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "zP1nTlht4RBx"
},
"source": [
"cost_func2 = cost + cost2 + Placeholder('a')*Constraint(constr, label='Kconstr')"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "CVPV6BkC5Yjq"
},
"source": [
"先ほどまで登場した項に対して、さらにもう一つcost2が追加されていますね。\n",
"あとはQUBO行列を出力するための準備を整えます。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "iMpP8A5E5U0I"
},
"source": [
"model2 = cost_func2.compile()"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "Fer6lVkG5iDC"
},
"source": [
"そしていつも通りQUBO行列の出力をしましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "HLlusSXl5XyX"
},
"source": [
"feed_dict = {'a': 2.0}\n",
"qubo2, offset = model2.to_qubo(feed_dict=feed_dict)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "AbG-4U7Q5wV6"
},
"source": [
"どんな答えが出てくるでしょうかね。さすがにわからないですよね。\n",
"そうした問題に対しては、量子アニーリングマシンを利用すると良いでしょう。\n",
"いつものようにソルバーを指定して、埋め込みをして準備をします。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "NM2tiZzpK36s"
},
"source": [
"dw_sampler = DWaveSampler(solver='Advantage_system1.1', token=token)\n",
"sampler = EmbeddingComposite(dw_sampler)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "eBB4IpgiK-Dt"
},
"source": [
"シミュレータを利用したい場合にはOpenJijでOKですね。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "X0FbICrVLBPm"
},
"source": [
"sampler = SQASampler()"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "os8HcbpdLFNH"
},
"source": [
"これで準備完了です。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "LUSRzg7d5sqv",
"outputId": "fb45b583-4aa0-4e92-b02b-bca4d6cce2d9"
},
"source": [
"sampleset = sampler.sample_qubo(qubo2,num_reads=10)\n",
"print(sampleset.record)"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"[([0, 0, 0, 1, 0, 1, 0, 0, 0, 1], -20.51305737, 1, 0.)\n",
" ([0, 0, 1, 1, 0, 0, 0, 0, 0, 1], -20.11880306, 2, 0.)\n",
" ([0, 0, 0, 0, 0, 1, 1, 0, 0, 1], -20.04456898, 1, 0.)\n",
" ([0, 1, 1, 0, 0, 0, 0, 0, 1, 0], -19.55954699, 1, 0.)\n",
" ([0, 1, 1, 0, 0, 0, 1, 0, 0, 0], -19.51930286, 3, 0.)\n",
" ([0, 0, 0, 0, 0, 0, 1, 1, 0, 1], -19.22660989, 1, 0.)\n",
" ([1, 0, 1, 1, 0, 0, 0, 0, 0, 1], -18.30759083, 1, 0.)]\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "Jfg2RK-x6ZJu"
},
"source": [
"先ほどとは異なる結果になったはずです。重さだけでいうと簡単だけど、荷物の組み合わせによって影響が出てくる。\n",
"ある荷物を持つと異なる荷物の影響が出てくるといった問題では、量子アニーリングが有効となります。\n",
"ポイントは二次関数を利用して、制約条件を課すことで現実的な問題を考える。\n",
"そして$x_i$が足し算されるようなただの一次関数だけではなく、\n",
"**$x_ix_j$が入った二次関数の問題で、互いの影響を考慮しなければならない問題を考える。**\n",
"そこで量子アニーリングの活路が見出されるということになります。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "InTp27KF6zbJ"
},
"source": [
"### ポートフォリオ最適化問題"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "XuWEIKQl633i"
},
"source": [
"株式銘柄や金融資産を保有する際に、基本的には期待収益と言いますが、\n",
"儲かりそうな金融資産を手に入れるというのが基本です。\n",
"ただそのときに、自分の予算の限界があり、全てを買うわけにはいかないため、制限があります。\n",
"こうした状況で、最適な金融資産の選択を行う問題を**ポートフォリオ最適化問題**と言います。\n",
"\n",
"実際の株式銘柄のデータを活用して、量子アニーリングマシンを利用してこの問題を解いてみましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "72qs9abK52hE"
},
"source": [
"import datetime\n",
"start = datetime.datetime(2018, 1, 1)\n",
"end = datetime.datetime(2020, 12, 31)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "mE_6bA4d8ZBd"
},
"source": [
"ここではまず株式銘柄のデータをどこからどこまでを利用するかの期間を指定しています。\n",
"次に東京証券取引所の証券コードをいくつか選んだリストを用意します。\n",
"pythonでは[]の中にカンマ区切りで指定して複数の要素を並べたものをリストと呼びます。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "Q_Gbz-By6X0q"
},
"source": [
"stockcodes=[\"2502\", \"3382\", \"4661\", \"6178\", \"6758\", \"7203\", \"8053\", \"8604\", \"9020\", \"9433\"]\n",
"## アサヒ, セブン&アイ, オリエンタルランド, 日本郵政, ソニー, トヨタ自動車, 住友商事, 野村ホールディングス, JR東日本, KDDI"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "w2Tsh32T8wQ6"
},
"source": [
"もちろんこれ以外にもご自身で興味のある証券コードで指定していただいても構いません。\n",
"この証券コードで指定された銘柄について、\n",
"start時点で購入した株式が、毎日の取引で変動するわけですが、その日の終値をみて、\n",
"どれだけの割合増減したのか収益率を調べることにします。\n",
"少々長めのfor文を利用したコードとなります。\n",
"コード内のコメントで説明しながらにしましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "3Xg4pO1q8vKm"
},
"source": [
"#証券データの読み込みに利用するライブラリ\n",
"import pandas as pd\n",
"import pandas_datareader.stooq as stooq\n",
"\n",
"rates = []\n",
"#for文で各証券コードについて処理を繰り返していく\n",
"for sc in stockcodes: \n",
" # f\"{sc}.jp\"で証券を指定して、start,endで始まりと終わりの時期を指定のうえ、.read()読み込みます。その中でも終値だけを知りたいので['Close']としました。\n",
" df = stooq.StooqDailyReader(f\"{sc}.jp\", start, end).read()[['Close']]\n",
"\n",
" # pandasのデータフレーム形式で結果を得ることができます。その結果はsort_valuesで日付'Date'で昇順（小さい順）に並べています。\n",
" df = df.sort_values(by='Date',ascending=True)\n",
" \n",
" #前日の終値と今日の終値を比較して、前日比を調べています。\n",
" return_rate = np.zeros(len(df.values))\n",
" for k in range(len(df.values)-1):\n",
" return_rate[k+1] = (df.values[k+1][0] - df.values[k][0])/df.values[k][0]\n",
" \n",
" rates.append(return_rate)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "PKZxmIxxC_uU"
},
"source": [
"この結果、ratesは、それぞれの証券コードの収益率（前日の終値と当日の終値の差を前日の終値終値との比をとったもの、要するに前の日よりも儲かるかどうか）が格納されたリストが出来上がります。\n",
"appendで、結果を次から次へと格納することができます。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "tag91Xk9DYTG"
},
"source": [
"これで下準備は終了しました。\n",
"その株式銘柄の良さは、平均的な収益率から決められそうですね。\n",
"平均的な収益率は大きい方が嬉しいので、先程の荷物の問題と同じです。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "ZUevGdPi_ToR"
},
"source": [
"N = len(stockcodes)\n",
"\n",
"w = np.zeros(N)\n",
"for k in range(N):\n",
" w[k] = rates[k].mean()"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "BsfKHyR6Epdx"
},
"source": [
"rates[k]から各株式銘柄の収益率の時系列データが出てきます。\n",
"それを平均するために、numpyの機能であるmean()を利用しました。\n",
"もっとわかりやすくnp.mean(rates[k])としても同じです。"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "hhqFxtruE6V9"
},
"source": [
"この株式銘柄のうちいくつかを選んで、平均的な収益率を高めるとしましょう。\n",
"これまでと同じようにまずは、文字式の準備をしましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "T-bHrLHy_eZm"
},
"source": [
"x = Array.create('x', shape=N, vartype='BINARY')"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "Pr3UoWDWFn3U"
},
"source": [
"株式銘柄をいくつ選ぶかということを単純にいつもの罰金法を利用して考えてみましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "xAwn2_BKFxIe"
},
"source": [
"K = N//2\n",
"constr = (sum(x)-K)**2"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "VUJcCPkVF44P"
},
"source": [
"ここでN//2としたのは、ただの割り算ではなく、偶数でも奇数でも整数の答えを返してくれる計算方法を利用しています。3//2=1、5//2=2といった具合です。割り切れる場合にはただの割り算と同じ（4//2=2）です。\n",
"次に収益率を最大化するように選択するために、コスト関数を次のようにします。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "O_sIEeYc_z-9"
},
"source": [
"cost = 0\n",
"for i in range(N):\n",
" cost = cost - w[i]*x[i]"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "bxixpf8aGUVL"
},
"source": [
"先ほどまでと同じですね。これで最適化をする準備は整いました。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "Da30NxsZGTuI"
},
"source": [
"cost_func = cost + Placeholder('a')*Constraint(constr, label='Kconstr')\n",
"model = cost_func.compile()"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "bhuqJIdkIAnq"
},
"source": [
"これでQUBO行列を書き出す準備が整いました。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "sYENxSm3H_-3"
},
"source": [
"max_coeff = np.max(abs(w))\n",
"\n",
"feed_dict = {'a': 10.0*max_coeff}\n",
"qubo, offset = model.to_qubo(feed_dict=feed_dict)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "uWmibeXELP8y"
},
"source": [
"いつもの通りD-Waveマシンを利用する場合には、その準備をしましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "DBDsH27gLS5b"
},
"source": [
"dw_sampler = DWaveSampler(solver='Advantage_system1.1', token=token)\n",
"sampler = EmbeddingComposite(dw_sampler)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "zik98XB_LTys"
},
"source": [
"シミュレータを活用する場合はコチラですね。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "uB28ugJhLUAO"
},
"source": [
"sampler = SQASampler()"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "Nqi2dxdQIKyj",
"outputId": "fc833595-fd19-47c2-d66f-0dbf73214161"
},
"source": [
"sampleset = sampler.sample_qubo(qubo,num_reads=10)\n",
"print(sampleset.record)"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"[([0, 0, 1, 0, 1, 0, 1, 1, 1, 0], -0.28074252, 1, 0. )\n",
" ([0, 0, 1, 1, 1, 0, 0, 1, 1, 0], -0.28050151, 1, 0. )\n",
" ([0, 0, 0, 1, 1, 1, 1, 1, 0, 0], -0.28027962, 1, 0. )\n",
" ([0, 1, 0, 0, 1, 1, 1, 0, 1, 0], -0.28018786, 1, 0. )\n",
" ([0, 0, 1, 1, 0, 1, 1, 1, 0, 0], -0.27998866, 1, 0. )\n",
" ([0, 1, 1, 1, 0, 0, 1, 1, 0, 0], -0.27955377, 3, 0. )\n",
" ([0, 1, 0, 1, 0, 0, 1, 1, 1, 0], -0.27822792, 1, 0. )\n",
" ([1, 0, 0, 1, 0, 0, 1, 1, 0, 1], -0.27910096, 1, 0.1)]\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "DuUH7GCOLncm"
},
"source": [
"バッチリとどの株式銘柄が良いかを選んでくれました。\n",
"しかしながら、実際の株式市場は複雑です。単純に平均的な指標だけで考えたポートフォリオが安定的に収益を得るかどうかはわかりません。\n",
"ある銘柄の収益率が上昇すると、別の銘柄も収益率が上昇するかもしれません。\n",
"または下降するかもしれません。そうした関係性を考慮して、もっと良いポートフォリオの組み方があるかもしれません。\n",
"そこで荷物のときと同じようにある銘柄を持つと他の銘柄とどのような関係にあるのかを考慮して、コスト関数に組み込んでみましょう。\n",
"\n",
"確率的に変動していく数に対して、その平均からのズレを見てみましょう。\n",
"そのズレが大きいほど変動幅は大きいということになります。\n",
"こうした平均からのズレの大きさ（二乗で計算する）を足し上げたものを**分散**と言います。\n",
"その変動幅が大きいものはリスクを伴う（儲かるかもしれないけど、暴落するかもしれない）ので、それはできることなら避けたいはずです。\n",
"\n",
"またある銘柄で平均からズレたときに、影響して他の銘柄でも平均からズレるかもしれません。\n",
"そのズレの積を取ってたしあげたものを**共分散**と言います。\n",
"一度下がり始めた場合に、同じように下がるようなものは共分散の値が正の値とります。\n",
"それは一緒に崖から落ちるようなものなので、避けたいです。選びたくないです。\n",
"そうした考えから次のようなコスト関数を追加することを思いつきます。\n",
"\\begin{equation}\n",
"\\sum_{i=1}^N \\sum_{j=1}^N W_{ij} x_i x_j\n",
"\\end{equation}\n",
"ここで$W_{ij}$は分散・共分散行列となります。\n",
"こうして定式化されたポートフォリオ最適化問題は、Markowitzの現代ポートフォリオ理論の基本的な手法の一つとなっています。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "5iO-4-joIkHD"
},
"source": [
"cost2 = 0\n",
"for i in range(N):\n",
" for j in range(N):\n",
" cost2 = cost2 +x[i]*x[j]*np.sum((rates[i]-w[i])*(rates[j]-w[j]))/len(rates[i])"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "nb9l3JlrOtM5"
},
"source": [
"これを足した新しいコスト関数でリスクの最小化も考慮したポートフォリオを設計することができます。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "wIzJqgDUOPrX"
},
"source": [
"cost_func2 = cost2 + Placeholder('a')*Constraint(constr, label='Kconstr')\n",
"model2 = cost_func2.compile()"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "QF4MRRE_O5nW"
},
"source": [
"QUBO行列を出力する準備が整いました。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "jVjCpAw1O4uq"
},
"source": [
"max_coeff = np.max(abs(w))\n",
"\n",
"feed_dict = {'a': 10.0*max_coeff}\n",
"qubo2, offset = model2.to_qubo(feed_dict=feed_dict)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "HowQDe4aPDLm"
},
"source": [
"この後の流れもお手の物ですね。D-Waveマシンのうちどれを使うかを選び、埋め込みの準備をする。"
]
},
{
"cell_type": "code",
"metadata": {
"id": "MaNdVJYSPCav"
},
"source": [
"dw_sampler = DWaveSampler(solver='Advantage_system1.1', token=token)\n",
"sampler = EmbeddingComposite(dw_sampler)"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "pzLF2na6PMQ0"
},
"source": [
"お値段のかからないシミュレータを使いたい場合はコチラ（マシンタイム使い切ってしまった場合もコチラ！）"
]
},
{
"cell_type": "code",
"metadata": {
"id": "FM4O7W28PQeK"
},
"source": [
"sampler = SQASampler()"
],
"execution_count": null,
"outputs": []
},
{
"cell_type": "markdown",
"metadata": {
"id": "WLQC1ovHPTG0"
},
"source": [
"QUBO行列を投じて、それでは問題を解いてみましょう。"
]
},
{
"cell_type": "code",
"metadata": {
"colab": {
"base_uri": "https://localhost:8080/"
},
"id": "DpT1FGS0PQ2q",
"outputId": "6fa21279-3724-4d6f-e56b-61518cda3757"
},
"source": [
"sampleset = sampler.sample_qubo(qubo2,num_reads=10)\n",
"print(sampleset.record)"
],
"execution_count": null,
"outputs": [
{
"output_type": "stream",
"text": [
"[([0, 1, 0, 1, 0, 1, 0, 0, 1, 1], -0.27644783, 1, 0.)\n",
" ([1, 1, 1, 1, 0, 0, 0, 0, 0, 1], -0.27628292, 1, 0.)\n",
" ([0, 1, 0, 0, 1, 0, 1, 0, 1, 1], -0.27626747, 1, 0.)\n",
" ([0, 0, 1, 0, 1, 0, 0, 1, 1, 1], -0.27596826, 1, 0.)\n",
" ([1, 1, 1, 1, 0, 0, 1, 0, 0, 0], -0.27593477, 1, 0.)\n",
" ([1, 0, 1, 1, 0, 0, 1, 0, 1, 0], -0.27574202, 1, 0.)\n",
" ([1, 0, 0, 0, 1, 0, 0, 1, 1, 1], -0.2757227 , 1, 0.)\n",
" ([1, 0, 1, 0, 1, 1, 1, 0, 0, 0], -0.27563545, 1, 0.)\n",
" ([1, 0, 1, 0, 1, 1, 0, 0, 1, 0], -0.27560824, 1, 0.)\n",
" ([1, 1, 1, 0, 0, 0, 1, 0, 1, 0], -0.27545233, 1, 0.)]\n"
],
"name": "stdout"
}
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "ngEMcg1APe5U"
},
"source": [
"先ほどから入れ替わったものもあれば、そのまま利用される株式銘柄もあり、味わい深い結果になったかと思います。\n",
"実際の運用に活かしてみてはいかがでしょうか？\n"
]
},
{
"cell_type": "markdown",
"metadata": {
"id": "pLbFcgdEPmp-"
},
"source": [
"コードサンプル提供（東北大学大学院情報科学研究科・井口大輔・渡辺大地）"
]
}
]
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment