Level 2 · Quantum Computing: Qubits and Gates

Level 1 built the toolbox: arrows with a weight, inner products with a conjugate, matrices that keep every length. This level gives those objects their quantum names. A unit vector of complex numbers is a state, the weights of its entries are the odds of what a measurement shows, and a unitary matrix is a gate. The first half is about states and reading them: the Born rule, outcome labels as bit strings, seeded sampling of many shots, the phases that counts can never show, measuring in other bases, and the Bloch sphere, a globe with a point for every qubit. The second half is about gates: the named gates, phase gates and rotations, when two gates count as the same, circuits that run left to right while their matrix multiplies right to left, undoing a circuit, and interference, where two paths to the same outcome cancel.

From this level on, numpy is loaded for every quantum problem: import numpy as np just works. Nothing needs it: the lists and loops of Level 1 solve every problem of this level, and the lesson shows both where it helps.

The paths meet again. Drawing many shots with a seeded random.Random is the simulation technique of the Statistics path at Level 2. The phase gates multiply by the turning arrows e^{i\lambda} of the Fourier transform in Signals, Systems & Control. Outcome labels are bit strings, read with the shifts and masks of the Algorithms path, and the running totals of the sampling procedure are a sorted list that binary search can scan.

numpy: complex arrays

A numpy array holds numbers of one type, and dtype=complex makes them complex. @ is the matrix product, .conj().T the dagger, and np.allclose compares with a tolerance. The list code of Level 1 still works side by side:

import numpy as np
X = [[0, 1], [1, 0]]
psi = [0.8, 0.6j]
print([sum(X[i][j] * psi[j] for j in range(2)) for i in range(2)])

Xa = np.array(X, dtype=complex)
psia = np.array(psi, dtype=complex)
print((Xa @ psia).tolist())
print((Xa * psia).tolist())
print(Xa.shape, psia.shape, np.allclose(Xa.conj().T @ Xa, np.eye(2)))

[0.6j, (0.8+0j)]
[0.6j, (0.8+0j)]
[[0j, 0.6j], [(0.8+0j), 0j]]
(2, 2) (2,) True

Xa @ psia is the matrix-vector product, the same as the list version. Xa * psia is not: * multiplies entry by entry (each row by psia), which is never what a gate does. np.eye(2) is the identity matrix.

Print arrays through .tolist(), which turns them into plain lists. Floats leave rounding noise behind, so round for printing; np.round rounds both parts of complex numbers, and adding 0 turns a -0.0 into 0.0. The examples below use a small helper show for that:

import numpy as np, math
H = np.array([[1, 1], [1, -1]], dtype=complex) / math.sqrt(2)
v = H @ H @ np.array([0.6, 0.8j])
print(v.tolist())

def show(a, digits=4):
    return (np.round(np.asarray(a, dtype=complex), digits) + 0).tolist()

print(show(v), show(H))

[(0.5999999999999999-1.7896914536605906e-17j), (-1.3422685902454429e-17+0.7999999999999998j)]
[(0.6+0j), 0.8j] [[(0.7071+0j), (0.7071+0j)], [(0.7071+0j), (-0.7071+0j)]]

A few more tools. np.vdot(u, v) is the inner product \langle u|v\rangle, conjugating the first vector as Level 1 did; np.linalg.norm(v) is the length; np.abs(v) ** 2 gives every weight at once. np.kron is the Kronecker product, which joins the states of several qubits into one; it does the heavy lifting next level.

import numpy as np, math
u = np.array([1, 1j]) / math.sqrt(2)
v = np.array([1, -1j]) / math.sqrt(2)
print(complex(np.vdot(u, v)), round(float(np.linalg.norm(u)), 6), np.round(np.abs(u) ** 2, 4).tolist())
print(np.allclose(u, [0.7071067811865476, 0.7071067811865476j]), np.allclose(u, v))
print(np.kron([1, 0], [0.6, 0.8]).tolist())

0j 1.0 [0.5, 0.5]
True False
[0.6, 0.8, 0.0, 0.0]

Wrap single numbers in complex(...) or float(...) before printing: numpy 2 prints its own numbers as np.complex128(...) inside lists and tuples. Answers you return may be numpy arrays or lists; the grader treats them the same.

Common mistakes: * where @ was meant. Forgetting dtype=complex: an array made from real numbers is real, and storing 1j into it raises an error. Printing raw arrays, whose layout changes between numpy versions. np.allclose on two lists of different lengths, which raises an error instead of saying False.

A qubit is a unit vector

A qubit is a pair of complex numbers [a, b] with |a|^2 + |b|^2 = 1, written

where |0\rangle = [1, 0] and |1\rangle = [0, 1] are the two basis states. The entries a and b are the amplitudes. n qubits have 2^n basis states, one per bit string of length n, and the state is a list of 2^n amplitudes: amplitude i belongs to basis state i. For two qubits, [a0, a1, a2, a3] means a_0|00\rangle + a_1|01\rangle + a_2|10\rangle + a_3|11\rangle.

A state must be a unit vector. A list that is not one, say from a log that dropped a scale factor, describes the same state as the list divided by its length:

