List of Primes

Computers are good at finding lists of prime numbers. And there are plenty of ways of doing so.

First a simplist approach:

#!/usr/bin/python3

count=0
i=2

while count<1000:
    prime=True
    for f in range (2,i):
        if (i%f==0):
            prime=False
            break
    if prime:
        count+=1
        print(count,i)
    i=i+1

This finds the first thousand primes by repeatedly trying to divide each integer by every smaller integer, and printing those for which no divisor is found. (The % operator gives the remainder when the left-hand side is divided by the right-hand side.)

It is correct, but a little slow. It takes about 1.25s to find the first thousand primes, 5.7s to find the first two thousand primes, and 26s to find the first four thousand.

But we know that, with the exception of two, all primes are odd, and no odd number has an even divisor, so both the outer and inner loop can run over odd numbers only.

#!/usr/bin/python3

print(1,2)
count=1
i=3

while count<1000:
    prime=True
    for f in range (3,i,2):
        if (i%f==0):
            prime=False
            break
    if prime:
        count+=1
        print(count,i)
    i=i+2

Times are now 0.7s to 1,000, 2.6s to 2,000 and 12.8s to 4,000. And, importantly, it gives the same result.

But we can do much better. As a number's factors will always appear in pairs, one of the pair must be no bigger than the square root of the number being factorised. So we can replace the start of the loop over f, our trial factors, with

      for f in range (3,int(sqrt(i))+1,2):

which also requires us to add as the second line of our code

from math import sqrt

Now the times are much better: under 0.2s to to 1,000, under 0.2s to 2,000, and about 0.3s to 4,000. This code can print the first ten thousand primes in under a second. (And it gives the 10,000th prime as 104,729, which is correct.)

If one has sufficient ambitious to wonder about the first 100,000 primes, then it starts to take some time -- about 22 seconds. But we can do even better.

We could keep a record of all primes found so far, and use just those as potential trial factors.

#!/usr/bin/python3
from math import sqrt

print(1,2)

count=1
i=3
primes=[]

while count<100000:
    prime=True
    for f in primes:
        if (f>sqrt(i)):
            break
        if (i%f==0):
            prime=False
            break
    if prime:
        count+=1
        primes.append(i)
        print(count,i)
    i=i+2

This runs in about 17s, so only about 20% or so faster.

Python vs C

Would it have been better to write this in a traditional compiled language, such as C? A line-by-line translation of the above into C might result in something like the following.

#include<stdio.h>
#include<stdlib.h>
#include<math.h>

int main(){
  int count,i,j,f,prime;
  int *primes;

  printf("1 2\n");
  
  count=1;
  i=3;
  primes=NULL;

  while(count<100000){
    prime=1;
    for(j=0;j<count-1;j++){
      f=primes[j];
      if (f>sqrt((double)i)) break;
      if (i%f==0){
	prime=0;
	break;
      }
    }
    if (prime){
      primes=realloc(primes,count*sizeof(int));
      primes[count-1]=i;
      count++;
      printf("%d %d\n",count,i);
    }
    i+=2;
  }
  return 0;
}

This we can compile and run with

  $ gcc -o primes primes.c -lm
  $ time ./primes

It runs in about 1.4s, so over ten times faster than the python version. One can control how hard the compiler attempts to optimise its own generation of code. If the above is compiled with

  $ gcc -o primes -O primes.c -lm

in order to turn on more optimisation, the run-time drops to under 0.9s.

So for this example C is about twenty times faster than python. For the example given on the Linpack page, their speeds are essentially identical. And for ease of use, particularly when adding simple graphics, such as in the logistic map example, python surely wins.

So neither C nor python is always superior to the other. Sometimes the use of one is more appropriate, and sometimes the other. And there are other possible choices too, such as Fortran...