# a practical prime counting function

Looking at Prime counting functions on Wikipedia, I only found formulas with no hint on how people got there. So, to better understand, I've decided to build one from scratch, starting from a naive kind of Sieves of Eratosthenes. Here is what I've found so far: $$\text{Let } k = \left\lfloor\frac {x - 1} {6} \right\rfloor \text{ and a function } p1(x) \text{ such that: }$$ $$p1(x)=k-\sum _{i=1}^{6i^2-2i\leqslant k} \left( \left\lfloor \frac{k+i}{6i-1}\right\rfloor -i+1 \right) - \sum _{i=1}^{6i^2+2i\leqslant k} \left( \left\lfloor \frac{k-i}{6i+1}\right\rfloor -i+1 \right) + \sum _{i=0}^{i\leqslant \frac{k-29}{25}} \left( \left\lfloor \frac{k+5i+6}{30i+35}\right\rfloor -1 \right) + \sum _{i=0}^{i\leqslant \frac{k-29}{35}} \left\lfloor \frac{k-5i-4}{30i+25}\right\rfloor + \sum _{i=0}^{i\leqslant \frac{k-64}{35}} \left( \left\lfloor \frac{k+7i+13}{42i+77}\right\rfloor -1 \right) + \sum _{i=0}^{i\leqslant \frac{k-57}{49}} \left( \left\lfloor \frac{k-7i-8}{42i+49}\right\rfloor -1 \right)$$ $$\text{and let } k=\left\lfloor \frac{x+1}{6}\right\rfloor \text{ and a function } p2(x) \text{such that}:$$ $$p2(x)=k-\sum _{i=1}^{i\leqslant \frac{k+1}{7}} \left\lfloor \frac{k-i}{6i-1}\right\rfloor + \sum _{i=0}^{i\leqslant \frac{k-21}{25}} \left( \left\lfloor \frac{k+5i+4}{30i+25}\right\rfloor -1 \right) + \sum _{i=0}^{i\leqslant \frac{k-41}{35}} \left\lfloor \frac{k-5i-6}{30i+35}\right\rfloor + \sum _{i=0}^{i\leqslant \frac{k-41}{35}} \left( \left\lfloor \frac{k+7i+8}{42i+49}\right\rfloor -1 \right)$$ $$\text{Then for any } x\leqslant 874, \text{we have:}$$ $$2+p1(x)+p2(x)= \pi (x)$$ Although the result looks like complicated, this is only very simple sums and programming an algorithm that matches it is really straightforward. Here for example is the corresponding Mathematica code:

Manipulate[
Module[{panel, p1, p2, pcount},

p1[x_] :=
Module[{k, i, sum, sum1, sum2, sum3, sum4, sum5, sum6},
k = Floor[(x - 1)/6];
sum1 = (i = 1; sum = 0;
While[6*i^2 - 2*i <= k, sum += Floor[(k + i)/(6*i - 1)] - i + 1;
i++]; sum);
sum2 = (i = 1; sum = 0;
While[6*i^2 + 2*i <= k, sum += Floor[(k - i)/(6*i + 1)] - i + 1;
i++]; sum);
sum3 = (i = 0; sum = 0;
While[i <= (k - 29)/25,
sum += Floor[(k + 5*i + 6)/(30*i + 35)] - 1; i++]; sum);
sum4 = (i = 0; sum = 0;
While[i <= (k - 29)/35, sum += Floor[(k - 5*i - 4)/(30*i + 25)];
i++]; sum);
sum5 = (i = 0; sum = 0;
While[i <= (k - 64)/35,
sum += Floor[(k + 7*i + 13)/(42*i + 77)] - 1; i++]; sum);
sum6 = (i = 0; sum = 0;
While[i <= (k - 57)/49,
sum += Floor[(k - 7*i - 8)/(42*i + 49)] - 1; i++]; sum);
k - sum1 - sum2 + sum3 + sum4 + sum5 + sum6
];

p2[x_] :=
Module[{km, i, sum, sum1, sum2, sum3, sum4, sum5, sum6},
km = Floor[(x + 1)/6];
sum1 = (i = 1; sum = 0;
While[i <= (km + 1)/7, sum += Floor[(km - i)/(6*i - 1)]; i++];
sum);
sum2 = (i = 0; sum = 0;
While[i <= (km - 21)/25,
sum += Floor[(km + 5*i + 4)/(30*i + 25)] - 1; i++]; sum);
sum3 = (i = 0; sum = 0;
While[i <= (km - 41)/35,
sum += Floor[(km - 5*i - 6)/(30*i + 35)]; i++]; sum);
sum4 = (i = 0; sum = 0;
While[i <= (km - 41)/35,
sum += Floor[(km + 7*i + 8)/(42*i + 49)] - 1; i++]; sum);
km - sum1 + sum2 + sum3 + sum4
];

pcount[x_] := 2 + p1[x] + p2[x];

panel =
Panel@Column[{
Grid[{
{"x", x}, {"\[Pi](x)", PrimePi[x]}, {"pcount(x)", pcount[x]}
}, Alignment -> {Left, Top}, Frame -> All]
}];

panel
]
, {{x, 301}, 1, 1500, 1, Appearance -> "Open"}
]


