Simpsons 1/3 rule function not working

I'm trying to do a simpsons 1/3 rule that you can break up an integration and run the rule on smaller segments to get a more accurate answer. It's fairly basic however when inputed as a function it only gives me the first x interation. When I run the code seperately it works perfect.

function [duck] = simp13(func, a, b, n)

%(b-a)/6 * f(x0) + 4f(x1) + f(x2)
if b <= a; error("b cannot be greater than a"); end
% n equals number of times rule is implimented
l = b-a;
c = l/(2*n);

xvals = a:c:b;
duck = 0;

count = 1;

for i = 1:n
    duck = duck + (xvals(:,count + 2) - xvals(:,count))/6*(func(xvals(:,count)) + 4*func(xvals(:,count + 1)) + func(xvals(:,count + 2)));        
    count = count + 2;    
end
end

Answer 1

[Alan Stevens]

Are you looking for something like this?

a = 1; b = 3;
f = @(x) exp(x);
n = 6;
Ysimp = simp13(f,a,b,n)

    3.0000    1.0754

    5.0000    2.5762

    7.0000    4.6708

    9.0000    7.5940

   11.0000   11.6737

   13.0000   17.3673

Ysimp = 17.3673

Yexact = exp(b)-exp(a)

Yexact = 17.3673

function [duck] = simp13(func, a, b, n)

%(b-a)/6 * f(x0) + 4f(x1) + f(x2)
if b <= a; error("b cannot be greater than a"); end
% n equals number of times rule is implimented
l = b-a;
c = l/(2*n);

xvals = a:c:b;
duck = 0;

count = 1;

for i = 1:n
    duck = duck + (xvals(:,count + 2) - xvals(:,count))/6*(func(xvals(:,count)) + 4*func(xvals(:,count + 1)) + func(xvals(:,count + 2)));        
    count = count + 2;
    disp([count, duck])
end
end