Derivation and approximation of fourth order approximation schemes ( Central difference, Right-skewed, Left-skewed) for second order derivative and their error comparison for function e^(x)*cos(x)

1) FOURTH ORDER APPROXIMATION OF THE SECOND ORDER DERIVATIVE CENTRAL DIFFERENCE

`trianglex^2*(del^2f)/(delx) approx a*f(i-2) + b*f(i-1) + c*f(i) + d(f(i+1) +e*f(i+2)` ...(equation 'I')`

Expanding the above equation using Taylor table method,

		Taylor series terms
		`f(i)`	`trianglex*f'(x)`	`trianglex^2 * f''(x)`	`trianglex^3*f'''(x)`	`trianglex^4*f^(IV)(x)`	`trianglex^5*f^V(x)`	`trianglex^6*f^(VI)(x)`
Numerical Stencil	`a*f(i-2)`	a	-2a	`(4a)/(2!)`	`-(8a)/(3!)`	`(16a)/(4!)`	`-(32a)/(5!)`	`(64a)/(6!)`
	`b*f(i-1)`	b	-b	`b/(2!)`	`-(b)/(3!)`	`(b)/(4!)`	`-(b)/(5!)`	`(b)/(6!)`
	`c*f(i)`	c	0	0	0	0	0	0
	`d*f(i+1)`	d	d	`d/(2!)`	`d/(3!)`	`d/(4!)`	`d/(5!)`	`(d)/(6!)`
	`e*f(i+2)`	e	2e	`(4e)/(2!)`	`(8e)/(3!)`	`(16e)/(4!)`	`(32e)/(5!)`	`(64e)/(6!)`
C0-efficent in equation		0	0	1	0	0	?	?

We can get 5 equations from the above Taylor table,

`a+b+c+d+e = 0`... equation (1)

`-2a+b+0+d+2e = 0`... equation (2)

`(4a)/2 + b/2 +0 + d/2 + (4e)/2 = 1`... equation (3)

`-(8a)/6 - b/6 +0 +d/6 +(8e)/6 = 0`... equation (4)

`(16a)/(24) + b/24 + 0 + d/24 + (16e)/24 = 0`... equation (5)

By solving the above 5 equations we get,

a= `-1/12`

b= `4/3`

c= `-5/2`

d= `4/3`

e= `-1/12`

If we substitute the values of a,b,d and e in the Taylor table than we get the co-efficient of equation in 6th column = 0

`trianglex^2(del^2f)/(delx) approx (2a +1/2b + 0 + 1/2d + 2e ) *f''(x)*trianglex^2 + ((64a)/720 + b/720 + 0 + d/720 + (64e)/720)*f^(VI) (x) *trianglex^6` ...(equation 'II')

In the above equation `trianglex^6 = 0` because the approximation is only till fourth order,

`equate 'I' and 'II' and dividing by `trianglex^2`,

`(del^2f)/(delx^2) = [(a*f(i-2) + b*f(i-1) +c*f(i) + d*f(i+1) + e*f(i+2))/(trianglex^2)]`

substituting the values of a,b,c,d, and e in the equation 'I' we get

`(del^2f)/(delx^2) = [(-f(i+2) + 16f(i+1) - 30f(i) + 16f(i-1) - f(i-2))/(12trianglex^2)]`

From equation 'II' we get,

`trianglex^2(del^2f)/(delx^2) approx trianglex^2f''(x) + trianglef^(VI)(x) (constants) = a*f(i-2) + b*f(i-1) + c*f(i) + d(f(i+1) +e*f(i+2)`

dividing throughout by `trianglex^2`,

`(del^2f)/(delx^2) approx f''(x) + trianglex^4f^(IV)(x) (constants) = (a*f(i-2) + b*f(i-1) + c*f(i) + df(i+1) +e*f(i+2))/(trianglex^2)`

Here the power of `trianglex` is 4 and hence this proves the approximation is of fourth order accurate.

2) FOURTH ORDER APPROXIMATION OF THE SECOND ORDER DERIVATIVE RIGHT-SKEWED

