Week 3 - Taylor table method and Matlab code

1. AIM:

• Derive the 4th order approximations of the Central difference, skewed right-sided difference, skewed left-sided difference of the second-order derivative.
• To write a program in MATLAB to evaluate the second-order derivative of the analytical function and compare it with the Central difference, skewed right-sided difference, skewed left-sided difference of the second-order derivative.
• Plot the graph that compares the absolute error between the Central difference, skewed right-sided difference, skewed left-sided difference of the second-order derivative.

2. OBJECTIVE:

• To Compute the derivative both numerically and analytically.
• Plot a graph that compare between absolute error between the Central difference, skewed right-sided difference, skewed left-sided difference of the second-order derivative.

3. BODY OF THE CONTENT:

• To write a program in MATLAB to evaluate of the Central difference, skewed right-sided difference, skewed left-sided difference of the second-order derivative, we want to first derive the 4th order approximations of the Central difference, skewed right-sided difference, skewed left-sided difference of the second-order derivative using the ‘Taylor table method’.
• Now we derive the approximations following

SECOND ORDER DERIVATIVE:

Here we have to decide that how many nodes we require to find the 4th order approximation, we have the formula i.e. p+q in which p is the order of the derivative & q is the order of approximation by using this we will get the nodes. By using this formula p+q we got nodes for the skewed left & right sided difference scheme. For central difference we are using the formula p+q-1 & got nodes for the central difference scheme.

1. Central difference Scheme = p+q-1 = 2+4-1 = 5
2. Skewed right sided scheme = p+q = 2+4 = 6
3. Skewed left sided scheme = p+q = 2+4 = 6

Firstly, set up the stencil for the scheme.

• Skewed right sided scheme (i, i+1, i+2, i+3, i+4, i+5)
• Skewed left sided scheme (i, i-1, i-2, i-3, i-4, i-5)
• Central difference Scheme (i-2, i-1, i, i+1, i+2)

Now, using Taylor series expansion, we can write

`f(i+1)=f(i)+f^' (i)dx+f^'' (x)((dx^2)/2!)+f^''' (i)((dx^3)/3!)+ f^'''' (i)((dx^4)/4!)+f'''''(i)((dx^5)/5!)`

`f(i-1)=f(i)-f^' (i)dx+f^'' (x)((dx^2)/2!)-f^''' (i)((dx^3)/3!)+ f^'''' (i)((dx^4)/4!)-f'''''(i)((dx^5)/5!)`

`f(i+2)=f(i)+f^' (i)2dx+f^'' (x)((4dx^2)/2!)+f'''(i)((8dx^3)/3!)+ f^'''' (i)((16dx^4)/4!)+f'''''(i)((32dx^5)/5!)`

`f(i-2)=f(i)-f^' (i)2dx+f^'' (x)((4dx^2)/2!)-f^''' (i)((8dx^3)/3!)+ f^'''' (i)((16dx^4)/4!)-f'''''(i)((32dx^5)/5!)`

`f(i+3)=f(i)+f^' (i)3dx+f^'' (x)((9dx^2)/2!)+f'''(i)((27dx^3)/3!)+ f^'''' (i)((81dx^4)/4!)+f'''''(i)((243dx^5)/5!)`

`f(i-3)=f(i)-f^' (i)3dx+f^'' (x)((9dx^2)/2!)-f^''' (i)((27dx^3)/3!)+ f^'''' (i)((81dx^4)/4!)-f'''''(i)((243dx^5)/5!)`

`f(i+4)=f(i)+f^' (i)4dx+f^'' (x)((16dx^2)/2!)+f'''(i)((64dx^3)/3!)+ f^'''' (i)((256dx^4)/4!)+f'''''(i)((1024dx^5)/5!)`

`f(i-4)=f(i)-f^' (i)4dx+f^'' (x)((16dx^2)/2!)-f^''' (i)((64dx^3)/3!)+ f^'''' (i)((256dx^4)/4!)-f'''''(i)((1024dx^5)/5!)`

`f(i+5)=f(i)+f^' (i)5dx+f^'' (x)((25dx^2)/2!)+f'''(i)((125dx^3)/3!)+ f^'''' (i)((625dx^4)/4!)+f'''''(i)((3125dx^5)/5!)`

`f(i+5)=f(i)-f^' (i)5dx+f^'' (x)((25dx^2)/2!)-f^''' (i)((125dx^3)/3!)+ f^'''' (i)((625dx^4)/4!)-f'''''(i)((3125dx^5)/5!)`

