Why is the Leibniz method for approximating pi so inefficient

Question

I've been playing around with algorithms for computing pi. One that I noticed is the leibniz algorithm.

It states that pi can be approximated like this

$n = 1$ (The first odd number)

$$\pi = 4 \left(\frac{1}{n} - \frac{1}{n+2} + \frac{1}{n+4} - \frac{1}{n+6} + \cdots \pm \frac{1}{n + x} \right)$$

The nature of this algorithm is to start with the first odd number 1. Add 1 divided by that odd number, then subtract 1 divided by the next odd number. Flip the signs on every iteration. Finally, multiply this result by 4.

I created an algorithm for this in JavaScript.

function leibnez(n) {
    var fraction = 0;
    var denominator = 1;
    var sign = "+";
    for (var i = 1; i <= n; i++) {
        sign = i % 2 == 0 ? "-" : "+";
        if (sign == "+") {
            fraction += (1 / denominator);
        }
        else {
            fraction -= (1 / denominator);
        }
        denominator += 2;
    }
    return fraction * 4;
}

When calling leibniz(10000) which will perform 10000 such fractional additions/subtractions, I get a value of 3.14149 which is accurate only to 4 decimal places.

Calling leibniz(100000) results in a value of 3.14158 which is accurate only to 5 decimal places.

Is this algorithm really that inefficient at approximating pi or is the discrepency an error in JavaScript's floating point addition.

Answer 1

Due to conditional convergence, the error term is comparable to the last term of the series:

$$\pi=4 \left[ 1-\frac{1}{3}+\ldots+\frac{(-1)^{n-1}}{2n-1} \right]+(-1)^{n} \left( \frac{1}{n}-\frac{1}{4n^{3}}+\ldots+\frac{E_{2k}}{2^{2k}n^{2k+1}} +\ldots \right)$$

Modern calculations

BBP algorithm $$\pi =\sum_{n=0}^{\infty} \frac{1}{16^{n}} \left( \frac{4}{8n+1}-\frac{2}{8n+4}-\frac{1}{8n+5}-\frac{1}{8n+6} \right)$$

Ramanujan series

$$\frac{1}{\pi} = \frac{\sqrt{8}}{9801} \sum_{n=0}^{\infty} \frac{(4n)!}{(n!)^{4}} \frac{(1103+26390n)}{396^{4n}}$$

See also Borwein's algorithms for further interests.

Answer 2

The Leibniz series is OK for calculating $\pi$ to reasonable precision, if Euler acceleration is applied to it. The idea is to replace $$a_{2n+1}-a_{2n+2}+a_{2n+3}-a_{2n+4}+\cdots$$ with $$\frac12a_{2n+1}+\frac12(a_{2n+1}-a_{2n+2})-\frac12(a_{2n+2}-a_{2n+3})+\frac12(a_{2n+3}-a_{2n+4})-\cdots$$ For a convergent series, this leaves the sum the same, but the terms are smaller because $$\frac1n-\frac1{n+1}=\frac1{n(n+1)}$$ So they decrease with $\frac1{n^2}$ instead of $\frac1n$. Now, since $$\frac1{n^2}-\frac1{(n+1)^2}=\frac{2n+1}{n^2(n+1)^2}$$ We can repeat the process, replacing $$\frac12(a_{2n+1}-a_{2n+2})-\frac12(a_{2n+2}-a_{2n+3})+\frac12(a_{2n+3}-a_{2n+4})-\frac12(a_{2n+4}-a_{2n+5)}+\cdots$$ with $$\frac12\cdot\frac12(a_{2n+1}-a_{2n+2})+\frac12\left(\frac12(a_{2n+1}-a_{2n+2})-\frac12(a_{2n+2}-a_{2n+3})\right)-\frac12\left(\frac12(a_{2n+2}-a_{2n+3})-\frac12(a_{2n+3}-a_{2n+4})\right)+\cdots$$ After $k$ repititions starting from the $n^{th}$ term, each term which started out of order $\frac1n$ is now of order $\frac{k!}{2^kn^{k+1}}$ To compute to quad precision, I chose $n=100$ and $k=24$ so as to be of order about $3.70\times10^{-34}$. Here is my program.

program leibniz
   implicit none
   integer, parameter :: qp = selected_real_kind(33,4931)
   integer, parameter :: Nsum = 100, Nacc = 24
   integer, parameter :: N = Nsum+Nacc
   real(qp), parameter :: one = 1
   real(qp) a(N)
   integer i

   a = [((-1)**(i-1)*one/(2*i-1),i=1,N)]
   do i = 1, Nacc
      a(Nsum+i:N) = [a(Nsum+i)/2,(a(Nsum+i:N-1)+a(Nsum+i+1:N))/2]
   end do
   write(*,*) 4*sum(a)
end program leibniz

And it printed out

  3.14159265358979323846264338327951

The last digit seems to be off by one.