Faßt man je zwei aufeinanderfolgende Glieder der Leibnizreihe zu einem zusammen, entsteht aus (2) die Reihe

 π/8 = 1/(1·3) + 1/(5·7) + 1/(9·11) + .... ;                       (4)

sie entspricht der Reihe des Inders NILAKANTHA SOMAYAJI (geb. 1444):

        1        1         1

π/8 = ------ + ------ + ------- + ... ,                            (4a) 

      2² - 1   6² - 1   10² - 1

der lange vor Leibniz (1646-1716) und Gregory (1638-1675) lebte und hundert Jahre alt wurde. (Die Reihe (4a) wird übrigens in der Formelsammlung von [3] falsch zitiert: dort steht sie mit alternierenden Vorzeichen.)

Sowohl die Leibnizreihe (2) wie die aus ihr durch Umformung gewonnenen Reihen (4) bzw. (4a) sind, so schön sie aussehen, für die praktische Berechnung von π leider gänzlich ungeeignet. Dies erkennt man z. B. an (2): Das 1000. Glied dieser Reihe ist ungefähr gleich 0,0005; d. h., um nur vier sichere Nachkommastellen von π zu erhalten, muß man mindestens 10000 Glieder der Leibnizreihe berücksichtigen!

Die Leibnizreihe liefert deshalb so schlechte Ergebnisse, weil die Arkustangensreihe (3) für x=1 nur noch äußerst langsam konvergiert. (Daß sie es überhaupt tut und nicht divergiert, ist nicht selbstverständlich und in der Hochschulmathematik Gegenstand eines eigenen Beweises.) 

Wählt man das Argument x in (3), absolut gesehen, kleiner als 1, konvergiert die Reihe schneller, und dies kann man auf andere Weise für die angenäherte Berechnung von π ausnutzen. Es gilt nämlich nicht nur π/4 = arctan 1, sondern nach dem Schweizer Mathematiker EULER (1707-1783) auch: 

π/4 = arctan1/2 + arctan1/3.

Hiermit läßt sich die Kreiszahl auf eine vorgegebene Anzahl von Stellen sehr viel leichter berechnen als mit der Leibnizreihe. Mit der von EULER ebenfalls gefundenen Formel

π/4 = 5·arctan1/7 + 2·arctan3/79 ,

in der die Argumente der Arkustangensterme noch kleiner sind als in der vorhergehenden, erhielt er nach [3] zwanzig Stellen von π in weniger als einer Stunde – ein gewaltiger Fortschritt gegenüber dem, was oben über Ludolf van Ceulen berichtet wurde!

Bei den π-Berechnern war bis zur Erfindung des Computers und noch eine Weile danach die Formel von MACHIN (1680-1751):

π/4 = 4·arctan1/5 – arctan1/239                               (5)

eine der beliebtesten. Bei ihr ist das Argument des zweiten Arkustangensterms besonders klein. Sie spielt auch bei einem der folgenden Programme eine Rolle. Im übrigen bestehen die meisten anderen, auf Arkustangensbasis zur angenäherten Berechnung von π geeigneten Formeln aus drei oder mehr Termen.

Im Jahre 1995 wandten sich die beiden Amerikaner STANLEY RABINOWITZ und STANLEY WAGON wieder der ursprünglichen Leibnizreihe zu, die sie mit Hilfe der sog. Euler-Transformation wie folgt umformten [3a],[5],[6]:

        1     2     3

π = 2 + -(2 + -(2 + -(2 + ...))) .                            (6)

        3     5     7         

Wie diese Umformung im einzelnen vor sich ging, soll hier nicht erklärt werden und ist für das Weitere auch nicht erforderlich. Die geschachtelte Klammerschreibweise in (6) bewirkt, daß nicht mehr nur addiert, sondern auch multipliziert wird. Dadurch erhöht sich die Konvergenzgeschwindigkeit beträchtlich, so daß es bei geeigneter Programmierung nicht schwerfällt, 1000 Stellen von π und mehr in kurzer Zeit zu berechnen.

