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.


Projects by deepak virkar

Title:  Simulation of quasi 1D Supersonic Nozzle flow simulation by use of Macormack method Objective; 1.To  implementing both the conservative and non-conservative forms of the governing equations. 2. To write separate functions for conservative and non-co Read more

Aim: To create an animation of pendulum by use second order diffrential equation of motion objective:1. To write a program that solves the  ODE                 2. Create an animation video by use movie(M)command. Conversio Read more

Title : Steady state analysis & Transient State Analysis Objective:   1.To Solve the 2D heat conduction equation by using the point iterative techniques for 1. Jacobi,  2. Gauss-seidel,   3. Successive over-relaxation.       Read more


Loading...

The End