{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Simple Example of a Simulation in astronomix" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Imports" ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "# ==== GPU selection ====\n", "from autocvd import autocvd\n", "autocvd(num_gpus = 1)\n", "# ruff: noqa: E402\n", "# =======================\n", "\n", "# general\n", "import jax.numpy as jnp\n", "\n", "# astronomix containers\n", "from astronomix import (\n", " SimulationConfig,\n", " SimulationParams,\n", ")\n", "\n", "# astronomix functions\n", "from astronomix import (\n", " get_helper_data,\n", " get_registered_variables,\n", " construct_primitive_state,\n", " finalize_config,\n", " time_integration\n", ")\n", "\n", "# plotting\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Simulation Setup" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Let us set up a very simple simulation, mostly with default parameters.\n", "\n", "First we get the configuration of the simulation, which contains parameters that typically do not change between simulations, changing which requires (just-in-time)-recompilation." ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "config = SimulationConfig(\n", " box_size = 1.0,\n", " num_cells = 100,\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next we setup the simulation parameters, things we might vary" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "params = SimulationParams(\n", " t_end = 0.2, # the typical value for a shock test\n", ")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "With this we generate some helper data, like the cell centers etc." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "helper_data = get_helper_data(config)\n", "registered_variables = get_registered_variables(config)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next we setup the shock initial conditions, namely\n", "\\begin{equation}\n", "\\left(\\begin{array}{l}\n", "\\rho \\\\\n", "u \\\\\n", "p\n", "\\end{array}\\right)_L=\\left(\\begin{array}{l}\n", "1 \\\\\n", "0 \\\\\n", "1\n", "\\end{array}\\right), \\quad\\left(\\begin{array}{l}\n", "\\rho \\\\\n", "u \\\\\n", "p\n", "\\end{array}\\right)_R=\\left(\\begin{array}{c}\n", "0.125 \\\\\n", "0 \\\\\n", "0.1\n", "\\end{array}\\right)\n", "\\end{equation}\n", "with seperation at $x=0.5$." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "# setup the shock initial fluid state in terms of rho, u, p\n", "shock_pos = 0.5\n", "r = helper_data.geometric_centers\n", "rho = jnp.where(r < shock_pos, 1.0, 0.125)\n", "u = jnp.zeros_like(r)\n", "p = jnp.where(r < shock_pos, 1.0, 0.1)\n", "\n", "# get initial state\n", "initial_state = construct_primitive_state(\n", " config = config,\n", " registered_variables = registered_variables,\n", " density = rho,\n", " velocity_x = u,\n", " gas_pressure = p,\n", ")" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "OPTIMAL_BACKEND: using the PALLAS backend (GPU compute capability >= 8.0).\n", "Setting time integrator to RK4_SSP for finite difference solver mode.\n", "Automatically setting open boundaries for Cartesian geometry.\n" ] } ], "source": [ "config = finalize_config(config, initial_state.shape)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Running the simulation" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "final_state = time_integration(initial_state, config, params, registered_variables)\n", "rho_final = final_state[registered_variables.density_index]\n", "u_final = final_state[registered_variables.velocity_index]\n", "p_final = final_state[registered_variables.pressure_index]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Visualization" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAABL4AAAHDCAYAAAAqZtO0AAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjYsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvq6yFwwAAAAlwSFlzAAAPYQAAD2EBqD+naQAAlRlJREFUeJzs3Xl4VOX5//HPzGSHJCyBBEIgLLLJHgQBrUujqQvWrlStKFa0Kv1a+bVVVKRqK2rVYlsqLS3Ftlqo1K2VohalVkVRFhc2ZTNsCYQtECAhM+f3x8lMEgiQZWaeOWfer+uansPkzMw96Zhnzn2e+348lmVZAgAAAAAAAFzGazoAAAAAAAAAIBJIfAEAAAAAAMCVSHwBAAAAAADAlUh8AQAAAAAAwJVIfAEAAAAAAMCVSHwBAAAAAADAlUh8AQAAAAAAwJVIfAEAAAAAAMCVSHwBAAAAAADAlUh8AYZ5PB799Kc/NR0GAMCQSI4DW7Zskcfj0dy5cyPy/AAAALGOxBdca+7cufJ4PKFbSkqKOnfurKKiIv3qV7/SwYMHTYfYoHfffVc//elPtX//ftOhAACOc8UVVygtLe2UY8g111yjpKQk7dmzJ4qRNd7ChQu54AIAMaih85fevXtr0qRJKi0tNR0e4FgkvuB6DzzwgP7yl7/oqaee0g9+8ANJ0g9/+EMNHDhQH3/8seHopCNHjujee+8N/fvdd9/V/fffT+ILAGLQNddcoyNHjuiFF15o8OeHDx/WSy+9pK985Stq3759lKM7Ubdu3XTkyBFde+21ofsWLlyo+++/32BUAIBTCZ6//OY3v9Ho0aP11FNPadSoUTp8+LDp0ABHSjAdABBpl1xyiYYPHx7695QpU/TGG2/o8ssv1xVXXKG1a9cqNTXVWHwpKSnGXhsA0DRXXHGF0tPT9eyzz2r8+PEn/Pyll15SRUWFrrnmGgPRnSg4YwAA4Bx1z19uvPFGtW/fXk888YReeuklXXXVVSccX1FRoVatWkU7zGazLEtHjx41eg6G+MKML8SlCy+8UFOnTtUXX3yhv/71r6H7161bp29+85tq166dUlJSNHz4cL388sv1HhucgvzOO+9o8uTJ6tChg1q1aqWvfe1r2r17d71jP/zwQxUVFSkrK0upqanq3r27brjhhnrH1O3t8tOf/lQ//vGPJUndu3cPTXPesmWLzjvvPA0ePLjB99OnTx8VFRW19NcCADiN1NRUff3rX9fixYu1a9euE37+7LPPKj09XVdccYX279+vH/7wh8rLy1NycrJ69eqlRx55RIFA4LSvs3LlSl1yySXKyMhQ69at9eUvf1nvvffeCcft379fd9xxh/Lz85WcnKwuXbpo/PjxKisrk3Rij6/rr79eM2fOlKR65TSWZSk/P19f/epXT3iNo0ePKjMzUzfffHNTflUAgDC58MILJUmbN2/W9ddfr9atW2vjxo269NJLlZ6eHrrYEggENGPGDJ155plKSUlRdna2br75Zu3bt6/e8zXmHGXevHkqKChQenq6MjIyNHDgQD355JOhn//0pz+Vx+M5IdbgudKWLVtC9+Xn5+vyyy/Xq6++quHDhys1NVW/+93vJKlFYyXQWMz4Qty69tprdffdd+u1117TxIkTtXr1ao0ZM0a5ubm666671KpVK/3973/XlVdeqX/84x/62te+Vu/xP/jBD9S2bVtNmzZNW7Zs0YwZMzRp0iTNnz9fkrRr1y5dfPHF6tChg+666y61adNGW7Zs0fPPP3/SmL7+9a/rs88+09/+9jf98pe/VFZWliSpQ4cOuvbaazVx4kR9+umnGjBgQOgxH3zwgT777LN65ZIAgMi55ppr9PTTT+vvf/+7Jk2aFLp/7969evXVV3XVVVfJsiydd9552r59u26++WZ17dpV7777rqZMmaKdO3dqxowZJ33+1atX69xzz1VGRoZ+8pOfKDExUb/73e90/vnn67///a9GjhwpSTp06JDOPfdcrV27VjfccIOGDRumsrIyvfzyy9q2bVtoDKnr5ptv1o4dO/T666/rL3/5S+h+j8ej7373u3r00Ue1d+9etWvXLvSzf/7znyovL9d3v/vdMPz2AABNtXHjRkkKldBXV1erqKhI55xzjh577DGlpaVJsv/Gz507VxMmTND//d//afPmzfrNb36jlStX6p133lFiYmKjzlFef/11XXXVVfryl7+sRx55RJK0du1avfPOO7r99tub9R7Wr1+vq666SjfffLMmTpyoPn366PDhw80eK4EmsQCX+tOf/mRJsj744IOTHpOZmWkNHTrUsizL+vKXv2wNHDjQOnr0aOjngUDAGj16tHXGGWec8LyFhYVWIBAI3X/HHXdYPp/P2r9/v2VZlvXCCy+c9vUty7IkWdOmTQv9+xe/+IUlydq8eXO94/bv32+lpKRYd955Z737/+///s9q1aqVdejQoVO+DgAgPKqrq61OnTpZo0aNqnf/rFmzLEnWq6++aj344INWq1atrM8++6zeMXfddZfl8/ms4uLi0H3HjwNXXnmllZSUZG3cuDF0344dO6z09HTrS1/6Uui+++67z5JkPf/88yfEGByfNm/ebEmy/vSnP4V+dtttt1kNfQVcv369Jcl66qmn6t1/xRVXWPn5+fXGPABA+AXPM/7zn/9Yu3fvtrZu3WrNmzfPat++vZWammpt27bNuu666yxJ1l133VXvsf/73/8sSdYzzzxT7/5FixbVu78x5yi33367lZGRYVVXV5/0mGnTpjU4lgTfQ91zmW7dulmSrEWLFtU7tiljJdASlDoirrVu3VoHDx7U3r179cYbb+jb3/62Dh48qLKyMpWVlWnPnj0qKirS559/ru3bt9d77E033VRveu+5554rv9+vL774QpLUpk0bSdK//vUvHTt2rMWxZmZm6qtf/ar+9re/ybIsSZLf79f8+fN15ZVXOqquHwCczOfz6Tvf+Y6WLl1ar5Tj2WefVXZ2tr785S/rueee07nnnqu2bduGxpSysjIVFhbK7/frrbfeavC5/X6/XnvtNV155ZXq0aNH6P5OnTrp6quv1ttvv63y8nJJ0j/+8Q8NHjz4hBnJkhosPzmd3r17a+TIkXrmmWdC9+3du1f//ve/dc011zTrOQEATVdYWKgOHTooLy9P3/nOd9S6dWu98MILys3NDR1zyy231HvMc889p8zMTF100UX1xp2CggK1bt1ab775pqTGnaO0adNGFRUVev3118P2nrp3735Ca5bmjpVAU5H4Qlw7dOiQ0tPTtWHDBlmWpalTp6pDhw71btOmTZOkE3q5dO3atd6/27ZtK0mhGvrzzjtP3/jGN3T//fcrKytLX/3qV/WnP/1JlZWVzY53/PjxKi4u1v/+9z9J0n/+8x+VlpbWW60LABB5wX4qzz77rCRp27Zt+t///qfvfOc78vl8+vzzz7Vo0aITxpTCwkJJJ44pQbt379bhw4fVp0+fE37Wr18/BQIBbd26VZJd+lK39D0cxo8fr3feeSd0Eee5557TsWPHGGcAIIpmzpyp119/XW+++abWrFmjTZs21UsaJSQkqEuXLvUe8/nnn+vAgQPq2LHjCWPPoUOHQuNOY85Rbr31VvXu3VuXXHKJunTpohtuuEGLFi1q0Xvq3r37Cfc1d6wEmooeX4hb27Zt04EDB9SrV69Q88Qf/ehHJ20S36tXr3r/9vl8DR4XnI3l8Xi0YMECvffee/rnP/+pV199VTfccIMef/xxvffee2rdunWTYy4qKlJ2drb++te/6ktf+pL++te/KicnJzQ4AACio6CgQH379tXf/vY33X333aHZuHUbDF900UX6yU9+0uDje/fuHc1wG+073/mO7rjjDj3zzDO6++679de//lXDhw9vMBEHAIiMESNG1FuV/njJycnyeuvPYQkEAurYsWO9Wbt1dejQQVLjzlE6duyoVatW6dVXX9W///1v/fvf/9af/vQnjR8/Xk8//XToeRri9/sbvL+hFRydOlbCeUh8IW4Fm/oWFRWFykkSExPDnkQ6++yzdfbZZ+vnP/+5nn32WV1zzTWaN2+ebrzxxgaPP1Upic/n09VXX625c+fqkUce0YsvvqiJEyeeNAkHAIica665RlOnTtXHH3+sZ599VmeccYbOOussSVLPnj116NChJo8pHTp0UFpamtavX3/Cz9atWyev16u8vLzQa3z66adNjvtU40y7du102WWX6ZlnntE111yjd955h+bCAOAAPXv21H/+8x+NGTOmwSTT8U53jpKUlKSxY8dq7NixCgQCuvXWW/W73/1OU6dOVa9evULVLvv37w+VT0oKzRhubMzNGSuBpqLUEXHpjTfe0IMPPqju3bvrmmuuUceOHXX++efrd7/7nXbu3HnC8bt3727ya+zbty80+ytoyJAhknTKcsdgr679+/c3+PNrr71W+/bt080336xDhw6xyhYAGBKc3XXfffdp1apVoX9L0re//W0tXbpUr7766gmP279/v6qrqxt8Tp/Pp4svvlgvvfRSvf5hpaWlevbZZ3XOOecoIyNDkvSNb3xDH330kV544YUTnuf48aeuxowza9as0Y9//ONQPzMAQGz79re/Lb/frwcffPCEn1VXV4f+5jfmHGXPnj31fu71ejVo0KB6x/Ts2VOS6vXhqqioCM0Ia2zMzRkrgaZixhdc79///rfWrVun6upqlZaW6o033tDrr7+ubt266eWXX1ZKSooku5b+nHPO0cCBAzVx4kT16NFDpaWlWrp0qbZt26aPPvqoSa/79NNP67e//a2+9rWvqWfPnjp48KBmz56tjIwMXXrppSd9XEFBgSTpnnvu0Xe+8x0lJiZq7NixoROVoUOHasCAAXruuefUr18/DRs2rJm/GQBAS3Tv3l2jR4/WSy+9JEn1El8//vGP9fLLL+vyyy/X9ddfr4KCAlVUVOiTTz7RggULtGXLFmVlZTX4vD/72c/0+uuv65xzztGtt96qhIQE/e53v1NlZaUeffTReq+xYMECfetb39INN9yggoIC7d27Vy+//LJmzZqlwYMHN/j8wXHm//7v/1RUVHRCcuuyyy5T+/bt9dxzz+mSSy5Rx44dW/y7AgBE1nnnnaebb75Z06dP16pVq3TxxRcrMTFRn3/+uZ577jk9+eST+uY3v9moc5Qbb7xRe/fu1YUXXqguXbroiy++0K9//WsNGTJE/fr1kyRdfPHF6tq1q773ve+FLpTMmTNHHTp0UHFxcaNibslYCTSJwRUlgYgKLqUbvCUlJVk5OTnWRRddZD355JNWeXn5CY/ZuHGjNX78eCsnJ8dKTEy0cnNzrcsvv9xasGDBCc97/BLAb775piXJevPNNy3LsqwVK1ZYV111ldW1a1crOTnZ6tixo3X55ZdbH374Yb3H6bhl7C3LXto3NzfX8nq9JywHbFmW9eijj1qSrIceeqj5vyAAQIvNnDnTkmSNGDHihJ8dPHjQmjJlitWrVy8rKSnJysrKskaPHm099thjVlVVVei4hsaBFStWWEVFRVbr1q2ttLQ064ILLrDefffdE15jz5491qRJk6zc3FwrKSnJ6tKli3XddddZZWVllmVZ1ubNmy1J1p/+9KfQY6qrq60f/OAHVocOHSyPx9PgcvS33nqrJcl69tlnm/mbAQA01cnOM+q67rrrrFatWp3057///e+tgoICKzU11UpPT7cGDhxo/eQnP7F27NhhWVbjzlEWLFhgXXzxxVbHjh2tpKQkq2vXrtbNN99s7dy5s95rLV++3Bo5cmTomCeeeCL0Huqev3Tr1s267LLLGoy3sWMl0BIeyzrFXHgAMenJJ5/UHXfcoS1btpywuiQAAC11xx136I9//KNKSkqUlpZmOhwAAIBmI/EFOIxlWRo8eLDat2+vN99803Q4AACXOXr0qPLy8nT55ZfrT3/6k+lwAAAAWoQeX4BDVFRU6OWXX9abb76pTz75JNRTBgCAcNi1a5f+85//aMGCBdqzZ49uv/120yEBAAC0GIkvwCF2796tq6++Wm3atNHdd9+tK664wnRIAAAXWbNmTWil41/96lehVb4AAACcjFJHAAAAAAAAuJLXdAAAAAAAAABAJJD4AgAAAAAAgCs5osdXIBDQjh07lJ6eLo/HYzocAHA8y7J08OBBde7cWV4v10AYZwAgvBhn6mOcAYDwaso444jE144dO5SXl2c6DABwna1bt6pLly6mwzCOcQYAIoNxxsY4AwCR0ZhxxhGJr/T0dEn2G8rIyDAcDQA4X3l5ufLy8kJ/X+Md4wwAhBfjTH2MMwAQXk0ZZxyR+ApOB87IyGCgAIAwotzCxjgDAJHBOGNjnAGAyGjMOEPBPQAAAAAAAFyJxBcAAAAAAABcicQXAAAAAAAAXMkRPb4AuJvf79exY8dMh+EqiYmJ8vl8psNwHT6r4ZeUlHTaJagBAABiHd8Twy9c3xNJfAEwxrIslZSUaP/+/aZDcaU2bdooJyeHxsJhwGc1crxer7p3766kpCTToQAAADQZ3xMjJ1zfE0l8ATAmOEB07NhRaWlpJGjCxLIsHT58WLt27ZIkderUyXBEzsdnNTICgYB27NihnTt3qmvXrvxeAQCA4/A9MTLC+T2RxBcAI/x+f2iAaN++velwXCc1NVWStGvXLnXs2JGyxxbgsxpZHTp00I4dO1RdXa3ExETT4QAAADQa3xMjK1zfE2mqAcCIYP17Wlqa4UjcK/i7pddAy/BZjazg1HW/3284EgAAgKbhe2Jkhet7IokvAEYxFThy+N2GF7/PyOD3CgAAnI7vM5ERrt8riS8AAAAAAAC4UpMTX2+99ZbGjh2rzp07y+Px6MUXXzztY5YsWaJhw4YpOTlZvXr10ty5c5sRKgDEhvPPP18//OEPG3Xsli1b5PF4tGrVqrA9p2T/XfV4PKweg1PiswoAJ+J8BgDi63tik5vbV1RUaPDgwbrhhhv09a9//bTHb968WZdddpm+//3v65lnntHixYt14403qlOnTioqKmpW0ABg0vPPP9/o5op5eXnauXOnsrKyJNl/3C+44ALt27dPbdq0adZzAo3FZxUATsT5DADE1/fEJie+LrnkEl1yySWNPn7WrFnq3r27Hn/8cUlSv3799Pbbb+uXv/wlAwUAR2rXrl2jj/X5fMrJyQnrcwKNxWcVAE7E+QwAxNf3xCYnvppq6dKlKiwsrHdfUVFRk6a/tcRnK/6rqor9UXmtqPF4Zckjy5tg73sSFPAly5+QqkBCqqqTMmT5kk1HGXe6Z7VS5zappsNAFJx//vkaMmSIZsyYofz8fN10003asGGDnnvuObVt21b33nuvbrrpJkn2tODu3btr5cqVatOmjS644AJJUtu2bSVJ1113nebOnVvvOSXpL3/5i5588kmtX79erVq10oUXXqgZM2aoY8eORt4znInPqssFAtLR/dKRfZJlSV6v5PFJXl/tVh5JVu1jPD7J46k5xmvf5LHvO+lWx+2r/n5LWFbD+/YdJ/lZQ/dbNfsn2VqBOvfJ/rcVqP1Z3Zus2p8FX8PySwG/vfUfkwLVkr/Kvi/4vAnJUkobKbWNlJYlJSS18JeDWGH6fOZQ+T5t+eitqLxW1Hi8sjxeWR6fLG+iAr4kWd5kVSelcy5jkM/r0bCubZWUQCvweBBP3xMjnvgqKSlRdnZ2vfuys7NVXl6uI0eOKDX1xERBZWWlKisrQ/8uLy9v9utbC3+iAdXrmv14p9prtVap1U6brBwtDZypdwNnapPVSfYXYERCcoJXy+4pVGZq7E3tdArLsnTkWMuWqm2u1ERfs1cNefzxx/Xggw/q7rvv1oIFC3TLLbfovPPOU58+feodl5eXp3/84x/6xje+ofXr1ysjI6PBv4GSvTTygw8+qD59+mjXrl2aPHmyrr/+ei1cuLBZMSK8TH1WW/I5lfisOp6/Wtr8X+nTf0gb/iNV7K5J1CDmJKVLQ6+RRtwkte9pOhq0kOnzmZLNazRg8fhmP96JjlhJ2mx10mdWrj4L5Om1QIE2WF1MhxUXvj28ix795mDTYTiWU89nJHd/T4x44qs5pk+frvvvvz8sz3UwrYs2HzoalueKDZY8NTefAvLKrwSrWomqVrJVqRQdlVeW2nkOqZ3nkPqpWJf5lkmSdng7aVbqjVqWOMLwe3Cfz0oPqrI6oLJDlSS+WuDIMb/63/eqkdde80CR0pKa9yfx0ksv1a233ipJuvPOO/XLX/5Sb7755gmDhM/nC03/7dixY716+OPdcMMNof0ePXroV7/6lc466ywdOnRIrVu3blacCB9Tn9WWfE4lPquO9skCadFddrLreEmt7RlcgTqzkoJbnF7d2XHBWXDHz3yrO5vOmyj5am7BGXTySNVHpKMHpCP7paqD0vuzpPd/J/UbK135lJTMfw/xJJznM4nJqdrkzQ/92yNLVpQvZnt0/IzMhjUmLq8Cdc5n/PJZfiWqWklWpVrpsLyylOqpUn/PF+qvLySf9BPN10rfIL2cfLneSxihgMfX0reE45QfPaadB45q694jpkNxNKeez0ju/p4Y8cRXTk6OSktL691XWlp6yqzglClTNHny5NC/y8vLlZeX16zXHz75H816nGNZll3ycHCnVL5D2rHKvjq89X119u/UAxUPSmd+TfrKI1J69mmfDo0z9IHXtO/wMQUCjftSAHcZNGhQaN/j8SgnJ0e7du1q0XMuX75cP/3pT/XRRx9p3759CgTsWR3FxcXq379/i54b8YvPqkN9NE964fuSLCm1ndT/q9KAr0tZve1/n6ykzrLqlyRadcr+LH/9Mr+TlgrWeWztP058nXo/O8mJ70mvQp/iRPn4x9QtvWzovhP2a7bHl3R6faeIp4UCAWnzEum9p6TPX5PWvix1Gix96UeReT1EnOnzmW59h0n3fdSsxzpOIGAnjivKpN3rpd1rpa0fSJ+/qqH+jzX08MdS11HSNc9Jyemmo3WVVz7eqdueXSH/CWXniBdu/p4Y8cTXqFGjTpjG9vrrr2vUqFEnfUxycrKSk6nrbhaPR0prZ9+yz5TOuEg678dS5SHpv49IS2dKq1+QNr4pXfuClDvMdMSu4K358kzeq2VSE31a84CZJrGpic2/cnj8yiUejyf0R705KioqVFRUpKKiIj3zzDPq0KGDiouLVVRUpKqqqmY/L8LH1Ge1JZ9Tic+qI9VNeg2/QbrkUXumUWM0lDTyeCR5FaOT/t3B65V6XmjfVvxFenmStHyudM4dNTPL4DScz0SR1yulZNq39j2lvpfa9+/fKn04R/rgD1LxUunZcXbyK6mV2XhdxFszZFgkvlrEqeczkru/Jzb5W8+hQ4e0YcOG0L83b96sVatWqV27durataumTJmi7du3689//rMk6fvf/75+85vf6Cc/+YluuOEGvfHGG/r73/+uV155JXzvAqeX3Fq6+EFp4DellyZJJR9Lz10n3fw/uwkrWsTrDSa+GChawuPxtGh6rhMkJdkzM/z+k5cgrVu3Tnv27NHDDz8cujr84YcfRiU+NA6fVRuf1Qj79B/1k16XPm6fFMI5Bn5Len2qdGCrPfurT+NXEkTkcD7jQG3ypMJpUr/LpT9fKX3xjvS3q6Sr50uJLC4VDh4u5IdFPHxHlJz3PbHJ354+/PBDDR06VEOHDpUkTZ48WUOHDtV9990nSdq5c6eKi4tDx3fv3l2vvPKKXn/9dQ0ePFiPP/64/vCHP7D0rymdBkvX/VNq003aXyy9dFsDqyihqYJXSPyMFDiNbt26yePx6F//+pd2796tQ4cOnXBM165dlZSUpF//+tfatGmTXn75ZT344IMGokU847NqWFWFtPDHkiypYAJJL6dKTJGGftfe/+APZmNBCOczDpZbIF2zwO5tuPm/0oLvmY7INXxcyEcTOO17YpO/QZ1//vmyLOuE29y5cyVJc+fO1ZIlS054zMqVK1VZWamNGzfq+uuvD0PoaLbUNtK3n5Z8SdK6f9mNV9EiwVJHxgmcTm5uru6//37dddddys7O1qRJk044pkOHDpo7d66ee+459e/fXw8//LAee+wxA9EinvFZNWzFn6XDe6S23aVLHyPp5WQFE+zthsXS3s1mY4Ekzmccr+tI6eq/2+cy61+RvlhqOiJXCF7Ip2cxGsNp3xM9lgOKeMvLy5WZmakDBw4oIyPDdDju8f7vpX//2F6d6MbXpc5DTUfkWGMefkPb9x/RS7eN0eC8NqbDcYSjR49q8+bN6t69u1JSUkyH40qn+h3zd7W+U/0++KxGFr/fBlRXSb8aIpVvly6fIQ2fYDoitNRfvi5tXCyNuV266AHT0UQF40x9/D4i4OX/k1Y8LfX+il3yiBZ5c90uTZj7gQbmZuqfPzjHdDiOwfeYyArX+QyXD+PZiIlSvyukwDHp9Wmmo3G04IV4pgYDAFrs4/l20qt1jjTkatPRIBzOqinHWvEX6dhRs7EAbjHmdkke6bNFUulq09E4Hj2L4WYkvuKZxyMVPWTP+Nr8X6n4fdMROVbtqo4MFACAFgj4pbd/ae+P/oGUwKpwrnBGkZSRKx3Za7eZANBy7XtK/b9q77/zpNlYXICexXAzEl/xrk2eNOQqe/+tR83G4mBeVkEBAITDmpekvRul1LZSwfWmo0G4+BKk3jWN0HevNxsL4Cbn/NDefrJA2veF0VCcjp7FcDMSX5DOmSx5fNKG/0jblpuOxpFoBgkACIvgrIWR35eSW5uNBeGV1t7eHtlnNg7ATToPlXpcIFl+aelvTEfjaFSwwM1IfEFq110aNM7ef+sXZmNxqOBA4WegAAA0156N0s5VkjdBOmui6WgQbqlt7S2JLyC8grO+VvxFOrLfZCSOFip15HwGLkTiC7Zz/5/k8Uqf/Vva+ZHpaByHqcEAgBZb+097m3+O1Kq92VgQfilt7O3R/SajANyn+3lSVm+p+ohdwYJmCTa353wGbkTiC7asXtKAb9j7b88wGooTeYKljowUQFjMnDlT+fn5SklJ0ciRI7Vs2bJTHj9jxgz16dNHqampysvL0x133KGjR1k5DQ4TbHre93KzcSAymPEFRIbHI/W9zN5f94rZWBzMy/kMXIzEF2qN/oG9XfcvqWKP2VgcxlczUrAKCtBy8+fP1+TJkzVt2jStWLFCgwcPVlFRkXbt2tXg8c8++6zuuusuTZs2TWvXrtUf//hHzZ8/X3fffXeUIwdaoHyntO0Dez94Agd3IfEFRE6fmr+bG/4jVVeZjcWh6PEFNyPxhVqdBkudhkj+KunjeaajcRRKHYHweeKJJzRx4kRNmDBB/fv316xZs5SWlqY5c+Y0ePy7776rMWPG6Oqrr1Z+fr4uvvhiXXXVVaedJQbElPU1sxRyh0sZnc3GgsgIJb72Gw0DcKXcAql1tlRZLn3xtuloHCmU+AoYDgSIABJfqK/gOnu7/GmyOE3A1OD4YlmWbrrpJrVr104ej0dt2rTRD3/4w7C+xk9/+lMNGTIkrM/pBFVVVVq+fLkKCwtD93m9XhUWFmrp0qUNPmb06NFavnx5KNG1adMmLVy4UJdeeulJX6eyslLl5eX1bm7EZ9VB1taUOfajzNG1UtvY26P7ObMEws3rlXp/xd5ft9BsLA7FjK/4E0/fExNMB4AYM+Cb0qv3SGXrpa3vS13PNh2RIwSbQVLpGB8WLVqkuXPnasmSJerRo4e8Xq9SU1NNh+UKZWVl8vv9ys7Ornd/dna21q1b1+Bjrr76apWVlemcc86RZVmqrq7W97///VOWOk6fPl33339/WGOPRXxWHeLIPmnL/+z9vmPNxoLICTa3twL2rJRgIgxAePS5VFrxtLT+39Klv6htwotG8dZMiSHxFT/i6XsiM75QX0qGdObX7f0VfzYbi4MEr5DQ4ys+bNy4UZ06ddLo0aOVk5Ojjh07Kj093XRYcWvJkiV66KGH9Nvf/lYrVqzQ888/r1deeUUPPvjgSR8zZcoUHThwIHTbunVrFCOOHj6rDvH561KgWurQ115sBu6UmCIlptn79PkCwq/HefZ/Y+XbWKW+GWrPZwwHgqiJp++JJL5womC546fPS0cPmI3FIYKljhZXSFzv+uuv1w9+8AMVFxfL4/EoPz9f559/fr1pwfn5+XrooYd0ww03KD09XV27dtXvf//7es9z5513qnfv3kpLS1OPHj00depUHTt2LMrvJvZkZWXJ5/OptLS03v2lpaXKyclp8DFTp07VtddeqxtvvFEDBw7U1772NT300EOaPn26AicpJ0pOTlZGRka9m9vwWXWQtf+0t6zm6H7BWV9H95uMAnCnxFSp54X2/nrKHZuqtmcx5zPxIN6+J5L4wom6nCV16CdVH5E+ec50NI5QWxNvOBCnsyypqsLMrZGD/JNPPqkHHnhAXbp00c6dO/XBBx80eNzjjz+u4cOHa+XKlbr11lt1yy23aP369aGfp6ena+7cuVqzZo2efPJJzZ49W7/85S/D8mt0sqSkJBUUFGjx4sWh+wKBgBYvXqxRo0Y1+JjDhw/L660/nPl8PkkR/PJm6rPahPfDZ9Uhjh21VyGT6O8VD1jZEYis4Kq49PlqMh+ljuHhgPMZKf6+J9LjCyfyeKRh46VXp0gr/yqddaPpiGJeaGowA0XLHDssPWRoNbO7d0hJrU57WGZmptLT0+Xz+U46A0mSLr30Ut16662S7Cshv/zlL/Xmm2+qT58+kqR77703dGx+fr5+9KMfad68efrJT37SwjfifJMnT9Z1112n4cOHa8SIEZoxY4YqKio0YcIESdL48eOVm5ur6dOnS5LGjh2rJ554QkOHDtXIkSO1YcMGTZ06VWPHjg0lwMLO1Ge1kZ9Tic+qY2x93/48tc62V1aGu5H4AiLrjCLJ45VKP5H2F0ttupqOyDE8tG4JDwecz0jx9z2RxBcaNmic9PpUacdKqWwDPUdOIzjZhKnBCBo0aFBo3+PxKCcnR7t27QrdN3/+fP3qV7/Sxo0bdejQIVVXV7uy3K45xo0bp927d+u+++5TSUmJhgwZokWLFoUa3hcXF9eb4XXvvffK4/Ho3nvv1fbt29WhQweNHTtWP//5z029BUfhs2rY5v/a2+7n0Yg5HgQb2h/ZbzIKwL1atberV7a+L236rzTsWtMROUZtqaPhQBBT3PI9kcQXGtaqvV0j//lrdrnjBVNMRxTTWP43TBLT7CsVpl47nE+XmFjv3x6PJ9RvaunSpbrmmmt0//33q6ioSJmZmZo3b54ef/zxsMbgZJMmTdKkSZMa/NmSJUvq/TshIUHTpk3TtGnTohBZDVOf1TB/TiU+q8Ztqkl89TjPbByIjlDiixlfQMR0/5Kd+NryNomvJgj2LOZ8poVcdD4jued7IokvnNzAb9Umvs6/iyvRpxBKfLEKSst4PI2enutk7777rrp166Z77rkndN8XX3xhMCI0GZ9VhMPRA9KOFfZ+dxJfcYFSRyDy8s+R3vqFnfiyLM5hGomexWESJ98RJWd9TyTxhZPrc6mUkCrt3WiXPOYOMx1RzApeIaHHFxrjjDPOUHFxsebNm6ezzjpLr7zyil544QXTYQEn4LMaYVvekayA1K6H1CbPdDSIhlDia7/RMABX6zJC8iZK5dukfZvtv7E4La+XnsVoGid9T2RVR5xccmup76X2/icLzMYS41j+F01xxRVX6I477tCkSZM0ZMgQvfvuu5o6darpsIAT8FmNsLr9vRAfUtrY26P7TUYBuFtSmtRluL2/5W2zsThI8EI+5zNoLCd9T/RYDvhkl5eXKzMzUwcOHIjJRmmutv7f0t++I7XOkSavkbwRWiHN4Sb++UO9vqZU078+UFeNYPWYxjh69Kg2b96s7t27KyUlxXQ4rnSq3zF/V+s71e+Dz2pkxfXvd+bZ0u610rfmSmd+zXQ0iIZPn5cWTJC6jZEmLDQdTUQxztTH7yPK3viZXe44aJz09d+bjsYRdpUf1YiHFsvn9WjjQ5eaDscx4vp7TBSE63yGGV84tZ5ftqflHyqRtvzPdDQxK1TqSFE8AKAxDpbaSS9Jyv+S2VgQPfT4AqIj/xx7G+zzhdPy1FSwcD4DNyLxhVNLSJL6X2nvf/Kc0VBiGaWOAIAm2fyWvc0ZaK+kjPhAjy8gOkJ9vrbbfb5wWt46awBwTgO3IfGF0xv4LXu75p9SdaXZWGJUsBkkF0gAAI2yeYm9pb9XfEltY2+Z8QVEFn2+msxbZ/VLzmngNiS+cHpdR0npnaTKA9LGN01HE5O8TA0GADSWZUmbamZ89TjfaCiIsuCMr+oj0rEjZmMB3C7/XHtL4qtRvN66iS/OaeAuJL5wel6v1P+r9v7q2Fye1LTgOMEgAQA4rX1bpAPFkjfBvriE+JGcIXlqFgqi3BGILPp8NUndUkcu5sNtSHyhcYKrTa1fSLljA2p7fBkOxIECgYDpEFyL32148fuMjLjsI/LFO/Y2t0BKbm02FkSXxyOlZNr7R/cbDQVwvS5nSb4k+nw1Ut1Sx3gcmluK74mREa7viQlheRa4X5cRUnpn6eAOaeMbUp9LTEcUU4IDBTO+Gi8pKUler1c7duxQhw4dlJSUFFpNBi1jWZaqqqq0e/dueb1eJSUlmQ7J0fisRo5lWdq9e7c8Ho8SExNNhxM9Xyy1t91Gm40DZqS2lY7spc8XEGlJaVLucKn4XWnLO1K7HqYjimk+Sh2bhe+JkRPO74kkvtA4Xq905pXSe7+1yx1JfNUTHCf8DBKN5vV61b17d+3cuVM7duwwHY4rpaWlqWvXrvJ6mdzbEnxWI8vj8ahLly7y+XymQ4me4IyvriS+4lJoZUcSX0DEdT3bTnxtfU8adq3paGJa3VwN5zSNx/fEyArX90QSX2i8/lfaia91C6VjR6XEFNMRxQxKHZsnKSlJXbt2VXV1tfx+v+lwXMXn8ykhIYErTmHCZzVyEhMT4yvpdbCkpuTGI3UdaToamBBKfO03GgYQF/Jq/s5uXWY2DgeoV+pI1V6T8D0xcsL1PZHEFxqvy1lSRq5dJ79xsdT3MtMRxYzgKigBGkE2WXDqalyVOcGR+KwiLL54197mDKjt9YT4ktrG3jLjC4i8vBH2tuwz6fBeKa2d2XhimM9DqWNL8D0xtlH/gsbzeu1ZX5K0+kWTkcQcSh0BAI0STHx1G2M2DphDqSMQPWntpKze9j6zvk6JUke4GYkvNE3d1R2PHTUbSwypbW5vOBAAQGwrrmls33WU2ThgDokvILqCs762vm82jhjn8XhCyS9mfMFtSHyhaboMlzK6SFWH7NUdIal2FZRwLbcKAHChI/uk0tX2Pis6xq+UNvb26H6TUQDxgz5fjUbfYrgViS80jccj9b/C3l/zotFQYknw6oifKV8AgJMpfl+SJbXvJbXuaDoamMKMLyC68s62t9uXS/5jZmOJcb5QFQvnNHAXEl9ouv5ftbfr/y1VV5qNJUZQ6ggAOK3imv5elDnGNxJfQHS172X/d1d9RCr52HQ0MY2L+XArEl9oui4jpNY5UmW5tGmJ6WhiQrC5PaWOAICTorE9pDqJr/1GwwDihtdrn79IlDueBqWOcCsSX2g6r7dOueNLZmOJEV4v04IBAKdQdVjasdLe78aMr7iW2sbeMuMLiB4a3DeKj3MauBSJLzRPsNxx3b+k6iqzscSA4NURf8BwIACA2LT9QylQLaV3ltp0Mx0NTArO+Dp6QAr4zcYCxIuuNX2+mPF1SpQ6wq1IfKF5uo6SWnW0v7Rtect0NMZ5WfoXAHAqxe/Z226jas8sEJ+CqzrKsr9HAYi8zsMkj08q3y7t32o6mphF32K4FYkvNI/XJ/Uba+9T7hhaAYUeXwCABgXLa4KriyF+JSRJia3s/aP7jYYCxI2kNKnTIHufcseTom8x3IrEF5ovWO649l+Sv9psLIZ5gqWODBIAgOMFAtLWD+z9YJ8ZxDdWdgSiL2+kvaXc8aRqe3wZDgQIMxJfaL5uY6S09tKRvdKW/5mOxiimBQMATmr3WqnygD3LJ3uA6WgQC0h8AdHX5Sx7u+0Ds3HEsNDFfE5q4DIkvtB8vgSp72X2/tqXzcZimK/mvySmBQMAThAsq+lSYI+dQGhlx/0mowDiSzDxVfKxdOyI2VhiFH2L4VYkvtAydcsd43hlIq6OAABOqpj+XjhOKPHFjC8gatp0lVpn2yvs7vzIdDQxqbZvseFAgDAj8YWW6X6evTpRxa7aFaviEKWOAICTCs746jrSbByIHaFSx/1GwwDiisdDueNp0LcYbkXiCy3jS5T6XGrvx/HqjkwLBgA06NAuad9mSXVOuAB6fAFmdBlub0l8Nchbkx3gnAZuQ+ILLRcqd3zZXrkqDgVXQGGMAADUE5wN3bG/lJJpNhbEjuQMe1t5wGwcQLwJXoDYSuKrIbWljpzUwF1IfKHlel4gJaVLB3dK2z80HY0R9PgCADQoWOaYN8JsHIgtCcn21n/MbBxAvOk8VPL4pIM7pAPbTUcTc7yhcxrDgQBhRuILLZeQLPX5ir0fp+WOlDoCABoU6u9FY3vU4U20tyS+gOhKaiVln2nvU+54Ag/nNHApEl8Ij2C545qX47Lej1JHAMAJjh2Rdqyy9/NobI86fAn21l9lNg4gHtHg/qRqF+zipAbuQuIL4dGrUEpsJR0olnasMB1N1FHqCAA4wY5VUuCY1Kqj1DbfdDSIJb4kexuoNhsHEI+Cpeckvk7AxXy4FYkvhEdiqnTGRfb+mpfNxmIApY4AgBNsrWlsnzeitn4EkCh1BEwKzvjasUqqZtZlXVzMh1uR+EL4hModX4q7ywS+0LRgw4EAAGJHcNUwyhxxPF8w8cVJNxB17XpIqW0lf6VU8onpaGIKF/PhViS+ED5nXCwlpEj7NsfdIEI9PACgHsuSti2z91nREccLJr4odQSiz+Ohz9dJUOoItyLxhfBJbm33+pLibnVHVkABANSzb4tUsdsuaes0xHQ0iDWUOgJmdQn2+VpmNo4YQ6kj3IrEF8Kr/5X2Ns7KHb2UOgIA6grOIug0SEpMMRsLYg+ljoBZXYbb220fmo0jxlDqCLci8YXw6l1kr1S053Np9zrT0URN7bRgBgkAgKStNbMIulDmiAZQ6giYlTtMkkfa/4V0aJfpaGIGfYvhViS+EF4pGVLPC+39OCp3DJY6Mi0YACCpTn+vs8zGgdhEqSNgVkqm1KGvvU+frxD6FsOtSHwh/Oqu7hgnGCQAACFVFVLJp/Y+M77QEEodAfPyaHB/PPoWw61IfCH8+lwieROkXWuk3Z+ZjiYqgqWOTPgCAGjHSsnyS+mdpMwupqNBLKLUETAvtLIjfb6C6FsMtyLxhfBLbSv1ON/ej5NZX6FGkIwSAIBQf6+zai+fA3VR6giYF0x8bV8h+UlCS/QthnuR+EJkBMsd18ZH4stDqSMAIChYNpNHmSNOglJHwLysPlJyhnSswq5UAX2L4VokvhAZfS6TPD6p5BNpz0bT0UQc04IBAJIky2JFR5wepY6AeV5vzeqOos9XDc5p4FYkvhAZrdpL3b9k78dBuaOv5r8kZnwBQJzbt1k6XGaXsnUabDoaxCpKHYHYELxAQZ8vSXX7FnNOA3ch8YXICa3u+KLRMKKBUkcAgCRpa82sgU6DpcQUs7EgdlHqCMSGLqzsWBd9i+FWJL4QOf3GSh6vtPMjae9m09FEVGhacMBwIAAAs7bVlDnS3wun4kuq2bGkgN9oKPFq5syZys/PV0pKikaOHKlly5ad8vgZM2aoT58+Sk1NVV5enu644w4dPXo0StEiYroMt7d7PpcO7zUbSwzwUOoIlyLxhchplSXln2Pvr33ZbCwR5mPGFwBAqp01EDyZAhriTajdp9wx6ubPn6/Jkydr2rRpWrFihQYPHqyioiLt2rWrweOfffZZ3XXXXZo2bZrWrl2rP/7xj5o/f77uvvvuKEeOsEtrJ7Xrae9vX242lhgQmvHFOQ1cplmJL66QoNH6X2lvV79oMoqIY5AAAKjqsFS62t6nsT1OJVjqKFHuaMATTzyhiRMnasKECerfv79mzZqltLQ0zZkzp8Hj3333XY0ZM0ZXX3218vPzdfHFF+uqq6467TkQHIJyx5Bgjy+Lcxq4TJMTX1whQZMEyx13rJD2F5uOJmKYFgwA0M6P7FX6WudImV1MR4NYFip1FCs7RllVVZWWL1+uwsLC0H1er1eFhYVaunRpg48ZPXq0li9fHkp0bdq0SQsXLtSll14alZgRYXk1ia+tJDKD5zR+TmrgMk1OfHGFBE3SuqPUbYy97+LVHVkBBQBQr8yx5uQBaJDXJ6nmM0KpY1SVlZXJ7/crOzu73v3Z2dkqKSlp8DFXX321HnjgAZ1zzjlKTExUz549df7555/yQn5lZaXKy8vr3RCjgjO+tq+I+4a9Xi7mw6WalPjiCgmaJbi6o4vLHVkBBQBQm/g6y2wccAZWdnSMJUuW6KGHHtJvf/tbrVixQs8//7xeeeUVPfjggyd9zPTp05WZmRm65eXlRTFiNEnHM6XEVlLlAalsvelojPLRvgUu1aTEF1dI0Cz9rpDkkbZ/6NpyR0odAQDa9qG9JfGFxgiWOwaY8RVNWVlZ8vl8Ki0trXd/aWmpcnJyGnzM1KlTde211+rGG2/UwIED9bWvfU0PPfSQpk+frsBJZghNmTJFBw4cCN22bt0a9veCMPElSLnD7P04L3f0smAXXCriqzpyhQRKz3Z9uSPN7QEgzh3YLh3cIXl8UuchpqOBEwRXdvTT4yuakpKSVFBQoMWLF4fuCwQCWrx4sUaNGtXgYw4fPiyvt/5pk8/nk3TyJuDJycnKyMiod0MMC67Euy2+E19czIdbNSnxxRUSNNuZV9pbl5Y71q6AYjgQAIAZwTLH7P5SUiuzscAZKHU0ZvLkyZo9e7aefvpprV27VrfccosqKio0YcIESdL48eM1ZcqU0PFjx47VU089pXnz5mnz5s16/fXXNXXqVI0dOzaUAIPDBVfi3RrfKzv6arIDXMyH2yQ05eC6V0iuvPJKSbVXSCZNmtTgY5p7hSQ5ObkpoSHW9btCWvjj2nLHNl1NRxRWXlZAAYD4Rn8vNBWljsaMGzdOu3fv1n333aeSkhINGTJEixYtCrVzKS4urnf+cu+998rj8ejee+/V9u3b1aFDB40dO1Y///nPTb0FhFvwb3fZeunIPim1rdl4DAmVOnJOA5dpUuJLsq+QXHfddRo+fLhGjBihGTNmnHCFJDc3V9OnT5dkXyF54oknNHToUI0cOVIbNmzgCkk8Ss+Wuo2WvnjHLncc/QPTEYWVh1JHAIhv9PdCU1HqaNSkSZNOeuF+yZIl9f6dkJCgadOmadq0aVGIDEa07iC17S7t2yxtXy71Kjz9Y1yIUke4VZMTX1whQbOd+TU78bX6RdclvoKljgwSABCH/MeknavsfRJfaCxKHYHYkjfCTnxt/SBuE1/0LYZbNTnxJXGFBM3k4nJHVkABgDhW+qlUfVRKaSO162k6GjgFpY5AbOlylvTx/LhucM/FfLhVxFd1BEKC5Y6S61Z35OoIAMSxUJnjcMnLVys0EqWOQGzJq2lwv225dJJF2NyOHl9wK76dIbrO/Jq9Xf2C2TjCjEECAOJYsLF97nCzccBZKHUEYkvHM6XEVlLlAbvJfRyibzHcisQXoqvfFZI8dtPIfV+YjiZsvDSCBID4FVrRkcQXmoBSRyC2+BKk3GH2/tb4LHf0cU4DlyLxhehKz5byz7H317xoNJRwoscXAMSpw3ulvZvs/dwCs7HAWUKljiS+gJgRXKAkTvt8eb2c08CdSHwh+oLljp8+bzaOMGJaMADEqe3L7W27nlJaO7OxwFlCpY4kvoCYEezztfUDs3EYEjqnYcoXXIbEF6Kv/1clj9de+n3PRtPRhAUroABAnKrb2B5oCkodgdgTnPFVtl46ss9sLAbQvgVuReIL0dcqS+r+JXvfJeWONLcHgDhFY3s0F6WOQOxplWXP4JVqL2zEER/tW+BSJL5gxplft7efumN1Ry+ljgAQfyyrttSxC/290ESUOgKxKW+kvd36vtk4DOCcBm5F4gtm9BtrX+ks/UQq+9x0NC3mpdQRAOLPno3S0f2SL1nKHmg6GjgNpY5AbOoav4kvDzO+4FIkvmBGWjupxwX2/mrnz/oKljpKlDsCQNzYXlMG02mwlJBkNhY4j5cZX0BMCs742rZc8lebjSXK6FsMtyLxBXNctLqjtzbvxRUSIAxmzpyp/Px8paSkaOTIkVq27NTLiu/fv1+33XabOnXqpOTkZPXu3VsLFy6MUrSIWzS2R0v46PEFxKSsPlJypnSsQir91HQ0UeVlVUe4FIkvmNP3Mnua/+61Uuka09G0iLdO5otxAmiZ+fPna/LkyZo2bZpWrFihwYMHq6ioSLt27Wrw+KqqKl100UXasmWLFixYoPXr12v27NnKzc2NcuSIO8EZX7n090IzUOoIxCavV8qrWd1x66kvvLkNpY5wKxJfMCe1jdSr0N5f7exZX/VKHRkogBZ54oknNHHiRE2YMEH9+/fXrFmzlJaWpjlz5jR4/Jw5c7R37169+OKLGjNmjPLz83Xeeedp8ODBUY4cceXYUamkZiYAM77QHJQ6ArErThvcU+oItyLxBbNCqzs+b6+O5VCUOgLhUVVVpeXLl6uwsDB0n9frVWFhoZYuXdrgY15++WWNGjVKt912m7KzszVgwAA99NBD8vv9J32dyspKlZeX17sBTVLysT1TJy1LatPNdDRwIkodgdiVN8LextmML0od4VYkvmBWn69ICSnS3o32SYRD1Z/xZTAQwOHKysrk9/uVnZ1d7/7s7GyVlJQ0+JhNmzZpwYIF8vv9WrhwoaZOnarHH39cP/vZz076OtOnT1dmZmbolpeXF9b3gThQt79XnTEAaDRKHYHYlVsgebzSgWKpfIfpaKLGS6kjXIrEF8xKTpfOuNjed3CTe0odAXMCgYA6duyo3//+9yooKNC4ceN0zz33aNasWSd9zJQpU3TgwIHQbevWrVGMGK4Q6u9FmSOaiVJHIHYlp0vZA+z9OCp3rO3xZTgQIMxIfMG8Ad+wt6udW+5Yr9SRkQJotqysLPl8PpWWlta7v7S0VDk5OQ0+plOnTurdu7d8Pl/ovn79+qmkpERVVVUNPiY5OVkZGRn1bkCThGZ80dgezUSpIxDbQn2+4qfc0RcsdXToORlwMiS+YN4ZF0uJraT9xdL25aajaRZKHYHwSEpKUkFBgRYvXhy6LxAIaPHixRo1alSDjxkzZow2bNigQCAQuu+zzz5Tp06dlJSUFPGYEYcqyqT9X9j7rOiI5qLUEYhtcdjg3uul1BHuROIL5iWlSX0usfcdWu7o9VLqCITL5MmTNXv2bD399NNau3atbrnlFlVUVGjChAmSpPHjx2vKlCmh42+55Rbt3btXt99+uz777DO98soreuihh3TbbbeZegtwu+BFmqzeUkqm2VjgXJQ6ArEt2OB+50dS1WGzsURJqNQxcJoDAYdJMB0AIEka8HXp0wXS6heki38meZ2Xk/V67NlelDoCLTNu3Djt3r1b9913n0pKSjRkyBAtWrQo1PC+uLhY3jp/I/Ly8vTqq6/qjjvu0KBBg5Sbm6vbb79dd955p6m3ALfbRn8vhEGo1LHhkmwAhrXpKqV3kg7ulHaskPLPMR1RxPlobg+XIvGF2NCrUErOlA7ukIqXSvljTEfUZF6PRwHLotQRCINJkyZp0qRJDf5syZIlJ9w3atQovffeexGOCqgRnPGVO8xsHHC2UKljtdk4ADTM47HLHde8KBW/FxeJLy89vuBSzptWA3dKSJb6j7X3P/m72ViaiZp4AIgDllWb+OrCjC+0AKWOQOzrWtNftDg+Lq55WdURLkXiC7Fj4Lfs7eoXpWrnTfsPXiHxM1IAgHvt3SQd3S/5kqWOZ5qOBk7mCya+nPedB4gbXc+2t1uXxUXjKy7kw61IfCF25J8rtc6xTyg2Lj7t4bEmeIWEcQIAXCzY36vTYCmBVUPRAsHEF6WOQOzKHmCvPl95QNq91nQ0EceFfLgViS/EDq9PGvANe/9j55U7emkGCQDuF+rvVWA2DjgfpY5A7PMlSHln2fvFS83GEgVcyIdbkfhCbBn4TXu7/t9S5UGzsTQRzSABIA5sr5nxRX8vtBSljoAzxFGfLw/nM3ApEl+ILZ2HSu16StVHpHULTUfTJNTEA4DLVVdKJZ/Y+8z4QktR6gg4Q95IexsHiS8f5zNwKRJfiC0ejzTo2/b+J8+ZjaWJWAUFAFyu5FN7dk5ae6ltvulo4HSUOgLO0GW45PFJB7ZKB7aZjiaiQucz7u/jjzhD4guxZ0BNuePGN6RDu83G0gT0+AIAl6vb3ytYDwI0F6WOgDMkp0s5A+19l8/6onUL3IrEF2JPVi+75NHyS2teNB1No7EKCgC4HI3tEU6UOgLOESd9vriQD7ci8YXYNLCm3NFBqzuyCgoAuFywsX0uje0RBpQ6As7RtabP19b4SHz5OZ+By5D4Qmwa8HVJHmnbMmnvZtPRNArNIAHAxY7sk/ZssPdzh5mNBe5AqSPgHHln29vS1dLRA2ZjiSBvTXbA4nwGLkPiC7EpPUfq/iV7/9N/mI2lkTyUOgKAe+1YaW/b9ZDS2pmNBe5AqSPgHBmd7EVNrIC09QPT0USMh1JHuBSJL8Suuqs7OuCPL6s6AoCL0d8L4UapI+AsXUfb2+J3zcYRQT5WdYRLkfhC7Oo3VvIlS7vXSaWfmo7mtILN7ZkaDAAutI3EF8KMUkfAWbrVNLj/YqnZOCKI5vZwKxJfiF0pmVLvInvfAU3uvV5mfAGAK1kWM74Qfr6kmh1LCviNhgKgEbqNsbfbP5SOHTUbS4QEL+ST+ILbkPhCbBv4LXv76T9ifs5taBUUMl8A4C4HtkkVuyRvgpQz0HQ0cAtvQu0+5Y5A7GvXQ2rV0Z6luWOF6Wgiggv5cCsSX4htZ1wsJWdK5dtjvp6eUkcAcKngbK/sM6XEVLOxwD2CpY6SFCDxBcQ8j0fqVtPn64t3zMYSIaFSRzJfcBkSX4htiSlS/yvs/Y/mmY3lNGhuDwAuRZkjIiFU6ihmfAFOESx3dGmfL0od4VYkvhD7Bl9lb1e/KB07YjSUUwmVOjJQAIC7bK8paSHxhXDy+iTVnGWS+AKcIdjgfuv7kr/abCwRQKkj3IrEF2Jf11FSZlep6qC07hXT0ZyUt+a/Jq6QAICLBPzSjpX2PokvhFuw3JFSR8AZOva3F+CqOiSVfmI6mrCjZzHcisQXYp/XKw0eZ+9/PN9sLKcQHCjo8QUALrJ7vXSsQkpqLWX1Nh0N3CZY7uivMhsHgMbx+qS8s+39L2K7/3Bz0LMYbkXiC84w6Dv2dsNi6WCp2VhOovYKieFAAADhs/1De9t5aE1pGhBGwZUdXVgyBbhWqMG9GxNflDrCnUh8wRmyekm5wyXLL326wHQ0DaIZJAC4UKix/TCzccCdKHUEnCeY+CpeKrnse39t4std7wsg8QXnGFwz6+ujv5mN4yQodQQAFwolvoabjQPuRKkj4DydhkgJqdLhPVLZZ6ajCSt6FsOtSHzBOQZ8Q/ImSiWfSKWrTUdzAlZBAQCXqTosla6x92lsj0ig1BFwnoQkqUvNxZAtb5uNJcwodYRbkfiCc6S1k3oX2fsxOOsrWOrIKigA4BIlH9sl9q2zpYzOpqOBG1HqCDhT/rn29ot3zMYRZpQ6wq1IfMFZBl9lbz/+e8xdHWWgAACXCZU5Fkg1f+OBsKLUEXCm/DH2dsvbrurzxYV8uBWJLzjLGRdLqe2kQ6XSxjdMR1OPzxvs8WU4EABAeNDYHpFGqSPgTLnDJV+yfU6yZ4PpaMKmtmex4UCAMCPxBWdJSJIGfdve/+hZs7Ecx1MzUHCFBABcYvsKe0t/L0QKpY6AMyWmSHkj7H0X9fmiggVuReILzhMsd1y3UDqyz2wsdQSnBjNQAIALHN4r7dts73ceajYWuBeljoBzdatT7ugSrOoItyLxBefpNFjq2F/yV0qfPm86mhCmBgOAiwRne7XvJaW2NRsL3CtU6siML8Bx8s+xty7q8xWa8RUwHAgQZiS+4DwejzTkans/hlZ3ZGowALhI3cb2QKSESh3p8QU4Tpezavp8lUh7N5mOJiw4n4FbkfiCMw38tuTxSds+kMo+Nx2NpDqroDBQAIDzkfhCNFDqCDhXYorUZbi9v+V/ZmMJE0od4VYkvuBM6dlSry/b+6tio8l97RUSw4EAAFrGsmoTX51Z0RERRKkj4Gx1yx1doO75jEXyCy5C4gvONeQae/vR36SA32wsknzeYI8vBgkAcLQDW6XDZXZSImeg6WjgZpQ6As4WSny944o+X8HEl+SKtwOEkPiCc/W5REptJx3cKW18w3Q0Co4TfqZ8AYCzBWd7ZQ+wS1mASKHUEXC2LmfZ/x0f3OGKPl++Ookvyh3hJiS+4FwJydKgcfb+yr+YjUWUOgKAa9DfC9HirZnxRakj4EyJqVJusM+X88sdPXWyA/QthpuQ+IKzDa0pd1y3UKrYYzQUSh0BwCW2r7C3JL4Qab6aHl+UOgLOFSp3dH6De0od4VYkvuBsOQOlTkOkwDHpk78bDYVSRwBwAX+1tGOlvU/iC5FGqSPgfN2/ZG83v+X4bJG3Nu9FqSNchcQXnG/od+3tir8YHWwodQQAFyhbLx07LCW1lrLOMB0N3I5SR8D5upwlJaRIh0qlss9MR9Mi3no9vgwGAoQZiS8438BvSr5kaddqaecqY2EEr5BwdQQAHCzY36vTEMnrMxoK4gCljoDzJaZIeSPt/c1vmY2lheomvqhigZuQ+ILzpbaV+o2191eYa3JPjy8AcIFgf68ulDkiCih1BNyhbrmjg9UtdeScBm5C4gvuECx3/GSBVHXYSAiemisk/oCRlwcAhENwxlfnYWbjQHyg1BFwh+7n2dst/5MCzj0Z8HkpdYQ7kfiCO3Q/T2rTTao8IK150UgIlDoCgMMdOyKVrrb3aWyPaAiVOpL4Ahyt8xC7N+SRfVLpp6ajaTYPpY5wKRJfcAevVxp2rb2//GkjIfg8lDoCgKOVfCJZfqlVRymzi+loEA9CpY4kvqJt5syZys/PV0pKikaOHKlly5ad8vj9+/frtttuU6dOnZScnKzevXtr4cKFUYoWMc+XKHUbbe+7pNyRcxq4CYkvuMeQ70oen7T1PWnXuqi/fKjUkUECAJwpWOaYWyDVueoNRAyljkbMnz9fkydP1rRp07RixQoNHjxYRUVF2rVrV4PHV1VV6aKLLtKWLVu0YMECrV+/XrNnz1Zubm6UI0dMc0mfr2C5IxO+4CYkvuAeGZ2k3kX2/oo/R/3lg6ugMEgAgEOFEl/090KUUOpoxBNPPKGJEydqwoQJ6t+/v2bNmqW0tDTNmTOnwePnzJmjvXv36sUXX9SYMWOUn5+v8847T4MHD45y5IhpwcTXF+84OpnNxXy4EYkvuMuw6+ztR3+Tqiuj+tK+mv+a6PEFAA4VXNGRxBeihVLHqKuqqtLy5ctVWFgYus/r9aqwsFBLly5t8DEvv/yyRo0apdtuu03Z2dkaMGCAHnroIfn9/pO+TmVlpcrLy+vd4HLZA6WUNlLVIWnHKtPRNFuobzFX8+EizUp8UROPmNWrUErvLB3ZK639Z1RfOjTji0ECAJzn8F5p70Z7nxUdES2UOkZdWVmZ/H6/srOz692fnZ2tkpKSBh+zadMmLViwQH6/XwsXLtTUqVP1+OOP62c/+9lJX2f69OnKzMwM3fLy8sL6PhCDvF6p+7n2/ub/mo2lBbyhvsWGAwHCqMmJL2riEdN8CdLQ79r7K6Lb5N5DqSMAONeOlfa2bXcprZ3ZWBA/fDWJL0odY1ogEFDHjh31+9//XgUFBRo3bpzuuecezZo166SPmTJlig4cOBC6bd26NYoRw5ju59nbTUuMhtESvtA5DSc1cI+Epj6gbk28JM2aNUuvvPKK5syZo7vuuuuE44M18e+++64SE+3BPT8/v2VRA6cy7FrprV/YjSXLNkhZvaLysqFpwQwSAOA8O4JljgVm40B88THjK9qysrLk8/lUWlpa7/7S0lLl5OQ0+JhOnTopMTFRPp8vdF+/fv1UUlKiqqoqJSUlnfCY5ORkJScnhzd4xL4e59vbre9LVYelpDSj4TRHcG0XenzBTZo04ytaNfFAi7TpKp1xsb2//E9Re9ngCiiMEQDgQNtJfMEASh2jLikpSQUFBVq8eHHovkAgoMWLF2vUqFENPmbMmDHasGGDAoFA6L7PPvtMnTp1ajDphTjWvpeU0UXyV0nFDZ8fxzpv6JyGkxq4R5MSX9GqiacZJFps+A32dtUz0rEjUXnJ0Aoo1DoCgLNYFis6wgxKHY2YPHmyZs+eraefflpr167VLbfcooqKilBFy/jx4zVlypTQ8bfccov27t2r22+/XZ999pleeeUVPfTQQ7rttttMvQXEKo+ndtbXpjeNhtJcPtq3wIWaXOrYVHVr4n0+nwoKCrR9+3b94he/0LRp0xp8zPTp03X//fdHOjS42RkXSZl50oGt0pqXpMHfifhLUuoIAA5VvkM6VCp5fFLOINPRIJ5Q6mjEuHHjtHv3bt13330qKSnRkCFDtGjRotDF/eLiYnm9tfMD8vLy9Oqrr+qOO+7QoEGDlJubq9tvv1133nmnqbeAWNbjfGnVXx3b54uL+XCjJiW+olUTP2XKFE2ePDn07/LyclZCQdN4fVLBddIbP5M++GNUEl9cHQEAhwrO9sru78h+LHAwSh2NmTRpkiZNmtTgz5YsWXLCfaNGjdJ7770X4ajgCj1qGtyXfCId2i217mA2nibiYj7cqEmljtGqiU9OTlZGRka9G9BkQ8dL3gRp2zJ74ImwYD18gMwXADhLMPHVmTJHRBmljoD7tO4oZQ+w9zf/12wszeD10LcY7tOkxJdETTwcJD1b6nu5vf9h5Jvce7g6AgDOxIqOMIVSR8CdHNznK7hgF+c0cJMm9/iiJh6OMvwGac2L0sfzpcKfSimRmz1IqSMAOFAgIO1YZe+T+EK0UeoIuFOPC6Slv5E2LrGnTgWvkDtAMFR6fMFNmtXcnpp4OEb3L0lZvaWyz+zk14iJEXspr4erIwDgOHs2SJXlUkKq1KGv6WgQbyh1BNyp2yjJlySVb5P2bJSyepmOqNG8XMyHCzW51BFwFI9HOqsm2bXs9xEtVqfUEQAcKNTfa4jki/hi10B9lDoC7pTUSsobae87rNwxWOpocU4DFyHxBfcb/B0pqbU96yuCDSa5OgIADhRMfFHmCBModQTcK9Tna4nJKJqMUke4EYkvuF9Khp38kqRlsyP2MjSCBAAHCja27zzUbByIT5Q6Au7V8wJ7u/ktRyW3uZgPNyLxhfgQLHdcv1DavzUiL+ENljoySgCAM1RXSiWf2PvM+IIJviR7awWkgN9sLADCq9MQKbWd3Udy2wemo2m04IJdlDrCTUh8IT507Gs3urcC0odzIvISHprbA4CzlH4q+avsE5O2+aajQTzy1ukr56AZIQAaweuTel5o729YbDaWJgiVOnJOAxch8YX4MeIme7viaenY0bA/fW2pY9ifGgAQCdtryhxzhzlqqXm4SLDUUaLcEXCjXl+2txv+YzaOJqDUEW5E4gvxo/clUmaedHiP9OmCsD89pY4A4DChxBdljjAkWOooMeMLcKPgjK+dH0kVZWZjaSRvTYaAKha4CYkvxA9fgjSiptfX0t9KYf5jTqkjADhMcEXHzsPMxoH45fVJCtYVkfgCXCc9R8oeKMmSNr5pOppGoccX3IjEF+LLsPFSYpq0a7W9wkoY+ZgWDADOcbRcKvvM3s8l8QWDWNkRcLdewT5fzih3DF7M9wcMBwKEEYkvxJfUttKQq+39954K61MzLRgAHGTnKkmWlNlVat3RdDSIZ8FyR3+V2TgAREbPmj5fG9+QArGfTQq1b+GcBi5C4gvxZ+Qt9vazRdKejWF7Wi+ljgDgHMEyx9yhZuMAgis7+qvNxgEgMrqeLSW2kip22asJx7jggl2UOsJNSHwh/mT1ks4okmRJ788K29OGenzF/oUcAECwsT39vWAapY6AuyUkS93PtfcdUO5IqSPciMQX4tPZNbO+Vj4jHdkflqf0MeMLAJwjmPjqMtxsHACljoD71S13jHGUOsKNSHwhPvU4X+p4pnSsQvpwTliekkECCJ+ZM2cqPz9fKSkpGjlypJYtW9aox82bN08ej0dXXnllZAOEsx0slcq3SfJInQabjgbxjlJHwP161SS+ipfai6vEMNq3wI1IfCE+eTzS6B/Y++89JR07GoanZFVHIBzmz5+vyZMna9q0aVqxYoUGDx6soqIi7dq165SP27Jli370ox/p3HPPjVKkcKwdNbO9OvSVktPNxgJQ6gi4X/ueUrueUqBa2rTEdDSnVNvjy3AgQBiR+EL8GvhNKaOL3Wjyo7+1+OmCgwRXR4CWeeKJJzRx4kRNmDBB/fv316xZs5SWlqY5c04+O9Pv9+uaa67R/fffrx49ekQxWjhSqLE9/b0QAyh1BOJD7yJ7+/mrZuM4jdoeX5zTwD1IfCF++RKlUbfZ++/+Sgr4W/R0oVJHBgmg2aqqqrR8+XIVFhaG7vN6vSosLNTSpUtP+rgHHnhAHTt21Pe+971GvU5lZaXKy8vr3RBHgv29SHwhFlDqCMSHMy6yt5+/HtPTqWjfAjci8YX4Nmy8lNJG2rtJWvtyi57KS6kj0GJlZWXy+/3Kzs6ud392drZKSkoafMzbb7+tP/7xj5o9e3ajX2f69OnKzMwM3fLy8loUNxzEsurM+CowGwsgUeoIxItuY6TEVtKhUmnnR6ajOanggl3kveAmJL4Q35JbSyNusvffntGiv/Aero4AUXfw4EFde+21mj17trKyshr9uClTpujAgQOh29atWyMYJWLK3k3S0f12eVnHM01HA1DqCMSLhGR7gS3JnvUVo0KljpzTwEUSTAcAGDfyZundX0s7V9nNJnte0Kynqe3xFb7QgHiTlZUln8+n0tLSeveXlpYqJyfnhOM3btyoLVu2aOzYsaH7AoGAJCkhIUHr169Xz549T3hccnKykpOTwxw9HGHHSnubM0hKSDIbCyDVKXVkxhfger0vlta/Yvf5Ou/HpqNpEKWOcCNmfAGtsqRh19r7/3u82U8TKnUk8wU0W1JSkgoKCrR48eLQfYFAQIsXL9aoUaNOOL5v37765JNPtGrVqtDtiiuu0AUXXKBVq1ZRwogTUeaIWBMqdaTHF+B6vWr6fG37UKrYYzaWk+BiPtyIGV+AJI25XfrwT9KW/0lfvCt1G93kp6DUEQiPyZMn67rrrtPw4cM1YsQIzZgxQxUVFZowYYIkafz48crNzdX06dOVkpKiAQMG1Ht8mzZtJOmE+wFJrOiI2EOpIxA/MnOl7IFS6SfShv9Ig8eZjugEXMyHGzHjC5CkzC7S0Gvs/f8+2qyn8IWa2zNIAC0xbtw4PfbYY7rvvvs0ZMgQrVq1SosWLQo1vC8uLtbOnTsNRwlH8h+Tdn5s7zPjC7GCUkcgvoRWd3zNbBwnwcV8uBEzvoCgcyZLK/8qbXpT2vqBlHdWkx7uZVowEDaTJk3SpEmTGvzZkiVLTvnYuXPnhj8guMOutVL1ESk5Q2p3Yu83wAhKHYH40rtIevsJe8aXv1ryxdYpOSvVw42Y8QUEte0mDf6Ovf9W02d90QgSAGLcjhX2tvNQyctXIMQISh2B+JI7XEpta68wvG2Z6WhOEOzxZXFOAxfhWx9Q1zmTJY/Xnnq8fUWTHuql1BEAYhuN7RGLvDUzvih1BOKDL0E6o8jeX7/QbCwNCJY6+pnyBRch8QXU1b6nNPDb9v6bP2/SQ2sbQYY7KABAWAQvaNDYHrEkWOZEqSMQP/pcYm/X/9tsHA2g1BFuROILON75d9pXXzf8R9r4ZqMfxowvAIhhVRXSrjX2PjO+EEsodQTiT68v2//t79kglX1uOpp6WLALbkTiCzheux7SWd+z91+/r9FTuFgBBQBi2M6PJSsgpXeSMjqbjgaoRakjEH+S06X8c+z9GCt3DLbADDDlCy5C4gtoyJd+Yq/6VfKx9MlzjXqIj1UdASB20d8LsYpSRyA+9bnU3sZYuaOHUke4EIkvoCGt2kvnTrb333hQOnb0tA+p7fHFKAEAMSeU+KK/F2IMpY5AfOr9FXu79X2posxsLHVQ6gg3IvEFnMzI70sZXaQDW6X3nzrt4V5KHQEgdjHjC7GKUkcgPrXJk3IG2mX4n79mOpoQzmngRiS+gJNJTJW+PNXeX/LwaRtPeil1BIDYVFEm7f/C3u80xGgowAlCpY4kvoC4Eyp3jJ0+Xx5mfMGFSHwBpzJonNTjAqn6qPTiLZL/5P03KHUEgBi1fYW9bX+GlNrGaCjACUKljiS+gLjT5xJ7u+GNRrVWiQYvPb7gQiS+gFPxeKSv/sZudL/tA+ndX530UKYFA0CM2lGT+KLMEbGIUkcgfnUaIqV3lo5VSJv/azoaSZIvuKoj5zRwERJfwOlkdpEuecTeXzJdKl3d4GFcHQGAGEV/L8QySh2B+OXxSH0vs/fX/tNsLDWoYoEbkfgCGmPwVVLvS+wVl567XjpYesIhwR5ffq6OAEDssCwSX4htlDoC8a3f5fZ2/cJTtlWJFg8X8+FCJL6AxvB4pLFP2lORyz6T5l4mle+sd0iw1NEi8QUAsWP/F9LhPXY5Wc4A09EAJ6LUEYhv3cZIqW3tsap4qeloKHWEK5H4AhorPVua8IqUmSft+Vyae6l0YFvox5Q6AkAMCs72yhkoJSSbjQVoiC+Y+KoyGwcAM3yJtas7xkC5I6WOcCMSX0BTtOshXf+K1KartHeT9LvzpLd+IR3ZVyfxxSABADEjuKJj7jCzcQAn42PGFxD3+o21t+v+ZZfoG0SpI9wowXQAgOO07SZdv1D669ftssc3fib975dK7/s1jff5tNHqLGtvf3kSU+3ZBcHeHXUHMY/XLp/0eCWPT/L67H8DAMKL/l6Idb6amYjM+ALiV48LpMRWUvl2eyVig2MWK9XDjUh8Ac3RJk+65V1p9QvS2zOkXauV8slf9UDNRVv9anoLX+C4JJjXV5skS0iyvyQnpEhJaVJSK/uWkimltLF7BKR3smPM7CJ16CslprYwHgBwIH+1tPMje78zM74QoxKCze1JfAFxKzFFOuMiac2LdrmjwcSXjxlfcCESX0Bz+RKlQd+WBn5L2rhYRz//r9569x319OxQj+QD8hw7Iqm5I8ZxjwvUWeGl+kgT40yS8kZK+edK/b8qdezbzJgAwGF2r5WOHZaS0qWs3qajARoWnPFVXWk2DgBm9Rtbm/j68jRj1SDBlerp8QU3IfEFtJTHI/Uq1LEu5+mm/74mSVo39StKSfDaCavqypqBy2NvLUuyAvVvAb9OmiSre3yg2r4iXF1p345VSFUVUuUh6eh+6ch+6cheqXyHdGCrtO8L6XCZtOV/9u2/j0ijfyCdfxezwAC4X6jMcajkpa0pYlRw0QU/iS8grp1xsX3Bes8Gafc6qWM/I2F4KHWEC5H4AsLEW+eqjGXJHjV8ibVNa02wLLsJ/+b/SusWShtel96ZYV9J+upvpG6jzcUGAJEWSnwNNxsHcCrBXqDVlDoCcS0lw+719fmr0pqXjSW+KHWEG3H5EwgTn7c28RUzV0g8Hql9T2n4DdJ3F0jfedbu/7V3o/T0FdKutaYjBIDICa3oSGN7xDBmfAEIOvNKe7v6BWMhsFI93IjEFxAmdcvw/bE6UPS9TLrtfanH+VLgmPTvnxhfMhkAIqKqQtq1xt4n8YVYxowvAEF9LrX/Juxea+wCNaWOcCMSX0CY1Ct1DBgM5HRSMqXLZ9jNdDe/Ja15yXREABB+Oz+y/xind5YyOpmOBjg5ZnwBCEptI/X8sr1vaNZXsIqFUke4CYkvIEzqJr5i/gpJu+7SOT+091+9R6o6bDQcAAi7UH+vYWbjAE4nOOMrUC0FYvnKGYCoOPNr9nb1C0YqM0KljmS+4CIkvoAwqdPiK/YTX5I05odSZp5Uvk16+wnT0QBAeG370N5S5ohYF0x8Scz6AiD1ucSuzCj7rLZkP4q8lDrChUh8AWHi8XhCNfEx2+OrrqQ0qejn9v47v5IObDcbDwCEE43t4RTBUkdJ8tPnC4h7KRlSr0J730C5o4fm9nAhEl9AGAWnBjtmnOh3hdRpiH2FefN/TUcDAOFxaJd0oFiSR+o81HQ0wKnVnfFFg3sAktFyR3p8wY1IfAFh5HPaFRKPR8o/x94PlgUBgNMFZ3t16GNfOQdimcdTm/yi1BGAJPX5il3uuGeDVPppVF86VOpI5gsuQuILCKNQqaOTBopgGdB2El8AXCLU2J4yRziEr6bcsZrEFwBJyenSGRfZ+5/+I6ovTakj3IjEFxBGjit1lKQuw+1t6Wrp2BGzsQBAOLCiI5wmITjji1JHADUGfMPefvKPqK74WlvBErWXBCKOxBcQRrU18Q4aKTLzpFYd7WXUd35kOhoAaBnLYsYXnIcZXwCO1+cSKSnd7lm59f2ovay3JkPgqPMZ4DRIfAFh5MhSR4+ndtYXfb4AON3eTdLR/XYioeOZpqMBGocZXwCOl5gq9Rtr73/y96i9rJdSR7gQiS8gjLxOnRpMny8AbhGc7dVpUG0yAYh1web2zPgCUNegb9nb1S9EbdXXUI+v6FVXAhFH4gsIo+AqKJbTrpCEZnwtNxsHALRUcOZq7nCzcQBNESx1ZFVHAHXlf8luSXJkn7Txjai8pONWqQcagcQXEEa1Pb4MB9JUnYdJ8tg9BA7tMh0NADQf/b3gRMHZiVGa0QHAIXwJdZrcR6fcMXghn8QX3ITEFxBGwanBjurxJUkpGVKHPvb+dmZ9AXCo6iqp5GN7vwuJLzgIM74AnEyw3HHdQqnyYMRfzuPU1i3AKZD4AsLI0VdIcmlwD8DhSj+xm4OntpPadjcdDdB4oeb2x8zGASD2dB4mtespVR+R1r0S8Zdz5Cr1wGmQ+ALCKFgT78hxogsN7gE43PYV9ja3oHaZXcAJgjO+aG4P4HgejzTo2/b+R/Mi/nKhC/lM+YKLkPgCwihU6ujEzFdwxtf2FSzjAsCZQo3tKXOEw4RmfJH4iqaZM2cqPz9fKSkpGjlypJYtW9aox82bN08ej0dXXnllZAMEggaNs7eblkgHtkf0pRy7Sj1wCiS+gDDy1vwX5cipwR37S4lpUmW5tOdz09EAQNMFexR2YUVHOExoxhfN7aNl/vz5mjx5sqZNm6YVK1Zo8ODBKioq0q5dp17kZ8uWLfrRj36kc889N0qRApLadZe6nSPJkj6O7KwvL6WOcCESX0AY1ZY6OnCg8CVInYbY+ztWGg0FAJrsyL7apH3nYWZjAZoqgeb20fbEE09o4sSJmjBhgvr3769Zs2YpLS1Nc+bMOelj/H6/rrnmGt1///3q0aNHFKMFJA252t6ufCaifVWCpY6OW6wLOAUSX0AYeUOrOhoOpLmyetnbfV+YjQMAmiqYsG/bXWrV3mwsQFP5akodmfEVFVVVVVq+fLkKCwtD93m9XhUWFmrp0qUnfdwDDzygjh076nvf+95pX6OyslLl5eX1bkCL9P+qlNhK2rtR2tq4stzm8Dq5ZzFwEiS+gDDyOHlVR0nKzLO3B7aajQMAmipY5kh/LziRjx5f0VRWVia/36/s7Ox692dnZ6ukpKTBx7z99tv64x//qNmzZzfqNaZPn67MzMzQLS8vr8VxI84lt7aTX5K06pmIvYzjz2eABpD4AsKothmkQweKzC729sA2s3EAQFNto78XHCzY3J5VHWPSwYMHde2112r27NnKyspq1GOmTJmiAwcOhG5bt3JREWEQLHdc/YJUdTgiL+Fz+vkM0IBmJb5YAQVomM/r8KnBJL4AOJFlMeMLzhZsbu+n1DEasrKy5PP5VFpaWu/+0tJS5eTknHD8xo0btWXLFo0dO1YJCQlKSEjQn//8Z7388stKSEjQxo0bT3hMcnKyMjIy6t2AFus2RmrT1V6Mat0rEXmJ2ub2EXl6wIgmJ75YAQU4OU+ox5dDR4pQqeM2B2fvAMSdA1ulil2SN0HKGWg6GqDpgs3tmfEVFUlJSSooKNDixYtD9wUCAS1evFijRo064fi+ffvqk08+0apVq0K3K664QhdccIFWrVpFGSOix+uVBgeb3P8lMi9BqSNcqMmJL1ZAAU7O8QNFRmdJHqn6iHR4j+loAKBxtn1ob7MHSImpZmMBmiPU4+uY2TjiyOTJkzV79mw9/fTTWrt2rW655RZVVFRowoQJkqTx48drypQpkqSUlBQNGDCg3q1NmzZKT0/XgAEDlJSUZPKtIN4MuVqSR9r8X2nvprA/veNbtwANaFLiKxoroABO5vhSx4RkqXVNo1ca3ANwiu3094LDBWd80dw+asaNG6fHHntM9913n4YMGaJVq1Zp0aJFoYb3xcXF2rlzp+EogQa07Sb1vNDeX/HnsD99KPHl1FXqgQYkNOXgU62Asm7dugYfE1wBZdWqVY1+ncrKSlVW1g78LP8Lp3B8qaNk9/k6VGKXO3YeajoaADi94IyvXBJfcCgfze1NmDRpkiZNmtTgz5YsWXLKx86dOzf8AQGNVXC9tHGxtPKv0vl31y6QEQbM+IIbRXRVx+asgCKx/C+cy/GljhIN7gE4i/+YtHOVvc+MLzhVAs3tATRBn0vsKo2K3dJn/w7rU3trMgSOPp8BjtOkxFc0VkCRWP4XzlV7hcRwIC3RpibRvJ//7gA4QOlqqfqolJIptetpOhqgeXw0twfQBL5Eaeh37f0P/xTWp/aGKljC+rSAUU1KfEVrBRSW/4VT+TzBHl8OznyFVnYk8QXAAbYHyxwLai9TA04TLFNixheAxho23t5uelPauzlsT+t1w/kMcJwm9fiS7BVQrrvuOg0fPlwjRozQjBkzTlgBJTc3V9OnTw+tgFJXmzZtJOmE+wE3qBkn5HfyQEGpIwAn2VbT2J7+XnAyenwBaKq2+XaT+41v2E3uC6eF5Wld0boFOE6TE1/jxo3T7t27dd9996mkpERDhgw5YQUUL1dcEadcUepI4guAkwRnfNHfC07mY8YXgGYomGAnvlb+RTr/rtp+gS3g9brgfAY4TpMTXxIroAAn4/O6YGpwsNSxYpd07KiUmGI2HgA4mSP7pbLP7P3cAqOhAC2SQI8vAM3Q5xIpvbN0cIe0+kVp8LgWP2XoQj6ZL7gIU7OAMAqVOjp5oEhtKyW2svfLt5uNBQBOZccKe9s2X2rV+NWjgZgTbG7vJ/EFoAl8idJZN9j77z8lheHiO6WOcCMSX0AYuaLU0eOpU+5Ig3sAMYz+XnCLUHP7Y2bjAOA8BRPs5PmOldK2D1v8dK44nwGOQ+ILCCNfqCbe4SMFfb4AOAH9veAWPkodATRTqyxp4Dft/WW/a/HTBXt8OXqxLuA4JL6AMApNDXb6JZI2NX2+SHwBiFWWVXtlmxlfcLoEmtsDaIERN9nb1S9I5Ttb9FTB8xlH9ywGjkPiCwgjj1umBlPqCCDW7f9COlwmeROlnIGmowFahhlfAFqi8xAp72wpUC0t/1OLnopSR7gRiS8gjFzTDDK4suN+El8AYlRwtlfOQFafhfMFV3UMHJMCAbOxAHCmkTfb2w/n2CuzN1Nt4svh5zNAHSS+gDAK9vhy/NRgenwBiHXbaxrb098LbuBLqt2n3BFAc/QbK2XkShW7pY/+1uynqS11dME5DVCDxBcQRsFSR7/T5wbXTXwx4AGIRds+sLddzjIbBxAO9RJflDsCaAZfojRqkr3/7q+kgL9ZTxOc8SVR7gj3IPEFhJFrauIzciV57C/fFWWmowGA+qorpZ0f2fvM+IIb1E18VTPjC0AzDRsvpbaV9m6S1r7crKcIruooUe4I9yDxBYSRzy09vnyJUnone/9AsdlYAOB4JZ/a5WBp7aW23U1HA7Sc12sv1CAx4wtA8yW3lkbU9Pp6+5fNqtyok/dyfhULUIPEFxBGrmoGSZ8vALEqWOaYO1yqU5IBOFoCKzsCCIMRN0mJafbM6E1LmvzwuqWObjilASQSX0BYedxS6iiR+AIQu+jvBTcKljv6j5mNA4CztWpvlzxK9qyvJvJR6ggXIvEFhJGv5r8oVwwSbfLs7f6tZuMAgONt/9DedikwGwcQTsEZX5Q6AmipUbdJ3gRp83+l4veb9NC6E6n9bjinAUTiCwirUKmjG6Z8ZdYkvsqZ8QUghhzaLe3bIskj5ZL4gosEZ3zR3B5AS7XpKg252t5ffH+TahbrlToGwh0YYAaJLyCMKHUEgAgLzvbK6i2lZJqNBQgnZnwBCKfz7pJ8ydIX70gbFjf6YXUTX66oYgFE4gsIK69bVnWUahNflDoCiCXbgmWO9PeCy/hobg8gjDJzpRET7f3FP5UCjZu+VXdVR1ec0wAi8QWEVbAZpKtmfB0uk44dMRsL4s7MmTOVn5+vlJQUjRw5UsuWLTvpsbNnz9a5556rtm3bqm3btiosLDzl8XC4UGP74WbjAMLNl2hv/ZQ6AgiTc/+flJwhlXwirX6+UQ/xeDyhPl/0+IJbkPgCwshVPb5S2khJre39A9uNhoL4Mn/+fE2ePFnTpk3TihUrNHjwYBUVFWnXrl0NHr9kyRJdddVVevPNN7V06VLl5eXp4osv1vbtfG5dJ+CXtq+w90l8wW0SmPEFIMzS2kmj/8/ef+NnjV41NnhOQ94LbkHiCwgjj5tKHT2e2gb3Byh3RPQ88cQTmjhxoiZMmKD+/ftr1qxZSktL05w5cxo8/plnntGtt96qIUOGqG/fvvrDH/6gQCCgxYsb388CDlH2mVR1UEpsJXXoZzoaILyCze2Z8QUgnM6+RWrVQdq3WVr6m0Y9xBfqW+yCcxpAJL6AsPK5qbm9RIN7RF1VVZWWL1+uwsLC0H1er1eFhYVaunRpo57j8OHDOnbsmNq1axepMGFKsMwxd5jkSzAbCxBuzPgCEAnJraWLHrD335wu7V5/2oeESh1dc1KDeEfiCwgjr9dlV0dCiS9mfCE6ysrK5Pf7lZ2dXe/+7OxslZSUNOo57rzzTnXu3Lle8ux4lZWVKi8vr3eDA4QSXwVm4wAiITTji8QXgDAbfJV0xsX235cXb7VbB5wCpY5wGxJfQBiFSh3dcnWEGV9wmIcffljz5s3TCy+8oJSUlJMeN336dGVmZoZueXl5UYwSzba1JvGVN8JsHEAkBGd8NbIHDwA0mscjXT7DbnS//UNp6cxTHu6qleoBkfgCwsp9pY70+EJ0ZWVlyefzqbS0tN79paWlysnJOeVjH3vsMT388MN67bXXNGjQoFMeO2XKFB04cCB027qVz3jMO3pA2r3O3u9yltlYgEjwUeoIIIIyc6Wih+z9N38ula4+6aFeN61UD4jEFxBWXrc1gmTGF6IsKSlJBQUF9RrTBxvVjxo16qSPe/TRR/Xggw9q0aJFGj789Kv9JScnKyMjo94NMW77ckmW1Kab1Lqj6WiA8EuguT2ACBv6XalXoVR9VHr6ipMmv4LnNPT4gluQ+ALCyHXTgtsEZ3xtlwIBs7EgbkyePFmzZ8/W008/rbVr1+qWW25RRUWFJkyYIEkaP368pkyZEjr+kUce0dSpUzVnzhzl5+erpKREJSUlOnTokKm3gEjY9qG9pcwRbsWMLwCR5vFIX58t5QySDpdJcy+Xdn58wmHBcxrLLec0iHskvoAw8rhtxld6J8njtRthHi4zHQ3ixLhx4/TYY4/pvvvu05AhQ7Rq1SotWrQo1PC+uLhYO3fuDB3/1FNPqaqqSt/85jfVqVOn0O2xxx4z9RYQCVuX2VvKHOFWCTS3BxAFae2k616WOg+TjuyVnh4rrX6xXsN7H6WOcBnWAgfCyHWDhC/RTn6Vb5f2b6W8CFEzadIkTZo0qcGfLVmypN6/t2zZEvmAYJZl1a7oSOILbhVc1bGaUkcAEZbaVhr/ovTMt6St70vPXSe16yGNmiT1vUw1E74odYRrkPgCwsjrtlUdJbvPV/l2u8F9lwLT0QCIR3s2SEf3SwkpUvYA09EAkREsdWTGF4BoSMmUvvu89M4Madlsae8m6ZXJ0iuT9YbStCGps7osypfS06WEVMnrlSxJVkCSZW9PqHKxau5rYBs6xDru2Eaod1xDz3XaJ2jiazTyeY6Pq+57O+EY69Tvve7vyQr+fgPHvaZH8iZIvgT7/5M2eVLb7lL7nlKPC6TEk69oHu9IfAFh5LpSR8lOfG19nwb3AMwJzvbqPLS2HAxwmwRmfAGIsuTW0oX3SmN+KK38q7R8rlS2XunWYQ31bpCKN5iOEI3Vsb/0raelDr1NRxKTSHwBYeS6UkeJlR0BmBfq73X6FTsBx2LGFwBTkltLZ3/fvlVX6ruP/k2tD23W3RfmqGu6Vzp2RLL8du9feewm+XX3gyzLvt/jOe5nxx1Xdz9YWFnvvuN+drLHnHD4Se4/1fOe9JCTHXPc/SeLK/g7qHfMaX4PdR/j9Z14nxWQAtX2reqQtO8Lad8WafNb0q410u/Pl8Y+KQ361unfX5wh8QWEkTtLHYMrO241GweA+BVc0bELKzrCxUIzvkh8ATAoIVlfJHTV1kAH3dR7tLp2bWs6IpzOwVLpH9+TtvxPev5Gu03NOT80HVVMYVVHIIy8rix1JPEFwKDKg9Ku1fY+je3hZqEZX8fMxgEg7gXPaSw3ndO4WXq2NP4laczt9r+XzTYbTwwi8QWEUW3iy3Ag4USpIwCTdqy0p/Zn5kkZnUxHA0ROAqWOAGKDz43nNG7n9Unn3Sl5fFL5Nmk/kxbqIvEFhFGw1NHvpqsjwcTX4T1S1WGzsQCIP/T3Qrzw0dweQGwItp7yk/lylqRWUqfB9n7xUrOxxBgSX0AYeb0unBackiklpdv75dvNxgIg/gRXdKS/F9yOGV8AYoQr27fEi26j7S2Jr3pIfAFh5AkOEgHDgYSTx1On3JEpswCiyLKkre/b+3kjzcYCRBozvgDECF/oYr7hQNB0Xc+2t8XvmY0jxpD4AsLI59arI/T5AmDCng3SkX1SQoqUM9B0NEBkBRNfzPgCYFjwYj6ljg6UV5P42rVGOrzXbCwxhMQXEEbBHl+uS3y1qVnZkSaJAKIpONur8zApIclsLECkBUsdq0l8ATDLtec08aB1B6l9L3s/2CcVJL6AcHLlqo4SpY4AzAiVOdLfC3EgNOOLUkcAZgXPach7OVTXUfaWPl8hJL6AMAo2t3fd1ZH2Z9jb0k/NxgEgvgSvVNLfC/GAGV8AYkTwnIZSR4ci8XUCEl9AGHnduvRv56H2dtda6dhRs7EAiA9H9ku719n7zPhCPGDGF4AYQamjwwUb3G9fIR07YjaWGEHiCwgj104LzuwipbWXAtXSrtWmowEQD7Z9aG/b9ZBaZZmNBYiG4IwvEl8ADHNt+5Z40a6H1DpbChyTdqw0HU1MIPEFhJFrSx09ntpZX/zxBBANof5elDkiTvjqJL7c9j0CgKP4Qhfz+VvkSB5P7ayvL941G0uMIPEFhJFrSx0lqdMQe7tjlckoAMQLGtsj3tRduZRZXwAM8gTPaUh8OVeoz9d7ZuOIESS+gDBybamjVGfG1yqjYQCIA/5qaftye58ZX4gXwRlfEg3uARhFqaMLBC8c7lxlNIxYQeILCCNXN4LsPMTe7l5Lk0QAkbVrjVR1SErOkDr0NR0NEB0+ZnwBiA0+L6WOjpfRxd5WlNkXFOMciS8gjGqvjrhwkMjIlVp1sBvcl9LgHkAEBcscuwyXvD6zsQDR4vVK3gR7nxlfAAzyuLl9S7xIay/JI8mSDu8xHY1xJL6AMAomvvxuHCM8njp9vmhwDyCCti6zt5Q5It6EGtyT+AJgDqWOLuBLqEl+SarYbTaWGEDiCwgjb81/Ua6dFhwsd6TPF4BI2lrTiLXLWWbjAKIt2OC+mlJHAOa4un1LPGnVwd5W7DIbRwwg8QWEkatLHaXaBvc0SQQQKQe2S/uLJY+XFR0Rf5jxBSAG0OPLJVrXJL4OMeOLxBcQRqFSx4DhQCIlWOq4iwb3ACKkeKm9zRkoJaebjQWINmZ8AYgBHref08SLVh3tLaWOJL6AcAomvlx7dSSjs/0H1PJLJZ+ajgaAGxXXlDl2HWU2DsCE0IwvEl8AzKHU0SVaBxNflDqS+ALCyPWDhMdT2+eLckcAkRBKfJ1tNg7AhARKHQGYR6mjS7TKsreUOpL4AsLJ642DFVCCfb5Y2RFAuB09IJXWzCZlxhfikY9SRwDm1ZY6uvmkJg60YsZXEIkvIIxCze3dPEh0HmZvP1skHdlnNhYA7rL1A0mW1La7lJ5jOhog+pjxBSAG1C7YZTgQtExrenwFkfgCwsj1pY6S1PNCKauPdHiP9MbPTEcDwE2Cje2Z7YV45Uu0t8z4AmCQLx7OaeJBK1Z1DCLxBYRRXJQ6JiRJlz1u73/wR2n7CrPxAHAP+nsh3vmY8QXAvNoZX24+qYkDwcRXxW4pzv+/JPEFhJE3Xurhu58rDfy2JEt65f9JAb/piAA4XXWVtP1De58ZX4hXwVLHahJfAMzxUOroDsHEV+CYdHS/0VBMI/EFhFGw1DEuVkC5+GdScoa0Y4X04RzT0QBwup0fSdVHpdR2UtYZpqMBzAg2t/dT6gjAnLho3xIPElOk5Ex7P87LHUl8AWEUV40g07OlC++19xf+SPr7ddLuz8zGBMC56vb3qvlbCsQdZnwBiAG+mswXeS8XaJVlb+N8ZccE0wEAbhIqdYyXUWL496TST6UVf5HWvCitfVk6o0jK7i+17yVldpES06SEFPsWPJn1eOyR1LIkKyCpZhu6WTX31WxPpe4xoUMbeoxH8nolj1dKbCWltpFSMmsbCQMwK5T4or8X4hgzvqJq5syZ+sUvfqGSkhINHjxYv/71rzVixIgGj509e7b+/Oc/69NPP5UkFRQU6KGHHjrp8YCTeeKlfUs8aN1R2rtROkTiC0CYeGvmUMZFqaMk+RKkK34tjbxFevPn0rp/SZ/92745Rbse0qBx9q1dd9PRAPEpEGBFR0CqnfFF4ivi5s+fr8mTJ2vWrFkaOXKkZsyYoaKiIq1fv14dO3Y84fglS5boqquu0ujRo5WSkqJHHnlEF198sVavXq3c3FwD7wCIHEodXSTU4L7MbByGkfgCwiiuSh3ryu4vfecZu0fPlnekPRvs28Gd0rGjUvWRmqXZ687i8tizrzyq2Xrr3Oep3ZfqlD2dpPypXlmUp+aw444Nzi4LVEvHDkuV5fb9ezdJS6bbtz6XSd/8o5SYGo7fCoDG2rVaOrLPno3ZeYjpaABzgjO+KHWMuCeeeEITJ07UhAkTJEmzZs3SK6+8ojlz5uiuu+464fhnnnmm3r//8Ic/6B//+IcWL16s8ePHRyVmIFp88bBSfbxoXZPIp9QRQLjE/dK/nQbbNyfwV9urm2xYLH30N2nTEmn9K9KqZ6Wzvmc6OiC+bHnb3nY9m/JjxDdmfEVFVVWVli9frilTpoTu83q9Kiws1NKlSxv1HIcPH9axY8fUrl27SIUJGBM6pyHz5XzBGV9xXupIc3sgjILTgqmHdwBfgt3scfA4afyL0lem2/cv/Y0U8BsNDYg7wcRX/jlm4wBM89HcPhrKysrk9/uVnZ1d7/7s7GyVlJQ06jnuvPNOde7cWYWFhQ3+vLKyUuXl5fVugFN4KHV0j1CpI6s6AgiT4NURxggHGnqtlNLGLn1cv9B0NED8CARqE1/dv2Q2FsC04IxHZnzFtIcffljz5s3TCy+8oJSUlAaPmT59ujIzM0O3vLy8KEcJNF/ctm9xo1CpI4kvAGFSWw/PKOE4ya1rSxzf/bXZWIB4UvqpXXac1No5pdJApCQw4ysasrKy5PP5VFpaWu/+0tJS5eTknPKxjz32mB5++GG99tprGjRo0EmPmzJlig4cOBC6bd26NSyxA9EQPKeJmwW73IxSR0kkvoCw8lDq6GwjbrIbC299Xyp+33Q0QHygvxdQK9jc3k/iK5KSkpJUUFCgxYsXh+4LBAJavHixRo06+cqyjz76qB588EEtWrRIw4cPP+VrJCcnKyMjo94NcArOaVyEUkdJJL6AsKLU0eHSc6RB37b33/2V2ViAeBHq73Wu2TiAWBCa8UWpY6RNnjxZs2fP1tNPP621a9fqlltuUUVFRWiVx/Hjx9drfv/II49o6tSpmjNnjvLz81VSUqKSkhIdOnTI1FsAIoZSRxcJljoeOyxVVZiNxaBmJb5mzpyp/Px8paSkaOTIkVq2bNlJj509e7bOPfdctW3bVm3btlVhYeEpjwecjFJHFxj1A3u77hVpz0azsQBuFwhIX7xj75P4Amqb2zPjK+LGjRunxx57TPfdd5+GDBmiVatWadGiRaGG98XFxdq5c2fo+KeeekpVVVX65je/qU6dOoVujz32mKm3AESML95XqneTpNZSQqq9H8fljglNfcD8+fM1efJkzZo1SyNHjtSMGTNUVFSk9evXq2PHjiccv2TJEl111VUaPXq0UlJS9Mgjj+jiiy/W6tWrlZubG5Y3AcSK0LRgBgnn6thX6naO9MXb9kyU9j1NRwS4V6i/Vzr9vQBJSqgpdWTGV1RMmjRJkyZNavBnS5YsqffvLVu2RD4gIEZ4WdXRPTweu9zxQLFd7tiuu+mIjGjyjK8nnnhCEydO1IQJE9S/f3/NmjVLaWlpmjNnToPHP/PMM7r11ls1ZMgQ9e3bV3/4wx9CNfSA29QtdaQZpIN1qmlWu3u92TgAt9vyP3vbbZTka/K1OMB9QjO+SHwBMMfDjC93aU2D+yYlvqqqqrR8+XIVFhbWPoHXq8LCQi1durRRz3H48GEdO3ZM7dq1a1qkgAMEE18Sfb4cLau3vS0j8QVEVKi/1zlm4wBiRQLN7QGYV9u+xXAgCI9WNZV5cdzgvkmXV8vKyuT3+0O170HZ2dlat25do57jzjvvVOfOneslz45XWVmpysraAb+8vLwpYQLG+OokvgKWJa88pzgaMatDX3vLjC8gcgL+Ov29SHwBkmpnfFHqCMCgUKkjmS93aM3KjlFd1fHhhx/WvHnz9MILLyglJeWkx02fPl2ZmZmhW15eXhSjBJrPU+e/KPp8OViHPvb2wFap8qDZWAC32vahdPSAlJIp5dDfC5BUu6ojM74AGESpo8u0otSxSYmvrKws+Xw+lZaW1ru/tLRUOTk5p3zsY489pocfflivvfaaBg0adMpjp0yZogMHDoRuW7dubUqYgDGUOrpEWrvaKcFln5mNBXCrDf+xtz0vpL8XEBRMfFUdNhsHgLjm9VDq6CqUOjYt8ZWUlKSCgoJ6jemDjepHjRp10sc9+uijevDBB7Vo0SINHz78tK+TnJysjIyMejfACY4vdYSDBWd97SbxBUTEhtftba+LzMYBxJKMmhXPK3ZJx46ajQVA3PLVZAk4n3EJSh2bXuo4efJkzZ49W08//bTWrl2rW265RRUVFZowYYIkafz48ZoyZUro+EceeURTp07VnDlzlJ+fr5KSEpWUlOjQoUPhexdAjKiT95KfSyTOFkp8Na5/IYAmOLRb2rHS3u918p6fQNxJay8ltbb39xebjQVA3ArN+OJ8xh0odWxac3tJGjdunHbv3q377rtPJSUlGjJkiBYtWhRqeF9cXCyvtzaf9tRTT6mqqkrf/OY36z3PtGnT9NOf/rRl0QMxxltvxpfBQNByNLgHImdjzczxnEFSevapjwXiiccjtc2XSj+V9n8hdehtOiIAcchDqaO7hEodSXw1yaRJkzRp0qQGf7ZkyZJ6/96yZUtzXgJwpODSv5JkMTXY2bJqTjaY8QWEX7C/F7O9gBO16WYnvvZtMR0JgDjlC67qyPmMO7SuSXwdPSBVV9b2k4wjUV3VEXA7L6WO7hGc8bX/C+nYEbOxAG4S8EsbamZ8nUF/L+AEbfPtLYkvAIZ4vazq6CqpbaWkdHt/72azsRhC4gsIIw+lju7RuqOU0kayAtKeDaajAdxjx0rpyF4pOVPqMsJ0NEDsIfEFwLBQqWPAcCAID49HyjrD3o/TFetJfAFhFpz1Ramjw3k89PkCIuHzmtUce5wn+ZrVcQFwt7bd7O2+L8zGASBueSl1dJ9gGxcSXwDCweelGaRrdKDPFxB2G2oSX5Q5Ag2rO+OLk04ABvhobu8+zPgCEE7BqcF+vqw6HzO+gPA6tFvavsLep7E90LA2Xe1t1UHpyD6zsQCIS14PPb5chxlfAMIpNDWYSyTO16GPvSXxBYTH6hckWVLnoVJGZ9PRALEpMVVqnWPv74vPJsQAzPJQ6ug+wfOass/jcjYxiS8gzIJTg+Pw74n7ZNUMEHs3StVVZmMB3OCTv9vbQePMxgHEOhrcAzCI1i0u1La75PFJVYekgztNRxN1JL6AMPNS6ugemV2kpNZSoFrau8l0NICz7dkobftA8nilM79uOhogtoUSXzS4BxB9oVJHMl/ukZAktetu78dhNQuJLyDMmBrsIh5PnXr4+BsggLD65Dl72+MCKT3bbCxArGPGFwCDOJ9xqdB5zedm4zCAxBcQZsGpwRYDhTt07G9vd6w0GwfgZJYlfTzf3qfMETi9tt3sLYkvAAbUljpyPuMqcdzgnsQXEGahUseA4UAQHvlj7O2mJUbDABxt+3K7XDgxTep7melogNgXnPG1n1JHANFXW+poOBCEVxxXspD4AsLMw/K/7tLjfHu7Y5V0eK/JSADn+rimqX3fy6Xk1mZjAZwglPjaKvmrjYYCIP54KXV0J0odAYQLA4XLZHSWOvSVZEmb/2s6GsB5/MekT/9h71PmCDRO6xzJlyxZfql8m+loAMQZLuS7VNYZ9vbgTuloudlYoozEFxBmtT2+DAeC8Olxgb3d+KbZOAAnWv2idLhMatWhdgYlgFPzeqU2Xe19VnYEEGW+UOLLcCAIr9Q2UuuaBYb2xNesLxJfQJjV9vhipHCNnhfa201vktEEmuLYUWnxA/b+yJslX4LZeAAnYWVHAIZ4a7IEzPhyoWC54+74anBP4gsIM5b/daH8MZI3UdpfbDfoBtA4y34nHSiW0jtLZ99mOhrAWVjZEYAhlDq6WLDcMc5WdiTxBYRZ7fK/hgNB+CS1kvJG2vubKHcEGqVij/TW4/b+l6dKSWlm4wGchpUdARjiY1VH9wo1uCfxBaAFvFwhcaee59tb+nwBjfPWo1LlASlnIE3tgeag1BGAIZzPuFicruxI4gsIs1CpI1O+3KVHTZ+vzf9jaXngdHZ+LH3wB3v/4p9JXp/ZeAAnakOpIwAzWKXexYKJr70b7ZW34wSJLyDMWAXFpToPkVLa2DNYdqw0HY3rzZw5U/n5+UpJSdHIkSO1bNmyUx7/3HPPqW/fvkpJSdHAgQO1cOHCKEWKE2xYLM29TApUS2dczEqOQHO16273lzy8R1o+13Q0AOKIl9Yt7pWRKyWm2d/T1v3LdDRRQ+ILCDOmBruU1yd1/5K9//4sqbrSbDwuNn/+fE2ePFnTpk3TihUrNHjwYBUVFWnXrl0NHv/uu+/qqquu0ve+9z2tXLlSV155pa688kp9+umnUY4c+nCO9My3pMpyqdsY6Wu/Mx0R4FzJ6dL5d9r7r/xI+mKp2XgAxI3Q+QyZL/fxeqWh19r7z99sV7PEARJfQJixqqOLDbna3n66QPpDoVS2wWw8LvXEE09o4sSJmjBhgvr3769Zs2YpLS1Nc+bMafD4J598Ul/5ylf04x//WP369dODDz6oYcOG6Te/+U2UI49Th3ZJ7/9O+sNF0r/ukCy/NOg70rUvSGntTEcHONu5P5L6XykFjknzvyvt32o6IgBxgFJHlyt6SOp7ueSvlP52VVxUsySYDgBwm+AVko+3HZBHHsPRIKy8Z6ntebPV+72fKLHkY/lnnauyrpfoaOs8VbbqrOrEDFneBFken+Spf13BitBnoVfBhUptlR6R5zahqqpKy5cv15QpU0L3eb1eFRYWaunShmc7LF26VJMnT653X1FRkV588cVIhhpi7d2syl1RSoI25mPU2O+oDT2XVfM/liVZAXkC1fYJt79aqj4iT+UhqeqgPBW75d2/RZ59m+2bZS/7ZHm8Onbunaoe8/+kgEeqoh8e0GKX/kopezbKW/qJAnMvl7/HBbLa9ZSVkSv5kmR5E+r30TvV34BmDEXJ2X3kadO16Q8E4FiemvOZiiq/3v68zHA0iATPkEd15r4ytSl9T9Vzr9D+nDE6nNlLR9K7K5CQooAn0R5f6pzTROp8JjWzvXoNPicizx1E4gsIs0Sf/QfhF6+uNxwJIqOVsvWgZiT+VqO0Rtmb/mE0mu1d3lNuj35GYwinsrIy+f1+ZWdn17s/Oztb69ata/AxJSUlDR5fUlJy0teprKxUZWVtuWp5eXmzYz62ar5S3pre7Me7wcpAL73kH61X/Gdr92ttpNdeMx0S4Cq5ulEvJk9Vh/1b5F3xp6i+dtVFDylpzG1RfU0AZgXPZ3YfrNR3//i+4WgQKa11o/6atEtDqjYpq/jfxuL4OKVAGvxGRF+DxBcQZhPGdNfv39rE1GBXS9eD1kM6+9gydQ0UKzuwS9mBXUqxjsonv3zyy2NZ8tRcdo/kvL+0pKQIPrt7TZ8+Xffff39Ynstq1VFrAt3C8lwmeGTVu4JnSQrII0se+eXTMflUbflUqUQdUqoOWak6oNbaYmXrCytbGwOdtUttzb0BIA5sVwddUvmILvCtVL6nRPmeEmV79itBfiXUjDuq+e823HqmtQ/7cwKIbf06ZejSgTnatLvCdCiIqHTdaz2mwdWfqJu/WN0CxcoJlCjROqYEVctn+SU1PLqc7Ey3OaPQkVZdmvGopiHxBYTZlUNzdeXQXNNhICouMB2A62RlZcnn86m0tLTe/aWlpcrJyWnwMTk5OU06XpKmTJlSrzyyvLxceXl5zYo5acQE5Q8d36zHAkDTfDvqr5iU6Dv9QQBcJdHn1W+vKTAdBqLmy0ZfvXsUXoPEFwAgZiQlJamgoECLFy/WlVdeKUkKBAJavHixJk2a1OBjRo0apcWLF+uHP/xh6L7XX39do0aNOunrJCcnKzk5OSwxezwepSUxnAIAAACxiG/qAICYMnnyZF133XUaPny4RowYoRkzZqiiokITJkyQJI0fP165ubmaPt3uq3X77bfrvPPO0+OPP67LLrtM8+bN04cffqjf//73Jt8GAAAAgBhA4gsAEFPGjRun3bt367777lNJSYmGDBmiRYsWhRrYFxcXy+utXWFm9OjRevbZZ3Xvvffq7rvv1hlnnKEXX3xRAwYMMPUWAAAAAMQIj2XFfgfu8vJyZWZm6sCBA8rIyDAdDgA4Hn9X6+P3AQDhxd/V+vh9AEB4NeXvqveUPwUAAAAAAAAcisQXAAAAAAAAXInEFwAAAAAAAFyJxBcAAAAAAABcicQXAAAAAAAAXInEFwAAAAAAAFyJxBcAAAAAAABcicQXAAAAAAAAXInEFwAAAAAAAFyJxBcAAAAAAABcKcF0AI1hWZYkqby83HAkAOAOwb+nwb+v8Y5xBgDCi3GmPsYZAAivpowzjkh8HTx4UJKUl5dnOBIAcJeDBw8qMzPTdBjGMc4AQGQwztgYZwAgMhozzngsB1yGCQQC2rFjh9LT0+XxeBr9uPLycuXl5Wnr1q3KyMiIYISxJV7ft8R7573z3hvLsiwdPHhQnTt3ltdL1TvjTNPx3uPvvcfr+5Z474wzLcc403S8d957PL33eH3fUnTGGUfM+PJ6verSpUuzH5+RkRF3Hx4pft+3xHvnvcef5rx3rsDXYpxpPt57/L33eH3fEu+dcab5GGeaj/fOe48n8fq+pciOM1x+AQAAAAAAgCuR+AIAAAAAAIAruTrxlZycrGnTpik5Odl0KFEVr+9b4r3z3nnviK54/v3z3uPvvcfr+5Z47/H63mNBPP/+ee+893gSr+9bis57d0RzewAAAAAAAKCpXD3jCwAAAAAAAPGLxBcAAAAAAABcicQXAAAAAAAAXInEFwAAAAAAAFzJ8YmvmTNnKj8/XykpKRo5cqSWLVt2yuOfe+459e3bVykpKRo4cKAWLlwYpUjDqynve/bs2Tr33HPVtm1btW3bVoWFhaf9PcWypv5/HjRv3jx5PB5deeWVkQ0wgpr63vfv36/bbrtNnTp1UnJysnr37h0Xn3lJmjFjhvr06aPU1FTl5eXpjjvu0NGjR6MUbXi89dZbGjt2rDp37iyPx6MXX3zxtI9ZsmSJhg0bpuTkZPXq1Utz586NeJxuF6/jjMRYE49jDeNMfI0zEmNNLGCcYZyJp3FGit+xhnHG4DhjOdi8efOspKQka86cOdbq1autiRMnWm3atLFKS0sbPP6dd96xfD6f9eijj1pr1qyx7r33XisxMdH65JNPohx5yzT1fV999dXWzJkzrZUrV1pr1661rr/+eiszM9Patm1blCNvuaa+96DNmzdbubm51rnnnmt99atfjU6wYdbU915ZWWkNHz7cuvTSS623337b2rx5s7VkyRJr1apVUY685Zr63p955hkrOTnZeuaZZ6zNmzdbr776qtWpUyfrjjvuiHLkLbNw4ULrnnvusZ5//nlLkvXCCy+c8vhNmzZZaWlp1uTJk601a9ZYv/71ry2fz2ctWrQoOgG7ULyOM5bFWBOPYw3jTPyNM5bFWGMa4wzjTDyNM5YVv2MN44zZccbRia8RI0ZYt912W+jffr/f6ty5szV9+vQGj//2t79tXXbZZfXuGzlypHXzzTdHNM5wa+r7Pl51dbWVnp5uPf3005EKMWKa896rq6ut0aNHW3/4wx+s6667zrGDRFPf+1NPPWX16NHDqqqqilaIEdPU937bbbdZF154Yb37Jk+ebI0ZMyaicUZSYwaJn/zkJ9aZZ55Z775x48ZZRUVFEYzM3eJ1nLEsxpp4HGsYZ+J7nLEsxhoTGGcYZywrfsYZy4rfsYZxxmZqnHFsqWNVVZWWL1+uwsLC0H1er1eFhYVaunRpg49ZunRpveMlqaio6KTHx6LmvO/jHT58WMeOHVO7du0iFWZENPe9P/DAA+rYsaO+973vRSPMiGjOe3/55Zc1atQo3XbbbcrOztaAAQP00EMPye/3RyvssGjOex89erSWL18emj68adMmLVy4UJdeemlUYjbFDX/jYkm8jjMSY008jjWMM4wzjeWWv3OxgHGGcSYoHsYZKX7HGsaZponE37mElgZlSllZmfx+v7Kzs+vdn52drXXr1jX4mJKSkgaPLykpiVic4dac9328O++8U507dz7hwxTrmvPe3377bf3xj3/UqlWrohBh5DTnvW/atElvvPGGrrnmGi1cuFAbNmzQrbfeqmPHjmnatGnRCDssmvPer776apWVlemcc86RZVmqrq7W97//fd19993RCNmYk/2NKy8v15EjR5SammooMmeK13FGYqyJx7GGcYZxprEYa8KHcYZxpi63jzNS/I41jDNNE4lxxrEzvtA8Dz/8sObNm6cXXnhBKSkppsOJqIMHD+raa6/V7NmzlZWVZTqcqAsEAurYsaN+//vfq6CgQOPGjdM999yjWbNmmQ4t4pYsWaKHHnpIv/3tb7VixQo9//zzeuWVV/Tggw+aDg2IC4w18YFxhnEGMIVxJn7E61jDOBNejp3xlZWVJZ/Pp9LS0nr3l5aWKicnp8HH5OTkNOn4WNSc9x302GOP6eGHH9Z//vMfDRo0KJJhRkRT3/vGjRu1ZcsWjR07NnRfIBCQJCUkJGj9+vXq2bNnZIMOk+b8/96pUyclJibK5/OF7uvXr59KSkpUVVWlpKSkiMYcLs1571OnTtW1116rG2+8UZI0cOBAVVRU6KabbtI999wjr9edOf+T/Y3LyMjgCnwzxOs4IzHWxONYwzjDONNYjDXhwzjDOFOX28cZKX7HGsaZponEOOPY31ZSUpIKCgq0ePHi0H2BQECLFy/WqFGjGnzMqFGj6h0vSa+//vpJj49FzXnfkvToo4/qwQcf1KJFizR8+PBohBp2TX3vffv21SeffKJVq1aFbldccYUuuOACrVq1Snl5edEMv0Wa8//7mDFjtGHDhtDAKEmfffaZOnXq5IgBIqg57/3w4cMnDAbBwdLuqehObvgbF0vidZyRGGvicaxhnGGcaSy3/J2LBYwzjDNB8TDOSPE71jDONE1E/s41uy1+DJg3b56VnJxszZ0711qzZo110003WW3atLFKSkosy7Ksa6+91rrrrrtCx7/zzjtWQkKC9dhjj1lr1661pk2b5sjlf5v6vh9++GErKSnJWrBggbVz587Q7eDBg6beQrM19b0fz8kroDT1vRcXF1vp6enWpEmTrPXr11v/+te/rI4dO1o/+9nPTL2FZmvqe582bZqVnp5u/e3/t3fHLKnFcRyHvcuxN1GCQktLU46+Czdxdw/aaml1kebGoL2lxU3cztgW+AocGgu+dyrupcsFb6HX/3kecPIIvz8H+cAP8dzd5fn5OY+Pj+n1ehkOh7s6wj95eXlJXdep6zqtVivT6TR1XWe1WiVJLi4uMhqNPq5/f/Tv+fl5np6ecnNz4xHzX9TUziRa08TW6EzzOpNoza7pjM40qTNJc1ujM7vtzF4vvpJkNpvl6OgoVVXl7Owsy+Xy473BYJDxePzb9ff39zk+Pk5VVTk5OcnDw8OWJ/4em5y70+mk1Wp9el1eXm5/8G+w6T3/1T5HItn87IvFIv1+P+12O91uN9fX13l7e9vy1N9jk7O/vr7m6uoqvV4vBwcHOTw8zGQyyXq93v7gXzCfz//43X0/63g8zmAw+PSZ09PTVFWVbreb29vbrc9dmqZ2JtGaJrZGZ5rVmURr/gc6ozNN6kzS3NbozO468yMp/HdyAAAAADTS3v7HFwAAAAD8jcUXAAAAAEWy+AIAAACgSBZfAAAAABTJ4gsAAACAIll8AQAAAFAkiy8AAAAAimTxBQAAAECRLL4AAAAAKJLFFwAAAABFsvgCAAAAoEgWXwAAAAAU6ScL/LocMCI8ZwAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fig, axs = plt.subplots(1, 3, figsize=(15, 5))\n", "\n", "axs[0].plot(r, rho, label='initial')\n", "axs[0].plot(r, rho_final, label='final')\n", "axs[0].set_title('Density')\n", "axs[0].legend()\n", "\n", "axs[1].plot(r, u, label='initial')\n", "axs[1].plot(r, u_final, label='final')\n", "axs[1].set_title('Velocity')\n", "axs[1].legend()\n", "\n", "axs[2].plot(r, p, label='initial')\n", "axs[2].plot(r, p_final, label='final')\n", "axs[2].set_title('Pressure')\n", "axs[2].legend()" ] } ], "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.10.18" } }, "nbformat": 4, "nbformat_minor": 2 }