{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"../../../images/qiskit_header.png\" alt=\"Note: In order for images to show up in this jupyter notebook you need to select File => Trusted Notebook\" align=\"middle\">"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Operators Overview\n",
    "\n",
    "The latest version of this notebook is available on https://github.com/Qiskit/qiskit-tutorial."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction\n",
    "\n",
    "This notebook shows how to use the `Operator` class from the *Quantum Information* module of Qiskit to create custom matrix operators and custom unitary gates, and to evaluate the unitary matrix for a quantum circuit."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-09-29T00:15:24.371649Z",
     "start_time": "2018-09-29T00:15:22.358409Z"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "from qiskit import QuantumCircuit, ClassicalRegister, QuantumRegister\n",
    "from qiskit import execute, BasicAer\n",
    "from qiskit.compiler import transpile\n",
    "from qiskit.quantum_info.operators import Operator, Pauli\n",
    "from qiskit.quantum_info import process_fidelity\n",
    "\n",
    "from qiskit.extensions import RXGate, CnotGate, XGate"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Operator Class\n",
    "\n",
    "The `Operator` class is used in Qiskit to represent matrix operators acting on a quantum system. It has several methods to build composite operators using tensor products of smaller operators, and to compose operators.\n",
    "\n",
    "### Creating Operators\n",
    "\n",
    "The easiest way to create an operator object is to initialize it with a matrix given as a list or a Numpy array. For example, to create a two-qubit Pauli-XX operator:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[0.+0.j 0.+0.j 0.+0.j 1.+0.j]\n",
       " [0.+0.j 0.+0.j 1.+0.j 0.+0.j]\n",
       " [0.+0.j 1.+0.j 0.+0.j 0.+0.j]\n",
       " [1.+0.j 0.+0.j 0.+0.j 0.+0.j]], input_dims=(2, 2), output_dims=(2, 2))"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "XX = Operator([[0, 0, 0, 1], [0, 0, 1, 0], [0, 1, 0, 0], [1, 0, 0, 0]])\n",
    "XX"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Operator Properties\n",
    "\n",
    "The operator object stores the underlying matrix, and the input and output dimension of subsystems. \n",
    "\n",
    "* `data`: To access the underly Numpy array, we may use the `Operator.data` property.\n",
    "* `dims`: To return the total input and output dimension of the operator, we may use the `Operator.dim` property. *Note: the output is returned as a tuple `(input_dim, output_dim)`, which is the reverse of the shape of the underlying matrix.*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0.+0.j, 0.+0.j, 0.+0.j, 1.+0.j],\n",
       "       [0.+0.j, 0.+0.j, 1.+0.j, 0.+0.j],\n",
       "       [0.+0.j, 1.+0.j, 0.+0.j, 0.+0.j],\n",
       "       [1.+0.j, 0.+0.j, 0.+0.j, 0.+0.j]])"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "XX.data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(4, 4)"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "input_dim, output_dim = XX.dim\n",
    "input_dim, output_dim"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Input and Output Dimensions\n",
    "\n",
    "The operator class also keeps track of subsystem dimensions, which can be used for composing operators together. These can be accessed using the `input_dims` and `output_dims` functions.\n",
    "\n",
    "For $2^N$ by $2^M$ operators, the input and output dimension will be automatically assumed to be M-qubit and N-qubit:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Input dimensions: (2, 2)\n",
      "Output dimensions: (2,)\n"
     ]
    }
   ],
   "source": [
    "op = Operator(np.random.rand(2 ** 1, 2 ** 2))\n",
    "print('Input dimensions:', op.input_dims())\n",
    "print('Output dimensions:', op.output_dims())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If the input matrix is not divisible into qubit subsystems, then it will be stored as a single-qubit operator. For example, if we have a $6\\times6$ matrix:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Input dimensions: (6,)\n",
      "Output dimensions: (6,)\n"
     ]
    }
   ],
   "source": [
    "op = Operator(np.random.rand(6, 6))\n",
    "print('Input dimensions:', op.input_dims())\n",
    "print('Output dimensions:', op.output_dims())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The input and output dimension can also be manually specified when initializing a new operator:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Input dimensions: (4,)\n",
      "Output dimensions: (2,)\n"
     ]
    }
   ],
   "source": [
    "# Force input dimension to be (4,) rather than (2, 2)\n",
    "op = Operator(np.random.rand(2 ** 1, 2 ** 2), input_dims=[4])\n",
    "print('Input dimensions:', op.input_dims())\n",
    "print('Output dimensions:', op.output_dims())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Input dimensions: (2, 3)\n",
      "Output dimensions: (2, 3)\n"
     ]
    }
   ],
   "source": [
    "# Specify system is a qubit and qutrit\n",
    "op = Operator(np.random.rand(6, 6),\n",
    "              input_dims=[2, 3], output_dims=[2, 3])\n",
    "print('Input dimensions:', op.input_dims())\n",
    "print('Output dimensions:', op.output_dims())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also extract just the input or output dimensions of a subset of subsystems using the `input_dims` and `output_dims` functions:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Dimension of input system 0: (2,)\n",
      "Dimension of input system 1: (3,)\n"
     ]
    }
   ],
   "source": [
    "print('Dimension of input system 0:', op.input_dims([0]))\n",
    "print('Dimension of input system 1:', op.input_dims([1]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Converting classes to Operators\n",
    "\n",
    "Several other classes in Qiskit can be directly converted to an `Operator` object using the operator initialization method. For example:\n",
    "\n",
    "* `Pauli` objects\n",
    "* `Gate` and `Instruction` objects\n",
    "* `QuantumCircuits` objects\n",
    "\n",
    "Note that the last point means we can use the `Operator` class as a unitary simulator to compute the final unitary matrix for a quantum circuit, without having to call a simulator backend. If the circuit contains any unsupported operations, an exception will be raised. Unsupported operations are: measure, reset, conditional operations, or a gate that does not have a matrix definition or decomposition in terms of gate with matrix definitions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[0.+0.j 0.+0.j 0.+0.j 1.+0.j]\n",
       " [0.+0.j 0.+0.j 1.+0.j 0.+0.j]\n",
       " [0.+0.j 1.+0.j 0.+0.j 0.+0.j]\n",
       " [1.+0.j 0.+0.j 0.+0.j 0.+0.j]], input_dims=(2, 2), output_dims=(2, 2))"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Create an Operator from a Pauli object\n",
    "\n",
    "pauliXX = Pauli(label='XX')\n",
    "Operator(pauliXX)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[1.+0.j 0.+0.j 0.+0.j 0.+0.j]\n",
       " [0.+0.j 0.+0.j 0.+0.j 1.+0.j]\n",
       " [0.+0.j 0.+0.j 1.+0.j 0.+0.j]\n",
       " [0.+0.j 1.+0.j 0.+0.j 0.+0.j]], input_dims=(2, 2), output_dims=(2, 2))"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Create an Operator for a Gate object\n",
    "Operator(CnotGate())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[ 7.07106781e-01+0.j         -4.32978028e-17-0.70710678j]\n",
       " [ 4.32978028e-17-0.70710678j  7.07106781e-01+0.j        ]], input_dims=(2,), output_dims=(2,))"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Create an operator from a parameterized Gate object\n",
    "Operator(RXGate(np.pi / 2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[ 0.70710678+0.j  0.70710678+0.j  0.        +0.j ...  0.        +0.j\n",
       "   0.        +0.j  0.        +0.j]\n",
       " [ 0.        +0.j  0.        +0.j  0.70710678+0.j ...  0.        +0.j\n",
       "   0.        +0.j  0.        +0.j]\n",
       " [ 0.        +0.j  0.        +0.j  0.        +0.j ...  0.        +0.j\n",
       "   0.        +0.j  0.        +0.j]\n",
       " ...\n",
       " [ 0.        +0.j  0.        +0.j  0.        +0.j ...  0.        +0.j\n",
       "   0.        +0.j  0.        +0.j]\n",
       " [ 0.        +0.j  0.        +0.j  0.70710678+0.j ...  0.        +0.j\n",
       "   0.        +0.j  0.        +0.j]\n",
       " [ 0.70710678+0.j -0.70710678+0.j  0.        +0.j ...  0.        +0.j\n",
       "   0.        +0.j  0.        +0.j]], input_dims=(2, 2, 2, 2, 2, 2, 2, 2, 2, 2), output_dims=(2, 2, 2, 2, 2, 2, 2, 2, 2, 2))"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Create an operator from a QuantumCircuit object\n",
    "qr = QuantumRegister(10)\n",
    "circ = QuantumCircuit(qr)\n",
    "circ.h(qr[0])\n",
    "for j in range(1, 10):\n",
    "    circ.cx(qr[j-1], qr[j])\n",
    "\n",
    "# Convert circuit to an operator by implicit unitary simulation\n",
    "Operator(circ)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using Operators in circuits\n",
    "\n",
    "Unitary `Operators` can be directly inserted into a `QuantumCircuit` using the `QuantumCircuit.append` method. This converts the `Operator` into a `UnitaryGate` object, which is added to the circuit.\n",
    "\n",
    "If the operator is not unitary, an exception will be raised. This can be checked using the `Operator.is_unitary()` function, which will return `True` if the operator is unitary and `False` otherwise."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Operatator is unitary: True\n",
      "         ┌──────────┐┌─┐   \n",
      "q1_0: |0>┤0         ├┤M├───\n",
      "         │  unitary │└╥┘┌─┐\n",
      "q1_1: |0>┤1         ├─╫─┤M├\n",
      "         └──────────┘ ║ └╥┘\n",
      " c0_0: 0 ═════════════╩══╬═\n",
      "                         ║ \n",
      " c0_1: 0 ════════════════╩═\n",
      "                           \n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "{'11': 1024}"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Create an operator\n",
    "XX = Operator(Pauli(label='XX'))\n",
    "\n",
    "# Check unitary\n",
    "print('Operatator is unitary:', XX.is_unitary())\n",
    "\n",
    "# Add to a circuit\n",
    "qr = QuantumRegister(2)\n",
    "cr = ClassicalRegister(2)\n",
    "circ = QuantumCircuit(qr, cr)\n",
    "circ.append(XX, [qr[0], qr[1]])\n",
    "circ.measure(qr, cr)\n",
    "print(circ)\n",
    "\n",
    "backend = BasicAer.get_backend('qasm_simulator')\n",
    "job = execute(circ, backend, basis_gates=['u1','u2','u3','cx'])\n",
    "job.result().get_counts(0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that in the above example we initialize the operator from a `Pauli` object. However, the `Pauli` object may also be directly inserted into the circuit itself and will be converted into a sequence of single-qubit Pauli gates:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "         ┌───────────┐┌─┐   \n",
      "q1_0: |0>┤0          ├┤M├───\n",
      "         │  Pauli:XX │└╥┘┌─┐\n",
      "q1_1: |0>┤1          ├─╫─┤M├\n",
      "         └───────────┘ ║ └╥┘\n",
      " c0_0: 0 ══════════════╩══╬═\n",
      "                          ║ \n",
      " c0_1: 0 ═════════════════╩═\n",
      "                            \n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "{'11': 1024}"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Add to a circuit\n",
    "circ2 = QuantumCircuit(qr, cr)\n",
    "circ2.append(Pauli(label='XX'), [qr[0], qr[1]])\n",
    "circ2.measure(qr, cr)\n",
    "print(circ2)\n",
    "\n",
    "# Simulate\n",
    "job2 = execute(circ2, backend)\n",
    "job2.result().get_counts(0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Combining Operators\n",
    "\n",
    "Operators my be combined using several methods. \n",
    "\n",
    "### Tensor Product\n",
    "\n",
    "Two operators $A$ and $B$ may be combined into a tensor product operator $A\\otimes B$ using the `Operator.tensor` function. Note that if both A and B are single-qubit operators, then `A.tensor(B)` = $A\\otimes B$ will have the subsystems indexed as matrix B  on subsystem 0, and matrix $A$ on subsystem 1."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[ 0.+0.j  0.+0.j  1.+0.j  0.+0.j]\n",
       " [ 0.+0.j -0.+0.j  0.+0.j -1.+0.j]\n",
       " [ 1.+0.j  0.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j -1.+0.j  0.+0.j -0.+0.j]], input_dims=(2, 2), output_dims=(2, 2))"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A = Operator(Pauli(label='X'))\n",
    "B = Operator(Pauli(label='Z'))\n",
    "A.tensor(B)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Tensor Expansion\n",
    "\n",
    "A closely related operation is `Operator.expand`, which acts like a tensor product but in the reverse order. Hence, for two operators $A$ and $B$ we have `A.expand(B)` = $B\\otimes A$ where the subsystems indexed as matrix A  on subsystem 0, and matrix $B$ on subsystem 1."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[ 0.+0.j  1.+0.j  0.+0.j  0.+0.j]\n",
       " [ 1.+0.j  0.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j  0.+0.j -0.+0.j -1.+0.j]\n",
       " [ 0.+0.j  0.+0.j -1.+0.j -0.+0.j]], input_dims=(2, 2), output_dims=(2, 2))"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A = Operator(Pauli(label='X'))\n",
    "B = Operator(Pauli(label='Z'))\n",
    "A.expand(B)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Composition\n",
    "\n",
    "We can also compose two operators $A$ and $B$ to implement matrix multiplication using the `Operator.compose` method. We have that `A.compose(B)` returns the operator with matrix $B.A$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[ 0.+0.j  1.+0.j]\n",
       " [-1.+0.j  0.+0.j]], input_dims=(2,), output_dims=(2,))"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A = Operator(Pauli(label='X'))\n",
    "B = Operator(Pauli(label='Z'))\n",
    "A.compose(B)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also compose in the reverse order by applying $B$ in front of $A$ using the `front` kwarg of `compose`:  `A.compose(B, front=True)` = $A.B$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[ 0.+0.j -1.+0.j]\n",
       " [ 1.+0.j  0.+0.j]], input_dims=(2,), output_dims=(2,))"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A = Operator(Pauli(label='X'))\n",
    "B = Operator(Pauli(label='Z'))\n",
    "A.compose(B, front=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Subsystem Composition\n",
    "\n",
    "Note that the previous compose requires that the total output dimension of the first operator $A$ is equal to total input dimension of the composed operator $B$ (and similarly, the output dimension of $B$ must be equal to the input dimension of $A$ when composing with `front=True`).\n",
    "\n",
    "We can also compose a smaller operator with a selection of subsystems on a larger operator using the `qargs` kwarg of `compose`, either with or without `front=True`. In this case, the relevant input and output dimenions of the subsystems being composed must match. *Note that the smaller operator must always be the argument of `compose` method.*\n",
    "\n",
    "For example, to compose a two-qubit gate with a three-qubit Operator:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  1.+0.j  0.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j -1.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  1.+0.j  0.+0.j]\n",
       " [ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j -1.+0.j]\n",
       " [ 1.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j -1.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j  0.+0.j  1.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j  0.+0.j  0.+0.j -1.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]], input_dims=(2, 2, 2), output_dims=(2, 2, 2))"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Compose XZ with an 3-qubit identity operator\n",
    "op = Operator(np.eye(2 ** 3))\n",
    "XZ = Operator(Pauli(label='XZ'))\n",
    "op.compose(XZ, qargs=[0, 2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  1.+0.j  0.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j -1.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  1.+0.j  0.+0.j]\n",
       " [ 0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j -1.+0.j]\n",
       " [ 1.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j -1.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j  0.+0.j  1.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]\n",
       " [ 0.+0.j  0.+0.j  0.+0.j -1.+0.j  0.+0.j  0.+0.j  0.+0.j  0.+0.j]], input_dims=(2, 2, 2), output_dims=(2, 2, 2))"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Compose YX in front of the previous operator\n",
    "op = Operator(np.eye(2 ** 3))\n",
    "YX = Operator(Pauli(label='YX'))\n",
    "op.compose(XZ, qargs=[0, 2], front=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Linear combinations\n",
    "\n",
    "Operators may also be combined using standard linear operators for addition, subtraction and scalar multiplication by complex numbers. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[-1.5+0.j  0. +0.j  0. +0.j  0. +0.j]\n",
       " [ 0. +0.j  1.5+0.j  1. +0.j  0. +0.j]\n",
       " [ 0. +0.j  1. +0.j  1.5+0.j  0. +0.j]\n",
       " [ 0. +0.j  0. +0.j  0. +0.j -1.5+0.j]], input_dims=(2, 2), output_dims=(2, 2))"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "XX = Operator(Pauli(label='XX'))\n",
    "YY = Operator(Pauli(label='YY'))\n",
    "ZZ = Operator(Pauli(label='ZZ'))\n",
    "\n",
    "op = 0.5 * (XX + YY - 3 * ZZ)\n",
    "op"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "An important point is that while `tensor`, `expand` and `compose` will preserve the unitarity of unitary operators, linear combinations will not; hence, adding two unitary operators will, in general, result in a non-unitary operator:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "op.is_unitary()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Implicit Conversion to Operators\n",
    "\n",
    "Note that for all the following methods, if the second object is not already an `Operator` object, it will be implicitly converted into one by the method. This means that matrices can be passed in directly without being explicitly converted to an `Operator` first. If the conversion is not possible, an exception will be raised."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Operator([[0.+0.j 1.+0.j]\n",
       " [1.+0.j 0.+0.j]], input_dims=(2,), output_dims=(2,))"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Compose with a matrix passed as a list\n",
    "Operator(np.eye(2)).compose([[0, 1], [1, 0]])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Comparison of Operators\n",
    "\n",
    "Operators implement an equality method that can be used to check if two operators are approximately equal. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Operator(Pauli(label='X')) == Operator(XGate())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that this checks that each matrix element of the operators is approximately equal; two unitaries that differ by a global phase will not be considered equal:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Operator(XGate()) == np.exp(1j * 0.5) * Operator(XGate())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Process Fidelity\n",
    "\n",
    "We may also compare operators using the `process_fidelity` function from the *Quantum Information* module. This is an information theoretic quantity for how close two quantum channels are to each other, and in the case of unitary operators it does not depend on global phase."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Process fidelity = (1+0j)\n"
     ]
    }
   ],
   "source": [
    "# Two operators which differ only by phase\n",
    "op_a = Operator(XGate()) \n",
    "op_b = np.exp(1j * 0.5) * Operator(XGate())\n",
    "\n",
    "# Compute process fidelity\n",
    "F = process_fidelity(op_a, op_b)\n",
    "print('Process fidelity =', F)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that process fidelity is generally only a valid measure of closeness if the input operators are unitary (or CP in the case of quantum channels), and an exception will be raised if the inputs are not CP."
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "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.3"
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}