Welcome, guest | Sign In | My Account | Store | Cart

I was looking for a fast way to list all primes below n, so far i came up to this with the numpy solution the fastest. It does primes up to 10e6 in 15ms in my old machine, and it is capable of reaching 10e9.

Python, 33 lines
 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36``` ``` #numpy solution import numpy def primesfrom2to(n): """ Input n>=6, Returns a array of primes, 2 <= p < n """ sieve = numpy.ones(n/3 + (n%6==2), dtype=numpy.bool) for i in xrange(1,int(n**0.5)/3+1): if sieve[i]: k=3*i+1|1 sieve[ k*k/3 ::2*k] = False sieve[k*(k-2*(i&1)+4)/3::2*k] = False return numpy.r_[2,3,((3*numpy.nonzero(sieve)[0][1:]+1)|1)] #pure-python solution def primes1(n): """ Returns a list of primes < n """ sieve = [True] * (n/2) for i in xrange(3,int(n**0.5)+1,2): if sieve[i/2]: sieve[i*i/2::i] = [False] * ((n-i*i-1)/(2*i)+1) return [2] + [2*i+1 for i in xrange(1,n/2) if sieve[i]] def primes2(n): """ Input n>=6, Returns a list of primes, 2 <= p < n """ n, correction = n-n%6+6, 2-(n%6>1) sieve = [True] * (n/3) for i in xrange(1,int(n**0.5)/3+1): if sieve[i]: k=3*i+1|1 sieve[ k*k/3 ::2*k] = [False] * ((n/6-k*k/6-1)/k+1) sieve[k*(k-2*(i&1)+4)/3::2*k] = [False] * ((n/6-k*(k-2*(i&1)+4)/6-1)/k+1) return [2,3] + [3*i+1|1 for i in xrange(1,n/3-correction) if sieve[i]] ```

I am looking for some feedback,

1. Does anyone know a faster way?
2. Any recommendations in how to improve this code?

This code was first posted here: Fastest way to list all primes below N in python.
There you can find some benchmarks.

Thomas Lehmann 13 years, 7 months ago

primes1 does not return primes only!

• The first output is by myself.
• the second output is from your primes1 function -> primes1(100).
• instead: primes2 is working fine.

[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97]

[2, 3, 5, 7, 11, 13, 17, 19, 23, 25, 29, 31, 35, 37, 41, 43, 47, 49, 53, 55, 59, 61, 65, 67, 71, 73, 77, 79, 83, 85, 89, 91, 95, 97]

robert.william.hanks (author) 13 years, 7 months ago

@Thomas Lehmann: Strange, well i am using python 2.6.

To test your claim i just copy/paste the code above in ideone.com.

here is the link -> http://ideone.com/MyojM it contains code & output.

As you can see it is working fine there.

Lorenz Fischer 7 years, 4 months ago
``````pypy -m timeit -r10 -s"from sympy import sieve" "primes = list(sieve.primerange(1, 10**6))"
10 loops, best of 10: 2.03 msec per loop
``````

So fast, so good. But if you are looking for the fastest isprime-function, you better use python with gmpy2.

``````python -m timeit -n10 -s"from gmpy2 import is_prime, mpz" "is_prime(mpz(2**521-1));is_prime(mpz(2**1279-1))"
10 loops, best of 3: 29.2 msec per loop
pypy -m timeit -n10 -s"from sympy import isprime" "isprime(2**521-1);isprime(2**1279-1)"
10 loops, best of 3: 214 msec per loop
``````

I couldn't install gmpy2 for pypy, hence i had to compare python with gmpy2 against pypy with sympy.

 Created by robert.william.hanks on Sat, 24 Jul 2010 (MIT)

### Required Modules

• (none specified)