Numerical solution to Exercise 6.2 in Kreps

Kreps6_2.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import scipy.optimize\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "f = lambda a: np.exp(-24*a) - 0.5 * np.exp(-50*a) - 0.5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAgIAAAEyCAYAAACI17j8AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xd8VuX9//HXJ3tAQhIChr2CCIiMMBQHKiC2VbA4sKK4\nigO1VTv012GrtbVqW6tWFEHFibiKOKqAOFBWkCGyEoYQCBAIECA7uX5/5OD3FoMBMk6S+/18PM7j\nvs91X+e6P+cYyTtnmnMOERERCU4hfhcgIiIi/lEQEBERCWIKAiIiIkFMQUBERCSIKQiIiIgEMQUB\nERGRIKYgICIiEsQUBERERIKYgoCIiEgQC/O7gLrQvHlz16FDB7/LEBERqTNLlizZ5ZxLrqpfUASB\nDh06kJ6e7ncZIiIidcbMvjmafjo0ICIiEsQUBERERIKYgoCIiEgQUxAQEREJYgoCIiIiQUxBQERE\nJIgpCIiIiAQxBQEREZEgpiAgIiISxILizoLyw4pKy9hXUMKBwlLyi8soKCkjv7iMwpIySsrKKS1z\nFJeVU1buvrdsWIgRHhpCWGjFa3R4KDERoUR5r02jwomLDiMyLNSHNRMRkaooCDRShSVl7MgrJHtf\nIdv3Vbzm7C9i14Eidh8sYveBYvbkF7OvoITCkvJarycqPIS4qHASYiJIahJBYmwESbERtIiLokXT\nSE6Ij6JlXBStmkXTJFI/liIidUX/4jZgB4tK2bjrIOtzDrAh5yBbcvPZnJvPlj357Mgr+l7/JpFh\nJDWp+AXcNjGGU9o0Iy46jPjocOKjw2kSFUZ0eBixkRV/zUeGhRIRFlLxF3+IERZqGPbteOXOUVbu\nKvYalDuKS8sp9PYmVOxVKGV/YSl5BSXsLyxlb34JufnF5B4s5utteew+UEReYen36mwWE07rZtG0\nSYimXWIMHZrH0jEplvbNY0mJiyIkxL63jIiIHB8FgQaguLScjJ37Wbu9YlqzfT/rduwne1/ht33M\noFV8xS/PM1OTaZsYQ6tm0aTER3FCfBQnxEURWw//0j6052JHXhHb8wrZtreArD35bN1TwPqcg8xd\nm0Nx6f/tsYgOD6Vzi1hSWzSlS4smdG3ZlJNSmtK6WTRmCggiIseq/v1mCHJl5Y51O/azbMtevtq6\nj5Vb97Emez/FZRW/DCNCQ+jcogkDOyaS2rIpnZrH0im5Ce2TYogKb3jH4aPCQ2mfFEv7pNhKPy8v\nd2TnFfLNroNs3H2Q9TsPkrFzPws27OatpVu/7RcXFUa3lDh6tIqjV5t4Tm7djE7NY7X3QESkCubc\n908Aa2zS0tJcfX0M8f7CEpZ8s4f0TXv4cvMelm/Zy8HiMgCaRoXRs1U8J7eJp0erOLqnxNGheSzh\nobrYAyCvsISMHftZlb2f1dl5rM7OY032fgpKKrZfbEQoPVvH07d9An3bJdC3XTOSmkT6XLWISN0w\nsyXOubSq+mmPQB07UFTKwg27mb9+Nws35vL1tn2UOwgNMU5Kacrofm3o064Zfdom0D4pRru7f0Bc\nVDj92ifSr33it22lZeWszznIiqyKPSrLt+zl6U83UOpd8dCxeSz9OyQwoGMSAzsm0iZBhxREJLhp\nj0AtKyt3LM/ay2frdjEvM4elm/dSWu6ICAuhT9tmDOyYyMBOSfRp14yYCOWy2lBYUsaKrH18ubli\nz0v6N7nszS8BICU+ilM7JzG4c3MGd2nOCfFRPlcrIlIzjnaPgIJALdiXX8LH63by8docPlmXQ+7B\nYszg5NbxDO7SnDO6NKdv+4QGeUy/MSgvd2TsPMCijbtZsCGX+Rt2k3uwGIBOybGcmZrMWScmM6hj\nEtER+m8kIg1TnQcBMxsB/BsIBSY75x447PNI4HmgH7AbuMw5t8n77G7gOqAMuM0598EPjWlmHYFp\nQBKwBLjSOVd8pNrqIghk7yvgw6938OGq7SzYkEtZuSMhJpwhJ7ZgyInJnJmaTEJsRK3WIMenvNyx\nenseX2TuZl7mLhZs2E1RaTkRYSEM7JjI2Se2YOhJLWmXFON3qSIiR61Og4CZhQLrgGFAFrAYuNw5\ntyqgz81AL+fcjWY2BrjIOXeZmXUHXgEGAK2A2UBXb7FKxzSz6cCbzrlpZvYksNw5N/FI9dVWENi6\nt4D3VmTzzlfZLN+yF4DOybEM73ECw7q35JQ2zQjVWesNTmFJGYs25vLJuhw+XruT9TkHAUht0YRz\nT2rJsO4t6dO2ma5IEJF67WiDQE2dfj4AyHTObfD+Mp8GjDysz0hgqvf+deBcqzhLayQwzTlX5Jzb\nCGR641U6prfMOd4YeGOOqqH1qFJ+cSlT5m3kxAmTGPzAR9z/3mrKyx2DYnbSatkUznfp/HZEN2IL\ndnLuOWczZMiQb5cdP348Q4YMYebMmQDMnDmTIUOGMH78+G/7DBkyhCFDhrB27VoAHn74YYYMGcLD\nDz8MwNq1a7/to3FrZ9wRw87lj+Mv4Q8/6c6cO4dw+t7ZJGz6CCvaz+TPNjB64hek/mo6Z/zyceZl\n7KKkrLxRbgeNq3E1bt2O65eaCgKtgS0B81leW6V9nHOlwD4qdu0fadkjtScBe70xjvRdmNl4M0s3\ns/ScnJzjXK3vCzHjX7PW4ULCaLb5U54d3Z6Zt55OWsxuIgpza+x7pP6IKT9I/PYl3NrTseQPw7gq\ntZzI/dvYGtmWsVMWMuD+2ezqOJyCuHaVPo9BRKQ+q6lDAxcDI5xz13vzVwIDnXO3BPRZ6fXJ8ubX\nAwOBPwELnHMveu1TgPe9xb43ZkD/Ll57W+B951zPI9VX04cGdh0oormuRw96BcVlfLIuh/e+ymb2\n6h3kF5fRvEkE5/dMYVSfVvRtl6BLE0XEN3V9H4GtQNuA+TZeW2V9sswsDIin4qTBH1q2svbdQDMz\nC/P2ClT2XbVKIUAAoiNCGdHzBEb0PIHCkjLmrtnJOyuymZ6+hRcWfEO7xBhG9m7FqD6t6ZzcxO9y\nRUQqVVNBYDGQ6p3NvxUYA/zssD5vA+OA+cDFwEfOOWdmbwMvm9k/qThZMBVYBFhlY3rLzPXGmOaN\nOaOG1kPkuESFh3L+ySmcf3IK+wtL+ODrHcxYtpX/zM3ksY8y6dOuGRf3a8NPerUiPjrc73JFRL5V\nk5cP/gh4hIpL/Z5xzt1vZvcC6c65t80sCngB6APkAmOccxu8ZX8HXAuUAr90zr1/pDG99k5UhIBE\nYCkw1jn3/cfteerzLYalcduZV8iMZdt4bckW1u04QERYCCN6nMCY/m0Z1ClJVx6ISK3RDYUCKAiI\n35xzrNyax2tLtjBj2Tb2FZTQPimGS9Packm/NrSI0x0NRaRmKQgEUBCQ+qSwpIz/rdzOK4s2s3Bj\nLqEhxvDuLblyUHtO7ZykEwxFpEYoCARQEJD6akPOAaYt3sL09C3szS+hU3IsYwe25+K0NsRF6VwC\nETl+CgIBFASkvissKePdFdm8uPAblm7eS0xEKKP7tmHcaR3o0kJXHIjIsVMQCKAgIA3Jyq37ePbz\nTcxcvo3isnLOSG3Odad35KyuyTpsICJHTUEggIKANES7DhTxysLNvLDgG3buLyK1RROuP6MjI3u3\n1pMrRaRKCgIBFASkISsuLWfm8m08/dkG1mzfT/MmEVx9WgeuHNSB+BidRyAilVMQCKAgII2Bc44v\n1u9m0qcb+GRdDrERoVw+oB3XndGRlPhov8sTkXpGQSCAgoA0Nqu25THp0/XMXJFNiMFP+7ThpiGd\n6dA81u/SRKSeUBAIoCAgjdWW3Hye/mwD0xZvobSsnAtOacWEs7vQtWVTv0sTEZ8pCARQEJDGbuf+\nQqZ8tpEXFnxDfnEZ5/c8gdvOTeWklDi/SxMRnygIBFAQkGCx52Axz3y+kec+38T+olIFApEgpiAQ\nQEFAgs3e/GKembeRZ71A8KOTT+D2oV1J1SEDkaChIBBAQUCC1b78EqbM28CUeRspKCljVO/W/GJo\nKu2TdFKhSGOnIBBAQUCCXe7BYp76ZD1T52+itMxxaf+2/OLcVFrqqYcijZaCQAAFAZEKO/MKeXxu\nJq8s2kxoiHHt4I7ccFZn4qN1YyKRxkZBIICCgMh3bd6dzz9mrWXGsm3ER4cz4ezOXHVqB926WKQR\nOdogEFIXxYhI/dIuKYZ/j+nDu7edTu+2zfjre2s49x+fMGPZVsrLG/8fByLyfxQERIJYj1bxTL12\nAC9dP5D46HB+MW0Zo574nAUbdvtdmojUEQUBEWFwl+a8c+vp/PPSU8jZX8SYSQu44YV0vtl90O/S\nRKSWKQiICAAhIcZP+7Zh7q+G8OvzTuSzjF0M/ecn/PW91eQVlvhdnojUEgUBEfmOqPBQJpzdhY9/\nNYRRvVvz9GcbOPuhj5m2aDNlOn9ApNFREBCRSrWIi+KhS05h5i2n0yk5lrve/IpR//mcJd/k+l2a\niNQgBQER+UE9W8cz/YZT+feY3uTsL2L0xPnc8eoydu4v9Ls0EakBCgIiUiUzY2Tv1sy58ywmnN2Z\nd1Zkc+7Dn/DMvI2UlpX7XZ6IVIOCgIgctdjIMH59Xjc+uP1M+rRP4N53VvGTx+axeJMOF4g0VAoC\nInLMOjaPZeo1/XlybF/yCkq45Mn5/Pq15eQeLPa7NBE5RtUKAmaWaGazzCzDe004Qr9xXp8MMxsX\n0N7PzL4ys0wze9TMzGt/yMzWmNkKM3vLzJp57R3MrMDMlnnTk9WpX0SOn5kxomcKs+88ixvP6sxb\nS7dy7j8+Znr6FoLh1uUijUV19wjcBcxxzqUCc7z57zCzROAeYCAwALgnIDBMBH4OpHrTCK99FtDT\nOdcLWAfcHTDkeudcb2+6sZr1i0g1xUSEcdf53Xj3tjPonNyE37y+gsueWkDmzv1+lyYiR6G6QWAk\nMNV7PxUYVUmf84BZzrlc59weKn7JjzCzFCDOObfAVfz58Pyh5Z1zHzrnSr3lFwBtqlmniNSyE09o\nyvQbTuXB0b1Yt3M/P/r3PB6ZvY6i0jK/SxORH1DdINDSOZftvd8OtKykT2tgS8B8ltfW2nt/ePvh\nrgXeD5jvaGZLzewTMzvjuCsXkRoXEmJc2r8ts+84ixE9T+CR2Rn8+FGdTChSn1UZBMxstpmtrGQa\nGdjP+6u+Rg8MmtnvgFLgJa8pG2jnnOsD3AG8bGZxR1h2vJmlm1l6Tk5OTZYlIlVo3iSSRy/vw7PX\n9KeguIxLnpzP7//7FQeKSqteWETqVJVBwDk31DnXs5JpBrDD28WP97qzkiG2Am0D5tt4bVv57i7/\nQ+14410N/AS4wgsZOOeKnHO7vfdLgPVA1yPUPck5l+acS0tOTq5qNUWkFpx9Ygs+vP1Mrh3ckZcW\nbmb4Pz/h47WV/TMhIn6p7qGBt4FDVwGMA2ZU0ucDYLiZJXgnCQ4HPvAOKeSZ2SDvaoGrDi1vZiOA\n3wAXOufyDw1kZslmFuq970TFCYYbqrkOIlKLYiPD+OMF3XnjptOIiQzj6mcXc8f0ZezN16WGIvVB\ndYPAA8AwM8sAhnrzmFmamU0GcM7lAvcBi73pXq8N4GZgMpBJxV/3h84FeBxoCsw67DLBM4EVZrYM\neB24MWAsEanH+rZL4N3bTue2c7rw9rJtDPvXp8xatcPvskSCngXD9b5paWkuPT3d7zJExPP1tn38\n6rUVrM7O46I+rbnngu40i4nwuyyRRsXMljjn0qrqpzsLikid69EqnhkTBvOLc1OZubxi78Bs7R0Q\n8YWCgIj4IiIshNuHdeW/EwaTFBvB9c+n8+vXlrO/sMTv0kSCioKAiPiqZ+t43r7ldCac3Zk3vsxi\nxCOf8cX6XX6XJRI0FARExHcRYSH8+rxuvH7TaUSEhfCzpxfy55lfU1iiuxKK1DYFARGpN/q2S+C9\n285g3KntefbzTVzw2DxWbt3nd1kijZqCgIjUK9ERofx5ZE+ev3YA+wpKuOiJz3ni40zKyhv/FU4i\nflAQEJF66cyuyXzwyzMZ1r0lD/5vLZdPWsDWvQV+lyXS6CgIiEi9lRAbwX9+1pd/XHIKq7LzGPHI\np8xcvs3vskQaFQUBEanXzIzR/drw3m1n0KVFE259ZSl3Tl+uBxiJ1BAFARFpENolxTD9hlO57Zwu\nvLU0ix8/+hkrsvb6XZZIg6cgICINRnhoCHcMP5Fp40+lpLSc0RO/4OlPN1CuEwlFjpuCgIg0OAM6\nJvLeL87g7BNbcP97q7nmucXsOlDkd1kiDZKCgIg0SM1iInjqyn7cN7IH8zfs5vx/646EIsdDQUBE\nGiwz48pTOzBjwmCaRoUxdvJCHpm9TvccEDkGCgIi0uCdlBLHzFtOZ1Tv1jwyO4Mrpyxk5/5Cv8sS\naRAUBESkUYiNDOMfl57Cgxf34svNe/jRv+fxRaYOFYhURUFARBoNM+PStLa8fcvpNIsJZ+yUhTw2\nJ0NXFYj8AAUBEWl0urZsyowJg7nglFb8Y9Y6rnluMbkHi/0uS6ReUhAQkUYpNjKMRy7rzV9G9WT+\n+t385NHPWLp5j99lidQ7CgIi0miZGWMHteeNm04jJMS49Kn5vLDgG5zToQKRQxQERKTRO7lNPO/c\nejqnd2nOH/67kjunL6eguMzvskTqBQUBEQkKzWIimDKuP3cM68pby7Zy0ROfs2nXQb/LEvGdgoCI\nBI2QEOO2c1N57poBbM8r5MLH5zF3zU6/yxLxlYKAiASds7omM/OW02mTEMO1UxfzqC4xlCCmICAi\nQaltYgxv3HQao3q35p+z1jH+hSXkFZb4XZZInVMQEJGgFR0Ryj8vPYU/XdCdj9fuZNR/Pmd9zgG/\nyxKpUwoCIhLUzIyrB3fkxesHsi+/hFGPf85Ha3b4XZZInal2EDCzRDObZWYZ3mvCEfqN8/pkmNm4\ngPZ+ZvaVmWWa2aNmZl77n8xsq5kt86YfBSxzt9d/rZmdV911EBEZ1CmJt289nXZJMVw3NZ3HP8rQ\n/QYkKNTEHoG7gDnOuVRgjjf/HWaWCNwDDAQGAPcEBIaJwM+BVG8aEbDov5xzvb3pPW+s7sAYoIfX\n9wkzC62B9RCRINe6WTSv33gaF/RqxcMfruOWl5eSX1zqd1kitaomgsBIYKr3fiowqpI+5wGznHO5\nzrk9wCxghJmlAHHOuQWuIno/f4TlD/++ac65IufcRiCTinAhIlJt0RGh/HtMb+4+vxvvrczmkifn\ns3Vvgd9lidSamggCLZ1z2d777UDLSvq0BrYEzGd5ba2994e3H3KLma0ws2cC9iAcaazvMLPxZpZu\nZuk5OTnHtEIiEtzMjBvO6swz4/qzeXc+Ix+fx5Jvcv0uS6RWHFUQMLPZZraykmlkYD/vr/qaOqg2\nEegM9AaygX8cy8LOuUnOuTTnXFpycnINlSQiweTsbi14a8JpxEaGcfmkhbyWvqXqhUQamKMKAs65\noc65npVMM4Ad3i5+vNfKbtO1FWgbMN/Ga9vqvT+8HefcDudcmXOuHHia/9v9f6SxRERqXJcWFY80\n7t8xgV+/voK/vbeaMt18SBqRmjg08DZw6CqAccCMSvp8AAw3swRvF/9w4APvkEKemQ3yrha46tDy\nh8KF5yJgZcD3jTGzSDPrSMUJhotqYD1ERCrVLCaC564ZwNhB7Xjq0w3c8MISDhTpJEJpHGoiCDwA\nDDOzDGCoN4+ZpZnZZADnXC5wH7DYm+712gBuBiZTcdLfeuB9r/1B77LCFcDZwO3eWF8D04FVwP+A\nCc45PUZMRGpVeGgIfxl1Mn++sAcfrdnBxRO/IGtPvt9liVSbBcN1smlpaS49Pd3vMkSkkfh0XQ4T\nXv6SyLAQnr4qjT7tKr19ioivzGyJcy6tqn66s6CIyDE6s2syb918GtERoYyZtIB3V2RXvZBIPaUg\nICJyHLq0aMp/bx5Mz9bxTHj5S/4zN1N3IpQGSUFAROQ4JTWJ5KXrB3LBKa146IO1/Pr1FRSXlvtd\nlsgxCfO7ABGRhiwqPJRHx/SmY1IMj36Uyba9BUwc24/46HC/SxM5KtojICJSTWbGHcNP5OFLTmHR\nxlxdUSANioKAiEgNubhfG56/dgDb8wq56IkvWJG11++SRKqkICAiUoNO69KcN286jYjQEC57agEf\nrdnhd0kiP0hBQESkhqW2bMpbE06jS4smXD81nZcXbva7JJEjUhAQEakFLZpGMW38IM7smsz/e+sr\nHv5grS4vlHpJQUBEpJbERoYx+ao0xvRvy+NzM7lz+nJdXij1ji4fFBGpRWGhIfztpyfTulk0/5i1\njpwDRUwc248mkfrnV+oH7REQEallZsat56by4OhefLF+N2MmzSdnf5HfZYkACgIiInXm0v5tefqq\nfqzfeZDRE79g466DfpckoiAgIlKXzunWklfGD+JAUSmjJ37B8i2614D4S0FARKSO9W7bjDduOo2Y\niFAuf3oBn2Xk+F2SBDEFARERH3RsHsubN51Gu8QYrn1uMTOWbfW7JAlSCgIiIj5pERfFqzecSp92\nCfxi2jKe/Xyj3yVJEFIQEBHxUXx0OM9fO4Dh3Vvy55mrdOMhqXMKAiIiPosKD+WJK/pyWVrFjYf+\nMGMlZeUKA1I3dEcLEZF6ICw0hAdGn0xCbARPfrKevfkl/PPS3kSE6e81qV0KAiIi9YSZcdf53UiI\nCedv768hr7CUJ8f2JSZC/1RL7VHUFBGpZ244qzMPju7FvIwcrpyyiH0FJX6XJI2YgoCISD10af+2\nPHFFX1Zk7WXMpAW6JbHUGgUBEZF6akTPFKaM68+mXQe59Kn5ZO3J97skaYQUBERE6rEzuybz4vUD\n2H2giEuenM/6nAN+lySNjIKAiEg91699ItPGn0pJWTmXPTWf1dl5fpckjUi1goCZJZrZLDPL8F4T\njtBvnNcnw8zGBbT3M7OvzCzTzB41M/PaXzWzZd60ycyWee0dzKwg4LMnq1O/iEhD0b1VHK/ecCrh\noSGMmbSAZXpYkdSQ6u4RuAuY45xLBeZ4899hZonAPcBAYABwT0BgmAj8HEj1phEAzrnLnHO9nXO9\ngTeANwOGXH/oM+fcjdWsX0Skweic3ITpN5xKfHQ4Vzy9gAUbdvtdkjQC1Q0CI4Gp3vupwKhK+pwH\nzHLO5Trn9gCzgBFmlgLEOecWuIr7aT5/+PLeHoJLgVeqWaeISKPQNjGG1248lZRm0Yx7ZhEfr93p\nd0nSwFU3CLR0zmV777cDLSvp0xrYEjCf5bW19t4f3h7oDGCHcy4joK2jmS01s0/M7IxqVS8i0gC1\njIvi1fGD6NKiCT9/Pp0Pv97ud0nSgFUZBMxstpmtrGQaGdjP+6u+pm+OfTnf3RuQDbRzzvUB7gBe\nNrO4I9Q93szSzSw9J0fP+haRxiWpSSQvXz+IHq3iufmlL3lnxTa/S5IGqsog4Jwb6pzrWck0A9jh\n7eLHe61sH9VWoG3AfBuvbav3/vB2vPHCgJ8CrwbUUuSc2+29XwKsB7oeoe5Jzrk051xacnJyVasp\nItLgxMeE88J1A+jbLoHbXlnKG0uyql5I5DDVPTTwNnDoKoBxwIxK+nwADDezBO8kweHAB94hhTwz\nG+SdC3DVYcsPBdY45779yTazZDML9d53ouIEww3VXAcRkQaraVQ4z13bn9M6N+fO15bz8sLNfpck\nDUx1g8ADwDAzy6DiF/cDAGaWZmaTAZxzucB9wGJvutdrA7gZmAxkUvHX/fsBY4/h+ycJngms8C4n\nfB24MWAsEZGgFBMRxuRxaZx9YjL/762veH7+Jr9LkgbEKg7tN25paWkuPT3d7zJERGpVUWkZt7y8\nlFmrdvD7H5/E9Wd08rsk8ZGZLXHOpVXVT3cWFBFpJCLDQnniir6c3/ME/vLuap78ZL3fJUkDoCAg\nItKIhIeG8NjlfbjglFY88P4aHv8oo+qFJKiF+V2AiIjUrLDQEB65rDfhIcbDH66jrBx+MTTV77Kk\nnlIQEBFphEJDjIcuOQUz41+z11HuHLcPq/RqawlyCgIiIo1UaIjx4MW9CDH495wMHHD70FS857uJ\nAAoCIiKNWmiI8ffRvQgx49E5GZSXO+4c3lVhQL6lICAi0siFhBh/++nJmMHjczMxgzuGKQxIBQUB\nEZEgEBJi/PWikwF47KNMzIw7dM6AoCAgIhI0DoWBcud4dE4GBjqBUBQERESCSUiI8cBPe+FcxQmE\nZvDLoQoDwUxBQEQkyIR4JxA64JHZGYSFGLeco/sMBCsFARGRIHQoDJSVOx7+cB1hoSHceFZnv8sS\nHygIiIgEqdAQ4+FLTqG03PHA+2sICzE9qCgIKQiIiASx0BDjX5eeQnm54y/vriY0xLhmcEe/y5I6\npCAgIhLkwkJDeGRMb0rLy/nzzFWEh4YwdlB7v8uSOqKnD4qIiPfUwr6c060Fv//vSqanb/G7JKkj\nCgIiIgJARFgIT1zRlzNSm/PbN1YwY9lWv0uSOqAgICIi34oKD2XSlWkM6JDIHdOX8/5X2X6XJLVM\nQUBERL4jOiKUZ67uT++2zbj1laXMWb3D75KkFikIiIjI98RGhvHsNf3p3iqOm176knkZu/wuSWqJ\ngoCIiFQqLiqcqdcMoGNSLD9/Pp30Tbl+lyS1QEFARESOKCE2ghevH0hKfBTXPLuYFVl7/S5JapiC\ngIiI/KDkppG89POBxMeEc9Uzi1i7fb/fJUkNUhAQEZEqpcRH8/L1g4gMC2HslIVs3HXQ75KkhigI\niIjIUWmXFMNL1w+krNwxdvJCtu0t8LskqQEKAiIictS6tGjK89cOIK+whLGTF5Kzv8jvkqSaFARE\nROSY9Gwdz7NX9yd7XyFXTlnIvvwSv0uSaqh2EDCzRDObZWYZ3mvCEfqN8/pkmNm4gPb7zWyLmR04\nrH+kmb1qZplmttDMOgR8drfXvtbMzqvuOoiIyLFJ65DIpKv6sSHnINc8t4j84lK/S5LjVBN7BO4C\n5jjnUoE53vx3mFkicA8wEBgA3BMQGGZ6bYe7DtjjnOsC/Av4uzdWd2AM0AMYATxhZqE1sB4iInIM\nzkhN5tHLe7Nsy15ueGEJRaVlfpckx6EmgsBIYKr3fiowqpI+5wGznHO5zrk9wCwqfonjnFvgnKvs\nZtaB476vdwwWAAAUmklEQVQOnGtm5rVPc84VOec2AplUHiRERKSWjeiZwt9H9+KzjF384pVllJaV\n+12SHKOaCAItA36RbwdaVtKnNRD4TMssr+2HfLuMc64U2AckHe1YZjbezNLNLD0nJ+do1kNERI7D\nJWlt+cNPuvO/r7dz95tfUV7u/C5JjkHY0XQys9nACZV89LvAGeecM7N68RPgnJsETAJIS0urFzWJ\niDRW153ekX0FJTw6J4O46HB+/+OTqNiJK/XdUQUB59zQI31mZjvMLMU5l21mKcDOSrptBYYEzLcB\nPq7ia7cCbYEsMwsD4oHdAe2BY+mh2SIiPrt9aCp5BSVMmbeRxNgIJpzdxe+S5CjUxKGBt4FDVwGM\nA2ZU0ucDYLiZJXgnCQ732o523IuBj5xzzmsf411V0BFIBRZVcx1ERKSazIw//qQ7o3q34qEP1vLS\nwm/8LkmOQk0EgQeAYWaWAQz15jGzNDObDOCcywXuAxZ7071eG2b2oJllATFmlmVmf/LGnQIkmVkm\ncAfe1QjOua+B6cAq4H/ABOecTlUVEakHQkKMhy45hXO6teD3/13JuysqOxdc6hOr+CO7cUtLS3Pp\n6el+lyEiEjQKisu4cspClmft5Zmr+3NGarLfJQUdM1vinEurqp/uLCgiIjUuOiKUKVf3p3NyE254\nYQnLt+jxxfWVgoCIiNSK+Ohwnr92AElNIrj62UWszzlQ9UJS5xQERESk1rSIi+L5awcSGmJcNWUR\n2/cV+l2SHEZBQEREalXH5rE8d80A9hWUcNUzC9mbX+x3SRJAQUBERGpdz9bxTLqyH5t25XPd1HQK\ninWxV32hICAiInXitC7NeWRMb77cvIdbX/lSzyWoJxQERESkzvzo5BT+fGEPZq/eye/eWkkwXMJe\n3x3VLYZFRERqylWndmBnXhGPz82kRVwkdw4/0e+SgpqCgIiI1Lk7h3clZ38Rj32USYumkVx5age/\nSwpaCgIiIlLnzIz7L+rJ7oNF/PHtr0luGsmInil+lxWUdI6AiIj4Iiw0hMcu70vvts24bdoyFm/K\n9bukoKQgICIivomOCGXKuP60aRbN9VPTydix3++Sgo6CgIiI+CoxNoKp1w4gIiyEcc/o7oN1TUFA\nRER81zYxhmev7k9eYSlXP7uIvMISv0sKGgoCIiJSL/RsHc/EsX3J3HmAm15cQnGpbjhUFxQERESk\n3jgjNZm/j+7F55m7+e0bK3TDoTqgywdFRKReGd2vDdn7Cnj4w3W0ahbFr8/r5ndJjZqCgIiI1DsT\nzu7C1r2F/GfuelLioxk7qL3fJTVaCgIiIlLvmBn3jezBjrxC/jhjJSnxUZx7Uku/y2qUdI6AiIjU\nSxU3HOpDj1bx3PLyUr7K2ud3SY2SgoCIiNRbsZFhTLk6jcTYCK6dupgtufl+l9ToKAiIiEi91qJp\nFM9d05+ikjKueW4x+/J1j4GapCAgIiL1XmrLpky6Ko1vdh/khhfTKSot87ukRkNBQEREGoRBnZJ4\n+JJTWLAhl7vf+Er3GKghumpAREQajJG9W7N5dz7/mLWOtokx3D6sq98lNXgKAiIi0qDcck4XvsnN\n599zMmiXGMPofm38LqlBq9ahATNLNLNZZpbhvSYcod84r0+GmY0LaL/fzLaY2YHD+t9hZqvMbIWZ\nzTGz9gGflZnZMm96uzr1i4hIw2Nm/PWikxncJYm73lzBF+t3+V1Sg1bdcwTuAuY451KBOd78d5hZ\nInAPMBAYANwTEBhmem2HWwqkOed6Aa8DDwZ8VuCc6+1NF1azfhERaYAiwkJ44op+dEiK5cYXlpC5\nc7/fJTVY1Q0CI4Gp3vupwKhK+pwHzHLO5Trn9gCzgBEAzrkFzrnswxdwzs11zh26WHQBoP0+IiLy\nHfHR4Tx7TX8iwkK55rnF7D5Q5HdJDVJ1g0DLgF/k24HK7v/YGtgSMJ/ltR2t64D3A+ajzCzdzBaY\nWWXBQ0REgkSbhBgmj0sjZ38RP38+ncISXVZ4rKoMAmY228xWVjKNDOznKq7jqNFrOcxsLJAGPBTQ\n3N45lwb8DHjEzDofYdnxXmBIz8nJqcmyRESkHundthn/urQ3X27ey69eW055uS4rPBZVBgHn3FDn\nXM9KphnADjNLAfBed1YyxFagbcB8G6/tB5nZUOB3wIXOuW/39zjntnqvG4CPgT5HqHuScy7NOZeW\nnJxc1deJiEgDdv7JKdx1fjfeWZHNP2et87ucBqW6hwbeBg5dBTAOmFFJnw+A4WaW4J0kONxrOyIz\n6wM8RUUI2BnQnmBmkd775sBgYFU110FERBqBG87sxJj+bXl8biavL8nyu5wGo7pB4AFgmJllAEO9\necwszcwmAzjncoH7gMXedK/Xhpk9aGZZQIyZZZnZn7xxHwKaAK8ddpngSUC6mS0H5gIPOOcUBERE\npOLRxaN6MrhLEne/uYJFG3P9LqlBsGC4RWNaWppLT0/3uwwREakD+/JLuGji5+w5WMxbNw+mQ/NY\nv0vyhZkt8c6p+0F61oCIiDQq8THhPDOuPw64dqqeVlgVBQEREWl0OjSP5amx/diSm8/NLy+hpKzc\n75LqLQUBERFplAZ2SuJvP+3F55m7ueftr/W0wiPQQ4dERKTRurhfGzJ3HuDJT9aT2qIJ1wzu6HdJ\n9Y6CgIiINGq/Oe9ENuQc4L53VtGheSxnn9jC75LqFR0aEBGRRi0kxPjXZb3pdkIct768lHU79ICi\nQAoCIiLS6MVGhjF5XBrREaFcN1UPKAqkICAiIkGhVbNonr4qjZ15Rdz04pcUl+pKAlAQEBGRINK7\nbTMevLgXizbl8vv/fqUrCdDJgiIiEmRG9m5N5s4DPPZRJl1bNuX6Mzr5XZKvtEdARESCzu1Du3Je\nj5b89b3VzF1b2YNzg4eCgIiIBJ1DVxKceEIct728lMydwXslgYKAiIgEpZiIiisJIsNDuG5qOnvz\ni/0uyRcKAiIiErRaN4vmqSv7kb23kAkvf0lpED6TQEFARESCWr/2ifzlop58nrmbv7y72u9y6pyu\nGhARkaB3aVpb1m7fz5R5GznxhKZcPqCd3yXVGe0REBERAe4+vxtndk3mjzNWsmhjrt/l1BkFARER\nESAsNITHLu9D24QYbnpxCVl78v0uqU4oCIiIiHjio8N5elwaxaXljH9+CQXFZX6XVOsUBERERAJ0\nTm7Co5f3YfX2PH79+vJGfxtiBQEREZHDnN2tBb8d0Y13VmTzxMfr/S6nVikIiIiIVOKGMzsxsncr\nHv5wLbNX7fC7nFqjICAiIlIJM+Pvo3vRs1U8v3x1GZk7D/hdUq1QEBARETmCqPBQnrqyH5FhIYx/\nPp19BSV+l1TjFARERER+QKtm0Uwc24/Nufn8ctpSysob18mDCgIiIiJVGNAxkXsu7MHctTn8c9Za\nv8upUbrFsIiIyFEYO7Adq7bt4z9z19M9JZ4f90rxu6QaUa09AmaWaGazzCzDe004Qr9xXp8MMxsX\n0H6/mW0xswOH9b/azHLMbJk3XV/VWCIiIrXJzPjzhT3p1z6BX722nDXb8/wuqUZU99DAXcAc51wq\nMMeb/w4zSwTuAQYCA4B7AgLDTK+tMq8653p70+SjGEtERKRWRYSFMPGKvsRFhzH++SXszS/2u6Rq\nq24QGAlM9d5PBUZV0uc8YJZzLtc5tweYBYwAcM4tcM5lH8P3HXEsERGRutAiLoqJY/uRva+AW19p\n+CcPVjcItAz4Rb4daFlJn9bAloD5LK+tKqPNbIWZvW5mbY91LDMbb2bpZpaek5NzFF8nIiJydPq2\nS+DekT35LGMXD33QsE8erDIImNlsM1tZyTQysJ+ruBlzTcWimUAH51wvKv7qn1pF/+9xzk1yzqU5\n59KSk5NrqCwREZEKlw9ox88GtuPJT9bzzoptfpdz3Kq8asA5N/RIn5nZDjNLcc5lm1kKsLOSbluB\nIQHzbYCPq/jO3QGzk4EHj3csERGR2vKnC3qwdvt+fv3aCrq0aEK3E+L8LumYVffQwNvAoTP3xwEz\nKunzATDczBK8E/uGe21H5IWKQy4EVh/vWCIiIrUlIiyEJ67oS5OoMG54YUmDvPNgdYPAA8AwM8sA\nhnrzmFmamU0GcM7lAvcBi73pXq8NM3vQzLKAGDPLMrM/eePeZmZfm9ly4Dbg6qrGEhER8UPLuCgm\nXtGXrXsK+OW0pZQ3sJMHrbE/ZxkgLS3Npaen+12GiIg0Yi/M38QfZnzNbeemcsewrn6Xg5ktcc6l\nVdVPtxgWERGpAWMHtefifm14dE4GsxrQY4sVBERERGqAmfGXUT05uXU8d7y6jA05DeOxxQoCIiIi\nNSQqPJSJY/sSFmrc8MISDhaV+l1SlRQEREREalCbhBgeu7wv63MO8Js3VlDfz8VTEBAREalhp6c2\n51fnnci7K7KZMm+j3+X8IAUBERGRWnDTWZ05r0dL/vb+Guav3131Aj5REBAREakFZsbDl5xCh6QY\nbn3lS7L3FfhdUqUUBERERGpJ06hwnrqyHwXFZdz80pcUl5b7XdL3KAiIiIjUoi4tmvLgxaewdPNe\n7n93ld/lfI+CgIiISC37ca8Urj+9I1Pnf8N/l271u5zvUBAQERGpA789vxsDOiZy15srWLM9z+9y\nvqUgICIiUgfCQ0N4/Gd9iIsK58YXlpBXWD+eVKggICIiUkdaNI3iP1f0JWtPAb+avrxe3GxIQUBE\nRKQO9e+QyN0/OokPV+1g0qcb/C5HQUBERKSuXTu4Az8+OYUHP1jLgg3+3mxIQUBERKSOmRl/v7gX\n7ZNiuOXlpezMK/StFgUBERERHzSJDOPJsf04WFTKLS8vpaTMn5sNKQiIiIj4pGvLpjww+mQWbcrl\n7WXbfKkhzJdvFREREQBG9m5NctNITu2U5Mv3KwiIiIj47LTOzX37bh0aEBERCWIKAiIiIkFMQUBE\nRCSIKQiIiIgEMQUBERGRIKYgICIiEsSqFQTMLNHMZplZhveacIR+47w+GWY2LqD9fjPbYmYHDuv/\nLzNb5k3rzGxvwGdlAZ+9XZ36RUREgl119wjcBcxxzqUCc7z57zCzROAeYCAwALgnIDDM9Nq+wzl3\nu3Out3OuN/AY8GbAxwWHPnPOXVjN+kVERIJadYPASGCq934qMKqSPucBs5xzuc65PcAsYASAc26B\ncy67iu+4HHilmnWKiIhIJaobBFoG/CLfDrSspE9rYEvAfJbXViUzaw90BD4KaI4ys3QzW2BmlQUP\nEREROUpV3mLYzGYDJ1Ty0e8CZ5xzzsxcTRXmGQO87pwrC2hr75zbamadgI/M7Cvn3PrDFzSz8cB4\nb/aAma2t4dqaA7tqeMxgo21YfdqG1adtWH3ahjWjprdj+6PpVGUQcM4NPdJnZrbDzFKcc9lmlgLs\nrKTbVmBIwHwb4OOjKY6KIDDhsHq2eq8bzOxjoA/wvSDgnJsETDrK7zlmZpbunEurrfGDgbZh9Wkb\nVp+2YfVpG9YMv7ZjdQ8NvA0cugpgHDCjkj4fAMPNLME7SXC41/aDzKwbkADMD2hLMLNI731zYDCw\nqlprICIiEsSqGwQeAIaZWQYw1JvHzNLMbDKAcy4XuA9Y7E33em2Y2YNmlgXEmFmWmf0pYOwxwDTn\nXODhhpOAdDNbDswFHnDOKQiIiIgcJ/vu71k5WmY23jv8IMdJ27D6tA2rT9uw+rQNa4Zf21FBQERE\nJIjpFsMiIiJBTEFAREQkiCkIAGY2wszWmlmmmVV2m+RIM3vV+3yhmXUI+Oxur32tmZ13tGM2NjW9\nDc2srZnNNbNVZva1mf2i7tbGP7Xxs+h9FmpmS83sndpfC3/V0v/PzczsdTNbY2arzezUulkbf9TS\nNrzd+395pZm9YmZRdbM2/jjebWhmSd6/fQfM7PHDlulnZl95yzxqZlYjxTrngnoCQqm4D0EnIAJY\nDnQ/rM/NwJPe+zHAq9777l7/SCrugLjeG6/KMRvTVEvbMAXo6/VpCqxrzNuwtrZjwHJ3AC8D7/i9\nng1xG1JxC/XrvfcRQDO/17UhbUMq7ia7EYj2+k0HrvZ7XevpNowFTgduBB4/bJlFwCDAgPeB82ui\nXu0RqHjoUaZzboNzrhiYRsUzFAIFPlPhdeBcL4mNpOISxyLn3EYg0xvvaMZsTGp8Gzrnsp1zXwI4\n5/YDqznKW1M3YLXxs4iZtQF+DEyug3XwW41vQzOLB84EpgA454qdc3tpvGrl55CKG9hFm1kYEANs\nq+X18NNxb0Pn3EHn3DygMLCzVdy0L85VPKPHAc9T+fN9jpmCwNE9C+HbPs65UmAfkPQDyx738xUa\nqNrYht/ydpn1ARbWYM31UW1tx0eA3wDlNV9yvVMb27AjkAM86x1emWxmsbVTfr1Q49vQVdwR9mFg\nM5AN7HPOfVgr1dcP1dmGPzRmVhVjHhcFAanXzKwJ8AbwS+dcnt/1NDRm9hNgp3Nuid+1NGBhQF9g\nonOuD3CQSh65LkdmFXeVHUlFqGoFxJrZWH+rkkMUBCqehdA2YL6N11ZpH2+3Vjyw+weWPZoxG5Pa\n2IaYWTgVIeAl59ybtVJ5/VIb23EwcKGZbaJi9+Q5ZvZibRRfT9TGNswCspxzh/ZIvU5FMGisamMb\nDgU2OudynHMlwJvAabVSff1QnW34Q2O2qWLM4+P3SRV+T1Sk/Q1UJNVDJ3X0OKzPBL57Usd0730P\nvntizAYqThKpcszGNNXSNjQqjoE94vf6NeTteNiyQ2j8JwvWyjYEPgNO9N7/CXjI73VtSNsQGAh8\nTcW5AUbFsfFb/V7X+rgNAz6/mqpPFvxRjdTr9warDxPwIyrOSl8P/M5ruxe40HsfBbxGxYkvi4BO\nAcv+zltuLQFncFY2ZmOeanobUnHWrANWAMu8qUZ+6OvzVBs/iwGfD6GRB4Ha2oZAbyDd+3n8L5Dg\n93o2wG34Z2ANsBJ4AYj0ez3r8TbcBOQCB6jYI9Xda0/ztt964HG8uwNXd9IthkVERIKYzhEQEREJ\nYgoCIiIiQUxBQEREJIgpCIiIiAQxBQEREZEgpiAgIiISxBQEREREgtj/B0vJOO0mU4ZbAAAAAElF\nTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10dfa9f60>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(figsize=(8,5))\n",
    "a_min, a_max = 0, 0.01\n",
    "a = np.linspace(0, 0.01, 100)\n",
    "ax.plot(a, f(a))\n",
    "ax.hlines(0, a_min, a_max, linestyles=':')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.0032034184388007483"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a_star = scipy.optimize.brentq(f, 0.002, 0.004)\n",
    "a_star"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "CE_a = 20"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def CE(x, p, a):\n",
    "    x = np.asarray(x)\n",
    "    p = np.asarray(p)\n",
    "    return -np.log(p @ np.exp(-a * x)) / a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "19.999999999999982"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "CE([20], [1], a_star)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "24.600322886059406"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "CE_b = CE([50, -10], [0.6, 0.4], a_star)\n",
    "CE_b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "22.204239619220147"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "CE_c = CE([50, 0, 10], [0.4, 0.3, 0.3], a_star)\n",
    "CE_c"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}