Skip to content

Instantly share code, notes, and snippets.

@santteegt
Created July 17, 2020 00:46
Show Gist options
  • Save santteegt/3e9462eff34179a57aaf4e5b3a52f8f7 to your computer and use it in GitHub Desktop.
Save santteegt/3e9462eff34179a57aaf4e5b3a52f8f7 to your computer and use it in GitHub Desktop.
robot-marbles-part-6.ipynb
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# cadCAD Tutorials: The Robot and the Marbles, part 6\n",
"In parts [1](../robot-marbles-part-1/robot-marbles-part-1.ipynb) and [2](../robot-marbles-part-2/robot-marbles-part-2.ipynb) we introduced the 'language' in which a system must be described in order for it to be interpretable by cadCAD and some of the basic concepts of the library:\n",
"* State Variables\n",
"* Timestep\n",
"* State Update Functions\n",
"* Partial State Update Blocks\n",
"* Simulation Configuration Parameters\n",
"* Policies\n",
"\n",
"In [part 3](../robot-marbles-part-3/robot-marbles-part-3.ipynb) we covered how to describe the presence of asynchronous subsystems within the system being modeled in cadCAD. And [part 4](../robot-marbles-part-4/robot-marbles-part-4.ipynb) introduced Monte Carlo simulations.\n",
"\n",
"In this notebook, we'll cover cadCAD's support for A/B testing, a useful feature when analyzing different design options for a system. Let's start by copying the base configuration with which we ended Part 4. Here's the description of that system:\n",
"\n",
"__The robot and the marbles__ \n",
"* Picture a box (`box_A`) with ten marbles in it; an empty box (`box_B`) next to the first one; and __two__ robot arms capable of taking a marble from any one of the boxes and dropping it into the other one. \n",
"* The robots are programmed to take one marble at a time from the box containing the largest number of marbles and drop it in the other box. They repeat that process until the boxes contain an equal number of marbles.\n",
"* The robots act __asynchronously__ and __non-deterministically__; at every timestep each robot acts with a probability P: 50% for robot 1 and 33.33% for robot 2."
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"%%capture\n",
"from numpy.random import rand\n",
"\n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
"# List of all the state variables in the system and their initial values\n",
"initial_conditions = {\n",
" 'box_A': 10, # as per the description of the example, box_A starts out with 10 marbles in it\n",
" 'box_B': 0 # as per the description of the example, box_B starts out empty\n",
"}\n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
"\n",
"\n",
"\n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
"# Settings of general simulation parameters, unrelated to the system itself\n",
"# `T` is a range with the number of discrete units of time the simulation will run for;\n",
"# `N` is the number of times the simulation will be run (Monte Carlo runs)\n",
"# In this example, we'll run the simulation once (N=1) and its duration will be of 10 timesteps\n",
"# We'll cover the `M` key in a future article. For now, let's leave it empty\n",
"simulation_parameters = {\n",
" 'T': range(10),\n",
" 'N': 50, # We'll run the same simulation 50 times; the random events in each simulation are independent\n",
" 'M': {}\n",
"}\n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
"\n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
"# We specify the robot arm's logic in a Policy Function\n",
"def robot_arm(params, step, sL, s):\n",
" add_to_A = 0\n",
" if (s['box_A'] > s['box_B']):\n",
" add_to_A = -1\n",
" elif (s['box_A'] < s['box_B']):\n",
" add_to_A = 1\n",
" return({'add_to_A': add_to_A, 'add_to_B': -add_to_A})\n",
" \n",
"robots_probabilities = [0.5,1/3] # Robot 1 acts with a 50% probability; Robot 2, 33.33%\n",
"\n",
"def robot_arm_1(params, step, sL, s):\n",
" _robotId = 1\n",
" if rand()<robots_probabilities[_robotId-1]: # draw a random number between 0 and 1; if it's smaller than the robot's parameter, it acts\n",
" return robot_arm(params, step, sL, s)\n",
" else:\n",
" return({'add_to_A': 0, 'add_to_B': 0}) # otherwise, the robot doesn't interfere with the system\n",
"\n",
"def robot_arm_2(params, step, sL, s):\n",
" _robotId = 2\n",
" if rand()<robots_probabilities[_robotId-1]: # draw a random number between 0 and 1; if it's smaller than the robot's parameter, it acts\n",
" return robot_arm(params, step, sL, s)\n",
" else:\n",
" return({'add_to_A': 0, 'add_to_B': 0}) # otherwise, the robot doesn't interfere with the system\n",
"\n",
"\n",
"\n",
" \n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
"# We make the state update functions less \"intelligent\",\n",
"# ie. they simply add the number of marbles specified in _input \n",
"# (which, per the policy function definition, may be negative)\n",
"def increment_A(params, step, sL, s, _input):\n",
" y = 'box_A'\n",
" x = s['box_A'] + _input['add_to_A']\n",
" return (y, x)\n",
"\n",
"def increment_B(params, step, sL, s, _input):\n",
" y = 'box_B'\n",
" x = s['box_B'] + _input['add_to_B']\n",
" return (y, x)\n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
"\n",
"\n",
"\n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
"# In the Partial State Update Blocks, \n",
"# the user specifies if state update functions will be run in series or in parallel\n",
"# and the policy functions that will be evaluated in that block\n",
"partial_state_update_blocks = [\n",
" { \n",
" 'policies': { # The following policy functions will be evaluated and their returns will be passed to the state update functions\n",
" 'robot_arm_1': robot_arm_1,\n",
" 'robot_arm_2': robot_arm_2\n",
" },\n",
" 'variables': { # The following state variables will be updated simultaneously\n",
" 'box_A': increment_A,\n",
" 'box_B': increment_B\n",
" }\n",
" }\n",
"]\n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # "
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"%%capture\n",
"from cadCAD import configs\n",
"from cadCAD.configuration import append_configs\n",
"\n",
"configs.clear()\n",
"\n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
"# The configurations above are then packaged into a `configs` object\n",
"append_configs(initial_state=initial_conditions, #dict containing variable names and initial values\n",
" partial_state_update_blocks=partial_state_update_blocks, #dict containing state update functions\n",
" sim_configs=simulation_parameters #dict containing simulation parameters\n",
" )\n",
"\n",
"from cadCAD.engine import ExecutionMode, ExecutionContext, Executor\n",
"exec_mode = ExecutionMode()\n",
"exec_context = ExecutionContext(exec_mode.multi_mode)\n",
"executor = Executor(exec_context, configs) # Pass the configuration object inside an array\n",
"raw_result, tensor, sessions = executor.execute() # The `execute()` method returns a tuple; its first elements contains the raw results\n",
"\n",
"%matplotlib inline\n",
"import pandas as pd\n",
"df = pd.DataFrame(raw_result)"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<div>\n",
"<style scoped>\n",
" .dataframe tbody tr th:only-of-type {\n",
" vertical-align: middle;\n",
" }\n",
"\n",
" .dataframe tbody tr th {\n",
" vertical-align: top;\n",
" }\n",
"\n",
" .dataframe thead th {\n",
" text-align: right;\n",
" }\n",
"</style>\n",
"<table border=\"1\" class=\"dataframe\">\n",
" <thead>\n",
" <tr style=\"text-align: right;\">\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" <th>box_A</th>\n",
" <th>box_B</th>\n",
" <th>simulation</th>\n",
" </tr>\n",
" <tr>\n",
" <th>run</th>\n",
" <th>timestep</th>\n",
" <th>substep</th>\n",
" <th></th>\n",
" <th></th>\n",
" <th></th>\n",
" </tr>\n",
" </thead>\n",
" <tbody>\n",
" <tr>\n",
" <th rowspan=\"5\" valign=\"top\">1</th>\n",
" <th>0</th>\n",
" <th>0</th>\n",
" <td>10</td>\n",
" <td>0</td>\n",
" <td>0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>1</th>\n",
" <th>1</th>\n",
" <td>10</td>\n",
" <td>0</td>\n",
" <td>0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>2</th>\n",
" <th>1</th>\n",
" <td>9</td>\n",
" <td>1</td>\n",
" <td>0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>3</th>\n",
" <th>1</th>\n",
" <td>8</td>\n",
" <td>2</td>\n",
" <td>0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>4</th>\n",
" <th>1</th>\n",
" <td>6</td>\n",
" <td>4</td>\n",
" <td>0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>...</th>\n",
" <th>...</th>\n",
" <th>...</th>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" <td>...</td>\n",
" </tr>\n",
" <tr>\n",
" <th rowspan=\"5\" valign=\"top\">50</th>\n",
" <th>6</th>\n",
" <th>1</th>\n",
" <td>5</td>\n",
" <td>5</td>\n",
" <td>0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>7</th>\n",
" <th>1</th>\n",
" <td>5</td>\n",
" <td>5</td>\n",
" <td>0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>8</th>\n",
" <th>1</th>\n",
" <td>5</td>\n",
" <td>5</td>\n",
" <td>0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>9</th>\n",
" <th>1</th>\n",
" <td>5</td>\n",
" <td>5</td>\n",
" <td>0</td>\n",
" </tr>\n",
" <tr>\n",
" <th>10</th>\n",
" <th>1</th>\n",
" <td>5</td>\n",
" <td>5</td>\n",
" <td>0</td>\n",
" </tr>\n",
" </tbody>\n",
"</table>\n",
"<p>550 rows × 3 columns</p>\n",
"</div>"
],
"text/plain": [
" box_A box_B simulation\n",
"run timestep substep \n",
"1 0 0 10 0 0\n",
" 1 1 10 0 0\n",
" 2 1 9 1 0\n",
" 3 1 8 2 0\n",
" 4 1 6 4 0\n",
"... ... ... ...\n",
"50 6 1 5 5 0\n",
" 7 1 5 5 0\n",
" 8 1 5 5 0\n",
" 9 1 5 5 0\n",
" 10 1 5 5 0\n",
"\n",
"[550 rows x 3 columns]"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"from IPython.display import display\n",
"tmp_rows = pd.options.display.max_rows\n",
"pd.options.display.max_rows = 10\n",
"display(df.set_index(['run', 'timestep', 'substep']))\n",
"pd.options.display.max_rows = tmp_rows"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXAAAAEGCAYAAAB8Ys7jAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nOzdd1gU19fA8e/QQRQVLCgo9gYCiiIICjbsbdXYNRp7bNHYey+xa0zsNWoE7L2AgmAFFHsv2MUGSOe+f+jPNyYWyi4rej/Pw6PsztxzRuTs7OydexQhBJIkSVLWo6PtBCRJkqT0kQVckiQpi5IFXJIkKYuSBVySJCmLkgVckiQpi9LLzGAWFhbCxsYmXfvGxMSQLVs29Sb0FcfVZmx5zN9HbHnMWSfumTNnngkh8vznCSFEpn1VrFhRpJefn1+6980IbcXVZmx5zN9HbHnMWScucFp8pKbKSyiSJElZlCzgkiRJWZQs4JIkSVlUpn6IKUmS9CmJiYlEREQQFxensRhmZmZcunRJY+NnNK6RkRFWVlbo6+unalxZwCVJ+ipERESQPXt2bGxsUBRFIzGioqLInj27RsbOaFwhBJGRkURERFCkSJFUjfvFSyiKoqxQFOWJoijn//FYbkVRDiiKcu3dn7lSFU2SJOkT4uLiMDc311jx/topioK5uXma3oGk5hr4KqDuvx4bBhwSQpQADr37XpIkKUO+1+L9P2k9/i9eQhFCHFUUxeZfDzcBPN79fTXgDwxNU+Q08Kvblccnz3KyZX2M8plrKsxHRT16SJxteYwscmdqXEmSpC9RRCrWA39XwHcKIWzfff9SCJHz3d8V4MX/vv/Ivt2B7gD58uWruHHjxjQn+cCzx6cSS/NYaSYEOjmzY9b3B4w8nTL1DCE6OhpTU9NMi6ftuNqMLY9Z+7HNzMwoXry4RuMmJyejq6ur0RgZjXv9+nVevXr1wWOenp5nhBBO/9n4Y3f3/PsLsAHO/+P7l/96/kVqxknvnZgpKSli/9I14uzouWJn2fpiPSXFekqKvVVaiYszl4mom3fTNW5q7F++TuyppBLrKSn8G/cUMfcfaSzWv2XVu8ayYmx5zNqPffHiRY3Hff369Wefv3XrlihXrpzG4s6ZM0cYGhqKly9ffnLbj/07oOY7MR8rimIJ8O7PJ+kcJ1UURUG/uDXlJ/SnwYVdNLi0m/KTBpASn0DorzPYXrQWeyo04/zkxby6fEOtsfWLFqRO8CYcfxvKowNB7CrbgOvLNv/vhUuSJCnVNmzYQKVKlfD19VXLeOmdRrgd6ARMe/fnNrVkk0pmpYthNrIXtiN7EX3zHnd99nHPZz/nRs3l3Ki5mJUrgbWqDtaqOuS0K5Xhyx46urqUGdQFqyY1OdFtFCe7jeLOhp04L52EaVFrNR2VJEn/c2bAZF6EXVbrmLkcSlNyYr8vbpeUlES7du0ICQmhXLlyrFmzhuDgYAYPHkxSUhKVKlVi8eLFnDt3jq5du3Ly5EmSk5OpXLkymzZtwtbW9qPj3rhxg+joaH7//XcmT57Mjz/+mOFjSs00wg1AMFBKUZQIRVG68rZw11YU5RpQ6933WmFa1Jqyv/6E1/G/aXLXn4rzRmJonpPzE39nj30TdpT0ImzYb0SeOpfhs+bsxQtT89BqKv85geenz7PLtiGX56wiJTlZTUcjSZK2Xblyhd69e3Pp0iVy5MjB7Nmz6dy5M5s2bSI8PJykpCQWL15MpUqVaNy4MaNGjWLIkCG0b9/+k8UbYOPGjbRu3Rp3d3euXLnC48ePM57sx66raOorM1cjfPPoqbj6xwZxqPaP4i/dMmI9JcXWwp7i9MAp4kngaZGSnJyhuDH3Hgq/hj3eXouv3EK8CL+SpvwyElvTvrZro99yXG3G/tqO+Wu5Bm5tbf3++0OHDgkPDw/h7u7+/rGDBw+KZs2aCSGEiI+PF+XLlxeVK1cWSUlJn41brlw5cfXqVSGEEAMHDhQLFiz46LaZcQ38q2ecz4ISPVpTY/8Kmj8+hvOKKZjZluDaovUccGvLVqtqnOoznkeHg0lJSkrz+CZW+am+fTGuf80i+uY99lZoTvj4hSQnJGjgaCRJyiz/vuSaM+dHJ9gBEBkZSXR0NFFRUZ+9AefChQtcu3aN2rVrY2Njw8aNG9mwYUOGc/1mC/g/GZrnotiPKjx2/knzJ8G4rv8NCxdHbq705XDNzmyxdOPETyN5sPdomgqwoijYtGlIg4u7sW7pRfi4BeytqCLy1DnNHYwkSRp19+5dgoODAfjrr79wcnLi9u3bXL9+HYC1a9dSvXp1AHr06MHEiRNp164dQ4d++lYYb29vxo0bx+3bt7l9+zYPHjzgwYMH3LlzJ0O5fhcF/J8MzLJj07YR7j4LUD0Nxs17PvlruXJn027863XDN68rQR2HELHtIEmxqbul1ShPbqqun0X1HX+Q8OIV+6v8QMjg6SS9idXw0UiSpG6lSpVi0aJFlClThhcvXjBw4EBWrlxJy5YtsbOzQ0dHh549e7JmzRr09fVp27Ytw4YN49SpUxw+fPijY/r4+NCsWbMPHmvWrBnpuS/mn77rxaz0splQSOVFIZUXyXHxPDxwjHs++4nYdojba7ehl80EvUpluPPkDQXqV0ff9PMtkQo29KTBBSfChs7k8qwVRGw9iPOySeTzcM6kI5IkKSNsbGy4fPm/s19q1qxJaGjoB4917NiRjh07AqCrq8uJEyc+Oe65c+f+s5jV7NmzM5zvd3cG/im6RoZYNaqBy6ppqJ4E4blvOTbtGpEQdpVjPwzEN48LR5v14da6bSS8ivrkOAZm2an8xwRq+q0B4JBnR072GPPZfSRJktLjuz4D/xQdfX0s67hhWceNmFYelNXNxj2f/dzz3U/E1oPo6OuTr5YLhVR1KNik5kfXScnn4Uz9c9sJHzufy7NXcX+XP5UWj8OqUQ0tHJEkSZkhPDycDh06fPCYoaEhBw8e1Eg8WcC/QNHVIZ+HM/k8nKk4byTPTpx9W8x99nHip1EoPcaSt3qltzcONauNsWXe9/vqmRjjOHMohVrV40TXkRxt3IvCbRpScd5IjPLIxbEk6VtjZ2dHWFjYfx6PitLMO3B5CSUNFB0d8rg4UuG3oTS+eYi6p30oM+Qn3kQ84nSfCWwpWI0D7m25PHcVMXcfvN/PvFJ5vE77YDe+L/e897GrTD1u/7VD3o4vSVKGyAKeToqikLuiLQ5TfqHh5b3UD9+B3difSXwVRcjAqWwr7Mneyi24OH0JUdfvoGtggN2Yn6kbugXT4oUJajeYI4178SbikbYPRZKkLEoWcDVQFIWctiWxG/sz9c/toOGVvdhPHQRCEDZsFjtK1GG3fWPCJy5CURRqH9tAhdnDeXwomJ1l63Ptz42IlBRtH4YkSVmMLOAakKNkEcoN607dUz40vnUIx1nD0DM1IXzMfHaVa8Bu24bER77AZe0MzCvZcarnWA7V7ETU9YxN6pck6fuSoQKuKEp/RVHOK4pyQVGUAepK6ltiamNFmV9+pM6xjTS9fxSnhWMwtszDxalLCGzRj+hbERjkMuOJ/0l22Tbk0m/L03VrvyRJGXf79u3PLkiVXnfu3MHY2BgHBwfs7e1xdXXlypUrGR433QVcURRboBtQGbAHGiqKotl2GlmcSYF8lOzTjpqH19Ds0TEqL5lIjpI2JEbFALxf33yjfjlehF3ScraSJKlTsWLFCAsL4+zZs3Tq1IkpU6ZkeMyMTCMsA5wQQrwBUBTlCNAcmJHhrL4DRnlyU7xbK4p3a0XCi1dE7PDjns8+7m9/eyvuHsemANzf+SeWdaqio6+vzXQlKVMN+HsOYRFX1Tqmg1VJJtb76YvbaWo98H96/fo1uXLlyvAxpaon5kd3VJQyvG3k4ALE8rY7/WkhRN9/bZfhnpjw/fQNTHkTR+z+47ya9+FKZcZ1XTGuVgHDiqVRDDRbzL+2Xonfclxtxv7ajvmfPTGHbv+d8Afq7a5lV6AYUxr0+Gxvyjt37mBnZ8f+/fupUqUKvXv3xsbGhpUrV7J9+3ZKlChB9+7dsbe3p0+fPkyYMIH4+HhiY2MpWLAggwYN+ui4t27dwtnZmRIlShAVFUVsbCyHDx/G2vq/DWHU3hPzU19AV+AMcBRYDMz93PaZuR64umhzzeTd42a97//5v6+/c1QQgW1/EXd99onEmDcaifu1rRP9LcfVZuyv7Zi/5fXAw8PDP+i1uXHjRuHl5fXRbTNtPXAhxHIhREUhRDXgBaDe9zzfOePqFWj56gwlerd9/5hBbjMe7g0kQNUXnzwuBLTsx+2Nu0iMitZippL07dDEeuD/1rhxY44ePZruHP8no7NQ8r77sxBvr3//leGMpA/o5zCl0qKx1DqyjuwlbIi5fR+rxjVw855PkY5NeBpwhqA2v+CTx4UjjXtyc/UWEl68+vLAkiR9lCbWA/+3wMBAihUrluFcM7oWio+iKOZAItBHCPEywxlJH5W3WiXqnd3G+fELufTbCh7uC6DS4nE4LRzDs6BQ7vns457vAe7v8EPR0yNfjSoUUtXBqmktjPKaazt9Scoy/rceeJcuXShbtizz58+nSpUqtGzZ8v2HmP9eDzw5ORlXV1cOHz5MjRofX7Duxo0bODg4IITAwMCAZcuWZTjXDBVwIYR7hjOQUk3P2AiHaYMp1LIux7uO5GjTPhRqVY+K80dRce5IKswZQeSp8LfF3HsfJ3uM4VSvceSp5oS1ygvrZrUxKZhP24chSV8tTa0HXrhwYWJj1d/gRd6JmQXlrmhL3VPelJ80gIitB9lVtgG31m0DwKJyeRyn/0qj6weoF7qVsiN6EPc4kjN9J7LVqhr7XVtzadYKom9HaPkoJEnKKFnAsygdfX1sR/aiXtg2cpQqQnCHIfg36P5+FURFUcjlUAb7iQNoeHE3DS7upvzE/iTHxhE6eDrbi9Rkr1NzLkz9k9dXb2n5aCTp2xAeHo6Dg8MHX87OmuvIJdcDz+LMyhSjVsB6rv3+F2eHz2ZXuQY4TB9MiZ5tUHR0PtjObFRvbEf1JurG3beXWXz2c3bEbM6OmI2ZbUmsVXUopKojl7mVtEYI8Z9ZIFlJRtcDT+vvnjwD/wbo6OpSqm8H6p/fgYWLI6f7TOCgR4dPnllnL1aIskO64XViM03u+lNh7ggMcuXg/IRF7C7fmCcdxxA2YjbPz5yXxVzKNEZGRkRGRn63/+eEEERGRmJkZJTqfeQZ+DfE1MYKz33LubV6C2cGTmV3+caUH9+X0oO6oKP38R91NmtLSvfvROn+nYh99JSIrQc5u2wTl2Ys4+LUP8lmUxDr5nWwVtXBoorDB2f1kqROVlZWRERE8PTpU43FiIuLS1OBzOy4RkZGWFlZpXpcWcC/MYqiULRzcyy93DjVZwJhw2Zx5++9VFkxhVz2pT+7r3H+PJTo2Yb7pS1xsbMnYtsh7vns5+qCdVyevRLjAnmxalabQqo65HF3+uSLgiSlh76+PkWKFNFoDH9/fxwdHTUaIzPjytOpb5SxZV6q+S7EzXs+sfcfs9dJxdlRc0iOi0/V/obmuSjWpQUeu5bQ/GkwLutmYu5sz83l3hyq0YktBdw50X00D/YFkJKYqOGjkSTpY+Qp1DeukMqLfJ7OhPwyjQuT/+Cez36cl08mj2uFVI9hYJadIu0aU6RdY5Ji3vBgz1Hueu/jzoad3Fj6N/o5c2DVuAbWqjpY1nFD18hQg0ckSdL/yDPw74Bh7py4rJqGx95lJL2J44BbW073m0RidEyax9LLZkKhFnVx2zgH1dPjVNv2O1aNaxCx/TBHm/TGJ08VAlsP5K73XpJi3mjgaCRJ+h95Bv4dKeDlToPzOzg7Yg5XF67j/vbDVF4yAcs6bukaT9fIEKvGNbFqXJPkhAQe+53gns9+IrYe5O6m3egaG2FZ1x1rVR0KNvTEwCy7mo9Ikr5v8gz8O6Of3RSnBaOpdXQdukYG+Hl15fiPw4l/nrFlbHQNDCjg5Y7zkok0exBAzcOrKdpFReTxMILb/4pvXhf8G3Tnxgpv4iNfqOloJOn7ltHVCAe+64d5XlGUDYqiZP78HCld8ro5US9sG+VG9OTW2m3sKtuAuz771DK2jp4e+TyrUGnhGJpGHKX2sQ2U/Lk9ry5c50TXkfjmq8rh2j9y7Y8NxD7S3JQxSfrWZaQnZkGgH+AkhLAFdIHW6kpM0jxdI0PsJw+k7mkfjAvkJbBFPwJa9CP5ufqWo1V0dMjjWoEKs4bR+NYh6p72ocyQn4i584BTvcaxpYA7B6q14/K81cTce6i2uJL0PcjoNXA9wFhRlETABHiQ8ZSkzJbLoQxeJ/7m0qyVhI9bAAcCeW1nT44SNmqNoygKuSvakruiLfaTB/Lq/FXu+uznns9+QgZMIWTAFPRL23Cxy3UKqbwwLfrfdlOSJP2/dPfEBFAUpT8wmbc9MfcLIdp9ZBvZEzMLxU6885BnfWegmyc3FouGomNkkClxk+49JvZoCDH+p0m5/nalRL3i1hhXc8SoWgX0C1tqNP739nPWZlxtxs6qcdXeExPIBRwG8gD6wFag/ef2kT0xs0bsPdMXiPVKKRHUaahISUnJ1Nh+fn4i6uZdcfG35WKfyw/ve4HuKFNPhI2aI56HXdJITt/jz1kec9aJiwZ6YtYCbgkhngohEgFfwDUD40lfCaPKttiO6cOt1Vu4sWxzpsc3LWJNmUFdqBO0kaYRR6m4YDRG+cy5OOVP9jg0YUeJOoQOncmzk+e+24WPJAkydg38LlBFURQT3l5CqQmcVktWktbZju7Ns+AwTvedSO6K5chdoZxW8jApmI9SP7en1M/tiXsS+W59ln1cnr2KSzOWYWJtiXXz2lirvLBwdURHV1creUqSNqT7DFwIcQLwBkKA8HdjLVFTXpKW6ejq4rp+JkZ5zQlQ9f0qGiUb5TWneLdWeO5djupJEFVWTSOXQ2mu/bGRg9XasdWqOqd6j+PRoWBSkpK0na4kaVyG5oELIcYKIUoLIWyFEB2EEKlbKUnKEowscuO2eR6x958Q1HEoIiVF2ym9Z5DLjKKdmlF9+x+ongbjumE2edwqcHP1Vg7X6syW/FU53nUE93cfITkhQdvpSpJGyDsxpc+ycLbHcfYwHuz04+L0pdpO56P0s5ti07oB7pvno3oajLvvQvJ7uXF3816ONOiObx4XgtoP5t7WgyTFxmk7XUlSG7kWivRFJfu049mxEM6NmotFFXvyeVbRdkqfpGdijHWz2lg3q01yfAKPDgZxz2cfEdsOc3v9DvSymVCgfjWsVV4UqF8N/ezamUYnSeogC7j0RYqiUHnpRF6cvcyx1r9QN2QLJgXzaTutL9I1NKBgAw8KNvAgJTGRJ0dOcddnHxG+B7i7eS86hgZYerlhrfIiJXfmzHeXJHWSl1CkVNE3zYa7zwKSYmI59sOALNfEQUdfn/y1XKm8eDxNHwRQ68g6ivf4gRchFzneaSiPmg3Gr95PXF+2mbinz7WdriSliizgUqqZlSlG5WWTeHoshLBhs7SdTrrp6OqSt1olnOaNoskdP+oc/5tsqppEXb3NyW6j2JK/KodqdOTqovW8efBY2+lK0ifJAi6liU3rBpT8uT2XZ69U2+qF2qTo6GDhbI9ZTxWNrh+gbsgWyg7vTuzDp5z+eQJbrapzwK0Nl+esIubOfW2nK0kfkAVcSjPHWUMxd7bn+I/DeX31lrbTURtFUcjtWBb7SQNpeGkPDS7swm58XxKj3xDyy1S22dRgbyUVF6Yt4fW129pOV5JkAZfSTtfAALfN89A10CdA1e+bbZ1mVrY4dqP7UD9sG42u7cdh+mAUHR3ODp/FzpJe7C7fiPDxC3l54Zq8pV/SClnApXTJZm2J64bZvLpwjZO9xn3zBSx78cKUHdINrxObaXLHjwpzhqNvlp3w8QvZbduQXWXqcXbkHJ6HXPjm/y2kr4cs4FK6Wdauit24n7m9dhvXl2zSdjqZJluhApQe0JnaAX/R7P5RKv0+FmOr/FycvpS9FZuzvVgtQgZP52lw6Fd196r07clIR55SiqKE/ePrtaIoA9SZnPT1sx3VG8u67pzpN4nI0+HaTifTGVvmpUSvttQ8uIpmjwJxXjaJHKWLcnX+Wg64tmZrIQ9O95vE4yMnSUlO1na60jcmI4tZXRFCOAghHICKwBtgi9oyk7IERUcH13UzMcqfh8AW/TPcHDkrM7LITbGuLfHcvZTmT4JwWTsD80p23Fj6N4c8OrC1gDsne4zh4f7ALDePXvo6qetOzJrADSHEHTWNJ2Uhhua5cNs8j4NubQnuMITqO/5A0fm+r84Z5MxBkfZNKNK+CYnRMTzcc5S73vu4vX4H15dswiCXGQUb18BaVQdhKM/MpfRRVwFvDWxQ01hSFmRRuTwV5g7ndJ8JXJj6J7Yje2k7pa+Gvmk2CrWsR6GW9UiKjePR/kDu+uwnYutBbq3egmJixLEmtbBW1aFAXXf0sploO2Upi8hQT0wARVEMeNvMuJwQ4j+3rcmemFkvdnrjCiF4OXkFsX6nMJ/RH8OKZTItdkZpI65ITCI+9DJRh06SfOICKa+iUQz1Maxsi1G1Chi52KGTzVhj8eX/7awTV+09Mf/3BTThbUPjL24re2JmjdgZiZsYHSN2lq0vvPNUETERjzI1dkZo++ecnJgoHh4KEid7jxO+llXFekqKDQblhF/DHuLGSh8RF/lCI3G15Xv7OX+NPTH/pw3y8on0jl42E9x85pMcG0dgq/7yw7pU0tHTI38NFyotGkvTiKPUDvyLEn3a8fLcFY7/OBzffFU5XKcL1/7cSOzjZ9pOV/pKZKiAK4qSDajN24bGkgSAWeliOC+fzLOgUEKHzNR2OlmOoqNDnqoVqTh7OE1uH8brlDdlBnch+lYEp3qOZWsBdw5Wb8+VBWt5c18utvU9y2hLtRghhLkQQvsNE6WvSuFW9SnVvyNX5q7m7uY92k4ny1IUBXMnOxymDqLR1X3UO7uNcqN6ER/5kjP9JrHVqhr7XH7g0m/Lib51T9vpSpns+57rJWmUw4xfsXBx5HiXEby6fEPb6WR5iqKQq3xpyo/vR4PzO2l4eQ/2kweSkpBI6K8z2F60FnsqNOP85MXy3/s7IQu4pDG6Bga4/T0XXSNDAlv0/2YXvdKWHKWKUm5ET+qd8aXxzYM4zhyCjqEB50bNZVeZ+uyybci5sfN5ce6yXJ/lGyULuKRRJlb5qbphNq8uXudkjzGykGiIaRFrygzuilfwJpreO0LF+aMwtMjFhUmL2WPfhB0lvQgb9huRp87Jn8E3RBZwSePy13Kl/IR+b+9C/ENOWNI0E6v8lOrbgVr+a2n2MJDKf07AtKg1l2atZF/llmwvUpMzv0wlPvy6XGwri5NNjaVMUW5ET54Fh3FmwBRyO9liXqm8tlP6LhjlNad49x8o3v0H4p+/5P72w9z12c+1RetJSUhk69RVWDWrjbWqDnmrVUJHT5aErESegUuZQtHRwWXtDIwt8xLQoj/xkS+0ndJ3xzB3Top2bo7Hjj9QPT1OzlFdsXCtwM2Vvhyu2Zktlm6c+GkkD/YeJTkhQdvpSqkgC7iUaQxz58TNex5xj54S1GGIfPuuRfo5TDGpWRl37/monh3H3WcB+WtX5c7fe/Cv1w3fvK4EdRxCxLaDJMXGaTtd6RPk+yUpU5k72VFx3khO9RrH+cmLsRvdR9spfff0TIyxbl4H6+Z1SI6L59HBIO757Cdi2yFur92GXjYTCjSo/naxrfrV0TfNpu2UpXdkAZcyXfEerXl6LITwsQuwqOKAZe2q2k5JekfXyJCCDT0p2NCTlMREHvuf5J7PPiK2HOTu33vQNTLEsq471qo6FGxUAwOz7NpO+bsmC7iU6RRFofIf43kReomgtoOoG7KFbNaW2k5L+hcdfX0sa1fFsnZVnBaN5WngGe757Oee79ulcHX09clXy4VCqjoUbFITI4vc2k75uyOvgUtaoZfNBHef+STHJxDYaoD80Owrp6OrS77qlXGaP4qmd/2pE7yJUv078PryTU78NIot+d04VLMTV39fT+zDJ9pO97uR0cWsciqK4q0oymVFUS4piuKirsSkb1+OUkWpsmIKkcfDCP11hrbTkVJJ0dHBoooDjjOH0vjGQeqGbKHssG7E3n/M6T4T2FKwGgfc23J57ipi7j7QdrrftIxeQpkH7BVCtHjX2EG2EpHSpFCLupQa2Jkrc1Zh4eII+eUHZFmJoijkdixLbseylJ84gFcXr3PPZx/3fPYTMnAqIQOnkruSHYVUdbBWeZG9eGFtp/xNSXcBVxTFDKgGdAYQQiQA8n2wlGaO0wfz/OQ5Tv40CpP2dbkcdjvTc4hLiAYPj0yNKYRgZ3ggSkJspsYFeH3tNgkXbqj1mBVFIWe5EuQsVwK7MT/z+tptInz3c9dnP2HDZhE2bBam9iW5XL8UJYqWy/R/729RuluqKYriACwBLgL2wBmgvxAi5l/byZZqWSy2NuImP33Bs74zSH78PFPj/lOu8T0wrlYh0+J5Xw1kUehOypvbMLtGd3QzqRF0cuQrnnafRMrrGCzm/YpB2SIaj5n0KJK4gFD+uBfAroLx5IpJ4adLBriVdsK4WgX0ilujKIrG88iqv1OfaqmWkQLuBBwHqgohTiiKMg94LYQY/al9nJycxOnTp9MVz9/fHw8tvGJrK642Y2srbkpiIv77DuDm5pa5cZOS2OHeBh48w+u0DzlK2Gg8ZtCNc1Sf3Yviea24/OgOI+p2YnITzTeCTklK4nDNzkSeCofsxhgaGlE3xDdTZpBsPnOIVstG0rKcO6cuhnFbRFHlRgIdA95QsIAV1s1rY92iLuaV7FA09GKWVX+nFEX5aAHPyDXwCCBCCHHi3ffewLAMjCd953T09dExNcEgZ45Mj51rbHde9plBoKovdY7/jZ6J5poJP3n9nFbLRlLY3JLgX5fRftEIpuxdjUtROxraafbF6+zIOTw5egqXtTO4FvOS5/1mEtz+V6rvWoKOrq7G4l55dIcuaydTpYgt63pO4ciRI5yIv8vE3Su4XDIHvR7nJGbeWi79tgITq/xYNa9NIddBy+IAACAASURBVJUXFlUraDSvrC7dL3NCiEfAPUVRSr17qCZvL6dIUpajl98c1/UzeXn+Gqd6j9fYkqvJKcm0XTGGyJjXeHebQk6T7PSr0BhH65J0WDmeW880N2sjYttBLs1YRvGerSnSvgkGpQrjtGA0D/cFcmHSYo3FjYmPRbVkOEb6BmzuNgUDPX30dfUYVb8LoSPXUrpQcaZaRLBmRg0KLx1J7orluP7nJg5Wb8/WgtU42XMMjw4Gyf6qH5HR9yl9gfWKopwDHIApGU9JkrSjQN1q2I7pw63VW7ixbLNGYozbuYxDV06zqPVgHKxLAmCgq493t6kAtFgynLjEeLXHjbpxl+BOw8jtZEvFuSPfP16sWyuKdGxK+PiFPNgXoPa4Qgh6/jWdi49u8VeX8VjlyvvB82UtixA4+E/mthzI0VvnqHdhFRcG1qTZkyCqbppD3uqVuL1uB4dr/4hvfjeOdxnO/V3+JMfL+RKQ8Z6YYUIIJyFEeSFEUyGEXGJOytJsR/cmfx03TvedyPOQC2ode/f5ICbtWUkX10Z0cW30wXNF8xRkTecxhNy7Qv+/56g1blJsHAGqvig6OrhtnoeuocH75xRFodLiceS0LUFQ28Fqn7f9Z8AW1p3cy/iG3ahdxvmj2+jq6NK/xg+cH/0Xzjbl6L1xJnWWDSHBwxa3TXNp/jQY9y2LKFC/Gvd8D3CkYQ9887pwrN0g7vnuJ+lN5s/i+VrIOzEl6R90dHVxXT8To7zmBKj6Ev/8pVrGvR35gPYrx+FgVZKFPwz66DaNyrszzKsjSwK3sjp4l1riApzuM56XZy/jsm4mpjZW/3lez8QYN+/5pCQmEtiyv9rObk/dvkj/zXOoV86FkXU7f3H7IhYF2N9vPss7jORsxDXKT+rAjP1rwUAP66a1cF07k+ZPgvDYvYRCLevyaF8gAaq++ORxIaBlP25v3EViVLRacs8qZAGXpH8xssiN2+Z5xN5/QnDHoRle9jYuMZ4WS0aQIlLw7j4FYwOjT247sVF3PEtWpOeGGZyLuJahuAA3lm/m5kpfyo3qRcH61T+5XY6SRXBZNY3Ik+cIGTQtw3Ejo1/RYukI8ucwZ23nceikclaJoih0cW3ExTEbqFu2CkO3LMJ5elfOvvu30DUwoEC96jgvm0yzR8eocWgVRTs15WlgCEFtfsEnjwtHGvfk5uotJLx4leHj+NrJAi5JH2HhbI/j7GE82OXPxWlLMjTWgM1zOXP3Mqs7jaFYnv+eAf+Tnq4eG7pOIJdJdlRLhvMqNv1nlM9DL3KqzwTy13LFblzfL25v3bwOpX/5kWuL1nP7rx3pjpuSkkKHVeN49DoS725TMDc1S/MYBXLmwbfHNDZ3m0LEyyc4Te3M6O1/Ep/4/+8OdPT0yF/DhUq/j6PZ/aPUClhPiV5teBF2meOdh+GT15XDXl25vmQTcU8i0308XzNZwCXpE0r2aUfh1g04N3oejw4Fp2uMtSf28GfAFobUaU8T+2qp2idfDnM2dZ3ErciHdFkzKV0zYhJeviawRT8MLXLh+tesVE/Fc5g2iDxuFTnRbTSvLl5Pc1yAyXtXsedCMHNbDKCSTdl0jQFvz8ZbVKjBxTEbaVvJi0l7VuI4pSPBN8P/u62ODnndnKg4ZwRN7vjhdXIzZQb9SPSNu5zsMYYtlm4c9OxAzBY/3tx/nO6cvjaygEvSJyiKQuWlE8leqgjH2vyS5l/88PvX6bF+GtVLODK5cc807etewoHpzXrjG+bP7EN/pWlfkZJCcKehxNx9iNvmeRjlSf1NOjr6+lTdNAd9UxMCVH3TfE35wKUTjN25lHaVvehZrXma9v0Uc1MzVncew56f5xCTEEfV37oz4O85RMe9+ej2iqJgXqk8DtMG0+jafuqFbaPcyJ7EP3nOq/kb2WpVjf2urbk0awXRtyPUkqO2yAIuSZ+hb5oNd58FJL+J49gPA1I9F/l1bAyqJcMxMzZlY9dJ6Omm/Z65X2q2pbmDB0O3/E7AtbBU73dp5jLubz9MhVlDyePimOa4JgXyUXXjbKKu3uZEt9Gpfgdw7/lj2q4YS9n8Rfiz7TC13xpft5wL50evp3c1FfP8NmE3qR0HLp347D6KopDLvjTlJ/SnwYVd5Fk9nvKTBpAcF0/o4OlsL1KTvU7NuTD1T15fvaXWfDODLOCS9AVmZYpRedkknh4LIXTob1/cXghBl7WTuPnsAZt+mkR+M/N0xVUUhZUdR1PUogA/LB/Fo1dfvo772P8EZ0fMoVCrepTs2yFdcQHyeVah/OSB3N20m6sL131x+4SkRFotG0l8UgI+3aeSzVAzd7JmN8rGwtaDOfrLHxjo6lNnfn+6rJnEi5jXqdpfv1B+bEf2ol7IFhrfOIjjzCEoenqcHTGbnaXqssuuEefGLeBl+BWN3cylTrKAS1Iq2LRuQMmf23Nlzirueu/97LZzD2/EJ9SPqU16Ua1E2s+A/ymHcTZ8uk/l5Zso2qwYTVJy0ie3jX34hGOtfyF7SRucl03K8Blw2SE/UbCRJ6GDpvPs+OffAfzqu4Djt86zosMoSuXX/JKx7iUcODtqLcO8OrLmxB7KTmjDljD/NI1hWtSaMoO74nX8b5reO0LFeSMxzG3G+QmL2F2+MTtL1SVsxGyenzn/1RZzWcAlKZUcZw3F3Nme411GfPLt9rEbZxniu5Cm9tUZXLudWuLaFSzOH22H4n81hNE7Pj4jJiUxkcBWA0iMisHNez762TO+4p6io4PL6ukYW+UjsGV/4p5+fKXITacPMN/vbwbUaE2LCjUyHDe1jPQNmdq0NyeHLid/DnOa/zmMlktHpOqdyr+ZWOWnVL+O1DqyjmYPA6n0x3iy2RTk0oxl7HVSsb1oTUIGTeNpUEiGp5WqkyzgkpRKugYGuP09F10DfQJU/UiK+fBDtCevn9Nq6ShszC1Z1Wm0Wq8Bd6xSn+5uTZm2bw3bzx79z/Nhw2fzNPAMzksnkrNcCbXFNchlhrvPAuKePieo7SBSkpM/eP7Sw1t0XTcF16J2zGj+s9ripkWFQqU5OWwFU5r0Yse5QMpOaMPq4F3pPms2zmdBiR6tqbF/Bc0fH8N5xRTMypXg6sJ1HKjahq3W1Tn18wQe+x0nJenT74gygyzgkpQG2QoVwPWvWby6cI2Tvca9LxLJKcm0WTGG529e4919KmbG6l9zel6rgVSwLkXH1RO48fT/Z0/c893P5VkrKNG7LTZtG31mhPTJ7VgWp4VvF5Q6P37h+8ej496gWjIcE30j/u42Gf10fFCrLvq6egyv24mwkWsok78wnddMpN7CgdyJfJihcQ3Nc1HsRxUeO/9E9fQ4rut/w6KKAzdX+HCoRie2FHDnRPfRPNgXoJXFtjLaE/O2oijhiqKEKYqSvoW+JSmLsazjht24n7m9dhvXl2wCYMyOJRy+cprFbX7F3kp9Z8D/ZKRviHf3qegoOrRYMoLYhDheX7vN8R+HY165PBVmD9dIXIBiXVtQtHNzzk/8nQd7jiCEoPv6aVx5fJcNXSdQMGfeLw+SCUrntyFg0J8s+GEQgTfOUm5iWxb6byZFDZc99HOYYtO2Ee4+C1A9DcbNez75a7pwZ8NO/Ov+hE9eV4I7DSVi+yGS49S/INnHqOMl01MI8UwN40hSlmE7qjfPgsM4028Sp/IkMeXAan6q2pjOLg01GreIRQHWdh5Lw98H0Wf9DFrMOYWip/f20s4/FqlSN0VRcFo0huchFwhqP4S7S7ux4fR+JjXuQc3SlTQWNz10dHT42aMlDe2q0mP9dPpumsXG0wdY3n7kl3dOJb1sJhRSeVFI5UVyXDwPDxzjns9+7m8/zK01W9EzNaFAAw8KtfCiQL3U3cCVrjw0NrIkfcMUHR1c181kZdWmDN69EIdCxVjwiUWq1K2BXVWGe3Vi6r7V6Ce+Yez6hWQrXFDjcfVMjHH3WcCCOirG7ltGfVsXhnt10njc9LIxL8DevnNZc3w3A73nUn58G+xfG2GyZrpmAuqAaJyb+GcKKQmJ8PwILDmCslyX3Lly41K5MoYm6u37nu6WagCKotwCXgAC+FMI8Z+PyGVPzKwXWx5z6iQkJ9Jn1zwevHzKzOv5KTOuf7pagaUn9usdRxlxdRvXChqyyOtniucqkClxX8XH0G3HLHgRxawYO6wHdExz3PTGTq+E8ze4PW81awpFcSe/IWi49+YLI0GM/oePFY3WZWGbERibZEvXmJ/qiYkQIt1fQMF3f+YFzgLVPrd9xYoVRXr5+fmle9+M0FZcbcaWx5w63dZNEfR0FounThTrKSnCJy7KlNiRZ86LDYa2wrtuR1FwWENRbLRKvIh5rfG4SclJwmt+f2Hws5tY98swsZ6S4ua6bWmOm57Y6ZEQFS1O9Zso1iulxJZCHuL+3qMai3vl0R0xZc8qUXFKJ0FPZ0FPZ2E/qb2YuGu5uPjgZobjAqfFR2pqhi6hCCHuv/vziaIoW4DKwH/nOEnSN2Z18C6WBm5jmFdHejTpRVD4C86NmY9FFQfy13LVWNyEF68IUPXFKK85ddfOwfJVBNVn9+LHNZPw7TFNo53dJ+1eyb6Lx/mz7TDauDTk8Ol7nOw+hlwOZdQ6dVEdHh44xsnuo4m5fZ+SP7fHfspA9LObctXfXy3jCyG48OAm3qGH8Qn15/yDGwBUtinL9GZ9UDl6frDy5OMrd9QS99/SXcAVRckG6Aghot79vQ4wQW2ZSdJX6lzENXpumIFnyYpMbNQdRVFwXjKRl2GXOdbmF+qFbsXEKr/a44qUFII6DiX2/hNqBazHyCI3rha5mdm8LwO95/LbgfX8Wqe92uMC7Lt4nPG7l9PRuT7d3JqgKApVN85mj2MzAlV98TrlrZabhzIq4cUrQgZN4+ZKX3KUKkKtgPXkdfvvlYf0EEIQcvcKPqF++IT6cfXJXRRFwa2YPXNbDqS5gwfWufOpJVZqZeQMPB+w5d0rvh7wlxDi8/cYS1IW9yo2GtWS4eQyyc6GrhPeL1Kll80EN5/57KvUgsBWA6jpvwZdA/XOCrk4fSkPdvpRccFoLJzt3z/ev8YPBN0MZ/i2xVS2KUv1khXUGvfu80e0WzEW2wJFWdx2yPuzfGPLvFTdNIfDNTtzoutIqm6aq9F3AF9yb8sBTvUeT/zT55Qd3gO7MX3QNTLM0JgpKSmcuH0B75DD+Ib5czvyIbo6uniUcGRAjR9o5uCR7rVu1CHdBVwIcROw/+KGkvSNEELQZc3bdbr9BiwiX44Pf3HNShfDeflkjv0wkLAhMz9oHpxRjw4Fc27UXAq3bkDJPh/eoq8oCsvaj+BsxDV+WD6K0BFrsDSzUEvc+MQEWi4dQUJyIt7dpmLyr25C+apXxn7KQMKG/oaF62pKD+islrhpEfvoKaf7TuSe9z5yOZTBY/cScjumfx3y5JRkAq6H4RPqh2/oER68eoq+rh61y1RmdP0uNC7vjoVpTjUeQfrJaYSSlEqzD/2Fb5g/v6n64l7C4aPbFG5Vn6fHQrgybw0Wro4UblU/w3Hf3H/MsTa/kL1UESovnfjRs9z/LXrlPKMrrZeP5lD/BelawvbfBvnM5+Tti/h0n0rJfIU+uk2ZX3/iWVAoob/OxLySHXmqVsxw3NQQQnBrzVZCBk4l6U0s9lN+oczgLujo6395539JTE7i8JXT+IT6sTXsKE+jX2Ckb0jdslVoUcGThnZuGrm7NqNkAZekVAi4FsbQLb/T3MGDX2q2/ey2jjOH8PzUeU50HUnO8qUwK10s3XFTEhM59sMAkt/E4e6zAH3TT09Dsy1YjD/bDqXDqvGM2LaYGc2/3Ebtc/46uY9FR7z5pWYbmjt6fnI7RVGosmoae51UBLYaQL3QrRjl1exlhZg79znZYwwP9wWSp2oFKi+blOZ/57jEeA5cOolPqB/bzgbwMjYKU0MTGti6onL0pF45F0yN1DtvW91kAZekL3j0KpIflo+iqEUBVnb88iJV/1v0ao9jUwJb9MfrxN/oZUtfIQgd+htPj4XgumE2ZmW+XKDaO9fj2I1zzDywHtei5Wnq8OlGxp9z8eEtuq2filsxe6Y16/PF7Q1y5sDdez77XX7gWNtBeO5bnuo2bmkhUlK4+vtfnB02C4CKC0ZTsnfbVM+/j01KwDvkMD6hfuwMP0Z0/BvMjE1pUr4aKkcP6pR1xkg/Y9fNM5Ms4JL0GUnJSbRZMZqXb6LY13cuOYxTdyOGiVV+qm6YzeE6XTjZYwwua2em+QO+u957uTJnFSV/bo9N6wap3m9uy4GcvnuZTqsncKbAKorntU5T3Ki4GFRLhmFqaMymnyalepGqXA5lcFo0hhNdRxI+dj72kwamKe6XvL5ykxNdR/L0WAiWXm5U/nNCqu5AfRUbzc7wQHxC/dgdHkR8ciIWpjlp7VQLlaMnNUo5YaCX9ssuXwNZwCXpM0bvWIL/1RBWdxqDXcHiado3fy1Xyk/ox7nR88hTtQIlen3+0ss/vb56i+NdRmDubI/jrKFpimuob8DmbpOpOLUzLZaOIPjXpRj/68PHTxFC0G3dVK4+vsfB/gsokDNPmmIX69KCZ0GhXJj8BxYujhRs4JGm/T8mJTGRS7+tIHz8QvRMjKmyejpFOjT57AtiZPQrtp8LwCfUjwOXT5KQlIilmQX1ijjRr2F73Ivbq+UzAm3L+kcgSRqy/exRpu1bQ3e3pnSskr4PI8uN6Pl20asBU8jtZIt5pfJf3Ccp5g0Bqn7oGhrgtnleuqYj2pgXYF3ncTT4fRB9Nv7Gio6jUrXfQv/NbDpzkKlNe+NZKn0fRlZcMJrnZy4Q3GEIdUN8MbWx+vJOn/A89CInuozgRdglrFt44bRwDMb5Pj7D5vHrSLaEHcEn1A+/qyEkpyRTKHd++lRX0cKxBlWK2HL06FE80nlcXyNZwCXpI248jaDj6glUsC7FvFbpvxSg6OjgsnYGeys0J6BFf+qF+GJonuuT2wshONlrHK8uXMNz33KyWVumO3Y9W1dG1fuRibtXULVYebpWbfzZ7Y/fPM8gn/k0snNjSO303xCkZ2yEu88C9lRoRmCLftQO3JDm+djJcfGEj1/IpZnLMcyTG3efBVg3r/Of7SJePME31A/vUD8Cb5xFCEGJvNb8WrsdKkdPKhYqrdW56ZomC7gk/UtsQhwtloxAR9HBu/vUDH+oZZg7J27e8zhQtQ1B7X/FY9eST37odn3JJm6v3Ybd+L5Y1q6aobgAYxt05fit8/TZ+BsVCpXC0brUR7d7GvWClstGYJUzL6s7jUEnHYty/ZNpUWtc1kznaJPenBkwmcp/pP4m7SeBpznRdSRRV29T9MfmVJg1DINcZu+fv/n0/vu7IU/cvgCAbYFijKnfBZWjJ7YFin3TRfufZAGXpH/pu2kWYRFX2dl7FkUs0r7K38eYO9lRcd5ITvUax/lJv2M35r/txyJPh3Om3yQs67pjO6q3WuLq6uiy/sfxVJjaiRZLRnB62EpyZcvxwTbJKcm0WzmWp1EvCR6y9D/Pp5dV45qUHdqNi9OXYuFagaIdm352+8SoaMKGz+baovVksymI5/4V71/ELj289a5o+xMWcRWACtalmNy4JypHz0xppPw1kgVckv5hZdBOlgftYETdTjSwy/gZ8D8V79Gap8dCCB+3EIsqDljWcXv/XHzkCwJb9Mcofx5c181M17K0n5Iney7+/mky1Wb1pNPqCWztOeODM+zxu5Zz4NJJlrYb/skz9PQqP2kAz06c5VTPseR2LENOu4+P/2BfACe7j+HNvYeU7NeB8pMGcOnVQxbvWIJPqB8XH75tIu1S1I7fVH1p7uCpthfXrEz2xJSkd8LuXaX3xpnUKOXEhEbd1T6+oihU/mM8ZmWLE9R2EDH33vZrFCkpBHUYQuyDJ7htnvfZa+Tp5VLUjlmqfuwID2TGgXXvH99zPoiJu1fQ2aXBF6+Rp4eOnh5VN8zGIGd2AlT9SHgV9cHz8c9fEtxpKP51f0LXxBDzrVPY7JGTsjM74zC5A5P3rCKPaU7mt/qFe1O2E/TrUgbVaieL9zsZPgNXFEUXOA3cF0Jotp+UJGlIdEIsPy0dQW6THGzoMgFdHfXfhAJvF71y95nP3kotCGzZn1pH1xG9fg9Re47itGgMFpW/PEslvfp6tuLYzXOM3PYHzjbleBzzgj47F1O+YHEWtf5VY9eNjfPnoeqmuRzy7MiJLiNw854PvJ3nfrLPeML1Y7gxwJkA41fc3TMDPR1dapRyYkjt9jS1r0beHLk1kte3QB2XUPoDlwD1XDiTpEwmhGD6yc3ciXyI/y+/a7xg5ChVlCorphDYsj9HGvUi6sAxCrdtmKZ54unxv0Wvzt2/TuvlozDTMyEpORmf7v9dpErd8ro74TBtEKG/ziBk9Bz8zxxhcswdzjQ04YWBCQZJt6hjVZkJjbrTqLwbubOZfXlQKWMFXFEUK6ABMBn4RS0ZSd+lqLgYmv0xlNA7VzDYMzNTYyenpPA0+gWzW/SnarHMWWCzUIu6lBrQiStzV6NX2BLnJR9fpErdshtlw7vbVCpP78KTqBf49piW5js10yMhKZFbdcuw/nwJjkZsJNpGB2PFlHr2brSoUIMGtlVTfZer9P8y2hPTG5gKZAcGf+wSiuyJmfViZ3ZcIQTjg/8iIOICNQvaYWiQ+WtR5Dcyo61tjUydfiaSkonefBBRqTQ5imfuLIqwJze4+/wxjUtrrntQfFIipx5d5WjEeYIeXCImMQ4TPQOc4nLgUtwRj4ruGOmpd830L8mqv1Nq74kJNAR+f/d3D2Dnl/aRPTGzRuzMjjv30EZBT2cxfd+a7+aYv4bYmogbFRsjNp7aL1ouGSGy9fcQ9HQWuX6pLTqvniB2nAsQcQnxGoudGlk1LhroiVkVaKwoSn3ACMihKMo6IYRmejpJ36SgG+cY7DOfJvbV+LV2e44cOaLtlKQ0evkmih3vFovae+E48UkJ5M2ei/aV66Jy9MCjZMVUL4glpU1GOvIMB4YDKIriwdtLKLJ4S6n25PVzWi0bSWFzS1alYplW6evxNOoF284exSfUj0NXTpOYnETBnHno4d4UlaMnVYuV19hMHun/yZdFSSuSU5Jpu2IMkTGvCf51KTlNsms7JekLHr56xpYwf7xD/DhyLZQUkUIR8wL09/wBlaMnlW3KZvgWfClt1FLAhRD+gL86xpK+D+N2LuPQldMs7zASB+uS2k5H+oQ7kQ/xDfPHJ9SPoJvhCCEonb8ww706oqrgiYNVSfnOSYvkGbiU6XafD2LSnpV0cW1EF9dG2k5H+pdrT+7iE/q2aJ++cwmA8gWLM67BT7SoUIOylkW0nKH0P7KAS5nqduQD2q8ch4NVSRb+MEjb6Ui8nYl28f1iUX6cu38dgEqFyzKtaW9Ujp6ZMldcSjtZwKVME5+YQMulI0kRKXh3n5LqLjGS+gkhuPr8Pge2LcYn1J8rj++gKAquRe2Y3aI/zR08KGye/rXIpcwhC7iUaQZsnsPpO5fY2nMGxfKkv0uLlD4pKSmcvH0R79DD+Ib6cyvyATqKDh4lK9DPsyXNHDywNPt4txvp6yQLuJQp1p7Ywx8BWxhapwNN7KtpO53vRnJKMoHXz+IT6odvmD/3Xz5FX1ePWqUroSrqwtBW3bAwzantNKV0kgVc0rjw+9fpsX4aHiUrMKlxD22n881LTE7C/+oZvEP82Hr2CE+iXmCkb4hXWWemNu1NIzs3cppkx9/fXxbvLE4WcEmjXsfGoFoynJwm2dnQZeI30Qn8axSfmMCByyfxCfVj29kAXrx5TTZDYxrYuqJy9KR+OVdMjUy0naakZvK3SdIYIQRd1k7i5rMH+A1YRH4zc22n9E15kxDHnvNB+IT6s/N8IFFxbzAzNqWRnRsqR0+8yjrLD4q/cbKASxoz9/BGfEL9mNm8L+4lHLSdzjfhdWwMu84fwyfUj93ng4hNjMc8mxmtKtRE5ehJzdKVMNDT13aaUiaRBVzSiGM3zjLEdyHNHKozqJZmGxV8657HvGL7uQB8Qv3Yf+kkCUmJ5M9hzo8uDVE5elKthIO8NPWdSvdPXVEUI+AoYPhuHG8hxFh1JSZlXU9eP6fV0lHYmFuyUi5SlS6PX0ey9d1iUX5XzpCUkox1rnz0rqZC5eiBa9Hyct0RKUNn4PFADSFEtKIo+kCgoih7hBDH1ZSblAUlpyTTZsUYnr95ze6fl2FmrJ1mGFnR/ZdP8H13C3vA9bOkiBSK57FiUK22qBw9cSpcRr4YSh/IyHKyAoh+963+u6/0t/eRvgljdizh8JXTrOw4CnurEtpO56t369kDNl0+yvCT6zh+6zwAZS2LMLJeZ1o4emJXsLgs2tInZbQnpi5wBigOLBJCnFBLVlKWtDM8kCl7V/NT1cZ0dvlPdz3pnSuP7uAT6od36GFC710FwNG6JJMa90Dl6Enp/DbaTVDKMjLUE/P9IIqSE9gC9BVCnP/Xc7InZhaLnZ64D6Of033/fCyz5WZhrV4Y6KZvJkRWOubUEkJw89Ujjkac5+i989x+/RiAMrmtqWZti1PuolpZLEr+3846cdXeE/PfX8AY3nblkT0xs3jstMaNTYgTFSZ3FDkH1hI3nkRkamx1UXfclJQUcfLWBTHUd6EoPlol6OkslF5VRLVZPcW8QxvF3chHGoudWvL/dtaJi7p7YiqKkgdIFEK8VBTFGKgNTE/veFLW1e/v2YTcu8KO3r9RNE9BbaejNSkpKQTfCn+77kioP3eeP0JXR5capSoyuHY7mtpXI18OeTOTpD4ZuQZuCax+dx1cB/hbCLFTPWlJWcXq4F0sDdzGcK9ONLRz03Y6mS4pOYmj18LwCfVjy9kjPHz1DAM9fWqXrszYBj/RuLw75qZm2k5T+kZlZBbKOcBRjblIWcy5iGv03DADz5IVmdCom7bTyTQJSYkcvnIan1A/mkw9QgAAFeZJREFUtp49yrPolxjrG1KvnAsqR08a2FWV0yelTCFv35LS5VVsNKolw8llkp0NXSd883cCxibEsf/S28Witp8L4FVsNKaGJjS0q0oLR0/qlnMhm6GxttOUvjPf9m+dpBFCCLqsmcStyIf4D1z0zV7XjY57w54LwXiHHmbX+SBi4mPJaZydpvbVUDl6UrtMZYz0DbWdpvQdkwVcSrPZh/7CN8yfWap+uBX/thapehUbzY5zAfiE+rP34nHiEuPJY5qLdpW8UDl64lmqIvrf+LsNKeuQ/xOlNAm4FsbQLb+jcvRkYM022k5HLZ5Fv2TXzVPMOL+Ng5dPkZicRAGzPPxUtREqR0/cizugq6Or7TQl6T9kAZdS7dGrSH74v/buPD6q+urj+Ock7PsSZN9FWYIwqICImIAIgiAQ9MGKYsEiD4JKrVZpXVrrglWrVqsPFdQiWJUBUVREdAKyCAqTQAi7bGERZAkQIGQ5zx9zQVD2zNzhJuf9euWVO0Pmfn8/CCc3d+49v3F/plFcLcbf/mdP3+K9PXMXU1OSmRwMMHtNkLz8PBpUrcm9ibeQ5EukXYMW1izKXPCsgJuzkpuXy63jH2Xvwf18MfIlKpQuG+0hnbNNu7cfaxY174elqCqXXFSPh7oOpH5+eYb2vc3TP5RM0WMF3JyVRz8ZS/LqJbwz6DFa1r442sM5a2t3bMYfDOAPJvPdxnQAWtZuzOM9hpDkS6RFrUaICMnJyVa8jedYATdn9HHqHJ794j8M7diHO9r3iPZwzih923qnaAdIzVgDwOX1mvJMn+H0a53AJdXrRXmExoSHFXBzWut2ZnDHO3/l8npNefmWUdEezkmpKikZq/EvCeBPCbBy+0YAOjRqyQtJ99LPl0CDqrWiPEpjws8KuDmlQ0cO03/saGIkhslDn76grnlWVRZtWH7s9MgPP20hRmK4tomPEdfeTN/W11KrUrVoD9OYiLICbk5p5PsvkJKxmunDX7ggjmDz8vOYt25pqFlUSjIZe3ZQLCaWLk2v5OFut9On1bVUK1852sM0xjUF6UZYF/gPUJ3QSjxjVfXlcA3MRNdb86czbv4n/Kn7nfRseXXUxpGTl8vs1UuYHAzwUepsfty3m5LFStCteTue6j2MXi07UrlshaiNz5hoKsgReC7wgKouEZHywGIR+VJV08M0NhMla/dsZWTg/+hy6RX8JQpNqrJzjjBr5Xf4gwGmLZ3D7qx9lClRip7xV5PkS6BHfAfKl/LeZYzGhFtBuhFuA7Y52/tFZAVQG7AC7mG7DmTy+PyJVC1bgUmD/+raHYgHjxxmxvIF/GvB+3w37Un2Hc6iQqmy9LqsI0m+RLo1b0+ZEqVcGYsxXhGuJdUaAHOAeFXd94s/syXVPJI9N2M5Ly3+iMzsLP7R+W7i4+pHNO9gTjbfblvJnM1pLNy2ksN5OZQvXpqOdVrQqU48bapfTAmX+o4UpX/naOdGM9uruRFbUg0oR2hh435n+lpbUu3CzN6e+ZPePHa0MqydtvrbQH3DPyFiWbsPZOrb86drr9ce0JIjrlGGtdPqD92gwyY+q7NWLNJZX82KWPbpFIV/5wslN5rZXs0l3EuqAYhIccAPTFTVKQXZl3GfqvLuohnc/+E/OJB9iKd6D+PB6wcy75u5Yc3ZsW83H6XOwR8M8PWq78nNz6NO5YsYdk1fknyJdGjc8tipmuTtyWHNNqYwK8hVKAKMA1ao6ovhG5Jxw6bd27l74rPMSP+Wqxq1ZNzA0TSr2TBs+9+ydwdTU2bjDwaYsyaFfM2nUVxtRnUZQP82nbmyfnO7dd2YAirIEfjVwO3AMhFJcZ4braqfFXxYJlLy8/N5fc4UHv7oXyjKK7f8nuHXJoXlzcoNu7Yeu7FmwQ/LAGhWowGjuw8iyZdIqzpNrGgbE0YFuQplLmD/Gz1k1faN3PXu08xdl0rXZm0Ze9vDBb5BZ/WPm471HVm8aSUAreo04cleQ0nyJYb1qN4YcyK7E7MIyM3L5flZk3hi+puULl6St+74M4Pa9zyvo2FVJW3rumNH2mlb1wHQtkFzxvS9hyRfIo2r1Qn3FIwxJ2EFvJBL2byaIROeYsnmVfRrncBrAx6kRsVzW8NSVVmyaRX+YIDJwa9Zs2MzIkLHxq146eZR9GudQN0q1SM0A2PMqVgBL6QO52Tz5GfjGTPzXeLKVWTy754mqU3ns359fn4+365PO9Z3ZMOubcTGxJLQxMeozgPo2zrhnH8QGGPCywp4ITRvXSpDJjzNqh83Mqh9D17sfx9VylY84+ty83L5Zm0qryyexm0zXmBr5k6Kxxaja7O2PNpjML0vu4a4cpVcmIEx5mxYAS9EDhw+yOhpr/Pq7MnUrVydGSNfolvz9qd9TU5eLl+v+h5/MMBHKXPYeWAPJWKL0SP+avq3SeTGlh2pWDo6d+sZY07PCnghMTN9IUMnPsOmPT8y4tr+PHXTsFM2fDqck83M9IX4g8l8vPQb9h7aT7mSZegZ34EkXyLldudyQ9duLs/AGHOurIB73O6sTB7wv8LbCz7l0ur1+eaBN7i6catffV1W9iE+X76AyUu+5tO0+RzIPkjF0uW46bJOJPkSuL55u2MLNiQnJ7s8C2PM+bAC7mH+JV9zz/vP89OBTEZ3H8SjPQafsGpO5qEDTF82F38wwIzl33IoJ5u4cpUYcMV1JPkS6XzpFZQoVjyKMzDGFIQVcA/anrmLEe8/jz8YwFf3EmaMeInWdS8BQu1gpy0N9R2ZtfI7juTmULNiHIM79CLJl8g1F7eimEsd/owxkWX/kz1EVXnn208Z9eHLHMrJ5pk+w3ngut+w60Amb8yZgj8YILB6CXn5edSrUoN7rk2iv68z7RvGExMTE+3hG2PCzAq4R2zYtZWhE5/lyxWL6Ni4FU/ceBdpW3+gy0sjmLsuFVWlyUV1ebDrbST5Erm8XlPrO2JMIVfQdrLjgRuBHaoaH54hmePl5+fz2uzJPDLtdbKyD9Ggak0O5x7hupdHAhBfqzGP9RhMki+R+FqNrWgbU4QU9Aj8beBVQosbmzDbuG8H1R7qzu6snxc52rBrG1XKVOCp3sNI8iVyaY3IrppjjLlwFaiAq+ocZzk1EyaqyuJNK2k3Zgj5mn/s+asatSTJl0C/1ok0jCtYB0FjTOFQ4DUxnQI+/VSnUGxNzDNTVVbuzmBORhr/XTn7hD+7vXlnejVuR7UyZ74VPlxsrcSikW1z9k5uJNfEbACknc3X2pqYP8vNy9U5q4N63/svat1HeivD2p3w8cacKYVuzhdyts25aGR7NZdIrIlpzk1uXi6z1wTxBwNMTZnN9n27KFmsBGVKlDr2NYM79OL5fiOpXLaC3RFpjDktK+ARlp1zhK+cZlHTUuewKyuTMiVKcUOLq7i+WTvmrktlwsLPaVC1Jv++7RGua9Y22kM2xnhEQS8jfA9IAOJEJAN4XFXHhWNgXnboyGFmpH+LPxjgk6Vz2Xc4i/KlytCrZUeSfIl0b3EVs1cv4e5JY8jYu4P7Ow/gb73vpmzJ0tEeujHGQwp6Fcqt4RqI1+0/nMVnafPxB5P5bPl8srIPUblMBfr5EkjyJdK1aVtKFi/BrgOZDJs0hgkLP6dZjQbM+8NYrmrUMtrDN8Z4kJ1CKYA9Wfv4xGkW9UX6QrJzj3BR+coMbNudJF8CCZdcTnGn74iq8uHirxjx/vPsztrHoz0G86fud1KyeIkoz8IY41VWwM/Rzv17mJb6c7Oo3Pw8aleqxt3X9CHJl8jVjS8jNib2hNds3buTe/77PB+lzubyek2Zee8rtKrTJEozMMYUFlbAz8LWvTuZmjIbfzDA7DVB8jWfhlVrcX/nAST5EmnboPlJm0WpKuPnf8ID/lfIzs3hub4jGNVlgHUDNMaEhVWSU9i4axv+YIDxydNI/2ATqkrTGvV5pNsdJLVJpHWdS07bd+SHnVsYOvEZvlr1PZ2a+Hhz4CM0uaieizMwxhR2VsCPs2bHJvzBZPzBAN9vXAFAo4o1eKLnXfRv05nmNRuecR95+Xn8M/Ahf/r4DWJjYnj91ocY2rGPtXM1xoRdkS7gqkr6tvVMXvI1/pQAy7asA+DK+s15ts9wknyJZKSvIyEh4az2l75tPUMmPMW369PoEd+BN279I3WrVI/gDIwxRVmRK+CqSnDzKvzBAP5gMqt+3IiI0KFRS17sfx/9WidQv2rNY1+fkb7ujPs8kpvDmJkT+Nvnb1G+ZBne/e0T/ObKbtba1RgTUUWigOfn57Nww3KnaAfYsGsbMRJDwiVtuDfxZvq2TqBmxbjz2vf3G1cwZMJTLN2ylgFXdOXlm0dxUYUqYZ6BMcb8WqEt4Hn5ecxdm4o/GGBKSjJb9u6keGwxrmt6JX++4bfc1KoTceUqnff+Dx45zBPT/80Ls96jRoWqTBv2HL1bdQrjDIwx5vQKVQHPycslsGox/mCAj1Jns2P/HkoVL0m35u14ps9werXsSKUy5QucM3v1Eu5692nW7szgdx1v4rm+I8KyX2OMOReeL+CHc7L5csUi/MEAHy+dy56D+yhbsjQ94zuQ5EukR4sOlCtVJixZ+w5l8cepr/LGN1NpFFebr+57lc5Nf92i1xhj3FDQZlbdgZeBWOBNVX02LKM6g6zsQ8xYvgB/MJnpaXPZf/ggFUuXO9YsqlvzdpQ+rkVrOHy6bB7D3hvD1r0/8fsut/Jk77tPaANrjDFuO+8CLiKxwGtAVyAD+E5EPlbV9HAN7nhZOYeZtOgL/MEAny9fwKGcbKqWrcgtbbqQ5EukS9MrKVGseNhzM7OzGPjW40xc9AUtajZi8oNP066hrd9sjIm+ghyBtwXWquoPACLyX+AmIOwFXP63/Umfr1q2IgvWp7FgfVq4I4/Z+NM2juTn8kTPu3ik+6CI/JAwxpjzcd5rYopIf6C7qt7lPL4daKeqI37xdQVeE/PGKY+TlZNNlVLlaRFXjxjcu6sxVmFgfBcaVqrhWuZRXl2/z4vZNueike3V3LCviQn0J3Te++jj24FXT/caWxPTG9k256KRbXP2Ti6nWBOzIIeyW4C6xz2u4zxnjDHGBQUp4N8BTUSkoYiUAAYAH4dnWMYYY87kvN/EVNVcERkBfEHoMsLxqro8bCMzxhhzWgVdE/Mz4LMwjcUYY8w5sCbVxhjjUVbAjTHGo6yAG2OMR1kBN8YYjzrvOzHPK0xkJ7DxPF8eB/wUxuFc6LnRzLY5F41sm7N3cuurarVfPulqAS8IEfleT3YraSHNjWa2zbloZNucvZ9rp1CMMcajrIAbY4xHeamAjy1iudHMtjkXjWybs8dzPXMO3BhjzIm8dARujDHmOFbAjTHGozxRwEWku4isEpG1IvKwS5njRWSHiERuvbaT59YVkYCIpIvIchG5z8XsUiKySERSney/uJXt5MeKSFBEprucu0FElolIioh872JuJRGZLCIrRWSFiFzlUu6lzlyPfuwTkftdyh7lfG+lich7IuLKyuAicp+TuTzScz1Z7RCRKiLypYiscT5XDkvYyVZ5uJA+CLWqXQc0AkoAqUBzF3I7AW2ANJfnWxNo42yXB1a7MV8nT4ByznZxYCHQ3sW5/x6YBEx3+e98AxDnZqaT+w5wl7NdAqgUhTHEAtsJ3SgS6azawHqgtPP4A+BOF3LjgTSgDKEOrLOAiyOY96vaATwHPOxsPwyMCUeWF47Ajy2erKpHgKOLJ0eUqs4Bdkc65yS521R1ibO9H1hB6BvfjWxV1QPOw+LOhyvvcotIHaAn8KYbedEmIhUJ/UcfB6CqR1R1bxSG0gVYp6rne4f0uSoGlBaRYoQK6lYXMpsBC1X1oKrmArOBfpEKO0XtuInQD2ycz33CkeWFAl4b2Hzc4wxcKmjRJiINAB+hI2G3MmNFJAXYAXypqm5lvwQ8BOS7lHc8BWaKyGJnEW43NAR2Am85p43eFJGyLmUfbwDwnhtBqroFeB7YBGwDMlV1pgvRacA1IlJVRMoAPThxOUg3VFfVbc72dqB6OHbqhQJeJIlIOcAP3K+q+9zKVdU8VW1NaI3TtiISH+lMEbkR2KGqiyOddQodVbUNcANwj4h0ciGzGKFfs19XVR+QRehXa9c4SyH2Bj50Ka8yoSPRhkAtoKyIDIx0rqquAMYAM4EZQAqQF+nc04xHCdNvtl4o4EVu8WQRKU6oeE9U1SnRGIPz63wA6O5C3NVAbxHZQOgUWWcRedeFXODYkSGqugOYSui0XaRlABnH/YYzmVBBd9MNwBJV/dGlvOuA9aq6U1VzgClABzeCVXWcql6uqp2APYTeW3LTjyJSE8D5vCMcO/VCAS9SiyeLiBA6L7pCVV90ObuaiFRytksDXYGVkc5V1UdUtY6qNiD07/u1qkb8yAxARMqKSPmj28D1hH7ljihV3Q5sFpFLnae6AOmRzv2FW3Hp9IljE9BeRMo43+ddCL3HE3EicpHzuR6h89+T3Mg9zsfAIGd7EDAtHDst0JqYbtAoLZ4sIu8BCUCciGQAj6vquEjnEjoavR1Y5pyLBhitofVHI60m8I6IxBL64f6Bqrp6SV8UVAemhuoJxYBJqjrDpeyRwETnwOQH4Lcu5R79YdUVuNutTFVdKCKTgSVALhDEvVvb/SJSFcgB7onkG8Ynqx3As8AHIjKEUEvtW8KS5VzWYowxxmO8cArFGGPMSVgBN8YYj7ICbowxHmUF3BhjPMoKuDHGeJQVcOMpTge/4c52LeeytEhltRaRHpHavzEFZQXceE0lYDiAqm5V1f4RzGpNqG+GMRckuw7ceIqIHO1GuQpYAzRT1XgRuZNQh7eyQBNCTZNKELopKhvooaq7RaQx8BpQDTgI/E5VV4rIzYRuuMgDMgnd9r0WKE2odcMzwHTgn4TakxYHnlDVaU52X6AioUZr76qqq73UTdF0wd+JacwvPAzEq2prp1vj8XeKxhPq3liKUPH9o6r6ROQfwB2EOh6OBYap6hoRaQf8C+gMPAZ0U9UtIlJJVY+IyGPAFao6AkBEniZ0m/9gp+XAIhGZ5WS3dfIPAt+JyKeq6triEKZosgJuCpOA00N9v4hkAp84zy8DLnM6PHYAPnRunQco6XyeB7wtIh8QarJ0MtcTarr1B+dxKaCes/2lqu4CEJEpQEfACriJKCvgpjDJPm47/7jH+YS+12OAvU673BOo6jDniLwnsFhELj/J/gVIUtVVJzwZet0vz0XauUkTcfYmpvGa/YSWmjtnTl/19c75biSklbPdWFUXqupjhBZaqHuSrC+AkU4nPUTEd9yfdXXWPSxN6Fz8vPMZozHnwgq48RTnNMU8Z8HYv5/HLm4DhohIKrCcn5fn+7uEFjZOA+YTWns1ADR3Fv79H+BJQm9eLhWR5c7joxYR6uG+FPDb+W/jBrsKxZgCcq5COfZmpzFusSNwY4zxKDsCN8YYj7IjcGOM8Sgr4MYY41FWwI0xxqOsgBtjjEdZATfGGI/6f1B24swQu+IMAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"ax = None\n",
"for i in range(simulation_parameters['N']):\n",
" ax = df[df['run']==i+1].plot('timestep', ['box_A', 'box_B'],\n",
" grid=True,\n",
" xticks=list(df['timestep'].drop_duplicates()), \n",
" yticks=list(range(1+max(df['box_A'].max(),df['box_B'].max()))),\n",
" legend = (ax == None),\n",
" colormap = 'RdYlGn',\n",
" ax = ax\n",
" )"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## A/B testing\n",
"In order to simulate two different versions of the same system, we create two `Configuration` objects to pass to the `Executor` instead of just one. For example, suppose we wanted to test the system under two different sets of initial conditions"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"initial_conditions_1 = {\n",
" 'box_A': 10, \n",
" 'box_B': 0 \n",
"}\n",
"\n",
"initial_conditions_2 = {\n",
" 'box_A': 10, \n",
" 'box_B': 4 \n",
"}"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Just like we did before, we package those initial conditions along with the partial state update blocks and the simulation parameters into `Configuration` objects"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [],
"source": [
"configs.clear()\n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
"# The configurations above are then packaged into a `configs` object\n",
"append_configs(initial_state=initial_conditions_1, #dict containing variable names and initial values\n",
" partial_state_update_blocks=partial_state_update_blocks, #dict containing state update functions\n",
" sim_configs=simulation_parameters #dict containing simulation parameters\n",
" )\n",
"\n",
"# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # \n",
"# The configurations above are then packaged into a `configs` object\n",
"append_configs(initial_state=initial_conditions_2, #dict containing variable names and initial values\n",
" partial_state_update_blocks=partial_state_update_blocks, #dict containing state update functions\n",
" sim_configs=simulation_parameters #dict containing simulation parameters\n",
" )"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"And now we can execute the simulation of those two different versions of the system in parallel. In order to accomplish that, we use an `ExecutionContext` of mode `multi_proc`"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"exec_mode = ExecutionMode()\n",
"exec_context = ExecutionContext(exec_mode.local_mode) # For A/B testing we use multi_proc instead of single_proc\n",
"executor = Executor(exec_context, configs) # Pass the configuration objects inside an array"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"And we are ready to execute the simulation. The `execute()` method will return a list of tuples - the first element of those tuples correspond to the datapoints of each one of the versions of the system being simulated."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"\n",
" ___________ ____\n",
" ________ __ ___/ / ____/ | / __ \\\n",
" / ___/ __` / __ / / / /| | / / / /\n",
"/ /__/ /_/ / /_/ / /___/ ___ |/ /_/ /\n",
"\\___/\\__,_/\\__,_/\\____/_/ |_/_____/\n",
"by cadCAD\n",
"\n",
"Execution Mode: local_proc\n",
"Configuration Count: 2\n",
"Dimensions of the first simulation: (Timesteps, Params, Runs, Vars) = (10, 0, 1, 2)\n",
"Execution Method: local_simulations\n",
"Execution Mode: parallelized\n",
"ValueError: sim_configs' N must > 0\n"
]
},
{
"ename": "TypeError",
"evalue": "'NoneType' object is not iterable",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mTypeError\u001b[0m Traceback (most recent call last)",
"\u001b[0;32m<ipython-input-8-0e2521e50f90>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m 1\u001b[0m \u001b[0;31m# %%capture\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 2\u001b[0m \u001b[0mresults\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0;32mfor\u001b[0m \u001b[0mraw_result\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtensor\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0msessions\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mexecutor\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mexecute\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 4\u001b[0m \u001b[0mdf\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mpd\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mDataFrame\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mraw_result\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 5\u001b[0m \u001b[0mresults\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mappend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdf\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/envs/cadcad/lib/python3.8/site-packages/cadCAD/engine/__init__.py\u001b[0m in \u001b[0;36mexecute\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 138\u001b[0m \u001b[0msim_executors\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mvar_dict_list\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mstates_lists\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mconfigs_structs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0menv_processes_list\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mTs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mSimIDs\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mRunIDs\u001b[0m \u001b[0;31m#Ns\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 139\u001b[0m )\n\u001b[0;32m--> 140\u001b[0;31m \u001b[0mfinal_result\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mget_final_results\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0msimulations_results\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpartial_state_updates\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0meps\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0msessions\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mremote_threshold\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 141\u001b[0m \u001b[0;32melif\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mexec_context\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0mExecutionMode\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdistributed\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 142\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"Execution Method: \"\u001b[0m \u001b[0;34m+\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mexec_method\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m__name__\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;32m~/anaconda3/envs/cadcad/lib/python3.8/site-packages/cadCAD/engine/__init__.py\u001b[0m in \u001b[0;36mget_final_results\u001b[0;34m(simulations, psus, eps, sessions, remote_threshold)\u001b[0m\n\u001b[1;32m 102\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mget_final_results\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0msimulations\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpsus\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0meps\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0msessions\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mremote_threshold\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 103\u001b[0m \u001b[0mflat_timesteps\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mtensor_fields\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 104\u001b[0;31m \u001b[0;32mfor\u001b[0m \u001b[0msim_result\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpsu\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mep\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mlist\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mzip\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0msimulations\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpsus\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0meps\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 105\u001b[0m \u001b[0mflat_timesteps\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mappend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mflatten\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0msim_result\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 106\u001b[0m \u001b[0mtensor_fields\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mappend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mcreate_tensor_field\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mpsu\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mep\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
"\u001b[0;31mTypeError\u001b[0m: 'NoneType' object is not iterable"
]
}
],
"source": [
"# %%capture\n",
"results = []\n",
"for raw_result, tensor, sessions in executor.execute():\n",
" df = pd.DataFrame(raw_result)\n",
" results.append(df)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"for df in results:\n",
" ax = None\n",
" for i in range(simulation_parameters['N']):\n",
" ax = df[df['run']==i+1].plot('timestep', ['box_A', 'box_B'],\n",
" grid=True,\n",
" xticks=list(df['timestep'].drop_duplicates()), \n",
" yticks=list(range(1+max(df['box_A'].max(),df['box_B'].max()))),\n",
" legend = (ax == None),\n",
" colormap = 'RdYlGn',\n",
" ax = ax\n",
" )"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The same thing can be done for any sort of variation in the system being modeled: wheter it's the inclusion of a new state update function, a change in the order in which partial state update blocks get executed, or completely different policy and state update functions, all we need to do is create a different `Configuration` object for each one of the variations and pass those to the `Executor`."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.8.3"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment