%simpson.m implements the Composite Simpson's Rule
%I = simpson(fun,a,b,npanel)
%
% Input:     fun    = (string) name of m-file that evaluates f(x)
%            a, b   = lower and upper limits of the integral
%            npanel = number of panels to use in the integration
%                     Total number of nodes = 2*npanel + 1
%
% Output:    I = approximate value of the integral from a to b of f(x)*dx

%Recreate example 11.7 from slide 20
func = inline('exp(-x.^2)','x');
%We want 2 panels, 5 points imply 2 applications of Simpson's Rule and
    %hence 2 Panels.
nPanels = 2;
fprintf('Example 11.7 recreated.\n');
fprintf('Composite Simpsons Rule with %d Panels gives us %f',nPanels, (2/sqrt(pi)) * simpson(func,0,1,nPanels));

%Watch the composite Simpson's Rule Converge to a value as nPanels
    %increases
fprintf('\n\nWatch the composite Simpsons Rule Converge to a value as nPanels increases.\n')
for i = 1:10
    nPanels = i;
    fprintf('Composite Simpsons Rule with %d Panels gives us %f\n',nPanels, (2/sqrt(pi)) * simpson(func,0,1,nPanels));
end

fprintf('\n\n\nRun NMMs demo.\n');
%To Run NMM's demo, just type
demoSimp


%Recreate the Nice Pictures of the curve we're integrating and the interpolant Simpson's Rule fits
	%try with different functions
%plotSimpInt(fun,0,5,3)