NTIC First Steps

Section 21.1 First Steps

It might seem at first there is very little we can say about this function; after all, thus far we've seen no particular pattern in the primes themselves (other than that they are nearly all odd). You may wish to see what the function looks like to confirm this sense. It is a not particularly smoothly increasing function with no upper bound (recall Theorem 6.2.1).

Figure 21.1.1. Plot of prime pi function (`plot(prime_pi,2,100)`)

🔗

Sage note 21.1.2. Syntax for counting primes.

The syntax for this function is prime_pi(n).

🔗

Subsection 21.1.1 A funky formula

🔗

Given the skepticism of the paragraphs so far this chapter, you may be surprised to learn there are exact formulas for this function, as well as for the

$n$ th prime. The following formula (for

$n > 3$ ) is one of my favorites (see the Appendix of the exhaustive Hardy and Wright, [C.2.2], and also Exercise 21.5.1):

$\begin{equation*} \pi(n)=-1+\sum_{j=3}^{n}\left((j-2)!-j\left\lfloor\frac{(j-2)!}{j}\right\rfloor\right)\text{.} \end{equation*}$

🔗

Can you see why this is not useful in practice? So there is plenty left for us to discuss.

🔗

On the other hand, it works! We can confirm this by using the following code (non-interactive).

xxxxxxxxxx
 
def primeish(n):
    if n==1:
        return 0
    elif n==2:
        return 1
    elif n==3:
        return 2
    else:
        result = -1
        fact = 1
        for j in range(3,n+1):
            fact = fact*(j-2)
            result += (fact - j*floor(fact/j))
        return result
​
import math
def plotprimeish(n):
    n = int(math.floor(n))
    return primeish(n)
​
pretty_print(html("The number of primes up to 20000 this formula gives is $%s$"%primeish(20000)))
pretty_print(html("The real function in Sage gives $%s$"%prime_pi(20000)))
pretty_print(html("And let's compare plots:"))
plot(lambda x:plotprimeish(x), (x,2,100)) + plot(prime_pi,2,100,color='black')

🔗

Sage note 21.1.3. Cython.

It's possible to significantly speed up many such computations by converting to Cython, a way to take Python/Sage and turn it into the much-faster compiled language C. For a project, try to speed this function up using Cython!

🔗

Sage note 21.1.4. Not all algorithms are equal.

Don't forget that just because an algorithm works, doesn't guarantee it will be useful in practice! However, it's often useful to get something correct first, and only then try to optimize.

🔗

Subsection 21.1.2 A very low bound

🔗

On a more computationally feasible note, one can find a very rudimentary (lower) bound on this function. Recall that unadorned logarithms are the natural log.

🔗

Fact 21.1.5.

There are at least

$\begin{equation*} \left\lfloor\frac{\log(\log(x)/\log(2))}{\log(2)}\right\rfloor+1=\left\lfloor\log_2(\log_2(x))\right\rfloor+1 \end{equation*}$

primes less than or equal to $x\text{.}$

🔗

Proof.

In Saidak's proof [C.7.22] of the infinitude of the primes, he constructs the sequence

\begin{equation*} 2,\; (2)+1,\; (2(2+1))+1,\; (2(2+1)(2(2+1)+1))+1, \ldots \end{equation*}

Then he shows, similarly to Euclid's proof, that there is at least one new prime divisor in each element of the sequence (even if not necessarily a larger one). So the $n$th prime can be no bigger than the $n$th term of this sequence. (See Exercise 21.5.3.)

By induction, we see that this term (and hence the $n$th prime) is less than or equal to $2^{2^{n-1}}\text{.}$

The case $n=1$ is clear, since the first prime is $2\text{.}$
The $n$th term is the previous terms multiplied together, plus 1, which by induction is less than

\begin{equation*} 2^{2^0}2^{2^1}\cdots 2^{2^{n-2}}+1=2^{1+2+4+\cdots +2^{n-2}}+1=2^{2^{n-1}-1}+1\leq 2^{2^{n-1}} \end{equation*}

(this uses the same type of technique as in Subsection 4.5.2).

So when $\pi(x)=n\text{,}$ the $n$th term in the sequence is $2^{2^{\pi(x)-1}}\text{,}$ which can't be less than $n$ itself (the $n$th prime is certainly at least $n$). Take two logs of this to get

\begin{equation*} \log(\log(2^{2^{n-1}}))=\log(2^{n-1}\log(2))=(n-1)\log(2)+\log(\log(2)) \end{equation*}

\begin{equation*} =(\pi(x)-1)\log(2)+\log(\log(2))\geq \log(\log(n))\text{.} \end{equation*}

This yields the given statement¹See also [C.2.1, Corollaries 2.7 and 2.8] for this proof, but connected more directly to Euclid's proof of the infinitude of primes., with the floor function accounting for the fact that $\pi$ takes only integer values.

🔗

As you can see below, this is not a very useful bound, considering there are actually 25 primes less than 100, not three! Each of the inequalities in the proof was in a sense ‘wasteful’. Note also that the floor function is only necessary for

$x\lt 5\text{.}$

Figure 21.1.6. Plot of prime pi versus log log

🔗

Subsection 21.1.3 Knowledge from nowhere

🔗

Finally, although it may not seem evident, you should know that it is not necessary to actually find all the first

$n$ primes (even of a particular type) to compute how many there are, at least not always.

🔗

Definition 21.1.7.

Let $\phi(n,a)$ to be the number of positive integers less than $n$ which are not divisible by any of the first $a$ primes

🔗

Try Exercise 21.5.2 to see how this function works.

🔗

It is possible to develop the recursive formula

$\begin{equation*} \phi(n,a)=\phi(n,a-1)-\phi\left(\left\lfloor\frac{n}{p_a}\right\rfloor,a-1\right)\text{,} \end{equation*}$

🔗

where

$p_a$ denotes the

$a$ th prime, which allows use of a type of inductive argument to compute

$\phi(n,a)$ without having to use many computational resources. It is then not too hard to use a counting argument to prove that

$\begin{equation*} \pi(n)=\pi(\sqrt{n})+\phi(n,\pi(\sqrt{n}))-1\text{.} \end{equation*}$

🔗

This is the typical way to calculate

$\pi$ in software without actually counting primes, and with some speedups it can be quite efficient.

🔗

Interestingly, this is also how one finds the

$n$ th prime. You use an approximation to the

$n$ th prime like

$n\log(n)$ and then check values of

$\pi(n)$ near that point to see where the value changes, which should lead you exactly to the prime you seek. (Recall Sage note 4.2.1 about %time when using the following cell.)

xxxxxxxxxx
 
%time nth_prime(10^7)