import math
import numpy as np
raw = [2, 1j, 0, -2]
length = math.sqrt(sum(abs(x) ** 2 for x in raw))
state = [x / length for x in raw]
print(length, [complex(round(x.real, 4), round(x.imag, 4)) for x in map(complex, state)])

s = np.array(raw, dtype=complex)
s = s / np.linalg.norm(s)
print(round(float(np.linalg.norm(s)), 12))

3.0 [(0.6667+0j), 0.3333j, 0j, (-0.6667+0j)]
1.0

The weights of raw are 4, 1, 0 and 4, so its length is 3. Two qubits, four amplitudes; basis state |10\rangle (index 2) has amplitude 0.

Common mistakes: a state of length 3 or 6: the length is always a power of 2. Normalising with x ** 2 instead of abs(x) ** 2. Dividing by the total weight instead of its square root.

The Born rule

Measuring all the qubits of a state gives one basis state, at random. Outcome i comes out with probability

the weight of its amplitude. This is the Born rule, and it is the only way the amplitudes ever reach the outside world. For a unit vector the probabilities add up to 1:

import numpy as np
state = [0.1, 0.7j, -0.7, 0.1]
p = [abs(a) ** 2 for a in state]
print([round(x, 4) for x in p], round(sum(p), 10))
print(np.round(np.abs(np.array(state)) ** 2, 4).tolist())
print([a ** 2 for a in state])
plot(["00", "01", "10", "11"], p, kind="bar", title="odds of each outcome")

[0.01, 0.49, 0.49, 0.01] 1.0
[0.01, 0.49, 0.49, 0.01]
[0.010000000000000002, (-0.48999999999999994+0j), 0.48999999999999994, 0.010000000000000002]

Press Run to see the odds as bars: |01\rangle and |10\rangle come out half the time each, |00\rangle and |11\rangle only one time in a hundred each. The amplitude 0.7j and the amplitude -0.7 give the same odds: a minus sign or a factor 1j changes nothing here. The last line is the classic slip: squaring the arrow itself turns 0.7j into -0.49, a "probability" below zero. Take the size first, then square; (1j) ** 2 is -1, but abs(1j) ** 2 is 1.

Aside (optional): what a measurement does. On a real machine the reading also changes the state: after showing outcome i, the qubits are in basis state i, and measuring again shows i again. So nobody can read the same state twice to get a second sample; every shot starts with a freshly prepared copy. Simulators skip this step and just draw from the odds.

Common mistakes: a ** 2 instead of abs(a) ** 2. a.real ** 2, which forgets the imaginary part. Applying the Born rule to a list that is not normalised.

Practice: Raw Amplitudes from the Chip Log

Outcome labels: bit strings

Outcomes are reported as bit strings: outcome i is i written in binary with n digits. Qubit 0 is the leftmost bit, the most significant one. In 3 qubits:

format(i, f"0{n}b") writes the label, int(bits, 2) reads it back, and the bit of qubit q is (i >> (n - 1 - q)) & 1: shift right until that bit is last, then keep only the last bit. bits[q] gives the same bit as a character.

n = 3
for i in [0, 1, 3, 6]:
    bits = format(i, f"0{n}b")
    print(i, bits, int(bits, 2), [(i >> (n - 1 - q)) & 1 for q in range(n)], [int(c) for c in bits])
print(format(5, "b"), format(5, "03b"), format(5, "06b"))

0 000 0 [0, 0, 0] [0, 0, 0]
1 001 1 [0, 0, 1] [0, 0, 1]
3 011 3 [0, 1, 1] [0, 1, 1]
6 110 6 [1, 1, 0] [1, 1, 0]
101 101 000101

Without the 0{n} part, format drops the leading zeros: outcome 1 of three qubits must read 001, not 1. With labels in hand you can ask about one qubit alone: qubit 0 reads 1 in every outcome whose leftmost bit is 1, so its probability is the sum of their probabilities. For the state of the last section:

state = [0.1, 0.7j, -0.7, 0.1]
n = 2
p1 = sum(abs(a) ** 2 for i, a in enumerate(state) if (i >> (n - 1 - 0)) & 1)
print("qubit 0 reads 1:", round(p1, 4), "from", [format(i, "02b") for i in range(4) if (i >> 1) & 1])

qubit 0 reads 1: 0.5 from ['10', '11']

Aside (optional): the other order. Some toolkits number the qubits the other way, with qubit 0 as the rightmost bit. Neither choice is wrong, but mixing them silently swaps qubits. On this path qubit 0 is always the leftmost bit.

Common mistakes: (i >> q) & 1, which counts from the right and swaps qubits. bin(i), which adds a 0b prefix and drops leading zeros. Reading the label "011" as index 6.

Many shots, with a seed

A quantum computer never shows its amplitudes. You prepare the same state many times, measure each copy (one shot), and get back counts: how often each label came out. A simulator imitates this with random numbers, and problems on this path fix the procedure exactly so that every correct program prints the same counts:

Use rng = random.Random(seed). For each shot draw r = rng.random() and pick the first outcome i whose running total p[0] + ... + p[i] is greater than r (if rounding leaves no such i, pick the last outcome with nonzero probability).

