diff --git a/docs/demos/global_rotations.ipynb b/docs/demos/global_rotations.ipynb new file mode 100644 index 00000000..957dae87 --- /dev/null +++ b/docs/demos/global_rotations.ipynb @@ -0,0 +1,746 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "2c84924c", + "metadata": {}, + "source": [ + "# Subjecting codes of various dimensionality to global rotations\n", + "\n", + "This tutorial reproduces **Fig. 4a** of [Bluvstein *et al.*, *Nature* **649**, 39 (2026)](https://www.nature.com/articles/s41586-025-09848-5), \"A fault-tolerant neutral-atom architecture for universal quantum computation\", using **Tsim**.\n", + "\n", + "The experiment is simple to state. Take a logical qubit encoded in a quantum error-correcting code, prepare it in the logical $|+_L\\rangle$ state, apply the *same* single-qubit $Z$-rotation $R_Z(\\varphi)$ to every physical qubit at once, and read out the logical $X$ operator and the code stabilizers as a function of $\\varphi$.\n", + "\n", + "A global physical rotation acts as a **transversal logical gate**, and which logical gate you get depends on the *geometry* of the code. We will see that:\n", + "\n", + "- a single unencoded qubit just traces out a smooth cosine,\n", + "- the **2D colour code** ([[7,1,3]] Steane code) shows flat *plateaus* at multiples of $90^\\circ$ - the logical Clifford angles $S, S^\\dagger, Z$,\n", + "- the **3D colour code** ([[15,1,3]] quantum Reed-Muller code) shows an *extra* plateau at multiples of $45^\\circ$ - the angles $T, T^\\dagger, TS, (TS)^\\dagger$ of its transversal **non-Clifford $T$ gate**.\n", + "\n", + "That $45^\\circ$ plateau is the whole point. It is the gate that lifts a stabilizer code out of the Clifford group and makes universal, fault-tolerant computation possible. Tsim is a natural tool for this: a global $R_Z(\\varphi)$ on the 15-qubit code is 15 non-Clifford rotations at once, which Tsim handles with ZX-calculus stabilizer-rank decomposition rather than a dense $2^{15}$ statevector." + ] + }, + { + "cell_type": "markdown", + "id": "b6366e2c", + "metadata": {}, + "source": [ + "## The physics\n", + "\n", + "### Global rotations as transversal gates\n", + "\n", + "The transversal operation here is a global $Z$-rotation,\n", + "\n", + "$$U(\\varphi) = \\bigotimes_{j} R_Z(\\varphi), \\qquad R_Z(\\varphi) = e^{-i\\varphi Z/2}.$$\n", + "\n", + "For a CSS code this commutes with every $Z$-type stabilizer, so the $Z$ checks are untouched; it does *not* commute with the $X$-type stabilizers or with $X_L$, and that is what we probe.\n", + "\n", + "At special angles $U(\\varphi)$ maps the code space onto itself and realizes a logical gate. At $\\varphi = \\pi/2$ this is always a logical Clifford. Whether $\\varphi = \\pi/4$ realizes a logical **$T$** depends on the code: it works exactly for codes whose $X$-stabilizers are *triorthogonal*. The [[15,1,3]] Reed-Muller code is built that way (it has a transversal $T$); the [[7,1,3]] Steane code is not (its $X$-checks have weight 4, only doubly-even), so it has no protected $45^\\circ$ point.\n", + "\n", + "The **Eastin-Knill theorem** forbids a transversal gate set that is both universal and unitary, which is why no single code plateaus everywhere. The paper circumvents this with logical measurement: a transversal CZ between a Reed-Muller $|T_L\\rangle$ and $|+_L\\rangle$ followed by measurement teleports an $H$, giving a universal $\\{H, T, \\text{CNOT}\\}$ set from fully transversal operations.\n", + "\n", + "### Why plateaus appear, and why stabilizer *signs* matter\n", + "\n", + "Both codes have distance 3, so they correct any single-qubit error. Near an angle where $U(\\varphi)$ is a clean logical gate, a small deviation $\\delta\\varphi$ looks like a low-weight coherent error the code absorbs, so the logical expectation stays *flat* - these are the plateaus - and the stabilizers *revive* toward $\\pm 1$. A key subtlety from the paper: for the non-Clifford $T$ point, the stabilizers must have the **correct ($+1$) signs**. Flipping signs (here, the four corners of the Reed-Muller tetrahedron) destroys the $45^\\circ$ plateau - the \"3D negative stabilizer\" control. This is the statement that **logical magic requires physical entanglement**: $|+_L\\rangle$ is an eigenstate of a tensor-product operator $X_L = X_1X_2X_3\\dots$, but $|T_L\\rangle$ is an eigenstate of $\\tfrac{1}{\\sqrt2}(X_1X_2X_3\\dots + Y_1Y_2Y_3\\dots)$, a superposition spanning the code that is necessarily entangled.\n", + "\n", + "### What Tsim computes\n", + "\n", + "For each $\\varphi$ Tsim builds the circuit as a ZX diagram, decomposes the non-Clifford rotations by stabilizer rank, and **samples** the readout. We use Tsim's sampler for the headline result (including the 15-rotation 3D code), and cross-check against an exact statevector reference computed from the Clifford encoder plus an analytic diagonal rotation." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "f7e6f001", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:08:57.098379Z", + "iopub.status.busy": "2026-06-10T23:08:57.098204Z", + "iopub.status.idle": "2026-06-10T23:08:58.567376Z", + "shell.execute_reply": "2026-06-10T23:08:58.566761Z" + } + }, + "outputs": [], + "source": [ + "import numpy as np\n", + "import stim\n", + "import matplotlib.pyplot as plt\n", + "\n", + "import tsim" + ] + }, + { + "cell_type": "markdown", + "id": "cd9f1768", + "metadata": {}, + "source": [ + "## Building the codes\n", + "\n", + "The [[7,1,3]], [[15,1,3]] and [[16,6,4]] codes in the paper are all members of the **quantum Reed-Muller family** built from the hypercube encoding circuit of Extended Data Fig. 10a (ref. 128). We will build that literal circuit below, but first we define the codes by their stabilizer groups - which is convenient and easy to verify:\n", + "\n", + "- **[[7,1,3]] Steane / 2D colour code** - the self-dual CSS code of the $[7,4,3]$ Hamming code.\n", + "- **[[15,1,3]] Reed-Muller / 3D colour code** - $X$-checks from the punctured $\\mathrm{RM}(1,4)$ code (weight-8), $Z$-checks from $\\mathrm{RM}(2,4)$ (weight-8 and weight-4). This triorthogonal structure is exactly what gives the transversal $T$.\n", + "\n", + "We ask `stim` for the Clifford that maps $|0\\dots0\\rangle$ onto the state stabilized by the generators together with $X_L$ (so the prepared state is the $+1$ eigenstate of $X_L$, i.e. $|+_L\\rangle$), and verify it with `stim.Circuit.flow_generators()` as the issue suggests." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "9f22dd96", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:08:58.569135Z", + "iopub.status.busy": "2026-06-10T23:08:58.568929Z", + "iopub.status.idle": "2026-06-10T23:08:58.573521Z", + "shell.execute_reply": "2026-06-10T23:08:58.572979Z" + } + }, + "outputs": [], + "source": [ + "def _pauli(kind, support, n):\n", + " \"\"\"Build a stim Pauli-string label, e.g. ('X', [0, 2], 4) -> 'X_X_'.\"\"\"\n", + " chars = [\"_\"] * n\n", + " for q in support:\n", + " chars[q] = kind\n", + " return \"\".join(chars)\n", + "\n", + "\n", + "def plus_logical_encoder(n, x_stabs, z_stabs, x_logical, corner_flips=None):\n", + " \"\"\"Clifford circuit preparing |+_L> from |0...0>.\n", + "\n", + " The prepared state is the simultaneous +1 eigenstate of every stabilizer\n", + " and of X_L. ``corner_flips`` appends a final X(pi) pulse to the listed\n", + " qubits: flipping the four corners of the Reed-Muller tetrahedron gives the\n", + " 'incorrect stabilizer signs' of the paper's '3D negative stabilizer'\n", + " control, which destroys the transversal-T plateau at 45 degrees.\n", + " \"\"\"\n", + " generators = [stim.PauliString(_pauli(\"X\", s, n)) for s in x_stabs]\n", + " generators += [stim.PauliString(_pauli(\"Z\", s, n)) for s in z_stabs]\n", + " generators.append(stim.PauliString(_pauli(\"X\", x_logical, n)))\n", + "\n", + " tableau = stim.Tableau.from_stabilizers(\n", + " generators, allow_redundant=False, allow_underconstrained=False\n", + " )\n", + " circuit = tableau.to_circuit(\"elimination\")\n", + " for q in corner_flips or []:\n", + " circuit.append(\"X\", [q])\n", + " return circuit" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "de269c13", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:08:58.574921Z", + "iopub.status.busy": "2026-06-10T23:08:58.574783Z", + "iopub.status.idle": "2026-06-10T23:08:58.579944Z", + "shell.execute_reply": "2026-06-10T23:08:58.579040Z" + } + }, + "outputs": [], + "source": [ + "# --- [[7,1,3]] Steane code (2D colour code) ----------------------------------\n", + "# Self-dual: X and Z checks share the [7,4,3] Hamming support.\n", + "steane = dict(\n", + " n=7,\n", + " x_stabs=[[3, 4, 5, 6], [1, 2, 5, 6], [0, 2, 4, 6]],\n", + " z_stabs=[[3, 4, 5, 6], [1, 2, 5, 6], [0, 2, 4, 6]],\n", + " x_logical=list(range(7)), # transversal all-X representative of X_L\n", + ")\n", + "\n", + "# --- [[15,1,3]] quantum Reed-Muller code (3D colour code) --------------------\n", + "# Qubit j (0-indexed) sits at the 4-bit point (j + 1) of the punctured tesseract.\n", + "# X-checks = punctured RM(1,4) (the four coordinate functions, weight 8).\n", + "# Z-checks = RM(2,4) (those four, plus the six degree-2 monomials, weight 4).\n", + "def _bit(j, k):\n", + " return (j >> k) & 1\n", + "\n", + "\n", + "_rm_x = [[j for j in range(15) if _bit(j + 1, k)] for k in range(4)]\n", + "_rm_z = [list(s) for s in _rm_x] + [\n", + " [j for j in range(15) if _bit(j + 1, a) and _bit(j + 1, b)]\n", + " for a in range(4)\n", + " for b in range(a + 1, 4)\n", + "]\n", + "# The four corners of the tetrahedron are the weight-1 points e_k (1,2,4,8),\n", + "# i.e. qubits 0, 1, 3, 7. Flipping them sets incorrect stabilizer signs.\n", + "reed_muller = dict(\n", + " n=15, x_stabs=_rm_x, z_stabs=_rm_z, x_logical=list(range(15)), corners=[0, 1, 3, 7]\n", + ")\n", + "\n", + "# --- Unentangled reference: 15 bare qubits \"analysed as a 3D code\" -----------\n", + "# As in the paper, 15 unentangled physical qubits in |+>^15 are read out with\n", + "# the *same* Reed-Muller operators (X_L and the weight-8 X-checks), so the\n", + "# only difference from the 3D curve is the missing entanglement.\n", + "unentangled = dict(\n", + " n=15, x_stabs=_rm_x, z_stabs=[], x_logical=list(range(15)), product=True\n", + ")\n", + "\n", + "print(\"Steane X-stab weights :\", [len(s) for s in steane[\"x_stabs\"]])\n", + "print(\"RM X-stab weights :\", [len(s) for s in reed_muller[\"x_stabs\"]])\n", + "print(\"RM Z-stab weights :\", [len(s) for s in reed_muller[\"z_stabs\"]])" + ] + }, + { + "cell_type": "markdown", + "id": "abfdb867", + "metadata": {}, + "source": [ + "`flow_generators()` reports how the Clifford encoder propagates the input Paulis. Below, the input $Z_{14}$ flows to the all-$X$ logical $X_L$ and the inputs $X_i$ flow to the code stabilizers - confirming we prepared the intended $|+_L\\rangle$ of the [[15,1,3]] code." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "f9ffd2a3", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:08:58.581446Z", + "iopub.status.busy": "2026-06-10T23:08:58.581212Z", + "iopub.status.idle": "2026-06-10T23:08:58.585815Z", + "shell.execute_reply": "2026-06-10T23:08:58.585103Z" + } + }, + "outputs": [], + "source": [ + "encoder = plus_logical_encoder(\n", + " reed_muller[\"n\"], reed_muller[\"x_stabs\"], reed_muller[\"z_stabs\"], reed_muller[\"x_logical\"]\n", + ")\n", + "print(f\"[[15,1,3]] encoder: {len(encoder)} gates, {len(encoder.flow_generators())} flow generators\")\n", + "for flow in encoder.flow_generators()[:6]:\n", + " print(\" \", flow)\n", + "print(\" ...\")" + ] + }, + { + "cell_type": "markdown", + "id": "995d6a53", + "metadata": {}, + "source": [ + "## The hypercube encoding circuit (Extended Data Fig. 10a)\n", + "\n", + "The paper prepares these codes with a single *hypercube encoding circuit* (ref. 128): place $2^m$ qubits at the corners of an $m$-cube, apply the transversal CNOT network $E = \\bigl(\\begin{smallmatrix}1&0\\\\1&1\\end{smallmatrix}\\bigr)^{\\otimes m}$ (at step $t$, a CNOT across every edge in direction $t$), with a programmable product input - $|+\\rangle$ on rows of low Hamming weight, $|0\\rangle$ on rows of high weight. A different input pattern (set by local $Y(\\pi/2)$ pulses) gives a different code from the *same* entangling structure:\n", + "\n", + "- **[[16,6,4]]** tesseract - full $m=4$ cube, no puncture;\n", + "- **[[15,1,3]]** - $m=4$ cube punctured at the $0\\dots0$ corner;\n", + "- **[[7,1,3]]** - $m=3$ cube punctured at the $0\\dots0$ corner.\n", + "\n", + "We build this literal circuit and check, by statevector overlap, that it prepares exactly the same $|+_L\\rangle$ as the stabilizer encoder above. (For the rotation sweep we then use the deterministic stabilizer encoder, since it prepares the identical state without the puncturing measurement.)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "a3f86c96", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:08:58.587067Z", + "iopub.status.busy": "2026-06-10T23:08:58.586933Z", + "iopub.status.idle": "2026-06-10T23:08:58.940650Z", + "shell.execute_reply": "2026-06-10T23:08:58.940060Z" + } + }, + "outputs": [], + "source": [ + "def hypercube_plus_logical(m, r_x):\n", + " \"\"\"Prepare |+_L> via the paper's hypercube encoder on an m-cube.\n", + "\n", + " |+> on rows of Hamming weight <= r_x (and on the logical corner 0...0),\n", + " |0> elsewhere, then the transversal CNOT network E; the 0...0 corner is\n", + " punctured by projecting it onto |+> (an X-basis logical measurement).\n", + " Returns the statevector of the remaining 2^m - 1 physical qubits.\n", + " \"\"\"\n", + " n = 1 << m\n", + " circuit = stim.Circuit()\n", + " for q in range(n):\n", + " if q == 0 or bin(q).count(\"1\") <= r_x:\n", + " circuit.append(\"H\", [q]) # |+> input; all others stay |0>\n", + " for t in range(m): # E = ((1,0),(1,1))^{otimes m}: CNOT across edge-direction t\n", + " for q in range(n):\n", + " if _bit(q, t) == 0:\n", + " circuit.append(\"CX\", [q, q | (1 << t)])\n", + " sim = stim.TableauSimulator()\n", + " sim.do_circuit(circuit)\n", + " psi = np.asarray(sim.state_vector(endian=\"big\"), dtype=np.complex128).reshape([2] * n)\n", + " branch0 = psi[(0, *([slice(None)] * (n - 1)))].reshape(-1)\n", + " branch1 = psi[(1, *([slice(None)] * (n - 1)))].reshape(-1)\n", + " sub = branch0 + branch1 # project the punctured corner onto |+>\n", + " return sub / np.linalg.norm(sub)\n", + "\n", + "\n", + "def _stabilizer_statevector(code):\n", + " sim = stim.TableauSimulator()\n", + " sim.do_circuit(plus_logical_encoder(code[\"n\"], code[\"x_stabs\"], code[\"z_stabs\"], code[\"x_logical\"]))\n", + " return np.asarray(sim.state_vector(endian=\"big\"), dtype=np.complex128)\n", + "\n", + "\n", + "overlap_steane = abs(np.vdot(hypercube_plus_logical(3, 1), _stabilizer_statevector(steane)))\n", + "overlap_rm = abs(np.vdot(hypercube_plus_logical(4, 1), _stabilizer_statevector(reed_muller)))\n", + "print(f\"hypercube vs stabilizer encoder, ||:\")\n", + "print(f\" [[7,1,3]] (3-cube): {overlap_steane:.6f}\")\n", + "print(f\" [[15,1,3]] (4-cube): {overlap_rm:.6f}\")\n", + "assert overlap_steane > 1 - 1e-6 and overlap_rm > 1 - 1e-6\n", + "print(\" -> the literal hypercube circuit prepares the same |+_L> state.\")" + ] + }, + { + "cell_type": "markdown", + "id": "5de691c6", + "metadata": {}, + "source": [ + "## The rotation experiment\n", + "\n", + "For a given angle $\\varphi$ the circuit is: encode $|+_L\\rangle$, apply $R_Z(\\varphi)$ to every physical qubit, then read out $X_L$ and the $X$-stabilizers. Tsim takes angles in units of $\\pi$, so $\\varphi$ radians is `R_Z(phi / pi)`. We read out with `MPP` (Pauli-product measurement) on the logical operator and each $X$-stabilizer, tagged with `OBSERVABLE_INCLUDE` and `DETECTOR`; measuring the products directly keeps the sampler's output count small, which matters at 15 qubits.\n", + "\n", + "Tsim also shows *why* the 3D case is the hard one: `tcount()` reports the number of non-Clifford gates." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "376e3f17", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:08:58.946772Z", + "iopub.status.busy": "2026-06-10T23:08:58.946478Z", + "iopub.status.idle": "2026-06-10T23:08:58.965232Z", + "shell.execute_reply": "2026-06-10T23:08:58.964108Z" + } + }, + "outputs": [], + "source": [ + "def rotated_circuit(code, phi, *, corner_flips=None):\n", + " \"\"\"|+_L>, a global R_Z(phi) (phi in radians), then MPP readout of X_L and\n", + " the X-stabilizers, annotated for Tsim's detector sampler.\n", + " \"\"\"\n", + " n = code[\"n\"]\n", + " enc = plus_logical_encoder(\n", + " code[\"n\"], code[\"x_stabs\"], code[\"z_stabs\"], code[\"x_logical\"], corner_flips\n", + " )\n", + " text = str(enc) + \"\\n\"\n", + " text += f\"R_Z({phi / np.pi}) \" + \" \".join(map(str, range(n))) + \"\\n\"\n", + " # rec[-(n_stab + 1)] is X_L; the last n_stab records are the stabilizers.\n", + " text += \"MPP \" + \"*\".join(f\"X{q}\" for q in code[\"x_logical\"]) + \"\\n\"\n", + " for sup in code[\"x_stabs\"]:\n", + " text += \"MPP \" + \"*\".join(f\"X{q}\" for q in sup) + \"\\n\"\n", + " n_stab = len(code[\"x_stabs\"])\n", + " text += f\"OBSERVABLE_INCLUDE(0) rec[{-(n_stab + 1)}]\\n\"\n", + " for i in range(n_stab):\n", + " text += f\"DETECTOR rec[{-(n_stab - i)}]\\n\"\n", + " return tsim.Circuit(text)\n", + "\n", + "\n", + "# At phi = pi/4 the global rotation is a transversal T: one T gate per qubit.\n", + "print(\"non-Clifford gate count (tcount) at phi = pi/4:\")\n", + "print(\" [[7,1,3]] :\", rotated_circuit(steane, np.pi / 4).tcount())\n", + "print(\" [[15,1,3]]:\", rotated_circuit(reed_muller, np.pi / 4).tcount())" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "2e94b5d2", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:08:58.967170Z", + "iopub.status.busy": "2026-06-10T23:08:58.966960Z", + "iopub.status.idle": "2026-06-10T23:08:58.983385Z", + "shell.execute_reply": "2026-06-10T23:08:58.979426Z" + } + }, + "outputs": [], + "source": [ + "# Tsim's ZX timeline view of the Steane circuit at phi = pi/4.\n", + "rotated_circuit(steane, np.pi / 4).diagram(\"timeline-svg\", height=260)" + ] + }, + { + "cell_type": "markdown", + "id": "6cc904e7", + "metadata": {}, + "source": [ + "### Exact reference values\n", + "\n", + "For ground truth we compute the expectation values exactly. The encoder is Clifford, so `stim`'s tableau simulator gives $|+_L\\rangle$ as a statevector; the global rotation is then applied analytically, multiplying the amplitude of a computational basis state $|b\\rangle$ by $e^{i\\varphi\\,\\mathrm{popcount}(b)}$ (up to a global phase). An $X$-type operator flips the bits on its support, so $\\langle X_S\\rangle = \\mathrm{Re}\\,\\langle\\psi|\\psi_{\\text{flip}(S)}\\rangle$.\n", + "\n", + "This is exact and cheap (one $2^n$ vector reused for every angle). It is *not* how Tsim or the paper work at scale - a dense $2^{15}$ contraction of the rotated diagram exhausts memory - but it is the right yardstick for these small codes." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "b52bf992", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:08:58.986554Z", + "iopub.status.busy": "2026-06-10T23:08:58.986340Z", + "iopub.status.idle": "2026-06-10T23:08:59.131746Z", + "shell.execute_reply": "2026-06-10T23:08:59.130104Z" + } + }, + "outputs": [], + "source": [ + "def _support_mask(support, n):\n", + " \"\"\"Bitmask of an operator support; qubit q occupies bit (n - 1 - q).\"\"\"\n", + " mask = 0\n", + " for q in support:\n", + " mask |= 1 << (n - 1 - q)\n", + " return mask\n", + "\n", + "\n", + "def _rotated_statevector(code, phi, corner_flips=None):\n", + " if code.get(\"product\"):\n", + " # Unentangled reference: |+>^n product state, no encoder.\n", + " prep = stim.Circuit(\"H \" + \" \".join(map(str, range(code[\"n\"]))))\n", + " else:\n", + " prep = plus_logical_encoder(\n", + " code[\"n\"], code[\"x_stabs\"], code[\"z_stabs\"], code[\"x_logical\"], corner_flips\n", + " )\n", + " sim = stim.TableauSimulator()\n", + " sim.do_circuit(prep)\n", + " psi0 = np.asarray(sim.state_vector(endian=\"big\"), dtype=np.complex128)\n", + " idx = np.arange(psi0.size)\n", + " popcount = np.array([int(b).bit_count() for b in idx])\n", + " return psi0 * np.exp(1j * phi * popcount), idx\n", + "\n", + "\n", + "def exact_curves(code, phis, *, corner_flips=None):\n", + " \"\"\"Raw (undecoded) (phi) and the mean ||(phi).\"\"\"\n", + " n = code[\"n\"]\n", + " mask_logical = _support_mask(code[\"x_logical\"], n)\n", + " mask_stabs = [_support_mask(s, n) for s in code[\"x_stabs\"]]\n", + " logical, stabilizer = [], []\n", + " for phi in phis:\n", + " psi, idx = _rotated_statevector(code, phi, corner_flips)\n", + " logical.append(float(np.real(np.vdot(psi, psi[idx ^ mask_logical]))))\n", + " if mask_stabs:\n", + " stabilizer.append(\n", + " float(np.mean([abs(np.real(np.vdot(psi, psi[idx ^ m]))) for m in mask_stabs]))\n", + " )\n", + " else:\n", + " stabilizer.append(np.nan)\n", + " return np.array(logical), np.array(stabilizer)\n", + "\n", + "\n", + "# Sanity check: at phi = 0 the state is |+_L>, so every expectation value is +1.\n", + "for name, code in [(\"steane\", steane), (\"reed_muller\", reed_muller)]:\n", + " xl0, s0 = exact_curves(code, [0.0])\n", + " assert np.isclose(xl0[0], 1.0, atol=1e-6), (name, xl0)\n", + " assert np.isclose(s0[0], 1.0, atol=1e-6), (name, s0)\n", + "print(\"phi = 0 self-check passed: = = +1 for both codes.\")" + ] + }, + { + "cell_type": "markdown", + "id": "f6481096", + "metadata": {}, + "source": [ + "### Sampling with Tsim\n", + "\n", + "Now the part that actually needs Tsim: turning a circuit with many non-Clifford rotations into samples. We compile a detector sampler and estimate each expectation value from Monte-Carlo shots via $\\langle P\\rangle = 1 - 2\\,\\Pr(\\text{outcome} = 1)$." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "b8912af6", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:08:59.138705Z", + "iopub.status.busy": "2026-06-10T23:08:59.138532Z", + "iopub.status.idle": "2026-06-10T23:08:59.144336Z", + "shell.execute_reply": "2026-06-10T23:08:59.143433Z" + } + }, + "outputs": [], + "source": [ + "def sampled_curves(code, phis, shots, *, corner_flips=None, seed=1, postselect=False):\n", + " \"\"\"Estimate (phi) and mean ||(phi) with Tsim's sampler.\n", + "\n", + " With ``postselect=True`` the logical estimate keeps only shots whose\n", + " syndrome is trivial (all detectors = 0), matching the decoded/postselected\n", + " top panel of Fig. 4a.\n", + " \"\"\"\n", + " logical, stabilizer = [], []\n", + " for phi in phis:\n", + " circuit = rotated_circuit(code, phi, corner_flips=corner_flips)\n", + " sampler = circuit.compile_detector_sampler(seed=seed)\n", + " det, obs = sampler.sample(shots=shots, separate_observables=True)\n", + " if postselect and code[\"x_stabs\"]:\n", + " accepted = obs[~det.any(axis=1)]\n", + " logical.append(1.0 - 2.0 * accepted.mean() if accepted.size else np.nan)\n", + " else:\n", + " logical.append(1.0 - 2.0 * obs.mean())\n", + " if code[\"x_stabs\"]:\n", + " stabilizer.append(float(np.mean(np.abs(1.0 - 2.0 * det.mean(axis=0)))))\n", + " else:\n", + " stabilizer.append(np.nan)\n", + " return np.array(logical), np.array(stabilizer)" + ] + }, + { + "cell_type": "markdown", + "id": "e7348512", + "metadata": {}, + "source": [ + "Before the full sweep, confirm the headline claim: **Tsim's sampler runs the 3D [[15,1,3]] code**, where a generic angle is 15 non-Clifford rotations at once. We check a few angles against the exact values, including a generic $30^\\circ$ (the genuinely hard case) and $45^\\circ$ (15 transversal $T$ gates)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "dfc1f11d", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:08:59.146604Z", + "iopub.status.busy": "2026-06-10T23:08:59.145927Z", + "iopub.status.idle": "2026-06-10T23:09:55.998706Z", + "shell.execute_reply": "2026-06-10T23:09:55.997983Z" + } + }, + "outputs": [], + "source": [ + "check_phis = np.deg2rad([0, 30, 45, 90])\n", + "xl_exact, _ = exact_curves(reed_muller, check_phis)\n", + "xl_sampled, _ = sampled_curves(reed_muller, check_phis, shots=8_000)\n", + "print(\" phi exact Tsim\")\n", + "for phi, e, s in zip(np.rad2deg(check_phis), xl_exact, xl_sampled):\n", + " print(f\" {phi:5.1f} {e:+.3f} {s:+.3f}\")" + ] + }, + { + "cell_type": "markdown", + "id": "199fae3b", + "metadata": {}, + "source": [ + "### Postselection and acceptance fraction\n", + "\n", + "The paper sharpens its plateaus with *postselection*, accepting only shots whose stabilizers are satisfied (it quotes acceptance fractions of 46% / 74% for the 3D / 2D codes). We implement the same idea: project the rotated state onto the $+1$ eigenspace of the $X$-stabilizers and report the accepted $\\langle X_L\\rangle$ and the acceptance fraction. On *noiseless* data the curves already lie almost entirely in the code space, so postselection barely moves them and the acceptance stays high; with physical noise it is what removes the off-codespace weight (alongside the paper's purity normalization)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "eefd7f2d", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:09:56.000107Z", + "iopub.status.busy": "2026-06-10T23:09:55.999966Z", + "iopub.status.idle": "2026-06-10T23:09:57.231789Z", + "shell.execute_reply": "2026-06-10T23:09:57.231392Z" + } + }, + "outputs": [], + "source": [ + "def postselected_logical(code, phi, corner_flips=None):\n", + " \"\"\" conditioned on all X-stabilizers = +1, and the acceptance fraction.\"\"\"\n", + " n = code[\"n\"]\n", + " psi, idx = _rotated_statevector(code, phi, corner_flips)\n", + " for sup in code[\"x_stabs\"]: # project onto +1 eigenspace of each X-stabilizer\n", + " psi = 0.5 * (psi + psi[idx ^ _support_mask(sup, n)])\n", + " accept = float(np.real(np.vdot(psi, psi)))\n", + " if accept < 1e-12:\n", + " return 0.0, 0.0\n", + " xl = float(np.real(np.vdot(psi, psi[idx ^ _support_mask(code[\"x_logical\"], n)])) / accept)\n", + " return xl, accept\n", + "\n", + "\n", + "print(\" phi raw postsel acceptance\")\n", + "for deg in (0, 22.5, 30, 45):\n", + " raw, _ = exact_curves(reed_muller, [np.deg2rad(deg)])\n", + " xl_ps, acc = postselected_logical(reed_muller, np.deg2rad(deg))\n", + " print(f\" {deg:5.1f} {raw[0]:+.3f} {xl_ps:+.3f} {acc:.3f}\")" + ] + }, + { + "cell_type": "markdown", + "id": "31ae191e", + "metadata": {}, + "source": [ + "## Sweeping the global angle\n", + "\n", + "Now we sweep $\\varphi$ from $-180^\\circ$ to $180^\\circ$ for all four configurations of Fig. 4a:\n", + "\n", + "- **3D colour** - the [[15,1,3]] Reed-Muller code,\n", + "- **2D colour** - the [[7,1,3]] Steane code,\n", + "- **Unentangled** - 15 bare qubits in $|+\\rangle^{\\otimes 15}$, read out with the same Reed-Muller operators,\n", + "- **3D negative stabilizer** - the Reed-Muller code with the four tetrahedron corners flipped (the paper's control: incorrect stabilizer signs, which should *remove* the $45^\\circ$ plateau).\n", + "\n", + "Following the paper's analysis, the two *encoded* curves show the **decoded** logical projection - $\\langle X_L\\rangle$ postselected on a satisfied syndrome - which turns the oscillating raw expectation of a many-qubit code into the flat transversal-gate plateaus. The unentangled and negative-stabilizer controls are shown **raw**: with no (correctly signed) code space to project onto, $\\langle X_L\\rangle = \\cos^{15}\\varphi$ collapses away from $0^\\circ$ and $\\pm180^\\circ$, and the two control curves hug each other - exactly the behaviour that distinguishes a genuine code state from its look-alikes. The bottom panel shows the raw $X$-stabilizer revival for all four. The open markers are Tsim samples (syndrome-postselected, on the 2D curve and on the 3D code at the multiples of $45^\\circ$). Raise `SAMPLE_SHOTS` for tighter error bars." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "290026e1", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:09:57.247216Z", + "iopub.status.busy": "2026-06-10T23:09:57.246966Z", + "iopub.status.idle": "2026-06-10T23:14:43.425883Z", + "shell.execute_reply": "2026-06-10T23:14:43.424506Z" + } + }, + "outputs": [], + "source": [ + "phi_deg = np.linspace(-180, 180, 181)\n", + "phi_rad = np.deg2rad(phi_deg)\n", + "phi_dots = np.linspace(-180, 180, 73) # marker spacing, as in the figure\n", + "phi_dots_rad = np.deg2rad(phi_dots)\n", + "\n", + "# Top panel, as in Fig. 4a:\n", + "# - 3D and 2D codes: the *decoded* logical projection (postselected on a\n", + "# satisfied syndrome). The raw of a many-qubit code is multi-frequency\n", + "# and oscillates; projecting onto the code space recovers the flat\n", + "# transversal-gate plateaus.\n", + "# - Unentangled and 3D-negative-stabilizer: the *raw* projection. With no\n", + "# (correct) code space to project onto, = cos(phi)^15 collapses to\n", + "# zero away from 0 and +-180 degrees -- these two curves hug each other.\n", + "xl_3d = np.array([postselected_logical(reed_muller, p)[0] for p in phi_rad])\n", + "xl_2d = np.array([postselected_logical(steane, p)[0] for p in phi_rad])\n", + "xl_un = exact_curves(unentangled, phi_rad)[0]\n", + "xl_neg = exact_curves(reed_muller, phi_rad, corner_flips=reed_muller[\"corners\"])[0]\n", + "\n", + "dots_3d = np.array([postselected_logical(reed_muller, p)[0] for p in phi_dots_rad])\n", + "dots_2d = np.array([postselected_logical(steane, p)[0] for p in phi_dots_rad])\n", + "dots_un = exact_curves(unentangled, phi_dots_rad)[0]\n", + "dots_neg = exact_curves(reed_muller, phi_dots_rad, corner_flips=reed_muller[\"corners\"])[0]\n", + "\n", + "# Bottom panel: the raw X-stabilizer revival for all four configurations.\n", + "st_3d = exact_curves(reed_muller, phi_rad)[1]\n", + "st_2d = exact_curves(steane, phi_rad)[1]\n", + "st_un = exact_curves(unentangled, phi_rad)[1]\n", + "st_neg = exact_curves(reed_muller, phi_rad, corner_flips=reed_muller[\"corners\"])[1]\n", + "st_dots_3d = exact_curves(reed_muller, phi_dots_rad)[1]\n", + "st_dots_2d = exact_curves(steane, phi_dots_rad)[1]\n", + "st_dots_un = exact_curves(unentangled, phi_dots_rad)[1]\n", + "st_dots_neg = exact_curves(reed_muller, phi_dots_rad, corner_flips=reed_muller[\"corners\"])[1]\n", + "\n", + "# Tsim sampled overlay, postselected on a trivial syndrome to match the top panel.\n", + "SAMPLE_SHOTS = 8_000\n", + "phi_marks_2d = np.deg2rad(np.linspace(-180, 180, 13))\n", + "phi_marks_3d = np.deg2rad(np.arange(-180, 181, 45)) # Clifford + T angles: cheap to compile\n", + "xl_2d_s, _ = sampled_curves(steane, phi_marks_2d, shots=SAMPLE_SHOTS, postselect=True)\n", + "xl_3d_s, _ = sampled_curves(reed_muller, phi_marks_3d, shots=SAMPLE_SHOTS, postselect=True)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "62e0b7dd", + "metadata": { + "execution": { + "iopub.execute_input": "2026-06-10T23:14:43.427731Z", + "iopub.status.busy": "2026-06-10T23:14:43.427536Z", + "iopub.status.idle": "2026-06-10T23:14:43.644710Z", + "shell.execute_reply": "2026-06-10T23:14:43.643872Z" + } + }, + "outputs": [], + "source": [ + "def _unit(curve):\n", + " \"\"\"Normalize to unity maximum expectation, as in Fig. 4a (Methods).\"\"\"\n", + " peak = np.nanmax(np.abs(curve))\n", + " return curve / peak if peak > 0 else curve\n", + "\n", + "\n", + "C_3D, C_2D, C_UN, C_NEG = \"#c0392b\", \"#1f3a6b\", \"#4f8a4f\", \"#b9b9b9\"\n", + "gate_labels = [\"Z\", \"(TS)†\", \"S†\", \"T†\", \"I\", \"T\", \"S\", \"TS\", \"Z\"]\n", + "gate_ticks = np.arange(-180, 181, 45)\n", + "\n", + "fig, (ax_top, ax_bot) = plt.subplots(\n", + " 2, 1, figsize=(8.2, 5.4), sharex=True, gridspec_kw=dict(height_ratios=[1.55, 1])\n", + ")\n", + "\n", + "\n", + "def _curve(ax, x_line, y_line, x_dots, y_dots, color, label=None):\n", + " ax.plot(x_line, _unit(y_line), color=color, lw=1.1, zorder=1)\n", + " ax.plot(\n", + " x_dots, _unit(y_dots), \"o\", color=color, ms=3.2,\n", + " mfc=color, mec=\"white\", mew=0.4, label=label, zorder=2,\n", + " )\n", + "\n", + "\n", + "# Top panel: logical X projection.\n", + "_curve(ax_top, phi_deg, xl_3d, phi_dots, dots_3d, C_3D, \"3D colour\")\n", + "_curve(ax_top, phi_deg, xl_2d, phi_dots, dots_2d, C_2D, \"2D colour\")\n", + "_curve(ax_top, phi_deg, xl_un, phi_dots, dots_un, C_UN, \"Unentangled\")\n", + "_curve(ax_top, phi_deg, xl_neg, phi_dots, dots_neg, C_NEG, \"3D negative\\nstabilizer\")\n", + "ax_top.plot(np.rad2deg(phi_marks_3d), xl_3d_s, \"o\", color=C_3D, ms=5.5, mfc=\"white\", label=\"3D (Tsim)\")\n", + "ax_top.plot(np.rad2deg(phi_marks_2d), xl_2d_s, \"s\", color=C_2D, ms=4, mfc=\"white\", label=\"2D (Tsim)\")\n", + "ax_top.set_ylabel(\"Logical X\\nprojection\")\n", + "ax_top.set_ylim(-1.06, 1.06)\n", + "ax_top.set_yticks([-1, 0, 1])\n", + "ax_top.legend(loc=\"center left\", bbox_to_anchor=(1.01, 0.5), fontsize=8.5, handlelength=1.1, frameon=False)\n", + "for tick, label in zip(gate_ticks, gate_labels):\n", + " ax_top.text(tick, 1.13, label, ha=\"center\", fontsize=9)\n", + "\n", + "# Bottom panel: normalized absolute stabilizer, all four configurations.\n", + "_curve(ax_bot, phi_deg, st_3d, phi_dots, st_dots_3d, C_3D)\n", + "_curve(ax_bot, phi_deg, st_2d, phi_dots, st_dots_2d, C_2D)\n", + "_curve(ax_bot, phi_deg, st_un, phi_dots, st_dots_un, C_UN)\n", + "_curve(ax_bot, phi_deg, st_neg, phi_dots, st_dots_neg, C_NEG)\n", + "ax_bot.set_ylabel(\"Normalized abs.\\nstabilizer\")\n", + "ax_bot.set_ylim(-0.03, 1.03)\n", + "ax_bot.set_yticks([0, 0.5, 1.0])\n", + "ax_bot.set_xlabel(r\"Global $\\varphi$\")\n", + "ax_bot.set_xticks(range(-180, 181, 90))\n", + "ax_bot.set_xticklabels([\"-180°\", \"-90°\", \"0\", \"90°\", \"180°\"])\n", + "\n", + "fig.subplots_adjust(hspace=0.08)\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "id": "0aa2e9cf", + "metadata": {}, + "source": [ + "## Reading the figure\n", + "\n", + "The top panel is the logical $X$ projection; the bottom panel is the (normalized) magnitude of the $X$-stabilizers. The labels along the top are the logical gate realized at each angle.\n", + "\n", + "- **Unentangled (green).** Fifteen bare qubits read out with the code's operators: $\\langle X_L\\rangle = \\cos^{15}\\varphi$ collapses to zero away from $0^\\circ$ and $\\pm180^\\circ$, and the weight-8 stabilizers $\\langle X_S\\rangle = \\cos^{8}\\varphi$ revive only at multiples of $180^\\circ$. No entanglement, no protection.\n", + "- **2D colour / Steane (blue).** Decoded plateaus at multiples of $90^\\circ$, where the global rotation realizes a logical Clifford ($S, S^\\dagger, Z$); the stabilizers revive only there (dipping to $1/2$ in between).\n", + "- **3D colour / Reed-Muller (red).** An *extra* plateau at $\\pm 45^\\circ$ and $\\pm 135^\\circ$ - the transversal **$T, T^\\dagger, TS, (TS)^\\dagger$** gates - with the stabilizers reviving there too. This is exactly the feature the experiment uses for universality.\n", + "- **3D negative stabilizer (grey).** The same code with the four tetrahedron corners flipped. Its raw projection hugs the *unentangled* curve: with incorrect stabilizer signs there is no $45^\\circ$ plateau and no protection, even though the connectivity and gate count are identical to the red curve. Correct stabilizer signs - i.e. the right entangled state - are essential for the non-Clifford gate.\n", + "\n", + "The open Tsim markers sit on the decoded curves, including the red circles at $\\pm 45^\\circ$ - the points where the 15-qubit code is running 15 non-Clifford gates at once. Reproducing those is the heart of Fig. 4a, and is exactly the regime where Tsim's stabilizer-rank sampler is needed (a dense statevector of the rotated diagram is out of reach).\n", + "\n", + "### Why this matters\n", + "\n", + "The $45^\\circ$ plateau is the fingerprint of a transversal non-Clifford gate, and it only exists when the code carries genuine in-block entanglement - $|T_L\\rangle$ is a code-spanning superposition, whereas a logical Pauli is a product operator. In the paper this is the basis for universal logic via transversal teleportation between 2D and 3D codes, and it underlies an error-corrected Bell (CHSH) test that reaches $1.99(3)\\times\\sqrt2$, linking quantum contextuality to the entanglement that error correction must protect.\n", + "\n", + "### Things to try\n", + "\n", + "- Increase `SAMPLE_SHOTS` and the number of sampled markers and watch the Monte-Carlo points tighten onto the exact curves.\n", + "- Flip a different set of qubits via `corner_flips`, or rotate about another axis (`R_X`, `R_Y`), and watch the plateaus move.\n", + "- Build the high-rate [[16,6,4]] tesseract code with the hypercube encoder *without* puncturing (the full 4-cube) and study one of its six logical qubits." + ] + } + ], + "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.12.7" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/mkdocs.yml b/mkdocs.yml index acc5a219..7d685945 100644 --- a/mkdocs.yml +++ b/mkdocs.yml @@ -17,6 +17,7 @@ nav: - demos/encoding_demo.ipynb - demos/magic_state_distillation.ipynb - demos/from_stim_to_tsim.ipynb + - demos/global_rotations.ipynb - API Reference: reference/ diff --git a/test/unit/demos/__init__.py b/test/unit/demos/__init__.py new file mode 100644 index 00000000..e69de29b diff --git a/test/unit/demos/test_global_rotations.py b/test/unit/demos/test_global_rotations.py new file mode 100644 index 00000000..5fe31626 --- /dev/null +++ b/test/unit/demos/test_global_rotations.py @@ -0,0 +1,187 @@ +"""Tests for the ``docs/demos/global_rotations.ipynb`` tutorial (issue #119). + +These mirror the constructions used in the notebook and guard the physics it +reproduces from Fig. 4a of Bluvstein et al. (Nature 649, 39, 2026): the +Reed-Muller code encoders, the literal hypercube encoding circuit, and the +transversal-T plateau that distinguishes the 3D code from the 2D code. +""" + +import numpy as np +import stim + +import tsim + + +def _pauli(kind: str, support: list[int], n: int) -> str: + chars = ["_"] * n + for q in support: + chars[q] = kind + return "".join(chars) + + +def _bit(j: int, t: int) -> int: + return (j >> t) & 1 + + +# [[7,1,3]] Steane and [[15,1,3]] Reed-Muller code definitions. +STEANE = { + "n": 7, + "x_stabs": [[3, 4, 5, 6], [1, 2, 5, 6], [0, 2, 4, 6]], + "z_stabs": [[3, 4, 5, 6], [1, 2, 5, 6], [0, 2, 4, 6]], + "x_logical": list(range(7)), +} +_RM_X = [[j for j in range(15) if _bit(j + 1, k)] for k in range(4)] +_RM_Z = [list(s) for s in _RM_X] + [ + [j for j in range(15) if _bit(j + 1, a) and _bit(j + 1, b)] + for a in range(4) + for b in range(a + 1, 4) +] +REED_MULLER = { + "n": 15, + "x_stabs": _RM_X, + "z_stabs": _RM_Z, + "x_logical": list(range(15)), + "corners": [0, 1, 3, 7], +} + + +def _encoder(code: dict, corner_flips: list[int] | None = None) -> stim.Circuit: + n = code["n"] + gens = [stim.PauliString(_pauli("X", s, n)) for s in code["x_stabs"]] + gens += [stim.PauliString(_pauli("Z", s, n)) for s in code["z_stabs"]] + gens.append(stim.PauliString(_pauli("X", code["x_logical"], n))) + circuit = stim.Tableau.from_stabilizers( + gens, allow_redundant=False, allow_underconstrained=False + ).to_circuit("elimination") + if corner_flips: + circuit += stim.Circuit("X " + " ".join(str(q) for q in corner_flips)) + return circuit + + +def _plus_logical_statevector( + code: dict, corner_flips: list[int] | None = None +) -> np.ndarray: + sim = stim.TableauSimulator() + sim.do_circuit(_encoder(code, corner_flips)) + return np.asarray(sim.state_vector(endian="big"), dtype=np.complex128) + + +def _support_mask(support: list[int], n: int) -> int: + mask = 0 + for q in support: + mask |= 1 << (n - 1 - q) + return mask + + +def _expectations(code: dict, phi: float, corner_flips: list[int] | None = None): + n = code["n"] + psi0 = _plus_logical_statevector(code, corner_flips) + idx = np.arange(psi0.size) + popcount = np.array([int(b).bit_count() for b in idx]) + psi = psi0 * np.exp(1j * phi * popcount) + logical = float( + np.real(np.vdot(psi, psi[idx ^ _support_mask(code["x_logical"], n)])) + ) + stab = float( + np.mean( + [ + abs(np.real(np.vdot(psi, psi[idx ^ _support_mask(s, n)]))) + for s in code["x_stabs"] + ] + ) + ) + return logical, stab + + +def _hypercube_plus_logical(m: int, r_x: int, r_z: int) -> np.ndarray: + """Literal hypercube encoder (Gong-Renes / Ext. Data Fig. 10a).""" + del r_z # the |0> rows are the complement of the |+> rows for these codes + n = 1 << m + lines = [f"H {q}" for q in range(n) if q == 0 or bin(q).count("1") <= r_x] + lines += [ + f"CX {q} {q | (1 << t)}" for t in range(m) for q in range(n) if _bit(q, t) == 0 + ] + circuit = stim.Circuit("\n".join(lines)) + sim = stim.TableauSimulator() + sim.do_circuit(circuit) + psi = np.asarray(sim.state_vector(endian="big"), dtype=np.complex128).reshape( + [2] * n + ) + # Puncture qubit 0 in the X basis (project onto |+>) and drop it. + branch0 = psi[(0, *([slice(None)] * (n - 1)))].reshape(-1) + branch1 = psi[(1, *([slice(None)] * (n - 1)))].reshape(-1) + sub = branch0 + branch1 + return sub / np.linalg.norm(sub) + + +def test_flow_generators_prepare_logical_x(): + # The encoder maps an input Z onto the all-X logical X_L (PauliString index 1 == X). + enc = _encoder(REED_MULLER) + all_x_flows = [ + f + for f in enc.flow_generators() + if all(f.output_copy()[q] == 1 for q in range(15)) + ] + assert all_x_flows + + +def test_reed_muller_stabilizer_weights(): + assert sorted(len(s) for s in REED_MULLER["x_stabs"]) == [8, 8, 8, 8] + assert sorted(len(s) for s in REED_MULLER["z_stabs"]) == [ + 4, + 4, + 4, + 4, + 4, + 4, + 8, + 8, + 8, + 8, + ] + + +def test_phi_zero_is_plus_logical(): + for code in (STEANE, REED_MULLER): + logical, stab = _expectations(code, 0.0) + assert np.isclose(logical, 1.0, atol=1e-6) + assert np.isclose(stab, 1.0, atol=1e-6) + + +def test_hypercube_matches_stabilizer_encoder(): + # The literal hypercube circuit prepares the same |+_L> state, exactly. + steane_sv = _plus_logical_statevector(STEANE) + rm_sv = _plus_logical_statevector(REED_MULLER) + assert abs(np.vdot(_hypercube_plus_logical(3, 1, 1), steane_sv)) > 1 - 1e-6 + assert abs(np.vdot(_hypercube_plus_logical(4, 1, 2), rm_sv)) > 1 - 1e-6 + + +def test_reed_muller_has_transversal_t_plateau(): + # 3D code: at 45 degrees the global rotation is a logical T and the + # stabilizers revive (transversal T). 2D Steane code: they do not. + rm_logical, rm_stab = _expectations(REED_MULLER, np.deg2rad(45)) + assert np.isclose(rm_logical, 1 / np.sqrt(2), atol=1e-6) + assert np.isclose(rm_stab, 1.0, atol=1e-6) + _, steane_stab = _expectations(STEANE, np.deg2rad(45)) + assert steane_stab < 0.9 + + +def test_negative_stabilizer_removes_plateau(): + # Flipping the four tetrahedron corners destroys the 45-degree plateau. + logical, stab = _expectations( + REED_MULLER, np.deg2rad(45), corner_flips=REED_MULLER["corners"] + ) + assert abs(logical) < 1e-6 + assert stab < 1e-6 + + +def test_tsim_sampler_matches_exact_on_3d_code(): + # tsim's stabilizer-rank sampler runs the 15-rotation 3D code. + n = REED_MULLER["n"] + text = str(_encoder(REED_MULLER)) + "\n" + text += "R_Z(0.25) " + " ".join(map(str, range(n))) + "\n" + text += "MPP " + "*".join(f"X{q}" for q in REED_MULLER["x_logical"]) + "\n" + text += "OBSERVABLE_INCLUDE(0) rec[-1]\n" + sampler = tsim.Circuit(text).compile_detector_sampler(seed=1) + obs = sampler.sample(shots=4000, separate_observables=True)[1] + assert abs((1.0 - 2.0 * obs.mean()) - 1 / np.sqrt(2)) < 0.05