🚀 go-pugleaf

In article <1993Oct14.183052.15240@ast.saic.com>, 
  hilljr@jupiter.saic.com (Jay R. Hill) writes:
|> 
|> LYNCH%XAVIER.BITNET@uga.cc.uga.edu (Chip Lynch) wrote on Wed, 13 
  Oct 1993 18:48:00 EST
|> 
|> |> The Feigenbaum number is hailed in the New Guiness book as the
|> |>         "Newest Mathematical Constant"
|> |> 
[...]
|> 
|> Here is a C++ program which illustrates the Feigenbaum Constant calculation
|> using the Mandelbrot Set iteration formula.  Most texts use the logistic
|> formula  x=r*x*(1-x).   This program shows the constant computed for two
|> period doubling cascades.  One starts with the period 1 cardioid and the
|> second starts with the period 3 cardioid.  At the end of the cascades are 
|> Feigenbaum points, for which I compute -1.40115518909205060052382678789... 
|> and -1.7798180758349913642620671...   Change the comment format and this 
|> program becomes C.  I prefer C++, a simple change of double to complex
|> here and there and you can explore the cascades for complex components.
|> When you do, you can try n=8, b=(-0.156520166833755062,-1.0322471089228318),
|> a=(-0.015608533223847222 -1.03664566584786215).

I am so sorry, I gave you the wrong period 8 root from my root file.  Use
for 'a' this one: a=(-.162867707067255135, 1.037313240810424765).   To make 
up for it, here is a complex math version of the program:

// C++ Compute Feigenbaum Universal Constant
// Complex version, (C) Jay R Hill, 1993

#include <stdlib.h>
#include <complex.h>
complex Newton(long n, complex c){
// find root of z=z^2 + c iterated n times
  complex z, z1, nc=c; double absz=1; long i,j;
  j=0;  // j limits iterations
  while(j<  7 && absz>(double)1e-13){ ++j;
    z=(complex)0;
    z1=(complex)0;
    for(i=0; i<n; ++i){
      z1=2*z1*z+1;        // first derivative
      z=z*z+nc;           // iterate
    }
    nc-=z/z1;         // Newton update
    absz=abs(z);      // for early out
  }
  return(nc);
}
void go(long n0, complex a, complex b){
// compute period doubling cascade
  double F=4.0; complex Tmp=a,NewC=a+(a-b)/F,OldC=b; long i,n=2*n0;
  for(i=0; i<10; ++i){
    NewC=Newton(n, NewC);
// make new estimate for Feigenbaum Universal Constant
    F=real((Tmp-OldC)/(NewC-Tmp));
    printf("%.16lf %.16lf %.16lf\n",real(NewC),imag(NewC),F);
    OldC=Tmp;
    Tmp=NewC;
// use estimate for Feigenbaum Universal Constant to find next
// period doubling cascade
    NewC+=(NewC-OldC)/F;
    n*=2;
  }
}
void main(void)
{
  complex a;  // center of period n component
  complex b;  // center of period n/2 component
  long n;     // period
  printf(
"      real(c)            imag(c)       F: (c2-c1)/(c3-c2)\n");
  n=2;
  b=(complex)0.0;
  a=(complex)-1.0;
  go(n,a,b);  // look at sequence starting at period 1 cardioid
  printf("\n");
  n=6;
  b=(complex)-1.754877666246692760049508896358;
  a=(complex)-1.772892903381623799434128230874;
  go(n,a,b);  // look at sequence starting at period 3 cardioid
  printf("\n");
  n=8;
  b=(complex)-1.940799806529484752232090979655;
  a=(complex)-1.941782090318198160140455136524;
  go(n,a,b);  // look at sequence starting at real period 4 cardioid
  printf("\n");
  b=complex(-.156520166833755061798986650788, 1.03224710892283180167160153243);
  a=complex(-.162867707067255135276520417258, 1.03731324081042476523420269342);
  go(n,a,b);  // look at sequence starting at complex period 4 cardioid
}
/*
      real(c)            imag(c)       F: (c2-c1)/(c3-c2)
-1.3107026413368328 0.0000000000000000 3.2185114220380888
-1.3815474844320612 0.0000000000000000 4.3856775985683472
-1.3969453597045607 0.0000000000000000 4.6009492765379756
-1.4002530812147831 0.0000000000000000 4.6551304953917194
-1.4009619629448422 0.0000000000000000 4.6661119478231585
-1.4011138049397749 0.0000000000000000 4.6685485815265810
-1.4011463258269472 0.0000000000000000 4.6690606602496461
-1.4011532908499256 0.0000000000000000 4.6691715552590010
-1.4011547825466131 0.0000000000000000 4.6691951766774649
-1.4011551020224637 0.0000000000000000 4.6692001427014445

-1.7782668211110817 0.0000000000000000 3.3523470290916921
-1.7794823875518131 0.0000000000000000 4.4209164957074565
-1.7797460275941450 0.0000000000000000 4.6107049216787335
-1.7798026381018630 0.0000000000000000 4.6570866957290749
-1.7798147692180581 0.0000000000000000 4.6665539104053879
-1.7798173676438036 0.0000000000000000 4.6686407015227296
-1.7798179241614378 0.0000000000000000 4.6690807007454609
-1.7798180433511361 0.0000000000000000 4.6691756249117882
-1.7798180688779419 0.0000000000000000 4.6691975235440442
-1.7798180743450029 0.0000000000000000 4.6692008116173804

-1.9420855403320916 0.0000000000000000 3.2370530358863174
-1.9421546694898462 0.0000000000000000 4.3896095909440485
-1.9421696903322636 0.0000000000000000 4.6022157634979770
-1.9421729168962865 0.0000000000000000 4.6553678498442297
-1.9421736083767791 0.0000000000000000 4.6661678198897549
-1.9421737564910764 0.0000000000000000 4.6685600605824327
-1.9421737882135650 0.0000000000000000 4.6690629783611408
-1.9421737950075979 0.0000000000000000 4.6691691284062635
-1.9421737964626664 0.0000000000000000 4.6692185588198223
-1.9421737967743034 0.0000000000000000 4.6691143778927930

-0.1647250612100108 1.0389057856466151 3.3173727898425440
-0.1651440512497371 1.0392690840074308 4.4117483369316961
-0.1652348396826442 1.0393480657843397 4.6084499102036540
-0.1652543311392327 1.0393650329859927 4.6566119207233827
-0.1652585078135890 1.0393686692822659 4.6664502295581052
-0.1652594024296227 1.0393694481759046 4.6686186913402787
-0.1652595940336001 1.0393696149961775 4.6690759451025707
-0.1652596350695145 1.0393696507241976 4.6691752814875516
-0.1652596438581618 1.0393696583760599 4.6691944637674432
-0.1652596457404239 1.0393696600148612 4.6691890367931270
*/

Warmly,
Jay "Not tonight honey, I'm in a complex period doubling cascade' Hill
-- 
{ hilljr@jupiter.saic.com } begin writeln(3*ln(640320)/sqrt(163):17:15) end.
void main(){double sqrt(), y=1/sqrt(2.), a=.5, m=1,z; int n=0;
for(;m*=2,z=sqrt(1-y*y),y=(1-z)/(1+z),a=a*(1+y)*(1+y)-m*y,n<4;n++);
printf("%17.15lf\n",1/a);}