Taylor series

The Taylor series of a function is an infinite sum of terms that are expressed in terms of the function's derivatives at a single point. For most common functions, the function and the sum of its Taylor series are equal near this point.

`f(x)=f(a)+(f'(a))/(1!)(x−a)+(f''(a))/(2!)(x−a)^2+(f'''(a))/(3!)(x−a)^3+...`

As the degree of the Taylor polynomial rises, it approaches the correct function. This image shows sin x and its Taylor approximations by polynomials of degree 1, 3, 5, 7, 9, 11, and 13 at x = 0.

For example, the below graph shows approximation at x=0, here we can see that as the number of terms increases the two curves becomes coincident.

The exponential function ex (in blue), and the sum of the first n + 1 terms of its Taylor series at 0 (in red).

The exponential function ex (in blue), and the sum of the first n + 1 terms of its Taylor series at 0 (in red).

Deriving the 4th order approximations of the second-order derivative

Let's consider a numerical stencile with a set of grid points as shown in the fig below

,

right-sided difference left -sided difference

Above shown both numerical stencils are Skewed 2nd order right & left -sided difference scheme for fourth order opproximation

central-difference

And this is symmetric central difference for 4th order opproximation

1) Skewed 2nd order right-sided difference scheme

The number of Nodes required for 4th order approximation of a second-order derivative can be given by the equation

N=P+Q where

N is the number of nodes

Q is the order of derivative

P is the order of approximation

So, here we need 6 nodes in order to find the 2nd order derivative with the 4th order of approximation.

Consider the value of `(∂^2y)/(∂x^2)` with the 4th order of approximation is given by the equation

`(∂^2y)/(∂x^2)=(af(i)+bf(i+1)+cf(i+2)+df(i+3)+ef(i+4)+gf(i+5))/(Deltax^2)`

Now, using taylor series to find the values of discrete poits with repect to point i

`af(i)=af(i)`

`bf(i+1)=bf(i)+bΔxf'(i)+b(Δx^2)/(2!)f''(i)+b(Δx^3)/(3!)f'''(i)+b(Δx^4)/(4!)f''''(i)+b(Δx^5)/(5!)f'''''(i)`

`cf(i+2)=cf(i)+c(2Δx)f'(i)+c((2Δx)^2)/(2!)f''(i)+c((2Δx)^3)/(3!)f'''(i)+c((2Δx)^4)/(4!)f''''(i)+c((2Δx)^5)/(5!)f'''''(i)`

Here we have `2Δx` because from i to i+2 the distance is 2 * `Δx` so `2Δx` as well as i+3,i+4,....

`df(i+3)=df(i)+d(3Δx)f'(i)+d((3Δx)^2)/(2!)f''(i)+d((3Δx)^3)/(3!)f'''(i)+d((3Δx)^4)/(4!)f''''(i)+d((3Δx)^5)/(5!)f'''''(i)`

`ef(i+4)=ef(i)+e(4Δx)f'(i)+e((4Δx)^2)/(2!)f''(i)+e((4Δx)^3)/(3!)f'''(i)+e((4Δx)^4)/(4!)f''''(i)+e((4Δx)^5)/(5!)f'''''(i)`

`gf(i+5)=gf(i)+g(5Δx)f'(i)+g.((5Δx)^2)/(2!)f''(i)+g.((5Δx)^3)/(3!)f'''(i)+g.((5Δx)^4)/(4!)f''''(i)+g.((5Δx)^5)/(5!)f'''''(i)`

In the above equations, we want only the coefficients of `Δx^2f''(i)`and we assume the coefficient to be one and for the remaining terms we assume the coefficients to be zero.

`a+b+c+d+e+g=0`

`b+2c+3d+4e+5g=0`

`b/2+(4c)/2+(9d)/2+(16e)/2+(25g)/2=1`

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

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

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

Taylor Table
Assumed Function		Taylor Series Expansion
		f(i)"	Δf′(i)	`Δx^2f''(i)`	`Δx^3f'''(i)`	`Δx^4f''''(i)`	`Δx^5f'''''(i)`
Stencil	af(i)	a	a	0	0	0	0
	bf(i+1)	b	b	`b/2`	`b/6`	`b/24`	`b/120`
	cf(i+2)	c	2c	`(4c)/2`	`(8c)/6`	`(16c)/24`	`(32c)/120`
	df(i+3)	d	3d	`(9d)/2`	`(27d)/6`	`(81d)/24`	`(243d)/120`
	ef(i+4)	e	4e	`(16e)/2`	`(64e)/6`	`(256e)/24`	`(1024e)/120`
	gf(i+5)	g	5g	`(25g)/2`	`(125g)/6`	`(625g)/24`	`(3125g)/120`
