[seqfan] Extension of A097486, Calculating Pi in the Mandelbrot Set
Robert Munafo
mrob27 at gmail.com
Fri Jan 8 04:52:58 CET 2010
I will (when able) be submitting a program, more terms, and web link for
A097486.
I discovered that the last term in the current entry, A(7)=31415928, is off
by one. When using the higher precision of PARI/GP, and also when using my
C++ quad-precision library (which is a whole lot faster) I get 31415927. I
suspect it is the result of round-off error in normal double-precision
floating-point (in IEEE double precision, this algorithm loops forever if
you try to calculate A(8)!)
Can anyone with Mathematica or Maple or Matlab verify? The algorithm is
pretty simple as seen in the PARI/GP program below. "I" is the square root
of -1. z, c, and the operation "z^2+c" are all complex quantities. "abs(z)"
is the modulus of z.
Program: (PARI/GP):
A097486(n) = local(a,c,z); c = 0.1^n*I - 0.75; z=c; a=0; while(abs(z)<4.0,
{ z=z^2+c; a=a+1 } ); a
I am also able to compute 3 more terms in about an hour.
First 11 terms: 3, 33, 315, 3143, 31417, 314160, 3141593, 31415927,
314159266, 3141592655, 31415926537
Link: http://www.mrob.com/pub/muency/seahorsevalley.html
--
Robert Munafo -- mrob.com
More information about the SeqFan
mailing list