Ann_manual_Project.ipynb

{
  "cells": [
    {
      "attachments": {},
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "##  Download & Import required pacages"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "!pip install -q --upgrade pip\n",
        "!pip install -q matplotlib\n",
        "!pip install -q math\n",
        "!pip install -q numpy\n",
        "!pip install -q copy\n",
        "!pip install -q keras"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 2,
      "metadata": {
        "id": "oAffhrNjKkJ3"
      },
      "outputs": [],
      "source": [
        "import matplotlib.pyplot as plt\n",
        "from keras.datasets import mnist\n",
        "import numpy as np\n",
        "import math\n",
        "import copy"
      ]
    },
    {
      "attachments": {},
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Load Database and normalize"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 3,
      "metadata": {
        "id": "VMZAMJsT2oVS"
      },
      "outputs": [],
      "source": [
        "def label_set(y, last_layer_lengh):\n",
        "  label_list = []\n",
        "  for label_value in y:\n",
        "    label_vector = np.zeros((last_layer_lengh,1))\n",
        "    label_vector[label_value, 0] = 1\n",
        "    label_list.append(label_vector)\n",
        "\n",
        "  return label_list"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 14,
      "metadata": {},
      "outputs": [],
      "source": [
        "def load_train_set():\n",
        "    # Reading The Train Set\n",
        "    train_images_file = open('train-images.idx3-ubyte', 'rb')\n",
        "    train_images_file.seek(4)\n",
        "    num_of_train_images = int.from_bytes(train_images_file.read(4), 'big')\n",
        "    train_images_file.seek(16)\n",
        "\n",
        "    train_labels_file = open('train-labels.idx1-ubyte', 'rb')\n",
        "    train_labels_file.seek(8)\n",
        "\n",
        "    train_set = []\n",
        "    for n in range(num_of_train_images):\n",
        "        image = np.zeros((784, 1))\n",
        "        for i in range(784):\n",
        "            image[i, 0] = int.from_bytes(train_images_file.read(1), 'big') / 256\n",
        "        \n",
        "        label_value = int.from_bytes(train_labels_file.read(1), 'big')\n",
        "        label = np.zeros((10, 1))\n",
        "        label[label_value, 0] = 1\n",
        "        \n",
        "        train_set.append((image, label))\n",
        "    return train_set\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 15,
      "metadata": {},
      "outputs": [],
      "source": [
        "def load_test_set(): \n",
        "    # Reading The Test Set\n",
        "    test_images_file = open('t10k-images.idx3-ubyte', 'rb')\n",
        "    test_images_file.seek(4)\n",
        "\n",
        "    test_labels_file = open('t10k-labels.idx1-ubyte', 'rb')\n",
        "    test_labels_file.seek(8)\n",
        "\n",
        "    num_of_test_images = int.from_bytes(test_images_file.read(4), 'big')\n",
        "    test_images_file.seek(16)\n",
        "\n",
        "    test_set = []\n",
        "    for n in range(num_of_test_images):\n",
        "        image = np.zeros((784, 1))\n",
        "        for i in range(784):\n",
        "            image[i] = int.from_bytes(test_images_file.read(1), 'big') / 256\n",
        "        \n",
        "        label_value = int.from_bytes(test_labels_file.read(1), 'big')\n",
        "        label = np.zeros((10, 1))\n",
        "        label[label_value, 0] = 1\n",
        "        \n",
        "        test_set.append((image, label))\n",
        "    \n",
        "    return test_set\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 16,
      "metadata": {},
      "outputs": [],
      "source": [
        "train_set = load_train_set()\n",
        "test_set = load_test_set()"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 9,
      "metadata": {},
      "outputs": [],
      "source": [
        "train_X = [item[0] for item in train_set]\n",
        "train_y = [item[1] for item in train_set]\n",
        "test_X = [item[0] for item in test_set]\n",
        "test_y = [item[1] for item in test_set]"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "def create_batch_list( train_X, train_y, batch_size):\n",
        "  train_X_batch_list = [train_X[i:i+batch_size] for i in range(0, len(train_X), batch_size)]\n",
        "  train_y_batch_list = [train_y[i:i+batch_size] for i in range(0, len(train_X), batch_size)]\n",
        "  batch_list = [zip(train_X_batch_list[i], train_y_batch_list[i]) for i in range(len(train_y_batch_list))]\n",
        "  return batch_list"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 18,
      "metadata": {
        "id": "70wKfaVXdcMf"
      },
      "outputs": [],
      "source": [
        "from collections import Counter\n",
        "def check_data_distribution():\n",
        "  print(\"train classes: \", set(train_y))\n",
        "  print(\"test classes: \", set(test_y))\n",
        "  train_dist = Counter(train_y)\n",
        "  print(train_dist)\n",
        "  test_dist = Counter(test_y)\n",
        "  print(test_dist)\n",
        "  return train_dist, test_dist"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 19,
      "metadata": {
        "id": "LEtqFV-xLCOB"
      },
      "outputs": [],
      "source": [
        "def barchart_ploter(dist_dict):\n",
        "  bars = list(dist_dict.keys())\n",
        "  hight = list(dist_dict.values())\n",
        "  plt.bar(bars, hight, width = 0.4)\n",
        "  plt.xlabal(\"numbers\")\n",
        "  plt.ylabal(\"distribution\")\n",
        "  plt.show()"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 10,
      "metadata": {
        "id": "BuDbFgjmLGFZ"
      },
      "outputs": [],
      "source": [
        "def check_data(index):\n",
        "  print(\"label: \",train_y[index])\n",
        "  plt.imshow(train_X[index])"
      ]
    },
    {
      "attachments": {},
      "cell_type": "markdown",
      "metadata": {
        "id": "lJw9kCHrLTci"
      },
      "source": [
        "## Forward Phase"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 32,
      "metadata": {
        "id": "UQvUkJU37ik5"
      },
      "outputs": [],
      "source": [
        "def W_init(row, col):\n",
        "    # random_w_list = np.random.uniform(low = 0, high = 1, size = row *col)\n",
        "    # random_w_matrix = random_w_list.reshape(row,col)\n",
        "    random_w_list = np.random.normal(0, 1, size = (row,col))\n",
        "    return random_w_list\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 33,
      "metadata": {
        "id": "PhhWVTDo7ik6"
      },
      "outputs": [],
      "source": [
        "def b_init(row):\n",
        "    # np.random.uniform(low = 0.0, high = 1.0, size = row)\n",
        "    random_b_list = np.zeros((row,1)) \n",
        "    random_b_matrix = random_b_list.reshape(row,1)\n",
        "    return random_b_matrix"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 34,
      "metadata": {
        "id": "MHmpEy8b7ik7"
      },
      "outputs": [],
      "source": [
        "def init_matrix(matrix_layer_list):\n",
        "    W_list = []\n",
        "    b_list = []\n",
        "    for layer_index in range(len(matrix_layer_list)-1):\n",
        "      W_list.append(W_init(matrix_layer_list[layer_index+1], matrix_layer_list[layer_index]))\n",
        "      b_list.append(b_init(matrix_layer_list[layer_index+1]))\n",
        "\n",
        "    return W_list, b_list"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 35,
      "metadata": {
        "id": "of6WqjHO_IF3"
      },
      "outputs": [],
      "source": [
        "def Sigmoid(x):\n",
        "  result = 1/(1+ np.exp(-x))\n",
        "  return result"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 36,
      "metadata": {
        "id": "wB5mf_u9KspI"
      },
      "outputs": [],
      "source": [
        "def next_leyar_activation_calculator(matrix_W, b_vector, a_vector):\n",
        "  resault_vector = Sigmoid((matrix_W @ a_vector) + b_vector)\n",
        "  return resault_vector "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 37,
      "metadata": {
        "id": "w1sCeBnkKySy"
      },
      "outputs": [],
      "source": [
        "def forward(input_data, W_list, b_list):\n",
        "  a_list = [input_data]\n",
        "  for i in range(len(W_list)):\n",
        "    res =next_leyar_activation_calculator(W_list[i], b_list[i], a_list[i])\n",
        "    a_list.append(res)\n",
        "  return a_list"
      ]
    },
    {
      "attachments": {},
      "cell_type": "markdown",
      "metadata": {
        "id": "Pwggts_pLir1"
      },
      "source": [
        "## Backward Phase"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 39,
      "metadata": {
        "id": "Zh7PPB6v7KaG"
      },
      "outputs": [],
      "source": [
        "def Sigmoid_derivate(x):\n",
        "  s = Sigmoid(x)\n",
        "  result = s * (1-s)\n",
        "  return result"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 40,
      "metadata": {
        "id": "B1VWmV_SAYm3"
      },
      "outputs": [],
      "source": [
        "def calculate_layer_dervatives(a_vector, b_vector, matrix_w, dcost_da_vector, layer_num):\n",
        "  z = (matrix_w @ a_vector) + b_vector\n",
        "  b_derivative =  dcost_da_vector * Sigmoid_derivate(z)\n",
        "  row, col = matrix_w.shape\n",
        "  a_derivative = np.zeros((col,1))\n",
        "  w_derivative = b_derivative @ np.transpose(a_vector)\n",
        "  if layer_num  != 0:\n",
        "    a_derivative = np.transpose(matrix_w) @ b_derivative\n",
        "\n",
        "  return b_derivative, a_derivative, w_derivative"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 41,
      "metadata": {
        "id": "14Zudi0WENL6"
      },
      "outputs": [],
      "source": [
        "def backpropagation_dev_calculator(a_list, b_list, w_list, label):\n",
        "  dcost_da_last_layer = 2*(a_list[-1] - label)\n",
        "  layers_num = len(b_list)\n",
        "  dcost_da_list = [0 for i in range(len(b_list))]\n",
        "  dcost_da_list[layers_num-1] = dcost_da_last_layer\n",
        "  \n",
        "  db_list = []\n",
        "  dw_list = []\n",
        "  for i in range(layers_num-1, -1, -1):\n",
        "    db, da, dw = calculate_layer_dervatives(a_list[i], b_list[i], w_list[i], dcost_da_list[i], i)\n",
        "    db_list.append(db)\n",
        "    dw_list.append(dw)\n",
        "    dcost_da_list[i-1] = da\n",
        "  \n",
        "  return db_list[::-1], dw_list[::-1]\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 42,
      "metadata": {
        "id": "1HNSp6kBENVv"
      },
      "outputs": [],
      "source": [
        "def backpropagation_grad_calculator(w_list, b_list, batch):\n",
        "\n",
        "  Gw_list = [np.zeros(mat.shape) for mat in w_list]\n",
        "  Gb_list = [np.zeros(mat.shape) for mat in b_list]\n",
        "  for image, label in batch:\n",
        "    # print(\"back_grad_Cal\", image.shape)\n",
        "    a_list = forward(image, w_list, b_list)\n",
        "    db_list, dw_list = backpropagation_dev_calculator(a_list, b_list, w_list, label)\n",
        "    for index in range(len(Gw_list)):\n",
        "      Gw_list[index] += dw_list[index]\n",
        "      Gb_list[index] += db_list[index]\n",
        "      \n",
        "  return Gw_list, Gb_list\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 44,
      "metadata": {
        "id": "Z7Zr8E9isYvn"
      },
      "outputs": [],
      "source": [
        "def cost(w_list, b_list, train_X, train_y):\n",
        "  cost_value = 0 \n",
        "  for sample, label in zip(train_X, train_y):\n",
        "    \n",
        "    a_list = forward(sample,w_list, b_list)\n",
        "    predicted_label = a_list[-1]\n",
        "    cost_value += ((predicted_label - label)**2) \n",
        "\n",
        "  return cost_value/len(train_X)\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 45,
      "metadata": {
        "id": "sQ0-Z9c_83I2"
      },
      "outputs": [],
      "source": [
        "def calculate_cost(batch, weights_matrix_list, bias_list):\n",
        "    cost = 0\n",
        "    for image, label in batch:\n",
        "        predicted_output = forward(image, weights_matrix_list, bias_list)\n",
        "        cost += (np.linalg.norm(predicted_output[-1] - label))\n",
        "    return cost/len(batch)\n",
        "    "
      ]
    },
    {
      "attachments": {},
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Train and evaluation"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 46,
      "metadata": {
        "id": "UwBOSN1Ha8IT"
      },
      "outputs": [],
      "source": [
        "from tqdm import tqdm\n",
        "def SGD(matrix_layer_list, learning_rate, batch_size, epoch, train_X, train_y):\n",
        "  w_list, b_list = init_matrix(matrix_layer_list)\n",
        "  cost_list = []\n",
        "  for i in tqdm(range(epoch)):\n",
        "    batch_list = create_batch_list( train_X, train_y, batch_size)\n",
        "    for batch in batch_list:\n",
        "      Gw_list, Gb_list = backpropagation_grad_calculator(w_list, b_list, batch)\n",
        "      for index in range(len(w_list)):\n",
        "        w_list[index] -= (learning_rate * Gw_list[index])/ batch_size\n",
        "        b_list[index] -= (learning_rate * Gb_list[index])/ batch_size\n",
        "    cost_list.append(cost(w_list, b_list,train_X, train_y))\n",
        "  return cost_list, w_list, b_list"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 47,
      "metadata": {
        "id": "hV4ff-G8odvE"
      },
      "outputs": [],
      "source": [
        "def cal_accuracy(test_X, test_y, w_list, b_list):\n",
        "  correct_prediction_num = 0\n",
        "  for img, label in zip(test_X, test_y):\n",
        "    a_list = forward(img, w_list, b_list )\n",
        "    predicted_label = a_list[-1]\n",
        "    if np.argmax(predicted_label, axis=0)== np.argmax(label, axis=0):\n",
        "      correct_prediction_num += 1\n",
        "  return (correct_prediction_num/len(test_X))*100   "
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 48,
      "metadata": {
        "id": "y6giI4KVvgrc"
      },
      "outputs": [],
      "source": [
        "def show_result(test_X, test_y, w_list, b_list, cost_list):\n",
        "  acc = cal_accuracy(test_X, test_y, w_list, b_list)\n",
        "  print(acc)\n",
        "  norm_list = [(np.linalg.norm(x)**2) for x in cost_list]\n",
        "  plt.plot(norm_list)\n",
        "  plt.show()\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 49,
      "metadata": {
        "id": "zvzfRqb6ENgy"
      },
      "outputs": [],
      "source": [
        "def main(epoch_num):\n",
        "  matrix_layer_list = [28*28, 8*8, 16, 10]\n",
        "  \n",
        "  cost_list, w_list1, b_list1 = SGD(matrix_layer_list, 2, 32, epoch_num, train_X, train_y)\n",
        "  show_result(test_X, test_y, w_list1, b_list1, cost_list)\n",
        "  \n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": 50,
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 1000
        },
        "id": "wVtD42-driFF",
        "outputId": "a58351f5-4ddb-4fd9-e167-a2f79c50a943"
      },
      "outputs": [
        {
          "name": "stderr",
          "output_type": "stream",
          "text": [
            "100%|██████████| 20/20 [10:17<00:00, 30.87s/it]\n"
          ]
        },
        {
          "name": "stdout",
          "output_type": "stream",
          "text": [
            "94.78\n"
          ]
        },
        {
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYYAAAD4CAYAAADo30HgAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAjRUlEQVR4nO3deXRcd3338fdX+zaWLFmS9yW25MRZHBvVhDSBUJPghAaTkoKBggtp05QEEpanTcpzKE9O6QlpSVraBBogJFAaZ4G2JgSy2AECIU5kx3G8W97t2JJsy7JlWfv3+WOunJE8skfWcqWZz+scnbnzu7975zuj0Xx07/3dO+buiIiIdEsLuwARERlZFAwiItKDgkFERHpQMIiISA8KBhER6SEj7AIGw7hx43z69OlhlyEiMqqsXr36kLuX9m5PimCYPn061dXVYZchIjKqmNnueO3alSQiIj0oGEREpAcFg4iI9KBgEBGRHhQMIiLSg4JBRER6UDCIiEgPKR0ML26u48Ff1YRdhojIiJLSwfC7mkN8a8U2urr0nRQiIt1SOhgqyyO0tHext6E57FJEREaMlA6GivICALYcPB5yJSIiI0eKB0MEgG11TSFXIiIycqR0MBRkZzCpKJettdpiEBHpltLBAFBZXsDWWm0xiIh0UzCUR9he30RHZ1fYpYiIjAgpHwwV5RHaOrrYfUQjk0REQMFAZTAyaZuOM4iIAAoGZpVFg0HHGUREohIKBjNbZGZbzKzGzO6MMz/bzB4P5q8ys+kx8+4K2reY2fuDtilm9qKZbTSzDWZ2e0z/r5nZfjNbG/xcNwjPs095WRlMKdbIJBGRbmf9zmczSwceAK4G9gGvmdlyd98Y0+0moMHdZ5nZEuAbwEfNbA6wBLgQmAi8YGaVQAfwJXdfY2YRYLWZPR+zzvvd/Z8H60mezezyCNu0xSAiAiS2xbAAqHH3He7eBiwDFvfqsxh4NJh+ClhoZha0L3P3VnffCdQAC9z9gLuvAXD348AmYNLAn865qSiPsONQE+0amSQiklAwTAL2xtzfx+kf4qf6uHsH0AiUJLJssNtpHrAqpvk2M1tnZg+b2dh4RZnZzWZWbWbV9fX1CTyNvlWWF9De6ew6dGJA6xERSQahHnw2swLgJ8Ad7n4saP42MBO4FDgAfDPesu7+kLtXuXtVaWnpgOqoKIteGkMHoEVEEguG/cCUmPuTg7a4fcwsAygEDp9pWTPLJBoKP3b3n3Z3cPdad+909y7gu0R3ZQ2pWWUFpBk6AC0iQmLB8BpQYWYzzCyL6MHk5b36LAeWBtM3Aivd3YP2JcGopRlABfBqcPzh+8Amd78vdkVmNiHm7g3A+v4+qf7KyUxnWkk+2+oUDCIiZx2V5O4dZnYb8CyQDjzs7hvM7G6g2t2XE/2Q/5GZ1QBHiIYHQb8ngI1ERyLd6u6dZnYF8EngTTNbGzzU37n7M8C9ZnYp4MAu4K8G7dmeQUWZrpkkIgIJBANA8IH9TK+2r8ZMtwB/2seyXwe+3qvtt4D10f+TidQ02CrLI6zYXEdrRyfZGelhlCAiMiKk/JnP3SrKC+jscnZqZJKIpDgFQ6CyXCOTRERAwXDKeaX5pKeZLqYnIilPwRDIzkhnekmehqyKSMpTMMSo1DWTREQUDLEqyiPsOnyClvbOsEsREQmNgiFGZXkBXQ7b67XVICKpS8EQo3tkknYniUgqUzDEmF6ST0aa6QC0iKQ0BUOMrIw0zivN17kMIpLSFAy9VJRHdDE9EUlpCoZeKssi7DnSzMk2jUwSkdSkYOilsrwAd6ip0+4kEUlNCoZeKk5dM0m7k0QkNSkYeplekkdWehpbdZxBRFKUgqGXjPToyCSdyyAiqUrBEEdleYQtB7XFICKpScEQR2V5AfuPnuREa0fYpYiIDDsFQxzdB6C3aWSSiKQgBUMclRqZJCIpTMEQx9TiPLIz0vRtbiKSkhQMcaSnGbPKCnTNJBFJSQqGPlSWR7QrSURSkoKhDxXlBRxobOFYS3vYpYiIDCsFQx8qy/SlPSKSmhQMfXj729y0O0lEUouCoQ+Tx+aSm5muA9AiknIUDH1ISzMqygv0pT0iknIUDGdQUaZrJolI6lEwnEFleQF1x1tpbNbIJBFJHQqGMzh1aQztThKRFJJQMJjZIjPbYmY1ZnZnnPnZZvZ4MH+VmU2PmXdX0L7FzN4ftE0xsxfNbKOZbTCz22P6F5vZ82a2LbgdOwjP85xUjtc1k0Qk9Zw1GMwsHXgAuBaYA3zMzOb06nYT0ODus4D7gW8Ey84BlgAXAouAB4P1dQBfcvc5wGXArTHrvBNY4e4VwIrgfigmFuZQkJ2hcxlEJKUkssWwAKhx9x3u3gYsAxb36rMYeDSYfgpYaGYWtC9z91Z33wnUAAvc/YC7rwFw9+PAJmBSnHU9CnzonJ7ZIDCLXjNJB6BFJJUkEgyTgL0x9/fx9of4aX3cvQNoBEoSWTbY7TQPWBU0lbv7gWD6IFAerygzu9nMqs2sur6+PoGncW4qNWRVRFJMqAefzawA+Alwh7sf6z3f3R3weMu6+0PuXuXuVaWlpUNWY2V5hENNbRw50TZkjyEiMpIkEgz7gSkx9ycHbXH7mFkGUAgcPtOyZpZJNBR+7O4/jelTa2YTgj4TgLpEn8xQqNCX9ohIikkkGF4DKsxshpllET2YvLxXn+XA0mD6RmBl8N/+cmBJMGppBlABvBocf/g+sMnd7zvDupYC/9vfJzWYZuuaSSKSYjLO1sHdO8zsNuBZIB142N03mNndQLW7Lyf6If8jM6sBjhAND4J+TwAbiY5EutXdO83sCuCTwJtmtjZ4qL9z92eAe4AnzOwmYDfwkUF8vv1WPiabSE6GrpkkIinDov/Yj25VVVVeXV09ZOv/8LdfJj3NeOKv3jVkjyEiMtzMbLW7V/Vu15nPCagsL2Bb7XGSIURFRM5GwZCAirIIDc3tHGrSyCQRSX4KhgToS3tEJJUoGBJQOb4A0JBVEUkNCoYElBZkU5SXydY6jUwSkeSnYEiAmVFZFmGrrpkkIilAwZCgivICtmpkkoikAAVDgirLIxxr6aDueGvYpYiIDCkFQ4Iqdc0kEUkRCoYEVZZ3j0zSAWgRSW4KhgSVFGRTkp+lA9AikvQUDP1QUV7AVn1pj4gkOQVDP1SWR6ipbdLIJBFJagqGfqgoj3C8tYMDjS1hlyIiMmQUDP0wWyOTRCQFKBj6oXtk0jaNTBKRJKZg6IeivCxKI9ls0RaDiCQxBUM/dX9pj4hIslIw9FNFWYRtdU10dWlkkogkJwVDP1WWR2hu62T/0ZNhlyIiMiQUDP00O/jSnm060U1EkpSCoZ9mlUWHrG45qJFJIpKcFAz9VJibyfgxOToALSJJS8FwDnTNJBFJZgqGc1BZHqFGI5NEJEkpGM7B7PIILe1d7G1oDrsUEZFBp2A4BxX60h4RSWIKhnNQoYvpiUgSUzCcg4LsDCYV5SoYRCQpKRjOUUV5gXYliUhSUjCco8ryCNvrm+jUyCQRSTIJBYOZLTKzLWZWY2Z3xpmfbWaPB/NXmdn0mHl3Be1bzOz9Me0Pm1mdma3vta6vmdl+M1sb/Fw3gOc3ZCrLI7R1dLH78ImwSxERGVRnDQYzSwceAK4F5gAfM7M5vbrdBDS4+yzgfuAbwbJzgCXAhcAi4MFgfQCPBG3x3O/ulwY/z/TvKQ2PSo1MEpEklcgWwwKgxt13uHsbsAxY3KvPYuDRYPopYKGZWdC+zN1b3X0nUBOsD3f/DXBkEJ5DKGaVdQeDDkCLSHJJJBgmAXtj7u8L2uL2cfcOoBEoSXDZeG4zs3XB7qax8TqY2c1mVm1m1fX19QmscnDlZWUwpVgjk0Qk+YzEg8/fBmYClwIHgG/G6+TuD7l7lbtXlZaWDmN5b6ssi+j7n0Uk6SQSDPuBKTH3JwdtcfuYWQZQCBxOcNke3L3W3TvdvQv4LsGup5Fo/rSxbKk9zqYDx8IuRURk0CQSDK8BFWY2w8yyiB5MXt6rz3JgaTB9I7DS3T1oXxKMWpoBVACvnunBzGxCzN0bgPV99Q3bn71zGpGcDO5/fmvYpYiIDJqzBkNwzOA24FlgE/CEu28ws7vN7INBt+8DJWZWA3wRuDNYdgPwBLAR+CVwq7t3ApjZY8Dvgdlmts/MbgrWda+ZvWlm64D3Al8YpOc66ArzMvmLK87juY21vLmvMexyREQGhUX/sR/dqqqqvLq6OpTHPt7SzpX3vsi8KUX84NMjdq+XiMhpzGy1u1f1bh+JB59HlUhOJn/17pm8uKWe1bsbwi5HRGTAFAyDYOnl0xhXkMV9z28JuxQRkQFTMAyCvKwMbnnPTH5Xc5jfbz8cdjkiIgOiYBgkf3bZNMrHZHPf81tIhuM2IpK6FAyDJCczndveO4vXdjXw0rZDYZcjInLOFAyD6CN/MIVJRbl88/mt2moQkVFLwTCIsjPS+dwfzeKNvUdZubku7HJERM6JgmGQffgdk5lWksc3n9tKl77ER0RGIQXDIMtMT+P2hRVsPHCMZzccDLscEZF+UzAMgcWXTmJmaT73v7BVX/0pIqOOgmEIpKcZd7yvkq21TTy97q2wyxER6RcFwxD5wMUTOH98hH95YRsdnV1hlyMikjAFwxBJSzO+cHUlOw+d4L9fP+NXUIiIjCgKhiF0zZxyLp5UyLdWbqOtQ1sNIjI6KBiGkJnxxasr2XvkJE+u3nv2BURERgAFwxC7anYp86cW8e8ra2hp7wy7HBGRs1IwDDEz40vXzOZAYwvLXt0TdjkiImelYBgGl88s4Z0zinngV9s52aatBhEZ2RQMw6B7q6H+eCs/emVX2OWIiJyRgmGYLJhRzJUV4/jOr3fQ1NoRdjkiIn1SMAyjL10zmyMn2nj05V1hlyIi0icFwzC6dEoRC88v4z9+vZ3Gk+1hlyMiEpeCYZh94epKjrV08P3f7gy7FBGRuBQMw+yiSYVce9F4Hv7tThpOtIVdjojIaRQMIfjC1ZWcaOvgoZd2hF2KiMhpFAwhqCyP8MG5E3nkd7uoP94adjkiIj0oGEJy+8IKWjs6+c6vt4ddiohIDwqGkJxXWsCH50/mh7/fxbp9R8MuR0TkFAVDiP7uugsoi+Tw2R+vobFZw1dFZGRQMIRobH4W//7xedQea+FLT76Bu74fWkTCp2AI2bypY7nr2gt4YVMt33tJ5zaISPgSCgYzW2RmW8ysxszujDM/28weD+avMrPpMfPuCtq3mNn7Y9ofNrM6M1vfa13FZva8mW0LbscO4PmNCp/+w+ksunA89/xyM9W7joRdjoikuLMGg5mlAw8A1wJzgI+Z2Zxe3W4CGtx9FnA/8I1g2TnAEuBCYBHwYLA+gEeCtt7uBFa4ewWwIrif1MyMe//0EiYV5XLbf73O4SYNYRWR8CSyxbAAqHH3He7eBiwDFvfqsxh4NJh+ClhoZha0L3P3VnffCdQE68PdfwPE+/c4dl2PAh9K/OmMXmNyMnnwE/M50tzGF554g64uHW8QkXAkEgyTgNgvLN4XtMXt4+4dQCNQkuCyvZW7+4Fg+iBQHq+Tmd1sZtVmVl1fX5/A0xj5LppUyN9fP4ffbK3nwV/VhF2OiKSoEX3w2aPDdOL+6+zuD7l7lbtXlZaWDnNlQ+fjC6ay+NKJ3Pf8Vl7efijsckQkBSUSDPuBKTH3JwdtcfuYWQZQCBxOcNneas1sQrCuCUBdAjUmDTPjH2+4mBnj8vn8Y2upO94SdkkikmISCYbXgAozm2FmWUQPJi/v1Wc5sDSYvhFYGfy3vxxYEoxamgFUAK+e5fFi17UU+N8Eakwq+dkZPPiJd9DU2s7nH3udTh1vEJFhdNZgCI4Z3AY8C2wCnnD3DWZ2t5l9MOj2faDEzGqALxKMJHL3DcATwEbgl8Ct7t4JYGaPAb8HZpvZPjO7KVjXPcDVZrYNeF9wP+XMHh/hHz50Ma/sOMK/vLA17HJEJIVYMpxtW1VV5dXV1WGXMST+5qk3eKJ6H498+g+4anZZ2OWISBIxs9XuXtW7fUQffBa4e/FFnD8+whceX8tbR0+GXY6IpAAFwwiXk5nOg5+YT1tHF5977HXaO7vCLklEkpyCYRQ4r7SAez58Cat3N3DvLzeHXY6IJDkFwyhx/dyJfOpd0/juSzt5bsPBsMsRkSSmYBhFvvKBC7h4UiFfevIN9hxuDrscEUlSCoZRJDsjerwB4Nb/WkNrR2fIFYlIMlIwjDJTivP45p/O5c39jXz955vCLkdEkpCCYRS65sLx/OWVM/jh73fzszfeCrscEUkyCoZR6m8Wnc/8qUV8+ck3eOzVPfpaUBEZNAqGUSozPY3vfqqKBTOKueunb/L5ZWs53tIedlkikgQUDKNYSUE2j356Af/n/bP5+bq3uP7ffsv6/Y1hlyUio5yCYZRLSzNufe8slt38Llrau/iTB1/m0Zd3adeSiJwzBUOSWDCjmGduv5IrKsbx98s38Nf/uYbGk9q1JCL9p2BIIsX5WXzvU1V85boLeGFTLR/41kus3Xs07LJEZJRRMCSZtDTjL999Hk/c8i7c4cZvv8z3XtqhXUsikjAFQ5KaP3Usz3z+ShZeUMY//HwTf/FoNQ0n2sIuS0RGAQVDEivMy+Q7f/YO/t8HL+SlbYe47lsvUb3rSNhlicgIp2BIcmbG0sun85O/vpysjDQ++tArPPBiDV36HmkR6YOCIUVcPLmQpz93BddeNJ5/enYLS3/wKoeaWsMuS0RGIAVDConkZPJvH5vHP95wMa/uPMK1//oSL9ccCrssERlhFAwpxsz4+Dun8j+3/iGRnAw+/r1V/OUPq9l88FjYpYnICKFgSFEXTBjD05+7gi9fU8krOw5z7b++xO3LXmfXoRNhlyYiIbNkGN9eVVXl1dXVYZcxah1tbuM/frODH/xuJ+2dzkeqpvD5hbOYUJgbdmkiMoTMbLW7V53WrmCQbnXHW3jwxe38eNVuzIxPXjaNz141k5KC7LBLE5EhoGCQhO1raOZbK7bx1Op95Gam85krZvAXV55HYW5m2KWJyCBSMEi/1dQ1cf8LW/n5ugMU5mZyy3tmsvTyaeRlZYRdmogMAgWDnLP1+xu57/mtrNxcR2kkm9veO4slC6aQnZEedmkiMgAKBhmw6l1HuPfZLby68wiTinK5430V3DBvEhnpGtwmMhopGGRQuDsvbTvEPz+3hXX7GplWkseHLp3E9XMnMKssEnZ5ItIPCgYZVO7OsxtqeeTlnazaeQR3OH98hOvnTuT6SyYytSQv7BJF5CwUDDJk6o618MybB/jZugOs3t0AwNzJhVw/dyIfuGSCzocQGaEGFAxmtgj4VyAd+J6739NrfjbwQ+AdwGHgo+6+K5h3F3AT0Al83t2fPdM6zewR4D1A97fa/7m7rz1TfQqGkWNfQzM/X3eAn617i/X7o5fZ+IPpY7l+7kSuvWgCpRGdEyEyUpxzMJhZOrAVuBrYB7wGfMzdN8b0+SxwibvfYmZLgBvc/aNmNgd4DFgATAReACqDxeKuMwiGp939qUSfnIJhZNp56ARPv/EWP1v3Fltrm0gzuHzmOP74kgksumg8RXlZYZcoktL6CoZEBqQvAGrcfUewomXAYmBjTJ/FwNeC6aeAfzczC9qXuXsrsNPMaoL1kcA6ZZSbMS6fzy2s4HMLK9hy8DhPr3uLn73xFnf+9E3+7/+s592VpXzg4gn80flljM1XSIiMFIkEwyRgb8z9fcA7++rj7h1m1giUBO2v9Fp2UjB9pnV+3cy+CqwA7gyCpQczuxm4GWDq1KkJPA0J0+zxEWaPn80Xr65k/f5jp0Ji5eY60gzeMW0sCy8o530XlDGztIDo/xUiEoaReArrXcBBIAt4CPhb4O7endz9oWA+VVVVo/8IeoowMy6eXMjFkwv520Xn8+b+RlZsquWFTXXc84vN3POLzUwryWPh+eUsvKCMBTOKydR5EiLDKpFg2A9Mibk/OWiL12efmWUAhUQPQp9p2bjt7n4gaGs1sx8AX06gRhmF0tKMuVOKmDuliC9eM5u3jp5kxeY6Vm6q5T9X7ebh3+0kkp3Bu2eX8r4LyriqUrucRIZDIsHwGlBhZjOIfngvAT7eq89yYCnwe+BGYKW7u5ktB/7LzO4jevC5AngVsL7WaWYT3P1AcIziQ8D6gT1FGS0mFuXyycum8cnLptHc1sFvtx1ixaY6Vmyu4+frDpBmUDWtmD+6oEy7nESGUKLDVa8D/oXo0NKH3f3rZnY3UO3uy80sB/gRMA84AiyJObD8FeAzQAdwh7v/oq91Bu0rgVKi4bEWuMXdm85Un0YlJbeuLmddzC6nTQeiw2CnleTxh7PGUVqQzdi8TMbmZ1GUlxWdzsuiKC+TguwMhYdIH3SCmySN7l1OKzbVsmZ3A8daOvrsm5luFOb2DIuxeVkU5WdSlJtF+ZhsLplcyHnjCkhLU4BIalEwSNLq6Ozi6Ml2jja3cbS5nYbmdhqa2zja3EZDc7S94UR329u3bZ1dp9YRyclg7uQi5k0t4tIp0R99QZEku4GcxyAyomWkpzGuIJtx/fggd3dOtneyv+Eka/ceZe3eo7y+5ygP/mo7nV3Rf5amFuedColLpxZx4cQxutS4pARtMYjEaG7rYP3+Y7y+p+FUYBxobAEgKz2NCyaOYd6Ut7csphbn6RiGjFralSRyjg42trB2bwOv7z3K2j1HWbevkZPtnQAU52cxb0oR86eNZd7UIuZOLiI/WxviMjpoV5LIORpfmMOiwgksumgCED2msbW2Kdj91MCaPQ2s2FwHQJrB7PFjmD+1iPlTxzJ/2liml2irQkYXbTGIDIKjzW28HhyneH1PA2v3HOV4a3S01Ni8TOZNHXsqLOZO0VaFjAzaYhAZQkV5Wbx3dhnvnV0GQGeXU1PXxJo9DazZHd0NtbLXVkV011MhE4tyGT8mh7IxOYzJ0XkXEj5tMYgMk8bmdl7f28CaOFsV3XIz0ykfk035mBzKx+QwvjCHskg24wuj98sjOZSNySYnU6OjZOC0xSASssK8TK6aXcZVwVZFV5ezt6GZg40t1B5vpbaxhdpjLRw81kLdsVbW7j3KwQ0ttHV0nbauorxMxgfhMbEol8ljc5lUlMuk4LZ8TA7pOmFPzpGCQSQkaWnGtJJ8ppXk99nH3Wk82U7tsVYOHosGR21jC7XHWzjY2MrBYyd5c38jR0609VguI82YUJQTDYuiPCaNzWVyTHBMKMrRORnSJwWDyAhmZhTlRa8BNXt8pM9+zW0dvHX0JPsaTrL/6En2x9y+vP0QB4+1ELvX2AxKC7KZUpzHeePymVlWwMzSAmaVFTBlbC4ZutR5SlMwiCSBvKwMZpVFmFUWPzzaO7s42NjSKzia2XOkmV9vrefJ1ftO9c1MN6aX5DMrCIuZZfnR29ICjaZKEfoti6SAzPQ0phTnMaU4L+78xpPt7KhvYnv9CWrqmthe38SW2uM8t7H21CVCACYU5pzasphZms+McQUU50cvTliUl0luZrpGVSUBBYOIUJgbPddi3tSxPdrbOrrYc+QENXUn2F7fxPa6Jmrqm3iyei8n2jpPW09WehqFeZkU5UaDojA3CI3u+3lZp6aLcrMYk5tBbmY62Znp5GSmkZWepmAZARQMItKnrIy0uLuo3J3aY63sOnzi1FVto1e4bafxZHC/uZ39R0+y8a1Gjp5spzlOkPRmBjkZ0ZDIDm5zuoMjIzodOy8/O4NJRbnBQfw8phbnaSjvIFAwiEi/mRnjC6PnWSSqtaOTxpPtNPYIkXZOtnfS2t5Ja0cXLe2dwU8wHdPW2t7F0ea26LyOaFtTS8dpWy7lY7KZVpzP1JI8phXnRW9L8plWnEdRXqa2SBKgYBCRYZGdkU5ZJJ2ySOJhcjbuTkNzO7sPn2DPkWZ2H47+7Dlygt9srafueGuP/pGcDKaV5J0KjklFuUD0+lcdXU5bZxcdnU5HZxftXcFtp9MetLd3BbdBOzjjCqInIE4szGV8YQ4TgsCM5GQO2vMcbgoGERm1zIzi/KzoVW57HR8BONnWGQRGTHAcaWbDW408u+EgHV19X/khM93ISEsjM93ITE8jI7iflZFGRpqdGtK7dm8jh5paT1u+IDvj7aAY0x0YuaeCY0JhDoW5I3MLRsEgIkkrNyud2eMjcc8B6ejs4lBTG2ac+uDPTOsOAOvXB3ZbR9eps9YPNLZwsPFkcBu9v632EHXHW+idQ1np0eBJM8MsetJjmhlpFg29NCM6j+B+GsH8oL8Z/3jDxSyYUTzAV6onBYOIpKSM9LR+HSM5k6yMMw8HhmgQ1Te19giMuuMtdHY6XQ5d7ri/Pd3lBPffbiNmXrQ/5GcP/sF2BYOIyDDISE9jQmEuEwpzwy7lrHTeu4iI9KBgEBGRHhQMIiLSg4JBRER6UDCIiEgPCgYREelBwSAiIj0oGEREpAdz7/taIaOFmdUDu89x8XHAoUEsZ7CpvoFRfQOj+gZuJNc4zd1LezcmRTAMhJlVu3tV2HX0RfUNjOobGNU3cKOhxt60K0lERHpQMIiISA8KBngo7ALOQvUNjOobGNU3cKOhxh5S/hiDiIj0pC0GERHpQcEgIiI9pEwwmNkiM9tiZjVmdmec+dlm9ngwf5WZTR/G2qaY2YtmttHMNpjZ7XH6XGVmjWa2Nvj56nDVFzz+LjN7M3js6jjzzcy+Fbx+68xs/jDWNjvmdVlrZsfM7I5efYb19TOzh82szszWx7QVm9nzZrYtuD39S4qj/ZYGfbaZ2dJhrO+fzGxz8Pv7bzMr6mPZM74XhrC+r5nZ/pjf4XV9LHvGv/UhrO/xmNp2mdnaPpYd8tdvwDz4Orlk/gHSge3AeUAW8AYwp1efzwLfCaaXAI8PY30TgPnBdATYGqe+q4CnQ3wNdwHjzjD/OuAXgAGXAatC/F0fJHriTmivH/BuYD6wPqbtXuDOYPpO4BtxlisGdgS3Y4PpscNU3zVARjD9jXj1JfJeGML6vgZ8OYHf/xn/1oeqvl7zvwl8NazXb6A/qbLFsACocfcd7t4GLAMW9+qzGHg0mH4KWGj9+TbwAXD3A+6+Jpg+DmwCJg3HYw+ixcAPPeoVoMjMJoRQx0Jgu7uf65nwg8LdfwMc6dUc+x57FPhQnEXfDzzv7kfcvQF4Hlg0HPW5+3Pu3hHcfQWYPNiPm6g+Xr9EJPK3PmBnqi/43PgI8NhgP+5wSZVgmATsjbm/j9M/eE/1Cf44GoGSYakuRrALax6wKs7sd5nZG2b2CzO7cHgrw4HnzGy1md0cZ34ir/FwWELff5Bhvn4A5e5+IJg+CJTH6TNSXsfPEN0CjOds74WhdFuwq+vhPnbFjYTX70qg1t239TE/zNcvIakSDKOCmRUAPwHucPdjvWavIbp7ZC7wb8D/DHN5V7j7fOBa4FYze/cwP/5ZmVkW8EHgyTizw379evDoPoUROVbczL4CdAA/7qNLWO+FbwMzgUuBA0R314xEH+PMWwsj/m8pVYJhPzAl5v7koC1uHzPLAAqBw8NSXfQxM4mGwo/d/ae957v7MXdvCqafATLNbNxw1efu+4PbOuC/iW6yx0rkNR5q1wJr3L2294ywX79AbffuteC2Lk6fUF9HM/tz4I+BTwThdZoE3gtDwt1r3b3T3buA7/bxuGG/fhnAnwCP99UnrNevP1IlGF4DKsxsRvBf5RJgea8+y4HuESA3Aiv7+sMYbME+ye8Dm9z9vj76jO8+5mFmC4j+7oYluMws38wi3dNED1Ku79VtOfCpYHTSZUBjzG6T4dLnf2phvn4xYt9jS4H/jdPnWeAaMxsb7Cq5Jmgbcma2CPgb4IPu3txHn0TeC0NVX+wxqxv6eNxE/taH0vuAze6+L97MMF+/fgn76Pdw/RAdNbOV6IiFrwRtdxP9IwDIIboLogZ4FThvGGu7guhuhXXA2uDnOuAW4Jagz23ABqKjLF4BLh/G+s4LHveNoIbu1y+2PgMeCF7fN4GqYf795hP9oC+MaQvt9SMaUAeAdqL7uW8iesxqBbANeAEoDvpWAd+LWfYzwfuwBvj0MNZXQ3T/fPd7sHuU3kTgmTO9F4apvh8F7611RD/sJ/SuL7h/2t/6cNQXtD/S/Z6L6Tvsr99Af3RJDBER6SFVdiWJiEiCFAwiItKDgkFERHpQMIiISA8KBhER6UHBICIiPSgYRESkh/8PnGEL6S+RhiMAAAAASUVORK5CYII=",
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "needs_background": "light"
          },
          "output_type": "display_data"
        }
      ],
      "source": [
        "main(20)\n"
      ]
    }
  ],
  "metadata": {
    "colab": {
      "collapsed_sections": [],
      "provenance": []
    },
    "gpuClass": "standard",
    "kernelspec": {
      "display_name": "Python 3.10.0 64-bit",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.10.0"
    },
    "vscode": {
      "interpreter": {
        "hash": "b38687464f26cc788bd6aab5d6ae24f3673c8f039be5f165e71e31106dd7d5ca"
      }
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}