beam: modular multiplicative inverse via extended Euclid, brute-force verified

wizard-work/beam/mod_inverse.py

@@ -0,0 +1,31 @@
+"""Modular multiplicative inverse via extended Euclid.
+
+modinv(a, m) returns x in [0, m) with a*x ≡ 1 (mod m), or None when
+gcd(a, m) != 1 (no inverse exists). Handles negative a, m>1 only.
+"""
+def egcd(a, b):
+    if b == 0:
+        return (a, 1, 0)
+    g, x1, y1 = egcd(b, a % b)
+    return (g, y1, x1 - (a // b) * y1)
+
+def modinv(a, m):
+    if not isinstance(m, int) or m <= 1:
+        raise ValueError("modulus must be int > 1")
+    a = a % m
+    g, x, _ = egcd(a, m)
+    if g != 1 and g != -1:
+        # normalize gcd sign
+        if g < 0:
+            g = -g
+            x = -x
+        if g != 1:
+            return None
+    x %= m
+    # verify
+    if (a * x) % m != 1:
+        return None
+    return x
+
+if __name__ == "__main__":
+    print(modinv(3, 11))

wizard-work/beam/test_mod_inverse.py

@@ -0,0 +1,45 @@
+import random
+from mod_inverse import modinv
+
+def test_basic():
+    assert modinv(3, 11) == 4
+    assert modinv(10, 17) == 12
+    assert modinv(1, 5) == 1
+
+def test_none():
+    assert modinv(2, 4) is None
+    assert modinv(6, 9) is None
+    assert modinv(0, 5) is None
+
+def test_negative_a():
+    assert modinv(-1, 7) == 6
+    assert ( (-1)*modinv(-1,7) ) % 7 == 1 % 7
+
+def test_brute():
+    for m in range(2, 200):
+        for a in range(0, m):
+            x = modinv(a, m)
+            g = None
+            # compute gcd
+            import math
+            if math.gcd(a % m, m) != 1:
+                assert x is None, (a, m, x)
+            else:
+                assert x is not None and 0 <= x < m
+                assert (a * x) % m == 1
+
+def test_random_large():
+    random.seed(1)
+    for _ in range(2000):
+        m = random.randint(2, 10**6)
+        a = random.randint(0, m-1)
+        x = modinv(a, m)
+        import math
+        if math.gcd(a, m) != 1:
+            assert x is None
+        else:
+            assert x is not None and (a*x)%m == 1
+
+if __name__ == "__main__":
+    test_basic(); test_none(); test_negative_a(); test_brute(); test_random_large()
+    print("ALL PASS")