mc.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "2cf93199",
   "metadata": {},
   "source": [
    "<header><font size=\"+2\"> <center> <b> Solving the Ising model with the Metropolis Algorithm </b></center></font></header>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "5615bbd0",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import multiprocessing as mp\n",
    "import itertools as it"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "9faaeeb6",
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "plt.rcParams[\"figure.figsize\"] = (6, 6)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "f7b3e5ca",
   "metadata": {},
   "outputs": [],
   "source": [
    "Tc = 2.269185\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1c9250cf",
   "metadata": {},
   "source": [
    "(I know these could be expressed as matrix ops, so far they're fast enough, though)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "a3b5f2dc",
   "metadata": {},
   "outputs": [],
   "source": [
    "def p(ind, N):\n",
    "    if ind == -1:\n",
    "        return N - 1\n",
    "    return ind % N"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "c6b8457e",
   "metadata": {},
   "outputs": [],
   "source": [
    "def energy(s, J=1):\n",
    "    E = 0\n",
    "    N = s.shape[0]\n",
    "\n",
    "    for i in range(N):\n",
    "        for j in range(N):\n",
    "            dE = s[i, j] * (\n",
    "                s[p(i+1, N), j] + s[p(i-1, N), j] \n",
    "                + s[i, p(j+1, N)] + s[i, p(j-1,N)]\n",
    "            )\n",
    "            E -= 0.5 * J * dE\n",
    "\n",
    "    return E\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "672840a6",
   "metadata": {},
   "source": [
    "# Use burn - in"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "3461f614",
   "metadata": {},
   "outputs": [],
   "source": [
    "def MCstep(s, N, E, T, J, nsteps):\n",
    "    for itr in range(nsteps):\n",
    "        indx = np.random.randint(0, N)\n",
    "        indy = np.random.randint(0, N)\n",
    "\n",
    "        dE = 2 * J * s[indy, indx] *  (  s[p(indy + 1, N), indx]\\\n",
    "                                       + s[p(indy - 1, N), indx]\\\n",
    "                                       + s[indy, p(indx - 1, N)]\n",
    "                                       + s[indy, p(indx + 1, N)])\n",
    "\n",
    "        # accept\n",
    "        if (dE < 0) or np.random.rand() <= np.exp(-1 * (dE / T)):\n",
    "            \n",
    "            s[indy, indx] = s[indy, indx] * -1\n",
    "            E = E + dE\n",
    "\n",
    "    return s, E\n",
    "\n",
    "def simulateBI(N, T, J = 1, burn_in = 3000, nsteps = 1000, n_outputs = 1000):\n",
    "    ss = []\n",
    "    Es = []\n",
    "    Ms = []\n",
    "    Mavgs = []\n",
    "    \n",
    "    # random init\n",
    "    s = 2 * (np.random.uniform(size=(N, N)) > 0.5) - 1\n",
    "    E = energy(s, J)\n",
    "\n",
    "    # burn-in\n",
    "    s, E = MCstep(s=s, N=N, E=E, T=T, J=J, nsteps = burn_in)\n",
    "    \n",
    "    #explore equlibrium state    \n",
    "    for itr in range(n_outputs):     \n",
    "        \n",
    "        s = np.copy(s)\n",
    "        ss.append(s)\n",
    "        \n",
    "        Es.append(E)\n",
    "        Ms.append(np.sum(s))\n",
    "        Mavgs.append(np.mean(s))\n",
    "        \n",
    "        s, E = MCstep(s=s, N=N, E=E, T=T, J=J, nsteps = nsteps)\n",
    "    \n",
    "    return ss, Es, Ms, Mavgs "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "id": "393730a3",
   "metadata": {},
   "outputs": [],
   "source": [
    "ss, Es, Ms, Mavgs = simulateBI(64, T=Tc, J=1, burn_in = 15000, nsteps = 1000, n_outputs=2000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "d31c479f",
   "metadata": {},
   "outputs": [],
   "source": [
    "np.array(ss).tofile(\"mcdata/new_64_\" + str(1) + \".dat\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "81eef1db",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAW4AAAFkCAYAAAAAOdHpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAAsTAAALEwEAmpwYAAAUJElEQVR4nO3dzXEbSZrG8efdWAOw1OqqjUB7QFIesD1gjwkcD2ZiTJA8EE2YoAekByT7tNemBRMaXveUe2CyBUIAUYVEZuVT+P8iGALxQWTWx6OMeiurIqUkAICP/5i6AQCAcQhuADBDcAOAGYIbAMwQ3ABghuAGADP/WfsLIqLofMOzs7NDNeUgHh8fNz5fo53bvmuI3pbbNiV9HKPV8ijtT2k7ayxPl2VXw5T70ePj479SSh83vRa1z+MuDe7ezjOPiI3P12jntu8aorfltk1JH8dotTxK+1PazhrL02XZ1TDlfhQRjyml802vcagEAMwQ3ABghuAGADPVi5NnZ2d6eHh489yYY1lD39vqWJTLMciWx+I36fF45Sab2umyLbXisi6PCSNuADBDcAOAGYIbAMwQ3ABgpnpxspUxBZQpi0I1ipsUj3YrLYiPKQBjPyzj4RhxA4AZghsAzBDcAGCG4AYAM91fHbCGoUWQHme29XY1PYpHh1U6C3bK7aPVCQJTb3MN9w2uDggAc7HzdMCIOJW0lKSU0k1+7lLSs6RlSum6ZgMBAG8NGXH/Iwf2SUQsc2grpXQnSRFxUbOBAIC33g3uiLiSdB8Ry5TSdUrpSdJnSU/5LU+STiu3EQCwYtehkl/yv98j4pukv0tarL3nw/qHcuBfFbfuAIYWEnosRA7FjLO3elweQ9tUWiTvse/HoPUyHjLl/Y+U0nNEPOoljJ8lnbz3gXzc+1rq86wSAHC2K7jv9SOkF3oJ7Sf9GHUvJd1WaBcAYIt3j3HnouTitQCZj3PfSFrm5xavRUoAQBuzn4DjfOx6k1ZXrnOegNPjcd6SNrnfLm/o54dyXpcjbZ2AM5vLum7jMiOyRI2NptVyq7ET9BhUJW2a+v6hm0wZntv63dv9aWti5iQAmCG4AcAMwQ0AZghuADDTfXFyDoWE3pUs4x7PeKiBdvajRvG6dDtuvR8w4gYAMwQ3AJghuAHADMENAGaqFyfPzs708PCw9+dbTd0e+t01iqVTF5RcZpe2WvZzm469ydTtLDHHmZNjZ+Yy4gYAMwQ3AJghuAHADMENAGYIbgAww5T3zOXMChdDl2ePy9hlW3A+M6REjRs29HjN8/cw4gYAMwQ3AJghuAHADMENAGa6v8u7yxRzl3Zu0msBpoW5TXl30Wpq/tCp5GOm0Tc8YWLrXd4ZcQOAGYIbAMwQ3ABghuAGADOTXI/bZZbS1O0Ze43eff/mMaMQeThTb1tDv39MO6fu0zaMuAHADMENAGYIbgAwQ3ADgJlJLus69QH/qb9/Xcsil/MNVYdqNdttzN+kkLmfKS/rvE0P+wYjbgAwQ3ADgBmCGwDMENwAYKary7pOedC/VUHLuUhVoxjXQ6Gnhd7W+5jLmJb+zUNzmXm9zYhlzGVdAWAuCG4AMENwA4AZghsAzEwyc3Kb3gpaU95vDp56K0JuU6OdpftLaZuOaX9lxA0AZnYGd0T8OyJuI+JvK89dRsRFRFzVbR4AYN2QEfdvKaVfU0pfpZfQlqSU0l3+/aJi+wAAa4YE9yIiliu/f5b0lB8/STo9eKsAAFsNCe4TSd8j4lv+fbH2+of1D0TEVUQ8RMTDp0+flFIa9DNURPz008qm7x7zM0bJMurRpv5MuS5R3yH2g7k5RB7uDO6U0nVK6VnScz5M8qyXMN/1mfOU0vnHjx+H9wgAsNO7wZ1HzuuHQu71Y9S9lHRboV0AgC12jbj/Kb0pSN6klG4kLXNRcvFapAQAtPHuBJx8iOT3/HOz8vzX/JDQBoDGupo5uemA/NBChts94/Y1tO1Tz4wbui57XBcu7XR27AXKUsycBAAzBDcAmCG4AcAMwQ0AZrq652QrJUXQGlwKXzWWUY/3KSy5vPDURbfe2lR6WVeXfaPUlnXEPScBYC4IbgAwQ3ADgBmCGwDMdDVzcpMaxZapC0ioa8z6LdkWXLajKQuWY4rCGI4RNwCYIbgBwAzBDQBmCG4AMENwA4CZ6meVnJ2d6eHh4c1zpRXt3qb1blOjcj7ltHOX5T6lbcty6rM4HBzzmSZj9y1G3ABghuAGADMENwCYIbgBwExX1+MuaYvztaJL9ViQmvo62+iby77VCtfjBoCZI7gBwAzBDQBmCG4AMNPV9bhLbhbqMmOt9OapLkraPmYZuczmPPRNibd93mV5tDLXGxAz4gYAMwQ3AJghuAHADMENAGa6Kk7WKBr0Vog45kJRKy4F4NJCosPleKfe/6b+/k0OsewZcQOAGYIbAMwQ3ABghuAGADNdFSdLZjnVKOrMUW+FqtLvrjFjdejnx7R96Lbtsm26tLNHh9gPGHEDgBmCGwDMENwAYIbgBgAzkxQnx1ym0uWyjEMLDqVtb1XMq8H5u0s/3+M2C1+MuAHAzODgjogvK48vI+IiIq7qNAsAsM2g4I6IC0nL/PhSklJKdyuvAQAa2RncEbGU9LTy1OeV358knVZoFwBgiyHFyWVK6W6lOLNYe/3D+gfyIZSth1FqFKmOpfjjfE/BGrNgW3Fe7s5cTk4Y4xDbzbvBHREXr4dEVjxLOnnvcymla0nX+W94L2UA6MyuEff3fAx7IWkZEaeS7vVj1L2UdFutdQCAn7x7jDul9HsecZ8oh3VK6UYvIX4habFhRA4AqChqHy+qcajE5RhXq+NzLsdaj/kYd2/brMuEoiM/xv2YUjrf9EJXl3UdasqVuW2hO1+is5XS9VYyO7VVQXzb98wtgFr1p8dl1MPAgpmTAGCG4AYAMwQ3AJghuAHAjGVxcpNWFf4xn+3xsq69Kb1vY+nfLPnuMQVp50sWDzW3/kj93o+VETcAmCG4AcAMwQ0AZghuADAzm+Lk1GoUztwLO+u4NOpuUxfZWR/TGLvcGXEDgBmCGwDMENwAYIbgBgAzR1mcLCnAtCwYOheKhi4nl4Jlj4Vih+3YvRjf6/bJiBsAzBDcAGCG4AYAMwQ3AJghuAHAzFGeVVKih4ry3JVOxy69ZnoNU56d0OOZEVOeQdLqptU1MeIGADMENwCYIbgBwAzBDQBmKE6iO61u9nssavS9t+KiVNYmt+2DETcAmCG4AcAMwQ0AZghuADBDcbJjPczQQj1Di2k9rvPS2Yclerxud2uMuAHADMENAGYIbgAwQ3ADgJnqwX12dqaU0t4/U5q6PRHx0w92Y7kd1tT7wbpN63fbT6ne+v6KETcAmCG4AcAMwQ0AZghuADDT/cxJZg96mts6Ku1PSVFrzD04S/eXoe3spUi3S41lVOoQM2YZcQOAmZ0j7oi4yA9/TSn9PT93KelZ0jKldF2veQCAde+OuCPiVNJpSulO0mlELHNoKz+3GuwAgAbeDe6U0u8ppa8RsZD0lFJ6kvRZ0lN+y5Ok07pNBACsGlqcPJf0R368WHvtw/qbI+JK0pUkffr0ad+2bTV1caFEq6KOy/IoNXRbGFPgO3R7ainZD1yKizWMWee97keDipP5sMgvK8e2T3a8/zqldJ5SOv/48WN5KwEAf9p1jPtLHj1LPwL7Xj9G3UtJt7UaBwD42a4R9zdJT7kAucgj6RtJy5Xn7qq3EgDwp3ePcedi5Gsh8m7l+a/rzwEA2uh+5mQNre71N2UBqEYxzrko7NLOUs7r6Fgc4n6dzJwEADMENwCYIbgBwAzBDQBmZlOcLC0Ezm3GWY2CVI9FrtL11mOfDs1lm92k1foZs4x6OGmBETcAmCG4AcAMwQ0AZghuADAzm+JkKZfLrZbeE5Bi3O73HmJmG8Zz2TZ72BYYcQOAGYIbAMwQ3ABghuAGADMUJyuaW7HFpT8u7cRuXKZ2M0bcAGCG4AYAMwQ3AJghuAHADMENAGY4q6RjNaZjH3q67rFMrZ/b9dqnNnQ7ntt2dCiMuAHADMENAGYIbgAwQ3ADgBmKk2bmNj29hpY3fh369yhavsXyKMOIGwDMENwAYIbgBgAzBDcAmLEsTrrczPUYZoKVFuOmXh7HsI6ORck2N2ad95A1jLgBwAzBDQBmCG4AMENwA4AZy+Jkj2oUtEqKsKUFlNL+lHy+9FKxFBfno0YhsIfi4rqx2ywjbgAwQ3ADgBmCGwDMENwAYIbi5IHMbQZeq5mPPRaKhmrV9h6Xscvs5aHcLsfLiBsAzLw74o6IhaRl/vmcUvp7fv5S0rOkZUrpunIbAQArdo24/yLpPKV0I0kRcZVDWymlu/zcRd0mAgBWvRvcKaXrlRH1UtKTpM/5X+V/T+s1DwCwblBxMiKWkr6nlO4i4re1lz9seP+VpKuV3wc1ptWsQOzWqrA6psjVWwHYufhcymUfnHIGcE1Dzyq5TCn9NT9+lnTy3pvzKP1akiLCYw0DgImdZ5VExGVK6Wt+fCrpXtIiv7yUdFutdQCAn7wb3Lnw+CUiHiPiUdJJLlQu82uL1yIlAKCNqH2sasyhEpfjZkPVOD7mMumjxJg+9noMsidz26/GMN8+HlNK55teYObkSD3em67HjbO3QuKxOOaQPibMnAQAMwQ3AJghuAHADMENAGa6Kk4eeoZlDaX3QyxFgW8+jqGQ2GMxfw4YcQOAGYIbAMwQ3ABghuAGADNdFSc3cSlYtLqv39w4z7B0bvs2vV0mt9X+P/VJB2Mx4gYAMwQ3AJghuAHADMENAGa6L06WcpiNWarH61f3WtQ5pBr3M3TZDlsVZp2XUU2MuAHADMENAGYIbgAwQ3ADgBmCGwDMTHJWicv00hoV7d76WKrluqxxJsPQvznH6e1D9Xj9+0NzW5eMuAHADMENAGYIbgAwQ3ADgJnqxcmzszM9PDzU/pqtSoobUxciS75/yu8u/XzLYtjQ5dSq2NrKmP4MLda2uib9mO8Z+l6KkwCAqghuADBDcAOAGYIbAMzM/nrcrbS6FrGLY7gO+hhzLZId0pR9n3oG8NjvZsQNAGYIbgAwQ3ADgBmCGwDMUJwcieLRW8dSXCzlvN1MOYN3zN+ssS2WXM635r7BiBsAzBDcAGCG4AYAMwQ3AJiZpDhZ45KSOKypLwuLeXC5v2yp1ts7I24AMLMzuCPiMiJuNzx3ERFX9ZoGANhkZ3CnlG5Wf4+Iy/z8Xf79ok7TAACb7HOo5LOkp/z4SdLp4ZoDANhln+LkYu33D+tvyIdQriTp06dPe3zF+8YULA9d3DyWYsuxOIbLrZbM/hvzXpeCtPO6fLXPiPtZ0sl7b0gpXaeUzlNK5x8/ftyrYQCAzfYJ7nv9GHUvJd1ufysA4NCGnFVyIel8pSh5I2mZn1+8FikBAG3sPMadg/m/1p77mh8S2gDQWNQuKEREk4pFq8LIHAobh1S63McsT4fi15hLjh7DtjRmnbWarTvlfSTHiIjHlNL5pteYOQkAZghuADBDcAOAGYIbAMzM5p6TpUWuuRWKSmeStloeDgXHbWrMPpybMcXaKWd4uhXJGXEDgBmCGwDMENwAYIbgBgAzsylOlprb5T25Z2RfjqEgXqr08rFDPz+H5c6IGwDMENwAYIbgBgAzBDcAmKlenDw7O9PDw0Ptrxk1Q6uE+z0nKToezpiCY2/bR+l2XGM7ajWjsbd1sQ9G3ABghuAGADMENwCYIbgBwMxsZk46F91aXoJ16N/ENJhh+VbLe5o6YcQNAGYIbgAwQ3ADgBmCGwDMENwAYKars0qmnG7bo6H9bLU8OCPlsFotuzHf09v6nOM2d4g+MeIGADMENwCYIbgBwAzBDQBmuipO9lZ0mOt02X3VWD9zKz6NuXwB29dbLjf7nfr7JUbcAGCH4AYAMwQ3AJghuAHATNQuBEXE4C+YsihVOmuzVcHCuXC3yZjipEMhs4fCVYmh1wPnOtlNPKaUzje9wIgbAMwQ3ABghuAGADMENwCY6Wrm5JTFp6FFmRpFlRqFnt6KdlKdmXEufXfRankO/Zs9FjFbzSB+DyNuADCz14g7Ii4lPUtappSuD9oiAMC7Ro+4c2grpXSXf784dKMAANvtc6jks6Sn/PhJ0unhmgMA2GWfQyWLtd8/rL8hIq4kXeVf/0/S/+7xPa9/a9+P1vLfkv51yD9Yo48j/+bB+zSliJhVfzTx+qm0D9quoy3Lo0Z//mfbC/sE97Okk/fekI97X0tSRDxsm7bpaG79kebXJ/rTv7n1qXV/9jlUcq8fo+6lpNuDtQYAsNPo4E4p3Uha5qLk4rVICQBoY6/TAVNKX/PDIaE9t9MF59YfaX59oj/9m1ufmvan+mVdAQCHxcxJADBTNbgj4jIiLvLpgZZyH243PGfZr4hYRMRp7sOXleed+3SRf2bRn1dz6k9E/DsibiPibyvP2fZpZR+6XHmuWX+qBfdcZljmYuyfZtCvv0g6f+1XRFw59ykiTiWd5rafRsTSuT+vcpuX+bF9fyT9llL69bU+NoM+/SPvQydTbHM1R9xznWFp3a+U0vXK9WWWeumDbZ9SSr+nlL5GxELSU0rJuj+SFBGv6+WVdX+yRe7XK9s+5RH1fUQs8/7UfJurGdyLtd9/mmFparH2u2W/8k70PY8QFmsvO/bpXNIf+fFi7TW3/ixzGLxarL3u1h/pZdLe94j4ln9frL3u1Kdf9NLe7xHxLQ8aFmvvqdqfmsH9rB0zLE09ax79ukwp/TU/fpZ5n/J/QL+sXLnSsj8RcbFhbsSzTPvzKo9MnyU9u6+j7I/cn0e9XN7jWQ37UzO45zrD0r5fEXG5cqzxVMZ9iogvK8WgZ73sPLb90cso7iKH29J9/Uh/1lHWDx049+l+5fFCL9td0/5UC+65zLDM7T9fKT5Y9yu3+0tEPEbEo6QT8z59k/S00vZr5/7kY/Z3evkPaJGfs+1P9k/pTUHyxrlPue2L1wLkFNscE3AAwAwTcADADMENAGYIbgAwQ3ADgBmCGwDMENwAYIbgBgAzBDcAmPl/yiHFJdFWrOEAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.pcolormesh(ss[100], cmap=\"binary\");"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "36d9da1e",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fc490790f70>]"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX8AAAFkCAYAAADIVTs/AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAAsTAAALEwEAmpwYAAAqIUlEQVR4nO3deXhTVf4/8PfJ0p0S2gKCUEoBRWXREkBABEZwgfE7o6KO4yxfHa3r+HVGFMV9R3TGcZf+ZkbHdQRc2aGURUGwCzuCQClLC3Tf1yTn90duQpKmbWra3rbn/Xoenqc5SW4+TS/vnJx77rlCSgkiIlKLQe8CiIio4zH8iYgUxPAnIlIQw5+ISEEMfyIiBTH8iYgUZNK7gEDFxcXJhIQEvcsgIupSMjMzC6WUvX3bu0z4JyQkICMjQ+8yiIi6FCHEUX/tHPYhIlIQw5+ISEEMfyIiBTH8iYgUxPAnIlIQw5+ISEEMfyIiBTH8iYgUxPAnIlIQw5+ISEEMfyIiBTH8iYgU1O3DP+tYCb7cfkLvMoiIOpVuH/5VdTbkltToXQYRUafS7cPfaBCwOaTeZRARdSrdPvxNBgPsDH8iIi/dPvzZ8yciaqzbh7/JINjzJyLy0e3D32gQsNkZ/kREnrp9+JuMAnaHQ+8yiIg6le4f/hzzJyJqpNuHv9Fg4LAPEZGPbh/+7PkTETXW/cOfY/5ERI10+/DnPH8iosa6ffjzDF8iosa6ffiz509E1Fi3D3+e4UtE1Fi3D3/2/ImIGgs6/IUQJUKItUKIh1pomy2EmC6ESG6ura05e/6c7UNE5MnUBtu4XkqZ2lybEGI2AEgpU4UQyUKI6QAsvm1+thM0ru1DRNRYWwz7WIQQiS20jQWQrf2cDSCpibY2J4QAo5+IyFtbhH8MgGIhxMJm2iw+z4ltos2L9o0gQwiRUVBQ0AalEhEREMCwjzZkE+PTnO0aopFSpmiPKxVCzJZSLvFtA1DqZxv+2rxo20kBAKvV+rM78OLnPpGIqJtqMfyllEuauk87UJshpcxqrg1AOs709BMBrNVu+7YREVEHCHbYZxHgdUB3ib82rT3RdaBXSpnqry3IWoiIKEBBzfaRUpYCyNL+LWmqTWtfoP2Y2lwbERG1v25/khcRETXG8CciUhDDn4hIQQx/IiIFKRH+PMOXiMibEuFPRETelAh/nuFLRORNifAnIiJvDH8iIgUpEf4GIeDg1byIiNyUCH+TkZdyJCLypET4G3kRdyIiL0qEv8kgYON1fImI3JQIf17Hl4jImxLhbzIa0MCePxGRmxLhb2bPn4jIixLhbzIaGP5ERB6UCH+zUXDYh4jIgyLhz54/EZEnJcLfZDCgwc6ePxGRixLhbzYKhj8RkQclwp/LOxAReVMj/DnsQ0TkRYnwDzEZ0MADvkREbkqEv8kgYGPPn4jITY3wN7LnT0TkSYnwNxu5qicRkSclwt9k4EleRESelAh/s1GgnmP+RERuioQ/e/5ERJ6UCH8Tx/yJiLyoEf4GzvYhIvKkRPibjZznT0TkSZHw5/IORESelAh/k1Fw2IeIyIMS4W82GnjAl4jIgxLhb+IF3ImIvCgR/mau7UNE5EWh8OewDxGRixLhz5O8iIi8KRH+Zp7kRUTkRYnwN/EkLyIiL8qEfwMv4E5E5KZE+IcYDWiwsedPROSiRPiHmoxcz5+IyIMS4R9iMqCugeFPROSiRPgbDQJ2yTF/IiIXJcKfiIi8KRP+Qu8CiIg6EWXCn4iIzjAFuwEhRBKARACQUi7R2mYDKAWQKKVMaU0bERG1v7bo+T+ihX6MECJRC3RIKVMBQAgxPdC2NqilSTzcS0R0RlDhL4RIBpAuhEiUUqZIKbMBjAWQrT0kG0BSK9qIiKgDBNvzHwIgFkCxEGKhEMICwOLzmNhWtHkRQiQLITKEEBkFBQVBFWoQgOR0TyIiAAGM+WvDMzE+zdmu4RoAh6WUpUKITADJcI7h+z4+0DYv2nGAFACwWq1BJXeIyYg6mwNhZmMwmyEi6hZaDH/XQdwmpONMgFvgDPRsnOnVJwJYq90OpK3dhJoMDH8iIk1Qwz7aB4PFdbBWG/dfAiBRa7NIKVMDbQvyd2mWM/zt7fkSRERdRtBTPaWUC7QfU9uirb2Emoxc34eISKPMSV6hZvb8iYhclAn/MJMRtez5ExEBUCj8nT1/hj8REaBS+POALxGRmzLhbzYakPZjvt5lEBF1CsqE/9Wj+qO4ql7vMoiIOgVlwj8sxIA6XseXiAiASuFvNqKugWP+RESASuHPqZ5ERG7KhL/ZKNDAYR8iIgAKhb8QvIovEZGLMuFPRERnMPyJiBTE8CciUpBS4c+LOBIROSkV/kRE5KRU+BsEYHew/09EpFT4h5mNXNmTiAiqhT/P8iUiAqBa+JsNqOX6PkREqoW/keFPRAQlw5/DPkRESoV/qNmAWh7wJSJSK/ydB3wZ/kRESoV/qNmAOg77EBGpFf7s+RMROakV/mYj6mzs+RMRKRb+nOdPRAQoF/4c9iEiApQLfwNqOexDRKRY+POALxERAMXCP5Rn+BIRAVAs/HnAl4jISbHw53r+RESAYuEfajJw2IeICIqFP6d6EhE5KRf+PMOXiEi18DfxgC8REaBY+JuMBtgdUu8yiIh0p1T4ExGRE8OfiEhBDH8iIgUpF/4c8SciUjD8iYiI4U9EpCTlwt9kEGiw80QvIlKbcuHPs3yJiJQMf57lS0SkXvjzal5EROqFP6/mRUQEmILdgBAiCUAiAEgpl2htJQAyAKyVUi7Q2mYDKAWQKKVMaaqtvXHYh4iobXr+j2ihHyOESNTarpdSzvAJfkgpU7Xb0/21tUEtLeLVvIiIggx/IUQygHQhRKKUMkVKma3dZfH4IACAsQBc92UDSGqird05x/w57ENEagu25z8EQCyAYiHEQiGERWuPcbVpty0+z4ttoq3dhXLYh4io5TF/bXgmxqc52zVcA+CwlLJUCJEJIBnAAo8x/VKPcX3fbfhr833tZG2biI+Pb6nUgITxOr5ERC2Hv+sgbhPScSbALQBKtcDOkFJm+TzOov2cCGCtdtu3zfe1UwCkAIDVam2TNdl4HV8ioiCHfbQPBovrYK0W1osAr4O8S7THJWqPs0gpU/21BVNLoHiGLxERIKTsGoscW61WmZGREfR21u47jZzCKtx+aWLLDyYi6uKEEJlSSqtvu3InecVEhqCwqk7vMoiIdKVc+PeOCkVRZb3eZRAR6Uq58I+NCkFRJXv+RKQ25cI/IsSIqnrO9iEitSkX/kIIvUsgItKdcuHv0lVmORERtQclwz8q1MShHyJSmpLhHxvJg75EpDYlwz863IzyGpveZRAR6UbJ8A8xGVBv57APEalLyfAPNRlQx5U9iUhhSoZ/iMmAOjvDn4jUpWb4Gw2o58qeRKQwJcM/1MTwJyK1KRr+XNOfiNSmZPiHsOdPRIpTOPw51ZOI1KVm+BsNHPYhIqUpGf6hZoY/EalNyfDvFRGC0uoGvcsgItKNkuHfu0coCriwGxEpTMnwj+GqnkSkOCXD32w0wObgxVyISF1Khj8RkeqUDX/ndE/O9SciNSkb/nFRISisrNe7DCIiXSgb/jGRoShm+BORopQN/57hZpTXcq4/EalJ2fA3mwTqeUEXIlKUsuEfYjTAZud0TyJSk7LhbzYa0MCePxEpiuFPRKQghcNf8IIuRKQsZcPfaBA4mF+pdxlERLpQNvzHJsTgdHmt3mUQEelC2fDv3SMUW7OLICVn/BCRepQN/zCzEafL61DTwPV9iEg9yoY/AFw40ILKOpveZRARdTilw/+8fj1QXsMlHohIPUqH/4UDLfjhSIneZRARdTilwz8pvhd255bpXQYRUYdTOvwHxUbiaFGV3mUQEXU4pcM/xGSAndfyJSIFKR3+RESqYvgTESlI+fA3CMGhHyJSjvLhHxFi5Fm+RKQc5cM/LMSImnqGPxGpRfnwjzAz/IlIPcqHfziHfYhIQQz/ECOq67m4GxGpxRTMk4UQSQAWAyjVmlKllHOFELO1tkQpZYr22IDaOlq4mT1/IlJPsD3/GCnlECnlGAC3A1ioBTqklKkAIISYHmhbkLX8LBEhRtQy/IlIMUGFvyu4NYlSymwAYwFka23ZAJJa0dbhws1GVPOALxEppk3G/IUQyVLKJdpNi8/dsa1oa7RdIUSGECKjoKCgDSptLDzExNk+RKScFsf8teGZGJ/mbJ9e/wwArjH7Uj+PD7TNi3YcIAUArFZru5yGG242oqiqrj02TUTUabUY/h49er+EEBafpnSc6dUnAlir3Q6krcNFhBhxooQ9fyJSS1sM+8QAKHbd0D4sErUDuBYpZWqgbW1QS6uFccyfiBQkpOwai5pZrVaZkZHR5ts9XlyNP3+6HV/dM6nNt01EpDchRKaU0urbrvxJXgNjIhBiUv5tICLFMPUAVNXZONefiJTC8AewN68c936SpXcZREQdhuGvSf0xX+8SiIg6DMMfwLcPTcPQPlF6l0FE1GEY/gAG9ApHbGSI3mUQEXUYhj8AIYTeJRARdSiGv8YgBOpsnPFDRGpg+Gv6W8KRX841fohIDQx/TVxUCAoqGf5EpAaGvyYuKhSFFQx/IlIDw18T1yMERVX1epdBRNQhGP6auKhQ7Mkt07sMIqIOwfDXJMX3wsfbjqGrrHJKRBQMhr8mMtR5XZusYyU6V0JE1P4Y/h4+uX083tlwGHYHe/9E1L0x/D1c0K8nUn/Mx0dbj+pdChFRu2L4e4gOdw79VNXbdK6EiKh9Mfw9uNb4qapj+BNR98bw9/HRn8Zj1wlO+SSi7o3h72PS0FgUVfJkLyLq3hj+PoQQiAozwcEZP0TUjTH8/YiL4lIPRNS9Mfz96NMjDLtOlOpdBhFRu2H4+zH9vL5YujNP7zKIiNoNw9+P8YkxOM0LuxBRN8bw98NsNHCJByLq1hj+TTCbeE1fIuq+GP5NGNgrAseKqvUug4ioXTD8mzBxaBxmvLoJN/9zq96lEBG1OYZ/Ey4aaAEAbD5UpG8hRETtgOHfhIExEe6fRz65WsdKiIjaHsO/GevnTAUAVNTZeHlHIupWGP7NGBwXiSV3TnD+/MgKnashImo7DP8WWBNiMC4hRu8yiIjaFMM/AM9dMwIAUFPPef9E1D0w/ANwTt8euME6ACfLavQuhYioTTD8AxQZakI1e/5E1E0w/AMkIJCRU6x3GUREbYLhH6B+PcPw1NJ9epdBRNQmGP4Bum3yYIwbzFk/RNQ9MPwDJITgdX2JqNtg+LdC3+gwzvghom6B4d8KIwf0xJ7ccr3LICIKGsO/FcYM6oXbP8jgOj9E1OUx/FshKb4XAKC81qZzJUREwWH4t4LRIAAAt3+QoXMlRETBYfi30q8u7I9QE982IuramGKt9Nis8/HtwUKs3Xda71KIiH42hn8r9e4RCsA59LNg1X7kFFbpXBERUesFFf5CiCQhxGEhRKb27yWtvUQIsVYI8ZDHY2cLIaYLIZKba+sKFlw3CgDw9obD2HG8VN9iiIh+hmB7/jFSyiFSyjEAbgewUGu/Xko5Q0q5AHCGPABIKVO129P9tQVZS4e5YexA989lNQ06VkJE9PMEFf6u4NYkSimztZ8tQohEj/vGAnDdlw0gqYm2LuPbh6bhX3+0oqiqXu9SiIharU3G/IUQyVLKJR5NMQCKhRCubwIWn6fENtHWaLtCiAwhREZBQUFblNpmBsZEYECvCBRV1uldChFRq5laeoA2POO7nGW2T69/BoAU1w0pZYr23FLt+aV+tuGvzYu2nRQAsFqtne602v6WMOSVNr3Wz7GiamQdK8GvLzq7A6siImpZi+Hv06NvRAhh8bmdDCBDSpnl0ZyOMz39RABrtdu+bV1KjzAzqurseOyr3bCEh2Ds4BhMOae3+/4Ve05i/sr9DH8i6nRaDP8AxADwvMTVIgCJHgd0lwCAEOIh7aCuxeMgb6O2ruaHnGL84LrC13og+4WZMBgEpJSYv3I/AKCkqh69IkN0rJKIyFvQY/5Symwp5R0et0ullFlSyiVSyrke7QuklKmuGUBNtXU1M0ee5f65Z7jZPfunos65/s99lw3Da+sO6lIbEVFTeJJXkG4eP8j982/Hx2Pu57uQ8PByjHpqDWaPGYAbxw7E+1tyUMgDw0TUiTD8g3RB/2hMGhqLjQ9OhcMhscZj2YeCijqcbQnHnVOG4GhRtY5VEhF5Y/gHyRIRgo9vuxiDYiPdB3tvHh8PAGiwOwAAYWYDKmp5MhgRdR4M/zY0cWgccubPwl1Th+Cs6DC8cv1oAMDx4hr873vp/AAgok6D4d8OBvSKwNZ5l6G/JRwAkF9RCwAY+dQa97cBIiI9Mfw7wNwrh7t/5lpARNQZMPw7wIize+K8ftEAnOHfYHdg6svrsTevTOfKiEhVbXGSFwVg5f9NxjsbDuOyv21ERIgR1fV2nCipwQX9eyKvtAbV9XYM7ROld5lEpAj2/DtQiHb5x+p6OwDnmb9SSjy3fB/+8K9tWLYrT8/yiEgh7Pl3oD9OGISSqnq8uf4QJg2NxcNf7EZRVT1W7D4FALj3k+345aj+OldJRCpgz78DmYwGzLniXOTMn4WXtKuBvbz6gM5VEZGKGP46GdArAvddNgwAsPfpK9ztDof0WgqirKYB9TZODyWitiWk7HTL5PtltVplRkaG3mW0KSklpAQMBoHc0hpMfikNl53XF2u1JSI+uW087v9sB6af3xcvXDNS52qJqCsSQmRKKa2+7ez560gIAYNBAADOtoTDOijGHfwAkPxhJvIr6vDJtmP44Pscnaokou6I4d+JLLpzAq68wLlE9JNXn49KbVloAEjbn69XWUTUDXG2Tyfz7u/HIL+8FtHhZjy9dB/O6RuF2MhQVHl8EBARBYvh3wn1iQ4DAOx88nL0DDcDAK567VvU1NsRHmIEANjsDgghYNSGjYiIWoPDPp2YK/gBYEJiLLYdKXLffvzrvRgybwWyCyr1KI2IujiGfxcxNqEX/ve9dNTU25FTWIVPfzgGAPjhSHELzyQiaozDPl3EVSP74e6pQ/DA4h1nzgieNhTHS3iFMCJqPfb8u5ARZ/d0B/+ep6/ADdaB2HK4CLUNdp0rI6KuhuHfhZyvLQud/cJMRIWaEBMVgu3HSrHppwL3Y2rq7diWXdTUJoiIADD8u5T4mAjMGtnPfWJYVKhz1M4gzsz4yTxaghtTtuK7g4W61EhEXQPDvwsxGATeujnJq23B7FEo9bg6mOvEsN/9a1uH1kZEXQvDv4vrGW52Xxryu4OFuPOjTPd9e3J5pTAi8o+zfbo4S7gZd3yYiRCjwONf7wVw5gNh14kyhJoMWLrrJEb0j8bl2tIRRETs+XdxPcKcJ4K5gv/eaUOx88nL8dpvLsSBU+WY8eomvL7uID7adqxNX/cXf9uAQ/md4wSzoso6lFbX610GUZfC8O/iYqNCEBVqwmOzzsORF2dizhXnAgCG9I7C0l0nAQB3TEnEpp8KcN07W/DQkp1IeHg5Dp6uwMUvrEOD3ftaAUcKqzD7nS348PsczF+53+9rOhwS2QVVmP73jXA49F0SfE9uGcY8l4pl2u9KRIFh+HdxfaPDsOfpK3Db5EQIj1k/5/TtgeIqZ294yjm9kRAbgcyjJVi1x3mewIxXN+FUeS2GPboSdTY7Eh5ejnU/nsa0VzYg42gJHv96L97deBj3/3d7o9c87LGkRHZhVYs1XvfOFry78TCOFTV/QtqqPadw23/SA/q9AaDB7sDzy38EABRU1LXw6M6tq1xXg7oPhn83FWIyYPPDv8BN4wbigv49sfy+yQCA8lob/me08zrBCbERAIAPthwFAPzpP86L5Ywa0NO9ne99zhlYlH4cJ0prAADTzu2NJ7/Z02ItmUdLMH/lfry1/hAA56J0/hzKr0Dqj/kBBeGBUxV4b/MRd30FlY3Dv6nX6WwOnq7A4EdW6F1Gq9gdEnadv/VRcBj+3djZlnC8eO0o9Aw3IzLUhE9uHw8AeP2miwAAGx6chnUPTMGC1fsRH+P8IEh7YApusA4EANw9dQiS4nvheHE1KutsOFVWi4c+34Vb3kvH/dOH4ZlfjcDmQ0Uoq2mA3SHx9NK97uC22R0or21ATb3z7OOoUBPyymqwas8pDH10Jf747x+8at2TW4bPMo4DAHJLa1BaXd/oQ+C+T7cjZdNhAMAV/9iEF1Y4h6UuirfgeHE1aurtqLPZ8dQ3e5Hw8HIMfXRlm7+nbaG2wVmny4xXNwGA7kNorfHYV3tw7ydZepehu4rahi77rY2zfRQycUgcDj5/FQBg/7NXAgAS4yLRYJdYMHsUcktqkNg7Cv0t4WiwO3DLpMH4eNtRfJZ+HG9qvXaXK0echYExEbhseB+MfnoNeoSaUFFng9lowITEWHy49ajXBWh2PDED17y9Bc+v2AcA2PhTAU6X12L8C+vw+V0T8f6WHBwvrsGDV5yLr7bn4pU1P2HezOG4NmkA4qJCcSi/At/szMM3O/PQu0coAGDdA1Nw2d824tPbL8Yf/v0DzntiVaPfubLO5j4Zzrd97ue78MZvLnKfNNecoso6LNyUjXkzz/N7/8HTFUjbn487pgzxe/+baQcxaWgc3kg7hNoGO8LMRiz8/RiYjQaMODsaJVUNeHnNAcy9cniLtTz+1R5cPbo/hvfrgegwc4uPb2tSSvfCgm0lv7zWvZR5oP6zJQdDekfhkmFxbVqLS1l1A3pGNP/+jnxqDf71RysuO6+vu23n8VIsyTyBZ3894me/dl5pDd5IO4gXrx31s7fREvb8FWM2Ov/kYWbndQGEEDjy4kxcnBiL68YMcN93y6TBAICrR/dvFPxv35yE4Wc5l5qYe5UzrCrqbDjbEo6UTdm45f30RlceMxkNsCb0wvHiGmx/fAZuHh+PW993ju9f984WLN2Zh/VzpuIPEwbhlTU/AQBeWLEf1udScf9/t2Nr9pnVS5fuPImv7pmEIb2jcPiFmQgzG93HN1zWz5mKP10yGEt35mHF7pO47p0t2HyoEHtyy9Bgd2DS/DQs33USy3Y7DxSfLKtp9n3LKapCyqZslNU0IL+ittGB8hmvbsKLPgfICyvr3L38V9b8hGve3oK0/fnYcrgIafvz8fTSvThdXovoMDNyS2vwzobDePLrPbjgiVVYknkCq/eeQo6fYyofbj2KGxZ+j8kvrcf1725ptu728AePb20VtQ1e950sq2nVcNCe3DLMWbwT415Y59VeVFnn/tbYlCe/2Yvf/Wub10wv1/s1+50tkFLi6je++1nHgwoq6jD6mTWN2u0OiYSHl+OV1Qfcba7hUsD5jfdXb23Gh1uPNrv9I80cKyupqsfE+Wn49Ifjra67NRj+5HWg2Fd0mBn/uXUcMh6bjvunD0OI0YCZI/u57z+nbw9sfeQyAMCX90xEXFQIAODzuybi87sm4JPbx2PHEzMAAA9ecS6W33cJekWG4PlrRmJvXrl7O1sfuQyD4yLRI8yMG60DMW/mcPSKMEMI4KsdeXjsqz344u6JePiq4Ujbn4+YCOfruC5m4zntdNs857aiw8x45IvduPvjLGQeLcHN/9yGX77xHT7aehRlNQ24eXw80o8U48CpCkx4Mc1v0OaX1yLh4eV4Zuk+hBgNGP30Gox7fh2GPbqy0QcAAKz78TQOF1Si3uaA9blUPP7VHhw8faa2nPmz3D8fLapGRk4Jxg2OwaxRzvd0SeYJVNXbMWfxTtzxYSamvrIB1769GQs3Ooe7nl++z/38spoGpOeU4J5PspDw8HI4HBLXvL25USD709z6T5MXpLmHMg7lV8Jmd+C6d7Zg2a48ZOQUY/+pCqT+dQoMAliw6oDXcye8mIbP0o9j008FWLH7JCbNT8Pp8tomXyu3tAZLMk+c+Z2qG3C8uBo3pmx1f0v0J7+8FiEmZ3zlldYi61gJAGDqKxsgpUTG0RLsO1mO3bllKPGZBiylxMgnV6PQz3Eilx3HSwEA27Xturg6Ca4O0Tl9oxARYsSavadQWWfDhgMF8HXt25uxbFceJs1PQ57reNkrG1Be24B9eeU4cKrC67ocuaVnOiKLM9rvA4DDPtSiKef0BgDceslg3Dh2YKP7+/QIRf+eYYiLDEXGYzOa3E5EiAkX9D9zMHnFfZMxOC7SfXUyl5dmO7/qJl86BMeKqnGqvBZvbziEpPheiI+JwKL04+hn8R4iSHtgCirrbLA5JPpqwwe3XJKAxN6R+POnzhlLowf0xM4TZXhr/WGMS4jB3KuGY9RTa9y9tKmvbMBjs87DuxuzkTZnCqLDzO4DyjtPlGFQbASOesxYGvboSrx03UjcODYeYxN6ITzE5O4F3nFpIgBgceYJLNbCbbI2PJH610sx/e+b8O3BQnx7sBDrHpiCIb2jsONYGvpGh+LT5Ivxqzc34+5pQ3Hfp9uRdawUWcdKcdP4eGw4UIAv754Ih3R+YwKA5do018R5zoPGX2Tlos5mx+8uHoSIEO//4uc8uhKZj0/HjSlb8cGt4xAVZsKos3vCZDSgtLoeY59PRYNdoqLOhugwM6b/fSO+umcSMo+WoEeYCRsOFGDKOb0xtE8Uvrn3Evzyje9w07h4nN8/2n3MYt6Xu71ec9/JcvSNDkN2QSX+8tkOLL5zIs55zHk85q3fnlmupKrOhjfSDuKf3x0BAAw/q0ejfaiwsg7Xv/s9hvWJwgBLOK4ccRZmvv4tAOBRbUhut3Zm+6zXvwMA3JSyFe/dMhZ1NgdG9O+J6nobKups+PFkOSYP642/fLYDs0b2w7ThfbDzRCmS4nu5PxiueXsLcubPwphn12LrvMtwyUvr3bV8vSMXAgLV9XYkf5gJf/afKkfWsVIM3HcauaU1OHCqAv16OvfPK1/dhLyyMx+MmY9Nx5jnUnGDdQAevOJcvLz6APbmleN6v1sOHsOfAhYdZvY7xmwwCGzRev+tcX7/6BYfEx8bgfjYCIwbPA4AEBcVirQ5Uxs9LrF3lN96rx7dH1drs5sAID2nGNe/+z3umjoE0WFm/GbsQPw3/ThiIkNQWWvDc9rU0VFPnfnKHxsZgjduugjf7MyD2WjA0nsvwcinVmP6eX3x8bZj7oX1Prh1HF5LPYhXU3/Cwk3ZGNYnCge1bySr778U52phNrRPD0weFocHLj8Xb68/hCFa7bdPHoz42AiEmoxYdf+l7lB/87cX4d5PtuPZpfvgkBIXxfcC4DxuM3F+GvpGh6Gosg752vDGk984T/izOSQuGRqH7w8XIflS51TgervDfaDcNXyz+v5Lse1IEfr0CEWD3RngGw8UYMb5znHsX7+1GQCw4UAB4qJCEar1uM/SQmzm699i91OXY+RT3sMkS+6cgB3HS3HLe+nY9OA0zPtyN3aeKHMHPwDc80kW4mMiMKR3JF5fdxBm05nBCJtWy6KM47h+zAB8uT0Xf120E4Bz2OTbh6Zh9d5T7sc/v8L5t/tH6kF32w3WAViUcQJvrT+E1XtPAwDCtSHP48XOHvaX23Px5fZcrJ8zFde+vQVf3TMJX2blureRX16Loqp6DNMmENxxaSIWbsrG//13BwDnN83xHsNWr990Ed7ZcBhZx0rw8qoDuHXSYPx7s/MD7ZudebigfzSG9I5Ez3AzahrsKKl2flMb81yq9vuewFf3TMKFAy1IzymGze6Aydj2gzSiqxyptlqtMiMjo+UHErXgvMdX4elfXeCe1ZTw8HJMSIzFp8kX47ll+7D5cBF+PHlmSMo1VFNdb4NDwusA8ux3tiDjaAneuOki94dMbmkNJs1Pw/eP/AJhJmfQWCLMzQ6v+WN3SJRU1yMuKhRf78h1h43n0NEXWSfQp0cYrAm9MO+L3fjl6H649X3n/5PkSxOxKOM4SrVwuTbpbHyhhdqcy89BYWU93t+S497WhQMt+PMvhmLlnlNeQzGAc+rw1/dMwuC4SNgdEpGhJjgcEuNfXIeCijrMHHkWVuw+hUeuGo4Z5/dFVKgJfaLDsPtEGa5+87tGv5srQAHgtd9ciKtG9MPv/7UNYxNi8Ob6Q5g0NBYOB/Bp8sVIeHi5+3mjBvTErhPOnv2Pz1yJjT/lY87iXfj8rom44h+bcG7fHjhwugK/HNUPN1gH4lB+JZ5Z1nj4aFxCDEYP7Imbxw/CPZ9kYW9eOf5ndH98szPP/ZhV90/G3R9nIbvAezgwZ/4srNx9Eo98uRtDe0dhyV0TcaSwCtNe2YDIECP2PnMlblj4PU6X1+JoUTV2PnE5Rj+zBh/9abx7wcVnf3UBfj8hAQBQXFWPy1/dhMLKOiRfmoiUTdnImT8LP54sx1WvfYvePUKR/uj0JvaSlgkhMqWU1kbtDH9STVWdDWFmo/t4QW2DHaEmg1c419Tb8f++zcZ3hwqx6I4JTW7rs/RjmPv5bq9AllJiwerAZu4EqriqHl/vyMVN4+LdB+ubsvlQIVbsPomPtx3DoNgIjEuIcQ89ufww7zKUVDfgv+nH8N7mHCTGRSK7sAo582ehut6G859YjcvP74teESH4LOM4hvaJQupfp/h9vb98tgMbDuQj6/EZfj/gFqzaj+MlNdibW4b3bxmHKa+sx8e3jUeS9g3G9d67Qn7k2T2x+M4JmPnat8ivqHOvVAsA6Y9Ox9jnU7HxwakYFBsJKc8MUdU22LHrRBluWPg9jrw4E0IInC6vxcurDzT6MNv04DRc+rJzCGf84Bhcbx2IOYt34tlfj0BFbQOuuehs9OsZjtoGO/6zJQcrdp/EzhNlePd3SbhyhPP4TL3NAYeU7r+H536U9Oxa9ySEnPmzUFJVj16RIfgi6wRW7D6FlN+P8ZpltmrPScRGheKigRbUNNjdy7a43hPP/au1mgp/SCm7xL8xY8ZIoo7kcDikze5o8XENNnsHVNM6DodD/vu7bFlSVSellPLdDYdkysbD8uDpCrliV57XYxtsdvnTqXL59Y5cd9tb6w/K7IJKeSi/Qg6au0zWN/M7LtuZJye+uC7g2pp6T99ef0gOmrtMbj9WIqWUctDcZXLQ3GXywcU73D839/ymtl9eUy8/2pojG2x2WVNvk//+Lls6HA45aO4yOe+LXXJfXpmUUsrpf9sgV+852Wh7rv2griHwv/N3BwvkA4t2yA+2HAn4Of58sOWIPF1WE9Q2AGRIP5nKnj8RNUtK2eKQVSCPaYndIbHppwJMG94HgHMZkbKaBiTF93LPPgr2NZrTFr9DZ9RUz58HfImoWYEEYluEptEg3MEPwH0gvK2235LuGPzN4Tx/IiIFMfyJiBTE8CciUhDDn4hIQQx/IiIFMfyJiBTE8CciUhDDn4hIQQx/IiIFBX2GrxBiNoBSAIlSypRg24iIqP0F1fMXQkwHkC2lTAWQLYRI0gIdWhuEENMDbQumFiIiClywwz4ZABYLIZLg7L1nARgLIFu7PxtAUivaiIioAwQV/lLKUgALASwGMEZrtvg8LLYVbV6EEMlCiAwhREZBQeNrYxIR0c/T4pi/NjwT49OcLaVM1YZqUqWUC4QQL3mM4fs+PtA2L9pxANfxgQIhxNGW6m1CHIDCn/nc9sS6Wod1tQ7rar3OWlswdQ3y19hi+EsplzRzd5KUcoH284sAbgCQjjO9+kQAa7XbgbQ1V0fvlmptihAiw9961npjXa3DulqHdbVeZ62tPeoKdrZPihAiGc4xe89ZPA9p3wosHgd0A2ojIqL2F1T4a2P+jaZoenwbSG1tGxERtT9VTvLqrOcQsK7WYV2tw7par7PW1uZ1dZlr+BIRUdtRpedPREQeun34CyFma2cUJ+v0+hbXmc9CiJc82kuEEGuFEA/pVWugNXRkXdp7dVgIkan9e6k1tbZDPbOFEGv9tLX4HrVnfb51dab9rIn3TPd9zc97pvu+1szfrd33sW4d/p1kCYkbAFhdU2Y9/kjXSylnuA5661RrizXoUFeMlHKIlHIMgNvhPIkwoFrboxjfqc6BvkftXZ+fKdidZj9rYnq47vuan7o6w77W6O/WUftYtw5/dIIlJKSUKR6L1iV61GMRQiR6PFSPWgOpoUPr8pnymyil7EzvV1Ov2xnet868nwVah3L7WhN/tw7Zx7p7+Ft8bjdaQqKjaDtTsccOFwOgWAjh6m1YfJ7SEbUGUoO/tnYnhEj26al1hverqdcNtK3dddL9LNA6/LW1u86wr/n83fy9XqBtAevu4V+KFpaQ6ECzpZR3uG5on/ilAEpF08titKsAa+jwujQzPG90hvdL4+91A23rCJ1uP2tFHf7aOkJn2Nc8/27+Xi/QtoAFvZ5/J+dvqYkOJ4SY7TF+mATACiBDWwXVpUNr1caEA6nB0pF1abVZfG4HWmtHaJflS9pCZ9zPtFq4rzVfg+/frUP2sW7d89e+yiUKHZeQ0F77JdeMAjg/qRdp97kO2CzRodaAatDpPYwBUNzaWtujEG37Vs/X9n1dPd4337o6037mW1ugdXT0e6bRdV/z93frqH2MJ3kRESmoW/f8iYjIP4Y/EZGCGP5ERApi+BMRKYjhT0SkIIY/EZGCGP5ERApi+BMRKej/A+rhduTxt6WdAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(Es)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "8bb82a79",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAFkCAYAAAA5XmCyAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAAsTAAALEwEAmpwYAABCkklEQVR4nO3dd5gb1bk/8O+oS9vXXrd1XRds09xNrzaEkALBkN5IYtLzu7k3hCQX0uAGuCGVewOEFAIpEOCGEFowYEowrmCMsY3tdW+79nbtqszM+f0xZUfSSKvdVZvV9/M8PEgzs9JZrfzq6D3nvEcSQoCIiJzHVewGEBHR0DCAExE5FAM4EZFDMYATETkUAzgRkUMxgBMROZSnkE82evRoMXXq1EI+JRGR423cuPG4EKIh+XhBA/jUqVOxYcOGQj4lEZHjSZK0z+44UyhERA7FAE5E5FAM4EREDsUATkTkUAzgREQOxQBORORQDOBERA7FAE5E5FAM4EREDsUATkTkUAzgREQOxQBORORQDOBERDny+fs3orU7WrDnYwAnIsqRp7ceRUt3pGDPxwBORORQDOBERDkUiSsFey4GcCKiHHp00yH0xQoTxBnAiYhy6I9r9+PSn71UkOdiACciyoGo3N/r3t/Wi62HO/GuPAdyBnAiohxoD8cT7h9o68X2o915fU4GcCKiHIgrKqoD/fvEf/6BTQCAxzcfzttzMoATEeWAogpUBbwpx7/y59fz9pwM4EREOfD9x7eiytIDN4R87rw9JwM4EVEOvLCjFVNHVaQc740pONLZl5fnZAAnIsqRoM+NPT96t3n/d59ejAWTa/M2L5wBnIgoR1bvaIEkSbjpPXMBANUBD06eUIO4IvLyfAzgRETDJCsqAMAlSQCAa8+ZBgAIeN3wul2I6+dzjQGciGiYIrIWoAPe/gHLx798DmaOqYLXI+UtgKcOmRIR0aAYOe7emGweO3ViDQDA63IxhUJEVKqMCoR288C9bpeZYsk19sCJiIYpHJPx4SWT8J3L56ac83okxJgDJyIqTduPdGPmmCpU+lP7xEyhEBGVsGe2HsWiqXW257xuKW8pFAZwIqJhEgKYWBeyPedxu5hCISIqVRFZQcBrH059bhe+9pc38vK8DOBERMPUG1UQ8NgXrXK7pLw9LwM4EdEwrdvbBleaQN3ZF7c9ngsM4EREw7R0Wn3ac+fNasDdH1+Yl+flPHAiomGQFTVjmmTGmErMGFOZl+dmD5yIaADr97ZBVe3nckdkNaEGSiExgBMRDeDqu9bgWHfE9lx3JI4KmwU8hTBgAJckaYUkScskSVo5wHW35a5ZRESlJS4LdPbGEY72F6zqjcn442v7sXBybVHalDGAS5K0AgCEEKv0+8vSXLcMQFPOW0dEVCJ6ojI+8bt1uOmxreax257ajjtf2IXxtcGitGmgHvhiAM367WYAC5IvkCSpyXINEZFjGVUF7bT3xrD5QAc6+2Jmfe+TG7WSsfncuDiTgQJ4bdL9UTbXNAkhGMCJyPFm3/g0wlEZQqQOWH703rUAgFXbWvDVP78OALj+4TcBANF4fpbKD2SgAN4BIO0ER0mSlhnplQzXrJQkaYMkSRtaW1uH0EQiosI5+bvP4JmtxzJes7+t17z9kaWTcdHsMflulq2BAvh69PfCmwA8m3S+TR/gXAGgSZKklBSLEOIeIcQiIcSihoaGYTeYiCjfrDvrAMDEuiAuPKk/fh2wBPDJ9aG0qzDzLWMAF0I8DC0wLwNQaxnMfFY/v0k/Vo/UdAsRkWP85pU95u2vP7Q54Vx1wIvemIJlc7SedsinTRtcPLUOK88t3vyNAScvCiFu12+ushxbnnTNPQDuyW3TiIgK5+9vHLI9HokrqAx4sLs1jFMaqwEAoyp9UFQBj8tVtN43wIU8REQZbT3chYYqP473RLF6hzaOV+H3ICar8KcpIVsoDOBERABi+rZnZzTVY4mlONW2I104e/po8/6W710CAOiLpy8hWygM4EREAKY3VOA3n1yEv6w8E9akyJaDnQnbpVUFvOjqi+MTv12Lp7ceLXxDLRjAiYgAqEJg3qTalONdkThGV/oBAHd9TJtod8FJY+B1u3DOjNEp1xcSAzgREYCYLOD1pIbEnqiMCr8bZzTVY854bRCzNuTF6/s78LEzJhe6mQlYD5yICEBUVuC3CeDRuAqf24W/rDzTPGZc52cOnIio+GKyah+QJUCSEqcK+vQA7rMJ+IXEAE5EZau5tQef+f16nOiJ4mhXf71vSdJ22nlowwHYzfJ+/7xGALDtsRcSAzgRla2dLT14bnsLFt68ChW+/ozy2OoAWrqjuP7hN2FT1wqVfg/mT67FlFEVBWxtKubAiahsHe+JmrdVS6SuCnjQHdHqoazb22b7s//3xbPz27gssAdORGXLCNKAlgM3VAW86OyLF6NJg8IATkRlK24J2tGEAO5Ba3cUF57UgL23Xl6MpmWFAZyIypaxs47P40rYxKHK70Frd6Rou81niwGciMqWUf+kaXRFwlTBqoAXrT1RBnAiolIlKyoe+9LZOHlCTcLxqoAHx7tjCBS52uBASrt1RER5FFfsS8JW+j1o7YkWfaXlQDiNkIjKVkxR4XW7cO7M0ZhQGzCP14S8eH57C86eYbePe+lgD5yIyspbhzrN//953QEEvW5cMb8R/37JSeY19RU+AMC/dp0oShuzxQBORGXlPb98BYC2IQMAhHypaZK6kBbAv3rRjMI1bAgYwImo7Dz79jFzCqHdTBOvWwuNX75oZkHbNVgM4ERUdh7ZeBBxRWDKqFDaglSbb7qk6NUGB1LarSMiyoOntx6FrKj46NLJKaViDTUhb4FbNXichUJEZWn93na4Hd6FdXjziYgGb/HUOtz14m488Nr+YjdlWBjAiaisNDVUYEy1Nuf75itOKXJrhocBnIjKis/twvo9Wo3vs6aX9kKdgTCAE1FZCfrcaOnWNnKoCpT+QGUmDOBEVFa8bheWzRkLoPibEg+Xs1tPRDQIQghI0DYtHl8TGPD6UscATkQjwo6j3VBVmx2ILXpjCgJeN4QQZr0TJ2MAJ6IR4dKfvYTDnX0Zrzn5u89AVlWoAnClWcDjJFzIQ0SOdrijD+v1nePbwjFMrAvZXtfRGwMA9EQVKKpwfP4bYA+ciBzu6beO4mt/eQMAcONjW9Ne9747/wUA2Hs8jO+972T8+OrTC9G8vGIAJyJHUy2bEW8+0AFZURPOf/GPGyGEQHtY64H3xRXMHleNaaMrCtrOfGAAJ6IRxZjjbXhyy1H0xRUYYT45wDsZc+BE5Gg3P7Et4X5nXxwTaoO44587ML4mCAAIRxX0RGUAwJ0fWVDwNuYLAzgROdbUG55IuF9f4TMD9eodreZUwbB+DADefer4wjUwz5hCIaIR49qzp+Lqu9bgoQ0HsOVQJ1z6TMFwTMayOWOxeGpdcRuYYwzgRDQi/OLD882pgdc//CYAwO/RtksLRxXEFBX3f2Zp0dqXDwzgRDQijK70wZe0Q4NXD+jdkTiicSXt9mlONbJ+GyIqW7Ii4PMkblDs1lMoe46HASDt9mlOxQBORI6kJNU9qQl6U3rYigCqAx7c/MQ2ZK6S4kychUJEjnS0K2Lefvr/nYvZ46ohJwX1zr44Gqr86IrION4TTX4Ix2MPnIgc6exbnzdvT6nXVlUunFKHDyxoNI8f6ehDhV/rpza3hgvbwAJgACciR7ti3gQEff2570+fNQ0AcON75mJnSw+icRXnzWooVvPyigGcStZn79tQ7CaQA3z78jkJ96sCWo/7inkTAABfvHA6vnXZbNxw2eyCty3fmAOnktQdiWPVtmPFbgaVuI8unYwxVYk760wdXYG9t16Ozr44AG0u+Jzx1ZgzvroYTcwr9sCpJB3r0gac/vHm4SK3hErVkmn1uOXKU9OeN+aE+70jN8yN3N+MHO3zD2wEAHz5T68XuSXkVF59EvhIW7xjNXJ/M3I0o3azsTT6RE8UO452F7NJVGB9MQVbD3cO+efdeiGUSv/IzRQzgFPJCUdldOvV467Sp4T9+uU9eO8vXylms6jA5tz0NC7/xdD/5saqy7qQ8zcvTocBnErOJ367DjFZK7rfHdEC+V0v7kZML8T/2fs24NXdx4vWPiqsnceG982rbgTsPp8OAziVnL6YAgC479olZgC3WrXtmO1xGlmWTqsHANsl8HFFNXPcA6nwuQe+yKEYwKnkGIsyJtYFEYkrttckV52j8tITkbPKbT/2pbNHXAErK/4roJIiKyo27msHAIQy9Jyisn1gp5HHSKdZdUdkVAW8A/7s6ZNq89Ci0jHgR5gkSSsAdABoEkLcY3N+mX5zuRDim7ltHpWbY5YNaYNeLYDHFRVnzxgFWRFmBbqozT9qGnm+cMF02w/rrkjcXHFZzjL2wPXgDSHEKv3+sqTzCwAs0M8vkCSpKV8NpfJwwlIxzkiltHRHzdV2sqoF7i7mwMuC3+NCNJ76Yd3ZF0d1Fj3wkW6gFMpiAM367WYACds5CyE2CSFulySpFkCzEKIZRMNgLflp5LnbwzHUhrR/rIoqUBvy4mBbb86fu6UrgneGOeOBckdAWwYfVVIDePPxMJoaKgrfqBIzUACvTbo/Ks11iwDstjshSdJKSZI2SJK0obW1dZDNo3Kz70R/YDYGn+KKCp/HBUnS8qGjKny4+6Vmc5eVXPnZczvxkV+vzelj0vD40vTAIzEFIR9TKAMF8A4A9QM9iJ5CmW6kXJLO3SOEWCSEWNTQMDJLOlLuHO3UivR/+939lePiioDP7ULA60ZPVEa9Pq9374ncBnAhBEbwhAXHkaCnUGxy4DH9Q73cDfQKrEd/L7wJwLPWk5Ik3SZJ0kr9bgeyCPZEmYRjMq5/10n43LnacIoQxpxfFwIeN8JRBbX6yjohcrtJlqwItHaPvF1bnEwL4Kk98JG4QfFQZHwFhBAPA2jSBy9rLYOZRiC/G0Cz5XzKLBWiweiNKXjvaRPM9Mm6vW346L1r4XW74Pe60BONmz0vNccTUbg4qHQYH84hnwe90dS/S5Q9cABZTCMUQtyu31xlObZc/38z+gc5V4FomPpiiu38b1UIBDxu/PaVvXhiyxEA/TNScmUk7pnoVHFFwOt2YUy1H28fSS1otXFvO65eOKkILSst/AijkhKJKwh4+wP4R5dOBgD0RGUEvC6MqwngPy+fg1s/cCp6Y7ldzDN1dAV7dSXCGLgeU+XHsa4ofvLPHZj1nafM826XhOmchcIdeai0ROJqQgB/aMMBAHoP3OvGoY4+LJlWj6isorM3ltPnDkdlTK4P5fQxaWhislbrZFSlHyd6onh440EAxkCzll4byUvks8XuBpUURQizjjMABDxaML/05HHwe1xo6YqiJuiFz+3CjY9tzelz90Rl7GrpwQOv7cvp49LgaT1wN3xuF+JK/2D1/662na1cthjAqaT9+pOLAAALJtfB43Zh6+FO1AS9mFgXzPlzGTU3jnVFcv7YNDg7W3oAaLvqxC0LefbmeO6/0zGAU0lJ/lJcE+xfLr3jWDfCMQW1IS9OaazBkmm5n7X60w+ezhobRTLj20+atz9671o8vvkwJElKKCcrqwLr9rRh7Z62wjewBDGAU0mb3lCJP1y7BIBWqRBIDOq5ngs+tjqAcLR/cPTiO1Zj6g1P4K1DQ9/ai7Ij64XKjIJl1nrfRm3wuKLif1fvKnzjShQDOJWMqKyk9Kx8HhfOm6Wt4PW4XAh4XWaVwuT8qNXGfe3oGMIg56S6EPZb6qzsbtW+sr/4DstAFMqa3ScAADt+eFnC8R++/2Q0t4axekcrrjufdfMABnAqIbtbMuc3/+sDp2LNDRebsw8CXhcOd/TZ9sKv+tWreOqto4NuQ0OVHyfCqYGfEx7y5/3/8y9M+9YTALRvWYr+93S5El/082Y14O0jXQAAPzf0AMAATiUk4M38dqwJehP2N/R73bjgx6vT5kMHejyrqKzA53Hp5Uu1FIr1g0FKyc5Trmw+0AHjpY7KauoHstAqE1rLx67gIh4ADOBUQhRV4GNnTM76eqMWxofuec32fMUgqtWFowoqfB5z0ExVBZb/9CXzvIvxuyAe2XQQsiJw9oz+wqchvxt9MSVhcHnyKM7XBxjAqYTIqoDHlf1bMuB1ozrgMasTGvbrJWmfTkqhPLXliDkl7Z9bjybstxmOygj5tdz6uj1t+OIfN2GfpdqhizmUgrjpsa3o7Itj2Zyx5rGTxlbB7ZLgYdokBV8RKhmykriIZyABjxsetwtTRoUSvnYf7uwDADz6+iE88eYR8/gX/rjJrDe+8v6NONIZwcZ9bXjzYAd6oomb5B7q6MNHl07BXR9bCACo5NTCvFk6rR4Lp9SZ4wwH2nsTgvVJ46rMDT0oEQM4lQxZVeFxZx/Aa4JetIVjqA540WfpTVsD8Zf+tAlCCNz5/E4AwJZDHea5C3+8Gve+vAfXP/wmemMyKiw/F5NV7D0RxsyxlfjB+0+23ViXcue+a5fg1g+cCgA40NYHr+WDfO6EavNbFmvVJOKrQSVDUQU8g+iBG1tqVQU86LGUHFWFSKhp0hWR8eN/vgNAS49YBfVNInqiSkLgb++NYfWOVgS9bpzaWIOD7bnfwo36Vfo9GFXhB6BVhbT2wGePq8YtV2jBfeqo0KC+pY10DOBUMuKKgHsQOfDLThkHQMuFd/TGEx7HWuy/tbt/afyf1x1IfBBJG7AMR2VUWMrYdvbFsWRaPSbUBhHyeRCx2daLhi8SV8xetZGmOtDWi9pgYsrE2OB6Yl0I9+rlFYgBnEpIRFYSvjoPxOilNVT58fbhLvO4rKjmP/gFk2tx14uJe223WAJ6OKqlTnqi/SmUV755IcbVBMzeerptvWj4jvdE0VCp9byNb0DNx8MYWx2wvb464OEccAu+ElQy9rSG0dRQOaif2XTjckwbVQFZFfjRk9ugqgKyKnBqYw3ObBqFqxZONEuRGpbc8px5e+O+dlT4PQhbBjEn1oXQ0tW/uUPA67bd1ouGz9i4AUBCGeGgzz40VQe98A9ifv9Ix1eCSkZUVm1348mkvsIHn8eFmKzi7pea0RtX8NahTjRU+fHnlWdg8dTMBa+O98Tg87jMnrhBscxq8afZGX2kkZXC/46yosLr0b51Wed5+z3274MPLZ6MWWOrCtI2J2AAp5IRk9UhbVRrTXGEozJ+9NR2bNjbDgBZTT9zSdp+mMmzUMzH9478FEokrmCGZcebQokpqjn3f2x1AHt+9G4ASPs+mDuhGlUBTik0MIBTSejojSVsWDwYRg8c6N+Y2OjNjanqz6U21gbNQlhWAa8bx3tiCbNQEh7fbb8z+kgy+8ani/K8siIS/uZGnRtOF8wOXyUqCRff8SJ+/fKetF+dM/F5XOjs02ahfOTX2rL65GC8eGodrpg/IWG++MvXXwgACPncON4TRYXf/rk9bhcUVeDFd1rx8s6RXZUw1+V5BxJX1JSpo6/fuBy1IV+anyArBnAqCUYFwHRBNBOfuz+At3Rrg4/JC4Lefep4fOPS2QnHJtWH8JeVZyDgdaO1O5oQ9Lf/8F04d+bohOs/9bt1+PXLe2zbMPWGJ3BiBOxqn648bz6fz5s0q6SugsE7WwzgVFIah7BVms/jwh/X7k86mhjAjV7e6RNrEo6f0TTK0gPvD+ABrxv3f2ZpwrVCAMEMMyCcmmZ5+q3+cgOymvl36I7EM54frIiscFbJMPCVc6DFt6wqdhNyytiBBdB604Nll3ZZOKUu4X5vTE+dSBLmTao1l20DgFuS0NIdTekJWrXrm0P4MqR4nBrA1+sDvgAQlzP3wE/93j9z+tzdEZmDksPACj0O1Nrt/K/qVuFY/zJ4aQhV/6wDXlcvnIi/vXEIKxZONI/9+XNnYOpobWn97VedBr/HhamjK8zz960ZeBd6Y5C03WazByNv7NSZKuMsi2a036FwAbWrL45qFgobMvbAHWYkFlVq7Y4OayNhY8rZVy6agdMn1abkcc+cPgrja7TUzEnjqhKCd7YCXjf+e8VpiNnMlTaOOXW5vSIETp9UCwD407rkVFSq4c4Xj8QVbD2s7THaFYmjOsge+FAxgDtMr95blRUVn71vA452Rgb4idL1jzcP4/X97bj5H29jxpjBrcC0MnrgC6fUJfQms/W1i2cOeI3f48KMMZUpGzvc8sTbWL9HS0FE487sgfdGZbz3tPEAgMMdfbbXbD/ahW/8dTMAIDLMTsTqHa24/BevAAC6+uSEnXZocPjdxWHCei73Z6t2YtW2Y7j27KkYVzP4oFUKvvyn183bF5zUMOTHMfLm1UEvhjILbuV5TQk73dsJ+ty2y/x//fIes27HcANbIZx7+/N4+fqLEo6FYwqWzx2H331qMXYc67b9uUPtffirXpIgElfSzpnPhvVDsCcaH9ZjlTv2wB2mVy+bulffLcYcnHO4wezEk8zogU+uDw0pGFT4Pbj2nGkZr7nrYwvTBnnVyIGXeA9cCIEDbX3mrkS7W3vw2fvW4zev7EHI79ZWnKZJA1nrlPTFFISjMv77me1DaocxZfTlna2IxAdfPoH6MYA7jNEDX79Xq5SnFHjhRa4kDwaqw/g9jBy4xyVhQm0A7zt9wrDaZsdYWCKEFgiFEJh6g7aTek9U+5uUeg/cWMRkzKi5+I4XsWpbCwBt/1C/x512INb694nKCo509uF/Xtg9pHZ869EtALTNjHtjcsKHAw0OA3iJ2rS/HQ+u358S6Iwe+DG9Wt7mAx3Yczyc8vOl7kRYa//c8dXwuV3D6oUZZWU9bheqAl784sPzc9JGOz6PCzFFxdGu/rGHXzy3E1V+T8n3wH/9krYIya6XHfK79Zoy9h9CxuD5lfMb0RdTkYu6Vz/+5zt4cstRBDgPfMiYfCpRH7x7DeKKgCRJuGbRJPP4p36/PuG6/129G1sPd+G+a5cUuonDYswUeftIF9Z86yK0hWP4h2X/yqEYzG4+QxXU66acfevzCce7o3LJ98CNXrTdTJqQ141AhqJdRgAfVxNAb0zG0S7Z9rpsLJ1Wj7WWnZHYAx86fvSVKCPAyUlT4uymER7rKs2ZKEaKwY7xe73v9AkYXxM0t9MajkJstRXyufHVP79ue67Ue+DVQS/G1wQQk1XsaulJOOdxu+D3uPHAa/vx1w0HUn7WCPr1IR+6IjLeOtQ5pDYoqoAraa5/pgVUlBlfuRKnJC1tNqrpjasOmD1O60IYp4jrv9d5s7TZJ7moPleIHnjI78HGfe2254ayEnPjvvaElaj51NIVwZhqLYAv+8mLKeeNsQTr9nQG43erCXnR0RtDOCoPKe0ViStMmeQQX8kSZUyNS/7Hfflp43Hl/EYsnFJnFmxK7tGUgo//Zi2A/nnryYweeEOV1vMezG706QxlFedghdJ83f/mu2YPqQd+1a9exZ7jPQNfmANrmk9g5phKcxYKgIRgapQksKtNEpNV/OLD81Hh8+DB9Qcgq2JI33i0AO7Gqq+fb1aDpKFjAC9B+06Eza+sspqaQvnmu2bj/JMazJV/Fb7SG8p4eedxAMDcm57BtiNdKedlRcXXLp6J8/UeeHXAi9e+dXFB2zgU6XqdZ04fZebAH3vjEN535yspaYp0CjWRKORzo7E2mJCG++jSKebtCr8b581qQE809UM3JqvwuV3oicaxYV874oqKoNed8GGQjYisIuB1Y8aYSlYdzAEG8BJkHcxL7oHHFRVet4Tp+qKSRVPqUB0svQBuZVc/RFYFvEm9bicsSApaPizv/Mh83HaVVhQr6HUjovfA3z7chTcPduKfbx/N+FhGgC/UVFAh9P09LUH3xvfMxebvXgJAy4N//rwm23GWmKLC55Hw/nmNmDepFrIiUBvyImwT7DOxplAq/R5s+8G7hvEbEQN4CbIW1U/ugUdlFX6v2wx+f/38mShQCnXIoja9tE/8dl1OB68e/vyZOXusTKw98HpLDzLodaNPn6NvpBYGyskbNcy7+go3hmHsXhTyufFVvYSAdYGSdXcjK60Hrk019HtckFWB2qDP3AEJAN77y1dSNpBOdqCt1/y9AW2FKw0dA3gJsq4mTO6B98ZkhLxu1IV8GF3pL0jed7CSd3XpS7NaNJdFjBYNsHlxrlgDjsflwhXzG/Hy9Rci4HOZKRQjcFsDlR1j0PCau9fkqbWJBPoD9CmNNfj68lkp1/g8LrSFY2hu7UkYrI3JKnweFyRJgoBWN3xcTQDHLZtYbD3cOWA+/75X9+LJLZm/mVD2GMBLUF2FDzdfcQoArQe+ekcLfvevPVjbfAJCAC6XhEn1IWz4z2VFbqm9I50RnDtzNL5y0QwAqcv9jbnGdQ7cNsvogZ8+sQbjawLwe9yYVB9CwNIDjykCt37g1AEDuPVzbleLfQ2SXPO7XdjZ0oN1lnnYVj6PC39ZfwDLf/oSPveHDQC0AldaCqU/XMiKwJgqv7mqEzDSSIm99+NJuxRdNGdsXhdalRsG8BJk9HaM29/9+1Z8//G38cF7XktYAGFnV0vPsMt9DtU7eiGks259Hm8c6MA5M0ZjXHUAfTEZXZE4bn1Kq50RjhoB3HlV6IwB48e+fA4m1YfM4359hSYAhKMyJtQG8cBrmUuzWhfUnOhJHSfIJSEEdrX0wOdxobU7/boB6+ynmKxi6g1P4KxbnzcHMQ1PbDmCyoAnIWAHLOMAxnMuulnbfGTnsW6oqoCsqPAWYLpnuWAAL0ExRYXf48LoSh8kCdh3ojfj9RJgzgZY9pMX8Vpz5iA/WIoqspptcMlPXzJv14V8WNo0CjdfcQp6YwrawzHc9eJunOiJYvOBDgDAqErn9cDT5e19bhdisoJIXDED+EAGO4NjOHpjCtrCMdQEvWjtjmLleU2211nbbZ2NErV0KoxvGmOqAgkBO+B1J2wabdTt2X+iF8t/+hIOdfRBVoRZ+oCGj69kCYrJWgC/5xOLUjbKveXKU1KuX7unDadZtrrKdWD42ap38LF712Z17YmeKE5prMaD150BQEs59MUVvLBdK5p0+9M78Gm9HMC00UOvAV4sMcU+ny9JEoTQCkQ9+vohjNIHONPl/wEtDWGMd+S7qmRvTMHscVWor/DhQFtf2sqK6Za1d0XiZs5+yqgQGmuDCPkSUyZ+rwuPbjqEDj2tYsw+OtCudUAUVSCuqjmZ808aBvASZKRQfG4XHtqQOKo/d3y17c9Yez6f/v36QU/vyuRoZ2TAbwGGvSfCCPk85g44QZ8bPREZ33v8bQBAr97Or140oyBL33OtsTaEc2aMtj23/Wg3DukbIoT8WiBcef+GtI8VV1R84YLp+MalJ+V9Ne3iW1ZBkiQtgLf3oiLD7I8bLpudcuyJN4+YPfAxVQF09MZQG/Li2/+3xbzGWCU87wfPAuj/5mi8F1/Y0QJZEfAOo3QwJeIrWYKMKVt2S44zFf6x1h4ZbkA40tmHXz63EwAwmIku3RE5YT/4kM+DJ7b0z2t/fPNhAMCL+kIfpzlpXBUe+OxS23PWQUtjVWOmOjUxRVvUMr4mgN5o/nrgRgDddqQLoyp96I0pCGWom36dnl65MGmTDSMHXhvyIhxTME3fmu6B17Q9RZOnTb59pBMBr8t8L37/8bfx6KaD7IHnEAN4CTLyjXZL5O0C+HP/fn7KsWzy4JEMS78PtPXhjmffAQCoAjjaFcHf9eBrZ69e0jZ5FV/I58YRm23f6h04gDkU7xzrSfs6xxUVPreEkM+Dzr543gafuyLaB8unzpqKkD4Iu7s1/XQ/Y2rqCztaE44bPXBjwZXxWN/9+1YASFlZebQzisVT68166XUhL/ae6E1ZwEVDxwBegu58YRdCPjdG63VCJtUH8ZlzpsHvcdn2yqfbbPWVrmKe1ewbn057zjr/3ChDalelzrBH3yGoJ5IYwNN9Y/j1JxYN2L6Rwth8I5m2qtaFCr8btzy5Dd/X00y5Fpe1v5+R9/74GVOy2gc0mRHAr1owEWu+dZE5pdLoeFvfM9f+fj3W7T2BqaMqEI7KuHj2GLTrRbKMbyc0fAzgJWrm2Epzs9cbL5+LG98zF1UBLwJp3vynTawZ0vM8+/YxnPrdZ1KOJ9SF1v9dZqqa53FJGF8TsO2BJ1v19fPLYibCxbPHAEg/cyUuC3jdLrMnK5CfJbUxRcXlp47HFy+cDgD44RWnmM+Zzp4fvdu8/V59hyMjReJ2SRhfEzRz3pKeNLOu4Hx+ewveOtSFsdV+hKMyRlf2lwtm/e/cGfn/ihxo6bT6hF6KkUOsDnpsK8UBwCVzxw7puf616zi6bQY8tx7uL0D16OuHAGTe9iwclXHyhJqUAB60/GMdW93/jWIke15PaX39Em2lY7oB5esfeROb9reb36p87vwEtriioqHKP6ier5FGeeQLZ+GX+sKb5FW/Lj2gG/PZ7TZErvR70BOVzaqTAFhONof4SpaQdLnSCfqMjuqAN+0/Qrsl9X9auz9lWXuyd9LsQv7fz+wAkDgNLlPNlS/8cROisoK+mJKwYs/4Rz5/ci2isopPnDllxH+FNlaYnjyhBjdfcYptdT+rWWOrAABt4WjG6wZy5/M78cKOlpTjRgG0wfr7l8/Gwil15m07X18+C4363HG7OuKyKvC7f+1NmNpaiuWPnYoBvES0hWNpc9LzJ2v/iH76wXmDmnr37f/bMmChpFd3n8h4/pVd/bNF1AwRXAigpSuKroic0Os2fGTJZHT0xrNa4OJU39Kn33ktH2BVAU/aAH7l/EZ86cIZZorlYHvfsJ7/x/98B5/+XeKWe+GojK8/tHlIhcNOm1hre9vqqxfPxMS6YNr3xpP6DCRjhs7XLp6J8Q6oOukUDOAlIK6oWPBDbe6sXY/Z6NEa07YGsvfWy83bPWmmE2a72Ofpt/oLDylC4N6Xm/E/L+yyvfbcmaPR2RdDZSAxv7r31stxtb6v52mNQ8vVO8F150/H6zcuTyhGVuHzpE2hxPRBTEB7jQS06Zu5tLOlRyuvkOeSlTFFxbkzU+fHf+lCrR7OLVdqZXdnja0qyQJsTsUAXgL+8Wb/9LyuiJzQW7r74wuH9djdkdSvtUKIhIU/xjErv8eFKr8nYUs3WRH4w5p9eOyNQymPuXRaPa5c0IiO3nhCAEvmxMU7g5E8la7C70mZmWNIri+ycV87zvzR87bXZmORnu6wcuvBMtOU0VzQpkS6cMW8CQnHjW+PbpeEqxZMxJzxVXltR7kp7Z0AyoQ1J9jaHUlY6HDpyeOyegxj67KT9HzquOoAjnZFUoLHnuNhXPjj1ZhjWdF5RlM9FFUkPO+8SbVYMKUO+9u01XRXLZiI3a092N/Wi+QYLCsq1u5pg9ftQmdfHBUZAni5LeLQUij2wfPZt4/lZC/Qzt44akJe1AS9aGpI/JZmDDDa1fjOFUnSHt/rduFnH5qPv73R3yGxptPuuOb0vLWhXLEHXgKsb/Kr71qD1UkLKLLxhQtmYN13LsYz/3YeAOCJr54DAAkF94H+XKR1mzOfx41HNx3CT/SFOwavS0J9yIcKnxt3XHM6fG4XGmuDmDW2KiEtYAR5j0vK2AOfVB/EmKryyn9W+NOnUAAMe1FLeziG03+g1cF5bntLyvZ6RqoslsfCWT6PG1/58+vYfjR167ygz52Q0qPcGjCAS5K0QpKkZZIkrbQ5VytJ0gL9mtvy08SRL+TzmPOl221G8rNR6fckBEdjJsTOpDrT1tz3HVefDq9bgs8t4dltx/DopsS6K10RGfe/tg/fuPQkAMDmgx041NGH7Ue78fkHNprXGSkfj0vvgaeps/Hy9RcllGAtB5V+D473RFPGDYQQWDqt3nY+/L0vN+NHT23L+LidfXF09MbMVJhRT9zoDKiqwMH2XvPvHVfylwMPed14dfcJ7M2yXg7lTsYALknSCgAQQqzS7yfvIHANgEVCiIf18ylBngYWV1VzKlauGNP3kqd2WVMqV85vxObvXgKPS9ul5WB7n1nXAtB2WAH6Z1VELV/DuyyPE5VVjK32w+OW0BOVURkoj2Xy2aj0e/Dc9hZzWqbh1qe2p9R2n1invQce2XQILw7wLezbj27Byvs3mrlto3CUkY1bt7cN59z2gpk6icn5y4EnL9Z65Atn4qmvnYurF07M23OSZqAe+GIAzfrtZgALrCeFEPcIIe7R7zZZrqVBkBWRUFzo++87OSePe/tVp6UMVnZZBjVdLq0Oh9fjMgc7v/+4VtdCoH/lpZGjH1PlxwV6gSOjpjeg5T8vO2W8uVKv0j+y53kPRrpFK9btygxfvGAGLj9tfEJ6Kx1ZVdHVF8dFd7wIAPjMfVrVw0hcgRDC/JsZPfCZY/I3eJi8r+XCKfWYM74a/301c975NlAAr026P8ruIkmSmgC0GT31pHMrJUnaIEnShtbWwed2y4GiqglB76zpti/zoC2aWpdSjzo5Jw5oedi39aAxd0L/NL/zZ2lLwdv0us5Pfu1c3PWx/lkx33z4TfzjzcPmdltGOiDTIGa5STdlzi4n/ZGlk3HhSdprfrijD3/dcCDt7BGv24XtR1MXYW0+2ImHLD8XUwRuvuIU/IeeBssHowf+riwH3Cl3BgrgHQCy2S12hRDiOrsTei99kRBiUUNDg90lZe9XLzabuUu/x4WZY3PTWwr5PCkbBdgGcJcLxmzBzQc68KMntfzrZ8+dBgBYoE8FG13pT6hj8eCGA/jyn15HNK7A73GZUwQZwBNt/+G7MG9SbcKxc2aMxpXzG1OuNTZN6IrI+MbDb+Jwh/28cOsHw+SkcYVvPrIFn/jtOgBAPGmqYj4EfR54XBL+96MLBr6Ycmqgv+x69PfCmwA8m3yBJEkrhBC367f5FxyCzQc6sHhqPT57zrSc7qYTTNriCtBSKB6XhK3fv9Q8tvqdloQe4au7T8DvcZmDk2fafCOwBoXemKKlYtxGCoUB3CrgdZuB2VAd9OK9p49PuTZ5Vkq6aYbWb1bGBgyPfvGslOviigqvJ79TN0M+N4JetznuQoWTMYDrg5NN+uBlrWUw81n9/8sA3CZJ0kZJkjYiu9466Yy520um1uNz5zahJujNWG9ksII+d0oKZcfRbqz7zrKEXvKxLq0Gxz++ok09FBCoDnjTTnELeF0JAb8rEkel3w2PvtMKA3iq5EzK7U9vt/021NKdWA8lXSmbqKzgvFnaN1qjEmVDpR+7/+vdCde9tLMVch5noABaAK8K8G9eDAN+txJC3C6EWGX0svVjy/X/rxJCTBdCLNT/S8mBk737X9uHuTfpZVwlbUDRm4NFHVZet5SSaw1HZdSm2Q/xFH2Z+/HuGGaMqYQkSfiSXoLUavsPL0u4f9+rexHSv0YDTKFkY+GUOiyZltrfOWv6aPzsg/PM++mWwMdkFdeePRUAMLEuhNqQFzUhL9wuyRw4dbskPLnlKNY0Z653M1xBrxvVad5TlF9cyFMkN/7trZRjAxQOHDS7ATQBpHzV/e2nFuGf+gKgx798DlyStnsKAHzj0tT9EQ2nNGqrOTft70BtyGs+bqb9FstV8t+2Juizrck9Y0wlrpjfaNkkIX1KzVimDgBv3HSJWT/+6oWT4HFp+18CwPyk/HuuhXwe9sCLhAG8BBjh9JGkhTSFctHssWZJ04DXhcOdEay3meaW7OqFk/DuU7WZB7Wh/hog5bBZw2AJAfzbg2+Y92U1c4lXo+OdrgeuCpGSVzf88IpTsOu/3o0ZDZW4ZO5YXH7aBNvrckVLobAHXgz8l1ZkUVkx86PRPCy2cEtSxp10khm1us+w+XqfzOdx4c4Pa+PWtWWyx+VQxVUV/6dvjPHqruNYvaPVHDOwY6Sj7PLXu1t7sH5v+4AlYj1uCX1xJe97UAaZAy8aBvAiO9jeZ/a2zpiWm/nfVm8e7MCyn7xo3h/on7KRP71ozsA7/HjdLjNtYgxcsu6FPaMGzTV3r8F9a/YCyFwHZe23L8ayOWNsP3yNgWm3S8KX9XKtdjwuCX0xZUi1wAej0u8x0zdUWPzYLIK3DnWatx/ZeBDr9CXVt1x5KpY25TaIVwY82HM8jIfWH0DQ507YKs2O0QNP9/XcyghAm25cbuZbyV5zq7bp87o9bZiuVwzMVBd7VKUfbeEYVm07hgm1wYQtyXqiMr56kRa4My3QMXaVz3cAnzu+Gv+ubx9HhcUeeBEcsizO6IrE8Z+XzwGgpSRW5Lh+hDGPeO+JML716JYBt/cK+LTrBwrglX6PGVQYvAcneWpnJr98fhcW35I4uas7ImeVcx6l/13yXYPd5ZISxkCocBjAi2CNZRuzPcfDed2l2+hRN9YFs1ok5Pe48aHFk2y3RbPa8r1LcNb01B1YyN6DK88wb9ttIm3nhsvmmLd/vmonVuv7Xf7t9UNZTdVc/53k2nM00jCA51B7OJYQJPtiCv6ybn/Kdb9/da95u7k1nFW6Yqi+vlz7atsbVRKqCWZy61WnDTiThNtiDY41NSYrAp86a+qAP1NhqY/z01Xv4OWdx/HwxoN4YsuRlE017HBl5MjHAJ5Dl/zsJdzzUn9Bxj3Hw7jh0S0p100ZFcJD152J1f9xAY50RvLaAzdW64XT7I1JhdcXV8y9IjNJ/hYUV1T8x183A8h+quZtV506+AaSYzCA55CiioTiQ+lmic0dX40ZYyoxZZRWhCifAbzS78HDnz8zpagVFd6qr59n3q7JYuVi8vviD2v6a7VnOy75wcWTs7uQHImzUHIolFR7xJUmzRBXBLxuyUxDpKsZnSu1IS96ojJ8bhcm1JbXlmalxLobUTZ7YWYah1g0hWWHiD3wnAp43VnlmeOKmjC1yxhozJeqgBfdERnzJtdi9TcuzOtzUXp+jxs//9C8rK+3bpTwxQumm7NKAJTd1nRkjwE8h/weV8JqynRV4JIDeL574DVBLzp6Y3l9DsrOoqn1WDZnTFbXWge3K/wenAjzb0iJGMBzaGdLD6Kyil0tPfjV6t2Q0xQiUlSRMDd3oCl7w2V8M+CchOJrrA3i3k8uzupa60yffM/lJmdiAM8RIQRisorNBzrw1w0HcNvT282dwJMXzyT3ywtRfjUmq1nlXak0eSwBfHSlP8OVVE74LzpHmo9rS6W7IjLu1qcSyvqc8OOWIv2qKtDVl7hTfCECeE9U5kYLDtRYG4THJSUUi9rwn1ygQxr+i86B7kgcF9/xInweF2KWQUyjEJF1Ct9jmw8lbEb78w/NQ1UBAuuRjj6cN5N7kjrNv264CF2ROCp8HkwZVYEP3fNasZtEJYQ98BwwUt1fSVqcETcDeH8Kxet24ab3zDXvv39eY0FWzIVjCqcQOlR1QNtpp471RigJA/gw/eK5neYqxxpLTWyvW4KsqKgOeBC29MC7I3LRij+x5KezDaauO5UHBvBhWrXtGI73aDlua56yOuBFXBGoCXnRZ+mBR+JKXmufpLN0Wj2mjq4o+PNS7njyvDEDOQ8D+DAd746a83NDPg/WfediAEB10AtZVVEd8CIc7e+BR+JqXpfOp/PgdWfabqJLzjFrbBU233RJsZtBJYQBfBiEEDjeE0NbjxbAFVWYwbk25IWiCtQEvQmFpCJxBf48L9yhkauGW9eRBSPJMHRFZMQU1dwia2y1HwF9WbzX7UJcEWisDeJIZ8T8mYisFKUHTkQjDwP4MBzviUKSgDcPduIz50zDwin18Lol1Ia8kKDNAx9T7cevVu9GJK6lUfYd7zWDPBHRcDCAD5IQ/TMBjndHMaEmCKB/gEmSJLxx0yXojsgIxxRz6pexeOfprUfzXvuEiMoDI8kgTfvWkwjrS+Pbe2OYWKcFcF9SgeYl0+rx3LZj5rLnDsvqy3xvMktE5YGRJEudlgDcFdFu98YUc053clBuaqhAS3fU3Pj3I79ea56zlgklIhoqBvAsvHOsG6d//5/mfWO5fG9MMXfjTp6jG/S60RaOodLvwRlN9TjeE8Wulh4ALEZERLnBAJ6F5BVwRgCPxBXUV2jTupJTKEGfFsBDlt52Z18M153flOfWElG5YADPgjFuebC9FwDMXXd6YwqWTNN2G/ck1TNpGl0JQAvkfXHt+qt+tYbL2YkoZxjAsxDTy8Kec9sL5v2/vX4IfXHFLNGanNeeMUYL4CGfBz5LeuWdY90gIsoFBvAsROOJO7rHZBX/78E30BuVzRRJyJdYEtbYPKEq4MG9n1ycdrYKEdFQsR54FiJJGxUbOfBjXVEzgNttlvC7Ty+G1+1CTdCFMVV+dPXF8ZWLZua/wURUFtgdzEBWVHzqd+vQF0vtgQPAiXDU3M/SbledC0/q37y2MuDFtNEVmDyKu4kTUW4wgGfwbw9txuodrdh7ImweO29Wg5kTX7+3HUGfGy/8xwVYNKUu42NV+T3ojsgZryEiGgymUDJobtXmbb99uMs81lDpT9g2Leh1Y1oWdbY9bsmcvUJElAvsgWdgTP/e39aLn39oHgBgTLUfv3x+p3mNJ8tByevfNRv3Xbsk100kojLGHngGxtTuNw50QJK0O+OqA9jdGs7wU/Yaa4O5bBoREXvg2VJUFS4JGFejbQx8wUnc4Z2IiosBPAOP22XOMlFVoPlHl6MmqK2kDHJTBiIqMgbwDIJeF978nrYHoaqvpzcCN3fVIaJiYwDPQAitTOxXLpqB+ZO1aYLGknmXxB3Ciai4OIiZhX+/5CTzNpfCE1GpYAAfJLc+NeXiOWNw6clji9waIipnDOCDNL4mgJ9/aB4uO2WcObWQiKgYGMDTUFRhm+f2uF14/7zGIrSIiChR2Sd0I3ElYWm8oTcmJ+ymQ0RUaso+gM++8Wl84YGNKcfDUcW2wiARUako+wAOAGuaT6QcC8dkBnAiKmkM4AAm1aXW6A5HZVQwhUJEJaysA7gQAo21QTTWpRaa6omyB05Epa2sA/gru47jUEcf4oqKlu4IFKN+LLQcuN02aUREpaKsA7ixVZqsCCy55Tk8ueWIea43JiPkZwqFiEpXWQfwlfdrs0+MQlUdfXHzXE9UZg+ciEpaWQdwg5E4ufFvb5nHtEFMBnAiKl0DBnBJklZIkrRMkqSVGc4/m/um5d915zfh9qtOw5sHOwAAF83u30U+HFWYQiGikpYxgEuStAIAhBCr9PvLkq8RQjycn6bln0uSMH1MJSJxbSVmb6x/1/gwUyhEVOIG6oEvBtCs324GsCC/zSmsuKwmlIcV/ZNQuJCHiEreQAG8Nun+qME+gSRJKyVJ2iBJ0obW1tbB/nhexRUVXo99RcFwVGEOnIhK2kABvANA/XCeQAhxjxBikRBiUUND6WwE/OD6/bhvzT5402zQ8PfNh1HBHDgRlbCBAvh69PfCmwA4crAy2YPr9+Obj2wB0L/Dzg/efzKE0KoT7m7tAQCE2AMnohKWMYDrA5RN+uBlrWUw0wzk+rlFxoCnE/ziuV3mbaMH7pIk1FV4MfvGp3HxHS8C6N99h4ioFA3YxRRC3K7fXGU5ttxyexWAutw3LX9Uy2il160F6eTNGyba1EchIiolZbeQZ23zCRzpjJj3vR4XfvepxbjgpAZzFsoH5jfiro8tLFILiYiyU3ZJ3qNdkYT7XpcLF+oLeCyzCFEb8hawVUREg1d2PfCgN3FmiZFCAfrngffGFAS8nIFCRKWt7AJ4cq7bbqCyoy/GAE5EJa/sAriRJhlXHcDeWy+HlBDQtbNt4RgCnrJ7aYjIYcouShk70I+q9KWcu3L+RABAd0SGJ80CHyKiUlF2USquaAG8viI1gF9+2nh86qyp6InIKeeIiEpN2QVwoweertJgbciL7igDOBGVvrIJ4DFZxdQbnsCOY904bWINFky2X3tUF0rtmRMRlaKyCeBG6uQ3r+zB+bMa8Lnzmmyv4/xvInKKsgngstK/TEeS0tc4YQ+ciJyifAK4qpq33RkCODdxICKnKKMA3t8Dv2bxxLTXjasJYNmcMWnPExGVirIJ4EYOHADG16SvNNhYG8S9n1xciCYREQ1L2QRwIwce8nGJPBGNDOUTwPUcuKKKAa4kInKGshmxiysCy+eOxWfPmVbsphAR5UTZ9MC7IzKaGiqwtGlUsZtCRJQTZRPAr7l7DR5Ys6/YzSAiypmyCeAAEI4pxW4CEVHOlE0OfPncsbjgpIZiN4OIKGfKpgfeWBvE/En2BayIiJyobAJ4VFbh95bNr0tEZaAsIlpUVvDndfvh4y47RDSClEVE+8Or2uwTFqoiopGkLAJ4VNZmn1T4uYyeiEaOsgjgQl897/cwgBPRyFEWAVxm/RMiGoHKIoBbS8kSEY0UDOBERA5VFgE86GXum4hGnrII4JIk4U+fXVrsZhAR5dSICeA7j3WnPXe0M4Lxtem3USMicqIRE8CX//Ql9ERlXHf/hpRzhzv7ML4mUIRWERHlz4gI4EKf6N0ejuGZrcfM+4aorCLAPDgRjTAjIoC39kQBAJ19cQBAS7d2/2B7Lz7+m7V461Bn0dpGRJQvIyKA/+DxtwH0B/DbntoOADjnthfw8s7j6OVGDkQ0Ao2IAL5kWj0AoKsvjoDXhTnjq7n7PBGNeCMigMdkFeOqA/jCHzdhVIUfPVGZi3eIaMQbEQG8L6aY6ZOQz424ojKAE9GI5+gAHldU7DzWjaisYpw+TbDC70FMViErTKEQ0cjm6AC+trkNy3/6Eu58YRfuuOZ0VPjcqPBrPfBDHX3Fbh4RUV45eosaj1syb/s9LsiqQMjnQUxRcfdLzQCA0ZU+fPGCGcVqIhFR3jg6gHuTAnhUViGEQEwWmDWmEgBw7ycXY96k2iK1kIgofxydQrGOU/rc2krLvSd6EVNUTKgN4varTmPwJqIRy9EBPCb3R3Cfx4VzZ46GogrEZRVRWYXf6+hfj4goI0dHuJjSv8LS53HB73FpKRRFRVRW4HM7+tcjIsrI0REuuQd+ydxxeN/pExBXVPTGFIT8jk7xExFl5OgAHolbArjbhWsWT8K/LZ+FmKyiOyKjKsAATkQjlyMC+NHOCNbtaUs5/v8efMO8bcxIkSQJAkBPNI5K9sCJaARzRAA/1hXBqm3H0p6//zNLIEn9UwrX7WnDrpYe7oVJRCOaIwJ4wOtGJJ5YEjYclc3b585sSPmZ15rbOAuFiEY0R0S4gNeVNoA/+dVz0/6c38MeOBGNXA4J4O6EAUtA27zh/fMmYO6E6gw/54hfj4hoSBwR4QKe1BTKd/72Fh5747Dt9d9771wA4DxwIhrRBpymIUnSCgAdAJqEEPcM9nwu+L0uROTEHviSqfWYNqrC9nqvx2W0LR/NISIqCRm7qHpwhhBilX5/2WDO54rf40I0rkC2FD+584VdmDWuyvZ6L3veRFQGBop0iwE067ebASwY5PmckCQJa/e0Ye5NzyQcnzvePv/t9zCAE9HIN1Ckq026P2qQ5yFJ0kpJkjZIkrShtbV1cK1LErP0wJfNGYuTG+0DeGNtcFjPQ0TkBAMF8A4A9cM4DyHEPUKIRUKIRQ0NqfO1sxXyuY3HAwAoqgqPyz7HPSVNbpyIaCQZKICvR38vuwnAs4M8nzMz9Q0a9rf1YuoNT0BWBVxpBikbqvzYectl+WoKEVFJyBjAhRAPA2jSBydrLYOVz2Y6nw+fObcJANAd0RbwqEKk7YEDHMgkopFvwGmEQojb9ZurLMeWZzqfD+fNHA1AW8ADALIi4M4QwImIRjrHdFONYL3tSBcAoCsic543EZU1xwRwIyVy8xPbAPQHciKicuWYAM50CRFRIscEcA5KEhElclRUnFATKHYTiIhKhqMC+IqFEwEA5+ozUoiIypmjNo289JRxGFsTwPSGSry883ixm0NEVFSOCuAnT6jByRNqoKqi2E0hIio6R6VQDC6XhJ99cF6xm0FEVFSODOAAcMX8xmI3gYioqBwbwImIyh0DOBGRQzGAExE5FAM4EZFDMYATETkUAzgRkUMxgBMRORQDOBGRQzGAExE5FAM4EZFDMYATETkUAzgRkUMxgBMROZQkROFqa0uS1Apg3xB/fDSAUtzFge0aHLZrcEq1XUDptm0ktmuKEKIh+WBBA/hwSJK0QQixqNjtSMZ2DQ7bNTil2i6gdNtWTu1iCoWIyKEYwImIHMpJAfyeYjcgDbZrcNiuwSnVdgGl27ayaZdjcuBERJTIST1wIiKycEQAlyRphSRJyyRJWlmk56+VJGmB3o7bLMfbJUl6VpKk64vV1mzbUMh26a/VbkmSNur/3TaYtuahPSskSXrW5tiAr1G+25fctlJ5r6V5zYr+XrN5vYr+XsvwN8v7e6zkA7gkSSsAQAixSr+/rAjNuAbAIiHEw3objBf6aiHEciHE7UVs64BtKEK76oUQ04UQCwF8DsDd2bY1H40x/m6GbF+jQrQvuW0okfeaTbuyakMR2lUK77WUv1mh3mMlH8ABLAbQrN9uBrCg0A0QQtwjhDAGIJos7amVJKnJcmkx2ppNGwraLuPNqGsSQpTS65XueYv+ugF8rw1WKbzX0vzNCvIec0IAr026P6oYjQAA/Q3RZnnT1ANokyTJ+NSvTfqRQrQ1mzbYHcs7SZJWJvWYSuH1Sve82R4rCL7XBqcU3mtJfzO758v2WNacEMA7oP0xSsEKIcR1xh39k7cDQIf+VagDBW5rlm0oeLt0y613SuH10tk9b7bHCoXvtcEphfea9W9m93zZHsuaZ6g/WEDr0f8p1QTg2fSX5o8kSSss+bQFABYB2CCE2GS5rKBt1fOj2bShtpDt0ttWm3Q/27YWQravkd2xvON7bdBtq026X/D3ms3frCDvsZLvgetfi5r05H5tUs6rIPTnvs0Y6Yb2ifmQfs4YhHi4CG3Nqg1Feg3rAbQNtq35aIj++Iusz538vMV63ZLbVirvteR2ZduGIrQLKPJ7ze5vVqj3GBfyEBE5VMn3wImIyB4DOBGRQzGAExE5FAM4EZFDMYATETkUAzgRkUMxgBMRORQDOBGRQ/1/BXSW4d2/FVkAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(Mavgs);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "4a6a8b87",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXkAAAFkCAYAAADFS0t4AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8/fFQqAAAACXBIWXMAAAsTAAALEwEAmpwYAAAPL0lEQVR4nO3dsXIbWXYG4HPsjZwYQw1zFRQ450JvAEdyyJpH4LyBd+YJVtIbDF9hlE5EvgGpSex0WdrAEbdU3MRVDlxnAzWHWCxFEE2QIA+/r0ql7otu9kUT/HlwbzeYVRUA9PRP2+4AAPdHyAM0JuQBGhPyAI0JeYDGhDxAY7/bxkG//fbbevny5TYODfBkffz48S9VtbvOPlsJ+ZcvX8bp6ek2Dg3wZGXmn9fdx3ANQGNCHqAxIQ/QmJAHaEzIAzQm5AEaE/IAjQl5gMaEPEBjQh6gMSEP0JiQB2hMyAM0tpVPoWS1lz/88tvyp7dvttgT4ClTyQM0JuQBGhPyAI0JeYDGhDxAY0IeoDEhD9CYkAdoTMgDNCbkARoT8gCNrQz5zNzPzKOltr2hfX9pu3lmHtxHRwFY38qQr6oP1zT/OLTvZOb0Muyr6jgiIjPnm+0mAGOsPVwzVOonmTmtqsOqOouI1xFxNmxyFhF7G+wjACONGZN/FREvIuJzZv6UmZOImCxt82J5p8w8yMzTzDw9Pz8fcVgA1jV24vVPVXURER8j4iAiLiJi56Ydhqp/VlWz3d3dkYcFYB1j/mjISVwF+iS+BPxZXFXz04g4Wt4JgId3m6tr5hExW5hc/RARk8vJ1aFC/xAR06FtcjkBC8B2razkh8D+Zqnt/bB4fFMbANvlZiiAxoQ8QGNCHqAxIQ/QmJAHaEzIAzQm5AEaE/IAjQl5gMaEPEBjQh6gMSEP0JiQB2hMyAM0JuQBGhPyAI0JeYDGhDxAY0IeoDEhD9CYkAdoTMgDNCbkARoT8gCNCXmAxoQ8QGNCHqAxIQ/QmJAHaGxlyGfmfmYefeWxd0vbzTPzYJMdBGC8lSFfVR+ua8/MeURMh+X9YdvjhccA2LJRwzWZOY2Is4Wm1wvrZxGxd8d+AbABY8fkp1W1GPKTpcdfLO+QmQeZeZqZp+fn5yMPC8A61g75zJxfDsssuIiInZv2q6rDqppV1Wx3d3fdwwIwwu9G7PN5GHOfRMQ0M/ci4iSuqvlpRFw7UQvAw7rN1TXziJgtTK7+OlTyOzEE+zA5O70M/2sqfQC2YGUlPwT2N9e0H0bE4cL6+2FRwAM8Em6GAmhMyAM0JuQBGhPyAI0JeYDGhDxAY0IeoDEhD9CYkAdoTMgDNCbkARoT8gCNjfmoYTbo5Q+//Lb86e2bLfYE6EglD9CYSv4RWazqATZBJQ/QmJAHaEzIAzQm5AEaE/IAjQl5gMaEPEBjQh6gMSEP0JiQB2hMyAM0JuQBGhPyAI0JeYDGVoZ8Zu5n5tHC+iQz94b2d0vbzTPz4L46C8B6VoZ8VX1YavouImaX7Zl5kJn7w7bHQ9t80x0FYH1rD9dU1WFVHQ6r04g4i4jXw/8x/L+3me4BcBejx+QzcxoRn4fqfbL08Itrtj/IzNPMPD0/Px97WADWcJeJ1/2q+n5YvoiInZs2Ht4BzKpqtru7e4fDAnBbo/7Ga2buV9X7YXkvIk7iqpqfRsTRV3bljhb/Duynt2+22BPgKbjN1TXziJhdTq4O6+8y82NmfoyInWESdjo8NrmcgAVgu1ZW8kNgf7O0/uqa7d4PiwIe4JFwMxRAY0IeoLFRE688rMXJVoB1qOQBGhPyAI0JeYDGhDxAY0IeoDEhD9CYkAdozHXyW+C6d+ChqOQBGhPyAI0JeYDGhDxAY0IeoDEhD9CYkAdoTMgDNCbkARoT8gCNCXmAxoQ8QGNCHqAxIQ/QmJAHaEzIAzQm5AEaE/IAja0M+czcz8yja9rmmXlwUxsA27Uy5Kvqw+J6Zu4P7cfD+vy6ts13FYB1jRmueR0RZ8PyWUTsfaUNgC0bE/KTpfUXX2n7O5l5kJmnmXl6fn4+4rAArGtMyF9ExM4t2v5OVR1W1ayqZru7uyMOC8C6fjdin5O4qtynEXE0rC+3AbBlt7m6Zh4Rs4XJ1Q8RMR3aJ1V1fF3bvfYagFtZWckPgf3NUtv7YfH4pjYAtsvNUACNCXmAxoQ8QGNCHqAxIQ/QmJAHaEzIAzQm5AEaE/IAjQl5gMaEPEBjQh6gMSEP0JiQB2hszB8N4ZF4+cMvvy1/evtmiz0BHish35xfBPC8Ga4BaEzIAzQm5AEaE/IAjQl5gMaEPEBjQh6gMSEP0JiQB2hMyAM0JuQBGhPyAI0JeYDGRn0KZWbuR8RFREyr6vBrbQBs19qVfGbOI+Ksqo4j4iwz94aAj6HtchsAtmzMcM1pRPycmXvxpWr/NSJeR8TZ8PhZROxtqH8A3MHaIV9VFxHxU0T8HBG/H5onS5u9WN4vMw8y8zQzT8/Pz9c9LAAjjB2uOa6qVxFxsTAWv3PTflV1WFWzqprt7u6O6iwA6xkzXLM3DNFERPwxvoT7SVxV89OIOLp71wC4qzFX1xxm5kF8GXtfvLrmP4cqf3I5AQvAdq0d8sOY/D9cIllV74dFAQ/wSLgZCqAxIQ/QmJAHaEzIAzQ26rNrWM/LH37ZdheAZ0olD9CYkAdoTMgDNCbkARoT8gCNCXmAxoQ8QGNCHqAxIQ/QmJAHaMzHGjSx+NEJn96+2WJPgMdEJQ/QmJAHaEzIAzQm5AEaE/IAjbm6piF/pAS4pJIHaEzIAzQm5AEaE/IAjQl5gMaEPEBjQh6gsVHXyWfmXkRMIyKq6sPQth8RFxExrarDTXUQgPHGVvI/DuG+k5nTIeCjqo4jIjJzvqkOAjDe2iGfmQcRcZKZ06o6rKqziHgdEWfDJmcRsbfBPgIw0phK/lVEvIiIz5n5U2ZOImKytM2L5Z0y8yAzTzPz9Pz8fMRhAVjX2OGaP1XVRUR8jIiD+DIWv3PTDkPVP6uq2e7u7sjDArCOMROvJ3EV6JP4EvBncVXNTyPi6I79AmAD1q7khwnXyeXk6lChf4iI6dA2uZyABWC7Rl1CWVXvh8Xjm9oA2C43QwE0JuQBGhPyAI0JeYDGhDxAY0IeoDEhD9CYkAdobNTNUDxNL3/45bflT2/fbLEnwENRyQM0ppK/J4tVM8C2qOQBGhPyAI0JeYDGhDxAY0IeoDEhD9CYkAdoTMgDNCbkARoT8gCNCXmAxoQ8QGNCHqAxIQ/QmJAHaEzIAzQm5AEaE/IAjd0p5DPz3cLyfmbOM/Pg7t0CYBNG/43XzJxHxHRY3o+IqKrjzDzIzHlVHW+oj0+Gv+sKPDajKvnMnEbE2ULT64X1s4jYu2O/ANiAscM106paDPnJ0uMvlncYKvzTzDw9Pz8feVgA1rF2yH9lKOYiInZu2q+qDqtqVlWz3d3ddQ8LwAhjxuQ/D+Pxk4iYZuZeRJzEVTU/jYijjfSOe7M4f/Dp7Zst9gS4T2tX8lX161DJ78QQ7FX1Ib4E/jwiJs9x0hXgMRp9dU1VHUbE4cL6+2FRwAM8Em6GAmhMyAM0Nnq4hufFRC08TSp5gMaEPEBjQh6gMSEP0JiQB2hs61fXPJWrNp5KP58S5xTun0oeoLGtV/JPUecKdFvPrfM5hW1SyQM0JuQBGhPyAI0JeYDGhDxAY0IeoDEhD9CYkAdoTMgDNCbkARrzsQZ81eJHDTw2PgYBbkclD9CYSv6OHnO1e19U0fB0qOQBGlPJ3+A5VulAL1sJ+f/6n7+uDFBDAg/HLzPoy3ANQGPParjGu4PtcN5he9YO+cycRMR0+Pe6qv4wtO9HxEVETKvqcIN9BGCkMcM130XErKo+RERk5sEQ8FFVx0PbfHNdBGCstUO+qg4XKvVpRJxFxOvh/xj+39tM9wC4i9ETr5k5jYjPQ/U+WXr4xTXbH2TmaWae/v///nXsYQFYw12urtmvqu+H5YuI2Llp4+EdwKyqZv/8L/96h8MCcFujQj4z96vq/bC8FxEncVXNTyPiaCO9A+BO1g75YVL1XWZ+zMyPEbEzTMJOh8cmlxOwAGzX2pdQDgH+6pr298OigOer3F0LD8sdrwCNCXmAxlp+rIHb6AG+UMkDNPaoKnmTcn34XsLjoJIHaOxRVfL34WsVpXF74DlQyQM01r6S5355RwSPm0oeoLEnUcmrFgHGUckDNCbkARp7EsM1EOEGKxhDJQ/QmEo+TOwCfankARpTybMx9zFmbhwe7kYlD9CYSn6JyhHoRCUP0JhKnkdn3XdTro6Cr1PJAzQm5AEaE/IAjT3pMXlXwgDcTCUP0JiQB2hMyAM0trEx+czcj4iLiJhW1eGmvi4A420k5IeAj6o6zsyDzJxX1fEmvvYyk60At7ep4ZrXEXE2LJ9FxN6Gvi4Ad7Cp4ZrJ0vqL5Q0y8yAiDobV//vzu//47w0d+6n7NiL+su1OPBJ3Phf5bkM92T6viyvOxZV/W3eHTYX8RUTs3LTBME5/GBGRmadVNdvQsZ805+KKc3HFubjiXFzJzNN199nUcM1JXFXz04g42tDXBeAONhLyVfUhIqaZOY+IyX1NugKwno1dQllV74fF2wS8SyyvOBdXnIsrzsUV5+LK2uciq+o+OgLAI+COV4DG7v1TKFfdCfuc7pS96blm5iS+TFpPI+J1Vf3hwTv4gG77fc/Md8/9XGTmXnx5XVzOf7UlL64Mz/X7qvr3Gx6/iBXn4l4r+cU7YYf1+TqPd3KL5/pdRMwuf4iH+wpauu33fWifPmDXHtwtz8WPw+tiJzPbno9b5MU8Is6Gx8+GX35t3fQLfZ3svO/hmlV3wj6nO2VvfK5Vdbjw23i6sG1HK7/vQ5h1PgeXbjwXwy/7k8ycDq+Rzudk1eviNCJ+vnxnU1W/PmTnHplbZ+d9h/xkaX35TthVj3cyWVq/9rkO4fa5+WWok6X1687FtHmgXZosrS+fi1dD2+fM/GkY1utqsrT+d+eiqi4i4qeI+Dkifv8wXXq0JkvrX83O+w75i7j5TthVj3dyEbd7rvtV9f0992XbLuKGc3GfH3D3CF3E6tfFn4aA+xhXHw3S0UWseF1ExHFVvYqIi8shi2fqIm6Znfcd8qvuhH1Od8qufK6ZuX95v0Hz8cZV5+JzZs6HH+LpMz8XJwvLk/jyw93VqnOxtzBE88d4PgXidW6dnfca8l+7EzYzj256vKNV52Jof5eZHzPzYzR+Ad/idfHr0LYT//i2tJVb/oxMLifWOl9RsupcRMTh5UeZR8R3nc9FxG+ZMFt8xzImO90MBdCYm6EAGhPyAI0JeYDGhDxAY0IeoDEhD9CYkAdoTMgDNPY3rUBNJ0spPVEAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.hist(np.abs(Mavgs), bins=50);\n",
    "plt.xlim(0, 1);"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.10"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}