Willans' Formula
Tags: algorithms, computerscience, mathematics, python
Published: 2022-12-01 19:15:14

In todays post we will be implementing [Willans’ Formula](https://www.jstor.org/stable/3611701). This is an algorithm for finding the nth prime number in the set of prime numbers and was first described by Willans in his 1964 paper _"On Formulae for the nth Prime Number"_ and the first version of the formula provided is as follows:

{% katex %}
{\displaystyle p_{n}=1+\sum _{i=1}^{2^{n}}\left\lfloor \left({\frac {n}{\sum _{j=1}^{i}\left\lfloor \left(\cos {\frac {(j-1)!+1}{j}}\pi \right)^{2}\right\rfloor }}\right)^{1/n}\right\rfloor }
{% endkatex %}

Don't panic at the sight of the formula, at its core it consists of basic arithmetic and trigonometric functions plus the numbers 1 and 2 and the constant [PI](https://en.wikipedia.org/wiki/Pi). Pseudo code for the formula could be presented as follows:

```
fn X'(i, n) -> floor (Y'(i, n) ^ (1 / n))

fn Y'(i, n) -> n / sum ( floor ( cos ( ((j - 1)! + 1 / j) * PI ^ 2 ) ) )

n = 5
result = sum ( map i in range(1, power(2, n)) => 1 + X'(i, n) )
```

> **Sidenote 1:**
>
> This algorithm is not efficient, as per [wikipedia](https://en.wikipedia.org/wiki/Formula_for_primes#Formulas_based_on_Wilson's_theorem):
>
> _In addition to the appearance of `(j-1)!`, it computes `p ^ n` by adding up `p ^ n` copies of 1; for example, `p ^ 5 = 1+1+1+1+1+1+1+1+1+1+1+0+0...+0 = 11`_

Willans also presented us alternative formulas, we will be implementing the version shown above but be aware that due to its inefficiency and terrible runtime complexity, you may not be able to generate primes higher than those in a double digit order. Some of the other formulas presented to us by Willans are more efficient but for now we shall focus purely on formula number 1.

> **Sidenote 2:**
>
> For reference, [you can see each version of Willans formula here](http://www.icodeguru.com/Embedded/Hacker's-Delight/109.htm).

## Tests

For brevity we will test the numbers from 0 through 10, below are a collection of [pytest](https://pytest.org) tests covering these cases. Thus, in `test_willans.py` we have the following:

```python
from main import willans
from pytest import raises

def test_0Throws():
    with raises(ZeroDivisionError):
       willans(0)

def test_1stPrimeIs2():
    assert willans(1) == 2

def test_2ndPrimeIs3():
    assert willans(2) == 3

def test_3rdPrimeIs5():
    assert willans(3) == 5

def test_4thPrimeIs7():
    assert willans(4) == 7

def test_5thPrimeIs11():
    assert willans(5) == 11

def test_6thPrimeIs13():
    assert willans(6) == 13

def test_7thPrimeIs17():
    assert willans(7) == 17

def test_8thPrimeIs19():
    assert willans(8) == 19

def test_9thPrimeIs23():
    assert willans(9) == 23

def test_10thPrimeIs29():
    assert willans(10) == 29
```

## Implementation

In `main.py` we have the following implementation of Willans formula:

```python
from math import cos, factorial, floor, pi


def X(i: int, n: int) -> int:
    y = Y(i, n)
    e = 1 / n
    r = pow(y, e)

    return floor(r)


def Y(i: int, n: int) -> float:
    def iterator(j: int) -> int:
        f = factorial(j - 1) + 1
        a = f / j
        b = a * pi
        c = cos(b)
        e = pow(c, 2)

        return floor(e)

    cs = [iterator(j) for j in range(1, i + 1)]

    return n / sum(cs)


def willans(n: int) -> int:
    xs = [X(i, n) for i in range(1, pow(2, n) + 1)]
    s = sum(xs)

    return 1 + floor(s)


lower_bound = 1
upper_bound = 6
result = [willans(n) for n in range(lower_bound, upper_bound + lower_bound)]

print(result)  # [2, 3, 5, 7, 11, 13]
```

> **Sidenote 3:**
>
> For the implementation of the algorithm I use [pythons built in math module](https://docs.python.org/3/library/math.html) but this has issues with values over a certain size when conducting division operations on large integers and floats.
>
> It would be worth using tools such as [pandas](https://pandas.pydata.org/) or [scikit learn](https://scikit-learn.org/stable/) if you wish to improve the efficiency and scale of the implementation.

## Conclusions

I originally believed, as I am sure most people do, that there was no uniform algorithm or formula for generating the nth prime number but Willans does indeed stand up to scrutiny thus far. That said, it is so inefficient and exponentially slow as an algorithm that practically it is useless even if correct.

I came across this algorithm a few weeks ago when watching this fantastic video by [Eric Rowland](https://ericrowland.github.io/) and I highly recommend giving it a watch.

{% embed https://www.youtube.com/watch?v=j5s0h42GfvM %}

After watching, I did further research as outlined in the links above and came up with the implementation of the formula seen in this article over that time. I find such algorithms and formulae to be utterly fascinating and perhaps in the future I will write more mathematics based algorithms in this series but for now, the article is at an end.

I hope that you found some value in todays post and if you have any questions, comments or suggestions, feel free to leave those in the comments area below the post!