Saturday, March 25, 2023

Simpson's 1/3 rule of numerical integration using MATLAB

 Simpson's 1/3 rule

 Simpson's 1/3 rule, also known as Simpson's rule, is a numerical method used to approximate the definite integral of a function. It is based on the idea of approximating the area under the curve by fitting parabolic segments to the curve.

To apply Simpson's 1/3 rule, we first divide the interval of integration [a, b] into an even number of subintervals of equal width, h = (b - a) / n, where n is an even number. The width of each subinterval is the same, so the distance between adjacent points where we evaluate the function is also the same.

Next, we approximate the area under the curve using parabolic segments. We fit a parabolic segment to the curve over every two adjacent subintervals. The parabolic segment passes through the endpoints of the two subintervals and the midpoint of the interval. The area under the parabolic segment represents an approximation of the area under the curve over the two subintervals.

The formula for Simpson's 1/3 rule is given by:

∫[a,b] f(x) dx ≈ h/3 * [f(a) + 4f(a+h) + 2f(a+2h) + 4f(a+3h) + ... + 2f(b-h) + 4f(b-2h) + f(b)]

where h is the width of each subinterval (h = (b - a) / n), and f(x) is the function being integrated.

In this formula, the term f(a) represents the value of the function at the left endpoint of the interval [a, b], and f(b) represents the value of the function at the right endpoint of the interval. The terms 4f(a+h), 4f(a+3h), 4f(a+5h), ..., 4f(b-2h) represent the values of the function at the endpoints of the even-numbered subintervals, multiplied by a factor of 4 because they contribute to parabolic segments. The terms 2f(a+2h), 2f(a+4h), 2f(a+6h), ..., 2f(b-h) represent the values of the function at the endpoints of the odd-numbered subintervals, multiplied by a factor of 2 because they do not contribute to parabolic segments.

The error of Simpson's 1/3 rule is of order O(h^4), which means that the error decreases as the width of the subintervals decreases. However, Simpson's 1/3 rule may not be accurate enough for some functions, and other numerical methods may be required. For example, Simpson's 3/8 rule is a modification of Simpson's 1/3 rule that uses three subintervals instead of two, and it has an error of order O(h^5). It can be summarized as:






Now let us solve using Matlab: Here is the code:

% Simpson 1/3 rule
% Define the function to integrate
f = @(x) 0.2+ 25*x-200*x.^2 + 675*x.^3-900*x.^4+400*x.^5;
% Define the integration limits
a = 0;
b = 0.8;
% Define the number of intervals (must be even)
n = 4;
% Calculate the interval width
h = (b-a)/n;
% Evaluate the function at the interval points
x = a:h:b;
y = f(x);
% Apply the Simpson's 1/3 rule
I = h/3 * (y(1) + 4*sum(y(2:2:end-1)) + 2*sum(y(3:2:end-2)) + y(end));
% Display the result
disp(['Approximate integral: ', num2str(I)]);


The above code gives the following result for the Simpsons 1.3 rule of integration
==========================================================

simpson_1_3
Approximate integral: 1.623

No comments:

Post a Comment