`trianglex^2*(del^2f)/(delx) approx a*f(i)+b*f(i+1)+c*f(i+2)+df(i+3)+ef(i+4)+gf(i+5)` ...(equation 'III')`

Expanding the above equation using Taylor table method,

		Taylor series terms
		`f(i)`	`trianglex*f'(x)`	`trianglex^2 * f''(x)`	`trianglex^3*f'''(x)`	`trianglex^4*f^(IV)(x)`	`trianglex^5*f^V(x)`	`trianglex^6*f^(VI)(x)`
Numerical Stencil	`a*f(i)`	a	`0`	`0`	`0`	`0`	`0`	`0`
	`b*f(i+1)`	b	b	`b/(2!)`	`(b)/(3!)`	`(b)/(4!)`	`(b)/(5!)`	`(b)/(6!)`
	`c*f(i+2)`	c	2c	`(2^2c)/(2!)`	`(2^3c)/(3!)`	`(2^4c)/(4!)`	`(2^5c)/(5!)`	`(2^6c)/(6!)`
	`d*f(i+3)`	d	3d	`(3^2d)/(2!)`	`(3^3d)/(3!)`	`(3^4d)/(4!)`	`(3^5d)/(5!)`	`(3^6d)/(6!)`
	`e*f(i+4)`	e	4e	`(4^2e)/(2!)`	`(4^3e)/(3!)`	`(4^4e)/(4!)`	`(4^5e)/(5!)`	`(4^6e)/(6!)`
	`g*f(i+5)`	g	5g	`(5^2g)/(2!)`	`(5^3g)/(3!)`	`(5^4g)/(4!)`	`(5^5g)/(5!)`	`(5^6g)/(6!)`
C0-efficent in equation		0	0	1	0	0	0	?

We can get 5 equations from the above Taylor table,

`a+b+c+d+e +g = 0`... equation (1)

`0 + b + 2c +3d +4e +5g = 0`... equation (2)

`0 + b/2 + (4*c)/2 + (9*d)/2 + (16*e)/2 + (25*g)/2 = 1`... equation (3)

`0 + b/6 + (8*c)/6 + (27*d)/6 + (64*e)/6 + (125*g)/6 = 0`... equation (4)

`0 + b/24 + (16*c)/24 + (81*d)/24 + (256*e)/24 + (625*g)/24 = 0`... equation (5)

`0 + b/120 + (32*c)/120 + (243*d)/120 + (1024*e)/120 + (3125*g)/120 = 0`... equation (6)

By solving the above 5 equations we get,

a= `15/4`

b= `-77/6`

c= `107/6`

d= `-13`

e= `61/12`

g=`-5/6`

If we substitute the values of a,b,d,e and g in the Taylor table than we get the co-efficient of equation in 6th column = 0

`trianglex^2(del^2f)/(delx) approx (0 + b/2 + 2c + (9c)/2+ 8e + (25g)/2 ) *f''(x)*trianglex^2 + (0a + b/720 + (64c)/720 + (729d)/720 + (4096*e)/720 + (15625*g)/720)*f^(VI) (x) *trianglex^6` ...(equation 'IV')`

In the above equation `trianglex^6 = 0` because the approximation is only till fourth order,

equate 'III' and 'IV' and dividing by `trianglex^2`,

`(del^2f)/(delx^2) = [(a*f(i) + b*f(i+1) +c*f(i+2) + d*f(i+3) + e*f(i+4)+g*f(i+5))/(trianglex^2)]`

Substituting the values of a,b,c,d,e and g in the equation 'III' we get,

`(del^2f)/(delx^2) = [(45f(i)-154f(i+1)+214f(i+2) - 156f(i+3) + 61f(i+4) -10f(i+5))/(12trianglex^2)]`

From equation 'IV' we get

`trianglex^2(del^2f)/(delx^2) approx trianglex^2f''(x) + trianglef^(VI)(x) (constants) = a*f(i) + b*f(i+1) +c*f(i+2) + d*f(i+3) + e*f(i+4)+g*f(i+5)`

dividing throughout by `trianglex^2`,

`(del^2f)/(delx^2) approx f''(x) + trianglex^4f^(IV)(x) (constants) = (a*f(i) + b*f(i+1) +c*f(i+2) + d*f(i+3) + e*f(i+4)+g*f(i+5))/(trianglex^2)`

Here the power of `trianglex` is 4 and hence this proves the approximation is of fourth order accurate.

3) FOURTH ORDER APPROXIMATION OF THE SECOND ORDER DERIVATIVE LEFT-SKEWED

