Poisson_random_process.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "provenance": [],
      "collapsed_sections": []
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "source": [
        "# Modeling of the Poisson process"
      ],
      "metadata": {
        "id": "-264yNOnP05v"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "To simulate the realisation of a Poisson process of length N with a given distribution law F(x) of intervals between jumps.\n",
        "The expectation M = 1 and the value of the coefficient of variation c = 2 are given for the distribution.\n",
        "Task:\n",
        "- Set the number of realisations K and their length T. According to the given M and c, select the parameters of the distribution F(x), simulate and present on the graph K realisations of a Poisson process of length T with a given distribution F(x).\n",
        "- Check the law of large numbers at t = T/2 and t = T\n",
        "- Check the central limit theorem for t = T/2 and t = T\n",
        "- According to the given realisations of the process, select the distribution function and test the hypothesis of agreement."
      ],
      "metadata": {
        "id": "pnCuJIWXPH_c"
      }
    },
    {
      "cell_type": "code",
      "execution_count": 1,
      "metadata": {
        "id": "Zwbl9Mo6GCUF"
      },
      "outputs": [],
      "source": [
        "import matplotlib.pyplot as plt\n",
        "import numpy as np\n",
        "from scipy import stats\n",
        "\n",
        "from scipy.stats import erlang\n",
        "import pandas as pd\n",
        "import statistics\n",
        "import math"
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "M = 1\n",
        "c = 2\n",
        "sigma = 0.5 #c = sigma/M\n",
        "lambd = 4\n",
        "k = 4"
      ],
      "metadata": {
        "id": "ixUi2zbGo1jh"
      },
      "execution_count": 3,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "def generate(m, a, c, x_start, number):\n",
        "    gen_list = [x_start]\n",
        "    x_prev = x_start\n",
        "    for elt in range(number-1):\n",
        "        x_prev = (a*x_prev+c)%m\n",
        "        gen_list.append(x_prev)\n",
        "    return gen_list"
      ],
      "metadata": {
        "id": "lUn_eMqkphEY"
      },
      "execution_count": 4,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "def get_recov_processes(K, T):\n",
        "  K_big = 20000\n",
        "  k = 4\n",
        "  sum_array_0 = [[] for i in range(K)]\n",
        "  nums_of_intervals = [[] for i in range(K)]\n",
        "  for num_proc in range(K):\n",
        "    u1 = list(np.array(generate(16807, 0x7fffffff, 0, num_proc+1, K_big*4))/16807)         \n",
        "    E = [0]*K_big*4\n",
        "    for i in range(K_big):\n",
        "        for j in range(k):\n",
        "          E[i] = E[i] - 1/lambd*math.log(u1[j*K_big + i])\n",
        "\n",
        "    S = 0\n",
        "    for i in range(K_big):\n",
        "      if S<T:\n",
        "        # print(S,T)\n",
        "        S += E[i]\n",
        "        sum_array_0[num_proc].append(S)\n",
        "      else: \n",
        "        nums_of_intervals[num_proc].append(i)\n",
        "        break\n",
        "  return E, sum_array_0, nums_of_intervals"
      ],
      "metadata": {
        "id": "mN1mPKL_cs9L"
      },
      "execution_count": 10,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "sum_array_0 = []\n",
        "S = 0\n",
        "nums_of_intervals =  0\n",
        "for i in range(K_big):\n",
        "  if S<T:\n",
        "    # print(S,T)\n",
        "    S += E[i]\n",
        "    sum_array_0.append(S)\n",
        "    nums_of_intervals +=1\n",
        "  else: break"
      ],
      "metadata": {
        "id": "jhnb1vZrU4Vr"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "K = 50\n",
        "T = 50"
      ],
      "metadata": {
        "id": "VkX3j8NUIMWg"
      },
      "execution_count": 2,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "E_0, sum_array_0, nums_of_intervals_0 = get_recov_processes(K, T)"
      ],
      "metadata": {
        "id": "eaDUIvVVSSAN"
      },
      "execution_count": 11,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "for i in range(K):\n",
        "  x = np.array(sum_array_0[i])\n",
        "  y = np.array([i for i in range(nums_of_intervals_0[i][0])])\n",
        "    \n",
        "  plt.step(x, y)\n",
        "  plt.xlabel('t');\n",
        "  plt.ylabel('N(t)');\n",
        "  plt.title('Poisson processes for Erlang-distributed intervals');\n",
        "plt.show()\n"
      ],
      "metadata": {
        "id": "Fo17nrhksyq_",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 295
        },
        "outputId": "8aa1d117-20f9-49c3-83e4-5b884073e30e"
      },
      "execution_count": 12,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEWCAYAAABhffzLAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOy9eZwcVbn//36qeu+eJbP2TDKThWxkA0KAGwVkEQQREdSL4hUIXgXvT/SlXBX9iixuXBX3DVyiqAiCIiiLohJjWDNZyEZC9nX2ySzd01tVnd8fVTPpTGayzpbkvF+vfk13Lec8darmqVOfes5zRCmFRqPRaE4ejJE2QKPRaDTDi3b8Go1Gc5KhHb9Go9GcZGjHr9FoNCcZ2vFrNBrNSYZ2/BqNRnOSoR3/YSIia0XkgpG2QwMicrWI7BSRhIicMQL1XyAiu4a73v4QkbtE5Dfe91qvTcxBKvsnInKH931Ej1lEfikiXx5g3QdE5G/DbdPRcLDjGE5OOscvIttEJOX9gzR6JyJ2qP2UUjOVUouGwUTNofkm8DGlVEwpteJYCxORRSKS9q6Jns+fB8HOYUUptcNrE/tg24nIjSKy5DDKu0Up9aXBsE1ElIhMHoyy+qKU+q1S6tLDtOOwjv1E56Rz/B5XKqViwFxgHvCFEbZnyBAR30jbMASMB9YezY4H6Q333Eh6PlcOsP+J2J4HMFhPDScaJ8r5P1kdPwBKqd3AM8AsABF5pyfptHu9wFN7tvWeFN7qfT9bROpEpNN7aviWtzwkIr8RkVavjKUiUumtqxaRJ0WkTUQ2iciH88q+S0R+LyIPikiXZ8O8gez2ek8fF5EtItIiIt8QEcNbd6OIvCAi3xaRVuAuESnyym4Wke0i8oWe7b19Piwir3t1rxORuXk2/8Hbb6uIfDxvn6NpgyIR+bmI1IvIbhH5co+DEZHJIvIvEenwjumRfo47KCIJwAReE5HN3vJTvfPV7rXdO/P2+aWI/FhEnhaRJHDh4V0dvftfICK7ROSzItIALOxnm9tFZHNe+12dt+5GEVkiIt8Ukb1eO16et36iiCz29v27iPxQPOlmAHsmeu3UJSLPAWV56yZ414Yvr+4t3rZbxZVETgV+AswX98mmfaB2kn5kCRH5vHd+tonIB/KWLxKR/+573N73xd7i17w6r/WWv0NEVnrn7UURmZO3/xkistyz/REgdJA22a8X77XBLSKy0Sv7h+Iy0LEHvfOzw7uWfyIiYW/dAedf3P+Vd+TV5xP3f6Tn/+ZREWnwruXFIjJzALvLROQvno1tIvJvyfu/HFKUUifVB9gGvNX7XoPbc/wSMBVIApcAfuAzwCYg0M9+LwEf9L7HgP/wvt8M/BmI4DqnM4FCb91i4Ee4F/DpQDNwkbfuLiANvN3b72vAywc5BgU8D5QAtcAbwH97624ELOBWwAeEgQeBJ4ACYIK3/Ye87d8L7AbOAgSYjNujNoBlwBeBADAJ2AK87Rja4HHgfiAKVACvAjd7634H/D+v3hBw7iGOf7L33e+dp897dl4EdAHTvPW/BDqAN/eU3U95i3rar591F3jt+X9A0GvPC4Bdedu8F6j2yr8W9zqqyjsfOeDDXnt8FNgDSF47ftOz/VygE/jNQY79JeBbni3ne8f6G2/dBK9tfF4bd+a1QxUwM8+mJX3KPaCdvGVf7tMOPXW/xTvOaf21Yd868s+Z9/sMoAk4x2uXG3D/x4JeW2wHPumd3/d4bfjlAdqkv7r+AhTj/n80A5cd5Ni/DTyJ+/9UgHv9fu0g5/+LwG/z9r8CeD3v901eOUHgO8DKPu3c06Zfw70R+b3PeT3XxZD7weGoZDR9vIsrAbR7F9ePvJN5B/D7vO0MXId4Qd5+PY5/MXA3UNan7JuAF4E5fZbXADZQkLfsa8Avve93AX/PWzcDSB3kGFTPhez9/h/gH3kX9o68dSaQBWbkLbsZWOR9/yvwiX7qOCe/HG/Z54CFR9kGlUAGCOctez/wvPf9QeABYNxhnMN8x38e0AAYeet/B9zlff8l8OAhylsEdHvXRM/nS966C7z2C+VtfwF5jr+f8lYCV+Wdj0156yKe/XFcp2QBkbz1v2EAx5+3fTRv2UMM7PjbgXfnt3meTf05/gf7WdbX8efX/Xvgjrw2PBLH/+OeNs5btgH3hnI+eTdHb92LHJnjP7ePnbcPsK3g3sBOyVs2H9h6kPM/GfeGG/F+/xb44gC2FXv2FPXTpvfgdsgm97fvUH5OVqnnXUqpYqXUeKXU/yilUrg9tu09GyilHGAnMLaf/T+E+4SwXlwpo+ex79e4jvRhEdkjIl8XEb9XdptSqiuvjO19ym7I+94NhOTgeuLOPmVVD7CuDLc3sb3P9j111wCb+yl/PFDtPYa2e4/Fn8d14HDkbTDes6M+r7z7cXv+4D5hCfCqJ9fcdJBjz6ca2Omdr/6OD/Zvj4H4uHdN9HzuyFvXrJRKD7SjiFyfJ1m040qHZXmb9J5bpVS39zXGvuuiO2/bXls9yaHnZfPnve33KqWSfY71ALxtrgVuwW3zp0Rk+sEagEO3U391Vw+08SEYD9zW5/qq8cqrBnYrzzvm1XUk9P1/GiiAoxz3Zrwsz45nveU97Hf+lVKbgNeBK0UkArwT9waMiJgicq+40l8nbocR9r8eevgG7tPq3zxJ7vYjPMaj5mR1/P2xB/diBEBEBPdC3N13Q6XURqXU+3Gd1v8Bj4lIVCmVU0rdrZSaAbwJeAdwvVd2iYgU5BVT21/ZR0BNn7L25JuY970F9zF5fJ/te+reCZzST/k7cXs9+c6wQCn1djiqNtiJ2+MvyyuvUCk10yuvQSn1YaVUNe4TyY/k8KJA9gA1fbTRvm2rODYG3F9ExgM/BT4GlCqlioE1uDexQ1GPe11E8pb1nlflRtX0vGz+qrf9GBGJ5m1fO6DRSv1VKXUJrsyz3rPzYMdzqHbqr+6e6y6J60B7iB+irJ3AV/pcXxGl1O9wj3Os9z+YX9dg0PcYW4AUrgzWY0eRcoM/BtoH3KfK9wNXAeu8mwHAdd6ytwJFuE9h0M/1oJTqUkrdppSahHvz+JSIXHyUx3VEaMe/j98DV4jIxV4P9TZcR/Vi3w1F5L9EpNzrZbZ7ix0RuVBEZov7wrIT1+E6SqmdXjlfE/fl5xzcHvOAL/EOg0+LyBgRqQE+ARzwMhRAuaF9vwe+IiIFnqP6VF7dPwP+V0TO9F6ATfa2eRXo8l5qhb2ezCwROeso26Ae+Btwn4gUioghIqeIyFu88t4rIuO8cvbi/rPl9+IH4hXcHt1nRMQv7liLK4GHD68Zj5korq3NACKyAC9Y4FAopbYDdbgv4AMiMh/X9kNtf7e3/bkDbS8ilSJyleeoM7jyZk97NgLjRCRwOHb2oafu83Bv6o96y1cC14hIxLthf6jPfo2474l6+Clwi4ic4113URG5wuscvYQrK33cO6fXAGcfha39sd+xe9fvT4Fvi0gFgIiMFZG3HaKch4FLcd/ZPJS3vAC3vVtxb4RfHagAcV9uT/ZucB24cvDhXPPHjHb8HkqpDcB/Ad/H7QVciRv2me1n88uAteJGmHwXeJ8nF8WBx3Ad3uvAv3ClD3B7BxNwe0iPA3cqpf5+DCY/gfvydSXwFPDzg2x7K26PbAuwBPdC/QWAUupR4Cvesi7gT0CJd8N4B+6L6K24bfIz3F7M0bbB9bgv7tbhOvfHcHuj4L5cfsUr70nc9w5bDtUI3vm5Erjcs/FHwPVKqfWH2rcPP5D94/iXHc5OSql1wH24zqoRmA28cAT1fgBXU24Fvox7A88cZPvrcN+/tAF34r4b6Q8D9wa/x9v2LbhOCuCfuEENDSLScgS2NuCetz24uvYtee38bVwtvBH4lbc+n7uAX3lyyn8qpepwX3j/wCtzE67+3nNOr/F+t+FKVn88AjsPRn/H/lmv/pc9eebvwLSDFeJ1ZF7CfarN73Q9iCtL7ca9zl8+SDFTvLoSXlk/Uko9f6QHdDT0RBZojiNERAFT8h4vNScI4oYurldK3TnStmhOXHSPX6MZQUTkLE/yMkTkMlx9+E8jbZfmxOaEGIWm0RzHxHFljFJgF/BRNQhpKDSag6GlHo1GoznJ0FKPRqPRnGQcF1JPWVmZmjBhwkibodFoNMcVy5Yta1FKlfddflw4/gkTJlBXVzfSZmg0Gs1xhYj0O+JZSz0ajUZzkqEdv0aj0ZxkaMev0Wg0Jxna8Ws0Gs1Jhnb8Go1Gc5IxpI5fRIpF5DERWS/udGXzRaRERJ4Td1q050RkzFDaoNFoNJr9Geoe/3eBZ5VS04HTcLM13o47W9QU4B/eb41Go9EME0MWxy8iRbhTqN0IvalWsyJyFe50ZuCmb12EmxZVo9FoTjgefeNRnt7y9H7LppdM57Nnj5zbG8oe/0TcySkWisgKEfmZNylEpZfLGtz83pX97SwiHxGROhGpa25uHkIzNRqNZuh4esvTbGjbMNJm7MdQjtz1AXOBW5VSr4jId+kj6yillJdb/gCUUg/gTr7NvHnzdCY5jUZz3DKtZBoLL1s40mb0MpQ9/l3ALqXUK97vx3BvBI0iUgXg/W0aQhs0Go1G04ch6/ErpRpEZKeITPOmNbwYdyqydcANwL3e3yeGygaNRqMZCfJ1/Q1tG5hWctCZHIedoU7SdivwW29i4y3AAtynjN+LyIdw56b8zyG2QaPRaIaVHl1/Wsk0ppVM4+2T3j7SJu3HkDp+pdRKYF4/qy4eyno1Go1mpBltun4+x0VaZo1Goxmt9BeuORrlnXx0ygaNRqM5BvoL1xyN8k4+usev0Wg0x8holnX6Qzt+jUajOUyOR1mnP7TUo9FoNIfJ8Sjr9Ifu8Ws0Gs0RcLzJOv2hHb9Go9EchNE+GOto0FKPRqPRHIR8eed4lHX6Q/f4NRqN5hCcCPJOPrrHr9FoNCcZusev0WhOKvoLyTwYR6vrP/TKDp5YufuI9plRXcidV8484rqOFN3j12g0JxVHOjHK0er6T6zczbr6ziPebzjQPX6NRnPSMVya/YyqQh65ef6Q13OkaMev0WhOOA4m5wxlSGa+vLOuvpMZVYVDUs+xoqUejUZzwnEwOWcoQzLz5Z0ZVYVcdfrYIannWNE9fo1Gc0IyUiGYo1XeyUc7fo1Gc9zTV9oZjhG2/UXtjGZ5Jx8t9Wg0muOevtLOcIyw7S9qZzTLO/noHr9GozkhGAlp53iQdfpD9/g1Go3mJEP3+DUazXHFSE6GcryEax4K3ePXaDTHFSM5GcrxEq55KHSPX6PRHHeMZLbM41XXz0c7fo1GM+oYqZG3fekbsnkk8s7af+/mjVcbe3/XBZ9nVehlAPxBk6Ly8CHLmF4ync+e/dkjtPrQDKnUIyLbRGS1iKwUkTpvWYmIPCciG72/Y4bSBo1Gc/wxUiNv+9I3ZPNI5J03Xm2kZVei9/eq0Ms0mDsG3cajYTh6/BcqpVryft8O/EMpda+I3O79HvxbmkajOa4ZLZOfHIu0UzYuxtW3zQXgyWdjlDFjVBzTSEg9VwEXeN9/BSxCO36N5qRkIElnJOWcfA4l7fSVc2CfpJMOd6PMHA/80gSgmWbKKWfhwv0dfzwe5/LLLz/GozgyhjqqRwF/E5FlIvIRb1mlUqre+94AVPa3o4h8RETqRKSuubl5iM3UaDQjwUCSzkjKOfkcStrpK+fAPklHmTkc7N7l5ZQzXU0fHKOPkaHu8Z+rlNotIhXAcyKyPn+lUkqJiOpvR6XUA8ADAPPmzet3G41Gc/wzGiSdwZJzYJ+kc379+QAsuHHBoNg4mAyp41dK7fb+NonI48DZQKOIVCml6kWkCmgaShs0Gs3IM9KSzrHIOT30J+u07EpQNi7We3xdXV1sS24jaFazqqSWQCDA7/75KqlUPYY/g0KhFIzL7eA9iUdwzABmcTnjS6MH1FcQO5WpU+84ugM+BEMm9YhIVEQKer4DlwJrgCeBG7zNbgCeGCobNBrN6GCkJZ1jkXN66E/WKRsXY+rZlb3Hl0wmKcoWkYnOpzVWRDQao7ujHTG7EdNCKVBq5AWMoezxVwKPi0hPPQ8ppZ4VkaXA70XkQ8B24D+H0AaNRjNKGGlJZzAGXvWVdXp51j2+8+vPBz88OfFcaoDHz5jCI3ffTtHsRZRPmMTXl94KwO9ufh/w6WOy5VgYMsevlNoCnNbP8lbg4qGqV6PRaDQHR4/c1Wg0g8JoHW0LR55Qra6ujleWLKO7K9u7LJex8QdNPrnw16zPi1PpNP10qgZMfxWbxkURhJ1tzdRY2/nTo9fin+mjJlFP/Pkt/FD9E4DNPwyxZlwFAiCCowqoqakAhlbb70EnadNoNIPCaB1tC0eeUG316tW0tjeTy+wLx/QHTSIFAdbLeprZF2KeMH1Y/lqCkRmYpgVKUZPbzjmZF3AMHxgmlY05CruzGIZ17Ac4COgev0ajGTRGWsfvYTD0/LBZxOSCNx+g6dc9+zJx4r3HefWKjQB8Qd0JbGfTn2tpbF3NmEuho/Oj5OqTFPIobdH1/GPOWK5763LKgFOOybpjQzt+jUZz2IwWOQcGDtE8mkRqrdlt7LX2lZWyO/BZMerGPM+jT9xBMl2P5Vg4StFkQchXzNd+8z9c1vQqn5UQQZVmPNuwbR+lVhFSDMHVaUzfcjoDRWRS7dRTw0+s2/neQ0t666lqS3HFK3t6f1sqxXMzpuIrDTOjupA7r5x5DC00MFrq0Wg0h81okXNg4BDNo0mkttfaTcru6F0eNouoiNaysXwZmzp3gZPBVjYKKDMhHZ3Phc2rmJzcRdBJU+h0Yts+ctkwIiAG+HxZTNOiyyhkbXQyfy+dj50ZHRO36B6/RqM5IkaLnAODF6LpK4kBMRYs2H+U7QvP/pnaRJhrirM8Z51G93Y380xkRhmxHUvpSIZ4ubGUs2+ayuLXLqFlV4Kx24s5a/efaP56kNyeBF+We+hq3sb1/+zkf+eUHRgO+rH9f95yTEdzeGjHr9FoDmCkR9r2x9Hmxq+rq2P16tV0d2T3i9KBfZE6630r2VO0h8XPLiaTaSKdasaybHY6Dvir+Z7cxAVNr3FV8jMAqEVCbWITmyvG8tBFN/G7VlBjBRkbwXeOzbolk5h7+xrsjMEHuYPw3L2Ezu8AFM/93K074wuSipVw2uzzhzyKpy9a6tFoNAcw0iNt++Noc+OvXr2ahoYGuruy+0XpwL5InT1Fe2g1WwHI5lpxVAocUP6xJKLn00ExV7U8z4zkZkQUBg4bouN5ruxsAJQSBIWJKwedtWwNtfUNODYYCoKzOjDiaTd8cxSge/wajaZfRpOk08PRSjvxeJzittMgQr8jb/MjdZYtv46WnV2sePX9bLqokJxTQWVnO/PHFLK1vYjG6gwFhR1E7OncuPtmPnPpHLeQhVewZvtevsg9VLcodhdU8ccPfoHvE2VTwadQTTt4szEeFjx1rM1wzOgev0aj0ZxkyGhIGHQo5s2bp+rq6kbaDI3mhKWvpt+j5Y9Ej/9QYZoD9fgPFpoZNosobjuNsnExrCs38/jrvySbc6WdNidKS7oDwx+noPRGTDNLzgmw067h2pZnubr5HwSsHHOSG+iIBlkxp5DAHofyX/iRFGTtQloDRSRMAwMDEaF05lY6zhFaosWEnRRLy85hiTqfzraCA+wu2dvCvLolmCicYIQd8VnMCLVw55xOuPzeY2pLEVmmlJrXd7nu8Ws0mgM0/dGk5fdwuJOi9BeaOcY3dr9Mmps6d2HbSQDaVZRcYDyB8GxEGdh2gA4Kyfl9XNP0d2YlNhHNpOiMBGiq9EFLgOhrBiFykFGY2RRdpo+sGG4oJ9B1jiDVFiEnheHkWOKcz05j/KC31dGiNX6NRgOMLk3/aLX8g4Vm9vIs1IbDfG7SJM6c+xBXr9jIhrYNXPPPTqaVTODq2+Zy7uN/o5hWaroa6CRKZNMcVlzk0FBSw9t2f42m4ltIX24Qe30SLbsSfP/cj1I2Lsb3cfPqbwzdAG1+XvzbqSTNELvfP43aFPz7/VcMYPkNAywfGrTj12hOUvLlndEy6haOPkyzJzTTyiaIx+NQt5BHV/+SpyVJxrTImhZKKbbafpS/mk+3vp3Anx7ijfAEQpkI1oRF1LR/n7pvwDdNH0HSVEgLnaEg69+6Fl+lxbatk3lg0SNM2DuGzggErXoyRUGuDH6H6oIGXnRsTMfmxfJzecE+n9QlQUQBVgpp284TF7sO3u+YFFiu7GMYBpFYjJKSkgOOMXjqdOKf//zRN/QAaKlHozlJyZd3RsuoWzj6MM2e0Mx4PM7s2bNh9WM8nWtiA1mypoUtyp0ExV9NMno+pgXiGMTb6rHF5OKuV5iW2gZAkDRFtNMZCtIYDSPBQoKJWl6pn8f0XSvJBLJk/ELGFyQRiFE9qYHCkgSmYyPK4QX7fHaY4xEFpuMwrnEHZ615aTCb8JjQPX6N5iRmJOWdwRh1e9AwzdVAIMq0+Gw+VpEB4AdNQdL+dzOzZBqPnzGFr/zvIgJND/Pwu95FwMkR9E3lnm0f5Na536Y7FCOW+APlQLn3yiA6NknAMFD+ENWzZnLtnfey/ZIzeD5UQFu6lJpHfMwqKaP73VUY/npmlv+FhZctZOHChXQmOkmf/3auvfPYXtgOBjqqR6M5SRjJyJ2BRt0eyvHnz3PbN1qnrqKIzZW1+HIx/EGTspoYHbsfo7PrBRwUoLByTZj+OEUV78eyDK7c+QJv6VpJ0ElT4HTRZbqSUkaCzExswvILaybMwCquZ8nmS9m9dhyntNYRzbhTLvpPSzJmajumssmZflpCZYDFxNBuFne9leXtbwK/YtO4U5jaleWHdeuwckEylukaLdBpN7Mp9Q9ChjecKxCF0klUjJ/EhTd+ZHAa3ENH9Wg0JzkjGblztKNu8+e57Ruts7mylpZoQa/EA9DV+SLZXAOusC6Y/koCwRlYlkE2bXJpx8vMSmyi0OmkyywkI0EAgipDU6iQluJyDNOPv72KVxrPJpZ8A5/TRdBynxjGTG0nVtqNJT5SvrBniY25x2CpM59NNeNRIYepXVku26MwlEHGMrGUgTfryjG35WCgpR6N5iTieJR2eua5XbjwNfKjdV5csbF3XtseFrxWT1cgzeemFnPm3Id45O7bAfhW7WMIMCm7g6wPPtL5QzrmgpgW3+n6POOenc34Xz/I7rteZEIDjL3rTUTvf4nA9r9TkLOZ0hlh8cWfIRL9HInuAtKr5iORJD5l8udxe/j6b+sJ3GAyZ1oFj1/+Zp6/72UAJnztI7025Es8b+a2o2vEQUL3+DUajeYkQ2v8Gs0JTH8hm0Pd4z/YnLeH2+Nf++/dPLD8h6wJLEUMwR80yWazBAIBpKSU5qxFt+0QMQ1mxcK07NyMozpoCtqM9SsWmBHWrLgY27Eor9nBOcFV1La3U5RJsz1cyTWn/Ij2Dh+nr36V2W2rMOygm0Tfc4cTJm6g/JRWlIDCwFFBorkUi8aex4vOeagmV7LxK4WJQnKwuWYi4wIx7lzZjdqdYlcmR8rsJJdOg98kUby/r+0u2svOOXVML5nOZ8/+7LE3fD9ojV+jOQkZiZDNwZjz9o1XG1kTWMreYCuG6TrZQCBANBrdz+mXB1y12upOYnpO/zTHpKlpIlYggCV+Sqt3MSHRQlE6Q0cwxKOll9MeKyC4J8HkljcoSHe5Tt9DYVNxSiuhsgyGErrNCJFciqCV5UXnPHb4JvRu60NhKoUEhFOzHZy7M0PV3hy7MjkSXiJQfyhEOmjTbaWOskUHH63xazQnOCOh6w9GqKYYQqVRyTM3PbPf8p45bvO1/V/8/jIqLtlBJFzNrMZv8VRsGbUxeDY7nQm8ipGB9F4/zyW+xauXnMWZQKDUpiDkZ+Kk6Zyr3kkuneaV9vu49hx4wfRjN6V5YdM3eMZy+PSL30YKOij5WIiC7t18qHWaGz668Ar3fe2NbsbN5+97mfoxflKdjZgm3Py99wCw4Fn3vcRoGRmtpR6N5gRiuEM2B0PWqaur47erf8t6WY9tOTi2QjmKvcFWKo1KfnbJ9SzcuYN/ZU8laztsVxXUWDv5VOd9AASCaXyBDKZpYSkfexJxSumiWCUARcbyYVtgY5B1/PgctysuKJSA6SgK0zkAbFHYYqKqLZz6IN2/mMDzZ53F2pnTyAUU9WNrmda1l+8ub8MwhJCVIu0Ls4NCAqkkaauS7dkcXZbCMlp47arniTXGkAYh4oscMDo6Ho9z+eWXH+NZGJgRk3pExBSRFSLyF+/3RBF5RUQ2icgjIhIYahs0mpOF4Q7ZHAxZZ/Xq1azKraKZ5l6nL4ZQ5pRzcfXFNDT+mX+mxrPFriBnO4y1d3JWdmmPHN/r9HvIOAGKVYIw7mxbWcsg65hY+PA5NoZyEMBQCtNxCFrujcARAIUtBlZDiNyqQhxRLD9tNptrxkPIYEZiC2+rz+KzqwBI+8K0B0sIpJL4cort2RwJ28AyWthbuQWAaEuU4mwxJeEDUzKMFMMh9XwCeB3oSb7xf8C3lVIPi8hPgA8BPx4GOzSak4LhlnYGQ9YJBALMjs/mnWtvBfYfhbts+XWYZoQ5BaUEXlVUtP2B8k6L3N5yrv3JUzz/3GxQAaxdQn2khL++fB0fcP4f6TFFvFDwW17d9gv8QZN7vvldtn/wegDG//pB3jj3SpSCl2tMSt5Tzz2Bu8kmHkC2v5PS7rF8LLiH8y77Od8um8Tp8RiPX/4mWHgF2xqvZXPZWC687dJeG+/5xNUAzJ50B4XA1be9rXfdwvqFUMLASeNGgCGVekRkHPAr4CvAp4ArgWYgrpSyRGQ+cJdS6m0HKUZLPRpNHw41J+5wRe70lXV6EqcNRHdHlpXOa2yNbupd5iiL9mAHMZmEXXoz/lAXjmnj2A4ApmGxnYnQ7lD16lY+VPsIU6o20BuCA2QtEwptsi0BGh+rwsRBKbAdk66gn3DWYWyTTXWimT2xcn58xn9iWBkcn5/1VQaRbJQ5ySImZf1YCAEFTRMDLJ9o0FAYojyR4L3LlrvyEAYOipzkENdEROHKRkYIf9BkTFEAJ+HKRy25Dsr8Rbyr9PwD2iNQHaX4ylOO9XQMyEhJPWoX7u0AACAASURBVN8BPgN4zUMp0K6U6nku2wX0+0woIh8RkToRqWtubh5iMzWa44uRnhM33+nnyzo9idMGorsry5bwRtoCrb3LDPFRkasmHZ1PfTSKmOlepw+QsQNEUp1E97QDiinxDaiYtV+5Wcsk1RIk+XoUsRRKgfKKKMjkiHQrUoafPbFylldOB8Dx+emKuKkUwrkIp2QClNgmfhQxbFaND9NQGKI00cE075jcRBAK23AQxxsgjOv0DdPsHUXsJHKorCshlfmLmBKuOap2HiqGTOoRkXcATUqpZSJywZHur5R6AHgA3B7/IJun0Rz3jHT+/IEknng8PqCs8fh9y4n4F1NbNe4A269esZGirte5N/ZnHn2+FoB7b7uXC277NcHcXm5f9xRV520nK1lUN6TagzTuvIpr77y3d3RsbN1mon6Y2JLFVDDmg1+i4uY5XHv/S6xNpCge186tz/+cL/74+1zzq9lks0X8d8sVFDtj8UmOokADf3zrGZTVxBizJ4mdbOPmDU+z4HPfOiAyp78RuT003b8KgIqb5xxN0w45Qyb1iMjXgA8CFhDC1fgfB96Glno0miNmJAZjwcEjdz59pr9X2lnFKlblVhEOZwiF9/XIHcvBttzut1JQ72SpDvmYUHgeKwzvxqFgp7+G8Wzje9u/QLTeR1cuiqDY1lZEriJFdKKbs8dXmiPXFmDJ02dRlElSkO0mExCCWcU5m/ewO1bOT2dfBb4gmTFjiTkwNeEDI4cZasVv5whbBfitKKIEnzIQZfLK9Axra8MoQ+FXQmO0kLJkBzeteYzNhe0ENuwl4ncjcyqtGoqTpfhDIcrHTzqgzXL1CfxVsRF3/MMu9SilPqeUGqeUmgC8D/inUuoDwPPAe7zNbgCeGCobNJoTiZHKn3+wyJ18aWe9rKcj0IHfn+6d1hDAtly9HUAEaoIhpoYclhnnsNNX0yuX1GZ3MDe5nOq2bipznWQdN+CvpStAYEIKX4mrmVutfpq2lJEx/BRkuwnYOYJZRUHSZnesnBWV0xEziBUqoDVoMK3bT7ktiNmFoRxACOai+O0AgqBEYRtZ1tZGaCko8G4EBlWJTk5v2opT3EHxlgwlXQFKQm5kTrk9luJABZHC4n7bzF8VI3J6+WCdgkFnJAZwfRZ4WES+DKwAfj4CNmg0xyUjJe8MJOssXLdP2ln87GLixHtz35859yHAlXdg/0idBc8uoLszQpXTyd2vj+XHVT6+vuTHPFNdhVT42EQFN9i385t3PYbv692YjiLSYPBmJsGCp5h13z/hNHhq5S8BGH9RK4+vuw47PpuPxnxk1r/OtfOChH3dnF3XxZpwhrWn/ZirMx8HhOK2Kja0reeFiY9zU8fnePCiAqzdexhrZXiiqQwaVlNR+xDc/BTwYVfWKdkn64x2KedQDIvjV0otAhZ537cAZw9HvRqNRqM5ED1yV6MZxQynrj/QPLj9jcStq6tjy9Zf4fcvJxAIEI9X8bUt7oClj5YkSXVWsnnR/6AsB59VgG12kSioY824U5gZfIPL2pdgmT5MJ4c/kQbczAf1RYKqTONgkFIBBDCVTbg0TUNHnIUbbuL0bRuZv/E1DEMxvW0XdlEV2970aVozDmATcvykjRybw81sGDeGHZXlIDYFppBTDhnJ7js4AVECAvWxKNVdCR571UBIkTVa6fYmT8+l0/vp+aNFwz8UOkmbRnMcMpy6/kDz4PY3Enf16tX4/cuJxfYSjcb2W5dOl9G25UyMTBSUwja7yAbqAXijajwXd73MtNQ2gqSRnN2r5QPYlVnMMRaWMnG8SUtsMWnbO4ZN9VNo8Y9h/tbXOKVjN2GlSMTGYdacTTIlZHICOKSNHO2+JKIUW6qKaSw2CYoQVH4yksUSN8xSyb5wTICqrgQXNSQQI4vhT9Nt+93MmriJ1vL1/NGu4R8KnaRNoxnlDKeufySjcAOBACUlp/dq+QVNbrhjU/35pOsrsWkkkvoXU2tv5MdVY3nk5lupW7GRgj0+Cmrm0h3Yw92bFrA7XcPvfa6TN0P3Yrd38osXPoQEC/nTt72w0IVXQHEDay6YQWxpmIKkn43T/4u2cYpzYyapnbA7sptXJv6ci9svJx6Pc9uCBSx+cgWxNPyhNUDLri6eCqyhsdhk87mX8PgZU/bX6hdeAWHgq27CtcUHCdc83tFSj0YzihiuJGsDhWi+qaiD/yhODrCXSzS2vLe3X+Q7lWefu4a/RFeQGvtXookJXPL6/8ek6n8z27eYHI4bygP4HJvdbWG6m2ycHHQmIzjKoHjmXsZM68AyfIRKu8k2B8n8pIaQk6MwGmJvwOyt+5+nz2PltDPJBQrB506bKIClwBEHy8jic3xutI5h01QYoDKRZMHKlVgKWs0EjcXFvCMc46wGP+lUBV1WE2sSD0E26c5/WzUbgOZtWymfMPG4dvxa6tFojgOGK8naQCGak8y2g468BYhG1xCL7cXvP4XC+vn8M2WTLHZtntI8lxLLx2zfYqZkNpM/x6zPsWlt9tGcjpDJ+UG5GTKLp3UQKM9iiUGmOYhvuY+Qk6NE2ewNmKTNfW5qxfRZbI9XYhuB3owNCkAUjtiYykDcFGw4Ro6yZCdzGhvJASkDcuEisuMncVq9SS5VSpfVxJbOtW5BgSjE9sk35RMmcuqbLzi6Bh7laKlHoxllDJe005+ss3DheigYeOQtwLLlzwFVnDn3IVcqkdfwYVOcK2F64kK6Y+APZFkvE/lG8i7+eNdbAdhw1tlIrU24UjFnTZotF8/j2jvvZdnTs8ilavnx2o/z4SU/YMb4Ksa//CAAU552h/n8/Xd/AKChKkdJVzfvf9EiURqgoP6npIrHUFkkLPjct4B94aPtJVvJ1SeZEpzH966eu1/+/qZVv4VwK2tYDQVw7Z1PDU6jHidoqUejGQGGK8naQJE6qfqNTAvsZUZV4X7LGxoaKI+W8M7gfB6J/oDnZRm2sy8dmlJg+DLYVhBfZ5igrwvHFHb4bGotk9vbChAccg2drOmsQtr35dwhq2iLRCicnaPM88E+xyZUvJtsW5w9z95AyMqwZFYVr06Pu/WJIErh92azaioyqWy3+eCiLhqLHWqtl8kE/d627pOBiWBg0GYkKXUKeFtmJl12I6bse/qIGuUknWaeb33yuJdzDoaWejSaUcRwJVkbKFJnWmAvBc6BWn48HucUO06uPsEiWcFOO+M6fTc7GQqwrSC5TAFRo5OgZHBEqLVM/iMVQHAwlM2azio6u/37Fx4Q2oNjKDslSaxoFwCGY2O1VpDcMguAVCDE8illtBW4YoQohZF376josJm+K0Nn2KE73Egm6KfTzO5XjYEr95Q4USZapWRVMk9wckk6zTTLjhNazjkYWurRaEaIEZd0KOxX0mm6fxVEQIImFTmDN625DTtTzaopQbZt3QrAK/deR/qzk9gwppT/fdMdPHf+O2m6fxUbKz6B07WH7g21+EJw0++eYtnTrlM/8+1ruPb+lzjTvAunw6Rp+Wlcvqee1fVr+cU129n7jksJBAIo5VCe6GbZu9/M1c/8EYDr100A9h/9C/C9O75BgRPi1i9/ep/9d/8WgIo7PzAo7Xcionv8Go1Gc5KhNX6NZggZ7glT+mr6A4VoNjQ0EI/HmVRazotbl7A0vAYAy9eF+BIogXqVIU6AS1/8KqHxzzPL9wIKwUARcWzS9V28lhjH1uhYJlp7AAPHl0YyDoGqHMVTkih/EYZtYNr76g6UNpFtq+C51eeydtIMMoEg26pqsUwDn60wFIxJWpT7trMtXMVXXktyWlsExB3FK7bNG2Y928xm9hpJxjgR3u2r7C0/lyrFH27VPX60xq/RjAjDPWFKX01/oBDNeDzO7NmzMdZ2szKwnl1+b70/gfgyKBSVBJm8dzK22cVp5ovMSG/BQOEDlGOzurOazpSfIifhleogWQNfRlE8JUmkxB31atrKy4rpkm2rIL2hlGXTT2fzuPF0RosA8NmKcNZhTNKiJNWBA0xI1XNGW5CINxxAbBtRim1mM21GgjFOhEn2/iOH/eFWIqeGBrVdTzS0xq/RDDEjPQfuwUI0//3cFnKmxcyq2Sy8bCF/feJKAH7R8HnW1XdyZv0TZKtW44hifegUynO3IEETc93TSLqBWGk3c1/by1kLitj+0B6Iz6blkxbbt0Cyqwh2fAvnT3eS9gn3febT3LDoL/yrajGh+TFsVU5xt8k7ljrMrZB9+v3CK1iT6OaOGQ/w+BlT2L38RQCm3vWm3slPQrXTqKaYBTzm7rPg5ArHPFa01KPRDBL9yTpDmVjtYBOkPHLz/N75b//V8S+2xrbiM0xiuXDvtu0BHwmfn6yzAzFqiRZ9hFCsCTsXJpkoIkWAq9r+xlXtizg1vYV/t0+meUuE6Kw2wlPTOALgw5SiA2x7oXwWLzhvwdljErIypP1Bto8dT2mXTU9waFuBj/EdNvfVJWkLrGar4clRjsJWgojPi9AxUDh02bvIWhDwgR2MEA90u44/Pls7/gHQUo9GM8T0J+sMZWK1g02QAvvmv90a20q7r51YLkzI2Rdi2ek3SfnAMGrwmefgD3YBkM2GMcUghY8r2//FjNQWtvjGsbOpgE6/QXRqgmBZGsSPSYD+eMF5CzuNcSCQ9gfpiBQSzjr0OH0BSrssLmnIUSqw1UjShhuWaSvx3iW4YZkKB0ulANfpRwIQD3QzO9rmOv3Z7+nXBs3AaKlHoxlERlrW6Us8HicbzBIhwie3uTLJBXedB7hz3AJc/0/X4Tt8GAJ+Ortu5orsmdxS30ShlcRI5TA2JohEqqEiRagkS65jLC3P3EB8y8Nc8eLfeOTu29muAjw/YSk+K0ayu4pcsVB++izOfLIJgJkxg3O7P0X5uILeHnrT/atYW97N5kAhM2NhFixY0CvnXIibM0e/pB18tNSj0RwhIz3qFqC4eydG+w6QLBVmhkI7woaCjWyJbu3dxhdI4g+kGRvOMTdqYPrc0bBFnV0UdCVQGAgOQZXFUA6GUqQNk4zh9geV+HBkX4I0hQ34CJemybTXsPfJd7MsvpPH3n4lhjJ7twLB8vLrVLXnuGNJEidST5NvB4ayccQkKwHEtvG7r4ppMxIUO0EuSI0jZ4HfB8U+HZ1zrGipR6MZJEZ61C1AYaqeMZIEM0GhHSHs+NkS3UpbYG/vNv5AGtO0mBs1qPLvm/y8oCtBMJNFcPDhxlkaSiEoMoYP20t94IjZmwZBIYAPQ4Jk2ypJbJoKyTf4x/zzsE0/huoZGyuA4LMVhSmbyQ1pxgWFJrORNrI4YmIZPi86x+ndp9gOUJsrcO325BwdnTN0aKlHozkKRnLULbiROssbdrI8vpW7s3cCsHT8amqo7bVr2fLrAFjZ2kWjBS2L5pDCpuKlpYwXPyXv+zobYzdwz6YF/J/9IJM6dvPItin4lOKiojjL4xcDsDG6mT9lTiXsK+NjxaU8OfP7rM6+QlXVeLLxOIGMzYdfTvGFO9/sGXcFAFef/h1yY4sY/2qKFa1+qqtqeiOLeuWc+HX75rfVL2iHDS31aDSHSY/EM1yROuvqOykvzjB++uPEGmP4G8PkbLeXXJgN0x5oZ0Z1K7XRvWz2NRIOdOETH0WpQgq63Pw0CjDKM5iNPioXur3yXSpGfUkMX0HPy1bB9DmkVICWTAGYZWSKr6TQCdFpZnisuIUuu5BCyXBlzuIfp1Wzs2IMIOR8QjCnuP7f7cwP/Y3tuWIM5WBhYCkTH4KJj1aji2InwAWpcQD75JxABX7ZQsXkJ7XjHwKOWeoREUNEzhCRK0TkIhGpGFwTNZrRTb7TH45InRlVhfgLV7KhbQPRlijhVBDHcTtqnYEUjdEUM4raMAt2UhRMMDbg4JiFxLpSBDPZ3vgZpylMwSpFWLlRM/UlMTqDQfJybpJSQZopArMMX+BUEKHTzLDdnyRjpjFC9ZyaNRjTXc7OijHkfK7r8FswttWiuPh1dmYL2UsORwwcZbjvDTAAxRgn1CvlQJ6cE2wiUrZdR+YMM4eUekTkFOCzwFuBjUAzEAKmikg3cD/wK6XyhuZpNCcowyHx5Ms7C559AJjGtMw01mU6aSk9ez/p5+W/vJOdyRqs0G468fGr5Je47Xc/YkwiBV9p5Bt1H8fXXsqPir6A9VYfT6lrKFv6AgVWFn9JjHS7nwv8hXxm0rkAvC82j5ZdXaSjK8gZGew3raentneuvZW1iRSO36TQNJg5xh0T8PilU2Dh3TyQOJ+yQCEf+fzHdGTOKOdwNP4vAz8GblZ9dCGv138d8EHgV4Nvnkaj0WgGm8PW+EUkqJTKHGrZUKA1fs1I8+gbj3LPS/cwr3LekGr76+o7CYWSjI88TVWymN2RrbQGmzEwUY6DWDbRcIY5JUmmRSwsMbHwMV66MfcI8QcjmK0WvgKhZF47PWLLTjvKuvZKMu1+9kai2L4KUsXvBDFRykdGhM0Tg2yrCWMb0FxksM8zuO8BDCU4ItimQ1myg/9+7Ql8TiE+FQPHpkVsilSImSmL4lytDskcBQxGOOdLh7msp8KQiLwqIq+JyFoRudtbPlFEXhGRTSLyiIj0P/RPoxlF9MTtD7W2P6OqkHSwjspkmMJsmKZwM12BLsSxEUfhNwUzmGFaNEe538bCh8LA3CMUvWwgSiFFBpF4DrxkBwaKde2VtGYjoBS2rwIzOANDFAoThRBEsaUmRGOxSTJkkDNdbd7FDdF0RFCGojTZyfSG7RgiGCqKcvxgmBQ5QcKJbspVLcWBch2SOYo5HI0/DowFwiJyBvtmTy4EIgfZNQNcpJRKiIgfWCIizwCfAr6tlHpYRH4CfAhXStJoRjXzKufx3qnvHbLye7T9cxZ+E6NhIpNrJ7GnahdN27dwffP83ukBX/zjFSTZTlfbWL4R+gKza+J88+vvQCwfTZ+ZSMt9acyOCMtnzAbl8P9OWUnqCZsIiujZZ7Gt8S0IJp+L3cIzzp0Qn83Vt83lr8/8m2onQ0ekA9PKUdj+O6ZO/cp+c9Wy8AoWNkyD+GwWfOE77qQtQMXNc1xdvwjc6cqbdS9/FHNIqUdEbgBuBOYBS9nn+DtxX+r+8ZCViESAJcBHgaeAuFLKEpH5wF1KqbcdbH8t9WiGi+EalduTQO0fXUvZFtiIQiEIJjnmFnYyuzANKCzDwmcZVDgR/GFvikEH/L4MRtJPwV+ihDam8bXlMGMO6wK1VLe1sbYqTlO8HMPxgw3K2YtjlrL0jP9iW20EC0VMusmpMDnDoqkwSc5fiqgcoDCzO4h1PMap077qOv66hbD6MV7YPpHnKKM6UMa7Ss8nV5/AkhyJzF5yqW78ptLyzijiqKUepdSvlFIXAjcqpS5SSl3ofa46lNMXEVNEVgJNwHPAZqBdKdUzjHAX7tNEf/t+RETqRKSuubn5UGZqNIPCcI3K7Umgti2wkUSgFUHcfPPKYXZhispgFkHhtw0COT/+UBbTm1s2lLOIJiyCrUFCG1MYbTZmzCEVDQLCruIyWsaEMbIJRJnYYmD5yskFp7CpNkxjsYmIDWKAkSPl78LyxQBBnASm1YA/8wax4vO4pnKMZ/Bj0LCatdQAMCXs/vVXxciku4lYMfymImJktLxzHHA4Us9/Ab9VSv1hgPWnAFVKqSV91ymlbOB0ESkGHgemH65hSqkHgAfA7fEf7n4azbEyXKNy4/E4HZkIZCMsueUJAO75xNUUl7YTLQgyp3QGNUtv5+Zt9Xzw9O8QDfkoqynAt3gzaTtI25+KKO7YQ646Rdv1MRZvckfMfu7OL/PI3bezZ2Mz/sj5PEc7poK//PRj/PnROmo7kixJ3gsLnuL2+27npaKXoODTWL5aZo4Z71k3a3+JB9xMmPVRqoly4Sev6l3c/Pkt5HwZZk55yF2gB2KNeg5H6vkEcBOwzPv0xPFPBt4CtAC3K6U2HqKcLwIp3DEBWurRjDjDmT//03+9nxVbniJtdKGUgR8HHIeOQAcFyTBV9tspOyXLW83nKAq0YKZNCjfty3Mv4qCUgZ0zmJrdwWbfOB5pvY5/nzEXW9LYkoW8fDlKZREJkPOV02Q6iCmYEYOcaeK3bULWRlJGCBBQGazAeGavX83sTauoyLZQnm0nahrYBWfQHipCcJOw7TW6GeNEOD+1b6rDiFFBt9PEzNLP6dz4o4xjkXq+C8wFfof73uZi7/du4INKqXf35/RFpNzr6SMiYeAS4HXgeaBnmN4NwBNHdUQazTEynPnzF+/5G2mzjaS/AxBwHAwrR0EyzLiWAjbWzGZGdCWFgRZiSYvA3v2D3ZQycJSJ6ThsNceyqb2KxWfMZXdZDIcMONb+FYofywiTNBViCuJzbwp+2yaQS2BJz8O+QkkQQ6WZuWkVZa31lGfbidhu/vu9oQLajTTivdob40QYbxXuV1W304Qvslbnxj+OOKwkbZ5k85z3OVyqgF+JiIl7g/m9UuovIrIOeFhEvgysAH5+hDZrNIPG8ObPNyiwSqmSe7im2e3v9ETqXL1iI8Fmh0jCpHZzhr9NHMtl6a+zIOE+6T51200ALLvwGlJWO6FPtkIGqlu7ee/Tf0I57YydNoUfTXLH2ca7K1jd1MqUUxI8fd2VzH9oPrXdHXynOcSz8Vv7nYbxkReeheLJXFu7GojCgqd44Ks/oDTt52pfhX5ZewJxOBr/Fw+yWimlvjTAilXAGf0s3wKcfdgWajQajWZQOZwef7KfZVHc+PtSoF/Hr9GMNvpq+j16/mCUWV/fRKaro3e5UgplGojh5bMPhukIpTD9hTglO3hpYpZZxir+9MxZjMl18nEMVDiHX7IEzAw2ih8ueojamhDLZ8/hP379BxCQBdfhiA+VMmksDVLZ2kzabqMlVMb9vBm11Y0M2pFpRRX4SGbbeN8DM3H8QsiyyWazvTYmXqmne6UbMbdy9wralJt38f6tFwECd36bVpKUED2mNtKMPg5H47+v54MbZRMGFgAPA5OG2D6NZtDoq+kPhp7fU2amqwPD3udUHUAMAzFcbTwVCGL7awlFTgVglrGKuNQzJtdJ2M5g4BDohsK9wuuxal6oP5Mzdr7GilmzaSivcp2+N5JWGSYKg8rWZqZvfp2WUBnro1MBEC9WIxC0KQl3UdP6AqKESRl4UzrArjHzmT3bTZzWvbKZXH0CgM2qlb1GNyLgViagFKUqzBQKdHjmCcZhafwiUoI74vYDuMnY5iql9h58L41m9DEUmv60kmnU/svNVHn3/T8C4Nr7X6K8YRGxWAf33nYvV6/YSKLrde6NPczXl16Dr9IiGjU4dWmZW8g4d/B6xc1z+OR9vwCgKlwHIpQk2vl06Gfu+u90smX8h0nH4jSG7gfg4UkfAmDDHe8AYMEv3SCOhTfWAe9l4UL3eBfccqCu76+KUXHzHJw7nqHI8XPLlz61b6U3oYqO0jnxOJxwzm8A1+D29n+olEoMh2H56HBOzbEwmBOo/OynD7F456PsLHMlko5ogqJkjEtfLqWrsoSq05JsLpvAssBZ9OS6mSdLmWKsI0iayamd+DMGmZiFavfT+EwFSyafx/Lxc7HFROXNcSsoGv//9u48Tq6qzv//63OXWrqqel/TnXRCEghZgISwLyIKsrkwIrhD1BFnxH1UcPiioqOMX8ffiPOdAUeJy4hGdBCQHQRZRCAbBAJZSNJJ79VrLV1Vt+695/fHrW46KxhIOuk+z8ejH6m6davqnKI4ffp9z1JTR0N/kn94qAdXWYQLDpl4CzuiFv9bE4y8KTgW4bDLqTOC1HWDyjPHt3mr8Vmah8qp7guTwCauDF4yutlsJFFqdCUfH58iQ4ZDpWfz6dmPv1rZ7nV6eOZh7o0s0vYlYBpwLdApIqnST1pE9rwhqKYdQt7MDVR2rHyCjqoOhmNB/6ciG2d6Xx1D5T7FinKmVexgbfQ42s1mCObjstBYSyvbSPgj2EUfW/JE0gaDm8pJ5mOsnnU8XXVNpUY/iIZGY52G/iTHbQ4afV8MnFCILDnWW95YmcJhj1ZzAzjB5bijJMJ8Zzrd3d00phLUqDhR30T5Nq8YSfolOE/hlzZQh0rPZqbbs3Nl9fDMSes1ox6llN6QXTvsvZkRj2+GOLZ5wU6vt+y+ZTS8GCUSjuB5EapTIR6/5H2sWv1BUn1dhDMz+UXyWq7r/zJzMWn7Uy0baz0ajp5HuM6jKpPnvKdswOQ9H7uFkVWriP8/mx01BlVOJT+4MRj1/P47f8qthUUMJNKcuvjPr5Zh+c3AEXBF0Dtfvnw5NEKj0xg83r0OgNCMKuy2Iaqt3rGhpMETdKwzleg9d7VJZX9n4+7peb2pAn3ZYLsJu5DGdkYwPIfBcodjmxezqLKFWOYZ8A3q3X6q3WEest/GE1awm5UohQC+qTBciA4rCphk3RDKFVwjWJatt7aFxoEcH30khZVP0jB4KygwHGjq7eWVyia+cv5nMVSBRL6NrDMNI9LFiY038/3uIgDR4vFk/LMZkGq2+Ek6rDQGxlicA+BTpN8cwchlmTv8EpedNK6yOtaZlN6M9fg17ZC3v7Nx9/S8vmyBkUIwI3a00ffNEPVGMxcccQFl2WeYJp0AVLvDRL0CT1qn0Wa2jm10jlKYLtiFoIOVdUO4voEvgukrDInSMORyzJZBzEKSsux6UDC6htorlU08Mm0JhioQcoNkNRbuYkZsLW/NeKCCqCbjn43vB+vs7DCHGDJyY42+wsOniEOeCIpap4+jYx07fwA61plSXteoHk07nOxvrLPr8y67+SkQWLHslLE9ZMfHIzft+AmdfjObV57Hpfb38ZQiO3MOx86tJvRM0Kf62v8E52//x+AXxMCdzeBk2XH0lzn5yX9l+/QPYRqd3PKROC+eGuZIp5fl590Pyy+kt20LZ/FZQuEo1zUOALBs2WdfLfDyC6EaWHb32Lr4x1x5Dvd94zoq/QR/13I29Vces3tFl18INOvejzMM9gAAIABJREFU/RSme/yapmlTjM74tUPa3jZGGW98Fp+XHUTUdCq8Exk2n9nj+cfHd7A4OoyhwDB8DAMcBYZYPBP6O56XJTudbykXe9dF0ABlgnIFZ8BGFCgROiqn09TbzuV/+BVFOzF2rl+ayOX5w5hGBWWhs4g4paGZZRar6rs4qpjE9nf+/zGPEFGKWcrHELO0g+7uxg/NHDQK1Phx3ms1UD/j1t1P1nn+lKEzfu2wtLeNUcYbn8WPb/TzsmOP5x8bG2Za2EUUGAYgipCAb1awVo6ni2k7nW/7LuO3Hh8lHrg5E9czUQK+Mmjo6+bYl1ZTtBN4ZhgIfiGMMo0KbGsG+ELOCpORHBviLrPdJBUqyOVHKSCiFBWlztnoxdo9GT80s8oPc4SKUFb+4p4/MJ3nT3k649cOea+V2Y/P4kctu28ZsGCPz7v1oSUUBlxiL11Aw1n3AbDm4Su49bjjGK6tZnEiygW/CJZfeGheiGs6v8isga1s2VjFwPE/wHBNbvzALAAuvP6LpGIRrLoWbjz3DEIhk8/8/gHWVHyCTU1zWHfpDG5f+3m6uru4r/EzPNb0GLCDkbbzAVgR+jYAy7mEYleW99ScuVMuP/7awvj9bcfoYZjaftBRjzbh9hXn7G0o5q1Pb+eOtcHIlPVdKZpa1tLc8jIAyVySFn8LS6KgSuPYn7NP5Hl7SbCvreXgYtGhplPKaPA8m6IZTKBq7i/yvb/cxrzwn7ByvRhJn9TWKPl0nPZpLXSU+bjKBVORtS3iXpHFI130xU+hM3YyVsEnE2+hJxziyQqHJjqoNYfG/mYQoIDQgNAAKIR+yVKjEpxXWEjWe3UiVdEF24K6BBQL9djh3p3jGx3baPugox7tkLWvOGdvQzHvWNvB+q5geOP8pnLs8ufGXqM/18/iMpdptofnBUsgPGcvocNsCZ6shBQVFCQMSvBVcI7tKaoyHrO255kX/hN1dhsIDG+L4aRsCMXoKPNJ2YISQQxFwndocDKY+PTGlpILTQM7QkHlecF2EeVRaw5RJsHyClL6iQLVGChMlIJqP8Ysr46inwH/1Vm5tgVlpT1Z7HDv7vGNjm20/aCjHu2QsD9DMOc3lbPiyiDeWXbfj6kneI1l9y0jYT3HYL6JP/ddy4orT+Fn9z5OdTrPhx7wOeqS73K9+TVmV80e21f29n9bTW/bFiprVvKL+JFY+V4G82VcecXD3PDXTyMJmP/0Azz8wY+B7/PCKcdx+pG340bTzH1umHnV83h6TS2dsTBnlM/iV5lBuuPdrG28jeXdR0Hjot02PxmNbh7pDnrw7/j6F3V0ox0UOurRJsT4eGdfM2vHRzoAg+ZjDJvPYKl+ykNpRIRF0QKzI3lCIjQZJo8YZ/OUeRoOYbpVEwAFCRNWBaa7nRiWQxuzqO1K87X/+iEKYVZzFy2NSUChgOp4hv6XKxnoPBIyveTrwrxwQYjeleWIPY1s5SUoBD/cA5FeAEwleKIwTIgFc3ZJUKDPcKkNVfOemjPJDg4wkhoCzyFGFVmvl0e6bqUuQTCTVkc32ptIRz3aIWV8vLOvmbXjIx1gbLRORTiDbRbwlOLoqENzyCcuAqL4q3EabcxiiMqx54VVgYSfBhS+G6Z5ZJgznv0LLcOdINDSlKQikUFQGMonm4mT6qiGTA8q3kDqRIdMux28mH00OSUoJGj0zWDRM08UBYE4QhSwxAfDpNaIMzc6HYCR1BDFfB68ItliF0nnReoScPToQCId3WgHgY56tAnzeuOdnSOdcmABV9UXWN85zB87rmXR/B+xfWsPLz53LjNm1WPPjmN6wrzKGbz3phvxHY/zZp5G/YxbuX3gW6VI51HM59ezrXYa/3HNFcxc++/0erD96VqOGO6k7PQv4YT/h0x5gc7zvoM965OoPxlYdjPp+DG8srCM7TNcTlhvsLi6hWX8Lijssrt5/Lu/BeCMay7drS6PfPNWCFPa1xbds9cmhI56tINib9sejjb8u0Y6ANudh1lQcx8nlA+zJvQWnjNPQRFcHA0LKGXgeTaW6eB5NtlMFdlQiFQ0DsD0nm7ef+dPmTFYJJZR5FDMmtZFS1MvgiLsFSmYNlY4wQ4Ws37wLNx8CBCKgONuYIQuHCu4umq4SYpVrUjNDCwDMIWcSpPA4njqQEwIxahPV5Knh5cHfrPb55BMByN0LmvVkY524OmoR5tQr7Xt4a6RDkDGfpallUmaw3meM0+hS2YigImglIBvIig8z6LoREBBzg62CKxKDbNg8yZMq47G4Qz1w+24KFqaklQmsqAUBdPGMSxsO8bG/NkMSTN4DkrABUbowlNDY+XxrTqcimpcM0O0lOHHiNBILGj0zSAK6k0M0Zd+geSwv9vnMBbr6EhHm0A66tEOmteKdsZHOgAnLbcwCdNUs5B4JsFc+nnwzA/B8gu5Yf17IBSjYU6BZNtWqjPPYYfey38vqQSzj5ufHaE818iVFxzBnNv6SafL+dbp/8CC0LeQ/jae3LqIjSedz5mP3YW38ON44Sbiw+tZuuN/eeCsc/i118I53T7TE2kWXxHl+CW3csq9j3PC+lUsLFgcLfXcsvjevdZnxad+Tl2FwWU36R69dujRPX5N07QpRmf82gF328bbuP6p61nasHSPPeT3/eoHrJI0xYpjsA2XcrNI1FQoBWFDIcpgu7TS7Hbw9sc3Y/kGL0xrZWNjA0o88D2sYhEx6+mstKnu6+fam25EKYXluRxfs5HorCIjiUZ6OuaxxTsdX4XGNjksImSjdcxUbUwrKydXyvQ7rSQd4QyuZeGXpl7lVJoaP8Zsu5quGT/jCj82Vo/nd8BLwerLJId93ePXJtxBz/hFZLqIPCIi60XkRRH5XOl4tYg8KCKbSv9WHagyaIeG0Yu6exuy+XL6zziV83DDTRjiUWYqQqIwBDwlKGXQ7HZwVF9b8ASl2NDYQDIRR/AwfAUSwxSfUKrIWc88yczBDpRSmL5HbGaBRHmOKoQt3umkwsF+uBAsl1AUsPJJmmKVlIfrx8rVHs6QMgrjlmdTxCRKdSjG5uZuLshmx7Y0hKDRT6aD23UVBkefsPMqn5p2qDiQGb8LfEkptVpEEsAqEXkQuAJ4WCl1g4hcDVwNfPUAlkM7BCxtWMr7jnzfXh8XFeLkmho+0XUNIdWDn5nNyOZvcvGXgsZz2c+W4gFf/9ZK/vuzFxHy5nF83XS+/Wgna1JlxGZVcvGXlnDZzU9x+upn6Z7WxH+e+gmuevgnYBkMchRDxvdQbCThtPPB//kU13/uYs4PvY9rKocRQzhuKEEfaU6Z/zMA7hhYhOM4LGps2W3WLQDLH4do/asjc7ZfTV3Tzpu1aNqh6KBFPSJyB/AfpZ+zlFJdItIEPKqUOmpfz9VRz6Hr9ayXnx5eRc6HwbK3kwyfCLDTIscRwydjtTJTtXEt1wGQ65vNH7ZfxYutRXzfBN8Ye5LCo6+6gdbOLr7w05+gEGobj6Wx5USsUAxDBaNp+gwLwxymUhXIuFX4voEybLaaO9gc6QfAVN7Y+SGxcHBpi7YDQnkxRg19HJdaT2e2YveKOVkIxaBpEQDJbVupmzlLN/zaIWNCh3OKyExgMfA00KCU6io91A007OU5nxSRlSKyMplMHoxiavvh9ayXn/Oht+iTDJ/IiD1jt5XtE6ailW2c7D+J51o4uXKG245kXUuY7oomUCZKvfpVFUyO6Gjn7U8/hiPB7rY19fMJR6owzRCKYFbtkCEkJMWIX4FHGFE+4hfZHkqSkwwAhvIRFAgUcUmZeYrigwi19LNUbWBHKj4W4ewkFIN43djdupmzOPq0s/bvg9S0g+iAD+cUkTjwe+DzSqmUjNuUQimlRGSPf3IopX4M/BiCHv+BLqe2/15rmOYlvzuZVhss1QyOx4ztrTs9/rHG7wDwjgEfmM2qYyt45aGXMaa/g4VulI93TmNk7ToWb13B/CeCvy7uvvh95IYK3PbRaym8/DKfs+uZCWyMG3xjaZSnzj+Dr9z8FF/u/hEbhj6B22dz5LafUPx8OVvazsGya/g/n/lHVnzzarrbXqF87nSWXfOD0jr+BPVZfiEQZ7ssog4d4WiTxwGNekTEBv4I3K+U+kHp2AZ01HPY2VukMzoD99uLzuU3L/2cJ/r68ErLCisV/LwSPZ1i7DSK1kymq+1c7X8DlUmgisGImEhlL32FafywcCWCoMTHzio665s5oqubb/34FuymhcQaF+CLAUoR8Rz8aDVuWTVGYYiwXU6mmGft6AXgEhOffiPNSGQ7xfIQhXAIMWyMQo7wjpcwPfBMyCQ86gwh7oxgiEGZVTYW5SQLcR3haIeliRjVI8BPgZdGG/2SO4HLS7cvB+44UGXQ3jx7i3RGZ+B299zFE/0dtDsOqrRvrFLBWpfF2GnkrVZa/B2cIY8QNp1So28DivxQHU9yIj0VwUplogzCrs+c9jbOXLMG0ysQa1xINN4IpUUb8mYYiZQTBYqEGPYUPSPDu5XPR/AKL1IoD1OMxkAJRiGHPdQHBI2+E/IxDIPyYp4yX2EbpcXYSlGOjnC0yeZARj2nAR8B1onI2tKxrwE3AL8VkY8DbcDuK1lph6R9RTqrVt+Oacao9av4O+9Sli1bFmyJCCxaUMv6zhSfUj9jeqKDRDbCykcvh1CMv14QjIgZarichZ1ZbhmcQ/2Vx/DMye8EINkUITm9kr6qYBOV/zr+Dk58+V1MN3M0V/pYoRDda9YDcGNrFVDJzaEqqp12Wub8nvsq+nnhzwsIRSLMnHUEy5Yt4/p/OBfXgOtufWDnSui18LUp4oA1/EqpJxgdLL27tx2o99U0TdP2Ta/Vo+3TaLY/muUDdHT8mt+89HOeGnp1ATPXzdJsuyypTRF3l/P397exbc4xCIruQWiMdFKltrE518ItQ/+Md1YY0yuQ2P4OfDHZ4vcxq6uTbT+5medvMWidey524yIqDRslYCmDrXi0vHQFjcMD2AVFb+crjHhtZEhSTIR5uxvHEMW9uQKg8DfOBY6kLCz4IgxuWcnyj/2MilyU8mju1R7+qNFNUDRtktNr9Wj7NL7RH515O5rnb8/lxs7zfYulUUVTyMM1QjxvHUu3BLtfNapOjvHX0OVN4/ncaXhGGNN3sIsZUD5KDGZ2dXPm2iARzJoGduMiTDOMkmB27VY87sOhqWDjSxW4wXWEnN+B8pI4FVGMsI2FCyiUAKVhnaapiJPCGuonVWr059aP7F5ZvWKmNkXoHr/2mvaU7ZtmjAV1R48dX758OXFnOf1+go+8627+6aE1mK7DJx78BYvPfoXje1ph2YP8+u9v5cT8S7TOL6Px4RVsCCvWnvglFsSjMKuVqxc8Tmf33/GHQhhThKO+exZn3LkGefFZwsUCF6cqKM9lOb79hzx5wQU0tAUXadNlBfI0cdUN3wJg2X3LeNp+L0dVHzW2r+6KT10IFSaX3fTng/jpadqhRy/SNom9nlm1e5LMJenPBTNbc26OqBXlvEqDWWZwDOVRUGAR5TnvPFaFgtFihuVSUGE63RYcO0S4WOTaFc8T9X3EsUC5pGPTsHJJflNeRIwIp0XqOMWKUIlQLkKEQUIyhK+m4Uonw5EfAor2zOkMDS0iWRbCCXcSKnZRCNu4vmAbHn44TiFUw2+XNAMwgs/RrwxwwsY1LHCC+YJ64TRtqtEbsUxBr2dW7Z705/rJuUGME7Wi1ERrmGX2U2sWgxPEJGqGsYohVoWW0m4GI24KKsJwaZ/bcDHPkZ1pPLMJKQDKBaWwcknaJIMPFELlnGpGmIlJBUIUISRDWORxpZOssRohWE5hsHAMWZpxwp0UQwUKtllq9BVRU6gSg+ebahgKBWUpw+DYjaupHHh11rdeOE3TAjrqmeRe77624+00e7Xk0T8fB0T46NvXjh378nd+SOaoSirdIR6/7FIuXrOJJuD2xXNZdc9CXtj4Jey8yXE7/kDWSDGQcPjq6f+IVzgSJ2Tw/kyYFtOktyzYcrHCqaQudCMNjU0c2flpGo0FfHHuKjasehdNNTESm9vxCn1EYj6xvs3QtGinSVU33ruCSqed9edfBsCKZ9ZATZ2eeKVpu9BRzyS060ic19Pwd3T8mnVbf8pAfoBwMU/M9DGNcevjSIGH1Hk8OBKMr8+YUXJ2lKJpYbsOobxDPhqjaaDAlfe8jFJCMdLAkWRotgp0VNUCQswP08QQERmiCCjVjCedZOyb6MotZZuayUDYYYQIBj5hfxjPc1AYiPLxIjaRYprI9k30107jT+e/g4HSXwVDoRZOWvcn3j+cBfSiaZqmo54pZE8jcV5Ld89dqMJ2RtwRYqaPZbo7Pa5UmCf8s+iLBtsnjDX6XhHLKZL3bOqHXY7Z5hAquigEO99DY8QmlJjG6JSOaoSIDGNQwEUoSieOsRpQ9OSPJRkSCqaLgY+Ni+87KOUhysfAI1TIEs4N4leU0dvawgA+I6WGv9JpZ+nWTSS3bQX0ommatjc66pmk9ifiGVRxnuRYLnG34blw4QWP7PT41x/6X5r8Dlq7mojWB2P4n7j4XC67+SmKbb28K2PjpJJUdt7Gf73jUu792t+z4wsrKGZ62BpuprYlQVXcItbxVbJGireGv0f1SI5l9RFOPW028UeOh75HEKWYMaeGZd3X8ZuNR5MfDvPWWBWtv/wFF6/ZBMDtaz8PPMsyCdbmGa3rimfWQLVeUE3T9kVHPYex11o4bdeGf/T8OUYXs4xeAKLiEMHBwudh/1ye4Cwsw0WpENU1x9JTKNIxkkVUEUfChHyXxuECXYkYsc1pElszvENZnG6EAMEg2Dmr0+oiafQSRvBQZP1Xx/wbvsL3BVcVUBTGpncrBD9sEc6ladq8DlAM2xHi+SIJ3+TrVy2j9ZU+FmxcTXlhGxgmIyKU2WVjk8t0vKNpr9JRzyT0Wgun7e38WUYvVaX16ING38PzLP4iZ9JhtqBUCMuqBqDPcXEwEBRh8sTyPsoLUz1UxOjKM2DA6UaIVjEBhSHB3rRJs5eMkcEFHF9hlOIYAN8XUKAoBKN9SgQfO58jNpgcO1JedDGNMA+ceCYACzaupnYgCYYJhk2ZXUZ1pHrsNXS8o2mvTUc9h7m/NdI5qvoojqouAHD+kltZtfqDADz/3DkUG5s5rrGO2xefOnb+xWs20ZXczDWD1xFWkFv1ZWhcxM0D7Wx0DULiMSPXgwJmrf8dhUWfwoiX8XQ0jxrOYPTNYWZ4iDlzV3CT/14q1AAX37eaeCHMmqNn4g3vIJS4lLfPuIk6ays0HvPqImmlJRXmz/gUAOvfdpEeqaNpbwLd49c0TZtidMZ/GHi9Wf748+YYXSxQHUSkyKPyNp4yXu3Fl3YrZEDVkDYSoATPtyiaFlE3x/TeIYpmAqWE7poy5qY8vv+XIe4zfP5sKk4x+4hZfZjiYuNiEAa/gJ/eQd6wKJgGhWgCs+AQax/A8guU59vxEYTSXC5bSEdDJKIuS+cN0pEMM5g0yZlhXiqbXippsPa+a4SxKHJCVY3O8DXtb6Az/sPY683yx583y+glZuWwTJenjNPYLjOB0qBKBUoZpKWcAhGUKmX4qkBVapiiJPAkjALmpjzO7w5y+MdMn23iE7P6GDGy2LiYwRxc/HwK5ULBNHDFwCw42Ok8ll8g5KZLTXhA2YKKGiSiHk01eQTo6i8jWYjRZyfwxnblFRDBokiVFXxVdYavaW+czvgPE683yx89b9XqDzIwsJaC00B10wl427byj8+v47IZ64ITl9396tDIxXOD8/tWk/hpmOdrv0iyoZVfN/l8uytLS6afe8r+Sp86lgYrSixqEyokuMT+PVvjDlcccy1/uP6rANx79AnEciES8XdhyjYuece9pff701gZz7hzDQCPv2vx2LH137yaGHDfST2AYvl5t7zxD03TtD3SUc8h7raNt3H9U9fTWt5KXbRur+fNMbpo8LcRc8u4P/ReNlhHEcKhoCJspxW76DCrtx3L9UuDNw06pzUzd9jh3x/dgmcl8K0YYcMgZIbwBYxCJ0YhiwI2RNNsieQAxaDtUe7aJF68n7xViWMlUCqIkBwyeOF6QolLSSUcHlwSBTEhFCNtDJE1UmRDDcScHk43/kDVhhyVWwpEBlzy1Rb3ndSzX3MQNE3bnY56DlNj2b5inwuuzTJ6abZ9QkWbDdY8ummiQJiMHyfkOpQXMoRdD3EVfpC0M3fY4fzuYJilb8XwDZuQaWEK+CjMQhZchQK2REYYsBwEKHdtylJpQHCsBK4RRgSUYeCF6xkpn08q4dPWYASNvhnsYZs1UjiSJ+b00Fh8GWCnRn/oiPDfNNtY07T9o6Oew8DShuAXdl1Z3T72vA2GZT78+BzUGUKTdPOXtlt22j+27SMfZX1Xiq+c/g/Mbyrnu+1teN4Iw1v+jZfrvoKqaKKmOgwovr/U5r//7QYKyuYf3v4vvN9fTQ0Gn/rOj1jxzauhFpLhBVTmLiIcjXLp94Nx9mfcuYYKfzMff7KDRd4iLr7k1dUwl913IwDLz18OvBuAFU9frWfaatpBpqOeQ8yuI3jW9a0jgcsHqkcI+SaNTsXYY4+GTuee8DkoE+IEPXClhB5poEl1cfedP2RT2zT6MwXCVidF7zi2NZxGeX0jp0ceotb6Cx4+PcUTEJawyUix1XBRXh4/2x709i2Dzlg1EcvHLjjE24cpen1YZh1+w6XgCMNxgweXmKTtLH2JCsryncze+gvsaIjaljgAVRtyhDYM7jTLFvRMW007kHTUc5jYdQSPbdgsiuZptn1iXninc58KnUifWUuackI4Y2NhGlU3JzmvkGqLYm/dRLbgUmZm2NF0CoVYI0faJrXmXwjJFgrKR1iCpaaxxczSL2lUPh28kCXkwmGifg4jn8NKBQurWWYd4chRDIdNuitNtjcGX6OsPYLpbKd6+DlsG6IJe6yslVsKVKdDO82yBT1KR9Mmgo56DkHjL24uu28ZS9y/YqZred/f/XWn825es4nQwDCt/nZ+UPwOx1/wwk6Pt936It11Zdz2vq/y7tmfpuOuWgSXLjPEPK9AlyrnwdmfovmlGVhGAS/UQ0XS4G2r19F66iv0dJ5O6y//M4h2gELseGgEb4YDwEM1Fah8hts3/DN1lz726jr+H9o9jtKRjqYdOnSPX9M0bYrRGf8EGZ/lFwq9ZPNduL5LT1FRb0Fr9ZmstU5DocgZVaRUJbYdJ+L1UaaCJZEVQg+NzFRb+PaObzB879vpiy7G8GxEGTRFo9QnynEtmyrzQfqLOTYMLaTRLKc2GuaVUA/toWEssXCVy4CkMfMO8fZe8iEbVZpyZTlJ3FAd0dh7cVU/nfZDAPhGGMN3KC+0YUTijBRHdsvwR+ksX9MOvoOe8YvILSLSKyIvjDtWLSIPisim0r9VB+r9D3Xjs3yn2A9+AU95NNjCvJDHWvNUuoxWBGGYSgoSAaBMDRGiAApE+TT77ZyaeprI0wYDZceSCzUjykAQ6hPlxO0IrjKoth6nrTCXlDRRH60gHorTbg8xbOZxlUvBy2EWitjpETzTZHx3wA3VUUgciav6UbkXS0cFw3cIuSkwgix/15Uyx9NZvqYdOg5kxv8z4D+AX4w7djXwsFLqBhG5unT/qwewDIe08bNsNwy8zIPusSw/bzkrvnk1zxxTx+KZ9dy++LSdZtjefc/nALiwpyx4kWV30/aRu0hnykkdHaEjF+XDDXPIFUfwE104bpFvdbl8I1KggIWR70csi6zhkzW6qIrk+fTMHZxR820gmE0bZPUe73rxMwBc/KVgSOaKT10IEfjyjXejadrh64BGPSIyE/ijUmph6f4G4CylVJeINAGPKqV2zwV2MRminvHRTjKXpD3dTtSKclFtBSdbm1hRPJeneQuVYqB8jx12Ky1GD5d7/01BDEKuS306RaR8gJFULUN3fBwbg2llCeqjZSAGnhmip/dlkirY/qQ6UqQhVgOl/rvCoM0eoDuSxVUO/ZLB811G0r14RgTDd4g5vTjkCRGhJtuMbXnUhoPN0JPDPnUVBpfdpBt+TTscHCrDORuUUl2l291Aw95OFJFPishKEVmZTCb3dtphY3y005/rB6AmUjO2E9bT/pmk7BmI6SKGySzpZr6/dqzRj+VzKAxGUnWkts7HQhAl1JfFiYciIBaG5ZBUETLhegSPhlg1FXYVwfJogqBoDw3Rb2QoeCMYhRxFJ8tobGO7wSblISLE/HLssEnUGAInOF5XYXD0CUt2rZqmaYeZCRvOqZRSIrLXPzeUUj8GfgxBj/+gFewAGo12xoY9lmKe5LYtJKSG+lgfN8T/yPFLbgXglHvreJ6385tb/xuA1cd9HoD/63bwvdkp5hYbiPphBlWONRVJFrzzbrwb30KisIN3zrqOLel/I6kG+Wh9Fb6v+Hiuj4Y59XQMbMDfHizpfNnX/23fhV5+IRDbaQawpmmHt4Pd8PeISNO4qKf3IL//QTUa7xQKvWxK7aDegl/et4hT/AJVAg8+MI9HzLN4rOILtEuCVq+f54bSXHXv4yinSGc0Tn1qiPXhU2iM13DE9o2I8rgVhRKTjAyQKb3XorwFv3gPCyuFttAAP/O/iBlP4iiXs3NtWLiM+EnatsYIucH6OFQDK5fDut/tvRLd66Bx0cH4uDRNO0gOdtRzJ3B56fblwB0H+f0PqtF4xyn2M83ymBMu4rkutiskTBDT5S+cSbvMoIFuFshm1hRPxHEclFukcaiPRVsHaIjXUmFXIcoDFSyx5os5NtyScWNwBEVHaIhh08FVRVLKxcIl7OcREfBymMUUZdPqg1E2634XNO5707gIFl1yID8mTdMOsgPW4xeRXwNnAbUi0g58HbgB+K2IfBxoAy49UO9/qDiq+iiuqi+MjdqJPh1c1qh753YAquQEOtM5piUW8p+LL2b58uUsYTtn/+kRXh54mZdmXAWRKClvGFb/mgh5Gt4+xN19z1UtAAAOxklEQVQD/0outo3a42/hiNUjFEyba0If5csL/oeXX/kAVUPDVPb0svDoj42Nytmj5T8KGncd5WjalHHAGn6l1Af28tDbDtR7apqmaa9Nz9x9E9y28TZuf+lnwUQsoOg5Y7Nwp1efyWDZKcFGJY7gGRauaREyCjgqQoeayeyUy3ef2k7I81AoTM9HMPAlSOIq7CqGi4Os7xlguCpHqiqHp8KYUiBmdRLFYYQwm/xmKv0UhmlhFEao6MzSetzfc/Fbntt7jj+a4esev6ZNOofKcM5J6Z4t97A51Y7nBcMeXd/FUx51piIZPZXtzEI8wfDBNS18Q3D8MGlVzqyMzzs7HOxisPCZIBjKx8cYy/CHi4P0ZPsxLEWqIodjFjGlgG0El3ZzhBhUcSIqmNHruy4UfKqmLebIExv2nePrDF/Tphy9OuebZEY0yjVHHMHxS24dG6553tMN/GedzexIPz9Zfj8An//i/6HYmeULz67jebWJtzrHoPJZfpV5jHtaL+Ih558o/LnIY63XkK6r5EehKIbp8Nemb+D0hLm+8rM0DA4y4jo8cuYVfHGT2neGD7AZ3avXNG2MbvjfgI6OX/Obl37Oyp4OYhVv4aqBU0nefRsR+xzKxeTl0yy6jUYanW5+tfTDnJCOsOz+HioyKRIjNmcykwrTZsgvsKBvE+/ovYM/tFxJ+uQRXHMI0+jl01YXguJ3I+/CiVtUZ1N4SpEpehQLHnv8o23XIZp6SKamaePoqOcN6O65iyf6g+UMsmWn0i4zAUgYYJseviE0+t3U53o5KR1lRg5sz6MslwHfRxDShW6SA8/TlMvQU7GUwQqPQtjFxME2MggKUQpHmaT9CL7n4riK3hnHcGrWCKKcXe0a7eg4R9O0cXSP/w0yzRhhv8BQ0WZpXR23Lz6VVas/yLrePv6n8+/54WabriV3MSc3h36/wHb3PjqHN5HONmFVvJ+zXvwnpqH4+fmfxM7MxbSfp9yqZoYV/EI5eXOKzs1Jnj/h43zlxr+h8dbRjqZpe6Eb/tepo+PXdPfcxaMDAzw1FKyH32GfQLr8IhKVCl+mMdCzkjt++0HKKrK0Zpr59ktCwStQ/eQVeLYBxRz5tgyeF8USRVW2jedarkJZQmVvE7bsYCQ0RD6dp3PbeuIFj9QrnahY894LtqeZtzra0TRtH3TU8zp199xFJrOep4aG2J7LAZCMnErGaiUkipn5bSx1nsUzTazUNKo6T8H0/dKzhSFnkO0jWxBfYVFO1GhGAcoS0tEyIoURnEiwfp2dzhIveLQMDJFJNJJsXsDMheV7LtieRuzoaEfTtH3QPf6/QTw+n8RImAXxYIG1+Q//kTp/G/+845vM+tUxPHL2WxluWkTT+gQFhMpwBe86rYoPPNrHgsxmzHAZ57WvAhS/v/BobsvWkIqUs/nac7n+cxfjG7Mpo4zWxhgnZ2qgvIbWX/7iNculYx1N0/4WuuHfh9FRO7fnFtITeicKKCgoNxWn3/8rctZM6mUYI1pN7piLObkP/D6TcrOC4eIgq7tf4eu/gXSkmWLe4titP0INu3iVJrenZ1EvKU5yt/HNa+7HCB0BRgRJpTli5bPkcw6RefN2L5QesaNp2huko559GB21s8M+nhFrBl5pknPCUFimYobTxtLhv1LefTIJKxhdo5TJUHGInkz/2NJpiXwHTUMrAQWVJmunH8mAKmO21U+NyiKlRMjyHI7o7CSWzhKZN4/yiy7avVB6xI6maW+Q7vG/BtOMEbMTwCAn8QAInOqsYm7vCK/c1YokTGpaPsSGijDri30YyUVUDrez4JX/5fRzX2C5dQmbYjOY5ieIhS/kM5/5Htfd/BTxgQ3UZXoos8toLQv2rL3s6zfQ9pGPwsJF+454dLSjadoboHv8mqZpU4zu8Y9z28bbuGf9jfTmhukLv4WBsotwYoKrqjky28UlmVXYZgKrykPKZrJ4+scwxMRwGmgaybMpb5CqWEs+ViBVvYAHWEKPUwUFly5VC8Dn/vk7VBs2b/GLiETwMil62jZS4RVpO2cx+V6HSH2otPPVHuhMX9O0N0j3+Me5Z8s9dI0M0Fv0GYiewojdSthXNKdTnOg/gl2dI2xECadnEOo4DhELMSJkigWS6SGyiSyelcF2cxA26VMVFDFBKcZvlhLyi8TcEUL5DPHBXiq8IjOcEQAi9SHK58f3Xkid6Wua9gbpHv8uTDFpLSujbjCBW9zBhxO/pOa5DyNHPUJhJE7VumvprjB5mKdIJNqpHZjPsat+RKO9ldpzz6CbOvzyYYpAZ/lbSbZt4fItt3GWXb57bj/aq1/24EGvp6ZpU9eUbfiffuhaBlMPkTbz3NX/YV4sPxHFMpS9jHg4R0dVOTO8ApXpNN4x/0Uk7uNmYqy01+BnHGZmOvEwcIx21i2ZSd6cy1CxkrRns3JgJqIU/b1d1DpZ/Gwa8m27xzc6ttE0bQJM2ahnMPUQRtkAWbPAi+Un0peIg4AAVjhLK8FMXLc01jKfqSHVcySe7+IWs4DgK0EBedMmHSoj7dl0O2WIUhjKp7GQZHHqJWaq/J7jGx3baJo2AaZsjx/AH6nmz2VHA1CTzvA15wssfewLbH7HrdhunqdWHsvayPt4y6N3EC4WuP06xde7f8o92bfiW8JPE++mLN3DF7MruGwxeoilpmmHhSnV8G954iY2PLuaO+0LeW7mdZihEQwFg1U2M7wOQukoz59zE5GyfoYzNRRCMaoHt1E9MMi28kY23n8aH+I0lAjiK5yRLBVeHpwsEJvo6mmapr0uUyrq6R26h97kcayub6QrWoVhOhgIs7yuYIG1cbHOQOd0Ypk0Dd1dbCtv5OGWV3e5EhQGPo1ukiXuRo6eGdORjaZph40p1eMfZZhFGtV2Plj8Lsc9fB1hX3jAz7KZJbxl8yPMO6qXdKKXFclFpJ0sf/joNwB48spTghcYG42jox1N0w4/U6rHr2mapk3yHv/n19zPs1vDFFFYoSxW6CNwFnSbBq1A00gFq6M78HyHoiqnIjuEdVwYMQ163Di/rz0Ps9pnsK2L+eE+WP7t4IX1MExN0w5jE9LjF5HzRGSDiGwWkasP1Ps8MFCkPR7DEx8rlMMwiwBM9zpYmnuOVPdsPN/BK6apyA7R4rQDsMlr5kFnMYZShEyYH+7j3eWbXn1hPQxT07TD2EHv8YuICfw/4BygHXhWRO5USq0/EO/XkMnw4TXPUXnyf1HTHaL9d838cuFH6Gz9I3f2fgbDn8HL37tgt+ctAb56IAqkaZo2wSaix38isFkptUUp5QC/Ad49AeXQNE2bkiai4W8Gdoy73146thMR+aSIrBSRlclkcr/eqDXkMT0zQgiTbC6CNxDCx6LOTtPshigzM9RH7P2rhaZp2mFKlFKvfdab+YYilwDnKaU+Ubr/EeAkpdRVe3vO0qVL1cqVKw9WETVN0yYFEVmllFq66/GJ6PF3ANPH3W8pHdM0TdMOgolo+J8F5orILBEJAe8H7pyAcmiapk1JB31Uj1LKFZGrgPsBE7hFKfXiwS6HpmnaVDUhE7iUUvcA90zEe2uapk11eskGTdO0KUY3/JqmaVOMbvg1TdOmGN3wa5qmTTEHfQLX/hCRJNA27lAt0DdBxZkous6T31SrL+g6H2itSqm6XQ8eFg3/rkRk5Z5mo01mus6T31SrL+g6TxQd9Wiapk0xuuHXNE2bYg7Xhv/HE12ACaDrPPlNtfqCrvOEOCwzfk3TNG3/Ha49fk3TNG0/6YZf0zRtijmsGv6DtUn7RBORW0SkV0ReGHesWkQeFJFNpX+rJrKMbyYRmS4ij4jIehF5UUQ+Vzo+mescEZFnROS5Up2/WTo+S0SeLn3HV5SWLp80RMQUkTUi8sfS/cle320isk5E1orIytKxCf9eHzYN/7hN2s8H5gMfEJH5E1uqA+ZnwHm7HLsaeFgpNRd4uHR/snCBLyml5gMnA58u/bedzHUuAGcrpY4FjgPOE5GTgX8F/j+l1BxgEPj4BJbxQPgc8NK4+5O9vgBvVUodN27s/oR/rw+bhp8ptEm7UuoxYGCXw+8Gfl66/XPgPQe1UAeQUqpLKbW6dDtN0DA0M7nrrJRSmdJdu/SjgLOB35WOT6o6i0gLcCHwk9J9YRLXdx8m/Ht9ODX8r2uT9kmsQSnVVbrdDTRMZGEOFBGZCSwGnmaS17kUe6wFeoEHgVeAIaWUWzplsn3H/x34CuCX7tcwuesLwS/zB0RklYh8snRswr/XE7IRi/bGKKWUiEy6cbgiEgd+D3xeKZUKOoSByVhnpZQHHCcilcDtwLwJLtIBIyIXAb1KqVUictZEl+cgOl0p1SEi9cCDIvLy+Acn6nt9OPX4p/om7T0i0gRQ+rd3gsvzphIRm6DR/5VS6n9Lhyd1nUcppYaAR4BTgEoRGe2QTabv+GnAu0RkG0FMezbwQyZvfQFQSnWU/u0l+OV+IofA9/pwavin+ibtdwKXl25fDtwxgWV5U5Wy3p8CLymlfjDuoclc57pSTx8RiQLnEFzbeAS4pHTapKmzUuoapVSLUmomwf+7f1JKfYhJWl8AEYmJSGL0NnAu8AKHwPf6sJq5KyIXEOSEo5u0/8sEF+mAEJFfA2cRLN/aA3wd+APwW2AGwRLVlyqldr0AfFgSkdOBx4F1vJr/fo0g55+sdT6G4MKeSdAB+61S6noROYKgR1wNrAE+rJQqTFxJ33ylqOeflFIXTeb6lup2e+muBdyqlPoXEalhgr/Xh1XDr2mapr1xh1PUo2mapr0JdMOvaZo2xeiGX9M0bYrRDb+madoUoxt+TdO0KUY3/Jq2n0SkUkT+caLLoWl/K93wa9r+qwR0w68ddnTDr2n77wZgdmmt9f870YXRtNdLT+DStP1UWkn0j0qphRNcFE37m+gev6Zp2hSjG35N07QpRjf8mrb/0kBioguhaX8r3fBr2n5SSvUDT4rIC/rirnY40Rd3NU3Tphjd49c0TZtidMOvaZo2xeiGX9M0bYrRDb+madoUoxt+TdO0KUY3/JqmaVOMbvg1TdOmmP8f3ZLcDxR0NeQAAAAASUVORK5CYII=\n"
          },
          "metadata": {
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "plt.hist(E_0, bins = 40);"
      ],
      "metadata": {
        "id": "CAAESCxgCkhe",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 265
        },
        "outputId": "e460c2b6-7109-4d13-c6ed-0d5d66211364"
      },
      "execution_count": 14,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAD4CAYAAAAO9oqkAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAATBElEQVR4nO3df6zddX3H8efLFpTMYUHuOtLWlcVmBslQaEqNi3GSlYLGkkwNZrOVMJsM3DRb4ur+WCPOBP/xB5uyNNLZOhUJ6uiw2DWAMfuj2IsiCOi4ooTbgL1SKDo2TfW9P86neryc23vu7b3n3JbnIzk53+/7+/l+v+/zbc993fM933NuqgpJ0vPbC4bdgCRp+AwDSZJhIEkyDCRJGAaSJGDxsBuYrbPOOqtWrlw57DYk6YRxzz33/LiqRnotO2HDYOXKlYyOjg67DUk6YSR5dKplniaSJBkGkiTDQJKEYSBJwjCQJGEYSJLoMwySLElyS5LvJnkoyWuSnJlkb5KH2/0ZbWySXJ9kLMl9SS7o2s6mNv7hJJu66hcmub+tc32SzP1DlSRNpd9XBh8HvlpVrwDOBx4CtgB3VNUq4I42D3ApsKrdNgM3ACQ5E9gKXASsAbYeDZA25l1d660/voclSZqJacMgyUuA1wE3AlTVz6vqaWADsKMN2wFc3qY3ADurYx+wJMnZwCXA3qo6VFVPAXuB9W3Z6VW1rzp/XGFn17YkSQPQzyeQzwEmgH9Ncj5wD/AeYGlVPd7GPAEsbdPLgMe61h9vtWPVx3vUnyPJZjqvNnjZy17WR+u9rdzylSmX/fC6N856u5J0ournNNFi4ALghqp6NfA//PqUEADtN/p5/5NpVbWtqlZX1eqRkZ5fryFJmoV+wmAcGK+qu9v8LXTC4UftFA/t/mBbfgBY0bX+8lY7Vn15j7okaUCmDYOqegJ4LMkftNLFwIPALuDoFUGbgFvb9C5gY7uqaC1wuJ1O2gOsS3JGe+N4HbCnLXsmydp2FdHGrm1Jkgag328t/Svgs0lOBR4BrqQTJDcnuQp4FHhbG7sbuAwYA55tY6mqQ0k+COxv466tqkNt+mrg08BpwO3tJkkakL7CoKruBVb3WHRxj7EFXDPFdrYD23vUR4Hz+ulFkjT3/ASyJMkwkCQZBpIkDANJEoaBJAnDQJKEYSBJwjCQJGEYSJIwDCRJGAaSJAwDSRKGgSQJw0CShGEgScIwkCRhGEiSMAwkSRgGkiQMA0kShoEkCcNAkoRhIEnCMJAkYRhIkugzDJL8MMn9Se5NMtpqZybZm+Thdn9GqyfJ9UnGktyX5IKu7Wxq4x9OsqmrfmHb/lhbN3P9QCVJU5vJK4M/rqpXVdXqNr8FuKOqVgF3tHmAS4FV7bYZuAE64QFsBS4C1gBbjwZIG/OurvXWz/oRSZJm7HhOE20AdrTpHcDlXfWd1bEPWJLkbOASYG9VHaqqp4C9wPq27PSq2ldVBezs2pYkaQD6DYMC/jPJPUk2t9rSqnq8TT8BLG3Ty4DHutYdb7Vj1cd71J8jyeYko0lGJyYm+mxdkjSdxX2O+6OqOpDkd4C9Sb7bvbCqKknNfXu/qaq2AdsAVq9ePe/7k6Tni75eGVTVgXZ/EPgynXP+P2qneGj3B9vwA8CKrtWXt9qx6st71CVJAzJtGCT5rSS/fXQaWAd8B9gFHL0iaBNwa5veBWxsVxWtBQ6300l7gHVJzmhvHK8D9rRlzyRZ264i2ti1LUnSAPRzmmgp8OV2tedi4HNV9dUk+4Gbk1wFPAq8rY3fDVwGjAHPAlcCVNWhJB8E9rdx11bVoTZ9NfBp4DTg9naTJA3ItGFQVY8A5/eoPwlc3KNewDVTbGs7sL1HfRQ4r49+JUnzwE8gS5IMA0mSYSBJwjCQJGEYSJIwDCRJGAaSJAwDSRKGgSQJw0CShGEgScIwkCRhGEiSMAwkSRgGkiQMA0kShoEkCcNAkoRhIEnCMJAkYRhIkjAMJEkYBpIkDANJEoaBJAnDQJLEDMIgyaIk30pyW5s/J8ndScaSfCHJqa3+wjY/1pav7NrG+1v9e0ku6aqvb7WxJFvm7uFJkvoxk1cG7wEe6pr/MPDRqno58BRwVatfBTzV6h9t40hyLnAF8EpgPfDJFjCLgE8AlwLnAm9vYyVJA9JXGCRZDrwR+FSbD/AG4JY2ZAdweZve0OZpyy9u4zcAN1XVz6rqB8AYsKbdxqrqkar6OXBTGytJGpB+Xxl8DHgf8Ms2/1Lg6ao60ubHgWVtehnwGEBbfriN/1V90jpT1Z8jyeYko0lGJyYm+mxdkjSdacMgyZuAg1V1zwD6Oaaq2lZVq6tq9cjIyLDbkaSTxuI+xrwWeHOSy4AXAacDHweWJFncfvtfDhxo4w8AK4DxJIuBlwBPdtWP6l5nqrokaQCmfWVQVe+vquVVtZLOG8B3VtWfAXcBb2nDNgG3tuldbZ62/M6qqla/ol1tdA6wCvgGsB9Y1a5OOrXtY9ecPDpJUl/6eWUwlb8Dbkryj8C3gBtb/UbgM0nGgEN0frhTVQ8kuRl4EDgCXFNVvwBI8m5gD7AI2F5VDxxHX5KkGZpRGFTV14CvtelH6FwJNHnM/wFvnWL9DwEf6lHfDeyeSS+SpLnjJ5AlSYaBJMkwkCRhGEiSMAwkSRgGkiQMA0kShoEkCcNAkoRhIEnCMJAkYRhIkjAMJEkYBpIkDANJEoaBJAnDQJKEYSBJwjCQJGEYSJIwDCRJGAaSJAwDSRKGgSQJw0CSRB9hkORFSb6R5NtJHkjygVY/J8ndScaSfCHJqa3+wjY/1pav7NrW+1v9e0ku6aqvb7WxJFvm/mFKko6ln1cGPwPeUFXnA68C1idZC3wY+GhVvRx4Criqjb8KeKrVP9rGkeRc4ArglcB64JNJFiVZBHwCuBQ4F3h7GytJGpBpw6A6ftpmT2m3At4A3NLqO4DL2/SGNk9bfnGStPpNVfWzqvoBMAasabexqnqkqn4O3NTGSpIGpK/3DNpv8PcCB4G9wPeBp6vqSBsyDixr08uAxwDa8sPAS7vrk9aZqt6rj81JRpOMTkxM9NO6JKkPfYVBVf2iql4FLKfzm/wr5rWrqfvYVlWrq2r1yMjIMFqQpJPSjK4mqqqngbuA1wBLkixui5YDB9r0AWAFQFv+EuDJ7vqkdaaqS5IGpJ+riUaSLGnTpwF/AjxEJxTe0oZtAm5t07vaPG35nVVVrX5Fu9roHGAV8A1gP7CqXZ10Kp03mXfNxYOTJPVn8fRDOBvY0a76eQFwc1XdluRB4KYk/wh8C7ixjb8R+EySMeAQnR/uVNUDSW4GHgSOANdU1S8Akrwb2AMsArZX1QNz9gglSdOaNgyq6j7g1T3qj9B5/2By/f+At06xrQ8BH+pR3w3s7qNfSdI88BPIkiTDQJJkGEiSMAwkSRgGkiQMA0kShoEkCcNAkoRhIEnCMJAkYRhIkjAMJEkYBpIkDANJEoaBJAnDQJKEYSBJwjCQJGEYSJIwDCRJGAaSJAwDSRKGgSQJw0CShGEgScIwkCTRRxgkWZHkriQPJnkgyXta/cwke5M83O7PaPUkuT7JWJL7klzQta1NbfzDSTZ11S9Mcn9b5/okmY8HK0nqrZ9XBkeAv62qc4G1wDVJzgW2AHdU1SrgjjYPcCmwqt02AzdAJzyArcBFwBpg69EAaWPe1bXe+uN/aJKkfk0bBlX1eFV9s03/BHgIWAZsAHa0YTuAy9v0BmBndewDliQ5G7gE2FtVh6rqKWAvsL4tO72q9lVVATu7tiVJGoAZvWeQZCXwauBuYGlVPd4WPQEsbdPLgMe6VhtvtWPVx3vUe+1/c5LRJKMTExMzaV2SdAx9h0GSFwNfBN5bVc90L2u/0dcc9/YcVbWtqlZX1eqRkZH53p0kPW/0FQZJTqETBJ+tqi+18o/aKR7a/cFWPwCs6Fp9easdq768R12SNCD9XE0U4Ebgoar6SNeiXcDRK4I2Abd21Te2q4rWAofb6aQ9wLokZ7Q3jtcBe9qyZ5Ksbfva2LUtSdIALO5jzGuBdwD3J7m31f4euA64OclVwKPA29qy3cBlwBjwLHAlQFUdSvJBYH8bd21VHWrTVwOfBk4Dbm83SdKATBsGVfVfwFTX/V/cY3wB10yxre3A9h71UeC86XqRJM0PP4EsSTIMJEmGgSQJw0CShGEgScIwkCRhGEiSMAwkSRgGkiQMA0kShoEkCcNAkoRhIEnCMJAkYRhIkjAMJEkYBpIkDANJEoaBJAnDQJKEYSBJwjCQJGEYSJIwDCRJGAaSJPoIgyTbkxxM8p2u2plJ9iZ5uN2f0epJcn2SsST3Jbmga51NbfzDSTZ11S9Mcn9b5/okmesHKUk6tn5eGXwaWD+ptgW4o6pWAXe0eYBLgVXtthm4ATrhAWwFLgLWAFuPBkgb866u9SbvS5I0z6YNg6r6OnBoUnkDsKNN7wAu76rvrI59wJIkZwOXAHur6lBVPQXsBda3ZadX1b6qKmBn17YkSQMy2/cMllbV4236CWBpm14GPNY1brzVjlUf71HvKcnmJKNJRicmJmbZuiRpsuN+A7n9Rl9z0Es/+9pWVauravXIyMggdilJzwuzDYMftVM8tPuDrX4AWNE1bnmrHau+vEddkjRAsw2DXcDRK4I2Abd21Te2q4rWAofb6aQ9wLokZ7Q3jtcBe9qyZ5KsbVcRbezaliRpQBZPNyDJ54HXA2clGadzVdB1wM1JrgIeBd7Whu8GLgPGgGeBKwGq6lCSDwL727hrq+rom9JX07li6TTg9naTJA3QtGFQVW+fYtHFPcYWcM0U29kObO9RHwXOm64PSdL88RPIkiTDQJJkGEiSMAwkSRgGkiQMA0kShoEkCcNAkoRhIEnCMJAkYRhIkjAMJEkYBpIkDANJEoaBJAnDQJKEYSBJwjCQJGEYSJIwDCRJGAaSJAwDSRKGgSQJWDzsBtSflVu+MrR9//C6Nw5t35IGwzAYoGH+QD8ex9O3QSKdGAwDzavpgsSwkBaGBRMGSdYDHwcWAZ+qquuG3NKsnKi//Q/LfB4vg0bq34IIgySLgE8AfwKMA/uT7KqqB4fb2XP5w/7Ecax/K4NC+k0LIgyANcBYVT0CkOQmYAMw8DDwh/3zw8n472zA6XgslDBYBjzWNT8OXDR5UJLNwOY2+9Mk35vl/s4CfjzLdeebvc3OQu1tYH3lwzNeZaEeM7C32Zqut9+basFCCYO+VNU2YNvxbifJaFWtnoOW5py9zc5C7W2h9gX2Nlsna28L5UNnB4AVXfPLW02SNAALJQz2A6uSnJPkVOAKYNeQe5Kk540FcZqoqo4keTewh86lpdur6oF53OVxn2qaR/Y2Owu1t4XaF9jbbJ2UvaWq5rIRSdIJaKGcJpIkDZFhIEk6ucMgyfok30sylmRLj+UvTPKFtvzuJCsXSF/vTDKR5N52+4tB9NX2vT3JwSTfmWJ5klzfer8vyQULqLfXJzncddz+YUB9rUhyV5IHkzyQ5D09xgzluPXZ27CO24uSfCPJt1tvH+gxZljP0X56G9rztO1/UZJvJbmtx7KZH7eqOilvdN6I/j7w+8CpwLeBcyeNuRr4lzZ9BfCFBdLXO4F/HtJxex1wAfCdKZZfBtwOBFgL3L2Aens9cNsQjtnZwAVt+reB/+7xbzqU49Znb8M6bgFe3KZPAe4G1k4aM/Dn6Ax6G9rztO3/b4DP9fq3m81xO5lfGfzqKy6q6ufA0a+46LYB2NGmbwEuTpIF0NfQVNXXgUPHGLIB2Fkd+4AlSc5eIL0NRVU9XlXfbNM/AR6i86n6bkM5bn32NhTtWPy0zZ7SbpOvaBnGc7Tf3oYmyXLgjcCnphgy4+N2ModBr6+4mPwk+NWYqjoCHAZeugD6AvjTdjrhliQreiwfln77H5bXtJf2tyd55aB33l6Ov5rOb5Ldhn7cjtEbDOm4tVMd9wIHgb1VNeVxG+BztN/eYHjP048B7wN+OcXyGR+3kzkMTmT/Aaysqj8E9vLrhNexfRP4vao6H/gn4N8HufMkLwa+CLy3qp4Z5L6nM01vQztuVfWLqnoVnW8dWJPkvEHtezp99DaU52mSNwEHq+qeudzuyRwG/XzFxa/GJFkMvAR4cth9VdWTVfWzNvsp4MJ57mkmFuxXh1TVM0df2lfVbuCUJGcNYt9JTqHzw/azVfWlHkOGdtym622Yx62rh6eBu4D1kxYN4znaV29DfJ6+Fnhzkh/SOc38hiT/NmnMjI/byRwG/XzFxS5gU5t+C3BntXdchtnXpHPJb6Zznneh2AVsbFfHrAUOV9Xjw24KIMnvHj0vmmQNnf/f8/6Do+3zRuChqvrIFMOGctz66W2Ix20kyZI2fRqdv2fy3UnDhvEc7au3YT1Pq+r9VbW8qlbS+flxZ1X9+aRhMz5uC+LrKOZDTfEVF0muBUarahedJ8lnkozReWPyigXS118neTNwpPX1zvnu66gkn6dzdclZScaBrXTePKOq/gXYTefKmDHgWeDKBdTbW4C/THIE+F/gikH84KDzm9o7gPvbOWaAvwde1tXbsI5bP70N67idDexI549bvQC4uapuG/ZzdAa9De152svxHje/jkKSdFKfJpIk9ckwkCQZBpIkw0CShGEgScIwkCRhGEiSgP8Hzb9hiRqgKosAAAAASUVORK5CYII=\n"
          },
          "metadata": {
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Checking the law of large numbers\n"
      ],
      "metadata": {
        "id": "u6XdJD2mgeN8"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "_, sum_array_1, nums_of_intervals_1 = get_recov_processes(K = 1,T = 500)\n",
        "N_t = list(range(nums_of_intervals_1[0][0]))\n",
        "t = sum_array_1[0]"
      ],
      "metadata": {
        "id": "isujDWfjUKwQ"
      },
      "execution_count": 93,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "num_avg = []\n",
        "for i in range(len(t)):\n",
        "    num_avg.append(N_t[i]/t[i])\n",
        "\n",
        "M = 1\n",
        "plt.axhline(y = M, color = 'r', linestyle = '-', label = 'Theoretical average')\n",
        "plt.plot(num_avg, label = 'Sample average');\n",
        "plt.xlabel('t');\n",
        "plt.ylabel('N(t)/t');\n",
        "plt.title('The law of large numbers with t = T');\n",
        "plt.legend();"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 295
        },
        "id": "aT9KPrhVhWcf",
        "outputId": "44d82a82-07b5-4458-8122-c4a467c8635f"
      },
      "execution_count": 99,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEWCAYAAABrDZDcAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deXxU5bnA8d8zk32HJOwEwr4vioCAFdxAq9BFC1S9Ur0uBa3a9trqVYvbra1Wq1ZrsSqKFrRaFRErLihWUdkFwg6RhC0LZN9mee8f52SYhCQkkCEk5/l+PvlkzjJnnncyOc+8y3mPGGNQSinlXK6WDkAppVTL0kSglFIOp4lAKaUcThOBUko5nCYCpZRyOE0ESinlcJoIWjERmSsirzTTseaLyIPNcawTfP2OIrJCRIpF5E91bG/R+FqCiPQUESMiYS0dS31EpEREejWwPVNELjiVMammO20/YMr6JwtajAEqAZ+9fOOpjyikbgDygASjF7e0GsaYuOrHIjIfyDbG3H0ixxKRWcB/G2MmNE90geM2+H9kjHm1OV+vNdIawWnMGBNX/QPsBS4LWtfWPrw9gIzmTgJicfzn/HSuVYSaw/6PTojj/0HagAgRedluUtksIqOqN4hIFxF5U0RyRWSPiPyiMQcUkXYissR+3hH7cTd72yQR2Ri074cisipo+XMR+UE9xx0nIqtEpND+Pc5ePx+4BrjDbmposCmhofjs7Z+KyEMi8gVQBvQSkYtEZJv92s+IyGci8t9Bz7lWRLbYx/tARHrU89rVzTXXiMheEckTkf8N2l6jCUtEJopIdtBypoj8j4h8KyKlIvK83Sz2vv03/EhE2tV62WtFZL+IHBCRXwcdyyUivxWRXSKSLyKvi0j7WnFeJyJ7gU9EJEpEXrH3LbD/Bh3rKOPPROTdoOUdIvLPoOUsERlhPzYi0kdEbgCuDPobvht0yBF2eQtF5DURiarjNQcCzwJn288vqOv9V6GhiaD1mwosApKAxcBfwDpJAO8CG4CuwPnAbSIyuRHHdAEvYn1LTwPKq48LfAX0FZEUEQkHhgFdRCReRKKBUcDntQ9on6DeA54EkoHHgPdEJNkYMwt4Ffij/S3to5OIr9rVWM1N8UAh8AZwp/3a24BxQbFNA+4CfgSk2vEvPE4ME4D+WO/rvfaJrLF+DFwI9AMuA963Xz/VLlvthD0J6AtcBPwmKFHeAvwAOBfoAhwBnq713HOBgcBkrGSbCHTHeh9uwnrvavsMOMdONF2ACOBsALH6A+KAb4OfYIyZR82/4WVBm38CTAHSsT4vs2q/oDFmix3PSvv5SXXEhZ3EC+r5+bau56jj00TQ+v3HGLPUGOMDFgDD7fVnAanGmPuNMVXGmN3Ac8CM4x3QGJNvjHnTGFNmjCkGHsI6oWCMKQdWAd8DzsRKNF8A44GxwA5jTH4dh/2+vW2BMcZrjFkIbMU6ETZJQ/EFmW+M2WyM8QIXA5uNMf+yl58EDgbtexPwe2PMFnv7/2F9i62zVmC7zxhTbozZYL8HwxvYt7anjDGHjDH7sJLO18aYdcaYCuAtYGQdr1VqjNmIlQBnBsX9v8aYbGNMJTAXuFxqNgPNtZ9bDniwEkAfY4zPGLPGGFNUOzj7s1IMjMD6O38A7BeRAVjv8+fGGH8TyvukMWa/MeYw1peTEU14bu3YZhtjkur5GXaix3U6x7YbtiHBJ7QyIMo+EfTA+qYeXMV2U8e39dpEJAZ4HOtbXHUzRbyIuO2E8xkwEci2Hx/BOkFU2st16QJ8V2vdd1i1lSZpRHwAWbVeO7BsjDHBzTVY79UTUnO0ktix1Y65Wu33Pa6e/epyKOhxeR3LtY8VXJbvgKH24x7AWyISfFL2AcHNPcHPXYBVG1gkIknAK1iJxFNHjNV/4z724wKsv/HZ1P83rk/t96pLE5+vQkxrBG1XFrCn1jemeGPMJY147q+wmj3GGGMSsL4VgnVyhKMnie/Zjz/DOkmcS/0nif1YJ65gacC+RpanKfEBBHc6HwCC+xAkeBnrvbqx1nsVbYz58gRiK8UamVKt0wkco7buQY/TsN5LsOK+uFbcUXZNo1rgfTDGeIwx9xljBmE1jV0K/Fc9r1n9Nz6Hxv+NT7aj/7jPF5Fn7T6Eun42n+TrO5YmgrbrG6BYRH4jItEi4haRISJyViOeG4/1zbTAbtv/Xa3tX2KdiEcD3xhjNmOd5McAK+o55lKgn4j8VETCRGQ6MAhY0vSiHTe+2t4DhorID+za0hxqnqCfBe4UkcEAIpIoIlecQFwA64FLRKS9iHQCbjvB4wS7R0Ri7Ph+Brxmr38WeKi6CUtEUu3+jjqJ1dE/VETcQBFWU1F9TTyfYfVNRBtjsrFqklOwmpbW1fOcQ0C91xQ0wiGgm4hE1LeDMeam4FFAtX4Gn8RrO5omgjbKbiK5FKs9dg/WGP2/Y3UWHs+fgWj7OV8B/6517FJgLVa7e5W9eiXwnTEmp5548u14fgXkA3cAlxpj8ppWsuPHV8dr5wFXAH+0X3sQsBqrKQtjzFvAH7CaTIqATVj9CidiAVafQSawjKMn7ZPxGbAT+Bh41BizzF7/BNYAgWUiUoz1Xoxp4DidsDrNi4At9nEX1LWjMWY7UILdlGj3JewGvghqfqvteWCQ3XH7duOLF/AJsBk4KCIn8rlQJ0j02h3lNPaIqmzgSmPM8paOR6mWpjUC5QgiMllEkkQkEmuopmB9g1bK8TQRKKc4G9iF1Zx0GfADe0ilUo6nTUNKKeVwWiNQSimHa3UXlKWkpJiePXu2dBhKKdWqrFmzJs8Yk1rXtlaXCHr27Mnq1atbOgyllGpVRKS+q+S1aUgppZxOE4FSSjmcJgKllHI4TQRKKeVwmgiUUsrhNBEopZTDaSJQSimH00TQCJv3F7Lmu8MtHYZSSoVEq7ugrCV8/8n/AJD58PdbOBKllGp+WiNQSimH00RwHMUVdd3XWyml2g5NBA34KOMQQ+cuO/6OSinVimkiaMDHWw+1dAhKKRVymggakFtcxYBO8Vw3IZ3YCHdLh6OUUiGhiaABeSWVpMRFEhnmotLrb+lwlFIqJDQRNCC/tJKUuAgiwlx4/Qa/X2/rCZBTXMFlT/2H8Q9/wpxX1+Kr433JKa5g4Td7qfD4WiBCpVRT6HUEDcgrriIlLpKIMCtfVvn8RLmc20RU4fGx/VAxtyxcx3f5ZcRHhfHexgOUVHrJOlLGlMGd+NVF/cktruTiJ1ZwpMxDWZWP6yakt3ToSqkGaCKoR2mll3KPj5T4SMJcAkCl109UuHMSQU5xBR9mHMItwt//s4edOSUAdE6M4s2fj+OMtCTuemsjC7/JAuCZT3fx0peZeP0GlwgRbhcLVmYSHxVGhcfHjLPSiAhzUVbl5ZevbaB/p3hmjO7ON3sOM3lwJ4oqPHSIj2rBEivlTJoI6pFXUglASlwk5XbzRqXXB4S3YFShlZlXyuMfbWdIl0Sum5DOdfNXs3FfIQB9OsTxk1HdyDpczkM/HEKv1DgA/u+HQ5l+Vho92sfw1e58vt5jTcVx5Zg0thws5hcL13HHG98C8MRHOxjRPYmPt+YA8O/NB3ni4x01YhjaNZGkmHByiyu5dkI6PxnVvcb2kkov976ziX4d47nxe70QkTrLcqS0isJyD0UVHoZ0SeRgUQU5xZW88tV3/OK8vqQlxzTfG6dUK6eJoB7V3357JsewO7cUgKo21mF8sLCCJd/u5821+yip9JB1uByAd9bv56GlWwC4dnw6keEuZk/sTXzUsUlQRBjRPQmAi4d25uKhnQPbeiTHMrZXezolRHHx0M68u2E/Ow6V4HYJcyb1oX1MOCWVXmIjw1i8YT+9U+PYerCILQeKAOGON77lQEEF4/ok0yM5hoOFFTz58U4+2mIN691+sJhN+wspLPcwqkd7/vST4eSXVvHlzjwefG8LheV1Xwz41e583r/1HJZ8e4Bv9hxmxlndGdMruTnfWqVaFTGmdXWAjho1ypyKm9f/5ZMdPLpsO5vum8xHGYe47bX1fPKrcwPfhFu7d9bv49ZF6wPL3dpFYww8Pn0E//PGBr7LL2NwlwTevXkCLlfd37pPVKl98m9IpdfHz19Zyyd27SHYXZcM4K11+9lyoIiUuAgSosPZnVtaY3RXdLibCwZ1xOP18/HWQ5w3oANul/DDkd246ZU1gQ7uMJfg9Rtmjk7jN1P6kxQT0axlVep0ISJrjDGj6tqmNYI6+P2Gz3fkkdY+hrjIMCLtzuLWOoS0uMLDom+ymDayCx3io1iwMpP7l2QQEebilkl9GNs7mVE92gHWN/zP/mcSeSWVRIW7mz0JAMdNAgCRYW6ev2YUew+XsTOnhHV7C+iZEsuATvEM6ZrIrHHpHCmrokN8JCLCO+v38frqLNJTYpnQJ5WLBnUMxF478Sy4djS3v76eXilxPDlzJJc/+yULv9nLur1HeO8X5+AOQZmVOp1pIqjD4g37+XrPYWaN6wlwdNTQaZgIjDFsyC6kqNzDGT3aUVTuISkmnMy8MgZ1ScAYww0vr2Hl7nwe+WAbVT6rDBP6pPD0lWeQGF13n0dKXOSpLEadRIQeybH0SI7l/IEda2yLCHPRMeFox/K0EV2ZNqJrncepnXjG9Unhi9+ch9sliAiL50zgrXXZzH03g6UbD3DZ8C7NXxilTmMhSwQi8gJwKZBjjBlSx3YBngAuAcqAWcaYtaGKpym+3JVHfGQY9146CKDG8NHTzV8+2cmfPtxe57Yrx6Tx6bZc9hWUM6BTPFsPFgPw4zO68cfLhzn6m2+Y++glNIkx4Vx9dk8WrcriloXrePC9DG6e1Ifp9iinhhwprSIpJrzeTmulWoNQ1gjmA38BXq5n+8VAX/tnDPBX+3eLW7e3gFE92wWaFiLDrCGjlZ7TKxG8v/EAj320nQsGdmDWuHQ+2ZqDwbBg5Xd4/YZXv95LSlwE101I544p/fn753sY3i2JCX1TWjr0047bJfzt6jN59rNd7Mwp4Z53NnP/kgy6JkUzOr095w/syEWDOiIiGGNY/d0RHv9wO1/uymdkWhLn9ElhbO9kxvXW97alFJZ7OFRUwfZDxfgNpMRFMCY9mZW78lm0ai8HCyuIiwqjrNLH+QM74PUbCsqq6JQYTdbhMvYXlBMZ7ibcLfTrGM/15/RyzJelkHYWi0hPYEk9NYK/AZ8aYxbay9uAicaYAw0dc1R8vFl95pkhiNZSKW4GjL6dW/at5JfZXwCwPrYTPxh6NS9sfZPzCnaH7LWbIisygR8NvpJOVSX8c/NCoow3sK3EFU6s38O+iARSPaVEGr26tykM8GlSOu+370dWZCJbYjpQEB7N1LwMulYW8VlSOhmxHelYVczA0lzWxHehOMxqplq59q90ripp2QK0UR5xsTM6mZUJ3fmoXR8Gl+Ywrug7DkTEszIhjX+374en1gWfEX4vVa4wkj2l9CnPp8QdSZE7kqyopBrbo31VpFUW4hE3Re5I8iJi6VeWy5nF+ykIi8IrLiYW7KFveR7tvBWklx/GKy6iWtH/lnz22WnZWdwVyApazrbXHZMIROQG4AaAYZGhbbs+GBGPEaFbZWFgXYT9x66Sox+yqUOuolf5Ef68672QxlOX1XFd+dmAH1PpcvPX7YtrJAGAOL81bLJbVdEpj60tEGBSwR4mFewBwIvwfz0m8kJn63+oa2UhD+5exg/yMojze/AiLOowjLt7XcRPBs3kmoNrmZq/lUi/l0Rf5UnFYoDc8FhcxpDiLTvJkrUu+WHRhBk/+yMTeKnjSN5KHUSly+rTSvBW8GViD57rchYAqVUlXJ67iTFFWXSvLCTRW8FXCd3ZEtuBsUV7ufDwrsD/iQEORsQR5feS5K3gcFg08b5KIszRGv87yQN4vvMolrXvQztPOQVh0XzYvm+N+NzGzxnF++lYVUyYMczM2YAB+pQfJsVbhh/YFpNKh6oSovxefCLE+6rYF5FAl6qi02p+n5asESwBHjbG/Mde/hj4jTGmwbGhoR4++tXufGbM+4pXrhsTaELZmVPCBY99xi3n9eGacT3xG8Pohz4GTv3tK3ccKmbGvK9IiA7n5WtH0729Xhh1qhwqqiAq3E18ZFido6k+257LUx/vYPV3RwAQgbHpycREuImJDOPmSX1IjY9kQ3YBS789wKHiSqYN78I3ew7j8flJiY9k6vAudEiI5PPteXyxM4/l23I4UmYl9uHdkxjXO5mcokpiItzceckAYiKs73L7Csrx+029n4eSSi9vrM6iXWwEkwZ0IKGOa0JCqaCsinfW7+fjrTlszC4gKSaC4d0S+eEZ3aj0+Cgo87A+u4CcogoSosLJLihnVeZhqk9PUeEuLh3WhQGd4rlwUEfS2sdwuLSKjfsKSWsfQ3pKbEj7abw+PztySth6sAi/H3bklLA7t4QDhRUUVXjIK66ktMr6whgd7mZI1wQOFFaQfcS6NsftspoUo8LdlFX5SIoJp11MBCWVXnx+Q7d20fRKieWMHu0Y1zuFnTnFZBwoZt3eI4Fj/GhkV34+sXeN/q2mOF2Hj+4Dgi8b7Wava1EHCq03vXPS0REp1cNHn/pkJ1/szONgYcUpiaXS62Pu4gwuHdaZ8X1SWLE9l5/NX0V0uJvn/utMTQKnWPAopbqc2y+Vc/ul8m12Aaszj/BtdgHfZhcS7nZxoLCcdzfsD+wbHxlGpc/Piu25hLmE2Mgwyqq8zFtxtOnR7RKmDu/CsG6JLN+Wy4rtuWzIKiA1PpL8kko27S/kgoEdWZ9VwIcZ1kV2KXERXDioE+1jwyks99AzOZZlmw+xPqsgMNgh3C30To3D5zeMTm/P2b2TCXO58BvDeQM61JhGxe837MotoV1sBAIkx0VijMHrN7hFcLmEsiovFR4/7WOtazCKKjx8uTOPd789QMb+IorKPRSUe/D5Db1TYzmzRztyS6r4MOMQb6+v+Z5EhlsTPHaIj2TOxD54/H46JUTxw5Fdj7nGIzkukon9O5zYH7OJwtwuBnZOYGDnhDq37ysoZ2N2IeFu4aMth9iTV8qATvHMntiHjAPWZyAuMoziCi8J0daV80dKq2gXG05JpY+13x1h7+GyGu+HS6Bfx3j6d4yntMrLnz7cTpjbxc8n9m7+8jX7ERtvMXCziCzC6iQuPF7/wKlwwD7Jd048+k8fHXQvgrV7CwBIT4llT14pheWeeodg1uer3fnc924Gd148gO/1S61znxXbc3nm0518tfswC7/Zy3UT0nlr3T46JUTx0rWj6dOhbVzY1hYN65bEsG5JNdbtLyjn6ue/pndqHBcM7MgFgzriFiEzv5S+HeOIiQgjt7iSf63Nxu0SxqQn06dDXOCzN2tcT/wGyj0+4iLDeGf9Pn775kbW7S2gU0IUN57bizCXsCunlDfXZNcY4dY1KZpZ43tyydDOeH1+Pth8kF321fJvrdvHq1/vDezbr2Mcae1j2ZVbQkFZFV6fobjSalJxCYzonsT2QyWUVnkJd7sY2jWRjfsKqfL6Gdg5gZJKDwcKKvD6De1iwhnXJ4V2MeG0j43k4iGdapxIiyo8rN9bQFJMOInR4XRvFxOS61ZOha5J0XRNigY4ZqhzUyzfmkNOcQUDOyfQt0N8jXPP0o0HmNi/7vPFyQpZ05CILAQmAinAIeB32BP1GGOetYeP/gWYgjV89GfHaxaC0DcN3f32Rt7dcIANv7uoxvqF3+y1/uEWrGHm6O6MTk9mzj/WsvQX5zCoS93fEuoz59W1vLfxAAM6xfPv2753zHaf39D7rqU11rldwplp7XjgB0Po3ym+6QVTbU6V14/X7w80D1Xz+w17D5fhdgkVHh99OsTV22xS4fGxM6cEn9+w/VAxf1m+k3C3iw7xkfRKjeVQUSVnpLUjv6SSMo+PTfsKGdgpgY6JUewvKGdnTgkj05JIio5g9XeHSYqJoHNiFJP6d+CMHkmBEXeq5bVI05AxZuZxthtgTqhe/0QYY/hsey7DuiUes23m6DQAtj4wBZdLWJ9l1Qz2FZTXmQiyDpdxzh+X88fLh3GgoILZk3oT7nbx6bYc3t9kVXy2HixmfVYB/TvG4zOGOPvCpzV2GzPA0z89g037C7nxe710+gNVQ0SYi4g6uhxdLqFnSmyjjhEV7mZIV+vzPrx7ElfUmuRPOYNeWRxky4Fisg6Xc8ukvvXuU1117WG3zy/6Zi+DuyTQxa4WVvvGnoWzeubN0ent2XawiP97fyt9O8Rz+4X9uOmVNfzg6S/o3zGevYfLuGpsGndMGcDb6/cRFe5i9d0XEhcZxveHdUYppUJFE0GQ6o7ifo1oemkXG0FMhJuPt+awfd5KPr/jvBrbd+TUHEv+079/hTHQv2M8f73qDHqlxnH9Oek89/keth0q5vwBHXju8z3sKyjn8+15XDK0c6CGoJRSoXQ6DWVtcdXD9NrFNK7zt8weLpZ1uJzyKh/Lt+Uw68Vv+C6/lIwDRUS4XTz4A2vkrDFw3YR0/n3bOYEZTG8+z6p5uF3C87PO4sZze7F040GKK73894RezV08pZSqk37lDFJQVgXQ6Lb4jgmRHCqyLhj6/pOfk1dSSVGFl3Mf+RSAy8/sxlVje3D325sAuP3CfjU67RKjw3n2qjPo38nqY/j1Rf3pnBBFfFR4kzuglVLqRGkiCHKkrAq3S0iIatzbsvD6scxbsZtFq7LYnWcNxxvXO5kvd+UDMMgeKvfPm86mtNJbZ1PPlCFH2//D3S5mjdf7+yqlTi1tGgpSUGZdE9DYKxR7pcbx+x8NZcO9R4ea3nJeXzrZFx4N6Gz1NZzVs/0pu/BFKaWaShNBkIIyay7/phAREmPCA9/+B3SK5+EfDyWtfcwxFxUppdTpSJuGghwpq6LdCY7VX3DdaNZnFdAuNoKJ/Tuw4g6tASilWgetEQQ5UuZp9Iih2pLjIk/q0nKllGopmgiC5JdUBibOUkopp9BEYPP5DXkllcedYVIppdoaTQS2/JJK/AY6aCJQSjmMJgJb9YVhHeNDewc0pZQ63WgisB0qsu5DoE1DSimn0URgO1SsiUAp5UyaCGwHCipwu4SUOB01pJRyFk0EtqwjZXRJijrhG0MrpVRrpWc9297DZXRvpzeDV0o5jyYCW9bhck0ESilH0kQAlFf5yCuppHv76OPvrJRSbYwmAuCgPXS09n2HlVLKCTQRAAcLdeioUsq5NBEAOYFrCPSqYqWU82gi4OhVxTrPkFLKiTQRYM0zFB3uJr6OeworpVRbp4kAq0bQMSGy0fcqVkqptkQTAZBXUklKnPYPKKWcKaSJQESmiMg2EdkpIr+tY3uaiCwXkXUi8q2IXBLKeOqTX1JFss4xpJRyqJAlAhFxA08DFwODgJkiMqjWbncDrxtjRgIzgGdCFU9DDpdW0T5WawRKKWcKZY1gNLDTGLPbGFMFLAKm1drHAAn240RgfwjjqZPfbzhSVqWzjiqlHCuUiaArkBW0nG2vCzYXuEpEsoGlwC11HUhEbhCR1SKyOjc3t1mDLCj34DfoTeuVUo7V0p3FM4H5xphuwCXAAhE5JiZjzDxjzChjzKjU1NRmDeBwqXWLSk0ESimnCmUi2Ad0D1ruZq8Ldh3wOoAxZiUQBaSEMKZj5JVUAeioIaWUY4UyEawC+opIuohEYHUGL661z17gfAARGYiVCJq37ec4jpRaiaBdjNYIlFLOFLJEYIzxAjcDHwBbsEYHbRaR+0Vkqr3br4DrRWQDsBCYZYwxoYqpLoXlHgASY8JP5csqpdRpI6RzKhhjlmJ1AgevuzfocQYwPpQxHE8gEURrIlBKOVNLdxa3uMJyD26XEBvhbulQlFKqRWgiKPeQGB2u8wwppRzL8YmgqMKrzUJKKUdz7LzLHp+fh97bws6cEhKiHPs2KKWUc2sEK7bnMv/LTLYcKCJBawRKKQdzbCLwBw1S1aYhpZSTOTgRHM0EmgiUUk7m2ERQUuENPNbpJZRSTubYRFBc4Qk8To3XRKCUci7HJoKioBqBJgKllJM5NhFojUAppSyOTQRF5UE1Au0jUEo5mGMTQXGl1giUUgocnAiqawSJ0eFEheuEc0op53JsIiip9HJO3xQ2/O6ilg5FKaValGMTQYXHR4xOPa2UUs5NBOUeH9HaJKSUUg5OBFU+orVGoJRSDk4EHp92EiulFA5OBBXaNKSUUoBDE4HH58fjM5oIlFIKhyaCCo8PQPsIlFIKhyaC8iorEWgfgVJKOTURVNcINBEopZTDE4E2DSmllEMTQZXWCJRSqpozE4FH+wiUUqpaSBOBiEwRkW0islNEflvPPj8RkQwR2Swi/whlPNV01JBSSh0V1tgdRSTdGLPneOuCtrmBp4ELgWxglYgsNsZkBO3TF7gTGG+MOSIiHU6kEE1VXuUHtGlIKaWgaTWCN+tY90YD+48GdhpjdhtjqoBFwLRa+1wPPG2MOQJgjMlpQjwn7GjTkCNbxpRSqobj1ghEZAAwGEgUkR8FbUoAohp4alcgK2g5GxhTa59+9mt8AbiBucaYf9cRww3ADQBpaWnHC/m4Suz7FcdFNrpCpJRSbVZjzoT9gUuBJOCyoPXFWN/oT/b1+wITgW7AChEZaowpCN7JGDMPmAcwatQoc5KvSXGFdXey+Kjwkz2UUkq1eo1JBDHAr4F5xpiVTTj2PqB70HI3e12wbOBrY4wH2CMi27ESw6omvE6TFVd6iQxzERGmTUNKKdWYM2Ea8E/gjyIyV0TGiIg04nmrgL4iki4iEcAMYHGtfd7Gqg0gIilYTUW7Gxv8iSqu8GhtQCmlbMdNBMaYPxhjzgMuATYA1wJrReQfIvJfItKxnud5gZuBD4AtwOvGmM0icr+ITLV3+wDIF5EMYDnwP8aY/JMvVsOKK7wkRGn/gFJKQROGjxpjioG37B9EZBBwMfAyMLme5ywFltZad2/QYwP80v45ZYorvMRrIlBKKaBp1xG4gOFAF6Ac2GSM+RPwpxDFFjLaNKSUUkc1Zvhob+A3wAXADiAXa9hoPxEpA/4GvGSM8Ycy0OZUXOGlQ3xDI1+VUso5GlMjeBD4K3Cj3ZQTYPcPzASuBl5q/vBCQ5uGlFLqqOOeDY0xMwFEJA7rKHoAABi+SURBVBKorLW5wBjz51AEFkraNKSUUkc1ZSB9XdcQNOW6gtPCqszDlFb56JgQ2dKhKKXUaaExfQSdsKaLiBaRkUD1NQQJWBebtSr/XJ1FYnQ4V47t0dKhKKXUaaExDeWTgVlYVwb/iaOJoAi4KzRhhU65x09ybITOM6SUUrbG9BG8BLwkIj82xtQ1A2mrUunx6dQSSikV5LhnRBG5SkSkviQgIr1FZELzhxYalV4/kXofAqWUCmhM+0gysF5E1gBrOHodQR/gXCAPqPPuY6ejSq+PSLfWCJRSqlpjmoaeEJG/AOcB44FhWFcWbwGuNsbsDW2IzavS69f+AaWUCtKoM6Ixxgd8aP+0alVeP5GxWiNQSqlqjRk+em8Dm40x5oFmjCfkKr1+IsO0j0Appao1pkZQWse6WOA6rP6DVpYIfETqqCGllApoTB9BYHZREYkHbgV+hnUz+lY382ilx0+k3rReKaUCGtVHICLtse4ZcCXW5HJnGGOOhDKwUKn0+onQUUNKKRXQmD6CR4AfYd08fqgxpiTkUYVQpden1xEopVSQxnw1/hXWzWjuBvaLSJH9UywiRaENr3kZY6xRQ9pHoJRSAY3pI2gzZ02v3+A3aCJQSqkgjjojVnqtm6jp8FGllDrKWYnA4wPQSeeUUiqIo86IR2sEjiq2Uko1yFFnxEAi0OsIlFIqwFFnxCrtI1BKqWM4KhFUeq0+Am0aUkqpoxx1RqxuGgrXK4uVUirAUWdEj50IdNSQUkodFdIzoohMEZFtIrJTROq9i5mI/FhEjIiMCmU8VT6tESilVG0hOyOKiBt4GrgYGATMFJFBdexXPaPp16GKpZrHZwDtI1BKqWChPCOOBnYaY3YbY6qwpq2eVsd+DwB/ACpCGAsAHq0RKKXUMUJ5RuwKZAUtZ9vrAkTkDKC7Mea9hg4kIjeIyGoRWZ2bm3vCAR1NBHLCx1BKqbamxb4ai4gLeAxrdtMGGWPmGWNGGWNGpaamnvBrVumoIaWUOkYoz4j7gO5By93sddXigSHApyKSCYwFFoeyw7i6j0BHDSml1FGhPCOuAvqKSLqIRAAzgMXVG40xhcaYFGNMT2NMT+ArYKoxZnWoAtI+AqWUOlbIzojGGC9wM/ABsAV43RizWUTuF5GpoXrdhmgfgVJKHatR9yw+UcaYpcDSWuvurWffiaGMBfQ6AqWUqoujzoger9VHoIlAKaWOctQZ0ePz43YJbpc2DSmlVDXHJQLtH1BKqZoclQiqfH5tFlJKqVocdVb0+PxEaCJQSqkaHHVW9HiN1giUUqoWR50VPT4/4WHaR6CUUsEclQi0j0AppY7lqLOi9hEopdSxHHVW9Pi0j0AppWpz1FlRryNQSqljOSoRVHm1j0AppWpz1FnR4/PrvQiUUqoWR50VtY9AKaWO5aizovYRKKXUsRyVCKp8fsK0RqCUUjU46qzo8xvCdApqpZSqwXGJwC2aCJRSKpijEoHfb3BpjUAppWpwVCLwGa0RKKVUbc5KBH60RqCUUrU4LBH4tbNYKaVqcVgiMHrjeqWUqsVRicBvwKV9BEopVYOjEoFVI2jpKJRS6vTiqNOiz+jwUaWUqs1RicCvF5QppdQxQpoIRGSKiGwTkZ0i8ts6tv9SRDJE5FsR+VhEeoQyHp/RzmKllKotZIlARNzA08DFwCBgpogMqrXbOmCUMWYY8Abwx1DFY4zBaGexUkodI5Q1gtHATmPMbmNMFbAImBa8gzFmuTGmzF78CugWqmB8fgOgNQKllKollImgK5AVtJxtr6vPdcD7dW0QkRtEZLWIrM7NzT2hYLyaCJRSqk6nRWexiFwFjAIeqWu7MWaeMWaUMWZUamrqCb2G32giUEqpuoSF8Nj7gO5By93sdTWIyAXA/wLnGmMqQxVMoGlI+wiUUqqGUNYIVgF9RSRdRCKAGcDi4B1EZCTwN2CqMSYnhLHg91u/9ToCpZSqKWSJwBjjBW4GPgC2AK8bYzaLyP0iMtXe7REgDviniKwXkcX1HO6k+aqbhjQPKKVUDaFsGsIYsxRYWmvdvUGPLwjl6wfTUUNKKVW306Kz+FSo7izWpiGllKrJMYlAO4uVUqpuIW0aOp1UJwKtEai2wOPxkJ2dTUVFRUuHok4zUVFRdOvWjfDw8EY/x3GJQO9QptqC7Oxs4uPj6dmzJ6K1XGUzxpCfn092djbp6emNfp5zmob0gjLVhlRUVJCcnKxJQNUgIiQnJze5puiYROCvbhrSfxzVRmgSUHU5kc+FYxKB1giUUqpuzkkEWiNQqtnk5+czYsQIRowYQadOnejatSsjRowgKSmJQYNqzzYfWm+//TYZGRmB5XvvvZePPvqoycfJzMxkyJAhzRlaq+GYRFA9xYTWCJQ6ecnJyaxfv57169dz0003cfvttweWXa7mP614vd56t9VOBPfffz8XXHDKrlU9IcYY/NUnpdOAc0YNBZqGWjgQpZrbbbfB+vXNe8wRI+DPfz6hp/p8Pq6//nq+/PJLunbtyjvvvEN0dDS7du1izpw55ObmEhMTw3PPPceAAQPIzMzk2muvJS8vj9TUVF588UXS0tKYNWsWUVFRrFu3jvHjxzNnzpxjnn/48GEWL17MZ599xoMPPsibb77JAw88wKWXXsrll1/OqlWruPXWWyktLSUyMpKPP/6Y/Px8rr76akpLSwH4y1/+wrhx4+otT0lJCdOmTePIkSN4PB4efPBBpk2bxm9/+1u6d+/OnDlzAJg7dy5xcXH8+te/5pFHHuH111+nsrKSH/7wh9x3331kZmYyefJkxowZw5o1a1i6dCkPP/wwq1atory8nMsvv5z77rsPgKVLl/LLX/6S2NhYxo8fz+7du1myZAmlpaXccsstbNq0CY/Hw9y5c5k2bVq9sTeWcxKBNg0pdUrs2LGDhQsX8txzz/GTn/yEN998k6uuuoobbriBZ599lr59+/L1118ze/ZsPvnkE2655RauueYarrnmGl544QV+8Ytf8PbbbwPWMNkvv/wSt9vN+eefX+fzp06dGjjxB6uqqmL69Om89tprnHXWWRQVFREdHU2HDh348MMPiYqKYseOHcycOZPVq1fXW56oqCjeeustEhISyMvLY+zYsUydOpXp06dz2223BRLB66+/zgcffMCyZcvYsWMH33zzDcYYpk6dyooVK0hLS2PHjh289NJLjB07FoCHHnqI9u3b4/P5OP/88/n222/p168fN954IytWrCA9PZ2ZM2cGYnnooYc477zzeOGFFygoKGD06NFccMEFxMbGntTfzDGJQO9HoNqsE/zmHirp6emMGDECgDPPPJPMzExKSkr48ssvueKKKwL7VVZas86vXLmSf/3rXwBcffXV3HHHHYF9rrjiCtxud4PPr8+2bdvo3LkzZ511FgAJCQkAlJaWcvPNN7N+/Xrcbjfbt29v8DjGGO666y5WrFiBy+Vi3759HDp0iJEjR5KTk8P+/fvJzc2lXbt2dO/enSeeeIJly5YxcuRIwKpR7Nixg7S0NHr06BFIAmAlj3nz5uH1ejlw4AAZGRn4/X569eoVuA5g5syZzJs3D4Bly5axePFiHn30UcAaRrx3714GDhzYYBmOxzGJwOvTRKDUqRAZGRl47Ha7KS8vx+/3k5SUxPomNmFVf9M90efX5fHHH6djx45s2LABv99PVFRUg/u/+uqr5ObmsmbNGsLDw+nZs2dgnP4VV1zBG2+8wcGDB5k+fTpgJY4777yTG2+8scZxMjMza3xz37NnD48++iirVq2iXbt2zJo167jj/40xvPnmm/Tv3/9Eil4vx7SYB2oE2jSk1CmXkJBAeno6//znPwHrhLZhwwYAxo0bx6JFiwDrpHvOOec06fnx8fEUFxcf85z+/ftz4MABVq1aBUBxcTFer5fCwkI6d+6My+ViwYIF+Hy+BmMvLCykQ4cOhIeHs3z5cr777rvAtunTp7No0SLeeOONQG1l8uTJvPDCC5SUlACwb98+cnKOvd1KUVERsbGxJCYmcujQId5///1A3Lt37yYzMxOA1157LfCcyZMn89RTT2Hs89m6desajL2xHJMIdBpqpVrWq6++yvPPP8/w4cMZPHgw77zzDgBPPfUUL774IsOGDWPBggU88cQTTXr+jBkzeOSRRxg5ciS7du0K7B8REcFrr73GLbfcwvDhw7nwwgupqKhg9uzZvPTSSwwfPpytW7cet339yiuvZPXq1QwdOpSXX36ZAQMGBLYNHjyY4uJiunbtSufOnQG46KKL+OlPf8rZZ5/N0KFDufzyy+tMVMOHD2fkyJEMGDCAn/70p4wfPx6A6OhonnnmGaZMmcKZZ55JfHw8iYmJANxzzz14PB6GDRvG4MGDueeeexr79jdIqjNLazFq1CjTUMdOfZZvy+FnL67iX7PHcUZauxBEptSps2XLlpNuF1anr5KSEuLi4jDGMGfOHPr27cvtt9/e6OfX9fkQkTXGmFF17e+YGoFfp6FWSrUSzz33HCNGjGDw4MEUFhYe09/Q3BzTWaxNQ0qp1uL2229vUg3gZDmnRmD0OgKllKqLYxKBT6eYUEqpOjknEegUE0opVSfHnBYDncUhmBBLKaVaM8ecFb06akipZvXQQw8xePBghg0bxogRI/j6669D+noTJ05scE4gdeIcM2oocIcyx6Q+pUJn5cqVLFmyhLVr1xIZGUleXh5VVVUtHVbIeb1ewsLa3mmz7ZWoHnqHMtVW3ffuZjL2FzXrMQd1SeB3lw2ud/uBAwdISUkJzCuUkpIS2Hb//ffz7rvvUl5ezrhx4/jb3/6GiDBx4kRGjhzJ559/TmlpKS+//DK///3v2bhxI9OnT+fBBx8kMzMzcEXt2rVrGTx4MC+//DIxMTE1Xn/ZsmX87ne/o7Kykt69e/Piiy8SFxdXY5/nnnuOefPmUVVVRZ8+fViwYEHgqtw9e/bgcrkoLS1lwIAB7N69m71799Y5TXbt6bBnzJjBrbfeSkVFBdHR0bz44ov079+fsrIyZs2axaZNm+jfvz/79+/n6aefZtSoUY2KtyU55vuxT5uGlGo2F110EVlZWfTr14/Zs2fz2WefBbbdfPPNrFq1ik2bNlFeXs6SJUsC2yIiIli9ejU33XQT06ZN4+mnn2bTpk3Mnz+f/Px8wJo1dPbs2WzZsoWEhASeeeaZGq+dl5fHgw8+yEcffcTatWsZNWoUjz322DEx/uhHP2LVqlVs2LCBgQMH8vzzz5OYmMiIESMC8S5ZsoTJkycTHh7ODTfcwFNPPcWaNWt49NFHmT17duBY1dNhP/bYYwwYMIDPP/+cdevWcf/993PXXXcB8Mwzz9CuXTsyMjJ44IEHWLNmTZPibUmOqREEriPQGoFqYxr65h4qcXFxrFmzhs8//5zly5czffp0Hn74YWbNmsXy5cv54x//SFlZGYcPH2bw4MFcdtllAEydOhWAoUOHMnjw4MD8PL169SIrK4ukpCS6d+8emHfnqquu4sknn+TXv/514LW/+uorMjIyAvtUVVVx9tlnHxPjpk2buPvuuykoKKCkpITJkycDBO5RMGnSJBYtWsTs2bOPO8119XTYYE1Cd80117Bjxw5EBI/HA8B//vMfbr31VgCGDBnCsGHDmhRvSwppIhCRKcATgBv4uzHm4VrbI4GXgTOBfGC6MSYzFLFojUCp5uV2u5k4cSITJ05k6NChvPTSS8yYMYPZs2ezevVqunfvzty5c2tMrVzdlORyuWpMV+1yuQK3o5Ra/6O1l40xXHjhhSxcuLDB+GbNmsXbb7/N8OHDmT9/Pp9++ilgJaO77rqLw4cPs2bNGs477zxKS0sbnOY6eGK6e+65h0mTJvHWW2+RmZnJxIkTG4yjsfG2pJA1DYmIG3gauBgYBMwUkdp3tb4OOGKM6QM8DvwhVPEE7lCmNQKlTtq2bdvYsWNHYHn9+vX06NEjcNJPSUmhpKSEN954o8nH3rt3LytXrgTgH//4BxMmTKixfezYsXzxxRfs3LkTsG40U9fNZYqLi+ncuTMej4dXX301sD4uLo6zzjqLW2+9lUsvvRS3293gNNe1FRYW0rVrVwDmz58fWD9+/Hhef/11ADIyMti4cWOT4m1JoewjGA3sNMbsNsZUAYuA2jfXnAa8ZD9+Azhfaqf/ZqJ3KFOq+ZSUlHDNNdcwaNAghg0bRkZGBnPnziUpKYnrr7+eIUOGMHny5MDdwZqif//+PP300wwcOJAjR47w85//vMb21NRU5s+fz8yZMxk2bBhnn302W7duPeY4DzzwAGPGjGH8+PE1po4Gq3nolVdeCdxMBuqf5rq2O+64gzvvvJORI0cGajEAs2fPJjc3l0GDBnH33XczePBgEhMTGx1vSwrZNNQicjkwxRjz3/by1cAYY8zNQftssvfJtpd32fvk1TrWDcANAGlpaWcG3xiisT7MOMRb67J5fPoIIsPcJ1ospU4LbXUa6szMTC699FI2bdrU0qE0mc/nw+PxEBUVxa5du7jgggvYtm0bERERpzyWpk5D3So6i40x84B5YN2P4ESOceGgjlw4qGOzxqWUUtXKysqYNGkSHo8HYwzPPPNMiySBExHKRLAP6B603M1eV9c+2SISBiRidRorpRyoZ8+erbI2ANYtM1vrlc+h7CNYBfQVkXQRiQBmAItr7bMYuMZ+fDnwiWltt0xTqoXov4qqy4l8LkKWCIwxXuBm4ANgC/C6MWaziNwvIlPt3Z4HkkVkJ/BL4LehikeptiQqKor8/HxNBqoGYwz5+flERUU16XmOuWexUm2Jx+MhOzu7xhh9pcD6ktCtWzfCw8NrrG/1ncVKqZrCw8NJT09v6TBUG+GYuYaUUkrVTROBUko5nCYCpZRyuFbXWSwiuUDTLy22pAB5x92r7XFiuZ1YZnBmuZ1YZmh6uXsYY1Lr2tDqEsHJEJHV9fWat2VOLLcTywzOLLcTywzNW25tGlJKKYfTRKCUUg7ntEQwr6UDaCFOLLcTywzOLLcTywzNWG5H9REopZQ6ltNqBEoppWrRRKCUUg7nmEQgIlNEZJuI7BSRNjPLqYi8ICI59t3eqte1F5EPRWSH/budvV5E5En7PfhWRM5ouchPjoh0F5HlIpIhIptF5FZ7fZstu4hEicg3IrLBLvN99vp0EfnaLttr9rTviEikvbzT3t6zJeM/GSLiFpF1IrLEXnZCmTNFZKOIrBeR1fa6kHy+HZEIRMQNPA1cDAwCZorIoJaNqtnMB6bUWvdb4GNjTF/gY45O730x0Nf+uQH46ymKMRS8wK+MMYOAscAc+2/alsteCZxnjBkOjACmiMhY4A/A48aYPsAR4Dp7/+uAI/b6x+39Wqtbsaazr+aEMgNMMsaMCLpeIDSfb2NMm/8BzgY+CFq+E7izpeNqxvL1BDYFLW8DOtuPOwPb7Md/A2bWtV9r/wHeAS50StmBGGAtMAbr6tIwe33gs451L5Cz7cdh9n7S0rGfQFm72Se984AlgLT1MtvxZwIptdaF5PPtiBoB0BXIClrOtte1VR2NMQfsxweB6ps1t8n3wa7+jwS+po2X3W4iWQ/kAB8Cu4ACY90ICmqWK1Bme3shkHxqI24WfwbuAPz2cjJtv8wABlgmImtE5AZ7XUg+33o/gjbOGGNEpM2OERaROOBN4DZjTJGIBLa1xbIbY3zACBFJAt4CBrRwSCElIpcCOcaYNSIysaXjOcUmGGP2iUgH4EMR2Rq8sTk/306pEewDugctd7PXtVWHRKQzgP07x17fpt4HEQnHSgKvGmP+Za92RNmNMQXAcqxmkSQRqf5SF1yuQJnt7YlA/ikO9WSNB6aKSCawCKt56AnadpkBMMbss3/nYCX90YTo8+2URLAK6GuPNIgAZgCLWzimUFoMXGM/vgar/bx6/X/ZIwzGAoVB1cxWRayv/s8DW4wxjwVtarNlF5FUuyaAiERj9YlswUoIl9u71S5z9XtxOfCJsRuQWwtjzJ3GmG7GmJ5Y/7efGGOupA2XGUBEYkUkvvoxcBGwiVB9vlu6Q+QUdrxcAmzHalP935aOpxnLtRA4AHiw2gWvw2oT/RjYAXwEtLf3FazRU7uAjcColo7/JMo9AasN9Vtgvf1zSVsuOzAMWGeXeRNwr72+F/ANsBP4JxBpr4+yl3fa23u1dBlOsvwTgSVOKLNdvg32z+bqc1aoPt86xYRSSjmcU5qGlFJK1UMTgVJKOZwmAqWUcjhNBEop5XCaCJRSyuE0ESjVDEQkSURmt3QcSp0ITQRKNY8kQBOBapU0ESjVPB4Gettzxz/S0sEo1RR6QZlSzcCeAXWJMWZIC4eiVJNpjUAppRxOE4FSSjmcJgKlmkcxEN/SQSh1IjQRKNUMjDH5wBciskk7i1Vro53FSinlcFojUEoph9NEoJRSDqeJQCmlHE4TgVJKOZwmAqWUcjhNBEop5XCaCJRSyuH+HzZlqFnZmVmRAAAAAElFTkSuQmCC\n"
          },
          "metadata": {
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "N_t = list(range(nums_of_intervals_1[0][0]))"
      ],
      "metadata": {
        "id": "p52tBYr9gU60"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "t_2 = [time for time in sum_array_1[0] if time < T/2]"
      ],
      "metadata": {
        "id": "6UgZUjxfeZfR"
      },
      "execution_count": 96,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "num_avg = []\n",
        "for i in range(len(t_2)):\n",
        "    num_avg.append(i/t_2[i])\n",
        "\n",
        "M = 1\n",
        "plt.axhline(y = M, color = 'r', linestyle = '-', label = 'Theoretical average')\n",
        "plt.plot(num_avg, label = 'Sample average');\n",
        "plt.xlabel('t');\n",
        "plt.ylabel('N(t)/t');\n",
        "plt.title('The law of large numbers with t = T/2');\n",
        "plt.legend();"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 295
        },
        "id": "qQrmJ5zfgx9_",
        "outputId": "f5fe0087-3db1-42e2-c4d2-1d184c202707"
      },
      "execution_count": 102,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEWCAYAAABrDZDcAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nO3deXxU1fn48c8zk2VCdpKwQwj7vsvqgguiVqF1Q1wKtdVa1K+19etXrVrcfrVqbdViFasC1orUFSlWXBCpohL2fQ+QECBk32c7vz/uJE5CEhLIMCTzvF+vvDJzt3nO3OQ+95xz77lijEEppVTosgU7AKWUUsGliUAppUKcJgKllApxmgiUUirEaSJQSqkQp4lAKaVCnCaCFkBEZovIP5ppW/NE5PHm2NZJfn57EflKRIpF5E91zA9qfMEgIt1FxIhIWLBjqY+IlIhIjwbmZ4jIRaczJtV8NBGcAXz/ZFU/XhEp93t/Q7Dja2a3AseAOGPMb4MdjGocY0yMMWYvnHqyFpGZIvLf5ouueruN+j8SkU9E5GIRmSEia0SkSEQyReSpMzkZB5ImgjOA758sxhgTAxwArvCb9maw42tmqcBW08x3Mool5P+eQ/VABo37PxKRaGAUsAJoA/waSAbGABcC9wQn+uAK+X+cFiRCRBb4mlS2iMioqhki0klE3hWRHBHZJyL/05gNikiiiCzxrZfve93FN+98Ednkt+ynIrLa7/1KEflxPdsdLyKrRaTQ93u8b/o8YAZwr+8srcGmhIbi883/UkSeEJGvgTKgh+9Mb4fvs18UkRUi8gu/dW4WkW2+7X0iIqn1fHZVc80METkgIsdE5Hd+82ucFYvIRBHJ9HufISL/KyIbRaRURF71NYt97NuHn4lIYq2PvVlEDolItojc47ctm4jcJyJ7RCRXRBaJSNtacf5cRA4AX4iIQ0T+4Vu2wLcP2tdRxp+JyEd+73eJyL/83h8UkWG+10ZEeonIrcANfvvwI79NDvOVt1BE3hYRRx2f2R94CRjnW7+gru8/gC4EvjbGVBpj/maMWWmMcRpjsoA3gQmnOZ4zgiaClmMKsBBIABYDfwXrIAF8BGwAOmP9of9aRCY3Yps24HWss/RuQHnVdoFvgd4ikiwi4cAQoJOIxIpIFNZZ1craG/QdoP4NPA8kAc8C/xaRJGPMTKx/tqd8Z2mfnUJ8VW7Cam6KBQqBd4D7fZ+9AxjvF9tU4AHgSiDFF/9bJ4jhbKAv1vf6sO9A1lhXAZOAPsAVwMe+z0/xla12wj4f6A1cDPyfX6K8E/gxcB7QCcgH5tRa9zygPzAZK9nGA12xvofbsL672lYA5/gSTScgAhgHIFZ/QAyw0X8FY8xcau7DK/xmXwtcAqRh/b3MrP2BxphtvnhW+dZPqCMufEm8oJ6fjXWt00iXYf191uVcYMspbLvF0kTQcvzXGLPUGOMB3gCG+qafBaQYYx71ndnsBV4BrjvRBo0xucaYd40xZcaYYuAJrAMKxphyYDXWP8dIrETzNdYZ01hglzEmt47N/sg37w1jjNsY8xawHetA2CQNxednnjFmizHGDVwKbDHGvOd7/zxw2G/Z24A/GGO2+eb/P6yz2DprBT6PGGPKjTEbfN/B0AaWre0FY8wR39nmSuA7Y8w6Y0wF8D4wvI7PKjXGbMJKgNP94v6dMSbTGFMJzAaulprNQLN965YDLqwE0MsY4zHGrDHGFNUOzve3UgwMw9rPnwCHRKQf1ve80hjjbUJ5nzfGHDLG5GGdnAxrwrq1Y5tljEmo52fIyW4XKxEsrT1RRG7GOrl55hS23WKFbHtiC+R/QCsDHL4DQSrWmbp/FdtOHWfrtYlIG+DPWGdxVc0UsSJi9yWcFcBEINP3Oh/rAFHpe1+XTsD+WtP2Y9VWmqQR8QEcrPXZ1e+NMca/uQbru3pOal6tJL7Yasdcpfb3HtOEIhzxe11ex/va2/Ivy35gsO91KvC+iPgflD2Af3OP/7pvYNUGFopIAvAPrETiqiPGqn3cy/e6AGsfj6P+fVyf2t9VpyauH1AiMhgoNMYcrDX9x8AfgIuMMceCElyQaY2g5TsI7Kt1xhRrjLmsEev+FqvZY4wxJg7rrBCsgyP8cJA41/d6BdZB4jzqP0gcwjpw+esGZDWyPE2JD8C/0zkb8O9DEP/3WN/VL2t9V1HGmG9OIrZSrM7GKh1OYhu1dfV73Q3ruwQr7ktrxe3w1TSqVH8PxhiXMeYRY8wArKaxy4Gf1vOZVfv4HBq/j0+1o/+E64vIS1LzKiD/n5NtvjmuNiAil2DVoK/w1cRCkiaClu97oFhE/k9EokTELiKDROSsRqwbi3VmWuBr2/99rfnfYB2IRwPfG2O2YB3kxwBf1bPNpUAfEbleRMJEZBowAFjS9KKdML7a/g0MFpEf+2pLt1PzAP0ScL+IDAQQkXgRueYk4gJYD1wmIm1FpAPW1Sen6iERaeOL72fA277pLwFPVDVhiUiKr7+jTmJ19A8WETtQhNVUVF8TzwqsvokoY0wmVk3yEqympXX1rHMEqPeegkY4AnQRkYj6FjDG3OZ/FVCtn4En+bk1+gdE5AKs/o6rjDHfn+Q2WwVNBC2cr4nkcqz22H1Y1+j/Hauz8ET+AkT51vkW+E+tbZcCa7Ha3Z2+yauA/caYo/XEk+uL57dALnAvcPlJVrkbjK+Ozz4GXAM85fvsAUA6VlMWxpj3gT9iNZkUAZux+hVOxhtYfQYZwDJ+OGifihXAbuBz4BljzDLf9OewLhBYJiLFWN/FmAa20wGr07wI2Obb7ht1LWiM2QmU4GtK9PUl7MW6ssZT1zrAq8AAX8ftB40vXrUvsDplD4vIaWmK8TWRDcA6uanyENb/yVK/2sbHpyOeM4008+XcSp0xfFdUZQI3GGOWBzseFTwici1wtTHm2mDHcibSGoFqVURksogkiEgk1qWagnUGrUJbAdaFB6oOetWQam3GAf/EuiZ+K/Bj3yWVKoT5NbOpOmjTkFJKhThtGlJKqRDX4pqGkpOTTffu3YMdhlJKtShr1qw5ZoxJqWtei0sE3bt3Jz09PdhhKKVUiyIi9d09r01DSikV6jQRKKVUiNNEoJRSIU4TgVJKhThNBEopFeI0ESilVIjTRKCUUiFOE0GA/WdzNocLK4IdhlJK1UsTQQBlFZRz2z/W8r/vbAh2KEopVS9NBAH0+TbrEbXWExOVUurMpIkggD7bZj3EKyUmMsiRKKVU/VrcWENnujdWZTDvmwwiwuzsOVoCQHGFK7hBKaVUA7RG0MzeTj9IudNDnCMMgyEyzEZJpTvYYSmlVL20RtCM3B4vu46U8NNxqfzuRwMwxnDzvNUcK3GeeGWllAoSrRE0o4zcMirdXvp2iAOsTuJYR7g2DSmlzmhaI2hGOw4XA9CvQ2z1tBhHmDYNqToZY1i56xgvfrmb7YeL6dMuluJKN9NHd+WmsamNvtqswuVh1d5csgsquHZUF8Lsen6nmkYTQTPafrgIu03o1S6melqsI4yiCk0ECkoq3dhF8BjDglUZfLQhm23ZRXSIc3Bhv/Zk5JYSbhce/nALy7YcYcb47uzPLWXBqv30TInmyhFduHRQh+oDvcvj5YvtR/l/S7exP7cMgJ1Hipk9ZWAQS6laIk0EzWj74WLSkqNxhNurp8U5wnG6vVS6PUSG2RtYW7V0xhg+3XqE7MIKBnSKIyk6ggWr9rNiZw75ZU4Kylw4wm3ERIZzrKSSoV0T+MOVg7lyROfqvw2v1/DGt/v5y2c7uWWB9SS+kamJ7DxSwp1vrUMEYiLDiI4Io6DcSYXLS/ekNsy9aSSr9uby+tcZ9EyJ5qZx3VmxM4cVO3I4p08yI1MTiXOEHxfzriPFPPf5LsqdHn73o/6kJUdTWO7iu315AJzbO4WoCDtZBeVEhtlIio7Q+2JaIU0EzcQYw9ZDRQzrllBjekyk9RWXVLiJjNFE0JK5PV6WbMzmky2HsdmE/7mgN319zYBuj5eHPtzMW98frLFORJiNc3ol0ykhmQ7xDrIKysnML+fXF/VmRLfE4z7DZhNmjO/O1SO7sC27iDC7jaFd4jEGPt9+lE1ZhRSVuyitdBMXFc64Hkmc2yeFiDAbF/Zvz8G8Mn6/eAvvr8ti7YECROC1r/cB0K1tG1KT2hAVbsduE9rHOXjzu/1EhtmxCVzwpxWE2wWXx1THExVup3f7GDZmFgLQKd7BwM7xhNuFttERrN6XjwiM75nMbef1oF2co1Hf5easQv725R6iIuzcODaVYV0TTrxSPQrLXbg8XuwirDuYT6+UWLoltQEgt6SS4go33ZOj2XeslLe+P4DT7WX66G7V+06BGGNOvNQZZNSoUeZMfGbxhoMFTJ3zNU9eOZjrRnernv7e2kx+s2gDX94zke7J0UGMUJ3I4cIK4qPC+WB9FqWVbjxew6L0g4gIgzrFse5gAftzy+icEEVJpRsR6BQfRUxkGJHhNlbuOsasiT2ZOaE7n2w5Qm5JJTeMSSUl9vTdUFha6eauhevIKXEyqX87ZozvztoDBWzOKmTLoUIOFVRQ7vRQ4fawP7eMSQPa8+SVg3F5DB+uzyKvzElSdARDuyTg9ho+3pzNxsxCJvVvT5vIMFbvy2PfsVJcHi+HiyoY2iWB8DAb3+w+RphduGFMKoM6x5HYJoLENhF0jHfQNjqC7YeLsduEcpeHD9Zl8eZ3B4hzhOE14HR7+e3FfYgMs7E1u5jPth0hKToCrzHYRPjFOT2YOqwT4XYbBWVOisrdRITZWLHzKLuPlvDmdwcoc3oIswlur3U8O69PCkUVLtYdKAAgNakN+3PLCLcLIoLb4+WuC/vwi3PSiI4MwxjDvzdl8/7aLFKTohnUOY4L+7Unvs3xtaiWSkTWGGNG1TlPE0Hz+P2Hm3lr9UHSH7yoRhV82ZbD3PrGGpbceTaDOscHMUJVn8OFFTy2ZCv/3pSNTcDr9y8xJq0tcVHhrN2fT5fEKGad34tJ/dtzML+MWxakExMZRmZ+OcdKKnl06iBuHJsavII0UWmlm+jI5mkU2J9byl8+28UH67OofUiJsNtwerzV78NswrSzunLv5H44PV5mvPY9W7OLAKsGckG/dpRUugmzCdmFFWzNLqJTvIMO8Q42ZRXWqLGIwMUD2jO0awJF5W7O6Z3MhswCXly+h1hHGD8d151wu/D5tqOM65nEdaO7Em6zMfujLXy4/hBR4XYuGdSBrPxyvs/Io2O8g/wyq8mtU7yDX03siYhw+ZCOJLSJaJbvKlg0EQSQMYZv9+Yx6801jO+VzJzrR9SYv2pPLtNf+ZZ/3jKG8T2TgxRl6Np5pJgvdxzl3D4p9PNd1gvWflu29YjV1LP5MAA3n52GMYYxPdrSIzmGogoXQ7qcuMmitNLNkaIKeqTEnHDZ1q6owsWx4kryy5zkl7o4kFdGVkE5Q7smEGYTHOE2+neMo2N8VPU6Hq/haHEFdhGSYyKx2X7ogzDGsHzHUd781jrrH9gpjh4p1r65oF87uidFExF2/FVSTrcXm1DvFVTGGNbsz+fdtVks2XiIOEc4t5yTVp3I1+zP5+6313PIN3JwVLidF28cwfl92zXn13VaaSIIoEWrD3LvuxuJCrfzj1+MZmRq2xrzN2cVcvkL/2XuTSO5eGCHIEXZ8rg8XjZlFVJY5sJmE7okRtE9yWpay8wvIyYyjKSYSI4WVZBb6qR/x7jjtvHKV3t5Yuk2wDorvfeSvgD87cs9JLQJZ09OKckxEVw6qCO3ntuDrm3bnL4CqjOGx2uwyfGDQ5Y7PRwpqqDM6eF/39nAnpwS/nbDSM7v1w5jTIvrNG8oEQSss1hEXgMuB44aYwbVMV+A54DLgDJgpjFmbaDiCQSv1/DyV3sY1DmORb8cR5uI47/OWIc1rVgvIW20HYeL+e2/1rM5q6jG9KhwO+F2qb4c99w+KWzKLCC/zDo7vHhAe4or3HywPosDeWUUV7i5bHAH7r6oD3/8zw4e/7eVFEamJiLA7CtSuWlcd+y2lvUPrZpXffs/KsJe3a83/+bRXP/Kt/xs3mriHGFUuL30TIlhUv92nNsnhcFd4pvlqkBjDHmlTpJO80CVgbxqaB7wV2BBPfMvBXr7fsYAf/P9bjFW7j7GnpxSnr12aJ1JACDW11+gN5WdmNPtZfZHW3jr+wMktongqauH0KtdDF6vYd+xUrZlF1PucjO0SwKHCit4deVeOiVEcePYVN76/iBfbLdGe+3XIZafDO9Mu9hIbjuvJ2F2G6/8dCSL0g+SkVvGbyf10ZuuVJMkx0Ty0Z1n8/rXGRzIKyM6ws7GzEJeWL6b57/YTUSYjYSocMqcHgZ0iqNDnIOENuFWshjQnk4JUXVut9LtYXNWEekZeazOyGfN/jzyy1yMTE2kb4dY0pKiGdApjvwyJ5MHdiA8QH+3AW0aEpHuwJJ6agQvA18aY97yvd8BTDTGZDe0zVGxsSZ95MgARNt0v+l5GcsT0vh27UtEGk+dyzjFRp8xv+WeAyu549C39W7LAKF6XuoUG1/Hp/Jah1GsTOjOzdnp3Jm1ikR3w092K7JH4PC6iTBevMA+RyJxnkqSXWUh+12q0ysvLIr02M6sju1CfpgDh9fN1uh25IW1IT/cQWFYFDbj5bLcHZxXuI+c8GhcYic/LIqv41PJcCTislk1ibTyPEYVZ9GlspDFyf0pCHOQG/7DlYZX5Wzm6T0fn/S4QLJixelvGmqEzoD/RdeZvmnHJQIRuRW4FWBI5Jkztv+a2E6MLs6sNwkARBgvkV4XxWH1x70xuj039b+G17a/x8iSQ4EI9YzkEhvzOozg5Y6jORYRTYTXzR/2fsL0oxsbtX6c54fB/GxAz4r8AEWqVN3ausu5OH83F+fvrnP+Xkcib7cbwj/bDWFJcv/q6RFeN+MLD3Bh/h6GlWQzsiSLFFdZ9fy7slYBsMeRyAFHAutiOvF8l/H0K8vhluzm7yNtETeUGWPmAnPB6izmyy+DGxBwrKSS/Y9/xvU/GQfn3d3gsrGPf0rxuKspmPw7thwqYnzPpOqOpgqXh3vnfE3h4WI+u+NhRl7S73SEH1QFZU62HCrisSVb2X64mHP7pDBzfCrjeybjCJ8a7PCUajY9gPuBe72GvTklpMRGEhVhxxhqjEBQn56+n4nGkLRqPz8efjFEneS9DQ10bgczEWQBXf3ed/FNaxGqblQZkXr83aG1OcLt/PO7A7y7JpNKt5enrhrClSM687v3N/OvNQfxGoiPCme177b+lqzM6eaFL3bjdHt58Ef9a1xZYYzhhS928+fPdmIMdIhz8PJNI5msV1OpVs5uE3q3P/k7mUWsO84DJZiJYDFwh4gsxOokLjxR/8CZZO2BfMJswuBG3CTm9t0Ac+mgDuzPK+OJpdt4Z20m3+/L4/ox3ZjUvz3f7svltf/uo8LladSZQnN5e/UBXv5qL+f2TuHeS/rW2+ndkMJyF++vzeSywR258dXv2HnEejJbWnI0Trd1B2peqZO1+/PZe6yUK4Z24keDO3B275TqITiUUsETyMtH3wImAskikgn8HggHMMa8BCzFunR0N9bloz8LVCzNzRjDqj25DOwU16iD9pwbRlBa6ebcPinsPlrClL/+l6z8ch778SBu8t3AYjC8vGIv6w8WMLZH0knF5fWaGjfjnMjBvDIe+WgriW0imPdNBknREdx5Ye8mf+4ji7fw3rosnvzPdlwew6szRvGXz3bx4AebAYgMs9E2OoJe7WKYdX4vrhrRucVdg61UaxawRGCMmX6C+Qa4PVCfH0iL0g+y/mABjzRyuN+Rfs1HvdrFsO7hSUTYbTUOhiNT2yIC3+7NPalEUO70MPGZ5cwcn8avJvZkf24p/96Uzc/PTjvu+uZKt4c5X+xmUXqmVZ7bxvH7D7cwd+VeRqQmMqBjHInR1u30uSWVxDrC67x7E+Dr3cd4b10Wkwa0Z83+fO66sDcX9m9PalI0//h2vw7upVQLoPXyJqpweXhsyTbG9UiqPptvqrpuPImPCmdolwSWbz/K9aO7sflQIRf0a9/obf5nSzZHiip57vOdxEeF88f/bKew3EW3tm24fEgnXB4vy7cfZcXOHNbsz2f74WIu6NeOX5yTRueEKH4zqQ8/emElN/z9O4Z2TWDhLWN58uNtzF+1n8gwG0O7JDCxXwo3jU2tvjfiaHEFd7+9nrTkaF6YPpzIsB+SW692MTouvlIthCaCJsrILaWk0s31Y7o1qRmmMSYNaM/Tn+zgV2+uZc3+/CaNWPqv9Ezax0WSX+bigfc30b9jHOF24d8bs0mOieR3729iT04psY4wOidE8dKNI7lk0A+dtAM6xfHh7RP4YvtR/vLZLi5/YSV7ckq5fkw3osLtpGfk8dR/dvD3lfv441VDGN4tgVvmp1NU4WL+zaNPa7+GUqp5aSJoon05pYDVEdrcLvYlgjX7revhF64+yC3npJHYJgKbTSh3enjxy91syCzkprGpTBpg1Rgy88v4Zk8uv5nUhz7tYyksd3LViC48tmQrb60+yOfbjtIh3sHfbhjBpAHt672rdkiXBAZ2imfppmx2Hinh8R/XHE1zY2YB9727iVsWpCNijd/zwvThdY7zo5RqOTQRNNG+XCsRBOLZAr3axdA9qQ25JU4GdY5nwaoMXlm5lwv6teOF6cN5ZeVeXvhiN+1iI7llQTptoyMY0iWe4V2tPogrR3SmS+IPA6f9aEgn5q/aT4+UaN771fhGDaNrtwl/u3Eku46U1KgxgJUo3ps1ng/XZ7H3WCk/GtyxUaNzKqXObJoImmhfTikpsZEBuexRRHjyqiE43V5E4Kevfc+Ensl8uvUIv/3XBrYdKmJsj7bM+9lonv98F5sPFfHljhy+35fH+J5JNZIAwKjURJ66aghn905u0ljqPVNi6FnPkMqOcDvTzupW5zylVMukiaCJMnJLSUsK3JPG/K8Y2vj7i4l1hPP857t49tOdAPzyvB44wu3ce0k/nG4vE59ezqHCCq4Z1eW4bdlswrVndT1uulJK+dMhGJto37GygPQP1KXq6pzbzutJr3YxRITZuGRQx+r5EWE2/u/SfvTrEKt35yqlTprWCJqguMLFsZLK0/7s4YgwG3//6Sgy88uJrzXOyNRhnZk6rPNpjUcp1bpoImiCjGPW6IBpyaf/SVbdk6NPewJSSoUGbRpqgm2HrSdm9Wqnd8oqpVoPTQRNsDmrkOgIOz30zFwp1YpoImiCjZmFDOwc3+x3FCulVDBpImgkl8fL1uwihjRi2GmllGpJNBE00s4jxTjdXgZ30USglGpdNBE00uasQgAdUkEp1epoImik9Ix84qPCSW17+i8dVUqpQNJE0AjGGL7Zk8u4HknaUayUanU0ETTCgbwysgrKmdDr5B4hqZRSZzJNBI3wzZ5cAMb1TA5yJEop1fw0ETTCN3tyaRcbSc8UvZFMKdX6aCJohHUH8jmre9saD5tXSqnWQhPBCeSVOsnML2eI3j+glGqlNBGcwCbf/QN6I5lSqrXSRHACmzILABikQ0sopVopTQQnsDGzkB7J0cQ5wk+8sFJKtUCaCE5gU1ahNgsppVo1TQQNKCxzkV1YQf+OccEORSmlAiagiUBELhGRHSKyW0Tuq2N+NxFZLiLrRGSjiFwWyHiaatfRYgD6tI8JciRKKRU4AUsEImIH5gCXAgOA6SIyoNZiDwKLjDHDgeuAFwMVz8nYdbQEgN76aEqlVCsWyBrBaGC3MWavMcYJLASm1lrGAFXtLvHAoQDG02Q7jxQTFW6nc0JUsENRSqmACWQi6Awc9Huf6ZvmbzZwo4hkAkuBO+vakIjcKiLpIpKek5MTiFjrtOtICb3bx+iIo0qpVi3YncXTgXnGmC7AZcAbInJcTMaYucaYUcaYUSkpKactuF1Hi+nVTvsHlFKtWyATQRbQ1e99F980fz8HFgEYY1YBDuCMGOKzsNzFkaJK+rTX/gGlVOsWyESwGugtImkiEoHVGby41jIHgAsBRKQ/ViI4fW0/DdibY3UU90zRGoFSqnULWCIwxriBO4BPgG1YVwdtEZFHRWSKb7HfAreIyAbgLWCmMcYEKqamOJBXBkBqkj6aUinVuoUFcuPGmKVYncD+0x72e70VmBDIGE7WQV8i6JqoiUAp1boFu7P4jHUwr5yU2EiiIuzBDkUppQJKE0E9DuSV0TVR7x9QSrV+mgjqcSCvjG5ttVlIKdX6aSLws2BVBit25uDyeMkuLNdEoJQKCQHtLG5p/rRsJ6NSE+me1AavgS6aCJRSIUBrBD4llW4Ky10cyCvjYF45gNYIlFIhQROBT1a+dfA/kFdGRm4pAF01ESilQoAmAp+sAuu+gUq3l6925hAbGUaneEeQo1JKqcDTROBTVSMAWLEzh/4d4xDRUUeVUq2fJgKfzIIfEkGl28uATvp4SqVUaNBE4JOVX06neAdVlYAB+pxipVSI0ETgk1VQTvfkaDrFW3cTa41AKRUqNBH4ZOWX0zkhiq5to7DbRB9Io5QKGXpDGVDp9nC0uJLOiVF0T44mPiocR7gONqeUCg2aCICc4koAOsY7mHZWtyBHo5RSp5c2DQH5pS4A2kZHBjkSpZQ6/TQRAHllTgDaRocHORKllDr9NBEA+aVWIkhsExHkSJRS6vTTRADklVbVCDQRKKVCjyYCIL/MiU0gzqFNQ0qp0KOJAKtGkNgmAptNxxZSSoUeTQRYNYJEbRZSSoUoTQRYNYK22lGslApRmgiw7iNI1EtHlVIhShMB1n0EesWQUipUhXwiMMaQ7+ssVkqpUBTyiaC40o3ba7RGoJQKWQFNBCJyiYjsEJHdInJfPctcKyJbRWSLiPwzkPHURe8qVkqFukaPPioiacaYfSea5jfPDswBJgGZwGoRWWyM2eq3TG/gfmCCMSZfRNqdTCFOhd5VrJQKdU2pEbxbx7R3Glh+NLDbGLPXGOMEFgJTay1zCzDHGJMPYIw52oR4mkVVItD7CJRSoeqENQIR6QcMBOJF5Eq/WXGAo4FVOwMH/d5nAmNqLdPH9xlfA3ZgtjHmP3XEcCtwK0C3bs37vIAjRdazCNrH6RDUSqnQ1Jimob7A5UACcMgLglMAABpGSURBVIXf9GKsM/pT/fzewESgC/CViAw2xhT4L2SMmQvMBRg1apQ5xc+s4VBBOXab0C62oZymlFKtV2MSQRvgHmCuMWZVE7adBXT1e9/FN81fJvCdMcYF7BORnViJYXUTPueUHCoop0OcA7uOM6SUClGN6SPoBvwLeEpEZovIGBFpzFFzNdBbRNJEJAK4Dlhca5kPsGoDiEgyVlPR3sYG3xyyCsrplKC1AaVU6DphIjDG/NEYcwFwGbABuBlYKyL/FJGfikj7etZzA3cAnwDbgEXGmC0i8qiITPEt9gmQKyJbgeXA/xpjck+9WI2XXVhBp4So0/mRSil1Rmn05aPGmGLgfd8PIjIAuBRYAEyuZ52lwNJa0x72e22A3/h+Tjuv15BdWM6PEjoG4+OVUuqM0JT7CGzAUKATUA5sNsb8CfhTgGILuGMllbg8RmsESqmQ1pjLR3sC/wdcBOwCcrAuG+0jImXAy8B8Y4w3kIEGQlZBOQCd4rWPQCkVuhpTI3gc+BvwS19TTjVf/8B04CZgfvOHF1iHCioAtEaglAppJ0wExpjpACISCVTWml1gjPlLIAI7HbILfTUCTQRKqRDWlCEm6rqHoCn3FZxRvtl9jJdW7KF9XCRxjkZ3lSilVKvTmD6CDljDRUSJyHCg6h6COKybzVqkpz7ZQWSYnddmnkXjbotQSqnWqTGnwpOBmVh3Bv+JHxJBEfBAYMIKvIIyJyNSE+nbITbYoSilVFA1po9gPjBfRK4yxtQ1AmmLVFzh1iYhpZSiEX0EInKjiEh9SUBEeorI2c0fWuAYYyiqcBEXpQ+sV0qpxpwSJwHrRWQNsIYf7iPoBZwHHAPqfPrYmarC5cXlMcRqjUAppRrVNPSciPwVuACYAAzBurN4G3CTMeZAYENsfsUVLgDiHFojUEqpRp0SG2M8wKe+nxavqCoRaNOQUko16vLRhxuYbYwxjzVjPKdFUYUbQJuGlFKKxtUISuuYFg38HKv/oOUlgnJtGlJKqSqN6SOoHl1URGKBu4CfYT2MvkWOPFpVI4iP0hqBUko16kgoIm2xnhlwA9bgciOMMfmBDCyQqjqLY7VGoJRSjeojeBq4Euvh8YONMSUBjyrAisqtGoE2DSmlVOMGnfst1sNoHgQOiUiR76dYRIoCG15gFFW4CLcLjvCmjLmnlFKtU2P6CFrd0bK4wkWcI1wHm1NKKZo2DHWrUVTu1ktHlVLKJzQTgY4zpJRS1UIyEVgjj2oiUEopCNFEUFTu0qYhpZTyCc1E4OssVkopFaKJoLjCTZzeVayUUkAIJgK3x0uZ00NMpNYIlFIKQjARlLk8AERH2oMciVJKnRkCmghE5BIR2SEiu0Wk3qeYichVImJEZFQg4wEoq6xKBNo0pJRSEMBEICJ2YA5wKTAAmC4iA+pYrmpE0+8CFYu/Uqc1zlCbCK0RKKUUBLZGMBrYbYzZa4xxYg1bPbWO5R4D/ghUBDCWatU1ggitESilFAQ2EXQGDvq9z/RNqyYiI4Cuxph/N7QhEblVRNJFJD0nJ+eUgqquEWgfgVJKAUHsLBYRG/As1uimDTLGzDXGjDLGjEpJSTmlzy3zJQKtESillCWQiSAL6Or3votvWpVYYBDwpYhkAGOBxYHuMC6t1KuGlFLKXyATwWqgt4ikiUgEcB2wuGqmMabQGJNsjOlujOkOfAtMMcakBzCm6hpBG60RKKUUEMBEYIxxA3cAnwDbgEXGmC0i8qiITAnU555IqXYWK6VUDQE9GhpjlgJLa017uJ5lJwYylipVNYIovXxUKaWAELyzuNTpIcJuIyIs5IqulFJ1CrmjYVmlWy8dVUopPyGXCEqdHu0fUEopPyGXCMqcbh1eQiml/IRcIiit9NBGB5xTSqlqIZcIypxuorVGoJRS1UIuEZRWevRmMqWU8hNyiaDM6dbhJZRSyk/IJYJSp9YIlFLKX8glgrJK7SNQSil/IZUIvF5DmUuvGlJKKX8hlQgq3B6MQWsESinlJ6QSQdXIo1ojUEqpH4RYIqh6OpnWCJRSqkpoJQJ9KI1SSh0npBJBudPXNKQ1AqWUqhZSicDp9gLoswiUUspPSB0RnR5NBEopVVtIHRFdHgNAhD2kiq2UUg0KqSNiVdNQuCYCpZSqFlJHRJenKhFIkCNRSqkzR0glAu0jUEqp44XUEbGqRqB9BEop9YOQOiJqH4FSSh0vpI6I1X0E2jSklFLVQuqIWHX5qHYWK6XUD0IqEVTfWaxNQ0opVS2gR0QRuUREdojIbhG5r475vxGRrSKyUUQ+F5HUQMbj8ngJtwsiWiNQSqkqAUsEImIH5gCXAgOA6SIyoNZi64BRxpghwDvAU4GKB6wagXYUK6VUTYE8Ko4Gdhtj9hpjnMBCYKr/AsaY5caYMt/bb4EuAYzHVyPQRKCUUv4CeVTsDBz0e5/pm1afnwMf1zVDRG4VkXQRSc/JyTnpgJweozeTKaVULWfEUVFEbgRGAU/XNd8YM9cYM8oYMyolJeWkP8fl8WpHsVJK1RLIR3VlAV393nfxTatBRC4CfgecZ4ypDGA8vj4C7ShWSil/gTw9Xg30FpE0EYkArgMW+y8gIsOBl4EpxpijAYwF0D4CpZSqS8COisYYN3AH8AmwDVhkjNkiIo+KyBTfYk8DMcC/RGS9iCyuZ3PNwuXxah+BUkrVEtCnuBtjlgJLa0172O/1RYH8/NqcHqM1AqWUqiWkjopOt0c7i5VSqpaQOiq6PIbwMO0sVkopfwFtGjrTuDxe4hwhVWTVSrlcLjIzM6moqAh2KOoM43A46NKlC+Hh4Y1eJ6SOijrEhGotMjMziY2NpXv37jp2lqpmjCE3N5fMzEzS0tIavV5IHRWdHq8+i0C1ChUVFSQlJWkSUDWICElJSU2uKYbUUVHvLFatiSYBVZeT+bsIqaOiy200ESilVC0hdVR0ebx61ZBSzSA3N5dhw4YxbNgwOnToQOfOnRk2bBgJCQkMGFB7tPnA+uCDD9i6dWv1+4cffpjPPvusydvJyMhg0KBBzRlaixFSicCpQ0wo1SySkpJYv34969ev57bbbuPuu++ufm+zNf//mNvtrnde7UTw6KOPctFFp/Ve1SYzxuD1eoMdRrWQu2pIm4ZUq/PrX8P69c27zWHD4C9/OalVPR4Pt9xyC9988w2dO3fmww8/JCoqij179nD77beTk5NDmzZteOWVV+jXrx8ZGRncfPPNHDt2jJSUFF5//XW6devGzJkzcTgcrFu3jgkTJnD77bcft35eXh6LFy9mxYoVPP7447z77rs89thjXH755Vx99dWsXr2au+66i9LSUiIjI/n888/Jzc3lpptuorS0FIC//vWvjB8/vt7ylJSUMHXqVPLz83G5XDz++ONMnTqV++67j65du3L77bcDMHv2bGJiYrjnnnt4+umnWbRoEZWVlfzkJz/hkUceISMjg8mTJzNmzBjWrFnD0qVLefLJJ1m9ejXl5eVcffXVPPLIIwAsXbqU3/zmN0RHRzNhwgT27t3LkiVLKC0t5c4772Tz5s24XC5mz57N1KlT6429sUIqEeigc0oF3q5du3jrrbd45ZVXuPbaa3n33Xe58cYbufXWW3nppZfo3bs33333HbNmzeKLL77gzjvvZMaMGcyYMYPXXnuN//mf/+GDDz4ArMtkv/nmG+x2OxdeeGGd60+ZMqX6wO/P6XQybdo03n77bc466yyKioqIioqiXbt2fPrppzgcDnbt2sX06dNJT0+vtzwOh4P333+fuLg4jh07xtixY5kyZQrTpk3j17/+dXUiWLRoEZ988gnLli1j165dfP/99xhjmDJlCl999RXdunVj165dzJ8/n7FjxwLwxBNP0LZtWzweDxdeeCEbN26kT58+/PKXv+Srr74iLS2N6dOnV8fyxBNPcMEFF/Daa69RUFDA6NGjueiii4iOjj6lfRYyicDjNXgNOuican1O8sw9UNLS0hg2bBgAI0eOJCMjg5KSEr755huuueaa6uUqK61R51etWsV7770HwE033cS9995bvcw111yD3W5vcP367Nixg44dO3LWWWcBEBcXB0BpaSl33HEH69evx263s3Pnzga3Y4zhgQce4KuvvsJms5GVlcWRI0cYPnw4R48e5dChQ+Tk5JCYmEjXrl157rnnWLZsGcOHDwesGsWuXbvo1q0bqamp1UkArOQxd+5c3G432dnZbN26Fa/XS48eParvA5g+fTpz584FYNmyZSxevJhnnnkGsC4jPnDgAP3792+wDCcSMonA5bHa47RGoFRgRUZGVr+22+2Ul5fj9XpJSEhgfRObsKrOdE92/br8+c9/pn379mzYsAGv14vD4Whw+TfffJOcnBzWrFlDeHg43bt3r75O/5prruGdd97h8OHDTJs2DbASx/33388vf/nLGtvJyMiocea+b98+nnnmGVavXk1iYiIzZ8484fX/xhjeffdd+vbtezJFr1fIHBUr3VWJQK8aUup0i4uLIy0tjX/961+AdUDbsGEDAOPHj2fhwoWAddA955xzmrR+bGwsxcXFx63Tt29fsrOzWb16NQDFxcW43W4KCwvp2LEjNpuNN954A4/H02DshYWFtGvXjvDwcJYvX87+/fur502bNo2FCxfyzjvvVNdWJk+ezGuvvUZJSQkAWVlZHD16/ONWioqKiI6OJj4+niNHjvDxxx9Xx713714yMjIAePvtt6vXmTx5Mi+88ALGGADWrVvXYOyNFTKJoKpGoE1DSgXHm2++yauvvsrQoUMZOHAgH374IQAvvPACr7/+OkOGDOGNN97gueeea9L61113HU8//TTDhw9nz5491ctHRETw9ttvc+eddzJ06FAmTZpERUUFs2bNYv78+QwdOpTt27efsH39hhtuID09ncGDB7NgwQL69etXPW/gwIEUFxfTuXNnOnbsCMDFF1/M9ddfz7hx4xg8eDBXX311nYlq6NChDB8+nH79+nH99dczYcIEAKKionjxxRe55JJLGDlyJLGxscTHxwPw0EMP4XK5GDJkCAMHDuShhx5q7NffIKnKLC3FqFGjTEMdO/XJLixn3B++4MkrB3Pd6G4BiEyp02fbtm2n3C6szlwlJSXExMRgjOH222+nd+/e3H333Y1ev66/DxFZY4wZVdfyIXN67HJbCU/7CJRSZ7pXXnmFYcOGMXDgQAoLC4/rb2huIdNZ7PS1A+qgc0qpM93dd9/dpBrAqQqZo6LTVyOI0M5ipZSqIWQSgXYWK6VU3ULmqKj3ESilVN1C5qjodGsiUEqpuoTMUdGpNQKlmtUTTzzBwIEDGTJkCMOGDeO7774L6OdNnDixwTGB1MkLmauGXB6rszhS+wiUOmWrVq1iyZIlrF27lsjISI4dO4bT6Qx2WAHndrsJC2t9h83WV6J6aB+Baq0e+WgLWw8VNes2B3SK4/dXDKx3fnZ2NsnJydXjCiUnJ1fPe/TRR/noo48oLy9n/PjxvPzyy4gIEydOZPjw4axcuZLS0lIWLFjAH/7wBzZt2sS0adN4/PHHycjIqL6jdu3atQwcOJAFCxbQpk2bGp+/bNkyfv/731NZWUnPnj15/fXXiYmJqbHMK6+8wty5c3E6nfTq1Ys33nij+q7cffv2YbPZKC0tpV+/fuzdu5cDBw7UOUx27eGwr7vuOu666y4qKiqIiori9ddfp2/fvpSVlTFz5kw2b95M3759OXToEHPmzGHUqFGNijeYQuao6NSxhpRqNhdffDEHDx6kT58+zJo1ixUrVlTPu+OOO1i9ejWbN2+mvLycJUuWVM+LiIggPT2d2267jalTpzJnzhw2b97MvHnzyM3NBaxRQ2fNmsW2bduIi4vjxRdfrPHZx44d4/HHH+ezzz5j7dq1jBo1imefffa4GK+88kpWr17Nhg0b6N+/P6+++irx8fEMGzasOt4lS5YwefJkwsPDufXWW3nhhRdYs2YNzzzzDLNmzareVtVw2M8++yz9+vVj5cqVrFu3jkcffZQHHngAgBdffJHExES2bt3KY489xpo1a5oUbzCFTI1A+whUa9XQmXugxMTEsGbNGlauXMny5cuZNm0aTz75JDNnzmT58uU89dRTlJWVkZeXx8CBA7niiisAmDJlCgCDBw9m4MCB1ePz9OjRg4MHD5KQkEDXrl2rx9258cYbef7557nnnnuqP/vbb79l69at1cs4nU7GjRt3XIybN2/mwQcfpKCggJKSEiZPngxQ/YyC888/n4ULFzJr1qwTDnNdNRw2WIPQzZgxg127diEiuFwuAP773/9y1113ATBo0CCGDBnSpHiDKaCJQEQuAZ4D7MDfjTFP1pofCSwARgK5wDRjTEYgYqlqGtI+AqWah91uZ+LEiUycOJHBgwczf/58rrvuOmbNmkV6ejpdu3Zl9uzZNYZWrmpKstlsNYarttls1Y+jFKlZa6/93hjDpEmTeOuttxqMb+bMmXzwwQcMHTqUefPm8eWXXwJWMnrggQfIy8tjzZo1XHDBBZSWljY4zLX/wHQPPfQQ559/Pu+//z4ZGRlMnDixwTgaG28wBeyoKCJ2YA5wKTAAmC4itZ9q/XMg3xjTC/gz8MdAxePSy0eVajY7duxg165d1e/Xr19Pampq9UE/OTmZkpIS3nnnnSZv+8CBA6xatQqAf/7zn5x99tk15o8dO5avv/6a3bt3A9aDZup6uExxcTEdO3bE5XLx5ptvVk+PiYnhrLPO4q677uLyyy/Hbrc3OMx1bYWFhXTu3BmAefPmVU+fMGECixYtAmDr1q1s2rSpSfEGUyCPiqOB3caYvcYYJ7AQqP1wzanAfN/rd4ALpXb6byZVVw3pWENKnbqSkhJmzJjBgAEDGDJkCFu3bmX27NkkJCRwyy23MGjQICZPnlz9dLCm6Nu3L3PmzKF///7k5+fzq1/9qsb8lJQU5s2bx/Tp0xkyZAjjxo1j+/btx23nscceY8yYMUyYMKHG0NFgNQ/94x//qH6YDNQ/zHVt9957L/fffz/Dhw+vrsUAzJo1i5ycHAYMGMCDDz7IwIEDiY+Pb3S8wRSwYahF5GrgEmPML3zvbwLGGGPu8Ftms2+ZTN/7Pb5ljtXa1q3ArQDdunUb6f9giMZatuUwH6zP4i/ThuswE6rFa63DUGdkZHD55ZezefPmYIfSZB6PB5fLhcPhYM+ePVx00UXs2LGDiIiI0x5LU4ehbhGdxcaYucBcsJ5HcDLbuHhgBy4e2KFZ41JKqSplZWWcf/75uFwujDG8+OKLQUkCJyOQiSAL6Or3votvWl3LZIpIGBCP1WmslApB3bt3b5G1AbAemdlS73wOZBvJaqC3iKSJSARwHbC41jKLgRm+11cDX5iW9sg0pYJE/1VUXU7m7yJgicAY4wbuAD4BtgGLjDFbRORREZniW+xVIElEdgO/Ae4LVDxKtSYOh4Pc3FxNBqoGYwy5ubk4HI4mrRcyzyxWqjVxuVxkZmbWuEZfKbBOErp06UJ4eHiN6S2+s1gpVVN4eDhpaWnBDkO1EnodpVJKhThNBEopFeI0ESilVIhrcZ3FIpIDNP3WYksycOyES7U+oVhuLXPoCMVyn0yZU40xKXXNaHGJ4FSISHp9veatWSiWW8scOkKx3M1dZm0aUkqpEKeJQCmlQlyoJYK5wQ4gSEKx3Frm0BGK5W7WModUH4FSSqnjhVqNQCmlVC2aCJRSKsSFTCIQkUtEZIeI7BaRVjvKqYhkiMgmEVkvIum+aW1F5FMR2eX7nRjsOE+ViLwmIkd9T7mrmlZnOcXyvG/fbxSREcGL/OTVU+bZIpLl29/rReQyv3n3+8q8Q0QmByfqUyMiXUVkuYhsFZEtInKXb3qr3dcNlDlw+9oY0+p/ADuwB+gBRAAbgAHBjitAZc0AkmtNewq4z/f6PuCPwY6zGcp5LjAC2HyicgKXAR8DAowFvgt2/M1Y5tnAPXUsO8D3dx4JpPn+/u3BLsNJlLkjMML3OhbY6Stbq93XDZQ5YPs6VGoEo4Hdxpi9xhgnsBCYGuSYTqepwHzf6/nAj4MYS7MwxnwF5NWaXF85pwILjOVbIEFEOp6eSJtPPWWuz1RgoTGm0hizD9iN9X/Qohhjso0xa32vi7GebdKZVryvGyhzfU55X4dKIugMHPR7n0nDX2xLZoBlIrJGRG71TWtvjMn2vT4MtA9OaAFXXzlb+/6/w9cM8ppfs1+rK7OIdAeGA98RIvu6VpkhQPs6VBJBKDnbGDMCuBS4XUTO9Z9prLpkq79mOFTKCfwN6AkMA7KBPwU3nMAQkRjgXeDXxpgi/3mtdV/XUeaA7etQSQRZQFe/911801odY0yW7/dR4H2sKuKRquqx7/fR4EUYUPWVs9Xuf2PMEWOMxxjjBV7hhyaBVlNmEQnHOiC+aYx5zze5Ve/rusocyH0dKolgNdBbRNJEJAK4Dlgc5JianYhEi0hs1WvgYmAzVlln+BabAXwYnAgDrr5yLgZ+6ruiZCxQ6Nes0KLVav/+Cdb+BqvM14lIpIikAb2B7093fKdKRATr2ebbjDHP+s1qtfu6vjIHdF8Hu4f8NPbEX4bV+74H+F2w4wlQGXtgXT2wAdhSVU4gCfgc2AV8BrQNdqzNUNa3sKrHLqw20Z/XV06sK0jm+Pb9JmBUsONvxjK/4SvTRt8BoaPf8r/zlXkHcGmw4z/JMp+N1eyzEVjv+7msNe/rBsocsH2tQ0wopVSIC5WmIaWUUvXQRKCUUiFOE4FSSoU4TQRKKRXiNBEopVSI00SgVDMQkQQRmRXsOJQ6GZoIlGoeCYAmAtUiaSJQqnk8CfT0jRP/dLCDUaop9IYypZqBb5TIJcaYQUEORakm0xqBUkqFOE0ESikV4jQRKNU8irEeK6hUi6OJQKlmYIzJBb4Wkc3aWaxaGu0sVkqpEKc1AqWUCnGaCJRSKsRpIlBKqRCniUAppUKcJgKllApxmgiUUirEaSJQSqkQ9/8BTiQYM5L6jJsAAAAASUVORK5CYII=\n"
          },
          "metadata": {
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Checking the central limit theorem"
      ],
      "metadata": {
        "id": "fWSI0470nso9"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "def get_num_norm(K, T):\n",
        "  t = T\n",
        "  N_t = [0]*K\n",
        "  max_N = [0]*K\n",
        "  for i in range(K):\n",
        "    for j, elem in enumerate(sum_array_2[i]):\n",
        "      if elem <= t:\n",
        "        # N_t[i] = j\n",
        "        max_N[i] += 1 \n",
        "\n",
        "  num_norm = []\n",
        "  for i in range(K):\n",
        "    val = (max_N[i] - statistics.mean(max_N))/(statistics.stdev(max_N))\n",
        "    num_norm.append(val)\n",
        "  return num_norm"
      ],
      "metadata": {
        "id": "POLB9FsPbGTf"
      },
      "execution_count": 59,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "_, sum_array_2, nums_of_intervals_2 = get_recov_processes(K = 100,T = 100)"
      ],
      "metadata": {
        "id": "4A01LIBCF9rN"
      },
      "execution_count": 55,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "num_norm = get_num_norm(K = 100,T = 100)"
      ],
      "metadata": {
        "id": "Zi_Bx3Znbdzj"
      },
      "execution_count": 60,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "plt.hist(num_norm, density=True);"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 265
        },
        "id": "HrBHEEV7au1G",
        "outputId": "d24142d5-a101-4b4e-d1fc-1705cc1aaeee"
      },
      "execution_count": 61,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD4CAYAAADiry33AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAQmElEQVR4nO3df6xfdX3H8edrZWCi06HcZFlbaNW6WKeB5Fq2GDFRfpSxtGzBWBaXmpE0LDRzYcuswUBSY4KYsH9WA01sZoysom7ZTajpmKDGOLQXYWDLOi6VwW1MrJToFglQeO+PezRfbm65p9zvvd/2c5+P5Kbn8+v7fZ/S+7qHc77n3FQVkqR2/caoC5AkLS6DXpIaZ9BLUuMMeklqnEEvSY07a9QFzHbeeefVmjVrRl2GJJ1RHnzwwZ9V1dhcY6dd0K9Zs4bJyclRlyFJZ5Qk/3OyMU/dSFLjDHpJapxBL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4067O2Ol+azZcc/I3vvJW68a2XtLr5VH9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BLUuN6BX2SjUkOJ5lKsmOO8euTPJrk4STfTbK+61+T5Lmu/+Ekdwx7ByRJr27eh5olWQHsAi4DpoEDSSaq6tDAtLuq6o5u/ibgdmBjN/ZEVV043LIlSX31OaLfAExV1ZGqegHYC2wenFBVvxhovh6o4ZUoSVqIPkG/Enh6oD3d9b1CkhuSPAHcBvzVwNDaJA8l+XaS9y+oWknSKRvaxdiq2lVVbwM+AXyq6/4JcH5VXQTcCNyV5I2z1ybZlmQyyeSxY8eGVZIkiX5BfxRYPdBe1fWdzF7gaoCqer6qnum2HwSeAN4xe0FV7a6q8aoaHxsb61u7JKmHPkF/AFiXZG2Ss4EtwMTghCTrBppXAY93/WPdxVySvBVYBxwZRuGSpH7m/dRNVZ1Ish3YD6wA9lTVwSQ7gcmqmgC2J7kUeBF4FtjaLb8E2JnkReBl4PqqOr4YOyJJmluv3xlbVfuAfbP6bh7Y/vhJ1n0d+PpCCpQkLYx3xkpS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BLUuMMeklqnEEvSY0z6CWpcQa9JDWu19MrNb81O+4Z2Xs/eetVI3tvSac/j+glqXEGvSQ1zqCXpMYZ9JLUOINekhrXK+iTbExyOMlUkh1zjF+f5NEkDyf5bpL1A2Of7NYdTnLFMIuXJM1v3qBPsgLYBVwJrAeuHQzyzl1V9e6quhC4Dbi9W7se2AK8C9gIfL57PUnSEulzRL8BmKqqI1X1ArAX2Dw4oap+MdB8PVDd9mZgb1U9X1U/Bqa615MkLZE+N0ytBJ4eaE8DF8+elOQG4EbgbOCDA2sfmLV25RxrtwHbAM4///w+dUuSehraxdiq2lVVbwM+AXzqFNfurqrxqhofGxsbVkmSJPoF/VFg9UB7Vdd3MnuBq1/jWknSkPU5dXMAWJdkLTMhvQX4s8EJSdZV1eNd8yrgV9sTwF1Jbgd+F1gH/GAYhUujMMpnGo2Kz1I6880b9FV1Isl2YD+wAthTVQeT7AQmq2oC2J7kUuBF4Flga7f2YJK7gUPACeCGqnppkfZFkjSHXk+vrKp9wL5ZfTcPbH/8VdZ+BvjMay1QkrQw3hkrSY0z6CWpcQa9JDXOoJekxhn0ktQ4g16SGmfQS1LjDHpJapxBL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxvYI+ycYkh5NMJdkxx/iNSQ4leSTJN5NcMDD2UpKHu6+JYRYvSZrfvL8zNskKYBdwGTANHEgyUVWHBqY9BIxX1S+T/CVwG/CRbuy5qrpwyHVLknrq88vBNwBTVXUEIMleYDPw66CvqvsH5j8AfHSYRerVrdlxz6hLkHQa63PqZiXw9EB7uus7meuAbwy0X5dkMskDSa6ea0GSbd2cyWPHjvUoSZLUV58j+t6SfBQYBz4w0H1BVR1N8lbgviSPVtUTg+uqajewG2B8fLyGWZMkLXd9juiPAqsH2qu6vldIcilwE7Cpqp7/VX9VHe3+PAJ8C7hoAfVKkk5Rn6A/AKxLsjbJ2cAW4BWfnklyEXAnMyH/04H+c5Oc022fB7yPgXP7kqTFN++pm6o6kWQ7sB9YAeypqoNJdgKTVTUBfA54A/DVJABPVdUm4J3AnUleZuaHyq2zPq0jSVpkvc7RV9U+YN+svpsHti89ybrvAe9eSIGSpIXxzlhJapxBL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BLUuMMeklqnEEvSY3rFfRJNiY5nGQqyY45xm9McijJI0m+meSCgbGtSR7vvrYOs3hJ0vzmDfokK4BdwJXAeuDaJOtnTXsIGK+q9wBfA27r1r4ZuAW4GNgA3JLk3OGVL0maT58j+g3AVFUdqaoXgL3A5sEJVXV/Vf2yaz4ArOq2rwDurarjVfUscC+wcTilS5L66BP0K4GnB9rTXd/JXAd84zWulSQN2VnDfLEkHwXGgQ+c4rptwDaA888/f5glSdKy1+eI/iiweqC9qut7hSSXAjcBm6rq+VNZW1W7q2q8qsbHxsb61i5J6qFP0B8A1iVZm+RsYAswMTghyUXAncyE/E8HhvYDlyc5t7sIe3nXJ0laIvOeuqmqE0m2MxPQK4A9VXUwyU5gsqomgM8BbwC+mgTgqaraVFXHk3yamR8WADur6vii7IkkaU69ztFX1T5g36y+mwe2L32VtXuAPa+1QEnSwnhnrCQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BLUuMMeklqnEEvSY0z6CWpcQa9JDXOoJekxhn0ktS4XkGfZGOSw0mmkuyYY/ySJD9MciLJNbPGXkrycPc1MazCJUn9zPvLwZOsAHYBlwHTwIEkE1V1aGDaU8DHgL+d4yWeq6oLh1CrJOk1mDfogQ3AVFUdAUiyF9gM/Droq+rJbuzlRahRkrQAfU7drASeHmhPd319vS7JZJIHklw914Qk27o5k8eOHTuFl5YkzafPEf1CXVBVR5O8FbgvyaNV9cTghKraDewGGB8fryWoSVJPa3bcM5L3ffLWq0byvi3qc0R/FFg90F7V9fVSVUe7P48A3wIuOoX6JEkL1CfoDwDrkqxNcjawBej16Zkk5yY5p9s+D3gfA+f2JUmLb96gr6oTwHZgP/AYcHdVHUyyM8kmgCTvTTINfBi4M8nBbvk7gckk/wncD9w669M6kqRF1uscfVXtA/bN6rt5YPsAM6d0Zq/7HvDuBdYoSVoA74yVpMYZ9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNc6gl6TGGfSS1DiDXpIaZ9BLUuMMeklqnEEvSY0z6CWpcQa9JDXOoJekxhn0ktS4XkGfZGOSw0mmkuyYY/ySJD9MciLJNbPGtiZ5vPvaOqzCJUn9zBv0SVYAu4ArgfXAtUnWz5r2FPAx4K5Za98M3AJcDGwAbkly7sLLliT11eeIfgMwVVVHquoFYC+weXBCVT1ZVY8AL89aewVwb1Udr6pngXuBjUOoW5LUU5+gXwk8PdCe7vr66LU2ybYkk0kmjx071vOlJUl9nBYXY6tqd1WNV9X42NjYqMuRpKb0CfqjwOqB9qqur4+FrJUkDUGfoD8ArEuyNsnZwBZgoufr7wcuT3JudxH28q5PkrRE5g36qjoBbGcmoB8D7q6qg0l2JtkEkOS9SaaBDwN3JjnYrT0OfJqZHxYHgJ1dnyRpiZzVZ1JV7QP2zeq7eWD7ADOnZeZauwfYs4AaJUkLcFpcjJUkLR6DXpIaZ9BLUuMMeklqnEEvSY0z6CWpcQa9JDXOoJekxhn0ktQ4g16SGmfQS1LjDHpJapxBL0mNM+glqXG9HlN8Jlmz455RlyBpCEb5vfzkrVeN7L0Xg0f0ktQ4g16SGmfQS1LjegV9ko1JDieZSrJjjvFzknylG/9+kjVd/5okzyV5uPu6Y7jlS5LmM+/F2CQrgF3AZcA0cCDJRFUdGph2HfBsVb09yRbgs8BHurEnqurCIdctSeqpzxH9BmCqqo5U1QvAXmDzrDmbgS92218DPpQkwytTkvRa9Qn6lcDTA+3prm/OOVV1Avg58JZubG2Sh5J8O8n753qDJNuSTCaZPHbs2CntgCTp1S32xdifAOdX1UXAjcBdSd44e1JV7a6q8aoaHxsbW+SSJGl56RP0R4HVA+1VXd+cc5KcBbwJeKaqnq+qZwCq6kHgCeAdCy1aktRfn6A/AKxLsjbJ2cAWYGLWnAlga7d9DXBfVVWSse5iLkneCqwDjgyndElSH/N+6qaqTiTZDuwHVgB7qupgkp3AZFVNAF8AvpRkCjjOzA8DgEuAnUleBF4Grq+q44uxI5KkuaWqRl3DK4yPj9fk5ORrXu+zbiSdqRbyjJ0kD1bV+Fxj3hkrSY0z6CWpcQa9JDXOoJekxhn0ktQ4g16SGmfQS1LjDHpJapxBL0mNM+glqXEGvSQ1zqCXpMYZ9JLUOINekhpn0EtS4wx6SWqcQS9JjTPoJalxBr0kNa5X0CfZmORwkqkkO+YYPyfJV7rx7ydZMzD2ya7/cJIrhle6JKmPeYM+yQpgF3AlsB64Nsn6WdOuA56tqrcDfw98tlu7HtgCvAvYCHy+ez1J0hLpc0S/AZiqqiNV9QKwF9g8a85m4Ivd9teADyVJ17+3qp6vqh8DU93rSZKWyFk95qwEnh5oTwMXn2xOVZ1I8nPgLV3/A7PWrpz9Bkm2Adu65v8lOdyr+lNzHvCzRXjdM4X77/4v5/2HM+DvIJ9d0PILTjbQJ+gXXVXtBnYv5nskmayq8cV8j9OZ++/+L+f9h+X9d9Dn1M1RYPVAe1XXN+ecJGcBbwKe6blWkrSI+gT9AWBdkrVJzmbm4urErDkTwNZu+xrgvqqqrn9L96mctcA64AfDKV2S1Me8p266c+7bgf3ACmBPVR1MshOYrKoJ4AvAl5JMAceZ+WFAN+9u4BBwArihql5apH2Zz6KeGjoDuP/L23Lff1jGfweZOfCWJLXKO2MlqXEGvSQ1blkFfZLPJfmvJI8k+Zckvz3qmpZSkg8nOZjk5STL5mNm8z3Co2VJ9iT5aZIfjbqWUUiyOsn9SQ51//Y/PuqaRmFZBT1wL/D7VfUe4L+BT464nqX2I+BPge+MupCl0vMRHi37R2YeP7JcnQD+pqrWA38A3LDM/vsDyyzoq+rfqupE13yAmc/1LxtV9VhVLcZdx6ezPo/waFZVfYeZT8ItS1X1k6r6Ybf9v8BjzHF3fuuWVdDP8hfAN0ZdhBbdXI/wWHbf6ILuqboXAd8fbSVL77R4BMIwJfl34HfmGLqpqv61m3MTM/9L9+WlrG0p9Nl/ablJ8gbg68BfV9UvRl3PUmsu6Kvq0lcbT/Ix4I+BD1WDNxHMt//LkI/hWOaS/CYzIf/lqvrnUdczCsvq1E2SjcDfAZuq6pejrkdLos8jPNSo7nHpXwAeq6rbR13PqCyroAf+Afgt4N4kDye5Y9QFLaUkf5JkGvhD4J4k+0dd02LrLr7/6hEejwF3V9XB0Va1dJL8E/AfwO8lmU5y3ahrWmLvA/4c+GD3Pf9wkj8adVFLzUcgSFLjltsRvSQtOwa9JDXOoJekxhn0ktQ4g16SGmfQS1LjDHpJatz/A/lZnXCO9y01AAAAAElFTkSuQmCC\n"
          },
          "metadata": {
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "stats.kstest(num_norm, 'norm')"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "Ze15ux2ZbBnF",
        "outputId": "d1db74ee-ff18-4ae5-e77b-b96fc2f44a87"
      },
      "execution_count": 58,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "KstestResult(statistic=0.08085802793543595, pvalue=0.5045661392249077)"
            ]
          },
          "metadata": {},
          "execution_count": 58
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "num_norm = get_num_norm(K = 100,T = 50)"
      ],
      "metadata": {
        "id": "my91q4OabnxK"
      },
      "execution_count": 62,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "plt.hist(num_norm, density=True);"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 265
        },
        "id": "ct29AlVUDF8_",
        "outputId": "5d62595b-02f2-4259-d884-397e14aa7600"
      },
      "execution_count": 63,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD4CAYAAAD8Zh1EAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAL3UlEQVR4nO3df6hf913H8edr6aLiioK9oDRht2gQLrV2cu0qgkpXIbWS+GuQgLJhIQgLVDbQSKVg/adboSoYcMEVRTZjrQ6DzciqVoZga2611qVZNJRoUoTe6dwcw9W4t3/kW/l6d2++J8v35tz77vMBge8558P9vg9pn5x8f5ybqkKStP29ZewBJEnzYdAlqQmDLklNGHRJasKgS1ITN431xLfcckstLi6O9fSStC298MILn62qhfWOjRb0xcVFVlZWxnp6SdqWkvzzRsd8yUWSmjDoktSEQZekJgy6JDVh0CWpCYMuSU0YdElqwqBLUhMGXZKaGO2boto+Fo88PcrzXnj0/lGeV9quvEKXpCYMuiQ1YdAlqQmDLklNGHRJasKgS1ITBl2SmjDoktSEQZekJgy6JDVh0CWpCYMuSU0YdElqwqBLUhMGXZKa8H7o18h7g0vaqrxCl6QmDLokNWHQJakJgy5JTRh0SWrCoEtSEwZdkpow6JLUhEGXpCYGBT3J3iTnkpxPcuQq634ySSVZnt+IkqQhZgY9yQ7gKHAfsAQcTLK0zrqbgQeB5+c9pCRptiFX6HcB56vqlap6HTgO7F9n3a8CHwT+a47zSZIGGhL0W4GLU9uXJvv+T5LvAXZX1VXvXJXkUJKVJCurq6vXPKwkaWPX/aZokrcAjwMfmLW2qo5V1XJVLS8sLFzvU0uSpgwJ+qvA7qntXZN9b7gZuB34yyQXgLuBE74xKkk31pCgnwb2JLktyU7gAHDijYNV9fmquqWqFqtqEXgO2FdVK5sysSRpXTODXlWXgcPAKeAs8GRVnUnySJJ9mz2gJGmYQb+xqKpOAifX7Ht4g7U/dP1jSZKuld8UlaQmDLokNWHQJakJgy5JTRh0SWrCoEtSEwZdkpow6JLUhEGXpCYMuiQ1YdAlqQmDLklNGHRJasKgS1ITBl2SmjDoktSEQZekJgy6JDVh0CWpCYMuSU0YdElqwqBLUhMGXZKaMOiS1IRBl6QmDLokNWHQJakJgy5JTRh0SWrCoEtSEwZdkpow6JLUhEGXpCYMuiQ1YdAlqYlBQU+yN8m5JOeTHFnn+M8l+YckLyb5qyRL8x9VknQ1M4OeZAdwFLgPWAIOrhPsj1XVd1XVncCHgMfnPqkk6aqGXKHfBZyvqleq6nXgOLB/ekFVfWFq8xuBmt+IkqQhbhqw5lbg4tT2JeCdaxcleR/wfmAncM9cppMkDTa3N0Wr6mhVfTvwi8Avr7cmyaEkK0lWVldX5/XUkiSGBf1VYPfU9q7Jvo0cB35svQNVdayqlqtqeWFhYfiUkqSZhgT9NLAnyW1JdgIHgBPTC5Lsmdq8H/in+Y0oSRpi5mvoVXU5yWHgFLADeKKqziR5BFipqhPA4ST3Av8NfA54z2YOLUn6akPeFKWqTgIn1+x7eOrxg3OeS5J0jQYFXeNbPPL02CNI2uL86r8kNWHQJakJgy5JTRh0SWrCoEtSEwZdkpow6JLUhEGXpCYMuiQ1YdAlqQmDLklNGHRJasKgS1ITBl2SmjDoktSEQZekJgy6JDVh0CWpCYMuSU0YdElqwqBLUhMGXZKaMOiS1IRBl6QmDLokNWHQJakJgy5JTRh0SWrCoEtSEwZdkpow6JLUhEGXpCYMuiQ1YdAlqQmDLklNDAp6kr1JziU5n+TIOsffn+TlJC8l+fMkb5//qJKkq5kZ9CQ7gKPAfcAScDDJ0pplfwcsV9UdwFPAh+Y9qCTp6oZcod8FnK+qV6rqdeA4sH96QVU9W1Vfmmw+B+ya75iSpFmGBP1W4OLU9qXJvo08AHxivQNJDiVZSbKyuro6fEpJ0kxzfVM0yU8Dy8Bj6x2vqmNVtVxVywsLC/N8akl607tpwJpXgd1T27sm+/6fJPcCDwE/WFVfns94kqShhlyhnwb2JLktyU7gAHBiekGSdwAfBvZV1WvzH1OSNMvMoFfVZeAwcAo4CzxZVWeSPJJk32TZY8DbgD9M8mKSExv8OEnSJhnykgtVdRI4uWbfw1OP753zXJKka+Q3RSWpCYMuSU0YdElqYtBr6NIYFo88PfYIN9yFR+8fewRtY16hS1ITBl2SmjDoktSEQZekJgy6JDVh0CWpCYMuSU0YdElqwqBLUhMGXZKaMOiS1IRBl6QmDLokNWHQJakJgy5JTRh0SWrCoEtSEwZdkpow6JLUhEGXpCYMuiQ1YdAlqQmDLklNGHRJasKgS1ITBl2SmjDoktSEQZekJgy6JDVh0CWpCYMuSU3cNGRRkr3AbwA7gN+uqkfXHP8B4NeBO4ADVfXUvAedtnjk6c388ZK0Lc28Qk+yAzgK3AcsAQeTLK1Z9i/Ae4GPzXtASdIwQ67Q7wLOV9UrAEmOA/uBl99YUFUXJse+sgkzSpIGGPIa+q3AxantS5N91yzJoSQrSVZWV1e/lh8hSdrADX1TtKqOVdVyVS0vLCzcyKeWpPaGBP1VYPfU9q7JPknSFjIk6KeBPUluS7ITOACc2NyxJEnXambQq+oycBg4BZwFnqyqM0keSbIPIMn3JrkEvBv4cJIzmzm0JOmrDfocelWdBE6u2ffw1OPTXHkpRpI0Er8pKklNGHRJasKgS1ITBl2SmjDoktSEQZekJgy6JDUx6HPoknob83cMXHj0/tGeuxuv0CWpCYMuSU0YdElqwqBLUhMGXZKaMOiS1IRBl6Qm/By6tIWM+XlwbX9eoUtSEwZdkpow6JLUhEGXpCYMuiQ1YdAlqQmDLklNGHRJasKgS1ITBl2SmjDoktSEQZekJgy6JDVh0CWpCYMuSU14P3RJb0pj3nv+wqP3b8rP9Qpdkpow6JLUhEGXpCYMuiQ1MSjoSfYmOZfkfJIj6xz/uiR/MDn+fJLFeQ8qSbq6mUFPsgM4CtwHLAEHkyytWfYA8Lmq+g7g14APzntQSdLVDblCvws4X1WvVNXrwHFg/5o1+4HfnTx+CnhXksxvTEnSLEM+h34rcHFq+xLwzo3WVNXlJJ8HvgX47PSiJIeAQ5PNLyY5N+O5b1n7Mxrpem5dzws8t02Rzf33/Jb8O7vOc377Rgdu6BeLquoYcGzo+iQrVbW8iSONpuu5dT0v8Ny2o67ntZEhL7m8Cuye2t412bfumiQ3Ad8E/Ns8BpQkDTMk6KeBPUluS7ITOACcWLPmBPCeyeOfAv6iqmp+Y0qSZpn5ksvkNfHDwClgB/BEVZ1J8giwUlUngI8Av5fkPPDvXIn+PAx+eWYb6npuXc8LPLftqOt5rSteSEtSD35TVJKaMOiS1MSWD3qSx5J8JslLST6e5JvHnmlekrw7yZkkX0my7T9aNesWEdtVkieSvJbk02PPMk9Jdid5NsnLk/8OHxx7pnlJ8vVJ/ibJ30/O7VfGnulG2PJBB54Bbq+qO4B/BH5p5Hnm6dPATwCfGnuQ6zXwFhHb1e8Ae8ceYhNcBj5QVUvA3cD7Gv2dfRm4p6q+G7gT2Jvk7pFn2nRbPuhV9cmqujzZfI4rn4NvoarOVtWsb8tuF0NuEbEtVdWnuPLprVaq6l+r6m8nj/8TOMuVb31ve3XFFyebb538af8JkC0f9DV+FvjE2ENoXevdIqJFHN4MJndIfQfw/LiTzE+SHUleBF4DnqmqNue2kS3xO0WT/Bnwresceqiq/mSy5iGu/BPxozdytus15NykMSV5G/BHwM9X1RfGnmdequp/gDsn77t9PMntVdXqfZC1tkTQq+reqx1P8l7gR4F3bbdvoM46t0aG3CJCW0ySt3Il5h+tqj8ee57NUFX/keRZrrwP0jroW/4llyR7gV8A9lXVl8aeRxsacosIbSGTW1x/BDhbVY+PPc88JVl44xNxSb4B+GHgM+NOtfm2fNCB3wRuBp5J8mKS3xp7oHlJ8uNJLgHfBzyd5NTYM32tJm9cv3GLiLPAk1V1Ztyp5iPJ7wN/DXxnkktJHhh7pjn5fuBngHsm/2+9mORHxh5qTr4NeDbJS1y52Himqv505Jk2nV/9l6QmtsMVuiRpAIMuSU0YdElqwqBLUhMGXZKaMOiS1IRBl6Qm/hePEKf2myVtGAAAAABJRU5ErkJggg==\n"
          },
          "metadata": {
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "stats.kstest(num_norm, 'norm')"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "CNSY_J7too9k",
        "outputId": "353df31d-0902-408b-9671-9a78ff3e2543"
      },
      "execution_count": 64,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "KstestResult(statistic=0.08810012705744186, pvalue=0.3965373977216896)"
            ]
          },
          "metadata": {},
          "execution_count": 64
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "# Checking the distribution?\n"
      ],
      "metadata": {
        "id": "4cLTd8LkwFJT"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "_, sum_array_0, nums_of_intervals = get_recov_processes(100, 100)"
      ],
      "metadata": {
        "id": "Lg-HW_khn1Z4"
      },
      "execution_count": 80,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "K = 100\n",
        "t = 100\n",
        "N_t = [0]*K\n",
        "max_N = [0]*K\n",
        "for i in range(K):\n",
        "  for j, elem in enumerate(sum_array_0[i]):\n",
        "    if elem <= t:\n",
        "      max_N[i] += 1 "
      ],
      "metadata": {
        "id": "ZOvT__0RlIoD"
      },
      "execution_count": 81,
      "outputs": []
    },
    {
      "cell_type": "code",
      "source": [
        "plt.hist(max_N, bins = 20);"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 265
        },
        "id": "e3W1kLpBdVr5",
        "outputId": "d3f5b310-731e-4239-8fb0-29b13e7399a9"
      },
      "execution_count": 82,
      "outputs": [
        {
          "output_type": "display_data",
          "data": {
            "text/plain": [
              "<Figure size 432x288 with 1 Axes>"
            ],
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAD4CAYAAAD1jb0+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAANZklEQVR4nO3df6xk5V3H8fdHtqXdlgiUK7bQ2yVK0JYIrddK21hNwcivuK3RZElQqsT7F5Y2xrqkiVUTExqr9TdmFbrYIpggWiK1BamVmhSUxS3dZaFgQQpdWCqmta3hR/v1jzkm18v9OefcmX32vl/JZGaeOTPn+9zn3s8988ycc1JVSJLa8x3TLkCSNB4DXJIaZYBLUqMMcElqlAEuSY3aMsmVnXDCCbVt27ZJrlKSmrdnz56vVNXM4vaJBvi2bdu4++67J7lKSWpekv9Yqt0pFElqlAEuSY0ywCWpUQa4JDXKAJekRhngktSoVQM8yTVJDiXZt6Dtt5Pcn+TeJH+T5NiNLVOStNhatsB3A+cuarsNOL2qfgD4AnDFwHVJklaxaoBX1R3A04vabq2q57u7dwInb0BtkqQVDLEn5i8Af7Xcg0nmgXmA2dnZAVanI9m2nbeM/dxHrrxgwEqkw1+vDzGTvA94HrhuuWWqaldVzVXV3MzMC3bllySNaewt8CTvBC4Ezi7PyyZJEzdWgCc5F3gv8KNV9c1hS5IkrcVavkZ4PfBZ4LQkjyW5FPgj4BjgtiR7k/zpBtcpSVpk1S3wqrpoiearN6AWSdI6uCemJDXKAJekRhngktQoA1ySGmWAS1KjDHBJapQBLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUqCFO6CAdFlo9GUSrdWv63AKXpEYZ4JLUKANckhplgEtSowxwSWqUAS5JjTLAJalRBrgkNcoAl6RGGeCS1CgDXJIaZYBLUqMMcElq1KoBnuSaJIeS7FvQdnyS25I82F0ft7FlSpIWW8sW+G7g3EVtO4Hbq+pU4PbuviRpglYN8Kq6A3h6UfN24Nru9rXA2weuS5K0inFP6HBiVR3sbj8BnLjcgknmgXmA2dnZMVe3OXmg/8np87MGf96ajt4fYlZVAbXC47uqaq6q5mZmZvquTpLUGTfAn0zySoDu+tBwJUmS1mLcAL8ZuKS7fQnwsWHKkSSt1Vq+Rng98FngtCSPJbkUuBL48SQPAud09yVJE7Tqh5hVddEyD509cC2SpHVwT0xJapQBLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUKANckhplgEtSowxwSWqUAS5JjTLAJalRBrgkNcoAl6RGGeCS1CgDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKAJekRhngktQoA1ySGmWAS1KjDHBJapQBLkmN6hXgSd6TZH+SfUmuT/KSoQqTJK1s7ABPchLwLmCuqk4HjgJ2DFWYJGllfadQtgAvTbIF2Ap8uX9JkqS12DLuE6vq8SQfBB4F/ge4tapuXbxcknlgHmB2dnbc1Wmdtu28ZeznPnLlBQNWImmj9JlCOQ7YDpwCvAp4WZKLFy9XVbuqaq6q5mZmZsavVJL0//SZQjkHeLiqnqqq54CbgDcPU5YkaTV9AvxR4KwkW5MEOBs4MExZkqTVjB3gVXUXcCNwD/D57rV2DVSXJGkVY3+ICVBV7wfeP1AtkqR1cE9MSWqUAS5JjTLAJalRBrgkNcoAl6RGGeCS1CgDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKAJekRhngktQoA1ySGmWAS1KjDHBJapQBLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUKANckhplgEtSowxwSWpUrwBPcmySG5Pcn+RAkjcNVZgkaWVbej7/94FPVNVPJ3kxsHWAmiRJazB2gCf5TuCtwDsBqupZ4NlhypIkrabPFvgpwFPAh5OcAewBLq+qbyxcKMk8MA8wOzvbY3WalG07b+n1/EeuvGCgSrSR+oyzY3x46DMHvgV4A3BVVb0e+Aawc/FCVbWrquaqam5mZqbH6iRJC/UJ8MeAx6rqru7+jYwCXZI0AWMHeFU9AXwpyWld09nAfYNUJUlaVd9vofwScF33DZQvAj/fvyRJ0lr0CvCq2gvMDVSLJGkd3BNTkhplgEtSowxwSWqUAS5JjTLAJalRBrgkNcoAl6RGGeCS1CgDXJIaZYBLUqMMcElqlAEuSY3qezTCTcEzl0g6HLkFLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUKANckhplgEtSowxwSWqUAS5JjTLAJalRBrgkNcoAl6RGGeCS1KjeAZ7kqCT/luTvhihIkrQ2Q2yBXw4cGOB1JEnr0CvAk5wMXAD8+TDlSJLWqu8JHX4PeC9wzHILJJkH5gFmZ2d7rq49fU4GoXY4zpqGsbfAk1wIHKqqPSstV1W7qmququZmZmbGXZ0kaZE+UyhvAX4yySPADcDbknx0kKokSasaO8Cr6oqqOrmqtgE7gE9V1cWDVSZJWpHfA5ekRg1yVvqq+jTw6SFeS5K0Nm6BS1KjDHBJapQBLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUKANckhplgEtSowxwSWqUAS5JjRrkWCiSpqPVE0n0qfuRKy8YsJK2uQUuSY0ywCWpUQa4JDXKAJekRhngktQoA1ySGmWAS1KjDHBJapQBLkmNMsAlqVEGuCQ1ygCXpEYZ4JLUqLEDPMmrk/xjkvuS7E9y+ZCFSZJW1udwss8Dv1xV9yQ5BtiT5Laqum+g2iRJKxh7C7yqDlbVPd3t/wYOACcNVZgkaWWDnNAhyTbg9cBdSzw2D8wDzM7Ojr0ODwAvaZr6njxjI3Ko94eYSV4O/DXw7qr62uLHq2pXVc1V1dzMzEzf1UmSOr0CPMmLGIX3dVV10zAlSZLWos+3UAJcDRyoqt8driRJ0lr02QJ/C/CzwNuS7O0u5w9UlyRpFWN/iFlV/wxkwFokSevgnpiS1CgDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKAJekRhngktQoA1ySGmWAS1KjDHBJatQgZ+Q53PU9k4bWx5/3kc8xPjy4BS5JjTLAJalRBrgkNcoAl6RGGeCS1CgDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXKAJekRhngktSoXgGe5NwkDyR5KMnOoYqSJK1u7ABPchTwx8B5wGuBi5K8dqjCJEkr67MF/kbgoar6YlU9C9wAbB+mLEnSavqc0OEk4EsL7j8G/PDihZLMA/Pd3a8neaC7fQLwlR7rb5F93hw2W58n2t98YFJrWtG6+9yz7tcs1bjhZ+Spql3ArsXtSe6uqrmNXv/hxD5vDputz5utv3D49LnPFMrjwKsX3D+5a5MkTUCfAP9X4NQkpyR5MbADuHmYsiRJqxl7CqWqnk9yGfBJ4Cjgmqrav46XeMG0yiZgnzeHzdbnzdZfOEz6nKqadg2SpDG4J6YkNcoAl6RGTSTAk7wnyf4k+5Jcn+QlSXYneTjJ3u5y5iRqmZQkl3f93Z/k3V3b8UluS/Jgd33ctOsc0jJ9/vUkjy8Y5/OnXWcfSa5JcijJvgVtS45rRv6gO9TEvUneML3Kx7fOPv9Ykq8uGO9fm17l41umzz/T/W5/O8ncouWv6Mb5gSQ/Mak6NzzAk5wEvAuYq6rTGX3guaN7+Feq6szusneja5mUJKcDv8hob9UzgAuTfC+wE7i9qk4Fbu/uHxFW6DPAhxaM88enVuQwdgPnLmpbblzPA07tLvPAVROqcWi7WXufAT6zYLx/c0I1Dm03L+zzPuCngDsWNnaHENkBvK57zp90hxrZcJOaQtkCvDTJFmAr8OUJrXdavh+4q6q+WVXPA//EaOC3A9d2y1wLvH1K9W2E5fp8RKmqO4CnFzUvN67bgb+okTuBY5O8cjKVDmedfT4iLNXnqjpQVQ8ssfh24IaqeqaqHgYeYrQhs+E2PMCr6nHgg8CjwEHgq1V1a/fwb3VvLT+U5OiNrmWC9gE/kuQVSbYC5zPa6enEqjrYLfMEcOK0CtwAy/UZ4LJunK850qaNOsuN61KHmzhpkoVtoJV+l9+U5HNJ/j7J66ZQ26RNbZwnMYVyHKP/UKcArwJeluRi4Arg+4AfAo4HfnWja5mUqjoAfAC4FfgEsBf41qJlCjhivsO5Qp+vAr4HOJPRP/DfmVaNk3CkjetaLOrzPcBrquoM4A+Bv51aYZvAJKZQzgEerqqnquo54CbgzVV1sHtr+QzwYSb0lmNSqurqqvrBqnor8F/AF4An/+8tdHd9aJo1Dm2pPlfVk1X1rar6NvBnHGHj3FluXI/kw00s2eeq+lpVfb27/XHgRUlOmF6ZEzG1cZ5EgD8KnJVka5IAZwMHFgx+GM2f7VvhNZqT5Lu661lGc8F/yehQA5d0i1wCfGw61W2Mpfq8aM73HRxh49xZblxvBn6u+zbKWYymDw8u9QINWrLPSb67+5smyRsZZcx/TqXCybkZ2JHk6CSnMPrQ+l8msuaq2vAL8BvA/Yz+eD8CHA18Cvh81/ZR4OWTqGVSF+AzwH3A54Czu7ZXMPrE/kHgH4Djp13nBPr8kW6c72X0i/7KadfZs4/XM5oKeo7RXOely40rEEYnPfn37mcwN+36J9Dny4D93e/AnYzebU+9DwP1+R3d7WeAJ4FPLlj+fd04PwCcN6k63ZVekhrlnpiS1CgDXJIaZYBLUqMMcElqlAEuSY0ywCWpUQa4JDXqfwGxKNHpaT6hwgAAAABJRU5ErkJggg==\n"
          },
          "metadata": {
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "stats.kstest(max_N,'norm')"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "aAG5RHG-lIq-",
        "outputId": "80bd5850-2a76-4589-e37d-8f3c988b495e"
      },
      "execution_count": 83,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "KstestResult(statistic=1.0, pvalue=0.0)"
            ]
          },
          "metadata": {},
          "execution_count": 83
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "from scipy.stats import poisson\n",
        "\n",
        "K = len(max_N)\n",
        "num_interv = 6\n",
        "\n",
        "mu=statistics.mean(max_N)\n",
        "std_1=statistics.stdev(max_N)\n",
        "pmax=max(max_N);\n",
        "pmin=min(sorted(max_N)[1:]);\n",
        "dis=(pmax-pmin)/(num_interv-1);\n",
        "num=[0]*num_interv\n",
        "\n",
        "i=0;\n",
        "while i<num_interv:\n",
        "  for j in range(i, K):\n",
        "    if (max_N[j]>=pmin) and (max_N[j]<pmin+dis):\n",
        "      num[i]=num[i]+1\n",
        "  pmin=pmin+dis\n",
        "  i=i+1\n",
        "\n",
        "pmin=min(max_N);\n",
        "\n",
        "F = [0]*num_interv\n",
        "F[0]=poisson.cdf(pmin,mu)\n",
        "F[1]=poisson.cdf(pmin+dis,mu)\n",
        "F[2]=poisson.cdf(pmin+2*dis,mu)\n",
        "F[3]=poisson.cdf(pmin+3*dis,mu)\n",
        "F[4]=poisson.cdf(pmin+4*dis,mu)\n",
        "F[5]=poisson.cdf(pmax,mu)\n",
        "\n",
        "\n",
        "p = [0]*num_interv\n",
        "p[0] = F[0]\n",
        "for i in range(1,num_interv-1):\n",
        "  p[i]=F[i]-F[i-1]\n",
        "p[5] = 1 - F[5]\n",
        "\n",
        "k = 0\n",
        "for i in range(num_interv):\n",
        "  k = k+((num[i]-K*p[i])**2)/(K*p[i])\n",
        "k"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "rL6WzWzelItd",
        "outputId": "87a99200-4986-4b40-9615-4974bfa1274a"
      },
      "execution_count": 84,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "48.36742150484645"
            ]
          },
          "metadata": {},
          "execution_count": 84
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "The deviation is too large"
      ],
      "metadata": {
        "id": "kO-1En2Xd6b1"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "num"
      ],
      "metadata": {
        "id": "zDGhJEUinTE2",
        "outputId": "ce4194c5-65e9-4a2b-df03-5f9b4ae6d6ce",
        "colab": {
          "base_uri": "https://localhost:8080/"
        }
      },
      "execution_count": 85,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "[17, 24, 36, 18, 2, 1]"
            ]
          },
          "metadata": {},
          "execution_count": 85
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "p"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "80y0c4eLmHLK",
        "outputId": "abf5a793-512c-432e-c6a7-b12945fa428e"
      },
      "execution_count": 86,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "[0.14119067277253525,\n",
              " 0.11457821804159454,\n",
              " 0.18785477052085547,\n",
              " 0.19843873434045461,\n",
              " 0.16342880825229,\n",
              " 0.07381032598558801]"
            ]
          },
          "metadata": {},
          "execution_count": 86
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "sum(p)"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "qHWk29SjlIwT",
        "outputId": "944c923c-65c0-4a2a-bb0f-4ab8a54585b8"
      },
      "execution_count": 87,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "0.8793015299133179"
            ]
          },
          "metadata": {},
          "execution_count": 87
        }
      ]
    }
  ]
}