function [I,m]=uniform_refinement(f,F,a,b)
% Run quadratures on uniform subdivisions and see how the error behaves
h = []; err = [];
for i=0:10,
    m  = 2^i; % number of subintervals
    xi = linspace(a,b,m+1); % endpoints of the intervals
    [I,E]=trapezoid(xi,f,F);
    h   = [h,xi(2)-xi(1)];
    err = [err,sum(E)];
end
%plot the error behaviour vs expected h^2
figure(2); clf;
fs = 24;
loglog(h,err, '-X',...
       h,h.^2,'-', ...
       'LineWidth',3,'MarkerSize',10);
xlabel('h','FontSize',fs); ylabel('error','FontSize',fs); 
handle = legend('Numerical error', 'h^2','Location','NorthWest');
set(handle,'FontSize',fs);
grid on;
% plot the error distribution for the last subdivision
figure(3); clf;
plot(0.5*(xi(1:end-1)+xi(2:end)), h(end)/(b-a)*E,'-x',...
     'LineWidth',3,'MarkerSize',10);
xlabel('x','FontSize',fs); ylabel('error','FontSize',fs); 
I = sum(I);