Havens' Haven

I have never done anything “useful”. No discovery of mine has made, or is likely to make, directly or indirectly, for good or ill, the least difference to the amenity of the world.

— G. H. Hardy

# Christopher Havens

If you were facing a 25 year stretch for murder and had landed in solitary for a while, you might think to yourself, “Maybe this would be a good time to catch up on some number theory”. Something like this occurred to Christopher Havens when he was convicted of murdering Randen Robinson and sent to the Washington State Monroe Correctional Complex northeast of Seattle. He had dropped out of high school, began using drugs, and killed Robinson during a drug-related dispute. After his conviction, bad behavior sent him to solitary confinement for a year.

While in solitary, Havens started working puzzles and then realized that the guards were handing out math problems to other prisoners. He did as much math as they could provide, quickly learning calculus and number theory, but eventually exceeded the guards’ ability to keep up with his level.

That’s when he sent this letter that eventually arrived on the desk of Matthew Cargo, editor for Mathematical Sciences Publishers, whose partner Marta Cerruti is Associate Professor of Materials Engineering at McGill University.

Her parents are mathematicians and her father, Umberto Cerruti, is a number theorist at the University of Torino in Italy. Umberto Cerruti sent Havens a problem to work on, and Christopher wrote back with a long, complicated solution that turned out to be correct.

Christopher has since published a paper jointly with Cerruti and two other mathematicians, in the journal Research in Number Theory, “Linear fractional transformations and nonlinear leaping convergents of some continued fractions”. He hopes to inspire other prisoners through mathematics and has started the Prison Mathematics Project (PMP) to connect inmates with mathematicians. They hold bi-weekly meetings to discuss math problems and work on larger projects. Beginning in 2016 the PMP has been holding Pi Day celebrations.

“You should never set aside your dreams and goals and dismiss them as impossible just because you are institutionally challenged,” Havens said. “Today is not just a day to celebrate pi, but mathematics as a whole and celebrate an opportunity we have to move our lives in a positive direction.”

Besides having a research paper published, Chris also submitted a problem to Math Horizons, a magazine published by the Mathematical Association of America (MAA) meant to be accessible to undergraduate mathematics majors. The problem is this, What is the smallest positive integer $y$ such that $1729y^2 + 1$ is a perfect square? This is a form of Pell’s Equation, and there’s a connection to the mathematicians Srinivasa Ramanujan and G. H. Hardy.

Ramanujan was a self-taught mathematician from India who produced thousands of new ideas in number theory. In early January 1913, he sent a letter with some of his results to Hardy who was at Cambridge at the time. Later, Hardy said that the letters were “certainly the most remarkable I have received” and that Ramanujan was “a mathematician of the highest quality, a man of altogether exceptional originality and power”. He invited Ramanujan to visit him in England, and Ramanujan arrived in April 1914.

Ramanujan (22 December 1887 – 26 April 1920) was ill most of his short life and Hardy visited him while he was in the hospital.

I remember once going to see him when he was ill at Putney. I had ridden in taxi cab number 1729 and remarked that the number seemed to me rather a dull one, and that I hoped it was not an unfavourable omen. “No,” he replied, “it is a very interesting number; it is the smallest number expressible as the sum of two cubes in two different ways.”

$$1729 = 1^3 + 12^3 = 9^3 + 10^3$$

The number $1729$ is now known as a taxicab number or the Ramanujan-Hardy number and is the number Havens chose for the problem. The drawing at the top of Hardy and Ramanujan was done by Chris Havens.

# Pell’s Equations

Pell’s equations are in the form

$$x^2 - ny^2 = 1$$

for integers $x,y$, and $n$. The problem that Chris Havens proposed is to find $x$ and $y$ such that $x^2 = 1729y^2 + 1$, which can be rearranged into a Pell equation,

$$x^2 - 1729y^2 = 1.$$

In the 7th century, Indian mathematician and astronomer Brahmagupta first studied equations in this form (he called it the Varga-Prakriti equation), and Pell had nothing to do with them. This is an example of Stigler’s Law of Eponymy, named for University of Chicago statistics professor Stephen Stigler, which says that no scientific discovery is named after its original discoverer. Stigler’s Law is itself an example, having been earlier described by sociologist Robert K. Merton. Leonard Euler misattributed the equations to Pell rather than to the several earlier mathematicians who worked on similar equations including Diophantus, Bhāskara II, Narayaṇa Paṇḍita, William Brouncker, and Pierre de Fermat.

The equations are hyperbolic functions and have solutions at points where both $x$ and $y$ are integers. This is an example of the equation $x^2 - 2y^2 = 1$ plotted in Desmos. Use x^2-2y^2=1 to plot the function, and plot solutions in the form $(x,y)$ as shown on the second line.

