calculate error of an integral using repeat simpsons rule in matlab

main.m
error = (1/15) * abs(S(h/2) - S(h))
36 chars
2 lines

main.m
function [I, error] = repeat_simpsons_rule(f, a, b, N)
% f: function to integrate
% a: lower limit of integration
% b: upper limit of integration
% N: number of steps

% Apply Simpson's rule with step size h
h = (b - a) / (2*N);
x = a:h:b;
y = f(x);

% Calculate the integral and the error
I = h/3 * (y(1) + 4*sum(y(2:2:end-1)) + 2*sum(y(3:2:end-2)) + y(end));
error = (1/15) * abs(h/2)^4 * abs(diff(f(x, 4)).*(b-a)).*(b-a);

end
430 chars
17 lines

main.m
% Example usage
f = @(x) exp(-x.^2); % function to integrate
a = 0;
b = 1;
N = 4; % Number of steps should be multiple of 2 for simpsons rule

[I, error] = repeat_simpsons_rule(f, a, b, N);
disp(I) % 0.7468
disp(error) % 1.2661e-06
232 chars
10 lines