• Skewed right sided scheme (i, i+1, i+2, i+3, i+4, i+5)

`(d^2 f)/(dx^2 )=((af(i)+bf(i+1)+cf(i+2)+df(i+3)+ef(i+4)+gf(i+5)))/((dx^2))+0(dx^4)`

Using Taylor table method to find the values of a, b, c, d, e, g.

	f(i)	f’(i)dx	f’’(i)dx^2	F’’’(i)dx^3	f’’’’(i)dx^4	f’’’’’(i)dx^5
a f(i)	a	0	0	0	0	0
b f(i+1)	b	b	b/2	b/6	b/24	b/120
c f(i+2)	c	2c	4c/2	8c/6	16c/24	32c/120
d f(i+3)	d	3d	9d/2	27d/6	81d/24	243d/120
e f(i+4)	e	4e	16e/2	64e/6	256e/24	1024e/120
g f(i+5)	g	5g	25g/2	125g/6	625g/24	3125g/120

Representing the above table in matrix form and we get the values of a, b, c, d, e, g.When we Solving the equation using MATLAB and we get

a = 3.7500

b = -12.8333

c = 17.8333

d = -13.0000

e = 5.0833

g = -0.0833

put the values in equation

`(d^2 f)/(dx^2 )=(((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.0833)f(i+5)))/((dx^2)) ……. Equation 1`

The above equation is used in MATLAB to solve the Skewed right sided 4^th order approximation of second order derivative.

• Skewed left sided scheme (i-5, i-4, i-3, i-2, i-1, i)

`(d^2 f)/(dx^2 )=((af(i-5)+bf(i-4)+cf(i-3)+df(i-2)+ef(i-1)+gf(i)))/((dx^2))+0(dx^4)`

Using Taylor table method to find the values of a, b, c, d, e, g.

	f(i)	f’(i)dx	f’’(i)dx^2	F’’’(i)dx^3	f’’’’(i)dx^4	f’’’’’(i)dx^5
a f(i-5)	a	-5a	25a/2	-125a/6	625a/24	-3125a/120
b f(i-4)	b	-4b	16b/2	-64b/6	256b/24	-1024b/120
c f(i-3)	c	-3c	9c/2	-27c/6	81c/24	-243c/120
d f(i-2)	d	-2d	4d/2	-8d/6	16d/24	-32d/120
e f(i-1)	e	-e	e/2	-e/6	e/24	-e/120
g f(i)	g	0	0	0	0	0

Representing the above table in matrix form and we get the values of a, b, c, d, e, g.When we Solving the equation using MATLAB and we get

a = -0.0833

b = 5.0833

c = -13.0000

d = 17.8333

e = -12.8333

g = 3.7500

put the values in equation

`(d^2 f)/(dx^2 )=(((-0.0833 )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)))/((dx^2)) ……. Equation 2`

The above equation is used in MATLAB to solve the Skewed left sided 4^th order approximation of second order derivative.

• Central difference Scheme (i-2, i-1, i, i+1, i+2)

`(d^2 f)/(dx^2 )=((af(i-2)+bf(i-1)+cf(i)+df(i+1)+ef(i+2)))/((dx^2))+0(dx^4)`

Using Taylor table method to find the values of a, b, c, d, e.

	f(i)	f’(i)dx	f’’(i)dx^2	F’’’(i)dx^3	f’’’’(i)dx^4
a f(i-2)	a	-2a	4a/2	-(8a/6)	16a/24
b f(i-1)	b	-b	b/2	-(b/6)	b/24
c f(i)	c	0	0	0	0
d f(i+1)	d	d	d/2	d/6	d/24
e f(i+1)	e	2e	b/2	8e/6	16e/24

Representing the above table in matrix form and we get the values of a, b, c, d, e.When we Solving the equation using MATLAB and we get

a = -0.0833

b = 1.3333

c = 2.5000

d = 1.3333

e = -0.0833

put the values in equation

`(d^2 f)/(dx^2 )=(((-0.0833)f(i-2)+( 1.3333)f(i-1)+(-2.5000)f(i)+( 1.3333)f(i+1)+(-0.0833)f(i+2)))/((dx^2)) ……. Equation 3`

The above equation is used in MATLAB to solve the Central difference 4^th order approximation of second order derivative.

Now, we got the equation for central difference, skewed right sided & skewed left sided for programming prospective. We are using these equations to solve the function numerically in MATLAB.

