function I = trapzc(a,b,k,f)

% trapzc is the composite trapezoidal formulas. We divide the interval
% by half k times. 

h = (b-a);  I = (feval(f,a) + feval(f,b) )*h/2;     % One step trapezoidal rule.

for l=1:k
  h1 = h/2;
  HM = 0;  M = 2^(l-1);                % Evaluate the added points
  for j=1:M
    HM = HM + feval(f,a+(2*j-1)*h1);
  end
  HM = HM*h;
  I = (I + HM )/2; 
  h = h1;                  % Go to the next level.
end