Put Shor's algorithm together for a small odd number N and a base a with \gcd(a, N) = 1, on the simulator.
The quantum part. t counting qubits get H gates, and controlled powers of the gate |w\rangle\mapsto|a\,w \bmod N\rangle act on a register that starts at |1\rangle, so the state becomes \frac{1}{\sqrt{Q}}\sum_{x=0}^{Q-1}|x\rangle\,|a^x \bmod N\rangle with Q = 2^t (qubit 0 is the leftmost, most significant bit of x). Then the inverse QFT, |x\rangle \mapsto \frac{1}{\sqrt Q}\sum_y e^{-2\pi i\,xy/Q}|y\rangle, acts on the counting register, which is measured. Outcome y comes up with probability
where v runs over the distinct values of a^x \bmod N.
Sampling. Use rng = random.Random(seed). For each shot draw d = rng.random() and pick the first outcome y whose running total P(0) + ... + P(y) is greater than d (if rounding leaves no such y, pick the last outcome with P(y) > 10^{-12}). Every shot uses exactly one draw.
The classical part, for each shot's outcome y.
- Expand y/Q as a continued fraction (start from z = y/Q; terms c_0 = \lfloor z\rfloor, then z\leftarrow 1/(z - c_0), until z is whole; use exact integers). Its convergents h_k/k_k follow k_k = c_k k_{k-1} + k_{k-2} with k_{-1} = 0, k_{-2} = 1. The candidate period r is the denominator of the last convergent whose denominator is less than N (r = 1 for y = 0).
- The shot fails if a^r \bmod N \ne 1, or r is odd, or x = a^{r/2} \bmod N is 1 or N - 1.
- Otherwise it succeeds with the factors \gcd(x - 1, N) and \gcd(x + 1, N).
Write shor_run(N, a, t, seed, shots) that runs up to shots shots and returns (k, (g1, g2)) for the first successful shot k (counting from 1), with g1 <= g2 the two gcds; or None if every shot fails.
Examples
Input: N = 15, a = 7, t = 8, seed = 1, shots = 5
Output: (2, (3, 5))
Explanation: the shots give y = 0 (candidate 1, fails), then y = 192: 192/256 = 3/4, r = 4,
x = 7**2 mod 15 = 4, gcd(3, 15) = 3, gcd(5, 15) = 5.
Input: N = 21, a = 2, t = 9, seed = 11, shots = 5
Output: (3, (3, 7))
Explanation: y = 171 (candidate 3: 2**3 mod 21 = 8, fails), y = 256 (candidate 2, fails),
y = 427 (convergent 5/6, r = 6, x = 8: gcd(7, 21) = 7, gcd(9, 21) = 3).
Input: N = 15, a = 14, t = 8, seed = 2, shots = 6
Output: None
Explanation: 14 has order 2 and 14**1 mod 15 = 14 = N - 1, so no shot can succeed.
Constraints
15 <= N <= 63, odd and not a prime power;2 <= a < N,gcd(a, N) == 1;1 <= t <= 12;1 <= shots <= 20- Q is up to 4096, so a double loop over all
xandy(16 million terms) is too slow in the browser: use numpy's FFT, or group the terms
Goals
- Compute the counting register's outcome distribution of order finding
- Sample outcomes with a prescribed seeded procedure
- Chain the classical steps (period candidate, check, gcds) and repeat until a shot succeeds