Nun aber gibt es außer der Arkustangensreihe noch eine andere Reihe, die ebenfalls zur π-Berechnung herangezogen werden kann: die Arkussinusreihe:

                1  x^3    1·3 x^5   1·3·5 x^7

arcsin x = x + ---·--- + ---·--- + -----·--- + ....

                2   3     2·4  5    2·4·6  7 

Setzt man in ihr x=1/2, so ergibt sich wegen sin π/6=1/2 nach leichter Umformung:

        1    3·3       3·3·5       3·3·5·7 

π = 3 + - + ------ + --------- + ----------- + ...  .          (7)

        8   4·32·5   4·6·128·7   4·6·8·512·9

Dies läßt sich, ohne Benutzung der Euler-Transformation, durch fortgesetztes, einfaches Ausklammern, auf die Form

           1·1       3·3       5·5           

π  =  3 + -----(3 + -----(3 + -----(3 + ...)))                 (8)

          8·1·3     8·2·5     8·3·7 

bringen. Sie ähnelt der Darstellung (6) von RABINOWITZ und WAGON, weist aber zwei Besonderheiten auf: Zum einen ist Gl.(8) allem Anschein nach neu; in den ausgedehnten Formelsammlungen von [3] und [4] wird sie nicht erwähnt. Zum andern – und das ist für das weitere Vorgehen von größerer Bedeutung – konvergiert (8) doppelt so schnell wie (6) ¹⁾. 

Auf (8) beruht das nächste Programm:

Program pitau; {1000 Stellen von Pi}

 uses crt,dos;const n=1000;

 Var  i,j,k:integer;c,d,q,u,x:word;

      a    :array[1..n+1] of word;

procedure divi(y:word);

begin

 c:=0;for j:=1 to n+1 do

 begin x:=a[j]+c;q:=x div y;a[j]:=q;

 d:=x-y*q;c:=10*d;end;

end;

procedure mult(y:word);

begin

 for j:=1 to n+1 do a[j]:=y*a[j];

 for j:=n+1 downto 2 do

 begin u:=a[j] div 10;a[j-1]:=a[j-1]+u;

 a[j]:=a[j] mod 10;end;

end;

Begin

 clrscr;k:=trunc(n*ln(10)/ln(4));

 for i:=k downto 1 do begin divi(8);

 divi(i);mult(2*i-1);divi(2*i+1);

 mult(2*i-1);a[1]:=a[1]+3;end;write(' ');

 for i:=1 to n+1 do begin write(a[i]);

 if i=1 then write('.');if (i mod 6=0) then

 write(' ');if wherex=80 then write('   ');

 end;write('... (1000 Stellen)');repeat

 until keypressed;

End.

Es liefert bei 500 MHz Taktfrequenz 1000 Stellen von π in einer Sekunde. Das Programm ist insofern besonders einfach, als es nur zwei Prozeduren enthält, in denen das schriftliche Dividieren und Multiplizieren nachgeahmt wird. Eine besondere Additionsprozedur ist nicht erforderlich: es genügt, bei jedem Schritt den Inhalt der ersten Speicherplatz-Zelle des verwendeten Arrays um 3 zu erhöhen.

Für einen möglichst schnellen Ablauf wurden hauptsächlich Word-Variable verwendet. Dies hat zur Folge, daß sich das Programm nur auf knapp 2000 Stellen erweitern läßt.

Die verschachtelte Klammerschreibweise wie in (6) und (8) wirkt sich auch günstig bei der Arkustangensreihe aus, wenn π mit deren Hilfe berechnet werden soll. Hier gilt, was anscheinend auch nirgends erwähnt wird:

                         1       1         3         5  