Lay the probabilities end to end on the stretch from 0 to 1: outcome 0 owns the first p_0 of it, outcome 1 the next p_1, and so on. rng.random() lands anywhere in [0, 1) with equal chance, so it falls in outcome i's piece with probability p_i. Here is the procedure on something that is not a qubit at all, a spinner with three coloured sectors:

import random
names = ["red", "blue", "gold"]
weights = [0.5, 0.3, 0.2]
rng = random.Random(2024)
for spin in range(6):
    r = rng.random()
    total = 0.0
    for i, w in enumerate(weights):
        total += w
        if total > r:
            break
    print(round(r, 3), "->", names[i])

0.47 -> red
0.728 -> blue
0.304 -> red
0.887 -> gold
0.41 -> red
0.717 -> blue

0.728 is past red's half but inside red plus blue (0.8), so it picks blue; 0.887 is past 0.8, so gold. The generator is made once, before the loop: a new random.Random(2024) per shot would draw the same r every time. The same seed always gives the same sequence, which is what makes counts checkable.

The fallback clause is for rounding. A running total of probabilities computed from amplitudes rarely ends at exactly 1:

total = 0.0
for p in [0.1] * 10:
    total += p
print(total, total == 1)

0.9999999999999999 False

If r happens to be above such a total, no outcome passes the test, and the rule picks the last outcome whose probability is not 0 (never one that cannot occur). For long lists, the running totals are sorted, so bisect.bisect_right(totals, r) finds the first total greater than r by binary search.

Common mistakes: seeding inside the loop. >= instead of >, or walking from the last outcome. Using random.random(), the shared generator, which other code may also draw from. Listing outcomes that never came out with a count of 0 when the task wants only those that came out.

Practice: Ten Shots at the Same Qubit

Counts never show phases

Multiply every amplitude by the same number of size 1, e^{i\varphi}. Every weight stays the same, so every probability stays the same, and (as the next section shows) that holds for every kind of measurement. Such a common factor is a global phase, and the two lists are the same state. A factor on only some amplitudes, a relative phase, is different: the odds of 0 and 1 stay the same, but the state changes.

import cmath
a = [0.28, 0.96]
glob = [cmath.exp(2.1j) * x for x in a]
rel = [0.28, -0.96]
for name, s in (("a", a), ("global", glob), ("relative", rel)):
    print(name, [round(abs(x) ** 2, 4) for x in s])

a [0.0784, 0.9216]
global [0.0784, 0.9216]
relative [0.0784, 0.9216]

All three give 0 with probability 0.0784 and 1 with probability 0.9216, so counts in this basis cannot tell them apart. Turned around, counts can estimate the sizes of the amplitudes (the share of shots that showed i estimates |a_i|^2, and its square root estimates |a_i|), but never their phases.

glob really is a: no experiment can ever tell them apart. rel is a different state, and the next section shows a measurement that sees the difference. That is why answers that are states (and, later, gates) are graded up to a global phase: [1j, 0] counts as [1, 0], but [0.28, -0.96] does not count as [0.28, 0.96].

Common mistakes: calling [0.28, -0.96] "the same up to phase" because the sizes match; the factor must be the same for every amplitude. Expecting the phases back from counts.

Practice: Working Backwards from a Tally Sheet, How Many Different States Did the Log Record?

Measuring in another basis

Measuring as 0 or 1 is one question. Any two orthogonal unit vectors e_1, e_2 make another one, a basis, and the Born rule becomes

with the conjugate on the first vector. For e_k = |0\rangle this is just |a|^2 again. The most used other basis is the X basis:

It sees exactly the relative phase that counts of 0 and 1 hide:

import cmath, math
r = 1 / math.sqrt(2)
plus, minus = [r, r], [r, -r]
def ip(u, v):
    return sum(x.conjugate() * y for x, y in zip(u, v))
for name, s in (("a", [0.28, 0.96]), ("global", [cmath.exp(2.1j) * 0.28, cmath.exp(2.1j) * 0.96]), ("relative", [0.28, -0.96])):
    print(name, round(abs(ip(plus, s)) ** 2, 4), round(abs(ip(minus, s)) ** 2, 4))

a 0.7688 0.2312
global 0.7688 0.2312
relative 0.2312 0.7688

The global-phase twin gives the same odds again: |\langle e|c\,\psi\rangle|^2 = |c|^2|\langle e|\psi\rangle|^2 and |c| = 1. The relative phase swaps the odds of + and −.

A third basis, the Y basis, uses |{+i}\rangle = \tfrac{1}{\sqrt2}(|0\rangle + i|1\rangle) and |{-i}\rangle = \tfrac{1}{\sqrt2}(|0\rangle - i|1\rangle). Its entries are complex, so here the conjugate matters:

