Week 3.5 - Deriving 4th order approximation of a 2nd order derivative using Taylor Table method

AIM:

To derive the fourth order approximation of a second order derivative using central differencing scheme skewed right sided difference and skewed left sided difference with the help of taylor table method and to compare the analytical derivative with obtained values for the given function.

SCHEME DERIVATIONS:

CENTRAL DIFFERENCING SCHEME:

NO OF NODES:

The numbers of nodes = p+q-1

Where p is order of derivative.

q is order of approximations.

The no of nodes = 2+4-1 = 5

DERIVATION:

To derive fourth order approximation for the second order derivative we start with these approximation,

Taylor series derivation for the above terms,

Similarly expand all the terms in the equation (1) and are represented in the form of table using taylor table method.

af(i-2) a -2a 2a -4/3 2a/3

bf(i-1) b -b -b/2 -b/6 b/24

cf(i) c 0 0 0 0

df(i+1) d d d/2 d/6 d/24

ef(i+2) e 2e 2e 4e/3 2e/3

sum of

co-efficient 0 0 1 0 0

The linear equation obtained from the above table are solved using matrix inversion method,

A*X = B

X = inv(A)*B

we get the coefficients as,

a = -0.0833, b = 1.3333, c = -2.5000, d = 1.3333, e = -0.0833

substituting these coefficient in equation (1) we get,

SKEWED RIGHT HAND DIFFERENCE SCHEME:

NO OF NODES:

The number of nodes = p+q

Where p is order of derivative.

q is order of approximations.

The no of nodes = 2+4 = 6

DERIVATION:

To derive fourth order approximation for the second order derivatives we start with these approximation,

Expand all the terms in the equation (2) using taylor series expansion and are represented in the form of table using taylor table method.

af(i) a 0 0 0 0 0

bf(i+1) b b b/2 b/6 b/24 b/120

cf(i+2) c 2c 4c/2 8c/6 16c/24 32c/120

df(i+3) d 3d 9d/2 27d/6 81d/24 243d/120

ef(i+4) e 4e 16e/2 64e/6 256e/24 1024e/120

gf(i+5) g 5g 25g/2 125g/6 625g/24 3125g/120

Sum of

co-efficient 0 0 1 0 0 0

The linear equation obtained from the above table are solved using matrix inversion method,

A.X = B

X = inv(A)*B

We have co-efficients as,

a = 3.7500, b = -12.8333 , c = 17.8333, d = -13.0000 e = 5.0833 g = -0.8333

SKEWED LEFT HAND DIFFERENCE SCHEME:

NO OF NODES:

The number of nodes = p+q

Where p is order of derivative.

q is order of approximations

The no of nodes = 2+4 = 6.

DERIVATION:

To derive a fourth order approximation for the second order derivative we start with these approximation,

Expand all terms in the equation (3) using taylor series expansion and are represented in the form of table using taylor table method.

af(i-5) a -5a 25a/2 -125a/6 625a/24 -3125a/120

bf(i-4) b -4b 16b/4 -64b/6 256b/24 -1024b/120

cf(i-3) c -3c 9c/2 -27c/6 81c/24 - 243c/120

df(i-2) d -2d 4d/2 -8d/6 16d/24 -32d/120

ef(i-1) e -e e/2 -e/6 e/24 -e/120

gf(i) g 0 0 0 0 0

Sum of

co-efficient 0 0 1 0 0 0

The linear equation obtained from the above table are solved using matrix inversion method,

A.X =B

X =inv(A)*B

We get the coefficients as

a = -0.83333. b = 5.0833, c = -13.0000 , d = 17.8333 , e = -12.8333, g = 3.7500

Substituting these coefficient in equation(3) we get

CENTRAL DIFFERENCING FUNCTION:

function Error_central_differencing = fourth_order_central_differencing35(x,dx)

% analytic function = exp(x)*cos(x)

% analytic derivative
% f''(x) = -2*exp(x)*sinx(x)

analytic_derivative = -2*exp(x)*sin(x);

% numerical derivative

b = 1.3333;
c = -2.5000;
d = 1.3333;
e = -0.0833;
a = -0.0833;


