Easy Proof of Three Recursive $\pi$-Algorithms -- Einfacher Beweis dreier rekursiver $\pi$-Algorithmen

TL;DR

This paper provides elementary algebra and calculus proofs demonstrating that Borwein's quartic algorithm and the Brent-Salamin algorithm both converge to pi, with the former being more iteration-efficient and having quartic convergence.

Contribution

It offers new elementary algebra and calculus proofs establishing the equivalence, convergence, and convergence rates of Borwein's and Brent-Salamin pi algorithms.

Findings

Borwein's quartic algorithm matches Brent-Salamin's output with fewer iterations.

Brent-Salamin algorithm converges to pi using integral calculus.

Borwein's algorithm exhibits quartic convergence, faster than quadratic.

Abstract

This paper consists of three independent parts: First we use only elementary algebra to prove that the quartic algorithm of the Borwein brothers has exactly the same output as the Brent-Salamin algorithm, but that the latter needs twice as many iterations. Second we use integral calculus to prove that the Brent-Salamin algorithm approximates $\pi$. Combining these results proves that the Borwein brothers' quartic algorithm also approximates $\pi$. Third, we prove the quadratic convergence of the Brent-Salamin algorithm, which also proves the quartic convergence of Borwein's algorithm. -- -- Dieses Paper besteht aus drei unabh\"angigen Teilen: Erstens beweisen wir mit elementarer Algebra, dass der Borwein-Algorithmus vierter Ordnung die gleichen Ergebnisse liefert wie der Brent-Salamin-Algorithmus, wobei letzterer doppelt so viele Iterationen ben\"otigt. Zweitens beweisen wir mit Integralrechnung, dass der Brent-Salamin-Algorithmus gegen $\pi$ konvergiert. Hieraus folgt, dass der Borwein-Algorithmus vierter Ordnung ebenfalls gegen $\pi$ konvergiert. Drittens beweisen wir die quadratische Konvergenz des Brent-Salamin-Algorithmus und somit auch die quartische Konvergenz des Borwein-Algorithmus.

Equations

$$\left\{\begin{array}{rl}{a}_0& :=1\\ b_0& :=\frac{1}{\sqrt{2}}\end{array}\right\}\text{leavevmode}\text{leavevmode}\text{and}\text{leavevmode}\text{leavevmode}\left\{\begin{array}{rlrl}{a}_n& :=\frac{a_{n-1}+b_{n-1}}{2}& & \text{(arithmetic mean)}\\ b_n& :=\sqrt{a_{n-1}\cdot b_{n-1}}& & \text{(geometic mean)}\\ c_n^2& :=a_n^2-b_n^2& & \end{array}\right\}$$

$$⟹\text{leavevmode}\text{output after N iterations:}\text{leavevmode}\text{leavevmode}\text{leavevmode}{p}_N:=\frac{(a_N+b_N{)}^2}{1-2\cdot {\sum }_{j=1}^N{2}^j\cdot c_j^2}$$

$$\left\{\begin{array}{rl}{k}_0& :=3-2\cdot \sqrt{2}\\ e_0& :=6-4\cdot \sqrt{2}\end{array}\right\}\text{leavevmode}\text{leavevmode}\text{and}\text{leavevmode}\text{leavevmode}\left\{\begin{array}{rl}{k}_n& :=\frac{1-\sqrt{1-k_{n-1}^2}}{1+\sqrt{1-k_{n-1}^2}}\\ e_n& :=e_{n-1}\cdot (1+k_n{)}^2-2^{n+1}\cdot k_n\end{array}\right\}$$

$$⟹\text{leavevmode}\text{output after N iterations:}\text{leavevmode}\text{leavevmode}\text{leavevmode}{\hat{\pi }}_N:=\frac{1}{e_N}$$

$$\left\{\begin{array}{rl}{y}_0& :=\sqrt{2}-1\\ z_0& :=6-4\cdot \sqrt{2}\end{array}\right\}\text{leavevmode}\text{leavevmode}\text{and}\text{leavevmode}\text{leavevmode}\left\{\begin{array}{rl}{y}_n& :=\frac{1-\sqrt[4]{1-y_{n-1}^4}}{1+\sqrt[4]{1-y_{n-1}^4}}\\ z_n& :=z_{n-1}\cdot (1+y_n{)}^4-2\cdot 4^n\cdot y_n\cdot (1+y_n+y_n^2)\end{array}\right\}$$

$$⟹\text{leavevmode}\text{output after N iterations:}\text{leavevmode}\text{leavevmode}\text{leavevmode}{\pi }_N:=\frac{1}{z_N}$$

$${\hat{\pi }}_N=p_N\text{and}{\pi }_N={\hat{\pi }}_{2N}=p_{2N}$$

$${\hat{\pi }}_N=p_N$$

$$K_0$$

$$=2\cdot (2-\sqrt{2})-1=3-2\cdot \sqrt{2}=k_0$$

$$E_0$$

$$=\frac{2}{3+2\cdot \sqrt{2}}=\frac{2\cdot (3-2\cdot \sqrt{2})}{9-4\cdot 2}=2\cdot (3-2\cdot \sqrt{2})=e_0$$

$$K_n$$

$$⟹K_n^2$$

$$⟹\sqrt{1-K_{n-1}^2}$$

$$K_n=\frac{a_n-b_n}{a_n+b_n}=\frac{1-b_n/a_n}{1+b_n/a_n}=\frac{1-\sqrt{1-K_{n-1}^2}}{1+\sqrt{1-K_{n-1}^2}}=\frac{1-\sqrt{1-k_{n-1}^2}}{1+\sqrt{1-k_{n-1}^2}}=k_n$$

$$E_n:=\frac{1}{p_n}=\frac{1-2\cdot {\sum }_{j=1}^n{2}^j{c}_j^2}{4\cdot a_{n+1}^2}$$

$$a_{n+1}^2\cdot E_n-a_n^2\cdot E_{n-1}$$

$$⟹E_n$$

$$E_n$$

$$E_n={\left(K_n+1\right)}^2\cdot E_{n-1}-2^{n+1}\cdot K_n$$

$$⟹E_n={\left(k_n+1\right)}^2\cdot e_{n-1}-2^{n+1}\cdot k_n$$

$${\pi }_N={\hat{\pi }}_{2N}$$

$$k_{n-1}=\sqrt{1-{\left(\frac{1-k_n}{1+k_n}\right)}^2}=\sqrt{\frac{(1+k_n{)}^2-(1-k_n{)}^2}{(1+k_n{)}^2}}=\frac{2\cdot \sqrt{k_n}}{1+k_n}$$

$$\frac{2\cdot Y_n}{1+Y_n^2}$$

$$⟹\sqrt{1-Y_{n-1}^4}$$

$$⟹\sqrt[4]{1-Y_{n-1}^4}$$

$$⟹Y_n$$

$$e_{2n}$$

$$\text{and}{e}_{2n-1}$$

$$e_{2n}$$

$$=e_{2n-2}\cdot \left[(1+k_{2n-1}{)}^2\cdot (1+k_{2n}{)}^2\right]-2^{2n}\cdot \left[k_{2n-1}\cdot (1+k_{2n}{)}^2+2\cdot k_{2n}\right]$$

$$e_{2n}$$

$$=e_{2n-2}\cdot \left[{\left(1+y_n^2+2\cdot y_n\right)}^2\right]-2^{2n}\cdot \left[2\cdot y_n\cdot (1+y_n^2)+2\cdot y_n^2\right]$$

$$=e_{2n-2}\cdot {\left(1+y_n\right)}^4-2^{2n+1}\cdot y_n\cdot \left(1+y_n+y_n^2\right)$$

$$Z_n:=e_{2n}$$

$$=z_{n-1}\cdot {\left(1+y_n\right)}^4-2^{2n+1}\cdot y_n\cdot \left(1+y_n+y_n^2\right)=z_n$$

$${\pi }_0={\hat{\pi }}_0=p_0$$

$${\hat{\pi }}_1=p_1$$

Full text

\setlength

-12.08mm\setlength-17.08mm

defi]Theorem

