{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Optimization with numpy Arrays" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:40.844567Z", "start_time": "2020-03-24T21:13:38.850353Z" }, "attributes": { "classes": [], "id": "", "n": "1" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:24.176583Z", "iopub.status.busy": "2021-11-07T04:51:24.171196Z", "iopub.status.idle": "2021-11-07T04:51:25.406771Z", "shell.execute_reply": "2021-11-07T04:51:25.407034Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Python implementation: CPython\n", "Python version : 3.8.1\n", "IPython version : 7.24.1\n", "\n", "krotov : 1.2.1+dev\n", "numpy : 1.20.3\n", "matplotlib: 3.4.2\n", "scipy : 1.6.3\n", "\n" ] } ], "source": [ "# NBVAL_IGNORE_OUTPUT\n", "%load_ext watermark\n", "import numpy as np\n", "import scipy\n", "import matplotlib\n", "import matplotlib.pylab as plt\n", "import krotov\n", "# note that qutip is NOT imported\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}}$" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `krotov` package heavily builds on QuTiP. However, in rare circumstances\n", "the overhead of `qutip.Qobj` objects might limit numerical efficiency, in\n", "particular when QuTiP's automatic sparse storage is inappropriate. If you know\n", "what you are doing, it is possible to replace `Qobj`s with low-level objects\n", "such as numpy arrays. This example revisits the [Optimization of a\n", "State-to-State Transfer in a Two-Level-System](01_example_simple_state_to_state.ipynb),\n", "but exclusively uses numpy objects for states and operators." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Two-level-Hamiltonian" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We consider again the standard Hamiltonian of a two-level system, but now we\n", "construct the drift Hamiltonian `H0` and the control Hamiltonian `H1` as numpy\n", "matrices:" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:40.862100Z", "start_time": "2020-03-24T21:13:40.847535Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.410921Z", "iopub.status.busy": "2021-11-07T04:51:25.410571Z", "iopub.status.idle": "2021-11-07T04:51:25.412255Z", "shell.execute_reply": "2021-11-07T04:51:25.411969Z" } }, "outputs": [], "source": [ "def hamiltonian(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 * np.array([[1, 0], [0, -1]], dtype=np.complex128)\n", " H1 = np.array([[0, 1], [1, 0]], dtype=np.complex128)\n", "\n", " def guess_control(t, args):\n", " return ampl0 * krotov.shapes.flattop(\n", " t, t_start=0, t_stop=5, t_rise=0.3, func=\"blackman\"\n", " )\n", "\n", " return [H0, [H1, guess_control]]\n" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:40.872139Z", "start_time": "2020-03-24T21:13:40.865354Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.414375Z", "iopub.status.busy": "2021-11-07T04:51:25.414045Z", "iopub.status.idle": "2021-11-07T04:51:25.415659Z", "shell.execute_reply": "2021-11-07T04:51:25.415324Z" }, "lines_to_next_cell": 2 }, "outputs": [], "source": [ "H = hamiltonian()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Optimization target\n", "\n", "By default, the `Objective` initializer checks that the objective is expressed with\n", "QuTiP objects. If we want to use low-level objects instead, we have to\n", "explicitly disable this:" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:40.880764Z", "start_time": "2020-03-24T21:13:40.876732Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.417718Z", "iopub.status.busy": "2021-11-07T04:51:25.417390Z", "iopub.status.idle": "2021-11-07T04:51:25.419081Z", "shell.execute_reply": "2021-11-07T04:51:25.418749Z" } }, "outputs": [], "source": [ "krotov.Objective.type_checking = False" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now, we initialize the initial and target states,\n" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:40.889029Z", "start_time": "2020-03-24T21:13:40.883324Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.421401Z", "iopub.status.busy": "2021-11-07T04:51:25.421076Z", "iopub.status.idle": "2021-11-07T04:51:25.422763Z", "shell.execute_reply": "2021-11-07T04:51:25.422425Z" } }, "outputs": [], "source": [ "ket0 = np.array([[1], [0]], dtype=np.complex128)\n", "ket1 = np.array([[0], [1]], dtype=np.complex128)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "and instantiate the `Objective` for the state-to-state transfer:" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:40.902604Z", "start_time": "2020-03-24T21:13:40.891777Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.426563Z", "iopub.status.busy": "2021-11-07T04:51:25.426113Z", "iopub.status.idle": "2021-11-07T04:51:25.428233Z", "shell.execute_reply": "2021-11-07T04:51:25.427899Z" } }, "outputs": [ { "data": { "text/plain": [ "[Objective[a₀[2,1] to a₁[2,1] via [a₂[2,2], [a₃[2,2], u₁(t)]]]]" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "objectives = [\n", " krotov.Objective(initial_state=ket0, target=ket1, H=H)\n", "]\n", "\n", "objectives" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Note how all objects are numpy arrays, as indicated by the symbol `a`." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Simulate dynamics under the guess field\n", "\n", "To simulate the dynamics under the guess pulse, we can use the objective's\n", "`propagator` method. However, the propagator we use must take into account the\n", "format of the states and operators. We define a simple propagator that solve\n", "the dynamics within a single time step my matrix exponentiation of the\n", "Hamiltonian:" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:40.913174Z", "start_time": "2020-03-24T21:13:40.905293Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.431096Z", "iopub.status.busy": "2021-11-07T04:51:25.430767Z", "iopub.status.idle": "2021-11-07T04:51:25.432460Z", "shell.execute_reply": "2021-11-07T04:51:25.432126Z" } }, "outputs": [], "source": [ "def expm(H, state, dt, c_ops=None, backwards=False, initialize=False):\n", " eqm_factor = -1j # factor in front of H on rhs of the equation of motion\n", " if backwards:\n", " eqm_factor = eqm_factor.conjugate()\n", " A = eqm_factor * H[0]\n", " for part in H[1:]:\n", " A += (eqm_factor * part[1]) * part[0]\n", " return scipy.linalg.expm(A * dt) @ state" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will want to analyze the population dynamics, and thus define the projectors\n", "on the ground and excited levels, again as numpy matrices:\n" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:40.923040Z", "start_time": "2020-03-24T21:13:40.917919Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.434867Z", "iopub.status.busy": "2021-11-07T04:51:25.434532Z", "iopub.status.idle": "2021-11-07T04:51:25.436229Z", "shell.execute_reply": "2021-11-07T04:51:25.435924Z" } }, "outputs": [], "source": [ "proj0 = np.array([[1, 0],[0, 0]], dtype=np.complex128)\n", "proj1 = np.array([[0, 0],[0, 1]], dtype=np.complex128)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will pass these as `e_ops` to the `propagate` method, but since `propagate`\n", "assumes that `e_ops` contains `Qobj` instances, we will have to teach it how to\n", "calculate expectation values:" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:40.930798Z", "start_time": "2020-03-24T21:13:40.926588Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.438416Z", "iopub.status.busy": "2021-11-07T04:51:25.438091Z", "iopub.status.idle": "2021-11-07T04:51:25.439784Z", "shell.execute_reply": "2021-11-07T04:51:25.439450Z" } }, "outputs": [], "source": [ "def expect(proj, state):\n", " return complex(state.conj().T @ (proj @ state)).real" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we can simulate the dynamics over a time grid from $t=0$ to $T=5$ and plot\n", "the resulting dynamics." ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:40.937342Z", "start_time": "2020-03-24T21:13:40.933147Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.441958Z", "iopub.status.busy": "2021-11-07T04:51:25.441636Z", "iopub.status.idle": "2021-11-07T04:51:25.443339Z", "shell.execute_reply": "2021-11-07T04:51:25.443009Z" } }, "outputs": [], "source": [ "tlist = np.linspace(0, 5, 500)" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:41.170271Z", "start_time": "2020-03-24T21:13:40.939839Z" }, "attributes": { "classes": [], "id": "", "n": "12" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.450164Z", "iopub.status.busy": "2021-11-07T04:51:25.449766Z", "iopub.status.idle": "2021-11-07T04:51:25.566157Z", "shell.execute_reply": "2021-11-07T04:51:25.565854Z" } }, "outputs": [], "source": [ "guess_dynamics = objectives[0].propagate(\n", " tlist, propagator=expm, e_ops=[proj0, proj1], expect=expect\n", ")" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:41.178342Z", "start_time": "2020-03-24T21:13:41.172114Z" }, "attributes": { "classes": [], "id": "", "n": "13" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.568718Z", "iopub.status.busy": "2021-11-07T04:51:25.568415Z", "iopub.status.idle": "2021-11-07T04:51:25.570163Z", "shell.execute_reply": "2021-11-07T04:51:25.569855Z" } }, "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": 13, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:41.395497Z", "start_time": "2020-03-24T21:13:41.181086Z" }, "attributes": { "classes": [], "id": "", "n": "14" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.621116Z", "iopub.status.busy": "2021-11-07T04:51:25.610173Z", "iopub.status.idle": "2021-11-07T04:51:25.693710Z", "shell.execute_reply": "2021-11-07T04:51:25.693967Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEGCAYAAABo25JHAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8rg+JYAAAACXBIWXMAAAsTAAALEwEAmpwYAAAkIElEQVR4nO3deXRc53nf8e8DYLDvCzcsBChSlChSEiVQkpfIbhrHsuJIiVPHUpM0iuTIp7VSp2ndo5ymTuy25zhx6pPUdtKoiZtjp7ZqO4tpWZGj1lITx6JEUKTERaJEcwW4ASSIfcfTP94LYAiCwJDAYADc3+ece+7cZe48Q2Luc++7XXN3REQkvrIyHYCIiGSWEoGISMwpEYiIxJwSgYhIzCkRiIjEXE6mA7hW1dXV3tjYmOkwRESWlT179nS4e81M25ZdImhsbKSlpSXTYYiILCtmduJq21Q0JCISc0oEIiIxp0QgIhJzSgQiIjGnRCAiEnNpSwRm9mUzO29mB66y3czsv5nZETN73czuSFcsIiJydem8I/hz4L5Ztn8A2BRNjwN/nMZYRETkKtLWj8Dd/97MGmfZ5UHgKx7Gwd5lZuVmttbdz6Qjnt3HL/IPb7WDGVkGRphnZRkAWWaYMbktvL58btF7E9lZFOXmUJibHU05FOaF18V5ORTn5WBm6fgaIiILLpMdymqBU0nLrdG6KxKBmT1OuGugoaHhuj7s1ROdfOGFIyzG4xdysozywlwqixJhXphLRVGCisJcVpXksaYsn1Wl+awpzaemJI9EtqpqRCRzlkXPYnd/CngKoLm5+bpO5R97zw187D03TByPcYdxdzyaw+XL4w5M7DNt2/DoOP3DY/QPj9I/PEbf0CgDI2P0DY3ROzRCZ/8Il/qHudg3TGf/CEc7euk8OUJn3zCj45eHbwZVRXmsKcujtryAhspC6qOpobKQuooC8nKyr/8fT0RkDplMBG1AfdJyXbQu7cyMbINsFrf4Znzc6ewf5mz3IOe6BznXPcTZrvD6bPcgR9v7ePFwO0Oj40mxwprSfOorC2mqKmLT6mJuWFXMxppiassLJou2RESuVyYTwU7gCTN7Grgb6EpX/cBSkZVlVBXnUVWcxy3rymbcx91p7xni5MX+y6cL/fyfN87xv1umStMKEtlsqCliY5QYNq0u5qY1pTRUFipBiEjK0pYIzOzrwHuBajNrBX4bSAC4+38HngXuB44A/cCvpCuW5cTMWFUa6hCaGyuv2N7ZN8yR9l7ePtfLkfO9HGnvpeV4J9/ed3pyn6LcbG5eW8qWdaVsWVvKzWtL2bymhPyEiphE5Eq23B5e39zc7Bp99Ep9Q6McOd/LG2e6eeNMN4fOdPPGmR56h0aB0BrqhppitqwrZVttGbfXl3PLujIKcpUcROLAzPa4e/NM25ZFZbHMrSgvh9vqy7mtvnxy3fi409o5wKEzXRw63c2hMz28cuzi5N1DdpZx4+oSbq8v49a6cm6rK+fG1cXkqBWTSKzojiCGzvcM8vqpLl5rvcRrrV28duoSXQMjAOQnsti6rozb6svZ3lDOnesrWFtWkOGIRWS+ZrsjUCIQ3J0TF/pDYogSxIG2rsnWS7XlBdyxvoLm9RXcub6Cm9aU6K5BZJlR0ZDMysxorC6isbqIB2+vBWBkbJw3znTTcryTPSc72X3sIt95LRQpFeVmc3tDOXc2VHBnYyXbG8opzU9k8iuILBuDI2Mcv9DHj873cbS9l/M9Q1zoG6Kjd5iB4TGGRscYHh1naHSc0XGPOsGG+ZMfuIkPN9fP9RHXTIlAZpTIzuLWunJurSvnUZpwd9ouDbDnRCd7TnTScryTL75whHEPfR02ry7hjvUV7GisoHl9JXUVBRpmQwTo6B2i5Xgnu49fpOX4RQ6c7mYsqWNpeWGCqqJcqorzqCnJIzc7i7xEFrnZWeRkG0RD3hjQUFmYlhhVNCTXrXdolH0nL9Fy4iJ7TnSy9+SlyVZKa0rzubOxgh3rK2hurOTmtaVkq2+DrHDuzvEL/ZMn/ZbjnRzt6AMgLyeL2+rLaV5fwU1rS7mhpoim6iIKcxfnelxFQ5IWxXk5vHtTNe/eVA3A2Lhz+GwPLScusvt4KE767utnJvfd3lDOjsZKmtdXcHtD+aL9AETSZWRsnEOnu6MTfyctJy7S0TsMhCv95vWVfGRHPc2NlWytLV2yw8XojkDSJrk4aeKHcvhcD+6h6erWdaU0N1ayo7GCO9dXUlOSl+mQRWbVOzTK3pOd7D7eScvxi+w9eYmBkTEgFNs0N1awI/qb3lBdvKR6+KvVkCwZXf0jvHqyc/KuYd+pSwxHrZMaqwonE0NzYyUbqotUzyAZ1d4zRMvxi7wSXcgcOhPK97MMtqwrpXl9ZbjLbaxgdWl+psOdlRKBLFlDo2McaOsO5aknwlVWZ3/o01BZlMudExXQjZVsXVdGbo6arUp6JJfv7z4W/h6PJZXvTxRt7ohaypUss5ZySgSybLg7P2rvuywxHL/QD4Qf4+314cd4Z2MFdzRUUFawvH6MsnSMjo3zxpmecOI/Hu5QO3qHgKny/R2NFexoWhkXIUoEsqyd7xlkz/FQLrvnxFTzu4lmqxO35s2NldSWqxe0zKyrf4S9p0LrtldPdvLqiU76hkP5fm15AXc1hb+juxoruaFmaZXvLwQlAllR+odDs9XdUSuN5B/0urL8y+oZblxdomarMTQ27rx1rmfypL/3ZCc/ag/FPFkGN64uSarYrWRdDC4g1HxUVpTC3BzeubGad24MzVZHx8Z582wPLccvsvtEJy8fu8DOqBd0SX4OdzSEoTG21ZWxrbaM6mK1TlppzncP8npr1+QV/2unLk1eHFQW5bK9vpwP3VHH9vpybq0vpzhPp75kuiOQFcc9jLo60TKp5fhF3jrXO7l9XVk+W2vLuLWujK21ITlUKTksG+e6B9nf2sX+ti4OtIX5+Z5Qtp+dZdy8toQ7GirY3lDO9voK1lcVqvUZuiOQmDGzyec+/+z2OgB6Bkc4eLp78gSyv62Lvzt0bvI9teUFbKstY1tdGTetKWHzmhJqyzVMRiZNDKP+5tluDpzunjzpt0cnfYuesfGujdWTCX1brZ6xcT2UCCQWSvIT3LOhins2VE2u6x4cCSeXpOTw3MGzU+/Jy+HGNSXcFE2b15SyeXUJZYVqqbTQLvQOcfhsD4fP9XD4bA9vnu3h7XM9k8U7Ew9W+rGJk35dGVvWllKkIp4FoaIhkSRdAyO8dS6ciA6f7Z48KfUMjk7us7YsnxtqimmqDmPFNNUUsaG6iNryAg3PPYuRsXFOXeznWEcfxzr6+FF7H8c6ejlyvm+y2SZARWGCzWtKuGlNeMTqjatDItZJf35UNCSSorKCxGRLkgnuzpmuwcmr1cNnezja0cff7Gu7LEEkskOR1IbqItZXhcRQW1EQ5uUFlBcmVnRRk7tzoW+Yts4BWjsHaLvUT1vnAKc6BzjW0cfJi/2XjbpZUZigqbqI926umSyO27ymhJrivBX977QUKRGIzMHMWFdewLryAv7J5lWT692di33DHOvo42h0lXusPcx/cKSDwZHxy45TmJs9eZy1pflUl+RSXZw3OdVEy2UFSythuDvdA6O09w5yvieMm9/eMzQ5ne8ZpO3SAKcvDVzxnUvycqitKODmtSXcv20NTdXFbKgpoqmqiIqi3Ax9I5lOiUDkOpkZVcV5VBXn0Zx0BwFTSeL0pcFwZXxpkNOXBmjrHOB01wBvnOnmYt/wZVfIExLZRkl+gtL8nDAvyKEkL5rnJyhIZJOXk0VuThZ5OVnkJbInx7DPSeozMb3Ud9zDkB5Do+MMjYwxODrO0Mg4g6NjDI2M0zs0QvfAKN2DI2GaeD0wwgxhksi2KIHlsXl1CT++edXkHVBdRSG1FQXq+b1MKBGIpEFykthWVzbjPuPjTmf/MB29w3T0DtHRG66wL/QN0z0wQs9gOBH3DI5yvrt3cnlgZOyKk/x8JLKNvJxsivNyKC3IoTQ/waqSfDbW5FBakKA0P0F5YYKaknDSr4lO/kvtzkWunxKBSIZkZU0li82UpPw+d2d03CcfZxjm4fGGI2Nh6I0Jl73GojuILPJzsslLZJGXk62e16JEILLcmBmJbCORnUWR+sHJAlBbNxGRmFMiEBGJOSUCEZGYUyIQEYk5JQIRkZhTIhARiTklAhGRmEtrIjCz+8zssJkdMbMnZ9jeYGYvmNleM3vdzO5PZzwiInKltCUCM8sGvgR8ANgCPGxmW6bt9lvAN9x9O/AQ8EfpikdERGaWzjuCu4Aj7n7U3YeBp4EHp+3jQGn0ugw4ncZ4RERkBulMBLXAqaTl1mhdst8BftHMWoFngV+b6UBm9riZtZhZS3t7ezpiFRGJrUxXFj8M/Lm71wH3A181syticven3L3Z3ZtramoWPUgRkZUsnYmgDahPWq6L1iV7DPgGgLu/BOQD1WmMSUREpklnItgNbDKzJjPLJVQG75y2z0ngnwKY2c2ERKCyHxGRRZS2RODuo8ATwPeANwitgw6a2WfM7IFot38L/KqZvQZ8HXjEfSEfuSEiInNJ6/MI3P1ZQiVw8rpPJb0+BLwrnTGIiMjsMl1ZLCIiGaZEICISc0oEIiIxp0QgIhJzSgQiIjGnRCAiEnNKBCIiMadEICISc0oEIiIxp0QgIhJzSgQiIjGnRCAiEnNKBCIiMadEICISc0oEIiIxp0QgIhJzSgQiIjGnRCAiEnNKBCIiMadEICISc0oEIiIxp0QgIhJzSgQiIjGnRCAiEnNKBCIiMadEICISc0oEIiIxp0QgIhJzSgQiIjGXk+qOZpYNrE5+j7ufTEdQIiKyeFK6IzCzXwPOAc8D342mZ1J4331mdtjMjpjZk1fZ5+fN7JCZHTSzr11D7CIisgBSvSP4BLDZ3S+keuDoDuJLwPuAVmC3me1090NJ+2wCfhN4l7t3mtmq1EMXEZGFkGodwSmg6xqPfRdwxN2Puvsw8DTw4LR9fhX4krt3Arj7+Wv8DBERmadU7wiOAi+a2XeBoYmV7v75Wd5TS0ggE1qBu6ftcyOAmf0jkA38jrs/N/1AZvY48DhAQ0NDiiGLiEgqUk0EJ6MpN5oW8vM3Ae8F6oC/N7Nt7n4peSd3fwp4CqC5udkX8PNFRGIvpUTg7p8GMLPiaLk3hbe1AfVJy3XRumStwMvuPgIcM7O3CIlhdypxiYjI/KXaamirme0FDgIHzWyPmd0yx9t2A5vMrMnMcoGHgJ3T9vkbwt0AZlZNKCo6mnr4IiIyX6kWDT0F/Ia7vwBgZu8F/gfwzqu9wd1HzewJ4HuE8v8vu/tBM/sM0OLuO6NtP2lmh4Ax4JPX0jJJRGQxjYyM0NrayuDgYKZDuar8/Hzq6upIJBIpv8fc5y5yN7PX3P22udYthubmZm9paVnsjxUR4dixY5SUlFBVVYWZZTqcK7g7Fy5coKenh6ampsu2mdked2+e6X2pNh89amb/0cwao+m3UBGOiMTM4ODgkk0CAGZGVVXVNd+xpJoIHgVqgL+KppponYhIrCzVJDDheuJLKRG4e6e7/2t3vyOaPjHRCUxERBbPc889x+bNm9m4cSOf/exnF+SYs1YWm9kfuPuvm9l3gCsqE9z9gQWJQkRE5jQ2NsbHP/5xnn/+eerq6tixYwcPPPAAW7Zsmddx52o19NVo/vvz+hQREZm3V155hY0bN7JhwwYAHnroIb797W+nNxG4+57o5e3u/ofJ28zsE8D/m9eni4gsU5/+zkEOne5e0GNuWVfKb//01btotbW1UV8/1U+3rq6Ol19+ed6fm2pl8S/PsO6ReX+6iIhk3Fx1BA8D/xxoMrPkXsElwMV0BiYispTNduWeLrW1tZw6NTWWZ2trK7W1tfM+7lx1BD8EzgDVwH9NWt8DvD7vTxcRkZTt2LGDt99+m2PHjlFbW8vTTz/N1742/+d5zVVHcAI4Abxj3p8kIiLzkpOTwxe/+EXe//73MzY2xqOPPsott8z/ziSlsYbM7B7gC8DNhGGos4E+dy+ddwQiIpKy+++/n/vvv39Bj5lqZfEXgYeBt4EC4KOEx1CKiMgyl2oiwN2PANnuPubu/xO4L31hiYjIYkl1GOr+6JkC+8zs9wgVyCknERERWbpSPZn/EqFe4Amgj/DksZ9LV1AiIrJ4Un1U5Yno5QDw6fSFIyIii22uDmX7mWGwuQnufuuCRyQiIotqrjuCDy5KFCIikpJHH32UZ555hlWrVnHgwIEFOeasdQTufmK2aUEiEBGRlD3yyCM899xzC3rMlCqLzazHzLqjadDMxsxsYYfdExGROd17771UVlYu6DFTrSwumXht4TloDwL3LGgkIiLLyd8+CWf3L+wx12yDDyzMU8euxTX3BfDgb4D3L3w4IiKy2FIda+hDSYtZQDMwmJaIRESWgwxcuadLqj2Lfzrp9ShwnFA8JCIiy1yqdQS/ku5ARERkbg8//DAvvvgiHR0d1NXV8elPf5rHHntsXsdMtWhoA/CHhApiB14C/o27H53Xp4uIyDX5+te/vuDHTLWy+GvAN4C1wDrgm8DCRyMiIosu1URQ6O5fdffRaPoLID+dgYmIyOJItbL4b83sSeBpQtHQR4BnzawSwN31IHsRkWUq1UTw89H8Y9PWP0RIDBsWLCIRkSXM3Qn9apcm96uOE3pVqbYaarrmI4uIrDD5+flcuHCBqqqqJZkM3J0LFy6Qn39tJfepthpKAP8SuDda9SLwJ+4+Msf77iO0NsoG/tTdZ+yBYWY/B3wL2OHuLamFLiKyuOrq6mhtbaW9vT3ToVxVfn4+dXV11/SeVIuG/hhIAH8ULf9StO6jV3uDmWUTHnD/PqAV2G1mO9390LT9SoBPAC9fU+QiIosskUjQ1LTyCkhSTQQ73P22pOXvm9lrc7znLuDIRF8DM3ua0Bv50LT9/hPwu8AnU4xFREQWUKrNR8fM7IaJhaiD2dgc76kFTiUtt0brJpnZHUC9u393tgOZ2eNm1mJmLUv5lkxEZDlK9Y7gk8ALZjbRk7gRmNewE2aWBXweeGSufd39KeApgObm5muvEhcRkatK9Y7gH4E/AcaBi9Hrl+Z4TxtQn7RcF62bUAJsBV40s+OE4St2mllzijGJiMgCSDURfAVoIpTnf4HQb+Crc7xnN7DJzJrMLJfQ52DnxEZ373L3andvdPdGYBfwgFoNiYgsrlSLhra6+5ak5RfMbHql72XcfdTMngC+R2g++mV3P2hmnwFa3H3nbO8XEZHFkWoieNXM7nH3XQBmdjcw55W7uz8LPDtt3aeusu97U4xFREQWUKqJ4E7gh2Z2MlpuAA6b2X7C0ytvTUt0IiKSdqkmgvvSGoWIiGRMqmMNnUh3ICIikhmpthoSEZEVSolARCTmlAhERGJOiUBEJOaUCEREYk6JQEQk5pQIRERiTolARCTmlAhERGJOiUBEJOaUCEREYk6JQEQk5pQIRERiTolARCTmlAhERGJOiUBEJOaUCEREYk6JQEQk5pQIRERiTolARCTmlAhERGJOiUBEJOaUCEREYk6JQEQk5pQIRERiTolARCTmlAhERGIurYnAzO4zs8NmdsTMnpxh+2+Y2SEze93M/q+ZrU9nPCIicqW0JQIzywa+BHwA2AI8bGZbpu22F2h291uBbwG/l654RERkZum8I7gLOOLuR919GHgaeDB5B3d/wd37o8VdQF0a4xERkRmkMxHUAqeSllujdVfzGPC3M20ws8fNrMXMWtrb2xcwRBERWRKVxWb2i0Az8LmZtrv7U+7e7O7NNTU1ixuciMgKl5PGY7cB9UnLddG6y5jZTwD/AXiPuw+lMR4REZlBOu8IdgObzKzJzHKBh4CdyTuY2XbgT4AH3P18GmMREZGrSFsicPdR4Ange8AbwDfc/aCZfcbMHoh2+xxQDHzTzPaZ2c6rHE5ERNIknUVDuPuzwLPT1n0q6fVPpPPzRURkbkuislhERDJHiUBEJOaUCEREYk6JQEQk5pQIRERiLq2thkRkEY2PhSmZWdLrLMjKXtyYZFlQIhDJtLER6GuH3vNh3tcOg10w2B3mQ12XL48ORtNwmI9F8/HRuT/LsiEnHxL5YZ6TBzkFYZ5XAvllSVP51OvCSiheBcWrobAasnXqWEn0vymSTsN9cOkUdLVC16loaoXu09GJ/zwMdF79/Ymi6GRcGuZFNZAoiE7iudE8H7Kj15dd8fvlx3KPEshQmI8MTiWVkQEY7oWOc1HS6YKRfmZmUFgVkkJxDRSvgbI6KK+HsnooXx+WE/nz/deTRaJEIDJf/RfhwpHLp4vHwkl/+knesqG0FkrXQc2N0PjucKVdVBPNV0FRNRRUQF5pZq+8R4dhKLoL6esISav3HPS2R/NoueMH0HMafPzy9xetCsmhvAGqNl4+FZRn5CvJzJQIRFIxPg5dJ+HcITh/6PKTfvLJ3rKhohEqN0Bdc7hCLquPrpbroGTt8imnz8mFnOqQmKpumH3fsdGQDC5Fdz2XToap6xSc3geHdoIn1V8UVk8lhepNsGoLrN4SkmRyvYYsCiUCken6LsD5g9FJP5q3vxmKTiaU1oaT4y0/e/mVbnkDZCcyF3umZOeE717eMPP20WHoPH7lndOR52HfX0ztl1cGq24OSWHVFlh9S1guqFiUrxFXSgQSbz1n4fTecNV6ei+ceQ16z05tL6gMJ6PbfyE6Od0Cq24KFauSupzcUBRWc+OV2wY64fwbcO5gmJ8/BPv/Eoa+PLVPWQOsuw3W3g7rboe126GoarGiX/GUCCQ+es9PO+nvg54z0UaDms2w4b3hxL96C6zeGipEVVSRXgUVsP6dYZrgDt1tU3dlZ14L/29vfGdqn7J6WHvbVGJYd3soxpJrpkQgK9PYKJx9HU7ugpMvQduecGIBwKD6Rmi6F9ZtD1eZa7ZBXnEmI5ZkZqFOpawObvzJqfUDl0JSOLMvJIYz++DNZ6a2VzRC3Q6ouyvU0azZFs+iumukRCArw1APtO6OTvy7oLUFRvrCtvKGcLU5cdJfe6uKdpargnLY8J4wTRi4FJL+6b3hb+DYP8D+b4ZtOfnh/72uOUoOO6B0bSYiX9LM3efeawlpbm72lpaWTIchmdZ9Jlzpn9wFp3bB2f2h+aJlhSKdhndAwz1hKl2X6WhlMU0UK516JVwQtL4S7iLGhsP20jqo3wH198D6d4S/l+XSkmsezGyPuzfPtE13BLL0jY9Dx+HoxP9ymF86EbYlCsPV3r2fDCf92ubQ+UriK7lYaeuHwrrRITjzerhjaN0dksTBvw7bckug/q6QFBreAbV3hk57MaJEIEvPyGC4zT+1a6qoZ/BS2Fa0Kpzw7/5YmK+5VWXAMrecvOguYMfUukunoouLl+DES/D9/xzWZyVCcdJEYqi/OwyxsYKpaEgyr/8inHp56qR/+tWp2/iqTVERT1TUU7lBrXgkPSb/DqPEcHovjI+EbTU3TyWGhneEDoLLzGxFQ0oEsrjcQ8eiiR/cyV2hsxZEV2K3T5346+9Wc0DJnJGB0NpsIjGcegWGe8K20rrwdzqRHGpuhqylPaq/6ggkc8ZG4dyBqWacJ3dNddjKKwtls9s+HJXN3hG7sllZwhIFYSyoxneH5fGxqb/lEz+E4z+AA98K2/LKoOHuqYuYdXcsq0H3lAhkYQ31QlvL1Im/tWVqaIayemj6sfBjqb8nDB0Qg9YaskJkZYcObGtvC3VU7qHRQvJFztt/F/bNzg31DMl3t0u4nkFFQzI/XW2hUvfUK+GHcHZ/NLiYRc047576ISzDclWRa5Jcz3ByF7S9mlTPcNPURdDqW8LYVLmFixaa6ghkYYyNhu7+J1+eOvl3nQrbcgpCs7vJK6AdYfx8kTgbGQiVzhOJ4eTL4UFDE8rqQz+Xopow5RZFDwuKnjGRlRM1jrAwb7o3JJHroDoCuT6DXVNtrk/uChVnE8U8JWvDVf47Ph7K+dWMU+RKiYLLx1Ga6BPT/iZ0vB2m3rNw8Wi4kxjuDw8K8rGZj/dTn7/uRDAbJQIJxkbC6I+nXw0n/La9YRRIPOqtewvc9nB0a3tXuJJRM06Ra5OVFerGVt08+35jo1MJwR3wME+kpyhJiSCO3MMVSNueUIbZtieM1TI6GLYXVIZini0PhKv+umaNzSOymLJzIHvxBkFUIljpxkbDA0DO7odz+0M3+9N7p3rqJgrDQGw7Phqab9beGZ45q6t9kdhQIlhJBrtD8c65A+EK/+z+8KCPiSv97NzQcuGWnwkn/HV3hOVMPhdXRDJOZ4DlqK8D2g9HlU5vTc27W6f2KagMY7Hv+GiYr9kWxuBXha6ITKNEsBS5h8f3dR4LwzFMTB1Hwkm//8LUvonC8PDv9e8MjwFcc2tov1+6TsU7IpKStCYCM7sP+EMgG/hTd//stO15wFeAO4ELwEfc/Xg6Y1oSRgah53R4Xm736fC4xO4z0HUyOumfgKHuy99TWB06oNz0wfBIxerN4cRfWrfkxzgRkaUtbYnAzLKBLwHvA1qB3Wa2090PJe32GNDp7hvN7CHgd4GPpCumBTU2Gp6ANTxtGuqBgYuhh2H/hanXA51huedMeD1dojCMn17RFDpkVTROTeXr9RhFEUmbdN4R3AUccfejAGb2NPAgkJwIHgR+J3r9LeCLZmaeju7Or34VfviF8BQrPMx9PBTDTLbTnVgen7afJy17eMjF2NDcn2nZYXyRwqpQZl+5IRThlKwNU+laKFkHJWtCL1wV5YhIBqQzEdQCp5KWW4G7r7aPu4+aWRdQBXQk72RmjwOPAzQ0NFxfNIWVoROHWeggRTS3rGnrbGr5sv2S9snJhdzicBWfWzQ1JQpDe/vCynDi18ldRJaBZVFZ7O5PAU9BGGvoug5y00+FSURELpPOWsY2IHm4ybpo3Yz7mFkOUEaoNBYRkUWSzkSwG9hkZk1mlgs8BOycts9O4Jej1/8M+H5a6gdEROSq0lY0FJX5PwF8j9B89MvuftDMPgO0uPtO4M+Ar5rZEeAiIVmIiMgiSmsdgbs/Czw7bd2nkl4PAh9OZwwiIjI79UQSEYk5JQIRkZhTIhARiTklAhGRmFt2D683s3bgxHW+vZppvZZjQN85HvSd42E+33m9u9fMtGHZJYL5MLMWd2/OdByLSd85HvSd4yFd31lFQyIiMadEICISc3FLBE9lOoAM0HeOB33neEjLd45VHYGIiFwpbncEIiIyjRKBiEjMxSYRmNl9ZnbYzI6Y2ZOZjifdzOzLZnbezA5kOpbFYmb1ZvaCmR0ys4Nm9olMx5RuZpZvZq+Y2WvRd/50pmNaDGaWbWZ7zeyZTMeyGMzsuJntN7N9Ztay4MePQx2BmWUDbwHvIzwyczfwsLsfmvWNy5iZ3Qv0Al9x962ZjmcxmNlaYK27v2pmJcAe4GdW+P+zAUXu3mtmCeAHwCfcfVeGQ0srM/sNoBkodfcPZjqedDOz40Czu6elA11c7gjuAo64+1F3HwaeBh7McExp5e5/T3jGQ2y4+xl3fzV63QO8QXgu9orlQW+0mIimFX11Z2Z1wE8Bf5rpWFaKuCSCWuBU0nIrK/wEEXdm1ghsB17OcChpFxWT7APOA8+7+0r/zn8A/HtgPMNxLCYH/s7M9pjZ4wt98LgkAokRMysG/hL4dXfvznQ86ebuY+5+O+G54HeZ2YotCjSzDwLn3X1PpmNZZO929zuADwAfj4p+F0xcEkEbUJ+0XBetkxUmKif/S+B/uftfZTqexeTul4AXgPsyHEo6vQt4ICozfxr4cTP7i8yGlH7u3hbNzwN/TSjuXjBxSQS7gU1m1mRmuYRnI+/McEyywKKK0z8D3nD3z2c6nsVgZjVmVh69LiA0iHgzo0Glkbv/prvXuXsj4Xf8fXf/xQyHlVZmVhQ1fsDMioCfBBa0NWAsEoG7jwJPAN8jVCB+w90PZjaq9DKzrwMvAZvNrNXMHst0TIvgXcAvEa4S90XT/ZkOKs3WAi+Y2euEC57n3T0WTSpjZDXwAzN7DXgF+K67P7eQHxCL5qMiInJ1sbgjEBGRq1MiEBGJOSUCEZGYUyIQEYk5JQIRkZhTIhCZhZmVm9m/il6vM7NvZTomkYWm5qMis4jGLHomLiO4SjzlZDoAkSXus8AN0aBubwM3u/tWM3sE+BmgCNgE/D6QS+jQNgTc7+4XzewG4EtADdAP/Kq7r9iev7I8qWhIZHZPAj+KBnX75LRtW4EPATuA/wL0u/t2Qo/ufxHt8xTwa+5+J/DvgD9ajKBFroXuCESu3wvRcw96zKwL+E60fj9wazQK6juBb4ZhkADIW/wwRWanRCBy/YaSXo8nLY8TfltZwKXobkJkyVLRkMjseoCS63lj9CyEY2b2YQijo5rZbQsZnMhCUCIQmYW7XwD+0cwOAJ+7jkP8AvBYNHLkQVb4I1JleVLzURGRmNMdgYhIzCkRiIjEnBKBiEjMKRGIiMScEoGISMwpEYiIxJwSgYhIzP1/3c4T64tgAqMAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "plot_population(guess_dynamics)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The result is the same as in the original example.\n", "\n", "## Optimize\n", "\n", "First, we define the update shape and step width as before:" ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:41.404539Z", "start_time": "2020-03-24T21:13:41.397827Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.696457Z", "iopub.status.busy": "2021-11-07T04:51:25.696157Z", "iopub.status.idle": "2021-11-07T04:51:25.697837Z", "shell.execute_reply": "2021-11-07T04:51:25.697530Z" } }, "outputs": [], "source": [ "def S(t):\n", " \"\"\"Shape function for the field update\"\"\"\n", " return krotov.shapes.flattop(\n", " t, t_start=0, t_stop=5, t_rise=0.3, t_fall=0.3, func='blackman'\n", " )" ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:41.412374Z", "start_time": "2020-03-24T21:13:41.406924Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.700189Z", "iopub.status.busy": "2021-11-07T04:51:25.699889Z", "iopub.status.idle": "2021-11-07T04:51:25.701588Z", "shell.execute_reply": "2021-11-07T04:51:25.701166Z" } }, "outputs": [], "source": [ "pulse_options = {\n", " H[1][1]: dict(lambda_a=5, update_shape=S)\n", "}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we can run the optimization with only small additional adjustments. This is\n", "because Krotov's method internally does very little with the states and\n", "operators: nearly all of the numerical effort is in the propagator, which we\n", "have already defined above for the specific use of numpy arrays.\n", "\n", "Beyond this, the optimization only needs to know three things: First, it must\n", "know how to calculate and apply the operator $\\partial H/\\partial \\epsilon$. We\n", "can easily teach it how to do this:" ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:41.423350Z", "start_time": "2020-03-24T21:13:41.415715Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.703907Z", "iopub.status.busy": "2021-11-07T04:51:25.703606Z", "iopub.status.idle": "2021-11-07T04:51:25.705319Z", "shell.execute_reply": "2021-11-07T04:51:25.705015Z" } }, "outputs": [], "source": [ "def mu(objectives, i_objective, pulses, pulses_mapping, i_pulse, time_index):\n", " def _mu(state):\n", " return H[1][0] @ state\n", " return _mu" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Second, the pulse updates are calculated from an overlap of states, and we\n", "define an appropriate function for numpy arrays:" ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:41.432372Z", "start_time": "2020-03-24T21:13:41.427434Z" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.707441Z", "iopub.status.busy": "2021-11-07T04:51:25.707142Z", "iopub.status.idle": "2021-11-07T04:51:25.708795Z", "shell.execute_reply": "2021-11-07T04:51:25.708495Z" } }, "outputs": [], "source": [ "def overlap(psi1, psi2):\n", " return complex(psi1.conj().T @ psi2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Third, it must know how to calculate the norm of states, for which we can use `np.linalg.norm`." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "By passing all these routines to `optimize_pulses`, we get the exact same\n", "results as in the original example, except much faster:" ] }, { "cell_type": "code", "execution_count": 18, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:45.052773Z", "start_time": "2020-03-24T21:13:41.436123Z" }, "attributes": { "classes": [], "id": "", "n": "15" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:25.711454Z", "iopub.status.busy": "2021-11-07T04:51:25.711153Z", "iopub.status.idle": "2021-11-07T04:51:27.478073Z", "shell.execute_reply": "2021-11-07T04:51:27.478319Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "iter. J_T ∫gₐ(t)dt J ΔJ_T ΔJ secs\n", "0 9.51e-01 0.00e+00 9.51e-01 n/a n/a 0\n", "1 9.24e-01 1.20e-02 9.36e-01 -2.71e-02 -1.50e-02 0\n", "2 8.83e-01 1.83e-02 9.02e-01 -4.11e-02 -2.28e-02 0\n", "3 8.23e-01 2.71e-02 8.50e-01 -6.06e-02 -3.35e-02 0\n", "4 7.37e-01 3.84e-02 7.76e-01 -8.52e-02 -4.68e-02 0\n", "5 6.26e-01 5.07e-02 6.77e-01 -1.11e-01 -6.05e-02 0\n", "6 4.96e-01 6.04e-02 5.56e-01 -1.31e-01 -7.02e-02 0\n", "7 3.62e-01 6.30e-02 4.25e-01 -1.34e-01 -7.09e-02 0\n", "8 2.44e-01 5.65e-02 3.00e-01 -1.18e-01 -6.15e-02 0\n", "9 1.53e-01 4.39e-02 1.97e-01 -9.03e-02 -4.64e-02 0\n", "10 9.20e-02 3.02e-02 1.22e-01 -6.14e-02 -3.12e-02 0\n" ] } ], "source": [ "opt_result = krotov.optimize_pulses(\n", " objectives,\n", " pulse_options=pulse_options,\n", " tlist=tlist,\n", " propagator=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.check_monotonic_error,\n", " iter_stop=10,\n", " norm=np.linalg.norm,\n", " mu=mu,\n", " overlap=overlap,\n", ")" ] }, { "cell_type": "code", "execution_count": 19, "metadata": { "ExecuteTime": { "end_time": "2020-03-24T21:13:45.060650Z", "start_time": "2020-03-24T21:13:45.054855Z" }, "attributes": { "classes": [], "id": "", "n": "16" }, "execution": { "iopub.execute_input": "2021-11-07T04:51:27.481172Z", "iopub.status.busy": "2021-11-07T04:51:27.480688Z", "iopub.status.idle": "2021-11-07T04:51:27.482856Z", "shell.execute_reply": "2021-11-07T04:51:27.482489Z" } }, "outputs": [ { "data": { "text/plain": [ "Krotov Optimization Result\n", "--------------------------\n", "- Started at 2021-11-07 05:51:25\n", "- Number of objectives: 1\n", "- Number of iterations: 10\n", "- Reason for termination: Reached 10 iterations\n", "- Ended at 2021-11-07 05:51:27 (0:00:02)" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "opt_result" ] } ], "metadata": { "hide_input": false, "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.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": false }, "varInspector": { "cols": { "lenName": 16, "lenType": 16, "lenVar": 40 }, "kernels_config": { "python": { "delete_cmd_postfix": "", "delete_cmd_prefix": "del ", "library": "var_list.py", "varRefreshCmd": "print(var_dic_list())" }, "r": { "delete_cmd_postfix": ") ", "delete_cmd_prefix": "rm(", "library": "var_list.r", "varRefreshCmd": "cat(var_dic_list()) " } }, "types_to_exclude": [ "module", "function", "builtin_function_or_method", "instance", "_Feature" ], "window_display": false } }, "nbformat": 4, "nbformat_minor": 2 }