Lecture-6B-HPC.ipynb

{
 "metadata": {
  "name": "Lecture-6B-HPC"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "# Lecture 6B - Tools for high-performance computing applications\n",
      "\n",
      "J.R. Johansson (robert@riken.jp) http://dml.riken.jp/~rob/\n",
      "\n",
      "The latest version of this [IPython notebook](http://ipython.org/notebook.html) lecture is available at [http://github.com/jrjohansson/scientific-python-lectures](http://github.com/jrjohansson/scientific-python-lectures).\n",
      "\n",
      "The other notebooks in this lecture series are indexed at [http://jrjohansson.github.com](http://jrjohansson.github.com)."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%pylab inline"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 1
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## multiprocessing\n",
      "\n",
      "Python has a built-in process-based library for concurrent computing, called `multiprocessing`. "
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import multiprocessing\n",
      "import os\n",
      "import numpy"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 2
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def task(args):\n",
      "    print \"PID =\", os.getpid(), \", args =\", args\n",
      "    \n",
      "    return os.getpid(), args"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 3
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "task(\"test\")"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "PID = 22959 , args = test\n"
       ]
      },
      {
       "output_type": "pyout",
       "prompt_number": 4,
       "text": [
        "(22959, 'test')"
       ]
      }
     ],
     "prompt_number": 4
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "pool = multiprocessing.Pool(processes=4)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 5
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "result = pool.map(task, [1,2,3,4,5,6,7,8])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 6
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "result"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "pyout",
       "prompt_number": 7,
       "text": [
        "[(22978, 1),\n",
        " (22978, 2),\n",
        " (22978, 3),\n",
        " (22978, 4),\n",
        " (22978, 5),\n",
        " (22978, 6),\n",
        " (22978, 7),\n",
        " (22978, 8)]"
       ]
      }
     ],
     "prompt_number": 7
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The multiprocessing package is very useful for highly parallel tasks that do not need to communicate with each other, other than when sending the initial data to the pool of processes and when and collecting the results. "
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## IPython parallel\n",
      "\n",
      "IPython includes a very interesting and versatile parallel computing environment, which is very easy to use. It builds on the concept of ipython engines and controllers, that one can connect to and submit tasks to. To get started using this framework for parallel computing, one first have to start up an IPython cluster of engines. The easiest way to do this is to use the `ipcluster` command,\n",
      "\n",
      "    $ ipcluster start -n 4\n",
      "\n",
      "Or, alternatively, from the \"Clusters\" tab on the IPython notebook dashboard page. This will start 4 IPython engines on the current host, which is useful for multicore systems. It is also possible to setup IPython clusters that spans over many nodes in a computing cluster. For more information about possible use cases, see the official documentation [Using IPython for parallel computing](http://ipython.org/ipython-doc/dev/parallel/).\n",
      "\n",
      "To use the IPython cluster in our Python programs or notebooks, we start by creating an instance of `IPython.parallel.Client`:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from IPython.parallel import Client"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 8
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "cli = Client()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 9
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Using the 'ids' attribute we can retreive a list of ids for the IPython engines in the cluster:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "cli.ids"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "pyout",
       "prompt_number": 10,
       "text": [
        "[0, 1, 2, 3]"
       ]
      }
     ],
     "prompt_number": 10
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Each of these engines are ready to execute tasks. We can selectively run code on individual engines:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def getpid():\n",
      "    \"\"\" return the unique ID of the current process \"\"\"\n",
      "    import os\n",
      "    return os.getpid()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 11
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# first try it on the notebook process\n",
      "getpid()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "pyout",
       "prompt_number": 12,
       "text": [
        "22959"
       ]
      }
     ],
     "prompt_number": 12
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# run it on one of the engines\n",
      "cli[0].apply_sync(getpid)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "pyout",
       "prompt_number": 13,
       "text": [
        "22163"
       ]
      }
     ],
     "prompt_number": 13
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# run it on ALL of the engines at the same time\n",
      "cli[:].apply_sync(getpid)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "pyout",
       "prompt_number": 14,
       "text": [
        "[22163, 22164, 22169, 22174]"
       ]
      }
     ],
     "prompt_number": 14
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We can use this cluster of IPython engines to execute tasks in parallel. The easiest way to dispatch a function to different engines is to define the function with the decorator:\n",
      "\n",
      "    @view.parallel(block=True)\n",
      "\n",
      "Here, `view` is supposed to be the engine pool which we want to dispatch the function (task). Once our function is defined this way we can dispatch it to the engine using the `map` method in the resulting class (in Python, a decorator is a language construct which automatically wraps the function into another function or a class).\n",
      "\n",
      "To see how all this works, lets look at an example:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "dview = cli[:]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 15
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "@dview.parallel(block=True)\n",
      "def dummy_task(delay):\n",
      "    \"\"\" a dummy task that takes 'delay' seconds to finish \"\"\"\n",
      "    import os, time\n",
      "\n",
      "    t0 = time.time()\n",
      "    pid = os.getpid()\n",
      "    time.sleep(delay)\n",
      "    t1 = time.time()\n",
      "    \n",
      "    return [pid, t0, t1]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 16
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# generate random delay times for dummy tasks\n",
      "delay_times = numpy.random.rand(4)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 17
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now, to map the function `dummy_task` to the random delay time data, we use the `map` method in `dummy_task`:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "dummy_task.map(delay_times)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "pyout",
       "prompt_number": 18,
       "text": [
        "[[22163, 1361240989.803429, 1361240989.890851],\n",
        " [22164, 1361240989.805938, 1361240990.295277],\n",
        " [22169, 1361240989.808163, 1361240990.723319],\n",
        " [22174, 1361240989.811059, 1361240990.506656]]"
       ]
      }
     ],
     "prompt_number": 18
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Let's do the same thing again with many more tasks and visualize how these tasks are executed on different IPython engines:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def visualize_tasks(results):\n",
      "    res = numpy.array(results)\n",
      "    fig, ax = subplots(figsize=(10, res.shape[1]))\n",
      "    \n",
      "    yticks = []\n",
      "    yticklabels = []\n",
      "    tmin = min(res[:,1])\n",
      "    for n, pid in enumerate(numpy.unique(res[:,0])):\n",
      "        yticks.append(n)\n",
      "        yticklabels.append(\"%d\" % pid)\n",
      "        for m in numpy.where(res[:,0] == pid)[0]:\n",
      "            ax.add_patch(Rectangle((res[m,1] - tmin, n-0.25),\n",
      "                         res[m,2] - res[m,1], 0.5, color=\"green\", alpha=0.5))\n",
      "        \n",
      "    ax.set_ylim(-.5, n+.5)\n",
      "    ax.set_xlim(0, max(res[:,2]) - tmin + 0.)\n",
      "    ax.set_yticks(yticks)\n",
      "    ax.set_yticklabels(yticklabels)\n",
      "    ax.set_ylabel(\"PID\")\n",
      "    ax.set_xlabel(\"seconds\")"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 19
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "delay_times = numpy.random.rand(64)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 20
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "result = dummy_task.map(delay_times)\n",
      "visualize_tasks(result)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAmsAAADOCAYAAABhAkWQAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHZlJREFUeJzt3X1wFIX9x/HPBQLxGkRNiKcgPX4o5BLIAw9GmigXBepD\nUUxikfAgkpbOOIUJSlqxaGuHgEx4MAjUaUcSLI5gGzVaIEUohwJNqSRGGwMi5UoC0oIYTAh6hOT3\nB+MNkQRDuLvdXN6vmczkNne7n+9uBj6zl721NDc3NwsAAACmFGJ0AAAAALSNsgYAAGBilDUAAAAT\no6wBAACYWHejA/iLxWIxOgIAAEC7tXXNZ1CfWWtubu6yX7/+9a8Nz8DszM/8zM/8zM/87fu6lKAu\nawAAAJ0dZQ0AAMDEKGtByul0Gh3BMF15don5md9pdARDMb/T6AiGCtb5Lc3f9UZpJ2WxWL7zPWAA\nAAAzuFRv4cwaAACAiVHWAAAATIyyBgAAYGKUNQAAABOjrAEAAJhY0N5uSpKmvzk9YNvad2KfoiOj\nA7a9y+GrbGaesT32ndgnSZ16hu/S2Y/RlTLqGF+43816DMya60oE40ytudw5g3G/BONMlyOoy5r9\nGnvAtrXz8M6Abu9y+CqbmWdsj52Hd0oK7O9FoHX2Y3SljDrGF+53sx4Ds+a6EsE4U2sud85g3C/B\nONPl4G1QAAAAE6OsAQAAmBhlDQAAwMQoawAAACZGWQMAADAxyhoAAICJ+bysVVdXKzU1VbGxsRoy\nZIhWrFghScrJyZHD4VB8fLzS0tJ06tQpSdLJkyeVmpqqXr16adasWd711NXVKTEx0fvVp08fzZkz\np8W2ioqKFBISorKyMl+PAQAAYAo+L2uhoaFavny5KisrVVpaqlWrVqmqqkrjxo1TZWWlKioqNGjQ\nIC1atEiSFBYWpgULFmjJkiUt1tOrVy+Vl5d7v77//e8rPT3d+/O6ujrl5+frtttu8/UIAAAApuHz\nD8W12Wyy2WySpPDwcDkcDh09elRjx471PicpKUlFRUWSJKvVquTkZB04cKDNdX7yySf63//+p5SU\nFO+yp59+Wk8++aTy8vLU3Nzc6utchS7v9/YEu+wJ9iuYDAAAwDfcH7jl/sDdruf69Q4Gbrdb5eXl\nSkpKarF8zZo1mjRpUotlFoulzfWsX79eDz/8sPdxWVmZjhw5onvvvVd5eXltvtY53dnx8AAAAH7y\n7ZNIO9buaPO5fitr9fX1ysjIUH5+vsLDw73Lc3Nz1aNHD2VmZrZ7XRs2bNC6deskSU1NTXr88ce1\ndu1a78/bOrMGAADQ2fmlrJ09e1bp6emaMmWKJkyY4F1eWFioTZs2adu2be1eV0VFhRobG5WYmCjp\n/N+qVVZWyul0SpKOHTum+++/X2+//baGDRvm0zkAAACM5vOy1tzcrKysLMXExCg7O9u7vKSkRHl5\nedqxY4fCwsJafV1rXn311RZn4Xr37q3jx497H6empmrp0qUUNQAAEJR8XtZ27dqldevWKS4uzns2\nbOHChZo9e7Y8Ho/3QoNRo0Zp9erVkiS73a66ujp5PB4VFxdry5Ytio6OliT96U9/0ubNm30dEwAA\noFPweVlLSUlRU1PTRcsvdbWn2+1u82cHDx685Pa2b9/e7mwAAACdDXcwAAAAMDHKGgAAgIlR1gAA\nAEyMsgYAAGBifr2DgdHcte6AbSvSGhnQ7V0OX2Uz84ztEWmNlBTY34tA6+zH6EoZdYwv3O9mPQZm\nzXUlgnGm1lzunMG4X4JxpsthaQ7Sj/+3WCzc2QAAAHQKl+otvA0KAABgYpQ1AAAAE6OsAQAAmBhl\nDQAAwMQoawAAACZGWQMAADAxyhoAAICJUdYAAABMjLIGAABgYpQ1AAAAE6OsAQAAmBhlDQAAwMQo\nawAAACZGWQMAADAxyhoAAICJUdYAAABMjLIGAABgYpQ1AAAAE6OsAQAAmFh3owP40/Q3pxsdoV32\nndin6Mhoo2NcpL25fJnfH/vCrPv3Qv7O2Bn2gS90tjkDldfo/XLh9o3O4gu2cJsk6Vj9MYOTBFYw\nHLuO2ndin5x2p54b85wh2w/qsma/xm50hHbZeXinKbO2N5cv8/tjX5h1/17I3xk7wz7whc42Z6Dy\nGr1fLty+0Vl8wV3rltR5/o/xlWA4dh218/BOQ8s5b4MCAACYGGUNAADAxChrAAAAJkZZAwAAMDHK\nGgAAgIn5vKxVV1crNTVVsbGxGjJkiFasWCFJysnJkcPhUHx8vNLS0nTq1ClJ0smTJ5WamqpevXpp\n1qxZLdbl8Xg0c+ZMDR48WA6HQ6+//rok6T//+Y/uuusuxcfHKzU1VUeOHPH1GAAAAKbg87IWGhqq\n5cuXq7KyUqWlpVq1apWqqqo0btw4VVZWqqKiQoMGDdKiRYskSWFhYVqwYIGWLFly0bpyc3Nls9m0\nf/9+VVVVafTo0ZKkuXPnavr06aqoqNAzzzyjefPm+XoMAAAAU/B5WbPZbEpISJAkhYeHy+Fw6OjR\noxo7dqxCQs5vLikpSTU1NZIkq9Wq5ORk9ezZ86J1FRQUtChiERERkqSqqirdeeedkiSn06ni4mJf\njwEAAGAKfv1QXLfbrfLyciUlJbVYvmbNGk2aNKnFMovF0uJxbW2tJGn+/PlyuVwaOHCgVq5cqaio\nKMXHx6uoqEizZ8/WG2+8obq6On3xxRe69tprW6zDVejyfm9PsMueYPfdcAAAAB3kcrnkcrna9Vy/\nlbX6+nplZGQoPz9f4eHh3uW5ubnq0aOHMjMzL/n6xsZG1dTUKDk5WUuXLtXy5cs1d+5cvfzyy1qy\nZIl+/vOfq7CwUHfccYf69u2rbt26XbQO53Snr8cCAAC4Yk6nU06n0/v42WefbfO5filrZ8+eVXp6\nuqZMmaIJEyZ4lxcWFmrTpk3atm3bd64jIiJCVqtVaWlpkqSMjAy99NJLkqQbbrhBRUVFks6XwqKi\nIl199dV+mAQAAMBYPv+btebmZmVlZSkmJkbZ2dne5SUlJcrLy1NxcbHCwsJafd2FLBaLxo8fr+3b\nt0uStm3bptjYWEnS559/rqamJknSokWLlJWV5esxAAAATMHnZ9Z27dqldevWKS4uTomJiZKkhQsX\navbs2fJ4PBo7dqwkadSoUVq9erUkyW63q66uTh6PR8XFxdqyZYuio6O1ePFiTZ06VdnZ2YqKilJB\nQYGk8+/zzps3TxaLRaNHj9aqVat8PQYAAIAp+LyspaSkeM96XejAgQNtvsbtdre6vH///tqxY8dF\ny9PT05Went7hjAAAAJ0FdzAAAAAwMcoaAACAiVHWAAAATIyyBgAAYGKUNQAAABPz6+2mjOaudRsd\noV0irZGmzNreXL7M7499Ydb9eyF/Z+wM+8AXOtucgcpr9H65cPtGZ/EFW7hNUuf5P8ZXguHYdVSk\nNdJ73I1gaf72p9EGCYvFctEH7QIAAJjRpXoLb4MCAACYGGUNAADAxChrAAAAJnbJslZYWKhhw4bJ\narXKarVqxIgRWrt2baCyAQAAdHltXg26du1a5efna9myZUpMTFRzc7PKy8uVk5Mji8WiadOmBTIn\nAABAl9Tm1aBJSUlav369BgwY0GK52+3WxIkT9Y9//CMgATuKq0EBAEBn0aGrQevq6i4qapJkt9tV\nV1fnu3QAAABoU5tlLSwsrM0XXepnAAAA8J023wa96qqrdPPNN7f6ooMHD6qhocGvwa4Ub4MCAIDO\n4lK9pc0LDKqqqvwWCAAAAO3D7aYAAAAM1qEza+Hh4bJYLG2u8Msvv/RNOgAAALSJM2sAAAAG69CZ\ntTNnzujFF1/UwYMHNXToUGVlZal79zafDgAAAD9o88zaj3/8Y/Xo0UMpKSnavHmz7Ha78vPzA52v\nwzizBgAAOotL9ZY2y9rQoUP10UcfSZIaGxs1cuRIlZeX+y+lj1HWAABAZ9GhOxhc+JYnb38CAAAY\no80za926dZPVavU+PnPmjK666qrzL+oEV4NyZg0AAHQWHbrA4Ny5c34LBAAAgPYJ6vc3p7853egI\nfrPvxD5FR0YbHaNVF2bbd2KfJLV4bNbckvnzXal9J/ap5ssajfm/MUZHCbhLHdvLPe7B8nti9BxG\nb78jOlNmW7hNz415zugY8IGgLmv2a+xGR/CbnYd3mna+C7PtPLxTklo8Nmtuyfz5rtTOwzvVcLYh\nqGdsy6WO7eUe92D5PTF6DqO33xGdKbO71m10BPhImxcYAAAAwHiUNQAAABOjrAEAAJgYZQ0AAMDE\nKGsAAAAm5vOyVl1drdTUVMXGxmrIkCFasWKFJCknJ0cOh0Px8fFKS0vTqVOnJEknT55UamqqevXq\npVmzZrVYl8fj0cyZMzV48GA5HA69/vrrLX5eVFSkkJAQlZWV+XoMAAAAU/B5WQsNDdXy5ctVWVmp\n0tJSrVq1SlVVVRo3bpwqKytVUVGhQYMGadGiRZKksLAwLViwQEuWLLloXbm5ubLZbNq/f7+qqqo0\nevRo78/q6uqUn5+v2267zdcjAAAAmIbPy5rNZlNCQoIkKTw8XA6HQ0ePHtXYsWMVEnJ+c0lJSaqp\nqZEkWa1WJScnq2fPnhetq6CgQPPmzfM+joiI8H7/9NNP68knn1TPnj25rRQAAAhafv1QXLfbrfLy\nciUlJbVYvmbNGk2aNKnFMovF0uJxbW2tJGn+/PlyuVwaOHCgVq5cqaioKJWVlenIkSO69957lZeX\nd9Frv+EqdHm/tyfYZU+wX/lQAAAAV8jlcsnlcrXruX4ra/X19crIyFB+fr7Cw8O9y3Nzc9WjRw9l\nZmZe8vWNjY2qqalRcnKyli5dquXLl2vu3Llau3atHn/8ca1du9b73LbOrDmnO30yCwAAgC85nU45\nnU7v42effbbN5/rlatCzZ88qPT1dU6ZM0YQJE7zLCwsLtWnTJr3yyivfuY6IiAhZrValpaVJkjIy\nMlRWVqa6ujpVVlbK6XRqwIABKi0t1f33389FBgAAICj5vKw1NzcrKytLMTExys7O9i4vKSlRXl6e\niouLFRYW1urrLmSxWDR+/Hht375dkrRt2zbFxsbq6quv1vHjx3Xo0CEdOnRIt912m95++20NGzbM\n16MAAAAYzudvg+7atUvr1q1TXFycEhMTJUkLFy7U7Nmz5fF4NHbsWEnSqFGjtHr1akmS3W5XXV2d\nPB6PiouLtWXLFkVHR2vx4sWaOnWqsrOzFRUVpYKCAl/HBQAAMDWfl7WUlBQ1NTVdtPzAgQNtvsbt\ndre6vH///tqxY8clt/fNmTcAAIBgxB0MAAAATIyyBgAAYGKUNQAAABOjrAEAAJiYX+9gYDR3rdvo\nCH4TaY007XwXZou0RkpSi8dmzS2ZP9+VirRG6qvGr4J6xrZc6the7nEPlt8To+cwevsd0Zky28Jt\nRkeAj1iag/TGmhaLhXuGAgCATuFSvYW3QQEAAEyMsgYAAGBilDUAAAATo6wBAACYGGUNAADAxChr\nAAAAJkZZAwAAMDHKGgAAgIlR1gAAAEyMsgYAAGBilDUAAAATo6wBAACYGGUNAADAxChrAAAAJkZZ\nAwAAMDHKGgAAgIlR1gAAAEyMsgYAAGBilDUAAAAT6250AH+a/uZ07TuxT9GR0UZHaVNH8/lzLl+t\nO1j3vdHrNruuMLu/Z/z2+q90e2Y9Jv7K1dZ6zbofvmELt+m5Mc8ZHQMmFNRlzX6NXTsP75T9GrvR\nUdrU0Xz+nMtX6w7WfW/0us2uK8zu7xm/vf4r3Z5Zj4m/crW1XrPuh2+4a91GR4BJ8TYoAACAiVHW\nAAAATIyyBgAAYGKUNQAAABOjrAEAAJgYZQ0AAMDEfF7WqqurlZqaqtjYWA0ZMkQrVqyQJOXk5Mjh\ncCg+Pl5paWk6deqUJOnkyZNKTU1Vr169NGvWrBbr8ng8mjlzpgYPHiyHw6E33nhDkvTiiy8qLi5O\niYmJGjVqlCoqKnw9BgAAgCn4vKyFhoZq+fLlqqysVGlpqVatWqWqqiqNGzdOlZWVqqio0KBBg7Ro\n0SJJUlhYmBYsWKAlS5ZctK7c3FzZbDbt379fVVVVuuOOOyRJkydP1ocffqjy8nI99dRTeuKJJ3w9\nBgAAgCn4/ENxbTabbDabJCk8PFwOh0NHjx7V2LFjvc9JSkpSUVGRJMlqtSo5OVkHDhy4aF0FBQXa\nv3+/93FERIQkqVevXt5l9fX1ioyM9PUYAAAApuDXOxi43W6Vl5crKSmpxfI1a9Zo0qRJLZZZLJYW\nj2trayVJ8+fPl8vl0sCBA7Vy5UpFRUVJklavXq1ly5bp9OnT2r17d6vbdxW6dPK/J+Uqc8meYJc9\nwe6jyQAAADrO5XLJ5XK167l+u8Cgvr5eGRkZys/PV3h4uHd5bm6uevTooczMzEu+vrGxUTU1NUpO\nTtbevXs1atQozZ071/vzxx57TJ9++qmWLVumGTNmtLoO53SnrrvnOjmnOylqAADANJxOp37zm994\nvy7FL2Xt7NmzSk9P15QpUzRhwgTv8sLCQm3atEmvvPLKd64jIiJCVqtVaWlpkqSMjAyVlZVd9LyJ\nEye2uhwAACAY+LysNTc3KysrSzExMcrOzvYuLykpUV5enoqLixUWFtbq6y5ksVg0fvx4bd++XZK0\nbds2xcbGSlKLv2/buHGj4uLifD0GAACAKfj8b9Z27dqldevWeT9aQ5IWLlyo2bNny+PxeC80GDVq\nlFavXi1Jstvtqqurk8fjUXFxsbZs2aLo6GgtXrxYU6dOVXZ2tqKiolRQUCBJWrVqlbZu3arQ0FD1\n6dPHuxwAACDY+LyspaSkqKmp6aLlrV3t+Q23293q8v79+2vHjh0XLX/++ec7nA8AAKAz4Q4GAAAA\nJkZZAwAAMDHKGgAAgIlR1gAAAEzMr3cwMJq71q1Ia6TctW6jo7Spo/n8OZev1h2s+97odZtdV5jd\n3zN+e/1Xuj2zHhN/5WprvWbdD9+whduMjgCTsjR/+wPOgoTFYrnos9sAAADM6FK9hbdBAQAATIyy\nBgAAYGKUNQAAABOjrAEAAJgYZQ0AAMDEKGtByuVyGR3BMF15don5md9ldARDMb/L6AiGCtb5KWtB\nKlh/YdujK88uMT/zu4yOYCjmdxkdwVDBOj9lDQAAwMQoawAAACYW1HcwAAAA6CzaqmRBe2/QIO2g\nAACgi+FtUAAAABOjrAEAAJgYZQ0AAMDEgrKslZSUKDo6WrfccosWL15sdJyAmjFjhq6//noNHTrU\n6CiGqK6uVmpqqmJjYzVkyBCtWLHC6EgB89VXXykpKUkJCQmKiYnRvHnzjI5kiHPnzikxMVHjx483\nOkrA2e12xcXFKTExUbfeeqvRcQKutrZWGRkZcjgciomJUWlpqdGRAmb//v1KTEz0fvXu3btL/fsn\nSYsWLVJsbKyGDh2qzMxMff3110ZH8pmguxr03LlzGjx4sLZu3aq+fftq5MiRevXVV+VwOIyOFhDv\nvfeewsPDNW3aNH300UdGxwm4Y8eO6dixY0pISFB9fb2GDx+uN998s8sc/4aGBlmtVjU2NiolJUVL\nlixRSkqK0bECatmyZdq7d6/q6ur01ltvGR0noAYMGKC9e/fquuuuMzqKIR555BGNHj1aM2bMUGNj\no06fPq3evXsbHSvgmpqa1LdvX+3Zs0c33XST0XECwu12684771RVVZV69uypiRMn6t5779Ujjzxi\ndDSfCLoza3v27NHNN98su92u0NBQPfzwwyouLjY6VsDcfvvtuvbaa42OYRibzaaEhARJUnh4uBwO\nh44ePWpwqsCxWq2SJI/Ho3PnznW5/7Rramq0adMm/eQnP+myV4R31blPnTql9957TzNmzJAkde/e\nvUsWNUnaunWrBg4c2GWKmiRdffXVCg0NVUNDgxobG9XQ0KC+ffsaHctngq6sHTlypMUvaL9+/XTk\nyBEDE8Eobrdb5eXlSkpKMjpKwDQ1NSkhIUHXX3+9UlNTFRMTY3SkgJozZ47y8vIUEhJ0/7S1i8Vi\n0ZgxYzRixAj94Q9/MDpOQB06dEh9+vTRo48+qmHDhumnP/2pGhoajI5liPXr1yszM9PoGAF13XXX\n6YknnlD//v1144036pprrtGYMWOMjuUzQfcvGh+GC0mqr69XRkaG8vPzFR4ebnScgAkJCdEHH3yg\nmpoavfvuu0F7n7zW/OUvf1FUVJQSExO77NmlXbt2qby8XJs3b9aqVav03nvvGR0pYBobG1VWVqbH\nHntMZWVl+t73vqfnnnvO6FgB5/F49Pbbb+uhhx4yOkpAHTx4UM8//7zcbreOHj2q+vp6vfLKK0bH\n8pmgK2t9+/ZVdXW193F1dbX69etnYCIE2tmzZ5Wenq4pU6ZowoQJRscxRO/evXXffffp/fffNzpK\nwOzevVtvvfWWBgwYoEmTJulvf/ubpk2bZnSsgLrhhhskSX369NGDDz6oPXv2GJwocPr166d+/fpp\n5MiRkqSMjAyVlZUZnCrwNm/erOHDh6tPnz5GRwmo999/Xz/4wQ8UERGh7t27Ky0tTbt37zY6ls8E\nXVkbMWKEDhw4ILfbLY/How0bNuj+++83OhYCpLm5WVlZWYqJiVF2drbRcQLqxIkTqq2tlSSdOXNG\n77zzjhITEw1OFTgLFy5UdXW1Dh06pPXr1+vOO+/Uyy+/bHSsgGloaFBdXZ0k6fTp09qyZUuXuirc\nZrPppptu0ieffCLp/N9txcbGGpwq8F599VVNmjTJ6BgBFx0drdLSUp05c0bNzc3aunVrUP0ZSNDd\nbqp79+5auXKlfvjDH+rcuXPKysrqMlcCStKkSZO0Y8cOff7557rpppv029/+Vo8++qjRsQJm165d\nWrdunffjC6Tzl3PffffdBifzv88++0yPPPKImpqa1NTUpKlTp+quu+4yOpZhutqfRPz3v//Vgw8+\nKOn8W4KTJ0/WuHHjDE4VWC+88IImT54sj8ejgQMHqqCgwOhIAXX69Glt3bq1y/29oiTFx8dr2rRp\nGjFihEJCQjRs2DDNnDnT6Fg+E3Qf3QEAABBMgu5tUAAAgGBCWQMAADAxyhoAAICJUdYAAABMjLIG\nAH7kcrm65E3lAfgOZQ0AAMDEKGsAgtrp06d13333KSEhQUOHDtVrr72mvXv3yul0asSIEbr77rt1\n7NgxSdKnn36qMWPGKCEhQcOHD9ehQ4ckSTk5ORo6dKji4uL02muvSTp/xszpdOqhhx6Sw+HQlClT\nvNssKSmRw+HQ8OHD9cYbb3iX79ixQ4mJiUpMTNSwYcNUX18fwD0BoLMKug/FBYALlZSUqG/fvtq4\ncaMk6csvv9Q999yjt956SxEREdqwYYN+9atf6aWXXtLkyZP11FNP6YEHHpDH49G5c+dUVFSkiooK\nffjhhzp+/LhGjhypO+64Q5L0wQcf6OOPP9YNN9yg5ORk7d692/thnNu3b9fAgQM1ceJE7wf0Ll26\nVKtXr9aoUaPU0NCgnj17GrZfAHQenFkDENTi4uL0zjvv6Mknn9TOnTt1+PBh/etf/9KYMWOUmJio\n3NxcHTlyRPX19Tp69KgeeOABSVKPHj101VVXadeuXcrMzJTFYlFUVJRGjx6tf/7zn7JYLLr11lt1\n4403ymKxKCEhQYcOHdK+ffs0YMAADRw4UJI0ZcoU743lk5OTNWfOHL3wwgv64osv1K1bN8P2C4DO\ngzNrAILaLbfcovLycm3cuFHz589XamqqYmNjL7rJ8zf31WzNt2/08s2ZsgvPjHXr1k2NjY0X3ebq\nwtf+8pe/1I9+9CNt3LhRycnJ+utf/6rBgwd3eDYAXQNn1gAEtc8++0xhYWGaPHmy5s6dqz179ujE\niRMqLS2VJJ09e1Yff/yxevXqpX79+qm4uFiS9PXXX+vMmTO6/fbbtWHDBjU1Nen48eN69913deut\nt15U4KTzJS46Olput1v//ve/JZ2/sfY3Dh48qNjYWP3iF7/QyJEjtX///gDsAQCdHWfWAAS1jz76\nSDk5OQoJCVGPHj30u9/9Tt26ddPs2bN16tQpNTY2as6cOYqJidEf//hH/exnP9Mzzzyj0NBQ/fnP\nf9aDDz6ov//974qPj5fFYlFeXp6ioqJUVVXV6s3ie/bsqd///ve67777ZLVadfvtt+v06dOSpPz8\nfG3fvl0hISEaMmSI7rnnnkDvDgCdEDdyBwAAMDHeBgUAADAxyhoAAICJUdYAAABMjLIGAABgYpQ1\nAAAAE6OsAQAAmNj/AzomRp+KrP2lAAAAAElFTkSuQmCC\n"
      }
     ],
     "prompt_number": 21
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "That's a nice and easy parallelization! We can see that we utilize all four engines quite well.\n",
      "\n",
      "But one short coming so far is that the tasks are not load balanced, so one engine might be idle while others still have more tasks to work on.\n",
      "\n",
      "However, the IPython parallel environment provides a number of alternative \"views\" of the engine cluster, and there is a view that provides load balancing as well (above we have used the \"direct view\", which is why we called it \"dview\").\n",
      "\n",
      "To obtain a load balanced view we simply use the `load_balanced_view` method in the engine cluster client instance `cli`:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "lbview = cli.load_balanced_view()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 22
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "@lbview.parallel(block=True)\n",
      "def dummy_task_load_balanced(delay):\n",
      "    \"\"\" a dummy task that takes 'delay' seconds to finish \"\"\"\n",
      "    import os, time\n",
      "\n",
      "    t0 = time.time()\n",
      "    pid = os.getpid()\n",
      "    time.sleep(delay)\n",
      "    t1 = time.time()\n",
      "    \n",
      "    return [pid, t0, t1]"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 23
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "result = dummy_task_load_balanced.map(delay_times)\n",
      "visualize_tasks(result)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAmsAAADOCAYAAABhAkWQAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHdJJREFUeJzt3XtwVPX9//HXBhLidiNqIK6CdBkUsgnkwsVIE2WjQL0U\nxSQWuYqkpTNOYYJCKxZt7RCQCReDQJ12hGBxBNuo0XIpQgkKNKWSGDUGRMrWBEwLYjAB6hKS7x/8\nsj9CEsB4lnN2eT5mMpP97J7PeX8+OZt9zTl7zrE1NTU1CQAAAJYUZnYBAAAAaB9hDQAAwMIIawAA\nABZGWAMAALCwzmYXECg2m83sEgAAAC5Ze+d8hvSetaamJn7O+fn1r39teg1W+2FOmBPmhDlhTpgT\nK8zLhYR0WAMAAAh2hDUAAAALI6xdQTwej9klWA5z0hpz0hpz0hpz0hpz0hpz0rZvOy+2posdKA1S\nNpvtoseAAQAArOBCuYU9awAAABZGWAMAALAwwhoAAICFEdYAAAAsjLAGAABgYSF7uylJmvzm5ID2\n73Q4JUk19TUBXU8w2Xt0r2K7xX7r56wo2Oo1QiiM+du+L4NxzFaqubkWs2qy0lxczN6je+VxeSSZ\n97kRDPMVDDUawelw6rnhz13Sa0M6rLmucQW0f2+t97KsJ5js+HxHu/NxoeesKNjqNUIojPnbvi+D\nccxWqrm5FrNqstJcXMyOz3f4Q5pZNQfDfAVDjUZo/l91KTgMCgAAYGGENQAAAAsjrAEAAFgYYQ0A\nAMDCCGsAAAAWRlgDAACwMMPDWlVVldLT0xUfH6/+/ftr6dKlkqRZs2bJ7XYrMTFRGRkZOn78uCTp\n2LFjSk9PV1RUlKZNm+bvp66uTsnJyf6f7t27a8aMGS3WVVhYqLCwMJWWlho9DAAAAEswPKyFh4dr\nyZIlqqioUElJiZYvX67KykqNHDlSFRUVKi8vV9++fTV//nxJUmRkpObOnauFCxe26CcqKkplZWX+\nn+9///vKzMz0P19XV6f8/HzddtttRg8BAADAMgy/KK7T6ZTTefYK4g6HQ263W4cPH9aIESP8r0lJ\nSVFhYaEkyW63KzU1Vfv372+3z08//VT//e9/lZaW5m97+umn9eSTTyovL09NTU1tLldcUOz/3ZXk\nkivJ9R1GBgAAYIyaj2v0mw9+c0mvDegdDLxer8rKypSSktKifeXKlRo7dmyLNpvN1m4/a9eu1cMP\nP+x/XFpaqkOHDunee+9VXl5eu8t6Jns6XjwAAECAOPs79ZvRv/E/fvbZZ9t9bcDCWn19vbKyspSf\nny+Hw+Fvz83NVUREhMaNG3fJfa1bt05r1qyRJDU2Nurxxx/X6tWr/c+3t2cNAAAg2AUkrJ0+fVqZ\nmZmaMGGCRo8e7W8vKCjQhg0btHXr1kvuq7y8XA0NDUpOTpZ09rtqFRUV8ng8kqSamhrdf//9evvt\ntzVw4EBDxwEAAGA2w8NaU1OTsrOzFRcXp5ycHH/7pk2blJeXp+3btysyMrLN5dry6quvttgL17Vr\nVx05csT/OD09XYsWLSKoAQCAkGR4WNu5c6fWrFmjhIQE/96wefPmafr06fL5fP4TDYYOHaoVK1ZI\nklwul+rq6uTz+VRUVKTNmzcrNjZWkvSnP/1JGzduNLpMAACAoGB4WEtLS1NjY2Or9gud7en1ett9\n7sCBAxdc37Zt2y65NgAAgGDDHQwAAAAsjLAGAABgYYQ1AAAACyOsAQAAWFhA72BgNm+tN6D9Ox3O\ny7KeYNLN3q3d+bjQc1YUbPUaIRTG/G3fl8E4ZivV3FyLWTVZaS4uppu9m+mfG8EwX8FQoxGat4VL\nYWsK0cv/22w27mwAAACCwoVyC4dBAQAALIywBgAAYGGENQAAAAsjrAEAAFgYYQ0AAMDCCGsAAAAW\nRlgDAACwMMIaAACAhRHWAAAALIywBgAAYGGENQAAAAsjrAEAAFgYYQ0AAMDCCGsAAAAWRlgDAACw\nMMIaAACAhRHWAAAALIywBgAAYGGENQAAAAvrbHYBgTT5zckdWs7pcEqSauprDKwmcPYe3avYbrHt\nPjaiT6MFun8j7T26V53DOiutV1rQbBOBFkx/v/MFsvZgnhcp+Os30t6jeyXpipyPK3U7cDqcem74\nc2aX0aaQDmuua1wdWs5b6/1Oy19uOz7f0aLW8x8b0afRAt2/kXZ8vkOOCIdq6muCpuZAC6a/3/kC\nWXswz4sU/PUbacfnOyQFz+eAka7U7aD5s9+KOAwKAABgYYQ1AAAACyOsAQAAWBhhDQAAwMIIawAA\nABZmeFirqqpSenq64uPj1b9/fy1dulSSNGvWLLndbiUmJiojI0PHjx+XJB07dkzp6emKiorStGnT\nWvTl8/k0depU9evXT263W6+//rok6d///rfuuusuJSYmKj09XYcOHTJ6GAAAAJZgeFgLDw/XkiVL\nVFFRoZKSEi1fvlyVlZUaOXKkKioqVF5err59+2r+/PmSpMjISM2dO1cLFy5s1Vdubq6cTqf27dun\nyspKDRs2TJI0c+ZMTZ48WeXl5XrmmWc0e/Zso4cBAABgCYaHNafTqaSkJEmSw+GQ2+3W4cOHNWLE\nCIWFnV1dSkqKqqurJUl2u12pqanq0qVLq75WrVrVIohFR0dLkiorK3XnnXdKkjwej4qKioweBgAA\ngCUE9KK4Xq9XZWVlSklJadG+cuVKjR07tkWbzWZr8bi2tlaSNGfOHBUXF6tPnz5atmyZYmJilJiY\nqMLCQk2fPl1vvPGG6urq9NVXX+naa69t0UdxQbH/d1eSS64kl3GDAwAA6KDi4mIVFxdf0msDFtbq\n6+uVlZWl/Px8ORwOf3tubq4iIiI0bty4Cy7f0NCg6upqpaamatGiRVqyZIlmzpypl19+WQsXLtTP\nf/5zFRQU6I477lCPHj3UqVOnVn14JnuMHhYAAMB35vF45PF4/I+fffbZdl8bkLB2+vRpZWZmasKE\nCRo9erS/vaCgQBs2bNDWrVsv2kd0dLTsdrsyMjIkSVlZWXrppZckSTfccIMKCwslnQ2FhYWFuvrq\nqwMwEgAAAHMZ/p21pqYmZWdnKy4uTjk5Of72TZs2KS8vT0VFRYqMjGxzuXPZbDaNGjVK27ZtkyRt\n3bpV8fHxkqQvv/xSjY2NkqT58+crOzvb6GEAAABYguF71nbu3Kk1a9YoISFBycnJkqR58+Zp+vTp\n8vl8GjFihCRp6NChWrFihSTJ5XKprq5OPp9PRUVF2rx5s2JjY7VgwQJNnDhROTk5iomJ0apVqySd\nPc47e/Zs2Ww2DRs2TMuXLzd6GAAAAJZgeFhLS0vz7/U61/79+9tdxuv1ttneq1cvbd++vVV7Zmam\nMjMzO1wjAABAsOAOBgAAABZGWAMAALAwwhoAAICFEdYAAAAsjLAGAABgYQG93ZTZvLXeDi3ndDi/\n0/KXWzd7txa1nv/YiD6NFuj+jdTN3k2dwzrL6XAGTc2BFkx/v/MFsvZgnhcp+Os3Ujd7N0nB8zlg\npCt1O2j+7LciW9P5V6MNETabrdWFdgEAAKzoQrmFw6AAAAAWRlgDAACwMMIaAACAhV0wrBUUFGjg\nwIGy2+2y2+0aPHiwVq9efblqAwAAuOK1ezbo6tWrlZ+fr8WLFys5OVlNTU0qKyvTrFmzZLPZNGnS\npMtZJwAAwBWp3bNBU1JStHbtWvXu3btFu9fr1ZgxY/SPf/zjshTYUZwNCgAAgkWHzgatq6trFdQk\nyeVyqa6uzrjqAAAA0K52w1pkZGS7C13oOQAAABin3cOgV111lW6++eY2Fzpw4IBOnjwZ0MK+Kw6D\nAgCAYHGh3NLuCQaVlZUBKwgAAACXhttNAQAAmKxDe9YcDodsNlu7HX799dfGVAcAAIB2sWcNAADA\nZB3as3bq1Cm9+OKLOnDggAYMGKDs7Gx17tzuywEAABAA7e5Z+/GPf6yIiAilpaVp48aNcrlcys/P\nv9z1dRh71gAAQLC4UG5pN6wNGDBAH330kSSpoaFBQ4YMUVlZWeCqNBhhDQAABIsO3cHg3EOeHP4E\nAAAwR7t71jp16iS73e5/fOrUKV111VVnFwqCs0HZswYAAIJFh04wOHPmTMAKAgAAwKUJ6eObk9+c\nbHYJl2Tv0b2K7RZrdhmXFWNuyelwSpJq6msuZ0kBFyzjulzboxnb/bnrDMX3ndHbWKjN0fnjCbXx\ntcXpcOq54c+ZXYahQjqsua5xmV3CJdnx+Y6gqdUojLklb61XUvBss5cqWMZ1ubZHM7b7c9cZiu87\no7exUJuj88cTauNrS/M2EUraPcEAAAAA5iOsAQAAWBhhDQAAwMIIawAAABZGWAMAALAww8NaVVWV\n0tPTFR8fr/79+2vp0qWSpFmzZsntdisxMVEZGRk6fvy4JOnYsWNKT09XVFSUpk2b1qIvn8+nqVOn\nql+/fnK73Xr99ddbPF9YWKiwsDCVlpYaPQwAAABLMDyshYeHa8mSJaqoqFBJSYmWL1+uyspKjRw5\nUhUVFSovL1ffvn01f/58SVJkZKTmzp2rhQsXtuorNzdXTqdT+/btU2VlpYYNG+Z/rq6uTvn5+brt\nttuMHgIAAIBlGB7WnE6nkpKSJEkOh0Nut1uHDx/WiBEjFBZ2dnUpKSmqrq6WJNntdqWmpqpLly6t\n+lq1apVmz57tfxwdHe3//emnn9aTTz6pLl26cFspAAAQsgJ6UVyv16uysjKlpKS0aF+5cqXGjh3b\nos1ms7V4XFtbK0maM2eOiouL1adPHy1btkwxMTEqLS3VoUOHdO+99yovL6/Vss2KC4r9v7uSXHIl\nub77oAAAAL6j4uJiFRcXX9JrAxbW6uvrlZWVpfz8fDkcDn97bm6uIiIiNG7cuAsu39DQoOrqaqWm\npmrRokVasmSJZs6cqdWrV+vxxx/X6tWr/a9tb8+aZ7LHkLEAAAAYyePxyOPx+B8/++yz7b42IGeD\nnj59WpmZmZowYYJGjx7tby8oKNCGDRv0yiuvXLSP6Oho2e12ZWRkSJKysrJUWlqquro6VVRUyOPx\nqHfv3iopKdH999/PSQYAACAkGR7WmpqalJ2drbi4OOXk5PjbN23apLy8PBUVFSkyMrLN5c5ls9k0\natQobdu2TZK0detWxcfH6+qrr9aRI0d08OBBHTx4ULfddpvefvttDRw40OihAAAAmM7ww6A7d+7U\nmjVrlJCQoOTkZEnSvHnzNH36dPl8Po0YMUKSNHToUK1YsUKS5HK5VFdXJ5/Pp6KiIm3evFmxsbFa\nsGCBJk6cqJycHMXExGjVqlVGlwsAAGBphoe1tLQ0NTY2tmrfv39/u8t4vd4223v16qXt27dfcH3N\ne94AAABCEXcwAAAAsDDCGgAAgIUR1gAAACyMsAYAAGBhAb2Dgdm8tV6zS7gk3ezdgqZWozDmlpwO\np6Tg2WYvVbCM63Jtj2Zs9+euMxTfd0ZvY6E2R+ePJ9TG15bmbSKU2JpC9MaaNpuNe4YCAICgcKHc\nwmFQAAAACyOsAQAAWBhhDQAAwMIIawAAABZGWAMAALAwwhoAAICFEdYAAAAsjLAGAABgYYQ1AAAA\nCyOsAQAAWBhhDQAAwMIIawAAABZGWAMAALAwwhoAAICFEdYAAAAsjLAGAABgYYQ1AAAACyOsAQAA\nWBhhDQAAwMI6m11AIE1+c7L2Ht2r2G6xhvVpdH9tcTqckqSa+pqArqfZ5RiTEZrrbJ6fZjX1NZd9\nzr6rjsy5lf5OVqrl2zK69mDb9ox0JY+9LcxHax15vzGPrYV0WHNd49KOz3fIdY3LsD6N7q8t3lqv\nJAV8Pc0ux5iM0Fxn8/w0O7ctGMYhdWzOrfR3slIt35bRtQfbtmekK3nsbWE+WuvI+415bI3DoAAA\nABZGWAMAALAwwhoAAICFEdYAAAAsjLAGAABgYYQ1AAAACzM8rFVVVSk9PV3x8fHq37+/li5dKkma\nNWuW3G63EhMTlZGRoePHj0uSjh07pvT0dEVFRWnatGkt+vL5fJo6dar69esnt9utN954Q5L04osv\nKiEhQcnJyRo6dKjKy8uNHgYAAIAlGB7WwsPDtWTJElVUVKikpETLly9XZWWlRo4cqYqKCpWXl6tv\n376aP3++JCkyMlJz587VwoULW/WVm5srp9Opffv2qbKyUnfccYckafz48frwww9VVlamp556Sk88\n8YTRwwAAALAEwy+K63Q65XSevfqww+GQ2+3W4cOHNWLECP9rUlJSVFhYKEmy2+1KTU3V/v37W/W1\natUq7du3z/84OjpakhQVFeVvq6+vV7du3YweBgAAgCUE9A4GXq9XZWVlSklJadG+cuVKjR07tkWb\nzWZr8bi2tlaSNGfOHBUXF6tPnz5atmyZYmJiJEkrVqzQ4sWLdeLECe3atavN9RcXFOvYf46puLRY\nriSXXEkug0YGAADQcd4PvPJ+4L2k1wbsBIP6+nplZWUpPz9fDofD356bm6uIiAiNGzfugss3NDSo\nurpaqamp2rNnj4YOHaqZM2f6n3/sscf02WefafHixZoyZUqbfXgme3TdPdfJM9lDUAMAAJbhSnLJ\nM9nj/7mQgIS106dPKzMzUxMmTNDo0aP97QUFBdqwYYNeeeWVi/YRHR0tu92ujIwMSVJWVpZKS0tb\nvW7MmDFttgMAAIQCw8NaU1OTsrOzFRcXp5ycHH/7pk2blJeXp6KiIkVGRra53LlsNptGjRqlbdu2\nSZK2bt2q+Ph4SWrx/bb169crISHB6GEAAABYguHfWdu5c6fWrFnjv7SGJM2bN0/Tp0+Xz+fzn2gw\ndOhQrVixQpLkcrlUV1cnn8+noqIibd68WbGxsVqwYIEmTpyonJwcxcTEaNWqVZKk5cuXa8uWLQoP\nD1f37t397QAAAKHG8LCWlpamxsbGVu1tne3ZzOv1ttneq1cvbd++vVX7888/3+H6AAAAggl3MAAA\nALAwwhoAAICFEdYAAAAsjLAGAABgYQG9g4HZvLVedbN3k7fWa1ifRvfXFqfj7O26Ar2eZpdjTEZo\nrrN5fpqd2xYM45A6NudW+jtZqZZvy+jag23bM9KVPPa2MB+tdeT9xjy2Zms6/wJnIcJms7W6dhsA\nAIAVXSi3cBgUAADAwghrAAAAFkZYAwAAsDDCGgAAgIUR1gAAACyMsHYFKS4uNrsEy2FOWmNOWmNO\nWmNOWmNOWmNO2vZt54WwdgXhTdMac9Iac9Iac9Iac9Iac9Iac9I2whoAAEAIIawBAABYWEjfwQAA\nACBYtBfJQvbeoCGaQQEAwBWGw6AAAAAWRlgDAACwMMIaAACAhYVkWNu0aZNiY2N1yy23aMGCBWaX\nY7opU6bo+uuv14ABA8wuxTKqqqqUnp6u+Ph49e/fX0uXLjW7JNP973//U0pKipKSkhQXF6fZs2eb\nXZJlnDlzRsnJyRo1apTZpViCy+VSQkKCkpOTdeutt5pdjiXU1tYqKytLbrdbcXFxKikpMbsk0+3b\nt0/Jycn+n65du/K/VtL8+fMVHx+vAQMGaNy4cfrmm28uukzInQ165swZ9evXT1u2bFGPHj00ZMgQ\nvfrqq3K73WaXZpr33ntPDodDkyZN0kcffWR2OZZQU1OjmpoaJSUlqb6+XoMGDdKbb755RW8nknTy\n5EnZ7XY1NDQoLS1NCxcuVFpamtllmW7x4sXas2eP6urq9NZbb5ldjul69+6tPXv26LrrrjO7FMt4\n5JFHNGzYME2ZMkUNDQ06ceKEunbtanZZltHY2KgePXpo9+7duummm8wuxzRer1d33nmnKisr1aVL\nF40ZM0b33nuvHnnkkQsuF3J71nbv3q2bb75ZLpdL4eHhevjhh1VUVGR2Waa6/fbbde2115pdhqU4\nnU4lJSVJkhwOh9xutw4fPmxyVeaz2+2SJJ/PpzNnzvBhLKm6ulobNmzQT37yE84yPwdz8f8dP35c\n7733nqZMmSJJ6ty5M0HtPFu2bFGfPn2u6KAmSVdffbXCw8N18uRJNTQ06OTJk+rRo8dFlwu5sHbo\n0KEWG0PPnj116NAhEyuC1Xm9XpWVlSklJcXsUkzX2NiopKQkXX/99UpPT1dcXJzZJZluxowZysvL\nU1hYyP277DCbzabhw4dr8ODB+sMf/mB2OaY7ePCgunfvrkcffVQDBw7UT3/6U508edLssixl7dq1\nGjdunNllmO66667TE088oV69eunGG2/UNddco+HDh190uZD778PFcPFt1NfXKysrS/n5+XI4HGaX\nY7qwsDB98MEHqq6u1rvvvnvF39fvL3/5i2JiYpScnMyepHPs3LlTZWVl2rhxo5YvX6733nvP7JJM\n1dDQoNLSUj322GMqLS3V9773PT333HNml2UZPp9Pb7/9th566CGzSzHdgQMH9Pzzz8vr9erw4cOq\nr6/XK6+8ctHlQi6s9ejRQ1VVVf7HVVVV6tmzp4kVwapOnz6tzMxMTZgwQaNHjza7HEvp2rWr7rvv\nPr3//vtml2KqXbt26a233lLv3r01duxY/e1vf9OkSZPMLst0N9xwgySpe/fuevDBB7V7926TKzJX\nz5491bNnTw0ZMkSSlJWVpdLSUpOrso6NGzdq0KBB6t69u9mlmO7999/XD37wA0VHR6tz587KyMjQ\nrl27LrpcyIW1wYMHa//+/fJ6vfL5fFq3bp3uv/9+s8uCxTQ1NSk7O1txcXHKyckxuxxLOHr0qGpr\nayVJp06d0jvvvKPk5GSTqzLXvHnzVFVVpYMHD2rt2rW688479fLLL5tdlqlOnjypuro6SdKJEye0\nefPmK/5Mc6fTqZtuukmffvqppLPfz4qPjze5Kut49dVXNXbsWLPLsITY2FiVlJTo1KlTampq0pYt\nWy7p6yYhd7upzp07a9myZfrhD3+oM2fOKDs7+4o/w2/s2LHavn27vvzyS91000367W9/q0cffdTs\nsky1c+dOrVmzxn/5Aens6dR33323yZWZ54svvtAjjzyixsZGNTY2auLEibrrrrvMLstS+JqF9J//\n/EcPPvigpLOH/8aPH6+RI0eaXJX5XnjhBY0fP14+n099+vTRqlWrzC7JEk6cOKEtW7bw3cb/JzEx\nUZMmTdLgwYMVFhamgQMHaurUqRddLuQu3QEAABBKQu4wKAAAQCghrAEAAFgYYQ0AAMDCCGsAAAAW\nRlgDgAAqLi7mBvAAvhPCGgAAgIUR1gCEtBMnTui+++5TUlKSBgwYoNdee0179uyRx+PR4MGDdffd\nd6umpkaS9Nlnn2n48OFKSkrSoEGDdPDgQUnSrFmzNGDAACUkJOi1116TdHaPmcfj0UMPPSS3260J\nEyb417lp0ya53W4NGjRIb7zxhr99+/btSk5OVnJysgYOHKj6+vrLOBMAglXIXRQXAM61adMm9ejR\nQ+vXr5ckff3117rnnnv01ltvKTo6WuvWrdOvfvUrvfTSSxo/fryeeuopPfDAA/L5fDpz5owKCwtV\nXl6uDz/8UEeOHNGQIUN0xx13SJI++OADffLJJ7rhhhuUmpqqXbt2+S9yuW3bNvXp00djxozxX0x3\n0aJFWrFihYYOHaqTJ0+qS5cups0LgODBnjUAIS0hIUHvvPOOnnzySe3YsUOff/65Pv74Yw0fPlzJ\nycnKzc3VoUOHVF9fr8OHD+uBBx6QJEVEROiqq67Szp07NW7cONlsNsXExGjYsGH65z//KZvNpltv\nvVU33nijbDabkpKSdPDgQe3du1e9e/dWnz59JEkTJkzw3wQ+NTVVM2bM0AsvvKCvvvpKnTp1Mm1e\nAAQP9qwBCGm33HKLysrKtH79es2ZM0fp6emKj49vdfPk5vtdtuX8G7007yk7d89Yp06d1NDQ0OqW\nVOcu+8tf/lI/+tGPtH79eqWmpuqvf/2r+vXr1+GxAbgysGcNQEj74osvFBkZqfHjx2vmzJnavXu3\njh49qpKSEknS6dOn9cknnygqKko9e/ZUUVGRJOmbb77RqVOndPvtt2vdunVqbGzUkSNH9O677+rW\nW29tFeCksyEuNjZWXq9X//rXvySdvYl1swMHDig+Pl6/+MUvNGTIEO3bt+8yzACAYMeeNQAh7aOP\nPtKsWbMUFhamiIgI/e53v1OnTp00ffp0HT9+XA0NDZoxY4bi4uL0xz/+UT/72c/0zDPPKDw8XH/+\n85/14IMP6u9//7sSExNls9mUl5enmJgYVVZWtnlj9y5duuj3v/+97rvvPtntdt1+++06ceKEJCk/\nP1/btm1TWFiY+vfvr3vuuedyTweAIMSN3AEAACyMw6AAAAAWRlgDAACwMMIaAACAhRHWAAAALIyw\nBgAAYGGENQAAAAv7P0QKTpCgqK3+AAAAAElFTkSuQmCC\n"
      }
     ],
     "prompt_number": 24
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "In the example above we can see that the engine cluster is a bit more efficiently used, and the time to completion is shorter than in the previous example."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Further reading\n",
      "\n",
      "There are many other ways to use the IPython parallel environment. The official documentation has a nice guide:\n",
      "\n",
      "* http://ipython.org/ipython-doc/dev/parallel/"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## MPI \n",
      "\n",
      "When more communication between processes is required, sophisticated solutions such as MPI and OpenMP are often needed. MPI is process based parallel processing library/protocol, and can be used in Python programs through the `mpi4py` package:\n",
      "\n",
      "http://mpi4py.scipy.org/\n",
      "\n",
      "To use the `mpi4py` package we include `MPI` from `mpi4py`:\n",
      "\n",
      "    from mpi4py import MPI\n",
      "\n",
      "A MPI python program must be started using the `mpirun -n N` command, where `N` is the number of processes that should be included in the process group.\n",
      "\n",
      "Note that the IPython parallel enviroment also has support for MPI, but to begin with we will use `mpi4py` and the `mpirun` in the follow examples."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Example 1"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%file mpitest.py\n",
      "\n",
      "from mpi4py import MPI\n",
      "\n",
      "comm = MPI.COMM_WORLD\n",
      "rank = comm.Get_rank()\n",
      "\n",
      "if rank == 0:\n",
      "   data = [1.0, 2.0, 3.0, 4.0]\n",
      "   comm.send(data, dest=1, tag=11)\n",
      "elif rank == 1:\n",
      "   data = comm.recv(source=0, tag=11)\n",
      "    \n",
      "print \"rank =\", rank, \", data =\", data"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Overwriting mpitest.py\n"
       ]
      }
     ],
     "prompt_number": 25
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "!mpirun -n 2 python mpitest.py"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "rank = 0 , data = [1.0, 2.0, 3.0, 4.0]\r\n",
        "rank = 1 , data = [1.0, 2.0, 3.0, 4.0]\r\n"
       ]
      }
     ],
     "prompt_number": 26
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Example 2\n",
      "\n",
      "Send a numpy array from one process to another:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%file mpi-numpy-array.py\n",
      "\n",
      "from mpi4py import MPI\n",
      "import numpy\n",
      "\n",
      "comm = MPI.COMM_WORLD\n",
      "rank = comm.Get_rank()\n",
      "\n",
      "if rank == 0:\n",
      "   data = numpy.random.rand(10)\n",
      "   comm.Send(data, dest=1, tag=13)\n",
      "elif rank == 1:\n",
      "   data = numpy.empty(10, dtype=numpy.float64)\n",
      "   comm.Recv(data, source=0, tag=13)\n",
      "    \n",
      "print \"rank =\", rank, \", data =\", data"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Overwriting mpi-numpy-array.py"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "\n"
       ]
      }
     ],
     "prompt_number": 27
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "!mpirun -n 2 python mpi-numpy-array.py"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "rank = 1 , data = [ 0.29017907  0.1301338   0.0335738   0.61340947  0.69981434  0.0110428\r\n",
        "  0.40157161  0.04627437  0.65849675  0.49130351]\r\n",
        "rank = 0 , data = [ 0.29017907  0.1301338   0.0335738   0.61340947  0.69981434  0.0110428\r\n",
        "  0.40157161  0.04627437  0.65849675  0.49130351]\r\n"
       ]
      }
     ],
     "prompt_number": 28
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Example 3"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%file mpi-matrix-vector.py\n",
      "# example from: http://mpi4py.scipy.org/docs/usrman/tutorial.html\n",
      "\n",
      "from mpi4py import MPI\n",
      "import numpy\n",
      "\n",
      "def matvec(comm, A, x):\n",
      "    m = A.shape[0] # local rows\n",
      "    p = comm.Get_size()\n",
      "    xg = numpy.zeros(m*p, dtype='d')\n",
      "    comm.Allgather([x,  MPI.DOUBLE], [xg, MPI.DOUBLE])\n",
      "    y = numpy.dot(A, xg)\n",
      "    return y\n",
      "\n",
      "comm = MPI.COMM_WORLD\n",
      "rank = comm.Get_rank()\n",
      "\n",
      "A = numpy.random.rand(5,5)\n",
      "x = numpy.random.rand(5)\n",
      "y_mpi = matvec(comm, A, x)\n",
      "    \n",
      "if rank == 0:\n",
      "\n",
      "    y = numpy.dot(A, x)\n",
      "\n",
      "    print \"y =\", y\n",
      "\n",
      "    print \"y_mpi =\", y_mpi"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Overwriting mpi-matrix-vector.py"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "\n"
       ]
      }
     ],
     "prompt_number": 29
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "!mpirun -n 1 python mpi-matrix-vector.py"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "y = [ 1.17583156  1.21950954  1.78492742  1.10535307  1.42770478]\r\n",
        "y_mpi = [ 1.17583156  1.21950954  1.78492742  1.10535307  1.42770478]\r\n"
       ]
      }
     ],
     "prompt_number": 30
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Example 4"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# prepare some random data\n",
      "a = numpy.random.rand(100)\n",
      "numpy.save(\"random-vector.npy\", a)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 31
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%file mpi-psum.py\n",
      "\n",
      "from mpi4py import MPI\n",
      "import numpy as np\n",
      "\n",
      "def psum(a):\n",
      "    r = MPI.COMM_WORLD.Get_rank()\n",
      "    size = MPI.COMM_WORLD.Get_size()\n",
      "    m = len(a) / size\n",
      "    locsum = np.sum(a[r*m:(r+1)*m])\n",
      "    rcvBuf = np.array(0.0, 'd')\n",
      "    MPI.COMM_WORLD.Allreduce([locsum, MPI.DOUBLE], [rcvBuf, MPI.DOUBLE], op=MPI.SUM)\n",
      "    return rcvBuf\n",
      "\n",
      "a = np.load(\"random-vector.npy\")\n",
      "s = psum(a)\n",
      "\n",
      "if MPI.COMM_WORLD.Get_rank() == 0:\n",
      "    print \"sum =\", s, \", numpy sum =\", a.sum()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Overwriting mpi-psum.py"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "\n"
       ]
      }
     ],
     "prompt_number": 32
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "!mpirun -n 4 python mpi-psum.py"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "sum = 51.8048933094 , numpy sum = 51.8048933094\r\n"
       ]
      }
     ],
     "prompt_number": 33
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Further reading\n",
      "\n",
      "* http://mpi4py.scipy.org\n",
      "\n",
      "* http://mpi4py.scipy.org/docs/usrman/tutorial.html\n",
      "\n",
      "* https://computing.llnl.gov/tutorials/mpi/"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## OpenMP\n",
      "\n",
      "What about OpenMP? OpenMP is a standard and widely used thread-based parallel API that unfortunaltely is **not** useful in Python. The reason is that the CPython implementation use a global interpreter lock, making it impossible to simultaneously run several Python threads. Threads are therefore not useful for parallel computing in Python, unless it is only used to wrap compiled code that do the OpenMP parallelization (Numpy does something like that). \n",
      "\n",
      "This is clearly a limitation in the Python interpreter, and as a consequence all parallelization in Python must use processes (not threads)."
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## OpenCL\n",
      "\n",
      "OpenCL is an API for heterogenous computing, for example using GPUs for numerical computations. There is a python package called `pyopencl` that allows OpenCL code to be compiled, loaded and executed on the compute units completely from within Python. This is a nice way to work with OpenCL, because the time-consuming computations should be done on the compute units in compiled code, and in this Python only server as a control language. "
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%%file opencl-dense-mv.py\n",
      "\n",
      "import pyopencl as cl\n",
      "import numpy\n",
      "import time\n",
      "\n",
      "# problem size\n",
      "n = 10000\n",
      "\n",
      "# platform\n",
      "platform_list = cl.get_platforms()\n",
      "platform = platform_list[0]\n",
      "\n",
      "# device\n",
      "device_list = platform.get_devices()\n",
      "device = device_list[0]\n",
      "\n",
      "if False:\n",
      "    print(\"Platform name:\" + platform.name)\n",
      "    print(\"Platform version:\" + platform.version)\n",
      "    print(\"Device name:\" + device.name)\n",
      "    print(\"Device type:\" + cl.device_type.to_string(device.type))\n",
      "    print(\"Device memory: \" + str(device.global_mem_size//1024//1024) + ' MB')\n",
      "    print(\"Device max clock speed:\" + str(device.max_clock_frequency) + ' MHz')\n",
      "    print(\"Device compute units:\" + str(device.max_compute_units))\n",
      "\n",
      "# context\n",
      "ctx = cl.Context([device]) # or we can use cl.create_some_context()\n",
      "\n",
      "# command queue\n",
      "queue = cl.CommandQueue(ctx)\n",
      "\n",
      "# kernel\n",
      "KERNEL_CODE = \"\"\"\n",
      "//\n",
      "// Matrix-vector multiplication: r = m * v\n",
      "//\n",
      "#define N %(mat_size)d\n",
      "__kernel\n",
      "void dmv_cl(__global float *m, __global float *v, __global float *r)\n",
      "{\n",
      "    int i, gid = get_global_id(0);\n",
      "    \n",
      "    r[gid] = 0;\n",
      "    for (i = 0; i < N; i++)\n",
      "    {\n",
      "        r[gid] += m[gid * N + i] * v[i];\n",
      "    }\n",
      "}\n",
      "\"\"\"\n",
      "\n",
      "kernel_params = {\"mat_size\": n}\n",
      "program = cl.Program(ctx, KERNEL_CODE % kernel_params).build()\n",
      "\n",
      "# data\n",
      "A = numpy.random.rand(n, n)\n",
      "x = numpy.random.rand(n, 1)\n",
      "\n",
      "# host buffers\n",
      "h_y = numpy.empty(numpy.shape(x)).astype(numpy.float32)\n",
      "h_A = numpy.real(A).astype(numpy.float32)\n",
      "h_x = numpy.real(x).astype(numpy.float32)\n",
      "\n",
      "# device buffers\n",
      "mf = cl.mem_flags\n",
      "d_A_buf = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=h_A)\n",
      "d_x_buf = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=h_x)\n",
      "d_y_buf = cl.Buffer(ctx, mf.WRITE_ONLY, size=h_y.nbytes)\n",
      "\n",
      "# execute OpenCL code\n",
      "t0 = time.time()\n",
      "event = program.dmv_cl(queue, h_y.shape, None, d_A_buf, d_x_buf, d_y_buf)\n",
      "event.wait()\n",
      "cl.enqueue_copy(queue, h_y, d_y_buf)\n",
      "t1 = time.time()\n",
      "\n",
      "print \"opencl elapsed time =\", (t1-t0)\n",
      "\n",
      "# Same calculation with numpy\n",
      "t0 = time.time()\n",
      "y = numpy.dot(h_A, h_x)\n",
      "t1 = time.time()\n",
      "\n",
      "print \"numpy elapsed time =\", (t1-t0)\n",
      "\n",
      "# see if the results are the same\n",
      "print \"max deviation =\", numpy.abs(y-h_y).max()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Overwriting opencl-dense-mv.py"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "\n"
       ]
      }
     ],
     "prompt_number": 38
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "!python opencl-dense-mv.py"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "/usr/local/lib/python2.7/dist-packages/pyopencl-2012.1-py2.7-linux-x86_64.egg/pyopencl/__init__.py:36: CompilerWarning: Non-empty compiler output encountered. Set the environment variable PYOPENCL_COMPILER_OUTPUT=1 to see more.\r\n",
        "  \"to see more.\", CompilerWarning)\r\n"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "opencl elapsed time = 0.110778093338\r\n"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "numpy elapsed time = 0.141111135483\r\n",
        "max deviation = 0.0119629\r\n"
       ]
      }
     ],
     "prompt_number": 39
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "### Further reading\n",
      "\n",
      "* http://mathema.tician.de/software/pyopencl"
     ]
    }
   ],
   "metadata": {}
  }
 ]
}