{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# 解微分方程" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Populating the interactive namespace from numpy and matplotlib\n" ] } ], "source": [ "%pylab inline" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## 积分求解" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 简单的例子" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "$$\\frac{dy}{dt} = sin(t)$$" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": true }, "outputs": [], "source": [ "def dy_dt(y, t):\n", " return np.sin(t)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "积分求解:" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [], "source": [ "from scipy.integrate import odeint\n", "\n", "t = np.linspace(0, 2*pi, 100)\n", "\n", "result = odeint(dy_dt, 0, t)" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAsgAAAEPCAYAAABfrjLnAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3Xd8zPcfwPHXN4m9zt5KrR9F7NEYQa1qa7VUjVKbqhkr\n1EoQlFaNqk0RatbeMWJFbFLzVBKbXCJLxn1+f0QjicRoxl0u7+fjcQ++3/vc997uzn3f9/m+P5+P\nppRCCCGEEEIIEcXK1AEIIYQQQghhTiRBFkIIIYQQIgZJkIUQQgghhIhBEmQhhBBCCCFikARZCCGE\nEEKIGCRBFkIIIYQQIoZEJciaphXVNO2QpmlXNE27rGnaDwm0m6Np2g1N0y5omlYlMc8phBBCCCFE\ncrJJ5OPDgSFKqfOapmUFPDVN26eU8vq3gaZpnwKllFKlNU2rBSwAaifyeYUQQgghhEgWiepBVko9\nUEqdf/n3QMALKBSn2RfAipdtTgE6TdPyJ+Z5hRBCCCGESC5JVoOsaVpxoApwKs5dhQHvGNs+QJGk\nel4hhBBCCCGSUpIkyC/LKzYAg172JL/WJM62rG8thBBCCCHMUmJrkNE0LR2wEfhDKbUlnia+QNEY\n20Ve7ot7HEmahRBCCCFEslNKxe28jSWxs1howBLgqlLq5wSa/QV0fdm+NmBQSj2Mr6FSSm4muI0f\nP97kMaTlm7z+73Dbvh3l5xd7n58fytUV1b8/ys+PkBDF6QMB/NZgDb2+DaVaNUWmjJGU5AYt7IMY\nNEgxb55i3z7FnaN3icAKpde/ev39/KKOpddHHzNJ9yuF0usxkJ0zf/myZo1iwgRFp06KGpVCyY6B\nooXCadVKMXGiYtva5/h+OxrjszjHdXWN/7XYvt3079N73uSzL69/Wr7J62+627tIbImFHdAZaKhp\n2rmXtxaapvXRNK3Py6R3J3Bb07SbwEKgfyKfUwiR1tjZgaMjGAxR2wYDj4dOZXP4Zww2zqJq6QBy\n5TTSq70/p4q0pVL1DMxxfs7jLsO4qbdhZ3kHfp5goH9/+KS6gQ/WTsNafwtmzIDQ0KjjOjqCszMU\nLx71p6Mj/PNP0uw3GKJuM2aQQ3+Barud6djCwPjx8MdcA6frDsXvlh+HGk6mU+sgQkLg12VZsd3u\nRMHCVnzZIog5n+/j/NfTiPyk2WuvBY6OUa+REEKIpGHqLD5GNq+EaYwfP97UIaRp8vrHsH27Un5+\nsff5+Sm1fbt6pjeodU0Xq76d/FX5nPdU9mxG1aKFUlOnKuW+4Z4KIYNSev2rx/Tv/+pY/27fufPa\n/vHVqyvl6hr/844fnzT7XV3fOZ6Y20ajUneO3lV/8I3q9XWAKltWKZ1Oqc+ahanpdTapCzt9lLFf\n/+jXKKHXzlzJZ9+05PU3LXn9TedlzvnmvPRtDVLqJgmy6Rw6dMjUIaRp8vrHECNBNBqVunDUX02t\nvVnVqxOusmVT6tOGQWoWg5XnNl8VHh7nMXr9q+QyoWQxngT20LZtyZ9Evkc8sZLaeP5tDx4o9eef\nSg3o4q8+5KYqXCBc9eql1OY/AlVAzyEJJtvmSD77piWvv2nJ628675Iga+odazGSm6ZpylxiEUIk\nsx07okoCdLpX+wwGIo4cx81Ynz8dz7PjWW0yBPvR8uvsfNomAw1sDWRycgQHh6jSCGfnqMf9W9Kg\n08UulYh57NQo7r8l5jaAoyNquAM3xq5gZ4UR7DyUiRMnFLV012jTKy/t9DMpMHskuLvH+1rj7g4t\nW5rm3yaEECakaRrqLYP0JEEWQqS8GMleRFYdh3cEsn7MeTY/+pgPilvxVeNntHapTenbe9FKFE84\nWaxfH5o1s8zkL4EfEezZA0eOxJs4B9ro2L/mERv67GVH9o7YVrHmq5bBtLs6OSpZtrQfEUIkg6j5\nB4SliC+3lARZCGFaCSR56pg7pzPWZ/ng82x88DEfWPnw1YB8fNU1EyVyvkzgYvYUSy/oKwklzv++\nRi9fu9Cps9nTwJk/d2Zlx3Yjtllu0nlgLtrfmkr2GePkNRUiAS+TJ1OHIZJAQu+lJMhCCNOK01t5\n/29//vjuIMuefUF4pDXdWvnR8adqfKg/GDXrw5vKCqS3883e8NqFZtSxe+UjVvRx51C2VnzR2oru\nXwXSYNcorKY4yWstRAySIFsOSZCFEGYr4omBbV3WsySsM+7u0K69Nd16Z8DuIwPaWOkpTjLv2LP8\neNICVpeewDLXTDz3j+TbfLv57udKFF09TZJjIZAE2ZJIgiyEMK14krOH1/1ZNPkBC93KUrxAKL3O\n9KbdlclkKf+B9BSnpARea+XkzDm9jqU/B7B2VTj2zTMywCELDYN3oNWVHyki7ZIE2XIkJkFO7EIh\nQggR3UOp/AycOAGdvgrjf5XScde6ONvXBHC05jC66ieRZd70V8lWzGRYp3vVgyySVgKvtXbcnaof\nGpibbTT/XA6kSfBf/DAgko+GN2feF3t47i0LkQgh0i7pQRZCJFpEBPy5PIiZY54SkKUg/Qtvpdsf\nTcipU9JTbK7i6VlWYxw5/KkL8xan58DucL7tFMGQ8BkUmzNc3i+RZkgPsuWQEgshRPKLp4wiyNfA\n0sm+zNrzEcWKgUPXh3zasyBW+ttRg+7eVBcrl+tN6y3vjc8Jb3752JWluqF8+pk1DrWOUKlzJXkv\nhcVLzQmyv78/9vb2ZMuWDVdXV7y8vDhw4AD58+enfPnyNGnSxNQhpigpsRBCJL9/B3oZDDx+DONH\nhlKitDVu90qzdi0c3mrgs7OTopLjGTOikqeWLV/vedTpJKEyB296bwwGivwxjRn6r7jVbiQVSobQ\n3LkuLWx9ObQtEKWQ0gshzJCLiwuNGzfGYDBw69YtGjduzNOnT+nbty+NGjUydXipivQgCyHe2YNr\n/szocIZld+z5qvAJhq2oRJnq2WXQnSVJ4L188aMzf6xPz4wfA9AVz8mEIotptrITWk55f4VlSa09\nyEajkSJFirB3714KFy6MtbU1ISEhzJ07FwcHBwCyZ89u4ihTlvQgCyGSzo4dUUlRDA+u+TOszW3K\n18lBRJUaXPYvysIdRaKSY5BBd5Ykgfcywxl3egzMzFXPUIae78Kwa72o00LH7kmnUX6xPy8YDFGf\nIyFEijl16hTW1tZUqFCBnDlz4uTkhJubG1myZOHgwYNpLjlOLOlBFkLEFqMH8UGojhlOoSxbHEmX\nbjaMHBRKobljYs9dLD3Eace/nw0HB4zTZ7KhhgsTp2ckm99dJszNS7N2WdH85eqBSN1Saw/ylClT\nuHjxIq6urqYOxWxID7IQIunodBhGTGF0o1N8VC6SiCMnuHwugl+mhUQlx87OUQPwnJ2ja5JFGhCz\n9KJ4caymONH+zAguuT9n6NR8DOv9nLo1Qjnafakkx0KYwOHDh6lRo4apw7AYkiALkVbFU0oRct/A\njO+8KFMjB49Lf8wFwwf88lcJCpXLIWUUaV0C77/VCXfad8/CRY8X9PXsRRePH/iss45L8468/uNJ\nSi+ESBaRkZGcOHGCKlWqmDoUiyEJshBpVYxZKSIiYOncYMqWMXLicUkObwtgcZ5RFNEfkxkpRJS3\nzHphPWsGXfSTufbZMJrUDeaTiXX5tvbf/HPRP6qtzHohRLI5f/48gYGB2NramuT59Xp9vPvv379P\ncHBwCkeTNCRBFiKtetkDuLPzGmzLh7Fiqi/rNqZj06pgyq0cLaUU4t3EKb3IMG0ig7wduHEmgA++\nqEzVmtYM6+mPYbiTlF4IkUxOnDhB/vz5yZ07d4o/9+3btzl58mS89+XNm5fp06encERJQxJkISxd\nPKUUGAxcXXCYFh11DLnaC5cbbXE7lo46TbNJKYV4Pwl8XrJfcmfS9IxcOerH8yXrKLvFhQVrdUT8\ntVNKL4RIYqdOnaJChQomee6FCxfSsWPHeO+zsbGhZcuWrFy5MoWjSjxJkIWwdDFKKQCe3vZnYP3z\nNBhXj2b1Q7jUdBif6eeizZRSCvEfvKX0osDyafyub8qextNZvyacKqOasb/ryldJspReCJFoJ0+e\nTNIEedSoUezdu/et7S5cuECRIkXe2KZGjRrs378/wfsPHTpEoUKF8Pb2fu84k5MkyEJYupc9euGj\nxjFn/FPKVbDCWLMOXqeeM9hnOOmnTZJSCpH04pReVF7Yj4MVBzNpTCh9Lg3giyp3uXHIR6aEE6lD\nAlfi3uvKR1IcIx7Pnj3j1q1bVKxYMVHHiWnatGk0bdr0re22b98e7wp9jRs3JiIiIno7b9683Lx5\nM95j1KtXj/z581O0aNF4H2sqkiALkQYcuaijqttPbJt0loObA5i3OAN5/j4mpRQi+cRTeqFNcaZN\nTjeu/m2NXYci1GmUkbGaM8HpJTkWZi7Olbj/dOUjKY4RD09PT4AkTZDflYeHB+XLl4+1z9fXF6UU\nNjY20ftsbW2j44zL09Mzenq6+B5rKpIgC2EpElgBr0sjHzp1NDK++Er23i5Nhb+mSCmFSH5v+Hxl\nCDEw8vk4LpwI4daem5T/n5EtY8/IinzCfP3bgeDoCHfu/LcrH0lxjHh4enpiZWX13iUWAQEBzJ07\nl507dzJr1iwgqjd63bp1tG/fPvrY8+fPZ+zYsWzZsoWNGzfy3XffRR8jODgYTXu13sa+ffsYMmQI\nBQoUYNWqVdH7c+bMiY+PT/T29evXGTduHLt372bq1Kk0btw4wceajFLKLG5RoQgh/jM/P6X691fK\nz0+Fhyv1y7QglSdjgBrRx6Ce9xoSdX+cdkKkuLifPz8/tf+LX9T/SoapT4tdUjfP+sffTogU8sZ8\nRK9XCqL+/K+S4hgxtG/fXpUrV+69H7dixQo1YsQI5efnp3r06KGUUmrfvn3Kz89PVa9eXSml1K5d\nu9SBAwdU69atlVJKGY1G9eGHH0Yfo1GjRq8dt2PHjurMmTOx9u3bt09NmTJFKaVUYGCgsrW1VX4v\n/283bNhQPXr0KMHHJkZC7+XL/W/MS6UHWQhL8bJ34mSPRVSvFMaWX7054qZw+fwYWaf/KKUUwjzE\nU3rReEVXLszcR4PvSlLrYysmDvHjxajxUpsszIvBEDUvvF7/an54UxwjjgsXLlCtWrX3flyLFi14\n8uQJFStWjH78J598wvLly+nWrRsAzZs3Z9++fXTp0gWImk4u5mp9cUshlFKcO3futXj8/f3JlSsX\nAJs2baJixYrodDpCQ0MJDAwkb968CT7WVCRBFsJCBATA92N1tDk6hJFe33LgaAbK1coupRTCvCTw\neUzf+lNGjM/EuQN+nPv5MJX3z+DIRUmOhZmIM+j0Pw1qTopjxBESEsLNmzffO6k8deoUjo6OLFmy\nBE9PTw4fPhx939q1a+ncuTM7XpY3HTp0iMaNGwOwcuVKevXqxe7duwEoUKAAgYGB0Y+9evUq5cqV\nA8DV1TV6//379ylVqhQAT548iV7QZN++fdSuXZvdu3fj5eUV72NNRRJkIVKbeGqNt64J5KOSIYQG\nvODK56PpqJ/6ato2IVILg4Giq6exRV+ZKaWX801HI72b38Xvjv9r7aQ2WaSopJgfPhnmmPfy8sJo\nNMbq1X0X+fLlo1q1avz111+sWbOGn376Kfq+Dz/8kO3bt1OrVi2Cg4PR6XTkyJEDgCxZsvDo0aPo\n3uAGDRpw+vTp6Mfmzp2bHDlysHbtWho0aBC9//z589i9HIzYsWNHfHx82LVrF0+ePCFdunQ8f/6c\nXLlyxftYk3lbDUZK3ZAaZCHeTYzaTF9fpdp+/kKV0T1Qh9bef622U2o4RaoRT22yocdQ1b9LgCqY\n2U+tWxqojMZ42gmRxFJTPrJy5UplY2OjgoODTfL8fn5+ytHR8Y1tQkJC1JAhQ1IootgSei9JiRpk\nTdOWapr2UNO0Swncb69pmr+maede3sYm9jmFSNN0OpSTM7+33oltxUg+ur+fC14ZsM/mKdO2idQr\nnt61HDPHMa/DETZssmbScH+++CQY30HTpTZZiJe8vLywtbUlU6ZMJnl+nU5Hnjx5ePLkSYJtXF1d\n6dOnTwpGlTS0qEQ6EQfQtHpAILBSKfXaJHyaptkDQ5VSX7zlOCqxsQiRFty5Az17QsDjFyy9WI0K\n+u1R9WxCWLCw63eYWnYZ83L/yLTp1nTvDjFmlxIiyWiaRmrJRz7//HNKlSrF7NmzTRaDUorFixfT\nq1ev1+7z9vbm7NmztGrVygSRJfxevtz/xm+QRPcgK6WOAn5vaSZfY0K8rzi1xkYjzJsRTHXbMJrW\nC+H4x8OjkuMkGgkthNkyGEj/ywzG67uzr7ELc3+JoEX1R9y9JLXJIm27dOlSdG2vqWiaFm9yDFC0\naFGTJceJlRKD9BTwsaZpFzRN26lpWvm3PkIIEWvVpVu3oFH9CP746QHHNj1mxKPh2EydLEtEC8sX\nZ+S/7cL+nKozhHqN0lGtphW//xyMUiTZqmRCpBbPnj3j7t27Jk+QLVWiSywANE0rDmxLoMQiGxCp\nlArWNK0F8ItSqkw87aTEQog4jM8MzGu1l4lX2jGm3BYG/dUY65PuUUlAzBpMgyGqhlOmbhOWZseO\nBD/vl3PWo3trP3Rl8rHkQ2eKzRkutcki0VJLicXevXsZNGgQXl5epg7FbCWmxCLZF7tWSj2P8fdd\nmqbN1zQtl1LqWdy2EyZMiP67vb099vb2yR2eEGbL2xu6d9cRFPQFx/3KUWb1XsidwPzFMq+xsFRv\n+LxXAE4cf8aM0pOo/rcTMxpb07Wr1CaLtMHNzY0mTZqYOoxUwc3NDTc3t/d6TEr0IOcHHimllKZp\nNYH1Sqni8bSTHmSRNsXpIVMKVv4WhMNoG4b8YMTh8QhsRg6LqjWW0ftCvPJvWYWDAxdGrqHL5ZGU\nzP6YhSszka90jtjt5AqLeEeppQe5XLlyLFiwQDoT38Ckg/Q0TVsLHAfKaprmrWnad5qm9dE07d85\nPb4ELmmadh74Gfg6sc8phEWJUWv86BG0/TyMn8Ya2Lf2KaOfSq2xEPGKpzbZo95Q/lclE7aVYfMf\nQbHbSZ2msACHDh2iY8eOLFmyhKxZs0Ynx8+ePcPFxYVly5bh6elp2iAtRJL0ICcF6UEWaZrBwNZv\n1tH3TA+6FTnAhF21yHBGao2FSNAbapOP29Sna/tQ7D7JxJyc48kxc5xceRHvzJx7kH19fbG3tydD\nhgxs2rSJMmWihnTNnj0bOzs7qlatyrfffsvq1atNHKl5SEwPsiTIQphYUBAMGQL7d4ezytseO/1q\nmddYiEQKvPIPwyvsYk+RHqxam466dU0dkUgtzDlBTsjAgQMZMWIERYsWpUWLFuzatcvUIZkFk5ZY\nCCHeUZx5jQE8Dj6nSplAwgLDON9sVFRyLPMaC5E4BgNZ50/nN31zfqm4hC/bGRk7FsK37nz9/5bM\nnSwsgNFoxNraGohK/kTiSYIsREqJUWscGQnOY0P47DOFs2MIy3MOIfuMcVJrLERixalN/mLN15z/\n1BHPk+HYTWrKjQE/v/q/JfXJwkKULVuWhw8fEhoaSvbs2U0djkWQEgshUpLBwJ2BP9H5b0fS3/+H\nlXsKUOTOMak1FiKpJFCbrI65M+9OSyaMNzK10lp6LrVDmykzw4jXpcYSi6dPn7J06VJy5MhBxYoV\nqVOnjqlDMgtSgyxEKrFuHQwcEInD01EMuzUAqw+LmzokIdKUK1fgmy/DKP33Xyw6X5OctsVMHZIw\nM6kxQRbxkxpkIcxcUBD06AFjx0Sys+FMHPQDsPpJao2FSGkfFTZwqr4Dhb5tSuW6WTm26/nbHySE\nSHMkQRYiqcUZjHf+PFSrEknETT1nGzlQfVEfqTUWwhRe1hxndJnInOXZmft7er5sp5jU/hKRT2Xw\nnhDiFUmQhUhqLwfjKT8Dc+ZAk0+MjPtgFSv6nybbjB9f1TvqdFFJsru7aeMVIq1wd49Vc/x5x6x4\nnlG4Xc1Ho0qP8b7sH9VOBu8JkeZJDbIQyeDpbX+6NbrLwxxlWFtxCiXnDpGBQEKYqchIcJkQwi8/\nhfP7zyG0ujBJBu+lYVKDbDlkkJ4QZuTYMfjmG+jQ3B/nRXlJr78uC38IkQqc2HSfr9uF0aZ7TlwW\nZCdDBlNHJExBEmTLIYP0hDCFOLXGRiNMGRfCl5+HsmBGIDPSjYlKjmXhDyHMn8FAnQNOnDunccdN\nj13tCG4tPiQLiwiRRkmCLMR/FWPhj4cPofkn4exe4suZXU9oeWRk9EIFMhhPCDMXY3GRXJWLsdnz\nA7pm3UydUfVZ32GjLCwiRBokJRZCJIbBwMFvV9Dl1AC+K7KP8bvrYHPKXRb+ECI1SWBxkTPLL/P1\nnDo0yXyc2RuLkXHOdKlNTgOkxMJySA2yECYQGRl1rvxtXgQrHzXnE/1iqTUWwsL4+0PvTkFc33Gd\nPw/lpZR9EVOHJJKZJMiWQ2qQhUhhjx5B8+ZwcG8Eni3GRSXHUmsshMXJoQy4FhtBr0nFqNMiB38u\nDzJ1SEKIFCAJshBvE2cw3pEjULWKkZpZr7K/4hAK/jxSao2FsEQva461Kc70H5eb3btg1KBgBja/\nwYuHMnhPpG56vf6N99+/f5/g4OAUiiZ+CcWYErFJgizE27wcjGd8ZmDaNGj/lZHFVRfg/PUlbKZO\nloU/hLBUcRYWqWafDc+L6fF5lom6Ff3RXwiIaieD90Qqc/v2bU6ePPnGNnnz5mX69OkpFNHr3hRj\nSsQmNchCvIOnt/3p2tAbQ95SuJafTNE5DjJQR4g0Sin4ZVoIUyZHsOjnYFlYxMKkhRrkkSNH4uLi\nEr29ZcsWrl69ipWVFYULF6ZLly4AeHh44OXlRdeuXU0eY1zvEpvUIAuRjE6fhmqNclCuaVHcPLNR\ndFIvOREKkYZpGgwenYm/VgfyQ59QRkZOISKrfCeI1OHChQsUKfJqsKm/vz+TJ09mzJgxjBo1ivnz\n5/PkyRMAatSowf79+00eY3ySOzZJkIX4V5xaY6Vg3oxgPmv6glmTg5iZfgzp9DdkMJ4QAgwGau93\nwtNT4/xOXxo3iOD+qv2ysIgwe9u3b6dRo0bR20eOHKF8+fLR27a2thw6dCh6O2/evNy8edOkMf6r\ncePGRERERG8nZ2ySIAvxrxgLfwQGwjdfhrHI5RnHtz6h7ckRsvCHECJKjIVF8lQtxs5zhWgcvpvq\nDva4dVsuC4sIs+bh4RErIfbx8UEX46qoTqfjxo0b0du2trZ4enqaNEYAX19flFLY2NikSGySIAvx\nr5eD7K72+5UalcPI4uXBiUtZKRV4PnZ9oQzGEyJtizN4zzq3jh/31mX59550PPkD01q4Ybx9JzqJ\nlpIsy6JpSXP7r3x9fZk0aRK7du2ievXqhIWF4evry+TJk9mxYwcTJkzg1q1bBAQEMHfuXHbu3Mms\nWbOiHx8cHIwWIwCDwUDGjBmjt9OnT09gYGD0ds6cOfHx8UmWGIF444wb4759+xgyZAgFChRg1apV\niYrtXdm8vYkQaceanToG7R3D9Gc96a4fDwV18a9+p0tgvxDC8iXwndBkbC08ukH7Vs1xL7mPledH\nklOSY4tjyvF7QUFBtGnThl27dpE7d27q169PeHh4rH1WVlbMnDmTOnXq4O3tTefOndm0aVP0MSIj\nI2MdM1u2bDx9+jR6OyQkhPz580dvZ8qUibCwMACmT59OSEhIvLF9++23FC9e/L1iXLBgAVu2bHkt\nzrgxNmnShGXLljFs2DCqVasWb2xJTRJkIYCwMBg2DHbtiGT/Jy7YuoyPqjWW3h8hxHsoktWAW40f\nGWHrRPX6L9iw7TlV6mczdVjCQqxbt47q1auTO3duALJkycLSpUtj7bt69SqZM2emRYsWHD58mIoV\nKzJmzJjoY8QsUQAoWbIkZ86cid5+8uQJVatWjd729/cnV65cAIwYMSJJYwRo3rx5dJyOjo7xxqiU\n4ty5c7GS47ixJTUpsRBpT5zBeN7eUN8ugrsnfDnTaCS2C/tLrbEQ4v29rDlOP20SPy/NzpTZmWna\nDJb0OS2D90SSiIiIoFSpUtHbJ0+eJDg4OHpfSEgIGzdupFWrVjg6OrJkyRI8PT05fPhw9GMKFCgQ\nq4Sifv36sep4z549S+PGjaO379+/H+s5kyrGoUOHcurUKcaOHRsdp5ubW7wxXr16lXLlygHg6ur6\nn2N7HzIPskh7Ygyw2X9GR5fORgaX+AuHH15g1aJZ7B5jgyGq3lDKKYQQb7NjR9SAvBjfIV6nAmjX\nTlEn2xXmHixPpoK6WN9BcoXK/JjzPMjPnz/H2dkZOzs7wsPDKVCgABUqVMDFxYU6depw/vx52rZt\nS6ZMmdi7dy8FCxbk9u3bfPXVVxQuXBiApUuXUrx48VizRKxatYp//vkHo9FIyZIl6dSpU/R9PXv2\nZO7cubHqlJMixvLly6PX6+ONM26MDx48YPTo0TRt2hR7e3sKFiz4TrElZh5kSZBFmmR8ZmBqiyPM\nu/Mpa+rMxX55NzlRCSGSRWAg9OwaxrWjD9m4UePDdVMlOTZj5pwgJwWDwcDMmTNxcnJ6a9vQ0FDG\njBkTa5BfSniXGN8lNlkoRIj34OcHX3TVscvYlDOPimH/c2s5UQkhkk3WrLB2Y3q6989M7Qbp2V71\nR/nOESaj0+nIkydP9GIgb+Lq6kqfPn1SIKrY3iXG5I5NEmSRppw7B9WqQelioRyq5kAh/XFZ+EMI\nkew0fwM/PPmRzX9G0ndQesY5hBJnoL4QKWbQoEFs3rz5jW28vb3JmTMnZcuWTaGoYntTjCkRW6JL\nLDRNWwq0BB4ppSom0GYO0AIIBroppc7F00ZKLETSiacWcNm8YEaMsWbu7Ag6eI54dYlT6gGFEMkp\nznfMw+v+fN3wIelyZ2fN5kzkKZkjdlsZ92BSll5ikZaYusRiGdA8oTs1TfsUKKWUKg30BhYkwXMK\n8WYxVsULDYXe3V4wfWwAh/e8oEN+N1n4QwiRcuIsLJK/TA72XchHleJ+VKscicfB51HtZOU9IcxG\nkgzS0zStOLAtvh5kTdN+Aw4ppda93P4baKCUehinnfQgi6RlMHBn4E98dXEsxYOusPTQh2QrKj3E\nQgjzsfneHY4tAAAgAElEQVSPIPr0NjJ57At6+4xHmyJXskxNepAtR2J6kFNioZDCgHeMbR+gCPAw\n/uZCJEwpRWBYIIZQQ6yb/wt/Al4EEGGMIMIYQaQxEq9rkax/nI56uobUcGrAivuFyfwkM5lsMpE5\nXWaypM9Cnsx5yJclH3kz5yWddTpT//OEEGlMm85Z+KiQD+0a+3GinQvz02cls6mDEkKk2Ep6cbP0\neH+aTZgwIfrv9vb22NvbJ19EwiwppXgQ+ACvJ17ceHoDnwAffJ774BPgg7e/Nz4BPigUuTLlIkeG\nHOgy6qJv2dJnI511Oqw0azxOWXPxnDVtSpzlgyZVeXx0N8F1qhFsFUlIeAjB4cEEhQfxJPgJj4Ie\n8ST4CdnSZyNflnzkz5qfD3J8QMmcJfkw54eUzBX1Z/4s+WOtDS+EEIlmMFBm41ROXhlB79ZefFyr\nKhu3WFOypKkDE8JyuLm5RS9C8q5SqsTCTSnl+nJbSiwEAH4hfpy5d4bzD87j9cQr6vbYi3TW6SiX\npxylc5WmWI5iFMlehKI5ilIkexGKZC9C9gzZXx0kzmC8Z8+gc4dwnnsbWF97FgV/HvlOA/GMyohf\niB+Pgh7xIPABdwx3uO13m1t+t6L/DAkPoXze8tjmt8W2gC22+W2plL8SOTLmeO14QgjxVnG+l5Sf\ngXmt9jL54hcsXhDB5x2zxm4rg/dShJRYWA6TLxTylgT5U+B7pdSnmqbVBn5WStWOp50kyBYsJDwE\nz/ueePh64HEv6vYg8AFVC1alSoEqlM9bnnJ5ylEubznyZM7z7geOcYI5e1vHl20jaa1zw8XhKela\nNk3SVfH8Q/25/OgyFx5e4MKDC1x4eIHLjy6TN0teahauiV1RO+yK2mFbwBYbq5S6OCOESLXimW0H\ng4ETv56hvUs1vu2biYkuGbF+LjPtpCRJkC2HSRNkTdPWAg2APETVFY8H0gEopRa+bDOXqJkugoDu\nSqmz8RxHEmQL8iLiBad9T3NQf5BDdw5x5t4ZyuUtR41CNahZuCY1CtXgf3n+h7WVdeKfzGBgadvt\njLzQkXnVl9N+XbsUO4lEGiO5+ewmJ31O4u7tjru3O97+3tQoXAO7onY0KtGIj4t+THrr9CkSjxDC\nMjy64c/XDR9gXbwYa8pMIO+s0ZIcpxBJkC2HyXuQk4IkyKmbUoprT6+x7do29t7ey0mfk5TNXZaG\nxRvSsERD6hWrR7YM2ZL8eUND4fvv4fjhMDbetKWcfhcUL57kz/M+noU847j3cY7dPcYB/QGuP71O\ngw8a0KxkM5qVakapXKVMGp8QInWIuHmHcaXXsqaQA+s32VCrlqkjShskQbYckiALkwiPDOfo3aNs\nu7aN7Te2ExIewmdlPqNFqRY0KN4AXcYk7O2I51Kk/kIAX7ZTlK6UicV5RpF1zA9Rq+KZ2WXIx0GP\n2X97P3tu7WHvrb1kSpeJz0p/Rrvy7bArapc0vehCCMvyb/mYgwNb+++hl0cvJnx1lX7ORdByJl3p\nmHhdakmQHR0dqVGjBq1btzZ1KGZLEmSRYsIiw9h3ax+uV1zZfn07pXKV4vMyn/N5mc+pXKBy8s3y\nEGcwy871gXTvZmSMQwQ/PB73au5QM18VTynFpUeX2Pr3VjZ6beRB4ANa/6817cq1w764vUw1J4R4\n/XvMYODm9z/TznMMlYznWXioDJkLmf/3XWqVWhLk0aNHM3z4cHLnzp1gG71eT4kSJV7bf//+fXLk\nyEHmzJY9qaAkyCJZRRojOfzPYVwvu7LJaxP/y/M/vq7wNW3LtaVQtkIpF4jBQOTosUxM58TSpYp1\nG2ywizwS7yCX1NKjcvPZTTZ5bWKT1yZuPLtBm/+1oattV+oWq4uVlhQLXQohUp0EBu8FHzhBvw2N\nObf3MRs3KEpvmCrJcTIw9wR57ty5lCxZknnz5jFhwgQ2btzI1KlTmTRpEsOGDSNLliwA3L59m1On\nTtGxY8fXjhEREYGTk1Os6XUtkSTIIllceXSFJeeWsPbyWgplK8TXH31N+4/a84HuA5PE8/gxdGob\nQvixk7ieLkn+GsVMEkdy8fb3xvWyKysurCAoPIgulbrQpVIXSucuberQhBBmQin4feoTxjka+W0B\ntO2bz9QhWRxzTpDXr1+PtbU1TZs2ZfTo0Tg6OjJ58mTmz59Pz549Wbx4cXTbkSNH4uLikuCxPDw8\n8PLyomvXrikRukkkJkGWLioRS2BYIEvOLqHOkjo0/aMpmWwy4fatG569PXGwczBZcnzyJFSrYqRa\n0BH23ShB/uUuUT3FFqRojqI42Dlwqd8lNrXfxPMXz6m7rC4fL/mYRZ6LCAwLNHWIQggT0/wN9PEd\nz86tEQwdYc3wgaGEh5s6KpFS3NzcsLe35/jx49SpU4fw8HBy585NeHg4Njavphe9cOECRYoUeeOx\natSowf79+5M75FRLEmQBwGnf0/T8qydFZxdl2/VtONZz5J/B/+Dc2JmyecqmXCA7dsRKfJWCX6cH\n80WzUOZW+p2pB2thU6p41GVFR0eLS5Ih6pdtlYJVmN18Nj5DfBhTbww7buyg2OxifL/ze648umLq\nEIUQphCj5rj6F4XwPG/D1b9u0ajiY+55+b/edscO08Qpkk2zZs3Yt28fV65c4cGDB+h0OiIjI5kx\nYwZVqlSJbrd9+3YaNWr02uMbN25MRERE9HbevHm5efNmisSe2kiJRRoWFhnGhqsb+OXULzwOekzv\nar3pVrkbBbIWMF1QMU4AgTY6enYN49rRh2wYe4GS39ZNtbXGScHb35tFZxex+OxiSuUqRb/q/Whb\nri0ZbDKYOjQhREqIpzbZ+MzAlH7ezN9elDXrbLD/LKsM3kskcy6xSIijoyP9+/encOHCALRu3ZrN\nmzfHGjjv6+tLly5dOHjwYPS+lStXkiFDBjp06JDiMacEqUEW7+Vh4EMWei7ktzO/UT5veX6o9QMt\nS7c0n+nGDAau9vuVdh4jsUvvwa8HPiJTQfmS/1d4ZDhbr21lvsd8/n7yN9/X/J6+1fuSK1MuU4cm\nhDCRfZsD6do5kh8GGBkZOA6rKU6SHP9HqSlBXrZsGVmzZkUpRfv27aP3N23alL1790Zv79u3j0WL\nFmFjY0OLFi3o0qULANu2beP69esMGzYsxWNPCZIgi3dy9fFVprtPZ+u1rbQv356BtQZSIV8FU4f1\nmj/+gCGDIpn+rCfd9eNNvvCHObv48CKzTsxi67WtdK7YmcG1B1MyV0lThyWEMAGfE950+PguORtV\nYeWfmcklv5n/k9SUICekcePGHDhwINa+b775hmHDhlGtWrXoffv378fDw4PRo0endIgpIjEJss2b\n7hSW4ZTPKaYem8oJnxP8UPMHbv1wyyx7G0NDYdAgcDsYycEm06g4bbxZLvxhTirlr8Ty1su59/we\nv576lVqLa2Ff3B6Hjx2oVUSW3RIizTAYKPLHNNyuOzDqqxNUrdyQ9RusqFnT1IFZHm1i0sz3r8a/\nfxJuZfVuQ8eaNm0a+7mU4ty5c7GSYwB/f39yyS+p+CmlzOIWFYpIKkajUe25uUfZL7dXH8z+QP16\n6lcVFBZk6rBe2b5dKT+/6M2bN5WqXDFcta/9j/LvOfTVfX5+SvXvH6utSNjzF8/VLyd/UcVmF1NN\nVjZRR+4cMXVIQojkFvd70s9PbWq+UOXNHqLmuAQpozFO2+3bTRJmapGa8hFHR0e1ZcsWNWvWLPXg\nwYPo/V27dlXPnz+P3r58+bJq06aNUkqptWvXRu//9ddf1f79+1Mu4BSW0Hv5cv8b81KZxcLCKKXY\ncX0HNRbVYMieIXxX+TtuDLzB9zW/J3M6M1oxx84uehaKTZugTm0jPXUbcR10guwzxr3qMdbponqQ\n3d1NG28qkTV9Vn6o9QM3Bt6g/Uft6ba1Gw1XNOSQ/lCqv2QohEiAu3vsK206HW3WtufEpP0sm/GE\nDm3DCAjg1eA9OzuThiuSTmRkJHXr1uX69evkz58/en+DBg04ffp09Hbu3LnJkSMHa9eupUGDBtH7\nz58/j518HuIlNcgWQinFvtv7+PHQjwSFBzGhwQTalGtj1quxhT0yMLLJWTY/rc/62rOoubi3lFIk\nsQhjBGsurcHpiBP5suRjfIPxfPLhJ8m3JLgQwqyEPjAwuPElDgTX5s+aM6m8sJ98z75FaqhBjrma\n3qhRo3B2dmb8+PHUrl0bAIPBwMyZM3FyckrwGKGhoYwZM4ZZs2alVNgpThYKSeMO6g9Sb1k9Bu0e\nxJDaQ7jQ9wLtyrcz6+RYr4e6n+m4na82Z33zUXNme/nSTgY2VjZ0te2K1wAvBtQYwMBdA2m0shEn\nfU6aOjQhRArIWEDHbzuKMvHOtzTZP4LfXHWYee4n3mL9+vUULFiQunXrUrx4cYoWLYq9vX10cgyg\n0+nIkycPT548SfA4rq6u9OnTJyVCTpXMN4MSb3Xa9zSNVjSi7/a+9Kvej8v9LtOhQgfzSozjLPwB\nsGV1ELWqvKBj6xC2lHYgl/5s1GA8C1z0w1xYW1nTsWJHLve/TOeKnWn/Z3taubbi8qPLpg5NCJGc\nDAaYMYNv9FM41syJBXMj6djgHgF3Da+3k4VFUoW4q+m5u7tjZ2fH3bt3Y7UbNGgQmzdvjvcY3t7e\n5MyZk7JlU3AhsFRGSixSoVvPbjHm4Bjc77ozwX4C3Sp3w8bKTCckiTFhfVhmHSMGvWDr6ue4rgyn\n1j6nV3VzMrF9igqNCGWBxwKmuU+jacmmTLSfyIc5PzR1WEKIpBT3e9VgIGTkBAY/n8zBnSGs/ysT\nVepnk+/fOMy9xGLr1q2EhIRw7949jEYjBQoUIEOGDNSoUYPiMi1qLDIPchrxOOgxk49MZs2lNQyt\nM5TBtQeb18C7hBgM6AfOosNFRwoGXGPZwQ/IdfXYaytCpbWV8cxBwIsAZp+Yza+nf6Vb5W6MrT8W\nXUY5QQphEeJZee/f79k19+wZ9IORSY5h9L33I9oUSY7/Ze4Jsnh3kiBbuJDwEGafnM2sE7P4puI3\njKs/jrxZ8po6rHf2558woF8ko58OZ/DtQWglips4IhHXg8AH/HjoR7Ze28q4+uPoU60P6azTmTos\nIUQyun7Qh/aNn1C6RWkWrcki+fFLkiBbDhmkZ6GUUmy8upHy88vjed+Tkz1PMqfFnFSTHAcHQ58+\nMHpkJDsbzmSIfhDaTKk1NkcFshbg989/Z1+XfWy9tpWKCyqy/fp2OUkIYakMBspsnMrJv3NSwMeD\nKraRnJSxu0JEkx5kM3Xx4UUG7x7Mk+An/NL8FxqWaGjqkN4szqW8K1egw5eRVMrlw2/l57ya21hq\n3cyeUopdN3cxbO8wCmcrzJwWcyift7ypwxJCJJV4apO3dFxHnxPfMnSwwuHHTEQv2JYGS9+kB9ly\nSImFBXka/JRxh8ax0WsjExpMoFe1XuY7AC+ml1+4ysmZJRt1jB5lxKXSGrr3TofWvJnUGqdC4ZHh\nLDizgMlHJtO1UlfG248ne4bspg5LCJFYCdQm3119lG+cypOlfDFWrklH/gxps0NDEmTLIQmyBTAq\nI4s8FzHu0Dg6fNSBiQ0nkitT6lof3e+OP30/uYmXTUXWVZlGuQU/pKkvVUv1MPAhow+MZs+tPbh8\n4kKnip1koREhLFTEEwMTm59giXcTltX8jWarOqe573FJkC2HJMip3Nn7Z+m3ox82VjYsaLmASvkr\nmTqk93b0KHTuDK0aBuCyIj+Z9F4g081YlJM+JxmwcwCZ02Vmbou52BawNXVIQojkcOcObiW60bXQ\nfr7sYMPUqZAhg6mDSjmSIFsOGaSXSvmH+vPDrh9osboFfar14Wj3o+afHMdZ+CM8HMY5hNK+VSjz\npwcyJ8voqORYFv6wOLWL1OZ0z9N0rtiZpn80Zfje4QSGBZo6LCFEUnq5sIi9fjnnm4/mn5th1Crn\nj9epgNfbycIiwoJJgmwCSinWXlpL+fnlCY0I5Wr/q3xX5TvzWgEvIXZ2UTVpBgO3b0N9uwg81t3m\n3J7HtDwyMqpWrXjxqD9fthOWw9rKmj7V+3C532UeBT2iwvwKbL++3dRhCSGSQszBe8WLk+snRzYU\nGcKA3hHUt9dYODs4apnqf9vZ2Zk64mSjaZrcLOCWqM+AuVxGSCslFncMd+i7vS/3A+/zW8vfqFO0\njqlDem/Kz8Dqr7Yw9HwXxvxvMz/89QlWJ9xl4Y80aP/t/fTb0Y/KBSrzS/NfKJStkKlDEkL8V29Y\nWOTvvPX4pqWBYlXysKjIRPLOGp3mapOF5ZAaZDMSaYzk19O/4nTEiWF1hjH84+GpciGGZ8+gb1+4\ncj6M1TdqUlm/RWqN07iQ8BCcjzrz25nfmGg/kX41+qWOqyFCiPfy4todxv1vPX/kG8qipTbS9yFS\nLalBNhOXHl7i46Ufs+XvLRzvcZzR9Uabf3Icp9YYYM/GQCqVCaFI3lA8GzpEJcdSa5zmZUqXCadG\nThzudpg1l9dQf1l9rj25ZuqwhBBJyWAgw5wZTNe3Z22dXxnQz0i/fhC0cffr5wCpTxYWINEJsqZp\nzTVN+1vTtBuapo2M5357TdP8NU079/I2NrHPmVq8iHjBuIPjaLSyET2r9OTgtwcpk7uMqcN6NzFq\njYODYWDvF/TqFsaKOQHMYhgZXSZKrbGI5aN8H3Gk2xHaf9Qeu6V2uBxzIcIYYeqwhBCJFac2ucHy\n7lxo6kCwIYzKI5pwqueiV+eANFCfLNKGRJVYaJpmDVwDPgF8AQ+go1LKK0Ybe2CoUuqLtxzLokos\nPHw96La1G2Vyl2Hep/NSZ22mwcCZXgvpfHYoVa3PM29vGXJeOSa1xuKt9H56em3rhSHUwNJWS81/\ndhYhRMLeUJu8IaQlA/ob6fvBbsauKU+6n2ekuYVFROqT7DXImqbVAcYrpZq/3B4FoJSaFqONPTBM\nKfX5W45lEQnyi4gXTDw8kSXnlvBzs5/5usLXiR5JaQphYTB5MixcEMkvTzvTUT9Vao3Fe1FKsfTc\nUkYdGEX/6v1xrO9Ieuv0pg5LCJHE7t2DHt8E8/Dw36zYlZ+KzQubOiQh3iglapALA94xtn1e7otJ\nAR9rmnZB07SdmqaVT+Rzmi0PXw+q/l6Vv5/8zcW+F+lYsaP5J8fx1BqfP/qcGmUDOO8RzoXPxkYl\nx1JrLN6Tpmn0qNqD833O43nfk5qLanLx4UVThyWESGKFMhvYWd6BAdOK0ahNdqaMCyHir51SmyxS\nNZtEPv5dunzPAkWVUsGaprUAtgDxFuJOmDAh+u/29vbY29snMryUkap7jf+tNXZ2JjyLjqnjQ/h1\nNsx0CqPrrdFoU15eKvu31lgunYn3VDh7YbZ13Mby88tpvLIxQ2sPxcHOARurxH79CCFM7mXNsTbF\nmR46HZ+08KdHy3/Ysq0xK8q5UG7BD1HnjJh1zEKkMDc3N9zc3N7rMYktsagNTIhRYjEaMCqlXN7w\nGD1QTSn1LM7+VFlicfHhRbps7kIJXQl+++w3CmQtYOqQ3p/BwOW+c+l21YE8hlss3lWYInek1lgk\nvbv+d+m+tTtBYUGsaL2CsnnKmjokIURixFOfrPwMLBzrzbh1HzGyzBaGrKqK9SypTRbmIyVqkG2I\nGqTXGLgHnOb1QXr5gUdKKaVpWk1gvVKqeDzHSlUJcqQxkp9O/MSM4zOY0WQG39p+m3p6jWMIC4Op\nU2HunEimPOtLz9uOaCWKmzgqYcmMysgCjwVMODyBsfXGMrDWQJk3WQgLdPs29OgUQvDJiyzZXYQK\nzaQ2WZiHZK9BVkpFAN8De4CrwDqllJemaX00TevzstmXwCVN084DPwNfJ+Y5zYHeT0/DFQ3ZeWMn\nHr086Fa5m/knx/HUGp/a/5yqpZ9z5kQ45z77kV56R7SZUmsskpeVZsWAmgM4/t1x1l9dT5NVTfAJ\n8DF1WEKIJPZhLgMHqjjQw7kkDVtnZ8KoUF5s2SW1ySJVkJX03kPMUfmj7EYxpM6Q1NPzFaP+Kyid\njnEjQlm7/AWzp4XRwWvCq1rjmHVicilMJLMIYwQux1yYc3oOv7b4lfYftTd1SEKIpBDnXOJ71Z9+\nLf/hdvr/saTiL9Ra3EvOOcJkZKnpJPQk+Am9tvVC76fnj7Z/UCFfBVOH9P4MBvZ3XUnv8/2wy3SO\n2TvLkudvqTUWpnfm3hk6bepErcK1+LXFr+TImMPUIQkhEiOB2uT1024zeEVlvs7vxuQ1Jck6f7ok\nxyLFSYKcRPbd2kf3rd3pWKEjTo2cyGCTwdQhvbeHD2HoUHA/HMF83y/4VD9f5jUWZiUoLAiHfQ7s\nurmLla1XUu+DeqYOSQiRDJ4+haG9A3Hb9JQ5CzPSqnd+U4ck0piUmAfZor2IeMGwPcPovrU7y1sv\nZ0bTGeafHMepNTYa4bdZwVQs+4KieUO58qlDVHIs8xoLM5MlfRbmt5zP3BZz6bChA44HHAmPDDd1\nWEKIJJbb2sCKAiNZviYDI0cqWrcM5+6yA1KbLMyKJMgJuPr4KrUW10Jv0HOh7wU++fATU4f0bv6d\n19hg4MIF+LhWBKtm3OfAuqdMCx9Glunjo3qO/53XWJJkYWZalmnJ+b7nOffgHPWW1eO2321ThySE\nSCoxao4bdizAhb8zUu3JbqoOs2fmpwcJf2yI3c7OzrTxijRLSiziUEqx0HMh4w6NY2rjqfSo0sP8\nZ6iII+CugYlfeLLKx54pFV35bmNLrE64S62xSFWMysicU3NwPurMz81+plOlTqYOSQiRWPHUJmMw\ncPPPc/RfW5eHFx8yb4E1dd2cpDZZJBupQX5PfiF+9NrWi1t+t3Bt52reixjE8yVjfGZgldMdRrtW\npnnd50z7syT59Kel1likaufun6Pjxo7ULFyTeZ/OI1uGbKYOSQiRDJSC9fMeMXzgCxq0yonLvKwU\nPh9/Qi2dOyIxpAb5PbjfdafKwioUyV6Ekz1OmndyDLFKKQDOHHqOXflnzD9SgS2rnrM076io5Fhq\njUUqV6VgFTx7e5LBOgNVFlbBw9fD1CEJIZKB5m+gg9dE/r5i5APvY9hWMuLi0YgXo8a/Oo9J6YVI\nIWm+BznSGMnUY1OZe3oui79YzGdlPkvxGP4zg4FHQ6YyJngsO7YbmTojHV2/DsNqXIw5JWWOSWFB\nNlzdQP8d/RlVdxRDag9JdeVPQogExD1XGQzc+n42Q56O5e8bVvz8v4V8OvfTqE4fOZ+JRJISi7e4\n9/wenTZ1QinF6rarKZw99SyDGRoKv/wCM6dH0vXZz/x44UtyVPogwfouuRwlLIXeT8/XG78mb+a8\nLG+9nDyZ85g6JCFEYr3h3LXLqiWDB4RTQn+AGbsqUrF56jlXC/MkJRZvsOfmHqr9Xo2GxRtyoOsB\n802O45m2bfXvQZQtFsypY2G4N3fiJ307ciycHtWuZcvXf1nrdJIcC4tRImcJjnY/Srk85ai6sCpH\n/jli6pCEEIn1hnNXizoGLjUdRsvxNfikbTZ6dn3BvZX7ZVo4kazSXA9yhDGCcQfHseriKla3XU2D\n4g2S/TkTJcZlp8MXdAwfEoHm68NPczNSz22ylFKING3XjV1039qdATUGMKbeGKytrE0dkhAiKcU5\ntxn+8Wdqm1Ms1jdmYKndDN9sR9Yicg4U70dKLOLwCfCh48aOZEmXhVVtVpE3S95kfb6kcsk9gLGd\n73AxohxTy62ivWtbmbZNiJd8A3zpvLkz1po1q9uuJn9WWZVLCIuRQOnFP5vP4rirLod2BjNudCTf\neU8k/bRJkhyLdyIJcgw7b+zku63fMbj2YEbYjcBKM6PqkgS+AG6sP8d4t4YcPAgjez2jn1MhMur/\nlmnbhIgj0hjJxMMTWXJuCX+0+YOGJRqaOiQhRAo489c9HFtd4kbRRkxwSkcn3Q6s60vnkXgzqUEG\nwiPDGblvJH2392VD+w2MqjvKvJJjeG3Ktn8u+tOjzlU+HtOAjz6Cm2cMDHk2Lio5lmnbhHiNtZU1\nkxpOYlmrZXyz6RsmH55MpDHS1GEJIZKTwUD1Pc7s0ZdlebW5LFoQQQWH5vzZYQPGZzItnEgci+5B\n9g3w5euNX5M1fVZWtVll3qPdDQZ8Bs3AhZGsWWdFv+9tGOaYkZza61PfSJ2VEAm79/weHTd2JIN1\nBv5o+wf5suQzdUhCiKQWz7Rwaowjexu74OiUkUif+0yckp7Pzk3GaoqTnC9FLGm6B3nvrb1UX1Sd\nFqVasOObHeaRHMeZkQKIWmJz0SF6OeiotHUSGVb+jtexZzjNzEjOnERdFoqZDOt0Udvu7ikevhCp\nQaFshTjQ9QA1CtWQWS6EsFTxnBu1Kc40y3gYj7M2/OiUgQm9fal8cBauu3VE/hX/+VdmvRAJsbge\n5EhjJJMOT2LxucWsbrsa++L2iQ8uqcT5xXv5eABTu11jz9Nq9O8RxqCnP5J7XH+ZCF2IJLL75m66\nbenG0DpDGf7xcPMrrxJCJL2X51o13IFdA3fi/Lg3j59qjCq2ls6un5E+n1yNTevS3CC9h4EP6bSp\nE0ZlZE27NRTIWiCJoks6ys/A8R5LmBnUlxNHIxnskI7+PV6Q3UXKKIRIDt7+3rTf0J68mfOyovUK\ncmbKaeqQhBDJJYHSi8OfuuD8UwaunXnOsKGK7+45k23Gj3KOTaPSVIJ89J+jdNzYkW6VuzHRfqJp\n50ONZ1aK8McGNszQM9utCs8ehTPon6H0uDqczOVk9TshkltYZBgj941ky7Ut/PnVn1QvVN3UIQkh\nksNbzqent9xjZptjHNC1o9t31gz86CDF21aV828akyZqkJVS/HT8J77880sWfb4Ip0ZOpl8sIMas\nFM+ewbTxIZQooVh4vCKOgwO51mIIA/XDyDxXVr8TIiWkt07P7OazmdFkBi1Wt2CBxwLMpXNACJGE\n3nQ+NRiouc+Z9fqanP18AlpYKNWG2/NV9dsc3/McpZBZL0S0VN2D7B/qT/et3fEJ8OHPr/7kA90H\nyQ7+jfQAAB/xSURBVBRdAhL4paqOueORqT6/D7zIRt/afFHQg8G/f0SVSpEyI4UQJnb96XW++vMr\nKuSrwO+f/U6W9FlMHZIQIrnFU3qBoyPPRzmzbE16fpkSRM4PstMn9wa+XtmSbBdlMS5LZtE9yBce\nXKDa79UomLUgR7sfTfnkGF6bv9j/HwPzW+2myujmdOydjVKfl+fvgEKs2F2AKvWzyYwUQpiBMrnL\ncKLHCdJZpaPW4lpce3LN1CEJIZJbAuffbBfd+WFkZq6fDcLpUit2pm/NB7Y6+m74hLO9f3s184X0\nLKc5qTJBXn5+OZ+s+oRJDScxr+U8MthkSN4nTGB6NtzdMU525kj3ZfTo8JziZdPjlqstM2dbc8PD\nwKjAseTXn3q1uIeUUghhFjKny8yyVssYVGsQ9ZbVY8PVDaYOSQiRnN5SemE9awbN9b+xudQILrv7\nU7RUBtqeHEH1Mv785vyUp8OmvOrQkuni0gallFncokJ5s5DwENXrr16q7K9l1eWHl9/aPsn4+SnV\nv3/Uny+3L3WYrEYNDlHFiilVoewL5YKDenD6nwTbx9oWQpiNM75nVImfS6jBuwarsIgwU4cjhEhJ\nbzhfR0QotXv5fdUeV5U9W6T64gul1i0NVMG9B8n5PZV7mXO+MS9NNT3Iej89dkvtMIQa8OjlwUf5\nPkr6J3lDTzHOzugHzmL6qGfYlgqkxdExGNNnZNvqAC41HsIIfX/yL3eJ1V5KKYQwf9UKVeNM7zPc\neHaDhisa4hvga+qQhBAp5Q3na+vnBpqdnsw6fS28OzjQtnkQi9dmodC6WXSr/Td7Vz4gbNSP0rNs\nqd6WQafUjTf0IO+4vkPlm5FPzT4xWxmNxsT+cEhYnF+Cxmd+6uxXU9SPI0JUpUpK5c0doXryuzq0\n9r6KjHy9vfySFCL1ijRGKucjzqrgzILq4O2Dpg5HCGFKbzi/37un1KyxT1UtTqicOSJUx45RPcv+\nPYdKPpBK8A49yGY9i0WkMZKJhyey9NxSXL90pW6xuknzZG+YJzG4ih1H+/zBzjxd2bIhgnR5c9C6\nrTWtGz+nzl+jsR45/NVKd+4yylUIS7P/9n66bO7C4FqDGWE3Ak1740BnIYQletN8yv8O0Hdw4N6E\n39lmO5YtezPj7q6om/MqX/QuQBOvOZScO0TyBDOVqhcKeRL8hE6bOhEWGcbadmv/26p4CX3A9+yB\nI0fA2ZmIrDo83Z6z3/EQ+9O3wONcOqqWD6WZx2Ra7+lP+SaF0fzjnx5GpmcTwjJ5+3vz1Z9fUTBb\nQZa3Wk6OjDlMHZIQwhwkMF0czs4EWOn+396dx1VZ5n0c/1yioiAqoikuaJqaK4pLGklqpeaoLeQ2\nZZup2WI9ZU9ZllozTjkz5fhULimWS7mbaY2aFpWauaG4gAtK7huKWyIC1/OHVGTgClzA+b5fr/Pi\nnMPtzdfzOuWX6/zu+2bhpMPMf3YRS2/oSTHfwtwZdo47D0yh7YcPULZGqd+3DwuD9u1VnB3JldO8\nGWM6GGNijTHbjTEvZ7HNqPTvbzDGNL7cPlfvW02TcU0ILh/M172+vnw5zmp2+PTpP5yGjcRETgx8\niyXef+Fv/v+iU6M9lAtIpU/3kyQ0bcdLg4pwMDaR75u9yGu7+lBv3vAL5VgzxSIepUqpKnz/2PdU\n8qtE04+aEn0o2nUkEckLLtEHSqYl0i1mGJN3tWJf+HMs+PQk9Rp7M9k8TI26RWlUN5mnbotmUoMR\nbKtxN/bV1zI/jVxWnUbzzLnrcjMYl7oBXsAOoBpQBFgP1Llom47AV+n3bwFWZrEvm5aWZsd8/64t\nN8THzl496Y8DI8ePWztkyJ/neY4ft3batExnhY7tSrTff3nSfhj2me3d7aStV2af9fVNs61aWfvS\nS9bOHn3IHqC8tbt2/eHPaYZIRH41ZcMUW3ZEWTt5w2TXUUQkr7pMf0jetsv+yC32vdcTbLdu1gYF\nWVvGP9XeHbTJDn3+mJ3TYazdtuaETUm5xL6mTcu8A2XVjbJ6fsGC7P7b5ztcwQzy9RbklsDCDI9f\nAV65aJsxQPcMj2OB8pnsyz4yvaetN7iM3brhm8zfHPHxf3r+l77P2dhVJ+yi2afs2Nun2hd6J9p2\nVTbbioGp1s/P2hYtrO3d7aT9P562a77YZ5N/PYvTr/vctev3fS5YoDeTiPxJ9MFoW3NUTfvUgqds\n0vkk13FEJK+5VH/IrG9Ya/fts3bOmEP2FYbbTm3P2GrVrPXxsTYkxNqHe5yz77ScY2e8f8iuCn/b\nHtqWaNOOZVGcM+lGl3xei345f5CeMeYBoL21tk/644eAW6y1z2bYZj7wD2vtivTHS4CXrbVrL9qX\n7fxCAwbfu4jCvoGcO3qKcx9O4ETHniRM+5qE1uEcO1uchAPnSPhhC3v96vLztnMkpvlRpYqhalWo\nGnCKWjP+Rv0JL1D/jvIEBfH7/PBLL/1+cB1oplhErsqJpBM8Ou9RDpw6wMyuM6lSqorrSCKS111i\nZhn4Uz855VWaLVtg0ybY8uMJ4icsIb5+J+L3e5OUBEGVU6l6divlmgQRsHMVAZ1aElCpOGW8zxAw\nL4LivR7Ae/okij7XH+9yJSmadJLEUa9T7+Xn8R71L/WcdDl+kJ4xJhzocAUF+W1r7fL0x0uA/7XW\nrrtoX9a/1GC8inhRuDD4+7embIkWlFr9NQHhbQioWoKAAAgIgDIph6n0zL1UWzmdCs2qUKgQv7/p\nrqQIazheRK6BtZYRy0cw8qeRTLlvCndUv8N1JBHJy67gZAGZLtRl0mlOeZXm559h9+pDHHn8f0l4\n7T2OUYaEBC7c9p0lacU6khs25Zz1JjkZEgMWceS2Rxg3pRy9l8+HatWcvRQuRUZGEhkZ+dvjYcOG\nXbYgX++IRQv+OGIxiAurwxePWPTI8DjLEYtMPwq46COJTJ+/2nkdjUyIyHVYunOprfCvCnb498Nt\nalqq6zgikt9cyUhGZqMRV9iNUo8l2GGRw2zFfwba756798/bezhyYQa5MBDHhYP0inL5g/RacImD\n9K55nkZFWERy2Z4Te2yL8S3sPZ/dY4+f1T86IpJNsirPWZyQ4OJulHBgp+04KMiGftjU7nvmEc0g\nZ+JKCvJ1nwfZGHM3MJILZ7SYYK39hzGmX/rq9Nj0bd4HOgBngMfsReMV6dtcyJyYCCNHwvPP//kj\niaye12iEiDiQnJrMi4teZGHcQmZ3m03D8g1dRxKRgiqrcY0M3SjqQBThM8K5p1oHRkQFUOT5F9WZ\nMpGvLxQiIpJfTI2eyvOLnufddu/SK7iX6zgi4oEioiJ4ecnLfNDxA7rV6+Y6Tp6mgiwikks2Hd5E\n+Ixw2lZry8gOI/Eu7O06koh4gKSUJJ796ll+2P0Dc7rPoW65uq4j5Xm5ciU9ERGB+jfUZ3Wf1Rw6\nc4iwj8PYfWK360giUsDtOr6L0IhQTiafZHWf1SrH2UgFWUQkm5T0LsnsbrPpWrcrzT9qzuK4xa4j\niUgB9eW2L2kxoQUPN3yYaeHT8PP2cx2pQNGIhYhIDvgu/jv+Ouev9GvSj8FhgylktB4hItcvNS2V\nYd8NY+L6iUwLn0ZoUKjrSPmOZpBFRBw6cOoA3Wd1x7eoL1Pum0KAT4DrSCKSjx05c4QH5zzI+bTz\nTAufRvkS5V1Hypc0gywi4lCgXyBLH15K/XL1aTKuCav3rXYdSUTyqR/3/EiTcU0ICQzh615fqxzn\nMK0gi4jkgrkxc+m3oB/DWg/jyaZPYsylr3IqIgIXLug26qdRDF82nPGdx9O5dmfXkfI9jViIiOQh\n2xO2Ez4jnAblGzC201hKFC3hOpKI5GEnz53kiS+eIO54HDO7zqS6f3XXkQoEjViIiOQhNQNqsvKJ\nlXh7edP8o+bEHIlxHUlE8qhNhzfR7KNmlC5WmuWPL1c5zmUqyCIiuciniA8R90Qw8NaBhH0cxqcb\nP3UdSUTymEkbJtHmkza8eturjOs8jmKFi7mO5HE0YiEi4siGgxt4YOYD3FX9Lt5r/56uvifi4c6e\nP8uz/32WZbuXMbPrTBqUb+A6UoGkEQsRkTwsuEIwa/qs4fCZw4RGhLLz+E7XkUTEke0J22kxoQVn\nzp9hdZ/VKseOqSCLiDhUqlgpZnadSa+GvWgxvgVzY+a6jiQiuWzm5pmERoTSv2l/Pr3/U10VLw/Q\niIWISB6xat8qus/qzj2172HEXSMo6lXUdSQRyUHnUs4xcPFAvtrxFTO7ziQkMMR1JI+gEQsRkXyk\neaXmrOu7jvjEeFpNbEV8YrzrSCKSQ+KOxREaEcq+U/tY23etynEeo4IsIpKH+Bf3Z273uXSv151b\nxt/CvNh5riOJSDabuXkmLSe05JHgR5jdbTali5V2HUkuohELEZE8auXelfSY1YN7b76Xd+58R2e5\nEMnnklKSeGHRCyyKW8SMB2bQpGIT15E8kkYsRETysRaVW7Cu3zp+PvEzt0bcyo5jO1xHEpFrtD1h\nOy0ntOToL0dZ13edynEep4IsIpKHlSlehjnd5vBYo8doOaEln238zHUkEblKU6OncmvErfRr0o/p\nD0ynVLFSriPJZWjEQkQkn4g6EEX3Wd0JqxrGqLtH4VPEx3UkEbmE08mneearZ1i5dyXTH5hOcIVg\n15EEjViIiBQojQMbs7bvWpJSkmj2UTM2HtroOpKIZCHqQBRNxjXBy3ixtu9aleN8RgVZRCQf8fP2\nY/J9k3np1pdoO6ktH6z6AH36JpJ3WGsZuXIk7aa0Y+jtQ5lwzwR8i/q6jiVXSSMWIiL51LaEbfx1\n9l+p6FeRiHsiKOtT1nUkEY925MwRHv/icQ6dPsRn4Z9Ro0wN15EkExqxEBEpwGoF1GJF7xXUDqhN\nozGNWLpzqetIIh5rcdxiGo1tRN2ydVn2+DKV43xOK8giIgXA4rjFPDbvMXo17MWbbd7UZapFcklS\nShKDlgxiVswsPrn3E9re2NZ1JLkMrSCLiHiIdjXaEdUvik2HN3HrhFuJPRrrOpJIgbf58GZuGX8L\nu0/uZn2/9SrHBYgKsohIAXGD7w3M7zmfJ0KeoNXEVoxePVoH8InkAGst7696n9s/vp0BzQcwq+ss\nAnwCXMeSbKQRCxGRAmjr0a08NPchbvC9gYguEZQvUd51JJECYf+p/fT+ojdHfznK1PunUiuglutI\ncpVydMTCGFPGGPO1MWabMWaxMaZ0FtvFG2OijTFRxphV1/rzRETkytUuW5sVj68gpEIIjcY24out\nX7iOJJLvzdg8g8ZjG9O8YnNWPL5C5bgAu+YVZGPMCOCotXaEMeZlwN9a+0om2+0Cmlhrj11mf1pB\nFhHJAct3L+fhzx+mddXWvNfhPUp6l3QdSSRfOX72OM/89xnW7F/D5Psm07xSc9eR5Drk9EF6XYBP\n0u9/Atx7qSzX8XNEROQ6hAaFsr7feop4FaHh6IZ8s+sb15FE8o0lO5cQPCaYMsXKENUvSuXYQ1zP\nCvJxa61/+n0DHPv18UXb7QROAKnAWGvtR1nsTyvIIiI5bOGOhfSZ34d7a9/L23e+rSt8iWThdPJp\nXlnyCp/Hfk7EPRG0q9HOdSTJJte9gpw+Y7wxk1uXjNulN9us2m2otbYxcDfwtDGm1dX8JUREJPt0\nuKkD0U9Gk3gukUZjG7FizwrXkUTynMj4SBqObsip5FNs7L9R5dgDFb7UN621d2X1PWPMIWNMBWvt\nQWNMIHA4i30cSP96xBgzF2gO/JDZtkOHDv3tfuvWrWnduvXl8ouIyFXyL+7P5PsmMydmDuEzwnmw\nwYO82eZNfIr4uI4m4tTp5NMMWjKIubFzGdNpDJ1qdXIdSbJBZGQkkZGRV/VnrvcgvQRr7TvGmFeA\n0hcfpGeM8QG8rLWnjDG+wGJgmLV2cSb704iFiEguO3LmCAMWDmDN/jVM6DKBsKphriOJOPFd/Hc8\n/sXj3BZ0GyPbj8S/+J+mRqWAuJIRi+spyGWAGUAQEA90s9YmGmMqAh9Za/9ijKkOzEn/I4WBqdba\nf2SxPxVkERFH5sXO46mvnvptNtnP2891JJFccfLcSQYtGcS8rfMY/ZfRdK7d2XUkyWE5WpCzmwqy\niIhbx88eZ+DigSzZtYRxncbR/qb2riOJ5Kj5W+fz9FdP075Ge0bcNUKrxh5CBVlERK7a4rjF9J3f\nl7CqYfy73b8p51vOdSSRbHXo9CEGLBzAugPrGNdpHG1ubOM6kuSinD4PsoiIFEDtarRj01ObKOdT\njvqj6zMxaiJawJCCwFrLxKiJNBjdgOqlqxP9ZLTKsWRKK8giIpKldQfW0W9BP3yL+DK201hql63t\nOpLINYk5EsNTXz3FqXOnGN9lPI0qNHIdSRzRCrKIiFyXkMAQVvZeyf117ic0IpShkUNJSklyHUvk\niv1y/hcGLRlE2MdhhNcJ56cnflI5lstSQRYRkUvyKuTFgFsGsP7J9Ww4tIEGoxvw3+3/dR1L5LK+\n2PoFdT+oy88nfib6yWieaf4MXoW8XMeSfEAjFiIiclW+2v4Vzy18jnrl6vFe+/e40f9G15FE/iA+\nMZ7nFj5H7NFYPuz4IXdUv8N1JMlDNGIhIiLZrmPNjmzqv4nmlZrT9KOmDI0cytnzZ13HEuFM8hne\n+PYNmoxrQrOKzYh+MlrlWK6JCrKIiFw178LevNrqVaL6RbH5yGbqfliXuTFzdbYLccJay6cbP+Xm\nD25mx7EdrO+3nsFhg/Eu7O06muRTGrEQEZHrtmTnEp5f+Dxlfcrybvt3CQkMcR1JPMTa/WsZsHAA\nSSlJ/KfDf7gt6DbXkSSP04VCREQk16SkpRARFcGQyCG0r9Gev7f9O5VKVnIdSwqovSf38vq3r7Nw\nx0L+1uZvPNroUR2AJ1dEM8giIpJrChcqTN8mfdn6zFYq+lWk4ZiGDI0cypnkM66jSQFyIukEg5YM\nInhMMIElAol9OpbeIb1VjiVbqSCLiEi2KuldkuF3DGdd33VsTdhKzf+ryQerPiA5Ndl1NMnHzqWc\n4z8r/0Ot92tx+MxhNjy5geF3DKdUsVKuo0kBpBELERHJUWv3r+W1b15jW8I23mzzJj3r99Rqn1yx\n1LRUpm+ezuBvBlOnXB3evuNtGpRv4DqW5GOaQRYRkTzju/jvGLR0EKeSTzG87XA61eqEMZf8N0o8\nWJpNY07MHIZEDsGvqB/D7xhO2xvbuo4lBYAKsoiI5CnWWuZvm89r37yGbxFf3rj9De6+6W4VZfmN\ntZZ5W+cxJHIIRb2K8mbrN+lwUwe9RyTbqCCLiEielJqWyuyY2bz1/Vt4e3nzetjrdKndRSXIg1lr\n+XL7lwyJHEJqWipvtnmTzrU66z0h2U4FWURE8rQ0m8bnsZ/z1vdvYa1lcNhg7q9zP4WMjiH3FClp\nKczYPIO3l72NMYY3wt7gvjr36T0gOUYFWURE8gVrLQu2LeCt79/idPJpXmz5Ig82fJBihYu5jiY5\nJCkliYlRE/nnin9SuWRlBt02SKMUkitUkEVEJF+x1rJ011L+/eO/WX9wPU83e5r+TfsT4BPgOppk\nk4RfEhi3dhyjVo2iacWmvBL6CqFBoa5jiQdRQRYRkXxr8+HNvPvju8yJnUPP+j35nxb/Q82Amq5j\nyTXaeGgjo34axayYWXSp3YWBLQfqdG3ihAqyiIjkewdPH+T9Ve8zdu1YQgJD6N+0P51qdaJwocKu\no8llpKalMn/bfEb9NIrYo7H0b9qffk37cYPvDa6jiQdTQRYRkQIjKSWJWVtmMXrNaHaf2E2fkD48\nEfIEFf0quo4mF9lzYg8fr/+YCVETCPQLZEDzAYTXDaeoV1HX0URUkEVEpGDacHADY9aMYfrm6bS5\nsQ2PBj9Kh5s6UMSriOtoHis5NZkF2xYwft14Vu5dSY/6PejduDdNKjZxHU3kD1SQRUSkQDt57iTT\nNk1j0oZJbD+2nZ71e/Jw8MM0rtBYZ0PIBdZa1h9cz9SNU5kcPZk6ZevQu3FvwuuG41PEx3U8kUyp\nIIuIiMfYcWwHU6KnMGnDJHyK+NCrYS+61utKdf/qrqMVOLFHY5m2aRrTNk3jfNp5etTrwSONHqFW\nQC3X0UQuSwVZREQ8jrWW5XuWMyV6CnNj5xJYIpDwOuGE1w2nTtk6Wlm+BtZatiZs5fPYz5m2aRpH\nfjlC93rd6VG/B80qNtNrKvmKCrKIiHi01LRUlu9Zzuwts5kTOwffIr7cX+d+OtbsSIvKLXQmjEs4\nn3qeH3b/wPyt81mwfQFJKUl0rtWZ7vW6c1vQbXgV8nIdUeSaqCCLiIikS7NprNm/hrkxc1kYt5D4\nxHjaVGtD+xrtaX9Te6qVruY6olPWWnYe38m38d+yZOcSFsct5qYyN9G5Vmc61+5McPlgrRRLgZCj\nBdkY0xUYCtwMNLPWrstiuw7ASMALGG+tfSeL7VSQRUQk1xw8fZCv475mUdwiFsctxr+4P2FBYYQG\nhRJaJZSbytxU4Avhz4k/8238txduu74lJS2FNje2oW21tnSs2ZFAv0DXEUWyXU4X5JuBNGAs8GJm\nBdkY4wVsBe4E9gGrgZ7W2phMtlVBdiQyMpLWrVu7juGx9Pq7pdffnbz02qfZNDYc3MCy3ctYvmc5\ny/csJzk1mdAqF8pySGAIwRWCKVO8jOuo1yzhlwTW7F/D6v2rWb1/Ncu+W4ZXdS9aV2tNm2ptaHtj\nW2oF1CrwvxTkFXnp/e9prqQgX/PwlbU29tcfcgnNgR3W2vj0bacB9wB/Ksjijv4jdUuvv1t6/d3J\nS699IVOIxoGNaRzYmGdveRaA3Sd2s3z3clbsWcHc2LlEH4qmVLFSBJcPvnCrEEztgNpU96+On7ef\n47/B706dO0Xs0VhijsYQcySGmKMxbDy8kSNnjhASGEKzis14sMGDVIuqxsiBI1WIHclL73/5s5w+\nOqESsCfD473ALTn8M0VERK5bUKkgghoE0bNBT+DCKnN8YjwbDm5gw6ENfLbpM7YnbGfn8Z2UKFqC\nGmVqUN2/OjX8a1DJrxI3+N7wh1tJ75LXVUbTbBonz50kMSmRI2eOsPfkXvae3Muek3t+u7/z+E6O\nJx2ndkBtbi57M3XK1uGhhg9Rt1xdagfU/sOBdVuKb1E5FsnCJQuyMeZroEIm33rVWjv/CvavmQkR\nESkQCplCVPevTnX/6txX577fnrfWcvD0QXYe30nc8TjijsWx9sBaDp85/IfbudRzlPIuRfEixSle\nuDg+RXwoXuTC10KmEKlpqaTa1D98PZtylsSkRBKTEjmdfJoSRUtQulhpAooHUKVUFaqUrELlkpUJ\nLh9M5ZKVqVq6KkGlgihkCjl8pUTyv+s+i4Ux5luynkFuAQy11nZIfzwISMvsQD1jjMq0iIiIiOS4\nHJtBvkhWP2QNUNMYUw3YD3QHema24eWCioiIiIjkhmv+DMYYc58xZg/QAvjSGPPf9OcrGmO+BLDW\npgDPAIuALcD0zM5gISIiIiKSV+SZC4WIiIiIiOQFzqf4jTEdjDGxxpjtxpiXXefxJMaYCGPMIWPM\nRtdZPJExpoox5ltjzGZjzCZjzADXmTyFMaaYMeYnY8x6Y8wWY8w/XGfyRMYYL2NMlDHmSg76lmxk\njIk3xkSnv/6rXOfxJMaY0saYWcaYmPT//7RwnclTGGNqp7/nf72dyOrfXqcryFdzIRHJfsaYVsBp\nYJK1toHrPJ7GGFMBqGCtXW+MKQGsBe7V+z93GGN8rLW/GGMKA8uAgdbaZa5zeRJjzAtAE8DPWtvF\ndR5PYozZBTSx1h5zncXTGGM+Ab6z1kak///H11p7wnUuT2OMKcSF7tncWrvn4u+7XkH+7UIi1trz\nwK8XEpFcYK39ATjuOoenstYetNauT79/mgsX0KnoNpXnsNb+kn63KOAFqCjkImNMZaAjMJ6sD/SW\nnKXXPZcZY0oBray1EXDhWC2VY2fuBOIyK8fgviBndiGRSo6yiDiTfqaXxsBPbpN4DmNMIWPMeuAQ\n8K21dovrTB7mPeAlIM11EA9lgSXGmDXGmD6uw3iQG4EjxpiJxph1xpiPjDE+rkN5qB7Ap1l903VB\n1hGC4vHSxytmAc+lryRLLrDWpllrGwGVgTBjTGvHkTyGMaYTcNhaG4VWMV0JtdY2Bu4Gnk4fuZOc\nVxgIAT601oYAZ4BX3EbyPMaYokBnYGZW27guyPuAKhkeV+HCKrKIRzDGFAFmA1OstZ+7zuOJ0j/e\n/BJo6jqLB7kV6JI+B/sZ0NYYM8lxJo9irT2Q/vUIMJcLI4+S8/YCe621q9Mfz+JCYZbcdTewNv39\nnynXBfm3C4mkt/nuwBeOM4nkCmOMASYAW6y1I13n8STGmLLGmNLp94sDdwFRblN5Dmvtq9baKtba\nG7nwMec31tqHXefyFMYYH2OMX/p9X6AdoLMZ5QJr7UFgjzGmVvpTdwKbHUbyVD258Mt5lrLrSnrX\nxFqbYoz59UIiXsAEHcGfe4wxnwG3AwHpF315w1o70XEsTxIKPAREG2N+LWeDrLULHWbyFIHAJ+lH\nMRcCJltrlzrO5Mk0bpe7ygNzL/yOTmFgqrV2sdtIHuVZYGr6wmAc8JjjPB4l/ZfCO4FLzt7rQiEi\nIiIiIhm4HrEQEREREclTVJBFRERERDJQQRYRERERyUAFWUREREQkAxVkEREREZEMVJBFRERERDJQ\nQRYRyYeMMaWMMf1d5xARKYhUkEVE8id/4CnXIURECiIVZBGR/OltoIYxJsoY847rMCIiBYmupCci\nkg8ZY6oCC6y1DVxnEREpaLSCLCKSPxnXAURECioVZBERERGRDFSQRUTyp1OAn+sQIiIFkQqyiEg+\nZK1NAJYbYzbqID0Rkeylg/RERERERDLQCrKIiIiISAYqyCIiIiIiGaggi4iIiIhkoIIsIiIiIpKB\nCrKIiIiISAYqyCIiIiIiGaggi4iIiIhkoIIsIiIiIpLB/wMvF4e3Wv7SRgAAAABJRU5ErkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "fig = figure(figsize=(12,4))\n", "p = plot(t, result, \"rx\", label=r\"$\\int_{0}^{x}sin(t) dt $\")\n", "p = plot(t, -cos(t) + cos(0), label=r\"$cos(0) - cos(t)$\")\n", "p = plot(t, dy_dt(0, t), \"g-\", label=r\"$\\frac{dy}{dt}(t)$\")\n", "l = legend(loc=\"upper right\")\n", "xl = xlabel(\"t\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### 高阶微分方程" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "抛物运动(竖直方向):\n", "\n", "$$\n", "\\frac{d^2x}{dt^2} = g - \\frac{D}{m}\\frac{dx}{dt}\n", "$$\n", "\n", "改写成如下形式:\n", "\n", "$$y = \\left[x, \\frac{dx}{dt}\\right] $$\n", "\n", "$$\\begin{aligned}\n", "\\frac{dy_0}{dt} &= y_1 \\\\\\\n", "\\frac{dy_1}{dt} &= -g - \\frac{D}{m} y_1 \\\\\\\n", "\\end{aligned}\n", "$$" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": true }, "outputs": [], "source": [ "def dy_dt(y, t):\n", " \"\"\"Governing equations for projectile motion with drag.\n", " y[0] = position\n", " y[1] = velocity\n", " g = gravity (m/s2)\n", " D = drag (1/s) = force/velocity\n", " m = mass (kg)\n", " \"\"\"\n", " g = -9.8\n", " D = 0.1\n", " m = 0.15\n", " dy1 = g - (D/m) * y[1]\n", " dy0 = y[1] if y[0] >= 0 else 0.\n", " return [dy0, dy1]" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": true }, "outputs": [], "source": [ "position_0 = 0.\n", "velocity_0 = 100\n", "t = linspace(0, 12, 100)\n", "y = odeint(dy_dt, [position_0, velocity_0], t)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEPCAYAAABY9lNGAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xnc1XPex/HXp6vJSJlsLVpEJEmFVErT1cI0lixjiSnZ\nhrnjzsxtEG40DHfWbLedxCTZJsnadjEG0aJSUm5CqDClspbrc//xPXG5nOpc13XO+Z7l/Xw8zsM5\nv+t3fr/PUV2f890+X3N3REREKqsVOwAREclNShAiIpKUEoSIiCSlBCEiIkkpQYiISFJKECIiklTG\nE4SZ3Wtmy81sXoVj15jZW2Y2x8weN7NfVfjZBWa22MwWmtlBmY5PRESSy0YLYhTQr9Kx54E93b0D\nsAi4AMDM2gLHAW0T77nVzNTKERGJIOO/fN39n8DKSscmuXt54uV0oFni+eHAWHdf5+5LgHeAzpmO\nUUREfi4Xvp2fAjydeL4jsLTCz5YCTbMekYiIxE0QZnYR8J27P7iJ01QLREQkgtqxbmxmJwEHA30q\nHP4IaF7hdbPEscrvVdIQEakGd7dUz43SgjCzfsC5wOHu/k2FH00ABphZHTPbGdgNeC3ZNdy9YB+X\nXnpp9Bj0+fT5ivHzFfJnc6/69+qMtyDMbCzQE9jezD4ELiXMWqoDTDIzgFfcfYi7LzCzh4EFwHpg\niFfnU4mISI1lPEG4+/FJDt+7ifOvBK7MXEQiIpKKXJjFJJWUlpbGDiGj9PnyWyF/vkL+bNVh+diD\nY2bqeRIRqSIzw3N9kFpERHKfEoSIiCSlBCEiIkkpQYiISFJKECIikpQShIiIJKUEISIiSSlBiIhI\nUkoQIiKSlBKEiIgkpQQhIiJJKUGIiEhSShAiIpKUEoSIiCSlBCEiIkkpQYiISFJKECIiklTG96SW\n3FJeDh99BO++C599BitXhse334JZeNSuDdtsA9tuGx7NmkGLFrDFFrGjF5Fs0pajBez772HePHjl\nlfCYNSskhgYNYJddoGHDkAgaNIBf/hLcw2PdOli1Cv79b/j8c/jww5BUGjaE3XaDDh3Co2NHaNcu\nJBQRyX1V3XJUCaLAfPUVPP88/OMfMHEi7LADdOsG++8P++0Hu+4K9epV/brr18PSpfD22zBnTnjM\nmgUffwxdu8IBB0Dv3tClixKGSK5SgihSs2bBbbfBI4/AvvvCkUdC//6hayiTPvsMXn4Z/vlPmDwZ\n3n8f+vaF3/423H+77TJ7fxFJXc4lCDO7FzgEWOHueyWObQuMA3YClgDHuvuqxM8uAE4BvgeGuvvz\nSa6pBEEYT3jsMbj2Wli2DM44A049FRo1ihfTJ5/Ac8/Bk0+GhNG5M/zud3D00bD99vHiEpHcTBA9\ngLXA/RUSxNXAZ+5+tZmdD2zj7sPMrC3wILAf0BSYDLR29/JK1yzqBOEeuo8uvjh051x8MRx6KJSU\nxI7sp778Ep59NrRqnnkGevWCQYNCrBrwFsm+nEsQAGbWEniyQoJYCPR09+Vm1hgoc/c2idZDubtf\nlTjvWWC4u79a6XpFmyDmzYM//hFWr4bLL4fDDw8zj3Ld6tWhtXP//TB/PgweDKefHga9RSQ7qpog\nYq2DaOTuyxPPlwMbOkV2BJZWOG8poSVR9L7+Gi68EPr0Cb9c58yBI47Ij+QAsPXWcPLJMG1aGLOo\nVQu6dw/jFRMmhBlXIpJbos83cXc3s001B5L+bPjw4T88Ly0tpbS0NL2B5ZBXX4WBA6FTJ5g7Fxo3\njh1Rzey6K1x1FVx2WWhV/O1v8Oc/w9ChcMopUL9+7AhFCkNZWRllZWXVfn/MLqZSd19mZk2AaYku\npmEA7j4icd6zwKXuPr3S9Yqii8kdbrgBRoyAO+4ILYZC5B6S4MiRoYVxxhkhWTRsGDsykcKSL11M\nE4DBieeDgfEVjg8wszpmtjOwG/BahPiiW7kSjjoKxo6F6dMLNzlA6Cbbf394+OGQKD7/HNq0gf/8\nz7D2QkTiyHiCMLOxwMvA7mb2oZmdDIwADjSzRUDvxGvcfQHwMLAAeAYYUhRNhUqWLAm/MJs1g5de\ngpYtY0eUPa1ahfUcCxbAlltC+/YwZEhYzS0i2aWFcjlm9uwwDXTYsPANutitWAHXXQd33w2//30Y\nqM/3MRiRWPKli0mSmDQJfvMbuOkmJYcNGjYMA9pvvRXWebRtCxdcELrgRCSzlCByxNNPh2/Ijz8e\nVh7LTzVsGAax33gDPv0UWrcOLYtvv40dmUjhUoLIAZMnw0knhfIUBxwQO5rc1qJF6G564YXwaNMG\nHnwwzIQSkfTSGERkL74YWgyPPw49esSOJv+88AL85S+h5MgNN4RqsiKSnMYg8sjrr4cidg89pORQ\nXT17hmnAf/xjmBY8cKCmxoqkixJEJEuXhrUNd94ZymdI9dWqFcqPvP12mBLcsWNYXKjxCZGaUYKI\n4Msvw14JZ59d2Avgsq1evVC2Y/r0UO9pr71CNVkRqR6NQWRZeXkYc9hmG7jnnvwptpePnn46lOzo\n2DHMgGrePHZEInFpDCLHXXJJKCVx++1KDpl28MGhPHq7drD33mFjpXXrYkclkj/Ugsii558PO77N\nmhX2ipbsWbwYzjwzrMy+886w051IsVELIkctWxbWOjzwgJJDDLvtFrZCPffcsMnS0KGwZk3sqERy\nmxJEFpSXh602//AHKOBtK3KeWVitPn9+mCiw557w1FOxoxLJXepiyoIRI8KA6dSpYUGX5IYpU8K2\np127hkV2atlJoVMXU46ZNSvMoBkzRskh1/TpEwaxmzQJU2LHjVPJDpGK1ILIoHXrwmDon/8MJ54Y\nOxrZlOnTw3anu+8Ot96qkuJSmNSCyCHXXhuqkA4aFDsS2ZwuXUJrb489oEOH0OLLg+8gIhmlFkSG\nvP02dO8OM2YU145whWDmzFC6o3XrsF5Fe2NLoVALIgeUl8Npp8HFFys55KN99w2JvXXrsOXpo4/G\njkgkDrUgMuDOO2HUqLCfdElJ7GikJl59NYwfde4MN98cSqSI5Cu1ICJbtSq0HG69VcmhEHTtGvYJ\n32ab0JqYNCl2RCLZoxZEmp1zDqxeDXfdFTsSSbdJk8JMpyOPDPtkb7ll7IhEqqaqLQgliDRatAi6\ndQsrdRs1ih2NZMLKlTBkSNgbe8wY2Gef2BGJpE5dTBGdey6cd56SQyHbZhsYOzZ0I/brF1bJf/99\n7KhEMiNqC8LMLgAGAuXAPOBkYCtgHLATsAQ41t1XVXpfzrUgJk+GM86ABQtgiy1iRyPZ8MEHPy6A\nvP9+aNEibjwim5M3LQgzawn8AdjH3fcCSoABwDBgkru3BqYkXue08nL4r/+Ca65RcigmLVqEek79\n+kGnTqFUh0ghidnFtBpYB9Q1s9pAXeBjoD8wOnHOaCDnN+V85JEwYHnkkbEjkWwrKYFhw0Ixxosv\nhpNPhrVrY0clkh7REoS7/xu4DviAkBhWufskoJG7L0+cthzI6R7977+H4cPhssu0Q1wx69QplOqo\nVSsMXM+YETsikZqLVl/UzFoBfwJaAl8Aj5jZwIrnuLubWdLBhuHDh//wvLS0lNJIGy089BBstx0c\ndFCU20sOqVcv7DP+8MNhu9Pzzgtdj7U0FUQiKSsro6ysrNrvjzZIbWbHAQe6+2mJ14OArkBvoJe7\nLzOzJsA0d29T6b05MUi9fj20bRvq9fTuHTsaySXvvw8nnAD168Po0ZrZJrkhbwapgYVAVzPb0swM\n6AssAJ4EBifOGQyMjxTfZv3977DjjtCrV+xIJNfstBO88ELoetp7b63AlvwUe5rreYQkUA7MAk4D\n6gMPAy3I4Wmu69aFvQNGjYKePaOGIjlu2rRQ8v3EE+Gvf4Vf/CJ2RFKstJI6S+6/H+67L2wjKrI5\nK1aEEuKrV8ODD4YWhki25VMXU95yD5sBnX9+7EgkXzRsCE89FaZCd+4MTzwROyKRzVOCqIZJk0KS\n0MwlqYpateAvfwnJ4eyzw1a0330XOyqRjVOCqIZrrgn/0LXuQaqja9ewZuLdd+GAA+C992JHJJKc\nEkQVvfFGqLd0/PGxI5F8tu22MH58mArbpUt4LpJrNEhdRYMGQbt2Gn+Q9Jk+HY477sd9JurUiR2R\nFCrNYsqgDz+EDh1C10CDBlm/vRSwf/871HFatiysxNYsJ8kEzWLKoJtvDlMVlRwk3TZ0OR1zTJjl\n9NRTsSMSUQsiZd98A82bh03sW7XK6q2lyPzrXzBgAAwcCJdfDrWjVUyTQqMWRIY8+mio0qnkIJnW\nvXuY5TRzJhx4YOh2EolBCSJFt98Of/xj7CikWOywAzzzTCjjsu++8OKLsSOSYqQuphTMmxd2DVuy\nRHV0JPuefRZOOimUDj/3XK2/kepTF1MG3HEHnHaakoPE0a8fvPYaPPZYmAq7atXm3yOSDkoQm7F2\nbSiudtppsSORYtaiRehmatYslBCfMyd2RFIMlCA246GH4Ne/DjOYRGLaYgu45ZZQMrxv37ARkUgm\naQxiMzp1gr/9LTTzRXLFm2/C734XNqu68caQPEQ2R2MQaTRnTqjjr6qtkmvatYPXX4dPP4UePeCD\nD2JHJIVICWITHngg1F7SpvOSi7beOqzPOeaYUPBvypTYEUmhURfTRqxfHwYGp06FNm0yeiuRGps2\nLVSGPfvsUEhSU2Elmap2MW1yEb+ZNQSOAX4NtAQceB94EXjE3VdUP9TcNmVKmDGi5CD5oFev0OV0\n9NGhOuzo0aGFIVITG+08MbN7gIeBesDtwGDgZOAOoD7wsJndnY0gY9jQvSSSL5o1gxdegMaNQ8G/\nBQtiRyT5bqNdTGbW3t3nbvLNKZyTCZnuYlqzJkxrXbw4lDwQyTf33RdWXd92W2hViID2g0iL++6D\nf/xDG8tLfps1K0yFPfZYuOIKVYWVDExzNbPDzGy2ma00szWJx+qahZnb1L0khWCffcK4xKxZYR3P\nZ5/FjkjyTSoTOG8gjD9s5+71E4+0DH+ZWQMze9TM3jKzBWbWxcy2NbNJZrbIzJ43s6xuz/PBB2Hf\n6UMPzeZdRTJj++1Dsb9OncJj5szYEUk+SSVBLAXmu3t5Bu5/I/C0u+8BtAcWAsOASe7eGpiSeJ01\n48aFZvkvf5nNu4pkTkkJjBgB114bWhIq0SGp2uwYhJl1BS4DpgHfJQ67u19foxub/QqY7e67VDq+\nEOjp7svNrDFQ5u5tKp2TsTGILl1Cf23fvhm5vEhU8+eHirAHHQTXXw916sSOSLIpE6U2LgfWAr8k\nTHmtR5jmWlM7A5+a2Sgzm2Vmd5nZVkAjd1+eOGc50CgN90rJ++/Du++GTVpECtGee4bS4e+/D336\naLc62bRU5jU0cfcDM3TvfYCz3P11M7uBSt1J7u5mlrSpMHz48B+el5aWUlpaWuOAHnsMDj9c+z5I\nYWvQIMzQu+wy2G+/UK6jS5fYUUkmlJWVUVZWVu33p9LFdDUwxd2fq/Zdkl+3MfCKu++ceH0AcAGw\nC9DL3ZeZWRNgWra6mLp1g0suUeVWKR4TJoS9Tv7nf+DUU2NHI5mW9nUQZrYWqEsYf1iXOOzpmMlk\nZi8Cp7n7IjMbnrgPwOfufpWZDQMauPuwSu9Le4JYuhQ6dAhNbrUgpJgsXAhHHAG9e8MNN2hcopDl\n1UI5M+sA3A3UAf6PUMqjhFDiowWwBDjW3VdVel/aE8RNN8Hs2TBqVFovK5IXvvgCTjwRPv88dDk1\nbhw7IsmEtCUIM2vl7v+3mZtt9pxMyESC6NEDhg2DQw5J62VF8kZ5OVx+Odx9dxiP69w5dkSSbulM\nEOOArYAJwAzgE8CAJkAnoD+wxt0H1DToqkp3gvj447AByyefaGcukSeegD/8Aa66Ck4+OXY0kk5p\n7WIys12BAUB3YKfE4feBl4Cx7v5uDWKttnQniP/931Ai+f7703ZJkbz21lthXOI3v4HrrtO4XKHI\nqzGI6kp3gujdO2y0cvjhabukSN5btQp+/3v48kt45BFVNi4E2pO6ilatghkz4MBMrPQQyWMNGoRp\nsN27h/USs2bFjkiyregTxHPPwa9/DXXrbv5ckWJTUhJKz1xzTehuevDB2BFJNhV9hfiJE1W5VWRz\njjkGdt89jEvMnh2K/5WUxI5KMi2V/SCmpHIsH33/PTzzjKa2iqSiffsf95c4+GBYuTJ2RJJpm9qT\neksz2w7YIbFHw4ZHS6BptgLMpFdfDfv4Nm8eOxKR/LDddqFbtm3bMC4xf37siCSTNtXFdAZwNrAj\nUHGbkTXALZkMKlsmToTDDosdhUh+qV0bRo6EvfeG0tKwsE4zAAtTKrWYhrr7TVmKJyXpmubarh3c\nc48qWYpU12uvhQ22Tj8d/vu/wVKeQCkxZGQdhJl1A1pSocXh7tGWlaUjQbz3HnTtGlZP1yr6uVwi\n1ffJJ3DUUaG7dtQoqFcvdkSyMWlfB2FmfweuBQ4A9qvwyGtPPRUG2pQcRGqmSRMoK4P69cOaiffe\nix2RpEsq01z3BdpmbI/PSCZODPVmRKTmttgidNfefDPsvz+MHQu9esWOSmoqle/PbxIK9BWMtWvh\nX//S6mmRdDKDoUNhzBgYMABuuQUK62tl8dloC8LMnkw8rQcsMLPXgG8Tx9zd+2c6uEwpKwtT9Lau\n8ZZHIlJZnz7w8sthZtOcOaEYpjYhyk+b6mK6LmtRZNnkyWo9iGRSq1bwyiswaFBIGI89Bg0bxo5K\nqqooq7nuuSeMHg2dOqUxKBH5mfJyGD48/HsbPz6snZB4MrEn9Zokh78AXgfOibEnRE0SxEcfhZIB\nK1aoloxItjzyCAwZEsYljjsudjTFq6oJIpVZTDcCHwJjE68HAK2A2cC9QGkVY4xq8uSw/4OSg0j2\nHHMM7LZbKPY3d27Y2lRTzHNfKi2Iue7evtKxN9y9o5nNcfcOGY0weUzVbkEMHBjKe59+epqDEpHN\nWrECjj4att0WHnggrJ2Q7MnEhkFfmdlxZlYr8TgW+Cbxs7wawHDXALVITA0bhn+DDRuG9RLvRtm0\nWFKVSoL4PTAIWJF4nAgMNLMtgbMyGFvazZsXygDsvHPsSESKV506cMcd8B//Ad26wdSpsSOSjSmq\nWUzXXQfvvAO33ZaBoESkyqZOhRNOgEsuCQlDxf4yK22D1GZ2vrtfZWY3J/mxu/vQakX48/uUADOA\npe5+mJltC4wDdgKWAMe6+6p03GvSJI09iOSS3r1DVYP+/cPg9U03aVFdLtlUF9OCxH9nVnjMqPA8\nXc5O3GtDk2AYMMndWwNTEq9r7Ntvw19E1YcRyS0bFtV9/HEYH/z009gRyQYpdzGZ2Vbu/mVab27W\nDLgPuAL4r0QLYiHQ092Xm1ljoMzd21R6X5W7mKZNg2HDYPr0NAUvImlVXh72lBg7Fp54IqxXkvTK\nRLnvbma2AFiYeN3RzG6tQYwVjQTOBcorHGvk7ssTz5cDjdJxI81eEslttWrBlVfCFVeE8hzjx8eO\nSFJZKHcD0A94AsDd3zCznjW9sZkdCqxw99lmVprsHHd3M0vaVBg+fPgPz0tLSyktTXqJH7zwQljy\nLyK57YQTwqK6o46CN9+Eiy7S4HV1lZWVUVZWVu33p7JQ7jV372xms91978SxGi+QM7MrCdNn1wO/\nBLYGHidsRlTq7svMrAkwraZdTF9/DTvsAMuXw1Zb1SRqEcmWjz+GI48M09LvvRfq1o0dUf7LxEK5\nD8yse+LidczsL8Bb1Q1wA3e/0N2bu/vOhPIdU919EDABGJw4bTBQ44bm9Omw115KDiL5ZMcdQ2n+\n2rVD9YOlS2NHVHxSSRD/AZwJNAU+AvZOvE63DU2CEcCBZrYI6J14XSMvvBD+golIftlyy1CS45hj\noEsXTTLJtqJYKNenD5xzTtiDWkTy05NPwimnwMiRoaaaVF3ayn1XWiDnQMWLpm2hXHVUJUF89x1s\nt11onv7qVxkOTEQy6s03w6K6Y48Ns51Ulblq0jkGUXFh3OH8dJFcOhfKZdSMGWFGhJKDSP5r1w5e\new1efTWUDl+9OnZEhS2lLqaKM5hyQVVaECNGhNlLI0dmOCgRyZrvvoOhQ+Gll2DCBNhll9gR5YdM\nzGLKaxqgFik8deqEopsbKsLWYKq/bEJBJ4j16+Hll6FHj9iRiEi6mcGZZ8KYMWEb0zvvjB1R4dlU\nNde1/Dj1dMtKe1O7u2+d0cjS4I03oHlz2H772JGISKb06RO6mvr3D4PY118f1k5IzW20BeHu9dy9\nfuJRu8Lz+vmQHABefFHdSyLFYLfdQkXYRYvgt7+FlStjR1QYCrqL6cUXoWeNq0aJSD5o0AAmTgxV\nE7p0gbffjh1R/ivYhXLuof7S3Llhyb6IFI977oELLwyrsA86KHY0uUOzmBIWLw77Tys5iBSfU0+F\nRx+FwYPhxhvDF0apuoJNEK+8Al27xo5CRGLp0SPMYrz7bjjjjLB2QqqmoBPE/vvHjkJEYtp555Ak\nli0LG4Z99lnsiPKLEoSIFLT69eEf/wi/D7p0gfnzY0eUPwpykHrNGmjcOEx1q1Mni4GJSE574IFQ\n2XnUKDjkkNjRZJ8GqYHXX4eOHZUcROSnBg2CJ56A00+Ha6/V4PXmFGSCUPeSiGzM/vuH3xFjxoTZ\nTt9+Gzui3FWQCeLVVzWDSUQ2rkWLUJ5j1Sro2xdWrIgdUW4quAThHhKEWhAisilbbRXWSvTqFQav\n586NHVHuKbgE8c47YR/bpk1jRyIiua5WLbjsMrjyylD074knYkeUWwqu5qEWyIlIVR1/PLRqBUcd\nBQsXwnnnhXLixa7gWhDqXhKR6ujcOfz+ePjhUKLjm29iRxRfwSUIzWASkepq1gz++U/4+mvo3Tts\nV1zMCipBrF0b6sHvnTO7Z4tIvqlbF8aNC6U5OneGOXNiRxRPtARhZs3NbJqZzTezN81saOL4tmY2\nycwWmdnzZtYg1WvOmAHt28MWW2QubhEpfLVqwV//CldfHabBjh8fO6I4YrYg1gF/dvc9ga7AmWa2\nBzAMmOTurYEpidcpmTED9tsvI7GKSBE67jh4+mk46ywYMaL4Vl5HSxDuvszd30g8Xwu8BTQF+gOj\nE6eNBo5I9ZozZkCnTumOVESK2X77wfTpP+4vUUyD1zkxBmFmLYG9gelAI3ffMDS0HGiU6nVmzoR9\n9017eCJS5Jo2DVsYf/NNcQ1eR18HYWb1gMeAs919jVWYfOzubmZJG3XDhw//4XlpaSkdO5aybBm0\naZPhgEWkKNWtCw89FBbWdekSFtV16BA7qk0rKyujrKys2u+PWu7bzH4BTASecfcbEscWAqXuvszM\nmgDT3L1Npff9rNz31KlwySWhvoqISCaNGxfGJe66C45IuRM8vrwp922hqXAPsGBDckiYAAxOPB8M\npDR/QOMPIpItxTJ4HXMMojswEOhlZrMTj37ACOBAM1sE9E683iyNP4hINhXD4HXB7CjXqhVMnAh7\n7BEpKBEpSl99BSedBEuXhq1NG6U8rSb78qaLKZ1Wrgz13Fu3jh2JiBSbDYPXBx4YBq8LaeV1QSSI\nmTNDeY2SktiRiEgx2rDyesSIsPK6UMqGR5/mmg4afxCRXDBgAOyyS+GUDS+IFoRmMIlIrqhYNvyk\nk/J7z+uCSBBqQYhILmnWLKy8/uqrsFNdvu55nfcJ4vPP4bPPNEAtIrllq63CgrrevcPg9bx5sSOq\nurxPEBsGqGvl/ScRkUJTec/rJ5+MHVHV5P0g9cyZGn8Qkdx2/PE/Dl6//Tacc05+DF7n/ffumTNh\nn31iRyEismlduoTB6zFj4NRT82PwOu8TxJw50LFj7ChERDavefNQUHTVqrCw7tNPY0e0aXmdINau\nhY8+gt13jx2JiEhqttoq1G/q0SO0Kt58M3ZEG5fXCWLevFB7qXbej6SISDGpVQuuuCKsvu7dO1SG\nzUV5nSDmzMn9DTtERDZm0CAYPx5OOw1Gjsy9suFKECIiEXXrBq+8AqNGwemnw3ffxY7oR0oQIiKR\n7bQT/OtfYa/rgw4KC4BzQd4miPLyMAahBCEihaB+/bCfROfOYfD6rbdiR5THCeK992CbbcJDRKQQ\nlJTA1VfDRRdBz57w3HNx48nbBKHuJREpVCefDI89FrYyvfnmeIPXShAiIjmoRw94+WW4/XY480xY\nty77MShBiIjkqF12CTOcliyB3/42bK+cTUoQIiI5bOutQxXYvfaCrl1h8eLs3TtvE8Snn0KrVrGj\nEBHJvJKSsJDunHPggANg6tTs3DdvE0S7duF/mohIsTj9dHjoITjhBLjjjszfLycThJn1M7OFZrbY\nzM5Pdk779tmOSkQkvl69QkXYkSPhT3+C9eszd6+cSxBmVgLcAvQD2gLHm9kelc/T+IOIFKtddw17\nSyxYAIcdBl98kZn75FyCADoD77j7EndfBzwEHF75JCUIESlmDRqEKrC77hq2PLj0Uli2LL33yMVC\n2U2BDyu8Xgp0qXySuphEpNjVrh0W0g0ZAjfdFLY/6NsXdtghTddPz2XSKqU1g9dfP/yH56WlpZSW\nlmYoHBGR3LbHHnDbbXDllaF8+Ndfh+OLFpWxeHFZta9rnmMFyM2sKzDc3fslXl8AlLv7VRXO8VyL\nW0Qk15kZ7m6pnp+LYxAzgN3MrKWZ1QGOAyZEjklEpOjkXBeTu683s7OA54AS4B53z4HCtyIixSXn\nuphSoS4mEZGqK4QuJhERyQFKECIikpQShIiIJKUEISIiSSlBiIhIUkoQIiKSlBKEiIgkpQQhIiJJ\nKUGIiEhSShAiIpKUEoSIiCSlBCEiIkkpQYiISFJKECIikpQShIiIJKUEISIiSSlBiIhIUkoQIiKS\nlBKEiIgkpQQhIiJJKUGIiEhSShAiIpJUlARhZteY2VtmNsfMHjezX1X42QVmttjMFprZQTHiExGR\neC2I54E93b0DsAi4AMDM2gLHAW2BfsCtZlZ0rZyysrLYIWSUPl9+K+TPV8ifrTqi/PJ190nuXp54\nOR1olnh+ODDW3de5+xLgHaBzhBCjKvS/pPp8+a2QP18hf7bqyIVv56cATyee7wgsrfCzpUDTrEck\nIiLUztQ/zu36AAAFzElEQVSFzWwS0DjJjy509ycT51wEfOfuD27iUp6J+EREZNPMPc7vXzM7CfgD\n0Mfdv0kcGwbg7iMSr58FLnX36ZXeq6QhIlIN7m6pnhslQZhZP+A6oKe7f1bheFvgQcK4Q1NgMrCr\nx8piIiJFLGNdTJtxM1AHmGRmAK+4+xB3X2BmDwMLgPXAECUHEZE4onUxiYhIbsuFWUxVYmb9Eovo\nFpvZ+bHjSScza25m08xsvpm9aWZDY8eUbmZWYmazzezJ2LGkm5k1MLNHE4tAF5hZ19gxpVNiEet8\nM5tnZg+a2RaxY6oJM7vXzJab2bwKx7Y1s0lmtsjMnjezBjFjrImNfL6NLlJOJq8ShJmVALcQFtG1\nBY43sz3iRpVW64A/u/ueQFfgzAL7fABnE7oQC7HpeiPwtLvvAbQH3oocT9qYWUvCpJJ93H0voAQY\nEDOmNBhF+F1S0TBgkru3BqYkXuerZJ8v6SLljcmrBEEYvH7H3Ze4+zrgIcLiuoLg7svc/Y3E87WE\nXzA7xo0qfcysGXAwcDeQ8kyKfJD4JtbD3e8FcPf17v5F5LDSaTXhC0xdM6sN1AU+ihtSzbj7P4GV\nlQ73B0Ynno8GjshqUGmU7PNtYpFyUvmWIJoCH1Z4XbAL6RLf2PYm/CEWipHAuUD55k7MQzsDn5rZ\nKDObZWZ3mVnd2EGli7v/mzDz8APgY2CVu0+OG1VGNHL35Ynny4FGMYPJsIqLlJPKtwRRiN0SP2Nm\n9YBHgbMTLYm8Z2aHAivcfTYF1npIqA3sA9zq7vsAX5Lf3RM/YWatgD8BLQmt2npm9vuoQWVYYgZl\nQf7OSXGRct4liI+A5hVeN+enpTnynpn9AngM+Lu7j48dTxp1A/qb2XvAWKC3md0fOaZ0WgosdffX\nE68fJSSMQtEJeNndP3f39cDjhD/TQrPczBoDmFkTYEXkeNIusUj5YGCzCT7fEsQMYDcza2lmdQiV\nXydEjiltLCwKuQdY4O43xI4nndz9Qndv7u47EwY3p7r7ibHjShd3XwZ8aGatE4f6AvMjhpRuC4Gu\nZrZl4u9pX8Jkg0IzARiceD4YKKQvaRsWKZ8LHL6hgsWm5FWCSHxzOQt4jvCXc5y7F8xMEaA7MBDo\nlZgKOjvxB1qICrHp/p/AGDObQ5jFdGXkeNLG3ecA9xO+pM1NHL4zXkQ1Z2ZjgZeB3c3sQzM7GRgB\nHGhmi4Deidd5KcnnO4WwSLkeYZHybDO7dZPX0EI5ERFJJq9aECIikj1KECIikpQShIiIJKUEISIi\nSSlBiIhIUkoQIiKSlBKEFD0z267CupNPzGxp4vkaM7slQ/c8K7GidWM/729mF2fi3iKp0joIkQrM\n7FJgjbtfn8F7GDAL2C+x+HNj58xOnLMuU7GIbIpaECI/ZwBmVrphYyMzG25mo83sRTNbYmZHmdm1\nZjbXzJ5JlMDGzPY1szIzm2Fmz26o61NJd2DhhuRgZkMTG/HMSax+3VAo7hXgoGx8YJFklCBEUrcz\n0IuwZ8DfCRvLtAe+Bg5JFFq8Gfidu3cibNhyRZLrHEAoWbHB+UDHxCYuZ1Q4/hrw67R/CpEU1Y4d\ngEiecOAZd//ezN4Earn7c4mfzSOUwW4N7AlMDj1ElBD2TqisBfBShddzgQfNbDw/LQ73MT/fEUwk\na5QgRFL3HYC7l5tZxXGBcsK/JQPmu3sqZbAr7olxCKGlcBhwkZm1S+z6VYvCLGooeUJdTCKpSWWT\no7eBHcysK4S9PcysbZLz3gc27DlgQAt3LyNsMPQrQrVNgCaJc0WiUIIQ+Tmv8N9kz+Hn3+w9Mdvo\naOAqM3uDMAtp/yTXf4mwAQ+ElscDZjaXMLPpRndfnfhZZ+DFmnwQkZrQNFeRLKswzbWLu3+3kXNq\nJc7ptLGpsCKZphaESJYlprDexaa3fDwUeFTJQWJSC0JERJJSC0JERJJSghARkaSUIEREJCklCBER\nSUoJQkREklKCEBGRpP4fPzntqp99LaMAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "p = plot(t, y[:,0])\n", "yl = ylabel(\"Height (m)\")\n", "xl = xlabel(\"Time (s)\")" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Integration successful.\n", "['hu', 'imxer', 'leniw', 'lenrw', 'message', 'mused', 'nfe', 'nje', 'nqu', 'nst', 'tcur', 'tolsf', 'tsw']\n", "cumulative number of function evaluations at each calculated point: [ 45 49 51 53 55 59 61 61 63 65 67 67 69 71 73 73 75 77\n", " 77 79 79 81 81 83 85 85 87 87 89 89 91 91 93 95 95 97\n", " 97 99 99 101 101 103 103 105 107 107 109 109 111 111 113 113 115 115\n", " 117 117 119 119 121 121 123 123 123 125 125 127 127 129 129 131 131 131\n", " 133 133 135 135 135 137 137 139 139 139 141 141 143 143 143 145 145 147\n", " 147 149 149 149 154 158 274 280 280]\n", "cumulative number of time steps [ 20 22 23 24 25 27 28 28 29 30 31 31 32 33 34 34 35 36\n", " 36 37 37 38 38 39 40 40 41 41 42 42 43 43 44 45 45 46\n", " 46 47 47 48 48 49 49 50 51 51 52 52 53 53 54 54 55 55\n", " 56 56 57 57 58 58 59 59 59 60 60 61 61 62 62 63 63 63\n", " 64 64 65 65 65 66 66 67 67 67 68 68 69 69 69 70 70 71\n", " 71 72 72 72 73 75 130 133 133]\n" ] } ], "source": [ "y, infodict = odeint(dy_dt, [position_0, velocity_0], t, full_output=True, printmessg=True, )\n", "print sorted(infodict.keys())\n", "print \"cumulative number of function evaluations at each calculated point:\", infodict['nfe']\n", "print \"cumulative number of time steps\", infodict['nst']" ] } ], "metadata": { "kernelspec": { "display_name": "Python 2", "language": "python", "name": "python2" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 2 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython2", "version": "2.7.10" } }, "nbformat": 4, "nbformat_minor": 0 }