Easy Proof of Three Recursive $\pi$-Algorithms

  

Einfacher Beweis dreier rekursiver $\pi$-Algorithmen

Lorenz Milla, August 2019

Abstract. This paper consists of three independent parts:

First we use only elementary algebra to prove that the quartic algorithm of the Borwein brothers has exactly the same output as the Brent-Salamin algorithm, but that the latter needs twice as many iterations.

Second we use integral calculus to prove that the Brent-Salamin algorithm approximates $\pi$. Combining these results proves that the Borwein brothers’ quartic algorithm also approximates $\pi$.

Third, we prove the quadratic convergence of the Brent-Salamin algorithm, which also proves the quartic convergence of Borwein’s algorithm.

English version: pp. Easy Proof of Three Recursive $\pi$-Algorithms    Einfacher Beweis dreier rekursiver $\pi$-Algorithmen–3.2

Zusammenfassung. Dieses Paper besteht aus drei unabhängigen Teilen:

Erstens beweisen wir mit elementarer Algebra, dass der Borwein-Algorithmus vierter Ordnung die gleichen Ergebnisse liefert wie der Brent-Salamin-Algorithmus, wobei letzterer doppelt so viele Iterationen benötigt.

Zweitens beweisen wir mit Integralrechnung, dass der Brent-Salamin-Algo-rithmus gegen $\pi$ konvergiert. Hieraus folgt, dass der Borwein-Algorithmus vierter Ordnung ebenfalls gegen $\pi$ konvergiert.

Drittens beweisen wir die quadratische Konvergenz des Brent-Salamin-Algo-rithmus und somit auch die quartische Konvergenz des Borwein-Algorithmus.

Deutsche Version: S. Easy Proof of Three Recursive $\pi$-Algorithms    Einfacher Beweis dreier rekursiver $\pi$-Algorithmen–3.2

Introduction: The Algorithms

This paper is about the following three recursive $\pi$-algorithms:

We prove that these three algorithms produce the same approximations of $\pi$, where the number of correct digits is being doubled or quadrupled with each iteration.

These results have been proven before, but we elaborate all intermediate calculations and we use only elementary algebra and integral calculus.

Our proof consists of three independent chapters:

1**.   Proof of Equivalence of the Algorithms    1**

We prove that the algorithms have the same outputs: $\widehat{\pi}_{N}=p_{N}$ and $\pi_{N}=p_{2N}$.

2**.   Proof of the Brent-Salamin Algorithm    2**

We prove that the output $p_{N}$ of Brent-Salamin converges to $\pi$ as $N\rightarrow\infty$.

3**.   Proof of Brent-Salamin’s Quadratic Convergence    3**

We prove the quadratic convergence of $p_{N}$:  $|\pi-p_{n+1}|<0.075\cdot|\pi-p_{n}|^{2}$.Here $\pi$ denotes the limit of $p_{N}$.

1. Proof of Equivalence of the Algorithms

We call two algorithms “equivalent” if they produce the same outputs. In this chapter we will prove that the three algorithms on p. 2 are equivalent. More precisely:

{theo}

For the outputs of the three algorithms on p. 1, where

• •

$p_{N}$ is the output of the Brent-Salamin Alg. 1,

• •

$\widehat{\pi}_{N}$ is the output of the Borweins’ quadratic Alg. 2,

• •

$\pi_{N}$ is the output of the Borweins’ fourth order Alg. 3,

it holds:

[TABLE]

thus these algorithms produce the same sequence of outputs if the outputs are calculated exactly.

Proof.

This has been proven by Brent [1] who used elliptic modular functions and by Guillera [2] who used a theorem of Gauss, but we will need only elementary algebra for the proof of $\widehat{\pi}_{N}=p_{N}$ in Prop. 1.2 and for the proof of $\pi_{N}=\widehat{\pi}_{2N}$ in Prop. 1.3. ∎

Remark 1.1**.**

When the computations are done using floating-point or interval arithmetic, the initial values and the iterations can only be done with finite precision. This produces rounding errors which propagate differently in the three algorithms. The outputs of the algorithms thus differ in the last decimals. To compute $D$ decimals of $\pi$ correctly, one has to compute all initial and intermediate values to some extended precision (e.g. to $D+50$ decimals), and the additional decimals have to be cut off in the end.

Proposition 1.2**.**

For the sequences defined in the Brent-Salamin Alg. 1 and the Borweins’ Alg. 2 on p. 1 it holds $e_{n}=1/p_{n}$ and $k_{n}=a_{n}/a_{n+1}-1$. In particular it holds

[TABLE]

thus these two algorithms produce the same sequence of outputs.

Proof.

We set $E_{n}:=1/p_{n}$ and $K_{n}:=a_{n}/a_{n+1}-1$ and prove by induction that it holds $E_{n}=e_{n}$ and $K_{n}=k_{n}$:

• •

First we prove $K_{0}=k_{0}$ and $E_{0}=e_{0}$:

[TABLE]

• •

Now we prove $K_{n}=k_{n}$ and $E_{n}=e_{n}$ using the induction hypothesis (which states $K_{n-1}=k_{n-1}$ and $E_{n-1}=e_{n-1}$):

[TABLE]

Using the induction hypothesis this yields $K_{n}=k_{n}$:

[TABLE]

From $p_{n}:=\frac{(a_{n}+b_{n})^{2}}{1-2\cdot\sum_{j=1}^{n}2^{j}c_{j}^{2}}=\frac{4\cdot a_{n+1}^{2}}{1-2\cdot\sum_{j=1}^{n}2^{j}c_{j}^{2}}$ we obtain

[TABLE]

This yields

[TABLE]

Using $\frac{c_{n}^{2}}{a_{n+1}^{2}}=\frac{a_{n}^{2}-b_{n}^{2}}{(a_{n}^{2}+b_{n})^{2}/4}=4\cdot\frac{a_{n}-b_{n}}{a_{n}+b_{n}}=4\cdot\left(\frac{2a_{n}}{a_{n}+b_{n}}-\frac{a_{n}+b_{n}}{a_{n}+b_{n}}\right)=4\cdot\left(\frac{a_{n}}{a_{n+1}}-1\right)$ we get:

[TABLE]

Here we replace $a_{n}/a_{n+1}$ by $K_{n}+1$ and obtain

[TABLE]

But we already proved $K_{n}=k_{n}$. Thus the induction hypothesis $E_{n-1}=e_{n-1}$ implies:

[TABLE]

Here we recognize the definition of $e_{n}$, thus we have proven $E_{n}=e_{n}$.

This proves $\widehat{\pi}_{N}=1/e_{N}=p_{N}$ for all $N\in\mathbb{N}$, thus the two algorithms produce the same sequence of outputs. ∎

Proposition 1.3**.**

For the sequences defined in the Borweins’ Alg. 2 and 3 on p. 2 it holds $y_{n}=\sqrt{k_{2n}}$ and $z_{n}=e_{2n}$. In particular it holds

[TABLE]

thus one iteration of Alg. 3 is equivalent to two iterations of Alg. 2.

Proof.

We set $Y_{n}:=\sqrt{k_{2n}}$ and $Z_{n}:=e_{2n}$ and prove by induction that it holds $Y_{n}=y_{n}$ and $Z_{n}=z_{n}$:

• •

First we observe that $Z_{0}:=e_{0}=6-4\cdot\sqrt{2}=z_{0}$. Then it holds $Y_{0}:=\sqrt{k_{0}}=y_{0}$, because $y_{0}^{2}=\left(\sqrt{2}-1\right)^{2}=3-2\cdot\sqrt{2}=k_{0}$.

• •

Now we prove $Y_{n}=y_{n}$ and $Z_{n}=z_{n}$ using the induction hypothesis (which states $Y_{n-1}=y_{n-1}$ and $Z_{n-1}=z_{n-1}$):

From $k_{n}:=\frac{1-\sqrt{1-k_{n-1}^{2}}}{1+\sqrt{1-k_{n-1}^{2}}}$ we get $\sqrt{1-k_{n-1}^{2}}=\frac{1-k_{n}}{1+k_{n}}$ and thus