arctan(1/m) = -(1 - -----(1 - -----(1 - -----(1 - ...)))), m>1.   (9)  

              m     3·m·m     5·m·m     7·m·m

Wird (9) auf die Formel (5) von MACHIN angewendet, kann man wie folgt programmieren:

 Program pimachin;{1000 Stellen, mit der Formel von Machin berechnet}

 uses crt,dos;{ pi/4 = 4 arctan(1/5) - arctan(1/239) }

  const n=1005;{5 Sicherheitsstellen}

  var i,j,k,m,nr,di:integer;

      c,d,q,u,x:word;

      a:array[1..2,1..n+1] of word;

      ta,te:real;h,mi,se,hs:word;

 procedure divi(y:word);

 begin

  c:=0;for j:=1 to n+1 do

  begin x:=a[nr,j]+c;q:=x div y;a[nr,j]:=q;

  d:=x-y*q;c:=10*d;end;

 end;

 procedure mult(y:word);

 begin

  for j:=1 to n+1 do a[nr,j]:=y*a[nr,j];

  for j:=n+1 downto 2 do

  begin u:=a[nr,j] div 10;a[nr,j-1]:=a[nr,j-1]+u;

  a[nr,j]:=a[nr,j] mod 10;end;

 end;

 procedure atn;

 begin

  a[nr,1]:=1;k:=trunc(n*ln(10)/ln(m)/2);

  for i:=k downto 1 do

  begin

   divi(2*i+1);mult(2*i-1);divi(m);divi(m);

   for j:=2 to n-1 do a[nr,j]:=9-a[nr,j];

   a[nr,n]:=10-a[nr,n];

  end;

  divi(m);

 end;

  Begin

   gettime(h,mi,se,hs);ta:=3600*h+60*mi+se+hs/100;

   clrscr;writeln;

   nr:=1;m:=5;atn;mult(4);

   nr:=2;m:=239;atn;

   for i:=n downto 2 do

   begin

    di:=a[1,i]-a[2,i];

    if di<0 then

    begin di:=di+10;a[1,i-1]:=a[1,i-1]-1;end;

    a[1,i]:=di;

   end;

   nr:=1;mult(4);

   writeln('    ¶=3.');write('    ');

   gettime(h,mi,se,hs);te:=3600*h+60*mi+se+hs/100;

   for i:=2 to n-4 do

   begin

    write(a[1,i]);

    if ((i-1) mod 10=0) then write(' ');

    if ((i-1) mod 50=0) then begin write('(',i-1,')');

    writeln;write('    ');end;

   end;

   write('Rechenzeit ',te-ta:3:2,' s');

   repeat until keypressed;

  End.

Das Programm ist nur wenig komplizierter als das erste, aber doppelt so schnell. Mit ihm erhält man bei derselben Taktfrequenz wie oben (die auch im folgenden weiter vorausgesetzt wird), 1000 Stellen in einer halben Sekunde. Es läßt sich bis auf 4600 Stellen erweitern.

Bereits 1994 veröffentlichte DAVID ADAMSON [7] ein Programm für 2150 Stellen, in dem er die Machin-Formel in folgender Form anwandte:

