Simpson’s 1/3 rule

In Simpson’s 1/3 rule the interval [a, b] is divided into two equal sub-intervals by the points x₀, x₁, x₂, where h = (b- a)/2 and x₂ = x1 +h.

In the previous article we generate Trapezoidal Rule from the general integration formula based on Newton’s forward interpolation formula

\[\int_{a}^{b}{f\left( x \right)dx}=nh\left[ {{y}_{0}}+\frac{n}{2}\Delta {{y}_{0}}+\frac{2{{n}^{2}}-3n}{12}{{\Delta }^{2}}{{y}_{0}}+\frac{{{n}^{3}}-4{{n}^{2}}+4n}{24}{{\Delta }^{3}}{{y}_{0}}+…. \right]……..(1)\]

Simpson’s 1/3 rule is obtained by putting n = 2 in (1). In this case, the third and higher order differences do not exist.

\[\int_{a}^{b}{f\left( x \right)dx}=\int_{{{x}_{0}}}^{{{x}_{n}}}{f\left( x \right)dx}=2h\left[ {{y}_{0}}+\Delta {{y}_{0}}+\frac{1}{6}{{\Delta }^{2}}{{y}_{0}} \right]\]

\[\Rightarrow \int_{a}^{b}{f\left( x \right)dx}=2h\left[ {{y}_{0}}+\left( {{y}_{1}}-{{y}_{0}} \right)+\frac{1}{6}\left( {{y}_{2}}-2{{y}_{1}}+{{y}_{0}} \right) \right]\]

\[\therefore \int_{a}^{b}{f\left( x \right)dx}=\frac{h}{3}\left[ {{y}_{0}}+4{{y}_{1}}+{{y}_{2}} \right]\]

The above rule is known as Simpson’s 1/3 rule or simply Simpson’s rule.

Composite Simpson’s 1/3 rule

Let the interval [a, b] be divided into n (an even number) equal sub-intervals by the points x₀, x₁, x₂, …, x_n, where x_i = x₀ + ih, i = 1, 2, …, n. Then

\[\int_{a}^{b}{f\left( x \right)dx}=\int_{{{x}_{0}}}^{{{x}_{2}}}{f\left( x \right)dx}+\int_{{{x}_{2}}}^{{{x}_{4}}}{f\left( x \right)dx}+….+\int_{{{x}_{n-2}}}^{{{x}_{n}}}{f\left( x \right)dx}\]

\[=\frac{h}{3}\left[ {{y}_{0}}+4{{y}_{1}}+{{y}_{2}} \right]+\frac{h}{3}\left[ {{y}_{2}}+4{{y}_{3}}+{{y}_{4}} \right]+….+\frac{h}{3}\left[ {{y}_{n-2}}+4{{y}_{n-1}}+{{y}_{n}} \right]\]

\[=\frac{h}{3}\left[ {{y}_{0}}+4\left( {{y}_{1}}+{{y}_{3}}+…+{{y}_{n-1}} \right)+2\left( {{y}_{2}}+{{y}_{4}}+…+{{y}_{n-2}} \right)+{{y}_{n}} \right]\]

This formula is known as Simpson’s 1/3 composite rule for numerical integration.

Error in Simpson’s 1/3 rule

The error committed in the formula is given by

\[{{E}_{1}}=\int_{{{x}_{0}}}^{{{x}_{2}}}{f\left( x \right)dx}-\frac{h}{3}\left[ {{y}_{0}}+4{{y}_{1}}+{{y}_{2}} \right]\simeq -\frac{1}{90}{{h}^{5}}{{f}^{iv}}\left( {{x}_{0}} \right)\]

which is the error in the interval [x₀, x₂].

The total error committed in composite Simpson’s 1/3 rd rule is given by

\[E\simeq -\frac{1}{90}{{h}^{5}}\frac{n}{2}{{f}^{iv}}\left( \xi  \right)=-\frac{b-a}{180}{{h}^{4}}{{f}^{iv}}\left( \xi  \right),{{x}_{0}}<\xi <{{x}_{n}},\left( \because b=a+nh \right)\]

where ξ is the point having largest value of the fourth order derivatives, and since the number of sub-intervals is n/2.

Algorithm of Simpson’s 1/3 rule

