A stupid “real-life” application of quadratic reciprocity

The application

During this year’s MOP, we used the following procedure to divide some of our students into two classes:

Let {p = 7075374838595186541578161} be prime. Take the letters in your name as it appears on the roster, convert them with A1Z26 and take the sum of cubes to get a number {s}. For example, EVANCHEN corresponds to {s = 5^3 + 22^3 + \dots + 14^3 = 16926}. Then you’re in Red 1 (room A155) if {s} is a quadratic residue modulo {p}, and Red 2 (room A133) otherwise.

The students were understandably a bit confused why the prime was chosen. It turned out to be a prank: if you ran the calculation on the 30-ish students in this class, it was actually just splitting them by last name. (That is, last names A-P in one class and R-Z in one class).

A few people have asked me how I managed to find a prime with this property, so I thought I’d post the source code here (without the names).

How I did it

After computing {s} for each student, you first start by taking a prime factorization of each {s} value. Let {Q} denote the set of primes dividing some {s}. Then for each {q \in Q} you want to commit to a value of the Legendre symbol {(\frac qp) = \pm 1}. This is an {\mathbb{F}_2} system of equations in {|Q|} variables, which one can solve easily, e.g. with z3. (Using sums of cubes was necessary to get enough different primes to show up. If you used say sum of squares, for this year’s students, it turned out there was no solution to the {\mathbb{F}_2} system, meaning no such prime {p} could be found.)

Then after choosing a solution {Q \to {\pm 1}}, quadratic reciprocity turns it into a condition on {p \pmod q}. To make things simpler I restricted to {p \equiv 1 \pmod 4}. The simplest thing to do would then be to pick a single residue class for each {q} that has the correct Legendre symbol, use the Chinese remainder theorem to put them all together, then search for a prime in that residue class.

However, if you do this the prime you get will have size like {\prod_{q \in Q} q}, which I thought was too large — it was something like 80 or 100 digits this way.

So instead, I arbitrarily picked a threshold on {Q} where I only used Chinese remainder theorem on the {q < 100}, and then just searched among that residue class for any prime that also passed a coin flip for each {q > 100}. That brought the size of the prime down to the {25}-digit one you see here, with only a minute or so of thinking in Python.

Code

Here’s the source code I used:


import string
from collections import defaultdict
from sympy import isprime, jacobi_symbol, legendre_symbol
from sympy.ntheory.modular import crt
from z3 import Bool, PbEq, Solver, sat
LARGE_PRIME_THRESHOLD = 100
# Read all lines from standard input
with open("rnames.txt") as f:
lines = [line.strip() for line in f if line.strip()]
# Sort by last name (everything after the first space)
sorted_lines = sorted(lines, key=lambda name: name.split()[1])
# Print the sorted lines
n = 0
nsums = []
for line in sorted_lines:
s = sum(
(string.ascii_uppercase.index(c) + 1) ** 3
for c in line.upper()
if c in string.ascii_uppercase
)
nsums.append(s)
print(line, s)
A = nsums[:14]
B = nsums[14:]
print(A, B)
# Factorization helper
def prime_factors(n):
i = 2
factors = defaultdict(int)
while i * i <= n:
while n % i == 0:
n //= i
factors[i] += 1
i += 1
if n > 1:
factors[n] += 1
return dict(factors)
def solve_completely_multiplicative_z3(A, B):
all_nums = A + B
primes = set()
factorizations = {}
# Get all primes involved and factor each number
for n in all_nums:
f = prime_factors(n)
factorizations[n] = f
primes.update(f)
primes = sorted(primes)
f_vars = {p: Bool(f"f_{p}") for p in primes} # True = +1, False = -1
solver = Solver()
def parity_constraint(factors, desired_sign):
# We only care about odd exponents, since (-1)^even = 1
relevant = [f_vars[p] for p, e in factors.items() if e % 2 == 1]
if not relevant:
# f(n) = 1 for empty product, so desired_sign must be +1
return desired_sign == 1
# Count how many of the relevant primes are assigned -1 (i.e., f_p = False)
# Count number of False values → parity must match desired_sign
# desired_sign = 1 ⇒ even number of Falses; -1 ⇒ odd number
k = len(relevant)
parity = 0 if desired_sign == 1 else 1
solver.add(PbEq([(v, 1) for v in relevant], k – parity))
for a in A:
parity_constraint(factorizations[a], +1)
for b in B:
parity_constraint(factorizations[b], -1)
if solver.check() == sat:
model = solver.model()
assignment = {p: 1 if model.eval(f_vars[p]) else -1 for p in primes}
return assignment
else:
return None
result = solve_completely_multiplicative_z3(A, B)
assert result is not None
print(result)
print(len(result))
residues = []
moduli = []
if 2 in result:
moduli.append(8)
if result.pop(2) == 1:
residues.append(1)
else:
residues.append(5)
for q in result:
if q > LARGE_PRIME_THRESHOLD:
continue
moduli.append(q)
if result[q] == 1:
residues.append(1)
else:
for t in range(2, q – 1):
if jacobi_symbol(t, q) == -1:
residues.append(t)
break
print(residues, moduli)
large_primes = [q for q in result.keys() if q > LARGE_PRIME_THRESHOLD]
r, M = crt(moduli, residues)
print(r, M)
print(f"{len(large_primes)} large prime conditions to deal with…")
for p in range(r, r + 10**8 * M, M):
if not isprime(p):
continue
if all(legendre_symbol(q, p) == result[q] for q in large_primes):
print(p)
break
else:
raise Exception("too bad so sad")
with open("prime.txt", "w") as f:
print(p, file=f)

view raw

qr.py

hosted with ❤ by GitHub



			

4 thoughts on “A stupid “real-life” application of quadratic reciprocity”

    1. [note: the above 2 Johns are the same the comment was sort of broken I was trying to see if Evan banned my email address (for posting a perhaps too much hint to a diamong in OTIS) and didn’t tell me but after using an alt email got through and main email also got through I found this wasn’t the case.]

      Like

Leave a comment