clear all clc
prompt = 'Enter n: '; n = input(prompt);
simprl(@(x)(exp(-x^2)),0,1,n);
function s=simprl(f,a,b,n)
array1 = ones(1,n); array2 = ones(1, 2*n); array3 = ones(1,n); array4 = ones(1,n);
fun = @(x)(exp(-x.^2)); true = integral(fun,a,b);
n = n * 2; even =0; odd = 0; y0 = f (a); yn = f (b); x = 0; h=(b-a)/( n / 2);
for i = 0: (( n) + 1)
x=0; even = 0; odd = 0; h=(b-a)/(i + 1); for j = 1:i
x = x + a + h;
r = mod(j,2);
if (r == 0)
even = even + f(x);
else
odd = odd + f(x); end
end
end
I = (h / 3) * ((y0 + yn) + 2 * even + 4 * odd);
fprintf('Value after %i', n/2 ) fprintf(' iterations: %.16f \n', I);
fprintf('True value for the integral: %.16f' , true);
even =0; odd = 0; y0 = f (a); yn = f (b); x = 0; h=(b-a)/( n );
for i = 0: ((2 * n) + 1)
x=0; even = 0; odd = 0; h=(b-a)/(i + 1); for j = 1:i
x = x + a + h;
r = mod(j,2);
if (r == 0)
even = even + f(x);
else
odd = odd + f(x); end
end
y = (h / 3) * ((y0 + yn) + 2 * even + 4 * odd); array2(1, i + 1) = y;
end
even =0; odd = 0; y0 = f (a); yn = f (b); x = 0; h=(b-a)/(n / 2);
for i = 0: n - 1
x=0; even = 0; odd = 0; h=(b-a)/(i + 1);
for j = 1:i
x = x + a + h;
r = mod(j,2);
if (r == 0)
even = even + f(x);
else
odd = odd + f(x); end
end
y = (h / 3) * ((y0 + yn) + 2 * even + 4 * odd); array1(1, i + 1) = y; array4(1, i + 1) = i; l = array1(i + 1); m = array2(2 * (i + 1)); rich = ( ((2 ^ 8) * m) - l ) / ( (2 ^ 8) -1); error = (rich - true ) / true;
error = abs(error);
array3(1, i + 1) = error;
end
a = array3(2:2:end); b = array4(2:1:(n/2) + 1);
loglog(b,a)
title('Composite Simpsons 1/3 rule error %'); xlabel('iterations') ; ylabel('error %');
end