HW1/Homework 1.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "\n",
    "## Homework 1\n",
    "\n",
    "\n",
    "Spencer Pease  \n",
    "ASTR 427\n",
    "\n",
    "--------------------------\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Imports\n",
    "\n",
    "%matplotlib inline\n",
    "\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Problem 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5.551115123125783e-17"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Problem 1a\n",
    "\n",
    "itr = 1 # iterator\n",
    "\n",
    "epsilon = 2**(itr) # offset from 1 (in terms of general floting point formula)\n",
    "result = 1.0 - epsilon # measured difference\n",
    "\n",
    "while result != 1:\n",
    "    itr -= 1\n",
    "    \n",
    "    epsilon = 2.0**(itr)\n",
    "    result = 1 - epsilon\n",
    "\n",
    "    \n",
    "epsilon # answer: smallest measuralbe offset below 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1.1102230246251565e-16"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Problem 1b\n",
    "\n",
    "itr = 1 # iterator\n",
    "\n",
    "epsilon = 2**(itr) # offset (in terms of general floting point formula)\n",
    "result = 1.0 + epsilon # measured difference\n",
    "\n",
    "while result != 1:\n",
    "    itr = itr - 1 \n",
    "    \n",
    "    epsilon = 2.0**(itr)\n",
    "    result = 1 + epsilon\n",
    "\n",
    "    \n",
    "epsilon # answer: smallest measurable offset above 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "8.98846567431158e+307"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Problem 1c\n",
    "value = 2.0 # current max value\n",
    "\n",
    "while value * 2.0 < np.inf:\n",
    "    value *= 2.0\n",
    "\n",
    "inc_fac = 1.0 # incrementing factor: value to multiple against current max\n",
    "\n",
    "while value * inc_fac < np.inf:\n",
    "    inc_fac += .000000000001\n",
    "    inc_fac *= inc_fac\n",
    "    \n",
    "maxnum = value \n",
    "maxnum # maximun representible positive number"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5e-324"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Problem 1d\n",
    "value = 2.0 # current min value\n",
    "\n",
    "while value / 2.0 > 0:\n",
    "    value /= 2.0\n",
    "\n",
    "minnum = value\n",
    "minnum # minimum representable positive number"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The smallest change above 1 is greater than the smallest change below one because the difference between representable numbers is smaller for lower values.\n",
    "\n",
    "It is easier to represent the minimun value than the maximun value becuase the minimum positive representable positive number has a non-fractional binary representation, allowing us to find the exact value. The maximim representable value, however, has a repeating factional representation in binary, so it is difficult to exactly find through emperical methods."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Problem 2\n",
    "\n",
    " Numerically:\n",
    " \n",
    " $\\lim _{x \\to 0} \\frac{1 - \\cos{x}}{x^{2}} = \\frac{1}{2}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false,
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAEICAYAAABPgw/pAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmYFNXVP/DvGWGUccsLKCJEXIiyRcAgKqg0bgOKGRRN\nQBIUI+CCC0aFuITBGOP6qhHjigvK4qtixFdF8NUxwRVRUVkElx8CCopRkdVhOL8/Trc2TS/VPXWr\npqu/n+fhYbq76tZtmD51+txbt0RVQURExa8s7A4QEZE/GNCJiCKCAZ2IKCIY0ImIIoIBnYgoIhjQ\niYgiolGQBxMRzpEkIiqAqkqubQLP0FU1sn/Gjh0beh/4/vj+Su29lcL784olFyKiiGBAJyKKCAZ0\nH8VisbC74BTfX/GK8nsDov/+vJJ86jP1PpiIBnk8IqIoEBFoQxwUJSIiNxjQiYgiggGdiCgiGNCJ\niCKCAZ2IKCIY0ImIIoIBnYgoIhjQiYgiImdAF5EJIrJKRN7Lss3fRWSJiLwrIl387SIREXnhJUN/\nAEBlphdFpC+A/VT1FwBGALjLp74REVEecgZ0VZ0N4Jssm1QBmBjf9g0Au4pIC3+6R0REXvlRQ28F\nYFnS4xXx54iIKECB3rEIAKqrq3/8ORaLRXKVNFWgthYoLw+7J0RUjGpqalBTU5P3fp5WWxSRNgCe\nVtUD07x2F4CXVPXR+ONFAHqp6qo025bEaovPPw/cfTcwbVrYPSGiKPB7tUWJ/0lnOoAh8YMeCuDb\ndMG8lKxbB7z9dti9IKJSk7PkIiKTAcQANBORzwCMBVAOQFX1HlV9VkSOF5GPAKwDMNRlh4uBKrB0\nKbB2LbDTTmH3hohKRc6ArqqnedhmpD/diZYFC4Du3cPuBRGVCl4p6kBimGD+/HD7QUSlhQHdIQZ0\nIgoSA7oDqsAOOzCgE1GwGNAdad+eAZ2IgsWA7oAqsM8+wNdfA2vWhN0bIioVDOiObLedZekLFoTd\nEyIqFQzoDnXsyLILEQWHAd0BVUCEAZ2IgsWA7lCnTsD774fdCyIqFQzoDiQy9M6dgXnzfrrQiIjI\nJQZ0h/bc04L5ypVh94SISgEDugOJDF0E6NIFePfdsHtERKWAAd2xzp0Z0IkoGAzoDiTXzJmhE1FQ\nGNAdkfjtQBjQiSgoDOgOJGfoBxwALFtmdzEiInKJAd2RRIbeuDHQoQPnoxORewzoDqTOO+fAKBEF\ngQHdEUm6pTbr6EQUBAZ0B1Iz9C5d7IpRIiKXGNAdSc7QDzzQauh1deH1h4iijwE9ALvuCrRsCSxa\nFHZPiCjKGNAdSFz6n+zgg4E5c8LpDxGVBgb0gHTvDrz5Zti9IKIoY0B3gBk6EYWBAT0gXbva3Ys2\nbQq7J0QUVQzoDqS7oUVFBbD//py+SETuMKA7klpyAVh2ISK3PAV0EekjIotEZLGIjE7z+i4iMl1E\n3hWR90XkDN97WkQy3XKOA6NE5FLOgC4iZQDGA6gE0BHAIBFpl7LZeQDmq2oXAL0B3CwijfzubDFh\nhk5EQfOSoXcHsERVl6pqLYCpAKpStlEAO8d/3hnA16q62b9uFpdMGXrHjsBnnwFr1gTbHyIqDV4C\neisAy5IeL48/l2w8gA4i8jmAeQAu9Kd7xStdht64sa28OHdu8P0houjzqyxSCeAdVT1KRPYDMEtE\nDlTVtakbVldX//hzLBZDLBbzqQsNR6YMHQAOOQR47TWgd+/g+kNExaWmpgY1NTV57+cloK8AsFfS\n49bx55INBfA3AFDVj0XkUwDtALyV2lhyQI+ydBk6APTsCdx/f7B9IaLikprsjhs3ztN+XkoucwC0\nFZE2IlIOYCCA6SnbLAVwDACISAsA+wP4xFMPIihbht6zp2XoW7YE1x8iKg05A7qq1gEYCWAmgPkA\npqrqQhEZISLD45tdA6CHiLwHYBaAy1T1P646XQwyZeh77AE0bQosXBhsf4go+jzV0FV1BoADUp67\nO+nnL2B1dPKgZ09g9myb9UJE5BdeKepAusW5kh1+OPDKK8H1h4hKAwN6CBIZOhGRnxjQHciVobdr\nB3z3HfD558H1iYiijwE9BGVlQI8eLLsQkb8Y0B3INm0xgXV0IvIbA7oj2UougAX0l18Opi9EVBoY\n0B3wkqEffDDw8cfA11+77w8RlQYGdEdyZejl5TbbhVk6EfmFAd0BLxk6ABx1FPDii277QkSlgwHd\nkVwZOsCATkT+YkB3wGuG3qULsHIl8MUXbvtDRKWBAd0RLxn6dtsBvXoBL73kvj9EFH0M6A54zdAB\nll2IyD8M6I54ydABBnQi8g8Desg6dADWrwc+KdnbgRCRXxjQHci1OFcyEaCyEpgxw22fiCj6GNAb\ngL59geeeC7sXRFTsGNAdyGdQFACOO86uGN240U1/iKg0MKA74rXkAtg9Rjt1Av79b3f9IaLoY0B3\nIN8MHQiv7LJpE/D998Efl4j8x4DuSD4ZOmABPYyB0YceAgYNCv64ROQ/BnQHCsnQDzoIWL0aWLrU\n//5ks2kT8MwzwJIlwR6XiPzHgO5Ivhl6WZlNX3z2WTf9yWTLFluC4I47gj0uEfmPAd2BQjJ0AKiq\nAp56yt++5KJqx504kbV0omLHgO5Ivhk6YBn6q68C333nf38y2bIFaNMGOPpoq6cTUfFiQHeg0Ax9\n552BI48MtuySuKr1gguA22+3AE9ExYkB3ZFCMnQAOOkk4J//9Lcv2SQC+uGHAxUVwMyZwR2biPzF\ngO5AoRk6AJx4IvD888FdNbpliw3IigAXXgjcfHMwxyUi/3kK6CLSR0QWichiERmdYZuYiLwjIh+I\nSMnfsqHQDH333YFf/jK4JXWTFxI77TRg0SJg7txgjk1E/soZ0EWkDMB4AJUAOgIYJCLtUrbZFcAd\nAPqpaicApzroa8k46STgySeDOVYiQweA8nLg4ouB668P5thE5C8vGXp3AEtUdamq1gKYCqAqZZvT\nADyhqisAQFVX+9vN4pLP8rnpJOrotbX+9SmT1L4OG2a3xOOFRkTFx0tAbwVgWdLj5fHnku0PoKmI\nvCQic0Tk9351sBTtsw+w337BlF1SA/pOOwHnngvcdJP7YxORv/waFG0E4CAAfQH0AXCViLT1qe2i\nU59B0YSBA4GpU+vfTi7JJZeE888HHnsM+OILd8d94w1g8WJ37ROVokYetlkBYK+kx63jzyVbDmC1\nqm4EsFFE/gWgM4CPUhurrq7+8edYLIZYLJZfj4tEfUouAHDqqcDVV9taK9tv70+f0klXHmreHPjd\n74D//m/gxhvdHHfCBOCdd4DXX7elB4joJzU1Naipqcl7P9Ec6aSIbAfgQwBHA/gCwJsABqnqwqRt\n2gG4HZadbw/gDQC/VdUFKW1pruNFwY03AqtW1b9sEYsBF10E9O/vS7fSSpxfk86zAIDly4HOnYEF\nC4AWLfw/7lln2ZWpt91mJR4iykxEoKo508ScJRdVrQMwEsBMAPMBTFXVhSIyQkSGx7dZBOB5AO8B\neB3APanBvNTUN0MHgim7ZBrAbd0aGDwYuOEGN8fdssWuTh07Fli50s0xiEpNzgzd14OVSIZ+ww3A\nV1/Vv1yxerUNjq5YYYOVLvz5z0CjRvZ3qs8/tzspzZ8PtGzp73FPPx3o3dvmvS9bBkya5G/7RFHi\nW4ZOhfEjQ2/eHDjiCOCJJ+rfVibpBkUT9twTOOMM4G9/c3Pc7bYDrroKeOUV4P/+z9/2N24ENmzw\nt02iho4B3QE/v4QMHQo88IB/7aXKNWd+9GjLnpcv9/e4iRPJjjsC48cDZ58NrF/vX/t33WUnw02b\n/GuTqKFjQHfEjwwdsLVdFiwAPv7Yn/ZSbdmSva8tWtgA5jXX+Hvcurqfvhn06wccfLBl635Ztw54\n7z1g1Cj/2iRq6BjQHfAzQy8vtzVWHnzQvzaTqWYuuSSMHm1ln4ULs2+Xj0TJJeHvfwcmT7b14P1Q\nVweccw4wa5a164clS4Bvv/WnLSIXGNAd8StDB6zs8tBDFqT85mWZgqZNgTFjLLD7JbV237y5lV7O\nPNOf2nddnfX78cdtFUk/TkaXXgocf7y/pSEiPzGgO+D3RJ7OnS3guVgKINugaLKRI4H33wcKuNYh\nreSSS8KAAcCBB9pUxvravNm+AXTubIuNDRgArF1bvzZra2320uDB9Tu5fv45cOWV1kciPzGgO+Jn\nhg5Y5nr//f62CXhfSGz77W22yyWX+HNXo9SSS8L48cDDD9vMl/qoq/up/TPPBA491MYC6nOyra21\nq2e//97m0Bfa1vz5wF//alM3Cw3qqja2QpSMAb1InHYaMGMG8OWX/raba1A02W9/a1n1lCn+HDfd\nN4PddwfuvtuWHqjPvVXr6mx+fcI//gF8+ilw7bWFt7l5s93Vado0O+Fcd11h7dTWAr16WbZ/+umF\nZfuLFgEdO9rYQyG++cbdQDuFhwHdARfXTjVtamWD++7zt10vg6IJIpahjhlT//JFupJLwq9/DfTt\nW78lARIll4QddrAlie+6q/Bb/G3ebCeJXXax+77efTcwcWL+7dTWWhtPPWUn6EKC+oYNwF57AXfc\nAVxxRf6/c1OmWHlr2rT89gOAZ54BXnst//04hdQ9BnRH/C65AMB55wF33ulv7TXftdsPP9yu8PzL\nX+p33Ewll4SbbrLFux55pLD2k0suCS1b2o1Dhg+3KY352rwZaNzYft5zT+C552ygON+bkdTWWjtN\nmlhQX7nSLuDKJ6j/8IO9n9mzbSbPsGH5/V5s2AD07GkDxtddl98J4b777BvG7bfnt1/btrY2UT6B\nfdw44F//8r49YN98ShUDugOuVjfo2tWysqef9q9Nr4OiyW64wer59Zk5kuu4FRU23XDUKOCTT/Jv\nP7XkktCtmy0IVlWV/we/tnbrNtu3t0z97LMtuOfTTuLEUFEBTJ9uSxX/7nfeb2ryww82pXW33Wyw\nfNky4JRTvM8Q2rQJ+NWvbLXLxx6zE4rXQLtxo31TmzABGDLE+6yf1auBjz4CevSwv72YMMH+r8aM\nsfecy/LlwB572J23vPRr0yabuTRrlrf+PPqora/k5TM+b54lJV58+qk/iRoDuiMuMnTAsvTx4/1r\nr5C7K+2xh10ENHJk4SevbCWXhC5dgD/9ycYPvHyYU9vP9A1g0CCbqVJVld8UyUTJJVnXrlbCOf10\n7zOAEsE4oaLCTtLr1tnKml4CUXIbO+1k+1dUAMcdB3z9de79E8syt2plGfDatcAxx3g7yW3caPX7\nV1+1353DDssdoFXtmNOn2zTcww7zNhazfr19C1m4EDjkEBtQzubbby3pWbnSfn9mz86+/TffAC+/\nbN9whgyxk042kyZZKfDYY4EPP8y+7S23AN27AyNG5G43FrP+zpyZfbtcGNAdcLn+2Cmn2C+1Xxf5\n5DMomuzcc+2X9LHHCj+ul3XQL7rIrla95JL82k+toae6+mpg333tZOG11JEuoAMWnP7nf2wNey+1\n5eQMPaFJE7t467/+y8YP1qzJ3kbqSaG83MpTPXrYjJ5cwSZ5nf0dd7T/xyOPtAD09tvZ99240cYk\nKirs+ojhw+242b45btxofSwrs0Rg5kybnnrWWdlPYOvXA3vvbSfN886zwHfbbZlnWq1bZ99aJk+2\n6aq/+Y19y8t0jLVrrXT1wQdAs2a2GN0jj2T+DK9dayeifv2sZHXVVZmTgjVrbMymSROgQwdLxDJl\n4f/5j5Xvzj3X2l60KPO/STYM6I64ytDLy+0DdNtt/rSXz6BoskaNbEDu4osLm43itdRTVmZB45ln\n8ptdky1DT7R7//32obvwQm8n4eQaeqpYzAZI+/cH3nwzezvpAjpgz02caEGld+/s2XJqQAfsPV1/\nvZUnjjzS7g2bSeqNU8rKbCrljTcClZXZ1w9KBHTAfs/PO++ngHvxxelLNxs2WGBL6NoVmDvXtu3W\nzX5OtWWLHatJEzvOWWfZCXPqVOvjZ59tu8+6dXaCAuzevO+/bwPPnTsD//73ttuvXWvb77STZdRP\nP23jN336WBkk1fff20n3oouspLJ4sX1beeaZbbdds8a+Ldx6q5XFpk0DDjrIvhEkq6uzE87gwZas\n9e5t6xBdcIG3b1vJAl8+F9WBHY6IKBqq4Wn5XKhqYH/scNF39dWqV17p9hhnn616xRX1b2f4cNU7\n7yx8/2+/VW3dWvXFF/Pb77DDVGfPzm+f++9XPeAA1TVrcm87aJDqI494a3f5ctW99lJ96KHs27Vp\no/rJJ7nbmzFDtXlz1ZdfTv/6TTepjhqVu51bblFt1Up17txtX3v4YdXBg7Pvv2iR6n77qV52merm\nzVu/duaZqvfem3nf775T7d9f9ZBDVJct2/q1vffO/O+wZYvq+PH2/pP//ZcsUd1338zHW7ZM9dhj\nVbt3V/3wQ3vuq69UmzbNvM8HH9j2vXurfvSRPTd1qupvfpN++2++UT3nHNWWLVWnTLG+Pvec6nHH\npd/+449VTzjBfudmzbLn9txz238PVdUNG1T/8hfr7xVXqK5dq9quneqCBdtuu26d6p//bNuOGaP6\n3nuqP/955vdYWakaj505YyxLLg4E8aXnkktsHnR954MXWnJJ2HVXqxOedZZ93fXKaw092dCh9lV0\nyJDcV6vmKrkka9XK2xTEbCWXZJWVVhYYMCD9Ou+ZSi6pLrrISmuVlTabJtkPP+S+1+wBB9gsljff\nBE44weq0CbnuVbvLLlbTr6qylTCTB+uSSy6pEiWYWbNsausZZ9jvaGrJJVXr1nbh3JAhVo+/4w4r\nQ1RUZN4nMTDbr58NmN56q5VEEiWXVD/7mV1g9sQTVl7q29emr2a6ecy++wL/+79Whho2zOrx334L\n7LzzttvusIMt5zBvnpVq2rWzmUfptq2osOmY8+bZ7KaePe1zlOk9zpiR+d8gFQO6I65q6An77Qcc\ndRRw7731a6fQQdFkJ5xgH8Irr8zvuIWcSMaPt8HYXEvtZpq2mEmHDlYHHTEi8wco06BoOkcfbYFj\n0KD0wTi1/p3JgAE2M+QPfwDuuSf/Npo3t+DasaMF5nnz7HkvNx8vK7NZRpMm2cn0iivs3yBbQE/o\n0gV46y373era1WYAZQvoieOdd55dhZsYj8gW0AE7aV98sdXWn3zS23z8ww6zgd/eve1eurnuBpZY\nwrp9ezspZDphAHZimjTJTuj9+9sAbbZtH3zQ/m3y+exk5SWN9+sPSqTkUl2tetVV7o/z1lv2VW3T\npsLbOPNM1fvuq39fVq+2r7JeyygHHaQ6Z05hx/ryS9V99sleUqmqUp02Lf+2X3lFdbfdVGtqtn2t\nWTMrA+TjtddUW7TYupxz5ZWq48bl186SJapt26qOHq1aV6d6222qF1yQXxuTJ1spZPJk1RNPVP3n\nP73vu2qVlSZ69lRt3Fh1/Xrv+z7+uJUXjjjC+z61tarXXJO7rJSsrs7KPTNmeN/n009VFy/2vn1Y\nwJJLuFxn6IBdGNKuXWGXnycUMg89nWbN7CrW3/8+95Q7oLCSS8Juu1nWOmoU8MYb6bfJp+SSrEcP\ny65OPXXbtlMvLPLi0ENttslVV9nsiUQ7Xkouydq2tSz0lVeAk0+22Q9es/yEQYOAF16wbPDpp/Pb\nf/fdrSzVr59lqLmy+2QDBti0wHzusduokX0jyOdK4USGX1npfZ+99wZ+8Qvv2zd0DOgOBHkf7Opq\nu5tQvhfeJPhRckmoqrILW84/39tx61O779TJph2efDKwdOm2rxca0AErZT34oH3VTp5X7rWGnqp9\ne7vA5YEHbE31H34orJ3mza0m37y51YDzDeiATd+bM8fWy2nXLr99y8psSuTq1fn/37VsaXVucosB\n3ZEgMnTAMsr27QtfWre+g6Kpbr7ZBuKmTs2+nZcrRXPp188CTGXltvN1862hpzr+ePvmU1X10/zl\nfGroqX7+c2vn1VdtvnOhystt3OSOO6wGXIimTW0NmX32KWz/Qk+U5B4DegSMG2cZ28aN+e/rZ4YO\n2NfxyZPtooh0F34kH9ePwHD++XYBSb9+W8+yyXWlqBd9+tjFTAMG2IUh9QnogAXSWbOsVNa2beHt\niNjg7THHFN4GRRMDugNBllwAu1y7S5fCltb1O0MHLGD98Y92WX2mxabqW3JJdu21wP7723rtiePV\np+SS7Oij7bL4gQP96XNFhc3+6N+//n0jSsWA7khQJZeEceMssOU7L92vQdFUl15qc2vHjEn/uh8l\nlwQRO5lt2WKZq2r9Sy7JevWyS9sPPzz4/1eifDCgOxB0hg7YGhG9e+c3kwDwv+SSUFZmt5KbNs1u\n1JzuuH7WYhs3tkx6wQK76MqPkkuyHj3SrwVC1JAwoDsSRiZ37bV24c2KFd73cVFySWja1ILsOefY\nIkbJ/Cy5JOy4o13E88ILNpODg3dUahjQHQgjQweANm3sSrlcV1Emc5WhJ3TrZgO2AwZsXQ7ys+SS\nLDHw2LOnrdtOVEo8faREpI+ILBKRxSIyOst2B4tIrYic7F8Xi1NYtdY//cmy1Hff9ba9yww9Ydgw\nm4M8ePBPa7D4XXJJtvvutkRpmzZu2idqqHJ+lEWkDMB4AJUAOgIYJCLbXJIQ3+46AM/73cliE1aG\nDthA5NixNm3QSz9cDYomE7FFkb77zk44gJuSC1Gp8/KR6g5giaouVdVaAFMBVKXZ7nwAjwP40sf+\nFa0wZ0MMH26r23lZEsB1ySWhvNwWq3riCbsK01XJhaiUeZnY1QrAsqTHy2FB/kcisieA/qraW0S2\neq0UhZmhA1bKuPNOu9jmxBOtrpxJECWXhGbNbA2RXr3sbjwctCTyl18f5VsBJNfWS362btjzlbt1\ns4HIyy/Pvl1QGXpC+/Y/LbhUyFokRJSZlwx9BYC9kh63jj+XrBuAqSIiAJoD6Csitao6PbWx6urq\nH3+OxWKIxWJ5drnhCztDT/jrX22d76FDMy+MFEQNPdVxx9mNAjIt6k9U6mpqalBTU5P3fjnvKSoi\n2wH4EMDRAL4A8CaAQaqa9r7zIvIAgKdVdVqa1zTX8aLg8stt0fxc2XEQpkyxwD53bvolT0880Wru\nJ54YfN+IyBsR8XRP0ZwlF1WtAzASwEwA8wFMVdWFIjJCRIan2yXv3kZMQzpnDRxo65wkfTHaStAl\nFyJyx9NqF6o6A8ABKc/dnWHbM33oF/lExAZIO3e2BaFSSy9BDooSkVv8KDvSkLLeFi2A228HTj/d\npjMmY4ZOFB0M6A40pJJLwqmnWpaeuLAnIYxBUSJygwHdkYYYJO+80+6MPj1p7hFLLkTRwY+yAw0x\nQwfsAqMpU2xtlWXxS8VYciGKDgZ0RxpqkOzRAxg1yu4Av3kzM3SiKOFH2YGGmqEnXHaZrR1+5ZXM\n0ImixKebdFGqhhwky8rs8vtu3YA1axp2X4nIO2boDjT0DB0AdtvNVj7csIGLZBFFBTN0R4oh6+3W\nDXj7bbuSlIiKHwO6A8WQoSd06BB2D4jILyy5OFIMGToRRQsDugPFlKETUXQwoBMRRQQDugNcH4WI\nwsCATkQUEQzojjBDJ6KgMaA7wEFRIgoDA7ojzNCJKGgM6A4wQyeiMDCgO8IMnYiCxoDuADN0IgoD\nA7ojzNCJKGgM6A4wQyeiMDCgO8IMnYiCxoDuADN0IgoDA7ojzNCJKGgM6A4wQyeiMDCgExFFBAO6\nA1w+l4jC4Cmgi0gfEVkkIotFZHSa108TkXnxP7NF5Jf+d5WIiLLJGdBFpAzAeACVADoCGCQi7VI2\n+wTAkaraGcA1AO71u6PFhhk6EQXNS4beHcASVV2qqrUApgKoSt5AVV9X1e/iD18H0MrfbhYXDooS\nURi8BPRWAJYlPV6O7AH7LADP1adTUcAMnYiC1sjPxkSkN4ChAA7PtE11dfWPP8diMcRiMT+70CAw\nQyei+qipqUFNTU3e+4nmiD4iciiAalXtE388BoCq6vUp2x0I4AkAfVT14wxtaa7jRcE55wAHHmh/\nExHVl4hAVXN+7/dScpkDoK2ItBGRcgADAUxPOdhesGD++0zBvJSUwDmLiBqgnCUXVa0TkZEAZsJO\nABNUdaGIjLCX9R4AVwFoCuAfIiIAalW1u8uON3SsoRNR0DzV0FV1BoADUp67O+nnYQCG+du14sUM\nnYjCwCtFHWGGTkRBY0B3gBk6EYWBAZ2IKCIY0B3g4lxEFAYGdCKiiGBAd4AZOhGFgQGdiCgiGNAd\nYYZOREFjQHeA0xaJKAwM6I4wQyeioDGgO8AMnYjCwIDuCDN0IgoaA7oDzNCJKAwM6I4wQyeioDGg\nO8AMnYjCwIDuCDN0IgoaA7oDzNCJKAwM6EREEcGA7gAX5yKiMDCgExFFBAO6A8zQiSgMDOhERBHB\ngO4IM3QiChoDugOctkhEYWBAd4QZOhEFjQHdAWboRBQGBnRHmKETUdAY0B1ghk5EYfAU0EWkj4gs\nEpHFIjI6wzZ/F5ElIvKuiHTxt5vFhxk6EQUtZ0AXkTIA4wFUAugIYJCItEvZpi+A/VT1FwBGALjL\nQV+LBjN0IgqDlwy9O4AlqrpUVWsBTAVQlbJNFYCJAKCqbwDYVURa+NpTIiLKqpGHbVoBWJb0eDks\nyGfbZkX8uVWpjT39dJ49LEIrVrDkQkTB8xLQffXHP1b/+HOzZjE0bx4LugvO/exnQKdOYfeCiIpV\nTU0Nampq8t5PNEfBV0QOBVCtqn3ij8cAUFW9PmmbuwC8pKqPxh8vAtBLVVeltKW5jkdERFsTEahq\nzu/9XmrocwC0FZE2IlIOYCCA6SnbTAcwJH7gQwF8mxrMiYjIrZwlF1WtE5GRAGbCTgATVHWhiIyw\nl/UeVX1WRI4XkY8ArAMw1G23iYgoVc6Si68HY8mFiChvfpZciIioCDCgExFFBAM6EVFEMKATEUUE\nAzoRUUQwoBMRRQQDOhFRRDCgExFFBAM6EVFEMKATEUUEA7qPClnuspjw/RWvKL83IPrvzysGdB9F\n/ZeK7694Rfm9AdF/f14xoBMRRQQDOhFRRAS+fG5gByMiihAvy+cGGtCJiMgdllyIiCKCAZ2IKCIC\nC+gi0keOKgkUAAAEYUlEQVREFonIYhEZHdRxgyAiE0RklYi8F3Zf/CYirUXkRRGZLyLvi8gFYffJ\nTyKyvYi8ISLvxN/f2LD75IKIlInI2yKSeoP3oici/09E5sX/D98Muz9+E5FdReQxEVkY/xweknHb\nIGroIlIGYDGAowF8DmAOgIGqusj5wQMgIocDWAtgoqoeGHZ//CQiewDYQ1XfFZGdAMwFUBWV/zsA\nEJEKVV0vItsBeAXABaoaqcAgIqMA/ArALqr667D74ycR+QTAr1T1m7D74oKIPAjgZVV9QEQaAahQ\n1TXptg0qQ+8OYImqLlXVWgBTAVQFdGznVHU2gEj+MqnqSlV9N/7zWgALAbQKt1f+UtX18R+3B9AI\nQKRmCohIawDHA7gv7L44Ioho+VhEdgFwhKo+AACqujlTMAeC+0doBWBZ0uPliFhQKAUisjeALgDe\nCLcn/oqXI94BsBLALFWdE3affHYLgEsRsRNVEgUwS0TmiMiwsDvjs30ArBaRB+Ils3tEpEmmjSN5\nViP/xcstjwO4MJ6pR4aqblHVrgBaAzhERDqE3Se/iMgJAFbFv2VJ/E/U9FTVg2DfQs6Ll0CjohGA\ngwDcEX+P6wGMybRxUAF9BYC9kh63jj9HRSBet3scwMOq+lTY/XEl/lX2JQB9wu6Lj3oC+HW8zjwF\nQG8RmRhyn3ylql/E//4KwJOwEm9ULAewTFXfij9+HBbg0woqoM8B0FZE2ohIOYCBAKI22h7V7AcA\n7gewQFVvC7sjfhOR5iKya/znJgCOBRCZAV9VvVxV91LVfWGfuxdVdUjY/fKLiFTEvz1CRHYEcByA\nD8LtlX9UdRWAZSKyf/ypowEsyLR9o4A6VSciIwHMhJ1EJqjqwiCOHQQRmQwgBqCZiHwGYGxiEKPY\niUhPAIMBvB+vMyuAy1V1Rrg9801LAA/FZ2KVAXhUVZ8NuU/kXQsAT8aXFWkEYJKqzgy5T367AMAk\nEWkM4BMAQzNtyEv/iYgigoOiREQRwYBORBQRDOhERBHBgE5EFBEM6EREBfJ7YT4RuV5EPogvwnVr\nvvszoBMRFe4BAJV+NCQihwHooaqdAHQC0F1EjsynDQZ0IqICpVuYT0T2FZHn4mvLvJx0UVDO5gDs\nICI7AGgCm1e/Kp/+BHJhERFRCbkHwAhV/VhEugO4E3aFZ1aq+rqI1AD4Iv7UeFX9MJ8DM6ATEfkk\nvvxADwCPiUhiKZDG8ddOAnA1tl71UgAsV9W+IrIfgHYA9ow//4KIzFDVV7wenwGdiMg/ZQC+ia+M\nuBVVfRK2eFgmJwF4XVU3AICIPAfgMNhNVzwfnIiICvfjwnyq+j2AT0XklB9fFPF6F7PPAPQSke3i\n67b0gt1QxjMGdCKiAsUX5nsVwP4i8pmIDIUtZvcHEXlXRD4A4PWWf4/DFt96H8A7AN5R1Wfy6g8X\n5yIiigZm6EREEcGATkQUEQzoREQRwYBORBQRDOhERBHBgE5EFBEM6EREEcGATkQUEf8fbSUrEU5E\nJf4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x843fdd8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "x = np.linspace(0, 6e-8, 1000) # domain to plot fucntion over\n",
    "func = ((1 - np.cos(x)) / (x**2)) # function to plot\n",
    "numerical = 1000*[.5] # line with the value of the numerical limit\n",
    "\n",
    "plt.xlim(0, 6e-8)\n",
    "plt.ylim(-0.1, 1.1)\n",
    "plt.plot(x, func)\n",
    "plt.plot(x, numerical)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the analytical solution, when x decreases past 1e-7, it becomes apparrent that the computer is unable to accurately represent continuous values, as seen in the sawtooth behavior of the function. Also, the function drops to 0 for x values lower than 1e-8 because the function involve x^2, and the computer can't distinguish between number changes less than ~1e-18, so it approximates the function to 0. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Problem 3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Problem 3a\n",
    "\n",
    "# Reads in a data file with the given file path\n",
    "# Returns an array of defined coordinates\n",
    "def read_dat(filename):\n",
    "    data = np.loadtxt(filename)\n",
    "    return data\n",
    "\n",
    "# finds the points defined in the given array of data that most closely bound the given interpolation point\n",
    "# Returns an array of two points\n",
    "def lin_interp_bounds(data, point):\n",
    "    index_x1 = 0\n",
    "    while data[index_x1 + 1][0] < point:\n",
    "        index_x1 += 1\n",
    "    \n",
    "    p1 = [data[index_x1][0], data[index_x1][1]]\n",
    "    p2 = [data[index_x1+1][0], data[index_x1+1][1]]\n",
    "    \n",
    "    return p1, p2\n",
    "\n",
    "# Finds the linear interpolation of a given point based on a given array of defined points\n",
    "# Returns the linear interpolation value\n",
    "def lin_interp_point(data, point):\n",
    "    \n",
    "    [p1, p2] = lin_interp_bounds(data, point)\n",
    "    x1 = p1[0]\n",
    "    y1 = p1[1]\n",
    "    x2 = p2[0]\n",
    "    y2 = p2[1]\n",
    "    \n",
    "    interpolation = y1 + (point - x1) * (y2 - y1)/(x2 - x1)\n",
    "    \n",
    "    return interpolation\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "4.5625"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Problem 3b\n",
    "\n",
    "file = open('hw1.dat') # location of file containing needed data\n",
    "data = read_dat(file) # array of values defined in the given file\n",
    "\n",
    "point = 4.75 # point to interpolate\n",
    "\n",
    "interpolation = lin_interp_point(data, point)\n",
    "interpolation # answer: linear interpolation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2.9516602099609388"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Problem 3c\n",
    "\n",
    "# given an array of defined points and a point to interpolate, construct the highest order polynomial possible\n",
    "# that fits the given data.\n",
    "# Return the interpolation of the given point using the constructed polynomial\n",
    "def neville(data, point):\n",
    "    \n",
    "    n = len(data)\n",
    "    p = n*[0]\n",
    "    \n",
    "    for k in range(n):\n",
    "        \n",
    "        for j in range(n-k):\n",
    "            \n",
    "            if k == 0:\n",
    "                p[j] = data[j][1]\n",
    "                \n",
    "            else:\n",
    "                p[j] = ((data[j+k][0] - point)*p[j] + (point - data[j][0])*p[j+1]) / (data[j][0] - data[j+k][0])\n",
    "          \n",
    "    return p[0]\n",
    "\n",
    "estimate = neville(data, point)\n",
    "estimate # answer: approximation given by Neville's algorithm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Problem 3d (setup)\n",
    "\n",
    "# Defines the actual function for the dataset\n",
    "# returns the the evaluation of the function at a given point\n",
    "# NOTE: MUST BE DEFINED BY THE USER\n",
    "def func_actual(x):\n",
    "    \n",
    "    result = 100 / x**2\n",
    "    return result\n",
    "\n",
    "# Defines a linear function determined by the two given point arrays\n",
    "# Returns the evaluation of the function at the given `x`\n",
    "def func_lin(p1, p2, x):\n",
    "    \n",
    "    x1 = p1[0]\n",
    "    x2 = p2[0]\n",
    "    y1 = p1[1]\n",
    "    y2 = p2[1]\n",
    "    \n",
    "    slope = (y2 - y1) / (x2 - x1)\n",
    "    \n",
    "    result = slope*(x - x1) + y1\n",
    "    return result"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the case of the given data, at **```x = 4.75```**, linear interpolation produces a value closer to the value given by the actual function than interpolation with Neville's algorithm.\n",
    "\n",
    "The polynomial given with Neville interpolation does not have enough data to closely match the actual function for values outside of the given data, as you can see when the polynomial increases rapidly around **```x = 5```**. Even though linear intepolation is generally less accurate, in this case the actual function lends itsel well to being approximated over small intervals with linear equations becuase the derivative is linear.\n",
    "\n",
    "The plot below shows the actual function, Neville interpolation, and linear interpolation about the point of interest, allowing you to compare wich interpolation is more accurate.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAfoAAAD/CAYAAAAZrxROAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XdclWX/B/DPhRsVGaKCIA7C3CQOwAWllolpw1yRo2GZ\nT5YNbahYT6ZNm6ZPWo4n87FfmpqauXJw3Io7J0txkAqKyjrf3x8XHAFBGQfO4PN+vc7rrPvc93V7\ny/mea30vJSIgIiIi++Rg6QIQERFR6WGgJyIismMM9ERERHaMgZ6IiMiOMdATERHZMQZ6IiIiO3bX\nQK+Umq2UOq+U2p/jtUlKqXil1J6s20OlW0wiIiIqjsLU6H8A8GA+r38mIm2zbqvNXC4iIiIyg7sG\nehHZAuByPm8p8xeHiIiIzKkkffSjlVL7lFLfK6Vqma1EREREZDaqMClwlVI+AJaLSOus5+4AEkVE\nlFL/BuAhIs8U8Fnm2CUiKgYRYcsplVjF4nxIRC7mePofAMvvtP2kSZNMj0NCQhASElKcw1qliIgI\nREREWLoYpcaez8+ezw0oH+e3Z08EBg8GBg60dGlKbuPGjdi4caPp+eTJky1XGLIrhQ30Cjn65JVS\n9UTkXNbTxwAcvNOH7fnLhogsIzMT2LQJ+P57S5fEPPJWghjoyVzuGuiVUj8BCAHgppSKBTAJQKhS\nyh+AEUA0gJGlWEYiotvExwONGwN16li6JETW7a6BXkQG5/PyD6VQFptkT90Q+bHn87PncwPs//zS\n00PwYH4Tf4kol0INxivRAZQSrnlPRObWvj3w8ceAvf6eUUpxMB6ZBQM9EdmcxESgSRPg4kWgcmVL\nl6Z0MNCTuTDXPRHZnD//1DV5ew3yRObEQE9ENuePP8D+eaJCYtM9EdkUEaB+fWDzZt18b6/YdE/m\nwho9EdmUAwcAR0f7DvJE5sRAT0Q2hc32REXDQE9ENoWBnqho2EdPRDYjJQWoVw84exaoWdPSpSld\n7KMnc2GNnohsxsaNQECA/Qd5InNioCcim/H770Dv3pYuBZFtKdYytUREZU0EWLFC99ETUeGxRk9E\nNuHAAaBSJeDeey1dEiLbwkBPRDYhu9lecXgaUZEw0BORTVixAggLs3QpiGwPp9cRkdXLXq3uwgWg\nShVLl6ZscHodmQtr9ERk9VavBu6/v/wEeSJzYqAnIqu3YgWn1REVF5vuiciqpacDdesCBw8Cnp6W\nLk3ZYdM9mQtr9ERk1SIjgcaNy1eQJzInBnoismpsticqGQZ6IrJaIsDSpcAjj1i6JES2i4GeiKzW\n4cNAWhrQtq2lS0JkuxjoichqLVkC9OvHbHhEJcFAT0RWKzvQE1HxcXodEVmlmBi99vy5c0DFcrjO\nJqfXkbmwRk9EVum334A+fcpnkCcyJwZ6IrJKS5YAjz5q6VIQ2T423ROR1clexObcOaBaNUuXxjLY\ndE/mwho9EVmd5cuB7t3Lb5AnMicGeiKyOkuXstmeyFzYdE9EVuXaNZ3XPiYGcHGxdGksh033ZC6s\n0RORVVm+HOjcuXwHeSJzYqAnIqvyv/8BTz5p6VIQ2Q823ROR1UhOBry9dbO9s7OlS2NZbLonc2GN\nnoisxrJlQNeuDPJE5sRAT0RWY9EiYMAAS5eCyL6w6Z6IrMKVK4CPDxAXBzg5Wbo0lsemezIX1uiJ\nyCr89htw//0M8kTmViaBfubrx3HpUlkciYhs1aJFHG1PVBrKpOk+qUptvFnhMzgMDce//gU0a1aq\nhyQiG3PpEtCoEXDmDFCjhqVLYx3YdE/mUiY1eqfta/G15xQ8vzkcYV2T8eCDwIoVgNFYFkcnImv3\n669Ajx4M8kSloWz66Nu0QcV9u+Af7Ijj1dvglfZbEREB3HMP8MknYLM+UTm3YAEwZIilS0Fkn8pu\nMF716sDMmXD4Yjp6ff84dvaaiIXz0rF/v16OcsQIYOfOMisNEVmJmBjg4EHg4YctXRIi+1T2o+77\n9gX27oXatRMdXgnGvLeO4O+/AT8/PRCnXTtg9mwgJaXMS0ZEFvDTT8ATTwBVqli6JET2yTLT6zw8\ngJUrgWeeAbp0QZ2FX2D8m0acPAm8/76eZtOgATB6NHDggEVKSERlQASYPx8ID7d0SYjsl+UT5pw4\nAQwdClSurKvyjRsD0Ekzvv9ev+TjAzz/vK7xV6tWqsUlojK0Zw/Qv7/+GlAcX54LR92Tudy1Rq+U\nmq2UOq+U2p/jNRel1Bql1N9KqT+UUrWKXQJfX2DTJqB3b6BDB+CbbwCjEd7ewOTJQHQ08Oabeo6t\nlxfw8su6P4+IbN/8+cBTTzHIE5Wmu9bolVKdAVwDME9EWme9Ng3APyLykVJqHAAXERlfwOcLnwL3\n6FFg+HCgalVdnW/SJNfbMTH65TlzctfyHR0Lt3sish4ZGfrH++bNegYO5cYaPZnLXWv0IrIFwOU8\nL/cFMDfr8VwA/cxSmnvvBbZs0bX7jh2BTz8FMjNNb/v46D78mBjgrbeAX37RS1qOHg1ERZmlBERU\nRtauBRo2ZJAnKm3FHYxXR0TOA4CInANQx2wlqlABeP11YNs2nVUnKAjYvz/XJhUrAn366Lf37QPc\n3YGwMP3bgCP2iWwDB+ERlQ1zjbo3/4g+X19g/XrgueeABx7QVfgbN27bzNsbmDRJ9+VPnKhH7Ht7\nA6NGsZZPZK2uXAF+/51L0hKVhYrF/Nx5pVRdETmvlKoH4MKdNo6IiDA9DgkJQUhISOGOopQO9H36\nAK+8ArRsCXz3nc6VmUeFCrrFv3dvID5e1+zDwoD69YGRI/UXCvvyiazDTz8BDz0E1K5t6ZJYj40b\nN2Ljxo2WLgbZoUJNr1NKNQSwXERaZT2fBuCSiEwz62C8u1m5EnjpJd1G/9lngKfnHTfPyABWrQJm\nzgQMBj2694UXuKgOkaW1bQt89BHQvbulS2K9OBiPzKUw0+t+AhAJwE8pFauUGg5gKoAeSqm/ATyQ\n9bz0PfwwcOiQHo3fujXwxRc6mhcgZ1/+nj1AzZp6veuQEOB//wPS0sqk1ESUw549wOXL+m+RiEqf\n5RPmFNeRI7p2/88/wNdfA126FOpjaWnA0qXAt98Cf/+tewaef15P8yGi0vfSS0C9esCECZYuiXVj\njZ7MxXYDPaDzZy5eDLz2GtCtm24LvEtzfk6HDumAv3ChHu83ejTQtSuTdxCVluvX9WDZffv0PRWM\ngZ7MxTK57s1FKZ0x58gR/a3RujUwdSpw82ahPt6ihU7EFx2tm/NffFHv4j//0V9IRGRe//d/eogN\ngzxR2bHtQJ+tRg3gww/13Ptt23QEX7JE1/gLwclJNyceOgR8/jmwfLlOzjNunM65T0Tm8f33wLPP\nWroUROWLbTfdF+TPP4FXX9WZdD79VA/xLaJTp4CvvgLmzdOz+V55BQgMLIWyEpUTx47poTRxcXoN\nK7ozNt2TudhHjT6vHj10J+DAgXpi/fDhwJkzRdpF48a6dn/6tE7ON3iwvv/ll1xZeYmokL79Vq9M\nzSBPVLbss0afU3Ky7refOVNPoh83TrfVF1Fmph6t/9lnQEKCruE/8wxQvXoplJnIzly7prvD9u4F\nGjSwdGlsA2v0ZC72WaPPyckJmDJF58M9exbw89Nt8qmpRdpNhQrA448DW7fqrF6bNukFOSZOBC7c\nMS8gEf33v3pGC4M8Udmz/0CfzcsL+OEHYM0aYPVqvVLe/PnFaocPDNRN+JGROsjfe68ezBcdbf5i\nE9k6ET27ZfRoS5eEqHwqP4E+W+vWejWNuXOBGTMAf3+9Ek4xuhfuuUen3j9yRDccBAToNLuHDpVC\nuYls1ObNQHo6M+ERWUr5C/TZunbV7fAffKDb3wMD9Wj9YgT8unX17L5Tp/S6Ow88ADz2GLB7dymU\nm8jGfPONXk2SiaiILMP+B+MVhtGoM+xNnKij9nvv6Qw6xXT9uk668/HHQKtWerdBQeYrLpGtOHtW\n//iNji7WGNhyjYPxyFzKb40+JwcHvY7toUM6m8ezzwKhocBffxVrd46OwJgxwMmTQL9+epZfz566\nAYGoPJk5U///Z5AnspwyqdE/veRpBHkFIdg7GC3cW6CCQ4VSPWaJZWQACxYA//63HsQ3caIO/MVs\ne0xL00MCpkzRg/4nT2byHbJ/16/rmSmbNwNNm1q6NLaHNXoylzIJ9LN2zUJkfCQi4yJx7to5dKjf\nwRT4A70C4VzVuVTLUGwZGXrFm3//W2fZe+cd4KGHShTwf/xR765lS91D0K6deYtMZC2++w5YtUqP\ndaWiY6AncynzPvrE64nYFr8NkXGRMMQbsOvsLjSo1cAU+IO9g+Hn5gcHZUW9CpmZegH7KVN0Wq+3\n3wYefVQ3+RdDaqrO+T1lil7g4/33dXp+InuRmamnnc6ZU+gVpCkPBnoyF4sPxsswZmD/+f0wxBlM\ntf7k1GQEegUi2CsYQd5B6FC/A2pUrlGq5SwUoxFYsUKP1E9KAt58ExgyBKhSpVi7u35dpwX96CPg\nwQd1k37jxmYuM5EFLFmiE1Ju28bR9sXFQE/mYvFAn5+EqwkwxBsQGacDf9T5KDR1a5qr1t/QuSGU\npb5BRIANG4Bp0/QAvldeAZ5/vtgjjpKTdWrdr77SOfXffVcP/ieyVZ066T+L/v0tXRLbxUBP5mKV\ngT6v1IxU7EnYkyv4G8WIYO9gU/AP8AxA1YpVzVTqItizR1fJ//xTJ79/+WU9gK8YLl7Uzfnz5uks\nYq+/DtSsaebyEpWyyEidOOr4cZ06moqHgZ7MxSYCfV4igtikWFM/vyHegMMXD6NVnVa5av31neqb\n9bh3FB0NTJ+uo3Tv3nqZ3GIsj5u9qwkTgLVr9f1zzwGVKpm1tESl5rHHdBY8prwtGQZ6MhebDPT5\nSUlLwa6zu3LV+h0rOeaq9fvX80elCqUcMS9f1tlyvvoKaNJEt1/26VOsqs3evXoYQGysbjR45BH2\nd5J1O3xY55o6fZorO5YUAz2Zi90E+rxEBCcunTAFfUO8Aacun0KAZ4Ap8Ad5BcG9unvpFCA9Xa98\nM326bpMfPRoYMQJwLvpUwtWrdTN+7drAJ59wSh5ZryFD9NTRt96ydElsX0GBvlq1audu3rzJUTyU\nS9WqVc/fuHGjXn7v2W2gz0/SzSRsP7PdFPi3x2+He3X3XIG/ZZ2W5k/os3078MUXOmIPHKiXuivi\nfLqMDL343qRJQI8eui+/fhn2TBDdzd9/A50764yQzIRXcgUFemv6TiXrcacWoHIV6PMyihGHLx6G\nIU7382+N24qEqwnoUL+DqZ/frAl9EhJ0TtCZM2+tbdu3b5E64K9e1QvozJyp0+y+/rpOuUtkaUOH\nAr6+elwJlRwDPRUFA30RZCf0yZ7Xn53QJ9hLB/4g7yA0dWtasql9aWnAr7/qSfQnTujc+s89B3h7\nF3oX0dG6/377dr14Tv/+7L8nyzlxQqd1PnGiWL1TlA8GeioKBvoSyE7oszV2q2mEf86EPsHewWhf\nv33xE/ocPAjMmKFT7XbuDIwcqdPsFnLw3qZNumZfs6buHbjvvuIVg6gknnlGzyqdPNnSJbEfDPRU\nFAz0Znb26tlcaXyjzkXBz80PQV5BCPLW/f2NnBsVrdZ/7RqwaJFukz9/Xg/cGz4caNDgrh/NzNSp\nRidM0Jl5//1vwM2tBCdIVASnT+sBosePA66uli6N/WCgz+3DDz/E6dOnMWvWLMTExKBRo0bIyMiA\ng4MDQkNDER4ejhEjRhRpn3FxcWjRogWSkpIsl4DNTBjoS1lqRir2nttrGuGfndAnyDvIVOsvUkKf\nfft0MvyFC4EOHXR16ZFHdJ79O7h8WQ/W+/lnICJCNw4wYQmVtqFDAR8fvUgTmY8tBvqGDRvixo0b\niI6ORrVq1QAAs2fPxoIFC7BhwwazHScmJgaNGzdGenp6iQJ9UY+Z88eFtblToLe+0tqgKhWrINAr\nEGODxuKXJ3/BmbFnsP3Z7RjQYgBik2IxZvUYuH3khsDvAzH2j7FYfGgxziSfKXiH/v7A118D8fE6\nJ+433+gh9q+8AkRFFfgxFxfgyy+Bdev0Gjzt2wMGQymcMFGWAwduTf8kUkrBaDRi+vTpt71u60Qk\nO5gW6/OZmZlmLlHhMdCXAqUUfJx9MLDlQHz18FfY9fwuXHj9AqZ2nwp3R3fM3z8f/jP94TPdBwN/\nGYgvt3+JXWd3IT0zPfeOqlUDwsN1Xv1t23RHfJ8+uiP+iy/0/Px8tGqlP/L668ATT+hegAI2JSqR\nt9/Wc+Y5nY6yvfHGG/j000+RnJyc7/tHjx5Fz5494ebmhmbNmmHx4sUAgB07dsDDwyNXIF2yZAn8\n/f0BAJMnT0Z4eHihyjBnzhw0b94cbm5u6NWrF2JjY/PdLiYmBg4ODjAajQCA0NBQTJw4EZ07d4aT\nkxMeeughXLp0CQDQrVs3AICzszOcnJywffv2ux7LwcEB3377Lfz8/ODn52d67auvvkKTJk1Qp04d\nvPnmm4U6p5JgoC8j1StXR0jDELzV5S0sG7QMF16/gLXha/HwPQ/j8MXDGPHbCLhMc0HXH7pi/Nrx\n+O3ob7iQcuHWDpo00evZRkfrrDm7dwP33KOb9BcvBm7ezHU8pXRjwJEjehR0ixZ6fXAL/qgkO7Nl\ni67Rv/iipUtC1qRdu3YICQnBxx9/fNt7169fR8+ePfHUU08hMTERP//8M0aNGoWjR4+iQ4cOqFGj\nBtavX2/afuHChRgyZIjpeWFaBn777TdMnToVS5cuxcWLF9GlSxcMGjSowO3z7nPhwoWYO3cuLl68\niNTUVHzyyScAgE2bNgEAkpOTkZycjI4dOxbqWL/99ht27NiBw4cPm15bunQp9uzZgz179uC3337D\nnDlz7npeJSIipXrTh6DCuHLjiqw5sUYiNkTIg/MflFof1hLfL30l/NdwmbFzhkSdi5KMzIxbH0hO\nFvnxR5H77xdxdRV59lmRDRtEMjNv2/f+/SKdOom0ayeyc2fZnRPZJ6NRJDhYZO5cS5fEfmV9dxb5\nO1Uvr1nyW3E0bNhQ1q1bJwcPHhRnZ2dJTEyU77//XkJDQ0VEZNGiRdK1a9dcnxk5cqS89957IiLy\n7rvvyogRI0REJDk5WapXry5xcXEiIhIRESHh4eEiIhIdHS0ODg6SmfVdFxISIrNnzxYRkV69esmc\nOXNM+8/MzBRHR0eJjY29rbz57eeDDz4wvf/tt99Kr1698t22MMdSSsnGjRtzHVMpJWvWrMl1jO7d\nu9/pn7VQCvr/IiKs0VuTWlVroUeTHpgUMgmrn1qNS+MuYcmAJejcoDO2n9mO/ov7w/UjV/SY3wOT\nNkzC6vNbcWVAX90pv2+fruGPGaNHRr35pn4tqxmsVSs9Fe+ll4CwMOBf/wKSkix8wmSzli/Xyyvn\nqGyRlTBXqC+JFi1aICwsDB9++GGu12NiYrBt2za4urrC1dUVLi4u+Omnn3Du3DkAwODBg7FkyRKk\np6fj119/RUBAALyKuBpoTEwMxowZYzqGm5sblFI4c+YO46JyqFfvVhZZR0dHXLt2rUTHyq/8OV/z\n8fHB2bNnC1W24qpYqnunEnFQDmhZpyVa1mmJ5wOeB5A7oc+0rdNMCX2CvIIQ3CMYwSMWwS8hDQ4L\nf9Zz7apW1Wl3Bw2Cg58fhg3Trf3jxwPNmwOffgoMGMBkO1R4aWnAuHG6B4mzOqggERERaNu2LV57\n7TXTa97e3ggJCcEff/yR72eaNWsGHx8frFy5EgsXLsTgwYOLfFxvb2+8++67d2yuL478ug0aNGhw\n12Pl97m4uDg0a9YMABAbGwtPT0/zFTQfrNHbmNqOtRHmF4YPHvgAG4ZuwOVxlzH/0fm4r959WHtq\nLXr9txdqrwxB75ZR+ODHZ7Drg5eQdvE80K2bHsT34YdwvXwSs2bprv0pU4AHH9T5yYkK4+uvgYYN\ngYcftnRJyJo1adIEAwYMwJdffml6LSwsDMeOHcOCBQuQkZGB9PR07Nq1C0ePHjVtM3jwYHzxxRfY\nvHkz+vfvX+D+pYBmhxdeeAFTpkwx9YknJSXhl19+KfJ+8nJ3d4eDgwNO5viyHDlyZJGOle3jjz/G\nlStXEBcXhy+++AIDBw4sVBmKi4HexlV0qIi2Hm3xUoeX8NPjP+H0mNM4NOoQRviPwOXUKxhzZSFc\nPOciYGJdfD7AB8ei1iMjqCMkIADBf32I3YtOoEcPoGNH4IMPdG2NqCDnz+sfh9OnsxWIbpe39jpx\n4kRcv37d9HqNGjWwZs0a/Pzzz/D09ISnpyfGjx+PtBxfPAMHDsSmTZvwwAMPwPUOGZhyHivn4379\n+mH8+PEYOHAgnJ2d0bp1a6xevbrI+8mrWrVqeOedd9CpUye4urpix44ddz1WQfvr27cvAgIC0LZt\nW/Tp08c0/3/Lli1wKoUpLEyYUw7kTeizLWYrOpxKxXOnnNFldyJUnTq40XMwJu0bgA3nm2PWfxQ6\nd7Z0qckaPfOMzn6Xz4BqMjNbTJhDd+fg4IATJ06gcePGZt0vM+NRLiKC2KRYU9C/vmkt/LeexONH\nHZBesTp+rdwVie2ewav/6Q3X2mz0IW3HDqBfP+DoUc6bLwsM9PaJgZ4sJiUtBbvO7MSptYtRbfka\ntDGchtNNYHurhlD9u8O7bzja+HRApQqFX1KX7IfRCAQHAy+8AAwbZunSlA8M9PapQoUKOH78OAM9\nWZ6IYPW8P7Bn6g/oefkvNE26iLW+CoeCfCEPPYg2ze9HkHcQ6lSvY+miUhn47jtg3jydJMcK03zb\nJQZ6KgoGeiq29HTgs8+AH6edx7QHFqPj5Z9Qa9seHPdyxOJGN7Dzvjpwb9cVwQ06Idg7GC3cW6CC\nA+dc2ZMzZ/TyCxs36gyLVDYY6KkoGOipxE6e1M22iYnA91/dQEDyBsiKFUhfvhSpxjTs8HfHzw2S\nsdzjKlo16qjn9XsHI9ArEM5VnS1dfComEd0vf999ekVEKjsM9FQUDPRkFiLAggXAG28Agwbp1Ps1\nqgtw8CCwahWwciWMe3Yj0d8PO1q5YlGDZCzNPIQGzj4I9grWy/Z6B8PPzQ8Oiu2/tmDxYh3g9+wB\nqlSxdGnKFwZ6KgoGejKrxETgtdeAv/4Cvv02T+KUK1d0St5Vq4DVqyGVKyOxS1vsaOmCJR5JWPfP\nTiTdTNJB3ysYwd7BaF+/PWpUrmGx86H8XboEtGwJ/N//AUFBli5N+cNAT0XBQE+lYt06YORIoF07\nnUAlR4poTQQ4dAj44w99MxgAf39c7RqI3S1dsdL5IrYmbMe+c/vQ1K0pgr2DTU3+DZ0b2sUa1rZs\nyBDAzQ3IkdiMyhADPRUFAz2Vmhs3gPfeA2bP1hnTRoy4w6js69f1sO01a4A//wRiY4Fu3ZAeGoKD\nretiQ5Wz2BofCUOcAUYxItg72BT8AzwDULVi1TI9t/Lsp59018zu3YCjo6VLUz7ZYqBv1KgRZs+e\njfvvvx8ffvghTp8+jVmzZlm6WHj44YcxaNCgQq9nb4sY6KnU7d8PPP88ULkyMGsWcO+9hfjQ+fPA\n+vXA2rX6lpEB3H8/JDQUZ9o1xWYVC0O8AZFxkTiSeASt6rQy1fiDvYNR36l+qZ9XeRQdDXTooBth\n7rvP0qUpv2w90Nuq0NBQhIeHm9LS2oo7BXquXkdm0bo1sHUrMGMG0KULMGoU8NZbevG8AtWtq0f1\nDRqkm/lPnQLWr4daswZeb72FQY6OGBQaCoSOwY2BHbHDIQGGeAPm75+PUStHoVrFarma+/3r+TOh\nTwllZgJPPw28/jqDPNkmo9EIBwsme7D08fNTotIopaKVUlFKqb1KqR3mKhTZpgoVgNGjgX37gAMH\ngDZtgA0bCvlhpYAmTYDnngN+/hk4dw5YsQJo2xZYuhTVOnZCt/uHY/ycY1iW9jguPL4d68LXopdv\nLxxJPIJnlj0Dl2ku6PZjN4xfOx7L/l6GiykXS/V87dG0afo65lhZlKhYJk+ebGoqj4mJgYODA+bN\nmwcfHx/UqVMHU6ZMMW0rIpg6dSp8fX3h7u6OgQMH4vLly6b3n3zySXh4eMDFxQUhISGm1eIAYPjw\n4Rg1ahR69+6NmjVrYuPGjbeVJTQ0FHPmzAEAzJ07F126dMEbb7wBV1dXNGnSxLRs7rvvvovNmzdj\n9OjRcHJywssvvwwAOHr0KHr27Ak3Nzc0a9YMixcvvuPxhw8fjhdffBE9e/aEk5MTQkNDERsba75/\n3KISkWLfAJwC4HKXbYTKp6VLRby9RYYOFbl4sYQ7y8wUOXhQ5OuvRZ58UqRuXb3zwYNFZswQOXhQ\nrqRckjUn1kjEhgjpOb+n1Pqwlvh+6Svhv4bLjJ0zJOpclGRkZpjj1OzS+vX6nzU21tIlIRGRrO9O\nm/pObdiwoaxbt05ERCIiIiQ8PFxERKKjo0UpJc8//7ykpqZKVFSUVKlSRY4ePSoiItOnT5egoCA5\ne/aspKWlyQsvvCCDBg0y7feHH36QlJQUSUtLk1dffVX8/f1N7w0bNkycnZ3FYDCIiEhqaupt5QoJ\nCZHZs2eLiMiPP/4olStXltmzZ4vRaJQZM2aIp6dnvtuKiKSkpIi3t7fMnTtXjEaj7Nu3T2rXri1H\njhzJ9/g3b96UYcOGiZOTk2zZskXS0tJkzJgx0rlz55L/A99BQf9fRKTETfcKXOqWCtC3L3D//cCk\nSTqj2tSpOk96sQbTOzjonbRoAbz0km7qP34c2LRJD/D7+GPUSkpCj+Bg9OjcGeg8CZmP++PI1VMw\nxBmwNW4rPt/2Oc5dO4cO9TswoU8ecXHA4ME6T4K3t6VLQyWhJptntopMMu84AKUUIiIiULlyZbRu\n3Rpt2rRBVFQUmjZtipkzZ+Kbb76Bh4cHAL28rY+PDxYsWAAHBwcMy7HAwsSJEzF9+nRcvXoVNWvW\nBKCXfQ0MDAQAVK5c+a5l8fHxMfXBDx06FKNGjcKFCxdQp87tKb1XrFiBRo0a4emnnwYAtGnTBo8/\n/jgWL16MCRMm3Hb8KlkJJ3r37o1OnToBAD744APUqlULZ86cQf36ZT+2qKSBXgD8qZTKBDBLRP5j\nhjKRHaln6H8sAAAbUklEQVRZU6fQfeopPRXvhx90P36JU6kqBfj56duzz+rXzp7VAwW2bAH+9S9U\nOHoULf390bJTJzwX3A/oPw2JNStgW/w2GOIMmLZ1Gnad3YUGtRqY5vQHeQehqVvTcjW1LzUVeOIJ\nYMwYoHt3S5eGSsrcAdqc6tata3rs6OiIa9euAdBN+48++qipb1tEUKlSJZw/fx5169bF22+/jV9+\n+QWJiYlQSkEphcTERFOg9y7ir9N6OeYCV6tWDQBw7dq1fAN9TEwMtm3bBldXV1PZMjMzTYG/oOPn\nfK169epwdXXF2bNnbTLQdxKRBKWUO3TAPyIiW8xRMLIvbdsC27bpxVFCQnRsnjDBzFO3PD2B/v31\nDQCuXdNrq27dCsycCQwfjtqurggLCkJYUBAQ+AkyBjbH/ktHYIgzYO3ptXh/0/tISk1CoFegKZtf\nh/od7Dqhz5gxQP36wLhxli4JlVcNGjTAnDlzEJRPZqYFCxZg+fLlWL9+PRo0aICkpCS4uLhkd2MA\ngFl/mOfdl7e3N0JCQkz9+IX5DADExcWZHl+7dg2XLl2Cp6en2cpZFCUK9CKSkHV/USm1BEAHALcF\n+ogcSbJDQkIQEhJSksOSjapQQbe6P/aYHuzVvLlOtNO3bzGb8++mRg3dd5A91cdo1IupGwymXx0V\nT59GW39/tO3YES917Av0m4IEl0rYdmY7IuMiMWHDBOw7tw9+bn65av2NnBvZRa1/xgy9WM2OHaV0\nDajQNm7cmO9AMnuRMzDnNXLkSLz99tuYO3cuGjRogIsXL8JgMOCRRx7B1atXUaVKFbi4uCAlJQVv\nvfVWqf7t1a1bF6dOnTI9DwsLw1tvvYUFCxZg4MCBEBFERUWhZs2aaNq0aYH7WblyJSIjI9GuXTtM\nmDABQUFBFqnNAyUI9EopRwAOInJNKVUdQE8Ak/PbNoKrYVAOHh46Icu6dXqU/n/+A3z1FWDm5Zlv\n5+Cgf100bw4884x+LSkJ2LUL2L5dd1D/61/wUAqPduiARzt0ANpPRGpYa+xNi8HW2K1YcnQJ3vjz\nDVNCn+y+fltM6LNsmU6Ks3kz4ORk6dJQ3krQ5Mn5fp1atTsF4Lzv5Xw+ZswYAEDPnj2RkJCAOnXq\nYMCAAXjkkUfw9NNP448//kD9+vXh5uaG999/HzNnzjRbufIry9ChQzFjxgyEh4dj+vTpWLNmDV59\n9VWMHTsWIoI2bdrgs88+u+M+Bw8ejIiICBgMBgQEBGDBggWm91q2bIl33nkHgwYNKtJ5FFexE+Yo\npRoBWALdT18RwH9FZGo+20lxj0H2Ly1N9+F//LEO+uPGWTgTm4gembZjh77t3KnTw7m7A+3bA+3a\nQQICENukNiKTDiIyLhKGeMNtCX2CvIPg5eRlwRO5s5079RoFv/+uk+OQ9bHFhDmkDR8+HN7e3njv\nvffK7JjMjEdWLzZWJ2nZsQP4/HO9NKrVNCUbjcDff98K+rt2AVFRekxAQAAQEIAbrZtjdz0jtiTf\nCv7ZCX2ya/7WktDn1Cmgc2c9XuKRRyxdGioIA73tYqAnuoP164GXX9bN+59/rldPs0oZGbq/f/du\nfduzRwd/d3egbVuIvz/O+tbFFrcUrL95BIZ4A05dPoW2Hm1NgT/IOwh1qt8+yrc0xcXpwZCvvw68\n+GKZHpqKiIHedo0YMQJeXl4M9EQFSU/Xtc333gMGDgQmTwayZrZYt8xMPbd/zx6dHnDvXn2rWBHw\n98fNls1w3KsaNjknYYU6BkPCTrhXd8+Vv7+FewtUcKhQKsWLj9dBftQoYOzYUjkEmREDPRUFAz3Z\npMREnWxn8WLg7bd1gCpELgzrIgKcOaMDf/YtKgo4exbSrBmu+DXAEc9K2OSchKWVTuKIXEB7z/am\nwG+uhD5nz+og/9xzwBtvlPy0qPQx0FNRMNCTTTt4UAenEyeAjz6ysv774rp2TS8IsH+/DvxRUcDB\ngzA6VkOib30c9aiEyFpJWFElBtd9fdC2cSdT8Pdz84ODKnxCypgYoGdPYPhwYPz4UjwnMisGeioK\nBnqyC2vW6L5lJye9+EpWdkn7IaJHJR44YLrJwQOQ48eR7O6E455VEFnrKvbXzkCV1m3h3e4BdGzc\n5Y4JfaKigN699b/bK6+U8flQiTDQU1Ew0JPdyMwE/vtfnVXP3x+YMsUM6XStXVqa7vs/dAg4eBA3\nonYj48A+VD1zAWfcKmGvaxou+rijYotWqNs+BC06Pwafek2xYYPCwIHA118DTz5p6ZOgomKgp6Jg\noCe7c/Mm8O23umbfs6fuy/f1tXSpylhqKnDsGNIPROHcjg24EbUL1U5Ew/1cMs46VcCRau6o1KYp\nGoUEoUHgg6jSsg3g4mLpUlMhMdBTUdwp0HPlObJJVavqkeMnTgBNmwKBgTp/fo7MlfavShWgVStU\nGvwUvKfPht+GKLgfT8IrQ1Mxsu4fyBg+CKnuqTj26ywcHPQAUjxr46pLdVxs3xwpw4YAn34KrFih\nWwsyMix9NkRm89dff91xoZvIyEj4+fnByckJy5YtK7NybdmyBc2aNSuz42VjjZ7swuXLOsPejBm6\nT/rtt/UPgPIkJgZ4/HGgUSNg9uzcaW1T0lKw68xOHIhag4SdG5B2+ACa/aPQNrkGGl1IQ81L16Aa\nNoRqeu+tVQHvuUffe3jYwehH22PLNfqQkBDs378f58+fR6VKd08SFRMTg0aNGiEjI8O0gl1J/PXX\nXwgPD0dsbGy+73fv3h39+vXD6NGjS3ysO3FwcMCJEyfQuNTze9+5Rl/S1euIrIKLi87b/tprwJdf\n6sxvoaF6tH779pYuXekSAebOBd58U6cQHjv29rhcvXJ1dGsUgm6NQoB+eoGR45eOwxBnwIy4SOw+\nvRXq5Ck8mFkRnS4nofkf++D5QzIqnYwGUlJ0v8g999y6z35crx5/BFAuMTEx2LJlC5ydnbFs2TI8\n/vjjd/2MiGQHqjIooS5j8+bNS/041rLwFZvuya44OwMTJ+om/KAgXcMNDdU53Y1GS5fO/M6cAcLC\ngC++AP78U//QKcx3i1IKfm5+GOo/FDP7zMSulw9i7bQEdH31c+x4KhQj+1WA+yN/w/edmhj5XRh+\nea0XYh9oB2O1qsBff+lfFG3aADVr6vvHHtO/qmbOBNauBU6fZndAOTVv3jwEBQVh2LBh+PHHH3O9\nd/PmTbz22mto2LAhnJ2d0bVrV9y8eRPdunUDADg7O8PJyQnbt2/H5MmTER4ebvpsTEwMHBwcYMz6\nQ/7xxx/RvHlzODk5wdfXF7NmzSpU+Xx9fXH69GmEhYXByckJaWlpaNSoEdavX2/aJuexs487b948\n+Pj4oE6dOpgyZYppW6PRiClTpsDX1xdOTk5o37494uPj0a1bN4gIWrduDScnJyxevPi2LoWjR48i\nNDQULi4uaNWqFZYvX256b/jw4Rg9erSpnEFBQTh9+nQhr0IeIlKqN30IIstISxOZP1/kvvtEfH1F\nPv9c5MoVS5eq5FJTRT77TMTdXSQiQp+nuWUaM+Xg+YMya9csGbZ0mPh95SdOHzpJ93ndZcL6CbLq\n+Cq5fC5aZPdukUWLRKZMEXnmGZGQEBFvb5HKlUWaNBHp0UNk5EiRadNEFi8W2bVL5NIl8xfYzmR9\nd9rcd6qvr6989913snv3bqlUqZJcuHDB9N6oUaMkNDRUEhISxGg0isFgkLS0NImOjhYHBwcxGo2m\nbSMiIiQ8PNz0PHubzMxMERFZuXKlnD59WkRENm3aJI6OjrJ3714REdm4caN4e3sXWMaGDRvK+vXr\ncz1ft25dvseOjo4WpZQ8//zzkpqaKlFRUVKlShU5evSoiIh89NFH0rp1azl+/LiIiOzfv18uZf3/\nVkrJqVOnTPvNWa709HTx9fWVqVOnSnp6uqxfv15q1qwpx44dExGRYcOGSe3atWXXrl2SmZkpQ4YM\nkUGDBhV4TgX9fxERNt2TfatUCXjqKWDIEL0M/Vdf6ZS6jz2mB+8FBtpWy7MIsHSpbqb389NryZdW\nC6SDckCLOi3Qok4LPBfwHAAg8XoitsVvQ2RcJKZtnYZdZ3ehQa0GCPYKRvBDwQj2fgx+bn66yTI1\nFYiO1s0r2TeDQdf2T53SywY3aqTXJ27USN8aNrx1q169dE6MSs2WLVsQGxuLJ598Ei4uLvD19cVP\nP/2EMWPGQETwww8/YMeOHahXrx4AIDAwMNfnJasJvzB69epletylSxf07NkTmzdvhr+/f6E+L0Xo\nJlBKISIiApUrV0br1q3Rpk0bREVFoWnTppg9ezY++eQT+GZN+2nVqlWhjmMwGJCSkoJx48YBAEJD\nQxEWFoaFCxdi4sSJAIBHH30UAQEBAIAhQ4bgtddeK3SZc2Kgp3JBKSA4WN8SEoB584ChQ3Ua+uHD\ndU79OwzStTijEfjtNz2dMCUF+OYbPa2wrNV2rI0wvzCE+YUBADKMGdh/fj8i4yKx9vRavLfpPSSn\nJiPQKxDBXnq53g4PdEGNHF/KAPQvlkuXbgX906d1kqAVK/TjmBjdLdCwIeDjo285HzdooPtp6Hbm\n+uVajP7yefPmoWfPnnDJmsY5aNAgzJ07F2PGjEFiYiJSU1PNNjBt1apVeO+993Ds2DEYjUbcuHED\nrVu3Nsu+81O3bl3TY0dHR1y7dg0AEBcXV6xzSkhIuG1mgI+PD86cOWN6nv2DKO8xi4qBnsodDw/d\nxfzmm8DmzcCCBcB99wHNmumA37cv4GUlS8knJwOLFumZcE5Outz9+gEVSmfdmyKr6FARbT3aoq1H\nW4zuoEcwJ1xNgCHeAEOcAe+ufxdR56Pg5+ana/3eOvg3cm4E5eYGuLkB7drdvmOjEbhwQbcIREfr\nwH/4MLBqlX4cE6NbBHIG/ux7b2997+Ghf8mVNxYakX/z5k3873//g9FohIeHBwAgLS0NV65cwYED\nB9CyZUtUrVoVJ0+evK3Wm18tvnr16rh+/brpeUJCgulxWloannjiCSxYsAB9+/aFg4MDHn300WIP\n5st7rHPnzhX6s97e3jh58mSRB/d5enoiLi4u12uxsbFoWgrThcrhXwGRphTQtau+ff21TrG7aJEe\nzNewIdCnD/DggzoOFWKGkNlkZAAbNuiR9CtW6MVoZszQ97bQzeBR0wOPNXsMjzV7DACQmpGKvef2\nIjIuEkuOLsEbf74BoxhNy/UGewcjwDMAVStWvbUTBwc9or9ePd2/kpeInlMZE6PTBsfG6se7dum1\neGNjgYsX9ee9vQu+ubvrY1GJLVmyBBUrVkRUVFSuKXX9+/fHvHnz8PHHH2P48OEYO3Ys5s2bh7p1\n62LHjh0ICAiAu7s7HBwccPLkSdxzzz0AAH9/f3z00UeIi4uDk5MTpk6datpnWloa0tLSULt2bTg4\nOGDVqlVYs2bNbT8gCsvf3x8///wzHnroIezbtw+//PJLrq6BO/2AePbZZzFhwgQ0a9YMvr6+OHDg\nALy8vODi4oJ69erh1KlT+db4O3bsCEdHR3z00UcYO3YstmzZghUrViAiIqJY53AnDPRE0KvihYXp\nW0YGsHUrsHy5XrP99GmdV79rVz1VLyDAvK3GIrrSun49sHq1HrR+zz1AeDgwfTpQu7b5jmUJVSpW\nQaBXIAK9AjE2aCxEBLFJsYiMi4Qh3oAxq8fgSOIRtKrTyrRwT5BXEOo71S94p0rptYtdXXVzTH7S\n0vS0hLi4W7e//wbWrdOP4+OBpCTA01M34WTf6tfPfV+vXtn+0rNR8+bNw4gRI1C/fu7rNnr0aIwZ\nMwbTpk3DJ598grfffhvt27dHSkoK2rRpgz/++APVqlXDO++8g06dOiEjIwOrV69G9+7dMWDAALRu\n3Rru7u4YN26caVR6jRo18OWXX6J///5IS0tDnz590Ldv30KXNW8Lwvvvv49BgwbB1dUV3bp1w5Ah\nQ3Dp0qUCt8/5fOzYsUhLS0PPnj3xzz//4N5778WSJUvg4uKCSZMm4emnn8bNmzcxa9YsuLu7mz5X\nqVIlLF++HC+++CKmTJkCLy8vzJ8/3/RDx5xT85gwh+guEhP1jLLNm4Hdu/US8x4euqk/O7dMgwZA\n3bpAnTq6Nbpy5VsVRRHgxg3dt37x4q3K56lTeun6PXt0pr+uXYGHHtKtCDm65sqF6+nXsfPMThji\nDaYfAI6VHBHkFWSq9fvX80elCmYOuDdv6h8D8fH6Fhenn2e/duaM7kJwc9OBP/vm6Xn7vaurWZtc\nbDlhDpU95ronMqPMTF0x/Ptv4NgxfYuL0/Hg/Hngn3+A9HTdj16pkq5YVq4M1KihY0HOsWX+/rqF\noLwF9rsREZy4dAKRcZGmwH/q8im09Wibq9bvXt397jsrqYwMfWGzfwCcPXvrcULCrcc3b+pfgB4e\nOvDnvM95c3MrVHcBAz0VBQM9URkT0fEhLU2npC+PY8LMLelmEnac2YGtcVthiDdge/x2uFd3z9XX\n38K9BSo4WGik4o0btwJ/QoK+nT1763H27epV3fxTr54O/NljEfI8V40bM9BToTHQE5HdyTRm4kji\nERjiDIiMj4QhzoCEawnoUL+DKfAHegXCuaqVTcNLTdUtBAkJwLlz+pb9I+D8edNrKiaGgZ4KjYGe\niMqFf67/Y0roExkfmSuhT5C3Dv5+bn5wUNY/0p5N91QUDPREVC7lTOiTPdDvtoQ+9TugRuUali7q\nbRjoqSgY6ImIsiRcTTAFfkO8AfvO7UNTt6am5n5TQh8LJy1goKeiYKAnIipAakYq9iTsyTW1L9OY\nqZv6s7L53ZbQpwww0FNRMNATERVSdkKfnIH/8MXDaFWnlanWH+wdfOeEPmZQ0Bd3tWrVzt28ebNu\nfp+h8qtq1arnb9y4ke9EXQZ6IqK7SElLwa6zu0zBPzIuEo6VHHPN6Td3Qp871dCIioKBnoioiPIm\n9ImMj8Tpy6cR4Blwq6+/hAl9GOjJXBjoiYjMIDuhT3Zz/7b4bajtWDtXrb9lnZaFTujDQE/mwkBP\nRFQKjGLEkYtHTJn8iprQh4GezIWBnoiojCReT8S2+G0wxBmwNW4rdifsNiX0yZ7a19StKZRSDPRk\nNgz0REQWkp6Zjv3n9+ca4Z+cmowgryD8PuR3BnoyCwZ6IiIrknA1AYZ4Ax5v/jgDPZkFAz0RkRVi\n0z2Zi/Wv7EBERETFxkBPRERkxxjoiYiI7BgDPRERkR1joCciIrJjDPRERER2jIGeiIjIjjHQExER\n2TEGeiIiIjvGQE9ERGTHGOiJiIjsGAM9ERGRHWOgJyIismMlCvRKqYeUUkeVUseUUuPMVSgiIiIy\nj2IvU6uUcgBwDMADAM4C2AlgoIgczbMdl6klIioiLlNL5lKSGn0HAMdFJEZE0gH8DKCveYpFRERE\n5lCSQF8fQFyO5/FZrxEREZGVqFgWB4mIiDA9DgkJQUhISFkclojIZmzcuBEbN260dDHIDpWkjz4Q\nQISIPJT1fDwAEZFpebZjHz0RURGxj57MpSRN9zsB+CqlfJRSlQEMBLDMPMUiIiIicyh2072IZCql\nRgNYA/2DYbaIHDFbyYiIiKjEit10X+gDsOmeiKjI2HRP5sLMeERERHaMgZ6IiMiOMdATERHZMQZ6\nIiIiO8ZAT0REZMcY6ImIiOwYAz0REZEdY6AnIiKyYwz0REREdoyBnoiIyI4x0BMREdkxBnoiIiI7\nxkBPRERkxxjoiYiI7BgDPRERkR1joCciIrJjDPRERER2jIGeiIjIjjHQExER2TEGeiIiIjvGQE9E\nRGTHGOiJiIjsGAM9ERGRHWOgJyIismMM9ERERHaMgZ6IiMiOMdATERHZMQZ6IiIiO8ZAT0REZMcY\n6ImIiOwYAz0REZEdY6AnIiKyYwz0REREdoyBnoiIyI4x0BMREdkxBnoiIiI7xkBPRERkxxjoiYiI\n7BgDPRERkR1joCciIrJjDPRERER2jIGeiIjIjjHQExER2TEGeiIiIjvGQE9ERGTHGOiJiIjsGAN9\nCW3cuNHSRShV9nx+9nxuAM+PiDQG+hKy9y8bez4/ez43gOdHRBoDPRERkR1joCciIrJjSkRK9wBK\nle4BiIjslIgoS5eBbF+pB3oiIiKyHDbdExER2TEGeiIiIjtmlkCvlKqilNqulNqrlDqglJpUwHZf\nKqWOK6X2KaX8zXHs0laYc1NKdVNKXVFK7cm6vWuJspaEUsohq+zLCnjf5q5dTnc6P1u/fkqpaKVU\nVNb/0R0FbGOz1+9u52fL108pVUsptVgpdUQpdUgp1TGfbWz22pF1qGiOnYhIqlIqVESuK6UqANiq\nlFolIqY/SqVULwBNROSerP/M3wEINMfxS1Nhzi3LJhF5xBJlNJMxAA4DcMr7hq1euzwKPL8stnz9\njABCRORyfm/awfW74/llsdXr9wWAlSLSXylVEYBjzjft4NqRFTBb072IXM96WAX6B0TeUX59AczL\n2nY7gFpKqbrmOn5pKsS5AYDNjo5VSnkBeBjA9wVsYrPXDijU+QE2fP2gy36nv2Wbvn64+/llb2NT\nlFJOALqIyA8AICIZIpKcZzNbv3ZkBcwW6LOaRvcCOAfgTxHZmWeT+gDicjw/k/Wa1SvEuQFAUFbT\n2u9KqeZlXMSS+hzAG8j/Bwxgw9cuy93OD7Dt6ycA/lRK7VRKPZfP+7Z+/e52foBtXr9GABKVUj9k\ndTnMUkpVy7ONrV87sgLmrNEbReQ+AF4AOtrQH9tdFeLcdgNoICL+AL4GsLSsy1hcSqneAM6LyD7o\nWpHN1YzupJDnZ7PXL0snEWkL3WrxklKqs6ULZGZ3Oz9bvX4VAbQF8E3W+V0HMN6yRSJ7ZPZR91lN\nTxsAPJTnrTMAvHM898p6zWYUdG4ici27eV9EVgGopJRytUARi6MTgEeUUqcALAQQqpSal2cbW752\ndz0/G79+EJGErPuLAJYA6JBnE1u+fnc9Pxu+fvEA4kRkV9bzX6ADf042fe3IOphr1H1tpVStrMfV\nAPQAcDTPZssAPJ21TSCAKyJy3hzHL02FObecfWZKqQ7QiYgulWlBi0lE3haRBiLSGMBAAOtF5Ok8\nm9nktQMKd362fP2UUo5KqRpZj6sD6AngYJ7NbPb6Feb8bPX6ZV2DOKWUX9ZLD0APGM3JZq8dWQ+z\njLoH4AFgrlLKAfrHwyIRWamUGglARGRW1vOHlVInAKQAGG6mY5e2u54bgCeUUi8CSAdwA8AAyxXX\nPOzk2hXIjq5fXQBLlE41XRHAf0VkjR1dv7ueH2z7+r0M4L9KqUoATgEYbkfXjqwEU+ASERHZMWbG\nIyIismMM9ERERHaMgZ6IiMiOMdATERHZMQZ6IiIiO8ZAT0REZMcY6ImIiOwYAz0REZEd+38JLpkY\nKuB/kAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x3fd2e80>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Problem 3d (plotting)\n",
    "\n",
    "[p1, p2] = lin_interp_bounds(data, point) # point arrays of the defined data points that bound the interpolation point\n",
    "y_center = (p1[1] + p2[1]) / 2 # average y-value of the bounding points\n",
    "\n",
    "x = np.linspace(p1[0]-1, p2[0]+1, 100) # domian to plot function over\n",
    "\n",
    "plt.xlim(p1[0]-1, p2[0]+1)\n",
    "plt.ylim(y_center - 10, y_center + 10)\n",
    "plt.plot(x, neville(data, x), label='Neville interp.')\n",
    "plt.plot(x, func_lin(p1, p2, x), label='linear interp.')\n",
    "plt.plot(x, func_actual(x), label='Actual function')\n",
    "plt.legend(loc='center left', bbox_to_anchor=(1, 0.5))\n",
    "plt.show()"
   ]
  }
 ],
 "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.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}