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

AIM:

To derive the fourth order approximations 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 number 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 a fourth order approximation for the second order derivative,we start with these approximation,

∂2u∂x2=af(i−2)+bf(i−1)+cf(i)+df(i+1)+ef(i+2)dx2∂2u∂x2=af(i-2)+bf(i-1)+cf(i)+df(i+1)+ef(i+2)dx2..........(1)

taylor series derivation for the above terms,

af(i−2)=af(i)−2af'(i)dx1!+4af''(i)dx22!+8af'''(i)dx33!+16af''''(i)dx44!af(i-2)=af(i)-2af′(i)dx1!+4af′′(i)dx22!+8af′′′(i)dx33!+16af′′′′(i)dx44!

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

dx2f''(i)dx2f′′(i)	f(i)	f'(i)dx	f''(i)dx2f′′(i)dx2	f'''(i)dx3f′′′(i)dx3	f''''(i)dx4f′′′′(i)dx4
af(i-2)	a	-2a	2a	−4a3-4a3	2a32a3
bf(i-1)	b	-b	b2b2	−b6-b6	b24b24
cf(i)	c	0	0	0	0
df(i+1)	d	d	d2d2	d6d6	d24d24
ef(i+2)	e	2e	2e	4e34e3	2e32e3
sum of co_efficients	0	0	1	0	0

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

A.X = B

⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣11111−2−10122120122−43−160164323124012423⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦.⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣abcde⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦=⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣00100⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦[11111-2-10122120122-43-160164323124012423].[abcde]=[00100]

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,

∂2u∂x2=(−0.0833)f(i−2)+(1.3333)f(i−1)+(−2.5000)f(i)+(1.3333)f(i+1)+(−0.0833)f(i+2)dx2∂2u∂x2=(-0.0833)f(i-2)+(1.3333)f(i-1)+(-2.5000)f(i)+(1.3333)f(i+1)+(-0.0833)f(i+2)dx2

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 a fourth order approximation for the second order derivative,we start with these approximation,

∂2u∂x2=af(i)+bf(i+1)+cf(i+2)+df(i+3)+ef(i+4)+gf(i+5)dx2∂2u∂x2=af(i)+bf(i+1)+cf(i+2)+df(i+3)+ef(i+4)+gf(i+5)dx2..........(2)

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

dx2f''(i)dx2f′′(i)	f(i)	f'(i)dx	f''(i)dx2f′′(i)dx2	f'''(i)dx3f′′′(i)dx3	f''''(i)dx4f′′′′(i)dx4	f'''''(i)dx5f′′′′′(i)dx5
af(i)	a	0	0	0	0	0
bf(i+1)	b	b	b2b2	b6b6	b24b24	b120b120
cf(i+2)	c	2c	4c24c2	8c68c6	16c2416c24	32c12032c120
df(i+3)	d	3d	9d29d2	27d627d6	81d2481d24	243d120243d120
ef(i+4)	e	4e	16e216e2	64e664e6	256e24256e24	1024e1201024e120
gf(i+5)	g	5g	25g225g2	125g6125g6	625g24625g24	3125g1203125g120
sum of co_efficients	0	0	1	0	0	0

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

A.X = B

⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣11111101234501242921622520168627664612560124162481242562462524011203212024312010241203125120⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦.⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣abcdeg⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦=⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣001000⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦[11111101234501242921622520168627664612560124162481242562462524011203212024312010241203125120].[abcdeg]=[001000]

X = inv(A)*B

we get the coefficients as,

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

substituting these coefficient in equation(2) we get,

∂2u∂x2=(3.7500)f(i)+(−12.8333)f(i+1)+(17.8333)f(i+2)+(−13.0000)f(i+3)+(5.0833)f(i+4)+(−0.8333)f(i+5)dx2∂2u∂x2=(3.7500)f(i)+(-12.8333)f(i+1)+(17.8333)f(i+2)+(-13.0000)f(i+3)+(5.0833)f(i+4)+(-0.8333)f(i+5)dx2

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,

∂2u∂x2=af(i−5)+bf(i−4)+cf(i−3)+df(i−2)+ef(i−1)+gf(i)dx2∂2u∂x2=af(i-5)+bf(i-4)+cf(i-3)+df(i-2)+ef(i-1)+gf(i)dx2..........(3)

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

dx2f''(i)dx2f′′(i)	f(i)	f'(i)dx	f''(i)dx2f′′(i)dx2	f'''(i)dx3f′′′(i)dx3	f''''(i)dx4f′′′′(i)dx4	f'''''(i)dx5f′′′′′(i)dx5
af(i-5)	a	-5a	25a225a2	−125a6-125a6	625a24625a24	−3125a120-3125a120
bf(i-4)	b	-4b	16b216b2	−64b6-64b6	256b24256b24	−1024b120-1024b120
cf(i-3)	c	-3c	9c29c2	−27c6-27c6	81c2481c24	−243c120-243c120
df(i-2)	d	-2d	4d24d2	−8d6-8d6	16d2416d24	−32d120-32d120
ef(i-1)	e	-e	e2e2	−e6-e6	e24e24	−e120-e120
gf(i)	g	0	0	0	0	0
sum of co_efficients	0	0	1	0	0	0

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

