From 9894c41ff072ad3b2b9651f393cfb429440fd2a9 Mon Sep 17 00:00:00 2001 From: Charles Weiss <69219836+weisscharlesj@users.noreply.github.com> Date: Sat, 24 Jun 2023 15:52:43 -0500 Subject: [PATCH] Added stochastic tutorial ipynb file Added tutorial-stochastic-simulation.ipynb file which demonstrates using random numbers to perform stochastic simulations of various processes. This tutorial was proposed in issue #184. --- tutorial-stochastic-simulations.ipynb | 777 ++++++++++++++++++++++++++ 1 file changed, 777 insertions(+) create mode 100644 tutorial-stochastic-simulations.ipynb diff --git a/tutorial-stochastic-simulations.ipynb b/tutorial-stochastic-simulations.ipynb new file mode 100644 index 00000000..2bfd1c8c --- /dev/null +++ b/tutorial-stochastic-simulations.ipynb @@ -0,0 +1,777 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "2931e270-7884-4521-928c-65fa5ca36d60", + "metadata": {}, + "source": [ + "# Simulating Using Random Variables\n", + "\n", + "Randomness, or apparent randomness, is a major part of the world we live from the movement of gas molecules in the air we breath to the outcome of a coin toss. While many events are not truely random, they are often close enough to random or seemingly random enough to be simulated as such. Despite many processes being random, their outcomes often conform to a statistical pattern if enough of these random events occur. For example, if we roll a six-sided die over a large number of rolls, we expect to roll a 5 about $\\frac{1}{6}$ of the time. As a general rule, the larger the number of random events observed, the closer the distribution of these events will *likely* be to the expected disbribution. In this tutorial, we will use a random number generator provided in NumPy to solve a number of problems through stochastic simulations. While the problems demonstrated below all have know solutions or equations, you may face other problems that do not have an analytical model or known solution. The methodology demonstrated below can be used to solve problems when no known solution or model exists. \n", + "\n", + "## What you'll do\n", + "- Generate random values using NumPy to power stochastic simulations\n", + "- Perform stochastic simulations to determine the\n", + " - Likelyhood of multiple children in a classroom of 25 having the same birthday\n", + " - Distribution of marbles in a Galton board\n", + " - Value of $\\pi$\n", + " - Amount of radioactive material left versus time\n", + "- Compare the above simulations to known values or an analytical model\n", + "\n", + "## What you'll learn:\n", + "- Generating large quantites of random values in a desired distribution\n", + "- Using random values to simulation a range of scenarios and solve problems\n", + "\n", + "## What you'll need:\n", + "- [NumPy](https://numpy.org/)\n", + "- [Matplotlib](https://matplotlib.org/)\n", + "- Python's `math` module\n", + "\n", + "The above two packages are imported using the following commands:" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "id": "50324dde-32fe-44b1-b15f-bce368ce16ca", + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np\n", + "import matplotlib.pyplot as plt\n", + "import math" + ] + }, + { + "cell_type": "markdown", + "id": "70b2a718-ba0c-4ed3-8a5d-4a7bd90f45c8", + "metadata": {}, + "source": [ + "## Generating random variables\n", + "\n", + "NumPy can generate random variables using a Generator from the `np.random` module. While these values are not truely random, they are close enough for the following applications. One of the advantages of using NumPy to generate variables over the standard `random` Python module is that NumPy can generate a large number of random variables efficiently without the use of a `for` loop. These variables can be generated in a variety of probability distribution including normal, poisson, and binomial; but for the following simulations, we will mostly use even distributions where all numbers in the range have an equal likelyhood of being generated. This is accomplished by first creating a Generator using the `np.random.default_rng()` function and then using `random(n)` or `np.integers(low, high=, size=n)` methods with the Generator to produce random values where `n` is the number of random values returned. Below are some common methods for generating random values.\n", + "\n", + "| Method | Description |\n", + "|-----------|-------------|\n", + "|`random(n)`| Generates `n` float values in the [0,1) range with an even distribution |\n", + "|`integers(low, high=, size=)`| Generates `size` integers in the requested range [low, high) in a uniform distribution|\n", + "|`normal(loc, scale=, size=)`| Generates `size` floats in a normal distribution centered on `loc` with a `scale` standard deviation|\n", + "|`binomial(t, p, size=)`| Generates `size` integers in a bionimial distribution with `t` trials for each value with a probability `p`|\n", + "|`choice(a, size=)`| Selects `size` elements from `a` list, tuple or ndarray with an even probability |\n", + "|`shuffle(a)`| Shuffles items in list, tubple, or ndarray `a` |" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "id": "8a3c0c28-f9bf-466d-8522-a38e11488499", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "array([0.39930577, 0.71741468, 0.28082333, 0.08272492, 0.96977132,\n", + " 0.56392675, 0.64432096, 0.57687146, 0.47536097, 0.12240199])" + ] + }, + "execution_count": 2, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "rng = np.random.default_rng(seed=18)\n", + "rng.random(10)" + ] + }, + { + "cell_type": "markdown", + "id": "a59b2279-3f9f-4874-9754-aa1e031634cf", + "metadata": {}, + "source": [ + "You may have guessed by looking at the returned values that `random()` produces values in the 0 $\\rightarrow$ 1 range, but what happens if we need values in a different range? We can modify these values by mutiplying them by a coefficient to increase the range and using addition or substration to shift the center of the range. For example, if we need random values from -10 $\\rightarrow$ +10, we can accomplish this by the following." + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "id": "f6cbfb8d-b153-494b-9532-c657b8fe461b", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "array([-3.72807758, 4.72417165, 8.14778066, 7.7720389 , 8.96641764,\n", + " -9.49192157, 4.7607717 , 3.45200153, 2.50375281, 2.75248376])" + ] + }, + "execution_count": 3, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "20 * rng.random(10) - 10" + ] + }, + { + "cell_type": "markdown", + "id": "cdc9a16c-6d94-449c-8a41-1dbbb30b649a", + "metadata": {}, + "source": [ + "It's worth noting that the `random()` method returns values from the [0,1) range which includes the lower end and excludes the upper end. This is close enough to [0,1] for our purposes." + ] + }, + { + "cell_type": "markdown", + "id": "b915dba8-48d1-4dde-b3bd-4f048502f78e", + "metadata": {}, + "source": [ + "## Solving the birthday problem\n", + "\n", + "The goal of the [birthday problem](https://en.wikipedia.org/wiki/Birthday_problem) is to determine the probability of multiple children in a classroom of 25 students having the same birthday. For simiplicity, we will ignore leap years (i.e., February 29 birthdays) and assume that the students are equilly likely to be born any day of the year. One method of determining the answer is to repeately simulate a classroom of 25 student birthdays by using a random number generator and see how often all students in the class have unique birthdays. The general steps for this processes are outlined below. \n", + "\n", + "1. Simulate 25 random birthdays with a random integer generator\n", + "2. Check to see if any of the above birthdays match\n", + "3. Record the outcome if all birthdays are unique\n", + "4. Repeat setps 1-3 a large number of times\n", + "\n", + "To make things easier, we will just generate integers 0 $\\rightarrow$ 364 representating days of the year instead of actual months and days." + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "id": "36e6f902-5218-4975-a625-7a54a05fe03b", + "metadata": {}, + "outputs": [], + "source": [ + "def all_unique(class_size, n_classrooms=100):\n", + " '''Randomly generates birthdays for groups of kids of a given group size and outputs the percent of time\n", + " all students have uniquie birthdays. All input and output values are ints or floats.\n", + " \n", + " (classroom size, number of classrooms) -> percentage classrooms with all unique birthdays\n", + " '''\n", + " \n", + " all_unique_class = 0 # number of classrooms with students NOT sharing birthdays\n", + " \n", + " for classroom in range(n_classrooms):\n", + " birthdays = rng.integers(0, high=365, size=class_size)\n", + " if class_size == np.unique(birthdays).size:\n", + " all_unique_class += 1\n", + " \n", + " return 100 * all_unique_class / n_classrooms" + ] + }, + { + "cell_type": "markdown", + "id": "45763a5f-887f-4a4b-a5d2-a6a04cb52fcc", + "metadata": {}, + "source": [ + "Below we simulate 50000 classrooms of 25 students and find that only around 43% of the time every student in a classroom has a unique birthday from all the other students." + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "id": "7c77811a-d825-4f52-8161-f037427c72f9", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "43.286" + ] + }, + "execution_count": 5, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "all_unique(25, n_classrooms=50000)" + ] + }, + { + "cell_type": "markdown", + "id": "dc951c53-4e08-41e0-b54b-1b88e0ad1061", + "metadata": {}, + "source": [ + "This problem also has an equation for the probability given below where **p** is the probability of all students having a unique birthday and **n** is the classroom size. Performing this calculation below, we find that 43.13% of the time, students should have all unique birthdays which is quite close to what was found in the simulation above.\n", + "\n", + "$$ p = \\frac {365!}{365^n(365-n)!} $$" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "id": "c95b4fcc-45a0-4925-91fb-39d3469615b9", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "0.4313002960305361" + ] + }, + "execution_count": 6, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "n = 25\n", + "math.factorial(365)/(365**n * math.factorial(365-n))" + ] + }, + { + "cell_type": "markdown", + "id": "5a18ddfb-e49f-44b5-ac98-c1e9512e15ad", + "metadata": {}, + "source": [ + "## Simulating a Galton board\n", + "\n", + "Next, we will simulate a [Galton board](https://en.wikipedia.org/wiki/Galton_board) which allows mables to fall down levels of staggered pegs. Each marble starts in the horizontal center of the board, and at each level, the marble hits a peg. To continue falling to the next level, the marble must move around the peg, and it has an equal likeliness of moving to the left as it does to the right. At the bottom of the board, the marbles collect in bins with the height of marbles yielding a normal distribution." + ] + }, + { + "cell_type": "markdown", + "id": "2a0bb53d-71bb-4d70-a308-2f49f1d099ea", + "metadata": {}, + "source": [ + "![](tutorial-stochastic/galton_board.png)" + ] + }, + { + "cell_type": "markdown", + "id": "745bf65b-30c0-4ffb-a4ba-5ccae89d9805", + "metadata": {}, + "source": [ + "**Figure 1.** An example of a Galton board with three marbles with three different potential paths (dotted lines) shown.\n", + "\n", + "The key operational detail of a Galton board is that every marble moves one increment to the left or right at each level, so our simulation centers around a series of randomly generated +1 and -1 to determine this movement by the following methodolog.\n", + "1. Generate a random series of +1 and -1 in an even distribution - one value for every level on the Galton board\n", + "2. Sum up the values in step 1 to determine the final location of the mable\n", + "3. Repeat steps 1-2 for each simulated marble\n", + "4. Visualize the result\n", + "\n", + "Below is the simulation for a single marble moving down a 10-layered Galton board. Note how the integers 0 and 1 generated by `integers()` and converted to a series of -1 and +1 values, respectively." + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "id": "32ce5473-354f-494a-85eb-68054974a2ea", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "4" + ] + }, + "execution_count": 7, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "moves = 2 * rng.integers(0, high=2, size=10) - 1\n", + "final_position = moves.sum()\n", + "final_position" + ] + }, + { + "cell_type": "markdown", + "id": "ec34bc52-8270-42bb-8503-cbddfe2d8b25", + "metadata": {}, + "source": [ + "To simulate additional marbles, we could use a `for` loop, but because we need the same number of integers for each marble, a more efficient way is to generate the random integers in a two-dimensional array with a row for every level on the Galton board and column for every marble simulated. We then sum up the values in each column to obtain the final positions of all marbles." + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "id": "5b9fa1dd-d9db-4511-98e3-faa61a0cbff0", + "metadata": {}, + "outputs": [], + "source": [ + "def galton(marbles, levels=10):\n", + " '''integer -> ndarray of integers\n", + "\n", + " This function takes in an integer number of marbles and number of levels and outputs \n", + " the horizontal position of these marbles at the bottom of a Galton board. \n", + " '''\n", + " moves = 2 * rng.integers(0, high=2, size=(levels, marbles)) - 1\n", + " final_position = moves.sum(axis=0)\n", + " \n", + " return final_position" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "id": "d2b772e2-61f8-46f6-90e7-94a553837ac6", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "array([ 0, -2, -2, 2, 4, -4, -4, -2, 2, 2, 0, -2, 4, 0, 4, -4, 0,\n", + " -8, 2, 2])" + ] + }, + "execution_count": 9, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "sim = galton(1200)\n", + "sim[:20]" + ] + }, + { + "cell_type": "markdown", + "id": "693e2121-9fe7-428b-b284-654654a990ac", + "metadata": {}, + "source": [ + "We visualize the results below using a histogram plot. Being that all values are integers, we set `align='left'` so that each bar resides directly over the integer it represents." + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "id": "8c5d4ff6-4238-4283-9b49-e2fdfc52dcab", + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjsAAAGwCAYAAABPSaTdAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAAA4T0lEQVR4nO3deVgW9eL//9cNKIKyiAqIC+KSey5YhFrHjiailaaf1LJySzuJiUuWlrunY5mZaaR1KrCOW37SOnkKI0w9Kq6Je7iE4gJYbrh8BYT5/dHP+9N9QOXG+/bGOc/Hdc11OTPve+Y1LPFqZu65LYZhGAIAADApN1cHAAAAcCbKDgAAMDXKDgAAMDXKDgAAMDXKDgAAMDXKDgAAMDXKDgAAMDUPVwcoCwoLC3Xq1Cn5+PjIYrG4Og4AACgBwzB08eJFhYSEyM3txudvKDuSTp06pVq1ark6BgAAKIXjx4+rZs2aN1xP2ZHk4+Mj6fcvlq+vr4vTAACAksjJyVGtWrWsf8dvhLIjWS9d+fr6UnYAALjL3OoWFG5QBgAApkbZAQAApkbZAQAApkbZAQAApkbZAQAApkbZAQAApkbZAQAApkbZAQAApkbZAQAApkbZAQAApkbZAQAApkbZAQAApkbZAQAApkbZAQAApkbZAQAApubh6gAAUBp1xv3L4ds8+mY3h28TgOtxZgcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJgaZQcAAJiah6sDAIDZ1Rn3L4dv8+ib3Ry+TcCsXHpmZ8aMGbrvvvvk4+OjwMBA9ejRQ2lpaTZjOnToIIvFYjP95S9/sRmTkZGhbt26ydvbW4GBgRo7dqyuXbt2Jw8FAACUUS49s7Nu3TrFxMTovvvu07Vr1/Taa6+pc+fO2r9/vypWrGgdN2TIEE2bNs067+3tbf13QUGBunXrpuDgYG3atEmZmZl67rnnVK5cOf3tb3+7o8cDAADKHpeWncTERJv5hIQEBQYGaseOHXrooYesy729vRUcHFzsNr7//nvt379fP/zwg4KCgtSyZUtNnz5dr776qqZMmaLy5cs79RgAAEDZVqZuUL5w4YIkKSAgwGb5okWLVLVqVTVr1kzjx4/XlStXrOtSUlLUvHlzBQUFWZdFRUUpJydH+/btK3Y/ubm5ysnJsZkAAIA5lZkblAsLCzVy5Ei1a9dOzZo1sy5/+umnFRoaqpCQEO3evVuvvvqq0tLStGLFCklSVlaWTdGRZJ3Pysoqdl8zZszQ1KlTnXQkAACgLCkzZScmJkZ79+7Vhg0bbJYPHTrU+u/mzZurevXq6tixo44cOaJ69eqVal/jx4/X6NGjrfM5OTmqVatW6YIDAIAyrUxcxho+fLhWrVqlH3/8UTVr1rzp2IiICEnS4cOHJUnBwcHKzs62GXN9/kb3+Xh6esrX19dmAgAA5uTSsmMYhoYPH66VK1dqzZo1CgsLu+VrUlNTJUnVq1eXJEVGRmrPnj06ffq0dUxSUpJ8fX3VpEkTp+QGAAB3D5dexoqJidHixYv19ddfy8fHx3qPjZ+fn7y8vHTkyBEtXrxYXbt2VZUqVbR7926NGjVKDz30kO69915JUufOndWkSRM9++yzmjlzprKysjRhwgTFxMTI09PTlYcHAADKAJeWnfnz50v6/cGBfxQfH68BAwaofPny+uGHHzRnzhxdvnxZtWrVUq9evTRhwgTrWHd3d61atUovvviiIiMjVbFiRfXv39/muTwASoYn/QIwI5eWHcMwbrq+Vq1aWrdu3S23Exoaqm+//dZRsQAAgImUiRuUAQAAnIWyAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATI2yAwAATM2lZWfGjBm677775OPjo8DAQPXo0UNpaWk2Y65evaqYmBhVqVJFlSpVUq9evZSdnW0zJiMjQ926dZO3t7cCAwM1duxYXbt27U4eCgAAKKMcUnbOnz9fqtetW7dOMTEx2rx5s5KSkpSfn6/OnTvr8uXL1jGjRo3SN998o+XLl2vdunU6deqUevbsaV1fUFCgbt26KS8vT5s2bdLChQuVkJCgSZMm3e5hAQAAE7C77Lz11ltatmyZdb53796qUqWKatSooV27dtm1rcTERA0YMEBNmzZVixYtlJCQoIyMDO3YsUOSdOHCBX3yySeaPXu2/vznPys8PFzx8fHatGmTNm/eLEn6/vvvtX//fv3jH/9Qy5YtFR0drenTpysuLk55eXn2Hh4AADAZu8vOggULVKtWLUlSUlKSkpKS9N133yk6Olpjx469rTAXLlyQJAUEBEiSduzYofz8fHXq1Mk6plGjRqpdu7ZSUlIkSSkpKWrevLmCgoKsY6KiopSTk6N9+/YVu5/c3Fzl5OTYTAAAwJw87H1BVlaWteysWrVKvXv3VufOnVWnTh1FRESUOkhhYaFGjhypdu3aqVmzZtZ9lS9fXv7+/jZjg4KClJWVZR3zx6Jzff31dcWZMWOGpk6dWuqsAADg7mH3mZ3KlSvr+PHjkn6/DHX9rIthGCooKCh1kJiYGO3du1dLly4t9TZKavz48bpw4YJ1un48AADAfOw+s9OzZ089/fTTatCggc6cOaPo6GhJ0s6dO1W/fv1ShRg+fLhWrVql9evXq2bNmtblwcHBysvL0/nz523O7mRnZys4ONg6ZuvWrTbbu/5uretj/pOnp6c8PT1LlRUAANxd7D6z8+6772r48OFq0qSJkpKSVKlSJUlSZmamhg0bZte2DMPQ8OHDtXLlSq1Zs0ZhYWE268PDw1WuXDklJydbl6WlpSkjI0ORkZGSpMjISO3Zs0enT5+2jklKSpKvr6+aNGli7+EBAACTsfvMTrly5fTyyy8XWT5q1Ci7dx4TE6PFixfr66+/lo+Pj/UeGz8/P3l5ecnPz0+DBw/W6NGjFRAQIF9fX7300kuKjIzUAw88IEnq3LmzmjRpomeffVYzZ85UVlaWJkyYoJiYGM7eAACA0j1n5/PPP1f79u0VEhKiY8eOSZLmzJmjr7/+2q7tzJ8/XxcuXFCHDh1UvXp16/THt7a/++67evTRR9WrVy899NBDCg4O1ooVK6zr3d3dtWrVKrm7uysyMlLPPPOMnnvuOU2bNq00hwYAAEzG7jM78+fP16RJkzRy5Ei98cYb1puS/f39NWfOHHXv3r3E2zIM45ZjKlSooLi4OMXFxd1wTGhoqL799tsS7xcAAPz3sPvMzrx58/T3v/9dr7/+utzd3a3L27Rpoz179jg0HAAAwO2yu+ykp6erVatWRZZ7enrafMwDAABAWWB32QkLC1NqamqR5YmJiWrcuLEjMgEAADiM3ffsjB49WjExMbp69aoMw9DWrVu1ZMkSzZgxQx9//LEzMgIAAJSa3WXn+eefl5eXlyZMmKArV67o6aefVkhIiN577z317dvXGRkBAABKze6yI0n9+vVTv379dOXKFV26dEmBgYGOzgUAAOAQpSo713l7e8vb29tRWQAAAByuRGWnVatWslgsJdrgTz/9dFuBAAAAHKlEZadHjx5OjgEAAOAcJSo7kydPdnYOAAAApyj1PTvbt2/XgQMHJElNmjRReHi4w0IBAAA4it1l58SJE3rqqae0ceNG+fv7S5LOnz+vtm3baunSpapZs6ajMwIAAJSa3U9Qfv7555Wfn68DBw7o7NmzOnv2rA4cOKDCwkI9//zzzsgIAABQanaf2Vm3bp02bdqkhg0bWpc1bNhQ8+bN04MPPujQcAAAALfL7jM7tWrVUn5+fpHlBQUFCgkJcUgoAAAAR7G77Lz99tt66aWXtH37duuy7du3KzY2VrNmzXJoOAAAgNtVostYlStXtnmo4OXLlxURESEPj99ffu3aNXl4eGjQoEE8kwcAAJQpJSo7c+bMcXIMAAAA5yhR2enfv7+k38/gLF68WFFRUQoKCnJqMAAAAEew654dDw8P/eUvf9HVq1edlQcAAMCh7L5B+f7779fOnTudkQUAAMDh7H7OzrBhwzRmzBidOHFC4eHhqlixos36e++912HhAAAAbpfdZadv376SpBEjRliXWSwWGYYhi8WigoICx6UDAAC4TXaXnfT0dGfkAAAAcAq7y05oaKgzcgAAADiF3WXnuv379ysjI0N5eXk2yx9//PHbDgUAAOAodpedX375RU888YT27NljvVdHkvUJy9yzAwAAyhK733oeGxursLAwnT59Wt7e3tq3b5/Wr1+vNm3aaO3atU6ICAAAUHp2n9lJSUnRmjVrVLVqVbm5ucnNzU3t27fXjBkzNGLECJ7BAwAAyhS7z+wUFBTIx8dHklS1alWdOnVK0u83LqelpTk2HQAAwG2y+8xOs2bNtGvXLoWFhSkiIkIzZ85U+fLl9dFHH6lu3brOyAgAAFBqdpedCRMm6PLly5KkadOm6dFHH9WDDz6oKlWqaNmyZQ4PCAAAcDvsLjtRUVHWf9evX18///yzzp49q8qVK1vfkQUAAFBWlPo5O38UEBDgiM0AAAA4XInLzqBBg0o07tNPPy11GAAAAEcrcdlJSEhQaGioWrVqZX2QIAAAQFlX4rLz4osvasmSJUpPT9fAgQP1zDPPcPkKAACUeSV+zk5cXJwyMzP1yiuv6JtvvlGtWrXUu3dvrV69mjM9AACgzLLroYKenp566qmnlJSUpP3796tp06YaNmyY6tSpo0uXLjkrIwAAQKnZ/QRl6wvd3KwfBMqHfwIAgLLKrrKTm5urJUuW6JFHHtE999yjPXv26P3331dGRoYqVarkrIwAAAClVuIblIcNG6alS5eqVq1aGjRokJYsWaKqVas6MxsAAMBtK3HZWbBggWrXrq26detq3bp1WrduXbHjVqxY4bBwAAAAt6vEZee5557j4yAAAMBdx66HCgIAANxtSv1uLAAAgLsBZQcAAJgaZQcAAJgaZQcAAJhaicpO69atde7cOUnStGnTdOXKFaeGAgAAcJQSlZ0DBw7o8uXLkqSpU6fyOVgAAOCuUaK3nrds2VIDBw5U+/btZRiGZs2adcOPh5g0aVKJd75+/Xq9/fbb2rFjhzIzM7Vy5Ur16NHDun7AgAFauHChzWuioqKUmJhonT979qxeeuklffPNN3Jzc1OvXr303nvv8fEVAABAUgnLTkJCgiZPnqxVq1bJYrHou+++k4dH0ZdaLBa7ys7ly5fVokULDRo0SD179ix2TJcuXRQfH2+d9/T0tFnfr18/ZWZmKikpSfn5+Ro4cKCGDh2qxYsXlzgHAAAwrxKVnYYNG2rp0qWSfv+08+TkZAUGBt72zqOjoxUdHX3TMZ6engoODi523YEDB5SYmKht27apTZs2kqR58+apa9eumjVrlkJCQop9XW5urnJzc63zOTk5pTwCAABQ1tn9bqzCwkKHFJ2SWrt2rQIDA9WwYUO9+OKLOnPmjHVdSkqK/P39rUVHkjp16iQ3Nzdt2bLlhtucMWOG/Pz8rFOtWrWcegwAAMB1SvXW8yNHjuill15Sp06d1KlTJ40YMUJHjhxxdDZ16dJFn332mZKTk/XWW29p3bp1io6OVkFBgSQpKyurSPHy8PBQQECAsrKybrjd8ePH68KFC9bp+PHjDs8OAADKhhJ/NtZ1q1ev1uOPP66WLVuqXbt2kqSNGzeqadOm+uabb/TII484LFzfvn2t/27evLnuvfde1atXT2vXrlXHjh1LvV1PT88i9/4AAABzsrvsjBs3TqNGjdKbb75ZZPmrr77q0LLzn+rWrauqVavq8OHD6tixo4KDg3X69GmbMdeuXdPZs2dveJ8PAAD472L3ZawDBw5o8ODBRZYPGjRI+/fvd0ioGzlx4oTOnDmj6tWrS5IiIyN1/vx57dixwzpmzZo1KiwsVEREhFOzAACAu4PdZadatWpKTU0tsjw1NdXuG5cvXbqk1NRU6/bS09OVmpqqjIwMXbp0SWPHjtXmzZt19OhRJScnq3v37qpfv76ioqIkSY0bN1aXLl00ZMgQbd26VRs3btTw4cPVt2/fG74TCwAA/Hex+zLWkCFDNHToUP3yyy9q27atpN/v2Xnrrbc0evRou7a1fft2Pfzww9b566/v37+/5s+fr927d2vhwoU6f/68QkJC1LlzZ02fPt3mfptFixZp+PDh6tixo/WhgnPnzrX3sAAAgEnZXXYmTpwoHx8fvfPOOxo/frwkKSQkRFOmTNGIESPs2laHDh1kGMYN169evfqW2wgICOABggAA4IbsLjsWi0WjRo3SqFGjdPHiRUmSj4+Pw4MBAAA4gt1l548oOQAAoKwr1UMFAQAA7haUHQAAYGqUHQAAYGp2lZ38/Hx17NhRhw4dclYeAAAAh7LrBuVy5cpp9+7dzsoCAHCxOuP+5fBtHn2zm8O3CdjD7stYzzzzjD755BNnZAEAAHA4u996fu3aNX366af64YcfFB4erooVK9qsnz17tsPCAQAA3C67y87evXvVunVrSdLBgwdt1lksFsekAgAAcBC7y86PP/7ojBwAAABOUeonKB8+fFhHjhzRQw89JC8vLxmGwZkdoBjOuOFT4qZPACgpu29QPnPmjDp27Kh77rlHXbt2VWZmpiRp8ODBGjNmjMMDAgAA3A67y86oUaNUrlw5ZWRkyNvb27q8T58+SkxMdGg4AACA22X3Zazvv/9eq1evVs2aNW2WN2jQQMeOHXNYMAAAAEew+8zO5cuXbc7oXHf27Fl5eno6JBQAAICj2F12HnzwQX322WfWeYvFosLCQs2cOVMPP/ywQ8MBAADcLrsvY82cOVMdO3bU9u3blZeXp1deeUX79u3T2bNntXHjRmdkBAAAKDW7z+w0a9ZMBw8eVPv27dW9e3ddvnxZPXv21M6dO1WvXj1nZAQAACi1Uj1nx8/PT6+//rqjswAAADhcqcrOuXPn9Mknn+jAgQOSpCZNmmjgwIEKCAhwaDgAAIDbZfdlrPXr16tOnTqaO3euzp07p3Pnzmnu3LkKCwvT+vXrnZERAACg1Ow+sxMTE6M+ffpo/vz5cnd3lyQVFBRo2LBhiomJ0Z49exweEgAAoLTsPrNz+PBhjRkzxlp0JMnd3V2jR4/W4cOHHRoOAADgdtlddlq3bm29V+ePDhw4oBYtWjgkFAAAgKOU6DLW7t27rf8eMWKEYmNjdfjwYT3wwAOSpM2bNysuLk5vvvmmc1ICAACUUonKTsuWLWWxWGQYhnXZK6+8UmTc008/rT59+jguHQAAwG0qUdlJT093dg4AAACnKFHZCQ0NdXYOAAAApyjVQwVPnTqlDRs26PTp0yosLLRZN2LECIcEAwAAcAS7y05CQoJeeOEFlS9fXlWqVJHFYrGus1gslB0AAFCm2F12Jk6cqEmTJmn8+PFyc7P7nesAAAB3lN1t5cqVK+rbty9FBwAA3BXsbiyDBw/W8uXLnZEFAADA4ey+jDVjxgw9+uijSkxMVPPmzVWuXDmb9bNnz3ZYOAAAgNtVqrKzevVqNWzYUJKK3KAMAABQlthddt555x19+umnGjBggBPiAAAAOJbd9+x4enqqXbt2zsgCAADgcHaXndjYWM2bN88ZWQAAABzO7stYW7du1Zo1a7Rq1So1bdq0yA3KK1ascFg4AACA22V32fH391fPnj2dkQUAAMDh7C478fHxzsgBAADgFDwGGQAAmJrdZ3bCwsJu+jydX3755bYCAQAAOJLdZWfkyJE28/n5+dq5c6cSExM1duxYR+UCAABwCLvLTmxsbLHL4+LitH379tsOBAAA4EgOu2cnOjpaX375paM2BwAA4BAOKzv/+7//q4CAAEdtDgAAwCHsvozVqlUrmxuUDcNQVlaWfv31V33wwQcODQcAAHC77C47PXr0sJl3c3NTtWrV1KFDBzVq1MhRuQAAABzC7rIzefJkZ+QAAABwCpc+VHD9+vV67LHHFBISIovFoq+++spmvWEYmjRpkqpXry4vLy916tRJhw4dshlz9uxZ9evXT76+vvL399fgwYN16dKlO3gUAACgLCtx2XFzc5O7u/tNJw8P+04UXb58WS1atFBcXFyx62fOnKm5c+dqwYIF2rJliypWrKioqChdvXrVOqZfv37at2+fkpKStGrVKq1fv15Dhw61KwcAADCvEreTlStX3nBdSkqK5s6dq8LCQrt2Hh0drejo6GLXGYahOXPmaMKECerevbsk6bPPPlNQUJC++uor9e3bVwcOHFBiYqK2bdumNm3aSJLmzZunrl27atasWQoJCSl227m5ucrNzbXO5+Tk2JUbAADcPUpcdq4Xjj9KS0vTuHHj9M0336hfv36aNm2aw4Klp6crKytLnTp1si7z8/NTRESEUlJS1LdvX6WkpMjf399adCSpU6dOcnNz05YtW/TEE08Uu+0ZM2Zo6tSpDssKAADKrlLds3Pq1CkNGTJEzZs317Vr15SamqqFCxcqNDTUYcGysrIkSUFBQTbLg4KCrOuysrIUGBhos97Dw0MBAQHWMcUZP368Lly4YJ2OHz/usNwAAKBssesmmwsXLuhvf/ub5s2bp5YtWyo5OVkPPvigs7I5jaenpzw9PV0dAwAA3AElPrMzc+ZM1a1bV6tWrdKSJUu0adMmpxad4OBgSVJ2drbN8uzsbOu64OBgnT592mb9tWvXdPbsWesYAADw363EZ3bGjRsnLy8v1a9fXwsXLtTChQuLHbdixQqHBAsLC1NwcLCSk5PVsmVLSb/fSLxlyxa9+OKLkqTIyEidP39eO3bsUHh4uCRpzZo1KiwsVEREhENyAACAu1uJy85zzz1n8zERjnDp0iUdPnzYOp+enq7U1FQFBASodu3aGjlypP7617+qQYMGCgsL08SJExUSEmJ9inPjxo3VpUsXDRkyRAsWLFB+fr6GDx+uvn373vCdWAAA4L9LictOQkKCw3e+fft2Pfzww9b50aNHS5L69++vhIQEvfLKK7p8+bKGDh2q8+fPq3379kpMTFSFChWsr1m0aJGGDx+ujh07ys3NTb169dLcuXMdnhUAANyd7P64CEfq0KGDDMO44XqLxaJp06bd9C3tAQEBWrx4sTPiAQAAE3Dpx0UAAAA4G2UHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYGmUHAACYmoerAwClUWfcv5yy3aNvdnPKdgEArsOZHQAAYGqUHQAAYGqUHQAAYGqUHQAAYGqUHQAAYGqUHQAAYGqUHQAAYGpluuxMmTJFFovFZmrUqJF1/dWrVxUTE6MqVaqoUqVK6tWrl7Kzs12YGAAAlDVl/qGCTZs21Q8//GCd9/D4v8ijRo3Sv/71Ly1fvlx+fn4aPny4evbsqY0bN7oiKgDgDnLGw0V5sKg5lfmy4+HhoeDg4CLLL1y4oE8++USLFy/Wn//8Z0lSfHy8GjdurM2bN+uBBx6401EBAEAZVKYvY0nSoUOHFBISorp166pfv37KyMiQJO3YsUP5+fnq1KmTdWyjRo1Uu3ZtpaSk3HSbubm5ysnJsZkAAIA5lemyExERoYSEBCUmJmr+/PlKT0/Xgw8+qIsXLyorK0vly5eXv7+/zWuCgoKUlZV10+3OmDFDfn5+1qlWrVpOPAoAAOBKZfoyVnR0tPXf9957ryIiIhQaGqovvvhCXl5epd7u+PHjNXr0aOt8Tk4OhQcAAJMq02d2/pO/v7/uueceHT58WMHBwcrLy9P58+dtxmRnZxd7j88feXp6ytfX12YCAADmdFeVnUuXLunIkSOqXr26wsPDVa5cOSUnJ1vXp6WlKSMjQ5GRkS5MCQAAypIyfRnr5Zdf1mOPPabQ0FCdOnVKkydPlru7u5566in5+flp8ODBGj16tAICAuTr66uXXnpJkZGRvBMLAABYlemyc+LECT311FM6c+aMqlWrpvbt22vz5s2qVq2aJOndd9+Vm5ubevXqpdzcXEVFRemDDz5wcWoAAFCWlOmys3Tp0puur1ChguLi4hQXF3eHEgEAgLvNXXXPDgAAgL0oOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQoOwAAwNQ8XB0AAACzqzPuXw7f5tE3uzl8m2bFmR0AAGBqlB0AAGBqlB0AAGBq3LMDK2dcU5a4rgwAcC3O7AAAAFOj7AAAAFOj7AAAAFOj7AAAAFOj7AAAAFOj7AAAAFOj7AAAAFOj7AAAAFOj7AAAAFMzTdmJi4tTnTp1VKFCBUVERGjr1q2ujgQAAMoAU3xcxLJlyzR69GgtWLBAERERmjNnjqKiopSWlqbAwECXZuMjGAAAdxNn/N1y9d8sU5zZmT17toYMGaKBAweqSZMmWrBggby9vfXpp5+6OhoAAHCxu/7MTl5ennbs2KHx48dbl7m5ualTp05KSUkp9jW5ubnKzc21zl+4cEGSlJOT4/B8hblXHL5NiaxkvXuyOiOnRFaykpWs/7ddwzBuPtC4y508edKQZGzatMlm+dixY43777+/2NdMnjzZkMTExMTExMRkgun48eM37Qp3/Zmd0hg/frxGjx5tnS8sLNTZs2dVpUoVWSwWl2TKyclRrVq1dPz4cfn6+rokQ0mR1TnI6hxkdY67JevdklMia2kYhqGLFy8qJCTkpuPu+rJTtWpVubu7Kzs722Z5dna2goODi32Np6enPD09bZb5+/s7K6JdfH19y/wP+XVkdQ6yOgdZneNuyXq35JTIai8/P79bjrnrb1AuX768wsPDlZycbF1WWFio5ORkRUZGujAZAAAoC+76MzuSNHr0aPXv319t2rTR/fffrzlz5ujy5csaOHCgq6MBAAAXM0XZ6dOnj3799VdNmjRJWVlZatmypRITExUUFOTqaCXm6empyZMnF7m8VhaR1TnI6hxkdY67JevdklMiqzNZDONW79cCAAC4e9319+wAAADcDGUHAACYGmUHAACYGmUHAACYGmWnDHjjjTfUtm1beXt73/DhhhkZGerWrZu8vb0VGBiosWPH6tq1a3c2aDEOHjyo7t27q2rVqvL19VX79u31448/ujrWDf3rX/9SRESEvLy8VLlyZfXo0cPVkW4qNzdXLVu2lMViUWpqqqvjFHH06FENHjxYYWFh8vLyUr169TR58mTl5eW5OpokKS4uTnXq1FGFChUUERGhrVu3ujpSETNmzNB9990nHx8fBQYGqkePHkpLS3N1rBJ58803ZbFYNHLkSFdHKdbJkyf1zDPPqEqVKvLy8lLz5s21fft2V8cqoqCgQBMnTrT5PZo+ffqtP+/pDli/fr0ee+wxhYSEyGKx6KuvvrJZbxiGJk2apOrVq8vLy0udOnXSoUOHXBP2Jig7ZUBeXp6efPJJvfjii8WuLygoULdu3ZSXl6dNmzZp4cKFSkhI0KRJk+5w0qIeffRRXbt2TWvWrNGOHTvUokULPfroo8rKynJ1tCK+/PJLPfvssxo4cKB27dqljRs36umnn3Z1rJt65ZVXbvkYdFf6+eefVVhYqA8//FD79u3Tu+++qwULFui1115zdTQtW7ZMo0eP1uTJk/XTTz+pRYsWioqK0unTp10dzca6desUExOjzZs3KykpSfn5+ercubMuX77s6mg3tW3bNn344Ye69957XR2lWOfOnVO7du1Urlw5fffdd9q/f7/eeecdVa5c2dXRinjrrbc0f/58vf/++zpw4IDeeustzZw5U/PmzXN1NF2+fFktWrRQXFxcsetnzpypuXPnasGCBdqyZYsqVqyoqKgoXb169Q4nvQVHfBgnHCM+Pt7w8/Mrsvzbb7813NzcjKysLOuy+fPnG76+vkZubu4dTGjr119/NSQZ69evty7LyckxJBlJSUkuy1Wc/Px8o0aNGsbHH3/s6igl9u233xqNGjUy9u3bZ0gydu7c6epIJTJz5kwjLCzM1TGM+++/34iJibHOFxQUGCEhIcaMGTNcmOrWTp8+bUgy1q1b5+ooN3Tx4kWjQYMGRlJSkvGnP/3JiI2NdXWkIl599VWjffv2ro5RIt26dTMGDRpks6xnz55Gv379XJSoeJKMlStXWucLCwuN4OBg4+2337YuO3/+vOHp6WksWbLEBQlvjDM7d4GUlBQ1b97c5iGJUVFRysnJ0b59+1yWq0qVKmrYsKE+++wzXb58WdeuXdOHH36owMBAhYeHuyxXcX766SedPHlSbm5uatWqlapXr67o6Gjt3bvX1dGKlZ2drSFDhujzzz+Xt7e3q+PY5cKFCwoICHBphry8PO3YsUOdOnWyLnNzc1OnTp2UkpLiwmS3duHCBUly+dfwZmJiYtStWzebr29Z889//lNt2rTRk08+qcDAQLVq1Up///vfXR2rWG3btlVycrIOHjwoSdq1a5c2bNig6OhoFye7ufT0dGVlZdn8HPj5+SkiIqLM/Z6Z4gnKZpeVlVXkadDX5115uchiseiHH35Qjx495OPjIzc3NwUGBioxMbHMnSr+5ZdfJElTpkzR7NmzVadOHb3zzjvq0KGDDh48WKb+sBiGoQEDBugvf/mL2rRpo6NHj7o6UokdPnxY8+bN06xZs1ya47ffflNBQUGxvzc///yzi1LdWmFhoUaOHKl27dqpWbNmro5TrKVLl+qnn37Stm3bXB3lpn755RfNnz9fo0eP1muvvaZt27ZpxIgRKl++vPr37+/qeDbGjRunnJwcNWrUSO7u7iooKNAbb7yhfv36uTraTV3/+1Pc71lZu5WBMztOMm7cOFkslptOZfU/uiXNbhiGYmJiFBgYqH//+9/aunWrevTooccee0yZmZllKmthYaEk6fXXX1evXr0UHh6u+Ph4WSwWLV++vExlnTdvni5evKjx48ffkVy3k/WPTp48qS5duujJJ5/UkCFDXJT87hYTE6O9e/dq6dKlro5SrOPHjys2NlaLFi1ShQoVXB3npgoLC9W6dWv97W9/U6tWrTR06FANGTJECxYscHW0Ir744gstWrRIixcv1k8//aSFCxdq1qxZWrhwoaujmQZndpxkzJgxGjBgwE3H1K1bt0TbCg4OLvIukuzsbOs6Rytp9jVr1mjVqlU6d+6cfH19JUkffPCBkpKStHDhQo0bN87h2Uqb9Xr5atKkiXW5p6en6tatq4yMDGdGtLLn65qSklLkM2fatGmjfv363ZH/ANr783vq1Ck9/PDDatu2rT766CMnp7u1qlWryt3d3fp7cl12drZTfmccYfjw4Vq1apXWr1+vmjVrujpOsXbs2KHTp0+rdevW1mUFBQVav3693n//feXm5srd3d2FCf9P9erVbX7fJalx48b68ssvXZToxsaOHatx48apb9++kqTmzZvr2LFjmjFjRpk7C/VH13+XsrOzVb16devy7OxstWzZ0kWpikfZcZJq1aqpWrVqDtlWZGSk3njjDZ0+fVqBgYGSpKSkJPn6+hb5ZXaEkma/cuWKpN/vhfgjNzc365kUZytp1vDwcHl6eiotLU3t27eXJOXn5+vo0aMKDQ11dkxJJc86d+5c/fWvf7XOnzp1SlFRUVq2bJkiIiKcGdHKnp/fkydP6uGHH7aeLfvPnwdXKF++vMLDw5WcnGx9vEBhYaGSk5M1fPhw14b7D4Zh6KWXXtLKlSu1du1ahYWFuTrSDXXs2FF79uyxWTZw4EA1atRIr776apkpOpLUrl27Im/hP3jw4B37fbfHlStXivzeuLu737H/jpZWWFiYgoODlZycbC03OTk52rJlyw3fXewyLr5BGoZhHDt2zNi5c6cxdepUo1KlSsbOnTuNnTt3GhcvXjQMwzCuXbtmNGvWzOjcubORmppqJCYmGtWqVTPGjx/v0ty//vqrUaVKFaNnz55GamqqkZaWZrz88stGuXLljNTUVJdmK05sbKxRo0YNY/Xq1cbPP/9sDB482AgMDDTOnj3r6mg3lZ6eXmbfjXXixAmjfv36RseOHY0TJ04YmZmZ1snVli5danh6ehoJCQnG/v37jaFDhxr+/v4272osC1588UXDz8/PWLt2rc3X78qVK66OViJl9d1YW7duNTw8PIw33njDOHTokLFo0SLD29vb+Mc//uHqaEX079/fqFGjhrFq1SojPT3dWLFihVG1alXjlVdecXU04+LFi9a/SZKM2bNnGzt37jSOHTtmGIZhvPnmm4a/v7/x9ddfG7t37za6d+9uhIWFGf/v//0/Fye3RdkpA/r3729IKjL9+OOP1jFHjx41oqOjDS8vL6Nq1arGmDFjjPz8fNeF/v9t27bN6Ny5sxEQEGD4+PgYDzzwgPHtt9+6Olax8vLyjDFjxhiBgYGGj4+P0alTJ2Pv3r2ujnVLZbnsxMfHF/uzW1b+P2revHlG7dq1jfLlyxv333+/sXnzZldHKuJGX7/4+HhXRyuRslp2DMMwvvnmG6NZs2aGp6en0ahRI+Ojjz5ydaRi5eTkGLGxsUbt2rWNChUqGHXr1jVef/11lz5a5Loff/yx2J/P/v37G4bx+9vPJ06caAQFBRmenp5Gx44djbS0NNeGLobFMMrAIxoBAACcxPUX1wEAAJyIsgMAAEyNsgMAAEyNsgMAAEyNsgMAAEyNsgMAAEyNsgMAAEyNsgMAAEyNsgOY2NGjR2WxWJSamurqKLc0ZcqUMvPhgXXq1NGcOXNuOsZiseirr766I3lKo0OHDho5cuRNxyQkJMjf3/+O5AFcibIDuNCAAQOsH1T5R2vXrpXFYtH58+dva/u1atVSZmammjVrdlvbsceNjul2XS9u7u7uOnnypM26zMxMeXh4yGKx6OjRow7ftzNcP57rU5UqVdS5c2ft3LnTIdtfsWKFpk+fbp0vrsD16dNHBw8edMj+gLKMsgOYVF5entzd3RUcHCwPDw9Xx3GYGjVq6LPPPrNZtnDhQtWoUeO2t52Xl3fb27DXDz/8oMzMTK1evVqXLl1SdHT0bZdcSQoICJCPj89Nx3h5eSkwMPC29wWUdZQd4C7x5ZdfqmnTpvL09FSdOnX0zjvv2KyvU6eOpk+frueee06+vr4aOnRokctYAwYMsDmbcH1au3atJOncuXN67rnnVLlyZXl7eys6OlqHDh2y7uP6ZY/Vq1ercePGqlSpkrp06aLMzExJv1+KWrhwob7++usi23711Vd1zz33yNvbW3Xr1tXEiROVn59v99ehf//+io+Pt1kWHx+v/v372ywrKCjQ4MGDFRYWJi8vLzVs2FDvvfeezZjrZ6HeeOMNhYSEqGHDhtZ1Fy9e1FNPPaWKFSuqRo0aiouLu2mu48ePq3fv3vL391dAQIC6d+9eorNMVapUUXBwsNq0aaNZs2YpOztbW7ZskXTr7/kHH3ygBg0aqEKFCgoKCtL//M//WNf98TJWhw4ddOzYMY0aNcr6fZGKv4w1f/581atXT+XLl1fDhg31+eef26y3WCz6+OOP9cQTT8jb21sNGjTQP//5z1seJ+BKlB3gLrBjxw717t1bffv21Z49ezRlyhRNnDhRCQkJNuNmzZqlFi1aaOfOnZo4cWKR7bz33nvKzMy0TrGxsQoMDFSjRo0k/f7Hf/v27frnP/+plJQUGYahrl272pSSK1euaNasWfr888+1fv16ZWRk6OWXX5Ykvfzyy+rdu7e1AGVmZqpt27aSJB8fHyUkJGj//v1677339Pe//13vvvuu3V+Lxx9/XOfOndOGDRskSRs2bNC5c+f02GOP2YwrLCxUzZo1tXz5cu3fv1+TJk3Sa6+9pi+++MJmXHJystLS0pSUlKRVq1ZZl7/99tvWr+W4ceMUGxurpKSkYjPl5+crKipKPj4++ve//62NGzdai6A9Z4u8vLwk/X6G6Vbf8+3bt2vEiBGaNm2a0tLSlJiYqIceeqjY7a5YsUI1a9bUtGnTrN+X4qxcuVKxsbEaM2aM9u7dqxdeeEEDBw7Ujz/+aDNu6tSp6t27t3bv3q2uXbuqX79+Onv2bImPE7jjXPyp68B/tf79+xvu7u5GxYoVbaYKFSoYkoxz584ZhmEYTz/9tPHII4/YvHbs2LFGkyZNrPOhoaFGjx49bMakp6cbkoydO3cW2feXX35pVKhQwdiwYYNhGIZx8OBBQ5KxceNG65jffvvN8PLyMr744gvDMAwjPj7ekGQcPnzYOiYuLs4ICgqyOabu3bvf8tjffvttIzw83Do/efJko0WLFjcc/8djGTlypDFw4EDDMAxj4MCBxqhRo4ydO3cakoz09PQbbiMmJsbo1auXTdagoCAjNzfXZlxoaKjRpUsXm2V9+vQxoqOjrfOSjJUrVxqGYRiff/650bBhQ6OwsNC6Pjc31/Dy8jJWr159y+MxDMM4d+6c8cQTTxiVKlUysrKybvk9//LLLw1fX18jJyen2O3/6U9/MmJjY22O6d1337UZEx8fb/j5+Vnn27ZtawwZMsRmzJNPPml07drV5rgnTJhgnb906ZIhyfjuu++KzQGUBZzZAVzs4YcfVmpqqs308ccf24w5cOCA2rVrZ7OsXbt2OnTokAoKCqzL2rRpU6J97ty5U88++6zef/9963YPHDggDw8PRUREWMdVqVJFDRs21IEDB6zLvL29Va9ePet89erVdfr06Vvuc9myZWrXrp2Cg4NVqVIlTZgwQRkZGSXK+58GDRqk5cuXKysrS8uXL9egQYOKHRcXF6fw8HBVq1ZNlSpV0kcffVRkn82bN1f58uWLvDYyMrLI/B+/Dn+0a9cuHT58WD4+PqpUqZIqVaqkgIAAXb16VUeOHLnpsbRt21aVKlVS5cqVtWvXLi1btkxBQUG3/J4/8sgjCg0NVd26dfXss89q0aJFunLlyk33dSs32ud/Hve9995r/XfFihXl6+tbop8BwFXMc9cicJeqWLGi6tevb7PsxIkTpd7WrWRlZenxxx/X888/r8GDB9u9j3LlytnMWywWGYZx09ekpKSoX79+mjp1qqKiouTn56elS5cWuQelpJo3b65GjRrpqaeeUuPGjdWsWbMib69funSpXn75Zb3zzjuKjIyUj4+P3n77bev9MNeV5Gt2K5cuXVJ4eLgWLVpUZF21atVu+tply5apSZMmqlKlil1vA/fx8dFPP/2ktWvX6vvvv9ekSZM0ZcoUbdu2zelvJy/uZ6CwsNCp+wRuB2d2gLtA48aNtXHjRptlGzdu1D333CN3d/cSb+fq1avq3r27GjVqpNmzZxfZx7Vr12zKwJkzZ5SWlqYmTZqUeB/ly5e3OdskSZs2bVJoaKhef/11tWnTRg0aNNCxY8dKvM3iDBo0SGvXrr3hWZ2NGzeqbdu2GjZsmFq1aqX69evf8izLH23evLnIfOPGjYsd27p1ax06dEiBgYGqX7++zeTn53fT/dSqVUv16tUrUlBK8j338PBQp06dNHPmTO3evVtHjx7VmjVrit1Pcd+X/3Sjfdrz/QfKIsoOcBcYM2aMkpOTNX36dB08eFALFy7U+++/b70xuKReeOEFHT9+XHPnztWvv/6qrKwsZWVlKS8vTw0aNFD37t01ZMgQbdiwQbt27dIzzzyjGjVqqHv37iXeR506dbR7926lpaXpt99+U35+vho0aKCMjAwtXbpUR44c0dy5c7Vy5Up7vww2hgwZol9//VXPP/98sesbNGig7du3a/Xq1Tp48KAmTpyobdu2lXj7Gzdu1MyZM3Xw4EHFxcVp+fLlio2NLXZsv379VLVqVXXv3l3//ve/lZ6errVr12rEiBGlPkt3q+/5qlWrNHfuXKWmpurYsWP67LPPVFhYaPOOsj+qU6eO1q9fr5MnT+q3334rdszYsWOVkJCg+fPn69ChQ5o9e7ZWrFhh988ZUNZQdoC7QOvWrfXFF19o6dKlatasmSZNmqRp06ZpwIABdm1n3bp1yszMVJMmTVS9enXrtGnTJkm/v4U7PDxcjz76qCIjI2UYhr799tsily1uZsiQIWrYsKHatGmjatWqaePGjXr88cc1atQoDR8+XC1bttSmTZuKfbeYPTw8PFS1atUbPkPohRdeUM+ePdWnTx9FRETozJkzGjZsWIm3P2bMGG3fvl2tWrXSX//6V82ePVtRUVHFjvX29tb69etVu3Zt9ezZU40bN9bgwYN19epV+fr6lur4bvU99/f314oVK/TnP/9ZjRs31oIFC7RkyRI1bdq02O1NmzZNR48eVb169W54aa1Hjx567733NGvWLDVt2lQffvih4uPj1aFDh1IdA1BWWIxbXWwHAAC4i3FmBwAAmBplBwAAmBplBwAAmBplBwAAmBplBwAAmBplBwAAmBplBwAAmBplBwAAmBplBwAAmBplBwAAmBplBwAAmNr/B3vHo8SMeMvyAAAAAElFTkSuQmCC\n", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "bins = np.arange(-10,12)\n", + "plt.hist(sim, bins=bins, align='left')\n", + "plt.xticks(bins[::2])\n", + "plt.xlabel('Horizontal Marble Position')\n", + "plt.ylabel('Number of Marbles');" + ] + }, + { + "cell_type": "markdown", + "id": "3b1ed737-e851-4272-b706-67fbb91f5a63", + "metadata": {}, + "source": [ + "One interesting trend that appears in the above histogram is that all locations are integers. This is because we used whole numbers for our left and right movements of marbles (i.e., +1 and -1) and used ten layers in the above simulation. There is no way to produce an odd integer by totaling a series of ten +1 or -1 integers. If we changed the number of layers to an odd number, we'd only get odd positions in the result, and if we used +1/2 and -1/2 for our horizontal movement, we'd get both even and odd integers." + ] + }, + { + "cell_type": "markdown", + "id": "7152df76-78b5-47fa-9f76-276c7e497cb7", + "metadata": {}, + "source": [ + "## Calculating pi\n", + "\n", + "As a second appliation of random number generaters, we will calclate a numerical value for $\\pi$. Imagine a circle with radius = $r$ inscribed inside a square of side length $2r$. " + ] + }, + { + "cell_type": "markdown", + "id": "65166f10-9009-4551-a14f-8205dc3b858c", + "metadata": {}, + "source": [ + "![](tutorial-stochastic/circle_square.png)" + ] + }, + { + "cell_type": "markdown", + "id": "88d705e8-b983-4aa2-b35b-ccb13c9aaa6a", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "**Figure 2.** The circle radius r inscribed in a square with side length of 2r.\n", + "\n", + "The areas of the circle and square are described by $A\n", + "_{circle} = \\pi r^2$ and $A_{square} = 4r^2$, respectively. If we divide the area of the circle by the area of the square, we get $\\pi / 4$.\n", + "\n", + "$$ \\frac { A_{circle} }{ A_{square} } = \\frac {\\pi r^2} {(2r)^2} = \\frac {\\pi r^2} {4r^2} = \\frac {\\pi}{4} $$\n", + "\n", + "This means that four times the ratio of the area of the circle and square equals $\\pi$.\n", + "\n", + "$$ \\frac { A_{circle} }{ A_{square}} \\times 4 = \\pi $$\n", + "\n", + "The challenge is now finding the ratio between the areas of the two shapes which can be done through a variety of creative means. Our method will be counting randomly placed dots. Imagine we painted the circle inscribed in a square on pavement outside and count the water droplets that fall in the circle and square during a rain storm." + ] + }, + { + "cell_type": "markdown", + "id": "7571b83c-e163-42b0-8329-f494045790a8", + "metadata": {}, + "source": [ + "![](tutorial-stochastic/circle_square_dots.png)" + ] + }, + { + "cell_type": "markdown", + "id": "5c12bceb-c2ad-4cc2-98b9-4148667dc1c7", + "metadata": {}, + "source": [ + "**Figure 3.** Random dots inside the circle and square used to estimate a value of $\\pi$.\n", + "\n", + "Because the rain drops fall randomly, the probability of a drop falling in one of the shapes is proportional to the area of that shape, so dividing the number of drops in the circle by the number of drop in the square provides a estimate of $\\frac { A_{circle} }{ A_{square} }$ which then can be used to calculate $\\pi$. Counting rain is a bit tedious and challenging, so we can simulate this using a random number generator to produce random ($x, y$) coordinates of dots inside the square. We then total the number of these that fall inside the circle versus the square and complete the calculation to obtain a value of $\\pi$.\n", + "\n", + "For smilicity, we will center our circle and square around the origin and give the circle a radius = 1. This requires the coordinates to fall in the [-1, 1) ranges along both the $x$- and $y$-axes. We have no random number generator that produces values in this range, but we can modify the `random()` method which generates floats in the [0, 1) range by multiplying by 2 and substracting 1 as demonstrated below.\n", + "\n", + "~~~python\n", + "rng = np.random.default_rng()\n", + "x = 2 * rng.random() - 1\n", + "~~~" + ] + }, + { + "cell_type": "markdown", + "id": "70c285c1-adae-453e-a895-fb515ea562ab", + "metadata": {}, + "source": [ + "Below is a function that performs this simulation and calculation of $\\pi$ with the only function argument `n_samples` as the number of dots generated. The general proceedure is ouline here.\n", + "1. Generate `n_samples` of ($x, y$) coordinates using the `random()` method\n", + "2. Calculate the distance from the origin using the `np.hypot()` function to determine if the dot is inside the circle (all dots are inside the square)\n", + "3. Calculate $\\pi$" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "id": "e70bfca5-3606-4c3c-ad76-fc715f1f5ef4", + "metadata": {}, + "outputs": [], + "source": [ + "def estimate_pi(n_samples):\n", + " '''(integer) -> pi (float)\n", + " \n", + " Estimates the value of pi by finding the ratio of the area of a unit circle/area of a unit square\n", + " and multiplying by 4.\n", + " '''\n", + "\n", + " n_samples = int(n_samples)\n", + " \n", + " # Step 1\n", + " rng = np.random.default_rng()\n", + " coords = rng.random(size=(n_samples, 2))\n", + "\n", + " # Step 2\n", + " dist_from_origin = np.hypot(coords[:,0], coords[:,1]) # distance from the origion\n", + " in_circle = dist_from_origin[dist_from_origin <= 1] # if <= 1, the point is inside the circle\n", + "\n", + " # Step 3\n", + " pi = 4 * (in_circle.size / n_samples)\n", + " \n", + " return pi" + ] + }, + { + "cell_type": "markdown", + "id": "d34736f4-21f4-4e15-a3b3-d68cbf90b782", + "metadata": {}, + "source": [ + "We can then call this function with 100 dots generated. Because we are using a random number generator, every time we call this function, we likley get a different value for $\\pi$." + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "id": "c3bac1dc-e420-4b2d-a0b7-8b3fce9f277b", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "3.08" + ] + }, + "execution_count": 12, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "estimate_pi(100)" + ] + }, + { + "cell_type": "code", + "execution_count": 13, + "id": "19575a49-311d-4c27-aa34-2cd8a3539265", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "3.16" + ] + }, + "execution_count": 13, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "estimate_pi(100)" + ] + }, + { + "cell_type": "markdown", + "id": "947f9bfd-9550-440b-820b-eeb746dad01d", + "metadata": {}, + "source": [ + "We can also generate a value for $\\pi$ using the `np.pi` to compare our answer to." + ] + }, + { + "cell_type": "code", + "execution_count": 14, + "id": "9d107cde-6f0a-48ae-9258-0392f41cc37d", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "3.141592653589793" + ] + }, + "execution_count": 14, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "np.pi" + ] + }, + { + "cell_type": "markdown", + "id": "dfdb0029-6ef8-49d9-997a-2b25ad65fe5e", + "metadata": {}, + "source": [ + "Odds are that the values calculated above for $\\pi$ above is close but not exactly 3.14. If we increase the number of dots, the value for $\\pi$ will on average become closer to the true value. While it may be tempting to just call the function with an extremely large number of dots, this can take substatially longer to computer." + ] + }, + { + "cell_type": "code", + "execution_count": 15, + "id": "d49d0955-f5e9-4214-abb2-2f0c9c4f470c", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "3.14157644" + ] + }, + "execution_count": 15, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "estimate_pi(1E8)" + ] + }, + { + "cell_type": "markdown", + "id": "cbc28bcb-b975-4d82-9f06-6a6378bf09c2", + "metadata": {}, + "source": [ + "## Simulating radioactive decay\n", + "\n", + "As a final example, we will simulate the radioactive decay of the radioactive nuclide $^{137m}$Ba to determine the quantity of undecayed material versus time. We will then compare it to an analytical model to see how well the two agree. Radioactive nuclei decay by first-order kinetics, which means that each radioactive nucleus has a fixed probability of decaying at any given time. This is analogous to how each six-sided die in a bucket has a fixed probability of rolling a one each time they are dumped out. This analogy is so good that we can simulate the decay process by generating a large number of random integers with fixed probabilities.\n", + "\n", + "The first thing we need to is to know the probability of a nucleus decaying in a single second. Because this is a first-order process, the rate (-dA/dt) is described by\n", + "\n", + "$$ \\frac{-dA}{dt} = kA $$\n", + "\n", + "where A is the quantity of $^{137m}$Ba, $t$ is time in seconds, and $k$ is the rate constant for this process. We also know that the half-life (t$_{1/2}$) of $^{137m}$Ba is about 151 seconds which is related to the rate constant by the following equation.\n", + "\n", + "$$ t_{1/2} = \\frac{ln(2)}{k} $$\n", + "\n", + "This means that $k$ = 4.59e$^{-3}$ s$^{-1}$. For first-order kinetics, this means that 4.59e$^{-3}$ fraction of the nuclei decays, so this is the probability we need to simulate this process using random numbers. The function below calculates the fraction of radioactive material remaining based on a rate constant `k` and allows the users to choose the number of simulated nuclei `n` and duraction of the simulation `n_final`." + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "id": "f5ec62ec-f980-43ea-b5ea-19d7a6de8450", + "metadata": {}, + "outputs": [], + "source": [ + "def decay(k, n=1000, t_final=1000):\n", + " '''floats and integers -> ndarray of integers\n", + "\n", + " Calculates the fraction of nuclei remaining at each second in the simulation given n initial \n", + " nuclei with a rate constant k of decaying each second and a final simulation time t_final.\n", + " '''\n", + "\n", + " n_undecayed = n\n", + " undecayed_array = np.full(t_final + 1, n)\n", + " \n", + " for second in range(1, t_final + 1):\n", + " decays = rng.binomial(1, p=k, size=n_undecayed).sum()\n", + " n_undecayed -= decays\n", + " undecayed_array[second] = n_undecayed\n", + " \n", + " fraction = undecayed_array / n \n", + " \n", + " \n", + " return fraction" + ] + }, + { + "cell_type": "markdown", + "id": "61e2bc98-4518-4b2b-a33d-b8303ee65a53", + "metadata": {}, + "source": [ + "Below we simulate the decay of $^{137m}$Ba using a hundred simulated nuclei." + ] + }, + { + "cell_type": "code", + "execution_count": 17, + "id": "bd9d6811-2248-4c12-b834-385658b40baa", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "array([1. , 1. , 1. , ..., 0.01, 0.01, 0.01])" + ] + }, + "execution_count": 17, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "decay(4.59e-3, n=100)" + ] + }, + { + "cell_type": "markdown", + "id": "8320f7c8-31ea-4272-8683-9cd477f943c6", + "metadata": {}, + "source": [ + "We can plot and compare this simulation to the analytical model for this simulation knowing that the fraction of undecayed nuclei ($\\frac{A}{A_0}$) is described by the following first-order equation where $A_0$ and $A_t$ are the initial quantities and quantities of nuclei at time $t$.\n", + "\n", + "$$ \\frac{A_t}{A_0} = e^{-kt} $$" + ] + }, + { + "cell_type": "code", + "execution_count": 18, + "id": "aff9c2ec-a032-4fc7-8ca0-46502ab1b31c", + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAGwCAYAAABVdURTAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABqx0lEQVR4nO3dd3gUVd/G8e9uekiFhIRAIIQiRHondI0iKIIVEQGxPRQVRSw8ishrAQuIHcWC2AAVERVQCNKU3pEOAQKkQUhCElJ33j9WFvNQzMImm3J/rmvN7pmy9w6Y/TFz5hyTYRgGIiIiIhWE2dkBRERERBxJxY2IiIhUKCpuREREpEJRcSMiIiIVioobERERqVBU3IiIiEiFouJGREREKhRXZwcobRaLhePHj+Pr64vJZHJ2HBERESkGwzA4ffo0YWFhmM2XPjdT6Yqb48ePEx4e7uwYIiIichni4+OpVavWJdepdMWNr68vYD04fn5+Tk4jIiIixZGRkUF4eLjte/xSKl1xc/ZSlJ+fn4obERGRcqY4XUrUoVhEREQqFBU3IiIiUqGouBEREZEKpdL1uRERkfLBYrGQl5fn7BhSitzd3f/1Nu/iUHEjIiJlTl5eHnFxcVgsFmdHkVJkNpupW7cu7u7uV7QfFTciIlKmGIZBQkICLi4uhIeHO+Rf8lL2nR1kNyEhgdq1a1/RQLsqbkREpEwpKCggOzubsLAwvL29nR1HSlFwcDDHjx+noKAANze3y96PymERESlTCgsLAa740oSUP2f/zM/+HbhcKm5ERKRM0vx/lY+j/sxV3IiIiEiF4tTiZsWKFfTp04ewsDBMJhPz5s37122WLVtGq1at8PDwoH79+syYMaPEc4qIiEj54dTiJisri+bNm/Pee+8Va/24uDhuvPFGevTowZYtW3jsscd44IEH+PXXX0s4qYiIyJVZtmwZJpOJtLQ0Z0exS3FPPhRXREQEU6dOddj+LsSpd0v16tWLXr16FXv9adOmUbduXSZPngxA48aNWbVqFW+++SY9e/YsqZjFdjIzl+TTuTSuoQk5RUQqk3/rKzJ+/Hi6d+9eOmEu0wsvvMC8efPYsmVLkfaEhAQCAwOdE+oylas+N6tXryYmJqZIW8+ePVm9evVFt8nNzSUjI6PIoyT8+lcinV5ewIez55XI/kVEpOxKSEiwPaZOnYqfn1+RtjFjxjgt25WO8hwaGoqHh4eD0pSOclXcJCYmEhISUqQtJCSEjIwMzpw5c8FtJk6ciL+/v+0RHh5eItlaVc1jm/v9vH7qMVJST5XIe4iIVEaGYZCdV+CUh2EYxcoYGhpqe/j7+2MymYq0+fj42NbduHEjbdq0wdvbm+joaPbs2VNkXz/++COtWrXC09OTyMhIJkyYQEFBgW35kSNH6Nu3Lz4+Pvj5+XHnnXeSlJRkW/7CCy/QokULPv74Y+rWrYunpycAaWlpPPDAAwQHB+Pn58c111zD1q1bAZgxYwYTJkxg69atmEwmTCaTrU/r/16WOnr0KAMGDKBq1apUqVKFNm3asHbtWgAOHDhA3759CQkJwcfHh7Zt27JkyZLi/2E7SIUfxG/s2LGMHj3a9jojI6NECpzg0HBOmv2pZqSyff1yrunZz+HvISJSGZ3JLyTqeef0rdz5fz3xdnfsV+Wzzz7L5MmTCQ4OZtiwYdx333388ccfAKxcuZLBgwfz9ttv06VLFw4cOMBDDz0EWC9tWSwWW2GzfPlyCgoKGDlyJP3792fZsmW299i/fz/ff/89c+fOxcXFBYA77rgDLy8vFi5ciL+/Px9++CHXXnste/fupX///uzYsYNFixbZihF/f//zsmdmZtKtWzdq1qzJ/PnzCQ0NZdOmTbZpMjIzM+nduzcvv/wyHh4ezJw5kz59+rBnzx5q167t0ON4KeWquAkNDS1SnQIkJSXh5+eHl5fXBbfx8PAondNpJhOpgc2olrqMtD1/gIobERG5gJdffplu3boB8Mwzz3DjjTeSk5ODp6cnEyZM4JlnnmHIkCEAREZG8uKLL/LUU08xfvx4YmNj2b59O3FxcbZ/qM+cOZOrr76a9evX07ZtW8B6KWrmzJkEBwcDsGrVKtatW0dycrLtO/GNN95g3rx5fPfddzz00EP4+Pjg6upKaGjoRbN//fXXpKSksH79eqpWrQpA/fr1bcubN29O8+bNba9ffPFFfvjhB+bPn8/DDz/sqEP4r8pVcdOxY0cWLFhQpG3x4sV07NjRSYmK8qnfEdYtw//kJnILCvFwdXF2JBGRcs/LzYWd/+ecm0a83Bz/e7xZs2a25zVq1AAgOTmZ2rVrs3XrVv744w9efvll2zqFhYXk5OSQnZ3Nrl27CA8PL3IFIioqioCAAHbt2mUrburUqWMrbAC2bt1KZmYm1apVK5LlzJkzHDhwoNjZt2zZQsuWLW2Fzf/KzMzkhRde4JdffiEhIYGCggLOnDnDkSNHiv0ejuDU4iYzM5P9+/fbXsfFxbFlyxaqVq1K7dq1GTt2LMeOHWPmzJkADBs2jHfffZennnqK++67j6VLlzJnzhx++eUXZ32EIkKaXQ/rJtKB7WzYd5xOjUumf4+ISGViMpkcfmnImf45Z9LZu6z+eVlnwoQJ3Hrrredtd7bvTHFUqVKlyOvMzExq1KhR5NLVWQEBAcXe78Wukpw1ZswYFi9ezBtvvEH9+vXx8vLi9ttvv+JOzfZy6t+WDRs20KNHD9vrs31jhgwZwowZM0hISChS7dWtW5dffvmFxx9/nLfeeotatWrx8ccfl4nbwAHMNVuS6laDqvkJHF03HxqPdHYkEREpR1q1asWePXuKXOr5p8aNGxMfH098fLzt7M3OnTtJS0sjKirqkvtNTEzE1dWViIiIC67j7u7+r3M6NWvWjI8//pjU1NQLnr35448/uPfee7nlllsAa1F16NChS+6zJDi1uOnevfsle6JfaPTh7t27s3nz5hJMdQVMJtLr9qLq3k+pdmQhhjFCc6OIiEixPf/889x0003Url2b22+/HbPZzNatW9mxYwcvvfQSMTExNG3alIEDBzJ16lQKCgoYMWIE3bp1o02bNhfdb0xMDB07dqRfv3689tprNGzYkOPHj/PLL79wyy230KZNGyIiImxXUGrVqoWvr+95fVYHDBjAK6+8Qr9+/Zg4cSI1atRg8+bNhIWF0bFjRxo0aMDcuXPp06cPJpOJcePG2c5KlaZydSt4eRDSZSivFw5gQvatHEjJcnYcEREpR3r27MnPP//Mb7/9Rtu2benQoQNvvvkmderUAayXsX788UcCAwPp2rUrMTExREZGMnv27Evu12QysWDBArp27crQoUNp2LAhd911F4cPH7YNsXLbbbdxww030KNHD4KDg/nmm2/O24+7uzu//fYb1atXp3fv3jRt2pRJkybZ7siaMmUKgYGBREdH06dPH3r27EmrVq0cfJT+ncko7k38FURGRgb+/v6kp6fj51cyIwkP/nQdK/am8N/ejXioa70SeQ8RkYoqJyeHuLi4ImO0SOVwqT97e76/deamBFxzlbWHeuyuZCcnERERqXwqTvfzMuT6GmfY5fI7bvEW0rPb4O/t9u8biYiIiEPozE0JCHNJ41W36fyfy6es2faXs+OIiIhUKipuSkKdaI75NMFsMsjYPM/ZaURERCoVFTclxKh/HQA+SesoKCz92+BEREQqKxU3JaTG1V0AaGPsYGtcgpPTiIiIVB4qbkqIS0Q0p1yDCTalc2LVZ86OIyIiUmmouCkpbl4cu+peAKrHL3JuFhERkUpExU0JqtPhVtKNKiTleXIw+bSz44iIiJOZTCbmzZtX4u8TERHB1KlTy8x+SpuKmxLkW6sxD9f6nmH5j7NYA/qJiFR4KSkpDB8+nNq1a+Ph4UFoaCg9e/bkjz/+ACAhIYFevXo5OeX5ZsyYccHZwdevX89DDz1U+oGukAbxK0kmE9c1qcHKA6n8tjOJ/3TTVAwiIhXZbbfdRl5eHp9//jmRkZEkJSURGxvLyZMnAQgNDXVyQvsEBwc7O8Jl0ZmbEhbT2Doh2Z4jx0k5nevkNCIiUlLS0tJYuXIlr776Kj169KBOnTq0a9eOsWPHcvPNNwNFL0sdOnQIk8nEnDlz6NKlC15eXrRt25a9e/eyfv162rRpg4+PD7169SIlJcX2Pt27d+exxx4r8t79+vXj3nvvvWi2KVOm0LRpU6pUqUJ4eDgjRowgMzMTgGXLljF06FDS09MxmUyYTCZeeOEF4PzLUkeOHKFv3774+Pjg5+fHnXfeSVJSkm35Cy+8QIsWLfjiiy+IiIjA39+fu+66i9OnS7drhoqbEhaWG8da78dY6PYMsbuS/n0DERG5sLysiz/yc+xY90zx1rWTj48PPj4+zJs3j9zc4v9jdvz48Tz33HNs2rQJV1dX7r77bp566ineeustVq5cyf79+3n++eftzvNPZrOZt99+m7/++ovPP/+cpUuX8tRTTwEQHR3N1KlT8fPzIyEhgYSEBMaMGXPePiwWC3379iU1NZXly5ezePFiDh48SP/+/Yusd+DAAebNm8fPP//Mzz//zPLly5k0adIV5beXLkuVNP9aVLekYDIbrNm2i7va1XZ2IhGR8umVsIsva3A9DPz23OvX60N+9oXXrdMZhv5y7vXUppB98vz1Xki3K56rqyszZszgwQcfZNq0abRq1Ypu3bpx11130axZs4tuN2bMGHr27AnAqFGjGDBgALGxsXTq1AmA+++/nxkzZtiV5X/980xPREQEL730EsOGDeP999/H3d0df39/TCbTJS+bxcbGsn37duLi4ggPDwdg5syZXH311axfv562bdsC1iJoxowZ+Pr6AjBo0CBiY2N5+eWXr+gz2ENnbkqapx95VRsCkHd4LZm5BU4OJCIiJeW2227j+PHjzJ8/nxtuuIFly5bRqlWrSxYn/yx8QkKsXRmaNm1apC05+cpuSlmyZAnXXnstNWvWxNfXl0GDBnHy5Emysy9SAF7Arl27CA8PtxU2AFFRUQQEBLBr1y5bW0REhK2wAahRo8YV57eXztyUAveIDpC6h4dM81i/uRc9OrR1diQRkfLnv8cvvszkUvT1k/svse7//Lv+se2Xn+kCPD09ue6667juuusYN24cDzzwAOPHj79onxg3N7dz0UymC7ZZLOem8TGbzRiGUWQf+fn5F81z6NAhbrrpJoYPH87LL79M1apVWbVqFffffz95eXl4e3tfzse8qH9mv1D+0qAzN6XA1NB6urGF+QBtFt8BBepYLCJiN/cqF3+4edqxrlfx1nWQqKgosrLs78NzMcHBwSQknJvWp7CwkB07dlx0/Y0bN2KxWJg8eTIdOnSgYcOGHD9etFB0d3ensLDwku/buHFj4uPjiY+Pt7Xt3LmTtLQ0oqKiLvPTlAwVN6Xhqt4kNB0BgG/hKSxz7nVuHhERcbiTJ09yzTXX8OWXX7Jt2zbi4uL49ttvee211+jbt6/D3ueaa67hl19+4ZdffmH37t0MHz6ctLS0i65fv3598vPzeeeddzh48CBffPEF06ZNK7JOREQEmZmZxMbGcuLEiQteroqJiaFp06YMHDiQTZs2sW7dOgYPHky3bt1o06aNwz6fI6i4KQ0mE9VveYX3TXcxNv9+Vrea7OxEIiLiYD4+PrRv354333yTrl270qRJE8aNG8eDDz7Iu+++67D3ue+++xgyZIitsIiMjKRHjx4XXb958+ZMmTKFV199lSZNmvDVV18xceLEIutER0czbNgw+vfvT3BwMK+99tp5+zGZTPz4448EBgbStWtXYmJiiIyMZPbs2Q77bI5iMv73wl0Fl5GRgb+/P+np6fj5+ZXqez/93TYWbdhFv45XM6Fvk1J9bxGR8iInJ4e4uDjq1q2Lp6fnv28gFcal/uzt+f7WmZtSdP3VIaTjw287k87rDCYiIiKOoeKmFHWqH0Rt99OMy55E9gfXOjuOiIhIhaTiphR5urnQon44N5jXUyV5Ixxe7exIIiIiFY6Km1LWo2kd4oy/R4D89l4oyHNqHhERkYpGxU0pu6ZRCM9Z/p4+PjMRjm9ybiARkTJKfRMrH0f9mau4KWX+Xm541uvMr4V/jwkwZzDof2ARERsXF+tow3l5OrNd2Zz9Mz/7d+ByafoFJ+jVtAbL9jenp8sGyEyCxO1Q4+KTqomIVCaurq54e3uTkpKCm5sbZrP+HV4ZWCwWUlJS8Pb2xtX1ysoTFTdOcH1UCM/OvYb7LQupVSMETw/ff99IRKSSMJlM1KhRg7i4OA4fPuzsOFKKzGYztWvXts2xdblU3DhBgLc7HeoFM3DffxlyVVtGVK3r7EgiImWKu7s7DRo00KWpSsbd3d0hZ+pU3DhJ76Y1WLnvBAv+SmHENY2cHUdEpMwxm80aoVguiy5kOsn1USGYTbDjWAanlrwJX9wCx3TnlIiIyJVSceMk1Xw86BBZDYBTe1bCgaWw6ycnpxIRESn/VNw4Ua+mNQCYn9PC2nD4T+eFERERqSBU3DhRz6tDMJngjxNVrA2ZSc4NJCIiUgGouHGi6r6etI2oSjIB1obMZKfmERERqQhU3DhZ7yahpBgB1hf5WZB72ql5REREyjsVN052Q5MaZONJpvH37Y6nDjk1j4iISHmn4sbJQv09aV0nkL+MCE571QTD4uxIIiIi5ZoG8SsDejUJ5Ykjw6nrG84XNZo7O46IiEi5pjM3ZcCNzWpw1Ahm5ZEcjqedcXYcERGRck3FTRlQw9+LdhFVAYPf12+1DuhXkOvsWCIiIuWSipsyok9z64B+ff+81ToVw4l9Tk4kIiJSPqm4KSN6Na2Bi9nMrsKa1oaU3c4NJCIiUk6puCkjgnw8iK5XjR2WutaGuOXODSQiIlJOqbgpQ25uHsZvljYAGNu/h/SjTk4kIiJS/qi4KUOuvzqUTaar2WqJxJSfBTvnOzuSiIhIuaPipgzx93Kj61Uh/FZoPXvD0XXODSQiIlIOqbgpY25uHsYGy1UAGEfWgEUjFouIiNhDxU0Zc23j6vzl2pjH84azo++vYNYfkYiIiD30zVnGeLu7ck1UGD9YujBvd5az44iIiJQ7Km7KoD7NwwD4edtxCi2Gk9OIiIiULypuyqCuDYMI9zzDnVnfcGrmPc6OIyIiUq6ouCmDPFxduLZRdUa5ziXo0M9w+E9nRxIRESk3VNyUUTGto/i2sBsAhSvfdHIaERGR8kPFTRnVsV41fvTqC4DpwO+Qp87FIiIixeH04ua9994jIiICT09P2rdvz7p1lx64burUqVx11VV4eXkRHh7O448/Tk5OTimlLT0uZhPNW3bgmFENs5EPxzY5O5KIiEi54NTiZvbs2YwePZrx48ezadMmmjdvTs+ePUlOTr7g+l9//TXPPPMM48ePZ9euXXzyySfMnj2b//73v6WcvHTc0qommy31Acjd8aOT04iIiJQPTi1upkyZwoMPPsjQoUOJiopi2rRpeHt78+mnn15w/T///JNOnTpx9913ExERwfXXX8+AAQP+9WxPedUo1I91/jcA4LbxE0jY6uREIiIiZZ/Tipu8vDw2btxITEzMuTBmMzExMaxevfqC20RHR7Nx40ZbMXPw4EEWLFhA7969L/o+ubm5ZGRkFHmUJ+Ht+jGt4Cbme98KLu7OjiMiIlLmuRZnpYyMDPz8/GzPL+Xsev/mxIkTFBYWEhISUqQ9JCSE3bt3X3Cbu+++mxMnTtC5c2cMw6CgoIBhw4Zd8rLUxIkTmTBhQrEylUV9W4TRYeHdWFKhpUtt6jg7kIiISBlXrDM3gYGBtn4wAQEBBAYGnvc4216Sli1bxiuvvML777/Ppk2bmDt3Lr/88gsvvvjiRbcZO3Ys6enptkd8fHyJZnS06n6edKofBMAPm485OY2IiEjZV6wzN0uXLqVq1aoA/P777w5546CgIFxcXEhKSirSnpSURGho6AW3GTduHIMGDeKBBx4AoGnTpmRlZfHQQw/x7LPPYr7AJJMeHh54eHg4JLOz3NqqJn/tO8jBDb9idA7F5Fm8s2MiIiKVUbGKm27dul3w+ZVwd3endevWxMbG0q9fPwAsFguxsbE8/PDDF9wmOzv7vALGxcUFAMOouHMw9bw6lJYeE4jISWDvpto0jO7r7EgiIiJl1mV1KF65ciX33HMP0dHRHDtmvVTyxRdfsGrVKrv2M3r0aKZPn87nn3/Orl27GD58OFlZWQwdOhSAwYMHM3bsWNv6ffr04YMPPmDWrFnExcWxePFixo0bR58+fWxFTkXk7e5KckBzAFI3/eDkNCIiImVbsc7c/NP333/PoEGDGDhwIJs2bSI3NxeA9PR0XnnlFRYsWFDsffXv35+UlBSef/55EhMTadGiBYsWLbJ1Mj5y5EiRMzXPPfccJpOJ5557jmPHjhEcHEyfPn14+eWX7f0Y5Y5X6/6wdBGtT8wn78Qh3IMinB1JRESkTDIZdl7PadmyJY8//jiDBw/G19eXrVu3EhkZyebNm+nVqxeJiYklldUhMjIy8Pf3Jz09vdh3dpUFhYUWtr8UTQtjF3ubPkHD2553diQREZFSY8/3t92Xpfbs2UPXrl3Pa/f39yctLc3e3UkxubiYSQ63Duh35oBmCRcREbkYu4ub0NBQ9u/ff177qlWriIyMdEgoubCoq1sA4JqZwMnMXOeGERERKaPsLm4efPBBRo0axdq1azGZTBw/fpyvvvqKMWPGMHz48JLIKH+rVcc6z1So6aTGvBEREbkIuzsUP/PMM1gsFq699lqys7Pp2rUrHh4ejBkzhkceeaQkMspZ1eqxpO1HvLIqA9cN8dzfuS4mk8nZqURERMoUuzsUn5WXl8f+/fvJzMwkKioKHx8fR2crEeW1Q/FZ6WfyaffyEnILLPw4shPNwwOcHUlERKTElWiH4rPc3d2JioqiXbt25aawqQj8vdy4oUkoZix8u7F8TSUhIiJSGop1WerWW28t9g7nzp172WGkGAoLGO35E295TuY/W/5Lzo1ReLpV3AEMRURE7FWs4sbf37+kc0hx5WZQe+eHAIwzPuLXHXfTt2W4k0OJiIiUHcUqbj777LOSziHF5V0V0+N/UfhaJLU4wcd//kbflvc7O5WIiEiZYXefm7i4OPbt23de+759+zh06JAjMsm/8Qogp6F18szwhMXEp2Y7OZCIiEjZYXdxc++99/Lnn+ePkLt27VruvfdeR2SSYqjS8nYAerms5bsNR5ycRkREpOywu7jZvHkznTp1Oq+9Q4cObNmyxRGZpDjqXUu+axXCTKkcWz8fi+Wy7ugXERGpcOwubkwmE6dPnz6vPT09ncLCQoeEkmJw88TU8h4Absz5mdUHTzo5kIiISNlgd3HTtWtXJk6cWKSQKSwsZOLEiXTu3Nmh4eTSXLs/zTtXfc79+U8ya73GvBEREYHLmH7h1VdfpWvXrlx11VV06dIFgJUrV5KRkcHSpUsdHlAuoUo1enTtweStq1i0I4GTmVFU8/FwdioRERGnsvvMTVRUFNu2bePOO+8kOTmZ06dPM3jwYHbv3k2TJk1KIqNcQpOa/jSt6Y+lsIDvNxx2dhwRERGns/vMDUBYWBivvPKKo7PIZXrX/V2CPJbyzJrxGN0aaDJNERGp1OwublasWHHJ5V27dr3sMHJ5avqacTXlUuP0dlYfOEl0/SBnRxIREXEau4ub7t27n9f2zzMFumOq9LnWaQ97f6GleT9frTui4kZERCo1u/vcnDp1qsgjOTmZRYsW0bZtW3777beSyCj/plZbANqZd7P8ryOcyMx1ciARERHnsfvMzYUm0bzuuutwd3dn9OjRbNy40SHBxA5hrcAnhGqZSdzESr7d0JTh3es5O5WIiIhT2H3m5mJCQkLYs2ePo3Yn9nDzhGb9AbjadIhv1h3RiMUiIlJp2X3mZtu2bUVeG4ZBQkICkyZNokWLFo7KJfYKuRqAKJdjHEnN5o8DJ+jSINjJoUREREqf3cVNixYtMJlMGEbRMwMdOnTg008/dVgwsVNoMwhpwrqA0bAVvll3RMWNiIhUSnYXN3FxcUVem81mgoOD8fT0dFgouQwhUfDQcnqkZPPq1pX89lcSyadzqO6rPxcREalc7C5u6tSpUxI5xBFcXGkU6ker2gFsOpLGtxuOMrJHfWenEhERKVXFLm5mzpxZrPUGDx582WHEAeLX80LQ7ww8chVfrz3Cf7pG4urisH7jIiIiZZ7J+N/OMxcRGBh48Z2YTGRlZVFQUFDmB/HLyMjA39+f9PR0/Pz8nB3H8aY2g7TD/Mc8nl+zr2LaPa25oUmos1OJiIhcEXu+v4v9T/r/Hbzv7GPnzp3ceeedGIbBddddd8Xh5Qr9fdfU7bWzAJi5+pATw4iIiJS+y75ecfr0aZ577jkaNmzIli1b+PXXX1m0aJEjs8nlCL4KgOsOvYGbqZA/D5xkX9JpJ4cSEREpPXYXN/n5+UyZMoW6devy7bff8tlnn7FmzRp69OhREvnEXnW72Z4OqXf27M1hZ6UREREpdcUubgzD4PPPP6d+/fpMnjyZV155hZ07d3L77beXZD6xV70e4FsDgNsjLQB8v+koGTn5zkwlIiJSaopd3DRr1owRI0YwYMAANm7cyF133UVWVhYZGRlFHlIG1OkEwFXuJ6hf3YfsvELmbjzq5FAiIiKlo9jFzV9//cWZM2d47bXXqFmzJoGBgUUeAQEBl7yjSkpRtXrgXxtTeDuGdLSOSzRz9WHNNyUiIpVCsce5+f3330syhzhSm/sh+hHw8OWW0AJeXbSHgyeyWLX/BF0bakoGERGp2Ipd3HTr1u3fV5KywTfE9tTHw5XbW9dixp+HmLn6kIobERGp8DR0bUWWfwZOHWLQ35emYncnE5+a7eRQIiIiJUvFTUW1bzG8UhPeak4993S6NAjCMDSon4iIVHwqbiqq6lHA3x2IfxzB0Gjr2ZtZ6+LJzC1wXi4REZESpuKmovKvCf2/tD4/uIzu7nuIDK7C6dwC5qyPd242ERGREqTipiJrdCM0vhkAc9J27utUF4DP/oyjULeFi4hIBVWsu6VuvfXWYu9w7ty5lx1GSkD1xrBrPiTv4rbew3njtz3Ep55h8c5EbmhSw9npREREHK5YZ278/f1tDz8/P2JjY9mwYYNt+caNG4mNjcXf37/EgsplCm1m/blrPl6FGdzT3tr35uOVcU4MJSIiUnKKdebms88+sz1/+umnufPOO5k2bRouLi4AFBYWMmLECPz8/EompVy+hjdArXZQJRiyUxncsQ4frjjAhsOn2BKfRovwAGcnFBERcSiTYRh2db4IDg5m1apVXHXVVUXa9+zZQ3R0NCdPnnRoQEfLyMjA39+f9PT0SluMPTFnK99vOkqf5mG8M6Cls+OIiIj8K3u+v+3uUFxQUMDu3bvPa9+9ezcWi8Xe3YkT3N/Z2rF4wfYEjqWdcXIaERERxyr29AtnDR06lPvvv58DBw7Qrl07ANauXcukSZMYOnSowwOKgxgGJGwBi4WoWq2JrleNPw+c5PM/D/Hf3o2dnU5ERMRh7C5u3njjDUJDQ5k8eTIJCQkA1KhRgyeffJInnnjC4QHFQdZ9BAufgsgeMHgeD3Spy58HTvLN2iM8em0DfDzs/qsgIiJSJtnd5+afMjIyAMpV35VK2+fm5AF4pxWYXGDMPixeVYl5czkHU7IYd1OU7VKViIhIWVSifW7A2u9myZIlfPPNN5hMJgCOHz9OZmbm5exOSkO1ehDaFIxC2PIlZrOJBzpHAvDJyoPkF6q/lIiIVAx2FzeHDx+madOm9O3bl5EjR5KSkgLAq6++ypgxYxweUByo7YPWn6vehMJ8bm1VkyAfD46n5zB/y3HnZhMREXEQu4ubUaNG0aZNG06dOoWXl5et/ZZbbiE2Ntah4cTBWt5jHe/mzCk4uBxPNxfb5ahpyw9g0ZQMIiJSAdhd3KxcuZLnnnsOd3f3Iu0REREcO3bMYcGkBJhdoHEf6/OdPwAwsENtfD1c2ZecSezuZCeGExERcQy7ixuLxUJhYeF57UePHsXX19choaQENbrR+jN+HQB+nm4M7GCdkuGDZfu5gv7lIiIiZYLdxc3111/P1KlTba9NJhOZmZmMHz+e3r17OzKblITa0TByPQz/09Z0X6cI3F3NbDqSxvpDp5wYTkRE5MrZXdxMnjyZP/74g6ioKHJycrj77rttl6ReffXVksgojuTuDcENwcXN1lTdz5PbW9cCrGdvREREyrPLGuemoKCAWbNmsW3bNjIzM2nVqhUDBw4s0sG4rKq049z8i0Mnsrhm8jIsBiwc1YXGNXRsRESk7CjRcW6ysrJwdXXlnnvu4bXXXuP999/ngQceuOzC5r333iMiIgJPT0/at2/PunXrLrl+WloaI0eOpEaNGnh4eNCwYUMWLFhwWe9daSX9Bd/dD/MfsTVFBFWhV9MaAHy4/ICzkomIiFwxu4ubkJAQ7rvvPlatWnXFbz579mxGjx7N+PHj2bRpE82bN6dnz54kJ1/4rp28vDyuu+46Dh06xHfffceePXuYPn06NWvWvOIslUphHuz4DnYXLQqHd6sHwE/bEohPzXZGMhERkStmd3Hz5ZdfkpqayjXXXEPDhg2ZNGkSx49f3gBwU6ZM4cEHH2To0KFERUUxbdo0vL29+fTTTy+4/qeffkpqairz5s2jU6dORERE0K1bN5o3b35Z719pBTW0/sw+AZkptuYmNf3p0iCIQovBNJ29ERGRcsru4qZfv37MmzePY8eOMWzYML7++mvq1KnDTTfdxNy5cykoKCjWfvLy8ti4cSMxMTHnwpjNxMTEsHr16gtuM3/+fDp27MjIkSMJCQmhSZMmvPLKKxe8Nf2s3NxcMjIyijwqPfcqUD3K+vzHEUUWPdyjPgDfbjhKQvqZ0k4mIiJyxS5rbimA4OBgRo8ezbZt25gyZQpLlizh9ttvJywsjOeff57s7Etf1jhx4gSFhYWEhIQUaQ8JCSExMfGC2xw8eJDvvvuOwsJCFixYwLhx45g8eTIvvfTSRd9n4sSJ+Pv72x7h4eH2f9iKqMmt1p/7foOMc2fe2kdWo33dquQVWvhw+UEnhRMREbl8l13cJCUl8dprrxEVFcUzzzzD7bffTmxsLJMnT2bu3Ln069fPgTGtLBYL1atX56OPPqJ169b079+fZ599lmnTpl10m7Fjx5Kenm57xMfHOzxXudTxXGdi4tcWWfTotQ0A+HrdEZIzckozlYiIyBVztXeDuXPn8tlnn/Hrr78SFRXFiBEjuOeeewgICLCtEx0dTePGjS+5n6CgIFxcXEhKSirSnpSURGho6AW3qVGjBm5ubri4uNjaGjduTGJiInl5eedNCQHg4eGBh4eHHZ+wknDzhHYPwbqP4K95cPUttkXR9arRuk4gGw+f4sMVBxl3U5TzcoqIiNjJ7jM3Q4cOJSwsjD/++IMtW7bw8MMPFylsAMLCwnj22WcvuR93d3dat25dZLJNi8VCbGwsHTt2vOA2nTp1Yv/+/VgsFlvb3r17qVGjxgULG/kX0Y9YZwq/+Z0izSaTyXb25qu1hzmRmeuMdCIiIpfF7kH8srOz8fb2dsibz549myFDhvDhhx/Srl07pk6dypw5c9i9ezchISEMHjyYmjVrMnHiRADi4+O5+uqrGTJkCI888gj79u3jvvvu49FHH/3XYuosDeJXPIZh0O/9P9kan8Z/ukUyttelz8SJiIiUJHu+v+2+LPXPwiYnJ4e8vLwiy+0pGPr3709KSgrPP/88iYmJtGjRgkWLFtk6GR85cgSz+dzJpfDwcH799Vcef/xxmjVrRs2aNRk1ahRPP/20vR9D/sligfR4sBRANetYNyaTiVHX1ue+GRv4YvVh/tO1HlWr6OyYiIiUfXafucnKyuLpp59mzpw5nDx58rzll7otuyzQmZsL+PNd+O1ZiOoHd35uazYMgz7vrmLHsQxG9qjHkz0bOS+jiIhUaiU6/cJTTz3F0qVL+eCDD/Dw8ODjjz9mwoQJhIWFMXPmzMsOLU4U2tT689AqKDw3TpHJZOLRa6x9bz7/8zBp2XkX2lpERKRMsbu4+emnn3j//fe57bbbcHV1pUuXLjz33HO88sorfPXVVyWRUUpanWjwCrSOWHz4jyKLrosKoXENPzJzC/h4ZZyTAoqIiBSf3cVNamoqkZGRgLV/TWpqKgCdO3dmxYoVjk0npcPFDRrdZH2+a36RRda+N9azN5/9EcdJ3TklIiJlnN3FTWRkJHFx1n/BN2rUiDlz5gDWMzr/e0u4lCMNe1p/Hj5/6oueV4fQtKY/WXmFmnNKRETKvMsa52br1q0APPPMM7z33nt4enry+OOP8+STTzo8oJSSWu2sP5P/gl0/F1lkMpl44nrrZJszVx8mSaMWi4hIGWb3reCPP/647XlMTAy7d+9m48aN1K9fn2bNmjk0nJQi3xAI7wDHN0Hq32dnCvIAA1w96NYwmDZ1Atlw+BTvLt3Pi/2aODWuiIjIxVz23FJn1alTh1tvvVWFTUVw64fQ7wPoNMr6+thGmFQbds7HZDIxpudVAMxaf4T41EtPjCoiIuIsdhc3jz76KG+//fZ57e+++y6PPfaYIzKJswRGQNPbz71OOwwFOfDtEDi2kQ6R1ehcP4j8QoO3Y/c5LaaIiMil2F3cfP/993Tq1Om89ujoaL777juHhJIyok4n8AkBwwLzRkJBrq3vzfebjnIwJdPJAUVERM5nd3Fz8uRJ/P39z2v38/PjxIkTDgklZURAOIxYA1WCIWUXLH+NlrUDiWlcHYsBby7R2RsRESl77C5u6tevz6JFi85rX7hwoW38G6lAvKtC7zesz9dOg53zeaJbDQB+2nqcXQkZTgwnIiJyPrvvlho9ejQPP/wwKSkpXHPNNQDExsYyefJkpk6d6uh8UhZE9QW/WpBxFJZNovHtn3Jjsxr8si2B1xbt5rOh7ZydUERExMbu4ua+++4jNzeXl19+mRdffBGAiIgIPvjgAwYPHuzwgFIGmEzQ61XY+Jn1LE7Vuoy5PotfdyTy+54U/jxwguh6Qc5OKSIiAlzGrOD/lJKSgpeXFz4+Po7MVKI0K7iDbPqCA7/PYOiJgQTUbMi8EZ0wm03OTiUiIhVUic4KDlBQUMCSJUuYO3cuZ2uj48ePk5mpu2cqjR3fU+/0Bka5z2fb0XR+2Z7g7EQiIiLAZRQ3hw8fpmnTpvTt25eRI0eSkpICwKuvvsqYMWMcHlDKqE6PAtDLfQsAr/+6h7wCixMDiYiIWNld3IwaNYo2bdpw6tQpvLy8bO233HILsbGxDg0nZVh4B8CEd/4prvbJ5EhqNl+tPezsVCIiIvYXNytXruS5557D3d29SHtERATHjh1zWDAp49y9IaA2AL8UPERtUxLvLN3P6Zx8JwcTEZHKzu7ixmKxUFhYeF770aNH8fX1dUgoKSc6Pmx7OsNrKmlZOXy4/KATA4mIiFxGcXP99dcXGc/GZDKRmZnJ+PHj6d27tyOzSVnX/iG4+1tw88ajdmt8yObjVQdJTM9xdjIREanE7L4V/OjRo/Ts2RPDMNi3bx9t2rRh3759BAUFsWLFCqpXr15SWR1Ct4KXAMPAAO6YtpoNh09xW6taTL6zubNTiYhIBWLP9/dljXNTUFDArFmz2LZtG5mZmbRq1YqBAwcW6WBcVqm4KTlb4tPo994fAMx/uBPNagU4N5CIiFQYJV7clGcqbkpITjose5X4LYvpnjaOlnWC+HZYR0wmDewnIiJXzp7v72JNvzB//vxiv/nNN99c7HWlAnHzhq3fEJ6TSk+3bSw43JKftyXQp3mYs5OJiEglU6wzN2Zz0X7HJpOJ/93s7L/QL3QnVVmiMzcl6Lfn4M93APiiIIaPqgxj8Zhr8HRzcXIwEREp7xw+/YLFYrE9fvvtN1q0aMHChQtJS0sjLS2NhQsX0qpVKxYtWuSQDyDlVLdnILAuAINcl9Ancw7TV+jWcBERKV1297lp0qQJ06ZNo3PnzkXaV65cyUMPPcSuXbscGtDRdOamhJ1Ogk9iIO0IAP0KX+XDJ+8lxM/TycFERKQ8K9GJMw8cOEBAQMB57f7+/hw6dMje3UlF4xsCQxdiVAkG4G5jAa8u2u3kUCIiUpnYXdy0bduW0aNHk5SUZGtLSkriySefpF27dg4NJ+WUfy1MI9aw6/ZlPFXwH+ZuOsaW+DRnpxIRkUrC7uLm008/JSEhgdq1a1O/fn3q169P7dq1OXbsGJ988klJZJTyqEoQjZu05LZWtQB4/scdFFoq1agDIiLiJMW6Ffyf6tevz7Zt21i8eDG7d1svNzRu3JiYmBiNaSLneaZXI/jrB6KTNrAwdhQ3XXedsyOJiEgFp0H8pGSt/RAWPgVALO1oMeZnqvl4ODmUiIiUNw4fxO9/xcbGEhsbS3JyMhaLpciyTz/99HJ2KRVVg+uwHOiNee8CIi2HeW3RHl69vZmzU4mISAVmd5+bCRMmcP311xMbG8uJEyc4depUkYdIEVUjMfeZCkBdcxLZm2az6Yj+noiISMmx+8zNtGnTmDFjBoMGDSqJPFIR+VS3Du53Ko4X3T7j4bnt+PzRm3Exq4+WiIg4nt1nbvLy8oiOji6JLFJRmUww/A8KA+oSYMpiwMl3+XrtYWenEhGRCsru4uaBBx7g66+/LoksUpG5V8Gl/0z21byFpYWteP3XPaScznV2KhERqYDsviyVk5PDRx99xJIlS2jWrBlubm5Flk+ZMsVh4aSCqdGMyPs/Y897q8g4lsFLv+zkrbtaOjuViIhUMHYXN9u2baNFixYA7Nixo8gyjXMj/8bFbOKVW5rS770/+HHLcW5pWZPuV1V3diwREalANM6NlC7DgOSdrJ77DmMOd8AUEM5vo7vh7X5ZoxKIiEglUaITZ4pckfxsmH4NHZO+4Q/PUdyT+SlTl+xzdioREalAiv3P5ZYtW17wspO/vz8NGzZk1KhRREVFOTScVEDuVSD6UVg1BSwFDHP9mc4rr+Pm5mE0qenv7HQiIlIBFPuy1IQJEy7YnpaWxqZNm1izZg1Lly6lU6dODg3oaLosVYZ8cj3EryXXcGNgtVnMGnkNri46mSgiIuez5/vbYX1unn32WdasWUNsbKwjdldiVNyUITvnwxzrYJAfFdyIuedLPNAl0smhRESkLHJKn5u7776b7du3O2p3UhlE3Qz9vwKgt8taJv+2h/jUbCeHEhGR8s5hxY2Li8t5k2iK/KuGN2AZ9idP1PicM/kWnvpuGxZLpbqBT0REHMxhxc3cuXPVoVjs5+KKOfRqXr29JZ5uZlYfPMFX6444O5WIiJRjxb5b6u23375ge3p6Ohs3buSXX35h4cKFDgsmlUuE1xnm1/qGwmObuX3Bq3RvGEx4VW9nxxIRkXKo2MXNm2++ecF2Pz8/rrrqKlasWEHHjh0dFkwqGXcfGqT8hsl8hpDcYzz9/Ta+vL89Zs0cLiIidip2cRMXF1eSOaSyc/PEVK0eJO2gnmsKiw+c5Ot1R7inQx1nJxMRkXJGg4pI2REYAcDQKOtfy4kLdunuKRERsZuKGyk7/i5uOvqdoF1EVbLyCnlmru6eEhER+6i4kbIjvD0Api1f8caNNfF0M/PH/pPMXH3IublERKRcUXEjZUeD68CrKuRnU3vuzfy3VyMAJi7czb6k004OJyIi5UWxiptbb72VjIwMAGbOnElubm6JhpJKys0L+r4LjW6CoQsZ1DGCbg2DyS2wMGrWFvIKNEikiIj8u2IVNz///DNZWVkADB06lPT09BINJZVYoxvhrq/ArwYmw8I7jXcS6O3GzoQM3lyy19npRESkHCjWreCNGjVi7Nix9OjRA8MwmDNnzkUnrRo8eLBDA0olVpiP36+j+KD719y1CKYtP0D3hsG0j6zm7GQiIlKGFWtW8D///JPRo0dz4MABUlNT8fX1xWQ6f3A1k8lEamqq3SHee+89Xn/9dRITE2nevDnvvPMO7dq1+9ftZs2axYABA+jbty/z5s0r1ntpVvByZNkkWDYRqkfxfLXXmbk5nZoBXix8rAt+nm7OTiciIqXI4bOCR0dHs2bNGlJSUjAMg71793Lq1KnzHpdT2MyePZvRo0czfvx4Nm3aRPPmzenZsyfJycmX3O7QoUOMGTOGLl262P2eUk5E9rD+TN7Jc/6/UruqN8fSzvDCj385N5eIiJRpdt8tFRcXR3BwsMMCTJkyhQcffJChQ4cSFRXFtGnT8Pb25tNPP73oNoWFhQwcOJAJEyYQGRnpsCxSxtRuDw17AeB+dA1v9m+B2QRzNx/jp63HnRxORETKKruLmzp16pCens7kyZN54IEHeOCBB5gyZcpldTLOy8tj48aNxMTEnAtkNhMTE8Pq1asvut3//d//Ub16de6///5/fY/c3FwyMjKKPKQc6fmy9efRdbQOMfNwj/oA/Hfudo6c1OjFIiJyPruLmw0bNlCvXj3efPNNUlNTSU1N5c0336RevXps2rTJrn2dOHGCwsJCQkJCirSHhISQmJh4wW1WrVrFJ598wvTp04v1HhMnTsTf39/2CA8PtyujOFnVSKhaz/r8g0482qMubeoEcjq3gIe/2aTbw0VE5Dx2FzePP/44N998M4cOHWLu3LnMnTuXuLg4brrpJh577LESiHjO6dOnGTRoENOnTycoKKhY24wdO5b09HTbIz4+vkQzioOZTNBplPV5wxtwpZC3B7QkwNuNbUfTeXXRbufmExGRMqfYs4KftWHDBqZPn46r67lNXV1deeqpp2jTpo1d+woKCsLFxYWkpKQi7UlJSYSGhp63/oEDBzh06BB9+vSxtVksFluGPXv2UK9evSLbeHh44OHhYVcuKWNaD7EO7FfFegt4WAC8cXtzHpi5gU9WxdEhshrXRYVceh8iIlJp2H3mxs/PjyNHjpzXHh8fj6+vr137cnd3p3Xr1sTGxtraLBYLsbGxdOzY8bz1GzVqxPbt29myZYvtcfPNN9OjRw+2bNmiS04V2d+FDUc3wtf9idnxFE+2sRbYY77dyrG0M04MJyIiZYndZ2769+/P/fffzxtvvEF0dDQAf/zxB08++SQDBgywO8Do0aMZMmQIbdq0oV27dkydOpWsrCyGDh0KWAcFrFmzJhMnTsTT05MmTZoU2T4gIADgvHapoLJSYO8iAEZU+ZOk6o8yMzmSR7/ZzKyHOuDmounSREQqO7uLmzfeeAOTycTgwYMpKCgAwM3NjeHDhzNp0iS7A/Tv35+UlBSef/55EhMTadGiBYsWLbJ1Mj5y5Ahms76w5G8hV0Oft2Hx85iykvk/nsPiMYIvD3dmyuK9PH1DI2cnFBERJyvWCMUXkp2dzYEDBwCoV68e3t7eDg1WUjRCcQWxbwnMfxhOJwDwl6UOd+WNY+qQrlzbWP1vREQqGoePUHwh3t7eNG3alKZNm5abwkYqkAYx8NgOqNMJgGrerrQ27+Gx2Vs4dCLLyeFERMSZdL1Hyi8XVxjyMzy+k6pPrOd0+DWczilg2Jcbyc4rcHY6ERFxEhU3Ur6ZzeBfE3c3F94f2IrgKm7sTjzN2LnbucwrriIiUs6puJEKIyR1Iyt8xlLfnMCPW47z+Z+HnB1JREScQMWNVByr3sQrfT9fVP+agS5LeOmXnWw4ZP9M9SIiUr5d1t1S+/bt4/fffyc5Odk2QvBZzz//vMPClQTdLVWB7ZwPcwbZXqYaPjzp+l8mjrqf6n6eTgwmIiJXyp7vb7uLm+nTpzN8+HCCgoIIDQ3FZDKd25nJZPfkmaVNxU0FZrHA+o9h+xw4uh6ARYVtmRY6gVkPdcDTzcXJAUVE5HKVaHFTp04dRowYwdNPP31FIZ1FxU0lYBiw43v4/n7W0JQhOU/Qu2VdptzZvEgxLiIi5UeJjnNz6tQp7rjjjssOJ1LiTCZoejtE9aVOw+YUmD34YfMxpi0/6OxkIiJSCuwubu644w5+++23ksgi4li3TqdGr6d4oU8UgWTw2q+7WLwz6d+3ExGRcs3uuaXq16/PuHHjWLNmDU2bNsXNza3I8kcffdRh4USuiKsHBNZhUEfosaI/c0834rFZrnw/IppGobokKSJSUdnd56Zu3boX35nJxMGDZfvUv/rcVEKpcRjTupBR4ELL7HcJC6zCjyM7Uc3Hw9nJRESkmOz5/rb7zE1cXNxlBxNxCv9amExm/C1pXBOQxJJTNfjPFxv58oH2uoNKRKQCuqJB/AzD0BD3Uva5uEGtNgC80iYbXw9XNhw+xRPfbsVi0d9fEZGK5rKKm5kzZ9K0aVO8vLzw8vKiWbNmfPHFF47OJuI4EdbZw6uve5UvbvLGzcXEL9sSePXX3U4OJiIijmZ3cTNlyhSGDx9O7969mTNnDnPmzOGGG25g2LBhvPnmmyWRUeTKtX0Q/GpBXiYtYgcy9cYwAD5cfpAvVh9ybjYREXGoy+pQPGHCBAYPHlyk/fPPP+eFF14o831y1KG4Eju6AT6OAQy4/iXeOXMDkxfvxWyCjwa1ISYqxNkJRUTkIkq0Q3FCQgLR0dHntUdHR5OQkGDv7kRKT6028MAS6/OarXkYOHrqDLM3xPPIN5uZ/Z8ONKsV4MyEIiLiAHZflqpfvz5z5sw5r3327Nk0aNDAIaFESkytNtaHyYTJZOKlW5rQtWEwZ/ILuW/Geg6fzHJ2QhERuUJ2n7mZMGEC/fv3Z8WKFXTqZO2k+ccffxAbG3vBokekTDIMWP0ubtknef/uZ7nzwzXsTMhg0Cfr+G5YR80iLiJSjtl95ua2225j7dq1BAUFMW/ePObNm0dQUBDr1q3jlltuKYmMIo63+l347TlY9SY+++cz47621KnmzZHUbAZ9so707HxnJxQRkctkd4fi8k4digWAvGz45i6IWw5VguHxncRnFHDbB3+SfDqX1nUC+eL+dni7231yU0RESoDDZwXPyMgo8vxSD5Fywd0bBn4HVapDVgp82JVwUzIz72+Hn6crGw+fYviXm8grsDg7qYiI2KlYxU1gYCDJyckABAQEEBgYeN7jbLtIueHqDi0HWp+n7IK3mtMoZzufDW2Hl5sLy/em8MS3WynUKMYiIuVKsc65L126lKpVqwLw+++/l2ggkVJ1zfNgcoGVbwAmSN5J63ad+OCeVjw4cwM/bT2Oj4crr9zSBJPJ5Oy0IiJSDHb3uTly5Ajh4eHn/aI3DIP4+Hhq167t0ICOpj43ckE5GXAqDmo0tzX9tOUYj87egmHAvdERjO8TpQJHRMRJHN7n5p/q1q1LSkrKee2pqanUrVvX3t2JlA2efucKm7xs+PVZ+uwfx+u3NcFkghl/HuLlX3ZpolgRkXLA7ltBDMO44L9eMzMz8fTU2CBSARzfZL1VHLg9YSsuvT7i8QWJfLwqDjdXM0/1vEpncEREyrBiFzejR48GwGQyMW7cOLy9vW3LCgsLWbt2LS1atHB4QJFSF9EZOj8Oq96Ek/u5JW4CmX3fZtyPf/HBsgO4u5h5/LqGzk4pIiIXUeziZvPmzYD1zM327dtxd3e3LXN3d6d58+aMGTPG8QlFnCHmBfDwg9gJELecQU2Wk3dTV178eSdvxe7DzcXEw9douhERkbLI7g7FQ4cO5a233iq3nXHVoViKzTDgreaQdhhCmsCwVXy44iATF+4G4LGYBoy6toEuUYmIlIIS7VA8depUCgoKzmtPTU3VIH5SsZhM8OBScHEH31DIP8N/utXj5a6eDHCJ5cDSz3l90V/qZCwiUsbYXdzcddddzJo167z2OXPmcNdddzkklEiZUSUIYiaAV6B1VGNgYFgSE90+4R33d3lgTU82vTcEIy/byUFFROQsu4ubtWvX0qNHj/Pau3fvztq1ax0SSqRM6TgCbp1+7rVfTQhtBkBVUyatT/zIqumjsWgkYxGRMsHu4iY3N/eCl6Xy8/M5c+aMQ0KJlDn/7FcT2Q2GrYQhP7Ot4SMARCfP4r0vZ2mqBhGRMsDu4qZdu3Z89NFH57VPmzaN1q1bOySUSLlQtwvN7n6JI7VuwsVkUG/fpzz6zWZyCwqdnUxEpFKzexC/l156iZiYGLZu3cq1114LQGxsLOvXr+e3335zeECRsq527ydJnxPHxpON+WV7Auln8pk2qDU+Hnb/7yUiIg5g95mbTp06sXr1asLDw5kzZw4//fQT9evXZ9u2bXTp0qUkMoqUbWEt8H/sT7oPHoe3uwur9p/g7ulrOJmZ6+xkIiKVkt3j3JR3GudGStK2uEQ2fv4Uq/Mi2V+1OzPvb0etQO9/31BERC7Jnu/vKzpvnpOTQ15eXpE2FQxSmTU7Podm/Eh/d09uO1md2z8oZOb97WgY4uvsaCIilYbdl6Wys7N5+OGHqV69OlWqVCEwMLDIQ6RSa/8fCG6ENzks9BhLt6yF3DFtNeviUp2dTESk0rC7uHnyySdZunQpH3zwAR4eHnz88cdMmDCBsLAwZs6cWRIZRcoPV48iY+K84vYJtXP2cM/Ha/lxyzEnBhMRqTzs7nNTu3ZtZs6cSffu3fHz82PTpk3Ur1+fL774gm+++YYFCxaUVFaHUJ8bKRWZKfBJDJw6BEDDnM/Jw40nrmvIw9fU13xUIiJ2KtG5pVJTU4mMjASs/WtSU62n2zt37syKFSsuI65IBeQTDH3etr18qFMtACYv3suT320jr8DirGQiIhWe3cVNZGQkcXFxADRq1Ig5c+YA8NNPPxEQEODQcCLlWt2ucONk6PgwY25owov9mmA2wfcbj3DvZ+tIP5Pv7IQiIhWS3cXN0KFD2bp1KwDPPPMM7733Hp6enjz++OM8+eSTDg8oUm6ZTND2Aej5Mrh5MujMVxz0uJuVHo8Tf3And7y/kkMnspydUkSkwrnicW4OHz7Mxo0bqV+/Ps2aNXNUrhKjPjfiNMsmwbKJtpc/F7bnvy6jeffu1nRtGOzEYCIiZV+J9bnJz8/n2muvZd++fba2OnXqcOutt5aLwkbEqbqMgXYPgVcgBiZuclnLG4WvMfSzNXy88iCVbDxNEZESY1dx4+bmxrZt20oqi0jF5uIKvV+Hpw9huuZZAK532Ugr9vLSL7t44tut5ORr0k0RkStld5+be+65h08++aQksohUHl3GQHAjAG7pGYOL2cTPmw4xZNpSEtNznBxORKR8s3v6hYKCAj799FOWLFlC69atqVKlSpHlU6ZMcVg4kQrLZIK7Z8Py17i7y9VE1KrJoq8m838n3+edt+6hzcD/o2O9as5OKSJSLtld3OzYsYNWrVoBsHfv3iLLNDCZiB0CI6Df+wBE1w+iUZeasAKGFX7D7Z80ove11/Bgl0jMJhO4uIPZ7hOtIiKVUrHvljp48CB169Yt9wWM7paSMstioXDWPbjs/eX8ZfcugIhOpZ9JRKSMKJG7pRo0aEBKSortdf/+/UlKSrr8lCJSlNmMy81TMfxrnbfoQEqmEwKJiJRPxT5zYzabSUxMpHr16gD4+vqydetW21QM5YXO3EiZZymEvCx2JWYwes5W4lOzKXTx4sWYYG4LPICpSnVoEOPslCIipapE55YSkRJmdgFPPxpH1GLWI9cTHVWXM4Umdi75AtO84fDVbTD/EduknCIiUlSxixuTyXRef5vy3v9GpKzz93Ljw0GtGXdTFIvpyGnDy7pg00z47j4oLHBuQBGRMqjYd0sZhsG9996Lh4cHADk5OQwbNuy8W8Hnzp1rd4j33nuP119/ncTERJo3b84777xDu3btLrju9OnTmTlzJjt27ACgdevWvPLKKxddX6S8M5lM3N+5Lh0iqzLmy1xqp6/nEdd5+B3biOWbAZhrtYbuzzg7pohImVHsMzdDhgyhevXq+Pv74+/vzz333ENYWJjt9dmHvWbPns3o0aMZP348mzZtonnz5vTs2ZPk5OQLrr9s2TIGDBjA77//zurVqwkPD+f666/n2LFjdr+3SHlydZg/Ux8bQlabEbxYcA8A5v2/kbfN/n9QiIhUZFc8ceaVat++PW3btuXdd98FwGKxEB4eziOPPMIzz/z7v0YLCwsJDAzk3XffZfDgwectz83NJTc31/Y6IyOD8PBwdSiWcm3R9gRWzH2PqnmJZLr4U7/3owyMcsO09RvY/QsMXQSu7s6OKSLiMOWmQ3FeXh4bN24kJubcnR9ms5mYmBhWr15drH1kZ2eTn59P1apVL7h84sSJRc4shYeHOyS7iDPd0LQGjzz2LBvrPMCMvGt4bt4O3vj8e4j9Pzi2EX56VP1xRKTScmpxc+LECQoLCwkJCSnSHhISQmJiYrH28fTTTxMWFlakQPqnsWPHkp6ebnvEx8dfcW6RsqCGvxdfPdCe52+KwtPNzHvH6jLbuM66cOs38G4bWPE65GQ4N6iISCkr17eCT5o0iVmzZvHDDz/g6el5wXU8PDzw8/Mr8hCpKMxmE/d1rsuCR7vQqnYgT+cO5YuCvwv9U3Gw9CWYFA5n0pyaU0SkNDm1uAkKCsLFxeW8kY6TkpIIDQ295LZvvPEGkyZN4rfffqNZs2YlGVOkzIsM9uHbYdGM7dWIF3mA/+Q9xq90OLeCYbH+TNkDy1+H45udE1REpBQ4tbhxd3endevWxMbG2tosFguxsbF07Njxotu99tprvPjiiyxatIg2bdqURlSRMs/FbOI/3erxyyOdSQi7jv/kPErn3Km8UO0NDmb8PSZVym74/SX4qDu8Hw3fPwgFeU7NLSLiaE6/LDV69GimT5/O559/zq5duxg+fDhZWVkMHToUgMGDBzN27Fjb+q+++irjxo3j008/JSIigsTERBITE8nM1Nw7IgANQnyZOzya//ZuxAnXUGYcC+OGd9fw7tJ95PvWBte/L+Em/wXb58C6D50bWETEwZxe3PTv35833niD559/nhYtWrBlyxYWLVpk62R85MgREhISbOt/8MEH5OXlcfvtt1OjRg3b44033nDWRxApc1xdzDzUtR6LH+9G14bB5BVYeOO3vdz4XQab+6+DQfOg+d3WlddNh/wcp+YVEXEkp49zU9o0caZUNoZhMH/rcf7vp52czMrDZIIB7WrzVJfqBExrBgU5cOvH0OwOZ0cVEbmocjPOjYiUPJPJRN8WNVkyuht3tK6FYcDXa4/Q/f2trGs4GsPFA9z+nrMqL8s6K7mISDmm4kakkgis4s7rdzRn9kMdaBTqS1p2PnduasJN/t+zwSvautJ391nHxzlzyrlhRUSugIobkUqmfWQ1fn6kMxNuvho/T1f+Op7B7dNWM3r2FvJTj0DqQfj1WWfHFBG5bCpuRCohVxczQ6Ij+H1Md+5qG47JBHM3H+PBlDutK2z5yjqNg4hIOaTiRqQSq+bjwaTbmjFvRCeahwewLLchWyz1rAunX4Mx6x74YbhzQ4qI2EnFjYjQPDyAH4ZHM+XO5nzkca+t3bT7Jwp2/uS8YCIil0G3gotIETn5hcyJ/ZMjq38gr6CAfFxJajCAsb0a0SDE19nxRKSSsuf7W8WNiFzQycxc3o7dx1drj1BgMWhhPsCEoFhqxDxC9agu4Orh7IgiUomouLkEFTci9jmYksma2a9x94m3bG2nXatijrqJKu6u1oYOIyGovvX54dVgdoXwtk5IKyIVlT3f366llElEyqnIYB8iR04g4SdXjG1zCCuIx7cgFbbNPLfS1beeK26SdsCCMeBXC1oNgsjuUCUYqtVzSn4RqXx05kZE7LL7j/lsXf0bx06dAcDNxYRLqwHcGdOZIG9XmPsg7Pj+/A07PgzXvwQmUyknFpGKQJelLkHFjciVMwyDVftP8MZve9kanwaAh6uZAe1q81B0GGH58bDkBTh1CFIPWDeqUh2GrQJf66S4zOxnvXx150xw93bCpxCR8kTFzSWouBFxHMMwiN2VzDu/77cVOW4uJm5rVYth3eoREVQF8s/AnMGQGgePbDi38dSmkHYEOoyAGyY65wOISLmh4uYSVNyIOJ5hGPyx/yTv/r6PNQdTATCboE/zMEZ0r89VIT7WszhV657baM0HsOgZwAT9v4SgBta+Od5VnfIZRKRsU3FzCSpuRErWxsOpvLt0P7/vSbG1db8qmAe7RBJdrxqmf/a5mTfCOtXDP938DrQaXEppRaS8UHFzCSpuRErHjmPpvL9sPwt3JHL2t0yjUF8e6BJJn+Y18HB1gTNp8NXtcHK/dSbyWm2hWX9o+4A6HotIESpuLkHFjUjpOnwyi8/+OMScDfFk5xUCEOzrwZCOdRjYvg6BVdytK2adhEMrIKqftbA5dQgyUzRejogAKm4uScWNiHOkZ+fz9bojzPgzjqSMXAA83cz0a1GTezrUoUlN/3MrH/oDvroDXNwguFHRHZldocvjUD+mFNOLiLOpuLkEFTcizpVXYGHB9gSmrzzIX8czbO0tawcwqEMdejetgaeRAx90glNxF95Jh5FwwyullFhEygIVN5eg4kakbDAMg/WHTvHlmsMs3JFAfqH1V1Ggtxt3tg1nUBNvamVug3/+irLkQ2EBNO//9+tCMCzWMzwiUqGpuLkEFTciZU/K6VzmbIjnqzWHOZ6eA1i73XRtEMxdbcO5tnEI7q7mohtZLLB2Gvw6Fro9DYER4OELV/UGs0vpfwgRKVEqbi5BxY1I2VVQaOH3PSl8seYwK/aeu5W8ahV3+rWoyZ1ta9Eo9O//b08egHdanb+T6lfDvT9rvByRCkbFzSWouBEpHw6dyGL2hni+33iU5NO5tvZmtfy5o004NzcPw7/gJCwaC7kZkH4UUnaf28EtH527fCUi5Z6Km0tQcSNSvhQUWli57wRzNsSzZFeSrW+Ou6uZ66NC6NeiJl0bBlsvW50d9bjHs9BlDJjN/7J3ESkvVNxcgoobkfLrZGYu87YcZ876ePYknba1B3i70btpDfo2D6NtLW/MWMC9ivXS1Y7vARO0/w946v95kfJKxc0lqLgRKf8Mw2D7sXTmbT7OT9uOk/KPy1Zh/p70aRFGvxY1abTp/zCtn25d4OoFMS9YixyNfixS7qi4uQQVNyIVS6HFYM3Bk8zbfIxFOxI5nVtgW9YiGCa7vk+9U6uKbuThbx0IsPPj1teGAfFroXaHUkwuIvaw5/tbF6RFpFxzMZvoVD+I1+9ozvrnYvhgYCtuuDoUdxczW1Lg2oQR9Ml9iS3mq89tlJsOBXnW54YB390Hn/aEL26BlD3O+SAi4jA6cyMiFVJGTj6xu5JYuD2R5XtTyC0oJIyTuJvyCfbxoE3jenRvcRVtIqri8vuLsHLyuY19Qs497z4W2gwt/Q8gIkXostQlqLgRqXyycgtYtieFhTsS+H13Mll/T+AJ1hGRr2lYlQdcF9Bo59uYLPlFN+75CnQcaX0e+39QmAcu7taZy/3CSvFTiFRuKm4uQcWNSOWWk1/Iqn0nWLgjkSW7kkg/c66Y8TfncEPNXNrXrUq7ulWpFegNvqFQJci6wsthkJ91bmchTaDLExDRGXyql/InEalcVNxcgoobETmroNDCpiNpxO5KInZ3MvuTM4ssrxtUhR5XVadLwyDa162K98qJkH0CNs4ouqPgRjBoHvjVsPbhWf4aNLsDqkaW2mcRqehU3FyCihsRuZjDJ7NYujuZpbuTWXPwpG3AQAB3FzOt6wTSpWEQPcLgqvy/MC+fBMk7Iaov3DnTuqJhwIQA6/OwltDxYevPwAjNeSVyBVTcXIKKGxEpjszcAlbuTWHFvhRW7D3BsbQzRZZXreJOdL1qxNRxo6vbX1Rtd5d1QUEuvHSBS1SBEdDrdWh4fcmHF6mAVNxcgoobEbGXYRgcOpnNyr8LnTUHT5L5j/F0AMKretG+bjU6RFajY4iFmqe3wdIXrXNe5f19uavlIOj7rvV5fg5kJYPJDG7emuhT5F+ouLkEFTcicqXyCy1siU9j5b4TrNyXwraj6RRaiv4qrRngRfvIqnSoW43urCP42FJMdaKh5UA4tMp651X82r/XNlkva0XdXPofRqScUHFzCSpuRMTRMnML2HAolbVxqaw5eJLtR9Mp+J9ip4a/J63rBNKqdiDRfilc9dsgTGdOgWEBSz54+EFoU3Dzgnu+d9InESm7VNxcgoobESlpWbkFbDx8irVxJ1lzMJVtR9OKdE4G8HA107xWAK3DqzBi7wP4pv89MrK7L/z3qPV54g44uv7ib9Sw57mxdpJ3w5HV4BUIjW/WjOhS4ai4uQQVNyJS2s7kFbL5yCk2Hj7FpiOn2HQkrcj4OoFk0M68m+o+7kQE+2KO6kOzWv40O/IFbrHPX3zHg+dDZDfr8/Ufwy9PWJ/7hELI39NNePjAdf9n7dAM1ru58s/83dfH0/EfVqSE2PP97VpKmUREKi0vdxei6wcRXd86GKDFYnDwRJa10Pm74Pkt2Q8jA8gADuwEoKfLKYZ4tcfPy8368HTF19MNl7OzmntXO/cmAXUgqCGc2AuZidbHWYYF+n9pfZ6ZDJMbWp9HdLF2cj67r/rXasZ0qRB05kZEpAxIP5PPlvg0Nh85xY5j6Ww9mk7K6dzz1nMxm2gY4kuzmv5cXdOPxjX8uCrUFz9PN+sKB5ZaCxiA9HhYNRVaDITer1nbTiedK27+V6228MCSc68XPw8nD1x4XfcqcOtHl/dhRS6DLktdgoobESkvkjJy2HY0ne1H09h+LJ1tR9M5mZV3wXVrBXrRKNSPqBq+NKphLXpqV/XGBQtYCsDVw7qixQIFZ2DDZ3Ag1tqWshcyjoJfTRi989xOp18LxzZcOJxXIDx9yPo8L9taVBmFULM1+NdyzAEQ+QcVN5eg4kZEyivDMEhI/7vgOZbGroTT7ErIICE954Lre7m5cFWoL41r+NKgui8NQnyoX92HUD9PTP+8/GQYsPtna5HSvP+59l0/QVbKhcO4eFhvawfY9TPMfwTOpFpfN7kN2j4IdTo64FOLWKm4uQQVNyJS0aRl57E70Vro7E44za7EDPYknia3wHLB9X08XKlX3YcG1a3FTv1gHxqE+FAr0BsXs519brJOwrxhkJsJR/4suszVEzo/Dt2fucxPJnKOOhSLiFQiAd7udIi0jo58VqHFIO5EFrsTM9iVkMH+5Ez2JWdy+GQ2mbkFbI1PY2t8WpH9eLiaiQy2Fjx1q3kTEVTF+qhWhUBvt6Jne86qUg0Gfmt9nn4U/poHq6ZA9kkoyLFeErOFKrDOr6VOy1LCdOZGRKQSySuwcOhkFvuTM20Fz/7kTA6kZJJ3kTM9AH6ertQNqkKdataCp26QNxHV/i58qrgXXbmwAE4nWJ97+IJXACTthNn3QOoB611a3Z4CTNaZ0/1rltjnlYpDZ25EROSC3F3NNAzxpWGIb5H2QovB0VPZtqLn0MksDp3I5tDJLBLSc8jIKWDrUetdXP/L38uN8KpehAd6UyvQi1qB3oRXtf6s5e2FN8DBZdbCBuDQSusDrCMzj9qqubXEoXTmRkRELulMXiGHU88VO4dOZBF3IovDJ7NJzLhwZ+Z/qlbFnVqBXjT1OU1MwXKapy3B3WzgdSYBc34WRv+vMDW+qRQ+iZRn6lB8CSpuREQcJzuvgCOp2cSnnuHoqWyOnrL+PPs6I6fgotve67KIF9xmEmN5D8O/Fi29T3JjwWKqeLhSxcMVn3/89HA1Y2rUG2p3sG586hBs+PTiwerHQN2ujv2w4lS6LCUiIqXC292VRqF+NAq98JdN+pl8jp0teE6dK4DiU7OZl34TK3KaEW/4kp+SRU3zHnq4f3PR93pnQxZ/BkF1Pw9aGbsZsuetiwfz9D9X3Kx4HbZ8ff46VYLhhklQs9W5tpS91olML8TFA4LqX/w9pcxQcSMiIiXG38sNfy83osIuXPycySskMSOHhPQzZB7zY/PeY2TlFpCVW0hWXgFZuQWcybd2dF6RG8b6jJMAbDTlkOty43n7M5vA292FXes9SDi4nmBfD25P2U/r1IPnv3nqQfhhGDy87lzbl7daR3a+kODGMHKN9XlGAmyccfEPHtHpXHGVdcJ6lqnJbVCt3sW3EYdRcSMiIk7j5e5C3aAq1A2qAvW6QNcu562Tk19IUkYOzdJzSEzPIeV0Lsmn67LzdAtSMnNJzsglJTOXtOy/z7jkA1lAsnUaipWmDlSnUZF9Xm0+xAOuCzh60osJby6nahV3qlZx59k8bwLcqmI2mTCZTJhN/P0cTF6B2OZaP50Ayydd4pM9fa64yT4Jv79sfVS/+vxV7/wcghpYn6cdgcL/OXNUJRg81Y3CHupzIyIiFUJuQSEnMvNIOZ1reySfzvnH81xOZeeRmpnH6dyL9wW6FB8PV/y93GjgfoqBhXNxdzHj7mrGw9X6093VBQ8XM5m1r8Fo2BN/LzcC8pPw/aDFxXc6bBWENrU+f6c1nNxfdLnJBTqOgEY3netzlHUCNv89GarZBaL6QUD4ZX2m8kIdii9BxY2IiOQWFJKWnc/JzDxSs/JIzc4jNTP33POsPNuyU3+/tlzBt2WQKYPWnsep4uFCFXdXvN1d8XRzwdvdhdTAJnhUCcDXw5XbNtyNb3a89UwRJlzyT9v2Yen2DOYeY60vknbCB/8zvUXNNuee3/rRuUtgm2bCxs8vHq7vu1C9sfX5tm9h7bRzy8wu0PwuaHgDmMzgG3puWdYJKLzAXGcuHtbBHR1MHYpFREQuwcPVhRA/F0L8PIu1vsVikH4mn7Qz+daf2Xmk//08Pfuf7flknMkn7Uye7XVugYUThh+/nvGDMxfae8rfD3iZcUWWhHKSka4/4kkei5fksXr5r/h6uFLXPZ37Pa/Dy1xIx+yl1pX/McnpD+sPkFfNDW93V64+uI/Ii02ACpCXde55ZuL5k6XGrwUeB98a8MTuc+3fDICj6zhPZHcY/OPF368U6MyNiIhICcrJL7QVQmnZ1p+Zuflk5hSQkVNAZm4Bp3Osr0/nFHA6t8D6PPdcW8ElThuFkEoTc1yRtrWWxmRah0+krimBSNPxC27rajbzl2sUFg9/vD1cqeeSTD3TUTxdXfBwhd6pXxCWcwBMkOUWzBcdf8HLzQUvdxeuWzOEqqe22fZl+vs/ObU6c6b/t1T935Grr5AuS12CihsRESlPDMMgt8BCxt8FkLUYOvvIJzvPemdZdm4h2XmFZOcVkJVXSHZuwXmvs/IKOZNXSF7hxafacITmtfz58eHODt1nubss9d577/H666+TmJhI8+bNeeedd2jXrt1F1//2228ZN24chw4dokGDBrz66qv07t27FBOLiIiUDpPJhKebC55uLlT3/ff1iyOvwMKZs0XR2QIo1/ozO6+QrNwCcvILyfl7vZz8Qs7kn/1pbcstsBZKZ/5elptvsT7PK8TL3cUxQS+T04ub2bNnM3r0aKZNm0b79u2ZOnUqPXv2ZM+ePVSvXv289f/8808GDBjAxIkTuemmm/j666/p168fmzZtokmTJk74BCIiIuWL+993d/l7u5XI/p19Ucjpl6Xat29P27ZteffddwGwWCyEh4fzyCOP8Mwzz5y3fv/+/cnKyuLnn3+2tXXo0IEWLVowbdq089bPzc0lNzfX9jojI4Pw8HBdlhIRESlH7LksZb7k0hKWl5fHxo0biYmJsbWZzWZiYmJYvXr1BbdZvXp1kfUBevbsedH1J06ciL+/v+0RHl6xxwEQERGp7Jxa3Jw4cYLCwkJCQkKKtIeEhJCYmHjBbRITE+1af+zYsaSnp9se8fEXGVZbREREKgSn97kpaR4eHnh4eDg7hoiIiJQSp565CQoKwsXFhaSkpCLtSUlJhIaGXnCb0NBQu9YXERGRysWpxY27uzutW7cmNjbW1maxWIiNjaVjx44X3KZjx45F1gdYvHjxRdcXERGRysXpl6VGjx7NkCFDaNOmDe3atWPq1KlkZWUxdOhQAAYPHkzNmjWZOHEiAKNGjaJbt25MnjyZG2+8kVmzZrFhwwY++ugjZ34MERERKSOcXtz079+flJQUnn/+eRITE2nRogWLFi2ydRo+cuQIZvO5E0zR0dF8/fXXPPfcc/z3v/+lQYMGzJs3T2PciIiICFAGxrkpbZp+QUREpPwpN+PciIiIiDiaihsRERGpUFTciIiISIWi4kZEREQqFBU3IiIiUqE4/Vbw0nb25rCMjAwnJxEREZHiOvu9XZybvCtdcXP69GkAzQ4uIiJSDp0+fRp/f/9LrlPpxrmxWCwcP34cX19fTCaTQ/edkZFBeHg48fHxGkOnBOk4lw4d59Kh41x6dKxLR0kdZ8MwOH36NGFhYUUG972QSnfmxmw2U6tWrRJ9Dz8/P/2PUwp0nEuHjnPp0HEuPTrWpaMkjvO/nbE5Sx2KRUREpEJRcSMiIiIVioobB/Lw8GD8+PF4eHg4O0qFpuNcOnScS4eOc+nRsS4dZeE4V7oOxSIiIlKx6cyNiIiIVCgqbkRERKRCUXEjIiIiFYqKGxEREalQVNw4yHvvvUdERASenp60b9+edevWOTtSuTJx4kTatm2Lr68v1atXp1+/fuzZs6fIOjk5OYwcOZJq1arh4+PDbbfdRlJSUpF1jhw5wo033oi3tzfVq1fnySefpKCgoDQ/SrkyadIkTCYTjz32mK1Nx9kxjh07xj333EO1atXw8vKiadOmbNiwwbbcMAyef/55atSogZeXFzExMezbt6/IPlJTUxk4cCB+fn4EBARw//33k5mZWdofpcwqLCxk3Lhx1K1bFy8vL+rVq8eLL75YZO4hHefLs2LFCvr06UNYWBgmk4l58+YVWe6o47pt2za6dOmCp6cn4eHhvPbaa475AIZcsVmzZhnu7u7Gp59+avz111/Ggw8+aAQEBBhJSUnOjlZu9OzZ0/jss8+MHTt2GFu2bDF69+5t1K5d28jMzLStM2zYMCM8PNyIjY01NmzYYHTo0MGIjo62LS8oKDCaNGlixMTEGJs3bzYWLFhgBAUFGWPHjnXGRyrz1q1bZ0RERBjNmjUzRo0aZWvXcb5yqampRp06dYx7773XWLt2rXHw4EHj119/Nfbv329bZ9KkSYa/v78xb948Y+vWrcbNN99s1K1b1zhz5oxtnRtuuMFo3ry5sWbNGmPlypVG/fr1jQEDBjjjI5VJL7/8slGtWjXj559/NuLi4oxvv/3W8PHxMd566y3bOjrOl2fBggXGs88+a8ydO9cAjB9++KHIckcc1/T0dCMkJMQYOHCgsWPHDuObb74xvLy8jA8//PCK86u4cYB27doZI0eOtL0uLCw0wsLCjIkTJzoxVfmWnJxsAMby5csNwzCMtLQ0w83Nzfj2229t6+zatcsAjNWrVxuGYf2f0Ww2G4mJibZ1PvjgA8PPz8/Izc0t3Q9Qxp0+fdpo0KCBsXjxYqNbt2624kbH2TGefvppo3PnzhddbrFYjNDQUOP111+3taWlpRkeHh7GN998YxiGYezcudMAjPXr19vWWbhwoWEymYxjx46VXPhy5MYbbzTuu+++Im233nqrMXDgQMMwdJwd5X+LG0cd1/fff98IDAws8nvj6aefNq666qorzqzLUlcoLy+PjRs3EhMTY2szm83ExMSwevVqJyYr39LT0wGoWrUqABs3biQ/P7/IcW7UqBG1a9e2HefVq1fTtGlTQkJCbOv07NmTjIwM/vrrr1JMX/aNHDmSG2+8scjxBB1nR5k/fz5t2rThjjvuoHr16rRs2ZLp06fblsfFxZGYmFjkOPv7+9O+ffsixzkgIIA2bdrY1omJicFsNrN27drS+zBlWHR0NLGxsezduxeArVu3smrVKnr16gXoOJcURx3X1atX07VrV9zd3W3r9OzZkz179nDq1KkryljpJs50tBMnTlBYWFjkFz1ASEgIu3fvdlKq8s1isfDYY4/RqVMnmjRpAkBiYiLu7u4EBAQUWTckJITExETbOhf6czi7TKxmzZrFpk2bWL9+/XnLdJwd4+DBg3zwwQeMHj2a//73v6xfv55HH30Ud3d3hgwZYjtOFzqO/zzO1atXL7Lc1dWVqlWr6jj/7ZlnniEjI4NGjRrh4uJCYWEhL7/8MgMHDgTQcS4hjjquiYmJ1K1b97x9nF0WGBh42RlV3EiZM3LkSHbs2MGqVaucHaXCiY+PZ9SoUSxevBhPT09nx6mwLBYLbdq04ZVXXgGgZcuW7Nixg2nTpjFkyBAnp6s45syZw1dffcXXX3/N1VdfzZYtW3jssccICwvTca7kdFnqCgUFBeHi4nLe3SRJSUmEhoY6KVX59fDDD/Pzzz/z+++/U6tWLVt7aGgoeXl5pKWlFVn/n8c5NDT0gn8OZ5eJ9bJTcnIyrVq1wtXVFVdXV5YvX87bb7+Nq6srISEhOs4OUKNGDaKiooq0NW7cmCNHjgDnjtOlfm+EhoaSnJxcZHlBQQGpqak6zn978skneeaZZ7jrrrto2rQpgwYN4vHHH2fixImAjnNJcdRxLcnfJSpurpC7uzutW7cmNjbW1maxWIiNjaVjx45OTFa+GIbBww8/zA8//MDSpUvPO1XZunVr3NzcihznPXv2cOTIEdtx7tixI9u3by/yP9TixYvx8/M774umsrr22mvZvn07W7ZssT3atGnDwIEDbc91nK9cp06dzhvKYO/evdSpUweAunXrEhoaWuQ4Z2RksHbt2iLHOS0tjY0bN9rWWbp0KRaLhfbt25fCpyj7srOzMZuLfo25uLhgsVgAHeeS4qjj2rFjR1asWEF+fr5tncWLF3PVVVdd0SUpQLeCO8KsWbMMDw8PY8aMGcbOnTuNhx56yAgICChyN4lc2vDhww1/f39j2bJlRkJCgu2RnZ1tW2fYsGFG7dq1jaVLlxobNmwwOnbsaHTs2NG2/Owtytdff72xZcsWY9GiRUZwcLBuUf4X/7xbyjB0nB1h3bp1hqurq/Hyyy8b+/btM7766ivD29vb+PLLL23rTJo0yQgICDB+/PFHY9u2bUbfvn0veCtty5YtjbVr1xqrVq0yGjRoUOlvUf6nIUOGGDVr1rTdCj537lwjKCjIeOqpp2zr6DhfntOnTxubN282Nm/ebADGlClTjM2bNxuHDx82DMMxxzUtLc0ICQkxBg0aZOzYscOYNWuW4e3trVvBy5J33nnHqF27tuHu7m60a9fOWLNmjbMjlSvABR+fffaZbZ0zZ84YI0aMMAIDAw1vb2/jlltuMRISEors59ChQ0avXr0MLy8vIygoyHjiiSeM/Pz8Uv405cv/Fjc6zo7x008/GU2aNDE8PDyMRo0aGR999FGR5RaLxRg3bpwREhJieHh4GNdee62xZ8+eIuucPHnSGDBggOHj42P4+fkZQ4cONU6fPl2aH6NMy8jIMEaNGmXUrl3b8PT0NCIjI41nn322yK3FOs6X5/fff7/g7+QhQ4YYhuG447p161ajc+fOhoeHh1GzZk1j0qRJDslvMox/DOUoIiIiUs6pz42IiIhUKCpuREREpEJRcSMiIiIVioobERERqVBU3IiIiEiFouJGREREKhQVNyIiIlKhqLgRERGRCkXFjYiUKffeey/9+vVzdgwRKcdcnR1ARCoPk8l0yeXjx4/nrbfeQgOni8iVUHEjIqUmISHB9nz27Nk8//zzRWbP9vHxwcfHxxnRRKQC0WUpESk1oaGhtoe/vz8mk6lIm4+Pz3mXpbp3784jjzzCY489RmBgICEhIUyfPp2srCyGDh2Kr68v9evXZ+HChUXea8eOHfTq1QsfHx9CQkIYNGgQJ06cKHZWwzB44YUXqF27Nh4eHoSFhfHoo4866lCISAlScSMiZd7nn39OUFAQ69at45FHHmH48OHccccdREdHs2nTJq6//noGDRpEdnY2AGlpaVxzzTW0bNmSDRs2sGjRIpKSkrjzzjuL/Z7ff/89b775Jh9++CH79u1j3rx5NG3atKQ+oog4kGYFFxGnmDFjBo899hhpaWlF2u+9917S0tKYN28eYD1zU1hYyMqVKwEoLCzE39+fW2+9lZkzZwKQmJhIjRo1WL16NR06dOCll15i5cqV/Prrr7b9Hj16lPDwcPbs2UPDhg3/Nd+UKVP48MMP2bFjB25ubo750CJSKnTmRkTKvGbNmtmeu7i4UK1atSJnUUJCQgBITk4GYOvWrfz++++2Pjw+Pj40atQIgAMHDhTrPe+44w7OnDlDZGQkDz74ID/88AMFBQWO+kgiUoLUoVhEyrz/PXNiMpmKtJ29C8tisQCQmZlJnz59ePXVV8/bV40aNYr1nmfP8ixZsoTFixczYsQIXn/9dZYvX64zOSJlnIobEalwWrVqxffff09ERASurpf/a87Ly4s+ffrQp08fRo4cSaNGjdi+fTutWrVyYFoRcTRdlhKRCmfkyJGkpqYyYMAA1q9fz4EDB/j1118ZOnQohYWFxdrHjBkz+OSTT9ixYwcHDx7kyy+/xMvLizp16pRwehG5UipuRKTCCQsL448//qCwsJDrr7+epk2b8thjjxEQEIDZbP21N2PGjEsOKhgQEMD06dPp1KkTzZo1Y8mSJfz0009Uq1attD6GiFwm3S0lIpXS+PHjWb58OcuWLXN2FBFxMPW5EZFKaeHChbz77rvOjiEiJUBnbkRERKRCUZ8bERERqVBU3IiIiEiFouJGREREKhQVNyIiIlKhqLgRERGRCkXFjYiIiFQoKm5ERESkQlFxIyIiIhWKihsRERGpUP4fvFkwR6DPw1UAAAAASUVORK5CYII=\n", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "seconds = np.arange(1001)\n", + "fraction = np.exp(-4.59e-3 * seconds)\n", + "\n", + "plt.plot(seconds, fraction, linestyle='-', label='Theoretical')\n", + "plt.plot(seconds, decay(4.59e-3, n=500), linestyle='--', label='Simulation')\n", + "plt.xlabel('Time, s')\n", + "plt.ylabel('Fraction of Undecayed Nuclei')\n", + "plt.legend();" + ] + }, + { + "cell_type": "markdown", + "id": "10d2f74f-06d1-418e-a9e0-14265c58171a", + "metadata": {}, + "source": [ + "The simulation is in good agreement with the theoretical model, but this is some discrepency. If we increase the number of simulated nuclei, the simulation should align better with the theoretical values." + ] + }, + { + "cell_type": "code", + "execution_count": 19, + "id": "bae004b9-9cf0-4975-aaa5-b9b410b69cde", + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjcAAAGwCAYAAABVdURTAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAABo/0lEQVR4nO3dd3gU5d7G8e/uppNGSEiBQAIBIdJ7ERCJICiKBRERFNsLooIcG0cBPRawgFhQFEXEBqiIBUQhgNJ7lQ5BWhqEJCQhbXfePyLRCMQsbLIp9+e65jI7Ozt77wTZH888xWQYhoGIiIhIJWF2dgARERERR1JxIyIiIpWKihsRERGpVFTciIiISKWi4kZEREQqFRU3IiIiUqmouBEREZFKxcXZAcqazWbjxIkT+Pj4YDKZnB1HRERESsAwDM6cOUNYWBhmc/FtM1WuuDlx4gTh4eHOjiEiIiKX4OjRo9SuXbvYY6pccePj4wMUXBxfX18npxEREZGSSE9PJzw8vPB7vDhVrrg5dyvK19dXxY2IiEgFU5IuJepQLCIiIpWKihsRERGpVFTciIiISKVS5frciIhIxWCz2cjNzXV2DClDbm5u/zrMuyRU3IiISLmTm5tLXFwcNpvN2VGkDJnNZiIjI3Fzc7us86i4ERGRcsUwDOLj47FYLISHhzvkX/JS/p2bZDc+Pp46depc1kS7Km5ERKRcyc/PJysri7CwMLy8vJwdR8pQUFAQJ06cID8/H1dX10s+j8phEREpV6xWK8Bl35qQiufc7/zcn4FLpeJGRETKJa3/V/U46neu4kZEREQqFacWN7/99ht9+/YlLCwMk8nE/Pnz//U1y5cvp1WrVri7uxMVFcXMmTNLPaeIiIhUHE4tbjIzM2nevDlTp04t0fFxcXFcf/31dO/ena1btzJq1Cjuv/9+fv7551JOKiIicnmWL1+OyWQiNTXV2VHsUtLGh5KKiIhgypQpDjvfhTh1tFTv3r3p3bt3iY+fNm0akZGRTJo0CYDGjRuzcuVK3njjDXr16lVaMUvsVEYOyRk5NArRgpwiIlXJv/UVGT9+PFdffXXZhLlEzz33HPPnz2fr1q1F9sfHx1O9enXnhLpEFarPzZo1a4iJiSmyr1evXqxZs+air8nJySE9Pb3IVhp+/j2BTi8t5OPZc0vl/CIiUn7Fx8cXblOmTMHX17fIvscff9xp2S53lueQkBDc3d0dlKZsVKjiJiEhgeDg4CL7goODSU9P5+zZsxd8zYQJE/Dz8yvcwsPDSyVbyyATW9z+jwmnH+dk4olSeQ8RkarIMAyycvOdshmGUaKMISEhhZufnx8mk6nIPm9v78JjN23aRJs2bfDy8qJTp07s3bu3yLm+++47WrVqhYeHB/Xq1eP5558nPz+/8PkjR45w00034e3tja+vL7fffjuJiYmFzz/33HO0aNGCDz/8kMjISDw8PABITU3l/vvvJygoCF9fX6655hq2bdsGwMyZM3n++efZtm0bJpMJk8lU2Kf1n7eljh07xsCBAwkICKBatWq0adOGdevWAXDw4EFuuukmgoOD8fb2pm3btixZsqTkv2wHqfST+I0ZM4bRo0cXPk5PTy+VAqdmzWDiXMKItMZxYO33BN40zOHvISJSFZ3NsxI9zjl9K3f9rxdebo79qnzmmWeYNGkSQUFBDBs2jHvvvZdVq1YBsGLFCoYMGcJbb71Fly5dOHjwIA8++CBQcGvLZrMVFja//vor+fn5jBgxggEDBrB8+fLC9zhw4ADffPMN8+bNw2KxANC/f388PT356aef8PPz4/3336dHjx7s27ePAQMGsHPnThYtWlRYjPj5+Z2XPSMjg27dulGrVi2+//57QkJC2Lx5c+EyGRkZGfTp04eXXnoJd3d3Zs2aRd++fdm7dy916tRx6HUsToUqbkJCQopUpwCJiYn4+vri6el5wde4u7uXWXPayZAuRB6PgwNLABU3IiJyvpdeeolu3boB8PTTT3P99deTnZ2Nh4cHzz//PE8//TR33303APXq1eOFF17gySefZPz48cTGxrJjxw7i4uIK/6E+a9YsrrzySjZs2EDbtm2BgltRs2bNIigoCICVK1eyfv16kpKSCr8TX3/9debPn8/XX3/Ngw8+iLe3Ny4uLoSEhFw0+xdffEFycjIbNmwgICAAgKioqMLnmzdvTvPmzQsfv/DCC3z77bd8//33PPzww466hP+qQhU3HTt2ZOHChUX2LV68mI4dOzopUVEBzfvA8Vk0SF9Hbl4+bq4V6vKKiJRLnq4Wdv3POYNGPF0tDj9ns2bNCn8ODQ0FICkpiTp16rBt2zZWrVrFSy+9VHiM1WolOzubrKwsdu/eTXh4eJE7ENHR0fj7+7N79+7C4qZu3bqFhQ3Atm3byMjIoEaNGkWynD17loMHD5Y4+9atW2nZsmVhYfNPGRkZPPfccyxYsID4+Hjy8/M5e/YsR44cKfF7OIJTv30zMjI4cOBA4eO4uDi2bt1KQEAAderUYcyYMRw/fpxZs2YBMGzYMN555x2efPJJ7r33XpYuXcrcuXNZsGCBsz5CEZEtryFzoQc1TOls2bSClh26OzuSiEiFZzKZHH5ryJn+vmbSuVFWf7+t8/zzz3PLLbec97pzfWdKolq1akUeZ2RkEBoaWuTW1Tn+/v4lPu/F7pKc8/jjj7N48WJef/11oqKi8PT05LbbbrvsTs32cuqflo0bN9K9+18FwLm+MXfffTczZ84kPj6+SLUXGRnJggULeOyxx3jzzTepXbs2H374YbkYBg5gdnUnzqctTc6s4PT2BaDiRkRE7NCqVSv27t1b5FbP3zVu3JijR49y9OjRwtabXbt2kZqaSnR0dLHnTUhIwMXFhYiIiAse4+bm9q9rOjVr1owPP/yQlJSUC7berFq1invuuYebb74ZKCiqDh8+XOw5S4NTi5urr7662J7oF5p9+Oqrr2bLli2lmOrymBrEwOYVBCWswDAMrY0iIiIlNm7cOG644Qbq1KnDbbfdhtlsZtu2bezcuZMXX3yRmJgYmjZtyqBBg5gyZQr5+fk89NBDdOvWjTZt2lz0vDExMXTs2JF+/frx6quv0rBhQ06cOMGCBQu4+eabadOmDREREYV3UGrXro2Pj895fVYHDhzIyy+/TL9+/ZgwYQKhoaFs2bKFsLAwOnbsSIMGDZg3bx59+/bFZDIxduzYwlapslShhoJXBJGdbuEV6yD+c3YoB5MznB1HREQqkF69evHjjz/yyy+/0LZtWzp06MAbb7xB3bp1gYLbWN999x3Vq1ena9euxMTEUK9ePebMmVPseU0mEwsXLqRr164MHTqUhg0bcscdd/DHH38UTrFy6623ct1119G9e3eCgoL48ssvzzuPm5sbv/zyCzVr1qRPnz40bdqUiRMnFo7Imjx5MtWrV6dTp0707duXXr160apVKwdfpX9nMko6iL+SSE9Px8/Pj7S0NHx9S2cm4SEz1vPbvmTG9G7E/3WrXyrvISJSWWVnZxMXF1dkjhapGor73dvz/a2Wm1LQo1FNAGL3JDk5iYiISNVTebqflyPXRPmy1byCtsf2kZr5Jf7VKta01SIiIhWZWm5KQXh1Lya4zeBOSyxbNq5ydhwREZEqRcVNaXD14Jh/Qa/1rB0/ODmMiIhI1aLippRYrrwJgMYnfybfWvbD4ERERKoqFTelJLzzAHJxoR7H2bVjk7PjiIiIVBkqbkqJi5c/B6u1BODkpm+dnEZERKTqUHFTivKirgMg6MRSJycRERGpOlTclKLIzv2xGSby8vI4mJDi7DgiIuJkJpOJ+fPnl/r7REREMGXKlHJznrKm4qYU+dSsy8O15nBL7v9YvPe0s+OIiEgpS05OZvjw4dSpUwd3d3dCQkLo1asXq1YVTAsSHx9P7969nZzyfDNnzrzg6uAbNmzgwQcfLPtAl0nFTSnr2KwRAL/8nuDkJCIiUtpuvfVWtmzZwieffMK+ffv4/vvvufrqqzl16hQAISEh5y1GWZ4FBQXh5eXl7Bh2U3FTyq5tXLAg2cGjx0hOSnRyGhERKS2pqamsWLGCV155he7du1O3bl3atWvHmDFjuPHGG4Git6UOHz6MyWRi7ty5dOnSBU9PT9q2bcu+ffvYsGEDbdq0wdvbm969e5OcnFz4PldffTWjRo0q8t79+vXjnnvuuWi2yZMn07RpU6pVq0Z4eDgPPfQQGRkFizsvX76coUOHkpaWhslkwmQy8dxzzwHn35Y6cuQIN910E97e3vj6+nL77beTmPjXd9tzzz1HixYt+PTTT4mIiMDPz4877riDM2fOXPqFvQQqbkpZiJ8Hk6t/wya3YRxd+oGz44iIVFy5mRff8rLtOPZsyY61k7e3N97e3syfP5+cnJwSv278+PE8++yzbN68GRcXF+68806efPJJ3nzzTVasWMGBAwcYN26c3Xn+zmw289Zbb/H777/zySefsHTpUp588kkAOnXqxJQpU/D19SU+Pp74+Hgef/zx885hs9m46aabSElJ4ddff2Xx4sUcOnSIAQMGFDnu4MGDzJ8/nx9//JEff/yRX3/9lYkTJ15WfntpbakyEBQWictBGx5xS4Cxzo4jIlIxvRx28eca9IRBX/31+LUoyMu68LF1r4KhC/56PKUpZJ06/7jn0uyK5+LiwsyZM3nggQeYNm0arVq1olu3btxxxx00a9bsoq97/PHH6dWrFwAjR45k4MCBxMbG0rlzZwDuu+8+Zs6caVeWf/p7S09ERAQvvvgiw4YN491338XNzQ0/Pz9MJhMhISEXPUdsbCw7duwgLi6O8PBwAGbNmsWVV17Jhg0baNu2LVBQBM2cORMfHx8ABg8eTGxsLC+99NJlfQZ7qOWmDIR36AdAg+wdnEk96dwwIiJSam699VZOnDjB999/z3XXXcfy5ctp1apVscXJ3wuf4OCCrgxNmzYtsi8pKemyci1ZsoQePXpQq1YtfHx8GDx4MKdOnSIr6yIF4AXs3r2b8PDwwsIGIDo6Gn9/f3bv3l24LyIiorCwAQgNDb3s/PZSy00ZqBvVhMOm2kRwjO2r59O6z/3OjiQiUvH898TFnzNZij5+4kAxx/7j3/Wjdlx6pgvw8PDg2muv5dprr2Xs2LHcf//9jB8//qJ9YlxdXf+KZjJdcJ/N9tcyPmazGcMwipwjLy/vonkOHz7MDTfcwPDhw3nppZcICAhg5cqV3HfffeTm5jq8w/Dfs18of1lQy00ZMJlMxIdcXfBg949OzSIiUmG5Vbv45uphx7GeJTvWQaKjo8nMtL8Pz8UEBQURHx9f+NhqtbJz586LHr9p0yZsNhuTJk2iQ4cONGzYkBMnihaKbm5uWK3WYt+3cePGHD16lKNHjxbu27VrF6mpqURHR1/ipykdKm7KiH/r2wBonL6a3LOO+0MuIiLlw6lTp7jmmmv47LPP2L59O3FxcXz11Ve8+uqr3HTTTQ57n2uuuYYFCxawYMEC9uzZw/Dhw0lNTb3o8VFRUeTl5fH2229z6NAhPv30U6ZNm1bkmIiICDIyMoiNjeXkyZMXvF0VExND06ZNGTRoEJs3b2b9+vUMGTKEbt260aZNG4d9PkdQcVNGrmjZjRME4WXKYf9qrTUlIlLZeHt70759e9544w26du1KkyZNGDt2LA888ADvvPOOw97n3nvv5e677y4sLOrVq0f37t0venzz5s2ZPHkyr7zyCk2aNOHzzz9nwoQJRY7p1KkTw4YNY8CAAQQFBfHqq6+edx6TycR3331H9erV6dq1KzExMdSrV485c+Y47LM5isn45427Si49PR0/Pz/S0tLw9fUt0/f+6sNX2BKXiF+LvjzV/+J/EEVEqrLs7Gzi4uKIjIzEw8Pj318glUZxv3t7vr/VclOGArsM5QtrD+btt2KzVamaUkREpMyouClDnerXoJqbhcT0HLYft2/+BBERESkZFTdlyN3FQu8GXgyx/EzGzy86O46IiEilpHluyli/0NNcdfATUo77YdhexWRWfSkiIuJI+mYtY806XEOu4UKAkcbh/Y6dOEpEpDKpYuNdBMf9zlXclDFfbx8OejYB4Niar/7laBGRqsdiKZhtODc318lJpKyd+52f+zNwqXRbygmyG9wAO7bie3S5s6OIiJQ7Li4ueHl5kZycjKurK2bdvq8SbDYbycnJeHl54eJyeeWJihsniGodAzteJDL/AAcS04kKLtv5dkREyjOTyURoaChxcXH88ccfzo4jZchsNlOnTp3CNbYulYobJ/AJb0KuyQ1fzvLjulVE3djb2ZFERMoVNzc3GjRooFtTVYybm5tDWupU3DiDxZWTQR0JSFzNkb2bABU3IiL/ZDabNUOxXBLdyHSSav3eoG3eB0w71ZK4k1pIU0RExFFU3DiJX1h9WtSvBcBPO+P/5WgREREpKRU3TtS7SSgAy7cfcnISERGRykN9bpyod2ASLdzGkHfSwtGUqwgP8HJ2JBERkQpPLTdOVD24Lo3MR2luPsTK9eudHUdERKRSUHHjTN5BJNZoB0DejnlODiMiIlI5qLhxsmqtbgOg9ZnlHDud5dwwIiIilYCKGyfzbXEzNsxcaf6DlRs2OzuOiIhIhafixtmqBZLk3wKAjO3fOzeLiIhIJaDiphyo1vwmABqnreRoim5NiYiIXA4VN+WAT/ObWOHZg1nWnvy4XRP6iYiIXA4VN+VBQCRHr57Cz7a2fL/thLPTiIiIVGgqbsqJ3k1CcDGb2B2fzoGkM86OIyIiUmGpuCknqnu5cmfdVB6xzGPBlsPOjiMiIlJhqbgpLwyDMaef4z+uXxO/5RcMw3B2IhERkQpJxU15YTZjib4egGYZK/n9RLqTA4mIiFRMKm7KEbfoGwC41rKJBVv+cHIaERGRiknFTXkS0YVsjyCCTGlkb/0am023pkREROyl4qY8cXHDpe1QADrlrGDzkdNODiQiIlLxqLgpZ1yuLJituKt5Bz9vPuDkNCIiIhWPipvyJvhKznqHk4k7u3/fTL7V5uxEIiIiFYqKm/LGZMJlyDx6maezMjOctYdSnJ1IRESkQlFxUw651mzItc3CAfhu63EnpxEREalYVNyUUzc2DwMMVu88QHae1dlxREREKgwVN+VUO9t21niM4iXbmyzelejsOCIiIhWG04ubqVOnEhERgYeHB+3bt2f9+vXFHj9lyhSuuOIKPD09CQ8P57HHHiM7O7uM0pYds384oSRzlXkHizfucnYcERGRCsOpxc2cOXMYPXo048ePZ/PmzTRv3pxevXqRlJR0weO/+OILnn76acaPH8/u3bv56KOPmDNnDv/973/LOHkZCIwiO6gpLiYbvnELOJmR4+xEIiIiFYJTi5vJkyfzwAMPMHToUKKjo5k2bRpeXl7MmDHjgsevXr2azp07c+eddxIREUHPnj0ZOHDgv7b2VFQeLW4HoL95OT9sO+HkNCIiIhWD04qb3NxcNm3aRExMzF9hzGZiYmJYs2bNBV/TqVMnNm3aVFjMHDp0iIULF9KnT5+Lvk9OTg7p6elFtgqjxZ1YTa40Nx9i5/qlzk4jIiJSIbiU5KD09HR8fX0Lfy7OueP+zcmTJ7FarQQHBxfZHxwczJ49ey74mjvvvJOTJ09y1VVXYRgG+fn5DBs2rNjbUhMmTOD5558vUaZyp1og+Y37Ydn1FR1T5nMgaQBRNb2dnUpERKRcK1HLTfXq1Qv7wfj7+1O9evXztnP7S9Py5ct5+eWXeffdd9m8eTPz5s1jwYIFvPDCCxd9zZgxY0hLSyvcjh49WqoZHc294/8BcIN5DT9s1HIMIiIi/6ZELTdLly4lICAAgGXLljnkjQMDA7FYLCQmFh3mnJiYSEhIyAVfM3bsWAYPHsz9998PQNOmTcnMzOTBBx/kmWeewWw+v1Zzd3fH3d3dIZmdonYb9jZ+mPHb/DmxLZmR1xmYzSZnpxIRESm3SlTcdOvW7YI/Xw43Nzdat25NbGws/fr1A8BmsxEbG8vDDz98wddkZWWdV8BYLBYADMNwSK5yx2Si7i3/4/ddSziTlsv6wyl0qFfD2alERETKrUvqULxixQruuusuOnXqxPHjBcsDfPrpp6xcudKu84wePZrp06fzySefsHv3boYPH05mZiZDhw4FYMiQIYwZM6bw+L59+/Lee+8xe/Zs4uLiWLx4MWPHjqVv376FRU5l5OFqoU/TUAC+3azlGERERIpTopabv/vmm28YPHgwgwYNYvPmzeTkFMy/kpaWxssvv8zChQtLfK4BAwaQnJzMuHHjSEhIoEWLFixatKiwk/GRI0eKtNQ8++yzmEwmnn32WY4fP05QUBB9+/blpZdesvdjVDh31U2h2daP2LCjBdk3XYmHa+Ut5kRERC6HybDzfk7Lli157LHHGDJkCD4+Pmzbto169eqxZcsWevfuTUJCQmlldYj09HT8/PxIS0sr8ciu8sBY+jKm317hN2tT0vvP5YZmYc6OJCIiUmbs+f62+7bU3r176dq163n7/fz8SE1Ntfd0UkKmFncA0Nm8k+XrNjk5jYiISPlld3ETEhLCgQPnD0leuXIl9erVc0gouYCAepytfRUWk0GrIzNITK9862mJiIg4gt3FzQMPPMDIkSNZt24dJpOJEydO8Pnnn/P4448zfPjw0sgof/KMKehcfaN5Nd+t3+fkNCIiIuWT3R2Kn376aWw2Gz169CArK4uuXbvi7u7O448/ziOPPFIaGeWcup0541UHn6wjpGz4CqNHU0wmzXkjIiLyd3Z3KD4nNzeXAwcOkJGRQXR0NN7eFWNZgIraofic3GWv4fbri6y3XYHp3kW0jQhwdiQREZFSZ8/3t90tN+e4ubkRHR19qS+XS+TWehCpq95nydlWnF5/RMWNiIjIP5SouLnllltKfMJ58+ZdchgpAd8wDty5mg8+2IDXzgTG35SPt/sl16giIiKVTom+Ff38/Eo7h9ihdWQQ9QKrcehkJgu3x3N723BnRxIRESk3SlTcfPzxx6WdQ+xgMpno37oW2xZ/hrF8IbSd6uxIIiIi5YbdQ8Hj4uLYv3//efv379/P4cOHHZFJSmBAnXSmuU2hX8Zs4o4ec3YcERGRcsPu4uaee+5h9erV5+1ft24d99xzjyMySQkE1G/NUddI3E357Fn6qbPjiIiIlBt2Fzdbtmyhc+fO5+3v0KEDW7dudUQmKaHMRrcBEHL4O/KtNienERERKR/sLm5MJhNnzpw5b39aWhpWq9UhoaRk6nW/BxsmWhq7WbdZ602JiIjAJRQ3Xbt2ZcKECUUKGavVyoQJE7jqqqscGk6K5xZQmzjftgCcWjXLyWlERETKB7snSHnllVfo2rUrV1xxBV26dAFgxYoVpKens3TpUocHlOJVazsYYtfT8vRPJKRmEeLv5exIIiIiTmV3y010dDTbt2/n9ttvJykpiTNnzjBkyBD27NlDkyZNSiOjFCOk/W1kmqqRZPjzw+rtzo4jIiLidJe8tlRFVdHXlrqQBWt3MGL+EcL8PFjx1DVYzFpMU0REKpdSXVvqt99+K/b5rl272ntKuUw9Wkfj/0s8J9KyWb43iR6Ng50dSURExGnsLm6uvvrq8/aZTH+1FGjEVNnzcLVwa6vafL1yB8tWrKBH49ucHUlERMRp7O5zc/r06SJbUlISixYtom3btvzyyy+lkVFK4P7A31nv/hA3H5vIidSzzo4jIiLiNHa33FxoEc1rr70WNzc3Ro8ezaZNmm/FGUKjryJvkYnWpv18vXget/Uf5OxIIiIiTmF3y83FBAcHs3fvXkedTuzlG8rRyP4A1N31nmYsFhGRKsvulpvt24sONzYMg/j4eCZOnEiLFi0clUsuQa3rn4J3vqC1bScrt+yga5vmzo4kIiJS5uwublq0aIHJZOKfI8g7dOjAjBkzHBZM7OceGMEfPi2pe2YLp1e8D23edXYkERGRMmd3cRMXF1fksdlsJigoCA8PD4eFkkvn0Xk4LHqQzqk/cCz5JWoHVXd2JBERkTJld5+bunXrFtnCw8NV2JQjwW1vJcVcAz8y+e3XJc6OIyIiUuZK3HIza1bJFmYcMmTIJYcRB7C4sKfrOzy0KB2XPf7cmm/F3cXi7FQiIiJlpsTLL1SvfvHbGyaTiczMTPLz88v9JH6VcfmFf8qz2ujyyjIS0rOZMqAF/VrWcnYkERGRy2LP93eJb0v9c/K+c9uuXbu4/fbbMQyDa6+99rLDy+VztZgZ1L4OAJ+v3u/kNCIiImXrkue5OXPmDM8++ywNGzZk69at/PzzzyxatMiR2eQyDA49yrdu4+kZ/z7bj6U6O46IiEiZsbu4ycvLY/LkyURGRvLVV1/x8ccfs3btWrp3714a+eQS+VtTaGneT2fz73yy+g9nxxERESkzJe5QbBgGs2bNYty4ceTn5/Pyyy9z3333YbGos2q5FNkNw+RCtPkPDmxfzak+jajh7e7sVCIiIqWuxC03zZo146GHHmLgwIFs2rSJO+64g8zMTNLT04tsUk54B0H0jQDcySLmbDzq5EAiIiJlo8Sjpczmv+ogk8l03vOGYWAymTRaqjw5shZm9CLbcOVmjw/54ckbcbE4bDkxERGRMmPP93eJb0stW7bssoNJGQtvjy24GR6J2+mW8RNLdnfguiYhzk4lIiJSqkpc3HTr1q00c0hpMJkwd/g/+G4Eg10W89TqISpuRESk0tM9isquyW2cbXAD4/OHsupQCvsTzzg7kYiISKlScVPZuXrgOehzLI16Y2Dm49WHnZ1IRESkVKm4qSKGdo4EYN7mY6Rk5jo5jYiISOlRcVNFtA81M7H6d7zOG3yxTpP6iYhI5aXipoownU3l9uyvucGyjpWrV5CTX76H7IuIiFyqEo2WuuWWW0p8wnnz5l1yGClFAZFwxfWw53v6ZC/kh209ua11bWenEhERcbgStdz4+fkVbr6+vsTGxrJx48bC5zdt2kRsbCx+fn6lFlQun7nd/QDcYlnBFyt+p4TzN4qIiFQoJWq5+fjjjwt/fuqpp7j99tuZNm1a4bpSVquVhx56qPLP+FvRRXbFGtAA75T9RCf/xJqDrekUFejsVCIiIg5ld5+bGTNm8PjjjxdZMNNisTB69GhmzJjh0HDiYCYTlvYPADDUsogZKw44OZCIiIjj2V3c5Ofns2fPnvP279mzB5vN5pBQUopa3InVozr1zfH47f+WA0kZzk4kIiLiUCVefuGcoUOHct9993Hw4EHatWsHwLp165g4cSJDhw51eEBxMHcfLN2e4LfVq1iTHI37qjhevrmps1OJiIg4jN3Fzeuvv05ISAiTJk0iPj4egNDQUJ544gn+85//ODyglIKOI3APvoMTH6zlm03HeLznFQRUc3N2KhEREYcwGZcxZCY9PR2gQnUktmfJ9MrMMAxufGcVO46nMrpHfR69trGzI4mIiFyUPd/flzSJX35+PkuWLOHLL7/EZDIBcOLECTIy1H+jojCZTDzSxpNPXF/BuvodsnLznR1JRETEIey+LfXHH39w3XXXceTIEXJycrj22mvx8fHhlVdeIScnh2nTppVGTikFPVx3YrFsJ8KawNz1f3DPVfWdHUlEROSy2d1yM3LkSNq0acPp06fx9PQs3H/zzTcTGxvr0HBSuizN+pPr4kNdcxJbf/2OPKtGu4mISMVnd3GzYsUKnn32WdzcinZAjYiI4Pjx4w4LJmXAzQtz8wEA9Mz+iR+3n3ByIBERkctnd3Fjs9mwWs9fdPHYsWP4+Pg4JJSUHZd29wJwrXkTs5duxGbTkgwiIlKx2V3c9OzZkylTphQ+NplMZGRkMH78ePr06ePIbFIWgq8kP6wNriYrrVMWsmxvkrMTiYiIXBa7i5tJkyaxatUqoqOjyc7O5s477yy8JfXKK6+URkYpZS5tC1pv7nSJ5f1l+5ycRkRE5PLYPVqqdu3abNu2jdmzZ7N9+3YyMjK47777GDRoUJEOxlKBNLmFvPUfMfPYlWw9cooNh1NoGxHg7FQiIiKXxO6Wm8zMTFxcXLjrrrt49dVXeffdd7n//vsvubCZOnUqEREReHh40L59e9avX1/s8ampqYwYMYLQ0FDc3d1p2LAhCxcuvKT3lj+5euL6f0vJbPkAubjy3vKDzk4kIiJyyewuboKDg7n33ntZuXLlZb/5nDlzGD16NOPHj2fz5s00b96cXr16kZR04X4fubm5XHvttRw+fJivv/6avXv3Mn36dGrVqnXZWQQe7FofkwmW7kliT0K6s+OIiIhcEruLm88++4yUlBSuueYaGjZsyMSJEzlx4tKGEE+ePJkHHniAoUOHEh0dzbRp0/Dy8mLGjBkXPH7GjBmkpKQwf/58OnfuTEREBN26daN58+aX9P5SVKTnWZ6rs53rzWuZukytNyIiUjHZXdz069eP+fPnc/z4cYYNG8YXX3xB3bp1ueGGG5g3bx75+SWbxj83N5dNmzYRExPzVxizmZiYGNasWXPB13z//fd07NiRESNGEBwcTJMmTXj55ZcvODT9nJycHNLT04tschF7fuTuxIlMdn2Pjdt3cCBJy2mIiEjFc0lrSwEEBQUxevRotm/fzuTJk1myZAm33XYbYWFhjBs3jqysrGJff/LkSaxWK8HBwUX2BwcHk5CQcMHXHDp0iK+//hqr1crChQsZO3YskyZN4sUXX7zo+0yYMAE/P7/CLTw83P4PW1W0HAJ1OuFuyuMRyzymLjvg7EQiIiJ2u+TiJjExkVdffZXo6GiefvppbrvtNmJjY5k0aRLz5s2jX79+DoxZwGazUbNmTT744ANat27NgAEDeOaZZ4pdz2rMmDGkpaUVbkePHnV4rkrDbIaY8QDcbvmVbds2Ency08mhRERE7GP3UPB58+bx8ccf8/PPPxMdHc1DDz3EXXfdhb+/f+ExnTp1onHjxsWeJzAwEIvFQmJiYpH9iYmJhISEXPA1oaGhuLq6YrFYCvc1btyYhIQEcnNzz1sSAsDd3R13d3c7PmEVV6cDNOiFy/6fGWX5mqnL2vJ6f/VpEhGRisPulpuhQ4cSFhbGqlWr2Lp1Kw8//HCRwgYgLCyMZ555ptjzuLm50bp16yKLbdpsNmJjY+nYseMFX9O5c2cOHDiAzfbXAo/79u0jNDT0goWNXKJrngXgRssadm9dzZFTxd9iFBERKU/sLm7i4+N5//33adu27UWP8fT0ZPz48f96rtGjRzN9+nQ++eQTdu/ezfDhw8nMzGTo0KEADBkyhDFjxhQeP3z4cFJSUhg5ciT79u1jwYIFvPzyy4wYMcLejyHFCW0GV94CwOuWd/lg6e9ODiQiIlJydt+W8vLyKvw5Ozub3NzcIs/7+vqW+FwDBgwgOTmZcePGkZCQQIsWLVi0aFFhJ+MjR45gNv9Vf4WHh/Pzzz/z2GOP0axZM2rVqsXIkSN56qmn7P0Y8m9ixpNo+PHUlnrs2pLEsJgsalf3+vfXiYiIOJnJMAy7loHOzMzkqaeeYu7cuZw6deq854sbll0epKen4+fnR1paml2FWFU16MO1rDpwikHt6/DSzU2dHUdERKooe76/7b4t9eSTT7J06VLee+893N3d+fDDD3n++ecJCwtj1qxZlxxayqdHr2kAwFcbjxGfdtbJaURERP6d3cXNDz/8wLvvvsutt96Ki4sLXbp04dlnn+Xll1/m888/L42M4kTtLfuYUX0m3Y11mvdGREQqBLuLm5SUFOrVqwcU9K9JSUkB4KqrruK3335zbDpxvkPLuebsLzzpMps56//gaIpGTomISPlmd3FTr1494uLiAGjUqBFz584FClp0/jkkXCqBjiPA3Y/65ng6sZ03Y/c7O5GIiEixLmmem23btgHw9NNPM3XqVDw8PHjsscd44oknHB5QnMzDF1rcCcCTLrOZv/mI1pwSEZFyze7RUv/0xx9/sGnTJqKiomjWrJmjcpUajZa6BJkn4e1WkJ3GI7kPYzS5lXfubOXsVCIiUoWU6mipf6pbty633HJLhShs5BJVC4QOBRMl3uPyMz9uj2fXCa2uLiIi5ZPdxc2jjz7KW2+9dd7+d955h1GjRjkik5RHre8BsyutzftpYjrE5MX7nJ1IRETkguwubr755hs6d+583v5OnTrx9ddfOySUlEM+wdDpYU5e9RxHCWHJ7kS2HDnt7FQiIiLnsbu4OXXqFH5+fuft9/X15eTJkw4JJeVUzHMExjzGta0aAqj1RkREyiW7i5uoqCgWLVp03v6ffvqpcP4bqdxG9miAqwVW7D/JmoPnL8EhIiLiTHYvnDl69GgefvhhkpOTueaaawCIjY1l0qRJTJkyxdH5pLwxDMKPLeA331cYlPp/TFy0h/kPdcJkMjk7mYiICHAJxc29995LTk4OL730Ei+88AIAERERvPfeewwZMsThAaWcMZlg5zeEnt3PE27fMPxoGAt3JHB9s1BnJxMREQEuc56b5ORkPD098fb2dmSmUqV5bhwgYSdMuwow6JPzMpkB0Sx+rBtuLpc9s4CIiMgFlfo8N/n5+SxZsoR58+ZxrjY6ceIEGRmaubZKCGkCV/YDYKjHr/xxKosv1x9xbiYREZE/2V3c/PHHHzRt2pSbbrqJESNGkJycDMArr7zC448/7vCAUk61uhuAG13W4EYeb8Xu50x2npNDiYiIXEJxM3LkSNq0acPp06fx9PQs3H/zzTcTGxvr0HBSjkV2Bd9auOelc5/fJk5l5vLBb4ecnUpERMT+4mbFihU8++yzuLm5FdkfERHB8ePHHRZMyjmzBdo9CMAjLvNwIZ8PV8SRmJ7t5GAiIlLV2V3c2Gw2rFbrefuPHTuGj4+PQ0JJBdHuAbiiD563TqV53SDO5lmZskQT+4mIiHPZXdz07NmzyHw2JpOJjIwMxo8fT58+fRyZTco7t2ow8EtM9a5mTO9GAMzZcJT9iWecHExERKoyu4ubSZMmsWrVKqKjo8nOzubOO+8svCX1yiuvlEZGqQDaRARwXXQgNgNeXLDb2XFERKQKu6R5bvLz85k9ezbbt28nIyODVq1aMWjQoCIdjMsrzXNTCnIzYcVk8nZ8S8vkcWRYXfl4aFu6X1HT2clERKSSsOf72+4ZigFcXFy46667LimcVEJmV9g+B9e0o7zW4HeG72nBiz/u4qqoQFwtmthPRETKVomKm++//77EJ7zxxhsvOYxUUC5u0H4Y/PIMPc8uoqZXaw4mZ/LZ2j8Y2jnS2elERKSKKdFtKbO56L++TSYT/3zZuYUTLzSSqjzRbalSciYR3moJeZnsrP8gN/x+NX6erix//GqqV3P799eLiIgUw+HLL9hstsLtl19+oUWLFvz000+kpqaSmprKTz/9RKtWrVi0aJFDPoBUQD7BcONbAFx59AtaBLuQdjZPQ8NFRKTM2d2huEmTJkybNo2rrrqqyP4VK1bw4IMPsnt3+R4po5abUmQY8E4bOHWApHr9aLerPxazmZ9GdqFhsOZAEhGRS1eqC2cePHgQf3//8/b7+flx+PBhe08nlYnJBL1fAZOZmseXMbhBPlabwQs/7jrvNqaIiEhpsbu4adu2LaNHjyYxMbFwX2JiIk888QTt2rVzaDipgKJiYNDXMGob9/e7FjeLmRX7T7Jkd5Kzk4mISBVhd3EzY8YM4uPjqVOnDlFRUURFRVGnTh2OHz/ORx99VBoZpaKJ6gGe1alboxr3XlUwWur5H34nO698dzYXEZHKwe55bqKioti+fTuLFy9mz549ADRu3JiYmJjCEVMiABgGj/mvIM/7MB+d7sS7yw8y+tqGzk4lIiKV3CXNUFyRqUNxGdr1HcwdQr7Fk2uzXuS4pRa/jOpKRGA1ZycTEZEKptRnKI6NjSU2NpakpCRsNluR52bMmHEpp5TKqFFfiOyKS9xvfOQ9jV5nxvLcD7/z8T1t1conIiKlxu4+N88//zw9e/YkNjaWkydPcvr06SKbSCGzGW5+HzyrUy9vP6Ndv2H53mR+2ZX4768VERG5RHbflgoNDeXVV19l8ODBpZWpVOm2lBPs+h7mDsaGha45kzD86rJkdDc83SzOTiYiIhVEqc5zk5ubS6dOnS45nFRB0TdCvasxY+Uhr2UcTz3L1GUHnJ1KREQqKbuLm/vvv58vvviiNLJIZdZ+GAD9zctxJ5cPfjvEweQMJ4cSEZHKyO4OxdnZ2XzwwQcsWbKEZs2a4erqWuT5yZMnOyycVCINekL74bg0vY1Oi20s25vMmHk7mP1AB8xmdS4WERHHsbvPTffu3S9+MpOJpUuXXnao0qQ+N853NCWLnm/8xtk8K6/c2pQBbes4O5KIiJRz9nx/a54bcYoPfzvIiwv34OvhwpL/dKOmj4ezI4mISDlWqh2KRS5L6lH4YRT3Hn6c5mFepGfn878fdjk7lYiIVCIl7nPTsmXLC0685ufnR8OGDRk5ciTR0dEODSeVkDUXts/BnJfFtA5X0Tm+ET9uj+eWVolc0yjY2elERKQSKHFx069fvwvuT01NZfPmzbRs2ZKlS5fSuXNnR2WTyqhGfej1Mvw4itCtb/Jwhzm8tSaFsfN/p/1jNajmfkmTZouIiBRyWJ+bZ555hrVr1xIbG+uI05Ua9bkpB2xWmNYFkn4nr80DdP+9D8dOn+XezpGM66vWPxEROZ9T+tzceeed7Nixw1Gnk8rMbIFeLwHgumkGk67xAmDm6ji2HNESHiIicnkcVtxYLJbzFtEUuaj63aFBLzCstN83mZubh2Az4Imvt5OdZ3V2OhERqcAcVtzMmzdPHYrFPtf+D0wWSNjBc119CfR250BSBlOW7Hd2MhERqcBK3HvzrbfeuuD+tLQ0Nm3axIIFC/jpp58cFkyqgJqNYNBXUKsVfp7VeflmHx78dBMf/HaQXlcG07JOdWcnFBGRCqjEHYojIyMvuN/X15crrriCxx57jI4dOzo0XGlQh+LybdTsLczfeoKomt78+MhVeLhq5XAREbHv+7vELTdxcXGXHUzkovKyYeNHvJb5I6urjSq8PfV070bOTiYiIhWMZiiW8sGaC7+9huvR1XzUuqCQ/uC3g2w9murcXCIiUuGouJHywcMXOj0KQNOtL3J/tA2bAY9/tU2jp0RExC4qbqT8aPcg+NaC3DM8afmSIJ+C0VOv/bzX2clERKQCUXEj5Ye7N9zxOQBu+xcwq9UBAD5aGcfK/SedmUxERCqQEhU3t9xyC+np6QDMmjWLnJycUg0lVVhYS+j2FACNNz/PIy0L+rz/56utpGblOjOZiIhUECUqbn788UcyMzMBGDp0KGlpaaUaSqq4bk9B9UjwDGBE9wbUC6pGYnoOz3y7EwcthSYiIpVYiYaCN2rUiDFjxtC9e3cMw2Du3LkXHWM+ZMgQhwaUKshsgTvngNkFjxr1mDIggFveXc2CHfH02FKTW1rVdnZCEREpx0o0id/q1asZPXo0Bw8eJCUlBR8fH0wm0/knM5lISUmxO8TUqVN57bXXSEhIoHnz5rz99tu0a9fuX183e/ZsBg4cyE033cT8+fNL9F6axK9imvXTb4z/NY1q7m78NLIL4QFezo4kIiJlyOGrgnfq1Im1a9eSnJyMYRjs27eP06dPn7ddSmEzZ84cRo8ezfjx49m8eTPNmzenV69eJCUlFfu6w4cP8/jjj9OlSxe731MqmFVvMnjjbTxRcyMZOfmMnrsVq023p0RE5MLsHi0VFxdHUFCQwwJMnjyZBx54gKFDhxIdHc20adPw8vJixowZF32N1Wpl0KBBPP/889SrV6/Y8+fk5JCenl5kkwrGZMFkzeH/cj+hlvtZNhw+zXvLDzg7lYiIlFN2Fzd169YlLS2NSZMmcf/993P//fczefLkS+pknJuby6ZNm4iJifkrkNlMTEwMa9asuejr/ve//1GzZk3uu+++f32PCRMm4OfnV7iFh4fbnVOcrP3/Qc0rsWSf5rO6CwGYvHgf6+PsbykUEZHKz+7iZuPGjdSvX5833niDlJQUUlJSeOONN6hfvz6bN2+261wnT57EarUSHBxcZH9wcDAJCQkXfM3KlSv56KOPmD59eoneY8yYMaSlpRVuR48etSujlAMWV7hhMgCRR77hqYYJ2Ax49MstpGRqeLiIiBRld3Hz2GOPceONN3L48GHmzZvHvHnziIuL44YbbmDUqFGlEPEvZ86cYfDgwUyfPp3AwMASvcbd3R1fX98im1RAdTpAy8EADIsfS6eADBLSs3n8q20aHi4iIkWUeFXwczZu3Mj06dNxcfnrpS4uLjz55JO0adPGrnMFBgZisVhITEwssj8xMZGQkJDzjj948CCHDx+mb9++hftsNlthhr1791K/fn27MkgF0uslOLEFU+JOZvi8QWuX/7J0TxIfrYzj/i7F970SEZGqw+6WG19fX44cOXLe/qNHj+Lj42PXudzc3GjdujWxsbGF+2w2G7GxsXTs2PG84xs1asSOHTvYunVr4XbjjTfSvXt3tm7dqv40lZ2HH9z6EQQ1xqPrSMZc3xSAiT/t0erhIiJSyO6WmwEDBnDffffx+uuv06lTJwBWrVrFE088wcCBA+0OMHr0aO6++27atGlDu3btmDJlCpmZmQwdOhQomBSwVq1aTJgwAQ8PD5o0aVLk9f7+/gDn7ZdKqmYjeGgNmEwMMgzWHDrNgh3xPPzFZhY82gU/T1dnJxQRESezu7h5/fXXMZlMDBkyhPz8fABcXV0ZPnw4EydOtDvAgAEDSE5OZty4cSQkJNCiRQsWLVpU2Mn4yJEjmM1a31P+5s8JJE0mExP6NWL78VSOppzlqa+3895drS44waSIiFQdJZqh+EKysrI4ePAgAPXr18fLq2LMGKsZiiuR7XNhyXPsu/o9rv8mkzyrwTN9GvNAV/W/ERGpbBw+Q/GFeHl50bRpU5o2bVphChupZA7EQvpxGm4Yy/g+DQCYuGgPaw6ecnIwERFxJt3vkYrr2ufBszrEb2NQ7lfc0rIWVpvBw19sJj7trLPTiYiIk6i4kYrLJwSunwSAadWbvNzdj8ahvpzKzGX4Z5vJybc6OaCIiDiDihup2K68BSK7Qn42Ht/cxfRbI/DzdGXr0VT+98MuZ6cTEREnUHEjFZvJBH3fBO9gSNxJ7R8H8vYt9TGZ4PN1R/hqo5bbEBGpauweCg6wf/9+li1bRlJSUuEMweeMGzfOIcFESiygHtz9I3x2K0R0oWv96ozs0YApS/bz7PydNArxpWltP2enFBGRMmL3UPDp06czfPhwAgMDCQkJKTKniMlksnvxzLKmoeCVmM0KZkvBjzaD+2dtZOmeJEL9PPhuRGdq+no4OaCIiFwqe76/7S5u6taty0MPPcRTTz11WSGdRcVNFXH2NOm5Bjd/uJ2DyZm0CPdn9oMd8HC1ODuZiIhcglKd5+b06dP079//ksOJlLqMZJjZF99vh/DRXc0KOxiPmbdDK4iLiFQBdhc3/fv355dffimNLCKOceYEnD4Mh1cQsfq/vHdnSyxmE99uOc60Xw85O52IiJQyuzsUR0VFMXbsWNauXUvTpk1xdS26UOGjjz7qsHAilyS0Odw+Ez6/HbZ9SafABjx3422Mnb+TV3/eQ1RNb66NDnZ2ShERKSV297mJjIy8+MlMJg4dKt//Mlafmypk/XRY+HjBz7fPYuy++ny69g+83Cx8M7wTjUP1+xcRqShKtUNxRafipopZ+CSsfx9cPMkf8gNDfray+uApavl78u1DnTSCSkSkgiiThTMBDMNQB00p33q9DFHXQv5ZXH5+infvbEm9wGocTz3L0JkbyMzJd3ZCERFxsEsqbmbNmkXTpk3x9PTE09OTZs2a8emnnzo6m8jls7jAbTOgyW3Qfyb+1dz5eGhbalRz4/cT6Yz4YjP5Vtu/n0dERCoMu4ubyZMnM3z4cPr06cPcuXOZO3cu1113HcOGDeONN94ojYwil8fDF277CPzrAFC3RjVm3h6Jh6uZ5XuTeXb+TrVAiohUIpfUofj5559nyJAhRfZ/8sknPPfcc8TFxTk0oKOpz00VZ82D2Odh7TQ2d5rKbbHe2Ax4otcVjOge5ex0IiJyEaXa5yY+Pp5OnTqdt79Tp07Ex8fbezqRsmV2gdwssOXRasMTvNIzEIDXft7Lt1uOOTmciIg4gt3FTVRUFHPnzj1v/5w5c2jQoIFDQomUGpMJer0EoS0gJ43+WXN5sGs9AJ78ejsr9ic7N5+IiFw2u29LffPNNwwYMICYmBg6d+4MwKpVq4iNjWXu3LncfPPNpRLUUXRbSgCIWwGf3ABmV2x3L+CRVa4s2B6Pl5uFz+9vT8s61Z2dUERE/qZUb0vdeuutrFu3jsDAQObPn8/8+fMJDAxk/fr15b6wESkU2QUa9wVbHuY5A5l8UyRdGgSSlWvlno83sC/xjLMTiojIJdIkflJ15ZyBaV3gdBw0vpHMfh9z10fr2HIklWBfd74e1onwAC9npxQREUqh5SY9Pb3Iz8VtIhWGuw/c+BZ4VocOw6nm7sLH97SlYbA3iek5DP5oHclncpydUkRE7FSilhuLxUJ8fDw1a9bEbDZjMpnOO8YwDEwmE1artVSCOopabuQ8OWcKCp0/JaZnc+t7qzl2+izRob7M/r8O+Hq4FnMCEREpbfZ8f5doVfClS5cSEBAAwLJlyy4/oUh58rfChuObCd75DZ/dO4bb3l/Lrvh0hn68gU/ubYe3e4n+dxERESezu8/NkSNHCA8PP6/1xjAMjh49Sp06dRwa0NHUciMXlZUCbzSBvEy48W12hfRj4PS1pJ3No11kADOHtsXLTQWOiIgzlOpoqcjISJKTz58LJCUlhcjISHtPJ1J+eAVAxxEFP//4GNG2fXx6Xzt83F1YH5fC/Z9s5Gxu+b7tKiIil1DcnOtb808ZGRl4eHg4JJSI01w9BqL7gS0fvh1Gs1BvPrmvHdXcLKw+eIoHP91Idp4KHBGR8qzEbeyjR48GwGQyMXbsWLy8/hoia7VaWbduHS1atHB4QJEyZTZD3zch7jc4tR+W/o9W1/6Pmfe24+4Z61mx/yTDP9vEtMGtcXexODutiIhcQIlbbrZs2cKWLVswDIMdO3YUPt6yZQt79uyhefPmzJw5sxSjipQRT3/o/WrBz6vehI0f0zYigI/ubouHq5lle5MZ8fkWcvLVgiMiUh6VuOXm3CipoUOH8uabb6ozrlRuzfrDqQPw60TIKZi/qWP9Gnw4pC33frKBJbsTGfbpJt67qzUermrBEREpT+weLZWWlobVai0cGn5OSkoKLi4u5b7o0WgpKTHDKGi56TgCLH/Nc7NifzIPzNpIdp6Nq6IC+WBIa42iEhEpZaU6WuqOO+5g9uzZ5+2fO3cud9xxh72nEym/TCa4atRfhY3NCgk76dIgiJlD2+HlZmHlgZPcM2MDGTn5To0qIiJ/sbu4WbduHd27dz9v/9VXX826descEkqk3MlOh89vg496woFYOtSr8dcw8cMpDP5oHWln85ydUkREuITiJicnh/z88/+VmpeXx9mzZx0SSqTccfUCw1Ywwd9nt8CiMbQO8+TzB9rj5+nKliOpDPpwLaczc52dVESkyrO7uGnXrh0ffPDBefunTZtG69atHRJKpNyxuMCAz6D1PQWP174LM66jmW0PXz7QgYBqbuw8nk7/99dwIlVFvoiIM9ndoXjVqlXExMTQtm1bevToAUBsbCwbNmzgl19+oUuXLqUS1FHUoVgu295F8M19kJsBZhfo+xb7w25kyIz1xKdlE+rnwaf3tSOqps+/n0tEREqkVDsUd+7cmTVr1hAeHs7cuXP54YcfiIqKYvv27eW+sBFxiCuugxHr/5rJeMUkGgS48vXwTtQLqkZ8Wja3TVvDliOnnZ1URKRKsrvlpqJTy404jM0GK16HprdBQD0AUjJzGfrxerYdS8PT1cK0wa3p1jDIyUFFRCq+Um25+bvs7GzS09OLbCJVhtkM3Z4sLGwAAvZ9xZf9w+jSIJCzeVbu/2QD32097sSQIiJVj93FTVZWFg8//DA1a9akWrVqVK9evcgmUmUdWALfjcBrzq18dGNNbmgWSp7VYNScrXzw20GqWCOpiIjT2F3cPPHEEyxdupT33nsPd3d3PvzwQ55//nnCwsKYNWtWaWQUqRiCGoN/OKQcwu2bwbzV/0ru6RSBYcDLC/fwzPyd5Fttzk4pIlLp2d3npk6dOsyaNYurr74aX19fNm/eTFRUFJ9++ilffvklCxcuLK2sDqE+N1Kq0o7De50gOxVqtca4aSof7XXnpYW7MQzo1jCId+5siY+H67+eSkRE/lKqfW5SUlKoV6+gj4Gvry8pKSkAXHXVVfz222+XEFekEvGrBf1ngrsfHN+EaUYv7m9iYdpdrfFwNfPrvmT6T9NcOCIipcnu4qZevXrExcUB0KhRI+bOnQvADz/8gL+/v0PDiVRI9bvDiLUQ1hKy02DOIHo1qsGcBzsS6O3OnoQz9Ju6ip3H05ydVESkUrK7uBk6dCjbtm0D4Omnn2bq1Kl4eHjw2GOP8cQTTzg8oEiF5BsG/T+BakHQYhBYXGke7s+3D3WiQU1vks7k0H/aGn7aEe/spCIilc5lz3Pzxx9/sGnTJqKiomjWrJmjcpUa9bmRMpWVAl4Bfz1eP50sq5kHf2/CygMnAXj0mihGxTTEbDY5KaSISPlXan1u8vLy6NGjB/v37y/cV7duXW655ZYKUdiIlLm/FzY2G2ycgdfPo/mk6Tbu7RwJwFtLD/Dgp5s4k61VxUVEHMGu4sbV1ZXt27eXVhaRys1kgoa9ALAseY5xXf14vX9z3FzMLNmdyC3vrubwyUwnhxQRqfjs7nNz11138dFHH5VGFpHKzWSCHuMhtDnkZcK7HbjNawtzHuxATR939idlcOM7K/l1X7Kzk4qIVGh297l55JFHmDVrFg0aNKB169ZUq1atyPOTJ092aEBHU58bcbrkvfDlQEg5WPC48ygS249h2Geb2HIkFZMJHotpyMPdo9QPR0TkT/Z8f9td3HTv3v3iJzOZWLp0qT2nK3MqbqRcsOZB7P9g9VsFj+/8ipx6PRj/3e/M3nAUgK4Ng5gyoAUB1dycGFREpHwoleLm0KFDREZGYjJV7H9JqriRcuXX1wADrhoNFhcAvtp4lLHf7SQ7z0aonwdTB7WiVR2t2yYiVVupjJZq0KABycl/9QUYMGAAiYmJl55SRKDbEwUri/9Z2JC8l/4NTMwf0ZnIwGrEp2Vz+7Q1zFgZp4U3RURKqMTFzT//Yl24cCGZmRrZIeIwNhvMewCmtqfRH7P54a5aXN80lHybwf9+3MVDn28mLUvDxUVE/o3do6VEpJRkJkNeNuRmwE9P4D2tLe8EzuP566NwtZj4aWcCvd/8jXWHTjk7qYhIuVbi4sZkMp3X36ai978RKVd8guHBZdBjHNTpBBiY1rzN3Tvu4ZfrzxJRw4sTadkMnL6WSb/sJc9qc3ZiEZFyqcQdis1mM71798bd3R0oWCjzmmuuOW8o+Lx58+wOMXXqVF577TUSEhJo3rw5b7/9Nu3atbvgsdOnT2fWrFns3LkTgNatW/Pyyy9f9Ph/UodiqTD2LITvH4GsgmUazt77K+PWmfhq0zEAWtbx580BLalTw8uZKUVEykSpdCi+++67qVmzJn5+fvj5+XHXXXcRFhZW+PjcZq85c+YwevRoxo8fz+bNm2nevDm9evUiKSnpgscvX76cgQMHsmzZMtasWUN4eDg9e/bk+PHjdr+3SLnWqA88tBbaD4e2D+BZpwWv9W/O2wOuJNgjjy1HUunz1grmbzmuzsYiIn9z2QtnXq727dvTtm1b3nnnHQBsNhvh4eE88sgjPP300//6eqvVSvXq1XnnnXcYMmTIec/n5OSQk5NT+Dg9PZ3w8HC13EjFYhgFMxwDHIjF+Lw/693a8VDaEE7hR5+mIbxwUxNqeLs7N6eISCkptYUzHS03N5dNmzYRExNTuM9sNhMTE8OaNWtKdI6srCzy8vIICAi44PMTJkwo0rIUHh7ukOwiZerv/dvit2EyrLTPWcNK76fpadnMwh0J9HzjNxbtjHdeRhGRcsKpxc3JkyexWq0EBwcX2R8cHExCQkKJzvHUU08RFhZWpED6uzFjxpCWlla4HT169LJzizhVl9EwbCXUjMYzP41pnlP5j/8KMjIzGPbZZh79cgunM3OdnVJExGkq9FDwiRMnMnv2bL799ls8PDwueIy7uzu+vr5FNpEKL6QpPPgrRHTBnH+WR7LfY1ngawSZ0vh+2wmufeM3Fu/SJJsiUjU5tbgJDAzEYrGcN9NxYmIiISEhxb729ddfZ+LEifzyyy80a9asNGOKlE8ubjDoa+j5Enj4E5bxO3P7WKgfVI2TGTk8MGsjj365heQzOf9+LhGRSsSpxY2bmxutW7cmNja2cJ/NZiM2NpaOHTte9HWvvvoqL7zwAosWLaJNmzZlEVWkfHL1gE4Pw/1LoO+bRHYZwIJHu/Bg13qYTfD9thPETP6VuRuOakSViFQZTr8tNXr0aKZPn84nn3zC7t27GT58OJmZmQwdOhSAIUOGMGbMmMLjX3nlFcaOHcuMGTOIiIggISGBhIQEMjIynPURRJwvsAG0vgcAD1cL/73Kj2V9UokO9SXtbB5PfrOdgdPXcihZ/5+ISOXn9OJmwIABvP7664wbN44WLVqwdetWFi1aVNjJ+MiRI8TH/zUC5L333iM3N5fbbruN0NDQwu3111931kcQKV9yMuDz26kb+xA/Vp/E1Pan8XA1s/ZQCte9uYK3Y/eTm6/ZjUWk8nL6PDdlTTMUS6Vns8Kyl2DVm2DLByD1qnGMPNKFX/clA1A/qBrP3XglXRoEOTOpiEiJ2fP9reJGpLI6dRBWvw2bPgbAqNGAw76tGXWkC9syqwNw3ZUhPHN9Y8IDtISDiJRvKm6KoeJGqhTDgBWvw/JXwJYHwJkHNzB5Ux6z1vyB1Wbg7mJm+NX1GdatPh6uFicHFhG5MBU3xVBxI1VS+gk4uh4SdkCPsQDsTTjD5s+f5btTtVlri6Z2dU/G3hBNz+hgTH+fEVlEpBxQcVMMFTcifzqwBD67FYDPLTfxfOYt5OJKp/o1+G+fxjSpZf9CuCIipaXCrC0lIk5UpyO0KlhsdpD1O9b7Pc2dLstZffAkN7y9ksfmbOV46lknhxQRsZ9abkSqut0/woL/QEbBem67qrXj0ZTbOGDUxs3FzL2dI3moe318PVydHFREqjLdliqGihuRC8jJgA0fwrKXwZpDtn8Ud3u8zbrDpwGo7uXKoz0acGf7Ori7qNOxiJQ9FTfFUHEjUozkfRD7PDTshdFyMLG7k3hr4QZ2njSwYaaWvycjezTglla1cLHorraIlB0VN8VQcSNiH+vyV4n/fQVPn+7LyowwwEREDS8eu7YhNzQLw2LWyCoRKX0qboqh4kbEDtY8+LAHxG8D4KT3FbyY1Y+fs67gLB40DPZm9LUN6XVliIaPi0ipUnFTDBU3InY6sq5gKYeDSyG/YPRUrqUaz1mH8kV2JwCiQ315+JoorrsyBLNackSkFKi4KYaKG5FLlHkKVkyC37+FMycAWNDwRZ7a25CMnII1rOoHVeOhq6O4sUUYruqTIyIOpOKmGCpuRC6TzQaxz8GmmfDoVk4b3ny8+jAzV8WRnl1Q5NSu7sn/datP/9a1taSDiDiEiptiqLgRcZCsFPAKKPjZZiN/Rm82uLXniT/aciyzoNWmpo87D3Spxx3twvHRPDkichlU3BRDxY1IKTi8Cmb2AcDmH8H64DuYEFePbeneAHi7u3BH23CGXhVJLX9PZyYVkQpKxU0xVNyIlIL8XNj2Bfz6GqQfA8Bw8WBnvXuZdCya5SkFLTwWs4neTUK4v0s9WoT7OzGwiFQ0Km6KoeJGpBRlp8P6D2D394XDxwFW913O1C05rDpwqnBf24jq3HdVPa6NDtZcOSLyr1TcFEPFjUgZsFlhx1ewfQ6YzHDXNwDsOpHO0XljmZsQzK/5V5KPC7X8PbmzfR1ubxNOkI+7k4OLSHml4qYYKm5Eypg1HywuBT9npcBr9cGwcdbFl5nW6/ghuwW7jAhcLSauaxLK4A51aRtRXZMCikgRKm6KoeJGxIkykmHp/2D7V4UTAhqY+N6zH2NP9yGdagBcEezDXR3rcnPLWni7uzgzsYiUEypuiqHiRqQcsObD1s9g13cFMx//6aOISbx2sBbZeTYAvNws3NAslNvbhNO6rlpzRKoyFTfFUHEjUs7sXQSLngKvQLhvMWk5Vr7deJiclW/zTXo0+4zagIl6gdXo3yacW1vVoqavh7NTi0gZU3FTDBU3IuVUfi64uP3184tBACS5hTMruwvv5fbGigWL2cTVDYPo36Y21zQKxs1FyzyIVAUqboqh4kakArDmw+w74dAysOYCkOIXzTzb1Xx+sj5/GMHYMOPv5UqfpqH0a1GLNnWra9FOkUpMxU0xVNyIVCDZ6bDza/hlLORmFO7+svFU3jgQQtKZHNzJJQdXwvw86dsijH4tatEoxEf9c0QqGRU3xVBxI1IBZSTDlk9h3yKI3w6jd2H1qM7aQ6fI+vkFopMW8LW1C4usbdlt1KFhsA83tahF32Zh1Knh5ez0IuIAKm6KoeJGpIKzWcH8t5XG3+8G8VsLH24yGvJeXl+W2FoDEB3qS+8mIfRuGkJUTZ8yDisijqLiphgqbkQqmbyzsGcBbPgQjq4Do2AY+UlLEJ/kdOft/H6Fh0bV9KZ3kxCuaxJCdKivbl2JVCAqboqh4kakEjv9B6x9DzZMB1s+Zzs+zvcBd/PTzgQ2HzjGo6a5zLRexzEjiLo1vLjuyhCujQ6mZZ3qWt9KpJxTcVMMFTciVcDZVEjaDd41oUZ9ALLWfozXolEA7Dbqss0aySajId9bO+HlVY2rr6jJNY1q0rVhEH6ers7LLiIXpOKmGCpuRKqoI2thyXNwZE2R3WdxY521MQ/nPUIGXljMJq4JN9OxUTjdmtSlXmA13b4SKQdU3BRDxY1IFZeRBIdXQsJ22PE1pB0F4JW2K1iy9zT7kzJ42eVDBliWsd+ozS8evTjVaBAdG4bRKaoGvh5q1RFxBhU3xVBxIyKFDAMSd8LxzdD6bgCOnMrC9bMbCT29ofCwo7YgZlp7ssDoTK3wSK6KCqRrw0Ca1/bHxaIZkkXKgoqbYqi4EZESOZNI9s7vMP36Cu7ZJwt3T8vvy8T8gQD4uZu5OtKT1ldE0KFeDRrU9NYtLJFSYs/3t0sZZRIRqVh8gvHo+CC0vgu2fQHrP4Tk3XS9oiY7LKGsOnCSGmcP83rc06w82ITZtmbsdb+SgLrNaBUVRod6Nbgi2EdLQog4gVpuRERKwjAg/Ti4eEK1GlhtBkfWzSfy53uKHGYzTBw2gtliNOBjl9sJjYymfWQAHerVoHGor4aci1wi3ZYqhoobEXGo5L2wbxG2g8uxHt+Ka05K4VP35D7JclsLAKJNh/F0NeNRuzmtImrQqk51Wtbxx9/LzUnBRSoWFTfFUHEjIqXGMCAzGRK2Y41bxb6ga1meGsy6uFNcf3gi/U2xnDE82WuEs97WiO22ehwK6ELzukG0rludVnWrExXkrVtZIheg4qYYKm5ExBlsi/4LGz7EbM0psj/N8OJHa0eezR+KgRkfDxea1fajSS0/mtXyp1ltP2pX91RHZanyVNwUQ8WNiDiNNQ9OHYBjG+DgMoy9CzHlZ3PYty1Pe7/AtqNpnM3LZ7rrZJINPw4aYRw0QvnD/Qpq165D01p+BVttP2r5q+CRqkXFTTFU3IhIuWHNg0O/grsP1GlPvtVG3K4NNPimZ5HDcgwXNtsast5oxI/WDuw3ahNQzY0rw3xpHOpLoxAfGof6Uj/IGzcXzbsjlZOKm2KouBGRci0rBQ4tL5hc8OR+bAk7MZ8+VPj0dpem3Jr1NHk2E95kkYMbeX/O6uFqMVE/yJvGob40DvWhUUhB8RPk4+6kDyPiOCpuiqHiRkQqFMOApF3wx+qC//Z8kWyTB3sSzuCy/EUaxn3KercOTMuOYWV2JHD+rapAbzca1PShQbA3UTW9iQryJirYmyBvd93akgpDxU0xVNyISKUx5y7Y/UPhQ5ubN5le4SS6hHHAGsKref05fCoLmwEBpGPDRAae5P/Z0uPr4UJUTW8a1PQpKHr+3Gr5e2rElpQ7Km6KoeJGRCoNw4ATm2HDR7B9Dtjy/3rOLxwe28nZXCv7Es8QMfda/M7sw4aJM3iz21abLbYo8rFw1Ajia2s3bBT01/FwNRNRo1rBFliNiBpeRARWIzKwGjV91NojzqHlF0REqgKTCWq1LtiunwypRyDlEJyO49ztKU83C83D/cEdOANmDPw4QwfzbjqYdwNw2qMOmbUHsj85g7iTmcwzPUXKKW8yTnmRt89CLq5sM3z5zajGGYs/GwL6Ehn4t8KnRjXCA7wI9vXQDMxSLqjlRkSkqrDmQ3YqnIkv6MNz6iAYVgi8Ato/CEB+vhXLizUwceGvht+sTRmSN6bw8Sr3R8gwPDlohJGCHzYPf8741MPb2xePGnVxr9OK2tU9qV3di5o+7rrdJZdMLTciInI+iwtUCyzYQppe8BAXswmGLixoAcrPKRiunn8WMpKwnk3nipod+bhGW+JOZhKflESt7afABFdwrOAEeUBKwfZdXCdGrnkYABM2fnAfi8niximPumT61sMaEIVbjQh8a9YhMLgWoX4eVHPX15JcPv0pEhGRv5jNULdTwfYPFiD4z607gK0OdN1YMGT95AEyz5zmbMoxSIkjPzcHT48r6OAawLHTZ7GmnaCJKQ5sQNZeyAIS/jr31PwbeS3/Dnw8XKjva2WEbTaGVw3wrYW7fwg+ATXxjmhNSIAvPu4u6vcjxVJxIyIil8ZsgcAGENgAM+Dz53ZOGHBuOsL87AyS9tQg9VQyuYl7cD21F6+sY/jlxONrS8PTAuTDmex80nNOcK37d5ABJP11vnTDk422K5hPd7b4dKWmjwdhXtDE5QjePn54+9XAJzCMQD9vavp4UMPbDVeLJjWsilTciIhIqXPx8KZmi97UvNCTGcnca9jo7xpAYno2JxOOsGv7g9gyknDLOIFr7mmq5yXiTzrXWLayNq8xR1POcjTlLHmmg7ztPrbI6bINV5INf7YQwLeW3mz27UFNX3dCqpmp6WFQzT+AAC83Aqq5UcPbjYBq7gR4ueHrqRahykLFjYiIOJd3EPBny4+HK1E1r4RmrxU9xmaDE5vJPb6NoX7N6OVZn6T0HMx/ZJK2IxSzNRvP/HRcsOJhyiPclEw4yfyc25a9iWfYm3iG9qbdvOb+AumGJzttkRw2Qlhj1OKU4UcyfuwjAnO1gsKnhpeFGt6uBFRzx9/Lgxo+7gRUc8Pf0w1/L1f8PF3x83LF281FnaTLIY2WEhGRysFmKxgNlpOOLS2ezFPHSK4WxXFLbZLP5OC3/1t67H72oi8fkfsoC2wdAOhvWc5rrh8UPnfWcOO4EUgeFvJw4fX82/nN1hyzCVp5nOBmlzUYrl5Y3DzB3QdXdy883N3wcHfnbHAr3ALC8fd0xc8d/Dxd8anmhZerRYWRHTRaSkREqh6zGbwCwCsAc/UIfCIKWoPqnXu+1SOQcw8k7YbjmyHrVMGSFtlp2M4k8HzPPjzkE01KZi7evx+ArX+d2tOUS5TpROFjdwtgA5sBkbn7GGR8UzBSLOv8WA9tfpSFfxZNvcwbeM91CmdxIw1Xck2upJgCOGmpSZ6LN7/53cRJvyvx8XChtimZK8+sxuoThpubB+6e1TAFNsDNPwwfT1e8PVzw9XDF3cWs22n/oOJGRESqDncfCG9XsP2NGQj8cwMgchj0GlLQGmRYCxY0zToJ1lyw5jG9dluyXf1IP5tHdpwXCXvyyM/JxJpzFnLSMfJzsOXnYbPmEx5chxY2f9LP5hGReQazYVCNHKqRA0CwcRryD0I+zE2/kkV/BABwg3kNI9zevuDHyDUsPJL3KD/b2mI2QS+3HQwzz+OUJYhsizdY3Ml19SHTPZg8d39O+jbF6hNKNTcLvuZsathScPOshpuXLx5urnh6eODl4YaXuzte7q54ubviajFV2KJJxY2IiMg/ubgXbOd4n98V2gPwcLVA8+4F20WM+fsDW1eMrKfIycogMyuTrKwM8lOOYUs9Sn5OJjfV6EY7l3AycvIJSEphZ3w3quUmYdiseOenUsN2Egs23ExWrKaCkWA2AzzzTtPcbS/k74X88zM8lPtX69EN5jW8c5GiCWBc3t3MsvbCYjbRzOUod1t+BrMFi9lMomttMl1rkOkeRJ6rL1leYZg8/PBwLSiaAox0LB5e+NcMp2/zsOKucKkqF8XN1KlTee2110hISKB58+a8/fbbtGvX7qLHf/XVV4wdO5bDhw/ToEEDXnnlFfr06VOGiUVERC6B2YzJuyYe3jXxAGr84+krijxqAAwuekBuFuRmgjWXDzz8yTTcyMq1kp1yBUeONsdIO44tJx1rzllM2am4ZyVgyU2nS/iVhHhEkpWbT4NT+8hM8Mbdlo3LBSqheKMgldVmEJZ/jH7mWLBSsOUVPfb/ckfxs63g+/pG82pGu73DSuuVTAp9rWoXN3PmzGH06NFMmzaN9u3bM2XKFHr16sXevXupWfP8Snn16tUMHDiQCRMmcMMNN/DFF1/Qr18/Nm/eTJMmTZzwCURERMqIm1fBBn/NLeThCr5REBF10ZcNLPKoGfBYwcKr1ryC2242a8HCqzYr77p4kWV1ITvfijUhlOT97uRZrRg5Gbif3g95WXhmHMVszeHm5vVp6teQs3lW6iUeJOewJ9W9fenWMKj0rkEJOH20VPv27Wnbti3vvPMOADabjfDwcB555BGefvrp844fMGAAmZmZ/Pjjj4X7OnToQIsWLZg2bdq/vp9GS4mIiFQ89nx/O3XqxtzcXDZt2kRMTEzhPrPZTExMDGvWrLnga9asWVPkeIBevXpd9PicnBzS09OLbCIiIlJ5ObW4OXnyJFarleDg4CL7g4ODSUhIuOBrEhIS7Dp+woQJ+Pn5FW7h4eGOCS8iIiLlUqVfdGPMmDGkpaUVbkePHnV2JBERESlFTu1QHBgYiMViITExscj+xMREQkJCLviakJAQu453d3fH3d39gs+JiIhI5ePUlhs3Nzdat25NbGxs4T6bzUZsbCwdO3a84Gs6duxY5HiAxYsXX/R4ERERqVqcPhR89OjR3H333bRp04Z27doxZcoUMjMzGTp0KABDhgyhVq1aTJgwAYCRI0fSrVs3Jk2axPXXX8/s2bPZuHEjH3zwQXFvIyIiIlWE04ubAQMGkJyczLhx40hISKBFixYsWrSosNPwkSNHMJv/amDq1KkTX3zxBc8++yz//e9/adCgAfPnz9ccNyIiIgKUg3luyprmuREREal4Ksw8NyIiIiKOpuJGREREKhUVNyIiIlKpqLgRERGRSkXFjYiIiFQqKm5ERESkUnH6PDdl7dzId60OLiIiUnGc+94uyQw2Va64OXPmDIBWBxcREamAzpw5g5+fX7HHVLlJ/Gw2GydOnMDHxweTyeTQc6enpxMeHs7Ro0c1QWAp0nUuG7rOZUPXuezoWpeN0rrOhmFw5swZwsLCiqxccCFVruXGbDZTu3btUn0PX19f/Y9TBnSdy4auc9nQdS47utZlozSu87+12JyjDsUiIiJSqai4ERERkUpFxY0Dubu7M378eNzd3Z0dpVLTdS4bus5lQ9e57Ohal43ycJ2rXIdiERERqdzUciMiIiKVioobERERqVRU3IiIiEilouJGREREKhUVNw4ydepUIiIi8PDwoH379qxfv97ZkSqUCRMm0LZtW3x8fKhZsyb9+vVj7969RY7Jzs5mxIgR1KhRA29vb2699VYSExOLHHPkyBGuv/56vLy8qFmzJk888QT5+fll+VEqlIkTJ2IymRg1alThPl1nxzh+/Dh33XUXNWrUwNPTk6ZNm7Jx48bC5w3DYNy4cYSGhuLp6UlMTAz79+8vco6UlBQGDRqEr68v/v7+3HfffWRkZJT1Rym3rFYrY8eOJTIyEk9PT+rXr88LL7xQZO0hXedL89tvv9G3b1/CwsIwmUzMnz+/yPOOuq7bt2+nS5cueHh4EB4ezquvvuqYD2DIZZs9e7bh5uZmzJgxw/j999+NBx54wPD39zcSExOdHa3C6NWrl/Hxxx8bO3fuNLZu3Wr06dPHqFOnjpGRkVF4zLBhw4zw8HAjNjbW2Lhxo9GhQwejU6dOhc/n5+cbTZo0MWJiYowtW7YYCxcuNAIDA40xY8Y44yOVe+vXrzciIiKMZs2aGSNHjizcr+t8+VJSUoy6desa99xzj7Fu3Trj0KFDxs8//2wcOHCg8JiJEycafn5+xvz5841t27YZN954oxEZGWmcPXu28JjrrrvOaN68ubF27VpjxYoVRlRUlDFw4EBnfKRy6aWXXjJq1Khh/Pjjj0ZcXJzx1VdfGd7e3sabb75ZeIyu86VZuHCh8cwzzxjz5s0zAOPbb78t8rwjrmtaWpoRHBxsDBo0yNi5c6fx5ZdfGp6ensb7779/2flV3DhAu3btjBEjRhQ+tlqtRlhYmDFhwgQnpqrYkpKSDMD49ddfDcMwjNTUVMPV1dX46quvCo/ZvXu3ARhr1qwxDKPgf0az2WwkJCQUHvPee+8Zvr6+Rk5OTtl+gHLuzJkzRoMGDYzFixcb3bp1KyxudJ0d46mnnjKuuuqqiz5vs9mMkJAQ47XXXivcl5qaari7uxtffvmlYRiGsWvXLgMwNmzYUHjMTz/9ZJhMJuP48eOlF74Cuf7664177723yL5bbrnFGDRokGEYus6O8s/ixlHX9d133zWqV69e5O+Np556yrjiiisuO7NuS12m3NxcNm3aRExMTOE+s9lMTEwMa9ascWKyii0tLQ2AgIAAADZt2kReXl6R69yoUSPq1KlTeJ3XrFlD06ZNCQ4OLjymV69epKen8/vvv5dh+vJvxIgRXH/99UWuJ+g6O8r3339PmzZt6N+/PzVr1qRly5ZMnz698Pm4uDgSEhKKXGc/Pz/at29f5Dr7+/vTpk2bwmNiYmIwm82sW7eu7D5MOdapUydiY2PZt28fANu2bWPlypX07t0b0HUuLY66rmvWrKFr1664ubkVHtOrVy/27t3L6dOnLytjlVs409FOnjyJ1Wot8hc9QHBwMHv27HFSqorNZrMxatQoOnfuTJMmTQBISEjAzc0Nf3//IscGBweTkJBQeMyFfg/nnpMCs2fPZvPmzWzYsOG853SdHePQoUO89957jB49mv/+979s2LCBRx99FDc3N+6+++7C63Sh6/j361yzZs0iz7u4uBAQEKDr/Kenn36a9PR0GjVqhMViwWq18tJLLzFo0CAAXedS4qjrmpCQQGRk5HnnOPdc9erVLzmjihspd0aMGMHOnTtZuXKls6NUOkePHmXkyJEsXrwYDw8PZ8eptGw2G23atOHll18GoGXLluzcuZNp06Zx9913Ozld5TF37lw+//xzvvjiC6688kq2bt3KqFGjCAsL03Wu4nRb6jIFBgZisVjOG02SmJhISEiIk1JVXA8//DA//vgjy5Yto3bt2oX7Q0JCyM3NJTU1tcjxf7/OISEhF/w9nHtOCm47JSUl0apVK1xcXHBxceHXX3/lrbfewsXFheDgYF1nBwgNDSU6OrrIvsaNG3PkyBHgr+tU3N8bISEhJCUlFXk+Pz+flJQUXec/PfHEEzz99NPccccdNG3alMGDB/PYY48xYcIEQNe5tDjqupbm3yUqbi6Tm5sbrVu3JjY2tnCfzWYjNjaWjh07OjFZxWIYBg8//DDffvstS5cuPa+psnXr1ri6uha5znv37uXIkSOF17ljx47s2LGjyP9QixcvxtfX97wvmqqqR48e7Nixg61btxZubdq0YdCgQYU/6zpfvs6dO583lcG+ffuoW7cuAJGRkYSEhBS5zunp6axbt67IdU5NTWXTpk2FxyxduhSbzUb79u3L4FOUf1lZWZjNRb/GLBYLNpsN0HUuLY66rh07duS3334jLy+v8JjFixdzxRVXXNYtKUBDwR1h9uzZhru7uzFz5kxj165dxoMPPmj4+/sXGU0ixRs+fLjh5+dnLF++3IiPjy/csrKyCo8ZNmyYUadOHWPp0qXGxo0bjY4dOxodO3YsfP7cEOWePXsaW7duNRYtWmQEBQVpiPK/+PtoKcPQdXaE9evXGy4uLsZLL71k7N+/3/j8888NLy8v47PPPis8ZuLEiYa/v7/x3XffGdu3bzduuummCw6lbdmypbFu3Tpj5cqVRoMGDar8EOW/u/vuu41atWoVDgWfN2+eERgYaDz55JOFx+g6X5ozZ84YW7ZsMbZs2WIAxuTJk40tW7YYf/zxh2EYjrmuqampRnBwsDF48GBj586dxuzZsw0vLy8NBS9P3n77baNOnTqGm5ub0a5dO2Pt2rXOjlShABfcPv7448Jjzp49azz00ENG9erVDS8vL+Pmm2824uPji5zn8OHDRu/evQ1PT08jMDDQ+M9//mPk5eWV8aepWP5Z3Og6O8YPP/xgNGnSxHB3dzcaNWpkfPDBB0Wet9lsxtixY43g4GDD3d3d6NGjh7F3794ix5w6dcoYOHCg4e3tbfj6+hpDhw41zpw5U5Yfo1xLT083Ro4cadSpU8fw8PAw6tWrZzzzzDNFhhbrOl+aZcuWXfDv5LvvvtswDMdd123bthlXXXWV4e7ubtSqVcuYOHGiQ/KbDONvUzmKiIiIVHDqcyMiIiKVioobERERqVRU3IiIiEilouJGREREKhUVNyIiIlKpqLgRERGRSkXFjYiIiFQqKm5ERESkUlFxIyLlyj333EO/fv2cHUNEKjAXZwcQkarDZDIV+/z48eN588030cTpInI5VNyISJmJj48v/HnOnDmMGzeuyOrZ3t7eeHt7OyOaiFQiui0lImUmJCSkcPPz88NkMhXZ5+3tfd5tqauvvppHHnmEUaNGUb16dYKDg5k+fTqZmZkMHToUHx8foqKi+Omnn4q8186dO+nduzfe3t4EBwczePBgTp48WeKshmHw3HPPUadOHdzd3QkLC+PRRx911KUQkVKk4kZEyr1PPvmEwMBA1q9fzyOPPMLw4cPp378/nTp1YvPmzfTs2ZPBgweTlZUFQGpqKtdccw0tW7Zk48aNLFq0iMTERG6//fYSv+c333zDG2+8wfvvv8/+/fuZP38+TZs2La2PKCIOpFXBRcQpZs6cyahRo0hNTS2y/5577iE1NZX58+cDBS03VquVFStWAGC1WvHz8+OWW25h1qxZACQkJBAaGsqaNWvo0KEDL774IitWrODnn38uPO+xY8cIDw9n7969NGzY8F/zTZ48mffff5+dO3fi6urqmA8tImVCLTciUu41a9as8GeLxUKNGjWKtKIEBwcDkJSUBMC2bdtYtmxZYR8eb29vGjVqBMDBgwdL9J79+/fn7Nmz1KtXjwceeIBvv/2W/Px8R30kESlF6lAsIuXeP1tOTCZTkX3nRmHZbDYAMjIy6Nu3L6+88sp55woNDS3Re55r5VmyZAmLFy/moYce4rXXXuPXX39VS45IOafiRkQqnVatWvHNN98QERGBi8ul/zXn6elJ37596du3LyNGjKBRo0bs2LGDVq1aOTCtiDiabkuJSKUzYsQIUlJSGDhwIBs2bODgwYP8/PPPDB06FKvVWqJzzJw5k48++oidO3dy6NAhPvvsMzw9Palbt24ppxeRy6XiRkQqnbCwMFatWoXVaqVnz540bdqUUaNG4e/vj9lc8NfezJkzi51U0N/fn+nTp9O5c2eaNWvGkiVL+OGHH6hRo0ZZfQwRuUQaLSUiVdL48eP59ddfWb58ubOjiIiDqc+NiFRJP/30E++8846zY4hIKVDLjYiIiFQq6nMjIiIilYqKGxEREalUVNyIiIhIpaLiRkRERCoVFTciIiJSqai4ERERkUpFxY2IiIhUKipuREREpFJRcSMiIiKVyv8DeEiQiGNnXfgAAAAASUVORK5CYII=\n", + "text/plain": [ + "
" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "seconds = np.arange(1001)\n", + "fraction = np.exp(-4.59e-3 * seconds)\n", + "\n", + "plt.plot(seconds, fraction, linestyle='-', label='Theoretical')\n", + "plt.plot(seconds, decay(4.59e-3, n=5000), linestyle='--', label='Simulation')\n", + "plt.xlabel('Time, s')\n", + "plt.ylabel('Fraction of Undecayed Nuclei')\n", + "plt.legend();" + ] + }, + { + "cell_type": "markdown", + "id": "7e193e9f-25ca-47fe-9809-16db9141e164", + "metadata": {}, + "source": [ + "## Further reading\n", + "\n", + "- To learn more about NumPy random number generators, see [Random Generator](https://numpy.org/doc/stable/reference/random/index.html)\n", + "- To learn more about random number distributions and their generation, see [Monte Carlo Simulation Part 5: Randomness & Random Number Generation](https://towardsdatascience.com/monte-carlo-simulation-421110b678c5?gi=f86d4ec9f452)\n", + "- For examples of using random variables for chemical simulations, see chapter 9 of [Scientific Computing for Chemists](https://weisscharlesj.github.io/SciCompforChemists/intro.html)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.10.8" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +}