lectures/lecture1/Performance-variation.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Performance variance in a simple problem\n",
    "\n",
    "A key thing to understand is that **computers are not a box that does arithmetic at a fixed rate**.  There is a **huge variance** in performance depending on **how the \"same\" algorithm is implemented**, both *within a single language* and when *comparing different languages*.\n",
    "\n",
    "First, we will look at a simple algorithm for a simple function — summation of *n* numbers — and see how the \"same\" calculation can vary in performance by **orders of magnitude**.\n",
    "\n",
    "The next task is to understand *why* this variation occurs, which is essential to making sure that your own code is fast and not slow!    More detail can be found in the accompanying [Boxes and registers](Boxes-and-registers.ipynb) notes."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# `sum`: An easy enough function to understand"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Consider the  **sum** function `sum(a)`, which computes\n",
    "$$\n",
    "\\mathrm{sum}(a) = \\sum_{i=1}^n a_i,\n",
    "$$\n",
    "where $n$ is the length of `a`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "10000000-element Array{Float64,1}:\n",
       " 0.271846  \n",
       " 0.358779  \n",
       " 0.890233  \n",
       " 0.0790911 \n",
       " 0.00939895\n",
       " 0.767968  \n",
       " 0.308546  \n",
       " 0.322207  \n",
       " 0.925622  \n",
       " 0.798759  \n",
       " 0.964225  \n",
       " 0.89269   \n",
       " 0.347821  \n",
       " ⋮         \n",
       " 0.273381  \n",
       " 0.121823  \n",
       " 0.573685  \n",
       " 0.960718  \n",
       " 0.598456  \n",
       " 0.0506493 \n",
       " 0.18194   \n",
       " 0.726292  \n",
       " 0.14169   \n",
       " 0.384121  \n",
       " 0.553164  \n",
       " 0.838971  "
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a = rand(10^7) # 1D vector of random numbers, uniform on [0,1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5.000180786929853e6"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sum(a)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The expected result is 0.5 * 10^7, since the mean of each entry is 0.5"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Benchmarking a few ways in a few languages"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Julia has a `BenchmarkTools.jl` package for easy and accurate benchmarking:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "#Pkg.add(\"BenchmarkTools\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true,
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "using BenchmarkTools, Compat"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#  1. The C language\n",
    "\n",
    "C is often considered the gold standard: difficult on the human, nice for the machine. Getting within a factor of 2 of C is often satisfying. Nonetheless, even within C, there are many kinds of optimizations possible that a naive C writer may or may not get the advantage of.\n",
    "\n",
    "The current author does not speak C, so he does not read the cell below, but is happy to know that you can put C code in a Julia session, compile it, and run it. Note that the `\"\"\"` wrap a multi-line string."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "c_sum (generic function with 1 method)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "C_code = \"\"\"\n",
    "#include <stddef.h>\n",
    "double c_sum(size_t n, double *X) {\n",
    "    double s = 0.0;\n",
    "    for (size_t i = 0; i < n; ++i) {\n",
    "        s += X[i];\n",
    "    }\n",
    "    return s;\n",
    "}\n",
    "\"\"\"\n",
    "\n",
    "const Clib = tempname()   # make a temporary file\n",
    "\n",
    "\n",
    "# compile to a shared library by piping C_code to gcc\n",
    "# (works only if you have gcc installed):\n",
    "\n",
    "open(`gcc -fPIC -O3 -msse3 -xc -shared -o $(Clib * \".\" * Libdl.dlext) -`, \"w\") do f\n",
    "    print(f, C_code) \n",
    "end\n",
    "\n",
    "# define a Julia function that calls the C function:\n",
    "c_sum(X::Array{Float64}) = ccall((\"c_sum\", Clib), Float64, (Csize_t, Ptr{Float64}), length(X), X)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5.000180786930093e6"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "c_sum(a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "c_sum(a) ≈ sum(a) # type \\approx and then <TAB> to get the ≈ symbolb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2.39349901676178e-7"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "c_sum(a) - sum(a)  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "isapprox (generic function with 3 methods)"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "≈  # alias for the `isapprox` function"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "search: \u001b[1mi\u001b[22m\u001b[1ms\u001b[22m\u001b[1ma\u001b[22m\u001b[1mp\u001b[22m\u001b[1mp\u001b[22m\u001b[1mr\u001b[22m\u001b[1mo\u001b[22m\u001b[1mx\u001b[22m\n",
      "\n"
     ]
    },
    {
     "data": {
      "text/markdown": [
       "```\n",
       "isapprox(x, y; rtol::Real=sqrt(eps), atol::Real=0, nans::Bool=false, norm::Function)\n",
       "```\n",
       "\n",
       "Inexact equality comparison: `true` if `norm(x-y) <= atol + rtol*max(norm(x), norm(y))`. The default `atol` is zero and the default `rtol` depends on the types of `x` and `y`. The keyword argument `nans` determines whether or not NaN values are considered equal (defaults to false).\n",
       "\n",
       "For real or complex floating-point values, `rtol` defaults to `sqrt(eps(typeof(real(x-y))))`. This corresponds to requiring equality of about half of the significand digits. For other types, `rtol` defaults to zero.\n",
       "\n",
       "`x` and `y` may also be arrays of numbers, in which case `norm` defaults to `vecnorm` but may be changed by passing a `norm::Function` keyword argument. (For numbers, `norm` is the same thing as `abs`.) When `x` and `y` are arrays, if `norm(x-y)` is not finite (i.e. `±Inf` or `NaN`), the comparison falls back to checking whether all elements of `x` and `y` are approximately equal component-wise.\n",
       "\n",
       "The binary operator `≈` is equivalent to `isapprox` with the default arguments, and `x ≉ y` is equivalent to `!isapprox(x,y)`.\n",
       "\n",
       "```jldoctest\n",
       "julia> 0.1 ≈ (0.1 - 1e-10)\n",
       "true\n",
       "\n",
       "julia> isapprox(10, 11; atol = 2)\n",
       "true\n",
       "\n",
       "julia> isapprox([10.0^9, 1.0], [10.0^9, 2.0])\n",
       "true\n",
       "```\n"
      ],
      "text/plain": [
       "```\n",
       "isapprox(x, y; rtol::Real=sqrt(eps), atol::Real=0, nans::Bool=false, norm::Function)\n",
       "```\n",
       "\n",
       "Inexact equality comparison: `true` if `norm(x-y) <= atol + rtol*max(norm(x), norm(y))`. The default `atol` is zero and the default `rtol` depends on the types of `x` and `y`. The keyword argument `nans` determines whether or not NaN values are considered equal (defaults to false).\n",
       "\n",
       "For real or complex floating-point values, `rtol` defaults to `sqrt(eps(typeof(real(x-y))))`. This corresponds to requiring equality of about half of the significand digits. For other types, `rtol` defaults to zero.\n",
       "\n",
       "`x` and `y` may also be arrays of numbers, in which case `norm` defaults to `vecnorm` but may be changed by passing a `norm::Function` keyword argument. (For numbers, `norm` is the same thing as `abs`.) When `x` and `y` are arrays, if `norm(x-y)` is not finite (i.e. `±Inf` or `NaN`), the comparison falls back to checking whether all elements of `x` and `y` are approximately equal component-wise.\n",
       "\n",
       "The binary operator `≈` is equivalent to `isapprox` with the default arguments, and `x ≉ y` is equivalent to `!isapprox(x,y)`.\n",
       "\n",
       "```jldoctest\n",
       "julia> 0.1 ≈ (0.1 - 1e-10)\n",
       "true\n",
       "\n",
       "julia> isapprox(10, 11; atol = 2)\n",
       "true\n",
       "\n",
       "julia> isapprox([10.0^9, 1.0], [10.0^9, 2.0])\n",
       "true\n",
       "```\n"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "?isapprox"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now benchmark the C code directly from Julia:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "BenchmarkTools.Trial: \n",
       "  memory estimate:  0 bytes\n",
       "  allocs estimate:  0\n",
       "  --------------\n",
       "  minimum time:     9.856 ms (0.00% GC)\n",
       "  median time:      10.049 ms (0.00% GC)\n",
       "  mean time:        10.233 ms (0.00% GC)\n",
       "  maximum time:     12.929 ms (0.00% GC)\n",
       "  --------------\n",
       "  samples:          488\n",
       "  evals/sample:     1"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "c_bench = @benchmark c_sum($a) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "C: Fastest time was 9.855969 msec\n"
     ]
    }
   ],
   "source": [
    "println(\"C: Fastest time was $(minimum(c_bench.times) / 1e6) msec\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Dict{Any,Any} with 1 entry:\n",
       "  \"C\" => 9.85597"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "d = Dict()  # a \"dictionary\", i.e. an associative array\n",
    "d[\"C\"] = minimum(c_bench.times) / 1e6  # in milliseconds\n",
    "d"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "#Pkg.add(\"PyPlot\")\n",
    "using PyPlot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjYAAAG0CAYAAAAhJm17AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4wLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvpW3flQAAIABJREFUeJzt3Xl0VFW6/vGnSEgRMYkGyCSzLTOG0SuDEtAGwyC264IoFxGEBYpgjI0kP0STq5Jgt14UGlBsgeu8WoTGRmnoloCIIENQVBoEg0SFGwdICGgYsn9/9KKWZQaSolJ1avP9rFVrcfbZ59T71uHAw64q4jLGGAEAAFigXrALAAAA8BeCDQAAsAbBBgAAWINgAwAArEGwAQAA1iDYAAAAaxBsAACANQg2AADAGgQbAABgDYINAACwBsEGAABYIzzYBdS18vJyffvtt4qKipLL5Qp2OQAAoAaMMTp+/LiSkpJUr17N12GsDzbffvutmjVrFuwyAACADwoLC9W0adMaz7c+2ERFRUn69wsTHR0d5GoAAEBNlJSUqFmzZp6/x2vK+mBz7u2n6Ohogg0AACGmth8j4cPDAADAGgQbAABgDYINAACwBsEGAABYg2ADAACsQbABAADWINgAAABrEGwAAIA1CDYAAMAaBBsAAGANgg0AALAGwQYAAFiDYAMAAKxBsAEAANYg2AAAAGuEB7sA2K1lxuoq9x3MHRLASgAAFwNWbAAAgDUINgAAwBoEGwAAYA2CDQAAsAbBBgAAWINgAwAArEGwAQAA1iDYAAAAaxBsAACANQg2AADAGgQbAABgDYINAACwBsEGAABYg2ADAACsQbABAADWCGqw2bhxo4YNG6akpCS5XC6tXLmywpw9e/bo5ptvVkxMjKKionTttdfq0KFDQagWAAA4XVCDzYkTJ5ScnKz58+dXuv/AgQPq27ev2rVrp7y8PH388ceaNWuWGjRoEOBKAQBAKAgP5pOnpqYqNTW1yv0zZ87U4MGD9eSTT3rGWrduHYjSAABACHLsZ2zKy8u1evVqtWnTRoMGDVJcXJz+4z/+o9K3q36prKxMJSUlXg8AAHBxCOqKTXWKiopUWlqq3NxcPf7445ozZ47WrFmjW2+9VevXr1e/fv0qPS4nJ0fZ2dkBrhYXo5YZq6vcdzB3SAArAQCc4+gVG0kaPny4HnjgAXXp0kUZGRkaOnSoFi1aVOVxmZmZKi4u9jwKCwsDVTIAAAgyx67YNG7cWOHh4erQoYPXePv27bVp06Yqj3O73XK73XVdHgAAcCDHrthERESoZ8+e2rt3r9f4vn371KJFiyBVBQAAnCyoKzalpaXav3+/Z7ugoEC7du1SbGysmjdvrunTp+u2227T9ddfr/79+2vNmjV6++23lZeXF7yiAQCAYwU12Gzfvl39+/f3bKenp0uSxo4dq6VLl+p3v/udFi1apJycHE2bNk1t27bV8uXL1bdv32CVDAAAHCyowSYlJUXGmGrnjB8/XuPHjw9QRQAAIJQ59jM2AAAAtUWwAQAA1iDYAAAAaxBsAACANQg2AADAGgQbAABgDYINAACwBsEGAABYg2ADAACsQbABAADWINgAAABrEGwAAIA1CDYAAMAaBBsAAGCN8GAXAODCtcxYXeW+g7lDAlgJAAQXKzYAAMAaBBsAAGANgg0AALAGwQYAAFiDYAMAAKxBsAEAANYg2AAAAGsQbAAAgDUINgAAwBoEGwAAYA2CDQAAsAbBBgAAWINgAwAArEGwAQAA1iDYAAAAaxBsAACANYIabDZu3Khhw4YpKSlJLpdLK1eurHLupEmT5HK5NHfu3ABWCAAAQklQg82JEyeUnJys+fPnVztv5cqV2rp1q5KSkgJUGQAACEXhwXzy1NRUpaamVjvnm2++0X333ae///3vGjJkSIAqAwAAoSioweZ8ysvLNWbMGE2fPl0dO3as0TFlZWUqKyvzbJeUlNRVeQAAwGEcHWzmzJmj8PBwTZs2rcbH5OTkKDs7uw6rAkJLy4zVVe47mMsqKAC7OPZbUTt27NAzzzyjpUuXyuVy1fi4zMxMFRcXex6FhYV1WCUAAHASxwab999/X0VFRWrevLnCw8MVHh6ur776Sg8++KBatmxZ5XFut1vR0dFeDwAAcHFw7FtRY8aM0Y033ug1NmjQII0ZM0bjxo0LUlUAAMDJghpsSktLtX//fs92QUGBdu3apdjYWDVv3lyNGjXyml+/fn0lJCSobdu2gS4VAACEgKAGm+3bt6t///6e7fT0dEnS2LFjtXTp0iBVBQAAQlVQg01KSoqMMTWef/DgwborBgAAhDzHfngYAACgtgg2AADAGgQbAABgDYINAACwBsEGAABYg2ADAACsQbABAADWINgAAABrEGwAAIA1CDYAAMAaBBsAAGANgg0AALAGwQYAAFiDYAMAAKwRHuwCgNpombG6yn0Hc4cEsBIAgBOxYgMAAKxBsAEAANYg2AAAAGsQbAAAgDUINgAAwBoEGwAAYA2CDQAAsAbBBgAAWINgAwAArEGwAQAA1iDYAAAAaxBsAACANQg2AADAGgQbAABgDYINAACwRlCDzcaNGzVs2DAlJSXJ5XJp5cqVnn2nT5/WjBkz1LlzZzVs2FBJSUm688479e233waxYgAA4GRBDTYnTpxQcnKy5s+fX2HfyZMntXPnTs2aNUs7d+7UW2+9pX379unmm28OQqUAACAUhAfzyVNTU5WamlrpvpiYGK1bt85rbN68ebrmmmt06NAhNW/ePBAlAgCAEBLUYFNbxcXFcrlcuuyyy6qcU1ZWprKyMs92SUlJIEoDAAAOEDLB5ueff1ZGRobuuOMORUdHVzkvJydH2dnZAawMuPi0zFhd5b6DuUMCWAkAeAuJb0WdPn1ao0aNUnl5uRYsWFDt3MzMTBUXF3sehYWFAaoSAAAEm+NXbE6fPq2RI0eqoKBA7733XrWrNZLkdrvldrsDVB0AAHASRwebc6Hmiy++0Pr169WoUaNglwQAABwsqMGmtLRU+/fv92wXFBRo165dio2NVVJSkv7zP/9TO3fu1N/+9jedPXtWR44ckSTFxsYqIiIiWGUDAACHCmqw2b59u/r37+/ZTk9PlySNHTtWWVlZWrVqlSSpS5cuXsetX79eKSkpAasTAACEhqAGm5SUFBljqtxf3T4AAIBfC4lvRQEAANQEwQYAAFiDYAMAAKxBsAEAANYg2AAAAGsQbAAAgDUINgAAwBoEGwAAYA2CDQAAsAbBBgAAWINgAwAArEGwAQAA1iDYAAAAaxBsAACANcKDXQBgo5YZq6vcdzB3SAArqV6o1AkANcWKDQAAsAbBBgAAWINgAwAArEGwAQAA1iDYAAAAaxBsAACANQg2AADAGgQbAABgDYINAACwBsEGAABYg2ADAACsQbABAADWINgAAABrEGwAAIA1CDYAAMAaBBsAAGCNoAabjRs3atiwYUpKSpLL5dLKlSu99htjlJWVpaSkJEVGRiolJUWfffZZkKoFAABOF9Rgc+LECSUnJ2v+/PmV7n/yySf19NNPa/78+dq2bZsSEhL029/+VsePHw9wpQAAIBSEB/PJU1NTlZqaWuk+Y4zmzp2rmTNn6tZbb5UkLVu2TPHx8Xr11Vc1adKkQJYKAABCgGM/Y1NQUKAjR45o4MCBnjG3261+/fpp8+bNVR5XVlamkpISrwcAALg4BHXFpjpHjhyRJMXHx3uNx8fH66uvvqryuJycHGVnZ9dpbU7QMmN1lfsO5g4JYCXOESqviQ11AoBTOXbF5hyXy+W1bYypMPZLmZmZKi4u9jwKCwvrukQAAOAQjl2xSUhIkPTvlZvExETPeFFRUYVVnF9yu91yu911Xh8AAHAex67YtGrVSgkJCVq3bp1n7NSpU9qwYYN69+4dxMoAAIBTBXXFprS0VPv37/dsFxQUaNeuXYqNjVXz5s2Vlpam2bNn66qrrtJVV12l2bNn65JLLtEdd9wRxKoBAIBT+bRiM2DAAB07dqzCeElJiQYMGFDj82zfvl1du3ZV165dJUnp6enq2rWrHnnkEUnSQw89pLS0NN17773q0aOHvvnmG61du1ZRUVG+lA0AACzn04pNXl6eTp06VWH8559/1vvvv1/j86SkpMgYU+V+l8ulrKwsZWVl+VImAAC4yNQq2HzyySeeX3/++eeer2RL0tmzZ7VmzRpdccUV/qsOAACgFmoVbLp06SKXyyWXy1XpW06RkZGaN2+e34oDAACojVoFm4KCAhlj1Lp1a3300Udq0qSJZ19ERITi4uIUFhbm9yIBAABqolbBpkWLFpKk8vLyOikGAADgQvj8de99+/YpLy9PRUVFFYLOuW81AQAABJJPwWbx4sW655571LhxYyUkJHj9iAOXy0WwAQAAQeFTsHn88cf1xBNPaMaMGf6uBwAAwGc+/Qd9R48e1YgRI/xdCwAAwAXxKdiMGDFCa9eu9XctAAAAF8Snt6J+85vfaNasWdqyZYs6d+6s+vXre+2fNm2aX4oDAACoDZ+CzfPPP69LL71UGzZs0IYNG7z2uVwugg0AAAgKn4JNQUGBv+sAAAC4YD59xgYAAMCJfFqxGT9+fLX7X3zxRZ+KwcWlZcbqKvcdzB0Sss9VV6rrwUl8rTNUrgMAZ/Mp2Bw9etRr+/Tp0/r000917NixSn84JgAAQCD4FGxWrFhRYay8vFz33nuvWrdufcFFAQAA+MJvn7GpV6+eHnjgAf3P//yPv04JAABQK3798PCBAwd05swZf54SAACgxnx6Kyo9Pd1r2xijw4cPa/Xq1Ro7dqxfCgMAAKgtn4JNfn6+13a9evXUpEkTPfXUU+f9xhQAAEBd8SnYrF+/3t91AAAAXDCfgs053333nfbu3SuXy6U2bdqoSZMm/qoLAACg1nz68PCJEyc0fvx4JSYm6vrrr9d1112npKQk3X333Tp58qS/awQAAKgRn4JNenq6NmzYoLffflvHjh3TsWPH9Ne//lUbNmzQgw8+6O8aAQAAasSnt6KWL1+uN998UykpKZ6xwYMHKzIyUiNHjtTChQv9VR8AAECN+bRic/LkScXHx1cYj4uL460oAAAQND4Fm169eunRRx/Vzz//7Bn76aeflJ2drV69evmtOAAAgNrw6a2ouXPnKjU1VU2bNlVycrJcLpd27dolt9uttWvX+rtGAACAGvEp2HTu3FlffPGFXn75Zf3rX/+SMUajRo3S6NGjFRkZ6e8aAQAAasSnYJOTk6P4+HhNnDjRa/zFF1/Ud999pxkzZvilOAAAgNrw6TM2zz33nNq1a1dhvGPHjlq0aNEFF3XOmTNn9PDDD6tVq1aKjIxU69at9d///d8qLy/323MAAAB7+LRic+TIESUmJlYYb9KkiQ4fPnzBRZ0zZ84cLVq0SMuWLVPHjh21fft2jRs3TjExMbr//vv99jwAAMAOPgWbZs2a6YMPPlCrVq28xj/44AMlJSX5pTBJ+vDDDzV8+HANGTJEktSyZUu99tpr2r59u9+eAwAA2MOnYDNhwgSlpaXp9OnTGjBggCTpn//8px566CG//s/Dffv21aJFi7Rv3z61adNGH3/8sTZt2qS5c+dWeUxZWZnKyso82yUlJX6rBwAAOJtPweahhx7Sjz/+qHvvvVenTp2SJDVo0EAzZsxQZmam34qbMWOGiouL1a5dO4WFhens2bN64okndPvtt1d5TE5OjrKzs/1WAwAACB0+fXjY5XJpzpw5+u6777RlyxZ9/PHH+vHHH/XII4/4tbg33nhDL7/8sl599VXt3LlTy5Yt0x//+EctW7asymMyMzNVXFzseRQWFvq1JgAA4Fw+rdicc+mll6pnz57+qqWC6dOnKyMjQ6NGjZL07/8/56uvvlJOTo7Gjh1b6TFut1tut7vOagIAAM7l04pNoJw8eVL16nmXGBYWxte9AQBApS5oxaauDRs2TE888YSaN2+ujh07Kj8/X08//bTGjx8f7NIAAIADOTrYzJs3T7NmzdK9996roqIiJSUladKkSX7/LA8AALCDo4NNVFSU5s6dW+3XuwEAAM5x9GdsAAAAaoNgAwAArEGwAQAA1iDYAAAAaxBsAACANQg2AADAGgQbAABgDYINAACwBsEGAABYg2ADAACsQbABAADWINgAAABrEGwAAIA1CDYAAMAa4cEuAKhMy4zVwS5BknPquJhVdw0O5g5xzDmd9HzVcVItQF1gxQYAAFiDYAMAAKxBsAEAANYg2AAAAGsQbAAAgDUINgAAwBoEGwAAYA2CDQAAsAbBBgAAWINgAwAArEGwAQAA1iDYAAAAaxBsAACANQg2AADAGgQbAABgDccHm2+++Ub/9V//pUaNGumSSy5Rly5dtGPHjmCXBQAAHCg82AVU5+jRo+rTp4/69++vd999V3FxcTpw4IAuu+yyYJcGAAAcyNHBZs6cOWrWrJmWLFniGWvZsmXwCgIAAI7m6LeiVq1apR49emjEiBGKi4tT165dtXjx4mqPKSsrU0lJidcDAABcHBy9YvPll19q4cKFSk9P1//7f/9PH330kaZNmya3260777yz0mNycnKUnZ0d4Ert1zJjdZX7DuYOCWAlgDNVd48ACBxHr9iUl5erW7dumj17trp27apJkyZp4sSJWrhwYZXHZGZmqri42PMoLCwMYMUAACCYHB1sEhMT1aFDB6+x9u3b69ChQ1Ue43a7FR0d7fUAAAAXB0cHmz59+mjv3r1eY/v27VOLFi2CVBEAAHAyRwebBx54QFu2bNHs2bO1f/9+vfrqq3r++ec1ZcqUYJcGAAAcyNHBpmfPnlqxYoVee+01derUSY899pjmzp2r0aNHB7s0AADgQI7+VpQkDR06VEOHDg12GQAAIAQ4esUGAACgNgg2AADAGgQbAABgDYINAACwBsEGAABYg2ADAACsQbABAADWINgAAABrEGwAAIA1CDYAAMAaBBsAAGANgg0AALAGwQYAAFiDYAMAAKwRHuwC4H8tM1ZXue9g7pAAVgLUXHW/b305xtff69w/Fx+uuV1YsQEAANYg2AAAAGsQbAAAgDUINgAAwBoEGwAAYA2CDQAAsAbBBgAAWINgAwAArEGwAQAA1iDYAAAAaxBsAACANQg2AADAGgQbAABgDYINAACwBsEGAABYg2ADAACsEVLBJicnRy6XS2lpacEuBQAAOFDIBJtt27bp+eef19VXXx3sUgAAgEOFRLApLS3V6NGjtXjxYl1++eXBLgcAADhUSASbKVOmaMiQIbrxxhvPO7esrEwlJSVeDwAAcHEID3YB5/P6669r586d2rZtW43m5+TkKDs7u46rqp2WGaur3Hcwd0gAK3FWLU5R3WsSaE6qBYFRF9f8Yr2XUXs2/p3g6BWbwsJC3X///Xr55ZfVoEGDGh2TmZmp4uJiz6OwsLCOqwQAAE7h6BWbHTt2qKioSN27d/eMnT17Vhs3btT8+fNVVlamsLAwr2PcbrfcbnegSwUAAA7g6GBzww03aPfu3V5j48aNU7t27TRjxowKoQYAAFzcHB1soqKi1KlTJ6+xhg0bqlGjRhXGAQAAHP0ZGwAAgNpw9IpNZfLy8oJdAgAAcChWbAAAgDUINgAAwBoEGwAAYA2CDQAAsAbBBgAAWINgAwAArEGwAQAA1iDYAAAAaxBsAACANQg2AADAGgQbAABgDYINAACwBsEGAABYg2ADAACsER7sAlC1lhmrg11CjYRKncCF4vd68FV3DQ7mDglgJdVzUp0X2+9bVmwAAIA1CDYAAMAaBBsAAGANgg0AALAGwQYAAFiDYAMAAKxBsAEAANYg2AAAAGsQbAAAgDUINgAAwBoEGwAAYA2CDQAAsAbBBgAAWINgAwAArEGwAQAA1nB8sMnJyVHPnj0VFRWluLg43XLLLdq7d2+wywIAAA7k+GCzYcMGTZkyRVu2bNG6det05swZDRw4UCdOnAh2aQAAwGHCg13A+axZs8Zre8mSJYqLi9OOHTt0/fXXB6kqAADgRI4PNr9WXFwsSYqNja10f1lZmcrKyjzbJSUlAakLAAAEX0gFG2OM0tPT1bdvX3Xq1KnSOTk5OcrOzvb5OVpmrK5y38HcIT6fF0BgVXcv287X3n3988/X56vqnIGu34Zz1gV/X9dAcfxnbH7pvvvu0yeffKLXXnutyjmZmZkqLi72PAoLCwNYIQAACKaQWbGZOnWqVq1apY0bN6pp06ZVznO73XK73QGsDAAAOIXjg40xRlOnTtWKFSuUl5enVq1aBbskAADgUI4PNlOmTNGrr76qv/71r4qKitKRI0ckSTExMYqMjAxydQAAwEkc/xmbhQsXqri4WCkpKUpMTPQ83njjjWCXBgAAHMbxKzbGmGCXAAAAQoTjV2wAAABqimADAACsQbABAADWINgAAABrEGwAAIA1CDYAAMAaBBsAAGANgg0AALAGwQYAAFiDYAMAAKxBsAEAANYg2AAAAGsQbAAAgDUINgAAwBrhwS7AFi0zVgf0uLrgpFoA2C/U/8ypi/p9PefB3CF+riSwKuu7vOykT+dixQYAAFiDYAMAAKxBsAEAANYg2AAAAGsQbAAAgDUINgAAwBoEGwAAYA2CDQAAsAbBBgAAWINgAwAArEGwAQAA1iDYAAAAaxBsAACANQg2AADAGgQbAABgDYINAACwRkgEmwULFqhVq1Zq0KCBunfvrvfffz/YJQEAAAdyfLB54403lJaWppkzZyo/P1/XXXedUlNTdejQoWCXBgAAHMbxwebpp5/W3XffrQkTJqh9+/aaO3eumjVrpoULFwa7NAAA4DDhwS6gOqdOndKOHTuUkZHhNT5w4EBt3ry50mPKyspUVlbm2S4uLpYklZSU1Og5y8tOVrmvunNUdxwAhLpQ//Mv1Os/Hyf1V9O/b3+pshrPjRljanUuRweb77//XmfPnlV8fLzXeHx8vI4cOVLpMTk5OcrOzq4w3qxZswuuJ2buBZ8CAEJSqP/5F+r1n4+T+vN3LcePH1dMTEyN5zs62Jzjcrm8to0xFcbOyczMVHp6ume7vLxcP/74oxo1auQ5pqSkRM2aNVNhYaGio6PrrvAgsb0/yf4ebe9Psr9H2/uT7O/R9v4kZ/dojNHx48eVlJRUq+McHWwaN26ssLCwCqszRUVFFVZxznG73XK73V5jl112WaVzo6OjHXch/cn2/iT7e7S9P8n+Hm3vT7K/R9v7k5zbY21Was5x9IeHIyIi1L17d61bt85rfN26derdu3eQqgIAAE7l6BUbSUpPT9eYMWPUo0cP9erVS88//7wOHTqkyZMnB7s0AADgMGFZWVlZwS6iOp06dVKjRo00e/Zs/fGPf9RPP/2kl156ScnJyRd03rCwMKWkpCg83PHZzie29yfZ36Pt/Un292h7f5L9Pdren2Rfjy5T2+9RAQAAOJSjP2MDAABQGwQbAABgDYINAACwBsEGAABYI+SDzfHjx5WWlqYWLVooMjJSvXv31rZt26o95pVXXlFycrIuueQSJSYmaty4cfrhhx88+5cuXSqXy1Xh8fPPP9d1O5Xypcc//elPat++vSIjI9W2bVv97//+b4U5y5cvV4cOHeR2u9WhQwetWLGirlqoVl30F8xruHHjRg0bNkxJSUlyuVxauXKl135jjLKyspSUlKTIyEilpKTos88+O+95FyxYoFatWqlBgwbq3r273n//fa/9ZWVlmjp1qho3bqyGDRvq5ptv1tdff+3X3qTg9ZeSklLheo4aNcqvvZ1TFz2e75y+ntcXwervrrvuqnANr732Wr/2VtN6fOkxJydHPXv2VFRUlOLi4nTLLbdo7969XnNC+T6sSX+BvA99FfLBZsKECVq3bp1eeukl7d69WwMHDtSNN96ob775ptL5mzZt0p133qm7775bn332mf7yl79o27ZtmjBhgte86OhoHT582OvRoEGDQLRUQW17XLhwoTIzM5WVlaXPPvtM2dnZmjJlit5++23PnA8//FC33XabxowZo48//lhjxozRyJEjtXXr1kC15VEX/UnBu4YnTpxQcnKy5s+fX+n+J598Uk8//bTmz5+vbdu2KSEhQb/97W91/PjxKs/5xhtvKC0tTTNnzlR+fr6uu+46paam6tChQ545aWlpWrFihV5//XVt2rRJpaWlGjp0qM6ePWtFf5I0ceJEr+v53HPP+bW3c+qix/Od09fz+iJY/UnSTTfd5HUN33nnnQvqxdd6fOlxw4YNmjJlirZs2aJ169bpzJkzGjhwoE6cOOGZE8r3YU36kwJ3H/rMhLCTJ0+asLAw87e//c1rPDk52cycObPSY/7whz+Y1q1be409++yzpmnTpp7tJUuWmJiYGP8X7ANfeuzVq5f5/e9/7zV2//33mz59+ni2R44caW666SavOYMGDTKjRo3yU+U1U1f9OeUaSjIrVqzwbJeXl5uEhASTm5vrGfv5559NTEyMWbRoUZXnueaaa8zkyZO9xtq1a2cyMjKMMcYcO3bM1K9f37z++uue/d98842pV6+eWbNmjb/aqSBQ/RljTL9+/cz999/vx+prxl89VndOf53XF4Hqzxhjxo4da4YPH37hRddSXfRojDFFRUVGktmwYYMxJvTvw1/7dX/GBO8+rI2QXrE5c+aMzp49W+Ff4ZGRkdq0aVOlx/Tu3Vtff/213nnnHRlj9H//93968803NWTIEK95paWlatGihZo2baqhQ4cqPz+/zvqoji89lpWVVTr/o48+0unTpyX9e8Vm4MCBXnMGDRqkzZs3+7H686ur/iTnXMNfKigo0JEjR7xee7fbrX79+lX52p86dUo7duyocL0GDhzoOWbHjh06ffq015ykpCR16tQpoNe0rvo755VXXlHjxo3VsWNH/f73v/f7SkZN+NJjMM/rtDry8vIUFxenNm3aaOLEiSoqKrrgc9aWv3osLi6WJMXGxkoK7fuwMr/u7xwn3IfVCelgExUVpV69eumxxx7Tt99+q7Nnz+rll1/W1q1bdfjw4UqP6d27t1555RXddtttioiIUEJCgi677DLNmzfPM6ddu3ZaunSpVq1apddee00NGjRQnz599MUXXwTi49vkAAAMnklEQVSqNQ9fehw0aJBeeOEF7dixQ8YYbd++XS+++KJOnz6t77//XpJ05MiRCj9IND4+vsIPHK1rddWfk67hL517fWvz2n///fc6e/ZstcccOXJEERERuvzyy2t83rpQV/1J0ujRo/Xaa68pLy9Ps2bN0vLly3Xrrbf6uYPz86XHYJ7XSXWkpqbqlVde0XvvvaennnpK27Zt04ABA1RWVnZB560tf/RojFF6err69u2rTp06ec4bqvfhr1XWn+Sc+7A6If//J7/00ksaP368rrjiCoWFhalbt2664447tHPnzkrnf/7555o2bZoeeeQRDRo0SIcPH9b06dM1efJk/fnPf5YkXXvttV4faOvTp4+6deumefPm6dlnnw1IX79U2x5nzZqlI0eO6Nprr5UxRvHx8brrrrv05JNPKiwszDPP5XJ5HWeMqTAWCHXRn9Ou4a/58tr7ckywrmld9Ddx4kTPrzt16qSrrrpKPXr00M6dO9WtWzc/VF07dXX/OOW+rIs6brvtNs+vO3XqpB49eqhFixZavXp1UP5yvJAe77vvPn3yySdVriz7el5/qov+nHYfViakV2wk6corr9SGDRtUWlqqwsJCz9sRrVq1qnR+Tk6O+vTpo+nTp+vqq6/WoEGDtGDBAr344otVrhDUq1dPPXv2DNq/9mvbY2RkpF588UWdPHlSBw8e1KFDh9SyZUtFRUWpcePGkqSEhIQKyb2oqKhCwg+Euujv14J9Dc9JSEiQpFq99o0bN1ZYWFi1xyQkJOjUqVM6evRojc9bF+qqv8p069ZN9evXD/g19aXHYJ7XyXUkJiaqRYsWIXcNp06dqlWrVmn9+vVq2rSp13lD9T78par6q0yw7sPqhHywOadhw4ZKTEzU0aNH9fe//13Dhw+vdN7JkydVr5532+f+lW+q+LFZxhjt2rVLiYmJ/i26lmra4zn169dX06ZNFRYWptdff11Dhw719N6rVy+tW7fOa/7atWvVu3fvOqv/fPzZ36855Rq2atVKCQkJXq/9qVOntGHDhipf+4iICHXv3r3C9Vq3bp3nmO7du6t+/fpecw4fPqxPP/00oNe0rvqrzGeffabTp08H/Jr60mMwz+vkOn744QcVFhaGzDU0xui+++7TW2+9pffee6/CP75C+T6Uzt9fZYJ1H1YrcJ9Trhtr1qwx7777rvnyyy/N2rVrTXJysrnmmmvMqVOnjDHGZGRkmDFjxnjmL1myxISHh5sFCxaYAwcOmE2bNpkePXqYa665xjMnKyvLrFmzxhw4cMDk5+ebcePGmfDwcLN169aA92dM7Xvcu3eveemll8y+ffvM1q1bzW233WZiY2NNQUGBZ84HH3xgwsLCTG5urtmzZ4/Jzc014eHhZsuWLYFur076C+Y1PH78uMnPzzf5+flGknn66adNfn6++eqrr4wxxuTm5pqYmBjz1ltvmd27d5vbb7/dJCYmmpKSEs85BgwYYObNm+fZfv311039+vXNn//8Z/P555+btLQ007BhQ3Pw4EHPnMmTJ5umTZuaf/zjH2bnzp1mwIABJjk52Zw5cybk+9u/f7/Jzs4227ZtMwUFBWb16tWmXbt2pmvXrn7vr656PN85a3reUO3v+PHj5sEHHzSbN282BQUFZv369aZXr17miiuu8Ht/ddXjPffcY2JiYkxeXp45fPiw53Hy5EnPnFC+D8/XX6DvQ1+FfLB54403TOvWrU1ERIRJSEgwU6ZMMceOHfPsHzt2rOnXr5/XMc8++6zp0KGDiYyMNImJiWb06NHm66+/9uxPS0szzZs3NxEREaZJkyZm4MCBZvPmzYFqqYLa9vj555+bLl26mMjISBMdHW2GDx9u/vWvf1U471/+8hfTtm1bU79+fdOuXTuzfPnyQLRTQV30F8xruH79eiOpwmPs2LHGmH9/FfPRRx81CQkJxu12m+uvv97s3r3b6xwtWrQwjz76qNfYn/70J9OiRQsTERFhunXr5vUVTGOM+emnn8x9991nYmNjTWRkpBk6dKg5dOiQFf0dOnTIXH/99SY2NtZERESYK6+80kybNs388MMPfu+vrno83zlret5Q7e/kyZNm4MCBpkmTJqZ+/fqmefPmZuzYsXXye7SueqzsfJLMkiVLPHNC+T48X3+Bvg995TKmivdfAAAAQow1n7EBAAAg2AAAAGsQbAAAgDUINgAAwBoEGwAAYA2CDQAAsAbBBgAAWINgA8Bv8vLy5HK5dOzYMUnS0qVLddlll3n2Z2VlqUuXLp7tu+66S7fcckvA6/SnX/cMILhC/qd7A3CO3r176/Dhw4qJianR/GeeeabKn9EGAL4g2ADwm4iICM9PFq6JmgYgAKgp3ooCUKWUlBRNnTpVaWlpuvzyyxUfH6/nn39eJ06c0Lhx4xQVFaUrr7xS7777rqTavy3z67ei3nzzTXXu3FmRkZFq1KiRbrzxRp04ccKzf8mSJWrfvr0aNGigdu3aacGCBV7n+/rrrzVq1CjFxsaqYcOG6tGjh7Zu3erZv3DhQl155ZWKiIhQ27Zt9dJLL3kd73K59MILL+h3v/udLrnkEl111VVatWqV15x33nlHbdq0UWRkpPr376+DBw967f/qq680bNgwXX755WrYsKE6duyod955p0avB4ALR7ABUK1ly5apcePG+uijjzR16lTdc889GjFihHr37q2dO3dq0KBBGjNmjE6ePHlBz3P48GHdfvvtGj9+vPbs2aO8vDzdeuutnreqFi9erJkzZ+qJJ57Qnj17NHv2bM2aNUvLli2TJJWWlqpfv3769ttvtWrVKn388cd66KGHVF5eLklasWKF7r//fj344IP69NNPNWnSJI0bN07r16/3qiM7O1sjR47UJ598osGDB2v06NH68ccfJUmFhYW69dZbNXjwYO3atUsTJkxQRkaG1/FTpkxRWVmZNm7cqN27d2vOnDm69NJLL+i1AVALwf0ZnACcrF+/fqZv376e7TNnzpiGDRuaMWPGeMYOHz5sJJkPP/zQ8xOHjx49aowxZsmSJSYmJsYz99FHHzXJycme7bFjx5rhw4cbY4zZsWOHkWQOHjxYaS3NmjUzr776qtfYY489Znr16mWMMea5554zUVFRVf6k4d69e5uJEyd6jY0YMcIMHjzYsy3JPPzww57t0tJS43K5zLvvvmuMMSYzM9O0b9/elJeXe+bMmDHDq+fOnTubrKysSmsAUPdYsQFQrauvvtrz67CwMDVq1EidO3f2jMXHx0uSioqKLuh5kpOTdcMNN6hz584aMWKEFi9erKNHj0qSvvvuOxUWFuruu+/WpZde6nk8/vjjOnDggCRp165d6tq1q2JjYys9/549e9SnTx+vsT59+mjPnj1V9tuwYUNFRUV5etuzZ4+uvfZauVwuz5xevXp5HT9t2jQ9/vjj6tOnjx599FF98sknPr4iAHxBsAFQrfr163ttu1wur7Fzf8mfe8vHV2FhYVq3bp3effdddejQQfPmzVPbtm1VUFDgOffixYu1a9cuz+PTTz/Vli1bJEmRkZHnfY5fBhJJMsZUGKus33PPb2rwDa4JEyboyy+/1JgxY7R792716NFD8+bNO+9xAPyDYAPAMVwul/r06aPs7Gzl5+crIiJCK1asUHx8vK644gp9+eWX+s1vfuP1aNWqlaR/r7Ts2rXL83mYX2vfvr02bdrkNbZ582a1b9++xvV16NDBE6TO+fW2JDVr1kyTJ0/WW2+9pQcffFCLFy+u8XMAuDB83RuAI2zdulX//Oc/NXDgQMXFxWnr1q367rvvPMEjKytL06ZNU3R0tFJTU1VWVqbt27fr6NGjSk9P1+23367Zs2frlltuUU5OjhITE5Wfn6+kpCT16tVL06dP18iRI9WtWzfdcMMNevvtt/XWW2/pH//4R41rnDx5sp566imlp6dr0qRJ2rFjh5YuXeo1Jy0tTampqWrTpo2OHj2q9957r1bhCcCFYcUGgCNER0dr48aNGjx4sNq0aaOHH35YTz31lFJTUyX9+y2eF154QUuXLlXnzp3Vr18/LV261LNiExERobVr1youLk6DBw9W586dlZubq7CwMEnSLbfcomeeeUZ/+MMf1LFjRz333HNasmSJUlJSalxj8+bNtXz5cr399ttKTk7WokWLNHv2bK85Z8+e1ZQpU9S+fXvddNNNatu2bYWvpQOoOy5TkzeNAQAAQgArNgAAwBoEGwAAYA2CDQAAsAbBBgAAWINgAwAArEGwAQAA1iDYAAAAaxBsAACANQg2AADAGgQbAABgDYINAACwBsEGAABY4/8DOJv7MRBb154AAAAASUVORK5CYII=",
      "text/plain": [
       "PyPlot.Figure(PyObject <matplotlib.figure.Figure object at 0x1378d4150>)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "PyObject Text(24,0.5,u'count')"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t = c_bench.times / 1e6 # times in milliseconds\n",
    "m, σ = minimum(t), std(t)\n",
    "\n",
    "plt[:hist](t, bins=500)\n",
    "xlim(m - 0.01, m + σ),\n",
    "xlabel(\"milliseconds\")\n",
    "ylabel(\"count\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 2. Python's built in `sum` "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The `PyCall` package provides a Julia interface to Python:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": true,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "#Pkg.add(\"PyCall\")\n",
    "using PyCall"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "PyObject <built-in function sum>"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Call a low-level PyCall function to get a Python list, because\n",
    "# by default PyCall will convert to a NumPy array instead (we benchmark NumPy below):\n",
    "\n",
    "apy_list = PyCall.array2py(a, 1, 1)\n",
    "\n",
    "# get the Python built-in \"sum\" function:\n",
    "pysum = pybuiltin(\"sum\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5.000180786930093e6"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pysum(a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pysum(a) ≈ sum(a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "BenchmarkTools.Trial: \n",
       "  memory estimate:  512 bytes\n",
       "  allocs estimate:  17\n",
       "  --------------\n",
       "  minimum time:     101.583 ms (0.00% GC)\n",
       "  median time:      102.946 ms (0.00% GC)\n",
       "  mean time:        103.999 ms (0.00% GC)\n",
       "  maximum time:     112.713 ms (0.00% GC)\n",
       "  --------------\n",
       "  samples:          49\n",
       "  evals/sample:     1"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "py_list_bench = @benchmark $pysum($apy_list)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Dict{Any,Any} with 2 entries:\n",
       "  \"C\"               => 9.85597\n",
       "  \"Python built-in\" => 101.583"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "d[\"Python built-in\"] = minimum(py_list_bench.times) / 1e6\n",
    "d"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 3. Python: `numpy` \n",
    "\n",
    "## Takes advantage of hardware \"SIMD\", but only works when it works.\n",
    "\n",
    "`numpy` is an optimized C library, callable from Python.\n",
    "It may be installed within Julia as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "BenchmarkTools.Trial: \n",
       "  memory estimate:  720 bytes\n",
       "  allocs estimate:  22\n",
       "  --------------\n",
       "  minimum time:     6.388 ms (0.00% GC)\n",
       "  median time:      6.545 ms (0.00% GC)\n",
       "  mean time:        6.642 ms (0.00% GC)\n",
       "  maximum time:     9.712 ms (0.00% GC)\n",
       "  --------------\n",
       "  samples:          751\n",
       "  evals/sample:     1"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "numpy_sum = pyimport(\"numpy\")[\"sum\"]\n",
    "apy_numpy = PyObject(a) # converts to a numpy array by default\n",
    "\n",
    "py_numpy_bench = @benchmark $numpy_sum($apy_numpy)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5.000180786929854e6"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "numpy_sum(apy_list) # python thing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "numpy_sum(apy_list) ≈ sum(a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Dict{Any,Any} with 3 entries:\n",
       "  \"C\"               => 9.85597\n",
       "  \"Python numpy\"    => 6.38756\n",
       "  \"Python built-in\" => 101.583"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "d[\"Python numpy\"] = minimum(py_numpy_bench.times) / 1e6\n",
    "d"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 4. Python, hand-written "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "PyObject <function py_sum at 0x13b615938>"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "py\"\"\"\n",
    "def py_sum(a):\n",
    "    s = 0.0\n",
    "    for x in a:\n",
    "        s = s + x\n",
    "    return s\n",
    "\"\"\"\n",
    "\n",
    "sum_py = py\"py_sum\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "BenchmarkTools.Trial: \n",
       "  memory estimate:  512 bytes\n",
       "  allocs estimate:  17\n",
       "  --------------\n",
       "  minimum time:     1.446 s (0.00% GC)\n",
       "  median time:      1.449 s (0.00% GC)\n",
       "  mean time:        1.452 s (0.00% GC)\n",
       "  maximum time:     1.466 s (0.00% GC)\n",
       "  --------------\n",
       "  samples:          4\n",
       "  evals/sample:     1"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "py_hand = @benchmark $sum_py($apy_list)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5.000180786930093e6"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sum_py(apy_list)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "true"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sum_py(apy_list) ≈ sum(a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Dict{Any,Any} with 4 entries:\n",
       "  \"C\"                   => 9.85597\n",
       "  \"Python numpy\"        => 6.38756\n",
       "  \"Python hand-written\" => 1445.53\n",
       "  \"Python built-in\"     => 101.583"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "d[\"Python hand-written\"] = minimum(py_hand.times) / 1e6\n",
    "d"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 5. Julia (built-in) \n",
    "\n",
    "## Written directly in Julia, not in C!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "sum(a) at <a href=\"https://github.com/JuliaLang/julia/tree/903644385b91ed8d95e5e3a5716c089dd1f1b08a/base/reduce.jl#L359\" target=\"_blank\">reduce.jl:359</a>"
      ],
      "text/plain": [
       "sum(a) in Base at reduce.jl:359"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "@which sum(a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "BenchmarkTools.Trial: \n",
       "  memory estimate:  0 bytes\n",
       "  allocs estimate:  0\n",
       "  --------------\n",
       "  minimum time:     6.441 ms (0.00% GC)\n",
       "  median time:      6.638 ms (0.00% GC)\n",
       "  mean time:        6.884 ms (0.00% GC)\n",
       "  maximum time:     10.581 ms (0.00% GC)\n",
       "  --------------\n",
       "  samples:          725\n",
       "  evals/sample:     1"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "j_bench = @benchmark sum($a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Dict{Any,Any} with 5 entries:\n",
       "  \"C\"                   => 9.85597\n",
       "  \"Python numpy\"        => 6.38756\n",
       "  \"Python hand-written\" => 1445.53\n",
       "  \"Python built-in\"     => 101.583\n",
       "  \"Julia built-in\"      => 6.44099"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "d[\"Julia built-in\"] = minimum(j_bench.times) / 1e6\n",
    "d"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 6. Julia (hand-written) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "mysum (generic function with 1 method)"
      ]
     },
     "execution_count": 34,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "function mysum(A)   \n",
    "    s = 0.0  # s = zero(eltype(A))\n",
    "    for a in A   \n",
    "        s += a   \n",
    "    end\n",
    "    s\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "BenchmarkTools.Trial: \n",
       "  memory estimate:  0 bytes\n",
       "  allocs estimate:  0\n",
       "  --------------\n",
       "  minimum time:     10.150 ms (0.00% GC)\n",
       "  median time:      10.365 ms (0.00% GC)\n",
       "  mean time:        10.546 ms (0.00% GC)\n",
       "  maximum time:     13.060 ms (0.00% GC)\n",
       "  --------------\n",
       "  samples:          474\n",
       "  evals/sample:     1"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "j_bench_hand = @benchmark mysum($a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Dict{Any,Any} with 6 entries:\n",
       "  \"C\"                   => 9.85597\n",
       "  \"Python numpy\"        => 6.38756\n",
       "  \"Julia hand-written\"  => 10.1501\n",
       "  \"Python hand-written\" => 1445.53\n",
       "  \"Python built-in\"     => 101.583\n",
       "  \"Julia built-in\"      => 6.44099"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "d[\"Julia hand-written\"] = minimum(j_bench_hand.times) / 1e6\n",
    "d"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Summary"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Python numpy..............6.39\n",
      "Julia built-in............6.44\n",
      "C.........................9.86\n",
      "Julia hand-written.......10.15\n",
      "Python built-in.........101.58\n",
      "Python hand-written....1445.53\n"
     ]
    }
   ],
   "source": [
    "for (key, value) in sort(collect(d), by=x->x[2])\n",
    "    println(rpad(key, 20, \".\"), lpad(round(value, 2), 10, \".\"))\n",
    "end"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Julia 0.6.0",
   "language": "julia",
   "name": "julia-0.6"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "0.6.0"
  },
  "toc": {
   "colors": {
    "hover_highlight": "#DAA520",
    "running_highlight": "#FF0000",
    "selected_highlight": "#FFD700"
   },
   "moveMenuLeft": true,
   "nav_menu": {
    "height": "212px",
    "width": "252px"
   },
   "navigate_menu": true,
   "number_sections": true,
   "sideBar": true,
   "threshold": "2",
   "toc_cell": false,
   "toc_section_display": "block",
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}