{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Optimization of a State-to-State Transfer in a Two-Level-System" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "attributes": { "classes": [], "id": "", "n": "1" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "qutip 4.3.1\n", "numpy 1.15.4\n", "scipy 1.1.0\n", "matplotlib 3.0.2\n", "matplotlib.pylab 1.15.4\n", "krotov 0.0.1\n", "CPython 3.6.7\n", "IPython 7.2.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", "%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{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", "The purpose of this example is not to solve an especially interesting physical\n", "problem but to give a rather simple example of how the package can be used in\n", "order to solve an optimization problem.\n", "\n", "## Define the Hamiltonian\n", "\n", "In the\n", "following the Hamiltonian, guess field and\n", "states are defined.\n", "\n", "The Hamiltonian\n", "$\\op{H}_{0} = - \\omega \\op{\\sigma}_{z}$\n", "represents a\n", "simple qubit with energy\n", "level splitting $\\omega$ in the basis\n", "$\\{\\ket{0},\\ket{1}\\}$. The control\n", "field\n", "$\\epsilon(t)$ is assumed to couple via\n", "the\n", "Hamiltonian $\\op{H}_{1}(t) =\n", "\\epsilon(t) \\op{\\sigma}_{x}$ to the qubit,\n", "i.e., the control\n", "field effectively\n", "drives\n", "transitions between both qubit\n", "states. For now, we initialize the control\n", "field as constant." ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "attributes": { "classes": [], "id": "", "n": "2" } }, "outputs": [], "source": [ "def ham_and_states(omega=1.0, ampl0=0.2):\n", " \"\"\"Two-level-system Hamiltonian\n", " \n", " Args:\n", " omega (float): energy separation of the qubit levels\n", " ampl0 (float): constant amplitude of the driving field\n", " \"\"\"\n", " H0 = - 0.5 * omega * qutip.operators.sigmaz()\n", " H1 = qutip.operators.sigmax()\n", " \n", " psi0 = qutip.Qobj(np.array([1,0]))\n", " psi1 = qutip.Qobj(np.array([0,1]))\n", "\n", " eps0 = lambda t, args: ampl0\n", " return ([H0, [H1, eps0]], psi0, psi1)\n", "\n", "H, psi0, psi1 = ham_and_states()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The projectors $\\op{P}_0 = \\ket{0}\\bra{0}$ and $\\op{P}_1 = \\ket{1}\\bra{1}$ are\n", "introduced since they allow for calculating the\n", "population in the respective\n", "states later on." ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "attributes": { "classes": [], "id": "", "n": "3" } }, "outputs": [], "source": [ "proj0 = psi0 * psi0.dag()\n", "proj1 = psi1 * psi1.dag()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Define the optimization target\n", "\n", "First we define the time grid of the\n", "dynamics, i.e., by taking the following\n", "values as an example, we define the\n", "initial state to be at time $t=0$ and\n", "consider a total propagation time of\n", "$T=5$. The entire time grid is divided into\n", "$n_{t}=500$ equidistant time steps." ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "attributes": { "classes": [], "id": "", "n": "4" } }, "outputs": [], "source": [ "tlist = np.linspace(0, 5, 500)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next, we define the optimization targets, which is technically a list of\n", "objectives, but here it has just one entry defining a simple state-to-state\n", "transfer\n", "from initial state $\\ket{\\Psi_{\\init}} = \\ket{0}$ to the target state\n", "$\\ket{\\Psi_{\\tgt}} = \\ket{1}$, which we want to reach at final time $T$. Note\n", "that we also have to pass the Hamiltonian $\\op{H}(t)$ that determines the\n", "dynamics of\n", "the system to the optimization objective." ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "attributes": { "classes": [], "id": "", "n": "5" } }, "outputs": [], "source": [ "objectives = [\n", " krotov.Objective(initial_state=psi0, target=psi1, H=H)\n", "]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In addition, we have to define and assign a shape function $S(t)$ for the update\n", "in each control iteration to each\n", "control field that will be updated. This shape\n", "usually takes care of\n", "experimental limits such as the necessity of finite ramps\n", "at the beginning and\n", "end of the control field or other conceivable limitations\n", "for field shapes: wherever $S(t)$ is zero, the optimization will not change the\n", "value of the control from the original guess." ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "attributes": { "classes": [], "id": "", "n": "6" } }, "outputs": [], "source": [ "def S(t):\n", " \"\"\"Shape function for the field update\"\"\"\n", " return krotov.shapes.flattop(t, t_start=0, t_stop=5, t_rise=0.3, t_fall=0.3, func='sinsq')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "At this point, we also change the initial control field $\\epsilon_{0}(t)$ from a\n", "constant to a shaped pulse that switches on smoothly from zero and again\n", "switches off at the final time $T$. We re-use the shape function $S(t)$ that we\n", "defined for the updates for this purpose (although generally, $S(t)$ for the\n", "updates has nothing to with the shape of the control field)." ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "attributes": { "classes": [], "id": "", "n": "7" } }, "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", "H[1][1] = shape_field(H[1][1])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Having defined the shape function $S(t)$ and having shaped the guess field, we\n", "now tell the optimization to also use $S(t)$ as the update-shape for\n", "$\\epsilon_0(t)$. In addition, we have to choose `lambda_a` for each control\n", "field. It controls the update magnitude of the respective field in each\n", "iteration." ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "attributes": { "classes": [], "id": "", "n": "8" } }, "outputs": [], "source": [ "pulse_options = {\n", " H[1][1]: krotov.PulseOptions(lambda_a=5, shape=S)\n", "}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "It is convenient to introduce the function `print_fidelity`, which can be passed\n", "to the optimization procedure and will be called after each iteration and thus\n", "provides additional feedback about the optimization progress." ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "attributes": { "classes": [], "id": "", "n": "9" } }, "outputs": [], "source": [ "def print_fidelity(**args):\n", " F_re = np.average(np.array(args['tau_vals']).real)\n", " print(\" F = %f\" % F_re)\n", " return F_re" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Simulate dynamics of the guess field\n", "\n", "Before heading towards the optimization\n", "procedure, we first simulate the\n", "dynamics under the guess field\n", "$\\epsilon_{0}(t)$." ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "attributes": { "classes": [], "id": "", "n": "10" } }, "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)$, which is, as chosen\n", "above, just a constant field (with a smooth switch-on and switch-off)" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "attributes": { "classes": [], "id": "", "n": "11" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEKCAYAAAA4t9PUAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3XmUXGd55/Hvr/eW1KW1pRgtlsAKYGzHYFkmITgnBhOZIZaTYGKzmTkGJQNOmOGQxGQxiQM5MJkTMkmcBAMGs8UYEwYliDieYMiwGbUXbMuO4vYqyUaS1bK6tfSmfuaPutUqN93qalXdure6f59z6nTVe5d6SmXfp97lvq8iAjMzs1PVlHUAZmbW2JxIzMysKk4kZmZWFScSMzOrihOJmZlVxYnEzMyq4kRiZmZVcSIxM7OqOJGYmVlVWrIOoB6WLVsWa9euzToMM7OGcvfddz8bEd3T7TcnEsnatWvp6enJOgwzs4Yi6clK9nPTlpmZVcWJxMzMquJEYmZmVXEiMTOzqjiRmJlZVVJNJJI2SdopqVfStZNsf5+khyTdL+nfJJ1etu0qSY8kj6vKys+T9EByzr+SpDQ/g5mZnVxqiURSM3ADcAlwJnClpDMn7HYvsCEizgFuA/5ncuwS4IPABcBG4IOSFifH/B3wLmB98tiU1mcwM7PppXkfyUagNyIeA5B0C7AZeKi0Q0TcWbb/D4C3Js9/CbgjIvqSY+8ANkn6FlCIiB8k5Z8FLgO+keLnqNrug0f5x3v2MHp8LOtQzKwBvPinCrz+7J+iURpc0kwkK4FdZa93U6xhTOVqTiSEyY5dmTx2T1L+EyRtAbYArFmzZiZx19SRoVEuu+F7PHt4iAb5b8LMMhRR/Ptnv3I2b74gu2vXTOTiznZJbwU2AL9Qq3NGxI3AjQAbNmyIWp13pr7+wDM8e3iIL77zAn7ujGVZhWFmDeL4WHDp33yHm7/3BFduXN0QtZI0O9v3AKvLXq9Kyp5H0muBPwAujYihaY7dkzw/6Tnz5PYHf8zpS+fxsy9amnUoZtYAmpvEFeevZufeAR5/9kjW4VQkzUSyHVgvaZ2kNuAKYGv5DpJeDnycYhLZV7bpduB1khYnneyvA26PiGeAfkmvTEZrvR34WoqfoSpjY0HPkwf52RcubYhfFWaWD6Ufnj1PHsw4ksqklkgiYhS4hmJSeBi4NSJ2SLpe0qXJbn8OLAC+LOk+SVuTY/uAP6WYjLYD15c63oF3A58EeoFHyXFH+6P7D3Po2Ajnnb54+p3NzBIvXLaARfNa6Xmib/qdcyDVPpKI2AZsm1B2Xdnz157k2JuAmyYp7wHOqmGYqXnomX4Azlm1KONIzKyRNDWJs1cu5OFnBrIOpSK+sz1Fj+w9THOTWLdsftahmFmDWb+8i959hxkby2ysUMWcSFL0n3sHWLt0Hm0t/mc2s5lZv2IBx0aOs+e5Y1mHMi1f4VLUu/8w65d3ZR2GmTWg9csXAPDIvvw3bzmRpGRsLNjdd4zTl83LOhQza0CnLy02iT914GjGkUzPiSQlewcGGT4+xpolTiRmNnPLFrTR2drMroNu2pqzSr8iVi92IjGzmZPE6iWdPNXnGsmcVfoV4RqJmZ2qNUvmscuJZO7akySSFyzqzDgSM2tUqxbPG7+W5JkTSUr2DgyydH6bh/6a2SlbXmhnYGiUo8OjWYdyUr7KpWRf/yDLCx1Zh2FmDWxFV/Easq9/aJo9s+VEkpJ9A0OsKLRnHYaZNbAVyY/RfQNOJHPS3v5Blnc5kZjZqSv9GN3bP5hxJCfnRJKC42PB/oGh8V8TZmanotQ87kQyBx04PMRY4D4SM6tKoaOF9pYmN23NRXuTjrEVbtoysypIYkWhY27XSCRtkrRTUq+kayfZfqGkeySNSnpjWfkvJgtdlR6Dki5Ltn1G0uNl285N8zOcitKX7hqJmVVrRaE994kktYWtJDUDNwAXA7uB7ZK2RsRDZbs9BbwDeH/5sRFxJ3Bucp4lFFdD/NeyXX4nIm5LK/Zq7R0ofuketWVm1Vpe6ODhp/uzDuOk0qyRbAR6I+KxiBgGbgE2l+8QEU9ExP3A2EnO80bgGxGR/3kCEvv6h5Bg2QInEjOrzoqujjndR7IS2FX2endSNlNXAP8woezDku6X9DFJubta7xsYZOn8dlqb3QVlZtVZUWjn8NAoh4fye3d7rq90kk4DzgZuLyv+APAS4HxgCfB7Uxy7RVKPpJ79+/enHmu5ff1DdLuj3cxqoHQt2Z/jWkmaiWQPsLrs9aqkbCbeBHw1IkZKBRHxTBQNAZ+m2IT2EyLixojYEBEburu7Z/i21ek7OszS+W11fU8zm52WJNeSviPDGUcytTQTyXZgvaR1ktooNlFtneE5rmRCs1ZSS0GSgMuAB2sQa031HRke//LNzKqxdH6xRjInE0lEjALXUGyWehi4NSJ2SLpe0qUAks6XtBu4HPi4pB2l4yWtpVij+faEU39B0gPAA8Ay4ENpfYZT5URiZrWyeH4rAAdznEhSG/4LEBHbgG0Tyq4re76dYpPXZMc+wSSd8xFxUW2jrK2R42MMDI6yeJ4TiZlVr/Sj9ECOE0muO9sbUelXw5IFTiRmVr15bS10tDZx8KgTyZzRl3zZS1wjMbMaWTKvbW72kcxVpS+71K5pZlatxfOdSOaU0pddGmlhZlatJU4kc8tB10jMrMacSOaYviPFeyc9asvMamXJ/LZcD/91IqmxviNDFDpaPM+WmdXMknltDAyNMjx6svlts+OrXY31HR3xzYhmVlOLk2tKXocAO5HU2EHf1W5mNbY05/NtOZHU2AEnEjOrscVOJHPLwSPDLHJHu5nVUN6nSXEiqbH+wREWdXror5nVTumacujYyDR7ZsOJpIZGjo9xdPg4BScSM6uh0jWl34lk9hsYLC6FWehIdVJlM5tjOlqbaW9pco1kLij9WnCNxMxqbWFnK4eOOpHMeidqJE4kZlZbCztb52aNRNImSTsl9Uq6dpLtF0q6R9KopDdO2HZc0n3JY2tZ+TpJdyXn/FKyjG8u9A+6RmJm6ZiTiURSM3ADcAlwJnClpDMn7PYU8A7gi5Oc4lhEnJs8Li0r/yjwsYg4AzgIXF3z4E9RqWmry30kZlZjczKRABuB3oh4LCKGgVuAzeU7RMQTEXE/UNEEMpIEXATclhTdDFxWu5Cr4xqJmaVlriaSlcCuste7mWQN9pPokNQj6QeSSsliKfBcRIye4jlT1X/Mo7bMLB2FztbcDv/N8xXv9IjYI+mFwDclPQAcqvRgSVuALQBr1qxJKcTn6x8coUkwvy3P/6xm1ogWdrYyMDTK8bGguUlZh/M8adZI9gCry16vSsoqEhF7kr+PAd8CXg4cABZJKl2ppzxnRNwYERsiYkN3d/fMoz8F/cdG6OpopSlnX7KZNb6FSZP5wGD+aiVpJpLtwPpklFUbcAWwdZpjAJC0WFJ78nwZ8CrgoYgI4E6gNMLrKuBrNY/8FA0MjlLodG3EzGpvYY6nSUktkST9GNcAtwMPA7dGxA5J10u6FEDS+ZJ2A5cDH5e0Izn8pUCPpB9RTBwfiYiHkm2/B7xPUi/FPpNPpfUZZqp/cISudne0m1nt5TmRpPrzOSK2AdsmlF1X9nw7xeapicd9Dzh7inM+RnFEWO70H3ONxMzSsXBefhOJ72yvof7BEd/VbmapyHONxImkhvqPjfgeEjNLhRPJHNE/OOoaiZmlwolkDjg+Fhwech+JmaWjo7WZtpxOJe9EUiOHk5l/u1wjMbOULMzp3e1OJDUyPs+Wp0cxs5Tkdb6taROJpJ+W9G+SHkxenyPpD9MPrbEc8qJWZpayhk0kwCeADwAjAMlsvVekGVQjOlEjcSIxs3Q0ciKZFxE/nFA2Oumec9j4zL/ubDezlDRyInlW0ouAAEhWMnwm1aga0IBrJGaWsryu217Jz+f3ADcCL5G0B3gceGuqUTWg/tJ67e4jMbOUFJKp5MfGIlezjE+bSJK5rV4raT7QFBED6YfVeEpD8ha0u2nLzNKxsLOViOJM46W5t/JgyquepPdNUQ5ARPxFSjE1pOLMvy25W3DGzGaP8rvbGyKRAF3J3xcD53NiLZFfBiZ2vs95xZl/8/PFmtnsk9dpUqZMJBHxJwCS/h14RalJS9IfA1+vS3QNpH9whC7fjGhmKSpdY/K2SmIlo7ZWAMNlr4eTMivjmX/NLG2lUaH9DZhIPgv8UNIfJ7WRu4CbKzm5pE2SdkrqlXTtJNsvlHSPpNFkWHGp/FxJ35e0Q9L9kn69bNtnJD0u6b7kcW4lsaRtwDP/mlnKSveple5by4tKRm19WNI3gFcnRf81Iu6d7jhJzcANwMXAbmC7pK1lS+YCPAW8A3j/hMOPAm+PiEckvQC4W9LtEfFcsv13IuK26WKop/7BEV7S0TX9jmZmp6jU6pG3Gsm0iUTSGuBZ4KvlZRHx1DSHbgR6k+HDSLoF2AyMJ5KIeCLZNlZ+YET8Z9nzpyXtA7qB58gpN22ZWdoWtLUgnbhvLS8q6R3+Osld7UAnsA7YCbxsmuNWArvKXu8GLphpgJI2Am3Ao2XFH5Z0HfBvwLURMTTJcVuALQBr1qyZ6dvOyNhYMDA06pl/zSxVTU2iq70ld1PJT9tHEhFnR8Q5yWM9xZrG99MPDSSdBnyOYnNaqdbyAeAlFIckLwF+b7JjI+LGiNgQERu6u7tTjfPw8CgRvqvdzNLX1dGau6atGa9HEhH3UFnNYg+wuuz1qqSsIpIKFGtDfxARPyh7/2eiaAj4NMXElqnSrwN3tptZ2gqdrY3X2T7hDvcm4BXA0xWcezuwXtI6ignkCuDNlQQlqY1in8xnJ3aqSzotIp5R8Rb7y4AHKzlnmgYGPfOvmdVHoaOlIWskXWWPdoq1hM3THRQRo8A1wO3Aw8CtEbFD0vWSLgWQdL6k3cDlwMcl7UgOfxNwIfCOSYb5fkHSA8ADwDLgQxV+1tSUaiReZtfM0lbobB3/8ZoXlfyEfigivlxeIOly4MtT7D8uIrYB2yaUXVf2fDvFJq+Jx30e+PwU57yogpjranzmXycSM0tZoaOV/mP9WYfxPJXUSD5QYdmcNd5H4qYtM0tZVw6btk42++8lwOuBlZL+qmxTAa+Q+DxeZtfM6qXQ2crhnK1JcrKf0E8DPcClwN1l5QPA/0gzqEZTGkHhSRvNLG2FjpbimiRDo+OzAWftZLP//gj4kaQvJB3nNoWBwRHmtzXT0jzj0dRmZjMyPk3KsZH8JxJJt0bEm4B7JcXE7RFxTqqRNZD+QU+PYmb1UWpCz9PIrZO1xbw3+fuGegTSyPqPjbpZy8zqYnwG4Bx1uJ+saeuZ5O+T9QunMfUPjrij3czqYnxNkhzNt3Wypq0BTkzWCKDktYCIiELKsTWM/sERlnd1ZB2Gmc0BJxa3aoCmrYjw4hoV6j82yhndbtoys/SdWNyqAWok5SS9Avh5ijWS71SysNVc4s52M6uXBe2lddvzUyOZdrxqsu7HzcBSinNbfUbSH6YdWKOICC+za2Z109LcxIL2fN3dXkmN5C3Az0TEIICkjwD3kYPJEvPg6PBxjo+FR22ZWd10deRrcatK7qB7GijvSW5nBuuKzHbj06O4acvM6qSQs8WtKvkZfQjYIekOin0kFwM/LM2/FRG/nWJ8uVeaHsVNW2ZWL4XOllwtblVJIvlq8ij5VjqhNKYTNRI3bZlZfRQ6Wvlx/2DWYYyb9uoXETef6sklbQL+N9AMfDIiPjJh+4XAXwLnAFeUr4Yo6Sqg1Kn/oVIcks4DPgN0Ulzr5L0R8RNTuNSLl9k1s3ordLbyyL7DWYcxrpJRW2+QdK+kPkn9kgYkTbuqiqRm4AbgEuBM4EpJZ07Y7SngHcAXJxy7BPggxbXhNwIflLQ42fx3wLuA9clj03SxpOnEMrtOJGZWH3lbk6SSzva/BK4ClkZEISK6KryrfSPQGxGPRcQwcAsTluiNiCci4n5gbMKxvwTcERF9EXEQuAPYJOk0oBARP0hqIZ+luG57ZkpfpkdtmVm9FFdJHCHDxpjnqSSR7AIePIXmo5XJsSW7k7Jqjl2ZPD+Vc6bixHrtTiRmVh+FzhbGAo4MH886FKCyzvbfBbZJ+jYwVCqMiL9ILaoakLQF2AKwZs2a1N6nf3CUjtYm2luaU3sPM7Ny5RM3lu50z1IlNZIPA0cp3kvSVfaYzh5gddnrVVR+/8lUx+5Jnk97zoi4MSI2RMSG7u7uCt925vqPeeZfM6uvUp9sXqZJqSSVvSAizjqFc28H1ktaR/FifwXw5gqPvR34s7IO9tcBH4iIUof/K4G7gLcDf30KsdWM59kys3o7MQNwPjrcK6mRbJP0upmeOFme9xqKSeFh4NaI2CHpekmXAkg6X9Ju4HLg45J2JMf2AX9KMRltB65PygDeDXwS6AUeBb4x09hqqTjPVvZVSzObO0p9snmZJqWSK+B/A94vaQgYYQbrkUTENor3epSXXVf2fDvPb6oq3+8m4KZJynuAU6khpaL/2AiL57dlHYaZzSHj67bnpEZSyQ2JXpfkJPoHR1mzdH7WYZjZHFIYr5E0Th8JSV/Fesomb4yIf08rqEZS7Gx305aZ1U9XzpbbnfYKKOmdwHspNkHdB7wS+D5wUbqh5V9EuLPdzOquraWJztZmBobyUSOppLP9vcD5wJMR8YvAy4HnUo2qQQyOjDFyPDz818zqLk9rklSSSAbLFrVqj4j/AF6cbliNwTP/mllWCp35WZOkkivgbkmLgP8D3CHpIPBkumE1hoFBz/xrZtkodORnTZJKRm39SvL0jyXdCSwE/iXVqBrEoeRL9DxbZlZvhc5W+o4MZx0GUOGorZKI+HZagTQiL7NrZlkpdLTy5IGjWYcBVNZHYlPwolZmlpVG62y3KfSPL2rlpi0zq69SZ3se1iSpKJFIOl3Sa5PnnZJ8tzuukZhZdgodrYwcDwZHJq4LWH+VLLX7LuA24ONJ0SqKI7jmvIHBUdpamuho9VokZlZfpZaQPAwBrqRG8h7gVUA/QEQ8AixPM6hG0T/otUjMLBula89AgySSoWTNdQAktQDZN8rlgOfZMrOslEaLHsrBvSSVJJJvS/p9oFPSxcCXgX9KN6zG0D84SpeH/ppZBsbXJGmQGsm1wH7gAeA3KK4v8odpBtUoXCMxs6wUcjQD8LSJJCLGIuITEXE5sAW4KyocbyZpk6SdknolXTvJ9nZJX0q23yVpbVL+Fkn3lT3GJJ2bbPtWcs7Stsz6azzzr5ll5URnewM0bSUX7oKkJcDdwCckfayC45qBG4BLgDOBKyWdOWG3q4GDEXEG8DHgowAR8YWIODcizgXeBjweEfeVHfeW0vaI2FfB50xFcZldJxIzq7+GqpEACyOiH/hV4LMRcQHwmgqO2wj0RsRjSWf9LcDmCftsBm5Ont8GvEaSJuxzZXJs7vQfG/HNiGaWiY7WZtpamhhohBoJ0CLpNOBNwD/P4NwrgV1lr3cnZZPuExGjwCFg6YR9fh34hwlln06atf5oksRTF4MjxxkaHXONxMwyU+hoaZjO9uuB2ynWLrZLeiHwSLphFUm6ADgaEQ+WFb8lIs4GXp083jbFsVsk9Ujq2b9/f81jK/0KcGe7mWWl0NHaGE1bEfHliDgnIt6dvH4sIn6tgnPvAVaXvV6VlE26T3J/ykLgQNn2K5hQG4mIPcnfAeCLFJvQJov7xojYEBEburu7Kwh3Zjzzr5llrauzNRed7VP+nJb015zkxsOI+O1pzr0dWC9pHcWEcQXw5gn7bAWuorgG/BuBb5ZGhElqotic9uqymFqARRHxrKRW4A3A/50mjlR4ni0zy1ohJzMAn6xdpqeaE0fEqKRrKDaLNQM3RcQOSdcDPRGxFfgU8DlJvUAfxWRTciGwKyIeKytrB25PkkgzxSTyiWriPFUDnvnXzDJW6Gzl6eeOZR3G1IkkIm6ealulImIbxRsYy8uuK3s+CFw+xbHfAl45oewIcF61cdVCv5fZNbOMFTpy3rRVkiyv+xNNXBFxUSoRNYhDSXWyy4nEzDLSCE1bJe8ve94B/BqQfQrMWCmRLHRnu5llpNDZytDoGIMjxzNdzmLaRBIRd08o+q6kH6YUT8M4dGyEtuYmOlq9yKSZZaN0+8HA4Gi+E0kyNUpJE8U+ioWpRdQgine1t5LR/ZBmZuO3H/QPjtDd1Z5ZHJU0bd1NsY9EFJu0Hqc4R9ac1n9slIUesWVmGTqxuFW2vQ2VNG2tq0cgjebQsRH3j5hZpsbXJMm4w72Spq0O4N3Az1Osmfw/4O+Tobtz1qFjIyxb0JZ1GGY2h5U3bWWpkp7izwIvA/4a+Jvk+efSDKoRuEZiZlk7MZV8zpu2gLMionwdkTslPZRWQI3CicTMsnZicav810jukTR+h3kyI29V06c0urGxoH/QicTMstXZ2kxLkxjIOJFUUiM5D/iepKeS12uAnZIeACIizkktupwaGBolwjP/mlm2JNHV0dIQTVubUo+iwYzP/OtEYmYZK3S2Zt60Vcnw3yfrEUgj8fQoZpYXeVjcyvN7nAInEjPLi0JnS+YzADuRnAInEjPLC9dIGpQTiZnlRaGjNfMpUlJNJJI2SdopqVfStZNsb5f0pWT7XZLWJuVrJR2TdF/y+PuyY86T9EByzF8pg1kTnUjMLC+6Oloy72xPLZFIagZuAC4BzgSulHTmhN2uBg5GxBnAx4CPlm17NCLOTR6/WVb+d8C7gPXJo+6jyg4dG6GlScxry27aZjMzKI7aOjp8nJHjY5nFkGaNZCPQGxGPRcQwcAuwecI+m4HSkr63Aa85WQ1D0mlAISJ+EBFBcfqWy2of+smV7mr3FPJmlrXyNUmykmYiWQnsKnu9OymbdJ+IGAUOAUuTbesk3Svp25JeXbb/7mnOCYCkLZJ6JPXs37+/uk8ygadHMbO8GJ+4McMO97x2tj8DrImIlwPvA74oqTCTE0TEjRGxISI2dHd31zS40qJWZmZZy8OaJGkmkj3A6rLXq5KySfeR1EJx5cUDETEUEQdgfKnfR4GfTvZfNc05U9fvGomZ5cT4miQZdrinmUi2A+slrZPUBlwBbJ2wz1bgquT5G4FvRkRI6k4665H0Qoqd6o9FxDNAv6RXJn0pbwe+luJnmJSbtswsL/LQtJXaWrERMSrpGuB2oBm4KSJ2SLoe6ImIrcCngM9J6gX6KCYbgAuB6yWNAGPAb0ZEX7Lt3cBngE7gG8mjrpxIzCwv8rC4VaqLjkfENmDbhLLryp4PApdPctxXgK9Mcc4e4KzaRlq5iKB/cHR8HQAzsywVxpfbnZ19JLPS4aFRjo+FayRmlgvz21poEpmuSeJEMkO+q93M8qSpSXR1tGY6caMTyQw5kZhZ3nR1tIxfm7LgRDJDh46WEklbxpGYmRUtmtfqRNJIDhwZBmDpAicSM8uHxfPaxq9NWXAimaGDR4tf1pL5TiRmlg9L57dx0ImkcfQlX9Yi95GYWU4sdiJpLAePDLOws5WWZv/TmVk+LJnXxsDQKMOj2Uwl76vhDB04MuxmLTPLlcXJNanU9F5vTiQzdPCoE4mZ5cvS5JrUl1HzlhPJDPUdGWHxPCcSM8uP8RqJE0ljOHhkmCXz3dFuZvlRaiXpc9NW/kUEfUeGWTK/PetQzMzGLXHTVuM4Mnyc4eNjrpGYWa6UbkdwImkApfZH95GYWZ60NDexsLPVfSSNoDQFgUdtmVneLJmf3TQpqSYSSZsk7ZTUK+naSba3S/pSsv0uSWuT8osl3S3pgeTvRWXHfCs5533JY3man6HcQScSM8upJfPbMruPJLVl/pI1128ALgZ2A9slbY2Ih8p2uxo4GBFnSLoC+Cjw68CzwC9HxNOSzqK4XO/KsuPekqyUWFd9TiRmllOL57Wx57ljmbx3mjWSjUBvRDwWEcPALcDmCftsBm5Ont8GvEaSIuLeiHg6Kd8BdErKfKjUgSNDwIkx22ZmebF0fhsHDg9l8t5pJpKVwK6y17t5fq3ieftExChwCFg6YZ9fA+6JiPJ/oU8nzVp/JEmTvbmkLZJ6JPXs37+/ms8xbl//EJ2tzXS1e712M8uX5YV2nj08xPGxqPt757qzXdLLKDZ3/UZZ8Vsi4mzg1cnjbZMdGxE3RsSGiNjQ3d1dk3j2DQyxvNDOFLnLzCwzywsdjMWJlpN6SjOR7AFWl71elZRNuo+kFmAhcCB5vQr4KvD2iHi0dEBE7En+DgBfpNiEVhd7+wdZ0dVRr7czM6vY8q5i6/++/tmVSLYD6yWtk9QGXAFsnbDPVuCq5PkbgW9GREhaBHwduDYivlvaWVKLpGXJ81bgDcCDKX6G5ynVSMzM8mZFofgjd2//YN3fO7VEkvR5XENxxNXDwK0RsUPS9ZIuTXb7FLBUUi/wPqA0RPga4AzgugnDfNuB2yXdD9xHsUbzibQ+w0T7+gdZ7hqJmeXQiuRH7r6B+tdIUu01johtwLYJZdeVPR8ELp/kuA8BH5ritOfVMsZKHR4a5cjw8fEvy8wsT5YtaEeaZTWS2ab05bhpy8zyqLW5iaXz29g7y/pIZpVSB5Y7280sr7q7Otg/4BpJbu0bcI3EzPJtRaHdNZI8K9VIlhdcIzGzfFrR1TH+o7eenEgqtLd/kI7WJt/Vbma5tbzQzv6B+t/d7kRSoR/3D7Ki0OG72s0st0p3tz9b5zm3nEgqtPvgMVYu6sw6DDOzKa1KrlG7Dx6t6/s6kVRo98GjrF48L+swzMymtHpJMZHs6qvvdPJOJBU4OjzKs4eHWbPUicTM8mtV8mN3V59rJLmz+2Axu69a7KYtM8uvjtZmlne185QTSf48eaD4paxZ4hqJmeXb6UvnjV+z6sWJpAKP7BsA4IzlCzKOxMzs5F7UvYDe/Yfr+p5OJBXo3XuY0xZ20NXRmnUoZmYndcbyBfQdGa7rsrtOJBV4ZN9h10bMrCGsX9EFFK9b9eJEMo3BkePs/PEALz2tkHUoZmbTeulpxUTy4J5DdXvPVBOJpE2SdkrqlXTtJNvbJX0p2X6XpLVl2z6QlO+U9EuVnrPWdjx9iOHjY5x3+uK038rMrGrLuzpYvaSgM2h4AAAF9UlEQVSTu588WLf3TC2RSGoGbgAuAc4ErpR05oTdrgYORsQZwMeAjybHnklxad6XAZuAv5XUXOE5a+qux/sAeMUaJxIzawznrVnM9if66jbnVpo1ko1Ab0Q8FhHDwC3A5gn7bAZuTp7fBrxGxcmsNgO3RMRQRDwO9Cbnq+ScNRMRfO3ep/mZ1Yvo7vL08WbWGF575gqePTzMd3ufrcv7pZlIVgK7yl7vTsom3SdZ4/0QsPQkx1ZyzpqICH7/qw+yc+8Ab964Oo23MDNLxWtfuoJlC9r43dvup7cOne6ztrNd0hZJPZJ69u/ffyrH86Lu+fzWRWfwpg1OJGbWODpam/n8Oy9g/YoFdC9IvzUlzcU19gDlV+BVSdlk++yW1AIsBA5Mc+x05wQgIm4EbgTYsGHDKTUUvvPVLzyVw8zMMveSnyrwuasvqMt7pVkj2Q6sl7ROUhvFzvOtE/bZClyVPH8j8M2IiKT8imRU1zpgPfDDCs9pZmZ1lFqNJCJGJV0D3A40AzdFxA5J1wM9EbEV+BTwOUm9QB/FxECy363AQ8Ao8J6IOA4w2TnT+gxmZjY9FSsAs9uGDRuip6cn6zDMzBqKpLsjYsN0+83aznYzM6sPJxIzM6uKE4mZmVXFicTMzKriRGJmZlWZE6O2JO0HnjzFw5cB9ZmwJj/8mecGf+bZr9rPe3pEdE+305xIJNWQ1FPJ8LfZxJ95bvBnnv3q9XndtGVmZlVxIjEzs6o4kUzvxqwDyIA/89zgzzz71eXzuo/EzMyq4hqJmZlVxYnkJCRtkrRTUq+ka7OOJ22SbpK0T9KDWcdSD5JWS7pT0kOSdkh6b9YxpU1Sh6QfSvpR8pn/JOuY6kVSs6R7Jf1z1rHUg6QnJD0g6T5Jqc5a66atKUhqBv4TuJjikr7bgSsj4qFMA0uRpAuBw8BnI+KsrONJm6TTgNMi4h5JXcDdwGWz/DsWMD8iDktqBb4DvDcifpBxaKmT9D5gA1CIiDdkHU/aJD0BbIiI1O+bcY1kahuB3oh4LCKGgVuAzRnHlKqI+HeK68LMCRHxTETckzwfAB4GVmYbVbqiqLSId2vymPW/JiWtAv4L8MmsY5mNnEimthLYVfZ6N7P8IjOXSVoLvBy4K9tI0pc08dwH7APuiIhZ/5mBvwR+FxjLOpA6CuBfJd0taUuab+REYnOepAXAV4D/HhH9WceTtog4HhHnAquAjZJmdTOmpDcA+yLi7qxjqbOfj4hXAJcA70marlPhRDK1PcDqsterkjKbRZJ+gq8AX4iIf8w6nnqKiOeAO4FNWceSslcBlyZ9BrcAF0n6fLYhpS8i9iR/9wFfpdhcnwonkqltB9ZLWiepjeJ68lszjslqKOl4/hTwcET8Rdbx1IOkbkmLkuedFAeT/Ee2UaUrIj4QEasiYi3F/4+/GRFvzTisVEmanwwgQdJ84HVAaqMxnUimEBGjwDXA7RQ7YW+NiB3ZRpUuSf8AfB94saTdkq7OOqaUvQp4G8VfqPclj9dnHVTKTgPulHQ/xR9Ld0TEnBgOO8esAL4j6UfAD4GvR8S/pPVmHv5rZmZVcY3EzMyq4kRiZmZVcSIxM7OqOJGYmVlVnEjMzKwqTiRmNSZpkaR3J89fIOm2rGMyS5OH/5rVWDJv1z/PhRmUzQBasg7AbBb6CPCiZGLER4CXRsRZkt4BXAbMB9YD/wtoo3hT5BDw+ojok/Qi4AagGzgKvCsiZvXd59bY3LRlVnvXAo8mEyP+zoRtZwG/CpwPfBg4GhEvpzijwNuTfW4EfisizgPeD/xtXaI2O0WukZjV153J2icDkg4B/5SUPwCck8xE/HPAl4tTgQHQXv8wzSrnRGJWX0Nlz8fKXo9R/P+xCXguqc2YNQQ3bZnV3gDQdSoHJuuhPC7pcijOUCzpZ2oZnFmtOZGY1VhEHAC+K+lB4M9P4RRvAa5OZm7dwSxf4tkan4f/mplZVVwjMTOzqjiRmJlZVZxIzMysKk4kZmZWFScSMzOrihOJmZlVxYnEzMyq4kRiZmZV+f/jAYLJqbvL0wAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plot_pulse(H[1][1], tlist)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The next line solves the equation of motion for the defined objective, which\n", "contains the initial state $\\ket{\\Psi_{\\init}}$ and the Hamiltonian $\\op{H}(t)$\n", "defining its evolution." ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "attributes": { "classes": [], "id": "", "n": "12" } }, "outputs": [], "source": [ "guess_dynamics = objectives[0].mesolve(tlist, e_ops=[proj0, proj1])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The plot of the population dynamics shows that the guess field does not transfer\n", "the initial state $\\ket{\\Psi_{\\init}} = \\ket{0}$ to the desired target state\n", "$\\ket{\\Psi_{\\tgt}} = \\ket{1}$." ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "attributes": { "classes": [], "id": "", "n": "13" } }, "outputs": [], "source": [ "def plot_population(result):\n", " fig, ax = plt.subplots()\n", " ax.plot(result.times, result.expect[0], label='0')\n", " ax.plot(result.times, result.expect[1], label='1')\n", " ax.legend()\n", " ax.set_xlabel('time')\n", " ax.set_ylabel('population')\n", " plt.show(fig)" ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "attributes": { "classes": [], "id": "", "n": "14" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEKCAYAAAD9xUlFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3XlwXed93vHvD+sFQOwANywEKdHUQlEiCVJy3DhuG9syY0uZbJayNI7kKNNaqdOk7ijT1Imc6YyzNJPUSxolcVMrtTROnMSMrMjRNFKcxFoIUpREUqJEkwQBcAOxg9gvfv3jPTgAIZK4pHBwsTyfmTPnnuXe+14S9z7nfd9z3mPujoiICEBOtgsgIiKLh0JBRERiCgUREYkpFEREJKZQEBGRmEJBRERiCgUREYkpFEREJKZQEBGRWF62C3CtampqvKmpKdvFEBFZUvbv33/B3Wvn2m/JhUJTUxMtLS3ZLoaIyJJiZq2Z7KfmIxERiSkUREQkplAQEZGYQkFERGIKBRERiSUWCmb2FTM7b2aHrrDdzOx/mtkxM3vNzHYkVRYREclMkjWFPwPuvsr2jwCbo+kh4A8TLIuIiGQgsesU3P07ZtZ0lV3uBb7q4X6gL5pZhZmtc/czSZRn38lu/umtTjAjx8AI85wcAyDHDDPibeHxpfO8nByKC3IpKsilpCCPooJciqPHxYW5lBflk5+rFjkRWbqyefFaHdA2Y7k9WveOUDCzhwi1CRobG6/rzQ609vCF546R9C2pSwvzqCjJp7K4gPKiMK8szqe2tJA1ZSnWlKVYW55iTWmKsqI8zCzZAomIXIMlcUWzuz8GPAbQ3Nx8XT/rv/ADN/ALP3DD1Osx6dPzySgpJt3xaNkBn5x+POnORNoZGptgaCwdTRMMR48vjk3QOzROz9AYPRfH6Bkap3dojNauIXqGxhgYmXhHmVL5OawpS7GuPEVDZTGNVcU0xFMRtasKFRoisqCyGQodQMOM5fpoXeLMjFwDWLgf3JHxNOf7RznbP8K5S6ZRTvcO8523OznXP3rJc4ryc2moKmJDdQk31K7ixtVhuqG2hNJU/oKVXURWjmyGwl7gYTN7ErgT6EuqP2ExSOXn0lhdTGN18RX3GRlP094zRFv3MKe6h+Lp5IWLPH/0POPp6UrSmrLCEBJRWNy0royb1pYqLETkXUksFMzsCeADQI2ZtQO/DuQDuPv/Ap4G9gDHgCHg55Iqy1KRys/lxtWl3Li69B3bxtOTnOoe4nvnBznWOcix84N8r/Mi3zjQweDodNNUY1UxN68r5ZZ15WG+voy6iiI1Q4lIRsyT7nmdZ83Nza5RUqe5O2f7R3jzzABHzvRz5Ew/b5zu50TXxbhTvSyVx83ryritrpxtDRXcUV9BQ5WCQmQlMbP97t48135LoqNZrszMWFdexLryIv71Tavj9UNjE7x5doA3zvRz5HQIi8dfbGX0n08AUFmcz7b6Cm6vLw/zhgpqSwuz9TFEZJFQKCxTxQV57GisZEdjZbxuPD3JW+cGeLWtj9faeznY1ssXn+tkMqpRrC9PcXtDBXc0VLBzQyVb68pJ5edm6ROISDao+WiFGxqb4PDpfl5t6+XV9j5ebevlVPcQAAW5OWytK2Pnhkp2bqhi54ZK1SZE5lnf8Dhn+obpHBilc2CUrsExhsfTjE1MMpaeZHQ8HR+4fez29ezeWHVd76PmI8lIcUEeu5qq2NU0/Yd2YXCUA6097I+m//NCK3/8T6HZaUN1MTsbK9nZVMnODZVsXl1Kbo76JkTm4u60dQ+z72Q3r3f08fb5Ad4+N8j5gdHL7p9jUJCXQ0FuDnnRSAnb6suvOxQypZqCzGl0Is2hjn4OtPbQ0trN/tYeLgyOAeEK7u0bKmneUElzUyV3NFRQXKBjDZGJ9CRvnBlg38luWlq7aTnZEwdASUEuN64pZXN07VFjVTE1qwqpLS2kelUBxfm5cRDMl0xrCgoFuWbuzqnuobgmsb+1h6PnBnCHvBzj1vVlNDdVsaspNDupyUlWgoujExxs6w0hcLKHA6d6GBpLA1BXUcSupsroe1HF5tWr4nHXFopCQRZU3/A4B0710HKym30ne3i1rZfRiUkAmqqL45BobqpiU02JToeVJe98/wgtrT1xCBw500960jGDm9eWxX/vzU2VrCsvynZxFQqSXWMTkxw63ReHRMvJbnqGxgGoKilg54bK+EuzdX05BXkaXVYWL3fn+IWLtJzs5uUToRm1tSuckJHKz2F7Q/h73tlUxfbGCsoW4cgCCgVZVGZ+qaZC4mT0pSrMy+H2hoo4JHY0VlJetPi+VLJyTKQnOXKmn5dPhFpAS2t33I9WHR3U7N5YRXNTFbeuL1sSQ+YrFGTR6xwYZX/rdEgcPt3PRFT93rKmlOamSnY1hS9eXUX2q9+yfA2PpXmlrYd9J0Jz0Mz+gMaqYpqbKtndVMWujUu3+VOhIEvO0FjoqGs5Gb6Yr5zqjcd1Wl+eYudUv8SGKras1amwcv16Lo5FZwX18PKJbg519MUHJDdF/QFTp2qvLU9lu7jzQqEgS1560nnzbH8cEi0nezjbPwKEU2F3xKfCVnFHQwVFBbr6Wt4pPekcOz/IgVM9HGjt4ZW2Xo6dHwTCBZq3N5THAbBjw/JtulQoyLLj7nT0Dl8SEkfPDQDhVNitdeVxSDQ3VVKzSqfCrkS9Q2O8cqqXA6d6eOVUGM5lqsZZWZzP9sbKuE/gthU0lItCQVaEvqFwKuxUSBxs72UsOhV2U00JOzdUsmNDJbfVlfOeNaU6y2mZGZuY5OjZAV5t7+WVU728cqqH4xcuAuGK4JvWlrFjQwU7GivZ3lhJU3XxkuwPmA8KBVmRpq6+njrLaX/r9KmwBbk53LSulK115dwWTQqKpWMqAF7v6OP1jj4OdfRx9OwAY+lwEFBdUsD2xkq2N4YQ2FZfTkmhrq6folAQYfrq69c7+ni9vS/+QZm6Z/bsoNiytpQta0r1Y5Jlg6MTvHVugDfPDFw2AMpSedxWXx7/v22r0z1C5qJQELkCd6e1aygOiNfb+zh0ejooABqqitiyJtzidMvaUm5aW8rGmpJ5H49mpRtPT3LiwkXePDvA0bP9HD07wJtnB2jvGY73mRkA2+oquK2uXAFwHRQKItdgctJp6xni6NmB8MN0LsxPXLhIOhq3uCA3h021JdxQu4qNNSU01ZSwsaaETTUlVJYUZPkTLG49F8c4fmGQ450XOXHhYjw/ceFifPSfm2NsqimJQ3jL2hDK9ZUKgPmgobNFrkFOjrGhuoQN1SV86Na18fqR8TTf6xwMTRlnB3jrbLjt6TOHz8ZhAVBelB8HRENVMXWVRdRVFLG+ooh15allf4bLyHiajt5hOnqG43l7zxCnuoc4fuEivVG/DoQzxRqri9lUs4oP3FQbAmBNGTesLqEwb3n/Oy0FCgWRq0jl53Lr+nJuXV9+yfrx9CRt3UOc7ApHvSe7wlHvi8e7+OuDHcyugNesKqSuIsX6iiLWlKWoWVVAzapCalaFoZKnhk1ebOExMp6m6+JYfAOYqenCYJif6R+ho2eYC4OX3hMgN8dYW5aioaqIPbetY1NNCZtqS9hUs4r6yiI1wy1iaj4SmWdjE5Oc6x+Jj5hP9w5zum+Y9ugourN/lIHRics+t7ggl9JUHmWpfEpTeZTOmJel8ijMz6UwuvFKQV5OeBxNuWZMfZunvtYerUlPOqPjk4xMpMN8PB0/Hh5PMzAyQf/IOH3D4/QPj9M/MkHf8Hh8eu9slcX51JYWsro0RV1FEfWVRXHtqL6qmDWlhfrhX2TUfCSSJQV5OTRUFdNQVXzFfaaOwC9ER91dg2N0Do7SfXGMgZFxBkYmGBiZoHdojLbuIfpHJhgYGY+HI58PORZqQqn8XMpSeZQV5VOWymd9eRFlRSGYyoryqS4poLa0MJ6qSwp1Gu8yplAQyYJUfi51FUXXPNCfuzOedkYnZt7DN8wno+qBETplp/pmDTAzUvk5pKKaRio/d0mM7CkLT6EgsoSYGQV5piN1SYz+skREJKZQEBGRmEJBRERiCgUREYkpFEREJKZQEBGRWKKhYGZ3m9lRMztmZo9cZnujmT1nZq+Y2WtmtifJ8oiIyNUlFgpmlgt8CfgIcAtwv5ndMmu3XwO+7u7bgfuALydVHhERmVuSNYXdwDF3P+7uY8CTwL2z9nGgLHpcDpxOsDwiIjKHJK9orgPaZiy3A3fO2uc3gL83s18ESoAfTLA8IiIyh2x3NN8P/Jm71wN7gMfN7B1lMrOHzKzFzFo6OzsXvJAiIitFkqHQATTMWK6P1s30IPB1AHd/AUgBNbNfyN0fc/dmd2+ura1NqLgiIpJkKOwDNpvZRjMrIHQk7521zyng3wKY2c2EUFBVQEQkSxILBXefAB4Gvg28QTjL6LCZfc7M7ol2+xXg583sVeAJ4BO+1O76IyKyjCQ6dLa7Pw08PWvdZ2c8PgK8L8kyiIhI5rLd0SwiIouIQkFERGIKBRERiSkUREQkplAQEZGYQkFERGIKBRERiSkUREQkplAQEZGYQkFERGIKBRERiSkUREQkplAQEZGYQkFERGIKBRERiSkUREQkplAQEZGYQkFERGIKBRERiSkUREQkplAQEZGYQkFERGIKBRERiSkUREQkplAQEZGYQkFERGIKBRERiSkUREQkplAQEZGYQkFERGJ5me5oZrnAmpnPcfdTSRRKRESyI6Oagpn9InAOeBb4VjQ9lcHz7jazo2Z2zMweucI+P2FmR8zssJl97RrKLiIi8yzTmsKngS3u3pXpC0c1iy8BHwTagX1mttfdj8zYZzPwq8D73L3HzFZnXnQREZlvmfYptAF91/jau4Fj7n7c3ceAJ4F7Z+3z88CX3L0HwN3PX+N7iIjIPMq0pnAceN7MvgWMTq1099+7ynPqCGEypR24c9Y+7wEws38BcoHfcPdnZr+QmT0EPATQ2NiYYZFFRORaZRoKp6KpIJrm8/03Ax8A6oHvmNlt7t47cyd3fwx4DKC5udnn8f1FRGSGjELB3R8FMLNV0fJgBk/rABpmLNdH62ZqB15y93HghJm9RQiJfZmUS0RE5lemZx9tNbNXgMPAYTPbb2a3zvG0fcBmM9toZgXAfcDeWfv8DaGWgJnVEJqTjl9D+UVEZB5l2nz0GPDL7v4cgJl9APhj4Puu9AR3nzCzh4FvE/oLvuLuh83sc0CLu++Ntn3IzI4AaeAz13KGk4jIQhofH6e9vZ2RkZFsF+WKUqkU9fX15OfnX9fzzX3uJnoze9Xdb59r3UJobm72lpaWhX5bERFOnDhBaWkp1dXVmFm2i/MO7k5XVxcDAwNs3Ljxkm1mtt/dm+d6jUxPST1uZv/NzJqi6ddQM4+IrDAjIyOLNhAAzIzq6up3VZPJNBQeAGqBv4qm2midiMiKslgDYcq7LV9GoeDuPe7+H919RzR9euqCMxERWTjPPPMMW7Zs4cYbb+Tzn//8vL/+VTuazez33f2XzOxvgXd0Prj7PfNeIhERuax0Os2nPvUpnn32Werr69m1axf33HMPt9xyy7y9x1xnHz0ezX933t5RRESuy8svv8yNN97Ipk2bALjvvvv45je/uXCh4O77o4d3uPsfzNxmZp8G/nHeSiIisoQ8+reHOXK6f15f85b1Zfz6x658CVhHRwcNDdPXBNfX1/PSSy/Naxky7Wj+2cus+8Q8lkNERBaBufoU7gd+EthoZjOvRi4FupMsmIjIYna1I/qk1NXV0dY2Pc5oe3s7dXV18/oec/UpfBc4A9QA/2PG+gHgtXktiYiIXNWuXbt4++23OXHiBHV1dTz55JN87Wvze2+yufoUWoFW4L3z+q4iInLN8vLy+OIXv8iHP/xh0uk0DzzwALfeOr81lozGPjKzu4AvADcThs7OBS66e9m8lkZERK5qz5497NmzJ7HXz7Sj+YvA/cDbQBHwScKtNkVEZBnJNBRw92NArrun3f1/A3cnVywREcmGTIfOHoruiXDQzH6b0PmccaCIiMjSkOkP+88Q+hEeBi4S7qj2o0kVSkREsiPT23G2Rg+HgUeTK46IiGTTXBevvc5lBsKb4u7b5r1EIiKSNXPVFD66IKUQEZGMPPDAAzz11FOsXr2aQ4cOzfvrX7VPwd1brzbNe2lEROSqPvGJT/DMM88k9voZdTSb2YCZ9UfTiJmlzWx+hwcUEZE5vf/976eqqiqx18+0o7l06rGFe73dC9yVVKFERBa9v3sEzr4+v6+59jb4yPzfTe1aXPO1Bh78DfDhBMojIiJZlOnYRz8yYzEHaAZGEimRiMhSkOUj+qRkekXzx2Y8ngBOEpqQRERkGcm0T+Hnki6IiIjM7f777+f555/nwoUL1NfX8+ijj/Lggw/O2+tn2ny0CfgDQueyAy8A/8ndj89bSUREZE5PPPFEoq+faUfz14CvA+uA9cBfAMmWTEREFlymoVDs7o+7+0Q0/TmQSrJgIiKy8DLtaP47M3sEeJLQfPRx4GkzqwJw9+6EyiciIgso01D4iWj+C7PW30cIiU3zViIRkUXM3QnX8C5O7lccwzQjGTUfufvGq0xXDAQzu9vMjprZsaimcaX9ftTM3Myar+dDiIgshFQqRVdX17v+4U2Ku9PV1UUqdf2t+5mefZQP/Hvg/dGq54E/cvfxqzwnl3Af5w8C7cA+M9vr7kdm7VcKfBp46ZpLLyKygOrr62lvb6ezszPbRbmiVCpFfX39dT8/0+ajPwTygS9Hyz8TrfvkVZ6zGzg2ddqqmT1JuODtyKz9fhP4LeAzGZZFRCQr8vPz2bhxY7aLkahMQ2GXu98+Y/kfzOzVOZ5TB7TNWG4H7py5g5ntABrc/VtmplAQEcmyTE9JTZvZDVML0cVs6XfzxmaWA/we8CsZ7PuQmbWYWctirraJiCx1mdYUPgM8Z2ZTVzA3AXMNfdEBNMxYro/WTSkFtgLPRz35a4G9ZnaPu7fMfCF3fwx4DKC5uXlx9vCIiCwDmdYU/gX4I2AS6I4evzDHc/YBm81so5kVEE5f3Tu10d373L3G3ZvcvQl4EXhHIIiIyMLJNBS+CmwkdAp/gXBdwuNXe4K7TwAPA98G3gC+7u6HzexzZnbP9RdZRESSkmnz0VZ3v2XG8nNmNvssondw96eBp2et++wV9v1AhmUREZGEZFpTOGBm8e03zexOQM08IiLLTKY1hZ3Ad83sVLTcCBw1s9cJd+jclkjpRERkQWUaCncnWgoREVkUMr3zWmvSBRERkezLtE9BRERWAIWCiIjEFAoiIhJTKIiISEyhICIiMYWCiIjEFAoiIhJTKIiISEyhICIiMYWCiIjEFAoiIhJTKIiISEyhICIiMYWCiIjEFAoiIhJTKIiISEyhICIiMYWCiIjEFAoiIhJTKIiISEyhICIiMYWCiIjEFAoiIhJTKIiISEyhICIiMYWCiIjEFAoiIhJLNBTM7G4zO2pmx8zskcts/2UzO2Jmr5nZ/zOzDUmWR0REri6xUDCzXOBLwEeAW4D7zeyWWbu9AjS7+zbgL4HfTqo8IiIytyRrCruBY+5+3N3HgCeBe2fu4O7PuftQtPgiUJ9geUREZA5JhkId0DZjuT1adyUPAn93uQ1m9pCZtZhZS2dn5zwWUUREZloUHc1m9tNAM/A7l9vu7o+5e7O7N9fW1i5s4UREVpC8BF+7A2iYsVwfrbuEmf0g8F+BH3D30QTLIyIic0iyprAP2GxmG82sALgP2DtzBzPbDvwRcI+7n0+wLCIikoHEQsHdJ4CHgW8DbwBfd/fDZvY5M7sn2u13gFXAX5jZQTPbe4WXExGRBZBk8xHu/jTw9Kx1n53x+AeTfH8REbk2i6KjWUREFgeFgoiIxBQKIiISUyiIiEhMoSAiIrFEzz4SkQUymQafjBYsmtn0stmMZZErUyiIZIs7jA3C4Hm4eAEunoehLhjpg5F+GO1/53xiBCZGw5QehYmxsM7Tc79fXmp6yk9BXhHkFUJ+MaTKIFUeTRXTj4sqoLgGVq2GVWugcFXy/y6SVQoFkSRMTsLQBehtg7426Gufng+cgYudMNgJE8OXf77lQGFZ+LEuLA/zsvWQXwS5hZBXEM2jKbcQcnLAp14geuDRfHIiCpRoGh8J7z0xCuPD0H8azr8BI70hgKZf6FL5JdMBsWo1lK6F8gaoaISKBihvhJIa1UqWMIWCyPVKT0DfKej6XjQdC1NvK/R1hCP5mQpWhR/Q0rVQfSOU1IYf1pJaKFkNq2qhuDocoResyt4P6+QkjA3AcG8IiYsXohA7F2o1g+fCdOEtOP6PMNp36fPziqKAaIDKpvBZazZD9Q1QsQFycrPysSQzCgWRuYwNQeebcP5IOJqeCoCekzA5Pr1fYXn44Vt3B9z8sfCjWF4fTQ3hx34pHEHn5Ew3H5HBzRCHe0MtaKpW1HtqeupoCc1hU3ILoHJjFBQ3Qs17YPUtUHsTFBQn9pEkcwoFkSnpCeg+Hv34H4Fzh8O8+wRxc0peCqpugNU3w80fDT9uU1Nx9dL40Z9vRRVhWnvbO7e5h36SrmNw4e3p2lTXMTj2LKTHoh0NqjbBmltg9a3T86qNqlksMIWCrEzp8XD0f/oVOH0QzhwMITAxErZbTvQjtRW2fTwcza65NTSH6Ecqc2ahj6GkBhrvunTbZDoE7vnDcO7I9PyNp5gO4SJYuzXUvtbfAeu3Q80WyNVPV1LM/QodSotUc3Ozt7S0ZLsYspSkJ0IAnDk4HQJnX59u8y8ohXW3h2lNdJRae1Po1JWFNzYEF46GgDh3CM68GqaxwbB9dlCsuyP8fykorsrM9rt785z7KRRk2RnqhraX4NQLcOrF8IMyVQOYCoCpo851d4QaQY6u41zUJidDk9OZg9M1u9lBsX47NOyC+l1QvxtK12S3zItMpqGgaJWlzT10+J56cToELhwN23Lyww9F84Nhvv6O0B+gAFh6cnKg9j1h2vYTYd3MoOg4AO374IUvT3f+VzROB0TDLlhzWziVV65KoSBLS3oCzr0+IwRegsGzYVuqHBruhNs/Do3vDUGgJqDl63JBMT4CZ1+Dtpeh/eXwd3LoG2FbXirUDBt2hb+ThrvCacByCTUfyeI2OhiOAKdCoL0Fxi+GbRWN4YvdeFcIgdqbVAuQd+rrCH9D7ftCWJw5OH3WU/Vm2PBeaPy+8HdU2bRszyBT85EsTf1noO3F6RA4eygM4WA5oRN4+0+FL2/DXVBel+3SylJQXhemW384LE+Mhn6JU98Nf2dHvgkHvhq2la4LBxiN7w1hsfqWFXe2mUJBsmdyMlwVO9UX0PZi6B+A0HFY3wzf/yshBOp3haEeRN6tvEJovDNMEP4OO98If4etL4T54b8K2wrLo33vCrWJuh3h+cuYQkEWzsRoOCU0DoGXYLgnbCupDV+83Q+F+dptkJuf3fLKypAT1ULX3Aq7PhlOXuhriwIiqk28/fdh39zCEAxTtYmG3eHCvWVEfQqSnKHu0IZ76oUQAB0Hpq8NqN483RfQeFc4LXSZtuXKMnCxK2rWjGoTZw6GQQax0MTUeNf0VN6wKP+WdZ2CLCz3MBDczLOCOt8I23Lyw+mgUyHQcGe4wlVkqRobgo7903/vbS+HQQQByuouPeBZJP0S6miWZKXHp0/9m2oKGjgTthWWh2r1bT8Wvhh1O3RqqCwvBcWw8fvDBGHIjnOHpy+abH1h+lTYwrLQJzYVEnU7F/Xgf6opSGaGuqdPDW17ORwlTd0LoLwxhMCGqJ219madGior21S/xMyLKs8fCdty8sKYWqtvhtot4VTqig1hGPWiqsS+O6opyPVzD1eKTtUA2l6ecZVwXugEbv65EAQNd4abv4jINLPoxkON0xfWDfdA274QEqcPhHtRvPrErOflhqbV/OLo5kkF4aK7qean9z4cRudNkEJBQifa6QPh6H9quIDh7rAtVTF9lXDDnbB+x6Ku+oosWkWV8J4PhWnKSB90vgX97dM3MLrYGd0ZbyRcZDcxMn3/7QXom1AorDSjg2EgsZkh0NsabbRQld2yJ5yb3XBnOEtITUEiyUiVh2E32JXtksQUCsvZcG8YevjsoTBU9OkDYQjpqaOOisZw5L/rk6EzeN3tUFia3TKLSFYpFJaDyUnoPRl+/M9FAXD2ULh/8JSS2jBA3M33hLMf1m/XYGAi8g4KhaUkPR6Ggeg8Gjp+L7wdPX57+hxpywlNPg27YdcD4RaJa27T2PIikpFEQ8HM7gb+AMgF/sTdPz9reyHwVWAn0AV83N1PJlmmRS89Dn3t4ce/52Ro7+86Fjqjuo9feqP40vVh2OA7fjJcor92a7hQRtcEiMh1SiwUzCwX+BLwQaAd2Gdme939yIzdHgR63P1GM7sP+C3g40mVKeumbmI+cAYGzkL/6TDvaws//j0nwzC/np5+Tk4+VG4I96W9aU+Y174n1AY0QJyIzLMkawq7gWPufhzAzJ4E7gVmhsK9wG9Ej/8S+KKZmS/mK+rcYXwYxodg7GI0HwrNN8M94SKvoe5wSufM+eD5EAYzj/SnlKwO47g33AXbNoTHFdG8bP2iuEReRFaGJEOhDmibsdwO3Hmlfdx9wsz6gGrgwryX5sDj8N0vRGfeeJj7jPnsdfHyjHWTEyEIyCCzClaF85KLKqG4CqpvCGO1l66DsnXTj1et0S0CRWTRWBIdzWb2EPAQQGNj4/W9SHFVuKzcLHTGEs0tZ9Y6m16+ZD8LTTkFxeFqw4JVsx6XTAdAUeWyH3NdRJanJEOhA2iYsVwfrbvcPu1mlgeUEzqcL+HujwGPQRj76LpKc9MPhUlERK4oyUtV9wGbzWyjmRUA9wF7Z+2zF/jZ6PGPAf+wqPsTRESWucRqClEfwcPAtwmnpH7F3Q+b2eeAFnffC/wp8LiZHQO6CcEhIiJZkmifgrs/DTw9a91nZzweAX48yTKIiEjmNNKZiIjEFAoiIhJTKIiISEyhICIiMYWCiIjEbKldFmBmnUDrnDteXg1JDKGxuOkzrwz6zCvDu/nMG9x9zpuoLLlQeDfMrMXdm7NdjoWkz7wy6DOvDAvxmdV8JCIiMYWCiIjEVlooPJbtAmSyl+ofAAAD1UlEQVSBPvPKoM+8MiT+mVdUn4KIiFzdSqspiIjIVayYUDCzu83sqJkdM7NHsl2epJnZV8zsvJkdynZZFoqZNZjZc2Z2xMwOm9mns12mpJlZysxeNrNXo8/8aLbLtBDMLNfMXjGzp7JdloVgZifN7HUzO2hmLYm+10poPjKzXOAt4IOE24LuA+539yNXfeISZmbvBwaBr7r71myXZyGY2TpgnbsfMLNSYD/ww8v8/9mAEncfNLN84J+BT7v7i1kuWqLM7JeBZqDM3T+a7fIkzcxOAs3unvh1GSulprAbOObux919DHgSuDfLZUqUu3+HcI+KFcPdz7j7gejxAPAG4T7gy5YHg9FifjQt6yM9M6sHfgj4k2yXZTlaKaFQB7TNWG5nmf9YrHRm1gRsB17KbkmSFzWlHATOA8+6+3L/zL8P/BdgMtsFWUAO/L2Z7Y/uWZ+YlRIKsoKY2SrgG8AvuXt/tsuTNHdPu/sdhPug7zazZdtcaGYfBc67+/5sl2WB/St33wF8BPhU1DyciJUSCh1Aw4zl+midLDNRu/o3gP/r7n+V7fIsJHfvBZ4D7s52WRL0PuCeqI39SeDfmNmfZ7dIyXP3jmh+HvhrQpN4IlZKKOwDNpvZRjMrINwLem+WyyTzLOp0/VPgDXf/vWyXZyGYWa2ZVUSPiwgnU7yZ3VIlx91/1d3r3b2J8D3+B3f/6SwXK1FmVhKdOIGZlQAfAhI7q3BFhIK7TwAPA98mdD5+3d0PZ7dUyTKzJ4AXgC1m1m5mD2a7TAvgfcDPEI4eD0bTnmwXKmHrgOfM7DXCwc+z7r4iTtNcQdYA/2xmrwIvA99y92eSerMVcUqqiIhkZkXUFEREJDMKBRERiSkUREQkplAQEZGYQkFERGIKBZGrMLMKM/sP0eP1ZvaX2S6TSJJ0SqrIVURjKD21UkaaFcnLdgFEFrnPAzdEA869Ddzs7lvN7BPADwMlwGbgd4ECwsVzo8Aed+82sxuALwG1wBDw8+6+bK84lqVPzUciV/cI8L1owLnPzNq2FfgRYBfw34Ehd99OuJL830X7PAb8orvvBP4z8OUFKbXIdVJNQeT6PRfdt2HAzPqAv43Wvw5si0Zr/T7gL8KwTAAULnwxRTKnUBC5fqMzHk/OWJ4kfLdygN6oliGyJKj5SOTqBoDS63lidC+HE2b24xBGcTWz2+ezcCLzTaEgchXu3gX8i5kdAn7nOl7ip4AHoxEuD7PMbwMrS59OSRURkZhqCiIiElMoiIhITKEgIiIxhYKIiMQUCiIiElMoiIhITKEgIiIxhYKIiMT+P2JTeMxj0mnmAAAAAElFTkSuQmCC\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plot_population(guess_dynamics)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Optimize\n", "\n", "In the following we optimize the guess field $\\epsilon_{0}(t)$ such\n", "that the intended state-to-state transfer $\\ket{\\Psi_{\\init}} \\rightarrow\n", "\\ket{\\Psi_{\\tgt}}$ is solved.\n", "\n", "The cell below carries out the optimization. It\n", "requires, besides the\n", "previously\n", "defined optimization `objectives`, information\n", "about the\n", "optimization functional\n", "$F$ (via `chi_constructor`) and the\n", "propagation method that should be used. In\n", "addition, the number of total\n", "iterations is required and, as an option, we pass\n", "an info-hook that after each\n", "iteration combines a complete printout of the state\n", "of the optimization with the\n", "`print_fidelity` function defined above.\n", "\n", "Here, we\n", "choose $F = F_{re}$ with\n", "\\begin{equation}\n", "F_{re}\n", "=\n", "\\Re\\Braket{\\Psi(T)}{\\Psi_{\\tgt}}\n", "\\end{equation}\n", "\n", "with\n", "$\\ket{\\Psi(T)}$ the\n", "forward propagated state of $\\ket{\\Psi_{\\init}}$." ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "attributes": { "classes": [], "id": "", "n": "15" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Iteration 0\n", " objectives:\n", " 1:|(2)⟩ - {[Herm[2,2], [Herm[2,2], u1(t)]]} - |(2)⟩\n", " adjoint objectives:\n", " 1:⟨(2)| - {[Herm[2,2], [Herm[2,2], u1(t)]]} - ⟨(2)|\n", " λₐ: 5.00e+00\n", " S(t) (ranges): [0.000000, 1.000155]\n", " duration: 0.6 secs (started at 2018-12-20 08:17:54)\n", " optimized pulses (ranges): [0.00, 0.20]\n", " backward states: None\n", " forward states: [1 * ndarray(500)]\n", " fw_states_T norm: 1.000000\n", " τ: (2.12e-01:0.50π)\n", " F = 0.000000\n", "Iteration 1\n", " duration: 1.3 secs (started at 2018-12-20 08:17:55)\n", " optimized pulses (ranges): [0.00, 0.29]\n", " backward states: [1 * ndarray(500)]\n", " forward states: [1 * ndarray(500)]\n", " fw_states_T norm: 1.000000\n", " τ: (3.16e-01:0.23π)\n", " F = 0.235157\n", "Iteration 2\n", " duration: 1.3 secs (started at 2018-12-20 08:17:56)\n", " optimized pulses (ranges): [0.00, 0.37]\n", " backward states: [1 * ndarray(500)]\n", " forward states: [1 * ndarray(500)]\n", " fw_states_T norm: 1.000000\n", " τ: (4.91e-01:0.14π)\n", " F = 0.444066\n", "Iteration 3\n", " duration: 1.3 secs (started at 2018-12-20 08:17:57)\n", " optimized pulses (ranges): [-0.07, 0.44]\n", " backward states: [1 * ndarray(500)]\n", " forward states: [1 * ndarray(500)]\n", " fw_states_T norm: 1.000000\n", " τ: (6.44e-01:0.10π)\n", " F = 0.611354\n", "Iteration 4\n", " duration: 1.3 secs (started at 2018-12-20 08:17:59)\n", " optimized pulses (ranges): [-0.14, 0.50]\n", " backward states: [1 * ndarray(500)]\n", " forward states: [1 * ndarray(500)]\n", " fw_states_T norm: 1.000000\n", " τ: (7.60e-01:0.08π)\n", " F = 0.735071\n", "Iteration 5\n", " duration: 1.3 secs (started at 2018-12-20 08:18:00)\n", " optimized pulses (ranges): [-0.20, 0.55]\n", " backward states: [1 * ndarray(500)]\n", " forward states: [1 * ndarray(500)]\n", " fw_states_T norm: 1.000000\n", " τ: (8.41e-01:0.07π)\n", " F = 0.821734\n", "Iteration 6\n", " duration: 1.4 secs (started at 2018-12-20 08:18:01)\n", " optimized pulses (ranges): [-0.26, 0.58]\n", " backward states: [1 * ndarray(500)]\n", " forward states: [1 * ndarray(500)]\n", " fw_states_T norm: 1.000000\n", " τ: (8.96e-01:0.06π)\n", " F = 0.880461\n", "Iteration 7\n", " duration: 1.3 secs (started at 2018-12-20 08:18:03)\n", " optimized pulses (ranges): [-0.30, 0.61]\n", " backward states: [1 * ndarray(500)]\n", " forward states: [1 * ndarray(500)]\n", " fw_states_T norm: 1.000000\n", " τ: (9.31e-01:0.05π)\n", " F = 0.919555\n", "Iteration 8\n", " duration: 1.3 secs (started at 2018-12-20 08:18:04)\n", " optimized pulses (ranges): [-0.34, 0.63]\n", " backward states: [1 * ndarray(500)]\n", " forward states: [1 * ndarray(500)]\n", " fw_states_T norm: 1.000000\n", " τ: (9.54e-01:0.04π)\n", " F = 0.945388\n", "Iteration 9\n", " duration: 1.6 secs (started at 2018-12-20 08:18:05)\n", " optimized pulses (ranges): [-0.36, 0.65]\n", " backward states: [1 * ndarray(500)]\n", " forward states: [1 * ndarray(500)]\n", " fw_states_T norm: 1.000000\n", " τ: (9.69e-01:0.04π)\n", " F = 0.962447\n", "Iteration 10\n", " duration: 1.3 secs (started at 2018-12-20 08:18:07)\n", " optimized pulses (ranges): [-0.39, 0.66]\n", " backward states: [1 * ndarray(500)]\n", " forward states: [1 * ndarray(500)]\n", " fw_states_T norm: 1.000000\n", " τ: (9.78e-01:0.03π)\n", " F = 0.973756\n" ] } ], "source": [ "oct_result = krotov.optimize_pulses(\n", " objectives,\n", " pulse_options=pulse_options,\n", " tlist=tlist,\n", " propagator=krotov.propagators.expm,\n", " chi_constructor=krotov.functionals.chis_re,\n", " info_hook=krotov.info_hooks.chain(\n", " krotov.info_hooks.print_debug_information, print_fidelity\n", " ),\n", " check_convergence=krotov.convergence.check_monotonic_fidelity,\n", " iter_stop=10,\n", ")" ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "attributes": { "classes": [], "id": "", "n": "16" } }, "outputs": [ { "data": { "text/plain": [ "Krotov Optimization Result\n", "--------------------------\n", "- Started at 2018-12-20 08:17:54\n", "- Number of objectives: 1\n", "- Number of iterations: 10\n", "- Reason for termination: Reached 10 iterations\n", "- Ended at 2018-12-20 08:18:08" ] }, "execution_count": 16, "metadata": {}, "output_type": "execute_result" } ], "source": [ "oct_result" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Simulate dynamics of the optimized field\n", "\n", "Having obtained the optimized\n", "control field, we can now\n", "plot it and calculate the\n", "population dynamics under\n", "this field." ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "attributes": { "classes": [], "id": "", "n": "17" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY0AAAEKCAYAAADuEgmxAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3Xd8lfX5//HXlb2BEHYCCXvPiAhaUdGyxI1atW5/X1edbam1rtbWUUetaEGrom0dWK0oCqLFgYgQZG9IgARCyB5k51y/P3K0ERkHknPuk3Ou5+NxHjn3yU3u94GQK/dniqpijDHGeCLE6QDGGGNaDysaxhhjPGZFwxhjjMesaBhjjPGYFQ1jjDEes6JhjDHGY1Y0jDHGeMyKhjHGGI9Z0TDGGOOxMKcDtLSkpCRNTU11OoYxxrQqK1euLFDVDkc7L+CKRmpqKhkZGU7HMMaYVkVEdnlynjVPGWOM8ZgVDWOMMR6zomGMMcZjVjSMMcZ4zIqGMcYYj1nRMMYY4zErGsYYYzwWcPM0jDHBq7qugR35FewsqKS4spbSqjoaXEpEWAgxEaF0aRNNcrtoenaIJTIs1Om4rZIVDWNMq1VT38BX2wv4clsBX+8oZGteOS49+p+LCA1hYNcETkhtx4QBnRjVox1hodbw4gkrGsaYVmfV7mLeyshm/tpcyqrriQwL4YTURM4a1Jm+neLomRRHUlwECdHhhIUItQ0uKqrr2VtaTXZRJev3lLJqdwlzlu7ihS+zSIyNYNqwrlx2Ynf6dIp3+u35NVH1oCy3Iunp6WrLiBgTeFSVjzfmMfuLTFbuKiYmIpSfDurMtOFdOalne6LCj725qaKmni+25jN/XS4fb9hHXYMytld7bj29Dyf1au+Fd+G/RGSlqqYf9TwrGsYYf7d0RwGPfrSZNTmlpCRGc824NKanpxAb2XKNJQUVNczNyOGlr7LIL69hdFoiv508gGEpbVvsGv7MioYxptXbV1rN/fPWs3BDHl3bRHHHmX05b0Q3r/Y/VNc18Pry3cxcvIPCAzVcNCqZX03sT1JcpNeu6Q+saBhjWi1V5R/f7ObRjzZT1+DiF2f04dqT046rCep4lVfX8df/buelJVnER4Xx8HlDmDyki8+u72tWNIwxrVLRgVp+OXcNn27ezyl9kvjDuYPp0T7WsTzb8sq5a+4a1uaUcu7wrjw4bTBtYsIdy+MtnhYNGz1ljPEbS3cUcMebqyk+UMcDZw/kyrGpiIijmfp0iuffN45l5uLtPPvf7azYWcysK0YxuFsbR3M5xdGBySIyUUS2iMh2EZlxmHOmi8hGEdkgIv/ydUZjjPepKi8tyeLyF78hNjKMd28ey1Xj0hwvGN8JDw3h9gl9efvGsagqFzy/lLdX5jgdyxGOFQ0RCQVmApOAgcClIjLwoHP6AL8BxqnqIOB2nwc1xnhVXYOLe95dz0MfbGTCgE68f8vJDOrqn7/FD09py/u3nsyoHu24e+4aHpi3gQZPZhMGECfvNEYD21U1U1VrgTeAcw4653pgpqoWA6jqfh9nNMZ4UVl1HVe+tJzXl+/mpvG9+Nvlo1p0GK03tI+L5NVrRnPdyWm8snQnN/5jJdV1DU7H8hkni0Y3ILvJcY77tab6An1F5CsRWSYiEw/1hUTkBhHJEJGM/Px8L8U1xrSkwooafvbCMpZnFfHERcP41cT+hIT4R3PU0YSFhnDv1IHcf/ZAFm3K47IXv6H4QK3TsXzC3xdbCQP6AOOBS4EXRORHM21UdbaqpqtqeocOHXwc0RhzrHJLq5g+62u25VXwws/TuWBUstORjsvV49KY+bORrNtTyvRZX7O/vNrpSF7nZNHYA6Q0OU52v9ZUDjBPVetUNQvYSmMRMca0UtlFlVz4/NfkldXw6jWjOa1/R6cjNcvkIV2Yc/Vo9pRUcensZewvC+zC4WTRWAH0EZE0EYkALgHmHXTOf2i8y0BEkmhsrsr0ZUhjTMvJLa3ishe/oby6jtevH8OJPQNjfaeTerXnlatHs6+0mktmL2NfaeAWDseKhqrWA7cAC4FNwFuqukFEHhKRae7TFgKFIrIRWAz8UlULnUlsjGmO/eXVXPbCNxQdqOW1a09kSLJ/jpA6XqPTEplzzWj2l9dw6QvLyC+vcTqSV9iMcGOM15VU1jJ91tfkFFfx6jWjSU9NdDqS12TsLOLyv39Dz6Q43vh/Y0iIah2zxz2dEe7vHeHGmFauuq6B6+ZksLOgkhevTA/oggGQnprI3y4fxbb95Vz3SkbADce1omGM8RqXS7njzdVk7CrmqYuHM7ZXktORfGJ8v448OX04K3YVcfM/v6W+weV0pBZjRcMY4zV/mL+Jj9bv494pA5gyNHBXiD2Us4d15aFzBvPp5v38Yf4mp+O0GP+eemmMabVe/DKTl77K4ppxaVx3Sk+n4zjiijE92F14gBe+zCK1fQxXjUtzOlKzWdEwxrS4z7bs5+EPNzFxUGfunTLA6TiOmjFpADsLK3nog410bx/D6f07OR2pWax5yhjTorIKDnDr66vo3zmBJy8e1mqWBvGW0BDhL5cMZ2DXBG791yo27i1zOlKzWNEwxrSY8uo6rn81g7AQYfYVo4iJsMYMgJiIMP5+5QnER4Vzw2sZrXqdKisaxpgW0ThSag1ZBQeYedlIUhJjnI7kVzolRPG3K0axv6yGX7yxqtUuqW5FwxjTIp5dvJ1PNuVx75QBQTO09lgNT2nLQ+cM4sttBTy5aIvTcY6LFQ1jTLMt3VHA059s5dzhXblqbKrTcfzaJaO7c+noFGYu3sGC9fucjnPMrGgYY5olv7yG295YTWpSLA+fN8Rvtmj1Zw9MG8SwlLbcPXcNO/IrnI5zTKxoGGOOW4N7xndZVR3PXTbS73fd8xeRYaE8f9lIIsJCuPmf31JV23qWGrGiYYw5bs8t3s6S7QU8OG0Q/TsnOB2nVenaNponpw9jS145T3zcevo3rGgYY47L1zsKeeqTrZwzvCsXn5By9D9gfmR8v45cnJ7CK0t3klta5XQcj1jRMMYcs5LKWu54czWp7a0fo7luHN+LepcyNyPH6SgesaJhjDkmqsq9/1lPQUUNf7lkBHHWj9EsPdrHcmJaIvPX5jodxSNWNIwxx2Temr18sDaX2yf0Cbjd95xy5sBObMkrJ7uo0ukoR2VFwxjjsb0lVdz7n/WM7N6W/zu1l9NxAsb4fh2Axvku/s6KhjHGIy6Xctdba3C5lKcuHk5YqP34aCk9k+JoEx3Ot7tKnI5yVPavbozxyEtfZfF1ZiH3nT2QHu1jnY4TUEJChBHd2/Lt7mKnoxyVFQ1jzFFt31/BYwu3MGFAJ6an2/BabxjVvR3b9ldQWlXndJQjsqJhjDmiBpfyq7fXEBMRyh/PH2zDa71kZI92AKzO9u8mKisaxpgjmrN0J9/uLuH+swfSMT7K6TgBa1hKW0IEVvl5E5UVDWPMYe0qPMBjCzdzev+OnDu8m9NxAlpcZBip7WPZnFvudJQjsqJhjDkkl0uZ8e91hIeE8PB51izlC307xbM1z4qGMaYVen3Fbr7OLOSeKQPo0iba6ThBoV/neHYWHqC6zn9XvbWiYYz5kb0lVfzpw82M692eS2wxQp/p1zkelzaOVvNXVjSMMT+gqtzz7joaXMoj5w+1Zikf6tspHoAt+/y3icrRoiEiE0Vki4hsF5EZRzjvAhFREUn3ZT5jgtH8dbl8tiWfX/60HymJMU7HCSqp7WOICAvx634Nx4qGiIQCM4FJwEDgUhEZeIjz4oHbgG98m9CY4FNWXceD729kSLc2XGl7fftcWGgIvTvEsdnuNA5pNLBdVTNVtRZ4AzjnEOf9HngUqPZlOGOC0Z8XbqGwooaHzxtMaIg1SzmhV8c4sgoOOB3jsJwsGt2A7CbHOe7XviciI4EUVZ3vy2DGBKPV2SW8tmwXPz8plaHJbZ2OE7TSkmLJKa6ktt7ldJRD8tuOcBEJAZ4E7vLg3BtEJENEMvLz870fzpgAU9/g4p531tExPpK7zurrdJyglpYUg0tht5/ureFk0dgDNB3Ll+x+7TvxwGDgMxHZCYwB5h2qM1xVZ6tquqqmd+jQwYuRjQlMryzdycbcMu4/exDxUeFOxwlqaUlxAH7bROVk0VgB9BGRNBGJAC4B5n33SVUtVdUkVU1V1VRgGTBNVTOciWtMYNpbUsWTi7Yyvl8HJg3u7HScoJfmXnY+q8A/52o4VjRUtR64BVgIbALeUtUNIvKQiExzKpcxwebB9zfgUuX359hSIf6gTUw4ibERZBX4Z/OUozvCq+qHwIcHvXbfYc4d74tMxgSTTzflsXBDHr+aaHMy/ElaUqzdaRhj/Et1XQMPvr+R3h3juO7knk7HMU00Fg3r0zDG+JEXvshkd1ElD5w9iIgw+1HgT1Lbx5BXVkNVrf8tXGjfKcYEoT0lVcz8bDuTBnfm5D5JTscxB/muqTCn2P/6NaxoGBOE/jh/EwC/nTLA4STmUJLbNRaNbCsaxhinfbW9gPnrcrlpfO/vfzgZ/9LdfaeRXVTlcJIfs6JhTBCpa3DxwLwNdE+M4YafWOe3v0qKiyA6PNQvZ4Vb0TAmiMxZupNt+yv43dSBRIWHOh3HHIaIkNwumuzWWDREpK+IfCoi693HQ0XkXu9HM8a0pPzyGv7yyTZO7duBCQM6Oh3HHEVKYgzZxa2zeeoF4DdAHYCqrqVxyQ9jTCvy6ILNVNc3cP/ZA23mdyvQPTGG7KJKVNXpKD/gSdGIUdXlB71W740wxhjvWLmrmLdX5nDtyT3p2SHO6TjGA8ntoqmoqaekss7pKD/gSdEoEJFegAKIyIVArldTGWNaTINLeWDeBjolRHLr6b2djmM89N1cDX8bdutJ0bgZmAX0F5E9wO3AjV5NZYxpMf/+Nod1e0q5Z/IAYiMdXW7OHAN/HXZ71O8gVc0EJohILBCiqv67ea0x5gcO1NTz54VbGNG9LdOGdXU6jjkG391p+Nuw28MWDRG58zCvA6CqT3opkzGmhcz6fAf7y2v42xWjrPO7lYmLDKNdTLjfNU8d6U4j3v2xH3AC/9sg6Wzg4I5xY4yf2VtSxewvMzl7WFdGdm/ndBxzHFLcI6j8yWGLhqo+CCAiXwAjv2uWEpEHgPk+SWeMOW6PL9yCS+HXE/s5HcUcp5TEGDbsKXU6xg940hHeCahtclzrfs0Y46fWZJfw7qo9XHdymq0v1YqltIthT0kVDS7/mavhyVCKV4HlIvKu+/hcYI73IhljmkNV+f0HG0mKi+Cm02yIbWuWkhhNXYOSV1ZN17bRTscBPLjTUNWHgauBYvfjalX9o7eDGWOOz0fr95Gxq5i7zupHnA2xbdW6++EIKk/WnuoOFADvuh+F7teMMX6mpr6BP320if6d45menuJ0HNNMKd/tq+FHRcOTX0Pm454NDkQDacAWYJC3Qhljjs8rX+0ku6iKf1x7IqEhNsS2tevaNhoR/GrhQk8m9w1peiwiI4GbvJbIGHNcCitqePa/2zmjf0fbwjVARISF0CUhihw/utM45v00VPVb4EQvZDHGNMNTn2ylsq6B30y2LVwDSXJijF9N8DvqncZBM8NDgJHAXq8lMsYcs6155fzrm91cMaYHvTvaKraBJKVdDF9tL3A6xvc86dOIb/K8nsY+jn97J44x5ng8PH8TcZFh3D6hr9NRTAtLSYwmr7yamvoGIsOc323Rk6KxUVXnNn1BRC4C5h7mfGOMD322ZT+fb83n3ikDaBcb4XQc08JS2sWgCnuKq/xiLxRP+jR+4+Frxhgfq29w8fD8TaS2j+HnJ6U6Hcd4wf/21fCPEVRHWuV2EjAZ6CYizzT5VAK2c58xfuGNFdls21/B3y4fRUTYMY9rMa1ASmLjTHB/matxpOapvUAGMA1Y2eT1cuAOb4YyxhxdWXUdTy3ayolpifx0kC0HF6g6xUcRERriNyOojrTK7RpgjYj8U1W9cmchIhOBvwChwIuq+shBn78TuI7GO5t84BpV3eWNLMa0NjMXb6eospbfTR1oe2UEsJAQoVu7aHL8ZAe/IzVPvaWq04FVIvKjJRZVdWhzLiwiocBM4EwgB1ghIvNUdWOT01YB6apaKSI3Ao8BFzfnusYEgt2Flby8ZCcXjExmcLc2TscxXpbcLtr/7zSA29wfp3rp2qOB7e7tZBGRN4BzgO+LhqoubnL+MuByL2UxplV5dMFmQkOEu8+yvTKCQffEGD5cl+t0DODIzVO57o/eag7qBmQ3Oc7hyDPNrwU+8lIWY1qNjJ1FzF+Xy+0T+tC5TZTTcYwPpCTGUFxZR3l1HfFR4Y5mOVLzVDn/W6gQQNzHAqiqJng5W9MslwPpwKmH+fwNwA0A3bvbArwmcLlcjXtldE6I4oaf9HQ6jvGR/612W8XArs4WjcOO0VPVeFVNaPKIb/qxBa69B2i6dnOy+7UfEJEJwG+Baapac5iss1U1XVXTO3To0ALRjPFP89bsZU1OKb/8aT9iImyvjGDx/bBbP+jX8Oi7zr2y7ck03mksUdVVLXDtFUAfEUmjsVhcAvzsoOuOAGYBE1V1fwtc05hWq6q2gUcXbGZItzacN6Kb03GMD/nTvhqebMJ0H43bu7YHkoBXROTe5l7YPYz3FmAhsAl4S1U3iMhDIjLNfdrjQBwwV0RWi8i85l7XmNbqxS8zyS2t5ndTBxJie2UElbYx4cRFhpHjB7PCPbnTuAwYpqrVACLyCLAa+ENzL66qHwIfHvTafU2eT2juNYwJBHll1Tz/+Q4mDe7M6LREp+MYHxORxmG3reFOg8aZ4U2HaERyiL4HY4z3PPHxFuoblBmT+jsdxTgkxU/21fCkaJQCG0TkFRF5GVgPlIjIMwetSWWM8YL1e0qZuzKHq8al0qN9rNNxjENS28ewq7ASl+tHc619ypPmqXfdj+985p0oxpiDqSoPz99Eu5gIbj6tt9NxjIPSkuKoqXeRW1ZNt7bRjuXwZI/wOb4IYoz5sUUb8/g6s5DfnzOINtHOjs83zkpLarzLzMo/4GjR8GT01FQRWSUiRSJSJiLlIlLmi3DGBLPaehd/+mgzvTvGcelom7Qa7L4vGoUHHM3hSfPU08D5wDpVdbYxzZgg8tqyXWQVHODlq04gLNT2ygh2nRIiiQ4PJSvf2aLhyXdiNrDeCoYxvlNSWcszn27jlD5JjO9nqxyYxmG3aUmxZBVUOJrDkzuNXwEfisjnwPfLeKjqk15LZUyQe/qTbZRX13HvFNsrw/xPWlIsG3Od7R3w5E7jYaCSxrka8U0exhgv2JFfwT+W7eKS0d3p19n+q5n/SUuKZXdRJXUNLscyeHKn0VVVB3s9iTEGgD/O30R0eCh3ntnX6SjGz6QlxdLgUrKLKunZIc6RDJ7caXwoImd5PYkxhiXbCvh0835uPr03SXGRTscxfibVPYJqp4MjqDwpGjcCC0SkyobcGuM9DS7lD/M3kpIYzdXjUp2OY/xQT3fRyHRwBJUnk/usUdUYH3hzRTab95Xz3GUjiQwLdTqO8UPtYiNoEx1OVoEfFw0AEWkH9KHJwoWq+oW3QhkTbMqr63hy0RZGpyYyaXBnp+MYP9arQyzb9zs37PaoRUNErgNuo3FnvdXAGOBr4HTvRjMmeMxcvIPCA7W8fJUNsTVH1q9zAh+tz0VVHfle8aRP4zbgBGCXqp4GjABKvJrKmCCyu7CSl5Zkcf6IZIYkt3E6jvFzA7rEU1JZR17ZIXe/9jpPikZ1kw2YIlV1M9DPu7GMCR6PLNhEaIjwq4n238ocXb9Ojd3Mm/c5Mx7Jk6KRIyJtgf8Ai0TkPWCXd2MZExyWZxXx4bp9/N+pveiUEHX0P2CCXv/OCQBs3lfuyPU9GT11nvvpAyKyGGgDLPBqKmOCgMul/P6DjXRpE8UNP+npdBzTSrSJCadLmyi2+GvRaEpVP/dWEGOCzbur9rBuTylPXzyc6AgbYms8169zPJscWoPK1ls2xgGVtfU8tnAzw1LaMm1YV6fjmFZmcNc2bN9fQVVtg8+vbUXDGAfM+jyTvLIa7ps6gJAQG2Jrjs2I7m2pdylrc3w/kNWjoiEiPURkgvt5tIjYLHFjjlNuaRWzvtjB1KFdGNUj0ek4phUa0b0dAN/u9sOiISLXA28Ds9wvJdM4ksoYcxweW7AFl8KMSf2djmJaqcTYCNKSYvl2d7HPr+3JncbNwDigDEBVtwEdvRnKmEC1clcx767aw/WnpJHcLsbpOKYVG9G9Lat2F+PrTVU9KRo1qlr73YGIhAG29asxx8jlUh6Yt4HOCVHcfFpvp+OYVm5Uj3YUVNSS6ePFCz0pGp+LyD1AtIicCcwF3vduLGMCz9yV2azbU8pvJvcnJuKYRrsb8yOn9m3cO/7TTXk+va4nRWMGkA+sA/4f8CFwrzdDGRNoSqvqeGzBFk5IbWdDbE2LSG4Xw4AuCSza6GdFQ1VdqvqCql4E3AB8o75uRDOmlfvLJ9sorqzlgWmDbBVb02LOHNiJlbuKKajw3eKFnoye+kxEEkQkEVgJvCAiT7XExUVkoohsEZHtIjLjEJ+PFJE33Z//RkRSW+K6xvjStrxyXv16J5eM7s6grraKrWk5Px3UCZfC/LW5PrumJ81TbVS1DDgfeFVVTwTOaO6FRSQUmAlMAgYCl4rIwINOuxYoVtXewFPAo829rjG+pKo8+P5GYiJCufssW8XWtKyBXRIYltKWl77Korbe5ZNrelI0wkSkCzAd+KAFrz0a2K6qme7RWW8A5xx0zjnAHPfzt4EzxEv39g0uJTO/gtzSKsqq67xxCROEPt6Yx5LtBdx1Vj8SYyOcjmMCjIhw+4Q+7Cqs5KEPNlDug59dngzheAhYCCxR1RUi0hPY1gLX7gZkNznOAU483DmqWi8ipUB7oKAFrv8DJZW1nP7E/9Zj7NomirOHdeXaU9LoGG9LVptjV13XwO8/2Ei/TvFcdmJ3p+OYAHVav45cPS6Vl7/ayfKsIj6+41SvXs+TpdHn0jjM9rvjTOACb4Y6ViJyA42d9HTvfnz/OWMjw3j64uFU1TVQWlXHyl3FvPBlJm+vzOHZn43kpF7tWzKyCQIvfJFJTnEV/7ruRMJCbZk34z33TR3I1KFdKKl08E5DRP7KESbxqeovmnntPUBKk+Nk92uHOifHPamwDVB4iCyzgdkA6enpxzWyKyo8lHNHdPvBa9vyyvm/f6zkypeX8/r1YxjVo93xfGkThPaWVPHcZzuYPKQzY3snOR3HBDgR8dk6Zkf69SeDxtFSh3s01wqgj4ikiUgEcAkw76Bz5gFXup9fCPzXl8N9+3SKZ+7/jaVLmyiufzWD/WXVvrq0aeUenr8Jlyr3TB7gdBRjWtRh7zRUdc7hPtcS3H0Ut9DYXxIKvKSqG0TkISBDVecBfwdeE5HtQBGNhcWnEmMj+PuV6Ux+ZgkPvL+B5y4b5esIppX5fGs+89flcvdZfW19KRNwjtqn4d7i9Ue/3avq6c29uKp+SOMM86av3dfkeTVwUXOv01y9O8Zz2xl9eHzhFhZv3s9p/W29RnNo1XUN3P/eenomxXK9beFqApAno6fubvI8isZO8HrvxPFfN/ykJ2+uyObxhVs4tW8H2zjHHNKszzPZWVjJa9eOJjLMtnA1gceTZURWNnl8pap3AuO9H82/hIeGcPuEPmzMLWPBhn1OxzF+aFfhAWZ+tp2pQ7twSp8OTscxxis8WUYksckjSUR+SuMopqBzzvBupCXFMvuLTKejGD+j2rjseURoCL+bevDCBsYEDk8Gj6/kfyOpvgbuonF5j6ATGiJceVIPVmeXsDrb99ssGv+1cEMei7fkc/uEPnRKsMmgJnB50jyVpqo93R/7qOpZqrrEF+H80QWjkomLDGPO0p1ORzF+4kBNPQ+9v4H+neO5amyq03GM8SpPmqeiROROEXlHRP4tIreLSND+KhUfFc6Fo5L5YO1e8st9txyx8V/P/Hcbe0ur+cO5g23mtwl4nnyHvwoMAv4KPOt+/po3Q/m7n5/Ug7oG5a2M7KOfbALatrxy/v5lFtPTk0lP9c2MXGOc5MmQ28Gq2rRnb7GIbPRWoNagZ4c4xvRM5K2MbG48tZcNvw1SLpfy23fXExsZxq8n9nc6jjE+4cmdxrciMua7AxE5kcaO8aB2yQnd2VVYybKsHy2FZYLE6yt2s3xnEb+dMoD2cZFOxzHGJzwpGqOApSKyU0R20jiC6gQRWScia72azo9NHNyZ+Kgw3lphTVTBaF9pNY98uJmxvdpz0ahkp+MY4zOeNE9N9HqKVigqPJTzRnTjjRXZPFhZR5uYcKcjGR+6f956ahtc/PG8Ibbntwkqngy53XWkhy9C+qvp6SnU1rt4b83BK7qbQLZg/T4Wbsjj9gl9SU2KdTqOMT5l4wObYXC3NgzulsDry7Px4YrtxkGlVXXc9956BnRJ4LpT0pyOY4zPWdFopovTU9iUW8b6PWVORzE+8OiCzRRU1PDoBUMItzkZJgjZd30zTRvejciwEN7M2O10FONly7OK+Nc3u7lmXBpDk9s6HccYR1jRaKY20eFMHtKF91btpaq2wek4xkuq6xr4zTtrSW4XzZ1n9XU6jjGOsaLRAi4+IYXymno+Wp/rdBTjJU99spUd+Qd4+LwhxER4MujQmMBkRaMFnJiWSGr7GN6wORsB6dvdxbzwRSaXjk7h1L62T4YJblY0WoCIMP2EFJZnFZGZX+F0HNOCqusauHvuGrq0ieaeyQOcjmOM46xotJALRyYTGiK8lZHjdBTTgp74eAuZ+Qd49IKhxEfZBE5jrGi0kI4JUZzWryNzM7KprA26LdQDUsbOIl5cksVlJ3bn5D5JTscxxi9Y0WhBN47vSeGBWuYsDeqJ8gGhqraxWapb22h+Y81SxnzPikYLGtUjkfH9OjDrix2UVNY6Hcc0w2MLN7OzsJLHLhxKXKSNljLmO1Y0WtivJ/anorqeB+ZtcDqKOU5LthXw8lc7ufKkHoztZc1SxjRlRaOFDeiSwC2n9+Y/q/fy8YZ9Tscxx6j4QC13vrWaXh1imTHJmqWMOZgVDS+4aXxvBnRJ4J5311NYYfuItxaqyox31lJN+U2NAAAO8klEQVRcWctfLhlBdESo05GM8TtWNLwgIiyEpy4eRll1Hb96e62tgNtKvJWRzcINefzyp/0Y3K2N03GM8UtWNLykf+cEZkzsz6eb9/OPb2wxQ3+XmV/BA/M2MrZXe647uafTcYzxW1Y0vOjqcamc2rcDf/hgI9v3lzsdxxxGXYOL299cTURYCE9MH0ZIiO3EZ8zhOFI0RCRRRBaJyDb3x3aHOGe4iHwtIhtEZK2IXOxE1uYQER6/aCixkWHc+vpqauptFVx/9NiCzazNKeVP5w+hS5top+MY49ecutOYAXyqqn2AT93HB6sEfq6qg2jcp/xpEWl1mxh0jI/i8QuHsim3jD8v3OJ0HHOQRRvzeOHLLC4f053JQ7o4HccYv+dU0TgHmON+Pgc49+ATVHWrqm5zP98L7Ada5RKjZwzoxBVjevDCl1ks3rLf6TjGLae4krvnrmFQ1wTunTLQ6TjGtApOFY1Oqvrd5hP7gE5HOllERgMRwI7DfP4GEckQkYz8/PyWTdpCfjtlAP07x3PHm6vZU1LldJygV1vv4pZ/rcLlUp67bCRR4Ta81hhPeK1oiMgnIrL+EI9zmp6njeNRDzsmVUS6AK8BV6uq61DnqOpsVU1X1fQOHfzzZiQqPJTnLx9FfYNy8z+/pbb+kG/F+MijCzazOruERy8cSo/2sU7HMabV8FrRUNUJqjr4EI/3gDx3MfiuKByyzUZEEoD5wG9VdZm3svpKWlIsj104lNXZJfzpo01OxwlaC9bv4+9LsrjypB7Wj2HMMXKqeWoecKX7+ZXAewefICIRwLvAq6r6tg+zedXkIV24elwqL3+1k/lrbXtYX9uWV85db61mWHIb7pliy4QYc6ycKhqPAGeKyDZggvsYEUkXkRfd50wHfgJcJSKr3Y/hzsRtWb+ZNIAR3dvy63+vtZ3+fKi0qo7rX80gOiKMv10xisgw68cw5lhJoC1xkZ6erhkZGU7HOKo9JVVMfeZL2sdF8u5NY21XOC9rcCnXvLKCpTsKeP36MaSnJjodyRi/IiIrVTX9aOfZjHCHdGsbzczLRpJVcIDb31hNgyuwire/+fPHW/h8az4PTBtkBcOYZrCi4aCxvZK4/+yBfLp5P098bBP/vOWDtXt5/rMd/OzE7lx2Yg+n4xjTqtmWZA67YkwPNuWW8dxnO+jfJYFpw7o6HSmgfLu7mLveWkN6j3Y8cPYgp+MY0+rZnYbDRIQHpw3mhNR2/HLuGlbtLnY6UsDYXVjJ9XMy6JQQxawrRhERZt/uxjSX/S/yAxFhITx/+Sg6JURx3ZwMdhUecDpSq1dSWctVryynQZVXrj6B9nGRTkcyJiBY0fATSXGRvHL1CbhUuerlFRQdqHU6UqtVU9/ADa+tJKeoitlXpNOzQ5zTkYwJGFY0/EjPDnG8eGU6e0uquG7OCqrrbCn1Y9XgUu54czXLs4p4/KKhjE6zkVLGtCQrGn5mVI9E/nLJcFZll3Dr66uob7A1qjylqtzzzjo+XLePe6cM4Jzh3ZyOZEzAsaLhhyYO7sKD0waxaGMed89dY3M4PKCqPDx/E29mZPOL03tz3Sm2Zasx3mBDbv3Uz09KpaKmnscWbCE6IpQ/njcEEduG9HCe/e92XlySxVVjU7njzL5OxzEmYFnR8GM3je9NZU0Dzy7eTnR4GL+bOsAKxyHMXLydJxZt5fwR3bhv6kD7OzLGi6xo+Lm7zurLgdp6Xvoqi/BQYcak/vZDsYlnPt3Gk4u2cs7wrjx24VBCQuzvxhhvsqLh50SE+6YOpMGlzPoik6q6Bh44e1DQ/3BUVZ76ZBvPfLqN80d24/ELhxEa5H8nxviCFY1WoHHW+CCiw0OZ9UUmlbUNPHrB0KD9IamqPLJgM7M+z2R6ejJ/Oj94/y6M8TUrGq2ESGPTVExEGE99spXqugaemD4s6PaEqK13MeOdtbzz7R4uH9Odh6YNDvq7LmN8yYpGKyIi3DahD9ERIfzxw83kl9cw64pRtI2JcDqaT1TU1HPjP1by5bYC7jqzL7ec3tv6d4zxMZun0Qrd8JNejRMAd5dw/nNLg2Ktqryyai6e9TVLdxTy2IVDufWMPlYwjHGAFY1W6pzh3fjn9SdSXFnLec8tJWNnkdORvGblriKm/nUJWQUHePHKdKanpzgdyZigZUWjFTshNZF3bhpHm+hwLpm9jJe/yiLQtu/91ze7uWT2MmIiQnn3pnGc1q+j05GMCWpWNFq5tKRY/nPTOMb368iD72/k1tdXUVFT73SsZquoqefOt1Zzz7vrGNsriXk3n0y/zvFOxzIm6FlHeABoExPO7CtGMeuLTB5fuJkNe8t4YvowRnZv53S047Imu4RfvLGK7KJKfnFGH247o48NqTXGT9idRoAICRFuHN+Lf10/htp6Fxc+v5Q/L9xCbX3rWSW3pr6Bpz/ZygXPL6Wu3sUbN5zEnWf2tYJhjB+xohFgxvRsz0e3n8IFI5N5dvF2pj27hOVZ/t9JnrGziCnPLOHpT7YxZWgXPrrtJ7YXhjF+SAKt4zQ9PV0zMjKcjuEXFm3M4/731rO3tJrzR3RjxuT+dIyPcjrWD+wrreaJj7cwd2UO3dpG84dzB3Naf+vsNsbXRGSlqqYf7Tzr0whgZw7sxLje7Zm5eDuzv8hk4YZ9XDUuletP6en4hMDy6jpmfZ7Ji0sycbng+lPSuH1CX2Ij7VvSGH9mdxpBIjO/gicWbWX+2lziI8O4+uQ0Lh/T3ed3HvvLqnnpq538c9kuymvqmTasK7/8aT9SEmN8msMY80Oe3mlY0Qgym3LLeGrRVj7emEd4qDB5SBeuGNODUT3aeW2GtculLMss5O2VOXywNpd6l4tJQ7pw46m9GNytjVeuaYw5NlY0zBFl5lfw2rJdvJ2RQ3lNPd3aRjNpcGcmDenM0OS2hIc2b4xEbb2LjF1FfL4lnw/W5rKnpIr4yDDOHdGNa09OIzUptoXeiTGmJVjRMB6pqKln4fp9fLguly+25VPXoMREhDKqRztG9WhH307x9OwQS2r7WKLCf7yirqpSXlNPXmk12/dXsH5vKev3lJGxs4gDtQ2Ehwon9UriwlHJnDWw0yG/hjHGeX5dNEQkEXgTSAV2AtNVtfgw5yYAG4H/qOotR/vaVjSOX2lVHUu2FbA8q5BvsorYvK/8B5+PCg8hPiqcuMgwGlxKfYOL4so6quoavj8nLETo3TGOkT3aMb5vB8b2TiLOOreN8Xv+PnpqBvCpqj4iIjPcx78+zLm/B77wWbIg1iY6nClDuzBlaBcAKmvrycw/QGbBAXYXHqCsup6yqjoqauoJCxHCQkNoGx1Op4QoOiZEkto+ln6d4+1uwpgA5lTROAcY734+B/iMQxQNERkFdAIWAEetgKZlxUSEMbhbG+usNsZ8z6kZ4Z1UNdf9fB+NheEHRCQEeAK4+2hfTERuEJEMEcnIz89v2aTGGGO+57U7DRH5BOh8iE/9tumBqqqIHKpj5SbgQ1XNOdpQUFWdDcyGxj6N40tsjDHmaLxWNFR1wuE+JyJ5ItJFVXNFpAuw/xCnnQScIiI3AXFAhIhUqOoML0U2xhhzFE71acwDrgQecX987+ATVPWy756LyFVAuhUMY4xxllN9Go8AZ4rINmCC+xgRSReRFx3KZIwx5ihscp8xxhiP52nYfhrGGGM8ZkXDGGOMxwKueUpE8oFdzfgSSUBBC8VpLYLtPQfb+wV7z8GiOe+5h6p2ONpJAVc0mktEMjxp1wskwfaeg+39gr3nYOGL92zNU8YYYzxmRcMYY4zHrGj82GynAzgg2N5zsL1fsPccLLz+nq1PwxhjjMfsTsMYY4zHrGi4ichEEdkiItvdG0MFNBF5SUT2i8h6p7P4ioikiMhiEdkoIhtE5DanM3mbiESJyHIRWeN+zw86nckXRCRURFaJyAdOZ/EVEdkpIutEZLWIeG1ZDGueovEbDNgKnAnkACuAS1V1o6PBvEhEfgJUAK+q6mCn8/iCe0XlLqr6rYjEAyuBcwP831mAWFWtEJFwYAlwm6oucziaV4nInTRu3JagqlOdzuMLIrKTxoVdvTo3xe40Go0GtqtqpqrWAm/QuLtgwFLVL4Aip3P4kqrmquq37uflwCagm7OpvEsbVbgPw92PgP5NUUSSgSmALX7qBVY0GnUDspsc5xDgP0yCnYikAiOAb5xN4n3upprVNO5bs0hVA/09Pw38CnA5HcTHFPhYRFaKyA3euogVDRN0RCQO+Ddwu6qWOZ3H21S1QVWHA8nAaBEJ2OZIEZkK7FfVlU5nccDJqjoSmATc7G6CbnFWNBrtAVKaHCe7XzMBxt2u/2/gn6r6jtN5fElVS4DFwESns3jROGCau33/DeB0EfmHs5F8Q1X3uD/uB96lsdm9xVnRaLQC6CMiaSISAVxC4+6CJoC4O4X/DmxS1SedzuMLItJBRNq6n0fTONhjs7OpvEdVf6OqyaqaSuP/4/+q6uUOx/I6EYl1D+5ARGKBswCvjIy0ogGoaj1wC7CQxs7Rt1R1g7OpvEtEXge+BvqJSI6IXOt0Jh8YB1xB42+fq92PyU6H8rIuwGIRWUvjL0eLVDVohqEGkU7AEhFZAywH5qvqAm9cyIbcGmOM8ZjdaRhjjPGYFQ1jjDEes6JhjDHGY1Y0jDHGeMyKhjHGGI9Z0TCmmUSkrYjc5H7eVUTedjqTMd5iQ26NaSb3OlYfBMtqwSa4hTkdwJgA8AjQy70o4DZggKoOFpGrgHOBWKAP8GcggsYJhjXAZFUtEpFewEygA1AJXK+qATtr27Ru1jxlTPPNAHa4FwX85UGfGwycD5wAPAxUquoIGmfj/9x9zmzgVlUdBdwNPOeT1MYcB7vTMMa7Frv37igXkVLgfffr64Ch7hV3xwJzG5fGAiDS9zGN8YwVDWO8q6bJc1eTYxeN//9CgBL3XYoxfs+ap4xpvnIg/nj+oHs/jywRuQgaV+IVkWEtGc6YlmRFw5hmUtVC4CsRWQ88fhxf4jLgWvcKpRsI8K2GTetmQ26NMcZ4zO40jDHGeMyKhjHGGI9Z0TDGGOMxKxrGGGM8ZkXDGGOMx6xoGGOM8ZgVDWOMMR6zomGMMcZj/x+/dehx/PdfHQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plot_pulse(oct_result.optimized_controls[0], tlist)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In contrast to the dynamics under the guess field, the optimized field indeed\n", "drives the initial state $\\ket{\\Psi_{\\init}} = \\ket{0}$ to the desired target\n", "state $\\ket{\\Psi_{\\tgt}} = \\ket{1}$." ] }, { "cell_type": "code", "execution_count": 18, "metadata": { "attributes": { "classes": [], "id": "", "n": "18" } }, "outputs": [], "source": [ "opt_dynamics = oct_result.optimized_objectives[0].mesolve(\n", " tlist, e_ops=[proj0, proj1])" ] }, { "cell_type": "code", "execution_count": 19, "metadata": { "attributes": { "classes": [], "id": "", "n": "19" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYUAAAEKCAYAAAD9xUlFAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4yLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvOIA7rQAAIABJREFUeJzt3Xd4VGX2wPHvyaQXSkgIJQkl9CoSuhQ7IqLYUdG1odh1f2tZdde+lrVgA7GsiiuIHVFAQKQoLRFpoYWEEjoJLYSQZOb9/XEn2YiQDCGTO5M5n+e5T+7cuXPnTJQ5OfdtYoxBKaWUAgiyOwCllFK+Q5OCUkqpMpoUlFJKldGkoJRSqowmBaWUUmU0KSillCqjSUEppVQZTQpKKaXKaFJQSilVJtjuAE5WXFycad68ud1hKKWUX0lPT99rjImv7Dy/SwrNmzcnLS3N7jCUUsqviMhmT87T20dKKaXKaFJQSilVRpOCUkqpMpoUlFJKldGkoJRSqozXkoKIfCAiu0Vk1QmeFxF5XUQyRWSFiJzurViUUkp5xpuVwofA4AqevwBo7d5GAWO9GItSSikPeG2cgjFmnog0r+CUi4GPjbUe6CIRqScijY0xO7wRz9JNecxbvwcRQYAgEUSw9oMEABH38XL71nEhyH1u6T7u80IdQUSEOogMdRAR6iAqNLhsv35kKJGhDsR9HaWU8nV2Dl5rCmwt9zjHfexPSUFERmFVEyQnJ1fpzX7bvI8352RS00tSh4cE0SAqjAbRoTSICiU+Jozk2EiSYiNJdm+xUaGaOJRSPsEvRjQbY8YD4wFSU1Or9LV+28AUbhuYUno9jAGXMRgo26fcvnGf5zKAAYO1b9zPuYwBA0VOF0eKnBwuclJQVMKRIicF7v39BcXkHi5ib/5RcvOL2JtfxMptB9mbf/QPsdUJD6Zd4zp0aFyH9o1j6NC4Lu0bxxDs0H4ASqmaZWdS2AYklXuc6D7mdeK+dRSEPX+dHylysnVfAVtyC9icV8DGPfms3XGQyWlbKShyAhAZ6qBbcj16NI+lZ/NYTm9Wn/AQhy3xKqUCh51JYQpwl4hMAnoBB7zVnuBrIkIdtEmIoU1CzB+Ou1yGLXkFrNx2gLRNeSzZtI8xszdgDESEOOjXKo6z2jXkzHbxNK4bYVP0SqnazGtJQUQmAoOAOBHJAf4JhAAYY8YBPwBDgEygALjRW7H4i6AgoXlcFM3jorioaxMADhYWk7Ypj5/X7eGntbuZtWYXAF0S6zKsaxMu6tqEhDrhdoatlKpFxNR0y+spSk1NNYE6S6oxhszd+cxeu5vvV+xg5bYDiECflg24pFtTLurShIhQvcWklPozEUk3xqRWep4mBf+1cU8+U37fzre/b2NTbgF1woO5vHsS1/VOpmV8tN3hKaV8iCaFAGKMYUl2HhMWbWb6qp2UuAxntIrj1gEtGdA6Tru7KuXLjAGXE6uboyn3s5zSf8NBwRBUtbsBniYFv+iSqiomIvRq2YBeLRuw+1Ahk5duZcKizdzwwRI6Na3D6IGtGNypEY4gTQ5K2eLIfti6GLalQ+5GyMuCw3vh6AE4egiMy7PrXPgK9LjZq6FqpVBLHS1x8s2ybYybm0X23sO0iIvirjNbcUm3ppoclKoJ+7fA6q9h9TewfRlgQIKgbhLEtoDoRhBeB8LqQHAY4J42AXFXBqX/Tst9R6ecDU1Oq1I4evtIAeB0GWas3smbP2WSseMgbRKi+et5bTmvQ4LeVlKquhkDWT/DorGwYYZ1rEk3aHMBNOsLTbtDaKQtoWlSUH/gchmmrdrJyzPXkbXnMF2T6vHw4Hb0SWlgd2hK1Q5Zc2HWP62qICoeUm+GrldbVYEP0KSgjqvE6eKr37bx2qz1bD9QyOCOjXj0wvYkxdrz14tSfi8vG6Y9CBt+tG4NDXwIulzpviXkOzQpqAoVFjt5b34Wb83ZiNMYRvVvyehBKUSFad8DpTzicsLid+Cnp0EcMPBB6DkKQnxzMKkmBeWRnQcKeWH6Wr5eto2EOmH8fUh7hnVtou0NSlUkfzd8eTNkz4PW58HQV6Fuot1RVcjTpKDTcAa4RnXDefWq0/hydF8S6oRz76Tfuf6DJWzOPWx3aEr5ps0LYVx/2LoEhr0B10z2+YRwMjQpKAC6N6vP13f048lhHVm2ZT/nvTqPt3/OpNjpYf9ppQLB8knw0VCrB9Ets+D06/83sKyW0KSgyjiChBv6NmfmAwM4s21DXpy+jqGvLyB9c57doSllL2Ng/svw9W2Q3AdunQONOtsdlVdoUlB/0rhuBONGdufd61M5WFjM5eMW8uR3qzniXutBqYBiDPz4GMx+CjpfCdd9BRH17I7KazQpqBM6t0MCMx8YyMjezfjPL5sYPGYeS7K1alABxBiY8SgsfBN63gbD34HgULuj8ipNCqpC0WHBPHVxJz69tRcuY7hqvFU1FBSV2B2aUt5VWiEsegt63Q4XvABBtf8rs/Z/QlUt+qbEMf3eAVzvrhouGDNfqwZVu/3ymrtCGAWDn691DconoklBeSwqLJgnL+7ExFt7l1UNT0zRqkHVQssnwawnoNPlMPiFgEkIoElBVUGflAZlVcOHv2rVoGqZ7Hnw7Z3QYgBc8nZA3DIqL7A+rao2x6sanvouQ3soKf+2bzNMvgEatIKrPvG5+YtqgiYFdUpKq4aRvZvxwS/ZDHl9PmmbtGpQfqioAD671prT6OpPIbyu3RHZQpOCOmVRpT2UbulFsdPFFe8s5JmpGRQWa9Wg/IQxMOVu2LkKLn8fGqTYHZFtNCmoatO3VRzT7xvANT2TeW9BNkPGzNfR0Mo/LBkPq76As/8Brc+1OxpbaVJQ1So6LJhnh3fmv7f04miJi8vHLeS5H9Zo1aB8186V1niENoPhjPvtjsZ2mhSUV/RrFcf0+/pzdY9kxs/LYsjr8/ltyz67w1Lqj4oK4IubICIWLn47oLqenogmBeU1MeEh/OvSznx8U08Ki5xcPvZX/jVNqwblQ2Y8Ans3wPBxEKVL04ImBVUDBrSJZ/r9A7gyNYl35mYx9I0F/L51v91hqUC3Ziqkfwj97oGUM+2OxmdoUlA1ok54CM9f1oWPburJ4aMlXPr2LzwzNYPDR3U0tLJBQR5Mvd+a/vrMx+yOxqdoUlA1amCbeGbcP4CreiTx3oJsznllLtNX7cTfloVVfm76w3Akz2pHqOWznp4sTQqqxtUJD+Ffl3bhy9F9qBsRwu2fpHPzR2lszSuwOzQVCNZNhxWfQf+/QuMudkfjczQpKNt0bxbL1LvP4LEL27MoK5dzX53LW3MytSFaec+R/TD1PmjYEfr/n93R+CSvJgURGSwi60QkU0QePs7zySIyR0SWicgKERnizXiU7wl2BHFL/5bMemAgg9o05KUZ6zj75blMWb5dbymp6jfrCcjfDZe8pbeNTsBrSUFEHMBbwAVAB2CEiHQ45rTHgMnGmG7A1cDb3opH+bYm9awlQD+9tRd1I0K4Z+Iyhr/9q46IVtUnJ93qbdTrdmjSze5ofJY3K4WeQKYxJssYUwRMAi4+5hwD1HHv1wW2ezEe5Qf6psTx3d1n8NLlXdi+/wiXjV3I6E/SWbfzkN2hKX/mcsL390N0Agz6000LVU6wF6/dFNha7nEO0OuYc54AfhSRu4Eo4BwvxqP8hCNIuCI1iQu7NOaduVm8vyCb6at3MqRzY+47uzWtE2LsDlH5m7QPYMdyuPwDCK9T+fkBzO6G5hHAh8aYRGAIMEFE/hSTiIwSkTQRSduzZ0+NB6nsERkazP3ntmH+g2cyemAKP6/dzXmvzePuictYte2A3eEpf5G/G2Y/DS0HQcdL7Y7G53kzKWwDkso9TnQfK+9mYDKAMWYhEA7EHXshY8x4Y0yqMSY1Pj7eS+EqX1U/KpQHB7dj/kNncfvAFH5as4uhbyzgyncWMmP1TpwubZBWFfjxcSgugCH/1rmNPODNpLAUaC0iLUQkFKshecox52wBzgYQkfZYSUFLAXVcsVGhPDS4HQv/fjaPDmnPtn1HuG1COme9/DPj5m5k98FCu0NUvmbLYlgxCfreDXGt7Y7GL4g3u/25u5i+BjiAD4wxz4rIU0CaMWaKuzfSu0A0VqPzg8aYHyu6ZmpqqklLS/NazMp/lDhdzFi9iw9/zWbppn04goRBbeK5vHsig9o2JCLUYXeIyk4uF7x/DhzYBnenQ1i03RHZSkTSjTGplZ7nb33BNSmo49m4J58v0nP4Mj2H3YeOEhHiYFDbeAZ3asSZ7RpSJzzE7hBVTVv5BXx5szWVRbdr7Y7GdpoUVEAqcbpYlJXH9NU7mLF6F3sOHSVIoGtSPc5oFUfflDhOS6qnVURtV3wE3uwBEfVh1FwIsrtPjf00KaiA53QZftuyj3nr9/BL5l6W5xzA6TI4goTWDaPpmliPTol1SYmPomVcNAl1whBtiKwd5r8Ms5+CG6ZCi/52R+MTNCkodYyDhcUszc5j+db9/J5zgBU5+9lfUFz2fGSog+YNomhUN5z46DDiY6wtJjyYyFAHEaHunyEOHEF/Th7GgMsY92YlJZcxuFwGpzG4XNbzztJjLus8YwxhIUFEhQYTHR5MTFgI9aJCiAkL1iRVFYd2wRunW11Qr/6v3dH4DE+TgjcHrynlU+qEh3B2+wTObp8AWF/GOw4Ukr33MFl7D5O1J5/NuQXsPlTI6u0H2JtfZGt316hQB43rRdC4bjgp8dG0SYihbaMYOjapQ3iI3v46oTnPQkkhnPuU3ZH4JU0KKmCJCE3qRdCkXgT9Wv1peAwulyGvoIj8whIKipwcKbZ+FhQ5TzhZn4jgEMERJIhYo7MdItbxIMER9MdzgkQICoKjxS7yj5ZYW2EJeYeL2H7gCDsPFLJ9/xE+T9vK4SJr9thQRxBdk+rSq0UDBraNp3tyfYKOU7kEpL0bYNkE6DkKGqTYHY1f0qSg1AkEBQlx0WHERYfZHQoul2Hb/iNk7DhI+uZ9LM7KZezcjbw5J5P4mDDO75jA8G5NOT25fmDfcvrpGQiO0GmxT4EmBaX8QFCQkBQbSVJsJOd3bARYbSRz1u5m+qqdfJGewyeLttChcR1G9mnGxac1ITI0wP55b18GGd/AwIcgWmc+qCptaFaqFjh8tIRvf9/Oxws3sXbnIRpEhTJqQEtG9mkWOMlhwnDY/jvcu1wnvTsOTxuatfOuUrVAVFgw1/RKZtq9/Zl8Wx86NKnDv6atZcCLc5i8dCuu2j4/VPY82PgT9H9AE8Ip0qSgVC0iIvRsEcuEm3vx5eg+NGsQxYNfrmD427+wMqeWzixrDMx6Euo0hR632B2N39OkoFQt1b1ZLF/c3odXr+rK9gOFXPL2L7wycz1FJS67Q6te636AbWlWW0JIhN3R+D1NCkrVYiLC8G6JzLp/IBd3bcLrszcw/O1f2LT3sN2hVQ+X0xq53KAVnKbzG1UHTQpKBYC6kSG8ctVpvDOyOzn7jnDRmwv4cfVOu8M6dSsmw561cNZj4AiQBnUv06SgVAA5v2Mjpt59Bi3iohg1IZ1XZ64/4UA8n+csgbkvQKMu0P7Y5d9VVWlSUCrAJMVG8vntfbiieyJjZm/gr5OX+2c7w4rPYF82DHpEZ0GtRlpvKRWAwoIdvHh5F5JjI3l55np2Hizk3etTiQrzk68EZwnMe8mqEtpeYHc0tYqmV6UClIhw99mtefmKrizOzuP6D5ZwsLC48hf6gpWT/1clBPK0Hl6gSUGpAHdZ90TeHNGN5Vv3c917izlQ4OOJwVkCc1/UKsFLNCkopbigc2PGXdedtTsOcdNHSykoKrE7pBMrqxIe1irBCzQpKKUAOKdDAmOuPo1lW/Yx+pPffLPxuawtoTO0HWJ3NLWSJgWlVJkLOjfmueGdmbt+D3/9fLmtiwwd18rPIS8LBmqV4C1+0tVAKVVTru6ZzP4jxTw/bS2N64bz9yHt7Q7J4iyBeS9aVUK7C+2OptbSpKCU+pPbB6awff8Rxs/Lok1CDJd3T7Q7JFj1hVUlXPVfrRK8SG8fKaWO6/GhHeib0oC/f7WS9M159gZT2uMoQasEb9OkoJQ6rhBHEG9fezqN64Vz24R0tu0/Yl8wq76AvI0w6CGtErxMk4JS6oTqRYby/g2pFBa7uOvT3yh22tAjqbTHUUJnaKtVgrdpUlBKVahVwxiev6wzy7bs56UZ62o+gFVfQm4mDHxQ5ziqAfobVkpVamiXJozs3Yzx87KYvWZXzb2xy2n1OEroBO2G1tz7BjBNCkopjzx6YXs6NqnDXz9fXnPtC1ol1Dj9LSulPBIe4uCta06nxGn46+TfcXl7YJvLafU4atgR2l3k3fdSZTxOCiLiEJEmIpJcunkzMKWU72keF8U/hnZgUVYe//l1k3ffbNVXkLvB6nGkVUKN8eg3LSJ3A7uAmcD37m2qB68bLCLrRCRTRB4+wTlXikiGiKwWkU9PInallA2uSE3knPYNeXH6WjJ3H/LOm7ic1qpqDTtolVDDPE2/9wJtjTEdjTGd3VuXil4gIg7gLeACoAMwQkQ6HHNOa+ARoJ8xpiNw30l/AqVUjRIRnru0M5GhDh6YvNw73VRLq4SBWiXUNE9/21uBAyd57Z5ApjEmyxhTBEwCjl1I9VbgLWPMPgBjzO6TfA+llA0axoTz3PDOrMg5wNifN1bvxUt7HDXsAO2HVe+1VaU8nfsoC/hZRL4HjpYeNMa8UsFrmmIlk1I5QK9jzmkDICK/AA7gCWPM9GMvJCKjgFEAycnalKGUL7igc2OGdmnMmz9lMqRzY1o1jK6eC6/+Gvauhys+1CrBBp7+xrdgtSeEAjHltlMVDLQGBgEjgHdFpN6xJxljxhtjUo0xqfHx8dXwtkqp6vDPizoSHhLE379aWT29kUp7HMW3h/bH3lhQNcGjSsEY8ySAiES7H+d78LJtQFK5x4nuY+XlAIuNMcVAtoisx0oSSz2JSyllr/iYMB69sD0PfbmSz9K2MqLnKVbyq7+Gvevg8v9olWATT3sfdRKRZcBqYLWIpItIx0pethRoLSItRCQUuBqYcsw532BVCYhIHNbtpKyTiF8pZbMrU5Po3TKW535Yw+6DhVW/kMtpzXEU3w46XFJ9AaqT4mmbwnjgAWPMHAARGQS8C/Q90QuMMSUichcwA6u94ANjzGoReQpIM8ZMcT93nohkAE7gb8aY3Cp/GqVUjRMRnhvemcFj5vP092t4Y0S3ql0o4xvYsxYu/8Bnq4Ti4mJycnIoLDyF5Odl4eHhJCYmEhISUqXXizGV3wcUkeXGmK6VHasJqampJi0trabfVilViVdmruf12RuYNKo3vVs2OLkXu1wwto+1P/pXCHJUf4DVIDs7m5iYGBo0aID44BTexhhyc3M5dOgQLVq0+MNzIpJujEmt7BqepuMsEXlcRJq7t8fQ2zxKqXJGD0yhab0InpiympKTHbtQWiUM+JvPJgSAwsJCn00IYFVtDRo0OKVKxtOkcBMQD3zl3uLdx5RSCoCIUAePD23P2p2H+GTRZs9f6HJZPY7i2kLH4d4LsJr4akIodarxeZQUjDH7jDH3GGNOd2/3lg44U0qpUud3bMQZreJ4ZeZ6cvOPVv4CgDXfwp417plQfbdK8BXTp0+nbdu2tGrViueff77ar19hUhCR19w/vxORKcdu1R6NUsqviQhPDOtAQZHTswV5yqqENn5RJdjN6XRy5513Mm3aNDIyMpg4cSIZGRnV+h6V9T6a4P7572p9V6VUrdWqYQw39mvOewuyubZXMzon1j3xyWumwO4MuOx9rRI8sGTJElq1akXLli0BuPrqq/n222/p0KFDJa/0XIVJwRiT7t49zRgzpvxzInIvMLfaIlFK1Rr3nN2aL3/bxrM/ZDDx1t7Hv8/t51XCk9+tJmP7wWq9ZocmdfjnRSceArZt2zaSkv43JjgxMZHFixdXawyeNjTfcJxjf6nGOJRStUhMeAj3ndOaRVl5/LT2BPNcZnwDu1fDAG1L8CUVVgoiMgK4BmhxTBtCDJDnzcCUUv5tRM9kPvxlE/+atpaBbeIJdpT7G9RZAnOes+Y46nSpfUGegor+oveWpk2bsnXr/+YZzcnJoWnTptX6HpW1KfwK7ADigJfLHT8ErKjWSJRStUqII4gHB7fj9k/SmZyWwzW9ys2LtGKStV7CVZ9olXASevTowYYNG8jOzqZp06ZMmjSJTz+t3rXJKmtT2AxsBvpU67sqpQLC+R0TSG1Wn1dmrmfYaU2IDguGkqPw8wvQpBu0G2p3iH4lODiYN998k/PPPx+n08lNN91Ex47VW7F4OiFebxFZKiL5IlIkIk4Rqd4WFqVUrSMi/P3C9uzNP8r4ee5JEH77GA5sgbMeAx8fCOaLhgwZwvr169m4cSOPPvpotV/f04bmN7HWO9gARAC3YC21qZRSFTo9uT4XdmnMu/Oy2JWbZ82EmtwXUs62OzR1HB5PRWiMyQQcxhinMeY/wGDvhaWUqk0eOr8dJS4Xv33xEuTvgrMf1yrBR3k6dXaBe02E30XkRazGZ9+c21Yp5XOSG0Tyl+4N6L38Y440G0REsxPOuq9s5ukX+0isNRHuAg5jrah2mbeCUkrVPvdGzaS+5DM2aITdoagKeLocZ+mUh0eAJ70XjlKqVjq8l+j0cayrP5A31sVw0a5DtE6ojmXeVXWrbEK8lSKy4kRbTQWplPJzc1+A4gIaDX+OqNBgXpm53u6I1AlUViloJ2Kl1KnZmwlpH0D3G6ib3ImbzwhlzOwNrMw5UPFkeeq4brrpJqZOnUrDhg1ZtWpVtV+/wkrBGLO5oq3ao1FK1T6zn4DgcBj0CAC39G9BvcgQ/v2jB1Nrqz/5y1/+wvTp0712fU8Hrx0SkYPurVAHrymlPLJlEaz5DvrdC9ENAWuyvNsHpjB3/R6WbtIp1E7WgAEDiI2N9dr1PW1oLmsREmsO3IuB3t4KSilVCxgDPz4G0Y2gz51/eOqGPs15f0E2L81Yx2ejTjC1tq+b9jDsXFm912zUGS6o/tXUTsZJjzUwlm+A870Qj1Kqtsj4BnKWwlmPQmjUH56KCHVw15mtWJKdx4LMvTYFqI7Ho0pBRMrPbRsEpAKFXolIKeX/Sopg1pPQsAOcdu1xT7m6ZxLvzN3Ia7M2cEarOP+rFmz+i95bPK0ULiq3nY81dfbF3gpKKeXnFo+Dfdlw7lMnnBo7LNjBHWe2In3zPq0WfIhHScEYc2O57VZjzLPGmBMsp6SUCmiHdlrjElqfD63PrfDUK1ITaVw3nDGzNmCMqaEA/duIESPo06cP69atIzExkffff79ar+/p7aOWwBisxmUDLATuN8ZkVWs0Sin/N/Mf4CyCwf+q9NTSauHxb1bxS2YuZ7SOq4EA/dvEiRO9en1Pbx99CkwGGgNNgM8B70amlPI/mxfCis+g793QIMWjl1zprhZem7VeqwUf4GlSiDTGTDDGlLi3T4BwbwamlPIzLidM+xvUaQr9/+rxy8KCHdwxKIW0zfv4dWOuFwNUnvA0KUwTkYdFpLmINBORB4EfRCRWRLw3ikIp5T/S/2P12z/vmT91Qa3MlT2SaFRHqwVf4Ol6Cle6f952zPGrsdoYWlZbREop/5O/G2Y/Dc37Q8fhJ/1yq20hhX98u5qFG3Pp28p32xaMMT7dffZUk6qnvY9aVLCdMCGIyGARWScimSLycAXnXSYiRkRSq/IhlFI2m/4IFBfAhS9XeUW1K1OTSKgTxms+3BMpPDyc3Nxcn43PGENubi7h4VW/u+9p76MQYDQwwH3oZ+AdY0xxBa9xYK3jfC6QAywVkSnGmIxjzosB7gUWn3T0Sin7bZgJq76wJryLb1vly4SHOLhjUCv+OWU1C7Ny6Zvie9VCYmIiOTk57Nmzx+5QTig8PJzExMQqv97T20djgRDgbffjke5jt1Twmp5AZmm3VRGZhDXgLeOY854GXgD+5mEsSilfcTQfpj4AcW3hjPtP+XJX9Uji7Z8zeW3WBp9MCiEhIbRo0cLuMLzK04bmHsaYG4wxP7m3G4EelbymKbC13OMc97EyInI6kGSM+d7jiJVSvmPOc3BgC1w0BoLDTvly4SEORg9MYUl2Hgu1J5ItPE0KThEp63TsHszmPJU3FpEg4BWg0r5rIjJKRNJEJM2XyzalAsqWRbB4LKTeBM36VNtlr+6ZTMOYMF6bpauz2cHTpPA3YI6I/CwiPwM/UfmX+TYgqdzjRPexUjFAJ+BnEdmENVp6yvEam40x440xqcaY1Pj4eA9DVkp5zdF8+Pp2qJtozW9UjcJDHIwelMJirRZs4WlS+AV4B3ABee79hZW8ZinQWkRaiEgoVvfVKaVPGmMOGGPijDHNjTHNgUXAMGNM2kl+BqVUTZv5OOzbBJeMg7CYSk8/WSPc1cKY2Vot1DRPk8LHQAusRuE3sMYlTKjoBcaYEuAuYAawBphsjFktIk+JyLCqh6yUslXmLGvN5T53QvN+XnmL8BAHtw9MYVGWVgs1TTzpbysiGcaYDpUdqwmpqakmLU2LCaVscTgXxvWD8Lowai6EeG+2m8JiJ/1fnENKfBSTRlVfm0WgEpF0Y0ylY8E8rRR+E5Gy5TdFpBeg38xKBRKXC76+DQpy4dLxXk0I8MdqYVGWVgs1xdOk0B34VUQ2uRuFFwI9RGSliKzwWnRKKd/x6xjInGlNid24a4285bW9komPCWPMrA018n7K88Frg70ahVLKt21eaM1t1OESSL25xt62tFp4emoGi7Ny6dWyQY29d6DydO6jzRVt3g5SKWWj/N3wxU1QLxmGvV7luY2q6tpeycRFhzFmtlYLNcHT20dKqUBUchQ+uw6O7IMrP7IamGuYVS205NeNuSzWtgWv06SglDo+Y2Dq/bB1MQwfW2PtCMdzba9mWi3UEE0KSqnjW/gW/P5fGPhwldZIqE4Rof+rFpZk59kaS22nSUEp9WcZ31qjltsPg4EP2R0NUFothOooZy/TpKCU+qONc+DLWyCxBwx/B4J842vCqhZS+CUzl6WbtFrwFt/4r62U8g05aTDpWmjQGq75DEIj7Y7oD8qqBR234DWaFJRSlh0r4L+XQ3Q8jPwKIuq9nZORAAATA0lEQVTbHdGfRIQ6uG1ACgsy92q14CWaFJRSkJMOHw2FkCgY+Q3ENLI7ohO6tneyVgtepElBqUC3eSF8fDGE14Mbf4BY315uMjI0mFEDWrIgcy9pWi1UO00KSgWyVV9ZCSEmAW6cBvWb2R2RR67r3YwGUaE6bsELNCkoFYiMgQWvwRc3QpNucPNMqNu08tf5iNJqYf6GvaRv1mqhOmlSUCrQHD1kJYNZ/4SOl8L130JkrN1RnbSRfZoRGxXKa9q2UK00KSgVSHaugvFnWoPTznkSLv/A6+sieEtkaDCjB6Ywf8NeXZ2tGmlSUCoQOEtg3r9h/CAoPADXT4Ez7qvxGU+r28g+zWhUJ5wXZ6zFk1UkVeU0KShV2+Wkw/vnwE9PQ7sL4Y5F0KK/3VFVi/AQB/ee05plW/YzM2OX3eHUCpoUlKqtDu6Ab++E986Cg9utW0VXfgRRtWuhmiu6J9IyLop//7gOp0urhVOlSUGp2ubgDpj2EIzpCssnQd974K406HSZ3ZF5RbAjiAfOa8P6Xfl8+/s2u8Pxe54ux6mU8mXGwOZfYel7sGaK9fi0EdD//3x+MFp1GNKpMR2bbOSVmesZ2qUJocH6925VaVJQyl8ZAztXQMYUqzdR7gZrVHKv26HHLQGRDEoFBQkPDm7HDR8sYeKSLdzQt7ndIfktTQpK+QuXE/KyYMtC2PQLbFoAB3NAgqD5GdDvXusWkY/NbFpTBrSOo1eLWN74KZPLuycSFaZfb1WhvzWlKmMMOIug6DAUF0DxETCuE5wsEOSwvqglyL3vKLd/zPHSYyWF1qCyo/nWz4K9cCDHaiA+kAN71sKedVByxHqbyDho1hcGPmj1KIqKq7Ffh68SsaqFy8b+yn9+yeaus1rbHZJf0qSgApsx1hdvbqZ1+yUvGw7thPxdkL8bDu+GwoNgnPbEJw6IaQxxrSH1JkjoAE1TIb6t348x8IbuzepzbocExs3N4qoeycTHhNkdkt/RpKAChzGwL9vqt78tHbb/Zo3wLT78v3OCw60v4egE64u3RX8IrwshkRAaZf0MiTzxamTGWLd5jMtKJH/Ydx1z3P2cy2WNKg6NhrAYCKtjrWVQt6kVR5CjZn4/tcTDF7Tj/Ffn8eqs9Tw3vLPd4fgdTQqqdsvfA1lzIGsuZP1s3YMHCI6Axl3h9JHWX+ENWllbTBOfWX5SVU1KfDTX9W7Gxws3cX2fZrRrVMfukPyKJgVV+xzYBmu+s7pmbv4VMNZf3i0GQIsHIKkXxLcDh/7vX1vdd05rvl62jWe/X8PHN/VE9Fabx/Rfhaodig5b3TKXfQKbf7GONexgNcS2GWxVBXobJmDUiwzlnrNb8/TUDH5ev4cz2za0OyS/oUlB+bddq2HxOGuxmKJ8iE2Bsx6DDsMhrpXd0SkbjezdjAkLN/Hs92vo3yqOYIfeFvSEV39LIjJYRNaJSKaIPHyc5x8QkQwRWSEis0XEP5Z9UvZyuWD9j/DRMBjbF1Z+AR0ugRunw93pMOBvmhAUocFBPDKkPZm78/lk0Wa7w/EbXqsURMQBvAWcC+QAS0VkijEmo9xpy4BUY0yBiIwGXgSu8lZMys85S2DlZFjwKuxdbzUKn/MEnH6DXy4So7zvvA4J9G8dx8s/rmdIl8Y0jPHPtSNqkjcrhZ5ApjEmyxhTBEwCLi5/gjFmjjGmwP1wEZDoxXiUv3KWwO8T4a0e8M1oCA6DS9+D+1bAGfdrQlAnJCI8OawjR0tc/OuHtXaH4xe8mRSaAlvLPc5xHzuRm4Fpx3tCREaJSJqIpO3Zs6caQ1Q+zRhY9SW83Qu+uR1CouCq/8Jt86HLFeAIsTtC5Qdaxkdz28CWfL1sm67Q5gGfaHkRkeuAVOCl4z1vjBlvjEk1xqTGx8fXbHDKHluXwPvnwhc3gSMMrpwAt82D9kN1JK86aXcMakVi/Qj+8e0qip0nmqJEgXeTwjYgqdzjRPexPxCRc4BHgWHGmKNejEf5g32b4PMbrYSwfytc/BbcPh86DNNBZarKIkIdPDmsIxt25/PBgmy7w/Fp3uySuhRoLSItsJLB1cA15U8QkW7AO8BgY8xuL8aifF1RAcz/N/z6hjXfz8CHoe/dEBZtd2Sqlji7fQLntE/gtVkbGNK5MUmxgTmbbGW89qeXMaYEuAuYAawBJhtjVovIUyIyzH3aS0A08LmI/C4iU7wVj/Jh62dY7QbzX4aOw+Ge3+DMRzQhqGr35MUdcQQJD3+1AmN06c7j8ergNWPMD8APxxz7R7n9c7z5/srHHdgG0x+ypqSIawM3TK01C8or39S0XgSPDGnHo1+vYuKSrVzTK9nukHyOjmhWNc/ltEYhz3kOXCVw1uPWOsLBoXZHpgLANT2T+X7FDp77YQ0D28bTtF6E3SH5FG25UzVr7wb4YDDM+Dsk94E7FsGA/9OEoGqMiPDCZV1wGcPDX+ptpGNpUlA1w+W0GpHHnWGNRr70Xbj284BaR1j5jqTYSB6+oB3zN+xlgk6B8Qd6+0h5395M+PYO2LoY2lwAF70GMY3sjkoFuJG9mzFn7W6e+X4NPVvE6roLblopKO9xOWHhWzCun7XG8PB3YMRETQjKJ4gIL13RlTrhIdwzcRmFxTYtuepjNCko78jdCB9eaLUdtBgIdyyGrlfraGTlU+Kiw3jlyq6s35XPM99nVP6CAKBJQVUvlwsWjYWx/WBXBlwyFq75DOo0tjsypY5rQJt4Rg1oySeLtjBl+Xa7w7Gdtimo6pOXBd/eZa181upcGPY61Glid1RKVer/zmvL71v28+AXy0mJj6Jjk7p2h2QbrRTUqXO5YPE7VnWwc6U1X9G1n2tCUH4jNDiIt649nfqRoYz6OJ28w0V2h2QbTQrq1ORlw0cXwbQH3eMOFkK367TtQPmd+Jgwxl3XnT35R7nr098CdjZVTQqqalwuWPKuVR3sWA7D3oDrvoS6uk6S8l9dk+rxr+Gd+XVjLo9+vTIgB7Zpm4I6eXnZMOVu2DQfUs6Ci16HekmVv04pP3BZ90Q25xXw+uwNJNQJ56/ntbU7pBqlSUF5zuWCtPdh5j8hyGFVB91G6q0iVevcf05rdh8s5I2fMomLDuOGvs3tDqnGaFJQnsnLdvcsWgApZ1s9i/RWkaqlRIRnLulE7uEi/jllNUECI/s0tzusGqFtCqpiLhcsHg9j+8LOFdp2oAJGsCOIt645nXPaN+Txb1fz8cJNdodUI7RSUCe2KwOm3mfNWaTVgQpAocFBvH1td+74bzr/+HY1+wuKufusVkgtvmWqlYL6s+IjMPspeKe/NdX1xW9rdaACVmliuPT0prwycz2PfLWSklrcXVUrBfVHG3+CqQ/Avmzoeg2c9wxENbA7KqVsFRocxMtXdKVJ3QjenJPJ5twCXh/RjfiYMLtDq3ZaKSjLwe3w5S0wYThIEFw/BYaP1YSglJuI8H/nt+XfV3Tlty37uPD1+SzOyrU7rGqnSSHQFR+BeS/BG90hYwoMeBBG/wotB9odmVI+6fLuiXxzZz8iQx2MeHcRz0zN4EhR7Zl2W28fBSpjYM0U+PEx2L8F2g+D856G+s3tjkwpn9e+cR2+u/sMnp+2lvcWZDNzzS7+MbQDZ7Vr6PeN0FopBKIti621DiZfD6ExcMN3cNUETQhKnYSY8BCeHd6Zibf2xiHCzR+lcdX4RaRvzrM7tFMi/ja3R2pqqklLS7M7DP+0cyXMfho2zICohjDoYTj9BnBowajUqSh2upi0ZAtjZm9gb34R3ZLrcWO/FpzXIYHwEIfd4QEgIunGmNRKz9OkEAB2rYb5L8OqLyG8LvS7D3rdBqFRdkemVK1y+GgJX6Tn8J9fstmUW0B0WDDndkjgnPYJ9GhRn4Yx4bbFpklBwZZFMP8VqzIIibISQb97IKK+3ZEpVau5XIYFmXv5fsUOpq/eyYEjxQAkxUaQEh9N8wZRxMeEERMeTEx4MKEOByIglE4lJmVTipX/iu7YpA5JsZFVikmTQqByFsO6H6wlMbcshIhY6D0aetwCkbF2R6dUwCl2uli9/SBLs/P4fet+NuUeZtPewxyuQo+lZy7pxHW9m1UpDk+Tgt5Mri0Obof0jyD9Q8jfCXWTYPDzcPr1eptIKRuFOII4LakepyXV+8PxwmInhwpLOFRYTLHTYDAYY1UGpful1YJg7TSq6/3bT5oU/FnRYVg3DVZ+DhtmgnFBq7Ohx2vQ+jxremullE8KD3EQHuLwuVHRmhT8TdFhyJ4Hq76Ctd9D8WGIaQJ974Luf4HYlnZHqJTyY5oUfJ0xsG+TVQlsmAHZ88F5FMLrQZcroPMVkNwXgnTIiVLq1Hk1KYjIYGAM4ADeM8Y8f8zzYcDHQHcgF7jKGLPJmzH5PGcJ7Flj9RzashA2L4RD263nGrSyGoxbnwvN+kFwqL2xKqVqHa8lBRFxAG8B5wI5wFIRmWKMySh32s3APmNMKxG5GngBuMpbMfmUkiI4sNWajXTPemsswa6VsHutVQmAdVuoWR9I7mOthdwgxd6YlVK1njcrhZ5ApjEmC0BEJgEXA+WTwsXAE+79L4A3RUSMv/WTdTmte/3FBeV+FsCRPDi8x73ttbZDO6zbQQe3WQ3DpaLiIaET9BoFjbpAUi+ol6zrHyulapQ3k0JTYGu5xzlArxOdY4wpEZEDQANgb7VH89sE+PUNwFhfxsb9k9I+YKbcvsuDfff5ziIoKaz8/UOjISoOohOgWV9rnqHSLTYFYhKq/SMrpdTJ8ouGZhEZBYwCSE5OrtpFImOhYTtrrQDE+ilygn08OMe97wixxgGEREJopPXlX7ofEWtVAFFxEBJRDb8JpZTyLm8mhW1AUrnHie5jxzsnR0SCgbpYDc5/YIwZD4wHa0RzlaJpd6G1KaWUOiFv9mNcCrQWkRYiEgpcDUw55pwpwA3u/cuBn/yuPUEppWoRr1UK7jaCu4AZWF1SPzDGrBaRp4A0Y8wU4H1ggohkAnlYiUMppZRNvNqmYIz5AfjhmGP/KLdfCFzhzRiUUkp5TofBKqWUKqNJQSmlVBlNCkoppcpoUlBKKVVGk4JSSqkyfrccp4jsATZX8eVxeGMKDd+mnzkw6GcODKfymZsZY+IrO8nvksKpEJE0T9YorU30MwcG/cyBoSY+s94+UkopVUaTglJKqTKBlhTG2x2ADfQzBwb9zIHB6585oNoUlFJKVSzQKgWllFIVCJikICKDRWSdiGSKyMN2x+NtIvKBiOwWkVV2x1JTRCRJROaISIaIrBaRe+2OydtEJFxElojIcvdnftLumGqCiDhEZJmITLU7lpogIptEZKWI/C4iaV59r0C4fSQiDmA9cC7WsqBLgRHGmIwKX+jHRGQAkA98bIzpZHc8NUFEGgONjTG/iUgMkA5cUsv/OwsQZYzJF5EQYAFwrzFmkc2heZWIPACkAnWMMUPtjsfbRGQTkGqM8fq4jECpFHoCmcaYLGNMETAJuNjmmLzKGDMPa42KgGGM2WGM+c29fwhYg7UOeK1lLPnuhyHurVb/pSciicCFwHt2x1IbBUpSaApsLfc4h1r+ZRHoRKQ50A1YbG8k3ue+lfI7sBuYaYyp7Z/5NeBBwGV3IDXIAD+KSLp7zXqvCZSkoAKIiEQDXwL3GWMO2h2PtxljnMaY07DWQe8pIrX2dqGIDAV2G2PS7Y6lhp1hjDkduAC403172CsCJSlsA5LKPU50H1O1jPu++pfAf40xX9kdT00yxuwH5gCD7Y7Fi/oBw9z32CcBZ4nIJ/aG5H3GmG3un7uBr7FuiXtFoCSFpUBrEWkhIqFYa0FPsTkmVc3cja7vA2uMMa/YHU9NEJF4Eann3o/A6kyx1t6ovMcY84gxJtEY0xzr3/FPxpjrbA7Lq0Qkyt1xAhGJAs4DvNarMCCSgjGmBLgLmIHV+DjZGLPa3qi8S0QmAguBtiKSIyI32x1TDegHjMT66/F39zbE7qC8rDEwR0RWYP3xM9MYExDdNANIArBARJYDS4DvjTHTvfVmAdElVSmllGcColJQSinlGU0KSimlymhSUEopVUaTglJKqTKaFJRSSpXRpKBUBUSknojc4d5vIiJf2B2TUt6kXVKVqoB7DqWpgTLTrFLBdgeglI97HkhxTzi3AWhvjOkkIn8BLgGigNbAv4FQrMFzR4Ehxpg8EUkB3gLigQLgVmNMrR1xrPyf3j5SqmIPAxvdE8797ZjnOgGXAj2AZ4ECY0w3rJHk17vPGQ/cbYzpDvwf8HaNRK1UFWmloFTVzXGv23BIRA4A37mPrwS6uGdr7Qt8bk3LBEBYzYeplOc0KShVdUfL7bvKPXZh/dsKAva7qwyl/ILePlKqYoeAmKq80L2WQ7aIXAHWLK4i0rU6g1OqumlSUKoCxphc4BcRWQW8VIVLXAvc7J7hcjW1fBlY5f+0S6pSSqkyWikopZQqo0lBKaVUGU0KSimlymhSUEopVUaTglJKqTKaFJRSSpXRpKCUUqqMJgWllFJl/h/62n5owzRT6AAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plot_population(opt_dynamics)" ] } ], "metadata": {}, "nbformat": 4, "nbformat_minor": 2 }