Twin Prime Algorithm Optimization

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.

• Python is surprisingly good, although of course C will be typically be faster. At this point concentrating on the algorithms will give you far larger gains. Look at one of rhw's Python SoEs (numpy or non-numpy) that Caridorc references, and modify to print just the twin primes. They are simple, small, and very fast. After doing that, if you still want more speed, look into modification to not sieve modular values that can't be twin primes. – DanaJ Nov 28 '14 at 3:24

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)

• Thanks!! but I must not use any libraries rather then built in functions – newhere Nov 27 '14 at 18:40
• I've added more comments. Also you can easily get rid of the dependency on itertools. – SylvainD Nov 28 '14 at 16:02

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

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.

• Using the Eratosthene sieve is not that convenient here as it is hard to determine a (precise enough) upperbound for the n-th pair of twin numbers. – SylvainD Nov 28 '14 at 10:08
• Maybe you can just create a huge sieve overnight and then just read it from a file. – Caridorc Nov 28 '14 at 13:17