% Simple program to implement the trapezoidal rule
clear
a = 0;
b = 2;
n = input('Input n ');
h = (b-a)/n; % Could also do h = x(2)-x(1), but this is normally only
% done when one of a, b or n is not known.
% Tabulate the function
x = linspace(a,b,n+1)'; % Note that we are generating n + 1 points.
y = x.^4.*exp(-x);
total = 0; % Set accumulation variable to 0
for i = 1:n+1
if(i == 1 | i == n+1)
total = total + y(i); % The endpoints get a weight of 1
else
total = total + 2*y(i); % The interior nodes get a weight of 2
end
end
approx = h*total/2; % Compute integral approximation
exact = 24 - 168/(exp(1)^2);
relerr = (approx-exact)/exact