Newton’s Algorithm

Helmut Knaust, Department of Mathematical Sciences, UTEP, El Paso TX 79968
hknaust@utep.edu
03/05/2014

Computing Approximations for (Sqrt[2])

Using exact numbers (integers or fractions) gives exact rational approximations:

f[x_]=x^2-2;
g[x_]=x-f[x]/f'[x];
NestList[g,2,8]

Using approximate numbers (decimals) gives “approximate” approximations :

f[x_]=x^2-2;
g[x_]=x-f[x]/f'[x];
NestList[g,2.,8]

Here is how you can get more accuracy:

f[x_]=x^2-2;
g[x_]=x-f[x]/f'[x];
NestList[g,N[2,100],8]

FixedPointList iterates until the result does not change anymore. 
	Do NOT run FixedPointList[g,2]! The computations will never stop (Why?).

FixedPointList[g,N[2,200]]

Running Newton’s Algorithm for f(x)=(x^2)+1

There is no real solution to (x^2)+1=0. The algorithm exhibits chaotic behavior:
1. Outputs are dense on the real line.
2. Outputs exhibit sensitive dependence on initial conditions (“Butterfly Effect”).

f[x_]=x^2+1;
g[x_]=x-f[x]/f'[x];

g[x_]=x-f[x]/f'[x];
NestList[g,2.0,100]

How many iterations does it take before the output lies in the interval [-1000,-999]?

iter=NestWhileList[g,3.0000,Function[s,s<-1000||s>-999]];
Length[iter]

Graphs of the first 500 approximations for various initial values:

Animate[ListPlot[NestList[g,z,500],PlotRange->{-10,10}],{z,-4.5,5.15,1.02}]

A visualization of  successive approximations:

Animate[Graphics[{Red,AbsolutePointSize[4],Map[(Point[{#,0}]&),iter[[1;;m]]]},Axes->{True,False},PlotRange->{{-10,10},{-1,1}},AspectRatio->1/10],{m,1,500,1}]

Animate (and Manipulate) can display arbitrary output, not just graphics commands:

Animate[iter[[m]],{m,1,10,1}]

Running the algorithm for initial conditions 2., and 2. + (1/(10^12)), respectively. The third colum shows the (absolute) difference between corresponding iterations. (This also means that these numerical approximations are unreliable.)

itera=NestList[g,2.0000,100];
iterb=NestList[g,2.0000+1/10^12,100];
Transpose[{itera,iterb,Abs[itera-iterb]}]

An Animation of the Algorithm for f(x)=cos(x)-x

f[x_]=Cos[x]-x;
g[x_]=x-f[x]/f'[x]
h[x_,z_]=f[z]+f'[z](x-z)

You can use a list as an iterator in an Animate or Manipulate command:

Animate[Plot[{f[x],h[x,z]},{x,-5,5},PlotStyle->{Blue,Red},PlotLabel->z,PlotRange->{-5,5},Epilog->{Red,AbsolutePointSize[5],Point[{z,f[z]}],Point[{g[z],0}]}],{z,NestList[g,3.,7]}]