Taylor table method and Matlab code

Fourth order approximations of the second order derivative using the following schemes are derived with the help of Programming

• Central difference
• Skewed right sided difference
• Skewed left sided difference

Wrote a Matlab program to evaluate the second order derivative of the analytical function exp(x)*cos(x) and compared it with the 3 numerical approximations derived.

• Plot which compares the absolute error between the above mentioned schemes is generated
• Compared with FDS and BDS approximations
• Advantanges of the skewed scheme over CD scheme is explained

1. Once after the formation of the taylor table, equations can solved in the form of matrix. Simple logical thinking can help us to write those equations, because all those equations follow a pattern, based on the scheme used. Following is the program to derive coefficients of the formula for 2nd derivative with 4th order approx for the above three schemes

clear all
close all
clc

% Equatins formed using taylor table is solved in octave in the form of matrix_type
%(del x)^2 will be the common denominator for all the coefficients of a second derivative algebraic equivalent equation
% The denominator is removed now to reduce the complexity, and it will be added in the final equation with the solved coefficients

% skewed right scheme acquires information from the rigth side points
% To obtain a 4th order approximated second derivative, we need information from 6 right side neighbouring points
a=[1 1 1 1 1 1];
b=[0 1 2 3 4 5];
c=(1/factorial(2))*[0 1 2^2 3^2 4^2 5^2];
d=(1/factorial(3))*[0 1 2^3 3^3 4^3 5^3];
e=(1/factorial(4))*[0 1 2^4 3^4 4^4 5^4];
f=(1/factorial(5))*[0 1 2^5 3^5 4^5 5^5];

p=[a;b;c;d;e;f];

% coefficients of f'' are made equal to 1
q=[0;0;1;0;0;0];
w = pq;
coefficients_of_second_derivative_with_4th_order_approx_skewed_right=rats(w)

% skewed left scheme acquires information from the left side points
% To obtain a 4th order approximated second derivative, we need information from 6 left side neighbouring points
g=[1 1 1 1 1 1];
h=[0 -1 -2 -3 -4 -5];
i=(1/factorial(2))*[0 1 (2^2) (3^2) (4^2) (5^2)];
j=(1/factorial(3))*[0 -1 -(2^3) -(3^3) -(4^3) -(5^3)];
k=(1/factorial(4))*[0 1 (2^4) (3^4) (4^4) (5^4)];
l=(1/factorial(5))*[0 -1 -(2^5) -(3^5) -(4^5) -(5^5)];

m=[g;h;i;j;k;l];

% coefficients of f'' are made equal to 1
n=[0;0;1;0;0;0];
t = mn;
coefficients_of_second_derivative_with_4th_order_approx_skewed_left=rats(t)

% central scheme requires information from points, equally from both the sides
% To obtain a 4th order approximated second derivative, we need information from 5 equally distribured neighbouring points
a_1=[1 1 1 1 1];
b_1=[-2 -1 0 1 2];
c_1=(1/factorial(2))*[(2^2) (1) (0) (1) (2^2)];
d_1=(1/factorial(3))*[-(2^3) -(1) (0) (1) (2^3)];
e_1=(1/factorial(4))*[(2^4) (1) (0) (1) (2^4)];

f_1=[a_1;b_1;c_1;d_1;e_1];

% coefficients of f'' are made equal to 1
g_1=[0;0;1;0;0];
h_1 = f_1g_1;
coefficients_of_second_derivative_with_4th_order_approx_central=rats(h_1)

% sum of the coefficients should be zero. This is a way to check our answers
check_1=w(1)+w(2)+w(3)+w(4)+w(5)+w(6);
Checking_forward=rats(check_1)
check_2=t(1)+t(2)+t(3)+t(4)+t(5)+t(6);
Checking_backward=rats(check_2)
check_3=h_1(1)+h_1(2)+h_1(3)+h_1(4)+h_1(5);
Checking_central=rats(check_3)



2. Program which results a plot that compares the absolute error between the above mentioned schemes is given below

clear all
close all
clc

%f(x) = exp(x)*cos(x);
%f''(x) = -2*exp(x)*sin(x);
x = 50 ;
dx = linspace(10^-5,10^-1,100);

funtion = exp(x)*cos(x);
exact_second_derivative = -2*exp(x)*sin(x)

for i=1:length(dx)

a(i) = (exp(x-(5*dx(i)))*cos(x-(5*dx(i))));
b(i) = (exp(x-(4*dx(i)))*cos(x-(4*dx(i))));
c(i) = (exp(x-(3*dx(i)))*cos(x-(3*dx(i))));
d(i) = (exp(x-(2*dx(i)))*cos(x-(2*dx(i))));
e(i) = (exp(x-dx(i))*cos(x-dx(i)));
f(i) = exp(x)*cos(x);
g(i) = (exp(x+dx(i))*cos(x+dx(i)));
h(i) = (exp(x+(2*dx(i)))*cos(x+(2*dx(i))));
j(i) = (exp(x+(3*dx(i)))*cos(x+(3*dx(i))));
k(i) = (exp(x+(4*dx(i)))*cos(x+(4*dx(i))));
l(i) = (exp(x+(5*dx(i)))*cos(x+(5*dx(i))));