[TABLE]

Next, $Y_{n}:=\sqrt{k_{2n}}$ yields $k_{2n}=Y_{n}^{2}$ and $k_{2n-1}=\frac{2\cdot\sqrt{k_{2n}}}{1+k_{2n}}=\frac{2\cdot Y_{n}}{1+Y_{n}^{2}}$. This implies:

[TABLE]

Thus we have proven $Y_{n}=y_{n}$ using the induction hypothesis $Y_{n-1}=y_{n-1}$ in the last step. It remains to prove $Z_{n}=z_{n}$:

From the definition of $e_{n}$ in Alg. 2 we obtain

[TABLE]

Putting this representation of $e_{2n-1}$ into the one of $e_{2n}$ yields:

[TABLE]

Using $k_{2n}=Y_{n}^{2}=y_{n}^{2}$ and $k_{2n-1}=\frac{2\cdot Y_{n}}{1+Y_{n}^{2}}=\frac{2\cdot y_{n}}{1+y_{n}^{2}}$ we obtain:

[TABLE]

Here we use the induction hypothesis $z_{n-1}=Z_{n-1}=e_{2n-2}$:

[TABLE]

Thus we have proven that $\pi_{N}=1/z_{N}=1/e_{2N}=\widehat{\pi}_{2N}$ holds for all $N\in\mathbb{N}$, thus that Alg. 3 produces every second output of Alg. 2. ∎

Proof of Thm. 1.

In Prop. 1.2 we proved $\widehat{\pi}_{N}=p_{N}$ and in Prop. 1.3 we proved $\pi_{N}=\widehat{\pi}_{2N}$ – thus both statements from Thm. 1 are proven, and the algorithms are equivalent. ∎

Remark 1.4**.**

The first outputs of the three equivalent algorithms are:

[TABLE]

2. Proof of the Brent-Salamin Algorithm

In this chapter we prove that the Brent Salamin algorithm converges to $\pi$. This proof elaborates [3] and uses only integral calculus like integration by parts or by substitution (also: two-dimensional substitution).

{theo}

It holds the following formula due to Gauß (1809), Brent (1976) and Salamin (1976):

[TABLE]

Here, $\operatorname{AGM}(1;1/\sqrt{2})$ denotes the arithmetic-geometric mean (i.e. the common limit of $a_{n}$ and $b_{n}$ from the Brent-Salamin algorithm on p. 1). In particular, the sequence

[TABLE]

of the Brent-Salamin algorithm on p. 1 converges to $\pi$.

Proof.

First we generalize initial values of the Brent-Salamin algorithm to

[TABLE]

Later (from Prop. 2.9 onwards) we will use $a=1$ and $b=1/\sqrt{2}$. On p. 2 we will continue the proof of Thm. 2, but first we proof some auxiliary propositions: ∎

Proposition 2.1**.**

The geometric mean $\sqrt{x\cdot y}$ and the arithmetic mean $\frac{x+y}{2}$ of two positive real numbers $x\neq y$ satisfy:

[TABLE]

Proof.

From $x\neq y$ we deduce:

[TABLE]

This yields $4xy<(x+y)^{2}$ and proves that the geometric mean $\sqrt{x\cdot y}$ is less than the arithmetic mean $\frac{x+y}{2}$. ∎

Proposition 2.2**.**

The sequences $a_{n}$ and $b_{n}$ of the Brent-Salamin algorithm 1 converge to a common limit which we call $\operatorname{AGM}(a,b)$. The convergence of $a_{n}\searrow\operatorname{AGM}(a,b)$ and of $b_{n}\nearrow\operatorname{AGM}(a,b)$ is strictly monotonic and it holds $c_{n+1}^{2}<\frac{1}{4}c_{n}^{2}$.

Proof.

Prop. 2.1 tells that $b_{n}<a_{n}$ holds for all $n$. This implies the strict monotonicity of $b_{n+1}=\sqrt{a_{n}\cdot b_{n}}>\sqrt{b_{n}\cdot b_{n}}=b_{n}$ and $a_{n+1}=\frac{a_{n}+b_{n}}{2}<\frac{a_{n}+a_{n}}{2}=a_{n}$. Both sequences are bounded by $b=b_{0}\leq b_{n}<a_{n}\leq a_{0}=a$ and thus convergent. For $c_{n+1}^{2}$ it holds:

[TABLE]

This proves that $c_{n}^{2}=a_{n}^{2}-b_{n}^{2}<4^{-n}\cdot c_{0}^{2}$ converges to zero and that $a_{n}$ and $b_{n}$ have the same limit. ∎

Proposition 2.3**.**

The value of

[TABLE]

is constant on the whole AGM sequence, i.e. it holds $I(a_{n},b_{n})=I(a_{0},b_{0})$ for all $n\in\mathbb{N}$.

Proof.

First we substitute $t=b\cdot\tan\Phi$. Then $1+\tan^{2}\Phi=\frac{1}{\cos^{2}\Phi}$ yields $\cos^{2}\Phi=\frac{b^{2}}{b^{2}+b^{2}\tan^{2}\Phi}=\frac{b^{2}}{b^{2}+t^{2}}$ and $\sin^{2}\Phi=1-\cos^{2}\Phi=\frac{b^{2}+t^{2}-b^{2}}{b^{2}+t^{2}}=\frac{t^{2}}{b^{2}+t^{2}}$. Further, it holds $\frac{\mathrm{d}t}{\mathrm{d}\Phi}=b\cdot(1+\tan^{2}\Phi)=b+\frac{t^{2}}{b}=\frac{t^{2}+b^{2}}{b}$, thus $\frac{\mathrm{d}\Phi}{\mathrm{d}t}=\frac{b}{t^{2}+b^{2}}$. This shows that the substitution yields the following representation of $I(a,b)$:

[TABLE]

Now we substitute $x=\frac{1}{2}\left(t-\frac{ab}{t}\right)$. This yields $2xt=t^{2}-ab$ and t=x\mathbin{\vbox{\hbox{\oalign{\hfil\scriptstyle+$\hfil\cr\kern-1.29167pt\cr$\scriptscriptstyle({-})\cr}}}}\sqrt{x^{2}+ab} (since $t>0$). Thus it holds $\frac{\mathrm{d}t}{\mathrm{d}x}=1+\frac{2x}{2\sqrt{x^{2}+ab}}=\frac{\sqrt{x^{2}+ab}+x}{\sqrt{x^{2}+ab}}=\frac{t}{\sqrt{x^{2}+ab}}$ and:

[TABLE]

Here we have denoted $\frac{(t^{2}+a^{2})\cdot(t^{2}+b^{2})}{t^{2}}$ by $f(x)$ (remember $t>0$). About $f(x)$ it holds:

[TABLE]

This yields:

[TABLE]

Here we use the fact that the integrand is even, thus $\frac{1}{2}\int\limits_{-\infty}^{\infty}$ yields $\int\limits_{0}^{\infty}$:

[TABLE]

Now we have proven that for any $a>b>0$ it holds $I\mathopen{}\left(\frac{a+b}{2},\sqrt{ab}\right)\mathclose{}=I(a,b)$. By induction this yields $I(a_{n},b_{n})=I(a_{0},b_{0})$ for all $n\in\mathbb{N}$. ∎

Proposition 2.4**.**

Let $I(a,b)$ be the integral from Prop. 2.3. Then it holds:

[TABLE]

Proof.

With $m:=\operatorname{AGM}(a,b)$, Prop. 2.2 tells that $a_{n}$ and $b_{n}$ converge to $m$. If we interchange the limit and the integration, Prop. 2.3 yields:

[TABLE]

∎

Proposition 2.5**.**

If we denote

[TABLE]

then it holds $L(b,a)+L(a,b)=I(a,b)$ and $L(b,a)-L(a,b)=\frac{a-b}{a+b}\cdot L(b_{1},a_{1})$.

Proof.

To prove the first equation, we substitute $\Phi^{\prime}=\frac{\pi}{2}-\Phi$. Then it holds $\cos(\Phi^{\prime})=\sin(\Phi)$ and $\sin(\Phi^{\prime})=\cos(\Phi)$, thus

