Untitled.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAEKCAYAAADuEgmxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3WmYFNX59/HvzcAw7CIIiCCg4AKCsogSNRJForgbNGoStxg0xpDl+ZtoEpeQGM1mjCYaNRo1UdRoVFRcAMG4RBAQ2QRFREGUHWTYZOA8L840VT3M9PT0VtM1v891zeWp6uqqc+yi7z5rmXMOERGRdDSKOgMiIlI8FDRERCRtChoiIpI2BQ0REUmbgoaIiKRNQUNERNKmoCEiImlT0BARkbQpaIiISNoaR52BXGvfvr3r3r17ja9v2rSJFi1aFC5DBRDHMkE8y6UyFYc4lglSl2vGjBmrnXN71XaO2AWN7t27M3369BpfnzJlCkOHDi1chgogjmWCeJZLZSoOcSwTpC6XmX2UzjnUPCUiImlT0BARkbQpaIiISNoUNEREJG0KGiIikjYFDRERSVvshtyKNHjOwYcfwsqVYAb9+0NpadS5krpIPFHVLNp8VEM1DZG4cA4eewwGDYL994chQ+DII2GvveCii2Dp0qhzKKns3AkPPQSnngotWvhA36kTjBoFCxdGnbtdFDRE4mD7dvj2t+HrX4eZM5Nf+/xzeOAB6N0b7rwz+BUr9cfcuXDUUfDNb8Kzz8KWLVBRAStWwD33+M/uD3+oF5+dgoZIsauogDPPhH/8I9hXVgaDB0O3bsG+8nK44gr48Y/9r1qpH8aP9zXCN9+s+ZidO+Gqq+BHP4o8cChoiBS766+H554Ltr/1Lfj4Y5g61fdtTJ7sf6km3Hqrr5XUg1+tDd7DD8Npp8GmTX67SRO45hpYtAi2bYMpU3xASfjzn+G22yLJaoKChkgxmzABbrop2L7qKt8UtVflunNmMHQozJgBX/tacNz99/svJ4nO+PFwwQWwY4ff7tbNNy3+5je+T6q0FI491gf98Gd31VWQYn29fFPQEClWW7fCpZcGNYYTToCbb65+xE1ZGTz6KFxySbDvt7+Fu+4qTF4l2YwZMHJkEDAOOQT+9z//36rKynwH+YABfnv7dt/3sX174fIboqAhUqxuv903QwG0bw///Cc0SvFPuqTEB4nTTgv2ff/7vhmr2Djn+2iKsW9m1SrfB7Vli9/u3h1efBH23rvm9zRt6oN+q1Z+e+FC30EeAQUNkWK0Zg3ceGOw/ctfQseOtb+vcWMYOxYOPdRvb9/umz7WrMlPPnPFOfjvf2H0aN8/06yZ/wItLfVfuhddBE8+6QcF1Gc7dsC55wbDn9u0gRdegM6da39vz55w7bXB9g03wMaNeclmKgoaIsXolltgwwaf7tULvvOd9N/bvDn85z/Qtq3f/uQTuOyy+tsx/sILcPjhvn3/9tvh3Xd9JzH4L+GPPvL9OGed5f9f3HNP/a2B/O538PLLPm3mm50OPDD993//+7Dvvj69apW/DwpMQUOk2GzZktwX8etf+1E3dbHffvDgg8H2E0/4zvH6ZOVKOPtsOOkk3wdQVdOmu+9bssRPhvvSl2D+/LxnsU7eeguuuy7Y/sUv4OST63aOsjL41a+C7b/+1fdtFZCChkixGTs2aE7q1i15ZE1dnHIKXH55sD16tB/qWR+89ppf/uTxx4N9ZWU+IEyY4Ccsbt3qA+jUqfCzn0G7dsGxU6f62kk4MEapvBzOPz9oPhsyJDmA1MX550PXrj69ahU88khu8pgmBQ2RYuJc8jj9K6/0HdyZ+uMfg+aR8vJIR+Xscv/98JWvwPLlwb4LL/S1iLvugmHDgg7hxCTGG2/0gwKuvTZYZ2vzZv++X/wi+qa3H/wgCMitWvlmqcYZLv3XuDF873vB9p//XNDyKWiIFJM334R33vHp5s39JL1sNG+e/AU2dapv7orKjTfCxRcHv8jbtfPzGe6/v/aO/ubNYcwYP4fhoIOSz3nZZdF1kj/+ONx3X7B9xx3Qo0d257z0Uj8YAGDWLHjjjezOVweRBg0zu8/MVprZ3BpeNzO7zcwWmdlsMxtQ6DyK1Cv/+leQPu+8oDM7GwMHJreT//rXBf0S2mXMGF8rSOjXz092O+mkup2nb1/ffzBiRLDvnnvgnHMK3v7P0qXJgxTOP9/X5rLVrp0/V8I//5n9OdMUdU3jfuDEFK+fBPSq/BsF3FmAPInUT9u3+1VsE3Lx5ZNw1VV+dBL4kUff/GYwOqsQfv1rvxxKwvHHw6uvBiOF6qplS3jqKb+kSsKTT/o5Kon5Efm2Y4e//vr1frtbN1/LyJULLwzS//43fPFF7s6dQqTP03DO/dfMuqc45HTgQeecA940sz3MbG/n3KcFyaBIjuzc6VeDeO65zixYkNk59p3zEiNWrwagvG0X/jX/y5DhuXZXQsvhD3L2tH403bIBPvyQD477DhMufbTWZzq8917mZQLo//xvOGJcMP9g6cHDeeH0p9nxcFnmJwWgCRx5P0d+2oHDJv7R75owgaUDz+CF7z7NjiY1nz/bMgEMemYMg155BYCd1ohxZz/EZ2PbZHfSsJ1H8Y0996XV2o9h7Vqe/8ELfHToaSnfkoty4ZyL9A/oDsyt4bVngaND25OAQanON3DgQJfK5MmTU75ejOJYJufiVa4773TO91Zm/vcQ5+3a+C1XZX2+6v5G8ljSjsu4My/XSfz9lJuSdrzEMFfG5pxf5+f8KmnHeE50TdmSt3IdzwS3A9u14zpuyMt1fsPVuzYe4Zwsz8d052r/zjbnXJZhJzuVNY1nnXO7LbpiZs8BNznnXqvcngT8xDk3o8pxo/DNV3Ts2HHgIymGoJWXl9OyZcuc5b8+iGOZIF7lGjOmN5Mnd8j4/U3Zymra0xK/GuqhzGI2h+Yqe0nu4Lt8l78BsJWmHMHUvFzr//g9v+cnu7YncjynMY4tNM/5tQCuZQxjuH7X9nhO4kye5Auqme+RhU58yiwOoyMrAZjEcQznJXaSxSi3GvRmHvPwX51bKKM9q9lMiwzPZjOcc4NqO6q+P+51GdA1tN0FWF71IOfc3cDdAIMGDXJDhw6t8YRTpkwh1evFKI5lgniVK9yUPWyYX8S0Lnp/PIWWz/uAsbJ1T4ac248heXoS6PyKW1j25Bt0WTubMrbxYptzuOms6WwrbVXt8cuXL6dzOstghAx754+MfDMIGAs6f4WnTxzHBU3yEzAAVnIdz07fySkzfgnACJ5n+r5ncdfw/1BRkhw4MikTQKOdFfzw2fPo+KkPGBuadeTFkQ/xnea5DxheHz55rA/7rJtHM7byh69O5J3up9d4dKpypb12ZTrVkXz+kbp56mTgecCAI4FptZ1PzVPxEadynXVW0Azw+OMZnOC73w1O8OMf5zx/u1mwwLkWLYJrnnaacxUV1R5a58/pt79NbhcZOtS58vLs85yu665Lvv7JJzu3dWvSIRnfe6NHB+dt1Mi5SZOyz29trrkmuOa3v53y0FTlIs3mqaiH3I4F/gccaGbLzOzbZna5mSWmqY4HFgOLgHuAKyLKqkhWEitgQ+qFaKvlHDzzTLB96qk5yVNKBx4If/tbsD1uHPzkJzUfnw7n/Dl++tNg35e/7B9v2iLTJpUM3HBD8tDe557zy5Uk1rPK1G23JU+8vO46OO647M6ZjvCqxc88k/d1t6IePXVeLa874HupjhEpBuF/x3WewD1rFixb5tNt2/pnSRfCN7/pJxL+4Q9++5ZbYI89kldaTVdFhV+y5N57g31Dh/ovuUIGDPCjwcaM8R/Kb37j9z3zDHz1q36Y7h571P2c997rZ30njByZ2f+nTAweDB06+LW6Vq6EadOSn/aXY1HP0xBpELKqaYRrGSNG1H1xwmzcfDOccUawfd11fj5FXQbQrF/vV6ANB4zTT4fnn/fzKaJg5ueGXH11sO+VV/yaUHVZ6NA5+Mtf/AzthCOO8Gte1fmDzlCjRn4dsYRx4/J7ubyeXUSALGsaEyYE6bquipqtkhK/QOIJJwT7xozxs9E//7z297/+up9xHg58F13kl9Yoy3YeRpbMfE3j5puDfQsWwOGH03ncuNqbeTZv9rWn738/2DdggF/2JLHER6GEmyxffDGvl1LQECmAjGsaGzf69aYShg3LWZ7SVlYGTz8Nw4cH+x59FA4+GB56CKtugcP58/2M5aOPhsWLg/0/+YmvcWS6WF+umfk+locfDr7oN2/mgD/9ydc6nn46+cMD3/fx8MP+YVB33x3sP/xwH+D33LNw+U847rjg18jbb0PlJNB8qCefnEi8hb936lTT+O9/g4X2DjsM9torp/lKW7NmvsP6Rz/yz3AAvwrtN7/JkLZtfTDbd1+/Uu7Uqb4fJqxVK79o38iRhc97Os47zz+f++yz/aNUwfcNnHGGDwJDhvj1nhJ9BmvX7v7+e+8tfA0joXVr3yz2xhu+yezll/1aW3mgoCFSABk3T02cGKSjqGWENWni2++PPhp++ENYsQKA0nXr/NpHNTnzTLj11szXkSqUvn39r/SbbmLnzTfTKFGDWrvWj7CqTrt2/ml8F19c63IreTdsWLDQ5MSJeQsaap4SKYCMm6fqU9BIOPdc/8jV66+v+dnWJSX+V/orr/hHy9b3gJHQrBmMGcPUhx7yzVYdapjF37WrH7q7cCFcckn0AQOS74/wfZNjqmmIFEBGNY3PPoO5lU8NKC31v/Dri7Ztd813mHHvvQxs3Ng/TbBZM+jZ0w8Lbt066lxmbNtee/kO8ptu8v0z8+bBpk1+tFffvnDAAYUbHZWuI47ww5c3bYIPP/R9Sfvtl/PLKGiIFEBGNY0pU4L0kCGFn8+QjsaN2XjggX7ORRyZQZ8+/q++Ky31y9uPH++3J03KS9CoZ6FSJJ4yqmm89lqQTjzrQiSVcPB+9dW8XEJBQ6QAMho9Ff5Hf8wxOc2PxFT4PlHQECledW6eWr8e5szx6ZKSvC4LITEyYEAw7HfJkmD5mRxS0BApgDo3T73+erBUR//+0S23IcWltDT5B0YeahsKGiIFUOeahpqmJFN5bqJS0BApgDrXNBQ0JFMKGiLFr041jW3bYPr0YLs+zc+Q+m/IkGBtr7lzff9YDiloiBRAnUZPzZoFX3zh0/vvH916U1KcWrSAfv2C7WnTcnp6BQ2RAqhT89TUqUFao6YkE+H7Jnw/5YCChkgB1Kl5KrwU+hFH5CU/EnPhoBG+n3JAQUOkAFTTkIIK/9iYOrVuT1qshYKGSAGkXdNYtSp4aFHTpnDooXnNl8RUr15+UUnwC0l+8EHOTq2gIVIAadc0wrWMAQP8ZC2RujLLWxOVgoZIAaQ9ekr9GZIr4ftHQUOkuKTdPDVjRpAePDhv+ZEGIBw0wvN+sqSgIVIAaTVPOZccNAYOzGueJObC98877wTPms+SgoZIAaRV01i+3HeEg1+gsGfPvOdLYmyvvaBLF5/euhUWLMjJaRU0RAogrZrG228H6f7969/jRKX4DBgQpGfOzMkpdVeKFEBaNY3wP+rwP3aRTCloiBSntEZPhf9R9++f1/xIA6GgIVKc0mqeUk1Dci18H82alXwjZkhBQyTPnEv+t1pt89SqVbB0qU+XlcHBBxckbxJznTtDhw4+vXEjzZYvz/qUChoieRZe9sfMYVbNQeFO8H79guchiGTDLKm20fK997I+pYKGSJ6F+zOqDRiQHDTUNCW5FLqfWr3/ftanU9AQybPk/owaVhtVJ7jkS7imUexBw8xONLOFZrbIzK6u5vWLzGyVmc2q/Ls0inyKZCN5uG0aQUM1DcmlqjWNLJdJjyxomFkJ8FfgJKA3cJ6Z9a7m0Eedc4dV/v29oJkUyYFam6c2bIBFi3y6cWM45JCC5EsaiO7dYY89AGjy+efw8cdZnS7KmsZgYJFzbrFz7gvgEeD0CPMjkhe1Nk/NmhWk+/Txo6dEcqVKZ3i28zWiDBr7AEtD28sq91X1NTObbWaPm1nXwmRNJHeSaxrVBA01TUm+5TBoRDmur7qKetV/Uc8AY51z28zscuAB4LjdTmQ2ChgF0LFjR6ZMmVLjRcvLy1O+XoziWCaIT7k2bGgCHAX4oFG1TAc9/zydKtPvt2zJJ0VW5rh8TmFxK1OHpk1JtP2vmTCBOccfn/nJnHOR/AFDgBdD29cA16Q4vgTYUNt5Bw4c6FKZPHlyyteLURzL5Fx8yvXZZ8753kfn9thj2+4H9OkTHPDaa4XPYJbi8jmFxa5M774b3GOdOlV7CDDdpfHdHWXz1FtALzPrYWalwLnAuPABZrZ3aPM04N0C5k8kJ1KOngovWW2mZ4JLfvTqBS1a+PRnn/m/DEUWNJxzFcCVwIv4YPCYc26emY0xs9MqDxttZvPM7B1gNHBRNLkVyVzKJUQWLgyiSo8e/jkaIrlWUpI8Km/OnIxPFelaBc658cD4KvuuC6WvwTdbiRStlDWNuXODtIbaSj716wdTp/r0nDlwwgkZnUYzwkXyLLmmoaAhEenbN0jPnp3xaRQ0RPIs5QOYFDSkUPr1C9JZNE8paIjkmWoaUi+Eaxrz50NFRUanUdAQybMa+zTKy2HJEp9u3BgOPLCg+ZIGZs892da+vU9v3RosXVNHChoieVbj2lPz5wfpAw6A0tKC5UkapvL99gs2MmyiUtAQybMa154KN0316VO4DEmDtalHj2Ajw85wBQ2RPKuxpqH+DCmwTfvvH2yopiFSP9XYEa6gIQVWHq5pKGiI1E81doQraEiBbd533+D584sXw8aNdT6HgoZInlU7T2PtWvj0U59u2hTCzQYieeJKS5NH6c2bV+dzKGiI5Fm1zVPhf6y9e/u1gUQKIcuZ4QoaInlWbfOUmqYkKlnODFfQEMmzale51XBbiUq4pqGgIVL/qKYh9Uq4pjF7tn80Ux0oaIjk2W59Gs4paEh0unaFNm18et06WL68Tm9X0BDJs91GT332mR89Bf6hS/vuG0m+pIEyS/6hUsfOcAUNkTzbrXmqai0jaZq4SAFk0RmuoCGSZ7s1T6lpSqIWvu/qOFdDQUMkz3Zbeyr8j1RBQ6IQHrGnoCFSv+y2yq1qGhK1cNCYPz/5Jq2FgoZIniX1abAz+Zed5mhIFNq3h44dfXrLFvjww7TfqqAhkmfhoNHpi6X+iX0A7doF/3BFCi1cyw3XfmuhoCGSZ+Ga//5b3g02NHJKohSu5SpoiNQf4ZpGj80Lgw31Z0iUMhxBpaAhkmfhmkaPTVVqGiJRUU1DpH4K1zS6l6umIfVEOGgsWJD22xQ0RPIsUdMooYKumxYFL2jklESpTRvo0sWnt29P+22N85QdEamUqGnszweU7tzmNzp3hrZto8uUCMApp8CaNf4HzA03pPUWBQ2RPEsEjUPQpD6pZ+68M0inGTTUPCWSZ4nmKQUNiYO0ahpm1gE4CugMbAHmAtOdc+nPPRdpoFTTkDhJGTTM7CvA1cCewNvASqAMOAPY38weB/7onPs83xkVKVaqaUic1FbTGAF8xzn3cdUXzKwxcApwAvBEJhc3sxOBPwMlwN+dczdXeb0p8CAwEFgDfN05tySTa4lEZccOKGUbvXg/2Nm7d3QZEslCyqDhnLsqxWsVwFOZXtjMSoC/4oPOMuAtMxvnnJsfOuzbwDrnXE8zOxf4LfD1TK8pEoUdO+BAFtKYynaq/faDFi2izZRIhtLqCDezf5pZm9B2dzOblOW1BwOLnHOLnXNfAI8Ap1c55nTggcr048DxZrUs1rNiBfzpT1lmTSR3du6s0jSl+RlSxNIdPfUaMNXMRpjZd4CXgFuzvPY+wNLQ9rLKfdUeU1mz2QC0S3nWZcvgttuyzJpI7uzYof4MiY+0Rk855+4ys3nAZGA10N8591mW166uxuAyOAYzGwWMAt/5wZIlvDp+PDuaN9/tzeXl5UyZMqWuea3X4lgmiE+5PvigG2eGgsb8Ro1YGYNyJcTlcwqLY5kgN+VKd8jtt4BrgQuAfsB4M7vYOfdOFtdeBnQNbXcBltdwzLLKjvc2wNqqJ3LO3Q3cDTDIzAEc064dHHHEbhedMmUKQ4cOzSLb9U8cywTxKdfkyck1jd7nnEPvfv0izFFuxeVzCotjmSA35Uq3eeprwNHOubHOuWuAy4H7s7oyvAX0MrMeZlYKnAuMq3LMOODCyvRI4GXn3G41jWrV8bm3IvlSsqWc/fBPRtthJXDggRHnSCRz6TZPnVFle5qZ7f4zvg6ccxVmdiXwIn7I7X3OuXlmNgY/cXAccC/wTzNbhK9hnJv2Beqw1K9IPrVbGSyHvqJNDzo3bRphbkSyU9vkvl8AdzjnqmsS+sLMjgOaO+eezeTizrnxwPgq+64LpbcCZ2dybgUNqS86rAzuxeVtD6BzhHkRyVZtNY05wDNmthWYCazCzwjvBRwGTAR+k9ccZkrNU1JPJAWNdr0izIlI9mrr0xjpnDsK34Q0D9+M9DnwL2Cwc+5HzrlVec5j3SSmcSxfDmt3qyCJFFzH1UHQ+GzPnhHmRCR7tdU0BppZN+AbwFeqvNYMv3hh/VJWBlsqszVvHhxzTLT5kQavUzhotD8gwpyIZK+2oPE34AVgP2B6aL/h50vsl6d8Za5ZsyBozJ2roCHRWruWPTb5keRbacratl1reYNI/Zayeco5d5tz7mD8yKb9Qn89nHP1L2CADxoJ6teQqIXuwXc5GFdSEmFmRLKX1jwN59x3852RnAkHDY2gkqiFgsZcDkExQ4pd/J7cVzVopDkXUCQvQj9c5nIIlQsWiBSt+AWN0lJo2dKn16zxq96KRKVK0GjUSEFDilv8ggYkLz2tfg2JinNJQWMefWgUz39x0oDE8xYOLz2tfg2JyooVvrYLbKQlH7OvahpS9BQ0RPKlSi3D0UhBQ4pePINGuHlKQUOiUqU/A1DzlBS9eN7C4ZrGvHkaQSXRqDLcFlBNQ4pePINGp06w554+vXEjLF2a+niRfFBNQ2Ionrewmfo1JFo7d6qmIbEUz6AB6teQaH30ka/lAp83accKOgIKGlL84hs0qvZriBTSnDm7kkta9cWv8anmKSl+8b2F1TwlUQoFjQ9b9t2VVk1Dil18g0a4eWr+fNixI7q8SMMTChqLWwRBQ2tPSbGLb9Bo186PogLYuhUWL442P9KwhILGBy367UprlVspdvENGqB+DYnGtm2wcOGuzQ/KglqvahpS7BpO0FC/hhTKu+8GzaH77ccma7nrJfVpSLGLd9DQsFuJQqhpir592bkz2FTzlBS7eAcN1TQkCuGg0a9f0hgMNU9JsYt30AjXNBYu9G3NIvk2e3aQrlLTUPOUFLt4B41WraBHD5+uqPBtzSL5VqV5KlzT0OQ+KXbxv4UPPTRIv/NOdPmQhmHtWli+3KebNoWePVXTkFiJf9A47LAgraAh+RauZfTpA40bV6lpKGhIcYt/0AjXNGbNii4f0jBU6c8A1DwlsRL/W7hq85QeyCT5VKU/A1DzlMRK/ING9+7QurVPr11L09WrI82OxFyV4bagmobES/xvYbOk2kbLRYsizIzE2s6dyfOBVNOQGIp/0ICkoNHigw8izIjE2pIlUF7u0+3bQ0f/4CVN7pM4iSRomNmeZjbBzN6v/G/bGo7bYWazKv/GZXzB0Agq1TQkb6r2Z5h/8FI4aGgZESl2UdU0rgYmOed6AZMqt6uzxTl3WOXfaRlfLdw8pZqG5Mvbbwfp0D2n5imJk6iCxunAA5XpB4Az8nq1Pn129UA2++QT2LQpr5eTBmrmzCA9YMCupJqnJE6iChodnXOfAlT+t0MNx5WZ2XQze9PMMg8szZrBQQcBYM4lNyOI5EoNQUOr3EqcNM7Xic1sItCpmpd+XofT7OucW25m+wEvm9kc59xu7UtmNgoYBdCxY0emTJmy24kO7tSJjvPnA/DeY4+xfOvWOmSjfisvL6+2zMWumMrVZN06jvrkEwB2NG3Ka599hlu1CoCtW78ElAKweXPxlCldxfQ5pSuOZYIclcs5V/A/YCGwd2V6b2BhGu+5HxhZ23EDBw501br5Zuf81D7nLr+8+mOK1OTJk6POQl4UVbleeCG4v448MumlPfcMXnrqqVcjymD+FNXnlKY4lsm51OUCprs0vr+jap4aB1xYmb4QeLrqAWbW1syaVqbbA0cB8zO+ohYulHyqoWkKNHpK4iWqoHEzcIKZvQ+cULmNmQ0ys79XHnMwMN3M3gEmAzc75zIPGuGFC2fPTm5oFslWmkFDHeFS7PLWp5GKc24NcHw1+6cDl1am3wD65uyinTpBhw6wcqUfPfXBB9CrV85OLw1ciqChIbcSJw1jRnhCuLahFW8lV9atg8WLfbpJk+QnRqK1pyReGtYtHP4FOGNGdPmQeAn/ADnkECgtTXpZNQ2Jk4YVNAYODNLTp0eXD4mX8EzwKk1ToJqGxEvDuoUHDQrSM2bo2RqSGyn6M5xTTUPipWEFjW7d2J54tsb69UE7tEg2agkaCWa71jAUKVoNK2iYsfGAA4Jt9WtItjZtggULfLpRo10PXkpQ05TETYO7jZOChvo1JFvhRwgfdBA0b570stadkrhpeEHjwAODDdU0JFspmqZANQ2JnwZ3G+8WNDQzXLJRS9BQTUPipsEFjW0dOsBee/mNDRvgvfeizZAUt2nTgnR4SHclrTslcdPgggZmcOSRwfabb0aXFylun38OlcvtU1JSa9BQ85TEQcO8jcNBY+rU6PIhxe2tt4JO8L59oUWL3Q5R85TETcMMGkccEaRV05BMhX9whO+pENU0JG4a5m18+OHBLKs5c/TMcMlM+AdHDUFDNQ2Jm4YZNFq3DlYi3bFDQ2+l7pxLrmmEmzxDVNOQuGm4t7GaqCQbS5b4Z7MAtGkD4aHcIRo9JXHTcING+Jfh//4XXT6kOIV/aAweXGM1Qs1TEjcNN2h86UtB+rXXtOKt1M1rrwXpIUNqPEzNUxI3Dfc2PvhgaNfOp1evDhadE0nHq68G6WOOqfEw1TQkbhpu0DCDo48OtsNfAiKprFsHc+f6dElJjZ3goJqGxE/Dvo3DvxAVNCRdr78eNGf27w8tW9Z4qGp4LiMiAAARYklEQVQaEjcKGgkKGpKucH9GiqYp0OgpiZ+GHTT69w+ef/DRR/Dxx9HmR4pDmv0ZoOYpiZ+GfRs3aZI88uWVV6LLixSHzZv9mlMJ4X6xaqh5SuKmYQcNgK98JUhPmhRdPmqzYYNfhnvaNL+ce5yeA1JRAbNmwQsvwDPP+KVdtm+POlfVe+21IG+9ewfL7NdANQ2Jm8ZRZyByw4bBL37h0xMn+g7OxLpUUVuzBu67D+6/P1iCO6FtWzjhBPje94p3jsmbb8Kf/wzPPQcbNya/1qwZBx1zjF859vDDo8lfdSZODNLDhtV6uGoaEjf67TNwoF8GAuCTT2DhwmjzA/7n6Z13Qs+e8JOf7B4wwA/7fOwxOPZY+o8e7Z9VXSwWL4aTT/ZNg488snvAANiyhU4vveRnW198sZ9LUx/UMWiopiFxo9u4cePkJqrwl0IU1q2DESPgiitg/fpgf5MmcMghPshVaRJpM3eu33/99cnfUvWNc3D33X6xyPHjk1/bZx9fcxo+HLp2TX7t/vuhX7/kvoQorF4Nb7/t0yUlcOyxtb5Fo6ckbhQ0wH9ZJUQZND780P+yfumlYN/++/smqnXrfFv/9OmwYoXvAxg1ygcT8N9OY8bAiSf6Zq36Zts2uOgiuOwy2LrV7zODb33LP2d76VJf7hdf9CPZpk1jVXhk0qefwpe/DP/5TyTZB+Dll4P0EUf41ZJroeYpiRsFDUhuZpg8OZpO2MWLYehQWLQo2Hf11X7m8cUXJz8VzgwOPRTuugtmz2Z9v37BaxMn+i/XTz4pWNZrtW4dfPWr8OCDwb6+fX2fxoMP+qHP4X4kMzj8cOaNGQPjxvn+G/DB5utfhyefLGz+EyZMCNJpNE2BmqckfnQbA/TqBd26+fTnnxd+ot+nn8JxxwXzRMrKfH/FTTf5dCoHHcSsW26Ba68Nvnjnz4ejjoL3389vvtOxZo1v/gsPZ774Yv8sisGDa3//qaf6EWMHHOC3KyrgnHOSa2OFsHNncpNauHZay9sSVNOQOFDQAP9le9ppwfa4cYW79qZN/ovxo4/8dlkZPP00nH12+ucoKfFNUw895PtowJ/v6KN9M1ZU1q71v8jDnfS/+x3cey80a5b+eXr29DXAXr38dkWF//8zb15u85vKzJmwfLlPt2uXcr2pMNU0JG4iuY3N7Gwzm2dmO81sUIrjTjSzhWa2yMyuzmumTj01SI8bV5hhrM75dv7EkwNLSnyb/fDhmZ3vvPN8wEl8Ia9c6Wswic7bQlq3zgeMRNAygwcegKuuymxIc+fOPnB06eK3P/8cTjnFX6cQnnkmSI8YEQTnWqimIXET1W+fucBZwH9rOsDMSoC/AicBvYHzzKx33nJ07LHQqpVPf/hhYX7F/uUv8PjjwfYdd8BJJ2V3zhEjfNNNYhhx1S/vQli/3jffJIKVGfzjH3DBBdmdd5994Nlng/6dJUvgkksKE+DDtc9wrbQWGj0lcRPJ5D7n3LsAlvoX52BgkXNuceWxjwCnA9VMWgisXOnni9Xk/ff3qWFKQykn7X8SB8x6DIA3rh7HWycckupSWen48Vuc/ef/R+J7ZNYxV/LKllGQIu812b1MR9Ph0kmceccwyrash7Vr2XLU8fznikms7nJYDnJfs9LN6znzb8Pp9HHw3PUJX/8789dfWOeyVf9ZHcr+5zzIKf/4mt986immjPwL73z5+1nlO5VW6z7mksqgu6OkCXd/MJwv0izLnDlBWs1TEgvOucj+gCnAoBpeGwn8PbT9LeAvtZ9zoPM/Pev+dz7/2rUxk8MyPk9tf21Y5xbTfdeOtxjoStma8+sMYLpbyx67dqxmT9ePWXkrV2vWu/9xRNLOb3NPXq51G1fu2thKqRvA9LyVazS37tp4geEZn+fMM52bPHmyixuVqXikKhcw3bnav7fzVtMws4lAp2pe+rlz7ul0TlHNPlfDtUYBo/zWwDRzuLtnOYUtlNGMrfRnFgczn3fJdYuY4z4uoQdLANhAa87hMb6gaY6vAzMZyAlMYAIn0Jb1tGMtkzie45nEbA7N6bVa8TkvcCJHMnXXvlHcxb1cmtPrJPwff+BLvMFAZtKUL3iUrzOAmWyk9rkTdfUNHtqVfoKvZXye1q2XUF5ezpQpU3KQq/pDZSoeuShX3oKGcy69gew1WwaEpwZ3AZbXcK27gbsBOnYc5M49N8VJly2jS6IzdTdtWDD+VPov+jcAtw56iOe+dGPdc57CUXPu4qzJwTyDJ0bcx6k998vqnKnLNIh7V0zgiieH0fyLDbRnDa+XHc9fz3qZ5e371fCeumm5eSWXP30SXVfN3LXvsaF30KzfKEZncd7U5WrK+PWP0mfsAMq2b6QnHzDxoCt5aPiDNRyfmb3Wv8/gB/1M9IpGTWh36UhG1zIKujr77w+XXNKd6dOXMHTo0JzmMWpTpkxRmYpETsqVTnUkX3+kbp5qDCwGegClwDtAn9rOOXDgwIyrZ8455556KmhP6NbNuR07Uh9fF/PnO9esWXD+730vJ6dNqyo9bZpzbdoE127f3rkZM7K/+JIlzh1wQHI7zF/+kv15XZrlGjs2+dpjx+bk2rtcf31w7tNPz/p0cWz2UJmKRy6ap6IacnummS0DhgDPmdmLlfs7m9l4AOdcBXAl8CLwLvCYcy7/Q5pOOimYgfzRR8lPacvGtm1+SOyWLX77kEPgD3/IzbnTcfjhflRVYumL1av9zPHnnsv8nPPm+UmE773ntxs1gnvu8SvvFsq55yaPyrr88tw9TMs5ePjhYPsb38jNeUWKWCRBwzn3pHOui3OuqXOuo3Puq5X7lzvnRoSOG++cO8A5t79zLrftRDUpLfVLVST89a+5Oe/PfhZMcmvaFMaOrX22d64NHuyXwthjD7+dmFj4y1/WfaHDsWP9BLfEciWlpfDvf8Ol+enDSOn226FHD5/esMEHkVws3Dh5cjCrvlUrPy9EpIHTIMDqfPe7QfqJJ2DZsuzO99JLcMstwfbvf+9rGlEYPBhefz1YNsU5uOEG/9jSdJZX//RT/+v+/POhvNzva9XKP0DprLPylu2UWreGf/4zGNP6yiu5qcWFx25fcEHdZrGLxJSCRnX69fOLB4L/xXrHHZmfa9UquPDCYHvECLjyyqyyl7Xevf3aT+El4f/3PxgwwK/rNHkyfPFF8NrOnX4ZjdGj/ZIejz4avNarl2/CC58rCkcdFTxMC/xaXDNn1nx8bRYvTp4F/v38zQMRKSYKGjUZHRr3c9ddvtmjrpzzzTWffea3O3TwM6Prw5MBO3b0TVW//GWwvPrOnb6J6bjjfL/OQQf51WjbtfPP67j9dv+M7IRvfcsv1d4vN6OwsnbttX7JcvArFZ9/fnJ+6+K22/znB365+QMPzE0eRYqcgkZNTjstaCdfu9YvtFdXv/998vIT99/vA0d9UVIC110Hs2f7L8awzZv9Uwznzk1+GBT4ZdknTfLLmqfxTImCadwY/vWvYJmRhQt9U6OrdnpPzT76CP72t2BbtQyRXRQ0alJSAr/6VbB9yy1169t46SW45ppge/To7NeVypeDDoLnn/d9GldcEQTLsPbt/S/3F1/0a0odd1zh85mOnj19jSjhwQfr3rx47bV+tBv4PqD6+rmJREBBI5XzzvMPCAL/AKAf/Si9X61vvw0jRwZLnB59tK911Hf9+vnRYosX+76YOXOCJcFXrPBLrw8fXj+a11K56CL/zI6EH/7Qd9SnY9o0X1tJ+N3v6n95RQpIQSOVRo2Sv+wff9w/4zqVBQv8L9ONG/12586+n6C0NH/5zIf27f0Ir/79Ye+9i2u1PTNfuxhYuaRMRYUP4lOnpn7fhg3+h0Lih8Gpp6b1HHCRhqSIvgkicvzx/rnWCT/8YfJjP8Nef92P4lmxwm/vsYf/hdupuiW4JK/Kynx/0r77+u1Nm/wS8TU98W/7dl87WbzYb7duDbfeWpi8ihQRBY10/OlPfhQR+GaqESP8F0pidveqVf7hQsce6zvNAZo390M2E++Twuvc2ffBtGvnt8vL4eSTfef/1q3BcWvX+oEP4WeP33MP7JfdmmAicaSgkY5mzfwkv8TieRUVvn+jY0ff8dqxo59MlpiF3KGDn2B29NHR5Vm8gw7yn0XXyrUvKyr8AIcuXXxT1Dnn+NpIuM9j9Gi/X0R2o6CRrl69/AS48EzujRvhgw+SO8ePPRbefBMG1fgUWym0Pn38ZzdkSLBvzRp45BHf37RpU7D/2mvVLCWSgoJGXXTp4jtT//SnYBkO8J3ERx3ll7KYPLn6IasSrX328TPX77jDN1tV1bevH+gwZoxGS4mkEMnjXota8+a+M/wHP/AzvTds8M1Re+4Zdc6kNo0a+cl+l13mZ7K/847/PHv08LUQBQuRWiloZMrMD0Xde++ocyJ11aiRn7Q3eHDUOREpOmqeEhGRtCloiIhI2hQ0REQkbQoaIiKSNgUNERFJm4KGiIikTUFDRETSZq6uTzWr58xsFfBRikPaA6sLlJ1CiWOZIJ7lUpmKQxzLBKnL1c05t1dtJ4hd0KiNmU13zsVqYag4lgniWS6VqTjEsUyQm3KpeUpERNKmoCEiImlriEGjlue1FqU4lgniWS6VqTjEsUyQg3I1uD4NERHJXEOsaYiISIYaZNAws1+Z2Wwzm2VmL5lZNU/lKS5m9nszW1BZrifNbI+o85QtMzvbzOaZ2U4zK+qRLGZ2opktNLNFZnZ11PnJBTO7z8xWmtncqPOSK2bW1cwmm9m7lffeD6LOU7bMrMzMppnZO5Vl+mVW52uIzVNm1to593llejTQ2zl3ecTZyoqZDQdeds5VmNlvAZxzP404W1kxs4OBncBdwP8556ZHnKWMmFkJ8B5wArAMeAs4zzk3P9KMZcnMvgyUAw865w6p7fhiYGZ7A3s752aaWStgBnBGMX9WZmZAC+dcuZk1AV4DfuCcezOT8zXImkYiYFRqARR95HTOveScq6jcfBPoEmV+csE5965zbmHU+ciBwcAi59xi59wXwCPA6RHnKWvOuf8Ca6PORy455z51zs2sTG8E3gX2iTZX2XFeeeVmk8q/jL/zGmTQADCzG81sKfAN4Lqo85NjlwDPR50J2WUfYGloexlF/kXUEJhZd6A/MDXanGTPzErMbBawEpjgnMu4TLENGmY20czmVvN3OoBz7ufOua7AQ8CV0eY2PbWVqfKYnwMV+HLVe+mUKQaqe/h40ddu48zMWgJPAD+s0jJRlJxzO5xzh+FbIAabWcbNibF9Rrhzbliahz4MPAdcn8fs5ERtZTKzC4FTgONdkXRW1eFzKmbLgK6h7S7A8ojyIrWobPd/AnjIOfefqPOTS8659WY2BTgRyGgAQ2xrGqmYWa/Q5mnAgqjykitmdiLwU+A059zmqPMjSd4CeplZDzMrBc4FxkWcJ6lGZafxvcC7zrlbos5PLpjZXonRlGbWDBhGFt95DXX01BPAgfiROR8BlzvnPok2V9kxs0VAU2BN5a43YzAi7EzgdmAvYD0wyzn31WhzlRkzGwHcCpQA9znnbow4S1kzs7HAUPzKqSuA651z90aaqSyZ2dHAq8Ac/PcDwM+cc+Ojy1V2zKwf8AD+3msEPOacG5Px+Rpi0BARkcw0yOYpERHJjIKGiIikTUFDRETSpqAhIiJpU9AQEZG0KWiIiEjaFDRERCRtChoieWZmh1c+56TMzFpUPtMgFkuJS8OjyX0iBWBmvwbKgGbAMufcTRFnSSQjChoiBVC55tRbwFbgS865HRFnSSQjap4SKYw9gZZAK3yNQ6QoqaYhUgBmNg7/xL4e+MeJFsUzXESqiu3zNETqCzO7AKhwzj1c+bzwN8zsOOfcy1HnTaSuVNMQEZG0qU9DRETSpqAhIiJpU9AQEZG0KWiIiEjaFDRERCRtChoiIpI2BQ0REUmbgoaIiKTt/wPaEijYzOPD9QAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x15202d8a20>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from sympy import *\n",
    "import mpmath\n",
    "# fourier 5回級数展開\n",
    "l = [-mpmath.pi, mpmath.pi]\n",
    "f0 = lambda x:sign(sin(x))\n",
    "cs = mpmath.fourier(f0, l, 5)\n",
    "f1 = lambda x: mpmath.fourierval(cs, l, x)\n",
    "mpmath.plot([f0, f1], xlim=l)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "∫_0^{¥pi/3} sin x dx = 0.5 ± 5.55112e-15\n"
     ]
    }
   ],
   "source": [
    "y, abserr = quad(sin, 0, pi/3)\n",
    "print(\"∫_0^{¥pi/3} sin x dx = %g ± %g\" % (y, abserr))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "def fourier(fun, n_max):\n",
    "    a = []\n",
    "    b = []\n",
    "    for n in range(0, 1+n_max):\n",
    "        res, err = quad(lambda x:fun(x)*cos(n*x), -pi, pi)\n",
    "        a.append(res/pi)\n",
    "        res, err = quad(lambda x:fun(x)*sin(n*x), pi, pi)\n",
    "        b.append(res/pi)\n",
    "        \n",
    "    def fn(x):\n",
    "        sum = a[0] / 2\n",
    "        for n in range(1, 1+n_max):\n",
    "            sum += a[n]*cos(n*x) + b[n]*sin(n*x)\n",
    "        return sum\n",
    "    return fn\n",
    "\n",
    "    def f(x):\n",
    "        x = (x+pi) % (pi * 2) - pi\n",
    "        if x>= 0:\n",
    "            return 1\n",
    "        else:\n",
    "            return 0\n",
    "    x_min = -pi - 0.5\n",
    "    x_max = pi*2 + 0.5\n",
    "    y_min = -0.25\n",
    "    y_max = 1.25\n",
    "    axis([x_min, x_max, y_min, y_max])\n",
    "    \n",
    "    f_fn = fourier(f, 10)\n",
    "    xs = linspace(x_min, x_max, 256)\n",
    "    plot(xs, amap(f, xs), 'b:', lw=1)\n",
    "    plot(xs, amap(f_fn, xs), 'r-', lw=1)\n",
    "    show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "'int' object is not callable",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-20-9921ded9bfb7>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mfourier\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0;36m2\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[0;32m<ipython-input-19-7971275df3f5>\u001b[0m in \u001b[0;36mfourier\u001b[0;34m(fun, n_max)\u001b[0m\n\u001b[1;32m      3\u001b[0m     \u001b[0mb\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m     \u001b[0;32mfor\u001b[0m \u001b[0mn\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mrange\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0mn_max\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 5\u001b[0;31m         \u001b[0mres\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0merr\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mquad\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;32mlambda\u001b[0m \u001b[0mx\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0mfun\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mcos\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mn\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m-\u001b[0m\u001b[0mpi\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpi\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      6\u001b[0m         \u001b[0ma\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mappend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mres\u001b[0m\u001b[0;34m/\u001b[0m\u001b[0mpi\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      7\u001b[0m         \u001b[0mres\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0merr\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mquad\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;32mlambda\u001b[0m \u001b[0mx\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0mfun\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0msin\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mn\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpi\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpi\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/.pyenv/versions/anaconda3-5.1.0/lib/python3.6/site-packages/scipy/integrate/quadpack.py\u001b[0m in \u001b[0;36mquad\u001b[0;34m(func, a, b, args, full_output, epsabs, epsrel, limit, points, weight, wvar, wopts, maxp1, limlst)\u001b[0m\n\u001b[1;32m    321\u001b[0m     \u001b[0;32mif\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mweight\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    322\u001b[0m         retval = _quad(func, a, b, args, full_output, epsabs, epsrel, limit,\n\u001b[0;32m--> 323\u001b[0;31m                        points)\n\u001b[0m\u001b[1;32m    324\u001b[0m     \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    325\u001b[0m         retval = _quad_weight(func, a, b, args, full_output, epsabs, epsrel,\n",
      "\u001b[0;32m~/.pyenv/versions/anaconda3-5.1.0/lib/python3.6/site-packages/scipy/integrate/quadpack.py\u001b[0m in \u001b[0;36m_quad\u001b[0;34m(func, a, b, args, full_output, epsabs, epsrel, limit, points)\u001b[0m\n\u001b[1;32m    386\u001b[0m     \u001b[0;32mif\u001b[0m \u001b[0mpoints\u001b[0m \u001b[0;32mis\u001b[0m \u001b[0;32mNone\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    387\u001b[0m         \u001b[0;32mif\u001b[0m \u001b[0minfbounds\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;36m0\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 388\u001b[0;31m             \u001b[0;32mreturn\u001b[0m \u001b[0m_quadpack\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_qagse\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfunc\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0ma\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mb\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mfull_output\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mepsabs\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mepsrel\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mlimit\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    389\u001b[0m         \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    390\u001b[0m             \u001b[0;32mreturn\u001b[0m \u001b[0m_quadpack\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0m_qagie\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mfunc\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mbound\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0minfbounds\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mfull_output\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mepsabs\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mepsrel\u001b[0m\u001b[0;34m,\u001b[0m\u001b[0mlimit\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m<ipython-input-19-7971275df3f5>\u001b[0m in \u001b[0;36m<lambda>\u001b[0;34m(x)\u001b[0m\n\u001b[1;32m      3\u001b[0m     \u001b[0mb\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m     \u001b[0;32mfor\u001b[0m \u001b[0mn\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mrange\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0mn_max\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 5\u001b[0;31m         \u001b[0mres\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0merr\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mquad\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;32mlambda\u001b[0m \u001b[0mx\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0mfun\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mcos\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mn\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m-\u001b[0m\u001b[0mpi\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpi\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      6\u001b[0m         \u001b[0ma\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mappend\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mres\u001b[0m\u001b[0;34m/\u001b[0m\u001b[0mpi\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      7\u001b[0m         \u001b[0mres\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0merr\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mquad\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;32mlambda\u001b[0m \u001b[0mx\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0mfun\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0msin\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mn\u001b[0m\u001b[0;34m*\u001b[0m\u001b[0mx\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpi\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mpi\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mTypeError\u001b[0m: 'int' object is not callable"
     ]
    }
   ],
   "source": [
    "fourier(1,2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/noriakioshita/.pyenv/versions/anaconda3-5.1.0/lib/python3.6/site-packages/ipykernel_launcher.py:11: RuntimeWarning: invalid value encountered in sqrt\n",
      "  # This is added back by InteractiveShellApp.init_path()\n",
      "/Users/noriakioshita/.pyenv/versions/anaconda3-5.1.0/lib/python3.6/site-packages/ipykernel_launcher.py:12: RuntimeWarning: invalid value encountered in sqrt\n",
      "  if sys.path[0] == '':\n",
      "/Users/noriakioshita/.pyenv/versions/anaconda3-5.1.0/lib/python3.6/site-packages/ipykernel_launcher.py:13: RuntimeWarning: invalid value encountered in sqrt\n",
      "  del sys.path[0]\n",
      "/Users/noriakioshita/.pyenv/versions/anaconda3-5.1.0/lib/python3.6/site-packages/ipykernel_launcher.py:14: RuntimeWarning: invalid value encountered in sqrt\n",
      "  \n",
      "/Users/noriakioshita/.pyenv/versions/anaconda3-5.1.0/lib/python3.6/site-packages/ipykernel_launcher.py:15: RuntimeWarning: invalid value encountered in sqrt\n",
      "  from ipykernel import kernelapp as app\n",
      "/Users/noriakioshita/.pyenv/versions/anaconda3-5.1.0/lib/python3.6/site-packages/ipykernel_launcher.py:16: RuntimeWarning: invalid value encountered in sqrt\n",
      "  app.launch_new_instance()\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAIABJREFUeJzt3XdUFFcbBvDnsjRBEZUiigrYC1bAglixK0bFHnuJxhJr7BqN5TOWaIwl9t5bNGqsWLCCoih2AQWlg4D0Zd/vD8SYxChld2d3eX/ncE6E3bkPZHmYvXNnRhARGGOM6Q49qQMwxhhTLi52xhjTMVzsjDGmY7jYGWNMx3CxM8aYjuFiZ4wxHcPFzhhjOoaLnTHGdAwXO2OM6Rh9KQa1sLAgOzs7KYZmjDGtdfv27WgisvzS4yQpdjs7O/j6+koxNGOMaS0hxMucPI6nYhhjTMdwsTPGmI7hYmeMMR3Dxc4YYzqGi50xxnQMFztjjOkYLnbGGNMxXOyMMaZjuNgZY0zHcLEzxpiO4WJnjDEdw8XOGGM6houdMcZ0DBc7Y4zpGC52xhjTMVzsjDGmY5RW7EIImRDCTwjxh7K2yRhjLPeUucf+HYBHStweY4yxPFBKsQshbAF0ALBRGdtjjDGWd8raY18B4HsACiVtjzHGWB7lu9iFEB0BRBLR7S88brgQwlcI4RsVFZXfYRljjP0HZeyxuwLwEEIEA9gLoIUQYuc/H0RE64nIiYicLC0tlTAsY4yxT8l3sRPRNCKyJSI7AL0AXCCir/OdjDHGWJ7wOnbGGNMx+srcGBFdBHBRmdtkjDGWO7zHzhhjOoaLnTHGdAwXO2OM6RgudsYY0zFc7IwxpmO42BljTMdwsTPGmI7hYmeMMR3Dxc4YYzqGi50xxnQMFztjjOkYLnbGGNMxXOyMMaZjuNgZY0zHcLEzxpiO4WJnjDEdw8XOdMK9iwG4dzFA6hjISM/A0V9P4W1UvNRRWAHGxc50wpaZuzHnq8V4eOOppDlun/HH6rGb8fjmc0lzsIKNi53phOg3cUhJSsP5nZclzXFhzxUYmxrhvvcjSXOwgo2LnWm9THkmokJiYFvJBpf2X4M8Qy5JjpR3Kbj+uy+sylrg1IZzkmRgDOBiZzogMiQaikwFareogfjoRNw+6y9Jjmu/+yI1OQ01m1RDYlwS3r1NkiQHY1zsTOuFBUYCABp2ckaRYqbw2ustSY7zu6/AqqwF6rRwfJ8rQpIcjHGxFyAJMYnwPXMPexYdwbwey3RmjzL8fYGWqVwKbt0a4NpRH6Qkpao1w9uoeNw+cw/Ne7miVIWSWbmCItWaQVVSk9OwYsRv2L/kd/hduK8zrxtdpi91AKYa6anpeHYnCI9uPMUTn+d4fOv534qmpL0Vol/HorC5qYQplSMsMAL6BjJY2BZHi75uOLnxPK4f80WL3o3VluHS/utQZCrQsq8brMpaAADevND+PXYiQmezfjAubIzkhJQPny9d0QZVXCqgsnMFVG1QEeVr28HA0EDCpOxjXOw6IjY8Dg+8HyPg6hM8vP4Ez/2CIM/IBABYlbVAZZcK6DC8FSo5lUfFuvYoUqywxImVJywoAlblLCGTyeDoVhWWtiVwYfcVtRb7hd1XYFejDOwdywEAihQvjPAg7S/2jLQMKBSEXlO6oP2wlnh2JwhPfV/gqe9z3PV6gPO7rgAADIwMULGeA6o3rITqrlVQ3bUyzC2LSpy+4OJi11KRr6Jw79JD+F8MgP+VR3jzPBwAYGhsgMouFdBtfEdUbVAJVRtURPGSxSROq1phgZGwcbAGAOjp6aF5L1ccWnEC8dEJKGphpvrxgyLw8PpTDFnY58PnbBys8UYH5thTk9MAAEYmhihqYQan1rXg1LrWh69Hhcbg0Y2neHTjGR5ef4Kjq07hwLLjAICyVUvD0a0aajbN+rAoVVyS76Eg4mLXEvHRCfA7fx9+5+/jrteDD2/zixQzRQ23qug4vBWqN66CinXtC9xb4rDACFSq5/Dh3y36umH/0mO4fOA6Oo1so/LxL+zOOljb/KN3CDYOVnh2J0jlY6taalJWsRubGH3y65a2JWDp2RBNPBsCyJoCfOr7Ag+8H8P/yiN47fXGifVnAWQdA6ndwhF1WjqidvPqOvWuUdNwsWsoeYYcD68/hc+fd3H7zF089wsGEcHErBBqNauOzqPaoVbz6rB3LAs9vYJ7DDwpPgmJse8+7LEDgEPNcihXzRYX9nirvNiJCF57vFGjcRVYl7P88Hkbe2t4H76FzMxMyGQylWZQpbQPe+yfLvZ/MjQ2RI3GVVGjcVX0mtoFmZmZeO4XDP+LAfDzeoCz2y/i+NrT0NMTqORcAU6ta8G5bW1Udqmg1T8nTcPFrkHiIuNx6+Qd3Dx5B7fP3ENyQgr0ZHqo1rAS+v/QA3Vb1URlp/KQ6ef/FyAhNhEv7gajiksFFCpcSAnppZG91PHjYhdCoEUfN2yZuQfhwZEoaWelsvFf3AvGy4ehGLtm2N8+b+NgjUx5JqJDY/9W+Nrmwx67ac6K/Z9kMhkqO5VHZafy6D7JA/IMOR7ffAbfM/dw+6w/di84hJ0/HkSR4oXh1KYW6revB+d2tWFWvIgyv40Ch4tdYqFP38D7yC1cP+aDRzeegYhQolQxNO3eCM7t6qBuyxowLar8lSsBV59gdufFWHVjIaq4VFT69tUle634x8UOAM17u2LLzD3w2nMVvad1Udn4F3ZdgUxfhqbdG/7t8yXf5wkLjNDqYn8bmXUxs6KWyjlWoW+g/2GPfuC8XkiITcSds/649acffE7dhdeeq9CT6aFG4ypo5OGMRl85w8be+ssbZn/Dxa5mRITgB69w6cB1eB++iZcPQwEAFevao9/s7mjQqR4q1LGHEEKlOYrbZB1QjQ17q9JxVO2vYv/7XrmNvTWqNaqMC3uuqKzYFQoFLu67Bue2tWFW4u97mNl53ryIQO3mNVQyvjrEvH99lLBRzQF4s+JF0KynK5r1dIVCocATnxe4cdwX1475YN3EbVg3cRvK17aDW9cGaNK9AcpULq2SHLqGi11NXj1+jYt7r+LivqsIefIGenoCNdyq4tsVg+DaxQVWZSzUmsfSNmuFQuSraLWOq2xhgREoUrzwJ9/VtOzjhlWjNyLQ/yUcapZT+tj3Lz9CVGgMhv3U719fsypjAZm+TOuXPEa9f32UKKX6lVV6enqoWr8iqtaviEHzeyMsMALeR27B+/ANbJ29F1tn74W9Y1k06+mK5r1c//Uujf2Fi12FYsLi4LXHG+d3XcFzvyAIIVCzaTV0+a4DGndxQTFrc8myFbM2h7GJkdaf9h4WFPmfv+BNezTEmnFbMLfbEljYllD62PKMTIz8eSAadKr3r6/J9GWwKmuBMC0/+zQsKAIWpYvD0NhQ7WPbOFij+8RO6D6xE6Jfx+DKoZu4uP8atszcgy0z96BK/Ypo2dcNzXo24jXz/8DFrmTpqem49rsPzmy7iNtn7kGhIFR2Lo+Ryweiac9GKntLm1tCCJSqWBIhT99IHSVfwgIjULGu/Se/VtTCDDuD12BR35UqGVvfQIau33X4z6/bOFh9uNyBtgp9+ga2lWykjgGL0iXQZWx7dBnbHpGvouC19xou7L6C1WM3Y92EbajfoS5aD2iG+h3qQt+Aa41/AkoS9OAVTm44h/M7LyMxLgmWZUqg55Sv0Kp/U42dFyxXzRYBV59IHSPPMjMzEfkyCm5d6//nYyxKFccyr7lqTPUXG3treB+5KcnYyqBQKPAyIBTu/ZpIHeVvrMpaouf3ndHz+84Iuv8SZ7dfwvldV3Dtdx+YWxVF6/5N0W5oS9hWKiV1VMnku9iFEGUAbAdQEoACwHoiUs0ukoZJT8vA5QPXcXzdGTy89gQGhvpw7eKCtoNboHaLGhq/Ltehph289lxFQmyiVi4viw6NhTwjU2PnWks6WCM+OhHJiSkwKaJ9S0ojgqOQnJiikuMTymLvWA7Dl/THkEV94fPnXfy5+TwOrTiB/UuPoXbz6ug4og1cv3IucHvxyvhu5QAmEtEdIUQRALeFEGeJ6KEStq2Rol/H4Nia0zi18TzeRiWgdEUbDF/SH60HNFXLKezKkn225lPfwL+dJq4t/mupo6aw+WjJY/ladtKGyYMnPlm396vkVF7iJF8m05ehQcd6aNCxHmLD4/DnZi+c3HAO83suR4lSxdDxm9bo8E0rFLMqGHPx+S52IgoDEPb+vxOFEI8AlAagc8X+9PYLHFx+HJcP3AApFKjfsR46j2qLOi0dtfLsz8rO5SGEwKPrT7W72MtrarFnLXnU1mIPuPYERoUMYe9YVuoouVK8ZDH0md4VPad0hs+pu/h99Slsm7MPuxccQos+bug2oSPsa2jX95RbSn1/IoSwA1AHwL8mFoUQwwEMB4CyZbXnh0pE8D1zD/t/Ooq7XgEwKVIIX41ui85j2mn9iROmRU3hUKsc7l0KQD90lzpOroUFRkBPpqf2paI5VdI+u9i1c2WM/+WHqNqwktZOY8hkf+3Fhzx5jSMrT+LMtos4vdULzu3qoOf3nVGzSTWVnzMiBaXtZgohCgM4BGAcESX88+tEtJ6InIjIydJS88/EUygUuHL4JkY5T8H0dgsQ+jQMw3/qh92v1mLE8oFaX+rZ6rRwxMNrT9R+YwplCA+OhFVZC6VcYkEVzIoXQWFzU61cUhoX8RaB915+uBuUtitTuTTGrhmG3a/WYeC8Xnjm+wKTmv+AcW6zcOuUH4hI6ohKpZRiF0IYIKvUdxHRYWVsUyoKhQJXDt3AiDqTMc9zKZISUjBhwwhsf/Eruk/yUMnp/VJyblsbGely3L3wQOoouRb2IkJj59ez2ThYaeVJSj5/3gWQ9frQJWYliqDvzG7YGbwGo1cNQXRoDGZ0WIgxDabh5sk7OlPw+S52kfU+ZhOAR0S0PP+RpEFEuHXKD6Ocp2Je92WQp8sxdcdYbH64Au2GtNTZS+E6NqkKE7NCuHrkltRRci0sMAI29qq7wJcylHSw1so99qtHb8GidHFUqPPpcwS0nVEhI3Qe1RZbn/6C8etHID4qATM7LsK4xjNx72KA1PHyTRl77K4A+gFoIYS4+/6jvRK2qzaPbz3DpBY/YEaHhXj3NgnfbxuNDQ+Wo2VfN419m68sBoYGaOjhhKtHbyE9LUPqODmWnJiC+OhEzd9jt7dGRHAUFAqF1FFyLCk+CT5/3oVbtwY6Of/8MQNDA7Qf2hJbnvyC79YOR+SraExq8QOmt1+AQP+XUsfLs3wXOxF5E5EgoppEVPv9x0llhFO1yFdRWNh3BcY0mI5XD0MxetUQbH60Aq36NdX4NejK1KK3G969TcLNE3ekjpJjmr7UMZuNgzUy0uWIfh0rdZQcu7T/OjLSMtCij/puLSg1fQN9dPymFbY+/QXDf+qHRzeeYUSdyVg+bB1iw+Okjpdr2rdGTwnSU9Ox88eDGFx1HK4euYU+07ti2/Nf0XlUW52dcvmceq1qwqJ0cZzaeE7qKDmWXewlNb7Y/1ryqC1ObTqPctVsUdm5gtRR1M6okBG6T/LAtuer0PW79ji7/SIGVf4OB5YdhzxDLnW8HCtwxe5z+i6GOU7Atjn74NKhLjY/WolB83tr5ZmByiLTl6H9UHf4/HkXoVpy7Zi/brCh4XPsWrbk8fGtZ3h86zk6ftNa56dhPseseBGMWD4QG+4vRw23Klg/eTtG1JmM+1ceSR0tRwpMscdFxmNh3xWY3m4BZPoyLD47G7P3T9TqmyAoU8cRrWBgqI+Dy/+QOkqOhAVGwLSoicZfCsG6nCX09ITWXAzswLLjMDErhFYDmkodRSPYViqFBX9Mx9yj3yM1KQ0Tms7G8mHr8O5tktTRPkvni52I4LX3KoZWH48rB2+g/5weWHd3Keq21I31ucpSzNocrQc2x5mtXogM0fxrtIcHaf5SRyBr7tayjAXCtGDJ48uHIbhy8AY8vm0LUzMTqeNolEYeztjwYDl6TPLA6a1eGFJ9PG78cVvqWP9Jp4s9ISYR83v9jIV9VqBUeWus81uCfnO6w9Co4M2j50T2nYZ2zD0gcZIvCw+K1PhpmGxZa9k1fypmy6y9KFTYGJ4TOkodRSMVMjXGsJ/6YdWNhShqUQSzPP6HpYPXICkhWepo/6Kzxe534T6G15qIa0dvYfCCPljhPR/lqpWROpZGsy5nCY9v2+DMVi88vxskdZz/pFAosopdS87+LWlvjTcvNHuP/d7FAFw9cgs9JnfWqgvZSaFSvfJY7fM/9J7WBWe3X8SIOpPx6OYzqWP9jc4Ve2ZmJrbN2YcprX6ESZFCWHVjEXpP66Lz69GVpe8sT5iVKIKVIzcgMzNT6jifFPMmDhnpco1fEZPNxsEabyPjNfayDelpGfhl1AaUtLOE50TeW88JA0MDDF7QB8suzQMpCOPdZuHA0mMac+aqThV7fHQCprdbgJ0/HoR7/yZY7btYZ8+cU5UixQpjxPKBeHzzGQ6v0MzTEbRlDXu27JyaegB157wDePXoNcasHgajQkZSx9EqNVyrYJ3fEjT0cML673dgbrclSIqX/sCqzhT7i3vBGOU8FfevPMaEDSPw/ZbRKGRqLHUsrdSiT2M06uyMLTN247mf5k3JZE9raNMcOwCNvP/pvUsB2Lf4KNoMbA6XdnWkjqOVCpubYvaBiRi5fCCuH7+NMQ1nSL5sWCeK/fpxX4xrPBOZ8kz8fHke2g1pKXUkrSaEwIQNI2BuVRRzPZciISZR6kh/Ex4UAT09oTVLVf/aY9esYo8KjcGCXj+jVIWS+HblIKnjaDUhBLqO64Cfzs1GQnQCxjSYjrte0l1YT+uL/fja0/ihy08oW9UWq24uKpBny6lCUQszzDowETGvY/FD1yVIT02XOtIH4UGRsCxjoTXXCTcrUQQmRQrhzYtwqaN8kJSQjNmdFyM1KQ1zDk0u0CfoKVOtptWx6uYilChVDNPazseF3VckyaG1xU5E2DHvAH4ZtREu7etiqdcPsChVXOpYOqVq/YqYvHU07l95hPm9fkZGumZcJCwsMEJrpmGArL25kg5WGrOWPTU5DXO++gmB/i8xc98E2FXn1WLKZGNvjRXe81HdtQoWff0Ljvyi/mNVWlnsRISNU3Zi+w/70WpAU/xweDLPp6tI816uGPPrUFw/5ou53ZYiNTlN6kgIC4xASS1Z6pjNxsFaI9ayJ8UnYUaHhfC/9BBTto3meXUVKWxuioUnp8O1iwvWjNuCPYuOqHV8rSv27FLfv/QYOo1sg0mbvuWljCrm8W0bjF0zDLdO+mFyy7mICZPuancpSamIi4jXmhUx2Wzss4pdysv3hgVFYHyT2Qi4+gRTd45Fiz5ukmUpCAyNDTFr3wS07OuGzTN2q7Xcta7Yd/548EOpj/l1iFbeRFobdRrRGrMPTkTw/Vf4tt73uHPOX5IcpCAMX9If9bTs5ts2DtZIT81AbPhbSca/9rsPRjlNQVRIDBaenI4WvQvOJXmlJNOXYfLWUR/K/eiqU2oZV6ta8Y/fzn6Yfhm9anCBvvqcFBp3qY9fbiyEqbkpprT+EStHrkdi3Du1Zgi6/wptBjZDZafyah03v7KPCdw5e0+t48ZFxmPxgFWY0+UnWNtZ4ddbi1DXvaZaMxR0MpkMk7eMgutXzlj93WZc3HdV9YMSkdo/6tWrR7l169Qdai3rTtM7LCB5hjzXz2fKk5KUSmvHb6HWsu7U1WIQHVh2jFKSUlU63umtXjS6/lRyF560c/5BlY2lKu/evqODPx8nd+FJ45vOogt7vCktJU1148Un0c75B8mjaD9qa9iTtszcQ+lp6Sobj31ZWkoajW8yi9oZ9aIHVx/naRsAfCkHHStIglNgnZycyNfXN8ePD336BqNcpqKkvRVWXPkRhQrz0ixN8OJeMNZP3o475+7DrEQRtB/aEm0GNYdtpVL53nZaShrunLuPK4du4OqRW0hOTEGZyqXgMaot2gxqrpUHy5Pik3By4wUcW/MnwoMiUaSYKRp3bQA3zwao1ay6Ui5OF/TgFU5vvoA/t3ghKT4ZDT2cMPR/X6NsldJK+A5YfiXEJGJMw+lITkjBGt/FsLQtkavnCyFuE5HTFx+n6cWempyGMQ2mITbsLdb4Ltaak1IKkgfej3Bw+XFcP+YLhYLgULMcnNvWhmOTaqhUzwHFrM0/+3wiQlzEWwT6v8Ljm89w/8pD3L/yGBlpGTAtaoLGXeqj1YCmqNmkmk5MvykUCvidv4+z2y/h2u8+SHmXCmMTIzg2rQbHxlVR2aUCHGqWRVELs89+v0SEmLA4PPV9gfuXH+HmyTsIefwa+gYyNO5aH54TPbRuyqogePkwBGMaTIe9Y1ksuzg3V+dj6EyxrxixHifWn8XCUzPg3Ka2ipOx/Ih+E4tL+67h6u+38PDaU2TKsy4iVqR4YViXs0RRSzMUKmwMmb4e5OlyJCWk4G1kPCJfRiM5MeXDduwdy6J28xpwaV8XtZpV0+nbFaanpsPv/H3cPOmHu14PEPL49YevFTY3hVVZCxS1NIOJWSHoG8iQKVcgOTHr5xYRHIWk+KxLxhoY6sOxSVU06uyCpj0awtyyqFTfEsuBi/uuYkHvFeg1tQuGLOyT4+fpRLHfOuWHGR0WovvEThi+pL8akjFlSUlKxVOfF3hxNxghT14jKjQG8VEJSHmXCkWmAvqG+ihUpBDMLc1gVcYCpSqUhL1jWVSsaw/ToqZSx5dMQmwint8JQvCDEIQ+C0P066yfW3JiCjIzMiHTl8G4sDHMLc1gWcYCZauURoU6dqhYz4Ev4KVllg1dizNbvbD88o+o3qhyjp6j9cWenJiCodXHw7SoCVb7LuabYzDGdEpyYgqG15wIw0KGWOe3JEcdl9Ni19jljjvnHUBUaAzGr/+GS50xpnNMihTC2DXDEPL4NQ4p+V7DGlnsb16E48gvJ9FmYHNUa5iztyiMMaZtXNrVgetXzti98BDiIpR38ppGFvv2ufsh05dh4PxeUkdhjDGVGrq4H9JTM7B7wWGlbVPjiv318zB47faGx7dt+GqNjDGdZ1vRBq0HNMOJDecQG66c6zBpXLEfXnECMn0ZPCd2kjoKY4ypRc8pX0GeLsexNaeVsj2NKvaUpFSc23EZTXs2QvGSxaSOwxhjamFb0QYu7evg1MbzkGfI8709jSr27FPH2w91lzoKY4ypVfuh7ogNf4vbZ/J/oTiNKvZLB67BskwJVHfllTCMsYLFuV1tFDY3xaWD1/O9LY0p9vTUdPidu49GHs58jXXGWIFjYGiA+h3qwuekX75vyKIxDfroxjOkpaRr3Q0UGGNMWeq1qoW3UQkIfhCSr+1oTLE/vP4UAHgahjFWYGX3X3Yf5pXGFHugfzBK2lnCrHgRqaMwxpgkbBysYVrUBIH3gvO1HaUUuxCirRDiiRDiuRBial628fpZGMrwzQAYYwWYEAJlKpfC6+dh+dpOvotdCCEDsBpAOwDVAPQWQlTL7XaiQmNhVcYiv3EYY0yrWZa1QFRITL62oYw9dhcAz4kokIjSAewF0Dk3GyAiJMQkokgJnoZhjBVsZsWLIDE2fzeJV0axlwbw8SHc0Pefy7FMeSYUmQoYm/CNAhhjBZuxiSHSUtLztQ1lFPunbsr4r7t3CCGGCyF8hRC+UVFRShiWMcbYpyij2EMBlPno37YA3vzzQUS0noiciMjJ0vLvN6SW6csg05chNSlVCXEYY0x7pSanwyifsxfKKHYfABWFEPZCCEMAvQAcy80GhBAwK1EY8dGJSojDGGPaKz46AWYlCudrG/r5DUFEciHEaACnAcgAbCaigNxux7KMBSJDovMbhzHGtFrkq2hY2JbI1zaUso6diE4SUSUiKk9EC/KyDdtKNnj1KFQZcRhjTCspFAqEPnkD24o2+dqOxpx56lDTDlEhMYiPTpA6CmOMSeLN83AkJ6agfG37fG1HY4q9eqNKAID7Vx5JnIQxxqTxwPsxAKBaw0r52o7GFHtllwowNjVSykXmGWNMG/meuYviJc1RrpptvrajMcVuYGgApza1ce2YLzIzM6WOwxhjapWemg6fU3fh0r4uhPjU6UE5pzHFDgBNuzdEbFgc7l18KHUUxhhTq+vHbyM5MQXNejbK97Y0qtgbejihsLkpTm08J3UUxhhTq5Mbz8HStgRqt6iR721pVLEbFTJCm0HNceXQTV7TzhgrMIIDQnDnrD86DG8FmUyW7+1pVLEDQJex7QEAB5bk6uRVxhjTWnsXH4GxiRE6jmillO1pXLFbl7NEq/5NcWL9WUS85IuFMcZ0W3BACC7s8obHt21Q1MJMKdvUuGIHgH5zukPoCWyavkvqKIwxpjJEhN8mbYOJWSH0nPKV0rarkcVuVcYC3Sd5wGvPVdz1eiB1HMYYUwnvI7fge/oe+s/pATMl3mhII4sdAHpP64JS5a3x8ze/ITU5Teo4jDGmVIlx7/DrmE1wqFUOnUe3Veq2NbbYjQoZYcKGkXjzPBzrJ++QOg5jjCkNEeGXURsRH5WASZu+hUw//ythPqaxxQ4AtZpVR7fxHXF87Wl4H7kpdRzGGFOK01u8cHHvVfSb0x0V6zooffsaXewAMHhhH1RyKo8lg1Yj9Om/bszEGGNa5dmdQKwavRG1W9RAr6nKO2D6MY0vdkMjA8w+MBEGhvqY3XkxEuPyd/duxhiTSmx4HOZ0+QlFLc0wffc4pZyM9CkaX+xA1tr22QcnISwwAnO7LUV6WobUkRhjLFdS3qVglsdiJMa8w7yjU1DMqqjKxtKKYgeAmk2qYdLmUbh3MQD/6/cLMuV8BUjGmHZIT8vAD92W4rlfEGbsHY8KdfJ3I40v0ZpiB4CWfd0wYtkAXDl4A8uGruXL+zLGNF5GegYW9PoZd876Y8KGEWjQsZ7Kx8z3zazVrdv4jkhNSsPW2XtBRCpZKsQYY8qQnpZV6td+98HoVUPQZmBztYyrdcUOAH1ndgMEsHXWXqQlp2Hqzu9gaGQgdSzGGPsg5V0Kfui2FHfO+mP0qiHoPEq5JyF9jlYWOwD0ndENhUyNsXbCViTELMAPhyejsLmp1LEYYwxxkfGY1WkRnt0OxMRN36L5ODhPAAAZmUlEQVTtIPXsqWfTqjn2f+o6rgOmbB+DgKuPMa7xTIQFRUgdSatkZmbyQWj2n4gImfJMKBQKqaNolZcPQzC24XQEPwjBnMOT1V7qgBbvsWdz/7oJLEoXx9xuSzGm/jTM2j8RtZpVlzqWRklNToP/pYe4f/khnt8Nwutn4YgNi0NaSjoAwMBQH+bWRVGqfEmUr2WHGo2roK67I0yL8juggiIhJhF3zvkj4OoTvLgXjLDACLyNjIc8I+sPv7GpEYrbFINtJRtUrOMAxyZV4ehWFYbGhhIn1yw3/riNRV+vhFEhQyz1+gFVXCpKkkMQkdoHdXJyIl9fX6VuM/TpG8z+6ie8fhaGYYu/RrfxHfN9Q1htRkS4c84fpzadx43jt5GWkg6Zvgx2NcqgTOVSKFGqOEzNTAABpCalITY8DqFPwxDk/xJpKenQN5DBqW1ttB/qjvod6kJPT6vf3LFPyJRn4urRWzi58Tz8zt+HIlMBY1MjONSyQ+mKJVHc2hzGpsYgBSEpPglRr2MR8vg1Xj4MhSJTgUKFjeHaxQXth7qjRuMqBfr3LTMzE7vnH8aOeQdQvrYd5h79HlZlLJQ+jhDiNhE5ffFxulLsAJCUkIwlg1bj6pFbcO3igokbR6JIscJKH0eTERG8D9/E9rn7EfwgBGYliqBpj0Zo1NkZNRpXgbGJ0WefL8+Q4/HNZ7h61AcX9ngjNiwOZaqUxoC5PdHEs0GB/uXVFQqFAud2XMaOeQcQHhQJq7IWaNHHDY06O6NSPYcvrjJLeZcC/0sPcfWoDy4fvI6k+GRUdi6PfnN6wKVdnQL3GomLeIv/9V+FO2f90ap/U4xdM+yLv2d5VSCLHcgqtkM//4GNU3ehRKlimLpjLBzdqqpkLE0T9OAVfvl2Ax54P0aZKqXRa8pXaNbLNc8rhuQZclw5dBO7FxxCcEAIajathvG/fQPbSqWUnJypS6D/SywfthZPfF6gYj0H9J3RDQ061cvzqe2pyWk4t+My9v10FOFBkajbqibGrh6K0hVslJxcM/n86Yclg1YjKT4Zo1YORruhLVX6hy2nxQ4iUvtHvXr1SNUe3XxK/SuMolZ63WnDlB2Ulpqu8jGlkpmZSft+OkptDXtSV4tBdGLDOZLL5Urbvlwupz9+O0NfFRtAHU370qnNF5S27U9RKBSUmpxKSQnJlJGeodKxpJaelk7v4pMoLSWNFAqFysZRKBR0eOUJamvYkzyth9C5nZeVOl56WjodWvEHeRTtRx1M+tCxtadV+v1ILfldCv0yagO5C08aWmM8Bd5/qZZxAfhSDjpW5/bYP5acmIJ1E7bh1KbzsKteBpM2f4vKzhVUPq46JSem4H/9fsH1Y75o3LU+xq0brrT7Jv5T9JtYLO73C+56BeCrMe0wYvkApVzEKDw4Ejf+uI0H3o8Q6P8KEcGRSE/963pApkVNUNLeCg41y6Faw8qo16ombBys8z2uOhERQp++we2z/nh04ymCH4QgPCgSyYkpHx5jbGIEaztLlK9tB0e3amjQsS4sSpfI99jyDDmWD1+Hs9suoUGnepi06VvVvUZex2DZ0LXwPX0PLfu6Yfz6b2BUSDXTElK5dykAy4euxZsXEej6XQcMWdRHbQeRC+xUzKfcPHkHK775DbFhcegytj0GzOuJQoULqW18VYmLeItp7RYg6P4rjFg2AF+Naafy+c1MeSbWT96BwytPoHlvV0zZNiZPZ/4SEa4f98XhFSdw72IAAMCqrAUq1nOAjb01zEoUgb6BDGkp6XgbGY83L8Lx3C8IcRHxAACHWuXg/nVTtBnYTKm3FFO22PA4/LnZC+d3XcarR68BACVKFUOFOvawcbCGuWVRGBobQJ6RifioeLwJjMBT3xeIeRMHIQTqta4JzwmdUNe9Zp7+36anpmNe92W4eeIO+s3ujq9ne6r8QLhCocCehUewdfZeVHetjPnHp+nEOSaJce+waeounNhwDjYO1pi4aSRqNVXvCrwCPxXzT+/evqMVI34jd+FJvct+Q5cP3dDqt4oxYbE0qMpY6mjal26duqP28ff87wi5C09aPHBVrn+Oz/wCaUzDaeQuPKmv3UjateAQhT5788XnKRQKCnnymg4uP06jG2Q9v4NJH1o1eiNFhUbn9VtRiTeB4bRs6Fpqa9iT3IUnjW86i47+eoreBIZ/8bkKhYKCA17Rtjn7qGfpYeQuPGlSyx/o1ePQXGWQZ8hp9leLyV140vF1Z/L6reTZxX1Xqa1hTxpZbzIlxCaqfXxlUSgUdHbHJfK0HkKtZd1p3cRtlJKUKkkW5HAqpsAUe7YHVx/TsJoTyF140pQ2P1LwwxDJsuTVu/gkGl5rInU07Uv+lx9KlmP73P3kLjxp+w/7c/R4hUJBB5YdozYGWfO8JzeeI3lG3o8FBPoH05JBq6mNQU/qYNKHts7eS6nJ0vzCZUtKSKZ1E7dRW8Oe1M64N638dgOFPHmd5+2lpabTkVUnqUvxAdS+UG86seFcjp/769hN5C486ciqk3keP79unLhNbQ170neNZ0j+/yYvnt0JpHFuM8ldeNLo+lPp2Z1ASfNwsX+GPENOh1eeoM7m/am1fg9aNXojvY2KlzRTTsnlcpreYQG1MehJPqfvSppFoVDQ4oGryF14fvFdg1wup6VD1pC78KQ5XX+i+JgEpeV4ExhO83v/TO7Ck/pXGEX3LgUobdu5cfPkHeplO5xa6XWnpYNXK/VdRExYLE1p8yO5C09aN3HbF98lXdjjTe7Ck9aM26K0DHl1cf81aqXXnRb0+Vlr3iVHvY6hpYNXUyu97tTNMmtBQmZmptSxuNhzIi7yLa0Y8Ru11u9BHkX70a4Fhyj5XYrUsT5rx7wD5C486diaP6WOQkREqcmpNNRxPHUvOeQ/y1qhUNCSQavJXXjSlpl7VPbL7Xfh/oeVUJtn7FbqyqDPSU9L/7B3PLTGeAq4/kQl48jlclo1eiO5C09aO2Hrfz4uKjSaOpv3p7GuMzRmVdHO+Qclf/eQE+/evqMtM/dQR9O+1NawJ62buI0S495JHesDLvZcCH4YQrM6/4/chSf1sBlKR1ad1MjlkY9uPqXW+j1oYd8VGrXn88wvkFrr96BlQ9d+8uu7Fhwid+FJ2+bsU3mW5HcptHTw6g9Tbaqe240Nj6OxjaaTu/CkX8duorSUNJWOp1Ao6NcxWX9Ejq09/cnH/NBtCXUw6UOvn4epNEtuZGZm0vQOC6idcW+NnP5MfpdC+346Sl1KDCR34Unzey2nNy++fDxE3dRS7ACWAHgMwB/AEQDmOXmephV7tvvej2h801lZB1jLfEO/r/5T5b+oOZWRnkHDak6gXrbD6d1bzdmDyLZ2wlZqpdednt8N+tvnH1x9TK1l6n8bfnLjOWpr2JMGVf2OwoMjVTLGq8eh1NduJHUw6UMX919TyRifIpfLaVr7BdTOqNe/1k/fuxRA7sKTds4/qLY8ORUTFktdSgykcW4zNWJagyir0PcvPUae1kPIXXjS1LY/0tPbL6SO9Z/UVeytAei//+/FABbn5HmaWuxEWXtEt8/eo7GuM8hdeFLP0sPowLJjlJyYLGmuI6tOkrvwpCuHb0ia478kxCZSZ/P+NKfrTx8+l5GeQUOqj6O+diMpKUH9P7+7Fx9QZ/P+1LvMNxTy9MurbnLj+d0g8rQaTJ7WQ+ixz3Olbjsn4iLfkqfVYBrbaPrfSnJCs9nUo9QwjT1QeXLjOXIXnnR+12VJcyTGvaNdCw5RN8tB5C48abL7XLrv/UjSTDmh9qkYAF0A7MrJYzW52LMpFAq6fc6fJjafQ+7Ck7oUH0BbZu6h2PA4tWdJfpdCntZDaGLzORo1BfNPW2buoVZ63T8sXfzjtzPkLjzp6tFbkmV6fjeIulkOol62w3O01DAngh68oq4Wg1TyByM3Tm2+QO7Cky4dyHq38NjnObkLTzq4/Lhkmb4kMzOTvqkzifqVHyXJ/H9kSDT9NmkbeZj1I3fhSdPaL9CKQs8mRbEfB/B1Th6rDcX+sYDrT2hOl8XUSq87tTPuTUsHr6YX94LVNv6hFX+Qu/DU+Bdg1OsYaq3fgzZM2UFyuZz6lR9Fo+tPlfyP0Yt7wdSlxEDqX2EUxUa8zde2Il5GUs/Sw6hHqWE5WnuvSnK5nAZVGUsj600mhUJBy4etpY6mfTVyqu5j1475kLvwpDPbL6ptzMe3ntGir1dSG4Oe1FrWneb3/pme+Um7dDEvlFbsAM4BePCJj84fPWbG+zl28ZntDAfgC8C3bNmyavtBKFPIk9e0YsRv1MGkD7kLT5rYfA5dOnBNpXsemZmZ1L/CKPqu8QyVjaFM0zssoD7lRpDP6bvkLjzJa6+31JGIKOuPcweTPjS20fQ8HzdJSkimoY7jqbN5fwr0V98f9s85tuZPcheeFHDtMXUpPoAW9VspdaQvyszMpMHVvqNvnaeodJy01HQ6t/Pyh5PhPMz60ZpxWygsKEKl46qS2vbYAQwAcB2ASU6fo2177P8UH5NAexcfpb52I7NW0pQaRltm7aGIl8o/SHf34gNyF550dsclpW9bFbLnUEe5TCGPov005uAzEdGlA9fIXXjS0iFrcv1chUJBP3RbQq31e9Dts/dUkC5vEuPeUVvDnjS6/lSNPgbzT4dXniB34amSP5Chz97Q+u93kKfVYHIXnjSg0hg6vPIEvYtPUvpY6pbTYs/XHZSEEG0BTAHQlIiS87MtbWJWvAh6ft8ZnhM74tZJPxxfdwa7FxzG7gWH4dSmFtoOboEGnZyUcoPtS/uvw9jECK5dXJSQXPXqta4FAHji8wJung006g47TTwbove0Ltiz6AhqNq2GVv2a5vi5R345Ce/DNzF8SX/Uda+pwpS5U9jcFI5NqsHv/H3o6QnUaVFD6kg50ry3K9ZN2IpL+6/D3rFcvreXlpIG78O38Ofm87jrFQA9mR4adqqHjiPaoK67Y4G7UUx+b433KwAjAGffX6DoBhGNyHcqLSGTydCwkxMadnJCeHAk/tx8AWe2XsSPPZajSPHCaN7LFa36N0Vl5wp5uoATEeHWyTuo26omCpkaq+A7UL6P7xpT062ahEk+bcC8nnjg/Ri/jt4ER7eqKGln9cXnBAeEYOOUnWjQqR48J3RUQ8rccWxcFX7n78OqrIXW3M7Q3LIoqrtWwc2TdzDwx1552gYRIeDqY5zdfgkX919DckIKStpbYeCPvdBmUHNYlCqu5NTaI1/FTkS6dQ3cfChpZ4WB83qh35zuuHPWH2e2XcSpTRdwbM1p2FayQYvebmje2zVXN6mIfBWNiJdR8JzYSYXJVce+ZlmpI/yLTCbDlO1jMLzmRPw8fB3+d3rWZ//oZmZmYung1TAtaoIJG0Zq5N2Bsn/OejLt2iut09IRO+YewLu3Sbm6+mNwQAgu7L4Crz3eCA+OgrGJERp3q482A5ujZtNqBW7v/FO0/mbWmkYmk8G5bR04t62Dd2+TcOXQDZzfdQU75h3A9rn7UaGOPZr2aIQmng1QqnzJz27rqe8LAEDVBpXUEV3pilmbSx3hk6zLWWLIor5YNXojvPZ4o0Uft/987PG1Z/DE5wWm7x6HYlZF1Zgy5z78nDXwj87nVG1QCUSE535BqN3881NIrx6/xuUD13Fp/zUEB4RkTTu510T/uT3RuIuLTlyGW5m42FWosLkp2g1piXZDWiL6dQwu7ruGS/uvYdO0Xdg0bRfK17ZD4y714drFBXbVy/xrb7Bc9TIYNL83ylWzleg7yB+zEpp7v9kO37jjzDYvrP9+Bxp2dv7kVFdCbCK2z9mHuu6OaNazkQQpc6ZwMe2Yfvmn7Nd1yJM3/yp2IsKLu8HwPnIT3odv4uXDUAghUN21Mkb9MhhNuzfU2B0HTcDFriYWpUvAc0IneE7ohIiXUbh88Aa8D9/Atjn7sG3OPpQqb42GnZzQoJMTajSuAn0DfZStUhp9pneVOnqutR7YDGe2XkSR4ppb7DKZDCOWDcD4JrNxZOXJT/6c9y46gqT4ZIxYPlAjp2CyZb+TaOLZUOIkuVPMOit3fFQCACA9LQP+lx7ixnFfXD/ui8hX0dDTE6jhVhWjVg5G464uSrmjVEFQIO6gpMliwuJw7XcfXD/mg7sXHiAjXQ7Toiao17oWnNvWQat+TfJ0hyIpyTPkSIpPVtnt15RpZqdFeHjtCXYGr4VJkb/ezsdFvEU/h1Fw82yAKdvGSJgwZ+Ii42FWorBSblWoTh5F+6GycwWYmhXC7bP+SE1Kg1EhQ9RtVRONPJzRoFM9mFtq5hSYFPjWeFoo5V0Kbp/1x80Td3DrlB+EAPaE/KbRe4va7vGtZxjTYDqGL+mP7h8dpN4ycw/2LDqCzY9W5OqAN8udjdN24cjKEzC3KgqXdnVQv0M91GlZQ+fuk6osXOxajogQ8yaW33qqwYRmsxH1Khrbnv8KPT09pKemo0/ZEajuWgVzj3wvdTydFxUaA4vSxXkHJgdyWuy8LkhDCSG41NXEY2QbhAdH4c65+wCAq0d9EB+dCI9v20icrGCwtC3Bpa5kXOyswGv0lQsKm5vi3M5LAIBzOy/BskwJ1GnpKHEyxvKGi50VeIZGBnD9ygXXj/kiPjoBd876o2n3RnyiC9Na/MplDEBDDyckJ6Rg57yDkGdkoqHHF6cxGdNYXOyMAajdogb09ASO/noKxqZGqNZQO8/2ZQzgYmcMAGBqZgK7GlnXXKniUgH6BnzuHtNeXOyMvZd9MS2HmnbSBmEsn7jYGXvP1MwEQNZFwhjTZlzsjP0TL6lmWo6LnbH3Ppwko/6TsRlTKi52xt5r0TfruuzO7WpLnISx/OFD/4y9V61BJZxVHJA6BmP5xnvsjDGmY7jYGWNMx3CxM8aYjuFiZ4wxHcPFzhhjOoaLnTHGdAwXO2OM6RgudsYY0zFc7IwxpmO42BljTMdwsTPGmI7hYmeMMR3Dxc4YYzqGi50xxnQMFztjjOkYLnbGGNMxSil2IcQkIQQJISyUsT3GGGN5l+9iF0KUAdAKwKv8x2GMMZZfythj/xnA9+BbADPGmEbIV7ELITwAvCaie0rKwxhjLJ++eDNrIcQ5ACU/8aUZAKYDaJ2TgYQQwwEMB4CyZcvmIiJjjLHcEER5m0ERQjgCOA8g+f2nbAG8AeBCROGfe66TkxP5+vrmaVzGGCuohBC3icjpS4/74h77fyGi+wCsPhowGIATEUXndZuMMcbyj9exM8aYjsnzHvs/EZGdsrbFGGMs73iPnTHGdAwXO2OM6RgudsYY0zFc7IwxpmO42BljTMdwsTPGmI7hYmeMMR3Dxc4YYzqGi50xxnQMFztjjOkYLnbGGNMxXOyMMaZjuNgZY0zHcLEzxpiO4WJnjDEdw8XOGGM6houdMcZ0TJ5vZp2vQYWIAvBSzcNaANDG+7FybvXS1tyA9mbn3DlXjogsv/QgSYpdCkII35zc3VvTcG710tbcgPZm59zKx1MxjDGmY7jYGWNMxxSkYl8vdYA84tzqpa25Ae3NzrmVrMDMsTPGWEFRkPbYGWOsQChwxS6EGCOEeCKECBBC/CR1ntwQQkwSQpAQwkLqLDkhhFgihHgshPAXQhwRQphLnelzhBBt3782ngshpkqdJyeEEGWEEF5CiEfvX9PfSZ0pN4QQMiGEnxDiD6mz5IYQwlwIcfD96/uREKKh1Jk+VqCKXQjRHEBnADWJqDqApRJHyjEhRBkArQC8kjpLLpwFUIOIagJ4CmCaxHn+kxBCBmA1gHYAqgHoLYSoJm2qHJEDmEhEVQE0ADBKS3Jn+w7AI6lD5MFKAH8SURUAtaBh30OBKnYAIwH8j4jSAICIIiXOkxs/A/gegNYcFCGiM0Qkf//PGwBspczzBS4AnhNRIBGlA9iLrJ0AjUZEYUR05/1/JyKrYEpLmypnhBC2ADoA2Ch1ltwQQpgBaAJgEwAQUToRvZU21d8VtGKvBMBNCHFTCHFJCOEsdaCcEEJ4AHhNRPekzpIPgwGckjrEZ5QGEPLRv0OhJQWZTQhhB6AOgJvSJsmxFcjaWVFIHSSXHABEAdjyfhppoxDCVOpQH9OXOoCyCSHOASj5iS/NQNb3WwxZb1mdAewXQjiQBiwN+kLu6QBaqzdRznwuNxH9/v4xM5A1ZbBLndlySXzic5K/LnJKCFEYwCEA44goQeo8XyKE6AggkohuCyGaSZ0nl/QB1AUwhohuCiFWApgKYJa0sf6ic8VORO7/9TUhxEgAh98X+S0hhAJZ13uIUle+//JfuYUQjgDsAdwTQgBZ0xl3hBAuRBSuxoif9LmfNwAIIQYA6AigpSb8Af2MUABlPvq3LYA3EmXJFSGEAbJKfRcRHZY6Tw65AvAQQrQHYAzATAixk4i+ljhXToQCCCWi7HdGB5FV7BqjoE3FHAXQAgCEEJUAGELDLz5ERPeJyIqI7IjIDlkvqrqaUOpfIoRoC2AKAA8iSpY6zxf4AKgohLAXQhgC6AXgmMSZvkhk/bXfBOARES2XOk9OEdE0IrJ9/5ruBeCClpQ63v/uhQghKr//VEsADyWM9C86t8f+BZsBbBZCPACQDmCAhu9FartfARgBOPv+3cYNIhohbaRPIyK5EGI0gNMAZAA2E1GAxLFywhVAPwD3hRB3339uOhGdlDBTQTAGwK73OwGBAAZJnOdv+MxTxhjTMQVtKoYxxnQeFztjjOkYLnbGGNMxXOyMMaZjuNgZY0zHcLEzxpiO4WJnjDEdw8XOGGM65v8pXdxhpWXuJAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x12561c5c0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from __future__ import division\n",
    "import matplotlib.pyplot as plt\n",
    "from numpy import sqrt #originally had from scipy import sqrt\n",
    "from numpy import meshgrid\n",
    "from numpy import arange\n",
    "\n",
    "xs = arange(-7.25, 7.25, 0.01)\n",
    "ys = arange(-5, 5, 0.01)\n",
    "x, y = meshgrid(xs, ys)\n",
    "\n",
    "eq1 = ((x/7)**2*sqrt(abs(abs(x)-3)/(abs(x)-3))+(y/3)**2*sqrt(abs(y+3/7*sqrt(33))/(y+3/7*sqrt(33)))-1)\n",
    "eq2 = (abs(x/2)-((3*sqrt(33)-7)/112)*x**2-3+sqrt(1-(abs(abs(x)-2)-1)**2)-y)\n",
    "eq3 = (9*sqrt(abs((abs(x)-1)*(abs(x)-.75))/((1-abs(x))*(abs(x)-.75)))-8*abs(x)-y)\n",
    "eq4 = (3*abs(x)+.75*sqrt(abs((abs(x)-.75)*(abs(x)-.5))/((.75-abs(x))*(abs(x)-.5)))-y)\n",
    "eq5 = (2.25*sqrt(abs((x-.5)*(x+.5))/((.5-x)*(.5+x)))-y)\n",
    "eq6 = (6*sqrt(10)/7+(1.5-.5*abs(x))*sqrt(abs(abs(x)-1)/(abs(x)-1))-(6*sqrt(10)/14)*sqrt(4-(abs(x)-1)**2)-y)\n",
    "\n",
    "for f in [eq1,eq2,eq3,eq4,eq5,eq6]:\n",
    "    plt.contour(x, y, f, [0])\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/noriakioshita/.pyenv/versions/anaconda3-5.1.0/lib/python3.6/site-packages/matplotlib/contour.py:1180: UserWarning: No contour levels were found within the data range.\n",
      "  warnings.warn(\"No contour levels were found\"\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAACz9JREFUeJzt3W+IpXd5h/Hr291JbOtKlChLMks3orZdNVYZgyX0XxIlalh9GUEJ9cWi1BBBiUmWvi9V/AMKZUlSCi6EovEPEqsb1EJfJLq7JqZxVZaQmjWGRGgxIDSz5PbFnOC6zs6ZOedxnpl7rg8szDnn4ffcZzl7zW+ePTOTqkKS1McfjD2AJGlYhl2SmjHsktSMYZekZgy7JDVj2CWpGcMuSc0YdklqxrBLUjO7xzjpHy78Ub32DQfGOLW2kNMPPQ7Aq/5i/6hzDKnjc9LWceLEiV9U1cunHZcxfqTA3j2X1VPPPrnp59XW8u6X3gTAl//330aeZDgdn5O2jiQnqmpp2nFeipGkZgy7JDVj2CWpGcMuSc0YdklqxrBLUjOGXZKaMeyS1Ixhl6RmDLskNWPYJakZwy5JzRh2SWrGsEtSM4Zdkpox7JLUzGBhT7IryfeTfG2oNSVJGzfkjv0W4NSA60mSZjBI2JMsAu8E7hxiPUnS7IbasX8auBV4fqD1JEkzmjvsSW4Anq6qE1OOO5TkeJLjy8vL855WknQBQ+zYrwYOJnkcuAe4Jsnnzz+oqo5U1VJVLS0sLAxwWknSauYOe1XdXlWLVbUfuBH4VlW9d+7JJEkz8X3sktTM7iEXq6rvAN8Zck1J0sa4Y5ekZgy7JDVj2CWpGcMuSc0YdklqxrBLUjOGXZKaMeyS1Ixhl6RmDLskNWPYJakZwy5JzRh2SWrGsEtSM4Zdkpox7JLUjGGXpGYMuyQ1Y9glqRnDLknNGHZJasawS1Izhl2SmjHsktSMYZekZgy7JDVj2CWpGcMuSc0YdklqxrBLUjOGXZKaMeyS1Ixhl6RmDLskNWPYJamZucOeZF+Sbyc5leTRJLcMMZgkaTa7B1jjLPCRqjqZZA9wIsmxqvrhAGtLkjZo7h17Vf28qk5OPn4WOAVcPu+6kqTZDHqNPcl+4I3Ag6s8dijJ8STHl5eXhzytJOkcg4U9yYuBLwIfrqpfnv94VR2pqqWqWlpYWBjqtJKk8wwS9iQLrET9aFXdO8SakqTZDPGumAB3Aaeq6pPzjyRJmscQO/argfcB1yR5aPLnHQOsK0mawdxvd6yq/wIywCySpAH4naeS1Ixhl6RmDLskNWPYJakZwy5JzRh2SWrGsEtSM4Zdkpox7JLUjGGXpGYMuyQ1Y9glqRnDLknNGHZJasawS1Izhl2SmjHsktSMYZekZgy7JDVj2CWpGcMuSc0YdklqxrBLUjOGXZKaMeyS1Ixhl6RmDLskNWPYJakZwy5JzRh2SWrGsEtSM4Zdkpox7JLUjGGXpGYGCXuS65P8OMnpJLcNsaYkaTZzhz3JLuBzwNuBA8B7khyYd11J0myG2LFfBZyuqseq6jngHuBdA6wrSZrBEGG/HHjinNtnJvdJkkYwRNizyn31Owclh5IcT3J8eXl5gNNKklYzRNjPAPvOub0IPHn+QVV1pKqWqmppYWFhgNNKklYzRNi/B7w6yRVJLgJuBL46wLqSpBnsnneBqjqb5EPAN4BdwN1V9ejck0mSZjJ32AGq6j7gviHWkiTNx+88laRmDLskNWPYJakZwy5JzRh2SWrGsEtSM4Zdkpox7JLUjGGXpGYMuyQ1Y9glqRnDLknNGHZJasawS1Izhl2SmjHsktSMYZekZgy7JDVj2CWpGcMuSc0YdklqxrBLUjOGXZKaMeyS1Ixhl6RmDLskNWPYJakZwy5JzRh2SWrGsEtSM4Zdkpox7JLUjGGXpGYMuyQ1Y9glqZm5wp7k40l+lOQHSb6U5JKhBpMkzWbeHfsx4HVVdSXwE+D2+UeSJM1jrrBX1Ter6uzk5gPA4vwjSZLmMeQ19vcDXx9wPUnSDHZPOyDJ/cDeVR46XFVfmRxzGDgLHF1jnUPAIYCXXXzpTMNKkqabGvaqum6tx5PcBNwAXFtVtcY6R4AjAHv3XHbB4yRJ85ka9rUkuR74GPA3VfWrYUaSJM1j3mvsnwX2AMeSPJTkXwaYSZI0h7l27FX1qqEGkSQNw+88laRmDLskNWPYJakZwy5JzRh2SWrGsEtSM4Zdkpox7JLUjGGXpGYMuyQ1Y9glqRnDLknNGHZJasawS1Izhl2SmjHsktSMYZekZgy7JDVj2CWpGcMuSc0YdklqxrBLUjOGXZKaMeyS1Ixhl6RmDLskNWPYJakZwy5JzRh2SWrGsEtSM4Zdkpox7JLUjGGXpGYMuyQ1M0jYk3w0SSW5dIj1JEmzmzvsSfYBbwV+Ov84kqR5DbFj/xRwK1ADrCVJmtNcYU9yEPhZVT080DySpDntnnZAkvuBvas8dBi4A3jbek6U5BBwCOBlF3spXpJ+X6aGvaquW+3+JK8HrgAeTgKwCJxMclVVPbXKOkeAIwB791zmZRtJ+j2ZGvYLqapHgFe8cDvJ48BSVf1igLkkSTPyfeyS1MzMO/bzVdX+odaSJM3OHbskNWPYJakZwy5JzRh2SWrGsEtSM4Zdkpox7JLUjGGXpGYMuyQ1Y9glqRnDLknNGHZJasawS1Izhl2SmjHsktSMYZekZgy7JDWTqs3/vdJJngH+Z5NPeymwHX8fq3Nvru06N2zf2Z17/f6kql4+7aBRwj6GJMeramnsOTbKuTfXdp0btu/szj08L8VIUjOGXZKa2UlhPzL2ADNy7s21XeeG7Tu7cw9sx1xjl6SdYift2CVpR9hxYU9yc5IfJ3k0yT+PPc9GJPlokkpy6dizrEeSjyf5UZIfJPlSkkvGnmktSa6fvDZOJ7lt7HnWI8m+JN9Ocmrymr5l7Jk2IsmuJN9P8rWxZ9mIJJck+cLk9X0qyV+OPdO5dlTYk/wd8C7gyqp6LfCJkUdatyT7gLcCPx17lg04Bryuqq4EfgLcPvI8F5RkF/A54O3AAeA9SQ6MO9W6nAU+UlV/DrwF+IdtMvcLbgFOjT3EDD4D/EdV/RnwBrbYc9hRYQc+CPxTVf0/QFU9PfI8G/Ep4FZg2/ynSFV9s6rOTm4+ACyOOc8UVwGnq+qxqnoOuIeVTcCWVlU/r6qTk4+fZSUwl4871fokWQTeCdw59iwbkeQlwF8DdwFU1XNV9X/jTvXbdlrYXwP8VZIHk/xnkjePPdB6JDkI/KyqHh57ljm8H/j62EOs4XLgiXNun2GbBPIFSfYDbwQeHHeSdfs0K5uV58ceZINeCTwD/OvkMtKdSf547KHOtXvsAYaW5H5g7yoPHWbl+b6UlS9Z3wz8e5JX1hZ4a9CUue8A3ra5E63PWnNX1Vcmxxxm5ZLB0c2cbYOyyn2jvy7WK8mLgS8CH66qX449zzRJbgCerqoTSf527Hk2aDfwJuDmqnowyWeA24B/HHes32gX9qq67kKPJfkgcO8k5N9N8jwrP+/hmc2a70IuNHeS1wNXAA8ngZXLGSeTXFVVT23iiKta6+8bIMlNwA3AtVvhE+gazgD7zrm9CDw50iwbkmSBlagfrap7x55nna4GDiZ5B/Ai4CVJPl9V7x15rvU4A5ypqhe+MvoCK2HfMnbapZgvA9cAJHkNcBFb/IcPVdUjVfWKqtpfVftZeVG9aStEfZok1wMfAw5W1a/GnmeK7wGvTnJFkouAG4GvjjzTVFn5bH8XcKqqPjn2POtVVbdX1eLkNX0j8K1tEnUm//aeSPKnk7uuBX444ki/o92OfYq7gbuT/DfwHHDTFt9FbnefBS4Gjk2+2nigqj4w7kirq6qzST4EfAPYBdxdVY+OPNZ6XA28D3gkyUOT++6oqvtGnGknuBk4OtkEPAb8/cjz/Ba/81SSmtlpl2IkqT3DLknNGHZJasawS1Izhl2SmjHsktSMYZekZgy7JDXza+zrashjCnWQAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x11a5b9c18>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from __future__ import division\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "xs = np.arange(-7.25, 7.25, 0.01)\n",
    "ys = np.arange(-5, 5, 0.01)\n",
    "x, y = np.meshgrid(xs, ys)\n",
    "\n",
    "eq1 = x**4-x**2+6\n",
    "eq2 = 12/abs(x)+1\n",
    "eq3 = 1/2*np.cos(2*x)+7/2\n",
    "\n",
    "for f in [eq1,eq2,eq3]:\n",
    "    plt.contour(x, y, f, [0])\n",
    "plt.show()"
   ]
  },
  {
   "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.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}