tutorial/108_trafficflow.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "obQW5jfvvg7x"
   },
   "source": [
    "# Traffic Flow Optimization with XYmixer on QAOA\n",
    "The traffic optimization problem is one of the most famous optimization problems in quantum computing."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "_EDJAEs4zaMM"
   },
   "source": [
    "## References\n",
    "Traffic flow optimization using a quantum annealer  \n",
    "Florian Neukart, Gabriele Compostella, Christian Seidel, David von Dollen, Sheir Yarkoni, Bob Parney   \n",
    "https://arxiv.org/abs/1708.01625"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "Fdh7XfRtz17m"
   },
   "source": [
    "## Example\n",
    "Consider a network of 9 intersections and 12 streets."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "colab_type": "text",
    "id": "Fdh7XfRtz17m"
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAb4AAAEuCAYAAADx63eqAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/Z1A+gAAAACXBIWXMAAAsTAAALEwEAmpwYAABF8UlEQVR4nO3deVxU9f4/8NcZZjgDw+Iusikqm4KImKmVJZZbaiZKmqVpV716M6ufZXt2u7fFa2nmvbfsppmlZS4V6c20cilb3AERcAFZFETUAQZmmDNzfn/0PXNB0GSYmc+Zc97Px4OHhc6Z9+gwr/PZOVEURRBCCCEqoWFdACGEEOJJFHyEEEJUhYKPEEKIqlDwEUIIURUKPkIIIapCwUcIIURVKPgIIYSoCgUfIYQQVaHgI4QQoioUfIQQQlSFgo8QQoiqUPARQghRFQo+QgghqqJlXQAhhMiNKIqw2WwQRRGiKILjOHAcBx8fH3Acx7o80koUfIQQ1ZOCzmKxQBAEXO+0No7joNVqwfM8BaGX4ug8PkKIWomiCIvFAovFct2wuxaO48DzPHiepwD0IhR8hBBVEgQBJpPJqcC7GsdxMBgM0GqpE80bUPARQlRFFEWYzWZYLBaXX5vneej1emr9yRwFHyFENURRhMlkgiAIbnsOrVYLg8FA4SdjtJyBEKIKngg9wLVdqMQ9KPgIIapgNpvdHnoSQRBgNps98lyk5Sj4CCGKJwiCW8b0rkdaGkHkh4KPEKJoUhcnC9TlKU8UfIQQRXN2jZ4rSOsEibxQ8BFCFEsOwcMyeEnzKPgIIYol7bfJkrQdGpEPCj5CiGK5orW3efNm3HzzzQgLC0NycjL279/PpA7iOrS/DiFEsVo7q/KHH37A4sWLsXr1aqSkpKCsrMyp61CLT15o5xZCiCKJogij0diqawwfPhwPPvggHnzwwVbXExwcTLu5yAR1dRJCFKm1rSybzYajR4/i4sWL6NevH3r37o0nn3wSdXV1TOohrkPBRwhRpNZ2Zl24cAFWqxVfffUVtm/fjr179yIzMxNLly5lUg9xHQo+QogitTZo/Pz8AACzZ89GSEgI2rdvj7/85S/YuXMnk3qI61DwEUIUqzVh06ZNG4SGhrpsXI7G9+SDZnUSQrxabW0t8vPzkZeXh7y8POTm5iIvLw/+/v5Yv349AgMDnb721KlTsWrVKgwbNgw6nQ7//ve/MWLECKeuRcEnHzSrkxAie6IoorS0tFGwSf994cIF9OjRA7GxsYiLi0NsbKzjq7WsViuefvppbNq0CXq9HuPHj8fLL78MvV7f4mvRrE75oOAjiiXtmCGKIkRRBMdx4DgOPj4+9AEkUw1bbw0DLi8vD4GBgY5Aaxhw3bp1g4+PT7PXMxqNshhb02g0CAoKYl0G+T8UfEQxpKCTjoO53lub4zhotVrwPK/aIJRuBiQXL15EQUEB2rVrhx49ejT5fVc+b0lJSZOuyby8PEfrTQq2hgEXHBzc4ucymUywWq0ufw0tpdPpYDAYWJdB/g8FH/F60kbEzm4GzHEceJ4Hz/OqCUCbzQYfHx/s378fjz/+ODiOQ1VVFS5evIiXX34Zc+fObXXwmUymZsfe8vPzHa23hsEWFxeHrl27XrP15gxBEFBTU+Oy6zkrICAAWi1NqZALCj7i1QRBcNmZZxzHwWAwqOoDqq6uDnl5eejbty/Ky8vxpz/9Ce+++y7CwsJu6PENW29SsEm/VlRUoGfPns2OvTnTenOGKIqoqqpi2t3JcRyCgoJUc1PlDSj4iFcSRRFms9ktm//yPA+9Xq/YD6rS0lIcPnwYJ0+exO23346UlBQAwNmzZ3HzzTc3ux+lzWZDZmZmk67J/Px8BAUFNTv25urWm7PMZjPMZjOz59fr9U5NhiHuo55bW6IY0onard2A+FosFgtsNhsMBoPiwq+urg4PP/wwOnXqhIiICNhsNtjtdmg0Gvz888+4++67UVtbC39//0aPq62txbRp0xzBdvfdd+OJJ55ATEyMx1pvzuJ5ntmZeFI3OpEXCj7iVdwdehKpC1Vp4ZeVlYW2bdvimWeeQXx8fKPf27lzJ5KTk5uEHvD7GFVWVpanynQpqQubxVif0t4/SkE7txCvYjab3R56EkEQmHaRtYbdbkdRUVGTbst+/fpBp9PhhRdewLx587By5UqcP38ewO9doNL4pt1ub/Q4b//wlmbwehLP86oaL/Ym9K9CvIYgCB4/0NNisUCn08n2A6ympqbJurfc3FycPHkSbdq0weLFizFjxgxH/VqtFh999BEA4PPPP8e6desQGBiI9PR0dOnSxTGpxduDrjl6vR42m80jN05arZbG9WRMnj/NhFxF6uJkwWQyMZ2VZ7fbUVxc3Oy6t8rKSkRHRzsmlIwdOxYLFy5ETExMswumbTYbCgoKUF5ejurqanTs2BG33HILjhw5goqKCvTv3x+AMoNP6vJ0d1e5VqulLk6Zo1mdxCuoYWZeTU1No51KpICTWm8NZ01Kv0ZGRkKjufERC1EUsXr1aqxbtw5xcXFITU1Fenq6Y4KLGtCMYELBR2RPSWuxpNbb1ftN5uXl4dKlS47WW8NlAddqvZHWoTWg6kXBR2TPG3ffqK6uRn5+fpOAO3nyJNq1a9fsureWtt5I67V2158LFy6gvLwct956K7XyvAgFH5G91uy3OHv2bOzduxe1tbXo1KkTHn30UUybNs2pa12936I0c/Lqlltubi4uX76MmJiYJgEXExPTqmNyiHs03OdVWtt4LRqNBj4+PuB5Ht988w1eeuklHDx4kILPi1DwEdlrzQ77J06cQPfu3cHzPPLz8zF27Fh89tln6Nu3b4uvVVtbi3fffbfR2JvUert6z8mIiAhqvXmxGz3Zw263Iy4uDqtXr8att97KsGLSEtQhTWRN+uBxVsNF2tKHV0FBgVPBJ63LGj9+PLXeFE46veOPaDQaLFiwAMuXL6fg8yLU4iOy5orxvf/3//4fNmzYgLq6OvTp0wfbtm1DQECAU9eiXfbJ1WpqatCtWzccOHAAUVFRrMshN4CCj8ia1Wp1yfo9m82G3377DT/99BMWLFgAnU7n1HUMBoPTjyXK9eSTT8Jut+PNN99kXQq5ATQIQWTNVfdlPj4+GDRoEEpLS7F69Wrm9RBlmT9/Pj788ENUV1ezLoXcAAo+Imuunikn7VziLJq5R5oTGRmJYcOGYc2aNaxLITeAgo/IliiKqKiocHp7qYqKCmzevBk1NTWw2Wz47rvvsHnzZgwZMsTpmij4yLU89thjWLFiBWw2G+tSyB+gMT4iG+Xl5Thw4AAOHjzo+PL19cXhw4edmlBy8eJFTJ8+HdnZ2RBFEeHh4ZgzZw6mT5/udI3BwcEUfqRZoihi4MCBeO655zBu3DjW5ZDroOAjTFy8eBGHDh1qFHImkwn9+/dv9BUREcF8uzKJRqOhrcPIdW3YsAGrVq3CDz/8wLoUch0UfMTtrly5gsOHDzsC7sCBA7h06RL69evnCLibbroJUVFRzbamWrNziytdvXMLIVezWq2IiorC119/7dRaUeIZFHzEpWpqanDkyJFGXZbnzp1DcnJyo5ZcdHT0De9s4o17dRL1eu2115Cfn08TXWSMgo84ra6uDkePHm3UXVlYWIjExMRGIRcXF9eqwFDS6QxE+SorK9GzZ0/k5uaic+fOrMshzaDgUyhpf0FXsVgsyMrKahRy+fn5iI+Pb9Rd2bt3b7cs8FbDeXxEOf785z8jJCQEixcvZl0KaQYFn4JcuHABhw4dgsViwfjx452+jtVqRU5OTqMxuZycHPTs2dMRcP3790diYqLHwoBlq49ae6SlcnJykJqaisLCQrphkiEKPi8m7SCv1Wrx0ksvYcuWLejZsydqa2sxbtw4PPDAAwgODr7uNWw2G3Jzcxu15DIzM9G1a9dG3ZV9+/aFv7+/h15Z81iN9dHYHnHGyJEjMXnyZDz00EOsSyFX8Zrgu9FjQtRq9+7diImJQWhoKPbs2YP169dj6tSpGDJkiOPvSxRFlJSUYN++fY6QO3LkCEJCQhqFXHJysmyn7dfV1cFisXjs+Xieh5+fn8eejyjHjh078NRTT+Ho0aP0GSUzsr2NbXgwpCAI1+3iko4Q4XleNUFot9tx+PBh7Ny5EzabDZMnT0ZoaCgAoE2bNrh48SIiIiIA/G+3EYvFgg8//BDHjh1D//798dJLL6Ffv35o27Yts9fRUnq9HjabzendXFpCq9VSNxVx2vDhw7FlyxZUVVX9Yc8L8SzZtfhEUYTFYoHFYnFqPIfjOPA8D57nFR2AP/74IxYtWoThw4cjISEBt912Gzp16gQAWLduHb788kts2rSJcZXuIYoiTCaTW8NPq9XCYDAo+j1E3K++vh6+vr6syyBXkVWLTxAEmEymVk1gEEURZrMZFosFBoNBsWMzp06dQo8ePfDII4+gffv2ju9fvnwZ77//Pl5//XWG1bkXx3EwGAyOf2dX43keer2eQo+0GoWePMmixdcwrFzN2z/EpK29UlNTGy0TKC8vx7Rp09CtWzcYjUb07dsXCxcuxDvvvIOamhq88MILDKv2HFfcLEmkQFXqzRKRB7vd7ti8wdXLjsiNYR581G31Pw239pJ2PpG29vr4448RFhbW7OMyMjLw4YcfYtCgQdiyZQtuvfVW+Pr6orS0FI888ghSUlI8/Eo8y9nucemGy9/fXxXd44S93NxcfPnllzAajXjllVfg4+PDuiRVYhp8ngg9idzCr6amptH+lS3d2stiseDUqVPIy8vDvn37EBgYCL1ejyVLlmDatGmIjIzE0KFDkZycfMNbg3m7hhOibDYb7Hb7Nf+sRqMBx3GYMmUK3nrrLcTExHiwUqImFRUVsFqtCA0NxahRoxAdHY2KigoMHjwY8+fPhyAI1MvgYUyDTy1T0290a6/4+PgbvgO02Wx49913sWvXLgwcOBDjxo1DXFycbIJdDm5kCczChQuh0+nw2muvMa6WKNWrr76KwMBAzJ8/H8899xzGjh2Ljh074v7778evv/7KujxVYhZ8Sl2MfPXWXgcOHMDJkyc9trUXaZnjx49j+PDhKCoqom4n4hZ79uzB4sWLsWPHDtx333149NFHMXToUKxduxZ33nkn8vPzodFocPvtt7MuVTWYtK+lLk4WTCaTy7afkrb2angSgbS1l7St15w5czy6tRdpmd69eyM8PBzffvstRo0axbocokC33347RFHE66+/jurqaseavilTpsDX1xfnzp3DjBkzkJ2dzbhS9WASfM6u0XMFaSJES4PoRrb2mjZtmiy29iItM2PGDKxZs4aCj7jNxo0bsW/fPkRFRSEhIQEAHKG3detW5OTk4IcffsDQoUMZV6oOHu/q9IYjZux2O06dOtUo5Jrb2qtfv34IDAz0cPXE1a5cuYJu3brh9OnTjdZEEuJKK1aswKOPPgoA2L9/P9555x388ssvGDZsGGbPno0BAwYwrlA9PN7ikyYbsNRwc2dRFFFYWNhoCcGhQ4fQtm1bR8B549Ze5Ma1adMGo0ePxvr16zF//nzW5RCFqqioQGJiIsLCwnD27FlMmjQJS5cudSxTojV9nuPxFp/JZILVanXqsZcvX8b8+fPxww8/oF27dnjxxRcxadKkFl9HFEWcOHECL7/8Mg4ePAg/P79GLbmUlBR07NjRqRqJd9q1axeeeuopHD58mHUpRKEuX76MkJAQ/OMf/3C0/AgbHm/xtWbN3sKFC+Hr64u8vDxkZWXhvvvuQ0JCAuLj41t0HY7jEBYW5ljcLW3uTNQrNTUVlZWVOHr0KPr27cu6HKJAbdu2xeTJk5vMHj579iz27duHjh07om3btkhISKB5Am7m0RafKIowGo1OPdZkMiEqKgr79+9Hz549AQBz5sxBly5dnD7lODg4mLoWiMOLL74Io9GIt99+m3UpRKFKSkpw/PhxjBgxAjk5Ofjggw9w4sQJ9O7dG5WVlSgtLcWAAQPwyiuvsC5V0TwafK1Zu5eZmYkRI0bg/Pnzju+98847+Omnn/Dpp586dU06YJQ0dObMGdx8880oLS2lzYWJW126dAl///vf0a5dO9x1112IjIxEhw4dYDabERMTg8LCQnoPupFH97JqTcbW1NQ0ORw1KCioVYvgWU+yIfLSvXt3JCQkICMjg3UpRMHq6+uxdOlSGI1G3H///bjpppsQEhICrVaL7Oxs3HPPPR7ZxlHNvCb4AgICUF1d3eh71dXVCAgIYFIPUaYZM2Zg9erVrMsgCubr64uDBw/i4YcfRlRUFDiOw9GjRzFv3jzcf//9GDBgAI3xuZlH+/laM57Wo0cPCIKA06dPo0ePHgCA7OxsxMXFMamHKFNaWhoWLFiAc+fO0aQn4jbDhw/HCy+8gJEjR2Lnzp3QaDS46aab8MMPP6Br166sy1M8rxnjA4CZM2eC4zisWLECWVlZSE9Px44dO1o8q1NCY3ykObNmzULPnj2xaNEi1qUQhRJFET/99BM2btyIAQMG4NZbb0XHjh1hMBhoPZ8HeM2sTuD3dTCPPPIIdu/ejbZt2+Kll15yah2fhGZ1kubs378fM2fOxIkTJ+j9QYgCeXwBu9FolMXYmkajaTJZhhDg9xu0+Ph4rF69GoMHD2ZdDlE4auF5nsdPKJVL1yIdQUOuheM4PPTQQ1izZg3rUogKcBwHQRBk0SBQC4+3+Fidw3c1Gt8j13Pu3DkkJCSguLgYBoOBdTlE4Xbv3g2z2YyRI0eyLkUVPN7ia3j6NSvSKdyEXEtoaCgGDRqELVu2sC6FqMCVK1ec3oGKtJzHg4/jOPA87+mnbYTneebhS+SP1vQRTxk7diwuXLiAX375hXUpquDx4APYBo8cgpd4h7FjxyI7OxtnzpxhXQpROB8fHzz66KNYtmwZ61JUgUnwcRzHbNzEYDBQa4/cEJ7nMWfOHJw6dYp1KUQFZs6ciZ07d6KoqIh1KYrn8cktDdXV1cFisXjs+Xieh5+fn8eej3g/u90Oi8VC7xviEY8//jh8fX3xxhtvsC5F0ZgGnyiKMJlMHtmQVavVUmuPECJrZ86cwYABA3D27FmaTexGTLo6JVKXp7uXFVDoEUK8Qffu3TFkyBCsXbuWdSmKxrTFJxFFEWaz2S3dnjzPQ6/XU+gRl6iursauXbsQHByM1NRU1uUQBdq7dy9mzZqFEydOQKNh2jZRLFms4OY4Dn5+ftDpdDCZTC7ZwcBTrUmifEajEfX19ejYsSOeeOIJ1NTUoLi4GPX19Rg5ciQEQaD3GXGZ2267DQaDAd988w1Gjx7NuhxFktXthFarRVBQkFMtNLvdDuD3wNPr9QgKCqIPI+ISK1euxHfffQcA0Ol0mD17Np555hksX74cgHy24SPKwHEcHnvsMVra4EayCj6gcXAFBARAp9PdUHP/yy+/BMdxTgcnIdfStWtXxw4uVVVVCA4OxogRIzBw4EAAwKlTp3DlyhWGFRKlue+++5CdnY3s7GzWpSiS7IJPwnGcY1JKUFAQgoODERAQAIPBAH9/fxgMBgQEBCA4OBht2rTBpk2bsG3bNgo84nJTp05FZmYmPvzwQ5SWliI4OBiCICAoKAgLFizAww8/jDfffJN1mURBeJ7HvHnz8Pbbb7MuRZFkG3xXk4JQp9PB19cXOp0OWq3WEXRpaWnYvHkz4yqJEnEch61bt6K0tBRjx45Fly5dsHv3bnz66acYNGgQ5s6diy+//JJ1mURh5syZg02bNqGiooJ1KYrjNcH3R+655x7s2rULJpOJdSlEgdq2bYvffvsNTzzxBPz9/TFy5EjU1tYiLi4OkydPRmhoKD777DPWZRIF6dSpE9LS0vDee++xLkVxFBN87du3x4ABA/DNN9+wLoUoUEhICCorK3H48GHH9wYMGICjR48CABYtWoSAgABG1RGlWrBgAf71r3+hvr6edSmKopjgA37v7ty0aRPrMohCjR8/HsuXL8fJkyexd+9eXLlyxTHBZejQoRg1ahTjConSJCYmolevXti4cSPrUhRFFgvYXaWsrAxxcXEoKyuDXq9nXQ5RmKqqKixatAhGoxHV1dUYMmQInnzySdZlEYX7+uuv8dJLL+HgwYONJu+JogibzQZRFCGKIjiOc5w1SpP8rk9RwQcAt99+OxYuXIixY8eyLoUoVG5uLjQaDWJiYliXQlTAbrcjLi4OH3zwAQYNGgSLxQJBEK670Yc0GZDneQrCZigu+FasWIHDhw/jww8/ZF0KUQHpTpsQdxFFEbt27UL37t3Rvn37Fj9eOoOUDuD+H8UFX0lJCZKSknD+/Hn4+vqyLocQQpwmCAJt4+gGiprcAgDh4eGIjo7GDz/8wLoUomDSkVrSVnmEuJIoiqirq0NNTY1LQk+6Zk1NDerq6lx2TW+luOADaDE7cT+O47Bw4UKabUdcTrqpctch3RaLxWWtSG+l2OD74osvYLPZWJdCFGzIkCFYs2YN6zKIgnjqcG5XdqF6I0UGX/fu3REeHo59+/axLoUo2Pjx43Hw4EEUFRWxLoUohNlsdnvoSQRBgNls9shzyY0igw+gxezE/fz8/JCeno6PPvqIdSlEAQRBcFv35rVISyPURnGzOiW5ublITU1FSUkJnWJM3ObAgQOYPHkyTp48Se8z4jRRFFFVVcWk61E6zk1NSx0U+5MaFxeHdu3a4ZdffmFdClGw/v37w8/Pj7rVSatYLBZm422iKHq8pcmaYoMPoNmdxP04jsOMGTNokgtxmhyCh2XwsqDYrk4AyMzMxLhx41BQUKCqZjzxrPLycsTGxqK4uBiBgYGsyyFeRhAE1NTUsC4DAQEBqlncrugWX2JiInQ6HQ4dOsS6FKJgnTt3xh133IHPP/+cdSnECznb2gsPD2/01b59ezz11FMer8MbKTr4OI5DWloa9uzZw7oUonAzZszA6tWrWZdBvJCzsypLSkocX3l5efDz88M999zjdB1qWves6K5O4Pe7GLvdDr1eT92dxG2sVqtj7Sid2kBulCiKMBqNrb7Ohg0b8MYbb+DIkSOt+pwLDg5Wxeekolt8AMDzPPz8/FTxj0nY0el0eOCBB+hUENIirmplbdiwAZMnT27155xaWn2Kb/ER4inZ2dkYOXIkzp49Cx8fH9blEC9gtVphMpladY3i4mL07dsXhw8fRteuXVt1LYPBAJ1O16preAPFt/gI8ZSEhASEhobi22+/ZV0K8RKuaHd89tlnGDhwYKtDz1X1eAPVBZ/dbkdubq4spg8T5aE1faQlXDEE8+mnn2LKlCkuqMY19XgDVXR11tbWwmKxoG3btnjrrbdw5MgRnDt3DitXrkR8fDydok1c5sqVK+jWrRtOnz7t1GnZRF1qa2tRV1fndNf4r7/+igkTJiA3N9cla0jVspZPFS2+FStW4LfffgMA/Pzzzxg5ciRGjRqF999/HwDoMFHiMm3atMHo0aOxfv161qUQGRFFESUlJdi+fTtef/11TJkyBb1790Z4eHiruhc3bNiAMWPGuGzjBLWMTSs/2vF783379u0YMWIE2rZti4EDB8JiseC1114DoJ5/bOIZM2bMwKJFizB//nzWpRAG6urqcPz4cWRmZiIzMxPHjh1DZmYmfHx8kJSUhKSkJIwaNQqLFi1CfHw8zGaz0+G3fPlyl9Wt0WhU0/OliuCbMmUKRo8eje3bt6OsrAwhISHw8fFB//79sWzZMowePRo9evRQRROfuF9qaiouXryIY8eOISkpiXU5xE1EUURxcXGTgCssLERMTAySkpLQp08fjB49GklJSejcuXOz1xEEAVar1cPVN6WmBoAqxvgAYP369di3bx/OnTuH4OBg5OTkoKamBoMGDYIgCAgJCcE//vEP1mUShXjxxRdRVVXl0jtywk5tbS2ys7MbBVxmZib0ej369OnjCLk+ffogLi4Ovr6+N3xt2qvT81QTfACwZcsWrFy5EgMGDMCUKVMcd+NFRUUYNWoUjh8/zrhCohRnzpzBzTffjNLS0hZ9CBK2RFHE2bNnm7TiiouLERsb2yjg+vTpg06dOrnkOVmdxSdR25l86oh3AFVVVVi9ejX+/ve/Y9CgQY7vm81mHDx4EGFhYY5uUEJaq3v37ujduzcyMjKQlpbGuhyvdvbsWeTk5KBv377o0qWLy65rMpmQlZXVpBUXEBDgaMXde++9WLx4MWJiYty2sJvjOPA8D7PZ7Jbr3wie51UTeoCKgi8oKAjFxcUQRRF1dXXw9fVFWVkZduzYgb179+LFF19ESEgILW0gLiOt6WsYfKIowmazQRRFx3uN4zj4+PjQ+w6///3Y7XZoNBrU1tZiwYIFOHLkCIKCgtC3b19MmzYNycnJLfo5tdvtKCwsbNKKKy0tRXx8vKMVl5aWhj59+qBDhw5ufpVN8TzP7Ew8KXjVRFVdnf/85z/x1VdfISYmBhUVFcjOzkbPnj0xc+ZMjB49mj58iEuZTCaEh4fjxIkTCAwMhCAI1/1g4zgOWq0WPM/TexG/j8tnZGRgw4YNAIC5c+eid+/eeOSRR64ZfNXV1U1acVlZWQgODm40FpeUlITo6GhZjWmxGutT09ieRFXBJwgCqqqqsG/fPthsNqSmpqJNmzasyyIKJJ2qff78eRgMhhaP80l34WrrgsrNzcXu3bshiiJqamrw+eefY82aNcjMzMS3336L5cuXIzg42PHnRVHEjz/+iLfeeguZmZkoKytDr169mozFtWvXjuGrunF1dXUePRdP2sRfbVQVfM2RFq9rNKpYy088QBAEmEwml3RbcRwHg8GgijvyM2fOYNq0aejfvz/i4+Nxxx13IC8vD9OnT0fv3r3h7++PWbNmYdKkSY0eV1BQgAMHDqBPnz6Ijo726mn5oijCZDI5fUZfS2i1WhgMBlXdWEmU/9N0DdKHEgUecRVRFGE2m116xy61fHieV/yZkllZWYiIiMBTTz2F0NBQAL8v0M7IyMCtt96KnJwcTJs2DTfddBO6devmeFxUVBSioqIYVe1a0o2Ou8NPzaEHqGTLsuaUl5fjyJEjrMsgCiHdqburm8pisbisFSkHV65cQX19faPvDRs2DLW1tXj22Wfx4IMPYv369fjll18cvx8XF4f27dvj1KlTni7Xo6Twc9eEE57nVR16gIpbfN9//z3Wr1+Pr7/+mnUpxMt5qntK6kL1pg8tm82GkydPNplReenSJezevRspKSmOPxsQEIAvv/wSALBx40bs3bsXJ06cQEZGBtq1a4eNGzciOjoaycnJrF6Ox3AcBz8/P+h0Ouo2dwPVjvFVVVUhPDwcxcXFjQbLCWkpmpDwu0uXLjUJuJycHISEhDSZURkVFdVkmMFqteL06dMoLCzEoUOHcOHCBUycOBHff/89tm3bhuTkZEydOhVDhgxh9ArZkCZKtXS5g7QsRK0Tpa5HtcEHAGPGjMGUKVMwdepU1qUQL6XGKeiCICA/P79JyBmNRiQmJjaaUZmYmHjDJwfY7Xa8//77+OKLL5CcnIyRI0eqLuSuR1oDarFYYLPZrnuqjEajwb59+6DRaDBy5EgKvKuoOvjWrFmDjIwMbNmyhXUpxAux3GrKU1tMXbx4sUnAnThxAmFhYU1acV27dqXJYh70R5shrFu3Dps2bXJ0H5P/UXXwXbp0CVFRUTh37hwMBgPrcoiXMZvNTLeZ0uv10Ov1LrmW1WpFXl5ek5CrqalpsglzQkICAgICXPK8xH2MRiMiIyNRVFREwzlXUXXwAcDw4cMxe/ZsTJw4kXUpxIt488bCFy5caBJweXl5iIiIaNKKi4yMpG4yL3bPPfdg4sSJePDBB1mXIiuqD7733nsPu3fvdmyLRMiN8IajZOrr65Gbm9tkE2az2dxsK87f39/D1RN3+/jjj/HZZ58hIyODdSmyovrgKy8vR2xsLMrKylzWbUSUz2Qytfrw0NOnT+OWW27BuHHjsGrVKqeuodPpYDAYUFZW1qQVl5+fj27dujUKuKSkJISHh1MrTiWk2etFRUW0PWMDql/Q0blzZyQlJeHbb7/FuHHjWJdDvIQr1uw9+eSTrV6TVlpaittuuw1WqxVJSUlISkrC0KFDsWDBAvTu3VuWyx6I5wQFBSE1NRVfffUVpk2bxroc2VB98AHAxIkTsXnzZgo+ckOkWXStsXnzZgQHB2PAgAE4c+aM09fp0KEDDh8+jNDQUGrFkWalp6dj/fr1FHwN0NxjABMmTEBGRkaTLZQIaY7NZmvV46uqqvDaa6/hlVdeaXUtGo0GnTt3ptAj1zR27Fjs27cPly9fZl2KbFDwAQgLC0NsbCy+//571qUQL9Da1t6rr76KBx54AOHh4bKohyhbYGAghg0bRuv5GqDg+z9paWnYvHkz6zKIF2hN0GRlZWHPnj2YN2+eLOoh6jBp0iRs3LiRdRmyofpZnZKCggLcfPPNOHfuHG3iSq7LarXCZDI59dh///vf+Nvf/uZYAG4ymWCz2RATE4M9e/Y4dU2DwQCdTufUY4k6VFdXIzw8HAUFBV5zKK87UfA1kJKSgqVLl2Lo0KGsSyEyI4oiiouLkZWVhfLycowaNcqpGZO1tbWorq52/P/KlStRVFSEN998Ex06dHCqNpb7dhLvkZaWhjFjxmDGjBmsS2GOfloakLo7KfjUzWg0Iisrq8mXn58fEhMT0a9fP4wfP96pa/v7+zdaKC6du+Zs6AHw6hPHieekp6fjww8/pOADtfgaycvLw9ChQ1FSUkKb7apAw/0pGwZcZWUlevfujcTERMcJA4mJiY3CyWg0ymJsTaPRICgoiHUZxAvU1NQgLCwMZ86cQfv27VmXwxQF31USEhLw3nvv4ZZbbmFdCnERURRRUlKCrKysRiF38uRJdO3a1RFsUtA1d1bc1Vyxc4srSDu3EHIjJk2ahJEjR+Lhhx9mXQpTFHxXWbx4MaqqqvDWW2+xLoU4wWg0Ijs72xFumZmZyM7OBs/zjVpviYmJ6NWrl9M7m3jDXp2EXG3jxo344IMPsGPHDtalMEXBd5WsrCyMGTMGhYWFtCjYSeXl5SgsLESvXr1u+BDSlpK6KRt2UWZmZqKyshK9evVqEnIdO3Z06fN78+kMRL1MJhNCQ0Nx+vTpVo0rezu6VbxKQkICeJ7HwYMHcdNNN/3hYY9qJ/39aLVafPnll1iyZAnMZjMCAgIwatQoPPjggwgLC3P83Tlz/dLS0ibjcPn5+YiMjHR0T86cOROJiYno3r27R8ZnOY4Dz/NMz+PjeZ7eg6RFDAYD5s2bh/LyclUHH7X4riKKIv79738jKSkJCQkJ172j5zgOWq0WPM9TEAK4fPky6urqEBoaisrKSsyfPx933303pk6dekPBV1VVhezs7CYhx/N8ozG41nZTuooaTmAnymO1Wh2fXWql3ld+FVEUYbFYYLFYMHnyZMf3/ugxVqvV8UbieV51d+H5+fnYv38/SkpKcM899yAxMRFmsxnff/89eJ7HsGHDAKDR34koisjPz8eRI0cahdzFixfRq1cvR8ilpaW5pZvSVTiOg8FgYDLWZzAYVPU+I65Dmx1Q8AH4faKCyWRyBJ0zXWWiKMJsNsNiscBgMKjibqq4uBgPP/wwkpOT0bNnT/j5+UEQBMyaNQtZWVno3LkzLl++jJCQkEaPMxqNWLp0KS5fvozExETMmDHD0U3pbWvSpBa/xWLx2HPyPK+K9xfxPGeHJLyNqrs6G4aVq/E8D71er+g30ddff42PP/4YS5cubXbD5U2bNiEjIwOvvvoqwsLCHN9X2g+XKIowmUwuOaPvj2i1WmrtEZcSBAHHjx9HaGiobHtXXE21q7SlDyt33albLJZGrUilaPh6UlNTYbFY8Oyzz2LmzJlYtmwZKisrHb+fkpKCqqoqFBUVNbqG0j60pS5Pd7fCKPSIq5hMJly6dAkAsGTJEixZsgRpaWnIzMwEANjtdpbluZ0q+0s8dYcudaF644eVKIo4d+5co0XfZWVl+OqrrxyTSvz9/bF161YAwNatW7F27VoEBARg+vTpuHjxIlatWoW4uDj07duX4SvxDCn8qAeBeIN33nkHffr0wejRo5GZmYm0tDSUl5dj1apVWLlypeJu2K+myuAzm80e6ZYCfg8/s9nMfAbi9VRXVzsWfTecbKLT6RwTTYYOHYrExETo9XrH46xWK86cOYOzZ88iPz8fUVFRAIDExES0b98evXr1wrx582T92l2J4zj4+flBp9O5rLXvqdYkURdfX19s27YNo0ePRlBQEJKTk8FxnONMUm8ba28p1Y3xsdpxQw47bAiCgPz8/EYLvrOysnDhwoVGsymlr86dO1/3ena7He+//z62bNmCpKQkjBkzBoMHD6alHWg8S7glP2J2ux319fXw8/NT5Sxh4hllZWW44447sHTpUixZsgSffvopQkJCcODAAXTp0gV79uxBUlIS+vTpw7pUt1BV8Kll3VXDbsqGX3l5eQgLC2uyq0mPHj0Uf4fHirTA32KxwGazXXfsRKPRwGaz4b777sOGDRtUv5Ewca/t27cjIyMDXbp0wWOPPYbAwEDH59M///lP/Pe//8XXX3/NuEr3UFXwmc1mpjtt6PX6Rl2FrtCwm7Lhl4+PT7OLvmlDY7ZuZCeg+++/H4MGDcL8+fMZV0uUrKKiAr/++ivGjBnj+F5hYSE2bNiATz75BBUVFTh//rwiT6pRTfB5+96KgiDg5MmTTXY1KS8vR3x8fJOQ+6NuSiJfu3btwpNPPokjR46wLoUoWGVlJYYNG4Yvv/wSn3zyCb766itUVlbipptuwpQpUzB06FAEBASwLtMtVBN83rKbviiKOH/+fJNxOKmb8uojdKibUnnsdjuioqLwxRdfIDk5mXU5RMEGDBiAnJwc3HfffRg9ejRGjhypil4h1QRfa85PW7VqFTZs2ICcnBykpaXhX//6l9N1NDw/raampskROllZWdBoNE3G4Xr37q2KNyT53eLFi3Hp0iWsWLGCdSlEwZYvX469e/diy5Ytju8pbYOJ5qgm+FpzYnZGRgY0Gg2+++47mM3mVgWf0WjEggULkJmZibKyMsTHxzcJuc6dOyv+jUeur7CwEP3790dJSYnLx4UJkVRWVuKZZ57BqlWrGn1fFEXU1NTg4MGD6NixIxISEhhV6B6qCD5RFGE0Glt9nb/97W84d+5cq4LPbrdj9+7d6N27N3r27EndlOSa7rzzTsyePRvp6emsSyEKJggCtFqt49eGHn/8cWRnZ2Pnzp2MqnMP5U3XaYbNZmNdgoNGo8G4ceMQGxtLoUeua+bMmVi9ejXrMojCaTQaZGdnO0LPaDRi//79eP/99yGKIq5cuSKrz1BXUMV2EHJr1MqtHiJP9957L+bPn4/i4mJERESwLoco1OXLlzF48GDs2LED2dnZOHPmDGpra+Hn54fo6GiMHTsWdrtdUTfqFHwMyK0eIk9+fn5IT0/H2rVr8fzzz7MuhyhU+/bt0bdvX7z00ksYMmQIQkJC0KtXL8THxzd76ooSqCL45DZRRG71EPmaOXMmJk+ejGeffVaRC4mJPNx3331YuXIlnnvuuSatOyXO8lTFT1Jr/9GkjaZtNhtsNlurN7lW2puIuE///v3h7++PvXv3si6FKNh9992H/v37O3YRknYWApT5eUWzOm/A66+/jjfeeKPR9xYtWoSnn37aqesFBwcr8s1E3GPZsmU4evQo1q5dy7oUQhRBFcEHtG4dnytpNBoEBQWxLoN4kYqKCkRHR6OoqIjeO4S4gCq6OgEwPxJIoqSZUcQzOnbsiGHDhuGzzz5jXQpRAbPZrPgT2FUTfDzPsy4BgHzqIN6F1vQRT9m+fbviu9VVE3xyOBxVGjgmpKVGjBiBs2fPIicnh3UpROH8/f3x/vvvsy7DrVQTfBzHMW9t0WnaxFlarRbTp0/HmjVrWJdCFG7YsGHIz89HUVER61LcRjXBB7ANHjkEL/FuM2bMwLp165w+ZYSQG6HT6TB+/Hhs2rSJdSluo6rg4ziO2dE+BoOBWnukVWJiYtCzZ09s376ddSlE4dLT07Fx40bWZbiNqoIP+L3LyNMtL57nZTOrlHg3muRCPGHo0KE4deoUCgsLWZfiFqoLPgDQ6/UeCyKtVkvnqRGXmTRpEvbu3YuysjLWpRAF0+l0uPfeexXb3anK4JO6PN0dflqtlro4iUsFBgZi7dq1tJCduF16ejo+//xz1mW4hWp2bmmOKIowm82wWCwuvzbP89Dr9RR6xOVsNhstiyFuJwgCunTpggMHDqBbt26sy3EpVbb4JBzHwc/PDwEBAS4LKI7jEBAQAD8/Pwo94hYUesQTtFotJkyYoMhWn6qDT6LVahEUFNSqFlppaSlMJhOCgoJoIgshRBGUOrtT1V2dzRFFETabDRaLBTab7bp71mk0Gvj4+IDneTz33HOor6/HsmXLPFgtIUBVVRV27tyJoKAg3HXXXazLIQoiCAJCQ0Pxyy+/oHv37qzLcRkKvj8gBaF0PhXHcY6txxq2DouKitC3b18UFhbSxAPidkajEfX19ejYsSP+9Kc/wWQyoaSkBM8++yxGjRoFu91OB9cSl/jzn/+MqKgoLFq0iHUpLkM/GX+A4zhotVrodDr4+vpCp9NBq9U26RKNjIzEXXfdhQ8//JBNoURVVq5ciV27dgEAfH19MWfOHDzzzDOOHgel765PPEeJszsp+Fzosccew4oVK2Cz2ViXQhSuW7du2Lp1K0RRRGVlJQwGA0aPHo3U1FQIgkDjzMRlhgwZguLiYpw+fZp1KS5DwedCAwcORPv27bFt2zbWpRCFmzp1KnJycrB27VpwHIeYmBiUlJQgLy8PTz/9NLZv346zZ8+yLpMogFarRVpamqJafRR8LsRxHB577DEsX76cdSlEBf7zn//g/PnzuHDhAoYNG4YePXqgrKwMt956K7Zs2YK5c+eyLpEohNJmd9LkFhezWq2IiorCtm3bkJSUxLoconDnz5/Hu+++i65du8JisaCwsBALFixAaGgohg8fjgceeADTpk1jXSbxcjabDWFhYdi3bx+io6NZl9Nq1OJzMZ1Oh7/85S94++23WZdCVGDp0qXw8/PDzJkzMXfuXNTX12PlypUAgLlz59KensQlfHx8mnR3iqIIQRBgtVpRX18Pq9UKQRDgDW0pGgF3g9mzZ6Nnz554/fXX0alTJ9blEAXr1q0bjhw54vh/Hx8fDBw4EACQmpqKQYMGsSqNKEx6ejrWrFkDk8n0hwEnzYbneb7J0i85oK5ON5k9ezbCw8Px4osvsi6FKJjRaER6ejoiIiLw7bffokePHvj888/RoUMH1qURhRBFERaLBRaLBTU1NfD3929RkEmHcLM8CPxqFHxucvz4cdx5550oLCykk9eJW5lMJvz3v/9FaGgoBg8ezLocoiCCIMBkMrmk+9JTp+LcUC0UfO4zfPhwPPjgg3jwwQdZl0JURNphiBBnKf3kGprc4kaPPfYYli1b5hWDvUQ5KPRIa4iiCJPJ5JbQAwCLxeKyVqSzKPjcaOTIkTCZTPjxxx9Zl0JUwG634/jx47BaraxLIV5KCj1BENz6PK7sQnUGBZ8baTQaLFiwgBa0E4+Q3m9bt25lXQrxUmaz2e2hJxEEAWaz2SPPdTUKPjebNm0a9uzZg4KCAtalEBWYOXMmVq9ezboM4oUEQXBb9+a1WCwWjwVtQzS5xQOeeuop2Gw2vPnmm6xLIQpXV1eH8PBwHD16FBEREazLIV5CFEVUVVUx6XrkOA5BQUEeHZum4POAs2fPol+/figsLERgYCDrcojCzZs3D6GhoXj++edZl0K8hNlsZtbtCAB6vR56vd5jz0ddnR7QtWtXpKam0ll9xCNmzpyJNWvW0Jl85IZIC9RZslgsHm1tUvB5iHRWH30YEXdLSUmBv78/9u7dy7oU4gVsNhvzJVeiKHr0HFMKPg8ZPHgwOnXqhMzMTNalEIXjOI4muZAb1prWXlFRESZNmoRu3bohNjYWTz75pNOTVTzZ6qQxPg+qr6+HIAjw9/dnXQpRuIqKCkRHR+Ps2bMIDg5mXQ6RMaPR6HSLb9KkSejYsSPeeustGI1G3HvvvZg+fTrmzJnT4mtpNBoEBQU5VUeLn8sjz0IAAL6+vhR6xCM6duyIYcOGKerwUOJ6oii2qpuzqKgI48ePh16vR+fOnTFs2DCcOHHCqWvZ7XaPdblS8BGiUNTdSf5Ia8fV5syZgy1btqC2thbnzp3Drl27cOeddzKr50ZR8DHWcLIL9ToTVxoxYgTOnj2LnJwc1qUQmWrtZ84tt9yC3NxcREZGonfv3khOTsbdd9/NrJ4bRcHHUF5eHpYuXYrnn38eNpuNNhcmLqXVajF9+nSsWbOGdSlEploTNHa7HWlpaRgzZgxKS0tx+vRpXLlyBS+99BKTelqCgs/DKioqcO7cOQDA448/jqKiIpw8eRL/+te/AIDJ9j1EuWbMmIF169bRxtWkWa252b58+TJKS0sxa9Ys8DyPdu3aYerUqdi5cyeTelqCgs/D3n//fccmwn379sXUqVPxt7/9DR9//DEAyOKQRqIcMTExiI6Oxvbt21mXQmREEARkZmZi27ZtTu/Y0r59e3Tt2hWrV6+GIAgwGo3YsGEDEhISnK6Lgk+hBg8ejM2bN6O+vh4nTpyAxWJBdHQ05s6di9LSUuzZswf79u1jXSZRkBkzZtAkFxWTQm7NmjV45JFHMGjQILRp0wbp6en44osv4Ovr6/S1161bh++++w49e/ZEv379oNVq8eqrrzp9PR8fH6cf2xK0jo+B22+/HXfeeSf27NmDJUuWoF+/frBYLOB5Hr/99htmzZqFY8eOsS6TKER1dTUiIiKQm5uLkJAQ1uUQNxIEASdOnMDBgwdx6NAhHDp0CFlZWQgPD0dKSorjKzk52bFmrjXr+FzJk+v4qF+NgY0bN2Lfvn3o2rUrevfuDQDgeR5lZWXIyMhAVlYWdu/ejTvuuINtoUQRAgMDMWHCBHz88cdYuHAh63KIi0ghJwXcwYMHkZmZifDwcPTv3x8pKSmYNGlSo5BrjlarlcUYsKdaewC1+JhZsWIFHn30UQDAr7/+ihUrVmD//v244447MHv2bAwaNIhxhURJfvzxR8yePRvHjx+n2cNe6OqQO3ToEDIzMxEWFoaUlBRH0P1RyF3r2jU1NW6q/MYFBAR4bI4DtfgYKS8vR9++fdGlSxcUFBRg4sSJeP311x1nqImiSB9QxGVuueUWCIKAX3/9FQMHDmRdjldz98/mH4VcSkoKJk6c6FTINcfHxwccxzHt7uQ4jlp8anDp0iWEhITgzTffxCOPPOL4QaLAI+7y+uuv48yZM1i1ahXrUrxOVVUVDh48iCtXrmDChAkuu64gCMjNzW00Jnd1yEktOXfuuaq28/ioxcdIu3btkJ6eDp7nwXEcbDab486LEHeYNm0aEhISsGzZMhgMBtblyJooirDb7fDx8UFJSQlmzZoFURRhMBhQXV2N8ePHtziIpJCTxuOkkAsNDXV0Vaalpbk95JrD87zHz8STcBwHnuc9+5zU4mOnqKgI+fn5rdrbjpCWGDNmDNLT0zFt2jTH96Sz0KQNizmOc3Q90Y3Y72dp+vn54bXXXkNRURFWrFiBvn374oEHHrhmD03DkGs48SQ0NLTJmJxcTs9gNdbnybE9CbX4GIqMjERkZCTrMoiKzJw5EytWrMD9998Pi8UCQRCue5fPcRy0Wi14nldNEIqiiCNHjmDHjh1o164dfHx8UFFRAQDo3Lkzzp07h6qqqibBV1JSgiVLluDQoUM4duyYI+RSUlIwYcIEWYVcc6R/Z0+ei8fzPJNNOyj4ZMJqtaK2tlbWPxjEu4miiBEjRiAuLg7V1dU3FGKiKMJqtcJqtTq6pKTueaX67bff8Pjjj2Po0KFITk5GREQEnnjiCXz33XcoKioCx3HIzc0F8PvaM0lAQAAiIyNx7733ol+/fl75s6zX62Gz2TyydaJWq/XouF5D1NUpE1arFd27d8c333zjWNtHiKsIggCTyeSSMRyO42AwGBS7vd5HH32EnTt3YtmyZejQoQMA4JdffsGSJUvQoUMHTJgwAcuWLcOmTZsQGBjIuFrXE0URJpPJreGn1WphMBiY3UBR8MnIK6+8guLiYpp1R1xGFEWYzWa3dF/xPA+9Xq+41l9FRQUeeOABdO3aFdXV1UhOTsacOXMcLbg1a9agsLAQzzzzDLMWi7sp/X1DwScjFy5cQGxsLE6ePOm40yTEWWq4c3eGzWZrNPFk0KBBmDRpUrPryDIyMvDBBx9gypQpuPXWW/Hrr7/i448/xp///GcMHz5c8cuPlNpTwL4C4tCpUydMmDABq1atwrPPPsu6HOLFPBF6wP8+GOUafjabDXl5eY3WyR07dgwhISGOiScJCQmNxuosFgtOnTqFvLw8/Pzzz4iKikL79u3xxBNPwGg0Ij09HUOGDAHgudMEWNFqtQgKCoLFYnF6uYMcx4apxSczmZmZGDVqFAoKClq1azpRt7q6Oo/PzvPz8/PY8zVHCrmG6+SuDrmUlBT069cPbdq0ue513nvvPXz77be4+eabMWbMGCQmJnruhciUtOzFYrHAZrPBbrdf889qNBr4+PjIdjYwBZ8Mpaam4k9/+hPuv/9+1qUQL6SG9VgNQ04KumPHjqFz586N1sn9UcgR53nz+k8KPhn66quv8Morr+C3336T/RuIyIsoiqiqqmK2A0dQUJDL37NXh9yhQ4dw9OhRR8hJQUchR24UBZ8M2Ww2xMbG4qOPPsLgwYNZl0O8iLfvuWiz2ZCfn99oTO7qkJNacm3btnVh5URNKPhkasWKFfjxxx+xceNG1qUQL8GytSdpSatPCrmGY3IUcsQTKPhkqrq6Gt26dcORI0doWzNyQ+R8rlrDkJOC7ujRo+jUqVOTMTkKOeJuFHwy9sQTT0Cn0+GNN95gXQrxAiaTyemTtMeMGYODBw86AqtLly44cOCAU9fSarUoLi5uNCZ35MgRR8g1HJOjkCMsUPDJWEFBAW666SacPXuWjpEhf8hoNDrdzdncqQ3OKi4uxvjx45t0V7Zr167V1ybEFWgBu4xFRUVhyJAh+OijjzB37lzW5RAZk6aUy0F4eDhOnTpFM5KJbGn++I8QlhYsWIC33377uotFCbHZbK2+xl//+lf06NEDI0aMwI8//uj0daSDlQmRKwo+mRsyZAj8/PywY8cO1qUQGWtta2/x4sU4cuQIcnJy8NBDD2HKlCkoKChgVg8h7kTBJ3Mcx+Gxxx7D8uXLWZdCZKy1QdO/f38EBgaC53lMmTIFN998M3bu3MmsHkLciYLPC0yePBnHjh3D8ePHWZdCZMod42mtCS8a3yNyRsHnBXiex9y5c7FixQrWpRCZMJlMOHz4MNavX48XXngBzz33nNNr+IxGI7777juYzWYIgoCNGzfi559/RmpqqtP1UfAROaPlDF6ivLwccXFxOHXqFNq3b8+6HOIBoijiwoULOHHiBHJzcx2/5ubmoqKiAtHR0YiPj0dcXBx69eqFu+66y6nAuXjxItLT03Hy5EloNBrExMTg2WefxdChQ52uPTg4mMKPyBYFnxeZMWMGYmNj8fTTT7MuhbiQIAgoKChoFGzSf/v4+DjCreGvkZGRTQ5Obc06PlfSaDQICgpiXQYh10TB50WOHj2KMWPGoKCgADqdzvF9bz4eRE1qamqQl5fXpPV2+vRpdOnSxRFsDUOuQ4cON3z91uzc4ko6nY42XCCyRsHnZYYOHYo5c+Zg4sSJsFgsEAThunf5HMdBq9XK9kBIFs6fPw+r1YrIyEjHjYKriKKI8vLyZltvFy9eRExMTKNgi4uLQ0xMjEsOcZXzXp2EyAkFnxcRRREHDhxAmzZt0KlTpxY/nuM48DwPnudVE4DSwn+pFZydnY0hQ4YgNTUVmzZtgiAITn1IC4KAM2fONGm9nThxAjqdrtnWW3Pdk67kbaczEMIK3ZZ5CUEQYDKZEBMT4/Q1RFGE2WyGxWKBwWBQxV25RtN44nJpaSlSUlJQWloK4I9nH0rdk1dPMDlz5gxCQ0MdwXbLLbfg4YcfbnH3pCtJNzYsz+NT000V8V7U4pO5hmHlajzPQ6/XK/aDqri4GAcPHsSpU6cwYsQIREdHY9KkSfjPf/6DuLg4XLlypcljamtr8fTTTztCrrKyEjExMU0mmERHR7uke9LVlHgCOyGupvxbfi8miiJMJhMEQXDL9S0WC2w2GwwGg+I+rGpqajBz5kyEh4cjIiICoihix44dGDBgAEJCQpq0BCU8z6Nnz54YNWqUY/bktf6sHHEcB4PBwGSsT4nvI6JM1OKTKXeHXkNarVZxH1r79+/HihUr8Ne//tXRPTx16lQ899xz6NWrFyIiInDPPffg+eefR0hICONqXa+urs4tvQTXwvO8LFvAhDSHgk+m6IPr+kRRxPnz5x2TSnr06IG77rrL0TqzWq146KGHUF9fj7Zt2+L222/HP//5T+Tn56Njx464fPkygoODsWPHDnTr1o3ti3EDunEi5Nqoq1OGBEHwaOgBv3d76nQ62U14sVqtOHPmTLPLA3ied4y5RUVFwW63O4JPp9Phk08+AQB8+umn+OqrrzB27FgsWrQIGo0GL7/8MgoKCtCtWzfU19fD19eX5ct0OanL093hR6FHvJG8PuWI406dBZPJxGxyQlVVVZPF3SdOnEBBQQHCw8MdE0tuu+02zJo1C3Fxcdfduk3aDaWsrAxVVVXw9/fHuHHjHME4aNAgx+xHpYWeRAo/mhxFSGPU1SkzZrOZ6XR0vV4PvV7vlmtL3ZPNtd4uX76M2NjYJmvfoqOjnapHFEWsXr0an3zyCeLj43HnnXfi3nvvdcOr8g7SchhX/LhLgSq33gFCbhQFn4woZQGy1WrF6dOnmyzuzs3NhV6vb7JzSXx8PCIiIrxq9qQ3EkURFosFFovFqfeYGjdAIMpEwScj3rblVFVVVaNQk0KuoKAAERERTVpvcXFxaNeunQdeAbkeaW9XaTmLtLtNczQaDXx8fGjLO6IoFHwy4uwmwxaLBQsXLsTu3btx5coVREVF4YUXXsBdd93lVB0NNxkWRRHnzp1rMvaWm5uLK1euIDY2tknrrWfPnm7rLiWuR5ucE7WhTnoZcXb2nSAICAsLw7Zt2xAeHo5vv/0WM2fOxE8//YTIyMgWX6+yshJz5851tOT8/Pwa7VwyduxYxMXFUfekQkgbmROiFtTikwlRFGE0Gl12vVtuuQWLFi3CuHHjWvxYm82GjIwMxMbGIjY2lronCSGKQrd5MmGz2Vx2rQsXLuD06dOIi4tz6vE+Pj544IEHqBVACFEk6qeSCVc1vK1WK2bPno3Jkye3+iQHQghRIgo+mXBF0NjtdsyZMwc6nQ7/+Mc/mNdDCCFyRH1ZMtHa2XOiKGL+/PmoqKjAxo0bodPpmNZDCCFyRcEnE60NmieeeAL5+fnYunWrSzabpuAjhCgVzeqUidbM6iwqKkJSUhJ4nm80IeWtt95Cenq6U9cMDg6m8COEKBIFn4wYjUZZjK1pNBoEBQWxLoMQQtyCJrfIiFyWD/j4+LAugRBC3IaCT0Z4nmddAgD51EEIIe5AwScjctgbUdqjkRBClIqCT0akY19YoiNnCCFKR8EnMyyDRw7BSwgh7kbBJzPS6dYsGAwGau0RQhSPgk+GtFqtx1teV68BJIQQpaLgkym9Xu+xINJqtXRwLCFENSj4ZErq8nR3+Gm1WuriJISoCu3cInOiKMJsNsNisbj82jzPQ6/XU+gRQlSFgs9LCIIAk8nkki3NPNWaJIQQOaLg8yKiKMJiscBisTgVgNJyBVqrRwhRMwo+LySKImw2GywWC2w2G+x2+zX/rEajgY+PD3iel8XOMIQQwhoFnwJIQSiKIkRRBMdxjq3HKOgIIaQxCj5CCCGqQssZCCGEqAoFHyGEEFWh4COEEKIqFHyEEEJUhYKPEEKIqlDwEUIIURUKPkIIIapCwUcIIURVKPgIIYSoCgUfIYQQVaHgI4QQoioUfIQQQlTl/wNsqZLABzzARAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import networkx as nx\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "G = nx.Graph()\n",
    "G.add_nodes_from([i for i in range(9)])\n",
    "G.add_edges_from([[0,1],(1,2),(0,3),(1,4),(2,5),(3,4),(4,5),(3,6),(4,7),(5,8),(6,7),(7,8)])\n",
    "options = {'node_size': 1200,'with_labels':'True'}\n",
    "pos = nx.spring_layout(G)\n",
    "nx.draw(G, pos, **options, node_color='#efefef')\n",
    "\n",
    "edge_labels={(0,1):'s0',(1,2):'s1',(0,3):'s2',(1,4):'s3',(2,5):'s4',(3,4):'s5',(4,5):'s6',(3,6):'s7',(4,7):'s8',(5,8):'s9',(6,7):'s10',(7,8):'s11'}\n",
    "nx.draw_networkx_edge_labels(G, pos, edge_labels=edge_labels)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "JHpcOxEKLhDv"
   },
   "source": [
    "## 2 Cars and 3 possible routes for each car\n",
    "Now there are two cars, going from point 0 to point 8. We propose three routes for each car to take, and consider the cost of each route to be counted as congestion if the cars take the same street. Then assign a qubit to each proposed route from Q0 to Q5.\n",
    "\n",
    "route for car1  \n",
    "Q0 : s0,s3,s8,s11  \n",
    "Q1 : s2,s7,s10,s11  \n",
    "Q2 : s0,s1,s4,s9  \n",
    "\n",
    "route for car2  \n",
    "Q3 : s0,s3,s8,s11  \n",
    "Q4 : s2,s7,s10,s11  \n",
    "Q5 : s0,s1,s4,s9  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "uyAnT5dwL4V-"
   },
   "source": [
    "## Calculation of congestion level\n",
    "The level of congestion is determined by whether the street appears within the proposed route or not.　We can calculate the congestion for the roads from s0 to s11 and add them together to get the overall congestion level. For example, s0 appears in Q0, Q2, Q3, and Q5 of the proposed routes, so the cost is $(Q0+Q2+Q3+Q5)^2$. The whole cost is,\n",
    "\n",
    "$$\n",
    "H = (Q0+Q2+Q3+Q5)^2 + (Q2+Q5)^2 + (Q1+Q4)^2 + (Q0+Q3)^2 + (Q2+Q5)^2 + (Q1*Q4)^2 + (Q0+Q3)^2 + (Q2+Q5)^2 + (Q1+Q4)^2 + (Q0+Q1+Q3+Q4)^2\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "from blueqat.pauli import qubo_bit as q\n",
    "\n",
    "hamiltonian = (q(0)+q(2)+q(3)+q(5))**2 + (q(2)+q(5))**2 + (q(1)+q(4))**2 + (q(0)+q(3))**2 + (q(2)+q(5))**2 + (q(1)*q(4))**2 + (q(0)+q(3))**2 + (q(2)+q(5))**2 + (q(1)+q(4))**2 + (q(0)+q(1)+q(3)+q(4))**2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Hard constraints\n",
    "There is one constraint that must be observed here. That is the constraint that each car should choose only one route. This can be easily solved by using XYmixer. The XYmixer is a Hamiltonian that determines how the problem is explored; the search is performed in such a way that the 01 and 10 state vectors are swapped. Here, car 1 and car 2 choose only one of the three routes, so only one of the three qubits, Q0, Q1, Q2, will be 1 and the rest will be 0. Likewise, only one of the qubits Q3, Q4, and Q5 will be 1. To maintain this state, a mechanism for exchanging qubits, called 2opt, is used. For the initial state, choose any state that satisfies the constraints."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Quantum Alternating Operator Ansatz"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3572.3314878940582\n",
      "sample:010001\n"
     ]
    }
   ],
   "source": [
    "from blueqat import vqe, Circuit\n",
    "from blueqat.pauli import X, Y\n",
    "from blueqat.utils import qaoa\n",
    "import time\n",
    "\n",
    "step = 4\n",
    "\n",
    "#mixer and init state\n",
    "mixer = 0.5*X[0]*X[1] + 0.5*Y[0]*Y[1] + 0.5*X[1]*X[2] + 0.5*Y[1]*Y[2] + 0.5*X[2]*X[0] + 0.5*Y[2]*Y[0]\n",
    "mixer += 0.5*X[3]*X[4] + 0.5*Y[3]*Y[4] + 0.5*X[4]*X[5] + 0.5*Y[4]*Y[5] + 0.5*X[5]*X[3] + 0.5*Y[5]*Y[3]\n",
    "init = Circuit(6).x[0].x[3]\n",
    "\n",
    "start = time.time()\n",
    "result = qaoa(hamiltonian, step, init, mixer)\n",
    "b = result.circuit.run(shots=10)\n",
    "sample = b.most_common(1)[0][0]\n",
    "\n",
    "print(time.time()-start)\n",
    "print(\"sample:\"+ str(sample))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the end, the combination of either Q2 and Q4 or Q1 and Q5 resulted in the least congestion. These solutions are considered correct because they do not go through the same road"
   ]
  }
 ],
 "metadata": {
  "colab": {
   "name": "tutorial017_traffic_flow_optimization_en.ipynb",
   "provenance": [],
   "version": "0.3.2"
  },
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}