A.X = B

⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣111111−5−4−3−2−102521629242120−1256−646−276−86−1606252425624812416241240−3125120−1024120−243120−32120−11200⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦.⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣abcdeg⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦=⎡⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢⎣001000⎤⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥⎦[111111-5-4-3-2-102521629242120-1256-646-276-86-1606252425624812416241240-3125120-1024120-243120-32120-11200].[abcdeg]=[001000]

X = inv(A)*B

we get the coefficients as,

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

substituting these coefficient in equation(3) we get

∂2u∂x2=(−0.8333)f(i−5)+(5.0833)f(i−4)+(−13.0000)f(i−3)+(17.8333)f(i−2)+(−12.8333)f(i−1)+(3.7500)f(i)dx2∂2u∂x2=(-0.8333)f(i-5)+(5.0833)f(i-4)+(-13.0000)f(i-3)+(17.8333)f(i-2)+(-12.8333)f(i-1)+(3.7500)f(i)dx2

MATLAB CODE:

MAIN CODE:

clear all
close all
clc

%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_differencing(x,dx(i));
    
    right_sided_skew(i) = fourth_order_right_sided_skew(x,dx(i));
    
    left_sided_skew(i) = fourth_order_left_sided_skew(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','location','northwest')

FUNCTION CODE FOR CENTRAL DIFFERENCING:

function Error_central_differencing = fourth_order_central_differencing(x,dx)

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

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

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

%numerical derivative

%central differencing
central_differencing = (((-0.0833)*((exp(x-(2*dx)))*(cos(x-(2*dx)))))+((1.3333)*((exp(x-dx))*(cos(x-dx))))-((2.5000)*((exp(x))*(cos(x))))+((1.3333)*((exp(x+dx))*(cos(x+dx))))-((0.0833)*((exp(x+(2*dx)))*(cos(x+(2*dx))))))/(dx^2);

%error
Error_central_differencing = abs(central_differencing- analytical_derivative);
end

FUNCTION CODE FOR SKEWED RIGHT-SIDED DIFFERENCE:

function Error_right_sided_skew = fourth_order_right_sided_skew(x,dx)

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

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

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

%numerical derivative

%central differencing
right_sided_skew = (((3.7500)*((exp(x))*(cos(x))))-((12.8333)*((exp(x+dx))*(cos(x+dx))))+((17.8333)*((exp(x+(2*dx)))*(cos(x+(2*dx)))))-((13.0000)*((exp(x+(3*dx)))*(cos(x+(3*dx)))))+((5.0833)*((exp(x+(4*dx)))*(cos(x+(4*dx)))))-((0.8333)*((exp(x+(5*dx)))*(cos(x+(5*dx))))))/(dx^2);

%error
Error_right_sided_skew = abs(right_sided_skew- analytical_derivative);
end

FUNCTION CODE FOR SKEWED LEFT-SIDED DIFFERENCE:

function Error_left_sided_skew = fourth_order_left_sided_skew(x,dx)

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

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

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

%numerical derivative

%central differencing
left_sided_skew = (((-0.8333)*((exp(x-(5*dx)))*(cos(x-(5*dx)))))+((5.0833)*((exp(x-(4*dx)))*(cos(x-(4*dx)))))-((13.0000)*((exp(x-(3*dx)))*(cos(x-(3*dx)))))+((17.8333)*((exp(x-(2*dx)))*(cos(x-(2*dx)))))-((12.8333)*((exp(x-dx))*(cos(x-dx))))+((3.7500)*((exp(x))*(cos(x)))))/(dx^2);

%error
Error_left_sided_skew = abs(left_sided_skew- analytical_derivative);
end

OUTPUT:

CONCLUSION:

As seen from the plot, the central differencing scheme has very less error compared to other schemes. Therefore, central differencing is better than other schemes because it uses information from both the sides from the point of consideration to give approximative and accurate values.

As the left and right skewed difference methods uses data only from one side of point of interest it gives less approximative and accurate values about the point.

WHY A SKEWED SCHEME IS USEFUL:

The system which doesn’t contain sufficient information on both side of the point of interest. In such cases we cannot use central differencing scheme so the skewed scheme is useful in such cases.

WHAT CAN A SKEWED SCHEME DO THAT A CD SCHEME CANNOT DO:

CDS can not be employed at the boundary nodes due to the unavailability of right and left nodes at the same time. In such cases, the skewed differencing scheme is used.