`trianglex^2*(del^2f)/(delx) approx a*f(i)+b*f(i-1)+c*f(i-2)+df(i-3)+ef(i-4)+gf(i-5)` ...(equation 'V')`

Expanding the above equation using Taylor table method,

		Taylor series terms
		`f(i)`	`trianglex*f'(x)`	`trianglex^2 * f''(x)`	`trianglex^3*f'''(x)`	`trianglex^4*f^(IV)(x)`	`trianglex^5*f^V(x)`	`trianglex^6*f^(VI)(x)`
Numerical Stencil	`a*f(i)`	a	`0`	`0`	`0`	`0`	`0`	`0`
	`b*f(i-1)`	b	-b	`b/(2!)`	`-(b)/(3!)`	`(b)/(4!)`	`-(b)/(5!)`	`(b)/(6!)`
	`c*f(i-2)`	c	-2c	`(2^2c)/(2!)`	`-(2^3c)/(3!)`	`(2^4c)/(4!)`	`-(2^5c)/(5!)`	`(2^6c)/(6!)`
	`d*f(i-3)`	d	-3d	`(3^2d)/(2!)`	`-(3^3d)/(3!)`	`(3^4d)/(4!)`	`-(3^5d)/(5!)`	`(3^6d)/(6!)`
	`e*f(i-4)`	e	-4e	`(4^2e)/(2!)`	`-(4^3e)/(3!)`	`(4^4e)/(4!)`	`-(4^5e)/(5!)`	`(4^6e)/(6!)`
	`g*f(i-5)`	g	-5g	`(5^2g)/(2!)`	`-(5^3g)/(3!)`	`(5^4g)/(4!)`	`-(5^5g)/(5!)`	`(5^6g)/(6!)`
C0-efficent in equation		0	0	1	0	0	0	?

We can get 5 equations from the above Taylor table,

`a+b+c+d+e +g = 0`... equation (1)

`0 - b - 2c -3d -4e -5g = 0`... equation (2)

`0 + b/2 + (4*c)/2 + (9*d)/2 + (16*e)/2 + (25*g)/2 = 1`... equation (3)

`0 - b/6 - (8*c)/6 - (27*d)/6 - (64*e)/6 - (125*g)/6 = 0`... equation (4)

`0 + b/24 + (16*c)/24 + (81*d)/24 + (256*e)/24 + (625*g)/24 = 0`... equation (5)

`0 - b/120 - (32*c)/120 - (243*d)/120 - (1024*e)/120 - (3125*g)/120 = 0`... equation (6)

By solving the above 5 equations we get,

a= `15/4`

b= `-77/6`

c= `107/6`

d= `-13`

e= `61/12`

g=`-5/6`

If we substitute the values of a,b,d,e and g in the Taylor table than we get the co-efficient of equation in 6th column = 0

`trianglex^2(del^2f)/(delx) approx (0 - b/2 - 2c - (9c)/2- 8e - (25g)/2 ) *f''(x)*trianglex^2 + (0a + b/720 + (64c)/720 + (729d)/720 + (4096*e)/720 + (15625*g)/720)*f^(VI) (x) *trianglex^6` ...(equation 'VI')`

In the above equation `trianglex^6 = 0` because the approximation is only till fourth order,

equate 'V' and 'VI' and dividing by `trianglex^2`,

`(del^2f)/(delx^2) = [(a*f(i) + b*f(i-1) +c*f(i-2) + d*f(i-3) + e*f(i-4)+g*f(i-5))/(trianglex^2)]`

Substituting the values of a,b,c,d,e and g in the equation 'V' we get,

`(del^2f)/(delx^2) = [(45f(i)-154f(i-1)+214f(i-2) - 156f(i-3) + 61f(i-4) -10f(i-5))/(12trianglex^2)]`

From equation 'VI' we get

`trianglex^2(del^2f)/(delx^2) approx trianglex^2f''(x) + trianglef^(VI)(x) (constants) = a*f(i) + b*f(i-1) +c*f(i-2) + d*f(i-3) + e*f(i-4)+g*f(i-5)`

dividing throughout by `trianglex^2`,

`(del^2f)/(delx^2) approx f''(x) + trianglex^4f^(IV)(x) (constants) = (a*f(i) + b*f(i-1) +c*f(i-2) + d*f(i-3) + e*f(i-4)+g*f(i-5))/(trianglex^2)`

Here the power of `trianglex` is 4 and hence this proves the approximation is of fourth order accurate.

a) The plots of comparison between absolute and fourth-order approximation for the above three schemes are shown below

Comparison plot on a log scale

Comparison plot on a normal scale

b) The skewed scheme is useful because

Skewed schemes are more useful than central differencing method when using initial or end boundary conditions because in the initial conditions the data or previous or before points are not available and while using end boundary conditions the data of the forward or after data points are not available and in such situations, the right-skewed or left-skewed schemes can be used as per the data of points available.

