Fractals/mcmullen

C programs by Curtis T McMullen

=How to run=

julia.tar

 * description

Steps:
 * run tcsh
 * make
 * copy ps2gif script to each subdirectory
 * add ./ to each run script
 * run

=Names=
 * "Eye" of elephant resting on internal angle p/q of main cardioid = principal Misiurewicz point

=Images=

=code=

Drawing dynamic external ray using inverse Boettcher map by Curtis McMullen


This method is based on C program by Curtis McMullen and its Pascal version by Matjaz Erat

There are 2 planes here :
 * w-plane ( or f0 plane )
 * z-plane ( dynamic plane of fc plane )

Method consist of 3 big steps :
 * compute some w-points of external ray of circle for angle $$\theta \,$$ and various radii (rasterisation)
 * $$w_0 = R * e^{2\pi i t} \,$$ where $$1 < R < \infty \,$$


 * map w-points to z-point using inverse Boettcher map $$\Psi_c = \Phi_{c}^{-1} \,$$
 * $$z_0 = \Psi_c (w_0) = f_c^{-n}(w_0^{2^n}) \,$$


 * draw z-points ( and connect them using segments ( line segment is a part of a line that is bounded by two distinct end points )

First and last steps are easy, but second is not so needs more explanation.

Rasterisation
For given external ray in $$f_0\,$$ plane each point of ray has :
 * constant value $$t\,$$ ( external angle in turns )
 * variable radius $$R\,$$

so $$w\,$$ points of ray  are parametrised by radius $$R\,$$ and can be computed using exponential form of complex numbers :


 * $$w = F_t(R) = R * e^{2*\Pi*i*t}\,$$

One can go along ray using linear scale : t:1/3; /* example value */ R_Max:4; R_Min:1.1; for R:R_Max step -0.5 thru R_Min do w:R*exp(2*%pi*%i*t); /* Maxima allows non-integer values in for statement */

It gives some w points with equal distance between them.

Another method is to use nonlinera scale. To do it we introduce floating point exponent $$r\,$$ such that :


 * $$R = 2^r\,$$

and
 * $$w = G_t(r) = 2^r * e^{2*\Pi*i*t}\,$$

To compute some w points of external ray in $$f_0\,$$ plane for angle $$t\,$$ use such Maxima code :

t:1/3; /* external angle in turns */ /* range for computing R ; as r tends to 0 R tends to 1 */ rMax:2; /* so Rmax=2^2=4 / rMin:0.1; /* rMin > 0  */ caution:0.93; /* positive number < 1 ; r:r*caution gives smaller r */ r:rMax; unless r<rMin do (  r:r*caution, /* new smaller r */   R:2^r, /* new smaller R */  w:R*exp(2*%pi*%i*t) /* new point w in f0 plane */ );

In this method distance between points is not equal but inversely proportional to distance to boundary of filled Julia set.

It is good because here ray has greater curvature so curve will be more smooth.

Mapping
Mapping points from $f_0\,$-plane to $f_c\,$-plane consist of 4 minor steps :
 * forward iteration in $$f_0\,$$ plane
 * $$w_n = w_0^{2^n} = R^{2^n}*e^{2\pi i 2^n t}\,$$

until $$w_n\,$$ is near infinity
 * $$abs(w_n) > R_{limit} \,$$


 * switching plane ( from $$f_0\,$$ to $$f_c\,$$)
 * $$z_n =w_n \,$$

( because here, near infinity : $$w = z \,$$)
 * backward iteration in $$f_c\,$$ plane the same ( $$n \,$$) number of iterations
 * last point $$ z_0 \,$$ is on our external ray

1,2 and 4 minor steps are easy. Third is not.
 * $$z_{n-1} = f_c^{-1} (z_n) = \sqrt{z_n - c}\,$$

Backward iteration uses square root of complex number. It is 2-valued functions so backward iteration gives binary tree.

One can't choose good path in such tree without extre informations. To solve it we will use 2 things :
 * equicontinuity of basin of attraction of infinity
 * conjugacy between $$f_0\,$$ and $$f_c\,$$ planes

Equicontinuity of basin of attraction of infinity
Basin of attraction of infinity ( complement of filled-in Julia set) contains all points which tends to infinity under forward iteration.

$$A_{f_c}(\infty) \ \overset{\underset{\mathrm{def}}{}}{=} \  \{ z \in  \mathbb{C}  : f^{(k)} _c (z)  \to  \infty\  as\  k \to \infty \}. $$

Infinity is superattracting fixed point and orbits of all points have similar behaviour. In other words orbits of 2 points are assumed to stay close if they are close at the beginning.

It is equicontinuity ( compare with normality).

In $$f_c\,$$ plane one can use forward orbit of previous point of ray for computing backward orbit of next point.

Detailed version of algorithm

 * compute first point of ray (start near infinity ang go toward Julia set )
 * $$w = 2^r * e^{2\Pi i t}\,$$ where $$r = r_{max}\,$$

here one can easily switch planes :
 * $$z = w\,$$

It is our first z-point of ray.
 * compute next z-point of ray
 * compute next w-point of ray for $$r = r * caution \,$$
 * compute forward iteration of 2 points : previous z-point and actual w-point. Save z-orbit and last w-point
 * switch planes and use last w-point as a starting point :$$z_n = w_{last}\,$$
 * backward iteration of new $$z_n\,$$ toward new $$z_0\,$$ using forward orbit of previous z point
 * $$z_0\,$$ is our next z point of our ray
 * and so on ( next points ) until $$ r < r_{Min}\,$$

Maxima CAS src code /* gives a list of z-points of external ray for angle t in turns and coefficient c */ GiveRay(t,c):= block( [r],  /* range for drawing  R=2^r ; as r tends to 0 R tends to 1 */  rMin:1E-20, /* 1E-4;  rMin > 0  ; if rMin=0 then program has infinity loop !!!!! */  rMax:2,    caution:0.9330329915368074, /* r:r*caution ; it gives smaller r */  /* upper limit for iteration */  R_max:300,  /* */  zz:[], /* array for z points of ray in fc plane */  /*  some w-points of external ray in f0 plane  */  r:rMax,  while 2^r=R_max) in f0 plane */  R:2^r,  w:rectform(ev(R*exp(2*%pi*%i*t))),  z:w, /* near infinity z=w */  zz:cons(z,zz),   unless r<rMin do  ( /* new smaller R */ r:r*caution, R:2^r, /* */  w:rectform(ev(R*exp(2*%pi*%i*t))), /* */  last_z:z, z:Psi_n(r,t,last_z,R_max), /* z=Psi_n(w) */ zz:cons(z,zz) ), return(zz) )$

colors
=References=