In the 12th century, Bhāskara II (1114–1185) continued the work of Brahmagupta (598–670) and developed an algorithm for finding solutions that he called the Chakravala Method.

# The Chakravala

Bhāskara named his method Chakra (meaning “wheel”) because of its iterative nature. Shailesh Shirali is a mathematician at the Sahyadri School in India who wrote “The Chakravala Method - Zeroing in on a Solution” in the mathematics journal At Right Angles. This description of the Chakravala is from his paper.

The first step to understanding the Chakravala comes from the Diophantus-Brahmagupta-Fibonacci identity,

$$\begin{aligned} (a^2+b^2)(c^2+d^2) &= (ad-bc)^2 + (ac+bd)^2 \\ &= (ad+bc)^2 + (ac-bd)^2. \end{aligned}$$

Brahmagupta extended this to

$$\begin{aligned} (a^2-nb^2)(c^2-nd^2) &= (ac+nbd)^2 - n(ad+bc)^2 \\ &= (ac - nbd)^2 - n(ad-nbc)^2. \end{aligned}$$

It’s not clear how these identities were discovered, but they are easy to verify by expanding and collecting terms. You can start to see how the Chakravala works for finding integer solutions to the equation $x^2-ny^2 = k$ if we already have two solutions,

$$\begin{aligned} a^2-nb^2 = k_1 \\ c^2-nd^2 = k_2 \end{aligned}$$

then

$$(ac+nbd)^2 - n(ad+bc)^2 = k_1k_2.$$

To complete the method, we need Bhāskarā’s lemma which says that if $a^2-nb^2 = k$, then

$$\left( \frac{ma+nb}{k} \right)^2 - n \left( \frac{a+bm}{k} \right)^2 = \frac{m^2-n}{k},$$

which can be easily verified with expansion and collection of terms.

Given a non-square positive integer $n$, we’re looking for a solution to $x^2 - ny^2 = k$. If we have two solutions $[a,b,k_1]$ and $[c,d,k_2]$ a third solution $[u,v,k_1 \times k_2]$ can be generated from the Brahmagupta identity using

$$\begin{aligned} u &= ac + nbd \\ v &= ad + bc. \end{aligned}$$

Suppose $n=10$ and we have two solutions $[7,2,9]$ and $[11,3,31]$ meaning that

$$7^2-10 \times 2^2 = 9 \text{ and } 11^2 - 10 \times 3^2 = 31.$$

Using the equations for $u$ and $v$,

$$\begin{aligned} u &= 7 \times 11 + 10(2 \times 3) = 77 + 60 = 137 \\ v &= 7 \times 3 + 2 \times 11 = 21 + 22 = 43, \end{aligned}$$

gives a new solution $[137,43,9 \times 31] = [137,43,279]$,

$$137^2 - 10 \times 43^2 = 18769 - 18490 = 279.$$

If we already have one triple $[a,b,k_1]$ that works, a second one can be derived by using $(m)^2 - n \times 1^2 = m^2-n$ giving the triple $[c,d,k_2] =[m,1,m^2-n]$ and calculating a new $u$ and $v$,

$$\begin{aligned} u &= ma + nb \\ v &= a + bm \end{aligned}$$

so

$$u^2 - nv^2 = (ma+nb)^2 - n(a+bm)^2 = k_1(m^2 - n).$$

The trick to finding a solution for $x^2 - ny^2 = 1$ is to start with an “almost” solution and iterate using the Chakravala until $k(m^2-n) = 1$. Finding $[a,b,k_1]$ is done by letting $b=1$ and then setting $a$ to the nearest integer greater than $\sqrt{n}$ so that $a^2 - n = k_1$ gives the smallest positive integer $k_1$. From Bhāskarā’s, lemma we need $a + bm$ to be evenly divisible by $k_1$ and $m^2$ to be as close to $n$ as possible. The steps in the algorithm are,

• Find an initial nearby solution to $a^2 - nb^2 = k_1$ where $a$ and $b$ have no common divisors
• Look for $a+bm$ divisible by $k$ and then select the $m$ that minimizes $|m^2-n|$
• Get new values $a',b'$ and $k'$ from

  $$\begin{aligned} \\ a' &= \frac{ma+nb}{|k|} \\[1.5em] b' &= \frac{a+bm}{|k|} \\[1.5em] k' &= \frac{m^2-n}{k}. \end{aligned}$$

Repeat the second and third steps until $k' = 1$.

Does the Chakravala ever end, or does the wheel keep on turning? Fortunately, the algorithm converges, that is, it eventually gets to $k=1$. The proof is given by A. Bauval in “An elementary proof of the halting property for the Chakravala Algorithm”.

But to make the Chakravala work requires understanding modular arithmetic.

# Modular Arithmetic

In number theory it’s often useful to describe a number by its remainder when divided by another number, so if

$$a = qn+r$$

then $n$ goes into $a$ $q$ times with a remainder of $r$. An example is the angle when going around a circle:

If the circle has a radius of $1$, then the angle $\theta$ ranges from $0$ to $2 \pi$, but once it gets back to the $x$-axis the remainder resets to $0$. That is, you can’t tell from looking at where the pointer is whether it’s gone around zero times, once, twice, or more. The sloped lines in this plot represent the remainder of the angle after dividing by $2 \pi$. Usually, the numbers $a,q,n$, and $r$ are integers and the equation is often written as

$$a = r \text{ mod}(n)$$

meaning that the remainder when $a$ is divided by $n$ is $r$, or “$a$ equals $r$ modulo $n$”. The $q$ is missing in the equation because it’s the integer number of times that $n$ divides $a$, and we’re only interested in the remainder $r$.

If $a = r_1 \text{ mod}(n)$ and $b = r_2 \text{ mod}(n)$ then $a+b = r_1 + r_2 \text{ mod}(n)$, since

$$\begin{aligned} a &= q_1 n + r_1 \\ b &= q_2n+r_2 \\ a+b &= (q_1 n + r_1) + (q_2 n + r_2) \\ &= (q_1 + q_2) n + (r_1 + r_2) \\ &= r_1 + r_2 \text{ mod}(n) \end{aligned}$$

For example, $14 = 4 \text{ mod}(5)$ and $17 = 2 \text{ mod}(5)$ so $14+17 = 31 = 1 \text{ mod}(5)$. Of course, adding $4+2 = 6$ so you have to apply the modulus again to get $6 = 1 \text{ mod}(5)$. Subtraction works the same way with $a - b = (q_1 n + r_1) - (q_2 n + r_2)$.

You can multiply numbers modulo $n$. Using the same notation from above,

$$\begin{aligned} a \times b &= (q_1 n + r_1) \times (q_2 n + r_2) \\ &= q_1 q_2 n^2 + q_1 r_2 n + q_2 r_1 n + r_1 r_2 \\ &= (q_1 q_2 n + q_1 r_2 + q_2 r_1)n + r_1 r_2 \\ &= r_1 r_2 \text{ mod}(n) \end{aligned}$$

So, $14 \times 17 = 238 = 3 \text{ mod} (5)$ and $4 \times 2 = 8 = 3 \text{ mod}(5)$.

The only arithmetic operation left is division. Since $10 = 4 \text{ mod}(6)$ it would be convenient if we could divide both sides by $2$ to get $\frac{10}{2} = 5 = \frac{4}{2} = 2 \text{ mod}(6)$, but clearly, that’s wrong, $5 \neq 2 \text{ mod}(6)$. The division operation that we’re used to means that if

$$\frac{x}{y} = z$$

then

$$x = yz.$$

If the solution to dividing $x$ by $y$ is $z$, then multiplication reverses this, so $y$ times $z$ must equal $x$. We need an equivalent for modulo division. Suppose there are two numbers $z_1$ and $z_2$ such that multiplication by $y$ gives the same answer modulo $n$,

$$z_1y = z_2y \text{ mod}(n).$$

This means that $z_1y = q_1 n + r_1$ and $z_2y = q_2n + r_2$ for some $q_1,q_2,r_1$ and $r_2$. Since $z_1y = z_2y \text{ mod}(n)$ then $r_1 = r_2$, so

$$\begin{aligned} z_1y &= q_1n = q_2n = z_2y \\ z_1y - z_2y &= (z_1 - z_2)y = kn \end{aligned}$$

for some integer $k$. We can find a solution to the division problem using the greatest common divisor of $n$ and $y$.

The greatest common divisor (gcd) of two integers is the largest integer that divides both of them. To find the gcd, find the prime factors of each, and then take the product of all the common factors,

$$\begin{aligned} 60 &= 2^2 \times 3 \times 5 \\ 18 &= 2 \times 3^2 \\ gcd(60,18) &= 2 \times 3 = 6 \end{aligned}$$

Suppose the $gcd(n,y) = d$ and that $d \neq y$. The second condition means that $y$ does not divide evenly into $n$. Divide $n$ by $d$ to get some integer $g = \frac{n}{d}$. Now $g$ is also a divisor of $n$ since $d = \frac{n}{g}$, but $g$ can’t be a divisor of $y$ because if it were, then both $d$ and $g$ would be divisors of $n$, as would the product $dg$ meaning that $d$ wasn’t really the greatest common divisor. So, if $g = \frac{n}{d}$ doesn’t go evenly into $y$, then it must divide $z_1 - z_2$, since

