questions.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Integration exercises\n",
    "\n",
    "1. Fackler and Miranda exercise 5.1: Demand for a commodity is given by $q(p)=2p^{-0.5}$. The price of the good falls from 4 to 1. Compute the change in consumer surplus **analytically** using calculus and **numerically** using three different methods: a Gauss-Legendre rule, Monte-Carlo, and Quasi Monte Carlo.\n",
    "\n",
    "    1. Start out by plotting the demand function. In that plot, label the axis properly ($p$ on the $x$ axis), and add 2 horizontal lines for the equilibrium quantities at both prices $p=1,p=4$. This should help you visualize the consumer surplus, and it should guide you in the analytical solution. The simplest way to solve this first part would be to edit this `IJulia` notebook, just adding a code cell below (for the plot), and then one for the analytic solution. You can get and see how to use `IJulia` [by clicking on this link](https://github.com/JuliaLang/IJulia.jl).\n",
    "    1. To implement the integration rules, write 3 functions, one for each sub question, and detailed below. Let each function take as argument the number of integration points. Use $n=10,100,1000$ points. Each of those functions should produce the result (of course), as well as a plot where we can see integration nodes vs function value. We want to get an idea where each method places the points, and how this might influence the different results. Ideally, the result would tell us how far the corresponding method is away from your analytic solution.\n",
    "        1. Gauss-legendre rule. Note that you have to change the function domain to $[-1,1]$ first. This is achieved with the following transformation\n",
    "    $$ \\int_a^b f(x)\\,dx = \\frac{b-a}{2} \\int_{-1}^1 f\\left(\\frac{b-a}{2}x + \\frac{a+b}{2}\\right)dx $$\n",
    "        1. monte-carlo, again taking $n$ as an argument.\n",
    "        1. pseudo monte-carlo. Use a Sobol Sequence.\n",
    "\n",
    "\n",
    "2. Fackler and Miranda exercise 5.5: A government stabilizes the supply of a commodity at $S=2$, but allows the price to be determined on the market. Domestic and export demand for the commodity are given by $D = \\tilde{\\theta_1} p^{-1}, X = \\tilde{\\theta_2} p^{-0.5}$, and where $\\log \\tilde{\\theta_1}$ and $\\log \\tilde{\\theta_1}$ are both normally distributed with means 0, variances 0.02 and 0.01 respectively, and covariance 0.01.\n",
    "\n",
    "    1. compute the expected price $E[p]$ and expected variance $V[p]$ using a kronecker product rule with gauss-hermite grids (10) in each dimension of the shock.\n",
    "    1. perform the same computation with a monte carlo integration scheme. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Here is the code used to plot the price/quantities function:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAArkAAAI6CAYAAAAjelDCAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAIABJREFUeJzs3Xd8FHX+x/H3bEgCIYWEFlAgoUoRQiiiGIqnIkqzIVFQmoge5fBOxYKACJ7KCeJPgVMEUQnY4CjneUhHBAVEFMGTGnpPoYcwvz+GhIQkkGw2zO7k9Xw85rFkZnb2s9lE33z5zPdrmKZpCgAAAHAQl90FAAAAAJ5GyAUAAIDjEHIBAADgOIRcAAAAOA4hFwAAAI5DyAUAAIDjEHIBAADgOIRcAAAAOA4hFwAAAI5DyAXgUS6XS7fddlu2fSNGjJDL5dLy5cvdPrcwpk2bJpfLpenTp3vsmrmZMGGC6tevr6CgILlcLk2YMKFIXw9FI7efSwC+h5AL+CiXy3XFzc/Pz6NBMb8Mw5BhGFfdV9BzPV2Tp82cOVN/+ctfVKpUKQ0ZMkQjRoxQixYtivQ1vUVycrLefPNNde/eXfXr15e/v79cLpcWL15sd2nFxrlz59SgQQO5XC5VrVrV7nIAr1DC7gIAuM8wDI0YMUKmaeZ6PCoq6toWJGnz5s0KCgry+LmFcd999+nmm29WpUqViuw1FixYIMMwtGDBAlWsWLHIXscb7dy5U88995wMw9D111+v8uXL6+DBg3aXVaw8//zz2r17d5H/ZQ7wJYRcwMcNGzbM7hKyqV27dpGcWxghISEKCQkp0tfYt2+fJBW7gCtJ1apV06JFi9S4cWOVKVNGvXr1KvLWEFyydOlSjR8/XhMnTlT//v3tLgfwGrQrAMXIoUOH1KdPH0VGRiooKEiNGzfW9OnTtWzZMrlcLr3yyivZzo+KilL16tVzvVZB+mzzcrVzP/roI8XGxiooKEgVK1ZUnz59ch0hbNOmjfz8/JSWlqZXXnlFN9xwg0qWLKnevXtLyrsnd+nSperXr5/q16+vsLAwBQUF6cYbb9Qrr7yis2fP5us9jBw5Ui6XS0uWLJFpmtnaRbJatGiR7rrrLpUtW1YlS5ZUnTp19PzzzyslJaXA7+dKMr6n+/fvV48ePVSxYkUFBQWpadOmSkhIyNd7KqgyZcqobdu2KlOmTKGvlZKSojfffFN/+tOfVKVKFQUGBqpChQrq3LmzVq9enetzMt7z0aNH1a9fP1WuXFklS5ZUgwYNNG3atFyfk5aWplGjRqlmzZoqWbKkqlevrmHDhuncuXMFqnfXrl1yuVzq3bu3fv/9d3Xp0kVly5ZVcHCw4uLitHDhwoJ+CwokNTVVPXv21B133KF+/foV6WsBvoaRXKCYOHr0qG6++Wbt3LlTcXFxatmypfbv368nn3xSd9xxR549s3kp6j7Xt956SwsXLtRDDz2k9u3ba+XKlZo6daqWLVumNWvWqGzZsjnqvP/++7V27Vq1b99e9957rypUqHDFWl9//XX9/vvvuuWWW9ShQwedOXNG3333nUaMGKFly5bp22+/vep7bNu2rQzD0NSpU5WYmJjZPpL1eZMnT9ZTTz2l4OBgPfjgg6pQoYKWLl2q119/XfPnz9d3332n0NDQfL+fqzl+/LhuueUWhYeHq3fv3kpKStJnn32mRx55RPv27dNf//rXfF3HDps3b9ZLL72k1q1bq0OHDgoPD1diYqLmzp2rr7/+WvPnz9edd96Z43lJSUlq2bKlAgMD9eCDD+rs2bP6/PPP1bt3b/n5+alHjx7Zzn/wwQc1d+5c1axZUwMHDtS5c+c0depU/fLLL27VvX37dt18881q2LCh+vfvr/3792vWrFlq3769EhIS9OCDD7p13asZOHCgkpOTNWXKlCK5PuDLCLmAjxs5cmSu+0uWLKnnnnsu8+vnn39eO3fu1JAhQzR27NjM/QMGDPC6G6RM09R//vMf/fDDD2rYsGHm/qefflrjx4/X0KFD9f777+d4TmJiojZt2qTw8PB8vc7EiRNz7VsePny4Xn31VX3xxRdXDSetWrVSq1attGTJEiUmJuZoH0lMTNTgwYMVEhKiH3/8UbVq1co89uc//1kTJ07Us88+q0mTJhX6/WTYuHGjunbtmm3kdujQoYqNjdWLL76o+++/P9v7HjFiRIH+wtKmTRu1bt26QDXlV7169bR//35FRERk279v3z41a9ZMQ4YM0aZNm3I87+eff1bfvn01adKkzPcyePBgNWzYUK+//nq2kDtjxgzNnTtXt9xyixYvXqyAgABJ1u9S06ZN3frL24oVK/TMM8/o73//e+a+jN+t/v37q3379goODpZk3ag3bty4Ar3Ovffem+13QZJmz56t6dOn68MPP9R1111X4JoBxzMB+CTDMEyXy5XnFhERkXluWlqaWbp0aTMsLMxMSUnJca2ePXuaLpfLHDlyZLb9UVFRZnR0dK6vP2LECNPlcpnLli3LUVfbtm0Lda5hGObjjz+e4zWTk5PNMmXKmEFBQea5c+cy97dp08Z0uVzmvHnzcq112rRppsvlMj/66KNcj1/u6NGjpmEYZp8+ffJ1ftYaLvfqq6+aLpfLfOmll3IcO378uBkaGlrg93MlhmGY/v7+5s6dO3Mcy/gcXnnllRzPudLP0uXb5T8nl8v4eVq0aFGB67+SQYMGmS6Xy9y9e3eO+oODg83U1NQcz2ndurXpcrnMkydPZu67/fbbc/15NE3rZyW3n8u87Ny50zQMwwwPDzdPnDiR43jG92L69Ok5nlOQ7/nlP7sHDx40y5cvb3bo0CHH96JKlSr5qh1wOnpyAR+Xnp6e63b06NHMc7Zs2aJTp04pJiYm1xuw2rRpk+cMDXYwDEOtWrXKsT80NFQxMTE6c+aMNm/enON4s2bNCvQ6p06d0pgxY9S8eXOVKVNGfn5+crlcKleunAzD0N69e91+Dxl++uknSVZbw+XKlCmjxo0b68yZM9qyZUuO4wV9PxmqVq2qatWq5dif8Tln1JThwoULef4c5ba9/PLLbtWVX9999526du2qqlWrqmTJkpl9zu+8844k5fq51KpVK3OkNKsqVapIslo4Mvz0009yuVxq2bJljvPbtGnjVs2xsbEqXbp0rte7/HterVq1An/PH3300WzX7du3ry5cuKAPPvjArXqB4oB2BaAYSE5OlpT3nf+RkZHXspx8uVqtGe8pP8/Jzfnz59W2bVv9+OOPuvHGG9WtWzeVL19e/v7+kqx/ws/vzWdXklFnXtOXZexPSkrKcczdmRrc+d55i9mzZ+vBBx9UqVKldMcdd6hGjRoqXbp05s19y5cvz/VzyeumtxIlrP/NpaenZ+5LTk5WREREjpsDJfd/F67l93z69OmaP3++pk+fXixn8wDyi5ALFANhYWGSlOfcpQcOHMh1v8vlUlpaWq7HcgtlnnS1WjPek7v+9a9/6ccff1Tv3r1zjIYdOHBAI0aMKNT1M2TUeeDAAdWtWzfH8f3792c7zxMK+r3zpp7cYcOGKTAwUOvWrcsxxdy+ffs8ssBJWFiYjh07pvT09BxBN6/fhaspyPfcnZ7cLl26qFGjRpIu/evAo48+mmOEV5L27Nkjl8slwzB0/PjxbDc1AsUJIRcoBm644QYFBQVpw4YNSk1NzdGysGTJklz/hxseHq5ffvkl1zDw448/Flm9pmlq2bJl6t69e7b9KSkp2rBhg0qWLJlrYCyIrVu3yjAM3XvvvTmOLV26tFDXzqpx48b66quvtHTp0hwtC8nJyR57P1klJiYqMTExx8pXS5YskWT903pWr7zySoECl2EYRRZyt23bpgYNGuQIuKZpasWKFR55jdjYWC1atEgrV67M8T4yvkcFtX79ep08eTJHy0LG71bjxo0z9yUlJRX4ex4dHZ0Zcm+++WadPHky1/M++OADlS5dWvHx8TIMQ4GBgW68G8AZ6MkFioESJUrokUceUUpKSo4RyrVr12rGjBm5Pq958+Y6f/68pk6dmm3/tGnTtGrVqqIqV5L08ccfa8OGDdn2DR8+XMnJyXr44Ycz2wrcFRUVJdM0cwTa7du3a+jQoR6bHq179+7y9/fXO++8o23btmU79tJLLyklJUU9evQo9PvJKj09Xc8991y2PusdO3ZowoQJ8vf31yOPPJLtfG/qyY2KitIff/yRY0R1+PDhufZhu6NXr14yTVMvvvhittaHY8eOafTo0W599snJyTlmOsn43SpTpky2v0wVtie3a9eu+uc//5nrJll/Of3nP/+pyZMnE3JRrDGSC/i4vKYQk7JPOzRmzBgtWrRI48eP148//qhbb71V+/bt02effaZ77rlHc+fOzfH8gQMHaurUqerfv7++/fZbValSRRs2bNDq1avVsWNHzZ8/v8jeV/v27dWyZUt17dpVlSpV0ooVK/Tdd9+pevXqeu211wp8vctvrOvYsaNq1qypt956Sxs3blTjxo21a9cuLViwQB06dNDMmTM98j6qVaum8ePHa8CAAYqNjVXXrl1Vvnx5LVu2TN9//73q1auXbdopT2jYsKHWrFmjJk2a6M4779Tx48f1+eefKzk5WW+++aaio6M9+nqS9Le//S3zZseVK1fKNE298cYb+vjjjyVZ/9zeuXPnq15nyJAhevLJJxUTE6P7779f/v7++u6777R582Z16tRJ8+bNK3St8fHxmjVrlubNm6cGDRqoc+fOSktL0xdffKHmzZvn+MtIfrRq1UpTpkzRmjVr1LJly8zfLdM0NXny5FxvigNQxGya1eGKTpw4Yb788svmXXfdZUZERJiGYVx16p+0tDSzbt26pmEY5j/+8Y9rVClgH3enHerTp49ZoUIFMygoyGzcuLE5ffp0c+nSpaZhGLlODfXdd9+ZrVu3zpyCrGPHjuYvv/yS57RgLpfLvO2227Ltc/fcjz76yGzcuLEZFBRkVqhQwezTp4954MCBHDW2adPG9PPzy/N7ldcUYnv27DG7d+9uXn/99WZQUJDZoEEDc+zYseb58+dzre1KrlbDwoULzXbt2pkRERFmyZIlzVq1aplDhw41k5OTC3ytK8mY/mr//v1mjx49zIoVK5qlSpUymzZtas6cOdOta+ZHVFRUoaYdyyrjcw8ODjbLly9v3n///eavv/5aoJ+jDD179jT9/PzMXbt2ZduflpZmjho1yqxRo4ZZsmRJMzo62hw2bJh57ty5An32GdOB9erVy9yyZYvZpUsXMyIiwixdurQZFxdnLly4MN/v2xMMwzCrVq16TV8T8FaGaXrRvEEX7dq1S9HR0apWrZqqV6+upUuXaurUqbk22Gd46623NHz4cJ06dUpvvvmmnn766WtYMeDbli1bprZt22rEiBFFPj0UipbL5VKbNm20ePFiu0spFjL+f9WzZ099+OGHdpcDIAuv7MmtXLmyDhw4oB07duiNN9646vydhw4d0qhRozR06FCvmusTAAAA9vDKkOvv75/vNdola7nKunXr5riZAgAAAMWTz9949sMPP2j69OlatWqVx+6GBoojwzD4HXIAPsdrj+854J18PuQOHDhQ8fHxat68uXbt2mV3OYBPat26dbYVoeC7+ByvrWrVqvE9B7yUT4fcqVOnatOmTZo9e7bdpQAAAMCL+GzITU1N1QsvvKBnn31WlStXLtBzjxw5om+++UZRUVEqVapUEVUIAAAAd50+fVo7d+5Uu3btVK5cuQI/32dD7ptvvqm0tDR17do1s01h9+7dkqTjx49r165dqly5cq6rCH3zzTc5lgsFAACA9/nkk0/cmlzAZ0Pu7t27dfz4cdWrVy/bfsMwNHr0aI0ZM0Y//fRT5mpPWUVFRUmyvmmeXC8e3mvIkCEaN26c3WXgGuHzLl74vIsXPu/iY/PmzerevXtmbisonw25gwcPzrYWuGTNl9uvXz/16tVLXbp0yXPpyowWhbp16yo2NrbIa4X9wsLC+KyLET7v4oXPu3jh8y5+3G0t9dqQ++677yopKUl79+6VJM2dOzezHWHQoEGKiYlRTExMtudktC3Ur19fHTt2vLYFAwAAwGt4bcgdO3asEhMTJVktCLNnz86cRaFHjx4KCQnJ9XnMVQgAAACvDbk7duwo8HOYrxAAAACSly7rC3hafHy83SXgGuLzLl74vIsXPm/kFyEXxQL/USxe+LyLFz7v4oXPG/lFyAUAAIDjEHIBAADgOIRcAAAAOA4hFwAAAI5DyAUAAIDjEHIBAADgOIRcAAAAOA4hFwAAAI5DyAUAAIDjEHIBAADgOIRcAAAAOA4hFwAAAI5DyAUAAIDjEHIBAADgOIRcAAAAOA4hFwAAAI5DyAUAAIDjEHIBAADgOIRcAAAAOA4hFwAAAI5DyAUAAIDjEHIBAADgOIRcAAAAOA4hFwAAAI5DyAUAAIDjEHIBAADgOIRcAAAAOA4hFwAAAI5DyAUAAIDjEHIBAADgOIRcAAAAOA4hFwAAAI5DyAUAAIDjEHIBAADgOIRcAAAAOA4hFwAAAI5DyAUAAIDjEHIBAADgOMU+5D70kPTWW3ZXAQAAAE8qYXcBdtu8WapQwe4qAAAA4EnFfiQ3MFA6d87uKgAAAOBJxT7kBgQQcgEAAJyGkBsgnT1rdxUAAADwpGIfcmlXAAAAcJ5iH3JpVwAAAHAeQi7tCgAAAI5DyGUkFwAAwHGKfcilJxcAAMB5in3IpV0BAADAeQi5tCsAAAA4TrEPubQrAAAAOE+xD7mM5AIAADgPIZeeXAAAAMfxypB78uRJDR8+XO3bt1fZsmXlcrk0ffr0bOeYpqlp06apc+fOqlq1qoKDg3XjjTdq9OjROluA1Eq7AgAAgPN4Zcg9cuSIRo0apS1btigmJkaGYeQ459SpU+rdu7eOHDmiJ598Um+//bZuuukmDR8+XHfffXe+X4t2BQAAAOcpYXcBualcubIOHDigChUqaN26dWrWrFmOcwICArRq1Sq1aNEic1+fPn1UrVo1jRgxQosXL9Ztt9121deiXQEAAMB5vHIk19/fXxUqVLjqOVkDboZ7771Xpmlq8+bN+XqtjHYF03SrVAAAAHghrwy5hbF//35JUrly5fJ1fkCAFXDT04uyKgAAAFxLjgu5b7zxhsLCwtS+fft8nR8QYD3SsgAAAOAcXtmT664xY8Zo8eLFmjhxokJDQ/P1nIyQe+6cVLp0ERYHAACAa8YxIXfWrFkaNmyY+vbtq379+uXrOUOGDNGZM2GSpG7drP7c+Ph4xcfHF2WpAAAAyCIhIUEJCQnZ9iUnJxfqmo4IuQsXLtRjjz2mjh07auLEifl+3rhx43T4cKzuukv64AOpSpUiLBIAAAC5ym2Qcf369WrSpInb1/T5ntw1a9bovvvuU/PmzTVr1iy5XAV7S/TkAgAAOI9Ph9zNmzerQ4cOql69uubNm6fAwMACXyPjKSwIAQAA4Bxe267w7rvvKikpSXv37pUkzZ07V7t375YkDRo0SIZhqF27dkpKStKzzz6r+fPnZ3t+jRo1cp1H93JZbzwDAACAM3htyB07dqwSExMlSYZhaPbs2Zo9e7YkqUePHjJNMzMADx06NMfzH3vssQKFXNoVAAAAnMNrQ+6OHTuuek66B1ZwoF0BAADAeXy6J9cTaFcAAABwHkIuIRcAAMBxCLn05AIAADhOsQ+59OQCAAA4T7EPubQrAAAAOA8hl3YFAAAAxyn2IdflkkqUYCQXAADASYp9yJWs0VxCLgAAgHMQcmWFXNoVAAAAnIOQK2uGBUZyAQAAnIOQK9oVAAAAnIaQK0IuAACA0xByRU8uAACA0xByRU8uAACA0xByRbsCAACA0xByRbsCAACA0xByRbsCAACA0xByRbsCAACA0xByRcgFAABwGkKurHYFenIBAACcg5ArRnIBAACchpArQi4AAIDTEHJFuwIAAIDTEHLFSC4AAIDTEHJFyAUAAHAaQq5Y8QwAAMBpCLlixTMAAACnIeSKdgUAAACnIeSKkAsAAOA0hFwxhRgAAIDTEHJljeSmpUmmaXclAAAA8ARCrqyQK1lBFwAAAL6PkCurXUGiZQEAAMApCLm6NJLLzWcAAADOQMgVIRcAAMBpCLm6FHJpVwAAAHAGQq4u9eQykgsAAOAMhFzRrgAAAOA0hFwRcgEAAJyGkCumEAMAAHAaQq4YyQUAAHAaQq4IuQAAAE5DyBXtCgAAAE5DyBUjuQAAAE5DyBUhFwAAwGkIuSLkAgAAOA0hV5K/v/VITy4AAIAzEHIlGYY1mstILgAAgDMQci8i5AIAADgHIfeiwEDaFQAAAJyCkHsRI7kAAADOQci9iJALAADgHITci2hXAAAAcA5C7kWM5AIAADiHV4bckydPavjw4Wrfvr3Kli0rl8ul6dOn53ruli1bdNdddykkJERly5bVo48+qiNHjhT4NQm5AAAAzuGVIffIkSMaNWqUtmzZopiYGBmGket5e/fuVVxcnLZv366///3veuaZZ7RgwQLdeeedOn/+fIFek5ALAADgHCXsLiA3lStX1oEDB1ShQgWtW7dOzZo1y/W80aNH6/Tp09qwYYOuu+46SVKzZs10xx13aNq0aerbt2++X5OeXAAAAOfwypFcf39/VahQ4arnffXVV+rQoUNmwJWkP/3pT6pdu7Y+++yzAr0mI7kAAADO4ZUhNz/27dunQ4cOqWnTpjmONW/eXD/99FOBrkfIBQAAcA6fDbn79++XJFWqVCnHsUqVKunYsWNKS0vL9/VoVwAAAHAOnw25p0+fliQFBgbmOFayZMls5+QHI7kAAADO4bMht1SpUpKks7kMv545cybbOflByAUAAHAOr5xdIT8y2hQy2hay2r9/vyIiIuTv73/FawwZMkRhYWGSpJ9/llJSpISEeMXHx3u+YAAAAOQqISFBCQkJ2fYlJycX6po+G3IrV66s8uXLa+3atTmO/fDDD4qJibnqNcaNG6fY2FhJ0oAB0sqVEvkWAADg2oqPzznIuH79ejVp0sTta/psu4Ik3X///Zo/f7727t2buW/RokX63//+p65duxboWrQrAAAAOIfXjuS+++67SkpKygywc+fO1e7duyVJgwYNUkhIiF544QV98cUXatOmjQYPHqzU1FSNHTtWjRo1Us+ePQv0eoGBhFwAAACn8NqQO3bsWCUmJkqSDMPQ7NmzNXv2bElSjx49FBISouuvv17Lli3T008/reeff14BAQHq0KGDxo4de9V+3MsFBDCFGAAAgFN4bcjdsWNHvs6rW7euvv7660K/Hu0KAAAAzuHTPbmeRMgFAABwDkLuRax4BgAA4ByE3IsYyQUAAHAOQu5FAQFSerq1AQAAwLcRci8KDLQe09LsrQMAAACFR8i9KCDAeqQvFwAAwPcRci/KCLn05QIAAPg+Qu5FGe0KhFwAAADfR8i9iHYFAAAA5yDkXkS7AgAAgHMQci8i5AIAADgHIfeijJ5c2hUAAAB8HyH3IkZyAQAAnIOQexEhFwAAwDkIuRcxhRgAAIBzEHIvYgoxAAAA5yDkXkS7AgAAgHMQci+iXQEAAMA5CLkX0a4AAADgHITci/z8JJeLkVwAAAAnIORmERBAyAUAAHACQm4WgYGEXAAAACcg5GYREEBPLgAAgBMQcrOgXQEAAMAZCLlZ0K4AAADgDITcLGhXAAAAcAZCbha0KwAAADgDITcL2hUAAACcgZCbBe0KAAAAzkDIzYJ2BQAAAGcg5GYRGMhILgAAgBMQcrMICZFSU+2uAgAAAIVFyM0iNFRKSbG7CgAAABQWITcLQi4AAIAzEHKzIOQCAAA4AyE3i9BQKTnZ7ioAAABQWITcLMLCpBMnpPR0uysBAABAYRByswgNtR5PnLC3DgAAABQOITeLjJBLXy4AAIBvI+RmQcgFAABwBkJuFhkhl5vPAAAAfBshN4uwMOuRkVwAAADfRsjNgnYFAAAAZyDkZhEcbD0ScgEAAHwbITcLl0sKCSHkAgAA+DpC7mVY9QwAAMD3EXIvExbGSC4AAICvI+ReJjSUkAsAAODrCLmXIeQCAAD4PkLuZQi5AAAAvo+QexluPAMAAPB9hNzLcOMZAACA7yPkXoZ2BQAAAN9HyL0MIRcAAMD3EXIvExoqpaZKFy7YXQkAAADc5fMhd+vWrerWrZuqVKmi0qVLq27duho1apROnz7t1vVCQ63H1FQPFgkAAIBrqoTdBRTGnj171KxZM4WHh2vgwIGKiIjQ999/r+HDh2v9+vWaPXt2ga8ZFmY9pqRc+jMAAAB8i0dD7qlTpzRz5kydPXtWd999t6pVq+bJy+cwffp0paSk6Pvvv9cNN9wgSerbt6/S09P18ccfKzk5WWEFTKoZI7n05QIAAPgut0Nunz59tGbNGv3666+SpHPnzqlFixaZX4eFhWnx4sVq3LixZyrNRerFnoIKFSpk2x8ZGSmXy6WAgIACX5OQCwAA4Pvc7sldsmSJ7rvvvsyvZ8yYoV9//VWffvqpfv31V0VGRmrkyJEeKTIvbdq0kWma6t27t37++Wft2bNHs2bN0qRJkzR48GCVKlWqwNck5AIAAPg+t0PugQMHFBUVlfn1nDlz1LRpU8XHx6tevXp6/PHHtWbNGk/UmKd27dpp1KhRWrhwoRo3bqyqVavq4Ycf1qBBgzR27Fi3rpkRcln1DAAAwHe53a5QunRpJSUlSZLOnz+vpUuXauDAgZnHQ0JClHwNkmJUVJRat26tBx54QBEREVqwYIFGjx6tyMhIPfXUUwW+XnCw9chILgAAgO9yO+TGxsbq/fffV9u2bTV37lylpqaqY8eOmce3bdumihUreqTIvMycOVP9+vXT1q1bValSJUlSly5dlJ6erueee07x8fEKDw/P8/lDhgzJcWNafHy8QkLiCbkAAADXSEJCghISErLtK+xgqWGapunOE9euXat27dopKSlJpmnqgQce0GeffZZ5vE6dOmrWrJk++eSTQhV4Ja1bt9aFCxe0YsWKbPvnzJmj+++/XwsXLtRtt92W43nr169XkyZNtG7dOsXGxuY4fv31Ut++0ogRRVU5AAAAruRqee1q3B7Jbdq0qbZs2aLFM1N1AAAgAElEQVRVq1apTJkyat26deaxpKQkPfXUU9n2FYWDBw8qIiIix/60tDRJVhuFO1jaFwAAwLcVap7c8uXLq3Pnzjn2lylTRoMHDy7MpfOldu3aWrhwobZu3aqaNWtm7p8xY4ZcLpcaNmzo1nVDQ7nxDAAAwJcValnf9PR0zZw5U0888YTuvfde/fLLL5KsHoqvvvpKBw8e9EiReXnmmWeUnp6uW2+9Va+++qomTpyou+++W3PnzlWfPn0UGRnp1nUZyQUAAPBtbofcpKQktWzZUg8//LASEhI0d+5cHT58WJIUHBysQYMG6e233/ZYobmJi4vTqlWr1LRpU02cOFFDhgzRjh07NGbMGL333ntuXzcsjJALAADgy9xuVxg6dKg2bdqkb775Ro0bN8626pifn58eeOAB/fvf/9aYMWM8UmhemjZtqvnz53v0mqGh0p49Hr0kAAAAriG3R3LnzJmjgQMH6o477pBhGDmO165dWzt37ixMbbahXQEAAMC3uR1yk5OTFR0dnefxtLQ0t2c3sBs3ngEAAPg2t0NujRo1tH79+jyP//e//1W9evXcvbytGMkFAADwbW6H3L59++rDDz/UrFmzlLGehGEYOnv2rF588UX95z//0RNPPOGxQq+lsDApNVW6cMHuSgAAAOAOt288Gzx4sDZt2qT4+HiVKVNGkvTwww/r6NGjOn/+vJ544gn16dPHY4VeS6Gh1uOJE5f+DAAAAN/hdsg1DEPvv/++HnvsMX3xxRf6448/dOHCBdWoUUNdu3ZVq1atPFnnNZURbFNSCLkAAAC+qFArnknSrbfeqltvvdUTtXiNjGCbnCxdf729tQAAAKDgCrXimVNlHckFAACA78n3SG50dLRcLpe2bNkif39/RUdH5zo/blaGYWjbtm2FLvJaCwuzHgm5AAAAvinfIbd169YyDEMulyvb107ESC4AAIBvy3fInTZt2hW/dpLgYOuRBSEAAAB8k9s9udOnT7/isr27du3S9OnT3b28rfz8rKDLSC4AAIBvcjvk9urVS6tWrcrz+OrVq9WrVy93L287Vj0DAADwXW6H3IxVzvJy8uRJlShR6BnKbBMWRsgFAADwVQVKoRs3btSGDRsyv16xYoXOnz+f47ykpCRNmjRJtWvXLnyFNmEkFwAAwHcVKOTOnj1bI0eOlGRNDzZ58mRNnjw513PLlCnjsz25khVyufEMAADANxUo5Pbr108dOnSQaZpq3ry5XnnlFbVv3z7bOYZhqHTp0qpRo4ZPtyswkgsAAOC7CpRCK1WqpEqVKkmSlixZorp166pChQpFUpjdQkOlPXvsrgIAAADucHuotXXr1p6sw+uULy8dPmx3FQAAAHBHvkNu27Zt5XK59M0336hEiRK67bbbrvocwzC0aNGiQhVol8hI6cAByTQlhy7sBgAA4Fj5nkLMNE1duHAh8+sLFy7INM0rblnP9zWRkdKpU9KJE3ZXAgAAgILK90ju0qVLr/i100RGWo8HD0ohIfbWAgAAgIJxezGI5cuX6/AVmlaPHDmi5cuXu3t522WE3AMH7K0DAAAABed2yG3btq0WLlyY5/FFixapbdu27l7edoRcAAAA31Vky/qePXtWfn5+7l7edmXKSAEBhFwAAABfVKApxBITE7Vz587Mr7ds2ZJrS0JSUpImT56satWqFbpAuxiGVLEiIRcAAMAXFSjkTp06VSNHjpRhGDIMQ6NHj9bo0aNznGeapvz8/PJc8tdXZEwjBgAAAN9SoJDbtWtXNWjQQKZpqmvXrho0aJDi4uKynZOxrG9MTIwqVqzo0WKvtchIa3YFAAAA+JYChdy6deuqbt26kqxR3VatWik6OrpICvMGkZHSTz/ZXQUAAAAKyu1lfR977DFP1uGVaFcAAADwTW6HXEnavHmzpk6dqu3bt+v48eM5Zlzw5WV9pUvtChcuSC6356EAAADAteZ2yP3444/Vq1cv+fv7q06dOgoPD89xztWmGfN2kZFSWpp0/LhUtqzd1QAAACC/3A65I0aMUOPGjfX111+rXLlynqzJa2TcN3fgACEXAADAl7j9j/D79u1T7969HRtwJVY9AwAA8FVuh9yGDRtq3759nqzF62SM5DKNGAAAgG9xO+S+9dZbmjJlilatWuXJerxKcLC1MZILAADgW9zuyX399dcVFhamuLg41atXT1WrVpWfn1+2cwzD0L/+9a9CF2knphEDAADwPW6H3I0bN8owDFWtWlUnTpzQb7/9luMcwzAKVZw3IOQCAAD4HrdD7s6dOz1Yhvci5AIAAPgelji4iooVCbkAAAC+plArnmVITU1VcnKyLly4kONY1apVPfEStslY9QwAAAC+o1Ahd+LEiXrrrbe0ffv2PM9JT08vzEvYLjJSOnxYOn9eKuGRvxIAAACgqLndrjBp0iT9+c9/Vs2aNfXqq6/KNE395S9/0dChQxUZGalGjRppypQpnqzVFpGRkmlaQRcAAAC+we2Q+84776hdu3b6+uuv1a9fP0nSPffco9GjR+u3335Tamqqjh496rFC7cKqZwAAAL7H7ZC7bds2dezYUZLk7+8vSTp37pwkKSwsTH379tV7773ngRLtRcgFAADwPW6H3LCwMJ0/f16SFBoaqqCgIO3evTvzeEhIiA44IBlWqGA9OuCtAAAAFBtuh9wGDRro559/zvy6RYsWmjhxovbu3avdu3dr8uTJql27tkeKtFNAgBQRwQwLAAAAvsTt+QK6d++uSZMm6ezZswoMDNTIkSN1++23Z04Z5u/vry+//NJjhdqJBSEAAAB8i9sht1evXurVq1fm1y1bttSmTZs0b948+fn56c4773TESK5EyAUAAPA1Hp35tXr16ho8eLAnL+kVIiOlvXvtrgIAAAD5xbK++cBILgAAgG9xeyTX5XLJMIyrnufrK55J0nXXSXv2WItC5OMtAwAAwGZuh9yXX345R8hNT0/Xzp07NWfOHNWpU0cdOnQodIHeoHp16eRJ6dAhqWJFu6sBAADA1bgdckeMGJHnsf3796tFixaOufGsenXrcds2Qi4AAIAvKJKe3EqVKql///4aNWpUUVw+h/Xr16tTp04qW7asSpcurRtvvFH/93//57HrZw25AAAA8H4enV0hq9KlS2vHjh1FdflM//3vf9WpUyfFxsbq5ZdfVnBwsLZt26Y9e/Z47DWCg60RXEIuAACAbyiSkPvrr79qwoQJRd6ukJqaqscee0wdO3bU559/XqSvVaMGIRcAAMBXuB1yo6Ojc51dISkpScnJyQoKCtKcOXMKVdzVfPrppzp06JBGjx4tSTp16pRKlSqVr1kfCqpGDWnrVo9fFgAAAEXA7ZDbunXrHGHSMAyFh4erRo0a6tatmyIiIgpd4JUsWrRIoaGh2r17tzp16qT//e9/Kl26tHr06KFx48YpMDDQY69Vo4b0zTceuxwAAACKkNshd9q0aZKkkydPKjU1VeXKlVOJEkXW4purP/74Q2lpaercubMef/xx/f3vf9fSpUs1YcIEJScn69NPP/XYa9WoYU0hlpoqhYR47LIAAAAoAm7NrrBr1y4NGDBA1apVU2hoqK677jqVLFlSUVFRevbZZ7Vr1y5P15mrEydO6PTp0+rZs6fGjRunLl26aPz48XriiSc0c+ZMbfNgE22NGtbj9u0euyQAAACKSIFD7rx589SwYUO999578vPzU8eOHfXwww+rQ4cOcrlcGjt2rGJiYrRgwYLM57z00kseLTpDqVKlJEndunXLtv/hhx+WaZr6/vvvPfZaGSGXm88AAAC8X4H6CzZv3qyuXbsqOjpakydPVlxcXI5zVqxYof79++uhhx7S2rVr9dprr+mTTz7Rq6++6rGiM1SuXFm//fabKl62QkOFChUkScePH7/i84cMGaKwsLBs++Lj4xUfH5/j3PLlranECLkAAACelZCQoISEhGz7kpOTC3XNAoXcMWPGqGzZslq5cmWeN5XFxcVpxYoVatiwoZo0aaKzZ8/qtddeK1SReWnSpIm+/fZb7d27V7Vq1crcv2/fPklS+fLlr/j8cePGKTY2Nl+vZRhMIwYAAFAUchtkXL9+vZo0aeL2NQvUrrB48WL16dPnqrMmREREqHfv3jp9+rSmTZumZ5991u0Cr6Rr164yTVNTpkzJtv/999+Xv7+/2rRp49HXq16dkAsAAOALCjSSe/ToUUVFReXr3OjoaPn5+al79+7u1JUvMTEx6t27t6ZOnaq0tDS1bt1aS5Ys0ZdffqkXXnhBkZGRHn29GjWkL7/06CUBAABQBAoUcsuVK5fvpXp37NiR2RtblCZPnqxq1app6tSpmjNnjqpVq6bx48dr4MCBHn+tGjWkxEQpLU3y9/f45QEAAOAhBWpXaNOmjaZMmaJjx45d8bxjx45pypQpuu222wpVXH74+flp2LBh2r59u86cOaPff/+9SAKuZIXc9HQr6AIAAMB7FSjkvvDCCzp69KhatWqlVatW5XrOqlWr1Lp1ax09elTPP/+8R4r0FkwjBgAA4BsK1K5Qr149zZgxQ48++qji4uIUFRWlRo0aKSQkRKmpqdq4caN27NihUqVKacaMGapXr15R1W2LqlWlEiUIuQAAAN6uwOvw3nfffYqJidEbb7yh+fPna86cOZnHKleurH79+ulvf/ubamQMezpIiRJStWqEXAAAAG9X4JArSdWrV9ekSZMkSSkpKUpNTVVISIhCQ0M9Wpw3Yq5cAAAA7+dWyM0qNDS0WITbDDVqSCtX2l0FAAAArqRAN57h0kiuadpdCQAAAPJCyC2gOnWkU6eYRgwAAMCbEXILqFEj6/Hnn+2tAwAAAHkj5BbQ9ddLERHShg12VwIAAIC8EHILyDCs0VxGcgEAALwXIdcNhFwAAADvRsh1Q0yMNcNCSordlQAAACA3hFw3ZNx89ssv9tYBAACA3BFy3VCvnuTvT8sCAACAtyLkuiEgQKpblxkWAAAAvBUh103cfAYAAOC9CLluiomxenLT0+2uBAAAAJcj5LqpUSPp9Glp61a7KwEAAMDlCLluyphhgb5cAAAA70PIdVO5ctJ119GXCwAA4I0IuYXAzWcAAADeiZBbCI0a0a4AAADgjQi5hRATI+3bJx0+bHclAAAAyIqQWwjNmlmP339vbx0AAADIjpBbCFFR1s1nK1faXQkAAACyIuQWgmFIcXHSihV2VwIAAICsCLmFFBcnrV0rnTpldyUAAADIQMgtpLg46fx5ac0auysBAABABkJuIdWvL5UpQ8sCAACANyHkFpLLJbVsScgFAADwJoRcD4iLs6YRO3/e7koAAAAgEXI9Ii5OOnmS1c8AAAC8BSHXA5o2lUqWpGUBAADAWxByPSAgQLrpJkIuAACAtyDkesitt1orn5mm3ZUAAACAkOshcXHS4cPS77/bXQkAAAAIuR7SsqXVtvDf/9pdCQAAAAi5HhIcLLVpI82fb3clAAAAIOR60D33SMuWSampdlcCAABQvBFyPeiee6Rz56Rvv7W7EgAAgOKNkOtBNWpIN9xAywIAAIDdCLke1qGD9O9/Sxcu2F0JAABA8UXI9bAOHaQDB6T16+2uBAAAoPgi5HrYLbdIYWHSggV2VwIAAFB8EXI9zN9fuusu+nIBAADsRMgtAvfcI61da7UtAAAA4Noj5BaB9u0ll0uaO9fuSgAAAIonQm4RKFdO+tOfpBkz7K4EAACgeCLkFpFHHrFWP9u92+5KAAAAih9CbhG5916pZEkpIcHuSgAAAIofQm4RCQ2VOnWSPv3U7koAAACKH0JuEereXdq4Ufr1V7srAQAAKF4IuUWoXTspIoLRXAAAgGuNkFuEAgKkrl2tkHvhgt3VAAAAFB+E3CL2yCPWDAsrV9pdCQAAQPHhqJA7evRouVwuNWzY0O5SMt1yixQdLU2ZYnclAAAAxYdjQu7evXv12muvKTg42O5SsnG5pCeflGbNkg4ftrsaAACA4sExIfevf/2rbr75ZjVp0sTuUnLo3VsyDEZzAQAArhVHhNzly5frq6++0vjx4+0uJVdly0rx8dLEiVJ6ut3VAAAAOJ/Ph9wLFy5o0KBBevzxx1W/fn27y8nTgAFSYqI0f77dlQAAADhfCbsLKKyJEycqMTFRixcvtruUK4qNlVq0kP7v/6TOne2uBgAAwNl8eiT32LFjGj58uF5++WVFRETYXc5VDRggffuttHmz3ZUAAAA4m0+P5L744osqW7asBgwY4NbzhwwZorCwsGz74uPjFR8f74nycnjgAenpp6UJE6z+XAAAAEgJCQlKSEjIti85OblQ1/TZkLt161a9//77evvtt7V3715JkmmaOnPmjNLS0rRr1y6FhoYqPDw8z2uMGzdOsbGx16pkBQZKgwZJo0ZJw4ZJlStfs5cGAADwWrkNMq5fv75Qs2b5bLvC3r17ZZqmBg0apOjoaEVHR6t69epas2aNfv/9d1WvXl2jRo2yu8wcBgyQSpWS3nzT7koAAACcy2dHchs0aKDZs2fn2P/iiy/qxIkTmjBhgqpXr25DZVcWFiYNHiy9/ro0dKhUsaLdFQEAADiPz4bcsmXLqlOnTjn2jxs3ToZhqGPHjjZUlT+DB0tvvSX94x/SG2/YXQ0AAIDz+Gy7wpUYhmF3CVcUHm61Lbz3nnTkiN3VAAAAOI/jQu6SJUv0888/213GVT39tPU4bpy9dQAAADiR40KuryhXzhrNHTdOujg5BAAAADyEkGuj55+XSpeWXnrJ7koAAACchZBro7AwaeRI6aOPpPXr7a4GAADAOQi5NuvXT7rhBumvf5VM0+5qAAAAnIGQa7MSJayFIZYulebNs7saAAAAZyDkeoG775Zuv92aceH0aburAQAA8H2EXC9gGNKECdLu3dKrr9pdDQAAgO8j5HqJunWlF16wVkDbuNHuagAAAHwbIdeLDB0q1aolPf64lJ5udzUAAAC+i5DrRQIDpfffl374QXr3XburAQAA8F2EXC/TsqX05JNW68Iff9hdDQAAgG8i5HqhN96QKlWSHnlESkuzuxoAAADfQ8j1QsHB0qefSj/9JI0YYXc1AAAAvoeQ66WaN7eW/H3tNWnZMrurAQAA8C2EXC/23HNSq1ZSjx7SkSN2VwMAAOA7CLlezM9P+vhjaxW0bt2k8+ftrggAAMA3EHK9XJUq0mefSUuXWjMuAAAA4OoIuT6gbVtrxoU337QCLwAAAK6MkOsjhgyR4uOlXr2sWRcAAACQN0KujzAM6YMPpHr1pLvvlnbtsrsiAAAA70XI9SFBQdL8+VKpUlL79tLx43ZXBAAA4J0IuT6mYkXp66+lgwele++Vzp61uyIAAADvQ8j1QXXqSPPmSWvWSA89xNK/AAAAlyPk+qhbbpFmz5b+/W+pe3fm0AUAAMiKkOvD7rrLmlLsyy+l3r2lCxfsrggAAMA7EHJ9XJcu0iefSJ9+ak0vxoguAACAVMLuAlB43bpJpin16CGdPCnNmCEFBNhdFQAAgH0YyXWI+Hjpq6+sG9K6dJFOn7a7IgAAAPsQch2kUydpwQJp2TLp9tulI0fsrggAAMAehFyHuf12ackS6Y8/rBkYtm2zuyIAAIBrj5DrQM2bS6tXW0sBt2ghff+93RUBAABcW4Rch6peXVq1SqpbV2rTRpoyxe6KAAAArh1CroOVLSt9+63Us6fUt680YACrowEAgOKBkOtwAQHS5MnSpEnSP/8p3XabtGeP3VUBAAAULUJuMfHEE9LSpdKOHVJMjPT113ZXBAAAUHQIucXILbdIGzZIN90k3X239Nxz0rlzdlcFAADgeYTcYqZcOWvBiNdfl956y5p94bff7K4KAADAswi5xZDLJT37rLRmjXTmjBQbK40fL6Wn210ZAACAZxByi7HYWGndOql/f2nIEKllS2njRrurAgAAKDxCbjFXqpQ1irtypZSaKjVpIg0dKp06ZXdlAAAA7iPkQpI1ivvTT9KIEVbobdBA+u9/7a4KAADAPYRcZAoIkF580WpZiIqS2rWTHnlE2rfP7soAAAAKhpCLHGrXlhYtkqZNs0Zza9WShg+XTpywuzIAAID8IeQiV4YhPfaYtHWrtRzw669bYfeDD5iFAQAAeD9CLq4oLMwKuFu2SG3bSo8/bq2Y9s03dlcGAACQN0Iu8iUqSpoxw5pbNzxcuusuK/QuXWp3ZQAAADkRclEgzZtLy5ZJ//qXlJxsBd3Wra0eXtO0uzoAAAALIRcFZhhSp07WQhJz50onT0q33y7FxVk3qhF2AQCA3Qi5cJthSB07Sj/+KC1YIKWlWdOO3XSTNGuWdP683RUCAIDiipCLQjMM6e67pdWrpf/8RwoNlbp1k2rUkMaOlZKS7K4QAAAUN4RceIxhWCO5334rbdhg9eu+8IJUpYr0l79I27fbXSEAACguCLkoEo0aWYtJ7NolDR4sffyxVLOmNSvD7Nm0MgAAgKJFyEWRqlRJevVVafduayGJ48el++6TqlaVhg2TEhPtrhAAADgRIRfXRFCQ1Lu3Nc/uTz9JnTtL48db8+/ec4/0xRfSmTN2VwkAAJzCp0Pu2rVrNWDAADVo0EDBwcGqVq2aHnroIf3xxx92l4YriImRJk6U9u+XJk+WjhyRHnzQGvXt31/67jumIQMAAIXj0yH39ddf1+zZs3X77bdrwoQJeuKJJ7R8+XLFxsbqt99+s7s8XEVwsLVM8Jo10ubN0lNPSV9/Ld16q1SrljRypPS//9ldJQAA8EWGafrumNnq1avVtGlTlShRInPf1q1bdeONN+rBBx/U9OnTc33e+vXr1aRJE61bt06xsbHXqlzkw4UL0vLl0vTpVgtDaqp1E1vXrtZob61adlcIAACuhcLmNZ8eyW3RokW2gCtJNWvWVP369bV582abqkJhuFxSmzbShx9KBw9KX30l1a0rjRkj1a4tNW4svfYaI7wAAODKfDrk5uXgwYMqV66c3WWgkEqVku69V0pIkA4dskZ2a9e2ZmuoU8fannnGGvllSjIAAJCV40LuJ598or1796pbt252lwIPCgqS7r/fWi748GHpX/+SWrWy5t9t3VqqWFHq0UP6/HMpJcXuagEAgN1KXP0U37FlyxYNGDBALVu21KOPPmp3OSgiQUFSp07WduGC9OOP0rx51vbJJ5K/v9Xy0KGDdOed1oivYdhdNQAAuJYcE3IPHjyoe+65R+Hh4fr8889l5CfV0Lfr81ySbvKXbrpPevU+ad8+acUKadky6ZO/StPOSxUrSDfddGmLiLC7agAAcFWFzGmOCLkpKSm66667lJKSopUrVyoyMjJfzxvSvbvCLtsXf3GDb6os6aGLW6ZDkuZd3AAAgNdJuLhllVzIa/p8yD179qw6dOigrVu3atGiRapTp06+nzvuk08UW7duEVYHb3L4sPTDD9Lq1dbj4SNSgL/UoIEUGys1aSI1bGjd8AYAAK6d3AYZ12/erCbdu7t9TZ8OuRcuXFDXrl21Zs0azZ07V82bNy/YBerWtdINioXyku5pJ90ja0W1TZukb7+1Zmd4eY50dIpUooQVdlu1sraWLaXwcLsrBwAABeXTIffpp5/WvHnz1KlTJx05ckSffvpptuOPPPKITZXB2xmGNYLboIH0l79YoXfzZivwLl8uzZghvfmmdV7DhlbYvflmqUULqUYNbmQDAMDb+XTI/fnnn2UYhubNm6d583I2XBJykV+GIdWrZ239+1uhd+dOK/AuWyYtXiy99551btmyVtjN2Jo1k8Iub+4GAAC28umQu2TJErtLgEMZhhQdbW2PPWbtO378Uk/v6tXSP/4hJSVdCsg33WS1OsTGWksR09sLAIB9fDrkAtdSeLjUrp21SdYcvf/736XQu2aNtThFWprk53ep5Tsj+MbESMHB9r4HAACKC0Iu4CaXS7rhBmvr2dPad/as9Ouv0vr11rZunbVK29mz1ohvnTrZg2+jRtzYBgBAUSDkAh4UGGgF2CZNLu1LS7Nualu37lL4nTNHOnXKOl65snTjjdm3unWlkiXteQ8AADgBIRcoYv7+1gwNDRtKvXpZ+9LTpd9/lzZulH75xRr9/fJLaexY67jLJdWqlTP8Vq9uHQMAAFdGyAVs4Od3aTaHbt0u7U9NtebvzQi+v/wivfOOdOSIdbxUKal2batFom7dS+0StWtzoxsAAFkRcgEvEhJyaWqyDKYpHTxoBd5Nm6QtW6xtyRLp0CHrHMOQoqIuhd6sIbhcOeb1BQAUP4RcwMsZhhQZaW133JH92LFjVtvDli1W3++WLdL8+dLbb1uzP0jWHL41a2bfatWyHitUIAADAJyJkAv4sIgIayW2m2/Ovv/sWWnrViv4bt16aVu5Utq799J5wcE5A3DGVqkS/b8AAN9FyAUcKDBQql/f2i536pS0fbu0bVv2ADxzprRrl9UeIVk9vtWrW20Ql2/R0VbAZhQYAOCtCLlAMRMUJDVoYG2XO3vWWs44a/jdtcta3vijj6QTJy6dGxycewDO2AjBAAA7EXIBZAoMtBasqFMn5zHTtJY23rHDCsJZt6VLrf0nT146PzhYqlpVuv56qUqVnI9Vqlg32gEAUBQIuQDyxTCs0dmIiOyLXWQwTetGuIzgu2OHtHu3tGePNR/wv/8tHThwqR1CkkJDcw/AWR8JwgAAdxByAXiEYUhly1pbbiFYks6dk/bvvxR+sz7+/LM1M8TBg9mDcFiYFXYrV7ZuhqtUKfc/M08wACArQi6AayYgQKpWzdrycu6ctG9fzhC8f7/VI7xihfXnM2eyPy8s7OpBuFIlq40CAOB8hFwAXiUg4NLNa3kxTSkpyQq7+/dboTjrnxMTpTVrrD+fOpX9ucHB1vzAFSteesz656yP4eHcPAcAvoqQC8DnGIYVQMPDraWR82Ka1lLJlwfhQ4es7eBBad066/HQIWsUOasSJXIPv+6RxUIAABLfSURBVLk9li8v+fsX7fsGAOQfIReAYxmGdXNbaGjuM0ZkZZpScvKl8Jvb47Zt0vffW1+npOS8RliYtYxyfrfwcMnPr2jeOwAUd4RcAJAViMuUsbbata9+/pkz2UeEDx+WjhzJvm3dKq1ebf352LHsN9RlvGZERN4hOONGvoxZLcLDrcfAwKL5HgCAkxByAcANJUta8wBXrZq/89PTrXmGLw/Cl2+//n979x4UZb3HcfyzCwKLkJQKIXkENc0kJC9ENaNJThqZzZRkdTxWGscuRmlnTlmZzWCjNpRO1+k2YIPjVKaT3e1iOR4VTxctj1ZSSQ0SBSK4CovCnj9WltuirO7uwz77fs38Zvf5PRe+8Hj5+PO3v2d36/vaWs/Xio7uHHy72m7bFxvLHGMAoYOQCwABEBbWOkLbXceOuYLxwYOtreN2S9/u3e23m5s91+ApCLeMYMfFuaZceNru08f1oUAACBaEXADooXr1cn2wLT7eu/Oam11zhrsKxG23y8pcaxQfOuQaOfY017iFzda9QNzVPpuNkWQAgUPIBQCTsVpbw2VKinfnNjW5gm5trSv4trSutisrpZ9+ar/P0yiy5ArtcXGtHwbsTouN7dwXHU1YBnBqhFwAgFtYWOvybKfD6ZSOHDl5QK6ra22HD7se+NG2r65Ocji6/hpWa/fCcFehOSamtUVFEZgBsyLkAgB8xmJpDZDnnXf613E4XAG4JQh3DMGe2qFDrgeBdAzRJxMW1j70tm0dA3F3+wjOQM9AyAUA9DiRka7mzQf1PGluluz29qHXbm/fOva1bP/5p/TLL+37Dh+Wjh8/+dfsKjh3DMPR0VLv3q2tO9vh/K0NdBu/XQAAptV2aoOvNDZ2HYxPFZ4rK10PFbHbXdM6jh51vdbXd+9rR0R4F4q93bZaffdzAoxGyAUAwAsREa0P6vCV5mZX4G0JvS2t7fbJ9rVsV1V53neyOc5tRUa6wq7N5mot733dZ7MxKg3/C+1fYnv3Gl0BAACySoo50dx6Sepzop2hpibXU/rq61tf6xukhrbb9a4w3NDQ2hyO1j5Hg9RwyPX+aINU03afw3W9U03laCs8XIqKdM1hjuzw2rGvbX/bvpYWEdH6GhHpum67vhPvGakOMmeY00I75M6caXQFAAD4XZik3idaj3H8RDtidCEwq9AOucXF0ogRRlcBAAC81Nzsmh/d2Ng64tzyvuOr+32j1HjifUND575uX8chOb2sN8zafuS540hzRIRrLem2r57eR0Z63n+yc9r1RUrhYUGyAsjevWc0IBnaIXfECGn0aKOrAAAAXrJKijrRAs3pdD12u2VaR8s0kI6tu/0Oh1TraA3SDofkqG0/ZcRTa2o6vfotls5TPrxpnqaMnKq1DfWegn6vXr4P3qEdcgEAALxksbSGNF+u3OGtpqb2ofdUodhT6845NTXdO8ebOdmedAzCTm+Hyzsg5AIAAAShsDDXihXR0UZX4tLc7Dn8djX9o+10E0+vZWVSYeHp10PIBQAAwBmzWluXiPOFb745s5DLYhoAAAAwHUIuAAAATIeQCwAAANMh5AIAAMB0CLkAAAAwHUIuAAAATIeQCwAAANMh5AIAAMB0CLkAAAAwHUIuAAAATIeQCwAAANMh5AIAAMB0CLkAAAAwHUIuAAAATIeQCwAAANMh5AIAAMB0CLkAAAAwHUIuAAAATIeQCwAAANMh5AIAAMB0CLkAAAAwnaAPuY2NjXrwwQeVlJSk6OhoZWZm6tNPPzW6LAAAABgo6EPurbfeqpUrV+of//iHnnnmGYWHhys7O1tbt241ujQAAAAYJKhD7o4dO/TGG29o2bJlWrZsme644w599tlnGjRokP79738bXR56kDVr1hhdAgKI+x1auN+hhfuN7grqkLt27VqFh4crNzfX3RcZGak5c+Zo27ZtKi8vN7A69CT8oRhauN+hhfsdWrjf6K6gDrk7d+7UsGHDFBMT064/IyPDvR8AAAChJ6hDbkVFhRITEzv1JyYmyul06sCBAwZUBQAAAKMFdcitr69XZGRkp/6oqCj3fgAAAISecKMLOBM2m00Oh6NTf0NDg3u/Jy3hd+/evf4rDj1KbW2tvvnmG6PLQIBwv0ML9zu0cL9DR0tOO91By6AOuYmJiR6nJFRUVEiSBgwY4PG8/fv3S5Jmzpzpt9rQ84wZM8boEhBA3O/Qwv0OLdzv0LJ//35dfvnlXp8X1CE3PT1dX3zxhex2e7sPn23fvl0Wi0Xp6ekez5s8ebKKi4uVnJzc5WgvAAAAjFNfX6/9+/dr8uTJp3W+xel0On1cU8Ds2LFDmZmZKigo0IIFCyS5noCWmpqq/v376z//+Y/BFQIAAMAIQT2Sm5GRoZycHC1cuFCVlZUaOnSoioqKVFZWpsLCQqPLAwAAgEGCeiRXco3cLlq0SMXFxaqpqVFaWpqWLFmiSZMmGV0aAAAADBL0IRcAAADoKKjXyQUAAAA8CamQ29jYqAcffFBJSUmKjo5WZmamPv30U6PLgh8cOXJEixcv1tVXX62+ffvKarXq9ddfN7os+MlXX32lefPmKTU1VTExMRo0aJBmzJihffv2GV0a/GDPnj268cYbNWTIEPXu3Vv9+/fXhAkT9N577xldGgLgiSeekNVqVVpamtGlwA++/PJLWa3WTi0sLEw7duzw6lpB/cEzb916661at26d5s+f7/6QWnZ2tr744gtddtllRpcHH6qqqlJ+fr4GDRrkXmoO5rV8+XJt3bpVOTk5SktL0x9//KFnn31Wo0ePVklJiS688EKjS4QPlZWVyW6367bbbtOAAQN09OhRvf3225o2bZpefvll3XHHHUaXCD8pLy/X0qVL2y0bCnO6//77NXbs2HZ9Q4cO9eoaITMnt2W5saeeekrz58+XJDkcDqWmpiohIUFbtmwxuEL40rFjx1RTU6P4+Hh9/fXXGjdunIqKijRr1iyjS4MfbN++XWPHjlV4eOu/20tLS3XRRRcpJyeHUfwQ4HQ6NXr0aDkcDu3Zs8focuAnN910k6qrq3X8+HFVV1fru+++M7ok+NiXX36piRMnau3atbr++uvP6FohM11h7dq1Cg8PV25urrsvMjJSc+bM0bZt21ReXm5gdfC1Xr16KT4+3ugyECCZmZntAq7k+hf/yJEjeXx3iLBYLBo4cKAOHTpkdCnwk82bN2vdunVauXKl0aUgQOx2u5qamk77/JAJuTt37tSwYcM6/RdHRkaGez8Ac6msrFS/fv2MLgN+cvToUVVXV+uXX37RihUr9OGHH7J8pEk1NzcrLy9Pubm5GjlypNHlIABuv/12nXXWWYqKilJWVpa+/vprr68RMnNyKyoqlJiY2Kk/MTFRTqdTBw4cMKAqAP5SXFys8vJyLVmyxOhS4CcPPPCAXnrpJUmS1WrVDTfcoGeffdbgquAPL774on777Td9/vnnRpcCP4uIiND06dOVnZ2tfv36ac+ePSooKND48eO1detWjRo1qtvXCpmQW19fr8jIyE79UVFR7v0AzOGHH37QvHnzdPnllzMP28Tmz5+vnJwcHThwQG+++aaamprkcDiMLgs+dvDgQS1evFiPPfaYzjnnHKPLgZ9deumluvTSS93bU6dO1Q033KC0tDQtXLhQH3zwQbevFTLTFWw2m8c//BoaGtz7AQS/yspKXXPNNTr77LP11ltvyWKxGF0S/GTYsGHKysrSzJkztWHDBh0+fFjTpk0zuiz42COPPKK+fftq3rx5RpcCgwwZMkTXXXedNm3aJG/WSwiZkJuYmKiKiopO/S19AwYMCHRJAHysrq5OU6ZMUV1dnT766COde+65RpeEAJo+fbr++9//sj6yiZSWluqVV15RXl6eysvLVVZWpv3796uhoUHHjh1TWVmZampqjC4TATBw4EA1NjbqyJEj3T4nZEJuenq6fvrpJ9nt9nb927dvl8ViUXp6ukGVAfAFh8OhqVOnqrS0VO+//76GDx9udEkIsJZpZ7W1tQZXAl8pLy+X0+lUXl6eUlJSlJKSosGDB6ukpEQ//vijBg8erPz8fKPLRAD8/PPPioqK8mqN5JCZkzt9+nQVFBTo5Zdf1oIFCyS5noBWVFSkzMxMJSUlGVwhgNPV3NysG2+8USUlJdqwYYN71RSY019//aX+/fu36zt+/LhWrVolm83Gwz9MJDU1VevXr+/U/8gjj8hut+uZZ57R4MGDDagM/lJVVdVpVZxdu3bp3Xff1TXXXOPVtUIm5GZkZCgnJ0cLFy5UZWWl+4lnZWVlKiwsNLo8+MHzzz+vQ4cOuddA3rBhg37//XdJUl5enmJjY40sDz60YMECvfvuu5o2bZqqqqq0evXqdvv//ve/G1QZ/GHu3Lmqq6vT+PHjlZSUpD/++EOrV6/Wjz/+qKefflrR0dFGlwgf6du3r8d51itWrJDFYtG1115rQFXwpxkzZshms+myyy5TfHy8/ve//+mVV15RTEyMli5d6tW1QuaJZ5Jr5HbRokUqLi5WTU2N0tLStGTJEtZVNKmUlBT99ttvHvf9+uuv+tvf/hbgiuAvEydO1ObNm7vcfyaLiaPnefPNN/Xaa6/p+++/V3V1tWJjYzVmzBjl5eV5PdKD4DRx4kQdPHhQu3btMroU+Nhzzz2n1atXq7S0VHV1derfv78mTZqkxx57zOtR+5AKuQAAAAgNIfPBMwAAAIQOQi4AAABMh5ALAAAA0yHkAgAAwHQIuQAAADAdQi4AAABMh5ALAAAA0yHkAgAAwHQIuQAAADAdQi4ABJGioiJZrdYuH1kNAHAh5AJAAK1atUpWq9XdbDabhg8frnvvvVd//vnnKc+3WCyyWCwBqBQAgpvF6XQ6jS4CAELFqlWrNHv2bOXn5ys5OVkNDQ3asmWLXn/9dSUnJ2v37t2Kiorq8nyn06ljx44pIiIigFUDQPAJN7oAAAhFU6ZM0ejRoyVJs2fP1jnnnKMVK1bonXfe0YwZMzodf/ToUUVHR8tisRBwAaAbmK4AAD1AVlaWnE6nfv31V/eUhs2bN+vuu+9WQkKCBg4cKKnrObkffvihJkyYoLPOOkt9+vRRRkaG1qxZ0+6YkpISTZkyRXFxcerdu7euuOIKbd26td0xdrtd999/v1JSUhQVFaWEhARdddVV2rlzp39/AADgY4zkAkAPUFpaKknq27evu+/uu+9WfHy8Fi9erCNHjkjyPCe3qKhIc+bMUWpqqh5++GHFxcXp22+/1ccff6ybb75ZkvT5558rOztbY8eO1eOPPy6r1arCwkJlZWVpy5YtGjt2rCRp7ty5Wrdune69916NGDFC1dXV2rJli/bu3av09PRA/CgAwCcIuQBggNraWlVXV7vn5Obn56t3796aOnWqNm7cKEnq16+fPvvss5N+0Kyurk733XefMjMztWnTpi6nMtx111268sor9f7777v75s6dqwsvvFCPPvqoPvroI0nSBx98oNzcXD355JPu4/71r3/54lsGgIAi5AJAgDmdTl155ZXubYvFouTkZK1Zs0aJiYnuvtzc3FOupPDJJ5/IbrfroYce6jLg7ty5U/v27dOiRYtUXV3dqY7i4mJ3X1xcnEpKSlRRUeGuBQCCESEXAALMYrHohRde0Pnnn6/w8HAlJCRo+PDhnY5LTk4+5bV+/vlnSdLIkSO7PGbfvn2SpFmzZnncb7VaVVtbqz59+ujJJ5/UbbfdpoEDB2rMmDHKzs7WrFmzlJKS0o3vDAB6DkIuABhg3Lhx7tUVumKz2XzytZqbmyVJTz31lEaNGuXxmJiYGElSTk6Oxo8fr/Xr12vjxo0qKCjQ8uXLtX79ek2ePNkn9QBAIBByASCIDRkyRE6nU7t379bgwYO7PEaSYmNjlZWVdcprJiQk6M4779Sdd96pqqoqXXzxxXriiScIuQCCCkuIAUAQu+qqqxQbG6ulS5fK4XB4PGbMmDEaMmSICgoK3Ks0tFVVVSXJNeJbV1fXbl+/fv00YMCALq8NAD0VI7kAEGDdedBkdx9GGRsbqxUrVig3N1fjxo3TLbfcorPPPlu7du1SfX29CgsLZbFY9Oqrryo7O1sjR47U7bffrqSkJJWXl2vTpk3q06eP3nnnHR0+fFjnnXeepk+frlGjRikmJkaffPKJvvrqKz399NNn+m0DQEARcgEgwE61YkJ3j2kxe/ZsJSQkaNmyZVqyZIl69eqlCy64QPPnz3cfM2HCBG3btk35+fl6/vnnZbfbde655+qSSy7R3LlzJUnR0dG65557tHHjRq1fv17Nzc0aOnSoXnzxRf3zn//0/hsFAANZnN0dLgAAAACCBHNyAQAAYDqEXAAAAJgOIRcAAACmQ8gFAACA6RByAQAAYDqEXAAAAJgOIRcAAACmQ8gFAACA6RByAQAAYDqEXAAAAJgOIRcAAACmQ8gFAACA6RByAQAAYDr/B1adplGcsHnlAAAAAElFTkSuQmCC",
      "text/plain": [
       "PyPlot.Figure(PyObject <matplotlib.figure.Figure object at 0x000000001B698908>)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "PyObject <matplotlib.text.Text object at 0x000000001B392710>"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "using PyPlot\n",
    "#Setting the demand function\n",
    "function demand(p)\n",
    "  return 2.0.*p.^(-0.5)\n",
    "end\n",
    "prices=linspace(0,5,200)\n",
    "dem=demand(prices)\n",
    "figure(1)\n",
    "plot(prices,dem, color=\"blue\")\n",
    "plot(prices, ones(200).*2, color=\"red\")\n",
    "plot(prices, ones(200), color=\"red\")\n",
    "title(\"Equilibria for p=1 and p=4\" )\n",
    "xlabel(\"Prices\")\n",
    "ylabel(\"Quantities\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Analytically, the change in consummer surplus (S) is computed as:\n",
    "$$\\Delta S(p_1,p_2)=S(p_1)-S(p_2)$$\n",
    "which is:\n",
    "$$\\Delta S(1,4)= S(1)-S(4)$$\n",
    "$$\\Delta S(1,4)= 2-\\int_1^\\infty 2p^{−0.5} dp - \\left(4 - \\int_4^\\infty 2p^{−0.5} dp \\right) $$\n",
    "$$\\Delta S(1,4)= -2 + \\int_1^4 2p^{−0.5} dp $$\n",
    "$$\\Delta S(1,4)= -2 + \\left[4 \\times p^{0.5} \\right]_1^4$$\n",
    "$$\\Delta S(1,4)= -2 + 8 - 4 = 2$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 0.4.3",
   "language": "julia",
   "name": "julia-0.4"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "0.4.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}