%Lab 2
%Written by Jeannette Spuhler, 2011-09-08
% Modified 2011-11-03


%Adaptive trapezoidal rule
%Pohl page 185
b= 10;
a1= -10;
a2= -5;
a3= -3.5;
n1=4;
n2=200;
n3=4;
h1 = (a2-a1)/n1;
h2= (a3-a2)/n2;
h3=(b-a3)/n3;

%first interval
x1=a1:h1:a2;
f1=peakfkt(x1);
T1 = h1*(sum(f1)-0.5*(f1(1)+f1(n1+1)));

%second interval
x2=a2:h2:a3;
f2=peakfkt(x2);
T2 = h2*(sum(f2)-0.5*(f2(1)+f2(n2+1)));

%third interval
x3=a3:h3:b;
f3=peakfkt(x3);
T3 = h3*(sum(f3)-0.5*(f3(1)+f3(n3+1)));

Ttotal= T1+T2+T3;

%plot function
plot(x1, f1); hold on;
plot(x2,f2); 
plot(x3, f3);
title('Lab2 - With Adaptivity'); 
ylabel('0.51+25*exp(-81*(11*x+45).^2)');
xlabel('x');
for i=1:n1
line([x1(i) x1(i)] ,[0 f1(i)]);
end

for i=1:n2
 line([x2(i) x2(i)] ,[0 f2(i)]);
end

for i=1:n3
 line([x3(i) x3(i)] ,[0 f3(i)]);
end


disp('Integral calculated with adaptivity: ');
disp(Ttotal);
