Compute π using integrals

In our previous post, we looked at

\begin{equation*} I = {\int_{0}^{1}}{\frac{x^{4}(1-x)^{4}}{x^2+1}}{\,dx}. \end{equation*}

Let's now consider the generalized integral

\begin{equation*} I_n = {\int_{0}^{1}}{\frac{x^{4n}(1-x)^{4n}}{x^2+1}}{\,dx}. \end{equation*}

How to solve

First, we expand the integrand. Let $c_k$ be the coefficients of ${x^{4n}(1-x)^{4n}}$ expanded,

\begin{equation*} {x^{4n}(1-x)^{4n}}= c_{8n} x^{8n} + c_{8n-1} x^{8n-1} + \ldots + c_1 x + c_0. \end{equation*}

Since ${x^{4n}(1-x)^{4n}}$ has no term of degree less than $4n$, for all $k < 4n$, $c_k = 0$. Hence, $c_k$ may be described as follows,

\begin{equation*} c_k = \begin{cases} (-1)^k \binom{4n}{k-4n} & k \geq 4n,\\ 0 & k < 4n. \end{cases} \end{equation*}

Next, we try to long divide this whole thing by $x^2 + 1$, and hopefully, we find an expression we can work with.

Let $q_n(x)$ be the polynomial obtained by dividing ${x^{4n}(1-x)^{4n}}$ with $x^2+1$ and let $r_n(x)$ be the remainder. That is,

\begin{equation*} {\frac{x^{4n}(1-x)^{4n}}{x^2+1}}= q_n(x) + r_n(x). \end{equation*}

It follows that

\begin{equation*} \label{eq:qx} x^2 q_n(x) + q_n(x) + (x^2+1)r_n(x) = {x^{4n}(1-x)^{4n}}. \tag{*} \end{equation*}

Since ${x^{4n}(1-x)^{4n}}$ is of degree $8n$ and $x^2 + 1$ is of degree 2, $q_n$ is of degree $8n-2$. Hence, we may write

\begin{equation*} q_n(x) = z_{8n-2} x^{8n-2} + z_{8n-3} x^{8n-3} + \ldots + z_1 x + z_0. \end{equation*}

Note that if we find an expression for $z_k$ and the remainder, we can compute the integral.

Working out the quotient polynomial

Here's how we determine $z_k$. We rewrite $\eqref{eq:qx}$ in such a way that the terms on the LHS align nicely with the terms on the RHS:

\begin{equation*} \begin{array}{l l l l} &c_{8n} x^{8n} &+ c_{8n-1} x^{8n-1} &+ c_{8n-2} x^{8n-2} + \ldots + c_3 x^3 + c_2 x^2 + c_1 x + c_0 \\ =&z_{8n-2} x^{8n} &+ z_{8n-3} x^{8n-1} &+ z_{8n-4} x^{8n-2} + \ldots + z_1 x^3 + z_0 x^2 \\ & & &+ z_{8n-2} x^{8n-2} + \ldots + z_3 x^3 + z_2 x^2 + z_1 x + z_0 \\ & & &+ (x^2+ 1) r_n(x) . \\ \end{array} \end{equation*}

The first line is simply our original integrand. The second line is the expression $x^2 q_n(x)$ and the third line is $q_n(x)$.

From inspecting the terms, we observe the following,

