Twin Prime Algorithm Optimization

Question

I have wrote the following algorithm to print out the Nth twin prime numbers:

def isprime(x):
    for i in range(2, (x + 2) / 2):
        if x % i == 0:
            return False
    return True

inputIndex = int(raw_input())        
tempIndex = 1
n=2
while True:
    if isprime(n) & isprime(n + 2):
        if inputIndex == tempIndex:
            print "%d,%d" % (n,n+2)
            break
        else:
            tempIndex += 1
    n += 1

I am looking to optimize the speed and cannot import libraries like math and so or store a list of primes.

Answer 1

A better isprime check

At the moment, your function returns that 0 and 1 are prime, this is not what I have learnt.

A faster isprime check

You check primality by testing divisibility of various candidates between 2 and (x + 2)/2 (non included). You'd better check for candidates between 2 and sqrt(x) (included) as the smallest divisisor will be no bigger than this.

Here's the corresponding code from my toolbox :

def is_prime(n):
"""Checks if a number is prime."""
if n < 2:
    return False
return all(n % i for i in range(2, int(math.sqrt(n)) + 1))

Itertools is good

Instead of going through n by incrementing it, you can use itertools.count

for n in itertools.count(2):
    if isprime(n) & isprime(n + 2):
        if inputIndex == tempIndex:
            print "%d,%d" % (n,n+2)
            break
        else:
            tempIndex += 1

A bit of math

You know that twin primes will be odd (if you are not convinced of this, just think about it). Therefore, you can limit yourself to the search of odd numbers : for n in itertools.count(3, 2):.

Don't check things twice

At the moment, the primality of the same number is checked twice. A way not to do so would be to generate prime numbers and to check this to generate twin pairs.

def yield_odd_primes():
    for n in itertools.count(3, 2):
        if isprime(n):
            yield n

def yield_twin_pairs():
    p1, p2 = 0, 0
    for p in yield_odd_primes():
        p1, p2 = p2, p
        if p1 + 2 == p2:
            yield p1, p2

inputIndex = 1500
for i, (p1, p2) in enumerate(yield_twin_pairs()):
    if i + 1 == inputIndex:
        print "%d,%d" % (p1, p2)
        break

Better user experience

When I ran your code for the first time, nothing happened and I thought Python was just crunching numbers. After a while, I had a look at the code and realised that the script was waiting for some input. By adding a prompt, you make this much clearer :

inputIndex = int(raw_input('Please enter an index: '))

Code organisation

At the moment, when one imports your code, he gets the prompt and everything. The usual way to do is to put the code "actually doing something" behind an if-main guard.

More code organisation

Your can split your code into smaller function. I found quite useful writing the following function to take the nth element of an iterable :

def get_nth_element(iterable, index):
    for i, e in enumerate(iterable):
        if i == index:
            return e

Tests

For such an algorithmic oriented problem, it is easy and a good habit to write test to ensure you do not break anything as you go through modifications.

The hard part is sometimes to rewrite your code in such a way that testing is made easy. It is now the case with my final version of the code.

Final code

import math
import itertools

def isprime(n):
    if n < 2:
        return False
    return all(n % i for i in range(2, int(math.sqrt(n)) + 1))

def yield_odd_primes():
    for n in itertools.count(3, 2):
        if isprime(n):
            yield n

def yield_twin_pairs():
    p1, p2 = 0, 0
    for p in yield_odd_primes():
        p1, p2 = p2, p
        if p1 + 2 == p2:
            yield p1, p2

def get_nth_element(iterable, index):
    assert index >= 0
    for i, e in enumerate(iterable):
        if i == index:
            return e

if __name__ == "__main__":
    inputIndex = int(raw_input('Please enter an index: ')) - 1
    print(get_nth_element(yield_twin_pairs(), inputIndex))

Some math to make the code faster

At the moment, we check the primality of each and every odd numbers. We can do better by considering division by 6 : all numbers can be written :

6k + 0 -> divisible by 2 & 3
6k + 1 -> potential prime
6k + 2 -> divisible by 2
6k + 3 -> divisible by 3
6k + 4 -> divisible by 2
6k + 5 -> potential prime

Thus, except for the obvious pair (3, 5), the only way two prime numbers can be separated by two are if they can be written (6*k + 5, 6*k + 7).

Once you have this, the code pretty much writes itself :

def yield_twin_pairs():
    yield (3, 5)
    for i in itertools.count(5, 6):
        if isprime(i) and isprime(i+2):
            yield (i, i+2)

Answer 2

A faster isprime check

Checking even numbers before and then only looping over only odd numbers makes this twice as fast.

def isprime(n):
    if n < 2:
        return False
    if n % 2 == 0 and n != 2:
        return False
    for i in range(3,int(math.sqrt(n)),2):
        if n % i == 0:
            return False
    return True

An even faster isprime check

All primes (past 2 and 3) are of the forms 6n+1 and 6n-1.

def isprime(n):
    # The same as above

primes = []
for i in range(6,MAX,6):
    if isprime(i+1):
        primes.append(i+1)
    elif isprime(i-1):
        primes.append(i-1)

An extremly fast prime checker

If you need to make several test in a row, the Sieve of Erastothenes is the best method

MAX = 1000 # depends on your needs
def eratosthenes(n):
    # Taken from Rosetta Code
    multiples = set()
    for i in range(2, n+1):
        if i not in multiples:
            print(i)
            multiples.update(range(i*i, n+1, i))

primes = eratosthenes(MAX)
def is_prime(n,primes):
    return n in primes

If you are really really desperate about performance look here for faster-than-you-can-imagine implementations of the sieve of Eratostenes.