 # Efficient prime k-tuple search

a guest
Jan 25th, 2023
46
1
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1.
2. def sieve_of_eratosthenes(num):
3.     prime = [True for i in range(num+1)]
4.     # boolean array
5.     p = 2
6.     while (p * p <= num):
7.
8.         # If prime[p] is not
9.         # changed, then it is a prime
10.         if (prime[p] == True):
11.
12.             # Updating all multiples of p
13.             for i in range(p * p, num+1, p):
14.                 prime[i] = False
15.         p += 1
16.
17.     primes = []
18.
19.     # Print all prime numbers
20.     for p in range(2, num+1):
21.         if prime[p]:
22.             primes.append(p)
23.
24.     return primes
25.
26. def inverses_of_eratosthenes(p_m, primes):
27.     inverses = []
28.
29.     for p in primes:
30.         try:
31.             inverse = pow(p_m, -1, p)
32.         except:
33.             inverse = 0 # if no multiplicative inverse
34.
35.         inverses.append(inverse)
36.
37.     return inverses
38.
39. def get_primorial(m):
40.
41.     primes = sieve_of_eratosthenes(m*m+1)
42.
43.     p_m = 1
44.
45.     for i in range(0,m):
46.         p_m *= primes[i]
47.
48.     return p_m
49.
50. def fermat_prime_p(n):
51.     if n == 2:
52.         return True
53.     if not n & 1:
54.         return False
55.     if pow(2, n-1, n) == 1:
56.         return True
57.     return False
58.
59. def fermat_is_valid_pow(n, constellation_pattern):
60.     for offset in constellation_pattern:
61.         if not fermat_prime_p(n+offset):
62.             return False
63.     return True
64.
65. def get_eliminated_factors(primes, inverses, offset, p_m, m, t, constellation_pattern):
66.
67.     eliminated_factors = {}
68.     k_max = 200
69.
70.     i = 0
71.     for p in primes:
72.         p_m_inverse = inverses[i]
73.
74.         if p_m_inverse != 0:
75.
76.             for c_i in constellation_pattern:
77.                 f_p = ((-offset - c_i) * p_m_inverse ) % p
78.                 eliminated_factors[f_p] = True
79.
80.                 for k in range(0, k_max):
81.                         f_p+=p
82.                         eliminated_factors[f_p] = True
83.
84.         i+=1
85.     return eliminated_factors
86.
87. def wheel_factorization(p_m, o, t, eliminated_factors, constellation_pattern):
88.
89.     t_prime = t + p_m - (t % p_m)
90.
91.     f_start = int(t_prime / p_m)
92.
93.     f_max = 20000000
94.
95.     primality_tests = 0
96.     primes = 0
97.
98.     for f in range(f_start, f_start+f_max):
99.         if f not in eliminated_factors:
100.             primality_tests +=1
101.             candidate = p_m*f + o
102.             if fermat_is_valid_pow(candidate, constellation_pattern):
103.                 # print(candidate)
104.                 primes+=1
105.
106.
107.     print("primality_tests: ", primality_tests)
108.     print("primes: ", primes)
109.
110. if __name__ == "__main__":
111.     # Pick primorial number
112.     m = 4
113.     # Pick offset
114.     o = 97
115.     # Pick Difficulty
116.     t = 10_000_000
117.     # Pick Constellation Pattern
118.     # constellation_pattern = [0, 6, 8, 14, 18, 20, 24, 26]
119.     constellation_pattern = [0, 4, 6, 10, 12, 16]
120.     # Pick prime_table_limit
121.     prime_table_limit = 1_000_000
122.
123.     p_m = get_primorial(m)
124.
125.     print("Calculating prime table up to ", prime_table_limit)
126.     primes = sieve_of_eratosthenes(prime_table_limit)
127.
128.     print("Calculating inverses...")
129.     inverses = inverses_of_eratosthenes(p_m, primes)
130.
131.     print("Sieving... ")
132.     eliminated_factors = {}
133.     eliminated_factors = get_eliminated_factors(primes, inverses, o, p_m, m, t, constellation_pattern)
134.
135.     print("Primality Testing candidates...")
136.     wheel_factorization(p_m, o, t, eliminated_factors, constellation_pattern)
137.
138.
139.
140.
141.
142.
143.
Tags: