|
|
A244644
|
|
Consider the method used by Archimedes to determine the value of Pi (A000796). This sequence denotes the number of iterations of his algorithm which would result in a difference of less than 1/10^n from that of Pi.
|
|
13
|
|
|
0, 1, 3, 5, 6, 8, 10, 11, 13, 15, 16, 18, 20, 21, 23, 25, 26, 28, 29, 31, 33, 34, 36, 38, 39, 41, 43, 44, 46, 48, 49, 51, 53, 54, 56, 58, 59, 61, 63, 64, 66, 68, 69, 71, 73, 74, 76, 78, 79, 81, 83, 84, 86, 88, 89, 91, 93, 94, 96, 98, 99, 101, 103, 104, 106, 108, 109, 111, 113, 114
(list;
graph;
refs;
listen;
history;
text;
internal format)
|
|
|
OFFSET
|
0,3
|
|
COMMENTS
|
It takes on average 5/3 iterations to yield another digit in the decimal expansion of Pi.
The side of a 96-gon inscribed in a unit circle is equal to sqrt(2-sqrt(2+sqrt(2+sqrt(2+sqrt(3))))). This is the size of one of the two polygons that Archimedes used to derive that 3 + 10/70 < Pi < 3 + 10/71.
In the Mathematica program, we started with an inscribed triangle and a circumscribed triangle of a unit circle and used decimal precision to just over a 1000 places.
The perimeter of the circumscribed 3*2^n-polygon exceeds Pi by more than the deficit of the perimeter of the inscribed 3*2^n-polygon. If we were to give twice the weight of the inscribed 3*2^n-polygon to that of the circumscribed 3*2^n-polygon, then the convergence would be twice as fast!
Archimedes's scheme: set upp(0) = 2*sqrt(3), low(0) = 3 (hexagons); upp(n+1) = 2*upp(n)*low(n)/(upp(n)+low(n)) (harmonic mean, i.e., 1/upp(n+1) = (1/upp(n) + 1/low(n))/2), low(n+1) = sqrt(upp(n+1)*low(n)) (geometric mean, i.e., log(low(n+1)) = (log(upp(n+1)) + log(low(n)))/2), for n >= 0. Invariant: low(n) < Pi < upp(n); variant function: upp(n)-low(n) tends to zero for n -> inf. The error of low(n) and upp(n) decreases by a factor of approximately 4 each iteration, i.e., approximately 2 bits are gained by each iteration.
From Archimedes's scheme, set r(n) = (2*low(n) + upp(n))/3, r(n) > Pi and the error decreases by a factor of approximately 16 for each iteration, i.e., approximately 4 bits are gained by each iteration. This is often called "Snellius acceleration".
Note that replacing "5/3" by "log(10)/log(4)" would be better in the first comment. (End)
|
|
REFERENCES
|
Petr Beckmann, A History of Pi, 5th Ed. Boulder, Colorado: The Golem Press (1982).
Jonathan Borwein and David Bailey, Mathematics by Experiment, Second Edition, A. K. Peters Ltd., Wellesley, Massachusetts 2008.
Jonathan Borwein & Keith Devlin, The Computer As Crucible, An Introduction To Experimental Mathematics, A. K. Peters, Ltd., Wellesley, MA, Chapter 7, 'Calculating [Pi]' pp. 71-79, 2009.
Eli Maor, The Pythagorean Theorem, Princeton Science Library, Table 4.1, page 55.
Daniel Zwillinger, Editor-in-Chief, CRC Standard Mathematical Tables and Formulae, 31st Edition, Chapman & Hall/CRC, Boca Raton, London, New York & Washington, D.C., 2003, §4.5 Polygons, page 324.
|
|
LINKS
|
Frits Beukers and Weia Reinboud, Snellius versneld, (text in English), NAW 5/3 no. 1, pp. 60-63 (2002).
|
|
FORMULA
|
Conjecture: There exists a c such that a(n) = floor(n*log(10)/log(4)+c); where c is in the range [0.08554,0.10264]. Critical values to narrow the range are believed to be at a(74), a(133), a(192), a(251), a(310), a(366), a(425), a(484). - A.H.M. Smeets, Jul 23 2018
|
|
EXAMPLE
|
Just averaging the initial two triangles (3.89711) would yield Pi to one place of accuracy, i.e., the single digit '3'. Therefore a(0) = 0.
The first iteration yields, as the perimeters of the two hexagons, 4*sqrt(3) and 6. Their average is ~ 3.2320508 which is within 1/10 of the true value of Pi. Therefore a(1) = 1.
a(3) = 5 since it takes 5 iterations of Archimedes's algorithm to drive the averaged value of the circumscribed 96-gon and the inscribed 96-gon to yield a value within 0.001 of the correct value of Pi.
a(4) = 6 since it takes 6 iterations of Archimedes's algorithm to drive the averaged value of the circumscribed 3*2^6-gon and the inscribed 3*2^6-gon to yield a value within 0.0001 of the correct value of Pi.
|
|
MATHEMATICA
|
a[n_] := a[n] = N[2 a[n - 1] b[n - 1]/(a[n - 1] + b[n - 1]), 2^10]; b[n_] := b[n] = N[ Sqrt[ b[n - 1] a[n]], 2^10]; a[-1] = 2Sqrt[27]; b[-1] = a[-1]/2; f[n_] := Block[{k = 0}, While[ 10^n*((a[k] + b[k])/4 -Pi) > 1, k++]; k]; Array[f, 70]
|
|
CROSSREFS
|
|
|
KEYWORD
|
nonn,base,easy
|
|
AUTHOR
|
|
|
STATUS
|
approved
|
|
|
|