Instantly share code, notes, and snippets.

# StevenClontz/pi_approx.md

Last active August 29, 2015 13:58
Finding optimal rational approximations of pi

## I've got a beef with Pi "Approximation" Day

March 14 (3.14), often known as Pi Day, is a cute idea. I mean, June 28 would certainly be better, but that's not what this is about.

The chip on my shoulder is the shaft that July 22 (22/7) gets. Pi Approximation Day.

It's an apt description, but it has an odd implication. If 22/7 is a pi approximation, then I guess 3.14 = 157/50 = pi? Yeah okay.

I guess that, if you must, calling July 22 "Pi Inverse Day" would work for me. I'm not a fan of using the European "little endian" DD-MM-YYYY (and don't get me started on the nonsensical US convention MM-DD-YYYY), but the "big endian" ordering YYYY-MM-DD makes me happiest, and 07/22 is pretty close to 1/pi as long as you're not going to build a bridge with it.

But that's not what we've got here. March 31 is "Pi", and July 22 is the "Pi Approximation".

### Comparing the approximations

Which is actually the better approximation? You could break out a cell phone calculator and make short work of this in a couple of minutes, but let's use this simple question to explain the PiApprox Ruby class I've defined in pi_approx.rb.

Each PiApprox object represents a rational approximation of pi, initialized based upon our desired denominator. PiApprox.new 1000 represents 3.142, since 3142/1000 is the best approximation of pi with a denominator of 1000. Upon initialization we set four variables:

• @den is the denominator we initialized with
• @num is the numerator for which @num/@den is the best approximation of pi
• @to_f gives a Float representation of @num/@den
• @err measures the error between @to_f and Math::PI

We also can #inspect to get "#{@num}/#{@den}", and PiApprox is Comparable based upon its @err from Math::PI.

\$ irb -I .
irb(main):001:0> require 'pi_approx'
=> true
irb(main):002:0> march = PiApprox::PiApprox.new 100
=> 314/100
irb(main):003:0> july = PiApprox::PiApprox.new 7
=> 22/7
irb(main):004:0> july < march
=> true
irb(main):005:0> march.err - july.err
=> 0.00032816432244331395

So July's approximation of 22/7 is better than March's 3.14, y'know, by about three ten-thousandths. But it's something.

### How great are those decimal approximations anyway?

I know way too many digits of pi. Repeated watchings of Hard 'n Phirm's Pi have permantly engrained their circular hook in my mind:

3.14159265358979...

Expressed as a fraction in lowest terms, this decimal is pretty long: 314_159_265_358_979/100_000_000_000_000. Can we get those fifteen significant digits using a fraction with a denominator smaller than 100 trillion?

The #puts_improvements_up_to method serves to figure this out. Starting with 3/1 (only off by 14 hundredths!), this method runs through the denominators 2 through whatever and writes down the corresponding PiApprox object if its error is less than the previous object it wrote down.

It takes about five minutes or so, but by checking through 100_000_000 I got my answer:

