Another Python prime-listing algorithm

Question

I wrote a small program in Python that finds primes (how original) up to some limit. I tried to produce very clean code, so I'd be vary grateful if you could point out any improvement I could do to make it more pythonic and more readable.

I'm also interested in knowing whether or not it can be sped up (with the same method I used here). I tried to use generator, but since I want to read my list of primes multiple times it didn't work... I'm also wondering if a set or some other structure could be better than a list here.

Thank you very much !

from __future__ import print_function
import time

def measure_time(func):
    def wrapper(*args, **kw):
        starting_time = time.clock()
        res = func(*args, **kw)
        ending_time = time.clock()
        print('Duration : %fs' % (ending_time - starting_time))
        return res
    return wrapper

def ask_user_input():
    """Relentlessly ask user to enter an integer >1."""
    while True:
        try:
            lim = int(raw_input('What\'s the limit?   '))
            if lim == int(lim) and lim > 2:
                break
            else:
                print('Must be an integer larger than 1. Try again.')
        except (NameError, ValueError):
            print('Not a number. Try again.')
    return lim

def print_primes(primes):
    for p in primes:
        print(p)

@measure_time
def construct_primes_list(lim):
    """Builds the list of primes smaller than lim."""
    prime_divisors = [2]

    for num in xrange(2, lim):
        for d1 in prime_divisors:
            if num % d1 == 0:
                break
        else:
            prime_divisors.append(num)

    return prime_divisors

def main():
    lim = ask_user_input()
    primes = construct_primes_list(lim)
    # print_primes(primes)


if __name__ == '__main__':
    main()

Here is the benchmarking on my computer:

lim      | time
-----------------------
100      | 0.000143s
1000     | 0.003231s
10000    | 0.071654s
100000   | 2.598098s
1000000  | 169.052574s

EDIT :

As ChatterOne's comment suggested, I tried to compute the square root of num because I know I couldn't find any divisors after that. But it gets way slower. Here is the code :

for num in xrange(3, lim):
    sqrt_num = int(sqrt(num)+1)
    for d1 in [p for p in prime_divisors if p < sqrt_num]:
        if num % d1 == 0:
            break
    else:
        prime_divisors.append(num)

Any thoughts?

Answer 1

Regarding the code:

for num in xrange(3, lim):
    sqrt_num = int(sqrt(num)+1)
    for d1 in [p for p in prime_divisors if p < sqrt_num]:
        if num % d1 == 0:
            break
    else:
        prime_divisors.append(num)

One small improvement would be to either change the squarebrackets in the for statement into parenthesis or just move the test into the loop. Something like:

for num in xrange(3, lim):
    for d1 in prime_divisors:
        if d1 > sqrt_num:
             break
        if num % d1 == 0:
             break
    else:
        prime_divisors.append(num)

As written, the code first makes loops over prime_divisors to make a temporary list of the ones < sqrt_num, and then loops over that temporary list to run the modulus test. Changing the square brackets to parenthesis would make it a generator instead of a list so that performance hit would go away. Actually, the whole thing could be rewritten as:

for num in xrange(3, lim):
    prime_divisors.extend([] if any(p for p in prime_divisors if p < sqrt_num and num % p == 0) else [num])

..but that may be taking things a bit far.

Answer 2

You will probably get better performance from the sieve of Eratosthenes. The key difference is that the sieve never performs modulo. Instead it assumes primes, and then crosses them off as you find multiples of primes.

def sieve(n):
   """Return a list of the primes below n."""
   prime = [True] * n
   result = [2]
   append = result.append
   sqrt_n = (int(n ** .5) + 1) | 1    # ensure it's odd
   for p in range(3, sqrt_n, 2):
       if prime[p]:
           append(p)
           prime[p*p::2*p] = [False] * ((n - p*p - 1) // (2*p) + 1)
   for p in range(sqrt_n, n, 2):
       if prime[p]:
           append(p)
   return result