{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Optimization in 1-D and N-D. Debugging Functions\n",
    "====\n",
    "\n",
    "## Unit 11, Lecture 1\n",
    "\n",
    "*Numerical Methods and Statistics*\n",
    "\n",
    "----\n",
    "\n",
    "#### Prof. Andrew White,  April 19, 2020"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Goals:\n",
    "---\n",
    "\n",
    "1. Learn the meaning of root-finding and minimization, the two types of optimization\n",
    "2. Understand the iterative nature of these methods\n",
    "3. Debug common problems when defining functions, which is essential for optimization\n",
    "4. Be able to identify convex problems and understand their complexties\n",
    "5. Learn the two standard methods for minimize and root-finding and how to call them in Python"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "hide_input": false,
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Numerical Optimization - Root Finding\n",
    "====\n",
    "\n",
    "What is $x$ in\n",
    "\n",
    "$$\\cos (x) = x$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "To use a root-finding method, we must make our equation have one side be 0.\n",
    "\n",
    "$$\\cos (x) - x = 0$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Newton's Method - For finding Roots\n",
    "===="
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.7390851332151606 1.1102230246251565e-16\n"
     ]
    }
   ],
   "source": [
    "from scipy.optimize import newton\n",
    "\n",
    "root =  newton(lambda x: np.cos(x) - x, x0=0)\n",
    "print(root, np.cos(root) - root)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "$$x_{i+1} = x_i - \\frac{f(x_i)}{f'(x_i)}$$"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "<div class=\"naked-img\">\n",
    "    <img src=\"http://upload.wikimedia.org/wikipedia/commons/e/e0/NewtonIteration_Ani.gif\"/>\n",
    "    <p>     Source: Wikipedia, Newton's Method </p>\n",
    "</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "This method, like all this unit are *iterative*. This is modifiable by you, either through choosing tolerance, or maximum iterations. The functions have sensible defaults"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.73908511214521 3.526292269295794e-08\n"
     ]
    }
   ],
   "source": [
    "root =  newton(lambda x: np.cos(x) - x, x0=0, tol=1e-3)\n",
    "print(root, np.cos(root) - root)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.73908511214521 3.526292269295794e-08\n"
     ]
    }
   ],
   "source": [
    "root =  newton(lambda x: np.cos(x) - x, x0=0, tol=1e-3, maxiter=5)\n",
    "print(root, np.cos(root) - root)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "Newton's Method\n",
    "====\n",
    "*Use `root` instead, this is just for instructional purposes. `root` is better*\n",
    "\n",
    "**Type:** Root Finding\n",
    "\n",
    "**Discrete/Continuous:** Continuous\n",
    "\n",
    "**Dimensions:** 1\n",
    "\n",
    "**Derivative:** optional\n",
    "\n",
    "**Non-Convex:** not recommended\n",
    "\n",
    "**Python:** `newton`\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Let's try to solve a more interesting equation with Newton's method!\n",
    "\n",
    "$$ \\int_0^x e^{-s^2}\\,ds = \\frac{1}{2}$$\n",
    "\n",
    "Find $x$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "bad operand type for abs(): 'tuple'",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-5-b5fda6a0acf0>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m      3\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mequation\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[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      4\u001b[0m     \u001b[0;32mreturn\u001b[0m \u001b[0mquad\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;32mlambda\u001b[0m \u001b[0ms\u001b[0m\u001b[0;34m:\u001b[0m \u001b[0mnp\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mexp\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m-\u001b[0m\u001b[0ms\u001b[0m\u001b[0;34m**\u001b[0m\u001b[0;36m2\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m0\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[0m\n\u001b[0;32m----> 5\u001b[0;31m \u001b[0mroot\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mnewton\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mequation\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mx0\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[0;32m~/miniconda3/lib/python3.7/site-packages/scipy/optimize/zeros.py\u001b[0m in \u001b[0;36mnewton\u001b[0;34m(func, x0, fprime, args, tol, maxiter, fprime2, x1, rtol, full_output, disp)\u001b[0m\n\u001b[1;32m    316\u001b[0m         \u001b[0mq1\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfunc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mp1\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m*\u001b[0m\u001b[0margs\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    317\u001b[0m         \u001b[0mfuncalls\u001b[0m \u001b[0;34m+=\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 318\u001b[0;31m         \u001b[0;32mif\u001b[0m \u001b[0mabs\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mq1\u001b[0m\u001b[0;34m)\u001b[0m \u001b[0;34m<\u001b[0m \u001b[0mabs\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mq0\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[0m\u001b[1;32m    319\u001b[0m             \u001b[0mp0\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mp1\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mq0\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mq1\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mp1\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mp0\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mq1\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mq0\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    320\u001b[0m         \u001b[0;32mfor\u001b[0m \u001b[0mitr\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mrange\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mmaxiter\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[0;31mTypeError\u001b[0m: bad operand type for abs(): 'tuple'"
     ]
    }
   ],
   "source": [
    "#Note: this code is intentionally really bad\n",
    "from scipy.integrate import quad\n",
    "def equation(x):\n",
    "    return quad(lambda s: np.exp(-s**2), 0, x)\n",
    "root = newton(equation, x0=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Steps for Debugging Code\n",
    "===="
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "1. Correct any Python errors\n",
    "2. Restart the Kernel\n",
    "3. Correct any Python errors\n",
    "4. Print to discover logical errors"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "Most Common Mistakes\n",
    "----\n",
    "\n",
    "1. Copy-Pasta leaving in old variables\n",
    "2. Mixing Numpy arrays and `for` loops\n",
    "3. Forgetting order of function args and return values"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Review of Functions\n",
    "====\n",
    "\n",
    "For functions like `quad`, `newton`, and `minimize` the function should take in 1 argument and return 1 argument."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "Correct Examples\n",
    "====\n",
    "$$x^2 - 3x + 2$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2.0\n",
      "2.0\n",
      "[2. 2. 2. 2. 2. 2. 2. 2. 2. 2.]\n"
     ]
    }
   ],
   "source": [
    "def version_1(x):\n",
    "    return x**2 - 3 *x + 2\n",
    "\n",
    "version_2 = lambda x: x ** 2 - 3 *x + 2\n",
    "\n",
    "np_version = np.vectorize(version_1)\n",
    "\n",
    "print(version_1(3.))\n",
    "print(version_2(3.))\n",
    "\n",
    "some_threes = np.zeros(10)\n",
    "some_threes[:] = 3.0 # -> Notice how we don't replace the numpy array with 3, but instead make all elements in it equal to three\n",
    "#some_threes = 3 -> This would delete the numpy array and now some_threes would be a single 3\n",
    "\n",
    "print(np_version(some_threes))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "$$\\int_{\\pi}^x \\sin^2(s)\\,ds$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.5707963267948966\n"
     ]
    }
   ],
   "source": [
    "from scipy.integrate import quad\n",
    "import numpy as np\n",
    "from math import pi\n",
    "\n",
    "def integrate_sin2(x):\n",
    "    ans, err = quad(lambda s: np.sin(s) ** 2, pi, x)\n",
    "    return ans\n",
    "    \n",
    "print(integrate_sin2(2 * pi))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "2 Vectors, $\\vec{x}$ and $\\vec{y}$, where one component of $\\vec{x}$ is changing and we want to know the distance between the two vectors."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "5.408326913195984\n"
     ]
    }
   ],
   "source": [
    "from math import sqrt\n",
    "\n",
    "def distance(x, y):\n",
    "    sum_sq = np.sum((x - y)**2)\n",
    "    return sqrt(sum_sq)\n",
    "\n",
    "def my_distance(s):\n",
    "    x = np.zeros(3)\n",
    "    y = np.zeros(3)\n",
    "    x[0] = 2.0\n",
    "    x[1] = s\n",
    "    x[2] = -3.5\n",
    "    y[0] = 1.0\n",
    "    y[1] = -3.0\n",
    "    y[2] = 0.0\n",
    "    return distance(x, y)\n",
    "\n",
    "print(my_distance(1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "Incorrect Examples\n",
    "----"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "### No Return Value"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "slideshow": {
     "slide_type": "-"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "None\n"
     ]
    }
   ],
   "source": [
    "def version_1(x):\n",
    "    x**2 - 3 *x + 2\n",
    "\n",
    "print(version_1(3.))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "### Bad Return Value"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(1.5707963267948966, 1.743934249004316e-14)\n"
     ]
    }
   ],
   "source": [
    "def integrate_sin2(x):\n",
    "    return quad(lambda s: np.sin(s) ** 2, pi, x)\n",
    "\n",
    "print(integrate_sin2(2 * pi))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "### Too many arguments"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "distance() missing 1 required positional argument: 'y'",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-11-6fd995bbf352>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m      3\u001b[0m     \u001b[0;32mreturn\u001b[0m \u001b[0msqrt\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0msum_sq\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;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 5\u001b[0;31m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdistance\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m1\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[0m",
      "\u001b[0;31mTypeError\u001b[0m: distance() missing 1 required positional argument: 'y'"
     ]
    }
   ],
   "source": [
    "def distance(x, y):\n",
    "    sum_sq = np.sum((x - y)**2)\n",
    "    return sqrt(sum_sq)\n",
    "\n",
    "print(distance(1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Let's return to our example from above:\n",
    "\n",
    "$$ \\int_0^x e^{-s^2}\\,ds = \\frac{1}{2}$$\n",
    "\n",
    "Find $x$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.0\n"
     ]
    }
   ],
   "source": [
    "#note still wrong\n",
    "def equation(x):\n",
    "    ans, err = quad(lambda s: np.exp(-s**2), 0, x)\n",
    "    return ans\n",
    "root = newton(equation, x0=0)\n",
    "print(root)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.0"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "equation(root)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "We forgot to rearrange the equation to be equal to $0$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.5510394276090265 0.4999999999999999\n"
     ]
    }
   ],
   "source": [
    "root = newton(lambda x: equation(x) - 0.5, x0=0)\n",
    "print(root, equation(root))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "Scope\n",
    "----\n",
    "\n",
    "Scope means the set of variables and functions which are defined and accessible in your code"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "x = 4\n",
    "y = 2\n",
    "\n",
    "#Right now, there is x,y and all other functions/variables defined or imported above in the scope"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "Scopes nest\n",
    "----"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [],
   "source": [
    "x = 4\n",
    "y = 2\n",
    "#Here, I have x and y in scope\n",
    "def scope_example():\n",
    "    z = 2\n",
    "    #Here, I have x,y and z in scope\n",
    "    \n",
    "#Here I again have only x and y in scope"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "scrolled": true,
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2 Before function\n",
      "25 Inside function\n",
      "2 After Function\n"
     ]
    }
   ],
   "source": [
    "x = 4 \n",
    "y = 2\n",
    "print(y,\"Before function\")\n",
    "def scope_example():\n",
    "    y = 25 #This is a new version of y that exists only in this scope\n",
    "    print(y, \"Inside function\")\n",
    "scope_example()\n",
    "print(y, \"After Function\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[2] Before function\n",
      "[25] Inside function\n",
      "[25] After Function\n"
     ]
    }
   ],
   "source": [
    "x = 4 \n",
    "y = [2]\n",
    "print(y,\"Before function\")\n",
    "def scope_example():\n",
    "    y[0] = 25 #Here I'm not creating a y, but modifying y\n",
    "    print(y, \"Inside function\")\n",
    "scope_example()\n",
    "print(y, \"After Function\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "Things to remember about scope:\n",
    "\n",
    "1. Scopes nest, so that you can access things above your scope\n",
    "2. You can modify variables from any scope you can see, but ones you create disappear outside of the scope"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Returning to Optimization\n",
    "===="
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Applications of Optimization\n",
    "===\n",
    "\n",
    "1. Solving non-linear equations\n",
    "2. Solving systems of equations\n",
    "3. Optimizing equations with or without constraints\n",
    "4. Fitting models to data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Choosing which method to use\n",
    "----\n",
    "\n",
    "There are five things to consider when doing optimization:\n",
    "\n",
    "1. Is it 1 or N dimensions?\n",
    "2. Are you minimizing or root-finding?\n",
    "3. Are there constraints?\n",
    "4. Are there bounds?\n",
    "5. Is it convex or non-convex?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Identifying Convexity\n",
    "====\n",
    "\n",
    "If a problem has more than one minimum (derivative is 0), then the problem is non-convex. **The opposite of convex is non-convex. A concave function can be made convex with a negative sign. A non-convex function cannot be made convex**."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "Knowing about convexity can come from:\n",
    "\n",
    "* Plots\n",
    "* Running convex optimization in two starting positions and getting different minimums\n",
    "* Knowing something specific about the problem"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "Consider a function with two minimums:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [],
   "source": [
    "def two_well(x):\n",
    "    if x < 0.125:\n",
    "        return (x + 2) ** 2\n",
    "    if x >= 0.125:\n",
    "        return (x - 2) ** 2 + 1\n",
    "    \n",
    "np_two_well = np.vectorize(two_well)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAW4AAAD8CAYAAABXe05zAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xd4m9d59/HvAThAcYokOMQhkqKovalJDUvylEccrzixHad26qRNmrhNm93Ezdukad2MNs1y7LixHTux470iydrWprY4RFEkRYp77wXgvH+QVBVZEkGJwIMHuD/XpcukBeL5SQZuH9zPGUprjRBCCPOwGB1ACCHE2EjhFkIIk5HCLYQQJiOFWwghTEYKtxBCmIwUbiGEMBkp3EIIYTJSuIUQwmSkcAshhMkEeeJJ4+PjdUZGhieeWggh/NKhQ4eatNZ2dx7rkcKdkZFBfn6+J55aCCH8klLqrLuPlVaJEEKYjBRuIYQwGSncQghhMlK4hRDCZKRwCyGEybg1q0QpVQF0Ak7AobXO9WQoIYQQlzeW6YBrtdZNHksihBDCLdIqEUKIcfBBYT1P7yrD6fL8cZDuFm4NbFJKHVJKPXapByilHlNK5Sul8hsbG8cvoRBCmMArh6r43d4KrBbl8Wu5W7jztNYLgVuALyilVl/8AK31U1rrXK11rt3u1qpNIYTwC1prDla0sjgj1ivXc6twa61rhv/ZALwOLPFkKCGEMJPTDV20dA+wLCvOK9cbtXArpcKVUpEjXwM3Aic9HUwIIcxiX1kzAMsyvVO43ZlVkgi8rpQaefyLWus/ezSVEEKYyP6yFiZF20iLDfPK9UYt3FrrMmCeF7IIIYTpaK3ZX97M6ql2hge4HifTAYUQ4hqcaeyiqWuApVneuTEJUriFEOKa7C1rAfDajUmQwi2EENdkf1kzSVE20mMneO2aUriFEOIqaa3ZV9bCsqxYr/W3QQq3EEJctbKmbpq6+lnqxTYJSOEWQoirdn7+thRuIYQwh/1lLSREhpIR573+NkjhFkKIqzLU325mWVacV/vbIIVbCCGuSnlTNw2d/V6dvz1CCrcQQlyF/eXen789Qgq3EEJchX1lzdgjQ8mKD/f6taVwCyHEGI30t5dmenf+9ggp3CKgNXT0eeWoKeFfzjR2U9/Rz4op8YZcXwq3CFhtPQMs+cEWfvh+kdFRhMnsOTN0bnpetvf72yCFWwSww5WtAPxmV7nBSYTZ7C5tIiUmzKv7k1xICrcIWIfOtp7/ur130MAkwkycLs3eM82szI43pL8NUrhFAMuv+L/CPbJ0WYjRnKxup6PPwQqD2iQghVsEqEGni2Pn2nhwWToTQqzsLm0yOpIwid3D/W2jbkyCe2dOCuF3Cms66Bt0sSwrjurWXj6Uwi3ctKe0mWmJkdgjQw3LICNuEZAOVgytesudHEtedjxljd3UtPUanEr4ur5BJwcrWgxtk4AUbhGg8itaSYsNIynaxsqpQx95ZdQtRnO4spV+h4s8A9skIIVbBCCtNQcrWlicMbQ50LTESOIjQqXPLUa1p7QZq0UZsrHUhaRwi4BT3tRNc/fA+cKtlGJldhy7S5vQWlZRisvbfaaJuanRRNqCDc0hhVsEnJH+9kjhBsjLjqepa4BT9Z1GxRI+rrNvkOPn2g1vk4AUbhGADla0EhsewhT7/+3qlpc93Oc+Le0ScWn7y1pwurThNyZBCrcIQAcrWsidPPEvVr1Nigkjyx4ufW5xWbvPNBEaZGFh+kSjo0jhFoGloaOPs809f9EmGbEyO5795S0MOFwGJBO+bndpE4szYrEFW42OIoVbBJaDw8vcF2d+tHDnZcfTM+DkSGXrR35PBLba9l5K6rtYNdX4/jZI4RYB5mBFC2HBVmZNivrI7y3LisOikHaJ+IhdJUOvidU5doOTDJHCLQLKwYoWFqTHEGz96Es/OiyYuakxshBHfMSO040kRIYyPSnS6CjAGAq3UsqqlDqilHrHk4GE8JTOvkGKajvIvUR/e8TK7HiOnWuXbV7FeU6X5sPTTazOsRu2jevFxjLi/jIgR4UI0zpc2YZLw5IrFO410+w4XVraJeK84+faaO8d9Jk2CbhZuJVSqcCtwNOeDONyafoGnZ68hAhg+RUtWC2KBekxl33MgrQYIm1B7DjV6MVkwpftLGlCKViV7Rs3JsH9EfdPga8CHpsn1TfoZMkPPuBXO8546hIiwB0ob2HWpCjCQy+/m3GQ1cLK7Hh2lDTK8ncBwM7TjcxNiWZieIjRUc4btXArpW4DGrTWh0Z53GNKqXylVH5j49hHK7ZgK0nRNvaXtYz5Z4UYTb/DydGqtkvO377YddPs1HX0UVLf5YVkwpe19wxypLLVp9ok4N6IOw+4QylVAfwBWKeUeuHiB2mtn9Ja52qtc+32q/tDLs2MG942UdolYnydrO6g3+Ficcboq95G3qTbTzV4OpbwcbvPNOHSvjMNcMSohVtr/Q2tdarWOgO4H9iqtX7QE2GWZsbS73BxrKrdE08vAtj5gxPcGHEnR4cxLTGSHSXS5w50O0saibQFsSDt8vdFjOBT87iXZMaiFOyXg1vFODtQ3kKWPZz4CPeOm7pump2DFS109zs8nEz4Kq01O0sayZsST9Al5v0baUxptNbbtda3eSpMzIQQpidFsa9cCrcYPw6ni4PlLSzLcn9XtzU5dgadmj1n5LUYqM40dlHT3udzbRLwsRE3DLVLDp1tlY1+xLgprO2gs98xpsK9KGMiE0Ks7CiRPneg2nF+mbvvTAMc4XOFe1lWHH2DLk5UtxkdRfiJfcOtt2WX2FjqckKDrKyYEs/2UzItMFDtLGkkyx5O6sQJRkf5CJ8r3EuG31z7ZFqgGCf7yob62wlRtjH93Jppds619lLe1O2hZMJX9Qw42FvWzHU5CUZHuSSfK9yx4SFMS4w8P0oS4lpcTX97xHXnpwXK7JJAs6e0mQGHi/UzpHC7bVnWUJ970Cl9bnFtrqa/PSItdgJZ9nCZFhiAthQ3EBEa5NaCLSP4ZOFemhVHz4CTE9Uyn1tcm5GVuGPpb19oTY6dfWXN9A7IorBAobVmW3EDq6bGExLkkyXSNwv3SJ9blr+La7WvrJms+LH3t0esm55Av8Mle3QHkMLaDuo6+lg73TfbJOCjhTs+IpSpCRHS5xbXxOnSHChvYelVtElGLM2MIyI0iC1F9eOYTPiybcVDU0DXTpPCPWZLs2LJr2jBIX1ucZUKa0b621ffpwwJsrA6J54txQ24XDItMBBsKW5gXmo09kj3VtkawXcLd2Yc3QNOCmo6jI4iTOr8/O1rGHEDrJ+eSGNnv9xzCQDNXf0crWpj3fREo6Ncke8W7qyR+dzSLhFXZ6S/nXiV/e0Ra6cnYFFIuyQADC24Grq34ct8tnAnRNrIsodL4RZXZTz62yNiw0NYmD6RD4pk+bu/21rcQEJkKLMmRRkd5Yp8tnADLM+K42BFq/S5xZgV1V57f/tC62ckUljbQU1b77g8n/A9g04XO0saWTstAYvFNw4FvhyfLtx52fF09Ts4dk56i2Js9p4Zn/72iOuHV9BtKZZRt786WNFCZ7+DdT66WvJCPl24l2fFoRTskTm0Yox2n2kiy37t/e0R2QkRpMdOkD63H9tS1EBI0NCZo77Opwv3xPAQZiZHyeIHMSYDDhcHylvG9Q2olGL9jAT2nGmmZ0AOV/A3Wms2FdaxMjv+iodJ+wqfLtww1C45UtkmS46F245WtdEz4GTFlPEdOV0/I5EBh4tdp2Ug4W+Kajupaunlxpm+PQ1whM8X7hVT4hhwus6fGSjEaHaXNmFRQ6228bQkM5ZIm6yi9EebCutQCq6Xwj0+lmTGEmxV7D4joxzhnt2lTcxJiSZ6QvC4Pm+w1cKaHDtbixtwyipKv7KpoJ7cyRPdPpPUaD5fuCeEBLEgbSJ7SmU+txhdV7+Do1Vt5HnoBtONs5Jo6hrgcGWrR55feF9VSw+FtR3cODPJ6Chu8/nCDUN97pM17bT1DBgdRfi4A+XNOFzaY4V77TQ7IVYLfz5Z55HnF963uXCo9XWDSdokYJrCHYfWsvxdjO7D082EBllYNHmiR54/0hbMyqnx/PlknZxF6Sc2FtQxLTGSjPhwo6O4zRSFe15aDOEhVnZLu0SMYs+ZJhZnxGILtnrsGjfPSqK6rVc2QPMDLd0DHKxo4aZZ5hltg0kKd7DVwpLMWLlBKa6oobOP4rpOVmSP72ySi10/MxGrRUm7xA9sKarHpYfuXZiJKQo3DPW5yxq7qWvvMzqK8FEjy9w9vfItNjyEpZmx/LlACrfZbSqsZ1K0zec3lbqYaQr3yGKK3bKKUlzG7tImomxBzJoU7fFr3Tw7idKGLkobOj1+LeEZPQMOdp1u5MZZSSjl25tKXcw0hXt6UiSx4SFSuMUlaa3ZXdrMiinxWL2ws9vI1DFpl5jXtuJG+gZd3GSyNgmYqHBbLIoVU+L4sLRJ7uaLjzjb3EN1Wy95Hu5vj0iKtrEgPUbaJSb27oka4iNCzx9ObiamKdwAq3PsNHT2U1wnH0/FXxrZiMxT87cv5eZZSZys7qCqpcdr1xTjo2fAwdbiBm6ZneSVT2jjzVyFe6odgJ0ljQYnEb5mZ0kjKTFhZHpxLu7IR+yNMuo2na3FDfQNurh1brLRUa6KqQp3UrSNaYmR7DwthVv8n0Gniz1nmlmdY/fqTaaM+HBmJEfx7olar11TjI93j9dijwxlcYb52iTgRuFWStmUUgeUUseUUgVKqX/xRrDLWZ0Tz8HyVtkTWZx3+GwrXf0O1uTYvX7t2+clc6SyTdolJtLd72DbKfO2ScC9EXc/sE5rPQ+YD9yslFrm2ViXtzrHzoDTxf4y2eZVDNlR0kiQRXl84c2l3D53EoCMuk3kfJtkjjnbJOBG4dZDuoa/DR7+Zdi0jqHlzBZ2SJ9bDNtR0sjC9IlE2cZ3G1d3pMVOYF5aDG8fq/H6tcXVGWmT5Jq0TQJu9riVUlal1FGgAdistd5/icc8ppTKV0rlNzZ6rqjagq0sy4qTG5QCgMbOfgpqOlgzzfttkhG3z02moKaDssau0R8sDDXSJtlg4jYJuFm4tdZOrfV8IBVYopSafYnHPKW1ztVa59rtnn0TrZ5qp6ypW/qKgl3DN6pHZhwZ4ba5k1AK3jku7RJft6W4gX6Hiw0mbpPAGGeVaK3bgO3AzR5J46bVwzehZHaJ2FHSSFx4iKF7TSRF21icEctbx2pkcZiPe/d4DQkmb5OAe7NK7EqpmOGvw4DrgWJPB7uSKfZwUmLCpF0S4Fwuza7TTazOsWMx+GPv7fMmUdrQxal6WRzmq9p7B9lW3MiGOcmmbpOAeyPuZGCbUuo4cJChHvc7no11ZUopVufEs6e0mUGny8gowkAna9pp6R5gdY73Vktezi2zk7Ao5CalD3v/RC0DThcfX5BidJRr5s6skuNa6wVa67la69la6+95I9hoVk+10zl8vqAITCOfuFYZ2N8eER8RSl52PO8cr5V2iY9642g1mfHhzE31/O6RnmaqlZMXWpE9tAuctEsC146SRuakRPvMydy3z53E2eYejp9rNzqKuEhNWy/7ylq4c36K6bZwvRTTFu7osGAWpMWw/ZQU7kDU0TfI4co2n2iTjLhpdhIhQRZeP1JtdBRxkbeGW1gfmz/J4CTjw7SFG2Dt9AROVLfT0CGn4gSaPaVNOF2aNTkJRkc5LzosmBtmJvLWsRoGHHLvxZe8caSaBekxpjoQ+EpMXbjXzxh602471WBwEuFtW4sbiLQFsSA9xugof+HuhSm0dA/Iyl4fUlzXQXFdJ3fON/9NyRGmLtzTEiOZFG1jS5EU7kDicmm2FjeyJsdOsNW3XsKrptqJjwjhtcPnjI4ihr1xpAarRZl2C9dL8a1X/RgppVg3I4EPS5vodziNjiO85ER1O01d/ec/cfmSYKuFO+alsKWogbaeAaPjBDyXS/PW0WpWT433mZvY48HUhRtg/fREegacsltgANlS3IBF4VP97QvdtTCFAadLlsD7gL1lzdS093GnH8zdvpDpC/fyKXHYgi1sLZZ2SaDYWlzPwvSJxIaHGB3lkmZNimJaYqS0S3zAK/lVRNqCTHkg8JWYvnDbgq3kTYlnS3G9LHwIAHXtfZys7mD9jESjo1yWUoq7FqZwuLKN8qZuo+MErPbeQd4/WcfH5k/CFmw1Os64Mn3hBlg3I4Gqll7OyLaafm9kBpEv9rcvdOeCFCwKGXUb6O1jNfQ7XNyXm2Z0lHHnF4V77bShN7HMLvF/W4oaSJ0YxtSECKOjXFFilI1VU+386dA5nC75JGiEV/KrmJ4UyZwU8y9xv5hfFO5JMWHMSI6SPref6xt0sru0ifXTE0yxbPn+xWnUtvexo0Rel952qq6TY+fauTc3zRSvlbHyi8INsH56AvlnW2nvGTQ6ivCQvWXN9A46WefD/e0LrZ+RSHxECH84UGV0lIDzSn4VwVbFnX6yxP1iflO4105PwOnS7JDDFfzW1qIGJoRYWZppjk3wQ4Is3L0olS3FDbItgxcNOFy8fqSa62ckEudHc7cv5DeFe35aDLHhIWwpqjc6ivAArTVbixtYmR1vqhkCn8hNw+nSvHJIblJ6y9biBpq7B7g3N9XoKB7jN4XbalGsm57AtuIGOVzBDxXWdlDd1uvzs0kulmWPYGlmLH88WIVLblJ6xe/3nyUpymboOaSe5jeFG+CmWUl09DnYV9ZsdBQxzjYW1GNRcL1J+tsXun9JGpUtPfK69IKKpm52nW7ik0vSCfKxfWzGk1/9yVZNjScs2MqmAmmX+JtNBXXkTo41Zc/yltnJRNmCeOmg3KT0tBcPVGK1KO5f4n9zty/kV4XbFmxlTY6dTYV18rHUj5xt7qa4rpMbZ5lvtA1Dr8u7Fqay8WQdTV39RsfxW32DTl7Jr+LGmYkkRtmMjuNRflW4AW6anUh9Rz/HzslZlP5i5BOUmfebeHDZZAacLv5woNLoKH7rvRO1tPYM8tCyyUZH8Ti/K9zrpiUSZFFsKpR2ib/YWFDHzOQo0mInGB3lqmUnRLAyO54X9lXikJvnHvHCvrNk2cNZPiXO6Cge53eFO3pCMMuy4thYUGd0FDEOGjv7OVTZaurR9oiHV2RQ19EngwoPKKhp53BlGw8sneyXKyUv5neFG+CmWYmUNXZT2tBpdBRxjTYX1qP1UAvM7NZNTyB1Yhi/21NhdBS/88K+SmzBFu5Z6L9zty/kl4X7hplDo7ONMrvE9DYW1JEeO4FpiZFGR7lmVovioWWT2V/eQlFth9Fx/EZbzwBvHKnmjnmTiJ4QbHQcr/DLwp0UbWNeWgybpF1iah19g+w508RNsxL95uPvfblphAZZeG7vWaOj+I0XD1TSO+jkkZWZRkfxGr8s3AA3zkzk2Ll2att7jY4irtLQKljtF/3tERPDQ7hzfgpvHKmWDdHGwYDDxe/2VLAyO57pSVFGx/Eavy3cI2/2P5+UUbdZbSqoJz4ilIXpE42OMq4+vWIyvYNOXpSpgdfsvRO11Hf08+iqwBltgx8X7uyECKYnRfLeCTmw1Yx6BhxsLW7gplmJWCz+0SYZMWtSNHnZcTy7u5x+h9PoOKaltebpD8uYYg9njR/vS3Ipflu4AW6dk8zBilbq2mVLTbPZVtxI76CT2+b6537Kn1s9hYbOft48WmN0FNM6UN7CyeoOHlmZ6Xf/cx+NXxfuDXOTAWTUbULvnqghPiKUJSbZe3usVk2NZ0ZyFL/ZWSbbM1ylZz4sZ+KEYO5aEBhTAC80auFWSqUppbYppYqUUgVKqS97I9h4mGKPYEZyFO9K4TaV7v6hNsmGOUlY/XQkpZTisdWZnG7oYrscbTZmp+s72VRYz0PLJhMWYp792ceLOyNuB/AVrfUMYBnwBaXUTM/GGj+3zU3m0NlWatpkdolZbCluoG/Q5bdtkhG3zZ3EpGgbv95RZnQU0/nljjOEBVv5TF5g3ZQcMWrh1lrXaq0PD3/dCRQBKZ4ONl42zJF2idm8e7yGhMhQcif712ySiwVbLTyyMpP95S0cqWw1Oo5pVLX08ObRGj61NJ3Y8BCj4xhiTD1upVQGsADY74kwnpAZH86sSdIuMYuufgfbTjWyYU5yQNxwun9JOjETgvmfraVGRzGNX+88g0XBX6/KMjqKYdwu3EqpCOBV4HGt9UfW6yqlHlNK5Sul8hsbfevA3lvnJnOkso1zrT1GRxGj2FJUz4DDxW3DN5b9XURoEJ9dmcmW4gZOnGs3Oo7Pa+jo4+X8c9yzKJWkaP/ec/tK3CrcSqlghor277XWr13qMVrrp7TWuVrrXLvdt+ZU3jrcLnn/hCzG8XVvH6slOdrmd4turuThFRlE2YL4ry2njY7i8575sByH08Xn10wxOoqh3JlVooBngCKt9Y89H2n8TY4LZ05KNO9Iu8SndfQNsrMkcNokIyJtwTy6MosPiuo5WS2j7stp6urn+X1nuX3eJCbHhRsdx1DujLjzgIeAdUqpo8O/Nng417i7bW4yx6raqGjqNjqKuIzNBfUMOF3cGiBtkgt9Ji+DSFsQP9sqo+7L+eX2M/QNOvnS+qlGRzGcO7NKPtRaK631XK31/OFf73kj3Hi6Y/4klII3jlYbHUVcxhtHq0mLDWNBWozRUbwuOiyYR/Iy2VhQT2GNbPl6sdr2Xp7fd5a7F6YyxR5hdBzD+fXKyQslR4exPCuON45Uo7WsVPM19R197C5t4uPzU/xmC9exeiQvkyhbEP+56ZTRUXzOz7aWorWW0fawgCncAHcuSKGiuYejVXKQsK95+1gNLj303yhQRU8I5m/XZrO1uIF9Zc1Gx/EZlc09vHywivsXp5v63NHxFFCF+5bZSYQGWXjjiLRLfM1rh6uZlxZDVoB/DP7MigySomz88P1i+WQ47KcflGC1KL64LtvoKD4joAp3pC2YG2Ym8vbxWgblpG2fcaquk8LaDj4+37+XuLvDFmzlH27I4WhVm+wlz9AhwK8frebhFRkkRgXuvO2LBVThBvj4ghRaugfYWeJbi4QC2RtHq7FaFLfNk8INcPeiVHISI3hy46mAHmBorfnXd4qICQvmC2tltH2hgCvcq3PsxIaH8Lq0S3yCy6V580g1q6fGEx8RanQcn2C1KL5603TKmrp5KYBPydlcWM/esmb+/oYcosMC4xBgdwVc4Q62Wrh9bjKbC+vp7JMz/4y2v7yFmva+gL4peSnrZySwYkocP9pUQkv3gNFxvG7A4eLf3i8mOyGCTy1JNzqOzwm4wg1DMxf6HS7elx6i4d44Uk14iJUbZ/rPgcDjQSnFE3fMoqvfEZDTA5/fd5bypm6+desMgqwBWaauKCD/RuanxZAZH86rh84ZHSWg9Qw4eOd4DbfMSQ7IzfBHk5MYyaeXT+alA5UBtRS+vqOPn2wuYU2OnbXTEoyO45MCsnArpbhnUSr7y1tkCbyB3j1eS/eAk08sTjM6is96/PocYieE8N23CgJmeuD33ilkwOniex+bZXQUnxWQhRvgnkWpWBS8nF9ldJSA9XJ+FVnx4X5/YMK1iA4L5mu3TOfQ2Vb+cND/X6vbTzXw7vFa/m5tdsBvJHUlAVu4E6NsrJ2WwJ8OncMRwFOujHKmsYuDFa3cm5sWsEvc3XXvolSWZcXyg3eLqGvvMzqOx/QNOvnOmwVk2cN5bE3gHpLgjoAt3AD3LU6jobOfHTKn2+tezq/CalHcvUhmk4xGKcUP75rLgNPFP7950m9bJj/5oITKlh7+9c7ZhAbJPY8rCejCvW56AvERofwxAD6C+pJBp4tXD1WzdloCCZGyGs4dGfHhfOXGHDYX1vOeHx4IcuhsC0/tLOOTS9JZMSXe6Dg+L6ALd7DVwt0LU9hS3EBDp/9+BPU120810tTVLzclx+iRvEzmpkbznTdP+tXrtWfAwT+8fIzUiWF869YZRscxhYAu3AD35qbhdGleOywrKb3ljwersEeGsnaabx1x5+uCrBZ+dO88uvod/NMrx3G5/KNl8sP3i6ls6eHJe+YRERpkdBxTCPjCnZ0QQe7kibx8sMpve4e+pKGjj22nGrhrYYosrLgKUxMj+fZtM9lR0sizeyqMjnPNthU38NzeszySl8myrDij45iGvHOATyxOo6ypm72yB7LHvXSgCqdLyzLma/Dg0nSun5HIv79fTEGNeRfmVLf18vcvH2VmchT/dNM0o+OYihRu4PZ5k4iZEMzze88aHcWvDTpdvHjgLGty7DJH9xoopfiPe+YSMyGYL/z+MO095ttzZ9Dp4u9ePIzDqfnFAwuxBcsskrGQws3QHsj35aaxqbCe2vZeo+P4rQ8K66nv6OehZZONjmJ6seEh/OKBhVS39fLlPx7BabJ+9w/eK+JwZRs/vHsOGfHyP/GxksI97MGlk3FpzUsHZGqgpzy/7ywpMWGsnS77T4yH3IxYnrhjFttPNfIjE21E9dKBSp7dXcFf5WVw21zZg/1qSOEelh43gety7Lx0oJIBh6ykHG+lDZ3sOdPMA8vSsVpkpeR4eWDpZD65JJ1fbD/Dm0d9f2bUnjNN/PMbJ1mTY+dbG2Tq39WSwn2Bh5ZPprGzn40F/rfAwWgv7KskxGrhvlyZuz3enrhjJkszY/nHV4759MlOp+o6+ZsXDpMZH87PPrVAZhVdA/mbu8CanATSYsN4fp/cpBxP3f0OXj10jg1zkuSUGw8IDbLy1KdzmWKP4PMvHOJYVZvRkT6ioqmbB5/Zjy3Ywm8/s5gom5xocy2kcF/AalE8uHQyB8pbKK7rMDqO33j9SDWd/Q4eWi43JT0lOiyY5x5ZQmx4CA8/e8Cn9u+uaevlgaf343C6eOHRpaTFTjA6kulJ4b7Ifblp2IItPPthhdFR/ILLpfnt7nLmpESzMF22b/WkhCgbL352GeEhQXzyN/s4UtlqdCTKm7q591d76egd5LlHljI1MdLoSH5BCvdFJoaHcPfCVF4/Uk1jZ7/RcUxv26kGyhq7+eyqTNm+1QvS4ybw8ueXExsewoNP72fXaeN63gU17dz7qz30DTp56bFlzEmNNiyLv5HCfQmPrsxk0OXi+b0VRkdi7KlcAAAPDElEQVQxvd/sKiM52saGOclGRwkYKTFhvPy55aTFTuAzzx7kub0VXs+wqaCOT/x6HyFWCy9/fjmzU6Rojycp3JeQZY9g/fREnt93lt4Bp9FxTOtkdTv7ylr4q7wMgmUGgVclRtn409+s4LocO995s4BvvHbcK69lp0vz480lPPb8IbLs4bz6tyuYYo/w+HUDjbybLuOvV2XS2jPIq4flQOGr9fSuMsJDrNwv+5IYIiI0iKc+ncvfXDeFlw5UcdvPdnn0puVQP3sP/73lNPcsSuXlzy0nOTrMY9cLZFK4L2NJZixzU6P57YflfrN9pjfVtvfyzvFaPrE4XaZ+GchqUXzt5um88OhSuvodfOznu3nirQLae8dvf5PeASf/veU0t/zXTkobuviv++fz5D1zZf8RDxq1cCulfquUalBKnfRGIF+hlOKzq7Ioa+rmg6J6o+OYzrO7K3BpzV/lZRgdRQArp8az6fE1fGpJOs/trWDNk9v47y2nr6mA9ww4eG5vBet/tJ0fby5h7bQENv39Gj42P0VuRHuYGm0PaqXUaqALeE5rPdudJ83NzdX5+fnjEM9YDqeLtT/aTmx4KG/87Qp5MbqptXuAvH/fyo0zE/np/QuMjiMucrK6nZ9sLmFLcQPhIVZunp3MnQsmsSQzdtSzHp0uzeHKVt49XssbR6tp6xlkQXoMX795OktlP+1ropQ6pLXOdeexox43obXeqZTKuNZQZhRktfA3a7L55usn+LC0iVVT5cQWdzy7u5yeASdfWJttdBRxCbNTonnmM4sprOngd3sqeO9ELa8ePkdokIX5aTFMTYwgdeIEIkKDsFoUHb2D1HX0UVLfybGqdrr6HYQGWVg/I4FH8jJZNHmiDGq8bNQRN8Bw4X7nSiNupdRjwGMA6enpi86e9Y9l4/0OJ9c9uZ20iUPzY8WVdfQNkvfDreRNiedXDy0yOo5wQ9+gk12nm9hX1kz+2VbONnfTdtEe3xNCrGTEhbNwcgxLMuNYNz1BjhkbZ+M64naX1vop4CkYapWM1/MaLTTIyudWZ/HE24XsL2uWj4OjeH7vWTr7HHxxnYy2zcIWbOWGmYncMDPx/L/r7Bukd9CJ06WJCA0iUm4w+xSZVeKG+5ekEx8Rws+2lhodxaf1DDh4elcZa6fZZcGFyUXagkmItJEcHSZF2wdJ4XaDLdjKX6/K4sPSJg77wP4PvurF/ZW09gzKaFsID3NnOuBLwF5gmlLqnFLqUc/H8j0PLptMbHgIP95UYnQUn9TV7+AX28+wMjueRZNjjY4jhF8btXBrrT+ptU7WWgdrrVO11s94I5ivCQ8N4m+vm8KHpU3sKW0yOo7PeWZXOS3dA3JatxBeIK2SMXhw2WQmRdv4942ncGc2TqBo6R7gN7vKuGlWIvPSYoyOI4Tfk8I9BrZgK1++firHqtrYVCirKUf8cnspPQMO/vFGGW0L4Q1SuMfo7oWpZNnD+c+Np3DKHibUtvfyu71n+fiCVNkkXwgvkcI9RkFWC1+5YRqnG7p4TXYO5KebT6O15vHrpxodRYiAIYX7KtwyO4l5aTE8ufEU3f0Oo+MY5mR1Oy8fquLh5RlyjqAQXiSF+ypYLIrv3j6Ths5+fr4tMBflaK353juFTJwQwt+tl9G2EN4khfsqLUyfyF0LUnh6VzmVzT1Gx/G690/WcaC8ha/cmEN0mKysE8KbpHBfg6/dMp0gq+L77xUaHcWr+gad/OC9IqYnRfKJ3DSj4wgRcKRwX4PEKBtfWJvNxoJ6PjwdOItyfrn9DOdae/nObTMJkrMkhfA6edddo0dXZjI5bgL//OZJ+gb9/2Dh0oYufrn9DB+bP4kV2fFGxxEiIEnhvka2YCvfv3MO5U3d/GzraaPjeJTWmm++fgJbsIVv3zrT6DhCBCwp3ONg5dR47l6Yyq93lFFU22F0HI95Jf8cB8pb+OaGGdgjQ42OI0TAksI9Tr596wyiw4L5+msn/HJFZUNHH99/r4glGbHcJzckhTCUFO5xMjE8hO/cPpNjVW08tbPM6DjjSmvNV189PjSb5K45WCxyvqAQRpLCPY7umDeJW2Yn8ePNpzhZ3W50nHHz4oFKtp9q5Bu3TCc7IcLoOEIEPCnc40gpxQ8+PofY8BAe/+NRv5hlUt7Uzb++U8SqqfF8enmG0XGEEEjhHncTw0P4z3vnUdrQxb+9V2R0nGsy4HDx+B+PEmxVPHnPPGmRCOEjpHB7wKqpdh7Jy+R3e8/y7vFao+NctR+8V8SxqjZ+ePdckqJtRscRQgyTwu0hX79lOgvSY/inPx3jdH2n0XHG7O1jNfzvngoeyctkw5xko+MIIS4ghdtDQoIs/PKBRUwICeJzzx+io2/Q6EhuK6nv5OuvHmfR5Il8Y8N0o+MIIS4ihduDkqJt/PxTCzjb0sPjfziKw+kyOtKoGjv7+atnDxIeGsT/fGoBwbIXiRA+R96VHrY0K44nbp/J1uIGvvNWgU8fMtw74OSzz+XT0j3AMw8vJjk6zOhIQohLCDI6QCB4aHkG1W19/GrHGVJiwvjC2myjI33EoNPFl/5whOPn2vjVg4uYkxptdCQhxGVI4faSr940jbr2Xp7ceIpIW5BPzYl2ujT/8PIxNhfW8y93zOKmWUlGRxJCXIEUbi+xWBT/cc88uvqdfOfNAlwuzWfyMo2OhdOl+eqfjvP2sRq+fst0Hl6RYXQkIcQopMftRSFBFn7xwEJumJnIE28X8pudZYb2vPsGnXzxxcO8evgcj18/lc+vmWJYFiGE+6Rwe1lIkIWff2oht85J5vvvFfHEWwWGzDbp6Bvk4d8e4P2TdXz71hk8fn2O1zMIIa6OtEoMEBJk4WefXEDKxDCe2llGVWsvP7lvPtETvHPobkl9J59//hBVrT381/3z+dj8FK9cVwgxPmTEbRCLRfHNDTP41ztns7OkkQ3/vYtDZ1s9ek2tNa8dPsedP99NR5+DFx5dKkVbCBOSwm2wB5dN5pXPL0cpuO/Xe/m394roGXCM+3XqO/r46+fy+YeXjzFrUhTvfmklS7Pixv06QgjPc6twK6VuVkqdUkqVKqW+7ulQgWZB+kTe/dIq7lmYyq93lnHDj3fy+pFz43KSTne/g59+UMK6/9zOh6VNfPvWGfzhseUkRsmmUUKYlRptVoNSygqUADcA54CDwCe11oWX+5nc3Fydn58/njkDxoHyFr77VgFFtR1k2cN5JC+TO+ZPIso2tv53bXsvv99XyYsHKmnpHuCW2Ul87ebpZMSHeyi5EOJaKKUOaa1z3XqsG4V7OfCE1vqm4e+/AaC1/rfL/YwU7mvjcmk2Fdbxs62lFNR0YAu2sCbHzsqpdnInTyQzPhxbsPUvfqar30FhTQeHK1v5oLCeQ5VD/fLrZyTy+TVTWDR5ohF/FCGEm8ZSuN2ZVZICVF3w/Tlg6dUEE+6xWBQ3z07mpllJnKhu55X8c2wtbmBjQT0ASkFceAi2YCtKQVv3IJ39/9cXn5kcxZfWTeXuhamkx00w6o8hhPAQdwr3pY49+cgwXSn1GPAYQHp6+jXGEjB0FNrc1BjmpsbwPa2paO7hZHU7ZY3d1HX00T/oxKU1MRNCsEeGMiM5ktkp0SRESv9aCH/mTuE+B6Rd8H0qUHPxg7TWTwFPwVCrZFzSifOUUmTGh5MpPWohAp47s0oOAlOVUplKqRDgfuAtz8YSQghxOaOOuLXWDqXUF4GNgBX4rda6wOPJhBBCXJJbS9611u8B73k4ixBCCDfIykkhhDAZKdxCCGEyUriFEMJkpHALIYTJSOEWQgiTGXWvkqt6UqUagbNX+ePxQNM4xhkvkmtsJNfYSK6x8cdck7XWdnce6JHCfS2UUvnubrTiTZJrbCTX2EiusQn0XNIqEUIIk5HCLYQQJuOLhfspowNchuQaG8k1NpJrbAI6l8/1uIUQQlyZL464hRBCXIFPF26l1D8qpbRSKt7oLABKqf+nlDqulDqqlNqklJpkdCYApdSTSqni4WyvK6VijM4EoJS6VylVoJRyKaUMnQHgqwdeK6V+q5RqUEqdNDrLhZRSaUqpbUqpouH/hl82OhOAUsqmlDqglDo2nOtfjM40QillVUodUUq94+lr+WzhVkqlMXRAcaXRWS7wpNZ6rtZ6PvAO8B2jAw3bDMzWWs9l6GDnbxicZ8RJ4C5gp5Ehhg+8/jlwCzAT+KRSaqaRmS7wv8DNRoe4BAfwFa31DGAZ8AUf+TvrB9ZprecB84GblVLLDM404stAkTcu5LOFG/gJ8FUucUyaUbTWHRd8G46PZNNab9Jajxw6uY+hU4oMp7Uu0lqfMjoHsAQo1VqXaa0HgD8AHzM4EwBa651Ai9E5Lqa1rtVaHx7+upOhgpRibCrQQ7qGvw0e/mX4+1AplQrcCjztjev5ZOFWSt0BVGutjxmd5WJKqe8rpaqAB/CdEfeFHgHeNzqEj7nUgdeGFyGzUEplAAuA/cYmGTLckjgKNACbtda+kOunDA00Xd64mFsHKXiCUuoDIOkSv/Ut4JvAjd5NNORKubTWb2qtvwV8Syn1DeCLwHd9IdfwY77F0Efc33sjk7u5fIBbB16Lj1JKRQCvAo9f9InTMFprJzB/+F7O60qp2Vprw+4RKKVuAxq01oeUUtd545qGFW6t9fWX+vdKqTlAJnBMKQVDH/sPK6WWaK3rjMp1CS8C7+Klwj1aLqXUw8BtwHrtxTmeY/j7MpJbB16Lv6SUCmaoaP9ea/2a0XkuprVuU0ptZ+gegZE3d/OAO5RSGwAbEKWUekFr/aCnLuhzrRKt9QmtdYLWOkNrncHQm26hN4r2aJRSUy/49g6g2KgsF1JK3Qx8DbhDa91jdB4fJAdej5EaGjU9AxRprX9sdJ4RSin7yKwppVQYcD0Gvw+11t/QWqcO16v7ga2eLNrgg4Xbx/1QKXVSKXWcoVaOT0yRAv4HiAQ2D09V/JXRgQCUUh9XSp0DlgPvKqU2GpFj+MbtyIHXRcDLvnLgtVLqJWAvME0pdU4p9ajRmYblAQ8B64ZfU0eHR5RGSwa2Db8HDzLU4/b49DtfIysnhRDCZGTELYQQJiOFWwghTEYKtxBCmIwUbiGEMBkp3EIIYTJSuIUQwmSkcAshhMlI4RZCCJP5/5YXPWQ2mNtoAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "x = np.linspace(-4, 4, 1000)\n",
    "plt.plot(x, np_two_well(x))\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "BFGS Optimization - Minimization\n",
    "====\n",
    "Broyden–Fletcher–Goldfarb–Shanno"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "BFGS\n",
    "====\n",
    "\n",
    "**Type:** Minimization\n",
    "\n",
    "**Discrete/Continuous:** Continuous\n",
    "\n",
    "**Dimensions:** N\n",
    "\n",
    "**Derivative:** optional\n",
    "\n",
    "**Non-Convex:** not recommended\n",
    "\n",
    "**Python:** `minimize` if no constraints or bounds are given"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "Nomenclature\n",
    "----\n",
    "\n",
    "In optimization, you have a function called the **objective function**. That's what you're minimizing. It always returns a single value. It is sometimes called the **fit**, the **error**, the **residual**, or the **penalty**."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAD8CAYAAACMwORRAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl4VOXd//H3NzsJSYAkrNkIhH03IJugYhXQggsqKIK7aLX2sa11aWtra1v1+bW2VVRwQVBRxAVUFLWCCwKSsIc1BEJCFhICISF7cv/+yOAzjQmZwEzOLN/XdXE5OXNm5uNJ+HBylvsWYwxKKaW8i5/VAZRSSjmflrtSSnkhLXellPJCWu5KKeWFtNyVUsoLabkrpZQX0nJXSikvpOWulFJeSMtdKaW8UIBVHxwdHW0SExOt+nillPJIaWlpRcaYmJbWs6zcExMTSU1NterjlVLKI4lIliPr6WEZpZTyQlruSinlhbTclVLKC2m5K6WUF9JyV0opL9RiuYvIKyJyVER2NvO8iMi/RCRDRLaLyAjnx1RKKdUajuy5LwImn+H5KUCy7c+dwPPnHksppdS5aLHcjTFfA8VnWGU6sNg02AB0EJFuzgrY2NbsEzz56R5Xvb1SSrmMMYYnPt5Fem6Jyz/LGcfcewDZdl/n2Jb9iIjcKSKpIpJaWFh4Vh+2I+cEz689wM4jrt84SinlTBsyi1n4zUH25pe6/LOcUe7SxLImZ902xiwwxqQYY1JiYlq8e7ZJ04b1IDjAj7c2HT6r1yullFXe3nSY8JAApgxy2cGNHzij3HOAOLuvY4FcJ7xvkyLbBTJ1cDdWbM2lorrOVR+jlFJOVVJewyc785k+rDvtgvxd/nnOKPeVwBzbVTOjgRJjTJ4T3rdZ16XEUVpZyyc7XfoxSinlNCu2HaGqtp6ZI+Pb5PMcuRRyKbAe6CsiOSJym4jME5F5tlVWAZlABrAQuMdlaW1GJ3UiMSqUtzdlt7yyUkq5gbc3ZTOgWwSDekS2yee1OCqkMWZWC88b4GdOS+QAEeHalDieXr2XzMIykmLat+XHK6VUq+w8UkJ67kkenz6wzT7TY+9QnXFeLP5+wrLUHKujKKXUGb216TDBAX5MH9rkhYQu4bHl3iUihIv6xrA8LYeaunqr4yilVJMqqutYsTWXKYO6Ehka2Gaf67HlDnD9yHiKyqr4cs9Rq6MopVSTVu3Io7SylutGxrW8shN5dLlf1DeGLhHBLP1er3lXSrmnpd8fpmd0GGOSotr0cz263AP8/bguJY6v9hVy5ESF1XGUUuq/7C8oJTXrODNHxiHS1P2eruPR5Q4N17wDelmkUsrtLP0+m0B/4ZrzYtv8sz2+3OM6hTIhOYZ3UrOp1ROrSik3UVlTx3tbcrh0YFei2we3+ed7fLkDzBoVR15JJV/tO7vByJRSytlWp+dzoryGWW10R2pjXlHuk/p3Ibp9MEu/10MzSin3sPT7w8R3CmVsr7Y9kXqaV5R7oL8f16XE8uWeAvJLKq2Oo5TycZmFZWzILGbmqDj8/Nr2ROppXlHuADNHxlNv9MSqUsp6S78/TICfMMOCE6mneU25x0eFMqFPDG9tOqwnVpVSlqmsqeOdtBwuG9iVzuEhluXwmnIHuPH8ePJKKlmzV0+sKqWs8cnOPE6U13DD+dacSD3Nq8p9Ur/OdIkI5o2NWVZHUUr5qDc2WHNHamNeVe4B/n5cPzKer/YVkl1cbnUcpZSP2ZN/ktSs49wwKt6yE6mneVW5Aw23+YKON6OUanNvbjxMkL+fJXekNuZ15d69Qzsu7teFZak5VNfqiVWlVNsor67l/c1HmDq4K53CgqyO433lDnDj6IahgFen51sdRSnlI1ZuzaW0qpYbzk+wOgrgpeU+ITmGuE7teH2DnlhVSrmeMYbF67Po2yWckYkdrY4DeGm5+/sJN56fwMaDxewrKLU6jlLKy23JPsGuvJPMHpPQ5kP7Nscryx0ahgIOCvDTvXellMu9vj6L9sEBXDW87eZIbYnXlnunsCCuGNyN9zYfoayq1uo4SikvVXyqmo+253H1iB60Dw6wOs4PvLbcAWaPSaCsqpYPthyxOopSykstS82muq6e2aPd40TqaV5d7sPjOjCwewSvb8jCGGN1HKWUl6mrN7yxMYvze3aiT5dwq+P8F68udxHhptEJ7MlvmMdQKaWc6et9hWQXV3DTGPfaawcvL3eAacO6ExESwGvfHbI6ilLKy7y2/hAx4cFcOqCr1VF+xOvLPTQogOtS4vh0Zz4FJ3UiD6WUcxwqOsXavYXceH48QQHuV6Xul8gFbhqTQJ0xvLFRx5tRSjnH4vVZBPqL5UP7Nscnyj0hKoyL+nbmzY2HdbwZpdQ5O1VVyztp2UwZ1M3SCTnOxCfKHWDOmASKyqr4ZGee1VGUUh7u/S1HKK2sZe5Y9zuReppD5S4ik0Vkr4hkiMhDTTwfLyJrRGSLiGwXkanOj3puJiTH0DM6TE+sKqXOScM4MocY1COCEfHuMY5MU1osdxHxB54DpgADgFkiMqDRar8FlhljhgMzgfnODnqu/PwaLovcfPgEO3JKrI6jlPJQ6zOPsa+gjDljEt1mHJmmOLLnPgrIMMZkGmOqgbeA6Y3WMUCE7XEkkOu8iM4zIyWW0CB/Xv3uoNVRlFIeatG6Q3QMDWTa0O5WRzkjR8q9B5Bt93WObZm9PwCzRSQHWAXc55R0ThYREsiM82L5aFsehaVVVsdRSnmY7OJyvthdwKxR8YQE+lsd54wcKfemfu9ofC//LGCRMSYWmAosEZEfvbeI3CkiqSKSWlhY2Pq0TjB3bCLVdfW8qZdFKqVaafH6Qw13vrvhHamNOVLuOUCc3dex/Piwy23AMgBjzHogBIhu/EbGmAXGmBRjTEpMTMzZJT5HvWLac2HfGF7fmKWXRSqlHHaqqpa3NmUzZVBXukW2szpOixwp901Asoj0FJEgGk6Yrmy0zmFgEoCI9Keh3K3ZNXfAzWMTKSyt4uMdbnlqQCnlht7bnENpZS23jEu0OopDWix3Y0wtcC+wGthNw1Ux6SLyuIhMs632S+AOEdkGLAVuNm48DOOE5BiSYsJ4dd0hHS1SKdWi+nrDq98dYkhspFtf/mjPoZHljTGraDhRar/s93aPdwHjnBvNdfz8hFvGJvK7FelsPnyc8xI6WR1JKeXGvskoIrPwFP+4fqhbX/5oz2fuUG3s6hGxRIQE8Mq6Q1ZHUUq5uVe+PUhMeDCXD3bvyx/t+Wy5hwUHMGtUPJ/syCO7uNzqOEopN7W/oJSv9hUyZ3SCW47+2BzPSeoCc8c23GGmQxIopZrzyrqDBAf4caObTaPXEp8u9+4d2jF1cDfe3pRNaWWN1XGUUm7mWFkV724+wtUjYukUFmR1nFbx6XIHuG18T0qralmWmmN1FKWUm3nDNkz4beMTrY7Saj5f7sPiOpCS0JFX1x2krl4vi1RKNaiqrWPx+iwu7BtD787uNfm1I3y+3KFh7z3neAWfpedbHUUp5SZWbs2lqKyK28b3tDrKWdFyBy4d2JW4Tu146VsdLVIp1TBm+8vfHqRvl3DG9/7RSCoeQcsd8PcTbhvXk7Ss46RlFVsdRyllsa/3F7Env5Q7JiR5zE1LjWm521ybEkdku0AWfJ1pdRSllMUWfp1Jl4hgtx+z/Uy03G3CggO4aXQCn+0q4GDRKavjKKUskp5bwrcZRdwyrqdH3bTUmOcmd4E5YxMI9PPj5W91710pX7Xw60zCgvyZNSre6ijnRMvdTufwEK4a3oN3UnM4VqYzNSnla3JPVPDh9jxmjoonsl2g1XHOiZZ7I7df0JOq2nqWbMiyOopSqo29uq7hijlPGbP9TLTcG0nuEs6kfp1ZvD6Liuo6q+MopdpISXkNb248zBVDuhHbMdTqOOdMy70J8y7sRfGpapalZre8slLKK7y+MYtT1XXcNaGX1VGcQsu9CSMTO3FeQkcWfpNJbZ3Os6qUt6usqePVdQeZ0CeGAd0jrI7jFFruzZg3sRc5xyv4eEee1VGUUi62PC2HorJq5k1MsjqK02i5N2NSv84kd27PC19l6jyrSnmxunrDwm8yGRobyZikKKvjOI2WezP8/IQ7JySxO+8kX+0rtDqOUspFPtmZR9axcuZN7OWxQw00Rcv9DKYP60G3yBCeX3vA6ihKKRcwxvDCVwfoGR3GpQO7Wh3HqbTczyAowI/bL0hi48Fi0rKOWx1HKeVkX+8vYueRk8ybmIS/n/fstYOWe4tmjYqjY2ggz6/NsDqKUsrJ5q/JoFtkCFcNj7U6itNpubcgNCiAW8f15IvdR9mdd9LqOEopJ0k9VMzGg8XccUGSRw8Q1hzv+z9ygTljEmkfHKDH3pXyIvPXHqBTWBAzR8VZHcUltNwdEBkayI2j4/loey6HdDhgpTxeem4JX+45yi1jEwkNCrA6jktouTvotvE9CfD344WvdO9dKU/3/NoDtA8OYM6YRKujuIyWu4M6h4cwc2Qc727O4ciJCqvjKKXOUsbRMj7ekcfs0QlEhnr2sL5nouXeCndNbBhQ6EXde1fKY81fk0FwgB+3X9DT6igupeXeCj06tOOaEbG8tSmboycrrY6jlGqlrGOnWLEtl9nnJxDdPtjqOC6l5d5K91zYm7p6w4s6kbZSHmf+mgP424YW8XYOlbuITBaRvSKSISIPNbPOdSKyS0TSReRN58Z0H/FRoUwf1p03NmZRpFPxKeUxco6X8+7mHGaNjKNzRIjVcVyuxXIXEX/gOWAKMACYJSIDGq2TDDwMjDPGDAR+4YKsbuNnF/Wmqraehd/o3rtSnuKFrw4g8n/nzrydI3vuo4AMY0ymMaYaeAuY3midO4DnjDHHAYwxR50b0730imnPFUO6s2R9FsWnqq2Oo5RqQV5JBcs25TDjvDi6d2hndZw24Ui59wDs55vLsS2z1wfoIyLrRGSDiExu6o1E5E4RSRWR1MJCzx5G9+cX96aipo4FeuxdKbc3f80B6o3hZxf5xl47OFbuTQ2V1nj2igAgGbgQmAW8JCIdfvQiYxYYY1KMMSkxMTGtzepWkruEc8WQ7ixef4hjeuxdKbeVe6KCtzdlc21KnFdMfO0oR8o9B7AffCEWyG1inRXGmBpjzEFgLw1l79Xun9Sw977wm4NWR1FKNWP+2gwMvrXXDo6V+yYgWUR6ikgQMBNY2WidD4CLAEQkmobDNF5/vKJ353B+qnvvSrktX91rBwfK3RhTC9wLrAZ2A8uMMeki8riITLOttho4JiK7gDXAr40xx1wV2p383Lb3vkCvnFHK7Ty3pmEehp9d1NviJG3PoeHQjDGrgFWNlv3e7rEBHrD98Sm9O4czbWh3Fn+Xxe3jk4gJ9+673pTyFDnHy1mW2rDX3sNHrpCxp3eoOsH9k5Kpqq3TESOVciP//k8GIsJ9F/veXjtouTtFUkx7rh4Ry5INWeSX6JgzSlntYNEplm/O4YZR8XSL9L29dtByd5r7JyVTX294ds1+q6Mo5fP++cU+Av2Fe3zsChl7Wu5OEtcplOtHxvH2pmyyi8utjqOUz9pXUMqKbbnMHZtI53DvH0OmOVruTnTvxb0REf79pe69K2WVZ77YR1hQAPMm+O5eO2i5O1W3yHbMPj+B5Wk5HCgsszqOUj5n55ESVu3I59ZxiXQMC7I6jqW03J3snot6ERLoz98/22d1FKV8zlOr99IhNJDbfWC89pZouTtZdPtgbr8giY935LEjp8TqOEr5jPUHjvH1vkJ+dmFvIkK8d25UR2m5u8AdF/SkY2ggT63eY3UUpXyCMYanVu+ha0QIN41JsDqOW9Byd4HwkEDuubA33+wv4rsDRVbHUcrrfb6rgC2HT3D/JcmEBPpbHcctaLm7yE1jEugWGcJTn+6lYXQGpZQr1NUb/vezvSRFh3HtebFWx3EbWu4uEhLozy8uSWZr9glWp+dbHUcpr/Xe5hz2FZTxwKV9CPDXSjtNt4QLXTMilt6d2/PUp3upqau3Oo5SXqeypo6/f76PobGRXD64m9Vx3IqWuwsF+Pvx0OR+ZBad4q1N2S2/QCnVKq+uO0ReSSUPTemPSFOTxvkuLXcXm9S/M6MSO/HPL/ZzqqrW6jhKeY3jp6qZvzaDi/t1ZkyvKKvjuB0tdxcTER6e2o+isioW6oQeSjnNc2syOFVVy28m97M6ilvScm8Dw+M7MnVwVxZ8ncnRUh0SWKlzlV1czuL1WVwzIpa+XcOtjuOWtNzbyIOX9aOmrp5/fK7DEih1rp78dA9+fvDApX2sjuK2tNzbSGJ0GDeNTuTtTdnsyT9pdRylPFZa1nE+2p7HnRck+exEHI7Qcm9DP5/Um/CQQJ74eLfe2KTUWTDG8OePdxETHsxdE317SN+WaLm3oQ6hQfx8UjLf7C9i7b5Cq+Mo5XE+3pHHlsMn+PWlfQkLDrA6jlvTcm9jN41OIDEqlL98vJtavbFJKYdV1tTx5Kd76Nc1nGt0mIEWabm3saAAPx6a0p/9R8tY+v1hq+Mo5TFeXXeI7OIKfnv5APz99Iallmi5W+CygV0YkxTF3z/fx4nyaqvjKOX2jpZW8uyX+/nJgC6MT462Oo5H0HK3gIjw+58OoKSihme+0PlWlWrJ05/upbqunken9rc6isfQcrdI/24RzBoVz5INWewvKLU6jlJua1v2Cd5Jy+HW8T1JjA6zOo7H0HK30AM/6UNYkD+Pf7RLL41UqgnGGP74YTrR7YO596LeVsfxKFruFopqH8wvLunDN/uL+GL3UavjKOV2VmzNZfPhEzw4uS/hOi9qq2i5W+ymMQkkd27P4x+lU1lTZ3UcpdxGaWUNT6zazdDYSGaM0EsfW0vL3WKB/n78cdpAsosrePErHTVSqdP+9Z/9FJVV8fj0QfjppY+t5lC5i8hkEdkrIhki8tAZ1pshIkZEUpwX0fuN7R3N5UO6MX9tBtnF5VbHUcpy+wtKeXXdIa5PiWNoXAer43ikFstdRPyB54ApwABglogMaGK9cODnwEZnh/QFv728P34i/PnjXVZHUcpSxhj+8GE6oUH+/PqyvlbH8ViO7LmPAjKMMZnGmGrgLWB6E+v9CXgK0AHLz0K3yHbcN6k3q9MLWLtXT64q37VqRz7rMo7x68v6EtU+2Oo4HsuRcu8B2E8AmmNb9gMRGQ7EGWM+cmI2n3P7+CSSYsJ4bKWeXFW+qbSyhsc/Smdg9whuOD/B6jgezZFyb+pMxg8XZYuIH/AP4JctvpHInSKSKiKphYU6KmJjQQF+/Hn6ILKOlTN/TYbVcZRqc3//fB9HS6t44qrBOn7MOXKk3HOAOLuvY4Fcu6/DgUHAWhE5BIwGVjZ1UtUYs8AYk2KMSYmJiTn71F5sbO9orhzWnRe+yuRAYZnVcZRqMzuPlPDad4e48fx4hulJ1HPmSLlvApJFpKeIBAEzgZWnnzTGlBhjoo0xicaYRGADMM0Yk+qSxD7g0csHEBzox+8+2Kl3riqfUF9v+O0HO+kUFsSvL9MJr52hxXI3xtQC9wKrgd3AMmNMuog8LiLTXB3QF8WEB/Pg5H58d+AYK7bmtvwCpTzcm98fZmv2CX57+QAi2+mdqM7g0FQmxphVwKpGy37fzLoXnnssdcOoeN5Ny+FPH+1iYp8YOoYFWR1JKZcoOFnJk5/sYVzvKKYP6251HK+hd6i6KX8/4W/XDKakouEWbKW81WMr0qmuq+eJKwcjoidRnUXL3Y316xrBXROTWJ6Ww7qMIqvjKOV0q9Pz+TQ9n/svSdbhfJ1My93N3XdxMolRoTzy/g699l15ldLKGh5bkU6/ruHccUGS1XG8jpa7mwsJ9OcvVw8m61g5//hin9VxlHKaJz/dQ0FpJX+7ZgiB/lpFzqZb1AOM7RXNzJFxLPw6k23ZJ6yOo9Q5W3/gGK9vOMwtY3vqNe0uouXuIR65vD+dw0N4cPl2qmvrrY6j1Fkrr67lN+9uJyEqVAcGcyEtdw8RERLIE1cNYm9BKc/q0ATKg/2/z/ZxuLicv109hHZB/lbH8Vpa7h5kUv8uXDW8B/PXZLAr96TVcZRqtbSs47yy7iCzR8czpleU1XG8mpa7h/n9FQPoEBrIr97ZpodnlEepqK7j1+9so3tkOx6a0t/qOF5Py93DdAwL4omrBrMr7yT//nK/1XGUcthTq/eQWXSKp2YMoX2wQzfHq3Og5e6BLhvYlatH9GD+2gNs1atnlAf47kARr647xNwxCYzrHW11HJ+g5e6hHvvpQDqHB/PLZVv15ibl1kora/j1O9tJjArlN1N0xMe2ouXuoSLbBfLkNUM4UHiKp1fvtTqOUs3680e7ySup4P9dN5TQID0c01a03D3YhD4x3DQ6gZe/Pci3+3XsGeV+Pt2Zz9up2dw1sRfnJXSyOo5P0XL3cI9M7U+vmDB++c5WTpRXWx1HqR8cPVnJw+9tZ1CPCP7nkj5Wx/E5Wu4erl2QP/+cOZziU9U88v4OnblJuYX6esOvlm+noqaOZ64fTlCAVk1b0y3uBQb1iOSBn/Rl1Y58lqflWB1HKV5bf4iv9xXy6OUD6N25vdVxfJKWu5e4c0IS5/fsxGMr08nUibWVhXblnuSvn+zh4n6dmX1+vNVxfJaWu5fw9xOemTmMoAA/7lu6hapavTxStb3y6lruXbqZDu0CeXrGEJ1ZyUJa7l6kW2Q7np4xlPTck/ztkz1Wx1E+6LEV6RwsOsUzM4cR1T7Y6jg+Tcvdy/xkQBduHpvIq+sO8cWuAqvjKB+yYusR3knL4d6LejO2l96FajUtdy/08NR+DOgWwa+WbyPneLnVcZQPOFBYxiPv7SAloSP3T0q2Oo5Cy90rBQf489yNI6itM9z75hYdPVK5VEV1Hfe8vpmgAD/+NWs4ATplnlvQ74KX6hkdxtMzhrA1+wR/WbXb6jjKi/1uxU72HS3lmZnD6d6hndVxlI2WuxebMrgbt4xLZNF3h/h4e57VcZQXWpaazfK0HO67qDcT+8RYHUfZ0XL3cg9P6c/w+A48uHwbGUdLrY6jvMjOIyX87oOdjO0Vxf06vIDb0XL3ckEBfsy/cQTtgvy5c0kapZU1VkdSXqD4VDV3LUmjU1gQ/5o1HH8/vZ7d3Wi5+4Buke149oYRZB0r54Fl26iv1/Fn1NmrravnvqWbKSyr4oXZ5xGt17O7JS13HzE6KYpHp/bn810FPLsmw+o4yoM9vXov6zKO8ecrBzE0roPVcVQztNx9yC3jErlqeA/+8cU+PkvPtzqO8kAfbDnCi19nMnt0PNelxFkdR52BlrsPERH+evVghvSI5Bdvb2V33kmrIykPsuXwcR58dzujkzrx2E8HWh1HtcChcheRySKyV0QyROShJp5/QER2ich2EfmPiCQ4P6pyhpBAfxbMSSE8JIDbX0ulqKzK6kjKA+SVVHDnkjS6RoTw/I3nEag3Krm9Fr9DIuIPPAdMAQYAs0RkQKPVtgApxpghwHLgKWcHVc7TJSKEhXNSOHaqinlL0nQESXVG5dW13LE4lYrqOl6am0LHsCCrIykHOPLP7yggwxiTaYypBt4CptuvYIxZY4w5PYjJBiDWuTGVsw2J7cD/XjuU1KzjPLh8u87gpJpUV2/4+dKt7Mo9yb9mDaNPl3CrIykHOVLuPYBsu69zbMuacxvwSVNPiMidIpIqIqmFhYWOp1QuccWQ7jw4uS8rtuby98/3WR1HuaE/fbSLL3YX8IdpA7m4Xxer46hWCHBgnabuTmhyN09EZgMpwMSmnjfGLAAWAKSkpOiuohu4e2IvDh8r599fZhDXMZTrRuoVEKrBK98eZNF3h7htfE/mjEm0Oo5qJUfKPQew/xsfC+Q2XklELgEeBSYaY/QsnYcQEf505SCOnKjgkfd30DkimAv7drY6lrLYJzvy+NPHu7hsYBcemdrf6jjqLDhyWGYTkCwiPUUkCJgJrLRfQUSGAy8C04wxR50fU7lSoH/DEAV9u4Zz9+ub2XL4uNWRlIW+O1DE/W9tZUR8R565XocW8FQtlrsxpha4F1gN7AaWGWPSReRxEZlmW+1poD3wjohsFZGVzbydclPhIYEsumUUnSOCuXXRJjKO6iTbvmjnkRLuXJxGQlQoL89NoV2Qv9WR1FkSq66SSElJMampqZZ8tmpe1rFTXPP8dwQH+PPOvDE6PrcPafjeryfIX3j3nrF0i9TvvTsSkTRjTEpL6+mdCOq/JESFseiWUZysqGH2SxspLNXTJ74g90QFNyzcSF19Pa/dOkqL3QtouasfGdQjkldvGUleSSU3vbyRE+XVVkdSLnS0tJIbX9rIyYoaltx2Psl6LbtX0HJXTUpJ7MTCOSlkFp5i7ivf6zjwXur4qWrmvPw9+SWVLLp1JIN6RFodSTmJlrtq1vjkaObfOIL03JPM0YL3OsdPVXPjSxvJLDrFS3NTOC+hk9WRlBNpuaszumRAF569YQQ7ckqY88r3nNSC9wrFp6q54aWNZBSW8dKcFMb1jrY6knIyLXfVosmDuvLcjbaCf1kL3tMVn95jtxX7BJ3Y2itpuSuHXDawq+0QTQk3LNzAMR0q2CMVnKzk+hfXk1lYxkItdq+m5a4cdunAriyYk8L+gjKuX7CBgpOVVkdSrZBdXM61L6wn90QFi24ZpcXu5bTcVatc1Lczr906ivySSma88B2Hj5W3/CJluYyjpVz7wnpKKmp4447RjOkVZXUk5WJa7qrVRidF8cbt51NaWcvVz69jR06J1ZHUGWw6VMw1z6+ntt7w9l2jGaaTWvsELXd1VobGdWD5vLEEB/hz/YL1fLVPx+d3R5/uzGf2SxuJCgvi/XvG0q9rhNWRVBvRcldnrXfn9rx3z1gSosK4bdEmlm3KbvlFqk0YY1i07iB3v5HGgO4RLL97LHGdQq2OpdqQlrs6J10iQlh2V8Mx3Aff3c5fVu2mrl7nYbFSTV09v1uxkz98uItJ/brw5u2j6aTznvocLXd1zsJDAnnl5pHcNDqBBV9ncteSNMqqaq2O5ZNKKmq4ddEmXt9wmLsmJPHiTefpsL0+SstdOUWgvx9/unIQf5w2kC/3FHD1/HVkFuqY8G1pb34p05/9lg2Zx3hzTgOUAAAMjUlEQVRqxhAentpfJ9rwYVruyqnmjk1k8a3nU1haxfRn1/H5rgKrI/mED7flcuVz6zhVXcebd4zmuhSdC9fXabkrpxufHM2H940nITqUOxan8vTqPdTW1VsdyytV1dbxxw/TuW/pFgZ2j+Dj+8YzMlEHAFNa7spFYjuGsnzeWK5PieO5NQeYuWADR05UWB3LqxwqOsWM59fz6rpD3Dw2kTfvGE3niBCrYyk3oeWuXCYk0J8nZwzhnzOHsSe/lKn//IZPd+ZZHcvjGWP4YMsRrvj3txwuLufFm87jD9MGEhSgf53V/9GfBuVy04f14OOfjychKpR5r2/mf97eSkmFjix5No6VVXHPG5v5xdtb6dc1nFX3X8BlA7taHUu5oQCrAyjfkBAVxrt3j+XZLzN4dk0G6w8c48kZQ5iog1c5bHV6Po++v4OTFbX8ZnI/7pyQpFfDqGbpnrtqM4H+fvzPT/rw/j1jaR8SwNxXvuf+t7ZQpMMHn1F+SSXzlqRx15I0YsJDWHnfOO6+sJcWuzojMcaauwlTUlJMamqqJZ+trFdZU8f8tQd4fm0GoUEBPDylH9elxOGnhfWD2rp63th4mKdX76Wmrp77L0nmjguSCPTXfTJfJiJpxpiUFtfTcldW2l9QyiPv72DToeMM7hHJYz8dQIpeyse6jCIe/3AXewtKGd87mieuGkRCVJjVsZQb0HJXHsMYw4qtufztkz3kn6zkiiHd+NWlfUmM9r0y219QytOr9/LZrgJiO7bj0an9mTyoKyL6G41qoOWuPE55dS0vrD3Awm8OUl1Xz3Upcdw/KZmukd5/7XZ2cTnPfLGf97fkEBoUwLyJSdx+QRIhgToujPpvWu7KYx0treS5LzN48/vDiAgzzovlrglJXnlYIuNoGS98dYAPthzBz0+YOyaBuy/sraM4qmZpuSuPl11czvy1B3g3LYfa+nqmDu7GLeN6MiK+g0cfpjDGsPFgMYvWHWL1rnyCA/yYOTKeuyYm0S2yndXxlJvTclde4+jJSl7+9iBvbjxMaVUtA7tHMGdMApcP6U77YM+5VeNkZQ0rt+ayZH0WewtKiWwXyE2jE7hlXCJR7YOtjqc8hJa78jqnqmp5f8sRFq8/xL6CMtoF+jN5UFeuGt6DMb2i3PISweraer7NKOS9zUf4bFcB1bX1DOgWwc1jE/np0O461rpqNaeWu4hMBv4J+AMvGWP+1uj5YGAxcB5wDLjeGHPoTO+p5a7OljGGtKzjvLflCB9ty+VkZS2R7QKZ1L8zlw7oyrjeUYSHBFqWr6S8hm8zilidns+aPUcpraqlY2gg04Z256oRsQyNjfTow0rKWk4rdxHxB/YBPwFygE3ALGPMLrt17gGGGGPmichM4CpjzPVnel8td+UMlTV1rN1byGe78vnP7qOUVNTg7ycMjY1kXO9oRiR0ZGhsB5eeoCwqq2Jb9gnSso6zLqOIHUdKqDfQKSyIS2z/4EzoE6MDeymncGa5jwH+YIy5zPb1wwDGmL/arbPats56EQkA8oEYc4Y313JXzlZTV0/qoYaCXXegiO05JT/M5xrXqR19u0TQu3N7esWEEdsxlK6RIXSNCHHo0MipqlryT1ZSUFJJzvEKDhSWcaCwjN15pT8MZRzgJwyP78DYXtGMT45meFwHAtzwUJHybI6WuyNno3oA9tPa5wDnN7eOMaZWREqAKKDIsbhKnbtAfz/G9IpiTK8ofkVfyqpq2XmkhG3ZJ9iWc4L9BWV8te8oNXX/vc8RHOBHeEgAYcEBBNnK2NBwvPxUVS2lVbVU1/73ZCNBAX4kRYcxPL4DN49NZFh8BwZ2jyA0yHNO8Crv5shPYlMHBxvvkTuyDiJyJ3AnQHx8vAMfrdTZax8cwOikKEYnRf2wrLaunuzjFeSeqCC/pJL8k5WcrKihtKqWsspaauv/r8QD/f1oHxxA+5AAOrQLomtkMF0iQujRoR2xHUN14C7l1hwp9xzAfkLGWCC3mXVybIdlIoHixm9kjFkALICGwzJnE1ipcxHg70fP6DB6+uDQBsq3OHJAcBOQLCI9RSQImAmsbLTOSmCu7fEM4MszHW9XSinlWi3uuduOod8LrKbhUshXjDHpIvI4kGqMWQm8DCwRkQwa9thnujK0UkqpM3Po7I8xZhWwqtGy39s9rgSudW40pZRSZ0uv01JKKS+k5a6UUl5Iy10ppbyQlrtSSnkhLXellPJClg35KyKFQNZZvjwa9xzaQHO1juZqPXfNprla51xyJRhjYlpaybJyPxcikurIwDltTXO1juZqPXfNprlapy1y6WEZpZTyQlruSinlhTy13BdYHaAZmqt1NFfruWs2zdU6Ls/lkcfclVJKnZmn7rkrpZQ6A7ctdxG5VkTSRaReRJo9qywik0Vkr4hkiMhDdst7ishGEdkvIm/bhit2Rq5OIvK57X0/F5GOTaxzkYhstftTKSJX2p5bJCIH7Z4b1la5bOvV2X32SrvlVm6vYSKy3vb93i4i19s959Tt1dzPi93zwbb//wzb9ki0e+5h2/K9InLZueQ4i1wPiMgu2/b5j4gk2D3X5Pe0jXLdLCKFdp9/u91zc23f9/0iMrfxa12c6x92mfaJyAm751y5vV4RkaMisrOZ50VE/mXLvV1ERtg959ztZYxxyz9Af6AvsBZIaWYdf+AAkAQEAduAAbbnlgEzbY9fAO52Uq6ngIdsjx8Cnmxh/U40DIMcavt6ETDDBdvLoVxAWTPLLdteQB8g2fa4O5AHdHD29jrTz4vdOvcAL9gezwTetj0eYFs/GOhpex//Nsx1kd3P0N2nc53pe9pGuW4Gnm3itZ2ATNt/O9oed2yrXI3Wv4+Gocpdur1s7z0BGAHsbOb5qcAnNMxeNxrY6Krt5bZ77saY3caYvS2sNgrIMMZkGmOqgbeA6SIiwMXActt6rwFXOinadNv7Ofq+M4BPjDHlTvr85rQ21w+s3l7GmH3GmP22x7nAUaDFmzTOQpM/L2fIuxyYZNs+04G3jDFVxpiDQIbt/doklzFmjd3P0AYaZkRzNUe2V3MuAz43xhQbY44DnwOTLco1C1jqpM8+I2PM1zQxC52d6cBi02AD0EFEuuGC7eW25e6gpibv7kHD5NwnjDG1jZY7QxdjTB6A7b+dW1h/Jj/+wXrC9ivZP0QkuI1zhYhIqohsOH2oCDfaXiIyioa9sQN2i521vZr7eWlyHdv2OD3ZuyOvdWUue7fRsPd3WlPf07bMdY3t+7NcRE5PyekW28t2+Kon8KXdYldtL0c0l93p28vSqdpF5AugaxNPPWqMWeHIWzSxzJxh+TnncvQ9bO/TDRhMwyxWpz0M5NNQYAuA3wCPt2GueGNMrogkAV+KyA7gZBPrWbW9lgBzjTGnZ6o+6+3V1Ec0sczRyd7P6WeqBQ6/t4jMBlKAiXaLf/Q9NcYcaOr1Lsj1IbDUGFMlIvNo+K3nYgdf68pcp80Elhtj6uyWuWp7OaLNfr4sLXdjzCXn+BbNTd5dRMOvOwG2va+mJvU+q1wiUiAi3YwxebYyOnqGt7oOeN8YU2P33nm2h1Ui8irwq7bMZTvsgTEmU0TWAsOBd7F4e4lIBPAx8Fvbr6un3/ust1cTzmWyd0de68pciMglNPyDOdEYU3V6eTPfU2eUVYu5jDHH7L5cCDxp99oLG712rRMyOZTLzkzgZ/YLXLi9HNFcdqdvL08/LNPk5N2m4QzFGhqOd0PD5N2O/CbgCPvJwFt63x8d67MV3Onj3FcCTZ5Vd0UuEel4+rCGiEQD44BdVm8v2/fufRqORb7T6Dlnbq9zmex9JTBTGq6m6QkkA9+fQ5ZW5RKR4cCLwDRjzFG75U1+T9swVze7L6cBu22PVwOX2vJ1BC7lv3+DdWkuW7a+NJycXG+3zJXbyxErgTm2q2ZGAyW2HRjnby9XnTU+1z/AVTT8a1YFFACrbcu7A6vs1psK7KPhX95H7ZYn0fCXLwN4Bwh2Uq4o4D/Aftt/O9mWpwAv2a2XCBwB/Bq9/ktgBw0l9TrQvq1yAWNtn73N9t/b3GF7AbOBGmCr3Z9hrtheTf280HCYZ5rtcYjt/z/Dtj2S7F77qO11e4EpTv55bynXF7a/B6e3z8qWvqdtlOuvQLrt89cA/exee6ttO2YAt7RlLtvXfwD+1uh1rt5eS2m42quGhv66DZgHzLM9L8Bzttw7sLsS0NnbS+9QVUopL+Tph2WUUko1QctdKaW8kJa7Ukp5IS13pZTyQlruSinlhbTclVLKC2m5K6WUF9JyV0opL/T/AS/4aG9J25KUAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "def obj(x):\n",
    "    return x**2\n",
    "\n",
    "x = np.linspace(-1,1,100)\n",
    "plt.plot(x, obj(x))\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "      fun: 2.5388964258140254e-16\n",
       " hess_inv: array([[0.5]])\n",
       "      jac: array([4.67689909e-08])\n",
       "  message: 'Optimization terminated successfully.'\n",
       "     nfev: 9\n",
       "      nit: 2\n",
       "     njev: 3\n",
       "   status: 0\n",
       "  success: True\n",
       "        x: array([1.59339149e-08])"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from scipy.optimize import minimize\n",
    "\n",
    "minimize(obj, x0=3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Minimize Return value\n",
    "----\n",
    "         fun: The value of the function at the minimum\n",
    "    hess_inv: The inverse of the the Hessian\n",
    "         jac: The value of the Jacobian\n",
    "     message: A string describing what happened\n",
    "        nfev: Number of function evaluations\n",
    "         nit: Number of iterations of the x point\n",
    "        njev: Number of times it computed the Jacobian\n",
    "      status: The single digit message (0 = success, != 0 some error)\n",
    "     success: Boolean indicating success\n",
    "           x: The minimum x"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Objective Functions\n",
    "===\n",
    "\n",
    "Be careful that your objective function is *convex* and it's minimum isn't at $\\infty$ or $-\\infty$. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "A good objective function\n",
    "----\n",
    "\n",
    "Minimize the following:\n",
    "\n",
    "$$f(x) = \\frac{(x - 4)^2}{2} + \\frac{(x - 2)^2}{4}$$ "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXoAAAD8CAYAAAB5Pm/hAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl4VOX9/vH3JwsJgQAJCRAIIewQEQXCIqJ1Q8UNXAtqRaVS61LXn9VaW2trq7W1Sl1RXKu4L2hxobhWWQRk3xL2ECABQiAEsj6/PzL4TTGQkJnJmUzu13Xlysw5ZzI3Z4Y7J8+cxZxziIhI+IrwOoCIiASXil5EJMyp6EVEwpyKXkQkzKnoRUTCnIpeRCTMqehFRMJcrUVvZs+ZWZ6ZLa1h3u1m5swsyXffzGySmWWb2WIzGxiM0CIiUnd12aJ/ATjz4Ilm1hkYCWysNnkU0NP3NRF40v+IIiLij6jaFnDOfWVm6TXM+gdwB/B+tWmjgZdc1eG2s82sjZmlOOe2HO45kpKSXHp6TU8hIiKHMn/+/O3OueTalqu16GtiZucBm51zi8ys+qxOwKZq93N8035U9GY2kaqtftLS0pg3b159ooiINFlmtqEuyx3xh7FmFgfcDfyuptk1TKvxZDrOucnOuUznXGZycq2/kEREpJ7qs0XfHegKHNiaTwUWmNkQqrbgO1dbNhXI9TekiIjU3xFv0Tvnljjn2jnn0p1z6VSV+0Dn3FZgGnCFb++bYUBhbePzIiISXHXZvXIqMAvobWY5ZjbhMItPB9YC2cAzwHUBSSkiIvVWl71uxtUyP73abQdc738sEREJFB0ZKyIS5lT0IiJhrlEX/YYde/nDB8soq6j0OoqISMhq1EWfnVfE89+s5+35OV5HEREJWY266E/p045jO7dh0swsSsorvI4jIhKSGnXRmxn/74ze5BbuZ+qcjbU/QESkCWrURQ8wvHtbhnVL5LHP11BcWu51HBGRkNPoi97MuP303mwvKuGlWXU6v4+ISJPS6IseIDM9kZN6J/PUl2vYvb/M6zgiIiElLIoe4PbTe7OruIxnv17ndRQRkZASNkXfr1Nrzjq6A1O+XsuOohKv44iIhIywKXqAW0f2Yl9ZBU98scbrKCIiISOsir5Hu3guHJjKy7M3kLtrn9dxRERCQlgVPcBNp/XEOcekmVleRxERCQlhV/SpCXFcNrQLb87PYW1+kddxREQ8F3ZFD3D9yT2IiYrg7zNWex1FRMRzYVn0yfExTBjRlX8v3sKSnEKv44iIeCosix5g4ondSIiL5q+frPQ6ioiIp8K26ONjo7n+5B58nbWdb7K3ex1HRMQzYVv0AJcP60LH1rE8+PFKqi5nKyLS9NRa9Gb2nJnlmdnSatMeMrOVZrbYzN41szbV5t1lZtlmtsrMzghW8LqIjY7klpG9WJxTyPQlW72MIiLimbps0b8AnHnQtBlAP+dcf2A1cBeAmWUAY4GjfI95wswiA5a2Hi4YmEqv9i3526erdMlBEWmSai1659xXwM6Dpn3qnDtw8vfZQKrv9mjgNedciXNuHZANDAlg3iMWGWH8+sw+rNu+l9fm6uIkItL0BGKM/mrgI9/tTsCmavNyfNM8dUqfdgzpmsijM7MoKtHFSUSkafGr6M3sbqAceOXApBoWq/FTUDObaGbzzGxefn6+PzFqZWbcNaoP24tKefbrtUF9LhGRUFPvojez8cA5wGXu/3ZpyQE6V1ssFcit6fHOucnOuUznXGZycnJ9Y9TZgLQEzjq6A5O/Wkv+Hp3GWESajnoVvZmdCfwaOM85V1xt1jRgrJnFmFlXoCcw1/+YgXH76b0pLa/k0Zk6NYKINB112b1yKjAL6G1mOWY2AXgMiAdmmNlCM3sKwDm3DHgDWA58DFzvnKsIWvoj1C25JZcOTWPq3E1k5+mEZyLSNFgoHEiUmZnp5s2b1yDPtb2ohJMe+oLjurflmSsyG+Q5RUSCwczmO+dqLbKwPjK2JkktY7j2J92YsXwbc9ftrP0BIiKNXJMreoAJI7rRoVUs909foVMjiEjYa5JF37xZJLee3otFm3bx4eItXscREQmqJln0ABcOTKVPh3j++slKSspD5vNiEZGAa7JFHxlh/PbsDDbt3MeL3673Oo6ISNA02aIHGNEziZN7J/PPz7LZubfU6zgiIkHRpIse4Ddn9aW4tIJH/6ODqEQkPDX5ou/ZPp5xQzrzypyNrMnXQVQiEn6afNED3HxaL5pHR/Lnf6/wOoqISMCp6Kk6iOr6U3owc2Ue/83S9WVFJLyo6H2uOj6dtMQ4/vjhcsp1JSoRCSMqep+YqEh+c1YfVm3bw+vzNtX+ABGRRkJFX80ZR3VgaNdE/v7panbvL/M6johIQKjoqzEz7jkng4LiUv45M8vrOCIiAaGiP0i/Tq25ZFBnXvh2PWu1u6WIhAEVfQ1uP6M3MVGR3K/dLUUkDKjoa5AcH8OvTq3a3fKLVXlexxER8YuK/hCuHN6VrkktuO/D5ZRpd0sRacRU9IfQLCqCe87py9r8vTq7pYg0air6wzi5dzt+0iuZR2dmsb2oxOs4IiL1oqI/DDPjd+dmsL+sgr9+vNLrOCIi9VJr0ZvZc2aWZ2ZLq01LNLMZZpbl+57gm25mNsnMss1ssZkNDGb4htA9uSVXH9+VN+blsHDTLq/jiIgcsbps0b8AnHnQtDuBmc65nsBM332AUUBP39dE4MnAxPTWDaf0IDk+ht+/v5TKSl1MXEQal1qL3jn3FbDzoMmjgRd9t18ExlSb/pKrMhtoY2YpgQrrlfjYaO4a1YdFOYW8tSDH6zgiIkekvmP07Z1zWwB839v5pncCqp8RLMc3rdE7f0AnBnVJ4MGPVlK4T+fBEZHGI9AfxloN02oc6zCziWY2z8zm5efnBzhG4JkZfzjvKAqKS3n401VexxERqbP6Fv22A0Myvu8HDh/NATpXWy4VyK3pBzjnJjvnMp1zmcnJyfWM0bD6dWrN5cO68PLsDSzLLfQ6johIndS36KcB4323xwPvV5t+hW/vm2FA4YEhnnBx28jeJMQ143fvL9MHsyLSKNRl98qpwCygt5nlmNkE4AFgpJllASN99wGmA2uBbOAZ4LqgpPZQ67hofj2qD/M3FPDO95u9jiMiUquo2hZwzo07xKxTa1jWAdf7GyrUXTQwldfmbuQv01cwsm97WsdFex1JROSQdGRsPUREGPeN7kdBcSkPfaojZkUktKno66lfp9ZccVw6r8zZyCIdMSsiIUxF74fbTu9FcssYfvveUir0wayIhCgVvR/iY6P57TkZLNlcyKtzNngdR0SkRip6P53bP4URPZL46yeryNuz3+s4IiI/oqL3k5lx3+ijKCmr5I8f6hqzIhJ6VPQB0C25Jded3J0PFuXy1erQP52DiDQtKvoA+eVJ3emW1ILfvreU/WUVXscREfmBij5AYqIi+dP5/di4s5jHPsv2Oo6IyA9U9AE0vHsSFwzoxNNfrWH1tj1exxERAVT0AXf32X1pGRPFXe8s0UnPRCQkqOgDrG3LGH57dgbzNxTw6tyNXscREVHRB8MFAztxfI+2PPjRSrbt1r71IuItFX0QmBn3jzma0opK7p22zOs4ItLEqeiDJD2pBTed1pOPlm7lk2VbvY4jIk2Yij6IrjmhGxkprbjnvaW6oLiIeEZFH0TRkRE8eGF/theV8MBHOj2CiHhDRR9kR6e25poTujF17ia+XbPd6zgi0gSp6BvAzaf1okvbOO56Zwn7SnV6BBFpWCr6BtC8WSR/ueBoNuwo5uEZq7yOIyJNjIq+gQzvnsSlQ9OY8t91LNhY4HUcEWlC/Cp6M7vFzJaZ2VIzm2pmsWbW1czmmFmWmb1uZs0CFbaxu2tUHzq0iuWOtxbrDJci0mDqXfRm1gn4FZDpnOsHRAJjgQeBfzjnegIFwIRABA0H8bHR/PmCo8nOK2LSzCyv44hIE+Hv0E0U0NzMooA4YAtwCvCWb/6LwBg/nyOsnNS7HRcNSuXpr9ayJKfQ6zgi0gTUu+idc5uBvwEbqSr4QmA+sMs5V+5bLAfoVNPjzWyimc0zs3n5+U3rqkz3nJ1BUstm3PbmQkrKNYQjIsHlz9BNAjAa6Ap0BFoAo2pYtMZz9TrnJjvnMp1zmcnJyfWN0Si1jovmgQv6s3qbhnBEJPj8Gbo5DVjnnMt3zpUB7wDDgTa+oRyAVCDXz4xh6eQ+7bh4UCpPfrGGRZt2eR1HRMKYP0W/ERhmZnFmZsCpwHLgc+Ai3zLjgff9ixi+fntOBu1bxXLbm4u0F46IBI0/Y/RzqPrQdQGwxPezJgO/Bm41s2ygLTAlADnDUuvm0TxwYX+y84r4x4zVXscRkTAVVfsih+ac+z3w+4MmrwWG+PNzm5Kf9Epm3JA0Jn+9ltMy2jM4PdHrSCISZnRkbAj47dl96ZwQx21vLGJvSXntDxAROQIq+hDQIiaKv118DJsKirl/uk5nLCKBpaIPEUO6JnLNCd14dc5GPl+V53UcEQkjKvoQcuvIXvRq35I73lrMzr2lXscRkTChog8hsdGRPPLTAewqLuU37yzBuRqPNRMROSIq+hCT0bEVt53em4+XbeXtBZu9jiMiYUBFH4KuOaEbQ7omcu+0ZWzaWex1HBFp5FT0ISgywnj4kmMw4JbXF1JeUel1JBFpxFT0ISo1IY4/junHvA0FPP75Gq/jiEgjpqIPYWMGdGLMsR2Z9FkW8zfo8oMiUj8q+hB335h+dGgVy82vf8+e/WVexxGRRkhFH+JaxUbz6Nhj2Vywj9+9v8zrOCLSCKnoG4HM9ERuPKUn736/mXcW5HgdR0QaGRV9I3HjKT0Ykp7IPe8tZd32vV7HEZFGREXfSERFRvDI2GOJiozgV1O/p7Rcu1yKSN2o6BuRjm2a8+CF/VmyuZC/frzS6zgi0kio6BuZM/t14GfDuvDsf9cxc8U2r+OISCOgom+E7j67LxkprbjtzUXk7trndRwRCXEq+kYoNjqSxy8bSFl5Jb+a+j1lOkWCiByGir6R6prUgj9fcDTzNhTwsC4sLiKH4VfRm1kbM3vLzFaa2QozO87MEs1shpll+b4nBCqs/K/Rx3Zi3JA0nvxiDZ+t1Hi9iNTM3y36R4GPnXN9gGOAFcCdwEznXE9gpu++BMnvz80gI6UVt7y+iJwCndJYRH6s3kVvZq2AE4EpAM65UufcLmA08KJvsReBMf6GlEOLjY7kycsHUlnpuP5V7V8vIj/mzxZ9NyAfeN7MvjezZ82sBdDeObcFwPe9XQByymF0aduChy7uz6JNu/jz9BVexxGREONP0UcBA4EnnXMDgL0cwTCNmU00s3lmNi8/P9+PGAJwZr8UJozoygvfruf9hboEoYj8H3+KPgfIcc7N8d1/i6ri32ZmKQC+73k1Pdg5N9k5l+mcy0xOTvYjhhxw56g+DE5P4M63l7Bq6x6v44hIiKh30TvntgKbzKy3b9KpwHJgGjDeN2088L5fCaXOoiMjePzSgbSMjeKX/5qv89eLCOD/Xjc3Aq+Y2WLgWODPwAPASDPLAkb67ksDadcqlscvHciGncXc/uYiKiud15FExGN+Fb1zbqFv+KW/c26Mc67AObfDOXeqc66n7/vOQIWVuhnSNZG7RvXhk2XbePJLXW9WpKnTkbFhasKIrow+tiN/+3QVn6+q8WMSEWkiVPRhysx44IL+9OnQipumfs96XaxEpMlS0Yex5s0imfyzQUREGBNfnkdRSbnXkUTEAyr6MNc5MY7Hxg1kTf5ebn19oT6cFWmCVPRNwIieSdx9Vl8+Xb6NR2ZmeR1HRBpYlNcBpGFcdXw6K7bsZtLMLPp2iGfU0SleRxKRBqIt+ibCzPjT+f0YmNaGW99YxNLNhV5HEpEGoqJvQmKiInnqZ4NIiIvmmpfmkbd7v9eRRKQBqOibmHbxsTwzPpNdxWVc8/J89pdVeB1JRIJMRd8EHdWxNY+MPZbFObv4f28txjntiSMSzlT0TdQZR3XgjjP68MGiXP7xH+2JIxLOtNdNE3btT7qxbnsRk2Zm0SUxjgsHpXodSUSCQEXfhJkZfxpzNJt27uPOdxbTKaE5w7q19TqWiASYhm6auGZRETx1+SDSEuP4xcvzWZNf5HUkEQkwFb3QOi6a568cQnSkceXzc8nfU+J1JBEJIBW9AJDWNo4p4wezfU8pV7/wHXt1AjSRsKGilx8c07kNj106gGW5hdzw6gLKKyq9jiQiAaCil/9xat/2/HFMPz5flc9v3l2ifexFwoD2upEfuWxoF7YV7mfSZ9m0i4/l9jN61/4gEQlZKnqp0S0je5FfVMJjn2eTHB/D+OHpXkcSkXpS0UuNzIw/ju7H9qJS7v1gGYktmnHuMR29jiUi9eD3GL2ZRZrZ92b2oe9+VzObY2ZZZva6mTXzP6Z4ISoygn+OG8DgLonc+sZCvtBFxkUapUB8GHsTsKLa/QeBfzjnegIFwIQAPId4JDY6kmevzKRnu3iu/dd85m/Y6XUkETlCfhW9maUCZwPP+u4bcArwlm+RF4Ex/jyHeK9VbDQvTRhCSuvmXPX8dyzP3e11JBE5Av5u0T8C3AEc2OG6LbDLOXfgaJscoJOfzyEhIKllDC9PGEKLmCh+NmUO2Xk6VYJIY1Hvojezc4A859z86pNrWLTGHbHNbKKZzTOzefn5+fWNIQ0oNSGOf/18KGZw+bNz2LSz2OtIIlIH/mzRHw+cZ2brgdeoGrJ5BGhjZgf25kkFcmt6sHNusnMu0zmXmZyc7EcMaUjdk1vy8oSh7Cur4NJnZ7OlcJ/XkUSkFvUueufcXc65VOdcOjAW+Mw5dxnwOXCRb7HxwPt+p5SQ0jelFS9dPYSCvWVc+swcXXtWJMQF4xQIvwZuNbNsqsbspwThOcRjx3RuwwtXDWbb7v2Me2a2zngpEsICUvTOuS+cc+f4bq91zg1xzvVwzl3snFMDhKnM9ESev3Iwubv2c9mzs9lRpJdaJBTppGbil6Hd2jLlykw27izm0mfmsF1lLxJyVPTit+Hdk5gyfjAbdu5l3GQN44iEGhW9BMTxPZJ47srB5BTsY9wzs8nbow9oRWrzbfZ2CveVBf15VPQSMMO7J/H8VYPJ3bWPsU9r10uRw/n34i1c8dxcHvx4ZdCfS0UvATWsW1teunoI+XtKuPipWWzcoYOqRA729vwcbpy6gGM7t+HOUX2C/nwqegm4zPREXrlmKEUl5Vzy9CydLkGkmlfmbOC2NxdxXPe2vDRhCK1io4P+nCp6CYr+qW14beIwyisr+enTs1i6udDrSCKee+rLNdz97lJO6dOOKeMHE9esYS4JoqKXoOnToRVv/OI4YqMjGTd5NnPW7vA6kognnHM88NFKHvhoJef0T+GpywcRGx3ZYM+vopeg6pbckrd+eRztWsVwxXNzmblim9eRRBpURaXjN+8u5akv13DZ0DQeHTuAZlENW70qegm6lNbNefPa4fTuEM/El+fzxrxNXkcSaRD7yyq4/pUFTJ27ketO6s6fxvQjMqKmk/wGl4peGkRii2a8es0whndvyx1vLebxz7NxrsYzWIuEhcJ9ZVwxZS4fL9vKPedkcMeZfai6NlPDU9FLg2kZE8WU8YMZfWxHHvpkFb+ftoyKSpW9hJ8thfu45KlZfL+pgEnjBjBhRFdP8zTMR74iPs2iIvjHJcfSLj6GZ75ex5bC/UwaO4DmzRrugymRYFqWW8jVL3zH3pIKXrhqCMf3SPI6krbopeFFRBh3n53BvedmMHPFNsZOnqXz40hY+HJ1Ppc8NYsIM9765XEhUfKgohcPXXl8V57+WSartu3h/Ce+YfW2PV5HEqm3V+Zs4OoXviOtbQveve54+nRo5XWkH6joxVMjM9rzxi+Oo6S8kguf+JYvVuV5HUnkiFRUOv7wwTLufncpJ/ZM4s1rj6ND61ivY/0PFb14rn9qG96//ng6J8Zx9Qvf8fw367RHjjQKe/aX8fMXv+P5b9Zz9fFdeXb8YFrGhN5Hnyp6CQkd2zTnzWuP49S+7fnDB8u58+0llJRXeB1L5JDWbd/L+U98y1dZ27n//H787twMT/aRrwsVvYSMFjFRPH35IG44uQevz9vEuMk6r72Epi9W5TH6sf+yo6iEf00YymVDu3gd6bBU9BJSIiKM28/ozeOXDmTFlj2c989vWLCxwOtYIkDVOWue/GINV7/wHR3bNGfaDSM4rntbr2PVSkUvIens/im8/cvhNIuK4KdPz+Ll2Rs0bi+e2rO/jGv/NZ8HP17JqKNTeOe64XROjPM6Vp3Uu+jNrLOZfW5mK8xsmZnd5JueaGYzzCzL9z0hcHGlKcno2IoPbhjBiB5J3PPeUm57YxHFpeVex5ImaNXWPYx+7Bv+syKPe87J4LFxAxrsFMOB4M8WfTlwm3OuLzAMuN7MMoA7gZnOuZ7ATN99kXppHRfNlPGDufm0nry7cDNjHv+GLO1vLw3ojXmbGP34f9m9v5xXfz6UCSO6enbOmvqqd9E757Y45xb4bu8BVgCdgNHAi77FXgTG+BtSmraICOPm03rx0tVD2FFUynmPfcPb83O8jiVhrri0nNveWMQdby1mQOcEpt80gqHdQn88viYWiHFPM0sHvgL6ARudc22qzStwzv1o+MbMJgITAdLS0gZt2LDB7xwS/rbt3s+NU79n7rqdXDCwE/eN7heS+y1L47Y8dzc3Tl3A2u17ufGUntx0as+Q3HXSzOY75zJrXc7fojezlsCXwP3OuXfMbFddir66zMxMN2/ePL9ySNNRXlHJpM+yeeyzLNIS45g0bgD9U9vU/kCRWjjneP6b9Tzw0UraxEXzj58eGzLnq6lJXYver71uzCwaeBt4xTn3jm/yNjNL8c1PAXRMuwRUVGQEt47sxWsTj6O0vJILnviWxz/P1imPxS95u/dz1Qvfcd+HyzmhZxIf33xiSJf8kfBnrxsDpgArnHMPV5s1DRjvuz0eeL/+8UQObUjXRKbfdAJnHNWBhz5ZxSVPz2LDjr1ex5JGaPqSLZz+yFfMWrODP5x3FM+OzySxRTOvYwVMvYduzGwE8DWwBKj0Tf4NMAd4A0gDNgIXO+d2Hu5naehG/OGc4/2Fudzz/tKq63Oe1ZfLhqY1uj0jpOHtKi7l3mnLeG9hLv1TW/PwJcfSo11Lr2PVWYON0QeCil4CIXfXPu54azH/zd7O8O5tefDC/o3mgBZpeJ8u28rd7y2lYG8p153cgxtP6UF0ZOM6hlRFL02Sc46pczdx/7+X44Bfn9mHy4d1Cck9JsQb24tKuO+D5UxblEvflFY8dFF/+nVq7XWselHRS5OWU1DMXe8s4eus7QxIa8NfLjg6pC4EIQ3POcdb83O4f/oK9paUc/3JPbjupB40i2pcW/HVqeilyTswdn/fh8vZva+Mn5/QjV+d2qNRHbougZGdt4d73lvGrLU7yOySwAMXHk2PdvFex/JbXYte73gJW2bGmAGdOLFXMn+evoKnvlzDtIWb+d25GZxxVAd9WNsE7C0pZ9JnWUz5eh1xzSL505h+XDokjYgmNpSnLXppMr5bv5N73lvKyq17OKFnEveck0Gv9o1/q05+zDnHtEW5PPDRSrYU7ueiQancOaoPSS1jvI4WUBq6EalBeUUlL83awCP/Wc3e0gouH5rGzaf1IiGM9plu6hbn7OIPHyxn/oYC+nVqxb3nHkVmeqLXsYJCRS9yGDv3lvLwjFW8OmcjLWOiuP7kHowfnk5sdKTX0aSeNu0s5m+fruL9hbkktWzGHWf04cJBqWG9x5WKXqQOVm3dwwMfreDzVfl0atOcW0b24vwBncK6HMLNjqISnvhiDS/P2kBEBPx8RDd+8ZNuxMdGex0t6FT0Ikfg2+zt/OWjlSzZXEj35BbcMrIXZ/VLaXIf2jUmhcVlPPP1Wp77Zh37yyq4aFAqt47sTYfWsV5HazAqepEj5Jzjk2Vb+funq8nKK6JPh3huOKUHo/qlaAs/hOwqLuW5b9bz/Dfr2LO/nLP7p3DLab0a1akLAkVFL1JPFZWODxblMumzLNbm76V7cguuP7kH5x7TsdEdIh9O8vbs57n/ruflWevZW1rB6Rntufm0XmR0bLoHwqnoRfxUUen4aOkWHvssm5Vb95DSOparj+/K2CGdm8T4b6jIzivi2a/X8s6CzZRVVnJO/45cf3J3HemMil4kYJxzfLEqn6e/WsPstTuJj4niosxUxh+XTnpSC6/jhaXKSsdXWfm88O16vliVT0xUBBdnpvLzEd20zqtR0YsEwaJNu3jum3X8e/EWKpzjpF7JXDq0Cyf3TiZKwzp+K9hbytsLcnh1zkbWbt9LcnwMlw1N42fDutA2zA52CgQVvUgQ5e3ezytzNvLq3I3k7ymhfasYfprZmQsHpdKlrbY4j0RFpWP22h28MW8THy3ZSmlFJQPS2nDl8HRG9Utp1CcdCzYVvUgDKKuo5LOVeUydu5EvV+fjHAxOT+CCgamM6teBNnE64vZQVm/bw3vfb+bd7zezpXA/8bFRXDCgE2OHpNE3RePvdaGiF2lgubv28d7Czbw9P4c1+XuJijBG9Ezi7KNTOK1v+yZ/mgXnHNl5RUxfspV/L8ll9bYiIiOME3smccHAVEZmtNeRyUdIRS/iEeccy3J388HiXD5ctIXNu/YRGWEMTk9gZEYHTunTjvS2cU3i7JllFZUs2FDAzJV5zFi+jXXb92IGg7skcs4xKYzql0JyvMbe60tFLxICnHMsyilkxvKtzFi+jdXbigBIS4zjxF5JHN89iaHd2obNhaidc6zJL2LW2p18vTqfb9fsoKiknOhIY1i3tpx+VAdG9m3fpI5eDSYVvUgI2rijmC9X5/Hl6u18u2Y7xaUVAPRuH09megID0xIYkNaGrkktGsUW//6yCpbl7ub7jQUs2FjA3HU72V5UCkCnNs05sVcyP+mVzPE92urYgyDwvOjN7EzgUSASeNY598ChllXRS1NUVlHJ4pxCZq/dwey1O1i4cRd7SsoBiI+NIiOlFUd1bE2flHh6tGtJj3YtaeVRWTrn2FK4n+y8IrLyilieu5tluYVk5xVRXlnVIakullqOAAAIM0lEQVQJzcnsksBx3dsyrFtb0hKbxvCUlzwtejOLBFYDI4Ec4DtgnHNueU3Lq+hFqnYzXJNfxIINBSzNLWTp5t2s3Lqb/WWVPyyT1DKGtMTmpCXG0SmhOR1axdKuVSzJ8TEkxDUjIS6a+NjoIzo3T0l5BYXFZezaV8bOvaVs272fvN0lbN29n007i9no+zrw18eBHEd1bMVRHVvRP7UNA9Pa0K6VhmMamteXEhwCZDvn1vrCvAaMBmosehGByAijV/v4/7nqVUWlY9POYrLyisjK28OG7VWl+936Aj5YvIWKypo31GKjI4hrFkXz6EgiI4yoCMMMKh2UV1ZSVu4oLi1nX1kFZRU1/4zm0ZF09v1SGdatLd3btaSn7y+LcLtSU7gLVtF3AjZVu58DDA3Sc4mErcgIIz2pBelJLRiZ0f5/5lVUOnbsLSFvdwn5e0ooKC6loLiM3fvK2FdWwd6ScvaXVVLpHOWVjspK90PpR0YYcc0iiYuJokWzSFr7/hpIiGtG+1YxtGsVS3xMlIZewkSwir6md8f/bDaY2URgIkBaWlqQYoiEr8gIo118LO3iNWQihxesY4tzgM7V7qcCudUXcM5Nds5lOucyk5OTgxRDRESCVfTfAT3NrKuZNQPGAtOC9FwiInIYQRm6cc6Vm9kNwCdU7V75nHNuWTCeS0REDi9YY/Q456YD04P180VEpG50/k8RkTCnohcRCXMqehGRMKeiFxEJcyFx9kozywc21PPhScD2AMYJlFDNBaGbTbmOjHIdmXDM1cU5V+uBSCFR9P4ws3l1OalPQwvVXBC62ZTryCjXkWnKuTR0IyIS5lT0IiJhLhyKfrLXAQ4hVHNB6GZTriOjXEemyeZq9GP0IiJyeOGwRS8iIofRKIrezC42s2VmVmlmmQfNu8vMss1slZmdcYjHdzWzOWaWZWav+86oGeiMr5vZQt/XejNbeIjl1pvZEt9yQb9+opnda2abq2U76xDLnelbh9lmdmcD5HrIzFaa2WIze9fM2hxiuQZZX7X9+80sxvcaZ/veS+nBylLtOTub2edmtsL3/r+phmVOMrPCaq/v74Kdq9pzH/a1sSqTfOtssZkNbIBMvauti4VmttvMbj5omQZZZ2b2nJnlmdnSatMSzWyGr4tmmFnCIR473rdMlpmN9zuMcy7kv4C+QG/gCyCz2vQMYBEQA3QF1gCRNTz+DWCs7/ZTwC+DnPfvwO8OMW89kNSA6+5e4PZalon0rbtuQDPfOs0Icq7TgSjf7QeBB71aX3X59wPXAU/5bo8FXm+A1y4FGOi7HU/VdZgPznUS8GFDvZ+O5LUBzgI+oupCRMOAOQ2cLxLYStW+5g2+zoATgYHA0mrT/grc6bt9Z03veyARWOv7nuC7neBPlkaxRe+cW+GcW1XDrNHAa865EufcOiCbquvV/sCqroV2CvCWb9KLwJhgZfU93yXA1GA9RxD8cI1f51wpcOAav0HjnPvUOVfuuzubqovTeKUu//7RVL13oOq9dKoF+Tp7zrktzrkFvtt7gBVUXaazsRgNvOSqzAbamFlKAz7/qcAa51x9D8b0i3PuK2DnQZOrv48O1UVnADOcczudcwXADOBMf7I0iqI/jJquTXvwf4S2wK5qpVLTMoF0ArDNOZd1iPkO+NTM5vsup9gQbvD96fzcIf5UrMt6DKarqdryq0lDrK+6/Pt/WMb3Xiqk6r3VIHxDRQOAOTXMPs7MFpnZR2Z2VENlovbXxuv31VgOvcHl1Tpr75zbAlW/yIF2NSwT8PUWtPPRHykz+w/QoYZZdzvn3j/Uw2qYdvBuRHVZpk7qmHEch9+aP945l2tm7YAZZrbS95u/3g6XC3gS+CNV/+Y/UjWsdPXBP6KGx/q9O1Zd1peZ3Q2UA68c4scEfH3VFLWGaUF7Hx0pM2sJvA3c7JzbfdDsBVQNTRT5Pn95D+jZELmo/bXxcp01A84D7qphtpfrrC4Cvt5Cpuidc6fV42G1XpuWqnNItDGzKN+WWE3LBCSjmUUBFwCDDvMzcn3f88zsXaqGDfwqrrquOzN7Bviwhll1WY8Bz+X7kOkc4FTnG5ys4WcEfH3VoC7//gPL5Phe59b8+M/ygDOzaKpK/hXn3DsHz69e/M656Wb2hJklOeeCfk6XOrw2QXlf1dEoYIFzbtvBM7xcZ8A2M0txzm3xDWPl1bBMDlWfIxyQStXnk/XW2IdupgFjfXtEdKXqt/Lc6gv4CuRz4CLfpPHAof5C8NdpwErnXE5NM82shZnFH7hN1QeSS2taNlAOGhM9/xDP1+DX+DWzM4FfA+c554oPsUxDra+6/PunUfXegar30meH+uUUKL7PAKYAK5xzDx9imQ4HPiswsyFU/Z/eEcxcvueqy2szDbjCt/fNMKDwwLBFAzjkX9ZerTOf6u+jQ3XRJ8DpZpbgG2o93Tet/oL9yXMgvqgqqBygBNgGfFJt3t1U7TGxChhVbfp0oKPvdjeqfgFkA28CMUHK+QJw7UHTOgLTq+VY5PtaRtUQRrDX3cvAEmCx702WcnAu3/2zqNqrY00D5cqmahxyoe/rqYNzNeT6qunfD9xH1S8igFjfeyfb917q1gDraARVf7IvrraezgKuPfA+A27wrZtFVH2oPTzYuQ732hyUzYDHfet0CdX2mAtytjiqirt1tWkNvs6o+kWzBSjz9dcEqj7XmQlk+b4n+pbNBJ6t9tirfe+1bOAqf7PoyFgRkTDX2IduRESkFip6EZEwp6IXEQlzKnoRkTCnohcRCXMqehGRMKeiFxEJcyp6EZEw9/8BHqTqzXLSzIcAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "def f(x):\n",
    "    return ((x-4)**2)/2+((x-2)**2)/4\n",
    "x = np.linspace(-10, 10, 100)\n",
    "plt.plot(x, f(x))\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "the minimum value occurs at [3.33333336] and is 0.6666666666666674\n"
     ]
    }
   ],
   "source": [
    "result = minimize(f, x0=0)\n",
    "print(f'the minimum value occurs at {result.x} and is {result.fun}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "A bad objective function\n",
    "----\n",
    "\n",
    "Minimize the following:\n",
    "\n",
    "$$\n",
    "f(x) = \\frac{(x - 4)} { 2 }\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "The minimum is at $-\\infty$!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "the minimum value occurs at [-1.7895697e+08] and is -89478487.25\n"
     ]
    }
   ],
   "source": [
    "result = minimize(lambda x: (x - 4) / 2, x0=0)\n",
    "print(f'the minimum value occurs at {result.x} and is {result.fun}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "A good objective function but bad `x0`\n",
    "---\n",
    "\n",
    "Minimize the following function\n",
    "\n",
    "$$4 \\left[ r^{-12} - r^{-6}\\right]$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "the minimum is at [0.]\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/whitead/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:1: RuntimeWarning: divide by zero encountered in power\n",
      "  \"\"\"Entry point for launching an IPython kernel.\n",
      "/home/whitead/miniconda3/lib/python3.7/site-packages/ipykernel_launcher.py:1: RuntimeWarning: invalid value encountered in subtract\n",
      "  \"\"\"Entry point for launching an IPython kernel.\n"
     ]
    }
   ],
   "source": [
    "result = minimize(lambda r: 4 * (r**(-12) - r**(-6)), x0=0)\n",
    "print(f'the minimum is at {result.x}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "*Our initial value was not in the domain!*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXYAAAD8CAYAAABjAo9vAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAHPFJREFUeJzt3XmQHOd53/HvM/ceswewF7AHLoEAwQMEuQQlk44k6iJpWYpsi6bOiqwKUrGsoqpsWVbsih07qSSVROWoFMmGJZlyRUdsiZRsnWYiStRJCQRAAiQgkAIJYHHtLrD3Mbsz8+SPmQUX4AJYYHp3MD2/T9XUzM403n6aAH588czb3ebuiIhIeETKXYCIiARLwS4iEjIKdhGRkFGwi4iEjIJdRCRkFOwiIiFTcrCb2SYz2zvvMWpmHwqiOBERuXIW5Dp2M4sCx4E73P1IYAOLiMiiBd2KeR3wS4W6iEj5xAIe7wHgi5fbqKWlxdeuXRvwrkVEwu3JJ58cdPfWy20XWCvGzBLACeAGdz+9wOc7gB0APT09tx05okm9iMiVMLMn3b33ctsF2Yq5F9i9UKgDuPtOd+91997W1sv+D0dERK5SkMH+DhbRhhERkaUVSLCbWS3wBuDhIMYTEZGrF8iXp+4+CawMYiwRESmNzjwVEQkZBbuISMgo2EVEQqaigv27B0/zye89X+4yRESuaRUV7D94bpBPPfbLcpchInJNq6hgX1mXYCyTJZPNlbsUEZFrVmUFe30SgLMTM2WuRETk2lVZwV6XAODMuIJdRORiKivY64vBrhm7iMhFVVaw1xVaMWfGM2WuRETk2lVRwb6iXq0YEZHLqahgTydjJKIRtWJERC6hooLdzFhZn1ArRkTkEioq2AFW1CU0YxcRuYSKC/aV9UkFu4jIJVRcsLfUqRUjInIpFRfsK+oSOvNUROQSgro1XpOZfdnMDprZATN7VRDjLmRlfZLJmRyTM9ml2oWISEULasb+P4Fvu/tmYCtwIKBxX2al1rKLiFxSycFuZg3AvwA+A+DuM+4+XOq4FzN3vRi1Y0REFhbEjH09MAD8rZntMbNPm1ldAOMuaO4Kj2cm9AWqiMhCggj2GHAr8Cl33wZMAH904UZmtsPMdpnZroGBgave2dyMfVCtGBGRBQUR7H1An7s/Ufz5yxSC/jzuvtPde929t7W19ap3NtdjVytGRGRhJQe7u58CjpnZpuJbrwOeLXXci6lNxKiJR7WWXUTkImIBjfNB4PNmlgAOA+8LaNwFrahLaFWMiMhFBBLs7r4X6A1irMVoqdf1YkRELqbizjyFuevFqBUjIrKQigz2FXUJzqoVIyKyoIoM9pX1CQYnZnD3cpciInLNqcxgr0swk80zntH1YkRELlShwV44+1Rr2UVEXq4yg71eZ5+KiFxMZQZ7ccauk5RERF6uMoNdlxUQEbmoigz2FcULgekkJRGRl6vIYE/Fo9QnYwyqFSMi8jIVGexQaMeoFSMi8nKVG+y6EJiIyIIqNthX1CXVYxcRWUDFBntLfULLHUVEFlCxwT7XY8/ndb0YEZH5KjbY29Ipsnnn7KTaMSIi81VssLc3FM4+PT06XeZKRESuLYHcQcnMXgTGgByQdfclv5tSW0MKgP7RDDesXuq9iYhUjqDueQrwWncfDHC8S2ovBrtm7CIi56vYVkxr/VwrRitjRETmCyrYHfhnM3vSzHYENOYlJWIRVtYlOD2mGbuIyHxBtWLudPcTZtYGPGpmB9398fkbFAN/B0BPT08gO21rSNGvVoyIyHkCmbG7+4nicz/wCLB9gW12unuvu/e2trYGsVs6GpJqxYiIXKDkYDezOjNLz70G3gjsL3XcxWhvSHFKM3YRkfME0YppBx4xs7nxvuDu3w5g3Mtqa0gxOJ4hm8sTi1bs98AiIoEqOdjd/TCwNYBarlh7QxL3wr1POxpT5ShBROSaU9HT3Pa01rKLiFyosoNdJymJiLxMhQd78SSlMa2MERGZU9HBvrI+ScTQWnYRkXkqOtijEaM1nVQrRkRknooOdphby65WjIjInFAEu1oxIiIvCUGwqxUjIjJf5Qd7OsXQ5CyZbK7cpYiIXBMqP9jn3UlJRERCEOxtxbXs/bouu4gIEIJgf+nsU83YRUQgRMF+akQzdhERCEGwN9fGiUdNt8gTESmq+GA3M9rSKX15KiJSVPHBDlrLLiIyXyiCvaMxpWAXESkKLNjNLGpme8zs60GNuVhqxYiIvCTIGfuDwIEAx1u09oYUY5ksE5lsOXYvInJNCSTYzawL+DXg00GMd6VWNxWWPJ4YnirH7kVErilBzdj/EvhDIH+xDcxsh5ntMrNdAwMDAe22oKu5BoC+IQW7iEjJwW5mbwb63f3JS23n7jvdvdfde1tbW0vd7Xk6m2oB6NOMXUQkkBn7ncBbzOxF4EvA3Wb2vwMYd9Ha0kniUaNvaHI5dysick0qOdjd/aPu3uXua4EHgO+6+7tLruwKRCJGZ1MNx9WKEREJxzp2gM7mGvXYRUQIONjd/Xvu/uYgx1ysrqZajqvHLiISrhn7wFiG6VndSUlEqltogn1uyaPWsotItQtNsHc2aS27iAiEKNi7VhTWsqvPLiLVLjTB3p5OEo1oLbuISGiCPRaNsKoxpbXsIlL1QhPsUPgCVT12Eal2oQr2Tq1lFxEJV7B3NddwanSamexFLzIpIhJ6oQr2zuYa3OHUiG6TJyLVK1TB/tJ12bUyRkSqV7iCXddlFxEJV7B3NKaImM4+FZHqFqpgT8QitDdoLbuIVLdQBTvMrWVXj11EqlcIg71WrRgRqWqhC/bOpsJa9mxOa9lFpDqVHOxmljKzn5nZU2b2jJn9hyAKu1pdzTXk8s5JrWUXkSoVxIw9A9zt7luBW4B7zOyVAYx7VdasrAPgyBn12UWkOpUc7F4wXvwxXnx4qeNerfWthWB/YXD8MluKiIRTID12M4ua2V6gH3jU3Z9YYJsdZrbLzHYNDAwEsdsFtaWT1CaiHB6cWLJ9iIhcywIJdnfPufstQBew3cxuXGCbne7e6+69ra2tQex2QWbGupY6XlCwi0iVCnRVjLsPA98D7gly3Cu1vrWewwMKdhGpTkGsimk1s6bi6xrg9cDBUsctxbqWOvqGJslkc+UsQ0SkLIKYsa8CHjOzp4GfU+ixfz2Aca/a+pY68g7HzmpljIhUn1ipA7j708C2AGoJzLqWwsqYwwMTvKItXeZqRESWV+jOPAVY2zK35FF9dhGpPqEM9saaOC31CQW7iFSlUAY7FNoxWhkjItUotMG+vqVeJymJSFUKbbCva61jcDzD6PRsuUsREVlW4Q324heoL2rWLiJVJrTBvl4rY0SkSoU22HtW1mKGvkAVkaoT2mBPxqJ0NdfoC1QRqTqhDXaAdS31ui67iFSdUAf7+pY6XhiYwL1s9/0QEVl24Q721jomZnIMjGXKXYqIyLIJd7C31APwXL/aMSJSPUId7Js6Cld2PHhqrMyViIgsn1AHe2s6SUt9koMnR8tdiojIsgl1sANcvyqtGbuIVJXQB/vmjjSHTo+RzeXLXYqIyLII4p6n3Wb2mJkdMLNnzOzBIAoLyuaOBjLZPC+e0YlKIlIdgpixZ4Hfd/frgVcCHzCzLQGMG4jNqwpfoB44qXaMiFSHkoPd3U+6++7i6zHgANBZ6rhBeUVbPbGIcfCUvkAVkeoQaI/dzNZSuLH1Ewt8tsPMdpnZroGBgSB3e0nJWJQNrfUc1IxdRKpEYMFuZvXAV4APufvLpsfuvtPde929t7W1NajdLsrmVWkOaMmjiFSJQILdzOIUQv3z7v5wEGMGaXNHAydGphmZ1N2URCT8glgVY8BngAPu/rHSSwre3Beo6rOLSDUIYsZ+J/Ae4G4z21t83BfAuIHZsqoB0KUFRKQ6xEodwN1/CFgAtSyZtnSS5tq4ZuwiUhVCf+YpgJmxuaOBZ7UyRkSqQFUEOxT67IdOjZHL66YbIhJuVRPs13c0MDWb4+jZyXKXIiKypKom2LesLnyBuv/4SJkrERFZWlUT7Ne1p0nGIuw9NlzuUkREllTVBHsiFuHmrkb2HB0qdykiIkuqaoIdYFtPM/tPjJLJ5spdiojIkqmuYO9uYiab1yV8RSTUqivYe5oB1I4RkVCrqmDvaEyxqjHFnqP6AlVEwquqgh1gW08Te45pxi4i4VV9wd7dzLGzUwyMZcpdiojIkqi+YO9pAtRnF5Hwqrpgv7GzkVjE2KMTlUQkpKou2FPxKDesbtCMXURCq+qCHQrLHp/uGyGby5e7FBGRwFVpsDcxOZPj0OnxcpciIhK4oG5m/Vkz6zez/UGMt9RuLZ6otOvI2TJXIiISvKBm7A8B9wQ01pLraq6hq7mGHz43WO5SREQCF0iwu/vjQMVMf82Mu17Rwk8On1GfXURCZ9l67Ga2w8x2mdmugYGB5drtRd21sYWx6Sz7dOMNEQmZZQt2d9/p7r3u3tva2rpcu72oX9nQAqB2jIiETlWuigFYUZfghtUN/PB5BbuIhEvVBjvAXa9oYffRISZnsuUuRUQkMEEtd/wi8BNgk5n1mdn7gxh3qd21sYXZnPPECxXzva+IyGXFghjE3d8RxDjL7fa1K0jEIvzouUFeu6mt3OWIiASiqlsxqXiU3jXN6rOLSKhUdbBDoR1z8NQY/WPT5S5FRCQQCvZXFJY9/kizdhEJiaoP9htWN9JSn+TRZ0+XuxQRkUBUfbBHI8abbmjnsYMDTM3kyl2OiEjJqj7YAe67aRVTszm+f6j8lzoQESmVgh24Y90KmmvjfGv/yXKXIiJSMgU7EItGeOOWDv7fgX4yWbVjRKSyKdiL7rmpg/FMVhcFE5GKp2AvunNDC+lUjG/tP1XuUkRESqJgL0rEIrzh+nYeffY0s7r5hohUMAX7PPfetIqRqVl+8ssz5S5FROSqKdjn+dWNLaSTMb6653i5SxERuWoK9nlS8Shv3baab+w7ycjkbLnLERG5Kgr2C7xjew+ZbJ6H9/SVuxQRkauiYL/ADasb2drdxBeeOIq7l7scEZErFtQdlO4xs1+Y2fNm9kdBjFlO79zezXP94zx5ZKjcpYiIXLGSg93MosD/Au4FtgDvMLMtpY5bTr++dTX1yRhfeOJouUsREbliQczYtwPPu/thd58BvgS8NYBxy6Y2EeNfblvN1/edZHhyptzliIhckSCCvRM4Nu/nvuJ7Fe2d29cwk83zld1a+igilSWIYLcF3nvZt45mtsPMdpnZroGBa//yuFtWN3D72mY+/YPDujCYiCxaPu/MZPNMzmQZnZ7l7MQM/WPTnBie4tjZyWW570MsgDH6gO55P3cBJy7cyN13AjsBent7K2K5yQfv3sh7P/szvvxkH++6Y025yxGperliaGayOTLZPJnZPDO54utsvvhZ4Xkmm2c2V3wvl2c2m2cmN+/9XJ7ZrDObm/dzzpmd93k2V/w8X3g/my9uU/xs7udscZtsLk/+Mun2ud/Zzquva13S/05BBPvPgY1mtg44DjwAvDOAccvuVze2cEt3E5987Je8/bZuEjGtDhW5kLuTyeaZyGSZnMkxNZtjaiZXfJ1laiZfeG82x/Tc57M5ps898kzN5JjOvvRzZi68ZwvP08Xn2Vwwc0IzSEQjxKMRErEI8agVXhffi0WNWDRCMhohGY9QF3lpm1g0Qjxi57YpvC78mnik+Gvn3otc8DoaYVN7OpBjuJSSg93ds2b2e8B3gCjwWXd/puTKrgFmxoOv38j7/vbnPLy7jwe295S7JJFAZHN5xqazjGeyjE1nGZueZTyTfekxnWUik2U8kys8z2SZzGSZyOSYmCkE+FyQT85kLztLvVA8aqRiUVKJKKl4pPA6HqUmHiWditESK7yfnPecLG6XjBcCOFl8PxGLkIwVAnouiBPRwvtzj3jUSBbfixcDOcyCmLHj7t8EvhnEWNea11zXys1djXzisef5zdu6iIf8D4RUBndnejbP0OQMw5OzjEzNMjI1U3wuPEansoXn6VlGp2YZm84WX2eZmr18n9cM6hIx6pJR6pIx6pMxahNRVjWmqE0UXtcWP69JRKmNF36uSUSpTRTeq5l7Lx4llYgUnuNR/T1aYoEEe5iZGQ++biPv/9wuHtl9nPtv7778LxK5QrO5PEMTMwyOz3BmIsPZiRnOjM8wNDnD2YmXnocnZxmanGFocpaZ7MUvLx2NGA2pGI01cRpq4qRTMToaU6SThdfp1NzzS6/rkzHqUzHSyRh1yUIYRyILrY2Qa52CfRHu3tzGzV2N/Pd//gX33tRBOhUvd0lSAfJ550xxRUT/WIaB+Y/xDINjGQbHMwyOF2baCzGD5toEzbVxVtQl6F5Ry81djTTXJmiqTdBUG6epJk5jbZymmgSNtXEaa+LUJaKYKZSrlYJ9EcyMP3/rjbztkz/iY48e4k9//YZylyRlNj2b4/ToNCdHpjk5MsXJkWlOj0xzanSaU6MZ+kenGRjLkF2g+ZxOxmhJJ2mtT7KpI82v1CVpqU+ysj7ByroEK+uTrKhLsKIuQWNNnKhmzXKFFOyLdEt3E++6o4fP/fhFfvPWLm7sbCx3SbKERqdnOXZ2kr6hKY4PTXF8uPB8YmSKE8NTDI6//IzkdCpGR0OKjsYUG9taaG9I0pZO0ZZO0lZ83ZpOkopHy3BEUk0U7Ffgw2/azLf3n+KPH9nHw797p2ZSFSyby3NieJojZyc4cmaSY2cnOVp8HDs7yeh09rzta+JROptrWN1Uww2rG1jdWENHY4rVTTWsaiyEeW1Cf53k2qA/iVegsSbOn/zaFj70f/byhSeO8J5XrS13SXIJ7s6p0WkOD0xweHCCFwYmePHMBC8OTnBsaPK8NdGJWITu5hp6VtRy25pmuptr6Wquoau5ls7mGppr4+pZS8VQsF+ht96ymq/s7uM/ffMAt69bweaOhnKXVPVmc3mOnJngudPjPN8/zvMDhecXBieYnHf6dk08ypqVtWxeleZNN3awdmUta1bWsWZlLe3plFaASGhYOW4m0dvb67t27Vr2/Qalf2yaN3/8h9QlY/zj792pVTLLJJd3jp6d5Benxjh0eoxfnB7j0KkxXhicOO9Lys6mGja01bOhtY71rfVsaCk8tzckNeuWimZmT7p77+W204z9KrSlU3zinbfyjr/5KR/+h6f51LtvVWAEbHR6loMnxzhwcvTc49Dp8fNOrOlZUct17fW8fks7G9vq2diWZkNbnXrdUvX0N+AqbV+3go/eu5n/+I0D7Hz8MP/m1RvKXVLFGhzPsP/4CM+cGD33fPTs5LnPm2vjbO5o4IHt3WzuSLOpo4GNbfXUJfXHV2Qh+ptRgvfftY49x4b5z986SGNNXNeSWYShiRmePj7Cvr5hnu4bYd/xEU6OTJ/7fM3KWm7sbOD+3i62rG5gy6pGtVBErpCCvQRmxsfu38pkJstHH9lHJGLc36tLDsyZmsmx/8QITx0bZu+xYZ7qG+bY2alzn69vqWP7uhXc1NnIjZ2NbFndQIO+rxApmYK9RMlYlE+9+zb+9d/t4iNfeZqIGb91W1e5y1p2+bzzwpkJ9hwdZs/RIfYeG+bgqTFyxS81O5tq2NrdyLvuWMPNXYUgV4iLLA0FewBS8Sh/895e3v+5n/MH//AUz/WP8eE3bgr1pUFHp2d56tgwu48Ms7sY5HPXO0knY2ztbuLfvnoDt3Q3cXN3I23pVJkrFqkeCvaApOJRPvuvbufP/+lZ/vr7h9l/fISPP7CNlfXJcpdWMnfn8OAEu48MsfvoMLuPDHGofwz3wkWqrmtLc99NHWzrbmZbTxMbWuu1JlykjLSOfQn8/a5j/MlX99NcWzhT9c03r6qoL/8mMtnCbPxoIcj3HB1iaLIwG29IxdjW08ytPc3cuqaJrd1NaqmILBOtYy+j+3u72bKqgT/88tN88It7+OLPjvJnb7mB65bhllhXKpd3Dg+Ms+fY8Ln++KHTY+fuiLOhtY43bGnntjWFMNdsXOTaV9KM3czeDvwZcD2w3d0XNQ0P+4x9Ti7vfOGJI/y37/yCsUyW11/fzu/cuY5Xrl9Rlhm8e+HMzX3HR9jXN3JuueF4pnDBq4ZUoTdemJE3sa27mcZazcZFrhXLNWPfD/wG8NcljhNK0Yjxnlet5b6bVvHQj1/k808c5dFnT7OpPc29N3Xwhi3tbFnVsCQhPzY9y3P94zx3eowDJ8d4tnj25ljxqoWJaITNq9L8xq2dbO0qtFTWt9RpNi4SAoH02M3se8AfaMZ+adOzOb629zh/v6uP3UeHcIeOhhS3dDdxY2cDW1Y30N1cS3tjinQydtHAd3fGM1mGJ2c5NTrNieHCjR6OnJnkxcHCFQznn/RTE49y/ar0uRN+bu5q5Lr2NIlYeFftiISReuzXoFQ8ym/f3sNv397D4HiG7x7o5/HnBth/fIRvP3PqvG3nbhScjBXuup7LOzPZPDPZPKPTs+ddcnZOc22ctS11vGr9Sta31nFde5pNHWm6mmt17XiRKnLZYDez/wt0LPDRH7v71xa7IzPbAewA6OnRqfct9Unuv7373M2xR6dnOXRqjBMj05wamaJ/NMPkbI7p2RyZbJ5YxEhEIyRiERpq4jTXxmmqTdDekGJ1Y4pVTTXU69opIsIigt3dXx/Ejtx9J7ATCq2YIMYMk4ZUnN61K8pdhoiEgJqsIiIhU1Kwm9nbzKwPeBXwDTP7TjBliYjI1SqpKevujwCPBFSLiIgEQK0YEZGQUbCLiISMgl1EJGQU7CIiIaNgFxEJmbJcj93MBoAjy77jq9MCDJa7iCUS5mODcB+fjq1ylXJ8a9y99XIblSXYK4mZ7VrMRXcqUZiPDcJ9fDq2yrUcx6dWjIhIyCjYRURCRsF+eTvLXcASCvOxQbiPT8dWuZb8+NRjFxEJGc3YRURCRsEOmNlnzazfzPZf5HMzs4+b2fNm9rSZ3brcNV6tRRzbu4rH9LSZ/djMti53jaW43PHN2+52M8uZ2W8tV22lWsyxmdlrzGyvmT1jZt9fzvpKsYg/l41m9k9m9lTx2N633DVeLTPrNrPHzOxAsfYHF9hmSTNFwV7wEHDPJT6/F9hYfOwAPrUMNQXlIS59bC8Ar3b3m4G/oPL6mw9x6ePDzKLAfwUq7bLSD3GJYzOzJuCTwFvc/Qbg7ctUVxAe4tK/bx8AnnX3rcBrgP9hZollqCsIWeD33f164JXAB8xsywXbLGmmKNgBd38cOHuJTd4K/J0X/BRoMrNVy1NdaS53bO7+Y3cfKv74U6BrWQoLyCJ+7wA+CHwF6F/6ioKziGN7J/Cwux8tbl8xx7eIY3MgbYU7utcXt80uR22lcveT7r67+HoMOAB0XrDZkmaKgn1xOoFj837u4+W/UWHwfuBb5S4iSGbWCbwN+Kty17IErgOazex7Zvakmb233AUF6BPA9cAJYB/woLvny1vSlTOztcA24IkLPlrSTNHdjxfHFngvVMuJzOy1FIL9rnLXErC/BD7i7rnC5C9UYsBtwOuAGuAnZvZTdz9U3rIC8SZgL3A3sAF41Mx+4O6j5S1r8cysnsK/FD+0QN1LmikK9sXpA7rn/dxFYSYRCmZ2M/Bp4F53P1PuegLWC3ypGOotwH1mlnX3r5a3rED0AYPuPgFMmNnjwFYgDMH+PuC/eGE99vNm9gKwGfhZectaHDOLUwj1z7v7wwtssqSZolbM4vwj8N7iN9mvBEbc/WS5iwqCmfUADwPvCclM7zzuvs7d17r7WuDLwO+GJNQBvgb8qpnFzKwWuINCPzcMjlL4lwhm1g5sAg6XtaJFKn4v8BnggLt/7CKbLWmmaMYOmNkXKXzz3lK8OfefAnEAd/8r4JvAfcDzwCSF2URFWMSx/XtgJfDJ4qw2W0kXYFrE8VWsyx2bux8ws28DTwN54NPufslln9eKRfy+/QXwkJnto9C2+Ii7V8oVH+8E3gPsM7O9xff+HdADy5MpOvNURCRk1IoREQkZBbuISMgo2EVEQkbBLiISMgp2EZGQUbCLiISMgl1EJGQU7CIiIfP/Adi0KPLsCSZOAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "r = np.linspace(0.9, 2, 100)\n",
    "y = 4 * (r**(-12) - r**(-6))\n",
    "plt.plot(r, y)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "the minimum is at [1.12246208] and its value is -0.9999999999999768\n"
     ]
    }
   ],
   "source": [
    "result = minimize(lambda r: 4 * (r**(-12) - r**(-6)), x0=1)\n",
    "print(f'the minimum is at {result.x} and its value is {result.fun}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "Maximizing\n",
    "---\n",
    "\n",
    "In order to maximzie a function, just add a minus sign"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "Maximize the following function:\n",
    "\n",
    "$$\n",
    "-\\left[x - \\cos(x)\\right]^2\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "the maximum is at [0.7390852]\n"
     ]
    }
   ],
   "source": [
    "#place - sign to make it a maxmimization problem\n",
    "result = minimize(lambda x: (x - np.cos(x))**2, x0=1)\n",
    "print(f'the maximum is at {result.x}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "Multiple Dimensions\n",
    "----\n",
    "\n",
    "Just indicate multiple dimensions by using a multidimensional `x0`. Note that your $x$ becomes a vector!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "Minimize the following:\n",
    "\n",
    "$$\n",
    "f(x, y) = \\frac{(x - 4)^2} { 2 } + \\frac{(y + 3)^2} { 5 }\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "maximum occurs when x = 3.999999680451366 and y = -3.0000083683518133\n"
     ]
    }
   ],
   "source": [
    "result = minimize(lambda x: (x[0] - 4)**2 / 2 + (x[1] + 3)**2 / 5, x0=[0,0])\n",
    "print(f'maximum occurs when x = {result.x[0]} and y = {result.x[1]}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Identifying Convexity for Minimization - Example\n",
    "----\n",
    "\n",
    "The best ways to identify convexity are:\n",
    "\n",
    "* plot it\n",
    "* try optimizing in multiple starting positions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "scrolled": false,
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYYAAAD8CAYAAABzTgP2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl8FfW9//HXJwkJO4Q9IQQIixB2iIiiopRNtICKiq0tWinqrbVXb/so1tve1t7eauvP7daqVBGsXsXdqCgCsqgIEmTfQ9hCEhJCQoCQ/fP74wx9HCAbnJMz5+R8no9HHjkz852ZD5PDeZ+Z+c6MqCrGGGPMGRFuF2CMMSa4WDAYY4w5iwWDMcaYs1gwGGOMOYsFgzHGmLNYMBhjjDmLBYMxxpizWDAYY4w5iwWDMcaYs0S5XcDF6NChg/bo0cPtMowxJqSsX7/+qKp2rKtdSAZDjx49SEtLc7sMY4wJKSJyoD7t7FCSMcaYs1gwGGOMOYsFgzHGmLNYMBhjjDmLBYMxxpiz+CUYRGSeiOSKyNYapouIPCsi6SKyWUSGe02bKSJ7nJ+Z/qjHGGPMxfPXHsN8YFIt068D+jg/s4HnAUSkHfBfwGXASOC/RCTWTzUZY4y5CH4JBlVdBRyrpclU4FX1WAO0FZE4YCKwRFWPqWoBsITaA8YYY8LSgfxTPLF4F7lFJQ2+rkCdY+gKHPIaznTG1TT+PCIyW0TSRCQtLy+vwQo1xphgtGhLDn9bnk55lTb4ugIVDFLNOK1l/PkjVeeqaoqqpnTsWOcV3cYY06h8vj2HQV3b0LVtswZfV6CCIRPo5jWcAGTVMt4YY4zjSFEJGw4WMiG5c0DWF6hgSAV+7PROGgUcV9VsYDEwQURinZPOE5xxxhhjHEu2HwFgwoAuAVmfX26iJyJvANcAHUQkE09PoyYAqvoCsAiYDKQDxcBdzrRjIvJHYJ2zqEdVtbaT2MYYE3Y+336EHu2b07dzy4Cszy/BoKq31zFdgZ/VMG0eMM8fdRhjTGNTVFLON3uPctfonohUd1rW/+zKZ2OMCWLLd+ZSXqlMHBCY8wtgwWCMMUHt821H6NAyhqHdAnftrwWDMcYEqZLySpbvymV8cmciIwJzGAksGIwxJmit2JVHcVkl1w+KC+h6LRiMMSZILdqSTWzzJoxKahfQ9VowGGNMECopr2TZjiNMHNCFqMjAflRbMBhjTBD6cs9RTpVVcl2ADyOBBYMxxgSlRVuyadOsCVf0ah/wdVswGGNMkCmtqGTp9iNMSO5MkwAfRgILBmOMCTpf7TnKidIKJrtwGAksGIwxJugs2pJDq6ZRjO7dwZX1WzAYY0wQKauoYsn2HMYndyY6yp2PaAsGY4wJIl/vPUpRSUXAL2rzZsFgjDFB5KNNWbRqGsWVfdw5jAQWDMYYEzRKyiv5fNsRrhvYhZioSNfq8EswiMgkEdklIukiMqea6U+JyEbnZ7eIFHpNq/SaluqPeowxJhR9sTOXk6UVTB3a1dU6fH5Qj4hEAs8B4/E8w3mdiKSq6vYzbVT1Qa/2PweGeS3itKoO9bUOY4wJdR9uPEzHVjGMSgr8RW3e/LHHMBJIV9UMVS0D3gSm1tL+duANP6zXGGMajeOny1m+M48bBscF9Bbb1fFHMHQFDnkNZzrjziMi3YGewBdeo5uKSJqIrBGRaX6oxxhjQs7ibTmUVVa5fhgJ/PPM5+qiTWtoOwN4R1UrvcYlqmqWiCQBX4jIFlXde95KRGYDswESExN9rdkYY4JK6sYsurdvzpCENm6X4pc9hkygm9dwApBVQ9sZnHMYSVWznN8ZwArOPv/g3W6uqqaoakrHjh19rdkYY4JG7okSVu89ypQh8Yi4exgJ/BMM64A+ItJTRKLxfPif17tIRC4BYoFvvMbFikiM87oDMBrYfu68xhjTmH2yOZsqhalD490uBfDDoSRVrRCR+4HFQCQwT1W3icijQJqqngmJ24E3VdX7MFN/4EURqcITUo9592Yyxphw8OHGLJLjWtO7Uyu3SwH8c44BVV0ELDpn3O/OGf59NfOtBgb5owZjjAlF+4+eYuOhQuZc18/tUv7Frnw2xhgXvftdJhEC04KgN9IZFgzGGOOSqirlve8Oc2WfjnRp09Ttcv7FgsEYY1yyJiOfw4WnmT4iwe1SzmLBYIwxLnlnfSatmkYxIbmz26WcxYLBGGNccLK0gk+35nDD4HiaNnHvTqrVsWAwxhgXLNqSzenyyqA7jAQWDMYY44p31meS1KEFwxPbul3KeSwYjDEmwA7mF/PtvmPcPCIhKG6BcS4LBmOMCbB3v8tEBG4cFjzXLnizYDDGmACqrFLeWZ/Jlb07EN+2mdvlVMuCwRhjAmjVnjwOF57m9pHB+/gACwZjjAmgN9YepEPLaMb1D65rF7xZMBhjTIAcKSph2c5cpo/oRnRU8H78Bm9lxhjTyLyddojKKmXGpd3qbuwiCwZjjAmAqirljW8PMbp3e3p0aOF2ObWyYDDGmAD4Mv1o0J90PsMvwSAik0Rkl4iki8icaqbfKSJ5IrLR+ZnlNW2miOxxfmb6ox5jjAk2b6w9SPsW0UxI7uJ2KXXy+QluIhIJPAeMBzKBdSKSWs0jOheq6v3nzNsO+C8gBVBgvTNvga91GWNMsMgtKmHpjiPcfWXPoD7pfIY/KhwJpKtqhqqWAW8CU+s570Rgiaoec8JgCTDJDzUZY0zQWLjuEBVVym1BftL5DH8EQ1fgkNdwpjPuXDeLyGYReUdEzmyd+s5rjDEhqbyyitfWHuCqPh1I6tjS7XLqxR/BUN0doPSc4Y+AHqo6GFgKLLiAeT0NRWaLSJqIpOXl5V10scYYE0ifbc3hSFEpd43u4XYp9eaPYMgEvPePEoAs7waqmq+qpc7gP4AR9Z3XaxlzVTVFVVM6duzoh7KNMabhLVi9n+7tm3NN305ul1Jv/giGdUAfEekpItHADCDVu4GIxHkNTgF2OK8XAxNEJFZEYoEJzjhjjAl5Ww8fJ+1AAT++vAcREcF3e+2a+NwrSVUrROR+PB/okcA8Vd0mIo8CaaqaCjwgIlOACuAYcKcz7zER+SOecAF4VFWP+VqTMcYEg/mr99M8OpJbUoLvKW218TkYAFR1EbDonHG/83r9MPBwDfPOA+b5ow5jjAkW+SdLSd2Uxa0pCbRu2sTtci5I8HeoNcaYEPTmukOUVVQx8/IebpdywSwYjDHGz8orq3htzQGu7N2BPp1buV3OBbNgMMYYP/tkczbZx0tCqouqNwsGY4zxI1XlxVUZ9O7UkmsvCZ0uqt4sGIwxxo++Sj/KjuwiZl+VFFJdVL1ZMBhjjB/NXZVBp1YxTB0W73YpF82CwRhj/GRb1nG+3HOUO0f3ICYq0u1yLpoFgzHG+MncVRm0iI7kh5d1d7sUn1gwGGOMH2QWFPPx5mxuH5lIm2ahdUHbuSwYjDHGD/6xKgMBfnJlT7dL8ZkFgzHG+Ci3qIQ31h3i5uEJxLdt5nY5PrNgMMYYH724KoPKKuXfru3ldil+YcFgjDE+OHqylNfXHmDq0Hi6t2/hdjl+YcFgjDE++MeXGZRWVPGza3u7XYrfWDAYY8xFOnaqjH9+c4DvD46nV4g8z7k+LBiMMeYizftqH8Vlldw/tvHsLYCfHtQjIpOAZ/A8we0lVX3snOkPAbPwPMEtD/iJqh5wplUCW5ymB1V1ij9qMg1LVck/VcahY8VkFpymsLiM46fLOVFSQWWVEhEhCBDTJJLY5k2IbR5NuxbRdGvXnITYZjSJtO8kJrTlnyzlla/3MXlQF/qG4K21a+NzMIhIJPAcMB7IBNaJSKqqbvdqtgFIUdViEbkP+AtwmzPttKoO9bUO07COnixl9d58Nh8qZPPh42zPKuJkacV57WKiIoiMEFShSpXSiqrz2kQIdI1txiWdWzOoaxsGJ7RhUEIbOrSMCcQ/xRi/+PuKvZwur+Sh8X3dLsXv/LHHMBJIV9UMABF5E5gK/CsYVHW5V/s1wB1+WK9pQKrK9uwiPtmczcrdeWzLKgI8H/zJ8a25cVhXkjq2ILFdcxJim9O+ZTStmkadd3+Yyirl+OlyCovLOHqyjIPHijmQf4p9R0+xPbuIZTuPoOpp26dTS0b37sAVvdpzea/2tAqxxyGa8JFVeJp/rjnAzcMT6N2pce0tgH+CoStwyGs4E7islvZ3A596DTcVkTQ8h5keU9UPqptJRGYDswESExN9KtjULLeohLfSDvHBxizSc08SFSEM7x7LryZewpW9O5Ac3/qCDgNFRgjtWngOIyV1hJE92501/WRpBdsOH2fDoUJW783nzXUHmb96P00ihdG9OzBpQBfGJXe2vQkTVJ5ZugcUfjGuj9ulNAh/BEN1NxzXahuK3AGkAGO8RieqapaIJAFfiMgWVd173gJV5wJzAVJSUqpdvrl4mzMLeeXr/Xy8OYvySmVkj3b86caBTB4YR2yL6AZbb8uYKC5Las9lSe25d0wvSisq2XCwkC925vLp1mzmvLeFiPe3cHmv9kwfkcCkAXE0iw7du1aa0Lc37yTvfJfJj0Z1JyG2udvlNAh/BEMm0M1rOAHIOreRiIwDHgHGqGrpmfGqmuX8zhCRFcAw4LxgMA1jc2YhT3y+m1W78/51V8g7r+hBjw7uXKgTExXJqKT2jEpqz8PX9WN7dhGLt+bw/sbDPLhwE7+L2cYNQ+K4NaUbQ7u1RSQ0H4RiQteTS3YTExXR6HoiefNHMKwD+ohIT+AwMAP4gXcDERkGvAhMUtVcr/GxQLGqlopIB2A0nhPTpoHtO3qKxz/dyWfbcmjbvAlzruvHDy5LpHUQHdcXEQbEt2FAfBv+fVxf1u47xtvrD/HBhize+PYQQxLa8JMre3LdwDiio6yXk2l4mw4V8snmbH4+tnejPrwpqr4flRGRycDTeLqrzlPVP4nIo0CaqqaKyFJgEJDtzHJQVaeIyBV4AqMKzzUVT6vqy3WtLyUlRdPS0nyuOxyVlFfy9+XpvLAyg+ioCGZd1ZO7r+wZUid6T5ZW8P53mbyyej8Zeafo1CqGH43qzg9HdaddAx72MuFNVbnlhW/Yn3+K5b+8JqT+z5whIutVNaXOdv4IhkCzYLg4X+05ysPvb+bQsdNMGxrPbyb3p1Prpm6XddGqqpRVe/KY9/V+Vu3Oo3l0JHeM6s6sq3rSqVXo/rtMcPp4cxb3/98G/nzTIG4fGZodYCwYzL+UlFfy2Kc7mb96P0kdW/Df0wZyRa8ObpflV3uOnOC55emkbsqiSWQEt49M5N4xvejSxgLC+K6kvJJxT66kZUwUnzxwFZERoXluq77B4Jcrn03w2pFdxM/f2EB67knuvKIHv57Ur1H26unTuRVPzxjGL8b15e/L03ltzQH+b+1BfnBZIvc38uPBpuHN+3ofmQWneX3WZSEbChfC9hgasQ83HubX726mddMmPHHLEK7u29HtkgLm0LFinlueztvrM2kaFcGsq5L46dVJtIyx70LmwuSeKGHsEysZldSOl2Ze6nY5PqnvHoN15WiEyiur+MNH2/jFmxsZnNCWTx64KqxCAaBbu+Y8dvNgPn/waq7u25Fnlu1hzF+W88rX+yir5jYdxtTk8U93UVJeyW8m93e7lICxYGhkTpSU85P563jl6/38ZHRPXp91GR1bhe9hlF4dW/L8HSN4/9+uoE/nlvzho+1MeGolS7cfIRT3lk1grc3I593vMpl1VRJJjei22nWxYGhEjhSVcNuLa/hmbz5/mT6Y330/2e5i6hiWGMsbPx3FK3ddSmSEMOvVNGa+so703BNul2aCVHllFf/5wVa6tm3GA99rvBezVcc+NRqJ9NyT3PT31RzIP8XLd17KrSnd6p4pzIgI117Sic/+/Wp+e0MyGw4WMPHpL/nDR9s4frrc7fJMkHn5q33syT3J76cMoHl0eJ2bsmBoBHblnGDG3G8orahi4T2XMybMzidcqCaREdx9ZU9W/PIabk3pxvzV+7n2iRW88e1Bqqrs8JKBzIJinlm6h3H9OzM+ubPb5QScBUOI25FdxO3/WENkhPDWPaMY2LWN2yWFjPYtY/jzTYP4+OdX0rtjSx5+bwvTX1jNjuwit0szLlJVfp+6DYDfT0l2uRp3WDCEsG1Zx7n9H2uIiYpg4ezLw+rkmD8NiG/DwntG8f9uGcL+/GJu+N+v+J9FOyguO/9BRKbxS92UxdIduTw4vk+jvXtqXSwYQlRG3kl+/PK3NG8SyZuzR7l2N9TGQkS4eUQCyx4awy0jEpi7KoPxT65iyfYjbpdmAij3RAn/lbqNYYltufvKJLfLcY0FQwjKOV7Cj17+FoDXZl1G9/YWCv4S2yKax24ezDv3Xk7LmCh++moaP301jcOFp90uzTQwVeWR97dSXFbJX6cPCYsrnGtiwRBiCovL+NHLazl+upz5d420w0cNJKVHOz5+4ErmXNePL/fkMf7Jlcz7ah+VdnK60UrdlMWS7Uf45YS+9O4U3v+vLBhCSGlFJT99NY0D+cXM/fEIBiXYieaG1CQygnvH9GLJg2MY2bMdj368nekvrGb3Ebv2obE5UuQ5hDQ8zA8hnWHBECJUld+8t5V1+wt44tYhje7uqMGsW7vmvHLnpTx12xD2Hz3F9c9+ydNLd9utNRqJyirlwYUbKS2v4q+3hPchpDP8EgwiMklEdolIuojMqWZ6jIgsdKavFZEeXtMedsbvEpGJ/qinMXpxVQbvfpfJL77XhylD4t0uJ+yICDcOS2DpQ2O4bmAcTy/dww3/+yXfHSxwuzTjoxdX7WX13nx+PyWZXnZoFvBDMIhIJPAccB2QDNwuIud2/r0bKFDV3sBTwOPOvMl4HgU6AJgE/N1ZnvHy+bYcHv9sJzcMjuPfx/Vxu5yw1r5lDM/ePox5d6ZwoqSCm59fzR8+2sapUuvaGoo2Hirkyc93c/3gOLtbgBd/7DGMBNJVNUNVy4A3ganntJkKLHBevwN8TzxPcZ8KvKmqpaq6D0h3lmcc+46e4qG3NjGoaxueuGUIns1m3Da2X2c+f/Bq7risO698vZ8JT61i5e48t8syF6CopJwH3thA59ZN+Z8bB9n/LS/+CIauwCGv4UxnXLVtVLUCOA60r+e8YaukvJL7XltPVKTw/B0jaNrEdqaCSaumTfjjtIG8fe/lxDSJYOa8b3lo4UYKTpW5XZqpQ1WV8h9vbSKr8DTPzBhKm2ah9/zmhuSPYKguZs/t01dTm/rM61mAyGwRSRORtLy88Phm9rsPt7LryAmevm0oXds2c7scU4NLe7Rj0QNXcf+1vUndlMW4J1eSuinLbusdxJ5fuZcl24/wm8n9SenRzu1ygo4/giET8D44lwBk1dRGRKKANsCxes4LgKrOVdUUVU3p2LHx3yTurXWHeCstk59f25trLunkdjmmDk2bRPLLiZeQev+VdI1txgNvbODuBXZhXDBauTuPJz7fxdSh8dw1uofb5QQlfwTDOqCPiPQUkWg8J5NTz2mTCsx0Xk8HvlDP16lUYIbTa6kn0Af41g81hbT03BP89sOtjO7dnl+M6+t2OeYCJMe35r37ruA/r+/PN3vzmfDkSuZ/bRfGBYsD+af4xZsbuKRzK/58k51XqInPweCcM7gfWAzsAN5S1W0i8qiITHGavQy0F5F04CFgjjPvNuAtYDvwGfAzVa30taZQVlZRxS/e3EiLmCieum2o9akOQVGRnmdMf/7g1Yzo0Y7ff7Sdm59fzc4cu2urmwqLy7hr/joAXvzRiLB7xsKFkFA8DpqSkqJpaWlul9Eg/rp4J88t38uLPxrBxAFd3C7H+EhV+XBjFo9+vJ2i0+XcO6YX94/tbR0JAqy0opIfv/wtGw4W8tqsyxjZMzzPK4jIelVNqaudXfkcRNL2H+P5FXu5NSXBQqGREBGmDevK0ofGMGVIPH9bns7kZ75kbUa+26WFDVXl4Xe3sHbfMf56y+CwDYULYcEQJE6UlPPgWxtJiG3O774/wO1yjJ+1axHNk7cN5dWfjKSssorb5q7h4fe22CNFG5iq8thnO3lvw2EeGt+XqUOtN3x9WDAEicc/28nhgtM8eesQWsbYsc/G6uq+Hfn8wav56VU9WbjuIOOfXMmnW7Kta2sD+fuKvby4MoMfXpbIz8f2druckGHBEAS+3XeM19Yc5K7RPa1PdRhoHh3FI9cn88HPRtOhZQz3vf4dP5m/jgP5p9wurVGZ//U+/rp4FzcO68ofpw60HkgXwILBZSXllcx5bzMJsc34jwnWNTWcDE5oy4f3j+aRyf35dt8xxj+1iieX7KakPKw75vnF62sP8PuPtjMhuTN/nT6YCOvdd0EsGFz23PJ0MvJO8T83DrLuc2GoSWQEP706iS9+eQ2TBnTh2WV7GP/USpbtsEeKXqyXvszgkfe3MrZfJ/73B8OIirSPuQtlW8xFO3OKeH7FXm4a3pWr+zb+q7lNzTq3bsqztw/j/356GTFRkdy9II1ZC9Zx6Fix26WFDFXlmaV7+O9PdnD9oDheuGMEMVHWLfhiWDC4pKpKefi9LbRp1oTfXn/uXcpNuLqiVwcWPXAVv5ncj9V78xn35Er+8tlOTpRY76XaVFYpf/hoO08t3c1Nw7vyzIyhREfZx9vFsi3nkne/y2TDwUIentyf2BbRbpdjgkh0VASzr+7Fsv8Yw3UDu/D3FXu59okVvLbmABWV9tS4c50qrWD2q2nMX72fu6/syRPTh9jhIx/Z1nNBUUk5j3+2k+GJbblpmPWrNtWLa9OMp2cM48OfjSapQ0v+84OtXPfMlyzfmWvdWx1Zhae55YVvWL4rlz9OG8hvb0i2E81+YMHggqeX7CH/VBmPTh1ob2JTpyHd2rLwnlG8cMcIyiuruGv+Ou54eS0bwvyxoit25XL9s19y8FgxL995KT8a1d3tkhoN6wYTYLtyTrDgm/38YGQiA7u2cbscEyJEhEkDuzC2XydeW3OAvy1P58a/r2Zsv048NL5vWL2XKquUp5fu5m/L07mkcyue++Fwe1azn9lN9AJIVfnBP9ayI6eI5f9xjZ1bMBftZGkFC1bvZ+6qDI6fLmfigM48OL4v/bq0dru0BpWRd5JfvbOZ9QcKuGVEAo9OHUizaOt5VF/1vYme7TEE0JLtR/gmI58/ThtooWB80jImip9d25sfXd6dl7/cx7yv9rF425eMT+7MPVcnNbor6CurlFecK5ljoiJ46rYh3Dgswe2yGi3bYwiQ8soqJj69iggRPvvFVdZrwvhVYXEZr3y9n1e/2U9BcTkjusdyz9VJjOvfOeTPY60/UMAfPtrG5szjjOvfif+5cRCdWjd1u6yQZHsMQWbhukNk5J3ipR+nWCgYv2vbPJoHx/flnjFJvJ2WyUtfZTD7n+vp0b45P7gskekjutEuxPZSswpP85fPdvLBxiw6t47hmRlDmTIk3u55FAA+7TGISDtgIdAD2A/cqqoF57QZCjwPtAYqgT+p6kJn2nxgDHDcaX6nqm6sa72htsdwsrSCa/66nKSOLVk4e5S9sU2Dq6is4tOtObz6zX7W7S8gOiqC6wfFcfvIRC7tERvU78FDx4p5fuVe3k47hIgw+6ok7rumFy3srsM+C9Qewxxgmao+JiJznOFfn9OmGPixqu4RkXhgvYgsVtVCZ/qvVPUdH+sIanNXZXD0ZBkvzewf1P8hTeMRFRnB94fE8/0h8ezKOcH/rT3Ae98d5v0Nh+nathnfHxLP1KHx9OvSKijek6rKdwcLeX3NAVI3ZREhwm2XduPeMb1IiG3udnlhx9c9hl3ANaqaLSJxwApVvaSOeTYB052gmA98fKHBEEp7DLlFJYz56wrG9u/Ecz8Y7nY5JowVl1Xw2dYcUjdl8eWeo1RWKX06tWRccmeuvaQTwxPbBvww55GiEj7dks3CtEx2ZBfRMiaK6SMSuGdMEnFtmgW0lnBQ3z0GX4OhUFXbeg0XqGpsLe1HAguAAapa5QTD5UApsAyYo6qlda03lILhtx9s5c11B1n60Bi6t2/hdjnGAJB/spRFW3P4ZHMWafsLqKhSWjeN4qo+HRnZsx0jusfSr0srvwdFZZWyLes43+zNZ+mOI6QdKEAVkuNa88NRiUwb2tUOGTUgvwWDiCwFqnsA8SPAgvoGw5k9CmCmqq7xGpcDRANzgb2q+mgN888GZgMkJiaOOHDgQO3/siCQWVDMtU+s4NaUbvzpxkFul2NMtYpKyvl6z1GW78pl1e6j5BSVANA8OpKB8W3o07klvTu1pE+nViTENqNT65g6bxFfVaXknyoj53gJe/NOsiOniJ3ZJ/juYAEnSioA6NelFZMHxTF5UBy9O9kFaoEQqD2Geh1KEpHWeELhz6r6dg3Lugb4pareUNd6Q2WP4eH3NvPu+sOs+NU1xLe13WIT/FSVrOMlrD9QwHcHCthy+DjpuSfPezZ1y5goYls0ISYqkpioCJpERlBaUUVJeSXFZRUUnCqnzOuGf9GREfTu1JIh3dowKqk9lye1ty6nLgjUyedUYCbwmPP7w2oKiQbeB149NxREJM4JFQGmAVt9rCdoHMwv5u20TO4Y1d1CwYQMEaFr22Z0bduMKUPiAU9Y5J0sJT33JNmFJeSeKCX3RAmFxeWUVVRRWlFJaUUVHaIiaR4dSdMmEbRrEUNcm6Z0adOUHu1bkNSxBU2sm3bI8DUYHgPeEpG7gYPALQAikgLcq6qzgFuBq4H2InKnM9+Zbqmvi0hHQICNwL0+1hM0nv1iD5ERwn3X9HK7FGN8IiJ0atWUTq3sG3648CkYVDUf+F4149OAWc7r14DXaph/rC/rD1YZeSd577tM7hrdk862u2yMCTG2b9cAnl22h5ioSO4dY3sLxpjQY8HgZwfyT5G6KYs7RiXSsVWM2+UYY8wFs2DwsxdWZhAVEcGsq5LcLsUYYy6KBYMfHSkq4d31mUxPSbBzC8aYkGXB4Ecvf7WPiqoq7rna9haMMaHLgsFPjheX8/qaA9wwON5ufWGMCWkWDH6y4Jv9nCqrtOsWjDEhz4LBD4rLKnjl6318r18n+sc17mfuGmMaPwsGP1i47hAFxeX827W2t2CMCX0WDD7yPKR8PyNA0J9JAAAOjklEQVS6xzKie+N6ALsxJjxZMPho6Y4jHDxWzN1X9nS7FGOM8QsLBh+9/NU+urZtxoTkzm6XYowxfmHB4IOth4/z7b5j3HlFj4A/EtEYYxqKfZr5YN5X+2gRHcltI7u5XYoxxviNBcNFyi0q4aPNWdyS0o3WTZu4XY4xxviNBcNF+ueaA1RUKXeN7uF2KcYY41c+BYOItBORJSKyx/kdW0O7ShHZ6Pykeo3vKSJrnfkXOo8BDXol5ZW8vvYg4/p3tttfGGMaHV/3GOYAy1S1D7DMGa7OaVUd6vxM8Rr/OPCUM38BcLeP9QTEx5uzOXaqzPYWjDGNkq/BMBVY4LxeAEyr74wiIsBY4J2Lmd9Nr605QK+OLbg8qb3bpRhjjN/5GgydVTUbwPndqYZ2TUUkTUTWiMiZD//2QKGqVjjDmUDXmlYkIrOdZaTl5eX5WPbF23r4OBsPFXLHqO54ss0YYxqXqLoaiMhSoEs1kx65gPUkqmqWiCQBX4jIFqComnZa0wJUdS4wFyAlJaXGdg3t9bUHaNYkkpuGJ7hVgjHGNKg6g0FVx9U0TUSOiEicqmaLSByQW8MyspzfGSKyAhgGvAu0FZEoZ68hAci6iH9DwBSVlPPBhiymDImnTTPromqMaZx8PZSUCsx0Xs8EPjy3gYjEikiM87oDMBrYrqoKLAem1zZ/MHlvfSanyyu5Y1R3t0sxxpgG42swPAaMF5E9wHhnGBFJEZGXnDb9gTQR2YQnCB5T1e3OtF8DD4lIOp5zDi/7WE+DUVVeW3uQIQltGJTQxu1yjDGmwdR5KKk2qpoPfK+a8WnALOf1amBQDfNnACN9qSFQ1u47RnruSf4yfbDbpRhjTIOyK5/r6fW1B2ndNIrvD453uxRjjGlQFgz1UHCqjMVbc7hpeALNoiPdLscYYxqUBUM9fLjxMGWVVdyaYndRNcY0fhYMdVBVFqZlMqhrG5LjW7tdjjHGNDgLhjpsyypiR3YRt6bYBW3GmPBgwVCHt9IOER0VwZQhNd6twxhjGhULhlqUlFfywYbDTBrQhTbN7UpnY0x4sGCoxefbj1BUUsFtl9pJZ2NM+LBgqMXbaYfo2raZ3V7bGBNWLBhqkFlQzFfpR7klJYGICLu9tjEmfFgw1OC97w6jCjfb7bWNMWHGgqEaqsoHGw4zKqkd3do1d7scY4wJKAuGamzOPE7G0VPcOMy6qBpjwo8FQzXe33CY6KgIJg2Mc7sUY4wJOAuGc5RXVvHRpizG9e9kT2kzxoQlC4ZzfLXnKPmnyrhxmJ10NsaEJ5+CQUTaicgSEdnj/I6tps21IrLR66dERKY50+aLyD6vaUN9qccf3t9wmLbNmzCmb0e3SzHGGFf4uscwB1imqn2AZc7wWVR1uaoOVdWhwFigGPjcq8mvzkxX1Y0+1uOTk6UVfL49hxsGxxEdZTtTxpjw5Oun31RggfN6ATCtjvbTgU9VtdjH9TaIxVtzKCmvst5Ixpiw5mswdFbVbADnd6c62s8A3jhn3J9EZLOIPCUiMTXNKCKzRSRNRNLy8vJ8q7oGH2w8TLd2zRieeN4RMWOMCRt1BoOILBWRrdX8TL2QFYlIHDAIWOw1+mGgH3Ap0A74dU3zq+pcVU1R1ZSOHf1//D+3qISv049y49CuiNgtMIwx4SuqrgaqOq6maSJyRETiVDXb+eDPrWVRtwLvq2q517KznZelIvIK8Mt61u13n27NoUphytB4t0owxpig4OuhpFRgpvN6JvBhLW1v55zDSE6YIJ6v6NOArT7Wc9E+3pxFvy6t6N2plVslGGNMUPA1GB4DxovIHmC8M4yIpIjIS2caiUgPoBuw8pz5XxeRLcAWoAPw3z7Wc1Gyj59m3f4Crh9kVzobY0ydh5Jqo6r5wPeqGZ8GzPIa3g+c19VHVcf6sn5/WbQlB4DrB1swGGOMddbHcxgpOa41SR1bul2KMca4LuyDIbOgmA0HC7lhiO0tGGMMWDCwaIunY9QNg6w3kjHGgAUDn2zOZnBCGxLb2wN5jDEGwjwYDuYXsynzODfYSWdjjPmXsA6Gj7dkATDZuqkaY8y/hHUwLNqSzdBubUmItcNIxhhzRtgGQ2ZBMVsPF3HdwC5ul2KMMUElbIPh821HAJg4wILBGGO8hW0wfLYth35dWtGjQwu3SzHGmKASlsFw9GQpafuPMcH2Fowx5jxhGQxLtx+hSmHigM5ul2KMMUEnLINh8bYcurVrRnJca7dLMcaYoBN2wXCipJyv0/OZmNzFntRmjDHVCLtgWL4rj7LKKiZaN1VjjKmWT8EgIreIyDYRqRKRlFraTRKRXSKSLiJzvMb3FJG1IrJHRBaKSLQv9dTH4m05dGgZw/DE2IZelTHGhCRf9xi2AjcBq2pqICKRwHPAdUAycLuIJDuTHweeUtU+QAFwt4/11KqkvJIVO3MZn9yZyAg7jGSMMdXxKRhUdYeq7qqj2UggXVUzVLUMeBOY6jzneSzwjtNuAZ7nPjeYr9OPcqqs0nojGWNMLQJxjqErcMhrONMZ1x4oVNWKc8Y3mMXbcmgVE8UVvTo05GqMMSak1fnMZxFZClR3pvYRVf2wHuuo7piN1jK+pjpmA7MBEhMT67Ha8/Xs0JI7Lu9OdFTYnXM3xph6qzMYVHWcj+vIBLp5DScAWcBRoK2IRDl7DWfG11THXGAuQEpKSo0BUpv7rul1MbMZY0xYCcRX53VAH6cHUjQwA0hVVQWWA9OddjOB+uyBGGOMaUC+dle9UUQygcuBT0RksTM+XkQWATh7A/cDi4EdwFuqus1ZxK+Bh0QkHc85h5d9qccYY4zvxPPFPbSkpKRoWlqa22UYY0xIEZH1qlrjNWdn2FlYY4wxZ7FgMMYYcxYLBmOMMWexYDDGGHMWCwZjjDFnCcleSSKSBxy4yNk74Lm4LthYXRfG6rowVteFaax1dVfVjnU1Cslg8IWIpNWnu1agWV0Xxuq6MFbXhQn3uuxQkjHGmLNYMBhjjDlLOAbDXLcLqIHVdWGsrgtjdV2YsK4r7M4xGGOMqV047jEYY4ypRaMMBhG5RUS2iUiViNR4Bl9EJonILhFJF5E5XuN7ishaEdkjIgud24X7o652IrLEWe4SEYmtps21IrLR66dERKY50+aLyD6vaUMDVZfTrtJr3ale493cXkNF5Bvn771ZRG7zmubX7VXT+8Vreozz7093tkcPr2kPO+N3ichEX+q4iLoeEpHtzvZZJiLdvaZV+zcNUF13ikie1/pneU2b6fzd94jIzADX9ZRXTbtFpNBrWoNsLxGZJyK5IrK1hukiIs86NW8WkeFe0/y/rVS10f0A/YFLgBVASg1tIoG9QBIQDWwCkp1pbwEznNcvAPf5qa6/AHOc13OAx+to3w44BjR3hucD0xtge9WrLuBkDeNd215AX6CP8zoeyAba+nt71fZ+8Wrzb8ALzusZwELndbLTPgbo6SwnMoB1Xev1HrrvTF21/U0DVNedwN+qmbcdkOH8jnVexwaqrnPa/xyYF4DtdTUwHNhaw/TJwKd4nnw5CljbkNuqUe4xqOoOVd1VR7ORQLqqZqhqGfAmMFVEBBgLvOO0WwBM81NpU53l1Xe504FPVbXYT+uvyYXW9S9uby9V3a2qe5zXWUAuUOcFPBeh2vdLLfW+A3zP2T5TgTdVtVRV9wHpzvICUpeqLvd6D63B87TEhlaf7VWTicASVT2mqgXAEmCSS3XdDrzhp3XXSFVX4fkSWJOpwKvqsQbP0y/jaKBt1SiDoZ66Aoe8hjOdce2BQvU8YMh7vD90VtVsAOd3pzraz+D8N+WfnF3Jp0QkJsB1NRWRNBFZc+bwFkG0vURkJJ5vgXu9Rvtre9X0fqm2jbM9juPZPvWZtyHr8nY3nm+eZ1T3Nw1kXTc7f593ROTMI4CDYns5h9x6Al94jW6o7VWXmupukG1V5zOfg5WILAW6VDPpEVWtzyNCpZpxWst4n+uq7zKc5cQBg/A8+e6Mh4EcPB9+c/E8Ae/RANaVqKpZIpIEfCEiW4Ciatq5tb3+CcxU1Spn9EVvr+pWUc24c/+dDfKeqkO9ly0idwApwBiv0ef9TVV1b3XzN0BdHwFvqGqpiNyLZ29rbD3nbci6zpgBvKOqlV7jGmp71SWg762QDQZVHefjIjKBbl7DCUAWnvuQtBWRKOdb35nxPtclIkdEJE5Vs50PstxaFnUr8L6qlnstO9t5WSoirwC/DGRdzqEaVDVDRFYAw4B3cXl7iUhr4BPgP53d7DPLvujtVY2a3i/VtckUkSigDZ7DA/WZtyHrQkTG4QnbMapaemZ8DX9Tf3zQ1VmXquZ7Df4DeNxr3mvOmXeFH2qqV11eZgA/8x7RgNurLjXV3SDbKpwPJa0D+oinR000njdBqnrO6CzHc3wfYCZQnz2Q+kh1llef5Z53bNP5cDxzXH8aUG0PhoaoS0RizxyKEZEOwGhgu9vby/nbvY/n+Ovb50zz5/aq9v1SS73TgS+c7ZMKzBBPr6WeQB/gWx9quaC6RGQY8CIwRVVzvcZX+zcNYF1xXoNT8DwTHjx7yROc+mKBCZy959ygdTm1XYLnZO43XuMacnvVJRX4sdM7aRRw3Pni0zDbqiHOsLv9A9yIJ0lLgSPAYmd8PLDIq91kYDeexH/Ea3wSnv+46cDbQIyf6moPLAP2OL/bOeNTgJe82vUADgMR58z/BbAFzwfca0DLQNUFXOGse5Pz++5g2F7AHUA5sNHrZ2hDbK/q3i94Dk1NcV43df796c72SPKa9xFnvl3AdX5+v9dV11Ln/8GZ7ZNa1980QHX9GdjmrH850M9r3p842zEduCuQdTnDvwceO2e+BtteeL4EZjvv5Uw854LuBe51pgvwnFPzFrx6WzbEtrIrn40xxpwlnA8lGWOMqYYFgzHGmLNYMBhjjDmLBYMxxpizWDAYY4w5iwWDMcaYs1gwGGOMOYsFgzHGmLP8fyaH21UeGpZNAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "x = np.linspace(-1, 1, 1000)\n",
    "plt.plot(x,  2 * x ** 3 - 0 * x **2  - x)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "scrolled": true,
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "      fun: -0.27216552697590846\n",
       " hess_inv: array([[0.20398928]])\n",
       "      jac: array([-1.11758709e-08])\n",
       "  message: 'Optimization terminated successfully.'\n",
       "     nfev: 24\n",
       "      nit: 5\n",
       "     njev: 8\n",
       "   status: 0\n",
       "  success: True\n",
       "        x: array([0.40824828])"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from scipy.optimize import minimize\n",
    "minimize(lambda x:  2 * x ** 3 - 0 * x **2  - x, x0=0.05)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {
    "scrolled": false,
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYYAAAD8CAYAAABzTgP2AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl8FfW9//HXJwkJO4Q9IQQIixB2iIiiopRNtICKiq0tWinqrbVXb/so1tve1t7eauvP7daqVBGsXsXdqCgCsqgIEmTfQ9hCEhJCQoCQ/fP74wx9HCAbnJMz5+R8no9HHjkz852ZD5PDeZ+Z+c6MqCrGGGPMGRFuF2CMMSa4WDAYY4w5iwWDMcaYs1gwGGOMOYsFgzHGmLNYMBhjjDmLBYMxxpizWDAYY4w5iwWDMcaYs0S5XcDF6NChg/bo0cPtMowxJqSsX7/+qKp2rKtdSAZDjx49SEtLc7sMY4wJKSJyoD7t7FCSMcaYs1gwGGOMOYsFgzHGmLNYMBhjjDmLBYMxxpiz+CUYRGSeiOSKyNYapouIPCsi6SKyWUSGe02bKSJ7nJ+Z/qjHGGPMxfPXHsN8YFIt068D+jg/s4HnAUSkHfBfwGXASOC/RCTWTzUZY4y5CH4JBlVdBRyrpclU4FX1WAO0FZE4YCKwRFWPqWoBsITaA8YYY8LSgfxTPLF4F7lFJQ2+rkCdY+gKHPIaznTG1TT+PCIyW0TSRCQtLy+vwQo1xphgtGhLDn9bnk55lTb4ugIVDFLNOK1l/PkjVeeqaoqqpnTsWOcV3cYY06h8vj2HQV3b0LVtswZfV6CCIRPo5jWcAGTVMt4YY4zjSFEJGw4WMiG5c0DWF6hgSAV+7PROGgUcV9VsYDEwQURinZPOE5xxxhhjHEu2HwFgwoAuAVmfX26iJyJvANcAHUQkE09PoyYAqvoCsAiYDKQDxcBdzrRjIvJHYJ2zqEdVtbaT2MYYE3Y+336EHu2b07dzy4Cszy/BoKq31zFdgZ/VMG0eMM8fdRhjTGNTVFLON3uPctfonohUd1rW/+zKZ2OMCWLLd+ZSXqlMHBCY8wtgwWCMMUHt821H6NAyhqHdAnftrwWDMcYEqZLySpbvymV8cmciIwJzGAksGIwxJmit2JVHcVkl1w+KC+h6LRiMMSZILdqSTWzzJoxKahfQ9VowGGNMECopr2TZjiNMHNCFqMjAflRbMBhjTBD6cs9RTpVVcl2ADyOBBYMxxgSlRVuyadOsCVf0ah/wdVswGGNMkCmtqGTp9iNMSO5MkwAfRgILBmOMCTpf7TnKidIKJrtwGAksGIwxJugs2pJDq6ZRjO7dwZX1WzAYY0wQKauoYsn2HMYndyY6yp2PaAsGY4wJIl/vPUpRSUXAL2rzZsFgjDFB5KNNWbRqGsWVfdw5jAQWDMYYEzRKyiv5fNsRrhvYhZioSNfq8EswiMgkEdklIukiMqea6U+JyEbnZ7eIFHpNq/SaluqPeowxJhR9sTOXk6UVTB3a1dU6fH5Qj4hEAs8B4/E8w3mdiKSq6vYzbVT1Qa/2PweGeS3itKoO9bUOY4wJdR9uPEzHVjGMSgr8RW3e/LHHMBJIV9UMVS0D3gSm1tL+duANP6zXGGMajeOny1m+M48bBscF9Bbb1fFHMHQFDnkNZzrjziMi3YGewBdeo5uKSJqIrBGRaX6oxxhjQs7ibTmUVVa5fhgJ/PPM5+qiTWtoOwN4R1UrvcYlqmqWiCQBX4jIFlXde95KRGYDswESExN9rdkYY4JK6sYsurdvzpCENm6X4pc9hkygm9dwApBVQ9sZnHMYSVWznN8ZwArOPv/g3W6uqqaoakrHjh19rdkYY4JG7okSVu89ypQh8Yi4exgJ/BMM64A+ItJTRKLxfPif17tIRC4BYoFvvMbFikiM87oDMBrYfu68xhjTmH2yOZsqhalD490uBfDDoSRVrRCR+4HFQCQwT1W3icijQJqqngmJ24E3VdX7MFN/4EURqcITUo9592Yyxphw8OHGLJLjWtO7Uyu3SwH8c44BVV0ELDpn3O/OGf59NfOtBgb5owZjjAlF+4+eYuOhQuZc18/tUv7Frnw2xhgXvftdJhEC04KgN9IZFgzGGOOSqirlve8Oc2WfjnRp09Ttcv7FgsEYY1yyJiOfw4WnmT4iwe1SzmLBYIwxLnlnfSatmkYxIbmz26WcxYLBGGNccLK0gk+35nDD4HiaNnHvTqrVsWAwxhgXLNqSzenyyqA7jAQWDMYY44p31meS1KEFwxPbul3KeSwYjDEmwA7mF/PtvmPcPCIhKG6BcS4LBmOMCbB3v8tEBG4cFjzXLnizYDDGmACqrFLeWZ/Jlb07EN+2mdvlVMuCwRhjAmjVnjwOF57m9pHB+/gACwZjjAmgN9YepEPLaMb1D65rF7xZMBhjTIAcKSph2c5cpo/oRnRU8H78Bm9lxhjTyLyddojKKmXGpd3qbuwiCwZjjAmAqirljW8PMbp3e3p0aOF2ObWyYDDGmAD4Mv1o0J90PsMvwSAik0Rkl4iki8icaqbfKSJ5IrLR+ZnlNW2miOxxfmb6ox5jjAk2b6w9SPsW0UxI7uJ2KXXy+QluIhIJPAeMBzKBdSKSWs0jOheq6v3nzNsO+C8gBVBgvTNvga91GWNMsMgtKmHpjiPcfWXPoD7pfIY/KhwJpKtqhqqWAW8CU+s570Rgiaoec8JgCTDJDzUZY0zQWLjuEBVVym1BftL5DH8EQ1fgkNdwpjPuXDeLyGYReUdEzmyd+s5rjDEhqbyyitfWHuCqPh1I6tjS7XLqxR/BUN0doPSc4Y+AHqo6GFgKLLiAeT0NRWaLSJqIpOXl5V10scYYE0ifbc3hSFEpd43u4XYp9eaPYMgEvPePEoAs7waqmq+qpc7gP4AR9Z3XaxlzVTVFVVM6duzoh7KNMabhLVi9n+7tm3NN305ul1Jv/giGdUAfEekpItHADCDVu4GIxHkNTgF2OK8XAxNEJFZEYoEJzjhjjAl5Ww8fJ+1AAT++vAcREcF3e+2a+NwrSVUrROR+PB/okcA8Vd0mIo8CaaqaCjwgIlOACuAYcKcz7zER+SOecAF4VFWP+VqTMcYEg/mr99M8OpJbUoLvKW218TkYAFR1EbDonHG/83r9MPBwDfPOA+b5ow5jjAkW+SdLSd2Uxa0pCbRu2sTtci5I8HeoNcaYEPTmukOUVVQx8/IebpdywSwYjDHGz8orq3htzQGu7N2BPp1buV3OBbNgMMYYP/tkczbZx0tCqouqNwsGY4zxI1XlxVUZ9O7UkmsvCZ0uqt4sGIwxxo++Sj/KjuwiZl+VFFJdVL1ZMBhjjB/NXZVBp1YxTB0W73YpF82CwRhj/GRb1nG+3HOUO0f3ICYq0u1yLpoFgzHG+MncVRm0iI7kh5d1d7sUn1gwGGOMH2QWFPPx5mxuH5lIm2ahdUHbuSwYjDHGD/6xKgMBfnJlT7dL8ZkFgzHG+Ci3qIQ31h3i5uEJxLdt5nY5PrNgMMYYH724KoPKKuXfru3ldil+YcFgjDE+OHqylNfXHmDq0Hi6t2/hdjl+YcFgjDE++MeXGZRWVPGza3u7XYrfWDAYY8xFOnaqjH9+c4DvD46nV4g8z7k+LBiMMeYizftqH8Vlldw/tvHsLYCfHtQjIpOAZ/A8we0lVX3snOkPAbPwPMEtD/iJqh5wplUCW5ymB1V1ij9qMg1LVck/VcahY8VkFpymsLiM46fLOVFSQWWVEhEhCBDTJJLY5k2IbR5NuxbRdGvXnITYZjSJtO8kJrTlnyzlla/3MXlQF/qG4K21a+NzMIhIJPAcMB7IBNaJSKqqbvdqtgFIUdViEbkP+AtwmzPttKoO9bUO07COnixl9d58Nh8qZPPh42zPKuJkacV57WKiIoiMEFShSpXSiqrz2kQIdI1txiWdWzOoaxsGJ7RhUEIbOrSMCcQ/xRi/+PuKvZwur+Sh8X3dLsXv/LHHMBJIV9UMABF5E5gK/CsYVHW5V/s1wB1+WK9pQKrK9uwiPtmczcrdeWzLKgI8H/zJ8a25cVhXkjq2ILFdcxJim9O+ZTStmkadd3+Yyirl+OlyCovLOHqyjIPHijmQf4p9R0+xPbuIZTuPoOpp26dTS0b37sAVvdpzea/2tAqxxyGa8JFVeJp/rjnAzcMT6N2pce0tgH+CoStwyGs4E7islvZ3A596DTcVkTQ8h5keU9UPqptJRGYDswESExN9KtjULLeohLfSDvHBxizSc08SFSEM7x7LryZewpW9O5Ac3/qCDgNFRgjtWngOIyV1hJE92501/WRpBdsOH2fDoUJW783nzXUHmb96P00ihdG9OzBpQBfGJXe2vQkTVJ5ZugcUfjGuj9ulNAh/BEN1NxzXahuK3AGkAGO8RieqapaIJAFfiMgWVd173gJV5wJzAVJSUqpdvrl4mzMLeeXr/Xy8OYvySmVkj3b86caBTB4YR2yL6AZbb8uYKC5Las9lSe25d0wvSisq2XCwkC925vLp1mzmvLeFiPe3cHmv9kwfkcCkAXE0iw7du1aa0Lc37yTvfJfJj0Z1JyG2udvlNAh/BEMm0M1rOAHIOreRiIwDHgHGqGrpmfGqmuX8zhCRFcAw4LxgMA1jc2YhT3y+m1W78/51V8g7r+hBjw7uXKgTExXJqKT2jEpqz8PX9WN7dhGLt+bw/sbDPLhwE7+L2cYNQ+K4NaUbQ7u1RSQ0H4RiQteTS3YTExXR6HoiefNHMKwD+ohIT+AwMAP4gXcDERkGvAhMUtVcr/GxQLGqlopIB2A0nhPTpoHtO3qKxz/dyWfbcmjbvAlzruvHDy5LpHUQHdcXEQbEt2FAfBv+fVxf1u47xtvrD/HBhize+PYQQxLa8JMre3LdwDiio6yXk2l4mw4V8snmbH4+tnejPrwpqr4flRGRycDTeLqrzlPVP4nIo0CaqqaKyFJgEJDtzHJQVaeIyBV4AqMKzzUVT6vqy3WtLyUlRdPS0nyuOxyVlFfy9+XpvLAyg+ioCGZd1ZO7r+wZUid6T5ZW8P53mbyyej8Zeafo1CqGH43qzg9HdaddAx72MuFNVbnlhW/Yn3+K5b+8JqT+z5whIutVNaXOdv4IhkCzYLg4X+05ysPvb+bQsdNMGxrPbyb3p1Prpm6XddGqqpRVe/KY9/V+Vu3Oo3l0JHeM6s6sq3rSqVXo/rtMcPp4cxb3/98G/nzTIG4fGZodYCwYzL+UlFfy2Kc7mb96P0kdW/Df0wZyRa8ObpflV3uOnOC55emkbsqiSWQEt49M5N4xvejSxgLC+K6kvJJxT66kZUwUnzxwFZERoXluq77B4Jcrn03w2pFdxM/f2EB67knuvKIHv57Ur1H26unTuRVPzxjGL8b15e/L03ltzQH+b+1BfnBZIvc38uPBpuHN+3ofmQWneX3WZSEbChfC9hgasQ83HubX726mddMmPHHLEK7u29HtkgLm0LFinlueztvrM2kaFcGsq5L46dVJtIyx70LmwuSeKGHsEysZldSOl2Ze6nY5PqnvHoN15WiEyiur+MNH2/jFmxsZnNCWTx64KqxCAaBbu+Y8dvNgPn/waq7u25Fnlu1hzF+W88rX+yir5jYdxtTk8U93UVJeyW8m93e7lICxYGhkTpSU85P563jl6/38ZHRPXp91GR1bhe9hlF4dW/L8HSN4/9+uoE/nlvzho+1MeGolS7cfIRT3lk1grc3I593vMpl1VRJJjei22nWxYGhEjhSVcNuLa/hmbz5/mT6Y330/2e5i6hiWGMsbPx3FK3ddSmSEMOvVNGa+so703BNul2aCVHllFf/5wVa6tm3GA99rvBezVcc+NRqJ9NyT3PT31RzIP8XLd17KrSnd6p4pzIgI117Sic/+/Wp+e0MyGw4WMPHpL/nDR9s4frrc7fJMkHn5q33syT3J76cMoHl0eJ2bsmBoBHblnGDG3G8orahi4T2XMybMzidcqCaREdx9ZU9W/PIabk3pxvzV+7n2iRW88e1Bqqrs8JKBzIJinlm6h3H9OzM+ubPb5QScBUOI25FdxO3/WENkhPDWPaMY2LWN2yWFjPYtY/jzTYP4+OdX0rtjSx5+bwvTX1jNjuwit0szLlJVfp+6DYDfT0l2uRp3WDCEsG1Zx7n9H2uIiYpg4ezLw+rkmD8NiG/DwntG8f9uGcL+/GJu+N+v+J9FOyguO/9BRKbxS92UxdIduTw4vk+jvXtqXSwYQlRG3kl+/PK3NG8SyZuzR7l2N9TGQkS4eUQCyx4awy0jEpi7KoPxT65iyfYjbpdmAij3RAn/lbqNYYltufvKJLfLcY0FQwjKOV7Cj17+FoDXZl1G9/YWCv4S2yKax24ezDv3Xk7LmCh++moaP301jcOFp90uzTQwVeWR97dSXFbJX6cPCYsrnGtiwRBiCovL+NHLazl+upz5d420w0cNJKVHOz5+4ErmXNePL/fkMf7Jlcz7ah+VdnK60UrdlMWS7Uf45YS+9O4U3v+vLBhCSGlFJT99NY0D+cXM/fEIBiXYieaG1CQygnvH9GLJg2MY2bMdj368nekvrGb3Ebv2obE5UuQ5hDQ8zA8hnWHBECJUld+8t5V1+wt44tYhje7uqMGsW7vmvHLnpTx12xD2Hz3F9c9+ydNLd9utNRqJyirlwYUbKS2v4q+3hPchpDP8EgwiMklEdolIuojMqWZ6jIgsdKavFZEeXtMedsbvEpGJ/qinMXpxVQbvfpfJL77XhylD4t0uJ+yICDcOS2DpQ2O4bmAcTy/dww3/+yXfHSxwuzTjoxdX7WX13nx+PyWZXnZoFvBDMIhIJPAccB2QDNwuIud2/r0bKFDV3sBTwOPOvMl4HgU6AJgE/N1ZnvHy+bYcHv9sJzcMjuPfx/Vxu5yw1r5lDM/ePox5d6ZwoqSCm59fzR8+2sapUuvaGoo2Hirkyc93c/3gOLtbgBd/7DGMBNJVNUNVy4A3ganntJkKLHBevwN8TzxPcZ8KvKmqpaq6D0h3lmcc+46e4qG3NjGoaxueuGUIns1m3Da2X2c+f/Bq7risO698vZ8JT61i5e48t8syF6CopJwH3thA59ZN+Z8bB9n/LS/+CIauwCGv4UxnXLVtVLUCOA60r+e8YaukvJL7XltPVKTw/B0jaNrEdqaCSaumTfjjtIG8fe/lxDSJYOa8b3lo4UYKTpW5XZqpQ1WV8h9vbSKr8DTPzBhKm2ah9/zmhuSPYKguZs/t01dTm/rM61mAyGwRSRORtLy88Phm9rsPt7LryAmevm0oXds2c7scU4NLe7Rj0QNXcf+1vUndlMW4J1eSuinLbusdxJ5fuZcl24/wm8n9SenRzu1ygo4/giET8D44lwBk1dRGRKKANsCxes4LgKrOVdUUVU3p2LHx3yTurXWHeCstk59f25trLunkdjmmDk2bRPLLiZeQev+VdI1txgNvbODuBXZhXDBauTuPJz7fxdSh8dw1uofb5QQlfwTDOqCPiPQUkWg8J5NTz2mTCsx0Xk8HvlDP16lUYIbTa6kn0Af41g81hbT03BP89sOtjO7dnl+M6+t2OeYCJMe35r37ruA/r+/PN3vzmfDkSuZ/bRfGBYsD+af4xZsbuKRzK/58k51XqInPweCcM7gfWAzsAN5S1W0i8qiITHGavQy0F5F04CFgjjPvNuAtYDvwGfAzVa30taZQVlZRxS/e3EiLmCieum2o9akOQVGRnmdMf/7g1Yzo0Y7ff7Sdm59fzc4cu2urmwqLy7hr/joAXvzRiLB7xsKFkFA8DpqSkqJpaWlul9Eg/rp4J88t38uLPxrBxAFd3C7H+EhV+XBjFo9+vJ2i0+XcO6YX94/tbR0JAqy0opIfv/wtGw4W8tqsyxjZMzzPK4jIelVNqaudXfkcRNL2H+P5FXu5NSXBQqGREBGmDevK0ofGMGVIPH9bns7kZ75kbUa+26WFDVXl4Xe3sHbfMf56y+CwDYULYcEQJE6UlPPgWxtJiG3O774/wO1yjJ+1axHNk7cN5dWfjKSssorb5q7h4fe22CNFG5iq8thnO3lvw2EeGt+XqUOtN3x9WDAEicc/28nhgtM8eesQWsbYsc/G6uq+Hfn8wav56VU9WbjuIOOfXMmnW7Kta2sD+fuKvby4MoMfXpbIz8f2druckGHBEAS+3XeM19Yc5K7RPa1PdRhoHh3FI9cn88HPRtOhZQz3vf4dP5m/jgP5p9wurVGZ//U+/rp4FzcO68ofpw60HkgXwILBZSXllcx5bzMJsc34jwnWNTWcDE5oy4f3j+aRyf35dt8xxj+1iieX7KakPKw75vnF62sP8PuPtjMhuTN/nT6YCOvdd0EsGFz23PJ0MvJO8T83DrLuc2GoSWQEP706iS9+eQ2TBnTh2WV7GP/USpbtsEeKXqyXvszgkfe3MrZfJ/73B8OIirSPuQtlW8xFO3OKeH7FXm4a3pWr+zb+q7lNzTq3bsqztw/j/356GTFRkdy9II1ZC9Zx6Fix26WFDFXlmaV7+O9PdnD9oDheuGMEMVHWLfhiWDC4pKpKefi9LbRp1oTfXn/uXcpNuLqiVwcWPXAVv5ncj9V78xn35Er+8tlOTpRY76XaVFYpf/hoO08t3c1Nw7vyzIyhREfZx9vFsi3nkne/y2TDwUIentyf2BbRbpdjgkh0VASzr+7Fsv8Yw3UDu/D3FXu59okVvLbmABWV9tS4c50qrWD2q2nMX72fu6/syRPTh9jhIx/Z1nNBUUk5j3+2k+GJbblpmPWrNtWLa9OMp2cM48OfjSapQ0v+84OtXPfMlyzfmWvdWx1Zhae55YVvWL4rlz9OG8hvb0i2E81+YMHggqeX7CH/VBmPTh1ob2JTpyHd2rLwnlG8cMcIyiuruGv+Ou54eS0bwvyxoit25XL9s19y8FgxL995KT8a1d3tkhoN6wYTYLtyTrDgm/38YGQiA7u2cbscEyJEhEkDuzC2XydeW3OAvy1P58a/r2Zsv048NL5vWL2XKquUp5fu5m/L07mkcyue++Fwe1azn9lN9AJIVfnBP9ayI6eI5f9xjZ1bMBftZGkFC1bvZ+6qDI6fLmfigM48OL4v/bq0dru0BpWRd5JfvbOZ9QcKuGVEAo9OHUizaOt5VF/1vYme7TEE0JLtR/gmI58/ThtooWB80jImip9d25sfXd6dl7/cx7yv9rF425eMT+7MPVcnNbor6CurlFecK5ljoiJ46rYh3Dgswe2yGi3bYwiQ8soqJj69iggRPvvFVdZrwvhVYXEZr3y9n1e/2U9BcTkjusdyz9VJjOvfOeTPY60/UMAfPtrG5szjjOvfif+5cRCdWjd1u6yQZHsMQWbhukNk5J3ipR+nWCgYv2vbPJoHx/flnjFJvJ2WyUtfZTD7n+vp0b45P7gskekjutEuxPZSswpP85fPdvLBxiw6t47hmRlDmTIk3u55FAA+7TGISDtgIdAD2A/cqqoF57QZCjwPtAYqgT+p6kJn2nxgDHDcaX6nqm6sa72htsdwsrSCa/66nKSOLVk4e5S9sU2Dq6is4tOtObz6zX7W7S8gOiqC6wfFcfvIRC7tERvU78FDx4p5fuVe3k47hIgw+6ok7rumFy3srsM+C9Qewxxgmao+JiJznOFfn9OmGPixqu4RkXhgvYgsVtVCZ/qvVPUdH+sIanNXZXD0ZBkvzewf1P8hTeMRFRnB94fE8/0h8ezKOcH/rT3Ae98d5v0Nh+nathnfHxLP1KHx9OvSKijek6rKdwcLeX3NAVI3ZREhwm2XduPeMb1IiG3udnlhx9c9hl3ANaqaLSJxwApVvaSOeTYB052gmA98fKHBEEp7DLlFJYz56wrG9u/Ecz8Y7nY5JowVl1Xw2dYcUjdl8eWeo1RWKX06tWRccmeuvaQTwxPbBvww55GiEj7dks3CtEx2ZBfRMiaK6SMSuGdMEnFtmgW0lnBQ3z0GX4OhUFXbeg0XqGpsLe1HAguAAapa5QTD5UApsAyYo6qlda03lILhtx9s5c11B1n60Bi6t2/hdjnGAJB/spRFW3P4ZHMWafsLqKhSWjeN4qo+HRnZsx0jusfSr0srvwdFZZWyLes43+zNZ+mOI6QdKEAVkuNa88NRiUwb2tUOGTUgvwWDiCwFqnsA8SPAgvoGw5k9CmCmqq7xGpcDRANzgb2q+mgN888GZgMkJiaOOHDgQO3/siCQWVDMtU+s4NaUbvzpxkFul2NMtYpKyvl6z1GW78pl1e6j5BSVANA8OpKB8W3o07klvTu1pE+nViTENqNT65g6bxFfVaXknyoj53gJe/NOsiOniJ3ZJ/juYAEnSioA6NelFZMHxTF5UBy9O9kFaoEQqD2Geh1KEpHWeELhz6r6dg3Lugb4pareUNd6Q2WP4eH3NvPu+sOs+NU1xLe13WIT/FSVrOMlrD9QwHcHCthy+DjpuSfPezZ1y5goYls0ISYqkpioCJpERlBaUUVJeSXFZRUUnCqnzOuGf9GREfTu1JIh3dowKqk9lye1ty6nLgjUyedUYCbwmPP7w2oKiQbeB149NxREJM4JFQGmAVt9rCdoHMwv5u20TO4Y1d1CwYQMEaFr22Z0bduMKUPiAU9Y5J0sJT33JNmFJeSeKCX3RAmFxeWUVVRRWlFJaUUVHaIiaR4dSdMmEbRrEUNcm6Z0adOUHu1bkNSxBU2sm3bI8DUYHgPeEpG7gYPALQAikgLcq6qzgFuBq4H2InKnM9+Zbqmvi0hHQICNwL0+1hM0nv1iD5ERwn3X9HK7FGN8IiJ0atWUTq3sG3648CkYVDUf+F4149OAWc7r14DXaph/rC/rD1YZeSd577tM7hrdk862u2yMCTG2b9cAnl22h5ioSO4dY3sLxpjQY8HgZwfyT5G6KYs7RiXSsVWM2+UYY8wFs2DwsxdWZhAVEcGsq5LcLsUYYy6KBYMfHSkq4d31mUxPSbBzC8aYkGXB4Ecvf7WPiqoq7rna9haMMaHLgsFPjheX8/qaA9wwON5ufWGMCWkWDH6y4Jv9nCqrtOsWjDEhz4LBD4rLKnjl6318r18n+sc17mfuGmMaPwsGP1i47hAFxeX827W2t2CMCX0WDD7yPKR8PyNA0J9JAAAOjklEQVS6xzKie+N6ALsxJjxZMPho6Y4jHDxWzN1X9nS7FGOM8QsLBh+9/NU+urZtxoTkzm6XYowxfmHB4IOth4/z7b5j3HlFj4A/EtEYYxqKfZr5YN5X+2gRHcltI7u5XYoxxviNBcNFyi0q4aPNWdyS0o3WTZu4XY4xxviNBcNF+ueaA1RUKXeN7uF2KcYY41c+BYOItBORJSKyx/kdW0O7ShHZ6Pykeo3vKSJrnfkXOo8BDXol5ZW8vvYg4/p3tttfGGMaHV/3GOYAy1S1D7DMGa7OaVUd6vxM8Rr/OPCUM38BcLeP9QTEx5uzOXaqzPYWjDGNkq/BMBVY4LxeAEyr74wiIsBY4J2Lmd9Nr605QK+OLbg8qb3bpRhjjN/5GgydVTUbwPndqYZ2TUUkTUTWiMiZD//2QKGqVjjDmUDXmlYkIrOdZaTl5eX5WPbF23r4OBsPFXLHqO54ss0YYxqXqLoaiMhSoEs1kx65gPUkqmqWiCQBX4jIFqComnZa0wJUdS4wFyAlJaXGdg3t9bUHaNYkkpuGJ7hVgjHGNKg6g0FVx9U0TUSOiEicqmaLSByQW8MyspzfGSKyAhgGvAu0FZEoZ68hAci6iH9DwBSVlPPBhiymDImnTTPromqMaZx8PZSUCsx0Xs8EPjy3gYjEikiM87oDMBrYrqoKLAem1zZ/MHlvfSanyyu5Y1R3t0sxxpgG42swPAaMF5E9wHhnGBFJEZGXnDb9gTQR2YQnCB5T1e3OtF8DD4lIOp5zDi/7WE+DUVVeW3uQIQltGJTQxu1yjDGmwdR5KKk2qpoPfK+a8WnALOf1amBQDfNnACN9qSFQ1u47RnruSf4yfbDbpRhjTIOyK5/r6fW1B2ndNIrvD453uxRjjGlQFgz1UHCqjMVbc7hpeALNoiPdLscYYxqUBUM9fLjxMGWVVdyaYndRNcY0fhYMdVBVFqZlMqhrG5LjW7tdjjHGNDgLhjpsyypiR3YRt6bYBW3GmPBgwVCHt9IOER0VwZQhNd6twxhjGhULhlqUlFfywYbDTBrQhTbN7UpnY0x4sGCoxefbj1BUUsFtl9pJZ2NM+LBgqMXbaYfo2raZ3V7bGBNWLBhqkFlQzFfpR7klJYGICLu9tjEmfFgw1OC97w6jCjfb7bWNMWHGgqEaqsoHGw4zKqkd3do1d7scY4wJKAuGamzOPE7G0VPcOMy6qBpjwo8FQzXe33CY6KgIJg2Mc7sUY4wJOAuGc5RXVvHRpizG9e9kT2kzxoQlC4ZzfLXnKPmnyrhxmJ10NsaEJ5+CQUTaicgSEdnj/I6tps21IrLR66dERKY50+aLyD6vaUN9qccf3t9wmLbNmzCmb0e3SzHGGFf4uscwB1imqn2AZc7wWVR1uaoOVdWhwFigGPjcq8mvzkxX1Y0+1uOTk6UVfL49hxsGxxEdZTtTxpjw5Oun31RggfN6ATCtjvbTgU9VtdjH9TaIxVtzKCmvst5Ixpiw5mswdFbVbADnd6c62s8A3jhn3J9EZLOIPCUiMTXNKCKzRSRNRNLy8vJ8q7oGH2w8TLd2zRieeN4RMWOMCRt1BoOILBWRrdX8TL2QFYlIHDAIWOw1+mGgH3Ap0A74dU3zq+pcVU1R1ZSOHf1//D+3qISv049y49CuiNgtMIwx4SuqrgaqOq6maSJyRETiVDXb+eDPrWVRtwLvq2q517KznZelIvIK8Mt61u13n27NoUphytB4t0owxpig4OuhpFRgpvN6JvBhLW1v55zDSE6YIJ6v6NOArT7Wc9E+3pxFvy6t6N2plVslGGNMUPA1GB4DxovIHmC8M4yIpIjIS2caiUgPoBuw8pz5XxeRLcAWoAPw3z7Wc1Gyj59m3f4Crh9kVzobY0ydh5Jqo6r5wPeqGZ8GzPIa3g+c19VHVcf6sn5/WbQlB4DrB1swGGOMddbHcxgpOa41SR1bul2KMca4LuyDIbOgmA0HC7lhiO0tGGMMWDCwaIunY9QNg6w3kjHGgAUDn2zOZnBCGxLb2wN5jDEGwjwYDuYXsynzODfYSWdjjPmXsA6Gj7dkATDZuqkaY8y/hHUwLNqSzdBubUmItcNIxhhzRtgGQ2ZBMVsPF3HdwC5ul2KMMUElbIPh821HAJg4wILBGGO8hW0wfLYth35dWtGjQwu3SzHGmKASlsFw9GQpafuPMcH2Fowx5jxhGQxLtx+hSmHigM5ul2KMMUEnLINh8bYcurVrRnJca7dLMcaYoBN2wXCipJyv0/OZmNzFntRmjDHVCLtgWL4rj7LKKiZaN1VjjKmWT8EgIreIyDYRqRKRlFraTRKRXSKSLiJzvMb3FJG1IrJHRBaKSLQv9dTH4m05dGgZw/DE2IZelTHGhCRf9xi2AjcBq2pqICKRwHPAdUAycLuIJDuTHweeUtU+QAFwt4/11KqkvJIVO3MZn9yZyAg7jGSMMdXxKRhUdYeq7qqj2UggXVUzVLUMeBOY6jzneSzwjtNuAZ7nPjeYr9OPcqqs0nojGWNMLQJxjqErcMhrONMZ1x4oVNWKc8Y3mMXbcmgVE8UVvTo05GqMMSak1fnMZxFZClR3pvYRVf2wHuuo7piN1jK+pjpmA7MBEhMT67Ha8/Xs0JI7Lu9OdFTYnXM3xph6qzMYVHWcj+vIBLp5DScAWcBRoK2IRDl7DWfG11THXGAuQEpKSo0BUpv7rul1MbMZY0xYCcRX53VAH6cHUjQwA0hVVQWWA9OddjOB+uyBGGOMaUC+dle9UUQygcuBT0RksTM+XkQWATh7A/cDi4EdwFuqus1ZxK+Bh0QkHc85h5d9qccYY4zvxPPFPbSkpKRoWlqa22UYY0xIEZH1qlrjNWdn2FlYY4wxZ7FgMMYYcxYLBmOMMWexYDDGGHMWCwZjjDFnCcleSSKSBxy4yNk74Lm4LthYXRfG6rowVteFaax1dVfVjnU1Cslg8IWIpNWnu1agWV0Xxuq6MFbXhQn3uuxQkjHGmLNYMBhjjDlLOAbDXLcLqIHVdWGsrgtjdV2YsK4r7M4xGGOMqV047jEYY4ypRaMMBhG5RUS2iUiViNR4Bl9EJonILhFJF5E5XuN7ishaEdkjIgud24X7o652IrLEWe4SEYmtps21IrLR66dERKY50+aLyD6vaUMDVZfTrtJr3ale493cXkNF5Bvn771ZRG7zmubX7VXT+8Vreozz7093tkcPr2kPO+N3ichEX+q4iLoeEpHtzvZZJiLdvaZV+zcNUF13ikie1/pneU2b6fzd94jIzADX9ZRXTbtFpNBrWoNsLxGZJyK5IrK1hukiIs86NW8WkeFe0/y/rVS10f0A/YFLgBVASg1tIoG9QBIQDWwCkp1pbwEznNcvAPf5qa6/AHOc13OAx+to3w44BjR3hucD0xtge9WrLuBkDeNd215AX6CP8zoeyAba+nt71fZ+8Wrzb8ALzusZwELndbLTPgbo6SwnMoB1Xev1HrrvTF21/U0DVNedwN+qmbcdkOH8jnVexwaqrnPa/xyYF4DtdTUwHNhaw/TJwKd4nnw5CljbkNuqUe4xqOoOVd1VR7ORQLqqZqhqGfAmMFVEBBgLvOO0WwBM81NpU53l1Xe504FPVbXYT+uvyYXW9S9uby9V3a2qe5zXWUAuUOcFPBeh2vdLLfW+A3zP2T5TgTdVtVRV9wHpzvICUpeqLvd6D63B87TEhlaf7VWTicASVT2mqgXAEmCSS3XdDrzhp3XXSFVX4fkSWJOpwKvqsQbP0y/jaKBt1SiDoZ66Aoe8hjOdce2BQvU8YMh7vD90VtVsAOd3pzraz+D8N+WfnF3Jp0QkJsB1NRWRNBFZc+bwFkG0vURkJJ5vgXu9Rvtre9X0fqm2jbM9juPZPvWZtyHr8nY3nm+eZ1T3Nw1kXTc7f593ROTMI4CDYns5h9x6Al94jW6o7VWXmupukG1V5zOfg5WILAW6VDPpEVWtzyNCpZpxWst4n+uq7zKc5cQBg/A8+e6Mh4EcPB9+c/E8Ae/RANaVqKpZIpIEfCEiW4Ciatq5tb3+CcxU1Spn9EVvr+pWUc24c/+dDfKeqkO9ly0idwApwBiv0ef9TVV1b3XzN0BdHwFvqGqpiNyLZ29rbD3nbci6zpgBvKOqlV7jGmp71SWg762QDQZVHefjIjKBbl7DCUAWnvuQtBWRKOdb35nxPtclIkdEJE5Vs50PstxaFnUr8L6qlnstO9t5WSoirwC/DGRdzqEaVDVDRFYAw4B3cXl7iUhr4BPgP53d7DPLvujtVY2a3i/VtckUkSigDZ7DA/WZtyHrQkTG4QnbMapaemZ8DX9Tf3zQ1VmXquZ7Df4DeNxr3mvOmXeFH2qqV11eZgA/8x7RgNurLjXV3SDbKpwPJa0D+oinR000njdBqnrO6CzHc3wfYCZQnz2Q+kh1llef5Z53bNP5cDxzXH8aUG0PhoaoS0RizxyKEZEOwGhgu9vby/nbvY/n+Ovb50zz5/aq9v1SS73TgS+c7ZMKzBBPr6WeQB/gWx9quaC6RGQY8CIwRVVzvcZX+zcNYF1xXoNT8DwTHjx7yROc+mKBCZy959ygdTm1XYLnZO43XuMacnvVJRX4sdM7aRRw3Pni0zDbqiHOsLv9A9yIJ0lLgSPAYmd8PLDIq91kYDeexH/Ea3wSnv+46cDbQIyf6moPLAP2OL/bOeNTgJe82vUADgMR58z/BbAFzwfca0DLQNUFXOGse5Pz++5g2F7AHUA5sNHrZ2hDbK/q3i94Dk1NcV43df796c72SPKa9xFnvl3AdX5+v9dV11Ln/8GZ7ZNa1980QHX9GdjmrH850M9r3p842zEduCuQdTnDvwceO2e+BtteeL4EZjvv5Uw854LuBe51pgvwnFPzFrx6WzbEtrIrn40xxpwlnA8lGWOMqYYFgzHGmLNYMBhjjDmLBYMxxpizWDAYY4w5iwWDMcaYs1gwGGOMOYsFgzHGmLP8fyaH21UeGpZNAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "x = np.linspace(-1, 1, 1000)\n",
    "plt.plot(x,  2 * x ** 3 - 0 * x **2  - x)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {
    "scrolled": false,
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "      fun: -1.146240476075504e+25\n",
       " hess_inv: array([[1]])\n",
       "      jac: array([0.])\n",
       "  message: 'Optimization terminated successfully.'\n",
       "     nfev: 48\n",
       "      nit: 1\n",
       "     njev: 16\n",
       "   status: 0\n",
       "  success: True\n",
       "        x: array([-1.78956955e+08])"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from scipy.optimize import minimize\n",
    "minimize(lambda x:  2 * x ** 3 - 0 * x **2  - x, x0=-0.5)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "Powell hybrid method\n",
    "====\n",
    "\n",
    "**Type:** Root finding\n",
    "\n",
    "**Discrete/Continuous:** Continuous\n",
    "\n",
    "**Dimensions:** N\n",
    "\n",
    "**Derivative:** optional\n",
    "\n",
    "**Non-Convex:** not recommended\n",
    "\n",
    "**Python:** `root` unless `method` argument specifies a different method"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "1D Example\n",
    "---\n",
    "\n",
    "This is exactly like `newton` from above. Solve this equation:\n",
    "\n",
    "$$\n",
    "\\cos x + \\sin x = x\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "    fjac: array([[-1.]])\n",
      "     fun: array([2.22044605e-16])\n",
      " message: 'The solution converged.'\n",
      "    nfev: 8\n",
      "     qtf: array([-2.42716958e-12])\n",
      "       r: array([1.64467312])\n",
      "  status: 1\n",
      " success: True\n",
      "       x: array([1.25872818])\n"
     ]
    }
   ],
   "source": [
    "from scipy.optimize import root\n",
    "#rearranged equation so all terms on one side\n",
    "result = root(lambda x: np.cos(x) + np.sin(x) - x, x0=1)\n",
    "print(result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "The result type is like what we saw for minimize. Similar terms are here, including the root and the value of the function at the root. Notice it's not exactly $0$ at the root."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[2.22044605e-16]\n"
     ]
    }
   ],
   "source": [
    "x = result.x \n",
    "print(np.cos(x) + np.sin(x) - x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "Solve the following system of equations:\n",
    "$$ 3 x^2 - 2 y = 4$$\n",
    "$$ x - 4 y ^ 2 = -2$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "    fjac: array([[-0.99129609, -0.13165126],\n",
       "       [ 0.13165126, -0.99129609]])\n",
       "     fun: array([2.87769808e-13, 3.16324744e-12])\n",
       " message: 'The solution converged.'\n",
       "    nfev: 10\n",
       "     qtf: array([4.70588464e-09, 2.08104163e-08])\n",
       "       r: array([-8.45721872,  2.95677276,  7.06016501])\n",
       "  status: 1\n",
       " success: True\n",
       "       x: array([1.39555277, 0.92135129])"
      ]
     },
     "execution_count": 38,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def sys(v):\n",
    "    #I'm using v here to distinguish from the x in the equations\n",
    "    #extract x and y\n",
    "    x = v[0]\n",
    "    y = v[1]\n",
    "    #compute equations\n",
    "    eq1 = 3 * x ** 2 - 2 * y - 4\n",
    "    eq2 = x - 4 * y**2 + 2\n",
    "    #pack into list\n",
    "    sys_eq = [eq1, eq2]\n",
    "    return sys_eq\n",
    "root(sys, x0=[1,1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "So the answer is $x = 1.40,\\: y = 0.921$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "You should not ignore the information in the output of root. Imagine this small modification:\n",
    "\n",
    "$$ 3 x^2 - 2 y^2 = 4$$\n",
    "$$ x^3 - 4 y = -2$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "    fjac: array([[-0.79311064, -0.60907759],\n",
       "       [ 0.60907759, -0.79311064]])\n",
       "     fun: array([-0.31652211,  0.4128055 ])\n",
       " message: 'The iteration is not making good progress, as measured by the \\n  improvement from the last ten iterations.'\n",
       "    nfev: 30\n",
       "     qtf: array([-3.93523929e-04, -5.20186955e-01])\n",
       "       r: array([-1.16194697e+01,  6.56882367e+00, -1.16411656e-03])\n",
       "  status: 5\n",
       " success: False\n",
       "       x: array([1.53592084, 1.30262823])"
      ]
     },
     "execution_count": 39,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def sys2(v):\n",
    "    #I'm using v here to distinguish from the x in the equations\n",
    "    #extract x and y\n",
    "    x = v[0]\n",
    "    y = v[1]\n",
    "    #compute equations\n",
    "    eq1 = 3 * x ** 2 - 2 * y**2 - 4\n",
    "    eq2 = x**3 - 4 * y + 2\n",
    "    #pack into list\n",
    "    sys_eq = [eq1, eq2]\n",
    "    return sys_eq\n",
    "root(sys2, x0=[1,1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "If you did not read the message or status, you might believe the method succeeded. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "More compact version:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "    fjac: array([[-0.99129609, -0.13165126],\n",
       "       [ 0.13165126, -0.99129609]])\n",
       "     fun: array([2.87769808e-13, 3.16324744e-12])\n",
       " message: 'The solution converged.'\n",
       "    nfev: 10\n",
       "     qtf: array([4.70588464e-09, 2.08104163e-08])\n",
       "       r: array([-8.45721872,  2.95677276,  7.06016501])\n",
       "  status: 1\n",
       " success: True\n",
       "       x: array([1.39555277, 0.92135129])"
      ]
     },
     "execution_count": 40,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "root(lambda x: [3 * x[0] ** 2 - 2 * x[1] - 4, x[0] - 4 * x[1] ** 2 + 2], x0=[1,1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "Importance of starting position\n",
    "----\n",
    "\n",
    "By the way, there are two roots to this function!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {
    "scrolled": true,
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "    fjac: array([[ 0.98296038,  0.18381755],\n",
       "       [-0.18381755,  0.98296038]])\n",
       "     fun: array([ 2.66453526e-15, -4.44089210e-16])\n",
       " message: 'The solution converged.'\n",
       "    nfev: 16\n",
       "     qtf: array([ 2.87224221e-12, -7.25968027e-13])\n",
       "       r: array([ 5.41242052, -0.44135572,  6.9600742 ])\n",
       "  status: 1\n",
       " success: True\n",
       "       x: array([ 0.87635927, -0.84799164])"
      ]
     },
     "execution_count": 41,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "root(lambda x: [3 * x[0] ** 2 - 2 * x[1] - 4, x[0] - 4 * x[1] ** 2 + 2], x0=[0, 0])"
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Slideshow",
  "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.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}