import math
r = 1 / math.sqrt(2)
plus_i, minus_i = [r, 1j * r], [r, -1j * r]
psi = [0.28, 0.96j]
with_conj = sum(x.conjugate() * y for x, y in zip(plus_i, psi))
without = sum(x * y for x, y in zip(plus_i, psi))
print(round(abs(with_conj) ** 2, 4), round(abs(without) ** 2, 4))
print("overlap of +i and -i:", abs(sum(x.conjugate() * y for x, y in zip(plus_i, minus_i))) ** 2)

0.7688 0.2312
overlap of +i and -i: 0.0

With the conjugate, \langle{+i}|\psi\rangle = \tfrac{1}{\sqrt2}(0.28 + (-i)(0.96i)) = \tfrac{1.24}{\sqrt2}; without it the two terms fight and give the odds of the other outcome. The same number |\langle u|v\rangle|^2 for any two states is their overlap: 1 when they are the same state up to a global phase, 0 when they are orthogonal, like |{+i}\rangle and |{-i}\rangle, which one measurement always tells apart.

Common mistakes: dropping the conjugate (harmless for |+\rangle, wrong for |{+i}\rangle). Squaring the inner product itself instead of its size. Forgetting the \frac{1}{\sqrt2} in the basis vectors, which doubles every probability.

Practice: Asking the Qubit a Different Question, Which Signals Can the Receiver Tell Apart?

The average reading

A readout often reports a qubit as +1 when it reads 0 and -1 when it reads 1. The average of many such readings tends to the expectation value of Z,

a number between -1 (always 1) and +1 (always 0). It is 0 for any state with even odds, such as |+\rangle. A simulation of 10000 shots agrees with the formula to about two digits:

import random
a, b = 0.28, 0.96
p0, p1 = abs(a) ** 2, abs(b) ** 2
print(round(p0 - p1, 4))
rng = random.Random(3)
readings = [1 if rng.random() < p0 else -1 for _ in range(10000)]
print(sum(readings) / len(readings))

-0.8432
-0.848

For a qubit inside a larger state, p_0 and p_1 are that one qubit's odds: add the probabilities of all outcomes whose bit for that qubit is 0, and of those whose bit is 1, as in the bit-string section.

Common mistakes: p_1 - p_0. Weighting by the amplitudes instead of their weights. Reading the wrong bit of the index for the qubit.

Practice: The Average Reading of One Qubit

The Bloch sphere

Up to a global phase, every qubit state can be written with two angles:

First turn the whole state so that a is real and not negative (a global phase, so allowed); then |a| = \cos\frac\theta2, |b| = \sin\frac\theta2, and \varphi is how far the angle of b is ahead of the angle of a. The angles are the latitude and longitude of a point on the unit sphere, the Bloch sphere:

|0\rangle is the north pole, |1\rangle the south pole, and the six states met so far sit on the axes:

import cmath, math
def from_angles(theta, phi):
    return [math.cos(theta / 2), cmath.exp(1j * phi) * math.sin(theta / 2)]
def point(theta, phi):
    return tuple(round(v, 4) + 0.0 for v in (math.sin(theta) * math.cos(phi), math.sin(theta) * math.sin(phi), math.cos(theta)))
half = math.pi / 2
for name, t, p in [("|0>", 0, 0), ("|1>", math.pi, 0), ("|+>", half, 0), ("|->", half, math.pi), ("|+i>", half, half), ("|-i>", half, 3 * half)]:
    s = from_angles(t, p)
    print(name.ljust(5), [complex(round(x.real, 4) + 0.0, round(x.imag, 4) + 0.0) for x in map(complex, s)], point(t, p))

|0>   [(1+0j), 0j] (0.0, 0.0, 1.0)
|1>   [0j, (1+0j)] (0.0, 0.0, -1.0)
|+>   [(0.7071+0j), (0.7071+0j)] (1.0, 0.0, 0.0)
|->   [(0.7071+0j), (-0.7071+0j)] (-1.0, 0.0, 0.0)
|+i>  [(0.7071+0j), 0.7071j] (0.0, 1.0, 0.0)
|-i>  [(0.7071+0j), -0.7071j] (0.0, -1.0, 0.0)

The half angles belong to the state and the whole angles to the point: |1\rangle has \theta = \pi but \cos\frac\pi2 = 0. Opposite points are orthogonal states. The equator (\theta = \pi/2) holds every state with even odds of 0 and 1; the relative phase \varphi is how far round the equator it sits. And z = \cos\theta = \cos^2\frac\theta2 - \sin^2\frac\theta2 = p_0 - p_1 is the average reading \langle Z\rangle of the last section.

At the poles one amplitude is 0, so \varphi means nothing there; the convention is \varphi = 0 whenever \theta is within 1e-9 of 0 or \pi.

Two tools read the angles back. A difference of angles is unchanged by a global phase, because a global phase turns both arrows by the same amount; cmath.phase answers in (-\pi, \pi], so a difference can land anywhere in (-2\pi, 2\pi), and % (2 * math.pi) brings it into [0, 2\pi):

import cmath, math
z1, z2 = 1j, -1 - 1j
for turn in (1, cmath.exp(2.5j)):
    w1, w2 = turn * z1, turn * z2
    diff = cmath.phase(w2) - cmath.phase(w1)
    print(round(cmath.phase(w1), 4), round(cmath.phase(w2), 4), round(diff, 4), round(diff % (2 * math.pi), 4))
