Taylor table method - Comparing the error between CDS, Skewed RSD, Skewed LSD, FDS and BDS approximation schemes over the range of "dx"

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

Once you have done this, write a program in Matlab to evaluate the second order derivative of the analytical function exp(x)*cos(x) and compare it with the 3 numerical approximations that you have derived.

a. Provide a plot that compares the absolute error between the above mentioned schemes

b. How do these errors compare with FDS and BDS approximations

Software used:- MATLAB

Solution:-

Matlab code for finding the coefficient a,b,c,d,e,h for central difference scheme

clear all
close all
clc

% finding the coefficients a,b,c,d,e for central difference scheme using
% taylor table method
% first we have to formulate the equations on paper and then we can solve
% the matrix in matlab
A = [1 1 1 1 1;
    -2 -1 0 1 2;
    2 1/2 0 1/2 2;
    (-2)^3/6 (-1)/6 0 1/6 2^3/6;
    2^4/24 1/24 0 1/24 2^4/24]

B = [0;0;1;0;0]
X = linsolve(A,B)

% Values of coefficients  
% a = -0.0833
% b =  1.3333
% c = -2.5000
% d =  1.3333
% e = -0.0833 

Matlab code for finding the coefficient a,b,c,d,e,h for skewed RSD scheme

clear all
close all
clc

% finding the coefficients a,b,c,d,e for Skewed right sided difference
% using taylor table method
% first we have to formulate the equations on paper and then we can solve
% the matrix in matlab

A = [1 1 1 1 1 1 ;
    0 1 2 3 4 5;
    0 1/2 2 3^2/2 4^2/2 5^2/2;
    0 1/6 2^3/6 3^3/6 4^3/6 5^3/6;
    0 1/24 2^4/24 3^4/24 4^4/24 5^4/24;
    0 1/120 2^5/120 3^5/120 4^5/120 5^5/120]
B = [0; 0; 1; 0; 0; 0]

X = linsolve(A,B)

% Values of Coefficients   
% a =  3.7500
% b = -12.8333
% c =  17.8333
% d = -13.0000
% e =  5.0833
% h = -0.8333

Matlab code for finding the coefficient a,b,c,d,e,h for skewed LSD scheme

clear all
close all
clc

% finding the coefficients a,b,c,d,e for Skewed left sided difference using
% taylor table method
% first we have to formulate the equations on paper and then we can solve
% the matrix in matlab

A = [1 1 1 1 1 1;
    -5 -4 -3 -2 -1 0;
    5^2/2 4^2/2 3^2/2 2^2/2 1/2 0;
    (-5)^3/6 (-4)^3/6 (-3)^3/6 (-2)^3/6 (-1)/6 0;
    5^4/24 4^4/24 3^4/24 2^4/24 1/24 0;
    (-5)^5/120 (-4)^5/120 (-3)^5/120 (-2)^5/120 (-1)/120 0]

B = [0; 0; 1; 0; 0; 0]

X = linsolve(A,B)

% values of coefficient   
% a = -0.8333
% b =  5.0833
% c = -13.0000
% d =  17.8333
% e = -12.8333
% h =  3.7500

Matlab code for comparing the error between different types of schemes over the range of dx

clear all
close all
clc

x = 10;
dx = linspace(pi/4,pi/400,30);
% Given function: f(x) = exp(x)*cos(x);
% using for loop to calculate the error for the range of dx
for i = 1:length(dx)
    
    central_difference_error(i) = CDS_error(x,dx(i))
    skewed_right_sided_difference_error(i) = Skewed_RSD_error(x,dx(i))
    skewed_left_sided_difference_error(i) = Skewed_LSD_error(x,dx(i))
    forward_difference_error(i) = FDS_error(x,dx(i))
    backward_difference_error(i) = BDS_error(x,dx(i))
    
end