[TABLE]

With $\sin^{2}+\cos^{2}=1$ we deduce the first equation:

[TABLE]

Next we prove an alternative representation of $L(a,b)$, similar to the one of $I(a,b)$ in eq. (2.2): again we substitute $t=b\cdot\tan\Phi$ and obtain

[TABLE]

Then we calculate $L(b,a)$ by interchanging $a$ and $b$:

[TABLE]

From $\frac{a^{2}}{a^{2}+t^{2}}-\frac{b^{2}}{b^{2}+t^{2}}=\frac{a^{2}(b^{2}+t^{2})-b^{2}(a^{2}+t^{2})}{(a^{2}+t^{2})(b^{2}+t^{2})}=\frac{a^{2}t^{2}-b^{2}t^{2}}{(a^{2}+t^{2})(b^{2}+t^{2})}=\frac{(a^{2}-b^{2})t^{2}}{(a^{2}+t^{2})(b^{2}+t^{2})}$ we deduce:

[TABLE]

And, as with $I(a,b)$, we substitute $x=\frac{1}{2}\left(t-\frac{ab}{t}\right)$:

[TABLE]

where again it holds $f(x):=\frac{(t^{2}+a^{2})\cdot(t^{2}+b^{2})}{t^{2}}=(2x)^{2}+(a+b)^{2}$ and thus:

[TABLE]

In this integral we recognize the representation (2.3) of $L(b_{1},a_{1})$. This yields:

[TABLE]

Thus we have also proven the second equation. ∎

Proposition 2.6**.**

Denoting $S:=\sum_{j=1}^{\infty}2^{j}\cdot c_{j}^{2}$ it holds:

[TABLE]

Proof.

From $4\cdot(a_{1}^{2}-b_{1}^{2})=4\cdot\left(\frac{a+b}{2}\right)^{2}-4ab=(a+b)^{2}-4ab=(a-b)^{2}$ we deduce (using both equations of Prop. 2.5):

[TABLE]

With the definition of $c_{n}$ this reads:

[TABLE]

Thus it holds for all $j\in\mathbb{N}$:

[TABLE]

Here we multiply with $2^{j}$ and use $I(a_{j},b_{j})=I(a,b)$ from Prop. 2.3:

[TABLE]

Now we add these equations for $0\leq j\leq n$ and obtain:

[TABLE]

In the second sum we shift the index $k=j+1$:

[TABLE]

This shows that the left side of (2.4) is a telescoping series, thus nearly all terms cancel each other out:

[TABLE]

Here we estimate $L(b_{n+1},a_{n+1})<I(b_{n+1},a_{n+1})=I(b,a)$ and use $c_{n+1}^{2}<4^{-n-1}c_{0}^{2}$ from Prop. 2.2:

[TABLE]

Thus the second term of (2.5) tends to [math] (for $n\rightarrow\infty$) and we obtain:

[TABLE]

Finally we use $L(b,a)=I(a,b)-L(a,b)$ from Prop. 2.5 and obtain:

[TABLE]

∎

Proposition 2.7**.**

The Gamma function $\Gamma(x):=\int_{0}^{\infty}t^{x-1}\cdot e^{-t}\mathrm{d}t$ satisfies for $\operatorname{Re}(x)>0$:

[TABLE]

Proof.

We prove the functional equation integrating by parts:

[TABLE]

To calculate $\Gamma(1/2)$ we substitute $s=\sqrt{t}$ and obtain $t=s^{2}$ and $\frac{\mathrm{d}t}{\mathrm{d}s}=2s$:

[TABLE]

We square this integral to obtain a twodimensional integral:

[TABLE]

Now we use polar coordinates $x=r\cdot\cos\Phi$ and $y=r\cdot\sin\Phi$:

[TABLE]

Thus the value of $\Gamma(1/2)=\sqrt{\pi}$ is also proven. ∎

Proposition 2.8**.**

The Beta function $B(u,v):=\int_{0}^{1}t^{u-1}\cdot(1-t)^{v-1}\mathrm{d}t$ satisfies for $\operatorname{Re}(u)>0$ and $\operatorname{Re}(v)>0$:

[TABLE]

Proof.

We start with $\Gamma(u)\cdot\Gamma(v)$:

[TABLE]

The substitution of $\left\{\begin{aligned} s&=x\cdot(1-y)\\ t&=x\cdot y\end{aligned}\right\}$ or $\left\{\begin{aligned} x&=s+t\\ y&=\frac{t}{x}=\frac{t}{s+t}\end{aligned}\right\}$ yields $0<x<\infty$ and $0<y<1$. The Jacobian matrix of this substitution is

[TABLE]

This yields:

[TABLE]

A final division by $\Gamma(u+v)$ proves the equation of Prop. 2.8. ∎

Proposition 2.9**.**

It holds

[TABLE]

Proof.

First, $\sin^{2}\Phi=1-\cos^{2}\Phi$ yields:

[TABLE]

Then we substitute $x=\cos\Phi$ with $\frac{\mathrm{d}x}{\mathrm{d}\Phi}=-\sin\Phi=-\sqrt{1-\cos^{2}\Phi}=-\sqrt{1-x^{2}}$ and thus $\frac{\mathrm{d}\Phi}{\mathrm{d}x}=\frac{-1}{\sqrt{1-x^{2}}}$:

[TABLE]

Next we substitute $x=t^{1/4}$ with $\frac{\mathrm{d}x}{\mathrm{d}t}=\frac{1}{4}\cdot t^{-3/4}$ to obtain the Beta function $B(x,y)$ from Prop. 2.8:

[TABLE]

Now we replace the Beta function by Gamma functions as in Prop. 2.8 and use the properties of the Gamma function from Prop. 2.7:

[TABLE]

∎

Proof of Thm. 2.

The integrals $I(a,b)$ and $L(a,b)$ satisfy:

[TABLE]

This implies $I(1,1/\sqrt{2})=\sqrt{2}\cdot I(\sqrt{2},1)$ and $L(1,1/\sqrt{2})=\sqrt{2}\cdot L(\sqrt{2},1)$ and then with Prop. 2.9:

[TABLE]

Prop. 2.6 yields

[TABLE]

Here we use $a=1$ and $b=1/\sqrt{2}$, thus $c_{0}^{2}=a^{2}-b^{2}=1/2$:

[TABLE]

And now Prop. 2.4 proves:

[TABLE]

Solving this for $\pi$ proves the Gaussian formula:

[TABLE]

The numerator of $p_{N}$ is $(a_{N}+b_{N})^{2}=4\cdot a_{N+1}^{2}$ which converges to $4\cdot\operatorname{AGM}(1;1/\sqrt{2})^{2}$ (cf. Prop. 2.2), thus the sequence $p_{N}$ converges to $\pi$. ∎

Remark 2.10**.**

Now we have proven that the Brent-Salamin algorithm approximates $\pi$. We have already proven the equivalence of the three algorithms, thus the two algorithms by the Borwein brothers also approximate $\pi$.

3. Proof of Brent-Salamin’s Quadratic Convergence

This chapter does not use the fact that the Brent-Salamin sequence $p_{n}$ converges to $\pi\approx 3.14159$. We only use the monotonic convergence of $a_{n}$ and $b_{n}$ from Prop. 2.2 and use the symbol $\pi$ as a placeholder for the limit of this sequence $p_{n}$.

{theo} The sequence $\displaystyle p_{n}:=\frac{(a_{n}+b_{n})^{2}}{1-2\cdot\sum_{j=1}^{n}2^{j}c_{j}^{2}}$ of the Brent-Salamin algorithm converges quadratically to its limit which we denote with the symbol $\pi$:

[TABLE]

In particular, the number of significant digits is approximately doubled with each iteration, but one has to do all calculations with the desired accuracy.

Proof.

First we denote the (de)nominator of the Brent-Salamin sequence with $X_{n}$ and $Y_{n}$. In Prop. 2.2 we proved that $a_{n}\rightarrow\operatorname{AGM}(1,1/\sqrt{2})$, thus it holds:

[TABLE]

