{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Ott-Antonsen reduction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "import jax; jax.config.update(\"jax_enable_x64\", True)\n",
    "import jax.numpy as jnp\n",
    "from jax import vmap, random\n",
    "from jaxkuramoto import theory, Kuramoto, odeint\n",
    "from jaxkuramoto.solver import runge_kutta\n",
    "from jaxkuramoto.distribution import Cauchy\n",
    "\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# number of oscillators\n",
    "n_oscillator = 10**4\n",
    "# define natural frequencies\n",
    "loc, gamma = 0.0, 1.0\n",
    "seed = 0\n",
    "dist = Cauchy(loc=loc, gamma=gamma)\n",
    "omegas = dist.sample(key=random.PRNGKey(seed), shape=(n_oscillator,))\n",
    "# define coupling strengths\n",
    "K = 3.0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "model = Kuramoto(omegas, K)\n",
    "# define initial conditions\n",
    "init_thetas = random.uniform(key=random.PRNGKey(seed), shape=(n_oscillator,)) * 2 * jnp.pi\n",
    "# define time points\n",
    "t0, t1 = 0.0, 50.0\n",
    "dt = 0.01\n",
    "# solve the Kuramoto model\n",
    "sol_full = odeint(model.vector_fn, runge_kutta, t0, t1, dt, init_thetas, model.orderparameter)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "model_oa = theory.OttAntonsen(dist, K)\n",
    "# compute the Ott Antonsen solution\n",
    "sol_oa = odeint(model_oa.vector_fn, runge_kutta, t0, t1, dt, 0.01*1j, model_oa.to_orderparam)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "orderparam_theory = theory.orderparam(K, dist)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.legend.Legend at 0x158823ac0>"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAgoAAAF9CAYAAABlDkCxAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAABnd0lEQVR4nO3dd3hUVfrA8e+b3kintwTpoCCCgogUQVAU7OvPjl12LatrwYIUu6trWRvrqthdFQuIoEhRFJQinVATSigJJKSRnvP7406GVDIh96a+n+eZZ5Jz7z3nzGXIvHOqGGNQSimllKqIV11XQCmllFL1lwYKSimllKqUBgpKKaWUqpQGCkoppZSqlAYKSimllKqUBgpKKaWUqpRPXVegNkVHR5uYmJi6roZSSilVK1atWnXIGNO8Jnk0qUAhJiaGlStX1nU1lFJKqVohIrtqmod2PSillFKqUhooKKWUUqpStR4oiEh7EflCRNJEJF1EZolIh2pc30NEPheRQyKSLSJbRORuJ+uslFJKNVW1OkZBRIKAhUAucD1ggCeARSJyijEmq4rr+7uuXwzcDKQBXYAQB6utlFJKNVm1PZjxFqAT0M0Ysx1ARNYB24DbgBcru1BEvICZwE/GmItLHFrkXHWVUkqppq22ux7GAcuLgwQAY0w88CswvoprhwE9OU4woZRSSil71Xag0AvYUEH6Rqwg4HjOcj0HiMhyEckXkSQReUVEAm2tpVJKKaWA2g8UIoHUCtJTgIgqrm3jev4M+AEYBTyHNVbhY7sqqJRSSqlj6mLBJVNBmnhwXXFQ86ExZrLr58Ui4g08IyI9jTGbymUscitwK0CHDh5PrlBKKaUUtd+ikIrVqlBWBBW3NJR02PX8Y5n0H1zPfSu6yBgzwxjT3xjTv3nzGq1iqZRSSjU5tR0obMQap1BWT6Bca0AF10L5Foni1oiiGtRLKaWUUhWo7UDhW2CgiHQqThCRGGCw69jxfI+1/sKYMumjXc+6iYNSStnoUGYuyRm5dV0NVcdqO1D4D5AAfCMi40VkHPANsAd4q/gkEekoIgUiUjwWAWPMYeBp4HYReUpERorIQ8BkYGbJKZdKKaVqbuBTPzHgyQXsTM6s66qoOlSrgYJr5cURwFbgA+AjIB4YYYwp+U4UwLuC+k0DHgCuAOYCdwDPYy3kpJRqoo7mFXDPp3+y70i2o+XkFxbx5+6qhlPVvn1Hshn8zEJSsvJsyy8pI4eCIqund8QLSzCmonHotS/+UBaZuQV1XY0mpdb3ejDG7DbGXGqMCTXGNDPGXGSMSShzToIxRowxU8qkG2PMi8aYzsYYP2NMR2PMZGNMfm2+BtW0HM7MZUdyJi/+sIWYh77j8jd/qxd/NDfvT+f1xfWzIS23oLBWy5s2exNfr9nHpFnrbf23+W37IS574zey86zXc8+na7j49d/4z887bSvjROUVFBHz0HfEPPQdZz6zkMQj2fSbXnas94k585mFnP7kT6XSNu1Pd/+cW1DI5yv3kFdg39CwNXuOcOGrS0nLrvzPuTGG4f9cTO/H5/OPz9dSWFT3/w+LzV67j7gD6eXS1+w54vjfi8Iiw+87rfH+xhiu/e/vPPndJmIe+s6W/KU+/MGrLf379zcrV+pQBlU95/5rCVsPlm96TXhmbJXX/rr9EAG+3pzWsaplQiq3N/Uo4UF+hPiXns1c/Edg2aQRtA5zZs2xDYlpnNQ8hEA/b4+v+e/SeKbP2cT8e86mW6tmHl9XVGSYuSyBS/q1IzTAB5GqZ03vPnyUyBA/ej8+35326Nge3HRWLLd/uIpzurfkigHtq8znUGYuOfmFtIsIcqcVFBbR+ZHv3b/7+3iR6/pgHNatOe9NOL1UHuk5+dz03gruPqcrHSKDOPt5a3X5L24fRP+YiiZ71cy4fy9l3d60cunv3jCA4d1bnFCe+YVF7E45yjkvLHGn9WkXxtoS5dx1Thde+Wmb+/fTYyJ5ZGwPurduhr+P5+8TgIRDWTz/wxaeu/QUepX4N6zo/9ZDX67j0xV7SqVFh/jTr0M4M67rX2kZv+88zE0zVzJx+En8pX97okL8q1VHTxQWGU56eC4Afzx8Di1CA5g6eyPv/poAQKCvN1EhfvzywHCP3teeys4rJDu/0B0gvnxlX2Kjgxn371/d5+x69oJVxpjKb5AHNFBQ6jj2p2Uz6OmFFR7b8dT5eHsd+08//J+LiT+UxaTzutM/JoJL31jmPuZJUFGRP3encvHrvwHw/GWn0DoskBUJKdx1Thf3H6aRPVrw9vUDPMpv3oYD3P7hKgDGntya167uV+4cYwxdH/2e/MJjfxtio4NJPJJNXkERP98/nPaRgaX+4OUXFrF2zxFu+2AVh0s0f981ojPXnRlDRJAfhUWGuev3M75vmwr/WP4Rn8IVb1n37C/92/PsZadU+BryC4vIzCngVA+/Pa98dCTRVXw4FAddL1/Zl/F92wKwaldKqX/Dsp6/7BQu6dcOL4GUrDw+W7mH5+ZtqfDcdycMYHg3zz+8F8Yd5Lfth3lkbI9S9yojJ58QfyuI6vLI3FL/Rn3ah7N2zxGg/Pttf1o2CYeOMuikqArLW77zMG3CApkyeyML45JKHZtz51lc8OpSj+o9dVwvrh3YES+v8v++iUeyCfH3ISzQ15028aNVzF1/ABEo+VFUtv5VfTM+3v+vstee6P/Fkq575w9+3prMExf15pqBHVm1K5VL3/jNfTzA14uc/PKtLRef2pZnLj2ZnLwiBj79E29eexpDuzbnUGYuIf4+BPh6HmilZefTZ+oPpdLO7tqcZgE+fLduvztNA4Vq0kBBVYcxhthJc0ulLbl/GK/8tJ0vV+8F4PPbBzEgJpKiIkOnh+dWlI1b3PQxFBlDz8nHvjk9d9kpnNw2jB6tQzHGuD8UjDFc8sZv/Ln7iEd17dMujDevPc3dsrDrcBZPfLeZ167qh5/PsR7G059cQFKJUezvTRjAsDIfYMkZuQx4csFxy7v7nC78fVRX9+/VaeK8/LR2PH95H/fvB9JymL12H+k5+by68FhXypYnxlT4DfXMp39iX1qOx+WB9WHXu22Y+/dftiWzfOdh7h/dnUVxSUx4b4X7WNz0Mew6fJTRL/0MwCv/dyp3ffJnhfm2jwzknO4tee+3hCrrsPbxc0t9SJaUlVvAj5sOMq5PG3ILiugxeR4AP98/nA5RQRzKzMUYGPDkAq4d2JF/nNuNPtN+4K4RnWnezJ+ebULp1SaM7o9Z161+bBSRwX7u/Hs8No/s/GPdQQnPjGVFQgqntAvD38fb/e8XGx1M/KFjm/he1LcNL1zR1x2UeuLZS0/mLwNKL263aV8657/yy7F7MflcwoJ8+evHq0t9qBV7/ep+dG4RwsH0HOauP8Anf+x2H3v3BisoLvlvNuXCntwwOLbC+pR9b94x7CQeHNPdo9eSW1CIr5dXqcAnO6/Q/e8D8PD53XlqbpxH+UH5IGLt5HPpM+0HhndrzrtlWqlKKlludYI3DRSqKbprtLnwlQtLpY2OGc2V3a8kuyCbiQsmlrtmfOfxXNT5IlJzUrl38b3ljv+l218YEzuGA1kHmPTLpHLHr+91PcPaDyM+LZ5py6aVO37rKbcyqM0g4lLiePaPZ8sdv7vf3fRt0Zc1SWt4efXL5Y4/ePqDdI/szrJ9y5ixbka545MHTSY2LJbFexYzc+PMcsefHvI0rYJbMS9+Hp9t+azc8ReHvUhEQARfb/+ab7Z/U+746yNfJ9AnkE/jPmV+wvxyx98d8y4A7214jyV7l5Q65u/jz5sj3wTgzbVv8vv+30sdD/cP51/D/wXAS6teYm3y2lLHWwa35JkhzwDw7B/PEpdS4j+rKaJjUCumnHIHFOQyZe1r7MraB6bI/ege2JIH240GU8SD8bM4kJtm/UEwBjCcJNEU/BmBN0UcPnU3aYXWh1NeQREH0nNodTSMU1La8/dRnZm4/0cSjmRRcpmPdlkR9Eq1Vh6f324jZcVmRNIjrSUFUsiPbbcS7O9NaIAv3l5CbkER4YkhdMuIJscrnwWtd5S7vkdaC07KjCTTJ4/FLXcS7OdDhOvDYW/qUU4+0oqOWeGMO7sl01JWUFBo1bvYqSltaJsdSkhrWBpV/IfYkJ5dQHpOPgMOt6NlTggHAzJZEbW3XPm3+PXmik6xLD26j6m7V5U6FhHkR6+ENoTnB7Ar+Ajrww+UOt42PJA7vPrQrNCfN3dvZnNYEj5e4h48B/BsxJn0axnF15k7+SYrHrBaEw66XsOYfV3wMd5sCktiZ0gKzZv5k19YxJGjVh/33/P7sy0pk3XhB0hvfpQQfx8ycvJJy87Hx3gxZl9Xrh/UkXu2LCMxsHTfckCRDyP3d7buyGlZ/JlzyD1QsmVoAJkpRQw/aM3yXha9m8P+R0tdH5YfwPOtBzF3/QEWRceT5mvVubhro7tfOA9GnGb1J29fyMHCY9e3aBaAJPsw4HA7urVsxmtef5LjVXrwXsfcME4+1JoLTmnN82YlucYKApIzcsktKKRDVjj/PW0IABMOLmBvaulBnp0yI+mZ1oIubYL5OGyj+54CiMB10d0YHxxLjlcBD6T+Rl5BId5eXuTkF3IkO4/uR1owsWMPUsjh8SOl/9/6entxf9s+DA9qR3x+OtNSVpCVW0Dq0WMtTaemtOG5M0/j7Y1b+cz7WCtMsJ8PWXkFlb73wgJ9aRbgy4MR/Vj422ES/NL4M3Kf694GYoz1HnmixUBifUP5KXMPz+xdQ2iAL/4+XiRnWkHyB51HsHz9ERIjj7Ay4ADGGA6k5xAR5EeArzfPRw5mxZYjzDkaT0JkCi2aBbjrkJyRy/DdnUq998q6INEKRJI6HGGjzyGySgy+LH7vAfwZua/ce69XVDgvNbf+7V46soa1uYfBGPa63n/BBX5VvveGJMUQ6OvN2o6JvD5hXo0DhbpYwlmpE1eYB/lHoSAXMg/D7HsgYz/kxIPJgaICKCq0goH8AlhkBSJERYJvmbd73lr4/X8ASPMovLxLf3sNyd3EPb5W3+zfk6LBy/pm7gd0EDjDeze3+/4Ki8GnZXPal2lOH+q1mxt8rah/k5Rvdj7HexdX+maSLcJmaQ55WA/AHzjXO4uLfLJI9fIiTqLLXX++dwJjfI5ywNub7RIF+cAR61g7gXFe8QzzySZ+hQ9EReLjSi92sfcOBvnkEpfmy1LvY2MoQoFQgSu8t9PXJ4813n7sl/By5Z9y8CvYk493gD/twsNKH8yG63w2EGsKWOwVSKqUGauQBh2T19CqsJDzg4PIkBArxipRv9it70BcEYQEWw/At8RrmOizmUBj+NQrhDwJgkzrvoW4jp+f8j74wHtezViSHQjZ0AxoJuCP4W6fDbACBoaH8rsEUFK4FHG3zzoAXtodhpe//7F7lwHdvQq522cNADle4cSJX6nrO0oBsRtX81cvSJZIdonrvef69yFvO6xbgADtmkfhW/K9lwl9vHK52+cPOAzzW0RzREqPO+/PXm73WWYtU9eyufXpDjTHuodDvPbAEldTeKsWpf7dAYZ77eZKn0yyk4XPvZqXOx6SsItmmVkUeHlBi2iKX12w63G+dwLdtyzmgLc37ZqX6c4oAtkwC7Kzrf9zUZHWdWXeeyxZwFl+vvwaWWL8Tj5EHO+9l+N67FzLRMlnmbc/yRLmvreC9f+TXf+B/AK8AwNpF9bMWoEnt8T7f/kbXFRYyLzCIFY2C0GA1gCuxhTvras5r6iI3JBgvikMPvbv5rrHE302lH7vlXG3zywA3stsxr6gQCJKvHb3ew9406v8e08OJcEm179dRBj4W91mxXVv6cF7726f1WBgSoY9Y2OaVIuCdj00MBkHYPdy2PM7JK6CpDjILTN4KygKQttAcHMICIeAsBKPUPAJBB9/8AlwPfyPPbz9QLwZ/uLPFOJFEV4UGi+KEJ69vC8f/bGXP3al88ejo/D39QXxsh4uz/2whf/+Yn3TNa5PuDWTRxHo50NuQREBPt7uP+BjXvqZbcnHmnSHdGnOz9sOlVpm1JTa8sT6uVN0MN/dNaTUYEJjDMt3pjAgJgIfby+WbjvE/I0H+GD5Lo9u69YnzsPPx4uftyZz3Tt/uOvdd1rpPv/K+nIPpOUw8OmfyqWvm3IuoQGlm9ZLNvt+d9dZbDuYyT2fram0buP6tOGV/zu11HXF9fj3wm3884etAGyeNsbjAZaedIu8f+PpnN21ucf92RV1NRWfu+9INi2a+ePjfey9UtHrKZv+8pV9ufvTNVXWtar63TxzBQs2J/HaVf3oGBXEZyv2ePzeADivdyveuOY0j84tKCzi7s/WcOPgWE7rGMGLP2zhFVf3UfG/UclBl/eM7MJLC6yBkO/feDo3v7+SUT1b8tpVx8bKDHr6J/aX6FoKC/Tl1f87lbO7VrwE//Kdh7lyxvLj1rO4G+u5eXG8vrh861x1PTq2BzcP6USvyfPIyissNV4pK7cAH2+pdGDn9e/8wZKtyQDcO6orL/64lRnXnsbnq/by46aDAHRuEcL2pEy+mngmB9Nz3eOK/n3Vqby8YBvbkjJpFuDD+imjKyyjJBHRrofq0EChAUiKg03fQNwcOGB9o8MnANqcCi16QvNu1iMiBpq1tj7wa6B4hP7xePphsfOp8yscxAVWc+WqXSmM6tmK6XM28feRXekzzRqIFP/0+Uyatb7ciO7jlV2Rsh9yvzwwnOvf+YOdJfqcT2oezE/3DXP//s7SeKaVef3Bft7MvXsIHaOCKy0rt6CQbo/OK5VWUV2TMnJ48YetTDq/B2GBvuXGfdw/uhs9Wjdj9tr9fPVnIhf1bcNLV57K9+v3c8dHq4Fj97X49X1660AGdqp4UF5Fnp67mbdKTGd8/8bTOT020t2fH+Lvw/op5yIirN1zhH1HsunSshmdooMr/fcE2J6UQXpOAe//lsDl/dszuHP5Vp9ihzJz6f+ENe6j+AP0kz92M2nWegBuH3oSD53XndcXb3cPiHziot48+rX1zfPGwbG886sVlL5weR/u+3wtn906kDMquA9PzNnE20vjy6X/8cg5BPn54O8aszJ9zibeX2YFEOf1bsXO5Cy2HMzgzhGdue/cbpW+lqoU/zsF+nqXGhdx4+BYJl/Ys9wgUX8fL7Y8cZ7797yCIvILi1i+8zC7U44yoZKxByVVNjOpWPF705OxRMVO7RDOpPN6uAfYFpt+UW+uHdjRozwqkpKVx/Xv/ME7NwygebPSf7+OF9QWj7PZkJjGR7/v5omLepcaTF0ZDRSqSQOFeqowHzZ+BSvfgd3LAIH2p0O38yDmbGh1Mvj4VZlNVTJy8vH2EoL8rGbgkh9aF/Vtw7SLenM0t7Dct+XjfVhn5xVy+pML+PDmM+jTPrxa9flsxW5aNAtgePcWFBQWkZyZS9z+DHq1CaVFaEDVGZTx245DXPUfq7+45IfIjuRMznlhCaEBPqwr8w1kzZ4jXPTasalUd53ThXtLDFI8nqIiw4xfdrIyIZXWYQFMv6i3R9fl5Bfi5116gFhyRi7/+HwtT11yMm3DrQGZL/64lVd+2saKR0by8Ffr3d+2qjtqPSu3wD317oYzY5gyrpe7HkXGuN8PTqvsQ+CuEZ25t8QH85KtyYQF+tKnXRhPfx/HNWdYH0q3frCST28dSHjQ8f8v7DqcxdDnF5dLL3vfjDG8sWQHz83bwqyJZxLo680DX6zjyzvOLDUAtroe+GIt/1tZfkzL8knn0CosgLSj+e4gGao3a6cye1OPMuXbjSzYbM3YKBmklB0Uu+VAhnugasIzY+nyyFw6t2jGVxPP5Lr//sEfCSmc3DaM2XeeBcCPmw5yy/srOSM2kgfP686p7cNtneJY9nWc9eyiCo+VnWXlKQ0UqkkDhXqmqAg2fAGLn4aUnRDZCU6bAKf8BZq1tK2YlKw8Pl2x2/1Nbe5dQ+jZJpT/rdzDA19YrRa/PjTC/QG1cV8aKxNSeXvpTu4f3Z1xfdrYVhenLd95mI5RQeXWVdi4L43urULL/aHJLyyiS4m1Aor/mNcH8zbs5/YPV5dKm3Red24belId1ahmFm9J4oZ3V5RLd+Kely1r4X1D6dQ8xNYyKmOModuj88grPDayv+yHXMmgKf7p82354E3PyeeCV5Zy/+hunNe7FXd+8idjerdyT3ctKa+giMIi4+6+KjnjqK5VFFAeb7pwVewIFHQwo6obSXEw+27YsxxangxXfgJdx7gHDNrp4Vnrmbfx2Kj781/5hc3TxriDhLLz7Hu1CaNXmzCuPzPG9ro4rbIm+V5twipM9/X2IuGZsaQdzedgRk69CRKAUosfFatuq019cnaX0n3sfdqF8drV/Ry558XrJbQND+Sbvw2uch0JO4kIz112inssSkWBQLeWzdhyMINlk0bY9gEdGuDLzw8Md/9+vHEWZVtM6kuQUNIFp7Rmzrr99GoTyg2DY+q0LhooqNplDKx4G+Y/DH7BcNEbcMqVtgYIadlWF0PZlQxLKjkPujb/iNZXYUG+hAVVPMe/rsRGlx4j8eFNZ3BGrP0rHNYWLy/h+7uH8NjXG/jk1oH4etsfFBfz9/Hmw5vOoEvLkDp5f4/oYc3yGdq1eYUfwrPvPKvUN3p1zLx7hnA0r5BT24fz4hV9a9QNZBcNFFTtyc+xWhHWfQpdRsP41yCk4pHMJyqvoIg+U38gKtiPVY+NAmBrUgbBft68dW1/Pl2xmzklFnh57IKetpav7BPs78Nf+renR+tmlS6m09D0aB3KF3ecWStlndWl8sGVTgsN8D3uWJL68OFXX3VvFer+2c+nfrR0aKCgakduJnx6FcQvgeGPwJB/ONLNsD7xCACHs/LYnpSBj5cXO13TEs/qEs1ZXaKZs87qAzy1Qzg3ndU4PoAaqxPtl1VK2UfDOuW8vKPw4aWQ8IvV1TD0AduChKSMHHo/Pp9Vu1I4mldQatrV64t2uKc23VjiG+mdI6wV9x4d28OWOiilVGOmLQrKWYUF8OVN1qJJl70DvS+xNfvirXAr2rxn1p+J7p+vG3Rs3vNtQ0+iX4cI+nU48R0dlVKqqdBAQTnrh0dgy1w4/5+2Bwm7Dx+t+iSXkiPLQ/x9TngbXqWUamq060E5Z8Ms+P1NGPhXOP0W27P/eo3VYtCpzOj4GwfH8t4EawGXtuGBfPPXwdXavlUppdQx2qKgnJGyE769C9oNgFFTHSnixR+tdf8X/mMYAL9uP0RYoK97O2E79p1XSqmmTgMFZT9jrCBBvKxxCd72z88vLCq/oujx1tpXSil1YjRQUPZb85E1w+GClyC8gyNF/PMHaznmZy452ZH8lVJKWTRQUPY6mgLzH4EOg6Df9bZnb4xh/sYDvOHaKnZoN3sXbFJKKVWaBgrKXj//E3LTYeyLti+olFtQyHkv/+JeQAmg1QnssqiUUspzOutB2Sd1F6z4D/S9ClravzTy7ztTSgUJdu06p5RSqnIaKCj7LHrKGsA47GFHst96MMP983OXnaJBglJK1QLtelD2SImH9f+DgRMhrPz+73bYtD+dqGA/frx3KJHBfo6UoZRSqjRtUVD2WPZvEG8Y9FfHiti0L52T24VpkKCUUrVIAwVVc5nJ8OeH0OcvENrGkSIKiww7kjPp1qqZI/krpZSqmAYKquZWvA0FOXDm3Y5kn5GTz/2fryW/0NAuIsiRMpRSSlVMxyiomiksgNUz4aRzoHlX+7MvMpw85Qf37zodUimlape2KKia2TYfMvZD/wmOZL83tfQOkcN0gSWllKpV2qKgambluxDSCrqOcST7oc8vBuDuc7pw7aCO+HprbKuUUrVJ/+qqE3dkD2xfAP2udWTjp7yCIvfPd53ThegQf9vLUEopdXwaKKgTt+ELwEDfqx3J/mB6DgDn9myJt5curqSUUnVBAwV14tZ/Ce0GQGSsI9kv23kYgOvPjHEkf6WUUlXTQEGdmKQ4OLgeel/mWBEPfLEOgA6ROiVSKaXqigYK6sRs+MLa16HXxY5kn5NfCECArxftNVBQSqk6o4GCqj5jYMOXEDMEmrV0pIgXf9wKwLCuLRzJXymllGc0UFDVl7QZUnY61poAMOPnnQD8bURnx8pQSilVtVoPFESkvYh8ISJpIpIuIrNEpIOH15pKHn0drrYqactc67nbeY4X1bttmONlKKWUqlytLrgkIkHAQiAXuB4wwBPAIhE5xRiT5UE27wFvlUnbamc9VRW2fA9tT4NmrRzJ3hiDv48X/3e6R/GjUkopB9X2yoy3AJ2AbsaY7QAisg7YBtwGvOhBHonGmOXOVVEdV8ZBSFwJIx51rIj0nAJyC4poGx7oWBlKKaU8U9tdD+OA5cVBAoAxJh74FRhfy3VRJ2Lr99Zzt/MdK+JAmrXQUqsw3QBKKaXqWm0HCr2ADRWkbwR6epjHHSKSKyJHRWShiAyxr3qqSlvmQXgHaOHpP1f1bdyXBkBsdLBjZSillPJMbQcKkUBqBekpQIQH138ITARGArcCUcBCERlW2QUicquIrBSRlcnJydWusCqhIA/if4bOo0CcW1L53v+tBaBn61DHylBKKeWZutg90lSQ5tGnjjHm2hK//iIi32C1UDwBnFXJNTOAGQD9+/evqGzlqcSVkJ8FJw13rIgFmw66f/bS/R2UUqrO1XaLQipWq0JZEVTc0nBcxpgM4DtgQA3rpTyxY5G1GmOMM709C+MOcvP7KwFrIyillFJ1r7YDhY1Y4xTK6glsOsE8hYpbKZTddi6GNv0gMNz2rNOy87nxvZXu35+99BTby1BKKVV9tR0ofAsMFJFOxQkiEgMMdh2rFhEJBcYCv9tVQVWJnDRIXOVYt8Pvrp0iAZZPOoeIYD9HylFKKVU9tR0o/AdIAL4RkfEiMg74BthDiUWURKSjiBSIyOQSaf8Qkf+IyFUiMkxErseaVtkKcG5Sv7IkLAVTCJ2cCRSemrsZgDWTR+m0SKWUqkdqdTCjMSZLREYA/wI+wOo2+Am4xxiTWeJUAbwpHchsAS52PcKAdKxA4SZjzB+1UP2mbedi8A2GdvYPB4k/lEXC4aMAhAdpS4JSStUntT7rwRizG7i0inMSKDMTwhgzG5jtXM3UcSX8Ch3OAB/7P8jX7T0CoEs2K6VUPaS7R6qqZadC0ibocKYj2X+5OhGAaeMrGueqlFKqLmmgoKq25w/AQIeBjmS/5UA6vt6Cr7e+HZVSqr7Rv8yqart+Ay9fa8dIm6Vl53MwPZe/Du9se95KKaVqTgMFVbXdy6FNX/ALsj3rbQczAOjaspnteSullKo5DRTU8eXnwL7V0GGQI9lf819rCYwWzfwdyV8ppVTNaKCgjm/faijMcyxQKN746eR2YY7kr5RSqmY0UFDHt+s369mhgYwRQX50b9UMfx9vR/JXSilVMxooqOPb8ztEd4OgivbyqrkdyZmc1DzEkbyVUkrVnAYKqnLGWPs7tHdmc87cgkJ2pxzlpObBjuSvlFKq5jRQUJU7sguOHnZkWiTA2j1pFBnopC0KSilVb2mgoCqXuMp6btPPkeyveGsZAK11EyillKq3NFBQlUtcDd7+0NLZpZX7xzgz/kEppVTN1fqmUKoBSVwFrfuAt68j2Q/r1pzDmXl4e0nVJyullKoT2qKgKlZYAPvWODY+AeBAWo4utKSUUvWcBgqqYslxUJDtWKBQWGTYeSiLk1roQEallKrPNFBQFSseyNjWmYGMianZ5BUU6dRIpZSq5zRQUBVLXAUB4RDZyZHsdxzKBHRqpFJK1XcaKKiKJa62uh3EmYGGuw5lAdAxyv4dKZVSStlHAwVVXn4OJG2CNqc6VkTC4aME+XnTPEQHMyqlVH2mgYIqL2kTmEJofYoj2RtjeO+3BLq0CEEcarFQSillDw0UVHkH1lvPrU52JPvVu48AsDc125H8lVJK2UcDBVXegXXgHwrhMc5kn5YDwDOXOtNioZRSyj4aKKjyDqyHlr3By5m3x74jVkvC6bG6dLNSStV3Giio0oqK4MAGx8YnACQeyaaZvw9hgc4sDa2UUso+Giio0lJ2Qn6WY+MTAPakHKVNeKBj+SullLKPBgqqtAPrrGeHAoXlOw/zU1ySrp+glFINhAYKqrQD68DLF5r3cCT7K2csByAlK8+R/JVSStlLAwVV2oH10Lw7+Pg5WsxfR3R2NH+llFL20EBBlbZ/nWMDGQuLDACndYxgeLcWjpShlFLKXlUGCiLiKyLjRSS2Niqk6lDGQchKcmx8woLNBwHo3SbUkfyVUkrZr8pAwRiTD/wPiHG8NqpuuQcyOtOikJyRC8BNZzmzI6VSSin7edr1sBPQtuLG7uAG67llL0ey33YwgxB/H9pH6tRIpZRqKDwNFJ4DHhGR5k5WRtWx5C3QrA0EhtuetTGGmct2kZlboBtBKaVUA+Lj4XkjgEggXkSWA/sBU+K4McZcb3flVC1LjoPm3RzJevqczY7kq5RSylmeBgpnAflAMnCS61GSKXeFaliKiqwWhX7OxHu+3lYrwrx7hjiSv1JKKWd4FCgYY3TGQ2OXtgfyj0KL7o5kX2QMgb7edG+lMx6UUqoh0XUUlCV5i/Xc3JlAISUrn4gg3QRKKaUaGo8DBREJFpG7ROQLEVkkIl1c6VeKiMefLiLS3pVHmoiki8gsEelQ3YqLyCQRMSKytLrXqgokx1nP0V0dyf63HYcI1d0ilVKqwfGo60FE2gOLgXZAHNAbaOY6PBwYCdzsQT5BwEIgF7gea2zDE8AiETnFGJPlYX06AY8ASZ6crzyQvAVCWkJQpO1ZH80rYH9aDvvTcmzPWymllLM8Hcz4AtaHexdgH1ByR58lwBQP87kF6AR0M8ZsBxCRdcA24DbgRQ/zeQP4COiG569BHY+DMx7iD1nx3/WDOjqSv1JKKed42vUwCnjcGLOb8jMcEoG2HuYzDlheHCQAGGPigV+B8Z5kICJXAf2ASR6WqapijNWi4ND4hB3JVqBw5enV7mFSSilVxzwNFPyAjEqOhWFNnfREL2BDBekbgZ5VXSwiEcC/gAeMMSkelqmqkr4P8jIca1G465M/AYiJCnYkf6WUUs7xNFBYB1xaybHzgFUe5hMJpFaQngJEeHD988BW4D0Py0NEbhWRlSKyMjk52dPLmpbigYwOtCjkFhS6fw7087Y9f6WUUs7ytH//eeAL19K7H7vSeorIeOAmrC4FT1W0OFOVa/qKyBDgOqCfMcbjBZ6MMTOAGQD9+/fXhaEq4mCgsCPJ6na4TscnKKVUg+TpgkuzRGQi8Axwoyv5fazuiL8ZY+Z5WF4qVqtCWRFU3NJQ0lvAf4G9IhLuSvMBvF2/Zxtjcj2shyopOQ6CoiE42vasN+9PBzRQUEqphsrjGQPGmDdF5ANgENZOkoeB34wxlY1dqMhGrHEKZfUENlVxbQ/X4/YKjqUCfwdeqkZdVDEHBzLGHUjHz8dLxycopVQD5ek6CtcB3xljDgMLyhyLBC4wxrzvQVbfAv8UkU7GmJ2u62OAwcBDVVw7vIK0lwBv4E5gewXHVVWMsVoUel/mSPZxBzLo2jIEH29dBFQppRoiT/96v0v5jaCKxbqOe+I/QALwjYiMF5FxwDfAHqyuBQBEpKOIFIjI5OI0Y8zisg/gCJDm+n2vh3VQJWUehJw0R1oUCgqLWLc3jZ6tdX8HpZRqqDwNFI432DAYKPAkE9fKiyOwZi58gLVoUjwwwhiTWaY872rUT50o90BG+6dGrkhIJS07n+HdWtiet1JKqdpRadeDiPTFWtio2IUi0rvMaYHAlVgrK3rEtWhTZVMti89JwIOZEMaYYZ6Wqyrh4GZQG/elAXBGpyjb81ZKKVU7jjdGYTzwuOtng7W3QkUOY02RVA1RchwEhEOI/d/64w9lER7kS2Swn+15K6WUqh3HCxRewlrYSICdwCXAn2XOyQUOVmddA1XPJMVZrQlSZQNOtcUfyiI2Wmc7KKVUQ1ZpoGCMSQPSAEQkFthvjMmr7HzVABkDyZuhR3XWy/Jc/KEsBmm3g1JKNWgeDRY0xuwC8kVknIj8U0TeFZGOACIyVETaOFpL5YysQ5CdCi162J51SlYe+9NytEVBKaUaOE/XUYgA5gJnAOlAM+BVYBfW1tEpwF0O1VE5xcEZD1sPWutw9WqrUyOVUqoh83T64fNAe6yFkaIpPSNhAXCOzfVStcHBPR72pmYDEBsdYnveSimlao+nSziPB/5hjFkmImW3ANyNFUSohiZ5C/iHQrPWtmed6AoU2oQH2J63Ukqp2uNpi0IIkFjJsQA8WPNA1UPJcVa3gwMzHvamHqVFM3/8fXRraaWUasg8DRS2AOdWcmwosN6e6qhalbzFkfEJAHtSj9IhMsiRvJVSStUeT7seXgNeE5E04GNXWriITAD+BtzqROWUg46mQFaSY7tG7knJ5vTYinYUV0op1ZB4FCgYY/4jIicBU4FpruQfgSLgOWPMRw7VTznFwYGMeQVF7EvLpr22KCilVIPnaYsCxpiHROQNrC6I5lhLN/9YvF20amAcnBq5aEsSxkC78EDb81ZKKVW7PA4UwL3w0n8cqouqTclbwDcYQtvZnvVtH6wCIMBPBzIqpVRDV61AQUTaY02FLDfnzRiz0K5KqVpQPOPBy7mdvM/v3cqxvJVSStUOT1dm7AR8BJxenOR6Nq6fDaBfHxuS5C3Qabjt2RpjCA3w4cI+bfDxdi4IUUopVTs8bVF4G+gA3APEAbo5VEOWfQQy9jsyPuHI0XzScwp0jwellGokPA0UBgA3GGO+dLIyqpYc2mo9OzDjIfGItSJjuwgdyKiUUo2Bp23De9FWhMbDwRkP+44UL92sgYJSSjUGngYKTwEPioi2JzcGyVvAJxDCO9ie9f60HEADBaWUaiw8XXDpAxHpDiSIyHIgtfwp5nrba6eckbQZoruAl/3jT/ekHMXPx4uoYD/b81ZKKVX7PJ31cAMwCSgE+lG+G8LYWy3lqOQt0PFMR7LetD+dHq2aIQ5sNKWUUqr2eTqYcSrwFXCTMeaIc9VRjstJh/S9jm0GtTc1m77twx3JWymlVO3zdIxCFPC6BgmNwKFt1nOLHrZnXVRk2J+WreMTlFKqEfE0UFgK2P/Jomqfg5tBHcrMJb/Q0Da83MKdSimlGihPux7uBv4nIqnAPMoPZsQYU2RnxZRDkuPA2x/CO9qe9YF0a8ZDy1ANFJRSqrHwNFDY7Hp+v5Ljphp5qbqUvMWa8eBt/z9XSpY1xjUqxN/2vJVSStUNTz8tpqEzGxqH5Dho19+RrFOPWoFCpE6NVEqpRsPTdRSmOFwPVRvysuDILjj1GkeyP5zpChSCNFBQSqnGQrf3a0rcezw4MzUyKSMXPx8vQgO1F0oppRoLj/+ii4gfcB7QDSg7Ws0YY6bbWTHlgOQt1rMDMx4AZvy8E0AXW1JKqUbE05UZ22BNkYzBGqtQ/ElQctyCBgr1XXIcePlCZKe6rolSSqkGwtOuh+eBZKADVpBwBtAJeBLY7vpZ1XfJWyCqM3j72p51dl4hAPeN6mp73koppeqOp4HCEOAFYJ/r9yJjTIIxZjLwBfCKE5VTNkuOc2x8wpo9RwBoFaZrKCilVGNSnSWc97kWVcoCIkocWwgMs7leym752ZCa4Nj4hAPp2QC6z4NSSjUyngYKe4Fo1887gHNLHDsdyLGzUsoBh7eDKXKsRWHW6kQA2kcGOZK/UkqpuuHprIdFwFDga+At4DUR6QvkA6Ndaao+c3jGQ/yhLEIDfAjw9XYkf6WUUnXD00DhUSASwBjzhoj4AH8BgoDnsFZuVPVZ0mYQb4g6yfasjTGkZOVx5YAOtuetlFKqbnna9ZAP7Cr+xRjzqjHmLGNMP2PMw8YYj7seRKS9iHwhImkiki4is0Skyk8YEekoIt+IyC4RyRaRQyKyWETO87TsJi05zpoW6WP/PgzpOQUczSuktQ5kVEqpRqfKQMHVenCY0uMSToiIBGENfuwOXA9cC3QBFolIcBWXhwCHsFo3zgduAjKBuSJySU3r1uglb3FsfMK+I9ZAxta6vbRSSjU6VXY9GGMKROQgUGhDebdgrbnQzRizHUBE1gHbgNuAF49Tj41YwYGbiHwHxAMTgFk21K9xKsiFlJ3Qc7wj2b+5ZAcAwX66dLNSSjU2nnY9fAjcbEN544DlxUECgDEmHvgVqPanmDGmAEjD6hpRlTm8A0whtOjhSPaHMnMBOKNTpCP5K6WUqjuefgVMAK4SkRXAN8B+ymw7bYx5x4N8ermuL2sjcLknFRERL6wAJxqrhaIrcLcn1zZZyXHWs0NdD4G+Pvj5eBGkLQpKKdXoePqX/TXXc1vgtAqOG8CTQCESSK0gPYXSizgdz3PAfa6fM4ErjTE/eXht05S8BcTLWr7ZAXtSjnJ2l+aO5K2UUqpueRooxNpYpqkgrTrbDb4EfAq0Aq4DPhaRy4wxcyo6WURuBW4F6NChiU7fS46DiBjwDbQ9a2MMu1OOMrhzdNUnK6WUanA8ChSMMbuqPssjqbjWYygjgopbGiqqy16slSIB5ojIYuCfQIWBgjFmBjADoH///hUFKY1f8hbHFlo6nJVHdn4h7SPtD0KUUkrVPU8HM9plI9Y4hbJ6AptOMM+VgDNt6o1BYb61fLND4xN2pxwFoIMu3ayUUo2Sx4GCiIwWka9EZJOI7Cz78DCbb4GBIuLellpEYoDBrmPV4hrYeBbW/hOqIod3QFG+Yy0KxWsotI3QFgWllGqMPAoUROR8YC7Wks3dgThgN9AeKAKWeFjef7BmUHwjIuNFZBzWLIg9lNgvwrUKY4GITC6RNkVEXhGRv4jIUBH5CzAPa1Oqxz0sv+lxz3hwJlA4lGFNjWweYv+Kj0oppeqepy0Kj2HNfDjf9fujxphhWN0I3sD3nmRijMkCRgBbgQ+Aj7AWTBphjMkscaq48i1Zv9VAb+BV4Aes2Q85wBBjzKcevo6mJzkOEIju6kz2mbl4ewkRQX6O5K+UUqpueTrroTswGav1wBRfZ4zZKiJTsAKJ/3mSkTFmN3BpFeckUGYmhDHmW06ge6LJS9pszXjwc2YMwaGMPKKC/fDyqs7EFaWUUg2Fpy0KRUCBMcYAyUDJeYb7APu3JFT2SI5zbEVGsFoUorXbQSmlGi1PA4UtQIzr55XAPSLSWkSaYy1+lGB/1VSNFeQ5OuMBrOWbmzfTQEEppRorT7sePgKKv5Y+Dizg2FoGhcBVNtdL2SFlBxQVQHPnWhQOZeTSpUUzx/JXSilVtzxdcOm1Ej+vEpGTgfOAQGCBMeZE10BQTkrabD23cGbGQ1GR4WCGtigopVRjdkK7+LhWR/yPzXVRdkuOs/Z4cGjGw57UoxQWGV1sSSmlGrFqBQoiMhwYhLU5VCLwmzFmsQP1UnYonvHgwB4PAIcy8wBoHR7gSP5KKaXqnkeBgohEAp8Dw7CmR6Zi7c8grr0WLjfGpDhUR3Wikrc4Oj4h/lAWoIstKaVUY+bprIdXgAHAtUCgMaY51viE64D+wMvOVE+dsII8azCjQ+MTADbuS8PHS+jROtSxMpRSStUtT7seLgQmGWM+Lk4wxuQDH7laG55wonKqBg5vd814cC5Q2Lw/nV5tw/DWxZaUUqrR8rRFoRDYVsmxLa7jqj5Jds14cDBQSDySTUyUDmRUSqnGzNNA4RvgL5UcuxL42pbaKPskOTvjoajIcDAtl1ZhOpBRKaUaM0+7HmYD/xKR77AGNR4EWgJXYG0MdbeIjCg+2Riz0O6KqmpK3gwRseDrzAf54aw88gqLaBOm20srpVRj5mmg8IXruT3WQktlfel6FqxZEd41rJeqqSRn93g4kJYDoC0KSinVyHkaKAx3tBbKXgW5kLITeo53rIg1e1IBtEVBKaUaOU+XcF7idEWUjQ5tA1PoaIvCY99sBLRFQSmlGjtPBzOqhiQ5znp2aNfI3IJjk1x0nwellGrcNFBojA5uBC8fx2Y8bDuYCUCzgBPaKkQppVQDooFCY3RwoxUk+DjzbX/rwQwAZt1xpiP5K6WUqj80UGiMDm6Elr0cy37rwUx8vYWY6GDHylBKKVU/eBQoiEiYiGhndEOQnQrpex0NFHYdzqJ9ZBC+3hpnKqVUY1flX3oR8QEOA+c6Xx1VY0mupZtb9nasiMQj2bQN12mRSinVFFQZKBhjCrBWYtT9HBqCg9a0RVr0dKyIxNRs2kVooKCUUk2Bp23HHwI3O1kRZZODGyAgHELbOJJ9SlYeh7PyiInS8QlKKdUUeDq/LQG4SkRWYG0QtR9rqWY3Y8w79lZNnZCDG61uB3Fm6+eN+9IA6N02zJH8lVJK1S+eBgqvuZ7bAqdVcNwAGijUtaIiOLgJTr3GsSI27ksHoFebUMfKUEopVX94GijEOloLZY8jCZCf5eiMh4370mkbHkh4kJ9jZSillKo/PN3rYZfTFVE2KB7I6OCMh/hDmXRuEeJY/koppeqXak2EF5FTRORvIvK4iLRypXUWkWbOVE9Vy8GNgECL7o5kb4xh16GjxEQFOZK/Ukqp+sejFgXXYksfApcAgjUmYTZwAHgO2Ao85FAdlacOboDITuDnzIyEQ5l5ZOQW0FFnPCilVJPhaYvCk8BI4FqgJVawUOx7YLTN9VIn4uAmR8cnbHDNeOipAxmVUqrJ8DRQ+D/gUWPMx0BKmWPxQIydlVInIC8LUnY6uyJjajYAsbrHg1JKNRmeBgpRwObj5KH7QNS1g5sA42iLwsH0HLy9hOgQ/edWSqmmwtNAIR4YVMmx04Et9lRHnbD9a6zn1n2cKyIth+Yh/nh7ObOYk1JKqfrH00DhfeAhEbkaKJ5Ab0RkOPB3dLGlurd/LQRGQlg7x4qYtXovEcG6foJSSjUlngYKzwHfAR9wbIzCUmABMM8Y86oDdVPVsX+t1Zrg0NLN2XmFFBnYciDdkfyVUkrVT54uuFQIXCkir2HNcGiBtfX0PGPMEgfrpzxRkGdtLz1oomNFJGXkADDpvB6OlaGUUqr+8XQJZwCMMb8AvzhUF3WikjdDUb6j4xOSMnIB6NpK19ZSSqmmpForM6p6av9a67l1X8eKmLV6LwBROkZBKaWalEoDBREpEpFCTx+eFigi7UXkCxFJE5F0EZklIh08uK6/iMwQkTgROSoiu0XkIxHRDav2rwO/ZhDh3K04nJkHQJeWus+DUko1JcfrepiGtVQzWCsx3ggEYi3dfBBoBVwAZAP/9aQwEQkCFgK5wPWu/J8AFonIKcaYrONcfiXQC3gF2Ii15fVjwEoR6WuM2eNJHRql/Wuh9Sng5UwDkTGGhMNZ9G4bir+PtyNlKKWUqp8qDRSMMVOKfxaRR4FdwGhjzNES6cHAfKDAw/JuAToB3Ywx2115rAO2AbcBLx7n2meNMcklE0TkV6w1Hm4BJntYh8alqNDa46Hf9Y4VkZ5dwNaDmTx0njObTSmllKq/PP0KehvwfMkgAcDVAvBP4HYP8xkHLC8OElx5xAO/AuOPd2HZIMGVtgtIxmpdaJoOb4f8o44OZEw8Yi3d3CFSd41USqmmxtNAIZpjCy2V5Ye1xLMnegEbKkjfCPT0MA83EemBNVWzsuWlGz/3QEbnA4U24YGOlaGUUqp+8jRQWAlMFZFS39xdv08BVniYTySQWkF6ChDhYR7FZfsAb2K1KFQ6RkJEbhWRlSKyMjm5XKNEw7d/LfgEQHRXx4r4dfshANpqoKCUUk2Op+so3IU1CHGHiCzHGszYEhgIHAWuqkaZpoK0E1lO8N/AmcBYY0xFwYdVmDEzgBkA/fv3r6jshi1xtbVjpHe1lsSolvd+SwAgOkSnRiqlVFPjUYuCMeZPoDPwAlAInOx6/ifQxRizxsPyUrFaFcqKoOKWhgqJyNPArcCNxpgfPL2u0SkssDaDatff8aJio4MRh5aHVkopVX9V+TVURPyAO4CfjDGP1LC8jVjjFMrqCWzyJAMReQR4CLjLGPNBDevTsCXHWQMZ2zoXKBQUFuHtJYw9ubVjZSillKq/qmxRMMbkAc9QcUtAdX0LDBSRTsUJIhIDDHYdOy4RuQtr3YVHdCMqIHGV9dy2n2NFHMzIpbDI6EBGpZRqojwdzLgZa/2DmvoPkAB8IyLjRWQc8A2wB3ir+CQR6SgiBSIyuUTalcBLwDxgoYgMLPGo9oyJRiFxJQRGQKQd/zQV2+ea8dA2QgMFpZRqijwdATcZeFlEVhlj1p9oYcaYLBEZAfwLa8tqAX4C7jHGZJY4VQBvSgcyY1zpY1yPkpYAw060Xg1W4mpoe5pjW0sDJKa6AgVtUVBKqSbJ00DhQSAE+FNEEoD9lJ69YIwxQz3JyBizG7i0inMSKDMTwhhzA3CDh/Vt/PKyIGkTdB/raDHFayhooKCUUk2Tp4FCIR4ONlS1ZN8aMEVWi4KD9qZmExXsR6Cf7vGglFJNkUeBgjFmmMP1UNXlHsjobKCQeCRbBzIqpVQT5sx2g8p5iasgvCMERztaTMKhLO12UEqpJszjQEFEWovIP0VkhYjsEJE/ROQ5EWnlZAVVJRJX1Uprwu6Uo/RqE+poOUoppeovjwIFEekKrMFayjkT+APIAu4G1ohIF6cqqCqQcQDS9ji+IuOsVXsBaK0tCkop1WR5OpjxWSAdOMM1IwGw1jsAfnAdv8T22qmK7V5mPXcY6FgRxhh+23EYgHN7tXSsHKWUUvWbp10Pw4HHSgYJAMaYXVi7Rw63t1rquHb/Dr5B0OoUx4r4+I/dLNtpBQqhAb6OlaOUUqp+8zRQ8AMyKjmW4TquasvuZdb4BG9nPsA37kvjhR+2AhAWqEGCUko1ZZ4GCmuAO0Wk1PlibSc40XVc1YbcDDiwDjoMcqyIsa8sJSUrD4DVj41yrByllFL1n6djFKYBc4DNIvIZ1sqMrYDLgS6As8sDqmP2rrQWWnJofMK/F24r9bu3l24trZRSTZmnCy7NE5ELcO3ciLW8sgFWARcYY35wroqqlN3LQbyg3QBHsv+nq8sBIOEZjf+UUqqp87RFAWPMPGCeiAQBEUCqMeaoYzVTFduzHFr2ggD71zYw5tj2HW9c7dzW1UoppRqOaq/MaIw5aoxJ1CChDhQWwJ4Vjo1P+HD5LvfP553c2pEylFJKNSy6hHNDcnA95Gc5Nj7hC9cCS3eO6OxI/koppRoeDRQakl2/Wc8OtCgUFRnW7k0D4N5RXW3PXymlVMOkgUJDEv8zRHWG0Da2Zz1zWYL7Z2vWq1JKKaWBQsNRWAAJv0Ls2Y5kv3r3EQA+uOl0R/JXSinVMGmg0FDsXwN5GY4FCrtTrLGpZ3V2dttqpZRSDYsGCg1F/BLrOWaI7VnnFxaxds8Rhndrrt0OSimlStFAoaGI/xla9IJg+7/xP/b1BgCSM3Ntz1sppVTDpoFCQ1CQa63I6FC3Q0GRtdDSG1ef5kj+SimlGi4NFBqCvSugIMexQOH3+MNEBvvRPjLIkfyVUko1XB4v4azq0M4l1v4OHc+0Peu4A+nsScm2PV+llFKNg7YoNAQ7foI2/SAw3Pasz3/5FwBuGRJre95KKaUaPg0U6rusQ5C4GrqMsj3rgsIiXMMTmHReD9vzV0op1fBp10N9t2MhYKCz/YFCwmFr7YR2EYF4eem0SKU8lZ6eTlJSEvn5+XVdFdVE+fr60qJFC0JD7d9JuCwNFOq7bT9CUBS0OdX2rLccyADgzWt0toNSnkpPT+fgwYO0bduWwMBAXXtE1TpjDNnZ2SQmJgI4Hixo10N9VlRkjU/oPBK87P+n2nIgHS+Bzi1CbM9bqcYqKSmJtm3bEhQUpEGCqhMiQlBQEG3btiUpKcnx8jRQqM/2/QlHDzvS7VBQWMSGfenERAcT4Otte/5KNVb5+fkEBgbWdTWUIjAwsFa6vzRQqM+2/wgInDTC9qxvfn8lC+OS6KhrJyhVbdqSoOqD2nofaqBQn22dD21Pg+AoW7P9afNBFm9JBqBbK+cHwiillGq4NFCor9ISYd9q6H6+7Vm/9fNOAPx8vHhwTDfb81dKKdV4aKBQX8V9Zz13v9D2rJv5W5Nd1k85V5tQlVK89957iAjh4eGkpqaWOlZQUICIMGXKlBqVMXv2bK666iq6du2Kl5cXw4YNq/TcpUuXcuaZZxIYGEirVq249957yc4uv4Lsxo0bOffccwkJCSEqKooJEyaQkpJS7rw9e/Zw2WWXERYWRmhoKJdccgm7d+8ud15qaio333wz0dHRBAcHM3LkSNavX1+j190YaKBQX8XNhuhu0Lyr7VnHH85iTK9W+PvoIEal1DFpaWk8++yzjuT99ddfs2bNGgYOHEi7du0qPW/dunWMGjWKFi1aMGfOHJ544gneffddbrjhhlLn7du3j2HDhpGdnc0XX3zBa6+9xoIFC7jgggsoKipyn3f06FFGjBhBXFwcM2fO5IMPPmDbtm0MHz6crKws93nGGMaNG8e8efN49dVX+fLLL8nPz2f48OHs3bvX9vvRoBhjmszjtNNOMw1C1mFjpkQYs2Cq7VkXFRWZro/MNU/M2Wh73ko1BZs2barrKtju3XffNYA599xzTVBQkNm/f7/7WH5+vgHM448/XqMyCgsL3T8PHjzYDB06tMLzLrroItO5c2eTl5fnTps5c6YBzKpVq9xp99xzjwkLCzOpqanutCVLlhjAfPnll+60l156yXh5eZlt27a503bu3Gm8vb3NCy+84E77+uuvDWAWLlzoTjty5IiJiIgwd9555wm95tpQ1fsRWGlq+NmpLQr10dZ5YAqh+wW2Z715fwa5BUX4+eg/vVKqtEcffRSAJ5980va8vTxYCyY/P5958+ZxxRVX4Ovr606/4oor8PPz45tvvnGnffvtt4wdO5bw8HB32tlnn02HDh3KnTdw4EA6d+7sTouNjWXw4MHlzmvTpg3Dhw93p4WFhXHhhReWOq8p0k+L+mjzHAht58hqjGv3HgGgf0yk7XkrpRq21q1b87e//Y0ZM2awa9euSs8rKCio8lFYWFjt8nfs2EFOTg69e/culR4QEMBJJ53Epk2bAMjOziY+Pr7ceQC9evVynwfWOIaanrd7924yMzOr/Xoai1pfwllE2gP/AkYBAiwA7jHGlB9ZUv7ap4D+wGlAJDDBGPOec7WtAznpsH0B9J8ANg80LCwyTJplDcwZ1MneKZdKNXVTZ29k0770Oq1DzzahPH5hrxrl8eCDD/LWW28xdepU3nnnnXLHExISiI2terfZjh07kpCQUK2yiwciRkRElDsWGRnpPp6amooxptLztmzZUirPys4rOXAzJSWFmJiYCs8rLjMkpGmuYlurgYKIBAELgVzgesAATwCLROQUY0zW8a4H7gTWAHOA6xysat2JmwOFuXDy5bZnPW/DAffPuhqjUqoikZGR3HfffUydOpUHH3yQk046qdTxNm3asGLFiirz8ff3r3bZVpd6xQsJFR+rznnFPDnPGONxfk1Nbbco3AJ0AroZY7YDiMg6YBtwG/BiFdeHGWOKRKQzjTVQWP85RMRYCy3Z7Jl5mwEY0iXa9ryVaupq+k2+Pvn73//Oq6++yuTJk/noo49KHfPz86Nv375V5nEiU6+Lv71XNMUxNTWVXr2sexwREYGIVHpecT7F51Z2XsmWhpItFmXPK86nqartMQrjgOXFQQKAMSYe+BUYX9XFxpiiqs5p0DIOws7FVmuCA+sbhPhbg4M+uOkM2/NWSjUeISEhTJo0ic8//5w1a9aUOpaQkICvr2+Vj7ItEZ446aST8Pf3Z+PGjaXSc3Jy2LlzJz179gQgKCiImJiYcucBbNq0yX0eWGMManpehw4dmmy3A9R+oNAL2FBB+kagZwXpTcvGr8AU2d7tUFhkOJSZy+b96Vx1Rgdb81ZKNU4TJ06kbdu27pkQxYq7Hqp6zJ49u9pl+vn5MWbMGP73v/9RUFDgTv/iiy/Izc1l3Lhx7rRx48bx3XffkZaW5k5bunQpu3btKnfe8uXL2blzpzstISGBX3/9tdx5iYmJLFmyxJ2Wnp7O7NmzS53XFNV210MkkFpBegrgSLuOiNwK3ArQoUM9/5Bc/zm0Ohma27us8tTZG3l/mTWCOTuv+iORlVJNj7+/P5MnT+bWW28tle7n50f//v2rnd+uXbvcYxsOHz6Ml5cXX3zxBQADBgygY8eOAEyZMoVBgwZxxRVX8Ne//pWEhATuv/9+LrvsMk477ViX7P3338+HH37IuHHjmDRpEmlpaTzwwAOcfvrpXHzxxe7zbrnlFv79738zfvx4nnjiCUSExx57jPbt23Pbbbe5zxs3bhyDBg3immuu4fnnnyciIoKnn34aYwwPPPBAtV9vY1IX0yMrGhni2DrCxpgZxpj+xpj+zZs3d6qYmkveCokrHRnEWBwkANqioJTy2IQJE+jSpYsteS1atIjLL7+cyy+/nLi4ODZt2uT+fdGiRe7z+vbty/z589m/fz9jx47l4Ycf5rrrrmPmzJml8mvbti2LFi3Cz8+PSy+9lNtvv53hw4czd+7cUms2BAcHs3DhQrp27cq1117L1VdfTWxsLAsXLizVneDl5cWcOXMYNWoUEydO5OKLL8bb25tFixbRvn17W+5BQyW1OaJTRA4CXxtjbiuT/jpwuTHGo09y12DGbVRzemT//v3NypUrq1HjWvTDo7D8Dbh3M4S0sC3boiJDp4fnAvDK/53KuD5tbMtbqaZo8+bN9OjRo66roRRQ9ftRRFYZY6rfBFRCbbcobMQap1BWT2BTBelNQ0EerPkEup1na5AAsCPZWiTkuctO0SBBKaVUtdV2oPAtMFBEOhUniEgMMNh1rGnaMheOHoJ+19uabWpWHqP+9TMAp3VsulN7lFJKnbjaHsz4H+BvwDci8ijWeIXpwB7greKTRKQjsAOYZoyZViJ9KNAcaOVK6i8imQDGmC9q5RU4YfVMCGsPJ42wJbs/4lO44q1lRAb7udM6RQfbkrdSSqmmpVYDBWNMloiMwFrC+QOsQYw/YS3hXHIhbQG8Kd/iMRUYWuL3v7oexdc0PKkJsGMRDHsIvGq+WmJKVh5XvLXM/TPAzBtPP6HFT5RSSqla3+vBtafDpVWck0AFH/zGmGHO1KoO/T7DChD62bPQ5PQ55Yd6nK0rMSqllDpBuntkXcpJh9XvQ6+LIdSegYZf/ZkIlF6mWVsTlFJKnahab1FQJaz5CPIyYOAdtmSXeCTb/fO7NwzgUGYezZtVf2MWpZRSqpgGCnWlqNBaN6H9QNs2gIrbb21xO/mCnvh4e9EqLMCWfJVSSjVd2vVQV+LmwJFdMPB227L8cvVeAC7r3862PJVSSjVtGijUBWNgyfMQ1Rl62LfZyNz1BwAIDfC1LU+llFJNmwYKdWHL93BwPQz5hy1TIgFiHvoOgLBADRKUUifmhx9+4LzzziMqKoqAgAC6du3Kgw8+SGrqsb38jhw5wpQpU1i9enW561966SVmzZpV7XKXLl2KiNCyZctSu0aeiClTprBw4cIa5aFK00ChthkDS56FiBjbNoBKy853//zWtfaMd1BKNS1PPfUUo0ePJiAggLfffpv58+dz++2389577zFgwAD27NkDWIHC1KlTbQ0Uijd8SkpK4vvvv6/R65g6daoGCjbTQKG2bfsB9q+xWhO87RlLejA9B4Dbh57EwE5RtuSplGo6Fi1axKOPPso999zDV199xcUXX8zQoUO59957Wb58OSkpKVx3nT1rvZSVnZ3N559/zrBhwwgKCiq3S6Sqexoo1KaiQlgwxWpN6HOlbdn+uv0QAOf0sHdDKaVU0/Dcc88RGRnJ008/Xe5YbGwsDz30EIsXL+b3338nNjYWgFtuuQURQUR47733iImJYdeuXXz00Ufu9BtuuKHKsr/++mvS0tLcWzvPmTOnVFcHQEJCAiLCW2+9xeTJk2ndujXh4eFceOGF7N27131e8ZoxTz75pLsOU6ZMcR//8MMP6dOnDwEBAURHR3Pttdeyf//+UmXFxMRwzTXX8Omnn9KjRw+Cg4Pp378/S5cuLXXeihUrGDVqFFFRUQQFBdGpUycmTpxY6pz4+Hiuvvpqmjdvjr+/P3379uWrr74qdc6UKVMQEbZt28bYsWMJCQmhY8eOTJs2jaKioirvX23QQKE2rfkYkjbByCngbd9YgqmzrdUY+7QLty1PpVTTUFBQwJIlSxg1ahQBARVPqR43zhp0PW/ePHfXwqRJk1i2bBnLli1j7NixfPXVV7Rq1YrRo0e70x977LEqy585cybh4eGMGzeO6667jtzcXD799NMKz3366afZvn0777zzDi+//DLLli3j6quvdh9ftsxavv6GG25w1+Hmm28GYMaMGVx77bX06NGDWbNm8cwzzzB//nyGDh1KZmZmqXJ++eUXXnjhBaZPn85nn31GYWEhF1xwAUeOHAEgMzOT0aNH4+3tzXvvvcfcuXOZPHlyqfEVe/bs4YwzzmDt2rX861//4ttvv6Vfv35ceumlfPtt+T0QL774YkaMGMHXX3/NRRddxOOPP15vWld0HYXakpcFi56EdgOg50W2Zbs9yXqDd2vZDD8fjfuUqjPfPwQH1tdtHVqdDOc9U61LDh8+THZ2NjExMZWeU3zs4MGDnHrqqQB06tSJgQMHus8p/tYcHR1dKv149u3bx4IFC7jpppvw9/dn5MiRtG3blpkzZ3LHHeUXouvYsSMff/yx+/fk5GTuv/9+9u3bR5s2bdzltm3btlQdCgsLeeyxxxg2bFipIKR79+4MGTKEd955h7vuusudnp6ezpo1a4iIsHbdbdWqFQMGDGDu3LlcddVVxMXFkZqaynPPPccpp5zivq5kC8qUKVMwxrBkyRKioqwu4dGjR7Nnzx4mT57sDr6K3XfffUyYMAGAkSNHsnDhQj755BN3Wl3ST5ba8tu/IWM/nPsE2LSk8rwNBxj54hIAJgyOsSVPpVTTYoxxPP+CgoJSj2IffvghhYWF7vEPXl5eXHPNNfz+++9s2bKlXF5jx44t9fvJJ58MwO7du49bhy1btpCUlFSq9QHgrLPOomPHjixZsqRU+qBBg9xBQkXldOnShfDwcG677TY+/PBD90DPkubNm8f5559PWFhYqdc+evRo1q5dS3p6+nFfW+/evat8XbVFWxRqw+Ed8MsL1p4OHTyLtKuSk1/I7R+ucv9+5ekdbMlXKXWCqvlNvr6Ijo4mMDCQhISESs8pPta+fftq579kyRKGDx9eKq04OHn//ffp0KEDvXr1cjfrjx8/nmeffZb333+fJ598stR1kZGRpX7397eWqM/JyTluHVJSUgBo3bp1uWOtWrVyH/e0nLCwMBYtWsT06dOZOHEiGRkZ9OrVi6lTp3Lppdaeh0lJSbz//vu8//77Fdbp8OHDhIaGHrfMql5XbdFAwWnGwNx/gI8/jC4/UOhEfPXnXv7+2Vr37zNvPN2WfJVSTY+Pjw9nn302P/74Izk5ORWOUyjuUx8xYkS18z/ttNNYsWJFufSVK1eyceNGgFLf3ot98MEHTJ8+HS+vmjd8F38IHzhwoNyxAwcO0L9//2rn2bdvX7788ksKCgpYuXIlTz/9NFdccQVr166ld+/eREVFMWTIEB588MEKr2/Txp6NAGuDBgpO2/Al7FgI5z0PoeWj2eoqKjKlgoRfHxpB2/DAGuerlGq67r//fkaOHMnDDz/Miy++WOpYfHw8zz77LGeffTZnnHGGe5ZAdnZ2uXz8/f3LpTdr1qzCD+KZM2ciInzxxRflvk3Pnz+fZ555hsWLF1c7OPHz8ytXh27dutGyZUs+/fRTbrrpJnf6b7/9xq5du7jvvvuqVUZJPj4+DBw4kOnTp/Ptt9+yefNmevfuzZgxY1i2bBm9evUiMLBh/43WQMFJGQfh+wegzakw4Kaqz/dA98nz3D/fP7qbBglKqRo755xzmDZtGpMnTyYhIYHrrruOiIgIVq9ezTPPPENYWBgffPABAC1btiQqKopPP/2UU045heDgYGJjY4mKiqJnz5788ssvzJkzh1atWhEdHV3hIMn8/Hw+/fRThg4dyiWXXFLueN++fXnppZeYOXNmtQOFnj178t133zFmzBgiIiJo06YNbdq0Ydq0adx2221cc801XHPNNSQmJvLII4/QpUuXag8YnDNnDjNmzOCiiy4iNjaWrKwsXnnlFZo1a8agQYMAmDZtGqeffjpnn302f/vb34iJiSE1NZUNGzawc+dO3nnnnWqVWZd0MKNTjIFv/2bNdrjoTduWas4rsObVLv7HMP46vLMteSql1GOPPcb3339PVlYWEyZM4Nxzz+X111/nuuuuY+XKlXToYI2D8vLy4u233yY1NZWRI0cyYMAAZs+eDVjTF7t168YVV1zBgAEDSq1hUNKcOXM4dOgQN954Y4XHw8PDueSSS/jyyy/LTV2syr///W+Cg4O58MILGTBgADNmzADg1ltv5YMPPmD9+vWMHz+eBx54gFGjRrFkyRJCQkKqVUaXLl0IDAxk+vTpnHfeeUyYMAEfHx9+/PFH2rWzNuXr0KEDK1eupE+fPjz88MOMGjWKO+64gyVLlpxQF05dEqdHvNYn/fv3NytXrqydwlb8F767F8Y8a9sOkdsOZjDqXz8DkPDM2CrOVko5YfPmzfTo0aOuq6EUUPX7UURWGWOqPwijBG1RcMKBDTD/Eeg0HE6/1ZYsCwqL3EHC/aO72ZKnUkopVRUNFOyWnQqfXQ2B4XDxW2DDiF2AT1ccm6c7cdhJtuSplFJKVUUHM9qpqBBm3QppiTBhLjRraVvWO5KtfrqlDw53r2eulFJKOU0DBbsYA/MmWbtDjn0B2tu7tsG7vyYQFexHu4ggW/NVSimljke7Huzy68vwx1sw6G8w4GZbs956MAOAk9uF2ZqvUkopVRUNFOzw54ew4HHofSmMmm5r1oczc3n31wQApo3rbWveSimlVFW066GmVr8P395lzXC46A3bBi8WO/u5RWTlFQLQIUq7HZRSStUubVGoiRVvw7d3Qudz4P8+tfZzsFFOfqE7SFBKKaXqgrYonIiiIlg4HZa+CF3HwBXv2x4kAKzZc8T986J/DLM9f6WUUqoqGihUV342fH0HbPwKTrsBzv8nePs6UtSVM5YDMOfOs4iNDnakDKWUUup4tOuhOg5tg7dHwsavrUGLF7zkWJCwPSnD/XNHHZuglHKAiFT5KN7U6YYbbnDvY6CaFm1R8IQxsO4zmHMv+AbA1Z9Dl1GOFvnJH8dWYmwW4EwwopRq2pYtW1bq94svvpg+ffqU2szJ39/+blXVsGigUJW0vfDdfbB1HnQcDJe+DaFtHCtuZ3ImzZv50z7C2j7641vOcKwspVTTNnDgwFK/+/v7Ex0dXS69vsnNzdUAphZp10NlCvJg+Rvw2kCI/xlGPwXXz3Y0SCgoLGLEC0sY/9qv7E/PAWBATKRj5SmlVHX9+eefDBkyhKCgILp06cKbb75Z7pz4+Hiuvvpqmjdvjr+/P3379uWrr74qd968efMYNGgQgYGBhIWFcdFFF7Fly5ZS5wwbNoyzzjqL2bNnc+qpp+Lv78/rr7/OySefzMUXX1wuz8WLFyMizJ8/374X3cRpoFCWMdZAxddOh3kPWUsxT1wGg/4KXt42ZG/YdySbvtN+4KArGCg2689EAHYmZ/HWkp0A+HrrP5FSqn5IT0/nqquu4pprruGbb75hwIAB3HHHHSxatMh9zp49ezjjjDNYu3Yt//rXv/j222/p168fl156Kd9++637vHnz5jF27FhCQkL47LPPeOONN9iwYQNnnXUWiYmJpcrdunUrd911F3feeSfz58/nnHPO4Y477mDOnDns27ev1LlvvfUWsbGxnHvuuc7ejCZEux6KFRbApq+tpZgPrIMWPeHqL6DzSLBpE6b7/reWL1fvdf9+watLWfHISPfvv20/ZEs5Sqm6MWHehHJpo2NGc2X3K8kuyGbigonljo/vPJ6LOl9Eak4q9y6+t9zxv3T7C2Nix3Ag6wCTfplU7vj1va5nWPthxKfFExsWa88LqURGRgavv/46w4cPB+Dss8/mhx9+4JNPPnGnTZkyBWMMS5YsISoqCoDRo0ezZ88eJk+ezLhx4wB49NFH6dSpE99//z0+PtZH0aBBg+jatSsvvPACL774orvcQ4cO8cMPP9C3b193WmxsLA899BD//e9/eeyxx9znzZo1i6lTp+rmeTZqUl9XjYHEI9nu3/+7NJ43v//D6mJ4tR98eZM1/XH863D7UmvAok1vtrd/2VkqSAA4ue2xvRtSsvL4eo0VGV9wSmtevrIvax/XiFgpVX8EBQW5AwKwxjR06dKF3bt3u9PmzZvH+eefT1hYGAUFBe7H6NGjWbt2Lenp6WRlZbF69Wr+8pe/uIMEsD78Bw8ezJIlS0qVGxMTUypIAGjWrBnXXHMNb7/9NkVFRQC8++67GGOYMKF8wKZOXJNqUdiwL43Bzyxk7aPDCTuwjBbzn+Ncr5UgBRS0GYDP6KfYGjGEVuFBhNrQzVCq7MS0cmkL45LIzC0gxN+Hh2etB+Dy09rx/OV9bC1bKVU73h3zbqXHAn0Cj3s8IiDiuMdbBbc67nGnWxMAIiIiyqX5+/uTk3OsGzUpKYn333+f999/v8I8Dh8+jK+vL8YYWrduXe54q1at2LVrV6m0is4DmDhxIm+88QZz585l7NixzJgxg4svvpiWLVtW52WpKjSpQCGUo/zT902CX7kD8tI4yyuEjwvP4bPC4fx98MWsik9lxsylAKydfC7/W7mHm4fEVtmElVdQxLakDD79Yw+3nt2Jyd9s4PTYKK4Z2IFlOw5z6wer3OfeODiWd36N5+S2YaxPTOOrPxMZ3bMl8zYeAOAfo7s5dwOUUsphUVFRDBkyhAcffLDC423atKGgoAAR4cCBA+WOHzhwwN1lUayyv8G9e/dmyJAhvPXWWwQEBLB9+3beeuutmr8IVUqtBwoi0h74FzAKEGABcI8xZvdxL7SuDQCmA9cA4cAa4EFjzM+elN1RDjDKayV7W4zkqZ0nsaSoD7n4AXBbiQ9zgKv/u5wNiek8OXczvz98Di1DAyrMMyk9h9Of+sn9+wfLrUh40ZZknp0XV+rczi1CeOyCHtw2tBPB/j70fnw+j329gce+3uA+p7JylFKqIRgzZgzLli2jV69eBAYGVniOv78/p512Gp9//jlTpkzB29tqwd21axe//fYbd955p8flTZw4kWuuuYbU1FS6du3KiBEjbHkd6phaHaMgIkHAQqA7cD1wLdAFWCQinqxR/F/gFmAycAGwH5gvIn09KT/etOa03DcZtv1KfigaQC5+PH3JyRWeuyEx3f3zGU/9REGh1Qf27Lw4vlu3332sZJBQlR/uORsRoWVoACH+PgT5le7e+Pn+4ZVcqZRSDcO0adNIS0vj7LPPZubMmSxZsoSvv/6aJ554ghtvvNF93vTp09m2bRsXXHABs2fP5pNPPmHUqFGEhYVx3333eVzepZdeSnR0NL/++iu33XabEy+pyavtFoVbgE5AN2PMdgARWQdsA24DXqzsQhHpA1wF3GiMedeVtgTYCEwDxlVVeCaBNCvxkied153/O70Df+5O5X8rrYGGP903lHNeWFLu2nNf+pmCQsPulKMAzFrdgkv6Vb2c6Sv/dypdWoTQo3VouWN/Th5Ft0fnAfDylX11G2mlVIPXoUMHVq5cyZQpU3j44YdJTk4mKiqK3r17c/3117vPGzNmDN999x1Tp07liiuuwM/Pj2HDhvHcc8/Rpo3n69X4+voyfvx4Zs6cWSp/ZR8xxtReYSI/AQHGmMFl0pcAGGOGHufax4DHgHBjzNES6VOBh4BQY0zu8cpv16W3mfLOtzzx3WYA4p8+3933VVBYRKEx+Pt4E/PQd9V6XZec2pZ2kUEE+3nj6+3F4M7R+Pt44eMttIs4/of/ocxcftmWzMWn6hrqSjUEmzdvpkePHnVdDeVSUFBA586dGTJkCB988EFdV6fWVfV+FJFVxpj+NSmjtlsUegHfVJC+Ebjcg2vjSwYJJa71Azq7fq5Uq7AALjq1LU98t5k+7cNLDZDx8fZy34zHL+zJ1NmbmHPnWfRuG1Zl4PD85X3w9jqxaZTRIf4aJCilVDWlp6ezYcMGPv74Y/bs2VOt7gpVPbUdKEQCqRWkpwDl5914fm3x8SpFh/jzx8PnEB1S+TrhEwbHcsOZMe5AYtmkEexIymLt3iPsT8tmyoW9iDuQwQWvLuXNa0474SBBKaXUiVm9ejXDhw+nRYsWvPzyy+XWWVD2qYvpkRX1dXjySSsncq2I3ArcClbfGUALD2YWlGxtaB0WSOuwQM7qEu1O6902jIRnxnpQbaWUUnYbNmwYtdl13pTV9sqMqVT8zT+CilsLSko5zrXFx8sxxswwxvQ3xvRv3ry5xxVVSimlVO0HChuxxhqU1RPY5MG1sa4plmWvzQO217x6SimllCqptgOFb4GBItKpOEFEYoDBrmNVXetLiUGPIuID/AX4oaoZD0opZRdt8lb1QW29D2s7UPgPkAB8IyLjRWQc1iyIPYB73U0R6SgiBSIyuTjNGLMG+Ax4SURuFpFzgE+BWODx2nsJSqmmzNfXl+zs7KpPVMph2dnZ+Pr6Ol5OrQYKxpgsYASwFfgA+AiIB0YYYzJLnCqAdwX1mwC8CzwBfAe0B8YYY1Y7XHWllAKgRYsWJCYmcvToUW1ZUHXCGMPRo0dJTEykRYsWjpdX67MeXHs6XFrFOQlUMJvBGJMN3Ot6KKVUrQsNtVZZ3bdvH/n5+XVcG9VU+fr60rJlS/f70UlNavdIpZSyQ2hoaK38gVaqPqjtMQpKKaWUakA0UFBKKaVUpTRQUEoppVSlNFBQSimlVKU0UFBKKaVUpTRQUEoppVSlpCktGCIiGcCWuq5HIxcNHKrrSjRyeo+dp/e4duh9dl43Y0yzmmTQ1NZR2GKM6V/XlWjMRGSl3mNn6T12nt7j2qH32XkisrKmeWjXg1JKKaUqpYGCUkoppSrV1AKFGXVdgSZA77Hz9B47T+9x7dD77Lwa3+MmNZhRKaWUUtXT1FoUlFJKKVUNjT5QEJH2IvKFiKSJSLqIzBKRDnVdr4ZIRNqJyKsiskxEjoqIEZGYCs6LEJG3ReSQiGSJyAIRObkOqtzgiMhlIvKliOwSkWwR2SIiT4tIszLn6T0+QSIyWkQWisgBEckVkb0i8j8R6VnmPL3HNhKRea6/GU+USdf7fIJEZJjrnpZ9HClzXo3ucaMOFEQkCFgIdAeuB64FugCLRCS4LuvWQHUGrgBSgV8qOkFEBPgWGAPcCVwK+GLd83a1VM+G7B9AIfAw1j18A7gD+FFEvEDvsQ0igVXA34BzgUlAL2C5iHQEvcd2E5H/A/pUkK732R53AYNKPEYWH7DlHhtjGu0DuBvrj27nEmmxQAFwb13Xr6E9AK8SP98MGCCmzDnjXenDS6SFASnAK3X9Gur7A2heQdp1rns6Qu+xY/e9m+ue3qf32PZ7Gw4cAP7PdU+fKHFM73PN7u0w1/0beZxzanyPG3WLAjAOWG6M2V6cYIyJB37FunmqGowxRR6cNg7YZ4xZVOK6NGA2es+rZIxJriB5heu5retZ77H9Drue813Peo/t8xyw0RjzSQXH9D47r8b3uLEHCr2ADRWkbwR6VpCuau5497yDiITUcn0ag6Gu582uZ73HNhARbxHxE5EuwFtY33o/dR3We2wDETkLq0VsYiWn6H22x0ciUigih0Xk4zLj8Gp8jxt7oBCJ1Z9eVgoQUct1aSqOd89B73u1iEhbYBqwwBhTvBSr3mN7/A7kAluBU7C6dpJcx/Qe15CI+GIFYP80xlS2x47e55pJA17A6goeAUzHGp+wTERauM6p8T1uCns9VLRQhNR6LZoOQe+5LVyR/jdYY2omlDyE3mM7XAuEAp2wBpH+KCJnGWMS0HtshweBQODJ45yj97kGjDF/An+WSFoiIj8Df2ANcHwUG+5xYw8UUrGiqbIiqDjCUjWXQuX3HPS+e0REArBGKncChhpj9pY4rPfYBsaY4q6c30XkeyABeAi4Hb3HNeJq+n4E65uuv4j4lzjsLyLhQAZ6n21njFktIluBAa6kGt/jxt71sBGrf6asnsCmWq5LU3G8e77bGJNZy/VpcFxNtl8CpwPnG2PWlzlF77HNjDFHgO1YU4BB73FNdQICgA+xPoiKH2C13qQCJ6P32SklWxFqfI8be6DwLTBQRDoVJ7gWCBrsOqbs9y3QVkSKB+AhIqHAheg9r5JrrYSPgHOA8caY5RWcpvfYZiLSEmu9lR2uJL3HNbMGGF7BA6zgYThWYKb32WYi0h/oijUGB2y4x416rwfXokprgWysvhqDNdijGXCKRqvVJyKXuX48B6uJdiKQDCQbY5a4PuiWAu2B+7G+OUzCGizWxxizp/Zr3XCIyBtY9/VJYE6Zw3uNMXv1HteMiHwFrAbWAelYf1T/DrQCTjfGbNV77AwRMcCTxphHXb/rfa4BEfkIiMd6Px8BTsW6f0eBfsaYQ7bc47peMKIWFqTogNWMm47VJ/Y1ZRYJ0ke17qep5LG4xDmRwDtYfWNHgZ9cb8g6r399f2D1k1d2j6foPbblHj+ItTLjEde924I1Oj+mzHl6j+2/96UWXNL7XOP7OQkr4E3DWgNkD9Zuka3tvMeNukVBKaWUUjXT2McoKKWUUqoGNFBQSimlVKU0UFBKKaVUpTRQUEoppVSlNFBQSimlVKU0UFBKKaVUpTRQUEodl4hcJCL3lkkbJiJGRIbVTa2UUrVF11FQSh2XiLwHjDTGtCuRFoprzxRjTHpd1U0p5bzGvnukUsoBruCgon0olFKNjHY9KKUq5WpNuB5rUxnjeiRU1PUgIotFZKmIjBGRNSKSLSJ/isgZIuIjIk+JyH4RSRGR91x7sZQsK0hEnhWReBHJcz0/4lqrXilVR7RFQSl1PNOB5lh7249zpeUCYZWc3xl4HmtTq0zgOawd6r7F+ntzA9DDdU4S8ACAiPgA87G6M6YD64GBwGNY69TfZ+urUkp5TAMFpVSljDE7RCQZyDMltrw+ziDGKOBMY8xO13lewDdArDFmpOuc+SJyNnA5rkAB+D/gLGCoMeZnV9pPIgLwuIg8a4xJsu+VKaU8pU16Sik7bS0OElziXM/zy5wXB7QTVyQAjAF2Ab+5uil8XK0MPwC+WK0LSqk6oC0KSik7pZb5Pe846T6AN1AAtAA6Ym2VW5EouyqolKoeDRSUUvXBYSAeuKKS4wm1VxWlVEkaKCilqpILBDpcxjzgUiDTGBNX1clKqdqjgYJSqiqbgEgRuQNYCeQ4UMZHwASsAYwvAGsBP+AkrNkWFxljjjpQrlKqChooKKWq8jbWYMKngHCsQYc32FmAMSZfREYDDwG3ArFAFrAD+I5jYx2UUrVMl3BWSimlVKV0eqRSSimlKqWBglJKKaUqpYGCUkoppSqlgYJSSimlKqWBglJKKaUqpYGCUkoppSqlgYJSSimlKqWBglJKKaUqpYGCUkoppSr1/9GBgCBMsuvsAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 576x432 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(8, 6))\n",
    "plt.rcParams.update({\"font.size\": 16})\n",
    "plt.xlabel(\"time\")\n",
    "plt.ylabel(\"order parameter\")\n",
    "plt.xlim(t0, t1)\n",
    "plt.plot(sol_full.ts, sol_full.observables, label=f\"N={n_oscillator}\")\n",
    "plt.plot(sol_oa.ts, sol_oa.observables, label=\"Ott-Antonsen\")\n",
    "plt.plot([t0, t1], [orderparam_theory, orderparam_theory], label=\"Theory\", ls=\"dashed\")\n",
    "plt.legend()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "base",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.7 | packaged by conda-forge | (default, Sep 29 2021, 19:24:02) \n[Clang 11.1.0 ]"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "7d3977cd516aab4e7ac34ec0977b1608119a96f0bf9dd8c065a30d2af58323e0"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
