CHAPTER 14(b),  Complex Dynamics Continued

14.6 Newton’s Method for Cubics

In fact, Mathematica has a built in command for graphing Julia sets: JuliaSetPlot[f[z],z]:

p[z_]:=z^3-1;
theRoots=z/.NSolve[p[z]==0,z]
cp=Compile[{{z,_Complex}},Evaluate[p[z]]];
n=Compile[{{z,_Complex}},Evaluate[Simplify[z-p[z]/p'[z]]]];
bail=150;
orbitData=Table[NestWhileList[n,x+I y,Abs[cp[#]]>0.01&,1,bail],{y,-1,1,0.01},{x,-1,1,0.01}];
numRoots=Length[Union[theRoots]];
sameRootFunc=Compile[{{z,_Complex}},Evaluate[Abs[3 p[z]/p'[z]]]];
whichRoot[orbit_]:=Module[{i,z},z=Last[orbit];i=1;
Scan[If[Abs[z-#]<sameRootFunc[z],Return[i],i++]&,theRoots];
If[i≤numRoots,{i,Length[orbit]},None]];
rootData=Map[whichRoot,orbitData,{2}];
colorList={{cc,0,0},{cc,cc,0},{0,0,cc}};
cols=rootData/.{{k_Integer,l_Integer}:>(colorList[[k]]/.cc->(1-l/(bail+1))^8),None->{0,0,0}};
Graphics[{Raster[cols]}]

g[z_]:=z^3;

h[z_]:=z-(g[z]-1)/g'[z];

ContourPlot[Arg[FixedPoint[h,x+I y,30]],{x,-1.5,1.5},{y,-1.5,1.5},ContourLines->False,PlotPoints->125,ColorFunction->Hue,Frame->False]

If we use Mathematica’s built in command for graphing Julia sets: JuliaSetPlot[f[z],z], we get a version of the above pictures. This is because the Julia set of the Newton function of p(z)=(z^3)-1, is just the boundary of the basin of attraction of the various zero’s of p(z).

g[z_]:=z^3-1;
h[z_]:=z-(g[z])/g'[z]

JuliaSetPlot[h[z],z]

Exercises 14.6.

1.     This exercise refers to Exercises 14.6 # 4(c) concerning the Relaxed Newton Method. Use Mathematica to show that if g(z)=(((z-a))^2)(z-b), then (Subscript[N, 2])(z)=z-2g(z)/g'(z), the relaxed Newton function, is conjugate to the quadratic p(z)=(z^2)-3/4, via h(z)=(3z+a-4b)/(2(z-a)) (so that p(h(z))=h((Subscript[N, 2])(z))  for all z).

14.7 Newton’s Method and Halley’s Method

We can use Mathematica to obtain formula related to Halley’s Method (see Section 14.7.12 of the text book and Exercise 14.7 # 11).

Clear[f]

P[z_]:=z-(f[z]/f'[z])/D[f[z]/f'[z],z]

This is the Newton function for g(z)=f(z)/f'(z), which is the start of Halley’s method. We use P in place of N because Mathematica reserves N for numerical approximation.

P[z]

Simplify[%]

P'[z]

Simplify[%]

P''[z]

Simplify[%]

P'''[z]

Simplify[%]

The following command substitutes f(a)=0 in P'''(a).

P'''[a]/.{f[a]->0}

Simplify[%]

P''[a]/.{f[a]->0}

Exercises 14.7.

1.     This exercise refers to Exercises 14.7 # 11 concerning the Halley’s Method. Use Mathematica to show as above that (Subscript[N, g])(α)=α, (Subscript[N', g])(α)=0, (Subscript[N'', g])(α)=-f''(α)/f'(α) and (Subscript[N''', g])(α)=-4{f'''(α)/f'(α)-3/((4[f''(α)/f'(α)])^2)}. Note that this is slightly different to what is requested in Exercise 14.7 # 11.

2.     Following on from the last exercise using Halley’s Method, use Mathematica to show that if g=f/(((f'))^(1/2)), then (Subscript[N, g])(α)=α, (Subscript[N', g])(α)=0, (Subscript[N'', g])(α)=0 and (Subscript[N''', g])(α)=-Sf(α),  where Sf is the Schwarzian derivative of f.