Now we denote the differences with

[TABLE]

In Prop. 2.2 we proved that $a_{n+1}>\operatorname{AGM}(a,b)>b_{n+1}$ holds for all $n$. This yields:

[TABLE]

Since the summation of $\delta_{n}$ contains only positive terms, it holds $\delta_{n}>2\cdot 2^{n+1}\cdot c_{n+1}^{2}$. Finally we have proven in Prop. 2.2 that $c_{n+1}^{2}<c_{n}^{2}/4$, thus we can use the geometric series to estimate $\delta_{n}$ (by setting the summation index to $j=n+1+k$):

[TABLE]

Thus we have proven:

[TABLE]

The difference between $p_{n}$ and its limit $\pi$ is:

[TABLE]

Here we use $1>Y_{n}>Y>0$ and $\delta_{n}>\varepsilon_{n}>0$ to obtain:

[TABLE]

Now, eq. (3.1) yields:

[TABLE]

In the proof of Prop. 2.2 we proved $c_{n+1}^{2}=(a_{n}-b_{n})^{2}/4$. Thus it holds:

[TABLE]

Here we see the quadratic convergence:

[TABLE]

Using first (3.2) and then (3.3) yields:

[TABLE]

Here we use $\operatorname{AGM}>b_{1}$ (Prop. 2.2) and $\operatorname{AGM}(1,1/\sqrt{2})>\sqrt{1\cdot 1/\sqrt{2}}=\sqrt[4]{1/2}$ and $2^{n}\geq 1$ to obtain:

[TABLE]

This proves the quadratic convergence of $p_{n}$. ∎

Remark 3.1**.**

Using $p_{1}>3.14057$, Thm. 3 yields $|\pi-p_{44}|<10^{-3.6\cdot 10^{13}}$. Thus $p_{44}=\widehat{\pi}_{44}=\pi_{22}$ is closer to $\pi$ than the current world record of $\pi\cdot 10^{13}$ digits (March 2019).

Remark 3.2**.**

For $p_{44}$, it even holds $|\pi-p_{44}|<10^{-4.8\cdot 10^{13}}$. This follows from the stronger error bound $0<\pi-p_{n}<\left(2^{n+4}\pi^{2}-8\pi\right)\cdot\exp\mathopen{}\left(-2^{n+1}\pi\right)\mathclose{}$ which is proven in [1, eq. (20)].

Einleitung: Die Algorithmen

Dieses Paper handelt von den folgenden drei rekursiven $\pi$-Algorithmen:

Wir beweisen, dass diese drei Algorithmen die gleichen Näherungen der Zahl $\pi$ berechnen, wobei die Anzahl korrekter $\pi$-Dezimalen mit jeder Iteration ungefähr verdoppelt bzw. vervierfacht wird.

Diese Resultate wurden bereits anderswo bewiesen, aber wir führen alle Rechnungen explizit aus und wir verwenden nur elementare Algebra und Integralrechnung.

Unser Beweis besteht aus drei unabhängigen Kapiteln:

1**.   Beweis der Äquivalenz der Algorithmen    1**

Wir beweisen, dass die drei Algorithmen die gleiche Ausgabe liefern, also dass $\widehat{\pi}_{N}=p_{N}$ und $\pi_{N}=p_{2N}$ gilt.

2**.   Beweis des Brent-Salamin-Algorithmus    2**

Wir beweisen, dass die Ausgabe $p_{N}$ des Brent-Salamin-Algorithmus gegen $\pi$ konvergiert.

3**.   Beweis der quadratischen Konvergenz    3**

Wir beweisen, dass $p_{N}$ quadratisch konvergiert: $|\pi-p_{n+1}|<0{,}075\cdot|\pi-p_{n}|^{2}$. Hierbei bezeichnet $\pi$ den Grenzwert von $p_{N}$.

1. Beweis der Äquivalenz der Algorithmen

Zwei Algorithmen, die exakt die gleichen Ergebnisse ausgeben, nennen wir äquivalent. Wir werden beweisen, dass die drei Algorithmen von S. 1 äquivalent sind. Genauer:

{theo}

Für die Ausgaben der drei Algorithmen auf S. 1, wobei

• •

$p_{N}$ die Ausgabe des Brent-Salamin Alg. 1 bezeichne,

• •

$\widehat{\pi}_{N}$ die Ausgabe des quadratischen Borwein-Alg. 2 bezeichne und

• •

$\pi_{N}$ die Ausgabe des quartischen Borwein-Alg. 3 bezeichne,

gilt:

[TABLE]

d.h. dass diese Algorithmen genau die gleichen Ergebnisse liefern, falls die Ergebnisse exakt berechnet werden.

Beweis.

Das wurde bereits von Brent [1] mit Hilfe elliptischer Modulfunktionen bewiesen und von Guillera [2], der eine Formel von Gauss verwendet. Wir benötigen nur elementare Algebra für den Beweis von $\widehat{\pi}_{N}=p_{N}$ in Satz 1.2 und von $\pi_{N}=\widehat{\pi}_{2N}$ in Satz 1.3. ∎

Bemerkung 1.1**.**

Bei einer tatsächlichen Implementierung der Algorithmen mit Hilfe von Gleitkomma- oder Intervallarithmetik kann man die Startwerte und die Iterationen nur mit einer endlichen Genauigkeit berechnen. Hier entstehen Rundungsfehler, die sich bei den verschiedenen Algorithmen unterschiedlich fortpflanzen. Die Ausgaben der Algorithmen unterscheiden sich also in den letzten Dezimalen. Um $D$ Dezimalen von $\pi$ korrekt zu berechnen, muss man von Anfang an alle Zwischenergebnisse auf einige zusätzliche Dezimalen berechnen (z.B. auf $D+50$ Dezimalen), die man am Ende wieder abschneidet.

Satz 1.2**.**

Für die Größen des Brent-Salamin-Algorithmus 1 und des Borwein-Alg. 2 auf S. 1 gilt $e_{n}=1/p_{n}$ und $k_{n}=a_{n}/a_{n+1}-1$. Insbesondere gilt

[TABLE]

d.h. diese beiden Algorithmen liefern genau die gleichen Ergebnisse.

Beweis.

Wir setzen $E_{n}:=1/p_{n}$ und $K_{n}:=a_{n}/a_{n+1}-1$ und beweisen dann per vollständiger Induktion, dass $E_{n}=e_{n}$ und $K_{n}=k_{n}$ gilt:

• •

Für den Induktionsanfang beweisen wir $K_{0}=k_{0}$ und $E_{0}=e_{0}$:

[TABLE]

• •

Beweise jetzt unter Verwendung der Induktionsvoraussetzungen ($K_{n-1}=k_{n-1}$ und $E_{n-1}=e_{n-1}$), dass $K_{n}=k_{n}$ und $E_{n}=e_{n}$ gilt:

[TABLE]

Hieraus folgt nun $K_{n}=k_{n}$ aufgrund der Induktionsvoraussetzung:

[TABLE]

Schließlich gilt $p_{n}:=\frac{(a_{n}+b_{n})^{2}}{1-2\cdot\sum_{j=1}^{n}2^{j}c_{j}^{2}}=\frac{4\cdot a_{n+1}^{2}}{1-2\cdot\sum_{j=1}^{n}2^{j}c_{j}^{2}}$ und somit

[TABLE]

Das liefert

[TABLE]

Mit $\frac{c_{n}^{2}}{a_{n+1}^{2}}=\frac{a_{n}^{2}-b_{n}^{2}}{(a_{n}^{2}+b_{n})^{2}/4}=4\cdot\frac{a_{n}-b_{n}}{a_{n}+b_{n}}=4\cdot\left(\frac{2a_{n}}{a_{n}+b_{n}}-\frac{a_{n}+b_{n}}{a_{n}+b_{n}}\right)=4\cdot\left(\frac{a_{n}}{a_{n+1}}-1\right)$ folgt:

[TABLE]

Hier ersetzen wir $a_{n}/a_{n+1}$ durch $K_{n}+1$ und erhalten

[TABLE]

