[This is the eighth in a series of posts on the *decadic numbers* (previous posts: A curiosity, An invitation to a funny number system, What does "close to" mean?, The decadic metric, Infinite decadic numbers, More fun with infinite decadic numbers, A self-square number).]

In my previous post, we found a decadic number

with the curious property that it is its own square.

We also discovered an algorithm for computing : starting with , we square each and take the remainder (that is, we keep only the last digits) to give us . Then is defined as the limit of the as .

(As an aside, notice why we are allowed to define as a limit: because of the funny decadic metric, the distance between successive is actually getting smaller by a factor of each time; they are *converging* on something, and that something is .)

So,

We can pretty easily turn this into a computer program to compute to as many digits as we like. As usual, I’m using my favorite programming language, Haskell. It doesn’t matter if you don’t know Haskell, I’ll try to explain things well enough that you can still follow along anyway.

First we just need to import something we’ll need later.

```
> module Decadic where
>
> import Control.Monad.State
```

Now, we define `u1`

to be 5, and the function `nextU`

takes and as input and produces .

```
> u1 :: Integer
> u1 = 5
>
> nextU :: Integer -> Integer -> Integer
> nextU n un = (un^2) `mod` (10^(n+1))
```

Now we can define the (infinite!) list of , which we’ll call `us`

. `generateUs`

takes the current values of and and generates a list with as its first element, and the remainder of the list computed by a recursive call to `generateUs`

(with and as its inputs). `us`

just gets things started by calling `generateUs`

with and .

```
> us :: [Integer]
> us = generateUs 1 u1
> where generateUs n un = un : generateUs (n+1) (nextU n un)
```

Let’s see if it works by asking for the first ten elements of the list `us`

:

```
*Decadic> take 10 us
[5,25,625,625,90625,890625,2890625,12890625,212890625,8212890625]
```

Seems good!

However, all is not roses. Let’s think about what our program had to do in order to calculate . Given , it first squared it to get … but wait a minute! ends in —of course it does, since that’s the defining property of the . But that means that *half* (more or less) of the work needed to square was completely wasted! There’s no point in computing the last half of , since we already know it is going to be . We only care about the rest of it.

Realizing that our program is doing too much work is one thing; turning this intuition into actual improvements is quite another! Initially it was far from obvious to me how to avoid the repeated work. But I finally figured it out.

The key is to remember at each step not just , but . However, since we know that ends with , we can remember along with the *prefix* (call it ) of that comes before the ending . For example, for we’ll remember and (since ). That is, we will always insist that

.

The key question now is: given and , how do we compute and ? Well, computing is easy: it’s just with an extra digit tacked on, namely, the final digit of . In particular, we can define to be the last digit of ; then .

It’s computing which is the tricky part. Since we want to have , the idea is to expand out , manipulate it into the desired form, and see what we get for :

Now, note that we can break into two parts: the last digit (which we’re calling ) and the rest. That is,

where denotes rounding down to the nearest integer, so means to divide by 10 and throw away the remainder—that is, without its final digit. Continuing,

and there we have it! We have expressed in the form . In particular,

Notice that is actually guaranteed to be an integer, since is divisible by 5—this is the same phenomenon we ran across in proving that our simple algorithm for computing works correctly. So we can factor out a division by ten (division is generally slow so we’d rather only do it once):

And now for some code! As one final optimization, instead of keeping track of , we keep track of , so we don’t have to keep recomputing every iteration.

First, the data structure we use to keep track of the current values of , , and :

```
> -- State for incrementally constructing u_n.
> -- Invariant: curT = 10^n; un^2 = pn*curT + un
> data UState = UState { pn :: Integer
> , un :: Integer
> , curT :: Integer
> }
> deriving Show
```

The initial state for has , , and .

```
> -- u_1 = 5; 5^2 = 25 = 2*10 + 5
> initUState = UState 2 5 10
```

The `uStep`

function is the heart of the algorithm: it takes the current values of , , and and updates them to , , and . It also returns as its result, which we can use to build up digit-by-digit.

```
> uStep :: State UState Int
> uStep = do
> u <- gets un
> p <- gets pn
> t <- gets curT
>
> let d = p `mod` 10 -- next digit
> u' = d * t + u -- prepend the next digit to u
> p' = (p + 2*d*u + d*d*t) `div` 10 -- see above proof
>
> put (UState p' u' (10*t)) -- record the new values
>
> return $ fromInteger d -- return the new digit
```

So, did we gain anything? As a concrete comparison, let’s see how long it takes to compute . Using our first, simple code, it takes 7.2 seconds:

```
*Decadic> :set +s
*Decadic> length . show $ us !! 10000
10001
(7.23 secs, 208746872 bytes)
```

(I just had it print the *number of digits* of to avoid wasting a ton of space printing out the entire number.) And using our new and hopefully improved code?

```
*Decadic> length . show . un . flip execState initUState
$ replicateM_ 10000 uStep
10001
(1.55 secs, 225857080 bytes)
```

Only 1.5 seconds! Nice!

The other nice thing about `uStep`

is that it spits out one digit of at a time, which we can use to define as an (infinite) list of digits—as if the digits were coming one at a time down a "tube", a -tube, one might say… get it, a … tube… heh… never mind.

```
> type TenAdic = [Int]
>
> u :: TenAdic
> u = 5 : evalState (sequence $ repeat uStep) initUState
```

```
*Decadic> reverse $ take 20 u
[9,2,2,5,6,2,5,9,9,1,8,2,1,2,8,9,0,6,2,5]
```

Nifty! Next time the *real* fun begins, as I show off some cool things we can do with our shiny new implementation of .

Pingback: Computing with decadic numbers | The Math Less Traveled

Superb analysis on how to work with decadic numbers, from a computing standpoint. I just passed the link to my colleagues. Nice work.