function fig6_6

%  this includes Fig 6.6 and Table 6.2

%  Integrate  f(x)  over  [a, b]  using Trap rule
%
%  Inputs:
%        a = left endpoint of interval
%        b = right endpoint of interval
%        n = number of subintervals

a=0; b=1;
exact=(exp(3)-1)/3;

%  interv = number of subinternals
interv=[1000 2000 4000 8000 16000 32000];
nk=length(interv);
for k=1:nk
    n=interv(k);
    T(k)=trap(a,b,n);
    err(k)=abs(exact-T(k));
end

% table 6.2
fprintf('\n RIE = relative iterative error \n\n')
fprintf(' k \t  n  \t\t I_M \t \t Error \t\t RIE \n')
fprintf(' %i \t %i\t %19.15f \t %5.1e  \n',1,interv(k),T(1),err(1))
%fprintf('\\, %i \\, & \\, %i \\, & \\,\\, %19.15f \\,\\, & \\,\\, %5.1e  \\,\\, & \\,\\,  \\\\ \n',1,interv(1),T(1),err(1))
for k=2:nk
    fprintf(' %i \t %i \t %19.15f \t %5.1e  \t %5.1e   \n',k,interv(k),T(k),err(k),abs((T(k)-T(k-1))/T(k)))
    %fprintf('\\, %i \\, & \\, %i \\, & \\,\\, %19.15f \\,\\, & \\,\\, %5.1e  \\,\\, & \\,\\, %5.1e  \\\\ \n',k,interv(k),T(k),err(k),abs((T(k)-T(k-1))/T(k)))
end
fprintf('\n')

%%% fig 6.6
figure(1)
clf
% get(gcf)
set(gcf,'Position', [25 1115 658 230])
subaxis(1,1,1,1,'MT',0.01,'MB',0.18,'MR',-0.02,'ML',0.03,'P',0.04)
co = [0 0 1;
    0 0.5 0;
    1 0 0;
    0 0.75 0.75;
    0.75 0 0.75;
    0.75 0.75 0;
    0.25 0.25 0.25];
set(groot,'defaultAxesColorOrder',co)

hold on
n=3;
xd=linspace(a,b,n+1);
h=xd(2)-xd(1);
for in=1:n
    plot([xd(in) xd(in+1)], [f(xd(in)) f(xd(in+1))],'--r','LineWidth',2)
    plot([xd(in) xd(in)], [0 f(xd(in))],'--r','LineWidth',2)
end
plot([0.999*b 0.999*b], [0 f(b)],'--r','LineWidth',2)

nx=100;
x=linspace(a,b,nx);
for ix=1:nx
    y(ix)=f(x(ix));
end
plot(x,y,'LineWidth',1.6)
box on
%grid on
xlabel('x-axis','FontSize',16)
ylabel('y-axis','FontSize',16)
axis([0 1 0 20])
set(gca,'XTick',[0   1/3   2/3   1])
set(gca,'XTickLabel',{'0','1/3','2/3','1'})
%set(gca,'YMinorGrid','off');
set(gca,'FontSize',16,'FontWeight','bold')

%exportgraphics(gca,'/Users/mark/Desktop/expT.eps')


function y=f(x)
y=exp(3*x);

%%% composite trap: n = number of subintervals
function s=trap(a,b,n)
xd=linspace(a,b,n+1);
h=xd(2)-xd(1);
s=0.5*f(xd(1));
for j=2:n
    s=s+f(xd(j));
end
s=h*(s+0.5*f(xd(n+1)));












