TAYLOR TABLE MODELLING FOR FOURTH-ORDER-APPROXIMATION OF SECOND ORDER DERIVATIVE

DESCRIPTION:

Using the central and skewed scheme methods initially, we have to develop the Taylor table, which should be converted into a matrix and solved for the constants. These constants help in solving the second order differential equation. Finally, the error results are plotted and studied.

 

1. CENTRAL DIFFERENCING SCHEME:

 In central differencing method, we will be solving with this format of stenciling the points.

`(∂^2 f)/(∂x^2) = a.f(i-2)+b.f(i-1)+c.f(i)+d.f(i+1) +e.f(i+2)`

 

 

This table can be converted into the matrix form as:

 

After solving the matrix for constants a, b, c, d, e we should substitute the values in the table for the fifth term of the taylor series and check for the value of the sum of the fifth terms.

By checking we get,

a = -0.0833

b = 1.33

c = -2.5

d = 1.33

e = -0.0833

and hence  `((-32*a)-(b)+(d)+(32e))/120` become zero.

Similar steps have to be followed for the skewed scheme also.

2. LEFT-SIDE SKEWED SCHEME:

CONSTANTS: 

a = 3.75

b = -12.833

c = 17.833

d = -13.0

e = 5.0833

g = -0.833

 

3. RIGHT-SIDE SKEWED SCHEME:

 

 

CONSTANTS: 

a = 3.75

b = -12.833

c = 17.833

d = -13.0

e = 5.0833

g = -0.833

 

 

MAIN CODE:

clear all
close all
clc

% Taylor table

x = pi/3;
dx = linspace(pi/4,pi/400,50);

% analytical function = exp(x)*cos(x)
% f'(x) = exp(x)*cos(x) - exp(x)*sin(x)
% f''(x) = -2*exp(x)*sin(x)

%1. central diff
%2. skew right scheme
%3. skew left scheme

for n = 1: length(dx)
error1(n) = taylorfunc(x,dx(n),1);
error2(n) = taylorfunc(x,dx(n),2);
error3(n) = taylorfunc(x,dx(n),3);
error4(n) = taylorfunc(x,dx(n),4);
error5(n) = taylorfunc(x,dx(n),5);
end

figure(1)
loglog(dx,error1,'linewidth',2)
hold on
loglog(dx,error2,'linewidth',2)
loglog(dx,error3,'linewidth',2)
legend('4th order central appx','4th order left skewed appx','4th order right skewed appx','Location','southeast')
xlabel('x_value')
ylabel('error')
title('COMPARISION OF CENTRAL AND SKEWED SCHEME ERRORS') 

figure(2)
loglog(dx,error1,'linewidth',1)
hold on
loglog(dx,error2,'linewidth',1)
loglog(dx,error3,'linewidth',1)
loglog(dx,error4,'linewidth',1,'LineStyle','- -')
loglog(dx,error5,'linewidth',1,'LineStyle','- .')
legend('4th order central appx','4th order left skewed appx','4th order right skewed appx','1st order forward appx','1st order backward appx','Location','southeast')
xlabel('x_value')
ylabel('error')
title('CENTRAL,FORWARD,BACKWARD AND SKEWED SCHEME ERRORS') 

FUNCTION CODE-1:

function error = taylorfunc(x,dx,i)


if i == 1
    
%central diff scheme
%calculating constants: a, b, c, d, e

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^-1*B;

% Extracting constants
a = X(1);
b = X(2);
c = X(3);
d = X(4);
e = X(5);

central_diff = ((a*f(x-(2*dx)))+(b*f(x-dx))+(c*f(x))+(d*f(x+dx))+(e*f(x+2.*dx)))/(dx^2);
error = abs(central_diff - f2(x));

end

if i == 2
    
%left skewed scheme
%calculating constants: a, b, c, d, e
A = [1 1 1 1 1 1; 0 -1 -2 -3 -4 -5; 0 1/2 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];
B = [0; 0; 1; 0; 0; 0];
X = A^-1*B;

