Capsule_Network_complete.ipynb

{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": 46,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "6_-rU3J1Kmhd",
        "outputId": "65db939f-3b83-45b5-d122-bf65a4b28550"
      },
      "outputs": [
        {
          "data": {
            "text/plain": [
              "<torch._C.Generator at 0x7f4968d6fab0>"
            ]
          },
          "execution_count": 46,
          "metadata": {},
          "output_type": "execute_result"
        }
      ],
      "source": [
        "# import resources\n",
        "import numpy as np\n",
        "import torch\n",
        "import time\n",
        "\n",
        "# random seed (for reproducibility)\n",
        "seed = 1\n",
        "# set random seed for numpy\n",
        "np.random.seed(seed)\n",
        "# set random seed for pytorch\n",
        "torch.manual_seed(seed)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 47,
      "metadata": {
        "collapsed": true,
        "id": "-0ZoCrZzKmhq"
      },
      "outputs": [],
      "source": [
        "from torchvision import datasets\n",
        "import torchvision.transforms as transforms\n",
        "\n",
        "# number of subprocesses to use for data loading\n",
        "num_workers = 0\n",
        "# how many samples per batch to load\n",
        "batch_size = 20\n",
        "\n",
        "# convert data to Tensors\n",
        "transform = transforms.ToTensor()\n",
        "\n",
        "# choose the training and test datasets\n",
        "train_data = datasets.MNIST(root='data', train=True,\n",
        "                            download=True, transform=transform)\n",
        "\n",
        "test_data = datasets.MNIST(root='data', train=False, \n",
        "                           download=True, transform=transform)\n",
        "\n",
        "# prepare data loaders\n",
        "train_loader = torch.utils.data.DataLoader(train_data, \n",
        "                                           batch_size=batch_size, \n",
        "                                           num_workers=num_workers)\n",
        "\n",
        "test_loader = torch.utils.data.DataLoader(test_data, \n",
        "                                          batch_size=batch_size, \n",
        "                                          num_workers=num_workers)\n",
        "\n",
        "\n",
        "primary_caps_total_time = 0\n",
        "digit_caps_total_time = 0\n",
        "dynamic_routing_total_time = 0\n",
        "decoder_total_time = 0\n",
        "total_capsule_net_time = 0\n",
        "capsule_loss_total_time = 0"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 48,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 268
        },
        "id": "A498vTDyKmhu",
        "outputId": "d38cf9c3-eee3-4621-b654-6d14e8cfbfa8"
      },
      "outputs": [
        {
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAABXEAAAD7CAYAAAAsAtcsAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3de9xNdfbA8fV1JyGUrujifk+JRjFCpUKKyKXUTKVfdCNdVCakpCuVpNSUSZJLpciESlJJmnENhYhckvs1+/fHY77zXd/pnM5znPPsfZ7n8369vH5rzTpn7/Wr3T77fO29jgmCQAAAAAAAAAAA0ZQv7AYAAAAAAAAAALGxiAsAAAAAAAAAEcYiLgAAAAAAAABEGIu4AAAAAAAAABBhLOICAAAAAAAAQISxiAsAAAAAAAAAEcYiLgAAAAAAAABEWJ5ZxDXGzDLG7DXG7Dz8Z1nYPSH6jDGljTETjTG7jDGrjTFXh90TMocxptLh887rYfeCaDPG3GKMmWeM2WeMeSXsfpA5jDHVjDEzjDHbjDErjDGXh90Tos0YU9gY89Lh65odxpgFxpiLw+4L0cbnFJJhjHndGLPeGLPdGPOdMeYvYfeE6ON8gyOR27+D55lF3MNuCYKg+OE/VcJuBhnhWRHZLyLlRKSziDxvjKkRbkvIIM+KyFdhN4GM8JOIDBSRl8NuBJnDGFNARCaLyHsiUlpEbhCR140xlUNtDFFXQER+FJEmIlJSRPqJyDhjTMUQe0L08TmFZAwWkYpBEJQQkdYiMtAYUz/knhB9nG9wJHL1d/C8togLJMwYc5SIXCEi9wdBsDMIgtki8o6IdA23M2QCY0xHEflVRD4KuxdEXxAEE4IgmCQiW8LuBRmlqoicKCJPBkHwWxAEM0TkM+FzCnEEQbArCIL+QRCsCoLgUBAE74nIDyLCwgpi4nMKyQiCYFEQBPv+kx7+c3qILSEDcL5BsvLCd/C8tog72Biz2RjzmTGmadjNIPIqi8jBIAi+c/63b0WEO3ERlzGmhIg8JCJ3hN0LgDzHiEjNsJtA5jDGlJOsa55FYfcCIPcxxjxnjNktIktFZL2IvB9ySwByobzyHTwvLeL2FZHTROQkERkpIu8aY/hbQMRTXES2e//bNhE5OoRekFkGiMhLQRCsDbsRALnaMhHZKCJ9jDEFjTEtJesR+WLhtoVMYYwpKCJjROTVIAiWht0PgNwnCIKbJev703kiMkFE9sV/BwAkJU98B88zi7hBEHwRBMGOIAj2BUHwqmQ9btgq7L4QaTtFpIT3v5UQkR0h9IIMYYypKyLNReTJsHsBkLsFQXBARNqKyCUiskFE7hSRcSKSqy9ekRrGmHwi8ppkzf6/JeR2AORih0f+zBaRk0WkR9j9AMhd8tJ38AJhNxCiQLIeOQRi+U5EChhjKgVBsPzw/1ZHeNwQ8TUVkYoissYYI5J1R3d+Y0z1IAjODLEvALlQEAT/kqy7b0VExBgzR0ReDa8jZAKT9QH1kmT9cGurw38hAADpVkCYiQsg9ZpKHvkOnifuxDXGlDLGXGiMKWKMKWCM6Swi54vI1LB7Q3QFQbBLsh75ecgYc5Qx5k8i0kay7loBYhkpWRendQ//GSEiU0TkwjCbQrQd/mwqIiL5JeuCo4gxJi//RSsSZIypffh4KWaM6S0iJ4jIKyG3heh7XkSqichlQRDsCbsZRB+fU8guY8xxxpiOxpjixpj8xpgLRaST5OIfHEJqcL5BEvLMd/A8sYgrIgVFZKCIbBKRzSLSU0Taej9YBfyem0WkqGTNHHxDRHoEQcCduIgpCILdQRBs+M8fyRrLsTcIgk1h94ZI6ycie0TkbhHpcjjuF2pHyBRdJeuHYjaKyAUi0sL5JXDgfxhjKojIjZL1JWeDMWbn4T+dQ24N0cbnFLIrkKzRCWtFZKuIDBWR24IgeCfUrpAJON8gW/LSd3ATBEHYPQAAAAAAAAAAYsgrd+ICAAAAAAAAQEZiERcAAAAAAAAAIoxFXAAAAAAAAACIMBZxAQAAAAAAACDCWMQFAAAAAAAAgAgrkJ0XG2OCdDWCbNscBMGxYTeRCI6b6AiCwITdQyI4ZiKFcw2SwXGDZHDcIBkcN0gGxw2SwXGDbOM7OJIQ81zDnbiZa3XYDQDIEzjXIBkcN0gGxw2SwXGDZHDcIBkcNwByQsxzDYu4AAAAAAAAABBhLOICAAAAAAAAQISxiAsAAAAAAAAAEcYiLgAAAAAAAABEGIu4AAAAAAAAABBhLOICAAAAAAAAQISxiAsAAAAAAAAAEcYiLgAAAAAAAABEGIu4AAAAAAAAABBhLOICAAAAAAAAQISxiAsAAAAAAAAAEcYiLgAAAAAAAABEWIGwGwAyUf369VV+yy232Lhbt26q9ve//93Gw4YNU7X58+enoTsAAAAgOU8//bTKe/XqZeOFCxeq2qWXXqry1atXp68xAAAi5qOPPrKxMUbVmjVrlvL9cScuAAAAAAAAAEQYi7gAAAAAAAAAEGEs4gIAAAAAAABAhOXKmbj58+dXecmSJRN+rzvbtFixYqpWpUoVG//f//2fqg0dOtTGnTp1UrW9e/fa+JFHHlG1v/3tbwn3hvDUrVtX5dOnT1d5iRIlbBwEgap17drVxq1bt1a1MmXKpKpF5BEXXHCBjceMGaNqTZo0sfGyZctyrCdEQ79+/Wzsf7bky/ffv7Nt2rSpqn388cdp7QtAZjj66KNVXrx4cRtfcsklqnbsscfa+IknnlC1ffv2paE7pFvFihVt3KVLF1U7dOiQjatVq6ZqVatWVTkzcfOWypUr27hgwYKqdv7559v4ueeeUzX3mDoSkydPtnHHjh1Vbf/+/SnZB9LLP27OPfdcGz/88MOq9qc//SlHegLiefLJJ1XuHrPu7yGlC3fiAgAAAAAAAECEsYgLAAAAAAAAABEW6XEK5cuXV3mhQoVs7N6yLCLSuHFjG5cqVUrVrrjiipT0s3btWhs/88wzqnb55ZfbeMeOHar27bff2pjHVjNHgwYNbPz222+rmj+iwx2h4P/7dx/l8ccnNGzY0Mbz58+P+T4kxn1sS0T/8544cWJOt5MWZ599to2/+uqrEDtB2K699lqV9+3b18bxHlP0R74AyDvcR+bdc4aISKNGjVRes2bNhLZ5wgknqLxXr17JNYdQbdq0ycaffPKJqvnjwJC31KhRw8b+tUf79u1t7I5uEhE58cQTbexfl6TqWsQ9NkeMGKFqt912m423b9+ekv0h9fzv1TNnzrTxhg0bVO34449XuV8H0sUdi3rTTTep2oEDB2z80Ucfpb0X7sQFAAAAAAAAgAhjERcAAAAAAAAAIoxFXAAAAAAAAACIsMjNxK1bt66NZ8yYoWr+vJR082f39OvXz8Y7d+5UtTFjxth4/fr1qrZ161YbL1u2LJUt4ggVK1bMxmeeeaaqvf766zb2573Fs3z5cpUPGTLExmPHjlW1zz77zMbu8SUiMnjw4IT3iSxNmzZVeaVKlWycqTNx/flip556qo0rVKigasaYHOkJ0eD/+y9SpEhInSAnnHPOOSrv0qWLjZs0aaJq7vxCX+/evVX+008/2dj9fQER/Tn4xRdfJN4sQlW1alUbuzMhRUQ6d+5s46JFi6qa/xny448/2tif91+tWjUbd+jQQdWee+45Gy9dujTRthGyXbt22Xj16tUhdoKocb+TtGrVKsRO4uvWrZvKX3rpJRu737mQOfwZuMzERVjc3zIqWLCgqs2ePdvG48aNS3sv3IkLAAAAAAAAABHGIi4AAAAAAAAARFjkximsWbPGxlu2bFG1VIxT8B8H/PXXX1X+5z//2cb79+9Xtddee+2I949oeeGFF2zcqVOnlGzTH8tQvHhxG3/88ceq5j7+X7t27ZTsPy/zH6P6/PPPQ+okdfxRHn/9619t7D7qLMJjq3lB8+bNbdyzZ8+Yr/OPhUsvvdTGP//8c+obQ1pcddVVNn766adVrWzZsjb2H4OfNWuWyo899lgbP/bYYzH352/HfV/Hjh3/uGHkGPea+NFHH1U197g5+uijE96mPw7qwgsvtLH/6KB7jnGPxd/LkRlKlSpl4zp16oTYCaJm+vTpNo43TmHjxo0qd8cZ+OPB/LGFrnPPPVfl/sgg5B2MisPvOf/881V+33332dhf0/nll1+S2oe/nZo1a9p45cqVquaPKks37sQFAAAAAAAAgAhjERcAAAAAAAAAIoxFXAAAAAAAAACIsMjNxHVnVvTp00fV3Jl+33zzjao988wzMbe5YMECG7do0ULVdu3apfIaNWrY+NZbb02gY2SS+vXrq/ySSy6xcbyZO/4s23fffVflQ4cOtfFPP/2kau6xunXrVlVr1qxZQvtHYvx5W7nBqFGjYtb8+YXIfRo3bqzy0aNH2zjenHh/7unq1atT2xhSpkCB/16KnXXWWar24osv2rhYsWKq9sknn9h4wIABqjZ79myVFy5c2Mbjxo1TtZYtW8bsbd68eTFrCNfll19u47/85S9JbcOf6eZfI//44482PuOMM5LaBzKHe44pX758wu87++yzVe7OS+azJ3d4/vnnbTxp0qSYrztw4IDKN2zYkNT+SpQoofKFCxfa+MQTT4z5Pr83PsMyXxAEKi9SpEhInSBKRo4cqfJKlSrZuHr16qrmXxMn6t5771V5mTJlbOz+Ro2IyLfffpvUPpKV+1Y8AAAAAAAAACAXYREXAAAAAAAAACIscuMUXP4jETNmzLDxjh07VK1OnTo2vv7661XNfdTdH5/gW7RokY1vuOGGxJtFZNWtW9fG06dPVzX3cR3/cY0PPvjAxp06dVK1Jk2aqLxfv3429h9/37Rpk439W+0PHTpkY3e0g4jImWeeaeP58+cLfl/t2rVtXK5cuRA7SY94j8z7xzNyn2uuuUbl8R4jnDVrlo3//ve/p6slpFiXLl1sHG98iv/f+1VXXWXj7du3x92H+9p44xPWrl2r8ldffTXudhGe9u3bJ/S6VatWqfyrr76ycd++fVXNHZ/gq1atWuLNISO548BeeeUVVevfv3/M9/m1X3/91cbDhw9PRWsI2cGDB20c7zyRKhdeeKHKjznmmITe53+G7du3L2U9IRr8sVNz584NqROEaffu3Sp313GOZOSGu25UoUIFVXPXbcIe68GduAAAAAAAAAAQYSziAgAAAAAAAECEsYgLAAAAAAAAABEW6Zm4vngz37Zt2xaz9te//tXGb775pqq5sy2QO1SuXFnlffr0sbE/X3Tz5s02Xr9+vaq5swB37typalOmTImbJ6No0aIqv/POO23cuXPnI95+btWqVSsb+/8MM5U72/fUU0+N+bp169blRDvIQWXLllX5ddddp3L3M8udOygiMnDgwPQ1hpQZMGCAyu+9914b+7PZn3vuORu7s9dF/ngOruu+++5L6HW9evVSuTvTHdHiXtv6v+Hw4Ycf2njFihWqtnHjxqT2lxtnziM2/zwVbyYukAodO3a0sXt+E0n8+v6BBx5IaU/IGe7MZRG9ruN/dz/99NNzpCdEj/u5VKtWLVVbsmSJjf3fIIrnqKOOUrn7WwHFihVTNXf+8vjx4xPeRzpwJy4AAAAAAAAARBiLuAAAAAAAAAAQYRk1TiEe9zGf+vXrq1qTJk1s3Lx5c1VzHzlD5ipcuLCNhw4dqmru4/Y7duxQtW7dutl43rx5qhb2o/nly5cPdf+ZokqVKjFrixYtysFOUsc9hv1HWL/77jsb+8czMlPFihVt/Pbbbyf8vmHDhql85syZqWoJKeY+4umOTxAR2b9/v42nTZumau5jXXv27Im5/SJFiqi8ZcuWKnc/T4wxquaO4Zg8eXLMfSBafvrpJxvnxKPujRo1Svs+EF358v33vh9G0SEZ/mi4u+++W+VnnHGGjQsWLJjwdhcsWGDjAwcOJNkdwuSPB/v0009tfOmll+Z0O4iIU045ReXumBV/BMctt9xi4+yMAnviiSdU3r59exu711kiIn/6058S3m66cScuAAAAAAAAAEQYi7gAAAAAAAAAEGEs4gIAAAAAAABAhOWambi7du2ysTsvQ0Rk/vz5Nn7xxRdVzZ8h6M5FffbZZ1UtCIIj7hPpUa9ePRu7M3B9bdq0UfnHH3+ctp4Qvq+++irsFqwSJUqo/KKLLrJxly5dVM2fZ+kaMGCAjf0ZUshM7rFQu3btuK/96KOPbPz000+nrSccmVKlSqn85ptvtrF/LeHOwW3btm3C+3DnB44ZM0bV/N8GcI0fP17lQ4YMSXifyHy9evWy8VFHHZXw+2rVqhWzNmfOHJV//vnn2W8MkebOweX7UN7jzu7v2rWrqvm/NxNL48aNVZ6d42j79u029mfpvv/++zaONzseQPTVrFnTxhMnTlS1smXL2tj/XZDsrOn07t3bxtdee23M1w0aNCjhbeY07sQFAAAAAAAAgAhjERcAAAAAAAAAIizXjFNwrVy5UuXubdKjR49WNf+REDf3HzP7+9//buP169cfaZtIoSeeeMLGxhhVc2+vj9r4hHz5/vv3KO6jakiN0qVLJ/W+OnXq2Ng/ntzHxk4++WRVK1SokI07d+6sau6/axH9yNcXX3yhavv27bNxgQL6NP3111/H7R3R5z8y/8gjj8R87ezZs1V+zTXX2Hjbtm2pbQwp454LRPQjYD738fbjjjtO1bp3727j1q1bq5r7yFnx4sVVzX9M1c1ff/11VXPHUSEzFStWTOXVq1e38YMPPqhq8UZO+Z9T8a5LfvrpJxu7x6mIyG+//Ra7WQCR536+iIi88847Ni5fvnxOtyOffvqpjUeOHJnj+0d0lClTJuwWcITc77b+SMGXXnrJxvGuSRo1aqRq99xzj43ddSGR/10PaN++vY397/nuet8LL7zw+/8PRAB34gIAAAAAAABAhLGICwAAAAAAAAARxiIuAAAAAAAAAERYrpyJ65s4caKNly9frmr+zIwLLrjAxg8//LCqVahQwcaDBg1StXXr1h1xn0jcpZdeqvK6deva2J8F6M5xihp3tovf94IFC3K6nYzkzpb1/xmOGDHCxvfee2/C26xdu7aN/Vk5Bw8etPHu3btVbfHixTZ++eWXVW3evHkqd+cz//zzz6q2du1aGxctWlTVli5dGrd3RFPFihVt/Pbbbyf8vu+//17l/rGCaNq/f7/KN23aZONjjz1W1X744Qcb++eweNyZpNu3b1e1E044QeWbN2+28bvvvpvwPhAdBQsWVHm9evVs7J9T3H//7mekiD5uPv/8c1W76KKLVO7P2nW5M+3atWunak8//bSN/f8WAGQe91rYvy5OVHZmbvvc730XX3yxqn3wwQdJ9YPM5P8+ADJPx44dbTxq1ChVc6+D/XPEihUrbHzWWWepmpu3adNG1U466SSVu9dI7vW5iMh1110Xt/eo4E5cAAAAAAAAAIgwFnEBAAAAAAAAIMJYxAUAAAAAAACACMsTM3FdCxcuVHmHDh1Uftlll9l49OjRqnbjjTfauFKlSqrWokWLVLWIBPhzQgsVKmTjjRs3qtqbb76ZIz3FUrhwYRv3798/5utmzJih8nvuuSddLeUqN998s41Xr16taueee25S21yzZo2NJ02apGpLliyx8dy5c5Pavu+GG25QuTsz05+JiszUt29fG2dnDtwjjzySjnaQZr/++qvK27Zta+P33ntP1UqXLm3jlStXqtrkyZNt/Morr6jaL7/8YuOxY8eqmj8T168jM7jXNv682gkTJsR839/+9jcb+9cWn332mY3dY+/3XluzZs2Y+3A/pwYPHqxq8T5D9+3bF3ObiC53nukffYadf/75Nh4+fHjaekL6+N+XmzZtauMuXbqo2rRp02y8d+/epPd5/fXX27hnz55JbweZb+bMmTb2fwcHmeeqq65SubvGduDAAVVzr5+vvvpqVdu6dauNH3/8cVVr0qSJjf15uf4cb3fubtmyZVXtxx9/tLF73hP532v0MHEnLgAAAAAAAABEGIu4AAAAAAAAABBheW6cgs9/5PG1116z8ahRo1StQIH//uNyHxUS0bdbz5o1K3UNItv8R/XWr1+fo/t3xyeIiPTr18/Gffr0UbW1a9fa2H8sYOfOnWnoLnd79NFHw24hKRdccEHM2ttvv52DnSBV6tatq/KWLVsm9D738XkRkWXLlqWsJ4Tniy++sLH7GPqRcK9D3MfIRP73cWfGsmSGggULqtwdi+BfP7g++OADlQ8bNszG/nWue/y9//77qlarVi2V79+/38ZDhgxRNXfUQps2bVRtzJgxNv7nP/+pau7ntPtopG/BggUxa8h57jnFfRT197Rr187G1atXV7XFixentjHkCHdc2aBBg9KyD3fkHOMU8jZ3JI/P/5ysUKGCjf2xeogGdySpiP73O3DgQFXzx5nG4p8jXnjhBRs3atQo4d78UQvuKI8ojU/wcScuAAAAAAAAAEQYi7gAAAAAAAAAEGEs4gIAAAAAAABAhOW5mbi1a9dW+ZVXXqnys88+28buDFyfP9Ppk08+SUF3SIV33nknx/fpzr/059ZdddVVNvbnXV5xxRXpbQwZb+LEiWG3gCR8+OGHKj/mmGNivnbu3Lk2vvbaa9PVEnKZokWL2tifgevPrBw7dmyO9ITsy58/v40HDBigar1797bxrl27VO3uu++2sf/v152De9ZZZ6na8OHDbVyvXj1VW758ucp79OhhY3dOnIhIiRIlbHzuueeqWufOnW3cunVrVZs+fbrE8uOPP9r41FNPjfk65LwRI0bY2J9vGM8NN9yg8ttuuy1lPSF3ufDCC8NuARFx8ODBmDV/hqn/WzSIHn/9Y8KECTZ2P/ezo2zZsip35/T7OnXqpPKFCxfGfK37e0VRxp24AAAAAAAAABBhLOICAAAAAAAAQITlynEKVapUUfktt9xi43bt2qna8ccfn/B2f/vtNxuvX79e1fxHGZFe/qMUbt62bVtVu/XWW1O+/9tvv13l999/v41LliypamPGjLFxt27dUt4LgOgpU6aMyuN9Rjz33HM23rlzZ9p6Qu4ybdq0sFtACriPm7vjE0REdu/ebWP/EXZ3ZEvDhg1VrXv37ja++OKLVc0dw/HQQw+p2ujRo1Ue7zHH7du323jq1Kmq5ub+Y4xXX311zG3611aIjqVLl4bdAlKsYMGCKm/ZsqWNZ8yYoWp79uxJ+f7d85SIyNNPP53yfSAzuY/f++eeqlWrqtwd0XLzzTentzEkJVX/bbtrLO3bt1c1d8TTypUrVW3cuHEp2X+UcCcuAAAAAAAAAEQYi7gAAAAAAAAAEGEs4gIAAAAAAABAhGXsTFx/lq07c8udgSsiUrFixaT2MW/ePJUPGjTIxu+8805S20RqBEEQM/ePjWeeecbGL7/8sqpt2bLFxv5Mua5du9q4Tp06qnbyySerfM2aNTb25xS68y6BRLgznitXrqxqc+fOzel2kCB3nmS+fIn/HemcOXPS0Q5yuQsvvDDsFpACDzzwQMxa/vz5bdynTx9V69+/v43POOOMhPfnvm/w4MGq5v72Q6q88cYbcXNkhmHDhtm4Z8+eqnb66afHfJ//uxTudvy5hUi/xo0b2/i+++5TtRYtWtj41FNPVbV487HjKV26tI1btWqlak888YTKixUrFnM77kzevXv3JtULMpM7/11E5KSTTlL5HXfckZPtIETuzOMePXqo2saNG23crFmzHOspLNyJCwAAAAAAAAARxiIuAAAAAAAAAERYpMcplCtXTuXVq1e38fDhw1WtatWqSe3jiy++UPljjz1m48mTJ6vaoUOHktoHcpb7+KGIvvX+iiuuULXt27fbuFKlSgnvw3/8eebMmTaO92gkkAh3PEh2HstHzqpbt67KmzdvbmP/82L//v02fvbZZ1Xt559/TkN3yO1OO+20sFtACmzYsMHGxx57rKoVLlzYxv5YJ9f777+v8k8++cTGkyZNUrVVq1bZOB3jE5D7LVq0SOXxzkV8d4oW9/tzzZo1Y77urrvuUvmOHTuS2p87ouHMM89UNX80nmvWrFkqf/75523sfudC3uMfN+71NXKXChUqqPwvf/mLjf3jYOTIkTZeu3ZtehuLAFYHAAAAAAAAACDCWMQFAAAAAAAAgAhjERcAAAAAAAAAIiz0mbilS5dW+QsvvGBjf95gsvPf3Pmljz/+uKpNmzZN5Xv27ElqH8hZn3/+ucq/+uorG5999tkx33f88cer3J+77NqyZYuNx44dq2q33nprQn0CR6pRo0Yqf+WVV8JpBP+jVKlSKvfPL65169bZuHfv3mnrCXnHp59+amN/djZzKDPH+eefb+O2bduqmjtDcuPGjar28ssv23jr1q2qxoxApJM7e1BE5LLLLgupE6RLjx490r4P/5z27rvv2tj/nrV3796094PMUKJECZW3adPGxhMnTszpdpBG06dPV7k7I/f1119XtQcffDBHeooK7sQFAAAAAAAAgAhjERcAAAAAAAAAIixHximcc845Ku/Tp4+NGzRooGonnXRSUvvYvXu3jZ955hlVe/jhh228a9eupLaPaFm7dq3K27VrZ+Mbb7xR1fr165fQNp9++mmVP//88zZesWJFdlsEkmaMCbsFABG3cOFCGy9fvlzV/PFTp59+uo03bdqU3saQLTt27LDxa6+9pmp+DkTB4sWLVb5kyRKVV6tWLSfbQTZce+21Nu7Zs6eqXXPNNUe8/ZUrV6rc/X7ujgAS+d+xHO5nGvAfHTp0UPm+fftU7p9/kHuMHj1a5QMGDLDx5MmTc7qdSOFOXAAAAAAAAACIMBZxAQAAAAAAACDCWMQFAAAAAAAAgAjLkZm4l19+edw8Fn/m0nvvvWfjgwcPqtrjjz9u419//TW7LSLDrV+/3sb9+/dXNT8HouaDDz5Qefv27UPqBNmxdOlSlc+ZM8fGjRs3zul2kIe5s/9FREaNGqXyQYMG2difg+hfawFAPKtXr1Z5rVq1QuoE2bVgwQIb33zzzar25Zdf2njgwIGqdswxx9h40qRJqjZ9+nQb+3MqN2zYkHyzgIh88sknKvdnbu/Zsycn20EOGjx4cNw8L+NOXAAAAAAAAACIMBZxAQAAAH0N9YgAACAASURBVAAAACDCTBAEib/YmMRfjHT7OgiCs8JuIhEcN9ERBIEJu4dEcMxECucaJIPjJgeVKFFC5ePGjVN58+bNbTxhwgRV6969u4137dqVhu6yheMGyeC4QTI4bpAMjhtkG9/BkYSY5xruxAUAAAAAAACACGMRFwAAAAAAAAAijEVcAAAAAAAAAIiwAmE3AAAAgORt375d5R06dFD5oEGDbNyjRw9V69+/v40XL16c+uYAAAAApAR34gIAAAAAAABAhLGICwAAAAAAAAARxjgFAACAXMQfr9CzZ8/fjQEAAABkDu7EBQAAAAAAAIAIYxEXAAAAAAAAACKMRVwAAAAAAAAAiLDszsTdLCKr09EIsq1C2A1kA8dNNHDMIBkcN0gGxw2SwXGDZHDcIBkcN0gGxw2yi2MGyYh53JggCHKyEQAAAAAAAABANjBOAQAAAAAAAAAijEVcAAAAAAAAAIgwFnEBAAAAAAAAIMLyzCKuMeZ1Y8x6Y8x2Y8x3xpi/hN0Tos8YM8sYs9cYs/Pwn2Vh94To43yDZBljOhpjlhhjdhljVhpjzgu7J0SXMeYWY8w8Y8w+Y8wrYfeDzOBc0/znz2/GmGFh94XoMsYUNsa8ZIxZbYzZYYxZYIy5OOy+EH3GmIrGmPeNMVuNMRuMMcONMdn9cXXkQVwTI7uMMdWMMTOMMduMMSuMMZeH3VM65JlFXBEZLCIVgyAoISKtRWSgMaZ+yD0hM9wSBEHxw3+qhN0MMgLnG2SbMaaFiDwqIt1F5GgROV9Evg+1KUTdTyIyUEReDrsRZA7nmqa4iBwvIntE5K2Q20K0FRCRH0WkiYiUFJF+IjLOGFMxxJ6QGZ4TkY0icoKI1JWsY+jmUDtC5HFNjOw6/JdDk0XkPREpLSI3iMjrxpjKoTaWBnlmETcIgkVBEOz7T3r4z+khtgQgl+J8gyT9TUQeCoJgbhAEh4IgWBcEwbqwm0J0BUEwIQiCSSKyJexekLGukKwFlk/DbgTRFQTBriAI+gdBsOrw59N7IvKDiPAX1Pgjp4rIuCAI9gZBsEFEpopIjZB7QvRxTYzsqioiJ4rIk0EQ/BYEwQwR+UxEuobbVurlmUVcERFjzHPGmN0islRE1ovI+yG3hMww2Biz2RjzmTGmadjNIDNwvkF2GGPyi8hZInLs4cd/1h5+5LBo2L0ByNWuEZG/B0EQhN0IMocxppyIVBaRRWH3gsh7SkQ6GmOKGWNOEpGLJWshF/hdXBMjhYyI1Ay7iVTLU4u4QRDcLFm3458nIhNEZF/8dwDSV0ROE5GTRGSkiLxrjOGOSvwhzjfIpnIiUlBErpSsY6auiNSTrEdWASDljDEVJOvR5lfD7gWZwxhTUETGiMirQRAsDbsfRN4nknXn7XYRWSsi80RkUqgdIeq4JkYylknWk0V9jDEFjTEtJesap1i4baVenlrEFRE5fGv1bBE5WUR6hN0Poi0Igi+CINgRBMG+IAhelaxb8luF3RcyA+cbZMOew/93WBAE64Mg2CwiTwjnGwDp01VEZgdB8EPYjSAzGGPyichrIrJfRG4JuR1E3OHjZapk3cxwlIiUFZFjJGvWKRAL18TItiAIDohIWxG5REQ2iMidIjJOsv7yKFfJc4u4jgLCjEpkXyBZt+UD2cH5BnEFQbBVsi4y3EeaebwZQDp1E+7CRYKMMUZEXpKsu+SuOPyFGYintIiUF5Hhh2+I2SIio4XFOMTBNTGSFQTBv4IgaBIEQZkgCC6UrCeqvwy7r1TLE4u4xpjjjDEdjTHFjTH5jTEXikgnEfko7N4QXcaYUsaYC40xRYwxBYwxnSXrlzGZ44SYON/gCIwWkZ6Hj6FjROR2yfqFVeB3Hf5sKiIi+UUk/38+r8LuC9FnjDlXskZFvRV2L8gYz4tINRG5LAiCPX/0YuDwHZQ/iEiPw59XpSRrDve/wu0MGYBrYmSbMab24WvhYsaY3iJygoi8EnJbKZcnFnEl629uekjW3+hsFZGhInJbEATvhNoVoq6giAwUkU0isllEeopI2yAIvgu1K0Qd5xska4CIfCUi34nIEhH5RkQGhdoRoq6fZD12eLeIdDkcMzMOibhGRCYEQbAj7EYQfYfnJ98oWbMpNxhjdh7+0znk1hB97UTkIsn6PrVCRA5I1oIcEA/XxEhGV8n6QfGNInKBiLQIgiDX/S6N4cdoAQAAAAAAACC68sqduAAAAAAAAACQkVjEBQAAAAAAAIAIYxEXAAAAAAAAACKMRVwAAAAAAAAAiLAC2XmxMYZfQYuOzUEQHBt2E4nguImOIAhM2D0kgmMmUjjXIBkcN0gGxw2SwXGDZHDcIBkcN8g2voMjCTHPNdyJm7lWh90AgDyBcw2SwXGDZHDcIBkcN0gGxw2SwXEDICfEPNewiAsAAAAAAAAAEcYiLgAAAAAAAABEGIu4AAAAAAAAABBhLOICAAAAAAAAQISxiAsAAAAAAAAAEVYg7AYAAEBiKleubOOpU6eqWv78+W1coUKFHOsJAAAAAJB+3IkLAAAAAAAAABHGIi4AAAAAAAAARBiLuAAAAAAAAAAQYczEBQAgooYNG6byq666ysalS5dWtffeey9HegIAAADCdtppp9l48ODBqnb55ZfbuHbt2qq2dOnS9DYGpBF34gIAAAAAAABAhLGICwAAAAAAAAARlmvGKVSvXt3Gl156qardcMMNNv7qq69U7Ztvvom5zaeeekrl+/fvP5IWAQD4H+XKlbPxhAkTVK1hw4YqD4LAxgsXLlS166+/Pg3dAQAAAOE799xzVT516lQbb9q0SdWeffZZG//888/pbQzIQdyJCwAAAAAAAAARxiIuAAAAAAAAAEQYi7gAAAAAAAAAEGEZOxP3xhtvVPnQoUNtXLx48ZjvO/3001XesWPHmK/15+fOnDkzOy0CyAHuf+9XXXWVqu3du9fG9evXV7Wjjz7axp07d1a1WbNm2XjdunVJ9bVhwwaVT548WeXz5s1LarvIfJUrV1a5+/l1zjnnxH3vPffcY2P/GNqyZUsKukOUGGNs/MYbb6haq1atbOz+LoCIyNq1a9PbGIBcp2vXrjZu2bKlqtWtW9fGVapUibuduXPn2viyyy5TtW3bth1Ji4AcddRRKnev2U888URV+9Of/mTjVatWpbMtpMkll1yi8vHjx6t8xIgRNr7vvvtUbffu3elrDAgRd+ICAAAAAAAAQISxiAsAAAAAAAAAEWaCIEj8xcYk/uI0K126tMqXLFli4+OOOy4l+/j1119V7j6q/eGHH6ZkH0fg6yAIzgq7iURE6bjJ64IgMH/8qvBl55gZMmSIjXv37p2WflLh0KFDKl+8eLGN/cek3TwCj39xrkmxhg0bqnz27NkxX+s+Ti8i0qVLFxv7x03EcNykQLFixWy8bNkyVTvppJNsfMMNN6jaqFGj0ttY+nDcIBkcNwkqW7asjf3zhDv6wP8ONGfOnJjbbNq0qcrdx92XLl2qav7ol5Bx3ITIH31w7LHHxnzt1q1bbfznP/9Z1UaPHm1j/3OyQYMGNt6xY0dSff4Ojps0O+OMM2z87bffqtqnn36qcne0lP9dK0py43dwpF3Mcw134gIAAAAAAABAhLGICwAAAAAAAAARxiIuAAAAAAAAAERYgbAbSNYvv/yi8gcffNDGjz/+uKq5M+XWrFmjauXLl4+5j1KlSqn8oosusnEEZuIiF6hQoYKNixYtqmqdOnWycY8ePWJuY8qUKSrv3r17irrLDO3atUvqfVu2bLHxv/71r6S24c/eqlKlio3980e9evVUXrNmTRsPGjRI1dx+IjATFylQuXJlG//jH/9QNX/urcs/vidPnpzaxhBpu3fvtvHy5ctVzZ2JG2+WIJCoO++8U+WFChWycbVq1VStc+fOMbfjzkGtUaNGirpDKkydOtXGFStWVDX3NwYee+wxVfO/d7mqVq2q8i+//NLG7mefiMgDDzxg44ceeuiPG0bkudezvXr1UjX3e47PPzbifSd/5JFHbOzPVXavodatW6dq7jkM0VWkSBGVu/O6//3vf6tahw4dVB7lObjIOe7vZbm/YyUicu+996rcn8ft6tevn40HDx6cou5SjztxAQAAAAAAACDCWMQFAAAAAAAAgAjL2HEKvhEjRtj4pptuUrU6derYePv27UnvY/jw4Um/F3lX8+bNbew/Gu2OTChZsqSqBUGQ0PYbNmx4BN1lvgsvvNDG/qNZ3333Xcz3uY8pr1+/PuV9HX300Sr3HweK99hY69atbeyPy0Bm6tq1q439f/fvv/++jf3PL//RQORdzz77rMqbNm1qY/9Rd+A/mjRponL30We/dvnll6s83qiXeNcolSpVsvHixYtVzX8UGunVokULlbujncaNG6dq99xzT1L7cMdniIg89dRTNnYfTRXRI78Yp5A7NGvWzMbXX399wu/bt2+fyl9//fXf3aaIyN133x1zO+656JVXXlE1d3QaomvAgAEqP+ecc2zsfp6IHNlaDnIPf/3jySeftHGDBg1Uzb9eiXf94h6L/rpClEZWcicuAAAAAAAAAEQYi7gAAAAAAAAAEGEs4gIAAAAAAABAhOWambiugQMHqvy+++6zcd26dZPebqFChZJ+L3K3UaNG2bhWrVqqdvbZZye0jR07dqh8zJgxNv7qq69U7Y033rDx3r17E+4zN1q5cuXvxmG79NJLVR5vBq4/F+zFF19MS0/IOXPmzFG5+9mzatUqVbv99tttzAxcxPLll1/GrHXo0EHlffv2VXk65n4jZ51wwgkqd68DTjvttJjv8+ftH3XUUTb2Z95+/fXXKj/zzDOz3aeISL58/71HxN0fcl6BAvqr3ooVK2w8duzYtOxz/PjxNvZn4hYpUsTGJUqUUDVmXWaG/v37q7xPnz4xX/vqq6/aeNOmTao2dOhQlbt1//v6tGnTbFy2bNmY73OPPURb4cKFbdylSxdVmzVrlo3Xrl2bUy0h4tz/9v3vyu5vQ/jnmkmTJql88uTJNu7WrZuqtW/f3sb+3F13LXD//v2Jtp0W3IkLAAAAAAAAABHGIi4AAAAAAAAARFiuHKfgP0oxe/ZsG3/44Yeq5j/6Ho87puHKK69MsjtkojJlyqh88ODBKr/uuuts/Msvv6ia+3jiI488omoLFy608Z49e1RtzZo1yTWLHOOPWHnmmWds7D+eEU+jRo1UvmDBgiNrDKFo06aNjc855xxVC4LAxm+99Zaq5fWRKEiO+yi8fy5q3bq1yl944YUc6Qmp1bx5cxv7jw6ecsopR7z96tWrq3zz5s0qdx9dPPHEE1Vt9OjRNj755JNj7mPx4sVH0iKO0MyZM1Ver149G+/evTst+/RHRLnKlStn46uvvlrVRowYkZZ+kFr+iJSiRYvaePXq1armjjT8o7E+Z5xxho3vvfdeVTv22GNtvGvXLlVzxztwPZU57rrrLhsXL15c1dzjBvgPdwyCOz5BRK/xtWrVKuFtLl++XOXudZd/bePu89tvv014H+nAnbgAAAAAAAAAEGEs4gIAAAAAAABAhLGICwAAAAAAAAARlitn4nbu3FnlderUsXHNmjWT3q47Wxd5y/3336/y66+/XuXDhg2zsT/HZ+fOnelrDDnuz3/+s427du2qatdee23M9x04cEDlvXr1svHSpUtT0xxyVKlSpVR+3nnnJfS+rVu3qnzt2rVJ7f/WW29VebwZmb17905qH4gud86yz5+Ri8zkzgzMzgxcdyZp3759VW3u3Lk2XrZsWdztbNmyxcb++SbeHNxVq1bZ2P+cRM4KY0bo999/b+NFixapWo0aNWxcqVKlHOsJqeP/9sxFF11kY3/OtvtbIDfffLOqlSxZUuVPPPGEjS+55BJVc39vZNCgQar2/PPPJ9I2IqZly5Y2/uyzz1Rt/vz5Od0OMoD/+0Eud15uqmzfvl3l/u8GhIk7cQEAAAAAAAAgwljEBQAAAAAAAIAIy9hxClWrVlX5xIkTbXzGGWeoWoECqfl/85133knJdhAdxYoVs7H/yKH7COBtt92majNnzlT5tGnTbBzGo2tInwYNGqj8ww8/tHH+/PkT3o7/6POaNWts/NtvvyXZHcLk/3urX7++jfPl039HeujQIRt/8sknCe/j9ttvj1nr2bOnyitUqBDztXfeeaeN/ceg161bl3A/ANLHfbxURKRhw4YJvc/9PBHR1y/+Y6rJijc+wec+1hilxw+RM9zxUQcPHgyxE6TDggULVO6OaPHHKTRr1szGLVq0ULUnn3xS5eXLl4+5z7/97W82dkfYIXM0btxY5e7nW61atZLebtOmTW28adMmVfPHuSCzGWN+NxbRo+qKFCmiaqeffrrK3fGH7nc3EZENGzbYuFOnTqoWpe9L3IkLAAAAAAAAABHGIi4AAAAAAAAARBiLuAAAAAAAAAAQYRk7E7datWoqP/XUU22cqhm4Pnc2oT+LEJmpX79+NvZn4o4bN87G7hxUEebe5iUdOnRQeXbm4LoKFSqk8ilTpth43rx5qvbuu+/a2J33LSKycOHCpPaP1GvSpInKzzvvPBu7M3BF9MzKeDMi69atG3ObIiKtW7eO+d5du3bZeO3atapWpUoVG48fP17VOnbsaOPVq1fH3D6A9HJnV4vouf2+OXPm2NidFymS/BzcY445RuUXXXSRjc8///yEehERef/995PaP3KHwoUL29ifTejasWNHTrSDFNu3b5/Kt2/fHvO1J554oo3ffvttVfNnWrq/HfHSSy+p2qRJk7LdJ6KlS5cuKl+yZImNf/jhh5jvc+eXiog8/vjjKnc/t/xjs3fv3jZ+9tlnE+4V0VSjRg0b+781c8cdd9jYv5by59663O9AIv/7HSmquBMXAAAAAAAAACKMRVwAAAAAAAAAiLCMHafgP2J811132fjRRx9VtXiP8mTHCSeckJLtIDruueceG/u35b/xxhs2ZnxC3jVhwgSVu6Nczj77bFUrW7ZsUvs466yzYuYPPvigqj311FM2HjJkiKpt3Lgxqf0jcUcffbSN3TE+vp9++knlr732mo1XrFihapUrV7Zxnz59VK1NmzYqd0cx+GNe3EfMSpYsqWozZsyIWUNmch9F9T+/kJlGjhypcvczZdu2bap29dVX23jDhg0p2f9NN92k8gEDBsR87aJFi2zsjx1KVT/ITBUrVrSxO8rHN3Xq1IS36f63UKdOHVVr1KiRjd966y1VW7ZsWcL7QHJSNYbJHcMydOhQVfvxxx9Tsg+E57rrrlO5+xnmj0FwR9D534NuvPFGlU+bNs3GrVq1UrXRo0fbeOXKlaqWnfMPomHLli02dr+PiejvzvFGtYiI7N6928aLFy9OZYs5hjtxAQAAAAAAACDCWMQFAAAAAAAAgAhjERcAAAAAAAAAIixjZ+L6nnnmGRsvX75c1UqVKhXzfQUK6H8Ew4cPt3GJEiVS1B2i6ssvv7SxP5fUPRb27NmjatOnT09vY4iMOXPmqPySSy6xcfny5VXNndlWrlw5VWvXrp3K3dlQ/uweV758+u/a7rjjDhvXr19f1S644AIbHzp0KOY2kbzGjRvb+Mknn4z5uhdffFHlDz30kI39Y8Od/ebP89qxY4fKx40bZ+PevXurWqVKlWw8YsSImNv56KOPVC1V8+yQs5iDm/u8/fbbcfNUu+yyy1T+wAMPxHztwYMHVe6eY5iBm7cULlxY5SeffLLKzz333IS2439Off311zY+88wzVa106dI2PuWUU1TN/Xw744wzVO3aa69NqBckLn/+/Co/77zzbBzvetY3ZcoUlfvnI2S+GjVq2Nhfc/E/U1zuf//+7Nrx48fHfN+bb76pcvea3f0dnN/bLqLPPZ4aNmyoau7nkH8c+Nzfu2EmLgAAAAAAAAAg5VjEBQAAAAAAAIAIYxEXAAAAAAAAACIs18zEdX3wwQcJv9af3ePOUvJng9WtW9fGFSpUUDVmCkbHOeeco/JvvvnGxvv371e1iy++2Ma9evVStfvvv9/G/vwdfx9Lly5NrllktDVr1sTNXf55adasWTbu2bOnqjVo0CCh/Tdp0kTl7ozUIUOGJLQNZE/t2rUTep07A9fnzmIS+d/ziatNmzYq//jjj23sz4OaPXt2zO089dRTNvZn6SL3+de//hV2C8gAkyZNUnm8Ocv+NdLIkSPT0hNSq2jRoio/7rjjbOzPnXU/U5o1axZzm0WKFFG5O6cwO/z3lSxZMuZrX375ZRv7s1Q3b95s41WrViXVCxI3duxYlbu/+ZCdWe3Mdc/9jj/++Ji1eN+dFy1aZON+/folvf/nn3/exv/+97+T3g6iZ+7cuSqvWbNmwu99+OGHU91OjuNOXAAAAAAAAACIMBZxAQAAAAAAACDCcuU4hewoVKiQyv0RCq4DBw7Y+LfffktbT/hjJ5xwgsrfe+89G5cvX17Vbr/9dhu//vrrqvbLL7/YePjw4armjlMoXry4qpUuXTqbHQPamDFjbPzmm2+q2j//+U8bn3/++Qlv0x0Hg/QoVaqUjf1xPJMnT475PnccT8WKFVXN3c6dd96pau74BBGRypUr2/gf//hHwttxxykg91u5cmXYLSCi3McI8+XT93IcOnQo5vv8cxGiwx+Z0L9/fxtfdtllqla1atWk9rF9+3Yb79ixQ9UOHjyo8gIFYn+9HDVqlI1HjBihavPnz0+qN6TeiSeeqPLu3bvb+IorrlA1dyyC/+/w22+//d1tiOjRHsh71q1bF7Pmn2OStXbt2pRsB9FXq1YtG2fn2iZTcScuAAAAAAAAAEQYi7gAAAAAAAAAEGEs4gIAAAAAAABAhOX5mbgDBw5M+LUvvfSSjZmxEi5/5lKJEiVs3LdvX1Xz5+DGcuutt8asuTNKRUQWLlyY0DaBRPjz5L7++msbZ2cm7nfffZeynvDH3Dlwv5fH4s9mct9Xu3ZtVVuzZo3KixQpYuMffvhB1c477zwbb9u2LaFeAORu/m8/1KtXz8bxzkUi+rpo+fLlaegOqTBp0iSVt2jRwsb79u1TtSlTptjY/wxx57r771u1apWN/e9AS5cuVbk7u/37779XtTvuuMPGO3fuFETTBRdcoPKHHnoo5mv79etnY//3Rdq2bWtjfybu4sWLj6RFZAD3txr835HICU2aNLFxqubsIpr27NljY//aZtasWSrfv39/TrSUVtyJCwAAAAAAAAARxiIuAAAAAAAAAERY6OMUypQpo/LRo0fb+I033lA1P0/GCSecoPIbbrgh4fdOmDDhiPeP1HjmmWdU7j7K49f83OU+HlipUiVVW716tY3vueceVdu+fXvizSLy/PPCX//6Vxv7jwmOGzcu5fvPnz+/yuvUqZPQ+/wxDHPnzk1ZT/h97uOmffr0UbU2bdrYuGHDhqpWt25dGx999NExt9+tWzeV+4+fbd682cb9+/dXtXXr1sXcLvKWwoULh90CQlSsWDEbd+nSRdXcR+19/nX2mDFjbOw/nojoaNmypcrdMQnt2rVTtQULFiS1jwIF/vuV8dFHH1W1k046SeUbN260cYcOHVSNEQrR1bRpUxvH++7UunVrlbsj544//nhVe+CBB2Juxx3RgdzJHdGT6MixI1GwYEGV33TTTTZ+7bXX0r5/5JyqVauq/Prrr7fxpk2bVO35559XeW4493AnLgAAAAAAAABEGIu4AAAAAAAAABBhLOICAAAAAAAAQISFPhPXn7lz2WWX2bhy5cqq9tNPP9nYn/23YsUKG9evX1/V3O3cddddqlaiRImYvT3++OMx949wDR48WOUHDhywcb169VStefPmMbdzzDHH2HjKlCmq1rt3bxu7xxdyB3du19SpU1WtVq1aNnaPkVQqV66cje+44w5Va9asWULbWLJkicpnz5595I0hLvdcs3v3blVz51B+9tlnqpbsLLAdO3ao3J3J/MEHHyS1TeR+rVq1UvmwYcNC6gQ5wZ+z/eKLL9r4yiuvjPm+22+/XeXDhw9XOXNwM4P/+fLrr7/aeOHChUlts0iRIip/6623bHzJJZeo2r59+1TesWNHG8+fPz+p/SPnufOyS5YsqWoff/yxjd977z1Vc+eQXnrpparmbsef8e/PrUTus3jxYhuvX79e1dx57f7M0uxwjz9/OxUrVrTxNddck/Q+EA3u+WTatGmq5s5m79u3r6qNHz8+vY2FgDtxAQAAAAAAACDCWMQFAAAAAAAAgAgLfZyC/4jfqaeeauNGjRqp2qxZs2y8atUqVXNv1z/vvPNUzX/MzOU/grR06VIbP/jgg6q2d+/emNtBuIYOHRp2C8gwTz31lI3d8Qk+95wkIrJs2TIb79mzJ+b7ihYtqnJ/lIs7QiHeOcp//Mx9vL5Xr14x34f0+Prrr23cqVMnVXP/nTZt2jThbb766qs2/ve//61q33zzjcrdRxqRt/z8888qX7RokY1r1KiR0+0gQtzHCEXij1BYuXKljf2RZshM3333ncrr1q1r45EjR6pamTJlbPztt9+q2vfff2/jPn36qFqVKlVs/MUXX6hajx49VL5gwYJE2kbEuONT/O/Hbu4+vi4i0rZtWxs//fTTqrZ161Ybjxo1StWO5BF6ZAZ3hMLDDz+sav7YSteYMWNsfNppp6lanTp1VH7vvffa2F+radmypY03b96cQMeIsiFDhtjYv+554403bBzv2MotuBMXAAAAAAAAACKMRVwAAAAAAAAAiDAWcQEAAAAAAAAgwkKfiTt37lyVf/755zZ+7bXXVO25556zccWKFVXNzxPlzuoREalevXpS2wGQWT766CMbd+jQIebr5s+fr3J3Rum2bdtivq9kyZIqr1evXnZbFBE9A1dE5PLLL7cx81HDNWXKlLg5kEr79+9Xebw5/S1atFC5//sDyHxVq1a18Z133hnzdf681IsvvjhtPSEc7rEgIjJgwAAb9+7dW9Xy5fvvHhok7QAAA+lJREFU/TsXXXRRzG2+8847KnePsalTpybVJ6LtuOOOi1nbtGmTjadPn65q/m/RuLp3727jd9999wi6Q6Z79tlnY9b8GabDhw+P+Vr/e5E7233gwIGq5l83IbM0b95c5V26dLGx/7s048ePz5GeooI7cQEAAAAAAAAgwljEBQAAAAAAAIAIC32cgs99XKdw4cKqVrx48Zjvcx9V7tSpU8zX+Y8/+48cAsgb3MfBxo4dq2odO3aM+b5kxyLEc/DgQZU/9dRTNn777bdV7Ysvvkj5/gFkngULFti4fv36qhbvegm5w/3332/jq666Kubr/FEaq1evTltPiAb32HBjIJ4lS5bErF155ZU2Nsao2i+//GJj/5H5f/7znynqDrmNe6zEG7WAvMUdkfrmm2/GfF23bt1UPnny5HS1FEnciQsAAAAAAAAAEcYiLgAAAAAAAABEGIu4AAAAAAAAABBhkZuJ69q3b5/KH3vssYTed/XVV6ejHQC5yKpVq2zcvXt3VXvnnXds3KxZM1X77rvvbNy6deuY21+6dGnc/c+YMSPma91ZlwDwewYNGmTjmjVrqtq4ceNyuh2kWY0aNVReokSJmK8dOXKkjd3PGgCI5dVXX7VxoUKFVM2drTxv3jxVc6+Zn3zyyTR1ByA3Klq0qMrd38cqWbKkqrm/EzNx4sT0NhZx3IkLAAAAAAAAABHGIi4AAAAAAAAARJgJgiDxFxuT+IuRbl8HQXBW2E0kguMmOoIgMGH3kAiOmUjhXINkcNwgGRw3MTz66KMqdx85XL16taq1atXKxsuWLUtvY9HAcYNkcNwgGRw3yDa+g/++Hj16qHz48OE2njNnjqo1b97cxv7Y1Vwq5rmGO3EBAAAAAAAAIMJYxAUAAAAAAACACGMRFwAAAAAAAAAirEDYDQAAAACI7cMPP1S5OxP3jjvuULU8MgcXAABkmAYNGtj43nvvVbWBAwfa+MUXX1S1PDIHNyHciQsAAAAAAAAAEcYiLgAAAAAAAABEGOMUAAAAgAj76KOPVF6gAJfwAAAgs3z55Zc2PuWUU0LsJHNxJy4AAAAAAAAARBiLuAAAAAAAAAAQYSziAgAAAAAAAECEZXeg1mYRWZ2ORpBtFcJuIBs4bqKBYwbJ4LhBMjhukAyOGySD4wbJ4LhBMjhukF0cM0hGzOPGBEGQk40AAAAAAAAAALKBcQoAAAAAAAAAEGEs4gIAAAAAAABAhLGICwAAAAAAAAARxiIuAAAAAAAAAEQYi7gAAAAAAAAAEGEs4gIAAAAAAABAhLGICwAAAAAAAAARxiIuAAAAAAAAAEQYi7gAAAAAAAAAEGH/DyO8YKtWba0ZAAAAAElFTkSuQmCC",
            "text/plain": [
              "<Figure size 1800x288 with 20 Axes>"
            ]
          },
          "metadata": {},
          "output_type": "display_data"
        }
      ],
      "source": [
        "import matplotlib.pyplot as plt\n",
        "%matplotlib inline\n",
        "    \n",
        "# obtain one batch of training images\n",
        "dataiter = iter(train_loader)\n",
        "images, labels = dataiter.next()\n",
        "images = images.numpy()\n",
        "\n",
        "# plot the images in the batch, along with the corresponding labels\n",
        "fig = plt.figure(figsize=(25, 4))\n",
        "for idx in np.arange(batch_size):\n",
        "    ax = fig.add_subplot(2, batch_size/2, idx+1, xticks=[], yticks=[])\n",
        "    ax.imshow(np.squeeze(images[idx]), cmap='gray')\n",
        "    # print out the correct label for each image\n",
        "    # .item() gets the value contained in a Tensor\n",
        "    ax.set_title(str(labels[idx].item()))"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 49,
      "metadata": {
        "collapsed": true,
        "id": "WbTwHZENKmhw"
      },
      "outputs": [],
      "source": [
        "import torch.nn as nn\n",
        "import torch.nn.functional as F\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 50,
      "metadata": {
        "collapsed": true,
        "id": "nEiWchRIKmhy"
      },
      "outputs": [],
      "source": [
        "class ConvLayer(nn.Module):\n",
        "    \n",
        "    def __init__(self, in_channels=1, out_channels=256):\n",
        "        '''Constructs the ConvLayer with a specified input and output size.\n",
        "           param in_channels: input depth of an image, default value = 1\n",
        "           param out_channels: output depth of the convolutional layer, default value = 256\n",
        "           '''\n",
        "        super(ConvLayer, self).__init__()\n",
        "\n",
        "        # defining a convolutional layer of the specified size\n",
        "        self.conv = nn.Conv2d(in_channels, out_channels, \n",
        "                              kernel_size=9, stride=1, padding=0)\n",
        "\n",
        "    def forward(self, x):\n",
        "        '''Defines the feedforward behavior.\n",
        "           param x: the input to the layer; an input image\n",
        "           return: a relu-activated, convolutional layer\n",
        "           '''\n",
        "        # applying a ReLu activation to the outputs of the conv layer\n",
        "        features = F.relu(self.conv(x)) # will have dimensions (batch_size, 20, 20, 256)\n",
        "        return features\n",
        "    "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 51,
      "metadata": {
        "collapsed": true,
        "id": "z_qooHsGKmh1"
      },
      "outputs": [],
      "source": [
        "from torchvision.transforms.functional import rotate, affine\n",
        "\n",
        "class PrimaryCaps(nn.Module):\n",
        "    \n",
        "    def __init__(self, num_capsules=8, in_channels=256, out_channels=32):\n",
        "        '''Constructs a list of convolutional layers to be used in \n",
        "           creating capsule output vectors.\n",
        "           param num_capsules: number of capsules to create\n",
        "           param in_channels: input depth of features, default value = 256\n",
        "           param out_channels: output depth of the convolutional layers, default value = 32\n",
        "           '''\n",
        "        super(PrimaryCaps, self).__init__()\n",
        "\n",
        "        # creating a list of convolutional layers for each capsule I want to create\n",
        "        # all capsules have a conv layer with the same parameters\n",
        "        self.capsules = nn.ModuleList([\n",
        "            nn.Conv2d(in_channels=in_channels, out_channels=out_channels, \n",
        "                      kernel_size=9, stride=2, padding=0)\n",
        "            for _ in range(num_capsules)])\n",
        "    \n",
        "    def forward(self, x):\n",
        "        global primary_caps_total_time\n",
        "        '''Defines the feedforward behavior.\n",
        "           param x: the input; features from a convolutional layer\n",
        "           return: a set of normalized, capsule output vectors\n",
        "           '''\n",
        "        s = time.time()\n",
        "        # get batch size of inputs\n",
        "        batch_size = x.size(0)\n",
        "        # reshape convolutional layer outputs to be (batch_size, vector_dim=1152, 1)\n",
        "        u = []\n",
        "\n",
        "        i = 0\n",
        "\n",
        "        # transformations: \n",
        "\n",
        "        # 0 -> rotate by 0 degrees\n",
        "        # 1 -> rotate by 60 degrees\n",
        "        # 2 -> rotate by 120 degrees\n",
        "        # 3 -> rotate by 180 degrees\n",
        "        # 4 -> rotate by 240 degrees\n",
        "        # 5 -> rotate by 300 degrees\n",
        "        # 6 -> shear, -22.5 degrees\n",
        "        # 7 -> shear, 22.5 degrees\n",
        "        # 8 -> scale, 0.5\n",
        "        # 9 -> scale, 1.5\n",
        "\n",
        "        for capsule in self.capsules:\n",
        "            transformed_x = None\n",
        "\n",
        "            if i == 1:\n",
        "                transformed_x = rotate(x, -60) \n",
        "            elif i == 2:\n",
        "                transformed_x = rotate(x, -30) \n",
        "            elif i == 3:\n",
        "                transformed_x = x \n",
        "            elif i == 4:\n",
        "                transformed_x = rotate(x, 30) \n",
        "            elif i == 5:\n",
        "                transformed_x = rotate(x, -60) \n",
        "            elif i == 6:\n",
        "                transformed_x = affine(x, 0, [0, 0], 1, [-22.5, 0])\n",
        "            elif i == 7:\n",
        "                transformed_x = affine(x, 0, [0, 0], 1, [22.5, 0])\n",
        "            elif i == 8:\n",
        "                transformed_x = affine(x, 0, [0, 0], 0.5, 0)\n",
        "            elif i == 9:\n",
        "                transformed_x = affine(x, 0, [0, 0], 1.5, 0)\n",
        "            else:\n",
        "                transformed_x = x\n",
        "                \n",
        "            output = capsule(transformed_x)\n",
        "            reshaped_output = output.view(batch_size, 32 * 6 * 6, 1)\n",
        "            u.append(reshaped_output)\n",
        "            \n",
        "\n",
        "            i += 1\n",
        "\n",
        "        # stack up output vectors, u, one for each capsule\n",
        "        u = torch.cat(u, dim=-1)\n",
        "        # squashing the stack of vectors\n",
        "        u_squash = self.squash(u)\n",
        "        primary_caps_total_time += time.time() - s\n",
        "        # print(\"primary caps forward time\", time.time() - s)\n",
        "        return u_squash\n",
        "    \n",
        "    def squash(self, input_tensor):\n",
        "        '''Squashes an input Tensor so it has a magnitude between 0-1.\n",
        "           param input_tensor: a stack of capsule inputs, s_j\n",
        "           return: a stack of normalized, capsule output vectors, v_j\n",
        "           '''\n",
        "        squared_norm = (input_tensor ** 2).sum(dim=-1, keepdim=True)\n",
        "        scale = squared_norm / (1 + squared_norm) # normalization coeff\n",
        "        output_tensor = scale * input_tensor / torch.sqrt(squared_norm)    \n",
        "        return output_tensor\n",
        "    "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 52,
      "metadata": {
        "collapsed": true,
        "id": "C2EeuUutKmh4"
      },
      "outputs": [],
      "source": [
        "import torch\n",
        "import torch.nn.functional as F\n",
        "\n",
        "def softmax(input_tensor, dim=1):\n",
        "    # transpose input\n",
        "    transposed_input = input_tensor.transpose(dim, len(input_tensor.size()) - 1)\n",
        "    # calculate softmax\n",
        "    softmaxed_output = F.softmax(transposed_input.contiguous().view(-1, transposed_input.size(-1)), dim=-1)\n",
        "    # un-transpose result\n",
        "    return softmaxed_output.view(*transposed_input.size()).transpose(dim, len(input_tensor.size()) - 1)\n",
        "\n",
        "# dynamic routing\n",
        "def dynamic_routing(b_ij, u_hat, squash, routing_iterations=3):\n",
        "    '''Performs dynamic routing between two capsule layers.\n",
        "       param b_ij: initial log probabilities that capsule i should be coupled to capsule j\n",
        "       param u_hat: input, weighted capsule vectors, W u\n",
        "       param squash: given, normalizing squash function\n",
        "       param routing_iterations: number of times to update coupling coefficients\n",
        "       return: v_j, output capsule vectors\n",
        "       '''    \n",
        "    global dynamic_routing_total_time\n",
        "    s = time.time()  \n",
        "    # update b_ij, c_ij for number of routing iterations\n",
        "    for iteration in range(routing_iterations):\n",
        "        # softmax calculation of coupling coefficients, c_ij\n",
        "        c_ij = softmax(b_ij, dim=2)\n",
        "\n",
        "        # calculating total capsule inputs, s_j = sum(c_ij*u_hat)\n",
        "        s_j = (c_ij * u_hat).sum(dim=2, keepdim=True)\n",
        "\n",
        "        # squashing to get a normalized vector output, v_j\n",
        "        v_j = squash(s_j)\n",
        "\n",
        "        # if not on the last iteration, calculate agreement and new b_ij\n",
        "        if iteration < routing_iterations - 1:\n",
        "            # agreement\n",
        "            a_ij = (u_hat * v_j).sum(dim=-1, keepdim=True)\n",
        "            \n",
        "            # new b_ij\n",
        "            b_ij = b_ij + a_ij\n",
        "\n",
        "    dynamic_routing_total_time += time.time() - s\n",
        "    \n",
        "    return v_j # return latest v_j\n",
        "    "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 53,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "8QVeRfnXKmh7",
        "outputId": "631b84f3-4fed-4d11-9ca8-07f2cf30431a"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Training on GPU!\n"
          ]
        }
      ],
      "source": [
        "# it will also be relevant, in this model, to see if I can train on gpu\n",
        "TRAIN_ON_GPU = torch.cuda.is_available()\n",
        "\n",
        "if(TRAIN_ON_GPU):\n",
        "    print('Training on GPU!')\n",
        "else:\n",
        "    print('Only CPU available')\n",
        "    "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 54,
      "metadata": {
        "collapsed": true,
        "id": "ALrLvUsqKmh-"
      },
      "outputs": [],
      "source": [
        "class DigitCaps(nn.Module):\n",
        "    \n",
        "    def __init__(self, num_capsules=10, previous_layer_nodes=32*6*6, \n",
        "                 in_channels=8, out_channels=16):\n",
        "        '''Constructs an initial weight matrix, W, and sets class variables.\n",
        "           param num_capsules: number of capsules to create\n",
        "           param previous_layer_nodes: dimension of input capsule vector, default value = 1152\n",
        "           param in_channels: number of capsules in previous layer, default value = 8\n",
        "           param out_channels: dimensions of output capsule vector, default value = 16\n",
        "           '''\n",
        "        super(DigitCaps, self).__init__()\n",
        "\n",
        "        # setting class variables\n",
        "        self.num_capsules = num_capsules\n",
        "        self.previous_layer_nodes = previous_layer_nodes # vector input (dim=1152)\n",
        "        self.in_channels = in_channels # previous layer's number of capsules\n",
        "\n",
        "        # starting out with a randomly initialized weight matrix, W\n",
        "        # these will be the weights connecting the PrimaryCaps and DigitCaps layers\n",
        "        self.W = nn.Parameter(torch.randn(num_capsules, previous_layer_nodes, \n",
        "                                          in_channels, out_channels))\n",
        "\n",
        "    def forward(self, u):\n",
        "        '''Defines the feedforward behavior.\n",
        "           param u: the input; vectors from the previous PrimaryCaps layer\n",
        "           return: a set of normalized, capsule output vectors\n",
        "           '''\n",
        "        global digit_caps_total_time\n",
        "\n",
        "        s = time.time()\n",
        "        \n",
        "        # adding batch_size dims and stacking all u vectors\n",
        "        u = u[None, :, :, None, :]\n",
        "        # 4D weight matrix\n",
        "        W = self.W[:, None, :, :, :]\n",
        "        \n",
        "        # calculating u_hat = W*u\n",
        "        u_hat = torch.matmul(u, W)\n",
        "\n",
        "        # getting the correct size of b_ij\n",
        "        # setting them all to 0, initially\n",
        "        b_ij = torch.zeros(*u_hat.size())\n",
        "        \n",
        "        # moving b_ij to GPU, if available\n",
        "        if TRAIN_ON_GPU:\n",
        "            b_ij = b_ij.cuda()\n",
        "\n",
        "        # update coupling coefficients and calculate v_j\n",
        "        v_j = dynamic_routing(b_ij, u_hat, self.squash, routing_iterations=3)\n",
        "        digit_caps_total_time += time.time() - s\n",
        "        # print('DigitCaps (dynamic routing included) forward time: ', time.time() - s)\n",
        "\n",
        "        return v_j # return final vector outputs\n",
        "    \n",
        "    \n",
        "    def squash(self, input_tensor):\n",
        "        '''Squashes an input Tensor so it has a magnitude between 0-1.\n",
        "           param input_tensor: a stack of capsule inputs, s_j\n",
        "           return: a stack of normalized, capsule output vectors, v_j\n",
        "           '''\n",
        "        # same squash function as before\n",
        "        squared_norm = (input_tensor ** 2).sum(dim=-1, keepdim=True)\n",
        "        scale = squared_norm / (1 + squared_norm) # normalization coeff\n",
        "        output_tensor = scale * input_tensor / torch.sqrt(squared_norm)    \n",
        "        return output_tensor\n",
        "    "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 55,
      "metadata": {
        "collapsed": true,
        "id": "kVMKhMZoKmiB"
      },
      "outputs": [],
      "source": [
        "class Decoder(nn.Module):\n",
        "    \n",
        "    def __init__(self, input_vector_length=16, input_capsules=10, hidden_dim=512):\n",
        "        '''Constructs an series of linear layers + activations.\n",
        "           param input_vector_length: dimension of input capsule vector, default value = 16\n",
        "           param input_capsules: number of capsules in previous layer, default value = 10\n",
        "           param hidden_dim: dimensions of hidden layers, default value = 512\n",
        "           '''\n",
        "        super(Decoder, self).__init__()\n",
        "        \n",
        "        # calculate input_dim\n",
        "        input_dim = input_vector_length * input_capsules\n",
        "        \n",
        "        # define linear layers + activations\n",
        "        self.linear_layers = nn.Sequential(\n",
        "            nn.Linear(input_dim, hidden_dim), # first hidden layer\n",
        "            nn.ReLU(inplace=True),\n",
        "            nn.Linear(hidden_dim, hidden_dim*2), # second, twice as deep\n",
        "            nn.ReLU(inplace=True),\n",
        "            nn.Linear(hidden_dim*2, 28*28), # can be reshaped into 28*28 image\n",
        "            nn.Sigmoid() # sigmoid activation to get output pixel values in a range from 0-1\n",
        "            )\n",
        "        \n",
        "    def forward(self, x):\n",
        "        '''Defines the feedforward behavior.\n",
        "           param x: the input; vectors from the previous DigitCaps layer\n",
        "           return: two things, reconstructed images and the class scores, y\n",
        "           '''\n",
        "        global decoder_total_time\n",
        "        s = time.time()\n",
        "        classes = (x ** 2).sum(dim=-1) ** 0.5\n",
        "        classes = F.softmax(classes, dim=-1)\n",
        "        \n",
        "        # find the capsule with the maximum vector length\n",
        "        # here, vector length indicates the probability of a class' existence\n",
        "        _, max_length_indices = classes.max(dim=1)\n",
        "        \n",
        "        # create a sparse class matrix\n",
        "        sparse_matrix = torch.eye(10) # 10 is the number of classes\n",
        "        if TRAIN_ON_GPU:\n",
        "            sparse_matrix = sparse_matrix.cuda()\n",
        "        # get the class scores from the \"correct\" capsule\n",
        "        y = sparse_matrix.index_select(dim=0, index=max_length_indices.data)\n",
        "        \n",
        "        # create reconstructed pixels\n",
        "        x = x * y[:, :, None]\n",
        "        # flatten image into a vector shape (batch_size, vector_dim)\n",
        "        flattened_x = x.contiguous().view(x.size(0), -1)\n",
        "        # create reconstructed image vectors\n",
        "        reconstructions = self.linear_layers(flattened_x)\n",
        "        decoder_total_time += time.time() - s\n",
        "        # print('Decoder forward time: {}'.format(time.time() - s), 'seconds')\n",
        "        \n",
        "        # return reconstructions and the class scores, y\n",
        "        return reconstructions, y"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 56,
      "metadata": {
        "collapsed": true,
        "id": "ATYXUK9PKmiE"
      },
      "outputs": [],
      "source": [
        "class CapsuleNetwork(nn.Module):\n",
        "    \n",
        "    def __init__(self):\n",
        "        '''Constructs a complete Capsule Network.'''\n",
        "        super(CapsuleNetwork, self).__init__()\n",
        "        self.conv_layer = ConvLayer()\n",
        "        self.primary_capsules = PrimaryCaps()\n",
        "        self.digit_capsules = DigitCaps()\n",
        "        self.decoder = Decoder()\n",
        "                \n",
        "    def forward(self, images):\n",
        "        global total_capsule_net_time\n",
        "        '''Defines the feedforward behavior.\n",
        "           param images: the original MNIST image input data\n",
        "           return: output of DigitCaps layer, reconstructed images, class scores\n",
        "           '''\n",
        "        s = time.time()\n",
        "        primary_caps_output = self.primary_capsules(self.conv_layer(images))\n",
        "        caps_output = self.digit_capsules(primary_caps_output).squeeze().transpose(0,1)\n",
        "        reconstructions, y = self.decoder(caps_output)\n",
        "        total_capsule_net_time += time.time() - s\n",
        "        # print(\"Time to run entire Capsule Network:\", time.time() - s, \"seconds\")\n",
        "        return caps_output, reconstructions, y\n",
        "    "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 57,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "Ri97KozmKmiG",
        "outputId": "0d1cb957-137a-403d-9558-06241f4e28c1"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "CapsuleNetwork(\n",
            "  (conv_layer): ConvLayer(\n",
            "    (conv): Conv2d(1, 256, kernel_size=(9, 9), stride=(1, 1))\n",
            "  )\n",
            "  (primary_capsules): PrimaryCaps(\n",
            "    (capsules): ModuleList(\n",
            "      (0): Conv2d(256, 32, kernel_size=(9, 9), stride=(2, 2))\n",
            "      (1): Conv2d(256, 32, kernel_size=(9, 9), stride=(2, 2))\n",
            "      (2): Conv2d(256, 32, kernel_size=(9, 9), stride=(2, 2))\n",
            "      (3): Conv2d(256, 32, kernel_size=(9, 9), stride=(2, 2))\n",
            "      (4): Conv2d(256, 32, kernel_size=(9, 9), stride=(2, 2))\n",
            "      (5): Conv2d(256, 32, kernel_size=(9, 9), stride=(2, 2))\n",
            "      (6): Conv2d(256, 32, kernel_size=(9, 9), stride=(2, 2))\n",
            "      (7): Conv2d(256, 32, kernel_size=(9, 9), stride=(2, 2))\n",
            "    )\n",
            "  )\n",
            "  (digit_capsules): DigitCaps()\n",
            "  (decoder): Decoder(\n",
            "    (linear_layers): Sequential(\n",
            "      (0): Linear(in_features=160, out_features=512, bias=True)\n",
            "      (1): ReLU(inplace=True)\n",
            "      (2): Linear(in_features=512, out_features=1024, bias=True)\n",
            "      (3): ReLU(inplace=True)\n",
            "      (4): Linear(in_features=1024, out_features=784, bias=True)\n",
            "      (5): Sigmoid()\n",
            "    )\n",
            "  )\n",
            ")\n"
          ]
        }
      ],
      "source": [
        "# instantiate and print net\n",
        "capsule_net = CapsuleNetwork()\n",
        "\n",
        "print(capsule_net)\n",
        "\n",
        "# move model to GPU, if available \n",
        "if TRAIN_ON_GPU:\n",
        "    capsule_net = capsule_net.cuda()"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 58,
      "metadata": {
        "collapsed": true,
        "id": "NQyj0NgQKmiH"
      },
      "outputs": [],
      "source": [
        "class CapsuleLoss(nn.Module):\n",
        "    \n",
        "    def __init__(self):\n",
        "        '''Constructs a CapsuleLoss module.'''\n",
        "        super(CapsuleLoss, self).__init__()\n",
        "        self.reconstruction_loss = nn.MSELoss(reduction='sum') # cumulative loss, equiv to size_average=False\n",
        "\n",
        "    def forward(self, x, labels, images, reconstructions):\n",
        "        global capsule_loss_total_time\n",
        "        s = time.time()\n",
        "        '''Defines how the loss compares inputs.\n",
        "           param x: digit capsule outputs\n",
        "           param labels: \n",
        "           param images: the original MNIST image input data\n",
        "           param reconstructions: reconstructed MNIST image data\n",
        "           return: weighted margin and reconstruction loss, averaged over a batch\n",
        "           '''\n",
        "        batch_size = x.size(0)\n",
        "\n",
        "        ##  calculate the margin loss   ##\n",
        "        \n",
        "        # get magnitude of digit capsule vectors, v_c\n",
        "        v_c = torch.sqrt((x**2).sum(dim=2, keepdim=True))\n",
        "\n",
        "        # calculate \"correct\" and incorrect loss\n",
        "        left = F.relu(0.9 - v_c).view(batch_size, -1)\n",
        "        right = F.relu(v_c - 0.1).view(batch_size, -1)\n",
        "        \n",
        "        # sum the losses, with a lambda = 0.5\n",
        "        margin_loss = labels * left + 0.5 * (1. - labels) * right\n",
        "        margin_loss = margin_loss.sum()\n",
        "\n",
        "        ##  calculate the reconstruction loss   ##\n",
        "        images = images.view(reconstructions.size()[0], -1)\n",
        "        reconstruction_loss = self.reconstruction_loss(reconstructions, images)\n",
        "\n",
        "        # return a weighted, summed loss, averaged over a batch size\n",
        "        ret = (margin_loss + 0.0005 * reconstruction_loss) / images.size(0)\n",
        "        capsule_loss_total_time += time.time() - s\n",
        "        # print(\"capsule loss iteration elapsed time:\", time.time() - s, \"seconds\")\n",
        "        return ret\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 59,
      "metadata": {
        "collapsed": true,
        "id": "JSvU-GfPKmiJ"
      },
      "outputs": [],
      "source": [
        "import torch.optim as optim\n",
        "\n",
        "# custom loss\n",
        "criterion = CapsuleLoss()\n",
        "\n",
        "# Adam optimizer with default params\n",
        "optimizer = optim.Adam(capsule_net.parameters())"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 60,
      "metadata": {
        "collapsed": true,
        "id": "L6nvFF_JKmiL"
      },
      "outputs": [],
      "source": [
        "def train(capsule_net, criterion, optimizer, \n",
        "          n_epochs, print_every=100):\n",
        "    '''Trains a capsule network and prints out training batch loss statistics.\n",
        "       Saves model parameters if *validation* loss has decreased.\n",
        "       param capsule_net: trained capsule network\n",
        "       param criterion: capsule loss function\n",
        "       param optimizer: optimizer for updating network weights\n",
        "       param n_epochs: number of epochs to train for\n",
        "       param print_every: batches to print and save training loss, default = 100\n",
        "       return: list of recorded training losses\n",
        "       '''\n",
        "\n",
        "    # track training loss over time\n",
        "    losses = []\n",
        "\n",
        "    # one epoch = one pass over all training data \n",
        "    for epoch in range(1, n_epochs+1):\n",
        "\n",
        "        # initialize training loss\n",
        "        train_loss = 0.0\n",
        "        \n",
        "        capsule_net.train() # set to train mode\n",
        "    \n",
        "        # get batches of training image data and targets\n",
        "        for batch_i, (images, target) in enumerate(train_loader):\n",
        "            # reshape and get target class\n",
        "            target = torch.eye(10).index_select(dim=0, index=target)\n",
        "\n",
        "            if TRAIN_ON_GPU:\n",
        "                images, target = images.cuda(), target.cuda()\n",
        "\n",
        "            # zero out gradients\n",
        "            optimizer.zero_grad()\n",
        "            # get model outputs\n",
        "            caps_output, reconstructions, y = capsule_net(images)\n",
        "            # calculate loss\n",
        "            loss = criterion(caps_output, target, images, reconstructions)\n",
        "            # perform backpropagation and optimization\n",
        "            loss.backward()\n",
        "            optimizer.step()\n",
        "\n",
        "            train_loss += loss.item() # accumulated training loss\n",
        "            \n",
        "            # print and record training stats\n",
        "            if batch_i != 0 and batch_i % print_every == 0:\n",
        "                avg_train_loss = train_loss/print_every\n",
        "                losses.append(avg_train_loss)\n",
        "                print('Epoch: {} \\tTraining Loss: {:.8f}'.format(epoch, avg_train_loss))\n",
        "                train_loss = 0 # reset accumulated training loss\n",
        "        \n",
        "    return losses\n",
        "    "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 61,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "Gi4BQu5WKmiM",
        "outputId": "24ed2876-bfac-44db-f3c1-58bac55bfb81"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Epoch: 1 \tTraining Loss: 0.45348430\n",
            "Epoch: 1 \tTraining Loss: 0.17312004\n",
            "Epoch: 1 \tTraining Loss: 0.12498576\n",
            "Epoch: 1 \tTraining Loss: 0.10535942\n",
            "Epoch: 1 \tTraining Loss: 0.10441681\n",
            "Epoch: 1 \tTraining Loss: 0.08628859\n",
            "Epoch: 1 \tTraining Loss: 0.08050859\n",
            "Epoch: 1 \tTraining Loss: 0.07868407\n",
            "Epoch: 1 \tTraining Loss: 0.07097200\n",
            "Epoch: 1 \tTraining Loss: 0.06067677\n",
            "Epoch: 1 \tTraining Loss: 0.06130396\n",
            "Epoch: 1 \tTraining Loss: 0.06282486\n",
            "Epoch: 1 \tTraining Loss: 0.05862658\n",
            "Epoch: 1 \tTraining Loss: 0.06048516\n",
            "Epoch: 1 \tTraining Loss: 0.06212459\n",
            "Epoch: 1 \tTraining Loss: 0.06078425\n",
            "Epoch: 1 \tTraining Loss: 0.05157406\n",
            "Epoch: 1 \tTraining Loss: 0.05177380\n",
            "Epoch: 1 \tTraining Loss: 0.05213749\n",
            "Epoch: 1 \tTraining Loss: 0.05032017\n",
            "Epoch: 1 \tTraining Loss: 0.05031982\n",
            "Epoch: 1 \tTraining Loss: 0.04958880\n",
            "Epoch: 1 \tTraining Loss: 0.04679158\n",
            "Epoch: 1 \tTraining Loss: 0.05749342\n",
            "Epoch: 1 \tTraining Loss: 0.05061885\n",
            "Epoch: 1 \tTraining Loss: 0.04643542\n",
            "Epoch: 1 \tTraining Loss: 0.04996473\n",
            "Epoch: 1 \tTraining Loss: 0.04024940\n",
            "Epoch: 1 \tTraining Loss: 0.04218979\n",
            "Epoch: 2 \tTraining Loss: 0.04852553\n",
            "Epoch: 2 \tTraining Loss: 0.03782086\n",
            "Epoch: 2 \tTraining Loss: 0.04273885\n",
            "Epoch: 2 \tTraining Loss: 0.04019849\n",
            "Epoch: 2 \tTraining Loss: 0.04233346\n",
            "Epoch: 2 \tTraining Loss: 0.04213146\n",
            "Epoch: 2 \tTraining Loss: 0.03594931\n",
            "Epoch: 2 \tTraining Loss: 0.03521657\n",
            "Epoch: 2 \tTraining Loss: 0.03840889\n",
            "Epoch: 2 \tTraining Loss: 0.03412043\n",
            "Epoch: 2 \tTraining Loss: 0.03558842\n",
            "Epoch: 2 \tTraining Loss: 0.03522476\n",
            "Epoch: 2 \tTraining Loss: 0.03343793\n",
            "Epoch: 2 \tTraining Loss: 0.03587388\n",
            "Epoch: 2 \tTraining Loss: 0.03585835\n",
            "Epoch: 2 \tTraining Loss: 0.03667655\n",
            "Epoch: 2 \tTraining Loss: 0.03308575\n",
            "Epoch: 2 \tTraining Loss: 0.03421831\n",
            "Epoch: 2 \tTraining Loss: 0.03449699\n",
            "Epoch: 2 \tTraining Loss: 0.03395573\n",
            "Epoch: 2 \tTraining Loss: 0.03412230\n",
            "Epoch: 2 \tTraining Loss: 0.03247327\n",
            "Epoch: 2 \tTraining Loss: 0.02989372\n",
            "Epoch: 2 \tTraining Loss: 0.03778269\n",
            "Epoch: 2 \tTraining Loss: 0.03434713\n",
            "Epoch: 2 \tTraining Loss: 0.03008195\n",
            "Epoch: 2 \tTraining Loss: 0.03459434\n",
            "Epoch: 2 \tTraining Loss: 0.02836240\n",
            "Epoch: 2 \tTraining Loss: 0.02818420\n",
            "Epoch: 3 \tTraining Loss: 0.03304808\n",
            "Epoch: 3 \tTraining Loss: 0.02748846\n",
            "Epoch: 3 \tTraining Loss: 0.03248798\n",
            "Epoch: 3 \tTraining Loss: 0.02823999\n",
            "Epoch: 3 \tTraining Loss: 0.02989688\n",
            "Epoch: 3 \tTraining Loss: 0.02710802\n",
            "Epoch: 3 \tTraining Loss: 0.02633414\n",
            "Epoch: 3 \tTraining Loss: 0.02669966\n",
            "Epoch: 3 \tTraining Loss: 0.02863755\n",
            "Epoch: 3 \tTraining Loss: 0.02365689\n",
            "Epoch: 3 \tTraining Loss: 0.02677263\n",
            "Epoch: 3 \tTraining Loss: 0.02847553\n",
            "Epoch: 3 \tTraining Loss: 0.02426153\n",
            "Epoch: 3 \tTraining Loss: 0.03009300\n",
            "Epoch: 3 \tTraining Loss: 0.02624928\n",
            "Epoch: 3 \tTraining Loss: 0.02542656\n",
            "Epoch: 3 \tTraining Loss: 0.02574847\n",
            "Epoch: 3 \tTraining Loss: 0.02598032\n",
            "Epoch: 3 \tTraining Loss: 0.02857920\n",
            "Epoch: 3 \tTraining Loss: 0.02619290\n",
            "Epoch: 3 \tTraining Loss: 0.02737601\n",
            "Epoch: 3 \tTraining Loss: 0.02580921\n",
            "Epoch: 3 \tTraining Loss: 0.02196284\n",
            "Epoch: 3 \tTraining Loss: 0.02994080\n",
            "Epoch: 3 \tTraining Loss: 0.02594417\n",
            "Epoch: 3 \tTraining Loss: 0.02365601\n",
            "Epoch: 3 \tTraining Loss: 0.02678345\n",
            "Epoch: 3 \tTraining Loss: 0.02214155\n",
            "Epoch: 3 \tTraining Loss: 0.02102673\n",
            "primary_caps_total_time 123.54049515724182\n",
            "digit_caps_total_time 63.51781153678894\n",
            "dynamic_routing_total_time 6.178102493286133\n",
            "decoder_total_time 19.780436515808105\n",
            "total_capsule_net_time 209.0126690864563\n",
            "capsule_loss_total_time 3.3438844680786133\n"
          ]
        }
      ],
      "source": [
        "# training for 3 epochs\n",
        "n_epochs = 3\n",
        "losses = train(capsule_net, criterion, optimizer, n_epochs=n_epochs)\n",
        "print(\"primary_caps_total_time\", primary_caps_total_time)\n",
        "print(\"digit_caps_total_time\", digit_caps_total_time)\n",
        "print(\"dynamic_routing_total_time\", dynamic_routing_total_time)\n",
        "print(\"decoder_total_time\", decoder_total_time)\n",
        "print(\"total_capsule_net_time\", total_capsule_net_time)\n",
        "print(\"capsule_loss_total_time\", capsule_loss_total_time)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 68,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 281
        },
        "id": "GWznbtWZKmiN",
        "outputId": "e6248d8e-c613-4a1f-a4bb-d77a37c39278"
      },
      "outputs": [
        {
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEICAYAAABPgw/pAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deZhcdZ3v8fe3qroq6S1Jp7uzdLpJh+wJSzAbKIgQFVCWGXUMg150RsG5MKLIM6J3Ls4weB8dvTrcAZc4jKCyGAEljiiDbCGyZQGE7J3O1lk7nfSSdHqpqu/9oyqd7nR30kk6qZzqz+t5eKyzVJ1fH08+/evv+Z1fmbsjIiLBF8p0A0REpH8o0EVEsoQCXUQkSyjQRUSyhAJdRCRLKNBFRLKEAl0Cz8x+b2Y39ve+IkFjGocumWBm+zst5gKtQCK9fLO7P3z6W3XizOxS4BfuPibTbZGBK5LpBsjA5O75h16b2Sbgc+7+xyP3M7OIu8dPZ9tEgkolFzmjmNmlZlZjZl81s53AT81smJn9l5nVmtm+9Osxnd7zopl9Lv36M2a2xMy+m953o5ldeYL7VprZYjNrMrM/mtn9ZvaLE/iZpqSPW29mK83smk7brjKzVeljbDOzO9Lri9M/Z72Z7TWzl81M/17lqHSByJloJFAEnAXcROo6/Wl6uQI4CNx3lPfPAdYCxcC/Ag+YmZ3Avo8AbwDDgX8CPn28P4iZ5QC/Bf4bKAX+HnjYzCald3mAVImpAJgOPJ9e/xWgBigBRgBfB1QflaNSoMuZKAl8w91b3f2gu9e5+xPu3uzuTcA3gfcf5f2b3f0n7p4AHgJGkQrFPu9rZhXALOAud29z9yXAohP4WeYC+cC30p/zPPBfwPXp7e3AVDMrdPd97r6i0/pRwFnu3u7uL7tueMkxKNDlTFTr7i2HFsws18x+bGabzawRWAwMNbNwL+/feeiFuzenX+Yf576jgb2d1gFsPc6fg/TnbHX3ZKd1m4Gy9OuPAVcBm83sJTO7ML3+O0AV8N9mVm1md57AsWWAUaDLmejInuhXgEnAHHcvBC5Jr++tjNIfdgBFZpbbaV35CXzOdqD8iPp3BbANwN2Xuvu1pMoxvwEWptc3uftX3H0ccA1wu5ldfgLHlwFEgS5BUECqbl5vZkXAN071Ad19M7AM+Cczi6Z7zlcf631mNqjzf6Rq8M3AP5hZTnp449XAY+nPvcHMhrh7O9BIqtyEmX3UzMan6/kNpIZ0Jns8qEiaAl2C4N+AwcAe4DXgD6fpuDcAFwJ1wD3AL0mNl+9NGalfPJ3/KycV4FeSav8PgP/h7mvS7/k0sCldSvpC+pgAE4A/AvuBV4EfuPsL/faTSVbSg0UifWRmvwTWuPsp/wtB5ESohy7SCzObZWZnm1nIzK4AriVV5xY5I+lJUZHejQSeJDUOvQb4O3d/M7NNEumdSi4iIllCJRcRkSyRsZJLcXGxjx07NlOHFxEJpOXLl+9x95KetmUs0MeOHcuyZcsydXgRkUAys829bVPJRUQkSyjQRUSyhAJdRCRLKNBFRLKEAl1EJEso0EVEsoQCXUQkSwQu0Jdu2st3n1lLPKGpoUVEOgtcoL+5ZR/3vVBFS1yBLiLSWeACPRpONblNgS4i0kXwAj2S+l5gBbqISFcBDHT10EVEehLYQG+NJzLcEhGRM0vwAj18KNDVQxcR6SxwgR47VHLRsEURkS4CF+iqoYuI9CxwgR5ToIuI9Chwga4euohIz4Ib6Kqhi4h0EbxA15OiIiI9Cl6gq+QiItKjwAZ6q0ouIiJdBC7QY2HN5SIi0pPABbpKLiIiPVOgi4hkicAFejhkhENGW0KTc4mIdBa4QIfU0EX10EVEugpmoEdCmm1RROQIgQ109dBFRLoKZqCr5CIi0k2fAt3MrjCztWZWZWZ3HmW/j5mZm9nM/mtid7FISA8WiYgc4ZiBbmZh4H7gSmAqcL2ZTe1hvwLgNuD1/m7kkVRyERHpri899NlAlbtXu3sb8BhwbQ/7/QvwbaClH9vXo5gCXUSkm74EehmwtdNyTXpdBzO7ACh3998d7YPM7CYzW2Zmy2pra4+7sYeohy4i0t1J3xQ1sxDwPeArx9rX3Re4+0x3n1lSUnLCx4xGQpoPXUTkCH0J9G1AeaflMel1hxQA04EXzWwTMBdYdCpvjGqUi4hId30J9KXABDOrNLMoMB9YdGijuze4e7G7j3X3scBrwDXuvuyUtBiVXEREenLMQHf3OHAr8AywGljo7ivN7G4zu+ZUN7An0UhYJRcRkSNE+rKTuz8NPH3Eurt62ffSk2/W0ankIiLSXTCfFNVcLiIi3QQy0FPj0DV9rohIZ4EMdPXQRUS6C2agh1Pj0N09000RETljBDPQIyHcIZ5UoIuIHBLYQAd9r6iISGfBDPSwAl1E5EiBDPRYTjrQ9XCRiEiHQAa6eugiIt0FM9DTNXQNXRQROSyQgR7TTVERkW4CGegdo1xUQxcR6RDMQA+HAfXQRUQ6C2agq+QiItJNsAM9oQm6REQOCWaga9iiiEg3wQx0DVsUEekmkIEeU6CLiHQTyEDXTVERke6CGeiqoYuIdBPMQNeDRSIi3QQ70NVDFxHpEMhAj4SMkCnQRUQ6C2SgmxnRSEglFxGRTgIZ6JD+omj10EVEOgQ30CNhjUMXEekksIEei6iHLiLSWWADXTV0EZGughvo4RBtcc22KCJySHADXSUXEZEuAh3ouikqInJYcANdwxZFRLoIbqDrpqiISBfBDnT10EVEOijQRUSyRGADPRbWTVERkc76FOhmdoWZrTWzKjO7s4ftXzCzd8zsLTNbYmZT+7+pXcVyVEMXEensmIFuZmHgfuBKYCpwfQ+B/Yi7n+Pu5wP/Cnyv31t6BI1yERHpqi899NlAlbtXu3sb8Bhwbecd3L2x02Ie4P3XxJ6phi4i0lWkD/uUAVs7LdcAc47cycxuAW4HosBlPX2Qmd0E3ARQUVFxvG3tQsMWRUS66rebou5+v7ufDXwV+Mde9lng7jPdfWZJSclJHS8aDpNIOonkKf9jQEQkEPoS6NuA8k7LY9LrevMYcN3JNKov9L2iIiJd9SXQlwITzKzSzKLAfGBR5x3MbEKnxY8A6/uviT1ToIuIdHXMGrq7x83sVuAZIAz8p7uvNLO7gWXuvgi41czmAe3APuDGU9loOBzorYkEkHOqDycicsbry01R3P1p4Okj1t3V6fVt/dyuY4qF1UMXEekssE+KdvTQFegiIkAWBLp66CIiKcENdJVcRES6CG6gH+qh6+EiEREgGwJdPXQREUCBLiKSNQIb6DGNchER6SLwga4auohISmADPRoOAyq5iIgcEtxAVw1dRKSLLAj0RIZbIiJyZgh+oKuGLiICBDnQ9aSoiEgXgQ30nLABCnQRkUMCG+hmRjQS0jh0EZG0wAY6pOZEV6CLiKQEOtCjkZBuioqIpAU/0NVDFxEBFOgiIlkj2IEeVqCLiBwS7EBXDV1EpEOgAz2mkouISIdAB7pq6CIihwU80MO0quQiIgIEPdB1U1REpEOgAz1VQ9f0uSIiEPBA1ygXEZHDgh3oKrmIiHQIdqBrtkURkQ6BD3T10EVEUhToIiJZItiBHg4RTzrJpGe6KSIiGRfsQNcXRYuIdAh0oMfSga4boyIiAQ/0jh66Al1EJOCBHlbJRUTkkEAHeixHPXQRkUP6FOhmdoWZrTWzKjO7s4ftt5vZKjP7s5k9Z2Zn9X9Tu4uGw4ACXUQE+hDoZhYG7geuBKYC15vZ1CN2exOY6e7nAo8D/9rfDe2JaugiIof1pYc+G6hy92p3bwMeA67tvIO7v+DuzenF14Ax/dvMnh0etqgZF0VE+hLoZcDWTss16XW9+Vvg9z1tMLObzGyZmS2rra3teyt7ceimqIYtioj0801RM/sUMBP4Tk/b3X2Bu89095klJSUnfTyVXEREDov0YZ9tQHmn5THpdV2Y2TzgfwHvd/fW/mne0enBIhGRw/rSQ18KTDCzSjOLAvOBRZ13MLMZwI+Ba9x9d/83s2fqoYuIHHbMQHf3OHAr8AywGljo7ivN7G4zuya923eAfOBXZvaWmS3q5eP6VceDRQp0EZE+lVxw96eBp49Yd1en1/P6uV19osm5REQOC/SToiq5iIgcpkAXEckSwQ50Tc4lItIhKwJdwxZFRAIe6KGQkRM2lVxERAh4oAPEImEFuogIWRDoudEwTS3tmW6GiEjGBT7QK4vzqN5zINPNEBHJuMAH+vjSfKp278fdM90UEZGMyopAbzjYTu3+0zIfmIjIGSsrAh2gavf+DLdERCSzsibQNyjQRWSAC3ygjywcRH4soh66iAx4gQ90M+Ps0nyqahXoIjKwBT7QAcaX5KuHLiIDXnYEemk+uxpbadQDRiIygGVNoINGuojIwKZAFxHJElkR6OXDBhONhDR0UUQGtKwI9Eg4xLjiPPXQRWRAy4pABzR0UUQGvKwJ9PEl+Wzd20xLeyLTTRERyYjsCfTSfJIOGzWVrogMUFkT6BNGpEa6rFcdXUQGqKwJ9MriPEKmoYsiMnBlTaDHImEqinI1dFFEBqysCXQ4/O1FIiIDUVYF+tml+Wzcc4B4IpnppoiInHZZFehTRhbSlkjydk1DppsiInLaZVWgz5s6gtxomIVLt2a6KSIip11WBXp+LMLV547mt3/eTpOm0hWRASarAh1g/uxymtsSLHp7e6abIiJyWmVdoJ9fPpTJIwt47A2VXURkYMm6QDcz5s8q551tDby7TTdHRWTgyLpAB7huRhnRSIjHlm7JdFNERE6brAz0oblRrpo+kqfe3E5zWzzTzREROS2yMtAB5s+uoKk1zg9e2MC72xrYd6ANd890s0RETpk+BbqZXWFma82syszu7GH7JWa2wsziZvbx/m/m8ZtTWcR5Y4Zw3wtVfPTflzDjX57lom89z5a65kw3TUTklIgcawczCwP3Ax8EaoClZrbI3Vd12m0L8BngjlPRyBNhZvzy5gup2r2fmn0HqdnXzLf/sIb//NNG/umaaZlunohIvztmoAOzgSp3rwYws8eAa4GOQHf3TeltZ9QkKoNywkwvG8L0siEAvLOtgceX13DHhyeRH+vLjy4iEhx9KbmUAZ0Hddek1x03M7vJzJaZ2bLa2toT+YiT8pmLxrK/Nc6TK2pO+7FFRE6103pT1N0XuPtMd59ZUlJyOg8NwIyKYZw3ZggPvbKJZFI3SEUku/Ql0LcB5Z2Wx6TXBdJn3juWDbUHWFK1J9NNERHpV30J9KXABDOrNLMoMB9YdGqbdepcdc4oivOjPPTKpkw3RUSkXx0z0N09DtwKPAOsBha6+0ozu9vMrgEws1lmVgN8Avixma08lY0+GbFImL+ecxbPr92tIYwiklX6VEN396fdfaK7n+3u30yvu8vdF6VfL3X3Me6e5+7D3f2MHhd4w5wKwmb8+/PrVUsXkayRtU+KHs2IwkHceNFYfrW8hs8+uJS6/a2ZbpKIyEkbkIEO8I8fmcI9103n1eo6rvp/L/N6dV2mmyQiclIG7NM1Zsan5p7FjIqh3PrIm8z/yWtMGlHA7MoiZo4t4uLxxQzLi2a6mSIifWaZmrBq5syZvmzZsowc+0j7W+M89MomXquuY8XmfRxoS1CcH+XRz89lwoiCTDdPRKSDmS1395k9blOgdxVPJFmxpZ5bHlmBu/PI5+cyUaEuImeIowX6gK2h9yYSDjG7sohHPz8XM+P6Ba+xbldTppslInJMCvRejC/N57Gb5hIOpUK9unZ/ppskInJUCvSjOLskn0dvmosD//PhFRxsS2S6SSIivVKgH8PZJfl876/OY83OJv75t2fsA7AiIgr0vrh0Uim3fOBsHlu6lV+/qal3ReTMNGDHoR+vL8+byNJN+/j6k+8yffQQivNj7G1uo7aplbU7m3h3WwPvbGsg6c4tHxjP1eeOJhSyLp/h7myua+atrfX8uaaBwdEQE0cUMHlkIZXFeUQj+v0qIidOwxaPw67GFq6692XqDrR12zY8L8r0siHsbmpl9Y5Gpo0u5B+umEx+LMLSTXt5Y+NeVmzZR31zOwCDckK0J5xEei6ZIYNzePhzczq+XUlEpCcah96PVu9o5Pfv7mTo4ByK8qIU5UWZMCKfkYWDMDOSSeept7fx3WfWsa3+YMf7xpXkMeusIs6vGMr55UOZUJpPwp3q2gOs3dnEd55ZSzyZ5Klb3sfIIYNO6c/g7vz6zW3MrixizLDcU3osEelfCvQMaI0n+O3bO8iLhpk5toiSgthR91+9o5GP//AVxpXks/DmCxkcDZ+ytv3s1U3c9dRKivOj/PQzszlnjP4qEAkKPViUAbFImI+/ZwxXnjPqmGEOMGVUIffOn8G72xu4feFbPU7r6+60tCc6yjQnYs3ORu753WrmVBYRi4T55IJXeXHt7hP+PBE5c6iHfob5yeJqvvn0asYV52EGbYkkbfEkza0JDrTFSTrkxyLMqBjKrLFFzKgYSm708L3tQTkhSgpiDM+LET7ipmxLe4Jr7lvC3gPt/OFLF5NMOp99cClrdjZxz3XTmT+rHDM7skkkkt7ts0QkM47WQ9colzPM5y6uJJ50VmzZRzQcIhoJEQ2HyI2FyYtGGBwNs7OhhaWb9vL9P66jt9/HIYPi/BiXTCzhL2eUMWfccO753SrW7drPz/5mNsX5qb8afnnzhfzdL5bztSffYeGyrXztyinMriwimXSeXb2LnyyuZt2uJn71hYuYNPLoc9q8vbWeIYNzGFuc19+nRUT6QD30AGtobmfljgbaE4f/P2xujbNnfyu1Ta1s3tvMc6t3s781TklBjNqmVm66ZBxfv2pKl8+JJ5I8saKG7z27jl2NrVwysYSte5vZuOcAY4YNpqU9QX4swlO3vo8hg3O6tcPd+fHiar79hzWMKBjE72+7WFMPi5wiuik6gLW0J/jj6l38esU2zOAHN7yn1/HuB9sSPPjKJhYs3kB5US43XTKOK6aN5K2t9cxf8BoXTyjmgRtndRlf39Ke4M4n/sxv3trOpZNK+FPVHi6dVMqCT7+nx/KNiJwcBbqctJ+/tpn//Zt3+eJl47n9Q5Oob27j9Y17+cELVbxd08AdH5rILR8YzwNLNnLP71bzL9dN59NzzwKgqaWd+16oYk9TGx+eNoJLJpYwKCeMu7N+934Wr6ulZt9BYpEQOeEQebEIH3tPGaUFp3b4pkgQKdDlpLk7X33izyxcVsOkEQWs292EOxQMivB/P3EeH5o2EoBk0vnMg0t5vbqORbe+j811B7jrqZXsamqhIBahsSVOXjTMrMoi1uxoYmdjC5C60dueSNKWSOIOIwpj/PBT7+GCimGZ/LFFzjgKdOkXLe0Jbn1kBc1tCS4cN5y5Zw/n3DFDiEW6jpmvbWrlynsX09KeZH9rnMkjC/jWx85l2uhCXquu4+l3dvJ6dR1TRhVy8YRiLp5YQtnQwR3vX72jkZt/vpwdDQe5+9rpXD+7oltbttQ18/Drm9ne0MLNl4zLmidsdza08OyqnVx1ziiG5x97uKsMPAp0Oe1eXl/LHb96mxsvGsvnLx5HTvj4Hnmob27ji4+9xeJ1tVw+uZQJIwooKYiRHwvz+3d38tK6WkJm5EXDNLbE+YsZZdzx4UmMLBzExj37Wbm9kZp9BxmUEyYvGiY3FqFyeB5TRhUQSbclmXTe2dbAs6t2dfylYKSGae5uamVnYwu7GloYPXQwf3txJdedX3ZK59vZ2dDCJxe8yua6ZmKREH95QRl/897KXr8G8dC/Xd2rGFgU6BJIiaRz7x/X8fjyGvbsb6MtkQRS5ZjrZ1cwf1YFubEwP3pxAw8s2YiTGq7Z0p7s9TNzo2HOLx9K2dDBvLx+DzsbWwiHjBGdHv4yM0oKYowsHMSIwhivb9zLmp1NjCwcxA1zKki4s2nPATbuOUBtUyst8SQH2xK0JZLkRcMMzY0yNDeHQTlh4okk7QknnnTGDs/lvPLU1A/nlA0hL3Z41PDuxhY+ueA1apta+dbHzuGVDXU8sbyG1niSD08bwVevmMy4knwgFeTPrNzF/3l6NaUFMe6/4QJGFPZ+vyGRdN7YuJcLzhra7a8pCR4FugSeu9N4ME7dgVbKi3K79fi31R/kJ4urCZkxbXQh08oKGTs8j9Z4kua2OPtb4qze2cTyTXtZtnkfW/Y2c9HZw/nwtJFcNrmUobm9D7N0d15aV8uPX6rm1eo6zGD0kMGMK8mjtGAQg6MhBkXCRCMhmtsS1De3UX+wnYNtCXLCIXLChpmxfncTW/em5veJRkJcMqGEj5w7kvPLh/G5h5ayo6GFn/3NbGaOLQJg74E2HnplE//xcjUt8SR/PbuCq88bzb3PreNPVXWML81ne/1B8mMRfvTpnu83rNzewNeffIe3axqYNjr1NPL40tQvhvZEkh+9uIEfvbSB2ZVF3HnllGM+ayCZp0AX6Se7G1soHJzqfZ+Iuv2tvF1Tz5L1dfz+3R3saEiVegbnhHnws7OYM254t/fUNrVy73PrePSNrSSSzpDBOdz+wYncMKeCqtr93PSz5exsaOF/Xz2VuZVF5IRDhEPGz1/bzANLNjIsN4dPzx3Lg69s5GB7gm9cPY2powr56hN/Zs3OJi6eUMxbW+s50BrnE+8p57oZZaze0chbW+tZvaORSSML+IsZZVwysaRPpbPWeIJ3tzVSUZR7zGkvmlra2VZ/kImlBd2mm35ieQ3f+sMaLp9cyi0fGE950YlNJLdnfytb9zZzTtmQjnLb8WppT3D7wrdobktw2+UTmJHBm/UKdJEzUDLpvLm1nmdX7eLyKaXMSvfMe7Ohdj9/qtrDR88dTVGnB7fqm9v4+0ff5OX1e7q95/rZ5dx5xRSG5Oawq7GFryx8myVVqf1GFMa457pz+ODUEew70MZ9L1Txs1c3dTyoNmrIICaNLODtrfXsa26nKC/KZZNLGTNsMCMKB1GaDusDbQmaW+PsaGjpmCa6NZ4kPxbhq1dO5obZFR1h3RZP8uyqXSypquXNLfWs3ZUaLXXumCHceeVkLjq7mINtCb6x6F0WLqth4oh8NtU1k0w6f3lBGTfMOYsxwwZTlBftdu8gmXTaEkla46mb8S+treV372zn1Q11JD01RfVlk0uZN2UEkbBRXXuA6tr97GtuY9LIAs4pG8q5Y4YwasigLp/dFk/yhV8s5/k1uxmam0N9czvzppTy5Q9OZNro7jfjdza08PL6Wt4/sYTSo5TCTpQCXSTLxRNJ/rShjqaWdtoTSdrjzqSRBZxXPrTLfsmk8+Arm9hef5AvzptA4aCuT/5u3dvMul1NTC8b0lGXb4snWbyull+/tY3XNtT1+H0AAGYwbXQhcyuHM6NiGI++sYUlVXuYPbaIL82bwEvra3l8WQ11B9ooHBRhRsUwZlQMZejgHBYsrmZ7QwuXTiphR30L63Y3ccul4/nSvAnUHWjjhy9u4NE3ttAaT90fiYZTcxYl3TnYnqC5LUFbvPu9k3HFeXzk3FGML83npbW1PL92d8d3EgCUFMQYOjiH6j0HOia9m15WyG2XT2TelFLiSefWR1ak7ln8xTlcc/5oHvzTRn68uJqmltQIrssml3L5lFKaWuI8/PoWnl+zm0TSyY9F+PIHJ3LjhWed8F8GPZ9nBbqI9JO2eJI9+1vZ3dSKAXmxMLnRCEMG53S50evuPL68hnt+t5qGg+2EQ8blk0v56zkVXDKhpNsTxw++son7X6giJxzi+588n/dPLOly3N1NLSzftI+djS3sbGxhd2Mr4ZCRGw0zOCfMoJwwsZwQsUiYWCTEBRXDmDKqoEtvO55I8nZNPZFQiMqSvI5faC3tCVbvaGTFlnp+9uomNtc1M210IaUFMV5YW8s3rp7KZ99b2fE5Dc3t/Gr5Vp5bvZulm/YST/8yKM6P8omZ5bx/Ygk/fHEDL62rZdKIAm6bN4FJIwsoH5Z70iOlFOgikjG1Ta28tK6WiycUH3U0DqRq6ofKI5kSTyT5zVvbue/59Wyqa+ZrV07m5vef3ev+DQfbWbJ+D+EQXDZ5REdguzv/vWoXd/92VceX3YQMyoYN5o4PTeLa88tOqH0KdBGR4xRPJNm8t5mz08NFT1RLe4JVOxrZtOdAarhrXTPzZ5Xz3vHFJ/R5mj5XROQ4RcKhkw5zgEE5YS6oGHZaprHQNxaJiGQJBbqISJZQoIuIZAkFuohIluhToJvZFWa21syqzOzOHrbHzOyX6e2vm9nY/m6oiIgc3TED3czCwP3AlcBU4Hozm3rEbn8L7HP38cD3gW/3d0NFROTo+tJDnw1UuXu1u7cBjwHXHrHPtcBD6dePA5ebJmkWETmt+hLoZcDWTss16XU97uPucaAB6D5tnIiInDKn9cEiM7sJuCm9uN/M1p7gRxUD3aeWE9C5ORqdm97p3PTsTDwvZ/W2oS+Bvg0o77Q8Jr2up31qzCwCDAHqjvwgd18ALOjDMY/KzJb19ujrQKdz0zudm97p3PQsaOelLyWXpcAEM6s0sygwH1h0xD6LgBvTrz8OPO+ZmiRGRGSAOmYP3d3jZnYr8AwQBv7T3Vea2d3AMndfBDwA/NzMqoC9pEJfREROoz7V0N39aeDpI9bd1el1C/CJ/m3aUZ102SaL6dz0Tuemdzo3PQvUecnY9LkiItK/9Oi/iEiWUKCLiGSJwAX6seaVGUjMrNzMXjCzVWa20sxuS68vMrNnzWx9+n9P/cz6ZyAzC5vZm2b2X+nlyvRcQ1XpuYeimW5jJpjZUDN73MzWmNlqM7tQ10yKmX05/W/pXTN71MwGBem6CVSg93FemYEkDnzF3acCc4Fb0ufjTuA5d58APJdeHohuA1Z3Wv428P30nEP7SM1BNBDdC/zB3ScD55E6RwP+mjGzMuCLwEx3n05qVN98AnTdBCrQ6du8MgOGu+9w9xXp102k/mGW0XVunYeA6zLTwswxszHAR4D/SC8bcBmpuYZg4J6XIcAlpIYa4+5t7l6PrplDIsDg9AOSucAOAnTdBC3Q+zKvzICUnrJ4BvA6MMLdd6Q37QRGZKhZmfRvwD8AyfTycKA+PdcQDNxrpxKoBRe+Yt4AAAGrSURBVH6aLkf9h5nloWsGd98GfBfYQirIG4DlBOi6CVqgSw/MLB94AviSuzd23pZ+YndAjU01s48Cu919eabbcgaKABcAP3T3GcABjiivDMRrBiB93+BaUr/0RgN5wBUZbdRxClqg92VemQHFzHJIhfnD7v5kevUuMxuV3j4K2J2p9mXIe4FrzGwTqbLcZaTqxkPTf0rDwL12aoAad389vfw4qYAf6NcMwDxgo7vXuns78CSpaykw103QAr0v88oMGOm68APAanf/XqdNnefWuRF46nS3LZPc/WvuPsbdx5K6Rp539xuAF0jNNQQD8LwAuPtOYKuZTUqvuhxYxQC/ZtK2AHPNLDf9b+vQuQnMdRO4J0XN7CpS9dFD88p8M8NNyhgzex/wMvAOh2vFXydVR18IVACbgb9y970ZaWSGmdmlwB3u/lEzG0eqx14EvAl8yt1bM9m+TDCz80ndLI4C1cBnSXXuBvw1Y2b/DHyS1AiyN4HPkaqZB+K6CVygi4hIz4JWchERkV4o0EVEsoQCXUQkSyjQRUSyhAJdRCRLKNBFRLKEAl1EJEv8f0HZ56P14xIEAAAAAElFTkSuQmCC",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "needs_background": "light"
          },
          "output_type": "display_data"
        }
      ],
      "source": [
        "plt.plot(losses)\n",
        "plt.title(\"Training Loss\")\n",
        "plt.show()"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 63,
      "metadata": {
        "collapsed": true,
        "id": "KAPp9ZQwKmiO"
      },
      "outputs": [],
      "source": [
        "def test(capsule_net, test_loader):\n",
        "    '''Prints out test statistics for a given capsule net.\n",
        "       param capsule_net: trained capsule network\n",
        "       param test_loader: test dataloader\n",
        "       return: returns last batch of test image data and corresponding reconstructions\n",
        "       '''\n",
        "    class_correct = list(0. for i in range(10))\n",
        "    class_total = list(0. for i in range(10))\n",
        "    \n",
        "    test_loss = 0 # loss tracking\n",
        "\n",
        "    capsule_net.eval() # eval mode\n",
        "\n",
        "    for batch_i, (images, target) in enumerate(test_loader):\n",
        "        target = torch.eye(10).index_select(dim=0, index=target)\n",
        "\n",
        "        batch_size = images.size(0)\n",
        "\n",
        "        if TRAIN_ON_GPU:\n",
        "            images, target = images.cuda(), target.cuda()\n",
        "\n",
        "        # forward pass: compute predicted outputs by passing inputs to the model\n",
        "        caps_output, reconstructions, y = capsule_net(images)\n",
        "        # calculate the loss\n",
        "        loss = criterion(caps_output, target, images, reconstructions)\n",
        "        # update average test loss \n",
        "        test_loss += loss.item()\n",
        "        # convert output probabilities to predicted class\n",
        "        _, pred = torch.max(y.data.cpu(), 1)\n",
        "        _, target_shape = torch.max(target.data.cpu(), 1)\n",
        "\n",
        "        # compare predictions to true label\n",
        "        correct = np.squeeze(pred.eq(target_shape.data.view_as(pred)))\n",
        "        # calculate test accuracy for each object class\n",
        "        for i in range(batch_size):\n",
        "            label = target_shape.data[i]\n",
        "            class_correct[label] += correct[i].item()\n",
        "            class_total[label] += 1\n",
        "\n",
        "    # avg test loss\n",
        "    avg_test_loss = test_loss/len(test_loader)\n",
        "    print('Test Loss: {:.8f}\\n'.format(avg_test_loss))\n",
        "\n",
        "    for i in range(10):\n",
        "        if class_total[i] > 0:\n",
        "            print('Test Accuracy of %5s: %2d%% (%2d/%2d)' % (\n",
        "                str(i), 100 * class_correct[i] / class_total[i],\n",
        "                np.sum(class_correct[i]), np.sum(class_total[i])))\n",
        "        else:\n",
        "            print('Test Accuracy of %5s: N/A (no training examples)' % (classes[i]))\n",
        "\n",
        "    print('\\nTest Accuracy (Overall): %2d%% (%2d/%2d)' % (\n",
        "        100. * np.sum(class_correct) / np.sum(class_total),\n",
        "        np.sum(class_correct), np.sum(class_total)))\n",
        "    \n",
        "    # return last batch of capsule vectors, images, reconstructions\n",
        "    return caps_output, images, reconstructions"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 64,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "XNiFrb9BKmiP",
        "outputId": "d4617ca1-d80d-42ed-a8a5-1059157cd994"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "Test Loss: 0.03245602\n",
            "\n",
            "Test Accuracy of     0: 99% (977/980)\n",
            "Test Accuracy of     1: 99% (1132/1135)\n",
            "Test Accuracy of     2: 99% (1026/1032)\n",
            "Test Accuracy of     3: 99% (1001/1010)\n",
            "Test Accuracy of     4: 99% (976/982)\n",
            "Test Accuracy of     5: 99% (890/892)\n",
            "Test Accuracy of     6: 98% (942/958)\n",
            "Test Accuracy of     7: 98% (1014/1028)\n",
            "Test Accuracy of     8: 99% (969/974)\n",
            "Test Accuracy of     9: 97% (985/1009)\n",
            "\n",
            "Test Accuracy (Overall): 99% (9912/10000)\n"
          ]
        }
      ],
      "source": [
        "# call test function and get reconstructed images\n",
        "caps_output, images, reconstructions = test(capsule_net, test_loader)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 65,
      "metadata": {
        "collapsed": true,
        "id": "q0U7P3xDKmiQ"
      },
      "outputs": [],
      "source": [
        "def display_images(images, reconstructions):\n",
        "    '''Plot one row of original MNIST images and another row (below) \n",
        "       of their reconstructions.'''\n",
        "    # convert to numpy images\n",
        "    images = images.data.cpu().numpy()\n",
        "    reconstructions = reconstructions.view(-1, 1, 28, 28)\n",
        "    reconstructions = reconstructions.data.cpu().numpy()\n",
        "    \n",
        "    # plot the first ten input images and then reconstructed images\n",
        "    fig, axes = plt.subplots(nrows=2, ncols=10, sharex=True, sharey=True, figsize=(26,5))\n",
        "\n",
        "    # input images on top row, reconstructions on bottom\n",
        "    for images, row in zip([images, reconstructions], axes):\n",
        "        for img, ax in zip(images, row):\n",
        "            ax.imshow(np.squeeze(img), cmap='gray')\n",
        "            ax.get_xaxis().set_visible(False)\n",
        "            ax.get_yaxis().set_visible(False)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 66,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 295
        },
        "id": "FoVaC7piKmiS",
        "outputId": "a8afa5cf-50f1-4a3e-ec54-805bd282c2d9"
      },
      "outputs": [
        {
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAABbwAAAEhCAYAAAC0kWbAAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOzdd7hV1bX38YGFJhGkKUoTUBCx0ATFCPYoSYxiidfyxJioSSxRo0mMphm9JrHEFCXRGB+jubZ4bVETG1ZAkSKIAtJBQBEEBCzoef94n7xvxm8Oz1p7s/cp63w//43lnGsv9pp7rrmm5xmjWU1NjQEAAAAAAAAA0NhtUd8XAAAAAAAAAABAJbDhDQAAAAAAAAAoBDa8AQAAAAAAAACFwIY3AAAAAAAAAKAQ2PAGAAAAAAAAABQCG94AAAAAAAAAgELYqpTGzZo1q6nWhaDBWVlTU9OpWidnLDUpVR1LZoynJoa5CZXC3IRKYm5CpTA3oZKYm1ApzE2oJMYTKqampqZZdJy/8MZnWVjfF4DCYCyhkhhPqBTGEiqJ8YRKYSyhkhhPqBTGEiqJ8YSqY8MbAAAAAAAAAFAIbHgDAAAAAAAAAAqBDW8AAAAAAAAAQCGw4Q0AAAAAAAAAKAQ2vAEAAAAAAAAAhcCGNwAAAAAAAACgENjwBgAAAAAAAAAUwlb1fQEAStOtWzcXX3PNNUmb4cOHu3jEiBFJm8WLF1f2wgAAAAAAAIB6xl94AwAAAAAAAAAKgQ1vAAAAAAAAAEAhsOENAAAAAAAAACgEcngDjYzm7D7uuOOSNuPHj3cx+bqbhiFDhrj45Zdfrsh577rrLhd/9atfrch5AQAAAABA3dhii/Tvns8880wX33DDDS5etGhR0mfkyJEuXrBgweZfXIXxF94AAAAAAAAAgEJgwxsAAAAAAAAAUAhseAMAAAAAAAAACoENbwAAAAAAAABAIVC0Emjgzj//fBdHRSrVCSecUK3LQQPRtWvX5Nhtt93m4pqamop81sEHH1yR8wAAAAAAgLrRoUMHF//yl79M2px22mku/vTTT10c7T3oMYpWAgAAAAAAAABQJWx4AwAAAAAAAAAKgQ1vAAAAAAAAAEAhNNoc3j179kyOPfHEEy7u1auXi9euXZv0ufzyy118zTXXbP7FAWXq1q1bckxzeKsLLrggObZ48eKKXRMapptvvjk51q9fv8x+q1evdvH111+fed6+ffuWeHWoL1tvvXVybNCgQS6+9NJLXXzkkUcmfbbYwv//cM3jZma2cOFCF+vz9NFHH036rFu3zsXr169P2qBxa9GiRXLshRdecPHAgQOTNkuXLnVx9+7dK3thqFPNmjVzcZT78cc//rGLTz/99JLPG9WqeP7551188cUXu3jChAmZn4OGLXrW9ejRw8U/+clPXHzyySdnnnfGjBnJsZ///Ocuvvfee11cqXopQETX6f3790/aHHrooXV1OQAauPbt2yfHnnrqKRcPGDCg5PNOnDgxOTZnzpySz1PX+AtvAAAAAAAAAEAhsOENAAAAAAAAACgENrwBAAAAAAAAAIXQIHN4R3loDznkEBefc845SZudd97ZxZpTTXOSmpm1bdvWxWeffbaLR44cmfQ57LDDkmPVkCfn+J/+9Kc6uRbUjX333Tc5pnm9NT/3ddddV9VrQsNw4YUXuviggw4q6zyaE/62227L7PPWW2+V9Vmoe9Ez67HHHqu1T5R/VHN2R200x/JNN92UeX1Tpkxx8cyZM5M2t956q4uffvrpzPOi4Yhywu+9994ujsYTeXAbt6OPPtrFo0ePdvFpp52WeY48YyBPmxEjRrj4hz/8oYuPOuqozHOgYWnVqpWL77vvvqRN1vtZVItCRfmR77zzThdrXZM333wz87xAuY499lgXR+smAE2X5uPW9ygzsz322MPFedZSmrP7mGOOSdqsWLEixxXWL/7CGwAAAAAAAABQCGx4AwAAAAAAAAAKgQ1vAAAAAAAAAEAhsOENAAAAAAAAACiEBlG08re//a2LTzrppKRNu3btNvtzWrZsmRz77ne/6+JtttnGxc2aNUv61FVhpTZt2iTHbrjhBhcPGzbMxaeffnpVrwmVpQUpr7766sw+FKlsGrp27eris846y8VbbZU9fa9atSo5RrGbYtFCxl/72tfq50JyGjRokIsHDhyYtOncubOLKVrZsG2//fYu/uMf/1hPV4K6cvzxxyfH/vznP7u4devWm/057733XnJM1/LR2l4NHTrUxVrgycxsxowZJV4dqqlXr14u/uc//1nrf69LF198sYvPOOOMeroSlGP48OEufuCBB5I2d9xxh4svuOCCql7Tv11xxRXJsS5duriYdXzDou9jRxxxRNJmzJgxLtZ7unz58qTPvffe6+JZs2YlbWbPnu1iLSa/1157BVecTQs/61p92rRpSR+dF1E9O+ywg4vHjRvn4u22266s82YVqYzGaWPAX3gDAAAAAAAAAAqBDW8AAAAAAAAAQCGw4Q0AAAAAAAAAKIR6yeE9atQoF3/zm990cfPmzavyuVtuuWVyTHN2a16sOXPmJH3y5PDWPMvvvPOOizV3s5nZRRdd5OJDDjkkaaM5xTVfKzm8GxfNIR+Ni/Hjx7uYHN5Ng+YY7d27d8nnWLJkSXJs0qRJZV8T6t/uu+/u4ksuucTFdVVjopo0t6Hmolu9enVdXg4y6Jqtffv29XQlqJZ99tnHxTfffHPSphI5ux966CEXX3rppZmf88gjjyRtdM7QPPO33XZb0kfrC6B6NOetvgeamZ1//vkurkTO7gkTJiTHpk+fnnktqlWrVi6O3i8/+eSTEq8OdeWHP/yhizt16pS00ffr66+/3sULFy6s+HWZxXVN0LD06dPHxVpL54QTTqjI55xyyikujmparFmzxsVa/61Dhw4VuRZ16KGHJsfI4V13NE98npzd+ky67777kjbnnnuui1esWFHG1TU8/IU3AAAAAAAAAKAQ2PAGAAAAAAAAABQCG94AAAAAAAAAgEJgwxsAAAAAAAAAUAj1UrSyElatWpUcmzJliovvuOMOF7/22muZ5126dKmLly1bVsbVZZs1a1ZyTAvM5bleNB7HH398cuyCCy5w8eLFizPbAHlpkR00fg888EB9X0JuN954Y3JMi98NHjw4aaMFW3faaScXU7SyYdPi2nn95je/qfCVoFJ0HaIF3/P44IMPkmNnnHGGi7Vo5dq1azPPO23atOTYqFGjau3Trl275JgWy63W+r+p0QKVZmbf+973XHzFFVdU5LNefvllF990000u1vFlZjZ69GgX5ylaqQXmW7RokbTZsGFD5nlQN7To++GHH57Z5/XXX3ex7g9USr9+/Vx88MEHJ20+/fRTF48dO7Yq14J8evbs6eI999yzTj43em5Fx+pCtA+H6tDnjVlaXDKPjRs3urhSxVUbA/7CGwAAAAAAAABQCGx4AwAAAAAAAAAKgQ1vAAAAAAAAAEAh1EsO73Hjxrn4Zz/7mYuPOuqopI/mdoxyYE+dOnXzL66OtGnTJjmmOe3y0LzfaDiuueYaF+fJxR2NgQkTJlTsmtC0nHjiicmxv/zlL/VwJchD84D+8Y9/TNrsvPPOLt5iC///rTXXYx5RrtF//OMfLp40aVLSZsGCBS6+9957Mz9rhx12cHGevJj6b0TDVlNTU9+XgArr1KlTyX0mTpzo4mOPPTZp89Zbb5V9Tf9WTm7dKAd527ZtXUwO78qIxk4lcnbrGtvM7LrrrnNxnnt49tlnl/zZs2fPdjH5uhu2oUOHuljXWtEzS/O9b9q0qfIXZmZ33nmni7feeuukzcqVK13897//vSrXgnyeeOIJF2ttmnXr1pV8zmiNHdX1KtX777+fHLv11ls3+7yVuDbkE92vrLzxr776anJszJgxlbqkRoe3SAAAAAAAAABAIbDhDQAAAAAAAAAoBDa8AQAAAAAAAACFUC85vNVVV11Va1xE/fv3T46ddtppmf00P+svfvGLil0TNs/555/v4jw5u7XN3XffXdFrQuM1fvx4F2vOvjy5uEaMGJEcO+uss1w8duzYMq4O1aA5bk8++eSkjeaa1GdClEt0zpw5Lr7ttttc/OyzzyZ9Jk+eXPvFVkiefM/l5A8GUDnnnHOOi0844YTMPr/61a9cvH79+ope078deuihJfeJ8n6/8cYblbgciPfeey85dscdd7h49OjRSRvNg/vnP//ZxTfffHPSJytn9+GHH54c69mzZ6190PiddNJJLtZ1x4cffpj0efTRR6tyLfr+37dv38w+t99+e1WuBZUxfPjwkvtobbpLL700aUNtgKbp+9//vov333//ks/x8MMPJ8fmzZtX9jU1dvyFNwAAAAAAAACgENjwBgAAAAAAAAAUAhveAAAAAAAAAIBCYMMbAAAAAAAAAFAIDaJoZUP2uc99Ljm27bbblnyetWvXulgLaOSlBVluuummss6DzXP88ccnx6699tpa+0T//brrrqvYNaFY1qxZ4+IrrrjCxQcffHDSp127di5u1apV0uaaa65xsRaLOuOMM5I+WYWgUBmnnHLKZp/j8ccfT44dc8wxm33eSunRo0fJffR7efLJJyt1OWhA6qpQKko3c+ZMF//kJz+ppytJn2tbbJH9tzta3Pfdd9+t6DXhs23cuDE5duqpp7p48ODBSZvly5e7OCo0mkXf17SQqlm6boqsWrXKxX/7299KvhbUn+hd/j998MEHybFp06Zt9udus802ybGf//znLm7RokXmeWbPnr3Z14LK0fnqwQcfLPkcBx54oIsPOeSQpE0550Xj0717dxefeeaZLt5qq+zt2kmTJrn4qquu2vwLM7Ojjz7axf369UvaaDHMu+66qyKfXUn8hTcAAAAAAAAAoBDY8AYAAAAAAAAAFAIb3gAAAAAAAACAQmjyObybNWvm4mOPPdbF559/ftJn2LBhJX/OhAkTXDx8+PDMPjU1Ncmxyy+/vOTPxubT+3X11Vdn9tGc3RdeeGFFrwlNy9SpU118++23J21OPvlkF0e5KTX/qebwvv/++5M+Rx11lIs1tyZKt+uuuybH8jwXlNYB0FzvDc2YMWPq+xJQYd/97nddrOuqyPvvv58ce+aZZyp2TSiGKPfuX//6Vxd37Ngx8zya//nQQw/dvAtDRb3yyisVOc+AAQNcfPPNN9f63/P605/+5OJx48aVdR5U3/e+973k2IgRI0o+T/PmzWv976NGjUqOaX7bc845J2nTu3fvWs87f/785Nj//M//1NoHdUtzLkf1krLstddeLo7e6e6++24Xf/vb307afPTRRyV/NurPlltumRzT51TPnj0zz7Np0yYXf/nLX3ZxtA6/8sorXRytgzQ/vZ4n2ptUUY2Xz3/+8y6u6zoq/IU3AAAAAAAAAKAQ2PAGAAAAAAAAABQCG94AAAAAAAAAgEJoUjm827Ztmxy77LLLXKw5u6McOHny16hycrO+9dZbybGbbrqp5POgdN26dXOx5tHS/25mNn78eBeXk7M7Oq/mR9U2xx13XOa17LfffiVfCxq2c889Nzmm9/2OO+4o+bxDhw5Njn3lK19x8dixY0s+LzzNZ2Zm1qZNm5LPo/nUV69eXfY11YVofKFY8qyRNCcuENG6FGZmX/rSlzL7ffrppy7+wx/+ULFrQsOlOZMr9bx56qmnKnIeVN8XvvCF5FiUN/c/RfsDH3zwQcWuqRS///3vk2Nr1qyphyvBZ/nGN75R8XNG6/+vf/3rLt52222TNvqMJKd3w7bvvvsmxw4++GAX6xo6mou0VsF5553n4mgPauutt671c6Jj5ex59u3bNzm23XbbuZgc3gAAAAAAAAAAlIENbwAAAAAAAABAIbDhDQAAAAAAAAAoBDa8AQAAAAAAAACF0GiLVvbs2TM5duKJJ7r46aefdvENN9yQ9Nl7771r/Zyrrroq81oGDBiQHPviF7+Y2U9pgcxLL7205HOgMrIKRS5evDjpc8IJJ5T8Occff7yLr7766qRNVMgyy5IlS0rug/J06NDBxb/5zW+SNvPnz3fxFVdc4eIPP/ywItfyzDPPuPiNN95I2vTr168in4XNM3LkyORYVCRZbdiwwcXr16+v2DXVhVGjRrlYi8tFdFwDKKajjz7axf/93/9d1nm0wPuvf/3rsq8JuP32212s744rVqyoy8tBgSxbtszFt9xySz1dCfJau3ZtyX1ef/11F++yyy4u3mqr7C25Y489Njn297//3cV33XVXydeG6mnevLmLL7nkkpLPsXz58uSYro0+97nPlXzeTz75JDmm+xUtWrRwcZ49qZdffjk59vbbb5d4dZXFX3gDAAAAAAAAAAqBDW8AAAAAAAAAQCGw4Q0AAAAAAAAAKIRGm8P75JNPTo797Gc/c/FHH33k4o8//jjpc80117j45z//uYvXrVuXeS2nnnpqcmz06NGZ/VRNTY2Lhw0blrS57bbbSj4vaqd5tM3MLrjgAhePHz/exVG+7iiv93+K8h5pzu4o9/Y999xT67VFn3vhhRfWei2onG9/+9suPumkkzL7bNq0ycU67+TRunXr5Nh+++3n4nLyda9evTo59vjjj5d8HtRO5/vPOqa0FsWNN95YsWuqC5qzO/o3z5gxw8X3339/Va8Jm+cb3/hGfV8CGql27dq5+LLLLnNxOXkpzcp7pqLx+93vfufiY445xsXt27cv67ydO3d28be+9S0X//SnPy3rvKi8Rx99NDnWpUuXzT7v9OnTXTx79uykzfe+9z0Xt2zZMvO8+j65Zs2aMq4Odek73/mOi6+77rrMPlk5vE855ZSkz3nnnVfytej72qpVqzLPgerZZ599XHz44YeXfI6obmGWqK7EokWLXPyLX/wiaaPj51e/+pWLzz777MzPnjhxYnKsnLz3lcRfeAMAAAAAAAAACoENbwAAAAAAAABAIbDhDQAAAAAAAAAoBDa8AQAAAAAAAACF0GiKVj7zzDMujgo6qhdffNHFv/71r5M2jz32WMnXooXgqlUc55vf/GZyTIvd5SloAE+LR2rhyMhvfvMbF0eFIvW8xx57rIuPO+64pE+eQhfnn3++i6+99tpar+2zrg/VUc5vcPfddy+5z4477ujiH//4x0mbM844o+Tz6pyiRRHNzObOnVvyeYFozstj2bJlLo4KqaL+jBw50sXbbruti7UwaaRZs2YVvSY0fNtss01y7H//939dvNdee5V83ui5FxVsQvFpweN9993XxYMGDUr6jB071sVt27bN/Bxdj6HhuOaaa3Id21wdOnRIjmnR+p133jlpo4XjJk+eXNkLQ9VpIciXXnqp5HPofd96662TNnneL/fff38Xd+/e3cUUraxfp512Wr18rhafNEuLgi9fvjxpo3tiWhQ1cvHFF7v4pptuynOJdYq/8AYAAAAAAAAAFAIb3gAAAAAAAACAQmDDGwAAAAAAAABQCA0yh/ePfvSj5Ng+++zj4ijX0fTp0138la98xcXr1q0r+VqiXLsPPPCAizV3cx4zZ85MjvXv39/FW265ZdLm61//uov/8Ic/uHj27NklX0tTo7ncovunebLvvvvuzD533XWXizV34Pjx45M+mp87Ou8999zj4gsvvDBpg/qjOUhPP/30ks/RqlWr5JjWKbjzzjtd3Llz55I/xyzNh6w5uzXHF5BXr169XHzFFVdk9onyc//yl7+s2DWh8jQPrubsrqmpyTzHP/7xj4peExq+KOftAQccsNnnffrpp5Njmi/8ww8/dPHHH3+82Z+Lhu/NN9+sNTYza9eunYtvvPHGzPOOGTPGxZdccknSZuXKlXkuEY3UiSeemByLcnYrrd20YcOGil0Tajd48GAXX3TRRUmbc845x8XvvPNOVa8JqIaTTz45OabrrSiHt+63qkmTJiXHtB7i2rVr81xineIvvAEAAAAAAAAAhcCGNwAAAAAAAACgENjwBgAAAAAAAAAUQoPI4X388ce7OMohqzm7V61albTRfpqzu2XLlkmf/fff38WHHXaYi6McOJ06dUqOZTnttNNc/Pjjjydtvv/977t44MCBSRu93rFjx7r4oIMOKvnaikzHlpnZcccdl9lvyZIltZ7n6quvTvpk5XLXnN5mZosXL3bxCSeckLTR/OFoWF588UUX62/dzGyLLfz/W9Q8cvfdd1/S5/DDDy/5WjZt2uTiWbNmJW20tsHcuXNL/hzUnyFDhtT3Jfw/Opdqzm7N6R2Jxv64ceM267rQ8M2ZM6e+LwFV1r59exdr/ZtKyTOWNGd8lH93/fr1FbsmlEbfbczS54c+FxYtWlSRzy4nR6/m/d5qqwbxOo06tPfee2e20Xc8M7O//vWv1bgc5KC5hqO6Ej/4wQ9cXFc5vLWOG1BpPXr0cHH37t0z+7z33nsuPvPMM5M2r7322uZdWB3gL7wBAAAAAAAAAIXAhjcAAAAAAAAAoBDY8AYAAAAAAAAAFAIb3gAAAAAAAACAQmgQVTa0OKMWqIw89NBDmcf69Onj4nPOOSfpc/bZZ9f6Oc2aNUuO1dTUuHjlypVJm3PPPdfFjz76qIvXrl2b9Pnud7/r4hYtWiRt9ttvPxfPmzcvaYP/b/z48ckxLSISFZu89tprN/uz9XOuu+66pE10DI3LX/7yFxf/4he/SNp06dLFxVqMKU9xvzwmTpzo4s9//vMVOS8qL3q2RMfUqFGjXPzTn/7Uxb/85S+TPhs3bnRxq1atXNy2bdukT+vWrV0cFZM+9dRTa7vU0Ntvv+3i66+/vuRzoH4NGjSovi8BjYDOTXvuuWf9XIiZjR492sVR4bivfe1rLo7W6aiMrl27uviRRx5J2myzzTYu1netW265JenzwgsvuPjhhx/OvJY8xQeVFgTXZyyKZ9ttt3VxVPhWRQUPV61aVbFrQmlatmyZ2Wbfffd18YIFC0r+nObNmyfHOnfu7OIbb7zRxV/4whdK/hwzs2uuucbFr776alnnAczMnn/+eRc/9dRTLp46dWpdXk7F8BfeAAAAAAAAAIBCYMMbAAAAAAAAAFAIbHgDAAAAAAAAAAqhXnJ4f+tb33LxBRdckNlH83P/9re/TdoceuihLtb8bjvuuGPeS/x/opyqt956q4t///vfJ20mT55c8mepDz/8MDn29NNPb/Z5mxLNo21mNmLECBdr/qvI8OHDXRzl/T7hhBNcfPfdd+e5RBTMJ598UpXzau2AN998M2lzyimnVOWzUXn33XdfcuzLX/6yi9u0aZN5Hs2tfdBBByVt5s+f7+Kdd97ZxVobwix99un4+6xj/ynKX3nwwQe7eObMmbWeAw3PyJEjXbzFFv5vJz799NOkj9Y6WLJkSeUvDPVKa+KcddZZ9XQl2Y466qjkmK7z/vWvf9XV5RSa5k83S9/PNF93pGPHji6++OKLkzYfffSRi1evXp153g4dOmS2UZrndM2aNSWfA43LHXfc4WKthRK54oorqnU5KIPu33znO99J2vzkJz9x8fLlyzPPO2TIEBdH+bgPPPDAWs+Rp2Zc9EzSOjjR+gv1Z9KkSS7WWiF16eOPP3bxr371q6TNlVde6eKi1KfgL7wBAAAAAAAAAIXAhjcAAAAAAAAAoBDY8AYAAAAAAAAAFAIb3gAAAAAAAACAQqiXopU33HCDi/Mk2M+TuL+cwiNq7dq1Ln744YeTNmeeeaaLNQk8GjYtZHn88cfX05WgiLTwoJnZj370IxePGTMm8zxz5sxx8eWXX+7i22+/vYyrQ0Nx//33J8cuuugiF994440lnzcqQBkdq4Y///nPLv71r3+dtImKraJx0fWYruGiYqY333xzVa8Jdevwww9Pjunvfcsttyz5vC+++KKLn3vuuaSNFvwdOnRo0uZ3v/udi6OCYKgbUXG/r371qy6OCsw/9NBDLt51110zP6t58+Yu3n777fNcYsleeumlqpwXDYfOX3369MnsM2PGDBc/+OCDFb0mbB69PxGdZ5588slqXY6zYsWK5NhVV13lYn2umZl98sknVbsmbL4//vGPLn7//feTNjrmTj31VBd37do183O0OGa0T6rjaf369ZnnLQr+whsAAAAAAAAAUAhseAMAAAAAAAAACoENbwAAAAAAAABAIdRLDu9yRHlxSxXllXz33XddfMQRR7h48uTJm/25AJqOqVOnJseOO+64ergSNDZ33nmni6OcpWeddZaLW7duXdVr+rcXXnghOXb11Ve7+IknnnDxxo0bq3pNAOrHD37wg+SY5k9WGzZsSI799re/dfGVV17p4jw5JjV3pZnZbrvt5uLvfOc7Lta1v5nZ0qVLMz8Lpbv33nvL6jdw4EAXn3jiiS6+5JJLkj69evUq67OynHvuuS6+5ZZbqvI5aDj2339/F/ft2zezz8SJE11MfuWG5e6773Zx//79kza639SjR4/M8+r+0sKFC5M2Og9qfvcJEyYkfTZt2pT52WjYtMbNX//618w+l112WbUup8niL7wBAAAAAAAAAIXAhjcAAAAAAAAAoBDY8AYAAAAAAAAAFEK95PDWvGs//vGPXfz0008nfQ4//PCSP+edd95x8eWXX560+cMf/lDyeQEAqLR169a5+KKLLspsM3r0aBcPHjy45M+NcuDq8/LJJ59M2nzwwQclfxYaP61t0rVrVxcvW7Ys6RMdQ+M1a9as5NgBBxzgYs1defTRRyd9NO9/pXz/+9938aOPPuriqCbB2rVrq3ItKI8+X/7yl7+4OMqF2rNnTxfr++VJJ52U9NGxfMghhyRtdP6KakKhaRk3blxy7MYbb6z7C0Fuq1evdvF5552XtDn//PNdfOSRR2aeV3NtP/bYY2VcHYBq4S+8AQAAAAAAAACFwIY3AAAAAAAAAKAQ2PAGAAAAAAAAABQCG94AAAAAAAAAgEJoVkrhjWbNmlGlo+l4paamZki1Ts5YalKqOpbMGE9NDHMTKoW5CZXE3IRKYW5CJTE3oVKYm1BJjCdUTE1NTbPoOH/hDQAAAAAAAAAoBDa8AQAAAAAAAACFwIY3AAAAAAAAAKAQ2PAGAAAAAAAAABQCG94AAAAAAAAAgEJgwxsAAAAAAAAAUAhseAMAAAAAAMV6kCsAACAASURBVAAACoENbwAAAAAAAABAIbDhDQAAAAAAAAAoBDa8AQAAAAAAAACFwIY3AAAAAAAAAKAQ2PAGAAAAAAAAABTCViW2X2lmC6txIWhwelT5/IylpqPaY8mM8dSUMDehUpibUEnMTagU5iZUEnMTKoW5CZXEeEKlfOZYalZTU1OXFwIAAAAAAAAAQFWQ0gQAAAAAAAAAUAhseAMAAAAAAAAACoENbwAAAAAAAABAIbDhDQAAAAAAAAAoBDa8AQAAAAAAAACFwIY3AAAAAAAAAKAQ2PAGAAAAAAAAABQCG94AAAAAAAAAgEJgwxsAAAAAAAAAUAhseAMAAAAAAAAACoENbwAAAAAAAABAIbDhDQAAAAAAAAAoBDa8AQAAAAAAAACFwIY3AAAAAAAAAKAQ2PAGAAAAAAAAABQCG94AAAAAAAAAgEJgwxsAAAAAAAAAUAhseAMAAAAAAAAACoENbwAAAAAAAABAIbDhDQAAAAAAAAAoBDa8AQAAAAAAAACFwIY3AAAAAAAAAKAQ2PAGAAAAAAAAABTCVqU0btasWU21LgQNzsqamppO1To5Y6lJqepYMmM8NTHMTagU5iZUEnMTKoW5CZXE3IRKYW5CJTGeUDE1NTXNouP8hTc+y8L6vgAUBmMJlcR4QqUwllBJjCdUCmMJlcR4QqUwllBJjCdUHRveAAAAAAAAAIBCYMMbAAAAAAAAAFAIJeXwBgAAAAAAAAA0bM2aNas1jnz66afVupw6xV94AwAAAAAAAAAKgQ1vAAAAAAAAAEAhsOENAAAAAAAAACgENrwBAAAAAAAAAIVQ6KKVeZKxV+tzampqao2BvLbYwv9/qVatWrm4devWSZ+ttvI/7ebNmydtNmzY4OJ33nmn3EtEA5VnDtS5qZw+UT/mPAD/Vs56jDmkWMp9tgAAAOD/yipAqftAZmadO3d2cY8ePTLbfPDBBy5evHhx0ufdd9918dq1a5M2H374oYs3bdqUtKkm/sIbAAAAAAAAAFAIbHgDAAAAAAAAAAqBDW8AAAAAAAAAQCE0mhzeeXL/ZeWzqVROb82pHOXJ2XLLLV2seQk/+eSTpM/HH3/s4k8//TRpEx1DcejYMjPbdtttXdy/f38Xn3jiiUmfPn36uHj69OlJm+eee87F48ePd7HmZDIjv2allDM35cnPpfnd27Rpk7Rp2bKli3VOieYYzb21bt26pM3GjRszz6MYT3VDx47OM5WqQ1HO/WQMNE156gCg4cgzZ2Tdv3LmmUqNkzznZS4CsDkqUacgz9xaqbkUQDHoGq1t27ZJmz322MPFRxxxROZ5ly1b5uL27dsnbebMmePiWbNmJW00F3hd13ThL7wBAAAAAAAAAIXAhjcAAAAAAAAAoBDY8AYAAAAAAAAAFAIb3gAAAAAAAACAQmiQRSvzFGzQopBmacL2PEV2tPCbFnlr3bp10kcTtm+//fZJm3bt2rl406ZNLl60aFHSZ+7cuS5es2ZN0uajjz5ycVT8Eo2XFhQ0MxsyZIiLzzzzTBdrEUszs/fff9/FUQHB+fPn19onKqDZ1MZbtYqo5ZnP9JgWpIzmHS1WGrXRIqj6OdG16Dy5YcOGpM369etdPHnyZBfreDMzW7VqlYspqlM6HUvNmzdP2nzuc59zcdeuXV3cvXv3pE+HDh1cnKfQst7PhQsXJn2WL1/u4mgs6fOSYs2puirwuOOOOybHdF7p3bu3iwcPHpz00fETzQdPPPGEi9966y0X6/oHlZFV1NYsfQZsvfXWSRs9pnE0Zlu0aOFiXTtHxZnzFFrW8abPp2h9rWsgnYfMeEZVS573Ph2X0diIxuV/iu6f3ufoWZc15vIUEcx7Pdg8ee5FniKQ1SjCa5Y9jrt06ZL00XV7NDfp9ehaK5rzovOg+Cox/vMUfmZ+qx69P/rs22mnnZI+Bx54oIt32WWXpI3OCRs3bnTxggULkj7vvfeei/XdMI9qjxX+whsAAAAAAAAAUAhseAMAAAAAAAAACoENbwAAAAAAAABAITSIHN558gdqbhrN+2eW5i7VXLRRPm7NXbrPPvu4eN999036aH6tKGec5pp8++23XRzlr5wwYUKtcXSeDz74wMXkOm1cdIzuuuuuSZvTTz/dxfvtt5+Lo3s+a9YsFz/33HNJG825pJpi7q1ycpaVc948Obw177Lmah82bFjSR9tEc1NWzmfNFW6WjrEoP9fSpUtdvG7dOhfPmzcv6aP/7ih3ZlMch7XR56M+CzX3tpnZgAEDXHzUUUe5eNCgQUmfTp06uTi653qPlyxZ4uIpU6YkfV555RUXz5w5M2mzbNkyF3/44YcuborPuTw5F1We346OJ82hHP0mdS3Vo0cPF0frpm222cbF0TNp6tSpLtZxgOrIs1bebrvtXBzNMx07dnSxPsP0HGbp3KT54HU8mqVzUZSbVtfYr732mot1jWSW1hx45513kja6tm+Kc1GWcnLCa45is3SM7bzzzi7WsWKW1hfo3Llz5udoztJorTJ79mwXa80lzZdslj4fNY989NnkVN580bMxq76A1uYyS9fCeWoH6Du5xlE/rZPRt2/fpI/W59E1kZnZihUrXKzPbq1REF0L81ndyao7F8mzpitn/o1q/+g6QMdGNJ/pszkap7zTVYbeZ3226XuemdnBBx/s4mh9r/W33nzzTRfrs8/MbOXKlZnnzVJuHYy8+AtvAAAAAAAAAEAhsOENAAAAAAAAACgENrwBAAAAAAAAAIVQLzm8s/LXRnlnNZdQ27ZtkzaaL7BXr14uPuCAA5I+Bx54YK19olzhmscoyou1evVqF2flXTVLc9WsWrUqaaOfRT7Bxm2XXXZx8aWXXpq0GTVqlIt1LEW5UP/0pz+5OMqT++6777pYc28VPc9WpXLglvNZOudFebN32mknFw8ZMsTFUd5lzVe5cePGzGvTuTXKuaZ5CKO5SXPC63iKcsTp91L0MVeqKPeefo86Tj7/+c8nff7rv/7LxTp2WrZsmfTJk9txhx12cHHXrl1drPUuzNLrj3KDa/63t956y8XROCn62MnK4R39+/PkhNT1lv7Wo/PqukNzM2seXbM0R6qOW7P0evPkLS/6fa8Gff5ofvXo3mh9kSjPrP7+99prLxdHOZe1VoCOxygXpD7XNH+kWZrPUuevaM7LyllqltbR0fHXFMdjOe90+r6maxczs379+rl46NChLt5jjz2SPprnW98V27Rpk/TR97EoN+2cOXNcrHlOX3zxxaSP5o1fvHhx0mbDhg0ubuo5lctZE+pzI1pr6hyntQGitYqO0TzrYL1fUU52/TfptfTs2TPps+eee2Z+to5bnSej9UCe9yCULvpes/aBNK921CfPZ+l8G815uh7T57BZuo7TWhnRfKbPR43N4rz2qF00nnT8HHTQQS4+5phjkj56n8ePH5+00WNadymqQajvBOW8o1V77cRfeAMAAAAAAAAACoENbwAAAAAAAABAIbDhDQAAAAAAAAAoBDa8AQAAAAAAAACFUC9FKzUJf54CJ1p4JCqYM2zYMBcPHjzYxVpAxyxNyq/FLqKCOVpkRIu1maUFJbTYZFSoS4u27Ljjjkmb2bNnu3jdunUujgpkoGHQIkpmZhdddJGLDznkkKSNJvJ/+umnXXzttdcmfbRYjo5Zs3SsNLViS9X69+YpVqLFSaIivAMHDnSxFq3s3r170keL1ETFQXS+0oJfWnzCLC1kGRV1mjt3rot1roqK+zJfeTp2ooI1WihS54zTTjst6dO/f38X6/hbu3Zt0keLzUSF4fQ3pMVXo0Jk+pzWZ5hZOo517ETP3OhZ3ViVU0iq3IJH2k/vaZ5CV3rf8xQN69ChQ9JGj+nasFpFK5taMUy95x07dnRxtL4eMWKEi/V5ZJYWttRxoWPALF3va9G3qIiyzhnR+kbX2HkKOut5Fi1alLTRtb3OO0UeN2blrW/yFAiN1kBagLJXr14ujoqr6pjTZ5sWXjPLvv7ovFqot0ePHkmf5cuXu1gLxZul45sigl6ewsX6nIjGmz5b9P06Gkta7E/XQNEcosei+UvHmxYQjOZfHV/RekeLCEbzIqpDx2U0h+jzT59J0XpZC6dGz1AdyzpWooKUeiyaf/XfsHTpUhfPmjUr6fP888+7OCqKSNHK0kXznha3/da3vuVifV6apWuaJ554Immjx1asWOHiPIWU86yD6nqtxF94AwAAAAAAAAAKgQ1vAAAAAAAAAEAhsOENAAAAAAAAACiEesnhrTQ3jebNMktzEg0dOjRps8cee7i4d+/eLtaccWZpfi3NDfrWW28lffLkydIcXfpvjPKUax/N3RQh31vDpfmvzjrrrKTNV7/6VRdH+VI1B9aVV17pYs3XbZbmr4xyJRU912R9yZPjUnMOdu3aNelz0EEHuXjQoEEu1hzLZuncFOU6XrhwoYs1p1eUz1nzekc52DQ/peY6jeoWMAY9HTtR3v/hw4e7+KSTTnKx5us2S+cizUX5+uuvJ30mTZrk4ihve1Yu4CiHnOar1JzeZmbvvPOOizV3YJQPvjHn8K7Ec7zcc2g+vqy1i1n6u23RooWLoxyqep6oTZ7cgOUo57vJ06cxzF/RvyMrB7vm4jYz23PPPV3crVu3zPPqcyPKZ6vPBX1mRc85/f1H6yadD3S8RflS9bOiMdoY7nk15RlPej+i9x39bqP3M73P8+fPd3E0NvS+67okeo5pfufouaW/Ec3VHH0vWssgyusbHcNni75n/Q6jPMft2rVzcZcuXVwc1czS55q+60frYJ3zonWJrv91Lu3Tp0/SR//d0dhfsmSJi6P5NktTq2dRLv2e9J5G+1j77befi4855hgXa80Cs3QfKKp9omtqnXei+6fjMlp76Tuc1j+IcoPrtTB2KiNai5x99tku3meffVwcPeseeughFz/66KNJG609oWMjzz1tiPMIf+ENAAAAAAAAACgENrwBAAAAAAAAAIXAhjcAAAAAAAAAoBDY8AYAAAAAAAAAFEK9VMvQZOZacCIq1KVFRKJk+VpUZN26dS6OCjho4YeZM2e6eMqUKZl9ouI3WjhMrz8q4qLXFyV4p8BJw6Xj+pBDDnHxZZddlvTRcaBFeczMLr/8chfrGI3GtV5LnmIB9V1QoLHKKthklhY00cIjRx55ZNJHixNqIZKouKQW2Z09e3bSZs6cOS5esWKFi6PikvpvjArxaEEfHZeNuahgNeQpcrXLLrskbY499lgX77XXXrWewywtVPrss8+6+Mknn0z66PMzmh90rGvRneha9Dm93XbbJW369u3r4hkzZrhYi26alVegqTHJKqIYFR3KU3ixnGeF3nedi6Lfun5ONH/pfS1nzqCQtxfdz6xnVlRAUAu4RfdPj2nRwajQmj6jNI4+R4vSDRgwIGmjhe71nSEaJ1poLM9319SUM56iNZHSQsVm6Vpk3rx5Lo7mPC0mqWNQx7FZWiQwKjinRQ61OGF0Lbom2rRpU9JG5zjW4V6e9bW+R0VFK9u2betiLVK5ww47JH30Hucp3q7HovW07nH07t3bxVqo0Cwd11Gh8aznZ57C8Yy/VJ53Or1no0ePTvp8+9vfdrHOO3kKskZrah0bOg60CKGZ2eLFi12sc2t0TM8b7VfoO6i+RyAfHV8jRoxI2ui7oPZ56aWXkj533323i/V+mVWuSGVDw194AwAAAAAAAAAKgQ1vAAAAAAAAAEAhsOENAAAAAAAAACiEBpEQWvPORHmqNS/WqlWrkjaad0ZzHS9btizpo200J5Hmf4uuN8ozqXlVNU9plO9N/43vvvtu0kZzwJFvq+HQez527FgXR3kxNX9SlOf7lVdecXFd5aplbOWTVZPAzKxdu3YuPvDAA10c5XvbfvvtXaxzXsuWLZM+eizK96bXq3kJo/Gl81WUCzAr32E05zVlUV5AzTP5la98JWmjudz0nmtOdjOzJ554wsX33HOPi1evXp300byAER0H+ryM6nHkyZ2pucs7duzo4mhcF13WfJwnb7auXfKcN+qj6xAdK9HY1j7vvfde0iZab2XJkztQ/43l5Bss0vNQ52L93qO8jlr7Ifpta75OzcetaxkzswULFrhYx5Lm/Dcz69mzp4t33333pI3OM+3bt3fxmjVrkj76PURtyslv2dToc2H9+vVJG10fRLm1dc2j80qUQ1k/O09djG7durlYx5dZmqNXn5lLlixJ+uhYjtZWeWo34f+L1teac13X22ZpHS3Nmx29n7355psu1mdW9LzKcz91ba/54aPfwqxZs1w8ffr0pI3+pvLkh2e8pXR9ENVc0zzxBxxwgIvPO++8pI/WpslT4+K1115z8dKlS5M2r776qot1HyvKo63jNJqjdX7V8RTNZzr/8t6Xj445nSOivSGtS6LPoFtuuSXpo3nZo9+/rvnLmSMaYt04/sIbAAAAAAAAAFAIbHgDAAAAAAAAAAqBDW8AAAAAAAAAQCHUSw5vzduiuR2j3NXTpk1zseb9M0vzCWm+rei8mmOtnHy8efKFaQ7SKGfqokWLXBzlUdQ8S+RHqhuaO1BzJ5mZXXvttS7u0aOHi6P7qXmZnn322aSNjmsdo1E+Ur3eKL+rYixli77rPHPGgAEDXHz00Ue7uGvXrkkfzc2s88yuu+6a9NG8qjoPmaX56HSsRPURNAdclBMuKx9XU88VmGec7Lnnni6OcrtrPQh9hj333HNJn9tvv93Fmitww4YNSR89bzQ/6NjRNjvttFPSR2sdRPkRNX+djv2oT9GVk4e6Er85XZ9Fx3TeifIu67VE+U+jcZhFv4fo31xOzu48522s9P7pevT5559P+sydO9fF0W9Q6+RovHbt2qSPrk00f22UD17zKWsOZjOz7t27u1jnpugZpjl781xvUxP9DnQ86W8nqvmh9zWqFaBjLM/aV8eP5moePHhw0uewww6rtY9Z+n6mv4co/+7KlStdHM1vTTmHd555WdtEtWt0faC/fbN0baV52qP9Ab1/mj85T50TrR1gltYc0Dom0fei9Q+iWi36OySHd3l0LorWM/vuu6+LL7zwQhfvscceSR8dL48//riLH3jggaTP66+/7uKo9onWmsiajyPROMiqV5FnPLGvkIruh+aEP/HEE1289957J300x7rWapo8eXLmZ0fzk45/7RPNe1k1vMzS+aiu11L8hTcAAAAAAAAAoBDY8AYAAAAAAAAAFAIb3gAAAAAAAACAQmDDGwAAAAAAAABQCA2iaKUW7YgKOmqi/jzJzvMk2Nfk7FokRYvjmKWFuIYMGZLZpnnz5i5+5513kj5aMEcLZpil35UWBIiS4VOUonRaUE6LCmqxSTOzUaNGuViLoFx11VVJnyeffNLF0bjWa9FiP3nub56CFMguVBb9vnQO0aKCZmnRSi3yFhV40N+6FqiICkdoAQothGGWFnXS4pdTpkxJ+mjxC70Ws7QgU56COUWWVWArKnY8ZswYF0fFl/Q8c+bMcfGtt96a9NGCW1qwTceaWXbBmogWy1m8eHHSZvny5S7u169f0kbHqP5e8hQmbEzjrRIFKfMU1M1TZEjjqOCc3ufWrVu7OCospvNDVMBNz6txpe5puc/MotDfthZnjAqK6jwTrR/yrLmVro31mRUVpNTnad++fZM2Oga1qKAWAzNL/43Rc67I46Jc5dz3cuY8nWei97NBgwa5eOTIkS4+4IADkj66vo/eQadOnepiLVoZFaXX80S/K8ZT7XR9Ha1pO3fu7GItjG1mtttuu7lYn2vRWmX+/Pku1vsXFe7VtUufPn2SNrpn0KVLFxcvWLAg6aPr8mj9n2c/AF5UFDnrHcnM7NRTT3Wxjrmo2PekSZNc/K9//cvFWpjULF0nRe+KOpbL2R8rpw1zV3l0zWOWFrLVopVRn3nz5rlYx1O0DtfP6dixY9JG51wdc1GBZi2iGxVX1TmrnHfMzcFfeAMAAAAAAAAACoENbwAAAAAAAABAIbDhDQAAAAAAAAAohAaRw1vzuES5j1SUm0pzMWkc5cDRfGCaK3T77bdP+uyxxx4ujvI76XlWrVrl4tmzZyd9NM9qlO8tKzdTObk/m7poXPTq1cvFmq9L83WbpXlwn3jiCRe/9NJLSZ8OHTq4WHMnmaVjVMdAlOtd8ydFuSij8QUvKw+zWXrPNP+bWXrP9PevubiiPpr7OLp/mo+rR48eSRs9duCBB7q4bdu2SR/NKRjlqyQnfO107Oyyyy5Jm7333tvF0dyk+XYfeeQRF0fPFp0j8uRGznpOm6X/Js0lqLnCozZRHkz9LP0eGGupcnMyZt3naM2hbbRuQXRPdb5atmxZ0kZzBea5z405d3t90e9I73G5v6+s52U0n2ktg5133tnFUT5ezUMZ5fXVfJEzZsxwseZTNYvnK5W1xmYNno9+T9GcoesZHRv77bdf0uewww5zsb6vRTVWdH0c5R/VYzqf5XnWkVPZy1OLQ9fT+s5klubJ1nzdZum6dunSpS6O8tLqvKjjMbp+vb5o/tIc3ppr99VXX0366Hhr1apV5vXqHkil6i4USfTOrbUBovE0bNgwF+szKMrDrrXbdDxFe0maGzl638/aM4vqsOixPGtDlEfrsGmNP7P0uaVtomeS5nzXNvq8NEvfO7V+hZlZz549Xay/h1mzZiV9nnrqKRfPnDkzabNkyRIX67jNk3t+c/AX3gAAAAAAAACAQmDDGwAAAAAAAABQCGx4AwAAAAAAAAAKoV5yeJdDc+BE+d40D5bmVIpyf3Xu3NnFrVu3dvGOO+6Y9Bk8eLCLu3XrlrTR3EeauybK2btmzRoXaz7LiH4PUS4nzdvV1PMy6VjSMWBmNmbMGBdrzrUoD+DixYtdPHnyZBdHuZI0h5fms4w+691333Wxjq3oWqJ8YnnyVTZ1eXJ467ENGzYkbaZNm+ZizW+VJ2ekzgfRHKg1BwYOHJi0+dKXvuRiHZeavz46r/57zJhXlI4LzbkY5QXU5000n8+dO9fFTz75pIujHH9Z+foqlb9P/81Rzkj9HnQ+Nktzqmr+yignYZHkyW1ard9bnvVCp06dXKxrreie6vMmTy5Kza+ZJ8ef5i01S/9N5P2uXZ7xFz0L9b7nyb+reUs1x2S/fv2SPppjMrqWhQsXunjixIku1tyoZum/O8rvGo2v2s4RHWO8pd+jvr+ZpWNBa+fou5hZun7RdZLWwDBL12zR/enevbuL33zzTRfrmtss/TdG82KeWhlNiX5nmj82eifXOSKqXaPjS2suRe9eAwYMcHGe+mL6vjZo0KCkjY4lzSeuewHRefPsD2gdg+j6q503t6GL5nPdB9KabGbpGlTzc0dr1Pbt27u4b9++Lo7mB52vtNabWboHsHr1ahdH7/p63/OMbWSL1iI6h+l9NzMbPny4i/V+aA0vM7P58+e7eIcddnBxNKfp3BPNlZo/XK9fP8cs/c1EY07XXFFtuWriL7wBAAAAAAAAAIXAhjcAAAAAAAAAoBDY8AYAAAAAAAAAFAIb3gAAAAAAAACAQqiXopVZxW+ixP2afF2Tqpulydh79+7t4ijRuhYj0Lht27ZJH/3sqKCBFjCZMWOGi996662kT1TsTjVv3tzFmtg+SpivxTmakuj70AJbUVGRPffc08Wa2F+LgZiZzZkzp9Y+UaGCNm3auDgaS1kFUKOiG/ob0uIGyEe/xzyFuqLiH1qkVn+TUfEGLXqSVXjNLC1WEtl7771d3L9/fxdrUUGztGBO9NlNvdiS0u9If9tR0Uqd37VYo5nZCy+84GIdW9FYyro35RZR03+jFnnTudYsLTqlxU7MzN5++20X628h+o0VvRBc1r8v+k1WQlQcV4sv6X2Oxtv69etdHD1D9VmW557m+XdTpHLz6W87+t3qc2LnnXd2cbQG0qLJWpCyS5cuSR8tQBcVedMilbNmzXJxVPRNx1+0ttI5Wot/RUXgmnphuOg3qt+troXN0oLyutaK1juvvfZarX2i9yG9p1EBTR2X+t6gxeTN0sK80WdXa95uDKJ/u84z+u6v7/lm6RwRrTv0Hut7fNRH54g8Rbk7duzoYl1fR5+l81d03jyFFHVdpOumpjzW/i3Pd6DfvxakNDN77rnnXKz3I1rfZBURjQqyasHA6P1M58433njDxdF7BGugysh6BzJLn2N9+vRJ2mhBby36uGzZsqSP3kMt2BzNETpWouKqCxYscLFef7du3ZI++j4btZk6dWpyrC7xF94AAAAAAAAAgEJgwxsAAAAAAAAAUAhseAMAAAAAAAAACqHqObzz5G7THGtRfm7NqTxs2LCkTb9+/Vysub+inIMtWrSotU+UL0lzgUX5AzUHaZRXVWneuCiXueaA0vNqnkyzNE9UlLupqPmcovyjnTp1cvHuu++etNE8cXovVq1alfTRvF86lnSsRcc0b5NZ+hvS8RflPtTcc1G+yqaeczmam3S86FiJ5gP9HqNc/GvXrnWx5nKM8mjpb1LjKL+o3ufoWjSfqM4ZUX7k9957L7ONfg9FnVPy0vuj+dSiPG46DqLce5o7UHOHRvliK3Evot+LHtPfx+DBg5M+WtsgujatUzB79mwXR/MZUtEcUaooJ6E+2zTPZLTe0TyYUV7MrGdStXKQRudtyvNXtPbUvLOan9vMbOjQoS7WvI6a3zY6prV28jxzX3/99aSN1snJs27S51qUy1nn1zxrLV23R+O8MY+3rN9lnt9t9L2tWLHCxdOnT3ex5uvOI8rVrM+kKO+y5jHdZZddXKzPKLM0l260HsvKs9yYx0U59Hmj6ya9D2bpu32U213nNP0d6zMsuhadH6I9hV133TXzvJrbXeeqKNe7Xn/0bNfr07ho80458tSiUIsWLUqO3X///S7Wexa9y+szSPM2R89Uzeut5zBL3y91fbZ06dKkT1O779Wi83W056T3MMrVrr9nncOi56POhfp8iX7vy5cvgha+tAAAE5xJREFUd3G0f6m/Ea37FeXn1uuP3kP1WVfXY5C/8AYAAAAAAAAAFAIb3gAAAAAAAACAQmDDGwAAAAAAAABQCGx4AwAAAAAAAAAKoV6KVmpCdE28fvDBByd9Ro4c6eKoqEhUaOc/RcnktWCOJvuPrl+LEWjBALO00I4W4okSxWtxjqgwnCal16JmWuzLLE12H523KPR+RWNCC0NokRGztJikJtePigxqgQAdJ1HhCy1u2rVr16SNjlEtjhEVwtEiFVEx06Yu+m1rMYZevXq5OCqYo2Nh4cKFSZu5c+e6WItY5SneoNcbFZPTAjlaNMwsLciq80FUTE7nFYoGZtP7owVqogJu+r1q0S6ztGBuVBxEZRXCylNULCqSpEXddC4dNWpU0qdt27Yujp6Fc+bMcbHOZ0V+hlWS3tfoPutzS5+ZUeG+9u3bu1jHdjRWdJxGhbn0WvQ8eQotV6oQTlMqHqffc1TYT4v0HX300UkbLQCuzyP97ZuZbb/99rV+dvRb12ds9DzS82qxr+hadC6K3hn0s/RZrms6s3T+iq63sYyvPM+KPPOOitaxWghy3rx5Lo7mA52/tE303NW5KU9xVS12p3Ng1CZas2XNcY1lXJQj+rfp96EFkJctW5b0mT9/vovLWTdFv/WscbvXXnslx/r27eviaIzq+8DMmTNrjc3SwpZRQXOdv3Ruaorrdv3N6XMguu/qvffeS47pnJ7nnU6fbdpH383M0udYtKeh413vc579ClRGNJ50b2+77bZL2uhYyLOXpWNZ1x5aoNIsfb+Pirb27NnTxdG4VPpZ0WfruKRoJQAAAAAAAAAAZWDDGwAAAAAAAABQCGx4AwAAAAAAAAAKoeo5vKOcZZrPZsSIES4+5phjkj477bSTi6M8OVGeov8U5aLU/DWaU0bz6pileXKiNpqrSfNGa45vszS3TpSjS/Neap6c1atXZ15vkWn+tyiXqOaEjNoovcdRjji9f5pvMBqfmmNZ83WZpfnD3n33XRe/+eabSR/NcRmN0aZGc2JF971Pnz4uHjRokIuj3Fv6m4vy/ulvWedFzVNolp1jLboWrXUQzaVa/0Dzu0+ZMiXpozkSo7Fc5FyT5dBxoPUhojy5mlcyyoOmbXQsReNP702eeVLPEz0/df765je/6eIox6V+lubrNjMbN26cizU3ZZ685Ujve57fqH63Ws/CLM0dqLmao7zLWuskqiuhYyOrLkvUJ5qb9N+kc2s01+bJO9xY6b9N76euV83MxowZ4+Jhw4YlbTQvtq5zdZxEx/R+Rrne9fp79OiRtNG1lM4hUe0AfQ5H84xej+bWjWibxvysjH4XWb+VPM+X6HereZfLyQ2unx310bkpymuqx3RslJPbPG+/pkTzvGrO7gkTJiR93njjDRdH+wP6DqTjLfqt67yoedq1roFZuh6LauLcc889Lv773//u4mhOKWe8aZ/GPO+US++H/o6jvRhdd0TPIK05oOM2z75CmzZtXBzVE9Oc8NEej16LzpvR+yUqI8++gs490XpTf5s6LqPz6hpan2PRHqK+h0Zrp8GDB7u4S5cuLl68eHHSZ+rUqS6eNWtW0kbnXHJ4AwAAAAAAAABQBja8AQAAAAAAAACFwIY3AAAAAAAAAKAQqp7DO8ql1b59exfvueeeLu7WrVvSR/NIRjkitY3mR4rysmmOK82PFOXfWrp0qYuj3HN6fVnXZpbm44lydOn1al7V6N+YJ5dUUWi+riiHvI5JzXtklua80nzJmtvNLL3Hmnsr+i1o/qQoZ+miRYtcrDlvJ0+enPTRHONNMXebypNrS++z5vTW+gNmZh07dsz8bP3dag6sKJ+tzhE6Lw4dOjTp88UvftHFWvvALM1Hp/kPH3zwwaSPzotZ+cWR0jkjT87I6Nmiz89OnTq5WHP8m6X3S38L0dykud6i8XbEEUe4eP/993dxNJ9pzYFHH300aaNjUnMUMp/lk+d70rGgz8wox6XmR9Y+ms/SLL2H0djQcajXHz3PdR6P/s26HmtK4ydaR+p3pvUEBgwYkPTp3bu3i6O1pp5H55Coj16fPp+isaRrrWgsad5M7RPVJNB6KdF3p2303xzV1mjMeZrzrJuyakLk6ZMnv7XG0e8469mmz0uzNBdzr169kjY6D2p9jSVLliR9oty/pcpTk6NIdK7W35s+R8zSd/LoOREd+0/Rmlbvec+ePV0c1SjRe/Pyyy8nbf75z3+6WHN2Z9UjQ3563/WZpO/pZumcoXnko/PquIzGm46fI4880sX6/maWvnP+61//StronsDKlStdTM2b6tHfe/Rd6x5TVFdP87drvYBoz1DXHtonmiu1zoqu68zSuio6/sePH5/0eeyxx1ysY9Cs/vcNms5OKAAAAAAAAACg0NjwBgAAAAAAAAAUAhveAAAAAAAAAIBCYMMbAAAAAAAAAFAIVS9aqYVizNKk6Vo8Jiomo8W8oiJb0bH/tG7duuSYFqWcNWuWi7UYgFlaEDC6Xi00oAnno0I877//votXrVqVtNHk91rQo6nJKmITFXfR71WLz5ilY0WLTWiBgehaosKqSosKaPESM7MXX3zRxRMnTnSxFrU0o0hFHlGhq6wiSVHBIy0UEc152kaLYUa0MNyQIUNcrMV+zdJiEzqnmKVz3NixY12sRQXN0uJLRS6aVClaoEMLk65Zsybpo4VJu3TpkrTp3r17refRMRDR55HOb2Zmw4cPd/HAgQOTNlrcS//Nr732WtLnkUcecfG4ceOSNlrUhfksWzkFKqNjuo6KilZq4Sddz0Sfo4W4tOCOWTp36n2PzqtjLir4lbUuKHJhuOjfpgW1dM6I1jfaJ0+Rtzx9dJ2ka6+owJOug/OcVwuRR4WitU9U9EmLAuv1RudtzLKK2pqlc0aeosj6The9Ryn9TUa/db0+XbPpOsrMbMSIES6Oin3rmJs9e7aLo3WTjqfoOcazrXZ6j6P3qjwFUHVcaByNUZ0XDznkEBdrAXGz9J38mWeeSdroO2c5RZUbcyHc+qT3PSrQrGvsBQsWJG1Wr15d6+dEY2O//fZzsRZ51/dEM7PXX3/dxU8++WTSZv78+S7W976irGUaIv1uo7093WecPn160kbf/fS5teOOOyZ9Onfu7GJdt+Upmhxd77x581yshXcfeOCBzD7ReSlaCQAAAAAAAABABbDhDQAAAAAAAAAoBDa8AQAAAAAAAACFUPEc3ppXKsqTqzS3nubIMzPr2LGjizVnl1maS2/p0qUu1ty1Zmavvvqqi+fOnZt5LZoXR/Ohmpltt912Lta85VGeXz2v5go3M9u4caOLNb+Z5goyS3PpFCmfk/5bNB9elO9t2bJlLtb8RGZpviQdf3o/zdKxr/m4Fy9enPTR3GBRbifN4aU5IqP7WaR7XC1RPim9R/pdR/lsNRdllOdbc3br7ziaDzSHt+b4iuZWzR84YcKEpM3111/vYs0RH+UgJc9k6fQ707x6mgPULM3Bpjn+zNJxoDkio9zgmi9Z57coh/cOO+zgYh3nZunz8ZVXXnHxP//5z6TP+PHjXRzlR9TnGvNZecrJBao5LqP7rsfy5CDVMRjl7M3K+Ryt+1Q0r2d9D0UaX3o/o+9Mnx26nom+D12fRs8EfXZo3mNd65ul63adU6I1rc6l0T3XOUTbRM9PHccrV67M/GyNo5owmgu8MY03vdboe8vK+x/NIfrepHUBzOJ10X+KcpTqWn3w4MEuPuigg5I++vyLap/oe8Kzzz7r4iinr54nut76zmva0JTz28gzRrNyXkdr+/79+7t4t912c3GUQ15rKk2bNi1pE9UGKFX0PWXVIGpM806lZOWAj3IN6/q4a9euSRt9lmm+92g+05o8+mzW9bOZ2cMPP+zip556Kmmjexq8r9WdPHUl9P5EdO2rz5OohpLWmsg6h1m6xxHVKdQ9Aa3FFK2LsvYmGwL+whsAAAAAAAAAUAhseAMAAAAAAAAACoENbwAAAAAAAABAIbDhDQAAAAAAAAAohIoXrdQE7lHxD014rkUdtLiRWZqcPSpAocViZs6c6WJNvG6WJnDXom9RkR0tehAVyNh2221drEUso3+jFjSJisdpgRNNDB8V2YwKNxaVjr+oSMyqVatcPGXKlKSNjtH58+e7WAtURH20eMkbb7yR9NFiBpr43yxfQTBky1NcYuHChS7WAiFanMnMrEePHi7W37pZOl60IGBUMEfnCB0bURHev/3tby5+8MEHkzZaTFV/IxRRqgz9HufNm+fiqFhur169XBwVKtlll11crIUuo/lBi9hom+iea8G56Pn53HPP1RpHBVG0mHQ0RzPH1Z2sgl/Rs07baBGuqFCirmeiYlFZa5XovDovllO0skjy/Fv1e9T5ISq8qIUto3Wv9tMCulHxJl2z6vwQvUPo3BS9D2hRVC2SGD1zdb2vn2OWrht07EdFeCtRpK6+6HiKxpd+11o4Uou1maWF4KI2WvhNx2BUPF7Pq2utaA2nY27SpElJG3226fP87bffTvro2M1TyCvP943aRd+ZHtNxoGPWLC1aqe/10dw0ceJEF+t6+7OurxoYO+mzTYsg6/0yS9chusY2S9dFWjBQY7N0ravvm48//njSR99Bo4KB7BE0HNH6U9eoujdklo5LLWAaPR91DGqR52j/UsdPVBRci13quI3e2RoD/sIbAAAAAAAAAFAIbHgDAAAAAAAAAAqBDW8AAAAAAAAAQCFUPIe3ypND5qWXXnJxlOevffv2Lo7yW2uePM2ppnlpovNoXsk8+WyjvGx6THN9NW/ePLNP9Nl6THM3RTkwm1I+pzy5afPkI9KxovngozyAOnY0d1ie3KJN6V7VtTxjQ3OQjh8/3sUrVqxI+mh+t27duiVtNPef5ryMaD0BzQGvOb7M0nyB5ISvPzq+9Ln37LPPJn0071/0zNptt91c3KFDBxdHzxal5507d27SZvbs2S5+9dVXkzZTp051seam0xy+Zuk8yfirO9F3rfmPNb/zunXrkj6vv/66i7XeQDQGNU9utIbLemZGc3Y5a7amJPo+9BmgdU2mT5+e9NF81tF59b7reaMaMzov6rjIk/c4yu2uOcazYrP0fUXHo1n679b1Y/R7KdKYjN4x9HvSOMqXruukaN2kx/RZF+Xw1rW5jqeo9smECRNcPHny5KSN1qPQcRvlac8zN/H8q12e70fbRGM0K8dy9E6nz0bdU4jeB/TZGK2BqnXPGUsp/U70dxrNBzqnL1myJGmjOd91zEU5vHUs6JyitcLM0r0j5pDGJ0/9MD2m64ioNkhUu6S2z42OFWltkoW/8AYAAAAAAAAAFAIb3gAAAAAAAACAQmDDGwAAAAAAAABQCFXP4R3lkNHcNJpjTfOAmqU5+vLkpsmTH7kSuY/y5FTSnHZRvrdy8jmTA7p25d7zrFzpUW56NC7RONA8WpqjL8pBqvlONZ+tmVnr1q1rvZZoPGkeOc3HrfkhzdJ5hvmg4dA5JMrJ9sADD7h4ypQpSZvevXu7uF27di6OclEqHW+a+90szS+o+SvN0jGpY5bx1/Dps07XZ5rL3SzN5/z888+7OE8Ob43N0rVgVr7kqE2U17Apj8M8ec/1t635ic3SfP1Z+SOjz85Tl0bvVfQ5ec6r7wyaUzXqU87zsynlwDSL1x06Z+j3Fv1uda2l9VPM0vWW5uyO5hnNkztv3jwXv/baa0kfnYu0fopZmktX/015xhM2X7nvdFm5a6OaEtOmTXOxrpt0TJiZTZw40cVRTn/W6fVHv3tdy5ilObujsdG5c+dazxuNDZ2btMZF9DlZz0cUU9YeIkrHX3gDAAAAAAAAAAqBDW8AAAAAAAAAQCGw4Q0AAAAAAAAAKAQ2vAEAAAAAAAAAhdCslAT4zZo1I1t+mfIU+MmjDgsWvFJTUzOkWidnLDUpVR1LZg1rPEW/9azffzlFd5ow5iZUSpOam/LQuSrPfLbVVr7+eTnFCT/rWNa1NLCigU16bqrUOrdUeZ6Nea6tgT1jCzM3aQFRs3TO0NgsLTSqcVSgWYthapxnbooKhGXNXw1s7ESazNwU/dZ1DOp4a9GiRdKnVatWLm7ZsqWLtRC5WVoEccOGDUmbAhQiLMzclGesRPNMmzZtam0TzSFZhW+bcGHCwown1L+amppwscdfeAMAAAAAAAAACoENbwAAAAAAAABAIbDhDQAAAAAAAAAohDRpGqpCc3Q1wnyCAHIoJzctADREeeYzzZOb5zzl5F1mHm3YGvL9acjXVnRR3uyPPvqo1hgoV/Rb1/zIGkfjb926dSV/DhqXcsaKmdkHH3xQ63mjPR7GC1B/+AtvAAAAAAAAAEAhsOENAAAAAAAAACgENrwBAAAAAAAAAIXAhjcAAAAAAAAAoBAoWllPKF4AAAAakkoUkyx3fZNVzJt1EwCg0ni2oJIYT0DDwl94AwAAAAAAAAAKgQ1vAAAAAAAAAEAhsOENAAAAAAAAACiEUnN4rzSzhdW4EDQ4Pap8fsZS01HtsWTGeGpKmJtQKcxNZahWfsoC5L1kbkKlMDehkpibUCnMTagkxhMq5TPHUrMCvGAAAAAAAAAAAEBKEwAAAAAAAABAMfyfduxABgAAAGCQv/U9vsJIeAMAAAAAsCC8AQAAAABYEN4AAAAAACwIbwAAAAAAFoQ3AAAAAAALwhsAAAAAgAXhDQAAAADAgvAGAAAAAGAhizvy4JAAwhAAAAAASUVORK5CYII=",
            "text/plain": [
              "<Figure size 1872x360 with 20 Axes>"
            ]
          },
          "metadata": {
            "needs_background": "light"
          },
          "output_type": "display_data"
        }
      ],
      "source": [
        "# display original and reconstructed images, in rows\n",
        "display_images(images, reconstructions)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 70,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "collapsed": true,
        "id": "v5LGx1LPKmiS",
        "outputId": "d4bdb79c-9792-49a6-c014-c037b9a95493"
      },
      "outputs": [
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "shear transformation\n",
            "15\n",
            "Test Loss: 0.03909738\n",
            "\n",
            "Test Accuracy of     0: 99% (974/980)\n",
            "Test Accuracy of     1: 99% (1132/1135)\n",
            "Test Accuracy of     2: 99% (1024/1032)\n",
            "Test Accuracy of     3: 98% (990/1010)\n",
            "Test Accuracy of     4: 98% (971/982)\n",
            "Test Accuracy of     5: 99% (888/892)\n",
            "Test Accuracy of     6: 97% (934/958)\n",
            "Test Accuracy of     7: 98% (1011/1028)\n",
            "Test Accuracy of     8: 98% (964/974)\n",
            "Test Accuracy of     9: 97% (982/1009)\n",
            "\n",
            "Test Accuracy (Overall): 98% (9870/10000)\n",
            "30\n",
            "Test Loss: 0.05525802\n",
            "\n",
            "Test Accuracy of     0: 99% (971/980)\n",
            "Test Accuracy of     1: 99% (1128/1135)\n",
            "Test Accuracy of     2: 98% (1021/1032)\n",
            "Test Accuracy of     3: 97% (983/1010)\n",
            "Test Accuracy of     4: 98% (970/982)\n",
            "Test Accuracy of     5: 99% (887/892)\n",
            "Test Accuracy of     6: 96% (928/958)\n",
            "Test Accuracy of     7: 98% (1009/1028)\n",
            "Test Accuracy of     8: 97% (948/974)\n",
            "Test Accuracy of     9: 96% (969/1009)\n",
            "\n",
            "Test Accuracy (Overall): 98% (9814/10000)\n",
            "45\n",
            "Test Loss: 0.12081570\n",
            "\n",
            "Test Accuracy of     0: 94% (927/980)\n",
            "Test Accuracy of     1: 97% (1108/1135)\n",
            "Test Accuracy of     2: 95% (987/1032)\n",
            "Test Accuracy of     3: 94% (952/1010)\n",
            "Test Accuracy of     4: 96% (946/982)\n",
            "Test Accuracy of     5: 98% (875/892)\n",
            "Test Accuracy of     6: 90% (863/958)\n",
            "Test Accuracy of     7: 93% (964/1028)\n",
            "Test Accuracy of     8: 90% (878/974)\n",
            "Test Accuracy of     9: 91% (919/1009)\n",
            "\n",
            "Test Accuracy (Overall): 94% (9419/10000)\n",
            "60\n",
            "Test Loss: 0.25856522\n",
            "\n",
            "Test Accuracy of     0: 83% (822/980)\n",
            "Test Accuracy of     1: 81% (921/1135)\n",
            "Test Accuracy of     2: 90% (933/1032)\n",
            "Test Accuracy of     3: 82% (829/1010)\n",
            "Test Accuracy of     4: 88% (867/982)\n",
            "Test Accuracy of     5: 94% (844/892)\n",
            "Test Accuracy of     6: 73% (702/958)\n",
            "Test Accuracy of     7: 86% (888/1028)\n",
            "Test Accuracy of     8: 76% (748/974)\n",
            "Test Accuracy of     9: 72% (736/1009)\n",
            "\n",
            "Test Accuracy (Overall): 82% (8290/10000)\n",
            "75\n",
            "Test Loss: 0.38287478\n",
            "\n",
            "Test Accuracy of     0: 69% (679/980)\n",
            "Test Accuracy of     1: 66% (760/1135)\n",
            "Test Accuracy of     2: 78% (807/1032)\n",
            "Test Accuracy of     3: 66% (667/1010)\n",
            "Test Accuracy of     4: 85% (835/982)\n",
            "Test Accuracy of     5: 86% (773/892)\n",
            "Test Accuracy of     6: 62% (600/958)\n",
            "Test Accuracy of     7: 77% (800/1028)\n",
            "Test Accuracy of     8: 62% (608/974)\n",
            "Test Accuracy of     9: 61% (621/1009)\n",
            "\n",
            "Test Accuracy (Overall): 71% (7150/10000)\n"
          ]
        }
      ],
      "source": [
        "from math import atan, degrees\n",
        "\n",
        "# shear_degrees = degrees(atan(0.2))\n",
        "\n",
        "# for i in range(1, 6):\n",
        "#   print(\"degree transformation\", i * 20)\n",
        "#   # convert data to Tensor *and* perform random affine transformation\n",
        "#   transform = transforms.Compose(\n",
        "#       [transforms.RandomAffine(degrees=i*20),\n",
        "#       transforms.ToTensor()]\n",
        "#       )\n",
        "\n",
        "#   # test dataset\n",
        "#   transformed_test_data = datasets.MNIST(root='data', train=False,\n",
        "#                                         download=True, transform=transform)\n",
        "\n",
        "#   # prepare data loader\n",
        "#   transformed_test_loader = torch.utils.data.DataLoader(transformed_test_data, \n",
        "#                                                         batch_size=batch_size,\n",
        "#                                                         num_workers=num_workers)\n",
        "\n",
        "#   # call test function and get reconstructed images\n",
        "#   _, images, reconstructions = test(capsule_net, transformed_test_loader)\n",
        "\n",
        "# print(\"scale transformation\")\n",
        "\n",
        "# transform = transforms.Compose(\n",
        "#     [transforms.RandomAffine(degrees=0, scale=[0.5, 2]),\n",
        "#     transforms.ToTensor()]\n",
        "#     )\n",
        "\n",
        "# # test dataset\n",
        "# transformed_test_data = datasets.MNIST(root='data', train=False,\n",
        "#                                       download=True, transform=transform)\n",
        "\n",
        "# # prepare data loader\n",
        "# transformed_test_loader = torch.utils.data.DataLoader(transformed_test_data, \n",
        "#                                                       batch_size=batch_size,\n",
        "#                                                       num_workers=num_workers)\n",
        "\n",
        "# # call test function and get reconstructed images\n",
        "# _, images, reconstructions = test(capsule_net, transformed_test_loader)\n",
        "\n",
        "print(\"shear transformation\")\n",
        "\n",
        "for i in range(1, 6):\n",
        "  print(i * 15)\n",
        "  transform = transforms.Compose(\n",
        "      [transforms.RandomAffine(degrees=0, shear=i * 15),\n",
        "      transforms.ToTensor()]\n",
        "      )\n",
        "\n",
        "  # test dataset\n",
        "  transformed_test_data = datasets.MNIST(root='data', train=False,\n",
        "                                        download=True, transform=transform)\n",
        "\n",
        "  # prepare data loader\n",
        "  transformed_test_loader = torch.utils.data.DataLoader(transformed_test_data, \n",
        "                                                        batch_size=batch_size,\n",
        "                                                        num_workers=num_workers)\n",
        "\n",
        "  # call test function and get reconstructed images\n",
        "  _, images, reconstructions = test(capsule_net, transformed_test_loader)\n",
        "\n",
        "\n",
        "# print(\"translate transformation\")\n",
        "\n",
        "# transform = transforms.Compose(\n",
        "#     [transforms.RandomAffine(degrees=0, translate=[0.2, 0.2]),\n",
        "#     transforms.ToTensor()]\n",
        "#     )\n",
        "\n",
        "# # test dataset\n",
        "# transformed_test_data = datasets.MNIST(root='data', train=False,\n",
        "#                                       download=True, transform=transform)\n",
        "\n",
        "# # prepare data loader\n",
        "# transformed_test_loader = torch.utils.data.DataLoader(transformed_test_data, \n",
        "#                                                       batch_size=batch_size,\n",
        "#                                                       num_workers=num_workers)\n",
        "\n",
        "# # call test function and get reconstructed images\n",
        "# _, images, reconstructions = test(capsule_net, transformed_test_loader)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "FvWOq_spKmiT"
      },
      "outputs": [],
      "source": [
        "# # obtain one batch of test images\n",
        "# dataiter = iter(transformed_test_loader)\n",
        "# images, labels = dataiter.next()\n",
        "# images = images.numpy()\n",
        "\n",
        "# # plot the images in the batch, along with the corresponding labels\n",
        "# fig = plt.figure(figsize=(25, 4))\n",
        "# for idx in np.arange(batch_size):\n",
        "#     ax = fig.add_subplot(2, batch_size/2, idx+1, xticks=[], yticks=[])\n",
        "#     ax.imshow(np.squeeze(images[idx]), cmap='gray')\n",
        "#     # print out the correct label for each image\n",
        "#     # .item() gets the value contained in a Tensor\n",
        "#     ax.set_title(str(labels[idx].item()))"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "bBGJlhw7KmiZ"
      },
      "outputs": [],
      "source": [
        "# original input images\n",
        "display_images(images, reconstructions)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": true,
        "id": "v6R-RfA9Kmia"
      },
      "outputs": [],
      "source": [
        "def vector_analysis(capsule_net, x, select_idx=1):\n",
        "    '''Generates perturbed iage reconstructions given some digit capsule outputs.\n",
        "       param capsule_net: trained capsule network\n",
        "       param x: a batch of digit capsule outputs\n",
        "       param select_idx: selects which image in a batch to analyze, default = 1 \n",
        "       return: list of perturbed, reconstructed images\n",
        "       '''\n",
        "    \n",
        "    classes = (x ** 2).sum(dim=-1) ** 0.5\n",
        "    classes = F.softmax(classes, dim=-1)\n",
        "\n",
        "    # find the capsule with the maximum vector length\n",
        "    # here, vector length indicates the probability of a class' existence\n",
        "    _, max_length_indices = classes.max(dim=1)\n",
        "\n",
        "    # create a sparse class matrix\n",
        "    sparse_matrix = torch.eye(10) # 10 is the number of classes\n",
        "    if TRAIN_ON_GPU:\n",
        "        sparse_matrix = sparse_matrix.cuda()\n",
        "    # get the class scores from the \"correct\" capsule\n",
        "    y = sparse_matrix.index_select(dim=0, index=max_length_indices.data)\n",
        "\n",
        "    # create reconstructed pixels\n",
        "    x = x * y[:, :, None]\n",
        "    \n",
        "    # flatten image into a vector shape (batch_size, vector_dim)\n",
        "    flattened_x = x.reshape(x.size(0), -1)\n",
        "    # select a single image from a batch to work with\n",
        "    flattened_x = flattened_x[select_idx]\n",
        "    \n",
        "    # track reconstructed images\n",
        "    reconstructed_ims = []\n",
        "    # values to change *one* vector dimension by\n",
        "    perturb_range = np.arange(-0.25, 0.30, 0.05)\n",
        "    \n",
        "    # iterate through 16 vector dims\n",
        "    for k in range(16):\n",
        "        # create a copy of flattened_x to modify\n",
        "        if TRAIN_ON_GPU:\n",
        "            transformed_x = torch.zeros(*flattened_x.size()).cuda()\n",
        "        else:\n",
        "            transformed_x = torch.zeros(*flattened_x.size())\n",
        "            \n",
        "        transformed_x[:] = flattened_x[:]\n",
        "        # iterate through each perturbation value\n",
        "        for j in range(len(perturb_range)):\n",
        "            # for each capsule output\n",
        "            for i in range(10):\n",
        "                transformed_x[k+(16*i)] = flattened_x[k+(16*i)]+perturb_range[j]\n",
        "\n",
        "            # create reconstructed images\n",
        "            reconstructions = capsule_net.decoder.linear_layers(transformed_x)\n",
        "            # reshape into 28x28 image, (batch_size, depth, x, y)\n",
        "            reconstructions = reconstructions.view(-1, 1, 28, 28)\n",
        "            reconstructed_ims.append(reconstructions)\n",
        "    \n",
        "    # return final list of reconstructed ims    \n",
        "    return reconstructed_ims\n",
        "        "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": true,
        "id": "D7A4gthwKmib"
      },
      "outputs": [],
      "source": [
        "# call function and get perturbed reconstructions\n",
        "reconstructed_ims = vector_analysis(capsule_net, caps_output, select_idx=1)"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "6Thuc2QGKmic"
      },
      "outputs": [],
      "source": [
        "fig = plt.figure(figsize=(10, 20)) # define figsize\n",
        "\n",
        "# display all ims\n",
        "for idx in range(len(reconstructed_ims)):\n",
        "    # convert to numpy images\n",
        "    image = reconstructed_ims[idx]\n",
        "    image = image.detach().cpu().numpy()\n",
        "    # display 16 rows of images\n",
        "    ax = fig.add_subplot(16, len(reconstructed_ims)/16, idx+1, xticks=[], yticks=[])\n",
        "    ax.imshow(image.squeeze(), cmap='gray')\n",
        "    "
      ]
    }
  ],
  "metadata": {
    "accelerator": "GPU",
    "colab": {
      "collapsed_sections": [],
      "name": "Capsule_Network_complete.ipynb",
      "provenance": []
    },
    "kernelspec": {
      "display_name": "Python 3",
      "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.9.1"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}