{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Problem Sheet 1\n",
    "\n",
    "_Quantum Software | Oxford Computer Science 2022_\n",
    "\n",
    "\n",
    "In this problem sheet, we will get the hang of working with the basic generators of circuits and ZX-diagrams concretely, using matrix calculations. For this, we will use the `sympy` library."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sympy import *\n",
    "from sympy.physics.quantum import TensorProduct\n",
    "from fractions import Fraction\n",
    "\n",
    "def T(*args):\n",
    "    if len(args) == 1: return args[0]\n",
    "    elif len(args) == 0: return Matrix([[1]])\n",
    "    else: return TensorProduct(args[0], T(*args[1:]))\n",
    "\n",
    "alpha = var(\"alpha\")\n",
    "beta = var(\"beta\")\n",
    "gamma = var(\"gamma\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With this library, we can construct matrices with `Matrix`. Then `*` is matrix multiplication, and `T` is tensor product.\n",
    "\n",
    "Note `T` takes any number of arguments, e.g. $A \\otimes B \\otimes C$ is `T(A,B,C)`, and if you want to be fancy: $A^{\\otimes n}$ is `T(*n*[A])`. (Python trivia: Why does that work??)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}1 & 2\\\\3 & 4\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "Matrix([\n",
       "[1, 2],\n",
       "[3, 4]])"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}37 & 54\\\\81 & 118\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "Matrix([\n",
       "[37,  54],\n",
       "[81, 118]])"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}1 & 2 & 2 & 4\\\\3 & 4 & 6 & 8\\\\3 & 6 & 4 & 8\\\\9 & 12 & 12 & 16\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "Matrix([\n",
       "[1,  2,  2,  4],\n",
       "[3,  4,  6,  8],\n",
       "[3,  6,  4,  8],\n",
       "[9, 12, 12, 16]])"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "M = Matrix([[1,  2],\n",
    "            [3,  4]])\n",
    "display(M)\n",
    "display(M * M * M)\n",
    "display(T(M,M))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also define some variables with `var` which can be used in mathmatical expressions, and substituted via `.subs(...)`. Variable names can be any string, but `sympy` knows how to pretty-print some variable names, e.g. Greek letters.\n",
    "\n",
    "NOTE: $\\sqrt{2}$ is `sqrt(2)`, $i$ is `I`, $\\pi$ is `pi`, and $e^x$ is `exp(x)`, so phases $e^{i \\alpha}$ are written `exp(i * alpha)`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\alpha$"
      ],
      "text/plain": [
       "alpha"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\displaystyle e^{i \\alpha}$"
      ],
      "text/plain": [
       "exp(I*alpha)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\displaystyle - i$"
      ],
      "text/plain": [
       "-I"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\frac{\\sqrt{2}}{2} + \\frac{\\sqrt{2} i}{2}$"
      ],
      "text/plain": [
       "sqrt(2)/2 + sqrt(2)*I/2"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alpha = var(\"alpha\")\n",
    "phase = exp(I * alpha)\n",
    "epi4 = exp(I * pi / 4)\n",
    "\n",
    "display(alpha)\n",
    "display(phase)\n",
    "display(phase.subs(alpha, -pi / 2))\n",
    "display(phase.subs(alpha, pi / 4) == epi4)\n",
    "display(re(epi4) + I * im(epi4))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "That should be enough to get going. If in doubt, [read the docs](https://docs.sympy.org/latest/index.html)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Question 0\n",
    "\n",
    "First some basics. Define matrices for:\n",
    " * `z0` := $|0\\rangle$, `z1` := $|1\\rangle$, `x0` := $|{+}\\rangle$, `x1` := $|{-}\\rangle$\n",
    " * `bz0` := $\\langle 0|$, `bz1` := $\\langle 1|$, `bx0` := $\\langle {+}|$, `bx1` := $\\langle {-}|$\n",
    " * `w` for the 2D identity matrix (\"wire\")\n",
    " * `swap` for the swap\n",
    "\n",
    "Compute various products and tensor products and show the results are sensible."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}1 & 0 & 0 & 1\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "Matrix([[1, 0, 0, 1]])"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}1\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "Matrix([[1]])"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}0\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "Matrix([[0]])"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "#SOLUTION\n",
    "\n",
    "z0 = Matrix([[1],[0]])\n",
    "z1 = Matrix([[0],[1]])\n",
    "x0 = (1/sqrt(2)) * Matrix([[1], [1]])\n",
    "x1 = (1/sqrt(2)) * Matrix([[1], [-1]])\n",
    "\n",
    "bz0, bz1, bx0, bx1 = (m.adjoint() for m in (z0,z1,x0,x1))\n",
    "w = eye(2)\n",
    "swap = Matrix([[1, 0, 0, 0],\n",
    "               [0, 0, 1, 0],\n",
    "               [0, 1, 0, 0],\n",
    "               [0, 0, 0, 1]])\n",
    "\n",
    "display(T(bz0,bz0) + T(bz1,bz1))\n",
    "display(bz0 * z0)\n",
    "display(bz0 * z1)\n",
    "display(swap * swap == T(w,w))\n",
    "display(T(swap, w) * T(w, swap) * T(swap, w) == T(w, swap) * T(swap, w) * T(w, swap))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Question 1\n",
    "\n",
    "Define a function that produces the matrix of a Z-spider. It should take 3 arguments: `m` for input legs, `n` for output legs, and `phase` for phase. The phase should have a default value of 0.\n",
    "\n",
    "Build this function in (at least) 2 different ways:\n",
    " 1. by building the $2^n \\times 2^m$ matrix of the spider explicitly (call this function `zs`)\n",
    " 2. by using sums, compositions, and tensor products of the generators from the previous question\n",
    "\n",
    "Test you implementation by comparing the matrices for various choices of inputs, outputs, and phases. (Don't forget to check no-legged spiders!)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[[True, True, True, True],\n",
       " [True, True, True, True],\n",
       " [True, True, True, True],\n",
       " [True, True, True, True]]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# SOLUTION\n",
    "\n",
    "def zs(m, n, phase=0):\n",
    "    return Matrix([[\n",
    "        (1 if i == 0 and j == 0 else 0) +\n",
    "        (exp(I * phase) if i == 2**m - 1 and j == 2**n - 1 else 0)\n",
    "      for i in range(2**m)] for j in range(2**n)])\n",
    "\n",
    "def zs_gens(m, n, phase=0):\n",
    "    return (\n",
    "        T(*n*[z0]) * T(*m*[bz0]) +\n",
    "        exp(I * phase) * T(*n*[z1]) * T(*m*[bz1])\n",
    "    )\n",
    "\n",
    "display([[zs_gens(i, j, alpha) == zs(i, j, alpha) for i in range(4)] for j in range(4)])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Question 2\n",
    "\n",
    "Similarly, define a function for X-spiders:\n",
    "\n",
    "1. by first defining the Hadamard gate `had` and using `zs` (call this function `xs`)\n",
    "2. using the generators from Question 0\n",
    "\n",
    "Again, check these two definitions agree for various numbers of legs and angles."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[[True, True, True, True],\n",
       " [True, True, True, True],\n",
       " [True, True, True, True],\n",
       " [True, True, True, True]]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# SOLUTION\n",
    "\n",
    "had = (1 / sqrt(2)) * Matrix([[1, 1], [1, -1]])\n",
    "\n",
    "def xs(m, n, phase=0):\n",
    "    return T(*n*[had]) * zs(m, n, phase) * T(*m*[had])\n",
    "\n",
    "def xs_gens(m, n, phase=0):\n",
    "    return (\n",
    "        T(*n*[x0]) * T(*m*[bx0]) +\n",
    "        exp(I * phase) * T(*n*[x1]) * T(*m*[bx1])\n",
    "    )\n",
    "\n",
    "display([[xs_gens(i, j, alpha) == xs(i, j, alpha) for i in range(4)] for j in range(4)])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Question 3\n",
    "\n",
    "Verify the following equations by matrix calculation:\n",
    "\n",
    "1. spider fusion for two Z-spiders, each having 2 inputs and 2 outputs, connected by a single leg\n",
    "2. the strong complementarity rule with 2 inputs and 2 outputs\n",
    "3. the Euler decomposition of the Hadamard\n",
    "\n",
    "In each case, find the right scalar factor so the matrices of the LHS and RHS agree exactly.\n",
    "\n",
    "**Hint:** If you have some variables, you may need to call the `simplify(...)` function before comparing matrices. If a constant phase is stuck and won't simplify any more, try splitting the real and imaginary parts and putting back together, i.e. `expr1 = re(expr) + I * im(expr)`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# SOLUTION\n",
    "\n",
    "sf_lhs = T(zs(2, 2, beta), w) * T(w, zs(2, 2, alpha))\n",
    "sf_rhs = zs(3, 3, alpha + beta)\n",
    "display(simplify(sf_lhs) == sf_rhs)\n",
    "\n",
    "sc_lhs = zs(1,2) * xs(2, 1)\n",
    "sc_rhs = T(xs(2, 1), xs(2, 1)) * T(w, swap, w) * T(zs(1,2), zs(1,2))\n",
    "display(sc_lhs == sqrt(2) * sc_rhs)\n",
    "\n",
    "euler_rhs = exp(-I * pi/4) * zs(1,1,pi/2) * xs(1,1,pi/2) * zs(1,1,pi/2)\n",
    "display(had == re(euler_rhs) + I * im(euler_rhs))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Question 4\n",
    "\n",
    "* Make an XOR using an X spider and a scalar multiplication, and show it acts like XOR on $\\{|0\\rangle, |1\\rangle\\}$.\n",
    "* Make a CNOT and TONC (a CNOT with a control on the second qubit and target on the first) using Z spiders, X spiders, swaps, and scalar multiplication. Check that CNOT acts as expected on basis states and `swap * CNOT * swap` is the same as `TONC`\n",
    "* Check XOR is associative and `CNOT * TONC * CNOT = swap`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "[True, True, True, True]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "[True, True, True, True]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# SOLUTION\n",
    "\n",
    "xor = sqrt(2) * xs(2, 1)\n",
    "cnot = T(w, xor) * T(zs(1,2), w)\n",
    "tonc = T(xor, w) * T(w, zs(1,2))\n",
    "\n",
    "display(tonc == swap * cnot * swap)\n",
    "display([\n",
    "    xor * T(z0, z0) == z0,\n",
    "    xor * T(z0, z1) == z1,\n",
    "    xor * T(z1, z0) == z1,\n",
    "    xor * T(z1, z1) == z0,\n",
    "])\n",
    "display([\n",
    "    cnot * T(z0, z0) == T(z0, z0),\n",
    "    cnot * T(z0, z1) == T(z0, z1),\n",
    "    cnot * T(z1, z0) == T(z1, z1),\n",
    "    cnot * T(z1, z1) == T(z1, z0),\n",
    "])\n",
    "display(xor * T(w, xor) == xor * T(xor, w))\n",
    "display(cnot * tonc * cnot == swap)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Question 5\n",
    "\n",
    "Build a 2-qubit parity-phase gate, i.e. a gate that acts as: $|x, y\\rangle \\mapsto e^{i (x \\oplus y) \\alpha} |x, y\\rangle$ in (at least) two ways:\n",
    "\n",
    "1. using CNOT gates and a Z phase gate\n",
    "2. using some Z-spiders and XOR\n",
    "\n",
    "Check the two definitions coincide and act correctly on Z basis states."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[True, True, True, True]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# SOLUTION\n",
    "\n",
    "pp1 = cnot * T(w, zs(1,1,alpha)) * cnot\n",
    "pp2 = T(w, zs(1, 0, alpha) * xor, w) * T(zs(1, 2), zs(1,2))\n",
    "\n",
    "display([\n",
    "    pp1 * T(z0, z0) == T(z0, z0),\n",
    "    pp1 * T(z0, z1) == exp(I * alpha) * T(z0, z1),\n",
    "    pp1 * T(z1, z0) == exp(I * alpha) * T(z1, z0),\n",
    "    pp1 * T(z1, z1) == T(z1, z1),\n",
    "])\n",
    "\n",
    "display(pp1 == pp2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Question 6\n",
    "\n",
    "Define a function that constructs an N-qubit parity phase gate $|x_1 ... x_n\\rangle \\mapsto e^{i (x_1 \\oplus ... \\oplus x_n) \\alpha} |x_1 ... x_n\\rangle$ either concretely or from gates and spiders. Do this however you like, and preferably in multiple ways. Check your answers, and try to come up with cool answers. :)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# SOLUTION\n",
    "\n",
    "def bits(k, n):\n",
    "    return [int(b) for b in bin(k)[2:].zfill(n)]\n",
    "\n",
    "def ppn(n, phase):\n",
    "    return diag(*[exp(I * phase) if sum(bits(i, n))%2 == 1 else 1 for i in range(2**n)])\n",
    "\n",
    "def ppn_rec(n, phase):\n",
    "    if n == 1:\n",
    "        return zs(1,1,phase)\n",
    "    else:\n",
    "        return T(cnot, eye(2**(n-2))) * T(w, ppn_rec(n-1, phase)) * T(cnot, eye(2**(n-2)))\n",
    "    \n",
    "ppn_rec(4, alpha) == ppn(4, alpha)"
   ]
  }
 ],
 "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}