Skip to content

Shor Algorithm Explained — Quantum Factoring Step by Step

DodaTech Updated 2026-06-23 11 min read

In this tutorial, you'll learn about Shor Algorithm Explained. We cover key concepts, practical examples, and best practices to help you understand and apply this topic effectively.

Shor's algorithm factors large integers in polynomial time using quantum period finding, breaking RSA and ECC encryption that classical computers cannot crack in any feasible timeframe.

What You'll Learn

By the end of this tutorial, you will understand how Shor's algorithm reduces factoring to period finding, how quantum Fourier transform extracts the period, how modular exponentiation is implemented, and what resources are needed to factor RSA-2048.

Why It Matters

Shor's algorithm threatens all public-key cryptography currently used on the internet (RSA, ECC, Diffie-Hellman). A large-scale quantum computer running Shor's algorithm could decrypt HTTPS traffic, forge digital signatures, and break Blockchain security. This is why NIST is standardizing post-Quantum Cryptography.

Real-World Use

If a fault-tolerant quantum computer with ~20 million qubits is built, Shor's algorithm could factor a 2048-bit RSA number in approximately 8 hours. Governments and standards bodies are actively preparing for this scenario through migration to post-quantum cryptographic algorithms.

Learning Path

flowchart LR
  A[Quantum Algorithms] --> B[Shor Algorithm]
  B --> C[Quantum Error Correction]
  B --> D[Quantum Cryptography]
  B --> E{You Are Here}
  style E fill:#f90,color:#fff
ℹ️ Info

Prerequisites: Understand QFT, phase estimation, and basic Python. Familiarity with RSA cryptography helps.

How Shor's Algorithm Works

Shor's algorithm has two parts: a classical reduction from factoring to period finding, and a quantum period-finding subroutine.

The Classical Reduction

Given an integer N to factor:

  1. Choose a random a coprime to N
  2. Find the period r of f(x) = a^x mod N
  3. If r is even and a^(r/2) ≠ ±1 mod N, then gcd(a^(r/2) ± 1, N) divides N
# classical_reduction.py
import math
import random