π = 3.2 + 1/25(-3.2/3 + 1/25(3.2/5 + 1/25(-3.2/7 + ...).).) 

    - 1/239(4 +1/239^2(-4/3 +1/239^2(4/5 +1/239^2(-4/7 +...).).).  

Sie ähnelt unserer geschachtelten Klammerschreibweise (9), ist aber aus leicht erkennbaren Gründen für das Programmieren weniger vorteilhaft als diese. Bei gleicher Schnelligkeit ist das Programm von Adamson recht lang und unübersichtlich.

10000 Stellen in einer Minute erhält man mit einer Formel von GAUSS (1777-1855):

 π=48*arctan(1/38)+80*arctan(1/57)+28*arctan(1/239)+96*arctan(1/268)

und dem Programm

 Program pigauss;{10000 Stellen von pi mit der Gauss-Formel}

 uses crt,dos;{pi=48*atn(1/38)+80*atn(1/57)+28*atn(1/239)+96*atn(1/268)}

 const n=10005;{5 Sicherheitsstellen}

 var i,j,k,m,nr,di:integer;

     c,d,q,u,x:word;

     a:array[1..2,1..n] of word;

     ch:char;

 procedure divi(y:word);

 begin

  c:=0;for j:=1 to n+1 do

  begin x:=a[nr,j]+c;q:=x div y;a[nr,j]:=q;

  d:=x-y*q;c:=10*d;end;

 end;

 procedure mult(y:word);

 begin

  for j:=1 to n+1 do a[nr,j]:=y*a[nr,j];

  for j:=n+1 downto 2 do

  begin u:=a[nr,j] div 10;a[nr,j-1]:=a[nr,j-1]+u;

  a[nr,j]:=a[nr,j] mod 10;end;

 end;

 procedure atn;

 begin

  a[nr,1]:=1;k:=trunc(n*ln(10)/ln(m)/2);

  for i:=k downto 1 do

  begin

   divi(2*i+1);mult(2*i-1);divi(m);divi(m);

   for j:=2 to n-1 do a[nr,j]:=9-a[nr,j];

   a[nr,n]:=10-a[nr,n];

  end;

  divi(m);

 end;

 procedure addi;

 begin

  for i:=n downto 2 do

   begin

    a[1,i]:=a[1,i]+a[2,i];u:=a[1,i]div 10;

    a[1,i-1]:=a[1,i-1]+u;

   end;

  for i:=1 to n do a[1,i]:=a[1,i] mod 10;

  end;

  procedure warten;

  begin

   writeln;writeln;

   write(i-1,' Stellen. ');

   if (i<n-4) then write('Weiter mit der Leertaste.')

   else write('Zurueck zum Programm mit der Leertaste.');

   repeat ch:=readkey until ch=' ';

   clrscr;writeln;

  end;

  Begin

   clrscr;writeln('Berechnet werden 10000 Stellen von Pi.');

   writeln('Bitte (bei 500 MHz Taktfrequenz) eine Minute warten!');

   nr:=1;m:=38;atn;mult(12);nr:=2;m:=57;atn;mult(20);addi;

   m:=239;atn;mult(7);addi;m:=268;atn;mult(24);addi;

   nr:=1;mult(4);clrscr;writeln('3.');

   for i:=2 to n-4 do

   begin

    write(a[1,i]);if i=1 then write('.');

    if (i-1) mod 4=0 then write(' ');

    if (i-1) mod 1000=0 then warten;

   end;

  End.

Angezeigt werden sie zu je 1000 beim Drücken der Leertaste. Ihre Richtigkeit überprüfen kann man bis zu 5000 mit [3], darüber hinaus mit [8].

 --------------------------
¹⁾ In (6) gehen die Faktoren k/(2k+1) vor den Klammern, k=1,2,3,..., für k→∞ gegen ½. Dies hat zur Folge, wie in [6] angemerkt wird, dass k=n/lg2 Faktoren berechnet werden müssen, um n Nachkommastellen von p zu erhalten.
In (8) gehen die Faktoren (2k-1)²/(8k(2k+1)) für k→∞ gegen ¼. Entsprechend sind von ihnen nur n/lg4 nötig, d. h. halb so viele wie bei (6). (Gemeint sind die ganzzahligen Anteile von n/lg2 und n/lg4, in den Programmen mit "trunc" bezeichnet.)
Ergänzung: mit dem Taschenrechner erhält man in zwölf Schritten mit (6) π≈3,141479649, d.h. drei richtige Nachkommastellen, und mit (8) π≈3,141592653; das sind acht. (Wie man dabei am besten vorgeht - vom Ende aus rückwärts -, wird hier erklärt.)

Literatur