{ "metadata": { "name": "" }, "nbformat": 3, "nbformat_minor": 0, "worksheets": [ { "cells": [ { "cell_type": "heading", "level": 1, "metadata": {}, "source": [ "Calculation of control fields for symplectic dynamics using L-BFGS-B algorithm" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Alexander Pitchford (agp1@aber.ac.uk)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Example to demonstrate using the control library to determine control\n", "pulses using the ctrlpulseoptim.optimize_pulse function.\n", "The (default) L-BFGS-B algorithm is used to optimise the pulse to\n", "minimise the fidelity error, which in this case is given by the\n", "'Trace difference' norm.\n", "\n", "This in a Symplectic quantum system example, with two coupled oscillators\n", "\n", "The user can experiment with the timeslicing, by means of changing the\n", "number of timeslots and/or total time for the evolution.\n", "Different initial (starting) pulse types can be tried.\n", "The initial and final pulses are displayed in a plot\n", "\n", "This example assumes that the example-control-pulseoptim-Hadamard has already been tried, and hence explanations in that notebook are not repeated here." ] }, { "cell_type": "code", "collapsed": false, "input": [ "%matplotlib inline\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "import datetime" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 1 }, { "cell_type": "code", "collapsed": false, "input": [ "from qutip import Qobj, identity, sigmax, sigmay, sigmaz, tensor\n", "from qutip.qip import hadamard_transform\n", "import qutip.logging as logging\n", "logger = logging.get_logger()\n", "#Set this to None or logging.WARN for 'quiet' execution\n", "log_level = logging.INFO\n", "#QuTiP control modules\n", "import qutip.control.pulseoptim as cpo\n", "import qutip.control.symplectic as sympl\n", "\n", "example_name = 'Symplectic'" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 2 }, { "cell_type": "heading", "level": 3, "metadata": {}, "source": [ "Defining the physics" ] }, { "cell_type": "code", "collapsed": false, "input": [ "#Drift\n", "w1 = 1\n", "w2 = 1\n", "g1 = 0.5\n", "A0 = Qobj(np.array([[w1, 0, g1, 0], \n", " [0, w1, 0, g1], \n", " [g1, 0, w2, 0], \n", " [0, g1, 0, w2]]))\n", "\n", "#Control\n", "Ac = Qobj(np.array([[1, 0, 0, 0,], \\\n", " [0, 1, 0, 0], \\\n", " [0, 0, 0, 0], \\\n", " [0, 0, 0, 0]]))\n", "ctrls = [Ac] \n", "n_ctrls = len(ctrls)\n", "\n", "initial = identity(4)\n", "\n", "# Target\n", "a = 1\n", "Ag = np.array([[0, 0, a, 0], \n", " [0, 0, 0, a], \n", " [a, 0, 0, 0], \n", " [0, a, 0, 0]])\n", " \n", "Sg = Qobj(sympl.calc_omega(2).dot(Ag)).expm()" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 3 }, { "cell_type": "heading", "level": 3, "metadata": {}, "source": [ "Defining the time evolution parameters" ] }, { "cell_type": "code", "collapsed": false, "input": [ "# Number of time slots\n", "n_ts = 1000\n", "# Time allowed for the evolution\n", "evo_time = 10" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 4 }, { "cell_type": "heading", "level": 3, "metadata": {}, "source": [ "Set the conditions which will cause the pulse optimisation to terminate" ] }, { "cell_type": "code", "collapsed": false, "input": [ "# Fidelity error target\n", "fid_err_targ = 1e-10\n", "# Maximum iterations for the optisation algorithm\n", "max_iter = 500\n", "# Maximum (elapsed) time allowed in seconds\n", "max_wall_time = 30\n", "# Minimum gradient (sum of gradients squared)\n", "# as this tends to 0 -> local minima has been found\n", "min_grad = 1e-20" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 5 }, { "cell_type": "heading", "level": 3, "metadata": {}, "source": [ "Set the initial pulse type" ] }, { "cell_type": "code", "collapsed": false, "input": [ "# pulse type alternatives: RND|ZERO|LIN|SINE|SQUARE|SAW|TRIANGLE|\n", "p_type = 'ZERO'" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 6 }, { "cell_type": "heading", "level": 3, "metadata": {}, "source": [ "Give an extension for output files" ] }, { "cell_type": "code", "collapsed": false, "input": [ "#Set to None to suppress output files\n", "f_ext = \"{}_n_ts{}_ptype{}.txt\".format(example_name, n_ts, p_type)" ], "language": "python", "metadata": {}, "outputs": [], "prompt_number": 7 }, { "cell_type": "heading", "level": 3, "metadata": {}, "source": [ "Run the optimisation" ] }, { "cell_type": "code", "collapsed": false, "input": [ "# Note that this call uses\n", "# dyn_type='SYMPL'\n", "# This means that matrices that describe the dynamics are assumed to be\n", "# Symplectic, i.e. the propagator can be calculated using \n", "# expm(combined_dynamics.omega*dt)\n", "# This has defaults for:\n", "# prop_type='FRECHET'\n", "# therefore the propagators and their gradients will be calculated using the\n", "# Frechet method, i.e. an exact gradient\n", "# fid_type='TRACEDIFF'\n", "# so that the fidelity error, i.e. distance from the target, is give\n", "# by the trace of the difference between the target and evolved operators \n", "result = cpo.optimize_pulse(A0, ctrls, initial, Sg, n_ts, evo_time, \n", " fid_err_targ=fid_err_targ, min_grad=min_grad, \n", " max_iter=max_iter, max_wall_time=max_wall_time, \n", " dyn_type='SYMPL', \n", " out_file_ext=f_ext, init_pulse_type=p_type, \n", " log_level=log_level, gen_stats=True)\n" ], "language": "python", "metadata": {}, "outputs": [ { "output_type": "stream", "stream": "stderr", "text": [ "INFO:qutip.control.pulseoptim:System configuration:\n", "Drift dynamics generator:\n", "[[ 1.0+0.j 0.0+0.j 0.5+0.j 0.0+0.j]\n", " [ 0.0+0.j 1.0+0.j 0.0+0.j 0.5+0.j]\n", " [ 0.5+0.j 0.0+0.j 1.0+0.j 0.0+0.j]\n", " [ 0.0+0.j 0.5+0.j 0.0+0.j 1.0+0.j]]\n", "Control 1 dynamics generator:\n", "[[ 1.+0.j 0.+0.j 0.+0.j 0.+0.j]\n", " [ 0.+0.j 1.+0.j 0.+0.j 0.+0.j]\n", " [ 0.+0.j 0.+0.j 0.+0.j 0.+0.j]\n", " [ 0.+0.j 0.+0.j 0.+0.j 0.+0.j]]\n", "Initial operator:\n", "[[ 1.+0.j 0.+0.j 0.+0.j 0.+0.j]\n", " [ 0.+0.j 1.+0.j 0.+0.j 0.+0.j]\n", " [ 0.+0.j 0.+0.j 1.+0.j 0.+0.j]\n", " [ 0.+0.j 0.+0.j 0.+0.j 1.+0.j]]\n", "Target operator:\n", "[[ 0.54030231+0.j 0.00000000+0.j 0.00000000+0.j 0.84147098+0.j]\n", " [ 0.00000000+0.j 0.54030231+0.j -0.84147098+0.j 0.00000000+0.j]\n", " [ 0.00000000+0.j 0.84147098+0.j 0.54030231+0.j 0.00000000+0.j]\n", " [-0.84147098+0.j 0.00000000+0.j 0.00000000+0.j 0.54030231+0.j]]\n", "INFO:qutip.control.pulseoptim:Initial amplitudes output to file: ctrl_amps_initial_Symplectic_n_ts1000_ptypeZERO.txt\n", "INFO:qutip.control.optimizer:Optimising pulse using L-BFGS-B\n", "INFO:qutip.control.pulseoptim:Final amplitudes output to file: ctrl_amps_final_Symplectic_n_ts1000_ptypeZERO.txt\n" ] } ], "prompt_number": 8 }, { "cell_type": "heading", "level": 3, "metadata": {}, "source": [ "Report the results" ] }, { "cell_type": "code", "collapsed": false, "input": [ "result.stats.report()\n", "print(\"Final evolution\\n{}\\n\".format(result.evo_full_final))\n", "print(\"********* Summary *****************\")\n", "print(\"Final fidelity error {}\".format(result.fid_err))\n", "print(\"Final gradient normal {}\".format(result.grad_norm_final))\n", "print(\"Terminated due to {}\".format(result.termination_reason))\n", "print(\"Number of iterations {}\".format(result.num_iter))\n", "print(\"Completed in {} HH:MM:SS.US\".format(datetime.timedelta(seconds=result.wall_time)))" ], "language": "python", "metadata": {}, "outputs": [ { "output_type": "stream", "stream": "stdout", "text": [ "\n", "------------------------------------\n", "---- Control optimisation stats ----\n", "**** Timings (HH:MM:SS.US) ****\n", "Total wall time elapsed during optimisation: 0:00:08.261914\n", "Wall time computing Hamiltonians: 0:00:00.208757 (2.53%)\n", "Wall time computing propagators: 0:00:07.326453 (88.68%)\n", "Wall time computing forward propagation: 0:00:00.052961 (0.64%)\n", "Wall time computing onward propagation: 0:00:00.051577 (0.62%)\n", "Wall time computing gradient: 0:00:00.611217 (7.40%)\n", "\n", "**** Iterations and function calls ****\n", "Number of iterations: 13\n", "Number of fidelity function calls: 18\n", "Number of times fidelity is computed: 18\n", "Number of gradient function calls: 17\n", "Number of times gradients are computed: 17\n", "Number of times timeslot evolution is recomputed: 18\n", "\n", "**** Control amplitudes ****\n", "Number of control amplitude updates: 17\n", "Mean number of updates per iteration: 1.3076923076923077\n", "Number of timeslot values changed: 17000\n", "Mean number of timeslot changes per update: 1000.0\n", "Number of amplitude values changed: 17000\n", "Mean number of amplitude changes per update: 1000.0\n", "------------------------------------\n", "Final evolution\n", "Quantum object: dims = [[4], [4]], shape = [4, 4], type = oper, isherm = False\n", "Qobj data =\n", "[[ 5.40302194e-01 -3.23386879e-07 9.59430701e-08 8.41471056e-01]\n", " [ 3.23386880e-07 5.40302194e-01 -8.41471056e-01 9.59430701e-08]\n", " [ 9.59430675e-08 8.41471056e-01 5.40302194e-01 2.00178239e-07]\n", " [ -8.41471056e-01 9.59430665e-08 -2.00178237e-07 5.40302194e-01]]\n", "\n", "********* Summary *****************\n", "Final fidelity error 4.955449758551298e-14\n", "Final gradient normal 2.5564261471023003e-06\n", "Terminated due to Goal achieved\n", "Number of iterations 13\n", "Completed in 0:00:08.261914 HH:MM:SS.US\n" ] } ], "prompt_number": 10 }, { "cell_type": "heading", "level": 3, "metadata": {}, "source": [ "Plot the initial and final amplitudes" ] }, { "cell_type": "code", "collapsed": false, "input": [ "t = result.time[:n_ts]\n", "\n", "fig1 = plt.figure()\n", "ax1 = fig1.add_subplot(2, 1, 1)\n", "ax1.set_title(\"Initial Control amps\")\n", "ax1.set_xlabel(\"Time\")\n", "ax1.set_ylabel(\"Control amplitude\")\n", "for j in range(n_ctrls):\n", " amps = result.initial_amps[:, j]\n", " ax1.plot(t, amps)\n", "\n", "ax2 = fig1.add_subplot(2, 1, 2)\n", "ax2.set_title(\"Optimised Control Amplitudes\")\n", "ax2.set_xlabel(\"Time\")\n", "ax2.set_ylabel(\"Control amplitude\")\n", "for j in range(n_ctrls):\n", " amps = result.final_amps[:, j]\n", " ax2.plot(t, amps)\n", "\n", "plt.show()" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "display_data", "png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEZCAYAAAC99aPhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xm8XPP9x/HXW1D7EtSWSNTSWopQBNVcW8VOKVVrq6EL\nUWItbYLWUmvxa+1aa2oXJCUi106E2LKIIIhEqCURGlnu5/fH59zcMbnL7Gdm7uf5eNxHZs6cOecz\nJ/fO53x3mRkhhBBCoRZJO4AQQgi1LRJJCCGEokQiCSGEUJRIJCGEEIoSiSSEEEJRIpGEEEIoSiSS\nUNMkDZV0WDuv/0PSmTkeq1HSUaWLrnIkDZJ0c9pxhM4pEkmoOpImS9opl33NbHczuzl535GSnsx6\n/Tdm9uccT23JT1txrS/pTkkfS/pc0iuSTpBU1N+RpAZJ7xdzDNqJO4Ryi0QSqlG7X+hpkLQO8Dzw\nLrCxma0A/BTYAli2Aufv0tEu5Y4hhLZEIglVLSllPCXpQkmfSnpbUt+M1xslHSXpe8BVwDaSvpD0\nafL6PyWdkzxeUdKDkj5KjvWApDVzDOUs4CkzO8nMpgOY2UQzO9TMZiTH31vSWEmfSRqZxNQc52RJ\nA5JSzOeSBkv6lqSlgWHAGkncMyWtnlRV3SXpZkkzgCMkrSFpiKRPJL0p6Vc5XsMV2vvcyTU8R9LT\nSQxDJK0s6VZJMySNktQjY/8mScdJeispnf1VkpLX1pX0ePIZP5Y0OMfrG2pYJJJQC7YCJgArAX8F\nrs94zQAzswnAMcCzZrasmXXNfD15rOS9ayU//wOuzDGGnYC72npR0vrAbUB/YGVgKPCApEUz4vgp\nsCuwNrAJcKSZfQn0BaYmcS9nZtOS9+wN3GlmyyfHHgy8B6wOHACcK2mHHGJfJIfPfRBwKLAmsA7w\nbPKersB4YGDW/vvipbHNgX2AXybbzwH+k5TY1gQuzyG+UOMikYRa8K6ZXW8+MdxNwOqSvt3Kfm1V\n7wjAzD41s3vNbLaZzQLOBfrkGMNKwLR2Xj8IeNDMRpjZfOAiYElg24x9LjezD83sM+ABYLMO4n7G\nzIYkj1dJjnWqmc0xs1eA64DDOwo8h89twI1m9o6ZzcRLSBPN7LHks9wJ9Mo67AVm9rmZvQ9cBhyc\nbJ8D9JS0ZhLnMx3FF2pfJJJQCz5sfmBmXyUPl8n3IJKWknR1Us00A3gcWL65WqYDnwBrtPP66nhp\noTlOA97H78qbfZjx+H90/BmmZDxeA/g0KcE0ey/r+K3K8XNPz3g8G/go63l2rJmdA96j5dqcgifG\nUZJel/SLjuILtS8SSagnbTXQN28fAKwPbJVUF/XBv/RySSSPAvu38/pUILMdQUB34IMcjt1a3Nkd\nDqYCXSVlfqGvxTeTTVvy/dy5dHRYK+vxBwBmNt3MjjazNfGqxr9L+k4Oxws1LBJJqCfTgW6SFsvY\nlvmFuQxeEpghqSsL1/tD21+uA4Ftk4blVWFBw/LNkpYD7gD2kLRjcv4B+J18LlU704GVkuO0GkdS\nhfQMcF7SSL8J3i5xSw7Hz/dz55JYT0oa8bvj7UL/BpD0U0ndkn0+x5NSUw7HCzUsEkmodq11BW7r\njnkEMBb4UNJHGfs2738Z3m7xX/xLeViuxzazt4FtgJ7AWEmf443vLwCzzGwi3lh9BfAxsAewl5nN\n6+hzJR0FbgfeTnpVrd7G5z44Of9U4B7gT2b2WCufM1u+nzuXa34/8CIwBniQlg4QPwCek/RFsk9/\nM5vcRlyhTijNha2SbpyXAV2A68zsglb2uRzYDfgK7+UyJtm+At7YuBH+S/5LM3uuUrGH0FlJagLW\nTZJrCOmVSJIBVlfiXR83BA6WtEHWPrvjv7DrAUcD/8h4+W/AUDPbAO9KOb4igYcQQviGNKu2tgIm\nmdlkM5uL95HfJ2ufvYF/AZjZ88AKklaVtDywvZndkLw2r3lQWAih7Kpq1oGQvjQTyZp8swvhFBbu\nytjaPt3wAV0fS7pR0kuSrpW0VFmjDSEAYGZdolorZEozkeR6V5Pdg8SARfERtX83s82BL4HTShhb\nCCGEHC3a8S5l8wHez75ZdxbuE5+9T7dkm4ApZvZCsv0uWkkkkqIIHkIIBTCznCcCTbNEMhpYT1JP\nSYvjU0wMydpnCMkUEJJ6A58nA54+BN5P5jcC2Bnv9rkQM4sfMwYOHJh6DNXyE9cirkVci/Z/8pVa\nicTM5kk6FngY7/57vZmNl3RM8vrVZjZU0u6SJuHVV5nTLRwH3JokobeyXgshhFAhaVZtYWbD8MFR\nmduuznp+bBvvfQXYsnzRhRBCyEWMbO8kGhoa0g6hasS1aBHXokVci8KlOrK93CRZPX++EEIoB0lY\njTS2hxBCqAMdJhJJ35U0QtLY5Pkmks4sf2ghhBBqQS4lkmuBP+ArnwG8RstqaCGEEDq5XBLJUubz\nXAELVn6bW76QQggh1JJcEsnHktZtfiLpANpfuzqEEEIn0mGvLUnrANcA2wKfAe8Ah1gNLFYTvbZC\nCCF/Je+1ZWZvmdlOwMrAd81su1IlEUl9JU2Q9KakU9vY5/Lk9Vck9cp6rYukMZIeKEU8IYQQ8tfm\nyHZJAzKeWsZ232B2STEnzljYamd8IsYXJA0xs/EZ+yxY2ErS1vjCVr0zDnM8MA5YtphYQgghFK69\nEsmywDLAFsBv8LVBugG/xqdwL1bBC1sBSOoG7I4vt5tzESyEEEJptVkiMbNBAJKeBDY3sy+S5wOB\noSU4d2uLVm2dwz5rAtOBS4GTgeVKEEsIIYQC5dJr69t8s7vv3GRbsQpd2EqS9gQ+MrMxrbweQgih\ngnKZ/fcmYJSke/Av7X1JqpuKVMzCVvsDeydtKEsAy0m6ycwOzz7JoEGDFjxuaGiIidlCCCFLY2Mj\njY2NBb8/p0kbJW0BbI+XIp5ISgJFkbQo8AawEzAVGAUc3Epj+7FmtnuysNVlZtY76zh9gJPMbK9W\nzhHdf0MIIU/5dv/tsEQiaS3gY+DeZJNJWsvM3iswRj9I8QtbfeNwxcQSQgihcLkMSHydli/qJYC1\ngTfMbKMyx1a0KJGEEEL+Sl4iMbONs06wOfC7AmILIYRQh/Jej8TMXmLhbrohhBA6qVzaSDJHuC+C\nD0b8oGwRhRBCqCm5dP9dlpY2knnAg8DdZYsohBBCTcklkYwzszsyN0j6KXBneUIKIYRQS3LptTXG\nzLJn3V1oWzWKXlshhJC/kvXakrQbPinimpIup2UqkmWJFRJDCCEk2qvamgq8iM/I+yItiWQmcEKZ\n4wohhFAjcqnaWiyZ5r30J5f6ApfhI9uvM7MLWtnncmA34CvgSDMbI6k7PgfYt/GOANeY2eWtvDeq\ntkIIIU+lrNq608x+CrzUvJhVBjOzTQqMsfn4xSxsNRc4wcxelrQM8KKk4ZnvDSGEUBntVW0dn/y7\n0GSIJbJgYSsASc0LW2Umg28sbCVpBUmrmtmHwIfJ9lmSxgNrZL03hBBCBbS3sNXU5N/JZTp3oQtb\ndcMXtgJAUk+gF/B8OYIMIYTQvvaqtmbR9qy6ZmbFrkxY6MJWmevHLwPcBRxvZrNae3OsRxJCCO2r\nyHok5ZCsLzLIzPomz08HmjIb3CVdBTSa2eDk+QSgj5lNl7QYPsp+mJld1sY5orE9hBDylG9je06T\nNkraXNLxko5LZv8thdHAepJ6SlocOAgYkrXPEODwJIbewOdJEhFwPT7qvtUkEkIIoTI6TCSS/oQ3\neHcFVgFulPTHYk9sZvOA5oWtxgH/bl7YKmNxq6HA28nCVlcDv03evh1wKLCDpDHJT99iYwohhJC/\nXMaRTAQ2MbPZyfMlgVfMbP0KxFeUqNoKIYT8laNq6wNgyYznS+C9p0IIIYScZv+dCYyV9EjyfBdg\nlKQr8N5b/csWXQghhKqXS9XWke28bGb2r5JGVEJRtRVCCPnLt2orte6/lRCJJIQQ8lfyNhJJeyW9\noj6T9EXyM7O4MEMIIdSLXKq23gL2A143s6aKRFUiUSIJIYT8laPX1hRgbK0lkRBCCJWRS6+tU4Fh\nkkYCc5JtZmaXlC+sEEIItSKXEsk5wCx8/Mgyyc+ypTi5pL6SJkh6U9KpbexzefL6K5J65fPeEEII\n5ZdLiWR1M9ul1CcuZmGrXN4bQgihMnIpkQyVtGsZzr1gYatkKd/mha0yfWNhK2AFSavl+N4QQggV\nkEsi+S3eRjK7xN1/W1u0as0c91kjh/eGEEKogA6rtsxsmTKdu9CFrfIiDcp41pD8hBBCaNGY/BQm\nlzYSJK0IrIc3uANgZk8UfFb3AdA943l3Fp4MMnufbsk+i+Xw3iTOQUWGGUII9a6BzJts6ay83p3L\nyPZ+wBPAI8BZ+Pohg/I6S+sKXtgqx/eGEEKogFzaSI7HG7cnm9kOQC9gRrEnLmZhq7beW2xMIYQQ\n8pfLFCmjzewHkl4GepvZbEnjzGzDyoRYuJgiJYQQ8pfvFCm5tJG8n7SR3AcMl/QZMLnA+EIIIdSZ\nvKaRl9QALAf8x8zmdLB76qJEEkII+Yv1SDJEIgkhhPyVY/bfEEIIoU2RSEIIIRQlEkkIIYSitNlr\nS9Is2p7GxMxsufKEFEIIoZa0mUjKOMdWCCGEOpJT1ZakTSUdJ+lYSZuW4sSSukoaLmmipEckrdDG\nfq0uYCXpQknjkwWv7pG0fCniCiGEkJ9c5to6HrgVWAVYFbhFUv8SnPs0YLiZrQ+MSJ5nn7t5Aau+\nwIbAwZI2SF5+BNjIzDYFJgKnlyCmEEIIecplipTX8KlRvkyeLw08Z2bfL+rE0gSgj5lNTxarajSz\n72Xtsw0w0Mz6Js9PAzCz87P22w/Y38wOzdoe40hCCCFP5RpH0tTG42KsmszkCzAdL+1ky2XxK4Bf\nAkNLFFcIIYQ85DLX1o3A85LuwReZ2he4IZeDSxoOrNbKS2dkPjEzk9Ra0aHD4oSkM4A5ZnZba68P\nGjRoweOGhgYaGho6OmQIIXQqjY2NNDY2Fvz+dqu2JC0CbAPMBn6If7E/aWZjCj5jy7EnAA1m9qGk\n1YGRrVRt9QYGZVRtnQ40mdkFyfMjgX7ATmY2u5VzRNVWCCHkqeRzbUl62cw2KzqyhY/7V+ATM7sg\naftYwcxOy9pnUeANYCdgKjAKODhZt6QvcDHezvLfNs4RiSSEEPJUjjaSRyUdIKmotdNbcT6wi6SJ\nwI7JcyStIekh6HABqyuAZfCp7cdI+nuJ4wshhJCDXEoks4ClgPl4FRfUyMj2KJGEEEL+Sr6wVYxw\nDyGE0J5cBiSOyGVbCCGEzqm9SRuXxKu0VpHUNeOl5Wh9LEcIIYROqL2qrWOA44E1gBcztn+BT1sS\nQggh5NTY3t/MLq9QPCUVje0hhJC/sqzZLmlboCcZJRgzu6mQACspEkkIIeSv5L22JN0CfAd4Ge8C\n3KzqE0kIIYTyy2WurS2ADePWPoQQQmtyGdn+OrB6KU9a7KJWGa8PkNSU1asshBBCBeWSSFYBxiVf\n+A8kP0OKPG+xi1ohqTuwC/BukbGEEEIoQi5VW4OSf5urtkQO07t3YG+gT/L4X0AjCyeTrYBJZjYZ\nQNJgYB+gea6tS4BTgPuLjCWEEEIROiyRmFkjMAEfiLgsMM7MHi/yvEUtaiVpH2CKmb1aZBwhhBCK\nlEuvrQOBC4Hm5HGlpJPN7M4O3leWRa2SEfd/wKu1FmxuK45Y2CqEENpX1oWtACS9CuxsZh8lz1cB\nRpjZJgWftIhFrYCH8HaVr5JduwEfAFs1x5hxjOhsFkIIeSrHeiQCPs54/gntlAByNAQ4Inl8BHBf\nK/uMBtaT1FPS4sBBwBAze93MVjWztc1sbbzKa/PsJBJCCKEycmls/w/wsKTb8ARyEDCsyPOeD9wh\n6ShgMnAg+KJWwLVmtoeZzZPUvKhVF+D6jEWtMkWRI4QQUpTrFCn7A9slT580s3vLGlWJRNVWCCHk\nr2RzbUlaD+9d9VTW9h8C08zsraIirYBIJCGEkL9StpFcBsxsZfvM5LUQQgih3USyamvjNJJta5cv\npBBCCLWkvUTS6vxXiSVKHUgIIYTa1F4iGS3p6OyNkvrxzRUTQwghdGLtNbavBtwLzKElcWwBfAvY\nz8ymVSTCIkRjewgh5K+kKyRKErADsDE+XmOsmT1WdJQVEokkhBDyV5aldmtVJJIQQshfOaZIKblS\nLGwl6ThJ4yW9LumCykReu4qZkK3exLVoEdeiRVyLwqWSSChyYStJO+BrmmxiZhsDF1Uq8FoVfyQt\n4lq0iGvRIq5F4dJKJHvjC1qR/LtvK/ssWNjKzOYCzQtbAfwGOC/Zjpl93Mr7QwghVEBaiaSoha2A\n9YAfSXpOUqOkH5Qv1BBCCO0pW2N7Bwtb/cvMVszY91Mz65r1/v2BvmbWL3l+KLC1mR0n6TXgMTM7\nXtKWwL/N7DutxBAt7SGEUIB8GttzmUa+0CB2aes1SdMlrZaxsFVra4l8AHTPeN4dL5WQ/HtPcp4X\nJDVJWsnMPsmKodh1U0IIIXQgraqtghe2Sl67D9gRQNL6wOLZSSSEEEJlpDKORFJX4A5gLZKFrczs\n88yFrZL9dsNnGm5e2Oq8ZPtiwA3AZvjI+wFm1ljpzxFCCKHOBySGEEIov7SqtsquvcGMnYmk7pJG\nShqbDN7sn3ZMaZLURdIYSQ+kHUvaJK0g6a5kYO84Sb3Tjiktkk5P/kZek3SbpG+lHVOlSLohabd+\nLWNbToPGm9VlImlvMGMnNBc4wcw2AnoDv0vjWkhaS9IXyfxthbz/C0k9SxDK8cA4wCT9U9I5JThm\nxUk6UtKTRR7mb8BQM9sA2AQYX2AsgyTdnDwu6v+5lWP3TDrTlO27Kvm96gdsbmbfx6vSf1au81Wh\nG/HvykwdDhrPVJeJhPYHM3YqZvahmb2cPJ6Ff1ms0dH7ki+q1yR9KWmapL9LWj7X80qaLGnHjDje\nM7NlC538LHnv5ELemxFTN2B34DpA+ESkbcYjaXVJ10uaKmlmcuc+SNJSRcZR9i/H5DzLSJolaWgr\nry0PbG9mNwCY2Twzm1HgqRZcw+z/52Sc11EFHrdSZuI3XEtJWhRYCu812imY2ZPAZ1mbcxk0vkC9\nJpL2BjN2WsmdVy/g+Q72GwCcDwwAlsNLMj2A4UlHh1wY/mVdTS4FTgaaMra1GmPSIeRZfNmE3ma2\nHLALsDywToniafP6JKXqYu0PvAc0SMoe9Ls28LGkGyW9JOnaIhJke//PVd8Ia2afAhfj12oq8LmZ\nPZpuVKnLZdD4AvWaSKr+l7fSJC0D3AUcn5RM2tpvOWAQcKyZPWJm883sXeBAoCdwaLLfoKR+fXBy\nt/6ipE2S127Ge+Q9kFRznJR9F57cqZ4j6elknyGSVpZ0q6QZkkZJ6pERV5Ok7ySPd0/qs2dKmpIk\nvub99pT0sqTPkmN/v3k7nkCuB4YCm9P+Sp8nAjPM7FAzew/AzKaY2Qlm9lpyzG0lvSDp8yTebTLi\naJR0tqSnkjgflrRS8vITyb+fJ6/1TkqAT0u6RNJ/gYGSlpN0k6SPkhLeGXlWGR2Bl76ebv5/yzAU\n+AHwQ3ymiG2AcyQNS67/8OZ68Yz/u36SPkhKaANoRca+XST9BdgeuDL5P768tdJYZqkled9Fkj6W\n9BawR9bxl88oJU5Jfoeaf6fWlfR48v/xsaTBuVwkSesAv8d/v9cAlpF0SC7v7QyS0mX736lmVnc/\n+B30fzKenw6cmnZcKV6PxYCHgd/nsG9fvJi/SCuv/RO4LXk8CO96/RO8TnkA8DbQJXn9HWDHjPf2\nxL/IF0meNwIT8Tvj5YCxwJv4+KAueHH6hoz3NwHfSR5PA7ZLHi8P9Eoe98LvnrbE75IPT+JYDC9h\nzQM+Sd4/G5gPnN3GdXgOGNjOdeqKVwccgt+Q/Qz4FFgx4/O9CayLJ6yR+Pxw4KW7Bdci2XZkct1/\nlxxvCeAmfHG5pZP3vAH8MmP/J9uJr0fyebvh9f+vZL3+XnINVsG/PD8FZgCb4qWwEcCfsv7vbgWW\nxNcn+gjYKeN34eY2/p9HNsfc2uvZ+wC/xqtf1wRWTF6bn3G8e4F/JHGsgpeuj05eux04PXm8OLBt\njn8fBwHXZTw/DPi/tP9uK/mT/L+8lvF8ArBa8nh1YEJ776/XEkl7gxk7leQO9npgnJldlsNbVgb+\na2ZNrbz2YfJ6s9Fmdo+ZzQcuwb/8cu35Y8CNZvaOmc0EhgETzeyx5Hh34omhNXOAjSQtZ2YzzGxM\nsv1o4Goze8HcTcDX+N32UGC6ma2Ef+kPx5NFW7riCactewBvmNmtZtZkZoPxP769sz7fJDObjY+b\n2ix5ra1SxVQz+7/k2s/Ff29PN7MvzUuFF+Nfcrk4DBhlZs2zQGwoabOM1+fjSXZFM5uK/9++ZWav\nmNnX+Bd29vU/y8z+Z2av4w20B+cYSz6lqAOBS83sAzP7DDi3+f1J9dxueOeR/5lP1noZLQ3jc4Ce\nktY0szlm9kyO55wA9Ja0ZPL3sjPeIaMzy2XQ+AJ1mUjMbB5wLH4XPg6fi6ugHil1YDu8WmMHebfX\nMZKye2hk+i+wslpvCF4dyJxpuXnKGsxvXaaQQ0N+hukZj2fzzalyZgPLtPG+/fFG88lJtUhz8uoB\nDEiqtT6T9Bl+R756EldmA6oB79L2l9wnHXyWNfC7+kzvZr3nw4zH/2vn8zTLbNdbGS9JvZux7T1y\nb+s7HE/GmM/60EjLF0Ozi4FbJb2SxDYs47XWrn9mfO+R+/91PlXNq7dynmY98GsyLeP/9yq8ZAJw\nCv7/OUre1f0XOQVn9gpe+hsNvJpsviaPmGuapNuBZ4DvSno/uW7nA7tImojXEpzf3jHKNtdW2sxs\nGN/8w+iUzOwp8rtheBa/i9+f5IsIFrSx9MWrCZt1z3h9EfxLe2rzqfMNNecdzUYD+8obpI+jZZaE\n94C/mNm52e+R1IfkS9jMHgcel/Q0MKmN0zwK7CfprCRJZvsAr9bL1IPcfufa+qyZ2/+Ll0p60tIt\ndy0ykndbJG2LV6mdKemUZPOywCaSBmSUNt82sy2T99yM39G3Zy28eq35cS49m7I/65fJv0sBzW11\nmZO7TkuOnXnOZu/jv5srtVZiNm8cPhpA0nbAo5IeN7O3OwzS7K/AXzvarx6ZWVsly51zPUZdlkhC\n4cy7gJ4FXCFpV0mLyXt73YH/Id+csfsWkvaTd5n8PX4X21xdNJ2Oezepjcdtv8HjOUTS8kkV2Bd4\nNQ3AtcCvJW0lt7SkPZIk+AwwT1L/5Bg/wdtS2nIJ3nbzL0lrJedeU9LF8gb8ocD6kg6WtKikg4Dv\nAQ/m8Jk+xtsJ2rw+yWe7A/iLvBtvD+AE4Jb2rk/iCOARYAO8zWNTvF1jSbwkV6gzk+qfjfA2mn/n\n8J5v/B4k1VEfAIclDeu/5JvX4Q6gf3KtVyRj/IKZTUs+1yWSlpW0iKR1JP0IQNJP5V28AT7Hk1hr\nVbShxCKRhIWY2YXAH/CVJ2fgyeFdvHF1bvNuwP14Pf6neKPzT5IvQIDz8C+ezySdmPGeb5wq63FH\nrzc7FHhH0gz8DvSQJO4X8YblK5OY3sSreEji/gn+BfgJXhd/dzvX4DNgW7xU8LykmXgp5XN8jNKn\nwJ54J4P/AicBeybb2/18ZvYV8BfgaUmfStq6jc9/HH4H/zbwJN7YfWP28TJJWgL4KXCFmX2U8TMZ\nvwk4vK3P3Fa8GR7HS3CPAhdaSxfZ7H0zH/8NOCD5nM1tdP3wbtj/xQcMP52x/7V4lfQreFXT3VnH\nOxxvSB+H/x/fSUuJ5gfAc5K+wH83+1uRY49Cbqpyri1JN+CNmR+ZjzRtbZ/L8Ya3r4AjMxpcQwVI\nGgisa2a5Nv6GGpWUSN8GFm2jE0bo5Kq1RNLakP0FJO2Of4mth9+R/qNSgYUFqm2wYQghJVWZSKz1\nIfuZFgzfN7PngRW08MjdUF4dD1IK9ST+r0ObarXXVmtToHTjm91JQxmZ2VlpxxAqI2lnKMWULaFO\n1WoigYWrVlpreIy7qBBCKIDlsVR5VVZt5SB7PfdutNGnvVJTDFT7z8CBA1OPoVp+4lrEtYhr0f5P\nvmo1kQwh6caYjGr+3FpmqgwhhFBBVVm1lQzZ74NP1fE+MBCfGgEzu9rMhspngJ2E97PPaSqEEEII\npVeVicTaHrKfuc+xlYilXjQ0NKQdQtWIa9EirkWLuBaFq8oBiaUiyer584UQQjlIwjpBY3sIIYQq\nEYkkhBBCUYpKJJK+K2mEpLHJ800knVma0EIIIdSCotpIJD2Bz+J5lZn1SlYXe93MNipVgMWINpLi\n/e9/8NRT8PzzMHo0vPsuTJsGX38NZtC1K6y6Kmy4IfTqBT/8IWy6KeS1sngIrZg/H15+GRob4ZVX\n4M03YcoU+OormDcPllsOVlkF1l8fvv99+NGPYOutYfHF04689uXbRlJsIhltZj+QNMbMeiXbXjaz\nzTp6byVEIinM/Pnw0ENw663w8MP+R7rNNvCDH8B3vgNrrAFLLOH7fvaZJ5bXX4cxY2DECE8y++4L\nv/qVJ5UQcmXmNy633AJ33QWrrw4NDbD55rDeerDWWrDMMtClC8ycCdOnwxtvtCSct97y371DD4Ud\ndoBFovK+IJVOJMPwNRPuTEokBwBHmdluBR+0hCKR5Gf2bLj2WrjsMlhpJU8E++3nd325MoOJE+Hf\n/4brrvOkc9ppsM8+UUoJbZs7139n/vpXL20ccQQcfLAnjnxMmwa33w7//Kf/Lp56Khx0ECy2WFnC\nrluVTiTr4Gsbb4vP1vsOcIhVyWIykUhy09QEgwfDGWfAxhvD6ad7CaTYL/7582HIEDj7bD/W+efD\nj39cmphDfTCD+++Hk0+Gbt38i3/XXYv/3TODRx6B887zUsull0LfNhemCNkqmkgyTro0sIiZfVH0\nwUooEknHJk3ykseXX8JFF0GfPqU/hxncey+cdBJssYX/UXfr1vH7Qn0bNw5+9zv4+GO45JLy3GSY\neTXtCSfABhvAVVd5KTm0ryLjSCQNyPg5ETgG6CfpxIxlVUMVM4MrroDevb3a6bnnypNEwO8uf/IT\nGDvWG+UgkpbJAAAcwElEQVQ328yr0CLHd07z58OFF3rj+P77e/tGuUqqEuy5p7fh9erlP4MHl+dc\nnVlBJRJJg/Bp278LbIlPoih8DetRZnZoCWMsWJRIWjdzJhx1FLz9ttdLr7tuZc8/fjz8/OfQo4e3\no6y8cmXPH9Lz3nve9rH44nDDDbD22pU9/wsvwGGHeRK74gr41rcqe/5aUZESiZkNMl/YqDuwuZkN\nMLMTgS2AHoUcM5ukvpImSHpT0qmtvN4gaYakMclPjF/Jwfjx3vtq5ZXh6acrn0TAqxiee8574Wy2\nmXctDvXvkUdgq628BDxiROWTCMCWW3oy+eQTL4F/0OriEyFvRc5Z/wawRMbzJYA3SjAXfhdgEtAT\nn/X3ZWCDrH0agCEdHMdCi8ZGs29/2+yGG9KOpMV995mtvLLZjTemHUkol6Ymsz//2Wz11c1Gjkw7\nGtfUZHbeeWZrrmk2Zkza0VSf5Lsz5+/sYmf/vQkYJekevGprX5K11Iu0FTDJkt5fkgYD+wDjs/aL\nDqU5GjwY+vf3rpE77ZR2NC322cdLJvvu63XlF1/sYwRCfZgzx6tRJ070Aa3V0tAtebf0ddf19pk7\n7vDxKqEwRQ3XMbO/4GuBfA58ChxpZueWIK7W1mRfM/v0wLaSXpE0VNKGJThvXbr2Wu8xNWJEdSWR\nZhtuCKNGwauvws9+5uNZQu37/HPvcjtrFowcWT1JJNMBB3g74YEHes/CUJhi59paC/gYuBe4D/gk\n2VasXFrIXwK6m9mmwBXJ+UOWa66Bc87xP+Tvfz/taNq2wgowbJiPRO7bF2bMSDuiUIypU326nE02\n8RHqSy2VdkRt22EHn8HhN7+B++JbpCDFVm0NpeVLfwlgbbzdpNi5trLXZO+Ol0oWsIwxK2Y2TNLf\nJXU1s08z9xs0aNCCxw0NDZ1q8ZqrroJzz4XHHkunUT1f3/qWV739/vfeEDp8eH6j6kN1ePddL/ke\ndZQPbq0FvXrB0KGw226w6KLeZbgzaWxspLGxsfAD5NOg0tEPsDlwfQmOsyjwFt7YvjitN7avSkv3\n5a2Aya0cpxTtTjXp5pvNunc3mzQp7Ujy19RkdsYZZhtvbDZ9etrRhHxMmmTWo4fZZZelHUlhRo0y\nW2UVs2HD0o4kXVS4sT07Kb0kaesSHGeepGOBh/EeXNeb2XhJxySvXw0cAPxG0jzgK+BnxZ63Xgwb\nBgMGeHXWOuukHU3+JK+OW3RRr3Z47DGfYThUtwkTYJdd4Mwz4Zhj0o6mMFtu6VO27LOPj4jfcsu0\nI6oNxc61NSDj6SJ4iaSrme1abGCl0BkHJD7/vBfLhwzx+bJq3dlne3XXY4/5TLChOk2a5NWR557r\nEy7WugcegKOPhiefrI1q4VKr9FK7ywLLJD+LAw/i3XRDCiZN8jupG2+sjyQC8Kc/wSGHwM47+yCy\nUH3ef9//fwYOrI8kArDXXjBokLeZfPRR2tFUv2JLJAea2R1Z235qZncWHVkJdKYSyYwZnjyOO857\nn9SbU0/1qroRI2DZZdOOJjSbPt2nGznmGDixDmfZ++Mf/Xdu5MjONZ1KpaeRX7CgVXvb0tJZEsn8\n+X4Htfba8H//l3Y05WEGv/61r5I3dGjLwlohPZ9+6oP49t/fSyP1qKnJx5gsv7zPC9dZ1tSpSCKR\ntBuwO3AQMJiWEebLAhua2VZ5H7QMOksiOekkHxU+bFh9L+Azf75Xc82e7WMTFi1pV5GQj1mzvIvv\n9tv7TL71/AU7axZst513Z+7fP+1oKqNSiWRToBdwNvBHWhLJTGCkmX2W90HLoDMkkptv9gbp55/3\n9dPr3Zw5Pp3Kyiv7KnixlGrlzZvnbXGrrgrXX1/fSaTZO+941fGtt1bn7BClVumqrcXMbG7BByiz\nek8kr70GO+7oa1VvVOwQ0Bry1VfeuNunj6+AFyqnuYrx3Xe9Z1M9l4CzNTb6sr2jRvkSCPUs30RS\nUOWApDvN7KfAS1r4dsTMbJNCjhtyN3OmzxN0ySWdK4mAT7cxZIhXN6y1Vn12LqhW553n07A//njn\nSiLg7UEnneTJ5IknfE2V4Aqt2lrDzKZK6tna6xZrtpeVmf8yd+3q06B0Vm+/7fM5XXUV7L132tHU\nv5tv9l5Mzz7becf0NDV51eq66/pNXL1KZc32alWvieTyy+Ff//KFqTp776UXXoA99vBqlq2LnlMh\ntOXRR72jw8iRPltzZ/bpp7D55nDppbDffmlHUx6VamyfRdsz9JqZLZf3QcugHhPJ6NGw++6+wuB3\nvpN2NNXhoYfgV7/qvKOQy+3VV71N6q67fMxI8M4te+1Vv3+HUSLJUG+JZNYsvxP685+9b3tocc01\n3g31mWdixuBSmjLFeytdeKGvFRNa/O1v3ovr6afrr72o4olE0ubA9kAT8LSZvVTUAUuo3hJJv34w\nd653ew0LO/NMr4J57LHqXv+iVsyY4W1QRxzhjczhm8y8WnWLLXyS0XpS0bm2JP0JX1q3K7AKcKOk\nPxZzzIxj95U0QdKbkk5tY5/Lk9dfkVQVo+nL5Z57vH76iivSjqR6nXMOrL8+/PznPngxFG7OHK//\nb2jwmaTDwiS44QYf8f7UU2lHk65ix5FMBDYxs9nJ8yWBV8xs/aKCkrrgC2TtjC9y9QJwsJmNz9hn\nd+BYM9s9mbr+b2bWO+s4dVEi+eADr9K6/37o3bvj/TuzOXO8Del73/Ok2xkGy5WaGRx2GHz5pbeL\ndOmSdkTV7YEHfMT7yy/7VCr1oNKz/34ALJnxfAmyVjIs0FbAJDObnAx4HMzCswrvjZeGMLPngRUk\n1d2qFU1NcPjhcOyxkURysfjicPfd3s//oovSjqY2nXEGvPWW1/9HEunYXnvBrrv6hKmdVbGJZCYw\nVtI/Jf0TeB2YIekKSZcXcdw1gfcznk9JtnW0T7cizlmVLrnE77L/8Ie0I6kdyy/vEztecQUMHpx2\nNLXlqqvgzjt9wGe0M+Xu4ou9J1dn/X0rdtq7e5OfZo0Zj4upU8r1vdlFr4XeV8trto8bB+ef711+\n484wP926ebfgnXbywXN9+qQdUfV74AE46yyv74+eb/lZemm47TZfv2T77WHN7NveKlfsmu1V2f1X\nUm9gkJn1TZ6fDjSZ2QUZ+1wFNJrZ4OT5BKCPmU3P2Kdm20jmzfMpQH7xC5/bKBRmxAhvfI+BdO0b\nNcp7ID30EGxVFXN316azz/axJQ89VNvtc5XutbWXpDGSPpP0RfIzs5hjJkYD60nqKWlxfLr6IVn7\nDAEOT+LoDXyemURq3aWXwjLL+HKfoXA77eRtJbvvDlOnph1NdcpcWTOSSHFOPx0+/NCvZWdSbK+t\nt4D9gNfNrKlkUbFgzZPLgC7A9WZ2nqRjAMzs6mSfK4G+wJfAL7LHsNRqiWTCBO+//8ILvlhVKN65\n53rd/xNPxAqLmT7+GLbd1seJHHNM2tHUh+ZZuV980ScVrUWVnkb+cWBHM6vKXvu1mEjmz/ckcthh\n8Nvfph1N/Wie/nzyZHjwwfobiVyIr76CHXaAXXbx2RJC6fzlLz5D8sMP12YVV6UTSW98cauRwJxk\ns5lZVcyLWYuJ5OKL/YtuxIhYtKnU5s3zmVu//e3OsyBTW+bN8yVyV1jBZ0rozNeiHObN8+76Rx9d\nm9XTlU4kw4EvgNfwKVIAMLOzCj5oCdVaInnjDW9gHzWqPieCqwazZvlo7d1394bRzqipyZeNnTrV\ne2rFuhrlMXas9xYcPRp69kw7mvxUOpG8bmYbF3yAMqulRDJ/vs+s+rOfde6BTZUwfbpf61//Gk44\nIe1oKssMTjzRb1YeecS7rYbyueACv87Dh9dWDUOlR7YPlbRrkccI+BojXbrA736XdiT1b9VV/Q/7\nb3/zeZI6k3PO8UktH3wwkkglDBjgU83U+wJ0xZZIZgFL4e0jzWu3x3okeXrzTZ+q+7nnYj2NSnrz\nTa/muuQSX3Gy3l1+uY/2f+opT6ahMmqxF2asR5KhFhJJU5PXox5wABx/fNrRdD6vv+6LNl13Hey5\nZ9rRlM+NN8LAgb74V48eaUfT+Vx4IQwb5ssc1EIVV6WrtpC0oqStJP2o+afYY3YmV17p/0a7SDo2\n3tgbnH/5S5+fqx5dfz386U9enRdJJB0nnujdra++Ou1IyqPYqq1+QH+gOzAG6A08a2Y7lia84lR7\niWTSJO8i+OyzsN56aUfTuT33HOy9N1x7rY/yrhfXXONjREaMiN+xtI0f7508Xnih+ntxVbpEcjw+\n5ftkM9sB6AXMKPKYnUJzF8wzzog/8GrQu7dXPRxzjI+ArwdXXeVJ5LHH4nesGmywAZx8sv/dV/H9\nbUGKTSSzzex/AJKWMLMJwHeLD6v+/f3vvmxu//5pRxKabbGFj0Tu3x9uuSXtaApnBn/9q88cPXJk\ndOCoJiee6GOZ6q2Kq9hp5N+XtCJwHzBc0mfA5KKjqnNvvw2DBsHTT8f08NVm0029QXS33WDaNJ+D\nqpZGfTc1eZfTRx/1369am8683i26qHd86NMH+vat/iquXJWs15akBmA54D9mNqeD3ds7Tlfg30AP\nPCkdaGaft7LfZHxhrfnAXDNbaN7SamwjaWryGWn32MO/pEJ1mjLFk0lDA1x2WW0k/K+/hiOP9KWZ\n778fVlwx7YhCWy64wDs/DB9enTcqFe+11czMGs1sSDFJJHEaMDxZ931E8rzVUwINZtartSRSra6+\nGmbP7nwjqmtNt24+3mLsWO+aPWtW2hG176OPfPLFr7/26rlIItVtwACYOdM7d9SDauzRvGAt9uTf\nfdvZtwpzedsmT/ZumDfcUBt3uJ3d8svDf/4DXbt6Y/zEiWlH1LoXX4Qtt/TS0113wZJLph1R6Ehz\nFdcZZ8C776YdTfGqMZGsmrFA1XSgrTG4BjwqaXTSDbmqNffSOukk770RasPii/tgxeOP99HJ99+f\ndkQtzPympG9fXwjt7LNrY7BbcBtt5I3vv/pV7ffiKraxvSDJrMGrtfLSGZlPzMwktXWJtzOzaZJW\nwRv6J5jZk9k7Vcua7Vdf7XPuDBiQyulDESTo188b4g880KuOLrww3bmqPv3Uuyq/8QY0NvqXUqg9\nJ58M99zjNyv9UrwdTmXN9mSOrbbeWNRcW8na6w1m9qGk1YGRZva9Dt4zEJhlZhdnba+KxvZ33vGq\nh6eegu+1+0lCtZsxw7sHP/20r+Pxwx9WPoYhQ+DYY73t5txzYYklKh9DKJ3XX/cFxqppRcWan2tL\n0l+BT8zsAkmnASuY2WlZ+ywFdDGzLyQtDTwCnGVmj2Ttl3oiae6ltfvufvcR6sO99/qX+S67+HiN\n1VorX5fYu+96EpswAf7xD1/ONdSHc8/1kmW1rKiYxlxbm0o6TtKxkjYt9njA+cAukiYCOybPkbSG\npIeSfVYDnpT0MvA88GB2EqkW//iH99I68cS0IwmltN9+/oW+6qo+X9c558DnC3VSL42PPvLfn169\nvGT76quRROrNKad4deX116cdSWGKnWvreKAfcA/eg2pf4Fozu7w04RUn7RLJ22/D1lt7ldZ3Y7x/\n3XrzTZ+K5KGHfFnVfv1KM134uHE+A8Ltt8PPfw5/+AOsvnrxxw3VqZqquCq9QuJrQG8z+zJ5vjTw\nnJl9v+CDllCaiaSpye8a99orGtg7i7fe8vU+br0VNtsMfvITH3ia65eCmd983H8/3H23P+7XzxvV\nY4R65/CXv/hU/8OGpVvFlUYi2Spjvq0lgVGRSPwLZfBgeOKJGDPS2cye7VPTDxni41CWWcarpTbe\n2NtSVlnFu+nOmweffOJtH2+8Ac8/7+/fYw/Yf39vW4v11DuXuXN9zNJvf+vDBdJS6URyInAk36za\n+qeZXVrwQUsorUTSPD38M8/A+utX/PShijQ1edXXmDE+jfj06fDxx/5aly4+Ar1HD1hnHa8G7dGj\nOhpbQ3pee81rM156Cbp3TyeGiiUSSYsA2wCzgR/i3YGfNLMxBR2wDNJIJHPnepfQQw6JmX1DCIX5\n85+9i/nQoencWFS6RPKymW1W8AHKLI1EMnCgV1GkXccZQqhdc+d6CfXYY331zkqrdCK5CHgOuDv1\nARutqHQiefZZ7xY6Zkz0rgkhFKe5iuu557zqs5IqnUhmAUvhU7nPTjYXNbK9lCqZSL74wnvqXHSR\nJ5MQQijW5Zd7L8CnnoLFFqvceWt+ZHspVTKR/PKX3hPnuusqcroQQidgBnvuCZtsAuedV7nzVnRk\nu6QRuWyrd3ff7d18L7ss7UhCCPVE8unmb7oJHnss7WjaVtDsv8l4kaWAVZIVDZstB3SqoVNvvw2/\n+Q08+KCPFwghhFL69rc9mRxxhLe/rrxy2hEtrNDZf38PHA+sAUzNeOkL4Bozu7I04RWn3FVbX3/t\nXX0PPdTXqwghhHI5+WRfXO2++8rfI7TSje39Sz2vlqSfAoOA7wFbmtlLbezXF7gM6AJcZ2YXtLJP\nWRNJ//6+tvfdd0dX3xBCec2ZA9ttBwcfXP5JYCve2C5pW6AnGdVkZnZTEcf7HtAEXA0MaC2RSOoC\nvAHsDHwAvAAcbGbjs/YrWyK5+26/Q3jpJVhhhbKcIoQQvuHdd318yR13wI9+VL7z5JtIilohUdIt\nwHeAl/EuwM0KTiRmNiE5dnu7bQVMMrPJyb6DgX2A8e29qVQmTfJ2kYceiiQSQqicHj284f3gg+GF\nF2CNNdKOyBW71O4WwIYpDEZcE3g/4/kUYOtKnHjmTNh7b18fe8stK3HGEEJo8eMfw69/7cs+jxxZ\n2fElbSl2YavXgbzHcEsaLum1Vn72yvEQqQx+aWqCww6DPn38PzKEENJwxhleG3LSSWlH4ootkawC\njJM0Cvg62WZmtnd7bzKzXYo87wdA5ryY3fFSyUIGDRq04HFDQwMNDQ0Fn3TgQPjsM7jzzoIPEUII\nRVtkEbj5ZthmG7jqquJvbBsbG2lsbCz4/cX22mpIHjYfRHgiebzgg7YceyRwkpm92Mpri+KN7Tvh\n3Y9HUebG9jvu8OUwR43yft0hhJC2SZN8CMJNN3mVV6lUdGS7mTUCE/CBiMsC44pNIpL2k/Q+0Bt4\nSNKwZPuCNdvNbB5wLPAwMA74d3YSKaUnnvBZOO+/P5JICKF6rLuu15AceiiMHZteHMWWSA4ELgSa\nk8ePgJPNrCoqf0pRIhk71mfgvO02X7EuhBCqzS23wB//6GuYlKInV6UHJL4K7GxmHyXPVwFGmNkm\nBR+0hIpNJFOmwLbbwvnnw89/XsLAQgihxM47z2cKfvxxWGml4o5V0aotvE3k44znnyTbat60abDz\nzj71SSSREEK1O/10nym4b18fplBJxZZILgQ2BW7DE8hBwKtmdkppwitOoSWS6dOhocGXyz3zzNLH\nFUII5WAGv/udV8kPHQpLL13YcdKYImV/YLvk6ZNmdm9RByyhQhLJ9OneJnLggd7dN4QQaklTE/Tr\nBxMmFD77RkUSiaT1gFXN7Kms7T8EppnZW3kftAzyTSSTJsGuu8KRR3rDVQgh1KKmJjjhBO9x+vDD\n+fc2rVQbyWVAa7VwM5PXas6LL/okaKecEkkkhFDbFlnEF9rbc08fZzJxYpnPV+D7VjWzV7M3JtvW\nLi6kyrvzTm+g+vvf4Zhj0o4mhBCKJ8E55/g0KttvD8OHl+9chSaS9mrdlijwmBU3b56XQE45xYt/\n++6bdkQhhFBaRx/tM3McdhhcdJFXe5VaoYlktKSjszdK6gcsNKVJNXrjDS/yvfoqjB4Nm2+edkQh\nhFAeffrAc8/Bvfd67cu0aaU9fqGJ5PfALyQ9LumS5Odx4Kjktao1Zw5cfLGvNHbYYd5FrtjBOyGE\nUO169vTBittuC5tuCldfDfPnd/i2nBTc/Ve+8tQOwMb4pI1jzeyx0oRVGpm9tsxgyBBf1XC99bwh\nar31Ug4whBBS8MorPt5kzhwfEb/jjt9cLrzi40hKLY812yfjvcTmA3PNbKtW9rGvvjL+/W8vhXTp\n4tOd9O1bvvhDCKEWNDXB7bd7g3zXrjBggPfy+ta3Kj9FSjm8BuwHPNHBfgY0mFmv1pJIszXW8It1\n8cUwZkznTSLFrDVQb+JatIhr0aKzXYtFFvHZO8aO9amgrrwSunWD3/62gGOVPrzimNkEM8u113OH\nGXP8eO+R9eMff7Po1tl0tj+S9sS1aBHXokVnvRZdusBBB/myvS+8AGutlf8xqi6R5MGARyWNTnqL\ntWq11SoYUQgh1LCePeG00/J/X7FL7RZE0nCgta/4P5jZAzkeZjszm5ZMXT9c0gQze7J0UYYQQshF\n1TW2N0uW2h3QVmN71r4DgVlmdnHW9ur8cCGEUOXyaWxPpUSSh1Y/iKSlgC5m9oWkpYEfA2dl75fP\nhQghhFCYqmsjyWXNdrxa7ElJLwPPAw+a2SPpRBxCCJ1b1VZthRBCqA1VVyIpFUl9JU2Q9KakU9OO\nJy2SuksaKWmspNcl9U87pjRJ6iJpjKRcO3XULUkrSLpL0nhJ4yT1TjumtEg6PfkbeU3SbZK+lXZM\nlSLpBknTJb2Wsa2rpOGSJkp6RFK7y2PVZSKR1AW4EugLbAgcLGmDdKNKzVzgBDPbCK8u/F0nvhYA\nxwPj8O7jnd3fgKFmtgGwCTA+5XhSIakn0A/Y3My+D3QBfpZmTBV2I/5dmek0YLiZrQ+MSJ63qS4T\nCbAVMMnMJpvZXGAwsE/KMaXCzD40s5eTx7PwL4s10o0qHZK6AbsD15HDYNZ6Jml5YHszuwHAzOaZ\n2YyUw0rLTPyGaylJiwJLAR+kG1LlJMMmPsvavDfwr+Txv4B2F9mo10SyJvB+xvMpybZOLbnz6oV3\nUOiMLgVOBsqwIkPNWRv4WNKNkl6SdG3SG7LTMbNPgYuB94CpwOdm9mi6UaVuVTObnjyeDqza3s71\nmkii2iKLpGWAu4Djk5JJpyJpT+AjMxtDJy+NJBYFNgf+bmabA1/SQfVFvZK0Dr78RU+8tL6MpENS\nDaqKJFOot/udWq+J5AOge8bz7nippFOStBhwN3CLmd2Xdjwp2RbYW9I7wO3AjpJuSjmmNE0BppjZ\nC8nzu/DE0hn9AHjGzD4xs3nAPfjvS2c2XdJqAJJWBz5qb+d6TSSjgfUk9ZS0OHAQMCTlmFKRrBtz\nPTDOzC5LO560mNkfzKy7ma2NN6Q+ZmaHpx1XWszsQ+B9Sesnm3YGxqYYUpomAL0lLZn8veyMd8jo\nzIYARySPjwDavQGt9pHtBTGzeZKOBR7Ge2Bcb2adskcKsB1wKPCqpDHJttPN7D8pxlQNovoTjgNu\nTW623gJ+kXI8qTCzV5LS6Wi8/ewl4Jp0o6ocSbcDfYCVk8HgfwLOB+6QdBQwGTiw3WPEgMQQQgjF\nqNeqrRBCCBUSiSSEEEJRIpGEEEIoSiSSEEIIRYlEEkIIoSiRSEIIIRQlEkkIJSJppWSK+jGSpkma\nkjz+QtKVaccXQrnEOJIQykDSQOALM7sk7VhCKLcokYRQPgKQ1NC8kJakQZL+JekJSZMl/UTSRZJe\nlTQsmcYcSVtIapQ0WtJ/muc9CqEaRSIJofLWBnbA13y4BV9AaBPgf8AeySSbVwD7m9kP8IWH/pJW\nsCF0pC7n2gqhihkwzMzmS3odWMTMHk5eew2fynx9YCPgUZ9DkC74OhkhVKVIJCFU3hwAM2uSNDdj\nexP+NylgrJl19qnMQ42Iqq0QKiuXRbXeAFaR1Bt8PRlJG5Y3rBAKF4kkhPKxjH9bewwLT2dvZjYX\nOAC4QNLLwBhgm3IGGkIxovtvCCGEokSJJIQQQlEikYQQQihKJJIQQghFiUQSQgihKJFIQgghFCUS\nSQghhKJEIgkhhFCUSCQhhBCK8v/tKsgSfEKC/wAAAABJRU5ErkJggg==\n", "text": [ "" ] } ], "prompt_number": 11 }, { "cell_type": "heading", "level": 3, "metadata": {}, "source": [ "Versions" ] }, { "cell_type": "code", "collapsed": false, "input": [ "from qutip.ipynbtools import version_table\n", "\n", "version_table()" ], "language": "python", "metadata": {}, "outputs": [ { "html": [ "
SoftwareVersion
OSposix [linux]
IPython2.3.1
matplotlib1.4.2
Cython0.21.2
QuTiP3.1.0
Python3.4.0 (default, Apr 11 2014, 13:05:11) \n", "[GCC 4.8.2]
SciPy0.14.1
Numpy1.9.1
Tue Jan 13 13:35:23 2015 JST
" ], "metadata": {}, "output_type": "pyout", "prompt_number": 12, "text": [ "" ] } ], "prompt_number": 12 } ], "metadata": {} } ] }