Math 4401 Numerical Methods

Notebook Author/Course Instructor: Barry McQuarrie
Date: Spring 2008
Lecture topic: Accelerating Convergence

Initializations

$TextStyle={FontFamily->"Times",FontSize->12};

Aitken's (Δ^2) Method

We want to improve the convergence of the sequence (Subscript[p, n]) which is given by the following recursion relation:

A Linearly Convergent Sequence

Clear[p,datap]

p[0]:=0.5;
p[n_]:=p[n]=Cos[p[n-1]]

len=20;
datap={};
datap=Table[{n,p[n]},{n,1,len}];

TableForm[SetPrecision[datap,10],TableHeadings->{None,{"n","p"}}]

Clear[p]

error=Table[{n,Abs[datap[[n,2]]-datap[[n-1,2]]]},{n,2,len}];
seqerrplot=ListLogPlot[error,PlotStyle->{Pink,PointSize[0.02]},AxesLabel->{"n","Logarithm of error"}]

Aitken's Method

Print["Aitken's Method"]

Aitkensdata=Table[{n,Part[datap,n,2]-( Part[datap,n+1,2]- Part[datap,n,2])^2/( Part[datap,n+2,2]- 2Part[datap,n+1,2]+ Part[datap,n,2])},{n,1,Length[datap]-2}];

TableForm[SetPrecision[Aitkensdata,15],TableHeadings->{None,{"n","Aitken's p"}}]

Error

error=Table[{n,Abs[Aitkensdata[[n,2]]-Aitkensdata[[n-1,2]]]},{n,2,Length[Aitkensdata]}];
aitkenserrplot=ListLogPlot[error,PlotStyle->{Red,PointSize[0.02]},AxesLabel->{"n","Logarithm of error"}]

Comaprison of Errors: Original Sequence vs Aitken's

Show[seqerrplot,aitkenserrplot,PlotRange->{Automatic,Automatic}]

Does Aitken's (Δ^2) Method Improve Convergence of Sequence Generated by Newton's Method?

We want to solve the equation f(x)=0for x, and then see if we can improve the convergence of the sequence.

f[x_] =Exp[6x]+3Log[2]^2Exp[2x]-Log[8]Exp[4x]-Log[2]^3
Plot[f[x],{x,-0.3,-0.1},PlotStyle->Blue,ImageSize->{500,Automatic}]

Note that the zero here does not appear to be a simple zero, so we might not acheive quadratic convergence with Newton's method.

Newton's Method

Clear[p0,p1,tol, max, i,p,q]

f[x_] =Exp[6x]+3Log[2]^2Exp[2x]-Log[8]Exp[4x]-Log[2]^3
p0=0.0;
tol=10.0^(-8);
max=300;

data={};
i=1;
Print["Newton's Method"]
Print["Input:  ","p0= ",p0," tol= ",tol]
While[i≤max ,
p=p0-f[p0]/f'[p0];
data=Append[data,{i,p,f[p]}];
If[Abs[f[p]]≤tol,{Print["Result: ","p= ",SetPrecision[p,10]," Error ",Abs[p-p0], " Steps= ",i],Break[]}];
p0=p;
i=i+1]
If[i>max,Print["Method Failed to Converge"]]

TableForm[SetPrecision[data,10],TableHeadings->{None,{"n","p","f[p]"}}]

Error in Newton's Method

error=Table[{n,Abs[data[[n,3]]]},{n,1,Length[data]}];
nmerrplot=ListLogPlot[error,PlotStyle->{Pink,PointSize[0.02]},AxesLabel->{"n","Logarithm of error"}]

Aitken's Algorthim

Print["Aitken's Method"]

Aitkensdata=Table[{n,Part[data,n,3]-( Part[data,n+1,3]- Part[data,n,3])^2/( Part[data,n+2,3]- 2Part[data,n+1,3]+ Part[data,n,3])},{n,1,Length[data]-2}];

TableForm[SetPrecision[Aitkensdata,15],TableHeadings->{None,{"n","Aitken's p"}}]

Error in Aitken's Method

error=Table[{n,Abs[Aitkensdata[[n,2]]-Aitkensdata[[n-1,2]]]},{n,2,Length[Aitkensdata]}];
aitkenserrplot=ListLogPlot[error,PlotStyle->{Red,PointSize[0.02]},AxesLabel->{"n","Logarithm of error"}]

Comaprison of Errors: Original Sequence vs Aitken's

Show[nmerrplot,aitkenserrplot,PlotRange->{Automatic,Automatic}]

Aitken's method certainly can improve the convergence of a sequence generated by Newton's Method. Rememeber, that Newton's method does not automatically have quadratic convergence. In this case, the convergence is linear, and so Aitken's method can help improve the convergence.

