I am trying to implement the algorithm described in this article Faster Sieve of Eratosthenes. I quite understand the idea but I cannot grasp how exactly could it be implemented via python code.
After some work I found a way to convert indexes from the sieve into numbers itself:
number = lambda i: 3 * (i + 2) - 1 - (i + 2) % 2
But the main problem is jumps i have to do after getting prime. Article explains it as:
6np ± p, where p is prime and n - some natural number.
Is there a way to describe jumps using index of the last found prime for such an idea?
Thanks in advance.
P.S. There is implementation in Objective-C I am quite new to programming and can only understand python and js code.
If you understand numpy as well as Python, look at this implementation of primesfrom2to, taken from this answer in StackOverflow.
def primesfrom2to(n):
# https://stackoverflow.com/questions/2068372/fastest-way-to-list-all-primes-below-n-in-python/3035188#3035188
""" Input n>=6, Returns a array of primes, 2 <= p < n """
sieve = np.ones(n/3 + (n%6==2), dtype=np.bool)
sieve[0] = False
for i in xrange(int(n**0.5)/3+1):
if sieve[i]:
k=3*i+1|1
sieve[ ((k*k)/3) ::2*k] = False
sieve[(k*k+4*k-2*k*(i&1))/3::2*k] = False
return np.r_[2,3,((3*np.nonzero(sieve)[0]+1)|1)]
In that answer I linked to, this routine was the fastest in building a list of primes. I use a variation of it in my own code exploring primes. Explaining this in detail would take much space, but it builds a sieve which leaves out multiples of 2 and of 3. Just two lines of code (the ones starting sieve[ and ending = False) mark out multiples of a newly-found prime. I think this is what you mean by "jumps... after getting prime". The code is tricky but is work studying. The code is for Python 2, legacy Python.
Here is some of my own Python 3 code that has more comments. You can use it for comparison.
def primesieve3rd(n):
"""Return 'a third of a prime sieve' for values less than n that
are in the form 6k-1 or 6k+1.
RETURN: A numpy boolean array. A True value means the associated
number is prime; False means 1 or composite.
NOTES: 1. If j is in that form then the index for its primality
test is j // 3.
2. The numbers 2 and 3 are not handled in the sieve.
3. This code is based on primesfrom2to in
<https://stackoverflow.com/questions/2068372/
fastest-way-to-list-all-primes-below-n-in-python/
3035188#3035188>
"""
if n <= 1: # values returning an empty array
return np.ones(0, dtype=np.bool_)
sieve = np.ones(n // 3 + (n % 6 == 2), dtype=np.bool_) # all True
sieve[0] = False # note 1 is not prime
for i in range(1, (math.ceil(math.sqrt(n)) + 1) // 3): # sometimes large
if sieve[i]:
k = 3 * i + 1 | 1 # the associated number for this cell
# Mark some of the stored multiples of the number as composite
sieve[k * k // 3 :: 2 * k] = False
# Mark the remaining stored multiples (k times next possible prime)
sieve[k * (k + 4 - 2*(i&1)) // 3 :: 2 * k] = False
return sieve
def primesfrom2to(n, sieve=None):
"""Return an array of prime numbers less than n.
RETURN: A numpty int64 (indices type) array.
NOTES: 1. This code is based on primesfrom2to in
<https://stackoverflow.com/questions/2068372/
fastest-way-to-list-all-primes-below-n-in-python/
3035188#3035188>
"""
if n <= 5:
return np.array([2, 3], dtype=np.intp)[:max(n - 2, 0)]
if sieve is None:
sieve = primesieve3rd(n)
elif n >= 3 * len(sieve) + 1 | 1: # the next number to note in the sieve
raise ValueError('Number out of range of the sieve in '
'primesfrom2to')
return np.r_[2, 3, 3 * np.nonzero(sieve)[0] + 1 | 1]
Ask me if there is something here you do not understand.