File:Dynamical plane with branched periodic external ray 0 for map f(z) = z*z + 0.35.svg

Summary[edit]

DescriptionDynamical plane with branched periodic external ray 0 for map f(z) = z*z + 0.35.svg	English: dynamical plane with branched periodic external ray 0 for map f(z) = z*z + 0.35. " Julia set for a polynomial ... which belongs to the ray R()(M) in parameter space. The right hand limit ray Ro bounces off infinitely many pre-critical points as it spirals in to the upper fixed point. Both fixed points have rotation number zero." [1] The countable family of R-analytic curves on which we cannot define the function “external argument”.[2]
Date	15 November 2020
Source	Own work
Author	Adam majewski
Other versions	BD BD DLD SAC LSM Papers Fixed points of polynomial maps. Part II. Fixed point portraits Lisa R. Goldberg; John Milnor Annales scientifiques de l'École Normale Supérieure (1993) Volume: 26, Issue: 1, page 51-98 ISSN: 0012-9593, see Fig 14 page 91 Atela, P. (1992). Bifurcations of dynamic rays in complex polynomials of degree two. Ergodic Theory and Dynamical Systems, 12(3), 401-423. doi:10.1017/S0143385700006854 Schleicher, D., “Rational parameter rays of the Mandelbrot set”, arXiv Mathematics e-prints, 1997. Exploring the Mandelbrot set. The Orsay Notes. Adrien Douady John H. Hubbard, page 64 fig 3 : " The countable family of R-analytic curves on which we cannot define the function “external argument”. Video: Binary decomposition around Julia set : Zn+1 = Zn^2 + c, c goes from .25 to 12
SVG development InfoField	The source code of this large SVG is valid.   This oversized plot was created with Gnuplot.    This plot uses embedded text that can be easily translated using a text editor.

Licensing[edit]

I, the copyright holder of this work, hereby publish it under the following license:

This file is licensed under the Creative Commons Attribution-Share Alike 4.0 International license.

You are free:

• to share – to copy, distribute and transmit the work
• to remix – to adapt the work

Under the following conditions:

• attribution – You must give appropriate credit, provide a link to the license, and indicate if changes were made. You may do so in any reasonable manner, but not in any way that suggests the licensor endorses you or your use.
• share alike – If you remix, transform, or build upon the material, you must distribute your contributions under the same or compatible license as the original.

https://creativecommons.org/licenses/by-sa/4.0CC BY-SA 4.0 Creative Commons Attribution-Share Alike 4.0 truetrue

Software[edit]

• Maxima CAS
  • draw package by Mario Rodríguez Riotorto^[3]
• gnuplot
• svgoptimizer: from 15.5 MB to 8.9 MB

Algorithm[edit]

Ray bifurcates at critical point and at it's preimages

Ray consist of segments ( curves between bifurcation points):

• r0: external ray 0 goes from (+) infinity along critical orbit towards critical point z=0 ( horizontal segment). Here bifurcates.
• segment r1 has two subsegments : r1a and r1b. These are two vertical segments from critical point z=0 towards it's two preimages : a(z) = f^-1(z) and b(z) = -a(z). So it is: [a(z), -a]
• each (sub)segment of r1 has 2 preimages
  • a(r1a) and b(r1a)
  • a(r1b) and b(r1b)
• repeat and infinity ( binary tree)

 Ray = r0 + r1 + r2 + ... + rn

Problems[edit]

inverse function[edit]

${\displaystyle f^{-1}(z)={\sqrt {z-c}}}$

the standard complex square-root ( csqrt) takes values in the right halfplane, but this is rotated by the squareroot. So some points from one list goes to preimage of anothe lists. Onbe can not use gnuplot option "points joined= true".

The soultion

• use lists with more points ( takes time) and gnuplot option "points joined= false"
• compute only preimages of one r1 subsegment and use symetry. It is 2-4 times faster and solves the problem of csqrt

length of the list[edit]

rn is a lis, but it consist of sublist ( list of lists). Each sublist rn[k] is a list of points describint segment (one of n=th preimage's of segment r0);

Notation

• n is a level of binary tree of preimages
• % denotes ratio of (sublists with size > PixelSize) / ( all sublists of rn)
• s denotes size ratio = size/PixelSize, where size = cabs( first - last)
• segments = number of segments in the n-th list = 2^n
• points = length(r[n][k]) = kMax + 2

List r11 has 4096 sublists ( segments) but only 77 of them are wider then Pixel size = 0.0024

r0 list: 1 / 1 =  1.0 segments has size > MinSize. MaxSize =  247 pixels 
r1 list: 2 / 2 =  1.0 segments has size > MinSize. MaxSize =  179 pixels 
r2 list: 4 / 4 =  1.0 segments has size > MinSize. MaxSize =  108 pixels 
r3 list: 8 / 8 =  1.0 segments has size > MinSize. MaxSize =  66 pixels 
r4 list: 16 / 16 =  1.0 segments has size > MinSize. MaxSize =  43 pixels 
r5 list: 32 / 32 =  1.0 segments has size > MinSize. MaxSize =  31 pixels 
r6 list: 64 / 64 =  1.0 segments has size > MinSize. MaxSize =  24 pixels 
r7 list: 128 / 128 =  1.0 segments has size > MinSize. MaxSize =  20 pixels 
r8 list: 220 / 256 =  0.859375 segments has size > MinSize. MaxSize =  19 pixels 
r9 list: 148 / 512 =  0.2890625 segments has size > MinSize. MaxSize =  18 pixels 
r10 list: 224 / 2048 =  0.109375 segments has size > MinSize. MaxSize =  17 pixels 
r11 list: 200 / 4096 =  0.048828125 segments has size > MinSize. MaxSize =  15 pixels 
r12 list: 196 / 8192 =  0.02392578125 segments has size > MinSize. MaxSize =  13 pixels 

Maxima CAS src code[edit]

/*

Batch file for Maxima CAS
save as a s.mac
run maxima : 
 maxima
and then : 
batch("s.mac");

or from bash : 
maxima -b f.mac >f.txt





(%i39) ray:GiveRay(depthMax,iMax,HalfPixelSize2)
depth d =  1 length(rd) = 2 
depth d =  2 length(rd) = 4 
depth d =  3 length(rd) = 8 
depth d =  4 length(rd) = 16 
depth d =  5 length(rd) = 32 
depth d =  6 length(rd) = 64 
depth d =  7 length(rd) = 128 
depth d =  8 length(rd) = 256 
depth d =  9 length(rd) = 512 
depth d =  10 length(rd) = 1024 
depth d =  11 length(rd) = 2048 
depth d =  12 length(rd) = 4096 
depth d =  13 length(rd) = 8192 
depth d =  14 length(rd) = 16384 
depthMax =  14 length(rA) = 32767 

==================
pth d =  1 length(rd) = 2 
depth d =  2 length(rd) = 4 
depth d =  3 length(rd) = 8 
depth d =  4 length(rd) = 16 
depth d =  5 length(rd) = 32 
depth d =  6 length(rd) = 64 
depth d =  7 length(rd) = 128 
depth d =  8 length(rd) = 256 
depth d =  9 length(rd) = 512 
depth d =  10 length(rd) = 832 
depth d =  11 length(rd) = 648 
depth d =  12 length(rd) = 496 
depth d =  13 length(rd) = 448 
depth d =  14 length(rd) = 440 
depth d =  15 length(rd) = 444 
depth d =  16 length(rd) = 444 
depth d =  17 length(rd) = 408 
depth d =  18 length(rd) = 312 
depth d =  19 length(rd) = 232 
depth d =  20 length(rd) = 156 
depthMax =  20 length(rA) = 5883 
           
           depth max = 20  depthMax_pre =  6  number of points in the ray =  11767 



*/




kill(all)$
remvalue(all)$
display2d:false$



/* ---------- functions ---------------------------------------------------- */

/* http://en.wikipedia.org/wiki/Complex_quadratic_polynomial  */

/* Forward iteration */
f(z):=float(rectform(z*z+c))$


/* 

backward iteration of complex quadratic polynomial 
when Im(z-c)  changes sign ( from quadrant (-,-) to (-,+)
then preimaga chages quadrant 

w = sqrt(z-c)

https://www.geogebra.org/m/YPGT79Zq

*/
  
b(z):=float(rectform(-sqrt(z-c)))$
a(z):=float(rectform(sqrt(z-c)))$


/* find fixed point alfa of function f(z,c)   */
GiveFixed_p(c):= float(rectform((1+sqrt(1-4*c))/2))$
GiveFixed_m(c):= float(rectform((1-sqrt(1-4*c))/2))$

/* 
converts complex number z = x*y*%i 
to the list in a draw format:  
[x,y] 
*/
d(z):=[float(realpart(z)), float(imagpart(z))]$

/* give Draw List from one point*/
dl(z):=[d(z)]$

ToPoints(myList):='points(map(d,myList))$


/* gives an orbit of z0 under fc where iMax is the length of the orbit */
 GiveForwardOrbit(z0,c, iMax):= block
(
  [i,z, zz],
   
   z:z0,
   zz :[z],
   i:1,
   while ( i<iMax ) 
    do
    ( 
      z:f(z),
      zz : endcons(z, zz),
      i:i+1
    ),
    
   return(zz)
)$


/* 
gives an orbit of z0 under fc where iMax is the length of the orbit 
note that backward orbit is a binary tree 
so this procedure gives only some branches from 2^iMax branches 

*/
 GiveBackwardOrbit(z0,c, iMax):= block
(
  [i,z,  orbit],
   
   
   orbit :[],
   z : z0,
   i:1,
   while ( i<iMax ) 
    do
    ( 
      z:a(z),
      orbit : cons( z, orbit),
      orbit : cons(conjugate(z), orbit), 
      orbit : cons( -z, orbit),
      orbit : cons(conjugate(-z), orbit), 
      
      
      i:i+1
    ),
    
    
    
   
   
   return(orbit)
)$






/* 
Line segment =  part of a line that is bounded by two end points 
input : 
- 2 endpoints : z1, z2 
- imax : number of points between endpoints
output: 
- list of (iMax+1) complex points


*/

GiveLineSegment(z1,z2,iMax):=block
(
 [l,z, dz, t , tmin, tmax,dt],
 dz : z2-z1,
 dt : 1/iMax, 
 tmin: dt, 
 tmax: 1,
 z:z1,
 l:[z1],
 for t:tmin thru tmax step dt 
 
  
    do
    (
     z : z1+dz*t,
     l : cons(z,l)
    ),
return(l)
)$


Give_r0(z1,z2,iMax):=[GiveLineSegment(z1,z2,iMax)]$


/* cabs2 is a fast version of cabs because it does not use sqrt */
cabs2(z) := realpart(z)*realpart(z) +imagpart(z)*imagpart(z)$



/* 
size of the segment 
= distance between fist  and last points of the list 

output = size : r*MinSize 
*/
CheckSegmentSize(L, MinSize ):=block(
	[z1,z2,d, r],
	z1:first(L),
	z2:last(L),
	d:cabs2(z1-z2),
	r: d/MinSize,
	
	return(r)



)$



CheckList(LL, name, MinSize):=block(
	[i,r,l, rMax],
	i:0,
	rMax:0,
	l : length(LL), 
	for L in LL do (
		r:CheckSegmentSize(L , MinSize),
		if (r>1.0) then i:i+1,
		if (r>rMax) then rMax:r
				
		),
	print(name ," list:", i, "/", l ,"= " , float(i/l) , "segments has size > MinSize. MaxSize = ", ?round(rMax) ,"pixels"),
	
	return([i/l, rMax]) 




)$











/* 

- read input z
- compute 2 preimages a and b=-a
- output list : [a,b]
 

Give2Preimages(z):= [a(z), b(z)]$


*/	
	
Give2Preimages(z):= block(
	[p],
	p: a(z),
	return([p,-p]) 
	
	)$


/*

input: list zz
output:  list of 2 sublists  

*/


Give2PreimagesOfList(zz, MinSize):=block(
	[La,Lb, LL],
	La:[],
	Lb:[],
	LL:[],
	
	
	
	for z in zz do (
		
	
	
		tmp : Give2Preimages(z),
		La:cons(tmp[1], La),
		Lb:cons(tmp[2], Lb)
				
	),
	/*
	optimisation : check the list ( r = segment size) 
	if (r<1) then list ={{Own based}}ly a one point
	if (r<2) then list = only 2 points
	
	*/
	 
	r: CheckSegmentSize(La, MinSize ),
	if ( r< 1.0) 
		then La : [first(La)]
		else (if ( r<2.0) then La : [first(La), last(La)]),
	
	r: CheckSegmentSize(Lb, MinSize ),
	if ( r< 1.0) 
		then Lb : [first(Lb)]
		else (if ( r<2.0) then Lb : [first(Lb), last(Lb)]),
	
	
	LL:[La , Lb],
	return(LL)



)$



	
/* 

regular complete binary tree depth 

- read input A
- compute full tree of preimages depth d
( number = sum(2^(d-1))
- output: unsorted list zz of complex points
 
lists 
 zz is output list 
 t1 input = previous level of the tree 
 t2 next level 
 temp  temporary list 
 
 length of precritical orbit :    4 094  for depth :  11 
 length of precritical orbit :  262 142  for depth :  17 
 
 if length is to big then save time is too long



*/

GiveBTree(z, depthMax):=block(

	[ zz, t1, t2, tmp],
	zz :[], 
	t1 :[z], 
	t2 :[], 
	tmp  :[], 
		
	for depth : 1 thru depthMax step 1 do (
		
		
		for z in t1 do (
			tmp : Give2Preimages(z),
			zz : append(tmp, zz),
			t2 : append(tmp, t2)
			
		),
		
		tmp:[],
		t1:t2,
		t2:[]
		
		
	),
	
	zz:flatten(zz), 
	return (zz)
)$




/* 

rA is a part of the ray
binary not complete tree
LL = list of lists
segments of the rA :  r0+r1+r2 + ... rd 
each list rd is a list of lists ( subsegments)


*/

Give_rA(depthMax, iMax, distanceMin):= block(

	[r0, rIn, rOut, rT, rA, d],
	
	
	rA : [], /* global output list */
	rIn : [], /* local input */
	rOut : [], /* local output */
	rT : [], /* local temporary list */
	
	/* r0  : initial part of binary tree */
	r0:[],
	r0: Give_r0(zcr, a(zcr) , iMax), /* compute  vertical segment decribed by iMax points*/
	/* rn : LL */
	for L in r0 do rA:cons(L, rA), /* save to the putput list */
	
	rIn : r0,
	
	
	/* binary tree of preimages  to max depth = dMax  */
	for d:1 thru depthMax step 1 do (
	
		
		
		/* compute 2 preimage of rIn */
		for L in rIn do(  
			if (length(L)>1) /* only if list L is not a single point */
				then ( 
	
					rT: Give2PreimagesOfList(L, distanceMin),
					rOut:cons(rT[1],rOut),
					rOut:cons(rT[2],rOut)
					)
			/* 	else print("segment to small for depth = ", d ) */
			), 
		for L in rOut do rA:cons(L,rA), /* copy preimages to global output */	
		
		print("depth d = ", d, "length(rd) =" , length(rOut)),
		
		/* swap */
		rIn: rOut,
		rOut:[],
		rT:[]
		
		
	
	
	
	
	), /* for d:1 ... */ 
	
	
	print("depthMax = ", depthMax, "length(rA) =" , length(rA)),
		
	
	

	return(rA)




)$



/*

ray is a list of 3 lists
ray : [rA, rB, rH]
each list consist of lists of points

ray = LLL = L(l(L))
list L = list ocf complex points  = zz

rH is a horizontal segment
rA is a binary tree
rB : is computed from rA : conjugate (rA) = symmetry


ray is a binary tree , but not complete


*/

GiveRay(MaxDepth, iMax, distanceMin):=block(

	[ray, rA, rH],
	
	ray:[],
	rA:[],
	rH:[],
	
	/* compute parts of the ray */
	rH : [[zcr,half]], /* horizontal part of the ray is a streight line along criticla orbit so it can be described by only 2 points */
	rA : Give_rA(MaxDepth, iMax, distanceMin), 
	
	
	
	/* convert and save to output list */
	for L in rH do (
		L:map(d, L), /* draw list */
 		L: points(L), /* draw command */
		ray:cons(L,ray) 
	),
	
	for L in rA do( 
		L:map(d, L), 
	 	L: points(L),
	 	ray:cons(L, ray)
	 	
	 	),

	/* symmetric part of the ray rB: conjugate(rA */
	for L in rA do( 
		L:map(conjugate,L),
		L:map(d, L), 
	 	L: points(L),
	 	ray:cons(L, ray)
	 	),
	 	
	 	
	 return(ray)

)$



GivePrecritical(zcr, depth_max):=block(
	[zz],
	zz:  GiveBTree(zcr, depth_max), /* trick : use smaller depth then for ray because GiveBTree has no optimisation */
	zz : append(zz, GiveBackwardOrbit(zcr,c,20)), /* trick  : add some important braches of binary tree */
	zz: ToPoints(zz),
	return(zz)

)$
   

compile(all)$



/* ===================== const ===================================================== */


half :1.2$
yMax: half$
yMin: - half$
iSide: 1000$
PixelSize: (yMax - yMin)/iSide$ 




HalfPixelSize : PixelSize/2$  /* 2 points are the same if the distance between them is smaller them  half pixel size */
HalfPixelSize2: HalfPixelSize * HalfPixelSize$ /* for cabs2 */


/* 
   integer constant values 
 proportional to the : 
 * quality of the image ( number of the detailes, but also precision of curves )
 * size of the svg file ( but not the png file )
 * time of computation 
*/


iMax: 10; /* number of points in the line segment   */
depthMax: 20$ /* depth of ray binary tree = number of preimages  */


/* precritical settings */
point_size_pre: 0.4$
color_pre: red$
key_pre:  "precritical points"$
depthMax_pre:15; /* number of the orbit points: forward and backward  */

/* parameter of the function : fc(z) = z^2 + c */
c:0.35;
zcr : 0; /* critical point */


 

/* ============================ computations =============================== */

zfp : GiveFixed_p(c) $ /* fixed point */
zfm : GiveFixed_m(c) $ /* fixed point */


/* precritical orbit = backward orbit of the critical point */

precritical: GivePrecritical(zcr, depthMax_pre)$

 	
ray: GiveRay(depthMax, iMax, HalfPixelSize2)$


 




/* draw , also some computations */

path:"~/Dokumenty/branched_ray/maxima/n12/"$ /*  pwd, if empty then file is in a home dir , path should end with "/" */
fileName: sconcat(path,  string(depthMax), "_", string(depthMax_pre),"_",  string(point_size_pre),"_", color_pre)$
/* draw it using draw package by */

 load(draw); 
/* if graphic  file is empty (= 0 bytes) then run draw2d command again */

 draw2d(
  user_preamble="set key top right; unset mouse",
  terminal  = 'png,
  file_name = fileName,
  
  title= " dynamical plane : branched periodic external ray 0/1 for map f(z) = z*z + 0.35 ",
  
  dimensions = [iSide, iSide],
 
  yrange = [-half, half],
  xrange = [-half, half],
 
  
  xlabel     = "zx ",
  ylabel     = "zy",
  
  
  
  point_type = filled_circle,
  
  
  
  /*  big points */ 
   points_joined = false,
  point_size    =  1.1,
  key= "fixed points",
  color = blue,
  points(dl(zfp)),
  key = "",
  points(dl(zfm)),
   
  key= "critical point",
  color = black,
  points(dl(zcr)),
    
 

 
  key = key_pre,
  point_size    =  point_size_pre,
  points_joined = false,
  color         = color_pre,
  precritical,
 
  
  
  points_joined = true,
  point_size    =  0.1,
  key="",
  color=black,
  ray 
  
 )$

draw2d(
  user_preamble="set key top right; unset mouse",
  terminal  = 'svg,
  file_name = fileName,
  
  title= " dynamical plane : branched periodic external ray 0/1 for map f(z) = z*z + 0.35 ",
  dimensions = [iSide, iSide],
 
  yrange = [-half, half],
  xrange = [-half, half],
 
  
  xlabel     = "zx ",
  ylabel     = "zy",
  
  
  
  point_type = filled_circle,
  
  
  
  /*  big points */ 
   points_joined = false,
  point_size    =  1.1,
  key= "fixed points",
  color = blue,
  points(dl(zfp)),
  key = "",
  points(dl(zfm)),
   
  key= "critical point",
  color = black,
  points(dl(zcr)),
  

 
  key =  key_pre, 
  point_size    =  point_size_pre,
  points_joined = false,
  color         = color_pre,
  precritical,
 
  
  
  points_joined = true,
  point_size    =  0.1,
  key="",
  color=black,
  ray 
    
 
 )$





  
print("depth max =", depthMax, " depthMax_pre = ", depthMax_pre, " number of points in the ray = ", length(ray))$
print("file", fileName, " saved to",path)$
print("file name: depthMax_depthMax_pre_ point_size_pre.png")$
print("Pixel size : ", PixelSize)$
/* print("length of precritical orbit : ", l_precritical, " for depth : ", depthMax)$ */


	

	

/*
CheckList(r0, "r0", HalfPixelSize)$
CheckList(r1, "r1", HalfPixelSize)$
CheckList(r2, "r2", HalfPixelSize)$
CheckList(r3, "r3", HalfPixelSize)$
CheckList(r4, "r4", HalfPixelSize)$
CheckList(r5, "r5", HalfPixelSize)$
CheckList(r6, "r6", HalfPixelSize)$

CheckList(r7, "r7", HalfPixelSize)$
CheckList(r8, "r8", HalfPixelSize)$
CheckList(r9, "r9", HalfPixelSize)$
CheckList(r10, "r10", HalfPixelSize)$

CheckList(r11, "r11", HalfPixelSize)$
CheckList(r12, "r12", HalfPixelSize)$
*/

References[edit]

1. ↑ Fixed points of polynomial maps. Part II. Fixed point portraits Lisa R. Goldberg; John Milnor Annales scientifiques de l'École Normale Supérieure (1993) Volume: 26, Issue: 1, page 51-98 ISSN: 0012-9593, see Fig 14 page 91
2. ↑ Exploring the Mandelbrot set. The Orsay Notes. by Adrien Douady John H. Hubbard. page 64 Fig 3
3. ↑ draw package ( Maxima-Gnuplot interface) by Mario Rodríguez Riotorto