Counting Primes

Math 242 Modern Computational Mathematics

Counting Primes

The prime-counting function π(x) counts the number of primes less than or equal to a real number x. For example π(8)=4 because there are four primes less than or equal to 8. We will compute values of π(x) using a list of primes found by the Sieve of Eratosthenes. Here is a fast version of the sieve:

sieveErat[n_Integer]:=Module[{nums=Range[n],m=Floor[Sqrt[n]]},
nums[[1]]=0;
Do[
If[nums[[i]]≠0,nums[[2i;;;;i]]=0],
{i,2,m}];
Select[nums,#≠0 &]
]

One way to count the primes less than or equal to n is simply to find the length of the list of primes returned by sieveErat[n]. For example:

Length[sieveErat[17]]

This works, since π(17)=7. However, we would like to compute lots of values of π(x). Doing the sieve of Eratosthenes for each value would be inefficient. Here is a more efficient method, implemented as a module countPrimes with the following specification.

Module: countPrimes
Arguments: positive integer n
Return: list containing values π(x), for x=1,2,3,...,n

The module will just do the sieve of Eratosthenes once, finding all primes up to n. Then the module will loop over the integers i from 2 through n, keeping track of how many primes it has found up to i.

countPrimes[n_Integer]:=Module[{primeList,j,i,counts=Table[0,n]},
primeList=sieveErat[n];
j=0; (* tracks position in primeList *)
Do[
If[j<Length[primeList],
If[i==primeList[[j+1]], (* then i has reached the next prime *)
 j++
];
];
counts[[i]]=j  (* there are j primes ≤ i *)
,{i,2,n}];
counts (* return the list of counts *)
]

This is quite fast:

Timing[countPrimes[100000];]

Exercise 1

Compute a few values of π(x) and verify that countPrimes agrees with your calculations. Then use ListPlot to plot π(x) up to x=m, for various positive integers m. What is the shape of the graph? Could you approximate the graph by a smooth function?

Density of Primes

In about 1793, Gauss conjectured that the density of primes near x is approximately 1/ln(x). Gauss had printed tables of prime numbers to work with. You have a computer and fast sieve algorithms. Can you provide evidence for or against Gauss’ conjecture?

To do this, you will need to approximate the density of primes using values of π(x). One way to do this is to fix an interval size d, and then compute the number of primes between x and x+d, divided by the length of the interval. That is:

((π(x + d) - π(x))/d)

Exercise 2

Plot an approximation of the prime density function. To do this, you will need to approximate the density of primes using values of π(x). One way to do this is to fix an interval size d, and then compute the number of primes between x and x+d, divided by the length of the interval. That is:

((π(x + d) - π(x))/d)

Note that our function primeCounts only accepts integer input. You might want to use primeCounts[Floor[x]] so that x can be a real number.

Can you provide evidence in support of Gauss’ conjecture about the density of primes? It may help to experiment with different interval sizes d. If you think Gauss’ conjecture is correct, then what continuous function would approximate the prime counting function π(x)?

Approximating the prime counting function

The prime counting function is not continuous. However, various continuous functions can be used to approximate the prime counting function. A simple approximation is given by:

(x/(ln(x)))

Another approximation is given by the logarithmic integral:

li(x)=(∫_0^x)(dt/(ln t))

In Mathematica, the logarithmic integral is implemented by LogIntegral.

Exercise 3

Plot the prime counting function, together with the approximations given above, on different intervals. Do you think these are good approximation? What is the percent error for each approximation at various values of x? Does the percent error decrease when x gets large?

Twin Primes

A pair of primes that differ by 2 are called twin primes. For example, 11 and 13 are twin primes, as are 29 and 31.

The twin prime counting function is denoted (Subscript[π, 2])(x). Specifically, (Subscript[π, 2])(x) is the number of primes p≤x such that p+2 is also prime.

Exercise 4

Write a function that computes (Subscript[π, 2])(x) according to the following specification.

Module: countTwinPrimes
Arguments: positive integer n
Return: list containing values (Subscript[π, 2])(x), for x=1,2,3,...,n

How does the graph of (Subscript[π, 2])(x) compare with the graph of π(x)?

Other Pairs of Primes

A pair of primes that differ by 4 are called cousin primes. A pair of primes that differ by 6 are called sexy primes (a pun on the Latin word for six).

Exercise 5

How do the counts of cousin primes and sexy primes compare to the counts of twin primes? Is one type of prime pair more common than the others? Investigate!

Counting primes mod 4

Interesting patterns emerge when considering the distribution of primes modulo an integer k. For example, consider the primes mod 4. Except for the prime 2, all primes are congruent to either 1 or 3 mod 4.

Exercise 6

Write a module to count the number of primes congruent to 1 mod 4. Then do the same for 3 mod 4. (You may choose to write a single module that does both, via an argument.)

You could also experiment with counting primes mod other integers, such as k=5 or k=6. What do you notice? Do you have any conjectures?