{ "cells": [ { "cell_type": "markdown", "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", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [], "source": [ "%matplotlib inline\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "import datetime" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": false }, "outputs": [], "source": [ "from qutip import Qobj, identity, sigmax, sigmay, sigmaz, tensor\n", "from qutip.qip import hadamard_transform\n", "import qutip.logging_utils 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'" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Defining the physics" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [], "source": [ "#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()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Defining the time evolution parameters" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# Number of time slots\n", "n_ts = 1000\n", "# Time allowed for the evolution\n", "evo_time = 10" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Set the conditions which will cause the pulse optimisation to terminate" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# 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" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Set the initial pulse type" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": false }, "outputs": [], "source": [ "# pulse type alternatives: RND|ZERO|LIN|SINE|SQUARE|SAW|TRIANGLE|\n", "p_type = 'ZERO'" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Give an extension for output files" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "collapsed": false }, "outputs": [], "source": [ "#Set to None to suppress output files\n", "f_ext = \"{}_n_ts{}_ptype{}.txt\".format(example_name, n_ts, p_type)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Run the optimisation" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "collapsed": false }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "INFO:qutip.control.dynamics:Setting memory optimisations for level 0\n", "INFO:qutip.control.dynamics:Internal operator data type choosen to be \n", "INFO:qutip.control.dynamics:phased dynamics generator caching True\n", "INFO:qutip.control.dynamics:propagator gradient caching True\n", "INFO:qutip.control.dynamics:eigenvector adjoint caching True\n", "INFO:qutip.control.dynamics:use sparse eigen decomp False\n", "INFO:qutip.control.pulseoptim:System configuration:\n", "Drift dynamics generator:\n", "Quantum object: dims = [[4], [4]], shape = (4, 4), type = oper, isherm = True\n", "Qobj data =\n", "[[ 1. 0. 0.5 0. ]\n", " [ 0. 1. 0. 0.5]\n", " [ 0.5 0. 1. 0. ]\n", " [ 0. 0.5 0. 1. ]]\n", "Control 1 dynamics generator:\n", "Quantum object: dims = [[4], [4]], shape = (4, 4), type = oper, isherm = True\n", "Qobj data =\n", "[[ 1. 0. 0. 0.]\n", " [ 0. 1. 0. 0.]\n", " [ 0. 0. 0. 0.]\n", " [ 0. 0. 0. 0.]]\n", "Initial state / operator:\n", "Quantum object: dims = [[4], [4]], shape = (4, 4), type = oper, isherm = True\n", "Qobj data =\n", "[[ 1. 0. 0. 0.]\n", " [ 0. 1. 0. 0.]\n", " [ 0. 0. 1. 0.]\n", " [ 0. 0. 0. 1.]]\n", "Target state / operator:\n", "Quantum object: dims = [[4], [4]], shape = (4, 4), type = oper, isherm = False\n", "Qobj data =\n", "[[ 0.54030231 0. 0. 0.84147098]\n", " [ 0. 0.54030231 -0.84147098 0. ]\n", " [ 0. 0.84147098 0.54030231 0. ]\n", " [-0.84147098 0. 0. 0.54030231]]\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(s) using GRAPE with 'fmin_l_bfgs_b' method\n", "INFO:qutip.control.pulseoptim:Final amplitudes output to file: ctrl_amps_final_Symplectic_n_ts1000_ptypeZERO.txt\n" ] } ], "source": [ "# 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" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Report the results" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", "------------------------------------\n", "---- Control optimisation stats ----\n", "**** Timings (HH:MM:SS.US) ****\n", "Total wall time elapsed during optimisation: 0:00:02.998303\n", "Wall time computing Hamiltonians: 0:00:00.112098 (3.74%)\n", "Wall time computing propagators: 0:00:02.615654 (87.24%)\n", "Wall time computing forward propagation: 0:00:00.022306 (0.74%)\n", "Wall time computing onward propagation: 0:00:00.019775 (0.66%)\n", "Wall time computing gradient: 0:00:00.222838 (7.43%)\n", "\n", "**** Iterations and function calls ****\n", "Number of iterations: 8\n", "Number of fidelity function calls: 14\n", "Number of times fidelity is computed: 14\n", "Number of gradient function calls: 13\n", "Number of times gradients are computed: 13\n", "Number of times timeslot evolution is recomputed: 14\n", "\n", "**** Control amplitudes ****\n", "Number of control amplitude updates: 13\n", "Mean number of updates per iteration: 1.625\n", "Number of timeslot values changed: 13000\n", "Mean number of timeslot changes per update: 1000.0\n", "Number of amplitude values changed: 13000\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.40298994e-01 -1.56965600e-06 -7.01444121e-06 8.41473111e-01]\n", " [ 1.56965600e-06 5.40298994e-01 -8.41473111e-01 -7.01444122e-06]\n", " [ -7.01444163e-06 8.41473111e-01 5.40298994e-01 1.05774201e-05]\n", " [ -8.41473111e-01 -7.01444163e-06 -1.05774201e-05 5.40298994e-01]]\n", "\n", "********* Summary *****************\n", "Final fidelity error 6.093073931117445e-11\n", "Final gradient normal 8.705152868582485e-06\n", "Terminated due to Goal achieved\n", "Number of iterations 8\n", "Completed in 0:00:02.998303 HH:MM:SS.US\n" ] } ], "source": [ "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)))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Plot the initial and final amplitudes" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XecVNX5x/HPd+kgRQGRKkUUAQUUpdlBA1gwdmNNYgz5\n2TUmdlFUiLFrYokaO0axKyqKgsZCU6QIKB0EBFF6331+f9y7y+xmywBz587uPO/X6752bj3PzMI+\nc8499xyZGc4551ymyYk7AOecc644nqCcc85lJE9QzjnnMpInKOeccxnJE5RzzrmM5AnKOedcRvIE\n5bKOpHclnVfK/kck3ZjktUZLuiB10aWPpCMkLYo7DudK4gnKVQiS5knqk8yxZtbPzJ4Ozztf0n+L\n7B9oZoNTFNfekl6W9JOkVZImS7pSUqWdvK4nF1fheYJyLiKS2gBjgYXAfmZWFzgVOBConYbyK0dd\nhnNR8gTlKpz8WpGkuyT9ImmupH4J+0dLukDSvsAjQA9JayWtDPc/Jem28PWukt6WtDy81tuSmiUZ\nyi3A52Z2pZktATCzmWZ2lpnll3WCpGmSVoZx7ZsQ5zxJfw5rXask/UdSdUm1gHeBJmHcayU1kTRI\n0nBJz0laDZwvqZqk+yQtDpf7JFVL8nO8X9JCSaslTZR0aMK+QWHN8DlJayRNCWuL10paFp53TJHP\nfIikceH13pC0W7ivenidFeHnMF5SoyQ/Y1eBeYJyFVU3YCbQALgTeEKSEg8ws+nAQOALM9vFzOoV\nc50c4N/AnkALYAPwUJIx9AGGl7RT0t7AMOByoCEwAnhLUtWEw04D+gKtgP2B881sHdAPWBzGvYuZ\nLQ6PHxCWWQ94Hrge6A50BjoBBwM3JBn/+PC83YAXgJclVU/YfzzwLLAr8DXwPsHn1RS4FXi0yPXO\nBX4HNAa2Ag+E288D6gLNgfoEv5MNScboKjBPUK6imm9m/zKzXOBpgj+K2/2t3MxWmNkrZrbezNYA\ntwOHJ3l6fWBJKftPB94xsw/MbAtwF1AD6JlwzANmttjMfgbeIkgYpfnCzF43szwz2wCcBdxqZsvM\nbDlBre6cZII3s+fC97/VzO4GqgH7JBzyqZm9b2ZbgZcJkuzQ8L28CLSUlJj0nzWzqWGCvRE4LbwX\nt4Xgs9rLzHLNbKKZrU4mRlexeYJyFdXS/Bdmtj58ucv2XkRSTUmPSpofNpt9AtRLspPDCoLEWJIm\nwPyEOPMI7lc1TThmacLr9ZT9HhaWVkb4ukkZ1wAgbF6cHjYvriSo5TRIOOTHhNcbgJ/CLwT56xSJ\nNzG2+UCV8HrPEtS+XgybIe+UVCWZGF3F5gnKZbuyhvO/iqDW0M3M6gCHhdtV8ikFPgROLmX/YoKm\nw+CCQRNkc+CHJK5dUtxFtxcqg6CZcjFlCO83/YWgiXHXsPlzFcm975I0LxLHFoKktsXMbjGz9gS1\nx+MImgNdlvME5bLdj0CzIvd9EtUmqA2sDG/q37wd174Z6Cnp75L2AJC0V9ghoB7wEnCspN5hjeEq\nYBPweZJx15dUt4zjhgE3SGooqQFwE/BcEtevTXCfaDlQWdJNQJ0kzivN2ZLaS6pJcI9quJnlSjpS\n0n5hrXQ1QeLK28myXAXgCcplu4+AacBSST8Vs/8+gvtCPwFfAu8le2Ezmw30AFoC0yStAl4BJgBr\nzGwmcDbwYHj944HjzWxzEteeQZB85oQ930pqtrstLG8yMAX4KtxWlvcJ3ut3BM1xG/nf5sPt9Szw\nFEGzZXXg0nD7HgQdO1YD04Ex4bEuy8knLHTORU3SaOA5M3s87lhc+eE1KOeccxkpqQQlqYakfco+\n0jnnnEuNMpv4JB1P8HxGVTNrJakzwXMVJ6QjQOecc9kpmRrUIIKnz1cCmNkkgqfanXPOucgkM5jk\nFjNbVXSUmIjiiVSDBg2sZcuWcYfhnHNZbeLEiT+ZWcOyjksmQU2T9BugkqS2BF1Dk3lOo0yS+gL3\nA5WAx81saJH9Cvf3J3iK/nwz+yrcNw9YA+QCW82sa1nltWzZkgkTJqQidOeccztI0vyyj0quie8S\noAPBA4TDCJ5VuHzHQwuED+X9g2DQy/bAmZLaFzmsH9A2XC4EHi6y/0gz65xMcnLOOVe+lFmDCscx\nuz5cUulgYJaZzQGQ9CLBSMzfJhwzAHjGgp4cX0qqJ6lx/tQFzjnnKq4SE5SktyjlXlMKevE1pfCT\n6YsIpkgo65imBCNEG/ChpFzgUTN7rLhCJF1IUPuiRYsWOxmyc865dCmtBnVX+PMkgqFI8sfvOpPC\noxjH5RAz+0HS7sAHkmaY2SdFDwoT12MAXbt2LZedO5xzLhuVmKDMbAyApLuL3ON5S1Iqehr8QOHR\njZvxv6M4l3iMmeX/XCbpNYImw/9JUM4558qnZDpJ1JLUOn9FUiugVgrKHg+0ldQqHEn6DODNIse8\nCZyrQHdglZktkVRLUu0wnlrAMcDUFMTknHMuQyTTzfwKYLSkOQRzwewJ/HFnCzazrZIuJhg1uRLw\npJlNkzQw3P8IwRTY/YFZBN3Mfxue3gh4LXw2qzLwgpklPcq0c865zJfUaOaSqgHtwtUZZrYp0qgi\n0rVrV/PnoJxzLl6SJibzeFCZNShJRWe27CQJM3tmh6NzzjnnypBME99BCa+rA70JJj3zBOWccy4y\nyTyoe0niejhV9YuRReScc86xYxMWrsNHM3fOORexZO5BJY4okUMwbt7LUQblnHPOJXMP6q6E11uB\n+Wa2KKJ4nHPOOSC5Jr7+ZjYmXD4zs0WS/hZ5ZM4557JaMgnq6GK29Ut1IM4551yi0kYz/xPwf0Br\nSZMTdtUGPos6MOecc9mttHtQLwDvAkOAaxK2rzGznyONyjnnXNYrLUGZmc2TdFHRHZJ28yTlnHMu\nSmXVoI4DJhJ0M1fCPgNaF3eSc845lwqlzQd1XPjTH8p1zjmXdqV1kjigtBPN7KvUh+Occ84FSmvi\nu7uUfQYcleJYnHPOuQKlNfEdmc5AnHPOuUTJjMVXneB5qEMIak6fAo+Y2caIY3POOZfFkhmL7xlg\nDfBguP4b4Fng1KiCcs4555JJUB3NrH3C+seSvo0qIOeccw6SG4vvK0nd81ckdQMmRBeSc845l1wN\n6kDgc0kLwvUWwExJUwhGm9g/suicc85lrWQSVN/Io3DOOeeKKDNBmdl8SbsCzROP9wd1nXPORSmZ\nbuaDgfOB2Wyb+t0f1HXOORepZJr4TgPamNnmqINxzjnn8iXTi28qUC/qQJxzzrlEydSghgBfS5oK\nbMrfaGYnRBaVc865rJdMgnoa+BswBciLNhznnHMukEwT33oze8DMPjazMflLKgqX1FfSTEmzJF1T\nzH5JeiDcPzlxCpCyznXOOVe+JVOD+lTSEOBNCjfx7VQ3c0mVgH8ARwOLgPGS3jSzxGGU+gFtw6Ub\n8DDQLclznXPOlWPJJKgu4c/uCdtS0c38YGCWmc0BkPQiMABITDIDgGfMzIAvJdWT1BhomcS5KXXV\nS9/w3tQlUV3eOefKjcP2bsjDZx8YeTnJPKgb1bxQTYGFCeuLCGpJZR3TNMlzAZB0IXAhQIsWLXY4\n2F571WfXmlV2+HznnKso9tp9l7SUk0wNCknHAh2A6vnbzOzWqIJKJTN7DHgMoGvXrlbG4SU66YBm\nnHRA2cc555xLjWRGkngEqAkcCTwOnAKMS0HZPxAMn5SvWbgtmWOqJHGuc865ciyZXnw9zexc4Bcz\nuwXoAeydgrLHA20ltZJUFTiDoCNGojeBc8PefN2BVWa2JMlznXPOlWPJNPFtCH+ul9QEWAE03tmC\nzWyrpIuB94FKwJNmNk3SwHD/I8AIoD8wC1gP/La0c3c2Juecc5kjmQT1tqR6wN+Brwh68P0rFYWb\n2QiCJJS47ZGE1wZclOy5zjnnKo5kevENDl++IultoLqZrYo2LOecc9kuqV58+cxsEwkP6zrnnHNR\nSaaThHPOOZd2nqCcc85lpBKb+BIHZi2OT/nunHMuSqXdg7q7lH0+5btzzrlIlZigIhyDzznnnCtT\nMkMdVQH+BBwWbhoNPGpmWyKMyznnXJZLppv5wwRj3/0zXD8n3HZBVEE555xzySSog8ysU8L6R5K+\niSog55xzDpLrZp4rqU3+iqTWQG50ITnnnHPJ1aCuBj6WNAcQsCfhoK3OOedcVEpNUJJyCEYzbwvs\nE26eGQ555JxzzkWm1ARlZnmS/mFmXYDJaYrJOeecS+oe1ChJJ0tS5NE455xzoWQS1B+Bl4FNklZL\nWiNpdcRxOeecy3LJzAdVOx2BOOecc4nKrEFJGpXMNueccy6VShvNvDpQE2ggaVeCLuYAdYCmaYjN\nOedcFiutie+PwOVAE2Ai2xLUauChiONyzjmX5Uobzfx+4H5Jl5jZg2mMyTnnnEuqk8SDknoCLROP\nN7NnIozLOedclktmuo1ngTbAJLaNwWeAJyjnnHORSWYsvq5AezOzqINxzjnn8iXzoO5UYI+oA3HO\nOecSJVODagB8K2kcUDBIrJmdEFlUzjnnsl4yCWpQ1EE455xzRZXZxGdmY4AZQO1wmR5u22GSdpP0\ngaTvw5+7lnBcX0kzJc2SdE3C9kGSfpA0KVz670w8zjnnMk8yQx2dBowDTgVOA8ZKOmUny70GGGVm\nbYFR4XrRcisB/wD6Ae2BMyW1TzjkXjPrHC4jdjIe55xzGSaZJr7rgYPMbBmApIbAh8DwnSh3AHBE\n+PppYDTw1yLHHAzMMrM5Ybkvhud9uxPlOuecKyeS6cWXk5+cQiuSPK80jcxsSfh6KdComGOaAgsT\n1hdReAzASyRNlvRkSU2Ezjnnyq9kalDvSXofGBaunw68W9ZJkj6k+O7p1yeumJlJ2t5nrB4GBhM8\nMDwYuBv4XQlxXAhcCNCiRYvtLMY551xckhnq6GpJJwGHhJseM7PXkjivT0n7JP0oqbGZLZHUGFhW\nzGE/AM0T1puF2zCzHxOu9S/g7VLieAx4DKBr167+sLFzzpUTJTbVSdpLUi8AM3vVzK40syuB5ZLa\n7GS5bwLnha/PA94o5pjxQFtJrSRVBc4IzyNMavl+TfAwsXPOuQpEJY1gJOlt4Fozm1Jk+37AHWZ2\n/A4XKtUHXgJaAPOB08zsZ0lNgMfNrH94XH/gPqAS8KSZ3R5ufxboTNDENw/4Y8I9rdLKXR6Wt6Ma\nAD/txPkVQbZ/Btn+/sE/g2x//7Dzn8GeZtawrINKS1DjzeygEvZNMbP9diK4cknSBDPrGnccccr2\nzyDb3z/4Z5Dt7x/S9xmU1huvXin7aqQ6EOeccy5RaQlqgqQ/FN0o6QKCGXadc865yJTWi+9y4DVJ\nZ7EtIXUFqhJ0TMhGj8UdQAbI9s8g298/+GeQ7e8f0vQZlHgPquAA6UigY7g6zcw+ijwq55xzWa/M\nBOWcc87FYWeHLHLOOeci4QkqSSVN/ZENJDWX9LGkbyVNk3RZ3DHFQVIlSV+HzwhmHUn1JA2XNEPS\ndEk94o4p3SRdEf4fmCppmKTqcccUtXC802WSpiZsS2rKpJ3lCSoJSUz9UdFtBa4ys/ZAd+CiLHv/\n+S4DpscdRIzuB94zs3ZAJ7Lss5DUFLgU6GpmHQkGEDgj3qjS4imgb5FtZU6ZlAqeoJJTMPWHmW0G\n8qf+yApmtsTMvgpfryH4w9S09LMqFknNgGOBx+OOJQ6S6gKHAU8AmNlmM1sZb1SxqAzUkFQZqAks\njjmeyJnZJ8DPRTYPIJgqifDniVGU7QkqOWVN/ZE1JLUEugBj440k7e4D/gLkxR1ITFoBy4F/h82c\nj0uqFXdQ6WRmPwB3AQuAJcAqMxsZb1SxSWbKpJ3mCcolTdIuwCvA5Wa2Ou540kXSccAyM8vmB9Qr\nAwcAD5tZF2AdETXrZKrwPssAgmTdBKgl6ex4o4qfBV3BI+kO7gkqOSVO/ZEtJFUhSE7Pm9mrcceT\nZr2AEyTNI2jePUrSc/GGlHaLgEVmll9zHk6QsLJJH2CumS03sy3Aq0DPmGOKy4/5s0qUMmXSTvME\nlZwSp/7IBpJEcO9hupndE3c86WZm15pZMzNrSfC7/8jMsuqbs5ktBRZK2ifc1Bv4NsaQ4rAA6C6p\nZvh/ojdZ1lEkQTJTJu20ZGbUzXpmtlXSxcD7bJv6Y1rMYaVTL+AcYIqkSeG268xsRIwxufS7BHg+\n/JI2B/htzPGklZmNlTQc+IqgZ+vXZMGwR5KGAUcADSQtAm4GhgIvSfo94ZRJkZTtI0k455zLRN7E\n55xzLiN5gnLOOZeRPEE555zLSJ6gnHPOZSRPUM455zKSdzN3LkaS6hMMtgmwB5BLMKQQwHozy9YH\nQZ3zbubOZQpJg4C1ZnZX3LE4lwm8ic+5DCVpbfjzCEljJL0haY6koZLOkjRO0hRJbcLjGkp6RdL4\ncOkV7ztwbud4gnKufOgEDAT2JRjVY28zO5hg+o9LwmPuB+41s4OAk8nSqUFcxeH3oJwrH8bnT28g\naTaQP83DFODI8HUfoH0wTBwAdSTtYmZr0xqpcyniCcq58mFTwuu8hPU8tv0/zgG6m9nGdAbmXFS8\nic+5imMk25r7kNQ5xlic22meoJyrOC4FukqaLOlbgntWzpVb3s3cOedcRvIalHPOuYzkCco551xG\n8gTlnHMuI3mCcs45l5E8QTnnnMtInqCcc85lJE9QzjnnMpInKOeccxnJE5RzzrmM5AnKOedcRvIE\n5ZxzLiN5gnLOOZeRPEG5jCaphaS1kirt4PlrJbVOcUxPSbotlddMF0ktJZmk2OeCC6eyX5SwPk3S\nESm8/jxJfVJ1PZd+nqBcSkk6X9IUSeslLZX0sKR623F+oT8qZrbAzHYxs9wdiSc8d86OnLujJDWW\n9ISkJZLWSJoh6RZJtXbyumlLLpJGS/pFUrWoy8pnZh3MbHRY/iBJz6WrbJeZPEG5lJF0FfA34Gqg\nLtAd2BP4QFLVOGNLF0m7AV8ANYAeZlYbOJrg82iThvJ3OnlJagkcChhwws5ez7kd5QnKpYSkOsAt\nwCVm9p6ZbTGzecBpQEvg7PC4QZKGS/pPWLv4SlKncN+zQAvgrbBp7i9Faw3hN/vbJH0eHvOWpPqS\nnpe0WtL48A9sflwmaa/wdX9J34bl/iDpzwnHHSdpkqSV4bX3T9jXJYxzjaT/ANVL+SiuBNYAZ4fv\nHzNbaGaXm9nk8Ho9wzhXhT97JpQ1WtJgSZ+F5Y2U1CDc/Un4c2X43nuENdbPJN0raQUwSFKOpBsk\nzZe0TNIzkupux6/zXOBL4CngvMQdYfPmPyW9G8bwmaQ9JN0X1rhmSOqScPw8SdeGn/svkv4tqdjP\nL7/2LKkvcB1weljGN4n7E44vVMuSdE74nldIur7ItXMkXSNpdrj/pfDLBJKqS3ou3L4y/J002o7P\ny0XEE5RLlZ4Ef7hfTdxoZmuBEQS1iHwDgJeB3YAXgNclVTGzc4AFwPFh09ydJZR1BnAO0JSgVvIF\n8O/wetOBm0s47wngj2GtpiPwEQQJCHgS+CNQH3gUeFNStbDm9zrwbHj9l4GTS/kc+gCvmllecTvD\nP4rvAA+EZd0DvCOpfsJhvwF+C+wOVAXyE+lh4c964efzRbjeDZgDNAJuB84PlyOB1sAuwEOlxFzU\nucDz4fKrYv5YnwbcADQANhF8/l+F68PD95ToLOBXBL+rvcNzS2Rm7wF3AP8J32ensgKW1B54mODf\nRROCz7ZZwiGXACcCh4f7fwH+Ee47j6CG2zw8byCwoawyXfQ8QblUaQD8ZGZbi9m3JNyfb6KZDTez\nLQR/zKoTNAcm699mNtvMVgHvArPN7MOw7JeBLiWctwVoL6mOmf1iZl+F2y8EHjWzsWaWa2ZPE/zh\n7R4uVYD7wlrhcGB8KbHVD99vSY4FvjezZ81sq5kNA2YAxxd5f9+Z2QbgJaBzKdcDWGxmD4bX20CQ\nEO4xsznhF4RrgTOSaf6TdAhBs+xLZjYRmE2QMBO9ZmYTzWwj8Bqw0cyeCe8T/of//fwfCmuRPxMk\n0DPLimMHnAK8bWafmNkm4EYg8UvCQOB6M1sU7h8EnBJ+JlsIfm97hb//iWa2OoIY3XbyBOVS5Seg\nQQl/BBuH+/MtzH8R1jQWEXyrTdaPCa83FLO+SwnnnQz0B+ZLGiOpR7h9T+CqsHlnpaSVBN+mm4TL\nD2ZmCdeZX0psKwjeb0maFHP+fILaYL6lCa/Xl/J+8i0ssl60jPlAZYIaVlnOA0aaWf7v6wWKNPOx\n/Z9/Ynzz2b7fdbKaUPjf1TqC30W+PYHXEn6/04Fcgs/kWeB94EVJiyXdKalKBDG67eQJyqXKFwS1\njpMSN0raBegHjErY3Dxhfw5BU8zicFNiIkgpMxtvZgMIms5eJ6idQPCH7XYzq5ew1AxrN0uAppKU\ncKkWpRTzIfDr8H0VZzHBH8tELYAfknkLSW4vWkYLYCuFE8n/kFSDoPnucAU9MJcCVwCdFN4n3EHN\nE163YNvvujTFvdd1QM2E9T0SXi+h8L+rmgS1onwLgX5FfsfVzeyHsGZ8i5m1J2iqPo6gmdPFzBOU\nS4mwue0W4EFJfSVVCTsrvERQQ3o24fADJZ0U1rYuJ0hsX4b7fiS4b5JSkqpKOktS3bBpcTXbmoD+\nBQyU1E2BWpKOlVSbIPFuBS4N39NJwMGlFHUPUAd4WtKeYdlNJd0TdrwYAewt6TeSKks6HWgPvJ3E\n21gexlzW5zMMuEJSq/ALQv79nOKaXxOdSFCraE/QrNgZ2Bf4lJ37g32RpGbh/bfrCZoBy/Ij0LJI\nop9E0FRZRVJXgma9fMOB4yQdEt43vJXCf98eAW5P+J00lDQgfH2kpP0UPGu3mqDJr9h7iC69PEG5\nlAk7NVwH3EXwH30swTfX3mG7f743gNMJblSfA5wUJg2AIcANYVPMn0mtc4B5klYT3JM4K4x7AvAH\ngo4EvwCzCDoZYGabCWqF5wM/h3G/SgnC+yw9Cf7IjZW0hqD2uAqYZWYrCL6hX0XQBPUX4LiEJrUS\nmdl6gns4n4WfT0n37Z4k+ELwCTAX2EjQSaAs5xHc/1pgZkvzF4LP5axk7mGV4AVgJEFHjtlAMg85\nvxz+XCEp/17hjQQdLX4h+DL0Qv7BZjYNuCjctiQ8puAhYOB+4E1gZPg7+ZKgcwkENbHhBP9mpwNj\nKPyFysVEhZvWnYuWpEEEN6PPjjsWFz1J84ALzOzDuGNx5Y/XoJxzzmWkyBOUpBqS9om6HOeccxVL\npE18ko4nuB9R1cxaSeoM3GpmPnyKc865UkVdgxpE0ONpJYCZTQJaRVymc865CiDqUZG3mNmqwo+Q\nRPecS1kaNGhgLVu2jKt455xzwMSJE38ys4ZlHRd1gpom6TdAJUltgUuBzyMus0QtW7ZkwoQJcRXv\nnHMOkFTaaCwFom7iuwToQPAg5jCC5wwuj7jMjLRxSy65ed6l3zkXn9w8ozw9WhRpDSp8sPD6cMkq\nI6ct5aUJC/l45vJCiUkCM2i2aw06NKnD4AEd2b1OabM3OOfcjlmxdhN3f/Ad4+b+zKxlawu2V62U\ngwRbcvM4fO+GnHRAM47vFMUQiTsnkl58kt6ilHtNcfXi69q1q0XdxPfsF/O48Y1phba1b1yHTs3r\n0ahONWYvX8f4uT+zdPXGQse8fckhdGy6PVP2OOdc8WYtW0Ofez4ptK3ZrjXo2KQuezaoiRCjZy7j\nux/XkNiw89e+7fjTEZHPq4mkiWbWtczjIkpQh4cvTyIYRiR/UrEzgR/N7IqUF5qEKBPUqvVb6HTr\nyELbhv2hOz3a1C/2+Lw847Z3pvPkZ3MLtlWvksPkm39F1cr+/LRzbvttzc2jx9CPWL5m28hi53Tf\nk0EndKBSjoo956sFv3DSPwt3DRh/fR8a1q4WWZyxJqiEICYUDaK4bekSVYIaOW0pFz47sWD9rYsP\nYb9mydeGnvpsLoPe+rZg/enfHczhe5fZwcU55wpMmPczpzzyRcH6n4/Zm4uPapv0+UVrXXef2omT\nD2xWyhk7LtkEFfVX9VqSCkZeltQKqFXWSZKeVDBV9dQS9kvSA5JmSZos6YAUxrxd7h45syA5HdCi\nHnOH9N+u5ARwfq9WzLmjPw12qQrAeU+O49pXJ6c8VudcxXTHiOkFyalyjph1e7/tSk4Ae+1em3lD\nj+WodrsDcNXL33DTG8X+CU6bqBPUFcBoSaMljQE+JrlefE8BfUvZ3w9oGy4XEkz1nHaD3pzGgx/N\nAmDwgA68+n+9KPLMV9JycsSEG45m0PHtARg2biG97x5drnrcOOfS7+SHP+exT+YAcOXRezPrjv5U\nrrTjf9qfPP8g7j09mP7rmS/mc9mLX6ckzh0R+WjmkqoB7cLVGUWmXSjtvJYEUzh3LGbfo8DocEI5\nJM0EjjCz0qbaTmkT3z0ffMcDo74HUt8k992Pazjm3qCqXa1yDjMG993hxOecq7h6DBnFklVBh6vX\n/q8nXVrsmrJrT5z/Myc/HNTKzum+J4NP/J8/xTssI5r4JJ1LMH9Op3A5Pdy2s5pSeBrpRRSeMjtS\nb36zuCA5Pfv71N8v2rtRbb668WgANm3No811I7wm5ZwrpOttHxYkpy+v7Z3S5ARw4J678fpFvQB4\n9sv5PPPFvJRePxlRN/EdlLAcSjA2X1q7mEu6UNIESROWL1++09ebvXwtlw4Lqrx3ndqJQ9tG05lh\nt1pVmTLoGADyDLoM/iCScpxz5c9Rd4/mp7VBY9TXNx7NHnWjeZayc/N6PH5uUNG56Y1pTF60MpJy\nShJpgjKzSxKWPwAHALuk4NI/AM0T1puF24qL4TEz62pmXRs23LlksnFLLr3vHgPA2d1bcEpEPVzy\n1a5ehclhklq5fgsDHvpvpOU55zLf754az5zl6wCYeEMfdq1VNdLy+rRvxKVH7QXACQ99xtpNWyMt\nL1G6H7hZR2pGM38TODfszdcdWFXW/adUOCCsxdStUYXbTtwv6uIAqFO9CuOv7wPAN4tWcePr8faq\ncc7F594PvuOjGcsA+PQvR1J/l+ieVUp05TH70KZh0AG7483vp6VMiP4e1FuS3gyXt4GZwOtJnDcM\n+ALYR9IGZ/F/AAAbsUlEQVQiSb+XNFDSwPCQEcAcYBbwL+D/InoLBa5/bQrrN+cCQZU6nRrWrsa7\nlx0KBG3Bb0wqtrLonKvAxny3nPvDe98vD+xB891qprX8D688vOD1H59Nz6DbUY9mflfC663AfDNb\nVNZJZnZmGfsNuGgnY9sum7fmAfDhlYeRU8IT2VHat3EdHjyzC5cM+5rLXpxEl+a70qJ+ev+BOufi\nsWzNRs57chwAg0/syEEtd0t7DJL4/Jqj6Dn0IyrnpKfxLeqRJP5mZn8ta1u6pGMsvqj9dfhk/jMh\n6MA487a+VKtcKeaInHNR2pKbR9vr3wWgz76NePy8WAbiSamM6GYOFNcW1i/iMiu0v52yPzWrBklp\nv5tHlnG0c668O+zOjwteV4TktD0iSVCS/iRpCsE9pMkJy1zAx/DZSVMH/QqAzbl53BzzUCTOuejc\n88F3Bc86fXdb9n23j6oG9QJwPEFvu+MTlgPN7OyIyswaOTnik6uPBODpL+an/dkE51z0Zi1bUzAg\nwHuXH5qVsxxE9Y7NzOYRdGRYk7AgKf139yqgFvVrcl3/YASpEx76jE1bc2OOyDmXKltz8wpGFh94\neBva7VEn5ojiEWUNCmAiMCH8OTFh3aXAhYe1YY9wNt5eQz8u42jnXHlx7APBQ/lVK+dwTb92ZRxd\ncUWSoMzsuPBnKzNrHf7MX1qXdb5L3qd/DZr6flq7iacSJj90zpVPb0z6gZk/rgFg8s3HxBxNvKLq\nJHFAaUsUZWarKpVyeCMc0HHQW9/y87rNMUfknNtRazZu4bIXJwHwwgXdqF4lux8jiepB3btL2WfA\nURGVm5U6Na/Hsfs35p3JSzhg8AfMG3ps3CE553bAfoOCR0e6t96Nnns1iDma+EWSoMzsyCiu60r2\n0JldeGdyMBzhkBHTubb/vjFH5JzbHo99Mrvg9QsXdI8xkswR9Vh81SVdKelVSa9IulxSNOPCZ7n8\nYUgAHv1kDj+s3BBzRM65ZK1Yu4k7RswA4IMr4hlOLRNF3bH+GaAD8CDwUPj62YjLzFpN6tXg94cE\ng8X3GvpRzNG4VNuam+cTV1ZQB972IQADOjehbaPaMUeTOaIeLLajmbVPWP9Y0rcRl5nVbjyuPU/8\nN+jN50195VNunvHxjGW8+vUiZi5dw6JfNrApHKwYoHKOUPgFu071KhzStgF99m1E3457UKVS9j3M\nWd4lNu3dd3rnGCPJPFEnqK8kdTezLwEkdcOfg4rcl9f2pvuQUTz6yRx+d0grGtXxVtXy4Pmx87n+\ntf8dumqfRrVp3bAWdWtUYfc61dmam8esZWuZvGgVS1dv5I1Ji3lj0uKC4284dl8uONSf5igPVq3f\nUtC0N+qqw5G8aS9R1KOZTwf2ARaEm1oQzAm1lWC0if0jK7wYFWE082Td/s63/OvToCblvfoyl5nx\n55cn88pXhWehueCQVpzXs2VSc/4sXrmBlyYs5L4Pvy+0vW+HPXj47AP8j14Ga3nNOwD8uktT7s2i\n2lOyo5lHnaD2LG2/mc2PrPBiZFOCgm3/+K/oszeX9WkbczSuqH+OnsWd780stO2dSw+hQ5O6O3zN\nWcvW0ueeMYW2/emINvy1b/aORpCpnvtyPjeEM2TPHdI/q75IZESCCgPZFWhOQnOimX0VaaElyLYE\ntWDFeg77ezAE0qSbjqZezaoxR+QAFv68nkPvLDw01cQb+qR0+u7VG7fQ7fZRbNiybYzGD644zG/A\nZ4h1m7bSIZw6/d3LDmXfxtk11l5GJChJg4HzgdkED+hC0LQXy4O62ZagAK5++Rtenhg0H3lTX/wu\neuGrgufVAEZecRh7R5g05q9Yx+F/H12wfmjbBjz7+26RleeS0+mWkazasIXe7XbnifMPijuctMuU\nCQtPA9qY2RFmdmS4+CgSaXTnKdtu8z37ZVpbVF2CVRu20PKadwqS0+96tWLe0GMjTU4Ae9avxbyh\nx3JZ76CJ99Pvf6LlNe+wYu2mSMt1JXvrm8Ws2rAFgMfOza4JCLdX1AlqKlAv4jJcKSTx7mWHAnDj\n61NZt2lrzBFln49m/EinW7bNfjzu+t7cdHz7Us5IvSuO3ptJN22b4PrA2z7kjUk/pDUGB5u25nLJ\nsK8BeHlgDyr5A7mlijpBDQG+lvS+pDfzl4jLdEXs27gOR+7TEIBD/uYP8KbTta9O4XdPBc3KHZrU\nYd7QY9m9djzd/uvVrMq8ocfSa6/6AFz24iQueiGW28FZq9/9nwLB+JkHtfSp8coSdYJ6GvgbMJRg\nANn8xaXZ4+cF7dy/rN9S6B6Ii4aZccjfPmLYuOAJi1tO6MA7lx4ac1SB5y/ozr2ndwLgnclL2O/m\n98nL8xEqovbf739izvJ1ALz8xx4xR1M+RJ2g1pvZA2b2sZmNyV8iLtMVo1KO+M+FwQCUF73wFZsT\nRiZwqZWbZ7S6dgSLfgnGQ3zv8kM5r2fLeIMq4tddmvHJ1cGYzms2baX1dSP830SEcvOMs58YC8AT\n53XNyunbd0TUn9KnkoZI6uHzQcWvW+v6dGwadGc99oFPY46mYtqwOZc2140oWJ886JiMna67Rf2a\nTL+1b8H63je8y+qNW2KMqOI66/EvAWharwa9920UczTlR9QJqgvQHbiDbc17d0VcpivFK3/qCcD3\ny9by2ayfYo6mYlm1fgv73vRewfr3t/ejTvUqMUZUthpVKzH7jv4F6/sPGsmy1RtjjKjimbJoFV/O\n+RmAD648LOZoypdIE1RC1/Ijt7ebuaS+kmZKmiXpmmL2HyFplaRJ4XJT6t9BxVOtciUeD7u2nvX4\nWHL93kNKLFuzkU63buupN3dI/3IzcGulHDF3SH/q1woe5D74jlEsWLE+5qgqBjPj+If+C8CdJ+9P\nzapRD39asUT+P0jSsZL+Iumm/CWJcyoB/wD6Ae2BMyUV1y/3UzPrHC63pjj0CqtP+0Y0rhv0JDvj\nsS9ijqb8W7xyAwffPgqAhrWrMW/oseVu2BpJTLzxaNrtETyXddjfP2b28rUxR1X+/em5oJdk1co5\nnHZQ85ijKX+inrDwEeB04BJAwKlAqePzhQ4GZpnZHDPbDLwIDIgs0Cw06qrDARg/7xe+XvBLzNGU\nXz+s3EDPcO6tNg1rMf76PjFHtHPeu/wwDg67P/e+ewyzlnmS2lGzlq3lvWlLARh/Xfn+dxGXqGtQ\nPc3sXOAXM7sF6AHsncR5TYGFCeuLwm3/c31JkyW9K6lDcReSdKGkCZImLF++fHvjr7BqVq3MXacG\nXY1//c/PfSK8HbBk1YaCiSH3bVyHUVcdEW9AKfLSwB4cslcDAPrcM4Y5XpPaIfmD9l7Xvx11a2b2\nvchMFXWCyp93fL2kJsAWoHGKrv0V0CKcsuNB4PXiDjKzx8ysq5l1bdiwYYqKrhhOObAZtasFbeK/\n+dfYmKMpX5at2UiPIUFy2qdR7YLROiqK5y7oRo/WwQO9R909hoU/+z2p7XHZi18XvL7wsDYxRlK+\nRZ2g3pZUD/g7QUKZB7yQxHk/EIyAnq9ZuK2Ama02s7Xh6xFAFUkNUhF0Nvn82qDPyhdzVjBp4cqY\noykfflm3ueCeU8v6NXn/iorZM2vYhd05oEUwUtmhd37Mj967Lylzlq8tmEBy4g3etLczou7FN9jM\nVprZKwT3ntqZWTK97cYDbSW1klQVOAMoNESSpD0U3omWdDDBe1mR2ndQ8dWuXoUhJ+0HwIn/+Mx7\n9ZVh7aatdBn8ARB0iPj4z0fEG1DEXv2/XuwTDmjb7Y5R/Lxuc8wRZba8POOou4Omvat/tU9Kp1DJ\nRmnrB2tmm8xsVZLHbgUuBt4HpgMvmdk0SQMlDQwPOwWYKukb4AHgDPMbKTvkzINbULNqJQBOe9R7\n9ZVk45ZcOoZz+FSrnMO463qXu956O+L9Kw5jjzpBr88DBn/AWh9wuEQDn5tY8PqiI/eKMZKKIfIJ\nCzNJNs4HlawNm3MLHjJ94YJu9NzLW0sT5eZZoREi5tzRn5wsGonazGh343tsCodD+u62fj5cTxHf\nLFzJgH98Fry++Rjq1vCOESXJlPmgXDlRo2olHvpNFwB+8/hYNibMxJrtzAonp9lZlpwgeE6q6LBI\nPsDsNltz8wqS0+ATO3pySpFIElTiuHvFLVGU6Xbecfs3oXWDWgAc/vePyzg6e7S6dlty+v72flk7\nh09OjgoNi9T6uhH+eELouAeD0SLq1azCOd2TedTTJSOqGtTdpSw+Fl8Gy++R9uPqTTzzxbxYY8kE\nB4QdIgBmDO5bboYvikqlHPHdbf0K1hOTd7Z6Y9IPzFi6BoAvr+0dczQVSyT/20oYg8+nfC8HqlTK\n4fWLegFw0xvTsnrg0GMf+LSg19qUQcdQvUqlmCPKDFUr5zBj8Lbmvu53jIoxmnitWr+Fy16cBMDz\nF3TzfyMpFvVQR1UkXSppeLhcLMkbZzNc5+b1OLFzEyAYODQbm3F+/9R4pi1eDQRTtNfO8FHJ0616\nlUp8c/MxACxdvbHg/ku2yR8g+NC2DejlHYtSLur2ioeBA4F/hsuB4TaX4e49vXPB63OfHBdjJOn3\n1+GTGTVjGQCj/3xEbFO0Z7q6Naow9rqgSeubhSu54Ons6iF7ecJoEc/87uAYI6m4ok5QB5nZeWb2\nUbj8Fjgo4jJdCkhi8qDgG/Kn3//EW98sjjmi9LhjxHT+MyEYBvLtSw6hZdhpxBWvUZ3qfBQOPPzh\n9B/5y/BvYo4oPT75bjmvh6NFjL++T1Y8DxeHqBNUrqSCgagktQa8/3I5Uad6FZ44L3hU4ZJhX7Ns\nTcW+H/XQR9/z2CdzgOBZsI5N68YcUfnQuuEuvHlxcN/ypQmLuP2db2OOKFqrNmwpaFW469RONKzt\no0VEJeoEdTXwsaTRksYAHwFXRVymS6He+zbimPbBFNUH315x70c9+d+53DXyOwAeOftAf1B5O+3f\nrB4vXNANgH99Opd7P/gu5oiiYWZ0uiW479SlRT1OObBZzBFVbJElKEk5BKOZtwUuJZgTah8z8wds\nypnHzt32wHe3Cthj6/mx87n17eBb/12ndqJvxz1ijqh86rlXAx4750AA7h/1Pf8cPSvmiFIv/3kn\ngNf+r1eMkWSHyBKUmeUB/wjH4JscLpuiKs9FK38UgWVrNnH9a1NijiZ1/jN+Ade/NhWA207s6N+I\nd9IxHfbggTODEUnufG8mj46ZHXNEqXPPyJkFPTvz78+6aEXdxDdK0snyO4jlXo2qlQpm4X1+7ALe\nnlz+O00MG7eAv74SJNtBx7fnbB8BICVO6NSEe04LJsMc8u4MHh5d/pPUp98v54GPghrhmxf3oo4/\ndpAWUSeoPwIvA5skrZa0RtLqiMt0EWnTcJeCPzwXv/A1M5aW31/lk/+dy7WvBsnp5uPbc36vVjFH\nVLGcdECzgn8rf3tvBveU43tS835axzlPBJ0ibjquPfs3qxdzRNkj6vmgaptZjplVNbM64XqdKMt0\n0TrpgGac1yOoafS979NyOdLEvR98V3DPafCADvzWk1MkTjqgWUFz3wOjvmfQm9Nijmj7rVy/mSPu\nGg3A8Z2a8LtD/N9KOkU9ksT/3FEvbpsrX24Z0JHOzYNvkQffMYo1G7fEHFHyrnttCveP+h6Au0/t\nxDk9WsYbUAV3QqcmPB52snnq83n8KWG+pEy3cUsunW8NxmLcs35NHgyTrUufqEYzry5pN6CBpF0l\n7RYuLYGmUZTp0uu1/+tJ7eqVAdhv0MhyMT3H6Y9+wQtjFwDw7/MP4mTvEJEWfdo34uWBPQB4d+pS\njr5nTMwRlW1Lbh7tbnyvYH10BZ85OVNFVYP6IzARaBf+zF/eAB6KqEyXRpKYfPO2nkztbnyPDZsz\nM0nl5Rn73PAuY+f+DAQ3uY9st3vMUWWXg1ruxodXBp1svl+2lpbXvENuhs4ntXlrHm2vf7dgffYd\n/X2kiJhENZr5/WbWCvizmbU2s1bh0snMPEFVEJKYkzA/0L43vcfqDGvuW7dpK62vG1EwE+yX1/b2\nm9wx2Wv3XfjqxqML1ttcN4JV6zPr38uGzbnsfcO25DQri+f/ygRRd5J4UFJPSb+RdG7+EmWZLr1y\ncgonqf0HjWTRL+tjjGibyYtW0uHm9wvWZwzuyx51feDXOO1Wqyrf375tPqlOt47ki9krYoxom+Vr\nNrHvTdua9Wbf0Z/KWT7/V9yi7iTxLMEEhYcQDBJ7EFDmPPSufMnJEXOH9KfBLlUBOORvHzN65rJY\nYxr89rec8FAwBUSTutWZO6S/z9WTIapUymHukP6026M2AGf+60v+/HK8g8xOmPczB93+IQASzB3S\n32tOGUBRjq0maTrQ3jJkALeuXbvahAnZNSVAup3/73GMnrkcgDMOas7Qk/dPa/kbt+QWurl90ZFt\nuPpX7dIag0veI2NmM/TdGQXr0275FbWqVU5rDEPenc6jY4JBgjs1r8cbF/kQRlGTNNHMyqysRJ2g\nXgYuNbMlkRWyHTxBpcfjn87htnemF6zPGNw3LbWX576czw2vTy1Y/+CKw2jbqHbk5bqds2DFeg77\n+7YhOq/p146Bh7cp5YzU2Lw1r9D9pkuP2osrj9kn8nJd5iSoj4HOwDigYBw+MzshskJL4Qkqfeb+\ntI4jwwccAa7vvy9/OKx1JGUtXrmBnkM/Klhv3bAWH15xODneRFNumBknPfw5Xy9YWbDt4z8fQauI\n5uN6YewCrksYU3LEpYfSvomPIZAumZKgDi9uu5nF8iCEJ6j0MjP63DOG2cvXFWx746JedGqeml50\nq9ZvoefQUaxL6N4+fGAPurbcLSXXd+k3bfEqjn1g24jhEoy7rk/K5lz67sc1HHPvJwXrDWtXY+y1\nvf3LTJplRIIKA2nEtll0x5lZbHfPPUHFY+oPqwpNUwDw6DkH8qsOOzatxdcLfuHX//y80Lbze7Zk\n0AkddjhGl1nuGTmzYHDWfMP+0J0eberv0PU++W55wSSD+V4e2IOD/MtMLDIiQUk6Dfg7MBoQcChw\ntZkNT+LcvsD9QCXgcTMbWmS/wv39gfXA+Wb2VWnX9AQVr6LNKgADOjfh7O57ckCLXUvsNZWbZ0yY\n9zNvfLO4YCSIfCd0asL9Z3T2BykrIDPjutemMmxc4d/5qQc248QuTenWarcSu4Hn5RlfL1zJsHEL\nGD5xUaF9Nxy7LxccGk1zs0tOpiSob4Cj82tNkhoCH5pZpzLOqwR8BxwNLALGA2ea2bcJx/QnmASx\nP9ANuN/MupV2XU9QmeGzWT9xybCv+Xnd5kLb69WsQsv6taheJYfqVSqxdNVGFv68vlATHkCjOtW4\npl87ft3FhyrKFu9NXcKtb33L4lWFByfepVplWjeshSSa1K3O7OVrWbZmEyuLPABct0YV7j29E0e1\na5TOsF0Jkk1QUffnzCnSpLeC5J69OhiYZWZzACS9CAwAvk04ZgDwTNiF/UtJ9SQ1zpQeg65kvfZq\nUDCiwLTFqxg+cRHzV6znl/WbWbpqI3VrVKFa5Vx2qVaZFvVr0ahONdo3rsMxHfagU7O6XlvKQn07\nNqZvx8aYGd8uWc27U5YyfclqFv2ygaWrNlK1cg7rN22lauUcdqtZlf2a1qV1g1qcdECzlN3zdOkX\ndYJ6T9L7wLBw/XTg3VKOz9cUWJiwvoigllTWMU2BQglK0oXAhQAtWrRIOnCXHh2a1KVDk7pxh+HK\nCUn+byaLRJqgzOxqSScRjCQB8JiZvRZlmcXE8BjwGARNfOks2znn3I6LJEFJ2gtoZGafmdmrwKvh\n9kMktTGzsuaA/gFonrDeLNy2vcc455wrp6KqQd0HXFvM9lXhvuPLOH880FZSK4KkcwbwmyLHvAlc\nHN6f6gasKuv+08SJE3+SND+J+EvSAPhpJ86vCLL9M8j29w/+GWT7+4ed/wz2TOagqBJUIzObUnSj\nmU0JJy0slZltlXQx8D5BN/MnzWyapIHh/keAEQQ9+GYRdDP/bRLXbbg9b6IoSROS6XlSkWX7Z5Dt\n7x/8M8j29w/p+wyiSlCldZupkcwFzGwEQRJK3PZIwmsDLtqh6JxzzmW8qKbbmCDpD0U3SrqAYGZd\n55xzrlRR1aAuB16TdBbbElJXoCrw64jKTIfH4g4gA2T7Z5Dt7x/8M8j29w9p+gyiHkniSKBjuDrN\nzD4q7XjnnHMuX+SDxTrnnHM7ItIp351zzrkd5QkqSZL6SpopaZaka+KOJ50kNZf0saRvJU2TdFnc\nMcVBUiVJX0t6O+5Y4hCOdzlc0gxJ0yX1iDumdJN0Rfh/YKqkYZKqxx1T1CQ9KWmZpKkJ23aT9IGk\n78Ofu0ZRtieoJISjq/8D6Ae0B86U1D7eqNJqK3CVmbUHugMXZdn7z3cZML3Moyqu+4H3zKwd0Iks\n+ywkNQUuBbqaWUeCZzTPiDeqtHgK6Ftk2zXAKDNrC4wK11POE1RyCkZXN7PNQP7o6lnBzJbkz7Vl\nZmsI/jA1jTeq9JLUDDgWeDzuWOIgqS5wGPAEgJltNrOVpZ9VIVUGakiqDNQEFsccT+TM7BPg5yKb\nBwBPh6+fBk6MomxPUMkpaeT0rBOOBNIFGBtvJGl3H/AXIC/uQGLSClgO/Dts5nxcUq24g0onM/sB\nuAtYQDBrwiozGxlvVLFplDC03FIgkom2PEG5pEnaBXgFuNzMVscdT7pIOg5YZmbZ/JB5ZeAA4GEz\n6wKsI6JmnUwV3mcZQJCsmwC1JJ0db1TxC0f1iaQ7uCeo5GT9yOmSqhAkp+fDEeqzSS/gBEnzCJp3\nj5L0XLwhpd0iYJGZ5dechxMkrGzSB5hrZsvNbAvBLA09Y44pLj9KagwQ/lxWxvE7xBNUcgpGV5dU\nleDG6Jsxx5Q2CqawfQKYbmb3xB1PupnZtWbWzMxaEvzuPzKzrPrmbGZLgYWS9gk39abwDNfZYAHQ\nXVLN8P9Eb7Kso0iCN4HzwtfnAW9EUUjUM+pWCCWNrh5zWOnUCzgHmCJpUrjtunBAX5c9LgGeD7+k\nzSGJGQQqEjMbK2k48BVBz9avyYJhjyQNA44AGkhaBNwMDAVekvR7YD5wWiRl+0gSzjnnMpE38Tnn\nnMtInqCcc85lJE9QzjnnMpInKOeccxnJE5RzzrmM5N3MnYuRpPoEg20C7AHkEgwpBLDezLL1QVDn\nvJu5c5lC0iBgrZndFXcszmUCb+JzLkNJWhv+PELSGElvSJojaaiksySNkzRFUpvwuIaSXpE0Plx6\nxfsOnNs5nqCcKx86AQOBfQlG9djbzA4mmP7jkvCY+4F7zewg4GSydGoQV3H4PSjnyofx+dMbSJoN\n5E/zMAU4MnzdB2gfDBMHQB1Ju5jZ2rRG6lyKeIJyrnzYlPA6L2E9j23/j3OA7ma2MZ2BORcVb+Jz\nruIYybbmPiR1jjEW53aaJyjnKo5Lga6SJkv6luCelXPllnczd845l5G8BuWccy4jeYJyzjmXkTxB\nOeecy0ieoJxzzmUkT1DOOecykico55xzGckTlHPOuYz0/ykfhd2XK5j7AAAAAElFTkSuQmCC\n", "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "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", " ax1.step(result.time, \n", " np.hstack((result.initial_amps[:, j], result.initial_amps[-1, j])), \n", " where='post')\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", " ax2.step(result.time, \n", " np.hstack((result.final_amps[:, j], result.final_amps[-1, j])), \n", " where='post')\n", " \n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Versions" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "
SoftwareVersion
QuTiP4.1.0
Numpy1.11.3
SciPy0.18.1
matplotlib2.0.0
Cython0.25.2
Number of CPUs4
BLAS InfoINTEL MKL
IPython5.1.0
Python3.6.0 |Anaconda 4.3.1 (64-bit)| (default, Dec 23 2016, 12:22:00) \n", "[GCC 4.4.7 20120313 (Red Hat 4.4.7-1)]
OSposix [linux]
Fri Jul 14 17:15:34 2017 BST
" ], "text/plain": [ "" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from qutip.ipynbtools import version_table\n", "\n", "version_table()" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.6.0" } }, "nbformat": 4, "nbformat_minor": 0 }