docs/simulate/noisy_simulation.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "DkA0Fobtb9dM"
   },
   "source": [
    "##### Copyright 2022 The Cirq Developers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "cellView": "form",
    "id": "tUshu7YfcAAW"
   },
   "outputs": [],
   "source": [
    "# @title Licensed under the Apache License, Version 2.0 (the \"License\");\n",
    "# you may not use this file except in compliance with the License.\n",
    "# You may obtain a copy of the License at\n",
    "#\n",
    "# https://www.apache.org/licenses/LICENSE-2.0\n",
    "#\n",
    "# Unless required by applicable law or agreed to in writing, software\n",
    "# distributed under the License is distributed on an \"AS IS\" BASIS,\n",
    "# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.\n",
    "# See the License for the specific language governing permissions and\n",
    "# limitations under the License."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "7zgataJVe0mU"
   },
   "source": [
    "# Noisy Simulation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "HYkRhx2pe2XX"
   },
   "source": [
    "<table class=\"tfo-notebook-buttons\" align=\"left\">\n",
    "  <td>\n",
    "    <a target=\"_blank\" href=\"https://quantumai.google/cirq/simulate/noisy_simulation\"><img src=\"https://quantumai.google/site-assets/images/buttons/quantumai_logo_1x.png\" />View on QuantumAI</a>\n",
    "  </td>\n",
    "  <td>\n",
    "    <a target=\"_blank\" href=\"https://colab.research.google.com/github/quantumlib/Cirq/blob/main/docs/simulate/noisy_simulation.ipynb\"><img src=\"https://quantumai.google/site-assets/images/buttons/colab_logo_1x.png\" />Run in Google Colab</a>\n",
    "  </td>\n",
    "  <td>\n",
    "    <a target=\"_blank\" href=\"https://github.com/quantumlib/Cirq/blob/main/docs/simulate/noisy_simulation.ipynb\"><img src=\"https://quantumai.google/site-assets/images/buttons/github_logo_1x.png\" />View source on GitHub</a>\n",
    "  </td>\n",
    "  <td>\n",
    "    <a href=\"https://storage.googleapis.com/tensorflow_docs/Cirq/docs/simulate/noisy_simulation.ipynb\"><img src=\"https://quantumai.google/site-assets/images/buttons/download_icon_1x.png\" />Download notebook</a>\n",
    "  </td>\n",
    "</table>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "846b32703c5c"
   },
   "outputs": [],
   "source": [
    "try:\n",
    "    import cirq\n",
    "except ImportError:\n",
    "    print(\"installing cirq...\")\n",
    "    !pip install --quiet cirq\n",
    "    print(\"installed cirq.\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "7d3daa944d6b"
   },
   "source": [
    "For simulation, it is useful to have `Gate` objects that enact noisy quantum evolution. Cirq supports modeling noise via *operator sum* representations of noise (these evolutions are also known as quantum operations or quantum dynamical maps). \n",
    "\n",
    "This formalism models evolution of the density matrix $\\rho$ via\n",
    "\n",
    "$$\n",
    "\\rho \\rightarrow \\sum_{k = 1}^{m} A_k \\rho A_k^\\dagger\n",
    "$$\n",
    "\n",
    "where $A_k$ are known as *Kraus operators*. These operators are not necessarily unitary but must satisfy the trace-preserving property\n",
    "\n",
    "$$\n",
    "\\sum_k A_k^\\dagger A_k = I .\n",
    "$$\n",
    "\n",
    "A channel with $m = 1$ unitary Kraus operator is called *coherent* (and is equivalent to a unitary gate operation), otherwise the channel is called *incoherent*. For a given noisy channel, Kraus operators are not necessarily unique. For more details on these operators, see [John Preskill's lecture notes](http://theory.caltech.edu/~preskill/ph219/chap3_15.pdf)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "b989869680d4"
   },
   "source": [
    "## Common channels"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "708b5b720a74"
   },
   "source": [
    "Cirq defines many commonly used quantum channels in [`ops/common_channels.py`](https://github.com/quantumlib/Cirq/blob/main/cirq-core/cirq/ops/common_channels.py). For example, the single-qubit bit-flip channel\n",
    "\n",
    "$$\n",
    "\\rho \\rightarrow (1 - p) \\rho + p X \\rho X\n",
    "$$\n",
    "\n",
    "with parameter $p = 0.1$ can be created as follows."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "386a49be9dc7"
   },
   "outputs": [],
   "source": [
    "import cirq\n",
    "\n",
    "\"\"\"Get a single-qubit bit-flip channel.\"\"\"\n",
    "bit_flip = cirq.bit_flip(p=0.1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "f1f501177e53"
   },
   "source": [
    "To see the Kraus operators of a channel, the `cirq.kraus` protocol can be used. (See the [protocols guide](../build/protocols.ipynb).)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "7635588faafe"
   },
   "outputs": [],
   "source": [
    "for i, kraus in enumerate(cirq.kraus(bit_flip)):\n",
    "    print(f\"Kraus operator {i + 1} is:\\n\", kraus, end=\"\\n\\n\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "3f9ff40f980f"
   },
   "source": [
    "As mentioned, all channels are subclasses of `cirq.Gate`s. As such, they can act on qubits and be used in circuits in the same manner as gates."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "91187990b395"
   },
   "outputs": [],
   "source": [
    "\"\"\"Example of using channels in a circuit.\"\"\"\n",
    "# See the number of qubits a channel acts on.\n",
    "nqubits = bit_flip.num_qubits()\n",
    "print(f\"Bit flip channel acts on {nqubits} qubit(s).\\n\")\n",
    "\n",
    "# Apply the channel to each qubit in a circuit.\n",
    "circuit = cirq.Circuit(bit_flip.on_each(cirq.LineQubit.range(3)))\n",
    "print(circuit)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "59651802c358"
   },
   "source": [
    "Channels can even be controlled."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "9be5383a7cd7"
   },
   "outputs": [],
   "source": [
    "\"\"\"Example of controlling a channel.\"\"\"\n",
    "# Get the controlled channel.\n",
    "controlled_bit_flip = bit_flip.controlled(num_controls=1)\n",
    "\n",
    "# Use it in a circuit.\n",
    "circuit = cirq.Circuit(controlled_bit_flip(*cirq.LineQubit.range(2)))\n",
    "print(circuit)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "755ba122d550"
   },
   "source": [
    "In addition to the bit-flip channel, other common channels predefined in Cirq are shown below. Definitions of these channels can be found in their docstrings - e.g., `help(cirq.depolarize)`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "6a9c9c76fb8f"
   },
   "source": [
    "* `cirq.phase_flip`\n",
    "* `cirq.phase_damp`\n",
    "* `cirq.amplitude_damp`\n",
    "* `cirq.depolarize`\n",
    "* `cirq.asymmetric_depolarize`\n",
    "* `cirq.reset`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "8b22bb323731"
   },
   "source": [
    "For example, the asymmetric depolarizing channel is defined by\n",
    "\n",
    "$$\n",
    "\\rho \\rightarrow (1-p_x-p_y-p_z) \\rho + p_x X \\rho X + p_y Y \\rho Y + p_z Z \\rho Z\n",
    "$$\n",
    "\n",
    "and can be instantiated as follows."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "b6a0966a29f6"
   },
   "outputs": [],
   "source": [
    "\"\"\"Get an asymmetric depolarizing channel.\"\"\"\n",
    "depo = cirq.asymmetric_depolarize(p_x=0.10, p_y=0.05, p_z=0.15)\n",
    "\n",
    "circuit = cirq.Circuit(depo(cirq.LineQubit(0)))\n",
    "print(circuit)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ddbc622c98da"
   },
   "source": [
    "## The `kraus` and `mixture` protocols"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "2dee355d2ae8"
   },
   "source": [
    "We have seen the `cirq.kraus` protocol which returns the Kraus operators of a channel. Some channels have the interpretation of randomly applying a single unitary Kraus operator $U_k$ with probability $p_k$, namely\n",
    "\n",
    "$$\n",
    "\\rho \\rightarrow \\sum_k p_k U_k \\rho U_k^\\dagger {\\rm ~where~} \\sum_k p_k =1 {\\rm ~and~ U_k U_k^\\dagger= I}.\n",
    "$$\n",
    "\n",
    "For example, the bit-flip channel from above\n",
    "\n",
    "$$\n",
    "\\rho \\rightarrow (1 - p) \\rho + p X \\rho X\n",
    "$$\n",
    "\n",
    "can be interpreted as doing nothing (applying identity) with probability $1 - p$ and flipping the bit (applying $X$) with probability $p$. Channels with these interpretations support the `cirq.mixture` protocol. This protocol returns the probabilities and unitary Kraus operators of the channel."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "d4d84c4d4fa7"
   },
   "outputs": [],
   "source": [
    "\"\"\"Example of using the mixture protocol.\"\"\"\n",
    "for prob, kraus in cirq.mixture(bit_flip):\n",
    "    print(f\"With probability {prob}, apply\\n\", kraus, end=\"\\n\\n\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "1ca339c3cbab"
   },
   "source": [
    "Channels that do not have this interpretation do not support the `cirq.mixture` protocol. Such channels apply Kraus operators with probabilities that depend on the state $\\rho$. \n",
    "\n",
    "An example of a channel which does not support the mixture protocol is the amplitude damping channel with parameter $\\gamma$ defined by Kraus operators\n",
    "\n",
    "$$\n",
    "M_0 = \\begin{bmatrix} 1 & 0  \\cr 0 & \\sqrt{1 - \\gamma} \\end{bmatrix} \n",
    "\\text{and }\n",
    "M_1 = \\begin{bmatrix} 0 & \\sqrt{\\gamma} \\cr 0 & 0 \\end{bmatrix} .\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "8e377332d2a7"
   },
   "outputs": [],
   "source": [
    "\"\"\"The amplitude damping channel is an example of a channel without a mixture.\"\"\"\n",
    "channel = cirq.amplitude_damp(0.1)\n",
    "\n",
    "if cirq.has_mixture(channel):\n",
    "    print(f\"Channel {channel} has a _mixture_ or _unitary_ method.\")\n",
    "else:\n",
    "    print(f\"Channel {channel} does not have a _mixture_ or _unitary_ method.\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "6f3ae4337219"
   },
   "source": [
    "To summarize:\n",
    "\n",
    "* Every `Gate` in Cirq supports the `cirq.kraus` protocol.\n",
    "  - If magic method `_kraus_` is not defined, `cirq.kraus` looks for `_mixture_` then for `_unitary_`.\n",
    "* A subset of channels which support `cirq.kraus` also support the `cirq.mixture` protocol.\n",
    "  - If magic method `_mixture_` is not defined, `cirq.mixture` looks for `_unitary_`.\n",
    "* A subset of channels which support `cirq.mixture` also support the `cirq.unitary` protocol."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "67c8ddbc10f3"
   },
   "source": [
    "For concrete examples, consider `cirq.X`, `cirq.BitFlipChannel`, and `cirq.AmplitudeDampingChannel` which are all subclasses of `cirq.Gate`.\n",
    "\n",
    "* `cirq.X` defines the `_unitary_` method. \n",
    "  - As a result, it supports the `cirq.unitary` protocol, the `cirq.mixture` protocol, and the `cirq.kraus` protocol.\n",
    "* `cirq.BitFlipChannel` defines the `_mixture_` method but not the `_unitary_` method.\n",
    "  - As a result, it only supports the `cirq.mixture` protocol and the `cirq.kraus` protocol.\n",
    "* `cirq.AmplitudeDampingChannel` defines the `_kraus_` method, but not the `_mixture_` method or the `_unitary_` method.\n",
    "  - As a result, it only supports the `cirq.kraus` protocol."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "0f5b77825a8a"
   },
   "source": [
    "## Custom channels"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "7614e0351462"
   },
   "source": [
    "There are two configurable channel types for channels not defined in `cirq.ops.common_channels`: `MixedUnitaryChannel` and `KrausChannel`.\n",
    "\n",
    "`MixedUnitaryChannel` takes a list of `(probability, unitary)` tuples and uses it to define the `_mixture_` method.\n",
    "\n",
    "`KrausChannel` takes a list of Kraus operators and uses it to define the `_channel` method.\n",
    "\n",
    "Both types also accept a measurement key as an optional parameter. This key will be used to store the index of the selected unitary or Kraus operator in the measurement results.\n",
    "\n",
    "An example of each type is shown below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "7cbd723fb567"
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "q0 = cirq.LineQubit(0)\n",
    "# This is equivalent to a bit-flip error with probability 0.1.\n",
    "mix = [\n",
    "    (0.9, np.array([[1, 0], [0, 1]], dtype=np.complex64)),\n",
    "    (0.1, np.array([[0, 1], [1, 0]], dtype=np.complex64)),\n",
    "]\n",
    "bit_flip_mix = cirq.MixedUnitaryChannel(mix)\n",
    "\n",
    "# This is equivalent to an X-basis measurement.\n",
    "ops = [np.array([[1, 1], [1, 1]]) * 0.5, np.array([[1, -1], [-1, 1]]) * 0.5]\n",
    "x_meas = cirq.KrausChannel(ops, key='x')\n",
    "\n",
    "# These circuits have the same behavior.\n",
    "circuit = cirq.Circuit(bit_flip_mix.on(q0), cirq.H(q0), x_meas.on(q0))\n",
    "equiv_circuit = cirq.Circuit(\n",
    "    cirq.bit_flip(0.1).on(q0),\n",
    "    cirq.H(q0),\n",
    "    # Measure in x-basis\n",
    "    cirq.H(q0),\n",
    "    cirq.measure(q0, key='x'),\n",
    "    cirq.H(q0),\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "360d5a316769"
   },
   "source": [
    "Alternatively, users can define their own channel types. Defining custom channels is similar to defining [custom gates](../build/custom_gates.ipynb).\n",
    "\n",
    "A minimal example for defining the channel\n",
    "\n",
    "$$\n",
    "\\rho \\mapsto (1 - p) \\rho + p Y \\rho Y\n",
    "$$\n",
    "\n",
    "is shown below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "55240879394e"
   },
   "outputs": [],
   "source": [
    "\"\"\"Minimal example of defining a custom channel.\"\"\"\n",
    "\n",
    "\n",
    "class BitAndPhaseFlipChannel(cirq.Gate):\n",
    "    def _num_qubits_(self) -> int:\n",
    "        return 1\n",
    "\n",
    "    def __init__(self, p: float) -> None:\n",
    "        self._p = p\n",
    "\n",
    "    def _mixture_(self):\n",
    "        ps = [1.0 - self._p, self._p]\n",
    "        ops = [cirq.unitary(cirq.I), cirq.unitary(cirq.Y)]\n",
    "        return tuple(zip(ps, ops))\n",
    "\n",
    "    def _has_mixture_(self) -> bool:\n",
    "        return True\n",
    "\n",
    "    def _circuit_diagram_info_(self, args) -> str:\n",
    "        return f\"BitAndPhaseFlip({self._p})\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "8ffff44ff2b2"
   },
   "source": [
    "Note: The `_has_mixture_` magic method is not strictly required but is recommended."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "2dea6bd62cfe"
   },
   "source": [
    "We can now instantiate this channel and get its mixture:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "0350d9e193ca"
   },
   "outputs": [],
   "source": [
    "\"\"\"Custom channels can be used like any other channels.\"\"\"\n",
    "bit_phase_flip = BitAndPhaseFlipChannel(p=0.05)\n",
    "\n",
    "for prob, kraus in cirq.mixture(bit_phase_flip):\n",
    "    print(f\"With probability {prob}, apply\\n\", kraus, end=\"\\n\\n\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "72486a155131"
   },
   "source": [
    "Note: Since `_mixture_` is defined, the `cirq.kraus` protocol can also be used."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ada775055c3f"
   },
   "source": [
    "The custom channel can be used in a circuit just like other predefined channels."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "1af5c9f60cab"
   },
   "outputs": [],
   "source": [
    "\"\"\"Example of using a custom channel in a circuit.\"\"\"\n",
    "circuit = cirq.Circuit(bit_phase_flip.on_each(*cirq.LineQubit.range(3)))\n",
    "circuit"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "46eb7a01e30d"
   },
   "source": [
    "Note: If a custom channel does not have a mixture, it should instead define the `_kraus_` magic method to return a sequence of Kraus operators (as `numpy.ndarray`s). Defining a `_has_kraus_` method which returns `True` is optional but recommended."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "b537c43e078c"
   },
   "source": [
    "This method of defining custom channels is the most general, but simple channels such as the custom `BitAndPhaseFlipChannel` can also be created directly from a `Gate` with the convenient `Gate.with_probability` method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "0937f8dad808"
   },
   "outputs": [],
   "source": [
    "\"\"\"Create a channel with Gate.with_probability.\"\"\"\n",
    "channel = cirq.Y.with_probability(probability=0.05)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "0461f377ac46"
   },
   "source": [
    "This produces the same mixture as the custom `BitAndPhaseFlip` channel above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "f0a85b33680a"
   },
   "outputs": [],
   "source": [
    "for prob, kraus in cirq.mixture(channel):\n",
    "    print(f\"With probability {prob}, apply\\n\", kraus, end=\"\\n\\n\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "c73ec3e63e23"
   },
   "source": [
    "Note that the order of Kraus operators is reversed from above, but this of course does not affect the action of the channel."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "96a696b543f7"
   },
   "source": [
    "## Simulating noisy circuits"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "d51d993ae869"
   },
   "source": [
    "### Density matrix simulation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "38dfcb60ef79"
   },
   "source": [
    "The `cirq.DensityMatrixSimulator` can simulate any noisy circuit (i.e., can apply any quantum channel) because it stores the full density matrix $\\rho$. This simulation strategy updates the state $\\rho$ by directly applying the Kraus operators of each quantum channel."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "a2973a8faef0"
   },
   "outputs": [],
   "source": [
    "\"\"\"Simulating a circuit with the density matrix simulator.\"\"\"\n",
    "# Get a circuit.\n",
    "qbit = cirq.GridQubit(0, 0)\n",
    "circuit = cirq.Circuit(cirq.X(qbit), cirq.amplitude_damp(0.1).on(qbit))\n",
    "\n",
    "# Display it.\n",
    "print(\"Simulating circuit:\")\n",
    "print(circuit)\n",
    "\n",
    "# Simulate with the density matrix simulator.\n",
    "dsim = cirq.DensityMatrixSimulator()\n",
    "rho = dsim.simulate(circuit).final_density_matrix\n",
    "\n",
    "# Display the final density matrix.\n",
    "print(\"\\nFinal density matrix:\")\n",
    "print(rho)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "ac8a991a426b"
   },
   "source": [
    "Note that the density matrix simulator supports the `run` method which only gives access to measurements as well as the `simulate` method (used above) which gives access to the full density matrix."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "0c659bb20098"
   },
   "source": [
    "### Monte Carlo wavefunction simulation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "6f3caf9cea92"
   },
   "source": [
    "Noisy circuits with arbitrary channels can also be simulated with the `cirq.Simulator`. When simulating such a channel, a single Kraus operator is randomly sampled (according to the probability distribution) and applied to the wavefunction. This method is known as \"Monte Carlo (wavefunction) simulation\" or \"quantum trajectories.\"\n",
    "\n",
    "Note: For channels which do not support the `cirq.mixture` protocol, the probability of applying each Kraus operator depends on the state. In contrast, for channels which do support the `cirq.mixture` protocol, the probability of applying each Kraus operator is independent of the state."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "d80a2bc7ce14"
   },
   "outputs": [],
   "source": [
    "\"\"\"Simulating a noisy circuit via Monte Carlo simulation.\"\"\"\n",
    "# Get a circuit.\n",
    "qbit = cirq.NamedQubit(\"Q\")\n",
    "circuit = cirq.Circuit(cirq.bit_flip(p=0.5).on(qbit))\n",
    "\n",
    "# Display it.\n",
    "print(\"Simulating circuit:\")\n",
    "print(circuit)\n",
    "\n",
    "# Simulate with the cirq.Simulator.\n",
    "sim = cirq.Simulator()\n",
    "psi = sim.simulate(circuit).dirac_notation()\n",
    "\n",
    "# Display the final wavefunction.\n",
    "print(\"\\nFinal wavefunction:\")\n",
    "print(psi)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "4c485bbc1dfd"
   },
   "source": [
    "To see that the output is stochastic, you can run the cell above multiple times. Since $p = 0.5$ in the bit-flip channel, you should get $|0\\rangle$ roughly half the time and $|1\\rangle$ roughly half the time. The `run` method with many repetitions can also be used to see this behavior."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "8143ae0c7a34"
   },
   "outputs": [],
   "source": [
    "\"\"\"Example of Monte Carlo wavefunction simulation with the `run` method.\"\"\"\n",
    "circuit = cirq.Circuit(cirq.bit_flip(p=0.5).on(qbit), cirq.measure(qbit))\n",
    "res = sim.run(circuit, repetitions=100)\n",
    "print(res.histogram(key=qbit))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "0f123f0e3c55"
   },
   "source": [
    "## Adding noise to circuits"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "473ac025a226"
   },
   "source": [
    "Often circuits are defined with just unitary operations, but we want to simulate them with noise. There are several methods for inserting noise in Cirq."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "d0de2c475f12"
   },
   "source": [
    "For any circuit, the `with_noise` method can be called to insert a channel after every moment."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "568530c8707b"
   },
   "outputs": [],
   "source": [
    "\"\"\"One method to insert noise in a circuit.\"\"\"\n",
    "# Define some noiseless circuit.\n",
    "circuit = cirq.testing.random_circuit(qubits=3, n_moments=3, op_density=1, random_state=11)\n",
    "\n",
    "# Display the noiseless circuit.\n",
    "print(\"Circuit without noise:\")\n",
    "print(circuit)\n",
    "\n",
    "# Add noise to the circuit.\n",
    "noisy = circuit.with_noise(cirq.depolarize(p=0.01))\n",
    "\n",
    "# Display it.\n",
    "print(\"\\nCircuit with noise:\")\n",
    "print(noisy)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "c2206a7b7dd4"
   },
   "source": [
    "This circuit can then be simulated using the methods described above."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "029509de9787"
   },
   "source": [
    "The `with_noise` method creates a `cirq.NoiseModel` from its input and adds noise to each moment. A `cirq.NoiseModel` can be explicitly created and used to add noise to a single operation, single moment, or series of moments as follows."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "10a2cd41bfe4"
   },
   "outputs": [],
   "source": [
    "\"\"\"Add noise to an operation, moment, or sequence of moments.\"\"\"\n",
    "# Create a noise model.\n",
    "noise_model = cirq.NoiseModel.from_noise_model_like(cirq.depolarize(p=0.01))\n",
    "\n",
    "# Get a qubit register.\n",
    "qreg = cirq.LineQubit.range(2)\n",
    "\n",
    "# Add noise to an operation.\n",
    "op = cirq.CNOT(*qreg)\n",
    "noisy_op = noise_model.noisy_operation(op)\n",
    "\n",
    "# Add noise to a moment.\n",
    "moment = cirq.Moment(cirq.H.on_each(qreg))\n",
    "noisy_moment = noise_model.noisy_moment(moment, system_qubits=qreg)\n",
    "\n",
    "# Add noise to a sequence of moments.\n",
    "circuit = cirq.Circuit(cirq.H(qreg[0]), cirq.CNOT(*qreg))\n",
    "noisy_circuit = noise_model.noisy_moments(circuit, system_qubits=qreg)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "a440a2aefdcf"
   },
   "source": [
    "Note: In the last two examples, the argument `system_qubits` can be a subset of the qubits in the moment(s)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "61cf807b3737"
   },
   "source": [
    "The output of each \"noisy method\" is a `cirq.OP_TREE` which can be converted to a circuit by passing it into the `cirq.Circuit` constructor. For example, we create a circuit from the `noisy_moment` below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "3e84e3e17421"
   },
   "outputs": [],
   "source": [
    "\"\"\"Creating a circuit from a noisy cirq.OP_TREE.\"\"\"\n",
    "cirq.Circuit(noisy_moment)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "11b18c640767"
   },
   "source": [
    "Another technique is to directly pass a `cirq.NoiseModel`, or a value that can be trivially converted into one, to the density matrix simulator as shown below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "3d16924842bf"
   },
   "outputs": [],
   "source": [
    "\"\"\"Define a density matrix simulator with a `cirq.NOISE_MODEL_LIKE` object.\"\"\"\n",
    "noisy_dsim = cirq.DensityMatrixSimulator(noise=cirq.generalized_amplitude_damp(p=0.1, gamma=0.5))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "0dcd4cf202be"
   },
   "source": [
    "This will not explicitly add channels to the circuit being simulated, but the circuit will be simulated as though these channels were present."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "99ac77a2c896"
   },
   "source": [
    "Other than these general methods, channels can be added to circuits at any moment just as gates are. The channels can be different, be correlated, act on a subset of qubits, be custom defined, etc."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "de8a285b926c"
   },
   "outputs": [],
   "source": [
    "\"\"\"Defining a circuit with multiple noisy channels.\"\"\"\n",
    "qreg = cirq.LineQubit.range(4)\n",
    "circ = cirq.Circuit(\n",
    "    cirq.H.on_each(qreg),\n",
    "    cirq.depolarize(p=0.01).on_each(qreg),\n",
    "    cirq.qft(*qreg),\n",
    "    bit_phase_flip.on_each(qreg[1::2]),\n",
    "    cirq.qft(*qreg, inverse=True),\n",
    "    cirq.reset(qreg[1]),\n",
    "    cirq.measure(*qreg),\n",
    "    cirq.bit_flip(p=0.07).controlled(1).on(*qreg[2:]),\n",
    ")\n",
    "\n",
    "print(\"Circuit with multiple channels:\\n\")\n",
    "print(circ)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "a9bb97a17256"
   },
   "source": [
    "Circuits can also be modified with standard methods like `insert` to add channels at any point in the circuit. For example, to model simple state preparation errors, one can add bit-flip channels to the start of the circuit as follows."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "524ded251565"
   },
   "outputs": [],
   "source": [
    "\"\"\"Example of inserting channels in circuits.\"\"\"\n",
    "circ.insert(0, cirq.bit_flip(p=0.1).on_each(qreg))\n",
    "print(circ)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "de5748df9323"
   },
   "source": [
    "### Simulation with realistic noise\n",
    "\n",
    "Cirq also provides a couple `NoiseModel`s which are designed to mimic the noise behavior seen on Google quantum hardware devices. As of July 19, 2022, models that mimic the Rainbow or Weber Google quantum processors are publicly available in Cirq. You can instantiate these noise models as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "id": "2d83326239b3"
   },
   "outputs": [],
   "source": [
    "import cirq_google\n",
    "\n",
    "processor_id = \"rainbow\"  # or \"weber\"\n",
    "# Load the calibration data\n",
    "cal = cirq_google.engine.load_median_device_calibration(processor_id)\n",
    "# Turn calibration data into a noise properties object\n",
    "noise_props = cirq_google.noise_properties_from_calibration(cal)\n",
    "# Build a noise model from the noise properties\n",
    "noise_model = cirq_google.NoiseModelFromGoogleNoiseProperties(noise_props)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "id": "640b41a79b49"
   },
   "source": [
    "While this `NoiseModel` can be used anywhere other noise models could be used, it is particularly useful in a [Quantum Virtual Machine](./quantum_virtual_machine.ipynb). A QVM combines a realistic noise model and a [Device](../hardware/devices.ipynb) object together and places them behind a `cirq.Engine`-style interface so that you can run circuits almost identically to how you would with a hardware device, and get results that approximate those a hardware device would produce. "
   ]
  }
 ],
 "metadata": {
  "colab": {
   "name": "noisy_simulation.ipynb",
   "toc_visible": true
  },
  "kernelspec": {
   "display_name": "Python 3",
   "name": "python3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}