Simpson's rule implementation

function I = simpsons()

f=@(x) (3+x.^1/2);
a = 1;
b = 6;
n = 20;
disp('n               I           h')
disp('________________________________')

while (n <= 140)

if numel(f)>1 % If the input provided is a vector
    n=numel(f)-1; 
    h=(b-a)/n;
    I = h/3*(f(1)+2*sum(f(3:2:end-2))+4*sum(f(2:2:end))+f(end));
else % If the input provided is an anonymous function
    h=(b-a)/n; 
    xi=a:h:b;
    I = h/3*(f(xi(1))+2*sum(f(xi(3:2:end-2)))+4*sum(f(xi(2:2:end)))+f(xi(end)));
end

fprintf('%f \t %f \t %f \n', n, I, h);
n=n+20;
end
end

Hey, I'm trying to implement Simpson's rule in matlab. I want to make n every loop n=n+20 and get different I values, but I do not get it. Where is the problem? I is numerical integration formula.

Answer 1

[Jan]

The number of steps does not matter, if the function is linear:

f=@(x) (3+x.^1/2)

x.^1/2 is (x .^ 1) / 2, which is the same as x/2. Operator precedence: power is higher than division. I guess you mean:

f = @(x) (3 + x.^(1/2))
% or: f=@(x) (3 + x.^0.5)
% or: f=@(x) (3 + sqrt(x))  % Fastest