From 6ef8a6758a155a4109f43ee58334caa02ecf6820 Mon Sep 17 00:00:00 2001
From: akiyamn
Date: Fri, 14 May 2021 14:10:56 +1000
Subject: Ass 3: almost done q1

---
 ass3/q1/twin_prime.py | 74 +++++++++++++++++++++------------------------------
 1 file changed, 31 insertions(+), 43 deletions(-)

(limited to 'ass3')

diff --git a/ass3/q1/twin_prime.py b/ass3/q1/twin_prime.py
index aaeeef8..f31f0fb 100644
--- a/ass3/q1/twin_prime.py
+++ b/ass3/q1/twin_prime.py
@@ -1,12 +1,13 @@
 import random
 import math
+import sys
 
 
 def even(n):
     return n % 2 == 0
 
 
-def is_really_prime(n):
+def naive_prime_test(n):
     if n == 1 or even(n):
         return False
     for i in range(n - 1, 2, -1):
@@ -14,10 +15,6 @@ def is_really_prime(n):
             return False
     return True
 
-def pi2(n):
-    c = 0.661618
-    return 2 * c * (n / (math.log(n) ** 2))
-
 
 def is_prime(x, k=20):
     if even(x):
@@ -29,7 +26,6 @@ def is_prime(x, k=20):
     while even(t):
         s += 1
         t //= 2
-    assert math.floor(t) == t
     witnesses = [random.randrange(2, x - 1) for _ in range(k)]
     for w in witnesses:
         if power_mod(w, x - 1, x) != 1:
@@ -38,7 +34,6 @@ def is_prime(x, k=20):
             if power_mod(w, 2 ** i * t, x) == 1 and power_mod(w, 2 ** (i - 1) * t, x) not in [1, x - 1]:
                 return False
     return True
-    # return all(map(lambda w: is_prime(x, w), [2, 3, 5, 7, 11, 13, 17]))
 
 
 def power_mod(base, power, mod):
@@ -52,47 +47,40 @@ def power_mod(base, power, mod):
     return result % mod
 
 
-# for i in range(3, 100):
-#     print(i, is_prime(i), is_really_prime(i))
-#     if is_prime(i) != is_really_prime(i):
-#         print("Liar", i, is_prime(i), is_really_prime(i))
-
-# is_really_prime(9)
-
-
-# for a in filter(is_really_prime, range(50000, 500000)):
-#     if not miller_rabins(a, max(4, math.ceil(math.log(a, 4)))):
-#         raise Exception(f"Didn't work for {a}")
-
-
 def generate_twin_prime(m):
     twin = False
     while not twin:
-        n = random.randrange(2 ** (m - 1), 2 ** m)
+        start, end = bit_range(m)
+        n = random.randrange(start + 6 - (start % 6), end, 6)
+        # n = random.randrange(2 ** (m - 1) + 6 - ((2 ** (m - 1)) % 6), 2 ** m, 6)
+        # n += 6 - (n % 6)
+        assert n % 6 == 0
         iterations = max(4, math.ceil(math.log(n, 4)))
         twin = is_prime(n - 1, iterations) and is_prime(n + 1, iterations)
-    if not (is_really_prime(n - 1) and is_really_prime(n + 1)):
-        raise Exception(f"{n-1}: {is_really_prime(n - 1)}, {n+1}: {is_really_prime(n + 1)}")
+    if not (naive_prime_test(n - 1) and naive_prime_test(n + 1)):
+        raise Exception(f"{n - 1}: {naive_prime_test(n - 1)}, {n + 1}: {naive_prime_test(n + 1)}")
+    assert n + 1 < 2 ** m
+    # print(f"{bin(2**m -1)[2:]}\n{bin(n-1)[2:]}\n{bin(n+1)[2:]}")
     return n - 1, n + 1
 
+def bit_range(m):
+    return 2 ** (m - 1), 2 ** m
+
+def main():
+    assert len(sys.argv) >= 2
+    twin = generate_twin_prime(int(sys.argv[1]))
+    write_twin(twin)
+    print(twin)
+
+
+def write_twin(twin):
+    with open("output_twin_prime.txt", "w") as file:
+        file.write(f"{twin[0]}\n{twin[1]}")
+
+
+for i in range(3, 32):
+    print(i)
+    generate_twin_prime(i)
 
-for bits in range(3, 64):
-    print(f"==== {bits} bits ==== {2 ** (bits - 1), 2 ** bits - 1}")
-    for _ in range(1000):
-        generate_twin_prime(bits)
-        if _ % 50 == 0:
-            print(".", end="")
-    print()
-
-# i = 0
-# while True:
-#     i += 1
-#     base, power, mod = random.randint(1, 1000), random.randint(1, 5000), random.randint(1, 1000)
-#     mine = power_mod(base, power, mod)
-#     real = (base ** power) % mod
-#     pass_test = mine == real
-#     if i % 10000 == 0:
-#         print(i)
-#     if not pass_test:
-#         # print(f"Failed for {base}^{power} mod {mod}. Got {mine}, expected {real}.")
-#         raise Exception(f"Failed for {base}^{power} mod {mod}. Got {mine}, expected {real}.")
+if __name__ == "__main__":
+    main()
-- 
cgit v1.2.3