{ "cells": [ { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "LzDhuJ4aSfhl" }, "source": [ "\n", "*This notebook contains course material from [CBE30338](https://jckantor.github.io/CBE30338)\n", "by Jeffrey Kantor (jeff at nd.edu); the content is available [on Github](https://github.com/jckantor/CBE30338.git).\n", "The text is released under the [CC-BY-NC-ND-4.0 license](https://creativecommons.org/licenses/by-nc-nd/4.0/legalcode),\n", "and code is released under the [MIT license](https://opensource.org/licenses/MIT).*" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "5j3CdtFbSfhn" }, "source": [ "\n", "< [Simulation and Optimal Control in Pharmacokinetics](http://nbviewer.jupyter.org/github/jckantor/CBE30338/blob/master/notebooks/07.01-Simulation-and-Optimal-Control-in-Pharmacokinetics.ipynb) | [Contents](toc.ipynb) | [Simulation in Pyomo](http://nbviewer.jupyter.org/github/jckantor/CBE30338/blob/master/notebooks/07.03-Simulation-in-Pyomo.ipynb) >
" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "8C_vr6C0Jf8I" }, "source": [ "# Soft Landing a Rocket" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "IJIPVYupiSrk" }, "source": [ "Landing a rocket on the surface of a planet was once a staple of science fiction, and then realized in the 1960's through multiple manned and unmanned landings on the moon. It's hard to overestimate the degree to which these missions inspired a new generation \n", " \n", "![Eagle In Lunar Orbit - GPN-2000-001210](https://upload.wikimedia.org/wikipedia/commons/thumb/c/cf/Eagle_In_Lunar_Orbit_-_GPN-2000-001210.jpg/256px-Eagle_In_Lunar_Orbit_-_GPN-2000-001210.jpg)\n", "\n", "[NASA Michael Collins [Public domain] via Wikimedia Commons](https://commons.wikimedia.org/wiki/File:Eagle_In_Lunar_Orbit_-_GPN-2000-001210.jpg)\n", "\n", "**Rocket Landing Videos** (these never get old):\n", "\n", "* [Apollo 11 Landing on the Moon, July 20, 1969](https://youtu.be/k_OD2V6fMLQ)\n", "* [SpaceX Falcon Heavy Side Boosters Landing at Kennedy Space Center, February 6, 2018 ](https://youtu.be/u0-pfzKbh2k)\n", "* [Blue Origin, November 24, 2014](https://youtu.be/9pillaOxGCo?t=103)\n", "\n", "Inspired by these examples, this notebook uses Pyomo and a simple model of a rocket to compute a control policy for a soft landing. The parameters used correspond to the descent of the Apollo 11 Lunar Module to the moon on July 20, 1969." ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "r1DkXPM5wif0" }, "source": [ "## Required Installations\n", "\n" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "WEEvIuiuqYXd" }, "source": [ "### Google Colab\n", "\n", "The following cell installs the necessary packages and solvers on Google Colab. This installation must be done for each Google Colab session." ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "colab": {}, "colab_type": "code", "id": "p00rSCoTwTd5" }, "outputs": [], "source": [ "!pip install -q pyomo==5.6.1\n", "!wget -N -q \"https://ampl.com/dl/open/ipopt/ipopt-linux64.zip\"\n", "!unzip -o -q ipopt-linux64\n", "ipopt_executable = '/content/ipopt'" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "-o3F-AGzL8Xi" }, "source": [ " ### MacOS\n", " \n", " On MacOS, replace the above cell with the following text to perform a one-time installation.\n", " \n", " !pip install -q pyomo\n", " !curl https://ampl.com/dl/open/ipopt/ipopt-osx.zip --output ipopt-osx.zip --silent\n", " !unzip -o -q ipopt-osx\n", " ipopt_executable = './ipopt'" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "2r-wdJM9qszE" }, "source": [ "### Windows PC\n", "\n", "Pyomo can be installed by executing\n", "\n", " !pip install -q pyomo\n", " \n", "into a Jupyter notebook cell. An ipopt binary executable is available from [ampl.com](https://ampl.com/products/solvers/open-source/). More specific instructions for a one-tme installation on Windows PC will be forthcoming. " ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "se0qY6Xfrlv5" }, "source": [ "## Python Initializations" ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "colab": {}, "colab_type": "code", "id": "ISnPZ93d5Zvr" }, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n", "from pyomo.environ import *\n", "from pyomo.dae import *\n", "import numpy as np" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "6VFVBpH9pt0M" }, "source": [ "## Version 1: Vertical Dynamics of a Rocket with Constant Mass\n", "\n", "For a rocket with a mass $m$ in vertical flight at altitude $h$, a momentum balance yields the model\n", "\n", "\\begin{align*}\n", "m\\frac{d^2h}{dt^2} & = - m g + v_eu \\\\\n", "\\end{align*}\n", "\n", "where $u$ is the mass flow of propellant and $v_e$ is the velocity of the exhaust relative to the rocket. In this first attempt at modeling and control we will neglect the change in rocket mass due to fuel burn.\n", "\n", "![LM illustration 02-IT](https://upload.wikimedia.org/wikipedia/commons/thumb/1/13/LM_illustration_02-IT.png/256px-LM_illustration_02-IT.png)\n", "\n", "[NASA Marshall Space Flight Center (NASA-MSFC)derivative work: Adert [Public domain]](https://commons.wikimedia.org/wiki/File:LM_illustration_02-IT.png)\n", "\n", "The complete Apollo lunar module was composed of descent and ascent stages, each containing a rocket engine and associated fuel tanks. The descent stage carried the entire assembly to the lunar surface. The total mass $m$ in the above model therefore consists of the dry and fuel masses of both stages. For the purpose of analyzing the descent of the lunar module to the lunar surface, the 'dry' mass consists of the total mass of the ascent stage plus the dry mass of the descent stage. \n", "\n", "The following data is for the [Apollo 11 Lunar Module](https://nssdc.gsfc.nasa.gov/nmc/spacecraft/display.action?id=1969-059C)." ] }, { "cell_type": "code", "execution_count": 0, "metadata": { "colab": {}, "colab_type": "code", "id": "r-GXDIYS4bao" }, "outputs": [], "source": [ "# lunar module\n", "m_ascent_dry = 2445.0 # kg mass of ascent stage without fuel\n", "m_ascent_fuel = 2376.0 # kg mass of ascent stage fuel\n", "m_descent_dry = 2034.0 # kg mass of descent stage without fuel\n", "m_descent_fuel = 8248.0 # kg mass of descent stage fuel\n", "\n", "m_fuel = m_descent_fuel\n", "m_dry = m_ascent_dry + m_ascent_fuel + m_descent_dry\n", "m_total = m_dry + m_fuel\n", "\n", "# descent engine characteristics\n", "v_exhaust = 3050.0 # m/s\n", "u_max = 45050.0/v_exhaust # 45050 newtons / exhaust velocity\n", "\n", "# landing mission specifications\n", "h_initial = 100000.0 # meters\n", "v_initial = 1520 # orbital velocity m/s\n", "g = 1.62 # m/s**2" ] }, { "cell_type": "markdown", "metadata": { "colab_type": "text", "id": "xCL3kAIJ0hDU" }, "source": [ "### First attempt at a solution\n", "\n", "For this first attempt at a solution, we will choose an arbitrary value for the length of the landing mission. The integration will start with the initial conditions, and we'll see what happens." ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "colab": { "base_uri": "https://localhost:8080/", "height": 226 }, "colab_type": "code", "executionInfo": { "elapsed": 6600, "status": "ok", "timestamp": 1556833347770, "user": { "displayName": "Jeffrey Kantor", "photoUrl": "https://lh5.googleusercontent.com/-8zK5aAW5RMQ/AAAAAAAAAAI/AAAAAAAAKB0/kssUQyz8DTQ/s64/photo.jpg", "userId": "09038942003589296665" }, "user_tz": 240 }, "id": "O10AnUPfr9GY", "outputId": "254fa7a5-87b7-45d7-b60c-f74d21d30888" }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAADRCAYAAACZ6CZ9AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl8ZFWd9/HPtypLJd3ppDcaeg1L\nC4OMbI0yMiDiBsqIzijK4Mag6KgjuIyPqM/gOqPiuKAOj4CIuKCCyODGIrKJ2NDN1t00CHQL3U3v\n6T2drer3/HFOJTdJJanupFKV5Pd+ve4rdc+9de65dSv1q3PuqXNkZjjnnHOVJlXuAjjnnHOFeIBy\nzjlXkTxAOeecq0geoJxzzlUkD1DOOecqkgco55xzFckDlHNlJuldkv44yPaTJT05gsc7VdLakcrP\nuVLxAOVchZFkkg7Lr5vZvWZ2eGL7XyW9sjylc270eIByzjlXkTxAOTdKJH1C0jOSdkl6XNIbC+xz\nT3z4qKTdkt6SbJKT9ENgPvCruP3jhZrskrUsSXWSrpG0TdLjwAl99p0t6ReSNktaLelDpTh/5/aV\nByjnRs8zwMlAI/BZ4EeSDkruYGanxIdHm9lkM/tZn+1vB54D/iFu/0oRx70EODQurwHemd8gKQX8\nCngUmAO8ArhI0mv24/ycG1EeoJwbJWZ2vZk9b2a5GHieAl48Coc+G/iimbWY2RrgssS2E4CZZvY5\nM+sws1XAlcBbR6Fczg2qqtwFcG6ikPQO4CNAc0yaDMwAsiU+9GxgTWL92cTjBcBsSdsTaWng3hKX\nybkheYBybhRIWkCombwCuN/MspIeAbQf2fWdgmAPUJ84VhqYmdi+HpgHrIjr8xPb1gCrzWzhfpTD\nuZLyJj7nRsckQmDZDCDpPOCoAfbdCBwySF59t/8FyEh6naRq4NNAbWL7z4GLJU2VNBf4t8S2B4Bd\nkv5P7EyRlnSUpF4dKZwrBw9Qzo0CM3sc+G/gfkKA+VvgvgF2/wzwA0nbJZ1dYPt/AZ+O2z9mZjuA\n9wNXAesINapkr77PEpr1VgO3AT9MlCsLnAkcE7dvifk07t+ZOjdy5BMWOuecq0Reg3LOOVeRPEA5\n55yrSB6gnHPOVaSSBShJV0vaJGl5Im2apNslPRX/To3pknSZpKclPSbpuMRz3hn3f0pS8hfwx0ta\nFp9zmSQNdgznnHNjSylrUNcAp/dJ+wRwR/zNxR1xHeAMYGFcLgAuhxBsCMO0vITwi/tLEgHncuA9\nieedPsQxnHPOjSEl7cUnqRn4tZkdFdefBE41s/VxDLK7zOxwSd+Nj69L7pdfzOy9Mf27wF1xudPM\njojp5+T3G+gYQ5V1xowZ1tzcPFKn7pxzbgBLly7dYmYzh9pvtEeSmGVm6+PjDcCs+HgOvYdiWRvT\nBktfWyB9sGMMqrm5mSVLlhR5Gj3au7J85uYVzGmqY+7UeuZMrWNOUx2zpmRIp/ZnkADnnBvfJD07\n9F5lHOrIzExSSX+ENdQxJF1AaFJk/vz5A+02qJY9Hdz++Ea27O7olV6VEgc1ZZjbFILW3Bi45k6t\nZ+7UOg5szFCd9j4qzjk3kNEOUBslHZRoftsU09cRxgrLmxvT1hGa+ZLpd8X0uQX2H+wY/ZjZFcAV\nAIsWLdqvYHlQYx1LPv0q9nZkWbd9L2u3tca/YVm3rZV7n9rMxp3tvZ6XUnhuCFp13UFs7tR65jTV\ncVBThtqq9P4UyTnnxoXRDlA3E+ai+VL8+7+J9A9K+imhQ8SOGGBuBf4z0THi1cDFZtYiaaekE4HF\nwDuAbw1xjJKqq0lz2AGTOeyAyQW3t3dlWb+9rTuIheC1l7Xb97J4dQvrH9lLLhEiJTigobY7YHUH\nr0RtLFPtAcw5N36VLEBJuo5Q+5kRZ/u8hBA0fi7pfMLYYPlxxn4LvBZ4GmgFzgOIgejzwINxv8+Z\nWUt8/H5CT8E64HdxYZBjlFVtVZrmGZNonjGp4PbObI4NO9pC4OoTxB5es43fLltPV653JW/G5Np+\nta+5iRpZfY0PVu+cG7t8LL5o0aJFtj+dJEZLNmds3JmogbXsTTQltvL89jY6srlez5k+qaZfzas7\nkHkAc86ViaSlZrZoqP2K+oSSdBLwiJntkfQ24Djgm2ZWVE8MN3zplJjdVMfspjpOaJ7Wb3suZ2ze\n3d5d88oHrrXb9rJy/U5uX7mRjq7eAWxadwDrCVrJjhyTaj2AOefKp9hPoMuBoyUdDXyUMBz/tcDL\nSlUwt29SKTFrSoZZUzIcv6D/9lzO2LKnnbXb9rKmpXdHjic27OL3Kzf1C2BT66u7g1UyiHkAc86N\nhmI/Ybpil+2zgG+b2ffiPR43RqRS4oCGDAc0ZDhufv/Rn5IBbF2fGthfNu7iD09sor3IADZvWujY\n4QHMOTccxX6C7JJ0MfA24BRJKaC6dMVyo22oAGZmbNnd0a8Jcc0gAaxQE2KyK70HMOfcYIr9hHgL\n8M/A+Wa2QdJ84NLSFctVGknMbKhlZkMtxxYZwNbExwM1IeYD2Lw+wWvetDrmNNVTV+Pd6J2byIbs\nxScpDfzezF4+OkUqj0rvxTfWDXQPbE1La/fvwfoGsBmTa3rd95o3rXdnDv8dmHNj04j14jOzrKSc\npEYz2zEyxXMTTTH3wPr2QlzTEh4vW7eDW1dsoDPb+8tU+CFz7+CVr43NbqqjpsqHknJuLCu2iW83\nsEzS7cCefKKZfagkpXITzlC9ELM5Y9Outp57Xy09fx9es43fLFtPNvFDZgkOnJLp3Xw4rb67SfGg\nxgxVPhaicxWt2AB1Y1ycK4t0ShzUWMdBjYV/B9aVzbFxV3t3rav777ZWFq9u4aY+Q0mF/DKJe2CJ\nWti0OmY1ZEj5aPTOlVVRAcrMfiCpDphvZk+WuEzO7bOqdIo5TeHeVCEdXWEoqTXbWruDV74X4j0F\nBvOtSaeY3ZRh3rT6Xt3n58UmxRmTa4iTODvnSqTYkST+AfgqUAMcLOkYwrh4ry9l4ZwbKTVVKeZP\nr2f+9PqC29s6s706buR7IK5taeXW5zfQsqf3dCqZ6lS859VT65o3tT4GsXoa6/1XGM4NV7FNfJ8h\nTLl+F4CZPSLpkBKVyblRl6lOc+jMyRw6s/Bo9HvauxL3v3qaD9e07GXJs9vY1dbVa/+GTFV3AJuX\nuPeVf+y/AXNuaMX+l3Sa2Y4+TRq5gXZ2bryZVFvF4Qc2cPiBDQW379jbGQNXCFr5GtjqLXu456nN\ntHUWGMh3Wk8Am9dnFA7vgehc8QFqhaR/BtKSFgIfAv5UumI5N7Y01lXTOKeRo+Y09ttmZmzd0xGb\nDvf2uge2vEAX+l49EHs1HYYANmtKhrR34HATQLEB6t+ATwHtwE+AW4HPl6pQzo0nkpgxuZYZkwuP\nwpGfSiUZwNbEKVXuf2Yrv9y5juTv6avTYk5TXXcHjnl9gti0Sd6Bw40PxQao15nZpwhBCgBJbwau\nL0mpnJtAklOpvKTA9vauLM9vb+sOXN1NiAN04KivSceAle/AEQLX/OmhKdHvf7mxoth36sX0D0aF\n0pxzI6y2Ks3BMyZx8ACzMe9u7+q595UIYmu3tXL/M1vZ05Httf+0STWh92G89zV/Wk8tzEfgcJVk\n0AAl6QzCVOxzJF2W2DQF6Cr8rMog6XTgm0AauMrMvlTmIjlXEpNrqzjiwCkcceCUftvMjJY9Hb2a\nDvPBa8W6HdzW5/5XSnBQY113h41k8Jo/rZ6ZDbXefOhGzVA1qOeBJcDrgaWJ9F3Ah0tVqOGKA9x+\nB3gVsBZ4UNLNZvZ4eUvm3OiSxPTJtUyfXMsx85r6bc/mjA35+1/Je2Atrdxb4AfMtVWp3sErfw9s\nWmhKnJLx33+5kTNogDKzR4FHJf0k7jtWRpJ4MfC0ma0CkPRT4CzAA5RzCemUukfgOPGQ6f22t3Vm\nE7/5yi97ea6llaUFfv/VWFfdU+vq1YTo3efdviv2HtTpjK2RJOYAaxLra6Hg/Wfn3CAy1WkOO2Ay\nhx1Q+AfMO1o7ea676bA1Pt7LE+t38fvHN9GR7fn9lwQHTcl017bmJ5sQp9Uzc7I3H7rehjOSxMEl\nKtOokXQBcAHA/Pnzy1wa58aexvpq/ra+kb+d2//3X7mcsXFXG89tDUHruZbWQZsPM9Wp7u7y8+OI\nG/na17xp9Uz23ocTznBGkhh8psPyWgfMS6zPjWm9mNkVwBUQJiwcnaI5NzGkEiPQF2q+6G4+zNe8\nYk3suZa9PLC6hd3tvZsP86NvhJpXXa/mQ58+ZXwaryNJPAgsjLW8dcBbCVPWO+cqxGDNh2bG9kTz\n4XOJ+1+PrtnOb/vM/1WVEnMS4x3O79OE2FhX7c2HY9C4HEnCzLokfZBQzjRwtZmtKHOxnHNFksTU\nSTVMnVTD0QV6H3Zlc6zf0dZd+8rf+3qupZVblq9nW2tnr/0bMlW9gtbcxGPvvFG5ZDZ0y5akRYQA\n1UxPUDMze1Hpija6Fi1aZEuWLCl3MZxzI2BXW2f3iBv5IPbs1p5pVDq6enfemN1Y1+v3XvOn99TE\npvvQUSNO0lIzWzTUfsXWoH4MfAxYjo9i7pyrcA2Zao6cXc2Rs/v/eDmXMzbtau+ueSU7b9z9l81s\n2tW780Z9Tbr7XlfvpsPQkSNTnR6t05pwig1Qm83sVyUtiXPOjYJUShzYmOHAxgwvPnhav+17O7Ks\n3da71rWmpZVnt+7h3gJTpxw4JdMrgC1I1L585uXhKTZAXSLpKuAOwn0oAMzsxpKUyjnnyqSuJs3C\nWQ0snNV/7i8zY/Pu9th1vpXntvZ0n7/v6S38Ymdb77yq+9a+wqC986dN8tpXEYoNUOcBRwDV9DTx\nGeAByjk3YUjigIYMBzRkWNTcv/bVt+t8sgnxvqe3sLez98C9XvsaXLEB6gQzO7ykJXHOuTFuqK7z\nW3Z3dAesfBPicy17Cta+Ct77mt7zI+baqvFf+yo2QP1J0pE+2Kpzzu0fScxsqGVmQy3HL+g/cWWo\nfcVa19bwg+XnWvYUvPeVHDYqWfNaMH0S86fVM7V+fPzuq9gAdSLwiKTVhHtQYpx1M3fOuXIKta8G\nDjtggHtfiZ6Hz27tqYXd9ZfNbO7T87ChtqpA8Arrs5vqqB4jo27sy2CxzjnnykASB0zJcMCUwve+\nWju6WNOyl2e37gk/WN66h+daWvnLpl384Yneg/aGGZwzLJg2qbvJcMG0niDWUEFTphQVoMzs2VIX\nxDnn3P6pr6ni8AMbOPzA/rWvXJzz69nunoexFtbSyi3LN9Cyp6PX/tMm1YRg1ee+14Lp9cxqyJBK\njV7ToQ8P7Jxz41gqJWY31TG7qY6/o/+cXzvbOkNzYSJwPbe1lYfXbOM3fcY8rKlKMW9qHQumT+Jb\n5xzLpBKPMO8ByjnnJrApmWpeOLuRF87uP2VKZzbH+u1tPNuyp/u+17NbW9mws436mtL3IvQA5Zxz\nrqDqdCo08U2v5+SFo3/8ogaLnQgkbQaGc69tBrBlhIoz1vi5T1wT+fz93PffAjObOdROHqBGiKQl\nxYzOOx75uU/Mc4eJff5+7qU/97HRGd4559yE4wHKOedcRfIANXKuKHcBysjPfeKayOfv515ifg/K\nOedcRfIalHPOuYrkAWqYJJ0u6UlJT0v6RLnLU0qS5km6U9LjklZIujCmT5N0u6Sn4t/+QzWPE5LS\nkh6W9Ou4frCkxfH6/0xSTbnLWCqSmiTdIOkJSSsl/d1EufaSPhzf88slXScpM56vvaSrJW2StDyR\nVvBaK7gsvg6PSTpupMrhAWoYJKWB7wBnAEcC50g6srylKqku4KNmdiRhhPsPxPP9BHCHmS0kzLo8\nngP1hcDKxPqXga+b2WHANuD8spRqdHwTuMXMjgCOJrwO4/7aS5oDfAhYZGZHAWngrYzva38N/QcJ\nH+hanwEsjMsFwOUjVQgPUMPzYuBpM1tlZh3AT4GzylymkjGz9Wb2UHy8i/ABNYdwzj+Iu/0AeEN5\nSlhakuYCrwOuiusCTgNuiLuM53NvBE4BvgdgZh1mtp0Jcu0Jo+7USaoC6oH1jONrb2b3AC19kge6\n1mcB11rwZ6BJ0kEjUQ4PUMMzB1iTWF8b08Y9Sc3AscBiYJaZrY+bNgCzylSsUvsG8HEgP3fBdGC7\nmXXF9fF8/Q8GNgPfj02cV0maxAS49ma2Dvgq8BwhMO0AljJxrn3eQNe6ZJ+DHqDcPpM0GfgFcJGZ\n7Uxus9AtdNx1DZV0JrDJzJaWuyxlUgUcB1xuZscCe+jTnDeOr/1UQi3hYGA2MIkJPkfeaF1rD1DD\nsw6Yl1ifG9PGLUnVhOD0YzO7MSZvzFfp499N5SpfCZ0EvF7SXwlNuacR7sk0xWYfGN/Xfy2w1swW\nx/UbCAFrIlz7VwKrzWyzmXUCNxLeDxPl2ucNdK1L9jnoAWp4HgQWxt48NYQbpzeXuUwlE++5fA9Y\naWZfS2y6GXhnfPxO4H9Hu2ylZmYXm9lcM2smXOc/mNm5wJ3Am+Ju4/LcAcxsA7BG0uEx6RXA40yA\na09o2jtRUn38H8if+4S49gkDXeubgXfE3nwnAjsSTYHD4j/UHSZJryXcm0gDV5vZF8tcpJKR9PfA\nvcAyeu7DfJJwH+rnwHzCiPBnm1nfG6zjhqRTgY+Z2ZmSDiHUqKYBDwNvM7P2cpavVCQdQ+ggUgOs\nAs4jfMkd99de0meBtxB6sj4MvJtwn2VcXntJ1wGnEkYt3whcAtxEgWsdg/a3Cc2ercB5ZrZkRMrh\nAco551wl8iY+55xzFckDlHPOuYrkAco551xF8gDlnHOuInmAcs45V5E8QDlXpDia9/sT67Ml3TDY\nc4ZxrGpJD5Ui7+GQ1Jwc4dq5UvIA5VzxmoDuAGVmz5vZmwbZfzj+HrivRHk7NyZ4gHKueF8CDpX0\niKRLk7UJSe+SdFOcJ+evkj4o6SNxYNU/S5oW9ztU0i2Slkq6V9IRAxzrdOB3yYQ4F9U1cU6iZZI+\nPFiekmZJ+qWkR+Py0pj+kZjHckkXxbRmhTmerozzHt0mqS5uOz6fB/CBRHleKOmB+Ho8JmnhSL7Y\nzmFmvvjiSxEL0AwsL7QOvAt4GmgAZhJGvH5f3PZ1wsC6EObRWRgfv4QwZFKhYz0A1PdJOx64PbHe\nNFiewM8Sx00DjTGPZYQBTycDKwij0jcTRkk4Ju7/c8LICACPAafEx5cmzvlbwLnxcQ1QV+5r5Mv4\nWvIDHTrnhu9OC/Nk7ZK0A/hVTF8GvCiOAv9S4PowOgwAtX0ziRPktZhZa59Nq4BDJH0L+A1w2xB5\nnga8A8DMssCOOFzVL81sTzzWjcDJhPHUVpvZI/G5S4FmSU2EQHhPTP8hYYI6gPuBT8V5sm40s6eK\nfaGcK4YHKOdGTnIctlxiPUf4X0sR5hA6Zoh8Tgdu7ZtoZtskHQ28BngfcDZwUZF5FiNZ/ixQN9jO\nZvYTSYsJkzj+VtJ7zewPI1AO5wC/B+XcvthFaMLbLxbmzlot6c0QRoePAaevfvef4v4zgJSZ/QL4\nNHDcEHneAfxrTE8rzIp7L/CGODL3JOCNMW2gMm8HtseaF8C5ifIcAqwys8sII1u/qNjXwrlieIBy\nrkhmthW4L3YuuHQ/szkXOD92OFhBmAivm6Q0cJiZPVHguXOAuyQ9AvwIuHiIPC8EXi5pGaHJ7kgz\newi4hnCPazFwlZk9PESZzwO+E4+rRPrZwPKYfhRw7VAn79y+8NHMnasgsabyNjN7X7nL4ly5eYBy\nzjlXkbyJzznnXEXyAOWcc64ieYByzjlXkTxAOeecq0geoJxzzlUkD1DOOecqkgco55xzFckDlHPO\nuYrkAco551xF8gDlnHOuInmAcs45V5E8QDk3AiSdLOnJcpdjNEh6o6Q1knZLOjZOcf/KcpfLjT8e\noJwbAWZ2r5kdXu5yjJKvAh80s8lFTNXh3H7zAOWc21cLCPNOOVdSHqDchCJptqRfSNosabWkDyW2\nfUbSzyVdK2mXpBWSFiW2Hyfp4bjtekk/k/SFuO1USWsT+/5V0sckPSZpR9w3k9h+pqRHJG2X9CdJ\nA85GK8kkvV/SU/HYn5d0aHzezljmmrjvVEm/jue3LT6em8jrXZJWxXxWSzo3ph8m6e5Y1i2Sflag\nHLWSdgNp4FFJzwywzzckPR+Xb0iqjdvulvRP8fFJ8bxeF9dfESc+dK6bByg3YUhKAb8CHiXMTvsK\n4CJJr0ns9nrgp0ATcDPw7fjcGuCXhNlopwHXEaZLH8zZhOnbDyZMh/6umNexwNXAe4HpwHeBm/Mf\n5AN4DXA8cCLwceAK4G3APMJstufE/VLA9wm1nPnA3sQ5TAIuA84wswbgpUA+KHweuA2YCswFvtW3\nAGbWbmaT4+rRZnZogXJ+KpbxGOBo4MWE6ekB7gZOjY9fBqwCTkms3z3I+bsJyANUgqSrJW2StHyE\n8rslfkP+dZ/070l6NH67vkHS5IHyKDbPPvt8RNLjMf87JC1IbJsv6TZJK+M+zTH93viN/pH4zfem\nmC5Jl0l6OuZ3XCKvr8Raxsq4jyTVS/qNpCfiti/1KdvZ8bgrJP0kkZ5NHP/mRPorJD0U0/8o6bCh\nznEQJwAzzexzZtZhZquAK4G3xu1HAl3Ar4FjgR8SPmQhfOhWAZeZWaeZ3UiYNn0wl5nZ82bWQgiM\nx8T0C4DvmtliM8ua2Q+A9niMgXzFzHaa2QpgOXCbma0ysx3A72J5MbOtZvYLM2s1s13AFwkf/nk5\n4ChJdWa2PuYH0EkIarPNrM3M/jjEuQ3kXOBzZrbJzDYDnwXeHrfdnSjLKcB/JdY9QLl+PED1dg3h\nG+9IuZSef86kD5vZ0Wb2IuA54IN9d5B0Vz54FJln0sPAopj/DcBXEtuuBS41s78hfLvdBGBmJ5vZ\nMWZ2DHA/cGPc/wxgYVwuAC6P5XspcBKhZnAU4cM//2HzVTM7gvCheZKkM+JzFgIXAyeZ2QuBixLl\n2ps/vpm9PpF+OXBuLNdP6Pk2Ptg5DmQBMDsG+O2StgOfBGbF7ZuAPwD3xPVWICOpCpgNrLPeU1Cv\nGeJ4GxKPW4H8F5EFwEf7lGNePMZANiYe7y2wPhkgfkH4rqRnJe2M59IkKW1me4C3AO8D1scvEkfE\nPD4OCHggfnn4lyHObSCzgWcT688mzut+4AWSZhGC9bXAPEkzCO/Fe3AuwQNUgpndA7Qk0xTa+m+R\ntDTWMo4Y4OmF8rsD2FUgfWfMW0AdYH332dc8++xzp5m1xtU/E5pskHQkUGVmt8f9dif2I+4zBTgN\nuCkmnQVca8GfCR92B8UyZ4AaoBaoBjbGb+53xvw7gIfyxwfeA3zHzLbF7ZuKOWVgSnzcCDw/2DnG\nc/h3SQ/G2tVnE3mtAVabWVNiaTCz18btW4CdA5RjPTAnXrO8eUWUv5A1wBf7lKPezK7bz/ySPgoc\nDrzEzKbQ04QmADO71cxeBRwEPEGoQWJmG8zsPWY2m9D0+D/52uo+ep4QgPPm03PNWoGlwIXA8vj+\n+BPwEeAZM9uyH8dz45gHqKFdAfybmR0PfAz4n5HIVNL3Cd+wj6BAe/8IOp/QBATwAmC7pBsVbvZf\nKindZ/83AHfkgyjhXk2yprAWmGNm9wN3Ej641wO3mtnKZEaSmoB/AO5IHP8Fku6T9GdJydpqRtKS\nmP6GRPq7gd8qdEB4O9CrybDvOUp6NaG292LCt/TjJeU/pB8Adkn6P5LqJKUlHSXphEIvXB/3A1ng\ng5KqJJ0Vj7E/rgTeJ+klsVl0kqTXSWrYz/ySGgg1qu2SpgGX5DdImiXpLIV7Ue3AbkKTH5LerJ7O\nFNsIXwxy+3H864BPS5oZa0b/Afwosf1uQotBvjnvrj7rznXzADUIhXtDLwWuV+hh9F3CN08k/aOk\n5QWWW4vJ28zOIzR9rCQ0uyDpvPx9GGAR4YP5EUm/3M/yvy3mc2lMqgJOJgTaE4BDiDfuE84hfMgM\nlfdhwN8Qai5zgNMknZzYXhXzuSze68kffyHhRvk5wJUxiAEsMLNFwD8D35CUvwH/YeC1ZjaXcPP/\na0Oc46vj8jCh9nZEPCbArUA9oZlxJ9BBCFpnDHW+8dv+PxKC4XZCB4VfEz7o94mZLSHUJr9NCAZP\n0/867K9vEGrlWwg1y1sS21KE2srzhJaClwH/GredACxW6KV3M3Bh4rrtiy8AS4DHgGWEa/CFxPa7\nCUH0ngHWnethZr4kFqCZ0PwAoWlp/TDzOxX49SDbTym0nfDNsnl/8oz7vJIQ/A5IpJ0I3J1Yfzuh\nyS2/PgPYCmQSad8FzkmsP0kI0v8O/N9E+n8AH0+sX00ITsky/T/gvMT6HcAJBcp+DfAmYCah6Sef\nPh94fIhz/G/gvcO8ZncR7m8Ntd/i5Pn44osvI7t4DWoQFpq5Vkt6M3T3aDt6iKcNKuaR74kmQrfm\nJ4Zd2N7HOJYQWF5vve/zPEi4hzQzrp8GPJ7Y/iZC4GtLpN0MvCOW+0Rgh5mtJ3TueFls7qomfBtf\nGY//BcL9omQnCAj3tU6N+8wgNPmtUvjtTm0i/aRYrm1Ao6QXxOe/KnGMgc7xVuBfYu0XSXMkHVDM\n6zYUSS+TdGA853cSOojcMtTznHP7qdwRspIWQpPUekKX27WE5pyDCR9CjxI+NP9jH/K7F9hMuCew\nlvBblhRwH6H5YznwY2BKgefeRYEaVKE8Y/rnCB/WAL8n9PJ6JC43J57/KnqaX64Bavoc8/Q+xxPw\nHeCZ+JxFMT1NCBAr4+vytZg+l3D/YmXi+O9O5PW1uP8y4K0x/aVx/dH49/zE8d+Y2HYXcEgR53hh\nfM4ywr2jQ4u8Xm+Mr2l7zPvWPtsviOm742v4unK/Z33xZTwvMiu6A5lzzjk3aryJzznnXEXyAOWc\nc64iVZW7AJVixowZ1tzcXO5iOOfcuLd06dItZjZzqP08QEXNzc0sWbKk3MVwzrlxT9KzQ+/lTXzO\nOecqlAco55xzFckDlHPOuYrkAco551xF8gDlnHOuInmAcs45V5E8QDnnnKtIFR+gJP1nYs4g4sjX\nXxjsOc4558a+ig9QwBlmtj2/YmG68NcOsr9zzrlxYCwEqHR+riAASXVA7SD7O+ecGwfGwlBHPwbu\nkPT9uH4e8IMylsc559woqPhTlhJYAAAOd0lEQVQAZWZflvQoYXpvgM+b2a1DPU/S1cCZwCYzO6qU\nZXTOOTfyKj5ARSuBLjP7vaR6SQ1mtmuI51wDfBu4tuSlc845N+Iq/h6UpPcANxCmFweYA9w01PPM\n7B6gpYRFc845V0IVH6CADwAnATsBzOwp4ICylsg551zJjYUA1W5mHfkVSVWAjUTGki6QtETSks2b\nN49Els4550bIWAhQd0v6JFAn6VXA9cCvRiJjM7vCzBaZ2aKZM4ec3NE559woGgsB6hPAZmAZ8F7g\nt8Cny1oi55xzJVfxAcrMcmZ2pZm9GbgAWGxmQzbxSboOuB84XNJaSeeXuqzOOedGTsV3M5d0F/B6\nQlmXApsk/cnMPjzY88zsnFEoHnvau7j8rmeoTqeorhI16RQ1Vamwnk6RTnwFGCysSr33yxmYGWZg\nGLn43JRACAkkhfWYNhIMI5cLN/ly4eDkzDBAQCol0hKpFKQk0ilRnU5RW5UiU52OS4pMVZqqtGjr\nzNHWmWVvZ5a9HeFvW2eWnW1d7Nzbyc69nezY28nOtvAXYHJtFZNqq2iIfydnqqivTpNKhbOUwvmn\nlF8vcB6h6L1fP7Pu1zUXzyuXfJwz9nZmae3I0tqepbUzy96OLva0Z5GgviZNfU0V9TVp6mrSTKqp\nIlOdQrEASl6bAco1mJxBNmd0ZY1sLkdnzrrX0ylIp1JUp8NrXpUS6VTv99dgr0EuvvnM+p9ja3sX\nrR1ZMtUpDmzMcOCUOg5srGXWlAwHNdbRVF9NW2eWPe1Z9naG16O1I0trRxe78texrZOde7vC37ZO\natIpZk3JcMCUDLOm1DKrIcOsKRmmTqomVeCFMQjvk46evPOPO7K57uuTzRlZg1zOMKzf/0Iqvjfy\n55w891yBf8D8/03+/ZROJRaJVErd/4fZ5PslZ935W3xde15zI5vLH7/344HUVKXCez2/ZMJ7v74m\nXfB/O2dGa0eW3e1d7G7rYk97F7vi42wuF18PkU7R/TilcA4Wy59N/D/sj4te+QLSqZH53BlIxQco\noNHMdkp6N3CtmV0i6bFyFypvd3sXl9/9DNnc/l3kia62KsWUumoa66qZkqlCElt3t7KrrSv887V3\njeprm06J+uo09bUhGNVVp6mvSQOwvbWT1o7wYb63I8ueji5KWbTqtKhKpUinRM5CoOrK5YZ9zL7n\nGAJvms27O1m2bgdbdncMnUkfddVpptRVMSVTTUOmih3ZTpY/v5Mtu9sH/WI20Qz0ZWqkj7Eveea/\nUO2rD71iIekR+mI8kLEQoKokHQScDXyq3IXpa9aUDM/852vJ5ozObI6ObI7OrhydWaOjK9fvW9tQ\nb1Cjby2p51thft/ub259alcjJf9NFEjUWsLxsrmeb5H5xx1dOdq7crR3ZmnrynbXmrqyRqYmTV2s\nVdV117B6Pswy1elBy2JmtHflaO3IJr6x9px/oW/FefkaTUpArxpX+HasfC0wUSOrTqu7RjQUM6Mz\na93lCWn0Wu/3HAb+MJBC8KhOpUgN8s00/+23K2sDnn/+OEpcy/zrMdQ5dnTl2LSrjQ072tiws43t\nrZ3dQayupopJsQZZX1PFlEwVDZlqaqoKV+W6sjm27O5g4842Nu5sY1vrwMEvU52OXwiqYv5hqalK\nxRpAsmYTrqPF2r3lev4Xcmbd1zoVP31TsZaVPOvkK9erhhZf32wutCZI4f8gJeJ7pXdtTYRj9NSc\ne7cy5GsvQ73me9p7vpTll70d2YL7C6irSdOQqWJybTWTM1VMrqliUm2aqnSqfytBrMWlUz3lScUa\nYiUbCwHqc8CtwB/N7EFJhwBPlblM/YRmgfSQH7hu30jqDmqVRhI1VaP/D55KiRSiVC9JTVWKuVPr\nmTu1fth5VaVjs2FjZgRKNn7VVKWoqaph6qSaEckvBElKXsMptYoNUJLOAW4zs+sJXcsBMLNVwD+V\nrWDOOedGRcUGKGA+cL2kauAO4HfAA8X04HPOOTf2VWw3czP7spmdRpic8FHgX4CHJP1E0jskzSpv\nCZ1zzpVSJdegAIijlv8yLkg6EjiDMEr5a8pYNOeccyVU8QFK0nEFkm8CvjnaZXHOOTd6KraJL+F/\ngD8DVwBXEkaHuB54UtKrB3qSpNMlPSnpaUmfGJ2iOuecGyljIUA9DxwbB3U9HjgWWAW8CvhKoSdI\nSgPfITQFHgmcE5sGnXPOjRFjIUC9wMxW5FfM7HHgiNjdfCAvBp42s1Vxqo6fAmeVuJzOOedGUMXf\ngwJWSLqcEGQA3gI8LqkW6BzgOXOANYn1tcBLSldE55xzI20sBKjvACcAF8X1+wi/ieoAXj6cjCVd\nQBghnfnz5w8nK+eccyNsLDTxfYMwosQbzeyNwDrg0xbsHuA564B5ifW5Ma0Xn7DQOecq11gIUG8C\nrpF0uKT3AO8HBuy9Fz0ILJR0sKQa4K3AzSUup3POuRFU8U18ZrYqjst3E/Ac8Goz2zvEc7okfZAw\nyGwauDrZ0cI551zlq9gAJWkZvUfEn0YINovjMPsvGuz5ZvZbwvTwzjnnxqCKDVDAmeUugHPOufKp\n2ABlZs+WuwzOOefKZyx0knDOOTcBeYByzjlXkTxAOeecq0geoJxzzlUk+QzqgaTNwHA6ZswAtoxQ\nccYaP/eJayKfv5/7/ltgZkMO3+MBaoRIWmJmi8pdjnLwc5+Y5w4T+/z93Et/7t7E55xzriJ5gHLO\nOVeRPECNnCvKXYAy8nOfuCby+fu5l5jfg3LOOVeRvAblnHOuInmAGiZJp0t6UtLTkj5R7vKUkqR5\nku6U9LikFZIujOnTJN0u6an4d2q5y1oqktKSHpb067h+sKTF8fr/LM4/Ni5JapJ0g6QnJK2U9HcT\n5dpL+nB8zy+XdJ2kzHi+9pKulrRJ0vJEWsFrreCy+Do8Jum4kSqHB6hhkJQmTEl/BnAkcI6kI8tb\nqpLqAj5qZkcCJwIfiOf7CeAOM1sI3BHXx6sLgZWJ9S8DXzezw4BtwPllKdXo+CZwi5kdARxNeB3G\n/bWXNAf4ELDIzI4iTPvzVsb3tb8GOL1P2kDX+gxgYVwuAC4fqUJ4gBqeFwNPm9kqM+sAfgqcVeYy\nlYyZrTezh+LjXYQPqDmEc/5B3O0HwBvKU8LSkjQXeB1wVVwXcBpwQ9xlPJ97I3AK8D0AM+sws+1M\nkGtPmPmhTlIVUA+sZxxfezO7B2jpkzzQtT4LuNaCPwNNkg4aiXJ4gBqeOcCaxPramDbuSWoGjgUW\nA7PMbH3ctAGYVaZildo3gI8Dubg+HdhuZl1xfTxf/4OBzcD3YxPnVZImMQGuvZmtA75KmNF7PbAD\nWMrEufZ5A13rkn0OeoBy+0zSZOAXwEVmtjO5zUK30HHXNVTSmcAmM1ta7rKUSRVwHHC5mR0L7KFP\nc944vvZTCbWEg4HZwCT6N39NKKN1rT1ADc86YF5ifW5MG7ckVROC04/N7MaYvDFfpY9/N5WrfCV0\nEvB6SX8lNOWeRrgn0xSbfWB8X/+1wFozWxzXbyAErIlw7V8JrDazzWbWCdxIeD9MlGufN9C1Ltnn\noAeo4XkQWBh789QQbpzeXOYylUy85/I9YKWZfS2x6WbgnfHxO4H/He2ylZqZXWxmc82smXCd/2Bm\n5wJ3Am+Ku43Lcwcwsw3AGkmHx6RXAI8zAa49oWnvREn18X8gf+4T4tonDHStbwbeEXvznQjsSDQF\nDov/UHeYJL2WcG8iDVxtZl8sc5FKRtLfA/cCy+i5D/NJwn2onwPzCSPCn21mfW+wjhuSTgU+ZmZn\nSjqEUKOaBjwMvM3M2stZvlKRdAyhg0gNsAo4j/Ald9xfe0mfBd5C6Mn6MPBuwn2WcXntJV0HnEoY\ntXwjcAlwEwWudQza3yY0e7YC55nZkhEphwco55xzlcib+JxzzlUkD1DOOecqkgco55xzFckDlHPO\nuYrkAco551xF8gDlXJHiaN7vT6zPlnTDYM8ZxrGqJT1UiryHQ1JzcoRr50rJA5RzxWsCugOUmT1v\nZm8aZP/h+HvgvhLl7dyY4AHKueJ9CThU0iOSLk3WJiS9S9JNcZ6cv0r6oKSPxIFV/yxpWtzvUEm3\nSFoq6V5JRwxwrNOB3yUT4lxU18Q5iZZJ+vBgeUqaJemXkh6Ny0tj+kdiHsslXRTTmhXmeLoyznt0\nm6S6uO34fB7ABxLleaGkB+Lr8ZikhSP5YjuHmfniiy9FLEAzsLzQOvAu4GmgAZhJGPH6fXHb1wkD\n60KYR2dhfPwSwpBJhY71AFDfJ+144PbEetNgeQI/Sxw3DTTGPJYRBjydDKwgjErfTBgl4Zi4/88J\nIyMAPAacEh9fmjjnbwHnxsc1QF25r5Ev42vJD3TonBu+Oy3Mk7VL0g7gVzF9GfCiOAr8S4Hrw+gw\nANT2zSROkNdiZq19Nq0CDpH0LeA3wG1D5Hka8A4AM8sCO+JwVb80sz3xWDcCJxPGU1ttZo/E5y4F\nmiU1EQLhPTH9h4QJ6gDuBz4V58m60cyeKvaFcq4YHqCcGznJcdhyifUc4X8tRZhD6Jgh8jkduLVv\nopltk3Q08BrgfcDZwEVF5lmMZPmzQN1gO5vZTyQtJkzi+FtJ7zWzP4xAOZwD/B6Uc/tiF6EJb79Y\nmDtrtaQ3QxgdPgacvvrdf4r7zwBSZvYL4NPAcUPkeQfwrzE9rTAr7r3AG+LI3JOAN8a0gcq8Hdge\na14A5ybKcwiwyswuI4xs/aJiXwvniuEByrkimdlW4L7YueDS/czmXOD82OFgBWEivG6S0sBhZvZE\ngefOAe6S9AjwI+DiIfK8EHi5pGWEJrsjzewh4BrCPa7FwFVm9vAQZT4P+E48rhLpZwPLY/pRwLVD\nnbxz+8JHM3eugsSaytvM7H3lLotz5eYByjnnXEXyJj7nnHMVyQOUc865iuQByjnnXEXyAOWcc64i\neYByzjlXkTxAOeecq0geoJxzzlWk/w/9tJNh5PpPYAAAAABJRU5ErkJggg==\n", "text/plain": [ "