try:
    math.acos(1.0000000000000002)
except ValueError:
    print("acos refuses", 2 * math.atan2(0.0, 1.0000000000000002))

1.5708 -2.3562 -3.927 2.3562
-2.2124 0.1438 2.3562 2.3562
acos refuses 0.0

Both rows end with the same 3\pi/4 = 2.3562: the second arrow is 135 degrees ahead of the first, whatever the common turn. For the sizes, \theta = 2\arccos|a| is fine in exact maths, but rounding can make abs(a) a hair above 1, and math.acos then raises an error; 2 * math.atan2(abs(b), abs(a)) gives the same angle and never fails.

Common mistakes: forgetting the factor 2 between the state's half angle and \theta. Taking the angle of b alone, which still contains the global phase. Leaving \varphi negative. Trusting cmath.phase(0) at a pole.

Practice: Pinning a Qubit on the Globe, Latitude and Longitude of a Qubit, Rebuilding a Qubit from Three Tally Sheets

Gates are unitary matrices

A one-qubit gate is a 2\times 2 unitary matrix U, a list of rows, and it changes a state by the matrix-vector product, new = U @ state, or with lists [U[0][0]*a + U[0][1]*b, U[1][0]*a + U[1][1]*b]. Unitary means U^\dagger U = I, so every gate keeps the length 1 of a state, and probabilities keep adding up to 1. The four most common ones:

import numpy as np, math
def show(a, digits=4):
    return (np.round(np.asarray(a, dtype=complex), digits) + 0).tolist()
X = np.array([[0, 1], [1, 0]], dtype=complex)
Y = np.array([[0, -1j], [1j, 0]])
Z = np.array([[1, 0], [0, -1]], dtype=complex)
H = np.array([[1, 1], [1, -1]], dtype=complex) / math.sqrt(2)
psi = np.array([0.28, 0.96j])
for name, U in (("X", X), ("Y", Y), ("Z", Z), ("H", H)):
    print(name, show(U @ psi), np.allclose(U.conj().T @ U, np.eye(2)))
print("H|0> =", show(H @ [1, 0]), " H|1> =", show(H @ [0, 1]))

X [0.96j, (0.28+0j)] True
Y [(0.96+0j), 0.28j] True
Z [(0.28+0j), -0.96j] True
H [(0.198+0.6788j), (0.198-0.6788j)] True
H|0> = [(0.7071+0j), (0.7071+0j)]  H|1> = [(0.7071+0j), (-0.7071+0j)]

X swaps the two amplitudes: on |0\rangle and |1\rangle it is the NOT of a classical bit. Z flips the sign of the |1\rangle amplitude, a relative phase: it changes no odds of 0 and 1 but swaps |+\rangle and |-\rangle. Y does both, with factors \pm i. H, the Hadamard gate, turns |0\rangle into |+\rangle and |1\rangle into |-\rangle: it swaps the Z basis and the X basis. Together with the identity I, which changes nothing, these are the everyday gates. The columns of a gate are what it does to |0\rangle and |1\rangle: U @ [1, 0] is column 0.

Common mistakes: j for 1j. Forgetting the \frac{1}{\sqrt2} in H, which doubles the length. Putting -1j in the wrong corner of Y (that gives -Y). Testing a gate only on |0\rangle and |1\rangle, where Z looks like the identity.

Phase gates

The phase gate turns the |1\rangle amplitude by an angle \lambda and leaves the |0\rangle amplitude alone:

Three of them have names: Z = P(\pi), the S gate S = P(\pi/2) = \mathrm{diag}(1, i), and the T gate T = P(\pi/4) = \mathrm{diag}(1, e^{i\pi/4}). Their daggers turn the other way: S^\dagger = \mathrm{diag}(1, -i) and T^\dagger = \mathrm{diag}(1, e^{-i\pi/4}), written "sdg" and "tdg" in gate lists. A phase gate never changes the odds of 0 and 1; on the Bloch sphere it turns the point about the vertical axis by \lambda. Starting from |+\rangle on the equator:

import numpy as np, cmath, math
def P(lam):
    return np.diag([1, cmath.exp(1j * lam)])
S, T = P(math.pi / 2), P(math.pi / 4)
r = 1 / math.sqrt(2)
plus = np.array([r, r])
print((np.round(S @ plus, 4) + 0).tolist(), (np.round(T @ plus, 4) + 0).tolist())
print((np.round(S.conj().T @ S, 4) + 0).tolist())
xs, ys = [], []
for k in range(12):
    a, b = P(k * math.pi / 6) @ plus
    phi = cmath.phase(b) - cmath.phase(a)
    xs.append(math.cos(phi))
    ys.append(math.sin(phi))
plot(xs, ys, kind="scatter", title="P(k pi/6)|+> seen from above the north pole")

[(0.7071+0j), 0.7071j] [(0.7071+0j), (0.5+0.5j)]
[[(1+0j), 0j], [0j, (1+0j)]]

S takes |+\rangle to |{+i}\rangle, a quarter turn round the equator; T goes half as far. Press Run to see twelve steps of \pi/6 walk once round the equator. np.diag builds a diagonal matrix from its diagonal.

Common mistakes: math.exp(1j * lam), which refuses complex numbers. Putting the phase on the |0\rangle amplitude. Treating "s" as its own inverse: S\,S = Z, not I.

Practice: One Pulse on the Test Chip, Stocking the Gate Drawer

Rotation gates

Three families of gates turn the Bloch point by any angle \theta about the x, y or z axis:

Note the half angle \theta/2 inside, the same half angle as in the Bloch angles: a turn of the point by \theta is a turn of the amplitudes by \theta/2. So R_y(\pi) takes the north pole |0\rangle to the south pole |1\rangle, and a full turn R_y(2\pi) brings every point back but multiplies the state by -1:

import numpy as np, cmath, math
def show(a, digits=4):
    return (np.round(np.asarray(a, dtype=complex), digits) + 0).tolist()
def ry(t):
    c, s = math.cos(t / 2), math.sin(t / 2)
    return np.array([[c, -s], [s, c]], dtype=complex)
def rz(t):
    return np.diag([cmath.exp(-0.5j * t), cmath.exp(0.5j * t)])
def P(lam):
    return np.diag([1, cmath.exp(1j * lam)])
print(show(ry(math.pi) @ [1, 0]), show(ry(2 * math.pi)))
print(show(rz(0.8)))
print(show(cmath.exp(-0.4j) * P(0.8)))

[0j, (1+0j)] [[(-1+0j), 0j], [0j, (-1+0j)]]
[[(0.9211-0.3894j), 0j], [0j, (0.9211+0.3894j)]]
[[(0.9211-0.3894j), 0j], [0j, (0.9211+0.3894j)]]

-I is the identity times a global phase, so it does nothing measurable. The last two lines show that R_z(\lambda) = e^{-i\lambda/2}P(\lambda): the same gate as P(\lambda) up to a global phase. A rotation by -\theta undoes a rotation by \theta, and so does a rotation by 4\pi - \theta (two full turns of the amplitudes).

Common mistakes: \theta instead of \theta/2 in the entries. Dropping a minus sign in R_x or R_y, which turns the other way. math.exp in R_z.

Practice: Turning the Rotation Dial

When two gates are the same

A gate U and the gate cU, with one number c of size 1, act the same on every state: (cU)\psi = c\,(U\psi), so the outputs differ by a global phase only, and no measurement can tell them apart. So, like states, gates are equal up to a global phase: A and B are the same gate when one c with |c| = 1 makes every entry of A equal to the matching entry of cB. Answers that are gates or circuit matrices are graded this way.

To test it, read c off one entry where B is far from 0 (the entry of largest size is a safe choice; never divide by an entry that is 0), then check that this one c fits all four entries. Here R_z(0.8) against P(0.8), which match, and against P(-0.8), which has entries of the same sizes but does not:

import numpy as np, cmath
def rz(t):
    return np.diag([cmath.exp(-0.5j * t), cmath.exp(0.5j * t)])
def P(lam):
    return np.diag([1, cmath.exp(1j * lam)])
A = rz(0.8)
for B in (P(0.8), P(-0.8)):
    a, b = A.ravel(), B.ravel()
    k = int(np.argmax(np.abs(b)))
    c = a[k] / b[k]
    print(k, round(abs(c), 6), round(cmath.phase(c), 4), np.allclose(a, c * b, rtol=0, atol=1e-9))

0 1.0 -0.4 True
0 1.0 -0.4 False

.ravel() flattens the 2\times 2 matrix into four numbers. In both cases c = e^{-0.4i} from the top-left entry; for P(-0.8) the bottom-right entry then fails, since e^{0.4i} \ne e^{-0.4i}e^{-0.8i}. Comparing only the sizes of the entries would call them equal, and so would comparing each entry with its own factor.

One more tool, for when a problem needs the global phase itself: the determinant \det\begin{pmatrix}p&q\\r&s\end{pmatrix} = ps - qr. Every rotation has determinant 1, a product has the product of the determinants, and multiplying a 2\times 2 matrix by c multiplies its determinant by c^2. So \det(cU) reveals c^2, and cmath.phase of it gives twice the angle of c (up to a whole turn):

import numpy as np, cmath, math
def ry(t):
    c, s = math.cos(t / 2), math.sin(t / 2)
    return np.array([[c, -s], [s, c]], dtype=complex)
def rz(t):
    return np.diag([cmath.exp(-0.5j * t), cmath.exp(0.5j * t)])
U = cmath.exp(0.3j) * rz(1.1) @ ry(0.7)
d = U[0, 0] * U[1, 1] - U[0, 1] * U[1, 0]
print(round(abs(d), 6), round(cmath.phase(d), 6), round(cmath.phase(complex(np.linalg.det(U))), 6))

1.0 0.6 0.6

The phase of the determinant is 0.6 = 2 \times 0.3. math.atan2(y, x) and cmath.phase are the safe ways to turn sizes and arrows into angles.

