Trapezoidal Rule

A brief introduction to the Trapezoidal rule and a uniform interval Composite Trapezoidal Rule implementation.

Trapezoidal Rule

The Trapezoidal Rule is used to evaluate a definite integral within some the limits and . Trapezoidal rule approximates with a polynomial of degree one ie. a line segment between the two points of the given function at with endpoints and , and then finds the integral of that line segment which is used to approximate the definite integral . The integral of
the line approximated function is simply the area of the trapezium formed by the straight line, with base width and bounded by the lines and in the positive side of the X axis. The approximate evaluation of the integral is found by finding the integral of the straight line with endpoints and can be. found with the below iterative formula:

The above iterative formula approximates with only one straight line. To represent the function more accurately with the straight lines using the Trapezoidal Rule, the function is divided in many intervals of the same uniformed width say , such that , and a line is fit within each interval to approximate the function and each such interval is evaluated with the above process and at last added to get the final integral. This is gives better results as the function is better approximated with the line segments within the small intervals, and the missed sections of a curves are better covered. This is known as the Trapezoidal Composite Rule. The iterative method for the Trapezoidal Composite rule could be derived from the above relation and is found to be:

The Process

The Trapezoidal Composite Rule can be implemented directly from the iterative formula presented above. The first and the last ordinate would be calculated separately and the other ordinates would be summed by a loop and then is multiplied by 2. Next these two parts are added together to get the integral. For demonstration purpose some sample functions are used.

1. 
2. 

FUNCTION trapezoidal(f, a, b, n)
  Find width of each interval h = (b-a)/n 
  Find y0= f(a), yn= f(b)
  Initialize x = a+h, i = 1, sum = 0
/*  Calculate the sum of the ordinates from (x + h) to (x + (n-1)h) */
  WHILE (i < n)
    DO
      sum = sum + f (x)
      x = x + h
  ENDWHILE
/*  Apply Trapezoidal Composite Formula to find the integral */
  Iz = (h/2) [(y0 + yn) + 2 * sum]
ENDFUNCTION

Sourcecode

#include <stdio.h>
#include <math.h>

double f1 (double x);
double f2 (double x);

double trapezoidal (double (*f) (double x), double a, double b, int n);

int
main (void)
{
  double a, b, n;
  double Iz;

  printf ("\nEnter a,b,n: ");
  scanf ("%lf %lf %lf", &a, &b, &n);

  printf ("\nf(x) = sin (2x) / (1+x)^5");
  /* Show integral computed with trapezoidal rule */
  Iz = trapezoidal (f1, a, b, n);
  printf ("\n\tI_Trapezoidal (f(x), %g, %g, %g) = %g", a, b, n,Iz);

  printf ("\n");

  printf ("\nf(x) = (1/x) + 5 + 10x^2;");
  /* Show integral computed with trapezoidal rule */
  Iz = trapezoidal (f2, a, b, n);
  printf ("\n\tI_Trapezoidal (f(x), %g, %g, %g) = %g", a, b, n,Iz);
 
  printf ("\n");

  return 0;
}

/* Takes a function pointer which should represent the single variable
 * function of which the integral is to be found. The function should 
 * receive a double and return a double. The limits a and b and the 
 * number of intervals n is provided. The value returned is the integral
 * of the passed function is the argument
 */
double
trapezoidal (double (*f) (double x), double a, double b, int n)
{
  double h;
  double y = 0, x, sum = 0, y0, yn;
  int i;

  /* Avoid calling NULL */  
  if (f == NULL)
    return 0;

  h = (b - a) / n;

  y0 = f (a);
  yn = f (b);

  for (i = 1, x = a + h; i < n; x = x + h, i++)
    sum = sum + f (x);

  y = (h / 2) * ((y0 + yn) + 2 * sum);

  return y;
}

/* Sample function 1 */
double
f1 (double x)
{
  return sin (2 * x) / pow ((1 + x), 5);
}

/* Sample function 2 */
double
f2 (double x)
{
  return (1 / x) + 5 + 10 * x * x;
}

Description

The two given function defined in the previous section is implemented with the functions f1 () and f2 (). These two functions accept a double type data, and returns double.

• double f1 (double x) : This function implements the equation sin (2 * x) / pow ((1 + x), 5) and returns the evaluated value.

• double f2 (double x) : This function implements the equation (1 / x) + 5 + 10 * x * x and returns the evaluated value.

• double trapezoidal (double (*f) (double x), double a, double b, int n) : This function evaluates a definite integral integral of a function with Trapezoidal method of a function f within limits a, and b with n equal intervals. This takes a pointer to a function which accepts one double type data and returns a double type data. The next argument is double a, which is interpreted as the lower limit of the integral, the next argument double b is considered as the upper limit of the integral. int n is the number of intervals to be taken. After finding the width of each interval, this first calculates the first and the last ordinate y0 = f(a) and y1 = f(b). Then it enters the for loop and computes the (n-2) values of the passed
  function starting from (x+h) up to (x+(n-1)h). Then it uses the composite formula to find the integral with Trapezoidal rule, and returns the value: (h / 2) * ((y0 + yn) + 2 * sum);

• int main (void) : The main function prompts the user to enter the upper and the lower limits and the number of intervals to be taken and calls the functions with the proper parameters. The trapezoidal() is called as Iz = trapezoidal (f1, a, b, n); which computes the integral of the function f1 within upper and lower limit a and b with the Trapezoidal Rule.

Error

When is continuous and M is the upper bound for the values of
on , then the error induced by the
Trapezoidal rule is found to be:

Although the theory tells us there will always be a smallest safe value of M, in practice we can hardly ever find it. Instead, we find the best value we can and go on from there to estimate . To make small we make small. It also can be graphically seen that a line fits less accurately for a certain interval in a higher degree function than a lower degree function. So increasing the interval h will result in better results, or implementing better fitting curve than a line needs to be selected.

Comments

The Trapezoidal evaluation function is implemented to be independent of the function of which the integral is to be found, by allowing passing of the function to be integrated as a parameter. This implementation could be incorporated into an algebraic equation parser, which can be used to evaluate some user entered functions to be evaluated like for example:
ans = trapezoidal (parse_expression (user_expression, x), a, b, n) ; the ‘user_expression’ is being considered a single variable function with the variable replaced with the value of ‘x’ .

Output

Some sample outputs of the above code is shown below

Run 1: input = lower limit a = 1, upper limit b = 2, intervals
n = 10

Enter a,b,n: 1 2 10

f(x) = sin (2x) / (1+x)^5
	I_Trapezoidal (f(x), 1, 2, 10) = 0.00513881

f(x) = (1/x) + 5 + 10x^2;
	I_Trapezoidal (f(x), 1, 2, 10) = 29.0438

Run 2: input = lower limit a = 5.3, upper limit b = 10.23,
intervals n = 100

Enter a,b,n: 5.3 10.23 100

f(x) = sin (2x) / (1+x)^5
	I_Trapezoidal (f(x), 5.3, 10.23, 100) = -3.20357e-05

f(x) = (1/x) + 5 + 10x^2;
	I_Trapezoidal (f(x), 5.3, 10.23, 100) = 3097.73

Links

1. Check for more: http://en.wikipedia.org/wiki/Trapezoidal_rule
2. Simpson’s 1/3^rd rule: https://phoxis.org/2011/03/11/simpson-13-rd-rule/

References

1. Calculus and Analytic Geometry : Thomas Finney
2. Images from : Wikipedia

: