function t = comptrap(f,a,b,n)
%
%  function t = comptrap(f,a,b,n)
%
%  Applies the composite trapezoidal rule to approximate the integral
%  of the function f over the interval [a,b], using n subintervals.
%
%  Example:
%     t = comptrap('sin',0,pi/2,128)
%
%  Note: compare this to the result of MATLAB built-in command quad:
%     quad('sin',0,pi/2)
%
h = (b-a)/n;
t = (feval(f,a) + feval(f,b))/2;
for i = 1:n-1
   t = t + feval(f,a+i*h);
end;
t = t*h;