% CLASS10OCT  Show commands for using Matlab to do trapezoid rule on integral
%    /5
%    |  sin(x) dx = 1 - cos(5) = 0.716337814536774
%    /0

format long

% one subinterval:
0.5 * (5-0) * (sin(0) + sin(5))

% ten subintervals, but done the dumb way with lots of typing
dx=5/10;
T10 = 0.5 * dx * ( sin(0) + 2*sin(0.5) + 2*sin(1) + 2*sin(1.5) + 2*sin(2) + 2*sin(2.5) + 2*sin(3) + 2*sin(3.5) + 2*sin(4) + 2*sin(4.5) + sin(5) )

% correct answer
I = 1 - cos(5)

% error
I - T10

% less dumb way
x = 0:dx:5;  %  equivalently:  x = linspace(0,5,11);
T10 = 0.5 * dx * [1 2 2 2 2 2 2 2 2 2 1]*sin(x)'
I - T10

% one-liners
n=10;  dx = 5/n;  x = linspace(0,5,n+1);  Tn = 0.5 * dx * [1 2*ones(1,n-1) 1] * sin(x)'

% we showed in class that n=103 subintervals would give error less than 10^{-3}
n=103;  dx = 5/n;  x = linspace(0,5,n+1);  Tn = 0.5 * dx * [1 2*ones(1,n-1) 1] * sin(x)'

% the error
I - Tn

% then we showed n=32275 would give error less than 10^{-8}
n=32275;  dx = 5/n;  x = linspace(0,5,n+1);  Tn = 0.5 * dx * [1 2*ones(1,n-1) 1] * sin(x)'
I - Tn

% now I switch to Simpson's rule
% it has a different pattern
format short
n=10; [1 repmat([4 2],1,(n-2)/2) 4 1]
format long

% far fewer points are needed for Simpson's to get the same accuracy as Trapezoid
n=3000;  dx = 5/n;  x = linspace(0,5,n+1);  Sn = (1/3) * dx * [1 repmat([4 2],1,(n-2)/2) 4 1] * sin(x)'
I-Sn
n=130;  dx = 5/n;  x = linspace(0,5,n+1);  Sn = (1/3) * dx * [1 repmat([4 2],1,(n-2)/2) 4 1] * sin(x)'
I-Sn

