{ "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 3.3.1\n", "qutip 4.5.0\n", "matplotlib.pylab 1.17.2\n", "scipy 1.3.1\n", "numpy 1.17.2\n", "krotov 1.2.0\n", "CPython 3.7.6\n", "IPython 7.17.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/+QKOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAixklEQVR4nO3de5Skd13n8fenu/p+mfuEYa4BohKRXGgDq1kNLokJBwmuwibegguOBxJFXVwD6yFsXM5hxUUPisKgc5J4MNGA0VGDIUcCURGYHhJyg5BhcpshyfRMz0xfZvpS3d/9o57qKTrVMzU9/Vy6+/M6p07X83vqqedbU8nzrd/l+f0UEZiZmc3WlHcAZmZWTE4QZmZWlxOEmZnV5QRhZmZ1OUGYmVldpbwDWEhr166Nbdu25R2GmdmisWfPnkMRsa7eviWVILZt20Z/f3/eYZiZLRqSnp5rn5uYzMysLicIMzOrywnCzMzqcoIwM7O6nCDMzKyu1BKEpM2S7pP0mKRHJb2nzmsk6WOS9kp6SNLFNfuuk/RE8rgurTjNzKy+NIe5loH/ERFfl9QD7JF0b0Q8VvOaq4DzksdrgT8DXitpNXAT0AdEcuyuiDiSYrxmZlYjtQQREc8BzyXPhyV9E9gI1CaIq4HbojLn+FckrZS0AbgMuDciBgEk3QtcCdyeVrzzdd/jB3ngmaNQlGnTpbwjAKAYUVQU5J8EFeRfpTj/HsVQlH+Pi7es4kdesTbvML5HJjfKSdoGXAR8ddaujcCzNdv7k7K5yuu993ZgO8CWLVsWJuAGfWdghP9+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+tatW7eQ8b3IyPgU3a5BmFkGlnQT0+lExFBEjCTP7wZaJK0FDgCba166KSnL3eh4eWb4mZlZmnrbSwwv1yYmSS+RpOT5JUksh4HdwHmSzpXUClwD7Morzqry1DQnJqfocoIwswwUYRRTalc7SbcDlwFrJe0HbgJaACLiE8DPAu+SVAZOANdEZfmksqQbgHuAZmBnRDyaVpyNGp2YAnANwswy0dPewtjkNJNT07Q05/NbPrWrXURce5r9fwL8yRz77gbuTiOu+RpNqnpOEGaWheo9VyNjZVZ1teYSQ96jmBaN6nAzNzGZWRaKMB+TE0SDRlyDMLMMVWd0zXPZUSeIBs00MXmqDTPLQG8B1oRwgmhQddKsrlYnCDNLX3cBVpVzgmiQm5jMLEvV/s7jEwVOEJK+T9K/VGdllfRqSb+bfmjF4iYmM8tStbWi6E1MnwLeB0wCRMRDVG5eW1ZOjmLyVBtmlr7qtWa04AmiMyK+Nqss3/u/czAyPkVLs2grOUGYWfpO1iCmcouhkQRxSNLLgQCQ9LPAc6lGVUCeh8nMstTUJDpbm3OtQTRyxbse2AH8gKQDwJPAL6QaVQGNjJd9k5yZZaqrrVTsBBER+4A3SOoCmiJiOP2wimfENQgzy1h3W76LBs15xZP0W3OUAxARH00ppkIaGXOCMLNsdbUVt4mpJ/n7/cAPc3LK7Z8CZndaL3mjE2VWdeYzYZaZLU9drSVGc+yknjNBRMT/BpB0P3BxtWlJ0geBf8okugIZGS+zeXVn3mGY2TLS3VbiuWNjuZ2/kVFM5wATNdsTSdmyMjJWptvTbJhZhrraSozmeCd1I1e824CvSbor2X4LcOvpDpK0E3gTcDAiXlVn/88DvwMIGAbeFRHfSPY9lZRNAeWI6GsgzlSNehSTmWVsMYxi+pCkzwH/OSn65Yh4oIH3voXKgkC3zbH/SeDHI+KIpKuoDKV9bc3+10fEoQbOk7rp6WB0YsrTbJhZprrbmos5iqlK0hbgEHBXbVlEPHOq4yLifknbTrH/yzWbXwE2nTbanFSreN2eZsPMMtTVVmJscpry1DSlHJYdbeQn8T+R3EUNdADnAo8DP7iAcbwD+FzNdgCflxTAJyNix1wHStoObAfYsmXLAoZ0UnUUgZuYzCxL1aH1oxNTrOgoYIKIiB+q3ZZ0MfDuhQpA0uupJIhLa4ovjYgDktYD90r6VkTcP0d8O6g0T9HX1xf1XnO2PNW3meWhs/XklN8rOloyP/8Zp6SI+Drf21cwb5JeDfw5cHVEHK45x4Hk70EqTVuXLMT55ssJwszykPeMro30QdTeUd0EXAx892xPnPRt/C3wixHx7ZrymSk9kudXADef7fnOxujMVN9OEGaWneqP0rxmdG3kitdT87xMpU/is6c7SNLtwGXAWkn7gZuAFoCI+ATwAWAN8KfJ9B3V4aznAHclZSXgryLinxv8PKlwDcLM8lD9UVrYGgTwWETcWVsg6a3AnXO8HoCIuPY0+98JvLNO+T7gggbiykx1PWonCDPL0skaRD4JopE+iPc1WLZkVYe5uonJzLJU2BpEcvPaG4GNkj5Ws6uXZbaiXDV79/hGOTPLUJE7qb8L9ANvBvbUlA8Dv5lmUEUzMlamuUm0lbIfh2xmy1dhO6mTeZG+IenTEbGsagyzjY6X6WptnlkLw8wsCx0tzTSpgDUISX8TEW8DHkjuaP4eEfHqVCMrkJHxKXras79JxcyWN0l0tea3qtypmpjek/x9UxaBFNnI+ORMW6CZWZbynNH1VE1MzyV/n84unGIaHZ/yEFczy0VXW3Nua0KcqolpmJOT9EFl3Yao/o2I3pRjK4zh8TK9HsFkZjnobisVspO6Z659y83oeJmXrmjPOwwzW4YK2cRUK5nB9VIqNYh/a3DBoCXjuFeTM7OcdLWVGBw9nsu5TzuwX9IHqCwxugZYC9wi6XfTDqxIRpJhrmZmWetqLWAfRI2fBy6IiDEASR8GHgT+T4pxFUZEZblR1yDMLA+VJqZ8+iAauTX4u0BtA3wbcCCdcIpnvDzN1HQ4QZhZLiqd1MWtQRwDHpV0L5U+iMuBr1XnZ4qIX08xvtyNeqpvM8tRV1uJifI0k1PTtGS8LnUjV727kkfVF9MJpZiqVbtO90GYWQ5qZ3Rd2dma6bkbWZP61iwCKSovFmRmeepOZnEYySFBNDKK6U2SHpA0KGlI0rCkoUbeXNJOSQclPTLHfkn6mKS9kh5KhtNW910n6YnkcV3jH2lheS0IM8vTyRpE9h3VjTRo/RFwHbAmInojoucM7qK+BbjyFPuvAs5LHtuBPwOQtJrKEqWvBS4BbpK0qsFzLiivR21meerKcVW5RhLEs8AjEfGiGV1PJyLuBwZP8ZKrgdui4ivASkkbgJ8E7o2IwYg4AtzLqRNNaqpZ201MZpaH7hxXlWvkqvc/gbslfQkYrxZGxEcX4PwbqSSgqv1J2VzlLyJpO5XaB1u2bFmAkL5X9UtxJ7WZ5aGrNb8E0UgN4kPAcSr3QvTUPAohInZERF9E9K1bt27B39+d1GaWp+4cm5gaueq9NCJeldL5DwCba7Y3JWUHgMtmlX8xpRhO6bg7qc0sR3muS91IDeJuSVekdP5dwC8lo5leBxxL1qG4B7hC0qqkc/qKpCxzI+NTtDY30er1qM0sBzOjmCayH8XUyM/idwHvlTQOTHIG60FIup1KTWCtpP1URia1UHmDTwB3A28E9lJpxvrlZN+gpN8DdidvdXNEnKqzOzWj42U6vZqcmeWkrdRES7OK2cR0NutCRMS1p9kfwPVz7NsJ7JzvuRfK6Hh5ppPIzCxrknJbE6LR9SBWUblXYWbSvmQI65I3OlF2B7WZ5aqrtcTIWAEThKR3Au+h0lH8IPA64D+An0g1soIYHZ+a6SQyM8tDXjO6NtLz+h7gh4GnI+L1wEXA0TSDKpIRryZnZjnrastn0aBGEsRYzWJBbRHxLeD70w2rONwHYWZ562orMZLDXEyNXPn2S1oJ/B1wr6QjwNNpBlUkx72anJnlrLutxHPHxjI/byOjmH46efpBSfcBK4B/TjWqAhkZL89Mt2tmlodCj2KqiogvpRVIEUVEpYnJNQgzy1GRO6mXrfHyNGWvR21mOetOahDzmFT7rDhBnMLx5Nb2Ls/kamY56morMR0wNjmd6XkbShCStkp6Q/K8Q1JhZnNNkxcLMrMiqF12NEuNLDn6K8BngE8mRZuojGha8jzVt5kVQVdOiwY1UoO4HvhRYAggIp4A1qcZVFG4BmFmRZDXsqONJIjxiJiobkgqAdn2lORkZCZBuA/CzPKT16JBjSSIL0l6P9Ah6XLgTuAf0g2rGGY6qV2DMLMcFbmJ6UZgAHgY+FUqazj8bppBFcVMDcJTbZhZjvLqpG7kTupp4FPApyStBjZF1oNxczLqTmozK4CTNYhs52NqZBTTFyX1JslhD5VE8YeNvLmkKyU9LmmvpBvr7P9DSQ8mj29LOlqzb6pm364z+EwLxp3UZlYE3Tk1MTVy5VsREUPJuhC3RcRNkh463UGSmoGPA5cD+4HdknZFxGPV10TEb9a8/teoTCVedSIiLmzwc6RidGKKlmZ5PWozy1W1mbuIndQlSRuAtwH/eAbvfQmwNyL2JaOg7gCuPsXrrwVuP4P3T53nYTKzImhqEp2tzYXspL4ZuIfKxX63pJcBTzRw3Ebg2Zrt/UnZi0jaCpwLfKGmuF1Sv6SvSHrLXCeRtD15Xf/AwEADYTVuxGtBmFlBdLWVMl80qJFO6jupDG2tbu8DfmaB47gG+ExE1PbAbI2IA0lC+oKkhyPiO3Xi2wHsAOjr61vQzvPRca9HbWbF0J3DokFzXv0k/TGnuCEuIn79NO99ANhcs70pKavnGip3bNe+/4Hk7z5JX6TSP/GiBJGm4xNTdPomOTMrgK627JuYTvXzuP8s33s3cJ6kc6kkhmuAn5v9Ikk/AKwC/qOmbBVwPCLGJa2lMtXH759lPGdsxDUIMyuIrtYSI2MFSRARcevZvHFElCXdQKX/ohnYGRGPSroZ6I+I6tDVa4A7Zt1b8Urgk5KmqfSTfLh29FNWRsfLnNPTnvVpzcxeJI9lR0/78zhZZvRFTU0R8ROnOzYi7qZy53Vt2QdmbX+wznFfBn7odO+fttFxr0dtZsVQyE5q4L01z9updFBnv/ZdDkYnvB61mRVDHutSNzKKac+son+X9LWU4imU0fEyna5BmFkB9LRnvy51I01Mq2s2m4DXACtSi6ggxstTTE6FO6nNrBC6WkuMTU5Tnpqm1JzN7A6NXP32UOmDEJWmpSeBd6QZVBFUJ8XyetRmVgTVdWlGJ6ZY0VGQBBER52YRSNF4oj4zK5LaCftWdLRkcs5GmpjagXcDl1KpSfwr8ImIyHa8VcaGk/HGPe3ZfBFmZqeSx6JBjfw8vg0YBv442f454C+Bt6YVVBFUO4N62l2DMLP85bHsaCNXv1dFxPk12/dJyvymtawNj00CXizIzIohj0WDGunp+Lqk11U3JL2Ws5+Go/BcgzCzIumaWXZ0MrNzNnL1ew3wZUnPJNtbgMclPQxERLw6tehyNJT0QXQ7QZhZAZxsYsquBtHI1e/K1KMooOqkWD1t7qQ2s/zlsexoI8Ncn84ikKIZGZ+k1CTaW7zcqJnlryuHTmpf/eYwPFamu72EpLxDMTOjrdREqUmZ1iCcIOYwMua1IMysOCRlPmGfE8QchsfLvknOzAol62VHU00Qkq6U9LikvZJurLP/7ZIGJD2YPN5Zs+86SU8kj+vSjLOe4bFJelyDMLMCyXrZ0dSugJKagY8DlwP7gd2SdtVZGe6vI+KGWceuBm4C+qhM77EnOfZIWvHONjJeZr1XkzOzAsl60aA0axCXAHsjYl9ETAB3AFc3eOxPAvdGxGCSFO4l4+G2I2Nl3yRnZoXS3VaamScuC2kmiI3AszXb+5Oy2X5G0kOSPiNp8xkei6Ttkvol9Q8MDCxE3EAyislNTGZWIL3tLctqmOs/ANuSu7HvBW490zeIiB0R0RcRfevWrVuwwNxJbWZF09NeYuhEdlNtpJkgDgCba7Y3JWUzIuJwRIwnm39OZVqPho5N03h5ionytJuYzKxQejtaGBpbGgliN3CepHMltQLXALtqXyBpQ83mm4FvJs/vAa6QtErSKuCKpCwT1Wk23MRkZkXS215ZdnSiPJ3J+VK7AkZEWdINVC7szcDOiHhU0s1Af0TsAn5d0pupLGU6CLw9OXZQ0u9RSTIAN0fEYFqxzuaZXM2siKrN3sNjk6zpbkv9fKleASPibuDuWWUfqHn+PuB9cxy7E9iZZnxzGXYNwswKqLejck0aGitnkiDy7qQupGFP9W1mBdRbU4PIghNEHdUmpl6PYjKzAqk2MQ2dyGaoqxNEHV5u1MyK6GQTk2sQuanWINzEZGZF0jtTg3CCyE31H9+jmMysSKrXpKym23CCqOPYiUnaW5poKzXnHYqZ2Yyu1hJNchNTro6dmGRFhzuozaxYmppET3uLm5jy5ARhZkXV25HdjK5OEHU4QZhZUfW0ZTcfkxNEHcdOlJ0gzKyQejtKvg8iT0MnJul1gjCzAuptdw0iV25iMrOi6mlvcR9EXspT04yMu4nJzIqp0sTkGkQuhpLM7ARhZkXU297C8HiZqelI/VxOELMcSzKzE4SZFdHKzsq16VgGtQgniFmcIMysyFZ3tQJw5PhE6udKNUFIulLS45L2Srqxzv7fkvSYpIck/YukrTX7piQ9mDx2zT42LU4QZlZkKzuTBDGafoJIbTY6Sc3Ax4HLgf3Abkm7IuKxmpc9APRFxHFJ7wJ+H/hvyb4TEXFhWvHNxQnCzIpsdTVBHF/cTUyXAHsjYl9ETAB3AFfXviAi7ouI48nmV4BNKcbTECcIMyuyah9EFjWINBPERuDZmu39Sdlc3gF8rma7XVK/pK9IestcB0nanryuf2Bg4KwChpP/6NVqnJlZkWTZB1GIBQ8k/QLQB/x4TfHWiDgg6WXAFyQ9HBHfmX1sROwAdgD09fWd9bivwdEJetpLtJbcf29mxdPZ2kxrc9Oib2I6AGyu2d6UlH0PSW8A/hfw5ogYr5ZHxIHk7z7gi8BFKcY64/DoBGu6XHsws2KSxMrOlkXfxLQbOE/SuZJagWuA7xmNJOki4JNUksPBmvJVktqS52uBHwVqO7dTMzg6PlOFMzMrotVdrYu7iSkiypJuAO4BmoGdEfGopJuB/ojYBXwE6AbulATwTES8GXgl8ElJ01SS2IdnjX5KzeGRCTat6sziVGZm87Kys2VxJwiAiLgbuHtW2Qdqnr9hjuO+DPxQmrHNZXB0ggs2rczj1GZmDVnd1cq3XxhJ/Tzuia0RERw5PsHqbjcxmVlxrexs5ehiv5N6sRkaKzM5Fe6kNrNCW93ZypHjk6lP2OcEUWMwGRXgTmozK7L1vW1MTcfMNSstThA1Bkcro2ydIMysyNb3tAFwcHgs1fM4QdQYGK4kiLXdbTlHYmY2t/W97QAcHBo/zSvPjhNEjeePVbLxS1a05xyJmdncqjWIF4Zcg8jM80PjtDY3zcyWaGZWROtmmphcg8jM88dOsL63jaYm5R2Kmdmc2krNrOpscQ0iS88PjfGSXjcvmVnxre9pdw0iS88fG3P/g5ktCut72zjoGkQ2IsI1CDNbNDau7ODA0ROpnsMJInFoZIKxyWk2rurIOxQzs9PauqaLQyMTDI+lty6EE0TiyUOjAJy7tivnSMzMTm/bmsqs008fPn6aV86fE0TiqSRBvGxtd86RmJmd3rbkx+xTh0dTO4cTRGLfoVFamsVLV7oPwsyKb2tSg6j+uE2DE0TiW88P8bK13ZSa/U9iZsXX2Vpi48oOvvn8cGrnSPVqKOlKSY9L2ivpxjr72yT9dbL/q5K21ex7X1L+uKSfTDPOiOCBZ45y0ZaVaZ7GzGxBXbhlJQ8+czS1908tQUhqBj4OXAWcD1wr6fxZL3sHcCQiXgH8IfB/k2PPp7KG9Q8CVwJ/mrxfKh797hDHTkxy8dZVaZ3CzGzB9W1dxYGjJ9h7MJ3V5dKsQVwC7I2IfRExAdwBXD3rNVcDtybPPwP8F1UWp74auCMixiPiSWBv8n4Lbmhskvfe+Q1aS01c/spz0jiFmVkq3vTql1JqEu/+9J5UFg9KM0FsBJ6t2d6flNV9TUSUgWPAmgaPBUDSdkn9kvoHBgbOOMiethI/8JIePnbNRazyOhBmtois62njD956Aa/ZuormFOaQKy34O2YsInYAOwD6+vrOOIVK4o+uuWjB4zIzy8JbLtrIWy6q+/v5rKVZgzgAbK7Z3pSU1X2NpBKwAjjc4LFmZpaiNBPEbuA8SedKaqXS6bxr1mt2Adclz38W+EJERFJ+TTLK6VzgPOBrKcZqZmazpNbEFBFlSTcA9wDNwM6IeFTSzUB/ROwC/gL4S0l7gUEqSYTkdX8DPAaUgesjYiqtWM3M7MVU+cG+NPT19UV/f3/eYZiZLRqS9kREX719vm3YzMzqcoIwM7O6nCDMzKwuJwgzM6trSXVSSxoAnp7n4WuBQwsYzmLgz7z0LbfPC/7MZ2prRKyrt2NJJYizIal/rp78pcqfeelbbp8X/JkXkpuYzMysLicIMzOrywnipB15B5ADf+alb7l9XvBnXjDugzAzs7pcgzAzs7qcIMzMrK5lnyAkXSnpcUl7Jd2YdzxZkPSUpIclPShpSc5uKGmnpIOSHqkpWy3pXklPJH+X1CLkc3zmD0o6kHzXD0p6Y54xLjRJmyXdJ+kxSY9Kek9SvmS/61N85gX/rpd1H4SkZuDbwOVUljXdDVwbEY/lGljKJD0F9EXEkr2ZSNKPASPAbRHxqqTs94HBiPhw8mNgVUT8Tp5xLqQ5PvMHgZGI+IM8Y0uLpA3Ahoj4uqQeYA/wFuDtLNHv+hSf+W0s8He93GsQlwB7I2JfREwAdwBX5xyTLYCIuJ/KGiO1rgZuTZ7fSuV/qiVjjs+8pEXEcxHx9eT5MPBNKuvXL9nv+hSfecEt9wSxEXi2Zns/Kf1DF0wAn5e0R9L2vIPJ0DkR8Vzy/HngnDyDydANkh5KmqCWTFPLbJK2ARcBX2WZfNezPjMs8He93BPEcnVpRFwMXAVcnzRNLCvJ0rbLoX31z4CXAxcCzwH/L9doUiKpG/gs8BsRMVS7b6l+13U+84J/18s9QRwANtdsb0rKlrSIOJD8PQjcRaWpbTl4IWm/rbbjHsw5ntRFxAsRMRUR08CnWILftaQWKhfKT0fE3ybFS/q7rveZ0/iul3uC2A2cJ+lcSa1U1sTelXNMqZLUlXRsIakLuAJ45NRHLRm7gOuS59cBf59jLJmoXiQTP80S+64licra9t+MiI/W7Fqy3/VcnzmN73pZj2ICSIaC/RHQDOyMiA/lG1G6JL2MSq0BoAT81VL8zJJuBy6jMg3yC8BNwN8BfwNsoTIt/NsiYsl06s7xmS+j0uQQwFPAr9a0zS96ki4F/hV4GJhOit9PpU1+SX7Xp/jM17LA3/WyTxBmZlbfcm9iMjOzOThBmJlZXU4QZmZWlxOEmZnV5QRhZmZ1OUGYzZOklZLenTx/qaTP5B2T2ULyMFezeUrmwfnH6sypZktNKe8AzBaxDwMvl/Qg8ATwyoh4laS3U5k9tAs4D/gDoBX4RWAceGNEDEp6OfBxYB1wHPiViPhW1h/CbC5uYjKbvxuB70TEhcBvz9r3KuC/Aj8MfAg4HhEXAf8B/FLymh3Ar0XEa4D3An+aRdBmjXINwiwd9yVz9Q9LOgb8Q1L+MPDqZCbOHwHurEytA0Bb9mGazc0Jwiwd4zXPp2u2p6n8f9cEHE1qH2aF5CYms/kbBnrmc2Ayf/+Tkt4KlRk6JV2wkMGZnS0nCLN5iojDwL9LegT4yDze4ueBd0j6BvAoXu7WCsbDXM3MrC7XIMzMrC4nCDMzq8sJwszM6nKCMDOzupwgzMysLicIMzOrywnCzMzq+v+zZZTMkShxNAAAAABJRU5ErkJggg==\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/+QKOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAA9aElEQVR4nO3dd3gU5fbA8e8hJIQECBA6ARI6SAkkIEURRaQpXQQEBQs2LCh61WtB77VcRcWKWOjSRZqCWFCRIiQQeu8h9EBIr+/vj1n8RQgQIJvZcj7Psw+7OzO7Z9gkZ+ct5xVjDEoppdT5itgdgFJKKdekCUIppVSeNEEopZTKkyYIpZRSedIEoZRSKk9F7Q6goJQrV86EhobaHYZSSrmV6Ojok8aY8nlt85gEERoaSlRUlN1hKKWUWxGRAxfbpk1MSiml8qQJQimlVJ40QSillMqTJgillFJ50gShlFIqT5oglFJK5UkThFJKqTw5dR6EiHQGPgR8gK+MMW+ft70GMB4oD8QDg4wxsY5t2cAmx64HjTHdnRJkZir8MRp8/cE3wLqVqAilqkCpqhBYDkSc8tZKuZus7BxOp2RyOiWDpPQssnMMWdmGHGPw9y1CgF9RShQrSlCALyWLFUX0d8etOS1BiIgP8CnQEYgF1orIAmPM1ly7jQYmG2MmicgtwFvAYMe2VGNMuLPi+1t6Ivz5PpicvLf7B0GFhtatUmOo0QbK1dWkoTzasbNpbDh0hm1HEtl/Kpl9J5M5cCqZ0ymZ+X6NQD8fKpcuTuUgf2qVL0G9SiWpV6kkdSuWpEQxj5mj69HEWQsGiUhrYJQxppPj8QsAxpi3cu2zBehsjDkk1leNBGNMKce2JGNMify+X2RkpLnqmdTGQHYmZKZARjIkHYWzcZAQCyd3wrGtcHwbpCdY+wcEW4mibmeo2wUCg6/ufZVyEbGnU1i+6yTLd50g+sBpjp1NB6zvQVWCilMjOIAawYFULFWM4EA/ygT6UaJYUYoWKYJPEaGIQFpWDsnpWSSlZ3EmJYMjCWkcOZNGXEIqe44nkZyRDUARgfqVStEitAyRoWVpVTOY8iWL2Xn6Xk1Eoo0xkXltc2YarwocyvU4Frj+vH02AL2xmqF6ASVFJNgYcwrwF5EoIAt42xgzz2mRikBRP+tWvDQEVYWqEf/cxxiI3wsHVlq3fb/DtoUgRaB6G2h6F1zXG4rlO6cpZRtjDFuPnGXhhiMs3XKUvSeTAagc5E+rmsE0DSlN02pBNKwcRHE/n2t+v5wcw+EzqWw/msjmwwlEHYhnVlQsk1ZZVR6ahgRxS/2KdGhQgeuqlNKmKRfhzCuIvlhXBw84Hg8GrjfGDM+1TxXgEyAM+APoAzQyxpwRkarGmMMiUhP4FehgjNlz3nsMA4YBVK9ePeLAgYuWFCl4xsCRGNj+PWyZB6d2gV8JaNQHWj5oNUcp5WKOJqQxK+oQ82IOs/dEMj5FhDa1gmlfrwLt6pSjdoUShfbHOTM7h61xZ1m+6wS/bD9OzKEzGAPVyhanR9Oq9GxWhdoVShZKLN7sUlcQtjYxnbd/CWC7MSYkj20TgUXGmDkXe79ramK6VsbAoTWwbjJsmWs1VdXqADc8BaE3an+FspUxhhW7TzF19QF+2naM7BxDq5pluaNpFbo0qkzZQD+7QwTgZFI6y7YfZ8GGOFbsPkmOgeuqlKJ/i2r0ah6i/RZOYleCKArsBDoAh4G1wEBjzJZc+5QD4o0xOSLyBpBtjHlFRMoAKcaYdMc+q4Ae53Vw/4OtCSK31NOw9mv463NIPgEhLaHja1afhVKFKDvH8MOmI3y6bDfbjyZSNtCPOyNDuLtlDaoHB9gd3iUdT0xj4YYjzF0Xy5a4swT6+dC7eQiDW9egbkW9qihItiQIxxt3BcZgDXMdb4x5Q0ReB6KMMQsczVBvAQariekxR1JoA4wDcrDmaowxxnx9qfdymQRxTmYqxHxjDaFNPGJ1aHd4FSo2tDsy5eGysnP4bv1hxv62h70nk6lVPpCHb6rFHU2r4O977f0JhckYQ8yhM0xZfYBFG4+QkZXDjXXK8Wj72rSqWVb7KgqAbQmiMLlcgjgnI8W6mvhzDGQkQosH4JaXrOGzShUgYwzLdhzn7cXb2XksiYaVSzH8ltp0vq4SRYq4/x/S+OQMpq85yIQV+zmZlE6z6qV5tH1tOtSv4BHnZxdNEK4gJR6WvQlrv4ISFaDzW9aoJ/0GpArA5sMJvPH9NlbtPUVYuUCe61SPzo0qeeQ37LTMbGZHxzLu9z3Enk6lYeVSPNu5Hu3rlvfI83U2TRCu5PA6WDTCGgFV+1bo/rE1a1upq5CUnsV7S3cwaeV+Sgf48WSHOgy8vjq+Pp5fRScrO4cFG+IY8/MuDsan0DK0LM91rkdkaFm7Q3MrmiBcTU62dSXx8yjw8YVu70PjvnZHpdyIMYYftxxl1IKtHEtMY9D1NRjZqR5BxX3tDq3QZWTlMHPtQT76dTcnEtO5tUEFXuzagJrldU5SfmiCcFWn9sB3D0HsWqu5qdt7EKDfftSlxSdn8OLcTSzZcpQGlUvxZq9GNKtexu6wbJeSkcWEFfsZ+9se0rOyue+GMB6/pY4Oj70MTRCuLDsLVn4Iy96CkpWh38QLZ3Er5fDr9mM8N2cTZ1Mzefq2ujxwQxhFvaA56UocT0zjnSU7mBMdS/mSxfhX5/r0blZVO7IvQhOEO4iNhtlDrCGxnd60ZmNrh5tySMnI4j+LtjF9zUHqVyrJB3eF06ByKbvDcmnrD55m1MKtbDh0hsgaZXird2Pq6ByKC2iCcBcp8TDvEdi5BK7rBd0/0dpOit3HE3l46jr2nEhiWLuaPN2xLsWKutd8Brvk5BjmrIvlzR+2kZyexWM31+aR9rX0/y8XTRDuJCfHanL65XWrxPiA6VC6ut1RKZvMjznMC3M3EeDnw4f9m9G2djm7Q3JLJ5PSeX3hVhZsiKN2hRL8r09jImpofx9cOkFo46WrKVIEbhgBd8+BM4fgi5vh4Gq7o1KFLD0rm5fnbebJGTFcV6UU3z9xoyaHa1CuRDE+GtCMCUNakJqRTd/PV/Hawi2kZWbbHZpL0wThqmp3gAd/sWZcT7oD1n9jd0SqkJxKSmfQV38xZfUBHmpXk2kPtqJiKX+7w/IIN9evwNIR7binVQ0mrNhPt4+Ws+HQGbvDclmaIFxZuTrwwM9QvTXMf9Saie0hTYIqbzuOJtLj0xVsjE3gowHNeKFrA6+Y9FaYAosV5bUejZh6//WkZGTTe+xKPvhpJ5nZF1lV0ovpT56rCygLg76FZoPg9//BwietobHK4yzbfpw+Y1eSnpXDzIda072pzrB3phvqlGPJU+3o0bQKH/6yi96frWT38US7w3IpmiDcgY+vNaKp3bOwbhLMGmwVAVQeY+rqA9w/aS01ggNYMLwt4dVK2x2SVwgq7sv7d4Uz9u7mxJ5O4faP/2TW2kN4yuCda6UJwl2IWFVgu46GHYthSk9rWKxya8YYPvhpJy/N28zN9Sow++HWVA4qbndYXqdL48r8+FQ7mlcvw3PfbuTJGTEkpmXaHZbtNEG4m5YPwp0TIW49TLwdkk7YHZG6Stk5hpfmbebDX3bRNyKEcYMjCPDTshB2qVDKnyn3X8/I2+ry/aYjdPvoTzbGnrE7LFtpgnBH1/WEgbMgfi9M7Apnj9gdkbpCaZnZDJ+2jm/+OsjDN9Xi3b5NtGSGC/ApIgy/pQ4zh7UiKzuHPmNX8tXyveTkeGeTk/5EuqtaN1ud12fjYEIXa86EcgupGdk8ODmKxZuP8lK3Bjzfpb6uY+BiIkPL8sOTN3JzvQr89/ttDJsSTUKq9zU5aYJwZ6FtYfA8qy9iQlfrikK5tJSMLO6buJY/d5/k3b5NeODGmnaHpC6idIAf4wZH8OodDfltx3F6fPIn246ctTusQqUJwt1VawH3LrCWM53QTZOEC0tKz2LIhLX8te8UH/QL587IanaHpC5DRBjaNowZw1qRkpFNr89W8N36WLvDKjSaIDxBlXC4dxFkpcKk7nDmoN0RqfMkpmUyZPwaog+c5sP+zejZrKrdIakrEBlalkVP3ECTkNKMmLmBV+dvJiPL8yfWaYLwFJUaWc1N6Wet0hwJh+2OSDkkpWdxz/g1xBw6wycDmnGHToBzSxVK+vPNA9fz4I1hTFp1gP5frOJIQqrdYTmVJghPUiUcBn0HyadgcndIPGZ3RF4vLTObByatZWNsAp8MbE6XxpXtDkldA1+fIvy7W0M+Hdic7UcTuePjP4na77nzkTRBeJqQCBg0xxr6Ork7JJ+0OyKvlZGVwyNTo/lrXzzv92tK50aV7A5JFZBuTSoz/7G2lChWlAFfrmbmWs9s1tUE4Ymqt4KBM+D0fmvGdVqC3RF5newcw4iZMSzbcYI3ejamR7j2OXiaOhVLMv+xG2hVM5h/fbuJUQu2kOVhBf80QXiqsHbQ/xs4vh2mD4BMz24rdSU5OYbnv93I95uO8O+uDRh4vS745KmCAnyZMKQFD9wQxsSV+7l3whpOJ2fYHVaB0QThyWrfCr0+hwMrYc79WgW2kLz5wzZmR8fyZIc6PNhO5zl4uqI+RXjp9oaMvrMpa/edpsenK9h5zDOqwmqC8HSN+0KXd2DH97DoSV1Pwsm+Wr6Xr/7cx5A2oTx1ax27w1GFqG9ECDMeakVqZja9Pl3BT1vdf5CIJghvcP0wuOlfsH4q/Pyq3dF4rIUb4vjv99vo2rgSr9zeUMtneKHm1cuwcPgN1KpQgmFTovh02W63Lh2uCcJbtH8BWjwAKz6EFR/ZHY3HWb33FM/M2kCL0DK83y+cIkU0OXirSkH+zHIs+PTujzsYOXuj206q09rC3kLEampKOQU/vQwlK0GTfnZH5RF2Hktk2OQoqgcH8OU9kfj7+tgdkrKZv68PY+4Kp2a5Enzw804OnU5h3KAIygT62R3aFdErCG9SxAd6jYPQG2Heo7Bvud0Rub2jCWkMGb8Gf18fJg5tQekA9/oDoJxHRHjy1jp8NKAZMYfO0OuzFew5kWR3WFdEE4S3KVoM7poCZWvCzLutYbDqqqRkZPHA5LUkpGYyYWgLQsoE2B2SckHdm1Zh+oOtSEzLotenK1i5230mr2qC8EbFy8Dds8GnGHxzp5bkuAo5OYZnZm1gS9xZPh7YjOuqBNkdknJhETXKMO+xtlQK8uee8WuYscY9Zl5rgvBWZWrAwJmQchKm9YN097r0tduYn3eyePNRXuzSgFvqV7Q7HOUGqpUNYM4jbWhTuxzPz93Emz9sI9vFV6rTBOHNqjaHvhPg6Eb4VifS5df8mMN89Otu+kWG8MCNYXaHo9xIKX9fxt8byT2ta/DFH3t5eGo0KRmu+3unCcLb1esMXd+FnUtg8XM6ke4yYg6d4dk5G2kZWpb/9myscx3UFSvqU4TXezRi1B0N+WXbMe78fBVHE9LsDitPTk0QItJZRHaIyG4ReT6P7TVE5BcR2Sgiv4lISK5t94rILsftXmfG6fVaPABtnoCor2HVp3ZH47KOJKTy4OQoKpQsxthBzfErqt+v1NUb0jaMr+9twf6TyfT49E82xbpeUU2n/YSLiA/wKdAFaAgMEJGG5+02GphsjGkCvA685Ti2LPAqcD3QEnhVRMo4K1YF3PoaNLjDmiOx80e7o3E5aZnZPDwlmtSMbL6+twXBJYrZHZLyADfXr8C3j7ahaJEi9Bu3iiWbj1z5iyQdh7iYAo8NnHsF0RLYbYzZa4zJAGYAPc7bpyHwq+P+slzbOwE/GWPijTGngZ+Azk6MVRUpYs2RqNjIKux3bKvdEbmUUQu2sCE2gff6NaVepZJ2h6M8SP1KpZj3WFvqVy7Jw1PX8dlv+SzPkZMDqz6DjyNg7jCnNA87M0FUBQ7lehzreC63DUBvx/1eQEkRCc7nsYjIMBGJEpGoEydOFFjgXssvEAbMAL8AmH6XLjbkMH3NQWasPcTwm2vT6Tpd9EcVvPIlizH9wVbc0bQK7yyxynOkZ2Vf/IDkk/BNX/jxBah2PfSfZlVLKGB2N6KOBG4SkfXATcBh4BL/K/9kjPnCGBNpjIksX768s2L0LkFVof9067J15iDISrc7IlvFHDrDq/O30K5ueUZ0rGt3OMqD+fv68FH/cJ66tQ7frotl8FdriM9rbYlTe+DLW2D/n3D7B9acpnK1nRKTMxPEYaBarschjuf+ZoyJM8b0NsY0A/7teO5Mfo5VThQSAT0/g4OrYOFTXjuy6WRSOo9MjaZiUDE+6h+OjxbgU04mIjx1a12rPEfsGXp+uoLdx3OtLXFkA4zvBOmJMHQxRN7nlCuHc5yZINYCdUQkTET8gP7Agtw7iEg5ETkXwwvAeMf9H4HbRKSMo3P6NsdzqrA06gM3PQ8bplkVYL1MVnYOw6etIz45g7F3R2iNJVWoujetwoxhrUjJyKLXZytZvusEHI6GibdbFRDu+9H6IudkTksQxpgsYDjWH/ZtwCxjzBYReV1Eujt2aw/sEJGdQEXgDcex8cB/sJLMWuB1x3OqMN30L7iuF/w8CrZ/b3c0hep/S7azem88b/VuTKOqWkZDFb7m1a3yHFVLF+etiXNIn9jTKpNz/49QvnCaO8WdF7PILTIy0kRFRdkdhufJSIGJXeHETusHs1JjuyNyuiWbj/Lw1GjuaV2D13s0sjsc5eWS47aS9VUXkrOLMKvRFwzv3YGiPgX33V5Eoo0xkXlts7uTWrk6vwCr09q/FEwfCMmn7I7IqQ7Fp/DsnA00DQnipW7nT9tRqpDF7yNwem9KFfdlXuPPGBOdwQOTozibllkob68JQl1eqcpw1zeQdAxm3wvZhfPDWdgysqx+B4BPBupMaWWzhFiY3B2y0pB75vNo3y682asxf+46SZcxy1m5x/nD0PU3QOVPSATc8SHsXw4//tvuaJzircXb2BCbwLt9m1KtrK7toGyUeAwm94DUMzD4O6h4HQADr6/OzIda41e0CAO//IsRM2PYfzLZaWFoglD5Fz4AWj0Ga8bBusl2R1OgftxylAkr9jO0bSidG+lkOGWjlHiY0hPOxllzHKo0+8fmiBpl+OGJG3n4plos3nyEm9/7jUemRjslFF2TWl2Zjq/D8S2w6GkoXx+qtbQ7omt2KD6FZ2db/Q4vdGlgdzjKm6WdhSm9rMlwd8+C6q3y3K24nw/Pd6nPfW1DmfrXQafNVdIrCHVlfIpaa0gEVbVmWp+Nszuia3Ku38Gg/Q7KZplpMGMgHNsM/SZDzfaXPaRCKX+e7liXp2+r55SQ9LdBXbmAstbIpvQkmHG39YPtpt5evN3R79BE+x2UfbKzrEW79i+HnmOtdVpcgCYIdXUqNoTe4yBuHSwa4ZblOH7dfozxK/YxpE0onRtVtjsc5a2MgUVPwfZF0Pl/0KSf3RH9TROEunoN7vj/chyrx9odzRU5npjGs7M30qByKV7oWt/ucJQ3+3kUrJ8C7Z6DVg/bHc0/aCe1ujY3/ctqM136ElRoALVutjuiy8rJMYycvZGk9Cxm9A+nWFEfu0NS3mrFh7BiDETeDze/aHc0F9ArCHVtihSBXp9DubowewjE77U7osuasHI/f+w8wUu3N6RORV38R9lk/VT46RW4rre1LrwLrm+uCUJdu2IlYcA06/6Mu63Oaxe1Ne4s/1u8nVsbVGTQ9dXtDkd5q22LYMHjUOsWayXHIq55FasJQhWMsjXhzglwYjt895C1HKKLScvM5skZ6wkK8OV/fRojLviNTXmBfcthzn1QpTn0mwJFXbeUvCYIVXBq3QId/2ONxvjjXbujucAb329j1/Ek3u/XlOASxewOR3mjuBiYPgDKhlmzpIuVsDuiS9IEoQpW68egyV3w25sutYbEz1uPMWX1AR64IYwb6+jytMoGJ3fD1D5QvDQMmmvNJ3JxmiBUwRKxivpVaQZzh8HxbXZHxPGzaTz37UYaVi7Fs52dM+NUqUs6G2eV0AAYPM+qROAGNEGogudb3CoP7htglQ5IPW1bKMYYRs7ZSEpGFh8N0CGtygYp8VZySD0Ng+ZAudp2R5RvmiCUcwRVhbumwJlDVodcdpYtYUz96yB/7DzBv7s2oHYFHdKqCllGMkzrB/H7YMD0CyqzujpNEMp5qreCbqNhz6/wy6hCf/t9J5N58/tttKtbnkGtahT6+ysvl5UBMwfD4WjoOx7CbrQ7oiumM6mVc0UMgaObYOXHUKlJodWZycrO4ZlZMfj6CO/0aaJDWlXhysmBeQ/Dnl+g+yfQ4Ha7I7oqegWhnK/z21CjrTUxKG59obzluD/2su7gGf7TsxGVgvwL5T2VAqzie0v+BZu/hVtHQfPBdkd01TRBKOfz8YU7J0FgeWumddJxp77dlrgExvy8k25NKtO9aRWnvpdSF/hjNKz5AloPh7ZP2R3NNcl3ghARHxGpIiLVz92cGZjyMCXKQ/9vrBEdMwdb7bNOkJ6VzdMzN1A6wI//9mikTUuqcEWNh2X/haYDrEmjbv7zl68EISKPA8eAn4DvHbdFToxLeaLKTaHHJ3BoNSx+zilv8f5PO9lxLJF3+jShTKDrljBQHmjLPGsp3jqdoPvHViFLN5ffTuongXrGmFPODEZ5gcZ9rU7rFWOgUmNocX+BvfTa/fF88cdeBrSszs31KxTY6yp1WXt/h7kPWmu03znRalb1APlNcYeABGcGorxIh1egdkfrKuLAygJ5yaT0LJ6eFUO1MgG81K1BgbymUvkSt96aEBpcGwbOBD/PWbo2vwliL/CbiLwgIk+fuzkzMOXBivhAn6+gTKjVH3Hm0DW/5BvfbyP2dCqj72xKYDEdva0Kyak9MLUvFC8Lg76F4mXsjqhA5TdBHMTqf/ADSua6KXV1ipeG/tMhOwNm3g0ZKVf9Usu2H2f6moMMu7EmLcNcvwCa8hBnj8CUntb9wd9BKc8bMZevr1rGmNcARKSE47Hrrgij3Ef5utD7S5je35oj0eerKx71kZCSyfNzN1KvYklGdKzrpECVOk/qaZja2xqVN2SRW9VXuhL5HcXUSETWA1uALSISLSLXOTc05RXqdYZbXoLNc2DlR1d8+H++38rJpAxG39kUf18txKcKQUYKTOsPp3ZD/2luV1/pSuS3iekL4GljTA1jTA3gGeBL54WlvMqNz0DDnvDzKNj1c74PW7b9OHOiY3nkplo0DglyWnhK/S07E+YMhUN/WVe/NW+yOyKnym+CCDTGLDv3wBjzGxDolIiU9xGBnp9BhYZW5ddTey57SEJqJi/M3UTdiiV4vINnXt4rF2MMLHgCdi6Bbu/BdT3tjsjp8j2KSUReFpFQx+0lrJFNShUMv0Drcr2Ij7UkY9rZS+7+xvdbOZGUzug7m+oaD6pw/PQKbJgGN/+7QOfvuLL8Joj7gPLAXMetvOM5pQpOmRrQbxLE77nkGhK/7TjOrKhYHmpXkyYhpQs3RuWdVnxo9ZG1HAbtnrU7mkKTrwRhjDltjHnCGNPccXvSGGPfMmHKc4W1g66jYfdPsPTfF2w+m2Y1LdWpUIInb61jQ4DK60RPtK4erusNnf/n9vWVrsQlh7mKyBhjzFMishAw5283xnR3WmTKe0UOhZO7YPWn1uzUlg/+venN77dx7GwaYx9tq01Lyvk2fwsLn7Jm/vca5xH1la7E5eZBTHH8O/pqXlxEOgMfAj7AV8aYt8/bXh2YBJR27PO8MeYHEQkFtgE7HLuuNsY8fDUxKDd123+spqbF/4KyYVD7Vv7YeYIZaw/x8E21CK9W2u4IlafbuRTmDoPqraHfZCjqfcUfL5kOjTHRjrvhxpjfc9+A8EsdKyI+wKdAF6AhMEBEGp6320vALGNMM6A/8FmubXuMMeGOmyYHb3OuHEeFBjB7KMmxm3n+243UKh/IU9q0pJxt/wqYNRgqXgcDZ3hUfaUrkd/rpXvzeG7IZY5pCew2xuw1xmQAM4Ae5+1jgFKO+0FAXD7jUd6gWEkYMAOK+pM+uS8ZZ4/zrk6IU84Wtx6m3QWlq8OgueDvvXNsLpkgRGSAo/8hTEQW5LotA+Iv89pVsarAnhPreC63UcAgEYkFfgAez7UtTETWi8jvIpLnat8iMkxEokQk6sSJE5cJR7ml0tWIueFzAtJP8l3wZzSv4p3f5FQhOb4dpvS2iu4NngeB5eyOyFaX64NYCRwBygHv5Xo+EdhYAO8/AJhojHlPRFoDU0SkkeM9qxtjTolIBDBPRK4zxvxjcLwx5gusWd5ERkZe0Imu3F9SehaP/S7cVvwpXk1616rZ1GucV40kUYXk9H6r+J6PL9wzD4LO/z7rfS6ZIIwxB4ADQOureO3DQLVcj0Mcz+V2P9DZ8V6rRMQfKGeMOQ6kO56PFpE9QF0g6iriUG7srR+2EZeQyu0PPwb7i1nLOQbXgZu8Zyy6KgSJR2FyD8hMhaE/QHAtuyNyCfkt1tdKRNaKSJKIZIhItohceqorrAXqiEiYiPhhdUIvOG+fg0AHx3s0APyBEyJS3tHJjYjUBOqgM7e9zordJ/nmr4M8cEMYETXKQLuR0KS/lSQ2zLA7POUpUuJhck9IOmGt6VBR65Cek9+VVT7B+gM/G4gE7sH6Rn9RxpgsERkO/Ig1hHW8MWaLiLwORBljFuAo+iciI7A6rIcYY4yItANeF5FMIAd42BhzuT4P5UGS0rN4bs5GapYL5Jnb6llPilhr/SbGwfzHoERFqHWzvYEq95Z6xmpWit8Ld8+CkEi7I3IpYszlm+5FJMoYEykiG40xTRzPrXcMT3UJkZGRJipKW6A8xcvzNjP1rwPMfqg1kaHnLQKUlgDju8CZg1ZzQOUm9gSp3FvaWSs5HNlo1QGre5vdEdlCRKKNMXlmxvwOc01xNBPFiMg7jm/83jWlUBWalXtOMmX1Ae5rG3ZhcgBr2OHds8G/FHxzZ4EsWaq8THoSfNMXjmyw6n95aXK4nPz+kR+M1Uw0HEjG6nzu46yglPdKTs/iX99uJDQ4gJHnmpbyElQV7p5jdSpO7WOt8KVUfmSkWKsYxkZBn6+hfje7I3JZ+S3Wd8AYk2qMOWuMec0Y87QxZrezg1Pe550l24k9nco7fZtS3O8yE+IqNoT+38DpfTB9IGSmFU6Qyn1lpsGMAbD/T2u4tBes6XAtLlesbxN5FOk751x/hFIFYfXeU0xadYChbUNpGZZH01Jewm6EnmPh2/vhu4eg7wSvK6im8ikrHWYOgr2/WwtUNbnT7ohc3uVGMd1eKFEor5eSYY1aqhEcwLOdLtG0lJfGfeFsHPz0MvwQbK32pRPpVG7ZmTB7qFVG/o4PIXyg3RG5hfxMlFPK6d5ZsoOD8SnMHNaKAL/8jr7Opc3jkHzCWtSleBno8HLBB6ncU1a6lRx2fG+tNRIxxO6I3Ea+fhNFJJH/b2ryA3yBZGNMqYsfpVT+rN0fz6RV+7m3dQ2urxl8dS8iAh1ft4bALh8NxUtbSUN5t8w0qyrrrqVWcsi1toi6vHwlCGNMyXP3RUSwqrK2clZQynukZmTz3JyNhJQpznOd61/bi4nA7R9YSWLpS9Zw2Ob3FEygyv1kpsKMgbDnV7h9jLUQlboiV9ybZyzzgE4FH47yNu8t3cG+k8n8r08TAotdRdPS+Yr4QO8voVYHWPgkbJl37a+p3E9GMkzrB3uWQY9PNTlcpfw2MfXO9bAIVrkNHVOorkn0gdN8vWIfd19fnTa1CrCsclE/uGuKVbb52wesdSVqdyi411euLT3JSg4HV0Gvz6Fpf7sjclv5vYK4I9etE1a57/MX/1Eq39Iys3l2zgaqBBXnha4NCv4N/AJh4EwoXx9m3G0NbVSeL+2sNXHy4GrrSlKTwzXJbx+EXp+pAvXBzzvZeyKZKfe3pERBNC3lpXhpq67/xNutFcLunm3Nm1CeKfmklRyObYa+43USXAHIb7nvmiKyUEROiMhxEZnvKMOt1BWLOXSGL//YS/8W1bixTnnnvllgObh3IZSpYTU77F/h3PdT9kiIhfGd4cR2q/CeJocCkd8mpmnALKAyUAWr7Pd0ZwWlPFd6VjbPzt5AxVL+vNjNCU1LeSlR3koSQSFWcb8DqwrnfVXhOLkLvu4EScdg8HdQV8fPFJT8JogAY8wUY0yW4zYVa3Efpa7IR7/sYtfxJN7s3ZhS/r6F98YlKlhJolRlq4rnwb8K772V88TFWFcO2ekwZBHUaGN3RB4lvwlisYg8LyKhIlJDRJ4DfhCRsiKSz6I5ytttik3g89/30jcihJvrVSj8AEpWgnsXWQsNTe1tFWxT7mv/Cqt/ybc4DF0ClZvaHZHHye+CQfsusdkYY2zvj9AFg1xbRlYO3T/5k/jkDH4acRNBAYV49XC+s3HW+sNnDsJd30CdW+2LRV2dLfNg7jCrb2nwPKv8u7oq17xgkDEm7BI325ODcn2fLNvN9qOJvNmrsb3JAaBUFRi6GMrVtdYF0Ml07sMYWPkJzB4CVcKtKwdNDk6T31FMviLyhIjMcdyGi4jNv+XKXWyJS+CzZbvp1awqtzasaHc4lnOjm6o2hzlDIWaa3RGpy8nJhsXPwdJ/Q8PucM98CLzK2l0qX/LbBzEWiAA+c9wiHM8pdUmZ2Tk8O3sjpQP8ePWOhnaH80/FS1ujXsLawbxH4K9xdkekLiYjBWYOhjVfQOvh0Hei1fegnCq/M5RaGGNy9wD9KiIbnBGQ8ixjf9vD1iNn+XxQBKUD/OwO50J+gTBgJsy5z/p2mngEbnlFFx1yJYlHraJ7h9dBl3fg+ofsjshr5Pe3IFtEap174Jgkl+2ckJSn2H70LB//uos7mlahc6NKdodzcb7+0G8yRN4Hf34A3w2z1hBQ9jscDV/cDMe3WcvLanIoVPm9gngWWCYiex2PQwEtv6EuKis7h+fmbKSUvy+vdb/O7nAuz6codHsfgqrBL69Z31rvmmo1Qyl7bJwFCx6HwApw/1Ko1NjuiLxOfq8gVgDjgBwg3nFfp6Oqixr72x42xibweo9GlA10waalvIjAjU9bRd4OrrYmYJ05ZHdU3icnG356BeY+CFUjYdgyTQ42yW+CmAyEAf8BPgZqAlOcFZRyb1viEvjI0bTUrUllu8O5ck36waBv4exh+PJmLc1RmFLircKKKz6EyPutYouBBVgKXl2R/CaIRsaYB4wxyxy3BwE3aDdQhS0jK4dnZm2gdIAfr7tD09LF1LwJHvjFWpVu0h0QNd7uiDxfbBSMawf7frdWBrz9ffDR0fR2ym+CWCcify8xKiLXAzptWV3go192sf1oIm/1akwZd2laupjyda0kUbM9LBph3bIy7I7K8xgDq8daTXoicN+P1oABZbv8dlJHACtF5KDjcXVgh4hswiq10cQp0Sm3sv7gaT77bTd3RoS4zoS4a1W8tLXw0C+vw4ox1miaPl/r7N2CknoGFgyHbQuhXlfo+RkUL2N3VMohvwmis1OjUG4vLTObZ2ZvoFIpf152tQlx16qID3R8zeooXfAEfH4D9BoHdW+zOzL3tm85fPewNffktv9aE+BE7I5K5ZLfFeUOODsQ5d5G/7iDvSeSmXr/9YVbxrswNe5rVQydPQSm3QltnoAOr2g7+ZXKSodf/2PVVCpb0xrCGpJnrThlM50uqq7Zmn3xfL1iH4Nb1eCGOh4+4qRcHXjgZ6uNfOVHMKELnNpjd1Tu49gW+PIWWPkxRAyBh5drcnBhmiDUNUlOz2Lk7A1UKxPA813q2x1O4fAtbo2y6TsBTu6EsW2tOk45OXZH5roy0+DXN2DcTdbKbwNnwR1jrFInymVpglDX5K3F2zh0OoXRdzYlsFh+u7Q8RKPe8OhqCL3BquM0uTuc1tbYC+xfYfXb/PHO//+f6bKgbkEThLpqy3edYOrqg9zfNoyWYV66sGCpKnD3bLjjI4hbD5+1tppPsjPtjsx+ySetUhkTu1pLgg76Fnp/oRPf3IgmCHVVElIyeW7ORmqVD2Rkp3p2h2MvEYi4Fx5ZCaFtYelL1oQvb52BnZVhdUB/1BzWf2ONTnp0NdTWlfvcjVMThIh0FpEdIrJbRJ7PY3t1EVkmIutFZKOIdM217QXHcTtERK9HXczL8zdzPDGd9/uF4+/rY3c4rqFMDatt/a5vIO0sTOhsDeNMiLU7ssJhDOxYDGNbW4v6VGsBj66CTm9oX4ObclqjsYj4AJ8CHYFYYK2ILDDGbM2120vALGPMWBFpCPwAhDru98cq51EF+FlE6hpjtMS4C5gfc5gFG+J4pmNdmlYrbXc4rkUEGtwOtW6GP96FVZ/C5rlWmeobn/bMSWDGwN7fYNkbELsWguvAwNk6T8QDOPMKoiWw2xiz1xiTAcwAepy3jwFKOe4HAXGO+z2AGcaYdGPMPmC34/WUzWJPp/DSvM1E1CjDI+1rXf4Ab+UXCLeOgsejrY7ZlR/Dh03hj9GQlmB3dAXDGGuy28RuMKUnnD0Ct4+xmto0OXgEZw47qQrkrpUcC1x/3j6jgKUi8jgQCJxrpKwKrD7v2AtqG4jIMGAYQPXq1QskaHVx2TmGZ2ZtICfH8EG/cIr6aBfWZZWuDr0+t9rhf3ndmiC24kNocT+0ehRKVLA7wiuXnQlb58OqT6yO+RKVoMu7Vj9M0WJ2R6cKkN3jEgcAE40x74lIa2CKiDTK78HGmC+ALwAiIyONk2JUDl8u38tf++J5t28TqgcH2B2Oe6nUCO6eBXEx1qp1f46BVZ9Zs7MjhlqTxVy9zMTZONgwHdaOh7OxVlPS7R9A0wG6PrSHcmaCOAxUy/U4xPFcbvfjqPNkjFklIv5AuXweqwrR5sMJvLd0B10aVaJvRIjd4bivKuHQb5I1+3rlx7BpNsR8AxUbW9/AG/aEEuXtjvL/ZabBrh9h/VTY/TOYHAhrZ5Xirt1R1+72cGKMc754i0hRYCfQAeuP+1pgoDFmS659FgMzjTETRaQB8AtWU1JDYBpWv0MVx/N1LtVJHRkZaaKitAK5M6RlZnP7x39yNjWTH59q5/5lvF1JeqKVJNaOh2ObQIpA6I1wXS+o1wVK2rCWd9pZ2P2TVWF151LITIaSVSB8oHUL1r4nTyIi0caYPOudOO0KwhiTJSLDgR8BH2C8MWaLiLwORBljFgDPAF+KyAisDushxspYW0RkFrAVyAIe0xFM9nl78XZ2H09i8n0tNTkUtGIlrbpOEUPh2GbYMg+2zIVFT1m3Cg2h5s3WAkZVmjvn6iL1DMStszqc9y+Hw+vAZFtrQTfpBw3usNbEKKLDmb2N064gCpteQTjHbzuOM2TCWoa2DeXVO9x4hTh3YoyVLHb/bA0fPbDKmokMEFTNqihbvh6UCYOyYRAUAsXLWskmr34MYyAjGVLjrTW2zxyEMweswnlHN8Lp/dZ+RYpaSSisnTWprVpLTQpe4FJXEJog1EUdT0yj64fLKRvox4LhN+iEOLtkplrLcR6JsUYNxcVYf9TPv6gWH2uJ1CJFrUQhRSArzWoyyusCvGxNqNQEKjexkk61VlCsRCGckHIltjQxKfeW4xjSmpSexbQHW2lysJNvcQi70bqdk51pzdA+vQ8SDkPaGUg9bTUXmWzrqgEDPn5W0vAPAv/SULoalK5hXXXokFR1GZogVJ7G/bGX5btO8lbvxtStWNLucNT5fHyt5qWyYXZHojyYjlFTF4g+cJrRS3fQrXFl+reodvkDlFIeSROE+oeE1EyemL6eykH+vNm7MeLqk7eUUk6jTUzqb8YYXpy7iWNn05j1cGuCiutay0p5M72CUH+bvuYQ3286wjO31aN5dQ+sOqqUuiKaIBQA24+e5bWFW7ixTjkealfT7nCUUi5AE4QiMS2TR6auo1RxX97r15QiRbTfQSmlfRBezxjDc3M2cjA+hWkPXE+Fkv52h6SUchF6BeHlvv5zH4s3H+W5TvW4vmaw3eEopVyIJggvFrU/nrcXb+e2hhUZpv0OSqnzaILwUieT0nls2jqqlinOu3c21fkOSqkLaB+EF8rOMTw5Yz1nUjL57tGWOt9BKZUnTRBeaPTSHazYfYp3+jahYZVSdoejlHJR2sTkZRZuiGPsb3sY0LI6/SK1zpJS6uI0QXiRzYcTeHbOBlqEluG17rr4j1Lq0jRBeImTSek8NCWasgF+fHZ3BH5F9aNXSl2a9kF4gYysHB6duo6TSel8+0gbypfUhWKUUpenCcILvL5oC2v2x/Nh/3AaVQ2yOxyllJvQdgYPN2XVfqauPsjDN9WiR3hVu8NRSrkRTRAebNn247y6YAu3NqjAs53q2R2OUsrNaILwUFviEhg+bR0Nq5Tiw/7N8NEKrUqpK6QJwgMdSUjlvolrCSruy/h7WxBYTLualFJXTv9yeJjEtEyGTlhLcno2cx5pTYVSWr5bKXV1NEF4kMzsHIZPW8+u40lMGNKC+pW0jIZS6uppE5OHyMkxjJy9gd93nuCNno1oV7e83SEppdycJggPYIzhtYVbmB8Tx3Od69G/ZXW7Q1JKeQBtYvIAH/2ym0mrDvDgjWE8clMtu8NRylaZmZnExsaSlpZmdyguxd/fn5CQEHx981/eXxOEm5u8aj8f/LyTvhEhvNi1gS78o7xebGwsJUuWJDQ0VH8fHIwxnDp1itjYWMLCwvJ9nDYxubE50bGOiXAVebt3Y/1lUApIS0sjODhYfx9yERGCg4Ov+KpKE4SbmrsulmfnbKBtrXJ8MrAZRX30o1TqHE0OF7qa/xP9q+KG5q0/zDOzN9CmVjBf3hOJv6+P3SEppTyQJgg3Mz/mME/PiqFVWDBf3dOC4n6aHJRydxMnTmT48OF5bmvTpg0A+/fvZ9q0aRd9jUmTJlGnTh3q1KnDpEmTCiQuTRBu5NvoWEbMjKFlWFm+HhKpyUEpL7By5Urg0gkiPj6e1157jb/++os1a9bw2muvcfr06Wt+bx3F5CbG/7mP1xdt5Yba5Rg3OIIAP/3olLqc1xZuYWvc2QJ9zYZVSvHqHZdesveNN95g0qRJVKhQgWrVqhEREcHIkSNp3749o0ePJjIykpMnTxIZGcn+/fsBOHToEO3bt+fw4cMMGjSIV199FYASJUqQlJTE888/z7Zt2wgPD+fee+9lxIgRf7/fjz/+SMeOHSlbtiwAHTt2ZMmSJQwYMOCaztWpf2VEpDPwIeADfGWMefu87R8ANzseBgAVjDGlHduygU2ObQeNMd2dGaurMsYw5uddfPjLLjpfV4kPB4RTrKheOSjlqqKjo5kxYwYxMTFkZWXRvHlzIiIiLnvcmjVr2Lx5MwEBAbRo0YJu3boRGRn59/a3336b0aNHs2jRoguOPXz4MNWqVfv7cUhICIcPH77mc3FaghARH+BToCMQC6wVkQXGmK3n9jHGjMi1/+NAs1wvkWqMCXdWfO4gJ8fw+qKtTFy5n74RIbzdu7GOVlLqClzum74zLF++nF69ehEQEABA9+75+27bsWNHgoODAejduzd//vnnPxKEHZz516YlsNsYs9cYkwHMAHpcYv8BwHQnxuNWUjOyefSbdUxcuZ/72obxTp8mmhyUcnNFixYlJycH4II5CecPQ72SYalVq1bl0KFDfz+OjY2latVrX0HSmX9xqgKHcj2OdTx3ARGpAYQBv+Z62l9EokRktYj0vMhxwxz7RJ04caKAwrbf8cQ0+n+xih+3HuWlbg14+fYGFNEFf5RyC+3atWPevHmkpqaSmJjIwoUL/94WGhpKdHQ0AHPmzPnHcT/99BPx8fGkpqYyb9482rZt+4/tJUuWJDExMc/37NSpE0uXLuX06dOcPn2apUuX0qlTp2s+F1f5StofmGOMyc71XA1jTCQwEBgjIhcUGTLGfGGMiTTGRJYv7xnVS3ceS6TXpyvZeSyJcYMieODGmjrpRyk30rx5c+666y6aNm1Kly5daNGixd/bRo4cydixY2nWrBknT578x3EtW7akT58+NGnShD59+lzQvNSkSRN8fHxo2rQpH3zwwT+2lS1blpdffpkWLVrQokULXnnllb87rK+FGGOu+UXyfGGR1sAoY0wnx+MXAIwxb+Wx73rgMWPMyou81kRgkTFmTl7bASIjI01UVFRBhG6b7zce4bk5GwgsVpSv721B45Agu0NSyu1s27aNBg0a2B3G30aNGkWJEiUYOXKk3aHk+X8jItGOL+MXcOYVxFqgjoiEiYgf1lXCgvN3EpH6QBlgVa7nyohIMcf9ckBbYOv5x3qKzOwc/rtoK49NW0e9SiVZMPwGTQ5KKds5bRSTMSZLRIYDP2INcx1vjNkiIq8DUcaYc8miPzDD/PNSpgEwTkRysJLY27lHP3mSowlpPDFjPWv2xTOkTSgvdm2AX1FXaflTSl2rUaNG2R3CVXPqPAhjzA/AD+c998p5j0flcdxKoLEzY3MFP2w6wgtzN5GRlcOYu8Lp2ezaRx0opVRB0em4NkhMy2TUgq18uy6WpiFBfHBXODXLl7A7LKWU+gdNEIVs6ZajjFqwhaNn03iiQx0ev6U2vjq/QSnlgjRBFJIjCam8On8LS7ceo36lknxyd3OaVy9jd1hKKXVR+tXVyZLTsxjz8046vPc7f+w6wb8612fh4zdoclDKQ506dYrw8HDCw8OpVKkSVatW/fvxufIbue3YsYP27dsTHh5OgwYNGDZsWJ6v64xy3pejVxBOkpGVw5zoWD74eScnEtPp2rgSz3duQPXgC39AlFKeIzg4mJiYGODCORAlSlzY1/jEE08wYsQIevSwKhFt2rTpgn3OlfOOiopCRIiIiKB79+6UKePcL5qaIApYSkYWM9Yc4svlezmSkEZkjTKMGxyhVwxK2WHx83D0wj+416RSY+jy9uX3y6cjR44QEhLy9+PGjS8cwOmsct6XowmigOw9kcTMtYeYFXWI0ymZtAwry1u9G3NT3fJaKkMpdVEjRozglltuoU2bNtx2220MHTqU0qVL/2MfZ5XzvhxNENfgRGI6P209xvyYw/y1Lx6fIsKtDSowrF1NImpcex0UpdQ1KsBv+s4ydOhQOnXqxJIlS5g/fz7jxo1jw4YNFCtWzO7QNEFcibTMbGIOnWHNvniW7zpB1IHTGAOhwQE826ked0aEUKGUv91hKqXcTJUqVbjvvvu47777aNSoEZs3b/7HIkNVq1blt99++/txbGws7du3d3pcmiDykJiWydGENOIS0th1LJEdRxPZcSyR7UcTycjKQQQaVCrFkx3q0LlRJepVLKnNSEqpq7JkyRI6dOiAr68vR48e5dSpUxes5dCpUydefPHFv9eZXrp0KW+9dUHd0wLn9QkiPjmDu8atIi0rm9SMHFIyskjJyP7HPuVK+FGvUkmGtAmlZWhZWoSWJSjA16aIlVLuKiUl5R8d0k8//TSxsbE8+eST+PtbrQ/vvvsulSpV+sdxuct5AwVWzvtynFbuu7Bdbbnv5PQsRs7egL+vD/6+PhT39aFCqWJUDvKnclBxapYPpFwJ+9sClVL542rlvl3JlZb79voriMBiRRk76PILiiullLfRmdRKKaXypAlCKeVxPKXpvCBdzf+JJgillEfx9/fn1KlTmiRyMcZw6tSpvzvC88vr+yCUUp4lJCSE2NhYTpw4YXcoLsXf3/8fI6jyQxOEUsqj+Pr6EhYWZncYHkGbmJRSSuVJE4RSSqk8aYJQSimVJ4+ZSS0iJ4AD1/AS5YCTBRSOu/C2c/a28wU9Z29xLedcwxhTPq8NHpMgrpWIRF1surmn8rZz9rbzBT1nb+Gsc9YmJqWUUnnSBKGUUipPmiD+3xd2B2ADbztnbztf0HP2Fk45Z+2DUEoplSe9glBKKZUnTRBKKaXy5PUJQkQ6i8gOEdktIs/bHU9hEJH9IrJJRGJE5MqX4XMDIjJeRI6LyOZcz5UVkZ9EZJfj3zJ2xljQLnLOo0TksOOzjhGRrnbGWNBEpJqILBORrSKyRUSedDzvkZ/1Jc7XKZ+zV/dBiIgPsBPoCMQCa4EBxpittgbmZCKyH4g0xnjsZCIRaQckAZONMY0cz70DxBtj3nZ8GShjjPmXnXEWpIuc8yggyRgz2s7YnEVEKgOVjTHrRKQkEA30BIbggZ/1Jc63H074nL39CqIlsNsYs9cYkwHMAHrYHJMqAMaYP4D4857uAUxy3J+E9YvlMS5yzh7NGHPEGLPOcT8R2AZUxUM/60ucr1N4e4KoChzK9TgWJ/5nuxADLBWRaBEZZncwhaiiMeaI4/5RoKKdwRSi4SKy0dEE5RFNLXkRkVCgGfAXXvBZn3e+4ITP2dsThLe6wRjTHOgCPOZomvAqxmpb9Yb21bFALSAcOAK8Z2s0TiIiJYBvgaeMMWdzb/PEzzqP83XK5+ztCeIwUC3X4xDHcx7NGHPY8e9x4DuspjZvcMzRhnuuLfe4zfE4nTHmmDEm2xiTA3yJB37WIuKL9cfyG2PMXMfTHvtZ53W+zvqcvT1BrAXqiEiYiPgB/YEFNsfkVCIS6OjcQkQCgduAzZc+ymMsAO513L8XmG9jLIXi3B9Jh1542GctIgJ8DWwzxryfa5NHftYXO19nfc5ePYoJwDEcbAzgA4w3xrxhb0TOJSI1sa4awFpydponnrOITAfaY5VBPga8CswDZgHVsUrD9zPGeEyn7kXOuT1Ws4MB9gMP5Wqbd3sicgOwHNgE5DiefhGrXd7jPutLnO8AnPA5e32CUEoplTdvb2JSSil1EZoglFJK5UkThFJKqTxpglBKKZUnTRBKKaXypAlCqaskIqVF5FHH/SoiMsfumJQqSDrMVamr5KiFs+hc5VSlPE1RuwNQyo29DdQSkRhgF9DAGNNIRIZgVQ8NBOoAowE/YDCQDnQ1xsSLSC3gU6A8kAI8aIzZXtgnodTFaBOTUlfveWCPMSYcePa8bY2A3kAL4A0gxRjTDFgF3OPY5wvgcWNMBDAS+KwwglYqv/QKQinnWOao158oIgnAQsfzm4AmjmqcbYDZVnkdAIoVfphKXZwmCKWcIz3X/Zxcj3Owfu+KAGccVx9KuSRtYlLq6iUCJa/mQEcN/30icidYVTpFpGlBBqfUtdIEodRVMsacAlaIyGbg3at4ibuB+0VkA7AFXe5WuRgd5qqUUipPegWhlFIqT5oglFJK5UkThFJKqTxpglBKKZUnTRBKKaXypAlCKaVUnjRBKKWUytP/AaWxgCgI0/J8AAAAAElFTkSuQmCC\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", " chi_constructor: chis_qubit\n", " mu: derivative_wrt_pulse\n", " S(t) (ranges): [0.000000, 1.000000]\n", " iter_start: 0\n", " iter_stop: 5\n", " duration: 0.3 secs (started at 2020-08-17 10:41:30)\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.5 secs (started at 2020-08-17 10:41:31)\n", " optimized pulses (ranges): [0.00, 2.06]\n", " āˆ«gā‚(t)dt: 7.64e-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.9 secs (started at 2020-08-17 10:41:33)\n", " optimized pulses (ranges): [0.00, 2.23]\n", " āˆ«gā‚(t)dt: 5.71e-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.9 secs (started at 2020-08-17 10:41:36)\n", " optimized pulses (ranges): [0.00, 2.33]\n", " āˆ«gā‚(t)dt: 8.19e-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.9 secs (started at 2020-08-17 10:41:39)\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: 3.0 secs (started at 2020-08-17 10:41:42)\n", " optimized pulses (ranges): [0.00, 2.15]\n", " āˆ«gā‚(t)dt: 6.40e-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", "opt_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": { "attributes": { "classes": [], "id": "", "n": "20" } }, "outputs": [ { "data": { "text/plain": [ "Krotov Optimization Result\n", "--------------------------\n", "- Started at 2020-08-17 10:41:30\n", "- Number of objectives: 1\n", "- Number of iterations: 5\n", "- Reason for termination: Reached 5 iterations\n", "- Ended at 2020-08-17 10:41:45 (0:00:15)" ] }, "execution_count": 20, "metadata": {}, "output_type": "execute_result" } ], "source": [ "opt_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": { "attributes": { "classes": [], "id": "", "n": "21" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEGCAYAAABo25JHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAArAElEQVR4nO3deXxcdb3/8dcnySSTfW3TpltaaGlLoaUtLfvFCyigggIiXERAARdQrqBewAXlp163iyKiWBWhiOxbQUDZRJaytKV7S+lC9zZp0uyZJJP5/v6YSQlttrY5cyaZ9/PxyCMzZ07mfE4Hzme+3+/5fr7mnENERJJXit8BiIiIv5QIRESSnBKBiEiSUyIQEUlySgQiIkkuze8A9ldJSYkrLy/3OwwRkQFl4cKFu5xzQ7p6bcAlgvLychYsWOB3GCIiA4qZbezuNXUNiYgkOSUCEZEkp0QgIpLklAhERJKcEoGISJJTIhARSXJKBCIiSW7AzSMQqKgPsWRzLTvqQlQ3tBJIM/IzA4wpymbCsByG5gb9DlFEBhAlggGiqTXMI4u2cv9bm1ixra7HfceVZDN7XBHHjCvm2HHFDM1TYhCR7ikRJLhIxPHQws38/Nl3qWps5YgR+XzrY4dxzLhiRhZmUpSdTnvEUd3YyvtVjazYWsebG6p4aul27ntrMwBTR+bz0cOH8bHDSzlkSA5m5vNZiUgisYG2QtnMmTNdspSYqAu1ce0DS3h+1U5mjinkf86YyMwxhX26kLdHHKu21/Hymkr+uXInSzbXADC2JJtTJw3l1EmlzBhTSFqqholEkoGZLXTOzezyNSWCxLSjNsTFf36TDbsa+c7HJ3HpceUH9U1+e20zz6/cyXOrKpi/bhdt7Y6CrAD/edhQTp1cyonjS8gNBvrxDEQkkSgRDDBba5r57B/mU9PUxpzPz+C4Q0r69f0bWsK8sqaS51bt5MXVFdQ0tRFINaaOLGDW2CKOHlvEjDGF5CkxiAwaSgQDSE1TK+fdMZ+ddSHuvXw2R44s8PR44fYIizbV8MLqnby5vprlW2sJRxwpBpOG53HEiHwml+VxeFkeE4flkZ2hYSWRgainRKD/qxNIW3uEK+cuZFNVE3O/OMvzJACQlprCrLFFzBpbBETvTnpnUw1vbahmwcZqnl2xg/vfjg46m0F5cTaThucyriSHMcVZjC3JZkxxNiU56RqEFhmglAgSyP8+vZq33q/m1gumccy4Yl9iyEpP4/hDSzj+0Gh3lHOO7bUhVm6rY+X2OlZuq2PFtjr+sWIn7ZEPWpPZ6amMKc6mrCDIsPwgw/KClOYFGZ6fybD8DEpyMsgLBkhJUbIQSTRKBAniqaXbuPO1DVx6XDlnTxvhdzh7mBllBZmUFWRy6uTSPdvb2iNs3d3MhqpGNu5q5P2qJjZWNbJldzMLNu6mpqltn/dKTTEKs9Ipzk6nKDudopzo4+LsDAqyAuRnRn/yMtNivwPkBQMEA6nxPGWRpKNEkAC21TRz/SPLmDGmkBvPnOR3OH0SSE2hvCSb8pJsOGzf10Nt7eyoDbGjLsSO2hC7GlqobmylurGVqtjvldvqqGpooS4U7vFYGWkpnZJE7HcwbU+iyMtMi/3e93luMI2AbpEV6ZESgc+cc9zw6DLaI45ff3Ya6WmD46IVDKR+kCh60dYeoba5jbrmNmpjP3Wh8Ie2dX6tsr6FdZUN1MX269xF1ZWs9NQPJYgPWhu9J5PcYJrmWsig51kiMLNRwFygFHDAHOfcrXvtY8CtwJlAE3Cpc26RVzElokcWbeXlNZX84JOTGVWU5Xc4vgikplCSEx1H2F/OORpb22NJoY265nCnx9FE8aHXQm3srA/xXkXDnn16ySNkp6d2kzDSOrVQAhRlp1OaF6Q0P4Pi7AxSNR4iA4SXLYIwcJ1zbpGZ5QILzew559zKTvucAYyP/cwGfh/7nRR21oW4+ckVHF1eyOePLfc7nAHJzMjJSCMnI40yMvf77yMRR2Nr+IOE0U3y6Px8R12INRX1e17r6g7s1BRjSE4GpflBSnMzGJYfpKwgkzFFWYwuzmJMcTY5uhVXEoRn/yU657YD22OP681sFTAC6JwIzgbmuuhkhjfMrMDMhsf+dlBzzvGdx5bTEo7w8/Om6m4an6SkGLnBALnBACMKDiyRNLSGqW1qo7qxlR11ISrqomMjO+ta2FkX4v2qRt5YX7XPWEhxdjqji7MYW5zNhGG5HFaay2HDchmeH9StuBJXcflKYmblwFHAm3u9NALY3On5lti2DyUCM7sSuBJg9OjRnsUZT/OWbOP5VTv5zpmTGNuHfnRJTCkpFu0qCgYYVZTF1B72rQu1samqiY1VTWysbtzz+PV1VTz6ztY9++UG0zisNJcpI/I5anQBR40qZFRRppKDeMbzRGBmOcAjwH8753qun9wN59wcYA5EZxb3Y3i+2NXQwg/mrWDaqAK+cMJYv8OROMkLBpgyIp8pI/L3ea22qY01FfWs3lHPmh31rN5RxwNvb+au198HoCg7nWmjCpg9tojjDy1h8vA8tSKl33iaCMwsQDQJ3Ouce7SLXbYCozo9HxnbNqj98MmVNLa084vzjtSAogCQnxXg6PIiji4v2rMt3B5hzc4G3tm8m8Wbali4aTcvrq4AoDArwHGHlHDShBJOmVR6QAPtIh28vGvIgD8Dq5xzt3Sz2zzgajO7n+ggcW0ijg9s2NVIY0u4y29y++u5lTt5csk2rjttAuNLc/shOhms0lJTmFyWx+SyPC6aPQaI3mDw+rpdvPpeFa+t3cXfl20nxZYxc0wRHz28lI8dPixp7z6TA+dZ0TkzOwF4BVgGRGKbbwRGAzjn7ogli98CpxO9ffQy51yPFeXiXXSuoi7ER375L5ra2nn4y8cxY0zhAb9XXaiN0255mcKsdOZdfcKgmTMg/nDOsXJ7tNzHP1fsYPWOegCOGVfEeTNGccaUYSoSKHuo+uhBmPPvdfzk6dUAHDW6gEe/ctwBD9rd8OgyHnh7E4999XimjiroxyhF4P1djTy5ZBsPL9rCxqomstJT+fgRw7nkuPJ+ac3KwNZTItBX0l68uLqCicNy+fm5R/LOphqeXrbjgN7ntbW7uO+tTVx+4jglAfFEeUk2XztlPP/65sk8/OVjOWtqGX9ftp1P3PYq5/9hPv9YsaPXWdiSnJQIehCJOJZuqWX22CLOnTGSicNy+dmzq2kJt+/X+1Q1tPCNBxZzyJBsvnHqBI+iFYkyM2aWF/HTc49k/g2n8J0zJ7F1dzNfumchp/zfv3h00RYlBPkQJYIebKpuoqm1nclleaSmGDeeOYlN1U3cM39jn9/DOcf/PLKUmqY2brtwOpnpqqQp8ZOfGeCKk8bx8rdO5ncXTSczPY1rH1zCab96mScWb1VCEECJoEcrt0enPUwangfASROGcNKEIdz24lqqGlr69B6/+9c6nl9VwfVnTGRyWZ5nsYr0JC01hTOPGM7fv3YCd3xuOoGUFK65fzGfuO1V3lhf5Xd44jMlgh6sr2wAYPzQD27z/O7HJ9Hc1s5/P7C4129Tzyzbzi/+8S5nTS3jsuPLvQxVpE9SUozTpwznmWtO5DcXHkVdcxsXzHmDq+5dxJbdTX6HJz5RIujB5upmSnIyPtSdM6E0lx+edTivvLeLHz65gu7uunp2+Xa+dt87TB9dwM/PO1LlASShpKQYZ00t44Xr/oNrT5vAC6t3csr/vcwfXl5HuD3S+xvIoKJE0IPNu5sYVbRvIbILjh7FlSeNY+78jVx93zvUNLXuea2xJczPnl3NV+5dxBEj87nrC7O0wpYkrGAgla+fMp4XrzuZkyYM4X+fWc25v3+dd2NzEiQ5aLZJD7bsbu7yVk8z44YzJpKfGeCW59bw0uoKjjukGDPjrQ3V1Da3cf7Mkdx89hQlARkQygoymXPxDJ5aup2b5q3gE7e9wjWnjOcrJx+qMihJQImgG+0Rx7aaZj5x5PAuXzczrvrIoZwyaSh/efV9Fm3ajQNOmTiUi44Zc1AzkEX8YGZ8cmoZxx1SzE3zVvDLf67h9XVV/Pqz0xiaF/Q7PPGQEkE3dtSFCEccIwt7rtsycVgePzvvyDhFJeK94pwMbrvwKE6aMISbnljBGbe+wi2fncZ/TBjid2jiEY0RdGNHbQiA4fn6JiTJx8w4f+Yo5l19PCU5GVxy51vc8twaIpp3MCgpEXSjsj46T2BIrsr7SvIaX5rLE1cfz3kzRvKbF97jq/cuorEl3PsfyoCiRNCNytiEsaFKBJLkgoFUfnHekXz345P458odnPv71zXnYJBRIuhGZX0LZtGVoUSSnZlx+Ynj+Mtls9ha08zZv32NZVtq/Q5L+okSQTcq61sozk4nLVX/RCId/mPCEB6/6niCgVQumDOf19bu8jsk6Qe6ynWjsj6k5f9EunDIkBwe/epxjCzM4rK/vM1TS7f5HZIcJCWCblTWt+jeaZFulOYFefBLxzJ1VD5fu+8dHnx7s98hyUFQIuhGZX0LQ9QiEOlWflaAe744mxMOLeHbjyzlb29u8jskOUBKBF1wzlHZ0EJJrgaKRXoSDKTyx8/P5COHDeHGx5Zxzxt9X6tDEocSQRfqW8K0tTtKstUiEOlNMJDKHRfP4NRJQ/ne48v5q5LBgKNE0IWaxjYACrICPkciMjBkpKXyu4tmcMrEoXzvieU8/s5Wv0OS/aBE0IWa5mhZ6YIsdQ2J9FV6Wgq3XzSd2WOLuO6hJTy/cqffIUkfKRF0YXdTtEVQqBaByH4JBlL50yVHM6Usj6/+bRGvr9M8g4FAiaALHQvNqEUgsv9yMtK467JZjCnK4oq7F2gG8gCgRNCF3Y3RRKAWgciBKcxO56+Xz6YgK50v3P22ahPFzJ3/Pjc9sZxQW7vfoXyIEkEXOrqG8jOVCEQOVGlekL9cdjShtna+cNfb1Da3+R2SrxZvruH7T6zg7vkbuev19/0O50OUCLpQ29xGbjBNdYZEDtKE0lz+8LkZbNjVyFf+upDWcMTvkHwRiThuemI5Q3IzOLwsj4cXbvE7pA/Rla4Lu5taKdT4gEi/OO7QEn56zpG8vq6KGx5dhnPJt7jNQws3s2RLLTeeOZGzp5WxtqKBivqQ32HtoUTQhd1NbRofEOlH584YyTdOncAji7bwmxfW+h1OXNU2t/HzZ9/l6PJCPjVtBEeXFwGw4P3dPkf2ASWCLtQ0tZKvFoFIv/r6KYdyzvQR/Or5NTyzbLvf4cTNr59fw+6mVn5w1uGYGZOG55GaYqzaXud3aHsoEXQh2jWkFoFIfzIzfvLpI5g2qoBrH1ySUBdCr7y3s5658zdywazRHF6WD0TnWowpzmLNznqfo/uAEkEXapraNEYg4oFgIJU5F88gLzONK+YuoDp2q/Zg5Jzj5qdWkpWeynWnTfjQa4eV5rJmZ4NPke1LiWAv4fYI9aGw6gyJeGRoXpA5F8+kor6Fr967kLb2wXkn0fOrKnjlvV1849QJFO9V0n780Bw2VjUmzF1USgR7qYnd61ygOQQinpk6qoCfnXsEb6yv5v89tdLvcPpdS7idH/19JYcOzeHiY8fs8/rIoiwiDrbXNvsQ3b6UCPZSHwoDkKdEIOKpTx81kitPGsfc+Ru5/63BtajNna++z8aqJr7/ickEupiPNKowC4Atu5UIElJDLBHkZKT5HInI4Pc/p0/kxPElfO+J5SzcWO13OP2ioi7Eb198j1MnlXLShCFd7jOyMBOAzdWJUXpDiWAv9S3RrqGcoBKBiNdSU4zfXjidEQWZfOmeRQnTVXIwbn5qJW0Rx3c/PqnbfYbnB0lNscHfIjCzO82swsyWd/P6yWZWa2aLYz/f9yqW/dHYEi0GlZuhriGReMjPCvDHz8+kuTXMl+5ZmHAF2fbHy2sqeWrpdq46+VDKS7K73S8tNYXh+UE2J0gxPi9bBHcBp/eyzyvOuWmxn5s9jKXPGtQiEIm78aW5/PqCo1i6pXbAlqEItbXzvceXM64kmy+fPK7X/csKMtlekxhlJjxLBM65fwMDrtNPYwQi/jhtcinXnjaBx97Zyp9f3eB3OPvt9pfWsqm6iR99egoZaam97l+aF0yYekN+jxEca2ZLzOwZMzu8u53M7EozW2BmCyorKz0NqL4lmghy1SIQiburP3IoZ0wZxk+eXsW/13j7/3p/WlvRwB0vr+Oco0Zw3CElffqbobkZ7KxrSYjWT6+JwMwmmNkLHX39ZnakmX23H469CBjjnJsK3AY83t2Ozrk5zrmZzrmZQ4Z0PQrfXxpCYdJSjIw0v3OkSPJJSTF++ZmpTCjN5Wv3vcP7uxr9DqlX7RHH9Y8sJSs9jRt7GCDe29DcDJrb2mmIffn0U1+udn8EbgDaAJxzS4ELDvbAzrk651xD7PHTQMDM+pZKPdTQEiYnmIaZ+R2KSFLKzkhjzsUzMYMr5i5IiAtlT/7y2gYWbNzND86aTMleM4h7UpoXBKCivsWr0PqsL4kgyzn31l7bDvqTMbNhFrvamtmsWCxVB/u+B6shFCY7Xd1CIn4aXZzF7f81nfW7Grn2gcVEIv53n3RlbUUDv/jHu5w6qZRPTRuxX387NDeaNCrqBkYi2GVmhwAOwMzOA3qtIWtm9wHzgcPMbIuZfdHMvmxmX47tch6w3MyWAL8BLnAJ0FlW3xLW+IBIAjj+0BK+c+Yk/rlyJ7e+8J7f4eyjPeL45kNLyExP5SfnTNnvXoShebFEkAADxn254l0FzAEmmtlWYAPwud7+yDl3YS+v/xb4bV+CjKeGUFh3DIkkiMuOL2fl9jpufeE9Jg3P5fQpw/0OaY/fvbSWxZtruPWCaQzNDe733w/t6BpKgBZBr1c859x64FQzywZSnHOJU0TbA42tYYqyVYJaJBGYGT/61BTWVjRw7YNLKC/JZuKwPL/D4q0N1fzq+TWcNbWMs6aWHdB75GakEQykJHaLwMyu7WY7AM65WzyKyVcNoTCji7L8DkNEYoKBVP5w8Qw+edurXDF3AfOuOoFCH7+s7W5s5Zr732FUURY//vT+dwl1MDOKszOoavB/TYaexghyYz8zga8AI2I/Xwamex+aPzRGIJJ4SvOC3HHxDHbWtnDV3xbREvanDEUk4vjWw0uoamjl9v+aTm7w4ErRFOekU92UwInAOfdD59wPgZHAdOfcdc6564AZwOh4BRhvGiMQSUzTRxfyv+ccwevrqrj2gSW0+3An0a+eX8Pzqyr4zscnMWVE/kG/X2FWekKs0taXK14p0DnS1ti2QSfcHqG5rZ0cFZwTSUjnzhjJ7qZWfvT3VeRlpvGTTx8Rtzk/Tyzeym0vruWCo0fx+S4WmzkQxdnprK3wf8nKviSCucBbZvZY7PmngLs9i8hHHZVHVXBOJHFdfuI4apra+O1La8kLBrj+jImeJ4OFG6v59sNLmVVexM1nH/i4wN6KstPZnQBdQ325a+jHZvYMcGJs02XOuXe8DcsfHWsR5KprSCShXffRCdQ2t/GHf68nHKv971UyWL61lkv/8jZlBZn87nPTSe/H8jOF2ek0tbYTamsnGOi9UJ1Xer3imdloYBfwWOdtzrnBtbYc7JnKrhaBSGIzM3541uGkphh/fnUDTa3t/OhTU0hN6d9ksGZnPZfc+Ra5GWn89fLZ+1VCoi+KY3c/VTe2UlaQ2a/vvT/6csX7O7FZxUAmMBZ4F+i2WuhA1VGCOlstApGEl5Ji3PTJyWRnpHL7S+uorG/h1xdM67ebPRZurOYLdy0gIy2Fv14+mxEeXKgLEyQR9NrGcc4d4Zw7MvYzHphFtHTEoLOnRaBEIDIgmBnf+thEfnjW4bz0bgXn/O61fqlY+sTirVz0pzcpyk7nka8cx7ghOf0Q7b46WgRVPt85tN+dXc65RcBsD2LxXYPWIhAZkC45rpy5X5hFRX0LZ9z6Cve8sfGA6vw3toT57uPLuOb+xUwpy+ehLx/LKA8nmHa0CHb7nAj6MkbQeYZxCtHJZNs8i8hHWp1MZOA6/tASnrnmRL798FK+9/hyHl20hetPn8jsccW9/m0k4nhm+Q5+/PeVbKsNccWJY/n26RMJpHq7LkmitAj6csXL7fQ4THTM4BFvwvGXBotFBrbh+ZnM/cIsHlq4hVv+uYbPznmDqaMK+MyMkZw0fgijijL33F0UiTg2VDXy/MqdPLhgM+sqGzmsNJeHLzyKmeVFcYk3LxggNcUSv0UArHTOPdR5g5l9Bniom/0HrPqOwWKtRyAyYJkZ588cxVlTy7jvrU387c1NfPfx5UC023dIbgaRiKOyvoXG1ujcoamjCvjNhUdx5pRhpHncCugsJcUozAoMiBbBDex70e9q24DX0BImOz21329BE5H4CwZSuez4sVx6XDlrKxp4Y30Vaysa2NXQSmqKUZSdzsRhuRx/aImn4wC9KcxKT9wWgZmdAZwJjDCz33R6KY9+WKEsETWEwuoWEhlkzIzxpbmML83tfWcfFGQFqG1u8zWGnq5624AFwFnAwk7b64FveBmUXxpaVHBOROIrPzPA1hp/1yTo9qrnnFsCLDGze51zg7IFsLf6ljA5B1lWVkRkf+RnprNqu7/rffXUNfSgc+584B0z2+eGXOfckZ5G5oPGljA5Gf7V+xCR5JOfGaDG58JzPfWDXBP7/Yl4BJIIGkJhSnK0OpmIxE9BVoDG1nba2iOez1voTk9dQ9tjvzfGLxx/RccI1DUkIvGTnxm95tQ1t1Hcz0Xt+qqnrqF6Pig2B2Cx5wY455z/K0j3s/pQm8pLiEhcdSSCmkRMBM65xLzXyiPOOd01JCJxl58VTQR+3kLap6uemU0HTiDaInh1MC5M09zWTsSpvISIxFdHi6C2yb9E0OvIhJl9n+jSlMVACXCXmX3X68DiTQXnRMQPBZkDo0VwETDVORcCMLOfAouBH3kYV9zVqwS1iPggPwESQV/uVdoGBDs9zwC2ehOOf9QiEBE/5HUMFvvYNdSXq14tsMLMniM6RnAa8FZH/SHn3Nc9jC9uGrU6mYj4IJCaQk5GWsJ3DT1Gp4XrgX95E4q/6rUWgYj4JD8zQE2zf7OLe73qOefujkcgflPXkIj4JT8zQF0ijxGY2SfM7B0zqzazOjOrN7O6eAQXT1q4XkT8kp/pbynqvgwW/xq4BCh2zuU553IH46xiLVMpIn6JFp5L7ESwGVjunNunAulgUh8Kk56aQkaaqo+KSHz5vThNX77+fht42sxeBlo6NjrnbvEsKh80tLSpNSAivogOFid2Ivgx0EB0LkG6t+H4pyGkOkMi4o/8rACt4QihtnaCgfj3SvTlylfmnJvieSQ+U8E5EfFLfqdJZcPy458I+jJG8LSZfdTzSHxWr4XrRcQne9YkCPnTPdSXRPAV4Fkza96f20fN7E4zqzCz5d28bmb2GzNba2ZLYxVOfdPYGiZXLQIR8YHf9YZ6TQSx20VTnHOZ+3n76F3A6T28fgYwPvZzJfD7vgTslQa1CETEJ3lBf0tR93U9gkKiF+w9xeecc//u6W+cc/82s/IedjkbmBu7LfUNMysws+EdS2TGW0NLmGy1CETEB353DfV65TOzy4kuZD+SaPnpY4D5wH8e5LFHEJ2j0GFLbNs+icDMriTaamD06NEHediu1YfUNSQi/kj4riGiSeBoYKNz7iPAUUCNl0HtzTk3xzk30zk3c8iQIf3+/q3hCC3hiO4aEhFfdKyDksiJINRpUZoM59xq4LB+OPZWYFSn5yPxaZ2DRpWXEBEfpcVKUdc1h305fl8SwRYzKwAeB54zsyeAjf1w7HnA52N3Dx0D1Po5PgAqOCci/vGz8FxfylB/OvbwB2b2EpAPPNvb35nZfcDJQImZbQFuAgKx97wDeBo4E1gLNAGXHUD8/aI+pGUqRcRfuUH/FqfZryufc+7l/dj3wl5ed8BV+3N8r3zQIgj4HImIJKv8zEBCTygb9DRGICJ+83NxGiUCOi1TqTECEfFJno9jBH1KBGY2xsxOjT3ONLNcb8OKrwaNEYiIzxK6RWBmVwAPA3+IbRpJ9A6iQaOhJfqPrxaBiPglPzNAY2s7be2RuB+7Ly2Cq4DjgToA59x7wFAvg4q3hlAYM8hK1+pkIuKPvFiPhB+tgr4kghbnXGvHEzNLAwbVspX1LWFy0tMwM79DEZEklZ/VUW8o/pPK+pIIXjazG4FMMzsNeAh40tuw4kuVR0XEb37WG+pLIrgeqASWAV8iOhHsu14GFW9anUxE/LanFLUPiaAvM4sjwB+BP5pZETAyNhls0GhoUYtARPy1pxR1IrYIzOxfZpYXSwILiSaEX3kfWvyoRSAifkv0rqF851wdcA7RhWRmA6d4G1Z8NYTCmkMgIr7KS/BEkGZmw4Hzgac8jscXahGIiN+CgVTS01J8qTfUl0RwM/APYK1z7m0zGwe8521Y8dUQCqvgnIj4zq/ZxX0ZLH6I6C2jHc/XA+d6GVQ8RSKOhlYNFouI//J8KkXd7dXPzG6jh4ljzrmvexJRnDW1teMc5GRoVrGI+CvaIoj/hLKevgYviFsUPuooOKeuIRHxW35mgF0Nrb3v2M+6TQTOubvjGYhf9hScU9eQiPgsLzPAusrGuB+316tfbHnKfbqInHP/6UlEcbZnmUrdNSQiPvNrlbK+XP2+2elxkOhAcfw7sTzS2NIOqEUgIv7ruGsoEnGkpMSvCGZf7hpauNem18zsLY/iiTutRSAiiSIvGCDioKE1vKf2UDz0pWuoqNPTFGAGkO9ZRHFWH9IylSKSGPaUmWhqS6xEQLS+kAOMaJfQBuCLXgYVTw0tWqZSRBJDR5mJeI8T9KVraGw8AvFLx+2j2WoRiIjP8jKj16F4TyrrS9dQEPgqcALRlsErwB3OuZDHscVFfUuYYCCFQGpfqm2IiHjHr1LUffkaPBeoB26LPf8v4B7gM14FFU/1qjMkIgnig0QQ3xsz+5IIpjjnJnd6/pKZrfQqoHhraAnvWTRaRMRPfpWi7kt/yCIzO6bjiZnNZhCVn6gPtWkOgYgkhJz0NFIsAccIiN4u+rqZbYo9Hw28a2bLAOecO9Kz6OJAi9KISKJISTHyfJhd3Jcr4OmeR+Gj+lCY4pwsv8MQEQGik8oSrkXgnNsYj0D8El2dTIPFIpIY8jPjnwiS/p7J+lCbuoZEJGH4sUpZUicC5xwNLRojEJHEkZcZ/1XKkjoRNLW2E3EqLyEiiSPaNRTfeQRJnQjqtTqZiCQYP+4aSupEoNXJRCTR5AUDtIYjhNra43bMpE4EdSFVHhWRxJLvw+zipE4EDVqmUkQSjB+F5zxNBGZ2upm9a2Zrzez6Ll6/1MwqzWxx7OdyL+PZ2wdrEWiMQEQSgx/1hjz7KmxmqcDtwGnAFuBtM5vnnNu7YN0DzrmrvYqjJ/UhjRGISGIZbF1Ds4C1zrn1zrlW4H7gbA+Pt9/qNUYgIgkm34dVyrxMBCOAzZ2eb4lt29u5ZrbUzB42s1FdvZGZXWlmC8xsQWVlZb8F2JEIstOVCEQkMXSUxa9tGhyJoC+eBMpjFUyfA+7uaifn3Bzn3Ezn3MwhQ4b028EbWsJkp6eSmmL99p4iIgfjgzGC+E0q8zIRbAU6f8MfGdu2h3OuyjnXEnv6J6Ilr+MmWmdIA8UikjgCqSlkp6cOmq6ht4HxZjbWzNKBC4B5nXcws+Gdnp4FrPIwnn00tIQ1UCwiCScvzhVIPbsKOufCZnY18A8gFbjTObfCzG4GFjjn5gFfN7OzgDBQDVzqVTxdqdeiNCKSgOJditrTq6Bz7mng6b22fb/T4xuAG7yMoSdKBCKSiPLiXIra78FiX2ktAhFJRPFepSypE0Ftc3jPPbsiIoki3ovTJG0icM5R19y251YtEZFEkZ8Z2FMUMx6SNhGE2iK0tkfUIhCRhJOfGaChJUy4PRKX4yVtIujof1MiEJFEk5cZHbuMV6tAiUCJQEQSTLxLUSsRKBGISIKJdwVSJQIlAhFJMPFek0CJQIlARBJMvEtRJ20iqFMiEJEEpa6hOOn4B1b1URFJNHlBJYK4qG2OlpfQWgQikmiCgRTSU1Ooi9OaBEmbCOqa29QtJCIJycziWoo6aRNBrRKBiCSwgqwAuxtb43KspE4EeRofEJEEVZKTTlVjS+879oOkTQS7m1opzFYiEJHEVJKTwa4GtQg8Vd3YSlF2ut9hiIh0qSQng131ahF4pj3iqGluoyg7w+9QRES6VJKTTn1LmFBbu+fHSspEsLupFeegWC0CEUlQJTnRL6pVcRgwTspEUB37h1XXkIgkquKORNDgffdQUiaCqtgAjFoEIpKoSnKi16ddSgTe2NMiyFEiEJHE1NE1tKteXUOeqI7dm6uuIRFJVHsSQRzmEiRlIugYfCnMUiIQkcSUmZ5KdnqqWgReqW5sJT8zQCA1KU9fRAaI4pwMjRF4paqhVQPFIpLwhuRmUBmHSWVJmQh21IUYlh/0OwwRkR4Nzw+yvbbZ8+MkZyKoDTEsT4lARBLbiIJMttWGiEScp8dJukQQiTh2qkUgIgNAWUEmreGI57OLky4R7GpsIRxxSgQikvDKCjIB2FbjbfdQ0iWCHbUhAHUNiUjCKyuIXqe8HidI2kQwPD/T50hERHo2ItYi2FoT8vQ4SZcItuyOZtbhBWoRiEhiy88MkJWeypbdTZ4eJ+kSwYZdjeQG0zSPQEQSnpkxtiSbdZWNnh4nKRPBuJJszMzvUEREenVYaS5rdtR7eoykSwTrKxsYNyTH7zBERPpkfGkuO+pC1Da3eXaMpEoEtc1tbKsNcciQbL9DERHpkwml0S+u7+30rlXgaSIws9PN7F0zW2tm13fxeoaZPRB7/U0zK/cynnc27QZg+uhCLw8jItJvjhiRD8DCjbs9O4ZnicDMUoHbgTOAycCFZjZ5r92+COx2zh0K/Ar4mVfxALz63i7SUoypowq8PIyISL8ZmhdkQmkO/3q30rNjeNkimAWsdc6td861AvcDZ++1z9nA3bHHDwOnmEejuC+u3slf39zIqZNKyc5I8+IQIiKe+OSRZcxfX8WfXlnvyft7eUUcAWzu9HwLMLu7fZxzYTOrBYqBXZ13MrMrgSsBRo8efUDBjC3J4dhxxXz/k3s3SkREEtsVJ41j/a5GDvHoRpcB8dXYOTcHmAMwc+bMAyrDN7Ykm79cNqtf4xIRiYdgIJVffXaaZ+/vZdfQVmBUp+cjY9u63MfM0oB8oMrDmEREZC9eJoK3gfFmNtbM0oELgHl77TMPuCT2+DzgReect4W3RUTkQzzrGor1+V8N/ANIBe50zq0ws5uBBc65ecCfgXvMbC1QTTRZiIhIHHk6RuCcexp4eq9t3+/0OAR8xssYRESkZ0k1s1hERPalRCAikuSUCEREkpwSgYhIkrOBdremmVUCGw/wz0vYa9ZyEtA5Jwedc3I4mHMe45wb0tULAy4RHAwzW+Ccm+l3HPGkc04OOufk4NU5q2tIRCTJKRGIiCS5ZEsEc/wOwAc65+Sgc04OnpxzUo0RiIjIvpKtRSAiIntRIhARSXJJkwjM7HQze9fM1prZ9X7HEw9m9r6ZLTOzxWa2wO94vGBmd5pZhZkt77StyMyeM7P3Yr8L/Yyxv3Vzzj8ws62xz3qxmZ3pZ4z9ycxGmdlLZrbSzFaY2TWx7YP2c+7hnD35nJNijMDMUoE1wGlEl8x8G7jQObfS18A8ZmbvAzOdc4N20o2ZnQQ0AHOdc1Ni234OVDvnfhpL+oXOuf/xM87+1M05/wBocM790s/YvGBmw4HhzrlFZpYLLAQ+BVzKIP2cezjn8/Hgc06WFsEsYK1zbr1zrhW4Hzjb55ikHzjn/k10LYvOzgbujj2+m+j/QINGN+c8aDnntjvnFsUe1wOriK53Pmg/5x7O2RPJkghGAJs7Pd+Ch/+oCcQB/zSzhWZ2pd/BxFGpc2577PEOoNTPYOLoajNbGus6GjTdJJ2ZWTlwFPAmSfI573XO4MHnnCyJIFmd4JybDpwBXBXrUkgqsaVPB3//J/weOASYBmwH/s/XaDxgZjnAI8B/O+fqOr82WD/nLs7Zk885WRLBVmBUp+cjY9sGNefc1tjvCuAxol1kyWBnrI+1o6+1wud4POec2+mca3fORYA/Msg+azMLEL0g3uucezS2eVB/zl2ds1efc7IkgreB8WY21szSia6NPM/nmDxlZtmxQSbMLBv4KLC8578aNOYBl8QeXwI84WMscdFxQYz5NIPoszYzI7q++Srn3C2dXhq0n3N35+zV55wUdw0BxG6z+jWQCtzpnPuxvxF5y8zGEW0FQHRt6r8NxnM2s/uAk4mW590J3AQ8DjwIjCZasvx859ygGVzt5pxPJtpd4ID3gS916j8f0MzsBOAVYBkQiW2+kWif+aD8nHs45wvx4HNOmkQgIiJdS5auIRER6YYSgYhIklMiEBFJckoEIiJJTolARCTJKRGI9MLMCszsq7HHZWb2sN8xifQn3T4q0otYrZenOip9igw2aX4HIDIA/BQ4xMwWA+8Bk5xzU8zsUqIVL7OB8cAvgXTgYqAFONM5V21mhwC3A0OAJuAK59zqeJ+ESHfUNSTSu+uBdc65acC39nptCnAOcDTwY6DJOXcUMB/4fGyfOcDXnHMzgG8Cv4tH0CJ9pRaByMF5KVYvvt7MaoEnY9uXAUfGqkceBzwULR8DQEb8wxTpnhKByMFp6fQ40ul5hOj/XylATaw1IZKQ1DUk0rt6IPdA/jBWQ36DmX0GolUlzWxqfwYncrCUCER64ZyrAl6LLRb/iwN4i4uAL5rZEmAFWiZVEoxuHxURSXJqEYiIJDklAhGRJKdEICKS5JQIRESSnBKBiEiSUyIQEUlySgQiIknu/wPd73priZWVrAAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plot_pulse(opt_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": { "attributes": { "classes": [], "id": "", "n": "22" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEGCAYAAAB/+QKOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAA3W0lEQVR4nO3dd3hUZfbA8e9JJyQQUmiBkNBBOqHaEEWKioKogCCIinVVFF30Z0FdV1awrmVFZAWUJrKIqIAFVBCBUELvNaEnARLSM+/vjztggAABMrmZmfN5nvvM3DZzLhPmzH2rGGNQSimlzuRjdwBKKaXKJk0QSimliqQJQimlVJE0QSillCqSJgillFJF8rM7gJISGRlpYmNj7Q5DKaXcyooVK44YY6KK2ucxCSI2NpaEhAS7w1BKKbciIrvPtU+LmJRSShVJE4RSSqkiaYJQSilVJE0QSimliqQJQimlVJE0QSillCqSJgillFJFcmk/CBHpBrwH+ALjjDGjzthfCxgPRAGpwABjTJJzXwGw1nnoHmNMT5cEmZcFv40B/yDwD7aWkCpQoTpUiIbykSDikrdWSqmyzGUJQkR8gQ+BLkASsFxEZhtjNhQ6bAww0RgzQUQ6A28AA537sowxLVwV3yk56bDobTCOovcHVYTKV0DlRlCtGdS6EiLqatJQSl2QMYbjWfkcOZFD6olcUjJySD2RR2ZuPlm5BWTlFZCZW0BOfgHGFP5aEUTAVwRfH8FHBF8f8PGRM7ZZJ0SFBtKvbUyJx+/KO4i2wDZjzA4AEZkK3AoUThCNgaeczxcAs1wYT9FCKsNLqVCQB3mZkHsCMg7A8X1wLAkOb4ZDG2HtV5DwmXVO+SgrUTToAfW7QrmwUg9bKWW/rNwCko9mkXw0i6S0TJLTTj7PYt/RLI5k5JBXcO5J2fx8hHIBvgT5+3IyNxjAmsfNUOCwFofBem4MDudj4bneWsaEuV2CiAb2FlpPAtqdcUwi0BurGKoXECoiEcaYFCBIRBKAfGCUMWaWyyIVAb8AaykXBhWjIbr16ccYAynbYPcfsHsx7PgVNswCHz+Iuwaa9YXGPcG/nMvCVEqVLofDcOB4NntSM9mTmsne1Ex2p1jPk9IyOZKRe9rxfj5CtbAgosPK0bFOJJUrBBJRPoDIkEDCywcQERJAePkAggP8KOfvS4DfpVcDG2MlDwA/X9dUJ9s9FtNw4AMRGQz8BiQDBc59tYwxySJSG/hFRNYaY7YXPllEhgJDAWJiSj57nkYEIutZS+tB4HBA8grYOBs2fAP/Gwo/PAPN7oI290NUA9fGo5S6bBk5+Rw4ls3B49kcOJbNAefj3jRnEkjNIrfgr+JnXx+helgQMeHB3NCoCjUqlaNGpWCiK5UjOqwcVSoEnSr2cTURwc/Xte8lrpqTWkQ6ACONMV2d688BGGPeOMfxIcAmY0yNIvZ9Dswxxsw41/vFx8cb2wbrczhg1++wcqKVMApyreKnK5+EmDNvmpTyDg6H4URuPidyCsjIyScjJ5/M3HwKHIZ8h6GgwPnoMOQ7HDiMweEAhzHOYharaMXhfOSMdeMsZnEU2v7Xtr/W8x2G49l5HMvK52hmLsey8jiamUfaiVzSc/LPirtCkB8xEcHEhAdTM9x6rBVenpjwYKqFBeHvol/rdhGRFcaY+KL2ufIOYjlQT0TisO4M+gL9zwgsEkg1xjiA57BaNCEilYBMY0yO85grgTddGOvl8fGB2tday4kjsGystWz+HmI6wvUvQa0OdkepVInJzXeQlJbJrpQT7DySyb6jWRxOz7GWDOvxWFae3WGequitUM6fsHL+VAz2J7x8ALUjyxMWHECVCkFUqxhElQpBVK0YRNUKQZQL8LU77DLDZXcQACLSA3gXq5nreGPM6yLyKpBgjJktIn2wWi4ZrCKmR51JoSPwCeDA6qvxrjHms/O9l613EEXJPQErJ8Gid6xK7/rdrERR5Qq7I1Pqohw8ns36fcdYn3yc9fuOs/HAcfamZuIo9NUR5O9D5dAgKocGEuVcwoIDCAn0JSTQn/KBvoQE+lEuwBd/Xx/8fAQ/Hx98faxiEh8R/Jwtc0SsL3YfEediFacU9XjymHOdoy7sfHcQLk0QpanMJYiTcjNh6X9g0buQcxxa3QM3jITgcLsjU+osxhh2p2SydGcKS3eksnRnKslHs07tj4ssT+NqFagTVZ7YyPLUiihPbEQw4eUD9AvZTWmCKAsyU60OeUv/A0EVrCTR8h6reEopG+XmO1i2M5WfNh7k500H2ZtqJYTw8gG0jQ2nTVw4TaMr0qhaKKFB/jZHq0qaJoiy5OAG+H641VQ2Oh5u/cDqhKdUKcovcLB4ewqzViXz04aDpOfkE+jnw5V1I7muQRTta0dQt3KI3hV4AbsqqVVRqjSGwd/Bmmkw73n45Bq47nno+Dj4aOWYcq2N+4/zVUISsxP3cSQjhwpBfnRvWpUujatyVd1IraBVp9EEYQcRaN4X6lwP3w2Dn0bCpu/gto+tfhZKlaC8Agdz1x1g0pLdLNuVSoCvD50bVua2ltFc1zCKQD9NCqpomiDsFBIFd06CdV/Dd0/Df66Gm8ZAi7t1rCd12Y5n5zHxj11MXLKbQ+k5xIQH83yPhtzRuiaVygfYHZ5yA5og7CYCTftA7FUw8wH45lHY+Tvc9BYEhtgdnXJDqSdyGb9oJxOW7CI9O59r6kfxr9tjubZ+FD6l1MtXeQZNEGVFaFUYOMtq6fTrKEhOgDs+h6pN7Y5MuYnj2Xl8vHA7ny/eRXZ+Ad2bVOWRTnVpEl3R7tCUm9IEUZb4+EKnv0OtjvD1/TCuC9z2ITS53e7IVBmWk1/AF3/u4YNftnI0K4+ezavzt851qVs51O7QlJvTBFEWxV0ND/0O0wbCjCFwcD1c94L2mVCnMcYwd90B/vnDRvamZnFV3UhGdG+odwyqxGiCKKtCKsOg2Vafid/fsvpP9B5rdbJTXm93ygle+mY9v245TMOqoUwc0pZr6kfZHZbyMJogyjK/QLjlfajSFOaOgPFd4e4Z1nwVyivl5Bfwn4U7+HDhNvx9hBdvbsygDrVcNh+A8m6aIMo6EWg31OofMW0gfNbFShJVGtsdmSpla5KO8vT0RLYeyuCmZtV48abGVK0YZHdYyoPpzw53Uec6GPKDNXf2+G5WU1jlFXLzHbz94xZ6ffQH6dn5/PfeNnzYv5UmB+VymiDcSdWmcN+PUKEafNEb1s20OyLlYpsPpHPbh4t5/+et3NqiOvOGXcN1DSrbHZbyElrE5G7CasKQuTClv9XCKfcEtBpod1SqhBljmLZ8Ly/PXk9okB+fDGxN1yuq2h2W8jKaINxRuUow4GuYNgBmPwZ5WVY9hfIIGTn5PD9zLbMT93FV3UjeuasFUaGBdoelvJAmCHcVEAz9plh3ET88A3mZcNWTdkelLtP6fcd4bPIqdqecYPiN9XmkU10dHkPZRhOEO/MLtIbj+N+D8NPLVpLo9JwO9OemZifu49kZiYSVC2Dq0A60jdNZB5W9NEG4O19/6P0p+JeDX/8F4gOdRtgdlboIBQ7DmPmb+XjhdtrGhvPRgFZEhmiRkrKfJghP4OMLt/wbjIGFb1hJ4+qn7Y5KFcPx7DyemLKKBZsP079dDCNvuYIAP21cqMoGTRCewscHev4bCvLg51fBNxA6PmZ3VOo89qZmMvi/y9idksk/bmvCgPa17A5JqdNogvAkPr7WrHQFuTD//8A3QFs3lVFrko4y5PPl5BUYvri/He1rR9gdklJn0QThaXz94PZx4Mi3Wjf5BUDrwXZHpQpZsOkQj05eSaXgAKYObUvdyjoxlCqbtLDTE/n6Q5/xULcLfPskrJ9ld0TKaeqyPdw/MYG4yPL879GOmhxUmaYJwlP5BcKdE6FmO2sq0+0L7I7Iqxlj+OCXrYyYuZYr60Yy7cEOVA7VsZRU2aYJwpMFBEP/qRBRD6beDckr7I7IKxljeHPeZsbM30KvltF8NiiekEAt3VVlnyYIT3dyWI7ykfBFHzi8xe6IvIrDYRg5ez0fL9xO/3YxvHVHc/x17gblJvQv1RtUqAYD/wc+fjDpNjiWZHdEXqHAYXj26zVMWLKbB66O4/XbmuiwGcqtaILwFhF1rDuJnHSY1Buy0uyOyKPlFzh4YuoqZqxI4onr6/F8j0aIDoGi3IwmCG9SrZk1wF/aTpg6APJz7I7IIxU4DE9NT2TOmv2M6N6QYV3qa3JQbkkThLeJvcrqTLd7Ecx6BBwOuyPyKA6H4ZkZidbAe90a8NC1dewOSalLpk0pvFHTPnB0D/z8CoTFwA0v2x2RR3A4DM/NXMvMlck81cUaqlspd6YJwltdNcxKEovetmapix9id0RuzRjDi9+sY1rCXv7WuS6PX1/P7pCUumyaILyVCPQYA8eT4bunoUINqH+j3VG5JWMMr3+3kS+X7uGha+vwVJf6doekVInQOghv5usHff4LVZvCV4Nh32q7I3JLHy3czrhFOxncMZa/d2ugFdLKY2iC8HaBIdB/OgSHw+Q7tY/ERZq6bA+j523m1hbVeenmxpoclEdxaYIQkW4isllEtonIWdOciUgtEflZRNaIyEIRqVFo3yAR2epcBrkyTq8XWhXu/grysmBKX8jJsDsitzB33X6e/99arq0fxeg+zbUTnPI4LksQIuILfAh0BxoD/USk8RmHjQEmGmOaAa8CbzjPDQdeBtoBbYGXRaSSq2JVQOVGVnHTwfUwc6g2f72AP7Yf4fEpq2lRM4yPB7TSWeCUR3LlX3VbYJsxZocxJheYCtx6xjGNgV+czxcU2t8V+NEYk2qMSQN+BLq5MFYFUO8G6PoGbP4Ofh5pdzRl1rrkYwyduILYyGDGD25DcIC29VCeyZUJIhrYW2g9ybmtsESgt/N5LyBURCKKeS4iMlREEkQk4fDhwyUWuFdr9yDE3weL34NVX9odTZmz72gWQz5fToUgPyYOaUdYcIDdISnlMnbfFw8HrhWRVcC1QDJQUNyTjTFjjTHxxpj4qKgoV8XoXUSg+7+gdif49gnYtdjuiMqM9Ow8hny+nKzcAv57b1uqVtT5HJRnc2WCSAZqFlqv4dx2ijFmnzGmtzGmJfB/zm1Hi3OuciFff7hjAlSKhWkDIHWH3RHZLq/AwaOTV7HtUAYfDWhFg6qhdoeklMu5MkEsB+qJSJyIBAB9gdmFDxCRSBE5GcNzwHjn83nAjSJSyVk5faNzmyot5cKg/zTAwOS7IOuozQHZxxjDS9+s57cth3m9VxOurqd3q8o7uCxBGGPygcewvtg3AtONMetF5FUR6ek8rBOwWUS2AFWA153npgKvYSWZ5cCrzm2qNEXUgTsnWXcQM+6Fgny7I7LFJ7/tYMqyPTzSqQ53tYmxOxylSo0YY+yOoUTEx8ebhIQEu8PwTCsnwuy/QZsH4KYxdkdTqn5Yu5+Hv1zJzc2q8X7fltrXQXkcEVlhjIkvap+2z1MX1uoeOLIF/vg3VG4Ibe63O6JSsX7fMZ6ankirmDDG3KEd4ZT3sbsVk3IXN7wC9bvB98/CjoV2R+NyRzJyGDpxBWHB/vxnYGuC/H3tDkmpUqcJQhWPjy/0/hQi68P0QZCy3e6IXCY338EjX6zkSEYOYwfGUzlUm7Mq76QJQhVfUAXoP9VKFpPv9Nh5rUd+u55lu1J5s08zmtaoaHc4StlGE4S6OJVi4a4vIG23NUS4h7VsmvTnbiY753W4tcVZnfeV8iqaINTFq9URbn7HqouY95zd0ZSYJdtTeGX2ejo3rMwzXRvYHY5SttNWTOrStBoIhzfBkg8gqiG0uc/uiC5L8tEsHp28kloRwbzbtwW+2mJJKb2DUJehy6tQ70b4/hnY8avd0VyynPwCHvlyJbn5DsbeE0+FIH+7Q1KqTNAEoS6djy/c/hlE1oPp97hty6bX5mwgce9RRvdpRp2oELvDUarM0AShLk9QBeg3FcTHLcdsmrkyiS/+3MPQa2rTvWk1u8NRqkzRBKEuX3ics2XTLrcas2nTgeM8/7+1tI0L51mtlFbqLJogVMmIvRJufhu2/wLz/8/uaC7oeHYeD01aQYUgfz7o3xI/X/2voNSZtBWTKjmt7oFDm+DPD62WTfH32h1RkYwxDJ+eyN60LKYOba89pZU6B/3ZpErWja9B3S7w/XDY+Zvd0RTpk992MH/DQZ7r3pA2seF2h6NUmaUJQpUsH1/o8xmE1ymTLZuW70pl9LzN3NS0GvddFWd3OEqVaZogVMkLqmiN2QQwpW+ZadmUdiKXx6esokalcoy6vSki2hlOqfPRBKFcI7x2odnohtjesskYw/CvEknJyOXD/q0I1c5wSl2QJgjlOnFXw01vwfafYf4LtoYyfvEuft50iOd6NKRJtI7QqlRxaCsm5VqtB8PhzfDnR9ZsdK0Hl3oIiXuPMuqHjXRpXIXBHWNL/f2Vcld6B6Fcr8trUPcG+O5p2Pl7qb718ew8HpuyksqhQYzu00zrHZS6CJoglOv5+kGf8c6WTQOteolSYIzhuZlr2Xc0m/f7tSAsOKBU3lcpT6EJQpWOwi2bJveF7GMuf8vJy/bw3Zr9DL+xAa1raX8HpS5WsROEiPiKSHURiTm5uDIw5YFOtWzaDjPuA0eBy95q4/7jvPLtBq6pH8WD19R22fso5cmKlSBE5G/AQeBH4DvnMseFcSlPFXc19BgD236E+S+65C2ycgv425RVVCznz9t3NsdHJ/9R6pIUtxXTE0ADY0yKK4NRXiL+Xms2uj8/tFo2tbqnRF/+jR82su1QBpPua0tkSGCJvrZS3qS4RUx7AdcXGivvcePrUOd6mPMU7FpUYi/788aDTFyym/uviuPqelEl9rpKeaPiJogdwEIReU5Enjq5uDIw5eFOtmyqFAvTBkLqzst+ycPpOTw7Yw0Nq4byTDed30Gpy1XcBLEHq/4hAAgttCh16cqFQf9pYBzWmE3Zxy/5pYwxPDMjkYycfN7v15JAP9+Si1MpL1WsOghjzCsAIhLiXM9wZVDKi0TUgbsmwaRe8PV91vSlPhf/5T5xyW4Wbj7Mq7deQf0q+ttFqZJQ3FZMTURkFbAeWC8iK0TkCteGprxG3DXQYzRsnQ8/vnTRp285mM7r32/kugZRDGxfywUBKuWdituKaSzwlDFmAYCIdAI+BTq6JizldeKHWLPRLfkAohoUu2VTdl4Bj09ZRYUgP97s01yH0lCqBBW3DqL8yeQAYIxZCJR3SUTKe3X9J9Tp7GzZtLhYp4yet5lNB9IZ3ac5UaHapFWpklTsVkwi8qKIxDqXF7BaNilVcnz9oM9/rZZN0wdC2q7zHv7blsN8tmgngzrU4rqGlUslRKW8SXETxBAgCpjpXKKc25QqWSdbNjkKnGM2Fd2yKfVELk9/lUi9yiE816NR6caolJcoVoIwxqQZYx43xrRyLk8YY9JcHZzyUhF14M6JcGSLcza6vNN2G2N4fuZajmXm8W7fFgT5a5NWpVzhvAlCRN51Pn4rIrPPXEolQuWdal8LN79tjdk050kw5tSuWauTmbv+AMO61OeK6jo7nFKucqFWTJOcj2Mu5cVFpBvwHuALjDPGjDpjfwwwAQhzHjPCGPO9iMQCG4HNzkP/NMY8dCkxKDfWejAc3w+/joLQatD5BfYfy+Klb9bTulYlhuoorUq51HkThDFmhfNpC2PMe4X3icgTwK/nOldEfIEPgS5AErBcRGYbYzYUOuwFYLox5mMRaQx8D8Q69203xrS4iGtRnqjTCEjfB7+NxoRW49k1TcgvMLx1R3N8dZRWpVyquJXUg4rYNvgC57QFthljdhhjcoGpwK1nHGOACs7nFYF9xYxHeQsRuOkdqN8N891wym3/ged7NCQ2UltZK+Vq572DEJF+QH8g7ow6h1Ag9QKvHY01CuxJSUC7M44ZCcx3zjdRHrih0L44Z+/t48ALxpizJjMWkaHAUICYGJ2/yGP5+rGn8wekbe7GB4Ef4l+9K3/daCqlXOVCdRB/APuBSOCtQtvTgTUl8P79gM+NMW+JSAdgkog0cb5njDEmRURaA7NE5ApjzGltHo0xY7F6eRMfH2/OfHHlGQochqf+t5WDMoIFYW8gU/rCkHnWXBJKKZc5bxGTMWa3MWahMaaDMebXQstKY0z+BV47GahZaL2Gc1th9wHTne+1BAgCIo0xOScnJ3LWg2wH6hf/spQnGff7DhJ2pzGsZwf87pkJvgHwxe1w7Mw/J6VUSSruYH3tRWS5iGSISK6IFIjIhcZmXg7UE5E4EQkA+gJnNo3dA1zvfI9GWAnisIhEOSu5EZHaQD2057ZX2nwgnbfmb6HrFVXo1TLa6mU9YAZkH4NJt8GJI3aHqJTHKm4l9QdYxUFbgXLA/VgtlM7JeYfxGDAPq8nqdGPMehF5VUR6Og97GnhARBKBKcBgY4wBrgHWiMhqYAbwkDHmQnUeysPk5jt4avpqQoP8+Gevpn8NxFetudXb+ugea5jwbJ3sUHkpRwFs+wnWfe2SlxdjLlx0LyIJxph4EVljjGnm3LbKGNPSJVFdgvj4eJOQkGB3GKoEvT1/M+//so1PBram6xVVzz5gy3yY2g9qtIEBMyEguPSDVMoOhzbC6smwZjpkHIAqTeDh4g1weSYRWWGMiS9qX3GH+850FhOtFpE3sSqRi3v3odRFS9x7lA8Xbqd3q+iikwNA/Ruh96fWREPTBkC/KeCnI7oqD3XiCKydAYmTYX8i+PhBvRuheV+o380lb1ncBDEQq6fzY8AwrMrn210SkfJ62XkFPDV9NZVDA3n5lgvMS9WkN+SegNmPwdf3W6PB+hb3z1qpMi4/F7bOg9VTrEdHPlRtBt1GQZM+EBLl0rcv7pSju51Ps4BXXBeOUvDm3M1sP3yCL+5rR8Vy/hc+odVAyEmHec/Bt49Dzw/AR29wlZsyBvattJLCuhmQlQYhVaD9w9C8H1Qpvck8L9RRbi1Wb+cinayPUKqkLNmewvjFO7mnQy2uqhdZ/BM7PAI5x2HhG+Drb/W+1iSh3Enablj7lVWvcGQz+AVBw5uspFD7OlvujC/0jjeXShRKAenZeQz/KpHYiGBGdL+ETnDX/h0KcuF3Z59OTRKqrDuRAutnWnULe/+0ttVsD7e8B41vs+ZHsdGFBuvbfb79SpWkf8zZyP5jWXz1UEeCAy7h15IIdH7RukVf9La1TZOEKmsyDsHm72HjHNixwKpXiGoE179k1StUqmV3hKcU63+hiKTzV1FTAOAPnDDGVDj3WUoV3y+bDjItYS8Pd6pD61qVLv2FRKz/aOBMEgI3va1JQtkrZTtsmgObvoO9ywADYTHQ/hFodqfVTFXK3ujExa2kDj35XKzeSrcC7V0VlPIuaSdy+fvXa2lYNZQnb6h3+S94KkkYWPSOtU2ThCpNOemwa7F1h7D9F2t2RLBaIHV6zqpbqHJFmUwKhV30fbyzp/MsEXkZGFHyISlv8+I36ziamcuEe9sS6FdC04eKwPUvW88XvWPdxt/yHvjo9KTKBQryIHkF7FhoLUnLrb85vyCo1RHa3A8Nult3DW6kuEVMvQut+gDxQLZLIlJe5dvEfcxZs59nujagcfUSLrE8mSR8/OC30ZCbAb3Ggl9Ayb6P8j7GWHcFOxbC9gWwaxHkpgMC1VtAx8ehdieo2Q78g+yN9TIU9w7ilkLP84FdnD35j1IX5dDxbF78Zh0taobxoKumDxWBzi9AYAX48UXIyYA7J+qwHOripR+Enb9aCWHHQmumQ4BKcdC0j5UQ4q6B4HA7oyxRxa2DuNfVgSjvYozh71+vITuvgLfubI6fr4vrB658HIIqwLdPWkOF959mrSt1LrmZsOcPKyFsXwCH1lvby1WCuGuhznVWUqgUa2eULlXcIqbawHtYFdMGWAIMM8boENzqkkxbvpcFmw/z8i2NqRMVUjpv2nowBITA/x6ECbfAgK+h/EV0xlOezRg4sAa2/WxVLu/50+pX4xsAMe2t4so610HV5l7T4KG4RUyTsYb37uVc74s1PPeZU4gqdUF7UzN5bc4GOtSOYFCH2NJ986Z9IDAUpt8Dn3WBu2dARJ3SjUGVHXnZsOt3q1/C5rl/FRtVaQJth1oJIaaj1xZJFjdBBBtjJhVa/0JEnnFFQMqzORyG4V8lIiKMvqMZPj42NPOr3xXumQ1T+lpJot80qNmm9ONQ9shJh03fw6ZvYdsvkHcC/MtD3c5Q/wWoewOEVrE7yjKhuAniBxEZAUzFKmK6C/heRMIBdDIfVVzjF+9k6c5U3uzTjBqVbPxVFtMO7vsRvrwdJtwMt4+DRrdc+DzlnvKyYOt8a2KdLfMgPxtCq0Hzu6BBD4i92q1bG7lKcScM2nme3cYY46ImKMWnEwaVfdsOpdPj/UVcUy+ST++J/2uGODudOAKT77LasHd7wxoxU3kGhwN2L4JVX1o9mHPToXwUXNHLGtKiRhuvqUs4n8ueMMgYE1eyISlvk1fg4KnpiZQP8OWfvZuWjeQAViX1oG9h5gMwdwQc3gTdR2tfCXd2fB+s/hJWfQFpuyCwIlxxm1X/VOsqnS/kIhS3FZM/8DDWXNEAC4FPjDF5LopLeZiPFmxnTdIxPrq7FZVDy9itfECw1Tfil39Y4zcd2gR3TYKQynZHpoqrIB+2zIWVE2Hbj2AcVrHRdf9nFR36l7M7QrdU3FT6MdYAfR851wc6t93viqCUZ1mbdIx//7KV21pUp0fTanaHUzQfX7jhZajaBGY9CmM7wV1fQHQruyNT55OZaiWF5ePg2F6rXuGqYdByAITbXvLt9oqbINoYY5oXWv9FRBJdEZDyLCenD40MCeSVnk3sDufCmtwOEfVg6t3w3+5w8zvQor/dUakzHVgHyz6xJtfJz7buFrq9AfW7axFSCSruv2SBiNQxxmyHUx3nClwXlvIUb83fzNZDGUwY0paKwcWYPrQsqNYMhi6ArwbDrIetcXZ6jIaA8nZH5t0K8q3+CsvGWn0X/MpB875Wf4VSnIbTmxQ3QTwDLBCRkz2nYwEdfkOd19IdKYxbtJMB7WO4tr5rJ1cvceUj4Z5v4Nd/wa9vQlIC3PE5VGlsd2TeJycdVk6CpR/D0T1QMQa6vAotB3rUuEdlUXETxGLgE+B64CgwD2u4DaWKlJGTz9NfJRITHszzPRrZHc6l8fGF656HWldarZw+vQ66jbKG7CgrrbA82bFkWPofWDEBco5ZPZq7/tPqt6DDtpeK4iaIicBx4DXnen9gEnCHK4JS7u/17zaQfDSLrx7scGnTh5Ylta+FhxbBzKEw50mrmOOW96FCGa1wd3f7E+GPD6y5mo2BxrdCx8cgurXdkXmd4v7PbWKMKXxvvUBENrgiIOX+Fmw6xJRle3no2jrEx3pIEUBIZRgw0yr//mkkfNQeeoyx2tbr3cTlczisns5LPrDqFwJCoe2D0P4ht5tkx5MUN0GsFJH2xpg/AUSkHaDdltVZ0k7k8uzXa2hYNZRhXUpg+tCyxMfH+sKqe4NVeT3zftj4DXR/EypUtzs695SXBYlT4c+PrAl4KkRDl9eg9SAIqmh3dF6vuAmiNfCHiOxxrscAm0VkLdZQG81cEp1yOyenD/383jYlN31oWRNZF4bMhT/+DQvfgO1trbqKtkO1iWVxZRy2+i4sHweZR6Bac+g9zurx7Osmrd28QHH/mru5NArlEWYXmj70iuoe/uvPxxeuetIqH//+GZj3HKyeDDe/DTXb2h1d2ZW8EpZ9ag2aV5AD9btBh8cg9iotqiuDijsW025XB6Lc28Hj2bw4ax0tY1w4fWhZFB4Hd38FG7+1xnL6rIs1GNz1L2lP3pPyc2HDN1bHtqTl1tDarQZadQxR9e2OTp2H3g+ry3Zy+tCc/ALeuqMUpg8ta0SgcU+o0xn+eN8qeto4B+KHwLXPeu+sdceSYeUESPgvnDgE4XWsZsIt+mv9gpvQBKEu29Tle1m4+TCv9LyC2qU1fWhZFBhi1UXED4GFo6zy9VVfQJshVjFKaFW7I3S9/Fxr0LxVk2DbT1Yz1Xo3Omdn66zDa7uZYs0H4Q50Pgh77EnJpNt7v9EyJoxJQ9rZM0NcWXV4M/w22ipv9/G3BpC78nHPnOT+0CYrKSROtSqdQ6tbdwotB1jFcKrMuuz5IJQqSn6BgyenrcLXR3izT3NNDmeKamDNVNfpOVj8rjXqaMJ46xd1m/uh7vXu3SM4dafVmW3dTDi4Dnz8oEF3aHmP+1+bAjRBqMvw4YLtrNxzlPf7tSQ6TMfbP6eIOtDz33DtCFjxubVMvgPCakGLu63OdhF17I7ywoyBI1utIqQNs6xZ+ABqtIVu/4ImvXUODQ/j0iImEekGvAf4AuOMMaPO2B8DTADCnMeMMMZ879z3HHAf1qixjxtj5p3vvbSIqXSt2pNGn/8s4ZZm1Xi3b0u7w3Ev+bmwaY51N7FrEWCgektrGswG3ctWssjNhD1LrF7OW+ZBmnP24arNrMR2RS/t6ezmzlfE5LIEISK+wBagC5AELAf6GWM2FDpmLLDKGPOxiDQGvjfGxDqfTwHaAtWBn4D6xphzDjGuCaL0nMjJ56b3fyevwPDDk1dTIUg7Nl2yY8lWMc3ar6wxiAAqxUG9LlD7OqtPRWm2gso6CskJsGsx7P7Duktw5IFfEMRdC/VvhHpdIaxm6cWkXMquOoi2wDZjzA5nEFOBW4HCYzgZoILzeUVgn/P5rcBUY0wOsFNEtjlfT0eQLQNem7OB3amZTH2gvSaHy1UxGjr+zVrSdsHWH61l5SRr3CewEkbNtlCliVWvEVnfKp66nBZBORnW3UDqTqsy/UAi7F8DR51dnsTXuqvp8Ig1j3PsVdbUrMqruDJBRAN7C60nAe3OOGYkMF9E/gaUB24odO6fZ5wbfeYbiMhQYChATIze5paGuesOMHX5Xh7pVId2tSPsDsezVIqFtg9YS1427FsFSctg7zLYsRDWTPvrWN8ACKkKoVWs5rPlKoF/sDX3sl85a07mglxryc+BrFTITLGW9ANw4vDp7x1e20oIrQdZjzXaWs12lVezu5K6H/C5MeYtEekATBKRYs9LaYwZC4wFq4jJRTEqp0PHs3lu5hqaRFfgyRu0B6xL+QdBrQ7WclJmqjWg3eHNkLrd+qJPP2BVHGcdhfwsa/C7glzreB8/K5H4BlgT6wRHQIUaVgKoFGc1P60UZ9V5BIbacpmqbHNlgkgGChdU1nBuK+w+nOM8GWOWiEgQEFnMc1UpcjgMw2esISuvgHfvakmAn3Z4KnXB4RDT3lrOx1EA4qNjG6nL5sr/5cuBeiISJyIBQF9g9hnH7MGapQ4RaQQEAYedx/UVkUARiQPqActcGKu6gIlLdvHblsP8302NqVtZix7KNB9fTQ6qRLjsDsIYky8ij2FNT+oLjDfGrBeRV4EEY8xs4GngUxEZhlVhPdhYzarWi8h0rArtfODR87VgUq615WA6b/ywic4NKzOgndb1KOUtdKgNdV7ZeQXc+sFijmTkMPfJa4gKDbQ7JKVUCdKhNtQle23OBjYfTOfze9toclDKy2hNozqnH9bu58ule3jwmtp0aqBDKCjlbTRBqCIlpWXy96/X0LxGRZ6+sYHd4SilbKAJQp0lr8DB41NWYQz8u18rbdKqlJfSOgh1lnd/2nJqlNaYCB1eQSlvpT8N1WkWbzvCRwu3c1d8TXo2r253OEopG2mCUKccycjhyWmrqR1Znpd7NrY7HKWUzbSISQFQ4DA8NT2RY1l5TBzSluAA/dNQytvpHYQC4INftvHblsO8fEtjGlWrcOETlFIeTxOE4veth3n35y30bhlN/7Y6lIZSyqIJwsvtO5rFE1NXU69yCP/o1QTRQd6UUk6aILxYbr6DxyavJCevgI8HtNZ6B6XUafQbwYuN+mETK/cc5YP+LakTpUN4K6VOp3cQXuq7NfsZv3gngzvGcnMz7e+glDqbJggvtO1QOn//eg0taobxfI9GdoejlCqjNEF4mWNZeTwwcQVB/j58dLeOs6SUOjetg/AiBQ7D41NWkZSWyeQH2lM9rJzdISmlyjBNEF5k9LzN/LrlMK/3akKb2HC7w1FKlXFavuAlZifu4z+/bqd/uxjublfL7nCUUm5AE4QXWJd8jGdnJNImthIjb7nC7nCUUm5CE4SHO5SezYOTVlApOICP7m6tldJKqWLTOggPlpVbwP0TEkg9kcv0BzsQFRpod0hKKTeiCcJDORyGJ6etYm3yMcYOjKdpjYp2h6SUcjNa3uChRs3dxLz1B3nxpsZ0aVzF7nCUUm5IE4QH+uLP3Yz9bQeDOtTi3itj7Q5HKeWmNEF4mF82HeTl2evp3LAyL97cWIfvVkpdMk0QHmT5rlQe/mIljatV4P1+LfHz1Y9XKXXp9BvEQ2zcf5whny8nOqwcn9/bhpBAbX+glLo8+i3iAfakZHLP+GWUD/Bj4n1tiQjR5qzKe+Xl5ZGUlER2drbdoZQpQUFB1KhRA39//2KfownCzR1Kz2bg+KXkFTiY/GAHalQKtjskpWyVlJREaGgosbGxWgfnZIwhJSWFpKQk4uLiin2eFjG5scPpOdz96VIOHc9h/OA21KsSandIStkuOzubiIgITQ6FiAgREREXfVelCcJNHcnI4e5xf7I3LZPxg9vQKqaS3SEpVWZocjjbpfybaIJwQ6knchkwbil7UjMZP6gNHepE2B2SUsoDaYJwM6kncun/6Z/sPHKCzwa1oWPdSLtDUkpdps8//5zHHnusyH0dO3YEYNeuXUyePPmcrzFhwgTq1atHvXr1mDBhQonEpQnCjew/lsWdnyxh55ETjBsUz5WaHJTyeH/88Qdw/gSRmprKK6+8wtKlS1m2bBmvvPIKaWlpl/3e2orJTew8coIB45ZyLCuPCUPa0r62FispdSGvfLueDfuOl+hrNq5egZcvMK/K66+/zoQJE6hcuTI1a9akdevWDB8+nE6dOjFmzBji4+M5cuQI8fHx7Nq1C4C9e/fSqVMnkpOTGTBgAC+//DIAISEhZGRkMGLECDZu3EiLFi0YNGgQw4YNO/V+8+bNo0uXLoSHWzNFdunShblz59KvX7/LulaXJggR6Qa8B/gC44wxo87Y/w5wnXM1GKhsjAlz7isA1jr37THG9HRlrGXZ+n3HGDR+GcbA1KHtaRKtI7MqVVatWLGCqVOnsnr1avLz82nVqhWtW7e+4HnLli1j3bp1BAcH06ZNG2666Sbi4+NP7R81ahRjxoxhzpw5Z52bnJxMzZo1T63XqFGD5OTky74WlyUIEfEFPgS6AEnAchGZbYzZcPIYY8ywQsf/DWhZ6CWyjDEtXBWfu1i09QgPf7GC0CA/Jt3fjjpRIXaHpJTbuNAvfVf4/fff6dWrF8HBVp+knj2L99u2S5cuRERYJQO9e/dm0aJFpyUIO7iyDqItsM0Ys8MYkwtMBW49z/H9gCkujMftTFm2h0H/XUb1sHLMeLijJgel3Jyfnx8OhwPgrD4JZzZDvZhmqdHR0ezdu/fUelJSEtHR0ZcRqcWVCSIa2FtoPcm57SwiUguIA34ptDlIRBJE5E8Rue0c5w11HpNw+PDhEgrbfg6H4Y3vN/LczLVcXS+SGQ93oHpYObvDUkoVwzXXXMOsWbPIysoiPT2db7/99tS+2NhYVqxYAcCMGTNOO+/HH38kNTWVrKwsZs2axZVXXnna/tDQUNLT04t8z65duzJ//nzS0tJIS0tj/vz5dO3a9bKvpay0YuoLzDDGFBTaVssYEw/0B94VkTpnnmSMGWuMiTfGxEdFRZVWrC51LCuPoZNW8IlzPodx98QTGlT8sVOUUvZq1aoVd911F82bN6d79+60adPm1L7hw4fz8ccf07JlS44cOXLaeW3btuX222+nWbNm3H777WcVLzVr1gxfX1+aN2/OO++8c9q+8PBwXnzxRdq0aUObNm146aWXTlVYXw4xxlz2ixT5wiIdgJHGmK7O9ecAjDFvFHHsKuBRY8wf53itz4E5xpgZRe0HiI+PNwkJCSURum027DvOw1+uIDktixdvbsygjrF2h6SU29m4cSONGjWyO4xTRo4cSUhICMOHD7c7lCL/bURkhfPH+FlceQexHKgnInEiEoB1lzD7zINEpCFQCVhSaFslEQl0Po8ErgQ2nHmuJ/kqYS+9PlpMdl4B0x5sr8lBKWU7l7ViMsbki8hjwDysZq7jjTHrReRVIMEYczJZ9AWmmtNvZRoBn4iIAyuJjSrc+smTHM/OY+Ts9cxcmUyH2hH8u39LInW4bqU8xsiRI+0O4ZK5tB+EMeZ74Psztr10xvrIIs77A2jqytjKgqU7UnhqeiL7j2XxeOe6PH59PZ0FTilVZmhPahtk5Rbw7k9bGPv7DmLCg5nxcEcdjVUpVeZogihlv245zIuz1rEnNZN+bWvywk2NKa/TgyqlyiD9ZiolB49n84/vNvJt4j5qR5VnygPtdZhupVSZpgXeLpaRk8/b8zfTafRC5q07wLAb6vPDE1drclDKQ6WkpNCiRQtatGhB1apViY6OPrV+cviNwjZv3kynTp1o0aIFjRo1YujQoUW+riuG874QvYNwkey8AqYn7OX9n7dyJCOXm5tV49muDYmJ0DmjlfJkERERrF69Gji7D0RIyNnD5Tz++OMMGzaMW2+1RiJau3btWcecHM47ISEBEaF169b07NmTSpVcW3epCaKEpWfn8eXSPYz7fSdHMnJoGxfOuEGNaFEzzO7QlPI+P4yAA2d/4V6Wqk2h+6gLH1dM+/fvp0aNGqfWmzY9uwGnq4bzvhBNECVk84F0pizbw9crk0jPzufqepE83KkFHWrr5OlKqXMbNmwYnTt3pmPHjtx4443ce++9hIWFnXaMq4bzvhBNEJfhwLFs5q0/wDerk1m55ygBvj50bVKVB66Oo1mNMLvDU0qV4C99V7n33nvp2rUrc+fO5ZtvvuGTTz4hMTGRwED7O8xqgrgI2XkFrNyTxp87Uvl962FW7TkKQP0qIbxwUyN6t6pBePkAe4NUSrmd6tWrM2TIEIYMGUKTJk1Yt27daZMMRUdHs3DhwlPrSUlJdOrUyeVxaYIoQnp2HgeOZbPvWDZbD6az+UA6mw5Yj7kFDnwEmkRX5JmuDeh6RVXqVtZ5GpRSl2bu3Llcf/31+Pv7c+DAAVJSUs6ay6Fr1648//zzp+aZnj9/Pm+8cda4pyXO6xNE6olc7vpkCdn5BWTlOsjMzSczt+C0YyJDAmlULZR7r4qlXVw48bHhVNAhuJVSFykzM/O0CumnnnqKpKQknnjiCYKCggAYPXo0VatWPe28wsN5AyU2nPeFuGy479J2qcN9n8jJZ/hXiQT5+xLk70s5f18qVwikWsUgqlUsR+2o8jp4nlJupKwN912WXOxw315/B1E+0I+PB1x4QnGllPI22pNaKaVUkTRBKKU8jqcUnZekS/k30QShlPIoQUFBpKSkaJIoxBhDSkrKqYrw4vL6OgillGepUaMGSUlJHD582O5QypSgoKDTWlAVhyYIpZRH8ff3Jy4uzu4wPIIWMSmllCqSJgillFJF0gShlFKqSB7Tk1pEDgO7L+MlIoEjJRSOu/C2a/a26wW9Zm9xOddcyxgTVdQOj0kQl0tEEs7V3dxTeds1e9v1gl6zt3DVNWsRk1JKqSJpglBKKVUkTRB/GWt3ADbwtmv2tusFvWZv4ZJr1joIpZRSRdI7CKWUUkXSBKGUUqpIXp8gRKSbiGwWkW0iMsLueEqDiOwSkbUislpELn4aPjcgIuNF5JCIrCu0LVxEfhSRrc7HSnbGWNLOcc0jRSTZ+VmvFpEedsZY0kSkpogsEJENIrJeRJ5wbvfIz/o81+uSz9mr6yBExBfYAnQBkoDlQD9jzAZbA3MxEdkFxBtjPLYzkYhcA2QAE40xTZzb3gRSjTGjnD8GKhlj/m5nnCXpHNc8EsgwxoyxMzZXEZFqQDVjzEoRCQVWALcBg/HAz/o813snLvicvf0Ooi2wzRizwxiTC0wFbrU5JlUCjDG/AalnbL4VmOB8PgHrP5bHOMc1ezRjzH5jzErn83RgIxCNh37W57lel/D2BBEN7C20noQL/7HLEAPMF5EVIjLU7mBKURVjzH7n8wNAFTuDKUWPicgaZxGURxS1FEVEYoGWwFK84LM+43rBBZ+ztycIb3WVMaYV0B141Fk04VWMVbbqDeWrHwN1gBbAfuAtW6NxEREJAb4GnjTGHC+8zxM/6yKu1yWfs7cniGSgZqH1Gs5tHs0Yk+x8PAT8D6uozRscdJbhnizLPWRzPC5njDlojCkwxjiAT/HAz1pE/LG+LL80xsx0bvbYz7qo63XV5+ztCWI5UE9E4kQkAOgLzLY5JpcSkfLOyi1EpDxwI7Du/Gd5jNnAIOfzQcA3NsZSKk5+STr1wsM+axER4DNgozHm7UK7PPKzPtf1uupz9upWTADO5mDvAr7AeGPM6/ZG5FoiUhvrrgGsKWcne+I1i8gUoBPWMMgHgZeBWcB0IAZraPg7jTEeU6l7jmvuhFXsYIBdwIOFyubdnohcBfwOrAUczs3PY5XLe9xnfZ7r7YcLPmevTxBKKaWK5u1FTEoppc5BE4RSSqkiaYJQSilVJE0QSimliqQJQimlVJE0QSh1iUQkTEQecT6vLiIz7I5JqZKkzVyVukTOsXDmnBw5VSlP42d3AEq5sVFAHRFZDWwFGhljmojIYKzRQ8sD9YAxQAAwEMgBehhjUkWkDvAhEAVkAg8YYzaV9kUodS5axKTUpRsBbDfGtACeOWNfE6A30AZ4Hcg0xrQElgD3OI8ZC/zNGNMaGA58VBpBK1VcegehlGsscI7Xny4ix4BvndvXAs2co3F2BL6yhtcBILD0w1Tq3DRBKOUaOYWeOwqtO7D+3/kAR513H0qVSVrEpNSlSwdCL+VE5xj+O0XkDrBG6RSR5iUZnFKXSxOEUpfIGJMCLBaRdcDoS3iJu4H7RCQRWI9Od6vKGG3mqpRSqkh6B6GUUqpImiCUUkoVSROEUkqpImmCUEopVSRNEEoppYqkCUIppVSRNEEopZQq0v8DYdDc7++AGqoAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "optimized_dynamics = opt_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.6" }, "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 }