$$(z_1-z_2)y = kn.$$

This means that

$$z_1y = z_2y \text{ mod}(n)$$

only when

$$z_1 = z_2 \text{ mod}\left( \frac{n}{d} \right).$$

So, a unique solution exists only when the $gcd(n,y) = 1$. In this case, $n$ and $y$ are said to be relatively prime or coprime. Let’s go back to the original problem of finding $y^{-1} = \frac{1}{y}$ such that $x y^{-1} = z \text{ mod}(n)$. One way to find $y^{-1}$ would be to check every element of the set ${1,2, \ldots, n}$ until we find a solution. But there’s a faster way using the greatest common divisor method discovered by Euclid.

# The Euclidian Algorithm

Euclid’s algorithm is a more practical method for finding the greatest common divisor of two numbers than listing their prime factors and then calculating the product of their common factors. Suppose you want to find $gcd(a,b)$ where $a > b$. Now, $a = q_1b + r_1$, and $b = q_2r_1 + r_2$ for some $q_1,q_2,r_1,r_2$. Write $r_1$ in terms of $r_2$: $r_1 = q_3r_2 + r_3$. Maybe you find that $r_3 = 0$, meaning that $r_2$ is a divisor of $r_1$. Then,

$$b = q_2(q_3r_2)+r_2 = (q_2q_3+1)r_2$$

and

$$\begin{aligned} a &= q_1b+r_1 = q_1(q_2q_3+1)r_2 + r_1 \\ &= q_1(q_2q_3+1)r_2 + q_3r_2 \\ &= (q_1(q_2q_3+1) + 1) r_2. \end{aligned}$$

This means that $r_2$ divides both $a$ and $b$ and if we let $w = q_2q_3+1$ then $b = w r_2$ and $a = (q_1w+1)r_2$, so $r_2$ is the largest divisor of both $a$ and $b$.

Using the example from above, we’ll find the $gcd(60,18)$ using the Euclidean Algorithm:

$$\begin{aligned} 60 &= 3 \times 18 + 6 \\ 18 &= 3 \times 6 + 0 \\ &\Rightarrow gcd(60,18) = 6. \end{aligned}$$

We get the answer in just two steps, and no factorization is needed. If $r_3$ isn’t zero, continue the process until some $r_n = 0$.

The Euclidean Algorithm can be extended to solve problems such as

$$ax+by = c.$$

As an example, let’s start by finding the greatest common divisor of $33$ and $27$,

$$\begin{aligned} 33 &= 1 \times 27 + 6 \\ 27 &= 4 \times 6 + 3 \\ 6 &= 2 \times 3 + 0 \\ &\Rightarrow \text{gcd}(33,27) = 3. \end{aligned}$$

Now, suppose we want to solve

$$3 = 33a + 27b.$$

From the gcd calculation above,

$$\begin{aligned} 6 &= 33 - 1 \times 27, \\ 3 &= 27 - 4 \times 6 \\ \end{aligned}$$

Replace the $6$ in the second equation with the one from the first to get

$$\begin{aligned} 3 &= 27 - 4 \times (33 - 1 \times 27) \\ &= -4 \times 33 + (1 + 4) \times 27 \\ &= -4 \times 33 + 5 \times 27. \end{aligned}$$

so $a = -4$ and $b = 5$.

# Modular Division

Going back to the division equation above,

$$\begin{aligned} \frac{x}{y} &= z \\ x &= y z \\ \Rightarrow x &= y z + k n \end{aligned}$$

In other words, if the solution to dividing $x$ by $y$ is $z$, then in modular arithmetic terms (modulo $n$), $x$ is $y$ times $z$ plus some constant multiple $k$ of $n$. This is exactly the extended Euclidean algorithm for the gcd that we just did. If $y$ and $n$ are coprime, meaning $gcd(n,y)=1$, then from the extended Euclidean algorithm there are integers $r$ and $s$ such that

$$1 = ry + sn$$

and

$$\begin{aligned} z &= rx + tn \\ k &= sz - ty \end{aligned}$$

for every integer $t$. Here’s the trick - when you calculate $1 = ry + sn \text{ mod} (n)$ the term with $sn$ is zero because $n$ is a divisor. So what’s left is

$$1 = ry \text{ mod}(n)$$

which means that

$$r = \frac{1}{y} \text{ mod}(n)$$

so $r$ is the inverse of $y$.