Aber wir haben bereits $K_{n}=k_{n}$ bewiesen und nach Induktionsvoraussetzung gilt $E_{n-1}=e_{n-1}$:

[TABLE]

Hier erkennen wir die Definition von $e_{n}$, also ist auch $E_{n}=e_{n}$ bewiesen.

Somit haben wir für alle $N\geq 0$ bewiesen, dass $\widehat{\pi}_{N}=1/e_{N}=p_{N}$ gilt, d.h. dass die beiden Algorithmen genau die gleichen Ergebnisse liefern. ∎

Satz 1.3**.**

Für die Größen der Algorithmen 2 und 3 auf Seite 2 gilt $y_{n}=\sqrt{k_{2n}}$ und $z_{n}=e_{2n}$. Insbesondere gilt

[TABLE]

d.h. eine Iteration des Algorithmus 3 entspricht genau zwei Iterationen des Algorithmus 2.

Beweis.

Wir setzen $Y_{n}:=\sqrt{k_{2n}}$ und $Z_{n}:=e_{2n}$ und beweisen dann per vollständiger Induktion, dass $Y_{n}=y_{n}$ und $Z_{n}=z_{n}$ gilt:

• •

Für den Induktionsanfang erkennen wir zunächst $Z_{0}:=e_{0}=6-4\cdot\sqrt{2}=z_{0}$. Außerdem gilt $y_{0}^{2}=\left(\sqrt{2}-1\right)^{2}=3-2\cdot\sqrt{2}=k_{0}$, also folgt $Y_{0}:=\sqrt{k_{0}}=y_{0}$.

• •

Beweise jetzt unter Verwendung der Induktionsvoraussetzungen ($Y_{n-1}=y_{n-1}$ und $Z_{n-1}=z_{n-1}$), dass $Y_{n}=y_{n}$ und $Z_{n}=z_{n}$ gilt:

Aus $k_{n}:=\frac{1-\sqrt{1-k_{n-1}^{2}}}{1+\sqrt{1-k_{n-1}^{2}}}$ folgt $\sqrt{1-k_{n-1}^{2}}=\frac{1-k_{n}}{1+k_{n}}$ und somit

[TABLE]

Aus $Y_{n}:=\sqrt{k_{2n}}$ folgt $k_{2n}=Y_{n}^{2}$ und $k_{2n-1}=\frac{2\cdot\sqrt{k_{2n}}}{1+k_{2n}}=\frac{2\cdot Y_{n}}{1+Y_{n}^{2}}$. Hieraus folgt:

[TABLE]

Somit ist bewiesen, dass $Y_{n}=y_{n}$ ist, wobei wir im letzten Schritt die Induktionsvoraussetzung $Y_{n-1}=y_{n-1}$ benutzt haben. Wir müssen nun noch $Z_{n}=z_{n}$ beweisen:

Aus der Definition von $e_{n}$ in Algorithmus 2 folgt

[TABLE]

Wenn wir diese Darstellung von $e_{2n-1}$ in die für $e_{2n}$ einsetzen erhalten wir:

[TABLE]

Mit $k_{2n}=Y_{n}^{2}=y_{n}^{2}$ und $k_{2n-1}=\frac{2\cdot Y_{n}}{1+Y_{n}^{2}}=\frac{2\cdot y_{n}}{1+y_{n}^{2}}$ folgt:

[TABLE]

Hier nutzen wir die Induktionsvoraussetzung $z_{n-1}=Z_{n-1}=e_{2n-2}$:

[TABLE]

Somit haben wir für alle $N\geq 0$ bewiesen, dass $\pi_{N}=1/z_{N}=1/e_{2N}=\widehat{\pi}_{2N}$ gilt, d.h. dass also Algorithmus 3 genau jedes zweite Ergebnis von Algorithmus 2 produziert. ∎

Beweis des Thm. 1.

In Satz 1.2 haben wir $\widehat{\pi}_{N}=p_{N}$ bewiesen und in Satz 1.3 haben wir $\pi_{N}=\widehat{\pi}_{2N}$ bewiesen – somit sind beide Aussagen des Thm. 1 bewiesen, und die drei Algorithmen sind äquivalent. ∎

Bemerkung 1.4**.**

Die ersten Ausgaben der drei äquivalenten Algorithmen sind:

[TABLE]

2. Beweis des Brent-Salamin-Algorithmus

In diesem Kapitel beweisen wir, dass der Brent-Salamin-Algorithmus gegen $\pi$ konvergiert. Der vorliegende Beweis arbeitet [3] aus und setzt nur Integrationstechniken wie partielle Integration und Integration durch Substitution (auch zweidimensional – also den Transformationssatz) voraus.

{theo}

Es gilt die Formel von Gauß (1809), Brent (1976) und Salamin (1976)

[TABLE]

wobei $\operatorname{AGM}(1;1/\sqrt{2})$ das arithmetisch-geometrische Mittel (also den gemeinsamen Grenzwert der Folgen $a_{n}$ und $b_{n}$ des Brent-Salamin-Algorithmus auf S. 1) bezeichnet. Insbesondere konvergiert die Folge

[TABLE]

des Brent-Salamin-Algorithmus auf S. 1 gegen $\pi$.

Beweis.

Wir verallgemeinern den Brent-Salamin-Algorithmus zunächst auf die Startwerte

[TABLE]

Später (ab Satz 2.9) werden wir $a=1$ und $b=1/\sqrt{2}$ setzen. Auf Seite 2 wird der Beweis von Thm. 2 fortgesetzt, zunächst beweisen wir einige Hilfssätze: ∎

Satz 2.1**.**

Für das geometrische Mittel $\sqrt{x\cdot y}$ und das arithmetische Mittel $\frac{x+y}{2}$ zweier positiver reeller Zahlen $x\neq y$ gilt:

[TABLE]

Beweis.

Zunächst gilt (weil n.V. $x\neq y$ ist):

[TABLE]

Hieraus folgt $4xy<(x+y)^{2}$ und somit, dass das geometrische Mittel $\sqrt{x\cdot y}$ kleiner als das arithmetische Mittel $\frac{x+y}{2}$ ist. ∎

Satz 2.2**.**

Die Folgen $a_{n}$ und $b_{n}$ des Brent-Salamin-Algorithmus 1 konvergieren gegen einen gemeinsamen Grenzwert, den wir $\operatorname{AGM}(a,b)$ nennen. Die Konvergenz von $a_{n}\searrow\operatorname{AGM}(a,b)$ und $b_{n}\nearrow\operatorname{AGM}(a,b)$ erfolgt streng monoton und es gilt $c_{n+1}^{2}<\frac{1}{4}c_{n}^{2}$.

Beweis.

Aus Satz 2.1 folgt, dass $b_{n}<a_{n}$ für alle $n$ gilt. Hieraus folgt die strenge Monotonie $b_{n+1}=\sqrt{a_{n}\cdot b_{n}}>\sqrt{b_{n}\cdot b_{n}}=b_{n}$ und $a_{n+1}=\frac{a_{n}+b_{n}}{2}<\frac{a_{n}+a_{n}}{2}=a_{n}$.Beide Folgen sind durch $b=b_{0}\leq b_{n}<a_{n}\leq a_{0}=a$ beschränkt und somit konvergent. Für ihre Abweichung gilt:

[TABLE]

Hiermit ist bewiesen, dass $c_{n}^{2}=a_{n}^{2}-b_{n}^{2}<4^{-n}\cdot c_{0}^{2}$ gegen Null konvergiert, und dass $a_{n}$ und $b_{n}$ also gegen den selben Grenzwert konvergieren. ∎

Satz 2.3**.**

Das Integral

[TABLE]

bleibt konstant über die ganze AGM-Folge, d.h. es gilt $I(a_{n},b_{n})=I(a_{0},b_{0})$ für alle $n\in\mathbb{N}$.

Beweis.

