Week 3.5 - Deriving 4th order approximation of a 2nd order derivative using Taylor Table method

Deriving Fourth-Order Approximations for Second-Order Derivative:

Let’s derive the following 4th-order approximations of the second-order derivative:

1. Central Difference:

The central difference approximation for the second-order derivative can be written as:

$$f''(x) \approx \frac{a_1 f(x-h) + a_2 f(x) + a_3 f(x+h)}{h^2}$$​

For a fourth-order accurate approximation, we'll use more points. The stencil will be $f(x-2h), f(x-h), f(x), f(x+h), f(x+2h)$, so the approximation becomes:

$$f''(x) \approx \frac{a_1 f(x-2h) + a_2 f(x-h) + a_3 f(x) + a_4 f(x+h) + a_5 f(x+2h)}{h^2}$$​

To determine the coefficients $a_1, a_2, a_3, a_4, a_5$ we expand $f(x)$ in a Taylor series about $x$x and solve for these coefficients.

2. Skewed Right-Sided Difference:

For a skewed right-sided difference, we use a different stencil. For fourth-order accuracy, we use five points:

$$f''(x) \approx \frac{b_1 f(x) + b_2 f(x+h) + b_3 f(x+2h) + b_4 f(x+3h) + b_5 f(x+4h)}{h^2}$$​

Here, we need to solve for $b_1, b_2, b_3, b_4, b_5$

3. Skewed Left-Sided Difference:

Similarly, for a skewed left-sided difference, we have:

$$f''(x) \approx \frac{c_1 f(x) + c_2 f(x-h) + c_3 f(x-2h) + c_4 f(x-3h) + c_5 f(x-4h)}{h^2}$$​

Again, the coefficients $c_1, c_2, c_3, c_4, c_5$ are determined by expanding the function in a Taylor series and solving a system of equations.

---

Step-by-Step Process to Derive the Coefficients Using Matrix Inversion

1. Set up Taylor series expansions: Expand $f(x+h),f(x+2h),\dots \; in\; Taylor\; series\; around$$f(x)$.

2. Match terms: To find the coefficients for fourth-order accuracy, match the derivatives up to the fourth derivative and solve the resulting linear system of equations.

3. Solve the system: Use matrix inversion to solve for the coefficients $a_i, b_i, c_i$.

---

MATLAB Code for Deriving the Coefficients

Here is an example MATLAB code that solves for the coefficients of the fourth-order central difference and skewed schemes using matrix inversion:

clc;

clear;

% Symbolic variables

syms h f0 f1 f2 f3 f4 f_1 f_2 f_3 f_4;

% Taylor series expansion up to h^4

T_f0 = f0; % f(x)

T_f1 = f0 + h*diff(f0) + (h^2/2)*diff(f0,2) + (h^3/6)*diff(f0,3) + (h^4/24)*diff(f0,4); % f(x+h)

T_f2 = f0 + 2*h*diff(f0) + (2*h^2/2)*diff(f0,2) + (8*h^3/6)*diff(f0,3) + (16*h^4/24)*diff(f0,4); % f(x+2h)

% Similarly for f(x-2h), f(x-3h), etc. using symbolic expansion

% System of equations for central difference

% For central difference: f(x), f(x+h), f(x-h), f(x+2h), f(x-2h)

% To solve the system of equations for coefficients, you can set up a

% matrix inversion problem. The idea is to match the Taylor expansion terms

% and solve the linear system.

% Matrix inversion method for central difference coefficients

A = [1, 1, 1, 1, 1;

0, 1, 2, 3, 4;

0, 1^2, 2^2, 3^2, 4^2;

0, 1^3, 2^3, 3^3, 4^3;

0, 1^4, 2^4, 3^4, 4^4];

b = [0; 0; 2; 0; 0]; % Right-hand side for second derivative

coeff = A\b;

% Display the coefficients

disp('Fourth-order central difference coefficients:');

disp(coeff');

MATLAB Code for Comparing Derivatives and Plotting Errors

% MATLAB Code for comparing numerical schemes for second-order derivatives

clc;

clear;

% Define the function and its analytical second derivative

f = @(x) exp(x) .* cos(x); % Function

f_analytical = @(x) -2*exp(x) .* sin(x); % Exact second derivative

x = pi/4; % Test point

h_values = logspace(-4, -1, 100); % Range of h values

% Initialize error arrays

central_error = zeros(size(h_values));

right_skewed_error = zeros(size(h_values));

left_skewed_error = zeros(size(h_values));

% Loop through different step sizes h

for i = 1:length(h_values)

h = h_values(i);

% Central difference 4th order approximation

central_diff = (-f(x - 2*h) + 16*f(x - h) - 30*f(x) + 16*f(x + h) - f(x + 2*h)) / (12*h^2);

central_error(i) = abs(central_diff - f_analytical(x));

% Right-skewed difference 4th order approximation

right_skewed_diff = (35*f(x) - 104*f(x + h) + 114*f(x + 2*h) - 56*f(x + 3*h) + 11*f(x + 4*h)) / (12*h^2);

right_skewed_error(i) = abs(right_skewed_diff - f_analytical(x));

% Left-skewed difference 4th order approximation

left_skewed_diff = (35*f(x) - 104*f(x - h) + 114*f(x - 2*h) - 56*f(x - 3*h) + 11*f(x - 4*h)) / (12*h^2);

left_skewed_error(i) = abs(left_skewed_diff - f_analytical(x));

end

% Plot the errors

figure;

loglog(h_values, central_error, 'r-', 'LineWidth', 1.5); hold on;

loglog(h_values, right_skewed_error, 'g--', 'LineWidth', 1.5);

loglog(h_values, left_skewed_error, 'b:', 'LineWidth', 1.5);

xlabel('Step size (h)');

ylabel('Absolute Error');

legend('Central Difference', 'Right-Skewed Difference', 'Left-Skewed Difference', 'Location', 'Best');

title('Error Comparison for 2nd Order Derivative Approximations');

grid on;

Explanation of Skewed Schemes:

Skewed schemes are useful in situations where boundary conditions limit the availability of data on one side of a point. For example, if we are near the boundary of a domain, we may not have access to data points on both sides of a point, making a central difference scheme unusable. In such cases, skewed (left-sided or right-sided) schemes are essential as they allow us to compute the derivative using points only on one side of the evaluation point.

CONCLUSION:

Boundary conditions: Skewed schemes are helpful near boundaries where data might not be available on both sides.

Flexibility: They provide flexibility in derivative approximations when the function or data is only available on one side of a point. Central difference schemes cannot be used in such cases.