% Extracting constants
a = X(1);
b = X(2);
c = X(3);
d = X(4);
e = X(5);
g = X(6);

skewleft_diff = ((a*f(x))+(b*f(x-dx))+(c*f(x-2*dx))+(d*f(x-3*dx))+(e*f(x-4*dx))+(g*f(x-5*dx)))/(dx^2);

error = abs(skewleft_diff - f2(x));

end

if i == 3
    
%right skewed scheme
%calculating constants: a, b, c, d, e
A = [1 1 1 1 1 1; 0 1 2 3 4 5; 0 1/2 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];
B = [0; 0; 1; 0; 0; 0];
X = A^-1*B;

% Extracting constants
a = X(1);
b = X(2);
c = X(3);
d = X(4);
e = X(5);
g = X(6);

skewright_diff = ((a*f(x))+(b*f(x+dx))+(c*f(x+(2*dx)))+(d*f(x+(3*dx)))+(e*f(x+(4*dx)))+(g*f(x+5*dx)))/(dx^2);
error = abs(skewright_diff - f2(x));

end

if i == 4
    %forward differencing
    fd_diff = (f(x+2*dx)-(2*f(x+dx))+f(x))/(2*(dx^2));
    error = abs(fd_diff - f2(x));
end

if i == 5
    %forward differencing
    back_diff = [f(x)-(2*f(x-dx))+f(x-2*dx)]/(2*(dx^2));
    error = abs(back_diff - f2(x));
end

end

FUNCTION CODE-2:

function [f_act] = f(x)

f_act = (exp(x)*cos(x));

end

FUNCTION CODE-3:

function out = f2(x)
out = -2*exp(x)*sin(x);
end

 

RESULT:

 

1.1 ERROR COMPARISION BETWEEN SKEWED AND CENTRAL DIFFERENCING(loglog plot):

 

 

 

1.2  ERROR COMPARISION BETWEEN SKEWED AND CENTRAL DIFFERENCING(normal plot):

 

 

2. COMPARISON ON THESE ERRORS WITH FORWARD AND BACKWARD DIFFERENCING:

 

 


Projects by Yokesh R

DEMO
Yokesh R · 2020-02-04 11:48:17

Hello Read more

Demo
Yokesh R · 2020-01-21 11:46:05

Hello Read more

Mixing tee project
Yokesh R · 2019-08-22 06:06:22

Case-1:       CASE-2:           CASE-3:         Read more

MAIN CODE: clear all close all clc A=[5 1 2 ; -3 9 4; 1 2 -7]; B=[10 ;-14; 33]; u = [0 1 2; 0 0 4; 0 0 0]; l = [0 0 0; -3 0 0; 1 2 0]; di = [5 0 0; 0 9 0; 0 0 -7]; mag = [0.25,0.5,0.75,1,1.5,2]; for i = 1 : length(mag) d = di*mag(i); Read more

DESCRIPTION: In this code, we are going to solve the two dimesnional convection equation defining a temperture distribution. The equation can be solved by 2 methods. 1. Implicit method 2. Explicit method In explicit method, we will be solving the equation sequencial Read more

DESCRIPTION: In this code, we are going to compare the results of the 1D convection distribution equation describing the velocity for a range of time-step values. As a result, we will be getting the plots on the final velocity profile on time marching and computational Read more

DESCRIPTION:      In this code, we are going to discretize the function for the first-order derivative using first, second and fourth order approximation techniques and studying the error bounced by every method by a bar chart.   CODE: clear all Read more

DESCRIPTION:      Here in this code, we are using a range of dx terms for discretizing the first order derivative and studying the results of the dx vs error plot  CODE: clear all close all clc x = pi/3; dx = linspace(pi/4,pi/4000,30); %y Read more

STOICHIOMETRIC COMBUSTION   DESCRIPTION:   Stoichiometric combustion defines the ideal combustion that needs to take place for a given chemical equation. This is actually the general form of equation for alkanes. The general form might differ for e Read more


Loading...

The End