Fractals/Iterations in the complex plane/fprays

Parabolic fixed point of period p is a landing point of p dynamic external rays. These rays divide neighborhood into curvilinear sectors.

On the main cardioid edit

Maxima CAS code :


kill(all);
remvalue(all);

DoublingMap(r):=
block([d,n],
 n:ratnumer(r),
 d:ratdenom(r),
 mod(2*n,d)/d)$

GivePeriod (r):=
block([rNew, rOld, period, pMax, p],
      pMax:100,
      period:0,
       
      p:1, 
      rNew:DoublingMap(r),
      while ((p<pMax) and notequal(rNew,r)) do
        (rOld:rNew,
         rNew:DoublingMap(rOld),
         p:p+1
        ),
      if equal(rNew,r) then period:p,
      period
);

/* f(z) is used as a global function
   I do not know how to put it as a argument */

GiveOrbit(r0,OrbitLength):=
block(
 [r,Orbit],
 r:r0,
 Orbit:[r], 
 for i:1 thru OrbitLength step 1 do
        ( r:DoublingMap(r),
          Orbit:endcons(r,Orbit)),
         
 return(sort(Orbit))

)$

compile(all);

R: 4985538889/17179869183;
p: GivePeriod(R);

orbit:GiveOrbit(R, p);

/* angles around critical point */
e1:first(orbit);
e2:last(orbit);

1/p edit

Video showing exernal dynamic rays landin on parabolic fixed point for internal angles 1/p of main cardioid

In the elephant valley[1][2] ( from parameter plane ) it is easy to find rays landing on the roots and dynamic external rays that land on the parabolic fixed point z.

  • first choose internal angle (= combinatorial rotation number) : 1/p
  • compute pair of parameter rays, angles of the wake :  
  • compute list of p external angles of dynamic rays :  
internal angle of main cardioid parameter c = root point parameter rays parabolic fixed point z dynamic rays
0/1 (0/1; 1/1) (1/1)
1/2 (1/3; 2/3) (1/3; 2/3)
1/3 (1/7; 2/7) (1/7; 2/7; 4/7)
1/4 (1/15; 2/15) (1/15; 2/15; 4/15; 8/15)
1/5 (1/31; 2/31) (1/31; 2/31; 4/31; 8/31; 16/31)
1/6 (1/63; 2/63) (1/63; 2/63; 4/63, 8/63, 16/33; 32/63)
1/7 (1/127; 2/127) (1/127; 2/127; 4/127, 8/127, 16/127; 32/127; 64/127)
1/p    

Note that :

  • internal ray 0/1 = 1/1
  • internal angle 1/p means that ray goes from period 1 component ( parent period = 1) to period p component ( child period = p)
  • as child period grows computations are harder
  • exponential growth[3] of  . One can easly create a numeric value that is too large to be represented within the available storage space ( integer overflow[4] ). For example   is to big for short ( 16 bit ) and long ( 32 bit) integer.

n/p edit

It is not so simple, as in 1/p case, to compute orbit portrait[7] of parabolic fixed point.

Algorithm :

  • choose child period p
  • compute internal angle ( rational number) = n/p ( where n<p and n/p is ... ( to do ))
  • compute angle of the wake
  • switch to dynamic plane : use one angle from pair of parameter rays ( rays with the same angles land on the parabolic fixed point) to compute orbit portrait of parabolic fixed point


Methods for computing anglew of the wake:

  • Combinatorial algorithm = Devaney's method
  • step method:
  • compute denominator of external angle =  
  • find parameter rays that land on the root point which is on the boundary of main cardioid :
    • compute all pairs for periods 1-p
    • remove pairs which land not on the main cardioid ( inside < 1/3; 2/3 > wake )
    • compute pairs of external angles which are not inside pairs of lower periods (see image on the right )
    • choose n-th pair of angles which land on the root point

Child period 5 edit

 
Parameter external rays for period 5. There are only 4 red arcs which are not inside black arcs .

From all 15 period five components only 4 components are directly connected to the main cardioid [8]

internal angle of main cardioid parameter c = root point parameter rays parabolic fixed point z dynamic rays
1/5 (1/31; 2/31) (1/31; 2/31; 4/31; 8/31; 16/31)
2/5 -0.504+0.568 i (9/31,10/31) (5/31 , 10/31 , 20/31 , 9/31 , 18/31)
3/5 (21/31,22/31) (11/31 , 22/31 , 13/31 , 26/31 , 21/31)
4/5 (29/31,30/31) (15/31 , 30/31 , 29/31 , 27/31 , 23/31)

Child period 7 edit

 
Parameter external rays for period 1-7. There are only 6 red arcs which are not inside black arcs .

From all 63 period seven components only 6 components are directly connected to the main cardioid [9]

internal angle of main cardioid parameter c = root point parameter rays parabolic fixed point z dynamic rays
1/7 (1/127; 2/127) (1/127; 2/127; 4/127; 8/127; 16/127, 32/127, 64/127)
2/7 (17/127; 18/127) (17/127; 34/127; 68/127; 9/127; 18/127, 36/127, 72/127)
3/7 (42/127; 43/127) (42/127, 84/127, 82/127, 37/127, 74/127, 21/127, 42/127)
4/7 (84/127; 85/127) (84/127, 41/127, 82/127, 37/127, 74/127, 21/127, 42/127)
5/7 (109/127; 110/127) (109/127, 91/127, 55/127, 110/127, 93/127, 59/127, 118/127)
6/7 (125/127; 126/127) (125/127, 123/127, 119/127, 111/127, 95/127, 63/127, 126/127)


Child period 11 edit

internal angle of main cardioid parameter c = root point parameter rays parabolic fixed point z dynamic rays
1/11
5/11 -0.6900598700150440+0.2760264827846140i (681/2047, 682/2047) -0.4797464868072486+0.1408662784207147i (341, 597, 661, 677, 681, 682, 1194, 1322, 1354, 1362, 1364)/2047

Maxima CAS code using doubling map:

(%i1) m(n,d):=mod(2*n,d)/d $

(%i2) m(681,2047);
                                     1362
(%o2)                                ----
                                     2047
(%i3) m(1362,2047);
                                     677
(%o3)                                ----
                                     2047
(%i4) m(677,2047);
                                     1354
(%o4)                                ----
                                     2047
(%i5) m(1354,2047);
                                     661
(%o5)                                ----
                                     2047
(%i6) m(661,2047);
                                     1322
(%o6)                                ----
                                     2047
(%i7) m(1322,2047);
                                     597
(%o7)                                ----
                                     2047
(%i8) m(597,2047);
                                     1194
(%o8)                                ----
                                     2047
(%i9) m(1194,2047);
                                     341
(%o9)                                ----
                                     2047
(%i10) m(341,2047);
                                     682
(%o10)                               ----
                                     2047
(%i11) m(682,2047);
                                     1364
(%o11)                               ----
                                     2047
(%i12) m(1364,2047);
                                     681
(%o12)                               ----
                                     2047
(%i13) m(681,2047);
                                     1362
(%o13)                               ----
                                     2047
                                     

Other version :

         
kill(all);
remvalue(all);

DoublingMap(r):=
block([d,n],
 n:ratnumer(r),
 d:ratdenom(r),
 mod(2*n,d)/d)$


GiveOrbit(r0,OrbitLength):=
block(
 [r,Orbit],
 r:r0,
 Orbit:[r], 
 for i:1 thru OrbitLength step 1 do
        ( r:DoublingMap(r),
          Orbit:endcons(r,Orbit)),
         
 return(sort(Orbit))

)$


r0: 681/2047$
 period : 11;
 
GiveOrbit(r0,period);                            
                                     
        341   597   661   677   681   681   682   1194  1322  1354  1362  1364
(%o6) [----, ----, ----, ----, ----, ----, ----, ----, ----, ----, ----, ----]
       2047  2047  2047  2047  2047  2047  2047  2047  2047  2047  2047  2047
float(%o6);
(%o7) [0.1665852467024914, 0.2916463116756229, 0.3229115779189057, 0.3307278944797264, 0.3326819736199316, 0.3326819736199316, 0.3331704934049829, 0.5832926233512458, 0.6458231558378115, 0.6614557889594529, 0.6653639472398633, 0.6663409868099658]
(%i8) 

child period 21 edit

See also wake 10/21


Maxima 5.45.1 https://maxima.sourceforge.io
using Lisp GNU Common Lisp (GCL) GCL 2.6.12
Distributed under the GNU Public License. See the file COPYING.
Dedicated to the memory of William Schelter.
The function bug_report() provides bug reporting information.
(%i1) bash("r.mac");
(%o1)                             bash(r.mac)
(%i2) load("r.mac");
(%o0)                                r.mac
(%i1) batch("r.mac");

read and interpret /home/a/Dokumenty/maxima/rotation/r.mac
(%i2) kill(all)
(%o0)                                done
(%i1) remvalue(all)
(%o1)                                 []
(%i2) DoublingMap(r):=block([d,n],n:ratnumer(r),d:ratdenom(r),mod(2*n,d)/d)
(%i3) GiveOrbit(r0,OrbitLength):=block([r,Orbit],r:r0,Orbit:[r],
                for i thru OrbitLength do
                    (r:DoublingMap(r),Orbit:endcons(r,Orbit)),
                return(sort(Orbit)))
(%i4) r0:699049/2097151
(%i5) period:21
(%i6) GiveOrbit(r0,period)

       349525   611669   677205   693589   697685   698709   698965   699029
(%o6) [-------, -------, -------, -------, -------, -------, -------, -------, 
       2097151  2097151  2097151  2097151  2097151  2097151  2097151  2097151
       
699045   699049   699049   699050   1223338  1354410  1387178  1395370
-------, -------, -------, -------, -------, -------, -------, -------, 
2097151  2097151  2097151  2097151  2097151  2097151  2097151  2097151


1397418  1397930  1398058  1398090  1398098  1398100
-------, -------, -------, -------, -------, -------]
2097151  2097151  2097151  2097151  2097151  2097151


(%o7)               /home/a/Dokumenty/maxima/rotation/r.mac

(%i9) float(%o6);
(%o9) [0.166666587193769, 0.2916666467984423, 0.3229166616996106, 0.3307291654249026, 0.3326822913562257, 0.3331705728390564, 
0.3332926432097641, 0.333323160802441, 0.3333307902006102, 0.3333326975501525, 0.3333326975501525, 0.3333331743875381, 0.5833332935968846, 
0.6458333233992212, 0.6614583308498053, 0.6653645827124514, 0.6663411456781129, 0.6665852864195282, 0.666646321604882, 0.6666615804012205, 
0.6666653951003051, 0.6666663487750763]
(%i10) 

result without sort(Orbit)

         699049   1398098  699045   1398090  699029   1398058  698965
(%o6)/R/ [-------, -------, -------, -------, -------, -------, -------, 
          2097151  2097151  2097151  2097151  2097151  2097151  2097151
1397930  698709   1397418  697685   1395370  693589   1387178  677205
-------, -------, -------, -------, -------, -------, -------, -------, 
2097151  2097151  2097151  2097151  2097151  2097151  2097151  2097151
1354410  611669   1223338  349525   699050   1398100  699049
-------, -------, -------, -------, -------, -------, -------]
2097151  2097151  2097151  2097151  2097151  2097151  2097151
(%o7)               /home/a/Dokumenty/maxima/rotation/r.mac
float(%o6);
(%o8) [0.3333326975501525, 0.6666653951003051, 0.3333307902006102, 0.6666615804012205, 0.333323160802441, 0.666646321604882, 0.3332926432097641, 
0.6665852864195282, 0.3331705728390564, 0.6663411456781129, 0.3326822913562257, 0.6653645827124514, 0.3307291654249026, 
0.6614583308498053, 0.3229166616996106, 0.6458333233992212, 0.2916666467984423, 0.5833332935968846, 0.166666587193769, 0.3333331743875381, 
0.6666663487750763, 0.3333326975501525]

Child period 27 edit

internal angle 13/27

  • c = -0.739880396515927 +0.115700424748225 i
  • fixed point alpha z = -0.496619178870972 +0.058046457062615 i
  • period 27


The 13/27-wake of the main cardioid is bounded by the parameter rays with the angles

  • 44739241/134217727 or p010101010101010101010101001
  • 44739242/134217727 or p010101010101010101010101010

Child period 34 edit

 
Root point with internal angle 13/34

See image by Arnaud Cheritat [10]

Angles of the wake external rays ( on parameter plane ) in different formats :[11]

 

 

There are 8 589 869 055 components of period 34.

External angle landing on the root points :

 

where denominator d is :

 

internal angle of main cardioid parameter c = root point external angles of the wake (decimal fractions) external angles of the wake ( binary fractions) parabolic fixed point z dynamic rays ( orbit portrait , only numerators)
1/34 (1/d; 2/d)
13/34 -0.392571548476155+0.585781365897037i (4985538889/d ; 4985538890/d) (p0100101001001010010100100101001001; p0100101001001010010100100101001010) -0.3695044586103295 ; 0.3368478218232787 [4985538889,9971077778,2762286373,5524572746,11049145492,4918421801,9836843602,

2493818021,4987636042,9975272084,2770674985,5541349970,11082699940,4985530697, 9971061394,2762253605,5524507210,11049014420,4918159657,9836319314,2492769445, 4985538890,9971077780,2762286377,5524572754,11049145508,4918421833,9836843666, 2493818149,4987636298,9975272596,2770676009,5541352018,11082704036]

Widest sector ( which incudes critical component ) is :

 ( 2492769445/17179869183; 11082704036/17179869183 )

Last componet = component to the left of component including critical point zcr = 0.0. This component is almost not changing when iPeriodChild is increasing , see this video

89 edit

34/89 edit

 
34/89 wake

See image by Arnaud Cheritat [12]

Let's find some info using program Mandel by Wolf Jung :

 
t = 34/89 = 0,382022472 // internal angle = rotational number  
c = -0.390837354761211  +0.586641524321638 i // Parameter c
z = -0.368804231870311  +0.337614334047815 i // fixed point alfa

denominator or external angle ( computed with this program ) :

 

Location of root point using book program by Claude Heiland-Allen:

1/3 = 0.3333
3/8 = 0,375
8/21 = 0,380952381
21/55 = 0,381818182
34/89 = 0,382022472
13/34 = 0,382352941
5/13 = 0,384615385
2/5 = 0.4
1/2 = 0.5

External angles of rays that land on the root point one can compute with book program by Claude Heiland-Allen :

 ./mandelbrot_external_angles 53 -3.9089629378291085e-01 5.8676031775674931e-01 89
.(01001010010010100101001001010010010100101001001010010100100101001001010010100100101001001)
.(01001010010010100101001001010010010100101001001010010100100101001001010010100100101001010)

Converting to other forms using gmp :

decimal fraction =  179622968672387565806504265 / 618970019642690137449562111 
decimal canonical form =  179622968672387565806504265/618970019642690137449562111
binary fraction  = 01001010010010100101001001010010010100101001001010010100100101001001010010100100101001001/11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111 
decimal floating point number : 0.290196557138708685358212600555 
decimal fraction =  179622968672387565806504266 / 618970019642690137449562111 
decimal canonical form =  179622968672387565806504266/618970019642690137449562111
binary fraction  = 01001010010010100101001001010010010100101001001010010100100101001001010010100100101001010/11111111111111111111111111111111111111111111111111111111111111111111111111111111111111111 
decimal floating point number : 0.290196557138708685358212602171 

Note that difference in floating form of external angles :

0.290 196 557 138 708 685 358 212 602 171 
0.290 196 557 138 708 685 358 212 600 555

is on 27-th decimal digit after decimal sign

0.(010 010 100 100 101 001 010 010 010 100 100 101 001 010 010 010 100 101 001 001 010 010 010 100 101 001 001 010 010 01)
0.(010 010 100 100 101 001 010 010 010 100 100 101 001 010 010 010 100 101 001 001 010 010 010 100 101 001 001 010 010 10)

and on 88-th binary digit after decimal sign.

Numerators of the orbit portrait ( external angles of rays landing on the fixed point alfa ) :

179622968672387565806504265
179622968672387565806504265
359245937344775131613008530
99521855046860125776454949
199043710093720251552909898
398087420187440503105819796
177204820732190868762077481
354409641464381737524154962
89849263286073337598747813
179698526572146675197495626
359397053144293350394991252
99824086645896563340420393
199648173291793126680840786
399296346583586253361681572
179622673524482369273801033
359245347048964738547602066
99520674455239339645642021
199041348910478679291284042
398082697820957358582568084
177195375999224579715574057
354390751998449159431148114
89811484354208181412734117
179622968708416362825468234
359245937416832725650936468
99521855190975313852310825
199043710381950627704621650
398087420763901255409243300
177204821885112373368924489
354409643770224746737848978
89849267897759356026135845
179698535795518712052271690
359397071591037424104543380
99824123539384710759524649
199648247078769421519049298
399296494157538843038098596
179622968672387548626635081
359245937344775097253270162
99521855046860057056978213
199043710093720114113956426
398087420187440228227912852
177204820732190319006263593
354409641464380638012527186
89849263286071138575492261
179698526572142277150984522
359397053144284554301969044
99824086645878971154375977
199648173291757942308751954
399296346583515884617503908
179622673524341631785445705
359245347048683263570891410
99520674454676389692220709
199041348909352779384441418
398082697818705558768882836
177195375994720980088203561
354390751989441960176407122
89811484336193782903252133
179622968672387565806504266
359245937344775131613008532
99521855046860125776454953
199043710093720251552909906
398087420187440503105819812
177204820732190868762077513
354409641464381737524155026
89849263286073337598747941
179698526572146675197495882
359397053144293350394991764
99824086645896563340421417
199648173291793126680842834
399296346583586253361685668
179622673524482369273809225
359245347048964738547618450
99520674455239339645674789
199041348910478679291349578
398082697820957358582699156
177195375999224579715836201
354390751998449159431672402
89811484354208181413782693
179622968708416362827565386
359245937416832725655130772
99521855190975313860699433
199043710381950627721398866
398087420763901255442797732
177204821885112373436033353
354409643770224746872066706
89849267897759356294571301
179698535795518712589142602
359397071591037425178285204
99824123539384712907008297
199648247078769425814016594
399296494157538851628033188

Code edit

References edit

  1. muency : elephant valley
  2. Visual Guide To Patterns In The Mandelbrot Set by Miqel
  3. integer number in wikipedia
  4. Integer overflow in wikipedia
  5. planetmath : San Marco fractal
  6. wikipedia : Douady rabbit
  7. wikipedia : orbit portrait
  8. Parameter rays of root points of period p components by A Majewski
  9. Parameter rays of root points of period p components by A Majewski
  10. Arnaud Cheritat - gallery
  11. knowledgedoor calculators: convert_a_ratio_of_integers
  12. Arnaud Cheritat - gallery