To find out 4th order Approximation error for second order derivative by use of CD SRS SLS

Objective

Derive the following 4th order approximations of the second order derivative

1. Central difference

2. Skewed right sided difference

3. Skewed left sided difference

4. plot the graph for diffrence of error for above scheme

Derivation of central diffrence

For derive 4th order approximation of second order derivative.folllowing node are consider for central diffrence method. nodes are i,(i-2),(i-1),(i+2),(i+1). the value of af(i-2),b(i-1)..etc find out by taylor series expansion. like this there is use of taylor series expansion for skewed left & right 4th order approx.

`af(i-2)= a[f_(i)-(2Δxf\'(i))/(1!)+f_(i)\'\'((2Δx)^2)/(2!) + f_(i)\'\'\'((2Δx)^2)/(2!) ] +......`

Taylor table for find out 4th order approximation

stensil	f(i)	dxf\'(i)	dx^2f\'\'(i)	dx^3f\'\'\'(i)	dx^4f^iv(i)	dx^5f^v(i)	dx^6f^vi(i)
af(i-2)	a	-2a	4a/2	-8a/2	16a/24	-32a/120	64a/720
bf(i-1)	b	-b	b/2	-b/6	b/24	-b/120	b/240
cf(i)	c	0	0	0	0	0	0
df(i+1)	d	d	d/2	d/6	d/6	-d/120	d/720
ef(i+2)	e	2e	4e/2	8e/6	8e/6	-32e/120	64e/720
sum	0	0	1	0	0	?	?

after that made 5 equation for finding the value a,b,c,d,e

equation1: a+b+c+d+e=0

equation2: (-2a)-b+0+d+2e=0

equation3: 2a+b/2+0+d/2+4e/2=1

equaion4: (-4a)+(-b/6)+0+d/6+8e/6=0

equation5: (16a/24)+b/24+0+d/6+8e/6=0

the value of a,b,c,d,e is find out by use of matrix in matlab

 A=[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];
B=[0;0;1;0;0];
x=A\\B;
disp(x)

command window shows the value for coefficent

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

after putting the value in in main equation we got final equation of 4th order approx.

`(d^2f)/dx^2= (-0.8333f_(i-2)+1.333f_(i-1)-2.5f_(i)+1.33f_(i+1)-0.0833f_(i+2))/dx^2`

Derivation for Skewed right side approximation

For derive 4th order approximation of second order derivative.folllowing node are consider for skewed right side approx. nodes are i,(i+1),(i+2),(i+3),(i+4).

	f(i)	dxf\'(i)	dx^2f\'\'(i)	dx^3f\'\'\'(i)	dx^4f^iv(i)	dx^5f^v(i)	dx^6f^vi(i)
af(i)	a	0	0	0	0	0	0
bf(i+1)	b	b	b/2	b/6	b/24	b/120	b/720
cf(i+2)	c	2c	4c/2	8c/2	16c/24	32c/120	64c/720
df(i+3)	d	3d	9d/2	27d/6	81d/24	243d/120	729d/720
ef(i+4)	e	4e	16e/2	64e/6	256e/24	1024e/120	4096e/720
sum	0	0	1	0	0	?	?

after that made 5 equation for finding the value a,b,c,d,e

equation1: a+b+c+d+e=0

equation2: 0+b+2c+3d+4e=0

equation3: 0+b/2+4c/2+9d/2+16e/2=1

equaion4: 0+(b/6)+(8c/2)+27d/6+64e/6=0

equation5: 0+(b/24)+16c/24+81d/24+256e/24=0

a,b,c,d,e value find out by use of matrix in matlab

A=[1       1         1     1           1;
    0       1         2     3           4;
    0       1/2       2     9/2         8;
    0       1/6      8/6     27/6       64/6;
    0       1/24      16/24  81/24      256/24];
B=[0;0;1;0;0];
x=A\\B;
disp(x)

command window shows the value for coefficent

a=2.9167 b=-8.6667 , c=9.5000 , d=-4.6667 , e=0.9167

after putting the value in in main equation we got final equation of 4th order approx.

`(d^2f)/dx^2= (2.9167f_(i)-8.6667f_(i+1)+9.5f_(i+2)-4.6667f_(i+3)+0.9167f_(i+4))/dx^2`

Derivation for Skewed left side approximation

For derive 4th order approximation of second order derivative.folllowing node are consider for skewed left side approx. nodes are i,(i-1),(i-2),(i-3),(i-4).

stensil	f(i)	dxf\'(i)	dx^2f\'\'(i)	dx^3f\'\'\'(i)	dx^4f^iv(i)	dx^5f^v(i)	dx^6f^vi(i)
af(i-4)	a	4a	16a/2	-64a/6	256a/24	-1024a/120	4096a/720
bf(i-3)	b	-3b	9b/2	-27b/6	81b/24	-243b/120	729b/720
cf(i-2)	c	-2c	4c/2	-8c/6	16c/24	-32c/120	64c/720
df(i-1)	d	-d	d/2	-d/6	d/24	-d/120	d/720
ef(i)	e	0	0	0	0	0	0
sum	0	0	1	0	0	?	?

