Case Study 1: The Rotation That Was Half as Large as It Should Be

"Every physicist has written this bug. The good ones recognize it in under a minute."

Executive Summary

A state-preparation routine is supposed to produce a qubit with a 25% chance of measuring 1. It produces 6.7% instead. The code is four lines long, there is no noise anywhere in the pipeline, and the routine is deterministic.

This case study is a complete debugging session on that failure — chosen because the bug is the single most common arithmetic error in quantum programming, because the method generalizes to problems far larger than four lines, and because it demonstrates the tool you will lean on for the rest of the book: a small circuit, simulated exactly, checked against a number you computed by hand.

Skills applied: Statevector inspection (§3.3); Bloch geometry (§3.4); rotation gates and the half-angle convention (§3.6); gate-to-matrix conversion (§3.8).

The Failure

The requirement: a routine prepare(p) that puts a qubit in a state measuring 1 with probability p. Here is the implementation, and it looks fine:

import numpy as np
from qiskit import QuantumCircuit

def prepare(p: float) -> QuantumCircuit:
    """Prepare a state with P(1) = p."""
    theta = np.arcsin(np.sqrt(p))     # sin(theta) = sqrt(p)  =>  sin^2(theta) = p
    qc = QuantumCircuit(1)
    qc.ry(theta, 0)
    return qc

The reasoning is stated in the comment and it is not obviously wrong: we want $\sin^2\theta = p$, so $\theta = \arcsin\sqrt{p}$.

The test:

from qiskit.quantum_info import Statevector

for p in (0.1, 0.25, 0.5, 0.75, 0.9):
    got = Statevector(prepare(p)).probabilities()[1]
    print(f"  requested {p:.2f}   got {got:.4f}   error {got - p:+.4f}")
  requested 0.10   got 0.0257   error -0.0743
  requested 0.25   got 0.0670   error -0.1830
  requested 0.50   got 0.1464   error -0.3536
  requested 0.75   got 0.2500   error -0.5000
  requested 0.90   got 0.3419   error -0.5581

Every value is wrong, and wrong by a different amount. That last detail matters: a constant offset suggests a bug in one place, while an error that grows with the input suggests a wrong functional form.

Step 1: Resist the Urge to Fix

The instinct is to start adjusting: try $2\theta$, try $\sqrt{p}$ without the arcsin, try a different gate. Two or three of those will get closer and one of them might even pass the test for the values you happened to check.

Do not. Guess-and-check on a four-line function is a bad habit that becomes unworkable on a forty-line one, and a fix you cannot explain is a bug you have relocated.

Instead: find out what the circuit actually does.

Step 2: Look at the State

qc = prepare(0.5)
sv = Statevector(qc)
theta = np.arcsin(np.sqrt(0.5))

print(f"theta      = {theta:.6f}  ({theta/np.pi:.4f} pi)")
print(f"amplitudes = {sv.data.round(6)}")
print(f"P          = {sv.probabilities().round(6)}")
theta      = 0.785398  (0.2500 pi)
amplitudes = [0.92388+0.j 0.38268+0.j]
P          = [0.853553 0.146447]

Now compare the amplitudes to what we expected. We asked for $P(1) = 0.5$, which needs an amplitude of $\sqrt{0.5} = 0.7071$. We got 0.38268.

And 0.38268 is $\sin(\pi/8)$, not $\sin(\pi/4)$.

There it is. The gate applied half the angle we passed it.

Step 3: Confirm Against the Matrix

Never conclude from one number. Read the definition:

$$R_y(\theta) = \begin{pmatrix} \cos\frac{\theta}{2} & -\sin\frac{\theta}{2} \\[4pt] \sin\frac{\theta}{2} & \cos\frac{\theta}{2}\end{pmatrix}$$

Applied to $|0\rangle = (1, 0)^\top$, this gives amplitudes $\bigl(\cos\frac{\theta}{2},\ \sin\frac{\theta}{2}\bigr)$, so

