% 4-23-2008
% Quad cumulative integration example

clear all
close all

% This is our sample grid size
dx = 0.1;
x=-10:dx:10;

% Trapz example
% Integrates sech through the given points
y=sech(x).^2;
trapz(x,y)

% Cumtrapz example
% Plots the cumulative integral (indefinite integral, kind of)
plot(x,cumtrapz(x,y),'k.')
hold on

% Cumulative quad example
% Loop through our grid integrating up to the current value of x(i) each
% time and storing it in the vector int_sech.  Note that we have to start
% at 2 instead of 1 since quad complains about an interval length of 0,
% i.e. since [-10, x(1)] = [-10,-10]
for i=2:length(x)
    int_sech(i-1) = quad(inline('sech(x).^2'),-10,x(i));
end

plot(x(2:end),int_sech,'ro')
title('Cumulative integral of sech(x)^2')
xlabel('x')
ylabel('\int sech(x)^2 dx')
legend('cumtrapz','quad')