Untitled.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Collecting krotov\n",
      "  Downloading krotov-1.2.1-py2.py3-none-any.whl (81 kB)\n",
      "\u001b[K     |████████████████████████████████| 81 kB 1.7 MB/s eta 0:00:01\n",
      "\u001b[?25hRequirement already satisfied: qutip>=4.3.1 in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from krotov) (4.6.2)\n",
      "Collecting loky\n",
      "  Downloading loky-3.0.0-py2.py3-none-any.whl (69 kB)\n",
      "\u001b[K     |████████████████████████████████| 69 kB 5.4 MB/s eta 0:00:011\n",
      "\u001b[?25hRequirement already satisfied: grapheme in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from krotov) (0.6.0)\n",
      "Requirement already satisfied: threadpoolctl in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from krotov) (2.1.0)\n",
      "Requirement already satisfied: scipy in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from krotov) (1.6.2)\n",
      "Requirement already satisfied: numpy in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from krotov) (1.20.1)\n",
      "Requirement already satisfied: glom in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from krotov) (20.11.0)\n",
      "Requirement already satisfied: packaging in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from qutip>=4.3.1->krotov) (20.9)\n",
      "Requirement already satisfied: boltons>=19.3.0 in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from glom->krotov) (21.0.0)\n",
      "Requirement already satisfied: attrs in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from glom->krotov) (20.3.0)\n",
      "Requirement already satisfied: face>=20.1.0 in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from glom->krotov) (20.1.1)\n",
      "Requirement already satisfied: cloudpickle in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from loky->krotov) (1.6.0)\n",
      "Requirement already satisfied: pyparsing>=2.0.2 in /home/mcditoos/anaconda3/lib/python3.8/site-packages (from packaging->qutip>=4.3.1->krotov) (2.4.7)\n",
      "Installing collected packages: loky, krotov\n",
      "Successfully installed krotov-1.2.1 loky-3.0.0\n"
     ]
    }
   ],
   "source": [
    "!pip install krotov\n",
    "import qutip\n",
    "import numpy as np\n",
    "import scipy\n",
    "import matplotlib\n",
    "import matplotlib.pylab as plt\n",
    "import krotov\n",
    "import math\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "import dill \n",
    "dill.settings['recurse'] = True\n",
    "\n",
    "H0=dill.load(open(\"H_0\", \"rb\"))\n",
    "He=dill.load(open(\"H_campo\",\"rb\"))\n",
    "Hde=dill.load(open(\"H_de\",\"rb\"))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "iter.      J_T    ∫gₐ(t)dt          J       ΔJ_T         ΔJ  secs\n",
      "0    -3.03e+00    0.00e+00  -3.03e+00        n/a        n/a     0\n",
      "1    -1.17e+04    2.12e+04   9.50e+03  -1.17e+04   9.50e+03     0 *\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "WARNING:matplotlib.legend:No handles with labels found to put in legend.\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAARQAAAEvCAYAAACEzT4XAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8QVMy6AAAACXBIWXMAAAsTAAALEwEAmpwYAAA4ZklEQVR4nO3dd3hUVfrA8e+bQigJNQhICyhVQISIIkU0lEBAFFnUdUVdXcTVXXF3LbvY8afi2tZ1LVixLIoGBQuoiAWkSChLlU4ggPQakpBk3t8fM8EQJiQkd+ZOMu/nee6TyS3nvLkzvNxz7p1zRFUxxhgnRLgdgDGm8rCEYoxxjCUUY4xjLKEYYxxjCcUY4xhLKMYYx0S5HUBZxcfHa0JCgtthGBN2Fi1atEdV6/vbVmETSkJCAmlpaW6HYUzYEZH04rZZk8cY4xhLKMYYx1hCMSEhPz+fnJwct8Mw5WQJxbhu9+7dnHvuufTt2xf7blnFZgnFuGrfvn3069ePVatWMWfOHKZOnep2SKYcLKEY1xw6dIjk5GRWr17N559/TuvWrbn//vvxeDxuhxbS7rrrLv7whz/w7LPP8uWXX7J169aQubKrsLeNTcWWmZlJSkoKS5YsYcqUKQwcOJCDBw9yzTXXMHnyZK6++mq3QwxJW7Zs4amnnqJq1apkZ2cfXx8bG0v79u1PWpo3b05ERPCuGyRUMtvpSkxMVHsOpWLKzs5m8ODBfPvtt0yaNIkRI0YA4PF46Ny5Mzk5OaxcuZKoKPv/rqj//ve/XHvttSxatIimTZuyevVqVq1adcKyY8eO4/vXqVOHJUuW0Lx5c8diEJFFqprod6OqVsila9euaiqenJwcTUlJUUAnTpx40vZPPvlEAX3jjTdciC70jR49WuPi4jQ3N7fYffbt26c//vijPv/88wrov/71L0djANK0mH+XrieGsi6WUCqe3NxcvfLKKxXQl19+2e8+Ho9HExMTtXnz5pqTkxPkCENfhw4dtH///qXev3Xr1jpw4EBHYzhVQrFOWRMU+fn53HDDDaSmpvLMM89wyy23+N1PRHj00UdJT0/n9ddfD3KUoW3fvn2sWLGCXr16lfqY5ORkvvvuO7KysgIY2a8soZiAU1VGjx7Ne++9x6OPPsqdd955yv379+9Pz549efTRR4P2D6Ei+PHHHwFOO6FkZWUxe/bsQIV1AksoJqBUlTFjxvDaa6/xj3/8g7Fjx5Z4TMFVyvbt23n55ZeDEGXFMGfOHKKjo+nWrVupj7n44ouJiYlhxowZAYyskOLaQqG+WB9K6PN4PHrPPfcooGPGjFGPx3Nax/fr10/r16+vhw8fDlCEFUv37t21e/fup31c//79tV27do7FgfWhGDc8/fTTjB8/nltuuYVnnnkGETmt48eNG8fu3bt5/vnnAxRhxZGVlUVaWtppNXcKFDw8mJ5e7KgDjrGEYgLik08+4e6772b48OG8+OKLp51MAC644AKGDBnCP//5Tw4cOOB8kBXITz/9RG5uLj179jztYwcOHAjAl19+6XRYJ7GEYhy3ePFirr32WhITE5k4cWK5ntR85JFHOHDgAM8884yDEVY8BZ2qPXr0OO1j27RpQ/PmzYPSj2IJxThq27ZtDBkyhHr16jFt2jSqV69ervI6d+7MiBEjePbZZ9mzZ49DUVY8c+bM4ZxzzqFu3bqnfayIkJyczMyZM8nNzQ1AdL8KakIRkaYi8q2IrBaRlSJyh2/9QyKyTUSW+pZBwYzLOCMzM5PLLruMQ4cO8emnn9KwYUNHyn3ooYc4evQo48ePd6S8iiY/P5+5c+eWqf+kQHJyMocPH2bevHkORnayYF+h5AF/VdV2wIXAbSLS3rftWVXt7Fu+CHJcppw8Hg8jR45kyZIlTJo0iXPPPdexstu1a8fvfvc7XnjhhRO+pxIuli1bxuHDh8uVUC699FKioqKYPn26g5GdLKgJRVV3qOpi3+vDwGqgcTBjMIExduxYpkyZwtNPP83gwYMdL/+BBx4gLy+Pxx57zPGyQ11B/0lZOmQL1KxZkx49egS+H6W4+8mBXoAEYAtQE3gI2AwsA94A6hRzzCggDUhr1qyZY/fVTfm88cYbCuioUaNO+1mT0zFq1CiNjo7WzZs3B6yOUDR8+HB14vP++OOPK6Dbt28vVzmE2pcDgVhgETDM93sDIBLvFdP/AW+UVIY92BYavvvuO42OjtakpCQ9duxYQOvaunWrxsTE6M033xzQekKJx+PRhg0b6rXXXlvuspYsWaKAvvXWW+Uq51QJJeh3eUQkGkgF3lPVKQCqulNV81XVA7wKlP7ZYuOa9evXM2zYMFq2bMmHH35IdHR0QOtr0qQJo0eP5s0332TdunUBrStUbNiwgV9++aVczZ0C5557Lg0bNgxosyfYd3kEeB1YrarPFFrfqNBuVwArghmXOX379+8nJSUFEeHzzz+nTp06Qan373//OxEREbz22mtBqc9tBf0n5emQLSAiDBgwgK+++or8/Pxyl+dPsK9QegDXAZcWuUX8pIgsF5FlwCXAqb+OalyVm5vL8OHD2bRpE1OmTOGss84KWt0NGjSgd+/efP7550Gr001z5syhbt26tGvXzpHykpOT2bdvX8Bm3Qz2XZ45qiqq2kkL3SJW1etUtaNv/WWqGn73BiuArKwsFi9ezE033cSsWbN49dVX6d27d9DjSElJYeXKlWzevDnodQfb7Nmz6dGjh2Pjwvbr1w8RCVizx56UNSfxeDxs3LiRqVOnMm7cOEaMGEHbtm2JjY2la9euvPPOO4wdO5brr7/elfhSUlIAKv1Vys6dO1m3bp0jzZ0C9erVo1u3bgFLKDYKsAHghx9+4N1332X58uWsWLGCI0eOHN/WsmVLOnXqxIgRI+jUqRPnnnsurVq1ci3W1q1b06pVKz7//HNuu+021+IItDlz5gDle/7En+TkZMaNG8fevXupV6+eo2VbQjFs2bKFQYMGERUVRefOnbnhhhvo1KkTHTt2pEOHDsTGxrod4klSUlJ46aWXyMzMpEaNGm6HExCzZ8+mWrVqdO3a1dFyk5OTefjhh5k5cyZXXXWVo2VbQglzqsptt92GqrJ06VISEhLcDqlUUlJSeO6555g1axZDhgxxO5yAmDNnDhdccAFVqlRxtNzzzz+fOnXqMGPGDMcTivWhhLnU1FQ+++wzxo0bV2GSCUDv3r2JjY2ttP0ohw8fZsmSJY43dwAiIyPp378/M2bMcHyWRksoYezAgQP86U9/okuXLvz5z392O5zTUqVKFfr378/nn39e8PR1pTJv3jw8Ho+jHbKFDRw4kF9++YVly5Y5Wq4llDB27733smvXLiZMmFAhZ+lLSUkhIyPD8X8UoWDOnDlERETQvXv3gJTfv39/AMfv9lhCCVNz5szhlVdeYcyYMY53+gXLoEHeYXMqY7Nn9uzZdO7cmbi4uICU36hRIzp37mwJxZRfTk4Oo0aNonnz5jz88MNuh1NmDRs2JDExkc8++8ztUBx17NgxFixYELDmToHk5GR+/PFHDh065FiZllDC0JNPPsnq1at58cUXQ/KW8OlISUlh/vz5lWp4yMWLF5OVlRWUhJKXl8esWbMcK9MSSphZs2YNjz76KFdfffXxJkNFNnjwYFQ1eBNZBYETAyqVRvfu3YmLi3P03FlCCSMej4dRo0ZRvXp1nnvuObfDcUSXLl1o0KBBpWr2zJkzh1atWtGgQYOA1lOlShWSkpKYMWOGY3fKLKGEkTfffJMffviBf/7znwH/sAZLREQEgwYN4ssvvwz4iO7B4PF4mDNnTsCvTgokJyeTnp7OmjVrHCnPEkqY2LlzJ3fddRe9e/fm97//vdvhOGrw4MEcOHCAuXPnuh1Kua1evZp9+/YFvP+kwIABAwDnbh9bQgkTd955J5mZmbzyyiuOfRU+VPTr14/o6OhKcfu44AuBwUooCQkJtG3btvIlFBFJFpE1IrJeRO51O57KZPr06UyaNImxY8fStm1bt8NxXFxcXKUZdGn27Nk0aNAgqINWJScn8/3335OVlVXuskIioYhIJPAfYCDQHrim0Hw9phwyMzO59dZbadeuHffcc4/b4QTM4MGDWbVqFZs2bXI7lHKZM2cOvXr1KtNc0GWVnJxMdnY233//fbnLComEgndQ6vWqulFVjwHvA0PLW2h2dja7du0qd3AV2YMPPkh6ejoTJkwgJibG7XACpjIMurR161bS09OD1twp0Lt3b6pWrerIJGCh8gWOxsDWQr9nABeUt9AHH3yQZ599llatWlG/fv3j60eMGMEf//hHjh496vdZjBtuuIEbbriBPXv2MHz48JO233rrrVx11VVs3bqV66677qTtf/3rXxkyZAhr1qzhlltuOWn7fffdR9++fVm6dCljxow5aftjjz3GRRddxNy5c/nHP/5x0vbnnnuOzp07M3PmTB599NGTtr/yyiu0adOG5557jqeffppGjRpx3333Hd/+zjvv0LRpUz744ANeeumlk47/6KOPiI+P56233uKtt946afsXX3xB9erVefHFF5k8efJJ27/77jsAnnrqqZNu51arVu34B3fcuHF88803J2yvV68eqampgHdA6qJTZzZp0oR3330XgDFjxrB06dITyh4/fjy33347AKNGjWLt2rUnHN+5c+fjt8x/97vfkZGRccL27t278/jjjwNw5ZVXsnfv3hO2JyUlcf/99wPeL9gVbSYMHjyYv/3tbwD06dOHok712du5cyfgff4k2J+9atWqOdKPEipXKP6u7066MS4io0QkTUTSdu/eXWKhI0eOJCYmhlWrVrF69epKcVvxdEycOJHo6GhatGjhdihBUbduXXbs2EFmZqbboZTJoUOHqFq1qqPTuJZW48aNufvuu8s/nEFxE/YEcwG6A18W+v3vwN9PdUxpJ/o6duyYPvzwwxoVFaWNGjXSzz//vFTHVXTff/+9AvrUU0+5HUrQzJw5UwGdOnWq26GUSYcOHXTAgAFuh1EiQmmir2IsBFqJSAsRqQJcDUxzouDo6GgeeOABFixYQL169UhJSeHmm2929AtRoUZVue+++2jUqBG33nqr2+EETa9evYiLi6uQ/Sj79u1jxYoVQXugLVBCIqGoah5wO/Al3gnUJ6vqSifr6NKlC2lpadx77728+eabdOzY0dEvRYWSmTNnMnv2bMaOHUv16tXdDidoKvKgSwUP5QW7Q9ZxxV26hPpSnrmN582bp61bt1ZAb7/9dj1y5EiZywo1Ho9Hu3Xrpk2bNtXs7Gy3wwm6gonblyxZ4nYop+Xuu+/W6OhoPXr0qNuhlIgK0OQJqgsvvJAlS5Zwxx138MILL9C5c+dK8dg2eG+b/vTTTzzwwAOV+jZxcQYOHAhUvNvH06dPp3v37lSrVs3tUMolLBMKcPwbt99++y15eXn07NmTp556yu2wysXj8XD//fdz1llnuTYJl9saNmzI+eefX6G+fbxu3TqWL1/OsGHD3A6l3MI2oRTo06cPy5YtY/jw4dx11128//77bodUZlOmTGHp0qU8+OCDREdHux2Oa1JSUliwYAGlebQgFBQ8d3PFFVe4HIkDimsLhfpSnj4Uf7Kzs7VXr14aExOjc+fOdbTsYMjLy9P27dtr27ZtNS8vz+1wXLVw4UIF9O2333Y7lFJJTEzU888/3+0wSg3rQylZTEwMU6ZMoUmTJgwdOrTCTcT9/vvvs2rVKh555BEiIyPdDsdVXbp0oWHDhhWi2ZOenk5aWhpXXnml26E4whJKIfHx8Xz22Wfk5uYyePBgDh486HZIpZKXl8dDDz1Ep06dKs0Hszwq0qBLU6ZMAag075sllCLatm3LRx99xJo1a7j66qvJy8tzO6QSvf3226xfv55x48ZVurFOyiolJYWDBw+G/N271NRUOnXqxNlnn+12KI6wT58fSUlJvPjii8yYMYM777zT7XBOKScnh4cffpjzzz+/0s7xWxYFgy6FcrNnx44dzJ07t9JcnYAllGL94Q9/4K9//SsvvPACL7zwgtvhFOv1119ny5YtjBs3LqhjaIS6uLg4Lr744pB+HuWTTz5BVStVQnH9bk1ZF6fv8viTl5enl112mUZEROgXX3wR8PpO19GjR/XMM8/Unj17qsfjcTuckPPss88qoJs2bXI7FL+SkpK0TZs2Fe69w+7ylE1kZCTvvfcenTp14qqrrmLFihVuh3SCl19+me3bt/Poo4/a1YkfBfP3zpw50+VITrZ3716+++47rrzyykr13llCKUFsbCyffvopsbGxDB48+PggOG47cuQITzzxBElJSVx88cVuhxOS2rVrR6NGjU4axCkUTJ06lfz8/ErxdGxhllBKoUmTJnz66afs2rWLoUOHOjKYb3m98MIL7Nq1i3HjxrkdSsgSEZKSkvjmm2/KP3CQw1JTU0lISKBLly5uh+IoSyil1LVrV959910WLFjAjTfe6OoH9ODBgzz55JOkpKTQvXt31+KoCJKSkti9e3dINVcPHjzI119/zbBhwypVcwcsoZyWYcOG8cQTT/DBBx8wcOBAtmzZ4koczz33HPv37+eRRx5xpf6KJCkpCQitfpSChycr1d2dAsX11jq9AP8EfgaWAR8DtX3rE4AsYKlvebk05QXjLo8/Ho9H//Of/2iNGjU0NjZWX3rpJc3Pzw9a/bNmzdK4uDgdNmxY0Oqs6Fq3bq2DBg1yO4zjrrjiCj3zzDOD+rlxEqe4yxPMhNIfiPK9Hg+M118TyorTLc+thFJg06ZNmpSUpIBecsklumHDhoDWd/DgQb3lllsU0LPPPlvXr18f0Poqkz/+8Y9ao0YNPXbsmNuh6JEjR7RatWp62223uR1KmZ0qoQStyaOqX6l3qEeA+UCTYNUdCAkJCXz99de88sorpKWl0bFjR/79738HpG9l+vTpnHPOObz66qvcddddLFu2LKgzy1V0SUlJZGZmsmDBArdDYcaMGWRlZVXO5g6482Ab8CnwO/31CiUTWAJ8D/QqTRluX6EUtmXLFk1OTlZAe/XqpWvXrnWk3L179+rIkSMV0HPOOUcXLFjgSLnhZt++fSoi+tBDD7kdil5zzTUaHx+vubm5bodSZgSryQPMBFb4WYYW2mcs3j4U8f0eA9Tzve6Kd8KvmsWUPwpIA9KaNWsW4NN2ejwej7755ptaq1YtrVq1qj799NPlGpckNTVVGzRooFFRUfrAAw+E5fiwTkpMTNSePXu6GkN2drbGxcXpzTff7Goc5RW0hFLSAlwPzAOqn2Kf74DEksoKpSuUwrZt26ZDhgxRQC+88EJdsmTJaSWWnTt36m9+8xsF9LzzztOlS5cGMNrwcc8992hUVJQePnzYtRg+/fRTBULyaxynIyQSCpAMrALqF1lfH4j0vW4JbAPqllReqCYUVe/Vyrvvvqt169ZVQKOiojQhIUEvvvhiHTlypN5///362muv6cyZM3XdunWanZ2tHo9H33vvPa1Xr55WqVJFH3vssZDoRKwsvv76awVcnejthhtu0Fq1amlOTo5rMTjhVAklmHMbv4C3efO172Ge+ao6GugNPCIieUA+MFpV9wUxLseJCNdeey19+/Zl6tSppKenk56ezubNm5k1axbbt28/ofNWRKhXrx579uyhe/fuvP7667Rr187Fv6Dy6dGjBzExMXzzzTd+57MOtNzcXKZOncqQIUOoUqVK0OsPlqAlFFX1O4KMqqYCqcGKI5gaNGjAqFGjTlqfm5tLRkbG8USTnp7Oli1b6NKlC7fcckvYD+EYCNWqVaNHjx6ufa/nu+++Y//+/ZX37o5PMK9QjE/BBObhMol5qEhKSmLs2LHs2rWLM844I6h1p6amUqNGDQYMGBDUeoPNHr03YaNv374AfPvtt0GtNz8/n48//phBgwZV+Im8SmIJxYSNrl27UqtWraB/r2fu3Lns2rWr0jd3wBKKCSORkZH06dMn6P0oqampxMTEuNIZHGyWUExY6du3L5s2bWLjxo1BqU9VmTJlCgMGDCAuLi4odbrJEooJKwXDGQTrKmXhwoVs3bq10o3MVhxLKCastG3bljPPPDNoCSU1NZWoqCguu+yyoNTnNksoJqwEc1hIVSU1NZVLL72UOnXqBLSuUGEJxYSdvn37smfPHpYvXx7QepYtW8aGDRvC4u5OAUsoJuwEa1jIjz76iIiICC6//PKA1hNKLKGYsNO4cWPatGkT0H4UVWXy5MlccsklQX8q102WUExY6tu3L99//z3Hjh0LSPnLli1j7dq1jBgxIiDlhypLKCYsJSUlcfTo0YANCzl58mQiIyO54oorAlJ+qLKEYsJSnz59iIiICEizp6C5c+mll1K/fn3Hyw9lllBMWKpTpw5du3YNSMfs0qVLWb9+Pb/5zW8cLzvUBS2hiMhDIrJNRJb6lkGFtv1dRNaLyBoRqdzf7zYho2/fvixYsIDDhw87Wm64Nncg+Fcoz6pqZ9/yBYCItAeuBs7BO0zkiyJiIwyZgEtKSiIvL48ffvjBsTILmjtJSUnEx8c7Vm5FEQpNnqHA+6qao6qbgPVAN5djMmHgoosuOj4spFMWL17Mxo0bw+7uToFgJ5TbRWSZiLwhIgXPIjfGO3VGgQzfOmMCqlq1avTs2dPRfpTJkycTFRUVVg+zFeZoQhGRmSKyws8yFHgJOAvoDOwAni44zE9RWkz5o0QkTUTSdu/e7WToJkwlJSWxfPlydu3aVe6yCpo7ffv2pV69eg5EV/E4mlBUta+qdvCzTFXVnaqar6oe4FV+bdZkAE0LFdME2F5M+RNUNVFVE8PtdpwJjIJhIWfNmlXustLS0ti8eXPYNncguHd5GhX69Qq8MwoCTAOuFpEYEWkBtAJ+ClZcJrx16dKF2rVrO9Ls+fDDD4mOjg7b5g4Ed9T7J0WkM97mzGbgFgBVXSkik/FOApYH3Kaq+UGMy4SxyMhILrnkEmbOnImq4psz6rQVNHf69esXNkMV+BO0KxRVvU5VO6pqJ1W9TFV3FNr2f6p6lqq2UdXpwYrJGPD2o6Snp5drWMiFCxeSnp4e1s0dCI3bxsa4yolhISdPnkx0dDRDhw51KqwKyRKKCXtt2rShcePGfPXVV2U6vqC5079/f2rXru1scBWMJRQT9kSEESNGMGXKlDLd7VmwYAFbt24N++YOWEIxBoBx48bRunVrrrvuOvbs2XNax06ePJkqVaqEzUDUp2IJxRigRo0aTJo0iT179nDzzTej6vfZypN4PB4+/PBDBgwYEPbNHbCEYsxx5513Hk888QRTp07l5ZdfLtUx8+fPJyMjw5o7PpZQjCnkjjvuIDk5mb/85S+sWLGixP0nT55MTEyMNXd8LKEYU0hERARvvfUWNWvW5JprriErK6vYfQuaO8nJydSsWTOIUYYuSyjGFNGgQQMmTpzIihUruPvuu4vdb+7cuWzfvt2aO4VYQjHGj+TkZO68805eeOEFPv30U7/7FDR3hgwZEuToQpclFGOK8fjjj9O5c2duvPFGtm8/8Qvw+fn5fPTRRwwaNIi4uDiXIgw9llCMKUZMTAyTJk0iKyuL66+//oS5kH/88Ud27NhhzZ0iLKEYcwpt27blX//6FzNnzuSpp546vn7y5MlUrVqVwYMHuxhd6LGEYkwJbrrpJoYPH87YsWNZuHAh+fn5pKamMmjQIGJjY90OL6RYQjGmBCLChAkTaNSoEb/97W+ZPn06v/zyizV3/JDSPmIcahITEzUtLc3tMEwYmT17Nn369KFatWp4PB527doVllcoIrJIVRP9bQvmEJAfFJrka7OILPWtTxCRrELbSvfMszFB1qtXL+677z4yMzNJSUkJy2RSkqANAamqVxW8FpGngYOFNm9Q1c7BisWYsrr//vvJzs7mqquuKnnnMBTMMWUBEO+gnSOAS4NdtzHlFRUVxfjx490OI2S50SnbC9ipqusKrWshIktE5HsR6eVCTMYYBzh6hSIiM4GGfjaNVdWpvtfXAJMKbdsBNFPVvSLSFfhERM5R1UN+yh8FjAJo1qyZk6EbYxzgaEJR1b6n2i4iUcAwoGuhY3KAHN/rRSKyAWgNnHQLR1UnABPAe5fHuciNMU4IdpOnL/CzqmYUrBCR+iIS6XvdEu9EX2Wfz8AY45pgd8pezYnNHYDewCMikgfkA6NVdV9JBS1atGiPiKSXos544PQGCXWe2zFY/fYZcLL+5sVtqLAPtpWWiKQV9xBOuMRg9dtnIFj126P3xhjHWEIxxjgmHBLKBLcDwP0YrH73uR1DUOqv9H0oxpjgCYcrFGNMkFSahCIiySKyRkTWi8i9fraLiDzv275MRLo4WHdTEflWRFaLyEoRucPPPn1E5GChb1U/4FT9herYLCLLfeWf9GBggM9Bm0J/21IROSQiY4rs4+g5EJE3RGSXiKwotK6uiHwtIut8P+sUc+wpPy/ljOGfIvKz7xx/LCK1izn2lO9XOep/SES2FTrPg4o51pFzcAJVrfALEAlsAFoCVYD/Ae2L7DMImA4IcCGwwMH6GwFdfK/jgLV+6u8DfBbg87AZiD/F9oCdAz/vxy9A80CeA7zPMHUBVhRa9yRwr+/1vcD4snxeyhlDfyDK93q8vxhK836Vo/6HgL+V4j1y5BwUXirLFUo3YL2qblTVY8D7wNAi+wwF3lav+UBtEWnkROWqukNVF/teHwZWA42dKNthATsHRSThHZKiNA8elpmq/gAUfQhyKDDR93oicLmfQ0vzeSlzDKr6larm+X6dDzQpS9llrb+UHDsHhVWWhNIY2Fro9wxO/gddmn3KTUQSgPOABX42dxeR/4nIdBE5x+m6AQW+EpFFvi9SFhWUc4D/J6ILBPocNFDVHeBN9MAZfvYJ1nkA+D3eq0J/Snq/yuN2X5PrjWKafQE5B5UloYifdUVvX5Vmn/IFIRILpAJj9ORvSy/G2wQ4F/g38ImTdfv0UNUuwEDgNhHpXTREP8c4fQ6qAJcBH/rZHIxzUBoBPw8AIjIWyAPeK2aXkt6vsnoJOAvojPfb/E/7C8/PunKfg8qSUDKApoV+bwJsL8M+ZSYi0XiTyXuqOqXodlU9pKpHfK+/AKJFJN6p+n3lbvf93AV8jPeytrCAngOfgcBiVd3pJ76AnwNgZ0Ezzvdzl599An4eROR6YDBwrfo6LYoqxftVJqq6U1XzVdUDvFpMuQE5B5UloSwEWolIC9//kFcD04rsMw0Y6bvTcSFwsODSuLxERIDXgdWq+kwx+zT07YeIdMN77vc6Ub+vzBoiElfwGm/H4IoiuwXsHBRSdLybwjEG9Bz4TAOu972+HpjqZ5/SfF7KTESSgXuAy1T1aDH7lOb9Kmv9hfvFriim3MCcg/L26obKgvcOxlq8PddjfetG4/32Mngv8f7j274cSHSw7p54LxeXAUt9y6Ai9d8OrMTbmz4fuMjhv7+lr+z/+eoJ6jnwlV8db4KoVWhdwM4B3sS1A8jF+z/uTUA94Btgne9nXd++ZwJfnOrz4mAM6/H2TxR8Fl4uGkNx75dD9b/je3+X4U0SjQJ5Dgov9qSsMcYxlaXJY4wJAZZQjDGOsYRijHGMJRRjjGMsoRhjHGMJxRjjGEsoxhjHWEIxxjgm6JOlOyU+Pl4TEhLcDsOYsLNo0aI9qlrf37YKm1ASEhJISyvTIFfGmHKQU0ywZ00eY4xjLKEYYxxjCcUY45gK24fiT25uLhkZGWRnZ5+0rWrVqjRp0oTo6GgXIjMm+D5btp1nvlpLTp6nVPs3qBnDlD/2KFedlSqhZGRkEBcXR0JCAr5xfADvmC979+4lIyODFi1auBihMYGX71GenPEzr/ywkQ6Na9Klud+ZRE5Sp3r5/7OtVAklOzv7pGQCICLUq1eP3bt3uxSZMcGxP/MYf35/CbPX7eF3FzbjgcHnUCUqeD0bJdZ0uhMJicjffRMHrRGRAYXWd/VNarRevJNNFQwFGCMiH/jWL/CNGl9mRZNJSeuNqSxWbT/EZf+Zw4KN+xh/ZUcevbxjUJMJlK5T9i0g2c/6Z1W1s2/5AkBE2uMdm/Ic3zEvikikb/+XgFFAK99SUOZNwH5VPRt4Fu/ESMaY0zDtf9sZ9tKP5OYpH9xyIVed38yVOEpMKHp6EwkNBd5X1RxV3YR3bM1uvkFza6rqPPWOOfk2v07AVHhipo+AJLHLCWNKJS/fw2NfrObPk5bQsXEtpv2pB+c1K12fSSCU53rI30RCxU0e1Nj3uuj6E45R72xrB/EONFwmxY2Ra2PnmspmX+Yxrn/zJyb8sJGR3Zvz3s0XckZcVVdjKmtCKW4ioeImDzrVpEKlnnBIREaJSJqIpPnrYK1atSp79+49KXkU3OWpWtXdk22MU1ZuP8iQf89h4eb9PDm8E48M7RD0/hJ/ynSXRwtN4iQirwKf+X4tbvKgDE6c37XwpEIFx2SISBRQi2KaWKo6AZgAkJiYeFLSadKkCRkZGX7v5hQ8h2JMRTd16TbuSV1GnepV+PCW7pzbtLbbIR1XpoQiIo301wmiCk8kNA34r4g8g3cOkFbAT6qaLyKHfZNLLQBG4p2KsuCY64F5wHBglpaxfRIdHW3PmZhKS1V5/pv1PDtzLd1a1OU/v+1C/bgYt8M6QYkJRUQmAX2AeBHJAB4E+ohIZ7xNk83ALQCqulJEJgOr8M7pepuq5vuKuhXvHaNqeCePLphA+nXgHRFZj/fK5GoH/i5jKpXcfA9/n7KcjxZlcGWXJjw+LPi3hEujwk70lZiYqDZ8gQkHh7Nz+eN7i5m9bg93JLViTN9Wrj5XJSKLVDXR37ZK9aSsMZXNjoNZ3PjmQtbvOsKTwzsxIrFpyQe5yBKKMSFq1fZD/P6thRzJyePNG8+nVyu/g6SFFEsoxoSg2et2c+u7i4mNieLD0d1p16im2yGViiUUY0LM5LSt/GPKcs4+I5Y3bzyfRrWquR1SqVlCMSZEqCrPzVzHv75ZR69W8bx4bRfiqlas8XssoRjjMlUlJ8/D2I9XkLo4g990bcJjwzoSHRl6t4VLYgnFmAD5YvkO/jVzHUdz88jPV/I8viXfQ77vdcHPAn/p15o/XXp2hR1uwxKKMQ7zeJTnZq7l+VnradeoJuc3r0tkhBAVKURFRHhfRwiRkd6fURERREUIHZrU4pI2Z7gdfrlYQjHGQUdy8vjLB0v5atVORiQ2YdzlHYiJiiz5wErCEooxDtmy9yh/eDuN9buP8OCQ9txw0cnDkVZ2llCMccDcDXv443uLUYWJN3ajZ6t4t0NyhSUUY8pBVXlnfjoPf7qKlvE1eHVkIgnxNdwOyzWWUIwpo2N5Hh6ctpJJP22hb7szePaqzhXuuRGnWUIxpgz2HMnhj+8u5qfN+7jtkrP4a782RESEV3+JP5ZQjDlNK7cfZNTbi9ibmcPz15zHZeee6XZIIcMSijGn4VB2LldPmE9sTBQfjb6IDo1ruR1SSLGEYsxpmL58B4ez85j4+26WTPyoeF8WMMZFUxZvo2V8Dc4LoYGhQ4klFGNKaeu+oyzYtI8rzmscdg+slZYlFGNKaerSbQBcfl7jEvYMX5ZQjCkFVWXK4m10a1GXpnWrux1OyLKEYkwpLN16gI17Mrmyi12dnIolFGNK4eMl24iJimBgx0ZuhxLSLKEYU4JjeR6m/W87/do3oGaYP1pfEksoxpTguzW7OHA0lyu72NzYJbGEYkwJpizeRnxsFXqF6ZAEp6PEhCIib4jILhFZUWhdXRH5WkTW+X7WKbTt7yKyXkTWiMiAQuu7ishy37bnxXcjX0RiROQD3/oFIpLg8N9oTJkdOHqMWT/v4rJzGxNVAQeNDrbSnKG3gOQi6+4FvlHVVsA3vt8RkfZ4Jzs/x3fMiyJSMP7dS8AooJVvKSjzJmC/qp4NPAuML+sfY4zTPlu2g2P5HobZ3Z1SKTGhqOoPwL4iq4cCE32vJwKXF1r/vqrmqOomYD3QTUQaATVVdZ56Z2d/u8gxBWV9BCSJPYZoQsSUxRm0bhDLOWdWjJn73FbWa7gGqroDwPezYKjuxsDWQvtl+NY19r0uuv6EY1Q1DzgI1CtjXMY4ZvOeTBZvOcCwLk3sUftScrpR6O+s6ynWn+qYkwsXGSUiaSKStnv37jKGaEzpTFmyDREY2tnGOymtsiaUnb5mDL6fu3zrM4CmhfZrAmz3rW/iZ/0Jx4hIFFCLk5tYAKjqBFVNVNXE+vVDfyZ6U3GpKh8vyaDHWfEVam5ht5U1oUwDrve9vh6YWmj91b47Ny3wdr7+5GsWHRaRC339IyOLHFNQ1nBglq+fxRjXpKXvZ+u+LK6wLwKelhIHWBKRSUAfIF5EMoAHgSeAySJyE7AF+A2Aqq4UkcnAKiAPuE1V831F3Yr3jlE1YLpvAXgdeEdE1uO9Mrnakb/MmHKYsngb1aIjSe7Q0O1QKpQSE4qqXlPMpqRi9v8/4P/8rE8DOvhZn40vIRkTCrJz8/ls2XaSOzSkRowNang67EkdY4r4ZvUuDmfn2bMnZWAJxZgiPl6SQYOaMVx0lj1qf7osoRhTyN4jOXy3ZjeXd25MpM2zc9osoRhTyKf/206eRxlm3ywuE0soxhQyZck22jeqSZuGcW6HUiFZQjHGZ/2uwyzLOGidseVgCcUYnymLtxEZIVxmj9qXmSUUYwCPR/lkyTZ6tYrnjLiqbodTYVlCMQaYv2kv2w9mW2dsOVlCMQb4YOFWYmOi6N++gduhVGiWUEzYW7r1AFOXbufaC5pRNTqy5ANMsSyhmLDm8SgPTVtJ/bgYbr/0bLfDqfAsoZiwlro4g6VbD3BvclvibM6dcrOEYsLWoexcxs9YQ5dmtW3cE4fYd7NN2Pr3N+vYm5nDGzckEmHf23GEXaGYsLR+1xHe/HEzVyU2pVOT2m6HU2lYQjFhR1V5+NOVVKsSyd8GtHE7nErFEooJO1+v2snsdXu4s29r4mNj3A6nUrGEYsJKdm4+4z5fResGsVzXvbnb4VQ61ilrwsprszeydV8W7918AdE2V7Hj7IyasLH9QBb/+XYDAzs0pMfZNrxjIFhCMWHjsS9W41HlH4PauR1KpWUJxYSF+Rv38tmyHdza5yya1q3udjiVliUUU+nl5Xt4aNpKGteuxuiLz3I7nErNEoqp9P770xZ+/uUw96W0s28TB5glFFOp7cs8xtNfreWis+rZtKJBYAnFVGpPf7WGIzl5PHTZOYjY93UCrVzPoYjIZuAwkA/kqWqiiNQFPgASgM3ACFXd79v/78BNvv3/rKpf+tZ35deJ1L8A7lBVLU9sJrzke5Sdh7LJ2J/FtgNH2bY/i637spi8aCs3XJRA6wY2LUYwOPFg2yWquqfQ7/cC36jqEyJyr+/3e0SkPXA1cA5wJjBTRFqraj7wEjAKmI83oSQD0x2IzVQy2w9k8dOmfWzak8m2A1lk7D/KtgNZ7DiQTZ7nxP+D4mOr0LtVfcb0be1StOEnEE/KDgX6+F5PBL4D7vGtf19Vc4BNIrIe6Oa7yqmpqvMARORt4HIsoRhg16Fs5m3cy7wNe5m3cS/pe48CIAIN4qrSuE41zmtahyGdqtG4TjWa1KlO49rVaFy7GtWqWAdssJU3oSjwlYgo8IqqTgAaqOoOAFXdISJn+PZtjPcKpECGb12u73XR9SYM7T2Sw/yN+5i3cQ/zNuxlw+5MAOKqRnFBi3qM7J5A95b1OPuMWKpEWRdgqClvQumhqtt9SeNrEfn5FPv66xHTU6w/uQCRUXibRjRr1ux0YzUh7J356bw7L501Ow8DUKNKJN1a1OWq85vSvWU87c+saZOXVwDlSiiqut33c5eIfAx0A3aKSCPf1UkjYJdv9wygaaHDmwDbfeub+Fnvr74JwASAxMRE67StJP41cx3PzlxLl2a1uWtAG7qfVY+OjWvZl/cqoDK/YyJSQ0TiCl4D/YEVwDTget9u1wNTfa+nAVeLSIyItABaAT/5mkeHReRC8d7XG1noGFOJqSrPfLWGZ2eu5couTfhw9EXcdsnZdGlWx5JJBVWeK5QGwMe+e/tRwH9VdYaILAQmi8hNwBbgNwCqulJEJgOrgDzgNt8dHoBb+fW28XSsQ7bSU1We+moN//l2A1clNuXxYR1tXNdKQCrq4x6JiYmalpbmdhimDFSVJ2b8zCvfb+Sabs34v8s7WDKpQERkkaom+ttmAyyZoFJVHv18Na/P2cR1Fzbn4cvOsWRSiVhCMUHjHRx6FW/N3cyNPRJ4YHB7exy+krGEYoLC41EemLaCd+dv4eaeLRib0s6SSSVkCcUEnMejjP1kBZN+2sItF7fk3uS2lkwqKUsoJqA8HuXeKcuYnJbB7ZeczV/7t7ZkUolZQjEBk+9R7v5oGamLM7gjqRVj+rayZFLJWUIxjlNVvl2ziydnrOHnXw7zl36t+XNSK7fDMkFgCcU4alH6PsZPX8NPm/fRvF51XvjteQzudKbbYZkgsYRiHLHml8P888s1zFy9k/pxMTx6eQeuOr+pPUIfZiyhmHLJ2H+UZ75ey8dLthEbE8VdA9pwY48Eqlexj1Y4snfdlMneIzm88O163pu/BREY1asloy8+izo1qrgdmnGRJRRTarsP57B6xyHmb9zLxLmbycrNZ0RiU+7o24pGtaq5HZ4JAZZQzEnyPcrmvZms2n6IVTsOHf+5+3DO8X0GdmjIX/u34ewzYl2M1IQaSyhhzONRdhzKZuPuI2zak8nPvxxm9Y5D/LzjMFm53pEloiOFs8+Io3er+rQ/sybtG9WkXaM4ale3po05mSWUMHDwaC4b9hxh0+5MNu3JZOOeI2zcncnmvZlk53qO71erWjTtGsVxTbdmx5OHjd1qTocllAoiY/9R5m/cx/yNe1myZf8JieBUjh7LY//R3OO/R0YIzepWp2V8DXq1iqdFfCwt69egZXwN6sfF2JOsplwsoYSorfuOMn/jXhZs8iaRjP1ZANSpHk1iQl1qVo0uVTkx0REk1KtOy/hYWtSvQbO61e3ZEBMwllBCxM5D2fywdvfxq5BtB35NIBe2rMfNPVtw4Vn1aH1GnA1IZEKWJRSXLd16gNfnbOKL5TvI9yh1a1ThghZ1GdW7JRe2rEerM2ItgZgKwxKKC/LyPXy1aievz9nEovT9xMVEceNFCQxPbGJXIKZCs4QSRIeyc5m8cCtv/riZbQeyaFq3Gg8Mbs+I85sSG2Nvhan47FMcBFv2HuXNuZv4MC2DIzl5dEuoy/2D29OvfQObDc9UKpZQAkRVmbdhLxPnbebrVTuJEGFwp0bc1LMlHZvUcjs8YwLCEorD9mUe46NFW5n001Y27cmkdvVoRl98FiO7J9CwVlW3wzMmoCyhOEBVmb9xH//9aQtfrviFY/keEpvX4U+Xns2gjo2oGh3pdojGBIUllHLYn3mM1MUZ/PenLWzcnUnNqlH89oJm/PaCZrRuEOd2eMYEXcgkFBFJBv4FRAKvqeoTLofkV05ePovTD/DBwi18seIXjuV56Nq8Dk/95mxSOjaiWhW7GjHhKyQSiohEAv8B+gEZwEIRmaaqq9yMS1XZdiCLJVsOeJet+1m57RDH8j3EVY3imvObcs0FzWjbsKabYRoTMkIioQDdgPWquhFARN4HhgJBTShHj+WxLOOgL4HsZ8nWA8fHAKkaHUGnxrW5sUcC5zWrzcWtz7CrEWOKCJWE0hjYWuj3DOCC8hb69aqdzF63m2N5HnLyPL6f+eQU+r3wuu0HsvCo99gW8TXodXY85zWrzXnN6tCmYZx9qc6YEoRKQvH3dJeetJPIKGAUQLNmzUosdMW2g0z733aqREYQEx1BlcgIqkRFEhMVQZWoCOKqRhETFUFMVCRVoiJoUqcaXZrVoXPT2jY2qjFlECoJJQNoWuj3JsD2ojup6gRgAkBiYuJJCaeoO/u15s5+rZ2K0RhTglC5hl8ItBKRFiJSBbgamOZyTMaY0xQSVyiqmicitwNf4r1t/IaqrnQ5LGPMaQqJhAKgql8AX7gdhzGm7EKlyWOMqQREtcS+zZAkIruB9FLsGg/sCXA4oR6D1W+fASfrb66q9f1tqLAJpbREJE1VE8M5BqvfPgPBqt+aPMYYx1hCMcY4JhwSygS3A8D9GKx+97kdQ1Dqr/R9KMaY4AmHKxRjTJBUmoQiIskiskZE1ovIvX62i4g879u+TES6OFh3UxH5VkRWi8hKEbnDzz59ROSgiCz1LQ84VX+hOjaLyHJf+Wl+tgfyHLQp9LctFZFDIjKmyD6OngMReUNEdonIikLr6orI1yKyzvezTjHHnvLzUs4Y/ikiP/vO8cciUruYY0/5fpWj/odEZFuh8zyomGMdOQcnUNUKv+B9XH8D0BKoAvwPaF9kn0HAdLzfbL4QWOBg/Y2ALr7XccBaP/X3AT4L8HnYDMSfYnvAzoGf9+MXvM8rBOwcAL2BLsCKQuueBO71vb4XGF+Wz0s5Y+gPRPlej/cXQ2ner3LU/xDwt1K8R46cg8JLZblCOT5Ak6oeAwoGaCpsKPC2es0HaotIIycqV9UdqrrY9/owsBrvGC+hJmDnoIgkYIOqlubBwzJT1R+AfUVWDwUm+l5PBC73c2hpPi9ljkFVv1LVPN+v8/F+ez4gijkHpeHYOSissiQUfwM0Ff0HXZp9yk1EEoDzgAV+NncXkf+JyHQROcfpuvGOIfOViCzyjR1TVFDOAd5vi08qZlugz0EDVd0B3kQPnOFnn2CdB4Df470q9Kek96s8bvc1ud4optkXkHNQWRJKaQZoKtUgTuUKQiQWSAXGqOqhIpsX420CnAv8G/jEybp9eqhqF2AgcJuI9C4aop9jnD4HVYDLgA/9bA7GOSiNgJ8HABEZC+QB7xWzS0nvV1m9BJwFdAZ2AE/7C8/PunKfg8qSUEozQFOpBnEqKxGJxptM3lPVKUW3q+ohVT3ie/0FEC0i8U7V7yt3u+/nLuBjvJe1hQX0HPgMBBar6k4/8QX8HAA7C5pxvp+7/OwT8PMgItcDg4Fr1ddpUVQp3q8yUdWdqpqvqh7g1WLKDcg5qCwJpTQDNE0DRvrudFwIHCy4NC4vERHgdWC1qj5TzD4NffshIt3wnvu9TtTvK7OGiMQVvMbbMbiiyG4BOweFXEMxzZ1AnwOfacD1vtfXA1P97BPQAb3EOyXMPcBlqnq0mH1K836Vtf7C/WJXFFNuYM5BeXt1Q2XBewdjLd6e67G+daOB0b7Xgneqjg3AciDRwbp74r1cXAYs9S2DitR/O7ASb2/6fOAih//+lr6y/+erJ6jnwFd+dbwJolahdQE7B3gT1w4gF+//uDcB9YBvgHW+n3V9+54JfHGqz4uDMazH2z9R8Fl4uWgMxb1fDtX/ju/9XYY3STQK5DkovNiTssYYx1SWJo8xJgRYQjHGOMYSijHGMZZQjDGOsYRijHGMJRRjjGMsoRhjHGMJxRjjmP8HA+D93XoP+j4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 288x360 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "\n",
    "proj0 = qutip.ket2dm(qutip.tensor(qutip.ket(\"0\"),qutip.ket(\"0\")))\n",
    "proj1 = qutip.ket2dm(qutip.tensor(qutip.ket(\"0\"),qutip.ket(\"1\")))\n",
    "proj2 = qutip.ket2dm(qutip.tensor(qutip.ket(\"1\"),qutip.ket(\"0\")))\n",
    "proj3 = qutip.ket2dm(qutip.tensor(qutip.ket(\"1\"),qutip.ket(\"1\")))\n",
    "omega=1\n",
    "g=0.2*omega\n",
    "ampl0=0.00001\n",
    "#T=8\n",
    "#ampl0=1\n",
    "T=math.pi/g\n",
    "nt=20\n",
    "tlist = np.linspace(0,T, nt)\n",
    "def hamiltonian(omega, ampl0, g):\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",
    "      HA = -1j*qutip.Qobj(H0(0.2))\n",
    "      H1 =-1j*qutip.Qobj(Hde(0.05,0.2))\n",
    "      H2= -1j*qutip.Qobj(He(0.05,0.2))\n",
    "      def guess_control(t, args):\n",
    "          return ampl0*krotov.shapes.flattop(\n",
    "              t, t_start=0, t_stop=T, t_rise=0.005, func=\"blackman\"\n",
    "          )\n",
    "      def guess_control2(t, args):\n",
    "          return ampl0*krotov.shapes.flattop(\n",
    "              t, t_start=0, t_stop=T, t_rise=0.005, func=\"blackman\"\n",
    "          )\n",
    "      \n",
    "      return [HA, [H1, guess_control]]\n",
    "\n",
    "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",
    "    )\n",
    "def plot_iterations(opt_result):\n",
    "    \"\"\"Plot the control fields in population dynamics over all iterations.\n",
    "\n",
    "    This depends on ``store_all_pulses=True`` in the call to\n",
    "    `optimize_pulses`.\n",
    "    \"\"\"\n",
    "\n",
    "    fig, [ax_ctr,ax] = plt.subplots(nrows=2, figsize=(4, 5))\n",
    "    n_iters = len(opt_result.iters)\n",
    "    EEnergy=np.zeros(nt)\n",
    "    for (iteration, pulses) in zip(opt_result.iters, opt_result.all_pulses):\n",
    "        controls = [\n",
    "            krotov.conversions.pulse_onto_tlist(pulse)\n",
    "            for pulse in pulses\n",
    "        ]\n",
    "        objectives = opt_result.objectives_with_controls(controls)\n",
    "        dynamics = objectives[0].mesolve(\n",
    "            opt_result.tlist, e_ops=[]\n",
    "        )\n",
    "        if iteration == 0:\n",
    "            ls = '--'  # dashed\n",
    "            alpha = 1  # full opacity\n",
    "            ctr_label = 'guess'\n",
    "            pop_labels = ['0 (guess)', '1 (guess)']\n",
    "        elif iteration == opt_result.iters[-1]:\n",
    "            ls = '-'  # solid\n",
    "            alpha = 1  # full opacity\n",
    "            ctr_label = 'optimized'\n",
    "            pop_labels = ['0 (optimized)', '1 (optimized)']\n",
    "        else:\n",
    "            ls = '-'  # solid\n",
    "            alpha = 0.5 * float(iteration) / float(n_iters)  # max 50%\n",
    "            ctr_label = None\n",
    "            pop_labels = [None, None]\n",
    "        ax_ctr.plot(\n",
    "            dynamics.times,\n",
    "            controls[0],\n",
    "            label=ctr_label,\n",
    "            color='black',\n",
    "            ls=ls,\n",
    "            alpha=alpha,\n",
    "        )\n",
    "    EField=np.transpose(np.array(opt_result.optimized_controls))\n",
    "    EEnergy[0]=(np.square(EField[0]))*(T/nt)\n",
    "    a=0\n",
    "    for i in range (1,nt):\n",
    "      a+=np.square(EField[i-1])\n",
    "      EEnergy[i]=(np.square(EField[i])+a)*(T/nt)\n",
    "      \n",
    "    \n",
    "    ax.plot(tlist,np.transpose(EEnergy))\n",
    "    plt.legend()\n",
    "    plt.show(fig)\n",
    "    \n",
    "    \n",
    "    \n",
    "\n",
    "\n",
    "\n",
    "H = hamiltonian(omega,ampl0,g)\n",
    "pulse_options = {\n",
    "    H[1][1]: dict(lambda_a=0.01, update_shape=S),   \n",
    "}\n",
    "objectives = [\n",
    "    krotov.Objective(\n",
    "        initial_state=qutip.Qobj(np.array([1,0,0,0,0,1/2,0,0,0,1/2,0,0,1/2,0,1/2])), target=qutip.Qobj(qutip.Qobj(np.array([1,0,0,2*np.sqrt(2),0,1/2,0,0,0,1/2,0,0,1/2,0,1/2]))), H=H\n",
    "      )\n",
    "  ]\n",
    "\n",
    "opt_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_ss,\n",
    "  info_hook=krotov.info_hooks.print_table(J_T=krotov.functionals.J_T_ss),\n",
    "  check_convergence=krotov.convergence.Or(\n",
    "      krotov.convergence.value_below('5e-3', name='J_T'),\n",
    "      krotov.convergence.check_monotonic_error,\n",
    "\n",
    "  ),\n",
    "    store_all_pulses=True,\n",
    ")\n",
    "plot_iterations(opt_result)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[Quantum object: dims = [[15], [1]], shape = (15, 1), type = ket\n",
       " Qobj data =\n",
       " [[ 1.00000000e+00]\n",
       "  [-1.09546587e+01]\n",
       "  [ 1.40097459e+01]\n",
       "  [ 1.11161768e+01]\n",
       "  [-2.47189390e+00]\n",
       "  [ 1.51497706e+02]\n",
       "  [-6.70922040e+01]\n",
       "  [ 5.68433539e+01]\n",
       "  [-3.44036440e+01]\n",
       "  [-1.35201769e+00]\n",
       "  [-6.17848523e+01]\n",
       "  [ 3.51302356e+01]\n",
       "  [ 4.43538888e-01]\n",
       "  [-2.02883894e-02]\n",
       "  [ 5.64576467e-01]]]"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "opt_result.states"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}