def classical_reduction(N, a=None):
    """Classical reduction from factoring to period finding."""
    if a is None:
        a = random.randrange(2, N - 1)

    g = math.gcd(a, N)
    if g != 1:
        return g  # Lucky guess

    print(f"Attempting to factor N={N} with a={a}")

    # Simulate finding the period r (normally done by quantum computer)
    # For small N, we can compute it classically
    r = find_period_classical(a, N)

    if r is None:
        print(f"  No period found (a may not be coprime)")
        return None

    print(f"  Found period r={r}")

    if r % 2 != 0:
        print(f"  r is odd, try different a")
        return None

    x = pow(a, r // 2, N)
    if x == N - 1 or x == 1:
        print(f"  a^(r/2) ≡ ±1 mod N, try different a")
        return None

    factor1 = math.gcd(x - 1, N)
    factor2 = math.gcd(x + 1, N)
    print(f"  Factors: {factor1} × {factor2}")
    return factor1, factor2

def find_period_classical(a, N, max_r=1000):
    """Classically find the period of a^x mod N (for simulation only)."""
    seen = {}
    val = 1
    for x in range(max_r):
        if val in seen:
            return x - seen[val]
        seen[val] = x
        val = (val * a) % N
    return None

print("=== Classical Reduction (Shor) ===")
for N in [15, 21, 35]:
    result = classical_reduction(N, a=random.randrange(2, N-1))
    print()

Expected output (varies with random a):

=== Classical Reduction (Shor) ===
Attempting to factor N=15 with a=7
  Found period r=4
  Factors: 3 × 5

Attempting to factor N=21 with a=8
  Found period r=2
  Factors: 7 × 3

Quantum Period Finding

The quantum part of Shor's algorithm finds the period r using two quantum registers and the Quantum Fourier Transform.

# quantum_period_finding.py
import numpy as np

class QuantumPeriodFinder:
    """Simulate the quantum period-finding subroutine of Shor's algorithm."""

    def __init__(self, n_counting=8):
        self.n = n_counting
        self.dim = 2 ** n_counting

    def qft(self, state):
        """Apply Quantum Fourier Transform."""
        n = int(np.log2(len(state)))
        dim = len(state)
        omega = np.exp(2 * np.pi * 1j / dim)
        QFT = np.zeros((dim, dim), dtype=complex)
        for i in range(dim):
            for j in range(dim):
                QFT[i, j] = omega ** (i * j) / np.sqrt(dim)
        return QFT @ state

    def inverse_qft(self, state):
        """Apply inverse QFT."""
        n = int(np.log2(len(state)))
        dim = len(state)
        omega = np.exp(-2 * np.pi * 1j / dim)
        IQFT = np.zeros((dim, dim), dtype=complex)
        for i in range(dim):
            for j in range(dim):
                IQFT[i, j] = omega ** (i * j) / np.sqrt(dim)
        return IQFT @ state

    def find_period(self, a, N):
        """
        Simulate quantum period finding.
        
        For simulation purposes, we compute the unitary
        U|y⟩ = |ay mod N⟩ and find its eigenvalues.
        """
        dim_N = N

        # Construct the modular multiplication unitary
        U = np.zeros((dim_N, dim_N), dtype=complex)
        for y in range(dim_N):
            U[(a * y) % dim_N, y] = 1.0

        # Initialize in equal superposition of first register
        counting_state = np.ones(self.dim, dtype=complex) / np.sqrt(self.dim)

        # Apply controlled-U^(2^k) operations
        # (simplified: directly apply QFT-based phase estimation)
        target_state = np.zeros(dim_N, dtype=complex)
        target_state[1] = 1.0  # Start in |1⟩

        # Build the full state
        full_dim = self.dim * dim_N
        full_state = np.kron(counting_state, target_state)

        # Apply controlled-U^(2^k) for each counting qubit (simplified)
        # In a real simulation, this would apply the modular exponentiation

        # For demonstration, we simulate the result:
        # The period r appears as peaks in the Fourier basis
        # Manually construct the expected output

        # Simulate measurement outcome for a known period
        # For a=7, N=15, the period is 4
        r = find_period_classical(a, N)

        # In the quantum algorithm, we'd measure |j/r⟩ in the Fourier basis
        peaks = self._simulate_peaks(r)

        return peaks, r

    def _simulate_peaks(self, r):
        """Simulate the peaks that appear after inverse QFT."""
        # The probability distribution peaks at multiples of 2^n / r
        peak_positions = [k * self.dim // r for k in range(r)]
        probs = np.zeros(self.dim)
        for p in peak_positions:
            probs[p] = 1.0 / len(peak_positions)
        return probs

# Helper function for period finding
def find_period_classical(a, N, max_r=1000):
    """Classical period finding for verification."""
    seen = {}
    val = 1
    for x in range(max_r):
        if val in seen:
            return x - seen[val]
        seen[val] = x
        val = (val * a) % N
    return None

print("=== Quantum Period Finding ===")
for a, N in [(7, 15), (2, 21), (4, 21)]:
    qpf = QuantumPeriodFinder(n_counting=8)
    r = find_period_classical(a, N)
    if r:
        print(f"a={a}, N={N}: period r={r}")
        print(f"  Continuing fraction: {qpf.dim}/{r} = {qpf.dim / r:.2f}")
        print(f"  Peak positions at multiples of 2^{qpf.n}/{r}")

Expected output:

=== Quantum Period Finding ===
a=7, N=15: period r=4
  Continuing fraction: 256/4 = 64.00
  Peak positions at multiples of 2^8/4
a=2, N=21: period r=6
  Continuing fraction: 256/6 = 42.67
  Peak positions at multiples of 2^8/6

Modular Exponentiation

The most resource-intensive part of Shor's algorithm is modular exponentiation: computing a^x mod N for a superposition of x values.

# modular_exponentiation.py
import numpy as np

def classical_mod_exp(a, x, N):
    """Compute a^x mod N classically (for verification)."""
    return pow(a, x, N)

class ModularExponentiationSimulator:
    """Simulate the modular exponentiation circuit."""

    def __init__(self, a, N, n_bits):
        self.a = a
        self.N = N
        self.n = n_bits
        self.dim = 2 ** n_bits

    def unitary(self):
        """Build the U matrix where U|x⟩ = |ax mod N⟩."""
        U = np.zeros((self.dim, self.dim), dtype=complex)
        for x in range(min(self.dim, self.N)):
            y = (self.a * x) % self.N
            U[y, x] = 1.0
        # Identity for |x⟩ where x >= N
        for x in range(self.N, self.dim):
            U[x, x] = 1.0
        return U

    def controlled_unitary(self, control_bit):
        """Simulate controlled-U operation."""
        U = self.unitary()
        # For controlled version, apply U when control is 1
        dim_full = self.dim * 2
        CU = np.eye(dim_full, dtype=complex)
        CU[self.dim:, self.dim:] = U
        return CU

    def verify(self):
        """Verify modular exponentiation."""
        U = self.unitary()
        print(f"=== Modular Exponentiation (a={self.a}, N={self.N}) ===")
        for x in range(min(5, self.N)):
            state = np.zeros(self.dim, dtype=complex)
            state[x] = 1.0
            result = U @ state
            y = np.argmax(np.abs(result))
            expected = classical_mod_exp(self.a, 1, self.N) if x == 1 else 0
            print(f"  U|{x}⟩ = |{y}⟩ (expected: |{(self.a * x) % self.N}⟩)")

# Test
ms = ModularExponentiationSimulator(a=7, N=15, n_bits=8)
ms.verify()

# Resource estimation
print(f"\n=== Resource Estimation ===")
print(f"To factor N=15: need ~8 qubits, ~100 gates")
print(f"To factor RSA-2048: need ~20M physical qubits (with error correction)")
print(f"To factor RSA-2048: ~2^13 logical qubits, ~10^11 T gates")
print(f"Estimated runtime: ~8 hours on fault-tolerant hardware")

Expected output:

=== Modular Exponentiation (a=7, N=15) ===
  U|0⟩ = |0⟩ (expected: |0⟩)
  U|1⟩ = |7⟩ (expected: |7⟩)
  U|2⟩ = |14⟩ (expected: |14⟩)
  U|3⟩ = |6⟩ (expected: |6⟩)
  U|4⟩ = |13⟩ (expected: |13⟩)

=== Resource Estimation ===
To factor N=15: need ~8 qubits, ~100 gates
To factor RSA-2048: need ~20M physical qubits (with error correction)
To factor RSA-2048: ~2^13 logical qubits, ~10^11 T gates
Estimated runtime: ~8 hours on fault-tolerant hardware

Complete Shor Simulation for N=15

Let us put it all together to factor 15 — the canonical example.

# shor_15.py
import math
import random

def shor_factor_15():
    """Complete Shor's algorithm simulation for N=15."""
    N = 15
    print("=== Shor's Algorithm: Factor 15 ===")

    attempts = 0
    while True:
        attempts += 1
        a = random.randrange(2, N - 1)
        print(f"\nAttempt {attempts}: a={a}")

        # Check gcd
        g = math.gcd(a, N)
        if g != 1:
            print(f"  Lucky! gcd({a}, {N}) = {g}")
            return g, N // g

        # Find period (quantum part simulated classically)
        r = None
        vals = {}
        val = 1
        for x in range(N * 2):
            if val in vals:
                r = x - vals[val]
                break
            vals[val] = x
            val = (val * a) % N

        if r is None:
            print(f"  No period found")
            continue

        print(f"  Period r = {r}")

        if r % 2 != 0:
            print(f"  r is odd, retry")
            continue

        x = pow(a, r // 2, N)
        if x == N - 1:
            print(f"  a^(r/2) ≡ -1 mod N, retry")
            continue

        factor1 = math.gcd(x - 1, N)
        factor2 = math.gcd(x + 1, N)

        if factor1 * factor2 == N:
            print(f"  Success! {N} = {factor1} × {factor2}")
            return factor1, factor2

factor1, factor2 = shor_factor_15()
print(f"\nFinal result: {factor1} × {factor2} = {factor1 * factor2}")

Expected output:

=== Shor's Algorithm: Factor 15 ===
Attempt 1: a=7
  Period r = 4
  Success! 15 = 3 × 5

Final result: 3 × 5 = 15

Common Mistakes

1. Thinking Shor's Algorithm is Only for Factoring

Shor's algorithm solves the hidden subgroup problem for abelian groups. This includes factoring (RSA), discrete logarithms (ECC, Diffie-Hellman), and certain problems in graph theory.

2. Underestimating Resource Requirements

Factoring RSA-2048 requires approximately 20 million physical qubits with surface code error correction, not a few thousand. The T-gate count is ~10^11, requiring days of computation.

3. Confusing Physical and Logical Qubits

When people say "Shor needs thousands of qubits," they mean logical qubits. Each logical qubit requires 100-1000 physical qubits for error correction, dramatically increasing the total.

4. Ignoring the Classical Post-Processing

After the quantum computation, Shor requires continued fraction expansion to extract r from the measured peaks, plus gcd computations. This classical part is straightforward but essential.

5. Assuming All Periods Work

Not every choice of a leads to a useful period. About 40% of random a values fail (odd period, or a^(r/2) ≡ -1 mod N). The algorithm must be repeated with different a values.

Practice Questions

1. How does Shor's algorithm reduce factoring to period finding?

By choosing a random a coprime to N and finding the period r of f(x) = a^x mod N. If r is even and a^(r/2) ≠ ±1 mod N, then gcd(a^(r/2) ± 1, N) yields a non-trivial factor of N.

2. What is the role of the Quantum Fourier Transform in Shor's algorithm?

QFT extracts the period r from the Superposition of modular exponentiation outputs. It creates interference peaks at multiples of 2^n / r, which are measured to determine r.

3. Why does Shor's algorithm threaten RSA?

RSA security relies on the difficulty of factoring large numbers. Shor's algorithm factors these numbers in polynomial time (O((log N)^3)), making RSA completely insecure against quantum adversaries.

4. How many qubits are needed to run Shor's algorithm for RSA-2048?

Approximately 2n logical qubits (where n = 2048), translating to about 20 million physical qubits with surface code error correction at current error rates.

Challenge: Factor N=21

Implement a complete simulation of Shor's algorithm for N=21. Determine which choices of a yield useful periods and which fail. Extend your implementation to try multiple a values automatically until a factorization succeeds.

Real-World Task: NIST Post-Quantum Cryptography Survey

Research the NIST post-Quantum Cryptography standardization finalists (CRYSTALS-Kyber, CRYSTALS-Dilithium, Falcon, SPHINCS+). For each algorithm, explain: what mathematical problem it relies on, why it resists Shor's algorithm, and its key sizes compared to RSA. Write a brief assessment of which should replace RSA in your systems.

FAQ

What is Shor's algorithm?

Shor's algorithm is a quantum algorithm that factors large integers in polynomial time, discovered by Peter Shor in 1994. It poses a direct threat to RSA and ECC cryptography.

How fast is Shor's algorithm?

Shor's algorithm runs in O((log N)^3) time, compared to the best classical algorithm (GNFS) which runs in sub-exponential time. For a 2048-bit number, this is the difference between 8 hours and billions of years.

Does Shor's algorithm affect symmetric cryptography?

No. Shor's algorithm only breaks public-key cryptography (RSA, ECC, Diffie-Hellman). Symmetric encryption (AES) and hash functions are only affected by Grover's algorithm, which provides a quadratic speedup.

Has Shor's algorithm been run on real hardware?

Yes, but only for tiny numbers like 15 and 21. Running Shor for cryptographically relevant sizes (2048-bit RSA) requires a fault-tolerant quantum computer not yet available.

When will Shor's algorithm become a practical threat?

Most estimates place large-scale fault-tolerant quantum computers capable of running Shor in the 2030s or later, but the timeline is uncertain.

What's Next

Quantum Error Correction Explained
Quantum Algorithms Overview
Quantum Cryptography Guide

You now understand how Shor's algorithm threatens modern cryptography. Next, you will learn how error correction enables fault-tolerant quantum computation.

Built by the developers of Doda Browser, DodaZIP, and Durga Antivirus Pro. Last updated: 2026-06-23.

Built by the developers of DodaTech

Doda Browser, DodaZIP & Durga Antivirus Pro