Forward Difference Approximation

Let u(x)=sin(x). We will compute the error in approximating u'(1). First define u(x):

u[x_]:=Sin[x]

Now define the forward difference approximation:

forward[x_,h_]:=((u[x+h]-u[x])/h)

Compute the error when h=0.1:

forward[1,0.1]-u'[1]

Compute the error when h=0.01:

forward[1,0.01]-u'[1]

Compute the error when h=0.001:

forward[1,0.001]-u'[1]

What do you observe?

We can use the Table function to compute the error for many step sizes at once:

Table[forward[1,(10^(-k))]-u'[1],{k,1,6}]//N

Error comparison for first derivative approximations

Next we will compare the errors in approximating the derivative of u(x) at x=1, using the forward difference, backward difference, and centered difference approximations.

We already defined the forward difference approximation. Now define the other two approximations

backward[x_, h_]:=(u[x]-u[x-h])/h;
centered[x_, h_]:=(u[x+h]-u[x-h])/(2h);

To make things easy, we will define three functions that return the error for each approximation of u'(x).

errForward[x_,h_]:=Abs[forward[x,h]-u'[x]];
errBackward[x_,h_]:=Abs[backward[x,h]-u'[x]];
errCentered[x_,h_]:=Abs[centered[x,h]-u'[x]];

We will also make a list of the values of h that we will use:

hVals=Table[(10^(-k)),{k,1,4}]//N

Now compute the errors for each approximation and print them in a table.

TableForm[Table[{errForward[1,h],errBackward[1,h],errCentered[1,h]},{h,hVals}],TableHeadings->{hVals,{"forward error","backward error","centered error"}}]

What do you observe?

Extend the table for a few more rows. What do you observe?

hVals=Table[(10^(-k)),{k,1,8}]//N

TableForm[Table[{errForward[1,h],errBackward[1,h],errCentered[1,h]},{h,hVals}],TableHeadings->{hVals,{"forward error","backward error","centered error"}}]

Something is going wrong for the the centered approximation!

The problem is round-off error. Round-off error is a big concern when doing computations with decimal numbers.

Precision[errCentered[1,0.000001]]

MachinePrecision means that Mathematica is using floating-point arithmetic, as implemented in the computer processor. For most computers, this means 64-bit “double precision” floating-point numbers, which use 53 bits of memory to store significant digits. This amounts to almost 16 decimal digits of precision, as the following commands show:

$MachinePrecision

N[53Log[10,2]]

If we avoid typing decimals, Mathematica will attempt to compute using arbitrary precision arithmetic. However, arbitrary precision calculations may take much more time than machine precision calculations. For example:

errCentered[1,(10^(-6))]

The formula above is exact, but probably not helpful. We can ask Mathematica to convert it to a decimal number with, say, 20 digits:

N[errCentered[1,(10^(-6))],20]

That' s better! Lastly, let’s see what Mathematica says is the precision of this calculation:

Precision[errCentered[1,10^(-6)]]

We can re-do our table using higher precision arithmetic:

hVals=Table[(10^(-k)),{k,1,8}]

N[TableForm[Table[{errForward[1,h],errBackward[1,h],errCentered[1,h]},{h,hVals}],TableHeadings->{hVals,{"forward error","backward error","centered error"}}],20]