Wir führen zunächst die Substitution $t=b\cdot\tan\Phi$ durch: Aus $1+\tan^{2}\Phi=\frac{1}{\cos^{2}\Phi}$ folgt dann $\cos^{2}\Phi=\frac{b^{2}}{b^{2}+b^{2}\tan^{2}\Phi}=\frac{b^{2}}{b^{2}+t^{2}}$ und $\sin^{2}\Phi=1-\cos^{2}\Phi=\frac{b^{2}+t^{2}-b^{2}}{b^{2}+t^{2}}=\frac{t^{2}}{b^{2}+t^{2}}$. Außerdem erhalten wir $\frac{\mathrm{d}t}{\mathrm{d}\Phi}=b\cdot(1+\tan^{2}\Phi)=b+\frac{t^{2}}{b}=\frac{t^{2}+b^{2}}{b}$, also $\frac{\mathrm{d}\Phi}{\mathrm{d}t}=\frac{b}{t^{2}+b^{2}}$. Die genannte Substitution liefert also folgende alternative Darstellung von $I(a,b)$:

[TABLE]

Nun substituieren wir $x=\frac{1}{2}\left(t-\frac{ab}{t}\right)$. Das führt auf $2xt=t^{2}-ab$ und (weil $t>0$) auf t=x\mathbin{\vbox{\hbox{\oalign{\hfil\scriptstyle+$\hfil\cr\kern-1.29167pt\cr$\scriptscriptstyle({-})\cr}}}}\sqrt{x^{2}+ab}. Somit gilt $\frac{\mathrm{d}t}{\mathrm{d}x}=1+\frac{2x}{2\sqrt{x^{2}+ab}}=\frac{\sqrt{x^{2}+ab}+x}{\sqrt{x^{2}+ab}}=\frac{t}{\sqrt{x^{2}+ab}}$, also:

[TABLE]

Hier haben wir $\frac{(t^{2}+a^{2})\cdot(t^{2}+b^{2})}{t^{2}}$ zu $f(x)$ zusammengefasst (beachte $t>0$). Für $f(x)$ gilt:

[TABLE]

Für $I(a,b)$ erhalten wir also:

[TABLE]

Hier nutzen wir, dass der Integrand eine gerade Funktion ist, weshalb $\frac{1}{2}\int\limits_{-\infty}^{\infty}$ in $\int\limits_{0}^{\infty}$ übergeht:

[TABLE]

Wir haben also für beliebige $a>b>0$ bewiesen, dass $I\mathopen{}\left(\frac{a+b}{2},\sqrt{ab}\right)\mathclose{}=I(a,b)$ gilt. Per vollständiger Induktion folgt hieraus $I(a_{n},b_{n})=I(a_{0},b_{0})$ für alle $n\in\mathbb{N}$. ∎

Satz 2.4**.**

Für das Integral $I(a,b)$ aus Satz 2.3 gilt:

[TABLE]

Beweis.

Mit $m:=\operatorname{AGM}(a,b)$ gilt nach Satz 2.2, dass $a_{n}\rightarrow m$ und $b_{n}\rightarrow m$ konvergieren. Vertauschen von Grenzwertbildung und Integration liefert dann mit Satz 2.3:

[TABLE]

∎

Satz 2.5**.**

Für das Integral

[TABLE]

gilt $L(b,a)+L(a,b)=I(a,b)$ und $L(b,a)-L(a,b)=\frac{a-b}{a+b}\cdot L(b_{1},a_{1})$.

Beweis.

Um die erste Gleichung zu beweisen, substituieren wir $\Phi^{\prime}=\frac{\pi}{2}-\Phi$. Dann gilt $\cos(\Phi^{\prime})=\sin(\Phi)$ und $\sin(\Phi^{\prime})=\cos(\Phi)$, also

[TABLE]

Schließlich folgt aus $\sin^{2}+\cos^{2}=1$, dass gilt:

[TABLE]

Als Nächstes beweisen wir eine alternative Darstellung von $L(a,b)$. Genau wie bei $I(a,b)$ in Gleichung (2.2) substituieren wir hierfür $t=b\cdot\tan\Phi$ und erhalten:

[TABLE]

Dann bilden wir die gesuchte Differenz $L(b,a)-L(a,b)$, wobei für den Ausdruck $L(b,a)$ die Variablen $a$ und $b$ vertauscht werden:

[TABLE]

Hier gilt $\frac{a^{2}}{a^{2}+t^{2}}-\frac{b^{2}}{b^{2}+t^{2}}=\frac{a^{2}(b^{2}+t^{2})-b^{2}(a^{2}+t^{2})}{(a^{2}+t^{2})(b^{2}+t^{2})}=\frac{a^{2}t^{2}-b^{2}t^{2}}{(a^{2}+t^{2})(b^{2}+t^{2})}=\frac{(a^{2}-b^{2})t^{2}}{(a^{2}+t^{2})(b^{2}+t^{2})}$, also:

[TABLE]

Und genau wie bei $I(a,b)$ substituieren wir nun $x=\frac{1}{2}\left(t-\frac{ab}{t}\right)$:

[TABLE]

wobei wieder $f(x):=\frac{(t^{2}+a^{2})\cdot(t^{2}+b^{2})}{t^{2}}=(2x)^{2}+(a+b)^{2}$ gilt, also:

[TABLE]

In diesem Integral erkennen wir die Darstellung (2.3) von $L(b_{1},a_{1})$. Es folgt:

[TABLE]

Somit ist auch die zweite Gleichung bewiesen. ∎

Satz 2.6**.**

Mit $S:=\sum_{j=1}^{\infty}2^{j}\cdot c_{j}^{2}$ gilt:

[TABLE]

Beweis.

Zunächst gilt $4\cdot(a_{1}^{2}-b_{1}^{2})=4\cdot\left(\frac{a+b}{2}\right)^{2}-4ab=(a+b)^{2}-4ab=(a-b)^{2}$. Daraus folgt (unter Nutzung beider Gleichungen aus Satz 2.5):

[TABLE]

Mit der Definition der $c_{n}$ können wir das wie folgt abkürzen:

[TABLE]

Also gilt für alle $j\in\mathbb{N}$:

[TABLE]

Hier multiplizieren wir noch mit $2^{j}$ und nutzen $I(a_{j},b_{j})=I(a,b)$ aus Satz 2.3:

[TABLE]

Nun summieren wir diese Gleichungen für $0\leq j\leq n$ und erhalten:

[TABLE]

In der zweiten Summe führen wir einen Indexshift $k=j+1$ durch:

[TABLE]

Somit erkennen wir, dass auf der linken Seite von (2.4) eine Teleskopsumme steht, in der sich fast alle Summanden gegenseitig auslöschen:

[TABLE]

Hier schätzen wir noch $L(b_{n+1},a_{n+1})<I(b_{n+1},a_{n+1})=I(b,a)$ und nutzen $c_{n+1}^{2}<4^{-n-1}c_{0}^{2}$ aus Satz 2.2:

[TABLE]

Also geht der zweite Ausdruck aus (2.5) gegen Null (für $n\rightarrow\infty$) und wir erhalten:

[TABLE]

Schließlich setzen wir noch $L(b,a)=I(a,b)-L(a,b)$ aus Satz 2.5 ein und erhalten

[TABLE]

∎

Satz 2.7**.**

Für die Gamma-Funktion $\Gamma(x):=\int_{0}^{\infty}t^{x-1}\cdot e^{-t}\mathrm{d}t$ gilt für $\operatorname{Re}(x)>0$:

[TABLE]

Beweis.

Die Funktionalgleichung folgt durch partielle Integration:

[TABLE]

Für $\Gamma(1/2)$ substituieren wir zunächst $s=\sqrt{t}$ und erhalten $t=s^{2}$ und $\frac{\mathrm{d}t}{\mathrm{d}s}=2s$:

[TABLE]

Dieses Integral quadrieren wir, um auf ein zweidimensionales Integral zu kommen:

[TABLE]

Hier bietet sich ein Übergang zu Polarkoordinaten an, also $x=r\cdot\cos\Phi$ und $y=r\cdot\sin\Phi$:

[TABLE]

Also ist auch der Wert $\Gamma(1/2)=\sqrt{\pi}$ bewiesen. ∎

Satz 2.8**.**