$$P(1) = \sin^2\!\frac{\theta}{2}, \quad \textbf{not} \quad \sin^2\theta.$$

Confirm it in code rather than trusting the reading:

from qiskit.quantum_info import Operator
from qiskit.circuit.library import RYGate

theta = np.pi / 4
print(Operator(RYGate(theta)).data.round(6))
print(f"cos(theta/2) = {np.cos(theta/2):.6f}   cos(theta) = {np.cos(theta):.6f}")
[[ 0.92388+0.j -0.38268+0.j]
 [ 0.38268+0.j  0.92388+0.j]]
cos(theta/2) = 0.923880   cos(theta) = 0.707107

The matrix entry is 0.92388, which is $\cos(\theta/2)$. Confirmed. The half-angle is in the gate definition, and it was in the documentation the whole time.

Step 4: The Fix, Derived

Now we can fix it by derivation rather than by guessing.

We want $P(1) = \sin^2\frac{\theta}{2} = p$. Therefore $\frac{\theta}{2} = \arcsin\sqrt{p}$, so $\theta = 2\arcsin\sqrt{p}$.

def prepare(p: float) -> QuantumCircuit:
    """Prepare a state with P(1) = p.

    Ry's matrix uses the HALF angle: applying Ry(theta) to |0> gives amplitudes
    (cos(theta/2), sin(theta/2)), so P(1) = sin^2(theta/2). To hit a target p we
    therefore need theta = 2*arcsin(sqrt(p)).
    """
    if not 0.0 <= p <= 1.0:
        raise ValueError(f"p must be in [0, 1], got {p}")
    theta = 2 * np.arcsin(np.sqrt(p))
    qc = QuantumCircuit(1)
    qc.ry(theta, 0)
    return qc
  requested 0.10   got 0.1000   error +0.0000
  requested 0.25   got 0.2500   error -0.0000
  requested 0.50   got 0.5000   error +0.0000
  requested 0.75   got 0.7500   error +0.0000
  requested 0.90   got 0.9000   error +0.0000

Note that the fix came with a docstring explaining why, and a bounds check. The comment is the actual deliverable: the next person to read this function — including you in three months — will not re-derive the half-angle from scratch.

Step 5: Why Half? (The Part Worth Understanding)

The factor is not an arbitrary convention. It is physics.

A qubit is a spin-½ system, and spin-½ systems have the property that rotating them through $2\pi$ does not return them to their original state — it returns them to $-1$ times it. You need $4\pi$ to get back exactly.

Verify it:

for turns in (0, 0.5, 1, 1.5, 2):
    qc = QuantumCircuit(1)
    qc.ry(turns * 2 * np.pi, 0)
    sv = Statevector(qc)
    print(f"  {turns:>4.1f} full turns: state {sv.data.round(4)}  P = {sv.probabilities().round(4)}")
   0.0 full turns: state [1.+0.j 0.+0.j]  P = [1. 0.]
   0.5 full turns: state [0.+0.j 1.+0.j]  P = [0. 1.]
   1.0 full turns: state [-1.+0.j  0.+0.j]  P = [1. 0.]
   1.5 full turns: state [-0.+0.j -1.+0.j]  P = [0. 1.]
   2.0 full turns: state [ 1.+0.j -0.+0.j]  P = [1. 0.]

After one full turn the state is $-|0\rangle$: correct probabilities, wrong sign. After two full turns it is genuinely back.

The Bloch sphere makes this visible. The Bloch sphere shows the state's direction, and one $2\pi$ rotation of the Bloch vector corresponds to only $\pi$ of "state rotation." The sphere is a two-to-one picture of the actual state space, and the half-angle in the gate matrix is the bookkeeping that connects them.

And here is why it usually does not bite: after one full turn the state is $-|0\rangle$, which differs from $|0\rangle$ by a global phase and is therefore unobservable (§3.7). The factor of two hides — until you either compute an angle from a target probability, as here, or control the operation, at which point the sign becomes relative and visible.

