extensions.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Extensions and compiled code "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "High-performance python typically means using extension modules (read: other languages)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Why extend?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Execution speed, obviously"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Vectorizing the algorithm is difficult, time-consuming, and inefficient"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Vectorization may be memory-limitiing (unless you use `memmap`, but that's another ball of wax)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Why should I not extend?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* It's not python, so you are typically dealing with the compile of code and linking to shared libraries. Yuck.  Installing and redistributing the code just got hard."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* It's not python.  Python is sooo much easier."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "What choices to I have?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Direct calls to compiled libraries: `ctypes`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Just-in-time compilation: `numba`, `PyPy`, and the now defunct `psyco`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Building python bindings: hand-binding, `SWIG`, `f2py`, and `cython`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using just-in-time compilation: `numba`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* `jit`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Overwriting numba_mandel.py\n"
     ]
    }
   ],
   "source": [
    "%%file numba_mandel.py\n",
    "\"\"\"\n",
    "Compute and plot the Mandelbrot set using matplotlib.\n",
    "\"\"\"\n",
    "\n",
    "import numpy as np\n",
    "import matplotlib.pylab as mpl\n",
    "\n",
    "from numba import jit\n",
    "\n",
    "@jit\n",
    "def mandel(x, y, max_iters):\n",
    "    \"\"\"\n",
    "    Given the real and imaginary parts of a complex number,\n",
    "    determine if it is a candidate for membership in the Mandelbrot\n",
    "    set given a fixed number of iterations.\n",
    "    \"\"\"\n",
    "    c = complex(x,y)\n",
    "    z = 0j\n",
    "    for i in range(max_iters):\n",
    "        z = z*z + c\n",
    "        if z.real * z.real + z.imag * z.imag >= 4:\n",
    "            return 255 * i // max_iters\n",
    "\n",
    "    return 255\n",
    "\n",
    "@jit(nopython=True)\n",
    "def create_fractal(min_x, max_x, min_y, max_y, image, iters):\n",
    "    height = image.shape[0]\n",
    "    width = image.shape[1]\n",
    "\n",
    "    pixel_size_x = (max_x - min_x) / width\n",
    "    pixel_size_y = (max_y - min_y) / height\n",
    "    for x in range(width):\n",
    "        real = min_x + x * pixel_size_x\n",
    "        for y in range(height):\n",
    "            imag = min_y + y * pixel_size_y\n",
    "            color = mandel(real, imag, iters)\n",
    "            image[y, x] = color\n",
    "\n",
    "    return image\n",
    "\n",
    "image = np.zeros((700, 1400), dtype=np.uint8)\n",
    "\n",
    "import time\n",
    "start = time.time()\n",
    "create_fractal(-2.0, 1.0, -1.0, 1.0, image, 20)\n",
    "print \"took %s\" % (time.time() - start)\n",
    "\n",
    "# mpl.imshow(image)\n",
    "# mpl.gray()\n",
    "# mpl.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "took 0.172127962112\r\n"
     ]
    }
   ],
   "source": [
    "!python2.7 numba_mandel.py"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Overwriting bubblesort.py\n"
     ]
    }
   ],
   "source": [
    "%%file bubblesort.py\n",
    "\"\"\"\n",
    "demonstrate compiled extensions using bubblesort\n",
    "\"\"\"\n",
    "\n",
    "import numpy as np\n",
    "from numba import jit\n",
    "\n",
    "def bubblesort(items):\n",
    "    length = len(items)\n",
    "    swapped = 1\n",
    "    for i in range(0, length):\n",
    "        if swapped:\n",
    "            swapped = 0\n",
    "            for ele in range(0, length-i-1):\n",
    "                if items[ele] > items[ele + 1]:\n",
    "                    temp = items[ele + 1]\n",
    "                    items[ele + 1] = items[ele]\n",
    "                    items[ele] = temp\n",
    "                    swapped = 1\n",
    "    return items\n",
    "\n",
    "jitbubblesort = jit(bubblesort)\n",
    "\n",
    "\n",
    "import time\n",
    "for N in (100,100,1000,1000,10000):\n",
    "\n",
    "    randoms = np.random.randint(0,1000,(N)).tolist()\n",
    "\n",
    "    print \"For N=%s\" % N\n",
    "    start = time.time()\n",
    "    x = bubblesort(randoms)\n",
    "    print \"%s: python\" % (time.time() - start)\n",
    "    assert np.all(sorted(randoms) == x)\n",
    "\n",
    "    start = time.time()\n",
    "    x = jitbubblesort(randoms)\n",
    "    print \"%s: numba\" % (time.time() - start)\n",
    "    assert np.all(sorted(randoms) == x)\n",
    "    \n",
    "    print ''\n",
    "\n",
    "\n",
    "# EOF"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "For N=100\n",
      "0.00119304656982: python\n",
      "0.191457986832: numba\n",
      "\n",
      "For N=100\n",
      "0.00129985809326: python\n",
      "8.01086425781e-05: numba\n",
      "\n",
      "For N=1000\n",
      "0.120245933533: python\n",
      "0.000295877456665: numba\n",
      "\n",
      "For N=1000\n",
      "0.120244026184: python\n",
      "0.000308990478516: numba\n",
      "\n",
      "For N=10000\n",
      "11.8938641548: python\n",
      "0.00219893455505: numba\n",
      "\n"
     ]
    }
   ],
   "source": [
    "!python2.7 bubblesort.py"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* `jit` with type annotations"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "See: http://numba.pydata.org/numba-doc/dev/user/examples.html"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* `autojit`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "For N=100\n",
      "7.39097595215e-05: python\n",
      "0.119683027267: numba\n",
      "\n",
      "For N=1000\n",
      "0.000465154647827: python\n",
      "1.4066696167e-05: numba\n",
      "\n",
      "For N=10000\n",
      "0.00473999977112: python\n",
      "4.91142272949e-05: numba\n",
      "\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/opt/local/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/IPython/kernel/__main__.py:14: RuntimeWarning: overflow encountered in double_scalars\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from numba import autojit\n",
    "import time\n",
    "\n",
    "def fib(N):\n",
    "    x = np.zeros(N, dtype=np.float64)\n",
    "    for i in range(N):\n",
    "        if i == 0:\n",
    "            x[i] = 0\n",
    "        elif i == 1:\n",
    "            x[i] = 1\n",
    "        else:\n",
    "            x[i] = x[i - 1] + x[i - 2]\n",
    "    return x\n",
    "\n",
    "jitfib = autojit(fib)\n",
    "\n",
    "\n",
    "for N in (100,1000,10000):\n",
    "\n",
    "    print \"For N=%s\" % N\n",
    "    start = time.time()\n",
    "    x = fib(N)\n",
    "    print \"%s: python\" % (time.time() - start)\n",
    "\n",
    "    start = time.time()\n",
    "    x = jitfib(N)\n",
    "    print \"%s: numba\" % (time.time() - start)\n",
    "    \n",
    "    print ''"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using python bindings: `cython`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Write python code, compile it as C -- from the commandline with `cython`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Writing pi.pyx\n"
     ]
    }
   ],
   "source": [
    "%%file pi.pyx\n",
    "\n",
    "def pied( int num ) :\n",
    "    return num * 3.14159265359\n",
    "\n",
    "def vec_pied( int r ) :\n",
    "    retList = []\n",
    "    cdef unsigned int i\n",
    "    for i in range(r):\n",
    "        retList.append( 3.14159 * i )\n",
    "    return retList\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Overwriting setup.py\n"
     ]
    }
   ],
   "source": [
    "%%file setup.py\n",
    "\n",
    "from distutils.core import setup\n",
    "from distutils.extension import Extension\n",
    "from Cython.Build import cythonize\n",
    "import numpy as np\n",
    "    \n",
    "extensions=[\n",
    "  Extension(\"pi\", [\"pi.pyx\"], include_dirs=[np.get_include()]),\n",
    "]   \n",
    "\n",
    "setup(\n",
    "  ext_modules=cythonize(extensions)\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Compiling pi.pyx because it changed.\n",
      "Cythonizing pi.pyx\n",
      "running build_ext\n",
      "building 'pi' extension\n",
      "creating build\n",
      "creating build/temp.macosx-10.8-x86_64-2.7\n",
      "/usr/bin/clang -fno-strict-aliasing -fno-common -dynamic -pipe -Os -fwrapv -DNDEBUG -g -fwrapv -O3 -Wall -Wstrict-prototypes -I/opt/local/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/numpy/core/include -I/opt/local/Library/Frameworks/Python.framework/Versions/2.7/include/python2.7 -c pi.c -o build/temp.macosx-10.8-x86_64-2.7/pi.o\n",
      "/usr/bin/clang -bundle -undefined dynamic_lookup -L/opt/local/lib -Wl,-headerpad_max_install_names -L/opt/local/lib/db48 build/temp.macosx-10.8-x86_64-2.7/pi.o -o /Users/mmckerns/notebooks/mystic_pathos_tutorials/parallel/pi.so\n"
     ]
    }
   ],
   "source": [
    "!python2.7 setup.py build_ext --inplace"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3.14159265359\n",
      "[0.0, 3.14159, 6.28318, 9.424769999999999]\n"
     ]
    }
   ],
   "source": [
    "import pi\n",
    "print pi.pied(1)\n",
    "print pi.vec_pied(4)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Write python code, compile it as C -- interactively in `IPython`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%load_ext Cython"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%%cython\n",
    "def cy_bubblesort(items):\n",
    "    length = len(items)\n",
    "    swapped = 1\n",
    "    for i in range(0, length):\n",
    "        if swapped:\n",
    "            swapped = 0\n",
    "            for ele in range(0, length-i-1):\n",
    "                if items[ele] > items[ele + 1]:\n",
    "                    temp = items[ele + 1]\n",
    "                    items[ele + 1] = items[ele]\n",
    "                    items[ele] = temp\n",
    "                    swapped = 1\n",
    "    return items"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "For N=100\n",
      "0.000560998916626: Cython\n",
      "\n",
      "For N=100\n",
      "0.000540971755981: Cython\n",
      "\n",
      "For N=1000\n",
      "0.0564270019531: Cython\n",
      "\n",
      "For N=1000\n",
      "0.0560281276703: Cython\n",
      "\n",
      "For N=10000\n",
      "5.55099987984: Cython\n",
      "\n"
     ]
    }
   ],
   "source": [
    "import time\n",
    "for N in (100,100,1000,1000,10000):\n",
    "\n",
    "    randoms = np.random.randint(0,1000,(N)).tolist()\n",
    "\n",
    "    print \"For N=%s\" % N\n",
    "    start = time.time()\n",
    "    x = cy_bubblesort(randoms)\n",
    "    print \"%s: Cython\" % (time.time() - start)\n",
    "    assert np.all(sorted(randoms) == x)\n",
    "\n",
    "    print ''"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Python code with type annotations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "100 loops, best of 3: 3.29 ms per loop\n"
     ]
    }
   ],
   "source": [
    "def py_fib(n):\n",
    "    a, b = 0, 1\n",
    "    for i in range(n):\n",
    "        b, a = a + b, b\n",
    "    return a\n",
    "\n",
    "%timeit py_fib(10000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%%cython\n",
    "def cy_fib(n):\n",
    "    a, b = 0, 1\n",
    "    for i in range(n):\n",
    "        b, a = a + b, b\n",
    "    return a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "100 loops, best of 3: 2.03 ms per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit cy_fib(10000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Overwriting fib.pyx\n"
     ]
    }
   ],
   "source": [
    "%%file fib.pyx\n",
    "\n",
    "def cy_fib(int n):\n",
    "    cdef int a = 0\n",
    "    cdef int b = 1\n",
    "    cdef int tmp\n",
    "    for i in range(n):\n",
    "        tmp = b\n",
    "        b = a + b\n",
    "        a = tmp\n",
    "    return a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Overwriting setup2.py\n"
     ]
    }
   ],
   "source": [
    "%%file setup2.py\n",
    "\n",
    "from distutils.core import setup\n",
    "from Cython.Build import cythonize\n",
    "\n",
    "setup(\n",
    "    name = \"cython fib\",\n",
    "    ext_modules = cythonize('fib.pyx'),  # accepts a glob pattern\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Compiling fib.pyx because it changed.\n",
      "Cythonizing fib.pyx\n",
      "running build_ext\n",
      "building 'fib' extension\n",
      "/usr/bin/clang -fno-strict-aliasing -fno-common -dynamic -pipe -Os -fwrapv -DNDEBUG -g -fwrapv -O3 -Wall -Wstrict-prototypes -I/opt/local/Library/Frameworks/Python.framework/Versions/2.7/include/python2.7 -c fib.c -o build/temp.macosx-10.8-x86_64-2.7/fib.o\n",
      "/usr/bin/clang -bundle -undefined dynamic_lookup -L/opt/local/lib -Wl,-headerpad_max_install_names -L/opt/local/lib/db48 build/temp.macosx-10.8-x86_64-2.7/fib.o -o /Users/mmckerns/notebooks/mystic_pathos_tutorials/parallel/fib.so\n"
     ]
    }
   ],
   "source": [
    "!python2.7 setup2.py build_ext --inplace"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "100000 loops, best of 3: 5.78 µs per loop\n"
     ]
    }
   ],
   "source": [
    "import fib\n",
    "%timeit fib.cy_fib(10000)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Aside: line profiling with `cython -a`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "!cython-2.7 -a fib.pyx"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "See: [result](fib.html)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Working with contiguous memory: `cython` with `numpy`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "%%cython\n",
    "import numpy as np\n",
    "from numpy cimport float_t\n",
    "\n",
    "def numcy_fib(int N):\n",
    "    cdef int i\n",
    "    cdef float_t[::1] x = np.zeros(N, dtype=float)\n",
    "    for i in range(N):\n",
    "        if i == 0:\n",
    "            x[i] = 0\n",
    "        elif i == 1:\n",
    "            x[i] = 1\n",
    "        else:\n",
    "            x[i] = x[i - 1] + x[i - 2]\n",
    "    return x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The slowest run took 4.08 times longer than the fastest. This could mean that an intermediate result is being cached \n",
      "10000 loops, best of 3: 31.4 µs per loop\n"
     ]
    }
   ],
   "source": [
    "%timeit numcy_fib(10000)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Not covered: working with native parallelism: `prange` and `nogil`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Not covered: leveraging external functions, structures, and classes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**EXERCISE:** See if `cython` or `numba` can help speed up our Monte Carlo example. Start with `roll.py`.  What are potential issues for using either within a parallel `map` that potentially distributes code and or objects to another process or computer?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using direct calls to native libraries: `ctypes`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Accessing shared memory (we've seen this already)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Direct access to shared libraries"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Writing C code and making direct calls to your C functions"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}