{ "cells": [ { "cell_type": "markdown", "id": "26bccb0d", "metadata": {}, "source": [ "# Simulating randomized benchmarking\n", "\n", "In this example, we will reproduce a randomized benchmarking experiment used in Figure 3a of [Piltz et. al.](https://www.nature.com/articles/ncomms5679?origin=ppub)\n", "\n", "Note: This example is quite computationally expensivem, hence for the full simulation we use [joblib](https://joblib.readthedocs.io/en/latest/) for parallel computing. However, you don't need this to run the demo.\n" ] }, { "cell_type": "code", "execution_count": 1, "id": "355292f1", "metadata": {}, "outputs": [], "source": [ "import matplotlib.pyplot as plt\n", "from scipy.optimize import curve_fit\n", "\n", "import numpy as np\n", "from qutip import sigmax, sigmay, sigmaz, basis, qeye, tensor, Qobj, fock_dm\n", "from qutip_qip.circuit import QubitCircuit, Gate\n", "from qutip_qip.device import ModelProcessor, Model\n", "from qutip_qip.compiler import GateCompiler, Instruction\n", "from qutip import Options\n", "from qutip_qip.noise import Noise" ] }, { "cell_type": "markdown", "id": "43a3e654", "metadata": {}, "source": [ "We build a two-qubit Processor, where the second qubit is detuned from the first one by $\\delta= 1.852$MHz. A sequence of $\\pi$-pulses with Rabi frequency of $\\Omega= 20$KHz and random phases are applied to the first qubit. We define noise such that the same pulse also applies to the second qubit. Because of the detuning, this pulse does not flip the second qubit but subjects it to a diffusive behaviour, so that the average fidelity of the second qubit with respect to the initial state decreases.\n", "\n", "Here, we reproduce these results with a two-qubit `Processor`.\n", "We start with an initial state of fidelity 0.975 and simulate the Hamiltonian\n", "\\begin{align}\n", "H=\\Omega(t)(\\sigma^x_0 + \\lambda \\sigma^x_1) + \\delta\\sigma^z_1\n", ",\n", "\\end{align}\n", "where $\\lambda$ is the ratio between the cross-talk pulse's amplitudes." ] }, { "cell_type": "markdown", "id": "3afdd083", "metadata": {}, "source": [ "In the cell below, we first build a Hamiltonian model called `MyModel`.\n", "For simplicity, we only include two single-qubit control Hamiltonians: $\\sigma_x$ and $\\sigma_y$.\n", "We then define the compiling routines for the two types of rotation gates RX and RY.\n", "In addition, we also define a rotation gate with mixed X and Y quadrature, parameterized by a phase $\\phi$, $\\cos(\\phi)\\sigma_x+\\sin(\\phi)\\sigma_y$.\n", "This will be used later in the example of custom noise.\n", "\n", "We then initialize a `ModelProcessor` with this model.\n", "In the `ModelProcessor`, the default simulation workflow is already defined, such as the `load_circuit` method.\n", "Since rotations around the $x$ and $y$ axes are the native gates of our hardware, we define them in the attribute \\texttt{native\\_gates}.\n", "Providing this native gates set, rotation around $z$ axis will be automatically decomposed into rotations around $x$ and $y$ axes.\n", "We define a circuit consisting of $\\pi/2$ rotation followed by a Z gate.\n", "The compiled pulses are shown in \\cref{fig:customize pulse}, where the Z gate is decomposed into rotations around $x$ and $y$ axes." ] }, { "cell_type": "code", "execution_count": 2, "id": "17cd84d0", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAksAAAFcCAYAAADYhvq8AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAABJ0AAASdAHeZh94AAAYiUlEQVR4nO3df6yd9X0f8PfXv01NnTrkh0nDAlTeUpeKteBobdeM1aRUsK5KhjZpYxqbhki0lImsDDVblmwaWwtMLa1om0gVyjZtQ2WaNLqKpVTqH1lXDxelxHJFaJrECU75YUywsY3j++yPc+/1OedenmP7xs/zvfd5vf6xOffecz/4PM/5vs/3+3m+T2maJgAALG9d3wUAANRMWAIAaCEsAQC0EJYAAFoISwAALYQlAIAWwhIAQAthCQCghbAEANBCWAIAaCEsAQC02NB3AfPcoA4A6FI51280swQA0EJYAgBoISwBALQQlgAAWtTS4H3RNE2T33/2xTz3wrG+SyHJu3dckr3vfUfWrzvnvjqATuz/6it5+muv9F0GSX5819uy6x2X9l3GotI0VVyIdtGK+KOvvZIPPvx/LtbTcwE+fdsP5wO739l3GQCL/vxbJ/Mj//73cmauijFx8LZv3Zh9H/+JbN6w/mL+GlfDLTh05PW+S2DK17wmQGW+cfSEoFSRV0+czrdOfLvvMhat+WW4cde/53tyycaLmlJ5E9+ea/L5P305SVLHZCbAWeOrLNe8a3t2XLKxx2qG6/lXT+ZL820zlax8JRlYWFpfSjasX/OTaZWa67sAgDc1Pqm0YZ2xoi/j/aw1TfQ5GuhcY8N2oDIVTWIwr6axQliiG2NtdN6UgNrMjb0xuVa3P+P/9maWGJwydgpUdPwDJJkMS9JSHeYqSkvCEp0zswTURlaqQ6l0FUJYonM1rUMDJFMzS1ShptdEWKJzFR3/AEkm+2NKMbfUl/GWDWGJwfHWA9SspoF50MYGi4paloQlOjKxDl3RGQCQTFx54sNdfyb/7esZK4QlOicrAbVxNVx9zCwxaBUd/wBJpnqW+itj8Cb3WapntBCW6MT4CVDR8Q+QpK6BedDGe5YqukuWsETnbB0A1Gayl9LcUl9cDQfzKjr+AZJMbx3QXx2cVdNYISwBMHg1zWIM2cQO3hWtQghLdGJ8k7d6Dn+AERfD1cfVcAybT3BAZcZnlizD9cfVcDCvnsMfYKSicXnYKt3AWFiicxUd/wBJpmaWLMT1aPxquB7LmCIs0ZmFU6Cmpj2ApK6BecgmluEqelGEJTpnZgmojdud1KeirCQs0b2Kjn+AJJP9MbJSf2wdAPMngZkloDbel+pT02siLNEZPUtArezgXQdbB8CCeo5/gCTTA7O01JfxDYz1LDFoFR3/AEn0LNXIzBKDVtNGYwBJXbMYjNQ0VghLdK6i4x8gidud1GL8n76msUJYojN2xQVqZWapEmPDRE2vibBEdxa2Dui3CoAl9CzVwdVwMK+i4x8gydT7krTUo7P/+HqWGDT7LAG1cSPdOhTLcAzd4qaUFZ0AAEldAzMjluEAoCJzepaqMNmz1FsZSwhLdK6mdWiAZOp9SVqqQk1jhbBE5+o5/AFG5mSlKoz3LFWUlYQlurNwEtR0AgAk3pdqpGeJQXM1HFCbyR28zS31ZfxKRD1LDFpFHxYAktTVH8OImSUGrZ7DH2BkYRbDnFK/JnuW6hkthCU6V9HxD5CkrlkMRmp6SYQlOmNXXKBWizNL3qZ6ZZ8lWFTRGQCQupZ8Bm3idif1vCbCEt2xdQBQqYW3JTPg/SpupMvQuTccUKs5Hd7VsQzHoNlnCaiNrFSHYhkORio6/gGS1DUwM2JmiUGr6PgHSFJXf8yQjc/s1fSaCEt0Rs8SUCtbB9SnprGiirB08ODB3HTTTdm2bVt27NiR2267LS+++GLfZXGR6FkCamMZrg7j9+Wr6TXZ0HcBpZTvveyyy7J9+/bcd999OXbsWB544IE888wz2bdvXzZt2tR3iXynLE4t9VoFwBK2DqhPTT1LvYelJD9//Pjx7N+/P1dccUWSZM+ePbnxxhvzyCOP5I477ui5PL7TKjr+AZKc7Y+xDNevyR286xktLngZrpRyaSnll0opXymlnCqlvFBK+Vwp5YdKKe8tpZwopXx26md+rJRyppTyC2MPf+iWW25ZDEpJsnfv3uzatSuPPvrohZZHhc72LNVzAgAkydxc3xWQZCIt1TRWrKRn6deTfDjJY0k+kuSBJCeSvLdpmoNJ/mWS20opP50kpZTvSvJIkj9J8on5x96V5O3XXXfdkiffs2dPnn766RWUR31GZ0E9hz/AyMIshomlelSUlVa0DHdzks80TfOxscd+cezv/yHJ30zy6VLK55N8KslfSPJXmqY5Nf89O5Nk586dS558586dOXLkSE6dOpXNmzevoExqU9MJAJDU1R8zZGvxRrpHk7yvlHL5cl9smmYuyT9Isi3J72Q0+/TvmqZ5auzbtiZZNgxt2bIlSXLixIkVlAgAsy0u+Zhaqsaa6FlKck+SH0hyqJSyr5TyyVLKVePf0DTNnyb5ZJLrkxxI8m+mnuNEkpw6dSrTTp48mSTZunXrCkqkRvUc/gAjZ5fhpKU+jW8dsCZ6lpqmeTTJVUk+muT5JD+X5EAp5aemvvUD839enuStU187nCSHDx9e8vyHDx/Ojh07LMGtIQvnQE0nAEDiQ1xNFuLSWlmGS9M0h5umebhpmp9JcmWSl5N8fOHrpZQ7k9w4/9imJL8x9fPfSPLiU0+Nr8yN7Nu3L9dee+1KyqNSFR3/AEns4F2jVb8MV0pZX0rZPv5Y0zQvZDTDtHn+e65Mcn+Sx5qmuS/JP0vy06WUvz/1dI89/vjjOXTo0OIDTz75ZJ599tnceuutF1Ietavn+AdIUtfAPHjzgbWmmaULvRru0iRfL6X8VpIvJDmWZG9GvUkfK6NFx9/MqCfpw0nSNM1vlFI+lOSXSym/2zTN8/PPdd8ll1xy5w033JC77rorx44dy/33359rrrkmt99++0r+36iU250AtWlsHVCNktFn6prGigsNS68neTijfqQPZjRD9VySjzRN82ullJ9N8teSfKhpmvGbvP2jJF9M8pmMth5I0zSHDhw4kLvvvjv33ntvNm3alJtvvjkPPvigfqU1xo10gVotbkopLVWjprHigsJS0zRvZHQ13D1v8vWHkjy0zOOHkmyffnz37t154oknLqQUVqGaTgCAxNVwNSmlJE2TuYrW4VbU4A3nZeFquIqmVgGSuvpjGKnpNanhRroMhGW4enzr5Ok88cVv5sTpM32XAlU4dOT1JFbhanB264B6BotBhaW3fffmbN+6se8yBmvD+nVJzuTE6TP5xlE7s/fpU//zQP73gT/vuwyozsYN67LzLVv6LmPQ1pVRi3dNe/INKiz90727cvXbtvVdxmDt/+rRHH39dA4e/lY+8T++2Hc5g/b/vnKk7xKgOutKcuf7r8odP35136UM2m//8eG88fqcZbgu/eTud2b/v9ibJGaVYN7Ce9CPft9b89Df+cu91gK12LRhXS7dYpzo2+fufn9Kki0b1/ddyqI1H5a2bFxf1T/4kC32AlT0aWGw5l+DLRvW563bbNEB1OOyCt+TXA1HZ9xGoB4LebV4UQBmEpbozOKNdPstgzHrZCWAmYQlOmOzt3os7HW1zswSwEzCEp0xLlfEHdYBzpmwRGeMy/VYWAo1swQwm7BEd+YHZj1L9ZCVAGYTluiMcbkeCxvjmlkCmE1YojPFzeGq42o4gNmEJTpjXK6Hq+EAzp2wRGeKnqV6LF4NJywBzCIs0RnDcj3O7uDdaxkAq4KwRGcMzPXRswQwm7BE5yzD9c8+SwDnTliiM253UhE9SwDnTFiiO4tbB/RaBRm/Gq7nQgBWAWGJzhiX62MZDmA2YYnOGJfrcXYH737rAFgNhCU6s9CzZBWuHnqWAGYTluiMcbke9lkCOHfCEp0xMNdHzxLAbMISnSkuh6tG07gaDuBcCUt0ZmESo5GVqmFmCWA2YQkGTIM3wGzCEp0xMNfD1gEA505YojPG5frIrwCzCUt0ZrFnqd8yiBvpApwPYQkGphnrsLc0CjCbsERnDMv10bMEMJuwRGeKdbjqWIYDmE1YojOG5TqMZ1UzSwCzCUt0xsRSJcZeAD1LALMJS3TIwFyD8bAqKwHMJizRGQNzffQsAcwmLNGZs8Oyhbg+jW8doGcJYDZhic7oWaqPmSWA2YQlOlMiLdVgsmdJWAKYRViiM8bl+liGA5hNWKIzluHqMNayZBkO4BwIS3Sm2DqgOrISwGzCEt0xMFehiRvpApwPYQmGZmIZrr8yAFYLYYnOGJfrMHlvOK8KwCzCEp1ZWPJpdHhXw8wSwGzCEp0xLtfBPksA50dYojPG5UrYOgDgvAhLdMa94eowcTVcj3UArBbCEp1Z7FnquQ7OWucdAGAmb5V0xixGHezgDXB+hCW6szAum1qqhgZvgNmEJTrjdid1mNxnqbcyAFYNYYnOuJFufSzDAcwmLNEZw3IdGrc7ATgvwhKdMYlRi4ltKXurAmC1EJZgYPQsAZwfYYnOLDR461nqma0DAM6LsERnjMt1mJhZ8g4AMJO3SjpT7LNUHfssAcwmLNEhA3NtLMMBzCYs0Zmz+yyZWupTM7Z3gAZvgNmEJTpjXK7D5MYBXhWAWYQlOmPFpz5mlgBmE5bojFmMOozv4K3BG2A2YYnOLPYsaVmqhpklgNk29F0Aw7EwLs81TZ4/eqLXWobs+KlvL/59nbQEMJOwRGcWlnzmmuTgN1/ruRoSM0sA58IyHJ254S+9PZs2OORq8e4dW7P78u19lwFQvdLU0UBSRRFcfK+dPJ3XTn579jdy0b390s3ZsF54BQbrnOfWhSUAYIjOOSz5WAkA0EJYAgBoISwBALQQlgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0MKNdGGA/u+XX85DT34pr79xpu9SoBpbN67PR3/i+/IjV1/WdylURliCAfrl3/1S/uDLL/ddBlRnrmmEJZYQlmCAjr8xupnxxvUll27Z2HM10L/XTp7O6TPN4rkB44QlGKCF+2d/95aNufbdb+m1FqjBF75+NC8deyN13Fue2mjwhgGamx8RyjnfcxvWtoVTYU5YYhnCEgyQAQGmjeJSY2qJZQhLMEALA4KJJZg3fzLMCUssQ1iCATo7IIhLkFiGo52wBAO0MCDoWYJJZpZYjrAEA6QvAyYtfnBwarAMYQkGaCErmViCkbPLcNISSwlLMEBz0hJMGZ0MepZYjrAEAzQnK8EEM0u0EZZggFwNB1PmTwVZieUISzBAVuFgeWaWWI6wBAOkZwkmWYajjbAEAyQrwRTLcLQQlmCAfHqGScXVcLQQlmCA7OANy7NhK8sRlmCAGlfDwQQ9S7QRlmCADAgwZT4tWYZjOcISDJBNKWGSmSXaCEswQLYOgOXJSixHWIIhkpVgQlncOkBaYilhCQZoYWZJWIIFtg7gzQlLMECLA4K9AyCJniXaCUswQGaWYIodvGkhLMEAGRBgkpkl2ghLMECuhoPlCUssR1iCAbIMB5POziz1WgaVEpZggIwHMGXsYgfbBzBNWIKBaZomVuFg0vi5ICsxTViCgZkYCGwdAEvoW2KasAQDMz4QiEowMv65Qd8S04QlGBgDAbQzs8Q0YQkGZmJmydQSJNGzRDthCQbGQADLORuXzCwxTViCgWmiZwmmjc+yikpME5ZgYCZ7lsQlmGZmiWnCEgyMniVYaqJnaa63MqiUsAQDYyCAZUxsHWBmiUnCEgyMfZagnbDENGEJBmZiIJCWIMnkqWAvMqYJSzAwc7ISLMONdHlzwhIMTDNxYbS4BImtA2gnLMHATKzCyUqQZHoZTlxikrAEA2MggHZ6lpgmLMHA6FmCpcZnWeekJaYISzAwBgJoZ/KVacISDIyeJViOG+ny5oQlGBgDASylwZs2whIMzOTGAaaWIMlEWhKVmCYswcDYwRuWmriRrpklplQRlg4ePJibbrop27Zty44dO3LbbbflxRdf7LssWJMa94aDVq6BYNqGvgsopXzvZZddlu3bt+e+++7LsWPH8sADD+SZZ57Jvn37smnTpr5LhDXFQABLTWwdYGaJKb2HpSQ/f/z48ezfvz9XXHFFkmTPnj258cYb88gjj+SOO+7ouTxYWwwEsJyxq+HmeiyDKl3wMlwp5dJSyi+VUr5SSjlVSnmhlPK5UsoPlVI+VUo5XUp52zI/9+lSytFSypb5hz50yy23LAalJNm7d2927dqVRx999ELLA97E+EBg6wAYcTUcbVbSs/TrST6c5LEkH0nyQJITSd6b5D9mNGv1t8d/oJSyKcnfSvJY0zQnSynvSvL26667bsmT79mzJ08//fQKygOWMzfRsyQtwTRZiWkrWYa7Oclnmqb52Nhjv7jwl1LKHyT5e0l+depnviejMJUkO5Nk586dS558586dOXLkSE6dOpXNmzevoEwAaFcmtg6Qlpi0kpmlo0neV0q5/E2+/tn5r1899tjfTXIoye/P//fWJMuGoS1bRqt0J06cWEGJwLSJmSUTS7CEiyCYtpKwdE+SH0hyqJSyr5TyyVLKVWNf/29JTmUUkFJK2Z7kliT/uTl77fKJJDl16tSSJz958mSSZOvWrSsoEZhmIICl9CzR5oLDUtM0jya5KslHkzyf5OeSHCil/NT8119J8njmw1JGvUqbk/ynsac5nCSHDx9e8vyHDx/Ojh07LMHBd9icfZZgqbFpVptSMm1Fm1I2TXO4aZqHm6b5mSRXJnk5ycfHvuWzSXaVUq7PKDQ93TTNgbGf/0aSF5966qklz71v375ce+21KykPWIaBAJaanFnqrQwqdUFhqZSyfn5ZbVHTNC9kNMM0PhX0O0leSvLPk7w/k7NKCx57/PHHc+jQocUHnnzyyTz77LO59dZbL6Q8oMXEQKBpCZaYk5aYcqFXw12a5OullN9K8oUkx5LsTXJ9ksWr45qmOV1K+a9J/kmSM0n+yzLPdd8ll1xy5w033JC77rorx44dy/33359rrrkmt99++wWWB7yZ8YFAVIIRM0u0udCw9HqSh5N8IMkHM5qhei7JR5qm+bWp7/1sRmHpyaZpljQnNU1z6MCBA7n77rtz7733ZtOmTbn55pvz4IMP6leCi8A4AMuwdQAtLigsNU3zRkZXw91zDt/+xvyfyy3BJUl2796dJ5544kJKAc6TrQNgqfFTQVsf01bU4H2O/nFGy3T/vYPfBcxgIIDljN0bzknClIt2I91Syt9I8v1J7kjyq03THL9Yv6vN5597KQ89+aU+fjVU6dUTpxf/bmIJRsZnWf/tbx/M9q3P9VcMSZKf3P3O/MMfu7LvMpJcxLCU5FeSvCPJ/0ryry7i72n10rFT+cM/O9LXr4eq7f3+d+Tyt9j4Fb756sk89dVXkiR/8s3Xeq6GJPmL77y07xIWXbSw1DTNey7Wc5+Py7Ztzvuu3NF3GVCdPVfuyJ3vvzpF4xKkaZocf+NM/vDLL/ddCvPe89bv6ruERaWSDeqqKAIAGIxz/qTYRYM3AMCqJSwBALQQlgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0EJYAgBoISwBALQQlgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0EJYAgBoISwBALQQlgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0EJYAgBoISwBALQQlgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0EJYAgBoISwBALQQlgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0EJYAgBoISwBALQQlgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0EJYAgBoISwBALQQlgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0EJYAgBoISwBALQQlgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0EJYAgBoISwBALQQlgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0GJD3wXMK30XAACwHDNLAAAthCUAgBbCEgBAC2EJAKCFsAQA0EJYAgBoISwBALQQloBVqYx8opTyV/uuBVjbhCVgtdqV5FNJdvZdCLC2CUvAavXD83/+Ua9VAGteaZqm7xoAzkspZV+S66cefrVpmrf0UA6wxtVybziA8/ELST6ZZHOSfz3/2NG+igHWNjNLwKpUSvlqkt9rmub2vmsB1jY9S8CqU0rZnuSKJH/cdy3A2icsAavRD87/KSwBF52wBKxGC2HpC71WAQyCsASsRj+Y5HDTNC/1XQiw9glLwGp0RZKv910EMAy2DgBWoz9L8tdLKfckeT7JwaZp9vdcE7BG2ToAWHVKKZcn+c0kP5pkW5KfbZrmV/qtClirhCUAgBZ6lgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0EJYAgBoISwBALQQlgAAWghLAAAthCUAgBbCEgBAC2EJAKCFsAQA0EJYAgBo8f8B5HopvMhoymwAAAAASUVORK5CYII=\n", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "class MyModel(Model):\n", " \"\"\"A custom Hamiltonian model with sigmax and sigmay control.\"\"\"\n", " def get_control(self, label):\n", " \"\"\"\n", " Get an available control Hamiltonian.\n", " For instance, sigmax control on the zeroth qubits is labeled \"sx0\".\n", "\n", " Args:\n", " label (str): The label of the Hamiltonian\n", "\n", " Returns:\n", " The Hamiltonian and target qubits as a tuple (qutip.Qobj, list).\n", " \"\"\"\n", " targets = int(label[2:])\n", " if label[:2] == \"sx\":\n", " return 2 * np.pi * sigmax() / 2, [targets]\n", " elif label[:2] == \"sy\":\n", " return 2 * np.pi * sigmay() / 2, [targets]\n", " else:\n", " raise NotImplementError(\"Unknown control.\")\n", "\n", "class MyCompiler(GateCompiler):\n", " \"\"\"Custom compiler for generating pulses from gates using the base class \n", " GateCompiler.\n", "\n", " Args:\n", " num_qubits (int): The number of qubits in the processor\n", " params (dict): A dictionary of parameters for gate pulses such as\n", " the pulse amplitude.\n", " \"\"\"\n", "\n", " def __init__(self, num_qubits, params):\n", " super().__init__(num_qubits, params=params)\n", " self.params = params\n", " self.gate_compiler = {\n", " \"ROT\": self.rotation_with_phase_compiler,\n", " \"RX\": self.single_qubit_gate_compiler,\n", " \"RY\": self.single_qubit_gate_compiler,\n", " }\n", "\n", " def generate_pulse(self, gate, tlist, coeff, phase=0.0):\n", " \"\"\"Generates the pulses.\n", "\n", " Args:\n", " gate (qutip_qip.circuit.Gate): A qutip Gate object.\n", " tlist (array): A list of times for the evolution.\n", " coeff (array): An array of coefficients for the gate pulses\n", " phase (float): The value of the phase for the gate.\n", "\n", " Returns:\n", " Instruction (qutip_qip.compiler.instruction.Instruction): An instruction\n", " to implement a gate containing the control pulses. \n", " \"\"\"\n", " pulse_info = [\n", " # (control label, coeff)\n", " (\"sx\" + str(gate.targets[0]), np.cos(phase) * coeff),\n", " (\"sy\" + str(gate.targets[0]), np.sin(phase) * coeff),\n", " ]\n", " return [Instruction(gate, tlist=tlist, pulse_info=pulse_info)]\n", "\n", " def single_qubit_gate_compiler(self, gate, args):\n", " \"\"\"Compiles single-qubit gates to pulses.\n", " \n", " Args:\n", " gate (qutip_qip.circuit.Gate): A qutip Gate object.\n", " \n", " Returns:\n", " Instruction (qutip_qip.compiler.instruction.Instruction): An instruction\n", " to implement a gate containing the control pulses.\n", " \"\"\"\n", " # gate.arg_value is the rotation angle\n", " tlist = np.abs(gate.arg_value) / self.params[\"pulse_amplitude\"]\n", " coeff = self.params[\"pulse_amplitude\"] * np.sign(gate.arg_value)\n", " if gate.name == \"RX\":\n", " return self.generate_pulse(gate, tlist, coeff, phase=0.0)\n", " elif gate.name == \"RY\":\n", " return self.generate_pulse(gate, tlist, coeff, phase=np.pi / 2)\n", "\n", " def rotation_with_phase_compiler(self, gate, args):\n", " \"\"\"Compiles gates with a phase term.\n", "\n", " Args:\n", " gate (qutip_qip.circuit.Gate): A qutip Gate object.\n", " \n", " Returns:\n", " Instruction (qutip_qip.compiler.instruction.Instruction): An instruction\n", " to implement a gate containing the control pulses.\n", " \"\"\"\n", " # gate.arg_value is the pulse phase\n", " tlist = self.params[\"duration\"]\n", " coeff = self.params[\"pulse_amplitude\"]\n", " return self.generate_pulse(gate, tlist, coeff, phase=gate.arg_value)\n", "\n", "\n", "# Define a circuit and run the simulation\n", "num_qubits = 1\n", "\n", "circuit = QubitCircuit(1)\n", "circuit.add_gate(\"RX\", targets=0, arg_value=np.pi / 2)\n", "circuit.add_gate(\"Z\", targets=0)\n", "\n", "myprocessor = ModelProcessor(model=MyModel(num_qubits))\n", "myprocessor.native_gates = [\"RX\", \"RY\"]\n", "\n", "mycompiler = MyCompiler(num_qubits, {\"pulse_amplitude\": 0.02})\n", "\n", "myprocessor.load_circuit(circuit, compiler=mycompiler)\n", "result = myprocessor.run_state(basis(2, 0))\n", "\n", "fig, ax = myprocessor.plot_pulses(\n", " figsize=(5, 3), dpi=120,\n", " use_control_latex=False\n", ")\n", "ax[-1].set_xlabel(\"$t$\")\n", "fig.tight_layout()" ] }, { "cell_type": "markdown", "id": "d856fe27", "metadata": {}, "source": [ "We now define a custom `ClassicalCrossTalk` noise object that uses the `Noise` class as the base.\n", "The `get_noisy_dynamics` method will be called during the simulation to generate the noisy Hamiltonian model.\n", "Here, we define a noise model that adds the same driving Hamiltonian to its neighbouring qubits, with a strength proportional to the control pulses strength applied on it.\n", "The detuning of the qubit transition frequency is simulated by adding a $\\sigma_z$ drift Hamiltonian to the processor, with a frequency of $1.852$ MHz.\n", "\n", "\n" ] }, { "cell_type": "code", "execution_count": 3, "id": "6ca2b2b6", "metadata": {}, "outputs": [], "source": [ "class ClassicalCrossTalk(Noise):\n", " def __init__(self, ratio):\n", " self.ratio = ratio\n", "\n", " def get_noisy_dynamics(self, dims=None, pulses=None, systematic_noise=None):\n", " \"\"\"Adds noise to the control pulses.\n", " \n", " Args:\n", " dims: Dimension of the system, e.g., [2,2,2,...] for qubits.\n", " pulses: A list of Pulse objects, representing the compiled pulses.\n", " systematic_noise: A Pulse object with no ideal control, used to represent\n", " pulse-independent noise such as decoherence (not used in this example).\n", " Returns:\n", " pulses: The list of modified pulses according to the noise model.\n", " systematic_noise: A Pulse object (not used in this example). \n", " \"\"\"\n", " for i, pulse in enumerate(pulses):\n", " if \"sx\" not in pulse.label and \"sy\" not in pulse.label:\n", " continue # filter out other pulses, e.g. drift\n", " target = pulse.targets[0]\n", " if target != 0: # add pulse to the left neighbour\n", " pulses[i].add_control_noise(\n", " self.ratio * pulse.qobj,\n", " targets=[target - 1],\n", " coeff=pulse.coeff,\n", " tlist=pulse.tlist,\n", " )\n", " if target != len(dims) - 1: # add pulse to the right neighbour\n", " pulses[i].add_control_noise(\n", " self.ratio * pulse.qobj,\n", " targets=[target + 1],\n", " coeff=pulse.coeff,\n", " tlist=pulse.tlist,\n", " )\n", " return pulses, systematic_noise" ] }, { "cell_type": "markdown", "id": "51e1c43b", "metadata": {}, "source": [ "Lastly, we define a random circuit consisting of a sequence of $\\pi$ rotation pulses with random phases.\n", "The driving pulse is a $\\pi$ pulse with a duration of $25 \\, \\mu\\rm{s}$ and Rabi frequency $20$ KHz.\n", "This randomized benchmarking protocol allows one to study the classical cross-talk induced decoherence on the neighbouring qubits.\n", "The two qubits are initialized in the $|00\\rangle$ state with a fidelity of 0.975.\n", "After the circuit, we measure the population of the second qubit.\n", "If there is no cross-talk, it will remain perfectly in the ground state.\n", "However, cross-talk induces a diffusive behaviour of the second qubit and the fidelity decreases.\n", "\n", "This simulation is repeated 1600 times to obtain the average fidelity. This may take several hours, therefore, the following code only takes two samples with $t=250$. The full simulation is in the commented lines." ] }, { "cell_type": "code", "execution_count": 4, "id": "4ef0af10", "metadata": {}, "outputs": [], "source": [ "\n", "\n", "\n", "\n", "def single_crosstalk_simulation(num_gates):\n", " \"\"\" A single simulation, with num_gates representing the number of rotations.\n", "\n", " Args:\n", " num_gates (int): The number of random gates to add in the simulation.\n", "\n", " Returns:\n", " result (qutip.solver.Result): A qutip Result object obtained from any of the\n", " solver methods such as mesolve.\n", " \"\"\"\n", " num_qubits = 2 # Qubit-0 is the target qubit. Qubit-1 suffers from crosstalk.\n", " myprocessor = ModelProcessor(model=MyModel(num_qubits))\n", " # Add qubit frequency detuning 1.852MHz for the second qubit.\n", " myprocessor.add_drift(2 * np.pi * (sigmaz() + 1) / 2 * 1.852, targets=1)\n", " myprocessor.native_gates = None # Remove the native gates\n", " mycompiler = MyCompiler(num_qubits, {\"pulse_amplitude\": 0.02, \"duration\": 25})\n", " myprocessor.add_noise(ClassicalCrossTalk(1.0))\n", " # Define a randome circuit.\n", " gates_set = [\n", " Gate(\"ROT\", 0, arg_value=0),\n", " Gate(\"ROT\", 0, arg_value=np.pi / 2),\n", " Gate(\"ROT\", 0, arg_value=np.pi),\n", " Gate(\"ROT\", 0, arg_value=np.pi / 2 * 3),\n", " ]\n", " circuit = QubitCircuit(num_qubits)\n", " for ind in np.random.randint(0, 4, num_gates):\n", " circuit.add_gate(gates_set[ind])\n", " # Simulate the circuit.\n", " myprocessor.load_circuit(circuit, compiler=mycompiler)\n", " init_state = tensor(\n", " [Qobj([[init_fid, 0], [0, 0.025]]), Qobj([[init_fid, 0], [0, 0.025]])]\n", " )\n", " options = Options(nsteps=10000) # increase the maximal allowed steps\n", " e_ops = [tensor([qeye(2), fock_dm(2)])] # observable\n", "\n", " # compute results of the run using a solver of choice with custom options\n", " result = myprocessor.run_state(init_state, solver=\"mesolve\",\n", " options=options, e_ops=e_ops)\n", " result = result.expect[0][-1] # measured expectation value at the end\n", " return result\n", "\n", "\n", "# The full simulation may take several hours\n", "# so we just choose num_sample=2 and num_gates=250 as a test\n", "num_sample = 2\n", "fidelity = []\n", "fidelity_error = []\n", "init_fid = 0.975\n", "num_gates_list = [250]\n", "\n", "# The full simulation is defined in the commented lines below.\n", "\n", "# from joblib import Parallel, delayed # for parallel simulations\n", "# num_sample = 1600\n", "# num_gates_list = [250, 500, 750, 1000, 1250, 1500]\n", "\n", "for num_gates in num_gates_list:\n", " # expect = Parallel(n_jobs=8)(delayed(single_crosstalk_simulation)(num_gates) for i in range(num_sample))\n", " expect = [single_crosstalk_simulation(num_gates) for i in range(num_sample)]\n", " fidelity.append(np.mean(expect))\n", " fidelity_error.append(np.std(expect)/np.sqrt(num_sample))\n" ] }, { "cell_type": "markdown", "id": "190e9545", "metadata": {}, "source": [ "We plot a recorded result as an illustration." ] }, { "cell_type": "code", "execution_count": 5, "id": "dc34fdf3", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAdMAAAEnCAYAAADy/amiAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAA9hAAAPYQGoP6dpAABA2ElEQVR4nO3dd3hc5Zn38e8t2ZZ77703jHu3bOOGTQmQQEI2gQ3JbtqGJLx5Nxt409g0SNmEBEhPCBCS7GYhdHC35SJL7hV3XHCvsixbliXd7x/nyIxlSZY0Gs1I+n2uay5rzjznzD3yaH5zznnO85i7IyIiIpWXFO8CREREajqFqYiISJQUpiIiIlFSmIqIiERJYSoiIhIlhamIiEiUFKYiIiJRUpiKiIhESWEqIiISJYWpiIhIlOIepmb2BTPbZ2a5ZpZhZmPLaFvfzL5lZnvC9hvNbE6xNslm9l0ze9fMLoZtv2lmFvtXIyIidVFcw9TM7gV+CvwnMBLYCMw1s/alrPI94LPAF4HBwK+Bf5jZiIg2XwM+DzwIDArv/0e4joiISJWzeA50b2YZwGp3fzC8nwQcBJ5098dLaH8Y+L67Px2x7EXgorvfF95/HTjm7v9SWpty1GVAZyC70i9ORERqi2bAYS8jMOtVYzFXMbMGwCjgsaJl7l5oZguACaWslgLkFlt2EUiNuL8S+IyZ9Xf3nWY2LHz8K2XUkhJuu0gnYHt5X4uIiNR6XYFDpT0YtzAF2gLJwLFiy48BA0tZZy7wFTNLA/YAM4APhdsp8jjQHNhuZgXhY1939xfKqOUR4NvFFx48eJDmzZuX46WIiEhtdO7cObp16wbXOVIZzzCtjC8DvyPYa3SCQH0G+FREm48AHwc+BmwFhgNPmNlhd3+2lO0+RnDutkgz4L3mzZsrTEVE5LriGaYngQKgQ7HlHYCjJa3g7ieAu8ysIdAGOEywJ7o3otmPgcfd/W/h/c1m1oNg77PEMHX3S8Clovvq+CsiIhURt9687p4HrCU4VAtc6YA0A0i/zrq57n6I4MvA3cArEQ83BgqLrVJAAlwGJCIitVO8D/P+FHjWzNYAmcBDQBOCQ7eY2XPAIXd/JLw/DugCbAj/fZQgJH8Usc3XgK+b2QGCw7wjCDof/THmr0ZEROqkuIapu/+3mbUDvgN0JAjJOe5e1CmpO1fvZTYkuNa0N3AeeBO4393PRrT5IvBd4JdAe4JDwb8Jn0NERKTKxfU600RlZs2BrKysLHVAEhGpw86dO0eLFi0AWrj7udLa6TyiiIhIlBSmZbicnx/vEkREpAZQmJbhid//N1t3vhvvMkREJMEpTMtw6sw5fv/XV/ndX1/l5JmseJcjIiIJSmFahqnjhpOUlMS2ne/yw6ef563F6eRd1qFfERG5mnrzliCyN+/FvHxeemsJO/ceBKB1y+bcNXsKQwb01khJIiK1XHl78ypMS1D80hh3Z9M7u3l5bhpnz50HYGDfHnxozlTatWkV32JFRCRmFKZRKO0600t5l1mwbDWL09dRUFBAcnIy0yaMYObksaQ0qB+/gkVEJCYUplG43qANJ06d4aW3l7J9934AWjZvyp2zpzBsUF8d+hURqUUUplEozwhI7s6WHXt5eW4ap88Gv9/+vbvxoTk30aFd62qsVkREYkVhGoWKDCeYdzmfhctXs2jFWvILCkhKSmLq+BHcPGUsDVMaVE/BIiISEwrTKFRmbN6Tp8/y8ty0K4M8tGjWhDtmTWbEkP469CsiUkMpTKMQzUD3W3e+yz/eXsqpcJCHvj278qFbptKpfdsYVCoiIrGkMI1CtLPGXM7PZ/GKtSxYvprL+QUkmTF53HBmTx1Ho4YpVV+wiIjEhMI0ClU1Bdvps+d4eW4am7fvAaBZk8Z8YFYqo4cO1KFfEZEaQGEahaqez3T77v289NYSTpw+C0Cv7p25+5ab6NKxXdTbFhGR2FGYRiEWk4Pn5+ezZNV65qdlknc5HzNj0uih3DJtPI0bNayS5xARkaqlMI1CLMK0yJmsbF6dv4wNW3cB0KRxI26fMZGxwweTlKR5B0REEonCNAqxDNMiO/ce4KW3lnLs5GkAunZqz9233ETPbp1i8nwiIlJxCtMoVEeYAhQUFLB89SbeXrKK3Et5AIweOpDbZ6bSolmTmD2viIiUj8I0CtUVpkWycy7wxsKVZK7figMpDepz85RxTBk/nHrJyTF/fhERKZnCNArVHaZFDhw6yktvLWX/oaMAtGvTkg/Omcqgvj2rrQYREXmfwjQK8QpTgEJ31mx8h9cXrCA75wIAN/TvxZ2zp9CudctqrUVEpK5TmEYhnmFa5GLuJealZZKWsYHCwsKIuVPHkNJAA+iLiFQHhWkUEiFMixw7eZp/vL2UHXsOANCiWVPumJWqAfRFRKqBwjQKiRSmEMydunXnu7w8N+3KAPoaRUlEJPYUplFItDAtcjk/nyUr17Fg+eoroyhNHDWEW6ZNoEnjRvEuT0Sk1ilvmMZ9yB0z+4KZ7TOzXDPLMLOxZbStb2bfMrM9YfuNZjanhHZdzOzPZnbKzC6a2WYzGx3bVxJ79evVY9aUsTz8hX9mxA39cXdWrNnMD558luWrN1JYWBjvEkVE6qS47pma2b3Ac8DngAzgIeDDwAB3P15C+x8C9wGfBrYDs4GfAhPdfX3YphWwHlgM/Ao4AfQD9rj7nnLWlZB7psXt3vceL721hCPHTwHQuUNbPnTLTfTp0SXOlYmI1A414jCvmWUAq939wfB+EnAQeNLdHy+h/WHg++7+dMSyF4GL7n5feP9xYJK7T46irhoRpgAFhYWkr9nMW4vTuZB7CYARQ/rzgZmptGrRLM7ViYjUbAl/mNfMGgCjgAVFy9y9MLw/oZTVUoDcYssuAqkR9+8A1pjZ383suJmtN7NPX6eWFDNrXnQDakwKJSclkTp2GI988RNMHHUjBqzfspPHn36O+WmZXM7Pj3eJIiK1XjzPmbYFkoFjxZYfAzqWss5c4Ctm1s/MksxsFvAhIHJ0+N7A54FdBIeBfwX8wsw+UUYtjwBZEbf3Kvha4q5p40Z8+PbpfOUz/0Sv7p3Ju5zPm4vT+eEv/8zm7XtQRzMRkdiJ22FeM+sMHCI435kesfxHwFR3H1fCOu2A3wEfABzYQ7An+yl3bxS2yQPWuPvEiPV+AYxx9xL3eM0shWCvt0gz4L2acJi3JO7Oui07eG3+crKycwDo37sbd82eSqf2beJcnYhIzZHwh3mBk0AB0KHY8g7A0ZJWcPcT7n4X0AToAQwEzgN7I5odAbYVW/UdoHtphbj7JXc/V3QDsivwOhKOmTHqxoE88uA/MzN1NMnJyezce5Cf/PoFXnxzCTkXix8pFxGRaMQtTN09D1gLzChaFnZAmgGkl7ZeuG6uux8C6gF3A69EPLwCGFBslf7A/ioou0ZJadCA22ZM4uEv3M/QQX0odGf56o384Mlnmbc0gx17D3D2XI3+3iAikhDi3Zv3XuBZ4LNAJsGlMR8BBrr7MTN7Djjk7o+E7ccBXYAN4b+PAr2Ake5+NmwzBlgJfBv4H2AswaHhz7j7C+Wsq8b05q2IXe8e5H/fXMzxk2euWj5rylgmjR6qOVRFRIqpCYd5cff/Bv4d+A5BQA4H5rh7Uaek7lzduagh8D2Cw7j/IDjnmloUpOE2VwMfBP4J2AJ8E3iovEFam/Xr1Y0BfXpcs3x+WiZvLFoRh4pERGqHevEuwN2fAp4q5bGbit1fCgwuxzZfB16vivpqm47tSu6AtHbjdpo1bsysKWNomJJSYhsRESlZ3IcTlOp1/OSpEpcXurNo5Vp+8ORzZKzfSqEupRERKTcNdF+C2nrOFCArO4cVazaxYNlq3B0zY2bqaNq2acWCtExOnD4LQNdO7fngnKn07t45vgWLiMRRjRhOMFHV5jAtcvZcNidPZ9G2dQtaNg8GfMovKGBZxkbmpWWQeykPgJFD+nO7hiYUkTpKYRqFuhCmZck+n8Obi9LJWL8VJ5itZkbqaKZNHEmD+vXjXZ6ISLVRmEahrodpkfeOHOcfby9l74HDALRs3pQ7Zk1m+A39MLM4VyciEnsK0ygoTN/n7mzYtovX5i/nTFYwwEOv7p354JypdOvUPs7ViYjElsI0CgrTa+VdvsziletYuHwNl/PzMWDciBu4dfoEmjXVYA8iUjspTKOgMC3dmaxsXl+4gnWbdwDBkIWzp45l8rjh1EtOjnN1IiJVS2EaBYXp9e09cJh/vL2U944cB6Bt6xbcMWsyQwb01vlUEak1FKZRUJiWT6E7qzds442FK8nOuQBAv15dufPmKXTp2C7O1YmIRE9hGgWFacXkXspj4fLVLElfT35Bgc6nikitoTCNgsK0ck6fPcdrC5azYesuIDifOnPyaKaOH0H9enEfBlpEpMIUplFQmEZn74FDvDw3jYOHg/OprVs25wMzUxk2uK/Op4pIjaIwjYLCNHqF7qzdtJ03Fq4kK/s8EFyfetfsKXTv3CHO1YmIlI/CNAoK06pzKe8yi1euZdGKtVzOzwdg9LBB3DZ9Ii2bN41zdSIiZYtZmJrZUuAPwN/d/WJUVSYohWnVO3sum9cXrmTtpu0ANKhfj+mTNN6viCS2WIbpE8DHgBTgf4A/uPuqypeaeBSmsbP/0FFenpvGvoNHgGC839tmTGLkjQNI0vlUEUkwMT3Ma2b1gDuATwC3ALuBPwLPu/uxSlWcQBSmseXubNi6i9cWvD/eb/cuHbhr9hR6ddP8qSKSOKrtnKmZtQc+A3wdSAbeBH7h7oui2nAcKUyrR97lfNJWrWfB8tVcyrsMwIgh/bl9xiRat9TvXUTir1rC1MzGAp8EPgqcA/4EdCE4DPxLd//3Sm88jhSm1etcOH9q5pX5U5OZOn4kM1JH0zClQbzLE5E6LJbnTNsD9xOEaD/gNeD3wFwPN2ZmqcDb7l4ju2sqTOPjvSPHeWXeMnbvew+AZk0bc9v0iYwZNoikpKQ4VycidVEswzQP2ENwjvRP7n6ihDbNgVfcfVqFNp4gFKbx4+5s2bGXV+ct4+SZLAC6dGzHXbOn0Ldn1zhXJyJ1TSzDdLK7L4uyvoSmMI2//Px8lmVuYl5aBrmX8gC4cWAfPjArlXatW8a3OBGpM2IZpouAD7n72WLLmwMvu/v0ipebWBSmieN8zgXeXrKKlWu34O4kJyUxacxQbp46jiaNGsa7PBGp5WIZpgVAJ3c/Xmx5e+CQu9f4K/AVponnyPFTvDp/Gdt37wegUcMUbp4yltQxQ6mnQfRFJEaqPEzNbGj44wZgOnA64uFkYA7wWXfvWYl6E4rCNHFt37OfV+ct48jxUwC0adWC22dOYtggDaIvIlUvFmFaCBQ1LulT6yLwRXf/YwVrxcy+AHwV6AhsDLeTWUrb+sAjBANGdAF2AF9z97dLaf8w8Bjwc3d/qJz1KEwTWGFhIZkbtvHm4nSyzweTkvfq1ok7b55Cj64d41ydiNQmsQjTHgQhuhcYC0T24s0Djrt7QUULNbN7geeAzwEZwEPAh4EBxQ8lh+1/CNwHfBrYDswGfgpMdPf1xdqOIRjy8BywWGFau1zKy2PRirUsXrnuyiD6GvRBRKpSjZk1xswygNXu/mB4Pwk4CDzp7o+X0P4w8H13fzpi2YvARXe/L2JZU2Ad8G/AN4ANCtPa6ey587y1eCWrN7yDA/WSk5kybjgzJ4+hUcOUeJcnIjVYecO0XD03zOwO4C13vxz+XCp3f7W8RZpZA2AUwWHYovULzWwBMKGU1VKA3GLLLgKpxZY9Dbzh7gvM7BvXqSMl3G6RZuUoXxJEy+ZN+ac7b2bKuBG8Mi+NXe++x6KVa0lft5XRQwcwdfwI2rRqEe8yRaQWK9eeaXi+tKO7Hw9/Lo27e3K5n9ysM3CI4BBtesTyHwFT3X1cCev8BRgG3EUweMQM4BUg2d1TwjYfJRgreIy755rZEsrYMzWzR4FvF1+uPdOax91Zs2k7r8xNI+fi+9+5hg/ux52zp2gOVRGpkPLumZZrjDZ3Tyo6fxn+XNqt3EEahS8DuwjOl+YBTwHPAIUAZtYN+DnwcXcvvgdbmseAFhE3DbVTQ5kZh44evypIATZs28UTv/8bB49ccxpeRCRq8R7w9CRQAHQotrwDcLSkFdz9hLvfBTQBegADgfMEHaMgOGzcHlhnZvlmlg9MBb4U3r8m8N39krufK7oB2dG/NImXnl1LnsYtKzuHn/32r7zw8jzOntN/sYhUnfKeM/1SeTfo7r+oQNs8M1tLcKj25fC5ksL7T11n3VzgUHipzN0EvXYBFgI3Fmv+DMGe7A8r0+NYapbWrVpgBpFnMMyMAb27s33PftZsfIeNW3dx04SRTJ80SjPTiEjUyjt0zP8pZzsHyh2moZ8Cz5rZGiCT4NKYJgQBiJk9RzCy0iPh/XEE15duCP99lGAP+0cA7p4NbIl8AjPLAU65+1XLpXbatnMvxbsCuDs9unbklmnjeWXeMvYeOMz8ZZmsWreFOdPGM27EDSRrZhoRqaS4XxoDYGYP8v6gDRuAL7l7RvjYEmCfuz8Q3p8K/AroTXB4903gYXc/XMb2l6BLY+qMrOwczp3PIfv8Bc6ey6Zl82Y0a9qY5k2b0KJZE9ydzdv38NqC5Zw8HcxM06l9Gz4wK5VBfXvGt3gRSSgxv840vKylF7DH3fMrtZEEpTCtG/ILClixehPz0jK5EHZYGtCnO3fMmkznDm3jXJ2IJIJYDnTfGHiSYDg/gP7uvtfMniQ4HHvNQAs1jcK0brlwMZf5y1azLGMDBYWFmBnjhg9mzrQJtGjWJN7liUgcxTJMfw5MIji3+TYwNAzTO4FH3X1EpatOEArTuunkmSxeX7CCjdt2AdCgfn2mTRzJtIkjSWmgTkoidVEsw3Q/cK+7rzKzbGBYGKZ9gXXuXuPTR2Fat7178DCvzFvG/veCq7OaNWnM7JvGMX7kEHVSEqljYhmmF4AhYYBGhukwIM3da/y4bQpTcXc2vbOb1xeuuNJJqX2bVtw+cxJDBvTWdG8idUQswzQN+Lu7PxmG6VB3fzc8Z9rP3edEU3giUJhKkYKCAlau3cLcpRnkXLgIQK/unbljVio9u3aKc3UiEmuxDNNU4C3gz8ADwG+AwcBEgvF011ay5oShMJXici9dYtGKtSxJX39lurdhg/ty24xJtGvdMr7FiUjMxPTSGDPrAzxMMOB80VRnP3T3zZUrN7EoTKU0Z89l8/biVWRu2IYDSUlJTBp9IzdPGUvTJo3jXZ6IVLEaM59pIlKYyvUcPnaS1xcs553d+wFomNKAGZNGM2X8cBrUrx/n6kSkqlRpmIbhUi5lPVlNoTCV8tq59wCvLVjBe+FsNC2aNeWWaeMZM2wQSer5K1LjVXWYFhKMu3td1TQNW0wpTKUiCt1Zv3kHbyxayZmsYDaaouEJB/bpoZ6/IjVYVYfp1Ii7PYHHgT8BRRN6TyAYEekRd3+2ciUnDoWpVMbl/HyWZ25k/rLVXMy9BEC/Xt24Y1YqXTu1j3N1IlIZsezNuxD4vbv/tdjyjwGfcfebKl5uYlGYSjRyLuayIC2TZas3UVAQzPg3auhAbp02gdYt9X4SqUliPWjDMHffVWx5f4KZWWp8l0aFqVSFU2eyeHNxOus27wCgXnIyk8cNZ2bqaBo3ahjn6kSkPGIZpjuAV9z9P4ot/xFwp7sPqES9CUVhKlXpwOFjvDZ/Obv3vQdA44YpzJoyltQxQ6lXr7xTCotIPMQyTG8FXgR2Axnh4rFAP+Bud3+zUhUnEIWpVDV3553d+3ht/gqOnjgFQOuWzblt+kSGD+lPkjopiSSkWA/a0A34PDAwXPQO8Gt3P1iJWhOOwlRipaCwkNUbtvHW4lWcO58DQNdO7bljVir9enWLc3UiUpwGbYiCwlRi7VLeZZauWs+iFWu4lHcZgIF9e3D7jEl06dguztWJSJGqvjRmKLDF3QvDn0vl7psqWmyiUZhKdcnOucC8pRmsXLuFwsJCDBipnr8iCSMWgzZ0dPfjEQM4lHSSxzVog0jFnTh9lrcWpbN+604AkpOTmTT6RmZNGUvTxo3iXJ1I3VXVYdoDOODuHv5cKnffX9FiE43CVOLl4OFjvLZgBbveDbofpDRowIxJo5gyfgQpDTTmr0h1q+owXQfMcPczZvYt4CfufqHKqk0wClOJtx179vPaghUcOnoCgGZNGzNn6njGjRhMcnKNP/gjUmNUdZheJJj4+z0zKwA6ufvxKqs2wShMJREUurNhy07eWLSS02eDv+F2rVty64yJDBvUV2P+ilSDqg7TdOA8sBz4NvCT8P413P07lSk4kShMJZHkFxSQvnYL85ZmcP7CRQC6dGzHxNE30q1zh6vaNm/ahBbNmsSjTJFaqarDdADwn0AfYCSwDcgvoam7+8hKVZxAFKaSiHIv5bEkfR2LV64j7/LlEtvMnjqOOTeNr+bKRGqvWI6AdKVnb3QlJi6FqSSy7PM5vLZgJWs2brtqXsTUscOYmTpGe6YiVai8YVrh2YvdPak2B6lIomvWtAm3Th9/zcVpyzM38tailVcOBYtI9dEo2yI10IlTZynpoFLGhm1sfGc30yfqchqR6lThPdNYMLMvmNk+M8s1swwzG1tG2/pm9i0z2xO232hmc4q1ecTMVptZtpkdN7OXw/O+IrVCSkoKxTvzGtCuTUtyL+Xx5uJ0vv/kn1i5ZvOVOVVFJHbiHqZmdi/wU4IOTiOBjcBcM2tfyirfAz4LfBEYDPwa+IeZjYhoMxV4GhgPzALqA/PMTCeTpFbYtnPvNXumDoy4oT/3f2gOrVs2J/v8Bf7+xiJ++Ks/s3HbLjQOt0jsxH2gezPLAFa7+4Ph/STgIPCkuz9eQvvDwPfd/emIZS8CF939vlKeox1wHJjq7mnlqEkdkCShZWXncO58DtnnL3D2XDYtmzejWdPGVy6NKelymu5dOvCBman07dk1ztWL1Bzl7YBUqXOmZtYSuIfgUpkfu/tpMxsJHHP3QxXYTgNgFPBY0bJwMP0FwIRSVksBcostuwiklvFULcJ/T5dSR0q43SLNytiWSNy1aFb29aT1kpOZPHYYY4YNunI5zYFDx3j62RcZ1LcHt2l2GpEqVZlLY4YCC4AsoCcwwN33mtn3gO7u/s8V2FZn4BAw0d3TI5b/iGAvclwJ6/wFGAbcBewBZgCvAMnunlJC+yTgVaClu5cYuGb2KMFgFFfRnqnUFtnnc5iXlnlldhqAEUP6c8u0CbRr3TK+xYkksJhdGkNwfvNP7t6Pq/cQ3wSmVGJ7FfVlYBewHcgDngKeAQpLaf80MAT4aBnbfIxg77XopuNgUqs0a9qEu2+dxsNfuJ+RQ/oDsH7LTh5/+nn+5/WFZGWXOKCZiJRTZfZMs4CR7r7HzLKBYeGeaQ9gh7s3rMC2GgAXgHvc/eWI5c8S7EneWca6DYE2wGHgceB2d7+hWJungDuBKe7+bgXq0jlTqdUOHT3BG4tW8s6ufQDUr5fM5LHDmZ46miaNyv0nLFLrxXLP9BJQUsL0B05UZEPungesJThUC1w5LDsDSC9tvXDd3PD8bD3gboJDvUXbsDBIPwhMr0iQitQFXTq24zMfu5MHH7iHXt06cTm/gEUr1/K9nz/D/GWZXMorebhCESlZZfZMf0+wR/gRgg49Q4EC4GUgzd0fquD27gWeJbjcJRN4KNz2QHc/ZmbPAYfc/ZGw/TigC7Ah/PdRoBfB3vLZsM0vgY8R7JXuiHi6LHe/7vAw2jOVusTd2bZrH28uWsnhYycBaNakMbOmjGXCqCHU05RvUofFcmzeFsD/AqMJer0eBjoS7Ene6u45FS3WzB4EvhpuZwPwJXfPCB9bAuxz9wfC+1OBXwG9CWaueRN42N0PR2yvtBf1SXf/UznqUZhKnVPozvotO3hr8SpOnckCoHXL5sy5aTyjbhxAUlLcL0sXqXYxC9MrK5qlEuyVNgXWufuCSm0oASlMpS7LLyggY91W5qZlkH3+AgCd2rfh1ukTuaF/L82jKnVKzMO0NlOYikDe5cssy9jIwhVruJh7CYCeXTtx24yJGvhB6oxYHub9UikPOcGlMrsJzp3W2AFBFaYi77twMZfFK9eydNUGLucH0xgP7NODW2dMpFun0kb9FKkdYhmm7wLtgMbAmXBxK4JLXM4D7YG9wDR3P1jx0uNPYSpyrazsHOanZZK+7v2BH4bf0I9bpk2gfZtWca5OJDZiGab/BHwG+Fd33xMu6wv8BvgtsAL4G3DU3e+pXPnxpTAVKd3J02d5e8kq1m3egQNJZowdcQOzp46lZXONxCm1SyzDdA9wt7tvKLZ8BPCiu/c2s4nhz50qXHkCUJiKXN/hYyd4Y1E623YGl3EXjQc8I3U0TRo3inN1IlUjlmF6gWBEoTXFlo8Blrp7YzPrCWxx96YVrjwBKExFym/vgcO8sXAFew8EV6c1TGnAtAkjmTphBCkNGsS5OpHoxDJM3yC4HvRf3X19uGwE8DuCQ7u3m9kHgB+4+42VfQHxpDAVqRh3Z/vu/byxaCWHjgYDoTVt0ohZk8cycdQQ6tWr1ARVInEXyzDtCDxPMORf0Zhj9YCFwP3hqEXTgPruPq8yxcebwlSkcgrd2bB1J28tTufk6WDgh1YtmnHzlLGMGT6YZA38IDVMdQzaMJBgPF4IBrjfUVb7mkRhKhKdgoICMtZvY15aBlnZwaBo7Vq3ZM5N4xk+pD9JGvhBaggN2hAFhalI1ci7nM/KNZtYsHwNOReCYbE7tW/DLdMmMGRAb42mJAkvpmFqZl2BO4DuwFU9DNz9KxXeYIJRmIpUrdxLeaRlbGDxyrXkXsoDoHuXDtw6bQL9e3dXqErCiuU50xnAqwQDMwwEtgA9ASMYo3d6JWtOGApTkdjIuZjLkpVrScvYQN7lYDSlPj26cOv0ifTu3jnO1YlcK5Zhmgm85e7fLpocHDgOvAC87e6/qnzZiUFhKhJb2edzWLB8DSvWbKagIBh5dFDfHtwyXUMUSmKJZZhmA8PdfY+ZnQFS3X2rmQ0DXnH3nlHUnRAUpiLV40xWNvPSMslcv5XC8LNo6KC+3DJtPI0aNuTc+Ryyz1/gzLlsWjVvRrOmjWnetAktmjWJc+VSV8QyTI8SjLv7jpltI5hL9NUwTFfU1IEaIilMRarXidNnmRsxRKGZ0al9myuTlUeaPXUcc24aX/1FSp1U3jCtzEVfq4DU8Oc3gf8ys68DfwwfExGpkHatW3Lfh+bw1c9/nBsH9sHdSwxSM2Nw/95xqFCkbJUZluQrBBOCA3w7/PleYFf4mIhIpXRq35ZP3Xs7Bw4f43/fWMzBw8euetzdycvLi1N1IqWr0J6pmSUDXYEDAO6e4+6fc/eh7n63u++PRZEiUrd079yBT917GyVdMLPpnd1XJisXSRQVCtNwwu95BPOXiojEjFkS0yaNvmb5ssyNfPfnzzB/WSaXtJcqCaIy50y3ADppISIxlb52M4tWrLlmedPGjbiYe4k3F6XzvV/8iaWr1nM5Pz8OFYq8rzK9eecAjwHfBNYCOZGPl9XbqaZQb16R+MvKzuHc+Zxrljdt3Ii9+w/x9pJVnDwTDKbfsnlTZk0Zy7jhg0lOTq7uUqUWi+WlMYURdyNXNsDdvca/kxWmIomvoKCAzA3vMC8tg7PnzgPQplULbp4yllFDB2qGGqkSsQzTqWU97u5LK7TBBKQwFak5Lufns3LNZhYsX835nGAw/XZtWjJn6niG39CPJIWqREGzxkRBYSpS81zKu8zy1RtZvGItORdzAejYrg2zbxrH0EF9Ne2bVEqsZ42ZDHyWoCPSh939kJndD7zr7ssrWXPCUJiK1Fy5l/JYlrGBxenrrlxC07lDW+bcNF7TvkmFxfIw793A8wQD298PDHb3vWb2IHCru99a+bITg8JUpOa7mHuJJenrWLpqw5VLaLp1bs8tN01gYN8eClUpl1iG6XrgZ+7+XNGsMWGYjiCYTaZjNIUnAoWpSO3x/rRvG8m7fBmAnl07MWfaePr36qZQlTLFcmzeAUBaCcuzgJaV2B5m9gUz22dmuWaWYWZjy2hb38y+ZWZ7wvYbw8t1Kr1NEam9mjRqyG0zJvHNLz/AtAkjqV+vHvveO8Kvn/8HTz/7Irv3vRfvEqUWqEyYHgX6lrA8lWDC8Aoxs3uBnwL/CYwENgJzzay0SQ2/R3C+9ovAYODXwD/CPePKblNEarmmTRpzx82T+fqXHmDy2OHUS05mz/5DPP3si/zyuZfYd/BIvEuUGqwyh3kfAe4DPgXMB24FegA/A77r7k9WcHsZwGp3fzC8nwQcBJ5098dLaH8Y+L67Px2x7EXgorvfV5ltlvAcOswrUsudPZfNgmWrWbVuKwWFweXzA/v24JabxtO9S40/WyVVpLyHeSsza8zjBHu0C4HGBId8LwE/qUSQNgBGEYyoBIC7F5rZAmBCKaulALnFll0knBauMts0s5Rwu0WaVeR1iEjN07J5M+65bTrTJ41mflommRu2sX33frbv3s+QAb2Zc9N4unRsF+8ypYao9HWmYWj1JZiCbZu7n6/ENjoDh4CJ7p4esfxHwFR3H1fCOn8BhgF3AXuAGcArQLK7p1Rym48STCd3Fe2ZitQdJ0+fZV5aJms2bafoc3HooL7MuWk8ndq3iXN1Ei8x64BkZveZWWN3z3P3be6eWZkgjcKXCeZO3Q7kAU8BzwCFZa10HY8BLSJuXaOsUURqmLatW/Kxu27m4X+7n5FD+mME0739+Fd/5rkX3+L4yTPxLlESWGU6IP0MOG5mfzGzW8M5TivrJFAAdCi2vANBR6druPsJd78LaEJwrnYgcJ73Oz9VZpuX3P1c0Q3IrvhLEZHaoH3bVtx/9y189fMfZ9jgvjiwfstOHv/l87zw8jxOnD4b7xIlAVUmTDsBHyUY5P5/gCNm9rSZTazohtw9j2DmmRlFy8LOQjOA9NLWC9fNdfdDBOd97yY41BvVNkVEinRq35YHPnwb//7ZjzFkQG/cnTUb3+Hxp57jb68u4PTZGj9BllShqMbmNbPGwAeBjwEzgffcvU8Ft3Ev8CzB5S6ZwEPAR4CB7n7MzJ4DDrn7I2H7cUAXYEP476NAL2Cku58tzzbLUZN684rIVQ4cPsbbi9N5Z/d+AJKTkhg7fDAzJ4+hdUt9TtRWsezNe4W7XzCzuUArgkOugyqxjf82s3bAd4COBCE5JyL0unP1+dCGBNea9iY4vPsmcH9RkJZzmyIiFdK9cwc+8/G72HfwCG8tSWfn3oOkr9tC5oZtjB0xmJmpCtW6rLID3RftkX6c4PDpQeCvwAvuvr1KK4wD7ZmKyPXsPXCIuUsz2Ln3IBDuqY4YzKzJY2jVonyfG6VNgN68aRNaNGtSpfVK5cRybN6/AbcDFwjOmb4QeQlKbaAwFZHy2rM/CNVd774fquNG3MDMyaOvG6pvL1nF3KUZ1yyfPXUcc24aH5N6pWJiGaYvEMwYM9fdC4o9NsTdt1Si3oSiMBWRigpCdRW73g3G+k1OSmL8yBuYkTqGVi1KHgcmKzuHU2eyePKZvwPwyY/cRquWzbVnmkCqbXJwM2sG/BPwr8Aod4/mUpmEoDAVkcras/8Qby9ZdWUA/euF6vLVG3nxzSUAmBkfuX0640cOqc6SpQwxD1MzmwL8C8FlKYeBl4AX3X11pTaYQBSmIhKt3fveY+7SjPdDNTmZ8eHh35bNg1A9ey6b7zzxDJGfw2bGtx765JU2El8x6c1rZh2BBwhCtDnBOdMU4C5331bpakVEapm+PbvSt2dXdu97j7eXrGLP/kOsWLOJVeu3MmHkDcxIHc2JU2cpvkPj7pw8naUwrWHKvWdqZq8BU4A3CM6Zvu3uBWZ2mWCC8FoTptozFZGqFhmqEOypjhzSjzWbdmjPNIFV+WFeM8sHfgH8yt13RSxXmIqIlNOudw/y9pJV7D1wGAjCs+hz2MyYOXkMk0YPVQekBBGLge5TCaYmW2tmGWb2oJm1jbJOEZE6pV+vbjz4wD382z9/iFYtml21V+ruzE/LJH3t5jhWKJVR7jB191Xu/mmCsXl/QzA+7+FwG7PCXr0iInIdZka/Xt340qc+wj23Tbtq3tTk5CTOZGWTlV2dk3FJtKIdm3cAQWek+4GWwHx3v6NqSosfHeYVkerk7lcO/7578AgA9ZKTmTDqRmakjqJFs6ZxrrDuqrbrTAHCadg+AHxKYSoiUjklhWr9esmMH6lQjZdqDdPaRmEqIvHk7uwMQ3VfxJ7q+JE3MH3S6FJHVJKqpzCNgsJURBKBu7Nz7wHmLs24sqdaNKC+ZqmpHgrTKChMRSSRuPuVEZWKrlNNSkpizLBBzJw8hratWsS5wtpLYRoFhamIJKrd+95jXlrmlVlqkswYPWwgM1PH0K5NqzhXV/soTKOgMBWRRLf3wGHmp2Wyfc9+ILjcZuSNA5g1eQwd2raOc3W1h8I0CgpTEakp9r93lHlpGWzbtQ8AA4YP6c+syWPp1L5NXGurDRSmUVCYikhNc/DwMealZbJlx94ry4YN7svNU8bSuUO7MtaUsihMo6AwFZGa6tDRE8xLy2TTO7uvLLtxYB9mTRlLt07t41hZzaQwjYLCVERquiPHTzI/bTUbtu6k6FN+cP9e3DxlLD26dIxrbTWJwjQKClMRqS2OnTjN/GWZrNuy88qg+gP79mD2lHH07NYpztUlPoVpFBSmIlLbHD91hgXLVrN203YKw8/9/r27M3vqWHp37xLn6hKXwjQKClMRqa1Onj7LguVrWL3xHQoLCwHo27Mrs6eOo0+PLphZnCtMLArTKChMRaS2O332HAuWryZz/TYKwlDt3b0zN08dR/9e3RSqIYVpFBSmIlJXnMnKZtGKNaSv20pBQQEAPbp2ZNbksQzu17POh6rCNAoKUxGpa86eO8/ilWtJX7uZy/lBqHbp2I6Zk8cwdFBfkupoqNaYMDWzLwBfBToCG4EvuntmGe0fAj4PdAdOAv8LPOLuueHjycCjwH3hNg8DfwK+5+V8sQpTEamrss/nsCR9PSvWbOJS3mUAOrRtzczJYxgxpD/JSUlxrrB61YgwNbN7geeAzwEZwEPAh4EB7n68hPYfA/4IfApYCfQnCMq/uftXwjb/D/gK8AlgKzAaeAb4urv/opx1KUxFpE7LuXCRtIwNLMvcyMXcSwC0adWCmamjGT1sEPWSk+NcYfWoKWGaAax29wfD+0nAQeBJd3+8hPZPAYPcfUbEsv8Cxrl7anj/deCYu/9LRJsXgYvufl8561KYiogAF3MvsWL1JpasWk/OhYsAtGzelOmTRjNuxA00qF8vzhXGVnnDNG7762bWABgFLCha5u6F4f0Jpay2EhhlZmPDbfQGbgXeLNZmhpn1D9sMA1KBt8qoJcXMmhfdAE1jLyICNGqYwszJY/jmlz/JXbOn0LxpE86eO89Lby3hez9/hkUr1pJ7KS/eZcZd3PZMzawzcAiY6O7pEct/BEx193GlrPcl4CcEkyPUA37t7p+PeDwJ+AHwH0ABkExwiPexMmp5FPh28eXaMxURudrl/HwyN2xj4fI1nMnKBqBxo4ZMHTecyeOG06hhSpwrrFoJf5i3MmFqZjcBfwO+QXCOtS/wc+B37v7dsM1HgR8TdGraCgwHngC+4u7PllJLChD5DmgGvKcwFREpWUFBAWs372DBstWcOH0WgJQGDZg8dihTx4+gaZPG8S2witSEMG0AXADucfeXI5Y/C7R09ztLWGcZsMrdvxqx7D7gt0BTdy80s4PA4+7+dESbbwD3ufvActamc6YiIuVQWFjIhm27WLBsNUeOnwKgQf16TBh1I9MmjqRFs6ZxrjA65Q3TuJ05dvc8M1sLzABehiuHaGcAT5WyWmOgsNiygvBfu06butWfW0SkGiQlJTFyyACG39CfrTv2Mn9ZJgcPH2fpqvUsX72J8SMGM33SaFq3rN07JvHuhvVT4FkzWwNkElwa04TgUhbM7DngkLs/ErZ/DfiKma3n/cO83wVec/eCiDZfN7MDBId5RxBcKvPHanlFIiJ1UJIZNw7sw5ABvdmx5wDzlmXy7oHDrFizmfS1WxjUrydjhw+mVRiqzZs2oUWzJnGuuuokwqAND/L+oA0bgC+5e0b42BJgn7s/EN6vB3wduB/oApwgDE93Pxu2aUYQsB8E2hMM2vBX4DvuXq4uZzrMKyISvT37D/HXV+Zx6sy1R0dnTx3HnJvGx6Gqikn4c6aJTGEqIlI1srJz2LpzL39/fdFVywf26cEt0yfQvXOHOFVWPgl/zlRERGq/Fs2aXJnqLdL2PfvZvmc/A/v0YObkMfTu3rlGD6qvMBURkZg5ey6bl95aetUyM+PGAb3ZsmPvlVDt2a0TM1PH1NiZahSmIiISMydOnaX46UR3Z/K44Xzg5sksXrGWzA3b2HfwCL//66t07tCWGamjGT64H0k1aFD9mlOpiIjUOO3atLxmT9PMaNu6BW1bteDDt0/nG1/+JNMmjiKlQX0OHzvJ8y++zWNPPUf62i3k5+fHqfKKUQekEqgDkohI1cjKziEtYwOLVqwBgiCdOXkMk0YPvebSmJyLuSzP3EhaxgYuXMwFgnOuN00YyYRRQ0hp0KDa61dv3igoTEVEqsbbS1Yxd2nGNcvLujTmUl4e6Wu3sCR9PVnZ54Fg/N8p44aTOnYYTRo1jGnNkRSmUVCYiohUjazsHM6dz7lmeXkGbcjPz2fNpu0sXLGGk6ezAEhpUJ+Jo25k6oSR1TLog8I0CgpTEZHEUVhYyMZtu1mwfDWHj50EIDk5mbHDBzF90mjatmoRs+dWmEZBYSoiknjcnXd272P+stXsO3gECM7BjhjSnxmTRtO5Q9sqf06FaRQUpiIiicvd2XvgMAuWrWb7nv1Xlg8Z0JsZqaPp2bVTlT2XwjQKClMRkZrh4JHjLFy+mk3bdlOUZn17dmXm5DH079Ut6gEgFKZRUJiKiNQsx0+eYeGKNazZtP3K8IXdOrdnZuoYhgzsQ1IlQ1VhGgWFqYhIzXQm6xyLV65n1botXA4HfOjQtjUzUkczckh/kpOTK7Q9hWkUFKYiIjVbds4F0lZtYPnqjeReCmbfbNWiGdMmjmLciME0qF+/XNtRmEZBYSoiUjtczL3EijWbWLpqPedzLgLQpHGjYACIMUNpfJ0BIBSmUVCYiojULnmX81m9YRuLVq7l9NkgE1Ma1GfCqBuZOn4ELZs3LXE9hWkUFKYiIrVTQWEhG7fuYuGKNe8PAJGUxOhhg5g+cRTt27a6qr3CNAoKUxGR2s3d2b57PwtXrGHP/kMAGHDjoL7MmDSK7l06AgrTqChMRUTqjncPHmbRirVs2bH3yrJ+vboxI3U0Hdu0oGXLlqAwrTiFqYhI3XPk+CkWrVzLus07rlyr2rZlU77x0L/CdcJUk4OLiIgAndq34eN33czXv/gJJo8dTv169Th09GS51tWeaQm0ZyoiIucvXGTuohXc84FZoD1TERGRimvauBGzpowtV1uFqYiISJQUpiIiIlFSmIqIiERJYSoiIhIlhamIiEiUFKYiIiJRUpiKiIhEqV68C0hk586Ven2uiIjUAeXNAY2AVAIz6wm8G+86REQkYXR190OlPag905KdDv/tCmTHs5AKaga8R82quybWDKq7OtXEmqFm1l0Ta4bY190MOFxWA4Vp2bLLGosx0ZhZ0Y81pu6aWDOo7upUE2uGmll3TawZqqXu625THZBERESipDAVERGJksK0ZJeA/wz/rUlqYt01sWZQ3dWpJtYMNbPumlgzJEDd6s0rIiISJe2ZioiIRElhKiIiEiWFqYiISJQUpiIiIlFSmJbAzL5gZvvMLNfMMsxsbBxrecTMVptZtpkdN7OXzWxAsTYNzexpMztlZufN7EUz61CsTXcze8PMLoTb+bGZVcugHWb2sJm5mT2R6DWbWRcz+3NY10Uz22xmoyMeNzP7jpkdCR9fYGb9im2jtZm9YGbnzOysmf3BzJrGqN5kM/uumb0b1rPHzL5pEVexJ0LNZjbFzF4zs8Phe+GuYo9XSY1mNtTMloV/uwfN7D9iVbeZ1TezH4bvkZywzXNm1jmedV/vd12s7a/DNg/Fs+by1m1mg8zsVTPLCn/nq82se8Tj8ftccXfdIm7AvQTdqz8JDAZ+C5wB2sepnreBB4AbgGHAG8B+oElEm18BB4DpwCggHVgR8XgysBmYDwwHbgFOAD+ohvrHEIxzvBF4IpFrBloB+4BngLFAL+BmoE9Em68BZ4E7gaHAK8BeoGFEm7eADcA4IBXYBfwlRjX/P+AkcBvQE7iHYDi1LyVSzeH/3/eADwIO3FXs8ahrBJoDR4E/h38vHwUuAJ+JRd1Ai/D9+RFgADAeyADWFNtGtdZ9vd91RLsPhnUdAh6KZ83lfI/0AU4BPwJGhPfvIOKzmTh+rlT5H3dNv4V/DE9F3E8K32wPx7u2sJ524RttSni/BZAH3BPRZmDYZnx4/xagAOgQ0eZzQBbQIIa1NgV2AjOBJYRhmqg1A48Dy8p43IAjwL9HLGsB5AIfDe8PCl/H6Ig2c4BCoHMMan4d+EOxZS8Cf07gmouHUpXUCHyeYFztBhFtHge2x6LuUtqMCdt1T4S6S6sZ6EIwlu0NBF8gH4p4LCF/18DfgOfLWCeunys6zBvBzBoQfJtZULTM3QvD+xPiVVcxLcJ/iwbjHwXU5+qatxN8OyuqeQKw2d2PRWxnLsG3yxtiWOvTwBvuvqDY8kSt+Q5gjZn9PTz8s97MPh3xeC+gY7G6swi+gEXWfdbd10Sst4Dgg2hcDGpeCcwws/4AZjaMYE/irQSuubiqqnECkObueRFt5gIDzKxVjGovrgXBh/fZiJoSqm4zSwKeB37s7ltLaJKoNd8G7DSzueHfZ0axQ8Fx/VxRmF6tLcFhgGPFlh8j+GOPq/AN9QTBYYst4eKOQJ67ny3WPLLmjpT8miBGr8vMPgqMBB4p4eGErBnoTfCNexcwm+CQ0S/M7BPFnres90dH4Hjkg+6eT/DlJxZ1P07wjX27mV0G1hMcAXghgWsurqpqjMd75gozawj8EPirvz/YeiLW/TUgH/hFKY8nYs3tCY50PUxw6utm4B/AS2Y2NeJ54/a5olljapangSEEex4Jy8y6AT8HZrl7brzrqYAkgvNd/y+8v97MhhAcBno2fmWV6SPAx4GPAVsJzgM9YWaH3T1Ra651zKw+8D8Eh6w/H+dySmVmo4AvAyM9PMZZQxTt+L3i7j8Lf95gZhMJ/j6Xxqes92nP9GonCY+nF1vegeBke9yY2VPA7cA0d38v4qGjQAMza1lslciaj1Lya4LYvK5RBN8k15lZvpnlA1OBL4U/H0vAmiE4b7et2LJ3gKLegkXPW9b74yjBa78i7CnYmtjU/WPgcXf/m7tvdvfngZ/x/hGBRKy5uKqqMR7vmcgg7UHwBTJyuq5Eq3tyWM+BiL/NHsB/mdm+BK0Zgs/mfK7/9xm3zxWFaYTw+P9aYEbRsvDQ6gyCXmHVzgJPEfRwm+7u7xZrsha4zNU1DyB4gxXVnA7caGaRfyCzCOboK/7mrAoLgRsJ9pKKbmuAFyJ+TrSaAVYQ9MqM1J+g9zQEvZKPcnXdzQnOI0XW3TLcAygyneBvLSMGNTcmOJcVqYD3/7YTsebiqqrGdGBKGG5FZgE73P1MLAqPCNJ+wEx3P1WsSaLV/TxBb+nhEbfDBF/KZidozUWfzasp++8zvp+FVdHzqjbdCC6NyQU+QdCr7TcEl8Z0iFM9vyTozDCV4Jh+0a1RRJtfhW+oaQR7hSuBlRGPF3UHn0twec1sgnMiMb80JqKGJVx7aUxC1UzQE/MyweUmfQkOneYAH49o87Xw/XAHwReGlyn5Eo51BJfXTCLo0RyrS2P+RNArs+jSmA8SdPX/YSLVTHC+a3h4c+D/hD93r6oaCTr/HAWeI+hMcm/4/xfN5Rql1k3Q2eUV4GD4Ho38+2wQr7qv97suof0+Sr40JmF+1+HjHyTorftpgr/PBwn2VlMjthG3z5Uq/+OuDbfwP2k/wfWmGcC4ONbipdweiGjTkOB86unwDf0S0LHYdnoAbxJcC3YC+AlQrxpfxxKuDtOErJngUPpmgi9U7wCfLva4Ad8JP0hyCXoO9i/WpjXwF4LrPbOAPwJNY1RvM4JOafuBi8Aegmv1GiRSzcBNpbyP/1SVNRLsdS0Lt/Ee8LVY1U3w5aW0v8+b4lX39X7XJbTfx7VhmlC/64g2nyLoIHiR4DrYO4ttI26fK5qCTUREJEo6ZyoiIhIlhamIiEiUFKYiIiJRUpiKiIhESWEqIiISJYWpiIhIlBSmIiIiUVKYioiIRElhKlLNzGyJmT0R7zqKhOM//9bMTpuZm9nweNdUlRLt9y21k8JU6hwz+1MYGg8XW36XmdXFIcHmAA8QDKXYCdhSZusYq2z4lbHeh4BvRluXSFkUplJX5QJfM7NW8S6kqphZg0qu2gc44u4r3f2oBxNBV7ko6ouKu5929+x4PLfUHQpTqasWEAyq/khpDcxsn5k9VGzZBjN7NOL+EjN70syeMLMzZnbMzD5tZk3M7Bkzyzaz3WZ2S7HN1zOzp8wsy8xOmtl3zcwitptkZo+Y2btmdtHMNprZPcVqWRJu4wkzO0kwE0ZJryPFzH5hZsfNLNfMlpvZmPCxPwFPAt3DvfV9pWzjnfDxkm4PlrJOifWVo56pwJcjtt8zfGxO2PasmZ0ys9fNrE851rtqj7Ws549o/wsz+1F46Pto5P952OYeM9sc/t+cMrMFZtakpN+D1A0KU6mrCgimWvuimXWNclufIJi8eCxBMP0K+DvB9E8jgXnA82bWuNg6+eE6Xwa+AvxrxOOPAP8MfI5giqufAX82s6klPHcewTRZnyulvh8Bd4dtRwK7gblm1jp87m8RzPrRiWAaupLcHf47I2zXk2Ae1Q8DvytlndLqu1496eE2O4W3g+F6TYCfAqPDOgqBf1gw53BZ61Xk9xFZdw7BnKr/AXzLzGYBmFkn4K8EM6kMIpjt5CWCmW+krqqqKZZ0062m3Aimz3o5/Dkd+EP4813Bn8SVdvu4dmqqDcCjEfeXAMsi7icD54HnIpZ1JJhKanzEOtsgmLUpXPY4sC38OYXgg3xCsef+PVfPKbkEWHed19qEIMw+FrGsPnAI+Gp4/yFg33W2M41grteU8P6o8DWVOs9vSfWVs54lREzXV8b224Y1DClrvcjlFXj+ZcW2kQk8Hv48MnzeHvF+L+uWODftmUpd9zXgE2Y2KIptbCr6wd0LgFMEc6IWORb+2z5i2Sp3j+zslA70M7NkgomPGwPzzex80Y1gT7VPsedee53a+hCExYqIGi8ThENFXvONwE53vxTeHwYcd/djZaxTUn2VrsfM+pnZX81sr5mdI/iyA8FE3eVV3uffVGy9I7z//7cRWAhsNrO/h4f1a825d6mcevEuQCSe3D3NzOYCjxHssUYq5NpDd/VL2Mzl4puNXObuHp4OLe+X16bhv7cR7DFFulTsfk45txmtoVz9BWFYsfulqcr6XiOYBP3TwGGC3+cWIBYdm0r6P02C4AtTeMh3InAz8EXg+2Y2zt3fjUEtUgNoz1QEHgY+AEwotvwEwbk3AMysOdCrip5zXLH744Fd4Z7tNoLQ7O7uu4vdSjsPWJo9vH/OEgAzq09wbnRbBbYzlKv31oZx7d5bVdWTR3C4nIg2bYABwPfcfaG7vwMU3xu8Zr1KPv91eWCFu38bGBFu84PlXV9qH+2ZSp3n7pvN7AXgS8UeWgQ8YGavAWeB7xB0XKoK3c3sp8BvCM7BfRH4v2E92Wb2E+BnYeea5UALggA45+7PlvdJ3D3HzH4F/NjMTgMHCDrUNAb+UJ5thDXcQPD6i/Qh6HRTIeWsZx8wLuyNex44DZwhOHz+GTM7QnBo9/Fim79mPXcvrMTzl8nMxhF0gJoHHCf4YtQOeKdcvwSplRSmIoFvAfcWW/YYwZ7o60AWwYX/VbVn+hzQiOBcXQHwc+C3EY9/k2DP+BGgN0GYrwN+UInnepjgKNTzQDNgDTDb3c+Uc/0+BGETuSe6GfhPM1vr7itKXq3S9fwEeJZgT7ER0Mvd95nZR4FfEBza3UHw5WdJxHavWY/3z6tW5Pmv5xwwhaDjVnOCQ8//193fKuf6UgvZ1X0gREREpKJ0zlRERCRKClMREZEoKUxFRESipDAVERGJksJUREQkSgpTERGRKClMRUREoqQwFRERiZLCVEREJEoKUxERkSgpTEVERKKkMBUREYnS/wcfQ4LALqwlyQAAAABJRU5ErkJggg==\n", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "# Recorded result of a full simulation\n", "num_gates_list = [250, 500, 750, 1000, 1250, 1500]\n", "fidelity = [\n", " 0.9566768747558925,\n", " 0.9388905075892828,\n", " 0.9229470389282218,\n", " 0.9075513000339529,\n", " 0.8941659320508855,\n", " 0.8756519016627652\n", "]\n", "\n", "fidelity_error = [\n", " 0.00042992029265330223,\n", " 0.0008339882813741004,\n", " 0.0012606632769758602,\n", " 0.0014643550337816722,\n", " 0.0017695604671714809,\n", " 0.0020964978542167617\n", "]\n", "\n", "\n", "def rb_curve(x, a):\n", " return (1 / 2 + np.exp(-2 * a * x) / 2) * 0.975\n", "\n", "\n", "pos, cov = curve_fit(rb_curve, num_gates_list, fidelity, p0=[0.001])\n", "\n", "xline = np.linspace(0, 1700, 200)\n", "yline = rb_curve(xline, *pos)\n", "\n", "fig, ax = plt.subplots(figsize=(5, 3), dpi=100)\n", "ax.errorbar(\n", " num_gates_list, fidelity, yerr=fidelity_error, fmt=\".\", capsize=2, color=\"slategrey\"\n", ")\n", "ax.plot(xline, yline, color=\"slategrey\")\n", "ax.set_ylabel(\"Average fidelity\")\n", "ax.set_xlabel(r\"Number of $\\pi$ rotations\")\n", "ax.set_xlim((0, 1700));" ] }, { "cell_type": "code", "execution_count": 6, "id": "5de0f329", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "qutip-qip version: 0.2.0\n" ] }, { "data": { "text/html": [ "
SoftwareVersion
QuTiP4.6.3
Numpy1.22.2
SciPy1.8.0
matplotlib3.5.1
Cython0.29.27
Number of CPUs12
BLAS InfoOPENBLAS
IPython8.0.1
Python3.9.0 | packaged by conda-forge | (default, Nov 26 2020, 07:53:15) [MSC v.1916 64 bit (AMD64)]
OSnt [win32]
Sun Feb 13 11:53:32 2022 W. Europe Standard Time
" ], "text/plain": [ "" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import qutip_qip\n", "print(\"qutip-qip version:\", qutip_qip.version.version)\n", "from qutip.ipynbtools import version_table\n", "version_table()" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.9.0" } }, "nbformat": 4, "nbformat_minor": 5 }