{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Getting started\n",
    "\n",
    "This tutorial shows you how to get started with the `povm_toolbox`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%load_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compute circuit output probabilities with `POVMSampler`\n",
    "Follow these instructions to get the probability distribution of a quantum circuit with the `povm_toolbox.sampler.POVMSampler` sampler."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Initialize QuantumCircuit\n",
    "The first step is to create the `qiskit.circuit.QuantumCircuits` from which you want to obtain the probability distribution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAPEAAACuCAYAAADnE+srAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjkuMCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy80BEi2AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAOgklEQVR4nO3df1DTd57H8VcCSPgRKj+0AUF+iCggP6zICVO7gwV7VPGcXt2651jvTkfrnatz65jp7t1ea3dvXGad3T3Xdg/25sbOdkrxdO1h2GvHWa4VPY/GIjdWgqzUWAL5rn4FK4YfNpD7w9GRI0gCyTf5fHk9ZpyOyTf5vJny5PvNN1+ixuVyuUBEwtIGegAimhlGTCQ4RkwkOEZMJDhGTCQ4RkwkOEZMJDhGTCQ4RkwkOEZMJDhGTCQ4RkwkOEZMJDhGTCQ4RkwkOEZMJDhGTCQ4RkwkOEZMJDhGTCQ4RkwkOEZMJDhGTCQ4RkwkOEZMJDhGTCQ4RkwkOEZMJDhGTCQ4RkwkOEZMJDhGTCQ4RkwkOEZMJDhGTCQ4RkwkOEZMJLjQQA9AE7lcLmBkJNBjeCc8HBqNJtBTzEqMOBiNjMD57W2BnsIrocffBXS6QI8xK/FwmkhwjJhIcIyYSHCMmEhwjJhIcIyYSHCMmEhwjJhIcIyYSHCMmEhwjJhIcIyYSHCMmEhwqo9YlmUYjUZkZmZCp9MhJSUF+/btg8PhwPbt26HRaHD06NFAj0l+dv+bUTR/LuHDJis+OmdDt3Qv0CP5jKp/FbGtrQ2VlZWQJAlRUVHIyclBb28vjhw5gq6uLvT19QEACgsLAzuon3wq30TFhU/wk5x8fG/RUrfbzDl9HC/OT8SHf7Ja4emUYb81iHfqLfj1yav44+2hR7drNMC61Sn47l/kYG1pcgAnnDnV7ollWUZVVRUkScL+/ftht9vR2toKSZJQXV2NxsZGmM1maDQa5OfnB3pc8oNLFhnPvPIhflzbNi5gAHC5ANPZbrzw2scw/uyzBx/EICjVRrx3717YbDbs2bMHhw8fhl6vf3Sf0WhEQUEBnE4n0tLSEBMTE8BJyR+ufXUXa1/7GJI8NOW2Pz12GW++c0mBqfxDlRFbLBbU19cjISEBhw4dcrvNihUrAAAFBQXjbr9+/To2bNgAvV6P2NhYvPrqq7h9+7bfZybf+v4/myH3D3u8/Y9qL8HaM+DHifxHlRHX1dVhbGwMW7ZsQXR0tNttIiIiAIyPeGBgAGVlZbDZbKirq0NtbS2am5uxfv16jI2NKTK7PwyOjkIeGXH7R416bzpwqumGV49xuYCaEx1+msi/VHliq6mpCQBQVlY26TY2mw3A+Ihra2vR09ODs2fPYuHChQCA5ORklJaWoqGhARs3bvTf0H701tUreOvqlUCPoZj3f/clRke9f4177D/+gEP7VvphIv9SZcQ3bjz4KZyamur2fqfTifPnzwMYH7HJZMKzzz77KGAAKCkpQUZGBk6fPj3tiIuKiiBJksfbR2i1aC8smdZa7uxYmIE/T0pxe1/l/3zqkzWysrIwFCRHK3ciKwHdKq8fJ8lDWJC8EBoo/3UYDAZcvHhxWo9VZcQOhwMAMDTk/qRGfX09ZFmGXq9Henr6o9vb29uxadOmCdvn5uaivb192vNIkoSenh6Pt48MCQEKp73cBJnR0Xh+3tO+e0I3ent7MTg66tc1PJY4AEzzgzd7e3qAAEQ8E6qM2GAwoL+/H62trSgpGb9Hs9vtOHDgAAAgPz9/3Gcl9/f3Y+7cuROeLy4uDlevXp3RPN6I0Ip3qiIpKSlo9sQDujHcncbjtGN3kbgg0efzeMLb75HHqTLi8vJyWCwWVFdXo6KiAllZWQAAs9mMrVu3QpZlAMpd5OHtYZJreFi4z53u7OyEJkg+d9p+axAL134Ap5evi1/fuRr/tPfv/DSV/4j3I98DRqMR8fHx6O7uRm5uLvLy8rB48WIUFxcjIyMDa9asATDx7aXY2FjcuXNnwvP19fUhLi5OidHJBxLnReKl8jSvHqPVarDz5SX+GcjPVBlxcnIympubsW7dOuh0OlitVsTFxaGmpgaNjY3o7OwEMDHi7Oxst69929vbkZ2drcjs5Bs/2bcS8+M8PzJ4c/dypCbpp94wCKkyYuBBkCaTCQMDAxgYGEBLSwt27twJh8MBq9UKrVaLZcuWjXvM+vXrce7cuUdvPwFAS0sLurq6UFVVpfSXQDOQnqzHmdpKJM2PnHLbH+wowD/sLPT/UH6icYl80eg0tLS0YNWqVViyZAk6Osa/uX/37l3k5eUhISEBBw8exPDwMIxGI+bNm4cLFy5Aq9AJJxFfE4cefzdoXhM/7ubtIdSc6EDNv3eg5+bguPteej4Ne76TjbLipABN5xuq3RNP5vLlywAmHkoDQExMDJqampCYmIjNmzdjx44dKC0thclkUixg8q358RH44a7lsH70Cv77N+sR/1Q4AMAQr8PJnz8vfMCASs9OP8mTIgaARYsWwWQyKTkSKSA0VIuSgqehCw8BAISEqOeHsnq+Eg9NFTGRaGbdnvjhddVEajHr9sREasOIiQTHiIkEx4iJBMeIiQTHiIkEx4iJBMeIiQTHiIkEx4iJBMeIiQQ3666dFkJ4OEKPvxvoKbwTHh7oCWYtRhyENBoNEIS/YE/BiYfTRIJjxESCY8REgmPERIJjxESCY8REgmPERIJjxESCY8REgmPERIJjxESCY8REgmPERIJjxESCY8REgmPERIJjxESCY8REgmPERIJjxESCY8REgmPERIJjxESCY8REguOHx5OqSfIgPm+X8Xn7bXxpG0Df1yMAgDsD9/FvpzqxIiceORmxCAsTd3+mcblcrkAPQeRLwyNOnDhjxTv1Flz435tTbh8bMwd/vTELu1/JxqKUGAUm9C1GTKrhcrnwnuka9h/+DLf6h6f1HJv/NANHXl+FeXERPp7OfxgxqYL91iB2HjwH09nuGT/XvFgd3vn7Ury8Nt0Hk/kfIybhWb68g4qd/4mem4M+fd43XluON3Yvf/AP3AUxRkxC67R+jdV/acLNvukdPk/lH3ctx8G/fcYvz+0rjJiEdW/wGxRuOoWu7gG/rvPeoW9hy7pMv64xE+KeV6dZ7/VfmL0O2Fy3Ad1nNsNct8Hjx3z30AXYb/n2UN2XZkXEsizDaDQiMzMTOp0OKSkp2LdvHxwOB7Zv3w6NRoOjR48GekzywqcX7Xj7A4vXjzMkRCL56SgYEiI9fkz/3ft47UfnvV5LKaq/2KOtrQ2VlZWQJAlRUVHIyclBb28vjhw5gq6uLvT19QEACgsLAzsoeeWtf7mk6HoNn3yFto7bKFwar+i6nlD1nliWZVRVVUGSJOzfvx92ux2tra2QJAnV1dVobGyE2WyGRqNBfn5+oMclD3Vcv4Omz+yKr/ur497v+ZWg6oj37t0Lm82GPXv24PDhw9Dr9Y/uMxqNKCgogNPpRFpaGmJixLtSZ7aqPdERkHXfM3VhwHE/IGs/iWojtlgsqK+vR0JCAg4dOuR2mxUrVgAACgoKHt32MPri4mKEh4cH/XuEs9F/mZXfCwPA4LAT5i/kgKz9JKqNuK6uDmNjY9iyZQuio6PdbhMR8eDSuscjvnbtGk6ePAmDwYCVK1cqMit5bnjEiS+u9Qds/c/bGbFimpqaAABlZWWTbmOz2QCMj/i5556D3W5HQ0MDysvL/Tskee2La/1wOgN3aUOr5XbA1p6Mas9O37hxAwCQmprq9n6n04nz5x+8bfB4xFqt73+uFRUVQZIknz/vbDQclgnot7q9z1y3Ycq3jgwJEY/+231m86TbSfIgVn6nYcLtp06fQfL77tefCYPBgIsXL07rsaqN2OFwAACGhobc3l9fXw9ZlqHX65Ge7t8L3SVJQk9Pj1/XmDX0CYDe/V0P3wP2RGiI1uNtHzcy4gy6/5eqjdhgMKC/vx+tra0oKSkZd5/dbseBAwcAAPn5+X4/eWUwGPz6/LPJcNhTmOyAVpKnvqrKkBCB0BAtnKNjkGT3P+Cf9Fzh4SFIWLDAk1G9MpPvEdVGXF5eDovFgurqalRUVCArKwsAYDabsXXrVsjygxMUSlzkMd3DJJqo4/odZP/ZSbf3uTv8/f+6z2xG8tNRkOQhpFR84PX6W195Eb9+0/27HYGi2hNbRqMR8fHx6O7uRm5uLvLy8rB48WIUFxcjIyMDa9asATD+9TAFv6zUpxAdGRaw9VfkJARs7cmoNuLk5GQ0Nzdj3bp10Ol0sFqtiIuLQ01NDRobG9HZ2QmAEYtGq9Vg+dK4gK0fjBGr9nAaALKzs2EymSbcfu/ePVitVmi1WixbtiwAk9FMVH1rIZpb/6j4uknzI1G4JPiunVZ1xJO5cuUKXC4XsrKyEBk58S2JEydOAADa29vH/T0tLQ1FRUXKDUpu/dXGLPzw7VaM3B9VdN1dLy8Nyk/FnJURX758GcDkh9KbNm1y+/dt27bh2LFjfp2NppYQq8O316bjN6Zriq0ZGqrBjpeyFFvPG4zYDX7YSfA7+DfP4Le/t8Ix5FRkvQPb8pE03/v3lZUQfMcGCpgqYgp+6cl6/PR7xYqslbNoLt7YvVyRtaZjVu6JH15XTWLbtWkpTGe/wu+abR4/5uFFHJ5cGAIAuvAQvPvj5xA+J2RaMyqBH5RHQnMMfoMXdn+M85d8f7Z6TpgWp35RjhdXp/j8uX1pVh5Ok3pERYbho1+9gLWlvr0UMjoyDI1vrw36gAHuiUklRkfH8Mv32/GDX17E0PDM3noqX5WEf33zWaQmTfKbFkGGEZOq/OHG1zD+3IyGT77C2Jh339oZyXp8f3sBtr+UJdQnujBiUqVu6R5qT1zFb39vRcf1rycNOn5uOFY/Y8Cul5dibekCaLXixPsQIybVcwx+g7arfejqvovh+6MIC9UiNmYOli+Nx8LEaKH2uu4wYiLB8ew0keAYMZHgGDGR4BgxkeAYMZHgGDGR4BgxkeAYMZHgGDGR4BgxkeAYMZHgGDGR4BgxkeAYMZHgGDGR4BgxkeAYMZHgGDGR4BgxkeAYMZHgGDGR4BgxkeAYMZHgGDGR4BgxkeAYMZHgGDGR4BgxkeAYMZHgGDGR4BgxkeD+DyBtQk9i1toDAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 287.294x200.667 with 1 Axes>"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from qiskit import QuantumCircuit\n",
    "\n",
    "qc = QuantumCircuit(2)\n",
    "qc.h(0)\n",
    "qc.cx(0, 1)\n",
    "qc.draw(\"mpl\", style=\"iqp\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "> **Note**: \n",
    ">\n",
    "> The `qiskit.circuit.QuantumCircuit` you pass to `povm_toolbox.sampler.POVMSampler` must not include measurements. These will later be specified by the chosen POVM."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Initialize `Sampler` and `POVMSampler`\n",
    "Next, create a `povm_toolbox.sampler.POVMSampler` instance, which will internally use a `qiskit.primitives.Sampler` instance."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "from numpy.random import default_rng\n",
    "from povm_toolbox.sampler import POVMSampler\n",
    "from qiskit.primitives import StatevectorSampler as Sampler\n",
    "\n",
    "rng = default_rng(22)\n",
    "\n",
    "sampler = Sampler(seed=rng)\n",
    "povm_sampler = POVMSampler(sampler=sampler)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Initialize measurement procedure\n",
    "Next, create a `povm_toolbox.library.POVMImplementation` instance."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "from povm_toolbox.library import ClassicalShadows\n",
    "\n",
    "measurement = ClassicalShadows(num_qubits=2, seed=12)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Run and get results\n",
    "Now that you have defined your `povm_sampler`, run it by calling the `povm_toolbox.sampler.POVMSampler.run` method, which returns an instance of `povm_toolbox.sampler.POVMSamplerJob`. You can get the results from the job (as a `qiskit.primitives.PrimitiveResult` object) with the `povm_toolbox.sampler.POVMSamplerJob.result` method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "PrimitiveResult([POVMPubResult(data=DataBin(povm_measurement_creg=BitArray(<shape=(), num_shots=256, num_bits=2>)), metadata=RPMMetadata(povm_implementation=ClassicalShadows(num_qubits=2), composed_circuit=<qiskit.circuit.quantumcircuit.QuantumCircuit object at 0x000001E6C4698F10>, pvm_keys=np.ndarray<256,2>))], metadata={'raw_results': PrimitiveResult([SamplerPubResult(data=DataBin(povm_measurement_creg=BitArray(<shape=(256,), num_shots=1, num_bits=2>), shape=(256,)), metadata={'shots': 1})], metadata={})})\n"
     ]
    }
   ],
   "source": [
    "job = povm_sampler.run([qc], shots=256, povm=measurement)\n",
    "result = job.result()\n",
    "print(result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "While this example only uses one `qiskit.circuit.QuantumCircuit`, you can sample multiple circuits by passing a list of `povm_toolbox.sampler.POVMSamplerPubLike` objects to the `povm_toolbox.sampler.POVMSampler.run` method.\n",
    "\n",
    "For each supplied `POVMSamplerPub`, the resulting `PrimitiveResult` object stores a corresponding `povm_toolbox.sampler.POVMPubResult` object."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Get the outcome counts\n",
    "\n",
    "From these results you can extract the outcome counts with the method `povm_toolbox.sampler.POVMPubResult.get_counts`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[Counter({(0, 0): 20, (4, 5): 18, (4, 0): 14, (5, 4): 13, (1, 1): 13, (4, 3): 12, (1, 5): 11, (2, 2): 10, (4, 2): 10, (0, 5): 9, (1, 4): 9, (3, 3): 8, (0, 2): 8, (1, 2): 7, (0, 4): 7, (2, 4): 7, (5, 0): 7, (4, 1): 7, (5, 2): 7, (2, 0): 7, (1, 3): 6, (0, 3): 6, (2, 1): 6, (5, 1): 6, (2, 5): 6, (5, 3): 5, (3, 5): 5, (3, 1): 5, (3, 0): 4, (3, 4): 3})] <class 'numpy.ndarray'>\n"
     ]
    }
   ],
   "source": [
    "pub_result = result[0]\n",
    "counts = pub_result.get_counts()\n",
    "print(counts, type(counts))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Instead of counts, you can also query the individual samples of the POVM measurement outcomes with the method `povm_toolbox.sampler.POVMPubResult.get_samples`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[(1, 2), (5, 4), (0, 4), (1, 1), (3, 0)] <class 'list'>\n"
     ]
    }
   ],
   "source": [
    "samples = pub_result.get_samples()\n",
    "print(samples[:5], type(samples))  # the first 5 samples"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compute an expectation value with the `POVMPostProcessor`\n",
    "\n",
    "After having sampled some POVM outcomes, follow these instructions to get the expected value of an observable for a given quantum circuit with the `povm_toolbox.post_processor.POVMPostProcessor`.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Initialize observables\n",
    "The first step is to define the observables whose expected value you want to compute. Each observable can be any `BaseOperator`, like the operators from `qiskit.quantum_info`. Among them it is preferable to use `qiskit.quantum_info.SparsePauliOp`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "from qiskit.quantum_info import SparsePauliOp\n",
    "\n",
    "observable = SparsePauliOp([\"II\", \"XX\", \"YY\", \"ZZ\"], coeffs=[1, 1, -1, 1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Initialize the `POVMPostProcessor`\n",
    "\n",
    "A `POVMPostProcessor` is initialized with an instance of a `POVMPubResult`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "from povm_toolbox.post_processor import POVMPostProcessor\n",
    "\n",
    "post_processor = POVMPostProcessor(pub_result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Get the expected value\n",
    "\n",
    "One can request the expectation value of an observable with the method `povm_toolbox.post_processor.POVMPostProcessor.get_expectation_value`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3.8828125\n"
     ]
    }
   ],
   "source": [
    "exp_value, std = post_processor.get_expectation_value(observable)\n",
    "print(exp_value)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For reference, we can compare our estimated expectation value to the exact value.\n",
    "\n",
    "The post-processor can also compute the (estimated) standard deviation of the estimator we built.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Exact value:     3.999999999999999\n",
      "Estimated value: 3.8828125\n",
      "\n",
      "Estimated standard deviation of the estimator: 0.26297472173936615\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from qiskit.quantum_info import Statevector\n",
    "\n",
    "exact_expectation_value = np.real_if_close(Statevector(qc).expectation_value(observable))\n",
    "print(f\"Exact value:     {exact_expectation_value}\")\n",
    "\n",
    "exp_value, std = post_processor.get_expectation_value(observable)\n",
    "print(f\"Estimated value: {exp_value}\")\n",
    "print(f\"\\nEstimated standard deviation of the estimator: {std}\")"
   ]
  }
 ],
 "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.10.11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}