Step 1. Start;
Step 2. Input function f(x);
Step 3. Read a,b,n; // the lower and upper limits and number of sub-intervals
Step 4. Compute h=(b-a)/n;
Step 5. Sum = [f(a)-f(a+nh)];
Step 6. for i=1 to n-1 step 2 do
			Compute sum = sum + 4*f(a+ih)+2*f(a+(i+1)h);
		endfor;
Step 7. Compute result = sum * h/3;
Step 8. Print result;
Step 9. Stop;

Simpson’s 1/3 rule implementation in C

/* Program Simpson’s 1/3 rule
   Program to find the value of integration of a function
   f(x) using Simpson’s 1/3 rule. Here we assume that 
   f(x) = x^3. */

#include<stdio.h>

void main()
{
	float f(float);
	float a,b,h,sum;
	int i,n;
	printf("Enter the values of a,b: ");
	scanf("%f%f",&a,&b);
	printf("Enter the value of n: ");
	scanf("%d",&n);
	if(n%2!=0)
	{
		printf("\nNumber of subdivision should be even");
		exit(0);
	}
	h=(b-a)/n;
	sum = f(a)-f(a+n*h);
	for(i=1;i<n;i++)
		sum += 4*f(a+ih)+2*f(a+(i+1)*h);
	sum *= h/3;
	printf("\nValue of the integration is %f",sum);
}

float f(float x)
{
	return(x*x*x);
}

Output

Enter the values of a, b: 0 1

Enter the value of n: 100

Value of the integration is 3.750000

Find the value of

\[\int_{1}^{5}{{{\log }_{10}}xdx}\]

Taking 8 sub-intervals, correct up to four decimal places by Simpson’s 1/3^rd rule.

Solution:

\[Let~~f\left( x \right)={{\log }_{10}}x\]

\[Here~~{{x}_{0}}=1,~~{{x}_{n}}=5,~~n=8\]

\[\therefore h=\frac{{{x}_{n}}-{{x}_{0}}}{n}=\frac{5-1}{8}=0.5\]

The tabulated values of f(x) for different values of x are given below,

By Simpson’s 1/3^rd rule, we have,

\[\int_{a}^{b}{f\left( x \right)dx}=\frac{h}{3}\left[ {{y}_{0}}+4\left( {{y}_{1}}+{{y}_{3}}+{{y}_{5}}+{{y}_{n7}} \right)+2\left( {{y}_{2}}+{{y}_{4}}+{{y}_{6}} \right)+{{y}_{8}} \right]\]

\[\therefore \int_{1}^{5}{{{\log }_{10}}xdx}=\frac{0.5}{3}[0+4\left( 0.17609+0.39794+0.54407+0.65321 \right)\]

\[+2\left( 0.30103+0.47712+0.60206 \right)+0.69897]\]

\[\therefore \int_{1}^{5}{{{\log }_{10}}xdx}=1.75744\simeq 1.7574\]

Correct upto four decimal places.

Find the value of

\[\int_{1.2}^{1.6}{\left( x+\frac{1}{x} \right)dx}\]

Taking 4 sub-intervals, correct up to four significant figures.

Solution:

\[Let~~f\left( x \right)=\left( x+\frac{1}{x} \right)\]

\[Here~~{{x}_{0}}=1.2,~~{{x}_{n}}=1.6,~~n=4\]

\[\therefore h=\frac{1.6-1.2}{4}=0.1\]

The tabulated values of f(x) for different values of x are given below:

By Simpson’s 1/3^rd rule, we have,

\[\int_{a}^{b}{f\left( x \right)dx}=\frac{h}{3}\left[ {{y}_{0}}+4\left( {{y}_{1}}+{{y}_{3}} \right)+2\left( {{y}_{2}} \right)+{{y}_{4}} \right]\]

\[\int_{1.2}^{1.6}{\left( x+\frac{1}{x} \right)dx}=\frac{0.1}{3}\left[ 2.033333+4\left( 2.069231+2.166667 \right)+2\left( 2.114284 \right)+2.225 \right]\]

\[\therefore \int_{1.2}^{1.6}{\left( x+\frac{1}{x} \right)dx}=0.847683\simeq 0.8477\]

<< Previous     Next>>

;