Program Syntax:

clear all 
close all
clc


% function = exp(x)*cos(x) ;
% analytical second order derivative 
%f'(x) = ((exp(x))*(-sin(x))) + ((cos(x))*(exp(x)));
%f''(x)= -2*exp(x)*sin(x)


% Skewed right sided scheme 4th order approximation of second order derivative
s = [1 1 1 1 1 1;0 1 2 3 4 5;0 1/2 4/2 9/2 16/2 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]
d = [0;0;1;0;0;0]
h = mldivide(s,d)
 
% Skewed left sided scheme 4th order approximation of second order derivative 
b = [1 1 1 1 1 1;-5 -4 -3 -2 -1 0;25/2 16/2 9/2 4/2 1/2 0;-125/6 -64/6 -27/6 -8/6  -1/6 0;625/24 256/24 81/24 16/24 1/24 0;-3125/120 -1024/120 -243/120 -32/120 -1/120 0]
n = [0;0;1;0;0;0]
k = mldivide(b,n)

% central difference 4th order approximation of second order derivative
y = [1 1 1 1 1;-2 -1 0  1 2;2 1/2 0 1/2 2;-8/6 -1/6 0 1/6 8/6;16/24 1/24 0 1/24 16/24]
z = [0;0;1;0;0]
a = mldivide(y,z)

x =pi/3;
analytical_second_derivative =  -(2*(exp(x)*sin(x)));
% numerical second order derivative
% Skewed right sided scheme = (d^2 f)/(dx^2 )=((af(i)+bf(i+1)+cf(i+2)+df(i+3)+ef(i+4)+gf(i+5)))/((dx^2))+0(dx^4)
% Skewed left sided scheme = (d^2f)/(dx^2 )=((af(i-5)+bf(i-4)+cf(i-3)+df(i-2)+ef(i-1)+gf(i)))/((dx^2))+0(dx^4)
% central_difference_scheme = (d^2 f)/(dx^2 )=((af(i-2)+bf(i-1)+cf(i)+df(i+1)+ef(i+2)))/((dx^2))+0(dx^4)

dx = linspace(10^-2,1,50);

for i = 1:length(dx)


% Skewed right sided scheme Fourth order approximatiom
Skewed_right_sided_scheme(i) = ((h(1)*(exp(x)*cos(x)) + h(2)*(exp(x+(dx(i)))*cos(x+(dx(i)))) + h(3)*(exp(x+(2*dx(i)))*cos(x+2*dx(i)))) + h(4)*(exp(x+3*(dx(i)))*cos(x+3*(dx(i)))) + h(5)*(exp(x+(4*dx(i)))*cos(x+(4*dx(i)))) + h(6)*(exp(x+(5*dx(i)))*cos(x+(5*dx(i)))))/(dx(i)).^2;

% Skewed left sided scheme Fourth order approximation
Skewed_left_sided_scheme(i)=  ((k(1)*(exp(x-(5*dx(i)))*cos(x-(5*dx(i)))) + k(2)*(exp(x-(4*dx(i)))*cos(x-(4*dx(i)))) + k(3)*(exp(x-(3*dx(i)))*cos(x-(3*dx(i)))) + k(4)*(exp(x-(2*dx(i)))*cos(x-(2*dx(i)))) + k(5)*(exp(x-(dx(i)))*cos(x-(dx(i)))) +(k(6)*(exp(x)*cos(x)))))/(dx(i)).^2;

% central difference scheme Fourth order approximatiom
central_difference_scheme(i) = (a(1)*(exp(x-(2*dx(i)))*cos(x-(2*dx(i)))) + a(2)*(exp(x-(dx(i)))*cos(x-(dx(i)))) + a(3)*(exp(x)*cos(x)) + a(4)*(exp(x+(dx(i)))*cos(x+(dx(i)))) + a(5)*(exp(x+(2*dx(i)))*cos(x+(2*dx(i)))))/(dx(i)).^2;

end

Skewed_right_sided_scheme_error  = abs(Skewed_right_sided_scheme - analytical_second_derivative)
Skewed_left_sided_scheme_error = abs(Skewed_left_sided_scheme - analytical_second_derivative)
central_difference_scheme_error  = abs(central_difference_scheme - analytical_second_derivative)
 