% plotting the error
loglog(dx,central_difference_error,'b','linewidth',2)
hold on
loglog(dx,skewed_right_sided_difference_error,'r','linewidth',2)
hold on
loglog(dx,skewed_left_sided_difference_error,'g','linewidth',2)
hold on
loglog(dx,forward_difference_error,'y','linewidth',2)
hold on
loglog(dx,backward_difference_error,'bla','linewidth',2)
legend('CDS','Skewed RSD','Skewed LSD','FDS','BDS','location','best')
xlabel('dx')
ylabel('Error')
title('Comparison of error between different approximation schemes over the range of dx')

Function for calculating the error between analytical solution and CDS

function out = CDS_error(x,dx)

analytical_sol = -2*exp(x)*sin(x)

numerical_sol_CDS = ((-0.0833)*(exp(x-(2*dx))*cos(x-(2*dx))) + 1.3333*(exp(x-dx)*cos(x-dx)) + (-2.5000)*(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(analytical_sol - numerical_sol_CDS)

end

Function for calculating the error between analytical solution and Skewed RSD

function out = Skewed_RSD_error(x,dx)

analytical_sol = -2*exp(x)*sin(x)

numerical_sol_Skewed_RS = (0.1023*(exp(x)*cos(x)) + (-0.3500)*(exp(x+dx)*cos(x+dx)) + 0.4864*(exp(x+(2*dx))*cos(x+(2*dx))) + (-0.3545)*(exp(x+(3*dx))*cos(x+(3*dx))) + 0.1386*(exp(x+(4*dx))*cos(x+(4*dx))) + (-0.0227)*(exp(x+(5*dx))*cos(x+(5*dx))))/((dx)^2);

out = abs(analytical_sol - numerical_sol_Skewed_RS)

end

Function for calculating the error between analytical solution and Skewed LSD

function out = Skewed_LSD_error(x,dx)

analytical_sol = -2*exp(x)*sin(x)

numerical_sol_Skewed_LS = ((-0.8333)*(exp(x-(5*dx))*cos(x-(5*dx))) + 5.0833*(exp(x-(4*dx))*cos(x-(4*dx))) + (-13.0000)*(exp(x-(3*dx))*cos(x-(3*dx))) + 17.8333*(exp(x-(2*dx))*cos(x-(2*dx))) + (-12.8333)*(exp(x-dx)*cos(x-dx)) + 3.7500*(exp(x)*cos(x)))/((dx)^2);

out = abs(analytical_sol - numerical_sol_Skewed_LS);

end

Function for calculating the error between analytical solution and FDS

function out = FDS_error(x,dx)

analytical_sol = -2*exp(x)*sin(x)

numerical_sol_FDS = (exp(x)*cos(x) -2*(exp(x+dx)*cos(x+dx)) + exp(x+(2*dx))*cos(x+(2*dx)))/(dx)^2;

out = abs(analytical_sol - numerical_sol_FDS);

end

Function for calculating the error between analytical solution and BDS

function out = BDS_error(x,dx)

analytical_sol = -2*exp(x)*sin(x)

numerical_sol_BDS = (exp(x-(2*dx))*cos(x-(2*dx)) -2*(exp(x-dx)*cos(x-dx)) + exp(x)*cos(x))/(dx)^2;

out = abs(analytical_sol - numerical_sol_BDS);

end

Report:

• First we have to derive the numerical solution equations of 4th order approximation for CDS, skewed right sided and skewed left sided scheme.
• We can do this by using Taylor table method. We have to derive the equations on paper and then those equations can be solved by finding the coefficients.
• The coefficients "a,b,c,d,e,h" can be calculated on matlab by solving the matrix formed by those equations.
• After finding the coefficients and formulating the final numerical solution equation for all the three approximation scheme, we can now write the main matlab programe.
• First we define the value of "x = 10" and range of dx from "dx=pi/4" to "dx=pi/400".
• Before writing the main program we define the function to calculate the error between analytical solution and numerical solution using different schemes.
• then we write the main prohram to calculate the error over the range of dx.
• Then we plot the result to compare the error between different types of schemes.