after that made 5 equation for finding the value a,b,c,d,e

equation1: a+b+c+d+e=0

equation2: 4a-3b-2c-d+0=0

equation3: 16a/2+9b/2+4c/2+d/2+0=1

equaion4: (-64a/6)+(-27b/6)+(-8c/6)+(-d/6)+0=0

equation5: (256a/24)+81b/24+16c/24+d/24+0=0

A=[1       1         1       1           1;
    -4     -3        -2      -1           0;
     8      9/2      4/2     1/2          0;
    -64/6   -27/6   -8/6     -1/6          0;
    256/24   81/24   16/24   1/24         0];
B=[0;0;1;0;0];
x=A\\B;
disp(x)

command window shows the value for coefficent.

a=2.9167 , b=-8.6667 , c=9.5000 , d=-4.6667 ,e=0.9167

after put the value in main equation we got final equation of 4th order approx.

`(d^2f)/dx^2= (2.9167f_(i)-8.6667f_(i-1)+9.5f_(i-2)-4.6667f_(i-3)+0.9167f_(i-4))/dx^2_`

Function out code for central diffrence

function out= central_fourth_order_approximation(x,dx)
% of analyatical function=exp(x)*cos(x)
% of analyatical derivative
% f\'\'(x)=-2*e(x)*sin(x)
second_order_analytical_derivative= -2*exp(x)*sin(x);
%numerical derivative
% use of cenrtal diffrencing
% fourth order approximation
central_diffrence=(((-0.08333)*((exp(x-(2*dx)))*(cos(x-(2*dx)))))+((1.3333)*((exp(x-dx))*(cos(x-dx))))-((2.5)*((exp(x))*(cos(x))))+((1.3333)*((exp(x+dx))*(cos(x+dx))))+((-0.0833)*((exp(x+(2*dx)))*(cos(x+(2*dx))))))/(dx^2);
out=abs(central_diffrence - second_order_analytical_derivative);
end

Function out code for skewed right side

function out= skewed_right_side_approximation(x,dx)
% of analyatical function=exp(x)*cos(x)
% of analyatical derivative
% f\'\'(x)=-2*e(x)*sin(x)
second_order_analytical_derivative= -2*exp(x)*sin(x);
%numerical derivative
% use of skewed right side 
% fourth order approximation
skewed_right_side =(((2.9167)*((exp(x))*(cos(x))))-((8.6667)*((exp(x+dx))*(cos(x+dx))))+((9.5)*((exp(x+(2*dx))*(cos(x+(2*dx))))))-((4.6667)*((exp(x+(3*dx)))*(cos(x+(3*dx)))))+((0.9167)*((exp(x+(4*dx)))*(cos(x+(4*dx))))))/(dx^2);
out=abs(skewed_right_side- second_order_analytical_derivative);
end 

Function out code for skewed left side

function out= skewed_left_side_approximation(x,dx)
% of analyatical function=exp(x)*cos(x)
% of analyatical derivative
% f\'\'(x)=-2*e(x)*sin(x)
second_order_analytical_derivative= -2*exp(x)*sin(x);
%numerical derivative
% use of skewed left side 
% fourth order approximation

skewed_left_side =(((0.9167)*((exp(x-(4*dx))*(cos(x-(4*dx))))))-((4.6667)*((exp(x-(3*dx))*(cos(x-(3*dx))))))+((9.5)*((exp(x-(2*dx))*(cos(x-(2*dx))))))-((8.6667)*((exp(x-dx))*(cos(x-dx))))+((2.9167)*((exp(x))*(cos(x)))))/(dx^2);
out=abs(skewed_left_side- second_order_analytical_derivative);
end

Main code

clear all
close all
clc
x=pi/2;
dx=linspace(pi/4,pi/40,10);
 second_order_analytical_derivative=-2*sin(x)*exp(x);

for i=1:length(dx)
    central_diffrence_error(i)= central_fourth_order_approximation(x,dx(i));
    skewd_right_side_error(i)=skewed_right_side_approximation(x,dx(i));
    skewd_left_side_error(i)=skewed_left_side_approximation(x,dx(i));
end

figure(1)
loglog(dx,central_diffrence_error,\'b\')
hold on
loglog(dx,skewd_right_side_error,\'r\')
hold on
loglog(dx,skewd_left_side_error,\'g\')
legend(\'central_diffrence_error\',\'skewd_right_side_error\',\'skewd_left_side_error\')
xlabel(\'dx value\')
ylabel(\'error value\')

Plot for comparison of absolute error between CDS, SRS,SLS

result is found at x=pi/2 & dx=linspace(pi/4,pi/40,10). so from graph the it is found that if we use central diffrence method for 4th order approximation there is error is less one compare with the other two method of skewed left & right side method. the central diffrence error gives less error because it take the data from symeterical node point from left & side right side.

why skewed system is useful?

in case of skewed for right side & left side system. it take inforamtion from current point & neighbouring point on left & right side so we can find out any order approximation. but in case of Cental diffrence system if we want to find out higher order approximation at boundry node, in that case there is no node is available on any side of node. , so it is not useful for Cental diffrence system.