{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Simulating Clifford randomized benchmarking using a generic noise model\n",
    "\n",
    "This tutorial demonstrates shows how to simulate Clifford RB sequences using arbitrary $n$-qubit process matrices.  In this example $n=2$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pygsti\n",
    "import numpy as np"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Get some CRB circuits\n",
    "\n",
    "First, we follow the [Clifford RB](../CliffordRB.ipynb) tutorial to generate a set of sequences.  If you want to perform Direct RB instead, just replace this cell with the contents of the [Direct RB](../DirectRB.ipynb) tutorial up until the point where it creates `circuitlist`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Specify the device to be benchmarked - in this case 2 qubits\n",
    "nQubits = 2\n",
    "qubit_labels = [0,1] \n",
    "gate_names = ['Gxpi2', 'Gypi2','Gcphase'] \n",
    "availability = {'Gcphase':[(0,1)]}\n",
    "pspec = pygsti.obj.ProcessorSpec(nQubits, gate_names, availability=availability, \n",
    "                                 qubit_labels=qubit_labels)\n",
    "\n",
    "#Specify RB parameters (k = number of repetitions at each length)\n",
    "lengths = [0,1,2,4,8,16]\n",
    "k = 10\n",
    "subsetQs = [0,1]\n",
    "randomizeout = False # ==> all circuits have the *same* ideal outcome (the all-zeros bitstring)\n",
    "\n",
    "#Generate clifford RB circuits\n",
    "exp_design = pygsti.protocols.CliffordRBDesign(pspec, lengths, k, qubit_labels=subsetQs, randomizeout=randomizeout)\n",
    "\n",
    "#Collect all the circuits into one list:\n",
    "circuitlist = exp_design.all_circuits_needing_data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Create a model to simulate these circuits\n",
    "Now we need to create a model that can simulate circuits like this.  Two things to note:\n",
    "\n",
    "1. RB circuits use our \"multi-qubit\" gate naming, so you have gates like `Gxpi2:0` and `Gcphase:0:1`.\n",
    "2. RB circuits do gates in parallel (this only matters for >1 qubits), so you have layers like `[Gypi2:0Gypi2:1]`\n",
    "\n",
    "In this example, we'll make a model with $n$-qubit process matrices, so this will be practically limited to small $n$.  We construct a model based on our standard 2-qubit X, Y, and CPHASE model, since this \n",
    "has all the appropriate gates.  To get a model with the multi-qubit labels, we'll use a standard multi-qubit \"model-pack\", which packages a `Model` object with relevant meta information needed by other protocols (like GST).  If you can't start with a standard model, then you'll need to create an `ExplicitOpModel` object of the appropriate dimension (see the [explicit models tutorial](../../objects/ExplicitModel.ipynb)) and assign to it gates with are, for instance `('Gxpi2',0)` rather than just `'Gxpi2'`.\n",
    "\n",
    "Here we import the `smq2Q_XYCPHASE` model pack:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "from pygsti.modelpacks import smq2Q_XYCPHASE"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We'll depolarize the target model and set one of the process matrices to a custom value as a demonstration.  Here is where you can set any 2-qubit process matrices you want to any of the gates:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "myModel = smq2Q_XYCPHASE.target_model().depolarize(op_noise=0.01, spam_noise=0.01)\n",
    "myModel[('Gx',0)] = np.kron( \n",
    "    np.array([[1, 0, 0, 0],\n",
    "              [0, 0.85, 0, 0],\n",
    "              [0, 0, 0, -0.85],\n",
    "              [0, 0, 0.85, 0]], 'd'),\n",
    "    np.array([[1, 0, 0, 0],\n",
    "              [0, 0.95, 0, 0],\n",
    "              [0, 0, 0.95, 0],\n",
    "              [0, 0, 0, 0.95]], 'd'))\n",
    "#print(myModel[('Gx',0)])\n",
    "myModel.operations.keys() #voila! you have gates like \"Gx:0\" rather than \"Gxi\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since, `ExplicitOpModel` objects (e.g., those in the model packs) don't know how to automatically simulate multiple gates in parallel (you'd need to add an operation for each layer explicitly), we'll just *serialize* the circuits so they don't contain any parallel gates.  This addresses point 2) above.  Then we can simulate our circuits using our `ExplicitOpModel`, creating a `DataSet`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "serial_circuits = [c.serialize() for c in circuitlist]\n",
    "ds = pygsti.construction.generate_fake_data(myModel, serial_circuits, 100, seed=1234)\n",
    "\n",
    "#See how the DataSet contains serialized circuits (just printing the first several layers for clarity)\n",
    "print(ds.keys()[10][0:7])\n",
    "print(circuitlist[10][0:5])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, we \"un-serialize\" the circuits in the resulting data-set (`ds`) using the `process_circuits` function.  This is needed because the RB experiment design calls for the original (parallel-gate) circuits, not the serialized ones.  The cell below updates the circuits for all the data we just simulated so the data counts are associated with the original circuits."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "#map circuits in dataset back to non-serialized RB circuits that we expect to have data for:\n",
    "unserialize_map = { serial_circuit: orig_circuit for (serial_circuit, orig_circuit) in zip(serial_circuits, circuitlist)}\n",
    "ds = ds.copy_nonstatic()\n",
    "ds.process_circuits(lambda c: unserialize_map[c])\n",
    "ds.done_adding_data()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Running RB on the simulated `DataSet`\n",
    "To run an RB analysis, we just package up the experiment design and data set into a `ProtocolData` object and give this to a `RB` protocol's `run` method.  This returns a `RandomizedBenchmarkingResults` object that can be used to plot the RB decay curve.  (See the [RB analysis tutorial](../RBAnalysis.ipynb) for more details.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = pygsti.protocols.ProtocolData(exp_design, ds)\n",
    "results = pygsti.protocols.RB().run(data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "results.plot()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}