% romberg2run
% Script to run the new version romberg2, for the Romberg method for integration.
% Author: Brenton LeMesurier
% Date: March 8, 2007
% Version for MATLAB7 only: function input as a function pointer.
format compact
format short
clc % clear the screen
% The integral $\int_1^3 f(x) dx$ has exact value 1; easy to check accuracy.
f=@(x) exp(x)/(exp(3)-exp(1))
a=1
b=3
errortolerance=1e-10
levelsmin=4 % always do at least four levels (16 intervals)
levelsmax=8
[Iabf,errorestimate,R]=romberg2(f,a,b,errortolerance,levelsmax,levelsmin);
fprintf('The final approximate integral is %g, with estimated absolute error %g\n',Iabf,errorestimate)
disp('Table of successive approximations')
[rows,columns]=size(R);
disp('Intervals Trapezoid      Simpson''s      other extrapolations')
for row=1:rows
    fprintf('%-10d',2^(row-1)), fprintf('%-15.8g', R(row,1:row)), fprintf('\n')
end

% Since we know the exact answer, here are the errors.
% (For functions with unknown integral, use Iabf in place of the exact value.)
Iabfexact=1; % f = exp(x)/(exp(3)-exp(1))
disp('Table of absolute errors')
for row=1:rows
    fprintf('%-10d',2^(row-1)), fprintf('%-15.3e', abs(Iabfexact-R(row,1:row))), fprintf('\n')
end