fourth_order_second_derivative_skewed_right(i) = (1/(12*(dx(i)^2)))*(45*(f(i))-154*(g(i))+214*(h(i))-156*(j(i))+61*(k(i))-10*(l(i)))
fourth_order_second_derivative_skewed_left(i) = (1/(12*(dx(i)^2)))*(45*(f(i))-154*(e(i))+214*(d(i))-156*(c(i))+61*(b(i))-10*(a(i)))
fourth_order_second_derivative_central(i) = (1/(12*(dx(i)^2)))*(-1*(d(i))+16*(e(i))-30*(f(i))+16*(g(i))-(h(i)))

fourth_order_skewed_right_error(i) = abs(fourth_order_second_derivative_skewed_right(i) - exact_second_derivative)
fourth_order_skewed_left_error(i) = abs(fourth_order_second_derivative_skewed_left(i) - exact_second_derivative)
fourth_order_central_error(i) = abs(fourth_order_second_derivative_central(i) - exact_second_derivative)

end

% Error comparison between the above three methods
figure(1)
loglog(dx,fourth_order_skewed_right_error)
hold on
loglog(dx,fourth_order_skewed_left_error,'r')
hold on
loglog(dx,fourth_order_central_error,'g')

h=legend('fourth order skewed right error','fourth order skewed left error','fourth order central error');
legend (h, "location", "northeastoutside");
xlabel('dx');
ylabel('Error');

Output Plot with comparatively less error when CD scheme is used

How do these errors compare with FDS and BDS approximations?

FDS and BDS acquires information from the current point and the right and left neighbouring points respectively. It results in 1st order approx for 1st order derivative, with two points information, where as central difference results in 2nd order approx for 1st order derivative, with two points information. With the same amount inofrmation being available, central scheme provides more accuracy due to the reduced truncation error.

These skewed schemes helps in gaining information from many available points, where as FDS and BDS expericiences this constrain, so that the skewed schemes can result in higher order approximation. Central difference schemes are far more better than skewed schemes, because due to the concept of symmetry, it gains from both the sides which increases its ability to find more accurate derivatives at a certain point.

Explanation for the shape of the curves obtained due to the effect of grid size, is shown below.

The absolute error curve of the skewed scheme and central scheme with many neighbouring points results in a U shaped curve (approximately), w.r.t mesh size(dx). This is due to the reason that, the influence of many neighbouring points and a very small value of dx, generates a large no.of arithmetic operations, which results in round-off error. Where as for the same range of dx values, the error curves obtained by FDS and BDS are not U shaped, due to less no.of arithmetic operations

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

skewed difference scheme works well whatever be the no.of points available on the grid. For example, if there is only node available at the left and 5 nodes at the right, and we are in need of 4th order approx for second derivative, then at that time only skewed schemes work out well. If we want to compute a single or higher order derivative at a corner node, then central scheme is no longer useful.

Solving equation of a simple pendulum with damping Krishnan Jothi · 2019-01-21 14:08:34

The objective of this project is to write a Matlab/Octave program that solves the following equation and Simulate the motion between 0-20 sec                                   &nb Read more

Shock tube simulation Krishnan Jothi · 2018-09-09 13:57:51

Link to the simulation video: Shock tube simulation video Objective: To understand the working of Shock Tube To understand the concept of \'Events\' in the Converge Studio Shock Tube setup and working:        The shock tube is a simple Read more

Prandtl Meyer Shock problem Krishnan Jothi · 2018-09-02 15:14:43

Objective: To understand the theory of formation of expansion fan due to infinite mach waves, originated at the sharp convex corner   To understand the reason for the usage of Neumann BC at the outlet for supersonic flow To understand the effect of SGS pa Read more

Conjugate Heat Transfer Simulation Krishnan Jothi · 2018-09-01 16:32:46

OBJECTIVE OF THIS PROJECT IS THE UNDERSTAND THE FOLLOWING : CHT simulation is used to simultaneously predict the heat transfer in the fluid and the solid portions of the domain. In fluid heat transfer takes place faster when compared to solids, therefore to overcome thi Read more

Transient simulation of flow over a throttle body Krishnan Jothi · 2018-08-19 12:15:22

1. Transient flow simulation of fluid over a throttle body is carried out using transient solver. No hydrodynamic mode can be used to check whether the mesh generation process is taking place correctly  2. Reynolds Averaged Navier-Stokes method is used to solve th Read more

Steady state simulation of flow over a throttle body Krishnan Jothi · 2018-08-19 10:45:59

1. Steady state flow simulation of fluid over a throttle body is carried out using Pressure-based steady solver  2. Reynolds Averaged Navier-Stokes method is used to solve the compressible fluid flow inside the bent pipe and RNG-K epsilon turbulence model is Read more

Flow over a backward facing step Krishnan Jothi · 2018-08-13 17:05:17

1. Backward step geometry is drawn in the converge studio. Dimensions are determined from the file provided by the skill lync. 2. Three different base mesh sizes are used (Converge studio resulted that the order of 10^-4 is too small and could result in data overflow, Read more

Channel flow simulation using CONVERGE CFD Krishnan Jothi · 2018-08-09 18:12:58

1. Channel flow case is set, with the help of tutorial videos provided by convergent science. 2. Three different base mesh sizes are used 2e-4m 1.5e-4m 1.2e-4m 3. Cygwin is used to run the simulation and in the converge studio the binary outputs are post proces Read more

BlockMesh Drill down challenge Krishnan Jothi · 2018-07-31 18:58:24

The main objective of this project is to simulate the flow through a backward facing step using icoFoam solver, with multiple meshes. But as an extensive study, flow through the backward facing step is simulated using two solvers, in this project. M Read more