{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Optimization of Dissipative Qubit Reset" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:26.374878Z", "start_time": "2019-02-12T04:47:25.200533Z" }, "attributes": { "classes": [], "id": "", "n": "1" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "matplotlib.pylab 1.16.4\n", "numpy 1.16.4\n", "krotov 0.3.0+dev\n", "scipy 1.2.1\n", "matplotlib 3.1.0\n", "qutip 4.3.1\n", "CPython 3.7.3\n", "IPython 7.5.0\n" ] } ], "source": [ "# NBVAL_IGNORE_OUTPUT\n", "%load_ext watermark\n", "import qutip\n", "import numpy as np\n", "import scipy\n", "import matplotlib\n", "import matplotlib.pylab as plt\n", "import krotov\n", "\n", "%watermark -v --iversions" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "$\\newcommand{tr}[0]{\\operatorname{tr}}\n", "\\newcommand{diag}[0]{\\operatorname{diag}}\n", "\\newcommand{abs}[0]{\\operatorname{abs}}\n", "\\newcommand{pop}[0]{\\operatorname{pop}}\n", "\\newcommand{aux}[0]{\\text{aux}}\n", "\\newcommand{int}[0]{\\text{int}}\n", "\\newcommand{opt}[0]{\\text{opt}}\n", "\\newcommand{tgt}[0]{\\text{tgt}}\n", "\\newcommand{init}[0]{\\text{init}}\n", "\\newcommand{lab}[0]{\\text{lab}}\n", "\\newcommand{rwa}[0]{\\text{rwa}}\n", "\\newcommand{bra}[1]{\\langle#1\\vert}\n", "\\newcommand{ket}[1]{\\vert#1\\rangle}\n", "\\newcommand{Bra}[1]{\\left\\langle#1\\right\\vert}\n", "\\newcommand{Ket}[1]{\\left\\vert#1\\right\\rangle}\n", "\\newcommand{Braket}[2]{\\left\\langle #1\\vphantom{#2} \\mid\n", "#2\\vphantom{#1}\\right\\rangle}\n", "\\newcommand{op}[1]{\\hat{#1}}\n", "\\newcommand{Op}[1]{\\hat{#1}}\n", "\\newcommand{dd}[0]{\\,\\text{d}}\n", "\\newcommand{Liouville}[0]{\\mathcal{L}}\n", "\\newcommand{DynMap}[0]{\\mathcal{E}}\n", "\\newcommand{identity}[0]{\\mathbf{1}}\n", "\\newcommand{Norm}[1]{\\lVert#1\\rVert}\n", "\\newcommand{Abs}[1]{\\left\\vert#1\\right\\vert}\n", "\\newcommand{avg}[1]{\\langle#1\\rangle}\n", "\\newcommand{Avg}[1]{\\left\\langle#1\\right\\rangle}\n", "\\newcommand{AbsSq}[1]{\\left\\vert#1\\right\\vert^2}\n", "\\newcommand{Re}[0]{\\operatorname{Re}}\n", "\\newcommand{Im}[0]{\\operatorname{Im}}$\n", "\n", "This example illustrates an optimization in an *open* quantum system,\n", "where the dynamics is governed by the Liouville-von Neumann equation. Hence,\n", "states are represented by density matrices $\\op{\\rho}(t)$ and the time-evolution\n", "operator is given by a general dynamical map $\\DynMap$." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Define parameters\n", "\n", "The system consists of a qubit with Hamiltonian\n", "$\\op{H}_{q}(t) = - \\frac{\\omega_{q}}{2} \\op{\\sigma}_{z} - \\frac{\\epsilon(t)}{2} \\op{\\sigma}_{z}$,\n", "where $\\omega_{q}$ is an energy level splitting that can be dynamically adjusted\n", "by the control $\\epsilon(t)$. This qubit couples strongly to another two-level\n", "system (TLS) with Hamiltonian $\\op{H}_{t} = - \\frac{\\omega_{t}}{2} \\op{\\sigma}_{z}$ with\n", "static energy level splitting $\\omega_{t}$. The coupling strength between both\n", "systems is given by $J$ with the interaction Hamiltonian given by $\\op{H}_{\\int}\n", "= J \\op{\\sigma}_{x} \\otimes \\op{\\sigma}_{x}$." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The Hamiltonian for the system of qubit and TLS is\n", "\n", "$$\n", " \\op{H}(t)\n", " = \\op{H}_{q}(t) \\otimes \\identity_{t}\n", " + \\identity_{q} \\otimes \\op{H}_{t} + \\op{H}_{\\int}.\n", "$$\n", "\n", "In addition, the TLS is embedded in a heat bath with inverse temperature\n", "$\\beta$. The TLS couples to the bath with rate $\\kappa$. In order to simulate\n", "the dissipation arising from this coupling, we consider the two Lindblad\n", "operators\n", "\n", "$$\n", "\\begin{split}\n", "\\op{L}_{1} &= \\sqrt{\\kappa (N_{th}+1)} \\identity_{q} \\otimes \\ket{0}\\bra{1} \\\\\n", "\\op{L}_{2} &= \\sqrt{\\kappa N_{th}} \\identity_{q} \\otimes \\ket{1}\\bra{0}\n", "\\end{split}\n", "$$\n", "\n", "with $N_{th} = 1/(e^{\\beta \\omega_{t}} - 1)$." ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:26.385786Z", "start_time": "2019-02-12T04:47:26.379865Z" }, "attributes": { "classes": [], "id": "", "n": "2" } }, "outputs": [], "source": [ "omega_q = 1.0 # qubit level splitting\n", "omega_T = 3.0 # TLS level splitting\n", "J = 0.1 # qubit-TLS coupling\n", "kappa = 0.04 # TLS decay rate\n", "beta = 1.0 # inverse bath temperature\n", "T = 25.0 # final time\n", "nt = 2500 # number of time steps" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Define the Liouvillian\n", "\n", "The dynamics of the qubit-TLS system state $\\op{\\rho}(t)$ is governed by the\n", "Liouville-von Neumann equation\n", "\n", "$$\n", "\\begin{split}\n", " \\frac{\\partial}{\\partial t} \\op{\\rho}(t)\n", " &= \\Liouville(t) \\op{\\rho}(t) \\\\\n", " &= - i \\left[\\op{H}(t), \\op{\\rho}(t)\\right]\n", " + \\sum_{k=1,2} \\left(\n", " \\op{L}_{k} \\op{\\rho}(t) \\op{L}_{k}^\\dagger\n", " - \\frac{1}{2}\n", " \\op{L}_{k}^\\dagger\n", " \\op{L}_{k} \\op{\\rho}(t)\n", " - \\frac{1}{2} \\op{\\rho}(t)\n", " \\op{L}_{k}^\\dagger\n", " \\op{L}_{k}\n", " \\right)\\,.\n", "\\end{split}\n", "$$" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "def liouvillian(omega_q, omega_T, J, kappa, beta):\n", " \"\"\"Liouvillian for the coupled system of qubit and TLS\"\"\"\n", "\n", " # drift qubit Hamiltonian\n", " H0_q = 0.5 * omega_q * np.diag([-1, 1])\n", " # drive qubit Hamiltonian\n", " H1_q = 0.5 * np.diag([-1, 1])\n", "\n", " # drift TLS Hamiltonian\n", " H0_T = 0.5 * omega_T * np.diag([-1, 1])\n", "\n", " # Lift Hamiltonians to joint system operators\n", " H0 = np.kron(H0_q, np.identity(2)) + np.kron(np.identity(2), H0_T)\n", " H1 = np.kron(H1_q, np.identity(2))\n", "\n", " # qubit-TLS interaction\n", " H_int = J * np.fliplr(np.diag([0, 1, 1, 0]))\n", "\n", " # convert Hamiltonians to QuTiP objects\n", " H0 = qutip.Qobj(H0 + H_int)\n", " H1 = qutip.Qobj(H1)\n", "\n", " # Define Lindblad operators\n", " N = 1.0 / (np.exp(beta * omega_T) - 1.0)\n", " # Cooling on TLS\n", " L1 = np.sqrt(kappa * (N + 1)) * np.kron(\n", " np.identity(2), np.array([[0, 1], [0, 0]])\n", " )\n", " # Heating on TLS\n", " L2 = np.sqrt(kappa * N) * np.kron(\n", " np.identity(2), np.array([[0, 0], [1, 0]])\n", " )\n", "\n", " # convert Lindblad operators to QuTiP objects\n", " L1 = qutip.Qobj(L1)\n", " L2 = qutip.Qobj(L2)\n", "\n", " # generate the Liouvillian\n", " L0 = qutip.liouvillian(H=H0, c_ops=[L1, L2])\n", " L1 = qutip.liouvillian(H=H1)\n", "\n", " # Shift the qubit and TLS into resonance by default\n", " eps0 = lambda t, args: omega_T - omega_q\n", "\n", " return [L0, [L1, eps0]]\n", "\n", "\n", "L = liouvillian(omega_q=omega_q, omega_T=omega_T, J=J, kappa=kappa, beta=beta)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Define the optimization target" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The initial state of qubit and TLS are assumed to be in thermal equilibrium with\n", "the heat bath (although only the TLS is directly interacting with the bath).\n", "Both states are given by\n", "\n", "$$\n", " \\op{\\rho}_{\\alpha}^{th} =\n", "\\frac{e^{x_{\\alpha}} \\ket{0}\\bra{0} + e^{-x_{\\alpha}} \\ket{1}\\bra{1}}{2\n", "\\cosh(x_{\\alpha})},\n", " \\qquad\n", " x_{\\alpha} = \\frac{\\omega_{\\alpha} \\beta}{2},\n", "$$\n", "\n", "with $\\alpha = q,t$. The initial state of the bipartite system\n", "of qubit and TLS is given by the thermal state\n", "$\\op{\\rho}_{th} = \\op{\\rho}_{q}^{th} \\otimes \\op{\\rho}_{t}^{th}$." ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:26.466604Z", "start_time": "2019-02-12T04:47:26.457181Z" }, "attributes": { "classes": [], "id": "", "n": "6" } }, "outputs": [], "source": [ "x_q = omega_q * beta / 2.0\n", "rho_q_th = np.diag([np.exp(x_q), np.exp(-x_q)]) / (2 * np.cosh(x_q))\n", "\n", "x_T = omega_T * beta / 2.0\n", "rho_T_th = np.diag([np.exp(x_T), np.exp(-x_T)]) / (2 * np.cosh(x_T))\n", "\n", "rho_th = qutip.Qobj(np.kron(rho_q_th, rho_T_th))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Since we are ultimately only interested in the state of the qubit, we define\n", "`trace_TLS`. It returns the reduced state of the qubit\n", "$\\op{\\rho}_{q} = \\tr_{t}\\{\\op{\\rho}\\}$ when passed\n", "the state $\\op{\\rho}$ of the bipartite system." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:26.482459Z", "start_time": "2019-02-12T04:47:26.472974Z" }, "attributes": { "classes": [], "id": "", "n": "7" } }, "outputs": [], "source": [ "def trace_TLS(rho):\n", " \"\"\"Partial trace over the TLS degrees of freedom\"\"\"\n", " rho_q = np.zeros(shape=(2, 2), dtype=np.complex_)\n", " rho_q[0, 0] = rho[0, 0] + rho[1, 1]\n", " rho_q[0, 1] = rho[0, 2] + rho[1, 3]\n", " rho_q[1, 0] = rho[2, 0] + rho[3, 1]\n", " rho_q[1, 1] = rho[2, 2] + rho[3, 3]\n", " return qutip.Qobj(rho_q)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The target state is (temporarily) the ground state of the bipartite system,\n", "i.e., $\\op{\\rho}_{\\tgt} = \\ket{00}\\bra{00}$. Note that in the end we will only\n", "optimize the reduced state of the qubit." ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:26.496869Z", "start_time": "2019-02-12T04:47:26.485940Z" }, "attributes": { "classes": [], "id": "", "n": "8" } }, "outputs": [], "source": [ "rho_q_trg = np.diag([1, 0])\n", "rho_T_trg = np.diag([1, 0])\n", "rho_trg = np.kron(rho_q_trg, rho_T_trg)\n", "rho_trg = qutip.Qobj(rho_trg)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next, the list of `objectives` is defined, which contains the initial and target\n", "state and the Liouvillian $\\Liouville(t)$ that determines the system dynamics." ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:26.506011Z", "start_time": "2019-02-12T04:47:26.501241Z" }, "attributes": { "classes": [], "id": "", "n": "9" } }, "outputs": [ { "data": { "text/plain": [ "[Objective[Ļā‚€[4,4] to Ļā‚[4,4] via [š“›ā‚€[[4,4],[4,4]], [š“›ā‚[[4,4],[4,4]], uā‚(t)]]]]" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "objectives = [krotov.Objective(initial_state=rho_th, target=rho_trg, H=L)]\n", "objectives" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In the following, we define the shape function $S(t)$, which we use in order to\n", "ensure a smooth switch on and off in the beginning and end. Note that at times\n", "$t$ where $S(t)$ vanishes, the updates of the field is suppressed." ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:26.514731Z", "start_time": "2019-02-12T04:47:26.508936Z" }, "attributes": { "classes": [], "id": "", "n": "10" } }, "outputs": [], "source": [ "def S(t):\n", " \"\"\"Shape function for the field update\"\"\"\n", " return krotov.shapes.flattop(\n", " t, t_start=0, t_stop=T, t_rise=0.05 * T, t_fall=0.05 * T, func='sinsq'\n", " )" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We re-use this function to also shape the guess control $\\epsilon_{0}(t)$ to be\n", "zero at $t=0$ and $t=T$. This is on top of the originally defined constant\n", "value shifting the qubit and TLS into resonance." ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:26.524627Z", "start_time": "2019-02-12T04:47:26.517577Z" }, "attributes": { "classes": [], "id": "", "n": "11" } }, "outputs": [], "source": [ "def shape_field(eps0):\n", " \"\"\"Applies the shape function S(t) to the guess field\"\"\"\n", " eps0_shaped = lambda t, args: eps0(t, args) * S(t)\n", " return eps0_shaped\n", "\n", "\n", "L[1][1] = shape_field(L[1][1])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "At last, before heading to the actual optimization below, we assign the shape\n", "function $S(t)$ to the OCT parameters of the control and choose `lambda_a`, a\n", "numerical parameter that controls the field update magnitude in each iteration." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "pulse_options = {L[1][1]: dict(lambda_a=0.01, update_shape=S)}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Simulate the dynamics of the guess field\n" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [], "source": [ "tlist = np.linspace(0, T, nt)" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:26.545108Z", "start_time": "2019-02-12T04:47:26.537953Z" }, "attributes": { "classes": [], "id": "", "n": "13" } }, "outputs": [], "source": [ "def plot_pulse(pulse, tlist):\n", " fig, ax = plt.subplots()\n", " if callable(pulse):\n", " pulse = np.array([pulse(t, args=None) for t in tlist])\n", " ax.plot(tlist, pulse)\n", " ax.set_xlabel('time')\n", " ax.set_ylabel('pulse amplitude')\n", " plt.show(fig)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The following plot shows the guess field $\\epsilon_{0}(t)$ as a constant that\n", "puts qubit and TLS into resonance, but with a smooth switch-on and switch-off." ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:26.773279Z", "start_time": "2019-02-12T04:47:26.547926Z" }, "attributes": { "classes": [], "id": "", "n": "14" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEGCAYAAAB/+QKOAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO3de5Skd13n8fenu/p+mfuEYa4BohKRXGgDq1kNLokJBwmuwibegguOBxJFXVwD6yFsXM5hxUUPisKgc5J4MNGA0VGDIUcCURGYHhJyg5BhcpshyfRMz0xfZvpS3d/9o57qKTrVMzU9/Vy6+/M6p07X83vqqedbU8nzrd/l+f0UEZiZmc3WlHcAZmZWTE4QZmZWlxOEmZnV5QRhZmZ1OUGYmVldpbwDWEhr166Nbdu25R2GmdmisWfPnkMRsa7eviWVILZt20Z/f3/eYZiZLRqSnp5rn5uYzMysLicIMzOrywnCzMzqcoIwM7O6nCDMzKyu1BKEpM2S7pP0mKRHJb2nzmsk6WOS9kp6SNLFNfuuk/RE8rgurTjNzKy+NIe5loH/ERFfl9QD7JF0b0Q8VvOaq4DzksdrgT8DXitpNXAT0AdEcuyuiDiSYrxmZlYjtQQREc8BzyXPhyV9E9gI1CaIq4HbojLn+FckrZS0AbgMuDciBgEk3QtcCdyeVrzzdd/jB3ngmaNQlGnTpbwjAKAYUVQU5J8EFeRfpTj/HsVQlH+Pi7es4kdesTbvML5HJjfKSdoGXAR8ddaujcCzNdv7k7K5yuu993ZgO8CWLVsWJN5GfWdghP9+y24iivEfWVFylJmduVKTuO+9l7F5dWfeocxIPUFI6gY+C/xGRAwt9PtHxA5gB0BfX1+ml8jP7NlPs8SX3/8TrO9pz/LU1qCiLIhVkDAoSBjF+V7yDiDx3NExfuwj97HrG9/l+te/Iu9wZqSaICS1UEkOn46Iv63zkgPA5prtTUnZASrNTLXlX0wnyvn72pODXLh5pZNDgakIVTuKUcMsFv+D1NqyppPvO6ebrz45yPWvzzuak9IcxSTgL4BvRsRH53jZLuCXktFMrwOOJX0X9wBXSFolaRVwRVJWGBPlaR4+cIyLt67KOxQzWwL6tq3mgWeOFKZ2BenWIH4U+EXgYUkPJmXvB7YARMQngLuBNwJ7gePALyf7BiX9HrA7Oe7maod1UTwzOMpEeZrzN/TmHYqZLQHft76b4bEyA8PjrO8tRqtEmqOY/o3T1COT0UvXz7FvJ7AzhdAWxL6BUQDOXduVcyRmthS8Yn0PAHsPjhQmQfhO6nl68lAlQWxzgjCzBfCK9d0A7B0YyTmSk5wg5unJQ6Os7W5lRUdL3qGY2RJwTm8bHS3NPH34eN6hzHCCmKf9R06waVVxxiub2eImiQ0r2nnu2Im8Q5nhBDFPzw+N8ZKCtBOa2dKwYWU73z06lncYM5wg5umFY2O8ZIUThJktnA0rOlyDWOxGx8sMj5c5xzUIM1tAG1a0c3B4nMmp6bxDAZwg5uWFoUoV8JzetpwjMbOlZMOKDiLg4PB43qEAThDz8nySINwHYWYLacPKyjXluaPFaGZygpiHmRqE+yDMbAGt6660Shwamcg5kgoniHl4/lil+uc+CDNbSGu6WwE4POompkVrYHicztZmutsyWU7DzJaJ1V2VBDHoGsTideT4xMwXaWa2UNpKzfS0lzg86gSxaA2OOkGYWTrWdrdxaMRNTIvWkeMTrOp0gjCzhbe6q5XDbmJavAZHJ1jjGoSZpWBNVyuDbmJavI6MTrDKCcLMUrCmu23pj2KStFPSQUmPzLH/tyU9mDwekTQlaXWy7ylJDyf7+tOKcT7GJqcYnZhyH4SZpWJtd6UGMT2d/9KjadYgbgGunGtnRHwkIi6MiAuB9wFfmrWs6OuT/X0pxnjGjh6fBHAfhJmlYnVXK9MBR09M5h1KegkiIu4HGl1H+lrg9rRiWUjVqt/qLi8UZGYLr/rj8+jx/Pshcu+DkNRJpabx2ZriAD4vaY+k7ac5frukfkn9AwMDaYYKwJFR1yDMLD3VVSqXdA3iDPwU8O+zmpcujYiLgauA6yX92FwHR8SOiOiLiL5169alHSuDSVZ3H4SZpWFFZyVBHHOCAOAaZjUvRcSB5O9B4C7gkhziqutIMvzMo5jMLA3VGsSx48s8QUhaAfw48Pc1ZV2SeqrPgSuAuiOh8lAdn7yyw30QZrbwqteWItQgUpttTtLtwGXAWkn7gZuAFoCI+ETysp8GPh8RozWHngPcJaka319FxD+nFeeZGhqbpKe9RKm5CJUvM1tqZvogClCDSC1BRMS1DbzmFirDYWvL9gEXpBPV2Rs6Uaa33bUHM0tHqbmJ7rZSIWoQ/hl8hobGJul185KZpWhFRwtHT3iY66IzdGKS3navA2Fm6VnR0cKQaxCLz9BY2TUIM0vVys6WQvRBOEGcoaETlU5qM7O0rOhocR/EYjQ0NulOajNL1crOFt9JvdhMTwcj425iMrN09SY1iIh8Z3R1gjgDw+NlInAntZmlamVHKxPlacYmp3ONwwniDFRHFbgGYWZp6u2o/AgdGsu3mckJ4gxUvyz3QZhZmrrbKglieKycaxxOEGdg6ETly6pmdzOzNFR/hA67BrF4uAZhZlmoDqV3DWIRqfZBrHAfhJmlqDtJECPjThCLRvXGFdcgzCxNPW5iWnyGkupet4e5mlmK3Em9CA2dmKSnrURzk/IOxcyWsGqCGHKCWDw81beZZaG5SXS3lRhZqglC0k5JByXVXS5U0mWSjkl6MHl8oGbflZIel7RX0o1pxXimhk6UPVGfmWWiu620pPsgbgGuPM1r/jUiLkweNwNIagY+DlwFnA9cK+n8FONs2LAn6jOzjPS0l5ZuH0RE3A8MzuPQS4C9EbEvIiaAO4CrFzS4eRqdKNPV1px3GGa2DPS0l5b9MNf/JOkbkj4n6QeTso3AszWv2Z+U1SVpu6R+Sf0DAwNpxsrIWJlu1yDMLAPd7S1LuonpdL4ObI2IC4A/Bv5uPm8SETsioi8i+tatW7egAc42Mj5Ft2sQZpaBJd3EdDoRMRQRI8nzu4EWSWuBA8DmmpduSspyNzpenhl+ZmaWpt72EsPLtYlJ0kskKXl+SRLLYWA3cJ6kcyW1AtcAu/KKs6o8Nc2JySm6nCDMLANFGMWU2tVO0u3AZcBaSfuBm4AWgIj4BPCzwLsklYETwDVRWT6pLOkG4B6gGdgZEY+mFWejRiemAFyDMLNM9LS3MDY5zeTUNC3N+fyWT+1qFxHXnmb/nwB/Mse+u4G704hrvkaTqp4ThJlloXrP1chYmVVdrbnEkPcopkWjOtzMTUxmloUizMfkBNGgEdcgzCxD1Rld81x21AmiQTNNTJ5qw8wy0FuANSGcIBpUnTSrq9UJwszS112AVeWcIBrkJiYzy1K1v/P4RIEThKTvk/Qv1VlZJb1a0u+mH1qxuInJzLJUba0oehPTp4D3AZMAEfEQlZvXlpWTo5g81YaZpa96rRkteILojIivzSrL9/7vHIyMT9HSLNpKThBmlr6TNYip3GJoJEEckvRyIAAk/SzwXKpRFZDnYTKzLDU1ic7W5lxrEI1c8a4HdgA/IOkA8CTwC6lGVUAj42XfJGdmmepqKxU7QUTEPuANkrqApogYTj+s4hlxDcLMMtbdlu+iQXNe8ST91hzlAETER1OKqZBGxpwgzCxbXW3FbWLqSf5+P/DDnJxy+6eA2Z3WS97oRJlVnflMmGVmy1NXa4nRHDup50wQEfG/ASTdD1xcbVqS9EHgnzKJrkBGxstsXt2Zdxhmtox0t5V47thYbudvZBTTOcBEzfZEUrasjIyV6fY0G2aWoa62EqM53kndyBXvNuBrku5Ktt8C3Hq6gyTtBN4EHIyIV9XZ//PA7wAChoF3RcQ3kn1PJWVTQDki+hqIM1WjHsVkZhlbDKOYPiTpc8B/Top+OSIeaOC9b6GyINBtc+x/EvjxiDgi6SoqQ2lfW7P/9RFxqIHzpG56OhidmPI0G2aWqe625mKOYqqStAU4BNxVWxYRz5zquIi4X9K2U+z/cs3mV4BNp4slL9UqXren2TCzDHW1lRibnKY8NU0ph2VHG/lJ/E8kd1EDHcC5wOPADy5gHO8APlezHcDnJQXwyYjYMdeBkrYD2wG2bNmygCGdVB1F4CYmM8tSdWj96MQUKzoKmCAi4odqtyVdDLx7oQKQ9HoqCeLSmuJLI+KApPXAvZK+FRH3zxHfDirNU/T19UW915wtT/VtZnnobD055feKjpbMz3/GKSkivs739hXMm6RXA38OXB0Rh2vOcSD5e5BK09YlC3G++XKCMLM85D2jayN9ELV3VDcBFwPfPdsTJ30bfwv8YkR8u6Z8ZkqP5PkVwM1ne76zMToz1bcThJllp/qjNK8ZXRu54vXUPC9T6ZP47OkOknQ7cBmwVtJ+4CagBSAiPgF8AFgD/GkyfUd1OOs5wF1JWQn4q4j45wY/TypcgzCzPFR/lBa2BgE8FhF31hZIeitw5xyvByAirj3N/ncC76xTvg+4oIG4MlNdj9oJwsyydLIGkU+CaKQP4n0Nli1Z1WGubmIysywVtgaR3Lz2RmCjpI/V7Oplma0oV83ePb5RzswyVORO6u8C/cCbgT015cPAb6YZVNGMjJVpbhJtpezHIZvZ8lXYTupkXqRvSPp0RCyrGsNso+NlulqbZ9bCMDPLQkdLM00qYA1C0t9ExNuAB5I7mr9HRLw61cgKZGR8ip727G9SMbPlTRJdrfmtKneqJqb3JH/flEUgRTYyPjnTFmhmlqU8Z3Q9VRPTc8nfp7MLp5hGx6c8xNXMctHV1pzbmhCnamIa5uQkfVBZtyGqfyOiN+XYCmN4vEyvRzCZWQ6620qF7KTumWvfcjM6XualK9rzDsPMlqFCNjHVSmZwvZRKDeLfGlwwaMk47tXkzCwnXW0lBkeP53Lu0w7sl/QBKkuMrgHWArdI+t20AyuSkWSYq5lZ1rpaC9gHUePngQsiYgxA0oeBB4H/k2ZgRRFRWW7UNQgzy0OliSmfPohGbg3+LlDbAN8GHEgnnOIZL08zNR1OEGaWi0ondXFrEMeARyXdS6UP4nLga9X5mSLi11OML3ejnurbzHLU1VZiojzN5NQ0LRmvS93IVe+u5FH1xXRCKaZq1a7TfRBmloPaGV1XdrZmeu5G1qS+NYtAisqLBZlZnrqTWRxGckgQjYxiepOkByQNShqSNCxpqJE3l7RT0kFJj8yxX5I+JmmvpIeS4bTVfddJeiJ5XNf4R1pYXgvCzPJ0sgaRfUd1Iw1afwRcB6yJiN6I6DmDu6hvAa48xf6rgPOSx3bgzwAkraayROlrgUuAmyStavCcC8rrUZtZnrpyXFWukQTxLPBIRLxoRtfTiYj7gcFTvORq4Lao+AqwUtIG4CeBeyNiMCKOAPdy6kSTmmrWdhOTmeWhO8dV5Rq56v1P4G5JXwLGq4UR8dEFOP9GKgmoan9SNlf5i0jaTqX2wZYtWxYgpO9V/VLcSW1meehqzS9BNFKD+BBwnMq9ED01j0KIiB0R0RcRfevWrVvw93cntZnlqTvHJqZGrnovjYhXpXT+A8Dmmu1NSdkB4LJZ5V9MKYZTOu5OajPLUZ7rUjdSg7hb0hUpnX8X8EvJaKbXAceSdSjuAa6QtCrpnL4iKcvcyPgUrc1NtHo9ajPLwcwoponsRzE18rP4XcB7JY0Dk5zBehCSbqdSE1graT+VkUktVN7gE8DdwBuBvVSasX452Tco6feA3clb3RwRp+rsTs3oeJlOryZnZjlpKzXR0qxiNjGdzboQEXHtafYHcP0c+3YCO+d77oUyOl6e6SQyM8uapNzWhGh0PYhVVO5VmJm0LxnCuuSNTpTdQW1muepqLTEyVsAEIemdwHuodBQ/CLwO+A/gJ9INrRhGx6dmOonMzPKQ14yujfS8vgf4YeDpiHg9cBFwNNWoCmTEq8mZWc662vJZNKiRBDFWs1hQW0R8C/j+dMMqDvdBmFneutpKjOQwF1MjV779klYCfwfcK+kI8HS6YRXHca8mZ2Y5624r8dyxsczP28gopp9Onn5Q0n3ACuCfU42qQEbGyzPT7ZqZ5aHQo5iqIuJLaQVSRBFRaWJyDcLMclTkTupla7w8TdnrUZtZzrqTGsQ8JtU+K04Qp3A8ubW9yzO5mlmOutpKTAeMTU5net6GEoSkrZLekDzvkFSY2VzT5MWCzKwIapcdzVIjS47+CvAZ4JNJ0SYqI5qWPE/1bWZF0JXTokGN1CCuB34UGAKIiCeA9WkGVRSuQZhZEeS17GgjCWI8IiaqG5JKQLY9JTkZmUkQ7oMws/zktWhQIwniS5LeD3RIuhy4E/iHdMMqhplOatcgzCxHRW5iuhEYAB4GfpXKGg6/m2ZQRTFTg/BUG2aWo7w6qRu5k3oa+BTwKUmrgU2R9WDcnIy6k9rMCuBkDSLb+ZgaGcX0RUm9SXLYQyVR/GEjby7pSkmPS9or6cY6+/9Q0oPJ49uSjtbsm6rZt+tMPtRCcSe1mRVBd05NTI1c+VZExFCyLsRtEXGTpIdOd5CkZuDjwOXAfmC3pF0R8Vj1NRHxmzWv/zUqU4lXnYiICxv9IGkYnZiipVlej9rMclVt5i5iJ3VJ0gbgbcA/nsF7XwLsjYh9ySioO4CrT/H6a4Hbz+D9U+d5mMysCJqaRGdrcyE7qW8G7qFysd8t6WXAEw0ctxF4tmZ7f1L2IpK2AucCX6gpbpfUL+krkt4y10kkbU9e1z8wMNBAWI0b8VoQZlYQXW2lzBcNaqST+k4qQ1ur2/uAn1ngOK4BPhMRtT0wWyPiQJKQviDp4Yj4Tp34dgA7APr6+ha083x03OtRm1kxdOewaNCcVz9Jf8wpboiLiF8/zXsfADbXbG9Kyuq5hsod27XvfyD5u0/SF6n0T7woQaTp+MQUnb5JzswKoKst+yamU/087j/L994NnCfpXCqJ4Rrg52a/SNIPAKuA/6gpWwUcj4hxSWupTPXx+2cZzxkbcQ3CzAqiq7XEyFhBEkRE3Ho2bxwRZUk3UOm/aAZ2RsSjkm4G+iOiOnT1GuCOWfdWvBL4pKRpKv0kH64d/ZSV0fEy5/S0Z31aM7MXyWPZ0dP+PE6WGX1RU1NE/MTpjo2Iu6nceV1b9oFZ2x+sc9yXgR863funbXTc61GbWTEUspMaeG/N83YqHdTZr32Xg9EJr0dtZsWQx7rUjYxi2jOr6N8lfS2leApldLxMp2sQZlYAPe3Zr0vdSBPT6prNJuA1wIrUIiqI8fIUk1PhTmozK4Su1hJjk9OUp6YpNWczu0MjV789VPogRKVp6UngHWkGVQTVSbG8HrWZFUF1XZrRiSlWdBQkQUTEuVkEUjSeqM/MiqR2wr4VHS2ZnLORJqZ24N3ApVRqEv8KfCIish1vlbHhZLxxT3s2X4SZ2anksWhQIz+PbwOGgT9Otn8O+EvgrWkFVQTVzqCedtcgzCx/eSw72sjV71URcX7N9n2SMr9pLWvDY5OAFwsys2LIY9GgRno6vi7pddUNSa/l7KfhKDzXIMysSLpmlh2dzOycjVz9XgN8WdIzyfYW4HFJDwMREa9OLbocDSV9EN1OEGZWACebmLKrQTRy9bsy9SgKqDopVk+bO6nNLH95LDvayDDXp7MIpGhGxicpNYn2Fi83amb568qhk9pXvzkMj5Xpbi8hKe9QzMxoKzVRalKmNQgniDmMjHktCDMrDkmZT9jnBDGH4fGyb5Izs0LJetnRVBOEpCslPS5pr6Qb6+x/u6QBSQ8mj3fW7LtO0hPJ47o046xneGySHtcgzKxAsl52NLUroKRm4OPA5cB+YLekXXVWhvvriLhh1rGrgZuAPirTe+xJjj2SVryzjYyXWe/V5MysQLJeNCjNGsQlwN6I2BcRE8AdwNUNHvuTwL0RMZgkhXvJeLjtyFjZN8mZWaF0t5Vm5onLQpoJYiPwbM32/qRstp+R9JCkz0jafIbHImm7pH5J/QMDAwsRN5CMYnITk5kVSG97y7Ia5voPwLbkbux7gVvP9A0iYkdE9EVE37p16xYsMHdSm1nR9LSXGDqR3VQbaSaIA8Dmmu1NSdmMiDgcEePJ5p9TmdajoWPTNF6eYqI87SYmMyuU3o4WhsaWRoLYDZwn6VxJrcA1wK7aF0jaULP5ZuCbyfN7gCskrZK0CrgiKctEdZoNNzGZWZH0tleWHZ0oT2dyvtSugBFRlnQDlQt7M7AzIh6VdDPQHxG7gF+X9GYqS5kOAm9Pjh2U9HtUkgzAzRExmFass3kmVzMromqz9/DYJGu621I/X6pXwIi4G7h7VtkHap6/D3jfHMfuBHamGd9chl2DMLMC6u2oXJOGxsqZJIi8O6kLadhTfZtZAfXW1CCy4ARRR7WJqdejmMysQKpNTEMnshnq6gRRh5cbNbMiOtnE5BpEbqo1CDcxmVmR9M7UIJwgclP9x/coJjMrkuo1KavpNpwg6jh2YpL2libaSs15h2JmNqOrtUST3MSUq2MnJlnR4Q5qMyuWpibR097iJqY8OUGYWVH1dmQ3o6sTRB1OEGZWVD1t2c3H5ARRx7ETZScIMyuk3o6S74PI09CJSXqdIMysgHrbXYPIlZuYzKyoetpb3AeRl/LUNCPjbmIys2KqNDG5BpGLoSQzO0GYWRH1trcwPF5majpSP5cTxCzHkszsBGFmRbSys3JtOpZBLcIJYhYnCDMrstVdrQAcOT6R+rlSTRCSrpT0uKS9km6ss/+3JD0m6SFJ/yJpa82+KUkPJo9ds49NixOEmRXZys4kQYymnyBSm41OUjPwceByYD+wW9KuiHis5mUPAH0RcVzSu4DfB/5bsu9ERFyYVnxzcYIwsyJbXU0Qxxd3E9MlwN6I2BcRE8AdwNW1L4iI+yLieLL5FWBTivE0xAnCzIqs2geRRQ0izQSxEXi2Znt/UjaXdwCfq9lul9Qv6SuS3jLXQZK2J6/rHxgYOLuIOfmPXq3GmZkVSZZ9EIVY8EDSLwB9wI/XFG+NiAOSXgZ8QdLDEfGd2cdGxA5gB0BfX99Zj/saHJ2gp71Ea8n992ZWPJ2tzbQ2Ny36JqYDwOaa7U1J2feQ9AbgfwFvjojxanlEHEj+7gO+CFyUYqwzDo9OsKbLtQczKyZJrOxsWfRNTLuB8ySdK6kVuAb4ntFIki4CPkklORysKV8lqS15vhb4UaC2czs1g6PjM1U4M7MiWt3VuribmCKiLOkG4B6gGdgZEY9Kuhnoj4hdwEeAbuBOSQDPRMSbgVcCn5Q0TSWJfXjW6KfUHB6ZYNOqzixOZWY2Lys7WxZ3ggCIiLuBu2eVfaDm+RvmOO7LwA+lGdtcBkcnuGDTyjxObWbWkNVdrXz7hZHUz+Oe2BoRwZHjE6zudhOTmRXXys5Wji72O6kXm6GxMpNT4U5qMyu01Z2tHDk+mfqEfU4QNQaTUQHupDazIlvf28bUdMxcs9LiBFFjcLQyytYJwsyKbH1PGwAHh8dSPY8TRI2B4UqCWNvdlnMkZmZzW9/bDsDBofHTvPLsOEHUeP5YJRu/ZEV7zpGYmc2tWoN4Ycg1iMw8PzROa3PTzGyJZmZFtG6mick1iMw8f+wE63vbaGpS3qGYmc2prdTMqs4W1yCy9PzQGC/pdfOSmRXf+p521yCy9PyxMfc/mNmisL63jYOuQWQjIlyDMLNFY+PKDg4cPZHqOZwgEodGJhibnGbjqo68QzEzO62ta7o4NDLB8Fh660I4QSSePDQKwLlru3KOxMzs9Latqcw6/fTh46d55fw5QSSeShLEy9Z25xyJmdnpbUt+zD51eDS1czhBJPYdGqWlWbx0pfsgzKz4tiY1iOqP2zQ4QSS+9fwQL1vbTanZ/yRmVnydrSU2ruzgm88Pp3aOVK+Gkq6U9LikvZJurLO/TdJfJ/u/Kmlbzb73JeWPS/rJNOOMCB545igXbfFCQWa2eFy4ZSUPPnM0tfdPLUFIagY+DlwFnA9cK+n8WS97B3AkIl4B/CHwf5Njz6eyhvUPAlcCf5q8Xyoe/e4Qx05McvHWVWmdwsxswfVtXcWBoyfYezCd1eXSrEFcAuyNiH0RMQHcAVw96zVXA7cmzz8D/BdVFqe+GrgjIsYj4klgb/J+C25obJL33vkNWktNXP7Kc9I4hZlZKt706pdSahLv/vSeVBYPSjNBbASerdnen5TVfU1ElIFjwJoGjwVA0nZJ/ZL6BwYGzjjInrYSP/CSHj52zUWs8joQZraIrOtp4w/eegGv2bqK5hTmkCst+DtmLCJ2ADsA+vr6zjiFSuKPrrloweMyM8vCWy7ayFsuqvv7+aylWYM4AGyu2d6UlNV9jaQSsAI43OCxZmaWojQTxG7gPEnnSmql0um8a9ZrdgHXJc9/FvhCRERSfk0yyulc4DzgaynGamZms6TWxBQRZUk3APcAzcDOiHhU0s1Af0TsAv4C+EtJe4FBKkmE5HV/AzwGlIHrI2IqrVjNzOzFVPnBvjT09fVFf39/3mGYmS0akvZERF+9fb5t2MzM6nKCMDOzupwgzMysLicIMzOra0l1UksaAJ6e5+FrgUMLGM5i4M+89C23zwv+zGdqa0Ssq7djSSWIsyGpf66e/KXKn3npW26fF/yZF5KbmMzMrC4nCDMzq8sJ4qQdeQeQA3/mpW+5fV7wZ14w7oMwM7O6XIMwM7O6nCDMzKyuZZ8gJF0p6XFJeyXdmHc8WZD0lKSHJT0oaUnObihpp6SDkh6pKVst6V5JTyR/l9Qi5HN85g9KOpB81w9KemOeMS40SZsl3SfpMUmPSnpPUr5kv+tTfOYF/66XdR+EpGbg28DlVJY13Q1cGxGP5RpYyiQ9BfRFxJK9mUjSjwEjwG0R8aqk7PeBwYj4cPJjYFVE/E6ecS6kOT7zB4GRiPiDPGNLi6QNwIaI+LqkHmAP8Bbg7SzR7/oUn/ltLPB3vdxrEJcAeyNiX0RMAHcAV+ccky2AiLifyhojta4Gbk2e30rlf6olY47PvKRFxHMR8fXk+TDwTSrr1y/Z7/oUn3nBLfcEsRF4tmZ7Pyn9QxdMAJ+XtO/9sPkAAALESURBVEfS9ryDydA5EfFc8vx54Jw8g8nQDZIeSpqglkxTy2yStgEXAV9lmXzXsz4zLPB3vdwTxHJ1aURcDFwFXJ80TSwrydK2y6F99c+AlwMXAs8B/y/fcNIhqRv4LPAbETFUu2+pftd1PvOCf9fLPUEcADbXbG9Kypa0iDiQ/D0I3EWlqW05eCFpv6224x7MOZ7URcQLETEVEdPAp1iC37WkFioXyk9HxN8mxUv6u673mdP4rpd7gtgNnCfpXEmtVNbE3pVzTKmS1JV0bCGpC7gCeOTURy0Zu4DrkufXAX+fYyyZqF4kEz/NEvuuJYnK2vbfjIiP1uxast/1XJ85je96WY9iAkiGgv0R0AzsjIgP5RxSqiS9jEqtAaAE/NVS/MySbgcuozIN8gvATcDfAX8DbKEyLfzbImLJdOrO8Zkvo9LkEMBTwK/WtM0vepIuBf4VeBiYTorfT6VNfkl+16f4zNeywN/1sk8QZmZW33JvYjIzszk4QZiZWV1OEGZmVpcThJmZ1eUEYWZmdTlBmM2TpJWS3p08f6mkz+Qdk9lC8jBXs3lK5sH5x+rMqWZLTSnvAMwWsQ8DL5f0IPAE8MqIeJWkt1OZPbQLOA/4A6AV+EVgHHhjRAxKejnwcWAdcBz4lYj4VvYfw6w+NzGZzd+NwHci4kLgt2ftexXwX4EfBj4EHI+Ii4D/AH4pec0O4Nci4jXAe4E/zSRqswa5BmGWjvuSufqHJR0D/iEpfxh4dTIT548Ad1am1gGgLfswzebmBGGWjvGa59M129NU/r9rAo4mtQ+zQnITk9n8DQM98zkwmb//SUlvhcoMnZIuWMjgzM6WE4TZPEXEYeDfJT0CfGQeb/HzwDskfQN4FC93awXjYa5mZlaXaxBmZlaXE4SZmdXlBGFmZnU5QZiZWV1OEGZmVpcThJmZ1eUEYWZmdf1/s2WUzMFdJBAAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plot_pulse(L[1][1], tlist)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We solve the equation of motion for this guess field, storing the expectation\n", "values for the population in the bipartite levels:\n" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [], "source": [ "psi00 = qutip.Qobj(np.kron(np.array([1,0]), np.array([1,0])))\n", "psi01 = qutip.Qobj(np.kron(np.array([1,0]), np.array([0,1])))\n", "psi10 = qutip.Qobj(np.kron(np.array([0,1]), np.array([1,0])))\n", "psi11 = qutip.Qobj(np.kron(np.array([0,1]), np.array([0,1])))\n", "proj_00 = qutip.ket2dm(psi00)\n", "proj_01 = qutip.ket2dm(psi01)\n", "proj_10 = qutip.ket2dm(psi10)\n", "proj_11 = qutip.ket2dm(psi11)" ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:26.872504Z", "start_time": "2019-02-12T04:47:26.775372Z" }, "attributes": { "classes": [], "id": "", "n": "15" } }, "outputs": [], "source": [ "guess_dynamics = objectives[0].mesolve(\n", " tlist, e_ops=[proj_00, proj_01, proj_10, proj_11]\n", ")" ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:27.116092Z", "start_time": "2019-02-12T04:47:26.874340Z" }, "attributes": { "classes": [], "id": "", "n": "16" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEGCAYAAAB/+QKOAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO3dd3gU5fbA8e9JSAgJECB0AiR0kBJIQIogikhTuggIChZsWFD0qteC3mu5iooVsdClizQFsaAiRUgg9N5D6IGQXt/fH7P4CxAgQDaT3T2f59mH3Z2Z3TNskrPzlvOKMQallFLqQl52B6CUUqpw0gShlFIqV5oglFJK5UoThFJKqVxpglBKKZWrInYHkF/Kli1rQkJC7A5DKaVcSlRU1EljTLnctrlNgggJCSEyMtLuMJRSyqWIyIFLbdMmJqWUUrnSBKGUUipXmiCUUkrlShOEUkqpXGmCUEoplStNEEoppXKlCUIppVSunDoPQkQ6Ax8B3sDXxph3LtheHRgPlAPigEHGmBjHtixgk2PXg8aY7k4JMiMF/hwNPn7g42/dileAkpWhZBUIKAsiTnlrpVxNZlY2p5MzOJ2cTmJaJlnZhswsQ7Yx+Pl44e9bhOJFixDo70OJokUQ/d1xaU5LECLiDXwGdARigLUissAYszXHbqOBycaYSSJyK/A2MNixLcUYE+as+P6RlgB/fQAmO/ftfoFQvoF1q9gIqreGsnU0aSi3duxsKhsOnWHbkQT2n0pi38kkDpxK4nRyRp5fI8DXm0qlilEp0I+a5YpTt2IJ6lYsQZ0KJShe1G3m6Lo1cdaCQSLSChhljOnkePwigDHm7Rz7bAE6G2MOifVVI94YU9KxLdEYUzyv7xcREWGueSa1MZCVARnJkJ4EiUfhbCzEx8DJnXBsKxzfBmnx1v7+QVaiqNMZ6nSBgKBre1+lComY08ks33WS5btOEHXgNMfOpgHW96DKgcWoHuRP9aAAKpQsSlCAL6UDfCletAhFvLzw9hK8BFIzs0lKyyQxLZMzyekciU/lyJlUYuNT2HM8kaT0LAC8BOpVLEnzkNJEhJShZY0gypUoaufpezQRiTLGROS2zZlpvApwKMfjGODGC/bZAPTGaobqBZQQkSBjzCnAT0QigUzgHWPMPKdFKgJFfK1bsVIQWAWqhJ+/jzEQtxcOrLRu+/6AbQtBvKBaa2hyN9zQG4rmOacpZRtjDFuPnGXhhiMs3XKUvSeTAKgU6EfLGkE0CS5Fk6qBNKgUSDFf7+t+v+xsw+EzKWw/msDmw/FEHohjVmQMk1ZZVR6aBAdya70KdKhfnhsql9SmqULCmVcQfbGuDh50PB4M3GiMGZ5jn8rAp0Ao8CfQB2hojDkjIlWMMYdFpAbwG9DBGLPngvcYBgwDqFatWviBA5csKZL/jIEj0bD9B9gyD07tAt/i0LAPtHjIao5SqpA5Gp/KrMhDzIs+zN4TSXh7Ca1rBtG+bnna1S5LrfLFC+yPc0ZWNltjz7J81wl+3X6c6ENnMAaqlilGjyZV6Nm0MrXKlyiQWDzZ5a4gbG1iumD/4sB2Y0xwLtsmAouMMXMu9X7X1cR0vYyBQ2tg3WTYMtdqqqrZAW56GkLaan+FspUxhhW7TzF19QF+3naMrGxDyxpluLNJZbo0rESZAF+7QwTgZGIay7YfZ8GGWFbsPkm2gRsql6R/86r0ahas/RZOYleCKALsBDoAh4G1wEBjzJYc+5QF4owx2SLyJpBljHlVREoDycaYNMc+q4AeF3Rwn8fWBJFTymlY+w38/QUknYDgFtDxdavPQqkClJVt+HHTET5btpvtRxMoE+DLXRHB3NOiOtWC/O0O77KOJ6SycMMR5q6LYUvsWQJ8vendLJjBrapTp4JeVeQnWxKE4427AmOwhrmON8a8KSJvAJHGmAWOZqi3AYPVxPS4Iym0BsYB2VhzNcYYY7653HsVmgRxTkYKRH9rDaFNOGJ1aHd4DSo0sDsy5eYys7L5fv1hxv6+h70nk6hZLoBHbq7JnU0q4+dz/f0JBckYQ/ShM0xZfYBFG4+QnplN29pleax9LVrWKKN9FfnAtgRRkApdgjgnPdm6mvhrDKQnQPMH4daXreGzSuUjYwzLdhznncXb2XkskQaVSjL81lp0vqEiXl6u/4c0Limd6WsOMmHFfk4mptG0Wikea1+LDvXKu8X52UUTRGGQHAfL3oK1X0Px8tD5bWvUk34DUvlg8+F43vxhG6v2niK0bADPd6pL54YV3fIbdmpGFrOjYhj3xx5iTqfQoFJJnutcl/Z1yrnl+TqbJojC5PA6WDTCGgFV6zbo/ok1a1upa5CYlsn7S3cwaeV+Svn78lSH2gy8sRo+3u5fRSczK5sFG2IZ88suDsYl0yKkDM93rktESBm7Q3MpmiAKm+ws60ril1Hg7QPdPoBGfe2OSrkQYww/bTnKqAVbOZaQyqAbqzOyU10Ci/nYHVqBS8/MZubag3z8225OJKRxW/3yvNS1PjXK6ZykvNAEUVid2gPfPwwxa63mpm7vg79++1GXF5eUzktzN7Fky1HqVyrJW70a0rRaabvDsl1yeiYTVuxn7O97SMvM4v6bQnni1to6PPYKNEEUZlmZsPIjWPY2lKgE/SZePItbKYffth/j+TmbOJuSwTO31+HBm0Ip4gHNSVfjeEIq7y7ZwZyoGMqVKMq/Otejd9Mq2pF9CZogXEFMFMweYg2J7fSWNRtbO9yUQ3J6Jv9ZtI3paw5Sr2IJPrw7jPqVStodVqG2/uBpRi3cyoZDZ4ioXpq3ezeits6huIgmCFeRHAfzHoWdS+CGXtD9U63tpNh9PIFHpq5jz4lEhrWrwTMd61C0iGvNZ7BLdrZhzroY3vpxG0lpmTx+Sy0ebV9T//9y0AThSrKzrSanX9+wSowPmA6lqtkdlbLJ/OjDvDh3E/6+3nzUvyltapW1OySXdDIxjTcWbmXBhlhqlS/O//o0Iry69vfB5ROENl4WNl5ecNMIuGcOnDkEX94CB1fbHZUqYGmZWbwybzNPzYjmhsol+eHJtpocrkPZ4kX5eEBTJgxpTkp6Fn2/WMXrC7eQmpFld2iFmiaIwqpWB3joV2vG9aQ7Yf23dkekCsipxDQGff03U1Yf4OF2NZj2UEsqlPSzOyy3cEu98iwd0Y57W1Znwor9dPt4ORsOnbE7rEJLE0RhVrY2PPgLVGsF8x+zZmK7SZOgyt2Oown0+GwFG2Pi+XhAU17sWt8jJr0VpICiRXi9R0OmPnAjyelZ9B67kg9/3klG1iVWlfRg+pNX2PmXgUHfQdNB8Mf/YOFT1tBY5XaWbT9On7ErScvMZubDrejeRGfYO9NNtcuy5Ol29GhSmY9+3UXvz1ey+3iC3WEVKpogXIG3jzWiqd1zsG4SzBpsFQFUbmPq6gM8MGkt1YP8WTC8DWFVS9kdkkcILObDB3eHMfaeZsScTuaOT/5i1tpDuMvgneulCcJViFhVYLuOhh2LYUpPa1iscmnGGD78eScvz9vMLXXLM/uRVlQKLGZ3WB6nS6NK/PR0O5pVK83z323kqRnRJKRm2B2W7TRBuJoWD8FdEyF2PUy8AxJP2B2RukZZ2YaX523mo1930Tc8mHGDw/H31bIQdilf0o8pD9zIyNvr8MOmI3T7+C82xnh2B7YmCFd0Q08YOAvi9sLErnD2iN0RqauUmpHF8Gnr+Pbvgzxyc03e69tYS2YUAt5ewvBbazNzWEsys7LpM3YlXy/fS3a2ZzY56U+kq6p5i9V5fTYWJnSx5kwol5CSnsVDkyNZvPkoL3erzwtd6uk6BoVMREgZfnyqLbfULc9/f9jGsClRxKd4XpOTJghXFtIGBs+z+iImdLWuKFShlpyeyf0T1/LX7pO817cxD7atYXdI6hJK+fsybnA4r93ZgN93HKfHp3+x7chZu8MqUJogXF3V5nDfAms50wndNEkUYolpmQyZsJa/953iw35h3BVR1e6Q1BWICEPbhDJjWEuS07Po9fkKvl8fY3dYBUYThDuoHAb3LYLMFJjUHc4ctDsidYGE1AyGjF9D1IHTfNS/KT2bVrE7JHUVIkLKsOjJm2gcXIoRMzfw2vzNpGe6/8Q6TRDuomJDq7kp7axVmiP+sN0RKYfEtEzuHb+G6ENn+HRAU+7UCXAuqXwJP7598EYeahvKpFUH6P/lKo7Ep9gdllNpgnAnlcNg0PeQdAomd4eEY3ZH5PFSM7J4cNJaNsbE8+nAZnRpVMnukNR18PH24t/dGvDZwGZsP5rAnZ/8ReR+952PpAnC3QSHw6A51tDXyd0h6aTdEXms9MxsHp0axd/74vigXxM6N6xod0gqn3RrXIn5j7eheNEiDPhqNTPXumezriYId1StJQycAaf3WzOuU+PtjsjjZGUbRsyMZtmOE7zZsxE9wrTPwd3UrlCC+Y/fRMsaQfzru02MWrCFTDcr+KcJwl2FtoP+38Lx7TB9AGS4d1tpYZKdbXjhu438sOkI/+5an4E36oJP7irQ34cJQ5rz4E2hTFy5n/smrOF0UrrdYeUbTRDurNZt0OsLOLAS5jygVWALyFs/bmN2VAxPdajNQ+10noO7K+Ltxct3NGD0XU1Yu+80PT5bwc5j7lEVVhOEu2vUF7q8Czt+gEVP6XoSTvb18r18/dc+hrQO4enbatsdjipAfcODmfFwS1Iysuj12Qp+3ur6g0Q0QXiCG4fBzf+C9VPhl9fsjsZtLdwQy39/2EbXRhV59Y4GWj7DAzWrVpqFw2+iZvniDJsSyWfLdrt06XBNEJ6i/YvQ/EFY8RGs+NjuaNzO6r2neHbWBpqHlOaDfmF4eWly8FQVA/2Y5Vjw6b2fdjBy9kaXnVSntYU9hYjV1JR8Cn5+BUpUhMb97I7KLew8lsCwyZFUC/Lnq3sj8PPxtjskZTM/H2/G3B1GjbLF+fCXnRw6ncy4QeGUDvC1O7SrolcQnsTLG3qNg5C2MO8x2Lfc7ohc3tH4VIaMX4OfjzcThzanlL9r/QFQziMiPHVbbT4e0JToQ2fo9fkK9pxItDusq6IJwtMUKQp3T4EyNWDmPdYwWHVNktMzeXDyWuJTMpgwtDnBpf3tDkkVQt2bVGb6Qy1JSM2k12crWLnbdSavaoLwRMVKwz2zwbsofHuXluS4BtnZhmdnbWBL7Fk+GdiUGyoH2h2SKsTCq5dm3uNtqBjox73j1zBjjWvMvNYE4alKV4eBMyH5JEzrB2mudelrtzG/7GTx5qO81KU+t9arYHc4ygVULePPnEdb07pWWV6Yu4m3ftxGViFfqU4ThCer0gz6ToCjG+E7nUiXV/OjD/Pxb7vpFxHMg21D7Q5HuZCSfj6Mvy+Ce1tV58s/9/LI1CiS0wvv750mCE9XtzN0fQ92LoHFz+tEuiuIPnSG5+ZspEVIGf7bs5HOdVBXrYi3F2/0aMioOxvw67Zj3PXFKo7Gp9odVq6cmiBEpLOI7BCR3SLyQi7bq4vIryKyUUR+F5HgHNvuE5Fdjtt9zozT4zV/EFo/CZHfwKrP7I6m0DoSn8JDkyMpX6IoYwc1w7eIfr9S125Im1C+ua85+08m0eOzv9gUU/iKajrtJ1xEvIHPgC5AA2CAiDS4YLfRwGRjTGPgDeBtx7FlgNeAG4EWwGsiUtpZsSrgtteh/p3WHImdP9kdTaGTmpHFI1OiSEnP4pv7mhNUvKjdISk3cEu98nz3WGuKeHnRb9wqlmw+cvUvkngcYqPzPzicewXRAthtjNlrjEkHZgA9LtinAfCb4/6yHNs7AT8bY+KMMaeBn4HOToxVeXlZcyQqNLQK+x3bandEhcqoBVvYEBPP+/2aULdiCbvDUW6kXsWSzHu8DfUqleCRqev4/Pc8lufIzoZVn8Mn4TB3mFOah52ZIKoAh3I8jnE8l9MGoLfjfi+ghIgE5fFYRGSYiESKSOSJEyfyLXCP5RsAA2aArz9Mv1sXG3KYvuYgM9YeYvgtteh0gy76o/JfuRJFmf5QS+5sUpl3l1jlOdIysy59QNJJ+LYv/PQiVL0R+k+zqiXkM7sbUUcCN4vIeuBm4DBwmf+V8xljvjTGRBhjIsqVK+esGD1LYBXoP926bJ05CDLT7I7IVtGHzvDa/C20q1OOER3r2B2OcmN+Pt583D+Mp2+rzXfrYhj89Rricltb4tQe+OpW2P8X3PGhNaepbC2nxOTMBHEYqJrjcbDjuX8YY2KNMb2NMU2BfzueO5OXY5UTBYdDz8/h4CpY+LTHjmw6mZjGo1OjqBBYlI/7h+GtBfiUk4kIT99WxyrPEXOGnp+tYPfxHGtLHNkA4ztBWgIMXQwR9zvlyuEcZyaItUBtEQkVEV+gP7Ag5w4iUlZEzsXwIjDecf8n4HYRKe3onL7d8ZwqKA37wM0vwIZpVgVYD5OZlc3waeuIS0pn7D3hWmNJFajuTSozY1hLktMz6fX5SpbvOgGHo2DiHVYFhPt/sr7IOZnTEoQxJhMYjvWHfRswyxizRUTeEJHujt3aAztEZCdQAXjTcWwc8B+sJLMWeMPxnCpIN/8LbugFv4yC7T/YHU2B+t+S7azeG8fbvRvRsIqW0VAFr1k1qzxHlVLFeHviHNIm9rTK5DzwE5QrmOZOceXFLHKKiIgwkZGRdofhftKTYWJXOLHT+sGs2MjuiJxuyeajPDI1intbVeeNHg3tDkd5uKTYrWR+3YWkLC9mNfyS4b07UMQ7/77bi0iUMSYit212d1Krws7X3+q09isJ0wdC0im7I3KqQ3HJPDdnA02CA3m524XTdpQqYHH7CJjem5LFfJjX6HPGRKXz4ORIzqZmFMjba4JQV1ayEtz9LSQeg9n3QVbB/HAWtPRMq98B4NOBOlNa2Sw+BiZ3h8xU5N75PNa3C2/1asRfu07SZcxyVu5x/jB0/Q1QeRMcDnd+BPuXw0//tjsap3h78TY2xMTzXt8mVC2jazsoGyUcg8k9IOUMDP4eKtwAwMAbqzHz4Vb4FvFi4Fd/M2JmNPtPJjktDE0QKu/CBkDLx2HNOFg32e5o8tVPW44yYcV+hrYJoXNDnQynbJQcB1N6wtlYa45D5abnbQ6vXpofn2zLIzfXZPHmI9zy/u88OjXKKaHomtTq6nR8A45vgUXPQLl6ULWF3RFdt0NxyTw32+p3eLFLfbvDUZ4s9SxM6WVNhrtnFlRrmetuxXy9eaFLPe5vE8LUvw86ba6SXkGoq+NdxFpDIrCKNdP6bKzdEV2Xc/0OBu13UDbLSIUZA+HYZug3GWq0v+Ih5Uv68UzHOjxze12nhKS/Derq+ZexRjalJcKMe6wfbBf1zuLtjn6HxtrvoOyTlWkt2rV/OfQca63TUghoglDXpkID6D0OYtfBohEuWY7jt+3HGL9iH0Nah9C5YSW7w1GeyhhY9DRsXwSd/weN+9kd0T80QahrV//O/y/HsXqs3dFcleMJqTw3eyP1K5Xkxa717A5HebJfRsH6KdDueWj5iN3RnEc7qdX1uflfVpvp0pehfH2oeYvdEV1RdrZh5OyNJKZlMqN/GEWLeNsdkvJUKz6CFWMg4gG45SW7o7mIXkGo6+PlBb2+gLJ1YPYQiNtrd0RXNGHlfv7ceYKX72hA7Qq6+I+yyfqp8POrcENva134Qri+uSYIdf2KloAB06z7M+6xOq8Lqa2xZ/nf4u3cVr8Cg26sZnc4ylNtWwQLnoCat1orOXoVzqtYTRAqf5SpAXdNgBPb4fuHreUQC5nUjCyemrGeQH8f/tenEVIIv7EpD7BvOcy5Hyo3g35ToEjhLSWvCULln5q3Qsf/WKMx/nzP7mgu8uYP29h1PJEP+jUhqHhRu8NRnig2GqYPgDKh1izposXtjuiyNEGo/NXqcWh8N/z+VqFaQ+KXrceYsvoAD94UStvaujytssHJ3TC1DxQrBYPmWvOJCjlNECp/iVhF/So3hbnD4Pg2uyPi+NlUnv9uIw0qleS5zs6ZcarUZZ2NtUpoAAyeZ1UicAGaIFT+8ylmlQf38bdKB6Scti0UYwwj52wkOT2TjwfokFZlg+Q4KzmknIZBc6BsLbsjyjNNEMo5AqvA3VPgzCGrQy4r05Ywpv59kD93nuDfXetTq7wOaVUFLD0JpvWDuH0wYPpFlVkLO00QynmqtYRuo2HPb/DrqAJ/+30nk3jrh220q1OOQS2rF/j7Kw+XmQ4zB8PhKOg7HkLb2h3RVdOZ1Mq5wofA0U2w8hOo2LjA6sxkZmXz7KxofLyFd/s01iGtqmBlZ8O8R2DPr9D9U6h/h90RXRO9glDO1/kdqN7GmhgUu75A3nLcn3tZd/AM/+nZkIqBfgXynkoBVvG9Jf+Czd/BbaOg2WC7I7pmmiCU83n7wF2TIKCcNdM68bhT325LbDxjftlJt8aV6N6kslPfS6mL/Dka1nwJrYZDm6ftjua65DlBiIi3iFQWkWrnbs4MTLmZ4uWg/7fWiI6Zg632WSdIy8zimZkbKOXvy397NNSmJVWwIsfDsv9CkwHWpFEX//nLU4IQkSeAY8DPwA+O2yInxqXcUaUm0ONTOLQaFj/vlLf44Oed7DiWwLt9GlM6oPCWMFBuaMs8ayne2p2g+ydWIUsXl9dO6qeAusaYU84MRnmARn2tTusVY6BiI2j+QL699Nr9cXz5514GtKjGLfXK59vrKnVFe/+AuQ9Za7TfNdFqVnUDeU1xh4B4ZwaiPEiHV6FWR+sq4sDKfHnJxLRMnpkVTdXS/rzcrX6+vKZSeRK73poQGlQLBs4EX/dZujavCWIv8LuIvCgiz5y7OTMw5ca8vKHP11A6xOqPOHPoul/yzR+2EXM6hdF3NSGgqI7eVgXk1B6Y2heKlYFB30Gx0nZHlK/ymiAOYvU/+AIlctyUujbFSkH/6ZCVDjPvgfTka36pZduPM33NQYa1rUGL0MJfAE25ibNHYEpP6/7g76Gk+42Yy9NXLWPM6wAiUtzxuPCuCKNcR7k60PsrmN7fmiPR5+urHvURn5zBC3M3UrdCCUZ0rOOkQJW6QMppmNrbGpU3ZJFL1Ve6GnkdxdRQRNYDW4AtIhIlIjc4NzTlEep2hltfhs1zYOXHV334f37YysnEdEbf1QQ/Hy3EpwpAejJM6w+ndkP/aS5XX+lq5LWJ6UvgGWNMdWNMdeBZ4CvnhaU8SttnoUFP+GUU7Polz4ct236cOVExPHpzTRoFBzovPqXOycqAOUPh0N/W1W+Nm+2OyKnymiACjDHLzj0wxvwOBDglIuV5RKDn51C+gVX59dSeKx4Sn5LBi3M3UadCcZ7o4J6X96qQMQYWPAk7l0C39+GGnnZH5HR5HsUkIq+ISIjj9jLWyCal8odvgHW57uVtLcmYevayu7/5w1ZOJKYx+q4musaDKhg/vwobpsEt/87X+TuFWV4TxP1AOWCu41bO8ZxS+ad0deg3CeL2XHYNid93HGdWZAwPt6tB4+BSBRyk8kgrPrL6yFoMg3bP2R1NgclTgjDGnDbGPGmMaea4PWWMsW+ZMOW+QttB19Gw+2dY+u+LNp9NtZqWapcvzlO31bYhQOVxoiZaVw839IbO/3P5+kpX47LDXEVkjDHmaRFZCJgLtxtjujstMuW5IobCyV2w+jNrdmqLh/7Z9NYP2zh2NpWxj7XRpiXlfJu/g4VPWzP/e41zi/pKV+NK8yCmOP4dfS0vLiKdgY8Ab+BrY8w7F2yvBkwCSjn2ecEY86OIhADbgB2OXVcbYx65lhiUi7r9P1ZT0+J/QZlQqHUbf+48wYy1h3jk5pqEVdWmJeVkO5fC3GFQrRX0mwxFPK/442XToTEmynE3zBjzR84bEHa5Y0XEG/gM6AI0AAaISIMLdnsZmGWMaQr0Bz7PsW2PMSbMcdPk4GnOleMoXx9mDyUpZjMvfLeRmuUCeFqblpSz7V8BswZDhRtg4Ay3qq90NfJ6vXRfLs8NucIxLYDdxpi9xph0YAbQ44J9DFDScT8QiM1jPMoTFC0BA2ZAET/SJvcl/exx3tMJccrZYtfDtLuhVDUYNBf8PHeOzWUThIgMcPQ/hIrIghy3ZUDcFV67ClYV2HNiHM/lNAoYJCIxwI/AEzm2hYrIehH5Q0RyXe1bRIaJSKSIRJ44ceIK4SiXVKoq0Td9gX/aSb4P+pxmlT3zm5wqIMe3w5TeVtG9wfMgoKzdEdnqSn0QK4EjQFng/RzPJwAb8+H9BwATjTHvi0grYIqINHS8ZzVjzCkRCQfmicgNxpjzBscbY77EmuVNRETERZ3oyvUlpmXy+B/C7cWe5rXE96yaTb3GedRIElVATu+3iu95+8C98yDwwu+znueyCcIYcwA4ALS6htc+DFTN8TjY8VxODwCdHe+1SkT8gLLGmONAmuP5KBHZA9QBIq8hDuXC3v5xG7HxKdzxyOOwv6i1nGNQbbjZc8aiqwKQcBQm94CMFBj6IwTVtDuiQiGvxfpaishaEUkUkXQRyRKRy091hbVAbREJFRFfrE7oBRfscxDo4HiP+oAfcEJEyjk6uRGRGkBtdOa2x1mx+yTf/n2QB28KJbx6aWg3Ehr3t5LEhhl2h6fcRXIcTO4JiSesNR0qaB3Sc/K6ssqnWH/gZwMRwL1Y3+gvyRiTKSLDgZ+whrCON8ZsEZE3gEhjzAIcRf9EZARWh/UQY4wRkXbAGyKSAWQDjxhjrtTnodxIYlomz8/ZSI2yATx7e13rSRFrrd+EWJj/OBSvADVvsTdQ5dpSzljNSnF74Z5ZEBxhd0SFihhz5aZ7EYk0xkSIyEZjTGPHc+sdw1MLhYiICBMZqS1Q7uKVeZuZ+vcBZj/cioiQCxYBSo2H8V3gzEGrOaBSY3uCVK4t9ayVHI5stOqA1bnd7ohsISJRxphcM2Neh7kmO5qJokXkXcc3fs+aUqgKzMo9J5my+gD3twm9ODmANezwntngVxK+vStflixVHiYtEb7tC0c2WPW/PDQ5XEle/8gPxmomGg4kYXU+93FWUMpzJaVl8q/vNhIS5M/Ic01LuQmsAvfMsToVp/axVvhSKi/Sk61VDGMioc83UK+b3REVWnkt1nfAGJNijDlrjHndGPOMMWa3s4NTnufdJduJOZ3Cu32bUMz3ChPiKgJS+GAAAB1kSURBVDSA/t/C6X0wfSBkpBZMkMp1ZaTCjAGw/y9ruLQHrOlwPa5UrG8TuRTpO+dcf4RS+WH13lNMWnWAoW1CaBGaS9NSbkLbQs+x8N0D8P3D0HeCxxVUU3mUmQYzB8HeP6wFqhrfZXdEhd6VRjHdUSBRKI+XnG6NWqoe5M9znS7TtJSbRn3hbCz8/Ar8GGSt9qUT6VROWRkwe6hVRv7OjyBsoN0RuYS8TJRTyuneXbKDg3HJzBzWEn/fvI6+zqH1E5B0wlrUpVhp6PBK/gepXFNmmpUcdvxgrTUSPsTuiFxGnn4TRSSB/29q8gV8gCRjTMlLH6VU3qzdH8ekVfu5r1V1bqwRdG0vIgId37CGwC4fDcVKWUlDebaMVKsq666lVnLIsbaIurI8JQhjTIlz90VEsKqytnRWUMpzpKRn8fycjQSXLsbznetd34uJwB0fWkli6cvWcNhm9+ZPoMr1ZKTAjIGw5ze4Y4y1EJW6Klfdm2cs84BOTohHeZj3l+5g38kk/tenMQFFr6Fp6UJe3tD7K6jZARY+BVvmXf9rKteTngTT+sGeZdDjM00O1yivTUy9czz0wiq3oWMK1XWJOnCab1bs454bq9G6Zj6WVS7iC3dPsco2f/egta5ErQ759/qqcEtLtJLDwVXQ6wto0t/uiFxWXq8g7sxx64RV7vvCxX+UyrPUjCyem7OByoHFeLFr/fx/A98AGDgTytWDGfdYQxuV+0s9a02cPLjaupLU5HBd8toHoddnKl99+MtO9p5IYsoDLSieH01LuSlWyqrrP/EOa4Wwe2Zb8yaUe0o6aSWHY5uh73idBJcP8lruu4aILBSREyJyXETmO8pwK3XVog+d4as/99K/eVXa1i7n3DcLKAv3LYTS1a1mh/0rnPt+yh7xMTC+M5zYbhXe0+SQL/LaxDQNmAVUAipjlf2e7qyglPtKy8ziudkbqFDSj5e6OaFpKTfFy1lJIjDYKu53YFXBvK8qGCd3wTedIPEYDP4e6uj4mfyS1wThb4yZYozJdNymYi3uo9RV+fjXXew6nshbvRtR0s+n4N64eHkrSZSsZFXxPPh3wb23cp7YaOvKISsNhiyC6q3tjsit5DVBLBaRF0QkRESqi8jzwI8iUkZE8lg0R3m6TTHxfPHHXvqGB3NL3fIFH0CJinDfImuhoam9rYJtynXtX2H1L/kUg6FLoFITuyNyO3ldMGjfZTYbY4zt/RG6YFDhlp6ZTfdP/yIuKZ2fR9xMoH8BXj1c6Gystf7wmYNw97dQ+zb7YlHXZss8mDvM6lsaPM8q/66uyXUvGGSMCb3MzfbkoAq/T5ftZvvRBN7q1cje5ABQsjIMXQxl61jrAuhkOtdhDKz8FGYPgcph1pWDJgenyesoJh8ReVJE5jhuw0XE5t9y5Sq2xMbz+bLd9GpahdsaVLA7HMu50U1VmsGcoRA9ze6I1JVkZ8Hi52Hpv6FBd7h3PgRcY+0ulSd57YMYC4QDnztu4Y7nlLqsjKxsnpu9kVL+vrx2ZwO7wzlfsVLWqJfQdjDvUfh7nN0RqUtJT4aZg2HNl9BqOPSdaPU9KKfK6wyl5saYnD1Av4nIBmcEpNzL2N/3sPXIWb4YFE4pf1+7w7mYbwAMmAlz7re+nSYcgVtf1UWHCpOEo1bRvcProMu7cOPDdkfkMfL6W5AlIjXPPXBMkstyTkjKXWw/epZPftvFnU0q07lhRbvDuTQfP+g3GSLuh78+hO+HWWsIKPsdjoIvb4Hj26zlZTU5FKi8XkE8BywTkb2OxyGAlt9Ql5SZlc3zczZS0s+H17vfYHc4V+ZdBLp9AIFV4dfXrW+td0+1mqGUPTbOggVPQEB5eGApVGxkd0QeJ69XECuAcUA2EOe4r9NR1SWN/X0PG2PieaNHQ8oEFMKmpdyIQNtnrCJvB1dbE7DOHLI7Ks+TnQU/vwpzH4IqETBsmSYHm+Q1QUwGQoH/AJ8ANYApzgpKubYtsfF87Gha6ta4kt3hXL3G/WDQd3D2MHx1i5bmKEjJcVZhxRUfQcQDVrHFgHwsBa+uSl4TRENjzIPGmGWO20OAC7QbqIKWnpnNs7M2UMrflzdcoWnpUmrcDA/+aq1KN+lOiBxvd0TuLyYSxrWDfX9YKwPe8QF462h6O+U1QawTkX+WGBWRGwGdtqwu8vGvu9h+NIG3ezWitKs0LV1KuTpWkqjRHhaNsG6Z6XZH5X6MgdVjrSY9Ebj/J2vAgLJdXjupw4GVInLQ8bgasENENmGV2mjslOiUS1l/8DSf/76bu8KDC8+EuOtVrJS18NCvb8CKMdZomj7f6Ozd/JJyBhYMh20LoW5X6Pk5FCttd1TKIa8JorNTo1AuLzUji2dnb6BiST9eKWwT4q6Xlzd0fN3qKF3wJHxxE/QaB3Vutzsy17ZvOXz/iDX35Pb/WhPgROyOSuWQ1xXlDjg7EOXaRv+0g70nkpj6wI0FW8a7IDXqa1UMnT0Ept0FrZ+EDq9qO/nVykyD3/5j1VQqU8Mawhqca604ZTOdLqqu25p9cXyzYh+DW1bnptpuPuKkbG148BerjXzlxzChC5zaY3dUruPYFvjqVlj5CYQPgUeWa3IoxDRBqOuSlJbJyNkbqFranxe61LM7nILhU8waZdN3ApzcCWPbWHWcsrPtjqzwykiF396EcTdbK78NnAV3jrFKnahCSxOEui5vL97GodPJjL6rCQFF89ql5SYa9obHVkPITVYdp8nd4bS2xl5k/wqr3+bPd////0yXBXUJmiDUNVu+6wRTVx/kgTahtAj10IUFS1aGe2bDnR9D7Hr4vJXVfJKVYXdk9ks6aZXKmNjVWhJ00HfQ+0ud+OZCNEGoaxKfnMHzczZSs1wAIzvVtTsce4lA+H3w6EoIaQNLX7YmfHnqDOzMdKsD+uNmsP5ba3TSY6uhlq7c52qcmiBEpLOI7BCR3SLyQi7bq4nIMhFZLyIbRaRrjm0vOo7bISJ6PVrIvDJ/M8cT0vigXxh+Pt52h1M4lK5uta3f/S2knoUJna1hnPExdkdWMIyBHYthbCtrUZ+qzeGxVdDpTe1rcFFOazQWEW/gM6AjEAOsFZEFxpitOXZ7GZhljBkrIg2AH4EQx/3+WOU8KgO/iEgdY4yWGC8E5kcfZsGGWJ7tWIcmVbXa6XlEoP4dUPMW+PM9WPUZbJ5rlalu+4x7TgIzBvb+DsvehJi1EFQbBs7WeSJuwJlXEC2A3caYvcaYdGAG0OOCfQxQ0nE/EIh13O8BzDDGpBlj9gG7Ha+nbBZzOpmX520mvHppHm1f88oHeCrfALhtFDwRZXXMrvwEPmoCf46G1Hi7o8sfxliT3SZ2gyk94ewRuGOM1dSmycEtOHPYSRUgZ63kGODGC/YZBSwVkSeAAOBcI2UVYPUFx15U20BEhgHDAKpVq5YvQatLy8o2PDtrA9nZhg/7hVHEW7uwrqhUNej1hdUO/+sb1gSxFR9B8weg5WNQvLzdEV69rAzYOh9WfWp1zBevCF3es/phihS1OzqVj+welzgAmGiMeV9EWgFTRKRhXg82xnwJfAkQERFhnBSjcvhq+V7+3hfHe30bUy3I3+5wXEvFhnDPLIiNtlat+2sMrPrcmp0dPtSaLFbYy0ycjYUN02HteDgbYzUl3fEhNBmg60O7KWcmiMNA1RyPgx3P5fQAjjpPxphVIuIHlM3jsaoAbT4cz/tLd9ClYUX6hgfbHY7rqhwG/SZZs69XfgKbZkP0t1ChkfUNvEFPKF7O7ij/X0Yq7PoJ1k+F3b+AyYbQdlYp7lodde1uNyfGOOeLt4gUAXYCHbD+uK8FBhpjtuTYZzEw0xgzUUTqA79iNSU1AKZh9TtUdjxf+3Kd1BERESYyUiuQO0NqRhZ3fPIXZ1My+Onpdq5fxrswSUuwksTa8XBsE4gXhLSFG3pB3S5Qwoa1vFPPwu6frQqrO5dCRhKUqAxhA61bkPY9uRMRiTLG5FrvxGlXEMaYTBEZDvwEeAPjjTFbROQNINIYswB4FvhKREZgdVgPMVbG2iIis4CtQCbwuI5gss87i7ez+3gik+9vockhvxUtYdV1Ch8KxzbDlnmwZS4setq6lW8ANW6xFjCq3Mw5VxcpZyB2ndXhvH85HF4HJstaC7pxP6h/p7UmhpcOZ/Y0TruCKGh6BeEcv+84zpAJaxnaJoTX7nThFeJciTFWstj9izV89MAqayYyQGBVq6JsubpQOhTKhEJgMBQrYyWb3PoxjIH0JEiJs9bYPnMQzhywCucd3Qin91v7eRWxklBoO2tSW9UWmhQ8wOWuIDRBqEs6npBK14+WUybAlwXDb9IJcXbJSLGW4zwSbY0aio22/qhfeFEt3tYSqV5FrEQhXpCZajUZ5XYBXqYGVGwMlRpbSadqSyhavEBOSRUetjQxKdeW7RjSmpiWybSHWmpysJNPMQhta93OycqwZmif3gfxhyH1DKSctpqLTJZ11YABb18rafgFgl8pKFUVSlW3rjp0SKq6Ak0QKlfj/tzL8l0nebt3I+pUKGF3OOpC3j5W81KZULsjUW5Mx6ipi0QdOM3opTvo1qgS/ZtXvfIBSim3pAlCnSc+JYMnp6+nUqAfb/VuhBT2yVtKKafRJib1D2MML83dxLGzqcx6pBWBxXStZaU8mV5BqH9MX3OIHzYd4dnb69KsmhtWHVVKXRVNEAqA7UfP8vrCLbStXZaH29WwOxylVCGgCUKRkJrBo1PXUbKYD+/3a4KXl/Y7KKW0D8LjGWN4fs5GDsYlM+3BGylfws/ukJRShYReQXi4b/7ax+LNR3m+U11urBFkdzhKqUJEE4QHi9wfxzuLt3N7gwoM034HpdQFNEF4qJOJaTw+bR1VShfjvbua6HwHpdRFtA/CA2VlG56asZ4zyRl8/1gLne+glMqVJggPNHrpDlbsPsW7fRvToHJJu8NRShVS2sTkYRZuiGXs73sY0KIa/SK0zpJS6tI0QXiQzYfjeW7OBpqHlOb17rr4j1Lq8jRBeIiTiWk8PCWKMv6+fH5POL5F9KNXSl2e9kF4gPTMbB6buo6TiWl892hrypXQhWKUUlemCcIDvLFoC2v2x/FR/zAaVgm0OxyllIvQdgY3N2XVfqauPsgjN9ekR1gVu8NRSrkQTRBubNn247y2YAu31S/Pc53q2h2OUsrFaIJwU1ti4xk+bR0NKpfko/5N8dYKrUqpq6QJwg0diU/h/olrCSzmw/j7mhNQVLualFJXT/9yuJmE1AyGTlhLUloWcx5tRfmSWr5bKXVtNEG4kYysbIZPW8+u44lMGNKcehW1jIZS6tppE5ObyM42jJy9gT92nuDNng1pV6ec3SEppVycJgg3YIzh9YVbmB8dy/Od69K/RTW7Q1JKuQFtYnIDH/+6m0mrDvBQ21Aevbmm3eEoZauMjAxiYmJITU21O5RCxc/Pj+DgYHx88l7eXxOEi5u8aj8f/rKTvuHBvNS1vi78ozxeTEwMJUqUICQkRH8fHIwxnDp1ipiYGEJDQ/N8nDYxubA5UTGOiXAVeKd3I/1lUApITU0lKChIfx9yEBGCgoKu+qpKE4SLmrsuhufmbKBNzbJ8OrApRbz1o1TqHE0OF7uW/xP9q+KC5q0/zLOzN9C6ZhBf3RuBn4+33SEppdyQJggXMz/6MM/MiqZlaBBf39ucYr6aHJRydRMnTmT48OG5bmvdujUA+/fvZ9q0aZd8jUmTJlG7dm1q167NpEmT8iUuTRAu5LuoGEbMjKZFaBm+GRKhyUEpD7By5Urg8gkiLi6O119/nb///ps1a9bw+uuvc/r06et+bx3F5CLG/7WPNxZt5aZaZRk3OBx/X/3olLqS1xduYWvs2Xx9zQaVS/LanZdfsvfNN99k0qRJlC9fnqpVqxIeHs7IkSNp3749o0ePJiIigpMnTxIREcH+/fsBOHToEO3bt+fw4cMMGjSI1157DYDixYuTmJjICy+8wLZt2wgLC+O+++5jxIgR/7zfTz/9RMeOHSlTpgwAHTt2ZMmSJQwYMOC6ztWpf2VEpDPwEeANfG2MeeeC7R8Ctzge+gPljTGlHNuygE2ObQeNMd2dGWthZYxhzC+7+OjXXXS+oSIfDQijaBG9clCqsIqKimLGjBlER0eTmZlJs2bNCA8Pv+Jxa9asYfPmzfj7+9O8eXO6detGRETEP9vfeecdRo8ezaJFiy469vDhw1StWvWfx8HBwRw+fPi6z8VpCUJEvIHPgI5ADLBWRBYYY7ae28cYMyLH/k8ATXO8RIoxJsxZ8bmC7GzDG4u2MnHlfvqGB/NO70Y6Wkmpq3Clb/rOsHz5cnr16oW/vz8A3bvn7bttx44dCQoKAqB379789ddf5yUIOzjzr00LYLcxZq8xJh2YAfS4zP4DgOlOjMelpKRn8di365i4cj/3twnl3T6NNTko5eKKFClCdnY2wEVzEi4chno1w1KrVKnCoUOH/nkcExNDlSrXv4KkM//iVAEO5Xgc43juIiJSHQgFfsvxtJ+IRIrIahHpeYnjhjn2iTxx4kR+xW274wmp9P9yFT9tPcrL3erzyh318dIFf5RyCe3atWPevHmkpKSQkJDAwoUL/9kWEhJCVFQUAHPmzDnvuJ9//pm4uDhSUlKYN28ebdq0OW97iRIlSEhIyPU9O3XqxNKlSzl9+jSnT59m6dKldOrU6brPpbB8Je0PzDHGZOV4rroxJgIYCIwRkYuKDBljvjTGRBhjIsqVc4/qpTuPJdDrs5XsPJbIuEHhPNi2hk76UcqFNGvWjLvvvpsmTZrQpUsXmjdv/s+2kSNHMnbsWJo2bcrJkyfPO65Fixb06dOHxo0b06dPn4ualxo3boy3tzdNmjThww8/PG9bmTJleOWVV2jevDnNmzfn1Vdf/afD+nqIMea6XyTXFxZpBYwyxnRyPH4RwBjzdi77rgceN8asvMRrTQQWGWPm5LYdICIiwkRGRuZH6Lb5YeMRnp+zgYCiRfjmvuY0Cg60OySlXM62bduoX7++3WH8Y9SoURQvXpyRI0faHUqu/zciEuX4Mn4RZ15BrAVqi0ioiPhiXSUsuHAnEakHlAZW5XiutIgUddwvC7QBtl54rLvIyMrmv4u28vi0ddStWIIFw2/S5KCUsp3TRjEZYzJFZDjwE9Yw1/HGmC0i8gYQaYw5lyz6AzPM+Zcy9YFxIpKNlcTeyTn6yZ0cjU/lyRnrWbMvjiGtQ3ipa318ixSWlj+l1PUaNWqU3SFcM6fOgzDG/Aj8eMFzr17weFQux60EGjkztsLgx01HeHHuJtIzsxlzdxg9m17/qAOllMovOh3XBgmpGYxasJXv1sXQJDiQD+8Oo0a54naHpZRS59EEUcCWbjnKqAVbOHo2lSc71OaJW2vho/MblFKFkCaIAnIkPoXX5m9h6dZj1KtYgk/vaUazaqXtDksppS5Jv7o6WVJaJmN+2UmH9//gz10n+Ffneix84iZNDkq5qVOnThEWFkZYWBgVK1akSpUq/zw+V34jpx07dtC+fXvCwsKoX78+w4YNy/V1nVHO+0r0CsJJ0jOzmRMVw4e/7OREQhpdG1Xkhc71qRZ08Q+IUsp9BAUFER0dDVw8B6J48Yv7Gp988klGjBhBjx5WJaJNmzZdtM+5ct6RkZGICOHh4XTv3p3SpZ37RVMTRD5LTs9kxppDfLV8L0fiU4moXppxg8P1ikEpOyx+AY5e/Af3ulRsBF3eufJ+eXTkyBGCg4P/edyo0cUDOJ1VzvtKNEHkk70nEpm59hCzIg9xOjmDFqFleLt3I26uU05LZSilLmnEiBHceuuttG7dmttvv52hQ4dSqlSp8/ZxVjnvK9EEcR1OJKTx89ZjzI8+zN/74vD2Em6rX55h7WoQXv3666Aopa5TPn7Td5ahQ4fSqVMnlixZwvz58xk3bhwbNmygaNGidoemCeJqpGZkEX3oDGv2xbF81wkiD5zGGAgJ8ue5TnW5KzyY8iX97A5TKeViKleuzP3338/9999Pw4YN2bx583mLDFWpUoXff//9n8cxMTG0b9/e6XFpgshFQmoGR+NTiY1PZdexBHYcTWDHsQS2H00gPTMbEahfsSRPdahN54YVqVuhhDYjKaWuyZIlS+jQoQM+Pj4cPXqUU6dOXbSWQ6dOnXjppZf+WWd66dKlvP32RXVP853HJ4i4pHTuHreK1MwsUtKzSU7PJDk967x9yhb3pW7FEgxpHUKLkDI0DylDoL+PTRErpVxVcnLyeR3SzzzzDDExMTz11FP4+VmtD++99x4VK1Y877ic5byBfCvnfSVOK/dd0K613HdSWiYjZ2/Az8cbPx9vivl4U75kUSoF+lEpsBg1ygVQtrj9bYFKqbwpbOW+C5OrLfft8VcQAUWLMHbQlRcUV0opT6MzqZVSSuVKE4RSyu24S9N5frqW/xNNEEopt+Ln58epU6c0SeRgjOHUqVP/dITnlcf3QSil3EtwcDAxMTGcOHHC7lAKFT8/v/NGUOWFJgillFvx8fEhNDTU7jDcgjYxKaWUypUmCKWUUrnSBKGUUipXbjOTWkROAAeu4yXKAifzKRxX4Wnn7GnnC3rOnuJ6zrm6MaZcbhvcJkFcLxGJvNR0c3flaefsaecLes6ewlnnrE1MSimlcqUJQimlVK40Qfy/L+0OwAaeds6edr6g5+wpnHLO2gehlFIqV3oFoZRSKleaIJRSSuXK4xOEiHQWkR0isltEXrA7noIgIvtFZJOIRIvI1S/D5wJEZLyIHBeRzTmeKyMiP4vILse/pe2MMb9d4pxHichhx2cdLSJd7Ywxv4lIVRFZJiJbRWSLiDzleN4tP+vLnK9TPmeP7oMQEW9gJ9ARiAHWAgOMMVttDczJRGQ/EGGMcdvJRCLSDkgEJhtjGjqeexeIM8a84/gyUNoY8y8748xPlzjnUUCiMWa0nbE5i4hUAioZY9aJSAkgCugJDMENP+vLnG8/nPA5e/oVRAtgtzFmrzEmHZgB9LA5JpUPjDF/AnEXPN0DmOS4PwnrF8ttXOKc3Zox5ogxZp3jfgKwDaiCm37Wlzlfp/D0BFEFOJTjcQxO/M8uRAywVESiRGSY3cEUoArGmCOO+0eBCnYGU4CGi8hGRxOUWzS15EZEQoCmwN94wGd9wfmCEz5nT08QnuomY0wzoAvwuKNpwqMYq23VE9pXxwI1gTDgCPC+veE4h4gUB74DnjbGnM25zR0/61zO1ymfs6cniMNA1RyPgx3PuTVjzGHHv8eB77Ga2jzBMUcb7rm23OM2x+N0xphjxpgsY0w28BVu+FmLiA/WH8tvjTFzHU+77Wed2/k663P29ASxFqgtIqEi4gv0BxbYHJNTiUiAo3MLEQkAbgc2X/4ot7EAuM9x/z5gvo2xFIhzfyQdeuFmn7WICPANsM0Y80GOTW75WV/qfJ31OXv0KCYAx3CwMYA3MN4Y86bNITmViNTAumoAa8nZae54ziIyHWiPVQb5GPAaMA+YBVTDKg3fzxjjNp26lzjn9ljNDgbYDzyco23e5YnITcByYBOQ7Xj6Jax2ebf7rC9zvgNwwufs8QlCKaVU7jy9iUkppdQlaIJQSimVK00QSimlcqUJQimlVK40QSillMqVJgilrpGIlBKRxxz3K4vIHLtjUio/6TBXpa6RoxbOonOVU5VyN0XsDkApF/YOUFNEooFdQH1jTEMRGYJVPTQAqA2MBnyBwUAa0NUYEyciNYHPgHJAMvCQMWZ7wZ+GUrnTJialrt0LwB5jTBjw3AXbGgK9gebAm0CyMaYpsAq417HPl8ATxphwYCTweYFErVQe6RWEUs6xzFGvP0FE4oGFjuc3AY0d1ThbA7Ot8joAFC34MJW6NE0QSjlHWo772TkeZ2P93nkBZxxXH0oVStrEpNS1SwBKXMuBjhr++0TkLrCqdIpIk/wMTqnrpQlCqWtkjDkFrBCRzcB71/AS9wAPiMgGYAu63K0qZHSYq1JKqVzpFYRSSqlcaYJQSimVK00QSimlcqUJQimlVK40QSillMqVJgillFK50gShlFIqV/8HpbGAKJyHeKQAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "def plot_population(result):\n", " fig, ax = plt.subplots()\n", " ax.plot(\n", " result.times,\n", " np.array(result.expect[0]) + np.array(result.expect[1]),\n", " label='qubit 0',\n", " )\n", " ax.plot(\n", " result.times,\n", " np.array(result.expect[0]) + np.array(result.expect[2]),\n", " label='TLS 0',\n", " )\n", " ax.legend()\n", " ax.set_xlabel('time')\n", " ax.set_ylabel('population')\n", " plt.show(fig)\n", "\n", "\n", "plot_population(guess_dynamics)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The population dynamics of qubit and TLS ground state show that\n", "both are oscillating and especially the qubit's ground state population reaches\n", "a maximal value at intermediate times $t < T$. This maximum is indeed the\n", "maximum that is physically possible. It corresponds to a perfect swap of\n", "the initial qubit and TLS purities. However, we want to reach this maximum at\n", "final time $T$ (not before), so the guess control is not yet working as desired." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Optimize" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Our optimization target is the ground state $\\ket{\\Psi_{q}^{\\tgt}}\n", "= \\ket{0}$ of the qubit, irrespective of the state of the TLS. Thus, our\n", "optimization functional reads\n", "\n", "$$\n", " F_{re} = 1 -\n", "\\Braket{\\Psi_{q}^{\\tgt}}{\\tr_{t}\\{\\op{\\rho}(T)\\} \\,|\\; \\Psi_{q}^{\\tgt}}\\,,\n", "$$\n", "\n", "and we first define `print_qubit_error`, which prints out the\n", "above functional after each iteration.\n" ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:27.124250Z", "start_time": "2019-02-12T04:47:27.118668Z" } }, "outputs": [], "source": [ "def print_qubit_error(**args):\n", " \"\"\"Utility function writing the qubit error to screen\"\"\"\n", " taus = []\n", " for state_T in args['fw_states_T']:\n", " state_q_T = trace_TLS(state_T)\n", " taus.append(state_q_T[0, 0].real)\n", " J_T_re = 1 - np.average(taus)\n", " print(\" qubit error: %.1e\" % J_T_re)\n", " return J_T_re" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In order to minimize the above functional, we need to provide the correct\n", "`chi_constructor` for the Krotov optimization. This is the only place where the\n", "functional (implicitly) enters the optimization.\n", "Given our bipartite system and choice of $F_{re}$, the equation for\n", "$\\op{\\chi}(T)$ reads\n", "\n", "$$\n", " \\op{\\chi}(T) =\n", " \\sum_{k=0,1} a_{k}\n", "\\op{\\rho}_{q}^{\\tgt} \\otimes \\ket{k}\\bra{k}\n", "$$\n", "\n", "with $\\{\\ket{k}\\}$ a\n", "basis for the TLS Hilbert space." ] }, { "cell_type": "code", "execution_count": 18, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:27.144420Z", "start_time": "2019-02-12T04:47:27.132277Z" }, "attributes": { "classes": [], "id": "", "n": "18" } }, "outputs": [], "source": [ "def TLS_onb_trg():\n", " \"\"\"Returns the tensor product of qubit target state\n", " and a basis for the TLS Hilbert space\"\"\"\n", " rho1 = qutip.Qobj(np.kron(rho_q_trg, np.diag([1, 0])))\n", " rho2 = qutip.Qobj(np.kron(rho_q_trg, np.diag([0, 1])))\n", " return [rho1, rho2]\n", "\n", "\n", "TLS_onb = TLS_onb_trg()\n", "\n", "\n", "def chis_qubit(fw_states_T, objectives, tau_vals):\n", " \"\"\"Calculate chis for the chosen functional\"\"\"\n", " chis = []\n", " for state_i_T in fw_states_T:\n", " chis_i = np.zeros(shape=(4, 4), dtype=np.complex_)\n", " for state_k in TLS_onb:\n", " a_i_k = krotov.optimize._overlap(state_i_T, state_k)\n", " chis_i += a_i_k * state_k\n", " chis.append(qutip.Qobj(chis_i))\n", " return chis" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We now carry out the optimization for five iterations." ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Iteration 0\n", " objectives:\n", " 1:Ļā‚€[4,4] to Ļā‚[4,4] via [š“›ā‚€[[4,4],[4,4]], [š“›ā‚[[4,4],[4,4]], uā‚‚(t)]]\n", " adjoint objectives:\n", " 1:Ļā‚‚[4,4] to Ļā‚ƒ[4,4] via [š“›ā‚‚[[4,4],[4,4]], [š“›ā‚ƒ[[4,4],[4,4]], uā‚‚(t)]]\n", " S(t) (ranges): [0.000000, 1.000000]\n", " duration: 0.5 secs (started at 2019-06-28 14:32:33)\n", " optimized pulses (ranges): [0.00, 2.00]\n", " āˆ«gā‚(t)dt: 0.00e+00\n", " Ī»ā‚: 1.00e-02\n", " storage (bw, fw, fw0): None, None, None\n", " fw_states_T norm: 1.000000\n", " Ļ„: (7.97e-01:0.00Ļ€)\n", " qubit error: 1.1e-01\n", "Iteration 1\n", " duration: 2.2 secs (started at 2019-06-28 14:32:34)\n", " optimized pulses (ranges): [0.00, 2.06]\n", " āˆ«gā‚(t)dt: 7.72e-02\n", " Ī»ā‚: 1.00e-02\n", " storage (bw, fw, fw0): [1 * ndarray(2500)] (1.3 MB), None, None\n", " fw_states_T norm: 1.000000\n", " Ļ„: (7.98e-01:0.00Ļ€)\n", " qubit error: 1.1e-01\n", "Iteration 2\n", " duration: 2.4 secs (started at 2019-06-28 14:32:36)\n", " optimized pulses (ranges): [0.00, 2.23]\n", " āˆ«gā‚(t)dt: 5.72e-01\n", " Ī»ā‚: 1.00e-02\n", " storage (bw, fw, fw0): [1 * ndarray(2500)] (1.3 MB), None, None\n", " fw_states_T norm: 1.000000\n", " Ļ„: (8.01e-01:0.00Ļ€)\n", " qubit error: 6.7e-02\n", "Iteration 3\n", " duration: 2.5 secs (started at 2019-06-28 14:32:38)\n", " optimized pulses (ranges): [0.00, 2.33]\n", " āˆ«gā‚(t)dt: 8.11e-02\n", " Ī»ā‚: 1.00e-02\n", " storage (bw, fw, fw0): [1 * ndarray(2500)] (1.3 MB), None, None\n", " fw_states_T norm: 1.000000\n", " Ļ„: (7.99e-01:0.00Ļ€)\n", " qubit error: 5.0e-02\n", "Iteration 4\n", " duration: 2.4 secs (started at 2019-06-28 14:32:41)\n", " optimized pulses (ranges): [0.00, 2.19]\n", " āˆ«gā‚(t)dt: 2.16e-01\n", " Ī»ā‚: 1.00e-02\n", " storage (bw, fw, fw0): [1 * ndarray(2500)] (1.3 MB), None, None\n", " fw_states_T norm: 1.000000\n", " Ļ„: (8.02e-01:0.00Ļ€)\n", " qubit error: 4.9e-02\n", "Iteration 5\n", " duration: 2.8 secs (started at 2019-06-28 14:32:43)\n", " optimized pulses (ranges): [0.00, 2.15]\n", " āˆ«gā‚(t)dt: 6.38e-02\n", " Ī»ā‚: 1.00e-02\n", " storage (bw, fw, fw0): [1 * ndarray(2500)] (1.3 MB), None, None\n", " fw_states_T norm: 1.000000\n", " Ļ„: (8.03e-01:0.00Ļ€)\n", " qubit error: 4.9e-02\n" ] } ], "source": [ "# NBVAL_IGNORE_OUTPUT\n", "# the DensityMatrixODEPropagator is not sufficiently exact to guarantee that\n", "# you won't get slightly different results in the optimization when\n", "# running this on different systems\n", "oct_result = krotov.optimize_pulses(\n", " objectives,\n", " pulse_options,\n", " tlist,\n", " propagator=krotov.propagators.DensityMatrixODEPropagator(\n", " atol=1e-10, rtol=1e-8\n", " ),\n", " chi_constructor=chis_qubit,\n", " info_hook=krotov.info_hooks.chain(\n", " krotov.info_hooks.print_debug_information, print_qubit_error\n", " ),\n", " check_convergence=krotov.convergence.check_monotonic_error,\n", " iter_stop=5,\n", ")" ] }, { "cell_type": "code", "execution_count": 20, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:39.603540Z", "start_time": "2019-02-12T04:47:39.598097Z" }, "attributes": { "classes": [], "id": "", "n": "20" } }, "outputs": [ { "data": { "text/plain": [ "Krotov Optimization Result\n", "--------------------------\n", "- Started at 2019-06-28 14:32:33\n", "- Number of objectives: 1\n", "- Number of iterations: 5\n", "- Reason for termination: Reached 5 iterations\n", "- Ended at 2019-06-28 14:32:46" ] }, "execution_count": 20, "metadata": {}, "output_type": "execute_result" } ], "source": [ "oct_result" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Simulate the dynamics of the optimized field" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The plot of the optimized field shows that the optimization slightly shifts\n", "the field such that qubit and TLS are no longer perfectly in resonance." ] }, { "cell_type": "code", "execution_count": 21, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:39.774187Z", "start_time": "2019-02-12T04:47:39.607074Z" }, "attributes": { "classes": [], "id": "", "n": "21" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEGCAYAAABo25JHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO3deXxcVf3/8ddnJpN9bZKu6UoXqKWUtoCsooggoKAswldF+LIooH5dEFERXH8uKIqKSFF2FNmKiGVTWQoW6AJ0o7U7bWibJs022Sc5vz9mEkKbZdrmzp1k3s/HYx6ZuXNn7mc6j97PnHPu+RxzziEiIqkr4HcAIiLiLyUCEZEUp0QgIpLilAhERFKcEoGISIpL8zuAfVVSUuImTJjgdxgiIoPK0qVLK51zpT09N+gSwYQJE1iyZInfYYiIDCpmtqW359Q1JCKS4pQIRERSnBKBiEiKUyIQEUlxSgQiIilOiUBEJMUpEYiIpLhBN49AYGddM29sraG8uon65ghpQSM/K8RBpTlMLs2lNC8DM/M7TBEZJJQIBonG1ggPvLaVR5ZtY9U7dX3uW5yTzpETh3HMQcWcdMgIRhdmJShKERmMbLAtTDN37lyXSjOLnXM8sqycG59ew866FmaWFXDGzFHMnTCMicU55GeFiHR0UBVuZeOuBtZX1LO8vJZXN+6mvKYJgBlj8jn5kJGcOmMkU0fkqrUgkoLMbKlzbm6PzykRJK/65jaueXg5T67cwayxhVx3+iHMnTAs7tevrwjz7OqdPLN6B6+/XQPApJIcPvK+kZw4rZTZ44pIT9MwkUgqUCIYhCrqmvncnYtZt7Oea06dxqXHTSIQ2P9f8hV1zTyzeidPrdzBoo1VtHc4ctKDHH1QMcdPKeX4KSVMLMlRa0FkiFIiGGR21DZz3m2LqAy38IfPzOGEqT0WDNxvdc1tLNpQxcJ1u1i4rpItVY0AjCnM4oSpJRw/pZRjDyqhIDs0oMcVEf8oEQwitY1tnHfbIsprmrj3kiM5fFyR58fcUtXAwnWVLFy3i/+sr6K+JYIZTC7NZWZZIYeNLWBmWSEHj8wjMxT0PB4RGXh9JQJdNZRE2to7uPzeJWyqbOCui49ISBIAGF+cw/jiHD7z/vFE2jt4c1sNL62r4s1tNbzw3woeWbYNgIDB2GHZTC7NZfKIXCaX5jJlRB4TirMpyAqpW0lkkFIiSCI/f2oNr27aza8+dRjHTC7xJYa0YIA544cxZ3x0UNo5xzu1zSzfWsOaHfWsrwizviLMwnWVtLZ3dL0uNyONsqKs2C2bEfmZlOZlUJqXwfDY32HZ6Qc0ziEi3lAiSBILVmzn9oWbuPDo8Xzi8DK/w+liZowpzGJMYRYfPXRU1/ZIewdv725kfUWYt3c3sq26KXZr5NWNu6lviez1XsGAUZKbTkluBkXZ6RRkhyjKDlGYlU5hdojC7HQKs0IU5YQoyEqnKDtEQVaItKCubBLxkhJBEiivaeKah5dz+LhCrjt9ut/hxCUtGGBSaS6TSnN7fL6xNcKu+pauW0XX32Yqw63UNLbyTm0TNY1t1DS20tHHUFVOepD8rGhSyM8MkZ+VFvsbu2Wm9fh8QXaI3PQ0tUJE+qFE4DPnHN9+dAUdzvGb8w8fMtf1Z6enMb44jfHFOf3u29HhqG+JUNvYRk1TK9Wx5BBNEm3UNbdR19RGbVP0/js1zaxprqe2qY365r1bHt2ZQV5GL4miWyIpyA69m1y6nksjKxTU2IcMeZ4lAjMbC9wDjAAcMM85d/Me+xhwM3Aa0Ahc5Jxb5lVMyejRZeW88N9dfO9j0xk7LNvvcHwRCBgFsRP1OPbt36C9wxFuibwnUdQ1RbqSR11TG3XNkehzsec3VzZ27dvY2t7n+6cFrFsSSetKFPlZ722ZdD5fkBViRH4mw/My1KUlg4aXLYII8HXn3DIzywOWmtmzzrnV3fb5KDAldjsKuDX2NyVU1DfzgydWM3d8ERcePcHvcAalYLckMnY/Xt/W3kH9HomirinSLal0JphI1/PlNU3RZNPU9p4B8+4CBqV5GYwsyGJkfgajCrIYWZDJmMIsJhTnMKEkm7xMzdOQ5OBZInDObQe2x+7Xm9lbwBigeyI4E7jHRSczvGJmhWY2KvbaIe/6x1bR1NbOz86ZqX5sn4SCAYblpDMsJ32/Xt/c1t6VIGqbItQ0trKzroUdtU3sqGtme20zG3c1dM3P6K44J53xxdlMKMlhUkkO00bmc8ioPMYUZqk7ShIqIWMEZjYBOBx4dY+nxgBbuz3eFtv2nkRgZpcDlwOMGzfOqzATasGK7Ty1agfXnDqNg3oZcJXklxkKkhkKMjw/s999wy0RtlU3srmykc1VDWypamBzZSOLNlTx6LLyrv3yMtI4eFQeB4/MZ9bYQuaML2J8cbaSg3jG80RgZrnAI8BXnHN910/uhXNuHjAPojOLBzA8X1Q3tHL931YyY0w+lx8/ye9wJEFyM9I4eGQ+B4/M3+u5cEuEtTvqWbOjjjXbo3/nv17Ova9sAWBYTjqzxxUxZ3wRx04uZsboArUiZcB4mgjMLEQ0CdzvnHu0h13K4T1du2WxbUPaD59YTU1jG/f871EaUBQgmiTmjI+e6Dt1dDjW7wqzdEs1S7dUs+ztav751k4AirJDHBcrFnjitFKG5/XfIhHpjZdXDRnwJ+At59xNvez2OPBFM3uA6CBxbTKOD1SFW2hoaWdc8YFf1fPc2goefb2cL31oMtNH7/3LUKRTIGBMHZHH1BF5XHBktEu0MtzCy+srefG/0dpQf3/zHczgiAnDOG3GSE6dMYqRBUoKsm88KzpnZscBC4EVQOelFd8GxgE45/4QSxa/A04levnoxc65PivKJbro3O6GVk765fOEWyLMv/JYZowp2O/3CrdE+MhNL5CTkcYTXz6OjDQVcJP955xjzY56nlm1kydXbmfNjnoA5o4v4ry5Yzl95ihyMjRVSKJUffQAzHtxA/9vwRoAjp1czH2XHLXfg3bX/20l976yhUeuOIbZCSooJ6ljw64wT67YzvzXy9mwq4Gc9CAfnzWaC44cx8yyQr/DE5/1lQjUQd2PhesqmToilxs+Np2X11fx4rrK/XqfVzdWcc+iLVx8zEQlAfHEQaW5fPFDU/jn1z7AI1cczWmHjuKx19/h4797mfNuW8Q/V++ko69aHpKylAj60NHhWLK5mqMnFfPpo8Yzblg2P1nwFu37+J+ptrGNr/71DcYXZ3P1KVM9ilYkysyYM34YN557GK995yS+e8Z0yqubuPSeJZz8qxd4cMlWIr1MhJPUpETQh63VjTS1tTN9dD7paQG+cco01uyoZ/7r8V/Y5Jzj2keXU1Hfwm8vOJzsdPXZSuLkZYa45LiJvPCNE7n5/FlkhoJc8/ByPvLrF/nH8u1qIQigRNCntbHBt6kj8gA4/dBRHFZWwI1Pr6G2sS2u97jz5c08uXIH3zhlmvppxTdpwQBnzhrDE186jts+O4egGVf9eRkfv+UlFm2o8js88ZkSQR/+uzOaCKbEEkEgYPzwrBlUhVv55iPL6W+g/d9rdvKjf6zm5OkjuEwTxyQJmBmnvG8kT33lBG467zCqG9q44PZX+NJfXmd7bZPf4YlPlAj6sLmqkRH5GeR2uwRvZlkh15w6jadW7eDGp9f2+trn11ZwxX3LmD46n5vPn6VZoJJUggHjk7PL+NfXP8D/nTSFp1ft4KRfvsBtL2zQ+EEKUiLoQ3l1E2MKs/baftnxk7jgyHH8/vkNXP+3lTR1K2XcEmnnlufWc+ndS5g8PJd7/vcojQtI0soMBfnqyVP551c/wDEHlfCTJ9dw9q3/YV2sNSypQWeoPpTXNHHY2L379c2MH501g5z0IH98aRMLVmznhKmlBMxYuG4XO+taOP3QUfzk7EPJV6lhGQTGFWfzx8/N5Ynl7/Ddx1Zy+m9e4qsnT+Wy4yeqDEoKUCLoRUeHY3ttE6fPHNXj88GAcd0Z0zllxkj+uHAjL62rxAGHlRVy4znjOX5KiapFyqBzxszRHDWxmOseW8HPnlrD82sr+M0FhzMijuqqMngpEfSior6FtnbXY9dQd0dMGMYRE4YlKCoR75XmZfCHz8zhkWXlfPexlZx280Ju+tQsPjC11O/QxCNq8/ViW3UjAGOK+k4EIkORmXHOnDL+/qVjKcnN4HN3vMZNz/5X8w6GKCWCXuysawFgpJrEksImD8/jsauO5Zw5ZfzmX+u48v5lNLZG+n+hDCpKBL2oDEcTQWlehs+RiPgrKz3IjefM5LrTD+GZ1Ts4+9ZFXS1mGRqUCHqxq76FYMAoyt6/tWxFhhIz49LjJ3HHRUewbXcjZ93yMiu21fodlgwQJYJeVIZbGJaTTlATwUS6nDhtOPOvOpaMtCDnz1vEf9bvXzVeSS5KBL3YVd9CSa66hUT2NHl4Lo9ccQxjirK46M7FPLki6RYVlH2kRNCLynCLxgdEejGyIJMHP380h5YVcOWfl/HXxW/7HZIcACWCXkRbBBofEOlNYXY6911yFCdMKeWbj6zgz68qGQxWSgQ9cM5RGW5Vi0CkH1npQW777Bw+OK2Ub89fwf2vbvE7JNkPSgQ9qGuK0NreQanGCET6lRkK8ofPzuFDBw/nO/Oj63LL4KJE0IPKhugcgmJ1DYnEJSMtyK2fmc2HDxnOdx9bycNLt/kdkuwDJYIe1MRWHyvUHAKRuGWkBbnl07M5bnIJ33xkOc+s2uF3SBInJYIe1Da1Amgymcg+ykiLjhkcOqaAL/7ldf6zQfMMBgMlgh5UN8RaBFlaS0BkX+VkpHHnRUcwoTiby+5ewvJtNX6HJP1QIuhBTVM0EahFILJ/inLSufeSoyjKSeeiOxezubLB75CSwv2vbuH7f19FS6S9/50TSImgBzWNrZhBXqaWaxDZXyPyM7nvkqNwznHxXYupbmj1OyRfrSyv5TvzV3Lny5u56+XNfofzHkoEPahpbKMgK6QF50UO0ISSHG6/cC7l1U18/r6lSfdLOFGcc9zw+CpKctM5eGRe0l1VpUTQg+rGVnULiQyQuROGceO5M3lt026ufWQFzqXe4jaPvVHO0i3VXHPKwXzi8DGsqwhTUd/sd1hdlAh6UNsUbRGIyMA4c9YYrv7IVOa/Xs6v/7nO73ASKtwS4ScL1jCzrIBz5pRx5MTo0rZLN1f7HNm7lAh6UNPYRlG2EoHIQLrqg5M5Z04ZN/9rHf9YnjoVS3/37/VU1LfwvY+/j0DAOGRUPsGA8db2Or9D66JE0IPqxlZNJhMZYGbGjz8xgznji7j6oTdZ/U7ynAi9sqmygT+9tJFPzh7D7HFFQLQkx/jibNburPc5uncpEfSgtlFdQyJe6CxFUZAV4rJ7lrB7iF9J9MMnVpORFuTaUw9+z/ZpI/JYtzPsU1R7UyLYQ1t7B/UtEQ0Wi3hkeF4m8y6cQ2W4hSvvX0pbe4ffIXniuTUV/HtNBV8+aTLD8zPf89yU4blsrmqgNZIcn12JYA+1TZ11htQiEPHKzLJCfnr2obyycTc/emK13+EMuNZIBz98YjWTSnK46JiJez1fNiybDgfba5t8iG5vSgR7qGmMNlWVCES89YnDy7j8hEncvWgLD7w2tBa1ues/m9hY2cB3Pzad9LS9T7Nji7IB2FatRJCUOlsEGiMQ8d43Tz2Y46eUcP3fVrF0S/JcTnkgKuqaufmf6zjp4OF8cNrwHvcpK8oCYOvuxkSG1islgj3UN0cAyMtUIhDxWjBg/PaCwxlZkMkV9y1lZ13yTLLaXz99ag1t7Y7vnjG9131GFWQSDNjQbxGY2R1mVmFmK3t5/kQzqzWzN2K3672KZV+EWzoTgeoMiSRCYXY68y6cQ7glwhcGeRmKRRuqeHRZOZedMJEJJTm97pcWDDCqIJOt1UO/RXAXcGo/+yx0zs2K3X7gYSxxC8daBLkZSgQiiXLwyHx+ce5hvP52DTf8bdWgLEPR3NbOd+avYNywbL70oSn97j+6MIvtNcnRAvIsETjnXgR2e/X+XulsEeQoEYgk1GmHjuKqDx7EA4u3cv+rg2/w+NbnN7CxsoEfnTWDzFCw3/1H5GcmTb0hv8cIjjazN83sSTN7X287mdnlZrbEzJbs2rXL04Dq1SIQ8c3XTp7GB6eV8r3HV/HapsHzO3J9RZhbn9/AmbNGc8LU0rheMzwvg511LUnR+uk3EZjZVDP7V2dfv5nNNLPrBuDYy4DxzrnDgN8Cj/W2o3NunnNurnNubmlpfP/I+6uhJUJ2epCgSlCLJFwwYPz6/MMZOyybK+9fmjTX2felo8PxnfkryAwFuO703geI9zQ8L4OmtvauXgg/xdMiuB34FtAG4JxbDpx/oAd2ztU558Kx+wuAkJmVHOj7HqhwS0StAREfFWSFmPfZOTS1tvOFe5fS3Jbcg8f3vbqFVzft5tunHUJpXkbcrxsRm21cUd/iVWhxiycRZDvnXttj2wGnMDMbaWYWu39kLJaqA33fA1XfEiFXVwyJ+GrKiDxu+tQs3twWXdUrGbpPerKpsoGfLFjDB6aW8qkjxu7Ta4fHkkZFnf+JIJ4zXqWZHQQ4ADM7B+i3hqyZ/QU4ESgxs23ADUAIwDn3B+Ac4AoziwBNwPkuCb7tcHOEPLUIRHx3yvtG8uWTpvCbf63j0DH5XHTs3qUa/NTe4bj6oTcJBY2fnT2T2O/auA3PjyWCJBgwjueMdxUwDzjYzMqBTcBn+nuRc+6Cfp7/HfC7eIJMpHBLRFcMiSSJr5w0hdXv1PLDf7zFtJH5HH1Qsd8hdbl94UaWbqnm15+axciCzP5fsIfOQnTJ0CLot2vIObfROfdhoBQ42Dl3nHNus+eR+STcrDECkWQRCBi/+tQsJhRnc9Wfl1FekxyDx29sreGXz6zlozNGcuas0fv1HnkZaWSGAsndIjCzr/WyHQDn3E0exeSrsMYIRJJKXmaIeRfO5azfvcwldy3mwS8cTb6PJWBqGlu56v5ljMjP5Kef3PcuoU5mRnFOBlVh/9dk6KtFkBe7zQWuAMbEbl8AZnsfmj/CLRojEEk2B5XmcsunZ7O+Iszn7/GvDIVz0XGBivpmbvmf2RQcYJXi4tx0djcmcSJwzn3fOfd9oAyY7Zz7unPu68AcYFyiAkwk55xaBCJJ6oSppfz8nJks2ljF1Q8tp6Mj8deW/P75DfzzrQq+fdohHDa28IDfryg7PSlWaYvnjDcC6B5pa2zbkNPc1kF7hyM3Q5VHRZLRJ2eXsbOuhZ89tYbheRlcd/oh+901s6/+sXw7Nz69ljNnjeaiYyYMyHsW56SzYZf/S1bGkwjuAV4zs/mxx2cBd3sXkn/qW6JrEeRm9F8nRET88YUPTGJnXTN/emkTGWkBvnHKNM+TwRtba/jag28wZ3zRfl0q2pthOYOkReCc+7GZPQkcH9t0sXPudW/D8kdDS7TfUV1DIsnLzLj+jOm0tnfw++c34IBrPEwGa3fUc/Gdr1Gal8Ftn50TV0G5eBXlpNPY2k5zW/uAvu++6veMZ2bjgEpgfvdtzrnBVx6wH++WoFbXkEgyCwSMH505g4BFq362tHVw3emHEBjgGmHrK8J8+o+vkJ4W4L5LjqIkN/4SEvEozkkHYHdDK6MLswb0vfdFPD99/0FsVjGQBUwE1gK9VgsdrN7tGlKLQCTZBQLGD8+cQSgY4I6XN1Fe08ivP3U4WekD88v6za01XHL3YsC4/9L397nQzP4qSpJEEM+EskOdczNjtynAkcAi70NLvHCzVicTGUzMjBs+9j5u+Nh0nlm9k3Nv+w9bqhoO+H2fXrWD8+e9QmYoyAOXv5/Jw3MHINq9dbYIqnweJ9jn9Qicc8uAozyIxXed5WDVIhAZXC4+diK3f3Yub1c1ctrNC3l46bb9KlTX3NbO9x5fxefvXcqUEbnMv/JYz5IAvNsiqPY5EcQzRtB9hnGA6GSydzyLyEdanUxk8Prw9BE8+ZUT+Opf3+Dqh97kwcVb+c7p8V3v75zj2dU7+fGCt9hS1cjFx07gm6ce7PkAbrK0COI54+V1ux8hOmbwiDfh+EsL14sMbmMKs/jLZe/nwSVb+cXTaznzlpc5csIwzplbxolTS7sKvXXauruR59ZW8OdX32bNjnomD8/l3kuO5Pgp3i6A1Sk/M0QwYMnfIgBWO+ce6r7BzM4FHupl/0Er3BwhLWBkpPm9gqeI7K9gwLjgyHGcPnMUDy7eyt2LNnPNw8uB6HX7nesAbK9tprYpeoHIIaPy+fnZM/nE7DGEgon7/x8IGEXZoUHRIvgWe5/0e9o26HWWl0jUTEUR8U5+ZohLj5/EJcdN5K3t9by8vpKNlWEqY0Xe5k4oYnJpLidMLWViSY5v/++LstOTt0VgZh8FTgPGmNlvuj2VzwCsUJaMVIJaZOgxM6aPzmf66Hy/Q+lRYXaoq2Xil77Oeu8AS4CPA0u7ba8HvuplUH6p13rFIpJgBVkhymv8XZOg17Oec+5N4E0zu985NyRbAHtSi0BEEq0gK523ttf7GkNfXUMPOufOA143s70uyHXOzfQ0Mh80tEYYFrucS0QkEQqyQtT4vCZBXz9//y/294xEBJIMws0Rxg3L9jsMEUkhhdkhGlrbaWvvSOgVS9311TW0PfZ3S+LC8Vd9S0RzCEQkoQqyokUu65raKB7gonbx6qtrqJ53i80BWOyxAc45l5xD8AdAYwQikmidiaAmGROBcy6vt+eGokh7B01t7SpBLSIJ1bnusZ+XkMb189fMZgPHEW0RvDQUF6ZpaI0uSpOj1clEJIE6WwS1jf4lgn5HJszseqJLUxYDJcBdZnad14ElmuoMiYgfCrMGR4vg08BhzrlmADP7KfAG8CMvA0s0rU4mIn4oSIJEEM+1Su8A3Uv2ZQDl3oTjn3Dn6mRqEYhIAuV3Dhb72DUUz1mvFlhlZs8SHSM4GXits/6Qc+7LHsaXMPXNWpRGRBIvFAyQm5GW9F1D8+m2cD3wvDeh+EtjBCLil4IsfwvP9XvWc87dnYhA/Nag1clExCfRROBfmYl4rho6w8xeN7PdZlZnZvVmVpeI4BJJXUMi4he/WwTxDBb/GvgcUOycy3fO5Q3JWcVauF5EfBItPJfciWArsNI5t1cF0qEk3BwhOz1IMKDVyUQksfxenCaen7/XAAvM7AWgpXOjc+4mz6LyQViL0oiIT/zuGornzPdjIEx0LsGQLdbfuV6xiEiiFWSHaIl00NzWTmYo8WVu4jnzjXbOzfA8Ep+pRSAifinoNqlsZEHiE0E8YwQLzOwjnkfiM5WgFhG/dK1J0OxP91A8ieAK4Ckza9qXy0fN7A4zqzCzlb08b2b2GzNbb2bLYxVOfaMWgYj4xe96Q/0mgtjlogHnXNY+Xj56F3BqH89/FJgSu10O3BpPwF6pb9YYgYj4Iz/T31LU8a5HUET0hN1VfM4592Jfr3HOvWhmE/rY5Uzgnthlqa+YWaGZjepcIjPRwi0R8tQiEBEf+N011O+Zz8wuJbqQfRnR8tPvBxYBHzrAY48hOkeh07bYtr0SgZldTrTVwLhx4w7wsHtzztGgq4ZExCdJ3zVENAkcAWxxzn0QOByo8TSqPTjn5jnn5jrn5paWlg74+7dEOoh0ONUZEhFfdBa7TOZE0NxtUZoM59waYNoAHLscGNvtcRk+rXPQWWdIXUMi4oe0WCnquqaIL8ePJxFsM7NC4DHgWTP7G7BlAI79OHBh7Oqh9wO1fo4PgBalERH/+Dm7OJ4y1J+I3f2emT0HFABP9fc6M/sLcCJQYmbbgBuAUOw9/wAsAE4D1gONwMX7Ef+A0DKVIuK3vEz/FqfZp5/AzrkX9mHfC/p53gFX7cvxvaLKoyLit4KsUFJPKBvytDqZiPitICtEXRIPFg95nQvX66ohEfFLvo9jBHElAjMbb2Yfjt3PMrM8b8NKrLBWJxMRnyV1i8DMLgMeBm6LbSojegXRkFGvriER8VlBVoiG1nba2jsSfux4WgRXAccCdQDOuXXAcC+DSrRwc4S0gJGRpp4yEfFHfuyHqB+tgnjOfC3OudbOB2aWBgypZSs7y0uYaZlKEfFHQXZnvaHETyqLJxG8YGbfBrLM7GTgIeDv3oaVWPUqQS0iPvOz3lA8ieBaYBewAvg80Ylg13kZVKJpURoR8VtXKWofEkE8M4s7gNuB281sGFAWmww2ZGhRGhHxW1cp6mRsEZjZ82aWH0sCS4kmhF95H1riaOF6EfFbsncNFTjn6oBPEl1I5ijgJG/DSiy1CETEb/lJngjSzGwUcB7whMfx+CLcHNEcAhHxVWYoSHpawJd6Q/Ekgh8ATwPrnXOLzWwSsM7bsBJLLQIRSQZ+zS6OZ7D4IaKXjHY+3gic7WVQidTe4WhsbVedIRHxXb5Ppah7PfuZ2W/pY+KYc+7LnkSUYCpBLSLJItoiSPyEsr7OfksSFoWPVIJaRJJFQVaIynBr/zsOsF7Pfs65uxMZiF8aWrQ6mYgkh/ysEBt2NST8uP3+DI4tT7lXF5Fz7kOeRJRgnQvXax6BiPjNr1XK4jn7Xd3tfibRgeLEd2J5RGMEIpIsOq8a6uhwBAKJK4IZz1VDS/fY9LKZveZRPAmnRWlEJFnkZ4bocBBujXTVHkqEeLqGhnV7GADmAAWeRZRgnctUqmtIRPzWvd5QUiUCovWFHGBEu4Q2AZd4GVQihVvaAbUIRMR/3ctMlBUl7rjxdA1NTEQgfunsGspJD/ociYikuvys6Ck50ZPK4ukaygSuBI4j2jJYCPzBOdfscWwJUd/cRnZ6kLSglqkUEX+92zWU2Otx4ukPuQeoB34be/w/wL3AuV4FlUjhFhWcE5Hk4NeaBPGcAWc456Z3e/ycma32KqBEq9fqZCKSJPwqRR1Pf8gyM3t/5wMzO4ohVH6iviVCbgJH50VEepObnkbASPiksnh+Cs8B/mNmb8cejwPWmtkKwDnnZnoWXQKEm9vIV9eQiCSBQMDIzwol32AxcKrnUfiovjnCiPxMv8MQEQGik8qSLhE457YkIhC/aFEaEUkmBT60CFL+msnoMpUaIxCR5ODHKmUpnQg6Ohzh1ojKS4hI0sjPSvwqZSmdCBpaIzmjhcEAAArESURBVDgHeeoaEpEkEe0aSuyEspROBFqdTESSTb4PaxKkdCLQojQikmzyM0O0RjpobmtP2DGVCECDxSKSNAp8mF2c4okgthaBxghEJEn4UW/I00RgZqea2VozW29m1/bw/EVmtsvM3ojdLvUynj1pjEBEko0f9YY8OwOaWRC4BTgZ2AYsNrPHnXN7Fqz7q3Pui17F0ZdwsxKBiCSXodY1dCSw3jm30TnXCjwAnOnh8fZZvdYrFpEk09U1lMArh7xMBGOArd0eb4tt29PZZrbczB42s7E9vZGZXW5mS8xsya5duwYswPqWCGaQk65EICLJobMIZm3j0EgE8fg7MCFWwfRZ4O6ednLOzXPOzXXOzS0tLR2wg4ebI9GyrwEbsPcUETkQ744RJG5SmZeJoBzo/gu/LLati3OuyjnXEnv4R6IlrxOmvrlNcwhEJKmEggFy0oNDpmtoMTDFzCaaWTpwPvB49x3MbFS3hx8H3vIwnr2o8qiIJKNEr0ng2VnQORcxsy8CTwNB4A7n3Coz+wGwxDn3OPBlM/s4EAF2Axd5FU9PtF6xiCSjRJei9vQs6JxbACzYY9v13e5/C/iWlzH0pa450jVCLyKSLPITXIra78FiX9U3t6nyqIgknUSvUpbSiaCuqa1rhF5EJFkkenGalE0Ezjlqm9rUNSQiSacgK0Rd89C4fDSpNbW109bulAhEJOnkZ6URbokQae9IyPFSNhF09r8pEYhIsnm3zERiWgVKBEoEIpJkEl2KOmUTQV1s+rYSgYgkm/zMxFYgTdlEoBaBiCSrgmwlgoRQIhCRZJXoUtRKBEoEIpJk1DWUIJ3/wKo+KiLJJtGrlKVsIqhraiMvM42g1iIQkSSTGQqQHgx0XdTitZRNBJpVLCLJyszIz0pTi8BrSgQikswKs9OpbmhNyLGUCEREklBJbjpVDS397zgAUjYRVDe2UpitRCAiyakkN4PKsFoEntrd0MqwnHS/wxAR6VE0EahF4JlIewc1jW0My8nwOxQRkR6V5KZT3xyhua3d82OlZCKoboyOxBerRSAiSao4N/pDtSoBA8Ypmgii/7DqGhKRZFXSmQgS0D2UkomgKjYAoxaBiCSrktzo+SkR4wQpmQh2x5paw3KVCEQkOXW2CCrr1TXkid2xa3PVNSQiyaorESRgLkFKJoLOwZeibCUCEUlOWelBctKDahF4ZXdDKwVZIULBlPz4IjJIFOdmJGR2cUqeCavCrRRrfEBEklxpXgYVdUoEnthR18zI/Ey/wxAR6dOogkx21DV7fpzUTAS1SgQikvzGFGZRXtOEc87T46RcIujocOysa2ZkgRKBiCS30YVZtEY6PJ9dnHKJoLKhhUiHY5QSgYgkudGFWQC8U9Pk6XFSLhHsqI32t41Q15CIJLnRhdHzlBLBAOtMBKMKsnyORESkb2NiLYLyGm8HjFMuEZTHMuuoQrUIRCS5FWSFyE4Psq260dPjpFwi2LirgbzMNBWcE5GkZ2ZMLMlhw64GT4+TcolgU2UDk0pyMDO/QxER6dfUEXms21nv6TFSLhFs3BVmUmmu32GIiMRl6og8ttc2U9vU5tkxUioR1Da18U5tMweV5vgdiohIXKaOiP5w9bJV4GkiMLNTzWytma03s2t7eD7DzP4ae/5VM5vgZTyvv10NwOxxRV4eRkRkwBw6pgCApVuqPTuGZ4nAzILALcBHgenABWY2fY/dLgGqnXOTgV8BP/MqHoCX1lWSFjBmji308jAiIgNmeH4mU4bn8vzaXZ4dw8sWwZHAeufcRudcK/AAcOYe+5wJ3B27/zBwknk0ivvvNTu579UtnHTIcHIz0rw4hIiIJz522GgWbazithc2ePL+Xp4RxwBbuz3eBhzV2z7OuYiZ1QLFQGX3nczscuBygHHjxu1XMBOKczh6UjE3fOx9+/V6ERG/XHr8RHbVtzCq0JuJsIPip7Fzbh4wD2Du3Ln7VYZvUmkud1585IDGJSKSCNnpafzwrBmevb+XXUPlwNhuj8ti23rcx8zSgAKgysOYRERkD14mgsXAFDObaGbpwPnA43vs8zjwudj9c4B/O68Lb4uIyHt41jUU6/P/IvA0EATucM6tMrMfAEucc48DfwLuNbP1wG6iyUJERBLI0zEC59wCYMEe267vdr8ZONfLGEREpG8pNbNYRET2pkQgIpLilAhERFKcEoGISIqzwXa1ppntArbs58tL2GPWcgrQZ04N+syp4UA+83jnXGlPTwy6RHAgzGyJc26u33Ekkj5zatBnTg1efWZ1DYmIpDglAhGRFJdqiWCe3wH4QJ85NegzpwZPPnNKjRGIiMjeUq1FICIie1AiEBFJcSmTCMzsVDNba2brzexav+NJBDPbbGYrzOwNM1vidzxeMLM7zKzCzFZ22zbMzJ41s3Wxv0V+xjjQevnM3zOz8th3/YaZneZnjAPJzMaa2XNmttrMVpnZ/8W2D9nvuY/P7Mn3nBJjBGYWBP4LnEx0yczFwAXOudW+BuYxM9sMzHXODdlJN2Z2AhAG7nHOzYht+zmw2zn301jSL3LOfdPPOAdSL5/5e0DYOfcLP2PzgpmNAkY555aZWR6wFDgLuIgh+j338ZnPw4PvOVVaBEcC651zG51zrcADwJk+xyQDwDn3ItG1LLo7E7g7dv9uov+BhoxePvOQ5Zzb7pxbFrtfD7xFdL3zIfs99/GZPZEqiWAMsLXb4214+I+aRBzwjJktNbPL/Q4mgUY457bH7u8ARvgZTAJ90cyWx7qOhkw3SXdmNgE4HHiVFPme9/jM4MH3nCqJIFUd55ybDXwUuCrWpZBSYkufDv3+T7gVOAiYBWwHfulvOAPPzHKBR4CvOOfquj83VL/nHj6zJ99zqiSCcmBst8dlsW1DmnOuPPa3AphPtIssFeyM9bF29rVW+ByP55xzO51z7c65DuB2hth3bWYhoifE+51zj8Y2D+nvuafP7NX3nCqJYDEwxcwmmlk60bWRH/c5Jk+ZWU5skAkzywE+Aqzs+1VDxuPA52L3Pwf8zcdYEqLzhBjzCYbQd21mRnR987ecczd1e2rIfs+9fWavvueUuGoIIHaZ1a+BIHCHc+7HPofkKTObRLQVANG1qf88FD+zmf0FOJFoed6dwA3AY8CDwDiiJcvPc84NmcHVXj7ziUS7CxywGfh8t/7zQc3MjgMWAiuAjtjmbxPtMx+S33Mfn/kCPPieUyYRiIhIz1Kla0hERHqhRCAikuKUCEREUpwSgYhIilMiEBFJcUoEIv0ws0IzuzJ2f7SZPex3TCIDSZePivQjVuvlic5KnyJDTZrfAYgMAj8FDjKzN4B1wCHOuRlmdhHRipc5wBTgF0A68FmgBTjNObfbzA4CbgFKgUbgMufcmsR/DJGeqWtIpH/XAhucc7OAb+zx3Azgk8ARwI+BRufc4cAi4MLYPvOALznn5gBXA79PSNQicVKLQOTAPBerF19vZrXA32PbVwAzY9UjjwEeipaPASAj8WGK9E6JQOTAtHS739HtcQfR/18BoCbWmhBJSuoaEulfPZC3Py+M1ZDfZGbnQrSqpJkdNpDBiRwoJQKRfjjnqoCXY4vF37gfb/Fp4BIzexNYhZZJlSSjy0dFRFKcWgQiIilOiUBEJMUpEYiIpDglAhGRFKdEICKS4pQIRERSnBKBiEiK+/+eOnUYo0zCZwAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plot_pulse(oct_result.optimized_controls[0], tlist)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "This slight shift of qubit and TLS out of resonance delays the population\n", "oscillations between qubit and TLS ground state such that the qubit ground\n", "state is maximally populated at final time $T$." ] }, { "cell_type": "code", "execution_count": 22, "metadata": { "ExecuteTime": { "end_time": "2019-02-12T04:47:40.231722Z", "start_time": "2019-02-12T04:47:39.777301Z" }, "attributes": { "classes": [], "id": "", "n": "22" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEGCAYAAAB/+QKOAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO3dd3hUZfbA8e9JT0ggpFFCCYTQpBOqilgoVhRRQUEEFcu6Koq76lpQ15WfXdeyWFiKCgIqIkrRVVQEgdA70kloIQkhIT3z/v64AwYIECCTO+V8nmeezC0zc4YJOXPfdsQYg1JKKXUiP7sDUEop5Z40QSillCqXJgillFLl0gShlFKqXJoglFJKlSvA7gAqS0xMjElISLA7DKWU8ijLli07aIyJLe+Y1ySIhIQEUlJS7A5DKaU8iojsPNUxbWJSSilVLk0QSimlyqUJQimlVLk0QSillCqXJgillFLl0gShlFKqXJoglFJKlcul8yBEpC/wFuAPfGSMGXPC8YbAOCAWyAQGG2NSncdKgTXOU3cZY65zSZDF+fDLqxAYAoFh1i28FlSvC9XjoVoMiLjkpZVSyp25LEGIiD/wLtALSAWWishMY8z6Mqe9Ckw0xkwQkcuAl4AhzmP5xph2rorvmMIcWPA6GEf5x0NqQNwFENcC6rSBhhdCdBNNGkqpMzLGkJ1fTMaRIjJyi8g8UkjmkWLyikrILyolr7iU/KJSCopLMcb6s2L9aRFEwF8Efz/BTwR/P/DzkxP2WX+HYiOCGdS5QaXH78oriM7AFmPMNgARmQL0A8omiJbAI877PwEzXBhP+cLj4JlMKC2G4jwoOgK5++DwHshOhfRNcGADrJkGKR9bj6kWBwkXQrOrIKk3hEZWedhKKfsVFJeSmpVP2qF8UrPySMvKP7adlpXPwdxCShynLsoW6C+EBvoTEuiPCBgDBusnGEod1s1hsO4bg8P5s2ytt/YNIj0uQcQDu8tspwJdTjhnFdAfqxnqBiBCRKKNMRlAiIikACXAGGOM65KHCAQEWbfQSKgRD/Edjz/HGMjYAjsXws7fYNvPsO4r8AuERj2g7UBocS0EhrosTKVU1XI4DPsOF7ArM49dmXnsdv607lsJoKwAP6FuZCjxkaFclBRDXEQw0eHBxIQHEVUtiOhqwURVCyIs2J/QQH8C/c+9G9gYK3kABJzH85yO3WsxjQLeEZE7gF+ANKDUeayhMSZNRBoDP4rIGmPM1rIPFpERwAiABg0qP3seRwRikqxbx6HgcEDaMtgwE9Z/DV/eDSGR0OYW6HQXxDZ1bTxKqfOWW1jCvuwC9h8uYF92AfucP3dnWUkgNTOfotI/m5/9/YS6kSE0iArjihZx1KsZSr2aYcTXDKVezVDiIkKONfu4mogQ4O/a1xJX1aQWkW7AaGNMH+f2EwDGmJdOcX44sNEYU6+cY+OBWcaY6ad6veTkZGPbYn0OB+z4FZZPtBJGaRE0uxouehjqd7YnJqVsVuowHCkq4UihdcspsNrdS5zNJtZPByUOQ0mpwWGsphTjbD45uu0wxtnsYjWvOI41w5R9jPO8crZLHIbDBcVk55dwKK+I7PxiDuUVk3WkiJzCkpPijggJoEFU2J+36D/v140MPa9v/e5IRJYZY5LLO+bKK4ilQJKINMK6MhgI3HpCYDFApjHGATyBNaIJEakJ5BljCp3nXAi87MJYz4+fHzS+xLodOQhLPrBuH39rdWpf/gw06Gp3lEpVmsKSUnZn5rMz4wjbDx5hz6EC0nMLSc8pID2nkPScQg4XnPzHt6r5ifWtPyIkkMjQQGqEBRJVLYjGMdWIDAuiVvUQatcItn5WD6F2jRDCguxuWHEfLvuXMMaUiMgDwFysYa7jjDHrROR5IMUYMxPoCbwkIgariekvzoe3AMaKiANrrsaYE0Y/ua9qMXDpk9D9QVgxCRa8AeP6QNMrrURRq6XdESp1VvZlF7BuTzZr0w6zbk82G/YdJi0rn7J9r2FB/sRFBBMbEUzz2tW5qEkQkWFBRIQEUC04gHDnLTTIn0B/wd/PjwA/axROgJ/gd/Snc3Sgn5/gJ+An1mge4c9tPxHED4Qy22XOtfZZTTDq/Lisiamq2drEdDpFR2Dxf2DBW1B42Oq/uPxZCIuyOzKlTmKMYUdGHou3ZbB4eyZLtmeSdigfsLrhGkVXo0Xd6iTGhtMoJoyG0dVoFF2NyLBA/YPsoU7XxKQJoqrkZcIvr8DisdbciitGQ/shVvOUUjYqKnGweHsGP6zfzw8bDhxLCDHhQXRuFEWnhChax9egRZ3qVAvW5hdvownCnexfD9+NsobKxidDv3esSXhKVaGSUge/bjnI1yvS+GHDAXILSwgJ9OOiJjH0bBZH18bRJMZW06sCH2BXJ7UqT62WcMe3sPpzmPskjO3xZ5+Fn7/d0Skvt25PNtNSUpm1eg8Hc4uoERrINW3q0KtlLbonxhAapL+D6k+aIOwgYk2sS7wcvh0JP4yGjd/C9e9b8yyUqkRFJQ5mr93LpEU7SdmZRVCAH1e0iOP6dvH0bBZHUIA2c6ryaYKwU3gs3DwJ1n4B3z4K/7kYrn4V2t2maz2p85adX8yEhTuYuGgnB3MLaRgdxlNXt2BAx3pEhgXZHZ7yAJog7CYCrQdAwkXWbOyv/wLbf4WrX4PgcLujUx4oI7eQjxdsZ+KineQWlnBJ01juuDCBS5Ji8auiWb7KO2iCcBcRtWHIDGvp8Z/HQFoK3DQeare2OzLlIbLzi3lv/hYmLNxBYYmDq1rV4f5LE7mgbg27Q1MeShOEO/Hzh55/h4bd4Yu74KNecP270OpGuyNTbqywpJRJi3byzk9byM4vpl/bujxwWRJN4vQKVJ0fTRDuqNHFcO+v8PkQmD4c9q+DS5/SORPqOMYYZq/dx7++20BqVj4XJ8Xw+JXN9YpBVRpNEO4qPA6GzrTmTPz6mjV/ov8HEFLd7siUG9hx8AjPzFzHL5vTaV47gkl3dubipFi7w1JeRhOEOwsIhmvfhlqtYc7j1ppOt0236lUon1RQXMrYn7fx7vwtBPn78ey1LRnStaHL6gEo36YJwt2JQJcR1vyIz4fAx72sJKGL/vmcVbsP8ei0VWw5kMs1berw9DUtqVU9xO6wlBfTrx2eIvFSGD7bqp09rq81FFb5hKISB6/N20T/9xeSW1DC+GGdeOfWDpoclMtpgvAktVvDnd9D9TrwSX9Y+6XdESkX27jvMP3e/Y1//7iF69vFM3dkD3o2i7M7LOUjtInJ00TWh+FzYPKt1ginoiPQYYjdUalKZoxh8pLdjP5mHdVDAvnw9mR6taxld1jKx2iC8EShNWHwF/D5bTDzASjOt/oplFfIKSjmya/W8s2qPVycFMMbt7QjJjzY7rCUD9IE4amCwmDQFJg2DGY/BsV5Vg1s5dHW7cnmL58uZ1dmHo/1acZ9lyTq8hjKNpogPFlAMNw8Ab66B3541rqS6Pm4LvTnob5emcbfpq+mZlgQU0Z0o3MjrTqo7KUJwtP5B0L/DyEg1FrDSfys5TqUxyh1GF6Zu4n//LyVzglRvDe4gzYpKbegCcIb+PnDdf8GDMz/l5U0Ln7E7qhUBWTnF/PQlBXM35TObV0a8Oy1F2h9BuU2NEF4Cz8/K0mUFsH/ngP/IOj+gN1RqdPYlZHHHeOXsCsjj39e34rBXRvaHZJSx9EE4U38/OH6/0BpMcz7h5UkdHSTW1qdeojh45dS4jB8elcXujSOtjskpU6iCcLb+AfAjR9ZSWL2YxAQBB3vsDsqVcZPGw9w/6fLiQ4PYvywzrost3Jb2tjpjfwD4ab/QpNe8M3DsG6G3REpp8lLdnHXxBQS46rx5f3dNTkot6YJwlsFBMPNE6F+F6uU6daf7I7Ipxlj+Pf//uCJL9dwUZMYPh/RjbgIXUtJuTdNEN4sKAxunQLRSTDlNkhbZndEPskYw5g5G3nt+83c0D6ej4YmUy1YW3eV+9ME4e2OLstRLRo+GQDpm+2OyKc4HIZnZ65j7M/buK1LA167qS2BWrtBeQj9TfUF1evAkBnWKKdJ10N2qt0R+YSSUgd/+2I1ExftZESPxvzz+la6bIbyKJogfEV0Igz+EgpzYFJ/yM+yOyKvVlLq4KHPVzJ9WSoPX5HEE1c2R3QJFOVhNEH4kjptYNBkyNoOUwZDSaHdEXmlUodh5NRVfLt6L09c2ZyHr2iqyUF5JE0QvibhIuj3HuxcADPuB4fD7oi8SqnD8Ni0VXyzag9/79ucey5JtDskpc6ZDqXwRW1uguzd1pIckQ3gimftjsgrOByGJ75czZcr0ni0V1Pu66nJQXk2TRC+6qKRcGgXLHjdqlKXPNzuiDyaMYanvl7L1JRUHrysCX+9PMnukJQ6b5ogfJUIXPUqHE6Dbx+F6vWgaW+7o/JIxhj++e0GPlu8i/t6JjKyV1O7Q1KqUmgfhC/zD4AB/4XarWHaHbBnpd0ReaT35m/l4wXbuaN7An/r00w7pJXX0ATh64LD4dapEBYFn92scyTO0uQlu3hl7iaub1eXZ65pqclBeRWXJggR6Ssim0Rki4g8Xs7xhiLyPxFZLSLzRaRemWNDReQP522oK+P0eRG14bZpVsnSyQOhMNfuiDzCnLV7+cdXa+jZLJZXbmqrk+CU13FZghARf+Bd4EqgJTBIRFqecNqrwERjTBvgeeAl52OjgGeBLkBn4FkRqemqWBUQ18Jqbtq/zqpxrcNfT2vhloM8OHkl7epH8t5tHXT5DOWVXPlb3RnYYozZZowpAqYA/U44pyXwo/P+T2WO9wG+N8ZkGmOygO+Bvi6MVQEkXQF9XoKNs6whsKpca9OyuXtiCgkxYYy7oxNhQTrWQ3knVyaIeGB3me1U576yVgH9nfdvACJEJLqCj0VERohIioikpKenV1rgPq3LPZB8J/z2Jqz41O5o3M6eQ/kMH7+UGqGBTBzehciwILtDUspl7L4uHgVcIiIrgEuANKC0og82xnxgjEk2xiTHxsa6KkbfIgJX/h807gnfPAQ7frM7IreRU1DM8PFLyS8q5b/DOlO7htZzUN7NlQkiDahfZruec98xxpg9xpj+xpj2wD+c+w5V5LHKhfwD4aYJUDMBPh8Mmdvsjsh2xaUO7v90OVsO5PLe4A40qx1hd0hKuZwrE8RSIElEGolIEDAQmFn2BBGJEZGjMTwBjHPenwv0FpGazs7p3s59qqqERsKtnwMGPhsIBdl2R2QbYwzPfL2OX/84yIs3tOLiJL1aVb7BZQnCGFMCPID1h30DMNUYs05EnheR65yn9QQ2ichmoBbwovOxmcALWElmKfC8c5+qStGJcPMkyNxqTaQrLbE7IluM/WUbk5fs4i+XJnJLpwZ2h6NUlRFjjN0xVIrk5GSTkpJidxjeaflEmPlX6DwCrnrF7miq1Ow1e7nv0+Vc27Yub93STuc6KK8jIsuMMcnlHdPxeerMOtwOBzfDwn9DbHPodKfdEVWJdXuyeWTqKjo2rMkrA9poclA+x+5RTMpTXPEcNO0L3z0G2362OxqXO5hbyIiJy4gMC+Q/gzsSEuhvd0hKVTlNEKpi/Pyh/4cQ0xSm3g4ZW+2OyGWKShzc/8lyDuYW8sGQZGIjgu0OSSlbaIJQFRdSHW6dYiWLz26B/EN2R+QSo79Zx5Idmbw8oA2t69WwOxylbKMJQp2dmglwyyeQtQOmD/O6kU2Tft95rK5Dv3YnTd5XyqdoglBnr2F3uOYN2PojzH3S7mgqzaKtGTw3cx2XNY9jVO9mdoejlO10FJM6Nx2GQPpGWPQOxDX3+JKlaYfy+ctny2kYHcabA9vhryOWlNIrCHUeej0PSb2tkU3bf7E7mnNWWFLK/Z8up6jEwQe3J1M9JNDukJRyC5og1Lnz84cbP4boJvD5EI8d2fTCrPWs2n2IV29qQ2JsuN3hKOU2NEGo8xNSHQZNAfGzqtF52MimL5en8snvu7inR2P6tqpjdzhKuRVNEOr8RTWCWyZZq7560MimDXsP8+RXa+jSKIrH+mintFIn0gShKkfCRX+ObJr3D7ujOaPDBcXc98kyqocE8u9b2xOgJUOVOomOYlKVp8PtcGAj/P6utWZT8jC7IyqXMYZRU1eRmpXP5BFdiYvQwj9KlUe/NqnK1fsFaNILvhvltiObxv6yjXnr9/PEVS3olBBldzhKuS1NEKpy+fnDgI8hKtEt12xasj2Tl+ds5Oo2dRh+YYLd4Sjl1jRBqMoXUsNaswlg8iC3qUaXdaSIh6asoEFUGGP6t0ZEJ8MpdTqaIJRrRDX+sxrd9OG2j2wyxjBq2ioycot459YOROhkOKXOSBOEcp1GF8PVr8GWH+D7p20N5eMF2/nfxgM8eVVzWsXrCq1KVYSOYlKu1fEO58im9yC2mbVdxVbtPsT/zdlI75a1GNo9ocpfXylPpVcQyvV6/xMSL4dvH4Xtv1bpSx8uKOaBycuJiwjh5QFttN9BqbOgCUK5nn8A3PRfq19i6hBrxnUVMMbw+Ber2XOogLcHtSMyLKhKXlcpb6EJQlWNkBrWmk0Anw2skpFNny7exXdr9jGqdzM6NtT5DkqdrQonCBHxF5G6ItLg6M2VgSkvFJ0IN090jmy6ExylLnupDXsP8/ys9fRoGss9PRq77HWU8mYVShAi8ldgP/A98K3zNsuFcSlv1agHXPUKbPke5rlmZFN+USl/nbyCGqGBvH5zW/y0+I9S56Sio5geApoZYzJcGYzyEcnDIX2TtWZTXHNrDadK9K/vNrDlQC6f3NmFmPDgSn1upXxJRZuYdgPuMR1WeYfeL1ojm2Y9Ajt+q7Sn/d+G/Uz6fSd3XdSIi5JiKu15lfJFFU0Q24D5IvKEiDxy9ObKwJSX8w+AAeOgZgJ8Phgyt5/3Ux7IKeBv01fTvHYEj/XV+g5Kna+KJohdWP0PQUBEmZtS5y40Em79HIzDqkZXcPicn8oYw2PTVpNbWMLbg9oTHOBfiYEq5Zsq1AdhjHkOQETCndu5rgxK+ZCjI5s+6Q9f3GkNhfU7+z/uExbu4OfN6Tzf7wKa1tLvLkpVhoqOYmolIiuAdcA6EVkmIhe4NjTlMxpfAle+DH/Mg++fOeuHb96fw79mb+TSZrEM6drQBQEq5ZsqOorpA+ARY8xPACLSE/gQ6O6iuJSv6XQnpG+ERe9YazZVcGRTQXEpD05eQfWQAF4e0FaX0lCqElW0D6La0eQAYIyZD1RzSUTKd/V5CRpfelYjm16Zu4mN+3J4ZUBbYiN0SKtSlanCo5hE5GkRSXDensIa2aRU5Tm6ZlPNhtaaTVk7Tnv6L5vT+XjBdoZ2a8ilzeOqJkalfEhFE8RwIBb40nmLde5TqnKF1oRBn4OjxLlmU/kjmzKPFPHotFUkxYXzxFUtqjhIpXxDhRKEMSbLGPOgMaaD8/aQMSbL1cEpHxXTBG6aAAc3O6vRFR932BjDk1+uITuvmLcGtickUIe0KuUKp00QIvKm8+c3IjLzxFvVhKh8UuKlcM3r1ppNsx4GY44dmrEyjTnr9vFI76a0rFvdxiCV8m5nGsU0yfnz1XN5chHpC7wF+AMfGWPGnHC8ATABiHSe87gx5jsRSQA2AJucp/5ujLn3XGJQHqzjHXB4D/z8fxBRBy57ir3Z+Tzz9TqSG9bk7ot1lValXOm0CcIYs8x5t50x5q2yx0TkIeDnUz1WRPyBd4FeQCqwVERmGmPWlzntKWCqMeZ9EWkJfAckOI9tNca0O5s3o7xQzycgZy/88gomog5/W92KklLDqze1xV9XaVXKpSraST20nH13nOExnYEtxphtxpgiYArQ74RzDHC0jaAGsKeC8ShfIQJXvwFN+2K+HUXo1tk8eXULEmJ0lLVSrnbaKwgRGQTcCjQ6oc8hAsg8w3PHY60Ce1Qq0OWEc0YD85z1JqoBV5Q51sg5e/sw8JQx5qRixiIyAhgB0KCB1i/yWv4B7Lz0HQ5t6ss7we8SWKcPoDOmlXK1M/VBLAT2AjHAa2X25wCrK+H1BwHjjTGviUg3YJKItHK+ZgNjTIaIdARmiMgFxpjjxjwaYz7AmuVNcnKyOfHJlXcodRgenfEHB+RxfqzxL2TyQBg+16oloZRymdM2MRljdhpj5htjuhljfi5zW26MKTnDc6cB9cts13PuK+tOYKrztRYBIUCMMabwaHEiZz/IVqBpxd+W8iYf/bqNlJ1ZjOzXjYChX4F/EHxyI2Sf+OuklKpMFV2sr6uILBWRXBEpEpFSETnT2sxLgSQRaSQiQcBA4MShsbuAy52v0QIrQaSLSKyzkxsRaQwkoTO3fdKmfTm8Nm8zfS+ozfXt4q36EYOnQ0E2TLoejhy0O0SlvFZFO6nfwWoO+gMIBe7CGqF0Ss4rjAeAuVhDVqcaY9aJyPMicp3ztEeBu0VkFTAZuMMYY4AewGoRWQlMB+41xpypz0N5maISB49MXUlESAAv3tDqz4X46rS16kgc2gWTbrCShVK+yFEKW36AtV+45OnFmDM33YtIijEmWURWG2PaOPetMMa0d0lU5yA5OdmkpKTYHYaqRK/P28TbP25h7JCO9Lmg9sknbJ4HUwZBvU4w+EsICqv6IJWyw4ENsPIzWD0VcvdBrVZw37mV7hWRZcaY5PKOVXS57zxnM9FKEXkZqxO5olcfSp21VbsP8e78rfTvEF9+cgBo2hv6f2gVGvp8MAyaDAG6oqvyUkcOwprpsOoz2LsK/AIgqTe0HQhN+7rkJSuaIIZgzXR+ABiJ1fl8o0siUj6voLiUR6auJC4imGevPUNdqlb9oSgXZv4VvrgLBvzXWhVWKW9QUgR/zIWVk62fjhKo3Qb6joFWAyA81qUvX9GSozudd/OB51wXjlLw8pxNbE0/wid3dqFGaOCZH9DhdijMgblPwjcPwnXvgJ9e4CoPZQzsWW4lhbXTIT8LwmtB1/ug7SCoVXXFPM80UW4N1mznch3tj1CqsizamsG437Zze7eGXJQUU/EHdvuLtTT4z2PAP9Cafa1JQnmSrJ2wZprVr3BwEwSEQPOrraTQ+FJbrozP9IrXVEkUSgE5BcWMmraKhOgwHr/yHCbB9XwcSotgwevWtiYJ5e6OZMC6L62+hd2/W/vqd4Vr34KW10NopK3hnWmxvp2nO65UZfrnrA3szc5n2r3dCQs6h29LInD5M9b9Ba8DAle/rklCuZfcA7DpO9gwC7b9ZPUrxLawfndbDbAqKrqJCv0vFJEc/mxqCgICgSPGGF2MX1WKHzfu5/OU3dzXM5GODWue+xMdSxIGFrxh7dMkoeyWsRU2zoKN38LuJYCByIZW02jrm6F2K7sjLFdFO6kjjt4Xa7ZSP6Crq4JSviXrSBF//2INzWtH8PAVSef/hCJw+bPW/QVvWNtXvaZJQlWdwhzY8Zt1hbD1R6s6IliTPC990upbiGtp/W66sbO+jnfOdJ4hIs8Cj1d+SMrXPPX1Wg7lFTFhWGeCAyqpfOiJScJRAte8CX5anlS5QGkxpC2DbfNh60+QlmL9zgWEQMPu0OkuaHYlRHrWqtMVbWLqX2bTD0gGClwSkfIpM1ft4dvVe3msT7PKLx96NEn4BcAvr0BhLtwwFgKCKvd1lO8xxroqOJoQdiyAohxAoG476P4gNO4J9btAYIi9sZ6Hil5BXFvmfgmwg5OL/yh1VvYfLuDpGWtp3yCSe3q4qHyoCFz2FARXh++ftibV3TRBl+VQZy9nP2z/2UoI2+ZDjrO+Wc1G0HqAVUc94WIIi7I1zMpU0T6IYa4ORPkWYwx//2I1hSWlvHZTWwL8Xdw/cOGDEFIdvnkYPh0Ag6ZY20qdSlEe7FpoJYStP8GBddb+0JrQ6BIrITTuaa0w7KUq2sTUGHgLq2PaAIuAkcYYXYJbnZMpS3czf1M6o69tSePY8Kp50Y53QFA4fHUPTLjWWuCvWnTVvLZyf8bAvtVWp/LWH2HX79a8Gv8gaNDVaq5MvBRqt/WZAQ8VbWL6DGt57xuc2wOxluc+sYSoUme0OzOPf85aT7fG0dzeLaFqX7z1AAiOgKm3w8e9rNoSUS5q3lLur6QQtv9izUvYNOfPZqO4C6DzCCshNOjus02SFU0QYcaYSWW2PxGRx1wRkPJuDofh0WmrEBFeuakNfn42DPNr2gdunwmTB8JHvazaEvXKXe1YeaOCw7BptjUvYeuPVr9UYDVochk0fQqaXAERteyO0i1UNEHMFpHHgSlYTUy3AN+JSBSAFvNRFTXut+0s2Z7JywPaUK+mjd/KGnSBO7+HT2+E8dfAjR9BC11ZxmsV58Mf86wlLf6YByUFEFEH2twMza6yOpc9eLSRq1S0YND20xw2xhjbr9G1YJD723Igh6veXkCPpBg+vD35zwpxdspNt64k0pZZSyh3vdfuiFRlcThgx69WYZ2N31rDUKvFwQU3QKsbrUJTPtKXcDrnXTDIGNOockNSvqa41MEjU1dRLciff/Vv7R7JAaz19Id+Y9WSmPN3SN8IV76scyU82eE9sPJTWPEJZO2AkBpwwfVW/1PDi7ReyFmo6CimQOA+rFrRAPOBscaYYhfFpbzMuz9tYXVqNu/d1oG4CDe7lA8Kg1smwY8vWLOu0zfCzRMhPM7uyFRFlZbA5jmwfCJs+R6MAxr1gEufspoOA0PtjtAjVTSVvo+1QN97zu0hzn13uSIo5V3WpGbzzo9buL5dXa5qXcfucMrn5w9XjLZq+379AHzQEwZ+CnXdpuy6Kk9eJiyfAEs+gsOpVr/CRSOh/WAdnVYJKpogOhlj2pbZ/lFEVrkiIOVdjpYPjQkP5rnr3HPFyuO0HgAxSTDlNhjX11q/qd0gu6NSJ9q3FpaMtYrrlBRYVwtX/p9Vm1mbkCpNRf8lS0Uk0RizFY5NnCt1XVjKW7w2bxN/HMhlwvDO1AirQPlQd1CnLYyYD9PugBn3Wh2dV70CQdVsDszHlZbA5tmweKz1mQSEQtuB0PkeqNXS7ui8UkUTxGPATyJydOZ0AqDLb6jT+n1bBh8t2M7grg24pKlri6tXumoxMGQG/Px/1kJ/acvgpvEQ18LuyHxPYQ4snwSL34dDu6BGA+j1gtWM5EXrHrmjiiaI34CxwIiLnZMAABcVSURBVOXAIWAu1nIbSpUrt7CEUdNW0SAqjCev8tA/qv4BcNk/rOWavxwBH1wKV46BDkPdfh1/r5CdBov/A8smQGG2NaO5z7+seQu6bHuVqGiCmAgcBl5wbt8KTAJuckVQyvO9+O160g7lM+2ebudWPtSdJF4K9y6Ar0bANw/Bxu/gurchorbdkXmnvatg4TtWrWZjoGU/6P4AxHe0OzKfU9H/ua2MMWUb+X4SkfWuCEh5vh837mfykt3ce0kiyQle0gQQUQsGfwVLPoAfRsO7XeCqV61Obb2aOH8OhzU8deG/rf6FoHCrb6HLPW5Vo9nXVDRBLBeRrsaY3wFEpAug05bVScqWDx3ZqxLKh7oTPz9rpnWTK6zO6y/vgg1fWxPrqte1OzrPVFwAq6fAonetAjzV463+hY5DrQluylYVTRAdgYUissu53QDYJCJrsJbaaOOS6JTHedpZPnT8sE6VVz7U3cQ0geFzYeHb8NNLsHW+VWe48wgdYllRRw7C0o9gyYeQdxBqt4H+H1rLYPh7yGg3H1DR3+a+Lo1CeYWZq/Ywy1k+9IK6Xv7tz8/fmpDVsh989xjMfcJa8+ea16F+Z7ujc197VlhJYc10KC2EpD5W/0LCxdpU54YquhbTTlcHojxblZQPdUdRjeG26bBhJsx+3KoxcUF/uPxpncl7VEmR9e+zeCykLrGW1m4/2OpfiG1md3TqNPR6WJ23Ki8f6m5ErCuJxMvgt7dh0Tuw4RvodCf0+JvvVq07vMcaorrsv5C730qYfcdA20EQGml3dKoCNEGo83a0fOhz111QdeVD3VFwhDVvInk4zH/JGvG04hNru9sDvlGEprS4zKJ5P1iL5jXpZV0tJF6uy2t7mArVg/AEWg/CHrsy8uj71i+0bxDJpOFd7KkQ567SN8HPL1vj+f0CocMQ6P5X7yxyn74JVkyCVVPgSLq1aF67W3XRPA9wunoQmiDUOSspdXDz2EX8cSCXOQ/3ID5Sl1QuV8ZWaxnxVVPAUWKVPO10l+d/o87aAWu/tG7714BfgLVYXofbrfemI7o8wnkXDFKqPO/+tJXluw7x9qD2mhxOJzoR+r0DPZ+w2uOXTYDNA6wriXa3WdXNohPtjvLMjIGMLVYT0roZkOb8Qlavk9W30OpGraHhZVx6BSEifYG3AH/gI2PMmBOONwAmAJHOcx43xnznPPYEcCfWqrEPGmPmnu619Aqiaq3YlcWA/yzi2jZ1eHOg1kw4KyVFsPEbWDoOdi6w9sV3tP7ANu3rXsmiOB92/W7Vcd48BzKd63XWbg2tBljzFnSms0ezpYlJRPyBzUAvIBVYCgwyxqwvc84HwApjzPsi0hL4zhiT4Lw/GegM1AV+AJoaY065xLgmiKpzpLCEq9/+leJSw+yHL6Z6iE5sOmfZqVYTzZppsG+1tS+qsdWxm3gp1OtctaOgCrIhNQV2LoSdv1mr2JYWgX+wVXOhaR/rFtmg6mJSLmVXE1NnYIsxZpsziClAP6DsGk4GqO68XwPY47zfD5hijCkEtovIFufz6QqybuCFWevZmZnHlLu7anI4XzXqwYUPWrfM7dbInz/mWaOAloy1zolKtCbf1boAYppBbFNryevz6b8ozLX6ELK2WyVW9662ElTWDuu4+Ft1MbrcY9VxbnSx1sPwQa5MEPHA7jLbqUCXE84ZDcwTkb8C1YAryjz29xMeG3/iC4jICGAEQIMG+o2mKsxZu48pS3dzf89EujT20fH9rhLVCDrfbd2K861Zx7uXQOpS2PI/WDX5z3P9g61hs+G1rVVlQ2tCYJhVezkw1BpeWlpk3UoKrdKceRnWLWcfHDlw/GvXbGQlhPZDoG47qN/FGrarfJrdndSDgPHGmNdEpBswSUQqXJfSGPMB8AFYTUwuilE5HThcwBNfrqZVfHUevqKp3eF4t8BQqw5Fw+5/7svLtIaTHtxkjYzK2Qe5+6x9BYeshe+K88BRbJ0v/uAfBAFBEBoFYdHWooJ12lrJqGYj62dUIoRULz8O5dNcmSDSgPpltus595V1J851nowxi0QkBIip4GNVFXI4DKOmrya/uJQ3b2lPUIAHD8/0VGFR0LCbdTud0hIQP88eQqvcgit/g5YCSSLSSESCgIHAzBPO2YVVpQ4RaQGEAOnO8waKSLCINAKSgCUujFWdwcRFO/hlczr/uLolTeJ8eLa0J/AP0OSgKoXLriCMMSUi8gBWeVJ/YJwxZp2IPA+kGGNmAo8CH4rISKwO6zuMNaxqnYhMxerQLgH+croRTMq1Nu/P4aXZG7mseRyDu2hfj1K+QmdSq9MqKC6l3zu/cTC3kDkP9yA2ItjukJRSlUhnUqtz9sKs9Wzan8P4YZ00OSjlY7ShUp3S7DV7+XTxLu7p0ZiezXQJBaV8jSYIVa7UrDz+/sVq2tarwaO9taiLUr5IE4Q6SXGpgwcnr8AY+PegDjqkVSkfpX0Q6iRv/rCZ5bsO8e9B7WkQHWZ3OEopm+hXQ3Wc37Yc5L35W7kluT7Xtq1rdzhKKRtpglDHHMwt5OHPV5IYG86z17W0OxyllM20iUkBUOowPDJ1Fdn5xUy6szNhQfqroZSv0ysIBcA7P27hl83pjL72AprX1oXblFKaIBTwy+Z03vzfZvp3iGdQ5/pnfoBSyidogvBxew7l89CUFTSNi+Cf17dCROwOSSnlJjRB+LCiEgcPfLacohIH7w3uoP0OSqnj6F8EH/bS7A0s33WId2/tQGKsLuGtlDqeXkH4qG9X7+W/v+1g2IUJXN2mjt3hKKXckCYIH/TH/hz+Nn0V7RtE8sSVLewORynlpjRB+JjsvGLunphCaJA/792m6ywppU5N+yB8SKnD8OCUFaQdyuezu7tSp0ao3SEppdyYJggf8vLcjfy8OZ1/3dCaTglRdoejlHJz2r7gI75emcbYn7dxW5cG3Kp1pZVSFaAJwgesTcvm71+splNCTZ699gK7w1FKeQhNEF7uwOECRkxMoWZYEO/d1lE7pZVSFaZ9EF4sr6iEuyamkJVXzLR7uxEbEWx3SEopD6IJwkuVOgwjP1/JmrRsPhySTKv4GnaHpJTyMNre4KXGzN7A3HX7eeaallzRspbd4SilPJAmCC806fedfPjrdoZ2a8iwCxvZHY5SykNpgvAy/9uwn9Ez13FZ8zievkbLhiqlzp0mCC+ydEcm93+6nJZ1qvP2oPYE+OvHq5Q6d/oXxEts2HuY4eOXEh8ZyvhhnQgP1vEHSqnzo39FvMCujDxuH7eEakEBTLyzM9HhOpxV+a7i4mJSU1MpKCiwOxS3EhISQr169QgMDKzwYzRBeLgDOQUMGbeY4lIHn93TjXo1w+wOSSlbpaamEhERQUJCgpbQdTLGkJGRQWpqKo0aVXzgijYxebD0nEJu+3AxBw4XMu6OTiTVirA7JKVsV1BQQHR0tCaHMkSE6Ojos76q0gThoQ7mFnLbR7+TmpXPf4d1okODmnaHpJTb0ORwsnP5N9EE4YEyjxQx+KPF7MrM4+M7kunaONrukJRSXkgThIfJPFLErR/+zvaDR/h4aCe6J8bYHZJS6jyNHz+eBx54oNxj3bt3B2DHjh189tlnp3yOCRMmkJSURFJSEhMmTKiUuDRBeJC92fnc9J+Fx5LDhU00OSjl7RYuXAicPkFkZmby3HPPsXjxYpYsWcJzzz1HVlbWeb+2jmLyENvScxny8RIO5xczcXhnumizklJn9Nw361i/53ClPmfLutXPWFflxRdfZMKECcTFxVG/fn06duzIqFGj6NmzJ6+++irJyckcPHiQ5ORkduzYAcDu3bvp2bMnaWlpDB48mGeffRaA8PBwcnNzefzxx9mwYQPt2rVj6NChjBw58tjrzZ07l169ehEVZVWK7NWrF3PmzGHQoEHn9V5dmiBEpC/wFuAPfGSMGXPC8TeAS52bYUCcMSbSeawUWOM8tssYc50rY3Vna9OyGTpuCQCTR3TVlVmVcmPLli1jypQprFy5kpKSEjp06EDHjh3P+LglS5awdu1awsLC6NSpE1dffTXJycnHjo8ZM4ZXX32VWbNmnfTYtLQ06tevf2y7Xr16pKWlnfd7cVmCEBF/4F2gF5AKLBWRmcaY9UfPMcaMLHP+X4H2ZZ4i3xjTzlXxeYoFfxzkvk+WUT00kEl3dqZxbLjdISnlMeyooPjrr79yww03EBZmzUm67rqKfbft1asX0dFWy0D//v1ZsGDBcQnCDq7sg+gMbDHGbDPGFAFTgH6nOX8QMNmF8XiczxbvYuh/lxBfM5Rp93bT5KCUhwsICMDhcACcNCfhxGGoZzMsNT4+nt27dx/bTk1NJT4+/jwitbgyQcQDu8tspzr3nUREGgKNgB/L7A4RkRQR+V1Erj/F40Y4z0lJT0+vrLhtV+owvPjtep78ag0XJ8Uw7d5u1I0MtTsspVQF9OjRgxkzZpCfn09OTg7ffPPNsWMJCQksW7YMgOnTpx/3uO+//57MzEzy8/OZMWMGF1544XHHIyIiyMnJKfc1+/Tpw7x588jKyiIrK4t58+bRp0+f834v7jKKaSAw3RhTWmZfQ2NMMnAr8KaIJJ74IGPMB8aYZGNMcmxsbFXF6lLZ+cXcMynlWD2Hj25PJiKk4munKKXs1aFDB2655Rbatm3LlVdeSadOnY4dGzVqFO+//z7t27fn4MGDxz2uc+fO3HjjjbRp04Ybb7zxpOalNm3a4O/vT9u2bXnjjTeOOxYVFcXTTz9Np06d6NSpE88888yxDuvzIcaY836Scp9YpBsw2hjTx7n9BIAx5qVyzl0B/MUYs/AUzzUemGWMmV7ecYDk5GSTkpJSGaHbZt2ebO7/dDlpWfk8fU1LhnZPsDskpTzOhg0baNGihd1hHDN69GjCw8MZNWqU3aGU+28jIsucX8ZP4soriKVAkog0EpEgrKuEmSeeJCLNgZrAojL7aopIsPN+DHAhsP7Ex3qTaSm76f/eQgqKS/n8nq6aHJRStnPZKCZjTImIPADMxRrmOs4Ys05EngdSjDFHk8VAYIo5/lKmBTBWRBxYSWxM2dFP3uRwQTGjZ67jy+VpdE+M5u1B7YnR5bqV8hqjR4+2O4Rz5tJ5EMaY74DvTtj3zAnbo8t53EKgtStjcweLt2XwyNRV7DtcwIOXJ/HgZU20CpxSym3oTGob5BeV8uYPm/ng1200iApj2r3ddDVWpZTb0QRRxeZvOsDTX69ld2Y+gzrX56mrW1JNy4MqpdyQ/mWqIvsPF/DCrPXMWr2XxrHVmHx3V7ol6npKSin3pQ3eLpZbWMJr8zZxySs/MW/dfkZe0ZTZD12syUEpL5WRkUG7du1o164dtWvXJj4+/tj20eU3ytq0aRM9e/akXbt2tGjRghEjRpT7vK5YzvtM9ArCRQqKS/l86W7e/t8fZBwp4tq2dXmsdzMaRGvNaKW8WXR0NCtXrgROngMRHn7ycjkPPvggI0eOpF8/ayWiNWvWnHTO0eW8U1JSEBE6duzIddddR82aru271ARRyXIKivnk9118vGA7B3ML6dwoio+vakG7+pF2h6aU75n9OOw7+Q/ueandGq4cc+bzKmjv3r3Uq1fv2Hbr1icP4HTVct5nogmikmzcd5gpS3bzxfJUcgpKuDgphvt7tqdr4yitj6uUOqWRI0dy2WWX0b17d3r37s2wYcOIjDz+C6WrlvM+E00Q52FfdgFz1+3j65VpLN91iCB/P/q2qs3dFzemdT2t2aCU7Srxm76rDBs2jD59+jBnzhy+/vprxo4dy6pVqwgOtn/CrCaIs1BQXMryXVn8vi2TX/9IZ8WuQwA0rRXOU1e34MYO9ahZLcjmKJVSnqZu3boMHz6c4cOH06pVK9auXXtckaH4+Hjmz59/bDs1NZWePXu6PC5NEOXIKShmX3YBe7IL+GN/Dpv25bBxn/WzqNSBn0Cr+Bo81qcZfS6oTZM4rdOglDo3c+bM4fLLLycwMJB9+/aRkZFxUi2HPn368OSTTx6rMz1v3jxeeumkdU8rnc8niMwjRdwydhEFJaXkFznIKyohr6j0uHNiwoNpUSeCYRcl0KVRFMkJUVTXJbiVUmcpLy/vuA7pRx55hNTUVB566CFCQkIAeOWVV6hdu/Zxjyu7nDdQact5n4nLlvuuaue63PeRwhJGTVtFSKA/IYH+hAb6E1c9mDo1QqhTI5TGsdV08TylPIi7LfftTs52uW+fv4KoFhzA+4PPXFBcKaV8jc6kVkopVS5NEEopr+MtTeeV6Vz+TTRBKKW8SkhICBkZGZokyjDGkJGRcawjvKJ8vg9CKeVd6tWrR2pqKunp6XaH4lZCQkKOG0FVEZoglFJeJTAwkEaNGtkdhlfQJiallFLl0gShlFKqXJoglFJKlctrZlKLSDqw8zyeIgY4WEnheApfe8++9n5B37OvOJ/33NAYE1veAa9JEOdLRFJONd3cW/nae/a19wv6nn2Fq96zNjEppZQqlyYIpZRS5dIE8acP7A7ABr72nn3t/YK+Z1/hkvesfRBKKaXKpVcQSimlyqUJQimlVLl8PkGISF8R2SQiW0TkcbvjqQoiskNE1ojIShE5+zJ8HkBExonIARFZW2ZflIh8LyJ/OH/WtDPGynaK9zxaRNKcn/VKEbnKzhgrm4jUF5GfRGS9iKwTkYec+73ysz7N+3XJ5+zTfRAi4g9sBnoBqcBSYJAxZr2tgbmYiOwAko0xXjuZSER6ALnARGNMK+e+l4FMY8wY55eBmsaYv9sZZ2U6xXseDeQaY161MzZXEZE6QB1jzHIRiQCWAdcDd+CFn/Vp3u/NuOBz9vUriM7AFmPMNmNMETAF6GdzTKoSGGN+ATJP2N0PmOC8PwHrP5bXOMV79mrGmL3GmOXO+znABiAeL/2sT/N+XcLXE0Q8sLvMdiou/Md2IwaYJyLLRGSE3cFUoVrGmL3O+/uAWnYGU4UeEJHVziYor2hqKY+IJADtgcX4wGd9wvsFF3zOvp4gfNVFxpgOwJXAX5xNEz7FWG2rvtC++j6QCLQD9gKv2RuOa4hIOPAF8LAx5nDZY974WZfzfl3yOft6gkgD6pfZrufc59WMMWnOnweAr7Ca2nzBfmcb7tG23AM2x+Nyxpj9xphSY4wD+BAv/KxFJBDrj+Wnxpgvnbu99rMu7/266nP29QSxFEgSkUYiEgQMBGbaHJNLiUg1Z+cWIlIN6A2sPf2jvMZMYKjz/lDgaxtjqRJH/0g63YCXfdYiIsDHwAZjzOtlDnnlZ32q9+uqz9mnRzEBOIeDvQn4A+OMMS/aHJJLiUhjrKsGsErOfuaN71lEJgM9sZZB3g88C8wApgINsJaGv9kY4zWduqd4zz2xmh0MsAO4p0zbvMcTkYuAX4E1gMO5+0msdnmv+6xP834H4YLP2ecThFJKqfL5ehOTUkqpU9AEoZRSqlyaIJRSSpVLE4RSSqlyaYJQSilVLk0QSp0jEYkUkfud9+uKyHS7Y1KqMukwV6XOkXMtnFlHV05VytsE2B2AUh5sDJAoIiuBP4AWxphWInIH1uqh1YAk4FUgCBgCFAJXGWMyRSQReBeIBfKAu40xG6v+bShVPm1iUurcPQ5sNca0Ax474VgroD/QCXgRyDPGtAcWAbc7z/kA+KsxpiMwCnivSqJWqoL0CkIp1/jJuV5/johkA984968B2jhX4+wOTLOW1wEguOrDVOrUNEEo5RqFZe47ymw7sP7f+QGHnFcfSrklbWJS6tzlABHn8kDnGv7bReQmsFbpFJG2lRmcUudLE4RS58gYkwH8JiJrgVfO4SluA+4UkVXAOrTcrXIzOsxVKaVUufQKQimlVLk0QSillCqXJgillFLl0gShlFKqXJoglFJKlUsThFJKqXJpglBKKVWu/weqitsvy+gxWwAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "optimized_dynamics = oct_result.optimized_objectives[0].mesolve(\n", " tlist, e_ops=[proj_00, proj_01, proj_10, proj_11]\n", ")\n", "\n", "plot_population(optimized_dynamics)" ] } ], "metadata": { "hide_input": false, "jupytext": { "formats": "" }, "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.7.3" }, "toc": { "base_numbering": 1, "nav_menu": {}, "number_sections": true, "sideBar": true, "skip_h1_title": false, "title_cell": "Table of Contents", "title_sidebar": "Contents", "toc_cell": false, "toc_position": {}, "toc_section_display": true, "toc_window_display": false }, "varInspector": { "cols": { "lenName": 16, "lenType": 16, "lenVar": 40 }, "kernels_config": { "python": { "delete_cmd_postfix": "", "delete_cmd_prefix": "del ", "library": "var_list.py", "varRefreshCmd": "print(var_dic_list())" }, "r": { "delete_cmd_postfix": ") ", "delete_cmd_prefix": "rm(", "library": "var_list.r", "varRefreshCmd": "cat(var_dic_list()) " } }, "types_to_exclude": [ "module", "function", "builtin_function_or_method", "instance", "_Feature" ], "window_display": false } }, "nbformat": 4, "nbformat_minor": 2 }