{ "metadata": { "name": "", "signature": "sha256:cb937b8f9d161aaaa0e521b57992dac0536c472afcd14ad83352bb92be77dc29" }, "nbformat": 3, "nbformat_minor": 0, "worksheets": [ { "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Integrating with SciPy and QuTiP\n", "\n", "By Thomas P Ogden (<t.p.ogden@durham.ac.uk>)\n", "\n", "_Here we introduce the Scipy ODE integration pack and the QuTiP solver for the specific problem of solving the Schr\u00f6dinger equation (for a closed quantum system) and Lindblad master equation (for an open quantum system)._" ] }, { "cell_type": "code", "collapsed": false, "input": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "%matplotlib inline" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 1 }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Integrating with SciPy\n", "\n", "Now we've gone through writing a bunch of ODE integrators, we're at the point where we find out why we'll (almost) never need to use them. Because like any good numerical library, SciPy has advanced solvers already available to us in its `scipy.integrate` toolbox. For example, the [`odeint`](odeintdoc) method allows us to:\n", "\n", "> Solve a system of ordinary differential equations using lsoda from the FORTRAN library odepack.\n", "\n", "[odeintdoc]:http://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.odeint.html\n", "\n", "These solvers use highly optimised methods with adaptive stepsizes, and can solve stiff or non-stiff problems. We always want to use them where we can. Let's test `odeint` with our `exp` ODE in the same way we did with the methods we wrote." ] }, { "cell_type": "code", "collapsed": false, "input": [ "def exp(y, t, args):\n", " \"\"\" An exponential function described as a first-order ODE. \"\"\"\n", " \n", " dydt = args['a']*y\n", " return dydt" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 6 }, { "cell_type": "code", "collapsed": false, "input": [ "y_0 = np.array([1.]) # Initial condition\n", "\n", "t = np.linspace(0., 5., 100+1) # Array of time steps\n", "\n", "solve_args = {'a':1}" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 24 }, { "cell_type": "code", "collapsed": false, "input": [ "from scipy.integrate import odeint\n", "\n", "y = odeint(exp, y_0, t, args=(solve_args,))\n", "\n", "plt.plot(t,y,'r', label='odeint')\n", "plt.plot(t, np.exp(solve_args['a']*t), 'k--', label='Known')\n", "\n", "plt.legend(loc=2)" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "pyout", "prompt_number": 25, "text": [ "<matplotlib.legend.Legend at 0x109f02b50>" ] }, { "metadata": {}, "output_type": "display_data", "png": "iVBORw0KGgoAAAANSUhEUgAAAaYAAAESCAYAAAC2KnFgAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XlcVXX+x/EXCIgIKriAAoqCIi4kKOqY6NVBxywpWiYc\nG0zUGpvKZeY3zUxNYrlMU5Oa47SYS2aZNVNpLpWUiGaJmqmZJSbY7RomLuyLwP39cY3RUuOynQu8\nn4/HeQjnnnPP+54H3U/fc77n+3WyWq1WREREHISz0QFEREQupcIkIiIORYVJREQcigqTiIg4FBUm\nERFxKCpMIiLiUFSYRETEoagwiYiIQ6lSYUpMTMTX15e+fftetr6goICJEycSERFBr1692L17NwCp\nqalERkYSHh7OkiVLaj+1iIg0Wk5VGflhx44deHp6kpCQwKFDhyrXT5w4keHDh5OYmEhZWRkFBQV4\nenoSGhpKcnIy/v7+REVFsXbtWsLCwur0g4iISONQpRZTdHQ03t7el63Lyclhx44dJCYmAuDi4kLr\n1q1JS0sjJCSEoKAgXF1diY+PZ/369bWfXEREGqVq32PKyMigffv23H333fTp04epU6dSVFSExWIh\nMDCwcruAgAAsFkuthBURkcbPpbo7lpWVsWfPHh555BGeffZZ7r33Xt544w1atmxp1/uEhITw9ddf\nVzeGiIg0EMHBwRw7duxnt6t2iykgIIC2bdsybtw4WrRowfjx49myZQv+/v6YzebK7cxmMwEBAVd9\nn6+//hqr1aqlisvs2bMNz9CQFp0vnS+dL8dZqtoIqXZh8vPzIyQkhN27d1NRUcGmTZuIiYkhKiqK\n9PR0MjMzKS0tZd26dcTGxlb3MCIi0sRUqTCNHz+eIUOGcPToUQIDA1m5ciUAL730EtOnT6dHjx5Y\nLBbi4+Np1qwZK1asIC4ujv79+5OYmKgeeSIiUmVV6i5epwGcnDA4QoOSkpKCyWQyOkaDofNlH50v\n++h82aeq3/cqTCIiUi+q+n2vIYlERMShVLu7eF3z8fHh3LlzRsdotLy9vTl79qzRMUREfsJhL+Xp\nEl/d0vkVkfqmS3kiItIgqTCJiIhDUWESERGHosJUz7799lucnZ355ptvqrT9K6+8Qr9+/eo4lYiI\n41BhcnATJkzgs88+q/L2q1atonv37nWYSESkbqkwiYiIQ1FhqgVnzpwhISGBjh070rFjR+6+++7K\nZ7CysrKIjY3Fx8eHAQMGsGvXrp/sv2zZMvr06UOrVq2IjIxk69atla/9uAVkMpl4+OGHSUhIoH37\n9oSEhLBhwwYAPv74Y6ZNm8bx48fx8vLCy8uL1NTUOv70IiK1y2EfsP1ZM2aAHZe4rqlfP1i0qNq7\nT5gwARcXF7788kusVivx8fH89re/ZePGjUyYMAFPT09OnDhBYWEht9xyy2X7Llu2jDlz5rBs2TKi\no6N57bXXiI2N5fPPPyc4OPiKx/vXv/7F8uXLWbp0KfPmzWPixIl89913/OIXv+C5555j7ty5pKen\nV/vziIjUutWrq7ypWkw1dPLkSd5//31mzZpF69atadOmDbNmzWLz5s2cPHmSbdu2MXPmTLy8vPD1\n9WXy5MmX7b948WKmT5/ODTfcgKenJ1OmTKFPnz689tprVz3myJEjuf322/Hy8uLBBx8kJyenshDp\noVkRcTSZO3dy+8SJVd6+4baYatDCqU0/TIo4cODAynWDBg0CbD3wAKKioipfu3Q7sE1R//jjjzN/\n/vzKdeXl5Zw8efKKx3Nycrrs/Tp16gRAXl5eTT6GiEidefWxx/ivHdurxVRDgYGBAOzevbty3Sef\nfAL8r2ikpaVVvnbpdgBBQUEkJSVx7ty5yiU3N5elS5dWK4+zs7NaTSLiMKwVFbySmsr1Xl5V3keF\nqYY6derE6NGjWbRoETk5OZw7d45FixYxduxYAgICMJlMPPPMM+Tl5XHq1ClWrVp12f4zZ87kmWee\nYdOmTRQUFFBUVMTOnTv56quvrnrMaxUePz8/vvnmG06dOlVbH1FEpNrK9+/nNyUlTL/zzirvo8JU\nC9asWUObNm0IDQ0lLCyMDh06sPrijb5XX32VsrIyunTpwo033sj06dNxcnKq3HfKlCk89thj/O1v\nf6NTp0506dKFefPmUVZWBtgu3V26/Q/rrmbkyJHExsbSs2dPvL292bFjRx18YhGRqnF5/XUednHh\njgULqryPRhdvonR+RaTOVVRAUBCEh8PGjRpdXEREDLZzJ5jNMGGCXbupMImISN145RVo2RJiY+3a\nTYVJRERqXUVxMdbXX4dbbrEVJztUqTAlJibi6+tL3759f/JaeXk5ERERjBs3rnJdamoqkZGRhIeH\ns2TJErsCiYhIw7fx8ccJOX+eYyaT3ftWqTBNmjSJd99994qvLV68mF69elX2FCsvLycxMZE333yT\nffv2sXz5co4cOWJ3MBERabheefllcp2c6PKb39i9b5UKU3R0NN7e3j9Z/+2337J582amTJlS2dMi\nLS2NkJAQgoKCcHV1JT4+nvXr19sdTEREGqacb75hg9nMnX364OrhYff+NbrHNHPmTJ588kmcnf/3\nNhaLpXI0BICAgAAsFktNDiMiIg3I6488QjEw8Q9/qNb+1R4rb+PGjXTo0IGIiAhSUlKq+zYAJCUl\nVf5sMpkwVeOapIiIOIZVb79NmJsb+YGBl32/V1W1C9OuXbvYsGEDmzdvpri4mNzcXBISErjvvvsq\nBzYF2yCnAQEB13yv6gQXERHHU3jwIKV5eUwcM4YRI0cyYuTIytfmzJlTpfeo9qW8+fPnYzabycjI\n4LXXXmPkyJGsXr2aAQMGkJ6eTmZmJqWlpaxbt45YO/uwNyQmk4l58+ZV/n7u3DmGDRvG8OHDycnJ\nMTCZiEj983jjDfY4O/PH55+v9ntUqTCNHz+eIUOGcPToUQIDA1m5cuVPtvmhV56LiwsrVqwgLi6O\n/v37k5iYSFhYWLUDOrpLx7Izm81ER0fj6+vL1q1bad26tcHpRETqUUWFbULAmBiade5c/fexGuxq\nERwgWpWYTCbrvHnzrIcOHbL6+/tbH3jggcrXVq5caQ0JCbGuXLnS2rNnT6u3t7f13nvvtZaXl1du\nc+DAAeuIESOs3t7e1m7dulnnzp1b+foDDzxgveeeeyq3jY6Otnbp0qXy9yeeeMI6duxYq9Vqtc6e\nPdsaExNjXbBggTUoKMjaoUMH6+zZs6+au6GcXxFpQD780GoFq/WVV674clW/dzTyQy3YtWsXw4YN\n4/777+eZZ5657LUTJ06we/dutm3bxsaNG1m9enXl7LQ5OTmMGjWK4cOHc+rUKTZt2sQLL7zA008/\nDcCoUaNITk4GID8/n/379wNUzla7detWRo0aVXms1NRUioqK2LdvH88//zyPPfYYu3btqvPPLyIC\nwEsvQatWttEeaqDhzmALV+29d7VegvZuXxVWq5WdO3fSokULxo8ff8XX582bh4+PD35+fgwdOpS9\ne/fym9/8hk2bNpGTk8Of/vQnXF1d6dmzJxMmTODFF1/kj3/8I8OHD6+8j/fFF18wcOBAunfvztat\nW+ncuTO7du1i4cKFlcdq3bo1s2fPxtnZmVtuuYWQkBD27dvHkCFDqv35RESqJD8f/vMfGD8eqvHs\n0qXUYqohJycn/vCHP/CrX/2K6Ohojh49etnrQUFB+Pj4VP7u7+9Pfn4+YLsn1bt3b1q0aFH5+qBB\ngyp7NbZq1YqoqCiSk5P54IMPGD16NDExMWzdupWdO3fi5eVFnz59Kve97rrrLnumzN/fX1Oui0i9\nePMvf+EPBQUU2jEh4NU06BaTvS2dmj5vdTXNmjVj1apV/P73v2fYsGFs3br1iuMK/lhgYCBffPEF\nhYWFeFz8P4xPPvmEzpfcNPyhEH355ZesWrWKoKAg7rnnHnr06EFMTEydfB4REXstWbMGs4sLT13S\nPby61GKqBdaLwzEtXbqUhIQETCYTe/bsueq2P2x/44030qpVK/75z39y4cIFvvrqK1577TUmT55c\nuX1MTAzvvvsuWVlZREZG4uPjQ9euXXn++eerVJismgxQROpY+tatpJw/z+SRI3FyrnlZUWGqBZdO\ndf6Pf/yD6dOnExMTg9lsvuK06D+sa926Ne+//z7btm3D19eXMWPGkJiYyMyZMyu3Hzx4MFarlZGX\n/F9ITEwMeXl5lxWmK03B/uNsIiJ1Yfmjj9IMmGjH9OnXoqnVmyidXxGpDRcKCwn08mJQhw6s/+67\na26rqdVFRKTO7V68mO8rKpgyZUqtvadaTE2Uzq+I1Ipx48hMSyPgxAlc3N2vualaTCIiUre+/RY2\nbyZo6tSfLUr2UGESEZHqWbXKNj5eYmKtvq0u5TVROr8iUiMVFRAcbFsuDp32c3QpT0RE6s6HH0Jm\nJtRip4cfqMXUROn8ikhN/KtfPwZlZhKVlQVVvL+kFpOIiNSJk59+yowDB3i9e/cqFyV7OOxYed7e\n3hq1oA55e3sbHUFEGqgX//hHyoF7//73Onl/h72UJyIijqesuJigli3p4+3Nu9nZdu2rS3kiIlLr\n3klKwlJRwbR77qmzY6jFJCIiVXZj+/YcPHeOjPx8ux+qrer3vcPeYxIREQeTns7L2dkcreWRHn5M\nl/JERKRqnn8eHxcXBs+ZU6eHqVJhSkxMxNfX97JZWc1mMyNGjKB3796YTCZWrVpV+VpqaiqRkZGE\nh4ezZMmSWg8tIiL1rKgIVq6EW26Bjh3r9FBVuse0Y8cOPD09SUhI4NChQwBkZWWRlZVFv379yM7O\npk+fPqSkpNC9e3dCQ0NJTk7G39+fqKgo1q5dS1hY2JUD6B6TiIjjW7nSNibeBx9ANadPr9V7TNHR\n0WRmZl62zs/PDz8/PwDatWtHVFQUFouFc+fOERISQlBQEADx8fGsX7/+qoVJREQcnNUKzzwDffrA\niBF1frhaucd07NgxDh8+zODBg7FYLAQGBla+FhAQgMViqY3DiIiIAfYvX847n31Gxf33Qz0MfFDj\nXnn5+fnEx8ezcOFCWrZsWa3RGpKSkip/NplMmEymmsYSEZFaMnf2bFKcnDDfdhseduyXkpJCSkqK\n3cerUWG6cOECt912G3fddRc333wzAP7+/pjN5sptzGYzAQEB13yfSwuTiIg4jhMffcTbJ0/yp8GD\n8WjXzq59f9zQmFPF3nzVvpRntVqZPHkyvXv3ZsaMGZXrBwwYQHp6OpmZmZSWlrJu3TpiY2OrexgR\nETHQ0pkzcQLuW7So3o5ZpV5548ePZ/v27WRnZ+Pr68tjjz1Gjx49iI6OJjw8vPLy3YIFCxgzZgzb\nt29nxowZlJWVMXXqVB588MGrB1CvPBERh1Tw/fcE+Pkxyt+f1y+5ElZdVf2+15BEIiJyRS/+9rdM\nXbOGHUuXMvS++2r8fipMIiJSfVYrxb16sbmkhLhjx3Byrnknbo2VJyIi1ffee7h/+SW3vvQS1EJR\nsodaTCIi8lMxMfDll3D8OLi51cpbaj4mERGpnv37bUMPTZ9ea0XJHmoxiYjI5SZMgHfeAbMZWreu\ntbdVi0lEROz2zccfs2ztWoonTarVomQPFSYREam0+IEHmGa18n18vGEZdClPREQAOH/iBIFBQdwc\nFMSajIxaf39dyhMREbu8MG0a+cAf5s83NIdaTCIiQkluLt28vQlr3Zrks2fr5Bh6wFZERKrszVmz\nOFlRwUt//rPRUdRiEhFp8srKsPboQUrz5pgOH66V4YeuRC0mERGpmtdfxykjgxFvv13vww9diVpM\nIiJNWUUFXHcdWK1w8GCdFia1mERE5Oe98w58/jmsWeMQrSVQi0lEpOmyWmHwYDh9Go4eBZe6bavo\nOSYREbmm1EWLeDwtjYIZM+q8KNlDLSYRkSZqhLc3X+XmcvzMGdzbtKnz4+kek4iIXFXKokWknD/P\n4ltvrZeiZA+1mEREmhqrFZO3N0fz8vj69Gla+PjUy2Fr9R5TYmIivr6+9O3b97L1qampREZGEh4e\nzpIlS352vYiIGG/bwoVsz8nhL3Fx9VaU7FGlFtOOHTvw9PQkISGBQ4cOAVBeXk5oaCjJycn4+/sT\nFRXF2rVr6dGjxxXXh4WFXTmAWkwiIvXHauUvnTuz+uRJvq6ne0s/qNUWU3R0NN7e3petS0tLIyQk\nhKCgIFxdXYmPj2f9+vXs2bPniutFRMQBbNvGgm+/5cDcuQ53b+kH1e4ubrFYCAwMrPw9ICAAi8Vy\n1fUiImIwqxVmz4ZOnWg3c6bRaa6q2r3ynJycai1EUlJS5c8mkwmTyVRr7y0iIhd98AHs3AlLloC7\ne50fLiUlhZSUFLv3q3Zh8vf3x2w2V/5uNpsJCAi46vprubQwiYhIHbBa4eGHISAApkypl0P+uKEx\nZ86cKu1X7Ut5AwYMID09nczMTEpLS1m3bh2xsbFXXS8iIsYpf/NNSEuzXcqrh9ZSTVSpMI0fP54h\nQ4Zw9OhRAgMDWblyJS4uLqxYsYK4uDj69+9PYmIiYWFhV10vIiLGKC8tZeBdd/FUu3Zw991Gx/lZ\nesBWRKSRW33PPUxctox1M2bw64ULDctR1e97FSYRkUasND+fUG9vvF1d2Zubi7OBg7VqrDwREeGF\nSZPILCvj2dmzDS1K9lCLSUSkkSr4/nuCO3akp5cX286excngiQDVYhIRaeKK/vlPfllRwe+feMLw\nomQPtZhERBqj776D7t1hzBj4z3+MTgNoBlsRkaZt9mwoKYG//93oJHZTYRIRaWw+/xyWL4ff/x5C\nQoxOYzddyhMRaWzGjoWPP4Zjx6BtW6PTVFLnBxGRJmjv0qW4b9lCn6eecqiiZA+1mEREGony0lIG\ntGlD/oULfJmTQzMPD6MjXUYtJhGRJmbFlCl8VlTEaw8+6HBFyR5qMYmINALnMjLoERxMr1atSHGA\nh2mvRN3FRUSakDm33spZq5XFzz/vkEXJHmoxiYg0cLlpaQQNGsSve/XiucOHjY5zVRpdXESkKbBa\nYfRostLScN2zh7Y9ehid6KrU+UFEpClYvx6Sk/F75hlw4KJkD7WYREQaqsJC6NMHPDzgs8/Awae1\nUItJRKSxmz8fMjJg2zaHL0r2UItJRKQBOrd7N22GDsVp/HhYvdroOFWizg8iIo2UtaKCX7Zti1d+\nPustFujQwehIVaLnmEREGqlX77+fbefPM+b22xtMUbJHjQvTsmXLGDJkCP3792fGjBmV61NTU4mM\njCQ8PJwlS5bU9DAiIgKcP3GCWc89R1TLltzz0ktGx6kTNSpMZ8+eZf78+WzdupU9e/Zw9OhR3n//\nfcrLy0lMTOTNN99k3759LF++nCNHjtRWZhGRJuuvN95IttXKc88/TzM3N6Pj1IkaFaYWLVpgtVrJ\nycmhqKiIwsJC2rRpQ1paGiEhIQQFBeHq6kp8fDzr16+vrcwiIk3SnmXLePbwYR7o14/ICROMjlNn\nalyYnn32WYKCgvDz8+P6669n4MCBWCwWAgMDK7cLCAjAYrHUOKyISJNVXEzEU0+x1MeHuVu2GJ2m\nTtWo4/vp06eZNm0aX3zxBd7e3txxxx1s2rQJJycnu94nKSmp8meTyYTJZKpJLBGRxmfuXFyOHuW+\n994DPz+j01RJSkoKKSkpdu9Xo8KUlpbG4MGDCbk4p/wdd9xBamoqcXFxmM3myu3MZjMBAQFXfZ9L\nC5OIiPzIwYPwxBOQkACjRxudpsp+3NCYM2dOlfar0aW86Oho9u7dy9mzZykpKWHLli2MHj2aqKgo\n0tPTyczMpLS0lHXr1hEbG1uTQ4mINE1lZTB5Mvj4wNNPG52mXtSoxdSqVSseeeQR4uLiKCwsZMyY\nMYwYMQJnZ2dWrFhBXFwcZWVlTJ06lbCwsNrKLCLSZJQ+9RRue/fCunXQtq3RceqFRn4QEXFQX2zY\nwOibb2ZtdDTR27eDnffvHY2GJBIRacAuFBUxpF07MouK+PzgQXz79DE6Uo1pdHERkQbs7zfdxN7C\nQt6YNatRFCV7qMUkIuJgPlu3jqj4eO7o0oVXMzONjlNrdClPRKQBqigupr+3N1mlpRw+ehSf4GCj\nI9UaXcoTEWmAnB99lH8XF1M4d26jKkr2UItJRMRRfPABxMTA734Hzz5rdJpap0t5IiINyZkzEB4O\nXl7w6afg4WF0olqnS3kiIg2F1QpTp8Lp0/DOO42yKNlDM9iKiBgs6+mn4a23YN48iIw0Oo7hdClP\nRMRARzZuJGrcOP4ZGsq9X3wBzo23vaB7TCIiDq4wO5tBAQFklZbyWVoa/gMGGB2pTukek4iIg3vw\n+uv5vKSEd+fObfRFyR6Nt80oIuLA1kybxvKjR/nrkCH86uGHjY7jUHQpT0SknpUfPkx437609fLi\nw1OncHF3NzpSvdClPBERR1RQQLP4eFK9vSl9//0mU5TsocIkIlJffnhe6fBh2m7ZAv37G53IIakw\niYjUl2eegbVrbc8r/epXRqdxWLrHJCJSH1JTYeRIuOkmePPNRv280tVU9fu+6Z0ZEZF6Ztm7l4fH\njuVCt27w0ktNsijZQy0mEZE6VHz+PKaAAD4vKODTTZvoMXas0ZEMo155IiIGs1ZUMDUykt0FBbz5\npz816aJkjxq3JwsKCpg4cSIRERH06tWL3bt3A5CamkpkZCTh4eEsWbKkxkFFRBqaf9x4I2syMpgb\nE0PcE08YHafBqPGlvIkTJzJ8+HASExMpKyujoKAAT09PQkNDSU5Oxt/fn6ioKNauXUtYWNhPA+hS\nnog0Qu8lJXHDnDnc2bkzr2Zk4KT7SvXT+SEnJ4cdO3aQmJgIgIuLC61btyYtLY2QkBCCgoJwdXUl\nPj6e9evX1+RQIiINx6FDDHrqKWb6+rJi/34VJTvV6GxlZGTQvn177r77bvr06cPUqVMpKirCYrEQ\nGBhYuV1AQAAWi6XGYUVEHN7Jk3DjjbRp3Zp/7ttHCx8foxM1ODXq/FBWVsaePXt45JFHePbZZ7n3\n3nt54403aNmypV3vk5SUVPmzyWTCZDLVJJaIiDHy8uDGG+HcOdixA/z9jU5kqJSUFFJSUuzer0b3\nmLKysujbty+nT58GYMuWLaxevZrp06eTlJTEu+++C8CCBQtwdnbmoYce+mkA3WMSkcagrAxuvhne\ne882PfoNNxidyOHUyz0mPz8/QkJC2L17NxUVFWzatImYmBiioqJIT08nMzOT0tJS1q1bR2xsbE0O\nJSLisKwVFawaNYqSzZvh3/9WUaqhGj/H9NJLL5GQkEB2djZ9+/bliSeeoFmzZqxYsYK4uDjKysqY\nOnXqFXvkiYg0BvN/9SseSUmhdMwY7rnnHqPjNHga+UFEpAZeuOsu7n3lFX7brRurvvoKZxeNW3A1\nGitPRKSO/ff//o9pr7zC2PbtWX7okIpSLVGLSUSkGnYvXcqw++9ngJcXW48fx6NdO6MjObyqft+r\nMImI2CstjcKRI3m4eXP+lpaGT3Cw0YkaBBUmEZG68NlnMGIEeHvrWSU7qTCJiNS2L76A4cPB3d1W\nlIKCjE7UoKgwiYjUpmPHYNgwsFpts9F27250ogZHvfJERGrJ1x9+yKTrrqOwpASSk1WU6pgKk4jI\nNRz74ANMo0fzTlERJ5Ytg969jY7U6KkwiYhcRfrWrZhGj6a4ooIPX3+dsFtvNTpSk6DCJCJyBV9t\n2cLwMWMotVrZ9t//En777UZHajL0mLKIyI8dPMhTt95KudXKtrfeovfNNxudqElRrzwRkUt9/DGM\nHUtJy5ZYVq2iW0yM0YkaDXUXFxGxV3Iy3HILdOxo+7lLF6MTNSrqLi4iYo+33rLNPtutm+3hWRUl\nw6gwiUiTt27yZPJuvRUiIyElBfz8jI7UpKkwiUiTZa2o4JGhQ4lfsYKFPXrYLt/5+Bgdq8lTrzwR\naZIuFBbyu379WJGezuQePfjrgQO2MfDEcGoxiUiTc/7ECW4MDGRFejqPDhvGsiNHcFFRchgqTCLS\ntGRksGjAALadPcuKSZOYs307Ts76KnQk6i4uIk3Hxx/DzTdzobSUzxYsIGraNKMTNSn12l28vLyc\niIgIxo0bV7kuNTWVyMhIwsPDWbJkSW0cRkSk+lavtk3w16oVrrt3qyg5sFopTIsXL6ZXr144OTkB\ntkKVmJjIm2++yb59+1i+fDlHjhypjUOJiNinrAxmzoSJE2HIEPjkEwgNNTqVXEONC9O3337L5s2b\nmTJlSmUTLS0tjZCQEIKCgnB1dSU+Pp7169fXOKyIiD3OpKczuXNnzi5aBA8+CO+9B+3aGR1LfkaN\nC9PMmTN58skncb7k5qHFYiEwMLDy94CAACwWS00PJSJSZZ+++ipRvXqx5rvvSJs5ExYvBldXo2NJ\nFdToOaaNGzfSoUMHIiIiSElJqfb7JCUlVf5sMpkwmUw1iSUiTZjVamX53Xdz/+rVtG/WjNQXX2TQ\n5MlGx2qSUlJSqlUbatQr769//Ssvv/wyLi4uFBcXk5uby2233cZ9991HUlIS7777LgALFizA2dmZ\nhx566KcB1CtPRGpJRUEBUyIiWJmezigfH17ZuZP2YWFGx5KL6n108e3bt/PUU0/xzjvvUFZWRmho\nKB988AGdOnVi4MCBrF27lrAr/IGoMIlIrThyBO68kz8eOkTLYcN4dOtWmrm5GZ1KLlHV7/taHZLo\nh155Li4urFixgri4OMrKypg6deoVi5KISI1ZrbB8ua1zg6cnT27ejNMNNxidSmpAD9iKSMOVkwP3\n3gvr1sEvfwkvv2ybS0kckuZjEpFGbfvixXwSGgr/+Q/Mn2/rCq6i1ChodHERaVBKcnP5W0wMT+3Z\nQ0yLFry/cycMHmx0LKlFajGJSINx4PXXGejry5N79jC1Z0/ePH5cRakRUmESEcdXWsqiMWMYcOed\nnCot5Z2//Y3njxzBUzPNNkq6lCcijm3/frj7boIOHuTXQUE88/77tO3e3ehUUofUYhIRx1RYCA89\nBFFR8P333LJ+Pa9kZKgoNQFqMYmIw7G+9x5O06ZBRgYkJsKTT4KPj9GxpJ6oxSQiDuO7zz5jQlAQ\n88eMsQ0XGZh/AAAPpElEQVS4um2b7eFZFaUmRYVJRAxXVlzMorg4QiMi+M+JEziNGAEHDoAGdG6S\nVJhExFDbFy8mok0bZr79NkPbt+dwcjJ//fBDcHc3OpoYRIVJRIxx/Djcdhv/mDGDvLIy3vrzn9mU\nlUXIL39pdDIxmMbKE5H6lZNjG0Jo0SJwcSHr/vtp/dBDtNB9pEav3qe9qC4VJpGmwVpSgtNzz8Hj\nj8OZMzBxoq1AdepkdDSpJxrEVUQcQkVZGa89+CC9vLzImDEDrrsO9u6FVatUlOSKVJhEpE5YKyrY\n+Oij9G/VivFLluDWrBnnliyB5GTo39/oeOLA9ICtiNQuq5WDy5Zxz6xZ7C4oINjFhTXTpjH+mWdw\ndtFXjvw8/ZWISO2wWm1zIj3+OK127SLbxYVlCQlMfPZZXD08jE4nDYg6P4hIzVRUwKZNtk4Ne/ZA\nYCD8+c9UTJqEc4sWRqcTB1LV73u1mESkWi4UFvLarFn03bqVfsePQ1AQvPCCrbedm5tuYEu16W9H\nROySd/Iki2+9lZBWrUh4/nlW5eTAyy/D0aMwdSq4uRkdURo4tZhEpErOHDjA36dM4YW9e8kFhrZq\nxdIZMxj76KPQrJnR8aQRqVGLyWw2M2LECHr37o3JZGLVqlWVr6WmphIZGUl4eDhLliypaU4RMYLV\nCqmp8Otf0ywykhf37uWGzp35ZPlyduTkcNOcOTirKEktq1Hnh6ysLLKysujXrx/Z2dn06dOHlJQU\nunfvTmhoKMnJyfj7+xMVFcXatWsJCwv7aQB1fhBxPPn5sHYt/OtfcPAgeHvD5MnkTZqEV69eRqeT\nBqpeOj/4+fnh5+cHQLt27YiKisJisXDu3DlCQkIICgoCID4+nvXr11+xMImI4zjw+uu88NhjjDl+\nnHFFRbZRGl58EcaPBw8PvIwOKE1CrXV+OHbsGIcPH2bw4MFYLBYCAwMrXwsICMBisdTWoUSkFuV+\n+y3LEhIY7OlJvzvvZPnhw3wVFgY7d8L+/TB5Mug5JKlHtdL5IT8/n/j4eBYuXEjLli1xcnKya/+k\npKTKn00mEyZNDiZStyoqIDWVvU8+ybDNmykCejVvzsJbbiHhqafwCQ42OqE0AikpKaSkpNi9X40f\nsL1w4QI33XQTN9xwAzNmzADgk08+ISkpiXfffReABQsW4OzszEMPPfTTALrHJFJ/vvzS1rV7zRr4\n5htKvbz4U5cu/GbWLKImTsTJWU+QSN2pl2kvrFYrEydOpF27djz99NOV68vKyggNDeWDDz6gU6dO\nDBw4UJ0fRAxi3r2b1+fOZZLZjM+BA+DsDKNHQ0IC3HyzLtNJvamXzg8fffQRa9asITw8nIiICMDW\nOhozZgwrVqwgLi6OsrIypk6dqo4PIvXo5Kef8tbf/866995jR24uAEHBwdz29NO2jgwXOy2JOCKN\nlSfSWJw4AW+/zd8XLuSvJ05gxXbfaPzQodz50EN0HzXK6ITSxGkGW5HGzmqFQ4dg/Xp46y1bDzpg\nV3AwH3Tpwm3Tp9MrNtbgkCL/o8Ik0ggVnz/P9qVLeWfdOnKOH+flggJwcoLBgyEuznbPqEcPo2OK\nXJFGFxdpJEqOHOHF2bN5d/t2Pvz+ewoBD+CGjh2p+Oc/cY6NhY4djY4pUmvUYhJxNOfPw7ZtsHUr\nbN1KxbFj+AKtXVy4ISyMG267jREPPEALHx+jk4rYRZfyRBqIwuxsPnrxRbZt2MA9RUUEHToE5eXQ\nsiWYTDB6NKcHDaL9oEFGRxWpEV3KE3FUeXl8snw5G954g9RDh0jLy+MCtv8Yo3r2JOgvf7E9ZzRo\nUOXcRu0NDSxSv9RiEqlrFgt89JFt7LmPPoIDB/i/8nIWAQM8PRnWqxcjYmMZOnkynnq+SBoxXcoT\nMUDR2bPsf+MNPtm8mY8//ZQBubk8dPEBVzw8bL3nrr+e7Ouuw33oUDx9fY0NLFKPdClPpK6VlsLh\nw7B3L7s2buT3W7dyqKiI8osvB7m40LdnT5gyBa6/3jaFhKsrAO2MSy3i8FSYRKqg4PvvOfj222Tv\n38+48nLbw6wHD9qKE+Dj5UX7Fi34c2QkA4cPJ+qOO+jYr5/BqUUaJl3KE7lURQVkZJCblsaiZcs4\nePQoB7//nmMXLmAFOgBZ3t44RUZC//7/W7p1sz3oKiJXpXtMItdQXlrKiV27+ColhRs8POCLL2yX\n5b74AgoLKQZaAV1cXQlv357wHj3o94tfEDFuHIGDBml6CJFqUGESAdvDqkePwtGjPP7iixz6+mu+\nzM7maHExJRc3yQbaduoEYWHQt69t6d2bom7daNFeHbVFaosKkzQJ1ooKTh85wvGPP+brTz/l2JEj\nPNChAz5mM6SnQ3Z25ba9gAuurvT09qZn58707N2bXkOG0P+WW3Dr0MG4DyHSRKgwSaNgrajAeuYM\nzmazbVqHEycgIwMyMrh9+3bey80l/5LtnYBPfH0Z2Ls3BAfbBjS9uJR37kwzTYonYhgVJmkYcnNt\nD6CazWzatIl9hw7xjcWC+cwZvsnLw1xaygZg5KX7eHpC167Mu3CB0+7uBAcH061PH7oNGEDXoUNx\nb9PGoA8jIteiwiTGsVohNxfz/v1kHDrEd8eOcfLECb47eZKTp09zn48PQ/Lz4eRJyP9fe+cO4D+A\nr7MznVu0ILB1awI7dGDq6NH0HjQIunSBoCDw8VEPOJEGSIVJaldpKWfS07EcOcLpjAxOm818b7Hw\n/alTjGvXjkFWK2RlwalTtqW4mLuBly55i+ZAJxcXngoO5tbwcPD3h06dICAAAgI45+WFR7duNG/V\nypjPKCJ1SoVJrqyiAnJy4Nw5Mj//nPTDhzmblcWZU6c4c/o0Z8+dY1zbtox0c7N1HDhzxvZvTg73\nAc/+6O2aAf/y8uJ3XbuCr+//Fj8/9peUcNbNDb8ePegUHk6bLl3UzVqkCVNhaoSsFRUUnT0Lubl4\nlJXZ7s/k5NiW3Fy279tHyuefk5ObS05eHucLCjhfVMSkVq24y9kZzp2zbXvxfP8JePJHx/AE5nt7\n80DXrtCunW1p2xbatWNPQQEnyspo37kzHYKD6RAainfXrji7aAAREfl5GivPaGVlUFhoWwoKoKCA\nb77+mqPp6RScP09BTg75OTkU5OcT1a4dQ9u2td1vycur/Pe5Y8dY8t135JWVkVdRQZ7VSjkwB3j0\nCof8EHgMW3Fp7exMa1dXvN3csHp7Q+/e4O1tW3x8wNubqcXFxJaX4xMYiE+XLvh064abp+dVP1LU\nxUVEpC7VWWFKTU1lxowZlJWVMXXqVB544IG6OtQ1lRUXU5KbS0luLsW5uZTk5+Pl6ko7Dw8oLoaS\nEtu/xcUcTk9n31dfUVxYSFFhIcVFRRQVFXG9ry+j/PygqOiyZfXx4/zbbKaovJzCHxarlRnOzswu\nL/9JltXA366Q8S/A0ObNbb3NPD3Bywu8vPD28qJnaSmtPDzw8vDAy9OT7LIyYoYNg4gI23atW9uW\nVq34q7s7f/P2xsXdvUrnpvvFpTFLSUnBZDIZHaPB0Pmyj85X3aiTwlReXk5iYiLJycn4+/sTFRVF\nTEwMYWFhV9x+3qhRXLhwoXIZ1rEjNwUEwIULtqW0FEpL+W9mJsszMyktL6e0vJyS8nJKKyr4batW\nzPL0tG1XUlL57xPFxfz5Cs3G/wP+cYUcm4CHrrD+T87OjPLyghYt/rd4eOBmtdLa3Z2Obm54uLnh\n4e5Oi+bN6d+jB4SH26Y58PCwzUTasiV35eUxLC8PTx8fWvr44Nm+PS3bt7dNfXCFYnLnxeVSSUlJ\nDElKuuJ5bH7FtU2bvjjso/NlH52vulEnhSktLY2QkBCCgoIAiI+PZ/369VctTI8kJwO2G+mugIub\nGzd5edmmCHBzsy2urhQWFJBdVIRbs2a4OTvT0s0NNxcXfAICoHt3aN7ctri5QfPmRJ8+zeNmM27N\nm9Pc3Z3m7u64e3jQt2tX2/Az7u62pXlzcHdnSkkJt5eW0tzTkxbe3ri3bo17mzZXvYcSf3GpqqCL\ni4iIXF2dFCaLxUJgYGDl7wEBAezevfuq25fk5eHi7v6zN9F/e3GpqiEXl6ryubiIiIhx6qQwOdnx\n8GNwcDDNvbzqIkajNWfOHKMjNCg6X/bR+bKPzlfVBQcHV2m7OilM/v7+mM3myt/NZjMBAQFX3PbY\nsWN1EUFERBqoOnnaccCAAaSnp5OZmUlpaSnr1q0jNja2Lg4lIiKNTJ20mFxcXFixYgVxcXGV3cWv\n1vFBRETkUoaP/CAiInIpDVwmIiIOxbDClJqaSmRkJOHh4SxZssSoGA1GYmIivr6+9O3b1+goDs9s\nNjNixAh69+6NyWRi1apVRkdyaMXFxQwaNIh+/foxePBgFi5caHSkBqG8vJyIiAjGjRtndBSHFxQU\nRHh4OBEREQwcOPBntzfkUl55eTmhoaGXjQyxdu1a3Ye6hh07duDp6UlCQgKHDh0yOo5Dy8rKIisr\ni379+pGdnU2fPn3Ytm2b/r6uobCwEA8PD0pKSujfvz9vv/02ISEhRsdyaE8//TT79u0jLy+PDRs2\nGB3HoXXt2pV9+/bh41O1J0UNaTFdOjKEq6tr5cgQcnXR0dF4e3sbHaNB8PPzo1+/fgC0a9eOqKgo\nTp48aXAqx+Zxccr5/Px8ysrKaN5cA1xdy7fffsvmzZuZMmWKZkeoInvOkyGF6UojQ1gsFiOiSCN3\n7NgxDh8+zODBg42O4tAqKiq47rrr8PX15f7777/sv0/5qZkzZ/Lkk0/irPnFqsTJyYmRI0cSERHB\nsmXLfnZ7Q86qPSNDiFRXfn4+8fHxLFy4kJYtWxodx6E5Oztz4MABjh07xr///W/2799vdCSHtXHj\nRjp06EBERIRaS1X00UcfceDAAV599VXmz5/Pjh07rrm9IYXJnpEhRKrjwoUL3Hbbbdx1113cfPPN\nRsdpMIKCghg7dizbt283OorD2rVrFxs2bKBr166MHz+eDz/8kISEBKNjObSOHTsCEBYWRlxcHGlp\nadfc3pDCpJEhpC5ZrVYmT55M7969mTFjhtFxHF52djbnz58H4MyZM2zZskW9P69h/vz5mM1mMjIy\neO211xg5ciSrV682OpbDKiwsJC8vD4DTp0+zefPmn/37MmQGW40MYb/x48ezfft2zpw5Q2BgII89\n9hiTJk0yOpZD+uijj1izZk1l91SABQsWMGbMGIOTOabvvvuOiRMnUl5ejp+fH7NmzeKXv/yl0bEa\nDN2auLZTp04RFxcHQNu2bZk5cyajR4++5j4a+UFERByKupSIiIhDUWESERGHosIkIiIORYVJREQc\nigqTiIg4FBUmERFxKCpMIiLiUFSYRETEofw/FNjvn+ICb7QAAAAASUVORK5CYII=\n", "text": [ "<matplotlib.figure.Figure at 0x109d92e10>" ] } ], "prompt_number": 25 }, { "cell_type": "markdown", "metadata": {}, "source": [ "So we see that the `odeint` integrator solved the system to a high level of accuracy \u2013 the two lines are indistinguishable at this scale. You might notice the basic methods we wrote had a similar interface to `odeint` \u2014 that was deliberate. But `odeint` can take a lot more arguments to optimise it for particular problems. I recommend reading through the documentation to understand the sorts of things you can do with it.\n", "\n", "### A Delayed Motivation\n", "\n", "So why did we go through Parts 1 to 3, if a better integration toolbox already exist? My hope is that by examining simple solvers yourself, you now:\n", "\n", "1. understand the principles of what `scipy.integrate` is doing. It's no longer a black box, it's just a more advanced version of what we just made;\n", "2. understand the principles of **accuracy**, **stability** and **computational cost** in ODE integrators.\n", "\n", "You many also find specific problems where the SciPy method isn't appropriate, so you'll need to write your own." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "##\u00a0Integrating the Schr\u00f6dinger Equation and Lindblad Master Equation with QuTiP\n", "\n", "Note: All of the below uses **QuTip v2.2.0**.\n", "\n", "In atom optics and quantum optics, the most common ODEs for us to be solving are those related to the evolution of a quantum system. **QuTiP** (The Quantum Toolbox in Python) is an excellent library for dealing with quantum systems, so we'll look at a basic example of solving evolution ODEs.\n", "\n", "We'll take a two-level atom with ground state $|g\\rangle$ and excited state $|e\\rangle$ and introduce an interaction with classical electromagnetic field.\n", "\n", "We'll write the vector representation of the system as\n", "\n", "$$\n", "\\Psi = \n", "\\left[\n", "\\begin{array}{c}\n", "c_g \\\\\n", "c_e\n", "\\end{array}\n", "\\right].\n", "$$\n", "\n", "Applying the electric dipole and rotating-wave approximations, the Hamiltonian is given by\n", "\n", "$$\n", "\\mathcal{H} = \n", "\\hbar\n", "\\left[\n", "\\begin{array}{cc}\n", "0 & \\Omega^*/2 \\\\\n", "\\Omega/2 & -\\Delta\n", "\\end{array}\n", "\\right]\n", "$$\n", "\n", "where $\\Omega$ is the Rabi frequency and $\\Delta$ the laser detuning.\n", "\n", "### Unitary Evolution with the Schr\u00f6dinger Equation\n", "\n", "To begin with we'll neglect spontaneous decay of the excited state, such that we have a closed quantum system. The dynamics of such a system are described by the time-dependent Schr\u00f6dinger equation,\n", "\n", "$$ \\mathrm{i} \\hbar \\frac{\\partial \\Psi}{\\partial t} = \\mathcal{H} \\Psi $$.\n", "\n", "\n", "\n", "\n", "Let's see how to solve this ODE \u2014 the time-dependent Schr\u00f6dinger equation \u2014 with QuTiP." ] }, { "cell_type": "code", "collapsed": false, "input": [ "import qutip as qu" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 4 }, { "cell_type": "markdown", "metadata": {}, "source": [ "First we'll specify the Hamiltonian." ] }, { "cell_type": "code", "collapsed": false, "input": [ "Omega = 2*np.pi*10. # [2\u03c0 MHz]\n", "Delta = 2*np.pi*0. # [2\u03c0 MHz]\n", "\n", "H = qu.Qobj([[ 0., Omega/2.],\n", " [Omega/2., -Delta]])" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 5 }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we'll specify the list of time points we want to solve the system over, and the intial condition \u2014 we'll start with all population in the ground state." ] }, { "cell_type": "code", "collapsed": false, "input": [ "t_max_in_pi = 5\n", "t_range = np.linspace(0, t_max_in_pi*np.pi/Omega, 100+1) # [\u00b5s]\n", "\n", "psi_0 = qu.basis(2,0) # Initial condition, all in ground state" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 6 }, { "cell_type": "markdown", "metadata": {}, "source": [ "To solve the system we use the `qu.mesolve` method. This has 5 required parameters:\n", "\n", "- the Hamiltonian to be solved (we'll pass `H`) \n", "- the initial condition (`psi_0`)\n", "- the list of time points we want to solve the system over (`t_range`)\n", "- a list of collapse operators. For a closed system we have no collapse operators, so this is blank (`[]`)\n", "- a list of operators whose expectation values we wish to find at each time point. If we leave this blank, we just get back the quanutm state $\\Psi$ which is all we want for now (so we pass `[]`).\n" ] }, { "cell_type": "code", "collapsed": false, "input": [ "ode_data = qu.mesolve(H, psi_0, t_range, [], []) # Solve with QuTiP solver\n", "\n", "print ode_data" ], "language": "python", "metadata": {}, "outputs": [ { "output_type": "stream", "stream": "stdout", "text": [ "Odedata object with mesolve data.\n", "---------------------------------\n", "states = True\n", "num_collapse = 0\n" ] } ], "prompt_number": 7 }, { "cell_type": "markdown", "metadata": {}, "source": [ "The return variable is an instance of the class `qu.Odedata`. As described in the [documentation][odedatadoc]: \n", "\n", "> This object is a qutip.Odedata class that stores all the crucial data needed for analyzing and plotting the results of a simulation.\n", "\n", "It's worth reading the rest of the docs for this class to see what you can do with it. As we specified an empty list of operators for which to find expectation values, the `states` property of `ode_data` contains a list of state vectors $[c_e; c_g]$ for each time point in `t_range`.\n", "\n", "We can then extract $c_g$ and $c_e$ and plot them over time.\n", "\n", "[odedatadoc]:https://qutip.googlecode.com/svn/doc/2.0.0/html/guide/dynamics/dynamics-data.html" ] }, { "cell_type": "code", "collapsed": false, "input": [ "c_g = np.zeros(len(t_range), dtype=np.complex)\n", "c_e = np.zeros(len(t_range), dtype=np.complex)\n", "\n", "for i, state in enumerate(ode_data.states):\n", " c_g[i] = state[0,0]\n", " c_e[i] = state[1,0] \n", "\n", "plt.plot(t_range, abs(c_g)**2, 'b', label=r'$c_g$')\n", "plt.plot(t_range, abs(c_e)**2, 'r', label=r'$c_e$')\n", "plt.xlabel(r't ($\\mu$s)')\n", "plt.legend()" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "pyout", "prompt_number": 8, "text": [ "<matplotlib.legend.Legend at 0x5311b10>" ] }, { "metadata": {}, "output_type": "display_data", "png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEiCAYAAACsmUZ+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXl8VNX5/5/JxhqWBKiAgiCrhCWQCQkkMSEJKlYNaBVk\nFWhrrVVb7Fet7RdsXaDV/kqr9msLLogiWjRiRZCEJQESMiGEJUBYBULYTMgKZJvz++PDJdtMZrvL\nufee9+s1xWZm7jnzmTP3Oed5zvMcC2OMBAKBQCDgDT+tOyAQCAQCgSOEgRIIBAIBlwgDJRAIBAIu\nEQZKIBAIBFwiDJRAIBAIuEQYKIFAIBBwiTBQAoFAIOASYaAEAoFAwCVuG6i1a9c++uKLL77e8u91\ndXWBs2bNWh0VFZU9ceLEnYWFhUPl7aJAIBAIzIhLA8UYsyQnJ2+eN2/eBxaLpVXZiVWrVs3p2bPn\n5ezs7KilS5e+sGjRojeV6apAIBAIzIRLA2WxWNjGjRvveeedd55kjFlaPp+enp44bdq0L4iIYmJi\nduTn549RoqMCgUAgMBcB7rzI39+/wc/Pz+7ouZKSktDQ0NASIhgzR6ssibaeEwgEAoGxcLSo8QSf\nN0mEhISUlpWVdZM648oIzZ7NqEMHRkSMXnqJEWPiwWpqiCUmEiMi1qsXsd/+llhBAS1+/HFi8+YR\n69gRzy1apH1fOXjU1jKaMgVjqEcPRr/+NaMDBxgtWLCYFixg1LkznnvySUZ2u/b91fxRX09s6lSM\noZAQYk8/TSw/nxYvXEjsZz8jFhyM5+bPJ2a3a99fjR8NDYymT8cY6taN0S9/yWjPHkY/+9lievJJ\nRl274rmZM/Farfur+cNuJzZ3LsZQly7EnniCmM3mq2kB7nbigw8+mPvCCy+83vLv//73vxf+5je/\neZMxRt9+++09M2fOXO3sGmiOsbIyxubNY4yIsXffZebGbmds9myI8Y9/MFZbe/OpxYsX4z8qKhj7\n+c/xmuXLteknJ9jtjC1cCCneeIOxmprG5yS9KisZe/ppvGbZMm36yQ12O2NPPQUxXnuNsevXbz51\nc3xVVzP23HN4zcsva9NPjvjtbyHF4sWMXbvW+HdJr6tXGXvpJbzmxRc16SJf/OEPEOP55zGWbnDj\nfu+T8fPIQL344ouvMcaoqKio7/Tp09cwxqi2tjZw+vTpayIiImzx8fFbi4qK+jpt7IaBYoyxujrG\n7rmHMX9/xjZskF0y/SB9uX/6U6unbt5AGGOsvp6xlBTGLBbGvvxSvf5xxquvQq7f/a71c031amhg\n7NFH8dpPP1Wvf9zx5psQ4Te/afVUs/FltzM2Zw5e+8EH6vWPM95+GxI8+SQkaUpTvex2xn72MzHJ\nZitXQoT581sJpqqBkuPR1EAxhoVBeDhjnToxlpfnq1I6RPpyFyxo/WtgjG3durX5H6qrGYuMZKxD\nB8ays9XpI0esXg25HnvMoVyt9Lp2jbGYGMaCghjLyFCnj1zx+eeY0Dz0ECx2C1qNr5oaxiZNYiwg\ngLG0NHX6yBHr1zPm58fY/fdjAt2SlnrV1TE2ZYqJJ9mbNuHDJyc38/xI6N5AMcbYuXOM9evHWO/e\njF244K1SOmTbNny5d9/t8Mt1ysWLjA0cyFjPnowVFSnXP87IymIsMJCx+PhmXiqXlJQwNnQoY927\nM3bqlGLd4489exhr146x6Gj4pNylrIyxsDDGunRh7OhR5frHGQcOMNaxI2MREYxVVbn/vspKxsaO\nxSS7oEC5/nHHkSOMBQczNmoUY+XlDl9iCAPFGGP792PS9tRTnqqkUxoaGBs9mrEBA7CM9JTDh7Es\nWLBA/r5xiN2OhWPfvoyVlnr+/hMnsOicOVP+vnGJ3c5YXBxjvXoxdvmy5+8/fZqxzp2x8jIJkycz\nFhLC2Pnznr+3uBgToPvuk79f3PLgg5jEnD3r9CWGMVCMMfbEE5ghnzjhrkI65uOPIf0nn3h/jWef\nhT/i0CH5+sUp69ZBrhUrvL/GCy/A25WfL1+/uGXDBgj21lveX2PJElxj9275+sUp6ens5qYbb1m2\nDNeQXMndu3fHVj8TPLp37+5QEzkMlAXXUQeLxcKctVdcTDRoENG0aUSrV6vWJfWprSUaPpwoOJgo\nL4/Iz8ud/pcvE91xB1FyMtG6dfL2kSPq64lGjsR/HzhAFOBW5l5rrlwhGjiQaMIEom++ka9/3GG3\nE40dS1RZSXT4MFFQkHfXqazE+AoLI0pPJ7L4lM7CLYwRRUXh/nPsGFH79t5d5+pVosGDiW6/nWjH\nDiI/PwupeW/VEovF8We98XffBo6vFs6TB7WxgmLMJLNcaZuQHFHVl182/CxX2keybp3v15Jmudu3\n+34tbpFW5x9/7Pu1li/Htb77zvdrcYq0Ol+50vdrvfsurvXVVzdXD6bA2WclI62giIjKyjDLjY42\n6Cy3qgrLxGHDiLZu9X1WavBZ7rVrREOGEPXpQ5Sd7fvHk2a5/fsT7dxpOLnkW51L1NRgrIaEENls\nvl+PM6TVucVCtH+/96vzptcbMYIoMJCooECsoORYQXE14rp1I3rhBaING4gyMrTujQIsX0508SLR\n66/Lc3cMDib6wx9g7DZv9v16nPHOO0RFRURLl8ojV8eOREuWEGVlEX39te/X445//5vo5EmMLzmM\nSbt2RH/8I4zdf/7j+/U448MPiY4cIXr1Vd+NExGu8corRAUFvl9LALhaQREZeJZbUoLlYUICUWqq\nfNeVZrnduxPt2WMYwcrLIZfVSrRxo3zXbTrL3bePyN9fvmtrSnU1VtNyrc4lGhqIxowhun6d6NAh\nCGcArl/HfaZvX0xY5JLLbieKjCTas0esoAy3giLCLPd3v8OgkaucExe8/z5RRQVmpHLSrh1WUXv3\nEmVmynttDfnoI6LSUqI//Une6wYEEC1ejFluerq819aUTz/F6vxPf5J3kuLvT/Tyy0THj8O1YRDW\nrcPq/I9/lFcuPz/5f+JmhrsVFBHu4717E82cSfSvf6nQMaVhDLGB0FAsC+Xm6lUI9uCDRKtWyX99\nlWEMk/aAACwK5aamBjPnSZOIPvtM/utrQnQ0fjgHD8q/iq6vJ+rXj2jcOMP4RhMSiM6eJTp6VP7Q\nmt1O5O8vVlByrKBk8LzKT5cuRI88QrRmDdFf/0rUubPWPfKRnTuJCguJ3ntPmet37Ahr/v77RH//\nO4J5OmbPHgSt33lHmeu3a0c0Zw7RW29ht37Pnsq0oxoHD2IXyV//qoyLNyCAaN48omXLiM6dg3XX\nMceOEW3bRvTaa8rs+3Dnms8+S5SfL097Y8YQ/e1v3r+/rq6Oli5dSj179qTq6mqaNGkShYeHy9M5\nH+HOxSexcCE2vRlihrtyJTY0/OQnyrWxcCEc6598olwbKrFyJVGHDkQzZijXxoIFRHV1cCXqnpUr\nERuaPVu5NubPx9Lggw+Ua0Ml3nsPRmTuXK17wgezZs2isLAweuKJJ+j48eOUk5OjdZduwqWLjwhu\nnjvvROx/1y6FO6Yk5eXYJz1rFtG77yrb1tix+DcvT9l2FKS6GnKlpGCXlZJMmIDUhoICHe8tqamB\nYElJRGvXKtvWpElEp09jCaLTLed1dfBWWq1E69cr144ztxdvZGVl0cKFC6ngxtbD4uJi6tWrFwV4\nsK3RVJskJCwWLAqysnS+bfPTTxEjWrhQ+bYWLsRmCR0bqP/8B6EUteQ6fBhjTLekpmI3iVqCnTwJ\n/5hO2bCB6MIFdeTSAzt27KBJkybd/P99+vTxyDgpDbcGiggei8BAeDB0y4oVRKNGEUVEKN/WY4+h\nVsuKFcq3pRArViA5NyZG+bYeeQTxTR3Lhc7370+UmKh8W9OmwaWhY8FWrMB+oilTtO4JH9xyyy3U\noUOHm/+/vLycdnHksuLaQPXq1bgxraZG6954QX4+UW4upmtq+JC6dSN6+GGijz/Gqk1nHDmCOmZq\nydW5M9H06fCMVVQo357snDpFlJaGgJoaLrf27eGqXrcOeX0649w5rKDmzZMnMdcIzJw5kywWC33w\nwQf00Ucf0fr16ykqKoq++OIL+uabb+j555+nTC3TV3ytleTJg7yoT7VxI+pbrV3r8Vu156mncCZP\nSYl6bW7bBsE+/FC9NmXiuedw7Iqa54JlZ0MuXZ6K+vvfo6L9mTPqtblvHwRbvly9NmXilVfQ9WPH\nlG/Lm3sdL5w7d44tXLiQMcbYM888w6pcHJDl7LOS0WrxOaKhARUFhg0j2rRJoY4pQU0N0S23wJfw\n8cfqtctYYwG77dvVa9dH6uvR5dhYdYuzMwYPbMeORLt3q9euz9jtcO2NGqV+4crISIzvffvUbdcH\nGEMZzH79UGhDafSyScIRX3zxBZWVldH8+fNp7ty59KGL3Uqm3CQh4e+P0Ep6OmLBuiE9HVvEZs1S\nt12LBW1mZqKygE7IyEBOklZy5eQQnTmjbts+kZ2NUghqC0aENvfvx24+nbB3L/Z3KLkT3yhERUXR\nhQsX6MMPP6TQ0FBN+8K9gSIimjoVKyldVThPTUXuU5MdMqoxdSqmjEruo5WZ1FTkPt19t/ptT52K\nf7/6Sv22vSY1FTuItIj2p6Q09kEnpKYiTHf//Vr3hH86dOhAv/vd7ygoKIji4+M17YsuDFREBNw/\nuvk9NDTgbjdlCsoWqM3IkUQDBuhGMMbQ1cmT4WpTmyFDUIlKJ3JBsC+/xOSna1f12+/XDzl3uhEM\nXY2JMUDVEBV45plnKDMzk06ePEn33Xefpn3RhYHy88OkbeNGnBHEPdnZRJcuNU7N1cZiQdtpabrY\nnpaXh7poWslFhLa3b9fJ5rRDh1C8VWvBsrKQVMQ5J07gNGZp4Sdom1WrVlFsbCy99NJL5K9xuX9d\nGCgiDK6rV3Vy7NGXX8L9cu+92vUhJQUH2Ml5VoVCfPklJiE//rF2fUhJ0ZEbWVq5PPCAdn1ISdGN\nG1mSSxgo/cH9Lj6J2lrkRU2dipqo3MIYDpoZPJjo22+160dDAzISExNRdZdjwsLgelFjd5Uz7PbG\nEjhffqldP9zCakUij5YlMHgZ524QG4vDp+UqzuoOet7F5ymm3sUnERSEGfbXX2NLMrccPAifgpbu\nFyJsf3zgASwJOM5yPnYMpay0lktyI2/axHmO89mzSP7WWjDJjZyezrUb+eJFHCagtVwC79CNgSLC\nDaSkRJkjlWQjNRU/Xi3dLxIpKZg6clw7TXK/PPigtv0gglzXrnHuRpa2GvLgr0pJQfVVjldQX3+N\nxR4Pcgk8R1cG6p57sCmOaxdMaipRVBSSdLUmKYmoUyeuBUtNJQoPR86p1tx1F6pFcSwXBBs+HFsP\ntSYqiuhHP+JasNRUottvRz6zQH/oykB17kyUnIxBx6V79/RpbEnjxZ/Qvj02anz1FYIsnHHhAsIo\nvMgVGMi5G7m0FKthXgST3MgbNnDpRq6sxGp46lQdH6dicnRloIiwVD99mtMqKzy5XyRSUmAJODqE\nTGL9ev7cLykpsAM7dmjdEwd88w02v/AmWGWltjtcnLBxIzZX8SSXwDN0Z6AeeAABbS69CqmpOGVx\n8GCte9LIffdhxxeHgqWmos5iWJjWPWnknnuw8ORQLgjWty/RuHFa96SRxES4NjgULDWVqEcPookT\nte6JwFt0Z6B69iQaP57D9J7KStS/462WSrdu2GfLmWDXr2PSff/9fLlfOnUiSkjgTi5sRti8GT5I\nnk6zbdcOJUA2buTK7263Y0fmlCnwRAr0CUcj3X2Sk7HT9soVrXvShO3bEbhITta6J61JTkZxT46K\nx+7YASPFq1xHj3JWPNZmwySIV8HOnOGqeGx+Pnb88iiXwH10a6Dsds7c3mlp8A3x6E+QfqXp6dr2\nowlpafA83nWX1j1pjSRXWpq2/WjG5s1YampRfNgVHAompQokJWnbD4Fv6KaSRFPq6ohCQlD1/5//\nlKFjcjBiBOID332ndU9a09CAMhwPPMBNGY5x4+BOy8jQuietYQzFiePjOSrCERODiD+Hm12IMQQT\nx4zhJhaVlASHwYED2rTvspLEs8/KV9pizBiiv/3N67fX1dXR0qVLqWfPnlRdXU2TJk2i8PBwt98v\nKkm0IDAQNw9uJmzFxSjgyas/wd8fM+/Nm7mIE/zwA87n4VUuiwU3uLQ0TnbnV1SgADHPgiUnE23Z\nwsX+/GvX4ELmVS7emDVrFoWFhdETTzxBx48fpxyOJkEBWnfAW5KTif77X6Lvv0cinqZIlpJnf0Jy\nMtF//kNUWIjjiTVkyxbYSd7lWr0aobsxYzTuzPbtWAXzLti//43gcFSUpl3ZsQNpWVwbKB9WPHKS\nlZVFBw8epKk3cuv+8Ic/UK9evTTuVSO6XEERNf5WuShLs3kztheOHq11T5wj/Vo5EGzzZhxjZLVq\n3RPncDe+OnYkmjBB6544Z9IkrKQ4EGzzZnhZ4uK07gn/7NixgyY1iWv26dOHAgL4Wbfo1kANH444\ngeZuPsbQicREvrb/tmTAAMQJNL6BMIYuJCRgkwSv9OmDlDYO7rcYX3Fx2hx+6S6hoTjEkAPB0tJg\nyzt10ron/HPLLbdQhw4dbv7/8vJy2rVrl4Y9ag7Hd9S2kdze6ekaxwkKClCpgWf3i0RyMkrl1NVp\n1oUTJ1AJRC9yZWZiO7xmFBURHT6sD8GSklC7qqpKsy5cvsx3fJM3Zs6cSRaLhT744AP66KOPaP36\n9RQVFUXvvfcepaen0+eff65p/3RroIgwCEtKMCA1Q5ox6uEXkZyMXBoNg6B6k+v6dY2r50suAr0I\nVl+PmJlGSJkUepCLB/z8/GjZsmU0b948mj17Ns2ePZtWrFhBlZWV1LFjR4928ynSP1cvqKurC5w1\na9bqqKio7IkTJ+4sLCwc2vR5xpjlF7/4xT/vuuuu7ePHj9+9bdu2eMV624LERPyrqZsvLQ2Vpfv1\n07ATbpKQoHmcIC0NUvFUDcoZcXFwQ2rqtUpLQ4rAyJEadsJNJk5ELqDG46tbN76qQekNadNEdHQ0\ndenSRdO+uDRQq1atmtOzZ8/L2dnZUUuXLn1h0aJFbzZ9Pi0tLenKlSvdt2/fftcnn3zy2DPPPLNc\nue4255Zb8LvV7PdQW4vZoh7cL0RIHouI0MyiNzRgB19SEl/ljZwRHEwUHa3hBEiKb+pFsPbtUVZL\nI8Gk+OakSaK8kS/Mnj2bNm/eTFu3bqWioiJN++LSQKWnpydOmzbtCyKimJiYHfn5+c023QYEBNRX\nVlYGM8YspaWlIV26dFH1eM3kZGwrvXZNzVZvkJVFVF2tL39CcjJyajQ4BTU3l6isTH9y5eXBlaw6\nBw4g21RvghUUIDdQZY4dQ8UlPcnFI1arlRYsWEAJCQk0duxYTfvich9VSUlJaGhoaAkRKkFYLJZm\nmZ4TJkzYdf78+d7Dhg07Ulxc3Gf58uXPtHW9JUuW3Pzv+Ph4io+P96rjEklJRH/9K4LZkyf7dCnP\nSUvDzj0fP4OqJCURvfYaNkuofOqvNLHmsVqPM5KSiP73fxHbeOQRlRvXQ35dS5qWPZozR9Wm9SiX\nkdi2bRttk/v0bsZYm4/p06ev2blz5wTGGNntdsttt912punzS5YsWfz73//+T4wxunTpUs/bb7/9\nVHl5eRdH10Jz8lJVxVhAAGMvvij7pV0zcSJjkZEaNOwD168z1r49Y7/+tepNJyYyNmqU6s36RF0d\nY8HBjP3iFxo0/uMfMzZkiAYN+0BDA2OhoYw9/rjqTT/8MGP9+jFmt6vedCuUuNfxirPPeuPvLm1M\nWw+XLr7ExMT0devWPUREtGnTprvj4uKaVU+rra0N6tmz52Uioq5du5a3b9/+estVlpJ06oSAaGam\nWi3e4No1VJjWWzZgu3Y4r0TlInh1dfCI6k2ugADk1KheM9Buh+9ab4L5+aFuoMqCMYZ7QFycPsJ1\nAvdwaaDmzp37YXFxcR+r1WpbtmzZ88uWLXv+3LlzfWfMmLGGiOi55557Y+vWrQkJCQlb4+Pjt730\n0kuvBgcHVyrf9UZiY7FzWtV8lZwcbJLQ2w2ECILt3Yst5yqRl0d09ap+5SooUDkOdfAgAnaxsSo2\nKhOxsUh4O39etSaPHUO4To9yCdrA1yWYJw9SaNm7fj1jRIxt367I5R3zpz+h0ZISFRuVie++Q983\nbVKtyb/8BU2eP69ak7KRkYG+f/WVio2+9RYaPXVKxUZlIicHfV+7VrUmV6xAk4cPq9Zkmyh1r+MR\nZ5+V1HDx6QHpCCZV3XyZmTirPCRExUZlIjoa+3BVFCwzk2jQIKQG6A2rlSgoSIPxdeutRP37q9io\nTISHq36WSmYmymEOHer6tQL9wHE1NPcJCYGtUO0GUl9PtGsX0ezZKjUoM5074yaikmBSOOXBB1Vp\nTnbatyeKjFRxfDGGm3t8vD4DKgEBmASpPAGKieFHru7du5OFl84oTPfu3RW7tiFWUETwPe/cqdJx\nNPn5qDemx4CKRGws8qFqahRv6tAhotJS/cu1Zw/S3hTn5EnEb/QcUImNRR5XWZniTZ07B8l4kqu0\ntNQjV9af/8yIiNH58yqEWy5fJkZE7LXXZLleaWmpYjoaykBVVRHt26dCY9LMkKdfhKfExsI45eYq\n3pRR5Kqvh01XHKMIxpgqhQyNIhcRPA2KIzWiA8EMZaCIVPIqZGbi+Iq+fVVoTCFiYvCvCoJlZhL1\n7o3TPvTKhAnYQa3a+AoJwXkfemX8eBzKpEIcKjMTXmvND5b0gbFjiTp0UHF8tWvH94FsNzCMgbr1\nVtgMxb9gxjAD0cHso0169sShWgoLJuWnxMbyEx/whq5dcR6lKjeQjAxMIHg+X8wVHTui7qNKE6AJ\nE/g+X8wVQUEqhu0yMzGB4Pl8sRvo+BfQmthYaM+UTBMuLMShM3oOqEhIgbuGBsWaOH0aRxoZRa6s\nLKS/KcaFC0THj+t/AkSEz5CbiwQ4hSgtRcqYUeTKzycqL1ewkaoqJCXqRDDDGajLl2FDFENyWejk\nC26T2Fj8Gg4cUKwJo8l17Rp+34phhICKRGwsSojs3q1YEzt3YkJqFLkYwwZhxcjKwoRUJ4IZykBJ\ns3RFl8mZmTifRw8HGrlCBcEyM3E+T1iYYk2ohipxzowMuMc0riItCxMnwq+r8PgKCkIagN6JioKb\nUvH7l58ffKI6wFAGavBg2A7Fv2CjFPzq1w8PhW8geg+nSPzoRzibUvHxFR2NDQZ6p3t3HNim8Piy\nWrHBQO+oUlc0MxM5kMHBCjYiHwa4bTRisTTGoRTh7FkEVXSyPHYLBQN3ly7B3Wo0uXbsQPKx7JSV\nEe3fbzzBsrLg6pOZq1cR4jKaXIrVFa2tRZ6EjgQzlIEigvbffw9bIjtGig9IxMY2BuZlRkfpFm4T\nG0t05QqKx8rOrl3GCahIxMYiu3nvXtkvvXs3ctOMJldtLYyU7OzZA8unI8EMZ6Ck9B5FAo07dyLh\nYtQoBS6uEQoKtnMndrKOGyf7pTVD8fHl748twEZB4fFFpJtwiltIdUUVG19NG9EBhjNQo0bBH61I\nxn92Nm4e/v4KXFwjhg8n6tJFEcGys5EKExQk+6U1Y+BApJApNr5Gj0Ywwij07Ut0222Kja8778Qm\nHKMQGoo4p2Lja+BABFN1guEMVGAgbopZWTJfuLoadZSiomS+sMb4+cHoyixYbS08CkaTy2LBZ5J9\nfDU0wK9jNMGIFBGMMdxvjSyXrGFhxnBRnQlmOANFhE1Qe/fKXAd1zx7cRKKjZbwoJ0RHIxeqqkq2\nS+bnQ3+jylVYiCRR2SgogP5GFezMGaLiYtkuefw4DpA0qlyXLhGdOiXjRYuKoL/OBDOkgYqKwgxe\n1oRKaQZopPiARFQUtqXZbLJdUpJLZxM2t5A+k6z5p2YQTEa/lZDLQ3QqmGENFJHMX3B2NhKtevSQ\n8aKcIBldGQXLzkboQc/1dJ1htcIzKvv46tGD6I47ZLwoJ4wdi0CkzOOrSxd919N1RlgYwpCyekWz\ns3Gw2ejRMl5UeQxpoHr3xkGksn3BOvXfuk1ICI4ilfEXYWS5OndG/qnsN5CoKGMkgLekXTskh8o8\nviIjjZEA3pKAAEyCZF9BRUToLgHcgF8viI6W8Qs+fZro4kXd+W89QhJMhsjs+fOQzOhy7d4tU8Ju\naSnRkSPGFyw3V5aE3epq5DMbXa78fNR+9JmaGsQ7dCiYYQ1UVBSSdc+dk+FiOvXfekR0NCrtnjzp\n86WkiYHR5aqoIDp8WIaLSVmZRhfs+nVZThS12TAxMLJcUVFIQt6zR4aL7d2LoLwOBTOsgZImC7J4\nFbKzUcBz5EgZLsYp0uCVQbDsbHgSwsN9vhS3yCgXLuLnp4sD5LxGxsCwdAkj7leSkHV86XjGaFgD\nNWYMXN+yuPmysnDz0POJaK4YMQLBFRkEy8pCXLx9exn6xSmDByN0J8v4ys5GZFwnBTy94rbbiPr0\nkW18DRmCpFaj0qsXcmplu3/16wf9dYZhDVRQEG6SPs9Arl3DElmHsw+P8PdH1NlHwerqEGowulyy\nJeza7Qhm6TA+4BEyCSYl6BpdLiJ8RlkSdnUsmGENFBG+kz17fDwBNS8PzmCdfsEeER2NGIEPJ6Du\n3w+bbha5Dh1CEXKvOXIEh0aaRbCTJ5GF6iWnTuHtZpArKgobjnwqfF1cjCRpnQpmaAMVFYUNLPn5\nPlxEx/5bj4mKQrWM3FyvL2E2uYh8rDxthg04EjLEocw0vmSJo+tcMEMbKOkL9smPm5VFNGCArgos\neo0MN5CsLOSh9esnU584JjISniufxld2Ng72GzJEtn5xy7hxiOP6OL46dTLGCc2uGDUKcVyf719S\nHpoOMbSBuvVWVDLweQai09mHx/ToQTRokE+CGTnftCVdumBviU/jS8poNoNgHTpg95KP4ysy0lgH\nCjhDlsLX2dmNlTx0iKENFJGPCbtSIpVO/bde4UPC7qVLRCdOmE8urxN2y8sRxDKbYDYb4roecu0a\n3PVmk8vrwte1tXDX61gwwxuoqCicsHvxohdvlqqBGjnhoiVRUThh98wZj98qxWLMJteVK0THjnnx\n5txcTARIlGf+AAAgAElEQVTMJlh1tVdHEkv7lcwmV22tlwcSHziA5GgdC2Z4AxUZiX+9KtSdk4Ol\nsc4KLPqED4Ll5CDf1Egn6LrC5/FFZOwE3Zb4OL6aXsIMyDK+hIHil7FjcdP0aqdVTk5jxq9ZGDUK\nRtkLwXJyEJMx0oGwrhg+HJ/X6/E1eDA2SZiFO+7A5/VyfN12G9EttyjQL07p2xebjrweXz176nrH\nkuENVKdOuGl6/AU3NCCJykyzWyIYpzFjPBaMMczyzDS7JUKwftw4H24gZhPMYsFvyksDJeTyAOkH\nqeMNOIY3UET4jmw2D+P+R47ghFOz/SKI8JmlE4Td5ORJFOU2q1z5+R4mhJ87hyRKswp28KBHCeEl\nJRhjZpXr6FEPE8IrK7EBR+eCmcZAlZZ6WKjbjA5vichIGOcjR9x+i9nlqqlBTNptpKCCWQVraPAo\n8m92uYg8zJ/fswczcp17gExhoKTvyKNlss2GRBczJFC2xAvBbDYkFY4YoVCfOMar8ZWTg6TVMWMU\n6RPXeCFYTg48VWbagCMREYF/Pb5/EQkDpQfCwnDz9GgnTE4ORoYRj+x0xZAhMM4e3kDGjtXdgZ2y\n0L8/YtEe30CkUgFm45ZbsNvBwwnQ8OHGLvjujO7dsZfG4/vXwIFIvtcxprj7Bgbi5un270E6WM2M\n/gSixrOJ3PxF1NcjR8WsclksjXFOt7DbzbmjpCkeCMaYOTdINCUy0osVus5XT0QmMVBE+IKlRD+X\n7NuHF5r9F7FvH4y1CwoKkOVvgN+D10RGIiZdWenGi48dQxUJswt24gR2P7jgzBlUKTG7XMXFbp4Q\nfvEiRDPA/culgaqrqwucNWvW6qioqOyJEyfuLCwsHNryNcuWLXs+PDx8b0RERO4333xznzJd9Q2r\nFTdRtxLYzZhA2RKrFUbajVLwZt4gIWG1Yqbv1hHdQrDG35Ybqyghl0dyGSb+ROSGgVq1atWcnj17\nXs7Ozo5aunTpC4sWLXqz6fM2m836+eef/8Rms1k3bNgw5fnnn1+mXHe9Rxrcbi2Tc3KQHde3r6J9\n4hoPUthzcuAnv+MOhfvEMR7fQDp1QlDFrIwbB9+oG4LZbEjPGzVKhX5xypgx2FPj9v3Lzw9xDZ3j\n0kClp6cnTps27QsiopiYmB35+fnNth1t2LBhyty5cz8MCAio79Wr16W1a9c+qlRnfUFKYHf7BqLz\nBDef6dsXR0S78Yuw2XCDNrNcPXogJu32DWTcOHOU5HZGly4w0G4IJhV00WlBblno0AEG2u37V1iY\nIUq6BLh6QUlJSWhoaGgJEZHFYmEWi6VZuuv58+d7X758ueeUKVM2XL16teNTTz311ogRI5w60pYs\nWXLzv+Pj4yk+Pt7rznuC2xnZZWVEhYVEs2er0i+ucUOw6mrkXN5/v0p94hir1Y2jEaTKn08/rUqf\nuMZqJfr2W/hGncxupPMzH39c5b5xiNVK9Omn2GPjdHOxtKNk6lRV+0ZEtG3bNtq2bZus13RpoEJC\nQkrLysq6ERExxiwtDVRwcHBlVVVV5w0bNkwpKyvrNnr06H2TJ0/+rkuXLhWOrtfUQKlNZCTR668j\ngb1jRycvkrLhzOzwloiMJPrqKxjtbt0cvmTvXtxEhFzQYO1axKidnm954ACMlBAMGnz4IQL6/fs7\nfMmRI5gECbmgwbvvYo/N0FY7AW6gYUmXlguOl19+2edrunTxJSYmpq9bt+4hIqJNmzbdHRcXl9H0\n+ejo6KyuXbuWExF16NDhWocOHa75+fl5czqO4riVwC6toaXsODPjRgq7geKxPuNW2E5swGnEDcGE\nXI24Nb4MVnLDpYGaO3fuh8XFxX2sVqtt2bJlzy9btuz5c+fO9Z0xY8YaIqJp06Z9UVtbGzRp0qQt\niYmJ6S+//PLizp07Vynfdc9xK4HdjBWmneFGCrsZK0w7Izzcjcr5UoVpJysGU+FG5fycHPMWdGmJ\nW5Xzc3IMVdLFwrw4OdXrxiwWpmZ7jujXj2jiRKI1a5y8oG9fooQEotWrVe0Xtwwdil9GaqrDp++4\nAzfm//xH5X5xyujR2AC6caOTF4wYQTRgANF//6tqv7hl/HjsAHASuxg3DnPFtDR1u8Urd92Fuo9O\nTwmPiUEcaudOVfvlCIvFQowxn7ZOmSZRV6LNuH9xMR7Cn9BIG4JJFaaFXI1IBTgczsMqK4kOHxaC\nNcVqdVo5//p1ov37hVxNsVrbqJwvlXQxkGCmNFBSHLEVUqzFQF+wz1itROfPw3C3QEpKFXI1YrVi\nbJ065eDJvXsNUWFaVqxWVM4vLGz11P79uOcKuRqxWrGCOnjQwZOHDxuupIspDRSRk4x/mw25KWas\nMO0MSTAHGyWkeKwB8gFlow25xAYcR7gxvoRcjZhtfJnOQEk3U4dfcG4uYgRO96CbkDFjEPl3IFhu\nLoLXTnagm5KwMMT9nY6vfv2IevVSvV/cMnQoIv9OxlevXtiEIwADBiAm53R8demCTV4GwXQGqnt3\nokGDHGzVlM4sN9DsQxY6doTRdrC3NTdXyNWSoCBslHC4FVgI1hp/f8wa2xhfZq5Q0hKLBZo4HF82\nG3aVGOiIION8Eg+wWh3MQE6fRtTfQP5b2ZAEaxL5v3CBqKhI3G8dIcX97U2zAa9cITp+XAjmCCny\nX1d380/V1agOL+RqjdWKGNS1a03+WFOD0wcMdv8ypYGKiCA6exYZ/zeRLJb4RbQmIoLohx+Q8X8D\nsZ/EORER2LB37FiTP4odJc6JiMCWvUOHbv5p714YeCFXayIisHlk//4mfzx4EAbeYPcv0xooohar\nKJsNJxuOHKlJn7hGEqyJXyE3F56E8HCN+sQxDuRqHGxmPLPcFU7GV9OnBI04HF8GLeliSgM1dix8\nuc0MVG4uggft2mnWL24ZNQrGu4lgNhvRnXcaomCy7AwfjtBdqwnQoEGiQokjBg0i6tq11fi69VZR\nocQRt96KWo+t7l+hoYarUGJKA9W5M24iN2cgdrsIYLdFu3YwUjd+EYwJudoiIAAry1Y3ECGYY6TI\nfxPBhFzOcSBX4wYvg+0oMaWBImr8ghkjBK8rKsQvoi0kwex2KirCEdxCLudERCCpv76eINaZM0Kw\ntoiIQFClpobKy4mOHhVytUVEBPJyq6oIxzMUFBjOvUdkYgNltWKTRFERiYi/O1itROXlRCdOGNXd\nLStWK3ZZHT5MYny5g9WKIP/+/WI/iRtYrXD87N1L2L3X0GBIi25aA9Vso4TNhgrAd96paZ+4polg\nublwY5n5CG5XNBtfublwvYgdJc5pMb6IxH6Stmh1/yIypEU3rYEaPRo32Zs3kPBw/EHgmBEjYMRt\nNrLZsNmxfXutO8Uvgwcjqd9mI/zPsGFEwcFad4tf+vXDMSQ3xteAAYj5Cxzzox+hwobNRrh/9e5N\n1KeP1t2SHdMaqA4dUJZmz27jVQBWhBuRf3Zjhivkahs/P6wAcm03dpQIwdqmSeRfyOUeNzdKGLgC\njmkNFBG+0wrbEQQZDfoFy0pEBLE9eVRR1iDkcoOICKLL+4pRdkMI5pqICGIFBXTx+6tCLjeIiCC6\ncKyCWGGhYS26qQ2U1Uo0pEJkBLqN1Up+V6tpGB0RcrmB1Uo0us54FaYVw2oli91O4bRXyOUGVitR\nOO0lC2OGHV+mNlAREURWslFd+86oqixomxs/guiAXAoL07gvOiAigiiCcsnuJ45wcYsbuyIiKFds\nkHCDceNw/yIiYaCMSFgYkdWSS2d6GKsCsGIMGUJX/TvT3aE2CgzUujP8c/vtRNGBuVTULQxBT0Hb\n9OlDJe37UGIXG3XponVn+CckhCi+Uy5d6tgfG0wMiKnvykFUS2Mon3LJmLMPuWkgf9pD4yiCHB1G\nI2iJhRhZySbGlwfYyCrGlwcYfXyZ2kDRwYMUxGppw2Vr86MRBA45epRod0ME9SttfjSCwAnff09d\n6kpp8xVr86MRBA4pLibacT2C+lQUorKLoG1KS+lH1Sdp+1UrXb6sdWeUwdwG6kZG4I6aCCos1Lgv\nOiA3lyiXIsi/rgbl/QVtc2N85bAIys/XuC86YM8ealwNSOUkBM65oVEuRTg+YdcAmNtA2WxU36U7\nnaSBhv2C5cRmIyrocGM7qxDMNTYbsaAgOkAjhVxuYLMR5VkcnYUjcMiNChJ5NM6wcpnbQOXmkl9k\nBHXsaDHsFywnublE3ccNJOrWTdxA3OHGES6htwQJudwgN5folrAe2F0iBHNNbi7R4MF0y7Bujo+A\nNwDmNVDXrhEdPEh+kVYaO5YM+wXLRX09ClNGWG9k/AvB2sZuJ9qzhywREUIuN2h2hIsQzD1uCGa1\nGteem9dA7duHu+6NG8jevTeORhA4pKAAp3JbrYT/OXAAfxA45tgxBPqtVrJaiY4cwTHwAsecOUN0\n+XKT8XXqFFFJidbd4peLF4nOniWyWikiguj8eWwyMRrmNVBNzpS2WnGvPXRI2y7xTLMjuCMiYM33\n79e0T1zTRLCICKwQ9u7Vtks802p8EYmNEm3R4v5FZMxFp3kNlM2GksC33nrz92DEL1gubDacyn3H\nHUSG/kXIhc2G5Nzhw8X4cgObjSgw8MYRLlIZCSGYc2w2FBcID6fRo4n8/Y3p5jOvgZIc3hYLDRqE\noxGM+AXLhSSXnx8R3XorUa9eQrC2yM0lGjuWKCCAevXCaRJCLufk5sI4tWtHmAkNGSIEa4vcXKLh\nw4k6d6aOHXEajhHtuTkNVFUVjjq9sRLw82tSul7QipoaePNulvtqcjSCwAE3d5Q0ZvgLuZzTbIOE\nhBDMOQ4EkzZKMKZhvxTAnAYqLw/fZIsbyL59uBkLmrN/PwpHNLuBWK0I2lVXa9Yvbjl8uNURLlYr\n0fHjRFeuaNgvTjl+nKi83MH4KirCUSWC5hQVYZNEi/tXSQnR6dMa9ksBzGmgmkVkgdWKm/CBAxr1\niWMkuZodORMRga3UIvLfGgeCibi/c5yOr6ZPChpxIJhRw8LmNVC33YZNEjcQvwfn5OYS9eiBOMpN\nROTfObm5ON598OCbfxJxf+fk5hK1b090551N/hgeDt+7EKw1ubk44XrUqJt/CgsjCgoy3v3LnAbK\nwRHJ/fsThYaK34MjJLksliZ/vOUWbJYw2i9CDmw2WKQmR7h07040aJCQyxE2G47LanaES6dOsFhC\nsNbYbLBITY5wadcO9spo9y/zGagrV+D0bnFEssVChs7I9parV5Gk6/BEaRHIbk1tLYKZDgQTcrWm\noQEh4TbHl9Ei/74gbZBwIJjVCheykU5mMJ+BysvDvw5OoIyIwM346lWV+8Qx+fkY8A4P7IyIwBkc\nZWWq94tbDh6EkXIyvs6cIbp0SYN+cUphIfbZOB1fly6hYoIAnDqFSbaT8VVRgSImRsF8BkpaAzs4\nU9pqxYxOHI3QiCSXwxuINIuTjL6gTcEkucQqqhG3xpcQrBGTjS9zGqiBA3FecgtE3L81NhtR795E\nffo4eFJE/ltjs2FsDRjQ6qnwcLiShVyN2GxEnTsTDR3q4MlRo7AZQAjWiM2GgFNYWKunhg9HWMpI\ncpnTQEVGOnyqTx88jPQF+0obcmFXycCBQrCm2GyYyjbbUQKCg3ETEXI1Iu0n8fd38GT79saM/PuC\ntKMkKKjVUwEBZLiTGcxloC5cgD/b6R0XT+XkqNgnjikrQ4ipDbmEYE2prkYMyo3xJeL+CNXl57sx\nvmw2Y0X+vaWhAbsgXIyvvDzkdBoBlwaqrq4ucNasWaujoqKyJ06cuLOwsNDRYpzsdrtfdHR01qZN\nm+6Wv5syIU0tXHzBx46JjH+iRl+2yxvI2bMi458ISct2u8vxdfmy8TL+vWH/fhgpl+OrogIzJbNz\n+DAmQS7G1/Xr2OxlBFwaqFWrVs3p2bPn5ezs7KilS5e+sGjRojcdve4f//jHrwoLC4daLBZ+54Y5\nOfAlhIc7fYkRA43eIi2MHAawJYyawu4NkmAO90wD6d4i5HJLLiFYU9wQTHrKKE4NlwYqPT09cdq0\naV8QEcXExOzIz88f0/I1Z86c6bdx48Z7Hnzwwa8YY62d77yQk4PgYseOTl8i3YyN8gX7Qk4Oikp3\n69bGi8LDYfSFYNCgX79mFUpaMnIkwgdCLmggVXp3yrBhSNoVgkGDrl2bVShpibT/yyhyBbh6QUlJ\nSWhoaGgJEZHFYmGOVkhPP/303//617/+5s9//vP/uFpBLVmy5OZ/x8fHU3x8vMed9grGMAt7+OE2\nX9atG3YUiQkbNJg0ycWLOnWC0TfKL8IX2txRAoKCYNOFXI1yOdhP0oi/P2aNQrDGDTh+ztcVFktj\n2E5ttm3bRtu2bZP1mi4NVEhISGlZWVk3IiLGmKWlAVq9evWskSNHHhg+fPhh6TVtXa+pgVKVEycQ\nWGrTnwCsVqK0NNi0Nn88BubcORwh7eJ+C6xWonXrzC3YDz8QnTxJ9MQTLl8aGUn03nuIeTvcvWYC\nKioQUpk+3Y0XR0YSLV+OgJWD3Wum4No1BO1++1uXL7VaiV59FeGqTp1U6NsNWi44Xn75ZZ+v6dLF\nl5iYmL5u3bqHiIg2bdp0d1xcXEbT53fs2BGzdevWhISEhK0bN26853/+53/+vGvXrgk+90xupBmY\nG3fcyEjE/M+dU7hPHONWfEAiMhLG/8QJRfvENdKU1c0JUHU1btBmZc8ezGfcGl9WK4zT/v2K94tb\n8vNxzpib9y+73Rj58y4N1Ny5cz8sLi7uY7VabcuWLXt+2bJlz587d67vjBkz1hAR/d///d8TkpG6\n5557Nv7lL3/57YQJE3Yp33UPkY7gHjHC5UtFXBafPSAAKRcukQQzsxvGZsPq0UGFkpYIuTyy50Iw\nIo8nQE3fomcsTMWEDIvFwtRsrxkTJ+IGsmOHy5dev46kyueeI3r9dRX6xiFJSciDcms3Y309UZcu\nRD//OdH/+3+K941Lfvxjou+/Rx6UC+x2BLJnzCD65z+V7xqPPPwwFgXHj7vxYsZQPX/KFKL331e8\nb1wyaxbR1q1uu3X69yeKjib69FOF+9UGFovFZcjHFeZI1K2rw3rXrYAKEthHjzbvhM1ub4zHuoWU\nwm5WwRjDZ3dTMD8/vNSschF5JFdj5N/sgrl5/yLSbqOE3JjDQB08iGWR278IfMG5ueZMYD96FEFs\nD34Pxkth94TTp5F964FgVitCKtevK9gvTnGjoEtrrFYE7SorFesXt1y5guoBHt6/Tp7E3h09Yw4D\n5UYFiZZYreZNYPdCLghmpBR2T/BCsMhIeEbNWDnfq/EVGYmV6p49ivSJa9wq6dIco8ShzGGgcnLg\n9B840O23mDkum5ODCtPDhnnwJrML1q4dsnDdxOxyuSjo0hqjlUjwBLdKujRn3DhjVM43j4FyUmHa\nGcOG4SZt1t+D0wrTzjBaCrsn5OQ4rTDtDKlyvlnlclHQpTVS5XyzCuaypEtzgoOJ7rxT/3IZ30BV\nV8Pt5JE/ATfnceP0PwPxFLcqTDvCYjFn5N+NCtPOMGPcXyro4oVc5hXMww0SEtLPUc+V841voPLy\nXFaYdkZkJG7WNTUK9ItT3Kow7YzISEwGqqtl7xe3uFFh2hlmrJwvFXTxenyZrXL+uXP4vF6Or8uX\nic6cUaBfKmF8A+VRSYTmREbiZr1vn8x94hgf5GpMYTdTINvH8UVkrlW6z+Or6UXMgAzja/duGfuj\nMsY3UNnZRLff3maFaWeMH49/9fwFe0p2NnIi26ww7Qwj/CI8JTubqHv3NitMOyMiAp5Rs8nVqZNb\nBV1aM3Yscu7MJphUYdhDRo5ETqee5TK+gcrKIoqK8uqtt91G1LcvLmEWJLm8qvnaqxfRHXeYT7Dx\n49usMO2Mrl0RyDabXFYr7IzHdOiAzShmEyw8HLtEPSQoCJMgPctlbANVVAQfbnS015eIisIkxgz8\n8ANKz/ggFwTLytJ3ZNZdKioQc/NBsOhojC8zyHXtGmK6Po+vnBxsTjE6dXXIgfLx/pWXp984urEN\nlDR18HIFRYSxceoU0cWLMvWJYyRD7INcEOzCBX1HZt1F2iLlg2BRUdg0YIaE8D17kJzs8/iqrnar\n5qHukUqN+Hj/qqnRb0K4sQ1UdjaWxm6V5HaMNDbMsIrKzm48H85rzCaYxdIYrPQCaXJsFrmIfDRQ\nZhpf0gTbxxUUkX7lMraByspCMpMPh5yNHUsUGKhvP667ZGWhSK5HCZQtGTUKsQKzCDZ8OIJJXjJs\nGN5uFrkGDkSo0msGDMAFzCBYdjZR794IhntJnz7Y8KRXuYxroGpq4Hz1yeHdGJfV6wzEXRoa4LHy\nUS5Y84gI4wvGGD6jj4L5+WEBptcbiLswhs/o8/iyWBrjnEZHEszHU6r1LJdxDdS+fTBSPvkTQHQ0\nclXq62XoF6cUFBBVVckiFwTLyzN2qe5jx4hKS2URLCoKIRUjF+o+e5bo/HkZx9fRo0QlJTJcjFMu\nXUI5cpnuX2fOEBUXy9AvlTGugZLBfysRFUV09SrRgQM+X4pbpAWPDHJBsLo6or17ZbgYp8goWHR0\n4xlcRkWW+JOEdBE9J/i4Qsbxpec4lHENVHZ2YyKTj5ghkJ2VRdSzp0cF352j51+Eu2Rl4RTh4cN9\nvpS0x8LockkHgfqM1QrfqJEFy85Gsti4cT5fKjwcYXg9ymVcA+VDgm5L+vdHIQq9+nHdwacE3Zb0\n7g3RjC6Ylwm6LeneHZsljC5XRARClD7TqRM24xhdsDFjEAT3kXbtsNlLj3IZ00CdP49TTmXxV+Gm\nLSVUGpHSUqLCQtnkAkYWrKoK/l4ZBZMSwo2YsFtTA2+v7ONr925jJuzW18PfK9MEmwhy5ebq78Br\nYxooWR3ejZc6dkz/Ryg7QnLlyygXLnb2LCp5GA2bDUEjmcfXDz+g2rfRyMtD0WXZx1dlJarJG42D\nB5GMLPME6Pp1/RW+Nq6BCgrCulYmpLFixLhsdjY8VV5VmHaGkQN30mfyIUG3JWaQS1YDJQTzCL3K\nZUwD5UOBRWdIJ8zq0Y/riqwsVD7u3FnGi44ZA/2NKtjQoThBWCZGjID+RpWrXz8kjcrGoEE4Zdeo\ngvXqhaRkmbj1VuivN7mMZ6BkKLDoiE6dsANJbzMQV9jtWBXKLBdWsOPGGU8wmRJ0W+Lvj9NKjCYX\nkSJyNSbsGlkwWXYsAb3G0Y1noPbvR9lkWf0JICrKeHHZw4dRlFsBuXDR3FwEIIzCyZM4plSh8bVv\nn7EOJD53DqFIxcbXoUNEZWUKXFwjSkqQhKzQ+Dp5Ul+Fr41noHbswL8TJsh+6YkTsYFr/37ZL60Z\nCsoFwWpqjHXCrsLjq6HBWHFOxccXEdGuXQpcXCN27sS/Co2vpk3oAeMZqMxMnKDrQ4FFZ8TGNjZh\nFDIzcYLuoEEKXDwmprERo5CZicQlr46EbZsJE+CKMZpcnTp5dSCsa8aPR2KV0QQLCmo8nVpGxo1D\nWpWe5DKWgWIM6kuWRGZuuw35p3r6gl0hySWju7uRXr2wmcBogsXEyJKg25Ju3ZB/ajS5oqNlStBt\nSceOuOsaTbDISJTdkJmgINh0PcllLAN17BiKLCpkoIhw6YwMYyRUnj6NIpIKyoWL79iB3Rh65+JF\nxAcUHl9ZWfpLqHTElSvIZ1Z8fNlsiDvrnepquMMVHl979+qnMLGxDFRGBv5V+Au+dAm2UO9IMynF\nbyBlZcY4AVUFweLiUJg4L0+xJlRj505M5BQfX7W1OCtG72Rno4qEwvcvu10/YTtjGajMTFQ8HTpU\nsSbi4hqb0juZmTgsb+RIBRsxmmAdO8qaAN4SI8U5MzPh2pMxn7k1MTHGCdxlZuKzKLKjBERHI6VB\nL3IZz0ApFlABQ4fCBkqLNT2TmYmdPf7+CjbSvz+yBI0iWFSUTyc0u0LasKKXG0hbZGaiQKxPJzS7\nont3orAw4wg2ZoxPJzS7onNnzK/0IpdxDFRREdGpUwr7E2D7YmL08wU74/Jl5EApLBcEi42FYHoO\n3JWXI0lJccGw6NR72O7qVaTAqSAXGtm1S98nitbWIviogmCxsUhlqKlRvCmfMY6BUiWg0tjEqVP6\nroMq5aeodgM5fx5Zgnpl1y5YDJXGV2kpclD1yu7d2OgheXgVJTYWCYr5+So0phB5edjoodL4qqnR\nxwGZxjJQnTvLdCJa2xghTpCZiVJ5EREqNGYUwQICFCqJ0ByjyGWxNCaHKopRBCNSxUDpKT3RWAZq\n4kTcRBRmzBjYQj2HVaRwioz1dJ1z550orKp3wcaNQ9apwgwciDMf9XADcUZmJjbfdOumQmN9+0I0\nPQuWkUE0eDBORlWYHj3wk9SDXMYwUKWl2Masir8KNnDCBH18wY6oqkIuhEpyIalVz4G769exjVkl\nwSwWuMb0Grarr1ctnNKIlG+nR8HsduzJV8UfCmJj0STvdUWNYaCk4lIq/iJiY2ETS0tVa1I2srIw\nMFW/gRw/TnThgoqNyoTNhiC2yuOrqAjJ1Hpj717knKp4v4Vgly/jaGi9UVCArGaVx1dFBRKpecYY\nBiojQ7H6Vc6QxpKeCi9KZGRgUSP7EQhtoec4geSaVCWgAiS59OgVVSFfvjV6FkzF+JOEXuRyaaDq\n6uoCZ82atToqKip74sSJOwsLC5tlwdbU1LR79NFH144fP353dHR01ubNm5OV664TFKxf5YzISNhE\n3r9gR2RmonhncLCKjY4di4QYvQo2YgQOyFOJsDDEb/RozzMzie64A3E01ZDiN3oULCMDpwnKeECh\nK/r100ddUZcGatWqVXN69ux5OTs7O2rp0qUvLFq06M2mz69Zs2ZGjx49fti9e/f4r7/++v5f/vKX\nbyvXXQeoUL/KER064Ih0vd1va2qwBVhV9wsRSgpER+tPsPp6bDFXWTApbKc3uex2hIJUH19Svp3e\nCmVKBa7j4hQtMOAIPaQnujRQ6enpidOmTfuCiCgmJmZHfn7+mKbP33777d8/8cQT/0dE1L59++tV\nVfKLskoAABtMSURBVFVyHhzumowM3ETi41VtlogoIQHJiHo6L23XLsT8NZALgu3fj2KGesFmQ2VN\njcbX0aOIRemF/ftx5p5m4+vMGaITJzRo3EsKC4mKizUbXxcvIgTGKy73ZJeUlISGhoaWEBFZLBZm\nsVia2dv4+PhtREQHDx4M+9nPfvav55577o22rrdkyZKm76V4X7+YtDTslVZ5BUVElJRE9MorRNu2\nEaWkqN68V6SlobSRJjeQpCSi3/+eaMsWounTNeiAF6SlYWabmKh600lJjV2YN0/15r1i82b8K/Vd\nVZoKpsgBZwogCZasfmSkqVxhYb5fb9u2bbRt2zbfL9QUxlibj+nTp6/ZuXPnBMYY2e12y2233Xam\n5Wtefvnl/x01atS+LVu2JLR1LTQnMyNHMpaYKP913aCmhrFOnRh78klNmvcKq5WxiRM1ary+nrFu\n3RibP1+jDnhBbCxj48Zp0rTdztiPfsTYY49p0rxXJCczNmKERo3b7Yz168fYtGkadcAL7r+fsYED\nNWt+yBDGpkxR5to37vcubUxbD5cuvsTExPR169Y9RES0adOmu+Pi4pp5xdesWTMjNzc3wmazWRMS\nErbKaz5dcOEC9klqMPsgwiaJu+7CDEQPlJbCJanJ7JYIS7dJkzBr5NnxLVFZiT35GglmsaDptDR9\nyHX9OmIamo0viwX3gi1b+E/wIUItqG3bNBQMcm3fjiwKHnFpoObOnfthcXFxH6vValu2bNnzy5Yt\ne/7cuXN9Z8yYsYaIaOPGjfecOnVqwN13370pISFh66RJk7Yo3+0bpKfjXw2/4KQkxAnOnNGsC26z\ndStudBrZc5CURHT2rD4O1JLimxoKlpSEkB3v+SpESLm4fp2D8VVWho1TvJOTg0mQxuOruhpHUXGJ\nr0swTx4kt4tv7lzGQkMZa2iQ97oecOAAY0SMrVihWRfc5uc/Zyw4mLHaWg07cewYBHvrLQ074SbP\nPMNY+/aMXbumWReKiiDXG29o1gW3ef55xgICGKus1LATly5BsFde0bATbrJ4MWMWC2MlJZp1oayM\nMX9/xn7/e/mvTWq4+LiFMbiKEhOxJ1cjRozAGT56cPNt3ozNEYGBGnbijjuIbr9dP4LFxqqaX9eS\nvn2Jhg/Xh1xpacgk6KzuPt7m9OyJJD9p8wHPpKWhvmNIiGZd6NoVOZ28yqVfA3X4MLZnaujeI2oe\nJ+D5/J6TJ/HQ1P1C1CjYli18n99TXIzzLjQXDHJt3873+T0lJTgxggO5INiuXSg6ySsVFfCrcSBY\nUhKyKa5c0bonrdGvgZKmlBx8wcnJRD/8gPPseIUjudCJigq+D6ThSLDkZBwVtGuX1j1xzpYtHMQ3\nJZKTsQGB5zIJ27djIwcHgiUnY3It9w5xOdCvgdq8GbkOt9+udU+a5RPwyubNcBcNHer6tYozaRJW\nUrwL1rMn0ahRWveE4uOxAZJ3ubp2Vel8MVfExCA3kle/FRH61qEDjkXQmKgouGV5lEufBoqD7ZlN\n6dMH56vw+AUTYaK2ZQtmSipXU3FMjx58xwkYgzXQOL4pERyMmwjPcm3ejMoEKhzH5poOHWCkeBWM\nCH2Li1PpQLa2CQxEugyPcmn/6/OG3bvhX+ZgeSyRnAyPwvXrWvekNXv3IgeKI7nQmawsPuMEBQXI\nseNIsORk5LDxGCc4eZLo+++5kgudOXiQz+NdioqIjhzhSrDkZJyG8/33WvekOfo0UJs3Y2abkKB1\nT26SlATjxOPxG9LMSINqPc5JSsImie3bte5JazSt1+OY5GSsVLaol2XoNhzK1Xjz59EvKvWJI8F4\nlUu/Bioigqh7d617cpO77oJ7g8dl8ubNCKWocJq0+8TEYPs2r4INGYIzCTjBaoWrj1e5+vXDiRfc\nMGYMjkfhVbBevYhGjtS6JzcZPhyhCt7k0p+BunQJLr577tG6J80IDkbKzNdfa92T5pSVwfXImVww\nTgkJEIynOj5VVVimcCZYYCBmuf/9L1/pDDU1RN99B7m4iG9K+PkR3X030YYNfKUz1NcTffst+sZB\nfFPCYsF3uGkTX2WP+FHIXb7+Gr/QqVO17kkrUlKQOnP0qNY9aeSbb/Cb4FAuCHbyJGIFvLBpE+66\nHAqWkkJ07hxfVXzS02HTOZQLgv3wA1/78zMyEEjkULCUFKLycr62m+vPQKWm4ijI0aO17kkrHnwQ\n/371lbb9aEpqKipdREZq3RMHPPAApm6pqVr3pJHUVGT2x8Ro3ZNW3HcftpvzJldwMFfh4EbuuQe7\n5HgTrH17osmTte5JK5KSiDp14ksufRmoyko4SVNSOPMngP79cbL5l19q3RNw/Tq8CQ8+yJU3oZFb\nbkFtHF4Eq6uDD+3++znZL92ckBDkRPEiV0MDJmNTpnCxW7o1wcG46375JR9uZMZw9588GZaAMzp0\ngE1PTeXHjczjbcs5HLtfJFJSUMHk/Hmte4IdOdXVXMsFwfbuJTp9WuueYEdhWRnXgqWkoMpXYaHW\nPcE4v3SJa7kg2Pff46hfrcnLQyV/jgVLScG9i5ciL/oyUKmp2JkzcaLWPXFKSgomSjxslkhNJerS\nhVP3i4R0FDEPftHUVEwjOcpPaQlPbuTUVGzeuPderXvSBvffz48bOTUVrowf/1jrnjjlvvvgPOBB\nLiIiC1Nx6WuxWJjX7dXWYmvm1KlE778vb8dkhDFstx08GO41rWhoIOrdG7lPa9Zo1w+3CAtDWaGt\n6p532Qy7HXulrVZ+fGhOsFpxE8nK0q4PvIxzt4iNRXggP1/bfvAwzt0gORkLvSNHfLuOxWIhxphP\nsRj9rKC2b8cWE46Xx0SYrKWkYHdTRYV2/di1i+jyZe7lAikp2N1UUqJdH/bswRY5HQjGgxu5oIDo\nxAldyAXB9u0jOnVKuz4cOwbRJI8Bx6SkwIXsq4GSA/0YqNRUoo4duXa/SEydini7ljPL1FQcSc9Z\nOo9jpk7FCua//9WuD6mp2CLHsftFQrrHrV+vXR9SUzEZe+AB7frgNpJgWvqtJJ+sDgyU5EbmwZGg\nDxef3U50222omLlunfwdk5mGBmRlJyQQffqp+u0zhnMBhw9HHhT3MNa4BVKrm8idd8Inmp6uTfse\nwBiq0g8cSLRxozZ9GDcOu6V5LO3lkFGjUHlGq9JaEyfizJS8PG3a95DISExAdu/2/hrmcfHl5uIA\nOR3MPogwEX/wQSSxa3HI3IED8GboRK5Gv+h33xFdvap++4WF2BqnE8EkubZsgddbbc6cwX1WJ3KB\nlBSiHTvg91abCxcQMNSRYCkpRDk58HpriT4M1Jdf4q5/331a98RtUlIQl9WiuOeXX+rI/SKRkoIZ\n5qZN6rctrdok34YOSEmBG3nDBvXbluTS0f0WnbXbtdleu349lr26CNgBXjbX8u/ia2iAL2PYMG1u\nXl5SU4M81ClTiD7+WL12GUOd0z59+CwU7pT6enQ6NlZdNy5jcP907OibP0Nl7HZ4RUeNUt+NGxmJ\n8c3zCdKtYAwHnPbrp/4uurg4oosXseuAwwIDjmAMmw67dPF+t6g5XHxpafApLFigdU88ol07olmz\ncK8tLVWv3YwMnOuiM7mwb3ruXMw2L15Ur92cHNQC1Jlgfn5E8+YhBnX2rHrt7tuHJE6dyQXDMH8+\nCs0dP65eu4WFqNa8YIFujBNRo1zZ2dh8qBX8G6gVK5CcqyP3i8TChZhpqrmCWrECs56HH1avTdlY\nsAArqVWr1GtzxQqsnqZPV69NmZg/HyupDz5Qr82VKxsnX7pj3jxY9vfeU6/NlSsx+ZozR702ZWL2\nbCRir1ypXR/4dvFdukR0661ETz1F9Ne/KtcxBbFaG90hSk+gysqwEe3xx4neeUfZthQjNhbfuxru\nkKoqCPaTn6h705KR5GSk2Jw8qXy9xWvX4IW9916iTz5Rti3FuP9+bLo6e1b5eou1tdh9PHEi0Rdf\nKNuWQvzkJ/CInjvneb1F47v4PvoIkWDd+RMaWbgQu+pyc5Vv65NPUCB24ULl21KMhQtxXsmOHcq3\n9dlnMFI6FmzhQpQxVGN3/JdfYhKkY7nQ+QsX1Nld8t//YrKlY8EWLkT+vFabJfhdQTGG3JTu3fk6\nz8VDKiowSZ81i+jdd5Vta+xY/KuTVAvHVFdjmp6SQvThh8q2NWEC7rgFBbqKDzSlpgZyJSURrV2r\nbFuTJsEYHjvGaXV8d6irayxppXSm85QpKFJ7+jR2IesQaY/a0KHIAvEEY6+gdu2Cm0fHsw8ixIMe\neQSrm6oq5drJy0NRcJ3LhWMIHnuM6PPPYTyUoqAA25MWLtStcSKC22XOHKxufvhBuXaOH4erZ/58\nHRsnIgRV5s3D1sfiYuXaOXsWO1gef1y3xokIXZ8/H6ccaVEpit+htmIFUefOuLvrnAULYJw+/1y5\nNlasQGb/Y48p14ZqLFiAgIeSVW5XrsTNavZs5dpQiQULsDD46CPl2njvvcadg7pH2l2i5Ar9/ffh\nBZo/X7k2VOLxxzGH06JGN58uvvJy+C1mziT617+U75jCMIayQ6GhypSGuXoVbsQHH1R3A5xiMEYU\nHo7pmxLnm9fUEPXtC5/VZ5/Jf30NiI6GO/ngQfkXhPX1iPVHRPBxjIwsJCRglXP0qPxLQskvNmQI\nlh4G4N57EUv3xFtpXBff6tW46+reXwUsFnyUXbuUOTft009xczKIXI2C5eUpc3LaunWI/BpGMHyU\nQ4eU2Vuyfj32FRhILnyYEyeQZyk3mzYhd9NAgi1ciJ18atdz5m8FdfUqDpnp3x/LDR3HB5pSWko0\nYAAmbnLWQ62pQZGN7t2x2DCIXFhFDxyIYLacFVHr64lGjIB7b98+XccHmlJdjQLBw4YhViTXOGho\nIBozBuPs0CHld2arxvXruM/07YtYpFyC2e0otfHDD0jS9XRvNqfU1WGjRJcumDe6s+iUYwVFjDHV\nHmjOBcuWMUbE2Pbtrl+rM155BR9t1y75rvn3v+OamzbJd01ueOMNfLgtW+S75r/+hWt+9ZV81+SE\nt97CR9uwQb5rrlqFa65dK981uWHFCny4L76Q75qffYZrfvihfNfkhNWr8dE+/ti919+43/tkM/ha\nQV25glnzhAk6OSfCM6RZ7tChqLji66StshLXCwtDHoxhVk8S0iy3Tx/UXPH1A0qr89tvhy/MYILV\n1iLWGRzs/iy3LWpqMFZDQ+Fp1fXuPUfU1xONHIlxsH+/78vDujqszoOCDLU6l7DbERquqkLx/6Cg\ntl9vvBXUCy8wZrEwlp/vnonWIW+/Ld8s9+WXca3du32/FresXIkPuW6d79eSVucZGb5fi1M+/tiz\nWW5bLF+Oa333ne/X4pZ16/AhV670/VrvvmvY1bnEN9/gI779tuvXkqFWUMXFqDY8bRo2SRgUuWa5\nly9jsTl5si7OcPQeaZZLhG1E3s5yDb46l7DbkbBdWeneLNcZ0up85EjsIzDYYrMRxnAQanExMpDb\nt/fuOgZfnUswRhQfj/Da8ePIBHKGsVZQTzzBWGAgYydOuDbNOkea5X7yiffXePZZxvz8GDt0SL5+\ncYs0y12xwvtrmGB1LrFhA+R66y3vr7FkiQlW5xJbtuDDvvGG99cwwepcYtcufNRXXmn7dSTDCooP\nA7V/P2P+/ow99ZR7CumchgbGRo9mbMAAxsrLPX//4cOMBQUxtmCB/H3jErudschIxvr2Zay01PP3\nnzjBWIcOjM2cKX/fOMRuZywujrFevRi7dMnz958+zVjnzow99JD8feOWyZMZCwlh7Px5z99bXMxY\nt26M3Xef/P3ilAcfZKxLF8bOnnX+GmMYqHPnGLvtNsZ692bswgXPVNIx27bBJt99N2O1te6/7+JF\nxgYOZKxnT8aKipTrH3dkZWGFfdddjF2/7v77fviBsaFDGevenbFTp5TqHXfs2cNYu3aMRUczdvWq\n+++7coWxESNw8zl6VLn+cceBA4x17MjYuHGMVVW5/77KSsbCwxnr1ImxggLl+scZR44wFhzM2KhR\nzifZ+jdQFRX4cjt3ZiwvzzuldIwU/1+wALNeV1RXYyHRoQNj2dnK9487pH2ujz3mnmDXrjEWE4Pl\npglcLy35z3/g1XzoIazaXVFTw9ikSYwFBDCWlqZ8/7hj/Xr4zX/8Y8bq6ly/vq6OsXvvxUxTzr39\nOmHTJnz05GTHk2x9G6i6Osbuuce0X67EH/6Ab+FPf2r93NatW2/+d309ltUWC2Opqer1jztefRWC\n/e53rZ5qqhdraGDskUfw2k8/Va9/nPHmm5DgN79p/VxTvex2xmbPZkZN4XEfaZvtk0+2mgQ1G192\nO2M//Sle++676vaRI6RJ9uOPt54zqmKgamtrA2fOnLl6/Pjx2RMmTNh55MiRoZ4836wxyUBducLY\nvHlo/l//klszXdH0xvD3vzefiSxevJgxhiX0z3/e+BpTY7cztnAhuxnUrqm5+ZSkF6usZOzpp/Ga\nZcu06Scn2O2M/epXkOLVV5t7RyW9qqoYe+45vObll7XpJ1f89rcQY/HiZv7Rm+Pr6lXGXnoJr3nx\nRU26yBPSJPv55+HlkVDFQK1YsWLBs88++/8YY5SRkRF73333/deT55s1RsTYrFmMtW+Ppl96SRnF\ndEZNDWNJSZCkVy/cLAoKGHv88cVs7ly4xokYW7RI655yQm0tY1OmQJQePRj79a8Z27+fLZ4/n7H5\n8xEPIGLsl790zxVocOrrGZs6FZKEhMBg5ecztnDhYvbTnyKW4Imr2fA0NDA2fTpE6dYNq6k9e9ji\nn/6UsV/8grGuXfHczJnu+U4Njt3euN7o0gWT6ZwclQzUjBkzPsnIyIhljJHdbrf07du3yJPnmzUm\nfYInnsAnENykro6xr7/GjSQgAN8M0WIWHAxPQlaWuHk0o74eruGHH8bmCSK2mAjGaf58xnbsEII1\noaEBMYNHH0VIThpfHTsyNncuKosJuZrQ0MBYejqM0I0J9WIi/PesWdiaLozTTex2bPyaMwcxcowv\nFQzU5MmTNxUUFNwp/f9bb731rCfPN2uMiImHeIiHeIiHOR6+GiiXafkhISGlZWVl3QitWSwWC/Pk\n+aYwX7OKBQKBQGAaXBbaSUxMTF+3bt1DRESbNm26Oy4uLsOT5wUCgUAg8AaXtfjq6uoC58yZs+r4\n8eODOnfuXLV69epZRETPPffcG2vWrJnh6Pm+ffueU6X3AoFAIDAsqhaLFQgEAoHAXYx2wotAIBAI\nDIIsBqquri5w1qxZq6OiorInTpy4s7CwcKir5129x8h4oxcRUVxcXEZCQsLWhISErU899dRb2vRe\nfdwdK2vXrn30xRdffN2T9xgRb/QiMu/4InKtWU1NTbtHH3107fjx43dHR0dnbd68Obm+vj7AjGPM\nG62IvBxfvm4DZMy7ZN6VK1fOdzfB12gPb/Sqrq7uaCaNPNHLbrdbkpKSNrdv3/7aiy+++Jo77zHy\nwxu9qqqqOplJI081e//99+c9+eSTbzPG6PLlyz0GDRp0zKz3ME+1Gjx48FFv71+ydNibZF5PEnyN\n9vBGr7y8vPARI0YcnDRpUnpycvJ3NpstQuvPwYtejDGqr6/3f++99x5/4YUXXnf3PUZ9eKPXnj17\nxpp1fLmj2datW+P3798/kjFGlZWVnXv37l1s1jHmjVbeji9ZXHwlJSWhoaGhJUQ4NbdlLlTL54mI\nSktLQ0JCQkqdvcfIeKqXxWJhgYGBdb/61a/+kZ6enrh8+fJnHn300bV2u90UMURXehER+fv7NzT9\nuzvvMSru6uXn52eX/r+ZxxeRa83i4+O3jRw58sDBgwfDJk+e/N1zzz33RklJSagZ72HeaOXt+PLy\n/OzmeJrM6+fnZw8JCSktLy/v6uw9Rsab5OcRI0YUhIWFHSQiGj58+OEePXr8cOHChVv69OlTrP4n\nUBd3k8Gb/t2TBHKj4c1nDwsLOzhy5MgDROYbX0TuafbHP/7xf9etW/fQ3/72t2cTEhK22mw2qxnv\nYd5oxRizeDO+ZJkheZrMGxsbmzlp0qQtZk3w9Uav119//cUlS5YsISI6f/5874qKii69e/c+r3rn\nNcCbZHAzJ5B789lfe+2135l1fBG51mzNmjUzcnNzI2w2mzUhIWGrO+8xKt5o5fX9Sw6fZG1tbeD0\n6dPXRERE2OLj47cWFRX1LSoq6jt9+vQ1zp539DetfatqPbzRq6KiIvj+++9fHxMTkxkXF7c9MzMz\nRuvPwYte0uODDz6YKwX9xfjyTC8zjy93NJszZ86HYWFhB+Lj47fGx8dvTUhI2GLWMeaNVt6OL5Go\nKxAIBAIuMU0QVCAQCAT6QhgogUAgEHCJMFACgUAg4BJhoAQCgUDAJcJACQQCgYBLhIESCAQCAZcI\nAyUQCAQCLhEGSiAQCARcIgyUQOAF5eXlXVeuXLmg5d+rq6s7zZ8//z1PrvXpp59O//rrr++Xr3cC\ngTEQBkog8IIrV650X7FixcKWf1+2bNnzc+bMWeXJtaZPn/7p3//+96fr6+tlKd4sEBgFYaAEAi94\n8cUXXz906NCdb7755iLpb3a73e+rr756MD4+fhsR0YYNG6Z89NFHs4mIXnrppVfPnDnT79ChQ3fG\nxsZmJiQkbL377rs3nT9/vjcR0fjx43d/++2392ryYQQCThEGSiDwgqVLl75w5513Hlq0aNGb0t+O\nHTs2uGmF5vT09MTw8PC9RER5eXlj+/Xrd2bdunUPJScnb966dWvCCy+8sPTMmTP9iIjCw8P3ZmRk\nxKn/SQQCfhEGSiDwAsaYpeXfSkpKQoODgyul/3/gwIGRYWFhB2tqatoFBQXVEhE988wzy+vr6wMe\nffTRtR988MG8rl27lhMRde3atfzy5cs91fsEAgH/CAMlEHiBo0PaevfufV46wO7q1asdr1692pGI\naPfu3ePHjBmTn5GREff+++8/npKSkrp27dpHJ0+e/N3y5cufISKqqKjoIp1SKhAIgAjKCgRe0KdP\nn+KSkpLQt99++5e//OUv3yYiGjBgwKmSkpJQIhil8vLyrt988819paWlITU1Ne0CAwPr4uLiMp58\n8sl3AgMD6xoaGvz/+c9//oKIaO/eveETJ07cqeVnEgh4Q5wHJRDIyB//+Mf/jYmJ2bFr164JsbGx\nmXfdddd2V++x2+1+kydP/u7bb7+9NzAwsE6NfgoEekAYKIFARq5fv97+qaeeeosxZnnnnXeebNeu\nXY2r93z22WePBAYG1k2dOvVLNfooEOgFYaAEAoFAwCVik4RAIBAIuEQYKIFAIBBwiTBQAoFAIOAS\nYaAEAoFAwCXCQAkEAoGAS4SBEggEAgGXCAMlEAgEAi4RBkogEAgEXPL/AZWxhxl6hbgAAAAAAElF\nTkSuQmCC\n", "text": [ "<matplotlib.figure.Figure at 0x52f08d0>" ] } ], "prompt_number": 8 }, { "cell_type": "markdown", "metadata": {}, "source": [ "So we see the probability of finding the atom in the excited state $|c_e|^2$ oscillates at the Rabi freqency $\\Omega$." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Non-Unitary Evolution with the Lindblad Master Equation\n", "\n", "If we wish to include spontanous decay to the environment, we need to follow non-unitary evolution of the system. To do this we need switch from following the state vector to the density operator $\\rho$, and instead of the Schr\u00f6dinger Equation we follow the evolution with the Lindblad master equation.\n", "\n", "\n", "$$\n", "\\dot{\\rho}(t) = \\frac{1}{\\mathrm{i}\\hbar} \\left[ \\mathcal{H}(t), \\rho(t) \\right] + \\mathcal{L}\n", "$$\n", "\n", "Where the Lindblad superoperator $\\mathcal{L}$ is given by\n", "\n", "$$\n", "\\mathcal{L} = \\sum_n \\left[ C_n \\rho C_n^\\dagger - \\tfrac{1}{2} \\left( \\rho C_n^\\dagger C_n + C_n^\\dagger C_n \\rho \\right) \\right]\n", "$$\n", "\n", "The collapse operators $C_n$ are the operators through which the system couples to the environment. In the example of including sponetanous decay from $\\left| e \\right>$ to $\\left| g \\right>$ in our two-level atom, we just have one collapse item\n", "\n", "$$\n", "C = \\sqrt{\\Gamma} \\left| e \\right> \\left< g \\right|\n", "$$\n", "\n", "where $\\Gamma$ is the natural decay rate.\n", "\n", "In QuTip, to solve the master equation we just use the [`qu.mesolve`](mesolvedoc) method as before, but this time pass our collapse operator instead of an empty list.\n", "\n", "[mesolvedoc]:http://qutip.googlecode.com/svn-history/r2510/doc/2.0.0/html/apidoc/functions.html#module-qutip.mesolve" ] }, { "cell_type": "code", "collapsed": false, "input": [ "Gamma = 2*np.pi*1. # [2\u03c0 MHz]\n", "\n", "# Collapse operator\n", "C_Gamma = qu.Qobj([[ 0., np.sqrt(Gamma)],\n", " [ 0., 0.]])\n", "\n", "# Solve with QuTiP solver\n", "ode_data = qu.mesolve(H, psi_0, t_range, [C_Gamma], [])\n", "\n", "for i, state in enumerate(ode_data.states):\n", " c_g[i] = state[0,0]\n", " c_e[i] = state[1,1] \n", "\n", "plt.plot(t_range, c_g.real, 'b', label=r'$c_g$')\n", "plt.plot(t_range, c_e.real, 'r', label=r'$c_e$')\n", "plt.xlabel(r't ($\\mu$s)')\n", "plt.legend()" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "pyout", "prompt_number": 48, "text": [ "<matplotlib.legend.Legend at 0x8b70e90>" ] }, { "metadata": {}, "output_type": "display_data", "png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEiCAYAAACsmUZ+AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJztnXl4jNcXx7+TSOxUglpqqaW1tihtIosslhI7rSCW0v4s\nVbTRWlu6IFpVtLaWEltQqaVCgixC7FRtraXW2CXEEkkmmfv744gGSSYz824zcz7PM8+QmffeM3fe\nec97zj2LTggBhmEYhtEaDmoLwDAMwzC5wQqKYRiG0SSsoBiGYRhNwgqKYRiG0SSsoBiGYRhNwgqK\nYRiG0SSsoBiGYRhNwgqKYRiG0SQFVlCrV6/uMXbs2KnP/l2v1zsFBQUtd3Nz2+vh4ZFw6tSpV6UV\nkWEYhrFHjCooIYSuVatW2/r3779Ep9M9V3Zi6dKlfcuVK3dr7969biEhIWOCg4O/l0dUhmEYxp4w\nqqB0Op2IjIx8e+7cuUOFELpnX4+Ojvbv2rXr7wDg6em568iRI43kEJRhGIaxLwoV5E2Ojo5ZDg4O\nhtxeS0pKcnV1dU0CSJnlZmVlk99rDMMwjG2Rm1FjChYHSbi4uCTfvXv3hWxhjCmhfv0EihUTAARG\njxYQgh95PSZOnKi6DNb04PXi9eL10s5DCixWUP7+/tHh4eHdACAqKqqNt7d3fH7vX7IEuHYNeP99\nYNo04McfLZWAYRiGsUVMUlDZ1tGVK1cq9+zZMwwA+vXrF3r16tVKzZo1OzBt2rTR06ZNG21snFKl\ngPnzgc6dgREjgPBw84RnGIZhbBedVKZYgSbT6UTO+R49Alq2BA4dArZtA7y8FBPFKoiLi4OPj4/a\nYlgNvF6mwetlGrxepqHT6SAs3INSVUEBQFIS4OEB3LgBHDkCVKummDgMwzCMTEihoFSvJOHqCmze\nDKSlARMmqC0NwzAMoxVUV1AAUKMGMHIksHw5cPiw2tIwDMNYhouLC3Q6nV08XFxcZFtH1V182aSk\nADVrAq+/DmzfDugsMgwZhmHU47F7S20xFCGvz2oTLr5sSpcGvvgCiIkBIiPVloZhGIZRG81YUACQ\nkQHUqwcULUoBE46OionGMAwjGWxB2ZgFBQDOzsDUqcDx40BoqNrSMAzDMGqiKQsKAIQA3N2By5eB\nc+eAwoUVEo5hGEYi2IKyQQsKoOCISZOAq1eBdevUloZhGIZRC81ZUABgMAC1agFVqwJxcfLLxTAM\nIyVsQUljQRWo3YbSODgAgwYBY8YAf/8N1K2rtkQMwzDSMXIkBYJJQaNGwMyZ5h+v1+sREhKCcuXK\n4eHDh/Dz80Pjxo2lEc5CNOfiy+a99wAnJ2DBArUlYRiGsV2CgoLQoEEDDB48GGfPnsX+/fvVFukJ\nmnTxZdOzJ+VEXbkCFCsmo2AMwzASYi0uvj179uD999/HiRMnAABXr15F+fLlUahQwZ1rdhUkkZMh\nQ4C7d4HVq9WWhGEYxvbYtWsX/Pz8nvy/UqVKJiknudG0gvLyov2n+fPVloRhGMb2qFChAooWLfrk\n/ykpKdi9e7eKEj2Npl18ADB7NjU1PHQIaNJEJsEYhmEkxFpcfAaDAWPHjkXdunXh+Lh0T+/evbF+\n/XoULlwY8fHxaN++PbzyadYnp4tP8wrqzh2gcmUgKAj4+WeZBGMYhpEQa1FQuXH16lVMnDgRv/zy\nC0aOHInJkyejePHieb7fbvegAKBMGeCdd2gfKj1dbWkYhmFsm71798Ld3R0AcOfOnXyVk9xoXkEB\nQI8ewL171BaeYRiGkQ83Nzdcv34doaGhcHV1VVUWq1BQLVtSO461a9WWhGEYxrYpWrQoxo0bB2dn\nZ/j4+Kgqi1UoKGdnoFMnYMMGaslh85w6BcyYQVnKiYlqS8MwjB0xYsQI7Ny5E+fOnUNAQICqsmg+\nSCKbTZuADh2AiAigXTuJBdMC584Bv/wCrF8P/PPP0681bgy0bw8MGwaUL6+OfAzDFBhrDpIwFbsO\nksimVSugVCkbdfPt2wc0awZMnw689BLw44/AxYvAiRNASAhQvDgweTLg5kbWFcMwjB1gNRYUAPTt\nS5bUjRtUp88m2LYN6NIFqFAB2LoVqFEj9/ft309WlMEAbNwING+urJwMwxQYtqDszIICKNz8zh0g\nOlptSSRi7VogIACoWRPYtStv5QQAb74J7NlDcff+/twsi2EYm8eqFFSrVkDJksBvv6ktiQT8/jvF\nz7/5JrBjB1lQxqhZE9i9m+rrd+tGlXQZhmFsFKty8QFUUWLLFuD6dSt28926BdSrB1SvTsrJ1FLt\nqamk2JKSgOPHAZVzFRiGeRp28dmhiw8gN19yMhAbq7YkFjB8OJCSAixZYl4fkWLFgOXLSUENGgTY\nyQ+BYRj7wuoUVJs25Oaz2mi+DRuAVauACROA+vXNH6dRI+Drr4HwcGDZMunkYxiG0QhW5+IDaPtl\n3z7g8mVAZ5EBqTB375Jrr1w54MABykC2hKwswNeXekcfPUouQ4ZhVIddfHbq4gMoUffKFdp+sSqC\ng4GbN4Fff7VcOQGAoyOwdCn9u18/dvUxDGNTWKWCevttet6yRV05TCIhgRTTqFHAG29IN2716sD3\n3wPx8Rx6zjCMTWGVLj4AeP11wMXFioIl2rUjt97Fi+YFRuRHVhbQoAFZVH/9Rc8Mw6gGu/ikcfFp\np/m8ibRtS4bDvXtUAknTHDlC5t4330ivnABSSF9+SXlVq1cDvXpJPwfDMNIxciRdF6SgUSNg5kyz\nD9fr9QgJCUG5cuXw8OFD+Pn5oXHjxtLIZiFW6eIDyCDJzAS2b1dbkgIQEkKhhx9+KN8c3bsDr70G\nTJxIC8MwDFMAgoKC0KBBAwwePBhnz57F/v371RbpCVbr4tPrgbJlgXffpSLgmuXMGaBOHdp7mjZN\n3rk2bqS+JAsXAgMHyjsXwzB5Yi0uvj179uD999/HiRMnAFC79/Lly6NQoYI71ziKLxecnKj00ZYt\nGg9e+/ZbEvbjj+Wfq0MHqor+1VdAerr88zEMY9Xs2rULfn5+T/5fqVIlk5ST3FitggJoH+rKFeDY\nMbUlyYMrV4DQUGDAgILV2rMUnY72uS5dIiuKYRgmHypUqICiRYs++X9KSgp2796tokRPY7UuPoCu\n/y+9RFs8o0dLNqx0BAcDs2aRm+/ll5WZUwigRQtqgHjhAqChuyE5efCAYkVy/NYYRjWsxcVnMBgw\nduxY1K1bF46Po3979+6NJUuWoFq1akhOTsY777yT7xhyuvisWkEBFMDywgtAXJykw1pOSgpQuTLQ\nuTPVzVOS7L2otWup7IYN8ugREBND3/uOHcChQ6Sg3ngD8PSkR5s2QJEiakvK2CPWoqBy4+eff8aj\nR4/w5ptvoly5cqhVq1a+71d1D0qv1zsFBQUtd3Nz2+vh4ZFw6tSpV3O+LoTQDRkyZF6LFi12vPXW\nW/vi4uJ8LBHIVNq2pRzYlBQlZy0Aa9YADx8CH32k/NwBAUC1asBPPyk/twIcPkx5cO3bA7NnkxIa\nNw745BNSUrNn031Bkyb0XoZhCs7x48fRpUsXuLu7o5TaOTxCiHwfCxcuHDhy5MgfhBCIj4/3CggI\n2JTz9a1bt7bq0aPHKiEEzp49W/O11177K6+xaDpp2bFDCECItWslH9oymjcXom5dIQwGdeYPCaGF\nOX5cnfllwGAQYuZMIZydhahcWYgNG4RITX3+fWlpQqxfL0SlSkIUKiTE5MlCZGYqLy9jv8hxrVOK\n/fv3i4ULF4qYmBhx6NAho+/P67M+/rtRHZPfw6gFFR0d7d+1a9ffAcDT03PXkSNHGuV8vVChQpn3\n798vKYTQJScnu5QqVeqeDHo0T9zdgRIlNNZl99Qpaiz43nvqVbMdOBAoXBiYM0ed+SXmzh3yWo4c\nSa67v/4COnbMfc+pcGF677FjQNeuwPjxgLc3lUFkGCZ/mjVrhoEDB8LX1xdNmjRRVRajO+hJSUmu\nrq6uSQDtIel0uqecjc2bN9997dq1inXq1Pnn6tWrlWbNmjUiv/EmTZr05N8+Pj7w8fExS/BsnJzo\n4hMTY9Ew0hIaCjg4UHdFtShbFggMpGKyU6cCpUurJ4uF6PWkaBISKObko48KpvddXKizSadOwPvv\nk+czLg4oXlx2kRnG7oiLi0Oc1MEAxkyswMDAsISEhOZCCBgMBl2VKlUu5Xx90qRJEydMmPC1EAI3\nb94sV7169fMpKSmlchsLMpm906eTNysxUZbhTSMzk/xP7dqpLYkQ+/fTwsyerbYkFjF4MH2MZcvM\nH2PjRiEcHIQICBBCr5dONobJDbmudVokr88KJVx8/v7+0eHh4d0AICoqqo23t3d8ztczMjKcy5Ur\ndwsASpcunVKkSJG0Z60sufH3p2dNWFHbt1P8+3vvqS0JJe2++SYwd67Gs5nzZu5cYP584LPPLDNI\nO3SgsSIigKFDrXY5GMa+MKbBMjIynAIDA8OaNm16wMfHJzYxMbFyYmJi5cDAwDAhBJKTk8t07tx5\nnY+PT6y7u/vuZcuWBeU1FmS6q8jKEsLFRYj+/WUZ3jQCA0mYtDS1JSFCQ8n82L5dbUlMJjpaCEdH\nsnqkCnIYN46W4+uvpRmPYXJDrmudFsnrs0ICC8rq86Cy6d6dullcuKBil927d6lixAcfAD/+qJIQ\nz5CWBlSpQsm7a9eqLU2BuXyZctwqVAD27JGuYr0Q1Ntx2TIKrMlR5YVhJMOa86BMhWvxFQA/P6rw\nc+6cikKsWkU18Pr3V1GIZyhSBOjdG/jjDwqFsxI++QRITQXWr5e2nYpOByxYANSuTYETDx9KNzbD\nMNJiUwoKUHkfKjSUGgeqHJr5HEFBQEYG8NtvaktSILZvJ2Nv/HhSJFJTtCiVKjx/nuZgGKkpU6YM\ndDqdXTzKlCkj2zrajItPCKos1KIFEBYmyxT5c/kyULUqMGUKMHasCgLkgxBAvXpAuXLUGl7DZGRQ\nW6vMTOD4cXlLFQ0bRoETu3YBzZvLNw/D2CPs4suBTkfRfDExKkVorV9Pz127qjC5EXQ6oE8fYOdO\n2qTTMD/8QHnO2SWM5GTqVNqeGziQtuoYhtEWNqOgAHLz3bwJnDypwuTr1gF16wKvvmr8vWqQ3QZe\n6cK1JpCYCHz9NVWIaNdO/vlKlqRml//8Qy20GIbRFjanoAAV9qFu3ybXmRatp2yqV6eSG8uXazYJ\nKDgYyMoCZs5Ubs7WrSmqb/p02pNiGEY72JSCqlYNqFFDBQX1xx90Ze3SReGJTaRPH/KfHTyotiTP\n8eefVAB+9GjlWmdlM3kyVUGfOFHZeRmGyR+bUlAAWVFxcaQvFGPdOgqQ0Fr03rN0706VVJctU1uS\n55g6lcLJP/5Y+bkrVwaGDyfjUrPdmaUiO/pk6VL6wEePUmQKw2gQm1RQd+8CR44oNOH9+8DWrWQ9\nqZYhXEBeeIFq/qxaRRVYNcKpUxRWPmyYejVtx4yhuW0y7Fyvp95g7u50F9CwIfk1+/ShxlrFi1Po\n5HffcbQIoylsTkG1aEHPikVTR0ZScq7W3XvZ9OkD3LpFSlUjTJtGEXsj8q2DLy9lypB78Y8/KOzc\nZti+nUpyfPQRWU+DB5MFffIkcOIE5WSMGkUL8NlnFOSzbBlgMKgtOcMYr8Un5QMK1aeqWVOIzp0V\nmUqInj2FKFfOejripacLUaaMEH36qC2JEEKIixepqeDw4WpLIsTDh0JUrCiEh4d6fSYl49o1Ibp0\nocKDNWpQB0djHyo6Wog33qBjGjUS4sQJZWRlbBIoUc3cGvH2ppQf2W8C09OBTZsoLtrRUebJJMLZ\nmRokbdyoib2H6dPpOThYXTkAoFgx4IsvqO/U5s1qS2MBly/TjyAykhLHT5yg79yYC9rPD9i/H1i5\nErh2DfDwsDFzkrE2bFZBJSUBf/8t80QxMbQHpeXw8tzo3h1ISVG9DfHNm1RyqE8fijHRAgMHAjVr\nUkSfRqPx8+fffwEvL1rc6GiqamJKxrODA9CzJ7B3L/Dii0DLlsDvv8snL8Pkg00qKC8vet65U+aJ\n1q2jbM/shlTWQsuWtFmucnXzWbNoT370aFXFeAonJ9qKOXSIokGtin/+obuzBw/o5snd3fyxqlcn\nU7JJE7qhmTNHMjEZpqDYTC2+nAgBvPQSBUysXCnjJFWrUkPA8HCZJpGRoCBgyxbg+nW6KivMo0cU\n3u3np7qefI60NMqpa9KElsgqOHeOFJJOR4ERDRpIM25qKllUGzdSMeS+faUZl7F5uBZfHuh0dCMZ\nHy+jm+bkSarN07atTBPITPfuQHIysGOHKtOvXUvdPz78UJXp86VIEcqLioykNCHNk54OvPsu7Snu\n2CGdcgJoY27tWrqT+OADYPdu6cZmGCPYpIICyM135YqM5Wuioui5TRuZJpCZNm0o/0Ul8yW7J5OP\njyrTG2XIEFqe7CAOTRMcTD7J0FB5akE6OVGrlipVKJ3i0iXp52CYXLBZBeXtTc+y7UNFRlILiypV\nZJpAZooWBdq3pw1wRctuUFBZQgLwv/9pN7fZxYUMhrAwCorTLKtX0/5QcDBFk8qFiwsliaWl0TwP\nHsg3F8M8xmYVVL169JuSJWE3NZUGfvttGQZXkG7dKGlX9miSp1mwgKLdtdR4ODdGjiQXsZLFa03i\nzBnSou7uVCtKburWJYV47Bh9eVYZ5mgaGRnkzf/9d3I2HDnCullJbFZBOTgAnp4yKagdO8jvb63u\nvWzatiVLSsEgj9RUKgPXrRtQtqxi05pFtWpAYCDw889UPktT6PW07+TsTEpDqUCXt98GQkLonFm6\nVJk5FebKFeDTT8lbWqwYUL8+na/vvAM0bkyBuxUrUkqCYiXV7BSbVVAAufnOnqWcQ0mJjKQLe7Yf\n0VopUYKUVHi4YqVt1qyhFKxBgxSZzmI+/ZTumH/+WW1JnmH2bLo6LlyovJs5OJg2eUeMoEAhG+HE\nCeC996ia/g8/AK+8QjUaly+nBgCHD9P5O2UKZZasWkUKq0ULRX9C9oWlpShMeUChUkfZ7N9PVVtW\nr5Z44FdeEaJtW4kHVYmVK2mRdu1SZDo3NyHq1LGuUkI+PkJUr66halaJiUKUKCFEhw7qyXD2rBDF\nignRpo11fZm5kJEhRHAw/QyKFhXio4+EOH/e+HHJyUJ8950Q1arRsb6+Qly6JLe01gMkKHVk0wpK\nrxeieHEhPvxQwkHPnaNlmzlTwkFV5O5dKob32WeyT/XXX7R0M2bIPpWkrFlDcm/apLYkjwkMFKJI\nEToX1eSnn2hhfv5ZXTks4PJlIZo3p48xZIgQt26ZPkZmJi1B8eJCvPCCEKtWSS+nNcIKqgC0aiVE\nw4YSDjhvHi3bP/9IOKjKtGxJZo3MfPSREIULC5GUJPtUkpKRIUSFCkK0a6e2JIIKugJCfPml2pII\nkZVFZkOJEkJcuKC2NCYTFSVE2bIkvhRK5cwZ8hAAQgQFCZGaavmY1owUCsqm96AAcpUfO0ZJoZIQ\nGUllYF55RaIBNUCnTlQm5/Rp2abQ6ylku1Mniq60JpycKCR+yxYq2KAaGRnUNKtGDarHpDYODsCv\nv9K/Bw2yqqi+FSso3qNiRdpf6tHD8jFr1aKA2EmTaPx27ahUJ2M+Nq+gPD3pec8eCQbLyKAaZ23a\naDeBxxw6dKDnP/6QbYpt24Dbt4HevWWbQlY++ICuxwsWqCjErFlUAXn2bNMKwMpJ9erA119T4rqM\n54+UrFtH/Rp9fKgmrpS5zYUKUaHh5ctJWbVuLeHNsT1iqQlmygMquPgePBDC0VGIsWMlGCwujuz3\ndeskGExjvP66EN7esg3fq5cQLi7Ujspa6dpVCFdXIR49UmHy5GQhSpUSon17FSY3QkaGEHXrUiO2\ntDS1pcmXLVuEcHISwt1diPv35Z1r3TohnJ2ptdbNm/LOpUXALj7jFC9ORT8TEiQYbOtW6vvk5yfB\nYBqjY0fq/XP7tuRDP3wIrF9PeSTOzpIPrxhDh1Ibl99+U2Hy778H7t2jGGet4eRE2cz//kvx2Rpl\nxw6q1FS/PvX7KlFC3vk6d6Yau6dOkbXGlpTp2LyCAqjv2v79EvTni42l6uWlSkkil6bo2JESOWTo\n1LdhAyXo9uol+dCK4udH7qC5cxWe+PZtcu+9+y7QsKHCkxeQ1q3pHPrmG+DqVbWleY4zZ0i8l1+m\n+8wXXlBm3jZtgIgImr9bN030CLUq7EZBpaVRop3Z3L9PWs7XVzK5NEWTJkClSnTLJzErVlAuafZ+\noLWi01ER2b17LTyXTGX6dDJDJ05UcFIzmDGDomHGjlVbkqd49IiK9xcqRDFO5copO7+vL+VTx8bS\n+WNFsSQFRwi6Mdm9m3ocSWXpW+ojNOUBFfaghBDi6lXaOpo+3YJBNm+mQbZvl0wuzTFoECVzSLjJ\ncvMm7QGOHi3ZkKqSnEwpSEOHKjThjRuUENurl0ITWsiYMfQ72bNHbUmeMHAgibR5s7pyfP45yTF1\nqrpySMrFi0JMnkxpKqSmnjzAeVAFp0YNIbp0sWCAUaNox9OWkxsiIuiU2LJFsiHnzKEhjx6VbEjV\n6dVLiNKlFToVgoOFcHCwnry7+/cpaczLSxMVJpYsofNv3Di1JaHl6NmT5PntN7WlsZD4eMqBy1ZI\nXl5C/PADXUNOnhTi4UNWUKbQp48Q5cpZ8Jt54w1Zo9w0waNHZEENGSLZkM2bC9GggWTDaYLsXNnl\ny2We6OpVMtf69pV5IonJviuR8EbHHI4do9JFPj5UVUYLPHpEEYTFiwtx6pTa0pjBtWt0MQWEqFJF\niK++yrOiCSsoE1iwgD7t6dNmHJycLIROJ8SkSZLLpTm6dhWicmVJ7n6zq0LZlEtDUAGFl1+mG0hZ\nGTmS/KNnzsg8kcSkp9MCNW5Mi6WSCPXrC/Hii6TntcTly5Ry0aiR5qPy/yMrS4jZsynVwdlZiPHj\nhXj4MN9DpFBQdhEkAVCgBGBmuHl273hbDZDISYcO1G9Agj4Cq1fTc8+eFg+lKRwcgAEDaNP7339l\nmuTuXdpZDwykEgXWhLMz8OWXwJ9/KtrKJSchIVSdfOFCqhahJV56CViyhH5in36qtjQF4OFDKrUx\nfDjg5kaleb75hnqRyI2lGs6UB1S0oLKyqJDj+++bcfCIEeQrsJrbHQu4fp3Mnq+/tniopk2FeOst\nCWTSIJcv09bQ+PEyTfDtt/Q9HD4s0wQyk5kpRL16Qrz6quL+tRMnKBk3MFDRaU1mxAj6itevV1uS\nfLh0iSxhnY6izEzwrEACC0pH4yiDTqcTSs73LO3b0x3v33+beOBrrwEvvkj1euyBN9+khGQL6kNd\nvEhVcKZN00bZODlo1w44epQ+q6OjhAPr9VRvr3ZtKq1lraxbB3TtSvX63ntPkSkNBkpnOHWKfufl\nyysyrVmkpwPNmwPnz5M1VbWq2hI9w969lG2cmkrNr9q1M+lwnU4HIYRFNeHsxsUHkJvvn3+oGkCB\nuXWLTFpbrB6RFwEBwL599NnN5Pff6blbN4lk0iADB5I3NCpK4oHXrKFGgMHBEg+sMJ07A82aUfXU\n9HRFppw7l+6rZs7UtnICgMKF6bqv1wP9+2us4WF8PF3zihcnRWWicpIKu1NQAOWSFZi4OHq2h/2n\nbNq3pz23LVvMHiI8HHj9daBmTQnl0hgdOlDS56JFEg4qBJU1qlOHuh1bMzodJWxeukSbQTJz6RLl\nCLduDQQFyT6dJNSuTV93bKyGujbv20c3qdWqkbavV081UYwqKL1e7xQUFLTczc1tr4eHR8KpU6ee\nq/07bdq00Y0bN/6zadOmByMiIgLkEdVymjWjsmG7dplwUGwsULIk0LSpbHJpjsaNgQoVqEaLGVy7\nRjcBtmw9ARQL0KcPFd+wwNh8mh07KLjgk08oGsPa8fenO8Np02Sv8/Pxx2SFLFhgXc0GPvgAaNmS\nAiYuXFBZmD//pD4kL74IREerb4Ya26RauHDhwJEjR/4ghEB8fLxXQEDAppyv79+/v9kbb7xxUK/X\nF7px40b5+vXrH89rLKgYJJGNm5sQHh4mHPDqqxrpVKcwAwZQNmpGhsmHzp1Lm7/Hj8sgl8Y4dkxI\n22C5fXtK2LOlhPDISCF3593YWJrim29km0JWLlygxon+/irmNx8/Th0cq1SRpAEllMiD6tmz58r4\n+HgvIQQMBoOucuXKiTlfnzRp0sTZs2d/lP3/48eP189zMg0oqOBg6upaoIC8K1eE5TWSrJTwcPrs\ncXEmH+rvT3pdA4UEFKFJEwp0sph//qE1nzhRgsE0hMEgxJtvClG9ulk3PMbIzKT1r1rVuvX6/Pn0\n9S9YoMLk168L8dJLVAVEorw7KRRUIWMWVlJSkqurq2sSQFF4Op3uqTC8a9euVbx161a5du3abU5N\nTS02bNiwn+rXr38ir/EmTZr05N8+Pj7w8fExw+4zHw8P8vkePgy4uxt5844d9KywjJqgVSvyh0ZE\nAC1aFPiwpCTaths92rrcLJbQrx8wYgTF0lhUbPynn8hvOHSoZLJpAp0O+Pxz2rRbsYIiAiRk6VLy\nTK1cCRQtKunQivK//1F8THAwedkUi+rLyKBquklJtP9hZt5dXFwc4rL37KXCmAYLDAwMS0hIaC4e\nW1BVqlS5lPP1UaNGfTd48OB5QgjcuXPnhapVq15MSUkpldtY0IAFlZ3m8913BXjzkCFClCxJt2j2\niL8/NaIzgUWLaH0PHZJJJg1y6xbl3QQHWzDIgweUpd+7t2RyaQqDgcyc2rUlzYvKLv3n5mYbFvv5\n81QbuEMHBT/PoEH0ow0Lk3RYKFFJwt/fPzo8PLwbAERFRbXx9vaOz/m6u7v7ntKlS6cAQNGiRR8V\nLVr0kYODg5YCJp/ixRcpsqxAkXzx8WRySZrkYkUEBFAyyfnzBT7k998p/6lxY/nE0hply9JSLV8O\nZGaaOUhYGDUkHDJEUtk0Q7YVdebMfyVGJCAkBLh+nfok2oLFXr06FeH44w/qoyY78+dTVMmYMVS1\nRGsY02AZGRlOgYGBYU2bNj3g4+MTm5iYWDkxMbFyYGBgWPZ7Pv744xm+vr4xHh4eu1atWtUjr7Gg\nAQtKCKoT5xDQAAAgAElEQVR1WL68kTuU27fprmLyZMXk0hynTtEa/Phjgd6ekkJluj75RGa5NMi6\ndbRUERFmHGww0EZWgwa2YQbkRVYWfca6dSXxSly8SLV0e/aUQDYNkZEhxGuv0ZbQvXsyTrRzpxCF\nClEQmAxeInCxWPPI3ozMdy9w/Xp6U3y8YnJpklq1hHj77QK9ddUqWrJdu2SWSYOkpwvh6irEu++a\ncfC+fbRwc+ZILpfmyD5JwsMtHuq99yjg6eJFCeTSGHv2UHWhjz+WaYKkJNKAtWsLcfeuLFNIoaBs\nINHCdJo3p+d83Xw7d1Kqd7NmisikWQICKBcsNdXoWzdupMRVNzcF5NIYzs7U0n7DBuDOHRMPnj+f\nMvatJbvUErp3p034qVMtai17+jQQGkoeUc2VCJIANzcKmpg1iwJAJEUIGvzGDXItly4t8QTSYZcK\nqn59oFQpI5XNd+6kmnRFiigmlyYJCKAyNbGx+b5Nrwc2b6YiFPa6ZdevHy3VmjUmHHTnDtW7CQqi\nk9LWcXSk4owHD1pUZ3DiRPppaqy7vKRMnUr7m4MGAVlZEg78669U6mXyZOCNNyQcWHrsUkE5OFCI\neZ4W1IMHwKFDgLe3onJpEm9vKqu/eXO+b9u1izpEdOigkFwapEkTuvkJDTXhoNBQ4NEjYPBg2eTS\nHH37Ug+MkBCzDj96lHT6iBHqFzqQkzJlgBkzgAMHJCyndeoUtc3w97eKWo92qaAACs47cYIuqs+x\ndy/dsnh5KS6X5ihcmOqwbN6cr0tm40Z6a6tWCsqmMXQ6sqL27KFgNaMIQe49NzegUSPZ5dMMhQtT\nXaLt28mSMpEvviCvlFX0UrKQXr0oDXHsWBOLXOdGRgYNWLQoJY9ZQSkt7UsoEx4edH3ItaNEfPx/\nZhZDlYwvXKBS8LkgBIXF+vsDJUooK5rW6N2bTp1lywrw5vh4uqO1J+spm0GDgBdeMNmKOnCA9vmC\ng8nCsHV0OuDHH4GUFGDCBAsH+/prqlCwaBFQqZIk8smN3Sqo7JZHubr5du6kRB572BMoCNlVtfMo\nHvv339Rnq2NHBWXSKJUqkcG5bFkB2icsWkTn2DvvKCKbpihVCvjwQ0qcy+PGJzcmTABcXYGRI2WU\nTWM0bAgMG0bpSocPmznIkSN0M9CvH9Cpk6TyyYndKqgSJagdxHOBEunp5OJj995/VK0KNGiQ5z7U\nxo303L69gjJpmL59yeDMt2p+Sgqwdi3Qs6cyrbO1yPDh5O777rsCvT0hAdi6lXJKS5aUWTaN8eWX\nFCE7bJgZfaP0emDAANLsM2bIIp9c2K2CAijcfN++Z7L/Dx0C0tI4QOJZAgLIsrx377mXNm6kYKDK\nlVWQS4N07kw3QPm6+VavpuCIgQMVk0tzlC9Pn3/ZMmrQaIRvvqGLtK2VKiwIpUsD335LWxJLl5p4\n8HffUaz6vHmAi4ss8smFXSsoDw9K7/nrrxx/jH9cycnTUxWZNEu7dqTJt29/6s83b5LBye69/yhe\nnHphrVlDOihXFi0iq9Se+ozlRnAwmQQzZ+b7tgMHgMhIeru9Gpx9+tBN9ejReQR35cbff5P59e67\nQJcussonB3atoLITdp9y8+3cCdStS7dqzH+4u9Nt3DNuvogICpJgBfU0ffuSsZnt/nyK48eB/fvJ\nerCFAnKW8PLLdPFcsCDfq+7kyRQUYaulCguCgwMFTNy6BXz1VQEOyMoi117JknSgFWLXCqpqVeCl\nl3IoqKws2jjg/afncXKiXtrPhJtv3AhUqUL7ecx/+PjQuuTqjlm8mNazd2+lxdImn31GuYfz5uX6\n8tGjFLk3YgTHLTVpQh14f/wROHnSyJsXLCD3xqxZVpswZtcKCiA3X0LC42vu8eN028vuvdxp1476\nuT/2iaal0aZ1+/ZsCDyLgwMVh4iKoooyT8jIIK3VsSNb6dk0akQ3P7Nm0Un1DN98Q0bA8OEqyKZB\nJk+mPc7hw/NJTbx+nZKnWrak3CcrhRWUB3DlCnDpEv4zpVhB5c7bb9PzYzdfXBzt4dlz9Yj86NOH\njPKwsBx/3LQJuH3bvoMjcmP0aNLkz5icf/9NwY7DhtlH3lNBKFuWUpqio4F16/J40yefkLKfM8eq\n7x51woKCjSZPptMJJecrCIcPUwTaihVAr4jeVHPuyhWr/lJlpWlTCg1OSMBHH9Fef1KSdXcylZNm\nzSi25EnBz4AAskAvXrTfooW5IQQlJ969S3lRj9emTx9KlbpwgQ3OnGRmUqrm/fukxJ/6/W3fTiVd\nJk4EcnQwVxqdTgchhEUXUru3oF57jczlhATQ/pOHByun/AgIAPbuhbidhIgIqh7Byilv+vWjHMmj\nR0Hu0chI+iMrp6fR6ciKOnv2iVlw/jxZn4MGsXJ6lkKFaB/q4kUKP39CWhrF4deqRQljVo7dK6hC\nhagU2umYRPLzsXsvfwICAIMBV3+NxPnz9F8mb3r2pHiI0FBQy12DgRQU8zxdutCFddo0QAhMn057\neVZQ01QVfHwoADIkhCxMALR2Z84Ac+faRCcGu1dQABlNrv8k/PcfJm+aNgXKl8e9MCp71K6dyvJo\nHFdXCiJZsVxAhIZSuP4rr6gtljZxdKQKsAcP4k54DH79lcL1OQE8b7KV+CefgOqNTZ0K9OhhM1Wb\nWUGBjCYP7EJm4WIcL20MBwegbVtUPh6JRg0ybbJZnNT06wdUvnkYuhMn2HoyRt++QIUKSPo0BOnp\n9lGx3BKqVAHGjyev6M2gj8lct7JyRvnBCgrAW28BHkjApYpv0RfM5MtD3wCUyryDQa/vVVsUq6Bt\nW2BQkVBkOBQmnwyTN0WKIG3Ix6h1YTs+9T2IV19VWyDtExwMvF8xAuX3/oHMcZ9bTaXygsAKCkBJ\n3Mfr+Au7BO8/FYStaA09CiFAl3t1c+ZpnJGBXliJjeiIO+BYaWPMx2DcRWmM1k1TWxSroLBIw0yM\nwCm8gpnCtsq8s4ICgH374AgD1l73gF6vtjDaZ0Nsaewt5ImX/mIFVSC2bEGJtCQsNvQzrR28HZKW\nBkybVwqbq38Il5hw6pfF5M+MGSh+7V+seOtHfDnVGVeuqC2QdLCCAoBduyB0OuxId3u6cCzzHAYD\nsGULcKlBAHTHjj3OcGbyJTQU4sUXkVivjWnt4O2QpUupCEKV70xrxWG3XL5MpSW6dEH/la2h19tW\n1CMrKABISIC+TkPcQ+nn+0MxT3HwIFUwL9XzcXx5Hj2imMckJQGbNkHXuzeC+hfCnj3A6dNqC6VN\nsrJIHzVtCnh2e5EKnS5dCpsyCaRm1Ci6a5wxAzVqUHWj1auBbdvUFkwaWEFlZgJ798LZ1xNVqxpp\nMscgIoIC+ZoPqEOVqPPosss8JiyMGsb164egIFo7tqJyZ906ytMdPfpxrnz2xff779UWTZvExlJP\nlzFjgOrVAdDa1a5Nubq5lDW0Ouy+1FHOWke9I3pxpSMj5Kh0BK51VACeqXXElY5yJ7vSUUoKle55\nsjZc6yh38ql1pJFKR1zqSBIS/kvQ9fCgajRPsrKZp7h+nRoOP6keERBAHfni4tQUS7v8/Tf5RHPk\nPg0cSDdAW7eqKJcGiY2lpRo16hnFPW4cnWM//KCabJpk/nzqvjBjxnM3hy1bUgWTqVOt353MCioh\ngVLVq1Z9UkSC3Xy5ExlJz0+qR/j4UHvTP/5QSyRts2wZXW179nzyp/btyRBYtEhFuTTItGnAiy9S\nnu5T1K0LdO8O/PQTcOeOKrJpjtu3gc8/J02UR5fcbL314Yf5tOSwAlhBJSQ8KRDboAE1jWUFlTub\nNwMVK+YotlGkCPXx2bTJun8FcmAwkIJq04auvI9xdiav1caN1BmVoWK6W7dSQ8Jcy8dNmECurNmz\nFZdNk4wfTw0eZ83Kcy+iQgVgyhRy961cqbB8EmLfCurSJSAx8Un9PUdH+ufOnSrLpUH0erqItGv3\nzG+iY0cKdeX4/KeJi6Nz6zmTAHjvPVrPFSuUF0uLfPstdRTIs537a68BnTrRBfnePUVl0xyHDwO/\n/EL7v/Xq5fvWQYOo9OOIERR5a43Yt4LavZuecxSI9fKirQO+u32aPXtoA/u54rABAaSxNm5URS7N\nsnQp9Sfv2PG5lxo0oICARYvY8Dx/nsKiBw0CXnghnzdOmEAuvjzawtsFBgN1bixXjiIgjODoSOfY\n/fukz6wR+1ZQCQm0h5KjQKy3Nz2zm+9pNm+m1iQtWz7zQvny1K+EFdR/PHwIhIdT3b08ohsHDqQ9\n7oMHFZZNY3z/PV1IP/7YyBubNqWOzt9/T+trjyxbRneK06bRXkQBqFuXdNmaNRQMaW3Yt4LavZsq\nxRYq9ORPTZuSHzw+XkW5NEhEBFmXpUrl8mLHjhTexwmVxPr1tEfQp0+eb+nRg3SXPQdL3LxJn79P\nnwK21JgwgVwb9mhFpaQAn31GPrtc3Mb58emnFJE+dCiQnCyTfDJhvwrqwQPaN3mm/5OzMxkEvA/1\nH5cu0d1+nr2fst1YmzYpJpOmWbqUEifzaX5ZujQFp4WFAampyommJWbPBtLT6bpbIDw8KCgnJIT8\nVvbExImknH/6ibK9TcDJCfj1V0pX/OQTmeSTCftVUPv3U22V5s2fe8nLi/Iq7X0/NpstW+g5z+65\ndesCNWuymw8Arl6l0Kk+fYxeSAYOpHPst98Ukk1D3LsHzJkDdO0K01pqfPMNXWlnzpRNNs1x7Bgp\npkGDgCZNzBqiUSMqOBEaCmzYILF8MmK/CiohgTb33d2fe8nbm/Yj9+xRQS4NsnkzGQR16uTxBp0O\n6NABiI4my9SeWbGCTp583HvZeHvTxXn+fAXk0hgLFgB379JF0ySaNaOIvunTrc9fZQ5CUGDECy9Q\nUVgL+Pxz0m8DB9J9lDVg3wqqfv1cQ4fc3Gjjlt185ILZvj2X8PJn6diR3mwrVSrNQQi6RXV3p4Jo\nRtDpgMGDgb17KRfIXkhLo0TSli1pz9dkvv6aXHzTp0sum+ZYuZI2xKdMAVxcLBrK2ZmGe/SItrEM\nBolklBH7VFDZ5lEu7j2AcjLeeIMDJQBag9TUfPafsvH0JGVvz1UlDh8GTGzr3q8fBUvYkxW1bBmV\nzTLZesqmYUOKMpk1y3oTfArCnTu0adSsGZk9EvDqq7T3Fx1tHTV47VNBnThBTvBnAiRy4uVF21S2\nUBHYEiIiqDisr6+RNzo5kRbbtIn29uyR0FBarB49CnxImTJAYCCwfLl97HlmZVFibtOmgJ+fBQN9\n+SX9OKdOlUw2zTF2LJU1WrBA0srCAwYA3bpRmcNDhyQbVhbsU0FlJ+jmYUEBpKDS04EDBxSSSaNE\nRJByKlasAG/u2JEijfbulV0uzZGRQf6TTp2MZJw+z+DBlNqzfLlMsmmItWuppcbYsRZ2DHjlFTI/\n582zzaaZe/aQYhoxgmLEJUSnA37+mcoh9exJe4FaxT4VVEICJZjWrJnnW7IjhO15H+r0abqY5Bm9\n9yxt25Kj2xozAi1l82aKLjPBvZdNs2a0eT1vnm1XljAYKAivbl3S4xYzaRJdbc32FWoUvZ7uWl56\nCfjqK1mmcHGhFIfz54HevbXr9DCqoPR6vVNQUNByNze3vR4eHgmnTp3KNSjUYDA4uLu774mKimoj\nvZgSs3v3kwKxeeHqSjEU9qygsnsRFlhBlSpFjWh+/922r7S5ERpKt6StW5t8qE5HdeiOH//PuLdF\nNmygzzh+vEQeq6pVqb95WJhtWe2zZgFHjwI//kgb4jLh6UnR65s3Uw60JhFC5PtYuHDhwJEjR/4g\nhEB8fLxXQEDAptzeN3PmzBFlypRJjoqKap3XWDSdyly/LgQgxHffGX3rkCFClCwpRGamAnJpEH9/\nIerWNfGgRYtofQ8dkkUmTXLzphCFCgkRHGz2EA8eCFGqlBC9e0sol4YwGIRo3FiI2rWF0OslHPj+\nfSEqVBDCzY0msXbOnxeiWDEhOnRQ7PMMGkQ/2ZUrpR338fXeqI7J72HUgoqOjvbv2rXr7wDg6em5\n68iRI42efc+lS5eqRkZGvt2pU6cNwsIOirKTS4HYvPDyomhWewoBzub+fYrga9/exAM7daLb4/Bw\nWeTSJGFh1OHUDPdeNsWLU+jvb7/ZZmBaRAQlv48b91RlMcspUYJCsPfuBVatknBgFRCCovUcHcm0\nUait9+zZZE0NHEiBqFrCqIJKSkpydXV1TQKoZbtOp3vOdzN8+PDZM2bM+CT7PfmNN2nSpCePODU6\nsSYkUKRVATKyfXzo2R4bxm7bRq7wArv3snF1pYULD7cfN19oKG1kN2xo0TDDhlGsxdy5EsmlEYSg\n1KXq1Wm/Q3L69aP1Hz2aknyslZ9/BmJiKL+ralXFpnV2puCVsmXp937mjHnjxMXFPXV9lwRjJlZg\nYGBYQkJCcyEEDAaDrkqVKpdyvr5s2bKgCRMmfC2EQP/+/RdHRka2yWssaMHF5+YmhIdHgd/+6qtC\ntGsnozwaZcAAIUqXFiIjw4yD584ln8Hx45LLpTmOHaPPOnOmJMO1by9E2bJCpKZKMpwmiIykJfr5\nZxkniYujSb75RsZJZOTCBSFKlBCiZUvVXJUnTtC5V6UKiWMpUMLF5+/vHx0eHt4NAKKiotp4e3s/\nlb66a9cuz9jYWF9fX9/YyMjItz/77LNvd+/enXf8tpo8ekSB/wVw72Xj50euLr1eRrk0hsFAG6dt\n2lB6k8l07kzuCXuI5lu8mHxWvXpJMlxwMKW+2ErIebb1VKWKRR5Q47RoQYX9pkwBLlyQcSIZEAJ4\n/3369y+/KObae5Z69agp6f37dN3TRHMCYxosIyPDKTAwMKxp06YHfHx8YhMTEysnJiZWDgwMDHv2\nvf3791+s6SCJHTvoLmvjxgIfsmYNHbJnj4xyaYyDB+kzh4ZaMIiHhxCvvy6ZTJokPZ1uObt2lWzI\n7GCCOnWEyMqSbFjV2LaNzqWfflJgsosXyQpp3dq6AiYWLKBFmjdPbUmEEELs3UvLWKeOEDdumD8O\nJLCgLDrY5MnUVlBTptBHvn27wIfcvEmHTJkio1waY9IkIXQ6+uxmM2MGLdzZs5LJpTnWrqXPGBEh\n6bDLl9OwmzZJOqziGAxCNGsmRNWqQqSlKTTpjz/S4i1dqtCEFnL6NGkDX19N3ZHs2CFE0aJC1Khh\nvqdeCgWlo3GUQafTCSXne4727YF//6We7ibw2mvAiy/aTx3UN9+kQCKLqrlfvEi74tOmmdDwx8po\n147yVS5elLQUjV4P1KhB9WZjYiQbVnHWrSOv26+/Au+9p9CkBgOFpJ06Rb/z8uUVmtgM0tOpms35\n8xQqrGBgREHYu5e89ampFKhqasCUTqeDsDCq234qSRgMFMGXTxO5vPDzo0PT02WQS2PcuEHlnUyO\n3nuWatWo4Jqt7kMlJgJRUUD//pIqJ4D2/YYPB2JjKTTbGsnKouTPV18tUOcR6XBwABYupLYvI0Yo\nOLEZjB5Ncd2LF2tOOQHU1eHAAaBWLeqmM3268oG59qOg/v6bik6ZECCRja8vxVfs2yeDXBrDaHNC\nU+jWjRbt4kUJBtMYS5bQTc+AAbIM/8EHlOJjDRWnc2PlSuDkSQqQkDTvqSDUq0flKlat0m6X540b\nqWLERx9JVPdJHqpUoWo63bpR6/g2bagEmmJY6iM05QE196CyNyJPnzb50ORk2pOZNEkGuTRG165C\nVK4s0R7zuXO05lOnSjCYhsjKEuLll2nfQEZGjhTC0VGIM2dknUZy0tNpeZo0UXFbJT1diPr1hXjx\nRSGuXlVJiDy4dEkIFxeKhlFsc84ysrJoe69UKSGcnIQYN06Ihw/zPwZKhJnbDLt2AeXKkb1qImXK\nUF6vNe8HFIS0NPJadewoUaTryy+Tj33FCgkG0xBxcbRvIFGPnrz47DNy9339tazTSM7ChbQ8kycb\n7XovH87OZEHdu0cluzMzVRLkGdLSqB1LejrJV7iw2hIVCAcHSiQ/fZqWc8oUct9+9RVw7pyME1uq\n4Ux5QE0LqkYNIbp0MfvwUaOEcHa2rQTKZ4mIIINnyxYJB50zhwY9elTCQVWmVy/KYlbgZAgOFsLB\nQYh//pF9Kkm4d49K43l5aSTSe8kSOv/GjVNbElqQwECS57ff1JbGIuLjhfDzI88SIISnpxA//ECR\npydOUG1JcJh5Abl6lT7q9OlmD7F5Mw2xfbuEcmmMQYMo4lVSr8OtW1RIdfRoCQdVkeRkIQoXFmLo\nUEWmu3GDaof26qXIdBYzZozQXt7gwIFCjnQAk/n8c5tzeV+8SCk4derQR3v6wQqqYPz2m8W/mnv3\naD9ACzdicpCVJUSlSkJ06ybD4O3aUf0UDeV5mM3MmULpau2jR9Od6okTik1pFmfOkJehb1+1JXmG\n1FRKGndxkaaGjzmEhtJ5M2CARkxLaTEYyA7YvZuqok+ezAqq4IwcKUSRIrRxagHu7vSwRQ4cEJZX\nj8iLFSto8B07ZBhcQQwGIV55heo5KsitW2TZvvuuotOaTMeOQhQvLsSVK2pLkgunT9MOf926Fmag\nm0FMDEUW+PpafA2yJqRQUPYRJJGQQNmnzs4WDePrC+zfT/uutsbGjbQR2q6dDIN36kQ941eulGFw\nBYmJoV3ioUMVnbZsWUrpWbMGOHZM0akLzNatdA5NmABUqqS2NLlQuzYJeP48NZVUqs95VBTlbNSu\nTRX+LbwG2Ru2r6AePqRkODPyn56ldWtKQLTFaL6NGymHuWxZGQYvXpxS0n/7jfpJWCtz51I7kXfe\nUXzq4GBqWDxunOJTG0WvB0aOBGrWBD7+WG1p8qFFCypvceIE0LYtJfPKybp1lOH66qsU+VmmjLzz\n2SC2r6D27yetYkYFiWdxdwdKlgQiIyWQS0NcvAj89ReFl8tG795AcrL1Ll5iIvUsHzgQKFJE8enL\nlCHrZNMmav6nJebMoTz4GTOsIGr67beB1aupREKHDlTHRw5WrKAbmTfeoJIg5crJM4+tY6mP0JQH\n1NiD+uor2mFOTpZkuE6dhKhe3bb2ObPra5qRw1xwMjKo8rfWN1Ly4osv6Dz691/VREhPpy2UGjWE\nePRINTGe4vx52h9r08bKfhMrVtD32aCBtDH8ev1/1ZZ9faklvZ0C3oMqADt3Ag0aSGZev/02tZtR\ntNyHzGzYANSpQ25y2XByogy/DRuApCQZJ5IBvZ66nbZtS1VcVcLZmTqBnzsHfPutamI8IWelpwUL\nVGtjZB69epE1f/061YxcvdryMc+eJU/NpEnkMYiIoHpVjNnYtoLKzAR27wa8vSUbsk0berZWT9Wz\npKSQe1xW9142H3xAGfRLlyowmYSsX08XMoWDI3LDzw8IDASmTpU5g78AzJtH3qsZM6g2sNXRujVV\n4339dVrUIUOoW6SpZGXRDUyjRlRFfdUqYNkyoGhR6WW2Nyw1wUx5QGkX3/795LtavVrSYV95RYi2\nbSUdUjVWrqQl2rVLoQnd3IR49VXr8gf5+JBfNzNTbUmEEEIkJpJbrUMH9WQ4e5YSiK3OtZcbGRlU\nKgagJkjDhpHv0hjJyUJ8+60Q1arRsX5+Qly+LLe0VgM4D8oI06fTR5S4WOTw4XQe20LZo65dhahY\nUcEc2sWL6TuJi1NoQgv580+Sd9o0tSV5iuxTe9065efOyqJSRqVL29j1+MQJId57j3KWHB2FCAig\nzPxlyyhR8NAhutmdPJlKexQrRl9CixZChIfbRiK6hLCCMkbHjkLUqiX5sNllj6KiJB9aUe7fp/zl\nYcMUnPThQ7qy9eyp4KQW0Ls3mSt37qgtyVNkZAjRqBEVR7h0Sdm5v/2Wzv8lS5SdVzESE4X49FOy\n9B0dRS41fOiubsAAIY4cUVtazSKFgrLdjroGA4V2du4MLFok6dCpqYCLC/Dhh9bbrwegtKR336V9\nBB8fBScePpx21a9ckSnxSiIuXqTknhEjNPlFnzlDUcwNGgA7dlAcitxs2UKNqbt0ofPHqgIjzCEj\n478u3AYDRRLVrMnBDwWAO+rmx8mTlHcjYYBENsWK0bDWHigRHk463MtL4YkHDaIf/pIlCk9sIjNn\n0hV45Ei1JcmV2rWptcWePcDYsfLP9/ffFEvQsCF9dTavnAAKnaxbl3rXd+9OARWsnBTDdhVUfDw9\ny3T1fftt0oGXL8syvOw8ekRJn127St6x3Dj161Nlj59/Vr6HdEFJTgZ++YVC46tUUVuaPHn33f8s\n+Q0b5JsnKYnyWosUoaojfI1mlMB2FdTOnUDlytQ0Twayw82jomQZXnaioqgKVPfuKgkwaBD5qOLi\nVBLACPPm0QJ9+qnakhjl++/J1de/P0U5S41eT4rw8mWKuK9aVfo5GCY3bFNBCUEWlLe3bH6IevWA\nl14in7w1snYt7aO1aKGSAN27U/L0nDkqCZAPaWnA7NlkJjdsqLY0RilcmArJOjvT93n8uHRjp6bS\nVxUTQwalu7t0YzOMMWxTQZ07B1y9Ksv+UzY6HRUW2LaNck+tifR04I8/KH5EiY31XClalBIjf/9d\nntt+S1i6FLh5k3quWwk1alCghKMjKanDhy0fMykJaNmSzpWffgL69rV8TIYxBdtUUDt30rPMu/9d\nugD37wPR0bJOIznbt1PLENXce9mMGEGbGtOmqSxIDvR6qiP0xhsKhzZaTp065DgoWZIqTuzZY/5Y\nFy7QNuHhwxRM8+GHkonJMAXGNhVUfDy1RahbV9Zp/PzoYvD777JOIzlr1wKlSwP+/ioLUr488P77\nVBbm0iWVhXnMokUUVvzll1YZplazJt2flS9P3++UKeSxLCgGA7XtcnMDbtygm5kuXeSTl2Hyw3YV\nlJcXdeCTkcKFKSdkwwYqx2UNZGSQvB07aqR32qhR9KyFPKPUVFJMnp4ydW5UhipV6CfQti0wfjwF\nTa5fbzxgMiYGaNaM6pxWrEh9PiXoUsMwZmN7CurqVboDlnH/KSddulB9yV27FJnOYiIjgTt3gB49\n1De6xYgAABKVSURBVJbkMVWrAn360A78zZvqyjJ7NhWFDQmxSuspJxUqkGtu+3ba7uvShZpKf/wx\nGawnT1LfvrAwYMwY+rn4+9O5vGwZcOgQBQIxjJrYXiWJsDAqpX/wIO0jyMyDB1QMYfBgyuvUOu+8\nQ5vpV66oGCDxLKdOkTt27Fhg8mR1ZLhzhyINvLwo0ceGyMwE5s8n192RI5QDlxMnJ1r+oCDgo49U\n6cfI2CBSVJKwPQX1wQe0yXL7tmIZqB07UkfaCxe0feN99y7dWf/vf2QsaIp336XkrEuXaINMaUaP\nBr77jr5IKwgtN5fMTLofOHyYztXXX6eO5Jpw9zI2BZc6yo2YGIq+UrA8QpcudF2VIrRXTtaupRDz\nPn3UliQXxo6l0MIfflB+7itXSGMHBdm0cgKAQoVoT6pPn/8+LisnRqvYloK6cIFyoPz8FJ22QwfS\nh+vWKTqtySxbRnfLTZuqLUkuNG5MVtS0acD588rOPX48Rbl8+aWy8zIMky+2paBiY+lZYQVVtixt\nMms53PzCBYrsCgrSsBvy++9J0ytZnHXrViA0lKIJZSqLxTCMediWgoqJoQQQFcKPunShas9aK4qQ\nzcqV9BwUpK4c+fLSS8Dnn1OQwubN8s93/z7tWdapA3zxhfzzMQxjErajoISgkg5+fqqYCJ0707MW\nrSghyL3n5QVUr662NEb4+GPyQw4fblqGqTmMHUsVUH/9lUPXGEaD2I6COnUKuHZNcfdeNlWqUPb9\nypXa6yBx6BDwzz8at56ycXamwm///gtMny7fPPHxVKh2+HCugMowGsV2FFRMDD2rpKAAoF8/qiSt\ntWi+5cvpuv/OO2pLUkBatqRCgZMnU0sOqUlNBQYOpD0ntfKuGIYxim0pqKpVKdlSJQIDqfyRlhrF\npqUBK1ZQpGGZMmpLYwIzZlDr4k6dKPxcKoSgXlRnz1I72uLFpRubYRhJsQ0FZTBQBJ9K+0/ZvPAC\nBUusXKmdFhxr1lDO8pAhaktiIlWqAL/9Bpw+TZVBpCp2OGECmZRff62qtc0wjHGMKii9Xu8UFBS0\n3M3Nba+Hh0fCqVOnXs35enp6euEePXqsfuutt/a5u7vv2bZtWyv5xM2Do0epRbcGLjj9+5Mof/yh\ntiTEnDkUpKaBpTEdPz9KoI2IAMaNs3y8+fOpvPf//ke5TwzDaBshRL6PhQsXDhw5cuQPQgjEx8d7\nBQQEbMr5+uLFi/sPHTp0jhACt27dKlu7du3TeY1F08nA9OlCAEIkJsozvglkZgpRubIQ7dqpLYkQ\n+/fTssyerbYkFjJ4MH2QpUvNH2PjRiEcHIQICBBCr5dONoZhcuXx9d6ojsnvYdSCio6O9u/atevv\nAODp6bnryJEjjXK+Xr169QuDBw+eDwBFihRJe/DgQQkZ9Gj+xMRQaHLlyopP/SyOjtR5NDKSggrV\nZM4coEQJCt6wambPpvJVAwcCs2aZFiYpBG3C9egBNGkCrF5N9X4YhtE8Rn+pSUlJrq6urkkAFXvV\n6XRPXR18fHziAOD48eMN/ve///08atSofGODJ02alPNY+FjatVSvp5BhDRWY698fmDqVtjo+/VQd\nGW7fBlatAgYMAEqVUkcGyXByogSzfv2oykR0NOUulS2b/3HJybT5tmYN0Lw51aLioAiGkYW4uDjE\nxcVJO6gxEyswMDAsISGhuRACBoNBV6VKlUvPvufLL7/84rXXXvsrJibGN7+xIIeLb8cOcv+sXSv9\n2BbQvLkQdesKYTCoM39ICC3L8ePqzC8LBoMQs2YJ4exMftQNG4RITX3+fY8eCbF+vRCVKglRqJAQ\nU6aQ75VhGMWABC4+oxaUv79/dHh4eLfmzZvvjoqKauPt7R2f8/WwsLCeBw8ebHrgwIFmzs7OGdKq\nzwKwZQu5bFq2VHzq/Ojfn/bi9+8H3npL2bmzsoB588grVr++snPLik5HibVeXuSy69SJErzefBNo\n0YKs6V27qBdYRgaVvNq0iQrRMgxjdRjtB6XX65369u279OzZs7VKlCjxYPny5UEAMGrUqOlhYWE9\n+/XrF3r48OEmZcuWvQ2QGzAmJibXmDFZ+kE1akTx3VKblhaSkkJbYp07k6tPSTZupGv32rVAt27K\nzq0Yjx5RakFcHHVgPHQIcHCgUu2enoCHB9CmDZcwYhiV4IaFV65QgdGQEGo4pzGCg2lP/8wZ5Qpl\nC0HGxLlzVMHcbuIBHjygCJWiRdWWhGEYcMNCCpUDgLZt1ZUjDz75hG7qv/tOuTm3bQN27qQ6qHaj\nnAAKV2TlxDA2hXVbUN27A3v3UkVqjTY5+uADqiR+4QK1W5cTIWi/68YNKsBQuLC88zEMw+SFfVtQ\nej2ZC23balY5AcBnn5GoSnQy/+MP4MABam3EyolhGGvHei2o+HjabAkPB7p2lWZMmQgMpP57ly5R\nPIccGAwUrJaaSo0T7cq9xzCM5rBvC2rzZk2Gl+fGmDHUvHXOHPnmWLuWShJ++SUrJ4ZhbAPrtaBe\nfx1wcaFQYyugXTtyv128SF0kpCQzE2jQgBTTX39RMBvDMIya2K8FdeUKmQsajd7LjfHjqfzQV19J\nP/bixdRQ+KuvWDkxDGM7WKeC0nh4eW54eFBdvOnTKadUKs6fp3wrb2/qRcUwDGMrWKeLr1s3YN8+\nTYeX58bdu1R9p1w5cvc5O1s2XlYWlTP66y/g2DGgWjVJxGQYhrEY+3TxpaUBW7fSpo4VKSeAIvjm\nzyfvZEiI5eNNn06l5376iZUTwzC2h/VZUBs2UIG7qCigdWtpBFOYnj0pOv7wYQpuMIcjR6hGaqdO\n1E3CynQ1wzA2jn3W4gsKogrm169TnyAr5NYtcvVVr051Tk2N6ktNJeWUnEyuPVdXWcRkGIYxG/tz\n8aWlUanuzp2tVjkBtAe1YAFZUK1b095UQbl1C/DzA06epJ59rJwYhrFVrEtBbdtGGa/vvKO2JBbT\ntSt1H9+/nwpiFKQ9/L//UmPYv/6iBrNvvy2/nAzDMGphXQrqt9+AMmUAf3+1JZGE7t2BiAhSPJ6e\n9JwX+/cD7u7AnTtATAwZkQzDMLaM9exBpacD5ctTiPmvv0ormMrs20dBiffukTXVqRM97t+nArCb\nNgF79lCkXmQk8MorakvMMAyTP/YVJLFpE9ChA5kc7dpJK5gGOHcOWLgQWL+eir3mpEkToH174MMP\nSUczDMNoHftSUP36UYDEjRuWZ7hqnNOnSQ8XLw4EBFDreIZhGGvCfhRURgaZDp07A0uWSC4XwzAM\nIy32E2a+fTuQkkJRBQzDMIxdYB0KavVqoFQpoFUrtSVhGIZhFEL7CurOHQov79GD+5gzDMPYEdpX\nUMuWAY8eAUOGqC0JwzAMoyDaDpIQAqhfHyhZkpKFGIZhGKtAiiCJQlIJIws7d1JSkI0l5jIMwzDG\n0bYF1bMnlU64csX0kt8MwzCMath2mPnNm9Q0qV8/Vk4MwzB2iHYV1OLFgF4PDBqktiQMwzCMCmjT\nxWcwALVqAVWrAnFxssvFMAzDSIvtuvi2bgXOnwcGD1ZbEoZhGEYltGdBCUGNjy5fphLfnJzLMAxj\nddhmmPnatZTztGgRKyeGYRg7RlsWVEYGUK8eULQocOQI4OiomGwMwzCMdNieBTV/PvU937yZlRPD\nMIydox0LKiUFqFkTaNQI2LYN0FmkeBmGYRgVsa0ovpAQICkJ+PZbVk4MwzCMRiyoc+eoKGz37lS9\nnGEYhrFqbKPle1IS4OEB3LhBgRHVqikmD8MwDCMP1h8k8egR0KEDcOEC7TuxcmIYhmEeo94eVFYW\nVSvfuxdYsQLw8lJNFK0Sx2WeTILXyzR4vUyD10t5jCoovV7vFBQUtNzNzW2vh4dHwqlTp1415fVc\nuXePyhht2ADMng1062bBR7Bd+AdhGrxepsHrZRq8XspjVEEtXbq0b7ly5W7t3bvXLSQkZExwcPD3\nprz+HP36ARUqAAsXAqNHA8OGWfgRGIZhGFvEqIKKjo7279q16+8A4OnpuevIkSONTHn9OdavB/r0\nAfbsodByhmEYhskNIUS+j9atW0edOHGiXvb/X3rppcumvJ7zAUDwgx/84Ac/7ONhTL8YexiN4nNx\ncUm+e/fuC6DZdDqdTpjyek4sDTlkGIZh7AejLj5/f//o8PDwbgAQFRXVxtvbO96U1xmGYRjGHIwm\n6ur1eqe+ffsuPXv2bK0SJUo8WL58eRAAjBo1anpYWFjP3F6vXLnyFUWkZxiGYWwWRStJMAzDMExB\n0U6xWIZhGIbJgSQKypxkXrMSfG0Ec5Ofvb294319fWN9fX1jhw0b9pM60itPQc+V1atX9xg7duxU\nU46xRcxZL8B+zy/A+Jqlp6cX7tGjx+q33nprn7u7+55t27a1yszMLGSP55g5awWYeX5ZGgYohMDC\nhQsHjhw58gchBOLj470CAgI2GXt90aJFA/I7xpYf5qzXw4cPi9nTGpmyXgaDQdeyZcttRYoUeTR2\n7NgpBTnGlh/mrNeDBw+K29Mambpmixcv7j906NA5QgjcunWrbK1atc7Y6zXM1LWqXbv2aXOvX5II\n3LNnz5Xx8fFeQtDJX7ly5URjrxs7xpYf5qzX4cOHG9evX/+4n59fdKtWrbYeOHCgqdqfQyvrJYRA\nZmam46+//vremDFjphb0GFt9mLNehw4damKv51dB1iw2Ntbn6NGjDYUQuH//fomKFStetddzzJy1\nMvf8ksTFl5SU5Orq6poEUEuNZ3Ohnn0dAJKTk11cXFyS8zrGljF1vXQ6nXByctJ/9NFHP0ZHR/vP\nmjVrRI8ePVYbDAa72EM0tl4A4OjomJXz7wU5xlYp6Ho5ODgYsv9vz+cXYHzNfHx84ho2bHjs+PHj\nDVq3br111KhR05OSklzt8RpmzlqZe35J0m7D1GReBwcHg4uLS3JKSkrpvI6xZcxJfq5fv/6JBg0a\nHAeAunXr/l22bNnb169fr1CpUqWryn8CZSloMnjOv5uSQG5rmPPZGzRocLxhw4bHAPs7v4CCrdlX\nX331RXh4eLeZM2eO9PX1jT1w4EAze7yGmbNWQgidOeeXJHdIpibzenl57fTz84ux1wRfc9Zr6tSp\nYydNmjQJAK5du1bx3r17pSpWrHhNceFVwJxkcHtOIDfns0+ZMmWcvZ5fgPE1CwsL63nw4MGmBw4c\naObr6xtbkGNsFXPWyuzrlxQ+yYyMDKfAwMCwpk2bHvDx8YlNTEysnJiYWDkwMDAsr9dz+5vavlWl\nHuas171790p26NBho6en505vb+8dO3fu9FT7c2hlvbIfS5Ys6Ze96c/nl2nrZc/nV0HWrG/fvqEN\nGjQ45uPjE+vj4xPr6+sbY6/nmDlrZe75xYm6DMMwjCaxm01QhmEYxrpgBcUwDMNoElZQDMMwjCZh\nBcUwDMNoElZQDMMwjCZhBcUwDMNoElZQDMMwjCZhBcUwDMNoElZQDGMGKSkppRctWjTw2b8/fPiw\n+IABA341ZaxVq1YF/vHHHx2kk45hbANWUAxjBnfu3CmzcOHC95/9+7Rp00b37dt3qSljBQYGrpo9\ne/bwzMxMSYo3M4ytwAqKYcxg7NixU0+ePFnv+++/D87+m8FgcNiwYUMnHx+fOADYvHlzu2XLlvUB\ngPHjx0++dOlS1ZMnT9bz8vLa6evrG9umTZuoa9euVQSAt956a9+WLVvaqvJhGEajsIJiGDMICQkZ\nU69evZPBwcHfZ//tzJkztXNWaI6OjvZv3LjxnwBw+PDhJlWrVr0UHh7erVWrVttiY2N9x4wZE3Lp\n0qWqANC4ceM/4+PjvZX/JAyjXVhBMYwZCCF0z/4tKSnJtWTJkvez/3/s2LGGDRo0OJ6enl7Y2dk5\nAwBGjBgxKzMzs1CPHj1WL1mypH/p0qVTAKB06dIpt27dKqfcJ2AY7cMKimHMILcmbRUrVryW3cAu\nNTW1WGpqajEA2Ldv31uNGjU6Eh8f77148eL3OnfuvH716tU9WrduvXXWrFkjAODevXulsruUMgxD\n8KYsw5hBpUqVriYlJbnOmTPnww8//HAOALz88svnk5KSXAFSSikpKaUjIiICkpOTXdLT0ws7OTnp\nvb2944cOHTrXyclJn5WV5Thv3rwhAPDnn3829vDwSFDzMzGM1uB+UAwjIV999dUXnp6eu3bv3t3c\ny8trZ4sWLXYYO8ZgMDi0bt1665YtW9o6OTnplZCTYawBVlAMIyFpaWlFhg0b9pMQQjd37tyhhQsX\nTjd2zJo1a951cnLSd+nSZZ0SMjKMtcAKimEYhtEkHCTBMAzDaBJWUAzDMIwmYQXFMAzDaBJWUAzD\nMIwmYQXFMAzDaBJWUAzDMIwmYQXFMAzDaBJWUAzDMIwm+T98O4nFeZuu7QAAAABJRU5ErkJggg==\n", "text": [ "<matplotlib.figure.Figure at 0x82dab50>" ] } ], "prompt_number": 48 }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we see the spontaneous decay acts as a damping on the oscillations and the system tends toward a steady state (similar to friction in the pendulum in Example 2).\n", "\n", "QuTip has many more useful features you'll find on its [website][qutip].\n", "\n", "[qutip]:http://www.qutip.org/" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Conclusion\n", "\n", "We've looked at some basic numerical integration methods: **Explicit Euler**, **Runge-Kutta**, introduced multistep methods with the **Adams-Bashforth Two-Step** method and solving stiff problems with the **Implicit Euler** method. We've looked at how to rewrite a higher-order ODE as a system of first order ODEs so that it can be solved. We've looked at the importance of **accuracy**, **stability** and **computational cost** in choosing numerical methods. Finally, we've introduced the **SciPy** ODE integration package and the **QuTiP** solver for the time-evolution of closed and open quantum systems.\n", "\n", "Many more familes of integration methods exist, some general and some useful for particular classes of problems. And there are many ways of optimising and tuning these methods, such as adaptive stepsize and use of the system Jacobian. _Numerical Recipes_ (see References) is a good place to start further reading." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## References\n", "\n", "- Jos Thijssen. _Computational Physics_. (Cambridge University Press, March 2007).\n", "- William H. Press et al. _Numerical Recipes 3rd Edition: The Art of Scientific Computing_ (Cambridge University Press, August 2007)\n", "- Lennart Edsberg. _An Introduction to Computation and Modeling for Differential Equations_. (Wiley-Blackwell, August 2008)\n", "- [SciPy](http://www.scipy.org/): Fundamental library for scientific computing\n", "- [QuTip](http://qutip.org/): Quantum Toolbox in Python" ] }, { "cell_type": "code", "collapsed": false, "input": [], "language": "python", "metadata": {}, "outputs": [] } ], "metadata": {} } ] }