hessian.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import tensorflow as tf\n",
    "import numpy as np\n",
    "from time import time\n",
    "from matplotlib import pyplot as plt\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "from tensorflow.python.ops.gradients_impl import _AsList\n",
    "from tensorflow.python.ops.gradients_impl import *\n",
    "from tensorflow.python.ops import check_ops\n",
    "from tensorflow.python.ops import tensor_array_ops\n",
    "\n",
    "\n",
    "def hessians(ys, xs, name=\"hessians\", colocate_gradients_with_ops=False,\n",
    "            gate_gradients=False, aggregation_method=None):\n",
    "  \"\"\"Constructs the Hessian of sum of `ys` with respect to `x` in `xs`.\n",
    "\n",
    "  `hessians()` adds ops to the graph to output the Hessian matrix of `ys`\n",
    "  with respect to `xs`.  It returns a list of `Tensor` of length `len(xs)`\n",
    "  where each tensor is the Hessian of `sum(ys)`. This function currently\n",
    "  only supports evaluating the Hessian with respect to (a list of) one-\n",
    "  dimensional tensors.\n",
    "\n",
    "  The Hessian is a matrix of second-order partial derivatives of a scalar\n",
    "  tensor (see https://en.wikipedia.org/wiki/Hessian_matrix for more details).\n",
    "\n",
    "  Args:\n",
    "    ys: A `Tensor` or list of tensors to be differentiated.\n",
    "    xs: A `Tensor` or list of tensors to be used for differentiation.\n",
    "    name: Optional name to use for grouping all the gradient ops together.\n",
    "      defaults to 'hessians'.\n",
    "    colocate_gradients_with_ops: See `gradients()` documentation for details.\n",
    "    gate_gradients: See `gradients()` documentation for details.\n",
    "    aggregation_method: See `gradients()` documentation for details.\n",
    "\n",
    "  Returns:\n",
    "    A list of Hessian matrices of `sum(y)` for each `x` in `xs`.\n",
    "\n",
    "  Raises:\n",
    "    LookupError: if one of the operations between `xs` and `ys` does not\n",
    "      have a registered gradient function.\n",
    "    ValueError: if the arguments are invalid or not supported. Currently,\n",
    "      this function only supports one-dimensional `x` in `xs`.\n",
    "  \"\"\"\n",
    "  xs = _AsList(xs)\n",
    "  kwargs = {\n",
    "      'colocate_gradients_with_ops': colocate_gradients_with_ops,\n",
    "      'gate_gradients': gate_gradients,\n",
    "      'aggregation_method': aggregation_method\n",
    "    }\n",
    "  # Compute a hessian matrix for each x in xs\n",
    "  hessians = []\n",
    "  _gradients = gradients(ys, xs, **kwargs)\n",
    "  for i, _gradient, x in zip(range(len(xs)), _gradients, xs):\n",
    "    check_rank = check_ops.assert_rank(x, 1, message='Cannot compute Hessian because '\n",
    "                                       'element %d of `xs` does not have rank one.' % i)\n",
    "    with ops.control_dependencies([check_rank]):\n",
    "      # Declare an iterator and tensor array loop variables to store the gradients\n",
    "      n = array_ops.size(x)\n",
    "      loop_vars = [array_ops.constant(0, dtypes.int32), tensor_array_ops.TensorArray(x.dtype, n)]\n",
    "      # Iterate over all elements of the gradient\n",
    "      _, hessian = control_flow_ops.while_loop(\n",
    "          lambda j, _: j < n,\n",
    "          lambda j, result: (j + 1, result.write(j, gradients(_gradient[j], x)[0])),\n",
    "          loop_vars\n",
    "      )\n",
    "        \n",
    "      hessians.append(hessian.stack())\n",
    "  return hessians"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2\n",
      "old .......... B 0.420970+-0.034222; E 0.137379+-0.011473\n",
      "new .......... B 0.380147+-0.064631; E 0.149488+-0.033953\n",
      "4\n",
      "old .......... B 0.954304+-0.142474; E 0.307502+-0.080685\n",
      "new .......... B 0.365630+-0.058224; E 0.121651+-0.005977\n",
      "8\n",
      "old .......... B 1.769587+-0.160420; E 0.570243+-0.072275\n",
      "new .......... B 0.348446+-0.051654; E 0.128656+-0.009134\n",
      "16\n",
      "old .......... B 3.444399+-0.189581; E 1.191949+-0.084632\n",
      "new .......... B 0.356584+-0.047020; E 0.132410+-0.004239\n",
      "32\n",
      "old .......... B 7.618027+-0.556610; E 2.424034+-0.265467\n",
      "new .......... B 0.406391+-0.152170; E 0.150345+-0.011243\n",
      "64\n",
      "old .......... B 17.382572+-1.143090; E 5.475313+-0.401230\n",
      "new .......... B 0.371568+-0.140755; E 0.153941+-0.010728\n"
     ]
    }
   ],
   "source": [
    "sizes = 2 ** np.arange(1, 7)\n",
    "num_runs = 10\n",
    "\n",
    "modes = {'old': tf.hessians, 'new': hessians}\n",
    "\n",
    "build_times = {mode: [] for mode in modes}\n",
    "eval_times = {mode: [] for mode in modes}\n",
    "\n",
    "\n",
    "for size in sizes:\n",
    "    print(size)\n",
    "    _x = np.random.normal(0, 1, size)\n",
    "    _A = np.random.normal(0, 1, (size, size))\n",
    "    \n",
    "    \n",
    "    for mode, method in modes.items():\n",
    "        print(mode, end=' ')\n",
    "        \n",
    "        eval_row = []\n",
    "        build_row = []\n",
    "        for run in range(num_runs):\n",
    "            print('.', end='', flush=True)\n",
    "            \n",
    "            with tf.Graph().as_default():\n",
    "                x = tf.constant(_x, dtype=np.float32)\n",
    "                A = tf.constant(_A, dtype=np.float32)\n",
    "                y = tf.reduce_sum(A * x * x[:, None])\n",
    "\n",
    "                start = time()\n",
    "                hess = method(y, [x])\n",
    "                build_row.append(time() - start)\n",
    "                session = tf.Session()\n",
    "\n",
    "            start = time()\n",
    "            value, = session.run(hess)\n",
    "            np.testing.assert_allclose(value, _A + _A.T, 1e-5, 1e-5)\n",
    "            eval_row.append(time() - start)\n",
    "        print(\" B %f+-%f; E %f+-%f\" % (np.mean(build_row), np.std(build_row),\n",
    "                                      np.mean(eval_row), np.std(eval_row)))\n",
    "            \n",
    "        build_times[mode].append(build_row)\n",
    "        eval_times[mode].append(eval_row)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xd4VGX2wPHvSSOU0Jv03juhuWIvWLCjoruK3XVdC7oq\n6+5adl137fUnNsRVBKkiir2hIAKhQ+gECIQWQhJCes7vjzuBgCGZJDO5U87neeaZmXvv3HtS3jn3\nLfe9oqoYY4wxgSbC7QCMMcaY0liCMsYYE5AsQRljjAlIlqCMMcYEJEtQxhhjApIlKGOMMQHJEpQx\nJuSIyA8icouf9v1XEXnbH/su5VjjReTv1XGsQGQJyhjjGhFJEpFsETlU4vGq23EVE5HTRSS55DJV\n/beq+jz5icgYEfn5uGPdoar/9PWxgkWU2wEYY8LeSFX9xu0gTOCxGpQxJuCISA0ROSgivUosa+Kp\nbTUVkQYi8qmI7BORNM/rVifY12Mi8kGJ9+1EREUkyvP+RhFJFJFMEdkiIrd7ltcGPgdalKjdtShl\nfxeLyBpPvD+ISPcS65JE5AERWSki6SLykYjElhJjd2A8MMxznIOe5RNF5F+e16eLSLKIPCgie0Uk\nRUQuFZELRGSDiBwQkb+W2GeEiDwsIptFJFVEpopIQ8+6WBH5wLP8oIgsFpFmlftr+Y8lKGNMwFHV\nXGAmMLrE4quAH1V1L85317tAW6ANkA1UtmlwL3ARUBe4EXhBRAaoahZwPrBLVet4HrtKflBEugCT\ngXuBJsBcYI6IxBwX9wigPdAHGFPKz5sI3AH84jlO/RPE2hyIBVoC/wDeAn4PDASGA/8QkQ6ebe8G\nLgVOA1oAacBrnnU3APWA1kAjz7Gzy/41VT9LUMYYt33sOYsvftzqWf4hxyaoaz3LUNVUVZ2hqodV\nNRN4EueLuMJU9TNV3ayOH4GvcL7svXE18Jmqfq2q+cCzQE3g5BLbvKyqu1T1ADAH6FeZOD3ygSc9\nx5oCNAZeUtVMVV0DrMFJggC3A4+oarIn4T8GXOmpOebjJKZOqlqoqgmqmlGFuPzC+qCMMW679AR9\nUN8BNUVkCLAb54t9FoCI1AJewKmZNPBsHycikapaWJGDi8j5wKNAF5yT9lrAKi8/3gLYVvxGVYtE\nZAdODafY7hKvD3s+U1mpJX6+4hrPnhLrs4E6ntdtgVkiUlRifSHQDHgfp/Y0RUTqAx/gJLP8KsTm\nc1aDMsYEJFUtAqbi1KKuBT711JYA7ge6AkNUtS5wqme5lLKrLJykU6x58QsRqQHMwKn5NPM0rc0t\nsZ/ybvewCycRFO9PcL74d5b385XC17eW2AGcr6r1SzxiVXWnquar6uOq2gOntncRcL2Pj19llqCM\nMYHsQ5xmtOs8r4vF4dQWDno6/h8tYx/LgVNFpI2I1APGlVgXA9QA9gEFntrUuSXW7wEaeT5XmqnA\nhSJylohE4yTOXGCBtz/gccdqdVz/VVWMB54UkbZwZJDJJZ7XZ4hIbxGJBDJwmvwqVPOsDpagjDFu\nm3PcdVCzileo6q84NaAWOCPqir2I09ezH1gIfHGinavq18BHwEogAfi0xLpMnMEEU3EGEVwLfFJi\n/TqcQRBbPP1jxzTPqep6nEEKr3hiGYkzbD6vor8EnCbNNcBuEdlfic8f7yWcn+UrEcnE+T0N8axr\nDkzHSU6JwI84zXwBReyGhcYYYwKR1aCMMcYEJEtQxhhjApIlKGOMMQHJEpQxxpiAFNQX6jZu3Fjb\ntWvndhgmxCUkJOxX1SZux1EdrEyZ6uBtmQrqBNWuXTuWLFnidhgmxInItvK3Cg1Wpkx18LZMWROf\nMcaYgGQJyhhjTEAKygQlIiNF5M309HS3QzHGGOMnQdkHpapzgDnx8fG3Hr8uPz+f5ORkcnJyXIjM\n/2JjY2nVqhXR0dFuh2LChJUp45agTFBlSU5OJi4ujnbt2uFMLBw6VJXU1FSSk5Np37692+GYMGFl\nyrglKJv4ypKTk0OjRo1CriABiAiNGjUK2TNZ456yms2tTBm3hFyCAkKyIBUL5Z/NuEdV56jqbfXq\nlX5XiVD+vwvlny3YhWSCqqir3/iFq9/4xe0wjEsStqXx2vebSNiW5nYoIcPKVHjzVZmyBFWN6tSp\nU+ryMWPGMH369GqOxgAkJB3gitcX8MyX67nu7YWWpIKMlanAk7AtjSvH+6ZMWYIyYe3Nn7YceZ1f\nUMTCLakuRmNM8Ju/aT/FtxmsapmyBAVk5uSz82C2T8+en3/+eXr16kWvXr148cUXj1mnqtx11130\n6NGDCy+8kL179/rsuMZ721Kz+GH9PiIEIgWioyIY2qGR22GFBCtT4Wv7gawjr6tapkJumHlJj89Z\nw9pdGWVuk5mTz9qUTABGjV9At+ZxxMWe+HqIHi3q8ujInmXuMyEhgXfffZdff/0VVWXIkCGcdtpp\nR9bPmjWL9evXs2rVKvbs2UOPHj246aabKvCTmaoqLFLun7qCmKgIXh7dn017DzG0QyMGtm3gdmgB\nzcqUKcvK5IPMWraL07s0YVD7hlUuUyGdoLyRkVNw5HWROu/LKkze+Pnnn7nsssuoXbs2AJdffjk/\n/fTTkfXz5s1j9OjRREZG0qJFC84888wqHc9U3BvzNrNkWxovXN2X83o257yyvx9NBViZCk+5BYXc\nP3UFTerU4KXR/alXs+oXPod0girvrAycDr1R4xdQpBAbHcFL1/Sv8lm0FjfAlsGGtrpnza50Xvh6\nAxf0bs6l/Vq6HU5QsTJlTuSFrzeyce8hJt44yCfJCawPioFtG9CteRytGtRk0i1DfdLEc+qpp/Lx\nxx9z+PBhsrKymDVrFsOHDz9m/ZQpUygsLCQlJYXvv/++ysc03snJL2TsRyuoXyuGJy/tbV9qHr6c\n39LKVPhZuj2NN+dt5ppBrTm9a1Of7Teka1DeiouNJi422mf9DwMGDGDMmDEMHjwYgFtuuYX+/fsf\nWX/ZZZfx3Xff0bt3b7p06XJMW7rxr+e/3sD6PZm8e+MgGtSOcTucgFHW/JaVYWUqfOTkF/LAtBU0\nrxvLIxd29+m+LUH5ydixYxk7duwxyw4dOgQ4TRGvvvqqG2GFtYVbUnnrpy1cO6QNZ/jwLM9UDytT\ngem5r9azZV8WH9w8pMp9jcezBAV8dPswt0MwfpaZk8/9U1fQpmEtHrnAt2d55resTIWHxUkHePvn\nrfx+aBtO6dzY5/u3BGXCwhNz1pKSns20O06mdg37tzemqg7nFfCXaStoWb8m4873z0lf2A+SMKHv\nqzW7mZaQzB9P72jXORnjI09/sZ6k1MM8c2Vfv530WYIyIW3/oVzGzVxFzxZ1ueesLm6HY0xIWLgl\nlYkLkhhzcjuGdfTf7CvW1mFClqry8IxVZOYWMPnqfsRE2fmYMVWVlVvAX6avoF2jWjw4oqtfj2Ul\nFuDdC52HCSnTliTzTeIeHjyvK12axbkdTnixMhWynvo8keS0bJ4Z1ZdaMf6t41iCMiFpx4HDPD5n\nDUM7NOSm39mtvI3xhZ837ueDhdu5+XftGdSuod+PZwnKhJziiWAjRHh2VF8iImy2CGOqKjMnn4dm\nrKRDk9o8cJ5/m/aKWYICyM2A9B2wY5FPdpeUlET37t259dZb6dmzJ+eeey7Z2dls3ryZESNGMHDg\nQIYPH866desoLCykQ4cOqCoHDx4kIiKCefPmATB8+HA2bdrkk5jCyds/bWFR0gEevbgnrRrUcjuc\noODLqY4AK1Mh6MnPEklJz+bZUX2JjY6slmOG9iCJzx+G3avK3iY3A3avdF5POA+a9YIadU+8ffPe\ncP5/yj30xo0bmTx5Mm+99RZXXXUVM2bM4N1332X8+PF07tyZX3/9lTvvvJPvvvuOLl26sHbtWrZu\n3crAgQP56aefGDJkCMnJyXTq1KkCP7BJTMngua82cF7PZlwxwCaC9ZbXUx1ZmQpLP6zfy5TFO7jj\ntI4MaFN9l2qEdoLyRk6JM0Ytct6XVZi81L59e/r16wfAwIEDSUpKYsGCBYwaNerINrm5uYBzVjdv\n3jy2bt3KuHHjeOuttzjttNMYNGhQleMIJ7kFhdz30XLq1ozm35fZRLCusTIVUtKz83l4xiq6NKvD\nfed0rtZjh3aC8uKsjB2LnLM8LYKomnDF29B6cJUPXaNGjSOvIyMj2bNnD/Xr12f58uW/2Xb48OGM\nHz+eXbt28cQTT/DMM8/www8/cOqpp1Y5jnDywtcbWbc7k3duiKdRnRrlf8BUnJWpsPPEnLXsO5TL\nm9cPpEZU9TTtFbM+qNaDnSaI+m3hhk98UpBKU7duXdq3b8+0adMA5xqdFStWADBkyBAWLFhAREQE\nsbGx9OvXjzfeeOOY2wmYsi1OOsAbnun+z+rezO1wwpuVqZDxzdo9zFiazJ2nd6RPq/rVfnxLUOA0\nP9Rr7beCVGzSpEm888479O3bl549ezJ79mzn8DVq0Lp1a4YOHQo4Z3+ZmZn07t3br/GEikO5BYyd\nupzWDWrxt4t6uB2OAStTIeDg4TzGzVpFt+Zx/PnM6m3aKybe3KkyUMXHx+uSJUuOWZaYmEj37hWc\nuLD4gsIbP/NRZP5VqZ8xhD08YyVTl+xg6u3DiPfDtRkikqCq8T7fcQCyMmWK3TtlGZ+uTGH2Xb+j\nZ4t6Pt23t2UqtPugvBUkhcj81jdr9zBl8Q7+eHpHvyQnU0lWpoLaF6t38/HyXdx3dhefJ6eKCJgm\nPhHpICLviMh0t2MxwSH1UC4Pz1xJ95Pqct/ZNhGsMb6QeiiXR2Y5EyzfeUZHV2Pxa4ISkQkisldE\nVh+3fISIrBeRTSLyMICqblHVm31x3GButixPKP9sFaGq/HXWKjKyC3jRJoL1u1D+vwvln60y/vHJ\nGjJy8nnuqr5ER7pbrvx99InAiJILRCQSeA04H+gBjBYRn/Vsx8bGkpqaGpL/dKpKamoqsbGxbofi\nuhlLd/Llmj08cF4Xuja3iWD9ycpU+Ph05S4+W5nCvWd3oVvzql+7VlV+7YNS1Xki0u64xYOBTaq6\nBUBEpgCXAGu92aeI3AbcBtCmTZvfrG/VqhXJycns27ev8oEHsNjYWFq1auV2GK5KTjvMY5+sYXD7\nhtx8Sge3wwl5VqbCw77MXP7+8Wr6tqrH7acGRrlyY5BES2BHiffJwBARaQQ8CfQXkXGq+lRpH1bV\nN4E3wRlxdPz66Oho2re32atDVVGR8sA051qX50b1JdImgvU7K1OhT1X528eryMor5Lmr+hLlctNe\nMTcSVGnfKKqqqcAd1R2MCS4T5m9l4ZYDPH1lH1o3tIlgjfGFT1bs4ss1e/jrBd3o1DRwmszdSJPJ\nQOsS71sBu1yIwwSZDXsyefrL9ZzToxmjBlqTjDG+sCcjh3/MXsOANvUDrsncjQS1GOgsIu1FJAa4\nBvjEhThMEMkrKOLeKcupGxvFU5fbRLC+5vPbbZigoKr8deYqcvILeTYAm8z9Pcx8MvAL0FVEkkXk\nZlUtAO4CvgQSgamquqaC+7XCFGZe/GYDa1MyeOryPjS2iWB9TlXnqOpt9eq5d1GmqX4zlu7k23V7\neXBENzo0qeN2OL/h71F8o0+wfC4wtwr79e7eNSYkJGw7wPgfN3NVfCvO6WETwRrjCynp2Tw+Zw2D\n2zXkxpPbuR1OqQJjqIYxJ5CVW8B9H62gRf2a/N0mgjXGJ1SVh2asoqBQeWZUHyICrGmvmCUoE9D+\n9VkiO9IO8/xV/YiLjXY7HGNCwkeLdzBvwz7GXdCNto1qux3OCVmCMgHru3V7mLxoO7cN78Dg9jYR\nrDG+kJx2mH99lsjJHRvx+yFt3Q6nTJagTEA6kJXHg9Ode9GMPdcmgjXGF4qKlIdmrERV+e8Vgdu0\nVywoE5SN4gttqsojs1aRnp3H81f1q/bbTBsTqiYt2s78Tak8cmGPoLjQPSgTlA2JDW2/+893fL56\nN1fHt6FHC/cnrDQmFGxPPcxTcxMZ3rkxowe3Lv8DAcBuWGgCyuerUtiVngPA9KU7uGxASwa2beBy\nVMYEt/NfnMfm/VlEifDfK/oEzYXuQVmDMqGpoLCIf312dFL7/IIiFm5JdTEiY4JfwrY01u3OJK+g\niPyiIlI8J4DBwBKUCRivfr+JnQdziI4UIgWioyIY2qGR22EZE9RmLk2m+LYPRUUaVCd91sRnAsKS\npAO8/O1GLu/fkuuGtmXhllSGdmhkzXvGVMH+Q7l8vmr3kffBdtJnCcq4Lj07n3umLKdVg1o8fklP\n4mKjLTEZU0WFRcrdk5eRlVfAs6P6sCcjN+hO+oIyQYnISGBkp06d3A7FVFHxkPLdGTlMv2OYzRZh\njI88//V6FmxO5Zkr+3DlwOAYtXe8oOyDsmHmoWPG0p18ujKFsed0oX+b4DmzMyaQfbN2D699v5lr\nBrVmVHxwJicI0gRlQkPS/iz+MXs1Q9o35I7TOrodjjEhYXvqYcZOXU6vlnV57OKebodTJZagjCvy\nCoq4e8oyoiMjeOHqfgF3ozRjglFOfiF/nJQAwOvXDSQ2OrhnYQnKPigT/F74ZgMrk9N5/boBtKhf\n0+1wjAkJj85ew5pdGUwYEx8UUxmVx2pQptot2Lyf8T9uZvTg1pzf+yS3wzEmJExdvIOPluzgz2d2\n4sxuoXFjT0tQplqlZeUx9qMVtG9c225AaIyPrN6Zzt9nr+aUTo259+zQmf0/KBOUzWYenFSVh2eu\nJDUrl5ev6U+tGGthDhRWpoJX+uF87py0lIa1Y3jpmtDqzw3KBGXDzIPT5EU7+HLNHh4a0Y1eLe1v\nF0isTAWnoiLl/mnLSUnP5rXrBtCoTg23Q/KpoExQJvhs2pvJE5+uYXjnxtz0u/Zuh2NMSHj9x818\nk7iXv13YgwEheB2hJSjjd7kFhfx58nJqx0Tx3FV9A/4unsYEg/mb9vPcV+sZ2bcF1w8L7Fu3V5Z1\nAhi/e/qL9SSmOENfm8bFuh2OMUFvd3oOd09eRocmdfjP5b2D5v5OFWU1KONXP6zfyzs/b+WGYW1D\nZuirMW7KKyjizkkJ5OQXMv73A6ldI3TrGaH7kxnX7cvM5YFpK+jaLI5xF3R3OxxjQsJTnyeydPtB\nXr22P52a1nE7HL+yBGX8QlX5y/QVZOQUMOmWoUE/5YoxgWDOil28Oz+Jm37Xnov6tHA7HL+zJj7j\nFxMXJPHD+n387cLudG0e53Y4xgS9TXszeWjGSga2bcC4C7q5HU61CMoEZRcVBrbElAyemruOs7o1\n5Q9DQ3N0kTHVKSu3gDs+WErN6Eheu3YA0ZFB+dVdYUH5U9pFhYErO6+Quycvo16taJ6+sk/Iji4y\npro4M7CsYsu+Q7wyuj/N64XPSNgy+6BE5BMv9nFAVcf4JhwT7J6cu5aNew/x/s2DQ+6qdmPc8N6C\nJOas2MWDI7pycqfGbodTrcobJNEduKWM9QK85rtwTDD7as1uPli4ndtO7cDwzk3cDseYoJewLY1/\nfZbI2d2bcsep4XdTz/IS1COq+mNZG4jI4z6MxwSp3ek5PDhjJb1a1uWBc7u6HU5AEJGxXmyWpapv\n+D0YE3T2H8rlT5OW0qJ+TZ4b1S8sZ2Apsw9KVacev0xEIkSkblnbmPBSVKSMnbqc3PwiXrqmPzFR\nQdm16Q9/AeoAcWU87nctOhOwCouUe6YsI+1wHq//fgD1akW7HZIrvLoOSkQ+BO4ACoEEoJ6IPK+q\nz/gzOBMc3vxpCws2p/LfK3rTsUloXzhYQe+r6hNlbSAitasrGBM8Xvh6A/M3pfL0lX3o2SJ8B4N5\ne6rbQ1UzgEuBuUAb4A9+i8oEjZXJB3n2y/Vc0Ls5V8W3djucgKKqD/piGxNevk3cw6vfb+KaQa3D\nvkx5m6CiRSQaJ0HNVtV8QP0XlgkGWbkF3DNlOU3javDUZTak/ERE5B4RqSuOd0RkqYic63ZcJvBs\nTz3MfR8tp2eLujx2cU+3w3GdtwnqDSAJqA3ME5G2QIa/gjLB4fE5a0hKzeKFq/uFbRu5l27ytECc\nCzQBbgT+425IJtDk5Bfyx0kJALx+3UCbHgwvE5SqvqyqLVX1AlVVYDtwhn9DM4Hs05W7mLokmbvO\n6MSQDo3cDifQFVctLwDeVdUVJZYZA8Bjn6xhza4MXri6H20a1XI7nIBQZoISkYtKW66OgrK28Seb\n6shdyWmHGTdzFf3b1Ofuszq7HU4wSBCRr3AS1JciEgcUuRyTCSBTF+9gyuId3HVGJ87qbrelKVbe\nKL5nRGQnZZ/t/Rv41HchlU9V5wBz4uPjb63O4xpn+Ot9Hy1HFV66un/YzAlWRTcD/YAtqnpYRBrh\nNPMZw+qd6fx99mp+16kR953Txe1wAkp5CWoP8Hw522z0USwmCLz2/SYWJ6XxojVDlEtEmqvqblUt\nApYWL1fVVCC15DZuxWjclZ6dz52TltKgVgwvXdOfyDC8GLcsZSYoVT29muIwQSBh2wFe+nYjl/Vv\nyaX9W7odTjCYCwzwwTYmBBUVKfdPXc6ug9l8dPswGtvclb9hNyw0XsnIyeeeKctpUT+WJy6x4a9e\n6isiZY12FWw0bFhK2JbGS99sYN7G/Tw2sgcD2zZwO6SAZAnKlEtV+dus1aSk5zDtjmHExdqQcm+o\nqo0TNr+RsC2NK19fgAIRAr1bhu9MEeWxHm5Trhe+3sAnK3ZxdXwrBrSxMz1jquLbxD1HZjkQYOHW\nA26GE9C8SlAiMsozNBYR+ZuIzBQRazcPA9OW7ODl7zYBMHPZThK2pbkckTHBq6hI+XnjfsCpPUVH\nRTDUriM8IW+b+P6uqtNE5BTgPOBZ4HVgiN8iM65LSc/m8TlrjrzPLyhi4ZZUay83ppJe/3EzK3em\nc8dpHYiLjWZoh0ZWnsrgbYIq9DxfCLyuqrNF5DH/hGQCQVZuATdPXEJhEdSIiqCgsMjO9ipJRDoC\nyaqaKyKnA32A/6nqQXcjM9Vp0dYDPPfVei7u24KHRnSzuSu94G2C2ikibwBnA/8VkRpY/1XIKixS\n7p68jPV7MnnnhnjiYqNZuCXVzvYqbwYQLyKdgHeAT4APcWaWMGHgQFYed09eRpuGtXjysl6WnLzk\nbYK6ChgBPKuqB0XkJJybsZkQ9M9P1/Ltur3889JenN61KYAlpqopUtUCEbkMeFFVXxGRZW4HZapH\n8fVOB7LymHnnyTYKtgK8nSz2MLAXOMWzqACbQSIkTZy/lYkLkrjllPb8YWhbt8MJFfkiMhq4gaPT\nglXLt5SI1BaR90TkLRG5rjqOaY719s9b+H79Pv52UXd62ZDyCvF2FN+jwEPAOM+iaOADfwVl3PHd\nuj088elazunRjHEXdHc7nFByIzAMeFJVt4pIe6pQfkRkgojsFZHVxy0fISLrRWSTiDzsWXw5MF1V\nbwUuruwxTeUs3Z7G01+sZ0TP5nbCVwne9iNdhvPPnQWgqruAOH8FVR6bzdz31uxK564Pl9GzRT1e\nuqafzQnmQ6q6FucEb6nn/VZVrcr9oCbiNLkfISKRwGvA+UAPYLSI9ABaATs8mxViqk364Xz+/OEy\nmteL5b9X2g09K8PbBJXnuQ+UgtNs4L+Qyqeqc1T1tnr1rLrsC7vTc7hp4mLq1Yzm7RviqRVjE4z4\nkoiMBJYDX3je9xORTyq7P1WdBxx/dedgYJOqblHVPGAKcAmQjJOk4ATlXURuE5ElIrJk3759lQ3L\nlKCq/GX6CvZk5PDqtQOoV9P6nSrD2wQ11TOKr76I3Ap8A7zlv7BMdcnKLeDm9xZzKKeACWMG0axu\nrNshhaLHcBLIQQBVXQ609/ExWnK0pgROYmoJzASuEJHXgTmlfVBV31TVeFWNb9KkiY/DCk/vLUji\nq7V7ePj8bvRrXd/tcIKWV6fKqvqsiJyDM7FlV+Afqvq1XyMzfldYpNwzZRmJKRm8M2YQ3U+q63ZI\noapAVdOPa+LRE21cSaW1H6mqZmH3nqpWq5LT+ffcdZzVrSk3n+Lr85Dw4nVbjqp+LSK/Fn9GRBqq\nqk0iFcT+9dlavkncyz8v6ckZnuHkxi9Wi8i1QKSIdAbuBhb4+BjJQOsS71sBu3x8DFOOjJx8/vTh\nUhrXieHZUX2t36mKvB3Fd7uI7AFWAkuABM+zCVLvLUji3flJ3PS79vxhWDu3wwl1fwZ6ArnAZJyW\niHt9fIzFQGcRaS8iMcA1OBcEm2qiqoybuYqdB7N5eXR/GtSOcTukoOdtDeoBoKeq7vdnMKZ6fLdu\nD4/PWcPZ3ZvxyIU2nNzfPNcRPuJ5VJmITAZOBxqLSDLwqKq+IyJ3AV8CkcAEVV1Txm6O3+dIYGSn\nTp18EWJY+nDRdj5bmcKDI7oS366h2+GEBG8T1GbgsD8DMdVjza50/vzhMnq0qMvLo204eXUQkXjg\nr0A7SpQ5Ve1Tmf2p6ugTLJ+Lc4feyuxzDjAnPj7+1sp8PtwlpmTw+Jy1nNqlCXec2tHtcEKGtwlq\nHLDA0weVW7xQVe/2S1TGL3an53DzxCXUrRnNOzcMsuHk1WcSztRgq4Ail2MxPpaVW8CfPlxK/ZrR\nPH9VXyLspM9nvP2GegP4DitgQat4OHlmTj7T7jjZhpNXr32qav1BIUhV+dvHq0nan8WkW4bSuE4N\nt0MKKd4mqAJVHevXSIzfHDOc/IZB9Ghhw8mr2aMi8jbwLce2QMx0LyTjC9MSkpm1bCf3nd2FYR3t\nVjS+5m2C+l5EbsO50K9kAbNh5kHgyc8S+SZxL09c0pMzutlwchfcCHTDmcOyuAVCcS6iNUFq455M\n/jF7NcM6NOKuM21wiT94m6Cu9TyPK7FMgQ6+Dcf42v9+SWLC/K3c+Lt2XG/Dyd3SV1V7ux1EWWwU\nX8Vk5xXypw+XUqdGlM1d6Ufe3m6jfSkPS04B7vt1e3nskzWc3b0pf7uwh9vhhLOFnolbA5bNb1kx\nj32yho17D/HC1f1oav25flNmDUpEzlTV70Tk8tLWWxt64Fq7K4O7PlxK95Pq8tI1/e0Mz12nADeI\nyFacJnLBmYaoUsPMjbs+XraTj5bs4E9ndGR4Z5u70J/Ka+I7DWf03shS1lkbeoDak5HDze8tJi7W\nGU5eu4bl1RjdAAAftElEQVQNJ3fZiPI3McFgy75DPDJrFYPaNeC+s7u4HU7IK/ObS1Uf9bx8QlW3\nllznuemaCTCH85zh5OnZ+Uy7YxjN61nzg1tEpK6qZgCZbsdiqi4nv5A/fbiMmKgIXh7dn6hIb28G\nYSrL21PrGcCA45ZNBwb6NhxTFYVFyt2Tl7N2VwZv3xBPzxbWn+CyD4GLcOauVI6dcdwGGQWZJz9L\nJDElgwlj4jmpXk23wwkL5fVBdcOZ5LLecf1QdQE7NQ8w/56byDeJe3j84p6c2a2Z2+GEPVW9yPNs\nrQ1Bbu6qFN5fuI3bTu1gZasalVeD6opzBlifY/uhMgGbsyuAvP9LEu/8vJUxJ7fjhpPbuR2OKUFE\nvlXVs8pb5iYbZn5i21MP89D0lfRrXZ8Hzu3qdjhhpbw+qNnAbBEZpqq/VFNM5bLCdKzv1+/l0U/W\ncFa3pvz9ooAezRxWRCQWqIUz63gDjjbx1QVauBZYKWyy2NLlFRRx1+SliMAro/sTE2X9TtXJ2+ug\nAiY5gV2zUdLaXRncNWkp3ZrX5eXRNpw8wNyO0//UzfNc/JgNvOZiXMZL//l8HSuT03n6yr60bljL\n7XDCjo0/DmIlh5NPGGPDyQONqr4EvCQif1bVV9yOx1TM12v3MGG+02w+oldzt8MJS/aNFqRsOHnw\nsOQUfJLTDvPAtBX0almXcRd0czucsFXeKL4yZzBX1ed9G47xxuKkA4ybuYrNew/ZcHJjfCy/sIi7\nJy+jsEh5dfQAakRFuh1S2CqvBhXnee4KDAKK72kzEpjnr6DMiSUkHWDUeKdLMCpCqF8rxuWIjAkt\nz321gaXbD/LK6P60a1zb7XDCWnmj+B4HEJGvgAGqmul5/xgwze/Rmd946dtNR16rKgu3pDKwbQMX\nIzInIiLHX9x+DFVdWl2xGO98v34v43/czLVD2jCyb0ANtAxL3vZBtQHySrzPA9r5PBpTpsmLtjNv\n474j76OjIhjawW6SFsCeK2OdAmdWVyDlsUs3YHd6DvdPXUG35nH8wy7XCAjeJqj3gUUiMgunYF0G\n/M9vUZnfmLsqhUdmreK0Lk248/SOLNmWxtAOjaz2FMBU9Qy3Y/BWuF8Hdf6L89i8PwsBXr19GLHR\n1u8UCLxKUKr6pIh8Dgz3LLpRVZf5LyxT0s8b93PvlOX0b9OA8b8fSM2YSIZYzSmoiEgvoAclpghT\nVTvJCwAJ29JI3O3M5xsdKaRn57sckSlW3ii+hiXeJnkeR9bZLd/9b/mOg9z2/hLaN67NhBsGUTPG\nzuyCjYg8CpyOk6DmAucDP2OtEAHh1e82HnldVGT9uoGkvBrU8bMwq+dZsNmY/W7T3kzGvLuIRnVi\neP/mwdSrFe12SKZyrgT6AstU9UYRaQa87XJMBpizYhffr7d+3UBV3ig+m4XZJclph/n924uIiojg\ng5uH2G2lg1u2qhaJSIGI1AX2Yid3rvtlcyr3T13BoHYNGHtOF5ZuP2j9ugGm3NttqOq6Ew2XtWGy\n/rH/UC7Xv7OIrLwCpt4+jLaN7FqMILdEROoDb+G0ShwCFrkbUnhbvzuT295fQptGtXjr+njq14ph\nWMfGbodljlNeE99Y4DZKHy4bUMNkQ0VmTj5j3l3ErvRs3r95CN1Pqut2SKaKVPVOz8vxIvIFUFdV\nV7oZUzhLSc9mzLuLqBkdycQbB9nF7gGsvCa+2zzPQTNcNpjl5Bdy6/+WsC4lk7euj2dQu4blf8gE\nPBGZDXwEzFbVJJfDCWvp2fmMmbCYzByndaJVA5uhPJB5NcxcRK4vbbkNk/WdgsIi/jx5GQu3HOCl\na/pxRrembodkfOd54GrgKRFZhJOsPlXVHHfDOiocLtTNLSjk9veXsHnfISbeOJgeLax1ItB5e/et\nQSUew4HHgIv9FFPYKSpSHp65iq/XOrdrv6RfS7dDMj6kqj96mvk6AG8CV+EMlAgYoX6PtaIi5YFp\nK1m45QDPjOrDKZ2tvykYeHuh7p9LvheRejizS5gqUlX+PTeR6QnJ3Ht2Z7tde4gSkZo4kyxfDQwA\n3nM3ovDy1OeJzFmxi4dGdOOy/q3cDsd4qbL3gzoMdPZlIOHq/37YzNs/b+WGYW255yz7lYYiEfkI\nGAJ8gXMn3R9UtcjdqMLHhJ+38tZPW7l+WFvuOM1G9wcTb/ug5nD0It0InCvip/orqHDx4a/beebL\n9VzSrwWPjuyJiN2uPUS9C1yrqoVuBxJu5q5K4Z+freW8ns2sjAUhb2tQz5Z4XQBsU9VkP8QTNuau\nSuGRj1dxetcmPDuqLxERVnBCjYg8qKpPq+oXIjKKEreoEZF/q+pfXQwv5P26JZV7P1rOgDYNeOma\n/kRaGQs6Xg2S8HTy/qiqPwLrgZ3+DSu0/bRxH/dMWcbANg14/bqBREd6O1bFBJlrSrwed9y6EdUZ\nSLjZsCeTW/+3hFYNavL29fE2O3mQKvObUUSGisgPIjJTRPqLyGpgNbBHRKyAVcKy7Wnc/n4CHZvU\n4Z0xNvlriJMTvC7tvfGR3ek5jJmwiBrRkbx342Aa1LYLcYNVeU18rwJ/BeoB3wHnq+pCEekGTMbp\n9DVe2rAnkxsnLqZxnRr876bB1Ktpk7+GOD3B69LeGx/I8MzEkp6dz0e3D6N1Q7sQN5iVl6CiVPUr\nABF5QlUXAnjm5/N7cKFkx4HD/OGdX4mOtMlfw0hfEcnAqS3V9LzG897+AXwsr6CIO95PYNPeQ7x7\n4yB6tQzNa7rCSXkJquRQ2Ozj1tkZoJf2ZeZy/YRFZOcVMvWOYbRpZGd14UBVrf22mhQVKX+ZvoIF\nm1N5blRfhndu4nZIxgfKS1B2BlhFxU0OKenZTLplCN2a2/Qqxvjaf79cx+zlu/jLeV25YqBdiBsq\nypss1s4AqyAnv5Bb3lvC+t2ZvHVDPAPb2uSvJjAF81x8E+dv5Y0ft/D7oW248/SObodjfMjGN/tJ\nQWERd324jMVJB3juqr6c0dUmfzWBK1jn4vt8VQqPf7qWc3o04/GLe9mFuCHGEpQfFBUpD81YxTeJ\ne3jCJn81xi8WJx3gno+W0691fV62C3FDUmXn4vM5EakN/B+QhzNX2SSXQ6oUVeXJuYnMWJrMfWd3\n4Q/D2rkdkjEhZ9PeTG55bwmt6tfknRvsesJQ5dcalIhMEJG9ngt8Sy4fISLrRWSTiDzsWXw5MF1V\nbyWIb+Xxfz9s5p2ftzLm5HbcfVbwtecbE+j2ZORww4TFREdG8N5Ng2loF+KGLH838U3kuCldRCQS\nZ0bn83EmnR0tIj2AVsAOz2ZBOanmpF+38cyX67msf0v+cVEPaw83xscyc/IZ8+5i0g7n8e6YQXYh\nbojzaxOfqs4TkXbHLR4MbFLVLQAiMgW4BEjGSVLLKSNxishtwG0Abdq08X3QlXDBS/PYnZHDgax8\nzuzWlKev7GOTvxrjY3kFRfzxg6Vs2JPJOzfE07tVcA3oMBXnxiCJlhytKYGTmFoCM4ErROR1YM6J\nPqyqb6pqvKrGN2ni/sV4CdvSSEzJ5EBWPiJwyyntbfJXY3xMVXloxkp+3rSf/1zem9NtVGxYcGOQ\nRGlVC1XVLODG6g6mqr5eu/uYG2Ut23GQkzvZ7aSN8aWnv1zPrGU7uf+cLoyKb+12OKaauHGqnwyU\n/A9rBexyIY4qKypS5m/af+R9dFQEQzs0cjEiY0LD1W/8wtVv/ALA/35J4vUfNjN6cBvuOtMGHoUT\nN2pQi4HOItIe575S1wDXuhBHlb27IIlVOzM4qW4NIiMjeOma/gxs28DtsIwJGV+s3s2jn6zh7O5N\n+ecldkfccOPXBCUik4HTgcYikgw8qqrviMhdwJdAJDBBVddUcL+uT8uyYU8m//1iHWd1a8rbN8Rb\nwTHGhzJz8tl/KJc/T15K31b1eWX0AKKsbzfs+HsU3+gTLJ8LzK3CfucAc+Lj42+t7D6qIq+giHun\nLCeuRhT/uaKPJSdjfChhWxqJuzNRdTqs7zmrk12IG6bslKQSXvhmA2tTMnjq8t40iavhdjjGhJQF\nm/ejnpFHIrA2JdPdgIxrLEFV0KKtBxj/42aujm/NuT2bux2OMT4hIiNF5M309HS3QyElPefI6xgb\neBTWLEFVQGZOPmOnLqd1g1r8fWQPt8MxxmcCZTbzxJQMpi3ZQVxsFK3qxzLplqE28CiMBcxkscHg\niTlr2XUwm2l3DKNODfvVGeNLeQVFjJ26gno1o2nTsBbRkRGWnMJcUNag3GiO+GJ1CtMSkrnz9E52\n40Fj/ODlbzeSmJLBU5f3sdlYDBCkNajqHsW3NzOHcTNX0atlXe4+q3N1HNKYsLJsexr/98MmrhjQ\ninN6NOOcHs3cDskEADtNKYeq8uD0lRzOK+TFq/sRE2W/MmN8KSe/kPunraB53Vgevdj6ds1R9m1b\njkm/bueH9fsYd343OjWNczscY0LOM1+uZ8u+LJ6+si91Y6PdDscEEEtQZdiy7xBPfpbI8M6Nud7u\njGuMzy3cksqE+Vv5w9C2nNLZJlk2x7IEdQL5hUXcN3UFMVERPHNlX7u/kzE+dii3gAemraBNw1qM\nu6Cb2+GYABSUgySqw2vfb2LFjoO8em1/mteLdTscY0LOk58lsvNgNtNuH0atGPsqMr8VlDUofw8z\nX7Y9jVe+28Rl/VtyUZ8WfjmGMeHsh/V7mbxoO7cN70B8O7tsw5QuKBOUP696P5xXwNipK2gWV4PH\nLu7p8/0bE+7SD+fz0IyVdG5ah/vO6eJ2OCaAWb36OP+em0hSahYf3jKUejVtRJExvvboJ6tJPZTH\n29cPIjbaZik3JxaUNSh/+X7dXj5YuJ1bTmnPsI42QaUxvvb5qhQ+Xr6Lu87sRO9W7s77ZwKfJSiP\nA1l5/GX6Sro1j+OB87q6HY4xIWf/oVwe+Xg1vVrW5U9n2K3bTfmsiQ9ntohxM1eSkZ3P+zcPpkaU\nNTsY40uqyl9nruJQTgHPX9XP5tozXrH/EmB6QjJfrtnD/ed2oftJdd0Ox5iQ8/HynXy11iljXZrZ\njCzGO0GZoHw5zHzHgcM8PmctQ9o35JbhHXwQnTGmpJT0bP4xew3xbRtYGTMVEpQJylfDzAuLlLFT\nlyPAc1f1JdJmizDGp4onWy4oVJ4dZWXMVExQJihfeXPeFhYnpfH4JT1p1aCW2+EYE3I+XLSdnzbu\n568XdKNd49puh2OCTNgmqNU703n+6/Vc0Ls5l/Vv6XY4xoSc7amHefKzRE7p1JjrhrR1OxwThMIy\nQeXkF3LfR8tpUCuGJy/tjYg1O5jw5uvpw4qKlAemrSBShKev7GOTLZtKCcsE9fQX69m49xDPjOpL\ng9oxbodjjOt8PX3YhPlbWZR0gEcv7kmL+jV9sk8TfsIuQc3ftJ8J87dyw7C2nNalidvhGBNyNu3N\n5Okv13N292ZcMcCaz03lhVWCSj+czwPTVtChSW0ePr+72+EYE3IKCosYO3UFtWMi+fflvaz53FRJ\nWM0k8ffZq9mXmcvMO0+mZozNFmGMr73+w2ZWJqfz2rUDaBpn91EzVROUNajKdOjOXr6TT1bs4p6z\nOtOnVX0/RmdMeFqzK52Xvt3IxX1bcGGfk9wOx4SAoExQFe3Q3XUwm79/vJr+berzx9M7+jk6Y8JP\nbkEhYz9aQYPaMTxxid1HzfhGyDfxFRUpf5m+goIi5YWr+hFlk1Qa43MvfrOR9XsymTAmnvq1bGSs\n8Y2Q/7aeuCCJ+ZtS+ftFPexKdmP8IGFbGm/8uJmr41tzZrdmbodjQkhIJ6gNezL5zxfrOLt7U64Z\n1NrtcIy33r3QeZiAl51XyAPTVnBSvZr87SIbGWt8K2Sb+PIKirh3ynLiakTx1OV9bLirMX7w3y/W\nsXV/Fh/eOoS42Gi3wzEhJiQTVMK2NJ7+Yh1rUzJ46/p4msTVcDskY4JbcY32xs+OLFqwaT8TFyQx\n5uR2nNyxsUuBmYC0YxEk/QTthkPrwZXeTcglqIRtaVz5+gIUiBShoU1lZIzPZebk85fpK+nQuDYP\njejmdjimuhUWQH4W5B2G/MOQl3X0OWUFfPsEoBBVE274pNJJKuQS1MIt+9Ej75SFW1IZ2LaBixGZ\nCsvNgJx05yysCmdfxody0iEn48jf5F+fJpKSns30P9pF70Hh/06Gwwfg1AegQdtjE0r+YU+iKS3h\nnGB5YZ53xy3Mc2pSlqAcQzs0RmQDqhAdFcHQDo3cDsl4K/sgLJsEu1cBChMvhCvehm4jIcKP43lK\nab4yJexYBHtWAwrvns+aQU/x0ZKW3Hl6Rwa0sZO/gHRwB2xfCDsWwqbvIG2Ls3zu/Sf+TFRNiKkF\n0bU9z7UgpjbUaV768uhapS8/sBU+/iOgEBnjNPNVUsglqIFtGzD9jpNZuCWVoR0aWe0pUKlC2lbY\n/qtTiHYsgr2JUKL+S2EeTL0eomKhUWdo3Bkad3Gem3SFRp0g2mbK9ruknyj+u2hRAT1//Qtf1e5E\n+/q3weGmUKuhu/GFu6JC5wRi+6+w/RfY8Stk7HTWxdSB2iUmxZYIGDgGBt1ybKKJruW7k8A2Q6FR\nx/DtgxKRkcDITp06lbp+4OcXMzAnHTq9DVgTUUAoyHXapnf86jmzWwRZe511NepB60HQ83KIrQtf\nPAxaBJE1YOidUJQP+zfArqWwZhZHk5hA/daepHXco3ZjsJGbvtFuuPPFpkXkSwwfFpzJNQ23Ef3F\ng/D136Dr+dDvOuh4FkQG5VdKcMk9BMmLj5al5MWQd8hZV7cltB7iJIk2Q6FpT6fcvHexc8IXGQN9\nR0MzP8/20XqwT5rng/K/SVXnAHPi4+Nv/c3KHYtKNBFdBFdOgA6nO2cK9oVVfbJSnQJUXDvauRQK\nc511DdpBxzOhzRBoPRSadDv27G3ZB06fxxVv//afPD8HDmx2Eta+Dc7z/g2wbYHTNl4str6TqJoc\nl7jqt7Uv0YpqPZh19U9lQ3okE7OHc+Y5FxF7ZmdIWQnLP4RVU2HtbKjTDPpc7SSrpjZwwmcydjk1\no+LWht2rQQsBcRJN32ucctRmqHPCdrzWg52BCj6o0VS30CupJZojKMyFj65zXkdEQ836ULOB8+VV\ns0Hp70tbF1XOSMDxw0/8hRoOVGH/RqfwbP/VSUypG511EdHQoh8MvtU5s2s9BOLKmW2gRl3nUdrv\nMjrWKZTHnwEWFTnNGvvXO7Hs3+A8b/zaSXjFImOgYccSzYVdIHO3c3ZpgzJKlbAtjdEpN5JHFCLC\nw+09TXon9XEe5zwBG790ktUvr8GCl6HFAOh/HfS6wilHxjtFhU5T946FTu1o+6+Qvt1ZF10LWg6E\n4WOdZNRqEMR6eYNJH9VoqlvoJah2wwHhSAfd0D85ySbnIGSnOR3x2WlwaDfsW+e8zy1nVvTo2scl\ns/pH3+ceOrZT/5L/c2oHNRv4t2O/Oh0/iCA/26kR7fj16CM7zVlXs6GThPpf5zy36F89/UQREc7Z\nY/3W0OnsY9dlp8H+TUdrW/s3wt61sO4zz5mox3sXV2lIbKhauCWVPJyLcAVYnJTG4PYlBh9FxUD3\nkc7j0F5YNc0Z7PLZ/fDFOOh2oVOr6nCG1V6LFV8n1DLeadkp2Rebm+FsU6e508ow9I9OQmreGyLD\n62Lo0PtvaT3Y+UNWpEZTVOhsX5zAckoksuyDv01uB7YcfV+QfXQ/hXkw8xbntURArUZOB2Xxc+0m\nTt9I7cZQq3GJZY2cxOdtE2R1jDorLHCazPKzIWufM0R12k1wcJvTl1SU72zXuIvzBVTcxNCoU+A1\npdZs4PRxtR507PKCPPjmUVj4f877Kg6JDVVDOzQiQqBIIaa8kbF1msKwPzl9hykrjjYBrpkFcSd5\nmgCvdQa5hApVp6zkZjonrLkZTp9Q7iFnWV7m0XV5h5zvj03fOP2sJTXp7tQ42wx1Tu4atAu8slTN\nQi9BAdzxU8W2j4h0RiJVZjTS1p/gfxd7OvVjYPj9TrI5vN/5Ys/a7zxSVjjPJ6qtRUSXSFzFjybH\nJTPP8uw05x9+w1dOQc/P9iSTw0df5x33/kTL8rOPW57lPJd2ncOaGdCsN5x8l1OAWg12kmuwioqB\nnpfBr+OP/v2qMCQ2VA1s24BpFR0ZK+I07bboB+f+EzZ8CcsnwYJXYP6LTs2h37XQ63J3mwBVnRPQ\nt891TlKH/hHqnuRJKJmeRFOcXEokmiPrPMuPTzalkQioEecc88j24gxaGPFvawothahq+VsFqPj4\neF2yZInbYVSsD6ogFw6nehKXJ4EdSWb7nMEFWfs8y/YfHZ1TVZExTlNb8ZDS6FrO+5gSr49ZXtt5\n3rYA1s919iGRcOYjThIOJeX8/UQkQVXjXYis2vm9TGXucWpUyybBvkRnpGa3C50m4Q5nOCeLvpKX\nBRkpkJni9DP+5nmX81yQc+J9SKSTVGrEOUO2a8RBDc9zTNxx7+uUsm2J99E1ncS9Y9Gxo+rCsFnZ\n2zIVmjWo6laRGltUDajbwnl4Iz/72GS2/ANnxBQA4rT797zU6ScrTjIxxQmneFnNyrddtx4CG74I\n7RpGWYMyjG/FNYOT/wzD7oKU5Z4mwGmwZqbTBNj3Gqe/KjvtxKPOCnI9CaY42ZSWhHYf7cspKbqW\nc5y4k5xaXFxzp19y41fOeolwmid/d4+TWKJifd/MFsSj6qqbJahAF13zaOc/OIMzEuc4CSMq1ins\n/vwHbz0YbvrSCpPxLRFnAE2L/nDuv2D9506ymv8S/PyC57ordWpUnc5y+kSLE1D2gd/uLyLaSTp1\nT4Km3Z2BSnHNPcmo+dGkVCPutwlnxyKnqb64RtPjEqcvzZ+CdFRddbMEFWzcSBhWmIw/RdVwWgF6\nXuokodl/cgYRABQVQNICaNzJmUOuzRCIa/Hb5FOrYeVrOlajCViWoIKRJQzfsjn4AkdcczjtIUia\nf7RG84eZ/v9/tzIVkCxBGWMCi9VojIclKGNM4LEajQFCZKoDY4wxoSYoE5SIjBSRN9PTy5miyBhj\nTNAKygSlqnNU9bZ69bycKNEYY0zQCcoEZYzxjoh0EJF3RGS627EYU1GWoIwJUCIyQUT2isjq45aP\nEJH1IrJJRB4uax+qukVVb/ZvpMb4h43iMyZwTQReBf5XvEBEIoHXgHOAZGCxiHwCRAJPHff5m1R1\nb/WEaozvWYIyJkCp6jwRaXfc4sHAJlXdAiAiU4BLVPUp4KLKHEdEbgNuA2jTpk2l4zXG14I6QSUk\nJOwXkW0nWF0PqMowv4p+3pvty9qmouu8WdYY2F9OTL4UKr/z45e3rUBM/tYS2FHifTIw5EQbi0gj\n4Emgv4iM8ySyY6jqm8Cbnu33+bFMVXQfVf37lrU+XMpURfdRXb9z78qUqobkA3izOj/vzfZlbVPR\ndd4sA5bY77ziv/Oq/hw+/p22A1aXeD8KeLvE+z8ArwTD37ei+6jq37eif+NQLFOB/jsv7xHKgyTm\nVPPnvdm+rG0qus7bZdUpVH7nbv8ey5IMtC7xvhWwq5qO7YvfS0X2UdW/b1nrw6VMVXQf1f07L1NQ\n37DQlE1ElmiY3GgvVHn6oD5V1V6e91HABuAsYCewGLhWVde4FWM4sTJVvUK5BmU8/QomOInIZOAX\noKuIJIvIzapaANwFfAkkAlMtOVUrK1PVyGpQxhhjApLVoIwxxgQkS1DGGGMCkiUoY4wxAckSlDHG\nmIBkCSqM2MzWxviWlSn/sgQV5Coy47XazNbGlMvKVOCwBBX8JgIjSi4oMeP1+UAPYLSI9Kj+0IwJ\nShOxMhUQLEEFOVWdBxw4bvGRGa9VNQ+YAlxS7cEZE4SsTAUOS1ChqbQZr1uKSCMRGY9nZmt3QjMm\nKFmZckFQ327DnJCUskxVNRW4o7qDMSYEWJlygdWgQpObM14bE4qsTLnAElRoWgx0FpH2IhIDXAN8\n4nJMxgQzK1MusAQV5GzGa2N8y8pU4LDZzI0xxgQkq0EZY4wJSJagjDHGBCRLUMYYYwKSJShjjDEB\nyRKUMcaYgGQJyhhjTECyBBUgROQpETldRC4tnsq/Ap8dIyKv+iu2co7dT0QucOPYxpTFylTwswQV\nOIYAvwKnAT+5HEtF9AOsMJlAZGUqyNlksS4TkWeA84D2OFevdwTOEpHpqvrEcds2AcYDbTyL7lXV\n+d5sIyKPeY5xEtAFGAsMxbm/zU5gpKrmi8hA4HmgDrAfGKOqKSLyA05hPwOoD9zsef8EUFNETgGe\nAnYDL3mOrcCpqppZpV+SMRVgZSqEqKo9XH7g3GvmFSAamF/Gdh8Cp3hetwESPa/HAK+Ws81jwM+e\nY/QFDgPne9bNAi71rFsANPEsvxqY4Hn9A/Cc5/UFwDfHH9vzfg7wO8/rOkCU279fe4Tfw8pUaDys\nBhUY+gPLgW7A2jK2OxvoIXJk5v+6IhJXgW0+V+eMbhUQCXzhWb4KaAd0BXoBX3s+HwmklNj3TM9z\ngmf70swHnheRScBMVU0u4+cxxl+sTIUAS1AuEpF+OLeXboVT9a/lLJblwDBVzT7uIxGlLS9RcMrb\nJhdAVYtEJF89p2RAEc7/ggBrVHXYCULO9TwXcoL/HVX9j4h8hnNGuFBEzlbVdSfYnzE+ZWUqtNgg\nCRep6nJV7QdsAHoA3wHnqWq/UgoSwFc4MyoDRwpjZbY5kfVAExEZ5vlstIj0LOczmcCRM04R6aiq\nq1T1v8ASnDNYY6qFlanQYgnKZZ4O2DRVLQK6qWpZzRF3A/EislJE1lL6nTy92aZUqpoHXAn8V0RW\n4DSRnFzOx77Haf5YLiJXA/eKyGrP57OBz709vjG+YGUqdNjtNowxxgQkq0EZY4wJSJagjDHGBCRL\nUMYYYwKSJShjjDEByRKUMcaYgGQJyhhjTECyBGWMMSYg/T/12Gpzl9IHFwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x12503a2b0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, (ax1, ax2) = plt.subplots(1, 2, True)\n",
    "\n",
    "# Build times\n",
    "for mode, times in build_times.items():\n",
    "    ax1.errorbar(sizes, np.mean(times, axis=1), np.std(times, axis=1) / np.sqrt(num_runs - 1), \n",
    "                 marker='.', label=mode)\n",
    "ax1.set_xscale('log')\n",
    "ax1.set_yscale('log')\n",
    "ax1.legend()\n",
    "ax1.set_ylabel('Build times [s]')\n",
    "ax1.set_xlabel('# elements')\n",
    "\n",
    "for mode, times in eval_times.items():\n",
    "    ax2.errorbar(sizes, np.mean(times, axis=1), np.std(times, axis=1) / np.sqrt(num_runs - 1), \n",
    "                 marker='.', label=mode)\n",
    "ax2.set_yscale('log')\n",
    "ax2.legend()\n",
    "ax2.set_title('Evaluation times')\n",
    "ax2.set_ylabel('Eval times [s]')\n",
    "ax2.set_xlabel('# elements')\n",
    "\n",
    "fig.tight_layout()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}