% The quarter circle (integrand circle.m) also exhibited slow convergence.
% We can improve that, to calculate pi faster!

format long  

O = '---- parabola ---- = O'
Nfor = 1e2;
pow = [1,2,3,4,5,6];
Nsteps = (10.^pow).*Nfor
Ntrials = numel(pow) 
for i = 1:Ntrials
   t0 = cputime;
   results(i) = midpoint2(@parabola, 0, 1, Nfor, 10^pow(i));
   execution_time = cputime - t0
end
errors = results(2:end) - results(1:end-1);
   % this error estimate presumes that there is much less error in
   % the next trial, which has 10 times the number of points evaluated. 

hold off
plot(log10(Nsteps(1:end-1)), log10(abs(errors)), 'o-k');
xlabel('log10 number of steps', 'fontsize', 16);
ylabel('log10 absolute error', 'fontsize', 16);
title('midpoint rule convergence: calculating \pi', 'fontsize', 18);
hold on


O = '---- circle ---- = O'
for i = 1:Ntrials
   t0 = cputime;
   results(i) = midpoint2(@circle, 0, 1, Nfor, 10^pow(i));
   execution_time = cputime - t0
end
errors = results(2:end) - results(1:end-1);
   % this error estimate presumes that there is much less error in
   % the next trial, which has 10 times the number of points evaluated. 
estimated_error = errors(end)
estimate_pi = 4*results(end)
actual___pi = pi
abs_error = estimate_pi - pi
plot(log10(Nsteps(1:end-1)), log10(abs(errors)), '+-b');


O = '---- circle2 (transformed from circle) ---- = O'
for i = 1:Ntrials
   t0 = cputime;
   results(i) = midpoint2(@circle2, 0, 1, Nfor, 10^pow(i));
   execution_time = cputime - t0
end
errors = results(2:end) - results(1:end-1);
   % this error estimate presumes that there is much less error in
   % the next trial, which has 10 times the number of points evaluated. 
estimated_error = errors(end)
estimate_pi = 4*(results(end) + 2*sqrt(2)/3)
actual___pi = pi
abs_error = estimate_pi - pi
O = '---- Calculations complete! ---- = O'

plot(log10(Nsteps(1:end-1)), log10(abs(errors)), '^-g');

legend('parabola','circle','circle2');

print('integration_test7.pdf','-dpdf')

hold off