{
"cells": [
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"# This cell is added by sphinx-gallery\n# It can be customized to whatever you like\n%matplotlib inline"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"------------------------------------------------------------------------\n\n::: {.meta}\n:property=\\\"og:description\\\": Investigation of the toric code degenerate\nground state and anyon excitations :property=\\\"og:image\\\":\n\n:::\n\n*Author: Christina Lee. Posted: 08 August 2022.*\n\nIntroduction\n============\n\nThe *toric code model* is a treasure trove of interesting physics and\nmathematics. The model sparked the development of the error-correcting\nsurface codes, an essential category of error correction models. But why\nis the model so useful for error correction?\n\nWe delve into mathematics and condensed matter physics to answer that\nvery question. Viewing the model as a description of spins in an exotic\nmagnet allows us to start analyzing the model as a material. What kind\nof material is it? The toric code is an example of a topological state\nof matter.\n\nA state of matter, or phase, cannot become a different phase without\nsome discontinuity in the physical properties as coefficients in the\nHamiltonian change. This discontinuity may exist in an arbitrary order\nderivative or non-local observable. For example, ice cannot become water\nwithout a discontinuity in density as the temperature changes. The\nground state of a **topological** state of matter cannot smoothly deform\nto a non-entangled state without a phase transition. Entanglement, and\nmore critically *long-range* entanglement, is a key hallmark of a\ntopological state of matter.\n\nLocal measurements cannot detect topological states of matter. We have\nto consider the entire system to determine a topological phase. To\nbetter consider this type of property, consider the parity of the number\nof dancers on a dance floor. Does everyone have a partner, or is there\nan odd person out? To measure that, we have to look at the entire\nsystem.\n\nTopology is the study of global properties that are preserved under\ncontinuous deformations. For example, a coffee cup is equivalent to a\ndonut because they both have a single hole. More technically, they both\nhave an [Euler\ncharacteristic](https://en.wikipedia.org/wiki/Euler_characteristic) of\nzero. When we zoom to a local patch, both a sphere and a torus look the\nsame. Only by considering the object as a whole can you detect the\nsingle hole.\n\n![A donut can be smoothly deformed into a\nmug.](../demonstrations/toric_code/torus_to_cup.png){.align-center\nwidth=\"70.0%\"}\n\nIn this demo, we will look at the degenerate ground state and the\nexcitations of the toric code model. The toric code was initially\nproposed in \"Fault-tolerant quantum computation by anyons\" by Kitaev.\nThis demo was inspired by \"Realizing topologically ordered states on a\nquantum processor\" by K. J. Satzinger et al. For further reading, I\nrecommend \"Quantum Spin Liquids\" by Lucile Savary and Leon Balents and\n\\\"A Pedagogical Overview on 2D and 3D Toric Codes and the Origin of\ntheir Topological Orders\\\".\n\nThe Model\n=========\n\nWhat is the source of all this fascinating physics? The Hamiltonian is:\n\n$$\\mathcal{H} = -\\sum_s S_s - \\sum_p P_p,$$\n\nwhere\n\n$$S_s = \\prod_{i \\in s} Z_i \\quad P_p = \\prod_{j \\in p} X_j.$$\n\nIn the literature, the $S_s$ terms are called the \"star\" operators, and\nthe $P_p$ terms are called the \"plaquette\" operators. Each star $s$ and\nplaquette $p$ is a group of 4 sites on a square lattice. You can compare\nthis model to something like the [Heisenberg\nmodel](https://en.wikipedia.org/wiki/Quantum_Heisenberg_model) used to\ndescribe spins interacting magnetically in a material.\n\nIn the most common formulation of the model, sites live on the edges of\na square lattice. In this formulation, the \"plaquette\" operators are\nproducts of Pauli X operators on all the sites in a square, and the\n\\\"star\\\" operators are products of Pauli Z operators on all the sites\nbordering a vertex.\n\nInstead, we can view the model as a checkerboard of alternating square\ntypes. In this formulation, all sites $i$ and $j$ are the vertices of a\nsquare lattice. Each square is a group of four sites, and adjacent\nsquares alternate between the two types of groups. Since the groups on\nthis checkerboard no longer look like stars and plaquettes, we will call\nthem the \"Z Group\" and \"X Group\" operators in this tutorial.\n\n![On the left, sites are grouped into stars around vertices and\nplaquettes on the faces. On the right, we view the lattice as a\ncheckerboard of alternating types of\ngroups.](../demonstrations/toric_code/stars_plaquettes2.png){.align-center\nwidth=\"70.0%\"}\n\nWe will be embedding the lattice on a torus via periodic boundary\nconditions. Periodic boundary conditions basically \"glue\" the bottom of\nthe lattice to the top of the lattice and the left to the right.\n\nModular arithmetic accomplishes this matching. Any site at `(x,y)` is\nthe same as a site at `(x+width, y+height)`.\n\n![By matching up the edges with periodic boundary conditions, we turn a\nsquare grid into a\ntorus.](../demonstrations/toric_code/converting_to_torus.png){.align-center\nwidth=\"70.0%\"}\n\nOn to some practical coding! First, let\\'s define the sites on a\n$4\\times 6$ lattice.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"import pennylane as qml\nimport matplotlib.pyplot as plt\nfrom matplotlib.patches import Polygon, Patch\nfrom itertools import product\nfrom dataclasses import dataclass\n\nimport numpy as np\n\nnp.set_printoptions(suppress=True)\n\nheight = 4\nwidth = 6\n\nall_sites = [(i, j) for i, j in product(range(width), range(height))]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We would like for our wire labels to match the sites. To do this, we\nwill be using an [Immutable Data\nClass](https://realpython.com/python-data-classes/#immutable-data-classes)\n. PennyLane allows wire labels to be any **hashable** object, but\niterable wire labels are currently not supported. Therefore we use a\nfrozen dataclass to represent individual wires by a row and column\nposition.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"@dataclass(frozen=True)\nclass Wire:\n i: int\n j: int\n\nexample_wire = Wire(0, 0)\nprint(\"Example wire: \", example_wire)\nprint(\"At coordinates: \", example_wire.i, example_wire.j)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Setting up Operators\n====================\n\nFor each type of group operator (X and Z), we will have two different\nlists: the \"sites\" and the \"ops\". The \"sites\" are tuples and will\ninclude virtual sites off the edge of the lattice that match up with\nlocations on the other side. For example, the site `(6, 1)` denotes the\nreal location `(0,1)`. We will use the `zgroup_sites` and `xgroup_sites`\nlists to help us view the measurements of the corresponding operators.\n\nThe \\\"ops\\\" list will contain the tensor observables. We will later take\nthe expectation value of each tensor.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"mod = lambda s: Wire(s[0] % width, s[1] % height)\n\nzgroup_sites = [] # list of sites in each group\nzgroup_ops = [] # list of operators for each group\n\nfor x, y in product(range(width // 2), range(height)):\n\n x0 = 2 * x + (y + 1) % 2 # x starting coordinate\n\n sites = [(x0, y), (x0 + 1, y), (x0 + 1, y + 1), (x0, y + 1)]\n\n op = qml.operation.Tensor(*(qml.PauliZ(mod(s)) for s in sites))\n\n zgroup_sites.append(sites)\n zgroup_ops.append(op)\n\nprint(\"First set of sites: \", zgroup_sites[0])\nprint(\"First operator: \", zgroup_ops[0])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We will later use the X Group operator sites to prepare the ground\nstate, so the order here is important. One group needs a slightly\ndifferent order due to interference with the periodic boundary\ncondition.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"xgroup_sites = []\nxgroup_ops = []\nfor x, y in product(range(width // 2), range(height)):\n x0 = 2 * x + y % 2 # lower x coordinate\n\n sites = [(x0 + 1, y + 1), (x0, y + 1), (x0, y), (x0 + 1, y)]\n\n if x == 2 and y == 1: # change order for state prep later\n sites = sites[1:] + sites[0:1]\n\n op = qml.operation.Tensor(*(qml.PauliX(mod(s)) for s in sites))\n\n xgroup_sites.append(sites)\n xgroup_ops.append(op)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"How can we best visualize these groups of four sites?\n\nWe use `matplotlib` to show each group of four sites as a Polygon patch,\ncoloured according to the type of group. The `misc_plot_formatting`\nfunction performs some minor styling improvements repeated throughout\nthis demo. The dotted horizontal lines and dashed vertical lines denote\nwhere we glue our boundaries together.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"def misc_plot_formatting(fig, ax):\n plt.hlines([-0.5, height - 0.5], -0.5, width - 0.5, linestyle=\"dotted\", color=\"black\")\n plt.vlines([-0.5, width - 0.5], -0.5, height - 0.5, linestyle=\"dashed\", color=\"black\")\n plt.xticks(range(width + 1), [str(i % width) for i in range(width + 1)])\n plt.yticks(range(height + 1), [str(i % height) for i in range(height + 1)])\n\n for direction in [\"top\", \"right\", \"bottom\", \"left\"]:\n ax.spines[direction].set_visible(False)\n\n return fig, ax\n\n\nfig, ax = plt.subplots()\nfig, ax = misc_plot_formatting(fig, ax)\n\nfor group in xgroup_sites:\n x_patch = ax.add_patch(Polygon(group, color=\"lavender\", zorder=0))\n\nfor group in zgroup_sites:\n z_patch = ax.add_patch(Polygon(group, color=\"mistyrose\", zorder=0))\n\nplt_sites = ax.scatter(*zip(*all_sites))\n\nplt.legend([x_patch, z_patch, plt_sites], [\"XGroup\", \"ZGroup\", \"Site\"], loc=\"upper left\")\n\nplt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The Ground State\n================\n\nWhile individual X and Z operators do not commute with each other, the X\nGroup and Z Group operators do:\n\n$$[S_s, P_p] = 0.$$\n\nSince they commute, the wavefunction can be an eigenstate of each group\noperator independently. To minimize the energy of the Hamiltonian on the\nsystem as a whole, we can minimize the contribution of each group\noperator. Due to the negative coefficients in the Hamiltonian, we need\nto maximize the expectation value of each operator. The maximum possible\nexpectation value for each operator is $+1$. We can turn this into a\nconstraint on our ground state:\n\n$$S_s |G \\rangle = +1 |G \\rangle \\qquad \\qquad P_p | G \\rangle = +1 |G\\rangle.$$\n\nThe wavefunction\n\n$$| G \\rangle = \\prod_{p} \\frac{\\mathbb{I} + P_p}{\\sqrt{2}} |00\\dots 0\\rangle = \\prod_{p} U_p |00\\dots 0 \\rangle,$$\n\nwhere $P_p$ (plaquette) denotes an X Group operator, is such a state.\n\n::: {.note}\n::: {.title}\nNote\n:::\n\nTo check your understanding, confirm that this ground state obeys the\nconstraints using pen and paper.\n:::\n\n$|G \\rangle$ contains a product of unitaries $U_p$. If we can figure out\nhow to apply a single $U_p$ using a quantum computer\\'s operations, we\ncan apply that decomposition for every $p$ in the product.\n\nTo better understand how to decompose $U_p$, let's write it concretely\nfor a single group of four qubits:\n\n$$U |0000 \\rangle =\n\\frac{\\left(\\mathbb{I} + X_1 X_2 X_3 X_4 \\right)}{\\sqrt{2}} |0000 \\rangle\n= \\frac{1}{\\sqrt{2}} \\left( |0000\\rangle + |1111\\rangle \\right).$$\n\nThis [generalized GHZ\nstate](https://en.wikipedia.org/wiki/Greenberger\u2013Horne\u2013Zeilinger_state)\ncan be prepared with a Hadamard and 3 CNOT gates:\n\n![](../demonstrations/toric_code/generalized_ghz_draw.png){.align-center\nwidth=\"50.0%\"}\n\nThis decomposition for $U_p$ holds only when the initial Hadamard qubit\nbegins in the $|0\\rangle$ state, so we need to be careful in the choice\nof the Hadamard\\'s qubit. This restriction is why we rotated the order\nfor a single X Group on the right border earlier.\n\nWe will also not need to prepare the final X Group that contains the\nfour edges of the lattice.\n\nNow let's actually put these together into a circuit!\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"dev = qml.device(\"lightning.qubit\", wires=[Wire(*s) for s in all_sites])\n\ndef state_prep():\n for op in xgroup_ops[0:-1]:\n qml.Hadamard(op.wires[0])\n for w in op.wires[1:]:\n qml.CNOT(wires=[op.wires[0], w])\n\n@qml.qnode(dev, diff_method=None)\ndef circuit():\n state_prep()\n return [qml.expval(op) for op in xgroup_ops + zgroup_ops]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"From this QNode, we can calculate the group operators\\' expectation\nvalues and the system\\'s total energy.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"n_xgroups = len(xgroup_ops)\nseparate_expvals = lambda expvals: (expvals[:n_xgroups], expvals[n_xgroups:])\n\nxgroup_expvals, zgroup_expvals = separate_expvals(circuit())\n\nE0 = -sum(xgroup_expvals) - sum(zgroup_expvals)\n\nprint(\"X Group expectation values\", xgroup_expvals)\nprint(\"Z Group expectation values\", zgroup_expvals)\nprint(\"Total energy: \", E0)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Excitations\n===========\n\nQuasiparticles allow physicists to describe complex systems as\ninteracting particles in a vacuum. Examples of quasiparticles include\nelectrons and holes in semiconductors, phonons, and magnons.\n\nImagine trying to describe the traffic on the road. We could either\nexplicitly enumerate the location of each vehicle or list the positions\nand severities of traffic jams.\n\nThe first option provides complete information about the system but is\nmuch more challenging to analyze. For most purposes, we can work with\ninformation about how the traffic deviates from a baseline. In\nsemiconductors, we don't write out the wave function for every single\nelectron. We instead use electrons and holes. Neither quasiparticle\nelectrons nor holes are fundamental particles like an electron or\npositron in a vacuum. Instead, they are useful descriptions of how the\nwave function differs from its ground state.\n\nWhile the electrons and holes of a metal behave just like electrons and\npositrons in a vacuum, some condensed matter systems contain\nquasiparticles that cannot or do not exist as fundamental particles. The\nexcitations of the toric code are one such example. To find these\nquasiparticles, we look at states that are *almost* the ground state,\nsuch as the ground state with a single operator applied to it.\n\nSuppose we apply a perturbation to the ground state in the form of a\nsingle X gate at location $i$ :\n\n$$| \\phi \\rangle = X_i | G \\rangle.$$\n\nTwo Z group operators $S_s$ contain individual Z operators at that same\nsite $i$. The noise term $X_i$ will anti-commute with both of these\ngroup operators:\n\n$$S_s X_i = \\left( Z_i Z_a Z_b Z_c \\right) X_i = - X_i S_s.$$\n\nUsing this relation, we can determine the expectation value of the Z\ngroup operators with the altered state:\n\n$$S_s |\\phi\\rangle = S_s X_i |G\\rangle = - X_i S_s |G\\rangle = - X_i |G\\rangle = - |\\phi\\rangle.$$\n\nThus,\n\n$$\\langle \\phi | S_s | \\phi \\rangle = -1.$$\n\n$S_s$ now has an expectation value of $-1$.\n\nApplying a single X operator noise term changes the expectation value of\n*two* Z group operators.\n\nThis analysis repeats for the effect of a Z operator on the X Group\nmeasurement. A single Z operator noise term changes the expectation\nvalues of *two* X group operators.\n\nEach group with a negative expectation value is considered an\nexcitation. In the literature, you will often see a Z Group excitation\n$\\langle S_s \\rangle = -1$ called an \"electric\" $e$ excitation and an X\nGroup excitation $\\langle P_p \\rangle = -1$ called a \"magnetic\" $m$\nexcitation. You may also see the inclusion of an identity $\\mathbb{I}$\nparticle for the ground state and the combination particle $\\Psi$\nconsisting of a single $e$ and a single $m$ excitation.\n\nLet's create a QNode where we can apply these perturbations:\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"@qml.qnode(dev, diff_method=None)\ndef excitations(x_sites, z_sites):\n state_prep()\n\n for s in x_sites:\n qml.PauliX(Wire(*s))\n\n for s in z_sites:\n qml.PauliZ(Wire(*s))\n\n return [qml.expval(op) for op in xgroup_ops + zgroup_ops]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"What are the expectation values when we apply a single X operation? Two\nZ group measurements have indeed changed signs.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"single_x = [(1, 2)]\n\nx_expvals, z_expvals = separate_expvals(excitations(single_x, []))\n\nprint(\"XGroup: \", x_expvals)\nprint(\"ZGroup: \", z_expvals)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Instead of interpreting the state via the expectation values of the\noperators, we can view the state as occupation numbers of the\ncorresponding quasiparticles. A group with an expectation value of $+1$\nis in the ground state and thus has an occupation number of $0$. If the\nexpectation value is $-1$, then a quasiparticle exists in that location.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"occupation_numbers = lambda expvals: 0.5 * (1 - expvals)\n\ndef print_info(x_expvals, z_expvals):\n E = -sum(x_expvals) - sum(z_expvals)\n\n print(\"Total energy: \", E)\n print(\"Energy above the ground state: \", E - E0)\n print(\"X Group occupation numbers: \", occupation_numbers(x_expvals))\n print(\"Z Group occupation numbers: \", occupation_numbers(z_expvals))\n\nprint_info(x_expvals, z_expvals)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Since we will plot the same thing many times, we group the following\ncode into a function to easily call later.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"def excitation_plot(x_excite, z_excite):\n x_color = lambda expval: \"navy\" if expval < 0 else \"lavender\"\n z_color = lambda expval: \"maroon\" if expval < 0 else \"mistyrose\"\n\n fig, ax = plt.subplots()\n fig, ax = misc_plot_formatting(fig, ax)\n\n for expval, sites in zip(x_excite, xgroup_sites):\n ax.add_patch(Polygon(sites, color=x_color(expval), zorder=0))\n\n for expval, sites in zip(z_excite, zgroup_sites):\n ax.add_patch(Polygon(sites, color=z_color(expval), zorder=0))\n\n handles = [\n Patch(color=\"navy\", label=\"X Group -1\"),\n Patch(color=\"lavender\", label=\"X Group +1\"),\n Patch(color=\"maroon\", label=\"Z Group -1\"),\n Patch(color=\"mistyrose\", label=\"Z Group +1\"),\n Patch(color=\"navy\", label=\"Z op\"),\n Patch(color=\"maroon\", label=\"X op\"),\n ]\n\n plt.legend(handles=handles, ncol=3, loc=\"lower left\")\n\n return fig, ax"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"With this function, we can quickly view the expectation values and the\nlocation of the additional X operation. The operation changes the\nexpectation values for the adjacent Z groups.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"fig, ax = excitation_plot(x_expvals, z_expvals)\nax.scatter(*zip(*single_x), color=\"maroon\", s=100)\nplt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"What if we apply a Z operation instead at the same site? We instead get\ntwo X Group excitations.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"single_z = [(1, 2)]\n\nexpvals = excitations([], single_z)\nx_expvals, z_expvals = separate_expvals(expvals)\nprint_info(x_expvals, z_expvals)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"fig, ax = excitation_plot(x_expvals, z_expvals)\nax.scatter(*zip(*single_z), color=\"navy\", s=100)\nplt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"What happens if we apply the same perturbation twice at one location? We\nregain the ground state.\n\nThe excitations of the toric code are Majorana particles: particles that\nare their own antiparticles. While postulated to exist in standard\nparticle physics, Majorana particles have only been experimentally seen\nas quasiparticle excitations in materials.\n\nWe can think of the second operation as creating another set of\nexcitations at the same location. This second pair then annihilates the\nexisting particles.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"single_z = [(1, 2)]\n\nexpvals = excitations([], single_z + single_z)\nx_expvals, z_expvals = separate_expvals(expvals)\n\nprint_info(x_expvals, z_expvals)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Moving Excitations and String Operators\n=======================================\n\nWhat if we create a second set of particles such that one of the new\nparticles overlaps with an existing particle? One old and one new\nparticle annihilate each other. One old and one new particle remain. We\nstill have two particles in total.\n\nAlternatively, we can view the second operation as moving a single\nexcitation.\n\nLet's see what that looks like in code:\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"two_z = [(1, 2), (2, 2)]\n\nexpvals = excitations([], two_z)\nx_expvals, z_expvals = separate_expvals(expvals)\n\nprint_info(x_expvals, z_expvals)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"fig, ax = excitation_plot(x_expvals, z_expvals)\nax.plot(*zip(*two_z), color=\"navy\", linewidth=10)\nplt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In that example, we just moved an excitation a little. How about we try\nmoving it even further?\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"long_string = [(1, 2), (2, 2), (3, 2), (4, 1)]\n\nexpvals = excitations([], long_string)\nx_expvals, z_expvals = separate_expvals(expvals)\n\nprint_info(x_expvals, z_expvals)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"fig, ax = excitation_plot(x_expvals, z_expvals)\nax.plot(*zip(*long_string), color=\"navy\", linewidth=10)\nplt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We end up with strings of operations that connect pairs of\nquasiparticles.\n\nWe can use a branch of topology called\n[Homotopy](https://en.wikipedia.org/wiki/Homotopy) to describe the\nrelationship between these strings and the wave function. Two paths\n$s_1$ and $s_2$ are **homotopy equivalent** or **homotopic** if they can\nbe continuously deformed into each other:\n\n$$s_1 \\sim s_2$$\n\nIn the following picture, assume the red \"X\" is a defect in space, like\na tear in a sheet. The two blue paths are equivalent because you can\nsmoothly move one into the other. You cannot move the blue path into the\ngreen path without going through the defect, so they are not equivalent\nto each other.\n\n![Paths are homotopy equivalent if they can be smoothly deformed into\neach other.](../demonstrations/toric_code/homotopy.png){.align-center\nwidth=\"40.0%\"}\n\nWe can divide the set of all possible paths into **homotopy classes**. A\nhomotopy class is an [equivalence\nclass](https://en.wikipedia.org/wiki/Equivalence_class) under homotopy.\nEvery member of the same homotopy class can be deformed into every other\nmember of the same class, and members of different homotopy classes\ncannot be deformed into each other. All the homotopy classes for a given\nspace $S$ form its [first homotopy\ngroup](https://en.wikipedia.org/wiki/Homotopy_group), denoted by\n$\\pi_1(S)$. The first homotopy group of a space is also called its\n*fundamental group*.\n\nHow do these mathematical concepts apply to the toric code model? To\nfind out, let\\'s look at a string of Z operations homotopic to the one\nabove.\n\nThis next string creates the same final state. Only the homotopy class\nof the path used to create the excitations influences the occupation\nnumbers and total energy. If the endpoints are the same and the path\ndoesn't wrap around the torus or other particles, the details do not\nimpact any observables.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"equivalent_string = [(1, 2), (2, 1), (3, 1), (4, 1)]\n\nexpvals = excitations([], equivalent_string)\nx_expvals, z_expvals = separate_expvals(expvals)\n\nprint_info(x_expvals, z_expvals)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"fig, ax = excitation_plot(x_expvals, z_expvals)\nax.plot(*zip(*equivalent_string), color=\"navy\", linewidth=10)\nplt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Contractible loops\n==================\n\nWe can also have a loop of operations that doesn't create any new\nexcitations. The loop forms a pair, moves one around in a circle, and\nthen annihilates the two particles again.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"contractible_loop = [(1, 1), (2, 1), (3, 1), (4, 1), (4, 2), (3, 3), (2, 3), (1, 2)]\n\nexpvals = excitations(contractible_loop, [])\nx_expvals, z_expvals = separate_expvals(expvals)\nprint_info(x_expvals, z_expvals)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"fig, ax = excitation_plot(x_expvals, z_expvals)\n\nax.plot(*zip(*contractible_loop), color=\"maroon\", linewidth=10)\n\nplt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The loop doesn't affect the positions of any excitations, but does it\naffect the state at all?\n\nWe will look at the probabilities instead of the expectation values to\nanswer that question.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"@qml.qnode(dev, diff_method=None)\ndef probs(x_sites, z_sites):\n state_prep()\n\n for s in x_sites:\n qml.PauliX(Wire(*s))\n\n for s in z_sites:\n qml.PauliZ(Wire(*s))\n\n return qml.probs(wires=[Wire(*s) for s in all_sites])\n\n\nnull_probs = probs([], [])\ncontractible_probs = probs(contractible_loop, [])\n\nprint(\"Are the probabilities equal? \", np.allclose(null_probs, contractible_probs))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The toric code\\'s dependence on the homotopy of the path explains this\nresult. All paths we can smoothly deform into each other will give the\nsame result. The contractible loop can be smoothly deformed to nothing,\nso the state with the contractible loop is the same as the state with no\nloop, our initial $|G\\rangle$.\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Looping the torus\n=================\n\nOn the torus, we have four types of unique paths:\n\n- The trivial path that contracts to nothing\n- A horizontal loop around the boundaries\n- A vertical loop around the boundaries\n- A loop around both the horizontal and vertical boundaries\n\n![The homotopy group of a torus can be generated by a horizontal loop\nand a vertical loop around the\nboundaries.](../demonstrations/toric_code/types_of_loops.png){.align-center\nwidth=\"50.0%\"}\n\nEach of these paths represents a member of the first homotopy group of\nthe torus: $\\pi_1(T) = \\mathbb{Z}^2$ modulo 2.\n\nNone of these loops of X operations create net excitations, so the wave\nfunction remains in the ground state.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"horizontal_loop = [(i, 1) for i in range(width)]\nvertical_loop = [(1, i) for i in range(height)]\n\nexpvals = excitations(horizontal_loop + vertical_loop, [])\nfig, ax = excitation_plot(*separate_expvals(expvals))\n\nax.plot(*zip(*horizontal_loop), color=\"maroon\", linewidth=10)\nax.plot(*zip(*vertical_loop), color=\"maroon\", linewidth=10)\n\nplt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We can compute the probabilities for each of these four types of loops:\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"null_probs = probs([], [])\nhorizontal_probs = probs(horizontal_loop, [])\nvertical_probs = probs(vertical_loop, [])\ncombo_probs = probs(horizontal_loop + vertical_loop, [])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"While X and Z operations can change the group operator expectation\nvalues and create quasiparticles, only X operators can change the\nprobability distribution. Applying a Z operator would only rotate the\nphase of the state and not change any amplitudes. Hence we only use\nloops of X operators in this section. I encourage you to try this\nanalysis with loops of Z operators to confirm that they do not change\nthe probability distribution.\n\nWe can compare the original state and one with a horizontal loop to see\nif the probability distributions are different:\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"print(\"Are the probabilities equal? \", qml.math.allclose(null_probs, horizontal_probs))\nprint(\"Is this significant?\")\nprint(\"Maximum difference in probabilities: \", max(abs(null_probs - horizontal_probs)))\nprint(\"Maximum probability: \", max(null_probs))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"These numbers seem small, but remember we have 24 qubits, and thus\n$2^{24}=16777216$ probability components. Since the maximum difference\nin probabilities is the same size as the maximum probability, we know\nthis isn't just random fluctuations and errors. The states are\ndifferent, but the energies are the same, so they are *degenerate*.\n\nThat compared a horizontal \"x\" loop with the initial ground state. How\nabout the other two types of loops? Let's iterate over all combinations\nof two probability distributions to see if any match.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"probs_type_labels = [\"null\", \"x\", \"y\", \"combo\"]\nall_probs = [null_probs, horizontal_probs, vertical_probs, combo_probs]\n\nprint(\"\\t\" + \"\\t\".join(probs_type_labels))\n\nfor name, probs1 in zip(probs_type_labels, all_probs):\n comparisons = (str(np.allclose(probs1, probs2)) for probs2 in all_probs)\n print(name, \"\\t\", \"\\t\".join(comparisons))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This table shows the model has four distinct ground states. More\nimportantly, these ground states are separated from each other by\nlong-range operations. We must perform a loop of operations across the\nentire lattice to switch between degenerate ground states.\n\nThis four-way degeneracy is the source of the error correction in the\ntoric code. Instead of 24 qubits, we work with two logical qubits (4\nstates) that are well separated from each other by topological\noperations.\n\n::: {.note}\n::: {.title}\nNote\n:::\n\nI encourage dedicated readers to explore what happens when a path loops\nthe same boundaries twice.\n:::\n\nIn this section, we\\'ve seen that the space of ground states is directly\nrelated to the first homotopy group of the lattice. The first homotopy\ngroup of a torus is $\\pi_1(T) = \\mathbb{Z}^2$, and the space of ground\nstates is that group modulo two, $\\mathbb{Z}_2^2$.\n\nWhat if we defined the model on a differently shaped lattice? Then the\nspace of the ground state would change to reflect the first homotopy\ngroup of that space. For example, if the model was defined on a sphere,\nthen only a single unique ground state would exist. Adding defects like\nmissing sites to the lattice also changes the topology. Error correction\nwith the toric code often uses missing sites to add additional logical\nqubits.\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Mutual and Exchange Statistics\n==============================\n\nThe hole in the center of the donut isn't the only thing that prevents\npaths from smoothly deforming into each other. We don't yet know if we\ncan distort paths past other particles.\n\nWhen one indistinguishable fermion of spin 1/2 orbits another fermion of\nthe same type, the combined wave function picks up a relative phase of\nnegative one. When fermions of different types orbit each other, the\nstate is unchanged. For example, if an electron goes around a proton and\nreturns to the same spot, the wave function is unchanged. If a boson\norbits around a different type of boson, again, the wave function is\nunchanged.\n\nWhat if a particle went around a different type of particle and\neverything picked up a phase? Would it be a boson or a fermion?\n\nIt would be something else entirely: an anyon. An anyon is anything that\ndoesn't cleanly fall into the boson/fermion categorization of particles.\n\nWhile the toric code is just an extremely useful mathematical model,\nanyons exist in physical materials. For example, fractional quantum Hall\nsystems have anyonic particles with spin $1/q$ for different integers\n$q$.\n\nThe **statistics** are described by the phase accumulated by moving one\nparticle around another. For example, if the particle picks up phases\nlike a fermion, then it obeys [Fermi-Dirac\nstatistics](https://en.wikipedia.org/wiki/Fermi\u2013Dirac_statistics).\n**Exchange statistics** are described by the phases that accumulate from\nexchanging the *same* type of particles. **Mutual** statistics are\ncharacterized by the phase acquired by moving one particle around a\nparticle of a *different* type.\n\nTo measure the mutual statistics of a Z Group excitation and an X group\nexcitation, we need to prepare at least one of each type of particle and\nthen orbit one around the other.\n\nThe following code rotates an X Group excitation around a Z Group\nexcitation.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"prep1 = [(1, 1), (2, 1)]\nprep2 = [(1, 3)]\nloop1 = [(2, 3), (2, 2), (2, 1), (3, 1), (3, 2), (2, 3)]\n\nexpvals = excitations(prep1, prep2 + loop1)\nx_expvals, z_expvals = separate_expvals(expvals)\n\nfig, ax = excitation_plot(x_expvals, z_expvals)\n\nax.plot(*zip(*prep1), color=\"maroon\", linewidth=10)\nax.plot(*zip(*(prep2 + loop1)), color=\"navy\", linewidth=10)\n\nplt.show()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"While we managed to loop one particle around the other, we did not\nextract the relative phase applied to the wave function. To procure this\ninformation, we will need the *Hadamard test*.\n"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Hadamard test\n=============\n\nThe [Hadamard\ntest](https://en.wikipedia.org/wiki/Hadamard_test_(quantum_computation))\nextracts the real component of a unitary operation\n$\\text{Re}\\left(\\langle \\psi | U | \\psi \\rangle \\right)$. If the unitary\noperation just applies a phase\n$U |\\psi\\rangle = e^{i \\phi} |\\psi \\rangle$, the measured quantity\nreduces to $\\cos (\\phi)$.\n\nThe steps in the Hadamard test are:\n\n1. Prepare the auxiliary qubit into a superposition with a Hadamard\n gate\n2. Apply a controlled version of the operation with the auxiliary qubit\n as the control\n3. Apply another Hadamard gate to the auxiliary qubit\n4. Measure the auxiliary qubit in the Z-basis\n\n![](../demonstrations/toric_code/Hadamard_test.png){.align-center\nwidth=\"50.0%\"}\n\n::: {.note}\n::: {.title}\nNote\n:::\n\nFor extra understanding, validate the Hadamard test algorithm using pen\nand paper.\n:::\n\nBelow we implement this algorithm in PennyLane and measure the mutual\nexchange statistics of an X Group excitation and a Z Group excitation.\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"dev_aux = qml.device(\"lightning.qubit\", wires=[Wire(*s) for s in all_sites] + [\"aux\"])\n\ndef loop(x_loop, z_loop):\n for s in x_loop:\n qml.PauliX(Wire(*s))\n for s in z_loop:\n qml.PauliZ(Wire(*s))\n\n@qml.qnode(dev_aux, diff_method=None)\ndef hadamard_test(x_prep, z_prep, x_loop, z_loop):\n state_prep()\n\n for s in x_prep:\n qml.PauliX(Wire(*s))\n\n for s in z_prep:\n qml.PauliZ(Wire(*s))\n\n qml.Hadamard(\"aux\")\n qml.ctrl(loop, control=\"aux\")(x_loop, z_loop)\n qml.Hadamard(\"aux\")\n return qml.expval(qml.PauliZ(\"aux\"))\n\nx_around_z = hadamard_test(prep1, prep2, [], loop1)\nprint(\"Move x excitation around z excitation: \", x_around_z)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We just moved two different types of particles around each other and\npicked up a phase. As neither bosons nor fermions behave like this, this\nresult demonstrates that the excitations of a toric code are anyons.\n\n::: {.note}\n::: {.title}\nNote\n:::\n\nI encourage dedicated readers to calculate the phase accumulated by\nexchanging:\n\n- A Z Group excitation and a Z Group Excitation\n- An X Group excitation and an X Group Excitation\n- A combination $\\Psi$ particle and an X Group excitation\n- A combination $\\Psi$ particle and a Z Group excitation\n- A combination $\\Psi$ particle with another $\\Psi$ particle\n\nThe combination particle should behave like a standard fermion. You can\ncreate and move combination particles by applying `PauliY` operations.\n:::\n\nIn this demo, we have demonstrated:\n\n1. How to prepare the ground state of the toric code model on a lattice\n of qubits\n2. How to create and move excitations\n3. The ground state degeneracy of the model on a toric lattice, arising\n from homotopically distinct loops of operations\n4. The excitations are anyons due to non-trivial mutual statistics Make\n sure to go and modify the code as suggested if you wish to gain more\n intuition.\n\nDo check the references below if you want to learn more!\n\nReferences\n==========\n\nAbout the author\n================\n"
]
}
],
"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.9.15"
}
},
"nbformat": 4,
"nbformat_minor": 0
}