The following integral expression provides an approximation to \(\pi\) for arbitrary positive $$n$$ and is exact in the limit \(n\rightarrow\infty\),

\begin{equation}
\pi = \lim_{n\rightarrow \infty} \frac{2n-1}{2} \left[\int_{-\infty}^\infty \frac{dx}{(1+x^2)^n}\right]^2.
\label{eq:pi}
\end{equation}

The proof relies on identifying the integral as a beta integral. First we note that the integral is an even function of $$x$$ and thus we can replace the limits with \((0,\infty)\) with the addition of a factor of 2. Then by employing the change of variable \(t = 1/(1+x^2)\) we arrive at

\begin{equation*} \int_{-\infty}^\infty \frac{dx}{(1+x^2)^n} = \int_0^1 t^{n-3/2} (1-t)^{-1/2}\,dt = B(n-1/2,1/2). \end{equation*}

Now the beta function can be defined in terms of gamma functions via \(B(x,y)=\Gamma(x)\Gamma(y)/\Gamma(x+y)\) and that combined with the value of \(\Gamma(1/2) = \sqrt{\pi}\), the right-hand-side of \eqref{eq:pi} becomes

\begin{equation*} \lim_{n\rightarrow \infty} \frac{2n-1}{2} \left[\frac{\Gamma(n-1/2)\sqrt{\pi}}{\Gamma(n)}\right]^2 = \pi \lim_{n\rightarrow \infty} \frac{\Gamma(n-1/2)\Gamma(n+1/2)}{\Gamma^2(n)}, \end{equation*}

where we have used the identity \(x\Gamma(x) = \Gamma(x+1)\). Now as we have a constant coefficient of \(\pi\) in this expression we now just require to show that what remains inside the limit goes to unity. To achieve this last step we use Stirling’s approximation for the gamma function which is given by

\begin{equation*} \Gamma(n) = \sqrt{\frac{2\pi}{n}}\left(\frac{n}{e}\right)^n\left(1+\mathcal{O}\left(\frac{1}{n}\right)\right). \end{equation*}

Substituting this for the four gamma functions above we can simplify the limit to

\begin{equation*} \lim_{n\rightarrow \infty} \frac{\Gamma(n-1/2)\Gamma(n+1/2)}{\Gamma^2(n)} = \lim_{n\rightarrow\infty} \frac{(1-1/4n^2)^n}{(1-1/2n)}. \end{equation*}

The denominator is clearly going to one. The numerator is a little trickier but it is actually going to one as well since \(e^x = \lim_{n\rightarrow\infty}(1+x/n)^n\) and so the numerator is going to \(e^{1/2}e^{-1/2}=1\).

Does it work?

Here is a plot for \(n\in[10,50]\)

And here are some value for the integral expression for values of \(n\) between 1 and 1,000,000.

So it looks to converge, albeit slowly (initially at around one significant digit for every digit in \(n\)).