MATLAB: Can someone tell me why the error is so large for the composite simpsons rule

composite simpson's approximationcomposite simpson's ruleerror bounds

My assignment was to code a composite simpsons rule where the exact value of the integral of xsin(x) from 0 to 1 is 0.301168678 but I keep getting 0.862149054988026 as the approximation and that is not within my error bound. Am I entering the composite simpsons rule wrong???
if true
% EV = sin(1)-cos(1); %exact value
a=0; %starting point
b=1; %endpoint
f=@(x) x.*sin(x);
%d=f(a);
for j=0:6
n=10^(j);
h=(b-a)/n;
Err_CS= 1/(36*(n)^4); %error bound
CS=(1/3)*(f(a) + 2.*sum(f((a + 2*h):2*h:(b - 2*h)))+4.*sum(f((a + h):2*h:(b -h)) +f(b)))*h; %composite simpsons rule
fprintf('%8.2e %14.8e %14.8e %14.8e\n', n, CS, abs(CS-EV), Err_CS)
end
end

Best Answer

  • Hi Briyahna,
    Yes the Simpson's rule expression is wrong, but only in the typo sense of having a misplaced parenthesis. Instead of
    CS=(1/3)*(f(a) + 2.*sum(f((a + 2*h):2*h:(b - 2*h)))+4.*sum(f((a + h):2*h:(b -h)) +f(b)))*h;
    it should be
    CS=(1/3)*(f(a) + 2.*sum(f((a + 2*h):2*h:(b - 2*h)))+4.*sum(f((a + h):2*h:(b -h))) +f(b ))*h;
    Even better would be to have (h/3) in front as per the usual Simpson expression rather than the h factor hiding all the way at the end.