diff --git a/maths/prime_sieve_eratosthenes.py b/maths/prime_sieve_eratosthenes.py index 32eef9165bba..40a8b9a13955 100644 --- a/maths/prime_sieve_eratosthenes.py +++ b/maths/prime_sieve_eratosthenes.py @@ -11,6 +11,10 @@ https://en.wikipedia.org/wiki/Sieve_of_Eratosthenes """ +from math import isqrt + +import numpy as np + def prime_sieve_eratosthenes(num: int) -> list[int]: """ @@ -45,10 +49,36 @@ def prime_sieve_eratosthenes(num: int) -> list[int]: return [prime for prime in range(2, num + 1) if primes[prime]] +def np_prime_sieve_eratosthenes(max_number: int) -> list[int]: + """ + Returns prime numbers below max_number. + See: https://en.wikipedia.org/wiki/Sieve_of_Eratosthenes + + >>> np_prime_sieve_eratosthenes(10) + [2, 3, 5, 7] + >>> np_prime_sieve_eratosthenes(2) + [] + """ + if max_number <= 2: + return [] + + # List containing a bool value for every odd number below max_number/2 + is_prime = np.ones(max_number // 2, dtype=bool) + + for i in range(3, isqrt(max_number - 1) + 1, 2): + if is_prime[i // 2]: + # Mark all multiple of i as not prime using list slicing + is_prime[i**2 // 2 :: i] = False + + primes = np.where(is_prime)[0] * 2 + 1 + primes[0] = 2 + return primes.tolist() + + if __name__ == "__main__": import doctest doctest.testmod() user_num = int(input("Enter a positive integer: ").strip()) - print(prime_sieve_eratosthenes(user_num)) + print(np_prime_sieve_eratosthenes(user_num)) diff --git a/project_euler/problem_070/sol1.py b/project_euler/problem_070/sol1.py index f1114a280a31..649b30a80693 100644 --- a/project_euler/problem_070/sol1.py +++ b/project_euler/problem_070/sol1.py @@ -26,30 +26,31 @@ References: Finding totients -https://en.wikipedia.org/wiki/Euler's_totient_function#Euler's_product_formula +https://en.wikipedia.org/wiki/Euler%27s_totient_function#Euler%27s_product_formula """ from __future__ import annotations import numpy as np +from maths.prime_sieve_eratosthenes import np_prime_sieve_eratosthenes -def get_totients(max_one: int) -> list[int]: + +def np_get_totients(limit) -> list[int]: """ Calculates a list of totients from 0 to max_one exclusive, using the definition of Euler's product formula. - >>> get_totients(5) + >>> np_get_totients(5) [0, 1, 1, 2, 2] - >>> get_totients(10) + >>> np_get_totients(10) [0, 1, 1, 2, 2, 4, 2, 6, 4, 6] """ - totients = np.arange(max_one) + totients = np.arange(limit) + primes = np_prime_sieve_eratosthenes(limit) - for i in range(2, max_one): - if totients[i] == i: - x = np.arange(i, max_one, i) # array of indexes to select - totients[x] -= totients[x] // i + for i in primes: + totients[i::i] -= totients[i::i] // i return totients.tolist() @@ -81,10 +82,10 @@ def solution(max_n: int = 10000000) -> int: >>> solution(10000) 4435 """ + totients = np_get_totients(max_n + 1) min_numerator = 1 # i min_denominator = 0 # φ(i) - totients = get_totients(max_n + 1) for i in range(2, max_n + 1): t = totients[i] @@ -97,4 +98,4 @@ def solution(max_n: int = 10000000) -> int: if __name__ == "__main__": - print(f"{solution() = }") + print(f"Solution : {solution()}")