Analysis: Why This Bug Is So Common

Three reinforcing reasons.

The convention is invisible at the call site. qc.ry(theta, 0) gives no hint that theta is halved internally. Nothing about the API surfaces it.

The wrong version is plausible and self-consistent. "I want $\sin^2\theta = p$, so $\theta = \arcsin\sqrt{p}$" is correct reasoning applied to a wrong premise, and it produces smooth, monotonic, plausible-looking output. Nothing crashes. Nothing warns.

It passes the endpoints. At $p = 0$ and $p = 1$ both versions give the same answer. A test suite that only checks the boundary cases — which is exactly what a hurried test suite checks — passes completely.

That last point is a general lesson about testing quantum code, and Chapter 27 returns to it: test the interior, not the endpoints. Endpoints are where the most bugs hide, because they are where the most functional forms agree.

The Method, Extracted

The debugging session was four steps and none of them involved changing the code:

  1. Look at the actual state, not the failing output. Statevector gave the amplitudes.
  2. Compare against a hand-computed expectation. We knew $\sqrt{0.5} = 0.7071$ and saw 0.38268.
  3. Confirm against the definition, in code. The matrix entry was $\cos(\theta/2)$.
  4. Derive the fix, and write down the derivation.

This scales. On a forty-qubit variational circuit you cannot inspect the statevector — but you can shrink to two qubits, where you can, and the bug will almost always survive the shrinking. Chapter 26 formalizes this as the first technique in its taxonomy.

Lessons

  1. Rotation gates use the half angle. $R_y(\theta)|0\rangle$ has $P(1) = \sin^2(\theta/2)$. When a rotation does half or twice what you expected, check this before anything else.
  2. An error that grows with the input indicates a wrong functional form, not a wrong constant. Read the shape of the error, not just its size.
  3. Do not guess-and-check. A fix you cannot derive is a bug you have relocated.
  4. Statevector before print. The amplitudes told the whole story in one line; the probabilities alone would not have.
  5. Confirm against the definition in code, not from memory. Operator(RYGate(θ)) takes three seconds and removes all doubt.
  6. Test the interior, not the endpoints. At $p = 0$ and $p = 1$ the buggy and correct versions agree exactly.
  7. The factor of two is spin-½ physics, not a convention someone chose. Understanding why makes it memorable in a way that memorizing it does not.
  8. Write the derivation into the docstring. That is the actual fix; the code change is a consequence.

Questions

  1. Derive the equivalent correction for $R_x$. Does the same factor of two apply? Verify.

  2. Write prepare_rz(phi) that applies a relative phase $\phi$ to a state already in $|+\rangle$. What half-angle trap does $R_z$ contain, and how would you detect it — given that the phase is invisible to a computational-basis measurement? (Hint: §3.5.)

  3. The buggy and correct versions agree at $p=0$ and $p=1$. Find a third value of $p$ where they agree, or prove none exists.

  4. Write a property-based test for prepare(p) that would have caught this bug, using 20 random values of $p$ in $(0,1)$ and a tolerance. What tolerance is appropriate for a statevector simulation, and why is it not the same tolerance you would use for a 1,024-shot sampled run?

  5. Suppose the bug had reached hardware and the routine was validated by sampling 1,000 shots. At $p = 0.25$ the true output is 7.3%. How many shots would you need for the discrepancy to be unambiguous? (Use the $1/\sqrt{N}$ result from Chapter 1.)

  6. prepare() now raises on $p$ outside $[0,1]$. Is arcsin of a value slightly above 1 — from floating-point error in an upstream calculation — a realistic risk? Write the defensive version and say what it should do.

  7. Hardest. Generalize prepare to two qubits: prepare2(p00, p01, p10, p11) producing a specified distribution over four outcomes, with all-real amplitudes. How many rotation angles do you need? Does an arbitrary distribution always have a solution with only $R_y$ gates and CNOTs, and if so, why? (Come back to this after Chapter 4 if it does not yield now.)