\$ irb -I .
irb(main):001:0> require 'pi_approx'; include PiApprox
=> Object
irb(main):002:0> puts_improvements_up_to 100_000_000
... checking 10
... checking 100
... checking 1000
... checking 10000
... checking 100000
... checking 1000000
... checking 10000000
... checking 100000000
3/1, error=0.14159265358979312
13/4, error=0.10840734641020688
16/5, error=0.05840734641020706
19/6, error=0.025074013076873403
22/7, error=0.0012644892673496777
179/57, error=0.0012417763968106676
201/64, error=0.000967653589793116
223/71, error=0.0007475831672580924
245/78, error=0.0005670125641521473
267/85, error=0.00041618300155787935
289/92, error=0.0002883057637061981
311/99, error=0.00017851217565167943
333/106, error=8.32196275291075e-05
355/113, error=2.667641894049666e-07
52163/16604, error=2.662132572162079e-07
52518/16717, error=2.6261055063869776e-07
52873/16830, error=2.5905622225153024e-07
53228/16943, error=2.555493043843171e-07
53583/17056, error=2.5208885512384427e-07
53938/17169, error=2.4867395653771496e-07
54293/17282, error=2.4530371511843896e-07
54648/17395, error=2.419772608952542e-07
55003/17508, error=2.386937456577698e-07
55358/17621, error=2.354523434000555e-07
55713/17734, error=2.3225224943246303e-07
56068/17847, error=2.2909267860526938e-07
56423/17960, error=2.259728666409444e-07
56778/18073, error=2.2289206702552633e-07
57133/18186, error=2.1984955322906785e-07
57488/18299, error=2.1684461559701163e-07
57843/18412, error=2.1387656268245792e-07
58198/18525, error=2.1094471902571854e-07
58553/18638, error=2.080484260424953e-07
58908/18751, error=2.0518704113570152e-07
59263/18864, error=2.023599372513729e-07
59618/18977, error=1.9956650154639988e-07
59973/19090, error=1.9680613672079517e-07
60328/19203, error=1.9407825835315862e-07
60683/19316, error=1.91382296677034e-07
61038/19429, error=1.8871769480455214e-07
61393/19542, error=1.8608390828234178e-07
61748/19655, error=1.8348040597970794e-07
62103/19768, error=1.809066687563643e-07
62458/19881, error=1.78362189018344e-07
62813/19994, error=1.7584647027391043e-07
63168/20107, error=1.7335902757764643e-07
63523/20220, error=1.708993875304543e-07
63878/20333, error=1.6846708605910976e-07
64233/20446, error=1.6606167019261875e-07
64588/20559, error=1.6368269628586063e-07
64943/20672, error=1.613297313518558e-07
65298/20785, error=1.5900235039723043e-07
65653/20898, error=1.5670013864266252e-07
66008/21011, error=1.5442269019061428e-07
66363/21124, error=1.5216960758124287e-07
66718/21237, error=1.4994050179240048e-07
67073/21350, error=1.4773499223963427e-07
67428/21463, error=1.4555270588800795e-07
67783/21576, error=1.4339327858436945e-07
68138/21689, error=1.4125635239281564e-07
68493/21802, error=1.3914157737104915e-07
68848/21915, error=1.370486111262892e-07
69203/22028, error=1.3497711837118231e-07
69558/22141, error=1.3292676959153482e-07
69913/22254, error=1.3089724326675878e-07
70268/22367, error=1.2888822364942598e-07
70623/22480, error=1.2689940165344638e-07
70978/22593, error=1.2493047396588963e-07
71333/22706, error=1.2298114349107436e-07
71688/22819, error=1.2105111935056811e-07
72043/22932, error=1.1914011599500895e-07
72398/23045, error=1.1724785364819468e-07
72753/23158, error=1.153740578629936e-07
73108/23271, error=1.1351845952134454e-07
73463/23384, error=1.1168079572243528e-07
73818/23497, error=1.0986080622998884e-07
74173/23610, error=1.0805823880133403e-07
74528/23723, error=1.062728434142457e-07
74883/23836, error=1.0450437626374764e-07
75238/23949, error=1.0275259754166655e-07
75593/24062, error=1.0101727232481039e-07
75948/24175, error=9.929816968679006e-08
76303/24288, error=9.759506314210853e-08
76658/24401, error=9.590773109025008e-08
77013/24514, error=9.423595459523426e-08
77368/24627, error=9.257952005015113e-08
77723/24740, error=9.093821651262601e-08
78078/24853, error=8.931183881344396e-08
78433/24966, error=8.770018311565764e-08
78788/25079, error=8.610305135547947e-08
79143/25192, error=8.452024724547869e-08
79498/25305, error=8.295157938320585e-08
79853/25418, error=8.139685903074678e-08
80208/25531, error=7.985590100290096e-08
80563/25644, error=7.832852322309236e-08
80918/25757, error=7.681454761154782e-08
81273/25870, error=7.531379786485104e-08
81628/25983, error=7.382610167638859e-08
81983/26096, error=7.23512889599931e-08
82338/26209, error=7.08891940703893e-08
82693/26322, error=6.943965269456953e-08
83048/26435, error=6.800250362815063e-08
83403/26548, error=6.65775887753739e-08
83758/26661, error=6.516475270501587e-08
84113/26774, error=6.376384265038837e-08
84468/26887, error=6.237470762116004e-08
84823/27000, error=6.099720062380243e-08
85178/27113, error=5.963117555296549e-08
85533/27226, error=5.827648985601286e-08
85888/27339, error=5.693300231257581e-08
86243/27452, error=5.5600575699088495e-08
86598/27565, error=5.427907279198507e-08
86953/27678, error=5.29683608085918e-08
87308/27791, error=5.166830785441334e-08
87663/27904, error=5.037878381131122e-08
88018/28017, error=4.909966211386063e-08
88373/28130, error=4.783081708481518e-08
88728/28243, error=4.657212482328532e-08
89083/28356, error=4.5323464981095185e-08
89438/28469, error=4.40847172100689e-08
89793/28582, error=4.285576471474428e-08
90148/28695, error=4.1636491143748344e-08
90503/28808, error=4.042678281024337e-08
90858/28921, error=3.922652780374847e-08
91213/29034, error=3.8035615101961184e-08
91568/29147, error=3.685393679120352e-08
91923/29260, error=3.568138584597591e-08
92278/29373, error=3.451785612895719e-08
92633/29486, error=3.336324505553989e-08
92988/29599, error=3.2217449597027326e-08
93343/29712, error=3.108036938925807e-08
93698/29825, error=2.995190540033832e-08
94053/29938, error=2.8831960374731125e-08
94408/30051, error=2.7720437500988737e-08
94763/30164, error=2.6617243076287878e-08
95118/30277, error=2.5522282953716058e-08
95473/30390, error=2.443546609498526e-08
95828/30503, error=2.335670146180746e-08
96183/30616, error=2.2285899792251485e-08
96538/30729, error=2.122297360074299e-08
96893/30842, error=2.0167835845796844e-08
97248/30955, error=1.9120402150463178e-08
97603/31068, error=1.808058769370291e-08
97958/31181, error=1.7048309430833797e-08
98313/31294, error=1.6023486537619647e-08
98668/31407, error=1.5006037745735057e-08
99023/31520, error=1.3995884451389884e-08
99378/31633, error=1.2992948050793984e-08
99733/31746, error=1.1997151716514054e-08
100088/31859, error=1.1008419065205999e-08
100443/31972, error=1.0026675489882564e-08
100798/32085, error=9.051846827645704e-09
101153/32198, error=8.083861136043424e-09
101508/32311, error=7.1226455844453085e-09
101863/32424, error=6.168130006756201e-09
102218/32537, error=5.220244680970154e-09
102573/32650, error=4.2789203291704325e-09
102928/32763, error=3.3440890057079287e-09
103283/32876, error=2.415684541290375e-09
103638/32989, error=1.493639878447084e-09
103993/33102, error=5.778906242426274e-10
104348/33215, error=3.3162805834763276e-10
208341/66317, error=1.2235634727630895e-10
312689/99532, error=2.914335439641036e-11
833719/265381, error=8.715250743307479e-12
1146408/364913, error=1.6107115641261771e-12
3126535/995207, error=1.1426415369442111e-12
4272943/1360120, error=4.04121180963557e-13
5419351/1725033, error=2.220446049250313e-14
42208400/13435351, error=2.0872192862952943e-14
47627751/15160384, error=1.5987211554602254e-14
53047102/16885417, error=1.199040866595169e-14
58466453/18610450, error=8.881784197001252e-15
63885804/20335483, error=6.217248937900877e-15
69305155/22060516, error=3.9968028886505635e-15
74724506/23785549, error=2.220446049250313e-15
80143857/25510582, error=4.440892098500626e-16
245850922/78256779, error=0.0
=> [3/1, 13/4, 16/5, 19/6, 22/7, 179/57, 201/64, 223/71, 245/78, 267/85, 289/92, 311/99, 333/106, 355/113, 52163/16604, 52518/16717, 52873/16830, 53228/16943, 53583/17056, 53938/17169, 54293/17282, 54648/17395, 55003/17508, 55358/17621, 55713/17734, 56068/17847, 56423/17960, 56778/18073, 57133/18186, 57488/18299, 57843/18412, 58198/18525, 58553/18638, 58908/18751, 59263/18864, 59618/18977, 59973/19090, 60328/19203, 60683/19316, 61038/19429, 61393/19542, 61748/19655, 62103/19768, 62458/19881, 62813/19994, 63168/20107, 63523/20220, 63878/20333, 64233/20446, 64588/20559, 64943/20672, 65298/20785, 65653/20898, 66008/21011, 66363/21124, 66718/21237, 67073/21350, 67428/21463, 67783/21576, 68138/21689, 68493/21802, 68848/21915, 69203/22028, 69558/22141, 69913/22254, 70268/22367, 70623/22480, 70978/22593, 71333/22706, 71688/22819, 72043/22932, 72398/23045, 72753/23158, 73108/23271, 73463/23384, 73818/23497, 74173/23610, 74528/23723, 74883/23836, 75238/23949, 75593/24062, 75948/24175, 76303/24288, 76658/24401, 77013/24514, 77368/24627, 77723/24740, 78078/24853, 78433/24966, 78788/25079, 79143/25192, 79498/25305, 79853/25418, 80208/25531, 80563/25644, 80918/25757, 81273/25870, 81628/25983, 81983/26096, 82338/26209, 82693/26322, 83048/26435, 83403/26548, 83758/26661, 84113/26774, 84468/26887, 84823/27000, 85178/27113, 85533/27226, 85888/27339, 86243/27452, 86598/27565, 86953/27678, 87308/27791, 87663/27904, 88018/28017, 88373/28130, 88728/28243, 89083/28356, 89438/28469, 89793/28582, 90148/28695, 90503/28808, 90858/28921, 91213/29034, 91568/29147, 91923/29260, 92278/29373, 92633/29486, 92988/29599, 93343/29712, 93698/29825, 94053/29938, 94408/30051, 94763/30164, 95118/30277, 95473/30390, 95828/30503, 96183/30616, 96538/30729, 96893/30842, 97248/30955, 97603/31068, 97958/31181, 98313/31294, 98668/31407, 99023/31520, 99378/31633, 99733/31746, 100088/31859, 100443/31972, 100798/32085, 101153/32198, 101508/32311, 101863/32424, 102218/32537, 102573/32650, 102928/32763, 103283/32876, 103638/32989, 103993/33102, 104348/33215, 208341/66317, 312689/99532, 833719/265381, 1146408/364913, 3126535/995207, 4272943/1360120, 5419351/1725033, 42208400/13435351, 47627751/15160384, 53047102/16885417, 58466453/18610450, 63885804/20335483, 69305155/22060516, 74724506/23785549, 80143857/25510582, 245850922/78256779]

