Numerical Integration Formulas

Here are the M-files to implement composite trapezoidal rule for equally spaced data and unequally spaced data. 

Composite Trapezoidal Rule (Equally Spaced Data) M-file

function I = trap(func,a,b,n,varargin)

% trap: composite trapezoidal rule quadrature

% I = trap(func,a,b,n,pl,p2,...):

% composite trapezoidal rule

% input:

% func = name of function to be integrated

% a, b = integration limits

% n = number of segments (default = 100)

% pl,p2,... = additional parameters used by func

% output:

% I = integral estimate

if nargin<3,error('at least 3 input arguments required'),end

if ~(b>a),error('upper bound must be greater than lower'),end

if nargin<4|isempty(n),n=100;end

x = a; h = (b - a)/n;

s=func(a,varargin{:});

for i = 1 : n-1

x = x + h;

s = s + 2*func(x,varargin{:});

end

s = s + func(b,varargin{:});

I = (b - a) * s/(2*n);

Composite Trapezoidal Rule (Unequally Spaced Data) M-file

function I = trapuneq(x,y)

% trapuneq: unequal spaced trapezoidal rule quadrature

% I = trapuneq(x,y):

% Applies the trapezoidal rule to determine the integral

% for n data points (x, y) where x and y must be of the

% same length and x must be monotonically ascending

% input:

% x = vector of independent variables

% y = vector of dependent variables

% output:

% I = integral estimate

if nargin<2,error('at least 2 input arguments required'),end

if any(diff(x)<0),error('x not monotonically ascending'),end

n = length(x);

if length(y)~=n,error('x and y must be same length'); end

s = 0;

for k = 1:n-1

s = s + (x(k+l)-x(k))*(y(k)+y(k+l))/2;

end

I = s;