Sum of Coefficients		0	0	1	0	0	0

Writing the coefficients in the form of a matrix

Solving using MATLAB to find the values of the coefficients

 A=[1 1 1 1 1 1;0 1 2 3 4 5;0 1/2 4/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_r=inv(A)*B

A =

    1.0000    1.0000    1.0000    1.0000    1.0000    1.0000
         0    1.0000    2.0000    3.0000    4.0000    5.0000
         0    0.5000    2.0000    4.5000    8.0000   12.5000
         0    0.1667    1.3333    4.5000   10.6667   20.8333
         0    0.0417    0.6667    3.3750   10.6667   26.0417
         0    0.0083    0.2667    2.0250    8.5333   26.0417


B =

     0
     0
     1
     0
     0
     0


x_r =

    3.7500
  -12.8333
   17.8333
  -13.0000
    5.0833
   -0.8333

2) Skewed 2nd order Left-sided difference scheme

The number of Nodes required for 4th order approximation of a second-order derivative can be given by the equation

N=P+Q where

N is the number of nodes

Q is the order of derivative

P is the order of approximation

So, here we need 6 nodes in order to find the 2nd order derivative with the 4th order of approximation.

Consider the value of `(∂^2y)/(∂x^2)` with the 4th order of approximation is given by the equation

`(∂^2y)/(∂x^2)=(af(i)+bf(i-1)+cf(i-2)+df(i-3)+ef(i-4)+gf(i-5))/(Deltax^2)`

Now, using taylor series to find the values of discrete poits with repect to point i

`af(i)=af(i)`

`bf(i-1)=bf(i)-bΔxf'(i)+b(Δx^2)/(2!)f''(i)-b(Δx^3)/(3!)f'''(i)+b(Δx^4)/(4!)f''''(i)-b(Δx^5)/(5!)f'''''(i)`

`cf(i-2)=cf(i)-c(2Δx)f'(i)+c((2Δx)^2)/(2!)f''(i)-c((2Δx)^3)/(3!)f'''(i)+c((2Δx)^4)/(4!)f''''(i)-c((2Δx)^5)/(5!)f'''''(i)`

Here we have `2Δx` because from i to i-2 the distance is 2 * `Δx` so `2Δx` as well as i-3,i-4,....

`df(i-3)=df(i)-d(3Δx)f'(i)+d((3Δx)^2)/(2!)f''(i)-d((3Δx)^3)/(3!)f'''(i)+d((3Δx)^4)/(4!)f''''(i)-d((3Δx)^5)/(5!)f'''''(i)`

`ef(i-4)=ef(i)-e(4Δx)f'(i)+e((4Δx)^2)/(2!)f''(i)-e((4Δx)^3)/(3!)f'''(i)+e((4Δx)^4)/(4!)f''''(i)-e((4Δx)^5)/(5!)f'''''(i)`

`gf(i-5)=gf(i)-g(5Δx)f'(i)+g.((5Δx)^2)/(2!)f''(i)-g.((5Δx)^3)/(3!)f'''(i)+g.((5Δx)^4)/(4!)f''''(i)-g.((5Δx)^5)/(5!)f'''''(i)`

In the above equations, we want only the coefficients of `Δx^2f''(i)`and we assume the coefficient to be one and for the remaining terms we assume the coefficients to be zero.

`a+b+c+d+e+g=0`

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

`b/2+(4c)/2+(9d)/2+(16e)/2+(25g)/2=1`

`-(b/6+(8c)/6+(27d)/6+(64e)/6+(125g)/6=0)`

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

`-(b/120+(32c)/120+(243d)/120+(1024e)/120+(3125g)/120=0)`

Taylor Table
Assumed Function		Taylor Series Expansion
		f(i)	Δf′(i)	`Δx^2f''(i)`	`Δx^3f'''(i)`	`Δx^4f''''(i)`	`Δx^5f'''''(i)`
Stencil	af(i)	a	-a	0	0	0	0
	bf(i+1)	b	-b	`b/2`	`-b/6`	`b/24`	`-b/120`
	cf(i+2)	c	-2c	`(4c)/2`	`-(8c)/6`	`(16c)/24`	`-(32c)/120`
	df(i+3)	d	-3d	`(9d)/2`	`-(27d)/6`	`(81d)/24`	`-(243d)/120`
	ef(i+4)	e	-4e	`(16e)/2`	`-(64e)/6`	`(256e)/24`	`-(1024e)/120`
	gf(i+5)	g	-5g	`(25g)/2`	`-(125g)/6`	`(625g)/24`	`-(3125g)/120`
Sum of Coefficients		0	0	1	0	0	0

Writing the coefficients in the form of a matrix

