Babylonians used floating point numbers

I am a history enthusiast. I listen to history podcasts, watch history videos and historical dramas. But I love when I get a piece of media that combines both history and STEM.

This is the case of Knuth’s “Ancient Babylonian Algorithms” which I found while reading the collection titled “Selected Papers on Computer Science” where it is presented in chapter 11.

The paper turned out to be interesting because it offered new interesting insights beside the already known facts about that ancient land. The big surprise came at the end of the first page: while it is known that Babylonians used a base 60 system – from which we get the number of minutes in an hour, the degrees in a circle, etc. –, I did not know that they also used floating point arithmetic (in base 60, nonetheless!).

What does this mean? The collection \(\{2, 20\}\) (that is, the juxtaposition of the symbols for \(2\) and \(20\) from the next image) can be interpreted as either meaning \(2 \times 60 + 20 = 140\) but also as \(2 + 20/60 = 2.3333..\). In fact, it can mean any number of the form \(140 \times 60^n\) where \(n\) is any integer, positive or negative.

How do we know this? From various tables of reciprocals going back as far as 2250 B.C. These tables transform the problem of dividing by a number into the equivalent one of multiplying by its inverse. A scribe needing to divide by, say, \(81\) (represented as \(\{1, 21\}\)) would instead multiply by the value \(\{44,26,40\}\) taken from such a reciprocal table². In fact, the scribe wouldn’t be doing the multiplication either, he would have used one of the multiplication tables laying around – the tables multiplying \(\{44,26,40\}\) with any integer \(k\) with \(1 \le k \le 20\) were commonplace. To get the right exponent of the powers of \(60\), the scribe would need to estimate the quantity, but this is easy. In fact, this is the algorithm that people used in the past to perform computations with slide rules – I had one of these laying around when I was young but thought it was very confusing, learned much later how it was supposed to work.

Let’s try to encode converting from a Babylonian number to ours. First, we can define a BabylonNumber to be a collection of integers:

type BabylonNumber = [Int]

If floating point were not being used, it would be easy to convert to base 10:

convert0 :: BabylonNumber -> Integer
convert0 = go 0
  where
    go n [] = n
    go n (x:xs) = go (n * 60 + x) xs

And indeed, this works for \(\{2,20\}\) but we won’t get \(1/81\) when we try \(\{44,26,40\}\):

ghci> convert0 [2, 20]
140
ghci> convert0 [44, 26, 40]
160000

This is to be expected, since we know that floating point numbers are being used. So, we can implement a convertAny function that would return an infinite list of all possible values. How will this work? We can generate the base value from convert0 and then construct two infinite lists where one contains numbers each being \(60\) times smaller than the previous and the other one contains numbers increasing by a factor of \(60\). Then, we just need to intercalate these 2 geometric series, the same way we prove that \(|\mathbb{N}| = |\mathbb{Z}|\):

convertAny :: BabylonNumber -> [Double]
convertAny ns = base : intercalate smaller larger
  where
    base = fromIntegral $ convert0 ns
    smaller = tail $ iterate (/60) base
    larger = tail $ iterate (*60) base

intercalate :: [a] -> [a] -> [a]
intercalate (a:as) (b:bs) = a:b:intercalate as bs
intercalate as bs = as ++ bs -- base case for finite lists

Now we can indeed test that our expected numbers show up:

ghci> take 5 $ convertAny [2,20]
[140.0,2.3333333333333335,8400.0,3.888888888888889e-2,504000.0]
ghci> 1/81
1.2345679012345678e-2
ghci> take 10 $ convertAny [44, 26, 40]
[160000.0,2666.6666666666665,9600000.0,44.44444444444444,5.76e8,0.7407407407407407,3.456e10,1.2345679012345678e-2,2.0736e12,2.0576131687242798e-4]

But, how do we know which number from the stream we actually mean? As said, this is based on estimating the quantity and using this to guess the exponent (i.e., search in the infinite stream). In our case, we could use the following:

estimate :: BabylonNumber -> Double -> Double
estimate ns ref = head . filter inRange $ convertAny ns
  where
    inRange x = let z = x / ref in 1 / 60 < z && z < 60

Which would produce desired results:

ghci> estimate [2,20] 100
140.0
ghci> 1/81
1.2345679012345678e-2
ghci> estimate [44, 26, 40] 0.01
1.2345679012345678e-2

Implementing addition, multiplication, etc. are left to the reader, together with the obvious improvements that would be needed. We’ll return back to the paper.

A similar proof that floating point was being used can be extracted from the YBC 7289 clay tablet, pictured below:

We see that the tablet represents a square and 3 different numbers. It is clear that \(\{30\}\) represents the side of the square. If we interpret the other numbers as being floating point in base 60, then \(\{1, 24, 51, 10\}\) could be

\[ 1 + \frac{24}{60} + \frac{51}{60^2} + \frac{10}{60^3} = 1.41421\ldots \approx \sqrt{2} \]

and the second one could be

\[ 42 + \frac{25}{60} + \frac{35}{60^2} \approx 30\sqrt{2} \]

That is, the first number gives us an approximation for \(\sqrt{2}\) whereas the second one gives the value of the length of the diagonal. Observe how the first number has an implicit factor of \(60^{-3}\) whereas the second one a factor of \(60^{-2}\) compared to the side of the square.

We can validate that with our code too, using corresponding estimates:

ghci> 30 * sqrt 2
42.42640687119285
ghci> estimate [42, 25, 35] 30
42.426388888888894
ghci> estimate [1, 24, 51, 10] 1
1.4142129629629632

I would stop here, however the original article has a few more topics of interest. So I’ll just present them from a very high level. First, the paper continues by giving several example of how “programming” was being used in those ancient times. There are tablets that describe procedures to solve simple problem where the scribe would just replace the input numbers, follow the text exactly – like a sermon, including when needing to, for example, multiply by 1 –, write down each partial result and end with the this is the procedure formula.

Still, these are straight-line programs – programs without any branching. There are no tablets containing if statements. When the algorithms are different (e.g., division by \(7\) needs to use a different procedure given the repeating nature of the expansion of \(1/7\)) the scribes would use different tablets. There are, however, a few cases where there seems to be some loop unrolling: there are a few tablets (e.g., VAT 8528) that solve problems related to compound interest. The same procedure is repeated on them several times until the answer is obtained.

Finally, the articles discusses how the impressive AO 6456 tablet containing reciprocals of numbers up to the 6th decimal place could have been computed. There is a theory that involves sorting many large numbers to interleave sequences, a task that now computers can do easily. You’ll probably enjoy reading this paper too.

Next time, we’ll jump in time to the most computationally advanced “brain” of our times.