I've found interesting the fact that you don't need to rely on Primes table nor need of trial division.

Although the range $x\leqslant 874$ (in fact, it also works for range $1001\leqslant x\leqslant 1224$) is quite small, I think that the formula can be extended to any range by adding necessary sums, but it then becomes difficult to avoid duplicates counts.

EDIT1: some people misunderstood "by adding necessary sums". I mean in fact that I'm quite sure that $p1(x)$ and $p2(x)$ valid for any $x$ can be written in the form of a closed single double-sum of $k$ (i.e: $pn(x)=\sum \sum f(k)$ ).

• With regards to your last statement there ("I think that the formula can be extended to any range by adding necessary sums"), I doubt that it can be finitely extended to cover an infinite range. This is because it is known that no non-constant polynomial function $P(n)$ with integer coefficients exists that evaluates to a prime number for all integers $n$. See here for more information on this. In addition to that, a similar question can be found here. Jun 22, 2015 at 14:51
• You've hidden a lot of hand-waving in "by adding necessary sums..." Sure, you can fix any formula by adding the parts that you've missed. Jun 22, 2015 at 14:53
• There is an approach substantially faster than Eratosthenes sieve for counting primes $\pi(x)$ up to a natural number $x$, due to Meissel with improvements by D.H. Lehmer (mentioned in the Wikipedia article). For a recent account of "record setting" in this regard, see Douglas Staple (2015), The Combinatorial Algorithm for Computing $\pi(x)$. Jun 22, 2015 at 15:05
• Hans Riesel, "Prime Numbers and Computer Methods for Factorization", 1994 has a very nice 26 page chapter on prime counting methods. I was able to quickly implement the Legendre, Meissel, and Lehmer methods from this, all of which handily beat sieving. Optimizing Lehmer is a larger task, as is LMO and extended LMO. Regardless, I recommend it as a source for a lengthier intro to the material as well as practically useful. Jun 22, 2015 at 17:14

I don't mean to be a pessimist, but the computational complexity of this "formula" is most likely greater than or equal to the complexity of just running the sieve itself. However, there are actually "formulas" for primes, they just have enormous computing complexity. For instance, look here. They also discuss the above issues.

I use "formulas" in quotation because a formula is supposed to vastly simplify the process of finding a solution to a problem. For instance, integration formulas vastly simplify the process of adding an (un)countable amount of things and indefinite amount of times! The formulas mentioned don't do that, so they generally aren't considered formulas.

Basically what you have, translates an algorithm into mathematical formalism.

Quite a few people have sought formulas for primes before; this effort falls into that body of work. You can convert this into a prime test via pcount(x+1)-pcount(x).

Unfortunately, this algorithm does not appear to be simpler than a direct sieve, and does not appear to be more compact than a list of primes (e.g. the storage required to store the algorithm is approximately equal to that required to simply store a list of all primes up to 874). It is not possible to have a single polynomial that yields all primes (this was proven by Legendre). If you allow floor functions (as the OP does), then a single formula suffices (this was proven by Mills).

If the formula were computer-generated, i.e. if 874 were a parameter, which you could turn to 2000 or whatever you like and automatically get the corresponding formula, then THAT part might be interesting. Unfortunately the actual formula presented is of limited interest itself.

• thanks you for your answer. In fact, I'm quite sure that these sums can be written in a close form (i.e: $pn(x)=\sum \sum f(k)$ ). But after reading answers and comments, I don't know if it's worth it. at least, I've learned a lot... Jun 24, 2015 at 13:51