Solving using MATLAB to find the values of the coefficients

 C=[1 1 1 1 1 1;0 -1 -2 -3 -4 -5;0 1/2 4/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]
  D=[0;0;1;0;0;0]
  x_l=inv(C)*D

C =

    1.0000    1.0000    1.0000    1.0000    1.0000    1.0000
         0   -1.0000   -2.0000   -3.0000   -4.0000   -5.0000
         0    0.5000    2.0000    4.5000    8.0000   12.5000
         0   -0.1667   -1.3333   -4.5000  -10.6667  -20.8333
         0    0.0417    0.6667    3.3750   10.6667   26.0417
         0   -0.0083   -0.2667   -2.0250   -8.5333  -26.0417


D =

     0
     0
     1
     0
     0
     0


x_l =

    3.7500
  -12.8333
   17.8333
  -13.0000
    5.0833
   -0.8333

3)symmetric central difference scheme

The number of Nodes required for 4th order approximation of a second-order derivative can be given by the equation

N=P+Q where

N is the number of nodes

Q is the order of derivative

P is the order of approximation

So, here we need 5 nodes in order to find the 2nd order derivative with the 4th order of approximation.

Consider the value of `(∂^2y)/(∂x^2)` with the 4th order of approximation is given by the equation

`(∂^2y)/(∂x^2)=(af(i-2)+bf(i-1)+cf(i)+df(i+1)+ef(i+2))/(Deltax^2)`

Now, using taylor series to find the values of discrete poits with repect to point i

`af(i-2)=af(i)-a2Δxf'(i)+a(2Δx^2)/(2!)f''(i)-a(2Δx^3)/(3!)f'''(i)+a(2Δx^4)/(4!)f''''(i)-a(2Δx^5)/(5!)f'''''(i)`

`bf(i-1)=bf(i)-b(Δx)f'(i)+b((Δx)^2)/(2!)f''(i)-b((Δx)^3)/(3!)f'''(i)+b((Δx)^4)/(4!)f''''(i)-b((Δx)^5)/(5!)f'''''(i)`

`cf(i)=cf(i)`

`df(i+1)=df(i)+d(Δx)f'(i)+d((Δx)^2)/(2!)f''(i)+d((Δx)^3)/(3!)f'''(i)+d((Δx)^4)/(4!)f''''(i)+d((Δx)^5)/(5!)f'''''(i)`

`ef(i+2)=ef(i)+e(2Δx)f'(i)+e((2Δx)^2)/(2!)f''(i)+e((2Δx)^3)/(3!)f'''(i)+e((2Δx)^4)/(4!)f''''(i)+e((2Δx)^5)/(5!)f'''''(i)`

Here we have `2Δx` because from i to i+2 the distance is 2 * `Δx` so `2Δx` as well as i+3,i+4,....

In the above equations, we want only the coefficients of `Δx^2f''(i)`and we assume the coefficient to be one and for the remaining terms we assume the coefficients to be zero.

Consider Taylor Table for central difference

Taylor Table
Assumed Function		Taylor Series Expansion
		f(i)	Δxf′(i)	`Δx^2f''(i)`	`Δx^3f'''(i)`	`Δx^4f''''(i)`
Stencil	af(i-2)	a	-2a	`( 4a)/2`	`−(8a)/6`	`(16a)/24`
	bf(i-1)	b	-b	`b/2`	`−b/6`	`b/24`
	cf(i)	c	0	0	0	0
	df(i+1)	d	d	`d/2`	`d/6`	`d/24`
	ef(i+2)	e	2e	`(4e)/2`	`−(8e)/6`	`(16e)/24`
Sum of Coefficients		0	0	1	0	0

Writing the coefficients in the form of a matrix

Solving using MATLAB to find the values of the coefficients

>> E=[1, 1, 1, 1, 1;-2, -1, 0, 1, 2;4/2, 1/2, 0, 1/2, 4/2;-8/6, -1/6, 0, 1/6, 8/6;16/24, 1/24, 0, 1/24, 16/24]
  F=[0;0;1;0;0]
  x_c=E^-1*F

E =

    1.0000    1.0000    1.0000    1.0000    1.0000
   -2.0000   -1.0000         0    1.0000    2.0000
    2.0000    0.5000         0    0.5000    2.0000
   -1.3333   -0.1667         0    0.1667    1.3333
    0.6667    0.0417         0    0.0417    0.6667


F =

     0
     0
     1
     0
     0


x_c =

   -0.0833
    1.3333
   -2.5000
    1.3333
   -0.0833

>> 