Für die Betafunktion $B(x,y):=\int_{0}^{1}t^{x-1}\cdot(1-t)^{y-1}\mathrm{d}t$ gilt im Bereich $\operatorname{Re}(u)>0$ und $\operatorname{Re}(v)>0$:

[TABLE]

Beweis.

Wir beginnen mit $\Gamma(u)\cdot\Gamma(v)$:

[TABLE]

Mit der Substitution $\left\{\begin{aligned} s&=x\cdot(1-y)\\ t&=x\cdot y\end{aligned}\right\}$ bzw. $\left\{\begin{aligned} x&=s+t\\ y&=\frac{t}{x}=\frac{t}{s+t}\end{aligned}\right\}$ gilt $0<x<\infty$ und $0<y<1$. Die Jacobi-Matrix der Substitution ist

[TABLE]

Also gilt:

[TABLE]

Eine abschließende Division durch $\Gamma(u+v)$ liefert die zu beweisende Gleichung. ∎

Satz 2.9**.**

Es gilt

[TABLE]

Beweis.

Zunächst gilt wegen $\sin^{2}\Phi=1-\cos^{2}\Phi$:

[TABLE]

Dann substituieren wir $x=\cos\Phi$, wobei $\frac{\mathrm{d}x}{\mathrm{d}\Phi}=-\sin\Phi=-\sqrt{1-\cos^{2}\Phi}=-\sqrt{1-x^{2}}$ und somit $\frac{\mathrm{d}\Phi}{\mathrm{d}x}=\frac{-1}{\sqrt{1-x^{2}}}$ gilt:

[TABLE]

Schließlich substituieren wir $x=t^{1/4}$ mit $\frac{\mathrm{d}x}{\mathrm{d}t}=\frac{1}{4}\cdot t^{-3/4}$, um auf die Betafunktion $B(x,y)$ aus Satz 2.8 zu kommen:

[TABLE]

Jetzt ersetzen wir die Betafunktionen durch Gammafunktionen mit Satz 2.8 und verwenden dann die Eigenschaften der Gammafunktion aus Satz 2.7:

[TABLE]

∎

Beweis des Thm. 2.

Zunächst gilt für $I(a,b)$ und ebenso für $L(a,b)$:

[TABLE]

Hieraus folgt $I(1,1/\sqrt{2})=\sqrt{2}\cdot I(\sqrt{2},1)$ und $L(1,1/\sqrt{2})=\sqrt{2}\cdot L(\sqrt{2},1)$ und dann mit Satz 2.9:

[TABLE]

Satz 2.6 liefert

[TABLE]

Hier setzen wir $a=1$ und $b=1/\sqrt{2}$ ein, also $c_{0}^{2}=a^{2}-b^{2}=1/2$:

[TABLE]

Nun verwenden wir Satz 2.4:

[TABLE]

Wenn wir diese Formel nach $\pi$ auflösen, erhalten wir die Gauß’sche Formel:

[TABLE]

Der Zähler der Folge $p_{N}$ ist $(a_{N}+b_{N})^{2}=4\cdot a_{N+1}^{2}$ und konvergiert (Satz 2.2) gegen$4\cdot\operatorname{AGM}(1;1/\sqrt{2})^{2}$, also konvergiert die Folge $p_{N}$ gegen $\pi$. ∎

Bemerkung 2.10**.**

Nun haben wir bewiesen, dass der Brent-Salamin-Algorithmus gegen $\pi$ konvergiert. Wir haben bereits die Äquivalenz der drei Algorithmen bewiesen, also folgt dass auch die beiden Algorithmen der Borwein-Brüder gegen $\pi$ konvergieren.

3. Beweis der quadratischen Konvergenz

Dieses Kapitel setzt nicht voraus, dass die Brent-Salamin-Folge gegen $\pi\approx 3{,}14159$ konvergiert. Wir setzen nur die monotone Konvergenz von $a_{n}$ und $b_{n}$ aus Satz 2.2 voraus und verwenden $\pi$ als Platzhalter für den Grenzwert von $p_{n}$.

{theo} Die Folge $\displaystyle p_{n}:=\frac{(a_{n}+b_{n})^{2}}{1-2\cdot\sum_{j=1}^{n}2^{j}c_{j}^{2}}$ des Brent-Salamin-Algorithmus konvergiert quadratisch gegen ihren Grenzwert, den wir mit $\pi$ bezeichnen:

[TABLE]

Insbesondere wird die Anzahl gültiger Stellen mit jeder Iteration ungefähr verdoppelt, wobei man von Anfang an mit der gewünschten Zielgenauigkeit rechnen muss.

Beweis.

Zunächst benennen wir Zähler und Nenner der Brent-Salamin-Folge mit $X_{n}$ und $Y_{n}$. In Satz 2.2 haben wir bewiesen, dass $a_{n}\rightarrow\operatorname{AGM}(1,1/\sqrt{2})$ gilt, also folgt:

[TABLE]

Wir bezeichnen jetzt die Abweichungen mit

[TABLE]

Weil für alle $n$ gilt $a_{n+1}>\operatorname{AGM}(a,b)>b_{n+1}$ (vgl. Satz 2.2) folgt

[TABLE]

Weiter werden bei $\delta_{n}$ nur positive Zahlen summiert, also gilt $\delta_{n}>2\cdot 2^{n+1}\cdot c_{n+1}^{2}$. Schließlich ist $c_{n+1}^{2}<c_{n}^{2}/4$ (vgl. Satz 2.2), also kann die Summe in $\delta_{n}$ mit der geometrischen Reihe abgeschätzt werden (setze hierfür den Summationsindex $j=n+1+k$):

[TABLE]

Insgesamt haben wir also bewiesen:

[TABLE]

Für die Abweichung zwischen dem Folgenglied $p_{n}$ und dem Grenzwert $\pi$ gilt:

[TABLE]

Weiter ist $1>Y_{n}>Y>0$ und $\delta_{n}>\varepsilon_{n}>0$, also gilt:

[TABLE]

Dann folgt mit (3.1):

[TABLE]

Im Beweis von Satz 2.2 haben wir $c_{n+1}^{2}=(a_{n}-b_{n})^{2}/4$ bewiesen. Hieraus folgt:

[TABLE]

In dieser Zeile erkennen wir die quadratische Konvergenz:

[TABLE]

Hieraus folgt mit (3.2):

[TABLE]

Hier nutzen wir $\operatorname{AGM}>b_{1}$ (Satz 2.2) bzw. $\operatorname{AGM}(1,1/\sqrt{2})>\sqrt{1\cdot 1/\sqrt{2}}=\sqrt[4]{1/2}$ und $2^{n}\geq 1$ und erhalten:

[TABLE]

Somit ist die quadratische Konvergenz von $p_{n}$ gegen $\pi$ bewiesen. ∎

Bemerkung 3.1**.**

Mit $p_{1}>3{,}14057$ folgt aus Thm. 3: $|\pi-p_{44}|<10^{-3{,}6\cdot 10^{13}}$. Insbesondere liegt $p_{44}=\widehat{\pi}_{44}=\pi_{22}$ näher an $\pi$ als der aktuelle Rekord von $\pi\cdot 10^{13}$ Dezimalen (März 2019).

Bemerkung 3.2**.**

Tatsächlich gilt sogar $|\pi-p_{44}|<10^{-4{,}8\cdot 10^{13}}$, was aus der in [1, Glg. (20)] bewiesenen Fehlerabschätzung von $0<\pi-p_{n}<\left(2^{n+4}\pi^{2}-8\pi\right)\cdot\exp\mathopen{}\left(-2^{n+1}\pi\right)\mathclose{}$ folgt.

Literatur

• [1]

Richard P. Brent.

The Borwein brothers, $\pi$ and the AGM, 2018.

https://arxiv.org/abs/1802.07558.

• [2]

Jesus Guillera.

Easy Proofs of Some Borwein Algorithms for $\pi$, 2008.

https://arxiv.org/abs/0803.0991.

• [3]

Nick Lord.

Recent Calculations of $\pi$: the Gauss-Salamin Algorithm.

The Mathematical Gazette, 76(476):231–242, 1992.

https://doi.org/10.2307/3619132.