{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Calculation of control fields for Lindbladian dynamics using L-BFGS-B algorithm" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Christian Arenz (christianarenz.ca@gmail.com), Alexander Pitchford (alex.pitchford@gmail.com)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Example to demonstrate using the control library to determine control pulses using the ctrlpulseoptim.optimize_pulse function. 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 'Trace difference' norm.\n", "\n", "This in an open quantum system example, with a single qubit subject to an amplitude damping channel. The target evolution is the Hadamard gate. For a $d$ dimensional quantum system in general we represent the Lindbladian\n", "as a $d^2 \\times d^2$ dimensional matrix by creating the Liouvillian superoperator. Here done for the Lindbladian that describes the amplitude damping channel. Similarly the control generators acting on the qubit are also converted to superoperators. The initial and target maps also need to be in superoperator form. \n", "\n", "The user can experiment with the strength of the amplitude damping by changing the gamma variable value. If the rate is sufficiently small then the target fidelity can be achieved within the given tolerence. The drift Hamiltonian and control generators can also be swapped and changed to experiment with controllable and uncontrollable setups.\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", "For more background on the pulse optimisation see:\n", "[QuTiP overview - Optimal Control](http://nbviewer.ipython.org/github/qutip/qutip-notebooks/blob/master/examples/example-optimal-control-overview.ipynb) " ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "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": {}, "outputs": [], "source": [ "from qutip import Qobj, identity, sigmax, sigmay, sigmaz, sigmam, tensor\n", "from qutip.superoperator import liouvillian, sprepost\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", "\n", "example_name = 'Lindblad'" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Defining the physics" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "Sx = sigmax()\n", "Sy = sigmay()\n", "Sz = sigmaz()\n", "Sm = sigmam()\n", "Si = identity(2)\n", "#Hadamard gate\n", "had_gate = hadamard_transform(1)\n", "\n", "# Hamiltonian\n", "Del = 0.1 # Tunnelling term\n", "wq = 1.0 # Energy of the 2-level system.\n", "H0 = 0.5*wq*sigmaz() + 0.5*Del*sigmax()\n", "\n", "#Amplitude damping#\n", "#Damping rate:\n", "gamma = 0.1\n", "L0 = liouvillian(H0, [np.sqrt(gamma)*Sm])\n", "\n", "#sigma X control\n", "LC_x = liouvillian(Sx)\n", "#sigma Y control\n", "LC_y = liouvillian(Sy)\n", "#sigma Z control\n", "LC_z = liouvillian(Sz)\n", "\n", "#Drift\n", "drift = L0\n", "#Controls - different combinations can be tried\n", "ctrls = [LC_z, LC_x]\n", "# Number of ctrls\n", "n_ctrls = len(ctrls)\n", "\n", "# start point for the map evolution\n", "E0 = sprepost(Si, Si)\n", "\n", "# target for map evolution\n", "E_targ = sprepost(had_gate, had_gate)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Defining the time evolution parameters" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "# Number of time slots\n", "n_ts = 10\n", "# Time allowed for the evolution\n", "evo_time = 2" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Set the conditions which will cause the pulse optimisation to terminate" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "# Fidelity error target\n", "fid_err_targ = 1e-3\n", "# Maximum iterations for the optisation algorithm\n", "max_iter = 200\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": {}, "outputs": [], "source": [ "# pulse type alternatives: RND|ZERO|LIN|SINE|SQUARE|SAW|TRIANGLE|\n", "p_type = 'RND'" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Give an extension for output files" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "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": {}, "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 = [[[2], [2]], [[2], [2]]], shape = (4, 4), type = super, isherm = False\n", "Qobj data =\n", "[[-0.10+0.j 0.00-0.05j 0.00+0.05j 0.00+0.j ]\n", " [ 0.00-0.05j -0.05+1.j 0.00+0.j 0.00+0.05j]\n", " [ 0.00+0.05j 0.00+0.j -0.05-1.j 0.00-0.05j]\n", " [ 0.10+0.j 0.00+0.05j 0.00-0.05j 0.00+0.j ]]\n", "Control 1 dynamics generator:\n", "Quantum object: dims = [[[2], [2]], [[2], [2]]], shape = (4, 4), type = super, isherm = False\n", "Qobj data =\n", "[[ 0.+0.j 0.+0.j 0.+0.j 0.+0.j]\n", " [ 0.+0.j 0.+2.j 0.+0.j 0.+0.j]\n", " [ 0.+0.j 0.+0.j 0.-2.j 0.+0.j]\n", " [ 0.+0.j 0.+0.j 0.+0.j 0.+0.j]]\n", "Control 2 dynamics generator:\n", "Quantum object: dims = [[[2], [2]], [[2], [2]]], shape = (4, 4), type = super, isherm = False\n", "Qobj data =\n", "[[ 0.+0.j 0.-1.j 0.+1.j 0.+0.j]\n", " [ 0.-1.j 0.+0.j 0.+0.j 0.+1.j]\n", " [ 0.+1.j 0.+0.j 0.+0.j 0.-1.j]\n", " [ 0.+0.j 0.+1.j 0.-1.j 0.+0.j]]\n", "Initial state / operator:\n", "Quantum object: dims = [[[2], [2]], [[2], [2]]], shape = (4, 4), type = super, 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 = [[[2], [2]], [[2], [2]]], shape = (4, 4), type = super, isherm = True\n", "Qobj data =\n", "[[ 0.5 0.5 0.5 0.5]\n", " [ 0.5 -0.5 0.5 -0.5]\n", " [ 0.5 0.5 -0.5 -0.5]\n", " [ 0.5 -0.5 -0.5 0.5]]\n", "INFO:qutip.control.pulseoptim:Initial amplitudes output to file: ctrl_amps_initial_Lindblad_n_ts10_ptypeRND.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_Lindblad_n_ts10_ptypeRND.txt\n" ] } ], "source": [ "# Note that this call will take the defaults\n", "# dyn_type='GEN_MAT'\n", "# This means that matrices that describe the dynamics are assumed to be\n", "# general, i.e. the propagator can be calculated using:\n", "# expm(combined_dynamics*dt)\n", "# prop_type='FRECHET'\n", "# and the propagators and their gradients will be calculated using the\n", "# Frechet method, i.e. an exact gradent\n", "# fid_type='TRACEDIFF'\n", "# and 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(drift, ctrls, E0, E_targ, 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", " out_file_ext=f_ext, init_pulse_type=p_type, \n", " log_level=log_level, gen_stats=True)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Report the results" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "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:01.340385\n", "Wall time computing Hamiltonians: 0:00:00.020057 (1.50%)\n", "Wall time computing propagators: 0:00:01.178558 (87.93%)\n", "Wall time computing forward propagation: 0:00:00.005692 (0.42%)\n", "Wall time computing onward propagation: 0:00:00.004413 (0.33%)\n", "Wall time computing gradient: 0:00:00.091266 (6.81%)\n", "\n", "**** Iterations and function calls ****\n", "Number of iterations: 201\n", "Number of fidelity function calls: 265\n", "Number of times fidelity is computed: 265\n", "Number of gradient function calls: 265\n", "Number of times gradients are computed: 265\n", "Number of times timeslot evolution is recomputed: 265\n", "\n", "**** Control amplitudes ****\n", "Number of control amplitude updates: 264\n", "Mean number of updates per iteration: 1.3134328358208955\n", "Number of timeslot values changed: 2640\n", "Mean number of timeslot changes per update: 10.0\n", "Number of amplitude values changed: 5280\n", "Mean number of amplitude changes per update: 20.0\n", "------------------------------------\n", "Final evolution\n", "Quantum object: dims = [[[2], [2]], [[2], [2]]], shape = (4, 4), type = super, isherm = False\n", "Qobj data =\n", "[[ 0.49398354 -9.35884799e-17j 0.43895309 +4.70500365e-04j\n", " 0.43895309 -4.70500365e-04j 0.50438676 +1.16197332e-16j]\n", " [ 0.43683752 -2.24436864e-03j -0.44201112 -3.02814020e-03j\n", " 0.43197985 +3.42846832e-03j -0.43695905 +4.27534789e-03j]\n", " [ 0.43683752 +2.24436864e-03j 0.43197985 -3.42846832e-03j\n", " -0.44201112 +3.02814020e-03j -0.43695905 -4.27534789e-03j]\n", " [ 0.50601646 +1.18358316e-16j -0.43895309 -4.70500365e-04j\n", " -0.43895309 +4.70500365e-04j 0.49561324 -1.33018845e-16j]]\n", "\n", "********* Summary *****************\n", "Initial fidelity error 0.7448841629400895\n", "Final fidelity error 0.005876671350497931\n", "Final gradient normal 0.00012944172414442896\n", "Terminated due to Iteration or fidelity function call limit reached\n", "Number of iterations 201\n", "Completed in 0:00:01.340385 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(\"Initial fidelity error {}\".format(result.initial_fid_err))\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": 10, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xm8HFWd9/HPNyEQIAkBLiYhEAIIaiBsXrbAKIw6DwQw\nL5GdwSHqAM6Ig6PO4/YoiAqOuIIKUUFAWUYWUTbZwQkQSJjIhQBOiAkQgoEAWUyAJPyeP+pc0rne\npW5udVfd9Pf9evXrdtdyzq+r6/av69SpU4oIzMzMqmZA2QGYmZl1xgnKzMwqyQnKzMwqyQnKzMwq\nyQnKzMwqyQnKzMwqyQnKmoakWyT9UzfzL5T0/3KWdY+kjxcXXf1JGispJG1QdixmeThBWb8maa6k\n9+dZNiIOjYhL03onS/rvDvNPi4iz6xFnXzm5WDNygjJbTzh52frGCcrWG+1HRZLOk/SKpD9LOrRm\n/j2SPi7pXcCFwP6Slkl6Nc3/haSvp+ebS7pR0ouprBslbZMzjoGSvijpaUlLJc2QtG2aN0HSw5IW\np78TOsR3tqSpab3bJLWk2felv6+mmPdP73eqpO9JWgScKWmApC9LmidpoaTLJG2WM+7P18Q8S9KH\nOmzb9rpelTQnvZeTJT2b6vqnmuV/kZpMb0/l3StpuzRPqZyFkpZIapO0a54Yrbk4Qdn6Zl/gKaAF\n+E/g55JUu0BEPAGcBjwQEUMiYngn5QwALgG2A8YAK4ALcsbw78DxwERgGPBRYLmkLYCbgB8CWwLf\nBW6StGXNuicAk4G3ARsCn03T35P+Dk8xP1DzfucAI4BvACenx8HADsCQXsT9NPB3wGbAWcAvJY2q\nmb8v8GiK/QrgKmBv4O3APwIXSBpSs/yJwNlkn8VM4Fdp+j+k97NzqusYYFHOGK2JOEHZ+mZeRPw0\nIlYDlwKjyL68eyUiFkXEtRGxPCKWkn35vzfn6h8HvhwRT0XmjxGxCDgM+N+IuDwiVkXElcCTwBE1\n614SEX+KiBXAfwF79FDX8xFxfipvBVlS+G5EzImIZcAXgOPyNP9FxK8j4vmIeDMirgb+F9inZpE/\nR8QladteDWwLfC0iXo+I24A3yJJVu5si4r6IeB34EtkR67bASmAo8E5AEfFERCzoKT5rPk5Qtr55\nof1JRCxPT4d0sWyXJG0i6aLUVLaErIltuKSBOVbfluxopKOtgXkdps0DRte8fqHm+XJ6jv3ZHuqY\nB2xAjiQt6SOSZqYmvFeBXcmOftr9peb5CoCI6DitNt63YkvJ8mVg64i4i+yo7kfAQklTJA3rKT5r\nPk5Q1qx6Gsb/M8A7gH0jYhhrmtjU9SpveRbYsZPpz5M1GdYaA8zPUWZX8Xac3rGOMcAq1k4ufyOd\nH/op8Elgy9Ts+Rj53m9Xtq0pfwiwRYqPiPhhRLwbGEfW1Pe5PtRj6yknKGtWfwG2kbRhF/OHkh0R\nvJrOHX21F2X/DDhb0k6pQ8Bu6TzTzcDOkk6QtIGkY8m+oG/MUeaLwJtk55W6cyXwaUnbp6TwTeDq\niFjVw3qbkiW7FwEkTSY7guqLiZIOTNv4bODBiHhW0t6S9pU0CPgr8BrZezNbixOUNau7gMeBFyS9\n1Mn87wMbAy8BDwK39qLs75KdP7oNWAL8HNg4nYc6nOzobBHwH8DhEdFZ/WtJzZXfAKamJrj9ulj0\nYuBysibJP5N9+Z+eo/xZwHeAB8iS93hgak/r9eAKssT+MvBuso4UkHUc+SnwClkT5CLg232sy9ZD\n8g0Lzaxokn4BPBcRXy47Fuu/fARlZmaVVGiCkrSxpHcUWaaZmTWnwpr4JB0BnAdsGBHbS9qD7BqJ\nDxZSgZmZNZUij6DOJLuo71WAiJgJbF9g+WZm1kSKHFxyZUQs7jiqTF8LlXQxWc+nhRGRq9trS0tL\njB07tq9Vm5lZHcyYMeOliNiqp+WKTFCPSzoBGChpJ+BTwP0FlPsLsqvOL8u7wtixY5k+fXoBVZuZ\nWdEkdRxRpVNFNvGdDuwCvE52seAS4Iy+FhoR95FdR2FmZk2ksCOodCHhl9Kj4SSdApwCMGbMmDJC\nKMwV057hhpl5Rr+pr0l7jOaEffv3tjSz/qvPCUrS7+jmXFOjevFFxBRgCkBra2u/vvr4hpnzmbVg\nCeNGlTd+5qwFSwCcoMysNEUcQZ2X/h4JjAR+mV4fTw8DVFrXxo0axtWn7l9a/cde9EDPC5k1iFsV\nmlOfE1RE3Asg6TsR0Voz63eS3FPBzPrMrQrNqchefJtK2iEi5gBI2p5shOQ+kXQlcBDQIuk54KsR\n8fO+lltl71t+MwesuBsuyXWn7rr4yqLFTN34YKC8ozizWm5VaD5FJqhPA/dImkN2D5ntgFP7WmhE\nHN/XMvqbA1bczdiVc4A9S4shq9/MrDxF9uK7NV3/9M406cl0q2dbB3MH7cAuk28qr/5vHlha3WZm\nUGCCkvSRDpN2l0RE5L7A1szMrF2RTXx71zwfDLwPeIRejABhZmbWrsgmvrXu2ilpOHBVUeVb4y1/\nY3XpJ4bdrdesedXzhoV/xaOZ91stQzZikw0HlhrDrAVLKnHti5mVo8hzULUjSgwAxgG/Lqp8a6wR\nQwczYuhgrp7sbr1mVo4iz0GdV/N8FTAvIp4rsHwzM2siRTbxTYyIe9NjakQ8J+lbBZZvZmZNpMgE\n9YFOph1aYPlmZtZEihjN/BPAvwA7SHq0ZtZQYGpfyzczs+ZUxDmoK4BbgHOAz9dMXxoRvtGgmZmt\nkyISVETEXEn/2nGGpC2cpMzMbF0UdQR1ODCDrJu5auYFsEMBdZiZWZMp4n5Qh6e/vijX1iu+SZ5Z\nuYroJLFXd/Mj4pG+1mFWBt8kz6xcRTTxfaebeQH8fQF1mJXCN8kzK08RTXwHFxGImZlZrSLH4htM\ndj3UgWRHTn8ALoyI14qqw6wZzVqwpPQjKZ8HszIUORbfZcBS4Pz0+gTgcuDoAuswayqT9hhddgg+\nD2alKTJB7RoR42pe3y1pVoHlN8RZv3ucWc8vKTWGz76xuvRbXVg1nLDvmNITQ9lHb9a8ihyL7xFJ\n+7W/kLQvML3A8pvGJhsOpGXIRmWHYWZWqiKPoN4N3C/pmfR6DPCUpDay0SZ2K7CuuvnqEbuUHQJc\nslnZEZiZla7IBHVIgWWZmVmTKyxBRcQ8SZsD29aW6wt1zcxsXRTZzfxs4GTgadbc+r3/Xah7y+fh\nhbZyY3ihDUaOLzcGswp53/KbOWDF3aU2f39l0WKmbnwwUN6F282myCa+Y4AdI+KNAstsTiPHw/ij\nyo7CrDIOWHE3Y1fOAfYsLYasfmukIhPUY8BwYGGBZTbeoeeWHYHZGtMvgbZrSg2hKkcOcwftwC6T\nbyqv/m8eWFrdzarIBHUO8D+SHgNeb58YER8ssA6z5tJ2TelNvj5ysLIUmaAuBb4FtAFvFliuWXMb\nOR585GBNqMgEtTwiflhgeW+RdAjwA2Ag8LOIcDucmdl6rsgE9QdJ5wC/Ze0mvj51M5c0EPgR8AHg\nOeBhSb+NiH43jJKZmeVXZIJq716zX820IrqZ7wPMjog5AJKuAiYBTlBmZuuxIi/Urdd9oUYDz9a8\nfg7Yt+NCkk4BTgEYM8ajLpuZ9XdFHkEh6TBgF2Bw+7SI+FqRdXQlIqYAUwBaW1ujh8XNzKziChvN\nXNKFwLHA6YDI7gO1XQFFzycbPqndNmmamZmtx4q83caEiPgI8EpEnEV2Vd/OBZT7MLCTpO0lbQgc\nR9YRw8zM1mNFNvGtSH+XS9oaWASM6muhEbFK0ieB35N1M784Ih7va7lmZlZtRSaoGyUNB74NPELW\ng++nRRQcETcDNxdRlplZf3XFtGe4YWb5ZzjGbT2sIffOK7IX39np6bWSbgQGR8Tioso3M2t2N8yc\nz6wFSxg3aljZoTREob342kXE69RcrGtmZsUYN2oYV5/aHLf8KLKThJmZWWHqcgRl64kX2uCSw0qr\nviq3eTCzcvQ5QUnaq7v5vuV7P1WBGyb6Ng9mza2II6jvdDOv/93y3TKtk7NHiXybB7Pm1ucEVccx\n+MzMrIkVdg5K0iDgE8B70qR7gIsiYmVRdZiZWfMospPET4BBwI/T65PStI8XWIc1meVvrObYix4o\npe5mut7ErIqKTFB7R8TuNa/vkvTHAsu3JtMyZCNeWlbe5XTjRg1j0h6jS6vfrNkVmaBWS9oxIp4G\nkLQDsLrA8q3JjBg6mBFDB3P1ZHczN2tGRSaozwF3S5pDdruN7YByu4GZWSHKbGoF+Owbq9lkw4Gl\n1W/lKCRBSRpANpr5TsA70uSn0pBHZtaPld3UCrDJhgNpGbJRqTFY4xWSoCLiTUk/iog9gUeLKNPM\nqqESTa2XbFZe3VaaIpv47pT0YeC6iPAt183WJyUPe8ULbTByfHn1J2U3dTZbz9IiE9SpwL8DqyS9\nRnYeKiKiebam2fqoAsNeMXJ86XFUoamz2XqWFnk/qKFFlWVmFVKBYa+qoBJNnU2msNttSLozzzQz\nM7M8ihjNfDCwCdAiaXOypj2AYUDzHIuamVmhimjiOxU4A9gamMGaBLUEuKCA8s3KMf0SaLum3Bgq\n0jnArAxFjGb+A+AHkk6PiPMLiMmsGtquKT9BVKBzgFlZiuwkcb6kCcDY2nIj4rKi6jBruJHjYfJN\nZUdh1pSKvN3G5cCOwEzWjMEXgBOUmZn1WpHXQbUC43yRrplZnVThvChkLQuHnlv3agrrZg48Bows\nsDwzM6vVfl60SRR5BNUCzJL0EPDW5dYR8cEC6zAza25NdF60yAR1ZoFlmZlZkyuyF9+9kkYAe6dJ\nD0XEwqLKNzOz5lLkUEfHAA8BRwPHANMk+QIOMzNbJ0U28X0J2Lv9qEnSVsAdwDp3OZF0NFnT4buA\nfSJiegFxWn9S5m0eyr5I16zJFZmgBnRo0ltE34/QHgOOBC7qYznWH5U9goJHcbCOfF+shioyQd0q\n6ffAlen1scAtfSkwIp4AkNTTorY+8m0erEqq8GOlyX40FdlJ4nOSjgQOTJOmRMT1RZXfE0mnAKcA\njBkzplHVmlmz8A+mhividhtvB0ZExNSIuA64Lk0/UNKOEfF0D+vfQecX+H4pIm7IG0dETAGmALS2\ntno0CzOzfq6II6jvA1/oZPriNO+I7laOiPcXEMNaZsyY8ZKkeX0oogV4qah46sQxFsMxFqM/xAj9\nI85miHG7PAsVkaBGRMTfjL0REW2SxhZQfq9FxFZ9WV/S9IhoLSqeenCMxXCMxegPMUL/iNMxrlHE\ndVDDu5m3cV8KlvQhSc8B+wM3pU4YZmbWBIpIUNMl/XPHiZI+TnaH3XUWEddHxDYRsVFEjIiI/9OX\n8szMrP8ooonvDOB6SSeyJiG1AhsCHyqg/DJMKTuAHBxjMRxjMfpDjNA/4nSMiYq6fZOkg4Fd08vH\nI+KuQgo2M7OmVFiCMjMzK1KRNyw0MzMrTFMlKEmHSHpK0mxJn+9kviT9MM1/VNJeeddtcJwnpvja\nJN0vafeaeXPT9JmS6ja4bo4YD5K0OMUxU9JX8q7bwBg/VxPfY5JWS9oizav7dpR0saSFkh7rYn7p\n+2OOGKuwL/YUY+n7Ys44y94ft5V0t6RZkh6X9G+dLNPYfTIimuIBDASeBnYg68DxR2Bch2Umko0f\nKGA/YFredRsc5wRg8/T80PY40+u5QEsFtuVBwI3rsm6jYuyw/BHAXQ3eju8B9gIe62J+FfbHnmIs\ndV/MGWOp+2LeOCuwP44C9krPhwJ/Kvs7spmOoPYBZkfEnIh4A7gKmNRhmUnAZZF5EBguaVTOdRsW\nZ0TcHxGvpJcPAtvUKZZ1jrFO69YzxuNZM9BxQ0TEfcDL3SxS+v7YU4wV2BfzbMeuNPL/urdxlrE/\nLoiIR9LzpcATwOgOizV0n2ymBDUaeLbm9XP87cbvapk86xalt3V9jLVHjQ/gDkkzlA2gWw95Y5yQ\nmgFukbRLL9dtVIxI2gQ4BLi2ZnIjtmNPqrA/9kYZ+2JeZe6LvVKF/VHZKEB7AtM6zGroPlnk7Tas\nwZR17f8Ya0aQBzgwIuZLehtwu6Qn0y+3RnsEGBMRyyRNBH4D7FRCHHkcAUyNiNpft1XZjv2C98VC\nlbo/ShpClhzPiIgl9agjr2Y6gpoPbFvzeps0Lc8yedYtSq66JO0G/AyYFBGL2qdHxPz0dyFwPdmh\nd8NjjIglEbEsPb8ZGCSpJc+6jYqxxnF0aE5p0HbsSRX2xx6VvC/2qAL7Ym+Vtj9KGkSWnH4V2d0p\nOmrsPlnPk25VepAdLc4BtmfNSbxdOixzGGufAHwo77oNjnMMMBuY0GH6psDQmuf3A4eUFONI1lxn\ntw/wTNquDdmWeesBNiM7L7Bpo7djKn8sXZ/cL31/zBFjqftizhhL3Rfzxln2/pi2yWXA97tZpqH7\nZNM08UXEKkmfBH5P1uPk4oh4XNJpaf6FwM1kvVRmA8uByd2tW2KcXwG2BH6s7G7DqyIbWXgE2bBT\nkO0wV0TErSXFeBTwCUmrgBXAcZHtyQ3ZljljhGw4rtsi4q81qzdkO0q6kqyHWYuyQZG/Cgyqia/0\n/TFHjKXuizljLHVf7EWcUOL+CBwAnAS0SZqZpn2R7EdIKfukR5IwM7NKaqZzUGZm1o84QZmZWSU5\nQZmZWSU5QZmZWSU5QZmZWSU1TTdzs6qQtCVwZ3o5ElgNvJheL4+ICaUEZlYx7mZuViJJZwLLIuK8\nsmMxqxo38ZlViKRl6e9Bku6VdIOkOZLOTfdeeijdF2jHtNxWkq6V9HB6HFDuOzArjhOUWXXtDpwG\nvIvsCv+dI2IfsnHvTk/L/AD4XkTsDXw4zTNbL/gclFl1PRwRCwAkPQ3clqa3AQen5+8HxqVhcACG\nSRoSaXBUs/7MCcqsul6vef5mzes3WfO/OwDYLyJea2RgZo3gJj6z/u021jT3IWmPEmMxK5QTlFn/\n9imgNd0tdhbZOSuz9YK7mZuZWSX5CMrMzCrJCcrMzCrJCcrMzCrJCcrMzCrJCcrMzCrJCcrMzCrJ\nCcrMzCrJCcrMzCrJCcrMzCrJCcrMzCrJCcrMzCrJCcrMzCrJCcoqSdIYScskDVzH9ZdJ2qHgmH4h\n6etFltkoksZKCkm+B5z1G05QVghJJ0tqk7Rc0guSfiJpeC/Wnyvp/e2vI+KZiBgSEavXJZ607px1\nWXddSRol6eeSFkhaKulJSWdJ2rSP5dY9uUg6UNL9khZLelnSVEl716s+szycoKzPJH0G+BbwOWAz\nYD9gO+B2SRuWGVujSNoCeADYGNg/IoYCHyDbHjs2oP51Tl6ShgE3AucDWwCjgbNY+46+Zo0XEX74\nsc4PYBiwDDimw/QhwIvAR9PrM4FrgKuBpcAjwO5p3uVktzFfkcr6D2AsEMAGaZl7gK8D96dlfgds\nCfwKWAI8DIytqT+At6fnE4FZqd75wGdrljscmAm8msrerWbeninOpSnuq4Cvd7Edvg60AQO62VYT\nUpyL098JNfPuAc4Gpqb6bgNa0rxn0vtZlh77AyenZb8HLEr1DwC+DMwDFgKXAZulMtbanh3iagVe\n7eFz/ijwBPAK8Htgu5p5HwCeTO/rAuBe4OM1n/sva5bt+LluBvwcWJA+m68DA9O8k4H/Bs5L9f4Z\nOLSmrC2AS4Dn0/zf5Pxc/2+qaynwFPC+sv+P/Ohivys7AD/69wM4BFjVxRffpcCV6fmZwErgKGAQ\n8Nn0hTMozZ8LvL9m3Y5fZPcAs8mORjYjSzh/At4PbJC+jC+pWb82QS0A/i493xzYKz3fM32R7wsM\nBP4pxbERsGH6ov90iveoFH9XCepB4KxuttMW6Uv0pBTv8en1ljXv72lgZ7KjsHuAczvbFmnayWm7\nn57K25gsicwGdiD7gXAdcHlXZdSUNYwsyV0KHAps3mH+pFTuu1JdXwbuT/Na0hd9++f66RRX3gR1\nPXARsCnwNuAh4NSa97gS+Of0+XyCLBm132j1JrIfDpunut+b43N9B/AssHVNPDuW/X/kR+cPN/FZ\nX7UAL0XEqk7mLUjz282IiGsiYiXwXWAwWXNgXpdExNMRsRi4BXg6Iu5Idf+a7IupMyuBcZKGRcQr\nEfFImn4KcFFETIuI1RFxKVmz1n7pMQj4fkSsjIhryI56urJler9dOQz434i4PCJWRcSVZEcdR3R4\nf3+KiBXAfwF7dFMewPMRcX4qbwVwIvDdiJgTEcuALwDH9dT8FxFLgAPJEsdPgRcl/VbSiLTIacA5\nEfFE2tbfBPaQtB3Z0enjNZ/r94EXeogbgFT+ROCMiPhrRCwkOyI8rmaxeRHx08jORV4KjAJGSBpF\nlkxPS5/pyoi4N63T3ee6mixRjZM0KCLmRsTTeeK1xnOCsr56CWjp4ktwVJrf7tn2JxHxJvAcsHUv\n6vpLzfMVnbwe0sV6Hyb7Ipwn6V5J+6fp2wGfkfRq+wPYNsW0NTA/IvuZnczrJrZFZO+3K1t3sv48\nsvM97Wq/2Jd3837aPdvhdcc65pEd8YygByn5nBwR2wC7prK+n2ZvB/ygZhu9DCjFvjVrf67RSVxd\n2Y7sR8CCmrIvIjuSavfWNomI5enpELLP6eWIeKWLcjv9XCNiNnAG2ZHdQklXSerNPmgN5ARlffUA\n2a/TI2snShpC9gv3zprJ29bMHwBsQ9ZkA9mv97qIiIcjYhLZF99vyI5OIPsi/UZEDK95bJKObhYA\noyWppqgx3VRzB/Ch9L468zzZF2etMWTnQnp8Czmnd6xjDFlz21/ohYh4EvgFWaKCbDud2mE7bRwR\n95Ntp9rPVbWvgb8Cm9S8Hlnz/FmyfaelptxhEbFLjjCfBbbooqdod58rEXFFRBxItq2CrIOPVZAT\nlPVJam47Czhf0iGSBkkaS5YEniPrANHu3ZKOTEdbZ5B9OT2Y5v2F7NxJoSRtKOlESZulJqglZB0y\nIGvOOk3SvspsKukwSUPJEu8q4FPpPR0J7NNNVd8lO5dzaWr6QtJoSd+VtBtwM7CzpBMkbSDpWGAc\nWe+5nryYYu5p+1wJfFrS9ukHwjeBq7tofn2LpHdK+oykbdLrbcnOkbV/NhcCX5C0S5q/maSj07yb\ngF1qPtdPsXYSmgm8J13XthlZsyMAEbGArDPIdyQNkzRA0o6S3tvD+2xf9xbgx5I2T5/Re9LsLj9X\nSe+Q9PeSNgJeIzvyfrOLaqxkTlDWZxHxn8AXyXpbLQGmkf2KfV9E1HZVvgE4ljWdBY5MSQPgHODL\nqUnmswWHeBIwV9ISsvMpJ6a4p5OdgL8gxTSb7MQ8EfEG2VHhyWRNWseSdTroVES8TNZLbyUwTdJS\nsqPHxcDsiFhE1rPsM2TNgf8BHB4RL3VRZG3Zy4FvAFPT9unqvN3FZD8I7iPrgPIaWSeKniwl61Aw\nTdJfyRLTYylWIuJ6sqOMq9I2fIzs6JgU/9HAuel97UTWu7A99tvJOjI8CszgbxPyR8g6pMwi+wyu\nofum0lonkW3vJ8k6RZyR6uzycyU7/3QuWdPzC2RH1V/AKqm9N4xZXUk6k6xX3T+WHYvVl6R7yHru\n/azsWKx/8xGUmZlVUq4EJWljSe+odzBmZmbtemzik3QE2bmFDSNie0l7AF+LiA82IkAzM2tOeY6g\nziTrvfQqQETMBLYvonJlA4S2SZopaXon8yXph5JmS3pU0l5F1GtmZtWXZ4DJlRGxeO3LQQq9ZuXg\nbnoyHUrWK2gnsl5GP0l/u9XS0hJjx44tLEAzMyvOjBkzXoqIrXpaLk+CelzSCcBASTuRXedwf18D\nzGkScFm6Ov1BScMljUrXQHRp7NixTJ/+NwdkZmZWAZK6G5XlLXma+E4HdiG7qPJKsutczlj30NYS\nwB2SZkg6pZP5o1l72JTnWHtomLdIOkXSdEnTX3zxxYLCMzOzsvR4BJUuEvxSehTtwIiYL+ltZPcO\nejIi7luXgiJiCjAFoLW1tX9f3DX9Emi7puwoYPxR0Dq57CjK48/BrFRdJihJv6Obc01F9OKLiPnp\n70JJ15N1xqhNUPNZe1yvbcg3dln/1nYNvNAGI8eXF8MLbdnfZv5i9OdgVqrujqDOS3+PJBtb65fp\n9fH0cvDJzii7DfaAiFianv8D8LUOi/0W+KSkq8g6Ryzu6fzTemPkeJh8U3n1X3JYeXVXiT8Hs9J0\nmaDa760i6TsR0Voz63eddQlfByOA61PvwA2AKyLiVkmnpfovJBtgcyLZWFrLAf+MNDNrEnl68W0q\naYeImAMgaXuyu1/2SSpv906mX1jzPIB/7WtdZmbW/+RJUJ8G7pE0h+wmZdsBp9Y1KjMza3p5evHd\nmq5/emea9GSHWyiYmZkVrscEJekjHSbtLomIuKxOMZmZmeVq4tu75vlg4H3AI4ATlJmZ1U2eJr61\n7sgpaThwVd0iMrNq8QXLVpJ1uWHhXyloNHMz6wfaL1gu0wtt1UiS1lB5zkHVjigxABgH/LqeQZlZ\nxfiCZStBnnNQ59U8XwXMi4jn6hSPmZkZkK+Jb2JE3JseUyPiOUnfqntkZmbW1PIkqA90Mu3QogMx\nMzOr1d1o5p8A/gXYQdKjNbOGAlPrHZiZ2VpeaCv/XJR7EjZUd+egrgBuAc4BPl8zfWlEvFzXqMzM\nao0/quwIfOuTEnSXoCIi5kr6m8FaJW3R1yQlaVuyi31HkPUSnBIRP+iwzEHADcCf06TrIqLjLTnM\nbH3XOrn8xFD20VsT6ukI6nBgBlkCUc28AHboY92rgM9ExCOShgIzJN0eEbM6LPeHiDi8j3WZmVk/\n0939oA5Pf+tyUW668eCC9HyppCeA0UDHBNVYt3y+GhcllnkXVzOzCuiuk8Re3a0YEY8UFYSkscCe\nwLROZk9InTTmA5+NiMeLqreyRo6vRpu7mVmJumvi+0438wL4+yICkDQEuBY4IyKWdJj9CDAmIpZJ\nmgj8Btipi3JOAU4BGDNmzLoHdOi5676umZkVprsmvoPrXbmkQWTJ6VcRcV0nMSypeX6zpB9LaomI\nlzpZdgowBaC1tTU6zjczs/4lz1h8g8muhzqQ7MjpD8CFEfFaXyqWJODnwBMR8d0ulhkJ/CUiQtI+\nZBcWL+q6eUSsAAANKElEQVRLvWZ5/WXpa7y07HW+dtEDpcXwlUWLaRmyESNKi8CsPHnG4rsMWAqc\nn16fAFwOHN3Hug8ATgLaJM1M074IjAGIiAuBo4BPSFoFrACOiwgfHVlDvLTsdZa/sbrUGJa/sZqX\nlr3uBGVNKU+C2jUixtW8vltSn3vaRcR/s3bX9c6WuQC4oK91ma2rTTYcyNWn7l9a/Y9/c2BpdVfJ\nFdOe4YaZ80uNwUezjZcnQT0iab+IeBBA0r7A9PqGZUbpN8obu3IOcwf19XI/K8INM+cza8ESxo0a\nVloMPpptvDwJ6t3A/ZKeSa/HAE9JaiMbbWK3ukVnza39RnklXRM2d9AOTN34YHYppXbraNyoYT6a\nbTJ5EtQhdY/CrCsl3iivvXPEKaXUbmY9JqiImCdpc2Db2uWLvFDXzMysozzdzM8GTgaeZs2t3wu7\nUNfMzKwzeZr4jgF2jIg36h2MmZlZuzwJ6jFgOLCwzrGYraXsC2XL7jVm1uzyJKhzgP+R9BjwevvE\niPhg3aIyo/wLZceNGsakPUaXVr9Zs8uToC4FvgW0AW/WNxyrirKPXgA++8bq0i+UrYKxK+eUe7M8\n3/7FSpInQS2PiB/WPRKrlLKPXiAbxaFlyEalxlC2qRtnYzaXei2Wb/9iJcmToP4g6Rzgt6zdxOdu\n5uu50o9eLtmsvLor4s5NJnLnJhO5enJzH0VWxfI3VnNsia0KAJP2GM0J+/bhlkL9SJ4EtWf6u1/N\nNHczN7Om0jJkI15a9nrPC9bRrAXZHYicoJJG3BfKzKzqRgwdzIihg0s9mi376K3R8hxBIekwsmbw\nwe3TIuJrfa1c0iHAD4CBwM8i4twO85XmTwSWAye7adHMrDkM6GkBSRcCxwKnk90e42hgu75WLGkg\n8CPgUGAccLykcR0WO5TsFu87kQ2J9pO+1mtmZv1DjwkKmBARHwFeiYizgP2BnQuoex9gdkTMSaNU\nXAVM6rDMJOCyyDwIDJc0qoC6zcys4vIkqBXp73JJWwMrgSKSxGjg2ZrXz6VpvV0GAEmnSJouafqL\nL75YQHhmZlamPAnqRknDgW8DjwBzgSvqGdS6iIgpEdEaEa1bbbVV2eGYmVkf5enFd3Z6eq2kG4HB\nEbG4gLrnk93Co902aVpvlzEzs/VQrl587SLidWou1u2jh4GdJG1PlnSOA07osMxvgU9KugrYF1gc\nEQsKqt/MrN+ZtWBJ6d3Nx209jK8eUf/xTXqVoIoUEaskfRL4PVk384sj4nFJp6X5FwI3k3Uxn03W\nzXxyWfE20hXTnuGGmeUeKLaPg2dm1dFsgxeXlqAAIuJmsiRUO+3CmucB/Guj4yrbDTPnl36rB4+D\nZ1Y9J+w7pmlGkYBuEpSkvbpb0RfM1te4UcM8Dp6ZNbXujqC+0808j8VnZmZ11WWC8hh8ZmZWph7P\nQUkaBHwCeE+adA9wUUSsrGNcZmbW5PJ0kvgJMAj4cXp9Upr28XoFZWZmlidB7R0Ru9e8vkvSH+sV\nkJmZGeQb6mi1pB3bX0jaASj3XuBmZrbey3ME9TngbklzyG63sR1NcsGsmZmVp9sEJWkA2WjmOwHv\nSJOfSkMerZfO+t3jzHp+SakxlH2RrplZFXSboCLiTUk/iog9gUcbFFPTGzdqWNMNaWJdq8LYa5P2\nGN1UIxhYNeRp4rtT0oeB69LQQ+u1RgyAaJZXFX6oTPvzy0z788uljg/pVoXmlCdBnQr8O7BK0mtk\n56EiIry3mNVZFcZeq8LgxW5VaE557gc1tOhKJX0bOAJ4A3gamBwRr3ay3FxgKVmvwVUR0Vp0LGbW\nvSokSWtOPXYzl3Rnnmm9dDuwa0TsBvwJ+EI3yx4cEXs4OZmZNZfuRjMfDGwCtEjanKxpD2AY0Kdj\n7Yi4reblg8BRfSnP6uSFNrjksHLrHzm+vPrNrFTdNfGdCpwBbA3MYE2CWgJcUGAMHwWu7mJeAHdI\nWk02/t+UrgqRdApwCsCYMW6O6LPxFfjNMHJ8NeIws1Kop455kk6PiPN7XbB0BzCyk1lfiogb0jJf\nAlqBIzvrIShpdETMl/Q2smbB0yPivp7qbm1tjenTp/c2ZDOzrrW3Jky+qdw41gOSZuQ5bZOnk8T5\nkiYAY2uXj4jLeljv/T0EeDJwOPC+rrqvR8T89HehpOuBfYAeE5SZmfV/eW63cTmwIzCTNWPwBdBt\nguqhzEOA/wDeGxHLu1hmU2BARCxNz/8B+Nq61mlmZv1LnuugWoFxBV+kewGwEXC7JIAHI+I0SVsD\nP4uIicAI4Po0fwPgioi4tcAYzMx6p+yOQ5Cdl21tjuFQ8ySox8jOJS0oqtKIeHsX058HJqbnc4Dd\nO1vOzKzhqtBh54W27K8T1FtagFmSHgLeGiQ2Ij5Yt6jMzKqmdXL5iaHso7cGy5Ogzqx3EGZmZh3l\n6cV3r6QRwN5p0kMRsbC+YZmZWbPLM9TRMcBDwNHAMcA0SRVojDUzs/VZnia+LwF7tx81SdoKuAO4\npp6BmZlZJ6rQk3DkeDj03LpXkydBDejQpLeIHEdeZmZWsCr0JGygPAnqVkm/B65Mr48FbqlfSGZm\n1qkq9CRsoDydJD4n6UjgwDRpSkRcX9+wzMys2XV3u423AyMiYmpEXAdcl6YfKGnHiHi6UUGamVnz\n6XI0c0k3Al+IiLYO08cD34yIIxoQ3zqR9CIwrw9FtAAvFRROvTjGYjjGYvSHGKF/xNkMMW4XEVv1\ntFB3TXwjOiYngIhokzS2D4HVXZ433h1J06t+B1/HWAzHWIz+ECP0jzgd4xrd9cYb3s28jYsOxMzM\nrFZ3CWq6pH/uOFHSx8nusGtmZlY33TXxnUF2u4sTWZOQWoENgQ/VO7CSdXlr+QpxjMVwjMXoDzFC\n/4jTMSZ5bvl+MLBrevl4RNxV96jMzKzp9ZigzMzMyuAhi8zMrJKaKkFJOkTSU5JmS/p8J/Ml6Ydp\n/qOS9sq7boPjPDHF1ybpfkm718ybm6bPlDS9xBgPkrQ4xTFT0lfyrtvAGD9XE99jklZL2iLNq/t2\nlHSxpIWSHutifun7Y44Yq7Av9hRj6ftizjjL3h+3lXS3pFmSHpf0b50s09h9MiKa4gEMBJ4GdiDr\n6PFHYFyHZSaSjTMoYD9gWt51GxznBGDz9PzQ9jjT67lASwW25UHAjeuybqNi7LD8EcBdDd6O7wH2\nAh7rYn4V9seeYix1X8wZY6n7Yt44K7A/jgL2Ss+HAn8q+zuymY6g9gFmR8SciHgDuAqY1GGZScBl\nkXkQGC5pVM51GxZnRNwfEa+klw8C29QplnWOsU7r1jPG41kzIHJDRMR9wMvdLFL6/thTjBXYF/Ns\nx6408v+6t3GWsT8uiIhH0vOlwBPA6A6LNXSfbKYENRp4tub1c/ztxu9qmTzrFqW3dX2MtUeXD+AO\nSTMknVKH+CB/jBNSM8Atknbp5bqNihFJmwCHANfWTG7EduxJFfbH3ihjX8yrzH2xV6qwPyobLWhP\nYFqHWQ3dJ/PcbsMqStklAB9jzUjzAAdGxHxJbwNul/Rk+uXWaI8AYyJimaSJwG+AnUqII48jgKkR\nUfvrtirbsV/wvlioUvdHSUPIkuMZEbGkHnXk1UxHUPOBbWteb5Om5Vkmz7pFyVWXpN2AnwGTImJR\n+/SImJ/+LgSuJzv0bniMEbEkIpal5zcDgyS15Fm3UTHWOI4OzSkN2o49qcL+2KOS98UeVWBf7K3S\n9kdJg8iS068iu4tFR43dJ+t50q1KD7KjxTnA9qw5ibdLh2UOY+0TgA/lXbfBcY4BZgMTOkzfFBha\n8/x+4JCSYhzJmuvs9gGeSdu1Idsybz3AZmTnBTZt9HZM5Y+l65P7pe+POWIsdV/MGWOp+2LeOMve\nH9M2uQz4fjfLNHSfbJomvohYJemTwO/JepxcHBGPSzotzb8QuJmsl8psYDkwubt1S4zzK8CWwI8l\nAayKbGThEWTDU0G2w1wREbeWFONRwCckrQJWAMdFtic3ZFvmjBGyYbtui4i/1qzekO0o6UqyHmYt\nkp4DvgoMqomv9P0xR4yl7os5Yyx1X+xFnFDi/ggcAJwEtEmamaZ9kexHSCn7pEeSMDOzSmqmc1Bm\nZtaPOEGZmVklOUGZmVklOUGZmVklOUGZmVklNU03c7OqkLQlcGd6ORJYDbyYXi+PiAmlBGZWMe5m\nblYiSWcCyyLivLJjMasaN/GZVYikZenvQZLulXSDpDmSzk33Xnoo3Rdox7TcVpKulfRwehxQ7jsw\nK44TlFl17Q6cBryL7Ar/nSNiH7Jx705Py/wA+F5E7A18OM0zWy/4HJRZdT0cEQsAJD0N3JamtwEH\np+fvB8alYXAAhkkaEmlwVLP+zAnKrLper3n+Zs3rN1nzvzsA2C8iXmtkYGaN4CY+s/7tNtY09yFp\njxJjMSuUE5RZ//YpoDXdLXYW2Tkrs/WCu5mbmVkl+QjKzMwqyQnKzMwqyQnKzMwqyQnKzMwqyQnK\nzMwqyQnKzMwqyQnKzMwq6f8D/LrX+ZmFOMwAAAAASUVORK5CYII=\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 Sequences\")\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", "fig1.tight_layout()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Versions" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
SoftwareVersion
QuTiP4.2.0
Numpy1.13.1
SciPy0.19.1
matplotlib2.0.2
Cython0.26
Number of CPUs4
BLAS InfoINTEL MKL
IPython6.1.0
Python3.6.2 |Continuum Analytics, Inc.| (default, Jul 20 2017, 13:51:32) \n", "[GCC 4.4.7 20120313 (Red Hat 4.4.7-1)]
OSposix [linux]
Thu Apr 05 11:50:14 2018 BST
" ], "text/plain": [ "" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from qutip.ipynbtools import version_table\n", "\n", "version_table()" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [] } ], "metadata": { "anaconda-cloud": {}, "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.2" } }, "nbformat": 4, "nbformat_minor": 1 }