Numerical discretization of Second order derivative up to fourth order approximation using Taylor table method and Matlab code
1. Central/symmetrical differencing scheme:
Scope:
To determine second order derivative value and compare by numerical approximation (up to fourth order) using Taylor’s table method.
Procedure:
Total No. of grip points considered for central differencing scheme = p + q-1=5
Where ‘p’- Order of derivative
‘q’- Order of approximation
Equation of Second order derivative by numerical discretization (fourth_order_approximation)
Figure.1 Representative diagram of central differencing consideration
Taylor’s Table:
|
|
f(i) |
dx*f’(i) |
dx2*f’’(i) |
dx3*f’’’(i) |
dx4*fIV(i) |
|
f(i-2) |
a |
-2a |
2a |
-4a/3 |
2a/3 |
|
f(i-1) |
b |
-b |
b/2 |
-b/24 |
-b/120 |
|
f(i) |
c |
0 |
0 |
0 |
0 |
|
f(i+1) |
d |
d |
d/2 |
d/6 |
d/24 |
|
f(i+2) |
e |
2e |
2e |
4e/3 |
2e/3 |
Considering & solving the unknowns by matrix form:
a=-0.0833
b=1.3333
c=-2.5
d=1.3333
e=-0.0833
Second order derivative equation:
2. Skewed right side differencing scheme:
Scope:
To determine second order derivative value and compare by numerical approximation (up to fourth order) using Taylor’s table method.
Procedure:
Total No. of grid points considered for skewed right side differencing scheme = p + q =6
Where ‘p’- Order of derivative
‘q’- Order of approximation
Equation of Second order derivative by numerical discretization(fourth_order_approximation)
Figure.2 Representative diagram of skewed right side differencing consideration
Taylor’s Table:
|
|
f(i) |
dx*f’(i) |
dx2*f’’(i) |
dx3*f’’’(i) |
dx4*fIV(i) |
dx5*fIV(i) |
|
f(i) |
a |
0 |
0 |
0 |
0 |
0 |
|
f(i+1) |
b |
b |
b/2 |
b/6 |
b/24 |
b/120 |
|
f(i+2) |
c |
2c |
4c/2 |
8c/6 |
16c/24 |
32c/120 |
|
f(i+3) |
d |
3d |
9d/2 |
27d/6 |
81d/24 |
243d/120 |
|
f(i+4) |
e |
4e |
16e/2 |
64e/6 |
256e/24 |
1024e/120 |
|
f(i+5) |
g |
5g |
25g/2 |
125g/6 |
625g/24 |
3125g/120 |
Considering & solving the unknowns by matrix form:
a= 3.75
b= -12.833
c= 17.833
d= -13.00
e= 5.0833
g= -0.833
Second order derivative:
3. Skewed left side differencing scheme:
Scope:
To determine second order derivative value and compare by numerical approximation (up to fourth order) using Taylor’s table method.
Procedure:
Total No. of points considered for skewed left side differencing scheme = p + q =6
Where ‘p’- Order of derivative
‘q’- Order of approximation
Equation of Second order derivative by numerical discretization(fourth_order_approximation)
Figure.2 Representative diagram of skewed left side differencing consideration
Taylor’s Table:
|
|
f(i) |
Dx*f’(i) |
Dx2*f’’(i) |
Dx3*f’’’(i) |
Dx4*fIV(i) |
Dx5*fIV(i) |
|
f(i) |
a |
0 |
0 |
0 |
0 |
0 |
|
f(i-1) |
b |
-b |
b/2 |
-b/6 |
b/24 |
-b/120 |
|
f(i-2) |
c |
-2c |
4c/2 |
-8c/6 |
16c/24 |
-32c/120 |
|
f(i-3) |
d |
-3d |
9d/2 |
-27d/6 |
81d/24 |
-243d/120 |
|
f(i-4) |
e |
-4e |
16e/2 |
-64e/6 |
256e/24 |
-1024e/120 |
|
f(i-5) |
g |
-5g |
25g/2 |
-125g/6 |
625g/24 |
-3125g/120 |
Considering & solving the unknowns by matrix form:
a= 3.75
b= -12.833
c= 17.833
d= -13.00
e= 5.0833
g= -0.833
Second order derivative:
4. Matlab Program:
Main Program:
- clear all
- close all
- clc
- %function f(i)=exp(x)*cos(x)
- %analytical_second_order_derivative
- %d2u/dx2=-2*(exp(x)*sin(x))
- x=2*pi
- dx=linspace(pi/4,pi/4000,20)
- analytical_second_order_derivative=-2*((exp(x))*sin(x))
- %Numerical_discritization(2nd_order_derivative&fourth_order_approximation)
- for i=1:length(dx)
- central_differencing_error(i)=central_differencing_fourth_order_approximation(x,dx(i));
- skewed_right_side_error(i)=skewed_right_side_fourth_order_approximation(x,dx(i));
- skewed_left_side_error(i)=skewed_left_side_fourth_order_approximation(x,dx(i));
- end
- figure(1)
- loglog(dx,central_differencing_error)
- hold on
- loglog(dx,skewed_right_side_error,\'black\')
- hold on
- loglog(dx,skewed_left_side_error,\'red\')
- legend(\'central_differencing_error\',\'skewed_right_side_error\',\'skewed_left_side_error\')
- xlabel(\'grid_spacing\')
- ylabel(\'fourth_order_approximation\')
- figure(2)
- plot(dx,central_differencing_error)
- hold on
- plot(dx,skewed_right_side_error,\'black\')
- hold on
- plot(dx,skewed_left_side_error,\'red\')
- legend(\'central_differencing_error\',\'skewed_right_side_error\',\'skewed_left_side_error\')
- xlabel(\'grid_spacing\')
- ylabel(\'fourth_order_approximation\')
Function Programs:
Central differencing:
function central_differencing_error = central_differencing_fourth_order_approximation(x,dx)
%function f(i)=(exp(x)*cos(x))
%analytical_second_order_derivative
%d2u/dx2=-2*(exp(x)*sin(x))
analytical_second_order_derivative=-2*(exp(x)*sin(x));
%Numerical_fourth_order_approximation
%system_of_equations_from_taylor\'s_table
%a+b+c+d+e=0
%-2a-b+d+2e=0
%4a/2-b/6+d/6-4e/3=0
%2a/3+b/24+d/24+2e=0
x1=[ 1 1 1 1 1
-2 -1 0 1 2
2 1/2 0 1/2 2
-4/3 -1/6 0 1/6 4/3
2/3 1/24 0 1/24 2/3 ]
x2=[0 0 1 0 0]\'
x3=inv(x1)*x2
%f(0i)=(exp(x)*cos(x))
%f(1i+1)=((exp(x+dx))*cos(x+dx))
%f(2i+2)=((exp(x+(2*dx)))*cos(x+(2*dx)))
%f(1i-1)=((exp(x-dx))*cos(x-dx))
%f(2i-2)=((exp(x-(2*dx)))*cos(x-(2*dx)))
central_differencing_fourth_order_approximation=(((x3(1))*((exp(x-(2*dx)))*cos(x-(2*dx))))+(x3(2)*((exp(x-dx))*cos(x-dx)))+((x3(3))*((exp(x))*cos(x)))+(x3(4)*((exp(x+dx))*cos(x+dx)))+((x3(5))*((exp(x+(2*dx)))*cos(x+(2*dx)))))/dx^2;
central_differencing_error=abs(analytical_second_order_derivative-central_differencing_fourth_order_approximation)
end
Skewed right side differencing:
function skewed_right_side_error = skewed_right_side_fourth_order_approximation(x,dx)
%function f(i)=(exp(x))*cos(x))
%analytical_second_order_derivative
%d2u/dx2=-2*(exp(x))*sin(x))
analytical_second_order_derivative=-2*((exp(x))*sin(x));
%Numerical_fourth_order_approximation
%system_of_equations_from_taylor\'s_table
%a+b+c+d+e+g=0
%b+2c+3d+4e+5g=0
%b/2+4c/2+9d/2+16e/2+25g/2=1
%b/6+8c/6+27d/6+64e/6+125g/6=0
%b/24+16c/24+81d/24+256e/24+625g/24=0
%b/120+32c/120+243d/120+1024e/120+3125g/120=0
x1=[ 1 1 1 1 1 1
0 1 2 3 4 5
0 1/2 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 ]
x2=[0 0 1 0 0 0]\'
x3=inv(x1)*x2
skewed_right_side_fourth_order_approximation=((x3(1)*((exp(x))*cos(x)))+((x3(2))*((exp(x+dx))*cos(x+dx)))+(x3(3)*((exp(x+(2*dx)))*cos(x+(2*dx))))+((x3(4))*((exp(x+(3*dx)))*cos(x+(3*dx))))+(x3(5)*((exp(x+(4*dx)))*cos(x+(4*dx))))+((x3(6))*((exp(x+(5*dx)))*cos(x+(5*dx)))))/dx^2
skewed_right_side_error=abs(analytical_second_order_derivative-skewed_right_side_fourth_order_approximation)
end
Skewed left side differencing:
function skewed_left_side_error = skewed_left_side_fourth_order_approximation(x,dx)
%function f(i)=(e^x)*cos(x)
%analytical_second_order_derivative
%d2u/dx2=-2*(e^x)*sin(x)
analytical_second_order_derivative=-2*((exp(x))*sin(x))
%Numerical_fourth_order_approximation
%system_of_equations_from_taylor\'s_table
%a+b+c+d+e+g=0
%-b-2c-3d-4e-5g=0
%b/2+4c/2+9d/2+16e/2+25g/2=1
%-b/6-8c/6-27d/6-64e/6-125g/6=0
%b/24+16c/24+81d/24+256e/24+625g/24=0
%-b/120-32c/120-243d/120-1024e/120-3125g/120=0
x1=[ 1 1 1 1 1 1
0 -1 -2 -3 -4 -5
0 1/2 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]
x2=[0 0 1 0 0 0]\'
x3=inv(x1)*x2
skewed_left_side_fourth_order_approximation=((x3(1)*((exp(x))*cos(x)))+((x3(2))*((exp(x-dx))*cos(x-dx)))+(x3(3)*((exp(x-(2*dx)))*cos(x-(2*dx))))+((x3(4))*((exp(x-(3*dx)))*cos(x-(3*dx))))+(x3(5)*((exp(x-(4*dx)))*cos(x-(4*dx))))+((x3(6))*((exp(x-(5*dx)))*cos(x-(5*dx)))))/dx^2
skewed_left_side_error=abs(analytical_second_order_derivative-skewed_left_side_fourth_order_approximation)
end
5. Results/Plots:
Loglog plot:
Usual plot:
Slope for various difference scheme:
Note: Kindly refer attached excel sheet also for reference.
6. Conclusion & Discussion:
On observation of the various plots obtained it is concluded that Central differencing scheme much effective scheme of discretization than any other skewed difference schemes.