# Continued fraction approximation for $\pi$

How to code this formula in Mathematica to approximate $\pi$? Thanks for your help!

$$\frac4{\pi}=1+\cfrac{1^2}{2+\cfrac{3^2}{2+\cfrac{5^2}{2+\cfrac{7^2}{2+\cdots}}}}$$

1 + ContinuedFractionK[(2 n - 1)^2, 2, {n, 1, Infinity}]
(*  4/π  *)


Pick a termination point less than Infinity to get an approximation.

Michael's method of using ContinuedFractionK[] is the canonical way. If you want to look at the forward recursion method manually, you can use repeated matrix multiplication to implement it:

al = Range[1, 45, 2]^2; (* partial numerators *)
bl = ConstantArray[2, Length[al]]; (* partial denominators *)

1 + Divide @@@ FoldList[{{0, 1}, #2}.#1 &, IdentityMatrix,
Transpose[{al, bl}]][[All, 2]]
{1, 3/2, 15/13, 105/76, 315/263, 3465/2578, 45045/36979, 45045/33976, 765765/622637,
14549535/11064338, 14549535/11757173, 334639305/255865444, 1673196525/1346255081,
5019589575/3852854518, 145568097675/116752370597, 4512611027925/3473755390832,
4512611027925/3610501179557, 4512611027925/3481569435902,
166966608033225/133330680156299, 166966608033225/129049485078524,
6845630929362225/5457995496252709, 294362129962575675/227848175409504262,
294362129962575675/234389556075339277, 13835020108241056725/10721947005578370344}


The other methods I mentioned in this answer can also be adapted to this case.

• Every time they add another 1000 functions, I feel it just takes some of the fun out of programming. :) – Michael E2 Jan 31 '17 at 20:01
• On the other hand, one usually hopes that the new built-ins outdo the stuff we can come up with, in either performance or generality. – J. M. is in limbo Feb 1 '17 at 8:08
pi[n_] := 1 + Fold[(2*#2 - 1)^2/(2 + #1) &, 2*n - 1, Reverse[Range[n]]]

pi /@ Range // N


{1.33333, 1.26316, 1.2766, 1.27176, 1.27401, 1.27279, 1.27353, 1.27305, 1.27338, 1.27314}

DiscretePlot[pi[n], {n, 10}, AxesOrigin -> {0, 4/Pi}, PlotRange -> All, FillingStyle -> Red] 1 + Fold[Inactivate[(2*#2 - 1)^2, Power]/(2 + #1) &, 2*4 - 1, Reverse[Range]] A quick hack :) (someone good in Mathematica can make this more functional )

rest[n_] := If[n < 1024, (n + 2)^2/(2 + rest[n + 2]), n]
1 + 1/(2 + rest) // N 4/Pi // N Showing speed of convergence

rest[n_, max_] := If[n < max, (n + 2)^2/(2 + rest[n + 2, max]), n];
data = Table[{i, 1 + 1/(2 + rest[1, i])}, {i, 3, 100}];
Show[ListLinePlot[data], Plot[4/Pi, {x, 0, 100}, PlotStyle -> Red],
PlotRange -> All] Lord Brouncker's Formula is the same as:

(Pi^2)/8 = 1/1^2 + 1/3^2 + 1/5^2 + ...

so for something simple in Mathematica

N[Sqrt[8 Sum[1/n^2, {n, 1, 99999, 2}]], 6]

3.14159