x = 0:pi/100:pi;
y = sin(x);
trapz(x,y)

ys = cumtrapz(x,y);
plot(x,ys)
xlabel('position x'); ylabel('integral of sin(x)')
grid on

Nrun = 10;
n = zeros(1,Nrun);
err = zeros(1,Nrun);

for k = 1:Nrun
  n(k) = 10*2^(k-1);
  x = linspace(0,pi,n(k)+1);
  err(k) = abs(2-trapz(x,sin(x)));
end%for
loglog(n,err);

F = [ones(size(n))' log(n')];
LinearRegression(F,log(err'))

errTrap = err;

for k = 1:Nrun
  n(k) = 10*2^(k-1);
  err(k) = abs(2-simpson(@sin,0,pi,n(k)));
end%for

loglog(n,err,'-+',n,errTrap,'-*');
xlabel('number of grid points'); ylabel('integration error')
grid on
legend('Simpson','trapezoidal')

F = [ones(size(n))' log(n')];
LinearRegression(F,log(err'))
