{
 "metadata": {
  "name": "trapezoid_rule"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Basic numerical integration: the trapezoid rule\n",
      "===============================================\n",
      "\n",
      "A simple illustration of the trapezoid rule for definite integration:\n",
      "\n",
      "$$\n",
      "\\int_{a}^{b} f(x)\\, dx \\approx \\frac{1}{2} \\sum_{k=1}^{N} \\left( x_{k} - x_{k-1} \\right) \\left( f(x_{k}) + f(x_{k-1}) \\right).\n",
      "$$\n",
      "<br>\n",
      "First, we define a simple function and sample it between 0 and 10 at 200 points"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%pylab inline"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "\n",
        "Welcome to pylab, a matplotlib-based Python environment [backend: module://IPython.zmq.pylab.backend_inline].\n",
        "For more information, type 'help(pylab)'.\n"
       ]
      }
     ],
     "prompt_number": 1
    },
    {
     "cell_type": "code",
     "collapsed": true,
     "input": [
      "def f(x):\n",
      "    return (x-3)*(x-5)*(x-7)+85\n",
      "\n",
      "x = linspace(0, 10, 200)\n",
      "y = f(x)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 2
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Choose a region to integrate over and take only a few points in that region"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": true,
     "input": [
      "a, b = 1, 9\n",
      "xint = x[logical_and(x>=a, x<=b)][::30]\n",
      "yint = y[logical_and(x>=a, x<=b)][::30]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 3
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Plot both the function and the area below it in the trapezoid approximation"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "plot(x, y, lw=2)\n",
      "axis([0, 10, 0, 140])\n",
      "fill_between(xint, 0, yint, facecolor='gray', alpha=0.4)\n",
      "text(0.5 * (a + b), 30,r\"$\\int_a^b f(x)dx$\", horizontalalignment='center', fontsize=20);"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD9CAYAAABHnDf0AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XlclOX+//HXAMMiiKLCoKCigrEI4p7+Tkl6sOykmXo0\ny+KrnXq0mubJrNOinQI8Vi6VnRYzTnVC2xSNSM0wy9SjaCmiKIuyu4zsOzO/P+4ACUQZB+4BPs/H\ng8fc3PfM3B9HfXNx3dd9XRqj0WhECCFEh2KldgFCCCHMT8JdCCE6IAl3IYTogCTchRCiA5JwF0KI\nDkjCXQghOqBmw33+/PnodDoCAwMbHXv99dexsrJCr9fX7YuIiMDHxwdfX1+2b99u/mqFEEJck2bD\nfd68ecTFxTXan5GRwY4dO+jfv3/dvuPHj7Nx40aOHz9OXFwcjz76KAaDwfwVCyGEuKpmw/2mm27C\nxcWl0f6nnnqKf/3rXw32bdmyhTlz5qDVavHy8sLb25sDBw6Yt1ohhBDXpMV97lu2bMHT05OgoKAG\n+7Ozs/H09Kz73tPTk6ysrOuvUAghRIvZtOTJpaWlhIeHs2PHjrp9zc1eoNFormmfEEKIq2vJbDEt\narmnpKSQnp7O0KFDGTBgAJmZmYwYMYK8vDw8PDzIyMioe25mZiYeHh5XLFC+jLz00kuq12ApX/JZ\nyGfRHj6LrCwjWq0RKysjp0+37blbqkXhHhgYSF5eHmlpaaSlpeHp6UlCQgI6nY6pU6cSHR1NZWUl\naWlpnDp1itGjR7e4ICGEsFRvvglVVTB9OgwapHY1zWs23OfMmcO4ceNITk6mb9++bNiwocHxy7tY\n/P39mTVrFv7+/kyePJl169ZJF4wQosMoKoJ//1vZ/vvf1a3lWmiMprT3r+eEGo1Jv2J0RPHx8YSE\nhKhdhkWQz6KefBb1LOmzWLMGFi6EP/0J9uxp+/O3NDsl3IUQ4iqqq8HbG86cgc2b4c47276Glman\nTD8ghBBX8eWXSrD7+MCUKWpXc20k3IUQohlGI9Tes/nUU2DVTlJTumWEEKIZcXEweTLodJCWBg4O\n6tQh3TJCCGFGr76qPC5erF6wm0Ja7kIIcQV79sDNN4OLi9Ln3rWrerVIy10IIcykttW+YIG6wW4K\nabkLIUQTDh2CkSPB0VFptffsqW490nIXQggzCA9XHh95RP1gN4W03IUQ4g+OHYOgILC1VUbI9O6t\ndkXSchdCiOu2fLkyvv3BBy0j2E0hLXchhLjMr79CcDDY2UFqKvTpo3ZFCmm5CyHEdVi2THl8+GHL\nCXZTSMtdCCF+l5AAI0YoNyulpoK7u9oV1ZOWuxBCmKi21f7II5YV7KaQlrsQQgAHD8KoUdCli9Jq\n1+nUrqghabkLIYQJnn1WeXzsMcsLdlNIy10I0ent3AmhodCtm9Jq79FD7Yoak5a7EEK0gMEAS5cq\n2888Y5nBbgppuQshOrVNm2D2bOVmpdOnlT53SyQtdyGEuEZVVfCPfyjbL71kucFuCgl3IUSn9cEH\nSmvdxwfmz1e7GvOSbhkhRKdUWKiE+rlzStfMX/+qdkXNM2u3zPz589HpdAQGBtbte/rpp/Hz82Po\n0KFMnz6dgoKCumMRERH4+Pjg6+vL9u3bTShfCCHaRni4Euxjx8LMmWpXY37Nhvu8efOIi4trsG/S\npEkkJiby66+/MnjwYCIiIgA4fvw4Gzdu5Pjx48TFxfHoo49iMBhar3IhhDBRWhqsWqVsr1oFGo26\n9bSGZsP9pptuwsXFpcG+0NBQrKyUl40ZM4bMzEwAtmzZwpw5c9BqtXh5eeHt7c2BAwdaqWwhhDDd\nkiVQWQlz58KYMWpX0zpsrufFH374IXPmzAEgOzubG2+8se6Yp6cnWVlZTb5uWe0EDkBISAghISHX\nU4YQQlyzPXvgiy+UycF+73iwSPHx8cTHx5v8epPD/dVXX8XW1pZ77rnnis/RXOF3ncvDXQgh2orB\nAIsWKdtLloCnp7r1NOePDd/ly5e36PUmhftHH31EbGws33//fd0+Dw8PMjIy6r7PzMzEw8PDlLcX\nQohWsX69svC1hwc8/bTa1bSuFo9zj4uLY+XKlWzZsgV7e/u6/VOnTiU6OprKykrS0tI4deoUo0eP\nNmuxQghhqosX66cZeO01cHRUt57W1mzLfc6cOezevZsLFy7Qt29fli9fTkREBJWVlYSGhgIwduxY\n1q1bh7+/P7NmzcLf3x8bGxvWrVt3xW4ZIYRoa88+C3o9TJigTDfQ0clNTEKIDm//fmU8u40N/PYb\n+PqqXVHLydwyQghxmZoaePRRMBph8eL2GeymkHAXQnRo69bVkJAAffvC88+rXU3bkW4ZIUSHlZJS\nxZAhUF6u5auv4K671K7IdNItI4QQQHl5BTNnnqO8XEtQ0Kl2HeymkHAXQnQ4paWlPPXUAY4c8aBr\n12qmT/9B7ZLanIS7EKJDKS4u5oMPNvPJJ6MAeOKJszg7l6hcVduTcBdCdBgFBQV89tlnfP75/6Oo\nyJ5RowqZMuW82mWp4romDhNCCEuh1+v5/PPPOX16KD/91B87OwP/+MeZDjmd77WQlrsQot07d+4c\nn332GXZ2/Xn//ZEAPPlkJp6elSpXph4JdyFEu5aTk8PGjRvp39+LDRtupKDAhtGjC5k5s3N2x9SS\nbhkhRLuVkZHBV199xQ033MDevb78+GN3nJyqefHFdKw6edNVwl0I0S6lpqYSExNDQEAApaXuvP56\nXwCWLMnA3b1K5erUJ+EuhGh3Tp48SWxsLEFBQTg6urBo0UBKS62ZMOESkyfr1S7PIki4CyHalWPH\njrFz506GDx+Ok5MTq1f3ITHREXf3Cp5/vvOOjvkjCXchRLuRkJDAjz/+yPDhw3F0dOTnn5355BN3\nrK2NvPpqGs7ONWqXaDEk3IUQ7cIvv/zCgQMHGDVqFPb29pw7p+Wll7wAeOSRLIYO7Xx3oTZHwl0I\nYdGMRiM//fQTR44cYdSoUdjZ2VFdDc89N4D8fC1jxhRy//15apdpcSTchRAWy2g0smvXLk6cOMGo\nUaPQarUArFrVlyNHuuLqWsnLL6d1+mGPTZFwF0JYJIPBwHfffUd6ejojR47ExkaJq2++6cHGjW7Y\n2Bj4179S6dmzWuVKLZOEuxDC4lRXVxMbG0tOTg4jRozA2toagBMnHAgP7w/A009nEBgo/exXIuEu\nhLAoVVVVxMTEcOnSJYYPH47V730uer0NTz89iIoKK+688wLTp19QuVLLJuEuhLAYFRUVfP3115SX\nlxMcHIzm90Hr5eUaFi8eRE6OHQEBJSxZclbGs19Fs5ch5s+fj06nIzAwsG6fXq8nNDSUwYMHM2nS\nJPLz8+uORURE4OPjg6+vL9u3b2+9qoUQHU5ZWRmff/45VVVVBAYG1gW7wQAvv+zF0aNOuLtX8Prr\np7Gzk3WYr6bZcJ83bx5xcXEN9kVGRhIaGkpycjITJ04kMjISgOPHj7Nx40aOHz9OXFwcjz76KAaD\nofUqF0J0GMXFxWzcuBErKyv8/f0bHHv33T5s394DR8caVq8+Ta9ecgH1WjQb7jfddBMuLi4N9sXE\nxBAWFgZAWFgYmzdvBmDLli3MmTMHrVaLl5cX3t7eHDhwoJXKFkJ0FAUFBURHR+Pg4MANN9zQ4NiW\nLT1Zv743VlZGIiJS8fYuV6nK9qfFo0Pz8vLQ6XQA6HQ68vKUmweys7Px9PSse56npydZWVlmKlMI\n0RHp9Xqio6NxcXHB29u7wbEffujOq68qI2OWLDnLuHGFapTYbl3XBVWNRlPXL3al401ZtmxZ3XZI\nSAghISHXU4YQoh06f/48n3/+OR4eHg0ahgAHDzrx3HMDMBg0PPRQNjNndr6RMfHx8cTHx5v8+haH\nu06nIzc3F3d3d3JycnBzcwPAw8ODjIyMuudlZmbi4eHR5HtcHu5CiM4nJyeHL7/8Ei8vL3r37t3g\nWFJSFxYv9qaqyopZs87x4IM5KlWprj82fJcvX96i17e4W2bq1KlERUUBEBUVxbRp0+r2R0dHU1lZ\nSVpaGqdOnWL06NEtfXshRAeXkZHBpk2bGDRoUKNgP3nSgccf96GkxJpbb9Xz979nyJBHEzXbcp8z\nZw67d+/mwoUL9O3bl5dffpmlS5cya9Ys1q9fj5eXF5s2bQLA39+fWbNm4e/vj42NDevWrWu2y0YI\n0fmkpaWxZcsW/P396dmzZ4NjJ0868OijgykosOHmm/NZtkyWyrseGqPR2KYDRjUaDW18SiGEBTh5\n8iTffvstgYGBdO/evcGx5GQHHnmkPthXrEhFqzVPTlRXV7N3714WLVpklvdTS0uzU+5QFUK0usTE\nRHbs2EFwcDBdu3ZtcOzYsS48+aRPXbBHRpov2DszCXchRKs6fPgwu3fvrls96XL793fl738fRFmZ\ndV2w29pKsJuDhLsQotXs27eP/fv3M3LkSBwcHBoc27HDhRde8KK62orbb7/Iiy+mYyOJZDbyUQoh\nzK6p1ZPqj0F0tBtvvOGJ0ajhnnvyWLgwUy6empmEuxDCrK60ehJAdTWsXNmPL790BZS1T+fPz5Xh\njq1Awl0IYTa1qyelpaU1WD0JoLDQmqVLB3LggDO2tgZefDGd2267pGK1HZuEuxDCLGpqaoiNjSU7\nO5uRI0fWrZ4EylDHZ54ZSEaGPT16VPHaaykEBckqSq1Jwl0Icd2utHoSwNatPYmM7EdFhRU+PqW8\n8UYKvXtXqlht5yDhLoS4LhUVFWzevJnS0tIGqyeVlVmxcmVfYmJ6ATB16gWWLDmLvb0MdWwLEu5C\nCJOVlZXx1VdfUVNTQ1BQUN3+o0cdefFFLzIy7LGzM7BkyVnuvPOiipV2PhLuQgiTlJSU8MUXX2Bj\nY0NAQACgjIb54IPefPhhbwwGDd7epfzzn+n4+JSpXG3nI+EuhGixwsJCNm3aRNeuXesW2fj1V0fC\nw/uTkuKARmPk/vtzefjhbLnjVCUS7kKIFrl06RKbNm2iV69eeHl5UVhozVtvefDVV8rYdU/Pcl58\n8QzDhxerXGnnJuEuzMJggHPnICMDsrKgoAAKC5WvoiLlq6am8escHMDJCRwdlUcnJ+jZE9zc6r/+\nMB2JUFHt6kl9+vShd+++bNzoynvv9aGgwAYbGwNhYXnMm5cjF00tgIS7aBG9Ho4cgcREOH4ckpLg\n7Fkl0CtbaXRbly5GPD01eHlR99W/f/12797IHY5tIDc3ly+++IL+/b04dcqXhQs9OXPGHoARI4p4\n5pmzDBwoC1hbCgl3cUVGI5w8Cd9/D7/8Avv3w+nTV35+z57Qty94ekKPHtC1Kzg7K19OTmBtbaCw\nsIhLl/RcuHCB/PwCNJougBPW1s5UVtpSUGDk0iUN+fm2FBTYUVjoQGmpDcnJkJzc9HmdnMDXF/z8\nlMfaL29vsLVtlY+m08nMzOTLL7/i4sX/x9tv+3LihPLrVL9+5SxYkMn48QXyA9bCSLiLBoqLIS4O\ntm2DnTuVFvnl7O0hOBiGDIGAAPD3V1rPnp7QpUvj99Pr9WRlZZGamkp6ejpdu2rp06crvXr1omfP\nnr/fxVgFND1MzmiEwkIjqalVpKcbyciwJjfXgfPnHbhwwQm93oXiYnsOHoSDBxu+1sZGCfygIBg6\ntP7R3d0cn1TnkZycziuvHOfnn/+P1FRnAHr2rGL+/BxmzDgvMzlaKFmJSVBSAl9/DZs2wfbtUFFR\nf8zVFSZOhJtughtvhMBAuGweqEaKi4vJysoiPT2d06dPU1NTg7OzMz169MDV1RVbMzaljUYjJSUl\nZGVVkJxsxZkz9mRkOJGd7cy5cz3R650xGhs3J93cGge+n5+08v/o3DkID7/Ahg32FBY6AeDqWklY\nWC7Tpl1oN/3qnXUlJgn3TspohD174KOP4PPPlRZ7rbFjYdo0uO02pYXe3FSsFRUV5OTkcPbsWVJS\nUsjPz8fZ2Znu3bvj6uraaHGGtlJRUcH586UkJVlx8qQ96enOnD3rQnZ2L8rL7Ro9X6s14u+vYehQ\nJeyDg5XHPyzz2WEYjUYqKyspLy+noqKCiooKysvLKS6uZMcOO2JiXNi3z4WaGuUvf9CgMmbPPsdf\n/nIRO7v29f9Xwr2NSLirq7QU/vtfWLsWjh6t3z92LMydC3fdpVygvJKamhry8vLIyMggJSWF3Nxc\nHB0d6datG66urnTr1s2iF0avrq4hNbWGY8esOXnSri70z5/v3uTz+/QxEBysITi4Pvi9veGyObFU\nZTQa68K5NqAv3y4rK6OkpISysrK67drnWFtbY21tjdFoR3Jyf379dSBHjvSjpES5SGplZWDcuEvc\nc89FRo0qard96hLubUTCXR16PaxZA2+9pWwD6HQwfz6EhcENNzT9OqPRiF6vJyMjg9TUVM6ePYut\nrS3dunWjV69euLi4NJj9r70qLbUiMdHq99B3IC2tKxkZLlRWNu6DcnAw4O9fw7Bh1gwbZoW/v3IB\nV6czfdTO5SH9x4CuDenS0lJKS0vrtsvKyqisrMTa2hobGxusra3RarV1oW1tbY2trS22trZotdrL\ntm3JzOzKwYPO7N/vzP/+15Wysvq/Q2/vUqZMuchtt+np2bPa1I/UYki4txEJ97Z14QKsWgVvvqmM\nNQcYNQqefBL++tem+5mLiorq+s1TUlIwGAwN+s21zXW6dyA1NZCZaceJE3YcO2ZDcrIS+np9011N\nzs41+PjU4Oenwdu7mv79K9DpSnF1LcXOrqxBS7q0tJTy8vJGIX15ONcGto2NTV04a7Va7OzsGnzf\n3G9KRiPo9TacPu1AYqIjx445kpjoyMWLDf8OfX1LuOWWfG65JZ8BA8rbbSu9KRLubUTCvW2UlcHq\n1RARUR/qoaHwwgvKxdHLVVRUkJ2dzdmzZzl16hRFRUV07doVFxcX3Nzc6NLUMJhOLD/fmlOnupCc\nbE9Ski3p6XacPetIaemVr8h26VJJr17FuLmV0bt3Oa6ulfTsaaBXLwNubuDqasDZ2YC9vaFFy81V\nVmooLbXi0iUtubm25OVpycuzJS/PlsxMO1JSHCgoaDycpWfPKkaPLmTUqCJGjy7E3b3KlI+iXeis\n4W7yIKaIiAg++eQTrKysCAwMZMOGDZSUlDB79mzOnDmDl5cXmzZtonv3pvsyReswGJQ+9eeeU+4W\nBbj1VnjpJaVfHZR+89zc3Lp+87y8PJycnOjWrRsDBw6kW7du6v0B2oHu3WsYNaqIUaOK6vbVtpDT\n0+3rvs6csSc315bsbFtKS205e7YHZ89e/f3t7Wvo0sWAg4MBrdbQ6HhVlRXFxdaUllpRVXX1nwSO\njjUMGlSGn18pQ4aUMGRICZ6eFR2qdS4aM6nlnp6ezoQJE0hKSsLOzo7Zs2dz++23k5iYSK9evViy\nZAkrVqzg0qVLREZGNjyhtNxbzfHj8PDDyigYUC7+vfYaTJxo5OLFiw36ze3t7XF2dqZXr1706NGj\nweIKwryMRigosCY7246cHFuysuy4eFHLpUs26PU2XLqkRa+3oajImvLyll2/sLY24uRUQ7du1eh0\nlb9/VaHTVdKnTwUDBpTj5lbVqYNcWu4t4OzsjFarpbS0FGtra0pLS+nTpw8RERHs3r0bgLCwMEJC\nQhqFuzC/sjIID4cVK6CqShnH/dJLZdx8czpnz6bx9tupGI1GnJ2d6dmzJ+PGjWuwtqVoXRqN0trv\n3r0Uf//SZp9rMEB5uRWlpVaUlVlTVaVpFMw2Nkqgd+lSg62tsVMHt7gyk/6H9+jRg8WLF9OvXz8c\nHBy49dZbCQ0NJS8vD51OB4BOpyMvL8+sxYrGDh6E++6DEyeU7++4I4uJE3dSUaHnwIFuuLi4EBwc\njIODg7qFimtiZQVduhjo0sUAtP+RKkI9JoV7SkoKq1evJj09nW7duvHXv/6VTz75pMFzNBrNFa/i\nL1u2rG47JCSEkJAQU8ro1KqrITzcyMsvG6mpscLN7SJhYT8xYkQ5bm79cHYeonaJQojrEB8fT3x8\nvMmvNyncDx48yLhx4+j5++1706dP55dffsHd3Z3c3Fzc3d3JycnBzc2tyddfHu6i5c6cgbvvhn37\nNICGadPO8NRT5+nSxVPt0oQQZvLHhu/y5ctb9HqTrqL5+vqyb98+ysrKMBqN7Ny5E39/f6ZMmUJU\nVBQAUVFRTJs2zZS3F8349lsYPhz27YPu3YtYteo3nn/+Al26SMerEKKeSS33oUOHcv/99zNy5Eis\nrKwYPnw4Dz30EEVFRcyaNYv169fXDYUU5lFTA8uXwyuvKKMvAgPPsnRpEjfc0Evt0oQQFkhuYmoH\niopgzhz45hvlgtu9955gwoQDBAYGqF2aEBZPhkIKi3TmDEyZokzy1aMHRESkUFm5k4CAMWqXJoSw\nYHLnigXbtw9Gj1aC/YYbIDb2IqWlWxk6dKjcdCSEaJYkhIWKjYUJE5QFE/78Z4iPr+Do0a/x9vaW\nuV6EEFcl4W6BPv0U7rxTufN0/nwl6A8d2omtrS29m5tsXQghfifhbmHefFNZNKO6GpYsgQ8+gKSk\n30hLS8PPz0/t8oQQ7YSEuwVZuRIWLFC2//UvZa6YCxfOs2vXLoKCgqSfXQhxzWS0jIV47TWlpQ7w\n3nvw4IPKPOtbtmxh0KBBqq1FKoRonyTcLcBrr8HTTyvbH3wADzygbO/cuROtVkufPn3UK04I0S7J\n7/kqW7Om6WA/evQoqamp0s8uhDCJhLuKPv4YFi5Utt97rz7Yz59X+tmHDh3aIRafFkK0PQl3lXzz\nDcybp2y//rrSxw5QWVlJTEwMAwYMkH52IYTJJNxV8NNPMHOmMhnYs8/CU0/VH9u5cyc2NjZ4eHio\nV6AQot2TcG9jJ08qc8WUlyut9VdfrT927NgxUlJSpJ9dCHHdJNzb0IUL8Je/QH6+cgfqO+9Qt/7l\nhQsX+P777wkKCpJ+diHEdZNwbyMVFXDXXZCSoiy28emnUJvhtf3sXl5eODk5qVuoEKJDkHBvA0aj\nMhLmp5/AwwO2boXLr5V+//33WFtb4+kpy+QJIcxDwr0N/OtfSkvd0RG2bYPL70lKTEzk9OnT0s8u\nhDArCfdWtmMHPPecsv3ppxAcXH/s4sWL0s8uhGgVEu6tKD0d7r4bDAZ44QXlImqtqqoqYmJi6N+/\nv/SzCyHMTsK9lZSVwfTpoNfD7bfDsmUNj3///fdoNBrpZxdCtAoJ91by+ONw+DAMGgSffKIsbF3r\n+PHjJCcn4+/vr16BQogOTcK9FXz6KXz4Idjbw1dfgYtL/bGLFy+yY8cOmTdGCNGqJNzN7PRpePhh\nZXvNGggKqj9W288u49mFEK3N5HDPz89n5syZ+Pn54e/vz/79+9Hr9YSGhjJ48GAmTZpEfn6+OWu1\neBUVMHs2FBfDrFn1k4HV2rVrF4D0swshWp3J4f7kk09y++23k5SUxG+//Yavry+RkZGEhoaSnJzM\nxIkTiYyMNGetFm/pUkhIgAEDlCl8a6cWAEhKSuLkyZMEBASoV6AQotMwKdwLCgrYs2cP8+fPB8DG\nxoZu3boRExNDWFgYAGFhYWzevNl8lVq47dth9WqwsYHoaOjWrf6YXq9n+/bt0s8uhGgzJoV7Wloa\nrq6uzJs3j+HDh/Pggw9SUlJCXl4eOp0OAJ1OR15enlmLtVSXLsHvP+dYvhxGj64/VlVVxdatW+nX\nr5/0swsh2oxJa6hWV1eTkJDAW2+9xahRo1i4cGGjLhiNRoPm8n6Jyyy7bNB3SEgIISEhppRhMZ54\nArKyYOzY+kWua8XHx1NTU0O/fv3UKU4I0S7Fx8cTHx9v8us1RqPR2NIX5ebmMnbsWNLS0gD46aef\niIiIIDU1lR9++AF3d3dycnK45ZZbOHHiRMMTajSYcEqL9fnnysXTLl3gyBHw8ak/duLECbZv386Y\nMWOwsZG1yIVQQ3V1NXv37mXRokVql3JdWpqdJnXLuLu707dvX5KTkwFl9aCAgACmTJlCVFQUAFFR\nUUybNs2Ut283cnLqhz2+9lrDYNfr9Xz33XcEBQVJsAsh2pzJqfPmm29y7733UllZyaBBg9iwYQM1\nNTXMmjWL9evX4+XlxaZNm8xZq8V57DFleoFbb60PeWjYz961a1f1ChRCdFomh/vQoUP53//+12j/\nzp07r6ug9uLLL+Hrr6FrV3j//YbDHnfv3k11dbX0swshVCN3qJrg0iVl7hiAyEjo27f+2IkTJ0hM\nTGTIkCHqFCeEEEi4m+Tvf4fcXPjTnxp2x1y6dEn62YUQFkHCvYW+/16ZFMzWVumOqZ3tsbq6mq1b\nt9K3b1+cnZ3VLVII0elJuLdAeXl9S/3FF8HXt/7Y7t27qaqqon///uoUJ4QQl5Fwb4GVK5VZH/38\n4Omn6/efPHmSY8eOST+7EMJiSLhfo7Q0CA9XttetU7plQPrZhRCWScL9Gi1YoHTL3HMP1M6WUF1d\nzbZt2/Dw8JB+diGERZFwvwZbt8K2beDsrNyJWuvHH3+koqICLy8v1WoTQoimSLhfRVmZ0moHePll\n6N1b2U5OTubo0aMEBgaqV5wQQlyBhPtVvPEGpKdDYKAy3QAoq1DFxcURGBgo/extaNOmTYwfP55j\nx46pXYoQFk/CvRk5ORARoWyvWaMsxFFTU8PWrVvx8PCg2+UrcohW95e//AU7OztZzUqIayDh3ox/\n/ANKSuDOO+GWW5R9e/bskX52lRw8eJBhw4ZdcZ0AIUQ9CfcrSEiAjz4CrVYZ3w5w+vRpfv31V+ln\nV8n+/fvRaDTExcURHh7O6dOn1S5JCIsl4d4EoxEWLVIen3hCmae9oKCAb7/9liFDhkg/exuIjo5m\n4sSJzJ07lzNnzgBKuN97773cdttt3Hzzzaxbt07lKoWwXBLuTdi8GX78EXr2hBdeUPrZt23bRp8+\nfejevbva5XV4Bw8eZNWqVaxevZqSkhL++c9/kpubi9ForPut6eLFi+Tn56tcqRCWS8L9D6qrYelS\nZXv5cujeXVlGsLS0VPrZ28ibb77J2LFjGTx4MEajEZ1OR1JSEsHBwXXP2bdvH+PGjVOxSiEsm4T7\nH3z4ISS6AwJPAAAU/UlEQVQng7c3PPQQpKSkcOTIEYKCgtQurVM4duwYx48fJzQ0FDs7OzZv3syr\nr76Ko6Nj3apWZ8+e5fTp08ydO1flaoWwXBLulykthWXLlO1XXoHS0gJiY2Oln70NxcbGAjRqlY8a\nNQorKyu2bdvGZ599xjvvvIO9vb0aJQrRLkhiXWbtWmVs+4gRMH16DZ9//g29e/eWfvY2tHv3bgYO\nHIiLi0uD/RqNhieffBKAO+64Q43ShGhXpOX+O71eWTIPlMe9e3+ipKSEAQMGqFtYJ3L27FnOnTvX\noG9dCGEaCfffRURAQQH8+c8wcGAqhw8flvHsbax2wXWZF1+I6yfhjtIV89ZbyvbzzxcTGxtLYGAg\nWq1W3cI6mUOHDgHg5+enciVCtH8S7ijdMOXlMG2akaysGHQ6nfSzq+DQoUPY2tpKV5gQZmByuNfU\n1DBs2DCmTJkCgF6vJzQ0lMGDBzNp0qR2c4NJVha8+66yfeedCZSUlDBw4EB1i+qEzpw5g16vx9vb\nG2tra7XLEaLdMznc16xZg7+/f90kTpGRkYSGhpKcnMzEiROJrL06aeEiIqCiAiZPLiE/f4/0s6vk\n8OHDAAwePFjlSoToGEwK98zMTGJjY/nb3/6G0WgEICYmhrCwMADCwsLYvHmz+apsJRkZ8P77oNEY\nCQ7ewpAhQ6SfXSUJCQkAeHt7q1yJEB2DSeG+aNEiVq5ciZVV/cvz8vLQ6XQA6HQ68vLyzFNhKwoP\nh8pKuPHGswQH2zQaWy3aztGjRwHLCPeamhqTX1tdXW3GSoQwXYtvYtq2bRtubm4MGzaM+Pj4Jp+j\n0WianXN7We1toEBISAghtStOt6GzZ2H9eqXVPmnSfgYNGtTmNQjFpUuXyMzMRKPRqP73sGvXLkpK\nSuquJbXUhg0bGD16NEOHDjVzZaKziY+Pv2LGXosWh/vevXuJiYkhNjaW8vJyCgsLue+++9DpdOTm\n5uLu7k5OTg5ubm5XfI/Lw10tK1dCVRUMH57Mbbf1U7ucTu23334DwMXFpU1GKWVkZPD6668zcOBA\nSkpKWLp0KRqNhkOHDnH48GEWL15s8nvPmzePxYsXs3Dhwmse9bNq1Sp27tzJuXPn+Pe//82IESNM\nPr/oOP7Y8F2+fHmLXt/ibpnw8HAyMjJIS0sjOjqaCRMm8PHHHzN16lSioqIAiIqKYtq0aS196zaT\nlwcffKBcK3j4Yb30s6usNtzbokumqqqKxx9/nIkTJ3Lx4kW2bNlCSUkJxcXFrF27lscff/y63t/G\nxoZnn32Wl1566Zq7aBYtWkRYWBi2trZyQV+YzXWPc6/tflm6dCk7duxg8ODB7Nq1i6W18+ZaoDfe\nMFJermHEiCyGD7dVu5xOr3bBax8fn1Y/1y+//EJ2djbDhw9n1qxZrF27FicnJzZs2MDkyZOxs7O7\n7nO4u7szaNAgtm3bds2vOXLkCP7+/tjayr9HYR7XNXHY+PHjGT9+PAA9evRg586dZimqNV26BG+/\nbQCseeKJQrXL6fRqamo4fvw40DbhfujQIVxcXPDw8MDDwwOAsrIyNm/ezNdff22288yePZtnn332\nmn+DPXz4MFOnTjXb+YXodHeovvUWlJRY4++fzZAhpWqX0+mlp6dTXl6ORqNpk3BPTEzE39+/wb6f\nfvqJPn364OzsbLbzDB48mIKCAk6cOHHV52ZmZnLhwgWGDx9utvML0amm/C0uhtWrle077vgNcFW1\nHkFdq93a2rpV7wxetmwZer2eX3/9FS8vLxYsWICHhwfPPPMM+/fvb3YxlqSkJGJjY7GysiInJ4fn\nn3+er776iqKiIs6fP89DDz2Ep6dng9dYWVkRHBzMvn378PX1bXDsf//7H19//TW9e/emqKiIQYMG\nYW1t3WiEjSnnFaJWpwr3995TpvYdOrQUX99cJNzVVxvuAwcObNUFUZYtW0ZWVhbTpk3jscceazAK\nITk5mbvuuqvJ12VmZrJ161aWLFlS9z7z5s1j2bJlGAwGHnzwQW644QbuvffeRq/t168fycnJDfZt\n2bKFdevW8cknn+Dq6kpubi4zZswgICCgweIj13NeIaATdctUVcGqVcr2Qw9dpJlh+KIN1Yb7DTfc\n0OrnOnnyJNB4ioPs7Oy6Jfz+6NNPP+WJJ56o+76srAxnZ2cCAwNxd3dn7ty5VxwT37VrV7Kzs+u+\nT05OJiIigsWLF+PqqjQs3N3dcXBwaNQlcz3nFQI6Ubhv2gSZmeDrC+PHF6tdjkC5mHr69Gmgbab5\nTU5OxsnJiT59+jTYX1xcfMVwv++++3BwcKj7/ujRo4wePRpQ7sResGDBFfvqu3fvTnFx/b+1devW\n4ejoyMSJE+v2paamUlBQ0Cjcr+e8QkAnCXejEV57TdlevBisOsWf2vKlp6dTWVmJRqNps3BvamIy\njUaDwWBo8jWX/yBIT0/n/PnzjBw58prOZzAY6uZeKioq4pdffmHMmDENZr08dOhQXf+8uc4rBHSS\ncN+1C44cATc3mDtX7WpErdr+aBsbmzbplklOTm7yPF27dqWw8OrDYg8ePIhWq21w8TUzM/OKzy8s\nLKz7jSAjIwODwdDowu3Bgwfx8/PDwcGBrKwss5xXCOgk4V7ban/iCbjsmpVQ2alTpwDlztTWvku4\noKCAvLy8Jodb9unTp8n1B8rLy1m7dm1d19H+/fvx8fGpu9HJYDDw8ccfN3vO2rH0jo6OgNLHfvn7\nJyQk1HXJREdHm+W8QkAnGC1z7BjExUGXLvDII2pXIy5XG15tsWZq7cXUpsI9ODiYtLS0Rvt//vln\nPv74Y3x9fbGxsSEjI6NB3/yHH37Y7EXNtLQ0xowZAygjZ3x8fOpa59XV1axYsYKqqio8PT3R6/X0\n6NHDLOcVAjpBuL/+uvI4fz707KluLaKh2nAPCAho9XOdOHGCrl27NtnnPnbsWN54441G+0eMGMGU\nKVM4ceIEJ0+e5KOPPiIyMpLw8HC0Wi3jx4+/4g+m6upqfvvtNxYsWAAo/fqRkZG88cYb5OXlYTAY\neOCBBxgxYgTbtm3jxIkTdaNjrue8QtTSGGuv+LTVCTUa2uqUeXnQrx9UV0NyMtTOJpuUlMT+/ftl\nkiYVFRUVMWHCBDQaDZs2bcLLy6tVz/fcc89RU1PDihUrGh2rrKxk8uTJREdH1w1RvF6//vor4eHh\nbNy40SzvJ0xXXV3N3r17WbRokdqlXJeWZmeH7nN/911lMY6pU+uDXViGlJQUAJydnVst2KOionjs\nsccAZTz95UMQL2dra8vs2bP57LPPzHbu//73v3KDkVBVhw33ykp45x1l+/ffjIUFSU1NBWg0BNCc\nYmNjsbW15dSpU2i12iuGO8D999/P3r17r2nUzNWkp6eTm5sr/eJCVR023L/4AnJzYcgQUGGhJ3EV\nteE+bNiwVjvHfffdh6urKxs2bGDlypUNxpf/kb29PS+88AKvvPLKdXUbVlRUsHLlSl599dVmVyMT\norV12Auqa9cqjwsWIFMNWKDaYZCt2XK/4447uOOOO675+QEBAcyYMYONGzdy9913m3TODRs28Nhj\nj8mEXkJ1HTLc9+9XvlxcQLo9LdOpU6dwcHBoNGOi2saMGVM3fNEUDz/8sBmrEcJ0HbJb5s03lccH\nH1TGtwvLkpOTQ1FREUOGDGm2q0QIYboOF+45OcokYVZW8OijalcjmpKUlAQgC0EL0Yo6XLi//74y\nve+dd0L//mpXI5qSmJgIUDfLoRDC/DpUuFdXK+EO0mq3ZMeOHcPR0bFN7kwVorPqUOH+zTfKnO0+\nPjBhgtrViKaUl5dz7NgxxowZg5XMvSxEq+lQ/7tqb1p6+GGZs91SHTx4kMrKSsaPH692KUJ0aB0m\nAlNS4LvvwM4O/u//1K5G1HrttdeYM2cO1dXVAMTFxeHs7Nzs3aJCiOtnUrhnZGRwyy23EBAQwJAh\nQ1j7+x1Der2e0NBQBg8ezKRJk5qcI7u1vPee8jh7Nvw+c6qwAAcOHKC8vByDwUBubi67du3innvu\nqZubXAjROkwKd61Wy6pVq0hMTGTfvn28/fbbJCUlERkZSWhoKMnJyUycOJHIyEhz19ukigr48ENl\nW+ZstyxDhw5l0qRJFBYW8vLLL9OvXz/CwsLULkuIDs+kcHd3d6+7bdzJyQk/Pz+ysrKIiYmp+48b\nFhbG5s2bzVdpM774Ai5cgOBguI6bC0UreOyxx0hMTGTatGnY2try5ptvYmPT9I3R1dXVvPPOO3z5\n5Zds3LiRRYsWyXJyQpjouqcfSE9P5/Dhw4wZM4a8vDx0Oh2grNCel5d33QVei3ffVR4ffljmkbE0\n3bt356233rqm50ZERODj48OMGTPIz8/n3XfflTlahDDRdYV7cXExM2bMYM2aNQ2WAQNlYvkrzYq3\nbNmyuu2QkBBCrmPaxuRk2LMHHB3hnntMfhuhslOnTrFjxw6eeeYZQFmlqXZtUSE6o/j4eOLj401+\nvcnhXlVVxYwZM7jvvvuYNm0aoLTWc3NzcXd3JycnBzc3tyZfe3m4X6/165XHWbPgDz9fRDty4MAB\ngoODsbW1rft+1KhRFBUVNWo4CNEZ/LHhu3z58ha93qQ+d6PRyAMPPIC/vz8LFy6s2z916lSioqIA\nZRWc2tBvLVVV8Pvp+NvfWvVUopU5OzvTq1cvAEpLS/nhhx8ICgri22+/VbkyIdonk1ruP//8M598\n8glBQUF1iy1ERESwdOlSZs2axfr16/Hy8mLTpk1mLfaPvvlGWSfVzw/Gjm3VU4lWduutt3LkyBG+\n++47Kisrue2229i7d6/FTQksRHthUrj/6U9/wmAwNHls586d11VQS9R2yTzwgFxIbe9sbW154YUX\n1C5DiA6j3d6hmpUFsbGg1cJ996ldjRBCWJZ2G+5RUWAwwNSpcIXrtkII0Wm1y3A3GOq7ZORCqhBC\nNNYuw/3nnyE1FTw9ITRU7WqEEMLytMtw/89/lMe5c0GW4BRCiMbaXbiXlSlrpIJcSBVCiCtpd+G+\ndSsUFsLIkeDvr3Y1QghhmdpduNd2yUirXQghrqxdhXteHsTFgY0N3H232tUIIYTlalfhHh0NNTUw\nebKMbRdCiOa0q3Cv7ZK5/3516xBCCEvXbsL92DFISIBu3eCOO9SuRgghLFu7CfePP1YeZ88Ge3t1\naxFCCEvXLsK9pgY+/VTZli4ZIYS4unYR7j/8oMwCOXAgjBundjVCCGH52kW4Xz62XeZtF0KIq7P4\ncC8pga++UrblxiUhhLg2Fh/uW7cqAX/jjTBokNrVCCFE+2Dx4b5xo/Iod6QKIcS1s+hwLyiAb79V\n+tn/+le1qxFCiPbDosN9yxaoqICbb4Y+fdSuRggh2g+LDnfpkhFCCNNYbLhfvAjbtysrLc2YoXY1\nQgjRvpg93OPi4vD19cXHx4cVK1aY/D5ffw3V1TBhAri6mrFAC3Lo0CG1S7AY8lnUk8+innwWpjNr\nuNfU1PD4448TFxfH8ePH+eyzz0hKSjLpvTpDl4z8w60nn0U9+SzqyWdhOrOG+4EDB/D29sbLywut\nVsvdd9/Nli1bWvw+eXmwaxdotXDXXeasUAghOgcbc75ZVlYWffv2rfve09OT/fv3t/h9vvwSDAZl\nUQ4XF3NWqNBoNOj1ehISEsz/5i2Qk5Ojeg2WQj6LevJZ1DPHZ2EwGLCxMWvUtQtm/RNrrnHil2t9\n3jffdPy5ZLZu3ap2CRZDPot68lnUM9dnsWDBArO8T3th1nD38PAgIyOj7vuMjAw8PT0bPMdoNJrz\nlEIIIZpg1j73kSNHcurUKdLT06msrGTjxo1MnTrVnKcQQghxDczacrexseGtt97i1ltvpaamhgce\neAA/Pz9znkIIIcQ1MPs498mTJ3Py5ElOnz7Ns88+W7ffXOPfO4KMjAxuueUWAgICGDJkCGvXrlW7\nJFXV1NQwbNgwpkyZonYpqsrPz2fmzJn4+fnh7+/Pvn371C5JNREREQQEBBAYGMg999xDRUWF2iW1\nmfnz56PT6QgMDKzbp9frCQ0NZfDgwUyaNIn8/Pyrvk+b3KFqzvHvHYFWq2XVqlUkJiayb98+3n77\n7U79eaxZswZ/f/9rvtDeUT355JPcfvvtJCUl8dtvv3Xa33rT09N5//33SUhI4OjRo9TU1BAdHa12\nWW1m3rx5xMXFNdgXGRlJaGgoycnJTJw4kcjIyKu+T5uEu7nGv3cU7u7uBAcHA+Dk5ISfnx/Z2dkq\nV6WOzMxMYmNj+dvf/tapL7YXFBSwZ88e5s+fDyhdnN26dVO5KnU4Ozuj1WopLS2lurqa0tJSPDw8\n1C6rzdx00024/GEMeExMDGFhYQCEhYWxefPmq75Pm4R7U+Pfs7Ky2uLUFi89PZ3Dhw8zZswYtUtR\nxaJFi1i5ciVWVhY7zVGbSEtLw9XVlXnz5jF8+HAefPBBSktL1S5LFT169GDx4sX069ePPn360L17\nd/785z+rXZaq8vLy0Ol0AOh0OvLy8q76mjb5H9XZf92+kuLiYmbOnMmaNWtwcnJSu5w2t23bNtzc\n3Bg2bFinbrUDVFdXk5CQwKOPPkpCQgKOjo7X9Kt3R5SSksLq1atJT08nOzub4uJiPv30U7XLshga\njeaaMrVNwv1axr93NlVVVcyYMYO5c+cybdo0tctRxd69e4mJiWHAgAHMmTOHXbt2cf/996tdlio8\nPT3x9PRk1KhRAMycObPT3qV68OBBxo0bR8+ePbGxsWH69Ons3btX7bJUpdPpyM3NBZS7dt3c3K76\nmjYJdxn/3pDRaOSBBx7A39+fhQsXql2OasLDw8nIyCAtLY3o6GgmTJjAf/7zH7XLUoW7uzt9+/Yl\nOTkZgJ07dxIQEKByVerw9fVl3759lJWVYTQa2blzJ/7+/mqXpaqpU6cSFRUFQFRU1LU1CI1tJDY2\n1jh48GDjoEGDjOHh4W11Wou0Z88eo0ajMQ4dOtQYHBxsDA4ONn777bdql6Wq+Ph445QpU9QuQ1VH\njhwxjhw50hgUFGS86667jPn5+WqXpJoVK1YY/f39jUOGDDHef//9xsrKSrVLajN33323sXfv3kat\nVmv09PQ0fvjhh8aLFy8aJ06caPTx8TGGhoYaL126dNX30RiNnbyzUwghOqDOPURBCCE6KAl3IYTo\ngCTchRCiA5JwF0KIDkjCXQghOiAJdyGE6ID+P77cvv/6VvTLAAAAAElFTkSuQmCC\n"
      }
     ],
     "prompt_number": 4
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Compute the integral both at high accuracy and with the trapezoid approximation"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from scipy.integrate import quad, trapz\n",
      "integral, error = quad(f, 1, 9)\n",
      "print \"The integral is:\", integral, \"+/-\", error\n",
      "print \"The trapezoid approximation with\", len(xint), \"points is:\", trapz(yint, xint)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "The integral is: 680.0 +/- 7.54951656745e-12\n",
        "The trapezoid approximation with 6 points is: 621.286411141\n"
       ]
      }
     ],
     "prompt_number": 5
    },
    {
     "cell_type": "code",
     "collapsed": true,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}