Capacitively-coupled-current-biased-JJ-phase-qubit.ipynb

{
 "metadata": {
  "name": "Capacitively-coupled-current-biased-JJ-phase-qubit"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "# Capacitively coupled current-biased Josephson junction phase qubits\n",
      "\n",
      "Robert Johansson (robert@riken.jp)"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "%pylab inline"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 1
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from scipy import *\n",
      "from scipy import optimize"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 2
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from wavefunction import *"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 3
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## Problem parameters"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "h = 6.626e-34\n",
      "h_ = h/(2*pi)\n",
      "e = 1.602e-19\n",
      "Phi0 = h / (2 * e)\n",
      "cf = 1\n",
      "\n",
      "H  = 6.626e-34\n",
      "H_ = H/(2*pi)\n",
      "\n",
      "Ic = 13.3e-6\n",
      "\n",
      "#Cj = 4.3e-12\n",
      "\n",
      "ceta = 0.005\n",
      "Cj = 4.3e-12\n",
      "j1 = 0.975\n",
      "\n",
      "dj2 = 0.001\n",
      "J2 = arange(0.9701, 0.98, dj2) \n",
      "\n",
      "# discretize the x coordinate\n",
      "xmin = -0.7 * pi\n",
      "xmax =  4 * pi\n",
      "N = 750\n",
      "D = arange(xmin, xmax, (xmax-xmin)/N)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 4
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## Potential"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def U_current_biased(x, args):\n",
      "    \"\"\"\n",
      "    Potential for a current-biased phase qubit potential\n",
      "    (the washboard potential)\n",
      "    \"\"\"\n",
      "    Ej = args['Ej']\n",
      "    Ic = args['Ic']\n",
      "    Ib = args['Ib']\n",
      "\n",
      "    u = - Ej * (cos(x) + Ib / Ic * x)\n",
      "\n",
      "    return u"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 5
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "def jjpq_current_biased_me(ceta, Cj, cf_p, J):\n",
      "\n",
      "    Ib = J * Ic \n",
      "    mm = Cj * (1 + ceta) * (Phi0/(2*pi))**2 * cf\n",
      "    Ej = Phi0/(2*pi) * Ic / cf\n",
      "    \n",
      "    U_args = {'Ej': Ej, 'Ic': Ic, 'Ib': Ib}\n",
      "    S_param = array([h_, mm])\n",
      "    \n",
      "    # Calculate the eigenfunctions\n",
      "    M = schrodinger_matrix(xmin, xmax, N, S_param, U_current_biased, U_args)\n",
      "    evals,evec = eigenvectors_sorted(M);\n",
      "    \n",
      "    # Find bound states\n",
      "    dmin = math.asin(Ib/Ic)\n",
      "    dmax = -math.asin(Ib/Ic) + pi\n",
      "    delta = dmax-dmin\n",
      "    \n",
      "    Umin = U_current_biased(dmin, U_args);\n",
      "    Umax = U_current_biased(dmax, U_args);\n",
      "    dU = Umax - Umin;\n",
      "    \n",
      "    boundidx = []\n",
      "    for i in range(0,N):\n",
      "        if evals[i] > Umin-0.5*dU and evals[i] < Umax+0.5*dU:\n",
      "            if inner(evec[:,i] * (dmin-delta < D) * (D < dmax), evec[:,i]) > 0.85:\n",
      "                boundidx.append(i)\n",
      "\n",
      "    # Evaluate matrix elements:\n",
      "    inner_prod = zeros((len(boundidx), len(boundidx))).astype(float)\n",
      "    expect_pos = zeros((len(boundidx), len(boundidx))).astype(float)\n",
      "    expect_kin = zeros((len(boundidx), len(boundidx))).astype(float)\n",
      "    energy_levels = zeros(len(boundidx)).astype(float)\n",
      "    \n",
      "    for i in range(0, len(boundidx)):\n",
      "        l = boundidx[i]\n",
      "        energy_levels[i] = evals[l]\n",
      "    \n",
      "        for j in range(0, len(boundidx)):\n",
      "            k = boundidx[j]\n",
      "    \n",
      "            inner_prod[i,j] = inner(evec[:,k], evec[:,l])\n",
      "            expect_pos[i,j] = inner(evec[:,k], evec[:,l] * D)\n",
      "            expect_kin[i,j] = inner(evec[:,k], derivative(evec[:,l], D))\n",
      "    \n",
      "    return energy_levels, expect_pos, expect_kin"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 6
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "E1, Ed1, Edd1 = jjpq_current_biased_me(ceta, Cj, cf, j1)\n",
      "P1 = h_ /1j  * Edd1 \n",
      "\n",
      "eval_idx = 0\n",
      "eval_store = zeros((len(J2),9)).astype(float)\n",
      "\n",
      "for j2 in J2:\n",
      "\n",
      "    E2, Ed2, Edd2 = jjpq_current_biased_me(ceta, Cj, cf, j2);\n",
      "    P2 = h_ / 1j * Edd2;\n",
      "\n",
      "    H = zeros((9,9)).astype(float)\n",
      "    Hint = zeros((9,9)).astype(float)\n",
      "    \n",
      "    for K1 in range(0,3):\n",
      "        for K2 in range(0,3):\n",
      "            for L1 in range(0,3):\n",
      "                for L2 in range(0,3):\n",
      "                    I1 = K1 + 3 * (L1)\n",
      "                    I2 = K2 + 3 * (L2)\n",
      "               \n",
      "                    H[I1,I2] =            (E1[K1] - (1*E1[0]+1*E1[1])/2) * (K1 == K2) * (L1 == L2)\n",
      "                    H[I1,I2] = H[I1,I2] + (E2[L1] - (1*E2[0]+1*E2[1])/2) * (K1 == K2) * (L1 == L2)\n",
      "            \n",
      "                    Hint[I1,I2] = (2*pi / Phi0)**2 * ceta / (Cj * ( 1 + ceta)) * P1[K1,K2] * P2[L1,L2] / cf\n",
      "\n",
      "    H = H + Hint\t\t\n",
      "    H = H / (2*pi * 1e9)\n",
      "\n",
      "    evals, evecs =  eigenvectors_sorted(H)\n",
      "    eval_store[eval_idx, :] = evals\n",
      "    eval_idx = eval_idx + 1"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stderr",
       "text": [
        "-c:25: ComplexWarning: Casting complex values to real discards the imaginary part\n"
       ]
      }
     ],
     "prompt_number": 7
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "fig, ax = subplots(1, 1, figsize=(12,6))\n",
      " \n",
      "for n in range(6):\n",
      "    ax.plot(J2, eval_store[:,n])\n",
      "    \n",
      "ax.axis('tight');\n",
      "ax.set_xlabel(r'$J_2$')\n",
      "ax.set_ylabel(r'$E_n$')\n",
      "ax.set_title(r'Energy levels of ...')"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "pyout",
       "prompt_number": 8,
       "text": [
        "<matplotlib.text.Text at 0x4713310>"
       ]
      },
      {
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAtAAAAGLCAYAAADu9qU1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XuUFNWBx/FfVXfPg4AIBuUxgAiigyCwCPhYY/vAVVcU\nEV0xRhMxm2iyWeNqiJ5kz5hz1MTHqsSsZk18rYnmmHUXjEpUsNWNIhB8rUYRFYQxCCKgCExPV939\nox90Tz/vzPR09/D9nFOnbt26VXVbyplf3a6qcYwxRgAAAABK4la6AwAAAEAtIUADAAAAFgjQAAAA\ngAUCNAAAAGCBAA0AAABYIEADAAAAFgjQAFDjXNfV+++/X9ZjhMNh/frXv+7Wfb7zzjuaNGmS9tln\nH91xxx3dum8AKKdgpTsAANXiwAMP1KZNmxQIBFJ13/jGN7RgwYIK9qo6OI4jx3G6dZ833nijTjzx\nRL366qvdul8AKDcCNAAkOI6jP/zhDzrhhBPKehzP8zJC+t5q3bp1OvrooyvdDQCwxi0cAFCC++67\nT3/7t3+rq666SgMHDtRBBx2kxYsXp9Zv375d8+bN09ChQ9XU1KQf//jH8n0/te0xxxyjK664Ql/+\n8pd17bXX6tNPP9XMmTPVv39/TZs2TT/60Y907LHHSpK+853v6Morr8w4/hlnnKHbbrutaD/b2tp0\n5ZVXauTIkRo8eLAuvfRS7d69W5LU3Nysxx9/PNU2Fotp0KBBqRHgZcuW6eijj9aAAQM0adIkPffc\nczmPsWbNGh133HHad999NWjQIJ133nl5+7No0SIddthhGjBggI4//ni9/fbbkqQTTjhBkUhE3/3u\nd7XPPvtozZo1RT8bAFSLmgrQF198sQ444ABNmDChy/tat26dpkyZosmTJ+uwww7T7bffntXme9/7\nnvr169flYwGoHcaYvOuWL1+uQw89VFu2bNEPfvADzZs3L7Xu61//uurq6vTee+/plVde0VNPPaVf\n/epXGduOHj1amzZt0jXXXKPLLrtM/fr108cff6z7779fDzzwQOoWia9//et66KGHUn355JNPtGTJ\nEn31q18t2v8f/vCHWrNmjV577TWtWbNGra2t+slPfiJJOv/88/XQQw+l2v7xj3/U/vvvr0mTJqm1\ntVWnn366/vVf/1Vbt27VzTffrLPPPltbtmzJOsaPf/xjnXLKKdq2bZtaW1v1ve99L2dfVq9erfPP\nP18LFizQJ598otNOO00zZ85ULBbT0qVLdeyxx+oXv/iFPvvsM40ZM6boZwOAalFTAfob3/hGxohP\nVwwdOlTLli3TK6+8ouXLl+vWW2/Vhg0bUutXrlypbdu2dfs9fwCqlzFGs2bN0oABA1JT+oNzI0eO\n1Lx58+Q4ji688EL99a9/1aZNm/Txxx/rySef1K233qrGxkYNGjRIl19+uR5++OHUtkOHDtV3vvMd\nua6rUCikRx99VNdee60aGhrU3Nysiy66KBWYp06dqv79+2vJkiWSpIcffljHH3+8Bg0aVLT/d999\nt/7t3/5N++67r/r27aurr7461Y+5c+dq0aJFqRHp3/72t5o7d64k6cEHH9Rpp52mU045RZJ00kkn\n6YgjjsgYsU6qq6vT2rVr1draqrq6ury3Yfzud7/T6aefrhNPPFGBQEBXXnmldu3apRdffDGjzwBQ\na2oqQB977LEaMGBARt17772nU089VUcccYS+8pWv6J133ilpX6FQSKFQSJK0a9cuhUIh9enTR1L8\n/sQf/OAHuvHGG/nhDuxFHMfRwoULtXXr1tSUPso8ePDgVDn582LHjh1at26d2tvbNWTIkFTw/va3\nv63Nmzen2g8fPjxV3rx5s2KxWEZdU1NTRl8uvPBCPfjgg5Li4fZrX/ta0f5v3rxZO3fu1JQpU1L9\nOPXUU/XJJ59IksaMGaPm5mYtWrRIO3fu1GOPPabzzz9fUvxbuUceeSTj4uFPf/qTNm7cmHWc5M/G\nadOmafz48br33ntz9uevf/2rRowYkVp2HEfDhw9Xa2trRh0A1Jqaf4jwH//xH/XLX/5SY8aM0csv\nv6zLLrssNWpTzIYNG3TaaadpzZo1uvnmmzVw4EBJ0h133KEzzzwz45clAOQzfPhw1dfXa8uWLXLd\n3OMS6UFx0KBBCgaDWr9+vQ4++GBJ0vr16zPaX3DBBZowYYJee+01vf3225o1a1bRfnz5y19WY2Oj\n3nrrLQ0ZMiRnm7lz5+qhhx6S53kaN26cDjroIEnSiBEj9LWvfU3/8R//UfQ4BxxwQKrdn/70J510\n0kk67rjjUvtKGjp0qN54443UsjFG69ev17Bhw4oeAwCqWU2NQHe0Y8cOvfTSSzrnnHM0efJkffvb\n306Nljz66KOaMGFC1nTqqaemtm9qatLrr7+u9957T7fddpvWrFmjjz76SL///e/13e9+l9FnYC/U\nmf/vhwwZopNPPllXXHGFPv/8c/m+r/fee0/PP/98zvaBQECzZ89WS0uLdu3apbffflv/+Z//mRGy\nm5qadMQRR+jCCy/UnDlzVF9fX7Qfruvqm9/8pi6//PLU6Hdra6ueeuqpVJvzzjtPf/zjH3XXXXdl\n3FN9wQUX6LHHHtNTTz0lz/O0e/duRSKRjNHipEceeSR1y9u+++4rx3FyXjice+65evzxx7V06VK1\nt7frlltuUUNDQ8YtH/ycBVCLajpA+76vfffdV6+88kpqevPNNyVJs2fP1htvvJE1Pfnkk1n7GTJk\niI499li9+uqrevXVV7VmzRqNGTNGBx10kHbu3KmxY8f29EcDUCEzZ85Uv379UtPZZ58tKfd7kNOX\nH3jgAUWjUY0bN04DBw7UOeeck7qgz7XtHXfcoe3bt2vw4MG66KKLNHfuXNXV1WW0ueiii/TGG28U\nvX0jfd8/+9nPNGbMGB155JHq37+/ZsyYodWrV6fWDx48WEcffbReeukl/cM//EOqvqmpSQsXLtT1\n11+v/fffXyNGjNAtt9ySM+CuXLlSRx55pPr166czzzxTCxYs0IEHHpjVbuzYsXrwwQf1T//0Txo0\naJAef/xxPfbYYwoG93z5WegWjksvvVSXXnppann8+PGphyA//PBD9evXL+PZFQDoKY6pscv/tWvX\naubMmamvBY855hh9//vf15w5c2SM0RtvvKHDDz+86H5aW1s1cOBANTY2auvWrTrqqKO0aNGirLDc\nr18/ff7552X5LACQNH/+fG3atCnjfuIXXnhBF1xwgdatW1fBngEAOqqpEei5c+fq6KOP1jvvvKPh\nw4fr3nvv1W9+8xv9+te/1qRJkzR+/HgtWrSopH395S9/0ZFHHqlJkybphBNO0DXXXJNzpJkHXACU\nwzvvvKPXX39dxhgtX75c99xzj84666zU+vb2dt1222365je/WcFeAgByqbkRaADoDVauXKm5c+fq\no48+0gEHHKBvfetbmj9/vqT4Bf7UqVM1adIkLV68WH379q1wbwEA6WoqQDMaDAAAgJ5QKCLX3Gvs\nIpE6GdMu162X6zYkpsa0cvZyIFB4vc32jlNPkK8hLS0tamlpqXQ30ItwTqG7cU6hu3FOdV2xrFdz\nAfq449pkjC/fb5Pv75Lv706bCi973q7EfKfa2z/t1PbZ4b3nAjzhHQAAoPJqLkBLkuO4CgQaFQg0\n9vix94T34oE7VwD3vF1qb99qvX08vEcT4d02kOefssN9/slxavJ0AQAA6FY1l4iWHbRMTsCRE3Ak\nVyWXbdpmlF1HChQqu3LcL0mBL2XWu4n9JMrBgKNgjvqscr7+hCS5RgpE5bu7ZdyojLtbxomXfe1O\nLLfFy85uGScqozYZv02ev0ueditmtslXm3yzW8bs7hDY80+etyvxxxIKhezSw3jxqeO+6uU4Aevz\nJRwOd/s5iL0b5xS6G+cUuhvnVPnV3EOEO9fslPGM5EvGM6WVfSN5nSj7if3Ylm36llbu7HZ5P2ty\nimXXyVdGQHeCaWE+T51CnpyGdjkNUZn6qJz6qFQflRra5dRFZeqiUl1UTmKuuqhMKCqFolKoTQpF\nZYLxsglGpWBUJtAmBdtkgm0ygTaZQFRyE2U3c5IJyjH1ck29HFMvxzTIVaKseNlV/FYX12mQ4yTq\nkoHc6RDSAw0KBOrlhEIK1Afl1gfl1ofk1gXlBoKSAnIcNxHcA3KcQKLspsp71pXaLiAp+49qAACA\n6uE4TsGHCGsuQNdQd6uaMTlCtSeZmCmtLhHMS6pLq+9UnWfkx3zJb4+Pnmu3jHbLN20yTlu87MTL\nvrNbctviy+7uePB2o/EQngzlyaAeSAR3eZLvyygmY3zJ8aSALwWNnKAvBX05iWW5vhRILAf81LKc\nRNn1ZJzEshPfl3F8SZ6M4yl+5WIkuXKSgdrZE8AzQ7ir9EDeMZR3X6jPbpf/WLn6lG+51Lb5jl/a\nZ8/ftqZecw8AqCIEaMCSMUamzcjf7cvb5cnf5acmb5cnf7efXbfLz6jvuF36vrxdnvy2dnlt7fKj\nMfltMZn2mFRvFOjjyO1j5PaRnEYp0EdyGo3cxsS8QXIbJbfByElO9ZLbIDn1vpx6I6feyK1TvFxn\n5NT5ibmRE4rPFfLjFweOL2N8GeNJ8mSMlyj7qXLmOj+jXcft9ixn7zP3sl9kP/Ztk+0kFQ3b+b4l\n6JkLh2CJ+6u2ZS5MAPR+BGigBiRDe0bw3l08jJcc4nPtq82XW+fKbXDlNu6ZAo2BeLkhR13H5T55\n1vXJvZ0b6rnwVTzEl3ox0N0XDrESLxiqc1kqfmFivxws0iZYwn6KtSnlgqVwm/iD1IU/V/E26Z+V\nW7mAalUsc9bcQ4RAb+Q4jpwGR26DKw3omWMaEx9l78yIurfLU/uW9tzrdno5t/N2enJcp8shPNAn\nkDPYp9en1gVDij+Bi+5S+rcLdsvGxIq0KeXCI71Nu3x/d4c2sSL9iRXsYyltin+OZBtf8echSr+g\nyA7nHdt1vIjIdUEQ7NZ95rpg6P59Zv/3iH97xAXI3ib+8yeW9v9h8v+3Pcsd551dXwwBGthLOY6j\nQGNAgcZAj4R2Y4xMu8kO1zvzjKjv3FP2vvDU/kl7/oCeZx+O65Q3oHdY5wR6/y/0+D3rrrgw6Zr4\nyFa+25IKXSxkBvjswN5x244XDbkuAPLtM3kRYrfPXBcRXe1nx+0lkzOU73n+odS5Y9k+e/uubFup\nY+/5Vq07g2f3hdd86yXJcUIdLsSSF1zBrPqurC+GWzgA9ErGGJmosRolLzaCnnOkPi28O0FnT7hu\ncPdMja7VsvX2DW78FZpAhZX6gHrW26E6PkDe8YH0jtt7vkwsJt+LyfiefC8m+bH427CML+N7kvFl\n/LTnPBJlGV++70nGxOuNH++3SV7QpNf72XPtaSfjyUjxdUrsM+88vg/Jj5cTy0Ymo86kHjb3UuvS\n643iFxByku0S65y0tsl1To5lx5eMk3izlSuZoCRXjglKxpWjoGQCchSUI1fKmCfrA5KTWHaC8RDq\nuGkhNB5yXTex7AbjD8y7QblOQE4gEYIDIbluYr0bkBsIxtsEQnICgfjcTXybGAjIDQZSbwjLNY93\nq3PrFcj864PcwgFgr+Q4jpx6R269q+C+5f9Rlx7Yk7fGpB46Tb9VpsBy+/b2jOWSt2/z5YSc0gI4\n4b2gjFeTlvr2oVLeWmTRJuerR0sNlrbBtLv32fEVqenBJd9rU4N72mYEnAJ1ufYp193z9xhcFZwH\nksuOcrcJFN+H3ETgKqVdrnm+Y5c67+Sx5Ugyyvy3jWX+26bPS15vs4/2ju3i672s7aO5t+9KP/Ot\nT567wcT5VAQj0ABQ43KF984E+NRDqru9xPbJ5fT1nvzdRiaxbBLh3Wlw41OjK6d+z7IaHKnBlVOf\nVm5wZOpdBR1HwcSbIAOeIzftPfZWobLUNnlek5m+LKPsd+IXCnaWbbKCn02bEkJkwWBaQljtlrDL\nvcmoQRnfnsSMgn2DBTMnARpA2fnGyDPxLx69ZNnEv5D08qzrbFu/yH6Sy34Jx7Rt66f1wSSWfSUe\n2ExfTiv7xmS2LcO61HKJ66TEQJUT/+LWTcyd9HJynaS69vjfTmqIOqqPxpfr26T6qKP6qFEoGi/X\ntRnVtcfnobb4F9RtMtrlGu12fLU7UjDoKBh0FQq6CgYdhYKuQqH4cl3QUSgUUF2iXB+Ml+tDrurr\nAmoIuWoIuqoPBdQYCqghFFBj0FVjXaJcF1CfUEChoCM36OYOl71oNB1A5/EaOyAPUyQg5QpspbbL\nF+RyBclS23W1T/lCYKl96kyQTK6TpIDjKJCcJyfFA1m+denLbvpyN7V1i+wnueyWcMz0tvlCp+s4\nuZfLuK5jEC5lXaVGEGPGaJfnaZfva6fva6fnlTTfVWK7nZ6nL3xfktTHddUnEMg7b3Tdom3S5405\n6oOMxAI1iwBd5UxixKdjaOlsuWPwydc+V8jx85RrYb8Ft82zr/jfAywcpnKFu1Lb5QtyucJYV9pV\nuu+lhE6XIIEq0u771iF9ZzLYW7QPOE7+wG0RztPnX0ovu64aXF7nBpRDrwvQP1u3rkuh0mq7EkJg\nrrJvcoe+XOX4m0CVChtunmBTrJwvCHYMQR3DUsdAVChQFtpXLe63kiNtAHo3Y4zajemWkP5FruXE\naHrU99XouvFgnQjVycCdrEuVO4TwXHW59hFyHH5WYq/T697Csbm9PSsQ1UlyXbfLATNXudPbFRiJ\n61jmBxMA9C6O46jOcVTnuto3WL5ftZ4x2pUWqjsG7i9y1G1tb1dr2i0tuYJ5ep2RyhLM0+sC/B5E\njam5Eega6i4AADWvPXG7S7FgnlEuIZinj7QHHafbgnlqSmtbzyg6LPW6WzhqqLsAAKAIY4zajOl0\nME+Wv+hQ/iKt7BmTN1x/qUDwzrUuffQ8Wc8Ieu9DgAYAAHu19FH0juE6o77DukKhPD3Yh1w3b7gu\nFLxLWceDopVBgAYAACgTY4x2+7518C51XdT3uz2Upy/zusXcCNAAAAA1Kpbj9hbbUP5F2u0wyXU7\nErfAhFxXfRPBum9asO6bPk+E7r4d1qVv03HbWr/vnAANAACALMnR8x2JYL0jLXTvSJ+nhe70efo2\nHbeNGZM3XHclrPfUW1sI0AAAAOhR7R1ua+lqIE+u2+l5augQsguNhNuE9bq0UfNe9x5oAAAAVLeQ\n62rfMrwH3STefV5qIN/uefooGs0byNPbGykVqIthBBoAAAB7vWhyxNz3NbyhgVs4AAAAgFIVy5xu\nD/YFAAAAqHkEaAAAAMACARoAAACwQIAGAAAALBCgAQAAAAsEaAAAAMACARoAAACwUFUBetu2bZoz\nZ46am5s1btw4LVu2rNJdAgAAADJU1Z/y/ud//meddtpp+v3vf69YLKYvvvii0l0CAAAAMlRNgN6+\nfbteeOEF3X///ZKkYDCo/v37Z7VraWlJlcPhsMLhcA/1EAAAAL1RJBJRJBIpuX3V/CnvV199Vd/6\n1rc0btw4vfbaa5oyZYpuv/129enTJ9WGP+UNAACAcquZP+Udi8W0atUqXXbZZVq1apW+9KUv6ac/\n/WmluwUAAABkqJoA3dTUpKamJk2dOlWSNGfOHK1atarCvQIAAAAyVU2AHjx4sIYPH67Vq1dLkp55\n5hkddthhFe4VAAAAkKlq7oGWpNdee02XXHKJotGoRo8erXvvvTfjQULugQYAAEC5FcucVRWgiyFA\nAwAAoNxq5iFCAAAAoBYQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEa\nAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAA\nsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECA\nBgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAA\nACwQoAEAAAALVRegPc/T5MmTNXPmzEp3BQAAAMhSdQH69ttv17hx4+Q4TqW7AgAAAGSpqgC9YcMG\nPfHEE7rkkktkjKl0dwAAAIAswUp3IN33v/993XTTTfrss8/ytmlpaUmVw+GwwuFw+TsGAACAXisS\niSgSiZTc3jFVMtT7hz/8QU8++aR+8YtfKBKJ6JZbbtFjjz2W0cZxHEamAQAAUFbFMmfV3MLx4osv\natGiRRo1apTmzp2rpUuX6sILL6x0twAAAIAMVTMCne65557TzTffzAg0AAAAelzNjEB3xFs4AAAA\nUI2qcgQ6H0agAQAAUG41OwINAAAAVCMCNAAAAGCBAA0AAABYIEADAAAAFgjQAAAAgAUCNAAAAGCB\nAA0AAABYIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0A\nAABYIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABY\nIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYIEAD\nAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYqKoAvX79eh1//PE67LDDNH78eC1YsKDSXQIAAAAyOMYY\nU+lOJG3cuFEbN27UpEmTtGPHDk2ZMkX/8z//o+bmZkmS4ziqou4CAACgFyqWOatqBHrw4MGaNGmS\nJKlv375qbm7WRx99VOFeAQAAAHsEK92BfNauXatXXnlF06dPz6hvaWlJlcPhsMLhcM92DAAAAL1K\nJBJRJBIpuX1V3cKRtGPHDoXDYf3oRz/SrFmzUvXcwgEAAIByq6lbOCSpvb1dZ599ti644IKM8AwA\nAABUg6oagTbG6KKLLtJ+++2nW2+9NWs9I9AAAAAot2KZs6oC9P/+7//qK1/5ig4//HA5jiNJuuGG\nG3TKKadIIkADAACg/GoqQBdDgAYAAEC51dw90AAAAEA1I0ADAAAAFgjQAAAAgAUCNAAAAGCBAA0A\nAABYIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABY\nIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYIEAD\nAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYIEADAAAA\nFgjQAAAAgAUCNAAAAGCBAA0AAABYCHZ2w9/85jfq06ePGhsbdcopp3Rnnwrbvl0KBOKT6+4pO07P\n9QEAAAB7rU4H6HHjxikUCunPf/5zd/anuOHDJc+TfD8+T5al7FBdrFzu9dV2LNtlm7aOw0UMAAC1\nwpg9WSoW6/q8O/ZRTccqotMB+r//+781cuRInXTSSZ3dRZbFixfr8ssvl+d5uuSSSzR//vzsRp99\nlntj388O1YXKXV3fnftKztvbpba27jtWofbd3daY0sN3uUJ8V9qWo67cx+GCBQD2yBUIc4Uj23K5\n2lbDMRxHCgbjv1O6Mu+OfSTnjY09d6xC88bGgqebY4wx+VbOmzdPM2fO1PTp0zVkyBB9/vnnamho\nUCgU0ksvvaShQ4fq5Zdf1rnnntvl897zPB1yyCF65plnNGzYME2dOlUPPfSQmpub93TWcVSgu6gk\nY8oTzMvRtpTlrtaVe5++H//BVw3hP99Uatvubted++QiBT0t/Wdprik9BHXHcjn22RPHyFeXHgjT\nw1CugGRbroa23XGM9GWXR+HyKZY5C45ADx48WLNmzUot19XVacmSJXr33Xc1ZswYHXXUURo5cmS3\ndHT58uUaM2aMDjzwQEnSeeedp4ULF2YEaElqaWlJlcPhsMLhcLccH12U/KEV7PSXGrCRPtLSU0G/\nUH2p7drbpd277ffXHce2bZd+kdKZ8J5r3pm67tpPT++7lOOl/xvkCkPFwpJtm2o4RrE2vp//fMsX\ngkpd7o59JEcIy32MzuwzeY4BnRCJRBSJREpuXzDtjBo1SpL0+OOP66233tK0adN00kkn6eSTT9YJ\nJ5ygU089tUudTdfa2qrhw4enlpuamvTyyy9ntUsP0MBeK330GSUzxsg3vjzjxee+J894ued+TF6s\nXV4sKuN58tqj8r12+e3t8r12ee3t8r2o/PaYTLLe92Ri7TK+Jz/WLuN5Ml4sPvc9+bFYqk6+l7ne\ni0m+L+N7MjEvsd6XvJiM70leW2J9vC7jAsr3pWiynFbnxedOsp0xiXJibowc35fj+fE6kygbIzcx\nd/xEG9/ITcydVH2yTnKMkZtY5ybXGe2pS5TdRH3AN/IdR54rea4jP2OSfNfdsxxw5LuujOvIDyTm\nrivjujIBR8Z15QcCMoFknSulyoF4OfH/S3I5HuQDUjBtXp8MZUE5GfM6OcGAnFBITiAYn4Jp82BQ\nTiAkNxiUG6qTEwjKDYbkBONzNxCSGwrJCYQUCNXJDYXkBkMKBOvi65PlUHweDIQUcAMKOIGMedAN\nZtQ5fEMCdJuOg7LXXnttwfYFA3Ry6Prv//7v9eyzz+qQQw6RMUau62rOnDld720afhBgb2eMSYW3\nZMjrGOxyhb6u1OU7ls3xS25X4rG72rbQsX3jy5GTCiCu4+YMKh3nruMWbdOxbXLfqeWAKzfoZrXJ\nbBeU69Rl1GXtTx23ybWf0upy9bWcdR3r08/5mB8r+O8a82NF/+2Tbfw8bYodI6NNwfVRef6u3H2I\nxeRFS+tn0T4UWOcbP+u8DLrBgoE7V12xbWzXp+q68Ri22xQru45L5kCXFQzQ11xzjSKRiI455hj1\n799f+++/v9zE1yONRW6utjVs2DCtX78+tbx+/Xo1NTV16zFQuvSRulzhJT2kdAwsFW9rua1NiOxq\n2CxUZ2QyfiHmCm35glyuINiVOpvjh9yQGgINJe+3s+G01P0Va8svTvQGyZ/RHcN2x3BebH1ntil1\nfTQWLdvwNJbYAAAQFElEQVQxim2Tq5xsb2RKDtvpQd2mbBPos7YvQ98KbV9s2XW4LSaXgg8R/vKX\nv9T06dO1bNkyrVixQq+//roCgYAmTpyoTz/9VI888ki3dSQWi+mQQw7RkiVLNHToUE2bNi3nQ4Sf\n7vy0YHjKVc4XBIuVbY6R83g91M98fS7UvlhbI5MxUpcr6KSHlI6BpeJti2xbajjsSp1tiGVUBADK\nr+PvvlKDdynlYvsqeIwy9CdX30pZl76czAKFgnapYbw7loNO+Y8VdIM6dNChBR8iLBigc/n888+1\nYsUK3X777Vq4cGGXT+R0Tz75ZOo1dvPmzdPVV1+d2VnHUf8b+hcNT7nK+cJTsbLNMXIer4r6WSyA\ndtwfYQ4AgL1b6puOPEG7UPju7uWePNbqf1rdvQE6acWKFZo6dWqn/0E6g9fYAQAAoNyKZc5OB+hK\nIEADAACg3IplTu4MBwAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADA\nAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEa\nAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAA\nsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECA\nBgAAACxUTYC+6qqr1NzcrIkTJ2r27Nnavn17pbsEAAAAZKmaAH3yySfrzTff1GuvvaaxY8fqhhtu\nqHSXAAAAgCxVE6BnzJgh1413Z/r06dqwYUOFewQAAABkC1a6A7ncc889mjt3bs51LS0tqXI4HFY4\nHO6ZTgEAAKBXikQiikQiJbd3jDGmfN3JNGPGDG3cuDGr/vrrr9fMmTMlSdddd51WrVql//qv/8pq\n5ziOerC7AAAA2AsVy5w9GqCLue+++3T33XdryZIlamhoyFpPgAYAAEC5FcucVXMLx+LFi3XTTTfp\nueeeyxmeAQAAgGpQNSPQBx98sKLRqAYOHChJOuqoo/Tv//7vGW0YgQYAAEC51dQtHMUQoAEAAFBu\nxTJn1bzGDgAAAKgFBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAA\nACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQ\noAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEA\nAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEAAAAL\nBGgAAADAAgEaAAAAsFBVAfqWW26R67r69NNPK90VAAAAIKeqCdDr16/X008/rZEjR1a6KwAAAEBe\nVROgr7jiCt14442V7gYAAABQULDSHZCkhQsXqqmpSYcffnjRti0tLalyOBxWOBwuX8cAAADQ60Ui\nEUUikZLbO8YYU77u7DFjxgxt3Lgxq/66667T9ddfr6eeekr77LOPRo0apZUrV2q//fbL7qzjqIe6\nCwAAgL1UsczZYwE6n//7v//TiSeeqD59+kiSNmzYoGHDhmn58uXaf//9M9oSoAEAAFBuVR+gOxo1\napT+/Oc/a+DAgVnrCNAAAAAot2KZs2oeIkxyHKfSXQAAAADyqroR6EIYgQYAAEC51dwINAAAAFDN\nCNAAAACABQI0AAAAYIEADQAAAFggQAMAAAAWCNAAAACABQI0AAAAYIEADQAAAFggQAMAAAAWCNAA\nAACABQI0AAAAYIEADQAAAFggQAMAAAAWCNAAAACABQI0AAAAYIEADQAAAFggQAMAAAAWCNAAAACA\nBQI0AAAAYIEADQAAAFggQAMAAAAWCNAAAACABQI0AAAAYIEADQAAAFggQAMAAAAWCNAAAACABQI0\nAAAAYIEADQAAAFggQAMAAAAWCNAAAACABQI0AAAAYIEADQAAAFggQAMAAAAWCNAAAACABQI0AAAA\nYKFqAvTPf/5zNTc3a/z48Zo/f36luwMAAADkFKx0ByTp2Wef1aJFi/T6668rFApp8+bNle4SAAAA\nkFNVBOg777xTV199tUKhkCRp0KBBedu2tLSkyuFwWOFwuMy9AwAAQG8WiUQUiURKbu8YY0z5ulOa\nyZMn68wzz9TixYvV0NCgm2++WUcccURWO8dxVAXdBQAAQC9WLHP22Aj0jBkztHHjxqz66667TrFY\nTFu3btWyZcu0YsUKnXvuuXr//fd7qmsAAABAyXosQD/99NN51915552aPXu2JGnq1KlyXVdbtmzR\nfvvt11PdAwAAAEpSFW/hmDVrlpYuXSpJWr16taLRKOEZAAAAVakqHiK8+OKLdfHFF2vChAmqq6vT\nAw88UOkuAQAAADlVxUOEpeIhQgAAAJRbscxZFbdwAAAAALWCAA0AAABYIEADAAAAFgjQAAAAgAUC\nNAAAAGCBAA0AAABYIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYIEADAAAAFgjQAAAAgAUCNAAA\nAGCBAA0AAABYIEADAAAAFgjQAAAAgAUCNAAAAGCBAA0AAABYIEADAAAAFgjQAAAAgIVgpTsAAAAA\nlIsxku9LsZjkeaXNiyFAAwAA1Cjfzwx/+cqF6nr73PMk15UCASkYLG1eDAEaAABUhY4jhbkCYb66\nzgTI3rDemNwBsFi51CDZ2Xl9fXn3bzt3HLtzsVh7AjQAAN0kGQA7hr18IdC2XbFgabtcjn125Ri+\nv2ekMFfwy7VcSmjsrvV1deXdf2faujzNVhEEaPRqkUhE4XC40t1AL9IbzyljSgt25ZxKCZfdEUDL\nfcxcAbBY+Gtri6hfv3DBNp0JkoXaJEcHO7ufrhy7WBvbkUJk640/p6oNARq9Gj9Eepd8I3vpUylt\nbNp1bPvooxG9+2640/ssNZD15GRM4dBX7ilfYMy1Pj34dXYf3dGPQpNtAGxpiailJVyW/2ewd+J3\nX/kRoIFukCs4dXddpfdZjnBqG2Kl0gJMqWGwM+0+/FB6+eXS91lX17Uw1hOT6zLqBwA2CNC9TPpX\nsR3DSfpyoXU2bat93erV0gsv2H3eztQ5Tu4wli+gldq23HXJ0byu7LMnQ2w13OvX0hKfAAB7L8cY\nYyrdiVI5DJEAAACgBxSKyDU1Al1DWR8AAAC9VBV8IQoAAADUDgI0AAAAYIEADQAAAFggQKOmLF68\nWIceeqgOPvhg/exnP8tav3XrVp111lmaOHGipk+frjfffFOS9M4772jy5MmpqX///lqwYIEk6aqr\nrlJzc7MmTpyo2bNna/v27T36mVBZ5Tinkm655Ra5rqtPP/20Rz4LKq9c59PPf/5zNTc3a/z48Zo/\nf36PfR5UXjnOqeXLl2vatGmaPHmypk6dqhUrVvToZ+oVDFAjYrGYGT16tPnggw9MNBo1EydONG+9\n9VZGmyuvvNL85Cc/McYY8/bbb5sTTzwxaz+e55nBgwebDz/80BhjzFNPPWU8zzPGGDN//nwzf/78\nMn8SVItynVPGGPPhhx+av/u7vzMHHnig2bJlS3k/CKpCuc6npUuXmpNOOslEo1FjjDGbNm0q8ydB\ntSjXOXXccceZxYsXG2OMeeKJJ0w4HC7zJ+l9GIFGzVi+fLnGjBmjAw88UKFQSOedd54WLlyY0eYv\nf/mLjj/+eEnSIYccorVr12rz5s0ZbZ555hmNHj1aw4cPlyTNmDFDbuIFw9OnT9eGDRt64NOgGpTr\nnJKkK664QjfeeGP5PwSqRrnOpzvvvFNXX321QqGQJGnQoEE98GlQDcp1Tg0ZMiT1beu2bds0bNiw\nHvg0vQsBGjWjtbU1I6A0NTWptbU1o83EiRP16KOPSor/4Fm3bl1WIH744Yd1/vnn5zzGPffco9NO\nO62be45qVa5zauHChWpqatLhhx9ext6j2pTrfHr33Xf1/PPP68gjj1Q4HNbKlSvL+ClQTcp1Tv30\npz/Vv/zLv2jEiBG66qqrdMMNN5TxU/ROBGjUjFL+kM4Pf/hDbdu2TZMnT9Ydd9yhyZMnKxAIpNZH\no1E99thjOuecc7K2ve6661RXV5c3XKP3Kcc5tXPnTl1//fW69tprU20M77DfK5TrZ1QsFtPWrVu1\nbNky3XTTTTr33HPL0n9Un3KdU/PmzdOCBQv04Ycf6tZbb9XFF19clv73ZjX1h1Swdxs2bJjWr1+f\nWl6/fr2ampoy2vTr10/33HNPannUqFE66KCDUstPPvmkpkyZkvUV6H333acnnnhCS5YsKVPvUY3K\ncU699957Wrt2rSZOnChJ2rBhg6ZMmaLly5dr//33L+fHQYWV62dUU1OTZs+eLUmaOnWqXNfVli1b\ntN9++5Xro6BKlOucWr58uZ555hlJ0pw5c3TJJZeU6yP0WoxAo2YcccQRevfdd7V27VpFo1H97ne/\n0xlnnJHRZvv27YpGo5Kku+++W8cdd5z69u2bWv/QQw9p7ty5GdssXrxYN910kxYuXKiGhobyfxBU\njXKcUxMmTNDHH3+sDz74QB988IGampq0atUqwvNeoFw/o2bNmqWlS5dKklavXq1oNEp43kuU65wa\nM2aMnnvuOUnS0qVLNXbs2DJ/kl6o0k8xAjaeeOIJM3bsWDN69Ghz/fXXG2OMueuuu8xdd91ljDHm\nxRdfNGPHjjWHHHKIOfvss822bdtS2+7YscPst99+5rPPPsvY55gxY8yIESPMpEmTzKRJk8yll17a\ncx8IFVeOcyrdqFGjeAvHXqQc51M0GjUXXHCBGT9+vPmbv/kb8+yzz/bY50HlleOcWrFihZk2bZqZ\nOHGiOfLII82qVat67gP1Eo4x3JwHAAAAlIpbOAAAAAALBGgAAADAAgEaAAAAsECABgAAACwQoAEA\nAAALBGgAAADAAgEaAHqx5cuX66tf/apOOOEEvf/++5XuDgD0CvwpbwDoxaZNm6YRI0Zozpw5GX/e\nFwDQeYxAA0Avt2zZMh111FGV7gYA9BoEaADoxTzP06ZNmzR48OBKdwUAeg0CNAD0Yq+//romTJhQ\n6W4AQK9CgAaAXuyll17SUUcdpZUrV+q5557TjTfeWOkuAUDNI0ADQC/z/vvvKxwOS5KefvppnX76\n6Vq5cqWmT5+uTz75RDt27KhsBwGgxvEWDgDoZQYMGKAzzjhDv/rVr3TWWWdp9OjRGj16tDzPUywW\nU9++fSvdRQCoaY4xxlS6EwCA8vvtb3+rk08+Wf3791coFKp0dwCgZhGgAWAvcP/99+v555+X67q6\n6667FAgEKt0lAKhZBGgAAADAAg8RAgAAABYI0AAAAIAFAjQAAABggQANAAAAWCBAAwAAABYI0AAA\nAIAFAjQAAABggQANAAAAWCBAAwAAABb+H8H5L1CEl2uCAAAAAElFTkSuQmCC\n"
      }
     ],
     "prompt_number": 8
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 8
    }
   ],
   "metadata": {}
  }
 ]
}