lisp or draw error

Hi,

Here is last version:
http://commons.wikimedia.org/wiki/Image:6furcation.gif

Regards

Adam



Adam Majewski pisze:
> Hi,
> I have changed code. It works better ( for biger iMax gets bind stack 
> overflow).
> Output images are slightly diffrent then  in :
> http://www.iec.csic.es/~miguel/Preprint4.ps
> I don't know why ( now ).
> 
> Regards
> 
> Adam
> 
> ================== code ====================================
> 
> 
> /* compute c for given radius, angle and l ( lambda=multiplier of fixed 
> point) */
> epsilon:0.001;
> radius:1+epsilon;
> angle:(3-sqrt(5))/2; /* in radians */
> l:radius * exp(%i*angle*2*%pi);
> c:l/2 - (l*l)/4;
> 
> trigsimp(%);
> rectform(%);
> c:float(%), numer;
> 
> 
> /* define function ( map) for dynamical system z(n+1)=f(zn,c)  */
> f(z,c):=z*z+c;
> 
> /* maximal number of iterations */
> iMax:200; /* to big couses bind stack overflow */
> EscapeRadius:10;
> 
> /* define z-plane ( dynamical ) */
> zxMin:-0.8;
> zxMax:0.2;
> zyMin:-0.2;
> zyMax:0.8;
> 
> /* resolution is proportional to number of details and time of drawing */
> iXmax:1000;
> iYmax:1000;
> 
> 
> /* compute critical point */
> zcr:rhs(solve(diff(f(z,c),z,1)));
> /* save critical point to 2 lists */
> xcr:makelist (realpart(zcr), i, 1, 1); /* list of re(z) */
> ycr:makelist (imagpart(zcr), i, 1, 1); /* list of im(z) */	
> 
> 
> /* ------------------- compute forward orbit of critical point ----------*/
> z:zcr; /* first point  */
> orbit:[z];
> for i:1 thru iMax step 1 do
> block
> (
>    z:f(z,c),
>    if abs(z)>EscapeRadius then return(i) else orbit:endcons(z,orbit)
>    );
> 
> 
> /*-------------- save orbit to draw it later on the screen 
> ----------------------------- */
> 
> /* save the z values to 2 lists */
> xx:makelist (realpart(zcr), i, 1, 1); /* list of re(z) */
> yy:makelist (imagpart(zcr), i, 1, 1); /* list of im(z) */
> for i:2 thru length(orbit) step 1 do
> block
> (
>    xx:cons(realpart(orbit[i]),xx),
>    yy:cons(imagpart(orbit[i]),yy)
> );		
> 
> 
> 
> 
> /* draw reversed orbit of beta  using draw package */
> load(draw);
> draw2d(
> 
>      terminal  = 'screen,
> 	pic_width  = iXmax,
>      pic_height = iYmax,
> 	yrange = [zyMin,zyMax],
>      xrange = [zxMin,zyMax],
> 	title= concat(" dynamical plane for f(z,c):=z*z+",string(c)),
> 	xlabel     = "Z.re ",
>      ylabel     = "Z.im",
>      point_type    = filled_circle,
> 	
>      points_joined = false,
> 	color		  =red,
> 	point_size    = 0.5,
> 	key = "critical orbit",
> 	points(xx,yy),
> 	key = "critical point",
> 	color		  =blue,
> 	points_joined = false,
> 	points(xcr,ycr)
>   );