Problem of the Week Special:

New Formulas for the New Year

The big news: David Bailey, Peter Borwein, and Simon Plouffe have found a new formula for Pi that allows one to compute, say, the millionth digit of Pi without computing any prior digits, and using NO high- precision software, and using essentially NO memory. Unfortunately it works only in base 16.

Here is their formula.

Pi = Sum[1/16^k * 
          (4/(8k+1) - 2/(8k+4) - 1/(8k+5) - 1/(8k+6)), 
             {k, 0, Infinity}]

Let's abbreviate this as

Pi = (16; (4, 0, 0, -2, -1, -1, 0, 0)).

The length of the vector (8 in this example), determines the denominators: 8k+1, 8k+2, etc.

Here is Mathematica code that implements the BBP algorithm for digit extraction: in essence it comes form looking at the fractional part of 16^d Pi.

p[i_, k_] := 8 k + If[i == 1, 1, i + 2]

PiHexExtractor[d_] := Module[{n = d - 1, mainSum = 0},
  Do[mainSum += Mod[{4., -2, -1, -1} . 
       Table[PowerMod[16, n - k, p[i, k]] / p[i, k], {i, 4}], 1],
    {k, 0, n}];
    
  rd = RealDigits[Mod[mainSum + 
     Sum[{4., -2, -1, -1} . 
       Table[16 ^ (n - k) / p[i, k], {i, 1, 4}], 
             {k, n + 1, n + 15}], 1], 16];

  Take[Join[Table[0, {-rd[[2]]}], rd[[1]]], 5]]
--------------------------------------------

PiHexExtractor[10 ^ 6]

Output =   {2, 6, 12, 6, 5}

These are the millionth and next 4 hex digits (in the BBP paper for a double-check). This takes 6-hours on my PowerMac 8100.

Pause and think about this: No high-precision, no memory consumption. This goes completely against conventional wisdom. Actually, Rabinowitz and Wagon had a paper in the Monthly in early 1995 that presented an algorithm that used no high-precision, and no floating-point arithmetic whatsoever. And it worked in base 10. But it would require a lot of memory for the millionth digit.

Copies of the BBP paper ("On the Rapid Computation of Various Polylogarithmic Constants") may be fetched from Peter Borwein's web site (www.cecm.sfu.ca/~pborwein). It will appear in Mathematics of Computation.

I learned of this work at the Canadian Math Society meeting in Vancouver in early December. Then Victor Adamchik (Wolfram Research, Inc.) and I, after using Mathematica to prove the BBP formula (which, as they show in their paper, reduces to some straightforward integrals), used Mathematica to discover some new formulas.

Many of the following (the logarithms) were known to BBP, but not all. To get used to the notation, look at the Leibniz series for Pi

Sum[(-1)^k 4/(2k+1)];

In current notation it is:

Pi = (-1; (4, 0))

We call this a one-term formula, because there is only one nonzero term in the coefficient vector. Of course, this one is in base 1, and so is useless.

Here is the series for log 2 that comes from the Maclaurin series for log(1 + x).

Log 2 = Sum[2^(-k) 1/k] = (2; (1))

Here is a new three-termer for Pi:

Pi = (-4; (2, 2,1, 0))

Here is one we found, but it turns out to be simply the representation one gets by using the Gregory series: the Maclaurin series for log((1+x)/(1-x)). Presumably the one-termers that one finds will be well-known!

log 3 = (4; (1, 0))

But here are two-termers for log 5 and log 7 and arctan 2.

log 5 = (-4; (2, 0, -1, 0))

log 7 = (8; (0, 3, 0, -3/2, 0, 0))

arctan 2 = (-4; (1, 0, 1/2, 0))

And a 3-termer for Pi * Sqrt[2]:

Pi * Sqrt[2] = (-8; (4, 0, 1, 0, 1, 0))

The following one-termer in base 10 is well-known and was used by BBP to extract far-out decimal digits. It also came out of our methods.

Log[9/10] = (10; (0, -1/5)) (well-known)

A nice feature of our methods is that the discovery comes along with a proof. We obtain huge ugly expressions as values of complicated integrals; these are identities that we then specialize by setting certain coefficients to special values (obtained by solving, say, 8x8 linear systems). Great fun!

Big questions: Is there such a representation for an irrational algebraic number? Sqrt[2]?

What about more logarithms: Log[23] is first unresolved case.

Our paper will appear in Mathematica in Education and Research (Springer) and copies are available via the WRI web site.

[More news about the Bailey-Borwein-Plouffe algorithm can be found at Mathsoft, as part of their Favorite Mathematical Constants Page. -Jeff]

© Copyright 1996 Stan Wagon, except for parenthetical comments. Reproduced with permission.