Code covered by the BSD License  

Highlights from
Simpsons Rule Demonstration

image thumbnail

Simpsons Rule Demonstration

by

 

Shows parabolas used in Simpson's rule

simpsons_rule(f,I,m)
function Plot = simpsons_rule(f,I,m)

% f is a function in one variable
% I is a 1x2 vector
% m is the number of points used to create m parabolas.
% For example
%
% simpsons_rule('exp(-x^2)',[-2,3],6)
%


close all
format long

h = (I(2)-I(1))/m;

for j = 1:1:m
    L(j) = [I(1)+h*(j-1)];
end

for j = 1:1:m
    R(j) = L(j)+h;
end

for j = 1:1:m
    c(j)=(L(j)+R(j))/2;
end

for j = 1:1:m
    YL(j) = subs(f,L(j));
    YR(j) = subs(f,R(j));
    Yc(j) = subs(f,c(j));
end

for j = 1:1:m
    polyfit([L(j) R(j) c(j)], [YL(j) YR(j) Yc(j)], 2);
end

for j = 1:1:m
    plot(L(j):.0001:R(j), polyval(polyfit([L(j) R(j) c(j)], [YL(j) YR(j) Yc(j)], 2),[L(j):.0001:R(j)]),'r');
    hold on
end

for j = 1:1:m
    plot(L(j),YL(j),'--ro','LineWidth',2,...
                'MarkerEdgeColor','b',...
                'MarkerFaceColor','g',...
                'MarkerSize',5);
    plot(R(j),YR(j),'--ro','LineWidth',2,...
                'MarkerEdgeColor','b',...
                'MarkerFaceColor','g',...
                'MarkerSize',5);
    plot(c(j),Yc(j),'or','LineWidth', 2);
end

Plot = ezplot(f, [I(1),I(2)]);
for j = 1:1:m
    jbfill(L(j):.0001:R(j),polyval(polyfit([L(j) R(j) c(j)],...
        [YL(j) YR(j) Yc(j)], 2),[L(j):.0001:R(j)]),...
        0*(polyval(polyfit([L(j) R(j) c(j)],[YL(j) YR(j) Yc(j)], 2),...
        [L(j):.0001:R(j)])),rand(1,3),rand(1,3),0,rand(1,1));
    grid on
end

Contact us