function int = trap_td(f,a,b,n)
% trapezoidal rule
h = (b-a)/n;
x=a:h:b;
int = 0;
for j=2:n
  int = int + f(x(j));
end  
int = 2*int;
int=int +f(a)+f(b);
int=0.5*h*int;
end

