Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

I don't have any experience with these types of formulae and am finding it difficult understanding how to use Bellards formula.

Say I wanted to get the $3^\text{rd}$ digit of $\pi$ (which equals 1), would I simply replace all the instances of $n$ in the formula for $1$, and then multiply by $\frac{1}{2^6}$ ?

If I wanted to get more than one digit, do I get them all without multiplying by $\frac{1}{ 2^6}$, sum and then multiply?

I have googled the problem extensively and the only resources that are available seem to explain how the formula was gotten and not its use.

share|cite|improve this question
I take it when you write "digit" you really mean "bit". We are considering binary, not decimal, expansions, right? – Gerry Myerson Feb 6 '12 at 23:15
up vote 7 down vote accepted

You should see how BBP formulas are used to extract the $p$-th digit.

The Bellard's formula, when truncated to $m$ terms will give you a rational number $r_m$, approximating $\pi$.

Suppose we are interested in finding $p^\text{th}$ binary digit of $\pi$. Choose $m$ is sufficiently large, so that $| \pi -r_m | \leqslant 2^{-p}$. In order to extract the digit, we need to compute $d_p = \lfloor 2^p r_m \rfloor \bmod 2$.

Consider $$ 2^{p-1} r_m = \sum_{n=0}^\infty (-1)^{n} 2^{p-1-10n} \Big( \frac{4}{10 n+1}-\frac{1}{10 n+3}-\frac{1}{2 (4 n+1)}-\frac{1}{2^4 (10 n+5)}-\frac{1}{2^4 (10 n+7)}+\frac{1}{2^6 (10 n+9)}-\frac{1}{2^6 (4 n+3)} \Big) = \Sigma_1 - \Sigma_2 - \Sigma_3 - \Sigma_4 - \Sigma_5 + \Sigma_6 - \Sigma_7 $$ Now the trick is not to compute the digits which precede the $p$-th one. To this end one splits the sum into 7 pieces. We consider just one: $$ \Sigma_1 = \sum_{n=0}^\infty \frac{(-1)^{n} 2^{1 +p-10n} }{10n+1} $$ Not computing the preceding digits is equivalent to only retaining the fractional part. This is done by computing the numerator modulo the numerator: $$ \sum_{n=0}^\infty \frac{(-1)^{n} 2^{1 +p-10n} \bmod (10 n+1) }{10n+1} $$

Here is a little Mathematica code for Bellard's algorithm:

Sigma[p_Integer, a_Integer, b_Integer, c_Integer, prec_] := 
 Module[{sum = 0, prev, den, n = 0, num, exp, 
   bprec = 1 + prec Log2[10.]},
   den = a n + b;
   exp = c + p - 1 - 10 n;
   If[exp > 0,
    num = Mod[PowerMod[2, exp, den], den, 1],
    num = 1; den = den*2^(-exp);
   prev = sum;
   sum += Divide[(-1)^n SetPrecision[num, prec], den];
   If[sum - prev == 0 && exp < -bprec, Break[]];

PthPiDigitBellard[p_Integer?NonNegative, prec_: MachinePrecision] := 
 Floor[2 Mod[Sigma[p, 10, 1, 2, prec] - Sigma[p, 10, 3, 0, prec] -
     Sigma[p, 4, 1, -1, prec] - Sigma[p, 10, 5, -4, prec] - 
     Sigma[p, 10, 7, -4, prec] + Sigma[p, 10, 9, -6, prec] - 
     Sigma[p, 4, 3, -6, prec], 1]]

In[166]:= {Table[PthPiDigitBellard[k, 20], {k, 1, 60}], 0} == 
 RealDigits[Pi, 2, 60, -1]

Out[166]= True
share|cite|improve this answer

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.