Skip to content

Instantly share code, notes, and snippets.

@AngelicosPhosphoros
Last active October 15, 2018 21:47
Show Gist options
  • Save AngelicosPhosphoros/09cb556cb24344717c889e8e05403180 to your computer and use it in GitHub Desktop.
Save AngelicosPhosphoros/09cb556cb24344717c889e8e05403180 to your computer and use it in GitHub Desktop.
HW from MIEM: Derivatives
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# © Timur Khuzin a.k.a Angelicos Phosphoros"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 1. Epsilon"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"My and np 2.220446049250313080847263336181640625e-16 == 2.220446049250313080847263336181640625e-16\n",
"This is True\n"
]
}
],
"source": [
"import numpy as np\n",
"\n",
"def machine_eps():\n",
" eps = 1.0\n",
" while(1.0+eps!=1.0):\n",
" old, eps = eps, eps/2.0\n",
" return old\n",
"\n",
"print(\"My and np {0:.55} == {1:.55}\\nThis is {2}\".format(machine_eps(),\n",
" np.finfo(float).eps,\n",
" machine_eps()==np.finfo(float).eps))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 3. First derivative"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"from matplotlib import pyplot as plt"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"def deriv(x, func, step=1e-10):\n",
" return func(x + 1j*step).imag/step"
]
},
{
"cell_type": "code",
"execution_count": 108,
"metadata": {},
"outputs": [],
"source": [
"# check for exp\n",
"\n",
"row = np.arange(-5,5,step=0.01, dtype=float)\n",
"lg = np.logspace(1e-10, 1e-1)\n",
"strict_deriv = np.exp(row)\n",
"lg_der = np.exp(lg)\n",
"my_deriv = deriv(row, np.exp)\n",
"lg_my_deriv = deriv(lg, np.exp)"
]
},
{
"cell_type": "code",
"execution_count": 109,
"metadata": {},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD8CAYAAAB5Pm/hAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3Xl8VPW9//HXJxtZ2EOQJUCQTcIeI3spigouReqvXtfWX6s/2p+1anutrXq7Wntt++vV3uqtRW1dLnWptwp1raKoIKABASFhkzWIEMMqELLM5/dHBhowEMjM5Ewm7+fjMY/JnHNmznvC4/HO4TtnvsfcHRERSVxJQQcQEZHYUtGLiCQ4Fb2ISIJT0YuIJDgVvYhIglPRi4gkOBW9iEiCU9GLiCQ4Fb2ISIJLCToAQKdOnTwvLy/oGCIizcrixYs/dfechraLi6LPy8ujqKgo6BgiIs2KmW06me00dCMikuBU9CIiCU5FLyKS4OJijL4+VVVVlJaWUlFREXSUhJKenk5ubi6pqalBRxGRJhK3RV9aWkqbNm3Iy8vDzIKOkxDcnfLyckpLS+ndu3fQcUSkicTt0E1FRQXZ2dkq+SgyM7Kzs/W/JJEWJm6LHlDJx4B+pyItT1wXvYhIIvvZ2XN5/ddLYr4fFf1xlJeXM3z4cIYPH06XLl3o3r37kceVlZUx2+/48eNZunTpKT3nzjvv5M0334xRIhGJhYrdFfxs7gTmvbQ35vuK2w9jg5adnX2kcH/605/SunVrbr311qO2cXfcnaSk4P5eVldXc/fddwe2fxFpnLVzNuP054yhaTHfl47oT9G6desYPHgw3/rWtygoKGDLli20b9/+yPqnnnqK66+/HoDt27dz6aWXUlhYyMiRI1m4cOHnXu/AgQNcdtllDB06lCuuuOKoD0pffvllxowZQ0FBAZdffjn79+8HIDc3l7vuuotx48bx3HPPcc011/D888/z97//nauuuurI819//XW+/OUvx+pXISIRWP1uOQBnjMuO+b6axxH9LbfAKQ5nNGj4cLjvvkY9tbi4mD//+c88+OCDVFdXH3e7m266idtuu43Ro0ezceNGLr74YlasWHHUNvfffz8dOnRg+fLlfPDBBxQWFgKwY8cO7rnnHubMmUNmZiZ33303v/vd77jjjjsAyMrKYv78+QDMmjULgMmTJ3PDDTdw8OBBMjIyePrpp7n88ssb9R5FJLZWLTsEQL+zc2O+r+ZR9HGmT58+nHXWWQ1u9/rrr7N69eojj3ft2nWkhA97++23ue222wAYMWIEgwYNAuDdd9+luLiYsWPHAlBZWcn48eOPPK++Ak9LS+O8887jxRdf5JJLLuHVV1/lvkb+MROR2CpZl0rP5FKyOqvoa8VZWWVlZR35OSkpCXc/8rju0Iu7895775GWduIxuPpOeXR3pkyZwhNPPNFghrouv/xyHnnkETIzMxkzZsxxtxORYBXvyGZQx21A7Iu+wTF6M/uTme0wsxX1rLvVzNzMOoUfm5n9p5mtM7PlZlYQi9DxJCkpiQ4dOrB27VpCoRDPPffckXXnnnsuDzzwwJHH9Z1NM2HCBGbOnAnAsmXLWLlyJQBjx47lrbfeYv369QDs37+ftWvXNphn0qRJLFq0iEceeUTDNiJxqqayhlUHe5Hfc3+T7O9kPox9FJhy7EIz6wGcB2yus/gCoF/4Nh34Q+QR49+vfvUrpkyZwqRJk8jN/edf5wceeID58+czdOhQ8vPzeeihhz733BtvvJHy8nKGDh3Kvffee2SM/rTTTjtS1sOGDWPs2LGsWbOmwSwpKSlccMEFvPbaa1x44YXRe5MiEjUb55VSQQb5Q5KbZH9Wd9jhuBuZ5QEvuPvgOsueBe4CZgGF7v6pmf0RmOvuT4a3WQ1MdPdtJ3r9wsJCP/bCIyUlJQwcOPDU3o2cFP1uRYL1wo/f40t3jWTBQysYff3ghp9wHGa22N0LG9quUadXmtlUYKu7LztmVXdgS53HpeFlIiISVlx0AICB5/dokv2d8oexZpYJ3AmcX9/qepbV+18GM5tO7fAOPXv2PNUYIiLNVvGaFLolbaNdz65Nsr/GHNH3AXoDy8xsI7UfGS8xsy7UHsHX/ROVC3xc34u4+wx3L3T3wpycBq9tKyKSMIo/6UB++3qrMSZOuejd/UN37+zuee6eR225F7j7J8Bs4Gvhs29GA3saGp8XEWlJPOQU7+9Ffo/PmmyfJ3N65ZPAAmCAmZWa2XUn2PwlYD2wDngIuCEqKUVEEsSWRR+zn9bkD2q6KcMbHKN39ysbWJ9X52cHvh15LBGRxFQ8ZxvQnfwx7Zpsn5rU7Dg0TbGIxELx+7VDNvmTm+aMG2guUyAEQNMUi0gsFK9KprOVkd2v6U5C0RH9KdI0xSISieKP25PfrrRJ99ksjujjbJZiTVMsIo3iIaf4s55cPTjKhdaAZlH08UbTFItIY2xd/Al76Ep+E89A0iyKPt66StMUi0hjfPhyKdCVoRPaN7htNGmMPkKaplhETtbyBbWfsw35Ul6T7ldFHwWaplhETsbyklR6JG+lfa+mO4ceTnKa4ljTNMVNS79bkWAMzVhDz7a7eWH7yKi8XkynKRYRkVNT+VklJRW9GdrnQJPvW0UvItIEVr2ykWpSGXLmiU/OiIW4Lvp4GFZKNPqdigTjwzfKABh6bucm33fcFn16ejrl5eUqpihyd8rLy0lPTw86ikiLs3xxFalU0v+8Xk2+77g9jz43N5fS0lLKysqCjpJQ0tPTjzozSESaxocbsshP30Bq5oAm33fcFn1qaiq9e/cOOoaISFQsL8/l7F4fAU1f9HE7dCMikih2frSLraGuDB14/LmxYklFLyISYx++sAmAIWNaB7J/Fb2ISIwtf2cPAEMvDObzsZO5ZuyfzGyHma2os+w3ZrbKzJab2XNm1r7OutvNbJ2ZrTazybEKLiLSXCxfYWRbOV2HnxbI/k/miP5RYMoxy14DBrv7UGANcDuAmeUDVwCDws/5LzNLjlpaEZFmaHlpR4a024wlNd0FwetqsOjd/W1g5zHL/uHuhz9VWAgc/v/IJcBT7n7I3TcA64DoTOogItIMVR2oYtn+vhT02RNYhmiM0X8DeDn8c3dgS511peFln2Nm082syMyKdK68iCSqVS9v4BDpFIxMDSxDREVvZncC1cDMw4vq2azer7a6+wx3L3T3wpycprtIrohIU1ryyg4ACi7sEliGRn9hysyuBS4GJvk/5ykoBXrU2SwX+Ljx8UREmrcPimrIZD/9z88LLEOjjujNbArwA2Cqu9edc3M2cIWZtTKz3kA/4L3IY4qINE9L1rdneJuPSE4L7ryUkzm98klgATDAzErN7DrgfqAN8JqZLTWzBwHcfSXwDFAMvAJ8291rYpZeRCSOhapDfLD3dArydgWao8GhG3e/sp7Fj5xg+7uBuyMJJSKSCNbN2cRn9GbEmcF+N1XfjBURiZElL24DoGBK089BX5eKXkQkRpYsrCSNQ+RfFOxMvCp6EZEYWbKuLUMyPyKtddNfPrAuFb2ISAx4yPlgdx4FPT8NOoqKXkQkFjYv2MpO70jBiKCTqOhFRGJiyexSAArOyw44iYpeRCQmFs+vIJlqhlxyetBRVPQiIrHwXkkbhmSsI6NjRtBRVPQiItEWqg7x3s6+jMrbEXQUQEUvIhJ1a/6xkT20Y+SoYC40ciwVvYhIlL33fO2kvaOmdQ04SS0VvYhIlC1aUEMb9nLGBcF+I/YwFb2ISJS9tz6HwvbBTk1cl4peRCSKKnZXsOxAX0adEdw1Yo+lohcRiaKlz66jijRGTkgPOsoRKnoRkSha9FI5AKMuzws2SB0qehGRKFq0JIXuSdvoVhDcxcCPpaIXEYmi97bmMqrLpqBjHOVkrhn7JzPbYWYr6izraGavmdna8H2H8HIzs/80s3VmttzMCmIZXkQknny6upyPqnsxcmhF0FGOcjJH9I8CU45Z9kNgjrv3A+aEHwNcAPQL36YDf4hOTBGR+LfoyfUAjJrcPuAkR2uw6N39bWDnMYsvAR4L//wYMK3O8se91kKgvZnFx1fDRERibP4/9pNCFSOv6R90lKM0doz+NHffBhC+P3zl2+7AljrblYaXfY6ZTTezIjMrKisra2QMEZH4MW9lBwqyVpPZKTPoKEeJ9oex9c3g4/Vt6O4z3L3Q3QtzcnKiHENEpGkd2nuI9/YOYPwZwV868FiNLfrth4dkwveH5+IsBXrU2S4X+Ljx8UREmofFT67hEOmMnxQ/X5Q6rLFFPxu4NvzztcCsOsu/Fj77ZjSw5/AQj4hIIps/u/aLUuO+1ifgJJ+X0tAGZvYkMBHoZGalwE+Ae4BnzOw6YDNwWXjzl4ALgXXAAeDrMcgsIhJ35i3JoF/qBjoPio8ZK+tqsOjd/crjrJpUz7YOfDvSUCIizYmHnPnb+zK1bzEQf0Wvb8aKiERo9SsbKPdsxo8LOkn9VPQiIhGa9/RWAMZfkRtwkvqp6EVEIjRvvpFjZfQ7Ly/oKPVS0YuIRGjeph6M6/IRlhQfFwM/lopeRCQCW4u28VF1L75wVnxNZFaXil5EJAJvPlI7kdnZV8TP/PPHUtGLiETgzTkhOtguhl0WXxOZ1aWiFxGJwBsb8pjYZTVJKfFbp/GbTEQkzm2cV8rG6h6cPfZQ0FFOSEUvItJIb/5pAwBnX90t4CQnpqIXEWmkN+bWnj8/6JK+QUc5IRW9iEgjeMh5c1Mfzs5dG7fnzx+mohcRaYR1czaxNdSVs8dXBx2lQSp6EZFGeOPRzQCcc22PBrYMnopeRKQR3nwnhW5J2+J2fpu6VPQiIqeoprKG10sHcG5e/M5vU5eKXkTkFC2euYpyz2bKhc2jQptHShGROPLqf5dhhDjvxgFBRzkpERW9mX3XzFaa2Qoze9LM0s2st5ktMrO1Zva0maVFK6yISDx45f2OFGaV0GlAdtBRTkqji97MugM3AYXuPhhIBq4AfgXc6+79gF3AddEIKiISD3Zt2M3CfYOYPKIs6CgnLdKhmxQgw8xSgExgG3AO8Gx4/WPAtAj3ISISN+bcX0KIZKZc1THoKCet0UXv7luB/wdsprbg9wCLgd3ufvgbBKVA9/qeb2bTzazIzIrKyprPX0YRadleeaGaduxh1Nfzg45y0iIZuukAXAL0BroBWcAF9Wzq9T3f3We4e6G7F+bk5DQ2hohIk/GQ8+pHfTi3ezEp6SlBxzlpkQzdnAtscPcyd68C/gaMBdqHh3IAcoGPI8woIhIXiv/+EaU13Zh8TlXQUU5JJEW/GRhtZplmZsAkoBh4E/hKeJtrgVmRRRQRiQ+vPFwKwOT/2yfgJKcmkjH6RdR+6LoE+DD8WjOAHwDfM7N1QDbwSBRyiogE7u/vtGdI+hp6jqn3o8e4FdEgk7v/BPjJMYvXAyMjeV0RkXhTvnYn8/YM5vbx84D4vT5sffTNWBGRk/DSb4upIYWp13cOOsopU9GLiJyE2S8k0zXpE868+oygo5wyFb2ISAMO7T3EK1sHM3XAGpJSml9tNr/EIiJNbO7vP+Qz2jD18oygozSKil5EpAGzn9pPJvs55+YhQUdpFBW9iMgJeMiZXdKPyV0/JL19etBxGkVFLyJyAkv+sorSmm5MvTD+LwJ+PCp6EZET+Ov920mhiqm3Dwo6SqOp6EVEjsNDzjOLT+fcTkvp2KdD0HEaTUUvInIcS/6yig3VPbns4oqgo0RERS8ichzP/L522Gbavw0OOkpEVPQiIvXwkPPMkj6cl9O8h21ARS8iUq+iJ0rYWN2j2Q/bgIpeRKRef31gB6lUMu1HzfNLUnWp6EVEjhGqDvHMkr6cm7OMDr3bBx0nYip6EZFjzH/wQzbV5HLFtENBR4kKFb2IyDGeeGAPWXzGpT8fHnSUqFDRi4jUcXDnQZ5ZNYxLT19K6y6tg44TFREVvZm1N7NnzWyVmZWY2Rgz62hmr5nZ2vB98z4vSURalL/f9QF7aMfXvpkZdJSoifSI/nfAK+5+BjAMKAF+CMxx937AnPBjEZFm4Ym/pNA9aRtn3zIs6ChR0+iiN7O2wATgEQB3r3T33cAlwGPhzR4DpkUaUkSkKexYWcbLOwq45qzVJKclBx0naiI5oj8dKAP+bGYfmNnDZpYFnObu2wDC9/VeSdfMpptZkZkVlZWVRRBDRCQ6nvxR7QXAv3p7btBRoiqSok8BCoA/uPsIYD+nMEzj7jPcvdDdC3NyciKIISISOQ85f365CwUZJQy6pG/QcaIqkqIvBUrdfVH48bPUFv92M+sKEL7fEVlEEZHYe/+xYpZVDOD/TEu8EYZGF727fwJsMbMB4UWTgGJgNnBteNm1wKyIEoqINIE//vtOsviMq36dGOfO15US4fO/A8w0szRgPfB1av94PGNm1wGbgcsi3IeISEzt2byHp9YWcNWAxbTNnRB0nKiLqOjdfSlQWM+qSZG8rohIU5p521IO8EW+eUenoKPEhL4ZKyItmoecP87qwoiMEs68ZmDQcWJCRS8iLdr7jxWzvGIA3/xyGZZkQceJCRW9iLRo9/9iF63Zx5W/SrwPYQ9T0YtIi7Vt6XaeWj+SbwxbQtvctkHHiRkVvYi0WP/1nRKqSeGm/8gLOkpMqehFpEU6uPMgD84fzNQu79PnnF5Bx4kpFb2ItEgzv1vEp96JW37QKugoMaeiF5EWx0POfU93YVj6ar54U+JMR3w8KnoRaXFe/eViVh7qxy1XJ+4plXWp6EWkxbn7N6n0SN7KVfeNDDpKk1DRi0iL8s79y5i3dxjfn7aOtNZpQcdpEip6EWlR7v5pFZ2tjOtntIyjeVDRi0gLUvR4Ma+WF/Ld81eS0TEj6DhNRkUvIi3GL3+4l/a2mxseLgg6SpNS0YtIi7BkZgnPbRvNTV9YmtDTHdRHRS8iLcK/ffczOtguvvfEiKCjNDkVvYgkvHfuX8bLZWfxwynLaNezXdBxmlzERW9myWb2gZm9EH7c28wWmdlaM3s6fJlBEZFAeMi5407omvQJNz7ecs60qSsaR/Q3AyV1Hv8KuNfd+wG7gOuisA8RkUZ55RdFzNs7jB/9yxoyO2UGHScQERW9meUCFwEPhx8bcA7wbHiTx4BpkexDRKSxqiuq+cG/t6N3ymaue2h00HECE+kR/X3AbUAo/Dgb2O3u1eHHpUD3CPchItIoD3/jXT6s6M+vb9raYr4FW59GF72ZXQzscPfFdRfXs6kf5/nTzazIzIrKysoaG0NEpF67N+3hR0/lM6HdUv7Xb1ru0TxEdkQ/DphqZhuBp6gdsrkPaG9mKeFtcoGP63uyu89w90J3L8zJyYkghojI5/380g8o947c94f0FjFD5Yk0uujd/XZ3z3X3POAK4A13vxp4E/hKeLNrgVkRpxQROQWrX17P75eM47oB8xhx5RlBxwlcLM6j/wHwPTNbR+2Y/SMx2IeISL085Nxw9W6yOMAvnh0YdJy4kNLwJg1z97nA3PDP64GWebKqiATuiW/N541d4/nDlW9z2uAJQceJC/pmrIgkjE9Xl/O9hwcypvWHTH98fNBx4oaKXkQSxvcvLmaPt2XG4+kkpajeDtNvQkQSwj/+fTGPrvsC3x8zn8Ff7hd0nLiioheRZm/Xht18/d+6kd9qHT9+qWWfM18fFb2INHs3nr2SHaFOPP5wFent04OOE3dU9CLSrP31ewv4y6Zx/Ojs+Zx5jU6nrI+KXkSarc0LtvKt+wZwVtZKbn9hXNBx4paKXkSapaoDVVx+/k6qPIW/zGpNamZq0JHilopeRJql2yfMZ+FnQ3j4lhX0ndQr6DhxTUUvIs3O7DsX8dvFE7lh8Fv8y71jg44T91T0ItKsrJy1jmt+OZCCjBJ++86ooOM0Cyp6EWk2ytfuZOpXUslMquD5N9vpVMqTpKIXkWah6kAVXxm5ia3Vp/H8H3fQY1S3oCM1Gyp6EYl7HnJuKFjA3N0jeOibRYy+fnDQkZoVFb2IxL0ff/EtHl49gTvGzuWrD2pWylOloheRuPb7r7zFL+ZN5Lr+7/CLd74YdJxmSUUvInHrye+8y83/8wWmdV3Ig8vGtPhrvzaWil5E4tJTN73LV+8fyYR2y3myeDgp6VG5IF6L1OiiN7MeZvammZWY2Uozuzm8vKOZvWZma8P3HaIXV0Ragpk3zOfq349iXLsVvLCqr06jjFAkR/TVwL+6+0BgNPBtM8sHfgjMcfd+wJzwYxGRk/L49Hl87Q+jmdB+OS+t6UfrLq2DjtTsNbro3X2buy8J/7wPKAG6A5cAj4U3ewyYFmlIEWkZfvuluVz70HjO7rCMF9cOIKtzVtCREkJUxujNLA8YASwCTnP3bVD7xwDoHI19iEjiClWH+N6Zc7n1hYlclruAFzcOIrNTZtCxEkbERW9mrYH/AW5x972n8LzpZlZkZkVlZWWRxhCRZurApwe48vSF3LtkIjcNe4unNoyiVdtWQcdKKBEVvZmlUlvyM939b+HF282sa3h9V2BHfc919xnuXujuhTk5OZHEEJFmauO8Usb13Mxft4zm1xfO5b4lE0hK0cmA0RbJWTcGPAKUuPt/1Fk1G7g2/PO1wKzGxxORRDX3vqWcNSGdDQe78uLPFvP9FyfqPPkYieTE1HHAV4EPzWxpeNkdwD3AM2Z2HbAZuCyyiCKSSGoqa7jnonf4yevj6Z+2iednJ9F/8llBx0pojS56d58HHO/P76TGvq6IJK7S97fx1fM/Ye7uiVzR813+OH8wbXPbBh0r4WkwTERizkPOM999l2GjWvH+7n48ev08/rJhjEq+ieg7xSISU1uLtnHDlzYz+5OxFGYWM/NvGfSfrBkom5KO6EUkJkLVIf549dvkn5XJa58M4TcXzWVBeX/6T+4ddLQWR0f0IhJ179y/jFtuS2PJwQmc02EJM57Nps85E4OO1WLpiF5EombjvFL+pccCJnxnGDsq2zHzhvm8/ukI+pzTK+hoLZqKXkQitnnBVr6V/zb9v9CZF0uH8rOz57L6k/Zc9cA4nRsfBzR0IyKNtmXRx9xz/VoeWjEGyOH6QQu589F+dC+cGHQ0qUNFLyKn7P3Hirn3J7v466aRGJ34xsCF3PFIH3qOmRB0NKmHil5ETkrlZ5XM+vFi7ns4i3f3DaUNe/lOwXxu/s++9Bqngo9nKnoROaGVs9bxp7tKeWLJIMp8DL1TNnPvtLf4xu9H0DZ3YtDx5CSo6EXkc0rf38azv1zDk69l897+waTSk6ndF/P16zcy5Y4CktN6Bh1RToGKXkQA2DS/lL/9Zh1/nZPNgs+GAF0Zlr6ae6e9xdW/HETOwDFBR5RGUtGLtFCH9h7inQdX8vIze3llRS7Fh/oCuQzPWMXd583lK//ai/6TBwADgo4qEVLRi7QQh/Yeomjmat5+fifvfJDFW2X5HKCANA7xxY4ruG7yXL50Yx79zjsDOCPouBJFKnqRBOQhp/T9bSyZtYX33z7IOyvas2jPGRxiKAD5rdbx9SFFXHBpJhNvyCer85kBJ5ZYUtGLNHPVFdWsf2sLK9/YzuL5FSxencXi8jzKvBvQjSRqKMhczQ0FC5lwfgbj/3dfOg3oC/QNOro0ERW9SDOxZ/MeNiz4hFXzyylZVknJ+lYUl3Vi7aGeVNIb6E0y1QxK/4iL+qzizOHFnHleR4ZOO52szvlAftBvQQKioheJA6HqEGWrytm2cidbS/ayYeUBNm5wNmxLZ8Oudmys6MIu7wC0AyCJGk5P3cLAjtu5qNdWBg5JIX9cB4ZccjoZHfUBqhwtZkVvZlOA3wHJwMPufk+s9iUSbzzkHPj0ADs37KF84z7KN++n/OMKyj+pZvu2ENu2J7FtVyu27WvDxxUd2B7KoYYcIOfIa6RzkLy0j+ndbidj+iwnr6fT+4xWDBibTb9JPUlvnwfkBfQOpTmJSdGbWTLwAHAeUAq8b2az3b04FvsTiZaqA1VU7K7g4K4KDuw6xL7tB9j36aHa284q9u6sZt+eEPv2Ovv2wb79Sew7kMS+g6nsPtiK8oosdla1pjzUgUNkAVn17ifHyujaaidds/YxpHs5XXNW07Wb0TUvjW79WtN79GmcNjgHS+oD9GnS34Eknlgd0Y8E1rn7egAzewq4BFDRNzMeckLVoePeaqrqPK7xz93XVNYcvbzuuqo6y8K3mqoQVYdCVFXU1N4fClFVGaK60qmq51ZdDVVVR9+qa4yqaqiqMiqrjYrKZA5WJlNRlczB6hQqqlM4WJNGRU0qB0OtqAilcdDTqSCdGlKBVKBNg7+bNA7Rxj6jTfIB2qQcpF3aQfp02MnINp+Q3b6Gjh0hu3MS2V3SyO6eTsfcTLJ7t6VT/46ktT766F0klmJV9N2BLXUelwKjor2TV+8u4rs/7wCAh5c5tXNfu/9zDuxj19Vd71BnOzv6vqHXONXtGrl9zF7jOOscI0QSIZJwkqgdfUsm3qVQRSpVpFBNqv3zlpF0iIzkStKTq8hIqaZDegXd0j4jPTVERqsa0tOcjPQQ6a0gIwPS0yEj08jISqJNhxTadEylbU4r2nRqRZvTMmtvXVuT1roV0ArIDvqti5xQrIq+visN+FEbmE0HpgP07Nm4eTPa5rRicOcddV7Tj9q51Ulx7Lq66w+vO2rZiV6j7rJGbnf09g3nqO+5sXqNpKTP35KT61meDElJVrsu/HPt8qPvk1Ps6GWHt022f66rc0ttlURqevJRt5S0JFIzUj53S0lPITUzlZT0FCzp8NG4iNQVq6IvBXrUeZwLfFx3A3efAcwAKCwsPOqPwMkaM30IY6Y3NqKISMsQq0sJvg/0M7PeZpYGXAHMjtG+RETkBGJyRO/u1WZ2I/AqtYO7f3L3lbHYl4iInFjMzqN395eAl2L1+iIicnJiNXQjIiJxQkUvIpLgVPQiIglORS8ikuBU9CIiCc7cG/VdpeiGMCsDNgWdoxE6AZ8GHaKJ6T23DHrPzUMvd29w0qS4KPrmysyK3L0w6BxNSe+5ZdB7TiwauhERSXAqehGRBKeij8yMoAMEQO+5ZdB7TiAaoxcRSXA6ohcRSXAq+igxs1vNzM2sU9BZYs3MfmNmq8xsuZlfDXQFAAACL0lEQVQ9Z2btg84UC2Y2xcxWm9k6M/th0Hlizcx6mNmbZlZiZivN7OagMzUVM0s2sw/M7IWgs8SCij4KzKwHtRdC3xx0libyGjDY3YcCa4DbA84TdXUucH8BkA9caWb5waaKuWrgX919IDAa+HYLeM+H3QyUBB0iVlT00XEvcBvHXC4xUbn7P9y9OvxwIbVXEEs0Ry5w7+6VwOEL3Ccsd9/m7kvCP++jtvi6B5sq9swsF7gIeDjoLLGioo+QmU0Ftrr7sqCzBOQbwMtBh4iB+i5wn/Cld5iZ5QEjgEXBJmkS91F7oBYKOkisxOzCI4nEzF4HutSz6k7gDuD8pk0Ueyd6z+4+K7zNndT+d39mU2ZrIg1e4D5RmVlr4H+AW9x9b9B5YsnMLgZ2uPtiM5sYdJ5YUdGfBHc/t77lZjYE6A0sMzOoHcJYYmYj3f2TJowYdcd7z4eZ2bXAxcAkT8xzdBu8wH0iMrNUakt+prv/Leg8TWAcMNXMLgTSgbZm9t/ufk3AuaJK59FHkZltBArdvblNjHRKzGwK8B/AF929LOg8sWBmKdR+0DwJ2ErtBe+vSuRrH1vt0cpjwE53vyXoPE0tfER/q7tfHHSWaNMYvTTG/UAb4DUzW2pmDwYdKNrCHzYfvsB9CfBMIpd82Djgq8A54X/XpeEjXWnmdEQvIpLgdEQvIpLgVPQiIglORS8ikuBU9CIiCU5FLyKS4FT0IiIJTkUvIpLgVPQiIgnu/wOa72w6UOMLIgAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEACAYAAABCl1qQAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAIABJREFUeJzt3XlcVXX+x/HXl1wQLc0lTdHcUUBExYVNbqmlTpNplmUuLZpZZss0Ttv8qimbrJnMymxybbHUSjMzM60uKIiCS6QoYKZJigqugIrA+f1hOeZYFxQ4l3vfz8ejxyOu53LfUo/75vM9536PsSwLERGRP+JjdwAREXF/KgsREXFJZSEiIi6pLERExCWVhYiIuKSyEBERl1QWIiLikspCRERcUlmIiIhLKgsREXGpit0Bykr9+vWt5s2b2x1DRKRSWb9+fbZlWQ1cHecxZdG8eXOSk5PtjiEiUqkYY3aV5DgtQ4mIiEsqCxERcUllISIiLnnMOYvzOXXqFJmZmZw4ccLuKB7F19cXf39/qlatancUEakgHl0WmZmZXHrppTRv3hxjjN1xPIJlWeTk5JCZmUmLFi3sjiMiFcSjl6FOnDhBvXr1VBRlyBhDvXr1NK2JeBmPLgtARVEO9DMVcQ9WscW88QnkZuWW+2t5fFnYKScnh9DQUEJDQ2nUqBFNmjQ583VBQUG5vW5UVBSbNm0q1XOefPJJvv3223JKJCJl7chPRxjaYg23vR7BtLvL/zNmHn3Owm716tU786b9zDPPUKtWLR599NHfHGNZFpZl4eNjX28XFhYyceJE215fREpnzdvfM/T+Ouwu7MbEPk4eWRRd7q+pycIG27dvJzg4mHvvvZfOnTuze/du6tSpc+bP582bx6hRowDYt28fgwYNIiwsjG7dupGYmPg/3y8/P5+bb76ZkJAQbr311t+cT1i2bBnh4eF07tyZIUOGkJeXB4C/vz/PPfcckZGRLFq0iGHDhvHpp5+yZMkShg4deub5K1euZODAgeX1oxCRUigqKOKFa51Ej2kPwOrp23jiKweXVLuk3F/beyaLhx6CUi7NuBQaCq++ekFPTU1NZfbs2bz11lsUFhb+7nHjx49nwoQJ9OjRg507d3L99dezefPm3xzzxhtvcPnll5OSksLGjRsJCwsDYP/+/bz44ot8/fXX+Pn5MXHiRKZMmcITTzwBQM2aNYmPjwdg8eLFAFx33XXcd999HD9+nBo1ajB//nyGDBlyQX9HESk7PyfvZXifLL497GBI0wT+szqI2s2aVtjre09ZuJlWrVrRtWtXl8etXLmStLS0M18fOnTozBv5r+Li4pgwYQIAnTp1IigoCICEhARSU1OJiIgAoKCggKioqDPPO18JVKtWjT59+rB06VIGDBjA8uXLefUCC1FEysZnT67lrn+25rjVlll3ruKOGVEYn4q90MR7ysLN3vBq1qx55t99fHywLOvM12cvI1mWxbp166hWrdoffr/zXaFkWRZ9+/blvffec5nhbEOGDGHmzJn4+fkRHh7+u8eJSPk6fvA4j8as483NMYTW2Ma8T6oR0K/8z0+cj85ZuAEfHx8uv/xyMjIyKC4uZtGiRWf+rHfv3kydOvXM1+e7yqlnz57MnTsXgO+++44tW7YAEBERQWxsLDt27AAgLy+PjIwMl3l69erF2rVrmTlzppagRGyyeVEG3Rpn8ubmGB7p4iQxqwUB/Vralkdl4SYmTZpE37596dWrF/7+/mcenzp1KvHx8YSEhBAYGMj06dP/57njxo0jJyeHkJAQJk+efOacRcOGDc+84Xfs2JGIiAjS09NdZqlSpQr9+vVjxYoV9O/fv+z+kiLiklVsMfWWWMIGNWV/QR2WPZfMv5MdVL+suq25zNnLH5VZWFiYde79LLZu3Ur79u1tSuTZ9LMVKXvZaTnc7djOZ1nd6Vs/iTnfNqdhsMv7El0UY8x6y7LCXB2nyUJExA18/fIGQgJP8WVWKJNvjGXp3i7lXhSlobIQEbFRQW4Bj/Vw0mdCKLWr5JH4wY88tCgGnyru9fbsPVdDiYi4mYwVOxl6Yz7J+Q7uaRfHK7FdqHmFe1596F7VJSLiBaxii9l3raLTtfXZcbwRCyck8p+tPd22KECThYhIhTr042HujUllwe5oHHU28t5XjfDv2sPuWC5pshARqSBxr39HxzZ5LNzdlReudbJyXwj+Xa+0O1aJqCzKkbYoFxGAU/mneCrKydXjg6nuc4qEORk8vrxiNgAsK1qGKkfaolxEfvhmF7ffcIy1eQ7ubLOK1+I6UatRLbtjlZomCxtoi3IRz2cVW7x7z2pCe9UlLd+fBQ+vYVZ6dKUsCvCiycLNdijXFuUiHuzQj4cZG5PK/N1R9Ky9ifeWNaBZeLjdsS6K15SFu9EW5SKeKXbKJob/pQF7i06fxJ6wJLpSnZv4PV5TFu72fqctykU8y6n8UzzTJ55/JvSkVdWfSJiZQdeRDrtjlRmds3AD2qJcpHLLWLGTyCvSeSHBwV1t49n4U326jgy0O1aZUlm4CW1RLlL5WMUWM+84/Uns7fmN+fjRNcxIq7wnsf+ItiiXC6KfrXi7nIyDjI5JZ9HeHlxdZyPvftWo0nzA7mzaolxEpJyseHE9Ie1O8vnezrzU38nKAx0rZVGUhspCRKSEThw+wSNdnFz7eBdqV8ln7Qc7+OtSh9ttJ14evOZqKBGRi7F5UQa3Dy0m5YSD+zvE8tI3XfGr72d3rArj8XXoKedk3Il+puJNiguLmTLo9D2x956sy+dPJ/FGSoxXFQV4+GTh6+tLTk4O9erVO+/nEKT0LMsiJycHX19fu6OIlLs9G7K489pMvsqJ4U9XrGPm1y1oGOz6w7SeyKPLwt/fn8zMTA4cOGB3FI/i6+v7m8t7RTzRor8lMvrlNuRbgUy7LY4x70djfLz3l06PLouqVavSokULu2OISCWSm5XLgz03Misjmi5+qcz92JeAfj3tjmU7jz9nISJSUmve/p7QpjnMzojkiQgnCftaE9Cvpd2x3ILKQkS83qn8Uzwd4yRqTCCF1iXEvv49E+MdVKv1x3uyeROPXoYSEXElY8VOhg3MZV2egxEtV/Patx2o3Uzn5M6lyUJEvJJVbPH2sDhCr21ARn4T5j+UwDs/RFG7WW27o7klTRYi4nX2bznAqF47WLKvJ73rrmfO8sY0CYuwO5Zb02QhIl5lyd/XEdwBvtrXkVcGOFm+rxNNwjx7X6eyoMlCRLxCblYuj1y9genbetLRN41vPjhM8ECH3bEqDU0WIuLxEmdsplOzbGZsi2JCNydr9zUneGAbu2NVKioLEfFYZy6JHd2OguIqfPtqCpPWOqh+WXW7o1U6WoYSEY+UtmwHwwafIDnfwfCWq3ldl8ReFE0WIuJRrGKLqbfE0ql/I3Ycb8RHj6zhXV0Se9E0WYiIx9izIYu7rs1keU4MfesnMWtFM64MDbc7lkfQZCEiHuGjR9bQIawacTmBvHlrHF/sC+PK0IZ2x/IYmixEpFI7vOsI4xybmbszkq41t/DeRzW0S2w50GQhIpXW1y9voEOrPObt7M7TMU7i97fVLrHlRGUhIpXO8YPHeahTLL0ndMbP5yRr5qTzjNNBVb+qdkfzWFqGEpFKZf37Wxl+dzW2FsQwrkMsk77p6nX3w7aDJgsRqRQKTxTyfG8nPYa35kihH8tfWM/rKTEqigqiyUJE3F768h8ZflMe6/Ic3NosganftKduqy52x/IqmixExG0VFxYz9ZZYQvs2JCO/CfPGJ/Dhrgjqtrrc7mheR5OFiLilzKS93NV3DysOxtCvQRIzvmxK486654RdNFmIiFuxii0+uD+eDt1rEH+wHW8NjWNpVhiNOzeyO5pX02QhIm4jOy2Hsb3S+fjnSCIuTeGdRbVp3UsfsHMHmixExC18/n/rCA4sYvHPXXixr5O47CBa97rK7ljyC00WImKro5lHeaTXd8xMjybEN42v3jtMyGCH3bHkHJosRMQ2sVM2EdL8KLPTI3g83Mm6fc0JGdzW7lhyHioLEalwxw8e5+HOsTgeCqWKKWLVW6m8kKA72LkzLUOJSIVKeieVEfdUZ1tBDPd3iGXSyjBqXqFzE+5Ok4WIVIiC3AL+Hu0k/I625BbVYMWkDbyREkPNK2raHU1KQJOFiJS7zYsyGHF7IRuPOxjZajWvft2BOld1tjuWlIImCxEpN0UFRUzq56TLoGb8fKI+ix5by5ztUdS5SvfDrmw0WYhIuUhf/iN3DM5lTa6Dm5qsYdqKNjRo393uWHKBNFmISJkqLizmtZtOb/63Lc+fD8Yl8NFPPWjQvr7d0eQiaLIQkTKzc3Umd/75AM7DMfRvkMR0bf7nMTRZiMhFs4otpo+Io0N0bdYfbsWMkav4XJv/eRRNFiJyUTKT9jK6fyZfZvfk6jobmbWkAc2jou2OJWVMk4WIXBCr2OLde1YT3M2PuOxA3rg5lpUHOtI8yt/uaFIONFmISKnt23yAMX12sDgrishLU5izqDate8XYHUvKkSYLESmVBQ8nEBTiw5dZHfnX9U5itZW4V9BkISIlkp2Ww3290/koM4KuNbfwzjxf2l/vsDuWVBBNFiLi0qK/JRLUvphPM7vwwrVOErIDaH99K7tjSQXSZCEiv+vgD4d4oFcqH+yKpFONrax87xAdbnLYHUtsoMlCRM5ryd/XEdS2gAW7uvHs1U7WZremw026MZG30mQhIr9x6MfDPNR7M+/uiCLEN41lcw4TOsRhdyyxmSYLETlj6TNJBLc+ztwdPfh7tJOkAy0IHRJgdyxxA5osRITDu47wcO8U5myPJrh6Bp/NOEyXYQ67Y4kb0WQh4uW+fD6Z4Jb5vLc9nCcjnSTvb0aXYe3tjiVuRpOFiJc6vOsIf+mTwqyMaAKrb+fTtw8RNsJhdyxxU5osRLzQsn8kEdwynzkZETwe7mTD/qaEjQi0O5a4MU0WIl7k8K4jPNInhdm/TBOL/pNG15EOu2NJJaDJQsRL/DpNvHPWNNF1pKYJKRlNFiIeTtOElAVNFiIebOkzSQS1zOfdjHCeiNA0IRdOk4WIBzr7U9j63ISUBU0WIh5myd/XEdT6xJlPYetzE1IWNFmIeIicjIM8eO1W5u6MJMQ3jc9nHKLz7Q67Y4mH0GQh4gEWTkgkqF0h83d24+mY03s6db5d04SUHU0WIpXYga3ZjLsunQW7I+hUYyvL5xyi4y0Ou2OJB9JkIVIJWcUWCx5OIDAIFu0O4/nep+830fEW7RAr5UOThUglk5Wyn/v77WDhntP3wp499zBBAxx2xxIPp8lCpJKwii3eHxtPUGgVlu4JZVK/0/fCDhrQ2u5o4gU0WYhUAplJe7n3+t0s3R9JxKUpzJpfi4B+DrtjiRfRZCHixqxiixkjVxHUzY9v9wcxZVAscdlBBPRraXc08TKaLETc1I9xu7ln4H5WHozGUWcjMz6pS6trYuyOJV5Kk4WImykuLOb1wbF0iLmcxINtmXZbHF8f6Eira66yO5p4MU0WIm4kbdkORt16jNVHY+hbP4n/fNaYZuE97Y4loslCxB0UnihkUj8nHfs3ZsuxZrwzejVf7AujWXgTu6OJAJosRGyX8nE6d48sJDnfwaDGiUxd1pJGIVF2xxL5DU0WIjY5efQk/9fTSZebW/DT8QYseHgNn/zcg0YhV9gdTeR/aLIQscHamZu5635fUk86GN5yNZO/DKRem3C7Y4n8Lk0WIhUob38ej3RxEj4qkKOFfix9Jol3f4iiXpu6dkcT+UMqC5EK8s2/NxLSJJvJGxzcG7SaLTtr0f/prnbHEikRlYVIOTu86wij28XR69FO+BgL56ubeHNzTy7zv8zuaCIlprIQKUeLn1hLYMvjzEqLZEI3JylZDYl5MNTuWCKlphPcIuVg3+YDjO+fwYLdEXT0TWPJ9EN0GeawO5bIBdNkIVKGrGKLd+9ZTWDIJXy6uwsT+zhJymlJl2G6xalUbposRMrIztWZjLkxi69yooi4NIWZ82rRrr/D7lgiZUKThchFKiooYsqgWIKj65CQE8AbN8ey6mAw7fprG3HxHJosRC7ClsXbufv2E6zNi6F/gySmLW5Ms3BtIy6eR5OFyAU4efQkz17tpNONzfghvxFz74vn8yxt/CeeS5OFSCmteft7Ro2vQepJB0OviufVZQE0aB9pdyyRcqXJQqSEju05xviOsUSOCeJYYQ2WPpPE3J2RNGhf3+5oIuVOZSFSAl88m0RQs6O8kRLNAx1XseWny7RVh3gVLUOJ/IH9Ww7w8J/S+WBXJIHVtxP/5hbC79EJbPE+mixEzuPXD9e173AJH+3qytMxTjbsb0r4PR3sjiZiC00WIufY4fyJe286wIqDpz9cN/19PwJvcNgdS8RWmixEflF4opB//9lJ8NX1STzYhqlDTn+4LvCG1nZHE7GdJgsRYOOH2xg9qpj1+Q7+3HAtby5phn9XnZsQ+ZVbTxbGmJrGmHeMMdONMbfbnUc8T352PhO6Oek6tDWZJ+oz/6EEFu/phn/XK+2OJuJWXJaFMcbXGLPOGPOdMWaLMebZC30xY8wsY8x+Y8zm8/xZX2NMmjFmuzHmsV8eHgR8bFnWaOCGC31dkfNZ8eJ6gq/M5uUkB3cGJLB1ezVumRyB8TF2RxNxOyWZLE4C11iW1REIBfoaY3qcfYAx5gpjzKXnPHa+hd45QN9zHzTGXAJMBfoBgcBtxphAwB/Y/cthRSXIKuJSdloOI1qt5trHu1DVFOF8dRPTt/Xk8hZ17I4m4rZcloV1Wu4vX1b95R/rnMNigMXGGF8AY8xo4LXzfK844OB5XqYbsN2yrB2WZRUA84ABQCanC6NEWUX+iFVs8f7YeNq3t/hwR3eeinLy3f4rdec6kRIo0RuwMeYSY8wmYD+wwrKstWf/uWVZHwFfAvN+ObdwF3BLKXI04b8TBJwuiSbAQuAmY8w0YMnvZPuzMebtI0eOlOLlxNv88M0urmuwgeFvRdK65l42LtzJc6sc+NbxtTuaSKVQorKwLKvIsqxQTv+W380YE3yeY14CTgDTgBvOmkZK4nyLxJZlWXmWZd1pWdZYy7Lm/k62JZZl3VO7du1SvJx4i1P5p5jUz0lwryvOXA4bfyiI4IFt7I4mUqmUamnHsqzDgJPzn3eIBoKBRcDTpcyRCTQ962t/YE8pv4fIb6ybvYWwejt47EsH/a78jq1Jedw3LwafKlrRFCmtklwN1cAYU+eXf68B9Aa2nXNMJ2A6p88z3AnUNcY8X4ocSUAbY0wLY0w14Fbgs1I8X+SMY3uO8WBoLD3uak92wWUsnJDIwj09aBKmy2FFLlRJfsW6EvjWGJPC6Tf1FZZlfX7OMX7AzZZl/WBZVjEwEth17jcyxnwIrAECjDGZxpi7ASzLKgTGAcuBrcACy7K2XOhfSrzXZ0+uJbDZMV7/LpqxwatI/dGPgZN6uH6iiPwhY1nnXthUOYWFhVnJycl2xxCb/Jy8l/EDdrFwTw+Cq2cw/Y2T9Bj1P6fWROQcxpj1lmWFuTpOi7dSqRUXFvPmrbG071qTL/Z05J/XOdlwsLmKQqSMaW8oqbS+/ySde+44SWJuDL3rruetj+rT6hqH3bFEPJImC6l08rPzeayHk86DW7A970reHbOarw50ptU1V9kdTcRjqSykUlk+MZngK7OZtNbBiDaJbEvzYfhbUdrPSaScqSykUti3+QBDm8fT96kwqvkU4nx1EzPTo6nXpq7d0US8gspC3FpxYTFvD4ujXUhVPtkVxrNXO/nuQBPt5yRSwXSCW9zW5kUZjBl5nIRjPXHU2chbH9QmoJ/D7lgiXkmThbid/Ox8Hg930mlQc9Jym/DO6NV8kxNKQL+WdkcT8VoqC3Erv57AfjHRwfA2iWzbCiPe1glsEbupLMQt7N20j1ubJZw5gf3t5E3MSo+mfkA9u6OJCCoLsVlRQRFTb4mlXSdfPt3d+cwJbMdDOoEt4k50gltss2l+GmPuPsW6vBh6Xb6BafPr0qaPw+5YInIemiykwuVm5fKXMCdht7ZiZ35D3h8bz4rsTrTp09zuaCLyO1QWUmGsYotPH19Le/+jvLLewd3tEtj2Q1VufzNSJ7BF3JyWoaRC7IrP5IGbfmbJvu508E1n/tTviRjT0+5YIlJCKgspV6fyTzF5cDzPLusKXM6/rncyfn4kVf2q2h1NREpBZSHlJn5aCvc+XIPNJx0MaLSW1xb60yzcYXcsEbkAOmchZS4n4yCjAlYRdV8IRwr9+PTxtXy6tzvNwpvYHU1ELpDKQsqMVWwxZ9RqAgIs5qSH89euTlIzazPghe52RxORi6RlKCkTWxZvZ+yIPFYdjSLy0hSmzc6hw00Ou2OJSBnRZCEXJW9/Ho+HOwm98Sq2HGvKjJGriDsYTIeb2todTUTKkMpCLthnT64lqPGhM5v+pW21uHtOND5V9L+ViKfRMpSU2q74TMYP/pnPsroTXD2DVa+lEHVftN2xRKQcqSykxApyC3hlcAL/WN4NH+rw8p+cPLhAn5kQ8QYqCymR2CmbuO9vtUg96WDglYlMWdSMpt0ddscSkQqixWX5Q/s2H2BEq9U4Hgolv6g6S/6+joV7etC0e2O7o4lIBVJZyHkVFRTx5q2xBHSoxrwd3XgqysmWvfW4/h/d7I4mIjbQMpT8j6R3Uhl7H6zPP32fialz6xDQz2F3LBGxkSYLOePgD4cYGxRH9zvasedEXeaNT2BFdicC+rW0O5qI2ExlIRQXFjP7rlUEtCliemoED3ZaxbZdfgyZEqH7TIgIoGUor5fycTr33XWC+GPRRFyawpszcuh4S4zdsUTEzWiy8FJHM4/ycOdYOt/ckrTcxsy8YxWrDgbT8ZYAu6OJiBtSWXgZq9jig/vjCbjqOFM2RjOqfQJpGZdw12xt0yEiv0/LUF4k9bPtjBt5jG8PR9LFL5XFb2TT7U7d2lREXFNZeIHcrFz+cUMyk5MiudTkMu22OEbPieSSapfYHU1EKgmtO3gwq9jio0fW0N7/KC8nORjRJpG0LUXc+0FPFYWIlIomCw+VtmwH44YdYuXBcDr6pjF/6vdEjNHOsCJyYVQWHiZvfx7PD0ji34kR+FGP1wfHcu97kVTx1X9qEblwWobyEFaxxSd/XUP7xod5MdHB7a3XkvZ9AeM+ilFRiMhF07uIB0hf/iMP3J7DVzmnl5w+fD2FyLFachKRsqOyqMTOXnKqQV2mDIrlvrlachKRsqdlqEro16uc2jU+wouJDoa2Wkf69wWM/0RLTiJSPvTOUsls+2IHDwz/71VO815PIXJslN2xRMTDqSwqiWN7jvHcjeuZnBRJTV3lJCIVTMtQbs4qtvjwgQTaNc3l5SQHw9sk6ionEalwerdxY5sXZfDAXbk4D0fQucZWPnkthx6jdJWTiFQ8lYUbOvLTEZ65cROvb4yktjnKW0PjGDVbezmJiH20DOVGiguLeWf0ato2LzizfXh6GoyZq72cRMRemizcxPr3tzJubCGJuVH0qPU9X0zLpsswbR8uIu5BZWGznIyDPDlgM29vjaKByWH23asZ8VaEbkQkIm5F70g2KSooYtptcbQNgBlbI3iw0yrSd1bjjhlRKgoRcTuaLGyw+s0UHni0GpuO98RRZyOvz6pF8MAYu2OJiPwu/QpbgfZsyGJYi3ii7w8hu+Ay5j+UwDc5oQQPbGN3NBGRP6SyqAAFuQW81N9JQJeafLQzjCcjnWzbU5tbJkdgfIzd8UREXFJZlLNl/0iiQ92f+dsyB1c3TCX16yyeX+2g5hU17Y4mIlJiKotysv3rXfy54Tr6P90VgC+eTeKzrO60uuYqm5OJiJSeTnCXsdysXCYOTOaVxHCqUZeX+jt5cH4E1Wq1sDuaiMgF02RRRqxiiw/ujyegyTFeTHQwpGUS6Rvz+etSB9VqVbM7nojIRdFkUQY2friNB8YUEH8sks41tvLR5GwixugeEyLiOVQWF+HA1myeGpTK9G1R1DMHmT5iFXdOj9A+TiLicbQMdQFO5Z/itZtiaRtUhZnbTn/6OuPHqox6J1pFISIeSWVRSl+/vIFOdXfy4MIYwi7/gZTFO5m8IYY6V9W2O5qISLlRWZTQj3G7GdQ4kd4TOpNfVJ1Fj63lqwOdCbyhtd3RRETKnc5ZuJCblcuLg5P5V3wPLqEuE/s4eWRBD3zrNLM7mohIhdFk8Tt+vRS2XZNjTIx3MLj5etKTjvLEVw586/jaHU9EpEJpsjiP5HdTeXBcIQnHIunil8qCV7KJGBNpdywREduoLM6SlbKfJwanMScjkgYmhxkjV3HnjEjdX0JEvJ7eBYGTR0/y8p+ctO3oy/sZ3Xm0axwZP1Xn7jnRKgoREVQWrHhxPcH19zLhCwcxV2xjy1d7eGmdg8v8L7M7moiI2/D6Zah9u05QxRSx7Llk+j7Vze44IiJuyViWZXeGMhEWFmYlJyeX+nnFhcUUFRRR1a9qOaQSEXFvxpj1lmWFuTrO6ycLnyo+Oi8hIuKC3iVFRMQllYWIiLikshAREZdUFiIi4pLKQkREXFJZiIiISyoLERFxyWM+lGeMOQDsusCn1weyyzCOiEhlcZVlWQ1cHeQxZXExjDHJJfkEo4iIt9IylIiIuKSyEBERl1QWp71tdwAREXemcxYiIuKSJgsREXFJZSEiIi6pLERExCWvv/nR+RhjagJvAgWA07KsuTZHEhGxlddMFsaYWcaY/caYzec83tcYk2aM2W6MeeyXhwcBH1uWNRq4ocLDioi4Ga8pC2AO0PfsB4wxlwBTgX5AIHCbMSYQ8Ad2/3JYUQVmFBFxS15TFpZlxQEHz3m4G7DdsqwdlmUVAPOAAUAmpwsDvOhnJCLye7z9jbAJ/50g4HRJNAEWAjcZY6YBS+wIJiLiTrz9BLc5z2OWZVl5wJ0VHUZExF15+2SRCTQ962t/YI9NWURE3Ja3l0US0MYY08IYUw24FfjM5kwiIm7Ha8rCGPMhsAYIMMZkGmPutiyrEBgHLAe2Agssy9piZ04REXekjQRFRMQlr5ksRETkwqksRETEJZWFiIi4pLIQERGXVBYiIuKSykJERFxSWYiIiEstm3VjAAAAEUlEQVQqCxERcUllISIiLv0/nTYtvPXqVsgAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"plt.plot(row, strict_deriv, color=\"red\", label='True deriv')\n",
"plt.plot(row, strict_deriv, color=\"blue\", label='True deriv')\n",
"plt.legend(loc='best')\n",
"plt.show()\n",
"\n",
"plt.loglog(lg, lg_der, color=\"red\", label='True deriv')\n",
"plt.loglog(lg, lg_my_deriv, color=\"blue\", label='True deriv')\n",
"plt.legend(loc='best')\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": 110,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Min, Max and Mean errors for range\n",
"0.0 2.842170943040401e-14 9.907543535581098e-16\n",
"Min, Max and Mean errors for logspace\n",
"0.0 8.881784197001252e-16 2.575717417130363e-16\n"
]
}
],
"source": [
"error = np.abs(strict_deriv-my_deriv)\n",
"print(\"Min, Max and Mean errors for range\")\n",
"print(np.min(error), np.max(error), np.mean(error))\n",
"error = np.abs(lg_der-lg_my_deriv)\n",
"print(\"Min, Max and Mean errors for logspace\")\n",
"print(np.min(error), np.max(error), np.mean(error))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 4. Richardson Check"
]
},
{
"cell_type": "code",
"execution_count": 111,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Maybe it is O(n^1)?\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD8CAYAAABn919SAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAHxJJREFUeJzt3Xl4VfW97/H3N/NASBgSCGEIg4ACIhAQtVUrWtEqarWO4ETFep+2ntZra9v7dD7nnntatdVqLSriVLUOrdjSqhUtOACGeZI5QMKUAAkJGXfyu3/sYBEDGfaw9vB5PU+ePa3s9VnZyYfFb03mnENERKJfgtcBREQkOFToIiIxQoUuIhIjVOgiIjFChS4iEiNU6CIiMUKFLiISI1ToIiIxQoUuIhIjksI5s969e7vCwsJwzlJEJOotW7aswjmX29507Ra6mc0BLgP2O+dGH/fa/wZ+BeQ65yrae6/CwkKKi4vbm0xERI5hZjs6Ml1HhlzmAlPbmMEA4CJgZ6eSiYhISLRb6M65hcDBNl56EPgeoLN7iYhEgC5tFDWzaUCZc25VB6adZWbFZlZcXl7eldmJiEgHdLrQzSwD+BHw445M75yb7Zwrcs4V5ea2O6YvIiJd1JU19KHAYGCVmZUA/YHlZtY3mMFERKRzOr3bonNuDZB39HFrqRd1ZC8XEREJnXbX0M3sBeAjYISZlZrZzNDHEhGRzmp3Dd05d0M7rxcGLY2ISIw50uDjgbc3MWPyIAp7Z4Z0Xjr0X0QkhN5ev48n399OeU1DyOelQhcRCaHXV5ZRkJPOhIE9Qj4vFbqISIgcqGlg4eYKLh/bj4QEC/n8VOgiIiEyf+1emlscV5zRLyzzU6GLiITIvJVlDO/TjZF9s8IyPxW6iEgIlB6q5eOSQ1xxRgFmoR9uARW6iEhIvLFqDwDTxoZnuAVU6CIiIfH6yjImDOrBgJ4ZYZunCl1EJMjW7a7ik73VXBmmjaFHqdBFRILsteVlJCcal4dxuAVU6CIiQeVrbuH1lWVMGdmHnIyUsM5bhS4iEkSLNldQUdPIV8cXhH3eKnQRkSB6ZXkpPTNTOH9EXvsTB5kKXUQkSKrqmnh7/T6mje1HSlL461WFLiISJPPX7KHR1+LJcAuo0EVEgubVZaUMy+vGmIJsT+avQhcRCYIdB45QvOMQV4/vH7ZD/Y+nQhcRCYLXlpdhBleOC+++58dSoYuIBMg5x2srSjlnaG/ys9M9y6FCFxEJ0OJtB9l1sI6rJ3izMfSodgvdzOaY2X4zW3vMc78ys0/MbLWZ/dnMckIbU0Qkcr308U6y0pK4ZHS+pzk6soY+F5h63HNvA6Odc6cDm4AfBDmXiEhUqKptYv7avVx5RgFpyYmeZmm30J1zC4GDxz33lnPO1/pwMdA/BNlERCLeX1aW0ehr4bqJA7yOEpQx9NuBv5/oRTObZWbFZlZcXl4ehNmJiEQG5xwvfryL0QXdGe3RvufHCqjQzexHgA94/kTTOOdmO+eKnHNFubm5gcxORCSirC07zIY9h7muyPu1c4Ckrn6jmd0CXAZMcc654EUSEYkOL368k9SkBKad4e3eLUd1qdDNbCrwfeA851xtcCOJiES+usZm5q3czVfG5JOdnux1HKBjuy2+AHwEjDCzUjObCfwOyALeNrOVZvZYiHOKiESU+Wv2UN3g49oI2Bh6VLtr6M65G9p4+skQZBERiRovfbyLwl4ZnDm4p9dRPqUjRUVEOmlreQ1LSw5y3cSBnp2Iqy0qdBGRTnp+8U6SE83zQ/2Pp0IXEemEusZmXlm2i4tH9SUvK83rOJ+hQhcR6YQ3Vu3mcL2PGZMHeR3lc1ToIiKd8OziHQzv041JEbQx9CgVuohIB63aVcmasiqmTx4UURtDj1Khi4h00LOLd5CRkshV4yJrY+hRKnQRkQ6orG3kjVW7uXJcAVlpkXFk6PFU6CIiHfDKslIafC1MPzPyNoYepUIXEWlHS4vj+SU7mTCoB6f16+51nBNSoYuItOODrRVsrzgSkbsqHkuFLiLSjrkflNC7WwpTR/f1OspJqdBFRE5ie8UR3vlkPzeeOcjza4a2R4UuInIScz/YTnKiMX3yQK+jtEuFLiJyAlV1Tby8rJTLx/aLuPO2tEWFLiJyAi8X76K2sZnbzxnsdZQOUaGLiLTB19zCUx+UMKmwJ6MLsr2O0yEqdBGRNvxzwz7KKuu4/QuFXkfpMBW6iEgb5nxQQkFOOhedFtm7Kh5LhS4icpy1ZVUs3X6QW88uJDEh8s6qeCLtFrqZzTGz/Wa29pjneprZ22a2ufW2R2hjioiEz5z3t5ORksi1Ewd4HaVTOrKGPheYetxz9wHvOOdOAd5pfSwiEvXKKuuYt2o3100cQHZ6ZJ5V8UTaLXTn3ELg4HFPXwE83Xr/aeDKIOcSEfHEnPe344Cvf3GI11E6ratj6H2cc3sAWm/zghdJRMQbVbVNvLB0J9PG9qMgJ93rOJ0W8o2iZjbLzIrNrLi8vDzUsxMR6bJnF5dQ29jMrHOjb+0cul7o+8wsH6D1dv+JJnTOzXbOFTnninJzc7s4OxGR0KpvambuhyWcNzyXU/Mj95znJ9PVQp8H3NJ6/xbg9eDEERHxxqvLS6moaeTO86Jz7Rw6ttviC8BHwAgzKzWzmcB/AxeZ2WbgotbHIiJRqbnF8fjCbYztn81ZQ3p5HafLktqbwDl3wwlemhLkLCIinnhr3V5KDtTy6E3jMYueA4mOpyNFRSSuOef4/b+2MqhXBhePip7D/NuiQheRuPavTeWsLq3iG+cNjarD/NuiQheRuOWc4+EFW+iXncbV4/t7HSdgKnQRiVsfbT3Ash2HuOv8oaQkRX8dRv8SiIh00UMLNpOXlcrXiqLrJFwnokIXkbi0dPtBFm87yJ3nDSUtOdHrOEGhQheRuPTwgs307pbCjZMGeh0laFToIhJ3Vuw8xKLNFdzxxSGkp8TG2jmo0EUkDj28YAs9MpKZPnmQ11GCSoUuInFl5a5KFnyyn5lfGExmarsHy0cVFbqIxJX739pIz8wUbj1nsNdRgk6FLiJxY/G2AyzaXMFd5w2lW4ytnYMKXUTihHOO+9/aSJ/uqcw4K7bGzo9SoYtIXPjXpnI+LjnENy84JWb2Oz+eCl1EYp5/7XwT/Xukc12MHBXaFhW6iMS8N9ftY01ZFXdPOSUmztlyIrG7ZCIi+K9G9MDbGxmSm8lV4wq8jhNSKnQRiWnzVpWxaV8N37lwOEmJsV15sb10IhLX6pua+fWbmxhTkM1XxuR7HSfkVOgiErPmflhCWWUdP7h0JAlRfjWijgio0M3sO2a2zszWmtkLZpYWrGAiIoE4dKSRR97dwgUj8zh7aG+v44RFlwvdzAqAbwNFzrnRQCJwfbCCiYgE4uEFWzjS4OO+S0Z6HSVsAh1ySQLSzSwJyAB2Bx5JRCQwOw4c4dnFJVxbNIDhfbK8jhM2XS5051wZ8GtgJ7AHqHLOvRWsYCIiXfWrNzeSlJDAdy8a7nWUsApkyKUHcAUwGOgHZJrZ9Damm2VmxWZWXF5e3vWkIiIdsHJXJX9dvYc7zh1CXvf42qwXyJDLhcB251y5c64JeA04+/iJnHOznXNFzrmi3NzcAGYnInJyzjl+8df19O6Wwqxzh3gdJ+wCKfSdwGQzyzAzA6YAG4ITS0Sk815fuZtlOw7xvYtHxuTpcdsTyBj6EuAVYDmwpvW9Zgcpl4hIpxxp8PF//76B0/tnc82E/l7H8URA/4Q5534C/CRIWUREuuzR97aw73ADj940IS4OImqLjhQVkai380Atjy/azlfHFTBhUA+v43hGhS4iUe+Xf1tPUoLx/Tg6iKgtKnQRiWqLNpfz1vp9fPOCYfSJs90Uj6dCF5Go1ehr4WdvrGdQrwxmfmGw13E8p0IXkaj1+KJtbNlfw08uP43UpNi8TmhnqNBFJCrtPFDLQ+9s5pLRfblgZB+v40QEFbqIRB3nHD+et5akBOPHl5/mdZyIoUIXkajz97V7eW9jOfd8eQT52elex4kYKnQRiSrV9U387I11jOrXnZvPGuR1nIgSfyc7EJGodv9bm9hf3cAfZhTF/EWfO0s/DRGJGqtLK3nmoxKmnzmIMwbkeB0n4qjQRSQqNPpauPfl1eRlpXHv1BFex4lIGnIRkajwyLtb2LivmidvKaJ7WrLXcSKS1tBFJOJt2HOYR97dwlXjCphyqvY5PxEVuohENF9zC/e+soqcjGR+fJn2OT8ZDbmISESbvWgba8sO8/ubxtMjM8XrOBFNa+giErG27K/hN//czKVj+nLJmHyv40Q8FbqIRKSm5hbu+dNKMlIS+dm00V7HiQoachGRiPTwgi2sKq3i0ZvGk5uV6nWcqKA1dBGJOMt3HuKRd7fw1fEFXKqhlg4LqNDNLMfMXjGzT8xsg5mdFaxgIhKfjjT4+O5LK+nbPY2fThvldZyoEuiQy2+BfzjnrjGzFCAjCJlEJI798m8b2HGwlhfvmKwDiDqpy4VuZt2Bc4FbAZxzjUBjcGKJSDz65/p9vLB0J3eeN4Qzh/TyOk7UCWTIZQhQDjxlZivM7AkzywxSLhGJM/sP13Pfa6s5Nb87371ouNdxolIghZ4EjAd+75wbBxwB7jt+IjObZWbFZlZcXl4ewOxEJFY1tzj+46WV1DT4+O31Z+j6oF0USKGXAqXOuSWtj1/BX/Cf4Zyb7Zwrcs4V5ebmBjA7EYlVj767hQ+3HuDn00YzvE+W13GiVpcL3Tm3F9hlZkfPYzkFWB+UVCISN5ZsO8CD/9zElWf042tF/b2OE9UC3cvlW8DzrXu4bANuCzySiMSLAzUNfPvFFQzqlckvrxqDmXkdKaoFVOjOuZVAUZCyiEgcaWlx3PPyKg4daWLOrRPplqoD1wOlI0VFxBOPLdzKexvL+T+Xncqoftlex4kJKnQRCbtFm8v59Zsbuez0fGZMHuR1nJihQheRsNp1sJZvv7CCU/Ky+J9rTte4eRCp0EUkbOqbmrnr+WX4WhyPzZhARorGzYNJP00RCQvnHD/681rWlh3myVuKGNxbB5YHm9bQRSQsnluyk1eXl3L3lFN0oecQUaGLSMh9uLWCn81bx5dG5HL3lFO8jhOzVOgiElLbymu467nlDO6dyW9vGEdCgjaChooKXURCpqq2ia8/XUxigjHn1ok6v3mIqdBFJCSamlu46/lllB6q4w8zJjCgp65/E2ray0VEgs45x0/mrePDrQe4/2tjmVjY0+tIcUFr6CISdI8v2sYfl+zkrvOHcvUEnUExXFToIhJUf1lRxn/N/4SvjMnn3i+PaP8bJGhU6CISNIs2l3PvK6uYPKQn9187Vnu0hJkKXUSCYm1ZFd94dhlDc7sx++Yi0pJ1GblwU6GLSMB2Hqjl1qc+Jicjhadvn6TdEz2ivVxEJCD7q+u5ec4SfC0tvHj7mfTpnuZ1pLilNXQR6bKDRxqZ/sQS9lc38OQtExmWpws8e0mFLiJdUlXXxM1zlrDjQC1P3FLEhEE9vI4U91ToItJpRxp83PbUUjbureaxGRM4e2hvryMJQSh0M0s0sxVm9tdgBBKRyFbf1MzXny5mVWkVD98wji+NyPM6krQKxhr63cCGILyPiES4+qZmZj27jMXb/Yf0Tx2d73UkOUZAhW5m/YGvAE8EJ46IRKraRh+3z/2YRZvL+X9Xn86V4wq8jiTHCXS3xd8A3wO0aVskhtU0+Mu8uOQgD1w7lqvG6fwskajLa+hmdhmw3zm3rJ3pZplZsZkVl5eXd3V2IuKR6vombpmzlGU7DvGb68epzCNYIEMu5wDTzKwEeBG4wMyeO34i59xs51yRc64oNzc3gNmJSLhV1TYx48mlrNpVycM3jGPa2H5eR5KT6HKhO+d+4Jzr75wrBK4HFjjnpgctmYh4at/heq79w0es213FozeN59Ix2gAa6XTov4h8zrbyGmY8uZTK2kbm3jaJc4ZpP/NoEJRCd869B7wXjPcSEW+tKa3i1qeWAvDirLMY0z/b40TSUVpDF5FPfbilgjueKSYnI4VnZ05iSG43ryNJJ6jQRQSA15aX8v1XVzOkdzeemTlJZ02MQip0kTjnnOPBtzfx0IItnDWkF49Nn0B2hs5nHo1U6CJxrL6pmXtfWc0bq3ZzbVF/fnnlGFKSdM6+aKVCF4lTB2oauOOZYpbvrOT7U0fyjfOGYKZrgEYzFbpIHFq3u4o7n11GeXUDv79pPJdoH/OYoEIXiTN/WVHGfa+tJic9hT/deRZjB+R4HUmCRIUuEieamlv4r/kbeOqDEiYN7skjN44nNyvV61gSRCp0kThQXt3AN/+4nCXbD3LbOYX88NJTSU7Uxs9Yo0IXiXGLtx3g7hdXUFXXxIPX6dS3sUyFLhKjmlscDy/YzEPvbKawVyZzbp3IqH46jD+WqdBFYtC+w/Xc/eIKFm87yFfHFfCLK0eTmao/91inT1gkxrz7yX7ueXkVdY3N/PprY7lmgoZY4oUKXSRGHGnw8Z/zN/DHJTsZ2TeL3904nmF5OrlWPFGhi8SAj0sOcs+fVrHrUC13njuE71w0nLTkRK9jSZip0EWiWH1TMw+8vYnHF21jQI8MXpp1FpMG9/Q6lnhEhS4SpZbtOMR9r65m8/4abjpzID+89FRt+Ixz+vRFoszh+iZ+9Y+NPLdkB/nd05h720TOH5HndSyJACp0kSjyj7V7+cm8tZRXN3Db2YO558vDtVYun9JvgkgUKD1Uy8/fWM9b6/dxan53Zs8o0km15HO6XOhmNgB4BugLtACznXO/DVYwEfFv9HzsX1v5/XtbMYP7LhnJzC8M1nlYpE2BrKH7gHucc8vNLAtYZmZvO+fWBymbSNxyzvH3tXv5z79toKyyjq+cns8PLz2Vgpx0r6NJBOtyoTvn9gB7Wu9Xm9kGoABQoYsEYMOew/z8jfV8tO0AI/tm8eKsyUwe0svrWBIFgjKGbmaFwDhgSTDeTyQe7TpYy4Nvb+LPK8vITk/mF1eO5oaJA0jS8Ip0UMCFbmbdgFeB/3DOHW7j9VnALICBAwcGOjuRmHPwSCO/W7CF5xbvwAzuPHcod503lOyMZK+jSZQJqNDNLBl/mT/vnHutrWmcc7OB2QBFRUUukPmJxJKaBh9Pvb+dPyzcRm2jj2uLBnD3haeQn61xcumaQPZyMeBJYINz7oHgRRKJbVV1TTz9YQlPvr+dqromLh7Vh3svHsGwvCyvo0mUC2QN/RxgBrDGzFa2PvdD59z8wGOJxJ6q2iae/GA7T32wnep6Hxee2odvTxnG6f21P7kERyB7ubwPWBCziMSk/YfrmfthCc98tIOaBh9TR/XlmxcMY3SBrh4kwaUjRUVC5JO9h3l84XbmrSrD1+K4dHQ+35oyjJF9u3sdTWKUCl0kiJxzLNxcwROLtrFocwXpyYncOGkgt39hMIN6ZXodT2KcCl0kCA7XN/HaslKeX7KTzftryM1K5d6LR3DTmQPJyUjxOp7ECRW6SADWlFbx3OIdzFu1m7qmZsb2z+ZX15zOtDP6kZqkKwZJeKnQRTrpcH0T81fv4YWlO1lVWkVacgJXjC1g+uRBjOmvDZ3iHRW6SAc0tzje31LBq8tKeXPdXhp8LQzL68ZPLz+Nq8b3JztdR3WK91ToIiexcW81ry0v5c8rythf3UB2ejLXFg3gq+MLOGNADv7j60Qigwpd5Dib9lXz19V7mL9mD1v215CUYJw/Io+rxxdwwal5GhuXiKVCl7jnnGPTvhr+tubfJW4Gkwp7cvMVo7h0TD69u6V6HVOkXSp0iUuNvhY+LjnIOxv2s+CTfZQcqMUMzhzck1vOGsXFo/uSl5XmdUyRTlGhS9yoqGngvY3lLPhkHws3VVDT4CMlKYGzh/Zi5heHcPGoPipxiWoqdIlZdY3NLC05yIdbKvhgawXrdh/GOcjLSuXysflcMLIP5wzrRUaK/gwkNug3WWJGU3MLq0ur+HBLBe9vqWDFzkoam1tITjTGDezBdy4czgUj8zgtvzsJCdo7RWKPCl2iVlVtE8t3HqJ4x0GKSw6xqrSS+qYWzOC0/O7cek4h5wzrzcTCHloLl7ig33KJCr7mFraWH2F1aSUrdlWyrOQQG/dVA5CYYIzq150bJg1kYmFPJg/pRc9MnT9F4o8KXSJOc4tja3kNq0urWFtWxerSStbvOUx9UwsAWWlJjB/Yg8tOz2dCYQ/OGJCjNXARVOjiIecc5dUNbNxXzca91WzaV83GfTVs2ltNXVMzABkpiYzul82NkwYxpn93xhTkMLh3JokaAxf5HBW6hJyvuYWyyjq2VxyhpOIIW8uPsHGfv8Ara5s+na53txSG98ni+kkDOL1/NmMKshncu5vKW6SDVOgSFA2+ZvZU1lN6qI7tB/zFXVJxhO0HjrDrYC1Nze7TabNSkxjeN4tLRuczok83hvfNYnifLB2NKRIgFbq0yznHodom9lbVU1ZZx+7KOsoq6yg71HpbWUd5dcNnvictOYHCXpkMz8viy6f1ZXDvDAp7ZTK4dya5Wak6qZVICARU6GY2FfgtkAg84Zz776CkkpBraXFU1/s4VNvIgSONlFc3UF7TQPnhev9tdQP7q/23FTUNn1nDBkhNSqAgJ51+OelcMCKPfjnp9MtJo6BHOoN7Z9InK037eouEWZcL3cwSgUeAi4BS4GMzm+ecWx+scHJyDb5maup9HGloprqhyX+/0Ud1vf+rsraRQ7VNVNY2td5v9N+v8z9ucZ9/TzPolZlKXlYquVmpDO+T9en9Pt3TKMhJp6BHOr0yU7SWLRJhAllDnwRscc5tAzCzF4ErgLgo9OYWR1NzC74Wh+/T238/19zSQlNz63MtLfiaHb7W2wZfC/VNzf4vXwsNR+83tT7v89//zHRNzVS3FnZNvY+aBt/n1prbkp6cSI+MZHIyUuiRmUx+djo5Gcn0yEj59LZnZgq5Wf4S75mZQlJiQhh+giISbIEUegGw65jHpcCZgcVp20PvbOb1lWU4AAcO/7iu/xYczn/b2m9tvsbR1499fMx0rfdP+v6tz/ta3KfzCqYEg7TkRP9XUgJpyYmkJieSlpxAenIiA3pmkJWaRGZqEt3SkuiW6v/KTP33ff/ziXRLTSYnI5m0ZJ27WyReBFLobf1/+3M1Z2azgFkAAwcO7NKM+nRPZWTf7mD+mZpZ6+1nH/tft2OeP+Zx6wRtvtb6HrT5/GffP8GM5AQjMSGBpEQjOdFISkjw3yYmkJhw3HOt0yUdM31qkr+k/bf++2nJiSQlmIYxRKTLAin0UmDAMY/7A7uPn8g5NxuYDVBUVNSl9drrJg7kuold+8dARCReBDJY+jFwipkNNrMU4HpgXnBiiYhIZ3V5Dd055zOzbwJv4t9tcY5zbl3QkomISKcEtB+6c24+MD9IWUREJADaP01EJEao0EVEYoQKXUQkRqjQRURihApdRCRGmAvFMewnmplZObCji9/eG6gIYpxooGWOD1rm+BDIMg9yzuW2N1FYCz0QZlbsnCvyOkc4aZnjg5Y5PoRjmTXkIiISI1ToIiIxIpoKfbbXATygZY4PWub4EPJljpoxdBERObloWkMXEZGTiLhCN7OpZrbRzLaY2X1tvJ5qZi+1vr7EzArDnzK4OrDM3zWz9Wa22szeMbNBXuQMpvaW+ZjprjEzZ2ZRvUdER5bXzK5t/ZzXmdkfw50x2Drwez3QzN41sxWtv9uXepEzmMxsjpntN7O1J3jdzOyh1p/JajMbH9QAzrmI+cJ/Gt6twBAgBVgFnHbcNP8LeKz1/vXAS17nDsMyfwnIaL1/Vzwsc+t0WcBCYDFQ5HXuEH/GpwArgB6tj/O8zh2GZZ4N3NV6/zSgxOvcQVjuc4HxwNoTvH4p8Hf8F0SbDCwJ5vwjbQ390wtPO+cagaMXnj7WFcDTrfdfAaZYdF+3rd1lds6965yrbX24GP/VoaJZRz5ngF8A/wPUhzNcCHRkee8AHnHOHQJwzu0Pc8Zg68gyO6B76/1s2rjiWbRxzi0EDp5kkiuAZ5zfYiDHzPKDNf9IK/S2LjxdcKJpnHM+oAroFZZ0odGRZT7WTPz/wkezdpfZzMYBA5xzfw1nsBDpyGc8HBhuZh+Y2WIzmxq2dKHRkWX+KTDdzErxX1fhW+GJ5qnO/r13SkAXuAiBjlx4ukMXp44iHV4eM5sOFAHnhTRR6J10mc0sAXgQuDVcgUKsI59xEv5hl/Px/w9skZmNds5VhjhbqHRkmW8A5jrn7jezs4BnW5e5JfTxPBPS/oq0NfSOXHj602nMLAn/f9VO9l+cSNehi22b2YXAj4BpzrmGMGULlfaWOQsYDbxnZiX4xxrnRfGG0Y7+Xr/unGtyzm0HNuIv+GjVkWWeCfwJwDn3EZCG/3wnsaxDf+9dFWmF3pELT88Dbmm9fw2wwLVubYhS7S5z6/DDH/CXebSPrUI7y+ycq3LO9XbOFTrnCvFvN5jmnCv2Jm7AOvJ7/Rf8G78xs974h2C2hTVlcHVkmXcCUwDM7FT8hV4e1pThNw+4uXVvl8lAlXNuT9De3eutwifYCrwJ/xbyH7U+93P8f9Dg/9BfBrYAS4EhXmcOwzL/E9gHrGz9mud15lAv83HTvkcU7+XSwc/YgAeA9cAa4HqvM4dhmU8DPsC/B8xK4MteZw7CMr8A7AGa8K+NzwS+AXzjmM/5kdafyZpg/17rSFERkRgRaUMuIiLSRSp0EZEYoUIXEYkRKnQRkRihQhcRiREqdBGRGKFCFxGJESp0EZEY8f8BcH5F9heBJcUAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Maybe it is O(n^2)?\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD8CAYAAABn919SAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAH/pJREFUeJzt3Xl81dWd//HXJwtJCElISNgTwr4IihAW9xWluHX664LW1gVltFPbTrW2tk7tdNqZLr9xabWjuOFSsdVapeMKblgVZN8h7CQhJIGQhZCE5N4zfyRaSoFccu/N9y7v5+ORx92+yfdzuDdvTs73fL/HnHOIiEj0S/C6ABERCQ0FuohIjFCgi4jECAW6iEiMUKCLiMQIBbqISIxQoIuIxAgFuohIjFCgi4jEiKSu3Flubq4rLCzsyl2KiES95cuX73PO5XW0XYeBbmZPAJcDlc65sUe9dgfwayDPObevo59VWFjIsmXLOtpMRESOYGa7AtkukCGXucD0Y+wgH5gG7D6pykREJCw6DHTn3CKg+hgv3QfcCejqXiIiEaBTB0XN7EqgzDm3OoBtZ5vZMjNbVlVV1ZndiYhIAE460M2sO/Aj4MeBbO+cm+OcK3LOFeXldTimLyIindSZHvpQYDCw2sx2AgOBFWbWN5SFiYjIyTnpaYvOubVA708ft4d6USCzXEREJHw67KGb2TzgY2CkmZWa2azwlyUiIierwx66c+7qDl4vDFk1IiIxpvGwj1+9uYkbzxpMfk73sO5Lp/6LiITRiytKefLDnZTXNoV9Xwp0EZEw8fsdT/x1B6cNzGJSYXbY96dAFxEJk4UbK9ixr4GbzhmCmYV9fwp0EZEwefSD7QzomcbnxnbNrG4FuohIGKzYfYClOw9ww1mFJCV2TdQq0EVEwuDh97aRlZbM1ZMLumyfCnQRkRDbWlnPWxsq+PoZg0hP6bplJxToIiIh9sj720lNTuD6Mwu7dL8KdBGRECqvbeTlVWV8uSifXj1SunTfCnQRkRB64q878Du4+ZwhXb5vBbqISIjUHmrhuSW7ufzUfmE/zf9YFOgiIiHyzOKdNBz28c/nDvVk/wp0EZEQaGrx8eSHOzlvRB5j+md6UoMCXUQkBF5YVsL+hsPccp43vXNQoIuIBK3F5+eRRds5Lb8nU4fkeFaHAl1EJEh/XllG6YFGvnXhsC65CNfxKNBFRILQ6vPzu3e3ckr/TC4c1bvjbwgjBbqISBD+smYPO/cf4rYLh3vaOwcFuohIp/n8jgff2cqovhlcMqaP1+Uo0EVEOuu1teVsq2rgtguHk5Dgbe8cAgh0M3vCzCrNbN0Rz/3azDaZ2Roz+7OZ9QxvmSIikcXf3jsf1rtHly1g0ZFAeuhzgelHPbcAGOucOxUoBu4KcV0iIhHtrQ172VxRz20XDouI3jkEEOjOuUVA9VHPveWca21/uBgYGIbaREQiknOO37y9lSG56Vx+an+vy/lMKMbQbwReP96LZjbbzJaZ2bKqqqoQ7E5ExFtvb6xkQ3kd37hgGIkR0juHIAPdzH4EtAK/P942zrk5zrki51xRXl5eMLsTEfGc3++4b2ExBTnduWp85PTOIYhAN7PrgMuBrzrnXOhKEhGJXG+u38v6PXV85+LhJHfR4s+B6tRid2Y2Hfg+cJ5z7lBoSxIRiUw+v+PeBcUMzUvnqvEDvC7nHwQybXEe8DEw0sxKzWwW8CCQASwws1Vm9nCY6xQR8dxfVu9hS+VBvjttZESNnX+qwx66c+7qYzz9eBhqERGJWC0+P/cvLGZ0v8yImXd+tMgaABIRiVAvrShl5/5D3D5tRMTMOz+aAl1EpAPNrT5+8/ZWTsvvyUWjvb2i4oko0EVEOvCHpSWU1TTyvUtGen5FxRNRoIuInEDjYR+/fWcrUwbncNawXl6Xc0IKdBGRE5j70U6q6pu5PcJ756BAFxE5rgMNh/nde1u5cFRvJg/2bq3QQCnQRUSO46F3t9LQ3Mr3p4/yupSAKNBFRI6hpPoQT3+8iy9OHMjIvhlelxMQBbqIyDHcu6AYM/jXaSO8LiVgCnQRkaOsK6vl5VVl3Hj2YPplpXldTsAU6CIiR/nlG5vISkvmlvOGel3KSVGgi4gc4YMtVXywZR/fvGAYWWnJXpdzUhToIiLt/H7HL17fxMDsNL52xiCvyzlpCnQRkXZ/XlnG+j113HHJSFKSEr0u56Qp0EVEgIbmVn75xiZOy+/JladF1tJygVKgi4gAD7+/jcr6Zu65YkzEXh63Iwp0EYl7pQcOMWfRdq4a358JBdlel9NpCnQRiXv/9fomzIiaU/yPR4EuInFt6c5qXl1Tzj+fO5T+PaPnJKJjUaCLSNzy+x0//csG+mWlRt1JRMfSYaCb2RNmVmlm6454LsfMFpjZlvbb6B10EpG49acVpawtq+UHnxtFWrfom6Z4tEB66HOB6Uc99wPgbefccODt9sciIlHjYHMrv3pzM6cXRO80xaN1GOjOuUVA9VFPXwU81X7/KeDzIa5LRCSs7l9QzL6DzdxzxSkRvxJRoDo7ht7HOVcO0H4buctgi4gcZfPeep78aCczJxUwPr+n1+WETNgPiprZbDNbZmbLqqqqwr07EZETcs7xb6+sIzM1iTsvHel1OSHV2UCvMLN+AO23lcfb0Dk3xzlX5JwrysvL6+TuRERC4+VVZXyyo5rvTx9Fdno3r8sJqc4G+nzguvb71wGvhKYcEZHwqWtq4eevtl2v5ctF+V6XE3KBTFucB3wMjDSzUjObBfwCmGZmW4Bp7Y9FRCLafQuK2d/QzM+uGhu112s5kaSONnDOXX2cly4KcS0iImGzYU8dT320k2unDGLcwCyvywkLnSkqIjHP73f8+JV19OzejTsuia0DoUdSoItIzJu3dDfLdh3grs+NIqt7dC0rdzIU6CIS0yrqmvjFa5s4c2gvvjhxoNflhJUCXURi2k/mr+ewz89//tO4mDkj9HgU6CISs95av5fX1+3l2xcPpzA33etywk6BLiIxqb6phR+/sp5RfTO4+ZwhXpfTJTqctigiEo1+9cZmKuqbePhrE0lOjI++a3y0UkTiyvJd1Ty7ZBfXn1kYUxff6ogCXURiSnOrjx/8aS39MlO5PYbnnB+LhlxEJKbcv3ALWyoP8uQNk+iREl8Rpx66iMSMVSU1PPL+Nr5cNJALRsbfMg0KdBGJCU0tPm7/4yr6ZKZy9+VjvC7HE/H194iIxKz7FhazraqBp26cTGZq7J7efyLqoYtI1Fux+wCPLtrO1ZPzOW9E/C6ko0AXkajW1OLjjhdW0y8rjR/OGO11OZ7SkIuIRLV7FxSzvaqBZ2dNISNOh1o+pR66iEStj7ft59EPtnPNlALOHp7rdTmeU6CLSFSqPdTC7X9cRWGvdO6+LL6HWj6lIRcRiTrOOe5+ZR2V9c386dYz6d5NUQbqoYtIFHpl1R7+snoP37l4OKfF0bVaOhJUoJvZv5rZejNbZ2bzzCw1VIWJiBxLSfUh/u3ldUwqzObW84d5XU5E6XSgm9kA4FtAkXNuLJAIzAxVYSIiR/P5Hbf/cTUA9355PIkJsb0C0ckKduApCUgzsxagO7An+JJERI7t4fe38cnOau77ymnk53T3upyI0+keunOuDPj/wG6gHKh1zr0VqsJERI60fNcB7l1QzBWn9efz4wd4XU5ECmbIJRu4ChgM9AfSzezaY2w328yWmdmyqqqqzlcqInGr5tBhvjVvJQN6pvHzfxob84s9d1YwB0UvBnY456qccy3AS8CZR2/knJvjnCtyzhXl5cXvNRZEpHOcc9zxwhoq65t48JrT4/bCW4EIJtB3A1PNrLu1/Xd5EbAxNGWJiLR5/K87WLixgh/OGM2pAzVF8USCGUNfArwIrADWtv+sOSGqS0SEVSU1/PKNTVwypg/Xn1nodTkRL6hZLs65e4B7QlSLiMhnahtb+OZzK+idkcqvv3iaxs0DoPNlRSTiOOe488XV7K1t4oVbziCru8bNA6FT/0Uk4jz6wXbeXF/B96eP4vSCbK/LiRoKdBGJKB9t3ccvXt/EjHF9uemcwV6XE1UU6CISMfbUNPLNeSsZkteDX2nc/KQp0EUkIjS1+Lj12eUcbvXz8LUT6ZGiQ3wnS/9iIhIR/v0v61ldWsvD105kWO8eXpcTldRDFxHPPf/JbuZ9UsI3zh/K9LF9vS4nainQRcRTK3cf4MevrOec4bncfslIr8uJagp0EfHMnppGbn56OX2zUnlg5um6vnmQNIYuIp44dLiVm55aRlOLj+dunkJOejevS4p66qGLSJfz+x3f/cNqNu2t47dXn86IPhlelxQTFOgi0uXuX1jMG+v38sMZo7lgVG+vy4kZCnQR6VLzV+/hN+9s5UsTBzLrbJ0JGkoKdBHpMqtKavjeC6uZVJjNz7TyUMgp0EWkS+za38CsuUvpnZnCw9dOJCUp0euSYo4CXUTCrrrhMNc/uRSfc8y9YTK9eqR4XVJM0rRFEQmrphYfNz21lLKaRp67aQpD83Raf7iohy4iYePzO77z/CpWltTwwFfGU1SY43VJMU2BLiJh87NXN/DG+r3cfdkYPjeun9flxDwFuoiExWMfbOfJD3dy41mDNT2xiwQV6GbW08xeNLNNZrbRzM4IVWEiEr1eXF7Kz17dyIxxfbn7stFelxM3gj0o+gDwhnPui2bWDegegppEJIq9uX4v3//TGs4elst9XxlPgi641WU6HehmlgmcC1wP4Jw7DBwOTVkiEo0+2rqP255bybgBWTzyNc0172rBDLkMAaqAJ81spZk9ZmbpIapLRKLM6pIabn56GYNz05l7wyTStYRclwsm0JOACcD/OOdOBxqAHxy9kZnNNrNlZrasqqoqiN2JSKTaWlnP9U9+Qk6Pbjw9azI9u+tSuF4IJtBLgVLn3JL2xy/SFvB/xzk3xzlX5JwrysvLC2J3IhKJdu1v4NrHPiEpMYFnZ02hT2aq1yXFrU4HunNuL1BiZp+uGXURsCEkVYlIVCipPsTVcxbT1Orj6RsnM6iXRl29FOwg123A79tnuGwHbgi+JBGJBqUHDnH1o4tpOOzj9zdNYXS/TK9LintBBbpzbhVQFKJaRCRK7Klp5JpHl1Db2MJzN01l7IAsr0sSdKaoiJykvbVNXPPoYg40HOaZWVMYN1BhHikU6CISsMq6tjCvqm9m7o2TGZ/f0+uS5AiaKCoiAdlT08hXH1tCRV0TT904mYmDsr0uSY6iQBeRDu3a38A1jy6hrrGFZ2ZNZuIgXQY3EinQReSEtlYe5KuPLaa51c9zN0/VmHkEU6CLyHFtLK/j2seWYGY8P3sqo/pqamIkU6CLyDGtLqnh6098QlpyIr+/WUvHRQMFuoj8g4+27mP2M8vp2T2ZeTdPJT9HV8aOBpq2KCJ/59U15Vz/5FL690zlhVvOUJhHEfXQReQzT3+8k3vmr2diQTaPXzeJrO7JXpckJ0GBLiI457h3QTG/fWcrF4/uw4PXnE5qshaniDYKdJE41+rzc/fL63h+aQkzJ+Xzs8+PJSlRo7HRSIEuEscamlv51ryVvL2pktsuHMZ3p43ATGuARisFukicKq9tZNbcZWyuqOdnnx/LtVMHeV2SBEmBLhKH1pbWMuuppRw67OOJ6ydx3gitJhYLFOgiceat9Xv59vOryEnvxp9uncLIvhlelyQhokAXiRPOOR7/6w5+/tpGTh3Yk0e/PpHeGVr/M5Yo0EXiQFOLj397eR0vLC9lxri+/PeXxpPWTdMSY40CXSTGldc2csszy1ldWsu3LhrOdy4aTkKCZrLEIgW6SAxburOaW59dTuNhH498bSKXntLX65IkjIIOdDNLBJYBZc65y4MvSUSC5Zzj2SW7+ff568nP6c68m6cyvI8Ofsa6UPTQvw1sBHShZJEI0NTi4yfz1/P80hIuGJnH/TNPJytN12SJB0Gd32tmA4HLgMdCU46IBGPnvga+8LuPeH5pCf9ywVAeu26SwjyOBNtDvx+4E9DfciIee21tOXe+uIbEBOPx64q4aHQfr0uSLtbpQDezy4FK59xyMzv/BNvNBmYDFBQUdHZ3InIcza0+/vPVjTz18S5OL+jJg9dMYEDPNK/LEg8E00M/C7jSzGYAqUCmmT3rnLv2yI2cc3OAOQBFRUUuiP2JyFFKqg/xL8+tYE1pLTedPZg7p4+iW5KulBivOh3ozrm7gLsA2nvodxwd5iISPq+sKuPul9cBaEqiAJqHLhJ16ppa+PHL63h51R4mDsrm/q+M1zJxAoQo0J1z7wHvheJnicjxfbKjmn/9wyr21jXx3Wkj+Mb5Q7UYhXxGPXSRKNDi83P/wmL+571t5Od054VbzmBCQbbXZUmEUaCLRLjiinrueGE1a0pr+dLEgdxz5Sn0SNGvrvwjfSpEIlSrz88ji7bzwMIt9EhN4ndfncCMcf28LksimAJdJAJt3lvP915s65XPGNeXn141ltweKV6XJRFOgS4SQY7slWekJvHQNRO47FT1yiUwCnSRCLGurJa7XlrL2rJaLju1Hz+98hR6qVcuJ0GBLuKxg82t3PtWMXM/2kFOeorGyqXTFOgiHnpz/V5+Mn89e+ua+OqUAr536ShdHVE6TYEu4oGymkbueWU9CzdWMKpvBg99dYLmlUvQFOgiXai51cfjf93Bg+9sxTn44YxR3HDWYJJ1tqeEgAJdpAs451iwoYKfv7aRXfsPMW1MH+65YgwDs3UNFgkdBbpImG2pqOen/7uBD7bsY1jvHjwzazLnDM/zuiyJQQp0kTCpPdTCfQuLeWbxLtK7JXLPFWO4duogDa9I2CjQRUKsudXHMx/v4sF3t1LX2MLVkwv47rQRmlMuYadAFwkRv9/x8qoy/vutYspqGjlneC53fW40Y/pnel2axAkFukiQnHO8X1zFL9/YzMbyOsYOyOSX/+9Uzh6e63VpEmcU6CJBWLn7AL9+czMfbdtPfk4aD8wczxWn9ichwbwuTeKQAl2kE1aV1HD/wmLe21xFTno37rliDNdMKSAlKdHr0iSOKdBFTsLa0lruW1jMO5sq6dk9mTunj+S6MwpJ14ITEgH0KRQJwLqyWu5fWMzCjZVkpSXzvUtHct2ZhVo5SCJKpz+NZpYPPA30BfzAHOfcA6EqTMRrzjmW7Kjm4fe38d7mKjJTk7h92giuP6uQjFRdQEsiTzDdi1bgdufcCjPLAJab2QLn3IYQ1SbiCb/fsWBjBQ+/v42Vu2vold6NOy4ZwdfPLCRTQS4RrNOB7pwrB8rb79eb2UZgAKBAl6h0uNXPK6vKePj9bWyraiA/J43/+PxYvjRxIKnJOtgpkS8kA4BmVgicDiwJxc8T6UoHGg4zb+lunvl4F+W1TYzul8kDM8dz2bh+JOk0fYkiQQe6mfUA/gR8xzlXd4zXZwOzAQoKCoLdnUjIbNpbx9wPd/LnlWU0t/o5c2gv/usL4zhvRB5mmkcu0SeoQDezZNrC/PfOuZeOtY1zbg4wB6CoqMgFsz+RYPn8jrc3VjD3o518tG0/qckJfGHCQK4/s5CRfTO8Lk8kKMHMcjHgcWCjc+7e0JUkEnpV9c28uLyU5z7ZRUl1I/2zUvn+9FHMnJRPdno3r8sTCYlgeuhnAV8D1prZqvbnfuicey34skSC5/c7Pty2j3mf7Oat9RW0+h2TB+fwg+mjufSUPhofl5gTzCyXvwIaaJSIU1nfxIvLS3n+kxJ2Vx8iu3sy159ZyMzJBQzr3cPr8kTCRqe5SUxobvXx7qZKXlpRxjubKmn1O6YOyeH2S0Zw6Sl9Ne1Q4oICXaKWc46VJTW8tKKUv6wup7axhbyMFG44q603PjRPvXGJLwp0iTol1Yd4eWUZL60sY8e+BlKSErj0lL58YcIAzh6Wq7FxiVsKdIkKe2oaeW1tOa+uLWfl7hoApgzO4dbzhvK5cX11bRURFOgSwfbWNn0W4st3HQBgTL9MvnfpSK48rT/5Od09rlAksijQJaKUVB9iwYYKXl9XztKdbSE+qm8Gd1wyghnj+jFE4+Iix6VAF0/5/Y61ZbUs3FjBgg0VbNpbD8DIPhncPm0EM07tp4ObIgFSoEuXa2rx8fG2/SzYWMHbGyuoqGsmwWBSYQ53Xzaai0f3oTA33esyRaKOAl3CzjnHtqoGFhVXsWhLFYu376epxU96t0TOG5nHxaP7cMHI3joFXyRICnQJi9rGFj7auo9FW6pYVLyPsppGAIbkpjNzUgHnj8zjjKG9tKiySAgp0CUkDh1uZfmuAyzevp/F26tZVVKDz+/ISEnizGG9+MYFQzl3eJ5mpoiEkQJdOqWhuS3Al+xoC/DVJTW0+h2JCca4AVl84/yhnDsij/H5PUnWiT4iXUKBLgGprG9ixa4aVu4+wNKd1awprf0swE8dmMXN5w5h6pBeFA3KJj1FHysRL+g3T/5Bi8/PxvI6Vuw6wIrdNazYfYDSA21j4MmJxtgBWcxuD/CJCnCRiKHfxDjn8zt27DvI2rJa1pXVsba0ljVlNTS1+AHok5nChIJsrjujkAmDsjmlf6auXCgSoRTocaTV52dr1UHWldWxrqyWdWW1bCiv49BhHwApSQmM7pfJ1ZMLmFCQzYRB2fTPStX6miJRQoEeg5xzlNc2sbminuK99W23FfVsqThIc2tbzzstOZFT+mfy5aJ8xg7IYuyATIbl9dCVCkWimAI9ivn9jvK6JnZUNbCt6uDfBXh9U+tn2/XJTGFEnwy+NnUQpwzIZNyALAbn9iAxQT1vkViiQI9wzjn2Nxxmx76Gv31Vtd3u3N/wWY8bICstmZF9M/j8+AGM6JvByD4ZjOjTg57ddQamSDxQoHvMOUdtYwulBxopPXCo/bbxs8dlBxqpb/5bbzs50cjP6c6Q3HTOHZHL4NweDM5NZ0heOr0zUjTeLRLHggp0M5sOPAAkAo85534RkqpihM/v2H+wmYq6Zirqmqiob6KirpnKuiYq6poor22i9EAjB48IbID0bokMzO7OwOw0Jg/OobBXOoPz0hmSm86Anmka5xaRY+p0oJtZIvAQMA0oBZaa2Xzn3IZQFRdp/H5HXVML1Q2HOXDoMNUNLVQ3NFPd0NL+uO1r38G2AK+qb8bv/v5nmEGv9BT6ZKYwMDuNqUN6MTA7rf2rLcSz0pLV0xaRkxZMD30ysNU5tx3AzJ4HrgIiItCdczS3+mlu8dPU6qOpxUdTi5+mFh+NLX973Nzqo76plYPNrdQ3tXCwqZX6plbqmlo52NxCffvjg82t1Da24Ds6odulJCXQK70b2endyEnvxsg+GfTJTKVPZgq9M1M/u5/bI0WnwotIWAQT6AOAkiMelwJTgivn2B56dyu/fnMzg3PT8fnd376cw99+6/O13/odfudo8R07eE8kMcHISE2iR0oSGanJZKQk0TczleG9k+iRmkTPtE8DO5mc9BRyuncjOz2ZnPRupCUnqlctIp4KJtCPlV7/kKJmNhuYDVBQUNCpHe072AzA2AFZJBokJiSQmNAWwAlmJCUYCQlGolnbcwlGcoKRkpxIanIiqckJpCYlktbtb/dTPn0+OZGMlLbAViiLSDQz506+JwtgZmcAP3HOXdr++C4A59x/He97ioqK3LJlyzq1PxGReGVmy51zRR1tF8xg7lJguJkNNrNuwExgfhA/T0REgtDpIRfnXKuZfRN4k7Zpi08459aHrDIRETkpQc1Dd869BrwWolpERCQImj8nIhIjFOgiIjFCgS4iEiMU6CIiMUKBLiISIzp9YlGndmZWBezq5LfnAvtCWE40UJvjg9ocH4Jp8yDnXF5HG3VpoAfDzJYFcqZULFGb44PaHB+6os0achERiREKdBGRGBFNgT7H6wI8oDbHB7U5PoS9zVEzhi4iIicWTT10ERE5gYgLdDObbmabzWyrmf3gGK+nmNkf2l9fYmaFXV9laAXQ5u+a2QYzW2Nmb5vZIC/qDKWO2nzEdl80M2dmUT0jIpD2mtmX29/n9Wb2XFfXGGoBfK4LzOxdM1vZ/tme4UWdoWRmT5hZpZmtO87rZma/af83WWNmE0JagHMuYr5ouwzvNmAI0A1YDYw5aptvAA+3358J/MHrurugzRcA3dvv3xoPbW7fLgNYBCwGiryuO8zv8XBgJZDd/ri313V3QZvnALe23x8D7PS67hC0+1xgArDuOK/PAF6nbcW3qcCSUO4/0nrony087Zw7DHy68PSRrgKear//InCRRfe6cR222Tn3rnPuUPvDxcDALq4x1AJ5nwH+A/gV0NSVxYVBIO29GXjIOXcAwDlX2cU1hlogbXZAZvv9LGBPF9YXFs65RUD1CTa5CnjatVkM9DSzfqHaf6QF+rEWnh5wvG2cc61ALdCrS6oLj0DafKRZtP0PH806bLOZnQ7kO+f+tysLC5NA3uMRwAgz+9DMFpvZ9C6rLjwCafNPgGvNrJS2dRVu65rSPHWyv+8nJagFLsIgkIWnA1qcOooE3B4zuxYoAs4La0Xhd8I2m1kCcB9wfVcVFGaBvMdJtA27nE/bX2AfmNlY51xNmGsLl0DafDUw1zn33+1rFD/T3mZ/+MvzTFjzK9J66KVA/hGPB/KPf4Z9to2ZJdH2p9qJ/sSJdIG0GTO7GPgRcKVzrrmLaguXjtqcAYwF3jOznbSNNc6P4gOjgX6uX3HOtTjndgCbaQv4aBVIm2cBfwRwzn0MpNJ2vZNYFtDve2dFWqAHsvD0fOC69vtfBN5x7UcbolSHbW4ffniEtjCP9rFV6KDNzrla51yuc67QOVdI23GDK51zy7wpN2iBfK5fpu3gN2aWS9sQzPYurTK0AmnzbuAiADMbTVugV3VplV1vPvD19tkuU4Fa51x5yH6610eFj3MUuJi2I+Q/an/up7T9QkPbm/4CsBX4BBjidc1d0OaFQAWwqv1rvtc1h7vNR237HlE8yyXA99iAe4ENwFpgptc1d0GbxwAf0jYDZhVwidc1h6DN84ByoIW23vgs4BbgliPe54fa/03WhvpzrTNFRURiRKQNuYiISCcp0EVEYoQCXUQkRijQRURihAJdRCRGKNBFRGKEAl1EJEYo0EVEYsT/AedTneQpC6XYAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Maybe it is O(n^3)?\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD8CAYAAAB5Pm/hAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAAErRJREFUeJzt3W2MXFd9x/Hvf2bXMQ8pCfEmiuwEU2FUEBJJtEpNkVoglIa0ilMpaYNKYyKrlmha0YLapu0L+vQCWpW0kRDUJQgHFUhKS2OhtBDlQWmrJs2GQMhDUUwakpUjvJDEbRVie3f+fTFn12Nn1jPrndnxHH8/0mrOOffMvefM7v58fWfunshMJEn1aox6AJKk4TLoJalyBr0kVc6gl6TKGfSSVDmDXpIqZ9BLUuUMekmqnEEvSZWbGPUAADZs2JCbN28e9TAkaaw8+OCDP8jMqV79Toqg37x5MzMzM6MehiSNlYj4Xj/9vHQjSZUz6CWpcga9JFXOoJekyhn0klQ5g16SKmfQS1LlxjroH3jqOT7x9e9waL416qFI0klrrIP+G997nhvv2st8y6CXpOWMddBLknoz6CWpcga9JFXOoJekyhn0klQ5g16SKmfQS1LlDHpJqlwVQZ856hFI0slrrIM+YtQjkKST31gHvSSpN4Nekipn0EtS5Qx6SapcX0EfEU9FxLcj4psRMVPaXhsRd0TEE+XxzNIeEXFjROyNiIcj4qJhTkCSdHwrOaN/Z2ZekJnTpX49cGdmbgHuLHWA9wJbytdO4FODGqwkaeVWc+lmG7C7lHcDV3S035xt9wFnRMS5qziOJGkV+g36BL4eEQ9GxM7Sdk5mPgtQHs8u7RuBZzqeO1vajhIROyNiJiJm5ubmTmz0kqSeJvrs9/bM3BcRZwN3RMR/Hadvt9uYXnbvambuAnYBTE9Pr+reVm+MlaTl9XVGn5n7yuN+4CvAxcD3Fy/JlMf9pfsscF7H0zcB+wY14E7R9d8USVKnnkEfEa+KiNMXy8B7gEeAPcD20m07cFsp7wGuKZ++2QocWLzEI0lae/1cujkH+Eq0/7DMBPCFzPyXiHgAuDUidgBPA1eV/rcDlwF7gReBawc+aklS33oGfWY+Cby1S/sPgUu6tCdw3UBGJ0laNe+MlaTKGfSSVDmDXpIqZ9BLUuUMekmqXBVBny4aK0nLGuugd81YSeptrINektSbQS9JlTPoJalyBr0kVc6gl6TKGfSSVDmDXpIqZ9BLUuWqCHrvi5Wk5VUR9JKk5Rn0klQ5g16SKmfQS1LlDHpJqpxBL0mVM+glqXIGvSRVzqCXpMpVEfQuGStJyxvroA8XjZWknvoO+ohoRsRDEfHVUn99RNwfEU9ExC0Rsa60n1bqe8v2zcMZuiSpHys5o/8Q8HhH/ePADZm5BXge2FHadwDPZ+YbgBtKP0nSiPQV9BGxCfh54DOlHsC7gC+XLruBK0p5W6lTtl8SXmORpJHp94z+r4DfBVqlfhbwQmbOl/ossLGUNwLPAJTtB0r/o0TEzoiYiYiZubm5Exy+JKmXnkEfEb8A7M/MBzubu3TNPrYdacjclZnTmTk9NTXV12AlSSs30UeftwOXR8RlwHrgx2if4Z8RERPlrH0TsK/0nwXOA2YjYgJ4DfDcwEcuSepLzzP6zPz9zNyUmZuBq4G7MvNXgLuBK0u37cBtpbyn1Cnb78r0k+6SNCqr+Rz97wEfjoi9tK/B31TabwLOKu0fBq5f3RAlSavRz6WbJZl5D3BPKT8JXNylz0vAVQMY2woGtqZHk6SxMt53xo56AJI0BsY66CVJvRn0klQ5g16SKmfQS1LlDHpJqpxBL0mVM+glqXJVBH16x5QkLWusg96/ci9JvY110EuSejPoJalyBr0kVc6gl6TKGfSSVDmDXpIqZ9BLUuUMekmqXBVB79LjkrS8sQ56b4yVpN7GOuglSb0Z9JJUOYNekipn0EtS5Qx6SaqcQS9JlTPoJalyPYM+ItZHxH9GxLci4tGI+OPS/vqIuD8inoiIWyJiXWk/rdT3lu2bhzsFSdLx9HNGfxB4V2a+FbgAuDQitgIfB27IzC3A88CO0n8H8HxmvgG4ofQbKm+MlaTl9Qz6bPu/Up0sXwm8C/hyad8NXFHK20qdsv2SiOGs7jqk3UpSVfq6Rh8RzYj4JrAfuAP4LvBCZs6XLrPAxlLeCDwDULYfAM7qss+dETETETNzc3Orm4UkaVl9BX1mLmTmBcAm4GLgTd26lcdup9kvu7qSmbsyczozp6empvodryRphVb0qZvMfAG4B9gKnBERE2XTJmBfKc8C5wGU7a8BnhvEYCVJK9fPp26mIuKMUn4F8G7gceBu4MrSbTtwWynvKXXK9rsy/UPCkjQqE727cC6wOyKatP9huDUzvxoRjwFfiog/Ax4Cbir9bwI+HxF7aZ/JXz2EcUuS+tQz6DPzYeDCLu1P0r5ef2z7S8BVAxmdJGnVvDNWkipn0EtS5aoIet/rlaTljXXQe2OsJPU21kEvSerNoJekyhn0klQ5g16SKmfQS1LlDHpJqpxBL0mVM+glqXJVBL33xUrS8sY66L0xVpJ6G+uglyT1ZtBLUuUMekmqnEEvSZUz6CWpcga9JFXOoJekyhn0klS5KoLeJWMlaXnjHfQuGitJPY130EuSejLoJalyPYM+Is6LiLsj4vGIeDQiPlTaXxsRd0TEE+XxzNIeEXFjROyNiIcj4qJhT0KStLx+zujngY9k5puArcB1EfFm4HrgzszcAtxZ6gDvBbaUr53ApwY+aklS33oGfWY+m5nfKOX/BR4HNgLbgN2l227gilLeBtycbfcBZ0TEuQMfuSSpLyu6Rh8Rm4ELgfuBczLzWWj/YwCcXbptBJ7peNpsaZMkjUDfQR8Rrwb+AfitzPyf43Xt0vayT7pHxM6ImImImbm5uX6HIUlaob6CPiImaYf832XmP5bm7y9ekimP+0v7LHBex9M3AfuO3Wdm7srM6cycnpqaOtHxS5J66OdTNwHcBDyemZ/o2LQH2F7K24HbOtqvKZ++2QocWLzEMyzpqrGStKyJPvq8HfhV4NsR8c3S9gfAx4BbI2IH8DRwVdl2O3AZsBd4Ebh2oCPu4H2xktRbz6DPzH9j+Uy9pEv/BK5b5bgkSQPinbGSVDmDXpIqZ9BLUuUMekmqnEEvSZUz6CWpcga9JFWujqD3xlhJWtZYB71LxkpSb2Md9JKk3gx6SaqcQS9JlTPoJalyBr0kVc6gl6TKGfSSVDmDXpIqV0XQe2OsJC1vrIM+XDVWknoa66CXJPVm0EtS5Qx6SaqcQS9JlTPoJalyBr0kVc6gl6TKGfSSVLmeQR8Rn42I/RHxSEfbayPijoh4ojyeWdojIm6MiL0R8XBEXDTMwS9Kb42VpGX1c0b/OeDSY9quB+7MzC3AnaUO8F5gS/naCXxqMMPszjVjJam3nkGfmfcCzx3TvA3YXcq7gSs62m/OtvuAMyLi3EENVpK0cid6jf6czHwWoDyeXdo3As909JstbZKkERn0m7HdLqZ0vYIeETsjYiYiZubm5gY8DEnSohMN+u8vXpIpj/tL+yxwXke/TcC+bjvIzF2ZOZ2Z01NTUyc4DElSLyca9HuA7aW8Hbito/2a8umbrcCBxUs8kqTRmOjVISK+CLwD2BARs8BHgY8Bt0bEDuBp4KrS/XbgMmAv8CJw7RDGLElagZ5Bn5nvW2bTJV36JnDdagclSRoc74yVpMpVEfTpqrGStKyxDnpvjJWk3sY66CVJvRn0klQ5g16SKmfQS1LlDHpJqpxBL0mVM+glqXIGvSRVroqgd81YSVreWAe9a8ZKUm9jHfSSpN4MekmqnEEvSZUz6CWpcga9JFXOoJekyhn0klS5KoLe+6UkaXljHfThYoKS1NNYB70kqTeDXpIqZ9BLUuUmRj0ASRoXC63k0HyLg/MLHJxvcfDwkfJLhxd46XCLHx1e4KXDC/zo8AIHy2Nn+1K/Qwu8NL/A9rdt5p0/cfZQx23QSzqptVrJ4VaLwwvJ4fkWhxdaHFpo1w911udL28JCCeDugdyud5Tn++8/3zrxz/hNNoP1E03Wr2uyfrLBKyabrJ9scnB+YYCvVndDCfqIuBT4a6AJfCYzPzaM40jqLTOZbyULreTwQov5hXZ9vnWkvFCCtF1vB9qx5aU+rc59JPML7bDtDN7F+sGj6i0Ozecx9RaHFjraOuqL+1pNuB4rAk6baHDaRLP9ONlRLu2vetXEcfo0S7379vWTzRLg7fL6ySavWNdk/USDieborpQPPOgjogl8EvhZYBZ4ICL2ZOZjgz6W1Kkz0Fql3Cr1hSyPraTVgvlWi1YmCy2W2jv7LO6js37UPkp56Rhdnnvk+LDQapU+pdyi5z6ODtSXh2u3UG6HedneOvJ6rKXJZjDZbCx9rWsG6yaO1Ccn2m2TzQavXDfBZLMdnEvPm2iwrtk4aj/rOtsmFve7uM9222nluZPNBuuPDfDJdnmiEcQpuJDFMM7oLwb2ZuaTABHxJWAbcNIHfWbSyvYvYCuTTJZ+aVvZ3t6uH913odXuu/i8Vud+WnRpO/o4nX26HbOVsJBZyu3+R+rd9kmpd9+e2Q6VI3WWgqZzHgstjhyz2z6PGvvR8zzyWh0zj1b2eO26vy7Hvv6L4+0MyZN9pbFmI2hGtB8bQSNgotmgEUGzQXtbs92n0QgmGw0mmsFEsx1QE41g/WSDiUapH7Wto63R0d7s3LbSPo0jfTvLXZ53JMhPzSA92Q0j6DcCz3TUZ4GfHMJxlvzcDfcydfppS8Gy4rDqCLxTQQQ0oh007ceOcuNIOUoALfaJoARULO2j2VFulLCKY/fdgEaj8bJ9NIL2MUqfKP2bHcdvlP5H7xMai6HZPDo8F0OyWYKx0bmtY/tEI47sY2k7NBuNpfE0ox1mx91HdBzrmG2L81zcJo3KMIK+20/0yyI0InYCOwHOP//8EzrQT73hLH7xwo0cWmi1f8H6DavSpzOsFkNn8exqsXy8oOkWVkvHPGqfXQKy4/jdjnlUSL4slJcJzEaPY5ayZ1zSqWUYQT8LnNdR3wTsO7ZTZu4CdgFMT0+f0Ln0pjNfyQ2/fMGJPFWSThnDeBv4AWBLRLw+ItYBVwN7hnAcSVIfBn5Gn5nzEfEbwNdof7zys5n56KCPI0nqz1A+R5+ZtwO3D2PfkqSV8W/dSFLlDHpJqpxBL0mVM+glqXIGvSRVLvIkuO8/IuaA753g0zcAPxjgcMaBcz41OOdTw2rm/LrMnOrV6aQI+tWIiJnMnB71ONaScz41OOdTw1rM2Us3klQ5g16SKldD0O8a9QBGwDmfGpzzqWHocx77a/SSpOOr4YxeknQcYxP0EXFpRHwnIvZGxPVdtp8WEbeU7fdHxOa1H+Vg9THnD0fEYxHxcETcGRGvG8U4B6nXnDv6XRkRGRFj/wmNfuYcEb9UvtePRsQX1nqMg9bHz/b5EXF3RDxUfr4vG8U4ByUiPhsR+yPikWW2R0TcWF6PhyPiooEOIMsaoifzF+0/d/xd4MeBdcC3gDcf0+fXgU+X8tXALaMe9xrM+Z3AK0v5g6fCnEu/04F7gfuA6VGPew2+z1uAh4AzS/3sUY97Dea8C/hgKb8ZeGrU417lnH8auAh4ZJntlwH/THuFvq3A/YM8/ric0S8tOJ6Zh4DFBcc7bQN2l/KXgUtivNfM6znnzLw7M18s1ftor+Y1zvr5PgP8KfDnwEtrObgh6WfOvwZ8MjOfB8jM/Ws8xkHrZ84J/Fgpv4Yuq9SNk8y8F3juOF22ATdn233AGRFx7qCOPy5B323B8Y3L9cnMeeAAcNaajG44+plzpx20zwjGWc85R8SFwHmZ+dW1HNgQ9fN9fiPwxoj494i4LyIuXbPRDUc/c/4j4P0RMUt7bYvfXJuhjcxKf99XZCgLjwxBPwuO97Uo+Rjpez4R8X5gGviZoY5o+I4754hoADcAH1irAa2Bfr7PE7Qv37yD9v/a/jUi3pKZLwx5bMPSz5zfB3wuM/8yIt4GfL7MuTX84Y3EUPNrXM7o+1lwfKlPREzQ/u/e8f6rdLLra5H1iHg38IfA5Zl5cI3GNiy95nw68Bbgnoh4iva1zD1j/oZsvz/bt2Xm4cz8b+A7tIN/XPUz5x3ArQCZ+R/Aetp/E6ZWff2+n6hxCfp+FhzfA2wv5SuBu7K8yzGmes65XMb4G9ohP+7XbaHHnDPzQGZuyMzNmbmZ9vsSl2fmzGiGOxD9/Gz/E+033omIDbQv5Ty5pqMcrH7m/DRwCUBEvIl20M+t6SjX1h7gmvLpm63Agcx8dlA7H4tLN7nMguMR8SfATGbuAW6i/d+7vbTP5K8e3YhXr885/wXwauDvy/vOT2fm5SMb9Cr1Oeeq9DnnrwHviYjHgAXgdzLzh6Mb9er0OeePAH8bEb9N+xLGB8b5xC0ivkj70tuG8r7DR4FJgMz8NO33IS4D9gIvAtcO9Phj/NpJkvowLpduJEknyKCXpMoZ9JJUOYNekipn0EtS5Qx6SaqcQS9JlTPoJaly/w8y9Gsmh0VvzgAAAABJRU5ErkJggg==\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"Maybe it is O(n^4)?\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAW4AAAEDCAYAAAAVyO4LAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvhp/UCwAADUtJREFUeJzt3X2MZXV9x/H3BxakVqykO20MT6NGaQmJQiaokFBFY5A28I9tlpSWpptuoNbYtEnTxn/68F+TamNCHzYtUdsKgq0tIWprhS1KADsrC/JQEQHrRtIdSkFJUyvw7R/3QrbLzD1nl3vuPb/l/Uom3Nk5M/P9McObu2fOmV+qCklSO45Z9gCSpMNjuCWpMYZbkhpjuCWpMYZbkhpjuCWpMYOFO8k1SQ4kubfHsR9Jsm/68mCSJ4eaS5Jal6Gu405yAfA08ImqOusw3u8DwNlV9cuDDCZJjRvsGXdV3Qo8cfCfJXlDks8n2ZvkS0l+YpN3vQy4dqi5JKl12xb8+XYDV1bVN5K8FfgT4MLn35jkdOB1wM0LnkuSmrGwcCd5FXAecEOS5//4FYcctgP4dFU9u6i5JKk1i3zGfQzwZFW9ZcYxO4D3L2geSWrSwi4HrKrvAo8k+VmATLz5+bcnOQM4Cbh9UTNJUouGvBzwWiYRPiPJ/iQ7gZ8Hdia5G7gPuPSgd7kMuK78dYWSNNNglwNKkobhnZOS1JhBfji5ffv2Wl1dHeJDS9JRae/evY9X1UqfYwcJ9+rqKuvr60N8aEk6KiX5Vt9jPVUiSY0x3JLUGMMtSY0x3JLUGMMtSY0x3JLUGMMtSY0ZVbg/+sVv8C8Pbix7DEkatVGF+0/3fJPbHnp82WNI0qiNKtySpG6GW5IaY7glqTGGW5IaY7glqTGGW5IaY7glqTGGW5IaM7pwu3mxJM02qnAny55AksZvVOGWJHUz3JLUGMMtSY0x3JLUmG19DkryKPA94FngmapaG3IoSdLWeoV76p1V5S/LlqQl81SJJDWmb7gL+Kcke5Ps2uyAJLuSrCdZ39hw+zFJGkrfcJ9fVecA7wXen+SCQw+oqt1VtVZVaysrK0c8kDdOStJsvcJdVd+Z/vMA8Bng3CGG8cZJSerWGe4kP5zkxOcfA+8B7h16MEnS5vpcVfLjwGcy+UUi24BPVtXnB51KkrSlznBX1cPAmxcwiySpBy8HlKTGGG5JaozhlqTGGG5JaozhlqTGjC7c3jgpSbONKtxx00lJ6jSqcEuSuhluSWqM4ZakxhhuSWqM4ZakxhhuSWqM4ZakxhhuSWrM6MLtnpOSNNuowu19k5LUbVThliR1M9yS1BjDLUmNMdyS1BjDLUmNMdyS1BjDLUmNMdyS1JjRhbvcdVKSZhpXuL11UpI69Q53kmOT3JXkpiEHkiTNdjjPuD8IPDDUIJKkfnqFO8kpwE8DfzHsOJKkLn2fcf8x8FvAc1sdkGRXkvUk6xsbG3MZTpL0Yp3hTvIzwIGq2jvruKraXVVrVbW2srIytwElSf9fn2fc5wOXJHkUuA64MMlfDzqVJGlLneGuqt+pqlOqahXYAdxcVZcPPpkkaVPjuo5bktRp2+EcXFV7gD2DTPLC5xjyo0tS+0b1jNsbJyWp26jCLUnqZrglqTGGW5IaY7glqTGGW5IaY7glqTGGW5IaY7glqTGjCnfiLTiS1GVU4ZYkdTPcktQYwy1JjTHcktQYwy1JjTHcktQYwy1JjTHcktSY0YW73LtMkmYaVbi9cVKSuo0q3JKkboZbkhpjuCWpMYZbkhpjuCWpMYZbkhpjuCWpMZ3hTnJCkq8kuTvJfUl+bxGDSZI2t63HMd8HLqyqp5McB3w5yeeq6o4hBvK+SUmarTPcNbkH/enpq8dNXwbpqzdOSlK3Xue4kxybZB9wAPhCVd25yTG7kqwnWd/Y2Jj3nJKkqV7hrqpnq+otwCnAuUnO2uSY3VW1VlVrKysr855TkjR1WFeVVNWTwB7gokGmkSR16nNVyUqS10wf/xDwbuDfhh5MkrS5PleVvBb4eJJjmYT++qq6adixJElb6XNVyT3A2QuYRZLUg3dOSlJjDLckNWZ04XbLSUmabVThjptOSlKnUYVbktTNcEtSYwy3JDXGcEtSYwy3JDXGcEtSYwy3JDXGcEtSY0YX7nLXSUmaaVTh9r5JSeo2qnBLkroZbklqjOGWpMYYbklqjOGWpMYYbklqjOGWpMYYbklqzOjC7Z6TkjTbqMLtlpOS1G1U4ZYkdTPcktQYwy1JjekMd5JTk9yS5IEk9yX54CIGkyRtbluPY54BfrOqvprkRGBvki9U1f0DzyZJ2kTnM+6qeqyqvjp9/D3gAeDkoQeTJG3usM5xJ1kFzgbu3ORtu5KsJ1nf2NiYz3SSpBfpHe4krwL+Fvj1qvruoW+vqt1VtVZVaysrK/OcUZJ0kF7hTnIck2j/TVX93ZADeeOkJM3W56qSAH8JPFBVHx52HG+dlKQufZ5xnw/8AnBhkn3Tl4sHnkuStIXOywGr6sv4VFiSRsM7JyWpMYZbkhpjuCWpMYZbkhpjuCWpMYZbkhozunC756QkzTaqcLvnpCR1G1W4JUndDLckNcZwS1JjDLckNcZwS1JjDLckNcZwS1JjDLckNWaE4fbWSUmaZVTh9sZJSeo2qnBLkroZbklqjOGWpMYYbklqjOGWpMYYbklqjOGWpMYYbklqzOjC7Z6TkjRbZ7iTXJPkQJJ7hx7GPSclqVufZ9wfAy4aeA5JUk+d4a6qW4EnFjCLJKmHuZ3jTrIryXqS9Y2NjXl9WEnSIeYW7qraXVVrVbW2srIyrw8rSTrE6K4qkSTNZrglqTF9Lge8FrgdOCPJ/iQ7hx9LkrSVbV0HVNVlixhEktTP6E6VeOekJM02qnDHXSclqdOowi1J6ma4JakxhluSGmO4JakxhluSGmO4JakxhluSGmO4Jakxowt34a2TkjTLqMLtnpOS1G1U4ZYkdTPcktQYwy1JjTHcktQYwy1JjTHcktQYwy1JjRlduN26TJJmG1W4vf9GkrqNKtySpG6GW5IaY7glqTGGW5IaY7glqTGGW5Ia0yvcSS5K8vUkDyX57aGHkiRtrTPcSY4FrgbeC5wJXJbkzKEHkyRtbluPY84FHqqqhwGSXAdcCtw/72GeK7hh7372ffvJeX9oSRrcSa88nuuvfPvgn6dPuE8Gvn3Q6/uBtx56UJJdwC6A00477YiGueodb+DOR/7ziN5Xkpbt1Scct5DP0yfcm92J/qLfKFJVu4HdAGtra0f0G0euOG+VK85bPZJ3laSXjT4/nNwPnHrQ66cA3xlmHElSlz7h/lfgjUlel+R4YAdw47BjSZK20nmqpKqeSfJrwD8CxwLXVNV9g08mSdpUn3PcVNVngc8OPIskqQfvnJSkxhhuSWqM4ZakxhhuSWpMaoDdeZNsAN86wnffDjw+x3Fa4JqPfi+39YJrPlynV9VKnwMHCfdLkWS9qtaWPcciueaj38ttveCah+SpEklqjOGWpMaMMdy7lz3AErjmo9/Lbb3gmgczunPckqTZxviMW5I0g+GWpMYsLdxdGxAneUWST03ffmeS1cVPOT891vsbSe5Pck+SLyY5fRlzzlPfTaaTvC9JJWn+0rE+a07yc9Ov9X1JPrnoGeetx/f2aUluSXLX9Pv74mXMOS9JrklyIMm9W7w9ST46/fdxT5Jz5j5EVS38hcmvh/0m8HrgeOBu4MxDjvlV4M+mj3cAn1rGrAtc7zuBV04fX9XyevuueXrcicCtwB3A2rLnXsDX+Y3AXcBJ09d/bNlzL2DNu4Grpo/PBB5d9twvcc0XAOcA927x9ouBzzHZPextwJ3znmFZz7hf2IC4qv4XeH4D4oNdCnx8+vjTwLuSbLaNWgs611tVt1TVf09fvYPJTkMt6/M1BvgD4A+B/1nkcAPps+ZfAa6uqv8CqKoDC55x3vqsuYBXTx//CI3voFVVtwJPzDjkUuATNXEH8Jokr53nDMsK92YbEJ+81TFV9QzwFPCjC5lu/vqs92A7mfwfu2Wda05yNnBqVd20yMEG1Ofr/CbgTUluS3JHkosWNt0w+qz5d4HLk+xn8nv9P7CY0ZbmcP97P2y9NlIYQJ8NiHttUtyI3mtJcjmwBvzUoBMNb+aakxwDfAT4pUUNtAB9vs7bmJwueQeTv1V9KclZVfXkwLMNpc+aLwM+VlV/lOTtwF9N1/zc8OMtxeDtWtYz7j4bEL9wTJJtTP6KNeuvJ2PWa8PlJO8GPgRcUlXfX9BsQ+la84nAWcCeJI8yORd4Y+M/oOz7ff0PVfWDqnoE+DqTkLeqz5p3AtcDVNXtwAlMfhnT0WrwDdaXFe4+GxDfCFwxffw+4OaanvlvUOd6p6cN/pxJtFs/7wkda66qp6pqe1WtVtUqk/P6l1TV+nLGnYs+39d/z+QH0STZzuTUycMLnXK++qz534F3AST5SSbh3ljolIt1I/CL06tL3gY8VVWPzfUzLPEnsxcDDzL5ifSHpn/2+0z+44XJF/cG4CHgK8Drl/3T5IHX+8/AfwD7pi83Lnvmodd8yLF7aPyqkp5f5wAfBu4HvgbsWPbMC1jzmcBtTK442Qe8Z9kzv8T1Xgs8BvyAybPrncCVwJUHfY2vnv77+NoQ39fe8i5JjfHOSUlqjOGWpMYYbklqjOGWpMYYbklqjOGWpMYYbklqzP8B6mo+FrvxyRQAAAAASUVORK5CYII=\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"f = lambda x: np.exp(5*x)\n",
"deriv_analyt = lambda x:5*np.exp(5*x)\n",
"def deriv_richardson(x, func, h):\n",
" Zh = (func(x+h)-func(x-h))/(h+h)\n",
" Zhq = (func(x+h/2)-func(x-h/2))/h\n",
" z_star = (4*Zhq-Zh)/3\n",
" return z_star\n",
"\n",
"hs = np.arange(1e-5, 1,1e-5)\n",
"richardson = np.vectorize(lambda h:deriv_richardson(0, f,h))(hs)\n",
"error = abs(richardson-deriv_analyt(0))\n",
"\n",
"for p in range(1,5):\n",
" print(\"Maybe it is O(n^%d)?\"%p)\n",
" plt.plot(hs, error/(hs**p))\n",
" plt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Прямая линия на последнем графике говорит о том, что ошибки являются функцией Const*n^4 = O(n^4)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 5. lnx не определена при x=0"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 6. Multidimension"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"(array([3.00000025, 3.00000025]),\n",
" array([0., 0.]),\n",
" array([0. , 3.00000025]))"
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def gradient_richardson(point, func, h, dims):\n",
" \"\"\"\n",
" point is list or ndarray\n",
" func is working with list and return float\n",
" \"\"\"\n",
" larg, rarg = np.array(point, dtype=float),np.array(point, dtype=float)\n",
" r = np.zeros(dims)\n",
" for i in range(dims):\n",
" larg[i]+=h; rarg[i]-=h\n",
" Zh = (func(larg)-func(rarg))/(h+h)\n",
" larg[i], rarg[i] = point[i]+h/2,point[i]-h/2\n",
" Zhq = (func(larg)-func(rarg))/h\n",
" z_star = (4*Zhq-Zh)/3\n",
" larg[i], rarg[i] = point[i],point[i]\n",
" r[i] = z_star\n",
" return r\n",
"\n",
"def f(point):\n",
" return point[0]**3 + point[1]**3\n",
"\n",
"gradient_richardson((1,1),f, 1e-10,2),\\\n",
"gradient_richardson((0,0),f, 1e-10,2),\\\n",
"gradient_richardson((0,1),f, 1e-10,2)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 7. Сделано в 4 и 6"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"1.0"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"1e-10*10**10"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"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.5"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment