Prime Number Theorem Python

Question

I am attempting to support the prime number theorem by using the attached code. I would like to somehow show that the average gap between prime numbers less than n is log(n). The code that I have right now can determine whether a certain number is prime, and then the second part calculates the prime gap for each consecutive prime in my range. Any ideas for how to proceed in python using the following code?

from pylab import *

def Is_Prime (n):
    if n==1:
        return False
    if n == 2 or n == 3:
        return True
    if n == 4:
        return False
    if n%2 == 0 or n%3 == 0:
        return False

    for i in range(5,int(n**0.5)+1,6):
        if n%i == 0 or n%(i+2) == 0:
            return False

    return True

k = 0
for i in range (1,100):
    if Is_Prime(i) == True:
        print(i)
        k+=1
print "Total number of prime numbers in [1,100] is", k


previous = 2
n = 0
for i in range(3,100000):
    if Is_Prime(i):
        n = n+1
        current  = i
        gn = current - previous
            print gn
        plot(n,gn,'rs')
        xlabel('n')
        ylabel('g(n)')
        previous = i
        if n == 100:
            break

show()

Answer 1

My suggestion of how to do this (in relation to your current code) is as follows:

previous = 2
n = 0
for i in range(3,10000000):
    if Is_Prime(i):
        n = n+1
        current  = i
        gn = current - previous
        previous = i
        if n % 1000 == 0:
            print "Now at the gap between prime",n,"and prime",n+1
            average_gap = (i-2)/float(n);
            # (n+1)st prime is i, 1st prime is 2. n gaps.
            # float(n) so that we get float, not int division!
            print "Average gap so far:", average_gap
            print "ln(p)/average_gap:", math.log(i) / average_gap
        if n % 100000 == 0:
            break

Note that by the 10000th prime, we are only at a ratio of 1.1 ish... and by the 100000th prime, we're still at 1.0831.

It takes a longgg time to get close to 1. This is partly because it actually more closely follows the reciprical of the Li(n) function: http://en.wikipedia.org/wiki/Logarithmic_integral_function (- as the gap length around n, is approximately log(n)... So the total density of primes below n is more appropriately the integral of this... But for large n, Li(n) ~ log(n), so it works out).

(See the first graph on the right of http://en.wikipedia.org/wiki/Prime_number_theorem - the blue line represents the ratio you'll get, the purple line represents the ratio you'd get using the logarithmic integral)

As another aside, (and as I suggested in my comment), it'd also be a lot faster to have a look at using a Sieve to generate your primes - eg see here: http://en.wikipedia.org/wiki/Sieve_of_Eratosthenes ...

But the above code works (at least for me)! Anyway, best wishes with your future prime endeavours!