f = @(x)exp(x)*cos(x);

central_differencing =( a*f(x-2*dx)+b*f(x-dx)+c*f(x)+d*f(x+dx)+e*f(x+2*dx))/dx^2;

%error

Error_central_differencing = abs(central_differencing - analytic_derivative );

end

LEFT SIDE SKEW FUCTION:

function Error_left_sided_skew = fourth_order_left_sided_skew35(x,dx)

% analytic function = exp(x)*cos(x)

% analytic derivative
% f''(x) = -2*exp(x)*sin(x)

analytic_derivative = -2*exp(x).*sin(x);

% numerical derivative
a = -0.8333;
b = 5.0833;
c = -13.0000;
d = 17.8333;
e = -12.8333;
g = 3.7500;

f = @(x) exp(x)*cos(x);
left_sided_skew = (a*f(x-5*dx)+b*f(x-4*dx)+c*f(x-3*dx)+d*f(x-2*dx)+e*f(x-dx)+g*f(x))/(dx)^2;

%error
Error_left_sided_skew = abs(left_sided_skew - analytic_derivative );

RIGHT SIDED SKEW FUNCTION:

function Error_right_sided_skew = fourth_order_right_sided_skew35(x,dx)

% analytic function = exp(x)*cos(x)

% analytic derivative
% f''(x) = -2*exp(x)*sin(x)

analytic_derivative = -2*exp(x)*sin(x);

% numerical derivative

a = 3.7500;
b = -12.8333;
c = 17.8333;
d = -13.0000;
e = 5.0833;
g = -0.8333;

f = @(x)exp(x)*cos(x);

right_sided_skew =(a*f(x)+b*f(x+dx)+c*f(x+2*dx)+d*f(x+3*dx)+e*f(x+4*dx)+g*f(x+5*dx))/(dx^2);


%error
Error_right_sided_skew = abs(right_sided_skew - analytic_derivative );

end

MAIN CODE:

clear all
close all
clc

%  Central differencing scheme

A = [1 1 1 1 1;(-2) (-1) 0 1 2;2 (0.5) 0 (0.5) 2; (-4/3) (-1/6) 0 (1/6) (4/3); (2/3) (1/24) 0 (1/24) (2/3) ];
B = [0;0;1;0;0];
X = inv(A)*B;

% Skewed right hand difference scheme

A2 = [1 1 1 1 1 1;0 1 2 3 4 5;0 (0.5) 2 (9/2) 8 (25/2);0 (1/6) (8/6) (27/6) (64/6) (125/6); 0 (1/24) (16/24) (81/24) (256/24) (625/24); 0 (1/120) (32/120) (243/120) (1024/120) (3125/120)];
B2 = [0;0;1;0;0;0];
X2 = inv(A2)*B2;

% value of x
x = pi/3;

% range of dx
dx = linspace(pi/4,pi/4000,20);

% for loop
for i = 1:length(dx)
central_differencing(i) = fourth_order_central_differencing35(x,dx(i)) ;

right_sided_skew(i) = fourth_order_right_sided_skew35(x,dx(i));

left_sided_skew(i) = fourth_order_left_sided_skew35(x,dx(i));

end

%plotting

figure(1)
loglog(dx,central_differencing ,'r')
hold on
loglog(dx,right_sided_skew,'g')
hold on
loglog(dx,left_sided_skew,'b')
xlabel('range of dx')
ylabel('truncation error')
title('range of dx vs error')
legend('central differencing','skewed right sided difference','skewed left sided difference')

CONCLUSION:

We assign the value of X and range of dx.

we define an anonymous function of analytical function and its derivative.

Use the matrices from above discussed section and solved it using inversion method .

Assigned the rows and column to the appropriate matrices.

Loop the order of approximation expression in a for loop under the range of dx. Here absolute error is also calculated from the derivative value calculated analytically.

The calculated valued are plotted to do a bettter comparison .

In first plot we can see that the value of all the approximation approached toward the analyticallly calculated value as he size of dx gets smaller.

On the second plot we can see that the error in the central difference is lowest with respect to other two.