CODE EXPLINATION

  clear all 
  close all 
  clc
  
  %inputs
  
  x=pi/3;
  
  %stepsize range
  dx=linspace(pi/4,pi/400,2000);
  
  %Analytical derivative
  
  %f(x) = exp(x)*cos(x)
  %f'(x)=(exp(x)*cos(x)-(exp(x)*sin(x)=exp(x)*(cosx-sinx)
  %Double derivative f" = -2*exp(x)*sin(x)
  
  analytical_derivative = -2*exp(x)*sin(x);
  
  %Skewed right-sided(4th order approimation)
  %Solving for the coefficients
  A=[1 1 1 1 1 1;0 1 2 3 4 5;0 1/2 4/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_r=inv(A)*B;
  
  %Skewed left-sided(4th order approimation)
  %Solving for the coefficients
  C=[1 1 1 1 1 1;0 -1 -2 -3 -4 -5;0 1/2 4/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];
  D=[0;0;1;0;0;0];
  x_l=inv(C)*D;
  
  %central differencing(4th order approimation)
  %Solving for the coefficients
  E=[1, 1, 1, 1, 1;-2, -1, 0, 1, 2;4/2, 1/2, 0, 1/2, 4/2;-8/6, -1/6, 0, 1/6, 8/6;16/24, 1/24, 0, 1/24, 16/24];
  F=[0;0;1;0;0];
  x_c=E^-1*F;
  
  %Anonymous Function
  
  f=@(x) exp(x).*cos(x);
  
  
  %numerical derivative
  
  right_difference=(x_r(1)*f(x)+x_r(2)*f(x+dx)+x_r(3)*f(x+2.*dx)+x_r(4)*f(x+3.*dx)+x_r(5)*f(x+4.*dx)+x_r(6)*f(x+5.*dx))./(dx.^2);
  
  left_difference=(x_r(1)*f(x)+x_r(2)*f(x-dx)+x_r(3)*f(x-2.*dx)+x_r(4)*f(x-3.*dx)+x_r(5)*f(x-4.*dx)+x_r(6)*f(x-5.*dx))./(dx.^2);
  
  cent_difference=(x_c(1)*f(x-2.*dx)+x_c(2)*f(x-1.*dx)+x_c(3)*f(x)+x_c(4)*f(x+dx)+x_c(5)*f(x+2.*dx))./(dx.^2);
  
  
  
  %absolute error
  
  error_r=abs(analytical_derivative-right_difference);
  error_l=abs(analytical_derivative-left_difference);
  error_c=abs(analytical_derivative-cent_difference);
  error=[error_r;error_l;error_c];
  
  
  
  figure(1)
  loglog(dx,error_r) 
  grid on
  hold on
  loglog(dx,error_l)
  loglog(dx,error_c)
  ylabel('Absolute Error (Log Scale)')
  xlabel('Grid Spacing (dx in Log Scale)')
  legend('Right Differencing','Left Differencing','Central Differencing')
     

inputs

x=pi/3;

stepsize range dx=linspace(pi/4,pi/400,2000);

Analytical derivative

`f(x) = e^x*cos(x)`

`f'(x)=e^x*(cosx-sinx)`

Double derivative `f'' = -2*e^x*sin(x)`
analytical_derivative = `-2*e^x*sin(x)`

Here we are creating the anonymous funtion with an handle with respect to x as shown

f=@(x) exp(x).*cos(x)

And then give the formulas and calculate the error

  %numerical derivative
  
  right_difference=(x_r(1)*f(x)+x_r(2)*f(x+dx)+x_r(3)*f(x+2.*dx)+x_r(4)*f(x+3.*dx)+x_r(5)*f(x+4.*dx)+x_r(6)*f(x+5.*dx))./(dx.^2);
  
  left_difference=(x_r(1)*f(x)+x_r(2)*f(x-dx)+x_r(3)*f(x-2.*dx)+x_r(4)*f(x-3.*dx)+x_r(5)*f(x-4.*dx)+x_r(6)*f(x-5.*dx))./(dx.^2);
  
  cent_difference=(x_c(1)*f(x-2.*dx)+x_c(2)*f(x-1.*dx)+x_c(3)*f(x)+x_c(4)*f(x+dx)+x_c(5)*f(x+2.*dx))./(dx.^2);

PLOTS

As we can see from the above figure the error for the central difference method is very less compared to the skewed methods.

Skewed right and left difference methods almost gives the same results in terms of error for a small value of time step.

• When data from the grid points is only available from either the left or right end i.e. at nodes closer to the boundaries,for computations that requires a high order approximations we make use of skewed schemes.Having a symmetric stencil for computation is not possioble in this scenarios like we do for central differencing schemes or having a low order accuracy.For better consistency of scheme to be maintained the values of dx are kept as low as possible for providing proper computational plots.
• A skewed scheme is useful when we have to predict using the boundary condition as we have only one direction to go. Hence the skewed scheme will be useful and we cannot use the central differencing scheme for such an application because it required information from both directions of the point.