Sieve of Eratosthenes without multiples of 2 and 3

Question

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.

You can find the Objective C source code on his GitHub website. — Barmar, Feb 13 '19 at 23:55
Thanks for an answer! I am quite newbie and can understand only python and js code, unfortunately... — Eric Rovell, Feb 13 '19 at 23:57
If you can understand JS, I think you'll understand the C code, they're very similar. — Barmar, Feb 13 '19 at 23:59
@Barmar You were right. After some time I've decided to give it a try and check the given solution. There were some reduntant code in here, but I got the working python solution really fast. After some refactoring I will post it here. Thank you! — Eric Rovell, Feb 16 '19 at 00:05

Rory Daulton · Accepted Answer · 2019-02-14T00:13:14.343

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.

Thanks a lot! I will need some time to consume this code :) I really want to understand the implementation of this non-multiples of 2;3 idea. — Eric Rovell, Feb 14 '19 at 00:17

Sieve of Eratosthenes without multiples of 2 and 3

1 Answers1