Number Theory for Programmers — Divisibility, GCD and Euclid's Algorithm, Modular Arithmetic, Primes and a Tiny RSA

Number theory used to be called the "purest" maths — useless in practice. Then it became the foundation of cryptography: HTTPS, UPI and digital signatures all depend on it. It's also everywhere in competitive programming ("print the answer modulo 10⁹ + 7").

1. Divisibility and primes

• a divides b (a | b) if b = a·k for some integer k.
• A prime has exactly two divisors: 1 and itself. 2 is the only even prime.
• Fundamental Theorem of Arithmetic: every integer > 1 is a product of primes in exactly one way (360 = 2³ · 3² · 5).
• To test if n is prime, trying divisors up to √n is enough. To find all primes up to N, use the Sieve of Eratosthenes.

def sieve(n):
    is_prime = [True] * (n + 1)
    is_prime[0] = is_prime[1] = False
    for p in range(2, int(n ** 0.5) + 1):
        if is_prime[p]:
            for multiple in range(p * p, n + 1, p):
                is_prime[multiple] = False
    return [i for i, ok in enumerate(is_prime) if ok]

def factorise(n):
    f, p = {}, 2
    while p * p <= n:
        while n % p == 0:
            f[p] = f.get(p, 0) + 1
            n //= p
        p += 1
    if n > 1:
        f[n] = f.get(n, 0) + 1
    return f

print("primes up to 50:", sieve(50))
print("360 =", factorise(360), "| 9973 =", factorise(9973))

primes up to 50: [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47]
360 = {2: 3, 3: 2, 5: 1} | 9973 = {9973: 1}

2. GCD and Euclid's algorithm

gcd(a, b) = gcd(b, a mod b), and gcd(a, 0) = a. It's one of the oldest algorithms (≈300 BCE) and runs in O(log min(a, b)) steps.

The extended Euclidean algorithm also finds integers x, y with a·x + b·y = gcd(a, b) (Bézout's identity) — the key to modular inverses.

def gcd_steps(a, b):
    while b:
        print(f"{a} = {a // b} x {b} + {a % b}")
        a, b = b, a % b
    return a

def ext_gcd(a, b):
    if b == 0:
        return a, 1, 0
    g, x, y = ext_gcd(b, a % b)
    return g, y, x - (a // b) * y

print("gcd =", gcd_steps(252, 105))
g, x, y = ext_gcd(252, 105)
print(f"252*({x}) + 105*({y}) = {252 * x + 105 * y}")
print("lcm(252, 105) =", 252 * 105 // g)

252 = 2 x 105 + 42
105 = 2 x 42 + 21
42 = 2 x 21 + 0
gcd = 21
252*(-2) + 105*(5) = 21
lcm(252, 105) = 1260

3. Modular arithmetic

a ≡ b (mod m) means m divides a − b — "same remainder". Clock arithmetic: 10 o'clock + 5 hours = 3 o'clock (mod 12).

Rules (why "answer mod 10⁹+7" works — you can reduce at every step):

• (a + b) mod m = ((a mod m) + (b mod m)) mod m
• (a × b) mod m = ((a mod m) × (b mod m)) mod m
• Division is not that simple — you multiply by a modular inverse instead (section 5).

Fast modular exponentiation computes aᵇ mod m in O(log b) multiplications by squaring:

def power_mod(a, b, m):
    result = 1
    a %= m
    while b:
        if b & 1:
            result = result * a % m
        a = a * a % m
        b >>= 1
    return result

MOD = 10 ** 9 + 7
print(power_mod(2, 10, 1000), pow(2, 10, 1000))
print(power_mod(3, 10 ** 18, MOD), pow(3, 10 ** 18, MOD))   # 10^18 multiplications would never finish

24 24
246336683 246336683

Python's built-in pow(a, b, m) does exactly this.

4. Fermat's little theorem and Euler's totient

• Fermat: if p is prime and p ∤ a, then aᵖ⁻¹ ≡ 1 (mod p).
• Euler's totient φ(n) = how many numbers in 1…n are coprime to n. φ(p) = p − 1; φ(p·q) = (p − 1)(q − 1) for distinct primes.
• Euler's theorem: if gcd(a, n) = 1, then a^φ(n) ≡ 1 (mod n).

from math import gcd

phi = lambda n: sum(1 for k in range(1, n + 1) if gcd(k, n) == 1)
print("phi(10) =", phi(10), "| phi(13) =", phi(13), "| phi(33) =", phi(33), "= (3-1)*(11-1)")
print("Fermat: 5^(13-1) mod 13 =", pow(5, 12, 13))
print("Euler : 7^phi(10) mod 10 =", pow(7, phi(10), 10))

phi(10) = 4 | phi(13) = 12 | phi(33) = 20 = (3-1)*(11-1)
Fermat: 5^(13-1) mod 13 = 1
Euler : 7^phi(10) mod 10 = 1

5. Modular inverse

The inverse of a modulo m is x with a·x ≡ 1 (mod m). It exists iff gcd(a, m) = 1. Find it with the extended Euclidean algorithm — or, when m is prime, as a^(m−2) mod m (from Fermat). Python 3.8+: pow(a, -1, m).

Example: inverse of 3 mod 11 is 4, since 3 × 4 = 12 ≡ 1.

6. A tiny RSA (how HTTPS keys work, in miniature)

1. Pick two primes p, q; n = p·q (public).
2. φ(n) = (p − 1)(q − 1) (secret).
3. Choose public exponent e with gcd(e, φ) = 1.
4. Private key d = e⁻¹ mod φ(n).
5. Encrypt: c = mᵉ mod n. Decrypt: m = cᵈ mod n.

p, q = 61, 53
n = p * q
phi_n = (p - 1) * (q - 1)
e = 17
d = pow(e, -1, phi_n)                  # modular inverse
print(f"public key (n={n}, e={e}), private d={d}, check e*d mod phi = {e * d % phi_n}")

message = [ord(ch) for ch in "UPI"]   # tiny numbers < n
cipher = [pow(m, e, n) for m in message]
plain = "".join(chr(pow(c, d, n)) for c in cipher)
print("cipher:", cipher, "-> decrypted:", plain)

public key (n=3233, e=17), private d=2753, check e*d mod phi = 1
cipher: [2310, 2933, 1486] -> decrypted: UPI

Why it's secure: anyone can see n and e, but finding d requires φ(n), which requires factoring n. Real RSA uses primes hundreds of digits long (n of 2048+ bits), and nobody knows how to factor such numbers in reasonable time with classical computers. (Real systems also add padding — never use textbook RSA like this for actual security.)

Common mistakes

• Dividing directly under a modulus — use the modular inverse.
• Forgetting to reduce intermediate results (overflow in C++/Java).
• Assuming an inverse exists when gcd(a, m) ≠ 1.
• Testing primality by checking divisors up to n instead of √n (slow).

Quick recap

• Primes: unique factorisation; test up to √n; sieve for many.
• gcd(a, b) = gcd(b, a mod b); extended Euclid gives ax + by = gcd.
• Mod rules for + and ×; fast power in O(log b); inverse exists iff gcd = 1.
• Fermat: aᵖ⁻¹ ≡ 1 (mod p); Euler: a^φ(n) ≡ 1 (mod n).
• RSA: n = pq, ed ≡ 1 (mod φ(n)), c = mᵉ, m = cᵈ (mod n).