Common mistakes: dividing by an entry that is 0 (the top-left entry of X is). Assuming c is always one of \pm 1, \pm i. Comparing with ==.

Practice: The Unlabelled Pulse, Any Pulse from Three Turns

Circuits run left to right

A circuit on one qubit is a list of gates, applied left to right: the gate at index 0 acts first. Running it is a loop that replaces the state by each gate times the state:

import numpy as np, math
def show(a, digits=4):
    return (np.round(np.asarray(a, dtype=complex), digits) + 0).tolist()
def ry(t):
    c, s = math.cos(t / 2), math.sin(t / 2)
    return np.array([[c, -s], [s, c]], dtype=complex)
S = np.diag([1, 1j])
H = np.array([[1, 1], [1, -1]], dtype=complex) / math.sqrt(2)
state = np.array([1, 0], dtype=complex)
for name, U in (("ry(0.6)", ry(0.6)), ("s", S), ("h", H)):
    state = U @ state
    print("after", name.ljust(7), show(state))
print("odds:", np.round(np.abs(state) ** 2, 4).tolist())
other = H @ (ry(0.6) @ (S @ np.array([1, 0])))
print("S first instead:", show(other), np.round(np.abs(other) ** 2, 4).tolist())

after ry(0.6) [(0.9553+0j), (0.2955+0j)]
after s       [(0.9553+0j), 0.2955j]
after h       [(0.6755+0.209j), (0.6755-0.209j)]
odds: [0.5, 0.5]
S first instead: [(0.8845+0j), (0.4666+0j)] [0.7823, 0.2177]

The same three gates in another order give another state with other odds: gates, like matrices, do not commute in general. With plain lists, the update of one step is a, b = U[0][0] * a + U[0][1] * b, U[1][0] * a + U[1][1] * b, in one assignment; updating a first and then computing b from the new a is a classic bug.

Gate lists in problems name the gates, such as ["h", "s"] or ("ry", 0.6); a dictionary from names to matrices turns names into matrices.

Common mistakes: running the list backwards. Updating a before computing b. Rounding inside the loop, which lets errors grow over a long circuit.

Practice: Reading the Pulse Schedule

One matrix for a whole circuit

Running U_0, then U_1, then U_2 gives U_2(U_1(U_0\psi)), and by the product rule of Level 1 that is

So the matrix of a circuit is the product of its gates in reverse order: the last gate stands on the left. One matrix then runs the whole circuit on any number of starting states with one product each. Check it on the circuit of the last section:

import numpy as np, math
def show(a, digits=4):
    return (np.round(np.asarray(a, dtype=complex), digits) + 0).tolist()
def ry(t):
    c, s = math.cos(t / 2), math.sin(t / 2)
    return np.array([[c, -s], [s, c]], dtype=complex)
S = np.diag([1, 1j])
H = np.array([[1, 1], [1, -1]], dtype=complex) / math.sqrt(2)
M = H @ S @ ry(0.6)
print(show(M))
print("column 0:", show(M[:, 0]), " run on |0>:", show(H @ (S @ (ry(0.6) @ [1, 0]))))
wrong = ry(0.6) @ S @ H
print("written in run order:", show(wrong), np.allclose(M, wrong))

[[(0.6755+0.209j), (-0.209+0.6755j)], [(0.6755-0.209j), (-0.209-0.6755j)]]
column 0: [(0.6755+0.209j), (0.6755-0.209j)]  run on |0>: [(0.6755+0.209j), (0.6755-0.209j)]
written in run order: [[(0.6755-0.209j), (0.6755+0.209j)], [(0.209+0.6755j), (0.209-0.6755j)]] False

Column 0 of M is what the circuit does to |0\rangle, and column 1 what it does to |1\rangle: a quick self-test for any circuit matrix. For a long list, start from the identity and multiply each new gate in from the left, M = U @ M; M = M @ U builds the circuit backwards. An empty circuit has the identity as its matrix.

Common mistakes: multiplying in run order. M = M @ U in the loop. Starting from a zero matrix instead of the identity.

Practice: One Matrix for the Whole Schedule, Different Pulses, Same Gate?

Undoing a circuit

The inverse of a unitary matrix is its dagger, U^{-1} = U^\dagger: swap rows and columns and conjugate every entry. X, Y, Z and H are their own daggers, so each undoes itself. S^\dagger and T^\dagger turn back by -\pi/2 and -\pi/4, and every rotation is undone by the same rotation with the negated angle, R(\theta)^\dagger = R(-\theta).

A circuit is undone by undoing its last gate first, like taking off shoes before socks. In matrices, (U_2U_1U_0)^\dagger = U_0^\dagger U_1^\dagger U_2^\dagger, which as a gate list means: first U_2^\dagger, then U_1^\dagger, then U_0^\dagger. The daggers, in reverse order:

import numpy as np, math
def show(a, digits=4):
    return (np.round(np.asarray(a, dtype=complex), digits) + 0).tolist()
def ry(t):
    c, s = math.cos(t / 2), math.sin(t / 2)
    return np.array([[c, -s], [s, c]], dtype=complex)
S = np.diag([1, 1j])
H = np.array([[1, 1], [1, -1]], dtype=complex) / math.sqrt(2)
circuit = [ry(0.6), S, H]
state = np.array([1, 0], dtype=complex)
for U in circuit:
    state = U @ state
undo = [U.conj().T for U in reversed(circuit)]
back = state
for U in undo:
    back = U @ back
print(show(back))
print(np.allclose(ry(0.6).conj().T, ry(-0.6)), show(S.conj().T))
same_order = state
for U in [U.conj().T for U in circuit]:
    same_order = U @ same_order
print(show(same_order))

[(1+0j), 0j]
True [[(1+0j), 0j], [0j, -1j]]
[(0.4127-0.2177j), (0.7823+0.4127j)]

The reversed daggers bring the state back to |0\rangle. The daggers in the original order land somewhere else.

Common mistakes: undoing in the original order. Keeping S as its own inverse (it needs S^\dagger). Transposing without conjugating, which turns S into itself.

Practice: Rewinding the Pulse Schedule

Interference

H is its own inverse, so two Hadamards in a row do nothing: HH = I. Follow the amplitudes to see why. The first H splits |0\rangle into \frac{1}{\sqrt2}|0\rangle + \frac{1}{\sqrt2}|1\rangle. The second H sends each part on to both outcomes, so there are two paths to each final outcome, and the amplitudes of the paths add:

import math
r = 1 / math.sqrt(2)
H = [[r, r], [r, -r]]
first = [H[0][0], H[1][0]]                       # H on |0>
for out in (0, 1):
    paths = [H[out][mid] * first[mid] for mid in (0, 1)]
    print("to", out, "via 0 and via 1:", [round(x, 4) for x in paths], "sum:", round(sum(paths), 4))

to 0 via 0 and via 1: [0.5, 0.5] sum: 1.0
to 1 via 0 and via 1: [0.5, -0.5] sum: 0.0

Both paths to 0 carry +\frac12 and reinforce; the paths to 1 carry +\frac12 and -\frac12 and cancel. That is interference: amplitudes add as arrows, so they can cancel, which probabilities never do (Level 1's weights of a sum against the sum of the weights).

A gate between the two Hadamards turns one path against the other. Z flips the sign of the path through |1\rangle, so the cancellation moves to the other outcome and HZH = X. A phase gate that turns less, like S, leaves the cancellation incomplete:

import numpy as np, cmath, math
H = np.array([[1, 1], [1, -1]], dtype=complex) / math.sqrt(2)
for name, lam in (("I", 0.0), ("S", math.pi / 2), ("Z", math.pi)):
    M = H @ np.diag([1, cmath.exp(1j * lam)]) @ H
    out = M @ [1, 0]
    print(name, (np.round(M, 4) + 0).tolist(), np.round(np.abs(out) ** 2, 4).tolist())

I [[(1+0j), 0j], [0j, (1+0j)]] [1.0, 0.0]
S [[(0.5+0.5j), (0.5-0.5j)], [(0.5-0.5j), (0.5+0.5j)]] [0.5, 0.5]
Z [[0j, (1+0j)], [(1+0j), 0j]] [0.0, 1.0]

A global phase can never do this, because it turns every path by the same angle. Only a relative phase between paths changes what interferes, which is why the phase gates, which change no odds on their own, matter so much once an H follows them.

Common mistakes: adding the probabilities of the paths instead of their amplitudes. Expecting a phase gate alone to change the odds of 0 and 1.

Practice: Two Hadamards and a Phase Dial

Practice

These problems are required on this path at this level; the path adds one problem from each of the level's topics, picked for you and kept.

• Raw Amplitudes from the Chip Log
• Ten Shots at the Same Qubit
• Latitude and Longitude of a Qubit
• Stocking the Gate Drawer
• Reading the Pulse Schedule
• One Matrix for the Whole Schedule

Checklist

You can now:

• Use numpy for quantum work: np.array(..., dtype=complex), @, .conj().T, np.vdot, np.allclose, and print results through np.round(...).tolist().
• Normalise a list of amplitudes and turn a state into probabilities with the Born rule p_i = |a_i|^2, taking the size before squaring.
• Label outcomes as bit strings with qubit 0 as the leftmost bit, read them back with int(bits, 2), and pick out one qubit's bit with (i >> (n - 1 - q)) & 1.
• Sample shots reproducibly with random.Random(seed) and running totals of the probabilities.
• Tell a global phase (the same state) from a relative phase (a different state), and explain why counts of 0 and 1 show neither.
• Measure in the X and Y bases with |\langle e|\psi\rangle|^2, conjugating the first vector, and compute \langle Z\rangle = p_0 - p_1.
• Place a qubit on the Bloch sphere with \theta and \varphi, and compute differences of angles with cmath.phase and % (2 * math.pi).
• Write the gates I, X, Y, Z, H, S, T, their daggers, P(\lambda) and the rotations with their half angles, and decide whether two gates are equal up to a global phase.
• Run a circuit left to right, build its matrix as the product in reverse order, and undo it with the daggers in reverse order.
• Explain interference with paths: why HH = I and how a phase between two Hadamards changes the odds.