### Some thoughts

• There is a huge jump between 333/106 and 355/113, with the error improving from about 8*10^(-5) to 3*10^(-7). The next improvement doesn't occur until 52163/16604, which is still just about 3*10^(-7).

• 245850922.to_f / 78256779 == Math::PI returns true, which makes sense: Ruby floats, roughly speaking, only have precision to around fifteen decimal places. So to go further than 10_000_000, I really need to use something like http://ruby-doc.org/stdlib-1.9.3/libdoc/bigdecimal/rdoc/BigMath.html#method-i-PI

So by the right choice of denominator, we can get more accurate digits of pi than digits in the denominator. But from another perspective, we always need to use more digits in the fraction than we get accurate digits of pi.

With one exception: 355/113 yields the first seven digits of pi accurately, but only needs six digits to express.

What I'm saying is, you can build me a bridge using 355/113, and I'll probably drive across it. But if you have time to compute 245850922/78256779, that's pretty good too.

This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters. Learn more about bidirectional Unicode characters
 module PiApprox class PiApprox include Comparable attr_reader :to_f, :num, :den, :err def initialize(den) raise ArgumentError, 'Argument is not an integer' unless den.is_a? Integer @den = den @num = (Math::PI*@den).round @to_f = @num.to_f / @den @err = (@to_f - Math::PI).abs end def inspect "#{@num}/#{@den}" end alias_method :i, :inspect def <=>(other) @err <=> other.err end end def puts_improvements_up_to(max) (2..max).each_with_object([PiApprox.new(1)]) do |n, arr| puts "... checking #{n}" if Math.log10(n) == Math.log10(n).to_i comp = PiApprox.new n arr << comp if comp < arr[-1] end.each {|i| puts "#{i.i}, error=#{i.err}"} end end