Shor's Algorithm in Qiskit: A Practical Tutorial

Shor's algorithm gets referenced constantly on this site, in our Q-Day piece and our QFT tutorial, because it's the algorithm that makes RSA a finite-lifetime cryptosystem. It rarely gets built. This tutorial builds it, using the standard small textbook example (factoring 15) so every piece stays visible instead of disappearing into a circuit too large to reason about.

A note on the code below: it targets Qiskit 2.5.0 and is written to run as shown. Where a step involves probabilistic measurement, we describe the theoretical distribution the math predicts rather than presenting a specific captured run as if it were the only possible outcome.

What Shor's algorithm reduces to

Factoring N reduces to a different problem: given some a coprime to N, find the smallest r such that a^r ≡ 1 (mod N). This is order finding, and once you have r, classical post-processing (a GCD computation) usually hands you a real factor. The 1994 paper is short precisely because almost all of it is number theory. The quantum part exists only to find r fast, since finding it classically takes exponential time.

For N = 15, pick a = 7 (coprime to 15, and a well-known example because the period is small enough to see clearly). The sequence 7¹, 7², 7³, 7⁴ mod 15 is 7, 4, 13, 1, so r = 4.

Step 1: Modular exponentiation as a quantum operator

The circuit needs a unitary U that maps |x⟩ to |a·x mod N⟩ for a work register, controlled by a counting register. Building a general modular exponentiation circuit is real engineering (see Qiskit's own algorithm library for the full construction), so this tutorial uses the standard shortcut for the textbook case: since 7 mod 15 has order 4, the controlled-U operations reduce to a small, explicit permutation matrix on 4 qubits (one for the counting register bit being tested, plus a 4-state work register encoded in 2 qubits for the cycle 1 → 7 → 4 → 13 → 1).

import numpy as np
from qiskit import QuantumCircuit, transpile
from qiskit.circuit.library import QFTGate
from qiskit_aer import AerSimulator
from fractions import Fraction

N, a = 15, 7

def c_amod15(power):
    """Controlled multiplication by 7^power mod 15, built from the known 4-cycle."""
    qc = QuantumCircuit(4)
    for _ in range(power):
        qc.swap(2, 3)
        qc.swap(1, 2)
        qc.swap(0, 1)
        for q in range(4):
            qc.x(q)
    U = qc.to_gate()
    U.name = f"7^{power} mod 15"
    return U.control()

This is the standard construction used in introductory Shor's tutorials for N = 15 specifically. It is not a general modular exponentiation circuit, and it would not extend to a larger N without a real arithmetic circuit behind it, worth being upfront about rather than implying this scales.

Step 2: The counting register and inverse QFT

Order finding is quantum phase estimation applied to the operator "multiply by a mod N." Three counting qubits are enough to resolve a period of 4 (since 2³ = 8 gives more than enough resolution for r = 4):

n_count = 3

qc = QuantumCircuit(n_count + 4, n_count)

for q in range(n_count):
    qc.h(q)

qc.x(n_count)  # work register starts in |1>

for q in range(n_count):
    qc.append(c_amod15(2 ** q), [q] + list(range(n_count, n_count + 4)))

qc.append(QFTGate(n_count).inverse(), range(n_count))
qc.measure(range(n_count), range(n_count))

The inverse QFT is doing exactly the job described in our QFT tutorial: pulling a phase that encodes the period out into a measurable computational-basis number.

Step 3: What the counting register measures

With 3 counting qubits and a true period of 4, quantum phase estimation concentrates the measurement outcomes on multiples of 8/4 = 2: the states 0, 2, 4, and 6, each with a theoretical probability of 1/4 in the ideal noiseless case. That's the mathematical prediction from the phase estimation formula, not a specific simulated run, and it's the number to check your own execution against if you run this yourself.

Expected measurement distribution (ideal, 3 counting qubits):
000 (0) -> 1/4
010 (2) -> 1/4
100 (4) -> 1/4
110 (6) -> 1/4

Step 4: Continued fractions turn a measurement into r

Each measured integer m encodes a phase s/r as m/2^n_count. Continued fractions recovers the fraction in lowest terms:

def phase_to_period(measured, n_count, N):
    phase = measured / (2 ** n_count)
    frac = Fraction(phase).limit_denominator(N)
    return frac.denominator

for measured in [0, 2, 4, 6]:
    r_candidate = phase_to_period(measured, n_count, N)
    print(f"measured={measured}: candidate r = {r_candidate}")

Feeding in 2 or 6 recovers r = 4 directly (measuring 0 gives no information and has to be discarded and retried, a known failure mode of the algorithm, not a bug in this implementation).

Step 5: From r to an actual factor

With r = 4 (even, which the algorithm requires) and a = 7:

r = 4
guess1 = np.gcd(a ** (r // 2) - 1, N)
guess2 = np.gcd(a ** (r // 2) + 1, N)
print(f"gcd(7^2 - 1, 15) = {guess1}")
print(f"gcd(7^2 + 1, 15) = {guess2}")

7² = 49. gcd(48, 15) = 3, and gcd(50, 15) = 5. Both nontrivial factors of 15, found without ever trying to divide 15 by anything directly.

Why this doesn't threaten RSA-2048 today

Everything above ran on 7 qubits for a 4-bit number. Factoring an RSA-2048 modulus needs a register sized to the number of bits in N, meaning thousands of logical qubits, and the logical qubit overhead to build even one of those on top of noisy physical hardware. Our Q-Day piece covers the published resource estimates for that gap in detail, and why they disagree with each other by orders of magnitude depending on the assumed error-correction overhead. The algorithm in this tutorial is exactly the algorithm that would eventually run at that scale. The distance between 15 and a 2048-bit modulus is the entire reason post-quantum cryptography migration has a runway measured in years, not months.

Try this next

• Run the counting-register measurement on AerSimulator yourself and compare your shot counts against the theoretical 1/4, 1/4, 1/4, 1/4 split above.
• Swap a = 7 for a = 2 or a = 4 against N = 15 and work out the new order by hand first, then check the circuit agrees.
• Read Qiskit's own textbook chapter on Shor's algorithm for the general (non-shortcut) modular exponentiation construction, and compare gate counts against the toy version here.