Shortest prime number generator we could come up with?

Question

I worked with a friend on the following prime number generator. Our goal was to make it as short and as fast as possible. This is what we came up with:

def p(n): # will print out primes between 0 and n
  s=[True]*(n/2)
  for i in range(int((n/2-1)/2) >> 1):
    for j in range((i*(i+3)<<1)+3,n/2,(i<<1)+3): s[j]=False
  return [2] + [((i<<1)+3) for i in range(n/2) if (s[i])]
print p(input()) # input number for n

Can it be made any faster? Also, it crashes due to memory errors if numbers larger than 10¹⁰ are inputted. Can this be fixed to make it more memory efficient?

Answer 1

0. Quick comments on your code

1. Clearer function names and variables in your code. What does p mean, and what is s? I can guess, but please.

2. Include a docstring or some indication of what your code does

3. Use the if __name__ == "__main__": module.

I did some quick speed comparisons of your function which you can see below.

magic1            21.9095528103 ms
magic2            26.5589244423 ms
Joseph Farah    1030.95786404   ms

I did 10 runs with primes up to \$2\$ million. Now note the big gap between the optimized functions and yours. Even if you were to speed up your code tenfold you would still be 10 times slower than the other functions.

You can add as many horses you want to your cart, but it will never beat a ferrari.

With this in mind i propose three solutions

1. Implement a better algorithm
2. Use Import
3. Give up

1. New algorithm

Wikipedia can shed some light on how to find larger primes.

For the large primes used in cryptography, it is usual to use a modified form of sieving: a randomly chosen range of odd numbers of the desired size is sieved against a number of relatively small primes (typically all primes less than \$65,000\$). The remaining candidate primes are tested in random order with a standard probabilistic primality test such as the Baillie-PSW primality test or the Miller-Rabin primality test for probable primes.

So in short you use for an example Wheel factorization to find possible primes, then use the Rabin_primality_test to test whether the candidates are primes. Without going too technical the Rabin_primality_test is a probabilistic test. If the function returns True p might be prime, if False we are guaranteed it is not. If n is composite then the Miller–Rabin primality test declares n probably prime with a probability at most \$4^{−k}\$. A friend told me banks used tested primes 80 times. The chances of a composite number passing \$80\$ tests is \$4^{-80}\$ which is in incredibly low number.

A comparison of many different sieves and ways to find primes can be found here https://stackoverflow.com/questions/2068372/fastest-way-to-list-all-primes-below-n. I used the fastest one on my machine primesfrom2to to generate the results above

import numpy as np
def magic1(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)]

2. Use Import

Python is magical in the sense that you can find libaries which pretty much do anything. The best part is that installing new packages is easy as \$\pi\$. There are many libaries for primes, I will mention a few below

• primesieve
• primefac
• pyprimesieve

primesieve] is Python bindings for the primesieve C++ library. It can generate list of primes at an incredible speed, and also has fast nthprime implementation. primefac is incredibly good at factoring primes, it can also generate primes, but is rather slow. pyprimesieve has functions such as prime_sum and has many nifty specialized functions.

So if you want a golfed answer the best you can get in terms of speed is

from primesieve import generate_primes

if __name__ == '__main__':

    print generate_primes(2*10**6)

3. Give up

Even if you are in your prime finding and generating large primes are incredibly hard. However this is a good thing as primes are incredibly important in internet security. If by some miracle you were to find an incredibly, incredibly fast prime generator it would break all the large banks and make payments over the net unsafe.

So perhaps it is a good thing your code is slow ;)

Answer 2

Shorter code != Better code

Honestly, if I saw code like this in a production code base, I would immediately reformat it at least for readability.

I understand that writing shortest possible code can be fun (and there's even a Code Golf site for that), but it is not necessarily something you want to do for serious code, and since you asked for a review of it, here goes.

1. Spacing: There is virtually no white space other than what is absolutely required by the Python compiler to work. Lines like for j in range((i*(i+3)<<1)+3,n/2,(i<<1)+3): s[j]=False just make the code very hard to read and understand.

2. Naming: As N3buchadnezzar mentioned, clearer naming would make this code much better. Meaningful naming is a key skill in writing programs that others (and your future self) can maintain effectively.

3. Documentation: Again as mentioned, short docstring / documentation to describe what this function is for is just a good habit to make in any programming language.

Edit: I removed the function argument type hint as this is a Python 3.5 addition (was originally def get_primes(n: int)) and would not compile with 2.7.

Here is your code, reformatted with that in mind:

def get_primes(n): 
    """
    Find all prime numbers between 0 and n
    Args:
        n (int): The input number to find primes up to.
    """
    primes = [True] * (n / 2)
    for i in range(int((n / 2 - 1) / 2) >> 1):
        for j in range((i * (i + 3) << 1) + 3, n / 2, (i << 1) + 3): 
            primes[j] = False
    return [2] + [((i << 1) + 3) for i in range(n / 2) if (primes[i])]

# test function
if __name__ == "__main__":
    print get_primes(input())

Answer 3

If your primary goal is the shortest program. Then you can use regular expressions for a much shorter program to find primes.

So in python it looks like:

def is_prime(n):
    import re
    for i in range(n):
       if re.match(r'1?$|(11+?)\1+$', '1' * i) is None:
           print i