/*gcc c.c -lm -Wall./a.outRoot point between period 1 component and period 987 component = c = 0.2500101310666710+0.0000000644946597Internal angle (c) = 1/987Internal radius (c) = 1.0000000000000000*/#include<stdio.h>#include<math.h> // M_PI; needs -lm also#include<complex.h>// numer of hyberolic components ( and it's centers ) of Mandelbrot set intGiveNumberOfCenters(intperiod){//int s = 0;intnum=1;intf;intfMax=period-1;//sqrt(period); // https://stackoverflow.com/questions/11699324/algorithm-to-find-all-the-exact-divisors-of-a-given-integerif(period<1){printf("input error: period should be positive integer \n");return-1;}if(period==1){return1;}num=pow(2,period-1);// minus sum of number of center of all it's divisors (factors)for(f=1;f<=fMax;++f){if(period%f==0){num=num-GiveNumberOfCenters(f);}}returnnum;}intListNumberOfCenters(intperiod){intp=1;intpMax=period;intnum=0;if(period==1||period==2){printf(" for period %d there is only one component\n",period);return0;}for(p=1;p<=pMax;++p){num=GiveNumberOfCenters(p);printf(" for period %d there are %d components\n",p,num);}return0;}intmain(){intperiod=15;ListNumberOfCenters(period);return0;}
"quick and dirty" algorithm : check if then colour c-point with colour n. Here n is a period of attracting orbit and eps is a radius of circle around attracting point = precision of numerical comutations
interval arithmetics and a divide and conquer strategy[12]
"Searching is performed using a quadtree subdivision of the grid which increases search resolution in areas of higher root density. Within a search area first a map of the function magnitude is made at moderate resolution, then points of local minimum are refined using 640 bit fixed point math. Root searching requires twice the precision of the desired result due to multiplying very small numbers. As roots are added the quadtree is subdivided further and the search continues. Only the upper half plane is searched due to symmetry.
Quadtree depth was purposely limited on higher periods to reduce run time since the number of roots becomes unmanageable." James Artherton
Centers of hyperbolic components for period 10Centers of hyperbolic components for period 1-10 computed using irreducible polynomials
Center of hyperbolic component is a point of parameter plain, for which periodic z-cycle is superattracting.
It gives a system of 2 equations defining centers of period hyperbolic components :
"These polynomials have integer coefficients. They can be obtained by a recurrent relation on the degree. Let Pd the polynomial of degree d. We have" [13]
GiveG[p]:=
block(
[f:divisors(p),
t:1], /* t is temporary variable = product of Gn for (divisors of p) other than p */
f:delete(p,f), /* delete p from list of divisors */
if p=1
then return(Fn(p,0,c)),
for i in f do t:t*GiveG[i],
g: Fn(p,0,c)/t,
return(ratsimp(g))
)$
Here is a table with degree of and for periods 1-10 and precision needed to compute coefficients of these functions (the roots can be calculated with lower precision, provided you work with the unexpanded form ).
period
1
1
1
16
2
2
1
16
3
4
3
16
4
8
6
16
5
16
15
16
6
32
27
16
7
64
63
32
8
128
120
64
9
256
252
128
10
512
495
300
11
1024
1023
600
Where:
fpprec is the number of significant decimal digits for arithmetic on bigfloat numbers[17]
Precision can be estimated as bigger than size of binary form of greatest coefficient :
period
1
1
1
16
2
1
1
16
3
3
2
16
4
6
2
16
5
15
7
16
6
27
13
16
7
63
34
32
8
120
67
64
9
252
144
128
10
495
286
300
11
1023
596
600
Here is Maxima code for
period_Max:11;
/* ----------------- definitions -------------------------------*/
/* basic function = monic and centered complex quadratic polynomial
http://en.wikipedia.org/wiki/Complex_quadratic_polynomial
*/
f(z,c):=z*z+c $
/* iterated function */
fn(n, z, c) :=
if n=1 then f(z,c)
else f(fn(n-1, z, c),c) $
/* roots of Fn are periodic point of fn function */
Fn(n,z,c):=fn(n, z, c)-z $
/* gives irreducible divisors of polynomial Fn[p,z=0,c] */
GiveG[p]:=
block(
[f:divisors(p),t:1],
g,
f:delete(p,f),
if p=1
then return(Fn(p,0,c)),
for i in f do t:t*GiveG[i],
g: Fn(p,0,c)/t,
return(ratsimp(g))
)$
/* degree of function */
GiveDegree(_function,_var):=hipow(expand(_function),_var);
log10(x) := log(x)/log(10);
/* ------------------------------*/
file_output_append:true; /* to append new string to file */
grind:true;
for period:1 thru period_Max step 1 do
(
g[period]:GiveG[period], /* function g */
d[period]:GiveDegree(g[period],c), /* degree of function g */
cf[period]:makelist(coeff(g[period],c,d[period]-i),i,0,d[period]), /* list of coefficients */
cf_max[period]:apply(max, cf[period]), /* max coefficient */
disp(cf_max[period]," ",ceiling(log10(cf_max[period]))),
s:concat("period:",string(period)," cf_max:",cf_max[period]),
stringout("max_coeff.txt",s)/* save output to file as Maxima expressions */
);
All these methods shows centers for period n and its divisors.
Animation showing colors proportional to absolute value of iterations 1-20Absolute value of first iterationColor shows in which quadrant 5-th iteration lands
color is proportional to magnitude of zn[edit | edit source]
This is radial nth-decomposition of exterior of Mandelbrot set ( compare it with n-th decomposition of LSM/M )
4 colors are used because there are 4 quadrants :
re(z_n) > 0 and im(z_n) > 0 ( 1-st quadrant )
re(z_n) < 0 and im(z_n) > 0 ( 2-nd quadrant )
re(z_n) < 0 and im(z_n) < 0 ( 3-rd quadrant )
re(z_n) > 0 and im(z_n) < 0 (4-th quadrant ).
"... when the four colors are meeting somewhere, you have a zero of q_n(c), i.e., a center of period dividing n. In addition, the light or dark color shows if c is in M." ( Wolf Jung )
Here is fragment of c code for computing 8-bit color for point c = cx + cy * i :
/* initial value of orbit = critical point Z= 0 */Zx=0.0;Zy=0.0;Zx2=Zx*Zx;Zy2=Zy*Zy;/* the same number of iterations for all points !!! */for(Iteration=0;Iteration<IterationMax;Iteration++){Zy=2*Zx*Zy+Cy;Zx=Zx2-Zy2+Cx;Zx2=Zx*Zx;Zy2=Zy*Zy;};/* compute pixel color (8 bit = 1 byte) */if((Zx>0)&&(Zy>0))color[0]=0;/* 1-st quadrant */if((Zx<0)&&(Zy>0))color[0]=10;/* 2-nd quadrant */if((Zx<0)&&(Zy<0))color[0]=20;/* 3-rd quadrant */if((Zx>0)&&(Zy<0))color[0]=30;/* 4-th quadrant */
One can also find cpp code in function quadrantqn from class mndlbrot in mndlbrot.cpp file from source code of program Mandel by Wolf Jung[21]
The differences from standard loop are :
there is no bailout test (for example for max value of abs(zn) ). It means that every point has the same number of iterations. For large n overflow is possible. Also one can't use test Iteration==IterationMax
Iteration limit is relatively small ( start for example IterationMax = 8 )
centers for periods 1-20 ( sum = 1 045 239 ) by James Artherton. Made with custom software using 320 bit fixed point math (approx. 92 decimal digits) and a quadtree subdivision of the grid to increase the search resolution in areas with higher root density. This adapts to very fine resolution in the tips of the antennae.[24]
a database of all islands up to period 16, found by tracing external rays (period, islandhood, angled internal address, lower external angle numerator, denominator, upper numerator, denominator, orientation, size, centre realpart, imagpart)[27]
real centers of components computed by Jay R. Hill[28]
It is a package for Ubuntu
A solver for Mandelbrot polynomials based on MPSolve.
A polynomial rootfinder that can determine arbitrary precision approximations with guaranteed inclusion radii.
It supports exploiting of polynomial structures such as sparsisty and allows for polynomial implicitly defined or in some non standard basis.
This binary, provided as an example of custom polynomial implementation in the MPSolve package, uses the particular structure of this family of polynomials to develop
an efficient solver that exhibit experimental O(n^2) complexity.
This package contains the specialization of mpsolve to Mandelbrot polynomials.
period 33 center
.gif)