4)Code for Matlab function

function [error] = Week_3_taylor_table_fun(x,dx)
%program to evaluate exp(x)*cos(x)

%First Finding analytical derivative for the given equation
%f(x) = (e^x)*cos(x)
%f'(x) = (e^x)*cos(x) - (e^x)*sin(x)
%f''(x) = -2*(e^x)*sin(x)
analytical_derivative = (-2)*(exp(x))*(sin(x));

%solving using numerical discretization of fourth order
%central differencing = f''(x) = (-f(i+2) + 16f(i+1) - 30f(i) + 16f(i-1) - f(i-2))/(12*dx^2)
central_difference = (-((exp(x+2*dx))*cos(x+2*dx)) + 16*((exp(x+dx))*cos(x+dx)) - 30* ((exp(x))*cos(x)) + 16*((exp(x-dx))*cos(x-dx)) - ((exp(x-2*dx))*cos(x-2*dx)))/(12*dx^2);

right_skewed = (45*((exp(x))*cos(x)) - 154*((exp(x+dx))*cos(x+dx)) + 214* ((exp(x+2*dx))*cos(x+2*dx)) - 156*((exp(x+3*dx))*cos(x+3*dx)) + 61* ((exp(x+4*dx))*cos(x+4*dx)) -10*((exp(x+5*dx))*cos(x+5*dx)))/(12*dx^2);

left_skewed = (45*((exp(x))*cos(x)) - 154*((exp(x-dx))*cos(x-dx)) + 214* ((exp(x-2*dx))*cos(x-2*dx)) - 156*((exp(x-3*dx))*cos(x-3*dx)) + 61* ((exp(x-4*dx))*cos(x-4*dx)) -10*((exp(x-5*dx))*cos(x-5*dx)))/(12*dx^2);

error(1) = abs(central_difference - analytical_derivative);

error(2) = abs(right_skewed - analytical_derivative);

error(3) = abs(left_skewed - analytical_derivative);

end

Code sample for Matlab plotting

close all
clear all
clc
%code to compare absolute error between three schemes and analytical error
x = pi/4;
dx = linspace(pi/4,pi/300,70);

%preallocating for speed
central_difference_error  = zeros(1,length(dx));
right_skewed_error  = zeros(1,length(dx));
left_skewed_error  = zeros(1,length(dx));

for i = 1:length(dx)
    error = Week_3_taylor_table_fun(x,dx(i));
    central_difference_error(i) = error(1);
    right_skewed_error(i) = error(2);
    left_skewed_error(i) = error(3);
end

figure(1)
loglog(dx,central_difference_error,'-db')
hold on
loglog(dx,right_skewed_error,'-sblack')
hold on
loglog(dx,left_skewed_error,'-or')
xlabel('dx')
ylabel('Absolute Error')
legend('Central Difference of 4th Order','Skewed Right Sided Difference of 4th Order','Skewed Left Sided Difference of 4th Order')
title('Comparison of 4th order approximations on 2nd order derivative on log scale')

figure(2)
plot(dx,central_difference_error,'-db')
hold on
plot(dx,right_skewed_error,'-sblack')
hold on
plot(dx,left_skewed_error,'-or')
xlabel('dx')
ylabel('Absolute Error')
legend('Central Difference of 4th Order','Skewed Right Sided Difference of 4th Order','Skewed Left Sided Difference of 4th Order')
title('Comparison of 4th order approximations on 2nd order derivative on normal scale')

5) A detailed explanation of the code and results/plot

1]Explanations of Matlab function and code

Matlab function was created with variable inputs for x and dx. 'x' was given a constant value of 'pi/4' and dx was given a range from pi/4 to pi/300 with a number of divisions of 70. values of x and dx can be changed as per the requirement. By lowering the number of divisions more accurate curve can be obtained. an analytical second-order derivative of the given function was calculated with the values of x and dx. Also, the fourth-order approximation for the second-order derivative of the given function was calculated and the difference between the numerical and analytical value was simulated in the function.

In the main Matlab code, a for loop was used to create a loop from 1st value of dx to the last value of dx and obtain results for every value of dx. Two graphs were plotted showing the comparison of error in all the three schemes ( central, right-skewed, and left-skewed). It can be seen from the graph that center difference approximation is most accurate than skewed approximation but as explained in one of the sections above it cannot be used in every situation. The right-skewed and left-skewed approximations have almost similar accuracy which can be seen from the plot. It can be seen as the value of dx goes on decreasing the size of error decreases and the accuracy of the scheme improves.