e05_distribution_of_primes.sagews

#prime_pi(x) = the number of primes up to and including x
prime_pi(10)

4

is_prime(12345678910987654321)

True

is_prime(1111111111111111111)

True

#primes up to 3,000,000
data = [x for x in [2..3000000] if is_prime(x)]

data

def prime_block_count(block_size, start):
    L = [start..start+block_size-1]
    #print L
    data = [x for x in L if is_prime(x)]
    return len(data)

prime_block_count(10,1)

4

prime_block_count(1000,2000000)

69

plot(prime_pi,2,3000000)

x

3000000

plot(lambda x: prime_pi(x)/x,2,3000000)

#primes in each of the 3000 blocks of length 1000 upto 3,000,000
data = [(x,prime_block_count(1000,1000*x-999)) for x in [1..3000]]

scatter_plot(data)

plot(log,2,3000000)

x

3000

plot(lambda x:1/log(x),2,3000000)

plot(prime_pi,2,3000000)+plot(lambda x: x/log(x),2,3000000, color="red")

plot(lambda x: prime_pi(x)*log(x)/x,2,3000000)+plot(1,2,3000000,color="red")

#try value sof C in the following expression to get a "better" fit:
#plot(prime_pi,2,3000000)+plot(lambda x: x/(log(x)-C),2,3000000, color="red")
plot(prime_pi,2,3000000)+plot(lambda x: x/(log(x)-.5),2,3000000, color="red")

plot(prime_pi,2,3000000)+plot(lambda x: x/(log(x)-.75),2,3000000, color="red")

plot(prime_pi,2,3000000)+plot(lambda x: x/(log(x)-.85),2,3000000, color="red")

plot(prime_pi,2,3000000)+plot(lambda x: x/(log(x)-.95),2,3000000, color="red")

plot(prime_pi,2,3000000)+plot(lambda x: x/(log(x)-1),2,3000000, color="red")

def points(y):
    return scatter_plot([(x,x+1) for x in [1..y]])

points(10)