figure(1)
loglog(dx,Skewed_right_sided_scheme_error,'linewidth',2) % loglog
hold on
loglog(dx,Skewed_left_sided_scheme_error,'r','linewidth',2) % loglog
hold on
loglog(dx,central_difference_scheme_error,'black','linewidth',2) % loglog
hold off
grid on
xlabel('dx')
ylabel('error')
legend('Skewed right sided error','Skewed left sided error','central difference scheme error') 

Syntax Explanation:

• We want to compute the derivative both analytically and numerically.
• Our computing function is and it is at line no. 4.
• We calculate the value of a, b, c, d, e & g for the right sided and left side scheme equation in the form of ’h’ and ‘k’ and the value a, b, c, d & e for central difference equation in the form of ‘a’. It shows in line no. 8 to 19.
• Compute the second order derivation at x=pi/3 and it is represented in line no. 20.
• We solve the analytical second order derivative in line no.21 and get below output screenshot

• Above the program syntax we derive the formula for skewed right sided, skewed left sided and central difference scheme of 4^th order approximation of second order derivative in the form of equation 1, 2 & 3.
• We are using this formula or equations to solve the function numerically.
• The formula for skewed right sided, skewed left sided and central difference scheme in line no.23, 24 & 25 respectively.
• We define the range of step size (dx) for the computing the function numerically using ‘linspace’ command in line no. 26 for first order derivative.
• y = linsapce (x1, x2, n)
• The linspace command generates ’n’ equally spaced elements between numbers x1 and x2 for y.
• Linspace command does same the above line for programming prospective.
• The line no. 26 shows ‘dx = linspace(10^-2,1,50)’ and we get below output screenshot for ‘dx’

• A ‘loop’, in general, refers to anything whose start and end points are connected. That is, once you begin from the start point and proceed further, you will reach the finish line just to find that it once again guides you to the same point where you started. Proceeding further means that you are repeating your same action again & again.
• Now, in programming point of view, a ‘loop’ is a set of instructions or commands that repeats itself for a specified number of times.
• We are using ‘i = 1: length(dx)’ for calculating the 4th order approximation of skewed right sided, skewed left sided and central difference scheme for second order derivative of each value of dx.
• We are using ‘for loop’ for getting multiple values of dx to calculate 4th order approximation of skewed right sided, skewed left sided and central difference scheme for second order derivative in line no. 29, 31 & 33 respectively and please find below output screenshots

Skewed right sided scheme

Skewed left sided scheme

Central difference scheme

• And we used the end key for ending the loop.
• After that we want to get the error for skewed right sided, skewed left sided and central difference scheme for second order derivative, we calculate the absolute difference in numerical and exact derivative for multiple values of dx in line no. 35, 36 & 37 respectively and we get the below output screenshots of the error

Skewed right sided error

Skewed left sided error

Central difference error

• When we got the skewed right sided, skewed left sided and central difference scheme error then we plot the graph to compare the step size(dx), skewed right sided, skewed left sided and central difference scheme error and we get below output screenshot of graph

Compares the error with FDS, BDS & CDS:

• From the above graph when we look carefully the error in Skewed right sided scheme maximum when the value of dx is greater but as the value of dx decreases the value of error also decreases and slightly increases then come approximately constant.
• In Skewed left sided scheme, the error maximum when the value of dx is greater but the value of dx decreases the value of error also decreases after that it again increases and come approximately constant.
• But, when we compare the FDS with BDS, we get the error in FDS is greater than BDS Because in FDS we are taking points from right side and in BDS from left side.
• In central difference the value of error is approximately constant at every value of dx, because in central difference we are taking only 2 neighbouring point from left and right side and it shows less error comparing with FDS & BDS. It is too accurate than FDS & BDS.

Describe why a skewed scheme is useful? What can a skewed scheme do that a CD scheme cannot do?

• Skewed is useful because we calculate either left or right sided point. If we want take more points and solve FDS or BDS we solve but in CDS we have to take same points from left and right considering current point.
• When we solve FDS or BDS we only take points from right or left, but advantage is that You can take many points you want.
• In CDS you have to takes point from both left and right.

4. Conclusion:

• We derived skewed right sided, skewed left sided and central difference scheme for fourth-order accurate of second order derivative.
• Hence, we had written a program to compute the derivative both numerically and analytically for skewed right sided, skewed left sided and central difference scheme for second order derivative.
• Also, we had calculated the error for absolute difference in numerical and exact derivative for multiple value for step size (dx).
• we plotted the graph to compare the step size(dx), skewed right sided, skewed left sided and central difference scheme error.