function y = traprule(f,a,b)
%
% y = traprule(f,a,b) returns the value
% of the trapezoidal rule to approximate
% the integral of f over the interval [a,b].
%
fa = feval(f,a);
fb = feval(f,b);
y = (b-a)*(fa+fb)/2;