{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "![](ChemE_logo.png \"ChemE Logo\")\n", "\n", "# ExempliPy\n", "by [Tony Saad](www.tonysaad.net)
\n", "Assistant Professor of [Chemical Engineering](www.che.utah.edu)
\n", "[University of Utah](www.utah.edu)\n", "\n", "\n", "A collection of example of problems solved numerically with Python. Applications span physics, chemical, mechanical, civil, and electrical engineering." ] }, { "cell_type": "markdown", "metadata": { "toc": "true" }, "source": [ "# Table of Contents\n", "

1  ExempliPy
2  Free Fall: ODE Time Integration
2.1  Method 1: Using Lists
2.2  Method 2: Using Numpy Arrays
3  Interpolation
" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Free Fall: ODE Time Integration" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "fragment" } }, "source": [ "Consider the free fall of an astronaut subject to drag. The governing equation is according to Newton's second law is:\n", "$$m \\frac{\\text{d}u}{\\text{d} t} = m g - c u$$\n", "or\n", "$$\\frac{\\text{d}u}{\\text{d} t} = g - \\frac{c}{m} u$$\n", "where $u$(m/s) is the (downward) speed of the astronaut, $g$(m/s/s) is the acceleration of gravity, and $c$(kg/s) is the drag coefficient. Here, the drag force acts in the direction opposite to the fall of the astronaut and is given by $F_\\text{D} = cu\\mathbf{j}$" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "## Method 1: Using Lists" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false, "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD8CAYAAABn919SAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHLFJREFUeJzt3XmYXHWd7/H3t5Ze0muWTqezNglJIJCEkCYEUAhCJIMg\nXMSMoEzkOuK4zKNe7ij6eK/Lo4+oM3PVGe7MBeEhI5tRkUVFgRAFhDR0VkL2dGfv9JpOb+mlqn73\nj6pIE5J0p7u6T52qz+t56jnn/OoU/T2c9Ccnv/qd3zHnHCIi4n8BrwsQEZHkUKCLiKQJBbqISJpQ\noIuIpAkFuohImlCgi4ikCQW6iEiaUKCLiKQJBbqISJoIjeQPGzdunCsvLx/JHyki4nvr1q1rdM6V\n9LffiAZ6eXk5VVVVI/kjRUR8z8z2DWQ/dbmIiKQJBbqISJpQoIuIpAkFuohImlCgi4ikCQW6iEia\nUKCLiKSJER2HLiJytpxzxBz0RmNEY45IzBHt84rEYu/ajjpHJOqIufh2fEmf9fg+rk+7c/G2mINY\nYr+YI76MvbPu+ra7E7U5nONd+5zYdryz3y0XT+accXnD+v9KgS4iQDzIuiMxunqjHO+N0tUbpas3\nRlckSndvjO5IfLs7EqU7EqM7EqMn8eqOROmJxOiNJtqi7q/bJ149UUfvibZYfD0SixGJOnqi8WUk\nEdDx9Xfa/M4MFk4brUAXkfdyLh6+rV29tHdFaOuK0N79zrKjO0JHT2LZHaWjO0JnT5TOnhPL+HpX\nb4zjve+sD0UoYGSFAoSDAbJCAbL6LMMhIxyMv5cdDpAXCBAOxttCwQDhgBEKWp/1AKFEWzAQbwsG\njaDF3wsaBBP7BM0IJvYNmBEKGIE+7SfWAwHe1RawPu2J7fiL+HsBw/quw18/FzAwe+f9E58xA+Od\n908sR4oCXcRDzjlauyIc7eihubOHox09HO3spaWzh5bOXo4d76XleHy7tStC2/FeWrt6aT0eoSfa\nfwAHA0ZeVpBRWSFGZQfJywqRmxVkXH4WuVm55IZD5GYFGJUVIiccJCccIDccJCccTCwDZIeD5ISC\nZIcDZIcCZIeC8WU4QHYwGA/tUIBgYOSCS05NgS6SZM45Wo9HqG/roqGtm/q27r+uN7b30NjeTVN7\nD00d8eXpuhQCBkW5YYpHZVGYG6Y4N8yU0bkU5oYpzAlTmBuiIDtEQU6YgpwQ+dkh8nNCFGSHycsO\nkpcdIjsUGNErRPGWAl3kLHV0RzjUcpxDR49zsOU4h1uOc+RYV/zVGl8e742+53M54QDj8rMZm5/N\nxOIc5k4qYkx+FmPzshg9KosxeVmMzsti9Kgwo/OyyM8KEdBVr5wFBbrISZxz1LV2U93YzoHmTvY1\ndbK/+Z1XS2fvu/YPBYzSwhwmFOUwZ2Ih15w3nglFOZQUZDO+IIfxhdmML8gmPzukq2UZVgp0yVjd\nkSjVDR3sqm9nd10b1Y0dVDd0sLepg86ed66wQwFj8uhcpowZxfVzy5g8OpdJxbmJ5ShKCrLVfywp\nQYEuac85x6GW42w93MrW2la217axs76NfU2dRBP91wGDKWNGcc64PBZPH8s5JXmcMzaPaWNHUVaU\nQyioe/Ak9SnQJa045zjQfJyNB1vYfKCFtw4dY1ttK61dESA+Hrh8bB6zSvP50NwyZpYWMKs0n3PG\n5ZEdCnpcvcjQKNDF1zq6I2w80MKbe5vZsL+FzQdbOJro484KBZhTVsgN8ycyp6yQORMLOW9CAaOy\n9Mde0pP+ZIuvtHT2sLa6ibXVzazbd5Stta1EYw4zmDW+gA/OmcC8KUXMn1zM7AkFhNVVIhlEgS4p\nrbMnwhs1zby2p4m/7G5ka20rzsWHAC6YMprPLZlBRfkYFkwtpjAn7HW5Ip5SoEvKqWnsYM32etbs\nqKeyupmeaIysYIAFU4v58rWzuHzGWOZPKdbVt8hJFOjiuVjMsW7/Uf6w5Qgvba+nprEDgHPH57Pi\n8mlcOauEimljyM3Sl5YiZ6JAF0/0RmNUVjfz3JZant9aR0NbN1mhAJfPGMudV5SzZNZ4po4d5XWZ\nIr6iQJcR45xj44EWntpwiGc319Lc0UNuOMjV55Ww7MIyPnDeePKz9UdSZLD02yPD7kBzJ7/ZcIin\nNhyiurGDrFCApXNKuXHeRK6aVaKuFJEkUaDLsOiNxli9rY5HK/fzyq5GABZPH8M/XDWDZXMnaESK\nyDBQoEtSHWo5zhNv7OeJNw/Q0NZNWVEOX752FrdWTGZSca7X5YmkNQW6JMWG/Uf52Ss1PLelFgdc\nPXs8H790Kktmj9fEVSIjRIEugxaNOV7cVscDL1dTte8ohTkh7rpyBndcNk1X4yIeUKDLWYtEYzy1\n8TD3rdlNTWMHk0fn8s0b57C8Ygp5GqUi4hn99smARaIxnt54mH97aRd7mzq5YGIh991+MdddUKrp\nZUVSgAJd+hWNOZ7ZdIifro5fkc8pK+T+OxaydE6pnsAjkkIU6HJGr+1u5Lu/28bW2lbOm1DAf35i\nIR+cU6pnXYqkIAW6nNKehna+//ttvLitnknFufzkYxdx47yJCnKRFDbgQDezIFAFHHLO3WBmY4Bf\nAOXAXmC5c+7ocBQpI+dYZy//58WdPLJ2HznhIF9ddh53XlFOTlh3c4qkurO5Qv8isA0oTGzfA6x2\nzt1rZvcktr+a5PpkhDjneHZzLd95divNHd3ctmgqX146i3H52V6XJiIDNKBAN7PJwIeA7wH/I9F8\nE7Aksb4S+BMKdF/a39TJN57ewss7G5g/uYiH77yECycVeV2WiJylgV6h/xj4ClDQp63UOVebWD8C\nlCazMBl+vdEYP3ulhp+s3kkoEOBbN87hjsvKdWeniE/1G+hmdgNQ75xbZ2ZLTrWPc86ZmTvN5+8C\n7gKYOnXqEEqVZKpuaOfLv9jIpoPHuO6CUr714QsoK9LdnSJ+NpAr9CuAD5vZ9UAOUGhmjwB1Zlbm\nnKs1szKg/lQfds7dD9wPUFFRccrQl5HjnOPRyv1873fbyA4HuO/2i/nQvDKvyxKRJOj39j7n3Nec\nc5Odc+XAx4CXnHOfAJ4BViR2WwE8PWxVSlLUt3XxqZVVfOOpLVSUj+aPX7pSYS6SRoYyDv1eYJWZ\nfQrYByxPTkkyHNbsqOfuVZvo6I7wrRvn8HeXlWtMuUiaOatAd879ifhoFpxzTcA1yS9JkikWc/z0\npV38ZPUuZpcW8G93LWZmaUH/HxQR39GdommspbOHL/9iI2t2NHDLxZP43s1z9bg3kTSmQE9TWw4d\n47OPruPIsS6+e/OFfPzSqZpISyTNKdDT0NMbD/GVX21mTF4Wqz5zGQumjva6JBEZAQr0NOKc4741\nu/nn53ey6Jwx/N+PX6xb90UyiAI9TfRGY/yvp7bwxJsHuPmiifzg1nlkh9RfLpJJFOhpoL07wuce\nXc/LOxv4wtXncvcHZ6m/XCQDKdB97sixLu58+E121rVx7y1z+dgiTa8gkqkU6D52oLmT23+2lub2\nHh765CVcNavE65JExEMKdJ/a19TB7Q9U0tbVy2OfXsz8KcVelyQiHlOg+9CehnZuf2AtPZEYj316\nseYuFxFAge47O+vauP2BSsDx+F2LOW9CYb+fEZHMoED3ka2HW/nEg5WEAsZjn17MueM1J4uIvEOB\n7hM1jR3c8WAl2aEAj316MeeMy/O6JBFJMf3Ohy7eO3Ksi0/8rBIHPPL3lyrMReSUFOgprqWzh797\nqJKWzh5W3rmIGSX5XpckIilKXS4prLMnwp0Pv8nexk4e/u+XMHeyRrOIyOnpCj1F9URifPaR9Ww6\n0MJPb1vA5TPGeV2SiKQ4XaGnIOccX/31Zv68s4EffGQuyy6c4HVJIuIDukJPQf/x5z38ZsMh7l46\ni7+9RHOziMjAKNBTzAtb6/jRH3fw4fkT+cIHzvW6HBHxEQV6Ctl+pJUvPbGBuZOK+OGt8zQFroic\nFQV6imju6OHvV1aRlx3i/jsqyAnr4RQicnb0pWgKiI9oWUd9WzerPnMZE4pyvC5JRHxIV+gp4Du/\nfZvKmmZ+dOs8LtI0uCIySAp0jz276TCPrN3PZ66czk0XTfK6HBHxMQW6h/Y1dfC1J9/i4qnF/M/r\nZntdjoj4nALdIz2RGP/4+AYCBj+9bQHhoE6FiAyNvhT1yA/+sJ3NB4/x/+5YyOTRo7wuR0TSgC4L\nPbB6Wx0PvlrDisumcd0Fuq1fRJJDgT7Cao8d5+5fbmJOWSFfu/58r8sRkTSiQB9BsZjjS09spCcS\n499vX6Cbh0QkqdSHPoJWvr6XyppmfnjrPKbrQRUikmS6Qh8h+5o6+OEfdrBkdgkfXTjZ63JEJA31\nG+hmlmNmb5jZJjN728y+nWgfY2YvmNmuxHL08JfrT7GY4yu/2kwoYHz/lrmadEtEhsVArtC7gQ84\n5+YDFwHLzGwxcA+w2jk3E1id2JZTeLRyH5U1zXzjhvMpK8r1uhwRSVP9BrqLa09shhMvB9wErEy0\nrwRuHpYKfe5Acyfff2477585juUVU7wuR0TS2ID60M0saGYbgXrgBedcJVDqnKtN7HIEKB2mGn3r\nxKPkAmbc+xHNby4iw2tAge6cizrnLgImA4vM7MKT3nfEr9rfw8zuMrMqM6tqaGgYcsF+8tgb+3lt\nTxNfv/58JhWrq0VEhtdZjXJxzrUAa4BlQJ2ZlQEklvWn+cz9zrkK51xFSUnJUOv1jfq2Lr7/++1c\nce5YblukrhYRGX4DGeVSYmbFifVcYCmwHXgGWJHYbQXw9HAV6Uf3PrednkiM796sUS0iMjIGcmNR\nGbDSzILE/wJY5Zz7rZm9Dqwys08B+4Dlw1inr1TtbebJ9Yf4/NUzOGdcntfliEiG6DfQnXObgQWn\naG8CrhmOovwsGnP876ffpqwoh89ffa7X5YhIBtGdokn22Bv72Vrbyjc+NIdRWZpZQURGjgI9iY52\n9PAvz+/gsuljuX6upsUVkZGlQE+iHz2/g7auCN++6QJ9ESoiI06BniRvHTzG42/s55OXlzOrtMDr\nckQkAynQk8A5xzef2cLYvGy+eO1Mr8sRkQylQE+CP2w5wvr9LXxl2WwKc8JelyMiGUqBPkSRaIwf\nPb+DmePz+cjFmudcRLyjQB+iX68/SHVDB/903WyCAX0RKiLeUaAPQVdvlB+/uIsFU4tZOkeTTYqI\ntxToQ/Dz1/dRe6yLry47T8MURcRzCvRBau3q5b4/7eaqWSUsnj7W63JERBTog3X/n6tp6ezln66b\n7XUpIiKAAn1Q6tu6ePDVGm6cP5ELJxV5XY6ICKBAH5R/f2k3vdEYdy+d5XUpIiJ/pUA/SwePdvL4\nG/tZfskUyjXXuYikEAX6WXrg5Wqcgy9ornMRSTEK9LPQ2N7NE28e4JaLJzFRD30WkRSjQD8LD71a\nQ080xj9cNcPrUkRE3kOBPkDHjvfy89f3cf3cMqaX5HtdjojIeyjQB+iRtfto647wWV2di0iKUqAP\nwPGeKA++WsOS2SUady4iKUuBPgBPvLmf5o4ePq+RLSKSwhTo/eiJxHjg5WoWlY/hkvIxXpcjInJa\nCvR+PLXxEIePdfG5q9V3LiKpTYF+BtGY4z//tIcLJhZy1awSr8sRETkjBfoZvLS9nurGDj67ZIbm\nOxeRlKdAP4OHX6uhrCiHZRdM8LoUEZF+KdBPY1ddG3/Z3cQnFk8jFNT/JhFJfUqq01j5+l6yQgFu\nWzTV61JERAZEgX4Kx4738uT6Q3x4/kTG5GV5XY6IyIAo0E/hl1UH6OyJ8snLy70uRURkwBToJ4nF\nHD9fu4+KaaN1m7+I+IoC/SR/2lnPvqZOVujqXER8pt9AN7MpZrbGzLaa2dtm9sVE+xgze8HMdiWW\no4e/3OH38Gv7KC3MZtmFGqooIv4ykCv0CHC3c24OsBj4vJnNAe4BVjvnZgKrE9u+tqehnZd3NvDx\nS6cR1lBFEfGZflPLOVfrnFufWG8DtgGTgJuAlYndVgI3D1eRI+W/XttLVlBDFUXEn87qMtTMyoEF\nQCVQ6pyrTbx1BChNamUjrK2rl1+tO8gN88ooKcj2uhwRkbM24EA3s3zg18CXnHOtfd9zzjnAneZz\nd5lZlZlVNTQ0DKnY4fTsplo6eqLccdk0r0sRERmUAQW6mYWJh/mjzrknE811ZlaWeL8MqD/VZ51z\n9zvnKpxzFSUlqTtj4aqqA8wqzeeiKcVelyIiMigDGeViwIPANufcv/Z56xlgRWJ9BfB08ssbGbvq\n2th4oIXlFVM0q6KI+FZoAPtcAdwBvGVmGxNtXwfuBVaZ2aeAfcDy4Slx+P1y3UFCAePmBZO8LkVE\nZND6DXTn3KvA6S5br0luOSOvNxrjyfUHueb88YzL15ehIuJfGT/Yes32ehrbe1heMcXrUkREhiTj\nA31V1UFKCrL1iDkR8b2MDvT6ti7W7Kjnlosn6SEWIuJ7GZ1iT204RDTm+OhCdbeIiP9lbKA751hV\ndZCF00Zz7vh8r8sRERmyjA30DQda2F3fzkcXTva6FBGRpMjYQP9l1UFyw0E+NK/M61JERJIiIwP9\neE+UZzcd5vq5ZRTkhL0uR0QkKTIy0J/feoT27ggfrVB3i4ikj4wM9Gc3HWZCYQ6Lysd4XYqISNJk\nXKAf6+zlzzsbuGFeGYGAJuISkfSRcYH+x61H6I06bpg/0etSRESSKuMC/dlNh5kyJpf5k4u8LkVE\nJKkyKtCb2rt5bU8TN86bqHnPRSTtZFSgP7flCNGY44Z56m4RkfSTUYH+282HmVGSx/llBV6XIiKS\ndBkT6HWtXVTWNHODultEJE1lTKD/bnMtzsGN83Wrv4ikp4wJ9N9uPsx5Ewo4d7y6W0QkPWVEoB9o\n7mT9/hZu1NhzEUljGRHov3urFoAbNbpFRNJYRgT6bzcfZv7kIqaOHeV1KSIiwybtA72msYMth1rV\n3SIiaS/tA/25LfHuluvnanSLiKS3tA/0F7bWMXdSEROLc70uRURkWKV1oNe3dbHxQAvXnl/qdSki\nIsMurQP9pW31OAdL5yjQRST9pXWgv7itjknFuZq7RUQyQtoGemdPhFd2NbJ0TqnmbhGRjJC2gf7q\nrka6IzH1n4tIxkjbQH9xWx0FOSEuna4HQYtIZkjLQI/GHKu31bNk9njCwbQ8RBGR90jLtNt44ChN\nHT0a3SIiGaXfQDezh8ys3sy29GkbY2YvmNmuxHL08JZ5dp7fWkcoYFw1q8TrUkRERsxArtAfBpad\n1HYPsNo5NxNYndhOGS9urWPx9LEU5Ya9LkVEZMT0G+jOuZeB5pOabwJWJtZXAjcnua5Bq25oZ09D\nB9eeP97rUkRERtRg+9BLnXO1ifUjQMp0Vr+4rQ6Aa9V/LiIZZshfijrnHOBO976Z3WVmVWZW1dDQ\nMNQf168XttZxflkhk0dr7nMRySyDDfQ6MysDSCzrT7ejc+5+51yFc66ipGR4v6Rsau9m3b6jGt0i\nIhlpsIH+DLAisb4CeDo55QzNmh0NxBws1d2hIpKBBjJs8XHgdWC2mR00s08B9wJLzWwXcG1i23N/\n3tnAuPxsLpxU6HUpIiIjLtTfDs65207z1jVJrmVIYjHHq7sauHr2eE3GJSIZKW3uFH37cCtHO3u5\nUjcTiUiGSptAf3lXfATNFeeO87gSERFvpE+g72xgTlkhJQXZXpciIuKJtAj09u4I6/cfVXeLiGS0\ntAj0tXua6I06rpyp7hYRyVxpEeiv7GogNxxkYXlKTfooIjKi0iTQG1k8fQzZoaDXpYiIeMb3gX6g\nuZPqxg7eP1P95yKS2Xwf6K/sagTgylnqPxeRzJYGgd7AxKIcZpTke12KiIinfB3okWiMv+xu5P0z\nS3S7v4hkPF8H+qaDx2jtivB+dbeIiPg70F/Z1YAZvE+3+4uI+D3QG5k3uZjiUVlelyIi4jnfBvqx\n471sPNCiu0NFRBJ8G+iv72kkGnOav0VEJMG3gf6X3U3kZQW5aEqx16WIiKQE3wb62uomKsrHEA76\n9hBERJLKl2nY1N7Nrvp2Lp0+xutSRERShi8D/Y2aZgAuPWesx5WIiKQOXwZ6ZU0zueEgcycVeV2K\niEjK8GWgr61u4uJpxWSFfFm+iMiw8F0itnT2sKOuTd0tIiIn8V2gv7n3KM7BpefoC1ERkb58F+iV\n1U1khQLM1/hzEZF38V+g1zSzYEoxOWE9bk5EpC9fBXprVy9vHz7GpdPVfy4icjJfBfq6vUeJOVis\n/nMRkffwVaCvrWkiHDQWTB3tdSkiIinHV4FeWd3M/MnF5Gap/1xE5GS+CfSO7ghvHTqm+VtERE7D\nN4G+bt9RojGnG4pERE7DN4FeWdNEMGAsnKb+cxGRUxlSoJvZMjPbYWa7zeyeZBV1KpXVzcydVERe\ndmg4f4yIiG8NOtDNLAjcB/wNMAe4zczmJKuwvo73RNl0sEX95yIiZzCUK/RFwG7nXLVzrgd4Argp\nOWW924b9R+mNOhar/1xE5LSGEuiTgAN9tg8m2pJubU0zAYOKcvWfi4iczrB/KWpmd5lZlZlVNTQ0\nDOq/Mak4h1sXTqYgJ5zk6kRE0sdQAv0QMKXP9uRE27s45+53zlU45ypKSkoG9YP+9pKp/PDW+YOr\nUkQkQwwl0N8EZprZOWaWBXwMeCY5ZYmIyNka9BhA51zEzL4A/BEIAg85595OWmUiInJWhjSo2zn3\ne+D3SapFRESGwDd3ioqIyJkp0EVE0oQCXUQkTSjQRUTShAJdRCRNmHNu5H6YWQOwb5AfHwc0JrEc\nL+lYUk+6HAfoWFLVUI5lmnOu3zszRzTQh8LMqpxzFV7XkQw6ltSTLscBOpZUNRLHoi4XEZE0oUAX\nEUkTfgr0+70uIIl0LKknXY4DdCypatiPxTd96CIicmZ+ukIXEZEz8EWgj+TDqIebme01s7fMbKOZ\nVXldz0CZ2UNmVm9mW/q0jTGzF8xsV2Lpi0dKneZYvmVmhxLnZaOZXe9ljQNhZlPMbI2ZbTWzt83s\ni4l2352XMxyLH89Ljpm9YWabEsfy7UT7sJ+XlO9ySTyMeiewlPhj7t4EbnPObfW0sEEys71AhXPO\nV2NrzexKoB34L+fchYm2HwLNzrl7E3/RjnbOfdXLOgfiNMfyLaDdOffPXtZ2NsysDChzzq03swJg\nHXAz8El8dl7OcCzL8d95MSDPOdduZmHgVeCLwC0M83nxwxX6iD2MWk7POfcy0HxS803AysT6SuK/\ngCnvNMfiO865Wufc+sR6G7CN+HN9fXdeznAsvuPi2hOb4cTLMQLnxQ+BPmIPox4hDnjRzNaZ2V1e\nFzNEpc652sT6EaDUy2KS4B/NbHOiSybluyn6MrNyYAFQic/Py0nHAj48L2YWNLONQD3wgnNuRM6L\nHwI93bzPOXcR8DfA5xP//Pc9F++7S+3+uzP7D2A6cBFQC/yLt+UMnJnlA78GvuSca+37nt/OyymO\nxZfnxTkXTfyeTwYWmdmFJ70/LOfFD4E+oIdR+4Vz7lBiWQ/8hniXkl/VJfo+T/SB1ntcz6A55+oS\nv4Qx4AF8cl4SfbS/Bh51zj2ZaPbleTnVsfj1vJzgnGsB1gDLGIHz4odAT5uHUZtZXuILH8wsD/gg\nsOXMn0ppzwArEusrgKc9rGVITvyiJfw3fHBeEl++PQhsc879a5+3fHdeTncsPj0vJWZWnFjPJT6g\nYzsjcF5SfpQLQGKo0o9552HU3/O4pEExs+nEr8oh/jzXx/xyLGb2OLCE+IxxdcA3gaeAVcBU4rNo\nLnfOpfyXjac5liXE/1nvgL3AZ/r0d6YkM3sf8ArwFhBLNH+deN+zr87LGY7lNvx3XuYR/9IzSPyi\neZVz7jtmNpZhPi++CHQREemfH7pcRERkABToIiJpQoEuIpImFOgiImlCgS4ikiYU6CIiaUKBLiKS\nJhToIiJp4v8DptIFXIJspDYAAAAASUVORK5CYII=\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "%matplotlib inline\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "dt = 0.5 # step size in s\n", "u=[0.0] # create a list for the velocity array. This contains the initial condition\n", "t=[0.0] # create a list for the time array. starts at t = 0.0\n", "tend = 30.0 # set end time\n", "\n", "c = 12.5 # drag coefficientkg/s\n", "m = 60 # object's mass, kg\n", "g = 9.81 # gravitational acceleration m/s/s\n", "\n", "# t[-1] returns the last element in the list\n", "while t[-1] < tend:\n", " unp1 = u[-1] + dt * (g - c/m*u[-1]) # time advance\n", " u.append(unp1)\n", " t.append(t[-1] + dt)\n", " \n", "# tplot = np.linspace(t0,tend,len(u)) # create an equally space array for time. This will be used for plotting.\n", "plt.plot(t,u)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "subslide" } }, "source": [ "## Method 2: Using Numpy Arrays" ] }, { "cell_type": "code", "execution_count": 138, "metadata": { "collapsed": false, "slideshow": { "slide_type": "fragment" } }, "outputs": [ { "data": { "text/plain": [ "[]" ] }, "execution_count": 138, "metadata": {}, "output_type": "execute_result" }, { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD8CAYAAABn919SAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHEBJREFUeJzt3XmYXHWd7/H3t6r3JUunu5PO2glZIAlJCE3YBNlBYRL0\nahQFo6Ogd1BxvHOVO/fOM+qoj+MMM947zjBGQYMDSBx0CIzKEkgybAmdjaxk3zrd6U46vaaXWn73\njyqSCOklvZ06pz6v56nnrEV9Dyf1yS+/Ouf8zDmHiIj4X8jrAkREZGAo0EVEAkKBLiISEAp0EZGA\nUKCLiASEAl1EJCAU6CIiAaFAFxEJCAW6iEhAZAzlhxUXF7vy8vKh/EgREd9bv379cedcSU/7DWmg\nl5eXU1lZOZQfKSLie2Z2sDf7qctFRCQgFOgiIgGhQBcRCQgFuohIQCjQRUQCQoEuIhIQCnQRkYAY\n0uvQRUT6Ix53ROOOWNwRjceTU3dmGnPEnCMWj59e/76XczjH6fl4cn3cQdy9O598xUnu74jFE9ud\nO7Nv3JFcPrPOuUSdjjP74BwfmT+eycX5g/r/R4EuIu8TjzvaozHaOmO0R+O0R2LJV5yOSIyOaJyO\naHIaOWs+GqczGqczlphGYu9fjsRccnpmvjOaCODoWeti8cQ0mgzraCyeCEcfMoP5k0Yq0EWke845\nTnXGaG6P0tIRSU6jtLRHae6I0toR5VRnjJaOKKc6orR0xDjVmVjX1hnjVOTMfFskMe2IxvtVU2bY\nyAqHyMoIkZmcZmWEyAonljPDRmY4RF5WmIxQRnJdiIywkRFKbH93PiNkZCTfEw4l3hcOGRmhxHJi\nPrFfKLn+9NTO7Hf2unDICBmEQkb49LIRCkHYEvuGLLHNDMKh5PSsbSEDI/GexHJiH7MzyyEDMxug\nM90zBbpIiojG4jS0RTjZ2snJUxHqWztpOJWYb2yL0NjWmZwmXk1tUZraEwEe60XTNSNk5GdnkJ8V\nJi87g7ysMHlZYUoLc8jNCpOXGSY3K0xuZpic5HxORigxzQyTnREmJzP0/mlm+HR4ZydDOxQauhCT\nMxToIoMoGotzvKWTuuYO6lraqWvuoLapg7qWDk60dHK8pYMTrZ2caOmgoS2C6yKXs8IhhuVmMiIv\nk+G5mZQW5jC1JIPCnEyG5SamhTnJaXYGBTkZFGSfeeVlJwJZgk2BLtJH8bijtrmDqoY2jja0nZ7W\nNLZzrKmdmqZEgJ+r8Tw8N5PigixGFWQzrbSAK6YUUZSfzaj8LEbmZzEyL5OReVkU5WcxIi+T3Mzw\nkP7TXfxJgS7SjdaOKAdOtHLwxCkO1Sdeh5OvqoY2IrE/TuvCnAzKhucwelgO00cXMiY5X1qYTemw\nHEoKsykuyFJrWQaFAl3SnnOOo43t7D7WzJ7aFvbWtbD/eCv7j7dyrKnjj/YdmZfJxKI8Zo8bzocu\nLmPciFzGjchl7Ihcxo7IoTAn06OjEFGgS5o52drJjpomdlQ3s6O66XSIt3bGTu9TlJ/F5OJ8rplW\nwuTifCYX5zNpVB4Ti/IU2JLSFOgSWHXNHWypamDz4Ua2VDWy/WgTNU3tp7cXF2Rz4ZhCPl4xgWmj\nC5hWWsi00gJG5md5WLVI3ynQJRA6ojG2VjWy/uBJNhxs4O0jDRxtTIS3GUwtSfzweFHZsNOvksJs\nj6sWGVgKdPGl5vYIbx2oZ+2+eioPnmTLkUY6Y4mbYSaNyqOivIg544czZ/wIZo0dRn62/qhL8OlP\nufhCW2eMyoP1vL73BG/sPcGWqkZicUdWOMTF44fz2avLuXTSSOZPHKmWt6QtBbqkJOcc+4+3suqd\nOlbtqmPtvhN0RONkhIy5E0bwZ9ddwJUXjGL+xJHkZOoSQBFQoEsKicbivHXgJM9vq+HlnbUcqj8F\nwJTifD51+USunV7CgvIidZ+IdEHfDPFURzTGa3uO84etNby0o5b61k6yM0JcPbWYL1wzmeumlzJx\nVJ7XZYr4ggJdhlws7nhz3wn+Y2MVf9haQ3NHlMLsDG64qJRbZ43hg9NL1AoX6QN9a2TIbK1q5Lcb\nq3h281FqmzsoyM7g1lljuGNuGVddMEq3w4v0kwJdBlVze4RnNh3lyXWH2Ha0icywcd2MUhbNG8tN\nF43WD5oiA0iBLgPOOcemww08ue4Qz26upi0S48IxhXxn0SwWzh3LiDzdiSkyGBToMmAisTi/31rD\nI6/uZ/PhBvKywiyaN5ZPLpjI3PHD9fhXkUGmQJd+a2yL8Kt1h/jF6weobmxncnE+31k0i4/OH0+B\nftwUGTL6tkmf1TV3sHTNXh5fe4hTnTGunDKK7945m+tnlGoIMhEPKNDlvNU2tfOTNft4fO1BOqNx\nFs4dy73XTmHW2OFelyaS1hTo0mu1Te38y6q9PLnuENG448554/jyDVOZXJzvdWkiggJdeqG1I8rS\nNftYumYfkVicj84fx/3XT2XSKAW5SCpRoEuXorE4yyuP8I8v7aKuuYPb55TxjVtnKMhFUlSvA93M\nwkAlUOWcu8PMioCngHLgALDYOXdyMIqUobd6Vx3ffW47u2tbqJg0kp/ccynzJ470uiwR6UboPPZ9\nANhx1vKDwErn3DRgZXJZfK6msZ0/e3w9Sx5dRyQW51/vns+vv3SlwlzEB3rVQjez8cDtwPeArydX\nLwKuS84vA1YB3xzY8mSoRGNxHnvjIA+98A7RuOMvbpnOvddO0fNVRHykt10uPwK+ARSetW60c646\nOV8DjB7IwmTobD7cwF/+dgvbjjZx3YwSvrNwth5ZK+JDPQa6md0B1Drn1pvZdefaxznnzMx18f77\ngPsAJk6c2I9SZaB1RuP835W7eHjVXkoKs3n40/O5bfYY3aIv4lO9aaFfDSw0sw8DOcAwM/s34JiZ\nlTnnqs2sDKg915udc0uBpQAVFRXnDH0Zejtrmvj6U5vZXt3E4orx/NUdMynMyfS6LBHphx5/FHXO\n/S/n3HjnXDnwSeBl59zdwApgSXK3JcAzg1alDJhY3PGT1XtZ+E+vUdvczk8/U8EPPzZXYS4SAP25\nDv0HwHIz+zxwEFg8MCXJYKlpbOerT25k3YF6bp01mu9/5GJGFWR7XZaIDJDzCnTn3CoSV7PgnDsB\n3DjwJclgeHX3cR741UbaIjEe+vhcPjp/nPrKRQJGd4oGXDzu+PEre/jHl3YxtaSAh+++lKmlBV6X\nJSKDQIEeYCdbO/naU5tYvauORfPG8v2PXKzBl0UCTN/ugNp+tIl7H6ukrrmDv7lzNndfPlFdLCIB\np0APoFd21vLlJzZQmJPJ8i9dybwJI7wuSUSGgAI9YB574wDfWrGNi8qG8ciSyxgzPMfrkkRkiCjQ\nAyIWd3zvP3fw6Gv7ufHCUv7fXZeov1wkzegbHwCnOqN89clNvLTjGJ+9qpy/umMmYY3pKZJ2FOg+\n19Qe4U9//hYbDp3k2wtnseSqcq9LEhGPKNB9rOFUJ595dB3bjzbxT3fN5/Y5ZV6XJCIeUqD7VF1z\nB/c8spZ9x1v5yT2XcuNFenqxSLpToPtQdWMbn/7pWqob23l0yWV8YFqx1yWJSApQoPvM4fpTfOpn\nb9LQGuGxzy/gsvIir0sSkRShQPeRY03tfOpnb9LUFuXxey9nznjdMCQiZ5zPINHioZOtndz9s7XU\nt3Ty2J8uUJiLyPuohe4DLR1RPvvzdRysP8Wyzy1grm7lF5FzUAs9xbVHYnxh2VtsO9rEw5+ez5UX\njPK6JBFJUWqhp7BILM6Xn9jA2v31/OgT83Rpooh0Sy30FOWc48Gnt/DSjlr+ZtFsFs0b53VJIpLi\nFOgp6l9W7eXpDUf485umc/cVk7wuR0R8QIGegv6wtZq/e/4dFs0by1dvnOp1OSLiEwr0FLO1qpE/\nf2ozl0wcwd/+tzkaZUhEek2BnkKONbXzhWWVFOVnsfSeCnIyw16XJCI+oqtcUkRbZ4x7H6ukqT3C\n0//9KkoKs70uSUR8RoGeApxzfOPpt9lS1chP76ngorJhXpckIj6kLpcU8Ms3D/Ls5qP8z1tncNNM\nXWsuIn2jQPfYliONfPe5HdxwYSlfuvYCr8sRER9ToHuoqT3C/U9soLggi4c+PpeQxgEVkX5QH7pH\nnHN889/f5mhDG0998UpG5md5XZKI+Jxa6B5Z9voBfr+1hm/cNoNLJ430uhwRCQAFugfePtLA9363\ngxsvLOULH5jidTkiEhAK9CHW0hHl/ic2UFqYw0OL1W8uIgNHfehD7Pu/28GRk238+otXMiJP/eYi\nMnDUQh9Ca3bV8cTaQ9x7zRQqNLiziAywHgPdzHLMbJ2ZbTazbWb27eT6IjN70cx2J6f6Za8bTe0R\nvvn021xQks/Xb57udTkiEkC9aaF3ADc45+YC84DbzOwK4EFgpXNuGrAyuSxd+O5z2znW1M5Di+fp\noVsiMih6DHSX0JJczEy+HLAIWJZcvwy4c1AqDIBXdtayvPIIX/rgBczTAM8iMkh61YduZmEz2wTU\nAi8659YCo51z1cldagA9hOQcGk9FePA3bzNjdCEP3DTN63JEJMB6FejOuZhzbh4wHlhgZrPfs92R\naLW/j5ndZ2aVZlZZV1fX74L95tvPbuN4SycPLZ5Ldoa6WkRk8JzXVS7OuQbgFeA24JiZlQEkp7Vd\nvGepc67COVdRUlLS33p9ZfWuOn6zsYr7r5/K7HHDvS5HRAKuN1e5lJjZiOR8LnAzsBNYASxJ7rYE\neGawivSj9kiMv35mK1OK87n/ej1FUUQGX29uLCoDlplZmMRfAMudc8+Z2RvAcjP7PHAQWDyIdfrO\nT9fs48CJU/zy8wvU1SIiQ6LHQHfOvQ1cco71J4AbB6Movztcf4ofv7KH2y8u45pp6dXNJCLe0Z2i\ng+Dbz24nHDL+zx0XeV2KiKQRBfoAW7njGC/tOMYDN06jbHiu1+WISBpRoA+g9kiMbz27jamlBXzu\n6slelyMiaUZPWxxAD6/ay+H6Np6493KyMvR3pYgMLaXOADl04hQPr97LwrljueqCYq/LEZE0pEAf\nID98fidhM/737fohVES8oUAfAG8faeC5t6v5wjWTGT0sx+tyRCRNKdD7yTnHD36/k6L8LO67VuOD\nioh3FOj9tGb3cV7fe4Kv3DCVwpxMr8sRkTSmQO+HeDzROp9QlMunLp/odTkikuYU6P3wzOYqdlQ3\n8Re3zNDzWkTEcwr0PuqIxvj753cxa+ww/mTOWK/LERFRoPfVv715iKqGNh780IWEQuZ1OSIiCvS+\naGqP8OOXd3PNtGI9TVFEUoYCvQ9+9l/7OXkqwjdvu9DrUkRETlOgn6em9gi/eG0/t80ao2HlRCSl\nKNDP0y/fOEhTe5Qv3zDV61JERP6IAv08nOqM8sir+7l+Rola5yKSchTo5+GJtYeob+1U61xEUpIC\nvZfaIzGWrtnHlVNGcemkIq/LERF5HwV6L/16/RFqmzv4ilrnIpKiFOi9EInF+ddVe5k/cQRXXjDK\n63JERM5Jgd4Lv91YRVVDG1+5YRpmuitURFKTAr0Hsbjj4VV7mTV2GNfN0F2hIpK6FOg9+M8t1ew/\n3spXbpiq1rmIpDQFejecc/xk9V6mlhZwy8wxXpcjItItBXo33jpwkm1Hm/jc1eV6oqKIpDwFejd+\n/tp+hudm8tFLxntdiohIjxToXahqaOP5bTV8csEEcrM0GpGIpD4Fehcee+MAZsZnriz3uhQRkV5R\noJ/Dqc4ov1p3mFtnjWbciFyvyxER6RUF+jn8dmMVjW0RPnf1ZK9LERHpNQX6ezjn+MVrB5g9bhgV\nk0Z6XY6ISK/1GOhmNsHMXjGz7Wa2zcweSK4vMrMXzWx3chqI9Ht1z3F217bwuasm60YiEfGV3rTQ\no8D/cM7NBK4A7jezmcCDwErn3DRgZXLZ937+2gGKC7K5Y26Z16WIiJyXHgPdOVftnNuQnG8GdgDj\ngEXAsuRuy4A7B6vIobL/eCsv76zl05dPJDtDlyqKiL+cVx+6mZUDlwBrgdHOuerkphpg9IBW5oFl\nrx8gM2x8+oqJXpciInLeeh3oZlYAPA18zTnXdPY255wDXBfvu8/MKs2ssq6url/FDqa2zhhPrz/C\n7ReXUVqY43U5IiLnrVeBbmaZJML8cefcb5Krj5lZWXJ7GVB7rvc655Y65yqccxUlJan7+Nnfb62m\nuSPKJxeodS4i/tSbq1wMeATY4Zz7h7M2rQCWJOeXAM8MfHlD51dvHaZ8VB6XT9Z4oSLiT71poV8N\n3APcYGabkq8PAz8Abjaz3cBNyWVf2lfXwrr99Sy+bIIuVRQR38roaQfn3KtAVyl348CW443llUcI\nh4yPzddTFUXEv9L+TtFILM7TG45w/YwSSofpx1AR8a+0D/RXdtZS19zBJy7Tj6Ei4m9pH+jLKw9T\nUpjN9RoAWkR8Lq0D/VhTOy/vrOVjl44nI5zW/ytEJADSOsX+ff0R4g4WV0zwuhQRkX5L20CPxx3L\nKw9z+eQiJhfne12OiEi/pW2gr91fz8ETp/jEZWqdi0gwpG2gP/XWIQpzMvjQbD0mV0SCIS0DvaUj\nyu+31rBw7lhys/SYXBEJhrQM9Be21dARjfORS8Z5XYqIyIBJy0Bfsfko40bkMn9iIEbNExEB0jDQ\n61s7eXX3ce6YW0YopAdxiUhwpF2g/25LNdG4Y+HcsV6XIiIyoNIu0FdsPsrU0gJmlg3zuhQRkQGV\nVoFe3djGWwfqWTh3rJ57LiKBk1aB/tzmapxD3S0iEkhpFegrNh9lzvjhlOtWfxEJoLQJ9P3HW9lS\n1ajWuYgEVtoE+opNRzGDO+Yo0EUkmNIi0J1zrNhcxYLyIsYM1zBzIhJMaRHo26ub2FvXysJ5ap2L\nSHClRaCv2HyUjJDxYT1ZUUQCLPCBHo87nttczbXTSxiZn+V1OSIigybwgb75SANVDW3cMUetcxEJ\ntsAH+gvbj5ERMm68cLTXpYiIDKrAB/qL249xxZRRDM/L9LoUEZFBFehA31vXwp7aFm6Zpda5iARf\noAP9xe3HALjpIgW6iARfoAP9hW01XDxuOGNH5HpdiojIoAtsoNc2t7PxcAO3zFTrXETSQ2ADfeWO\nWpyDW2aN8boUEZEhEdhAf2FbDZNG5TF9dIHXpYiIDIlABnpLR5TX9pzg5otGa2QiEUkbPQa6mT1q\nZrVmtvWsdUVm9qKZ7U5ORw5umedn9Tt1dMbi6m4RkbTSmxb6L4Db3rPuQWClc24asDK5nDJe2F5D\nUX4Wl05Kqb9nREQGVY+B7pxbA9S/Z/UiYFlyfhlw5wDX1WeRWJyXd9Zy00WlhEPqbhGR9NHXPvTR\nzrnq5HwNkDLXBq7dV09ze5RbZqq7RUTSS79/FHXOOcB1td3M7jOzSjOrrKur6+/H9eiF7TXkZob5\nwLTiQf8sEZFU0tdAP2ZmZQDJaW1XOzrnljrnKpxzFSUlJX38uN5xzvHCtmNcO72YnMzwoH6WiEiq\n6WugrwCWJOeXAM8MTDn9s7WqiZqmdnW3iEha6s1li08CbwAzzOyImX0e+AFws5ntBm5KLntu1Tu1\nmMF1Mwb3XwIiIqkoo6cdnHN3dbHpxgGupd9W76rj4nHDGVWQ7XUpIiJDLjB3ija2Rdh4uIEPTlfr\nXETSU2AC/fU9x4nFHdcq0EUkTQUm0NfsrqMwJ4NLJozwuhQREU8EItCdc6x+p46rLygmIxyIQxIR\nOW+BSL+9dS0cbWzng7q6RUTSWCACfdU7iTtQ1X8uIuksEIG+ZvdxppYWME5jh4pIGvN9oLdHYqzd\nd4Jrp6l1LiLpzfeBvnZ/PR3RuPrPRSTt+T7QV79TR3ZGiMsnF3ldioiIp3wf6Gt217FgcpGerigi\nac/XgV7V0Mae2hbd7i8igs8Dfc2uxOWKCnQREZ8H+up36hg7PIeppQVelyIi4jnfBnokFue1Pce5\ndnoJZhoMWkTEt4G+6XADzR1RdbeIiCT5NtBf3X2ckMFVUzUYtIgI+DjQ39x3gtnjhjM8N9PrUkRE\nUoIvA709EmPj4QaumDLK61JERFKGLwN946EGOqNxrpiiu0NFRN7ly0B/c98JQgYV5Qp0EZF3+TbQ\nZ48bzrAc9Z+LiLzLd4Gu/nMRkXPzXaCr/1xE5Nx8F+jqPxcROTdfBrr6z0VE3s9Xgd4eibHxkPrP\nRUTOxVeBvuHQSTpj6j8XETkXXwX6m/vq1X8uItIFnwW6+s9FRLrim0Bvj8TYpP5zEZEu+SbQ1X8u\nItI93wS6+s9FRLrXr0A3s9vM7B0z22NmDw5UUeei/nMRke71OdDNLAz8M/AhYCZwl5nNHKjCzqb+\ncxGRnvWnhb4A2OOc2+ec6wR+BSwamLL+mPrPRUR61p9AHwccPmv5SHLdgFP/uYhIzwb9R1Ezu8/M\nKs2ssq6urk//jXEjcvj4pRPUfy4i0o3+BHoVMOGs5fHJdX/EObfUOVfhnKsoKSnp0wd94rKJ/O3H\n5vStShGRNNGfQH8LmGZmk80sC/gksGJgyhIRkfOV0dc3OueiZvZl4HkgDDzqnNs2YJWJiMh56XOg\nAzjnfgf8boBqERGRfvDNnaIiItI9BbqISEAo0EVEAkKBLiISEAp0EZGAMOfc0H2YWR1wsI9vLwaO\nD2A5XtKxpJ6gHAfoWFJVf45lknOuxzszhzTQ+8PMKp1zFV7XMRB0LKknKMcBOpZUNRTHoi4XEZGA\nUKCLiASEnwJ9qdcFDCAdS+oJynGAjiVVDfqx+KYPXUREuuenFrqIiHTDF4E+lINRDzYzO2BmW8xs\nk5lVel1Pb5nZo2ZWa2Zbz1pXZGYvmtnu5HSklzX2VhfH8i0zq0qel01m9mEva+wNM5tgZq+Y2XYz\n22ZmDyTX++68dHMsfjwvOWa2zsw2J4/l28n1g35eUr7LJTkY9S7gZhLD3L0F3OWc2+5pYX1kZgeA\nCuecr66tNbNrgRbgMefc7OS6HwL1zrkfJP+iHemc+6aXdfZGF8fyLaDFOff3XtZ2PsysDChzzm0w\ns0JgPXAn8Fl8dl66OZbF+O+8GJDvnGsxs0zgVeAB4KMM8nnxQwt9yAajlq4559YA9e9ZvQhYlpxf\nRuILmPK6OBbfcc5VO+c2JOebgR0kxvX13Xnp5lh8xyW0JBczky/HEJwXPwT6kA1GPUQc8JKZrTez\n+7wupp9GO+eqk/M1wGgvixkAXzGzt5NdMinfTXE2MysHLgHW4vPz8p5jAR+eFzMLm9kmoBZ40Tk3\nJOfFD4EeNB9wzs0DPgTcn/znv++5RN9davffde9hYAowD6gGHvK2nN4zswLgaeBrzrmms7f57byc\n41h8eV6cc7Hk93w8sMDMZr9n+6CcFz8Eeq8Go/YL51xVcloL/JZEl5JfHUv2fb7bB1rrcT195pw7\nlvwSxoGf4pPzkuyjfRp43Dn3m+RqX56Xcx2LX8/Lu5xzDcArwG0MwXnxQ6AHZjBqM8tP/uCDmeUD\ntwBbu39XSlsBLEnOLwGe8bCWfnn3i5b0EXxwXpI/vj0C7HDO/cNZm3x3Xro6Fp+elxIzG5GczyVx\nQcdOhuC8pPxVLgDJS5V+xJnBqL/ncUl9YmZTSLTKITGe6xN+ORYzexK4jsQT444Bfw38B7AcmEji\nKZqLnXMp/2NjF8dyHYl/1jvgAPDFs/o7U5KZfQD4L2ALEE+u/ksSfc++Oi/dHMtd+O+8zCHxo2eY\nRKN5uXPuO2Y2ikE+L74IdBER6ZkfulxERKQXFOgiIgGhQBcRCQgFuohIQCjQRUQCQoEuIhIQCnQR\nkYBQoIuIBMT/BxSNHAFIHVXIAAAAAElFTkSuQmCC\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "%matplotlib inline\n", "import matplotlib.pyplot as plt\n", "import numpy as np\n", "dt = 0.5 # step size in s\n", "\n", "t0 = 0.0 # set initial time\n", "tend = 30.0 # set end time\n", "\n", "ndt = int( (tend-t0)/dt ) # number of timesteps that we will take\n", "t = np.linspace(t0,tend,ndt) # create an equally space array for time. This will be used for plotting.\n", "u= np.zeros(ndt) # allocate a numpy array of the same size as the number of timesteps\n", "\n", "\n", "c = 12.5 # drag coefficientkg/s\n", "m = 60 # object's mass, kg\n", "g = 9.81 # gravitational acceleration m/s/s\n", "\n", "n = 0 # just a counter\n", "while n < ndt-1:\n", " u[n+1] = u[n] + dt *(g - c/m*u[n]) # time advance\n", " n += 1\n", " \n", "plt.plot(t,u)" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "# Interpolation" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "fragment" } }, "source": [ "Use linear, polynomial, and cubic spline interpolants to interpolate the function\n", "$$ f(x) = e^{-x^2/\\sigma^2}$$\n", "on the interval $[-1,1]$. Start with $n=10$ samples and experiment with different values of the standard deviation, $\\sigma$." ] }, { "cell_type": "code", "execution_count": 151, "metadata": { "collapsed": false, "slideshow": { "slide_type": "subslide" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAecAAAFpCAYAAACmt+D8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzs3Xl8TFcfx/HPyWRfBBFbIvY9iSAoqvati61F0ZKEtlpL\nLbW01NZqtVS1pVS11hJKhdqrpbRaBLHEEoSQ2CKIyD6Z+/wx6qG1RCS5k+T3fr28ns7kzsw3npFv\nzp17zlGapiGEEEIIy2GldwAhhBBC3EvKWQghhLAwUs5CCCGEhZFyFkIIISyMlLMQQghhYaSchRBC\nCAsj5SyEEEJYGClnIYQQwsJIOQshhBAWRspZCCGEsDDWer1wsWLFtHLlyun18kIIIUSu27dv31VN\n09wfdZxu5VyuXDlCQ0P1enkhhBAi1ymlojJznJzWFkIIISyMlLMQQghhYaSchRBCCAuj22fOQggh\nckZ6ejrR0dGkpKToHaXAsre3x9PTExsbmyw9XspZCCHymejoaFxcXChXrhxKKb3jFDiaphEXF0d0\ndDTly5fP0nPIaW0hhMhnUlJScHNzk2LWiVIKNze3JzpzIeUshBD5kBSzvp7071/KWQghRLYzGAz4\n+fnd+TNlypRse+6wsDA2bNiQbc9niR75mbNS6nvgeeCKpmne9/m6Ar4AngWSgABN0/Znd1AhhBB5\nh4ODA2FhYTny3GFhYYSGhvLss8/myPNbgsyMnBcA7R7y9fZA5dt/XgdmP3ksIYQQ+U18fDxVq1bl\nxIkTAPTo0YNvv/0WgDfffBN/f39q1qzJ+PHj7zxm7969NGrUiFq1alG/fn3i4+MZN24cy5cvx8/P\nj+XLl+vyveS0R46cNU3boZQq95BDOgKLNE3TgL+VUoWVUqU0TbuYTRmFEEJk0ZAhQ7J9BOvn58eM\nGTMeekxycjJ+fn53br/77rt0796dmTNnEhAQwNtvv83169d57bXXAJg8eTJFixYlIyODli1bcujQ\nIapVq0b37t1Zvnw59erV4+bNmzg6OjJp0iRCQ0OZOXNmtn5fliQ7plJ5AOfvuh19+z4pZyEsREZG\nBpGRkURGRmIymbC1taVq1ap4eHjIhUMiRzzotHbr1q358ccfGTBgAAcPHrxz/4oVK5g7dy5Go5GL\nFy9y9OhRlFKUKlWKevXqAVCoUKFcy6+3XJ3nrJR6HfOpb7y8vHLzpYUocFJSUli9ejXBwcFs2bKX\nlJRyQDnMn2alAycpUuQKHTu2pkePHrRq1QorK7lGNL951Ag3t5lMJo4dO4ajoyPXr1/H09OTM2fO\nMG3aNPbu3UuRIkUICAgo8AuoZMe/xBigzF23PW/f9x+aps3VNM1f0zR/d/dH7pglhMiC9PR0vvnm\nG8qXb0jPngfZ8Gd/UqpOgZbNoIsVdNGgkwM07cT1UtNZ9GNH2rZdiZ9fY9atW4f5Eyohcsbnn39O\n9erVWbp0KYGBgaSnp3Pz5k2cnJxwdXXl8uXLbNy4EYCqVaty8eJF9u7dC0BCQgJGoxEXFxcSEhL0\n/DZyXHaMnNcCA5VSwUADIF4+bxZCH0eOHKFHjyEcOdEWaveBjgsxlvoEAIOypqSjJ9ZW1qQYU7iS\nvAgNDZPJAFHPcCT0DV7osIPnn1vEvHlfUaJECZ2/G5GX/fsz53bt2hEYGMi8efPYs2cPLi4uPPPM\nM3z44YdMnDiR2rVrU61aNcqUKUPjxo0BsLW1Zfny5QwaNIjk5GQcHBzYunUrzZs3Z8qUKfj5+d35\nLDu/UY/6LVkptQxoBhQDLgPjARsATdPm3J5KNRPzFd1JQKCmaY/cqNnf31+T/ZyFyD5ffz2bQW/v\nwVSvEjT5HBzj8Claj1frdOXZys9Sxa0KNob/r/ObnJ7MwcsH+fnEOhbsW8aF5Ei4VhF+GUqhC5v5\nadXbtGzZUsfvSGTVsWPHqF69ut4xCrz7/f+glNqnaZr/ox77yHLOKVLOQmSPjIwMBg4czZw1xaDj\nUih5iKdLtmVK+/dp7NU4c89hyiDkeAgj1k/kTNJhONkefvZhzqeVeOON13L4OxDZTcrZMjxJOcvV\nH0LkYenp6Tz3/CDmhNpB33EUKhnHyq4/seP1jZkuZgCDlYEXa7xIxPD9TGk2HeuKO+CN+fT/9BJj\nx36ag9+BEOJ+pJyFyKOMRiPPvfAOm22vw/OTqefeksgRB3mxRucsT4+ytrJmVNOhHBoYSqlCxeHV\n8UzeepmxY6dnc3ohxMNIOQuRB2maRscXx/BLkRNQJ5g3qo3l70HrcHN0y5bnr+5ejZMjd9OgcEdo\nP53Jfx1j8uS52fLcQohHk3IWIg96e/jnbHD9G6puYVzdmczp/gFWKnv/OTvZOvHn4JW0c+8Hz8xj\n7I7fWbEif282IISlkHIWIo+Z/c1KvorZDBV38GG9eUx8fkCOvZbBysCGN+fSrsib0GgpL3+9nNDQ\ng49+oBDiiUg5C5GHHDgQzlvrF0ONLQyoOI0xzwXl+GsqpVg/aCZ1rV9Ga76IZ4Z8SXx8fI6/rsjb\nnJ2d/3PfnDlzWLRokQ5p8p5cXb5TCJF1iYmJNB76KTRfSxvHt5n5yvBce20rZcXf7y7Ga9QlLrZY\nRINuhTi2abqsyy0eS//+/XP0+TVNQ9O0fLEMbd7/DoQoIBr2Gkdy0x8odbMVG9/J/aunra2sOTR+\nNTa3vDhRawkD3pPdYcXjmTBhAtOmTQOgWbNmjBo1ivr161OlShV27twJmOftjxgxgnr16uHr68s3\n33wDwK1bt2jZsiV16tTBx8eHNWvWAHD27FmqVq1K79698fb25vz58/d/8TxGRs5C5AFTZ6/kcOXF\n2MRX4MjkVdl+8VdmFXMuzPY31tN4QQNmX51L37AW1PWrpksWkUlDhkA2bxmJnx9kw4YaRqORPXv2\nsGHDBiZOnMjWrVv57rvvcHV1Ze/evaSmptK4cWPatGlDmTJlWL16NYUKFeLq1as89dRTdOjQAYCT\nJ0+ycOFCnnrqqSfOZClk5CyEhbtwMZZRf30GdgmseWU5RZ303TavUdVqDKk4HTwP0nT8WDIyMnTN\nI/KuLl26AFC3bl3Onj0LwJYtW1i0aBF+fn40aNCAuLg4Tp48iaZpvPfee/j6+tKqVStiYmK4fPky\nAGXLls1XxQwychbC4tUd+A6a7990sZlMe//aescB4PN+fVk1cCPna//E84Mms/HrcXpHEg9iYVtG\n3s3Ozg4Ag8GA0WgEzJ8bf/XVV7Rt2/aeYxcsWEBsbCz79u3DxsaGcuXK3dlW0snJKXeD5wIpZyEs\nSFRcIt/ujCTkwAUSU40knjxAbI3lOEc3Y+Xcd/WOd4+wjxfgPnY/m+y+JegLf/ZcM5CYasTJzppO\ntUvzWpMKlHXLfz80Rc5q27Yts2fPpkWLFtjY2BAREYGHhwfx8fEUL14cGxsbtm3bRlRUlN5Rc5Sc\n1hbCQmw7cYV2M3YSvOc8t1KNGNNSiC3yI6Q54ubcj+0RsXpHvEdRF2f6VpgErjEs2DaLW6lGNOBW\nqpHgPedpN2Mn205c0Tum0ElSUhKenp53/kyfnrmLGPv160eNGjWoU6cO3t7evPHGGxiNRnr16kVo\naCg+Pj4sWrSIatXy97UOsiuVEBYgKi6RdjN2kpz+/89vYw4FY2ywBPtDgylRuQ0ONgY2DWliMaPR\nfzKfPD6HDO81OO0bRbEaTe45xtIyFxSyK5VlkF2phMjjvt0ZSXqG6c7tpEvnMPqFQORTuFdqDUB6\nhol5O8/oFfE//slcwuMVuO5JYplVZKSm3HOMpWUWIq+Qz5yFsAAhBy5gNP3/LNbVlJVgSOfD/c74\nGv+/ZaP1WisILqlHxP9ofOQSDW7/QrHI2ZeVz2/gcugaStfsfucYo0lj9YEYPujkrVdMIfIkKWch\nLEBiqvHOf8ed3I3m+xtNfn+GwKgwkm3s7z3YeCmX091fldjEO/89Oyqe7VX9ueoTQnJsUxzc/v8L\nRGKa8X4PF0I8hJSzEBbAyc6aW6lGMjKMJJVYjuGGB6NPJdD0jW9JN9jcOc7ZzpojE9s+5JlyT8fx\nm7l1+5cK91vX+SLkE3oFpnI9fjkObm/fOc7JVn7MCPG45DNnISxAp9qlsbZSXI34BVOJCHptq8y7\nnUfeU8zWVorOtT10THmvfzIDxDoXYWHb/tTZ3ZD0Glu5GX0CsLzMQuQVUs5CWIDXmlQAYxopFX/C\nKaYGrqX9iHdwuecYG4MV/ZqU1ynhf73WpAI2hv//CDnhXo6uV1xQKa7cMPwIWF5mIfIKKWchLEBZ\nNycyrm4G14uM22rPlupP3/matZXCwcbA16/UsagpSWXdnPj6lTo42BjujKB31G7D8zt80Sr8TfzZ\nfRaXWeQeg8GAn58f3t7edO3alaSkpIcef78tJnNDaGgogwcPfugx27dv5/nnn8+lRGZSzkJYgIhz\nVzhbYinuEf40r1Qea2dnlDJ/xtyjvhebhjShedXiesf8j+ZVi7NpSBN61PfC2c6aYyUrMObMRQw3\nPLhZZCVNK7vrHVHoxMHBgbCwMI4cOYKtrS1z5szRO9J9+fv78+WXX+od4z+knIWwAM99OBbsbzLz\ntxvUmzicIxPbcubj5zgysS0fdPK26NFnWTcnPujkfSdzg/cG0fv3sphKHuaN6d/rHU9YgCZNmnDq\n1CkApk+fjre3N97e3sy4z7rfvXv3JiQk5M7tXr16sWbNGhYsWECXLl1o164dlStXZuTIkXeOWbZs\nGT4+Pnh7ezNq1Kg79zs7OzNixAhq1qxJq1at2LNnD82aNaNChQqsXbsWuHdUvGfPHho2bEjt2rVp\n1KgRJ06cyJG/j8yQyyiF0NmRyEucKhZM8WMN6Fo4HvL67jo9ezJzxEgWNfHi+yvTmZMRhMGg9E5V\nYOm9Y6TRaGTjxo20a9eOffv2MX/+fHbv3o2maTRo0ICmTZtSu/b/N3Tp27cvn3/+OZ06dSI+Pp5d\nu3axcOFClixZQlhYGAcOHMDOzo6qVasyaNAgDAYDo0aNYt++fRQpUoQ2bdoQEhJCp06dSExMpEWL\nFkydOpXOnTszduxYfvnlF44ePUqfPn3ubDn5j2rVqrFz506sra3ZunUr7733HqtWrcrOv7pMk5Gz\nEDrr8PFYsL3Ft9suowIDQeXxInNzw7FjBwL+KIOp+FH6fSaj54IoOTkZPz8//P398fLyom/fvvzx\nxx907twZJycnnJ2d6dKlCzt37rzncU2bNuXkyZPExsaybNkyXnzxRaytzePIli1b4urqir29PTVq\n1CAqKoq9e/fSrFkz3N3dsba2plevXuzYsQMAW1tb2rVrB4CPjw9NmzbFxsYGHx+fO1tU3i0+Pp6u\nXbvi7e3N0KFDCQ8Pz9m/pIeQkbMQOjp65gpn3JfjdvwpXojbA6++qnek7BEYyNfPtWd+o3IsNn3O\nPBk960avHSP/+cw5K3r37s2SJUsIDg5m/vz5d+7/Z4tJuHebyQexsbFB3f5l18rK6s7jrays7vvY\n999/n+bNm7N69WrOnj1Ls2bNspQ/O8jIWQgddZ06CexuMSv0FqpdOyhVSu9I2aN1a2xLedD5YEUy\nioczYvaPeicSFqBJkyaEhISQlJREYmIiq1evpkmTJv85LiAg4M7n0TVq1Hjoc9avX5/ff/+dq1ev\nkpGRwbJly2jatGmW8sXHx+PhYZ6Xv2DBgiw9R3aRchZCJxeuJnDUZSkOJxvQ/fRhCArSO1L2MRig\nTx8W/PU7JBRn1sEv0GkDPGFB6tSpQ0BAAPXr16dBgwb069fvns+b/1GiRAmqV69OYGDgI5+zVKlS\nTJkyhebNm1OrVi3q1q1Lx44ds5Rv5MiRvPvuu9SuXfuRo/KcJltGCqGTxu+8xy6Xj/loS3vejdgL\nMTFga6t3rOxz6hRUrkyD59qzp95GPq20lRG9WuqdqkDI61tGJiUl4ePjw/79+3F1ddU7TpbJlpFC\n5DEJSWnsYj5WZ2szOvRX6NUrfxUzQKVK0KQJPx4/CSkuTPhlmt6JRB6wdetWqlevzqBBg/J0MT8p\nKWchdNBn6hxwucSrV6qg0tIgE6fv8qSgILxOn6LsqadJKruF5VsO6Z1IWLhWrVoRFRXFkCFD9I6i\nKylnIXKZ0aix5vI8iC3HN+fCoU4dqFVL71g546WXwMmJpam2YLJm4NIpeicSIk+QchYil42evR5T\nicPUj/XH7siR/DtqBnB2hu7dabTtV5xPNuSq52p2HbyodyohLJ6UsxC5SNNgZugsSCrMKgcH8+fM\nPXvqHStnBQbCrVtMd64C1qm8OnOy3omEsHhSzkLkooVrw0ktv5lS0U/huW4DdOoERYvqHStnNW4M\nlSvTL+I4hlN1iSwazIUrKXqnEsKiSTkLkYtGrZoGJgPzK9eFuLj8fUr7H0pBYCBq5056ONcFxzgC\nP/tG71Qih126dImXX36ZihUrUrduXZ599lkiIiIe+phy5cpx9erV/9w/Z84cFi1alOnXnjx5MjVr\n1sTX1xc/Pz9279790OMnTJjAtGnm2QTjxo1j69atmX6tnCLLdwqRS8JPxXPFcxV2p+vR5sYB8PCA\n1q31jpU7eveGsWP5xsmRJbFl2WpcSEbG2xgMegcTOUHTNDp37kyfPn0IDg4G4ODBg1y+fJkqVao8\n9vP1798/08f+9ddfrFu3jv3792NnZ8fVq1dJS0vL9OMnTZr02PlygoychcglAV9MB7sEBpZvjtq0\nyVxYBaWdPDygbVscl6+g+s0GmEodYNJ3v+mdSuSQbdu2YWNjc0+p1qpViyZNmtyzRSPAwIED71kq\n89NPP8XHx4f69evf2Wby7pHtqVOnaNWqFbVq1aJOnTqcPn36nte+ePEixYoVu7OOdrFixShdujRg\nHpmPHDnyP89/t4CAAFauXHnn+PHjx1OnTh18fHw4fvw4AImJiQQFBVG/fn1q167NmjVrnvSv7D9k\n5CxELkhMyiDUsBB1vgYflbUHk6lgnNK+W1AQdO3K0gaNqR2znuk7ZzHx9RZ6p8r3hmwaQtil7N0z\n0q+kHzPaPXhHjSNHjlC3bt0sPberqyuHDx9m0aJFDBkyhHXr1t3z9V69ejF69Gg6d+5MSkoKJpPp\nnq+3adOGSZMmUaVKFVq1akX37t3vWWv7Uc//b8WKFWP//v18/fXXTJs2jXnz5jF58mRatGjB999/\nz40bN6hfvz6tWrXCySn79l2XkbMQuWDYrJ+gSBT1taew/eEHePppqFxZ71i564UXoGhR/Hb8iWtU\nfW6VW89vu2ValbhXjx497vzvX3/9dc/XEhISiImJoXPnzgDY29vj6Oh4zzHOzs7s27ePuXPn4u7u\nTvfu3e8ZmT/s+e+nS5cuANStW/fONpNbtmxhypQp+Pn50axZM1JSUjh37lyWvt8HkZGzEDlM02Dh\n0XlQrCg/tG4J338PI0fqHSv32dnBK6/AnDmMmfoRI69vo/83M4ho8IneyfK1h41wc0rNmjXvnBr+\nN2tr63tGuykp9165r+7az1xlcW9zg8FAs2bNaNasGT4+PixcuJCAgIAsPf8/p8fv3qJS0zRWrVpF\n1apVs5QvM2TkLEQOW/XLWVK9tuJ2vg4Vf90Gjo7QtavesfQRGAhpaQzTrLGKqsFJ12VcjcvQO5XI\nZi1atCA1NZW5c+feue/QoUPs3LmTsmXLcvToUVJTU7lx4wa//vrrPY9dvnz5nf9t2LDhPV9zcXHB\n09OTkJAQAFJTU0lKSrrnmBMnTnDy5Mk7t8PCwihbtmymnj+z2rZty1dffcU/G0cdOHAgS8/zMFLO\nQuSwYYu/BCsTk9p0huXLoVs3cHHRO5Y+/Pygdm0MixbRwrkRFD7Pm58H651KZDOlFKtXr2br1q1U\nrFiRmjVr8u6771KyZEnKlClDt27d8Pb2plu3bv/ZMvL69ev4+vryxRdf8Pnnn//nuRcvXsyXX36J\nr68vjRo14tKlS/d8/datW/Tp04caNWrg6+vL0aNHmTBhQqafPzPef/990tPT8fX1pWbNmrz//vtZ\nep6HkS0jhchB0ReMlJlWHsP1QqQ1G4lVQAD8/js884ze0fTz1VcweDDnN6zD65dArG/UJHXeNqxk\nqJBt8vqWkTmlXLlyhIaGUqxYsVx5PdkyUggLNfCLZeAazXMlmmC1YAFUrAhNmugdS189e4KtLWU2\n/0LJC34Yy25n0dqTj36cEAWIlLMQOcRkgvWXFkFCUeZ06g3bt0NAgHnFrILMzQ06doQlS5j4XDfQ\nFGN/mqV3KlEAnD17NtdGzU9KylmIHLL451MYy/9KqUt1KbVxs7mU+/TRO5ZlCAqCuDj62hXCcLYm\nMcVWcemyXBgmxD8yVc5KqXZKqRNKqVNKqdH3+bqrUupnpdRBpVS4UqqAra4gxH+N/ekrACY83xMW\nLDAv1VmmjL6hLEXr1uDhgWHhQloVbQCu0Qz+8v5Tb0TW6HU9kTB70r//R5azUsoAzALaAzWAHkqp\nGv86bABwVNO0WkAz4DOllO0TJRMiD7t42Ui0+48YIr3pW6QUnDtX8FYEexiDwXwWYdMmvunRDxIL\nE3JuMdIn2cPe3p64uDgpaJ1omkZcXBz29vZZfo7MLEJSHzilaVokgFIqGOgIHL07C+CizDO6nYFr\ngDHLqYTI44bOWgkuF2mZ/hyGRYugcGHz9pDi/wIC4KOPKPvrdopfrM2V8ptZtSmGl9p76J0sz/P0\n9CQ6OprY2Fi9oxRY9vb2eHp6ZvnxmSlnD+D8XbejgQb/OmYmsBa4ALgA3TVNMyFEAaRpsPrsEihZ\nmFkBb0CjJubPWJ/gt+h8qXJl85Xr8+czeuwQhkVuY9TS2bzU/kO9k+V5NjY2lC9fXu8Y4glk1wVh\nbYEwoDTgB8xUShX690FKqdeVUqFKqVD5jU7kV2u3XiGt3BaKRvtQac8+SEmRU9oPEhgIEREM8qqB\nVXQlIgv/SFycnIoVIjPlHAPcfRWL5+377hYI/KSZnQLOANX+/USaps3VNM1f0zR/d3f3rGYWwqKN\nXPwNGNIZ0qyDeR1tb2/I4g49+V7XruDkhPXixTSwrQvFIhg9+xe9Uwmhu8yU816gslKq/O2LvF7G\nfAr7bueAlgBKqRJAVSAyO4MKkRdcvw4RTqtQF8sysl5L2LPHPDos6HObH8TZ2byc6fLlfPHKm5Dq\nyA/h8+XCMFHgPbKcNU0zAgOBzcAxYIWmaeFKqf5KqX920v4AaKSUOgz8CozSNO1qToUWwlJNmrsL\nSh6ktqqN3dKlYG1t3olJPFhQENy6Rb2IszifrUlyhbX8/vdNvVMJoStZW1uIbOTa/TVuVlnIb23X\n0/ylV6FhQ1i9Wu9Ylk3ToGpVKFWKgc82ZFbKJ9S/8CW7vxmkdzIhsp2srS1ELjt8NI2bZUOwi6pO\n82vJcPmyXAiWGUqZ/5527ODDli/BtRLsTVtFaqrewYTQj5SzENlk+DfLwOkqL1YwTw+ieHFo317v\nWHlD795gZUXhkBDK3vBBK7uDeT+e1juVELqRchYiG2RkwLZryyGhKFO79IN168yFY2Ojd7S8wcMD\n2raFBQt4/4UeoDQ+XjdP71RC6EbKWYhssHLjRYwVtlD8ojelf90GRqOc0n5cgYEQE0OAY0mszlci\nplgIV67IZduiYJJyFiIbjFsxD6wyGNL89tzm+vWhxr+XoBcP1aEDFC2KYeFCGjnVBffjjJuzXe9U\nQuhCylmIJ3TzJkTYr0Vd8mJotQZw5IiMmrPCzg569YKQEKZ3CQKjLUuOLNE7lRC6kHIW4gl9tuAg\neIRSw+iN/bJl5jW0X35Z71h5U1AQpKVR73AEDlHVSCy/hn1haXqnEiLXSTkL8YRm7VwIJismd+4D\nS5dCly7mXajE4/PzM/+ZP58XKzUBxzhGzl2hdyohcp2UsxBPIDJSI84jBOvoSnSIz4AbN+SU9pMK\nCoL9+/m04QtwqzC/x/+IUTagFQWMlLMQT2DsN1ugyBlaFvNHLVgAXl7QooXesfK2nj3B1pZS6zdS\n7FI1MipuZPWmOL1TCZGrpJyFyCJNg5DI5ZBux7R2veGXXyAgAKzkn9UTcXODjh1hyRLeatAWDOlM\nXLFQ71RC5Cr5KSJEFu0LSyO5whqczlfDe1eoua0DAvSOlT8EBkJcHKNcK0JcSY4a1pKcrHcoIXKP\nlLMQWfTudz+C4zU6V2hkXq6zWTMoX17vWPlDmzbg4YFjcDDlEqqjld3B/JVReqcSItdIOQuRBSYT\nbL+2EhJd+aRmUzh92nwhk8geBoN5+dNNmxjTqB0ojU/WL9A7lRC5RspZiCzYuO0GxoobKRJTjdLr\nN4GLC7z4ot6x8pfAQDCZCLySjrpUhnOFfub6db1DCZE7pJyFyIIJwYvBOpW+3k/Djz9C9+7g6Kh3\nrPylcmV4+mkMixbhiw947GP6wkN6pxIiV0g5C/GY0tJgf+rPcN2d9x08ITFR5jbnlKAgiIhgms8z\noCnm/LFY70RC5AopZyEe04p1VzGV20aJq5UptGoVVK0KDRvqHSt/6toVnJxouS8C65gKXC31M9HR\nslOVyP+knIV4TB+tXgwGI0Or14M//jCPmpXSO1b+5OwM3bqhVqygqXNtKHaCSd/u1DuVEDlOylmI\nx3DrFhw3rIWrJRh83WBecOTVV/WOlb8FBsKtW3xZojpkGFh2OFjvRELkOClnIR7D/B8vopXdgdfN\nijgEL4d27aB0ab1j5W9PPw2VKlFj4+/Yx1TlVrmfOXrMpHcqIXKUlLMQj+GzjYvBysQEj5oQEyNz\nm3ODUubR844ddHXzBddo3p39s96phMhRUs5CZFJsLEQ5r4UrJXj1RJx5DegXXtA7VsHQuzdYWfFp\nWhFIs2dzTAiaXBcm8jEpZyEyafbSKCj7J9WTKmK9bh306gW2tnrHKhg8PaFtW0qu+plCF6uRWmkt\nu3an6p1KiBwj5SxEJs3atgSAT+zLmCc7y9zm3BUYCNHR9C9UERyv8d53K/VOJESOkXIWIhPOnoUr\n7j+jLpXiud0RULs2+PnpHatg6dABihZl7HkjJDvzZ/xajEa9QwmRM6SchciErxZHgOduGqaUwerA\nARk168GeXI4mAAAgAElEQVTODnr1wmXdRkpcqkpG5fVs/DVB71RC5AgpZyEyYcHepQB8mlzE/Dlz\nz546JyqggoIgLY3R1iXANpFJy5bpnUiIHCHlLMQjnD4N10qvxXDBg0a/7jOfXnVz0ztWweTnB35+\nvBUaDQmF2Z+2Tk5ti3xJylmIR5ix5AiUOkC7xFKoq1dlbrPegoKwPXiI8perYKq4hTWbZB9Jkf9I\nOQvxCIv3LwVNMeWCZl4NrE0bvSMVbD17gq0tk5LtwDqVD1bITlUi/5FyFuIhTp6EeM+12MaUoeaO\nMPNiGAaD3rEKNjc36NiRntvD4UYxDmvrSU/XO5QQ2UvKWYiHmLpwLxQPp1t8UVRGhlylbSkCA7GK\nu4bvpYqYyv/GinWX9E4kRLaSchbiIYKPLAeTFR8cjoPGjaFKFb0jCTB/tODhwceXksFgZPJPC/VO\nJES2knIW4gGOH9dI8FqLc3RZyh0/L6NmS2IwQO/etNtxGKu4Uhy33khamt6hhMg+Us5CPMCUBX+C\n20kCYh3A0RG6ddM7krhbYCBWJo2GMWXQvHayeHWU3omEyDZSzkI8wKoTKyDDwPu7zsJLL4GLi96R\nxN0qV4ann2bqiRiwMjHlZzm1LfIPKWch7uPoUY1bZX/G7Vx5isclydxmSxUURMPwGGyulOGU/SZS\nUvQOJET2kHIW4j4+mP8rFDnLW+c1qFABnnlG70jifrp2BScnWp0rAWX+4ruVJ/ROJES2kHIW4j7W\nnl4JRhuG/3kaAgJAKb0jiftxdoZu3ZgaegyAqRvk1LbIH6SchfiXw0dMJJVfh0dUOVzTFPTpo3ck\n8TCBgdS8lIjDhfJEuWwiOVnvQEI8OSlnIf5l/HfroVAMw0/cglatwMtL70jiYZ5+GipVotNpVyh9\ngJnLwvROJMQTk3IW4i6aBhvPhUC6Ha+FXZS5zXmBUhAYyOS9hwD4YusinQMJ8eSknIW4S9ihdFIq\nrKNiZFmcHVyhUye9I4nM6N2bcreg0PnKxBTeQlKS3oGEeDJSzkLc5f15IeB8hZGHr0KPHuDgoHck\nkRmenqg2behxzAZKhDN14S69EwnxRKSchbhN02DrpbWoVAdePX5N5jbnNUFBjDt4DExWzN75g95p\nhHgi1noHEMJS7AtLI7XCBnxOeeFQ1Rr8/fWOJB5Hhw6UsitCsbMluFxsMwkJGi4uMgVO5E2ZGjkr\npdoppU4opU4ppUY/4JhmSqkwpVS4Uur37I0pRM5799vl4HiNEYcumC8Ek7nNeYudHapXLwKOmMDt\nNJO//03vREJk2SPLWSllAGYB7YEaQA+lVI1/HVMY+BrooGlaTaBrDmQVIsdoGvwetx6rFCe6nkmE\nV17RO5LIisBARh6LgAwbvtuzTO80QmRZZkbO9YFTmqZFapqWBgQDHf91TE/gJ03TzgFomnYle2MK\nkbP+3ptEeoVN+J/wxL79C1CihN6RRFbUro171Vp4nK7C1eJbiL9p0juREFmSmXL2AM7fdTv69n13\nqwIUUUptV0rtU0r1zq6AQuSG0fOWgn08ww5fkLnNeV1gIEFH0qDwed7/ZqPeaYTIkuy6WtsaqAs8\nB7QF3ldKVfn3QUqp15VSoUqp0NjY2Gx6aSGejKbBroQNGJIK0TnBDp59Vu9I4kn06sXQ01GQbscP\nB5frnUaILMlMOccAZe667Xn7vrtFA5s1TUvUNO0qsAOo9e8n0jRtrqZp/pqm+bu7u2c1sxDZasdf\n8RgrbKHR0VLYvtIHbGz0jiSehJsbRdp3oMLJalwrtYXrNzL0TiTEY8tMOe8FKiulyiulbIGXgbX/\nOmYN8LRSylop5Qg0AI5lb1Qhcsbo7xaDbSJDw+WUdr4RFETQkSRwvszoOT/pnUaIx/bIctY0zQgM\nBDZjLtwVmqaFK6X6K6X63z7mGLAJOATsAeZpmnYk52ILkT00DfambMbmVhGed6sCNWvqHUlkhzZt\nGHjjJirVieXhP+qdRojHlqlFSDRN2wBs+Nd9c/51eyowNfuiCZHzfvk9lozyW2m2vxw2QX31jiOy\ni8GAa58gqpzYzIlKvxIbl4a7m63eqYTINFm+UxRo7y1cCDYpDDl+ybyWtsg/AgLod+QWOF7jndnB\neqcR4rFIOYsCS9MgLH0rdjfdaNegLRQurHckkZ2qVOGNIq6o5EKsjpDPnUXeIuUsCqx1W6PJqPAb\nzcKLYy2ntPMll9f7U/NYVRI8f+NibLLecYTINClnUWCNWbwQDOm8fT4BWrTQO47ICV270vd4Etgl\nMGTWQr3TCJFpUs6iQNI0COdXHK6XoG2XPmAw6B1J5AQXF/rV8cYq0Y31Z0P0TiNEpkk5iwLpx/WR\nmMrvoEV4MawCZd/m/Mz5rQH4hlch0XMHURdv6h1HiEyRchYF0vjgBWCVwaB0O6hQQe84Iic9/TQB\n5zLAJpm3v/5O7zRCZIqUsyhwTCaIsPkNp6ulad1noN5xRE5Tin6d2mC4WZKt0f9e3FAIyyTlLAqc\nhT8dwVR2F63Ci2HVrZvecUQucHrjDfzCK5HouYtT0XF6xxHikaScRYEz+af5oDTeLlYWnJz0jiNy\ng6cnvW86gnUag2bOefTxQuhMylkUKBkZcMZ5G64Xy9P87VF6xxG56LU3emF9vQw7r6zTO4oQjyTl\nLAqUmT/8hcnjAK0j3KFRI73jiFzk8HJ36oRXIrHMXo6cuaB3HCEeSspZFCgz1i0AYHiteqCUvmFE\n7rKzo5d1SbDKYMhXX+qdRoiHknIWBYbRCOeK7cTtXFWeGvGe3nGEDl4fNxib2Ir8ffMXvaMI8VBS\nzqLA+GjuRkwljtHqfEkoXVrvOEIH9g2fovbxSiR6HiD0RKTecYR4IClnUWB899sSMBl474VOekcR\nOnqlRDVQGiM+n6Z3FCEeSMpZFAipqSbOe+yiRGRNfF9/U+84QkdvTB2F7cXq7DHt0DuKEA8k5SwK\nhDHTFqMVPUurG+XAzk7vOEJHtqVL4Xu6Ckke4ew8cFDvOELcl7XeAYTIKVFxiXy7M5KQAxe4uH8l\n+NlSuW57ouISKesmi48UZF0r1SOUNfQePwFq9MPJzppOtUvzWpMK8t4QFkFGziJf2nbiCu1m7CR4\nz3niE5K4XjWUEqd8WZzmRbsZO9l24oreEYVOtp24wgJXH+zO1+Jy6XA04FaqkeA95+W9ISyGlLPI\nd6LiEnlryX6S0zMwmjQcDm9Hc7mE1yUvjCaN5PQM3lqyn6i4RL2jilz2z3sjCQNeJyuTXOokdhfN\np7blvSEsiZSzyHe+3RlJeobpzu2bTuGQ6kx8tefv3JeeYWLezjN6xBM6uvu9YePgC5qC2O33HCPv\nDWEJpJxFvhNy4AJGkwaAQ/w1Llfbh/vJWqS6Fr9zjNGksfpAjF4RhU7ufm8k+NbGLsqfcxVOoGn/\n/2VO3hvCEkg5i3wnMdV457+LHNkIDjdwTqj53+PSjP+5T+Rvd783rKwUJSKrkuJ+Ds9T964YJu8N\noTcpZ5HvONmZJyEozURUiXOopCKk12jz3+NsZbJCQfPPe+MfGe6NwWQgMWn/vcfJe0PoTMpZ5Dud\napfG2krR8Oh2LlTbT6GIehhs7e85xtpK0bm2h04JhV7+eW/8w7pSGWwjG3K8aiSl4y+b75P3hrAA\nUs4i33mtSQVsDFakXwkFm2Q07b9bQ9oYrOjXpLwO6YSe/nlv3M3uog+pRS5S+8AyQN4bwjJIOYt8\np6ybE0tqZrCjWhwqrjyFKte58zVrK4WDjYGvX6kji00UQGXdnPj6lTo42BjujKAdPBuD0ZaDhW9Q\n0pgk7w1hEaScRb5k+8NXXCl/FOdzTXB1tEUpcLazpkd9LzYNaULzqsUf/SQiX2petTibhjShR30v\nnO2scSxVGMPpxuypeY51tofkvSEsglz1IPKfiAg+vnEBlMZ77QMY/UZzvRMJC1PWzYkPOnnzQSdv\nABr3D2NXoW2sW7GKgHFjwN7+Ec8gRM6SkbPId7Tpn7HO9zrqXH3e6dtM7zgiD5g1uB+kOTKzgiMs\nXqx3HCGknEU+c/ky2zasIKF4NOUTWmJtrR79GFHg+dUogd2ZFhyoeZKkzz4Fk+nRDxIiB0k5i/xl\n5kw+rlkajLZ82KOP3mlEHtKu5AuYHG7ypbURfv5Z7ziigFOapunywv7+/lpoaKgury3yqcRE0sqV\nwSlQwxTdCOMP61EycBaZFHMhGc9pVfCKdSHqTFH44w+9I4l8SCm1T9M0/0cdJyNnkX98/z3BRcDo\ndAN/69ZSzOKxeJR2oNDZdpyrcIKYg3/Crl16RxIFmJSzyB+MRpg+nSl1S0OSG18O7q13IpEHverb\nBaxMTKpTBaZO1TuOKMCknEX+sGoV8RfPcqzSSexOPksD/6J6JxJ50MfD28L5+vzgk4S2JgQiIvSO\nJAooKWeR92kaTJ3KZ09VBus0nvV4/tGPEeI+XFysKHWlNYnu0ezxtIfPPtM7kiigpJxF3rd9O+zb\nx8yqwNUqzBjRRe9EIg97p31PMNoxqnFVWLgQLl/WO5IogKScRd736aecLl+U66VP4hr1HF5esvCd\nyLoBQTVQJ55jp1ckqcZUmDlT70iiAJJyFnnb4cOwaRPvNK4GQEDdl3QOJPI6OzuontYWk2MCS17w\nha+/hsREvWOJAkbKWeRt06ahOTmyrlgknH2aiUP+uz2kEI/rw8BukFCSiSUz4No1+P57vSOJAkbK\nWeRd0dGwdCnrXn4GY+FLeMY+j6ur3qFEftDh+cIYjr7IeffjXG5SB6ZPN0/XEyKXSDmLvOuLL0DT\nGO5wHVKdebdjL70TiXzCYIBmRZ4HQwbjG7jD2bOwcqXesUQBIuUs8qb4ePjmG25268jJQmGooy8S\n9Kqn3qlEPvLJsNZwsTYL049ClduLkui03LEoeKScRd40dy4kJPBRHRewTaWW6TnZgldkqzp1DDid\n7EJKkfP8FfQ87N8P27bpHUsUEJkqZ6VUO6XUCaXUKaXU6IccV08pZVRKySWzIuekpcGMGdCyJbMv\n7ITYanz4hiw8IrKXUtC7TmfIsGbozUNQvLgs6SlyzSPLWSllAGYB7YEaQA+lVI0HHPcJsCW7Qwpx\nj6VL4cIFDr7WmZuukVgfeZl27Rz0TiXyoTFDakLE8+wx/k3qgP6waZN5+p4QOSwzI+f6wClN0yI1\nTUsDgoGO9zluELAKuJKN+YS4l6bBtGng68vwM5vAZOA5zxcwGPQOJvIjDw/wvNIJzfEWM6tZgaOj\n+f0nRA7LTDl7AOfvuh19+747lFIeQGdgdvZFE+I+Nm6E8HDShw9h+43fIeJ5Jo2srXcqkY+NerET\nxHsydU8w9O1rPnMTHa13LJHPZdcFYTOAUZqmmR52kFLqdaVUqFIqNDY2NpteWhQoU6eCpyfflUoi\nwyGBolGd8fWVjZtFzgno44o62JvLzic4GdAJTCbzND4hclBmyjkGKHPXbc/b993NHwhWSp0FXgK+\nVkp1+vcTaZo2V9M0f03T/N3d3bMYWRRYoaHmTS6GDOHjrXPgVgkGtXtO71Qin3N2hvq25gsO3/vr\nW+jWDb75xjydT4gckply3gtUVkqVV0rZAi8Da+8+QNO08pqmldM0rRywEnhL07SQbE8rCrapU6FQ\nIc51bcc5+6Nw8BUGvFlM71SiABg3qB6cbs3a8xvJeGcYJCSYC1qIHPLIctY0zQgMBDYDx4AVmqaF\nK6X6K6X653RAIQCIjDSv0NS/P+9t/ASsTPimd0JOwIjc0LatNXbhvUhziCfYGAEtWphPbael6R1N\n5FOZ+sxZ07QNmqZV0TStoqZpk2/fN0fTtDn3OTZA0zRZ505kr88/B4OBjEEDWXnmZ4hswfsDntI7\nlSggDAboWbcNJBbjg/WfwogRcOGC+eIwIXKArBAmLN/Vq/Ddd/DKK6y4+AepDjewOdyHDh1k32aR\ne0YOLwkH+3CCI1xs6As+PuZpVbKkp8gBUs7C8n39NSQnwzvv8MGGKZDozks1W2Nrq3cwUZBUqwZl\n47qAwcS41RPgnXcgPNw8vU+IbCblLCxbcjLMnAnPPUeMpyvHTIfhQCDvDC2ldzJRAI3uVw+inuaH\n4ysxde9mXqVElvQUOUDKWVi2hQshNhZGjGDCmglgpVEipiu1Zd0RoYNevWywOtCXZIfrrAr/GYYM\nMU/vCw3VO5rIZ6ScheXKyIDPPoN69ch4ujE/HFsBkS0Z0tsHJeuOCB24uED7su0g0Y0J6z+E11+H\nQoVk9CyynZSzsFxr1sCpUzBiBCsPriLZ7iZqX1/69bPTO5kowEYOKwkH+nLUeJjzWjz072+e5hcZ\nqXc0kY9IOQvLpGnm0UiFCtClC5PWfwi33Gnt1ZJisu6I0FGTJuB2ticoGL92PLz9tnmu1fTpekcT\n+YiUs7BMf/4Jf/8Nw4dzLiGGo8ZwCAtixLDieicTBZxSMKBnNTjVlmUnVpBewh169YLvvzdP+xMi\nG0g5C8v06adQrBgEBPDeT++BgqKRr9Cihd7BhIDXXrOD0DdIsUnkh30/mKdVJSebp/0JkQ2knIXl\nOXYMfv4ZBgwg1dbAj5E/wYkXGPRqVazkHSssgKcnNCrWFOJLM3nzR1CzJjz7rHnaX3Ky3vFEPiA/\n6oTl+ewzsLeHAQP47u/vSLNJhr1v8vrrNnonE+KOkSOKwL43OWU6yfHY4+YlPWNjzdP/hHhCUs7C\nsly8CIsXQ2AguLsz5bdP4GoFmnk1pHRpvcMJ8X/PPQeup3tChoExq8dA06bg72/+5TIjQ+94Io+T\nchaW5auvID0dhg1jb8xezmvnYO9ghg111TuZEPewtoaBAV5wrAtrz68jMT0JRo40T/9bs0bveCKP\nk3IWliMhAWbPhi5doFIlxoSMgTR7XM92pX17vcMJ8V/9+1vDnkEYrdP4YvsX5vduhQrmCxplQwzx\nBKScheX47ju4cQNGjCAuKY5fL/8GB/vwZmAJrGUDKmGBPD2hZZU6cMGX6X9+jmZlBcOGwe7d8Mcf\nescTeZiUs7AM6enmPZufeQYaNOCTLZ9gMmTA3rd4802D3umEeKDhw5xg93DirK6yMWKj+XoJNzdZ\n0lM8ESlnYRlWrIBz52DECDJMGczZ9w2cbUT7ulXw8tI7nBAP1qYNFL/yAtwqzNifx4KjIwwYYJ4O\neOyY3vFEHiXlLPT3z1Kd1avDs88SHBZMgvVN2D2MYcPs9U4nxEMZDPD2QFcIHcyBWwc4de0UDBxo\nng742Wd6xxN5lJSz0N/WrXDwoHmVJSsrxm8cD9dL4ZXclpYt9Q4nxKMFBVmh9r8GJivGrRsH7u4Q\nEGCeFnjxot7xRB4k5Sz0N3UqlCoFvXqxJ3oPp42nYfcI3hnmLFtDijyhZEno1LIEhHdh5amVJKQm\nmC8MS083Tw8U4jFJOQt9hYXBL7/A4MFgZ8eokFGQao9teG9699Y7nBCZN2yYDeweQbohnS9+/wIq\nVzZPrZo92zxNUIjHIOUs9DVtGjg7Q//+RN+M5vfY32H/a/R5uTCusu6IyEMaNwbvIrUg2pfpf07H\npJnMS3reuGGeJijEY5ByFvo5dw6Cg+H116FwYSZtmoQGsHsYgwfL9CmRtygFI0fawV9juG51nRUH\nV0CDBuYNoKdPN5/iFiKTpJyFfmbMMP9EGzKEW2m3WBi+EI61o2H10nh76x1OiMfXrRsUufQ8XC/B\nmA1jzHeOGAHnz5unCwqRSVLOQh/Xr8PcufDyy1CmDF/t/Io0qzT4eyzDhtnqnU6ILLGzg8EDHeCv\nMUSmR/JH1B/mHTKqVTNf+ChLeopMknIW+pgzBxIT4Z13MGkmpu2cBtHV8dDq0rmz3uGEyLo331QY\nDveBZEdGhowEKyvzNMGDB83TBoXIBClnkftSU+HLL81LK9WqRXBYMNfUNfhrPKNG2mGQj5tFHlai\nBPR40Qn2DuCv638RERcBr7xinm8lS3qKTJJyFrlvyRK4dAlGjEDTNN7b8B5cK4Hz+Q4EBuodTogn\nN3SoAXYPhwwD76599/b57sHmaYNhYXrHE3mAlLPIXSaTefqUnx+0bMnW01uJMkbBnxMY+JY9zs56\nBxTiydWpAw193eBgN0LOhnAl8Qr072+eNjhtmt7xRB4g5Sxy1/r1cPy4+QpWpXgn5B1IKIThSG8G\nD5blwET+MWqUNfw1DpOViUmbJkGRIvDaa+bpg1FRescTFk7KWeSuqVOhbFno2pXQmFAOJR6Cv0fw\nag97SpXSO5wQ2eeFF6BK0UpwvBnzDs7jVtotGDLE/MUZM/QNJyyelLPIPX//DTt3wtChYGPD8JDh\nkGIHoYMZPlzeiiJ/sbKCMWOs4Y+PSDWk8umvn4KXl3n64LffmqcTCvEA8hNR5J6pU82n9vr2JSIu\ngh2xO2Bvf9o0dZJFR0S+1KMHlDbVg8i6TP97OsnpyeaPdBITzdMJhXgAKWeRO06ehNWr4c03wdmZ\nd356B4wG+PtdRo2SuVMif7KxgdGjrWHHNBKtEvnqj6+gVi1o3do8nTA1Ve+IwkJJOYvcMX26+SfV\noEHE3IxhXfR6COvFUz7FaN5c73BC5Jy+faHIzcZwrjof7fiItIw08+j50iXztEIh7kPKWeS8K1dg\nwQLo3RtKlmToT0PRNA3+nMiECQbZs1nka46OMHyYDeyYTjzxfLv7W2jVyjydcNo08/RCIf5Fylnk\nvFmzICUFhg8n5mYMq86sgrCXqVuxDG3a6B1OiJz31lvgcKEVXCjPuF/GYdQyzEt6Hj9unl4oxL9I\nOYuclZRkLucOHaBaNYaHDMekabDzIyZOlFGzKBiKFIGBA8yfPV/jGvP3zjdvYeXlJUt6ivuSchY5\na/58iIuDkSOJvhnNj6d/hLBu1CrrybPP6h1OiNwzYgTYR70Al8swZtMYMgxW5mmFO3eapxkKcRcp\nZ5FzjEbzhWANG0LjxoxYM+L2qHkKEyday6hZFCju7v989jyNWGL5bs930K8fFC4so2fxH1LOIuf8\n9BNERsKIEUTfjGbFqRUQ1g1vT086dNA7nBC5b/hwcIrqDJfK8e6md0l3sDNPL1y92jzdUIjbpJxF\nztA082igcmXo0IG3V75tHjXv+IQJE2TULAqmIkVg9Cgb2PYF19Q1vt71NQwaZJ5mOH263vGEBZFy\nFjnj998hNBSGD+dcQgyro1ZDWE/qVPSkc2e9wwmhn7ffBtfL7SG6KuN+HUeqe1F49VXzdMMrV/SO\nJyyElLPIGVOnmj9k692b15e+jmaygh2TmTbNgJW860QB5uICY8fYwG8zualuMm37NPP57pQU88wG\nIZByFjkhPBw2bIBBgzhy8zSbr2yG0Ddo+5SnrAYmBOZ5z8USmsIZPybvmExiRS/zNlazZpmnH4oC\nT8pZZL9p08zLIr31FoE/BEKqPeyYwKefyhraQoD5n8eE8ebRc7IhmfHrx8PIkeZph/Pn6x1PWAAp\nZ5G9YmLghx+gb192JIQTmhAKu0bT+yU3fH31DieE5Xj9dahs3wBONuXL/V9xo25NeOop84VhRqPe\n8YTOMlXOSql2SqkTSqlTSqnR9/l6L6XUIaXUYaXULqVUreyPKvKEL7+EjAy0IUMIWhoECa5Yhw7l\ngw/k90Ah7mZjA198YQ2/ziDdkMZbwQPMK5VERpqnIYoC7ZE/MZVSBmAW0B6oAfRQStX412FngKaa\npvkAHwBzszuoyANu3jTvUdu1K8tv7OF0+mn4/SOGvOWCl5fe4YSwPO3bQxtfHwjrSfCZYCKf8YFK\nlcwXVGqa3vGEjjIznKkPnNI0LVLTtDQgGOh49wGapu3SNO367Zt/A57ZG1PkCXPnws2bZLwzjIEh\nAyHOE9fIAN57TyY1C/EgM2YYUNunoGXY8MqiPuYrt0NDzdMRRYGVmXL2AM7fdTv69n0P0hfY+CSh\nRB6UlgYzZkDz5ky58gtxhjj4dTqfT3OkSBG9wwlhuapXh4G9PWDXSP5K+IsdTSuapyHKkp4FWrZ+\nEKiUao65nEc94OuvK6VClVKhsbGx2fnSQm/BwRATQ/zQt5j05wdwrj7+Th3o00fvYEJYvgkTrHA+\nNBQSivLKin5oAwaYpyOGh+sdTegkM+UcA5S567bn7fvuoZTyBeYBHTVNi7vfE2maNlfTNH9N0/zd\n3d2zkldYIk0zT5/y9ibowlLSbNJg4yzmfWsnC44IkQlFi8LHE4vCb1M5zznm+DuAg4P535UokDLz\no3MvUFkpVV4pZQu8DKy9+wCllBfwE/CqpmkR2R9TWLTNm+HwYY4MfJmfLoTA/iAGvlibWnLNvhCZ\n9uab4Gt6GS7VYPgfH5DSt495WmLMf8ZCogB4ZDlrmmYEBgKbgWPACk3TwpVS/ZVS/W8fNg5wA75W\nSoUppUJzLLGwPFOnonmU5qXLiyDNCdd94/jwQ1lwRIjHYTDA4kWOqF8+J9n+FgMqXoaMDPP0RFHg\nKE2ny/X9/f210FDp8Dxv3z7w92fZmK70tPkRNk1nycAh9OolV2gLkRWjRxv55HQvqLaKiNOtqbx+\nF5w/D4UK6R1NZAOl1P/au/f4HOv/geOvz33vfG9jB3ZwGmPIec2ZTCiHEkUlJBH6KulIfX9JxVd9\nVZJSISk6F5KvU4mQcsycaYwYG8bOB9t9f35/XJMJ27Dtvnfv/Xw87sfu+zrd78+uw/u+Ptfnc13b\ntNZRRU0nVwTFjZk6lWw/H4blroDTDbjFcxgPPCCJWYjrNXGiCyE7J4HVg7tr7kenphrdFEWFIslZ\nXL+4OPjmG57rV50srzTcfp7KFwt85VnNQtwADw/4ak49WPMquz2P8M2d9YxuiufP2zs0UYYkOYvr\nN20aB6qYmBH0J+zpz0cv9CA01N5BCVH+dewID9QdColNeLj+CTJOxRvdFUWFIclZXJ+kJGwfzaFP\nn0qQ602HtFcYOFAagQlRUt5/rzKVNkwlwzuDZ3r4G92q5JaeFYYkZ3F93n+fWQ3Osz8kCbdfXmHR\n/PP2jcUAACAASURBVAZSnS1ECfL1hRWzboM/hvBBsxR2Je4yui2KCsHF3gGI8uFoUgaz1x9m8R8n\nyM3IZOGn/2XMw25wpBWfPDGMwEB7RyiE82nTRvF4gxeYkf0/evVR/DhuAvf8psjIycPi7kKfFqE8\n0rEOtQIs9g5VlDA5cxZFWnPgFN3fXs+Xm4+RnpNHv10/8Ui0L7muNrx2jyWoRbq9QxTCab39nwgC\ntkzkWOhp3q90hDpH9qKB9Jw8vtx8jO5vr2fNgVP2DlOUMEnOolBHkzL414LtZOVaybNpPHKzscSv\n4LdG8ajfHiWwkzv/WrCdo0kZ9g5VCKd07FwGlZs2gD13MT36HHdsn/33uDybJivXKvugE5LkLAo1\ne/1hcq22vz/32vI1E3omQ0Jjgut0QbnYyLXamLM+zo5RCuG8Zq8/DJYcKmcNgCx/nm+XRetDmy+Z\nRvZB5yPJWRRq8R8nyLMZLUQtmal8VuNPrG6Z+CaOwM3faAGWZ9Ms+kPu/ytEabiwD1Zq5I3H9mGc\nDT3CkaQNl7Tcln3Q+UhyFoXKyMkDwC8zBf9tX5AY8Qd+W/rhF1Hn0unO59kjPCGc3oV9EKBqVHu8\n997C5g7rab38O1ytuRenk33QqUhyFoWyuLsQeXwfd377HRs6/4jnoUgsTftdPp2bNPwXojRY3C/u\nW0qBT9WHMGdU5pvo9bww/0OqpyQa08k+6FQkOYur05qpf/3IsC+X82rvw3DenUo+wzGbLz0IuJgU\nfVtUs1OQQji3Pi1CcTFdvImAm6Uq3klD0H5xjG7vzuyPJ3Pboc2yDzoZSc7iys6dgz59CJ33I326\n3AQhO/A90x93n5qXTepqNjG8Y207BCmE83ukYx1czZceqisHR+P2Z2eymn5Pp8ZDmPLtuzyz+iPI\nzb3KUkR5I8lZXG7LFoiM5OiyPXRoPoa8du/icSySKsF9L5nMxaTwdDUzc1Ck3ARBiFJSK8DCzEGR\neLqaLzmDDqoxBtPxGpzuPoWosE/Q734M0dFw/Lj9ghUlRpKzuEhrmDED2rfnSG41WtVeQHrPMVRO\nDSVm8o8MaFUTb3cXlAJvdxcGtKrJirEd6Vy/qr0jF8Kpda5flRVjO16yD/p6uPPoLR9iyjVxpNcz\n3Fp3C+kxh6BFC1ixwt4hixuktJ1upB4VFaW3bt1ql+8WV5CSAsOHw7ffcujWR+gQ+yoJd7XD5HmC\nfU/EEBEaYe8IhRBX8NnGzxi0cjDs70vHw3NZmdUVz73b4IUXYOJEcJGGYo5EKbVNax1V1HRy5ixg\nxw6IioJFizjw7Bza732PhA5DoHIcn981TxKzEA5sYLuBjAofCTctZL3PVHr4/U7ywNEweTJ06wYn\nT9o7RHEdJDlXZFrDhx9CmzaQlcXuj7fQ6dOHOd3oBYhYyVM3Pcl97e6zd5RCiCLMHDyTtm5todNk\n1qV+RqedMzjx9tewebNRzf3zz/YOUVwjSc4VVXo6DBoEo0ZBdDQ/TttNh8ebczbsHWwd3yDaN5o3\n7nvD3lEKIYpBKcWaZ9cQej4U3Xsoe9PW0W5afw5+EwP+/tC1K7zyClit9g5VFJMk54po925o2RK+\n/BImT+aD3svoMaAyRHxG7u1jqafqsWrMKpQ8oFmIcsPdxZ0/xv+BxepJ3j1dOJN7gPZD6vL7e9uM\nH+IvvQQ9esApeYJVeSDJuaKZNw9atYLkZKyrVjP21As8OtpEtVbfkdJtMIF5gWwbtw1Xs6u9IxVC\nXKOqPlXZ+NhGXDwg4+4mKJ/DdOruyaz2n6BnzYb1641q7vXr7R2qKIIk54oiMxOGDjVebdty+qcY\nek2NZvp0aNX9e/7q0A8LFmKei8HH3cfe0QohrlPT0KYsG7gMU2UrZ3rUp0nLOEaOUgz/fTjZv2wC\niwU6d4bXXwebregFCruQ5FwR7N8PrVvDJ5/AhAmsfWEVzbpVZe1auGvIl2xu0Ac3dzc2j9lMaKVQ\ne0crhLhB3Rp048s+X0KAlT8a1WXAkD3MnQsdRzfl6MJt0K8fjB8Pd94JSUn2DldcgSRnZ/f550Y3\nqcRErMtWMlG9TJfbzPj6wqjnP+R7vwG4eruy7pF13BR0k72jFUKUkP4392dOjznYgm18ZWrK/03c\nzIED0KyDD5/d8QX63ffgp5+Mau7ff7d3uOIfJDk7q+xsGDkSBg6EyEgOLdrJrVO68fLLRtuQOx98\niekpo3CxuPDLsF9oXbO1vSMWQpSwh9s+zLtd38VW08Z/4try+rTlNGoEgwYrBqz/F+dWbDJuUtKx\nI0ybdskzooV9SXJ2RrGx0LYtzJqFbdzzvNN3DU1vC2bHDpg3T+MROpI3zryCm8WNDY9soG1YW3tH\nLIQoJaM7jmZG1xnYwmyM/r0Xj4xZwKRJ8N130GRwc1a9sRPuuAOeegruvhuSk+0dskCSs/P59luI\njIS//iJ21s9Eb/wPTzxlJjoaYmLyWPhbH2blzsLdx52NozbSupacMQvh7B7r8Bif3vEphMDQXwZj\n9XqFjRs1Pj5w+z3eDPRayKlXPoClS43jh9xa2e4kOTuLnBwYMwb69ye7YQtefegQTcZ0ZudOo/fU\n55+nct//tWdJwBJ8vX3Z8q8t3FztZntHLYQoI4OjBvPDgB8wB5p56chLvDP/QbZsyePFF+GbbxQN\npo3ko2f3Y8u1Qvv28N57Us1tR5KcnUFcnHHNaMYMVvb9gCZJa5nwVmV694Y9e6Bt24M0GFSfzeGb\nqe5ZnX1P76NJUBN7Ry2EKGO9GvZi/SPr8fDzYIHHAroO78Czz6YREwONG8PwKeFEV49lX5uh8Nhj\ncP/9kJpq77ArJEnO5d3330NkJMf3p9O/7XG6LxqJMilWrYKvvoLtO36gydNNSIhKoHVga/Y9t49Q\nH+kuJURF1bZGW2IejyHYN5hNEZuIGBSBzbaHtWthzhzYfcCVZr+9z4u3/krWt/8zenvExNg77ApH\nknN5lZsLzzxDbp9+TLVMpIF1N0v/qMakSbBrF9x6q5Xxk8fT+9venI86z6D6g9gwegPebt72jlwI\nYWcRAREcePYAbQLbkBCZQLMXmzFvwVyGDTNui3DvvYpJP7ejUdAZliW1Nu6TMHu2VHOXIUnO5dGx\nY9CpE+ve3EwL/6M8F/8Et3YxsXcv/PvfkJj4Fzf3u5nXk1/HXNPM7J6zmX//fFxM8lxXIYTB192X\nXx//lTHNx2BtZmXYhmH0G9kPb+9MFiwwHmTlXsmDXmfn07fSav4a8So8+KDx0BxR6iQ5lzfLl3Oq\nWTeGbHmMTqwjwzeUJUtgyRKoXRvmfz6fiBERxDSLIdQvlJjRMQxvOdzeUQshHJBJmZh+13S+vudr\nPEI9+C7wO8LvDee3336jc2ejNvu112BVejsausby+mfVON+yvdGYRZQqSc7lRV4etuf/zayei2iQ\nuokv1ABeeMHYR+68ExITE7n9odt5cP2D5LTNoW94Xw4+c5BGVRvZO3IhhIPr37g/sU/F0qJKCxJa\nJtBuRjsefeZR8vIyGTcO9u5V3NbLjfH6NZr/+Q1rbn7GuB2wKDWSnMuDEyeIaT2CDq/1YiSzaNbe\nm5gYxeTJ4O5uZfrM6dQaVotVtVbhE+zDt/2+ZeHghVjcLPaOXAhRTlTzrcbWJ7YypdMUzBFmPjB/\nQNg9Yaz+eTW1asGiRUY36Oxqdbg1ZzmDHjKTMOBJ46E6osQpbacL/FFRUXqrdHQvUsYPPzPhvgNM\nz3oEf5883nzPg0GDQCnYunUrAyYOILZOLATA3bXvZk7/Ofh5+tk7bCFEObbn1B4e+PwBdqbshHiI\nzoxmzstzCA8PJysLpky28fprNjysGUwKmcno1Xdjaljf3mGXC0qpbVrrqKKmkzNnR2W1sm7oxzTt\nXYu3sh5l+L2pHDjqweDBcOhQLLcPvZ2W77UktmUsIcEhrBy4ku8e/E4SsxDihjWq2ogdT+zgo14f\n4V3Nm7V11xLxdASjxo3i/PkUXplkYtdeF1pH5jLm5PN0apxE7LQf7B22U5Hk7IAyDifyRNj3dJo3\nFLx9+GVFFh985U9W1gnuf+x+Ip6PYFXYKrxqefFm5zc5Ou4ot9W9zd5hCyGciFKKh6Me5vi44wxv\nMhyaw4duHxI8NJhxk8YRFJTCyq2BzHvrLLtUE5o+1YXp7b/Glplt79CdglRr28HRpAxmrz/M4j9O\nkJGTh8XdhT4tQnmkYx0Svj7AwMf9OGStzeNd9jBl8U3En/iTcTPGseTsEmwRNly0C4+2eJRJ3Sfh\n6+5r7+IIISqAo8lHeXLRkyw+shht1bjtdmNEsxG8/MTLZKX6MKJLLMsON+QW723MXxZIzY61Cj3W\n1QqomG1iilutLcm5jK05cIp/LdhOrtVGnu3i/94FTYNvslh5uA/VXBKZ90E2eeFHGP/JeLaxDcLA\n3erOiBYjmHD7BAK9Au1XCCFEhXX43GHGLhrL0r+WotGYD5m53e92/jNsCjs+tvL49LqYsTJh2HZm\nhZy//FhnUriaTcwcFEnn+lXtWBL7kOTsgI4mZdD97fVk5VovGe59OgePL0LZltWG9j6rqT9qCYtP\nLOBs9bPgCZV1ZZ5u/zRjO42VO3wJIRzCibQTvLLsFT7Z8wnZ5mw4A7VSa/FI2GgWTu3J9qxGRFdZ\nxl8P5GH1MF82v6ermRVjO1a4M2hJzg7o/xbv4svNxy75FVl7Swq71nQhLeQgIc3GcTT8NwgEZVNE\neUUx4c4J9GzYE5OS5gFCCMeTk5fD3M1zmbp6KnG2OABM8WbC99/Nn7unUDctF797YjgVdumJhYtJ\nMaBVTV7t09geYduNJGcH1PillaTn5OGTk0HAmZ1kH4pjcxV3TOHLsVZKABuYkwLxc+vK/knvEuAV\nYO+QhRCi2I6cPcLkJZP5OOZLrJWN23yazoTjEtuVtmnJ2OpFkuhflxxXdwC83V3Y/fLt9gy5zEly\ntqfsbDhyxHiUY1wc5+Ni2XlyB98f+4Mjfqlsq6bZV9X4v5uzLJjP+uKu2uFbqQ8uJj+Ugrgpvexb\nBiGEuE61x/+P7LwjpCUvxspOcoJT0K45KJuieaKmRbyZ8LP+VMoLY1S3rphrhxv3H65dG2rUABfn\nfQ5AcZNzuf4P2K0lYF4eHD/+d/K98Er/K5aDyYfYzxkOBMD+QON1IBBy6gH1wDNTYTvRFnb0Icjq\niWvHapj93S5ZvMWtXK8WIUQFZ3F3QROGe+BYAHJTM8ncmEiy7ykO1lzEnqa7Oe92GjjNc+e30OR3\naPo91D8D4Skmwt2CqRNQF0utuheT9oVXcLBxF6ZS4igtzMvtmfNVWz2XREtArSExEeLisB0+xJm4\nPSTEHyDxVBwJKcdJyDlLopeNBG9I8IZEC5z0VSR5XoxD2cA9GfQZqEQVmge1IC2lM5u+fgyd40Hg\nHTvwiki87Ksr6nUYIYTzuFL7Gq0h9fdwktfVxzUkHq/WL+HtsZ006yHOe6WRFwR5XpcuJzhdUTdJ\nUz0VQtIhJA2Cc1wJ8Q4mJKAWwSER+NdqgKpT52Ly9rv+GzGVal7JV6LV2kqp7sB0wAzM0Vq/9o/x\nKn98TyATeEhrvb2wZd5Icr5aq+eCCmsJqLUmOfEoCQe3kXhkDwknDpKQdJSE1BMk5iSRoNNI8LSR\n6A2nLGC9Qlss1zyFe5YZa4qV7DSNTgeVpghxC6FpSFPaN2xPy8iWtGnThkqVKrF2LfS+S5Nly6FK\nvy24BaVec9xCCFEeFHaMztgfzJmlzXGtlM3qVYqON3uRkJDA5s2b+WXzL2zYu4GDpw+SrJLBH5Qf\nuFVS2Lwh1/XyfOWWB0EZUDUDgtIh6LwrVd38CLIEEeRXnapBdQiq0ZCg8KYERDTH7O1zzTFfUBLH\n5xJLzkopM3AQ6AYcB7YAA7TWewtM0xN4HCM5twama61bF7bcG0nOV/pVVpDGhpvtBN2CE2nsn8Cx\nU7HEp8ZzPOc0x0kl3j2HTNfL53O1QkC2Cz45LrhnmclLh9Sz5zl9LpfcdCD/5Z7nTniNcOrVrUd4\neDj169enefPmNG7cGC8vr8uW+8UX8NBDEB4OL757hlfWbJW+f0IIp1bYWWhuvD/J37fCw93E0qXQ\nsuXl8589e5a9e/eyZ88edu/ezd69e4n9K5bjKcexednAG/ABdx/w93fH3ceEzctGpoeVFPc8ci/v\nvYXJBoHZJoKsHgSZfAj1qEqdymHUCbmJTacDWRFfhVwVgOLK1eYlUbNZksm5LTBRa317/ufnAbTW\nUwpM8yGwVmv9Rf7nA0C01vrk1ZZ7I8n5QqvnCyx/biEtcxtB/r9wrnIG8X42sgskXxcrhKZBQBr4\npCrMqRpbKmSmwbl0SMyAtHQgC1xdXQkODiY0NJRq1aoRGhr69/uaNWtSt25dQkNDMZmK17Vp2jR4\n6im45RZYvNiocTmalMGc9XEs+iOejPN5WNxc6NuiGsM71pYzZiGE0yjsWJd12kKPHnDqFCxcCLcX\ns9F2Xl4ex48f58iRI8TFxREXF8eJEydISEi4+EpMwOpiBW+obIFAC/h6g4cFzBawekOWN5z2gXhf\n0AVysUcuhCab8TpblcykFhDWD6vfxROmG21hXpINwqoBxwp8Po5xdlzUNNWAqybnG5FRIDEDmNMS\nONZ2OVnnvKibpAmPc8E12ZXcNC+qV66O1TsYa0AgFm9vLDUs+Pv74+/vT0BAwGV/LRYLqgQaG2gN\nkyfDiy/CPffAZ5+Bu9F7gFoBFl7t01iuKwshnFqhx7oA+O036N7deCb9V19B375FL9PFxYWwsDDC\nwsKIjo6+4jQ2m42kpCSSkpI4d+4cycnJnDt37rL36UeTcUs+jSnnFIlpCbj65JDnm0uKn5VdVXNI\nr7eSsMT7Lll2xvm8K35nSSvTZsFKqRHACICaNWte93Is7i6XnDmnNO2G+Z3XSffP4WT/LZysYQz3\ndndhpR360GkNL7wAr70GgwfD3LlO3TNACCGuS3AwrFkDPXtC//4wbx4MGnTjyzWZTFSpUoUqVaoU\ne56CNbLaBpkzu+BV4zi23vsvqeQuq940xambjQdqFPhcPX/YtU6D1nqW1jpKax11Lf+0f+rTIhQX\n08V/l3LxwFL/DNlHArFmGfXZLiZF3xbVrvs7rpfNBmPHGol55EhjY5PELIQQV+bnB6tWGZf+HnwQ\nPvzQPnEUzCvZxwKwZXjgVT/lkprUsswrxUnOW4B6SqnaSik34H5gyT+mWQI8qAxtgJTCrjffqEc6\n1sHVfGnoXg1OgM1E5sFgAFzNJoZ3rF1aIVyR1QojRsA778CTT8L770MxL00LIUSF5eMD//ufcQY9\nahS8+WbZx1Awr2TuC0G55uEZfuqSacoyrxSZOrTWecBjwEpgH/C11nqPUmqUUmpU/mTLgMNALDAb\n+FcpxQsY1zFmDorE09X89y8dt6BUXPzSydofiqermZmDIsu0cVVurvGr76OPjOvMb75Zqv3khRDC\nqXh6Gg3D+veHZ56Bl182LhGWlQt5xcPkQubBEDzrJmJytQHGGXNZ55VyexMSuLwlYPqvDUj6tQ6b\nd2cS1bDsEnNODgwYAIsWwZQpMH58mX21EEI4FasVhg83Lgk+/TRMnVq2JzqffpPNkHs9qHHfdsy1\nT5Z4b5oKcfvOf7YE3LcPbroJfvreQlTDsokhKwvuvhtWrIDp02HMmLL5XiGEcEZms1EDabEYNZCZ\nmfDuu2V3iXDJVx74+cHBeZF4eJTNd16JU10RbdgQunaFGTPg/PnS/760NOMaycqVMHu2JGYhhCgJ\nJpNxHH/uOaPtzrBhxhl1aYuLM2pAR47ErokZnCw5g9EQ68QJ+Oab0v2e5GS47TZYvx4WLDCqYYQQ\nQpQMpYxeLy+/bFRxP/CA0banNL3zjvHD4LHHSvd7isPpknP37lC/Prz1Vuk1JkhIgM6dYds240fA\nAw+UzvcIIURFphRMmABvvAFff23c0Ck7u3S+KyXFqE6/916oVva9cC/jdMnZZDLOnrdvh59/Lvnl\nHzoE7dvDwYPwww/Fu6ONEEKI6/f00zBzpnHM7d3buA5d0j74wLhU+eSTJb/s6+F0yRmMLk21ahk3\nAynJapAdO4zEnJJiJP7i3gtWCCHEjXn0UaN6e/Vqo+YyIaHkln3ypHG75Z49IarIdtRlwymTs6en\n8cCJ3buNX1sl4ZdfoFMncHMzrjO3LvSZW0IIIUrakCFGX+jdu6FNG+NvSXjuOaNL7PTpJbO8kuCU\nyRmgTx/jzHbCBONX0Y2YPRu6dTOuQ/z6q9EqXAghRNm76y5Yt87okdO+vdFb5kasW2c06n32Wahb\nt2RiLAlOm5yVMlre5eYa/ZCzsq59GVlZRlXKiBHQpQts3Ag1ahQ9nxBCiNJz882waROEhUGvXvD2\n29fXAPjYMbj/fqhd23hYkSNx2uQMEBEB8+cbK/Ghh4yHUhTXrl3QqpXRSODZZ2HpUqhcudRCFUII\ncQ1q1IANG4zk/OSTcMcdkJhY/PnT0ox5MjJgyRLw8iq9WK+HUydnMJrev/660Qy/f384d67w6VNS\njJaBkZFw+jQsXw7//a9x1xohhBCOw8cHFi827iC2erXRjXb69KIbAsfGQnQ07NljdIdtfIXHTdub\n0ydnMG6iPnWq8euoeXP4+GNITb043mYzVtLTTxvVJNOmGWfau3YZ/aaFEEI4JqVg9GijN03r1kYv\nnXr1jJOyY8cure4+edI42WrR4uLdwG67zX6xF6ZcP/jiWm3aBEOHGvfg9vCAmjWN+7ceOmQkaxcX\n4/r0c88Z1zSEEEKUH1obj5586y1Ys8YYFhwMISFGNfbhw8bJWNeuMHeufdoQFffBFxUqOYOx8n7/\n3ajKiI+H9HQjSbdsCT16GCtRCCFE+bZvH/z0E2zZYlzO9PExzqgHDjTaI9lLhXgq1fVQCtq2NV5C\nCCGcU8OG5bvba4W45iyEEEKUJ5KchRBCCAcjyVkIIYRwMJKchRBCCAcjyVkIIYRwMJKchRBCCAcj\nyVkIIYRwMJKchRBCCAcjyVkIIYRwMJKchRBCCAcjyVkIIYRwMJKchRBCCAcjyVkIIYRwMHZ7ZKRS\n6jRwtAQXGQicKcHl2ZOUxTE5S1mcpRwgZXFUzlKW0ihHLa11laImsltyLmlKqa3FeUZmeSBlcUzO\nUhZnKQdIWRyVs5TFnuWQam0hhBDCwUhyFkIIIRyMMyXnWfYOoARJWRyTs5TFWcoBUhZH5SxlsVs5\nnOaasxBCCOEsnOnMWQghhHAK5So5K6X6K6X2KKVsSqmrtqBTSnVXSh1QSsUqpcYXGO6vlPpRKfVn\n/l+/son8ijEWGYtSqr5SakeBV6pSamz+uIlKqfgC43qWfSn+jrNY/1el1BGl1K78eLde6/ylrZjr\npIZSao1Sam/+tvhEgXF2XydX2/YLjFdKqXfyx+9USkUWd96yVoyyDMwvwy6l1EalVLMC4664rdlD\nMcoRrZRKKbDdTCjuvGWtGGV5tkA5diulrEop//xxjrRO5iqlTimldl9lvP33E611uXkBDYH6wFog\n6irTmIFDQB3ADYgBbsof919gfP778cDrdizLNcWSX64EjD5yABOBZ+y9Tq6lLMARIPBG/xf2LAcQ\nAkTmv/cBDhbYvuy6Tgrb9gtM0xNYDiigDbCpuPM6YFnaAX7573tcKEth25qDliMaWHo98zpaWf4x\n/Z3Az462TvJjuQWIBHZfZbzd95Nydeastd6ntT5QxGStgFit9WGt9XngS+Cu/HF3AZ/kv/8E6FM6\nkRbLtcbSBTiktS7JG7eUlBv9vzrKeikyDq31Sa319vz3acA+oFqZRVi4wrb9C+4CPtWG34HKSqmQ\nYs5bloqMR2u9UWt9Lv/j70D1Mo6xOG7k/1ru1sk/DAC+KJPIrpHWeh1wtpBJ7L6flKvkXEzVgGMF\nPh/n4sEzSGt9Mv99AhBUloH9w7XGcj+Xb+iP51e5zLVnFT3FL4sGflJKbVNKjbiO+UvbNcWhlAoD\nWgCbCgy25zopbNsvaprizFuWrjWeYRhnOhdcbVsra8UtR7v87Wa5UqrRNc5bVoodj1LKC+gOfFdg\nsKOsk+Kw+37iUhoLvRFKqZ+A4CuM+rfW+vuS+h6ttVZKlWpT9cLKci2xKKXcgN7A8wUGvw+8irHB\nvwq8CTx8ozEXEkNJlKWD1jpeKVUV+FEptT//F2xx579hJbhOvDEOPGO11qn5g8t0nQiDUqozRnLu\nUGBwkduaA9kO1NRap+e3U1gM1LNzTDfqTuBXrXXBs9PytE7szuGSs9a66w0uIh6oUeBz9fxhAIlK\nqRCt9cn8KopTN/hdhSqsLEqpa4mlB7Bda51YYNl/v1dKzQaWlkTMV1MSZdFax+f/PaWUWoRRRbSO\nMlwvJVEOpZQrRmL+TGu9sMCyy3SdXEFh235R07gWY96yVJyyoJRqCswBemitky4ML2RbK2tFlqPA\njzu01suUUjOVUoHFmbeMXUs8l9X0OdA6KQ677yfOWK29BainlKqdf8Z5P7Akf9wSYEj++yFAiZ2J\nX4drieWyazf5yeOCvsAVWx2WkSLLopSyKKV8LrwHbuNizI6yXopTDgV8BOzTWr/1j3H2XieFbfsX\nLAEezG+N2gZIya/KL868ZanIeJRSNYGFwGCt9cECwwvb1spaccoRnL9doZRqhXFcTirOvGWsWPEo\npSoBnSiw/zjYOikO++8npdHKrLReGAe840AOkAiszB8eCiwrMF1PjFa0hzCqwy8MDwBWA38CPwH+\ndizLFWO5QlksGDtqpX/MPx/YBezM3zhCHLksGK0bY/JfexxxvRSzHB0wqq13AjvyXz0dZZ1cadsH\nRgGj8t8r4L388bso0OvhavuNHberosoyBzhXYD1sLWpbc9ByPJYfZwxGw7Z25XWd5H9+CPjyH/M5\n2jr5AjgJ5GLklGGOtp/IHcKEEEIIB+OM1dpCCCFEuSbJWQghhHAwkpyFEEIIByPJWQghhHAw2H0p\nGgAAACRJREFUkpyFEEIIByPJWQghhHAwkpyFEEIIByPJWQghhHAw/w9Xe4BDh/LLmgAAAABJRU5E\nrkJggg==\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import numpy as np\n", "from numpy import interp\n", "from numpy import polyfit, polyval, poly1d\n", "from scipy.interpolate import CubicSpline\n", "\n", "%matplotlib inline\n", "import matplotlib.pyplot as plt\n", "\n", "n = 10 # sampling points - we will use this many samples\n", "\n", "# we want to interpolate this gaussian data\n", "σ = 0.4\n", "y = lambda x: np.exp(-x**2/σ**2)\n", "\n", "# exact solution\n", "xe = np.linspace(-1,1,200) # create equally spaced points between -1 and 1\n", "ye = y(xe)\n", "\n", "plt.figure(figsize=(8, 6))\n", "\n", "# sampling points\n", "xi = np.linspace(-1,1,n)\n", "yi = y(xi)\n", "\n", "# plot sample point locations\n", "plt.plot(xi,yi,'o',markersize=10)\n", "\n", "plt.plot(xe,ye,'k-',label='Exact')\n", "\n", "# linear interpolation. Interpolate to to xe using sampling points xi\n", "ylin = interp(xe, xi, yi)\n", "plt.plot(xe,ylin,'r-',label='Linear')\n", "\n", "# polynomial interpolation. Interpolate to to xe using sampling points xi\n", "p = np.polyfit(xi, yi, n-1)\n", "ypoly =polyval(p,xe)\n", "plt.plot(xe,ypoly,'b-', label='Polynomial')\n", "\n", "# cubic spline interpolation. Interpolate to to xe using sampling points xi\n", "cs = CubicSpline(xi,yi)\n", "ycs = cs(xe)\n", "plt.plot(xe,ycs,'g-', label='Cubic Spline')\n", "\n", "# finalize plot\n", "plt.legend()\n", "plt.draw()" ] }, { "cell_type": "markdown", "metadata": { "slideshow": { "slide_type": "slide" } }, "source": [ "More examples coming soon!" ] } ], "metadata": { "celltoolbar": "Slideshow", "kernelspec": { "display_name": "Python [default]", "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.6.0" }, "toc": { "colors": { "hover_highlight": "#DAA520", "navigate_num": "#000000", "navigate_text": "#333333", "running_highlight": "#FF0000", "selected_highlight": "#FFD700", "sidebar_border": "#EEEEEE", "wrapper_background": "#FFFFFF" }, "moveMenuLeft": true, "nav_menu": { "height": "84px", "width": "252px" }, "navigate_menu": true, "number_sections": true, "sideBar": true, "threshold": 4, "toc_cell": true, "toc_section_display": "block", "toc_window_display": true, "widenNotebook": false } }, "nbformat": 4, "nbformat_minor": 2 }