CodeCogs

Computes the definite integral of a function using the trapezoidal rule.

You're viewing an older version of this page (#385). View the current version.

Interface

#include <codecogs/maths/calculus/quadrature/trapezoidal.h>

using namespace Maths::Calculus::Quadrature;

double trapezoidal(int n, double (*f)(double), double a, double b)
double trapezoidal(int n, const double *values, double a, double b)

Overview

This module computes the area beneath either a user supplied function or a set of discrete points, using an approximation which assumes the function is linear between each two consecutive points.

Consider a two times differentiable function $f:I \subset \mathbb{R} \rightarrow \mathbb{R}$ , two distinct abscissas $a < b \in I$ and a positive integer $n$ . Then the following approximation holds:

$\int_a^b f(x) \mathrm{d}x \approx \frac{b-a}{n} \left[ \frac{f(a) + f(b)}{2} + \sum_{k=1}^{n-1} f \left(a + \frac{b-a}{n} k\right) \right]$

(1)

with an error bound of:

$\epsilon \leq \frac{(b-a)^3}{12n^2} \sup_{x \in [a,b]} \left|\frac{d^2}{dx^2}f(x)\right|.$

(2)

This is a quadrature formula known as the trapezoidal rule. From the formula of the error bound you may notice that e.g. if one doubles the number of points $n$ , the approximation error is decreased four times. To better understand how this rule works, observe the following picture - the area between abscissas $x_1$ and $x_2$ is approximated by a trapezoid.

Below you will find examples that show how to use each version of the function and the error estimates in each case.

References

Mihai Postolache - "Metode Numerice", Editura Sirius
Mathworld, http://mathworld.wolfram.com/TrapezoidalRule.html

FUNCTION

trapezoidal

This version is to be used when you wish to pass the function to integrate as an argument, instead of passing an array with its values at distinct points.

Example 1

In what follows an approximation is found for the definite integral

$\int_1^3 \sin x \,dx$

(3)

and the absolute error from its actual value is estimated.

#include <codecogs/maths/calculus/quadrature/trapezoidal.h>
#include <stdio.h>
#include <math.h>

// number of points
#define N 100

// function to integrate
double f(double x)
{
  return sin(x);
}

// the primitive of f, to estimate errors
double pf(double x)
{
  return -cos(x);
}

int main()
{
  // compute the approximate area
  double fi = Maths::Calculus::Quadrature::trapezoidal(N, f, 1, 3),

  // use the Leibniz-Newton formula to find a more precise estimate
  realfi = pf(3) - pf(1);

  // display problem data
  printf("      f(x) = sin(x)\n");
  printf("    points = %d\n\n", N);

  // display the result and error estimate
  printf("   I(1, 3) = %.15lf\n", fi);
  printf("real value = %.15lf\n", realfi);
  printf("     error = %.15lf\n\n", fabs(fi - realfi));

  return 0;
}

Output

f(x) = sin(x)
    points = 100

   I(1, 3) = 1.530243792301766
real value = 1.530294802468585
     error = 0.000051010166819

Parameters

the number of sample points of the function f, from which to approximate

the function to integrate

the inferior limit of integration

the superior limit of integration

Returns

The definite integral of the given function from a to b.

FUNCTION

trapezoidal

This version is to be used when the function is not known analytically, but a table with its values at equally spaced abscissas is available.

Example 1

The code below gives an approximation for the definite integral

$\int_1^{1.3} \sqrt{x} \,dx$

(4)

and estimates the absolute error from its actual value.

#include <codecogs/maths/calculus/quadrature/trapezoidal.h>
#include <stdio.h>
#include <math.h>

// the primitive of f, to estimate errors
double pf(double x)
{
  return 2*sqrt(x*x*x)/3;
}

int main()
{
  // values of the function at equally spaced abscissas
  double ordinates[7] = {1, 1.0247, 1.04881, 1.07238, 
  1.09544, 1.11803, 1.14017},

  // compute the approximate area
  fi = Maths::Calculus::Quadrature::trapezoidal(7, ordinates, 1, 1.3),

  // use the Leibniz-Newton formula to find a more precise estimate
  realfi = pf(1.3) - pf(1);

  // display problem data
  printf("      f(x) = sqrt(x)\n");
  printf("    points = 7\n\n");

  // display the result and error estimate
  printf(" I(1, 1.3) = %.15lf\n", fi);
  printf("real value = %.15lf\n", realfi);
  printf("     error = %.15lf\n\n", fabs(fi - realfi));

  return 0;
}