www.pi314.net Boris Gourévitch The world of Pi - V2.57 modif. 13/04/2013 Pi-Day in
Home

Simon Plouffe

## Finding the th digit of a number a without needing to know the previous one thanks to the BBP formula

### 1 The BBP formula

Little reminder: On the 19 septembre 1995 at 0h29 (!), after month on research in the dark, Simon Plouffe, David Bailey and Peter Borwein de Vancouver discover the apperently simple and innocent formula
 (1)

now called the BBP formula (see Plouffe's page). Our three canadians noticed that this expression is very close to a decomposition in base 16 of . In the 1996 article [1] that is now famous, they show that you can use 1 to find the th digit of in base (all the more in base 16) without having worked out the previous digist, all of this in roughly a linear time ) and in space !! The space required being so small, they apply straight away this result by calculating the billionth digit of in base 2.

Let's stop now working on the calculation of the th digit of by considering the general case of a number and an associated BBP formula:

### 2 Finding the -th digit of a number without knowing the previous one

We will follow here, with a few adaption, the very clear explenation of Xavier Gourdon and Pascal Sebah on the subject [2] which already simplified the original article [1]. X. Gourdon knows a lot on claculation algorithm since he is the current holder (for quite a few years) of the speed record in calculating the decimals of with the program PiFast (even faster than the program made by Kanada's team).

To simplify the reasoning, we will restrain ourself to base , and the th bit of defines the th bit of the fraction part of , that is the part after the dot. The fraction part of a number shall be written as .

#### 2.1 Shift

The first idea relies on the following theorem

Theorem 2.1 The th bit of is obtained by calculating the th bit of i.e. the th bit of .

Proof. Decomposing in base

 (2)

where . The th bit of is therefore which as wanted is the th bit of .  __

Example 2.2 Let's find the 13th bit of . Since we have , in particular, . All that is left now is to work out the bit of . Changing into binary, we get

 (3)

hence the bit of is , and the bit is .

Finding the th bit of is therefore the same as finding the first bit of . Using the BBP formula 1, this correspond to the first bit of the serie

 (4)

where and are the sums

 (5)

Taking into account the fast decrease of for , we just need to evaluate the first few terms of (in floating point form!), plus exacly enough for the rest of the sum is less than the presicion needed (here bit after the dot).

Concerning , we can put its general term in the form

 (6)

with and in . By writting the euclidien division of by , we can see that , . This can also be written as

 (7)

#### 2.2 Exponentiation

Imagine that we want to calculate the bit of . Thanks to the previous principles, this comes down to calculate the bit of and so the bit of . Using 7, we need to find such that hence evaluate .

Calculating is done by a method called binary exponentiation (or binary powering), whose idea seems to have been already used before 200 BC [3], and even considered by the Egyptians to carry out multiplication! This simply comes down to using modulo arithmetic, in other words working in , in three steps :

1. We decompose in powers of :
2. We calculate each steps the power of untill we reach :
3. We use step 1 et 2 : .

More generally to obtain , we set to 1 and we take the biggest power of less than . Then we apply the following algorithm that is equivalent to our example :

1. If , Then ;  ; EndIf
2. If Then ; Go to step 1. ; EnfIf.

Finnaly we end up with to take into account . Due to the decomposition of in powers of , we use operations for this calculation, which is very weak !

#### 2.3 Final Steps

All that is left is to obtain and then to sum evrything on to obtain in 5. We need to note that if the calculation of every term of is done in floating point form and that, say, the last bit is wrong, then in the worst case scenario the of the remainders on the whole sum will produce an error on bits. For the values of non-excessive, we can in all calmness carry out the calculations of with an arithmetic precision of only bits.

For this reason and by noting that for , we sum termes each needing operations due to the exponentiation, the algorithm requires operations and in space, where is the complexity of the multiplication of two integers of size bits.

### References

[1]   D.H. Bailey, P.B. Borwein, S. Plouffe, “On the Rapid Computation of Various Polylogarithmic Constants”, Mathematics of Computation, 1997, vol. 66, pp. 903-913.

[2]   X. Gourdon, P. Sebah, “N-th digit computation”, preprint, 2003, http://numbers.computation.free.fr/Constants/Algorithms/nthdigit.pdf.

[3]   D. E. Knuth, “The Art of Computer Programming. Vol. 2: Seminumerical Algorithms”, Addison-Wesley, Reading, MA, 1981.