Using Steffensen's Method to Improve Convergence of Fixed Point Sequence from Linear to Quadratic

We want to solve the equation f(x)=0for x, using a basic fixed point scheme

Fixed Point Algorithm

g[x_]=Cos[x]

p0=0.2;
data =Table[{n,SetPrecision[Nest[g,p0,n],10],If[n>0,Abs[Nest[g,p0,n]-Nest[g,p0,n-1]]]},{n,0,25}];

TableForm[data,TableHeadings->{None,{"i","(Subscript[x, i])","Error"}}]

Error in Fixed Point

error=Table[{n-1,Abs[data[[n,3]]]},{n,2,Length[data]}];
fperrplot=ListLogPlot[error,PlotStyle->{Pink,PointSize[0.02]},AxesLabel->{"n","Logarithm of error"}]

Notice the convergence appears to be linear, as we expect from a fixed point scheme.

Steffensen's Method

Clear[p0,p1,tol, max, i,p,q,data]

g[x_]:=SetPrecision[Cos[x],25]

p0=SetPrecision[0,25];
tol=10.0^(-8);
max=300;

data={};
i=1;
Print["Steffensen's Method"]
Print["Input:  ","p0= ",p0," tol= ",tol]
While[i≤max ,
p1=g[p0];
p2=g[p1];
p=p0-(p1-p0)^2/(p2-2p1+p0);

data=Append[data,{p1,Abs[p1-p0]}];
data=Append[data,{p2,Abs[p2-p1]}];
data=Append[data,{p,Abs[p-p2]}];

If[Abs[p-p2]≤tol,{Print["Result: ","p= ",SetPrecision[p,10]," Error ",Abs[p-p0], " Steps= ",3*i],Break[]}];
p0=p;
i=i+1]
If[i>max,Print["Method Failed to Converge"]]

TableForm[SetPrecision[data,10],TableHeadings->{None,{"p","Error"}}]

Error in Steffensen's Method

error=Table[{n,Abs[data[[n,2]]]},{n,1,Length[data]}];
smerrplot=ListLogPlot[error,PlotStyle->{Red,PointSize[0.02]},AxesLabel->{"n","Logarithm of error"}]

Comaprison of Errors: Fixed Point vs Steffensen's

Show[fperrplot,smerrplot,PlotRange->{Automatic,Automatic}]

Steffensen's method is a fixed point method which does not evaluate derivatives and has quadratic converegence.

Meuller's Method

We want to solve the equation f(x)=0for x.

Algorthim

Clear[p0,p1,p2,d1,d2,dd,data,tol, max, i,p,q]

f[x_] =16x^4-40x^3+5x^2+20x+6

p0=0.5;
p1=-0.5;
p2=0.0;
tol=10.0^(-5);
max=300;

data={};
h1=p1-p0;
h2=p2-p1;
d1=(f[p1]-f[p0])/h1;
d2=(f[p2]-f[p1])/h2;
dd = (d2-d1)/(h2+h1);
i=3;
Print["Meuller's Method"]
Print["Input:  ","p0= ",p0," p1= ",p1," p2= ",p2," tol= ",tol]
While[i≤max ,
b = d2+h2*dd;
cd = N[Sqrt[b^2-4f[p2]*dd]];

If[ Abs[b-cd]<Abs[b+cd],ee=b+cd,ee=b-cd];

h=-2f[p2]/ee;
p=p2+h;

data=Append[data,{p,f[p]}];

If[Abs[h]≤tol,{Print["Result: ","p= ",SetPrecision[p,10]," Error ",Abs[p-p0], " Steps= ",i],Break[]}];

p0=p1;
p1=p2;
p2=p;
h1=p1-p0;
h2=p2-p1;
d1=(f[p1]-f[p0])/h1;
d2=(f[p2]-f[p1])/h2;
dd=(d2-d1)/(h2+h1);
i=i+1]
If[i>max,Print["Method Failed to Converge"]]

TableForm[SetPrecision[data,10],TableHeadings->{None,{"p","f[p]"}}]

Error

error=Table[{n,Abs[data[[n,2]]]},{n,1,Length[data]}];
mmerrplot=ListLogPlot[error,PlotStyle->{Red,PointSize[0.02]},AxesLabel->{"n","Logarithm of error"}]

References:

Heikki Ruskeepaa, Mathematica Navigator, Graphics and Methods of Applied Mathematics (1999) San Diego,Academic Press.
Richard Burden & J. Douglas Faires, Numerical Analysis, 7th ed. (2001) Pacific Grove, Brooks/Cole.