Skip to content

Instantly share code, notes, and snippets.

@pv
Created June 6, 2016 19:43
Show Gist options
  • Save pv/3a704c54b6b628762fc3afdfbfb55cc3 to your computer and use it in GitHub Desktop.
Save pv/3a704c54b6b628762fc3afdfbfb55cc3 to your computer and use it in GitHub Desktop.
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"[<matplotlib.lines.Line2D at 0x7fe6f876cb50>,\n",
" <matplotlib.lines.Line2D at 0x7fe6f86e1190>]"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXIAAAECCAYAAADjBlzIAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3X901PWd7/HnmwRRzABCQmn4pbRaMKmwdW04KhpqrYLL\nsndrrd6r26IFt16q9vQeLe0eCWzPsbB7Dq3tWkTXdhd1W731VPaUXlu6RlNUsO3FS/jRdBsRcNaS\ngJKxUUzgc//4zpDJZGYyvzLf+c68HufkMJl85ztvvoR33nl/39/P15xziIhIcI3yOwAREcmPErmI\nSMApkYuIBJwSuYhIwCmRi4gEnBK5iEjAKZGLiAScErmISMAVJZGb2dVmNq4Y7yUiUmnyTuRm1mhm\nlubr1cBCYEK+7yUiIkPllcjNrAl4GRhtZlVmttbMlprZqtg2zrl+4L084xQRkRTySuTOuR3Akein\ny4HDzrlngGNmdkPcpikrdhERyU8he+TzgV3Rx68C18Hp1soM4NwCvpeIiERVF3BfU4BI9HEEmAyn\nWyu3pXqRmWn5RRGRKOdc1h2MQlbkR4Ga6OMaoDvTFzrnSuZj9erVvsegeBSP4qnMeHJViEQe++nx\nLDA3+vgi4OeZ7qClpYXW1tYChCIiEjytra20tLTk/Pp8p1YuBmqBTwGbgRlm9hlgGvBYpvtpaWmh\nubk5n1BERAKrubk5r0SeV4/cOfcbIBT31H3RP5/KZj+r/m4V13zympJI5qUQQzzFk57iSU/xpFcq\n8bS2tubVlbB8+jKFYGZu7vfm0rasjdCY0PAvEBEpU2aG8/lkZ872du1lT9cev8MQEQmkQo4f5mzi\njokcazrmddZFRCpMWbRWet7rUVtFRCperq2VkkjkfscgIlIKAt0j1xy5iFSyfOfIVZGLiJSIQFfk\nIiKSu5JI5GqtiEglU2tFRKRMqLUiIlKhSiKRRyLDbyMiIsmVRCJfsEDJXEQkVyWRyHfvbmHz5la/\nwxAR8UVZnOycO9fR1gYhXaUvIhUs0Jfo9/Q4JXERqXiBTuR+xyAiUgo0figiUqFKIpGv+rtVurJT\nRCpWeZzs1K3eRESC3VrRrd5ERHJXEon8wroLaahr8DsMEZFAKonWim71JiKi8UMRkcALdI9cRERy\np0QuIhJwJZHIdYcgEalkZTFH7ncMIiKlQD1yEZEKpUQuIhJwSuQiIgGnRC4iEnBK5CIiAadELiIS\ncErkIiIBl1ciN7NGM0s582hmf2lml5vZnfm8j4iIpJZzIjezJuBlYLSZVZnZWjNbamar4jZb5Jz7\nFXDSzLROrYjICMg5kTvndgBHop8uBw47554BjpnZZxI3j36IiEiBFapHPh/YFX38KnBd9PE2M5sP\njHbO7U314vVPbWP9U9sIH40UKBwRkcpRXaD9TAFiWTgCfADAOffj6HMvp3vxvXuuBuC+lxrp/PqL\n1E/STSZERDJVqER+FKiJPq4BurN6dav3xwm3hw3/vJl/uOeOAoUlIlK6WltbC7Lya16rH5rZa8BH\ngBuBM5xzj5jZcuA959zmDPfhWO09HtOjilxEKleuqx/mXJGb2cVALfApYDOwJnqScxrQks2+mttv\n5OwPj+f6m5p54b9+ysH9B6k7q46ud7uYMW4GGBw8nvy5+MezJ81m4XkLdf9PEQmUfCvzkliPvPGf\nGmnvai/I/iaPncy6T67jzNFnpkz+yZ7TDwER8Vugb75szYab6eA8X0NhWmgaKz62gmnjpmX820Dt\n2bU0TW3SDwARyVmsIl+zZk1wE/mYuxs5MaEwFbkfpoems6Z5TdLkH+mLcP2c66kfV+93mCJS4gJd\nkVef3cOGH+1k5od76e3r9ZLg2Dq6eruYMX4GQMrnYo93/3E367avo9/1+/Z3SaXaqnlg0QNE3o+c\nruzVyhGRRIFO5B/4wGoefbSZxYub89pXuCfM0/ufpm5sHZA6+Sc+13G0g02/3cThnsN5vX+2Els5\nSu4ilaksWis9PY6Qz7krciJC64FW9nfv9zX5x5L7+ZPOV/9dpMIEuiL3O4ZcRU5E2PnGTnr7BlpC\n8cm/5owa1r6wljffeTPn95hy9hTuu/I++k71qdcuUuYCnchXr15Nc3Mzzc3NvsYyEmLJvru3+3Rl\nf6jnUE7VfLVV8/CSh5k2fpoqdZEyUhatFb9j8ENiKyfb5D5z3ExevO1FVegiZSTQFbnfMZSK+ORe\nc0YNd/3sLvpcX8rtJ545kQeve5DF5y9WdS5SBgKdyMu5tZKP2BROlVWl7bWfO/5ctt+6XdW5SECp\ntVIhYtX6yp+t5ODxg0O+Xju2lu3LtnNB7QU+RCcihRDoitzvGIIkciLC1t9v5ZanbxnSdqkeVc2e\nL+5RMhcJKCXyChPuCfNE+xPc/6v7OfbusdPPTxgzgR1f2KFkLhJAgU7kX/3qamprmzn//GZ6e+Hg\nQairg64umOGNZQ95rrYWmprw/UIiv3V0d9DwvQb6Tw0sTaDKXCRYyqJH3tjoaM9hzawpU+D+++HM\nM+HIEbj+eqivwPN9Hd0dND3SxNsn3j793LRx09h7x15Ns4gESKAr8qoqx8mT+e9r9Ghob4cLKrAQ\nTazMq6yKZ29+lqtmXeVzZCKSqUAn8lwr8mTGj4dvfcur0lO1aMq1NdPR3cHCf11IOBIGYNaEWbQt\na9NYokhABDqR9/Q4du6E3l4y6pHX1MBdd0Ff6mtlMjJ9OqxfX15JfVvnNq597FpOOu9XnFnnzGLX\n7bvUYhEJgEAn8lxiCIfh6ae95P7uu3D33XD8eO5xzJkDO3YEP5lHTkSYt3EenW93nn5u2y3b1GIR\nCYBcE/mokQgmWy0tLVnfeLS+HlauhM9+Fj7/edi50+uR52rfPm8fQRcaE2LDNRsGPdfb1+tTNCKS\nidbWVlpaWnJ+fWAr8mTiq3RI3qLp6IBNm+BwkrWpvv51GDcueVtn9mxYuDAYFXvkRISmh5vYd3Qf\noF65SFBUXGslH5GIV313d8M993iJOhNTpsDzzwdjKiaxV/7hiR/mtyt+q165SAlTIs/Rtm1w7bVk\nPP5YVQXf+Y73wyBWuZditR45EeFjD32M/3zrPwGNI4oEgRJ5jiIRuPRS8h5/nDYNVqyAuXNLJ6mH\ne8Is+MECOt/yTnw21jXy4m0vqioXKVFK5HmIb7Uk66vv3g3f/GbmVfu558L27aVxlWl8i6Xaqmm7\ntY350+b7HZaIJKFEPsI6OrxKOxzObPuJE+Eb3/B+SPh50jTcE+ZD3/kQ7/W/x5nVZ/KHL/1BJz1F\nSlSuibx6JIIpRxdcAPv3D63cDx1KPgVz7BjccUfyfV14Ibz8cnGS+evHXz992X7fyT72de9TIhcp\nMyVRkQf9DkGRCLS2wq5d3vIAx46l394MnnkGliwpQmwnIlz6z5fS3uWdBFCfXKT0lMXqh37HUEgd\nHdDYOPzyAdXV8MADA62XkVwmYFvnNhY9voj+U/2MHjWaF5a9oD65SAlSj7yExC5MqqkZfNL0qafg\nxz9O/brp02HNGi+5F3JJ3siJCAu+v4A9R/Zw7oRzef7zz6u9IlKClMgDIByG886D998ffttYxd7X\nV5ikHu4Jc+W/XMmBtw/QUNdA27I2tVdESkyg11qpFPX13ihjJkm5v987WXrXXV5F//d/n/nETDKv\nH3+dA28foP9UP3u79rKna0/uOxORkqKK3Afxc+vp1n5JVF0Njz0Gixdn30uPnfTc172PObVzdMJT\npARp/DBAQiG4Ku5K+bvv5vR67EePwqpV8OabQ1/X3w833uhV9M89l9uaL2ZZf4+ISIlTRV6CYhX7\noUNw773e/UgTVVXB449nXp2/dOglrvjBFfSf6te6KyIlqmRPdprZ1cAO51xPiq8rkacRicDWrV7V\nnqxKz/SGGJonFyl9RTnZaWaNlsXv5mZWDSwEJmQbmHhCIe/mGR0d8MMfen3yePv2wcaNXsJPu58x\nITZcu4HqUd4Ofnf0dzrhKVImMk7kZtYEvAyMNrMqM1trZkvNbFWq1zjn+oH3ChBnxYsl9D17vIuH\n4t1zD8ybN/xUS9PUJhrqGqi2amaOn8mMcTNGLmARKZqME7lzbgcQ69YuBw47554BjpnZDQBmdoOZ\nbY5+/KuZXQbo7FoBXXABvPrq0BHGzk5YsCB9ZR4aE2Lrf9/Kueecy4HjB1j8xGIiJ4Yp5UWk5OU6\nRz4f2BV9/CpwHYBz7knn3C3Rj78BdgAzgHPzDVQG1NfDK6/ArFmDn+/s9NZ8SUfz5CLlJ9fxwylA\nrJSLAJOTbRRtrdw23M7ibzoa5MWziqm+3luka+NGr7USc/fd0Nyc+uRn4+RGGuoa2HNkj9orIj6L\nLZaVr6ymVsysE5gNfB/4lnPulWjvfKVz7pacAtDUSl4iEa8/3undBIiqKnj22cFz6ol0ub5IaSrW\nJfqxN3gWmBt9fBHw82zfOF5LS0tBfipVolAI2toG2iwnT8KXvpS+V672ikhpaW1tHdSZyFbGFbmZ\nXQy0AjcBPwXWALuBRqAl17JaFXlhbNkCS5cOfL5tW+qqXJfri5SmEb9E3zn3GyD+f/t90T+fyvZN\npfDGjh38eW/v8K/R5foi5aEkVj9UayV/TU3eVZ4xX/ta6vZK+5F29h/dT/+pfl0YJFICitZaGSlq\nrRTOtm2waJG3uFa6k5660YRIaQr0euSqyAujqQlmz/YenzzpjSImq8p1YZBIaVFFLoNkWpXHr4ao\n+3iKlIZAV+RSOJlW5bELg7TuikjwlUQiV2ulcEIh2LBhYJXE3/3OW2hryHZqr4iUDLVWZIhIBC69\n1Fvids4cePHF5Jfsv3ToJRZ8fwEn3UmqrZq2W9vUXhHxkVorMsRwY+Izx89kdNVoAKqrqtVeEQko\nJfIy1N4O+/d7JzxTtVbAu1S//1Q/ACdPneRgz8EiRikihVISiVw98sJqbISGBhg9Gi680HucdDud\n8BQpCeqRS1KxGzg7502ypFrWVishipQO9chliK98xZspT3fnoPiVENuPtLPzjZ3FDVJE8qZEXqba\n273eeH8/7N2buk/eOLmR2ZO8wfOT7iR3/5+7NYYoEjAlkcjVIy+8TPvkoTEhNly7gepR3uD5vu59\nqspFikw9ckkp0z55bH3y9q52ABrrGrVGuYgP1COXpDLpk8eq8iqrAmB/934tbSsSIErkZSzTPjnA\nhbUX6uIgkYBSIi9jmfbJYfDFQX0n+9jXva9IUYpIvtQjL3ORiFeJNzSk7pGD+uQipSDQPXJNrYyc\nUAjmz0+fxEHTKyJ+0tSKDCsS8frljY2qykVKWaArchk5kYg3sXLFFeknV0DTKyJBpURe5rKZXIHB\n0yujRo1i4pkTixCliORDibzMZTO5At70St/JPgDeP/k+1/3bdbpkX6TEKZGXuVAI2trghRe8P4c7\n6dk4uZHzJpx3+vPX3npNJz1FSpxOdsoQ4Z4wC36wgM63OgGYNWEWbcvaqB9X73NkIuUt0Cc7NX5Y\nWurH1fPQXzx0+qRn59udLPjBArVYREaIxg8lL6lGEyMnIszbOI/OtztPP7ftlm1cNesqH6IUqQyB\nrshlZEUi8NJLg0cPIxHYsgUuugguvxzmzYNweODroTEhNlyzYdB+evt6ixSxiGRDFXmZi82Rxy7T\nb2vznr/0Uq8Sj1dfD6+84v0JXlXe9HAT+456666oVy4yslSRS1LJ5sjb22FfkjWxwmH4+McHKvPQ\nmBAPLH5gUK/84498nHBPeOiLRcQ3SuRlLtkc+cyZUF2dfPs33hh8BWjT1KZB44hvRN7QiU+REqNE\nXuaSzZG//rpXoQOMGgVr18KUKQOv6eyEjRu9ZB4aE+L5zz/PB2s+OPD1tzppPdBa3L+IiKSkHnkF\nivXN9+71qvS2NnjuOVi6dPB2c+bAjh1e8t+yfwtLfzSwQX1NPa8sf0X9cpECyrVHrkRegSIRL0Gb\neT3xUMh7bt48rxqPt2ULLFmSfBxx1jmz2HX7Lq2OKFIgI36y08wazSzjNzCzK8xsrZk9lm1QMnJi\n1fiiRd79PGNiLZj6hAL7b//WO/kZGhOibVmbWiwiJSijRG5mTcDLwGgzq4om6KVmtirNy3Y65+4D\nXi9EoFIY6VZDTBw/hMGTLPXj6tl43cZB+1vx7yvYsn+LTn6K+Cjj1oqZdQKzgVuBU865TWZ2O/CW\nc+5JM7sBWBLd3AEPAX8E+p1zB9LsV62VIkrWH09cSCuWvN94Y+C5WbNg1y7gjKEtFoA5k+awY/kO\ntVlE8lDMOfL5wK7o41eB6wCcc086526JfvwNMAW4C7jVzKbn8D4yAjJZDbG+HnbuhA8OdFHo7ITW\n1oEWy9TQ1EGv2Xd0n9osIj7JJZFPAWK/R0eAyck2cs792Dn3Jefcfc65Q7kGKIWXyX086+u9EcR4\nK1Z4Jz9DVs/OL+ykPjS4oX7rM7fS0d0xAhGLSDopLgtJ6yhQE31cA3TnG0T8ql/Nzc00Nzfnu0tJ\nEFsca+ZMb458uPt3Aixc6LVUYpMsb77pjSh6Y4n1vPKFV7jk4UsIv+Nd6dn9bjcN32tgzxf3cEHt\nBSP8NxIJvtbW1oKs/JpNj/w14CPAjcAZzrlHzGw58J5zbnPOAahHPuJiffH2du8Kz/7+gXVXhkvm\nyfrlAOvXexMtERdm7kNz6e4d+Hk+6axJPPqXj7LwvIXqmYtkYUR75GZ2MVALfArYDMwws88A04C8\nxwu1HvnIik2qnDwJ772X+f07YaBfnjiWeM890NTktVm2L9tO9aiBX+6OvnuUpT9aykcf/KjWZRHJ\ngNYjl2HFr4BYXe0l9FQTK6mEw3DJJYOXuoWByvy/TnRw+fcvp6u3a9DXa8+qZfut29VqEclAoK/s\nXL16tXrjIywS8RL5jBlw8KDXWsk0icekSuazZkWXx60Jc8kjlxCODN6gelS1+uYiacR65WvWrAlu\nIvc7BslcqmQeu5iImjCXPXoZB44fGPR19c1FhhfoitzvGCQ74bDXqklcl6W2FrZvhw/OjLD191u5\n+emb6Xf9g7aZHprOdxd/VwldJIlA31hCJzuDpb7eu8pz/frBz3d3e733rT8JsXjmZ9lzxx7qxtYN\n2uZQ5BBLf7SU2d+drZlzkSid7BTfRCLe5Eqyuw3V13tL49ZMSd43B6iyKh7/68dZfP5iVeciqLUi\nPolEYOtWuPnmgZtVxFRVweOPw59dEeaaHw3tm8eo3SLiCXQi19RK8HV0eFeCJp4EBZg+HdZ/O8LZ\n5+/kaN8hlv/78iG9c/BuVvHc557TdItUHE2tSMlIV53DQLuFiR1c9v3LBl0NGjOKUbQ0t3Dbn92m\nuw9JxQl0Re53DFJYHR1w2WXeyc9EVVXwyCNwZl2Yezsu42DPgaT7qLZqHl7yMNPGT6NpapNaLlIR\nAp3I1VopP+Gwl8wPHEi9zdRZEW5fu5PaWYe48xfJ2y2gHrqUP7VWpGRFIt46LYcOwfLlydst4K3I\nuHlrBzf87ErefOfNlPubHprO+qvXU3t2rap0KUuBrsj9jqFSxJayzWQJ20JLdzIUYPJkWP4/I5yY\nuZXNR+7mj72pEzqoSpfypEQuacUvnJXpErYjEUNrK6xc6a33kkrt1AjLWlr5t56VHI6k2RCYcvYU\n7r/qfvXSpSwokUtaL70EV1zhtTdGj/Zu9TZ/vj+xZJrQ68+L8MkvtPLLM1byxp/SJ3RQ60WCL9CJ\nXCc7R14mN132I6ZMEvqosyLc9JWdfPBD3fzw2D3DVukAU2umcmfTndx80c0aY5SSp5OdkrHYUra5\nLGE7kmIJfccOWLcu9UlR8Kr0qz+3k/HTD/Hg4RX0u760+662ah5Y9AB9p/q4fs71SupS0gJdkfsd\ng5SOcBieeAK+/W04fDj9tlPOD3PR/3iC31R/m6P9w2wMVFHFVy//Ko2TG9V+kZKkRC5lJXaV6J13\nwpEjw2x8RoQJDTu5aVk3T/XcQ3f/8K0XGDhROnHsRMaOHqvELr5TIpeyFGu77NoFmzYNX6VzRgSm\n7mTC9EP0f2IV75B+jDHe1Jqp3H7x7cydMldjjeKLQCdyneyUTOST1Ps+cS9/YrjSfsC00DRWfGwF\n08ZNI9IXUX9dRpROdkpFyuYEKeAl9ZmtMGk/E8bW0HPFXZwi/YnSeLGTppH3I8wYN0M9dhkRga7I\n/Y5Bgi0chqefhpoaWL06/SjjaTVhmPM09NUw8b+t5pjLrK8eb3poOmua19D1bhezJ81WO0bypkQu\nwsD6Lt3d3rIA2bRgGNMNkzqwpk24ccNPwSSKb8couUsulMhFksi6rw6D2jCcqIGr1kJN5idN4ym5\nSzaUyEWGEavWe3vh6FFYtQrezCQ/J1Ts/PkmmHAYsv7v5lFyl1SUyEWylFMbBgYSe1UvVPcOJPdz\nsm/HxMSS+/mTzgeDg8cPKsFXICVyGZafy9gGQXzF3tsLu3dnOBEDg9sxf6qDcYfyrtwBJo+dzLpP\nruPM0Wdy8PhB6s6q0zhkGQt0Itcc+cgrhWVsgyg2EVNX531esOQOeSX4+HHIurPq6Hq3S2ORAaY5\ncslIKS1jG3Txyf3dd7PotcPwyR3ySvAwsPSAqvjgCXRF7ncMlaAUl7EtF/G99oMHvQR/6FCOUzLH\nZ4ADJu+GK9ZBVVzpn2eCT1XFjz1jrNaaKRFK5DKsUl3GtlzFRh/37/eSe1eXd9HSXXdBXyYXlcYu\nWnqnzkvuEw5645DNa2Fcwq8AeSZ5GFhrJv6EqxJ9cSmRiwRE/JWoXV1eks/4ilQYPA454aDXojm7\ny0vy190FVQk/JQqQ5CF5olfbprCUyEUCLHFiJqcWDWRXxUPBknyqto1GKbOjRC5SppK1aApWxf+p\nDhauhnPS7KhAyT7xQqj4RK+JG48SuUiFSTzJOmOG9/zBg1n24hOTfOyEayzhFynRw+CFyJJV9uWe\n8JXIRWSQxF58LNFnNQsP6RN9urZNvAIme/BGLO+78r7TywrH9+yDPFNflERuZo3Ankwzr5ldCFwG\nTHTOrUuxjRK5SJElXugU68kXtG0TS/jx69OkU+BkD6kTfqn270c8kZtZE/BLYCJwElgN/Aa40Dl3\nf5rXXQ40Ouc2pvi6ErlIiUl18jWnRA9DL4TKJdHHjEDCT7YUQmJbpxgTOsWqyDuB2cCtwCnn3CYz\nux14yzn3pJndACyJbu6Ah4BXgI3OuVtT7FOJXCRg0vXnc5qZT1yILFllH2vlJBuxTGUEkn7i3aKS\ntXVy7ekXO5FvAh50zu00s/nAF51zn0uy/RLgLWCcc25rin0qkYuUqXQtnBkzvCUO7r0XjmR+O9WB\nEcv3aoYm+uFm6lMZgYQfb86kOexYvmPYZF7sRL4F+LJzbp+ZNQD/6JxblO2bR/epRC5SweLHK5NV\n9rGEn9WaNpA+4adbCmE4OSb9bbds46pZV6XfdY6JvDq3kDgK1EQf1wDdOe4HgJaWltOPtQqiSGUJ\nhWDJEu8jnU9/euiaNvHTOEMvoqqHV1am3+n+z8Jv7xh6EVWqtk4mEzqJDHDezUyYNfhLsVUP85Vt\nRf4a8BHgRuAM59wjZrYceM85tzmnAFSRi0iBJbuIKj7hxx5ndUMRGH7mfshJ3EdgXBi6Z/OdeS+y\ncoXPrRUzuxhoBW4CfgqsAXYDjUBLrtlY65GLiJ/STegU5CTu5D2ccbyB1/aHqE8x7KL1yEVEimy4\nk7jxz0UiXlsoVRKPV+weeUG1tLSoIheRwKivh5XDtN+zkW+vXBW5iEiJyLUiHzUSwWSrpaWlIGdu\nRUSCqLW1ddD0XrZUkYuIlIhAV+QiIpK7kkjkaq2ISCVTa0VEpEyotSIiUqGUyEVEAq4kErl65CJS\nydQjFxEpE+qRi4hUKCVyEZGAK4lErh65iFQy9chFRMqEeuQiIhVKiVxEJOCUyEVEAk6JXEQk4Eoi\nkWtqRUQqmaZWRETKhKZWREQqlBK5iEjAKZGLiAScErmISMApkYuIBJwSuYhIwJVEItccuYhUMs2R\ni4iUCc2Ri4hUKCVyEZGAUyIXEQk4JXIRkYBTIhcRCTglchGRgFMiFxEJuIwTuZk1mllW841m1mBm\nX80+LP+U2oVJiic9xZOe4kmv1OLJVUaJ3MyagJeB0WZWZWZrzWypma1K85oLgHeAMYUJtThK7R9W\n8aSneNJTPOmVWjy5qs5kI+fcDjM7Ev10OXDYOfeMmU0xsxucc0+a2Q3AEiB2mebvgSPAxWZW65zr\nLnj0IiKSWSJPMB94MPr4VeCLwJPOuSeBJxM3NrMPKImLiIycjNdaMbNOYDawBfiyc26fmTUA/+ic\nW5RzAGZaaEVEJCqXtVZyqciPAjXRxzVAXtV2LkGLiMiAbMYPYwn3WWBu9PFFwM8LGpGIiGQl06mV\ni4Fa4FPAZmCGmX0GmAY8NnLhiYjIcHxfj7yUmFmNc+4dn967EdhTKouzp4rHz2MUBDo+6en4pBa9\nTucs51xvtq/15crOTGfRixTLfWb2ezPbA5wTjeuvihnXcHP6xT5e8fFEP48/RjXFjMfMQmb2uJn9\nwcweNbNqP49PYjzR5/w8PuPN7Ftm9nMz+19+f/8kxhN9zrfjExfXx8zse9F/v29E3/vL0a8Nea5Y\n8UQfP2RmHcBOoD+XePy6RP/0LDpwLDqDXnRmdjZwFtDgnGsArovG9ZNixuWc24E3cw/Jj01Rj1d8\nPInHyDn3ZpHj+RRwK/AR4M+BryV5bz/imY13jcQC/D0+5znn7gauicbm9/fPoHhK4PsHMxsPLMS7\nOPHrwAvR964zs48nPDfZzC4pVjxmNgN4DZjtnLvEOfd+LvH4lcjnA7uij1/FS6B+uACYB7xhZstK\nJK5kMTT5GFfiMUoV40h5xjl3wjnXD+wFzk/y3sU8PrF4+qLxvIuPx8c5F3ufS4FHUrx30Y5PQjwP\n4//3D8Cngf+NN7ARfyz+H/78/4rFA16Bcg1w2MyujT6X9fHJZfywEKYAkejjCDDZjyCcc/8XWGRm\nHwH+A+8f1u+4Eo/NB/Culi12XAaDjtFsYJuZ/SxJjCMWTzSBY2ZjgMNAAz4en8R4nHO/xsfjE43l\nPCBWiBzxgZv+AAAB20lEQVTC5++fuHiagIudc74dHzP7NPATIBR9Ktn/r2TPjXQ84/DOUf4C+IWZ\nXQ780MzOzSUevyrybgo4i54v59zv8H5Cjsb/uJLN6Rd0dj9Dg05yOuf24x2jmfjz73cDcB+lc3xu\nAFbHPvHz+DjnXnPOfQHYAZxK8t5FPT4J8Xw0+pxfx2cZ3m8qm4BPADOSvHcxj08snoeAhbEeuHPu\nV8CvgUnkcHz8SuQlMYserapixgD/gn9xpZrTfzbJc8WIyyDpMdpb7HjMbDHws+jZ/J/j8/GJjyfa\n44zx5fjEeRv4IaXx/ROLJxz3edGPj3PuL5xzf43Xl/8P4Jv4eHzi4lkRjeef4r78jnPuj7nE48v4\noZkZsAbYDTQCLX6M3ZnZP+D9hN4CvA5s9yMu8+b0W4GbgJ8mxhDdrGhxJcRzJXHHyDn3q2L++5nZ\njcA64DhQBXwHqMen4xONZz1ekqrGawv8Ev+OTwve9Rw/Bk4AzyW+d3RTv+JZhI/fP3FxzcT7DeoO\nvGPyG2CWc26dmZ2F99veb2PPjWQsCfGchfe9/TLwa+dcey7xaI5cRCTgdIcgEZGAUyIXEQk4JXIR\nkYBTIhcRCTglchGRgFMiFxEJOCVyEZGAUyIXEQm4/w8eJWt/PrALwgAAAABJRU5ErkJggg==\n",
"text/plain": [
"<matplotlib.figure.Figure at 0x7fe6fb898c90>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"%matplotlib inline\n",
"\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"from scipy.optimize import nonlin\n",
"\n",
"def F4_powell(x):\n",
" A = 1e4\n",
" return [A*x[0]*x[1] - 1, np.exp(-x[0]) + np.exp(-x[1]) - (1 + 1/A)]\n",
"\n",
"xin = [-1, -2]\n",
"\n",
"iterates = []\n",
"\n",
"def callback(x, f):\n",
" iterates.append(np.array(x))\n",
"\n",
"try:\n",
" nonlin.newton_krylov(F4_powell, xin, f_tol=1e-3, maxiter=20000, verbose=0, callback=callback)\n",
"except nonlin.NoConvergence:\n",
" pass\n",
"\n",
"iterates = np.asarray(iterates)\n",
"\n",
"k = np.arange(len(iterates))\n",
"\n",
"plt.semilogy(k, abs(np.asarray(F4_powell(iterates.T)).T), '.')"
]
}
],
"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
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment