# Plotting fractals in Mathematica for S iteration

Dears, I'm trying to plot the attached picture figures in Mathematica software.

Let $$T(z)=z-\frac{P(z)}{P^{\prime}(z)}$$. The S iterative scheme is for $$z_{0}$$ a starting point and $$y_{n}=(1-\beta_{n})z_{n}+\beta_{n}Tz_{n}$$, $$z_{n+1}=(1-\alpha_{n})Tz_{n}+\alpha_{n}Ty_{n}$$, n=0,1,2,....

Now if $$P(z)=z^{3}-1$$ and considered the square domains centered at the origin in the complex plane $$D=[-2,2]\times[-2,2]$$ and To generate the basins of attraction and to study the behaviors of Mann iteration procedure, we take the corresponding square domain and divide it into $$250\times 250$$ grids. By using Newton’s operator N(z), we generate the sequence $$\{z_{n}\}$$; n = 1; 2; 3; corresponding to an iteration procedure, starting at each grid point z0. If the sequence $$\{z_{n}\}$$; n = 1; 2; 3; attempt a root of a polynomial with the accuracy of $$10^{-4}$$ in a $$k \leq 13$$ number of iterations, then the converging point $$z_0$$ is colored in a color assigned to k, otherwise the point is colored in white. The fig (b) obtain fig if $$\alpha_{n}=\beta_{n}=0.9$$. The complete paper is free to download at https://www.researchgate.net/publication/307016155_A_comparison_of_some_fixed_point_iteration_procedures_by_using_the_basins_of_attraction

I have tried the code in Mathematica but can not find any images like this.

P[z_] := z^3 - 1;
T[z_] := z - (T[z]/T'[z]);
α[n_] := 0.9;
β[n_] := 0.9;
z[0] = AnyTrue;
y[n_] := (1 - β[n])*z[n] + β[n]*T[z[n]];
z[n_] := (1 - α[n - 1])*T[z[n - 1]] + α[n - 1]*T[y[n - 1]];


## 1 Answer

Implementing this in Mathematica is not difficult. Lets start by implementing some iteration operators (the interested reader might add additional operators from the paper):

ClearAll[Tnewton,Tmann,Tkahn,TM]
Tnewton[P_,dPdz_]:=Module[{x,n},Inactive[Function][{x,n},x-P[x]/dPdz[x]]//Activate]
Tmann[alpha_][P_,dPdz_]:=Module[{x,n},Inactive[Function][{x,n},(1-alpha[n])x+alpha[n] Tnewton[P,dPdz][x,0]]//Activate]
Tkahn[alpha_][P_,dPdz_]:=Module[{x,n},Inactive[Function][{x},With[{y=(1-alpha[n])x+alpha[n] Tnewton[P,dPdz][x,0]},Simplify@Tnewton[P,dPdz][y,0]]]//Activate]
TM[alpha_][P_,dPdz_]:=Module[{x,n},Inactive[Function][{x},With[{y=(1-alpha[n])x+alpha[n] Tnewton[P,dPdz][x,0]},Simplify@Tnewton[P,dPdz][Tnewton[P,dPdz][y,0],0]]]//Activate]


These functors generate iteration functions after providing the polynomial, its derivative and optional parameters.

We continue by defining some test function, its derivative, resulting iterators and its roots:

P3=Function[{z},z^3-1];
dP3dz=Module[{z},Inactive[Function][{z},D[P3[z],z]]//Activate];
P3Tnewton=Tnewton[P3,dP3dz];
P3Tmann=Tmann[0.9&][P3,dP3dz];
P3Tkahn=Tkahn[0.9&][P3,dP3dz];
P3TM=TM[0.9&][P3,dP3dz];
P3roots=Solve[P3[z]==0,z][[All,1,2]];

P5=Function[{z},z^5+z^4+z^3+z^2+z];
dP5dz=Module[{z},Inactive[Function][{z},D[P5[z],z]]//Activate];
P5Tnewton=Tnewton[P5,dP5dz];
P5Tmann=Tmann[0.9&][P5,dP5dz];
P5Tkahn=Tkahn[0.9&][P5,dP5dz];
P5TM=TM[0.9&][P5,dP5dz];
P5roots=Solve[P5[z]==0,z][[All,1,2]];


The main iteration routine to compute the basins of attraction is

iteration[T_,Proots_,nmax_:13,epsilon_:10^-4,printQ_:False][z_]:=Module[{it,n,itn},
it[0]=z;
it[n_]:=it[n]=T[it[n-1],n];
n=0;
While[n<=nmax,
itn=it[n];
If[printQ,Print[{n,itn,Abs[itn-#]&/@Proots}]];
If[Min[Abs[itn-#]&/@Proots]<epsilon,Return[n];];
n++;
];
Return[-1];
]


This function returns $$-1$$ if the iteration fails to converge (up to a absolute precision epsilon) to one of the specified roots in the specified number of maximal iterations nmax. If the iteration converges the number of iteration steps n is returned.

Now we can generate some plots using

m=250;
aRange=Subdivide [-2.,2.,m];
bRange=aRange;

P3newton=Quiet@ParallelTable[{a,b,Quiet@iteration[P3Tnewton,P3roots,13][a+I b]},{a,aRange},{b,bRange},DistributedContexts->All];
P3mann=Quiet@ParallelTable[{a,b,Quiet@iteration[P3Tmann,P3roots,13][a+I b]},{a,aRange},{b,bRange},DistributedContexts->All];
P3khan=Quiet@ParallelTable[{a,b,Quiet@iteration[P3Tkahn,P3roots,13][a+I b]},{a,aRange},{b,bRange},DistributedContexts->All];
P3M=Quiet@ParallelTable[{a,b,Quiet@iteration[P3TM,P3roots,13][a+I b]},{a,aRange},{b,bRange},DistributedContexts->All];

ClearAll[m,aRange,bRange];

colorFuntion=Function[{k},If[k<0||k>13,White,{RGBColor[0.001,0.001,0.001],RGBColor[0.001,0.001,0.001],RGBColor[0.494,0.001,0.494],RGBColor[0.994,0.007,0.001],RGBColor[0.992,0.993,0.001],RGBColor[0.012,0.995,0.002],RGBColor[0.001,0.993,0.992],RGBColor[0.013,0.006,0.995],RGBColor[0.991,0.001,0.994],RGBColor[0.497,0.5,0.995],RGBColor[0.496,0.495,0.494],RGBColor[0.499,0.003,0.002],RGBColor[0.5,0.5,0.001],RGBColor[1.,1.,1.]}[[Floor[k]+1]]]];

plot=ListDensityPlot[Flatten[#1,1],PlotLegends->Placed[BarLegend[{colorFuntion,{1,14}}],Right],InterpolationOrder->0,PlotRange->All,ColorFunction->colorFuntion,ColorFunctionScaling->False,FrameLabel->{"a","b"},PlotLabel->#2<>" iteration for Subscript[P, 3](z)=z^3-1, where z=a+I b",ImageSize->350]&;
Grid[{{plot[P3newton,"Newton"],plot[P3mann,"Mann"]},{plot[P3khan,"Khan"],plot[P3M,"M"]}}]
Clear[plot]


resulting in

Similar results can be produced for P5 (following the reference the plot range is changed here):

m=250;
aRange=Subdivide [-10.,10.,m];
bRange=aRange;

P5newton=Quiet@ParallelTable[{a,b,Quiet@iteration[P5Tnewton,P5roots,13][a+I b]},{a,aRange},{b,bRange},DistributedContexts->All];
P5mann=Quiet@ParallelTable[{a,b,Quiet@iteration[P5Tmann,P5roots,13][a+I b]},{a,aRange},{b,bRange},DistributedContexts->All];
P5khan=Quiet@ParallelTable[{a,b,Quiet@iteration[P5Tkahn,P5roots,13][a+I b]},{a,aRange},{b,bRange},DistributedContexts->All];
P5M=Quiet@ParallelTable[{a,b,Quiet@iteration[P5TM,P5roots,13][a+I b]},{a,aRange},{b,bRange},DistributedContexts->All];

ClearAll[m,aRange,bRange];

colorFuntion=Function[{k},If[k<0||k>13,White,{RGBColor[0.001,0.001,0.001],RGBColor[0.001,0.001,0.001],RGBColor[0.494,0.001,0.494],RGBColor[0.994,0.007,0.001],RGBColor[0.992,0.993,0.001],RGBColor[0.012,0.995,0.002],RGBColor[0.001,0.993,0.992],RGBColor[0.013,0.006,0.995],RGBColor[0.991,0.001,0.994],RGBColor[0.497,0.5,0.995],RGBColor[0.496,0.495,0.494],RGBColor[0.499,0.003,0.002],RGBColor[0.5,0.5,0.001],RGBColor[1.,1.,1.]}[[Floor[k]+1]]]];

plot=ListDensityPlot[Flatten[#1,1],PlotLegends->Placed[BarLegend[{colorFuntion,{1,14}}],Right],InterpolationOrder->0,PlotRange->All,ColorFunction->colorFuntion,ColorFunctionScaling->False,FrameLabel->{"a","b"},PlotLabel->#2<>" iteration for Subscript[P, 5](z)=z^5+z^4+z^3+z^2+z, where z=a+I b",ImageSize->350]&;
Grid[{{plot[P5newton,"Newton"],plot[P5mann,"Mann"]},{plot[P5khan,"Khan"],plot[P5M,"M"]}}]
Clear[plot]


resulting in

EDIT: Following the first three comments on this answer I modified it to include M-iteration with fixed $$\alpha_n=0.9$$ and $$T_\mathrm{Newton}$$ as operator $$T$$. Compared to the original answer I changed the operators $$T$$ as well as the iteration routine to work with non-constant/functional $$\alpha$$: instead of the constant function 0.9& one could try out any function of n, e.g. 1.+#&. Concerning the M-iteration from the comments: I tried it only with $$T=T_\mathrm{Newton}$$ and constant $$\alpha_n$$ which seems to work rather poorly for $$P_5(z)$$. Other choices for $$T$$ and $$\alpha_n$$ might work better. To speed up the data generation I switched to using ParallelTable but most of the execution time is spent generating the plots.

• Dear @N0va I will thank you very very very much for your remarkable help. I will be more thankfull If you can do the following. Apr 10 at 17:26
• (1) Can I get exact colors as I attached the Picure? (2) For P5(z)=z^5+z^4+z^3+z^2+z then the Mann iteration figure is different from the paper. (3) How I set the following M-iteration which more general than the Khan-iteration $s_{n}=(1-\alpha_{n})z_{n}+\alpha_{n} Tz_{n}$, $y_{n}=Ts_{n}$ and $z_{n+1}=Ty_{n}$. Apr 10 at 17:30
• The M-iteration Paper is in the link researchgate.net/publication/… Apr 10 at 17:34
• I edited the answer to accommodate your wishes, which are all somewhat trivial modifications of the existing code. (1) Yes, if one bothers to extract the colors from the reference plots. (2) Have you adjusted the plot range to match the one given in the reference?! (3) Just add one additional stage to the Khan-iteration?! I have not looked into choices for $T$ and $\alpha_n$ in this case beyond the ones mentioned in the answer. Feel free to adapt/experiment.
– N0va
Apr 11 at 13:48
• Dear @N0va thanks for your kind help. I think you are too much expert in mathematica complicated figures and graphs. I have upload a similar question at mathematica.stackexchange.com/questions/266827/… If you can provide your kind answar Apr 24 at 17:30