function s = comp_trap(f,a,b,n)

%    A MATLAB CODE PERFORMING COMPOSITE TRAPEZOIDAL INTEGRATION
%
% INPUT:     f = the integrand input as a string; 
%                if the integrand is, say, f(x) = sin(x), 
%                give as first argument of comp_simp the expression
%                        inline('sin(x)')
%            a = lower limit of integration
%            b = upper limit of integration
%            n = number of intervals in which [a,b] is divided 
%
% OUTPUT:    The value of the definite integral

format long

h = (b-a)/n      %  Stepsize
sum  = 0.0;

for k = 1:(n-1)
  x = a + k*h;
  sum = sum + feval(f,x);
end

s = h * ( feval(f,a) + 2 * sum + feval(f,b) ) / 2.0;
