Barrett reduction in Python
Reduces x mod n using a precomputed reciprocal approximation instead of division.
- File
barrett_reduction.py, 59 lines- Level
- I0 not checked yet
- Toolchain
- Python 3.14, pytest 9
- Test framework
- pytest
- Run the tests
python3 -m pytest
# SPDX-License-Identifier: MIT
# Part of https://algorithm.tamnd.com/barrett-reduction/. The authors are listed in prov/MXGY63.jsonl.
"""Barrett reduction: x mod n with two multiplications and shifts, no division.
The modulus n >= 1 has k binary digits, so 2^(k-1) <= n < 2^k. The setup
computes mu = floor(4^k / n) once. Then every x with 0 <= x < 4^k, which
includes every product of two residues in [0, n), is reduced with two
multiplications, two shifts, one subtraction and at most two corrections.
Python integers have arbitrary size, so n may have any number of digits.
"""
from typing import NamedTuple
class Barrett(NamedTuple):
n: int # the modulus, n >= 1
k: int # the number of binary digits of n
mu: int # floor(4^k / n)
def barrett_setup(n: int) -> Barrett:
"""Precompute the constants for the modulus n >= 1."""
if n < 1:
raise ValueError("the modulus must be a positive integer")
k = n.bit_length() # @line L1
mu = (1 << (2 * k)) // n # @line L2
return Barrett(n, k, mu)
def barrett_reduce(x: int, ctx: Barrett) -> int:
"""Return x mod n for 0 <= x < 4^k."""
n, k, mu = ctx
if x < 0 or x >> (2 * k):
raise ValueError("x must satisfy 0 <= x < 4^k")
q = ((x >> (k - 1)) * mu) >> (k + 1) # @line L3
r = x - q * n # @line L4
while r >= n: # @line L5
r -= n
return r # @line L6
def barrett_mulmod(a: int, b: int, ctx: Barrett) -> int:
"""Return a * b mod n for residues 0 <= a, b < n."""
return barrett_reduce(a * b, ctx) # @line L7
def barrett_powmod(base: int, exponent: int, ctx: Barrett) -> int:
"""Return base^exponent mod n by binary exponentiation with Barrett products."""
if exponent < 0:
raise ValueError("the exponent must be non-negative")
result = barrett_reduce(1 % ctx.n, ctx)
b = base % ctx.n
while exponent:
if exponent & 1:
result = barrett_mulmod(result, b, ctx)
b = barrett_mulmod(b, b, ctx)
exponent >>= 1
return result
test_barrett_reduction.py 83 lines
# SPDX-License-Identifier: MIT
# Part of https://algorithm.tamnd.com/barrett-reduction/. The authors are listed in prov/MXGY63.jsonl.
import random
import pytest
from barrett_reduction import barrett_mulmod, barrett_powmod, barrett_reduce, barrett_setup
def corrections(x, ctx):
"""The number of subtractions of n made by line L5, from the quotient estimate."""
q = ((x >> (ctx.k - 1)) * ctx.mu) >> (ctx.k + 1)
return x // ctx.n - q
def moduli(rng):
"""Moduli 1, 2, powers of two, small and large primes, composites, several hundred digits."""
ms = [1, 2, 3, 4, 5, 7, 8, 9, 15, 16, 47, 50, 97, 255, 256, 257, 65535, 65536, 65537]
ms += [1 << e for e in (10, 31, 32, 63, 64, 127, 1000)]
ms += [(1 << e) - 1 for e in (61, 89, 127)] + [(1 << e) + 1 for e in (32, 64)]
ms += [10**9 + 7, 998244353, 2**61 - 1, 10**18 + 9, 10**300 + 1, 3**700, 10**400 - 3]
ms += [rng.randrange(1, 1 << rng.randrange(1, 1500)) for _ in range(120)]
return ms
def test_setup_constants():
ctx = barrett_setup(47)
assert (ctx.n, ctx.k, ctx.mu) == (47, 6, 87)
with pytest.raises(ValueError):
barrett_setup(0)
def test_worked_examples():
assert barrett_reduce(3561, barrett_setup(47)) == 36
assert barrett_reduce(46 * 46, barrett_setup(47)) == 1
ctx = barrett_setup(50)
assert barrett_reduce(47 * 49, ctx) == 3
assert corrections(47 * 49, ctx) == 2
def test_exhaustive_small_moduli():
for n in range(1, 130):
ctx = barrett_setup(n)
for x in range(1 << (2 * ctx.k)):
assert barrett_reduce(x, ctx) == x % n
assert 0 <= corrections(x, ctx) <= 2
def test_random_against_python_mod():
rng = random.Random(20261106)
for n in moduli(rng):
ctx = barrett_setup(n)
top = 1 << (2 * ctx.k)
xs = [0, 1, n - 1, n, (n - 1) * (n - 1), top - 1, max(top - n, 0)]
xs += [rng.randrange(top) for _ in range(40)]
for x in xs:
assert barrett_reduce(x, ctx) == x % n
assert 0 <= corrections(x, ctx) <= 2
for _ in range(30):
a, b = rng.randrange(n), rng.randrange(n)
assert barrett_mulmod(a, b, ctx) == a * b % n
for a in (0, 1, n - 1):
assert barrett_mulmod(a, n - 1, ctx) == a * (n - 1) % n
def test_powmod_matches_pow():
rng = random.Random(7)
for n in moduli(rng)[:60]:
ctx = barrett_setup(n)
for _ in range(5):
base = rng.randrange(n)
e = rng.randrange(1 << rng.randrange(1, 200))
assert barrett_powmod(base, e, ctx) == pow(base, e, n)
assert barrett_powmod(0, 0, ctx) == pow(0, 0, n)
def test_rejects_out_of_range_inputs():
ctx = barrett_setup(47)
with pytest.raises(ValueError):
barrett_reduce(-1, ctx)
with pytest.raises(ValueError):
barrett_reduce(1 << 12, ctx)
What the next level needs
- the tests did not run in this build