If $gcd(n,y) = d > 1$ and $x$ is an integer multiple of $d$, then

$$\begin{aligned} z &= r + t \frac{n}{d} \\[1.5em] k &= s - t \frac{y}{d}, \end{aligned}$$

but there won’t be a unique solution modulo $n$ (but there will be modulo $n/d$).

We now have methods for calculating addition, subtraction, multiplication, and under special conditions, division modulo $n$, and with the Chakravala algorithm we can find the solution to any Pell’s equation.

# The Algorithm

To solve the equation

$$a^2-nb^2 = 1$$

for a given positive integer $n$, start with values for $[a,b,k]$ and $gcd(a,b)=1$. An easy choice for $b$ is $1$ since setting $a$ to the next larger integer to $\sqrt{n}$ makes $k$ small.

• Find values of $m$ such that $k$ divides $a + bm$, and select the $m$ that minimizes $|m^2 - n|$.
• Recompute

  $$\begin{aligned} a' &\leftarrow \frac{ma+nb}{|k|} \\[1.5em] b' &\leftarrow \frac{a+bm}{|k|} \\[1.5em] k' &\leftarrow \frac{m^2-n}{k} \end{aligned}$$

If $k' = 1$ we have a solution, otherwise find a new $m$ satisfying the constraint and iterate until $k' = 1$.

Finding $m$ such that $k$ divides $a+bm$ means that there is some $j$ such that

$$a + bm = jk$$

or

$$bm = -a \; \text{mod}(k)$$

Using the extended Euclidean algorithm, we can find $\frac{1}{b} \text{ mod } k$ and an integer $t$ such that

$$t = \left(-a \frac{1}{b} \right) \text{ mod}(k).$$

The $m$ that satisfies $a + bm = 0 \text{ mod}(k)$ and minimizes $|m^2 - n|$ will have a remainder of $t$ when divided by $k$, so it’s just a matter of finding the multiple $j$,

Let $n = 17$ and start with $[a,b,k] = [5,1,8]$, so $5^2 - 17 \times 1^2 = 8$. Then we want $5 + 1 \times m = 8j$. This is equivalent to

$$1 \times m = -5 = 3 \; \text{mod}(8).$$

When $b=1$ the problem is trivial because $m$ has to be in the set $\{3,11,19,27, \ldots \} = \{8-5,16-5,24-5,32-5, \ldots \}$. For $n=17$, $|m^2-n|$ is minimized when $m=3$. Update $a,b,$ and $k$,

$$\begin{aligned} a' &\leftarrow \frac{ma+nb}{|k|} = \frac{3 \times 5 + 17 \times 1}{|8|} = 4 \\[1.5em] b' &\leftarrow \frac{a+bm}{|k|} = \frac{5 + 1 \times 3}{|8|} = 1 \\[1.5em] k' &\leftarrow \frac{m^2-n}{k} = \frac{3^2 - 17}{8} = -1 \end{aligned}$$

Calculating a new value for $m$, we need $4 + 1 \times m$ to be divisible by $k = -1$, and so any value of $m$ will do, and the value that minimizes $|m^2 - n|$ is $4$. Updating the values for $a,b$, and $k$ give $a = 33$, $b = 8$ and $k = 1$, so the algorithm stops with

$$\begin{aligned} 33^2 - 17 \times 8^2 &= 1 \\ 1089 - 17 \times 64 &= 1089 - 1088 = 1. \end{aligned}$$

Using the PARI/GP language, the code chakravala.gp finds the smallest values of $a$ and $b$ for a given $n$. This is a screenshot of the solution to Chris Havens’ problem with the Ramanujan-Hardy Taxicab number 1729:

# The Prison Mathematics Project

The Prison Mathematics Projects likes to hear from people willing to volunteer time to help inmates improve their math skills. From the website,

In early 2012, a prisoner by the name Christopher Havens began studying mathematics for the very first time within the prison inside of prison. Christopher was in a form of isolation that we call solitary confinement… prisoners know it as “the hole”. Having no other form of stimulation, mathematics occupied every hour of Christopher’s days, lasting the better part of a year until he was released back into the general population. In so many cases, the negative effects of prolonged isolation manifests itself through the behaviors in both the short and the long term. However, this case is one where a mix of isolation and the transformative powers of mathematics caused Christopher to undergo a steady chain of personal growth while igniting within him a passion for mathematics.

If you have some ability in mathematics and would like to volunteer, you can submit your name and email address, and the PMP will connect you with an inmate to become their “pen pal”. Donations are also greatly appreciated.

# Reach out!

If you have any questions or comments regarding this article you can visit this thread on GitHub discussions or email us directly