Chapter 43 — Code-Based Cryptography: McEliece and Goppa Codes#
“There are two types of cryptography in this world: cryptography that will stop your kid sister from reading your files, and cryptography that will stop major governments from reading your files.” — Bruce Schneier, Applied Cryptography, 1996
In this chapter we study code-based cryptography, one of the oldest families of public-key cryptosystems still considered secure against quantum computers. We implement McEliece encryption based on binary Goppa codes, explore error correction, and analyse the scheme’s security properties.
43.1 Historical Context#
In 1978, just one year after RSA was published, Robert McEliece proposed a radically different public-key cryptosystem based on error-correcting codes. Instead of relying on the difficulty of factoring integers, McEliece’s scheme relies on the hardness of decoding a general linear code — a problem known to be NP-hard.
McEliece chose binary Goppa codes (introduced by Valerii Denisovich Goppa in 1970) as the underlying code family. These codes possess an efficient algebraic decoding algorithm but, when disguised by a random scrambling matrix and permutation, appear indistinguishable from random linear codes.
In 1986, Harald Niederreiter proposed a dual variant that encrypts via syndromes rather than codewords. The two schemes were later shown to be equivalent in security.
Tip
Unlike RSA and elliptic curve cryptography, the McEliece cryptosystem has resisted all known quantum attacks. Shor’s algorithm does not apply because the underlying problem is not based on a hidden subgroup structure.
In the NIST Post-Quantum Cryptography standardization process (begun 2017), Classic McEliece advanced to the fourth round as a candidate for key encapsulation. Its conservative design — essentially McEliece’s original 1978 proposal with modern parameters — reflects decades of cryptanalytic confidence.
Year |
Event |
|---|---|
1970 |
Goppa introduces algebraic-geometry codes |
1978 |
McEliece proposes code-based public-key encryption |
1986 |
Niederreiter publishes the dual (syndrome-based) variant |
1988 |
Li, Deng, Wang show equivalence of McEliece and Niederreiter |
2008 |
Bernstein, Lange, Peters propose efficient implementations |
2017 |
Classic McEliece submitted to NIST PQC competition |
2022 |
Classic McEliece advances to NIST Round 4 |
43.2 Formal Definitions#
Definition 43.1 — Linear Code
A binary linear code \(\mathcal{C}\) with parameters \([n, k, d]\) is a \(k\)-dimensional subspace of \(\mathbb{F}_2^n\) with minimum Hamming distance \(d\). It can correct up to \(t = \lfloor (d-1)/2 \rfloor\) errors.
Definition 43.2 — Generator and Parity-Check Matrices
A generator matrix \(G \in \mathbb{F}_2^{k \times n}\) for code \(\mathcal{C}\) satisfies \(\mathcal{C} = \{\mathbf{m}G : \mathbf{m} \in \mathbb{F}_2^k\}\).
A parity-check matrix \(H \in \mathbb{F}_2^{(n-k) \times n}\) satisfies \(\mathcal{C} = \{\mathbf{c} \in \mathbb{F}_2^n : H\mathbf{c}^T = \mathbf{0}\}\).
We have \(GH^T = \mathbf{0}\) (mod 2).
Definition 43.3 — Syndrome Decoding
Given a received word \(\mathbf{r} = \mathbf{c} + \mathbf{e}\) (codeword plus error), the syndrome is \(\mathbf{s} = H\mathbf{r}^T = H\mathbf{e}^T\). The decoding problem is: given \(\mathbf{s}\), find the minimum-weight \(\mathbf{e}\) such that \(H\mathbf{e}^T = \mathbf{s}\).
Definition 43.4 — Binary Goppa Code
Let \(g(x) \in \mathbb{F}_{2^m}[x]\) be a polynomial of degree \(t\) with no repeated roots, and let \(L = (\alpha_0, \ldots, \alpha_{n-1})\) be a sequence of distinct elements from \(\mathbb{F}_{2^m}\) such that \(g(\alpha_i) \neq 0\) for all \(i\). The binary Goppa code \(\Gamma(L, g)\) is:
This code has parameters \([n, k \geq n - mt, d \geq 2t+1]\).
Definition 43.5 — McEliece Cryptosystem
Key Generation: Choose a binary Goppa code with generator matrix \(G_0\). Select a random invertible \(k \times k\) matrix \(S\) and a random \(n \times n\) permutation matrix \(P\). Set the public key \(G_{\mathrm{pub}} = S G_0 P\).
Encryption: To encrypt message \(\mathbf{m} \in \mathbb{F}_2^k\), compute \(\mathbf{c} = \mathbf{m} G_{\mathrm{pub}} + \mathbf{e}\) where \(\mathbf{e}\) is a random error vector of weight \(t\).
Decryption: Compute \(\mathbf{c}' = \mathbf{c} P^{-1}\), decode \(\mathbf{c}'\) using the Goppa decoder to remove errors, then recover \(\mathbf{m} = \mathbf{m}' S^{-1}\).
Theorem 43.1 — Hardness of Generic Decoding
The Computational Syndrome Decoding problem is NP-complete: given a random parity-check matrix \(H \in \mathbb{F}_2^{r \times n}\), a syndrome \(\mathbf{s} \in \mathbb{F}_2^r\), and an integer \(w\), decide whether there exists \(\mathbf{e}\) of Hamming weight \(\leq w\) with \(H\mathbf{e}^T = \mathbf{s}\).
— Berlekamp, McEliece, and van Tilborg (1978)
43.3 Implementation#
We implement the core components step by step: linear codes over \(\mathbb{F}_2\), a simplified Goppa code construction, and the full McEliece cryptosystem.
43.3.1 Linear Code Class#
import numpy as np
import math
class LinearCode:
"""
Binary linear code [n, k, d] over GF(2).
Represented by a generator matrix G (k x n) and parity-check
matrix H ((n-k) x n). All arithmetic is modulo 2.
"""
def __init__(self, G):
"""Initialise from a generator matrix G (k x n) over GF(2)."""
self.G = np.array(G, dtype=int) % 2
self.k, self.n = self.G.shape
self.H = self._parity_check_matrix()
self.r = self.n - self.k # redundancy
def _rref_gf2(self, M):
"""Reduced row echelon form over GF(2)."""
A = np.array(M, dtype=int) % 2
rows, cols = A.shape
pivot_row = 0
pivot_cols = []
for col in range(cols):
# Find pivot
found = -1
for row in range(pivot_row, rows):
if A[row, col] == 1:
found = row
break
if found == -1:
continue
# Swap
A[[pivot_row, found]] = A[[found, pivot_row]]
pivot_cols.append(col)
# Eliminate
for row in range(rows):
if row != pivot_row and A[row, col] == 1:
A[row] = (A[row] + A[pivot_row]) % 2
pivot_row += 1
return A, pivot_cols
def _parity_check_matrix(self):
"""Compute parity-check matrix H such that G H^T = 0 (mod 2)."""
G_rref, pivot_cols = self._rref_gf2(self.G)
n, k = self.n, self.k
r = n - k
# Identify non-pivot columns
non_pivot = [c for c in range(n) if c not in pivot_cols]
# Build H: for systematic form [I_k | P], H = [-P^T | I_r] = [P^T | I_r] over GF(2)
H = np.zeros((r, n), dtype=int)
for i, np_col in enumerate(non_pivot):
H[i, np_col] = 1
for j, p_col in enumerate(pivot_cols):
if j < len(pivot_cols):
H[i, p_col] = G_rref[j, np_col]
# Verify: G H^T = 0 mod 2
check = (self.G @ H.T) % 2
assert np.all(check == 0), "Parity check verification failed"
return H
def encode(self, message):
"""Encode a k-bit message to an n-bit codeword."""
m = np.array(message, dtype=int) % 2
assert len(m) == self.k, f"Message must be {self.k} bits"
return (m @ self.G) % 2
def syndrome(self, received):
"""Compute the syndrome of a received word."""
r = np.array(received, dtype=int) % 2
return (self.H @ r) % 2
def decode_syndrome(self, received, max_weight=None):
"""
Syndrome-based decoding: find the lowest-weight error pattern
matching the syndrome (brute-force for small codes).
"""
r = np.array(received, dtype=int) % 2
s = self.syndrome(r)
if np.all(s == 0):
return r, np.zeros(self.n, dtype=int)
if max_weight is None:
max_weight = self.n
# Try error patterns of increasing weight
from itertools import combinations
for w in range(1, max_weight + 1):
for positions in combinations(range(self.n), w):
e = np.zeros(self.n, dtype=int)
for p in positions:
e[p] = 1
if np.array_equal((self.H @ e) % 2, s):
corrected = (r + e) % 2
return corrected, e
return None, None # Decoding failure
def minimum_distance(self):
"""Compute minimum distance by enumerating all codewords."""
min_d = self.n + 1
for i in range(1, 2**self.k):
m = np.array([(i >> b) & 1 for b in range(self.k)], dtype=int)
cw = self.encode(m)
w = int(np.sum(cw))
if 0 < w < min_d:
min_d = w
return min_d
def __repr__(self):
d = self.minimum_distance()
return f"LinearCode[{self.n}, {self.k}, {d}]"
# --- Demo: Hamming [7, 4, 3] code ---
G_hamming = np.array([
[1, 0, 0, 0, 1, 1, 0],
[0, 1, 0, 0, 1, 0, 1],
[0, 0, 1, 0, 0, 1, 1],
[0, 0, 0, 1, 1, 1, 1],
], dtype=int)
hamming = LinearCode(G_hamming)
print(f"Code: {hamming}")
print(f"Generator matrix G ({hamming.k}x{hamming.n}):")
print(hamming.G)
print(f"\nParity-check matrix H ({hamming.r}x{hamming.n}):")
print(hamming.H)
# Encode and decode
msg = [1, 0, 1, 1]
codeword = hamming.encode(msg)
print(f"\nMessage: {msg}")
print(f"Codeword: {codeword}")
print(f"Syndrome: {hamming.syndrome(codeword)}")
# Introduce a single-bit error
received = codeword.copy()
received[2] ^= 1 # flip bit 2
print(f"\nReceived (error at pos 2): {received}")
print(f"Syndrome: {hamming.syndrome(received)}")
corrected, error = hamming.decode_syndrome(received, max_weight=1)
print(f"Decoded: {corrected}")
print(f"Error: {error}")
print(f"Match: {np.array_equal(corrected, codeword)}")
Code: LinearCode[7, 4, 3]
Generator matrix G (4x7):
[[1 0 0 0 1 1 0]
[0 1 0 0 1 0 1]
[0 0 1 0 0 1 1]
[0 0 0 1 1 1 1]]
Parity-check matrix H (3x7):
[[1 1 0 1 1 0 0]
[1 0 1 1 0 1 0]
[0 1 1 1 0 0 1]]
Message: [1, 0, 1, 1]
Codeword: [1 0 1 1 0 1 0]
Syndrome: [0 0 0]
Received (error at pos 2): [1 0 0 1 0 1 0]
Syndrome: [0 1 1]
Decoded: [1 0 1 1 0 1 0]
Error: [0 0 1 0 0 0 0]
Match: True
Implementation Note
The LinearCode class uses brute-force syndrome decoding by searching over
error patterns of increasing weight. This is \(O\binom{n}{t}\) for correcting
\(t\) errors — acceptable for small codes but infeasible for cryptographic
parameters. The Goppa decoder below uses algebraic structure for efficient
decoding.
43.3.2 Simplified Goppa Code#
import numpy as np
import math
from itertools import combinations
class GF2m:
"""
Arithmetic in GF(2^m) using polynomial representation.
Elements are integers whose bits represent polynomial coefficients.
"""
# Irreducible polynomials for small m (as integers)
IRREDUCIBLES = {
2: 0b111, # x^2 + x + 1
3: 0b1011, # x^3 + x + 1
4: 0b10011, # x^4 + x + 1
5: 0b100101, # x^5 + x^2 + 1
6: 0b1000011, # x^6 + x + 1
7: 0b10000011, # x^7 + x + 1
8: 0b100011011, # x^8 + x^4 + x^3 + x + 1
}
def __init__(self, m):
self.m = m
self.order = 2**m
self.mod_poly = self.IRREDUCIBLES[m]
# Precompute log and exp tables
self.exp_table = [0] * (2 * self.order)
self.log_table = [0] * self.order
self._build_tables()
def _build_tables(self):
"""Build logarithm and exponentiation tables."""
x = 1
for i in range(self.order - 1):
self.exp_table[i] = x
self.log_table[x] = i
x <<= 1
if x >= self.order:
x ^= self.mod_poly
# Extend exp_table for convenience
for i in range(self.order - 1, 2 * self.order):
self.exp_table[i] = self.exp_table[i % (self.order - 1)]
def mul(self, a, b):
"""Multiply two elements."""
if a == 0 or b == 0:
return 0
return self.exp_table[(self.log_table[a] + self.log_table[b]) % (self.order - 1)]
def inv(self, a):
"""Multiplicative inverse."""
assert a != 0, "Cannot invert zero"
return self.exp_table[(self.order - 1 - self.log_table[a]) % (self.order - 1)]
def div(self, a, b):
"""Division."""
return self.mul(a, self.inv(b))
def add(self, a, b):
"""Addition = XOR."""
return a ^ b
def pow(self, a, e):
"""Exponentiation by squaring."""
if e == 0:
return 1
if a == 0:
return 0
return self.exp_table[(self.log_table[a] * e) % (self.order - 1)]
def poly_eval(self, coeffs, x):
"""Evaluate polynomial at x. coeffs[i] = coefficient of x^i."""
result = 0
for i in range(len(coeffs) - 1, -1, -1):
result = self.add(self.mul(result, x), coeffs[i])
return result
def poly_mul(self, a, b):
"""Multiply two polynomials over GF(2^m)."""
if not a or not b:
return [0]
result = [0] * (len(a) + len(b) - 1)
for i, ai in enumerate(a):
for j, bj in enumerate(b):
result[i + j] = self.add(result[i + j], self.mul(ai, bj))
return result
def poly_mod(self, a, b):
"""Compute a mod b for polynomials over GF(2^m)."""
a = list(a)
# Strip trailing zeros
while len(a) > 1 and a[-1] == 0:
a.pop()
while len(b) > 1 and b[-1] == 0:
b.pop()
if len(a) < len(b):
return a
inv_lead = self.inv(b[-1])
while len(a) >= len(b):
if a[-1] == 0:
a.pop()
continue
coeff = self.mul(a[-1], inv_lead)
shift = len(a) - len(b)
for i in range(len(b)):
a[shift + i] = self.add(a[shift + i], self.mul(coeff, b[i]))
while len(a) > 1 and a[-1] == 0:
a.pop()
return a if a != [0] else [0]
def poly_gcd(self, a, b):
"""GCD of two polynomials."""
while b != [0] and any(x != 0 for x in b):
a, b = b, self.poly_mod(a, b)
return a
# --- Demo ---
gf = GF2m(4) # GF(2^4) = GF(16)
print(f"GF(2^4): {gf.order} elements")
print(f"Irreducible polynomial: {bin(gf.mod_poly)}")
# Show all elements and their inverses
print(f"\n{'Element':>10} {'Inverse':>10} {'Check':>10}")
for a in range(1, gf.order):
a_inv = gf.inv(a)
check = gf.mul(a, a_inv)
if a <= 5:
print(f"{a:>10} {a_inv:>10} {check:>10}")
print(f" ... ({gf.order - 1} nonzero elements total)")
GF(2^4): 16 elements
Irreducible polynomial: 0b10011
Element Inverse Check
1 1 1
2 9 1
3 14 1
4 13 1
5 11 1
... (15 nonzero elements total)
import numpy as np
import math
from itertools import combinations
class GoppaCode:
"""
Simplified binary Goppa code Gamma(L, g).
Parameters:
m: extension degree (field is GF(2^m))
t: error-correction capability (degree of Goppa polynomial)
goppa_poly: coefficients of g(x) in GF(2^m)[x] (degree t)
support: list of n elements from GF(2^m) with g(alpha_i) != 0
"""
def __init__(self, m, t, goppa_poly=None, support=None, seed=None):
self.field = GF2m(m)
self.m = m
self.t = t
rng = np.random.default_rng(seed)
if goppa_poly is None:
goppa_poly = self._random_goppa_poly(t, rng)
self.goppa_poly = goppa_poly
if support is None:
support = self._build_support()
self.support = support
self.n = len(support)
# Build parity-check matrix H over GF(2^m)
self.H_ext = self._parity_check_ext()
# Convert to binary parity-check matrix
self.H_bin = self._to_binary_matrix(self.H_ext)
# Build generator matrix
self.G, self.k = self._generator_matrix()
def _random_goppa_poly(self, t, rng):
"""Generate a random irreducible-like Goppa polynomial of degree t."""
# For simplicity: random monic polynomial of degree t
# (true McEliece requires irreducible g(x), but for demonstration
# any square-free g(x) with g(alpha_i) != 0 works)
coeffs = [int(rng.integers(1, self.field.order)) for _ in range(t)]
coeffs.append(1) # monic
return coeffs
def _build_support(self):
"""Build support set: all field elements where g(alpha) != 0."""
support = []
for alpha in range(self.field.order):
if self.field.poly_eval(self.goppa_poly, alpha) != 0:
support.append(alpha)
return support
def _parity_check_ext(self):
"""Build parity-check matrix over GF(2^m).
H[i][j] = alpha_j^i / g(alpha_j) for i=0,...,t-1, j=0,...,n-1
"""
n = self.n
t = self.t
H = np.zeros((t, n), dtype=int)
for j in range(n):
alpha_j = self.support[j]
g_val = self.field.poly_eval(self.goppa_poly, alpha_j)
g_inv = self.field.inv(g_val)
for i in range(t):
alpha_pow = self.field.pow(alpha_j, i)
H[i, j] = self.field.mul(alpha_pow, g_inv)
return H
def _to_binary_matrix(self, H_ext):
"""Expand GF(2^m) parity-check matrix to binary.
Each GF(2^m) entry becomes an m-bit column.
Result: (m*t) x n binary matrix.
"""
t, n = H_ext.shape
m = self.m
H_bin = np.zeros((m * t, n), dtype=int)
for i in range(t):
for j in range(n):
val = H_ext[i, j]
for b in range(m):
H_bin[i * m + b, j] = (val >> b) & 1
return H_bin
def _generator_matrix(self):
"""Compute generator matrix from binary parity-check matrix.
Use the kernel of H_bin to find G.
"""
H = self.H_bin.copy()
mt, n = H.shape
# Compute RREF of H
A = H.copy()
pivot_row = 0
pivot_cols = []
for col in range(n):
found = -1
for row in range(pivot_row, mt):
if A[row, col] == 1:
found = row
break
if found == -1:
continue
A[[pivot_row, found]] = A[[found, pivot_row]]
pivot_cols.append(col)
for row in range(mt):
if row != pivot_row and A[row, col] == 1:
A[row] = (A[row] + A[pivot_row]) % 2
pivot_row += 1
rank = len(pivot_cols)
k = n - rank
free_cols = [c for c in range(n) if c not in pivot_cols]
# Each free column gives a basis vector of the kernel
G = np.zeros((k, n), dtype=int)
for idx, fc in enumerate(free_cols):
G[idx, fc] = 1
for j, pc in enumerate(pivot_cols):
G[idx, pc] = A[j, fc]
# Verify G H^T = 0
assert np.all((G @ self.H_bin.T) % 2 == 0), "Generator check failed"
return G, k
def encode(self, message):
"""Encode a k-bit message."""
m = np.array(message, dtype=int) % 2
assert len(m) == self.k
return (m @ self.G) % 2
def syndrome(self, r):
"""Compute binary syndrome."""
return (self.H_bin @ np.array(r, dtype=int)) % 2
def syndrome_poly(self, r):
"""Compute syndrome polynomial S(x) = sum r_i / (x - alpha_i) mod g(x)."""
r = np.array(r, dtype=int) % 2
# S(x) = sum_{i: r_i=1} 1/(x - alpha_i) mod g(x)
# We compute this as a polynomial in GF(2^m)[x]
f = self.field
n = self.n
# Build S(x) = sum r_i * (product_{j != i} (x - alpha_j)) / (product_{j != i} (alpha_i - alpha_j))
# More directly: evaluate via partial fractions
# S(x) mod g(x): accumulate r_i / (x - alpha_i) mod g(x)
result = [0] # zero polynomial
for i in range(n):
if r[i] == 0:
continue
# 1 / (x - alpha_i) mod g(x)
# Compute as: g_inv(alpha_i) * ... actually, simpler:
# Just accumulate the syndrome values
pass
# Simpler approach: use the extended parity-check matrix
# S(x) = s_0 + s_1 x + ... + s_{t-1} x^{t-1}
# where s_i = sum_j r_j * alpha_j^i / g(alpha_j)
s = np.zeros(self.t, dtype=int)
for j in range(n):
if r[j] == 0:
continue
alpha_j = self.support[j]
g_inv = f.inv(f.poly_eval(self.goppa_poly, alpha_j))
for i in range(self.t):
s[i] = f.add(s[i], f.mul(f.pow(alpha_j, i), g_inv))
return list(s)
def decode(self, received):
"""
Decode using Patterson's algorithm (simplified).
1. Compute syndrome polynomial S(x)
2. Find error locator polynomial sigma(x)
3. Find roots of sigma(x) to locate errors
"""
r = np.array(received, dtype=int) % 2
f = self.field
# Step 1: Compute syndrome
s = self.syndrome_poly(r)
if all(si == 0 for si in s):
return r, np.zeros(self.n, dtype=int)
# Step 2: For small t, use brute-force error location
# Try all combinations of up to t error positions
bin_syndrome = self.syndrome(r)
for w in range(1, self.t + 1):
for positions in combinations(range(self.n), w):
e = np.zeros(self.n, dtype=int)
for p in positions:
e[p] = 1
if np.array_equal((self.H_bin @ e) % 2, bin_syndrome):
corrected = (r + e) % 2
return corrected, e
return None, None # Decoding failure
# --- Demo: small Goppa code with m=3, t=2 ---
goppa = GoppaCode(m=3, t=2, seed=42)
print(f"Goppa code: n={goppa.n}, k={goppa.k}, t={goppa.t}")
print(f"Field: GF(2^{goppa.m}) = GF({goppa.field.order})")
print(f"Goppa polynomial g(x) coefficients: {goppa.goppa_poly}")
print(f"Support set L: {goppa.support}")
print(f"Generator matrix shape: {goppa.G.shape}")
print(f"Parity-check matrix shape: {goppa.H_bin.shape}")
# Encode a random message
rng = np.random.default_rng(123)
msg = rng.integers(0, 2, size=goppa.k)
codeword = goppa.encode(msg)
print(f"\nMessage ({goppa.k} bits): {msg}")
print(f"Codeword ({goppa.n} bits): {codeword}")
print(f"Syndrome: {goppa.syndrome(codeword)}")
# Add errors and decode
error = np.zeros(goppa.n, dtype=int)
err_pos = rng.choice(goppa.n, size=min(2, goppa.t), replace=False)
for p in err_pos:
error[p] = 1
received = (codeword + error) % 2
print(f"\nError positions: {err_pos}")
print(f"Received: {received}")
decoded, found_error = goppa.decode(received)
if decoded is not None:
print(f"Decoded: {decoded}")
print(f"Error found: {found_error}")
print(f"Correct: {np.array_equal(decoded, codeword)}")
else:
print("Decoding failed!")
Goppa code: n=8, k=2, t=2
Field: GF(2^3) = GF(8)
Goppa polynomial g(x) coefficients: [1, 6, 1]
Support set L: [0, 1, 2, 3, 4, 5, 6, 7]
Generator matrix shape: (2, 8)
Parity-check matrix shape: (6, 8)
Message (2 bits): [0 1]
Codeword (8 bits): [1 0 1 0 1 1 0 1]
Syndrome: [0 0 0 0 0 0]
Error positions: [4 0]
Received: [0 0 1 0 0 1 0 1]
Decoded: [1 0 1 0 1 1 0 1]
Error found: [1 0 0 0 1 0 0 0]
Correct: True
Patterson’s Algorithm — Conceptual Overview
Patterson’s algorithm (1975) is the efficient decoding algorithm for binary Goppa codes. Given received word \(\mathbf{r}\) with syndrome polynomial \(S(x)\):
Compute \(T(x) = S(x)^{-1} + x \pmod{g(x)}\)
Find \(\tau(x) = \sqrt{T(x)} \pmod{g(x)}\) (possible because \(\mathrm{char} = 2\))
Use the extended Euclidean algorithm to find \(a(x), b(x)\) with \(a(x) = b(x) \tau(x) \pmod{g(x)}\) and \(\deg a \leq \lfloor t/2 \rfloor\), \(\deg b \leq \lfloor (t-1)/2 \rfloor\)
Error locator polynomial: \(\sigma(x) = a(x)^2 + x \, b(x)^2\)
Find roots of \(\sigma(x)\) among the support elements; these are the error positions.
The key insight is that over \(\mathbb{F}_2\), square roots exist and are unique, making step 2 possible. The algorithm runs in \(O(n t)\) time.
43.3.3 McEliece Cryptosystem#
Pedagogical Simplification
This implementation demonstrates the McEliece framework but differs from Classic McEliece (NIST submission) in several ways: the Goppa polynomial is not required to be irreducible, and decoding uses exhaustive error search rather than Patterson’s algorithm. See the Classic McEliece specification for the full construction.
import numpy as np
import math
from itertools import combinations
class McEliece:
"""
McEliece public-key cryptosystem based on binary Goppa codes.
Key generation: disguise a Goppa code with scrambling + permutation.
Encryption: encode message + add random errors.
Decryption: undo permutation, decode Goppa code, undo scrambling.
"""
def __init__(self, goppa_code, seed=None):
self.code = goppa_code
self.n = goppa_code.n
self.k = goppa_code.k
self.t = goppa_code.t
rng = np.random.default_rng(seed)
# Private key components
self.S = self._random_invertible_matrix(self.k, rng) # k x k scrambler
self.S_inv = self._invert_gf2(self.S)
self.perm = rng.permutation(self.n) # permutation
self.perm_inv = np.argsort(self.perm) # inverse permutation
# Public key: G_pub = S * G * P
P = np.eye(self.n, dtype=int)[self.perm] # permutation matrix
self.G_pub = (self.S @ goppa_code.G @ P.T) % 2
# Verify dimensions
assert self.G_pub.shape == (self.k, self.n)
def _random_invertible_matrix(self, size, rng):
"""Generate a random invertible binary matrix."""
while True:
M = rng.integers(0, 2, size=(size, size))
if self._gf2_rank(M) == size:
return M
def _gf2_rank(self, M):
"""Compute rank of a binary matrix."""
A = np.array(M, dtype=int) % 2
rows, cols = A.shape
rank = 0
for col in range(cols):
found = -1
for row in range(rank, rows):
if A[row, col] == 1:
found = row
break
if found == -1:
continue
A[[rank, found]] = A[[found, rank]]
for row in range(rows):
if row != rank and A[row, col] == 1:
A[row] = (A[row] + A[rank]) % 2
rank += 1
return rank
def _invert_gf2(self, M):
"""Invert a binary matrix using Gauss-Jordan elimination."""
size = M.shape[0]
aug = np.hstack([M.copy(), np.eye(size, dtype=int)])
for col in range(size):
# Find pivot
found = -1
for row in range(col, size):
if aug[row, col] == 1:
found = row
break
assert found != -1, "Matrix is not invertible"
aug[[col, found]] = aug[[found, col]]
for row in range(size):
if row != col and aug[row, col] == 1:
aug[row] = (aug[row] + aug[col]) % 2
return aug[:, size:]
def encrypt(self, message, seed=None):
"""
Encrypt a k-bit message.
c = m * G_pub + e, where wt(e) = t.
"""
rng = np.random.default_rng(seed)
m = np.array(message, dtype=int) % 2
assert len(m) == self.k, f"Message must be {self.k} bits"
# Encode
c = (m @ self.G_pub) % 2
# Add random error of weight t
e = np.zeros(self.n, dtype=int)
error_positions = rng.choice(self.n, size=self.t, replace=False)
e[error_positions] = 1
ciphertext = (c + e) % 2
return ciphertext, e
def decrypt(self, ciphertext):
"""
Decrypt a ciphertext.
1. Undo permutation: c' = c * P^{-1}
2. Decode Goppa code to remove errors
3. Undo scrambling: m = m' * S^{-1}
"""
c = np.array(ciphertext, dtype=int) % 2
# Step 1: Undo permutation
c_perm = c[self.perm_inv]
# Step 2: Decode using Goppa code
decoded, error = self.code.decode(c_perm)
if decoded is None:
return None # Decoding failure
# Step 3: Extract message bits and undo scrambling
# decoded = m' * G where m' = m * S
# We need to extract m' from the codeword
# Since G is in systematic-like form, solve m' G = decoded
m_prime = self._extract_message(decoded)
# Undo scrambling
m = (m_prime @ self.S_inv) % 2
return m
def _extract_message(self, codeword):
"""Extract message from codeword by solving m*G = codeword over GF(2)."""
G = self.code.G
k, n = G.shape
# Solve m*G = codeword, i.e., G^T m^T = codeword^T
aug = np.hstack([G.T, codeword.reshape(-1, 1)])
rows, cols = aug.shape
pivot_row = 0
for col in range(cols - 1):
found = -1
for row in range(pivot_row, rows):
if aug[row, col] == 1:
found = row
break
if found == -1:
continue
aug[[pivot_row, found]] = aug[[found, pivot_row]]
for row in range(rows):
if row != pivot_row and aug[row, col] == 1:
aug[row] = (aug[row] + aug[pivot_row]) % 2
pivot_row += 1
return aug[:k, -1]
def public_key_size(self):
"""Public key size in bits."""
return self.k * self.n
def private_key_size(self):
"""Approximate private key size in bits."""
# S (k*k) + permutation (n*log2(n)) + Goppa polynomial
return self.k * self.k + self.n * int(math.ceil(math.log2(self.n + 1))) + self.t * self.code.m
# --- Demo ---
print("Setting up McEliece with small Goppa code (m=3, t=2)...")
goppa = GoppaCode(m=3, t=2, seed=42)
mce = McEliece(goppa, seed=99)
print(f"Code parameters: n={mce.n}, k={mce.k}, t={mce.t}")
print(f"Public key size: {mce.public_key_size()} bits ({mce.public_key_size()//8} bytes)")
print(f"Private key size: ~{mce.private_key_size()} bits (~{mce.private_key_size()//8} bytes)")
print(f"G_pub shape: {mce.G_pub.shape}")
Setting up McEliece with small Goppa code (m=3, t=2)...
Code parameters: n=8, k=2, t=2
Public key size: 16 bits (2 bytes)
Private key size: ~42 bits (~5 bytes)
G_pub shape: (2, 8)
43.4 Experiments#
Experiment 1: McEliece Encrypt/Decrypt Demo#
We demonstrate the full McEliece encryption and decryption cycle with multiple messages, verifying correctness.
import numpy as np
import math
from itertools import combinations
# Use the McEliece instance from above
print(f"McEliece parameters: n={mce.n}, k={mce.k}, t={mce.t}")
print(f"Can correct up to {mce.t} errors in {mce.n}-bit codewords")
print("=" * 60)
rng = np.random.default_rng(2024)
n_trials = 10
successes = 0
for trial in range(n_trials):
# Random message
msg = rng.integers(0, 2, size=mce.k)
# Encrypt
ciphertext, error = mce.encrypt(msg, seed=trial)
# Decrypt
recovered = mce.decrypt(ciphertext)
if recovered is not None and np.array_equal(msg, recovered):
successes += 1
status = "OK"
else:
status = "FAIL"
if trial < 5: # Print first 5 trials in detail
print(f"\nTrial {trial + 1}:")
print(f" Message: {msg}")
print(f" Ciphertext: {ciphertext}")
print(f" Error (wt={int(np.sum(error))}): {error}")
if recovered is not None:
print(f" Recovered: {recovered}")
print(f" Status: {status}")
print(f"\n{'=' * 60}")
print(f"Results: {successes}/{n_trials} successful decryptions")
McEliece parameters: n=8, k=2, t=2
Can correct up to 2 errors in 8-bit codewords
============================================================
Trial 1:
Message: [0 1]
Ciphertext: [1 1 1 1 1 1 0 0]
Error (wt=2): [0 0 0 0 0 1 0 1]
Recovered: [0 1]
Status: OK
Trial 2:
Message: [0 0]
Ciphertext: [0 0 0 1 1 0 0 0]
Error (wt=2): [0 0 0 1 1 0 0 0]
Recovered: [0 0]
Status: OK
Trial 3:
Message: [0 0]
Ciphertext: [0 0 1 0 0 1 0 0]
Error (wt=2): [0 0 1 0 0 1 0 0]
Recovered: [0 0]
Status: OK
Trial 4:
Message: [1 1]
Ciphertext: [0 0 0 1 0 0 1 1]
Error (wt=2): [1 0 0 0 0 1 0 0]
Recovered: [1 1]
Status: OK
Trial 5:
Message: [1 1]
Ciphertext: [1 0 0 1 0 0 1 0]
Error (wt=2): [0 0 0 0 0 1 0 1]
Recovered: [1 1]
Status: OK
============================================================
Results: 10/10 successful decryptions
Experiment 2: Key Size Comparison#
One of the main drawbacks of code-based cryptography is the large public key size. We compare McEliece key sizes across different parameters and against other cryptosystems.
Show code cell source
import numpy as np
import matplotlib.pyplot as plt
import math
# McEliece key sizes for various parameter sets
# Using actual Classic McEliece submission parameters where applicable.
# Public key in systematic form: k * (n - k) bits.
# Format: (label, n, k, t, m)
param_sets = [
("Low (m=10)", 1024, 524, 50, 10),
("Medium (m=11)", 2048, 1498, 50, 11),
("mceliece348864 (L1)", 3488, 2720, 64, 12),
("mceliece460896 (L3)", 4608, 3360, 96, 13),
("mceliece6960119 (L5)", 6960, 5413, 119, 13),
("mceliece8192128 (L5)", 8192, 6528, 128, 13),
]
names = []
pub_key_kb = []
n_values = []
k_values = []
security_bits = []
print(f"{'Parameters':<26} {'n':>6} {'k':>6} {'t':>4} {'Pub Key (KB)':>14} {'~Sec (bits)':>12}")
print("=" * 72)
for label, n, k, t, m in param_sets:
# Systematic public key: k rows of (n - k) bits
pub_size_bits = k * (n - k)
pub_size_kb = pub_size_bits / 8 / 1024
# Rough security approximation (NOT a rigorous ISD estimate).
# See cell below for Prange/Stern work factors.
sec = min(t * m // 2, 256)
names.append(f"n={n}, t={t}")
pub_key_kb.append(pub_size_kb)
n_values.append(n)
k_values.append(k)
security_bits.append(sec)
print(f"{label:<26} {n:>6} {k:>6} {t:>4} {float(pub_size_kb):>12.1f} {sec:>10}")
print("\nNote: security column is a rough approximation (t*m/2).")
print("For precise estimates, see the ISD analysis below and consult")
print("the ISD literature (e.g., Bernstein et al., 2008).")
# Comparison with other systems
other_systems = {
"RSA-2048": 0.256,
"RSA-4096": 0.512,
"ECC P-256": 0.032,
"ECC P-384": 0.048,
"Kyber-768 (lattice)": 1.184,
"Kyber-1024 (lattice)": 1.568,
}
fig, ax = plt.subplots(figsize=(12, 6))
# McEliece bars
x_pos = np.arange(len(names))
bars1 = ax.bar(x_pos, pub_key_kb, color='#e74c3c', alpha=0.8, label='McEliece')
# Other systems
other_names = list(other_systems.keys())
other_sizes = list(other_systems.values())
x_pos2 = np.arange(len(other_names)) + len(names) + 1
colors = ['#3498db'] * 2 + ['#2ecc71'] * 2 + ['#9b59b6'] * 2
bars2 = ax.bar(x_pos2, other_sizes, color=colors, alpha=0.8)
# Labels
all_labels = names + [''] + other_names
ax.set_xticks(list(x_pos) + [len(names)] + list(x_pos2))
ax.set_xticklabels(all_labels, rotation=45, ha='right', fontsize=8)
ax.set_ylabel('Public Key Size (KB)')
ax.set_title('Public Key Sizes: McEliece vs Other Cryptosystems')
ax.set_yscale('log')
ax.grid(True, alpha=0.3, axis='y')
# Add value labels on bars
for bar, val in zip(list(bars1) + list(bars2), pub_key_kb + other_sizes):
ax.text(bar.get_x() + bar.get_width()/2, bar.get_height() * 1.1,
f'{float(val):.1f}' if val >= 1 else f'{float(val):.3f}',
ha='center', va='bottom', fontsize=7)
# Legend
from matplotlib.patches import Patch
legend_elements = [
Patch(facecolor='#e74c3c', alpha=0.8, label='McEliece (code-based)'),
Patch(facecolor='#3498db', alpha=0.8, label='RSA (factoring)'),
Patch(facecolor='#2ecc71', alpha=0.8, label='ECC (discrete log)'),
Patch(facecolor='#9b59b6', alpha=0.8, label='Kyber (lattice-based)'),
]
ax.legend(handles=legend_elements, loc='upper left')
plt.tight_layout()
plt.savefig('mceliece_key_sizes.png', dpi=150, bbox_inches='tight')
plt.show()
print(f"\nMcEliece NIST-1 public key is ~{float(pub_key_kb[2]/other_systems['Kyber-768 (lattice)']):.0f}x larger than Kyber-768")
print(f"McEliece NIST-1 public key is ~{float(pub_key_kb[2]/other_systems['RSA-2048']):.0f}x larger than RSA-2048")
Parameters n k t Pub Key (KB) ~Sec (bits)
========================================================================
Low (m=10) 1024 524 50 32.0 250
Medium (m=11) 2048 1498 50 100.6 256
mceliece348864 (L1) 3488 2720 64 255.0 256
mceliece460896 (L3) 4608 3360 96 511.9 256
mceliece6960119 (L5) 6960 5413 119 1022.2 256
mceliece8192128 (L5) 8192 6528 128 1326.0 256
Note: security column is a rough approximation (t*m/2).
For precise estimates, see the ISD analysis below and consult
the ISD literature (e.g., Bernstein et al., 2008).
McEliece NIST-1 public key is ~215x larger than Kyber-768
McEliece NIST-1 public key is ~996x larger than RSA-2048
Key Size Trade-Off
McEliece’s main practical drawback is its enormous public key. The Classic McEliece mceliece348864 parameter set (NIST Level 1, \(n=3488\), \(k=2720\), \(t=64\)) has a public key of approximately 255 KB — over 1,000 times larger than RSA-2048 and over 200 times larger than Kyber-768. However, its encryption and decryption operations are extremely fast compared to lattice-based schemes.
Experiment 3: Error Correction Visualisation#
We visualise how the Goppa decoder corrects errors by showing the codeword, error pattern, received word, and corrected output.
Show code cell source
import numpy as np
import matplotlib.pyplot as plt
from itertools import combinations
# Use the Goppa code from earlier
rng = np.random.default_rng(7)
# Generate a codeword
msg = rng.integers(0, 2, size=goppa.k)
codeword = goppa.encode(msg)
# Create error patterns of different weights
fig, axes = plt.subplots(3, 1, figsize=(12, 8), sharex=True)
error_weights = [1, 2, min(3, goppa.t + 1)] # up to t and one beyond
titles = []
for idx, w in enumerate(error_weights):
ax = axes[idx]
n = goppa.n
# Create error
error = np.zeros(n, dtype=int)
if w <= n:
err_pos = rng.choice(n, size=min(w, n), replace=False)
for p in err_pos:
error[p] = 1
received = (codeword + error) % 2
# Try to decode
decoded, found_error = goppa.decode(received)
success = decoded is not None and np.array_equal(decoded, codeword)
# Plot
positions = np.arange(n)
width = 0.25
# Codeword bits
ax.bar(positions - width, codeword, width, color='#3498db', alpha=0.7, label='Codeword')
# Error bits
ax.bar(positions, error * 0.8, width, color='#e74c3c', alpha=0.8, label='Error')
# Received bits
ax.bar(positions + width, received, width, color='#f39c12', alpha=0.6, label='Received')
# Mark error positions
for p in range(n):
if error[p] == 1:
ax.axvline(x=p, color='red', linestyle='--', alpha=0.3)
status = "CORRECTED" if success else "FAILED"
color = '#27ae60' if success else '#c0392b'
ax.set_title(f'Error weight = {w}: {status}', fontsize=11, color=color, fontweight='bold')
ax.set_ylabel('Bit value')
ax.set_ylim(-0.1, 1.3)
ax.legend(loc='upper right', fontsize=8)
ax.grid(True, alpha=0.2)
axes[-1].set_xlabel('Bit position')
fig.suptitle(f'Error Correction in Goppa Code [n={goppa.n}, k={goppa.k}, t={goppa.t}]',
fontsize=13, fontweight='bold', y=1.01)
plt.tight_layout()
plt.savefig('error_correction_viz.png', dpi=150, bbox_inches='tight')
plt.show()
Experiment 4: Security#
Analysis — Information Set Decoding
The best known generic attack against the McEliece cryptosystem is information set decoding (ISD). We estimate the work factor for different parameter choices.
Show code cell source
import numpy as np
import matplotlib.pyplot as plt
import math
def isd_work_factor(n, k, t):
"""
Estimate the work factor of basic information set decoding (Prange's algorithm).
Pick k random columns of H. If the error vector has weight 0 in those positions,
we can solve for the error directly. Probability of success:
p = C(n-t, k) / C(n, k)
Work factor ~ 1/p * k^2 (for Gaussian elimination)
"""
# Use log to avoid overflow: log2(C(n,k)) = sum log2(n-i+1) - sum log2(i)
def log2_comb(n_val, k_val):
if k_val > n_val or k_val < 0:
return -float('inf')
if k_val == 0 or k_val == n_val:
return 0
result = 0.0
for i in range(1, k_val + 1):
result += math.log2(n_val - i + 1) - math.log2(i)
return result
# log2(1/p) = log2(C(n,k)) - log2(C(n-t,k))
log2_inv_p = log2_comb(n, k) - log2_comb(n - t, k)
# Total work: 1/p * k^2 (Gaussian elimination)
log2_work = log2_inv_p + 2 * math.log2(k)
return log2_work
def stern_work_factor(n, k, t, p_param=2):
"""
Estimate work factor of Stern's algorithm (improved ISD).
Splits information set, looks for collisions. Reduces work by roughly
a factor of C(k/2, p)^2 / 2^l for optimal parameters.
Simplified estimate: ~ Prange work / C(k/2, p)
"""
prange = isd_work_factor(n, k, t)
# Improvement factor (simplified)
improvement = 0.0
half_k = k // 2
for i in range(1, p_param + 1):
improvement += math.log2(half_k - i + 1) - math.log2(i)
return prange - improvement
# Compute for various parameter sets
param_sets_security = [
("Original McEliece", 1024, 524, 50),
("Moderate", 2048, 1608, 40),
("NIST Level 1", 3488, 2720, 64),
("NIST Level 3", 4608, 3360, 96),
("NIST Level 5", 6688, 5024, 128),
("Conservative", 6960, 5413, 119),
("Maximum", 8192, 6528, 128),
]
labels = []
prange_wf = []
stern_wf = []
print(f"{'Name':<22} {'n':>6} {'k':>6} {'t':>4} {'Prange':>10} {'Stern':>10}")
print("=" * 62)
for name, n, k, t in param_sets_security:
pw = isd_work_factor(n, k, t)
sw = stern_work_factor(n, k, t, p_param=2)
labels.append(name)
prange_wf.append(pw)
stern_wf.append(sw)
print(f"{name:<22} {n:>6} {k:>6} {t:>4} {float(pw):>10.1f} {float(sw):>10.1f}")
fig, ax = plt.subplots(figsize=(10, 6))
x = np.arange(len(labels))
width = 0.35
bars1 = ax.bar(x - width/2, prange_wf, width, color='#e74c3c', alpha=0.8, label="Prange's ISD")
bars2 = ax.bar(x + width/2, stern_wf, width, color='#3498db', alpha=0.8, label="Stern's ISD")
# Security thresholds
ax.axhline(y=128, color='green', linestyle='--', alpha=0.5, label='128-bit security')
ax.axhline(y=192, color='orange', linestyle='--', alpha=0.5, label='192-bit security')
ax.axhline(y=256, color='red', linestyle='--', alpha=0.5, label='256-bit security')
ax.set_ylabel('Work Factor (log$_2$ operations)')
ax.set_title('McEliece Security: Information Set Decoding Work Factors')
ax.set_xticks(x)
ax.set_xticklabels(labels, rotation=30, ha='right', fontsize=9)
ax.legend(loc='upper left', fontsize=9)
ax.grid(True, alpha=0.3, axis='y')
ax.set_ylim(0, max(prange_wf) * 1.15)
plt.tight_layout()
plt.savefig('mceliece_security.png', dpi=150, bbox_inches='tight')
plt.show()
Name n k t Prange Stern
==============================================================
Original McEliece 1024 524 50 71.7 56.6
Moderate 2048 1608 40 112.1 93.8
NIST Level 1 3488 2720 64 165.6 145.8
NIST Level 3 4608 3360 96 208.3 187.9
NIST Level 5 6688 5024 128 286.9 265.4
Conservative 6960 5413 119 288.2 266.4
Maximum 8192 6528 128 325.5 303.1
Security Observation
The NIST Level 1 parameters \((n=3488, k=2720, t=64)\) achieve at least 128 bits of security against the best known classical and quantum attacks. Stern’s algorithm (and its successors like BJMM) improve on Prange but do not fundamentally change the exponential nature of the attack. Grover’s quantum search provides at most a quadratic speedup, roughly halving the security level — which the NIST parameters are designed to withstand.
Experiment 5: Structural Analysis of the Public Key#
We verify that the disguised public key \(G_{\mathrm{pub}} = S G_0 P\) looks random compared to the structured Goppa generator matrix \(G_0\).
Show code cell source
import numpy as np
import matplotlib.pyplot as plt
# Compare structure of original G vs public G_pub
G_original = goppa.G
G_public = mce.G_pub
fig, axes = plt.subplots(1, 3, figsize=(14, 4))
# Original generator matrix
ax0 = axes[0]
ax0.imshow(G_original, cmap='Blues', aspect='auto', interpolation='nearest')
ax0.set_title('Original Goppa G$_0$', fontsize=11)
ax0.set_xlabel('Column (n)')
ax0.set_ylabel('Row (k)')
# Public key
ax1 = axes[1]
ax1.imshow(G_public, cmap='Reds', aspect='auto', interpolation='nearest')
ax1.set_title('Public Key G$_{pub}$ = S G$_0$ P', fontsize=11)
ax1.set_xlabel('Column (n)')
ax1.set_ylabel('Row (k)')
# Random matrix for comparison
rng_comp = np.random.default_rng(555)
G_random = rng_comp.integers(0, 2, size=G_original.shape)
ax2 = axes[2]
ax2.imshow(G_random, cmap='Greens', aspect='auto', interpolation='nearest')
ax2.set_title('Random Binary Matrix', fontsize=11)
ax2.set_xlabel('Column (n)')
ax2.set_ylabel('Row (k)')
plt.suptitle('Structural Comparison: Goppa Code vs Disguised Public Key vs Random',
fontsize=12, fontweight='bold')
plt.tight_layout()
plt.savefig('key_structure_comparison.png', dpi=150, bbox_inches='tight')
plt.show()
# Statistical analysis: column weight distribution
def col_weight_stats(M):
weights = np.sum(M, axis=0)
return np.mean(weights), np.std(weights)
def row_weight_stats(M):
weights = np.sum(M, axis=1)
return np.mean(weights), np.std(weights)
print(f"\nStatistical Comparison:")
print(f"{'Matrix':<20} {'Col Mean':>10} {'Col Std':>10} {'Row Mean':>10} {'Row Std':>10} {'Density':>10}")
print("=" * 72)
for name, M in [("Original G0", G_original), ("Public G_pub", G_public), ("Random", G_random)]:
cm, cs = col_weight_stats(M)
rm, rs = row_weight_stats(M)
density = np.mean(M)
print(f"{name:<20} {float(cm):>10.2f} {float(cs):>10.2f} {float(rm):>10.2f} {float(rs):>10.2f} {float(density):>10.3f}")
Statistical Comparison:
Matrix Col Mean Col Std Row Mean Row Std Density
========================================================================
Original G0 1.25 0.43 5.00 0.00 0.625
Public G_pub 1.38 0.48 5.50 0.50 0.688
Random 0.88 0.78 3.50 0.50 0.438
Indistinguishability
The security of McEliece relies on the public key \(G_{\mathrm{pub}} = S G_0 P\) being computationally indistinguishable from a random generator matrix. The scrambling matrix \(S\) destroys any systematic structure in the rows, while the permutation \(P\) hides the column ordering that encodes the Goppa code’s algebraic structure. An attacker who could recover \(S\), \(P\), or \(G_0\) from \(G_{\mathrm{pub}}\) could use the efficient Goppa decoder to decrypt messages.
Experiment 6: Error Weight and Decoding Failure Rate#
We explore what happens when the error weight exceeds the correction capability \(t\) of the code.
Show code cell source
import numpy as np
import matplotlib.pyplot as plt
from itertools import combinations
# Test decoding success rate for various error weights
rng = np.random.default_rng(314)
n_tests = 30
max_weight = min(goppa.t + 3, goppa.n // 2)
weights = list(range(0, max_weight + 1))
success_rates = []
print(f"Goppa code: n={goppa.n}, k={goppa.k}, t={goppa.t}")
print(f"{'Weight':>8} {'Successes':>10} {'Trials':>8} {'Rate':>8}")
print("=" * 38)
for w in weights:
successes = 0
for trial in range(n_tests):
msg = rng.integers(0, 2, size=goppa.k)
codeword = goppa.encode(msg)
# Random error of weight w
error = np.zeros(goppa.n, dtype=int)
if w > 0:
err_pos = rng.choice(goppa.n, size=min(w, goppa.n), replace=False)
for p in err_pos:
error[p] = 1
received = (codeword + error) % 2
decoded, _ = goppa.decode(received)
if decoded is not None and np.array_equal(decoded, codeword):
successes += 1
rate = successes / n_tests
success_rates.append(rate)
print(f"{w:>8} {successes:>10} {n_tests:>8} {float(rate):>8.2f}")
fig, ax = plt.subplots(figsize=(9, 5))
colors = ['#27ae60' if w <= goppa.t else '#e74c3c' for w in weights]
ax.bar(weights, success_rates, color=colors, alpha=0.8, edgecolor='black', linewidth=0.5)
# Mark the correction capability
ax.axvline(x=goppa.t + 0.5, color='black', linestyle='--', linewidth=2, label=f't = {goppa.t} (correction limit)')
ax.set_xlabel('Error Weight', fontsize=11)
ax.set_ylabel('Decoding Success Rate', fontsize=11)
ax.set_title(f'Goppa Code [{goppa.n},{goppa.k}] Decoding vs Error Weight', fontsize=12)
ax.set_ylim(-0.05, 1.1)
ax.legend(fontsize=10)
ax.grid(True, alpha=0.3, axis='y')
plt.tight_layout()
plt.savefig('decoding_failure_rate.png', dpi=150, bbox_inches='tight')
plt.show()
Goppa code: n=8, k=2, t=2
Weight Successes Trials Rate
======================================
0 30 30 1.00
1 30 30 1.00
2 30 30 1.00
3 0 30 0.00
4 0 30 0.00
Danger
Using an error weight \(w > t\) in the McEliece scheme does not increase security — it causes decryption failure. The error weight must be exactly \(t\) for correct operation. An attacker who can induce the sender to use the wrong error weight could mount a denial-of-service attack.
Experiment 7: Syndrome Weight Distribution#
We analyse the distribution of syndrome weights for random error patterns to understand the information leakage from syndromes.
Show code cell source
import numpy as np
import matplotlib.pyplot as plt
# Analyse syndrome weights for the Hamming code
# (easier to enumerate all possibilities)
hamming_code = LinearCode(G_hamming)
# All possible error patterns and their syndromes
syndrome_weights_by_error_weight = {}
n = hamming_code.n
for val in range(2**n):
e = np.array([(val >> b) & 1 for b in range(n)], dtype=int)
ew = int(np.sum(e))
s = hamming_code.syndrome(e)
sw = int(np.sum(s))
if ew not in syndrome_weights_by_error_weight:
syndrome_weights_by_error_weight[ew] = []
syndrome_weights_by_error_weight[ew].append(sw)
fig, axes = plt.subplots(2, 2, figsize=(10, 8))
for idx, ew in enumerate([0, 1, 2, 3]):
ax = axes[idx // 2][idx % 2]
sw_list = syndrome_weights_by_error_weight.get(ew, [])
if sw_list:
max_sw = max(sw_list)
bins = np.arange(-0.5, max_sw + 1.5, 1)
ax.hist(sw_list, bins=bins, color='#3498db', alpha=0.8, edgecolor='black', linewidth=0.5)
ax.set_title(f'Error weight = {ew} ({len(sw_list)} patterns)', fontsize=10)
ax.set_xlabel('Syndrome weight')
ax.set_ylabel('Count')
ax.grid(True, alpha=0.3)
fig.suptitle('Syndrome Weight Distribution for Hamming [7,4,3] Code',
fontsize=12, fontweight='bold')
plt.tight_layout()
plt.savefig('syndrome_distribution.png', dpi=150, bbox_inches='tight')
plt.show()
print("Syndrome weight statistics:")
print(f"{'Error wt':>10} {'Count':>8} {'Mean S wt':>10} {'Max S wt':>10}")
print("=" * 42)
for ew in sorted(syndrome_weights_by_error_weight.keys()):
sw_list = syndrome_weights_by_error_weight[ew]
if ew <= 4:
print(f"{ew:>10} {len(sw_list):>8} {float(np.mean(sw_list)):>10.2f} {max(sw_list):>10}")
Syndrome weight statistics:
Error wt Count Mean S wt Max S wt
==========================================
0 1 0.00 0
1 7 1.71 3
2 21 1.71 3
3 35 1.37 3
4 35 1.37 3
Experiment 8: McEliece vs Niederreiter — Structural Comparison#
The Niederreiter variant encrypts by computing the syndrome of an error vector, rather than adding errors to an encoded message. We compare the two approaches.
Show code cell source
import numpy as np
import matplotlib.pyplot as plt
import math
def niederreiter_key_sizes(m, t):
"""
Niederreiter public key size.
Public key: (n-k) x n parity-check matrix in systematic form -> (n-k) x k bits.
But the message space is the set of weight-t vectors, so
message bits = floor(log2(C(n, t))).
"""
n = 2**m
k = n - m * t
r = n - k # = m*t
# Public key: systematic H -> r x k bits
pub_bits = r * k
# Message space
# log2(C(n, t))
msg_bits = 0.0
for i in range(1, t + 1):
msg_bits += math.log2(n - i + 1) - math.log2(i)
return n, k, r, pub_bits, int(msg_bits)
def mceliece_key_sizes(m, t):
"""McEliece public key size."""
n = 2**m
k = n - m * t
pub_bits = k * n
msg_bits = k
return n, k, n - k, pub_bits, msg_bits
params = [(10, 50), (11, 50), (12, 62), (12, 96), (13, 119)]
print(f"{'m':>3} {'t':>4} {'n':>6} | {'McEliece PK (KB)':>18} {'msg bits':>10} | {'Niederreiter PK (KB)':>22} {'msg bits':>10}")
print("=" * 100)
mc_sizes = []
ni_sizes = []
labels = []
for m, t in params:
mn, mk, mr, mc_pub, mc_msg = mceliece_key_sizes(m, t)
nn, nk, nr, ni_pub, ni_msg = niederreiter_key_sizes(m, t)
mc_kb = mc_pub / 8 / 1024
ni_kb = ni_pub / 8 / 1024
mc_sizes.append(mc_kb)
ni_sizes.append(ni_kb)
labels.append(f"m={m},t={t}")
print(f"{m:>3} {t:>4} {mn:>6} | {float(mc_kb):>16.1f} {mc_msg:>10} | {float(ni_kb):>20.1f} {ni_msg:>10}")
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 5))
x = np.arange(len(labels))
width = 0.35
ax1.bar(x - width/2, mc_sizes, width, color='#e74c3c', alpha=0.8, label='McEliece')
ax1.bar(x + width/2, ni_sizes, width, color='#3498db', alpha=0.8, label='Niederreiter')
ax1.set_xlabel('Parameter Set')
ax1.set_ylabel('Public Key Size (KB)')
ax1.set_title('Public Key Size Comparison')
ax1.set_xticks(x)
ax1.set_xticklabels(labels, rotation=30, ha='right')
ax1.legend()
ax1.grid(True, alpha=0.3, axis='y')
# Ratio plot
ratios = [mc / ni for mc, ni in zip(mc_sizes, ni_sizes)]
ax2.bar(x, ratios, color='#2ecc71', alpha=0.8)
ax2.set_xlabel('Parameter Set')
ax2.set_ylabel('McEliece / Niederreiter Size Ratio')
ax2.set_title('Relative Key Size (McEliece / Niederreiter)')
ax2.set_xticks(x)
ax2.set_xticklabels(labels, rotation=30, ha='right')
ax2.axhline(y=1, color='black', linestyle='--', alpha=0.5)
ax2.grid(True, alpha=0.3, axis='y')
for i, r in enumerate(ratios):
ax2.text(i, r + 0.02, f'{float(r):.2f}', ha='center', va='bottom', fontsize=9)
plt.tight_layout()
plt.savefig('mceliece_vs_niederreiter.png', dpi=150, bbox_inches='tight')
plt.show()
m t n | McEliece PK (KB) msg bits | Niederreiter PK (KB) msg bits
====================================================================================================
10 50 1024 | 65.5 524 | 32.0 284
11 50 2048 | 374.5 1498 | 100.6 334
12 62 4096 | 1676.0 3352 | 304.4 459
12 96 4096 | 1472.0 2944 | 414.0 652
13 119 8192 | 6645.0 6645 | 1254.9 892
McEliece vs Niederreiter
The Niederreiter variant uses the parity-check matrix \(H\) as the public key instead of the generator matrix \(G\). Since \(H\) is \((n-k) \times n\) vs \(G\) being \(k \times n\), and typically \(k > n-k\), Niederreiter has smaller public keys. However, encryption in Niederreiter requires generating a weight-\(t\) error vector that encodes the message (via a constant-weight encoding), making the message space \(\lfloor \log_2 \binom{n}{t} \rfloor\) bits instead of \(k\) bits.
Li, Deng, and Wang (1994) proved that the two schemes are equivalent in terms of security: breaking one implies breaking the other.
43.5 Exercises#
Exercise 43.1 (Warm-up) Construct the Hamming \([7, 4, 3]\) code’s syndrome decoding table: for each of the 7 possible single-bit error patterns, compute the syndrome and verify that all syndromes are distinct. Explain why this guarantees single-error correction.
Hint
The syndrome of error \(\mathbf{e}_i\) (error in position \(i\)) equals the \(i\)-th column of \(H\). Since \(H\) has \(2^r - 1 = 7\) distinct nonzero columns (one for each nonzero 3-bit vector), each single-bit error produces a unique syndrome.
Exercise 43.2 (Applied) Implement a brute-force attack on the McEliece cryptosystem for the small parameters used in this chapter. Given only the public key \(G_{\mathrm{pub}}\) and a ciphertext \(\mathbf{c}\), try all \(\binom{n}{t}\) possible error patterns to find the message. Measure how long this takes and compare with legitimate decryption.
Hint
For each candidate error \(\mathbf{e}\) of weight \(t\), compute
\(\mathbf{c}' = \mathbf{c} - \mathbf{e}\) and check whether \(\mathbf{c}'\) lies
in the row space of \(G_{\mathrm{pub}}\) (i.e., solve \(\mathbf{m} G_{\mathrm{pub}} = \mathbf{c}'\)
over \(\mathbb{F}_2\)). Use itertools.combinations for the enumeration.
Exercise 43.3 (Analysis) Compute the exact number of codewords in the Goppa code constructed in Section 43.3.2. Verify that it equals \(2^k\). Then compute the covering radius — the maximum Hamming distance from any \(n\)-bit vector to the nearest codeword.
Hint
Enumerate all \(2^k\) messages and their codewords. For the covering radius, you need to find the \(n\)-bit vector that is farthest from any codeword. For small \(n\), try all \(2^n\) vectors.
Exercise 43.4 (Theory) Prove that a binary Goppa code \(\Gamma(L, g)\) with \(\deg g = t\) has minimum distance \(d \geq 2t + 1\) (not just \(t + 1\) as for general alternant codes). The key step is showing that if \(g(x)\) is square-free over \(\mathbb{F}_{2^m}\), then \(\Gamma(L, g) = \Gamma(L, g^2) \cap \mathbb{F}_2^n\).
Hint
Show that for binary codes, the parity-check condition \(\sum c_i / (x - \alpha_i) \equiv 0 \pmod{g(x)}\) can be strengthened to \(\sum c_i / (x - \alpha_i) \equiv 0 \pmod{g(x)^2}\) because over \(\mathbb{F}_2\), the formal derivative of \(\sum c_i / (x - \alpha_i)\) equals \(\sum c_i / (x - \alpha_i)^2\), which is the square of the original sum in characteristic 2.
Exercise 43.5 (Challenge) Implement Patterson’s decoding algorithm (as described in the admonition in Section 43.3.2) for the small Goppa codes used in this chapter. Compare its performance with the brute-force decoder. For which parameter sizes does the algebraic decoder become faster?
Hint
The critical step is computing the square root \(\tau(x) = \sqrt{T(x)} \pmod{g(x)}\). Over \(\mathbb{F}_{2^m}\), the square root of \(a\) is \(a^{2^{m-1}}\). For polynomials modulo \(g(x)\), compute the matrix representation of the Frobenius map \(\phi: f(x) \mapsto f(x)^2 \pmod{g(x)}\) and invert it to get the square root.
43.6 Summary#
Concept |
Key Point |
|---|---|
Linear codes |
\([n, k, d]\) subspaces of \(\mathbb{F}_2^n\); correct \(\lfloor(d-1)/2\rfloor\) errors |
Goppa codes |
Algebraic codes with efficient decoding; \(d \geq 2t+1\) for degree-\(t\) polynomial |
McEliece encryption |
Disguise Goppa code with \(G_{\mathrm{pub}} = S G_0 P\); add weight-\(t\) errors |
Niederreiter variant |
Encrypt via syndromes; equivalent security; smaller keys |
Patterson’s algorithm |
\(O(nt)\) decoding for binary Goppa codes using square roots in char 2 |
Security basis |
NP-hardness of generic syndrome decoding; resists quantum attacks |
Key size trade-off |
Public keys are 100-1000x larger than RSA/ECC; fast encryption/decryption |
NIST standardisation |
Classic McEliece: conservative, 45+ years of cryptanalytic confidence |
Tip
Code-based cryptography provides one of the strongest known defences against quantum computers. While the large key sizes pose practical challenges, the McEliece system’s speed and its decades of resistance to cryptanalysis make it a compelling choice for high-security applications. In the next chapter, we will explore lattice-based cryptography, which offers a different balance of key sizes and computational efficiency.
References#
McEliece, R. J. (1978). A Public-Key Cryptosystem Based on Algebraic Coding Theory. DSN Progress Report, 42-44, 114–116.
Goppa, V. D. (1970). A New Class of Linear Correcting Codes. Problems of Information Transmission, 6(3), 207–212.
Niederreiter, H. (1986). Knapsack-Type Cryptosystems and Algebraic Coding Theory. Problems of Control and Information Theory, 15(2), 159–166.
Patterson, N. J. (1975). The Algebraic Decoding of Goppa Codes. IEEE Transactions on Information Theory, 21(2), 203–207.
Berlekamp, E. R., McEliece, R. J., and van Tilborg, H. C. A. (1978). On the Inherent Intractability of Certain Coding Problems. IEEE Transactions on Information Theory, 24(3), 384–386.
Bernstein, D. J., Lange, T., and Peters, C. (2008). Attacking and Defending the McEliece Cryptosystem. PQCrypto 2008, LNCS 5299, 31–46.
Bernstein, D. J. et al. (2017). Classic McEliece: Conservative Code-Based Cryptography. NIST PQC Submission.