\begin{equation*} \left \{ \begin{array}{l} z_{8n-2} = c_{8n} \\ z_{8n-3} = c_{8n-1} \\ z_{k} + z_{k-2} = c_{k} \end{array} \right . \end{equation*}

which gives us a recursive definition of $z_k$.

Zero coefficients

Some interesting patterns emerge if we inspect this sequence. Take a look at the coefficients $z_k$ for $n=2$ (largest index left),

[1, -8, 27, -48, 43, -8, -15, 0, 16, 0, -16, 0, 16, 0, -16]


and $n=3$,

[1, -12, 65, -208, 430, -584, 494, -208, 1, -12, 65, 0, -64, 0, 64, 0, -64, 0, 64, 0, -64, 0, 64]


Notice that for $k < 4n$, the coefficients seem to be either zero or some power of two. In fact, it can be proven that the element in the middle ($z_{4n-1}$) is always zero [1]. Because for $k < 4n$, $c_k = 0$, therefore all other odd elements to the right are also zero. Hence, we now know that $z_1 = 0$. Similarly, it can be shown using the properties of the binomial coefficient that $z_{4n-2} = (-4)^n$ [2], which means that $z_0 = -(-4)^n$.

 [1] To prove: $\sum_{k=0}^{2n-1} (-1)^k \binom{4n}{2k+1} = 0$.
 [2] To prove: $\sum_{k=0}^{2n} (-1)^k \binom{4n}{2k} = (-4)^n$.

Remainder term

Now consider the term $(x^2 + 1) r_n(x)$ in $\eqref{eq:qx}$. We see that this term is of degree less than 2. Since $x^{4n}(1-x)^{4n}$ has no term of degree lower than $4n$, it follows from the division that

\begin{equation*} (x^2 + 1) r_n(x) + z_0 + z_1 x = 0. \end{equation*}

Since $z_1 = 0$,

\begin{equation*} r_n(x) = \frac{z_0}{x^2+1} = \frac{(-4)^n}{x^2+1}. \end{equation*}

Thus, for our integrand we have,

\begin{equation*} \frac{x^{4n}(1-x)^{4n}}{x^2+1} = q_n(x) + r_n(x) = \sum_{k=0}^{8n-2} z_k x^k + \frac{(-4)^n}{x^2+1}. \end{equation*}

\begin{equation*} I_n = \sum_{k=0}^{8n-2} \frac{z_k}{k+1} - (-4)^{n-1} \pi. \end{equation*}

Notice that the quotient polynomial gives us a fraction and the remainder gives us an expression with π. Since the value of this integral is very small, we may use it to compute π:

\begin{equation*} \label{eq:pi} \pi \approx -\frac{1}{(-4)^{n-1}} \sum_{k=0}^{8n-2} \frac{z_k}{k+1}. \tag{**} \end{equation*}

Error bounds

Exactly how accurate is this estimate of π? Like in our previous post, we see that in the interval [0, 1],

\begin{equation*} \frac{x^{4n}(1-x)^{4n}}{x^2+1} \leq x^{4n}(1-x)^{4n}. \end{equation*}

So

\begin{equation*} \begin{aligned} I_n < {\int_{0}^{1}}{x^{4n}(1-x)^{4n}}{\,dx} = \sum_{k=0}^{8n} \frac{c_k}{k+1}. \end{aligned} \end{equation*}

Now we are ready to compute rational approximations to π.

Table

With $\eqref{eq:pi}$, we may program a simple routine using a symbolic library to find rational approximations to $\pi$. If we look at the error, we see that for every step increase in $n$, we obtain about two or three extra digits. The costs for this is that we need to compute four additional rows in Pascal's triangle. Certainly, it's possible to optimize the code (or perhaps find a better expression for the integral), however, I doubt it would be able to compete with any of the current record-breaking algorithms. Still, it remains an interesting way to compute $\pi$.

n rational approximation error
1 22/7 1.59e-03
2 47171/15015 4.57e-06
3 431302721/137287920 1.48e-08
4 741269838109/235953517800 5.04e-11
5 26856502742629699/8548690331301120 1.77e-13
6 211938730834003723543/67462193289708771840 6.33e-16
7 506119433541064524255449/161102819285860855603200 2.29e-18
8 19809071774292917047896724979/6305423381881718760060595200 8.39e-21
9 23285066731814401580687501596643/7411866941185812791748757094400 3.10e-23
10 89293478252053341114758995682016773/28422996899365886608045972478361600 1.15e-25
11 77997614500552946262295151821844788921/24827411794278189209115835981312819200 4.28e-28
12 53385652953509417766719017367204093188361/16993181115478931076058403260828562227200 1.60e-30
13 2504492442679095809835466877865587690003305358231/797204704377346869713059903864272972590717337600 6.02e-33
14 27866240172453844886325941887505217163994611104897858617/8870099737663958700556100446465640958076561638241075200 2.27e-35
15 3737271881299351053666441876167156192604573934552162661/1189610587174277674549768534026398549974217301531033600 8.55e-38
16 2199496269931868203800149756151730543362676840844586435379483/700121407343685092214243286289756273641462048387187186073600 3.24e-40
17 132570645054700093606473440370784391628703685360477563188116477/42198546938673298440864239484207883524155412032193942257664000 1.23e-42
18 29832337484492518217048916227659410964608663939042391473810429697/9495927949285245581337958144499514323433226651874145194402119680 4.66e-45
19 842755908791452603343975605728390570004752202303444739958124097773741955147/268257537408124350364026727600009545561127583382712382710348191778209792000 1.77e-47
20 66156338840129029362502085049678659745373047880820412086708394806039520515257/21058216686537761503576098116600749326548515295542922042762333054589468672000 6.75e-50
21 19208978128716291288284472076131229833907197187880665999064419495161046612530817093/6114407641858606996147679646207882639125089438506388063211023659509119218352128000 2.57e-52
22 194161760950485119583603495277389374815538366182218507109644401986701079185812018101/61803608029393290491398343929534957552325061335734232602973448135892175829545779200 9.82e-55
23 11376531846179208486779499075860814093243136394380595776565666024223647991552013566339271/3621262557123573816777089765864308356345692454910121468879503589693158929508344004608000 3.75e-57
24 33072589649906601294503039619501915430104748241262345867203946584680606098687119416797781273/10527332247264983629950069046266828718284283446809533417225181736459261497727838593220608000 1.44e-59
25 25422588845948287917840198990310955266893220891255517493881289524342496083962053135792429098585663/8092261362051105765248285301685603166696861085875287087233705766880133617989065376859671442227200 5.49e-62
26 46264531052986974409042524288673990665877753333636166970847031386641299179820533992829502611931261/14726457613822733014235618296761187744799693147268540465055933017205288205709830685816699381350400 2.10e-64
27 15118337974655890138153621692618786255874458647572745177428650890392802385633565461429689136646412227/4812316440000796768001141835184294518755566891467916410309369717790009488585371350904855292228403200 8.07e-67
28 2069928648601327427352439845133642469756287671058225941488213320736310619376143943839757508672264691186711/658878752544856169210929998187171572841546289913934734725669374808405315457981429539219541105191367475200 3.10e-69
29 1945427904999589890101598798890045028016567124243525244824182818377188279480704038637893031935822360758977292080733/619248935019189789053179275651408416868079582257816933219620037680829229923216655623409469247302624426776854528000 1.19e-71
30 231732713894571469968102101205967780695526360316041048788287205576407451806292974972467047521040580190265939101075713/73762813784841971728351645611122952962850797213242891037863061044561149361998215832515086874811852933591278498611200 4.56e-74
31 3276391597517381156722339574917443128660482033055065708467996703909632825271907608810720762683053048665415846495441582581/1042907836499285690943254466747130390957426471532303488753999865755368570446172107583986641627339851343735889599030886400 1.75e-76
32 14744914045249217409563438284635300992572088849462541609644118447684565944013718894666965951091297717732625917014942453716139/4693451911533054889370367338260552547336490092139263774838264473976963188770821966509512382626726442713299535344326148096000 6.75e-79
33 1809472427088685124865456761821112518606869891493260981446727916023734462133253943937697399103930900314412247151442840111929725959/575972962319306826489522834628755076951186522449311732752509519996508448485144742197277869306776632264505425508690335184715776000 2.59e-81
34 1145654685711889249769129784082023880664204709662470649074156123525837248838112856911829904235801637027201428011632652339566213475573/364673212614877939976015352622244473559869477243189812857522683251704159872045205282758919254541505428468725430181426036706639872000 9.99e-84
35 239974813888444249852898360727627543443755211013555059862043452073563515999052998410397597193481954205079096238042093742080440181303577497/76386355695807038084739543945241075049275490414689467699023939649313171834113902960821522439577785262595061187538615894965340914843648000 3.85e-86
36 38334433522726510759306283525375225056461436428102500322971275873573088154030933111430886139038191111696188535614112011156191123248109114121841/12202229171539165329134369214675841899577802410110905908955711336628981423701850335790621699407216015729407672254465241804476284875302764544000 1.48e-88
37 823469961672525301223173116537926628023498499469537047770337561775286582393238608924167522668085624424371918731695423872069637824319508091579/262118629775751995211437007297193775745344922693137328508794761264213674109249005725877301871080393456839150170381790884403688867245326336000 5.71e-91
38 28834831775394475155689761741074359292575367196542822593393645190239848906875523778477454385726190784244922494186830650413815442505610003720392059/9178412020554566087515748357688444107072715126868953788508119241459957495048314255675233689032508472509920902319913255605981519885260812189696000 2.20e-93
39 122298415725961133236587172274488877313768314029105795817058477350540613214149467972067973243342549201698492646394496552458046338144142727481148733477/38928794790188603515248339183013118605093418459858607126471663991864981609878489115047585504548399919796453257144725355266030530810245413326553088000 8.48e-96
40 48342614467886085771155676200886759726850439375211328930118361432995143143171994316294513521193274028793144221482388346320690176908109238458110886946384670429/15387932109099692199746965253984974932125901210095307633202237616558734891392454997644739224625470055554866306804514274588496291793420876440185453754187776000 3.27e-98