Skip to content

Instantly share code, notes, and snippets.

@juntyr
Last active February 20, 2024 09:01
Show Gist options
  • Save juntyr/c04f231ba8063a336744f1e1359f40d8 to your computer and use it in GitHub Desktop.
Save juntyr/c04f231ba8063a336744f1e1359f40d8 to your computer and use it in GitHub Desktop.
Display the source blob
Display the rendered blob
Raw
{
"metadata": {
"kernelspec": {
"name": "python",
"display_name": "Python (Pyodide)",
"language": "python"
},
"language_info": {
"codemirror_mode": {
"name": "python",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8"
}
},
"nbformat_minor": 4,
"nbformat": 4,
"cells": [
{
"cell_type": "code",
"source": "import numpy as np\nfrom matplotlib import pyplot as plt",
"metadata": {
"trusted": true
},
"outputs": [],
"execution_count": 1
},
{
"cell_type": "code",
"source": "def pdf(r, p, u):\n return (2.0 * np.pi * r) * p / (np.pi * u * (1 + r*r / u)**(p+1))",
"metadata": {
"trusted": true
},
"outputs": [],
"execution_count": 2
},
{
"cell_type": "code",
"source": "def sample_dispersal_distance(x, u):\n return np.sqrt(u * ((1.0 / x) - 1.0))",
"metadata": {
"trusted": true
},
"outputs": [],
"execution_count": 3
},
{
"cell_type": "code",
"source": "def sample_dispersal(u, size):\n angle = np.random.random(size) * np.pi * 2.0\n dist = sample_dispersal_distance(np.random.random(size), u)\n\n return (np.cos(angle) * dist, np.sin(angle) * dist)",
"metadata": {
"trusted": true
},
"outputs": [],
"execution_count": 4
},
{
"cell_type": "code",
"source": "x, y = sample_dispersal(5.0, 10000)\nplt.scatter(x, y, s=1)\nplt.axis(\"equal\")\nplt.show()",
"metadata": {
"trusted": true
},
"outputs": [
{
"output_type": "display_data",
"data": {
"text/plain": "<Figure size 640x480 with 1 Axes>",
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAjMAAAGdCAYAAADnrPLBAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAslklEQVR4nO3df3DU9Z3H8deGJJsQ2d0gISEYIsQWxR/gMTUExAnXDLGTttp63KmHwpWBNo06CiikKIgWwwCnJxwNenMFbrgTzpsprdqijApOjwAtBSO0UOEgQcIGWskupiE/v/cH3a+7+QGbkM3uJ3k+Zr4D+91PvvvZL8j35efz/ny/DsuyLAEAABgqLtodAAAAuBaEGQAAYDTCDAAAMBphBgAAGI0wAwAAjEaYAQAARiPMAAAAoxFmAACA0eKj3YG+0NbWppqaGg0ZMkQOhyPa3QEAAGGwLEsXL15UZmam4uK6Hn8ZEGGmpqZGWVlZ0e4GAADogdOnT+uGG27o8v0BEWaGDBki6fLJcLlcUe4NAAAIh9/vV1ZWln0d78qACDOBqSWXy0WYAQDAMFcrEaEAGAAAGI0wAwAAjEaYAQAARiPMAAAAoxFmAACA0QgzAADAaIQZAABgNMIMAAAwGmEGAAAYjTADAACMRpgBAABGI8wAAACjEWYAA23ZW6UpKz/Qlr1V0e4KAEQdYQYwUPmuEzpT16DyXSei3RUAiDrCDGCg4vwcjfQkqzg/J9pdAYCoc1iWZUW7E5Hm9/vldrvl8/nkcrmi3R0AABCGcK/fjMwAAACjEWYAAIDRCDMAAMBohBkAAGA0wgwAADAaYQYAABiNMAMAAIxGmAEAAEYjzAAAAKMRZgAAgNEIMwAAwGiEGQAAYDTCDAAAMBphBgAAGC2iYeajjz7St771LWVmZsrhcGj79u0h71uWpaVLl2rEiBFKTk5WQUGBPv3005A2ly5dUklJia6//npdd911euCBB1RbWxvJbgMAAINENMzU19dr/PjxWr9+fafvr1q1SmvXrtWGDRu0b98+paSkqLCwUJcuXbLbPPXUU3rrrbf05ptvavfu3aqpqdF3v/vdSHYbAAAYxGFZltUnH+Rw6Gc/+5nuv/9+SZdHZTIzM7VgwQItXLhQkuTz+ZSenq5NmzbpwQcflM/nU1pamv7rv/5Lf/d3fydJOnr0qG655RZVVFRo0qRJYX223++X2+2Wz+eTy+WKyPcDAAC9K9zrd9RqZk6ePCmv16uCggJ7n9vtVm5urioqKiRJBw4cUHNzc0ibm2++WaNGjbLbdKaxsVF+vz9kg7m27K3SlJUfaMveqmh3BQAQg6IWZrxeryQpPT09ZH96err9ntfrVWJiojweT5dtOlNWVia3221vWVlZvdt59KnyXSd0pq5B5btORLsrAIAY1C9XM5WWlsrn89nb6dOno90lXIPi/ByN9CSrOD8n2l0BAMSg+Gh9cEZGhiSptrZWI0aMsPfX1tZqwoQJdpumpibV1dWFjM7U1tbaP98Zp9Mpp9MZkX6j782clK2Zk7Kj3Q0AQIyK2sjM6NGjlZGRoffff9/e5/f7tW/fPuXl5UmSJk6cqISEhJA2x44dU3V1td0GAAAMbBEdmfniiy90/Phx+/XJkyd16NAhDR06VKNGjdKTTz6pH//4x/rKV76i0aNH67nnnlNmZqa94sntdmvOnDmaP3++hg4dKpfLpccff1x5eXlhr2QCAAD9W0TDzG9/+1tNmzbNfj1//nxJ0qxZs7Rp0yY988wzqq+v17x581RXV6e7775bO3bsUFJSkv0zr7zyiuLi4vTAAw+osbFRhYWF+slPfhLJbgMAAIP02X1moon7zAAAYJ6Yv88MAKD7uO8S0BFhBgAMwn2XgI4IMwBgEO67BHREzQwwQGzZW6XyXSdUnJ/DfXsAGIGaGQAhmJ4A0F8RZoABgukJAP0V00wAACAmMc0EAAAGBMIMAAAwGmEGAAAYjTADAACMRpgBAABGI8wAAACjEWYAAIDRCDMAAMBohBkAAGA0wgwAADAaYQYAABiNMAMAAIxGmAEAAEYjzAAAAKMRZgAAgNEIMwAAwGiEGQAAYDTCDAAAMBphBgAAGI0wAwAAjEaYAQAARiPMAIg5W/ZWacrKD7Rlb1W0uwLAAIQZADGnfNcJnalrUPmuE9HuCgADEGYAxJzi/ByN9CSrOD8n2l0BYACHZVlWtDsRaX6/X263Wz6fTy6XK9rdAQAAYQj3+s3IDAAAMBphBgAAGI0wAwAAjEaYAQAARiPMAL2A+6IAQPREPcw8//zzcjgcIdvNN99sv29ZlpYuXaoRI0YoOTlZBQUF+vTTT6PYY6Aj7osCANET9TAjSbfeeqvOnj1rb7/+9a/t91atWqW1a9dqw4YN2rdvn1JSUlRYWKhLly5FscdAKO6LAgDREx/tDkhSfHy8MjIyOuy3LEv/8i//omeffVb33XefJOk//uM/lJ6eru3bt+vBBx/s664CnZo5KVszJ2VH/HO27K1S+a4TKs7P6ZPPAwATxMTIzKeffqrMzEyNGTNG//iP/6jq6mpJ0smTJ+X1elVQUGC3dbvdys3NVUVFRZfHa2xslN/vD9mA/oDpLADoKOphJjc3V5s2bdKOHTtUXl6ukydPaurUqbp48aK8Xq8kKT09PeRn0tPT7fc6U1ZWJrfbbW9ZWVkR/Q5AX2E6CwA6irnHGdTV1Sk7O1svv/yybrnlFk2ZMkU1NTUaMWKE3ebv//7v5XA4tG3btk6P0djYqMbGRvu13+9XVlYWjzMAAMAgxj7OwOPx6Ktf/aqOHz9u19HU1taGtKmtre20xibA6XTK5XKFbAAGFpbLAwNHzIWZL774QidOnNCIESM0evRoZWRk6P3337ff9/v92rdvn/Ly8qLYSwCxjvoiYOCIephZuHChdu/erVOnTmnPnj36zne+o0GDBumhhx6Sw+HQk08+qR//+Mf6xS9+oU8++USPPvqoMjMzdf/990e76wBiGPVFwMAR9aXZn332mR566CH9+c9/Vlpamu6++27t3btXaWlpkqRnnnlG9fX1mjdvnurq6nT33Xdrx44dSkpKinLPAcSyvlouDyD6Yq4AOBLCLSACAACxw9gCYAAAgO4gzACGYZUOAIQizACGYZUOAIQizACGYZUOAISiABgY4Hh4JYBYRQEwgLAwbQXAdIQZoAeupQg31gp4mbZCfxZr/70hMggzQA9cy2hGrI2EzJyUrf9d/LdMMaFfirX/3hAZhBmgB65lNIOREKDv8N/bwEABMICrokgYQDRQAAyg1zBUDyCWEWaAAS6cAslYH6rvbpEnRaFA/0KYAQaY9hfycEZdrrVIONLhIfAdlv38cFifYdJIE8ELuDrCDDDAtL+Qd2fUpTsX1uC2kQ4Pxfk5GuSQWi2F9RmxPtIUzKTgBUQLBcDAAHMtxbxTVn6gM3UN8iQnKMUZf8VjBNoGQkOkC4j7a5Fyf/1eQDjCvX4TZgCELXBhrW9sUV1Ds0Z6kvW/i//2im25CAPoKVYzAehgy94qTVj+niYsf8+eKurO1FGgdmZh4dgrTtMQZAD0JcIMMICU7zqhuoZm1TU02zUYPanJuFpBcE+OSaErgJ4izAD9QLhBoDg/R8kJg+SQNDE71d7X2ShL4JhPvHGw28ue6xtb5ElO6FaBLYWuAHqKMAP0A+EGgZmTsjU0JVGWpANVF+x9gVGWzlYgvVNZozN1DVrz7rGwQk1g9Cfw+3BDkEkrjADEFsIM0A90Jwh01jYQYta8e8wORYF2RXdkypOcIF9Dc0hg6qz+Jvj4kq4YsNqPJvHASwA9RZgBDHGlqaRrDQKBURhJdtCZOSlbxfk59giOJWmQQ3YI6qz+JrgvCwvHypOcoPrGlk77zLQSgN5CmAEM0VsX/86OExhNWVg4NmTKadnPD4eEnOX33WYHpuL8HHmSE7qsjZk5KVspzvgOYaf9ZzKtBOBaxUe7AwDCE3zjud48TlfLqNe8e0ytluSQtLBwbIdRn5mTsq86EnS1z2JKCUBvYGQGMESkako6G6nZsrdKvr8W8bqTEzr9zPbTXp1Ng7Xvc7ijSyzTBtAdhBkgBkTi4t3VMQOBIrA6aWJ2aofpnvJdJ2Tp8qhM4FjtC37bB5PA62e3H9boxe/oiTcOduhTuFNL1NMA6A7CDBADInHx7mrEJXAPmMaWVp2pa9BHfzzfYcTn8v1o4mRJds1LcMHvmnePdQgm16ck2j9vSXrr45oOfQp3dIl6GgDdQZgBYkAkLt6dHTMQSFKc8XLGD+ryZ2dOytal5jb79cTsVBXn59gjNcGBJhBMjtT4Qo4RH+ewf9/VKFFX+1mmDaA7eNAkMEBs2VulNe8ek3S5oFdSSHFu8HszJ2XrKz/6pZrbvvznwSHJnRyvuoaWDse+Y6RbklR5JjTQeJITtLBwrD1K1P7BlMFP1u7qgZUABi6emh2EMANcOTgE3pOk5IQ4OeMH2XfxDVdyQpwa/jqa49DlqSYpNNC0XzF1rQ+k5IGWQP/GU7MBA0WyEPj6lEQNcnz5TKZgwVNIl5rbuh1kJNlBpjNdTRt1dzqp/fmhUBiARJgBYkokC4GP1PjUan35TCZJeuKNg8opfUf7T36uF++/TSM9yfrW+MuPLwjmaH/Qqwge7r3nq2m65blfdbnCqTsBrv35oVAYgESYAWJK8MW5t0Zpgp+x1P7C/05ljVqty7/OnJStidmp+sXHNapraJYn+ct7avZ0LvqOkW599Mfzamhuk/XXz2uvfUC50vduH14oFAYgUTMDxKy+KI799rpfq/KMT3eMdOsXj9+tGxe/06vHv7wEvE0Nza2SLoebP9c3qTg/R/tPfq53Kmt0a+aX+2ZOyqYoGICNmhnAcL01hRKYSmo/xbNlb5W9nPrP9U0RudtufWOLHWQS4hw6UuOzR2ECo0KVZ3whBbxMHQHoLkZmgH6iq5U9OaXvqNW6/MTrE2VF9v7ACMggh1R0R6be+rimx9NJ4QiscBrkkJbfd5v2n/xcv/jrjfUYhQHQGUZmgAEmuPYkuO6k6I5MO7AEC66liXSQ8SQn6PaRbrsfMydl667RQ5WccPnGfZ/XN3Z7ZIjnNwEIMCbMrF+/XjfeeKOSkpKUm5ur/fv3R7tLQEwJnp4JDjZrH7pTJ8qKtPahO0PaB4pnD1RdiGiQCaj+/C8hq6nWvHvMnoJqaG7r9goulmUDCDAizGzbtk3z58/XsmXL9Lvf/U7jx49XYWGhzp07F+2uARHVndGH4JU9V1sVFbyvOD9Hcd1de91NgfvWBNfCNLa02u97khO6XSNDbQ2AACNqZnJzc/W1r31N//qv/ypJamtrU1ZWlh5//HEtXrz4qj9PzQxM1Rsre9ofY8veKi37+WG1Wpfv2tvU0qbWPvhX4NvjM7X2oTvtxyoEAo4nOUGHlk2PfAcAGKff1Mw0NTXpwIEDKigosPfFxcWpoKBAFRUVnf5MY2Oj/H5/yAaYqDdGH9ofo3zXCTu8NDT3TZCRLj9Fe8Ly90KCzCDHl8+JAoCeivkw86c//Umtra1KT08P2Z+eni6v19vpz5SVlcntdttbVlZWX3QV6HXh3BTualNR7Y/RPhgNivAUU4AlhYzGeJITtPy+2ySJQl4A1yTmw0xPlJaWyufz2dvp06ej3SUgYsIphA0OPDMnZduriJITBoUs1+4Lo4YOVooz3n46N4W8AK5VzIeZYcOGadCgQaqtrQ3ZX1tbq4yMjE5/xul0yuVyhWxAfxXOVFT7wLCk6BaN9CRrSdEtkhTxAuBglWcu3zhvzbvHJF1+8GVXD8AEgHDEfJhJTEzUxIkT9f7779v72tra9P777ysvLy+KPQNiQ1dTUe1XLI30JGtidqqmrPxAkkJ+JprLAD7643m1Wpd/BYCeiPkwI0nz58/Xv/3bv2nz5s36wx/+oOLiYtXX1+uf/umfot01xLCBflO1wGjMsp8flnR5BOedypoOUzpPvHGwT+4zEyzcwt+B/mcIIDzxV28Sff/wD/+g8+fPa+nSpfJ6vZowYYJ27NjRoSgYCBY8tTJQnqocWPYsSfd8NU1eX4NaLdnhJbByKXhK562POz7JOpIcunwX4ECfFhaOtR/D0F776bHOHtcAAEaMzEjSY489pqqqKjU2Nmrfvn3Kzc2NdpcQ4wbiTdXKd51QXUOz6hqadaDqgpbfd5t9Dorzc+yVS4G78EpS0l+LgSPJoctPzB7pSdaL99+mA1UXQoJmVyu2urqrMQAEM+KmedeKm+ahP+jqQZLt2wSPzByouhByf5mJ2ak6UHXB/jXwXuAmesFGepJ0pu5Sr/S9/U3/2n+XcL8bIzPAwBLu9ZswAxiiu3cDDm4vKeRnr3RX4GCe5HjVNbRcU78dkl68/7ZOC5QD4SQw6sLTswEE6zd3AAZwWXemzbbsrVJ9Y4v9zKP2z2oKvBe8uunWTLeky+Ej4FqDjKQui4uDp40G4pQggN7DyAxgmHCmW640itPViE2gWDgSPMkJIYW+4U4tARjYGJkB+qlwCmGL83PkSU5QfWNLh2XNwaMgwTesy3An9Wm/w3lUAwCEgzADGCacKZmZk7KV4oxXXUNzh9ATHCIOVF1Qq3V5dZPX92Wx77fHZ4Y88qCnPMkJGulJ1sLCsVcMWABwLQgzgGHCHdFof9ff4ADxxBsHlVP6jq5PSbSDUdEdmRrkuBxk1j50p5YU3SJPcoKc8XFK6MHzDhySHWICgaqrgAUA14KaGcAgPakz6ax+ZvTid2TpcuA4ubLzB00Gfk5/bef+66hKc9vV/8kY5JCW33d5BVPw5weCDXUyAMJBzQzQD/XkxnGdTUsFbpR3pRvmBaaFHLq8IinFGa9l375VIz3JumOk2x7FCR61+fb4TI30JNtBpv3nz5yUbQcappoA9BZGZgCD9NYKoO4c52ptu9un7t4vB2iPlXADBzfNC0KYAa5N8MVDkn2X4YWFY694M7zeCD9AewTigYNpJmAA6+7TpoPbB4qDn3jjoP1+8PRW8POfOpvuutpUGEuyca24ySLaI8wA/VBngaKzkNJZ+3cqa9RqSe9Ufvk07eCLR6CWJnB34fbBKZYuNN0NdTADgRjtxUe7AwB638TsVHl9DZqYnWrvCw4pax+6M6R98Cqj/Sc/1zuVNSq6I9N+f+ak7JALR/DvA0P+gSdgB7eN9pRScEjjwgf0X4zMAP1Q8M3wAgL3kQkOKQHB/6e79qE7daKsqEPgCRY84nGlkZierL7q6nN6IpZGiQBEDgXAQD/QfgQk0iMigdEYT3KCUpzxYRX7SuqVe+QAGDhYzRSEMIP+rq8u+oFwMjE7VQeqLqi+sUV1Dc1hfW5P+hjtaSoA0cVqJmAA6avplMC00YGqC/rfxX+rhYVjQz73StNCPekjhZ4AwsHIDNCHTB9puFr/+3qEyNTzCCA8jMwAMehaC2Kj7WojJe1HXyK1NNr08wigdxFmgD7U31fXtA87kQodXZ1H7isDDExMMwG4Jlea8unr6SBWPwH9C9NMAPrElUZf+rqAt7+PfAHoHHcABnBNgu8eHG3t71QMYGBgmglAVLAiCcDVMM0EIKaxIglAbyHMALD15Wqg3q5vYSUTMHARZgDY+nK0pLeLgxnpAQYuwgwAm8mrgUzuO4BrQwEwAACISRQAAwCAAYEwAwxAFMsC6E8IM8AARLEsgP6EMAMMQBTLAuhPKAAGAAAxiQJgAAAwIBBmAACA0aIaZm688UY5HI6QbeXKlSFtqqurVVRUpMGDB2v48OF6+umn1dLSEqUeA4hFrM4CBrb4aHfghRde0Ny5c+3XQ4YMsX/f2tqqoqIiZWRkaM+ePTp79qweffRRJSQk6KWXXopGdwHEoODVWTyBGxh4oj7NNGTIEGVkZNhbSkqK/d57772n3//+99qyZYsmTJigb3zjG3rxxRe1fv16NTU1RbHXAGIJq7OAgS3qYWblypW6/vrrdeedd2r16tUhU0gVFRW6/fbblZ6ebu8rLCyU3+/XkSNHujxmY2Oj/H5/yAag/+rth1YCMEtUp5meeOIJ/c3f/I2GDh2qPXv2qLS0VGfPntXLL78sSfJ6vSFBRpL92uv1dnncsrIyLV++PHIdBwAAMaPXR2YWL17coai3/Xb06FFJ0vz585Wfn6877rhDP/jBD/TP//zPWrdunRobG6+pD6WlpfL5fPZ2+vTp3vhqAAAgBvX6yMyCBQs0e/bsK7YZM2ZMp/tzc3PV0tKiU6dOaezYscrIyND+/ftD2tTW1kqSMjIyujy+0+mU0+nsXscBAICRej3MpKWlKS0trUc/e+jQIcXFxWn48OGSpLy8PK1YsULnzp2z9+3cuVMul0vjxo3rtT4DAABzRa1mpqKiQvv27dO0adM0ZMgQVVRU6KmnntLMmTOVmpoqSZo+fbrGjRunRx55RKtWrZLX69Wzzz6rkpISRl4AAICkKIYZp9OprVu36vnnn1djY6NGjx6tp556SvPnz7fbDBo0SG+//baKi4uVl5enlJQUzZo1Sy+88EK0ug0AAGIMD5oEAAAxiQdNAgCAAYEwAwAAjEaYAQAARiPMAAAAoxFmAACA0QgzAADAaIQZAABgNMIMAAAwGmEGAAAYjTADAACMRpgBAABGI8wAAACjEWYAAIDRCDMAAMBohBkAAGA0wgwAADAaYQYAABiNMAMAAIxGmAEAAEYjzAAAAKMRZgAAgNEIMwAAwGiEGQAAYDTCDAAAMBphBgAAGI0wAwAAjEaYAQyzZW+Vpqz8QFv2VkW7KwAQEwgzgGHKd53QmboGle86Ee2uAEBMIMwAhinOz9FIT7KK83Oi3RUAiAkOy7KsaHci0vx+v9xut3w+n1wuV7S7AwAAwhDu9ZuRGQAAYDTCDAAAMBphBgAAGI0wAwAAjEaYAQAARiPMAAAAoxFmAACA0SIWZlasWKHJkydr8ODB8ng8nbaprq5WUVGRBg8erOHDh+vpp59WS0tLSJvKykpNnTpVSUlJysrK0qpVqyLVZQAAYKCIhZmmpibNmDFDxcXFnb7f2tqqoqIiNTU1ac+ePdq8ebM2bdqkpUuX2m38fr+mT5+u7OxsHThwQKtXr9bzzz+v119/PVLdBgAAhon4HYA3bdqkJ598UnV1dSH7f/WrX+mb3/ymampqlJ6eLknasGGDFi1apPPnzysxMVHl5eVasmSJvF6vEhMTJUmLFy/W9u3bdfTo0bD7wB2AAQAwT8zfAbiiokK33367HWQkqbCwUH6/X0eOHLHb3HPPPXaQCbQ5duyYLly40OWxGxsb5ff7QzYAANA/RS3MeL3ekCAjyX7t9XrDbtOZsrIyud1ue8vKyurNrgMAgBjSrTCzePFiORyOK27dmf6JlNLSUvl8Pns7ffp0tLsEAAAiJL47jRcsWKDZs2dfsc2YMWPCOlZGRob2798fsq+2ttZ+L/BrYF9XbTrjdDrldDrD6gcAADBbt8JMWlqa0tLSeuWD8/LytGLFCp07d07Dhw+XJO3cuVMul0vjxo2z2yxZskTNzc1KSEiw24wdO1apqam90g8AAGC2iNXMVFdX69ChQ6qurlZra6sOHTqkQ4cO6YsvvpAkTZ8+XePGjdMjjzyijz/+WO+++66effZZlZSU2KMqDz/8sBITEzVnzhwdOXJE27Zt06uvvqr58+dHqtsAAMAwEVuaPXv2bG3evLnD/g8//FD5+fmSpKqqKhUXF2vXrl1KSUnRrFmztHLlSsXHfzlgVFlZqZKSEv3mN7/RsGHD9Pjjj2vRokXd6gtLswEAME+41++I32cmFhBmAAAwT8zfZwYAAKA3EGYAAIDRCDMAAMBohBkAAGA0wgwAADAaYQYAABiNMAMAAIxGmAEAAEYjzAAAAKMRZgAAgNEIMwAAwGiEGQAAYDTCDAAAMBphBgAAGI0wAwAAjEaYAQAARiPMAAAAoxFmAACA0QgzAADAaIQZAABgNMIMAAAwGmEGAAAYjTAD9JIte6s0ZeUH2rK3KtpdAYABhTAD9JLyXSd0pq5B5btORLsrADCgEGaAXlKcn6ORnmQV5+dEuysAMKA4LMuyot2JSPP7/XK73fL5fHK5XNHuDgAACEO4129GZgAAgNEIMwAAwGiEGQAAYDTCDAAAMBphBgAAGI0wAwAAjEaYAQAARiPMAAAAoxFmAACA0QgzAADAaIQZAABgtIiFmRUrVmjy5MkaPHiwPB5Pp20cDkeHbevWrSFtKisrNXXqVCUlJSkrK0urVq2KVJcBAICB4iN14KamJs2YMUN5eXn693//9y7bbdy4Uffee6/9Ojj4+P1+TZ8+XQUFBdqwYYM++eQTfe9735PH49G8efMi1XUAAGCQiIWZ5cuXS5I2bdp0xXYej0cZGRmdvvef//mfampq0k9/+lMlJibq1ltv1aFDh/Tyyy8TZgAAgKQYqJkpKSnRsGHDdNddd+mnP/2pLMuy36uoqNA999yjxMREe19hYaGOHTumCxcudHnMxsZG+f3+kA0AAPRPUQ0zL7zwgv77v/9bO3fu1AMPPKAf/vCHWrdunf2+1+tVenp6yM8EXnu93i6PW1ZWJrfbbW9ZWVmR+QIAACDquhVmFi9e3GnRbvB29OjRsI/33HPPacqUKbrzzju1aNEiPfPMM1q9enW3v0R7paWl8vl89nb69OlrPiYAAIhN3aqZWbBggWbPnn3FNmPGjOlxZ3Jzc/Xiiy+qsbFRTqdTGRkZqq2tDWkTeN1VnY0kOZ1OOZ3OHvcDAACYo1thJi0tTWlpaZHqiw4dOqTU1FQ7iOTl5WnJkiVqbm5WQkKCJGnnzp0aO3asUlNTI9YPAABgjoitZqqurtbnn3+u6upqtba26tChQ5Kkm266Sdddd53eeust1dbWatKkSUpKStLOnTv10ksvaeHChfYxHn74YS1fvlxz5szRokWLdPjwYb366qt65ZVXItVtAABgGIcVvHyoF82ePVubN2/usP/DDz9Ufn6+duzYodLSUh0/flyWZemmm25ScXGx5s6dq7i4L0t5KisrVVJSot/85jcaNmyYHn/8cS1atKhbffH7/XK73fL5fHK5XNf83QAAQOSFe/2OWJiJJYQZAADME+71O+r3mQEAALgWhBkAAGA0wgwAADAaYQYAABiNMAMAAIxGmAEAAEYjzAAAAKMRZgAAgNEIMwAAwGiEGQAAYDTCDAAAMBphBgAAGI0wAwAAjEaYAQAARiPMAAAAoxFmAACA0QgzAADAaIQZAABgNMIMAAAwGmEGAAAYjTADAACMRpgBAABGI8wAAACjEWYAAIDRCDMAAMBohBkAAGA0wgwAADAaYQYAABiNMAMAAIxGmAEAAEYjzAAAAKMRZgAAgNEIMwAAwGiEGQAAYDTCDAAAMBphBgAAGI0wAwAAjBaRMHPq1CnNmTNHo0ePVnJysnJycrRs2TI1NTWFtKuurlZRUZEGDx6s4cOH6+mnn1ZLS0tIm8rKSk2dOlVJSUnKysrSqlWrItFlAABgqPhIHPTo0aNqa2vTa6+9pptuukmHDx/W3LlzVV9frzVr1kiSWltbVVRUpIyMDO3Zs0dnz57Vo48+qoSEBL300kuSJL/fr+nTp6ugoEAbNmzQJ598ou9973vyeDyaN29eJLoOAAAM47Asy+qLD1q9erXKy8v1f//3f5KkX/3qV/rmN7+pmpoapaenS5I2bNigRYsW6fz580pMTFR5ebmWLFkir9erxMRESdLixYu1fft2HT16NOzP9vv9crvd8vl8crlcvf/lAABArwv3+t1nNTM+n09Dhw61X1dUVOj222+3g4wkFRYWyu/368iRI3abe+65xw4ygTbHjh3ThQsXuvysxsZG+f3+kA0AAPRPfRJmjh8/rnXr1un73/++vc/r9YYEGUn2a6/XG3abzpSVlcntdttbVlZWr3wPAAAQe7oVZhYvXiyHw3HFrf30z5kzZ3TvvfdqxowZmjt3bq92viulpaXy+Xz2dvr06T75XAAA0Pe6VQC8YMECzZ49+4ptxowZY/++pqZG06ZN0+TJk/X666+HtMvIyND+/ftD9tXW1trvBX4N7OuqTWecTqecTueVvwwAAOgXuhVm0tLSlJaWFlbbM2fOaNq0aZo4caI2btyouLjQQaC8vDytWLFC586d0/DhwyVJO3fulMvl0rhx4+w2S5YsUXNzsxISEuw2Y8eOVWpqane6DgAA+qmI1MycOXNG+fn5GjVqlNasWaPz58/L6/WG1LlMnz5d48aN0yOPPKKPP/5Y7777rp599lmVlJTYoyoPP/ywEhMTNWfOHB05ckTbtm3Tq6++qvnz50ei2wAAwEARuc/Mzp07dfz4cR0/flw33HBDyHuBleCDBg3S22+/reLiYuXl5SklJUWzZs3SCy+8YLd1u9167733VFJSookTJ2rYsGFaunQp95gBAAC2PrvPTDRxnxkAAMwTc/eZAQAAiATCDAAAMBphBgAAGI0wAwAAjEaYAQAARiPMAAAAoxFmAACA0QgzAADAaIQZAABgNMIMAAAwGmEGAAAYjTADAACMRpgBAABGI8wAAACjEWYAAIDRCDMAAMBohBkAAGA0wgwAADAaYQYAABiNMAMAiFlb9lZpysoPtGVvVbS7ghhGmAEAxKzyXSd0pq5B5btORLsriGGEGQBAzCrOz9FIT7KK83Oi3RXEMIdlWVa0OxFpfr9fbrdbPp9PLpcr2t0BAABhCPf6zcgMAAAwGmEGAAAYjTADAACMRpgBAABGI8wAAACjEWYAAIDRCDMAAMBohBkAAGA0wgwAADAaYQYAABiNMAMAAIxGmAEAAEYjzAAAAKPFR7sDfSHwYHC/3x/lngAAgHAFrtuB63hXBkSYuXjxoiQpKysryj0BAADddfHiRbnd7i7fd1hXizv9QFtbm2pqajRkyBA5HI5rPp7f71dWVpZOnz4tl8vVCz1Ed3D+o4vzH12c/+ji/Pcty7J08eJFZWZmKi6u68qYATEyExcXpxtuuKHXj+tyufjLHEWc/+ji/EcX5z+6OP9950ojMgEUAAMAAKMRZgAAgNEIMz3gdDq1bNkyOZ3OaHdlQOL8RxfnP7o4/9HF+Y9NA6IAGAAA9F+MzAAAAKMRZgAAgNEIMwAAwGiEGQAAYDTCTBdOnTqlOXPmaPTo0UpOTlZOTo6WLVumpqamkHbV1dUqKirS4MGDNXz4cD399NNqaWkJaVNZWampU6cqKSlJWVlZWrVqVV9+FWOtWLFCkydP1uDBg+XxeDptw/nvW+vXr9eNN96opKQk5ebmav/+/dHuUr/w0Ucf6Vvf+pYyMzPlcDi0ffv2kPcty9LSpUs1YsQIJScnq6CgQJ9++mlIm0uXLqmkpETXX3+9rrvuOj3wwAOqra3tw29hprKyMn3ta1/TkCFDNHz4cN1///06duxYSBvOf+wjzHTh6NGjamtr02uvvaYjR47olVde0YYNG/SjH/3IbtPa2qqioiI1NTVpz5492rx5szZt2qSlS5fabfx+v6ZPn67s7GwdOHBAq1ev1vPPP6/XX389Gl/LKE1NTZoxY4aKi4s7fZ/z37e2bdum+fPna9myZfrd736n8ePHq7CwUOfOnYt214xXX1+v8ePHa/369Z2+v2rVKq1du1YbNmzQvn37lJKSosLCQl26dMlu89RTT+mtt97Sm2++qd27d6umpkbf/e53++orGGv37t0qKSnR3r17tXPnTjU3N2v69Omqr6+323D+DWAhbKtWrbJGjx5tv/7lL39pxcXFWV6v195XXl5uuVwuq7Gx0bIsy/rJT35ipaam2q8ty7IWLVpkjR07tu86briNGzdabre7w37Of9+66667rJKSEvt1a2urlZmZaZWVlUWxV/2PJOtnP/uZ/bqtrc3KyMiwVq9ebe+rq6uznE6n9cYbb9ivExISrDfffNNu84c//MGSZFVUVPRZ3/uDc+fOWZKs3bt3W5bF+TcFIzPd4PP5NHToUPt1RUWFbr/9dqWnp9v7CgsL5ff7deTIEbvNPffco8TExJA2x44d04ULF/qu8/0Q57/vNDU16cCBAyooKLD3xcXFqaCgQBUVFVHsWf938uRJeb3ekHPvdruVm5trn/sDBw6oubk5pM3NN9+sUaNG8efTTT6fT5Lsf+s5/2YgzITp+PHjWrdunb7//e/b+7xeb8iFVJL92uv1ht0GPcP57zt/+tOf1Nra2um55DxGVuD8Xunce71eJSYmdqgt48+ne9ra2vTkk09qypQpuu222yRx/k0x4MLM4sWL5XA4rrgdPXo05GfOnDmje++9VzNmzNDcuXOj1PP+oSfnHwD6QklJiQ4fPqytW7dGuyvopvhod6CvLViwQLNnz75imzFjxti/r6mp0bRp0zR58uQORaMZGRkdVnMEqtczMjLsX9tXtLdvM5B09/xfCee/7wwbNkyDBg3q9FxyHiMrcH5ra2s1YsQIe39tba0mTJhgt2lqalJdXV3I6AB/PuF77LHH9Pbbb+ujjz7SDTfcYO/n/Bsi2kU7seyzzz6zvvKVr1gPPvig1dLS0uH9QAFqbW2tve+1116zXC6XdenSJcuyvixAbWpqstuUlpZSgNoNVysA5vz3jbvuust67LHH7Netra3WyJEjKQDuZeqiAHjNmjX2Pp/P12kB6v/8z//YbY4ePUoBahja2tqskpISKzMz0/rjH//Y6fuc/9hHmOnCZ599Zt10003W17/+deuzzz6zzp49a28BLS0t1m233WZNnz7dOnTokLVjxw4rLS3NKi0ttdvU1dVZ6enp1iOPPGIdPnzY2rp1qzV48GDrtddei8bXMkpVVZV18OBBa/ny5dZ1111nHTx40Dp48KB18eJFy7I4/31t69atltPptDZt2mT9/ve/t+bNm2d5PJ6Q1WTomYsXL9p/vyVZL7/8snXw4EGrqqrKsizLWrlypeXxeKyf//znVmVlpXXfffdZo0ePthoaGuxj/OAHP7BGjRplffDBB9Zvf/tbKy8vz8rLy4vWVzJGcXGx5Xa7rV27doX8O/+Xv/zFbsP5j32EmS5s3LjRktTpFuzUqVPWN77xDSs5OdkaNmyYtWDBAqu5uTmkzccff2zdfffdltPptEaOHGmtXLmyL7+KsWbNmtXp+f/www/tNpz/vrVu3Tpr1KhRVmJionXXXXdZe/fujXaX+oUPP/yw07/rs2bNsizr8ujAc889Z6Wnp1tOp9P6+te/bh07dizkGA0NDdYPf/hDKzU11Ro8eLD1ne98J+R/vtC5rv6d37hxo92G8x/7HJZlWX01pQUAANDbBtxqJgAA0L8QZgAAgNEIMwAAwGiEGQAAYDTCDAAAMBphBgAAGI0wAwAAjEaYAQAARiPMAAAAoxFmAACA0QgzAADAaIQZAABgtP8HEoAq6pJ4WmwAAAAASUVORK5CYII="
},
"metadata": {}
}
],
"execution_count": 5
},
{
"cell_type": "code",
"source": "u = 10.0\n\ndists = sample_dispersal_distance(np.random.random(100000), u=u)\n\nplt.hist(dists, bins=250, range=(0, 20), density=True)\n\nplt.plot(np.linspace(0, 20, 100), pdf(r=np.linspace(0, 20, 100), p=1, u=u))\n\nplt.show()",
"metadata": {
"trusted": true
},
"outputs": [
{
"output_type": "display_data",
"data": {
"text/plain": "<Figure size 640x480 with 1 Axes>",
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAiwAAAGdCAYAAAAxCSikAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAA9hAAAPYQGoP6dpAABJzklEQVR4nO3deXxU1f3/8dedyQ4kgEACEgmLgCKCIkZwVzRYVOgmWitK1VZ/0q9KW0W/Fly+LVr3VhSrIliLotWiBQtCFDdQlEUWAUF2SUJYspCELHPP749JBiIJZJKZ3DuT9/PxmEdmJvfe+dzehnl7zrnnWMYYg4iIiIiLeZwuQERERORYFFhERETE9RRYRERExPUUWERERMT1FFhERETE9RRYRERExPUUWERERMT1FFhERETE9WKcLiBUbNtm165dtGnTBsuynC5HREREGsAYQ3FxMV26dMHjqb8dJWoCy65du0hPT3e6DBEREWmEHTt20LVr13p/HzWBpU2bNoD/hJOTkx2uRkRERBqiqKiI9PT0wPd4faImsNR0AyUnJyuwiIiIRJhjDefQoFsRERFxPQUWERERcT0FFhEREXE9BRYRERFxPQUWERERcT0FFhEREXE9BRYRERFxPQUWERERcT0FFhEREXE9BRYRERFxPQUWERERcT0FFhEREXE9BRYRERFxPQUWERERcT0FlgiWMWEuGRPmOl2GiIhI2CmwiIiIiOspsIiIiIjrKbBEKXUXiYhINFFgEREREddTYBERERHXU2ARERER11NgEREREddTYBERERHXU2ARERER11NgEREREddTYBERERHXU2ARERER11NgEREREddTYBERERHXU2ARERER11NgEREREddTYIkwWoFZRERaIgUWERERcT0FliikVhgREYk2CiwiIiLiegos0WD/VvjmHaiqcLoSERGRsIhxugBpmss9Syh5aiytrHI4rhdk/RkwgOV0aSIiIiGjFpZIVVXBpJgZPBP3N1pZ5VQYL+zdBDOv4pXYhznR2lnnbhkT5mqMi4iIRBwFlkhUuBOm/4ixMfMBmFJ1JYPLn4OzbwdvHOd5V/PfuAnc6n3X4UJFRERCQ4ElwqRbefD8ebDzSwpNEr+q+D2PVl1NIa3hkgfhti+Y5xtMjGXzh5hZsHOZ0yWLiIg0mQJLhPm1dy6U7oVO/RhR8Wc+sE+vvUH7HtxSeSdv+c7FYxnW/P1X9JyglhYREYlsjQosU6ZMISMjg4SEBDIzM1m6dGm927799ttccskldOzYkeTkZIYMGcL8+fOP2O7NN9+kb9++JCQk0L9/f957773GlBbdDhbxY++n/ufD/8xO06neTf9c+QsKTCtO8WxljPd9jV0REZGIFnRgmTVrFuPHj2fSpEksX76cAQMGkJWVxe7du+vc/uOPP+aSSy7hvffeY9myZVx44YVcccUVrFixIrDN4sWLueaaa7jxxhtZsWIFo0aNYtSoUaxZs6bxZxaNVs2itXUQOvSG7ucfddO9pPBI1dUA/C7mTVLZ1xwVioiIhIVljDHB7JCZmcngwYN55plnALBtm/T0dH77298yYcKEBh2jX79+jB49mokTJwIwevRoSkpKmDNnTmCbs846i4EDBzJ16tQGHbOoqIiUlBQKCwtJTk4O5pQigzEwJRP2bIDL/gKZvzlmi4mFzVtx93O6ZxNzfJmMq7w98LutD48Id8UiIiLH1NDv76BaWCoqKli2bBnDhg07dACPh2HDhrFkyZIGHcO2bYqLi2nfvn3gvSVLltQ6JkBWVtZRj1leXk5RUVGtR1Tb+gns2UCJiYcBVzdoF4OH/628kSrj4XLvF1zgWRneGkVERMIkqMCyZ88efD4fqamptd5PTU0lNze3Qcd47LHHOHDgAFdddVXgvdzc3KCPOXnyZFJSUgKP9PT0IM4kAi19AYB/+86BhJQG77bOdONl33AAHox5mXg0G66IiESeZr1LaObMmTzwwAO88cYbdOpU/4DRhrjnnnsoLCwMPHbs2BGiKl2o8HtY7+/+ecV3adC7P1X1U3aZ9pzgyeeX3gWhrk5ERCTsggosHTp0wOv1kpeXV+v9vLw80tLSjrrv66+/zk033cQbb7xxRPdPWlpa0MeMj48nOTm51iOa1LqrZ9l0MD7odjbfmuBbkkpI5G9VPwbgOu9CLOwQVioiIhJ+QQWWuLg4Bg0aRHZ2duA927bJzs5myJAh9e732muvMXbsWF577TVGjDhysOeQIUNqHRNgwYIFRz1mi1FV4Q8sAINvavRhZvvOpsgkkuHJ41zP6tDUJiIi0kyC7hIaP348L7zwAjNmzGDdunXceuutlJSUMHbsWMDfVTNmzJjA9jNnzmTMmDE8/vjjZGZmkpubS25uLoWFhYFtbr/9dubNm8fjjz/O+vXruf/++/nqq68YN25cCE4xwq3/D5Tshtap0PfyRh+mjATe8p0H+FtZ6qK5WkRExK2CDiyjR4/mscceY+LEiQwcOJCVK1cyb968wKDZnJwctm/fHtj+73//O1VVVdx222107tw58Lj99kO32A4dOpSZM2fy97//nQEDBvCvf/2L2bNnc8opp4TgFCPb52886n8y6AaIiWvSsf7huwSAizzLoWC7AoqIiESMmMbsNG7cuHpbP6ZPn17r9aJFixp0zJ///Of8/Oc/b0w5USuJg5xhbfC/GHhtk4+32XThU18/zvGuha9eBs5o8jFFRESaQ6MCizSPgZ5NxFg2O00HurbrFni/Ka0i//Bd4g8sy18hjgFUEKtWFhERcT0tfuhig6tbV760+4TsmAvtQeSY9lC6h8s8X4TsuCIiIuGkwOJiZ3j8geWrEAYWH15mVl0EwJgYzckiIiKRQYHFpbz4ON2zEQhtCwvA676LqDReBnk20s/aGtJji4iIhIMCi0udZG2jlVVOoUliozk+pMfOpy3z7MEAXKuZb0VEJAIosLjU4MO6gwyekN+C/EqVf4r/Ud7FJHIwZMcVEREJBwUWl6oZv7LM7h2W439p+rDN7kSSVc5FWsVZRERcToHFZfwtKXMY7PkWCP34lUMs5tpnAXC5d0mYPkNERCQ0FFhc6ARrN52sAspNDKtMj7B9zhyfP7Bc6FlJK8oC72sGXBERcRsFFheqmX9ltelBOU2bjv9ovjHd2GynkWBVcrFnedg+R0REpKkUWFwoHPOv1M1iTnW30BXez8P8WSIiIo2nwOJCZwTGr4RnwO3h5viGAHCe52uSKQn754mIiDSGAovLtKOIEz3fA+G7Q+hw35qubLSPJ96qYphnWdg/T0REpDEUWFxmUPXstt/ax1NAm2b4RCsw+PZydQuJiIhLKbC4wOF35TTf+JVDasaxnOtZTQoHmu1zRUREGkqBxWVqZrgN3/wrR/rOHM86+wRiLR+Xer9qts8VERFpKAUWN6kso7+1GfDPRNucarqFrvBoEjkREXEfBRY3+X45cZaPXNOOnaZjs350TbfQUM9a2lPUrJ8tIiJyLAosbvK9vztmmX0iYDXrR28zaay2M4ixbIZ7v2zWzxYRETkWBRY3yfePX1lvn+DIx//XlwnAJR6NYxEREXdRYHGTPf4J474zXRz5+PftQQAM9XxDEgcdqUFERKQuCiyuYSDfH1g2meMdqWCTOZ5tdifirUrO9awGtBCiiIi4gwKLS3SkEMoL8RmLrSbNoSosFla3smjWWxERcRMFFpfoVT0d/3bTiQpiHatjoX06ABd6V4Ltc6wOERGRwymwuERPaxfg3PiVGl/afSgySXSwimCnBt+KiIg7KLC4RE1gcWr8So0qYvjQHuh/8e1/Ha1FRESkhgKLS7ilhQUg2+fvFmKDAouIiLiDAotL1Ixh+c52PrAssk+l0nghfz0nWHlOlyMiIqLA4gatKKOLtQ9wRwtLEa1ZavcFYJhnucPViIiIKLC4Qg8rB4B8k0whrR2uxi+7+m4h3d4sIiJuoMDiAofGrzg74PZwNbc3n+lZTzIHHK5GRERaOgUWF3DT+JUa200qG+yuxFg2F3hWOV2OiIi0cAosLnDolmb3BBY41C10sVfjWERExFkKLC7Qy0W3NB9uYfXtzRd6VoKv0tliRESkRVNgcZqvkm5WLgCbbPeMYQFYaXqxxySTbJXCji+cLkdERFowBRan7d9KnOWj1MSTQ3unq6nFxsPH9qn+FxsX1PqdVnEWEZHmpMDitD3fAvCd6Yxx4eVY5Bvgf7JpobOFiIhIi+a+b8iWJn8D4L7xKzU+tk/FNhbkrYHC750uR0REWigFFqft2Qi4b/xKjQLasNL09L9QK4uIiDhEgcVpe9zdwgKwyDfQ/2Tj+47WISIiLZcCi5OMCbSwuDmwfGgP9D/Z/BFUVThai4iItEwKLE4qzoXyInzGYqtJc7qaeq0xGdCqI1QU6/ZmERFxhAKLk6rvENpuOlFBrMPF1M/ggV7D/C/ULSQiIg5QYHFSdWDZ5KJFD+tVE1g08FZERBygwOKkwBws7h2/EtDzIrA8sPsbKNzpdDUiItLCKLA4yeVzsNSS1B66Dgbgnr884XAxIiLS0iiwOKnmDiHb/YElY8JcHtt8AgAXeL52uBoREWlpFFicUlEKxf5Vmjebzg4X0zCLbP80/Wd71hBLlcPViIhIS6LA4pSi6mnu41pTQGtna2mgtSaDfJNCa+sgZ3g2OF2OiIi0IAosTinYDsCGg20By9FSGsrg4aPqVpYLPSudLUZERFoUBRanVN9ps8sc53AhwalZvfl8jWMREZFmpMDigIwJc/nrvz8A4HvTweFqgvOJ3R+fsejj2Uln9jpdjoiItBAKLA453vJ/2e+KsMBSSGtWml4AnOdd5XA1IiLSUiiwOKRLdevE9xHWJQTwkbqFRESkmSmwOOR4Kx+IvC4hgI/sUwE4x7MafJUOVyMiIi2BAosDPNikWfuAyOsSAlhterDPtCbZKoOdXzldjoiItAAKLA7oQCFxlo8q4yGPdk6XEzQbD59Ut7L87YWpDlcjIiItgQKLA7pWdwfl0h4fXoeraZyPfP7AonEsIiLSHBRYHNAlcIdQ5A24rVHTwnKqZwscyHe4GhERiXYKLA7oYu0BIjuw5NOWNXaG/8V3Hzhai4iIRD8FFgccXx1YIvEOocPV3C3EpoXOFiIiIlFPgcUBXSJ00rgfqpmPhe+ywbadLUZERKKaAosDama5jcRJ4w633JxIsUmE0r2Qs9LpckREJIopsDjg0KRxHR2upGmqiOEz+xT/i03ZzhYjIiJRTYGluR0sIsUqBSDHtHe4mKbTOBYREWkOCizNreh7AApNEgdIcriYpvu4ej4Wdi6Fsv3OFiMiIlFLgaW5FewAIr87qMb3dIQOfcDYsPkjp8sREZEopcDS3AprAktkD7itpdcw/89NC5ytQ0REopYCS3OrDiyRPGncD/3y42T/k03ZYAwZE+aSMWGus0WJiEhUaVRgmTJlChkZGSQkJJCZmcnSpUvr3TYnJ4df/OIX9O7dG4/Hwx133HHENtOnT8eyrFqPhISExpTmfoU7gcifNO5wX9p9KDXxUJwDeWudLkdERKJQ0IFl1qxZjB8/nkmTJrF8+XIGDBhAVlYWu3fvrnP78vJyOnbsyH333ceAAQPqPW5ycjI5OTmBx7Zt24ItLTJUB5ZInzTucOXEscQ+2f9C3UIiIhIGQQeWJ554gptvvpmxY8dy8sknM3XqVJKSkpg2bVqd22dkZPD0008zZswYUlJS6j2uZVmkpaUFHqmpqcGWFhkCg26jJ7DA4bc3az4WEREJvaACS0VFBcuWLWPYsGGHDuDxMGzYMJYsWdKkQg4cOEC3bt1IT09n5MiRrF0bhV0Lvioo3gVEX2BZZA8EoHLLYlpT6mwxIiISdYIKLHv27MHn8x3R+pGamkpubm6ji+jTpw/Tpk3jnXfe4dVXX8W2bYYOHcrOnTvr3ae8vJyioqJaD7cb+r//BGNTYbzkU39rUyTablLZbKcRa/kY6onCsCkiIo5yxV1CQ4YMYcyYMQwcOJDzzz+ft99+m44dO/L888/Xu8/kyZNJSUkJPNLT05ux4sapWaU5xxyHccf/9CH1ke0fo3SB52uHKxERkWgT1Ldmhw4d8Hq95OXl1Xo/Ly+PtLS0kBUVGxvLaaedxqZNm+rd5p577qGwsDDw2LFjR8g+P1y6VAeWaBpwe7iawHK+92vAOFuMiIhElaACS1xcHIMGDSI7+9DAStu2yc7OZsiQISEryufzsXr1ajp37lzvNvHx8SQnJ9d6uF3NKs27iJ45WA73uX0S5SaW46299LK+13wsIiISMjHB7jB+/Hiuv/56zjjjDM4880yeeuopSkpKGDt2LOBv+fj+++955ZVXAvusXLkS8A+szc/PZ+XKlcTFxXHyyf5bYR988EHOOussevXqRUFBAY8++ijbtm3jpptuCsEpukdNl9DOKG1hOUg8n9sncb53FRd4vmaTr6vTJYmISJQIOrCMHj2a/Px8Jk6cSG5uLgMHDmTevHmBgbg5OTls37691j6nnXZa4PmyZcuYOXMm3bp1Y+vWrQDs37+fm2++mdzcXNq1a8egQYNYvHhxINBEi2jvEgJ/t9D53lWc7/maF30jnC5HRESihGWMiYrBBkVFRaSkpFBYWOja7qENE/vRx7OT6yom8EnNvCVRpoe1iw/if0+5ieG08r9TSgJbH1ZwERGRujX0+zv6blVxK2MCXULRNgfL4TabzuywOxJvVXGW5xunyxERkSihwNJcDhbQ2joIRNfCh0eyWKTbm0VEJMQUWJpL9ZT8e00bDhLvcDHhVRNYLvSsRLc3i4hIKCiwNJcoXKW5PovtfpSbGNI9+fS0djldjoiIRAEFluYShas016eMBL6wTwLgIs8Kh6sREZFooMDSXAr9XULRPX7lkA9s/63s/m6hQzSZnIiINIYCS3Mp9i8OmWvaOVxI8/iwevXmwZ4NcLDQ2WJERCTiKbA0lwP+9ZfyTVtn62gm20wa39mdibV8sHmR0+WIiEiEU2BpLtWBZTdtna2jGdW0srwx8yVnCxERkYinwNJcWlgLCxwKLBd4vwbbdrYYERGJaAoszaGqHMr2A7C7BQWWL+2+HDAJdLIKIFeTyImISOMpsDSHA7sBqDBeCmjtcDHNp4JYPrNP8b/YuMDZYkREJKIpsDSHmu4g2gKWo6U0t5rbm/l2vrOFiIhIRFNgaQ6B8SspDhfS/Bb5/NP08/0yKNnjbDEiIhKxFFiaQ/UcLPktZA6Ww+XRnrV2N8DApoVOlyMiIhFKgaU5VI9haUkDbg9Xc7cQG993tA4REYlcCizN4UB1Cwstr0sI4ANf9TiWTQvx4nO2GBERiUgKLM0h0MLS8rqEAFaaXpDYDg4Wcrq10elyREQkAimwNIfAGJaW2cJi44FelwBwsXe5w9WIiEgkUmBpDi18DAsAfS4DYJhHgUVERIIX43QBUc+YFjkt/xF6XQyeGHqxi+5WjtPViIhIhFELS7iV7Qe7EoA9LXTQLQAJKZBxDgAXq5VFRESCpMASbtXjV/ab1lQQ63AxDuvzIwCGaRyLiIgESYEl3Kq7g1r0+JUavYcDcIa1AUr3OVyMiIhEEgWWcGvB0/IfoV031tnpxFi2Zr0VEZGgKLCEW00LCy1zDpYfWmgP8j/Z8J6zhYiISETRXULhVqwWFoCMCXMBGGidzm9jZsPGhVBVATFxzhYmIiIRQS0s4aYxLLV8bXr4w1tFMddOfMLpckREJEIosITZkq+/AVr4HCyHMXjIrl5bSJPIiYhIQymwhFlHqwCAfNo6Woeb1IxjucS7zD+xnoiIyDEosIRZp+rAoi6hQz61T+GgiaWrtQd2f+N0OSIiEgEUWMKpsoxkqxRQl9DhDhLPp/YpADz616cCA3JFRETqo8ASTtUDbstNLEUkOVyMu9TqFhIRETkGBZZwqrVKs+VoKW6T7Tsd21gM9HxHGnudLkdERFxOgSWcqtcR2q0Bt0fIpy3LzIkAXOr9yuFqRETE7RRYwikwLX9bZ+twqfm+wQBkeRRYRETk6BRYwiRjwlz++u6ngO4Qqs98+wwAMj3roETdQiIiUj8FljDqSCGgafnrs8OkstbuRoxl84c/P+x0OSIi4mIKLGEUmINFCx/Wa16gW+hLhysRERE3U2AJo8Ast2phqdd82x9YzvWshvJih6sRERG3UmAJI81ye2zfmq5sttOIt6pg44Jav8uYMFeTyomICKDAEjYWNh0CY1jaOluMq1mBVpb/zHpeAUVEROqkwBIm7SkmxrKxjcVekp0ux9VqxrFc6FlJPBUOVyMiIm6kwBImHS1/68o+2lBFjMPVuNsq04Mc057W1kHO9qxxuhwREXEhBZYw6WTtB9Qd1BAGD/N9/jlZNImciIjURYElTDQHS3BqxrFc4v0KfFUOVyMiIm6jvoow0RwswVlq92WfaU176wDX/PEJltj9nC5JRERcRC0sYaI5WILjw8uC6m6hEZ7PHa5GRETcRoElTDQHS/Dm2pkADPd+iRefw9WIiIibKLCEyaEWlraO1hFJFtv92Gda08Eq4izPN06XIyIiLqLAEiYdKQDUwhKMKmKY5zsTgMvVLSQiIodRYAmTmnlY8mnrbCER5j/2EMDfLRSD7hYSERE/BZZwqCihjVUGqIUlWEvtvuSbZNpZBzhHk8iJiEg1BZZwOLAbgDITRwkJDhcTWXx4+a/PP/hWdwuJiEgNBZZwKNkDUL2GkOVsLRFoju8sALK8XxFHpcPViIiIGyiwhENpdWAxWvSwMb40fcg17Ui2SjnXs8rpckRExAUUWMKhuoVln2njcCGRyeDhvepuocu96hYSEREFlvCobmHZh1pYGqumW+gSzzKoLHO4GhERcZoCSziUqEuoqVaYXnxvjqO1dRA2LXS6HBERcZgCSziU7gXUJdQUBg9zq1tZWPO2s8WIiIjjFFjCodZdQtJYNd1CfDsPyoudLUZERBylwBIOpRp0GwqrTA8222lQWQrr5jhdjoiIOEiBJRxKarqE1MLSNBazfef4n66a5WwpIiLiKAWWcChVl1Co/Ns+2/9ky0dQlONsMSIi4hgFllCrKPV3YaAuoVDYYVIh/SwwNqx+0+lyRETEIQosoVbdulJuYjhAosPFRIkBo/0/1S0kItJiKbCEWsnhk8ZpHaGQ6Pdj8MZB3hrI1QrOIiItkQJLqGkOltBLbAcnXup/rlYWEZEWSYEl1DTLbchlTJjLb1ad6H+x+k2wfc4WJCIizU6BJdQC6wiphSWUPrQHQkJbKM6BLR87XY6IiDSzRgWWKVOmkJGRQUJCApmZmSxdurTebXNycvjFL35B79698Xg83HHHHXVu9+abb9K3b18SEhLo378/7733XmNKc15gpWa1sIRSBbG8emCQ/4W6hUREWpygA8usWbMYP348kyZNYvny5QwYMICsrCx2795d5/bl5eV07NiR++67jwEDBtS5zeLFi7nmmmu48cYbWbFiBaNGjWLUqFGsWRN5AyxnfbQSUJdQOLztO9f/5Jt3oaLE2WJERKRZBR1YnnjiCW6++WbGjh3LySefzNSpU0lKSmLatGl1bp+RkcHTTz/NmDFjSElJqXObp59+muHDh/OHP/yBk046iYceeojTTz+dZ555JtjyHNfeKgLUJRQOy82JbLM7QWUJrJ/rdDkiItKMggosFRUVLFu2jGHDhh06gMfDsGHDWLJkSaOLWLJkSa1jAmRlZR31mOXl5RQVFdV6uMFx1YFFLSzhYPFvu3qq/hWvOluKiIg0q6ACy549e/D5fKSmptZ6PzU1ldzc3EYXkZubG/QxJ0+eTEpKSuCRnp7e6M8Ppfb4VxVWYAmPN6vOxzaWf6r+fVvImDCXjAlqbRERiXYRe5fQPffcQ2FhYeCxY8cOp0sCDu8SUmAJh+/pyKf2Kf4XK//pbDEiItJsggosHTp0wOv1kpeXV+v9vLw80tLSGl1EWlpa0MeMj48nOTm51sNxVeUkW2UA7NXEcWHzuu9C/5MV/8SL5mQREWkJggoscXFxDBo0iOzs7MB7tm2TnZ3NkCFDGl3EkCFDah0TYMGCBU06piOqZ7mtNF6KaOVwMdFroT2IfaY1FO/iPM8qp8sREZFmEHSX0Pjx43nhhReYMWMG69at49Zbb6WkpISxY8cC/q6aMWPG1Npn5cqVrFy5kgMHDpCfn8/KlSv55ptvAr+//fbbmTdvHo8//jjr16/n/vvv56uvvmLcuHFNPL1mVj0Hy37aoHWEwqeC2MAtzqO9iwA0lkVEJMrFBLvD6NGjyc/PZ+LEieTm5jJw4EDmzZsXGDSbk5PD9u3ba+1z2mmnBZ4vW7aMmTNn0q1bN7Zu3QrA0KFDmTlzJvfddx/33nsvJ554IrNnz+aUU05pwqk5oLRmWn51B4XbLN+F3BTzXy72LKcjBeTT1umSREQkjCxjjHG6iFAoKioiJSWFwsJC58azrHoT3r6Jz3z9uLbyf52poQV5O24ip3s2MbnyGp73XQHA1odHOFyViIgEo6Hf3xF7l5AraR2hZlUz+PYq7yIgKnK3iIjUQ4EllLRSc7Oa4xvCAZNAT08Og60NTpcjIiJhpMASSqVa+LA5lZLAHN9ZAIyOWeRsMSIiElYKLKFUoi6h5jarultohOdzkinRnUIiIlFKgSWUqudhUZdQ81lherHOTifRquBn3o+dLkdERMJEgSWUSvIBdQk1L4t/+C4F4Drv+1jYDtcjIiLhoMASSjWDbtUl1Kxm+86myCTR3ZPHuZ7VTpcjIiJhoMASKr5KOFgAqIWluZWSwL985wEwxvu+w9WIiEg4KLCESuk+AGxjUUBrh4tpef7huwSAizwrYf9WR2sREZHQU2AJldKadYRaY+t/1ma3xXTmY19/PJaBL19yuhwREQkxfbOGSonmYHHajOrBt6z4B1SWOVuMiIiElAJLqGhafsd9aJ/GDrsjlO2HNW85XY6IiISQAkuolGgOFqfZeHjVNwyA1f9+FKJjXU8REUGBJXQC0/KrhcVJs3wXcNDE0t+zFXZ+5XQ5IiISIgosoRKYg0UtLE4qoA3/8Q3xv/jiOWeLERGRkFFgCRUtfOgaL/uG+5+snQ0F2x2tRUREQkOBJUQ+X/MtoC4hN/jGZPCJ7xQwPvhcrSwiItFAgSVE2lMEqEvILV7wjQCgZMlL/ruGREQkoimwhEh7qxhQl5BbfGyfyno7nVZWOSyb7nQ5IiLSRAosoWDbtMMfWPaqS8glLF6o8rey8PlUqKpwthwREWkSBZZQKNuP1/LP+bFfE8e5xrv2UHJNOziQy+8m/tHpckREpAkUWEKh+g6hAtOKKmIcLkZqVBLDy1X+O4ZujpmrieRERCKYAkso1MzBovErrvOa7yIOmAT6enbApmynyxERkUZSYAkFrSPkWkW04nXfhf4Xi//qbDEiItJoCiyhoJWaXe3lquFUGQ9s+Qi+X+50OSIi0ggKLKFQWrPwoVpY3Oh7OvKOPdT/4uNHnS1GREQaRYElFGpaWDRpnGtNqRoFWLDhPcj52ulyREQkSAosoaB1hFxvs+kCp/zU/0KtLCIiEUeBJRRK9wGw37R2uBA5qvN+D1iw7j+Qt9bpakREJAgKLKFQVh1YUGBxtU4nwckj/c/VyiIiElEUWEKh1L+4XoEG3bpaxoS5DF9xlv/F2tmwe72j9YiISMMpsIRC9V1CmofF/dabE6Dv5YCBTx5zuhwREWkgBZamqiqHyhIACjSGJTKcfxcAvlX/gj0bHS5GREQaQoGlqaoH3FYZD0UkOVyMNEjnAdD7Mv+CldVjWTImzCVjwlyHCxMRkfoosDRV9YDbAloDlrO1SINkTJjLFavP9r9Y9YbuGBIRiQAKLE109dPvAeoOijSrTQ/m+s4EDGQ/6HQ5IiJyDAosTdSWA4BuaY5Ej1ddBZYXvp3HGZbuGBIRcTMFliZqZxUDuqU5Em02XeC0XwJwd+zrgHG2IBERqZcCSxO1q2lhUZdQRMpcPJiDJpbBnm+5yLPC6XJERKQeCixN1NZSl1Aky6M9L/uGA3BXzCywfQ5XJCIidVFgaaJ2qEso0j1XdQWFJom+nh3+u4ZERMR1FFiaSC0ska+I1jxXdaX/xYd/9k8GKCIirqLA0kTtLI1hiQbTfVnkmnZQuB2+mOp0OSIi8gMKLE2kLqHocJB4Hq0c7X/x0aNQnOdsQSIiUosCSxOpSyh6vG2fA8cPgopiTSYnIuIyCixNYduHJo5Tl1DEM3hg+CP+Fytf5cp7/upsQSIiEqDA0hTlhf4F9IAC1CUUDTKm7OYt3zkA3B87g4wJcxyuSEREQIGlaapXaj5gEqgkxuFiJFQeqbyGEhPP6Z5NjPJ85nQ5IiKCAkvTlO0HalZqlmixm3Y8U/VjACbEvsbJE94iY8Jch6sSEWnZFFiaorqFReNXos8033C22Z1Is/ZzW8xsp8sREWnxFFiaoqwmsGj8SrQpJ47/q/IvjHizdy69rJ1kTJirlhYREYcosDRFdQuLuoSi0wJ7EAt9pxFn+fhz7EtY2E6XJCLSYimwNEWZuoSim8XEyrGUmHjO9GxgtHeR0wWJiLRYCixNEWhhUZdQtNpFBx6vugqAe2Nm0pECZwsSEWmhFFiaQi0sLcJ0Xxar7O4kW6VMjH0l8L7GtIiINB8FlqaobmHZp0G3Uc3Gwz2VN1NlPFzh/Ry+fd/pkkREWhwFlqYo06DblmKtyWCa7zIAdv7zVk6e8JbDFYmItCwKLE2wK+d7QLc1txRPVv2UnaYDXa093B3zWq3fqXtIRCS8FFiaoB1aqbklKSOBCZU3A3B9zALO8awGUFAREWkGCiyNVVlGolUBQIEG3bYYn9r9mVF1CQCPxj5PcnVoFRGR8FJgaazqAbeVxssBEh0uRprTw1XXsMVOpbO1j/sPu2tIRETCR4GlsWoNuLWcrUWaVRkJ/K7yVnzG4ifeTxnuWep0SSIiUU+BpbG08GGLttz0ZqrvCgD+HPuiJpQTEQkzBZbGqpk0TrPctlhPVf2Mb+xutLcOMDn2BcA4XZKISNRSYGmsmmn51cLSYlUSw/jKWyk3MQzzruBX3nlOlyQiErUUWBpL0/ILsN6cwJ+qrgVgQsxM+H6ZwxWJiEQnBZbGKt0PaOFDgVd8l/Jf32DiLB+8ORbKCpwuSUQk6iiwNJZaWCTA4u7KX7Pd7ggF2+Dd34LReBYRkVBSYGmsmruENMutAEW0Ylzl/4AnFta9C1++6HRJIiJRpVGBZcqUKWRkZJCQkEBmZiZLlx59HopFixZx+umnEx8fT69evZg+fXqt30+fPh3Lsmo9EhISGlNa8ynToFupbZXpyYPlV/tfzL8Xdq0IrDGk6ftFRJom6MAya9Ysxo8fz6RJk1i+fDkDBgwgKyuL3bt317n9li1bGDFiBBdeeCErV67kjjvu4KabbmL+/Pm1tktOTiYnJyfw2LZtW+POqLkE5mHRGBY5ZJpvOO/7BoGvAl7/JcdR6HRJIiJRIejA8sQTT3DzzTczduxYTj75ZKZOnUpSUhLTpk2rc/upU6fSvXt3Hn/8cU466STGjRvHz372M5588sla21mWRVpaWuCRmprauDNqLmXqEpK6WPy+8hY4rhcU7eTZuKeJpcrpokREIl5QgaWiooJly5YxbNiwQwfweBg2bBhLliypc58lS5bU2h4gKyvriO0PHDhAt27dSE9PZ+TIkaxdu/aotZSXl1NUVFTr0WxsX+BOkAK1sMgPFNEKrn4N4pPJ9KxnYozWGxIRaaqgAsuePXvw+XxHtH6kpqaSm5tb5z65ubl1bl9UVERZWRkAffr0Ydq0abzzzju8+uqr2LbN0KFD2blzZ721TJ48mZSUlMAjPT09mFNpmoOF1MxqWkCr5vtciRwde8NPXsA2FtfFLOQX3mynKxIRiWiuuEtoyJAhjBkzhoEDB3L++efz9ttv07FjR55//vl697nnnnsoLCwMPHbs2NF8BVePXykyiVQR03yfK5Glz3Aeq/o5AA/ETIdtdbdCiojIsQX1bduhQwe8Xi95eXm13s/LyyMtLa3OfdLS0urcPjk5mcTExDr3iY2N5bTTTmPTpk311hIfH098fHww5YeO7hCSBnrWN5KTPdu43PsFe6ZdxY8rHmCH8bc4bn14hMPViYhEjqBaWOLi4hg0aBDZ2Yeat23bJjs7myFDhtS5z5AhQ2ptD7BgwYJ6twfw+XysXr2azp07B1Ne8ynVwofSUBZ/qPwNa+1udLCKmB77F9pS7HRRIiIRJ+guofHjx/PCCy8wY8YM1q1bx6233kpJSQljx44F/F01Y8aMCWx/yy23sHnzZu666y7Wr1/Ps88+yxtvvMGdd94Z2ObBBx/k/fffZ/PmzSxfvpxf/vKXbNu2jZtuuikEpxgGamGRYzh87pUyEhhbcRffm+Po6cnhxbjHiafC4QpFRCJL0AMwRo8eTX5+PhMnTiQ3N5eBAwcyb968wMDanJwctm/fHti+e/fuzJ07lzvvvJOnn36arl278uKLL5KVlRXYZv/+/dx8883k5ubSrl07Bg0axOLFizn55JNDcIphoFluJUi7accNFXfzr7j7OcPzLU/FTgH7SvB4nS5NRCQiWMZEx6InRUVFpKSkUFhYSHJycng/LPtB+ORxXq7K4oGq68P7WRJVMq11vBI3mXirCjJvgeEPg2U5XZaIiGMa+v3tiruEIo5muZVG+sKc5J9YDuCLqfDZ084WJCISIRRYGkOz3EoT/Mceyp8qf+F/sXASLH3B2YJERCKAAktjlGrQrTTNC77LeaZqpP/Fe7+HFa86W5CIiMspsDSGbmuWEHis6iqmVQ33v3j3t7DmrWPuo1WfRaSlUmBpjJouIbWwSJNYPFh1HZx+PRgb3v41rH/P6aJERFxJgSVYxqhLSELIgsufhP5XgV0Fb16v0CIiUgcFlmBVloKvHFCXkISIxwujnoOTR4KvAt64Dta87XRVIiKuopX7glXduoI3jlIcWstIokrNuJStf5oG3ltg9Zvw1o1QVQ4Dr3G4OhERd1ALS7Cqx6+Q2B7QhF8SQt4Y+PHzcNp1/jEts2+Br6Y5XZWIiCuohSVYZfv9PxPbOVuHRJ2alhaLLCbF5HFDzPsw506oKIGhv3W4OhERZ6mFJVg1gSWpvbN1SNQyeLi/6nqmVl3hf+P9+2DePWDbzhYmIuIgtbAEqzqwzN9c7nAhEt0sHq66mn2mNffGvgafPwtFu4hnlNOFiYg4Qi0swaoOLLqlWcLP4u++K/ifitvAEwvfzOaVuIcPtfKJiLQgCizBqgkstHK4EGkp3rXPhuvehvhkMj3r4aUs2Lcl8PuMCXM1A66IRD0FlmAFWlg0B4s0o+7nwa/mkWPaw54N7H/6HNjyca1NFFxEJJppDEuwStXCIs2vJoik8iDPxz3BQM9mqqaP5IGqMcAl6BZ7EYl2amEJlsawiIPyaM/oiom87TuHGMvmodjp/DnmRWKpcro0EZGwUmAJVmAMiwKLOKOcOMZX3sqfK6/BNha/iPmQ1+MeojN7Aa3oLCLRSYElWNWBpdCoS0ic5L+D6FeVf6DIJDHIs5G58fdwgWel04WJiISFAkswjFGXkLjKInsgIyr+xGo7g/bWAabH/YW7Yl4Hn7qIRCS6KLAEo7IssFKzuoTELXaYVH5a8QAzqi4B4P/FvAszroDCnQ5XJiISOgoswahuXakwXq3ULK5SQSyTqsZyW8X/UGwSYftieHYorHrT6dJEREJCgSUY1Ss1F9Ia3UYqbjTXPosrKv4Pjh8E5YXw9k3wr1/VOTuu5m0RkUiiwBKM6n/092v8irjYVtMZfjUfLrgHLC+secvf2vLdB06XJiLSaAoswdAtzRIhMv73fTLm9WfUwUlsttOgeBf848cw+zatRSQiEUmBJRiBW5oVWCQyrDS9GFHx58CAXFa+Cs+cCd+842xhIiJBUmAJRuCWZs3BIpGjjAQmVY3lp+WToENvKNkNb4xhauyTpLH3iLEsGtsiIm6kwBIMdQlJBFtm+tBn5738tWoUeGIY7v2S7Pjfc4v3XU3tLyKup8ASDE0aJxGunDieqLoKfv0RX9q9aWWVMyH2debF3Q2bsp0uT0SkXgoswSitua1ZXUIS4dJO4ecVkxhfcQv5Jpmenhx49Sfw+rV0s3Kdrk5E5AgKLMEoKwDUwiLRwuJt+zwuLn+caVXDqTIeWD+HhXF/YFLMDCjZG9hS41pExGkKLMGomYdFY1gkwh0ePopoxYNVY/hRxWQ+9A0g1vIxNmY+/HUgfPqkf0kKERGHKbAEQ2NYJIp9a9IZW3k311bcw1q7G5QXwcL74a+ncZ33feKodLpEEWnBLGOMcbqIUCgqKiIlJYXCwkKSk5PD8yH/lwZVZZxT/jQ7TcfwfIaIC1jYbLn6AHz4JyjcAcBO04G/Vf2Yt3znUkUMAFsfHuFkmSISBRr6/a0WloaqLIMqf9O45mGRaGfwwMBr6J33f/yx8gbyTFu6Wnt4JPYFsuN+zzXebLW4iEizUmBpqOoBt1XGwwESna1FpJlUEMs/fJdyXvlTPFT5S/aYZLp5djM59iU+jr8DFj8D5QecLlNEWgAFloaqXqm5QCs1SwtUThwv+X7EueVP8UDldeSY9qRZ++H9/2X/n/vw9H03QHGe02WKSBSLcbqAiBFYR0jdQdJylZHAy77L+KdvGKO8n3Kr9126e/K4Pebf8NRcOOVnXLa0P+tMN41vEZGQUgtLQ2lafpGACmJ5w3chF1c8zv+r+B+W2SeCrwK+nsl/4+/htdj/g7WzwVep+VtEJCTUwtJQuqVZ5Ag2Ht6zz+K9irPYOi4Vlkyhas1shni/gTevJ8+05c6Yi6BwIKQc73S5IhLB1MLSUGphETmqjGfyyFj2M84rf4pnqkaSb5JJtQq4PeZteKo/vHYNrJ8LPt1dJCLBUwtLQ6mFRVqYxnbl7KIDj1WN5umqn5Ll+ZLrYhaQ6VkPG97zP1p15PnCM3nTdz4LJ/+m1mdp3IuI1EeBpaECgUWDbkUaopIY5thDmFMxhK2/6wUrXoWvX4OSfH4TM5ffxMxl9cRnme07h44MIZ92gMKLiNRNgaWh1CUk0ngd+8ClD8HFE2HjAhb88zEu8Kykv2cr/T1buTfmnyy2+8HyfaQQR2Edf2cKMiItmwJLQ5X652EpVJeQSON5Y6Hvj7i50tCOIkZ4v2CU9zPO8HzLud418O44vor38pl9CizfB31HQFJ7p6sWERdQYGmo6pluC1CXkEiw6mod2U8yr/ou4VXfJXS1djPSs5jLvZ9zkmc7F3i/hnfHwZw7oNvZ0HcEXUhgFx0cOgMRcZoWP2yoJ/pB0U6uLH+IVaZn6I8vIgD0sHbxI88X/L7rOshbU+t3a+wMTrnwajjxUuhyGnh0o6NIpGvo97cCS0P9qQtUlnBe+ZNsN6mhP76I1LL14RGw9zv/nUXr38O3bQle67B/rpI6wImX+B89Lqy360hjX0TcraHf3+oSaoiqcqgsAXRbs0izOq4nGe/2AMbRnjFc7F3OhZ6VnOtZTZvSPf67jr5+DbCgy2k8sz2dT+3+vP7g7RAT53T1IhJCamFpiOI8eLw3PmPRq/wfGM23J+KoWKoY5PmWCzwrudCzkj6enT/YIAnSM6H7ufz4PS+rTA++e/hKZ4oVkaNSl1Ao7V4Pz2ayz7Tm9PK/h/bYItJkqezjHM8azvGu5hzPajpaRbV+X2wSWW6fyFK7L1/affja9GTDwz92qFoROZy6hEKpzH9Ls7qDRNwpj/a8ZZ/HW/Z5WNj0tnYyf6QFWz+hYN2HtLVKON+7ivO9qwAoNzHw0nPQdTCknwldz4TkzoHjadyLiPsosDRE9aRxdU1mJSLuYvCwwZwAZ42As27htAn/4SRrO4M9GxjsWc+Zng10sgpgxxf+xxL/fjtNB1baPfna7kmm1ZPVpruj5yEitSmwNISm5ReJOIfWQvLwjcngG18GM3xZgKGblccg61ueGFoBO77El7uGrtYeunr3cLn3CwB8xoIpj/JWbgfW2BlM+vW1kNafjEkfAYdaXzImzFVLjEgzUGBpCE3LLxJFLLaZNLaZNN7+FGAYrSjjVM9mBljfMcDjf3Sx9kH+On7qhZ96P4GX/4FtLLLj0lhnToCPv4HU/nRmLxgDluX0iYlENQWWhtBKzSJRrYREltj9WEI/8Pnf68h+TvVs5hRrK6d4tnCKZyudrX30tHLoSQ584G+JWZIAPHIfdDqZf2xuxQaTzv/d/HPo2BdaHefcSYlEGQWWhgiMYVGXkEhLkU87su1BZDMoEGI6UMhJnm2cZG3j3tN9kLeGqrz1xBwshO1LuK7mX9TpL/t/JnXwB5eOvaFDb+hwIhx3IqSkk3Hvf4HGD+zVwGBpaRRYGqI6sOw3bRwuRESctIcUPrFP5RNO5d6f+oPCyRNm08PKoY+1nT6enfSxdnDxcfugYBuU7oFtn/ofhzloYvlvXBpbTRos/BLa94TjejL4uc3k0xbwdy8dPk7m8NciLZECS0OU1tzWrBYWEfE7NKg3lvXmBNabE8CufisXEjnIuv/pCfkbePbNOfSwcuhh7SLDyiXBquQkawcnsQM+/TJwzC8ToNTEs910YrvpxAv3vcoO05ELPJ3YYTpBZRnEJob1fBSKxK0UWBpCtzWLSJDKSCDjr98DrYGrA+978ZFu7aa7lUt3K5eJQ+P9aybt+w7f/h0kWeX0tXbQlx1HHvRPf4BWnaBtOs/ExrDTdISluyClKyQf7/+Z2E4DgCUqKbA0wM6cXXS1NOhWRJrOh5etpjNbTWc+BKZ9AnARWx8eQd8J73C8lU83azfdrFxOsHaTbuWTbuXT1dpNslUGJbuhZDeXe6sP+N6c2h8Q2wpSjofkLv4Q06Zz9fMu0CYN2nSBVh210rVEHAWWBkiheuFDDboVkTDxd8nEBMIMDPjBFoatE4cy4qGZ/jljrHyOt/ZwvLWHztZeulh76WAV+Rdq3fOt/1EfTww5vmR2m7bVj3bc7m3rHz+zHmidytnPrGEPKY1ewuCHXUyar0aaSoHlWHyVtLHKALWwiIiTLDIeXAJ0Z209s/DGU0EXay9p1j7S2Ednay+drX1cd3IsFOf4Hwd2g11FZ2sfna19Rx7k9ZcA+Cyh+vXk/2FzWRJ7SWZwvz7QqiNPf76fvSaZB685399ak3Sc/46opPbgjQ3P6UuLp8ByLGUFgadFamERERcrJ44tpjNbTOda71/3i0OtHF58dKSAVGs/nayCwKMj++loFdHRKqCDVUhHCoi3qqC8iB6eInqQC+v8rTa313xz/GvGETUUmST2mTbMjmvNPtMG3n4XktozzpvvH2+T1N4/zqbmkdCW7g98gsGjFhg5Kq3WfCz538KUwRSaJAaUvxi644qIuJohmVI6WIUcRxEdrEL/c6uI9hRznFXIiB5xUJLPnvwc2lGM12rc14nPWBTRinbtO/L1Xosik8S5/XtBYlumfrGHIpPEXSPPhIS2kJAM8cmQkMzQJ7+imCRW//mn9Y7J0d1P7tfQ728FlmPZ/jlMy2KrncoFFU+G7rgiIlHEwiaFEo6zimjLAdpZB2hnFdOOYtpbB0jhAG0t//ttOUCKVUJbDpBoVYTk84tNIsUkcsAk0vuE4yG+NcS3YdaqQkpIoJhESkwC944aDPFtuOn19bx48wUQ1wri2lT/bAVxrcF7qPNBgSf8Gvr9rS6hYwmsI6TuIBGR+hg8FNCGgpoJNhv4n8LxVJBMCW2tEpIpIcUqIZnS6p8lJFulJFNKslVCCiUM7RoLB4ugvIiKkgLiLP80xG2sMtpQ5p9zb+f3geOP/uG33NzXAHgxDpjxeN1FxST4w0tsK+bH2ZSSADOmQmwr3l1XQKmJ5+qz+0JsEsQl+e/Mik2s3ifJ/zw2icunfkUZ8WTfPdz/fkyC/3ceb92f20AtdQCzAsux1MzBogG3IiIhV04c+cSRb9r53zhW0NlyqLWj94Q5xFNJG8poY5XS+rCff7+qD5QX85d3v6SNVUYrDtLKKqM1B2lFGa2tgyRxkFbWQY5P9EHFAbCr/J9RddD/YC99anqatmwC4MqarPHFomOe25z46idP/aH2OZtY4hOTaoeYmAS+2FnKQRPH+f1O8L8fkwCxCUc8v867CZbvrX4/3v/TG1f7dUz8D34X7/8ZwXP0KLAci1ZqFhFxlUOzDFuUE0c5cewxKf63qgNPxusAycDxxzze1gdGgDH0vucdf4ipDjJJlJNolVe/d+h5IuUkWeUkUc6vzuzE7C83kUQ5CVSQZJWTWP080SonkQraxVSCrzzwefFWJRws9D8Ok1kTjtavOmq9D8UC704/5nnVyRt/KLx44yAmrvq96p/V72VvLKSCGC4bcEL1+7H+3w35f9C+R+M+u4kaFVimTJnCo48+Sm5uLgMGDOBvf/sbZ555Zr3bL1q0iPHjx7N27VrS09O57777uOGGG2pt8+abb/LHP/6RrVu3cuKJJ/LII4/wox/9qDHlhVZgpWZ1CYmIRKPDl1moIJYC2hxq6TlGi8+vrhzBHYvnHn2jcv8Yn3gqScQfYhKsChKpYO6tZ0BlKVQd5LYZi0mggnirkgQqAs/jqSCeyjpfx1lVxFNZ/fD/Pg7/e3FU+u/0OpyvvFZ4qs/FNS1Ja76s/YtTR0dOYJk1axbjx49n6tSpZGZm8tRTT5GVlcWGDRvo1KnTEdtv2bKFESNGcMstt/DPf/6T7OxsbrrpJjp37kxWVhYAixcv5pprrmHy5MlcfvnlzJw5k1GjRrF8+XJOOeWUpp9lU6iFRURE6nEo7BydwcNB4jlIPPv9b/j3f3bPYVudFerysLCJo6pWiImrDjXv/zaTn/xtEXFWlf99/D9jqQoEoViq/K+p4vcXd/fPouyQoO8SyszMZPDgwTzzzDMA2LZNeno6v/3tb5kwYcIR2999993MnTuXNWvWBN67+uqrKSgoYN68eQCMHj2akpIS5sw5NMX0WWedxcCBA5k6dWqD6grbXUL/+hWseYsHK69jmu+y0B1XREQkgoRroG9Y7hKqqKhg2bJl3HPPPYH3PB4Pw4YNY8mSJXXus2TJEoYNG1brvaysLO64445a24wfP/6IbWbPnl1vLeXl5ZSXH2rWKiz09wUWFRU19HQaZt9uKDfkVcRg26WhPbaIiEiECPn36w+Oe6z2k6ACy549e/D5fKSmptZ6PzU1lfXr19e5T25ubp3bFxUVUVZWRmJiYr3b5Obm1lvL5MmTeeCBB454Pz09vaGnE6Snqx8iIiItT8pT4T1+cXExKSkp9f4+Yu8Suueee2q1yti2zb59+zjuuOOwQnjbVlFREenp6ezYsSO0XU0uEu3nqPOLfNF+jjq/yBft5xjO8zPGUFxcTJcuXY66XVCBpUOHDni9XvLy8mq9n5eXR1paWp37pKWl1bl9cnIyiYmJR92mvmMCxMfHEx8fX+u9tm3bNvRUgpacnByV/yc8XLSfo84v8kX7Oer8Il+0n2O4zu9oLSs16l58oR5xcXEMGjSI7OzswHu2bZOdnc2QIUPq3GfIkCG1tgdYsGBBre0bso2IiIi0XEEFFoDx48fzwgsvMGPGDNatW8ett95KSUkJY8eOBfxdNWPGjAlsf8stt7B582buuusu1q9fz7PPPssbb7zBnXfeGdjm9ttvZ968eTz++OOsX7+e+++/n6+++opx48aF4BRFREQk0gU9hmX06NHk5+czceJEcnNzGThwIPPmzQsMms3JyWH79u2B7bt3787cuXO58847efrpp+natSsvvvhiYA4WgKFDhzJz5kzuu+8+7r33Xk488URmz57t/Bws+LueJk2adET3UzSJ9nPU+UW+aD9HnV/ki/ZzdMP5Rc1qzSIiIhK9gu4SEhEREWluCiwiIiLiegosIiIi4noKLCIiIuJ6CizAlClTyMjIICEhgczMTJYuXXrU7RctWsTpp59OfHw8vXr1Yvr06c1TaCNMnjyZwYMH06ZNGzp16sSoUaPYsGHDUfdZtGgRlmUd8TjaUglOuf/++4+os2/fvkfdJ5KuX0ZGRp3X4rbbbqtz+0i4dh9//DFXXHEFXbp0wbKsI9YMM8YwceJEOnfuTGJiIsOGDWPjxo3HPO6bb75J3759SUhIoH///rz33nthOoOjO9r5VVZWcvfdd9O/f39atWpFly5dGDNmDLt27TrqMadPn37ENU1ISAjzmdTvWNfwhhtuOKLe4cOHH/O4kXANgTr/xizL4tFHH633mG66hg35XnDj32GLDyyzZs1i/PjxTJo0ieXLlzNgwACysrLYvXt3ndtv2bKFESNGcOGFF7Jy5UruuOMObrrpJubPn9/MlTfMRx99xG233cbnn3/OggULqKys5NJLL6WkpOSY+27YsIGcnJzAo1OnTs1QcfD69etXq85PP/203m0j7fp9+eWXtc5twYIFAPz85z8/6n5uvnYlJSUMGDCAKVOm1Pn7v/zlL/z1r39l6tSpfPHFF7Rq1YqsrCwOHjxY7zEXL17MNddcw4033siKFSsYNWoUo0aNqrVKfHM52vmVlpayfPly/vjHP7J8+XLefvttNmzYwJVXXnnM4yYnJ9e6ptu2bQtH+Q1yrGsIMHz48Fr1vvbaa0c9ZqRcQ6DWeeXk5DBt2jQsy+KnP/3pUY/rlmvYkO8FV/4dmhbuzDPPNLfddlvgtc/nM126dDGTJ0+uc/u77rrL9OvXr9Z7o0ePNllZWWGtM1R2795tAPPRRx/Vu82HH35oALN///7mK6yRJk2aZAYMGNDg7SP9+t1+++2mZ8+exrbtOn8fSdfOGGMA8+9//zvw2rZtk5aWZh599NHAewUFBSY+Pt689tpr9R7nqquuMiNGjKj1XmZmpvnNb34T8pqD8cPzq8vSpUsNYLZt21bvNi+//LJJSUkJbXEhUtc5Xn/99WbkyJFBHSeSr+HIkSPNRRdddNRt3HwNf/i94Na/wxbdwlJRUcGyZcsYNmxY4D2Px8OwYcNYsmRJnfssWbKk1vYAWVlZ9W7vNoWFhQC0b9/+mNsOHDiQzp07c8kll/DZZ5+Fu7RG27hxI126dKFHjx5ce+21tSYu/KFIvn4VFRW8+uqr/OpXvzrmAp+Rcu1+aMuWLeTm5ta6RikpKWRmZh71GkXydS0sLMSyrGOuhXbgwAG6detGeno6I0eOZO3atc1TYCMtWrSITp060adPH2699Vb27t171O0j9Rrm5eUxd+5cbrzxxmNu69Zr+MPvBbf+HbbowLJnzx58Pl9glt4aqamp9fb55+bm1rl9UVERZWVlYas1FGzb5o477uDss88+6izCnTt3ZurUqbz11lu89dZbpKenc8EFF7B8+fJmrLZhMjMzmT59OvPmzeO5555jy5YtnHvuuRQXF9e5fSRfv9mzZ1NQUMANN9xQ7zaRdO3qUvN3F8zfZM1+we7jBgcPHuTuu+/mmmuuOeqCcn369GHatGm88847vPrqq9i2zdChQ9m5c2czVttww4cP55VXXiE7O5tHHnmEjz76iMsuuwyfz1fvPpF6DWfMmEGbNm34yU9+ctTt3HoN6/pecOvfYdBT80vkuu2221izZs1Rx3iA/w+rT58+gddDhw7lu+++48knn+Qf//hHuMsMymWXXRZ4fuqpp5KZmUm3bt144403GvRfPJHkpZde4rLLLjvqEuyRdO1ausrKSq666iqMMTz33HNH3XbIkCG1FoMdOnQoJ510Es8//zwPPfRQuEsN2tVXXx143r9/f0499VR69uzJokWLuPjiix2sLPSmTZvGtddee8wBtG69hg39XnCDFt3C0qFDB7xeL3l5ebXez8vLIy0trc590tLS6tw+OTmZxMTEsNXaVOPGjWPOnDl8+OGHdO3aNej9zzzzTDZt2hSGykKrbdu29O7du95aI/X6bdu2jYULF3LTTTcFvW+kXDsg8HcXzN9kzX7B7uOkmrCybds2FixYcNTWlbrExsZy2mmnRcx17dGjBx06dDhqvZF2DQE++eQTNmzY0Ki/Szdcw/q+F9z6d9iiA0tcXByDBg0iOzs78J5t22RnZ9dKwocbMmRIre0BFixYUO/2TjPGMG7cOP7973/zwQcf0L1790YdZ+XKlXTu3DnE1YXegQMH+O677+qtNdKuX42XX36ZTp06MWLEiKD3jZRrB/7FUtPS0mpdo6KiIr744oujXqNIuq41YWXjxo0sXLiQ4447Luhj+Hw+Vq9eHTHXdefOnezdu/eo9UbSNazx0ksvMWjQIAYMGBD0vk5ew2N9L7j27zAkQ3cj2Ouvv27i4+PN9OnTzTfffGN+/etfm7Zt25rc3FxjjDETJkww1113XWD7zZs3m6SkJPOHP/zBrFu3zkyZMsV4vV4zb948p07hqG699VaTkpJiFi1aZHJycgKP0tLSwDY/PMcnn3zSzJ4922zcuNGsXr3a3H777cbj8ZiFCxc6cQpH9bvf/c4sWrTIbNmyxXz22Wdm2LBhpkOHDmb37t3GmMi/fsb471w74YQTzN13333E7yLx2hUXF5sVK1aYFStWGMA88cQTZsWKFYG7ZB5++GHTtm1b884775hVq1aZkSNHmu7du5uysrLAMa677jozYcKEwOvPPvvMxMTEmMcee8ysW7fOTJo0ycTGxprVq1e76vwqKirMlVdeabp27WpWrlxZ62+yvLy83vN74IEHzPz58813331nli1bZq6++mqTkJBg1q5d2+znZ8zRz7G4uNj8/ve/N0uWLDFbtmwxCxcuNKeffro58cQTzcGDBwPHiNRrWKOwsNAkJSWZ5557rs5juPkaNuR7wY1/hy0+sBhjzN/+9jdzwgknmLi4OHPmmWeazz//PPC766+/3px//vm1tv/www/NwIEDTVxcnOnRo4d5+eWXm7fgIAB1Pg6v+Yfn+Mgjj5iePXuahIQE0759e3PBBReYDz74oPmLb4DRo0ebzp07m7i4OHP88ceb0aNHm02bNgV+H+nXzxhj5s+fbwCzYcOGI34Xideu5tbrHz6uv/56Y4z/lso//vGPJjU11cTHx5uLL774iHM///zzA9vXeOONN0zv3r1NXFyc6devn5k7d24znVFtRzu/LVu21Ps3+eGHHwaO8cPzu+OOOwL/RqWmppof/ehHZvny5c1/ctWOdo6lpaXm0ksvNR07djSxsbGmW7du5uabbw78R2CNSL2GNZ5//nmTmJhoCgoK6jyGm69hQ74X3Ph3aFUXLyIiIuJaLXoMi4iIiEQGBRYRERFxPQUWERERcT0FFhEREXE9BRYRERFxPQUWERERcT0FFhEREXE9BRYRERFxPQUWERERcT0FFhEREXE9BRYRERFxPQUWERERcb3/D2JCp3oy3xWlAAAAAElFTkSuQmCC"
},
"metadata": {}
}
],
"execution_count": 6
},
{
"cell_type": "code",
"source": "def cdf(r, u, p):\n return 1.0 - (1.0 / (1.0 + r*r / u)**p)",
"metadata": {
"trusted": true
},
"outputs": [],
"execution_count": 7
},
{
"cell_type": "code",
"source": "plt.hist(dists, bins=100, range=(0, 20), density=True, cumulative=True)\n\nplt.plot(np.linspace(0, 20, 100), cdf(r=np.linspace(0, 20, 100), p=1, u=u))\n\nplt.show()",
"metadata": {
"trusted": true
},
"outputs": [
{
"output_type": "display_data",
"data": {
"text/plain": "<Figure size 640x480 with 1 Axes>",
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAiMAAAGdCAYAAADAAnMpAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8qNh9FAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA9Z0lEQVR4nO3deXhU5cH+8fucmWSSAAkCkoAEAggIqKyS4lK3aAS02vZXcUdUqpa2Iq0K7Vsob/sWl9allYpVAa3WolVxQUGMgooRLItSRQTZIQmLkMlGljnP74+EKZEkZEImZ2by/VzXXMmcnDO5D4dhbp6zWcYYIwAAAJfYbgcAAACtG2UEAAC4ijICAABcRRkBAACuoowAAABXUUYAAICrKCMAAMBVlBEAAOAqr9sBGsNxHO3evVvt2rWTZVluxwEAAI1gjFFRUZG6du0q265//CMqysju3buVnp7udgwAANAEO3bsULdu3er9eVSUkXbt2kmqXpnk5GSX0wAAgMbw+/1KT08Pfo7XJyrKyOFdM8nJyZQRAACizLEOseAAVgAA4CrKCAAAcBVlBAAAuIoyAgAAXEUZAQAArqKMAAAAV1FGAACAqygjAADAVZQRAADgKsoIAABwVchl5P3339dll12mrl27yrIsLViw4JjLLF26VEOHDpXP59PJJ5+sefPmNSEqAACIRSGXkZKSEg0aNEizZs1q1PxbtmzRmDFjdP7552vt2rWaNGmSbrnlFi1evDjksAAAIPaEfKO8UaNGadSoUY2ef/bs2erZs6f+9Kc/SZL69++vDz/8UA899JCys7ND/fUAACDGhP2YkdzcXGVlZdWalp2drdzc3HqXKS8vl9/vr/UAAACxKeSRkVDl5+crNTW11rTU1FT5/X6VlZUpMTHxqGVmzpypGTNmhDsaAAAxL2PKwmPOs/XeMS2QpH5hLyNNMXXqVE2ePDn43O/3Kz093cVEAACEV2NKQ6wKexlJS0tTQUFBrWkFBQVKTk6uc1REknw+n3w+X7ijAQAQktZcGMIp7GVk5MiRevPNN2tNW7JkiUaOHBnuXw0AQL0oFpEj5DJSXFysTZs2BZ9v2bJFa9euVYcOHdS9e3dNnTpVu3bt0jPPPCNJuu222/Too4/q7rvv1k033aR3331XL7zwghYu5C8BACB0lIjYE3IZ+fe//63zzz8/+PzwsR3jxo3TvHnzlJeXp+3btwd/3rNnTy1cuFB33nmnHnnkEXXr1k1PPvkkp/UCQCtHqcBhljHGuB3iWPx+v1JSUlRYWKjk5GS34wAAjoGiEV3CdTZNYz+/I/JsGgBA5KBYxBZLjnyqVKLKlagKJVgVUt6nUvvuUuIJrmSijABAK0bRiDyWHCWpvPphHVKSypWociVZ5UrSoeoSYVUcMf2QElURnCdR5UpQhRKt8uD0BKsiWD4SrYqjf+njkn40Txr4/ZZeXUmUEQCIWRSNlmCUoAq11SG1scqqv6pMbaxDaqNDamMdUluVKanm+zZHfE3SoWDBaBP8vrpQtJRyE6dDilNKcopkeVrs934bZQQAohBF43gZJalcbVWmdlapklWqdlap2qlMba0ytVOp2lllaquaxxE/a6NDalvzszY6JK/lhCWhYyyVyqcy+VRqfLW+L1OCSoPf+1SmeJWaBB1SfPV8pnraIfl0yNRMq3leZuJVJp8OKV5OzV1htv6CK7ACAGpQMhrLKFHlSlGJUqwSJatUKVaJUlSi5JrnyVapklWiZKtU7Y543q6mbDRniXCMpRIlVD9MgoqVqFKToGIlqPTIafKpxPx3Wql8KlGCyowvuHxpzfRDipdkNVvGSEYZAQAXUT6M2qpM7a1inaBinWAVqb1K1D74tVgpVrFSVKL2Vonaq7imdBQr3goc92+vMraKlKRik6giJalIiSqq+b7EJFRPM9XTDxeKwz8rVqKKTaKKlagyxcuE/96zMYsyAgBh0hqLhk8V6ii/Olh+dbSK1EF+dbCKdIJVpA6q+WoV6QQV6QSrWO1VrLjjKBWVxqNCtZHfJMmvNio0bWo995skFSkp+LzIJKpQbYIFo0w+tZbRh0hGGQGAZhCrxcOSoxNUrE5WYfVDhepk+dXRKlQn1Xy1/OogvzpafrVp4sGXZSZeB9RWhaatDpi2Nd+30UG100HTRgdrnhce/lozrZQyERMoIwBwDLFYNBJUrs7WQXXWAXW2DupEq1CdrQM6UYU6seb5idZBdZQ/5GMryo1X3yhZB0w77TftdEDttN9UP/9G7XTAtNMBtdVB0zb4/SFxc9TWjDICAN8SzeXDlqNOKlSa9Y3SrG+Uah2o+f6AOutAzfMDSrZKQ3rdb0xb7TMp1Q+laL9J1j6Tov1K1n5T86j5vliJYrQCoaCMAGjVoqt4GHWUX12t/TWPfepifaOu1v7/lg8daPQxGGUmXntMe+1R++qv5gTtNSnaq/baW1M89pr22q9kVfFxgTDibxeAViWSy4ctR511QN2svepm7av5ulcnWfvU1dqvk6x9SrAqj/k6AWNpj05QgTlBeaaD8k0HFZjq5wU10/ea9vIrSYxgIBJQRgDEjEguGoe1U6m6WwXqYRUo3dqr7tYepdc8TrL2HfN0VcdY2qsU7TadtMt0VF7w0SH4dZ9SFJB7V9MEQkUZARCVIrl4pKhYPa18ZVj5yrCrv/aw9qiHla8OVnGDy1YZW7tNR+0wnbXLdNJOc6J2q2PN952Ubzqqkn+6EWP4Gw0gKkRa+YhTlXpY+ept7VZvK0+97Dz1tPLUy8rTCccoHHtNsnaYztpmUrXDnKgdprO2O6nVZUMdGNVAq0MZARBxIql4JOqQelu71dfaqT72Lp1s7VYva7d6WAUNnvKab07QVpOmLU6atplUbTXVX7ebzipRYguuARD5KCMAXBcJ5SNelTWlY4f62Turv1o7lW7vrXeZIpOozaaLvjZdtdnpoi2mi7aYNG0xaSpTQgumB6IbZQRAq3OiDmigvU2nWNt1ir1dp1g71NvaXe8psXtNsjY53bTRnKRNpqs2mZO0yTlJe9RenI0CHD/KCIAW1bKjIEbp1h6dZm3RQHurBlrbNNDeqhOtwjrnLjRJ2mDStcFJ1waTro1ON31lTtIBJbdgZqD1oYwACKuWKx9GXbVfg+yvNcjerFOtzTrV3qr2VslRczrG0mbTRetNd613emi96a4vne7KUwcx0gG0PMoIgGbTkqMebVWqQfbXGmptDBaQukY8yo1XX5ru+tzpoS9Mhj53MvSlSeeYDiCCUEYARAGjDCtfZ9gbNNTaqCH2JvW1dsq2TK25Ko1HX5p0feb01meml/7j9NRXphvX5QAiHO9QAE0WrpEQjwIaYG3TCHu9zrC/0jB7g060/EfNt905UWtMH33q9NZap7c+NxkqV3xYMgEIH8oIgEYJ5y4YjwI63dqs79jrlWmv1zD7K7WzymrNU268+tT01iqnr1Y7fbTWOVl71T5smQC0HMoIgDqFs3xYctTf2q6R9uc60/5CI+wvjyoffpOkT5x+Wumcok+cfvqP6akKxYUtEwD3UEYAtIhUfaNzPOt0tr1OZ9v/Uadv7XY5aNroY2eAVjinaIXTX1+a7nJku5QWQEuijAAIyyiIV1Uabn+l8+y1Os/+VKfYO2r9vMT4tMLpr4+cgcp1BuoL012G8gG0SpQRAM2mg/y6wLNGF9hrdLa9TslH7HpxjKXPTC994JymDwKnaY3pw1kuACRRRoBWqTlHQnpbu3SRvUpZntUaam2sdbrtPpOsZc7pWhoYrA+c03RQ7Zrt9wKIHZQRIMY1/y4Yo4HWVl3i+USj7JU62d5d66frnAzlOEP1bmCI1pme7HoBcEyUEQCNYHSqtUWXej7WGHtFrTvZVhiPcp2BWuIMU05gqPLU0cWcAKIRZQRAvfpZ23WZJ1eX2h8rwy4ITi8z8VrqDNJbgRF6zxmiIiW5mBJAtKOMADHmeHfLpGm/Lvd8pCs8H6r/EWfAlJl45ThD9EZgpJY6g3RIvuONCgCSKCMAJCXqkEbbK/UDzwcaaX8RPAi13Hi11BmsNwLfUY4zVKXcXA5AGFBGgCh2fKMgRkOtjfqRZ5ku9Xxc6wqoK5xT9ErgbL0ZGCG/2h5/UABoAGUEaGVSVKwfej7Q1Z531cfeFZy+zemsFwPnaoFztnaaE11MCKC1oYwAUaTpIyFGQ6xNutabo0vtXCVYlZKkUuPTm06mXqw6VytNP07DBeAKyggQw3yq0GWeXN3oWaxT7a3B6V84PfRc4EK9GjhTxZwJA8BllBEgQh3P8SCddUDXeZfoGs+7wRvSHTJxesMZqeeqLtQac7Ikq5mSAsDxoYwAMaSftV0/9i7U9+yPFGcFJEm7TEf9veoi/TNwPpdjBxCRKCNA1DP6jr1et3pe1/meT4NTVzr9NLfqEr3tDFdAHhfzAUDDKCNAhAh9t4zRhfZq/cy7QIPtryVJAWPpLWeE/lZ1qT4zvZs/JACEAWUEiDKWHF1if6KfeRdogL1NUvXxIC8EztOTgdHablJdTggAoaGMAC5oysGplhyNslfqDu/L6mfvlCQVmwT9PXCRnqwarf1Kae6YANAiKCNAxDM6z16ru7wvaGDNSIjfJGluIFtzqy7hoFQAUY8yAkSwTGu97oqbr+H2V5KkIpOopwKjNKdqlPxq43I6AGgelBGgBYS6W6a3tUtTvM/rIs9qSdXHhMwLZOvxqkt1QMnhiAgArqGMABGkowp1h/dlXePJkddyVGVsPR+4QH+p+r726AS34wFAWFBGgAgQpyqN97yln3kXBO+euyQwTDOrrtZm09XldAAQXpQRoJmFukvmbHudZnjnqbedJ0n6zOmpP1Rdq4+dAeGIBwARhzICuOQk7dX/xD2rUZ5PJEl7TbLuq7paLwXO4e65AFoVygjQwjwK6GbPm7rT+5ISrQpVGVtPB7L1cNUPVcQddAG0QpQR4DiFslvmNGuz7o17Ini9kI+d/ppWeaO+MunhigcAEY8yArSAJB3SZO+LGu9ZJI9ldMC01f9VXat/Bb4ryXI7HgC4ijIChFmmtV4PxM1Wd3uvJOmVwFn6feV1XL4dAGpQRoAw8alCd3vna7xnkWzLaKfppF9X3qxlziC3owFARKGMACFo7PEhg61N+lPcY8HTdZ+vOl//V3WtijlAFQCOQhkBmpEtRxM9CzTJ+5I8llGBaa97KidoqTPE7WgAELEoI0AzSdN+PRI/S5n2l5KkVwNnalrljSpUW5eTAUBko4wADWjsbpmL7U90X9wTOsEqVrFJ0G8qx+sV55wwpwOA2NCkyzzOmjVLGRkZSkhIUGZmplauXNng/M8995wGDRqkpKQkdenSRTfddJP279/fpMBAJIlTlX7rnae/xT+kE6xifer00piKP1BEACAEIZeR+fPna/LkyZo+fbpWr16tQYMGKTs7W3v27Klz/uXLl+uGG27QzTffrM8//1wvvviiVq5cqQkTJhx3eMBNadqvF+L/Vzd635YkPV41Rv+v4rfaZtJcTgYA0SXkMvLggw9qwoQJGj9+vAYMGKDZs2crKSlJc+bMqXP+3NxcZWRk6Oc//7l69uyps88+W7feeusxR1OASDbS/lxv+H6tIfYmFZokja+4SzOrrlUlez4BIGQh/ctZUVGhVatWaerUqcFptm0rKytLubm5dS4zcuRI/epXv9Kbb76pUaNGac+ePfrXv/6l0aNH1/t7ysvLVV5eHnzu9/tDiQk0SeOODzG6zfO67vLOl8cy+tzpodsqJ2mHSQ17PgCIVSGNjOzbt0+BQECpqbX/4U1NTVV+fn6dy5x11ll67rnnNHbsWMXHxystLU0pKSmaNWtWvb9n5syZSklJCT7S07lvB9znU4X+EvcXTYn7pzyW0YtV39UPKmZQRADgOIX9PuVffPGF7rjjDk2bNk2rVq3SokWLtHXrVt122231LjN16lQVFhYGHzt27Ah3TKBBnXVA8+P/V5d5Plal8ehXlTfrrqpbVa54t6MBQNQLaTdNp06d5PF4VFBQUGt6QUGB0tLqPmhv5syZOuuss3TXXXdJkk4//XS1adNG55xzjn7/+9+rS5cuRy3j8/nk8/lCiQaEzUBri56M/5O6WN/ogGmr2yru1ArT3+1YABAzQhoZiY+P17Bhw5STkxOc5jiOcnJyNHLkyDqXKS0tlW3X/jUej0eSZIwJNS/QorLtlfpX/Ax1sb7RJqerLq/4HUUEAJpZyIf+T548WePGjdPw4cM1YsQIPfzwwyopKdH48eMlVe9i2bVrl5555hlJ0mWXXaYJEyboscceU3Z2tvLy8jRp0iSNGDFCXbt2bd61AUJwrANWx3kWa7r3GdmW0bLA6fpp5c9VxL1lAKDZhVxGxo4dq71792ratGnKz8/X4MGDtWjRouBBrXl5edq+fXtw/htvvFFFRUV69NFH9Ytf/ELt27fXBRdcoPvuu6/51gJoVkZ3e+frJ97XJEl/r8rSb6vGKSCPy7kAIDZZJgr2lfj9fqWkpKiwsFDJyclux0GMqGtkxKsq3Rf3hH7o+UCS9EDllZoVuFyS1cLpAKDlbL13TFhet7Gf31yhCaiRpEP6a9wjOs/zqaqMralVt+jFwHluxwKAmEcZQatwrOND2qlUc+Pv13D7K5WZeP2k8g695wxpoXQA0LpRRtDqpahYz8Tfq0H2ZhWaJN1YcY/WmD5uxwKAVoMyglatg/x6Nn6mBtjbtN+00/UVU/WFyXA7FgC0KpQRtFon6oD+Ef8H9bF3aY9pr2srfqWNppvbsQCg1aGMoFXqrAP6Z/zv1MvO127TQddW/FpbzNFXAwYAhB9lBDGpoQNWO6pQz8X/Qb3sfO1wTtTVlb/WTtO5BdMBAI4U9hvlAZEkRcX6e/y96mPv0m7TgSICABGAMoJWo61K9XT8vRpgb9Nek6JrKygiABAJKCNoFRJ1SHPiH9Bge7O+MW11bcWvOEYEACIEZQQxL05VejzuIY2wN8hvknR9xa/0lUl3OxYAoAYHsCLqNXx1VaN74/6m73rWqdT4NK7iHn3OdUQAIKIwMoKYdrd3vn7o+VBVxtZPKu/gyqoAEIEoI4hZN3gW6yfe1yRJU6omaKkz2N1AAIA6UUYQk0bZK/Rb7zOSpAcqr9S/Aue6nAgAUB/KCGLOUOsrPRz3V9mW0TNVF2lW4HK3IwEAGkAZQUzpqn16PP5B+axKvR0Ypt9WjZNkuR0LANAAzqZB1Knv7JkkHdKT8X/SiZZfXzg9NKlyohz6NgBEPP6lRkyw5OhPcY/VXF01WbdU/EKlSnA7FgCgESgjiAl3ev+lUZ5PVG68urVisnark9uRAACNRBlB1LvM/kg/9y6QJP2q8hatNn3dDQQACAllBFGtr7VD98U9IUmaXXWpXnK+63IiAECoKCOIWm1VqsfiHlaSVa73A6fp/qqr3I4EAGgCzqZBRKv/vjNG98X9Tb3tPO02HThzBgCiGP96Iyrd5FmkMZ6VqjAeTay4Q98o2e1IAIAmoowg6gyzNmiq9x+SpP+ruo6b3wFAlKOMIKp0kF+z4v+sOCug1wIj9XTgYrcjAQCOE2UEUcTo/rjHlWYd0Canq6ZUThCXegeA6EcZQdS43rNEWZ41Kjdxmlj5c66wCgAxgrNpEFHqO3umj7VTv/Y+J0maWXW1NpjuLRkLABBGjIwg4vlUoT/HPaoEq1JLA4M0L5DtdiQAQDOijCDi3e2dr/72du0zyfpl5W3iOBEAiC2UEUS0c+1PdbP3LUnSXZW3ap9SXE4EAGhulBFErBQV64G4xyVJ86ou1nvOEJcTAQDCgTKCiDUt7u/qbB3U104Xzay6xu04AIAwoYwgIp1vr9EPPR/IMZbuqrxV5Yp3OxIAIEw4tReuqe803mSVaGbck5KkJwOjtdr0bclYAIAWxsgIIs7/eJ9VmnVAm500/anqR27HAQCEGWUEEeVc+1Nd6V0mx1i6u/LH7J4BgFaAMoKI0U6lmhn3hCRpXiBb/zanuJwIANASKCOIGHd556ur9Y22Oql6oOpKt+MAAFoIZQQR4XTra13neUeSNLXqFpVxEzwAaDUoI3CdLUe/j5sj2zJ6JXCWcp2BbkcCALQgyghcd51niU63t8hvkvR/lde5HQcA0MK4zghaTF3XFTlRB/RL7wuSpPurxnLvGQBohRgZgat+E/eskq0yrXV66R+BC92OAwBwAWUErjnbXqfveXIVMJZ+XXmzHP46AkCrxL/+cEWcqjTDO0+S9EzgYn1uerobCADgGsoIXHGDZ7F623naa1L0IJd8B4BWjTKCFtdBft3hfUVS9UGrRUpyOREAwE2UEbS4O73/UrJVqv84GXop8F234wAAXEYZQYvqa+3QNZ4cSdLvKq/noFUAANcZQXjUdU0Ryeh/vM/KYxm9GRihFaZ/i+cCAEQe/luKFnO+vVbf9axTufFqZtXVbscBAEQIyghahFdV+h/vs5KkOYFR2mFSXU4EAIgUlBG0iOs879ScypusWVWXux0HABBBKCMIuzYq089qTuV9sOpHKuZUXgDAESgjCLubPW+po1Wkr50ueiFwnttxAAARhjKCsDpBfk3wVp9Z82DVjxSQx+VEAIBI06QyMmvWLGVkZCghIUGZmZlauXJlg/OXl5fr17/+tXr06CGfz6eMjAzNmTOnSYERXX7ifU3trDKtczL0pjPC7TgAgAgU8nVG5s+fr8mTJ2v27NnKzMzUww8/rOzsbG3YsEGdO3euc5krr7xSBQUFeuqpp3TyyScrLy9PjuMcd3hEti7arxs8SyRJD1SNlWEgDgBQh5DLyIMPPqgJEyZo/PjxkqTZs2dr4cKFmjNnjqZMmXLU/IsWLdKyZcu0efNmdejQQZKUkZFxfKkRceq6yNkd3pfksyqVGxig953TXUgFAIgGIf1XtaKiQqtWrVJWVtZ/X8C2lZWVpdzc3DqXee211zR8+HDdf//9Oumkk9S3b1/98pe/VFlZWb2/p7y8XH6/v9YD0aW3tUs/8iyTVH0zPMlyNxAAIGKFNDKyb98+BQIBpabWvmBVamqqvvzyyzqX2bx5sz788EMlJCTolVde0b59+/STn/xE+/fv19y5c+tcZubMmZoxY0Yo0RBhfuF9UR7L6O3AMK0xfdyOAwCIYGHfie84jizL0nPPPacRI0Zo9OjRevDBB/X000/XOzoydepUFRYWBh87duwId0w0owHWVo32rJRjLP2x6kq34wAAIlxIIyOdOnWSx+NRQUFBrekFBQVKS0urc5kuXbropJNOUkpKSnBa//79ZYzRzp071afP0f9r9vl88vl8oURDBDl8gbPXnZH6yqS7nAYAEOlCGhmJj4/XsGHDlJOTE5zmOI5ycnI0cuTIOpc566yztHv3bhUXFwenffXVV7JtW926dWtibESqftZ2jfJ8IsdY+kvVFW7HAQBEgZB300yePFlPPPGEnn76aa1fv1633367SkpKgmfXTJ06VTfccENw/muuuUYdO3bU+PHj9cUXX+j999/XXXfdpZtuukmJiYnNtyaICD/zLpAkvelkapOhbAIAji3kU3vHjh2rvXv3atq0acrPz9fgwYO1aNGi4EGteXl52r59e3D+tm3basmSJfrZz36m4cOHq2PHjrryyiv1+9//vvnWAhHhZGunRtsrJIlREQBAo1nGGON2iGPx+/1KSUlRYWGhkpOT3Y6DOmRMWahH4h7V5Z6P9FbgDN1eeafbkQAAjbT13jFhed3Gfn5zSUw0i97WLl1mV19r5i9V33c5DQAgmoS8mwao62qrf4p7VbZltCQwTF+YjJYPBQCIWoyM4Lj1sPJ1hb1ckvRnRkUAACGijOC4/cTzmjyW0buBwVpnerkdBwAQZSgjOC6ddUDf93wgSXqUM2gAAE1AGcFxGe9dpHgroE+cvlpt+rodBwAQhSgjaLK2KtW1nuqr8T5edZnLaQAA0Yoygia7yvOekq1SbXK6KscZ4nYcAECUooygSbyq0k3etyRJTwTGyPBXCQDQRHyCoEkus3PV1fpGe0x7LQic5XYcAEAUo4ygCYx+7H1DkjSvKlvlinc5DwAgmnEFVhzTt6+4eq79mfrbO1RifHo2cKFLqQAAsYKREYTsVs/rkqTnAxfIr7YupwEARDvKCEIy0NqqMz1fqMrYmlM1yu04AIAYQBlBSG70LJIkLXS+o93q5HIaAEAsoIyg0TrIr+95ciVVH7gKAEBzoIyg0a7yvCufVam1Ti+tMSe7HQcAECMoI2gUr6p0vfcdSdLTVdmSLHcDAQBiBmUEjZJt/1tdrG+01yRrofMdt+MAAGIIZQSNMs67WJL0j8CFqlCcy2kAALGEMoJjGmht1Qh7gyqNR89VZbkdBwAQYygjOKZxnupRkbecEdqjE1xOAwCINVwOHrV8+9LvHeRXru8jSZzOCwAID0ZG0KCrPO/JZ1XqU6eXVps+bscBAMQgygjqZcvRNd4cSdLTVReL03kBAOFAGUG9vmt/qm7WPh0wbTmdFwAQNpQR1Otaz7uSpJcC56hc8S6nAQDEKsoI6pSm/brAXi1Jej5wgctpAACxjDKCOl3pWSaPZbTCOUVfm5PcjgMAiGGUERzFlqOx3vckSc9VXehyGgBArKOM4Cjn2Wt1krVf35i2Wuyc4XYcAECMo4zgKNd4qk/nfSnwXQ5cBQCEHWUEtXTRfp1vr5XEgasAgJZBGUEtY73vyWMZ5QYGaLPp6nYcAEArwL1pWrkj70XjUUAf+pZKYlQEANByGBlB0Ln2p+pifaP9pp0WceAqAKCFUEYQ9CPPMknSK4GzVaE4l9MAAFoLyggkSSfIrwtrrrj6YuBcl9MAAFoTyggkSVd4liveCugzp6c2mO5uxwEAtCKUEUiS/p/nfUmMigAAWh5lBBpgbdVAe5vKjVevBc50Ow4AoJWhjCA4KrLEGaZCtXU5DQCgtaGMtHJxqtLlnuWSpH8FvutyGgBAa0QZaeUusFero1WkAtNeHzinux0HANAKUUZaucO7aF4OnKOAPC6nAQC0RpSR1qyoIHhTPHbRAADcwr1pWpEj70MjSRM8b+jXcY5WOyfra3OSS6kAAK0dIyOtlgnuovkX1xYBALiIMtJKDbC2qZ+9U+XGqzcC33E7DgCgFaOMtFKHT+fNcYbKrzYupwEAtGaUkVbIlqPveXIlSa8GznI5DQCgtaOMtEKZ9np1sb5RoUnSe85gt+MAAFo5ykgrdIVdvYtmYSBTFYpzOQ0AoLWjjLQyPlVolGeFJGlB4GyX0wAAQBlpdS6w1yjZKtMu01GfmH5uxwEAgDLS2lxRcxbNq4GzZNj8AIAIwKdRK5KiYp1vr5EkLeAsGgBAhKCMtCKjPSsUbwX0hdNDX5l0t+MAACCJMtKqHN5FsyBwpstJAAD4ryaVkVmzZikjI0MJCQnKzMzUypUrG7Xc8uXL5fV6NXjw4Kb8WoQoY8rC4OOsKfOUaX8px1h6jTICAIggIZeR+fPna/LkyZo+fbpWr16tQYMGKTs7W3v27GlwuYMHD+qGG27QhRde2OSwaLrLaq64+rHTX/nq6HIaAAD+K+Qy8uCDD2rChAkaP368BgwYoNmzZyspKUlz5sxpcLnbbrtN11xzjUaOHNnksGi6w2XkNYdREQBAZAmpjFRUVGjVqlXKysr67wvYtrKyspSbm1vvcnPnztXmzZs1ffr0Rv2e8vJy+f3+Wg80XU8rTwPtbao0Hi0KnOF2HAAAagmpjOzbt0+BQECpqam1pqempio/P7/OZTZu3KgpU6bo2WefldfrbdTvmTlzplJSUoKP9HTO/Dgel9rVRXG5c6oOqp3LaQAAqC2sZ9MEAgFdc801mjFjhvr27dvo5aZOnarCwsLgY8eOHWFMGfsu9XwsSXrD+Y7LSQAAOFrjhipqdOrUSR6PRwUFBbWmFxQUKC0t7aj5i4qK9O9//1tr1qzRT3/6U0mS4zgyxsjr9ertt9/WBRdccNRyPp9PPp8vlGioRx9rp/rZO1VuvHo7MNztOAAAHCWkkZH4+HgNGzZMOTk5wWmO4ygnJ6fOA1OTk5O1bt06rV27Nvi47bbb1K9fP61du1aZmZnHvwZo0KU1B66+7wySX21cTgMAwNFCGhmRpMmTJ2vcuHEaPny4RowYoYcfflglJSUaP368pOpdLLt27dIzzzwj27Z16qmn1lq+c+fOSkhIOGo6wsHosprjRV4PsIsGABCZQi4jY8eO1d69ezVt2jTl5+dr8ODBWrRoUfCg1ry8PG3fvr3ZgyJ0A61t6mXn65CJU44z1O04AADUyTLGGLdDHIvf71dKSooKCwuVnJzsdpyo8dj/XKfbva/rzcAI/aRykttxAAARauu9Y8Lyuo39/ObeNLHKGF1q15xFwy4aAEAEo4zEql2rlW7vVYnx6V1niNtpAACoF2UkVn3+siQpxxmqQ+I0aQBA5Ar5AFZEpowpC4PfW3K03PcPdbXYRQMAiHyMjMSgIdYmdbW+UZFJ1DJnkNtxAABoEGUkBl3i+USSlOMMUbniXU4DAEDDKCMxx2i0Z4Uk6a3ACJezAABwbJSRGHOqtUXdrH0qNT520QAAogJlJMaM9qyUJL3rDOYsGgBAVKCMxBSjUXb1LppF7KIBAEQJykgMOcXaoZ52gQ6ZOL3nDHY7DgAAjUIZiSGjag5cXeYMUokSXU4DAEDjUEZiyGi7+ngRzqIBAEQTykiMONnaqT72LlUYj3KcoW7HAQCg0SgjMWJUzajIB87pKlKSy2kAAGg8ykiMOHxK71sOu2gAANGFMhIL9n+t/vZ2VRqPlgSGuZ0GAICQcNfeKHTkHXol6XbPa7onTsp1BqhQbV1KBQBA0zAyEgOya3bRLGIXDQAgClFGolya9muwvVmOsfR2YLjbcQAACBllJMpd7Pm3JGmV6aN9SnE5DQAAoaOMRLlL7E8kSYsCZ7icBACApqGMRLH2KtII+0tJ0mKHMgIAiE6UkSiW5Vktr+XoC6eHdprObscBAKBJKCNRLLtmF81iDlwFAEQxykiUStIhfddeJ0laxC4aAEAUo4xEqXPtT+WzKrXVSdUGk+52HAAAmowyEqUu8dTsonGGS7LcDQMAwHGgjEShOFXpfHuNJGkxp/QCAKIcZSQKnWl/rmSrTHtMe60xJ7sdBwCA40IZiULZdvW9aN4ODJNhEwIAohyfZNHGCegizypJXOgMABAbvG4HwLFlTFkY/H6YtUEv+fzymyR97AxwMRUAAM2DkZEoc/jGeDnOEFXSJQEAMYAyElWMsu3qMvI2V10FAMQIykgU6WPtUoZdoHITp/ed092OAwBAs6CMRJGLa0ZFPnROVYkSXU4DAEDzoIxEkcPHi7ztsIsGABA7KCNRoov2a5C9WY6xlBMY6nYcAACaDWUkSmTVXFtklemjfUpxOQ0AAM2HMhIlLuYsGgBAjKKMRIFkFes79npJ0hJnmMtpAABoXpSRKHC+vVZxVkAbnG7aarq4HQcAgGZFGYkCnEUDAIhllJFIV1mm8+xPJXG8CAAgNlFGIt3mZWpjlWu36aB1pqfbaQAAaHbcaS3CHHmHXkm61/s3XeU9PCpiuRMKAIAwYmQkgtlylOVZLYmzaAAAsYsyEsGGWl+pk+VXoUnSCqe/23EAAAgLykgEu7jmqqs5zlBVsUcNABCjKCMRy3DVVQBAq0AZiVB9rZ3KsAtUbuL0vnO623EAAAgbykiEOjwq8oFzqkqV4HIaAADChzISobI9n0jiqqsAgNhHGYlAXbVPp9lb5RhLOYGhbscBACCsKCMR6KKas2j+bfpqv1JcTgMAQHhRRiLQ4eNFFnMWDQCgFaCMRJgUFSvTXi9JWsLxIgCAVoAyEmEusNfIazla76Rru0l1Ow4AAGFHGYkwF3tqLnTGqAgAoJVoUhmZNWuWMjIylJCQoMzMTK1cubLeeV9++WVddNFFOvHEE5WcnKyRI0dq8eLFTQ4c0yrLdK79mSSuugoAaD1CvuHJ/PnzNXnyZM2ePVuZmZl6+OGHlZ2drQ0bNqhz585Hzf/+++/roosu0h/+8Ae1b99ec+fO1WWXXaYVK1ZoyJAhzbIS0SxjysLg91n2Kj0ZX65dpqM+NxnuhQIAoAWFPDLy4IMPasKECRo/frwGDBig2bNnKykpSXPmzKlz/ocfflh33323zjjjDPXp00d/+MMf1KdPH73++uvHHT7WXFJzobPFgTMkWe6GAQCghYRURioqKrRq1SplZWX99wVsW1lZWcrNzW3UaziOo6KiInXo0KHeecrLy+X3+2s9Yp1XVcqyq68vsihwhstpAABoOSGVkX379ikQCCg1tfZZHqmpqcrPz2/Ua/zxj39UcXGxrrzyynrnmTlzplJSUoKP9PT0UGJGpUx7vdpbJdpnkvVv08/tOAAAtJgWPZvmH//4h2bMmKEXXnihzuNLDps6daoKCwuDjx07drRgSndcYlfvolkSGCaHk5wAAK1ISAewdurUSR6PRwUFBbWmFxQUKC0trcFl//nPf+qWW27Riy++WGs3T118Pp98Pl8o0aKaJUfZNaf0LnbYRQMAaF1C+i94fHy8hg0bppycnOA0x3GUk5OjkSNH1rvc888/r/Hjx+v555/XmDFjmp42Rg2xNqmzdVB+k6iPnIFuxwEAoEWFfGrv5MmTNW7cOA0fPlwjRozQww8/rJKSEo0fP15S9S6WXbt26ZlnnpFUvWtm3LhxeuSRR5SZmRk8tiQxMVEpKdwETpKya86iedcZogrFuZwGAICWFXIZGTt2rPbu3atp06YpPz9fgwcP1qJFi4IHtebl5Wn79u3B+f/2t7+pqqpKEydO1MSJE4PTx40bp3nz5h3/GkQ9EzxeZFFghMtZAABoeZYxxrgd4lj8fr9SUlJUWFio5ORkt+M0q1FT/6q3fFN1yMRpSPnjKlOC25EAAK3M1nvDcwhFYz+/OW3DZYcvdLbMGUQRAQC0SpQRl2UHd9FwFg0AoHWijLhp3yadYu9QpfEox+E+PQCA1oky4qYvq+/Pk+sMkF9tXQ4DAIA7Qj6bBk135B16JWlB/DMabEtvO8NdSgQAgPsYGXFJN2uvBtubFTAWp/QCAFo1yohLRtsfS5JWOP21T1z8DQDQelFGXDLGs0KStND5jstJAABwF2XEBelWgQYFd9FwSi8AoHWjjLhgjF09KpLrDNB+dtEAAFo5yogLxniqjxdhFw0AAJSRFtfDytdp9lZVGVuL2UUDAABlpKUd3kXzkTNQ3yi2bvoHAEBTUEZaGLtoAACojTLSgjKsPA20t9XsouGqqwAASJSRFjW6ZhfNcudUHVQ7l9MAABAZKCMt6NKaC529wS4aAACCuFFeGB15Y7xe1m6969umSuPR2+yiAQAgiJGRFnK55yNJ0ofOqSpUW5fTAAAQOSgjLcLo+/YHkqRXAme7nAUAgMhCGWkBw60N6m7vVbFJ0NsOu2gAADgSZaQF/MDzoSTprcAIHZLP5TQAAEQWykiY+VQRvNDZy845LqcBACDyUEbC7AJ7jVKsUu0yHfWx09/tOAAARBzKSJgd3kXzauAsGf64AQA4Cp+OYXSC/DrPXitJepmzaAAAqBNlJIwu9XysOCugdU6GNplubscBACAiUUbC6Iee6muLvBzgwFUAAOpDGQmXfRs12P5aVcbWa4Ez3U4DAEDEooyEy2fzJUnvO6drv1JcDgMAQOTiRnnN5Mib4tly9L5vrrpZXP4dAIBjYWQkDM6z16qbtU8HTFsu/w4AwDFQRsLges8SSdILgXNVrniX0wAAENkoI82su1Wgc+3PJEnPBbJcTgMAQOSjjDSzazw5si2jpYFB2m5S3Y4DAEDEo4w0I58qNNazVJL0d0ZFAABoFMpIMxpjf6wTrGLtNJ30njPE7TgAAEQFykgzus77jiTpH1UXyuGPFgCARuETs5kMtLZoqL1JFcaj+YHz3I4DAEDUoIw0k+s81aMibzqZXHEVAIAQUEaaQ9lBXeFZLkn6e9VFLocBACC6cDn4Jjry8u8/9ryuX8VVaL2TrlWmr4upAACIPoyMHKdEHdKPvdXF5KnAaEmWu4EAAIgylJHjdJ3nHXWy/NrqpHJTPAAAmoAychwSdUi3et+QJD0auEIBeVxOBABA9KGMHIfrPUsYFQEA4DhRRpqIUREAAJoHZaSJrvcsUUeriFERAACOE2WkKSpKGBUBAKCZUEaa4pMnGRUBAKCZcNGzRjjyAmcd5Nc7vvvVwWJUBACA5sDISIj+N26uOljF+tJJZ1QEAIBmQBkJwWj7Y13qWaEqY+uXlbcyKgIAQDOgjDRSRxXqd3FzJUmzApfrP6aXy4kAAIgNlJFG+t+4uepoFWm9012PVn3f7TgAAMQMykgjjLE/1hjPSlUaj35ZeZsqOe4XAIBmQxk5luI9+l3cHEnSXwOX63OT4W4eAABiDGWkIUUF0t9/oA5Wsb5weujRqivcTgQAQMxhf0MdMqYsVDdrj56Nm6kMu0B7TYp+XjmR3TMAAIRBk0ZGZs2apYyMDCUkJCgzM1MrV65scP6lS5dq6NCh8vl8OvnkkzVv3rym/NoW08faqX/Fz1CGXaAdzon6fxXTtcl0czsWAAAxKeQyMn/+fE2ePFnTp0/X6tWrNWjQIGVnZ2vPnj11zr9lyxaNGTNG559/vtauXatJkybplltu0eLFi487fFjs+EQvxP+v0qwD2uB00/+rmK5tJs3tVAAAxCzLGGNCWSAzM1NnnHGGHn30UUmS4zhKT0/Xz372M02ZMuWo+e+55x4tXLhQ//nPf4LTrrrqKh08eFCLFi1q1O/0+/1KSUlRYWGhkpOTQ4l7bMZIezdIX74urX9DylsrSVrjnKwbK+5Wodo27+8DACDCbL13TFhet7Gf3yEdBFFRUaFVq1Zp6tSpwWm2bSsrK0u5ubl1LpObm6usrKxa07KzszVp0qR6f095ebnKy8uDzwsLCyVVr1SzWv5nbVk6Tz3t/47qOMbSO85Q/bryJpXJllTavL8TAIAI0+yfr9963WONe4RURvbt26dAIKDU1NRa01NTU/Xll1/WuUx+fn6d8/v9fpWVlSkxMfGoZWbOnKkZM2YcNT09PT2UuMdhac0DAIDYl/JweF+/qKhIKSkp9f48Ik8PmTp1qiZPnhx87jiOvvnmG3Xs2FGWZTXb7/H7/UpPT9eOHTuaf/dPhIj1dWT9ol+sr2Osr58U++vI+jWdMUZFRUXq2rVrg/OFVEY6deokj8ejgoKCWtMLCgqUllb3QZ5paWl1zp+cnFznqIgk+Xw++Xy+WtPat28fStSQJCcnx+RfsCPF+jqyftEv1tcx1tdPiv11ZP2apqERkcNCOpsmPj5ew4YNU05OTnCa4zjKycnRyJEj61xm5MiRteaXpCVLltQ7PwAAaF1CPrV38uTJeuKJJ/T0009r/fr1uv3221VSUqLx48dLqt7FcsMNNwTnv+2227R582bdfffd+vLLL/XXv/5VL7zwgu68887mWwsAABC1Qj5mZOzYsdq7d6+mTZum/Px8DR48WIsWLQoepJqXl6ft27cH5+/Zs6cWLlyoO++8U4888oi6deumJ598UtnZ2c23Fk3k8/k0ffr0o3YJxZJYX0fWL/rF+jrG+vpJsb+OrF/4hXydEQAAgObEjfIAAICrKCMAAMBVlBEAAOAqyggAAHBVzJeRWbNmKSMjQwkJCcrMzNTKlSsbnH/p0qUaOnSofD6fTj75ZM2bN69lgjbBzJkzdcYZZ6hdu3bq3LmzrrjiCm3YsKHBZZYuXSrLso565Ofnt1Dqxvvtb397VM5TTjmlwWWiaftJUkZGRp3bY+LEiXXOH+nb7/3339dll12mrl27yrIsLViwoNbPjTGaNm2aunTposTERGVlZWnjxo3HfN0XX3xRp5xyihISEnTaaafpzTffDNMaNKyh9ausrNQ999yj0047TW3atFHXrl11ww03aPfu3Q2+5rx5847angkJCWFek/odaxveeOONR+W95JJLjvm60bANJdX5/rIsSw888EC9rxlJ27AxnwuR+D6M6TIyf/58TZ48WdOnT9fq1as1aNAgZWdna8+ePXXOv2XLFo0ZM0bnn3++1q5dq0mTJumWW27R4sWLWzh54yxbtkwTJ07Uxx9/rCVLlqiyslIXX3yxSkpKjrnshg0blJeXF3x07ty5BRKHbuDAgbVyfvjhh/XOG23bT5I++eSTWuu3ZMkSSdKPfvSjBpeL1O1XUlKiQYMGadasWXX+/P7779ef//xnzZ49WytWrFCbNm2UnZ2tQ4cO1fuaH330ka6++mrdfPPNWrNmja644gpdccUVte4E3lIaWr/S0lKtXr1av/nNb7R69Wq9/PLL2rBhg773ve8d83WTk5Nrbc9t27aFI36jHGsbStIll1xSK+/zzz/f4GtGyzaUVGu98vLyNGfOHFmWpR/+8IcNvm6kbMPGfC5E5PvQxLARI0aYiRMnBp8HAgHTtWtXM3PmzDrnv/vuu83AgQNrTRs7dqzJzs4Oa87msmfPHiPJLFu2rN553nvvPSPJHDhwoOWCNdH06dPNoEGDGj1/tG8/Y4y54447TO/evY3jOHX+PJq2nyTzyiuvBJ87jmPS0tLMAw88EJx28OBB4/P5zPPPP1/v61x55ZVmzJgxtaZlZmaaW2+9tdkzh+Lb61eXlStXGklm27Zt9c4zd+5ck5KS0rzhmkld6zhu3Dhz+eWXh/Q60bwNL7/8cnPBBRc0OE8kb8Nvfy5E6vswZkdGKioqtGrVKmVlZQWn2batrKws5ebm1rlMbm5urfklKTs7u975I01hYaEkqUOHDsecd/DgwerSpYsuuugiLV++PNzRmmzjxo3q2rWrevXqpWuvvbbWBfW+Ldq3X0VFhZ599lnddNNNx7whZLRsvyNt2bJF+fn5tbZRSkqKMjMzG9xG0bxdCwsLZVnWMe+tVVxcrB49eig9PV2XX365Pv/885YJ2ERLly5V586d1a9fP91+++3av39/g/NH6zYsKCjQwoULdfPNNx9z3kjdht/+XIjU92HMlpF9+/YpEAgErwx7WGpqar371/Pz8+uc3+/3q6ysLGxZm4PjOJo0aZLOOussnXrqqfXO16VLF82ePVsvvfSSXnrpJaWnp+u8887T6tWrWzBt42RmZmrevHlatGiRHnvsMW3ZskXnnHOOioqK6pw/mrefJC1YsEAHDx7UjTfeWO880bT9vu3w+y6U9+Th5UJdJhIcOnRI99xzj66++uoGbz7Wr18/zZkzR6+++qqeffZZOY6jM888Uzt37mzBtI13ySWX6JlnnlFOTo7uu+8+LVu2TKNGjVIgEKh3mWjdhk8//bTatWunH/zgBw3OF6nbsK7PhUh9H4Z8OXhEpokTJ+o///lPg8dUSNVvmn79+gWfn3nmmfr666/10EMP6e9//3u4Y4Zk1KhRwe9PP/10ZWZmqkePHnrhhRca9T+VaPPUU09p1KhRDd5qO5q2X2tWWVmpK6+8UsYYPfbYYw3OO3LkyFo3Dj3zzDPVv39/Pf744/rd734X7qghu+qqq4Lfn3baaTr99NPVu3dvLV26VBdeeKGLyZrfnDlzdO211x7zYNRI3YaN/VyIBDE7MtKpUyd5PB4VFBTUml5QUKC0tLQ6l0lLS6tz/uTkZCUmJoYt6/H66U9/qjfeeEPvvfeeunXrFvLyI0aM0KZNm8KQrHm1b99effv2rTdrtG4/Sdq2bZveeecd3XLLLSEvGy3b7/D7LpT35OHlQl3GTYeLyLZt27RkyZKQb8keFxenIUOGRMU2laRevXqpU6dODeaNtm0oSR988IE2bNjQpPdkJGzD+j4XIvV9GLNlJD4+XsOGDVNOTk5wmuM4ysnJqdVgjzRy5Mha80vSkiVL6p3fbcYY/fSnP9Urr7yid999Vz179mzS66xdu1ZdunRp5nTNr7i4WF9//XW9WaNt+x1p7ty56ty5s8aMGRPystGy/Xr27Km0tLRa28jv92vFihUNbqNo2q6Hi8jGjRv1zjvvqGPHjiG/RiAQ0Lp166Jim0rSzp07tX///gbzRtM2POypp57SsGHDNGjQoJCXdXMbHutzIWLfh81yGGyE+uc//2l8Pp+ZN2+e+eKLL8yPf/xj0759e5Ofn2+MMWbKlCnm+uuvD86/efNmk5SUZO666y6zfv16M2vWLOPxeMyiRYvcWoUG3X777SYlJcUsXbrU5OXlBR+lpaXBeb69jg899JBZsGCB2bhxo1m3bp254447jG3b5p133nFjFRr0i1/8wixdutRs2bLFLF++3GRlZZlOnTqZPXv2GGOif/sdFggETPfu3c0999xz1M+ibfsVFRWZNWvWmDVr1hhJ5sEHHzRr1qwJnk1y7733mvbt25tXX33VfPbZZ+byyy83PXv2NGVlZcHXuP76682UKVOCz5cvX268Xq/54x//aNavX2+mT59u4uLizLp16yJq/SoqKsz3vvc9061bN7N27dpa78ny8vJ612/GjBlm8eLF5uuvvzarVq0yV111lUlISDCff/55i6+fMQ2vY1FRkfnlL39pcnNzzZYtW8w777xjhg4davr06WMOHToUfI1o3YaHFRYWmqSkJPPYY4/V+RqRvA0b87kQie/DmC4jxhjzl7/8xXTv3t3Ex8ebESNGmI8//jj4s3Hjxplzzz231vzvvfeeGTx4sImPjze9evUyc+fObdnAIZBU5+PIzN9ex/vuu8/07t3bJCQkmA4dOpjzzjvPvPvuuy0fvhHGjh1runTpYuLj481JJ51kxo4dazZt2hT8ebRvv8MWL15sJJkNGzYc9bNo236HTz3+9mPcuHHGmOrTCn/zm9+Y1NRU4/P5zIUXXnjUep977rnB+Q974YUXTN++fU18fLwZOHCgWbhwYQutUW0Nrd+WLVvqfU++9957wdf49vpNmjQp+G9UamqqGT16tFm9enXLr1yNhtaxtLTUXHzxxebEE080cXFxpkePHmbChAnB/+AdFq3b8LDHH3/cJCYmmoMHD9b5GpG8DRvzuRCJ70OrJjwAAIArYvaYEQAAEB0oIwAAwFWUEQAA4CrKCAAAcBVlBAAAuIoyAgAAXEUZAQAArqKMAAAAV1FGAACAqygjAADAVZQRAADgKsoIAABw1f8HkkSWQ6ngs2UAAAAASUVORK5CYII="
},
"metadata": {}
}
],
"execution_count": 8
},
{
"cell_type": "code",
"source": "N = 1000000\n\na = np.linspace(0.0, np.pi*0.25, N)\nr = 0.5 / np.cos(a) # the radius of a cube of unit length\n\nself_dispersal = np.mean(cdf(r, u=u, p=1.0))\n\nx, y = sample_dispersal(u, N)\n\nexperimental = np.mean((np.abs(x) <= 0.5) & (np.abs(y) <= 0.5))\n\nself_dispersal, experimental",
"metadata": {
"trusted": true
},
"outputs": [
{
"execution_count": 9,
"output_type": "execute_result",
"data": {
"text/plain": "(0.03080578323898754, 0.030975)"
},
"metadata": {}
}
],
"execution_count": 9
},
{
"cell_type": "code",
"source": "",
"metadata": {
"trusted": true
},
"outputs": [],
"execution_count": null
}
]
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment