Skip to content

Instantly share code, notes, and snippets.

@goerz

goerz/06_example_3states_benchmark.ipynb

Last active June 13, 2021 17:57
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save goerz/34af142b78d7e344417d838bbea78aaf to your computer and use it in GitHub Desktop.
Save goerz/34af142b78d7e344417d838bbea78aaf to your computer and use it in GitHub Desktop.
Download ZIP
Benchmarking of propagators. Notebook derived from https://qucontrol.github.io/krotov/v1.2.1/notebooks/06_example_3states.html
Display the source blob
Display the rendered blob
Raw
06_example_3states_benchmark.ipynb
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Optimization of a Dissipative Quantum Gate - Propagator Benchmarks"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.304910Z",
"start_time": "2021-06-13T15:05:12.110832Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:56.348706Z",
"iopub.status.busy": "2021-01-13T19:03:56.347322Z",
"iopub.status.idle": "2021-01-13T19:03:57.862517Z",
"shell.execute_reply": "2021-01-13T19:03:57.862988Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Python implementation: CPython\n",
"Python version : 3.8.1\n",
"IPython version : 7.24.1\n",
"\n",
"sys : 3.8.1 (default, Aug 12 2020, 19:33:59) \n",
"[GCC 8.3.0]\n",
"scipy : 1.6.3\n",
"krotov : 1.2.1+dev\n",
"qutip : 4.6.1\n",
"numpy : 1.20.3\n",
"matplotlib: 3.4.2\n",
"\n"
]
}
],
"source": [
"# NBVAL_IGNORE_OUTPUT\n",
"%load_ext watermark\n",
"import sys\n",
"import os\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 copy\n",
"from functools import partial\n",
"from itertools import product\n",
"%watermark -v --iversions"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"$\\newcommand{tr}[0]{\\operatorname{tr}}\n",
"\\newcommand{diag}[0]{\\operatorname{diag}}\n",
"\\newcommand{abs}[0]{\\operatorname{abs}}\n",
"\\newcommand{pop}[0]{\\operatorname{pop}}\n",
"\\newcommand{aux}[0]{\\text{aux}}\n",
"\\newcommand{int}[0]{\\text{int}}\n",
"\\newcommand{opt}[0]{\\text{opt}}\n",
"\\newcommand{tgt}[0]{\\text{tgt}}\n",
"\\newcommand{init}[0]{\\text{init}}\n",
"\\newcommand{lab}[0]{\\text{lab}}\n",
"\\newcommand{rwa}[0]{\\text{rwa}}\n",
"\\newcommand{bra}[1]{\\langle#1\\vert}\n",
"\\newcommand{ket}[1]{\\vert#1\\rangle}\n",
"\\newcommand{Bra}[1]{\\left\\langle#1\\right\\vert}\n",
"\\newcommand{Ket}[1]{\\left\\vert#1\\right\\rangle}\n",
"\\newcommand{Braket}[2]{\\left\\langle #1\\vphantom{#2}\\mid{#2}\\vphantom{#1}\\right\\rangle}\n",
"\\newcommand{ketbra}[2]{\\vert#1\\rangle\\!\\langle#2\\vert}\n",
"\\newcommand{op}[1]{\\hat{#1}}\n",
"\\newcommand{Op}[1]{\\hat{#1}}\n",
"\\newcommand{dd}[0]{\\,\\text{d}}\n",
"\\newcommand{Liouville}[0]{\\mathcal{L}}\n",
"\\newcommand{DynMap}[0]{\\mathcal{E}}\n",
"\\newcommand{identity}[0]{\\mathbf{1}}\n",
"\\newcommand{Norm}[1]{\\lVert#1\\rVert}\n",
"\\newcommand{Abs}[1]{\\left\\vert#1\\right\\vert}\n",
"\\newcommand{avg}[1]{\\langle#1\\rangle}\n",
"\\newcommand{Avg}[1]{\\left\\langle#1\\right\\rangle}\n",
"\\newcommand{AbsSq}[1]{\\left\\vert#1\\right\\vert^2}\n",
"\\newcommand{Re}[0]{\\operatorname{Re}}\n",
"\\newcommand{Im}[0]{\\operatorname{Im}}$\n",
"\n",
"This example illustrates the optimization for a quantum gate in an open quantum system, where the dynamics is governed by the Liouville-von Neumann equation. A naive extension of a gate optimization to Liouville space would seem to imply that it is necessary to optimize over the full basis of Liouville space (16 matrices, for a two-qubit gate). However, [Goerz et al., New J. Phys. 16, 055012 (2014)][1] showed that is not necessary, but that a set of 3 density matrices is sufficient to track the optimization.\n",
"\n",
"This example reproduces the \"Example II\" from that paper, considering the optimization towards a $\\sqrt{\\text{iSWAP}}$ two-qubit gate on a system of two transmons with a shared transmission line resonator.\n",
"\n",
"[1]: https://michaelgoerz.net/research/Goerz_NJP2014.pdf"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## The two-transmon system"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We consider the Hamiltonian from Eq (17) in the paper, in the rotating wave approximation, together with spontaneous decay and dephasing of each qubit. Alltogether, we define the Liouvillian as follows:"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.310453Z",
"start_time": "2021-06-13T15:05:13.306100Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:57.881921Z",
"iopub.status.busy": "2021-01-13T19:03:57.881050Z",
"iopub.status.idle": "2021-01-13T19:03:57.883779Z",
"shell.execute_reply": "2021-01-13T19:03:57.883270Z"
},
"lines_to_next_cell": 2
},
"outputs": [],
"source": [
"def two_qubit_transmon_liouvillian(\n",
" ω1, ω2, ωd, δ1, δ2, J, q1T1, q2T1, q1T2, q2T2, T, Omega, n_qubit\n",
"):\n",
" from qutip import tensor, identity, destroy\n",
"\n",
" b1 = tensor(identity(n_qubit), destroy(n_qubit))\n",
" b2 = tensor(destroy(n_qubit), identity(n_qubit))\n",
"\n",
" H0 = (\n",
" (ω1 - ωd - δ1 / 2) * b1.dag() * b1\n",
" + (δ1 / 2) * b1.dag() * b1 * b1.dag() * b1\n",
" + (ω2 - ωd - δ2 / 2) * b2.dag() * b2\n",
" + (δ2 / 2) * b2.dag() * b2 * b2.dag() * b2\n",
" + J * (b1.dag() * b2 + b1 * b2.dag())\n",
" )\n",
"\n",
" H1_re = 0.5 * (b1 + b1.dag() + b2 + b2.dag()) # 0.5 is due to RWA\n",
" H1_im = 0.5j * (b1.dag() - b1 + b2.dag() - b2)\n",
"\n",
" H = [H0, [H1_re, Omega], [H1_im, ZeroPulse]]\n",
"\n",
" A1 = np.sqrt(1 / q1T1) * b1 # decay of qubit 1\n",
" A2 = np.sqrt(1 / q2T1) * b2 # decay of qubit 2\n",
" A3 = np.sqrt(1 / q1T2) * b1.dag() * b1 # dephasing of qubit 1\n",
" A4 = np.sqrt(1 / q2T2) * b2.dag() * b2 # dephasing of qubit 2\n",
"\n",
" L = krotov.objectives.liouvillian(H, c_ops=[A1, A2, A3, A4])\n",
" return L"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We will use internal units GHz and ns. Values in GHz contain an implicit factor 2π, and MHz and μs are converted to GHz and ns, respectively:"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.320502Z",
"start_time": "2021-06-13T15:05:13.311846Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:57.888423Z",
"iopub.status.busy": "2021-01-13T19:03:57.887670Z",
"iopub.status.idle": "2021-01-13T19:03:57.889938Z",
"shell.execute_reply": "2021-01-13T19:03:57.890439Z"
}
},
"outputs": [],
"source": [
"GHz = 2 * np.pi\n",
"MHz = 1e-3 * GHz\n",
"ns = 1\n",
"μs = 1000 * ns"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This implicit factor $2 \\pi$ is because frequencies ($\\nu$) convert to energies as $E = h \\nu$, but our propagation routines assume a unit $\\hbar = 1$ for energies. Thus, the factor $h / \\hbar = 2 \\pi$."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We will use the same parameters as those given in Table 2 of the paper:"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.329291Z",
"start_time": "2021-06-13T15:05:13.321825Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:57.896571Z",
"iopub.status.busy": "2021-01-13T19:03:57.895859Z",
"iopub.status.idle": "2021-01-13T19:03:57.898107Z",
"shell.execute_reply": "2021-01-13T19:03:57.898605Z"
}
},
"outputs": [],
"source": [
"ω1 = 4.3796 * GHz # qubit frequency 1\n",
"ω2 = 4.6137 * GHz # qubit frequency 2\n",
"ωd = 4.4985 * GHz # drive frequency\n",
"δ1 = -239.3 * MHz # anharmonicity 1\n",
"δ2 = -242.8 * MHz # anharmonicity 2\n",
"J = -2.3 * MHz # effective qubit-qubit coupling\n",
"q1T1 = 38.0 * μs # decay time for qubit 1\n",
"q2T1 = 32.0 * μs # decay time for qubit 2\n",
"q1T2 = 29.5 * μs # dephasing time for qubit 1\n",
"q2T2 = 16.0 * μs # dephasing time for qubit 2\n",
"T = 400 * ns # gate duration"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.339658Z",
"start_time": "2021-06-13T15:05:13.330219Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:57.903116Z",
"iopub.status.busy": "2021-01-13T19:03:57.902378Z",
"iopub.status.idle": "2021-01-13T19:03:57.904592Z",
"shell.execute_reply": "2021-01-13T19:03:57.905083Z"
}
},
"outputs": [],
"source": [
"tlist = np.linspace(0, T, 2000)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"While in the original paper, each transmon was cut off at 6 levels, here we truncate at 5 levels. This makes the propagation faster, while potentially introducing a slightly larger truncation error."
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.347692Z",
"start_time": "2021-06-13T15:05:13.340588Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:57.910272Z",
"iopub.status.busy": "2021-01-13T19:03:57.909252Z",
"iopub.status.idle": "2021-01-13T19:03:57.912081Z",
"shell.execute_reply": "2021-01-13T19:03:57.911569Z"
}
},
"outputs": [],
"source": [
"n_qubit = 5 # number of transmon levels to consider"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In the Liouvillian, note the control being split up into a separate real and imaginary part. As a guess control we use a real-valued constant pulse with an amplitude of 35 MHz, acting over 400 ns, with a switch-on and switch-off in the first 20 ns (see plot below)"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.355982Z",
"start_time": "2021-06-13T15:05:13.348622Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:57.916979Z",
"iopub.status.busy": "2021-01-13T19:03:57.916288Z",
"iopub.status.idle": "2021-01-13T19:03:57.918422Z",
"shell.execute_reply": "2021-01-13T19:03:57.918912Z"
}
},
"outputs": [],
"source": [
"def Omega(t, args):\n",
" E0 = 35.0 * MHz\n",
" return E0 * krotov.shapes.flattop(t, 0, T, t_rise=(20 * ns), func='sinsq')"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The imaginary part start out as zero:"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.364157Z",
"start_time": "2021-06-13T15:05:13.357783Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:57.923419Z",
"iopub.status.busy": "2021-01-13T19:03:57.922575Z",
"iopub.status.idle": "2021-01-13T19:03:57.924587Z",
"shell.execute_reply": "2021-01-13T19:03:57.925359Z"
}
},
"outputs": [],
"source": [
"def ZeroPulse(t, args):\n",
" return 0.0"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We can now instantiate the Liouvillian:"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.378361Z",
"start_time": "2021-06-13T15:05:13.365457Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:57.937665Z",
"iopub.status.busy": "2021-01-13T19:03:57.936822Z",
"iopub.status.idle": "2021-01-13T19:03:57.941402Z",
"shell.execute_reply": "2021-01-13T19:03:57.940864Z"
}
},
"outputs": [],
"source": [
"L = two_qubit_transmon_liouvillian(\n",
" ω1, ω2, ωd, δ1, δ2, J, q1T1, q2T1, q1T2, q2T2, T, Omega, n_qubit\n",
")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The guess pulse looks as follows:"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.381676Z",
"start_time": "2021-06-13T15:05:13.379165Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:57.947423Z",
"iopub.status.busy": "2021-01-13T19:03:57.946666Z",
"iopub.status.idle": "2021-01-13T19:03:57.948745Z",
"shell.execute_reply": "2021-01-13T19:03:57.949290Z"
}
},
"outputs": [],
"source": [
"def plot_pulse(pulse, tlist, xlimit=None):\n",
" fig, ax = plt.subplots()\n",
" if callable(pulse):\n",
" pulse = np.array([pulse(t, None) for t in tlist])\n",
" ax.plot(tlist, pulse/MHz)\n",
" ax.set_xlabel('time (ns)')\n",
" ax.set_ylabel('pulse amplitude (MHz)')\n",
" if xlimit is not None:\n",
" ax.set_xlim(xlimit)\n",
" plt.show(fig)"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.512446Z",
"start_time": "2021-06-13T15:05:13.382436Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:58.020401Z",
"iopub.status.busy": "2021-01-13T19:03:57.977601Z",
"iopub.status.idle": "2021-01-13T19:03:58.187322Z",
"shell.execute_reply": "2021-01-13T19:03:58.187868Z"
}
},
"outputs": [
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAX4AAAEGCAYAAABiq/5QAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAAi8klEQVR4nO3de5xdZX3v8c93bnuSuSQTEkK4xADSUsrRgPHSI7UUq0VLvbRotT2W9tCmx1vtsbZF67Ho69CXPdVq29PTHjgi2FpbRVGk1IqIWmsLJsgliAhIEEJIJte5ZLLn9jt/rLUnm2Fmsmcy65LZ3/frNa/Za+291/PLmslvnv2sZ/0eRQRmZtY8WooOwMzM8uXEb2bWZJz4zcyajBO/mVmTceI3M2sybUUH0IjVq1fHhg0big7DzOy4snXr1j0RsWb6/uMi8W/YsIEtW7YUHYaZ2XFF0mMz7fdQj5lZk3HiNzNrMk78ZmZNxonfzKzJOPGbmTWZzBK/pE5Jd0q6R9L9kt6f7r9O0qOS7k6/NmYVg5mZPVOW0zmrwEURMSSpHfimpH9On/u9iLghw7bNzGwWmSX+SOo9D6Wb7enXcVMD+skDI9yw9QnGJyaza0TK7tBZHTe7kFFmUWcdd4bHzjLuLA+ekePx92/l8nZ+8fxT6Wgrz8h6pjdwSWoFtgLPBv4qIu6Q9GbgKknvA24DroiI6gzv3QxsBli/fn2WYc7ozX+3lXueOJjZL5qXQTBrHnsGq7z9pWcVHcYU5bEQi6SVwI3A24G9wFNAB3A18EhEfGCu92/atCnyvHP30T3D/PSHvsb/uOQcLr/g9NzaLbssf1ey/DXM8jc803OS2ZGzO9+RYdTHa2fp1z/+bZ4aOMzt77ow97YlbY2ITdP351KyISIOSLoduDgiPpTurkr6OPCuPGKYj28+vAeAl559YsGRlEuWQwPH4ahD6rgN3HJy0dknctUtD7Br4DBrezuLDgfIdlbPmrSnj6RlwMuA70lal+4T8BpgW1YxLNR3nxxgxbJ2nnXC8qJDMbPj3KYNfQDc/fiBYgOpk2WPfx1wfTrO3wJ8OiJulvRVSWtIukp3A/8twxgW5MGnBjj7pJ7j8uKXmZXLWWt7AHh49xA/++MFB5PKclbPvcB5M+y/KKs2F0NE8P1dQ/zi+acUHYqZLQHdlTZO6u3kkd1DR39xTsozv6gk+oeqDFXHOfPE7qJDMbMl4tkndvNIvxN/ae3YPwLAKSuXFRyJmS0Vp6/uYvveQ0WHMcWJf5odB5LEf7ITv5ktkpNXLuPgyBjD1fGiQwGc+J/hyTTxn9LnxG9mi+Pklck0zp0HRwqOJOHEP82O/SP0VNro7WwvOhQzWyJqIwhPHjhccCQJJ/5pdhw47N6+mS2qdSuSHn9tRKFoTvzT7Dw4MvVDMjNbDGt7O5HgyYPu8ZfSnqEqa3oqRYdhZktIe2sLfcs72DP0jHqUhXDirxMR7B0a5YRuJ34zW1yruzvYM+jEXzoHR8YYnwxWO/Gb2SJb3V1xj7+M9gyNAslfZjOzxbS6u8Le4dGiwwCc+J9mb/rX+IQu9/jNbHGt7q54qKeMpnr8Pe7xm9niWt3TwfDoBCOjE0WH4sRfb++we/xmlo3VaV4pwzi/E3+dPUOjSLCqyz1+M1tctZGEfif+ctkzVGXV8g5aW7wAi5ktrlVpj39/CS7wOvHX2TtU5QTP6DGzDKxcltT/OnBorOBInPifZv+hMVYud+I3s8XXl+aWAyNO/KUyMDI29VfZzGwx9XS20SI4cGgJD/VI6pR0p6R7JN0v6f3p/tMl3SHpYUn/KKk0XewDh8ZY4cRvZhloaRErlrUv+aGeKnBRRDwX2AhcLOlFwJ8AH4mIZwP7gcszjGFeDo448ZtZdlYu72D/Uu7xR6K2unB7+hXARcAN6f7rgddkFcN8VMcnGBmbYOVyJ34zy8aKZe0cXOpj/JJaJd0N7AZuBR4BDkREbeHJJ4BTZnnvZklbJG3p7+/PMkyAqR+Ge/xmlpW+5Ut/qIeImIiIjcCpwAuAs+fx3qsjYlNEbFqzZk1WIU4ZSBN/rxO/mWVkyQ/11IuIA8DtwE8AKyW1pU+dCuzII4ajqfX4PZ3TzLKycnk7B5dyj1/SGkkr08fLgJcBD5D8Abg0fdllwBeyimE+PNRjZllbuayDweo4YxOThcaRZY9/HXC7pHuBbwO3RsTNwB8A75T0MHAC8LEMY2hYbdzNid/MstLXleSXoi/wth39JQsTEfcC582w/wck4/2l4h6/mWVtxVTZhtFCV/rznbupWuLv7czsb6GZNbna5JGBw+NHeWW2nPhTB0fG6Km00dbqU2Jm2ah1LAed+Mvh4KExT+U0s0z1dCY5ZvBwsWP8TvypgyNjvmvXzDLV4x5/ubhOj5llrTft8Q8UPKvHiT91wInfzDK2vKOV1ha5x18W7vGbWdYk0V1p8xh/WQwdHp8afzMzy0pPZ5t7/GUwPjHJyNgE3RX3+M0sWz2d7Qy4x1+84dEJALoqrQVHYmZLXW9nm2/gKoPhavJD6K54qMfMstXT2V74UM+cmU7SqcAbgJ8ETgZGgG3APwH/HBHFlphbJLXE3+XEb2YZ6+1s43sFD/XMmukkfZxkdaybSdbJ3Q10Aj8CXAz8oaQrIuIbeQSapSH3+M0sJz2dbYXP458r0304IrbNsH8b8DlJHcD6bMLK11Ti96weM8tY77J2hqrjRASSColh1jH+WtKX9CZJPfXPSbokIkYj4uGsA8zD1FBPhxO/mWWrp7ONyTgyqaQIjVzc/UvgXyX9WN2+D2QUTyGGqskPwEM9Zpa1MhRqayTxPwr8V+AGSa9L9xXz+SQjRy7uejqnmWWrDIXaGuniRkTcJemngE9JeiGwpDLkkGf1mFlOylCorZEe/06AiNgD/CwQwLlZBpW3oeo47a2i0ubbGswsW2Xo8R8100XEz9U9noyI34uIJZUhh6vjdFXaCrvCbmbNozbGX2TZhrnm8X+RpHc/o4h41VwHlnQa8AlgbXqcqyPizyVdCfwm0J++9D0Rccs8415UQ9Vxz+gxs1yUYfnFubLdh9LvAq4BfmOexx4Hfje9PtADbJV0a/rcRyLiQ3O8N1fD1XHP6DGzXNSuJdYmlRRh1mwXEV+vPZY0VL/diIjYyZHrA4OSHiC5E7h0hqrjntFjZrlY3tGKVGzib3SsftYhn0ZI2gCcB9yR7nqbpHslXSupb5b3bJa0RdKW/v7+mV6yaIaqE3R3uiSzmWVPEl0dbVP3DxVh1sQvaVXtC2iV1DdtX0MkdQOfBX4nIgaAvwbOBDaSfCL48Ezvi4irI2JTRGxas2bNPP5J85cM9bjHb2b56Kq0lnOoB9hK0tOvTXW5q+65AM442sEltZMk/U9GxOcAImJX3fPXkBSBK9SwL+6aWY66K20MjZYw8UfE6cdyYCVzIz8GPBARf1a3f106/g/wWpKib4UaSqdzmpnlobvSxlAZZ/VIOn+uN0bEXXM9D7wYeBNwn6S7033vAd4oaSPJp4btwG81GGsmIsKzeswsV12VttIO9Wwh6Y3vSbfr724K4KK5DhwR32Tmmj6FztmfbmRsgslwSWYzy09XpY19w4cKa3+ubPdO4FKSVbf+AbgxIoZyiSpHrtNjZnnrrrQxXOAY/1z1+D8aERcAbwdOA26T9Ol0mGbJGJ4qyexZPWaWj6LH+Bup1fMD4AvAl4EXkCy9uGR4ERYzy1syxl/cPP65Lu6eQbLQ+quBx0mGe/44IkZyii0XtXoZvrhrZnnprrQyOjHJ6PgkHQVUBZ4r2z0M3EvS2x8gWV/3zbUKlvVTNI9nwx7jN7Oc1dfr6WjryL39ubLdBzhSqqE7h1gKUbvA4lk9ZpaXWuIfqo7T11WixB8RV+YYR2Fqs3o81GNmeempS/xFmKtWz3tnK6CWPn+RpEuyCSs/Huoxs7wVXZp5rmx3H3CzpMMkdXr6gU7gLJICa18B/jjrALNWm1K1vN3TOc0sH10F9/jnGur5AvAFSWeRlF9YR3KR9++AzUtlds9QdYKujlZaWrzsopnlo3uqx1/MlM6jjm9ExEPAQznEUojh6rgv7JpZrmoLPxU11LOkFk1fiKFRV+Y0s3z1VJKFnwad+Ivhypxmljf3+As2dNiLsJhZvtpaW6i0tZQ38Uv6EUm3SdqWbj9H0nuzDy0fXoTFzIrQXWkr3zz+OtcA7wbGACLiXpIaPkvC8Og4Pb64a2Y56+4sbjGWRhL/8oi4c9q+4uqJLrLh6sTUeJuZWV66Osrd498j6UzSuj2SLgV2zv2W44eHesysCEUO9TSS8d4KXA2cLWkH8CjwXzKNKiej40lZ1G5f3DWznHVVWtkzNFpI2w0txBIRPwOsAc6OiAsiYvvR3ifpNEm3S/qupPslvSPdv0rSrZIeSr/PWg8oa67TY2ZFKXLB9bkWYnnnLPuBhurxjwO/GxF3SeoBtkq6Ffg14LaI+KCkK4ArgD9YQOzHbKoypy/umlnOejrLOdTTk37/UeD5wE3p9s8D0y/2PkNE7CS9FhARg5IeAE4hWdHrwvRl1wNfo6DEP1WL3z1+M8tZkRd35yrS9n4ASd8Azo+IwXT7SuCf5tOIpA3AecAdwNr0jwLAU8DaWd6zGdgMsH79+vk01zAP9ZhZUboqbRwanWByMnIvEtnIrJ61QP0ViFFmSdYzkdQNfBb4nYgYqH8uIoIjq3wx7bmrI2JTRGxas2ZNo83Ny5H1dj2d08zyNVWhczT/Xn8jXd1PAHdKujHdfg3JEM1RSWonSfqfjIjPpbt3SVoXETslrQN2zzPmRVMrieoev5nlrauuNHNPZ3uubTcyq+cq4NeB/enXr0fEURdgUXIV+GPAA9MuBN8EXJY+voxkMfdCDHvZRTMrSO3G0SLG+Y+a8SStB/YAN9bvi4gfHuWtLwbeBNwn6e5033uADwKflnQ58Bjw+gXEvSi83q6ZFaW7wOUXG8l4/8SRcfhlwOnAg8CPz/WmiPgmMNsVi5c2GmCWfHHXzIpS5Lq7jazA9Z/qtyWdD7wls4hyNFQdp6OthfbWpq9ObWY56y5w3d15Z7yIuAt4YQax5G7Ii7CYWUG6yjyrZ9odvC3A+cCTmUWUI6++ZWZFOXJxN/8F1xvJej11j8dJxvw/m004+RqqTnh838wKUfaLu9+NiM/U75D0OuAzs7z+uJH0+H3zlpnlb1l7Ky1Kln/NWyNj/O9ucN9xx7X4zawokgqr1zNXdc5XAK8ETpH0F3VP9bJEVuAaro6z/oTlRYdhZk2qqNLMc3V3nwS2AK8CttbtHwT+e5ZB5WWoOu5FWMysMN2dbeWa1RMR9wD3SPpkRCyJHv50w9Vx1+I3s8J0VdrKNatH0qcj4vXAdyQ9o4JmRDwn08gyNjkZDI96Vo+ZFae70lq6oZ53pN8vySOQvB1ZhMWzesysGF0dbewdOpR7u3MN9dRWz3osv3Dy45LMZla07kr5ZvUM8vRFUpRui2QNld6MY8uUK3OaWdFKN6snInpme24pcC1+MytakvhLdHG3XlqR8wKSHv83I+I7mUaVA5dkNrOidVdaGZ2YZHR8ko62/KoEH7UlSe8jWWrxBGA1cJ2k92YdWNYG3eM3s4IVVZO/kaz3K8BzI+IwgKQPAncD/zPDuDLnHr+ZFa2rriZ/X1dHbu028tniSaCzbrsC7MgmnPwcSfyezmlmxeguqCZ/I93dg8D9km4lGeN/GXBnrX5PRPx2hvFlpna3XE8l39XtzcxqyjzUcyN1C60DX2vkwJKuJbn5a3dEnJvuuxL4TaA/fdl7IuKWRoNdTMPVcVoEne1edtHMitFd0GIsjay5e/0Cj30d8L+BT0zb/5GI+NACj7loaiWZpdnWgzczy1ZRPf5GZvVcIuk7kvZJGpA0KGngaO+LiG8A+xYlygx4vV0zK1pXRzELrjcyzvFR4DLghIjojYieY7xr922S7pV0raS+2V4kabOkLZK29Pf3z/ayBRv2IixmVrBa5zPvVbgaSfyPA9si4hkVOhfgr4EzgY3ATuDDs70wIq6OiE0RsWnNmjWL0PTTucdvZkUr88Xd3wdukfR1oFrbGRF/Nt/GImJX7bGka4Cb53uMxTLsxG9mBetoa6GjtYWhnKdzNtLjvwo4RDKXv6fua94kravbfC2wbSHHWQzJxV3P4TezYnUVUJO/kS7vybXpmPMh6VPAhcBqSU8AfwRcKGkjyf0A24Hfmu9xF8tw1YuwmFnxiijU1kjmu0XSyyPiy/M5cES8cYbdH5vPMbLkMX4zK4MiavI3MtTzZuBLkkbmM52zzCLCY/xmVgpF1ORv5AauJVeXvzo+yfhkeKjHzArXVWnj4KHRXNtstB5/H3AWdcXa0hu0jktefcvMyqK70sqO/SXr8Uv6DZKF108lKcf8IuDfgYsyjSxDLslsZmXR1ZH/xd1GxvjfATwfeCwifho4DziQZVBZO9Lj93ROMytWEWP8jST+w3WLsFQi4nvAj2YbVrZqt0d3uySzmRWsp7ON4dFxFqc4QmMaGet4QtJK4PPArZL2A49lGVTWaosedHd6qMfMitVVaWMyYGRsguUd+eSkRmb1vDZ9eKWk24EVwJcyjSpjg4d9cdfMyqF++cXSJP56EfH1rALJk2f1mFlZ1K41DlcnFlgMZ/6acvmp2oUUD/WYWdFqNfnzvMDblIl/6PA4Eixv96weMytWd91QT14aSvySniXpZ9LHyyQd13fzDlbH6e5oo6XFyy6aWbGKqMnfyNKLvwncAPzfdNepJDN8jltefcvMyqKrpD3+twIvBgYAIuIh4MQsg8raUHXc4/tmVgplHeqpRsRUBSFJbST19I9bg4ddmdPMyqFralZPuRL/1yW9B1gm6WXAZ4AvZhtWtoaq4/S4x29mJVCb1TOUY72eRhL/FUA/cB/Jilm3AO/NMqisDVfHp062mVmRWlrE8o58l19s5M7dSeAa4BpJq4BTI8+iEhkYOuwxfjMrj7wLtTUyq+drknrTpL+V5A/AR7IPLTuDXn3LzEok7+UXGxnqWRERA8AvAJ+IiBcCLz3amyRdK2m3pG11+1ZJulXSQ+n3voWHvjBedtHMyqarku9QTyOJv03SOuD1wM3zOPZ1wMXT9l0B3BYRZwG3pdu5GhmbYDJcrsHMyiPvxVgaSfwfAP4FeDgivi3pDOCho70pXZpx37TdrwauTx9fD7ym8VAXx5Arc5pZyeQ91NPIxd3PkEzhrG3/APjFBba3NiJ2po+fAtbO9kJJm4HNAOvXr19gc89UO7mezmlmZdFVaZtaJyQPs2Y/SX/JHDdqRcRvH0vDERGS5jr+1cDVAJs2bVq0WUS1xO/pnGZWFnnP6pkr+23JoL1dktZFxM70usHuDNqY09RQj3v8ZlYS3ZXWcgz1RMT1sz13DG4CLgM+mH7/QgZtzGnQi7CYWcl0Vdo4PDbJ+MQkba3ZV8s/avZLl1t8xlBLRFx0lPd9CrgQWC3pCeCPSBL+pyVdTrJu7+sXEPMxGfYYv5mVTK0jOjw6wYplJUj8wLvqHneSXNg96meSiHjjLE8d9R6ALE2N8bvHb2YlUV+Tf8Wy9szba2RWz9Zpu/5N0p0ZxZM5L7RuZmWT92IsjQz1rKrbbAGeB6zILKKMDVfHaW8VlbamXHXSzEqotuB6Xhd4G+n2biUZ4xfJEM+jwOVZBpWlobRcg+RlF82sHLoryfBObUQia40M9ZyeRyB5GTrsZRfNrFx6lyU5qTSJX1In8BbgApKe/78CfxMRhzOOLRNDLtBmZiXT21nr8Y/l0l4jGfATwCDwl+n2LwN/C7wuq6CyNHjYq2+ZWbn0pjN5BkqU+M+NiHPqtm+X9N2sAsrawOEx1q3oLDoMM7MpXR2ttAgGRvIZ6mlkastdkl5U25D0QrIp55CLgyNjUx+rzMzKQBK9y9pL1eN/HvAtST9Mt9cDD0q6j6TW2nMyiy4DAyNjUx+rzMzKoreznYGR8iT+6YupHLcmJ4PB6ji9HuM3s5LpXdbGQFlm9UTEY3kEkoeh0XEicI/fzEonzx5/U92+WjupHuM3s7Lp6WzLbYy/yRJ/8jHKPX4zK5vezvbcbuBqrsSf/jWt3SVnZlYWvcs81JMJD/WYWVn1drYzPDrB+MRk5m01VeI/mCb+POpdm5nNR571epoq8demSrnHb2ZlU8tLeVzgba7En/b4vdC6mZXNVL2eHMo2NFfiPzxGT6WN1hbX4jezcqndWOoe/yIbGBn3VE4zK6UjPf7sE38hYx6StpOUep4AxiNiUx7tDhwec0lmMyulPEszF5kFfzoi9uTZoAu0mVlZTQ31eIx/cR04NMZKJ34zK6Hu9PrjgZHRzNsqKvEH8GVJWyVtnukFkjZL2iJpS39//6I0uv/QKKu6OhblWGZmi0kSfcs72De8dC/uXhAR5wOvAN4q6SXTXxARV0fEpojYtGbNmmNuMCLYf2iUPid+MyupVV3t7B9eoj3+iNiRft8N3Ai8IOs2h6rjjE0Efcs91GNm5dS3vIN9h5Zg4pfUJamn9hh4ObAt63b3px+f+pa7x29m5bSqqyOXHn8Rs3rWAjdKqrX/9xHxpawb3Z/+FfUYv5mVVV9XB/uWYuKPiB8Az8273drHp5Xu8ZtZSa1a3sH+Q6NMTgYtGVYYaJrpnLWPT+7xm1lZ9XV1MBnZ38TVNIm/9vFplXv8ZlZSq7qSySdZD/c0TeI/cGiM1ha5ZIOZlVZt8sn+jGf2NE3i33dolJXL2jMdNzMzOxYndFUAMr+Jq2kS//7hUVZ6Dr+ZlVhfOtST9ZTOpkn8/YNV1vRUig7DzGxWtckne534F8fuwSprezuLDsPMbFbL2luptLWwb7iaaTtNkfgjgl0DhznRPX4zKzFJrO3tZNeAE/8xGzg8TnV8khN73OM3s3I7qbeTpwYOZ9pGUyT+/sHkJJ7Y6x6/mZXb2hWd7HLiP3a7049N7vGbWdmt7anw1MHDRERmbTRF4t+V9vjXusdvZiV30opOquOTmS7B2ByJv9bj96weMyu52uzDLMf5myLx/3DfIfqWt9NdcbkGMyu3k1YkiX/nwZHM2miKxP/4vkOsX7W86DDMzI7qtL4kVz2+71BmbTRF4v/hvkOc5sRvZseBtb0VlrW38ugeJ/4FG5+YZMf+Eff4zey4IIlnnbCc7XuHM2tjySf+7XsPMT4ZnL66q+hQzMwacvrqLrbvceJfsAd2DgBwzsm9BUdiZtaYHz2ph0f3DjOY0UpchSR+SRdLelDSw5KuyLKt+3YcpL1VPPvE7iybMTNbNM97Vh8RcPfjBzI5fu6JX1Ir8FfAK4BzgDdKOier9r7+YD/P37CKSltrVk2YmS2qjaetRIJvPbI3k+MX0eN/AfBwRPwgIkaBfwBenUVDH7n1+zy4a5BXnHtSFoc3M8tET2c7Lz37RK77t+38ewbJv4jEfwrweN32E+m+p5G0WdIWSVv6+/sX1NAZa7r45Reu55eev35hkZqZFeSPfv7H2bShL5OJKcqyENCMDUqXAhdHxG+k228CXhgRb5vtPZs2bYotW7bkFaKZ2ZIgaWtEbJq+v4ge/w7gtLrtU9N9ZmaWgyIS/7eBsySdLqkDeANwUwFxmJk1pdyrlkXEuKS3Af8CtALXRsT9ecdhZtasCilXGRG3ALcU0baZWbNb8nfumpnZ0znxm5k1GSd+M7Mm48RvZtZkcr+BayEk9QOPLfDtq4E9ixjOYnFc8+O45sdxzU9Z44Jji+1ZEbFm+s7jIvEfC0lbZrpzrWiOa34c1/w4rvkpa1yQTWwe6jEzazJO/GZmTaYZEv/VRQcwC8c1P45rfhzX/JQ1LsggtiU/xm9mZk/XDD1+MzOr48RvZtZklnTiz3NR9wZi2S7pPkl3S9qS7lsl6VZJD6Xf+3KI41pJuyVtq9s3YxxK/EV6/u6VdH7OcV0paUd6zu6W9Mq6596dxvWgpJ/NMK7TJN0u6buS7pf0jnR/oedsjrgKPWeSOiXdKemeNK73p/tPl3RH2v4/piXZkVRJtx9On9+Qc1zXSXq07nxtTPfn9rufttcq6TuSbk63sz1fEbEkv0hKPj8CnAF0APcA5xQYz3Zg9bR9/wu4In18BfAnOcTxEuB8YNvR4gBeCfwzIOBFwB05x3Ul8K4ZXntO+vOsAKenP+fWjOJaB5yfPu4Bvp+2X+g5myOuQs9Z+u/uTh+3A3ek5+HTwBvS/X8DvDl9/Bbgb9LHbwD+MaPzNVtc1wGXzvD63H730/beCfw9cHO6nen5Wso9/twWdT8GrwauTx9fD7wm6wYj4hvAvgbjeDXwiUj8B7BS0roc45rNq4F/iIhqRDwKPEzy884irp0RcVf6eBB4gGSN6ELP2RxxzSaXc5b+u4fSzfb0K4CLgBvS/dPPV+083gC8VJJyjGs2uf3uSzoV+Dng/6XbIuPztZQTf0OLuucogC9L2ippc7pvbUTsTB8/BawtJrRZ4yjDOXxb+lH72rqhsELiSj9Wn0fSWyzNOZsWFxR8ztJhi7uB3cCtJJ8uDkTE+AxtT8WVPn8QOCGPuCKidr6uSs/XRyRVpsc1Q8yL7aPA7wOT6fYJZHy+lnLiL5sLIuJ84BXAWyW9pP7JSD67FT63tixxpP4aOBPYCOwEPlxUIJK6gc8CvxMRA/XPFXnOZoir8HMWERMRsZFkPe0XAGfnHcNMpscl6Vzg3STxPR9YBfxBnjFJugTYHRFb82x3KSf+Ui3qHhE70u+7gRtJ/kPsqn18TL/vLii82eIo9BxGxK70P+skcA1HhiZyjUtSO0ly/WREfC7dXfg5mymuspyzNJYDwO3AT5AMldRW/Ktveyqu9PkVwN6c4ro4HTKLiKgCHyf/8/Vi4FWStpMMR18E/DkZn6+lnPhLs6i7pC5JPbXHwMuBbWk8l6Uvuwz4QhHxzRHHTcCvpjMcXgQcrBveyNy0MdXXkpyzWlxvSGc4nA6cBdyZUQwCPgY8EBF/VvdUoedstriKPmeS1khamT5eBryM5PrD7cCl6cumn6/aebwU+Gr6CSqPuL5X98dbJOPo9ecr859jRLw7Ik6NiA0kOeqrEfErZH2+FvPKdNm+SK7Mf59kjPEPC4zjDJIZFfcA99diIRmbuw14CPgKsCqHWD5FMgQwRjJ2ePlscZDMaPir9PzdB2zKOa6/Tdu9N/2FX1f3+j9M43oQeEWGcV1AMoxzL3B3+vXKos/ZHHEVes6A5wDfSdvfBryv7v/AnSQXlT8DVNL9nen2w+nzZ+Qc11fT87UN+DuOzPzJ7Xe/LsYLOTKrJ9Pz5ZINZmZNZikP9ZiZ2Qyc+M3MmowTv5lZk3HiNzNrMk78ZmZNxonfmoKklZLeUrd9sqQb5nrPMbT1GknvW+B7v6IcqrRac/N0TmsKaT2bmyPi3Bza+hbwqojYs4D3XgacGhFXLX5kZgn3+K1ZfBA4M625/qeSNiit/S/p1yR9Xkld/e2S3ibpnWl99P+QtCp93ZmSvpQW2vtXSc+oQSPpR4BqLekrqff+F5K+JekHki5N96+T9I00nm2SfjI9xE3AG/M4Ida8nPitWVwBPBIRGyPi92Z4/lzgF0iKdV0FHIqI84B/B341fc3VwNsj4nnAu4D/M8NxXgzcNW3fOpI7bS8h+QME8MvAv0RSNOy5JHfeEhH7gYqkTCpUmgG0Hf0lZk3h9kjq2g9KOgh8Md1/H/CctArmfwY+U1f+vPLMw7AO6J+27/ORFE37rqRa+eZvA9emhdY+HxF3171+N3AyGRcrs+blHr9Zolr3eLJue5Kkg9RCUiN9Y93Xj81wnBGSeiqzHVswtfDMS0iqLV4n6VfrXtOZHscsE0781iwGSZYoXJBIat0/Kul1MLUm63NneOkDwLOPdjxJzwJ2RcQ1JCsvnV87LnASyVKdZplw4remEBF7gX9LL6T+6QIP8yvA5ZJqVVZnWsrzG8B50lGXw7sQuEfSd4BfIqnBDvA84D/iyOpLZovO0znNFpmkPwe+GBFfWeB7b4qI2xY/MrOEe/xmi++PgeULfO82J33Lmnv8ZmZNxj1+M7Mm48RvZtZknPjNzJqME7+ZWZNx4jczazL/HzzoRp2LqJjZAAAAAElFTkSuQmCC\n",
"text/plain": [
"<Figure size 432x288 with 1 Axes>"
]
},
"metadata": {
"needs_background": "light"
},
"output_type": "display_data"
}
],
"source": [
"plot_pulse(L[1][1], tlist)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Optimization objectives"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Our target gate is $\\Op{O} = \\sqrt{\\text{iSWAP}}$:"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.515905Z",
"start_time": "2021-06-13T15:05:13.513303Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:58.194672Z",
"iopub.status.busy": "2021-01-13T19:03:58.193879Z",
"iopub.status.idle": "2021-01-13T19:03:58.196231Z",
"shell.execute_reply": "2021-01-13T19:03:58.196734Z"
}
},
"outputs": [],
"source": [
"SQRTISWAP = qutip.Qobj(np.array(\n",
" [[1, 0, 0, 0],\n",
" [0, 1 / np.sqrt(2), 1j / np.sqrt(2), 0],\n",
" [0, 1j / np.sqrt(2), 1 / np.sqrt(2), 0],\n",
" [0, 0, 0, 1]]),\n",
" dims=[[2, 2], [2, 2]]\n",
")"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The key idea explored in the paper is that a set of three density matrices is sufficient to track the optimization\n",
"\n",
"$$\n",
"\\begin{align}\n",
"\\Op{\\rho}_1\n",
" &= \\sum_{i=1}^{d} \\frac{2 (d-i+1)}{d (d+1)} \\ketbra{i}{i} \\\\\n",
"\\Op{\\rho}_2\n",
" &= \\sum_{i,j=1}^{d} \\frac{1}{d} \\ketbra{i}{j} \\\\\n",
"\\Op{\\rho}_3\n",
" &= \\sum_{i=1}^{d} \\frac{1}{d} \\ketbra{i}{i}\n",
"\\end{align}\n",
"$$"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In our case, $d=4$ for a two qubit-gate, and the $\\ket{i}$, $\\ket{j}$ are the canonical basis states $\\ket{00}$, $\\ket{01}$, $\\ket{10}$, $\\ket{11}$"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.528122Z",
"start_time": "2021-06-13T15:05:13.516717Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:58.204101Z",
"iopub.status.busy": "2021-01-13T19:03:58.203376Z",
"iopub.status.idle": "2021-01-13T19:03:58.205630Z",
"shell.execute_reply": "2021-01-13T19:03:58.206096Z"
}
},
"outputs": [],
"source": [
"ket00 = qutip.ket((0, 0), dim=(n_qubit, n_qubit))\n",
"ket01 = qutip.ket((0, 1), dim=(n_qubit, n_qubit))\n",
"ket10 = qutip.ket((1, 0), dim=(n_qubit, n_qubit))\n",
"ket11 = qutip.ket((1, 1), dim=(n_qubit, n_qubit))\n",
"basis = [ket00, ket01, ket10, ket11]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The three density matrices play different roles in the optimization, and, as shown in the paper, convergence may improve significantly by weighing the states relatively to each other. For this example, we place a strong emphasis on the optimization $\\Op{\\rho}_1 \\rightarrow \\Op{O}^\\dagger \\Op{\\rho}_1 \\Op{O}$, by a factor of 20. This reflects that the hardest part of the optimization is identifying the basis in which the gate is diagonal. We will be using the real-part functional ($J_{T,\\text{re}}$) to evaluate the success of $\\Op{\\rho}_i \\rightarrow \\Op{O}\\Op{\\rho}_i\\Op{O}^\\dagger$. Because $\\Op{\\rho}_1$ and $\\Op{\\rho}_3$ are mixed states, the Hilbert-Schmidt overlap will take values smaller than one in the optimal case. To compensate, we divide the weights by the purity of the respective states."
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.538412Z",
"start_time": "2021-06-13T15:05:13.529657Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:58.211673Z",
"iopub.status.busy": "2021-01-13T19:03:58.210867Z",
"iopub.status.idle": "2021-01-13T19:03:58.213209Z",
"shell.execute_reply": "2021-01-13T19:03:58.213731Z"
}
},
"outputs": [],
"source": [
"weights = np.array([20, 1, 1], dtype=np.float64)\n",
"weights *= len(weights) / np.sum(weights) # manual normalization\n",
"weights /= np.array([0.3, 1.0, 0.25]) # purities"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The `krotov.gate_objectives` routine can initialize the density matrices $\\Op{\\rho}_1$, $\\Op{\\rho}_2$, $\\Op{\\rho}_3$ automatically, via the parameter `liouville_states_set`. Alternatively, we could also use the `'full'` basis of 16 matrices or the extended set of $d+1 = 5$ pure-state density matrices."
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.563694Z",
"start_time": "2021-06-13T15:05:13.539176Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:03:58.243251Z",
"iopub.status.busy": "2021-01-13T19:03:58.242407Z",
"iopub.status.idle": "2021-01-13T19:03:58.245490Z",
"shell.execute_reply": "2021-01-13T19:03:58.245982Z"
}
},
"outputs": [
{
"data": {
"text/plain": [
"[Objective[ρ₀[5⊗5,5⊗5] to ρ₁[5⊗5,5⊗5] via [𝓛₀[[5⊗5,5⊗5],[5⊗5,5⊗5]], [𝓛₁[[5⊗5,5⊗5],[5⊗5,5⊗5]], u₁(t)], [𝓛₂[[5⊗5,5⊗5],[5⊗5,5⊗5]], u₂(t)]]],\n",
" Objective[ρ₂[5⊗5,5⊗5] to ρ₃[5⊗5,5⊗5] via [𝓛₀[[5⊗5,5⊗5],[5⊗5,5⊗5]], [𝓛₁[[5⊗5,5⊗5],[5⊗5,5⊗5]], u₁(t)], [𝓛₂[[5⊗5,5⊗5],[5⊗5,5⊗5]], u₂(t)]]],\n",
" Objective[ρ₄[5⊗5,5⊗5] to ρ₅[5⊗5,5⊗5] via [𝓛₀[[5⊗5,5⊗5],[5⊗5,5⊗5]], [𝓛₁[[5⊗5,5⊗5],[5⊗5,5⊗5]], u₁(t)], [𝓛₂[[5⊗5,5⊗5],[5⊗5,5⊗5]], u₂(t)]]]]"
]
},
"execution_count": 15,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"objectives = krotov.gate_objectives(\n",
" basis,\n",
" SQRTISWAP,\n",
" L,\n",
" liouville_states_set='3states',\n",
" weights=weights,\n",
" normalize_weights=False,\n",
")\n",
"objectives"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The use of `normalize_weights=False` is because we have included the purities in the weights, as discussed above."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Propagator Benchmarks"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In the following, we benchmark different propagators in a context comparable to that internal propagations inside an optimizatin with Krotov's method. most importantly, we'll use discretized controls (numpy arrays, not functions):"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.576108Z",
"start_time": "2021-06-13T15:05:13.564617Z"
}
},
"outputs": [],
"source": [
"Liouvillian = [\n",
" L[0],\n",
" [L[1][0], krotov.conversions.discretize(L[1][1], tlist)],\n",
" [L[2][0], krotov.conversions.discretize(L[2][1], tlist)],\n",
"]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### QuTiP's mesolve"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We just call `qutip.mesolve` directly, propagating over the entire time grid in one go. This is our \"baseline\": the fastest way to propagate, but not an option for `krotov`, since we need to propagate individual time steps"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:05:13.588772Z",
"start_time": "2021-06-13T15:05:13.576869Z"
}
},
"outputs": [],
"source": [
"def prop_all_with_mesolve():\n",
" for obj in objectives:\n",
" obj.mesolve(H=Liouvillian, tlist=tlist)"
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:06:02.977669Z",
"start_time": "2021-06-13T15:05:13.589718Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"6.17 s ± 8.04 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n"
]
}
],
"source": [
"%%timeit\n",
"prop_all_with_mesolve()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### QuTiP's mesolve in a loop"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We can call `qutip.mesolve` for each time step."
]
},
{
"cell_type": "code",
"execution_count": 19,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:06:02.981659Z",
"start_time": "2021-06-13T15:06:02.978691Z"
}
},
"outputs": [],
"source": [
"def mesolve_loop(Liouvillian, rho0, tlist):\n",
" \"\"\"Call mesolve for each individual time step.\"\"\"\n",
" options = qutip.Options(store_states=False, store_final_state=True)\n",
" L0, (L1, u1), (L2, u2) = Liouvillian\n",
" rho = rho0\n",
" for i in range(len(tlist) - 1):\n",
" dt = tlist[i + 1] - tlist[i]\n",
" rho = qutip.mesolve(\n",
" [L0, [L1, u1[i : i + 2]], [L2, u2[i : i + 2]]],\n",
" rho,\n",
" tlist[i : i + 2],\n",
" options=options,\n",
" ).final_state\n",
" return rho"
]
},
{
"cell_type": "code",
"execution_count": 20,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:06:02.994991Z",
"start_time": "2021-06-13T15:06:02.982754Z"
}
},
"outputs": [],
"source": [
"def prop_all_with_mesolve_loop():\n",
" for obj in objectives:\n",
" mesolve_loop(Liouvillian, obj.initial_state, tlist)"
]
},
{
"cell_type": "code",
"execution_count": 21,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:08:03.838290Z",
"start_time": "2021-06-13T15:06:02.996035Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"15.1 s ± 28.6 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n"
]
}
],
"source": [
"%%timeit\n",
"prop_all_with_mesolve_loop()"
]
},
{
"cell_type": "markdown",
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T02:36:37.764288Z",
"start_time": "2021-06-13T02:36:37.758700Z"
}
},
"source": [
"We observe a rather large overhead from `mesolve` having to re-initialize its internal ODE solver in each time step"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### QuTiP's mesolve in a loop with a pre-summed Liouvillian in each time step"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"When we propagate by time steps, we can exploit the fact that the Liouvillian in each step is constant: Instead of passing a nested-list Liouvillian where the time-dependency just happens to be constant, we calculate the constant Liouvillian for that time step and pass that directly."
]
},
{
"cell_type": "code",
"execution_count": 22,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:08:03.842213Z",
"start_time": "2021-06-13T15:08:03.839280Z"
}
},
"outputs": [],
"source": [
"def mesolve_loop_pwc(Liouvillian, rho0, tlist):\n",
" \"\"\"Like mesolve_loop, but sum the Liouvillian into a single constant term for each time step.\"\"\"\n",
" L0, (L1, u1), (L2, u2) = Liouvillian\n",
" u1 = krotov.conversions.control_onto_interval(u1)\n",
" u2 = krotov.conversions.control_onto_interval(u2)\n",
" options = qutip.Options(store_states=False, store_final_state=True)\n",
" rho = rho0\n",
" for i in range(len(tlist) - 1):\n",
" dt = tlist[i + 1] - tlist[i]\n",
" rho = qutip.mesolve(\n",
" qutip.QobjEvo(L0 + u1[i] * L1 + u2[i] * L2), rho, [0, dt], options=options\n",
" ).final_state\n",
" # whether or not we wrap the Liouvillian in qutip.QobjEvo does not seem\n",
" # to make a difference\n",
" return rho"
]
},
{
"cell_type": "code",
"execution_count": 23,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:08:03.852382Z",
"start_time": "2021-06-13T15:08:03.844879Z"
}
},
"outputs": [],
"source": [
"def prop_all_with_mesolve_loop_pwc():\n",
" for obj in objectives:\n",
" mesolve_loop_pwc(Liouvillian, obj.initial_state, tlist)"
]
},
{
"cell_type": "code",
"execution_count": 24,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:09:25.863667Z",
"start_time": "2021-06-13T15:08:03.853844Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"10.3 s ± 13 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n"
]
}
],
"source": [
"%%timeit\n",
"prop_all_with_mesolve_loop_pwc()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We see a rather substantial speedup."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We can do the exact same propagation by wrapping it in Krotov's `Propagator` class (which is suitable to being passed to the `optimize_pulses` routine)"
]
},
{
"cell_type": "code",
"execution_count": 25,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:09:25.867276Z",
"start_time": "2021-06-13T15:09:25.864622Z"
}
},
"outputs": [],
"source": [
"class MesolvePWCStepPropagator(krotov.propagators.Propagator):\n",
" \"\"\"Like mesolve_loop_pwc, but as a krotov Propagator.\"\"\"\n",
"\n",
" options = qutip.Options(store_states=False, store_final_state=True)\n",
"\n",
" def __call__(\n",
" self, H, state, dt, c_ops=None, backwards=False, initialize=False\n",
" ):\n",
" L = H[0]\n",
" for L_ctrl in H[1:]:\n",
" L += L_ctrl[1] * L_ctrl[0]\n",
" return qutip.mesolve(\n",
" L, state, [0, dt], options=self.options\n",
" ).final_state"
]
},
{
"cell_type": "code",
"execution_count": 26,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:09:25.879050Z",
"start_time": "2021-06-13T15:09:25.868438Z"
}
},
"outputs": [],
"source": [
"def prop_all_with_mesolve_pwc_step_propagator():\n",
" for obj in objectives:\n",
" obj.propagate(tlist, H=Liouvillian, propagator=MesolvePWCStepPropagator())"
]
},
{
"cell_type": "code",
"execution_count": 27,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:10:48.706953Z",
"start_time": "2021-06-13T15:09:25.879934Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"10.4 s ± 10.4 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n"
]
}
],
"source": [
"%%timeit\n",
"prop_all_with_mesolve_pwc_step_propagator()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We observe no additional overhead."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### DensityMatrixODEPropagator"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The `krotov` package contains a `DensityMatrixODEPropagator`, which is was written by extracting the code from `mesolve` that handles the specific case of propagating a Liouville-von-Neuman equation."
]
},
{
"cell_type": "code",
"execution_count": 28,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:10:48.710164Z",
"start_time": "2021-06-13T15:10:48.708079Z"
}
},
"outputs": [],
"source": [
"def prop_all_with_density_matrix_ode_propagator():\n",
" for obj in objectives:\n",
" obj.propagate(\n",
" tlist,\n",
" H=Liouvillian,\n",
" propagator=krotov.propagators.DensityMatrixODEPropagator(\n",
" reentrant=True\n",
" ),\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": 29,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:11:54.853691Z",
"start_time": "2021-06-13T15:10:48.711303Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"8.26 s ± 15 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n"
]
}
],
"source": [
"%%timeit\n",
"prop_all_with_density_matrix_ode_propagator()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This is actually remarkably fast! The `reentrant=True` parameter means that the internal ODE solver gets re-initialized in each time step. Thus, the propagator is comparable to calling `mesolve` in a loop for each time step (the slowest propagator). It would seem that the overhead of `mesolve` figuring out that it's handling the special case of a Liouville-von-Neumann equation adds up to 50% of the runtime!"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### DensityMatrixODEPropagator with a pre-summed Liouvillian in each time step"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We might try to see if we can again get a speedup from summing the parts of the Liouvillian into a single constant object for each time step."
]
},
{
"cell_type": "code",
"execution_count": 30,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:11:54.858499Z",
"start_time": "2021-06-13T15:11:54.854699Z"
}
},
"outputs": [],
"source": [
"from qutip.superoperator import mat2vec, vec2mat\n",
"from qutip.cy.spconvert import dense2D_to_fastcsr_fmode\n",
"\n",
"class DensityMatrixODEPropagatorPWC(krotov.propagators.DensityMatrixODEPropagator):\n",
" \"\"\"Propagator where we sum drift and control terms.\"\"\"\n",
" def __call__(\n",
" self, H, state, dt, c_ops=None, backwards=False, initialize=False\n",
" ):\n",
" L = H[0]\n",
" for (L_c, coeff) in H[1:]:\n",
" L += coeff * L_c\n",
" if initialize or self.reentrant:\n",
" self._initialize([L], state, dt, c_ops, backwards)\n",
" else:\n",
" if self.reentrant:\n",
" self._initialize_integrator(self._y)\n",
" # only update the control values\n",
" for i in self._control_indices:\n",
" self._L_list[i][1] = H[i][1]\n",
" self._t += dt\n",
" self._r.integrate(self._t)\n",
" self._y = self._r.y\n",
" return qutip.Qobj(\n",
" dense2D_to_fastcsr_fmode(\n",
" vec2mat(self._y), state.shape[0], state.shape[1]\n",
" ),\n",
" dims=state.dims,\n",
" isherm=True,\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": 31,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:11:54.872127Z",
"start_time": "2021-06-13T15:11:54.859655Z"
}
},
"outputs": [],
"source": [
"def prop_all_with_density_matrix_ode_propagator_pwc():\n",
" for obj in objectives:\n",
" obj.propagate(\n",
" tlist,\n",
" H=Liouvillian,\n",
" propagator=DensityMatrixODEPropagatorPWC(\n",
" reentrant=True\n",
" ),\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": 32,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:13:05.195127Z",
"start_time": "2021-06-13T15:11:54.873214Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"8.79 s ± 24 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n"
]
}
],
"source": [
"%%timeit\n",
"prop_all_with_density_matrix_ode_propagator_pwc()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"It is somewhat surprizing that this is actually slower! Than again, it's not clear a-priori whether the overhead of summing is more or less than the overhead of applying the nested-list Liouvillian."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Optimization"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We now define the optimization parameters for the controls, the Krotov step size $\\lambda_a$ and the update-shape that will ensure that the pulse switch-on and switch-off stays intact."
]
},
{
"cell_type": "code",
"execution_count": 33,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:13:05.198442Z",
"start_time": "2021-06-13T15:13:05.196260Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:04:44.128878Z",
"iopub.status.busy": "2021-01-13T19:04:44.127809Z",
"iopub.status.idle": "2021-01-13T19:04:44.130192Z",
"shell.execute_reply": "2021-01-13T19:04:44.131054Z"
}
},
"outputs": [],
"source": [
"pulse_options = {\n",
" L[i][1]: dict(\n",
" lambda_a=1.0,\n",
" update_shape=partial(\n",
" krotov.shapes.flattop, t_start=0, t_stop=T, t_rise=(20 * ns))\n",
" )\n",
" for i in [1, 2]\n",
"}"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We now compare doing a few optimal control iterations, with different propagators. The \"0\" iteration consists of a single forward propagation. All iterations >= 1 consist of a backward and a forward propagation. We expect the iteration to take twice as long as iteration \"0\" plus the additional overhead of calculating pulse updates"
]
},
{
"cell_type": "code",
"execution_count": 34,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:14:36.520148Z",
"start_time": "2021-06-13T15:13:05.199304Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:04:44.180324Z",
"iopub.status.busy": "2021-01-13T19:04:44.138592Z",
"iopub.status.idle": "2021-01-13T19:07:01.833927Z",
"shell.execute_reply": "2021-01-13T19:07:01.833039Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"iter. J_T ∑∫gₐ(t)dt J ΔJ_T ΔJ secs\n",
"0 1.22e-01 0.00e+00 1.22e-01 n/a n/a 10\n",
"1 7.49e-02 2.26e-02 9.75e-02 -4.67e-02 -2.41e-02 26\n",
"2 7.41e-02 3.98e-04 7.45e-02 -8.12e-04 -4.14e-04 26\n",
"3 7.33e-02 3.70e-04 7.37e-02 -7.55e-04 -3.85e-04 27\n"
]
}
],
"source": [
"opt_result = krotov.optimize_pulses(\n",
" objectives,\n",
" pulse_options,\n",
" tlist,\n",
" propagator=MesolvePWCStepPropagator(),\n",
" chi_constructor=krotov.functionals.chis_re,\n",
" info_hook=krotov.info_hooks.print_table(J_T=krotov.functionals.J_T_re),\n",
" iter_stop=3,\n",
")"
]
},
{
"cell_type": "code",
"execution_count": 35,
"metadata": {
"ExecuteTime": {
"end_time": "2021-06-13T15:15:54.899821Z",
"start_time": "2021-06-13T15:14:36.521191Z"
},
"execution": {
"iopub.execute_input": "2021-01-13T19:04:44.180324Z",
"iopub.status.busy": "2021-01-13T19:04:44.138592Z",
"iopub.status.idle": "2021-01-13T19:07:01.833927Z",
"shell.execute_reply": "2021-01-13T19:07:01.833039Z"
}
},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"iter. J_T ∑∫gₐ(t)dt J ΔJ_T ΔJ secs\n",
"0 1.22e-01 0.00e+00 1.22e-01 n/a n/a 8\n",
"1 7.49e-02 2.26e-02 9.75e-02 -4.67e-02 -2.41e-02 23\n",
"2 7.41e-02 3.98e-04 7.45e-02 -8.12e-04 -4.14e-04 23\n",
"3 7.33e-02 3.70e-04 7.37e-02 -7.55e-04 -3.85e-04 23\n"
]
}
],
"source": [
"opt_result = krotov.optimize_pulses(\n",
" objectives,\n",
" pulse_options,\n",
" tlist,\n",
" propagator=krotov.propagators.DensityMatrixODEPropagator(reentrant=True),\n",
" chi_constructor=krotov.functionals.chis_re,\n",
" info_hook=krotov.info_hooks.print_table(J_T=krotov.functionals.J_T_re),\n",
" iter_stop=3,\n",
")"
]
}
],
"metadata": {
"hide_input": false,
"jupytext": {
"encoding": "# -*- coding: utf-8 -*-"
},
"kernelspec": {
"display_name": "Python (krotov)",
"language": "python",
"name": "krotov"
},
"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.1"
},
"toc": {
"base_numbering": 1,
"nav_menu": {},
"number_sections": true,
"sideBar": true,
"skip_h1_title": false,
"title_cell": "Table of Contents",
"title_sidebar": "Contents",
"toc_cell": false,
"toc_position": {},
"toc_section_display": true,
"toc_window_display": true
},
"varInspector": {
"cols": {
"lenName": 16,
"lenType": 16,
"lenVar": 40
},
"kernels_config": {
"python": {
"delete_cmd_postfix": "",
"delete_cmd_prefix": "del ",
"library": "var_list.py",
"varRefreshCmd": "print(var_dic_list())"
},
"r": {
"delete_cmd_postfix": ") ",
"delete_cmd_prefix": "rm(",
"library": "var_list.r",
"varRefreshCmd": "cat(var_dic_list()) "
}
},
"types_to_exclude": [
"module",
"function",
"builtin_function_or_method",
"instance",
"_Feature"
],
"window_display": false
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Display the source blob
Display the rendered blob
Raw
profile.svg
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
@goerz
Copy link
Author

goerz commented Jun 13, 2021

The above profile.svg (click on the link to load the file in its own browser tab!) contains a profiling analysis of just the final optimization from In [35]

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment