primes.py 2.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113
  1. # BUG!! There is a memory leak in this code
  2. UPPER_BOUND = 5000000
  3. PREFIX = 32338
  4. # exit(0)
  5. class Node:
  6. def __init__(self):
  7. self.children = {}
  8. self.terminal = False
  9. class Sieve:
  10. def __init__(self, limit):
  11. self.limit = limit
  12. self.prime = [False] * (limit + 1)
  13. def to_list(self):
  14. result = [2, 3]
  15. for p in range(5, self.limit + 1):
  16. if self.prime[p]:
  17. result.append(p)
  18. return result
  19. def omit_squares(self):
  20. r = 5
  21. while r * r < self.limit:
  22. if self.prime[r]:
  23. i = r * r
  24. while i < self.limit:
  25. self.prime[i] = False
  26. i = i + r * r
  27. r += 1
  28. return self
  29. def step1(self, x, y):
  30. n = (4 * x * x) + (y * y)
  31. if n <= self.limit and (n % 12 == 1 or n % 12 == 5):
  32. self.prime[n] = not self.prime[n]
  33. def step2(self, x, y):
  34. n = (3 * x * x) + (y * y)
  35. if n <= self.limit and n % 12 == 7:
  36. self.prime[n] = not self.prime[n]
  37. def step3(self, x, y):
  38. n = (3 * x * x) - (y * y)
  39. if x > y and n <= self.limit and n % 12 == 11:
  40. self.prime[n] = not self.prime[n]
  41. def loop_y(self, x):
  42. y = 1
  43. while y * y < self.limit:
  44. self.step1(x, y)
  45. self.step2(x, y)
  46. self.step3(x, y)
  47. y += 1
  48. def loop_x(self):
  49. x = 1
  50. while x * x < self.limit:
  51. self.loop_y(x)
  52. x += 1
  53. def calc(self):
  54. self.loop_x()
  55. return self.omit_squares()
  56. def generate_trie(l):
  57. root = Node()
  58. for el in l:
  59. head = root
  60. for ch in str(el):
  61. if ch not in head.children:
  62. head.children[ch] = Node()
  63. head = head.children[ch]
  64. head.terminal = True
  65. return root
  66. def find(upper_bound, prefix):
  67. primes = Sieve(upper_bound).calc()
  68. str_prefix = str(prefix)
  69. head = generate_trie(primes.to_list())
  70. for ch in str_prefix:
  71. head = head.children.get(ch)
  72. if head is None: # either ch does not exist or the value is None
  73. return None
  74. queue, result = [(head, str_prefix)], []
  75. while queue:
  76. top, prefix = queue.pop()
  77. if top.terminal:
  78. result.append(int(prefix))
  79. for ch, v in top.children.items():
  80. queue.insert(0, (v, prefix + ch))
  81. result.sort()
  82. return result
  83. def verify():
  84. left = [2, 23, 29]
  85. right = find(100, 2)
  86. if left != right:
  87. print(f"{left} != {right}")
  88. exit(1)
  89. verify()
  90. results = find(UPPER_BOUND, PREFIX)
  91. assert results == [323381, 323383, 3233803, 3233809, 3233851, 3233863, 3233873, 3233887, 3233897]