Real chips can read a qubit in the middle of a schedule and carry on, and the simulator you are building should too. A register of n qubits starts in |0\dots0\rangle (amplitude 1 at index 0) and the operations in ops are applied left to right. Each operation is a tuple (or list) whose first entry is a lower-case name:
("x", q),("z", q),("h", q),("s", q),("t", q): a gate on qubitq, with X=\begin{pmatrix}0&1\\1&0\end{pmatrix}, Z=\begin{pmatrix}1&0\\0&-1\end{pmatrix}, H=\tfrac{1}{\sqrt2}\begin{pmatrix}1&1\\1&-1\end{pmatrix}, S=\begin{pmatrix}1&0\\0&i\end{pmatrix}, T=\begin{pmatrix}1&0\\0&e^{i\pi/4}\end{pmatrix};("ry", q, theta),("rz", q, theta): R_y(\theta)=\begin{pmatrix}\cos\frac\theta2&-\sin\frac\theta2\\\sin\frac\theta2&\cos\frac\theta2\end{pmatrix} or R_z(\theta)=\begin{pmatrix}e^{-i\theta/2}&0\\0&e^{i\theta/2}\end{pmatrix} on qubitq;("cx", c, t): CNOT, which flips qubittin every basis state where qubitcreads 1 (c != t);("measure", q): read qubitq, as described below.
A gate G on qubit q acts as the full matrix I\otimes\dots\otimes G\otimes\dots\otimes I (G in place q, qubit 0's factor first). Bit order: basis state |b_0 b_1 \dots b_{n-1}\rangle has index int("b0b1...", 2), so qubit 0 is the leftmost bit, the most significant. In 3 qubits, |011\rangle is index 3.
Measuring qubit q. Let p0 and p1 be the totals of |a_i|^2 over the basis states where qubit q reads 0 and where it reads 1. Use rng = random.Random(seed), created once before the first operation. Each measure draws exactly one r = rng.random(), even when its outcome is certain, and picks the first outcome o (0, then 1) whose running total (p0, then p0 + p1) is greater than r (if rounding leaves no such outcome, pick the last outcome with nonzero probability). Then the register collapses: the amplitudes of basis states where qubit q does not read o become 0, and the others are divided by \sqrt{p_o}. The register keeps all n qubits.
Write run_with_readouts(n, ops, seed) that returns the pair (bits, state): bits is a string with one character "0" or "1" per measure, in the order they happened ("" if there are none), and state is the final list of 2^n amplitudes. Any global phase is accepted in the state: [1j, 0] is the same state as [1, 0].
Examples
Input: n = 1, ops = [("h", 0), ("measure", 0)], seed = 1
Output: ("0", [1, 0])
Explanation: p0 = p1 = 0.5 and the first draw is 0.1344 (rounded): 0.5 > 0.1344, so outcome 0.
Input: n = 2, ops = [("h", 0), ("cx", 0, 1), ("measure", 0), ("measure", 1)], seed = 0
Output: ("11", [0, 0, 0, 1])
Explanation: the first draw, 0.8444, gives qubit 0 the outcome 1 and the state collapses to |11>,
so qubit 1 reads 1 for certain (the second draw, 0.7580, is still used up).
Input: n = 3, ops = [("x", 0), ("measure", 0), ("h", 1), ("measure", 1)], seed = 1
Output: ("11", [0, 0, 0, 0, 0, 0, 1, 0])
Explanation: the first draw (0.1344) goes to the certain readout of qubit 0; the second
draw, 0.8474, is not below p0 = 0.5, so qubit 1 reads 1. The register ends in |110>, index 6.
Constraints
1 <= n <= 8,0 <= len(ops) <= 200- amplitudes are compared with a tolerance of
1e-6; the bits must match exactly
Goals
- Simulate a circuit that measures qubits part-way through, with reproducible random outcomes
- Collapse and renormalise the register after each readout and keep simulating
- See that measuring one qubit of an entangled pair fixes the other